diff --git a/.gitignore b/.gitignore index 1a136df..1f25ad4 100644 --- a/.gitignore +++ b/.gitignore @@ -12,3 +12,5 @@ docs/_build/ docs/_static/ docs/_templates/ docs/xml +.cache +.eggs \ No newline at end of file diff --git a/.travis.yml b/.travis.yml index 7de155a..85d72b8 100644 --- a/.travis.yml +++ b/.travis.yml @@ -1,6 +1,6 @@ sudo: false dist: trusty -language: cpp +language: generic env: global: @@ -35,28 +35,35 @@ before_install: jobs: include: - stage: build + env: + - CC=gcc-6 + - CXX=g++-6 script: - mkdir build && cd build - - export CC=gcc-6 && export CXX=g++-6 - cmake -DBUILD_TESTS=ON -DBUILD_PYTHON=ON .. - make VERBOSE=1 - ./test/TESTS - # same stage + env: + - CC=gcc-5 + - CXX=g++-5 script: - mkdir build && cd build - - export CC=gcc-5 && export CXX=g++-5 - cmake -DBUILD_TESTS=ON -DBUILD_PYTHON=ON .. - make VERBOSE=1 - ./test/TESTS - - stage: docs + + - stage: python_n_docs env: - global: - - secure: "q2oqY0BnXInnLb9Vuq4pI7rJkXzn9HnCZygqpng/9oTdeFEGA5V5j0eQRca5ui7/oL2grhr/bhujTUNb5t3XEWc3dgPx9Z2HDgiR+oyepUEZn3fZA8asErMlgTAYSdycwgAGzUvTbET2ojy0B8ebxEAR3VGWilSkOcbG5EOeIWeHNxdKBrPCmZZztmlVP8z4fM8C9KYEmkG0u86NdtoKFRGU8Lu0W4qhQ3axX/UDFiWxU7tQ2Y8YKSemBxuOY4POIdBhjcv37RFs/zBSs2p0ZVX06fuCqR/v0Z1Y0N/5W2dQNBzhL/BfWJMuf15/HrP3/CFAz6QS0sM+Oy3xwXbVD7V8Hoc4pagq/GNcQxcd07NIsR7L5tWD9GYKYT17wdgN21W5A8IuQEmfyKIyr9A0DSX/Fo43dUFrohuieJr/S0/FrPLL132iCJdAFeZRqJ+jjYY0xnHG1uSNgdnLfAUBvamnM/FwV+AsbrZ2U4HhvqwS+4dsXiSRpHAKG2+vraeUGm5I5zwxOCGrDbPt1aFYce7e3DVoX36ZX58tI/SxLm4FusWNGunEftaEsMqsUxqzF3q0lhnH1XQiaipw4XUo+6bJfP6aZirIQHjbLyd5tbobaUQZQZ7mFO3F0nph+PExmlcVmWnQbv65ytia0Dm21n/zPbHnykQlHkPioKu+1sk=" + - CC=gcc-6 + - CXX=g++-6 + - secure: "q2oqY0BnXInnLb9Vuq4pI7rJkXzn9HnCZygqpng/9oTdeFEGA5V5j0eQRca5ui7/oL2grhr/bhujTUNb5t3XEWc3dgPx9Z2HDgiR+oyepUEZn3fZA8asErMlgTAYSdycwgAGzUvTbET2ojy0B8ebxEAR3VGWilSkOcbG5EOeIWeHNxdKBrPCmZZztmlVP8z4fM8C9KYEmkG0u86NdtoKFRGU8Lu0W4qhQ3axX/UDFiWxU7tQ2Y8YKSemBxuOY4POIdBhjcv37RFs/zBSs2p0ZVX06fuCqR/v0Z1Y0N/5W2dQNBzhL/BfWJMuf15/HrP3/CFAz6QS0sM+Oy3xwXbVD7V8Hoc4pagq/GNcQxcd07NIsR7L5tWD9GYKYT17wdgN21W5A8IuQEmfyKIyr9A0DSX/Fo43dUFrohuieJr/S0/FrPLL132iCJdAFeZRqJ+jjYY0xnHG1uSNgdnLfAUBvamnM/FwV+AsbrZ2U4HhvqwS+4dsXiSRpHAKG2+vraeUGm5I5zwxOCGrDbPt1aFYce7e3DVoX36ZX58tI/SxLm4FusWNGunEftaEsMqsUxqzF3q0lhnH1XQiaipw4XUo+6bJfP6aZirIQHjbLyd5tbobaUQZQZ7mFO3F0nph+PExmlcVmWnQbv65ytia0Dm21n/zPbHnykQlHkPioKu+1sk=" script: - set -e - - export CC=gcc-6 && export CXX=g++-6 - - pip install doctr sphinx sphinx_rtd_theme breathe - - pip install -vvv . + - pip install numpy doctr sphinx sphinx_rtd_theme breathe + - python setup.py build + - python setup.py test + - python setup.py install - cd docs - doxygen Doxyfile - make html @@ -64,21 +71,16 @@ jobs: - doctr deploy . --built-docs docs/_build/html --no-require-master # TODO: remove no-require-master - - stage: deploy - script: skip - before_deploy: - - deactivate - - pip install --upgrade setuptools wheel twine + - stage: pypi_deploy + env: + - CC=gcc-6 + - CXX=g++-6 + - secure: "cEGT3jj7NFUKhDmTCEqLtJvR5jBXD9E/MfjN+0Qqn/uuaib2iD9IUKfr2iuimG8AwtN8Ihf/NwOuICdKb0e+kUnpBOq710S+UeyXLjzp6ezsGIZTtYOXO8EpDcdD2aULHA4laGJ9EsqCY3jFQI60jRuCEJnh00UX9JcaGymmGcwuY7yyGY121igFDYPbseApUIlTA9ISRdg2fHFercWmGHnHa6oxi9eAbLNISvaxIbbdVWiOPBp1PX6p90VI3Okdk3+MLdblFoWuZD5vF7jYxa72+/RoOrnEJ6oAnIVlEFs/XlSqjNmomINVEkEo8QQbshbxOMVYHW1otFbbJiwaRqnkvbtKC00uAn0P1zAjeekxLw+PxMfTPn64HyxQoUN8sjXX04w3FhwFZbsZPoVc0/RZXwcOto9mTTPygCL7pvfkvkGU+XDwkhbtWBrqgUuIsPDuvGtAUoeHWTZD5grP9QFaxZuG8J5if3p2CpSKg9BW6itA1qu5pdW6xDaUX0stS4pga2s4yb9kt8RybvbZx7aeMzbj4L21YHMbcJZiojvaG9eR3l2f8ORx23fu4b1uhtY8nZGQ89jmzCxpgqE5pb1jiCqI3sfTJYOCj+NYI0CHWlTdX8Xw8eyD7cXjfPaqvEwArmu//NukHvqta9ic3oNWnOokBggtiJgLjUkrbf4=" + - TWINE_USERNAME=mirthlessmarmot + script: + - pip install --upgrade numpy setuptools wheel twine - pwd - ls - - export CC=gcc-6 && export CXX=g++-6 - deploy: - - provider: pypi - user: mirthlessmarmot - password: - secure: "q/JwK15wRzViATJ4RtUmvV6YlVfv6yqPyKzH1W8lrzBvtpygi8zeSsI94dhxrR0vla0rEsjTUafN10UiflGLAwgTBOv80PLfch6miJTpFAazFVQb4bv3Qu59dQm3EYoHtgiMGgtO/gJJE1/jDGp9SiXS9EJMerL4l7kL7NERPFDGV3vwIYtmTsbzw8YKdwR53n6oUus5JY49a2AJH181yMy8iAhHmYpgqoyCnLlJqnvpKJw7svFlAZiaxv+9NoO9jkP/icdKCxnLJCLXhgwD6okXdFi/GrNj3pLG5FYL3Lz3dMv2aqJMkhkUnSrtWnSC1vvc5ZRaF5KcFSQ4QYFx7aRyxIYs7Df0ZsmeUOIGA++nk3YmkEF1UXM/VOa5SkmNzo8OKK88jRj8jxjE8kBiZw13TgjKr0WN4IjRoo8wmJfB3E5o51IUDR3L7YNdgGXglUOIJvkmSCWWw16wAazUX2ndAxMfyBoxP2jtHJHmvuyLHy1WcibXmhWoEbP+uwH0309ajH18LEOxmoqDBUYol4BWKczZQoF85WySLAlkVbpTmCWuIoMqhfISaatnBuVR7fM84q+aHCDUMytWR+55WDcqrM4s/v1QYc2SdXqBqYbu4lhlgIL9nfOEYn4oyQksboP1LBVULpGFuPAlYX3YrLyVZw9nzGubRXXnKdRpBDc=" - on: - all_branches: true - tags: true - server: "https://test.pypi.org/legacy/" - distributions: sdist \ No newline at end of file + - python setup.py sdist + - if [ -n "$TRAVIS_TAG" ]; then twine upload -r "https://test.pypi.org/legacy/" dist/*; else echo "package _not_ uploaded to testpypi"; fi + - if [ "$TRAVIS_BRANCH" = "master" -a -n "$TRAVIS_TAG" ]; then twine upload dist/*; else echo "package _not_ uploaded to pypi"; fi \ No newline at end of file diff --git a/catana/include/catana/points/Point.hpp b/catana/include/catana/points/Point.hpp index d625c64..a15d1b4 100644 --- a/catana/include/catana/points/Point.hpp +++ b/catana/include/catana/points/Point.hpp @@ -9,7 +9,7 @@ namespace catana { //! Representation of a particle in 3-dimensional space in spherical coordinates struct Point { //! Default constructor - Point() = default; + Point(); //! Construct from cartesian coordinates Point(const double& x, const double& y, const double& z); diff --git a/catana/src/points/Point.cpp b/catana/src/points/Point.cpp index a0b61ff..e873e8a 100644 --- a/catana/src/points/Point.cpp +++ b/catana/src/points/Point.cpp @@ -9,6 +9,9 @@ namespace catana { p.normalize(); } + Point::Point() + : r(0), p(0, 0) {} + Point point_from_box_position(const double& pos_x, const double& pos_y, const double& pos_z, const double& shift, const double& hubble_param) { Point obj( diff --git a/docs/python/catana.rst b/docs/python/catana.rst index bf22d55..9467c04 100644 --- a/docs/python/catana.rst +++ b/docs/python/catana.rst @@ -20,6 +20,7 @@ PointContainer .. autoclass:: PointContainer :members: + :special-members: __getitem__ PixelizedPointContainer diff --git a/python/py_basictypes.cpp b/python/py_basictypes.cpp index f82981e..bc00578 100644 --- a/python/py_basictypes.cpp +++ b/python/py_basictypes.cpp @@ -22,89 +22,89 @@ PYBIND11_PLUGIN(basictypes) // Binding for the "Point" class, with cartesian and spherical output methods py::class_(m, "Point", py::buffer_protocol(), R"pbdoc( - Representation of a particle in 3-dimensional space in spherical coordinates. - - **Convention**: - - longitude phi in (0, 2pi) - - lattitude theta in (0, theta) measured from z-axis. - - (can be converted to numpy array (r, theta, phi) by calling ``np.array(point)``) - )pbdoc") - .def(py::init<>(), R"pbdoc( - Default constructor (0,0,0) - )pbdoc") - .def(py::init(), R"pbdoc( - Constructor from spherical coordinates - )pbdoc", py::arg("r"), py::arg("theta"), py::arg("phi")) - .def("spherical", [](const Point& obj) { - return std::make_tuple(obj.r, obj.p.theta, obj.p.phi); - }, R"pbdoc( - get spherical point coordinates - - Returns: - tuple[float]: spherical coordinates (r, theta, phi) - )pbdoc") - .def("cartesian", [](const Point& obj) { - vec3 v = obj.p.to_vec3(); - return std::make_tuple(obj.r*v.x, obj.r*v.y, obj.r*v.z); - }, R"pbdoc( - get cartesian point coordinates - - Returns: - tuple[float]: cartesian coordinates (x, y, z) - )pbdoc") - .def_buffer([](Point& point) -> py::buffer_info { - return py::buffer_info( - &point.r, - sizeof(double), - py::format_descriptor::format(), - 1, - {3}, - {sizeof(double)} - ); - }); + Representation of a particle in 3-dimensional space in spherical coordinates. + + **Convention**: + - longitude phi in (0, 2pi) + - lattitude theta in (0, theta) measured from z-axis. + + (can be converted to numpy array (r, theta, phi) by calling ``np.array(point)``) + )pbdoc") + .def(py::init<>(), R"pbdoc( + Default constructor (0,0,0) + )pbdoc") + .def(py::init(), R"pbdoc( + Constructor from cartesian coordinates + )pbdoc", py::arg("x"), py::arg("y"), py::arg("z")) + .def("spherical", [](const Point& obj) { + return std::make_tuple(obj.r, obj.p.theta, obj.p.phi); + }, R"pbdoc( + get spherical point coordinates + + Returns: + tuple[float]: spherical coordinates (r, theta, phi) + )pbdoc") + .def("cartesian", [](const Point& obj) { + vec3 v = obj.p.to_vec3(); + return std::make_tuple(obj.r*v.x, obj.r*v.y, obj.r*v.z); + }, R"pbdoc( + get cartesian point coordinates + + Returns: + tuple[float]: cartesian coordinates (x, y, z) + )pbdoc") + .def_buffer([](Point& point) -> py::buffer_info { + return py::buffer_info( + &point.r, + sizeof(double), + py::format_descriptor::format(), + 1, + {3}, + {sizeof(double)} + ); + }); // Binding of the PointContainer class py::class_(m, "PointContainer", py::buffer_protocol(), R"pbdoc( - Collection of :class:`Point` objects + Collection of :class:`Point` objects - (can be converted to numpy array with shape (n_points, 3), where each point is of the form - (r, theta, phi) by calling ``np.array(point_container)``) - )pbdoc") - .def(py::init<>(), R"pbdoc( + (can be converted to numpy array with shape (n_points, 3), where each point is of the form + (r, theta, phi) by calling ``np.array(point_container)``) + )pbdoc") + .def(py::init<>(), R"pbdoc( Default constructor (empty container) - )pbdoc") - .def("__init__", [](PointContainer& oc, py::array_t array, std::string coord_type) - { - bool spherical; - if(coord_type==std::string("cartesian")) - spherical=false; - else if(coord_type==std::string("spherical")) - spherical= true; - else - throw std::runtime_error("coord_type either 'cartesian' or 'spherical'"); - - py::buffer_info info = array.request(); - if(info.ndim != 2 ) - throw std::runtime_error("Incompatible number of dimensions of array (needs 2)"); - if(info.shape[1] != 3) - throw std::runtime_error("Incompatible number of columns of array (needs 3)"); - - new (&oc) PointContainer(); // - double* data_ptr = (double*) info.ptr; // A pointer to the storage - double* data_ptr_end = data_ptr + info.shape[0]*3; // A pointer behind the storage - - if (info.strides[1] == sizeof(double)){ // Correct data ordering - for(double* current_ptr = data_ptr; current_ptr != data_ptr_end; current_ptr+=3){ - if(spherical) - oc.add_point(point_from_spherical_position(current_ptr, 1)); - else - oc.add_point(point_from_box_position(current_ptr, 0, 1)); - } - } else { // Wrong data ordering - throw std::runtime_error("Data ordering of array is ColumnMajor, but RowMajor needed."); + )pbdoc") + .def("__init__", [](PointContainer& oc, py::array_t array, std::string coord_type) + { + bool spherical; + if(coord_type==std::string("cartesian")) + spherical=false; + else if(coord_type==std::string("spherical")) + spherical= true; + else + throw std::runtime_error("coord_type either 'cartesian' or 'spherical'"); + + py::buffer_info info = array.request(); + if(info.ndim != 2 ) + throw std::runtime_error("Incompatible number of dimensions of array (needs 2)"); + if(info.shape[1] != 3) + throw std::runtime_error("Incompatible number of columns of array (needs 3)"); + + new (&oc) PointContainer(); // + double* data_ptr = (double*) info.ptr; // A pointer to the storage + double* data_ptr_end = data_ptr + info.shape[0]*3; // A pointer behind the storage + + if (info.strides[1] == sizeof(double)){ // Correct data ordering + for(double* current_ptr = data_ptr; current_ptr != data_ptr_end; current_ptr+=3){ + if(spherical) + oc.add_point(point_from_spherical_position(current_ptr, 1)); + else + oc.add_point(point_from_box_position(current_ptr, 0, 1)); } - }, R"pbdoc( + } else { // Wrong data ordering + throw std::runtime_error("Data ordering of array is ColumnMajor, but RowMajor needed."); + } + }, R"pbdoc( Create a PointContainer from a (numpy) array containing the coordinates Parameters: @@ -113,104 +113,108 @@ Create a PointContainer from a (numpy) array containing the coordinates Note: coordinate array has to be in row-major ordering. - )pbdoc", py::arg("coordinates"), py::arg("coordinate_type")="spherical") - - .def("__getitem__", [](const PointContainer& oc, size_t i) - { - return oc[i]; - }) - - // Buffer definition, so that we can call numpy.array(point_container) in python - .def_buffer([](PointContainer& oc) -> py::buffer_info { - return py::buffer_info( - (double*) oc.data(), - sizeof(double), - py::format_descriptor::format(), - 2, - { oc.size(), 3}, - { sizeof(double) * 3, sizeof(double)} - ); - }) - .def("add_point", &PointContainer::add_point, R"pbdoc( - Add a point to the collection - - Parameters: - point (Point): the point to add - )pbdoc", py::arg("point")); + )pbdoc", py::arg("coordinates"), py::arg("coordinate_type")="spherical") + + .def("__getitem__", [](const PointContainer& oc, size_t i) + { + if(i >= oc.size()) + throw pybind11::index_error("index out of bounds!"); + return oc[i]; + }) + + // Buffer definition, so that we can call numpy.array(point_container) in python + .def_buffer([](PointContainer& oc) -> py::buffer_info { + return py::buffer_info( + (double*) oc.data(), + sizeof(double), + py::format_descriptor::format(), + 2, + { oc.size(), 3}, + { sizeof(double) * 3, sizeof(double)} + ); + }) + .def("add_point", &PointContainer::add_point, R"pbdoc( + Add a point to the collection + + Parameters: + point (Point): the point to add + )pbdoc", py::arg("point")) + .def("__len__", [](const PointContainer& oc){ return oc.size(); }); py::class_(m, "PixelizedPointContainer", R"pbdoc( - Pixelized collection of :class:`Point` objects + Pixelized collection of :class:`Point` objects - Pixelization occurs along the angular axis according to the HEALPix scheme. - Radial coordinates remain the same. - )pbdoc") - .def(py::init(), R"pbdoc( + Pixelization occurs along the angular axis according to the HEALPix scheme. + Radial coordinates remain the same. + )pbdoc") + .def(py::init(), R"pbdoc( Default constructor (empty container) Parameters: nside (int): NSide resolution parameter of the HEALPix scheme. Must be a power of 2. - )pbdoc", py::arg("nside")) - .def(py::init(), R"pbdoc( + )pbdoc", py::arg("nside")) + .def(py::init(), R"pbdoc( Constructor from a :class:`PointContainer` Parameters: nside (int): NSide resolution parameter of the HEALPix scheme. Must be a power of 2. point_container (PointContainer): A container of :class:`Point` objects - )pbdoc", py::arg("nside"), py::arg("point_container")) - .def("get_countmap", &PixelizedPointContainer::get_countmap, R"pbdoc( - Obtain a HEALPix count map which contains the number of points which are stored in each pixel. + )pbdoc", py::arg("nside"), py::arg("point_container")) + .def("get_countmap", &PixelizedPointContainer::get_countmap, R"pbdoc( + Obtain a HEALPix count map which contains the number of points which are stored in each pixel. - Returns: - numpy.ndarray[int]: HEALPix map with point number counts + Returns: + numpy.ndarray[int]: HEALPix map with point number counts - Note: - The resulting array can be viewed with healpy. - )pbdoc"); + Note: + The resulting array can be viewed with healpy. + )pbdoc"); py::class_(m, "FunctionInterpolator", R"pbdoc( - A function sampler and interpolator - - Samples the given function with ``interpolation_points`` number of equidistantly distributed points between - ``x_min`` and ``x_max``. When called at a given position, interpolates linearly from the two closest sample - points. - - The FunctionInterpolator mainly serves to speed up computations since we do not need to do expensive - Python object calls once it is initialized. - )pbdoc") - .def("__init__", [](FunctionInterpolator& fi, std::function fct, - unsigned int interpolation_points, double x_min, double x_max) { - new (&fi) FunctionInterpolator(fct, interpolation_points, x_min, x_max, false); - }, R"pbdoc( - Construct a FunctionInterpolator from a function - - Parameters: - function (callable object): the function to be sampled and interpolated - interpolation points (int): number of points at which the function needs to be sampled - x_min (float): lower boundary of sample points - x_max (float): upper boundary of sample points - )pbdoc", py::arg("function"), py::arg("interpolation_points"), py::arg("x_min"), py::arg("x_max")) - .def("__call__", &FunctionInterpolator::operator(), R"pbdoc( + A function sampler and interpolator + + Samples the given function with ``interpolation_points`` number of equidistantly distributed points between + ``x_min`` and ``x_max``. When called at a given position, interpolates linearly from the two closest sample + points. Note that evaluating the function for values ``x >= x_max`` and ``x < x_min`` raises a + ValueError. + + The FunctionInterpolator mainly serves to speed up computations since we do not need to do expensive + Python object calls once it is initialized. + )pbdoc") + .def("__init__", [](FunctionInterpolator& fi, std::function fct, + unsigned int interpolation_points, double x_min, double x_max) { + new (&fi) FunctionInterpolator(fct, interpolation_points, x_min, x_max, false); + }, R"pbdoc( + Construct a FunctionInterpolator from a function + + Parameters: + function (callable object): the function to be sampled and interpolated + interpolation points (int): number of points at which the function needs to be sampled + x_min (float): lower boundary of sample points + x_max (float): upper boundary of sample points + )pbdoc", py::arg("function"), py::arg("interpolation_points"), py::arg("x_min"), py::arg("x_max")) + .def("__call__", &FunctionInterpolator::operator(), R"pbdoc( Evaluate the FunctionInterpolator at the given value x. Warning: - x must be within the boundaries [x_min, x_max] for which the FunctionInterpolator was set up. + x must be within the boundaries [x_min, x_max) for which the FunctionInterpolator was set up. Parameters: x (float): evaluation point - )pbdoc", py::arg("x")) - .def("__call__", [](FunctionInterpolator& fi, py::array_t x) { - auto stateful_closure = [&fi](double x){return fi(x);}; - return py::vectorize(stateful_closure)(x); - }, R"pbdoc( + )pbdoc", py::arg("x")) + .def("__call__", [](FunctionInterpolator& fi, py::array_t x) { + auto stateful_closure = [&fi](double x){return fi(x);}; + return py::vectorize(stateful_closure)(x); + }, R"pbdoc( Evaluate the FunctionInterpolator at the given values x. Warning: - each element in x must be within the boundaries [x_min, x_max] for which the FunctionInterpolator + each element in x must be within the boundaries [x_min, x_max) for which the FunctionInterpolator was set up. Parameters: x (numpy.ndarray[float]): evaluation points - )pbdoc", py::arg("x")); + )pbdoc", py::arg("x")); return m.ptr(); } \ No newline at end of file diff --git a/python/py_io.cpp b/python/py_io.cpp index cf2d7b6..f223f86 100644 --- a/python/py_io.cpp +++ b/python/py_io.cpp @@ -343,7 +343,8 @@ Constructor from a HEALPix map array sink (Sink): to where the filtered points are stored buffer_size (int): size of the buffer (in number of points). Faster if larger, but consumes more memroy. verbose (bool): if ``True``, print additional logging information to stdout - )pbdoc", py::arg("source"), py::arg("sink"), py::arg("buffer_size")=1000000, py::arg("verbose")=true) + )pbdoc", py::arg("source"), py::arg("sink"), py::arg("buffer_size")=1000000, py::arg("verbose")=true, + py::keep_alive<1, 2>(), py::keep_alive<1, 3>()) .def("add_filter", &io::FilterStream::add_filter, R"pbdoc( Add a Filter to the pipeline. @@ -351,7 +352,7 @@ Constructor from a HEALPix map array Parameters: filter (Filter): A radial or angular Filter instance - )pbdoc", py::arg("filter")) + )pbdoc", py::arg("filter"), py::keep_alive<1, 2>()) .def("set_source", &io::FilterStream::set_source, R"pbdoc( Reset the source to a new source @@ -360,7 +361,7 @@ Constructor from a HEALPix map array Parameters: source (Source): a source instance - )pbdoc", py::arg("source")) + )pbdoc", py::arg("source"), py::keep_alive<1, 2>()) .def("run", &io::FilterStream::run, R"pbdoc( Run the pipeline diff --git a/setup.cfg b/setup.cfg index e70847f..7d1ad97 100644 --- a/setup.cfg +++ b/setup.cfg @@ -7,4 +7,10 @@ tag_prefix = parentdir_prefix = catana- [metadata] -description-file = README.rst \ No newline at end of file +description-file = README.rst + +[aliases] +test = pytest + +[tool:pytest] +testpaths = test/python \ No newline at end of file diff --git a/setup.py b/setup.py index b1e391c..39c598a 100644 --- a/setup.py +++ b/setup.py @@ -76,5 +76,7 @@ def get_cmdclass(): ext_modules=[CMakeExtension('basictypes'), CMakeExtension('besseltools'), CMakeExtension('decomposition'), CMakeExtension('io')], cmdclass=get_cmdclass(), zip_safe=False, - install_requires=['numpy'] + install_requires=['numpy'], + setup_requires=['pytest-runner'], + tests_require=['pytest'] ) \ No newline at end of file diff --git a/test/python/besseltools_test.py b/test/python/besseltools_test.py new file mode 100644 index 0000000..d04abbf --- /dev/null +++ b/test/python/besseltools_test.py @@ -0,0 +1,58 @@ +import numpy as np +import catana +import pytest + + +class TestSphericalBesselZeros(object): + def test_default_constructor(self): + with pytest.raises(TypeError): + sfb = catana.besseltools.SphericalBesselZeros() + + def test_zeros(self): + sfb = catana.besseltools.SphericalBesselZeros(0) + + assert sfb[0] == pytest.approx(np.pi) + assert sfb[1] == pytest.approx(6.283185307179586) + assert sfb[2] == pytest.approx(9.42477796076938) + + sfb = catana.besseltools.SphericalBesselZeros(5) + assert sfb[9] == pytest.approx(38.8836) + + def test_computeupto(self): + sfb = catana.besseltools.SphericalBesselZeros(2) + sfb.compute_up_to(100) + + +class TestBesselIntegrator(object): + def test_known_values(self): + w = lambda r: 1 + assert catana.besseltools.double_sbessel_integrator(w, 0, 100, 1, 1) == pytest.approx(1.56582, abs=1e-5) + assert catana.besseltools.double_sbessel_integrator(w, 0, 100, 1, 2) == pytest.approx(0.785393, abs=1e-5) + assert catana.besseltools.double_sbessel_integrator(w, 4, 100, 1, 4.2) == pytest.approx(0.000133933, abs=1e-9) + + def test_k_symmetric(self): + w = lambda r: np.exp(-(r/40)**2) * r**2 + res1 = catana.besseltools.double_sbessel_integrator(w, 6, 100, 2.3, 4.5) + res2 = catana.besseltools.double_sbessel_integrator(w, 6, 100, 4.5, 2.3) + assert res1 == pytest.approx(res2) + + def test_function_interpolator(self): + w = lambda r: np.exp(-(r/40)**2) * r**2 + fi = catana.FunctionInterpolator(w, 100000, 0, 100) + + assert fi(0) == w(0) + assert fi(99.9) == pytest.approx(w(99.9)) + assert fi(20) == pytest.approx(fi(20)) + + with pytest.raises(ValueError): + fi(100) + + with pytest.raises(ValueError): + fi(-1) + + with pytest.raises(ValueError): + fi(120) + + res1 = catana.besseltools.double_sbessel_integrator(w, 6, 100, 2.3, 4.5) + res2 = catana.besseltools.double_sbessel_integrator(fi, 6, 100, 2.3, 4.5) + assert res1 == pytest.approx(res2, rel=1e-5) diff --git a/test/python/decomposition_test.py b/test/python/decomposition_test.py new file mode 100644 index 0000000..5c1ad49 --- /dev/null +++ b/test/python/decomposition_test.py @@ -0,0 +1,73 @@ +import numpy as np +import catana +import pytest +import os + + +@pytest.fixture +def datadir(): + return os.path.join(os.path.dirname(__file__), '..', 'test_data') + + +@pytest.fixture +def gaussian_source(datadir): + return { + 'correction_factor': 5535184878.03/(2500**3 * 4/3 * np.pi), + 'source': catana.io.CartesianTextSource(os.path.join(datadir, "gaussian_catalog.txt")), + 'c_ln': np.loadtxt(os.path.join(datadir, "gaussian_catalog.python.c_ln")), + 'c_ln_rev': np.loadtxt(os.path.join(datadir, "gaussian_catalog.pythonrev.c_ln")), + 'k_ln': np.loadtxt(os.path.join(datadir, "gaussian_catalog.python.k_ln")) + } + + +class TestSFB(object): + def test_raw(self, gaussian_source): + pc = gaussian_source['source'].get_point_container() + assert(np.array(pc).shape == (445, 3)) + + kclkk = catana.sfb_decomposition(pc, 10, 10, 2500, False, False) + + # correct volume (gaussian filter was applied on data) + k_ln = kclkk.k_ln + c_ln = kclkk.c_ln * gaussian_source['correction_factor']**2 + + assert np.allclose(k_ln, gaussian_source['k_ln']) + assert np.allclose(c_ln, gaussian_source['c_ln'], rtol=1e-2) + + def test_rev(self, gaussian_source): + pc = gaussian_source['source'].get_point_container() + assert(np.array(pc).shape == (445, 3)) + ppc = catana.PixelizedPointContainer(64, pc) + + kclkk = catana.sfb_decomposition(ppc, 10, 10, 2500, False, False) + + # correct volume (gaussian filter was applied on data) + k_ln = kclkk.k_ln + c_ln = kclkk.c_ln * gaussian_source['correction_factor']**2 + + assert np.allclose(k_ln, gaussian_source['k_ln']) + assert np.allclose(c_ln, gaussian_source['c_ln_rev'], rtol=1e-2) + + +class TestAnalyzer(object): + def test_raw(self, gaussian_source): + a = catana.Analyzer(gaussian_source['source']) + kclkk = a.compute_sfb(10, 10, 2500, False, False) + + # correct volume (gaussian filter was applied on data) + k_ln = kclkk.k_ln + c_ln = kclkk.c_ln * gaussian_source['correction_factor']**2 + + assert np.allclose(k_ln, gaussian_source['k_ln']) + assert np.allclose(c_ln, gaussian_source['c_ln'], rtol=1e-2) + + def test_rev(self, gaussian_source): + a = catana.Analyzer(gaussian_source['source']) + kclkk = a.compute_sfb_pixelized(10, 10, 2500, 64, False, False) + + # correct volume (gaussian filter was applied on data) + k_ln = kclkk.k_ln + c_ln = kclkk.c_ln * gaussian_source['correction_factor']**2 + + assert np.allclose(k_ln, gaussian_source['k_ln']) + assert np.allclose(c_ln, gaussian_source['c_ln_rev'], rtol=1e-2) diff --git a/test/python/io_test.py b/test/python/io_test.py new file mode 100644 index 0000000..bd68e53 --- /dev/null +++ b/test/python/io_test.py @@ -0,0 +1,124 @@ +import numpy as np +import catana +import pytest +import os + +@pytest.fixture +def datadir(): + return os.path.join(os.path.dirname(__file__), '..', 'test_data') + +class TestText(object): + def test_source(self, datadir): + s = catana.io.CartesianTextSource(os.path.join(datadir, "mock_data.txt")) + pc = s.get_point_container() + pc2 = s.get_point_container() + + assert len(pc) == len(pc2) == 20 + + def test_sinks(self, datadir): + pass + + +class TestGadget(object): + pass + + +class TestFilter(object): + def test_tophat(self): + data = np.zeros((5, 3)) + data[:, 0] = np.arange(5, 0, -1) + pc = catana.PointContainer(data) + tf = catana.io.TophatRadialWindowFunctionFilter(3.1) + tf(pc) + + assert np.array(pc).shape == (3, 3) + + def test_gaussian(self): + N = 1<<24 + data = np.zeros((N, 3)) + data[:, 0] = 1000 * np.random.uniform(0, 1, N)**(1/3) + pc = catana.PointContainer(data) + gf = catana.io.GaussianRadialWindowFunctionFilter(100) + gf(pc) + + assert np.array(pc).shape[0]/N == pytest.approx(0.001329340388, rel=0.1) + + def test_angular(self): + pc = catana.PointContainer(np.array([[1, 0, -1], [1, 0, 1]])) + af = catana.io.AngularMaskFilter(np.array([0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1])) + + af(pc) + assert np.array(pc).shape[0] == 1 + assert pc[0].spherical()[2] == -1 + + def test_angular_file(self, datadir): + pc = catana.PointContainer(np.array([[1, 0, -1], [1, 0, 1]])) + af = catana.io.AngularMaskFilter(os.path.join(datadir, 'testmap.fits')) + + af(pc) + assert np.array(pc).shape[0] == 1 + assert pc[0].spherical()[2] == -1 + + def test_angular_invalid_map(self): + with pytest.raises(RuntimeError): + af = catana.io.AngularMaskFilter(np.zeros(11)) + + +class TestFilterStream(object): + def test_tophat_stream(self): + N = 1000 + data = np.zeros((N, 3)) + data[:, 0] = np.random.permutation(np.arange(N)) + pc = catana.PointContainer(data) + source = catana.io.PointContainerSource(pc) + + pc2 = catana.PixelizedPointContainer(16) + sink = catana.io.PixelizedPointContainerSink(pc2) + + f = catana.io.TophatRadialWindowFunctionFilter(222.4) + + fs = catana.io.FilterStream(source, sink) + fs.add_filter(f) + fs.run() + + assert np.sum(pc2.get_countmap()) == 223 + + def test_multisource_stream(self): + N = 1000 + data = np.zeros((N, 3)) + data[:, 0] = np.random.permutation(np.arange(N)) + pc = catana.PointContainer(data) + source = catana.io.PointContainerSource(pc) + + pc2 = catana.PixelizedPointContainer(16) + sink = catana.io.PixelizedPointContainerSink(pc2) + + f = catana.io.TophatRadialWindowFunctionFilter(222.4) + + fs = catana.io.FilterStream(source, sink) + fs.add_filter(f) + + fs.run_totemp() + source.reset() + fs.run_totemp() + + fs.run_fromtemp() + assert np.sum(pc2.get_countmap()) == 446 + + def test_filterscope(self): + # This test was added because if the filter was not saved in a variable, this caused SIGSEGV + N = 1000 + data = np.zeros((N, 3)) + data[:, 0] = np.random.permutation(np.arange(N)) + pc = catana.PointContainer(data) + source = catana.io.PointContainerSource(pc) + + pc2 = catana.PixelizedPointContainer(16) + sink = catana.io.PixelizedPointContainerSink(pc2) + + fs = catana.io.FilterStream(source, sink) + + fs.add_filter(catana.io.TophatRadialWindowFunctionFilter(222.4)) + fs.run() + + assert np.sum(pc2.get_countmap()) == 223 \ No newline at end of file diff --git a/test/python/point_test.py b/test/python/point_test.py new file mode 100644 index 0000000..c20e7c3 --- /dev/null +++ b/test/python/point_test.py @@ -0,0 +1,62 @@ +import numpy as np +import catana +import pytest + +class TestPoint(object): + def test_default_constructor(self): + p = catana.Point() + assert p.cartesian() == (0.0, 0.0, 0.0) + assert p.spherical() == (0.0, 0.0, 0.0) + assert np.allclose(p, np.zeros(3)) + + def test_value_constructor(self): + p = catana.Point(0, 0, 1) + assert p.spherical() == (1, 0, 0) + assert p.cartesian() == (0, 0, 1) + + +class TestPointContainer(object): + def test_default_constructor(self): + pc = catana.PointContainer() + assert np.array(pc).shape == (0, 3) + + def test_add_point(self): + pc = catana.PointContainer() + pc.add_point(catana.Point()) + assert np.array(pc).shape == (1, 3) + assert np.allclose(pc, np.zeros((1, 3))) + + def test_numpy_constructor(self): + data = np.arange(9).reshape(3, 3) + pc = catana.PointContainer(data) + assert(np.allclose(pc, data)) + + def test_getitem(self): + pc = catana.PointContainer(np.ones((2, 3))) + point0 = pc[0] + point1 = pc[1] + assert np.allclose(point0, point1) + assert np.allclose(point0, [1, 1, 1]) + + with pytest.raises(IndexError): + point2 = pc[2] + + +class TestPixelizedPointContainer(object): + def test_default_constructor(self): + with pytest.raises(TypeError): + ppc = catana.PixelizedPointContainer() + + def test_empty_constructor(self): + ppc = catana.PixelizedPointContainer(4) + cm = ppc.get_countmap() + assert np.all(cm == np.zeros(4**2 * 12)) + + def test_pc_constructor(self): + data = np.array([[1, 0, 0], [1, 0, 0], [1, np.pi, 0]]) + pc = catana.PointContainer(data) + ppc = catana.PixelizedPointContainer(32, pc) + cm = ppc.get_countmap() + assert cm.shape == (32**2 * 12, ) + assert np.sum(cm) == 3 + assert np.count_nonzero(cm) == 2