diff --git a/SConstruct b/SConstruct index d2d5e8a0c6c..894aa166ad0 100644 --- a/SConstruct +++ b/SConstruct @@ -1106,6 +1106,7 @@ PyMODINIT_FUNC initcheck_tmv(void) def CheckNumPy(config): numpy_source_file = """ +#define NPY_NO_DEPRECATED_API NPY_1_7_API_VERSION #include "Python.h" #include "numpy/arrayobject.h" diff --git a/pysrc/CppShear.cpp b/pysrc/CppShear.cpp index a90af9e9a68..d83a5739129 100644 --- a/pysrc/CppShear.cpp +++ b/pysrc/CppShear.cpp @@ -20,10 +20,7 @@ */ #include "boost/python.hpp" #include "CppShear.h" - -#define PY_ARRAY_UNIQUE_SYMBOL GALSIM_ARRAY_API -#define NO_IMPORT_ARRAY -#include "numpy/arrayobject.h" +#include "NumpyHelper.h" namespace bp = boost::python; @@ -32,22 +29,19 @@ namespace { struct PyCppShear { - static bp::handle<> getMatrix(CppShear const & self) { + static bp::handle<> getMatrix(const CppShear& self) { static npy_intp dim[2] = {2, 2}; // Because the C++ version sets references that are passed in, and that's not possible in // Python, we wrap this instead, which returns a numpy array. double a=0., b=0., c=0.; self.getMatrix(a, b, c); - bp::handle<> r(PyArray_SimpleNew(2, dim, NPY_DOUBLE)); - *reinterpret_cast(PyArray_GETPTR2(r.get(), 0, 0)) = a; - *reinterpret_cast(PyArray_GETPTR2(r.get(), 1, 1)) = b; - *reinterpret_cast(PyArray_GETPTR2(r.get(), 0, 1)) = c; - *reinterpret_cast(PyArray_GETPTR2(r.get(), 1, 0)) = c; - return r; + double ar[4] = { a, c, c, b }; + PyObject* r = PyArray_SimpleNewFromData(2, dim, NPY_DOUBLE, ar); + return bp::handle<>(r); } static void wrap() { - static char const * doc = + static const char* doc = "CppShear is represented internally by e1 and e2, which are the second-moment\n" "definitions: ellipse with axes a & b has e=(a^2-b^2)/(a^2+b^2).\n" "But can get/set the ellipticity by other measures:\n" @@ -114,21 +108,17 @@ struct PyCppShear { struct PyCppEllipse { - static bp::handle<> getMatrix(CppEllipse const & self) { + static bp::handle<> getMatrix(const CppEllipse& self) { static npy_intp dim[2] = {2, 2}; // Because the C++ version sets references that are passed in, and that's not possible in // Python, we wrap this instead, which returns a numpy array. tmv::Matrix m = self.getMatrix(); - bp::handle<> r(PyArray_SimpleNew(2, dim, NPY_DOUBLE)); - *reinterpret_cast(PyArray_GETPTR2(r.get(), 0, 0)) = m(0,0); - *reinterpret_cast(PyArray_GETPTR2(r.get(), 1, 1)) = m(1,1); - *reinterpret_cast(PyArray_GETPTR2(r.get(), 0, 1)) = m(0,1); - *reinterpret_cast(PyArray_GETPTR2(r.get(), 1, 0)) = m(1,0); - return r; + PyObject* r = PyArray_SimpleNewFromData(2, dim, NPY_DOUBLE, m.ptr()); + return bp::handle<>(r); } static void wrap() { - static char const * doc = + static const char* doc = "Class to describe transformation from an ellipse\n" "with center x0, size exp(mu), and shape s to the unit circle.\n" "Map from source plane to image plane is defined as\n" @@ -152,7 +142,7 @@ struct PyCppEllipse { .def( "reset", ( void (CppEllipse::*)( - const CppShear &, double, const Position))&CppEllipse::reset, + const CppShear&, double, const Position))&CppEllipse::reset, bp::args("s", "mu", "p")) .def("fwd", &CppEllipse::fwd, "FIXME: needs documentation!") .def("inv", &CppEllipse::inv, "FIXME: needs documentation!") diff --git a/pysrc/Image.cpp b/pysrc/Image.cpp index 128d3be2cf1..a006a3b6088 100644 --- a/pysrc/Image.cpp +++ b/pysrc/Image.cpp @@ -79,8 +79,8 @@ struct PyImage { { CheckNumpyArray(array,2,isConst,data,owner,stride); bounds = Bounds( - xmin, xmin + PyArray_DIM(array.ptr(), 1) - 1, - ymin, ymin + PyArray_DIM(array.ptr(), 0) - 1 + xmin, xmin + GetNumpyArrayDim(array.ptr(), 1) - 1, + ymin, ymin + GetNumpyArrayDim(array.ptr(), 0) - 1 ); } diff --git a/pysrc/Noise.cpp b/pysrc/Noise.cpp index 3685170babc..2ac68cdaaac 100644 --- a/pysrc/Noise.cpp +++ b/pysrc/Noise.cpp @@ -20,10 +20,7 @@ */ #include "boost/python.hpp" #include "Noise.h" - -#define PY_ARRAY_UNIQUE_SYMBOL SBPROFILE_ARRAY_API -#define NO_IMPORT_ARRAY -#include "numpy/arrayobject.h" +#include "NumpyHelper.h" namespace bp = boost::python; diff --git a/pysrc/NumpyHelper.h b/pysrc/NumpyHelper.h index aee18a0f081..58fe473a03c 100644 --- a/pysrc/NumpyHelper.h +++ b/pysrc/NumpyHelper.h @@ -18,6 +18,8 @@ * You should have received a copy of the GNU General Public License * along with GalSim. If not, see */ +#ifndef NumpyHelper_H +#define NumpyHelper_H #include "boost/python.hpp" // header that includes Python.h always needs to come first @@ -31,10 +33,17 @@ #pragma warning (default : 47) #endif +#define NPY_NO_DEPRECATED_API NPY_1_7_API_VERSION #define PY_ARRAY_UNIQUE_SYMBOL GALSIM_ARRAY_API #define NO_IMPORT_ARRAY #include "numpy/arrayobject.h" +#if !defined(NPY_API_VERSION) || NPY_API_VERSION < 7 + #define NPY_OLD_API + #define NPY_ARRAY_ALIGNED NPY_ALIGNED + #define NPY_ARRAY_WRITEABLE NPY_WRITEABLE +#endif + namespace bp = boost::python; namespace galsim { @@ -46,8 +55,10 @@ template <> struct NumPyTraits { static int getCode() { return NPY_INT3 template <> struct NumPyTraits { static int getCode() { return NPY_FLOAT32; } }; template <> struct NumPyTraits { static int getCode() { return NPY_FLOAT64; } }; -static int Normalize(int code) +inline int GetNumpyArrayTypeCode(PyObject* array) { + PyArrayObject* numpy_array = reinterpret_cast(array); + int code = PyArray_TYPE(numpy_array); // Normally the return of PyArray_TYPE is a code that indicates what type the data is. However, // this gets confusing for integer types, since different integer types may be equivalent. In // particular int and long might be the same thing (typically on 32 bit machines, they can both @@ -62,23 +73,61 @@ static int Normalize(int code) // return the NumPy type for a C++ class (e.g. float -> numpy.float32) template -static bp::object GetNumPyType() +inline bp::object GetNumPyType() { bp::handle<> h(reinterpret_cast(PyArray_DescrFromType(NumPyTraits::getCode()))); return bp::object(h).attr("type"); } +inline int GetNumpyArrayNDim(PyObject* array) +{ + PyArrayObject* numpy_array = reinterpret_cast(array); + return PyArray_NDIM(numpy_array); +} + +inline int GetNumpyArrayDim(PyObject* array, int i) +{ + PyArrayObject* numpy_array = reinterpret_cast(array); + return PyArray_DIM(numpy_array,i); +} + +inline int GetNumpyArrayFlags(PyObject* array) +{ + PyArrayObject* numpy_array = reinterpret_cast(array); + return PyArray_FLAGS(numpy_array); +} + +inline PyObject* GetNumpyArrayBase(PyObject* array) +{ + PyArrayObject* numpy_array = reinterpret_cast(array); + return PyArray_BASE(numpy_array); +} + template -static void DestroyCObjectOwner(T* p) +inline int GetNumpyArrayStride(PyObject* array, int i) { - boost::shared_ptr * owner = reinterpret_cast< boost::shared_ptr *>(p); + PyArrayObject* numpy_array = reinterpret_cast(array); + return PyArray_STRIDE(numpy_array,i) / sizeof(T); +} + +template +inline T* GetNumpyArrayData(PyObject* array) +{ + PyArrayObject* numpy_array = reinterpret_cast(array); + return reinterpret_cast(PyArray_DATA(numpy_array)); +} + +template +inline void DestroyCObjectOwner(T* p) +{ + boost::shared_ptr* owner = reinterpret_cast*>(p); delete owner; } template struct PythonDeleter { - void operator()(T * p) { owner.reset(); } - explicit PythonDeleter(PyObject * o) : owner(bp::borrowed(o)) {} + void operator()(T* p) { owner.reset(); } + explicit PythonDeleter(PyObject* o) : owner(bp::borrowed(o)) {} bp::handle<> owner; }; @@ -88,18 +137,13 @@ static bp::object MakeNumpyArray( boost::shared_ptr owner = boost::shared_ptr()) { // --- Create array --- - int flags = NPY_ALIGNED; - if (!isConst) flags |= NPY_WRITEABLE; + int flags = NPY_ARRAY_ALIGNED; + if (!isConst) flags |= NPY_ARRAY_WRITEABLE; npy_intp shape[2] = { n1, n2 }; - npy_intp strides[2] = { stride * int(sizeof(T)), int(sizeof(T)) }; - bp::object result( - bp::handle<>( - PyArray_New( - &PyArray_Type, 2, shape, NumPyTraits::getCode(), strides, - const_cast(data), sizeof(T), flags, NULL - ) - ) - ); + npy_intp strides[2] = { stride* int(sizeof(T)), int(sizeof(T)) }; + PyObject* result = PyArray_New( + &PyArray_Type, 2, shape, NumPyTraits::getCode(), strides, + const_cast(data), sizeof(T), flags, NULL); // --- Manage ownership --- PythonDeleter* pyDeleter = boost::get_deleter >(owner); @@ -113,9 +157,13 @@ static bp::object MakeNumpyArray( PyCObject_FromVoidPtr(new boost::shared_ptr(owner), &DestroyCObjectOwner) ); } - reinterpret_cast(result.ptr())->base = pyOwner.release(); +#ifdef NPY_OLD_API + reinterpret_cast(result)->base = pyOwner.release(); +#else + PyArray_SetBaseObject(reinterpret_cast(result),pyOwner.release()); +#endif - return result; + return bp::object(bp::handle<>(result)); } template @@ -124,18 +172,13 @@ static bp::object MakeNumpyArray( boost::shared_ptr owner = boost::shared_ptr()) { // --- Create array --- - int flags = NPY_ALIGNED; - if (!isConst) flags |= NPY_WRITEABLE; + int flags = NPY_ARRAY_ALIGNED; + if (!isConst) flags |= NPY_ARRAY_WRITEABLE; npy_intp shape[1] = { n1 }; - npy_intp strides[1] = { stride * int(sizeof(T)) }; - bp::object result( - bp::handle<>( - PyArray_New( - &PyArray_Type, 1, shape, NumPyTraits::getCode(), strides, - const_cast(data), sizeof(T), flags, NULL - ) - ) - ); + npy_intp strides[1] = { stride* int(sizeof(T)) }; + PyObject* result = PyArray_New( + &PyArray_Type, 1, shape, NumPyTraits::getCode(), strides, + const_cast(data), sizeof(T), flags, NULL); // --- Manage ownership --- PythonDeleter* pyDeleter = boost::get_deleter >(owner); @@ -149,9 +192,13 @@ static bp::object MakeNumpyArray( PyCObject_FromVoidPtr(new boost::shared_ptr(owner), &DestroyCObjectOwner) ); } - reinterpret_cast(result.ptr())->base = pyOwner.release(); +#ifdef NPY_OLD_API + reinterpret_cast(result)->base = pyOwner.release(); +#else + PyArray_SetBaseObject(reinterpret_cast(result),pyOwner.release()); +#endif - return result; + return bp::object(bp::handle<>(result)); } // Check the type of the numpy array, input as array. @@ -168,7 +215,7 @@ static void CheckNumpyArray(const bp::object& array, int ndim, bool isConst, PyErr_SetString(PyExc_TypeError, "numpy.ndarray argument required"); bp::throw_error_already_set(); } - int actualType = Normalize(PyArray_TYPE(array.ptr())); + int actualType = GetNumpyArrayTypeCode(array.ptr()); int requiredType = NumPyTraits::getCode(); if (actualType != requiredType) { std::ostringstream oss; @@ -197,44 +244,40 @@ static void CheckNumpyArray(const bp::object& array, int ndim, bool isConst, PyErr_SetString(PyExc_ValueError, oss.str().c_str()); bp::throw_error_already_set(); } - if (PyArray_NDIM(array.ptr()) != ndim) { + if (GetNumpyArrayNDim(array.ptr()) != ndim) { PyErr_SetString(PyExc_ValueError, "numpy.ndarray argument has must be 2-d"); bp::throw_error_already_set(); } - if (!isConst && !(PyArray_FLAGS(array.ptr()) & NPY_WRITEABLE)) { + if (!isConst && !(GetNumpyArrayFlags(array.ptr()) & NPY_ARRAY_WRITEABLE)) { PyErr_SetString(PyExc_TypeError, "numpy.ndarray argument must be writeable"); bp::throw_error_already_set(); } - if (ndim == 2 && PyArray_STRIDE(array.ptr(), 1) != sizeof(T)) { + if (ndim == 2 && GetNumpyArrayStride(array.ptr(), 1) != 1) { PyErr_SetString(PyExc_ValueError, "numpy.ndarray argument must have contiguous rows"); bp::throw_error_already_set(); } - stride = PyArray_STRIDE(array.ptr(), 0) / sizeof(T); - data = reinterpret_cast(PyArray_DATA(array.ptr())); - PyObject * pyOwner = PyArray_BASE(array.ptr()); + stride = GetNumpyArrayStride(array.ptr(), 0); + data = GetNumpyArrayData(array.ptr()); + PyObject* pyOwner = GetNumpyArrayBase(array.ptr()); if (pyOwner) { - if (PyArray_Check(pyOwner) && PyArray_TYPE(pyOwner) == requiredType) { + if (PyArray_Check(pyOwner) && GetNumpyArrayTypeCode(pyOwner) == requiredType) { // Not really important, but we try to use the full array for // the owner pointer if this is a subarray, just to be consistent // with how it works for subimages. // The deleter is really all that matters. - owner = boost::shared_ptr( - reinterpret_cast(PyArray_DATA(pyOwner)), - PythonDeleter(pyOwner) - ); + owner = boost::shared_ptr(GetNumpyArrayData(pyOwner), + PythonDeleter(pyOwner)); } else { - owner = boost::shared_ptr( - reinterpret_cast(PyArray_DATA(array.ptr())), - PythonDeleter(pyOwner) - ); + owner = boost::shared_ptr(GetNumpyArrayData(array.ptr()), + PythonDeleter(pyOwner)); } } else { - owner = boost::shared_ptr( - reinterpret_cast(PyArray_DATA(array.ptr())), - PythonDeleter(array.ptr()) - ); + owner = boost::shared_ptr(GetNumpyArrayData(array.ptr()), + PythonDeleter(array.ptr())); } } } // namespace galsim + +#endif diff --git a/pysrc/SBShapelet.cpp b/pysrc/SBShapelet.cpp index 1ef84c5b9d0..6ecc72f7579 100644 --- a/pysrc/SBShapelet.cpp +++ b/pysrc/SBShapelet.cpp @@ -51,7 +51,7 @@ namespace galsim { boost::shared_ptr owner; int stride = 0; CheckNumpyArray(array,1,true,data,owner,stride); - int size = PyArray_DIM(array.ptr(), 0); + int size = GetNumpyArrayDim(array.ptr(), 0); if (size != PQIndex::size(order)) { PyErr_SetString(PyExc_ValueError, "Array for LVector is the wrong size"); bp::throw_error_already_set(); diff --git a/pysrc/module.cpp b/pysrc/module.cpp index 4dd2df3c771..c4c19dd0f6d 100644 --- a/pysrc/module.cpp +++ b/pysrc/module.cpp @@ -20,7 +20,9 @@ */ #include "boost/python.hpp" +#define NPY_NO_DEPRECATED_API NPY_1_7_API_VERSION #define PY_ARRAY_UNIQUE_SYMBOL GALSIM_ARRAY_API +// This is the only one that doesn't have NO_IMPORT_ARRAY. #include "numpy/arrayobject.h" namespace galsim { diff --git a/tests/test_InterpolatedImage.py b/tests/test_InterpolatedImage.py index 2620bed7e91..b11d4249183 100644 --- a/tests/test_InterpolatedImage.py +++ b/tests/test_InterpolatedImage.py @@ -243,22 +243,25 @@ def test_exceptions(): import time t1 = time.time() - # What if it receives as input something that is not an Image? Give it a GSObject to check. - g = galsim.Gaussian(sigma=1.) - np.testing.assert_raises(ValueError, galsim.InterpolatedImage, g) - # What if Image does not have a scale set, but dx keyword is not specified? - im = galsim.ImageF(5, 5) - np.testing.assert_raises(ValueError, galsim.InterpolatedImage, im) - # Image must have bounds defined - im = galsim.ImageF() - im.setScale(1.) - np.testing.assert_raises(ValueError, galsim.InterpolatedImage, im) - # Weird flux normalization - im = galsim.ImageF(5, 5) - im.setScale(1.) - np.testing.assert_raises(ValueError, galsim.InterpolatedImage, im, normalization = 'foo') - # Weird interpolant - give it something random like a GSObject - np.testing.assert_raises(RuntimeError, galsim.InterpolatedImage, im, x_interpolant = g) + try: + # What if it receives as input something that is not an Image? Give it a GSObject to check. + g = galsim.Gaussian(sigma=1.) + np.testing.assert_raises(ValueError, galsim.InterpolatedImage, g) + # What if Image does not have a scale set, but dx keyword is not specified? + im = galsim.ImageF(5, 5) + np.testing.assert_raises(ValueError, galsim.InterpolatedImage, im) + # Image must have bounds defined + im = galsim.ImageF() + im.setScale(1.) + np.testing.assert_raises(ValueError, galsim.InterpolatedImage, im) + # Weird flux normalization + im = galsim.ImageF(5, 5) + im.setScale(1.) + np.testing.assert_raises(ValueError, galsim.InterpolatedImage, im, normalization = 'foo') + # Weird interpolant - give it something random like a GSObject + np.testing.assert_raises(RuntimeError, galsim.InterpolatedImage, im, x_interpolant = g) + except ImportError: + print 'The assert_raises tests require nose' t2 = time.time() print 'time for %s = %.2f'%(funcname(),t2-t1)