[ VIGRA Homepage | Function Index | Class Index | Namespaces | File List | Main Page ]

vigra/numpy_array.hxx |

    00001 /************************************************************************/
    00002 /*                                                                      */
    00003 /*       Copyright 2009 by Ullrich Koethe and Hans Meine                */
    00004 /*                                                                      */
    00005 /*    This file is part of the VIGRA computer vision library.           */
    00006 /*    The VIGRA Website is                                              */
    00007 /*        http://hci.iwr.uni-heidelberg.de/vigra/                       */
    00008 /*    Please direct questions, bug reports, and contributions to        */
    00009 /*        ullrich.koethe@iwr.uni-heidelberg.de    or                    */
    00010 /*        vigra@informatik.uni-hamburg.de                               */
    00011 /*                                                                      */
    00012 /*    Permission is hereby granted, free of charge, to any person       */
    00013 /*    obtaining a copy of this software and associated documentation    */
    00014 /*    files (the "Software"), to deal in the Software without           */
    00015 /*    restriction, including without limitation the rights to use,      */
    00016 /*    copy, modify, merge, publish, distribute, sublicense, and/or      */
    00017 /*    sell copies of the Software, and to permit persons to whom the    */
    00018 /*    Software is furnished to do so, subject to the following          */
    00019 /*    conditions:                                                       */
    00020 /*                                                                      */
    00021 /*    The above copyright notice and this permission notice shall be    */
    00022 /*    included in all copies or substantial portions of the             */
    00023 /*    Software.                                                         */
    00024 /*                                                                      */
    00025 /*    THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND    */
    00026 /*    EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES   */
    00027 /*    OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND          */
    00028 /*    NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT       */
    00029 /*    HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,      */
    00030 /*    WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING      */
    00031 /*    FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR     */
    00032 /*    OTHER DEALINGS IN THE SOFTWARE.                                   */
    00033 /*                                                                      */
    00034 /************************************************************************/
    00035
    00036 #ifndef VIGRA_NUMPY_ARRAY_HXX
    00037 #define VIGRA_NUMPY_ARRAY_HXX
    00038
    00039 #include <Python.h>
    00040 #include <iostream>
    00041 #include <algorithm>
    00042 #include <complex>
    00043 #include <string>
    00044 #include <sstream>
    00045 #include <map>
    00046 #include <vigra/multi_array.hxx>
    00047 #include <vigra/array_vector.hxx>
    00048 #include <vigra/sized_int.hxx>
    00049 #include <vigra/python_utility.hxx>
    00050 #include <numpy/arrayobject.h>
    00051
    00052 int _import_array();
    00053
    00054 namespace vigra {
    00055
    00056 /********************************************************/
    00057 /*                                                      */
    00058 /*              Singleband and Multiband                */
    00059 /*                                                      */
    00060 /********************************************************/
    00061
    00062 typedef float NumpyValueType;
    00063
    00064 template <class T>
    00065 struct Singleband  // the last array dimension is not to be interpreted as a channel dimension
    00066 {
    00067     typedef T value_type;
    00068 };
    00069
    00070 template <class T>
    00071 struct Multiband  // the last array dimension is a channel dimension
    00072 {
    00073     typedef T value_type;
    00074 };
    00075
    00076 template<class T>
    00077 struct NumericTraits<Singleband<T> >
    00078 : public NumericTraits<T>
    00079 {};
    00080
    00081 template<class T>
    00082 struct NumericTraits<Multiband<T> >
    00083 {
    00084     typedef Multiband<T> Type;
    00085 /*
    00086     typedef int Promote;
    00087     typedef unsigned int UnsignedPromote;
    00088     typedef double RealPromote;
    00089     typedef std::complex<RealPromote> ComplexPromote;
    00090 */
    00091     typedef Type ValueType;
    00092
    00093     typedef typename NumericTraits<T>::isIntegral isIntegral;
    00094     typedef VigraFalseType isScalar;
    00095     typedef typename NumericTraits<T>::isSigned isSigned;
    00096     typedef typename NumericTraits<T>::isSigned isOrdered;
    00097     typedef typename NumericTraits<T>::isSigned isComplex;
    00098 /*
    00099     static signed char zero() { return 0; }
    00100     static signed char one() { return 1; }
    00101     static signed char nonZero() { return 1; }
    00102     static signed char min() { return SCHAR_MIN; }
    00103     static signed char max() { return SCHAR_MAX; }
    00104
    00105 #ifdef NO_INLINE_STATIC_CONST_DEFINITION
    00106     enum { minConst = SCHAR_MIN, maxConst = SCHAR_MIN };
    00107 #else
    00108     static const signed char minConst = SCHAR_MIN;
    00109     static const signed char maxConst = SCHAR_MIN;
    00110 #endif
    00111
    00112     static Promote toPromote(signed char v) { return v; }
    00113     static RealPromote toRealPromote(signed char v) { return v; }
    00114     static signed char fromPromote(Promote v) {
    00115         return ((v < SCHAR_MIN) ? SCHAR_MIN : (v > SCHAR_MAX) ? SCHAR_MAX : v);
    00116     }
    00117     static signed char fromRealPromote(RealPromote v) {
    00118         return ((v < 0.0)
    00119                    ? ((v < (RealPromote)SCHAR_MIN)
    00120                        ? SCHAR_MIN
    00121                        : static_cast<signed char>(v - 0.5))
    00122                    : (v > (RealPromote)SCHAR_MAX)
    00123                        ? SCHAR_MAX
    00124                        : static_cast<signed char>(v + 0.5));
    00125     }
    00126 */
    00127 };
    00128
    00129 template <class T>
    00130 class MultibandVectorAccessor
    00131 {
    00132     MultiArrayIndex size_, stride_;
    00133
    00134   public:
    00135     MultibandVectorAccessor(MultiArrayIndex size, MultiArrayIndex stride)
    00136     : size_(size),
    00137       stride_(stride)
    00138     {}
    00139
    00140
    00141     typedef Multiband<T> value_type;
    00142
    00143         /** the vector's value_type
    00144         */
    00145     typedef T component_type;
    00146
    00147     typedef VectorElementAccessor<MultibandVectorAccessor<T> > ElementAccessor;
    00148
    00149         /** Read the component data at given vector index
    00150             at given iterator position
    00151         */
    00152     template <class ITERATOR>
    00153     component_type const & getComponent(ITERATOR const & i, int idx) const
    00154     {
    00155         return *(&*i+idx*stride_);
    00156     }
    00157
    00158         /** Set the component data at given vector index
    00159             at given iterator position. The type <TT>V</TT> of the passed
    00160             in <TT>value</TT> is automatically converted to <TT>component_type</TT>.
    00161             In case of a conversion floating point -> intergral this includes rounding and clipping.
    00162         */
    00163     template <class V, class ITERATOR>
    00164     void setComponent(V const & value, ITERATOR const & i, int idx) const
    00165     {
    00166         *(&*i+idx*stride_) = detail::RequiresExplicitCast<component_type>::cast(value);
    00167     }
    00168
    00169         /** Read the component data at given vector index
    00170             at an offset of given iterator position
    00171         */
    00172     template <class ITERATOR, class DIFFERENCE>
    00173     component_type const & getComponent(ITERATOR const & i, DIFFERENCE const & diff, int idx) const
    00174     {
    00175         return *(&i[diff]+idx*stride_);
    00176     }
    00177
    00178     /** Set the component data at given vector index
    00179         at an offset of given iterator position. The type <TT>V</TT> of the passed
    00180         in <TT>value</TT> is automatically converted to <TT>component_type</TT>.
    00181             In case of a conversion floating point -> intergral this includes rounding and clipping.
    00182     */
    00183     template <class V, class ITERATOR, class DIFFERENCE>
    00184     void
    00185     setComponent(V const & value, ITERATOR const & i, DIFFERENCE const & diff, int idx) const
    00186     {
    00187         *(&i[diff]+idx*stride_) = detail::RequiresExplicitCast<component_type>::cast(value);
    00188     }
    00189
    00190     template <class U>
    00191     MultiArrayIndex size(U) const
    00192     {
    00193         return size_;
    00194     }
    00195 };
    00196
    00197 /********************************************************/
    00198 /*                                                      */
    00199 /*                a few Python utilities                */
    00200 /*                                                      */
    00201 /********************************************************/
    00202
    00203 namespace detail {
    00204
    00205 inline long spatialDimensions(PyObject * obj)
    00206 {
    00207     static python_ptr key(PyString_FromString("spatialDimensions"), python_ptr::keep_count);
    00208     python_ptr pres(PyObject_GetAttr(obj, key), python_ptr::keep_count);
    00209     long res = pres && PyInt_Check(pres)
    00210                  ? PyInt_AsLong(pres)
    00211                  : -1;
    00212     return res;
    00213 }
    00214
    00215 /*
    00216  * The registry is used to optionally map specific C++ types to
    00217  * specific python sub-classes of numpy.ndarray (for example,
    00218  * MultiArray<2, Singleband<int> > to a user-defined Python class 'ScalarImage').
    00219  *
    00220  * One needs to use NUMPY_ARRAY_INITIALIZE_REGISTRY once in a python
    00221  * extension module using this technique, in order to actually provide
    00222  * the registry (this is done by vigranumpycmodule and will then be
    00223  * available for other modules, too).  Alternatively,
    00224  * NUMPY_ARRAY_DUMMY_REGISTRY may be used to disable this feature
    00225  * completely.  In both cases, the macro must not be enclosed by any
    00226  * namespace, so it is best put right at the beginning of the file
    00227  * (e.g. below the #includes).
    00228  */
    00229
    00230 typedef std::map<std::string, std::pair<python_ptr, python_ptr> > ArrayTypeMap;
    00231
    00232 VIGRA_EXPORT ArrayTypeMap * getArrayTypeMap();
    00233
    00234 #define NUMPY_ARRAY_INITIALIZE_REGISTRY                                 \
    00235     namespace vigra { namespace detail {                                \
    00236     ArrayTypeMap * getArrayTypeMap()                                    \
    00237     {                                                                   \
    00238         static ArrayTypeMap arrayTypeMap;                               \
    00239         return &arrayTypeMap;                                           \
    00240     }                                                                   \
    00241     }} // namespace vigra::detail
    00242
    00243 #define NUMPY_ARRAY_DUMMY_REGISTRY                      \
    00244     namespace vigra { namespace detail {                \
    00245     ArrayTypeMap * getArrayTypeMap()                    \
    00246     {                                                   \
    00247         return NULL;                                    \
    00248     }                                                   \
    00249     }} // namespace vigra::detail
    00250
    00251 inline
    00252 void registerPythonArrayType(std::string const & name, PyObject * obj, PyObject * typecheck)
    00253 {
    00254     ArrayTypeMap *types = getArrayTypeMap();
    00255     vigra_precondition(
    00256         types != NULL,
    00257         "registerPythonArrayType(): module was compiled without array type registry.");
    00258     vigra_precondition(
    00259         obj && PyType_Check(obj) && PyType_IsSubtype((PyTypeObject *)obj, &PyArray_Type),
    00260         "registerPythonArrayType(obj): obj is not a subtype of numpy.ndarray.");
    00261     if(typecheck && PyCallable_Check(typecheck))
    00262         (*types)[name] = std::make_pair(python_ptr(obj), python_ptr(typecheck));
    00263     else
    00264         (*types)[name] = std::make_pair(python_ptr(obj), python_ptr());
    00265 //    std::cerr << "Registering " << ((PyTypeObject *)obj)->tp_name << " for " << name << "\n";
    00266 }
    00267
    00268 inline
    00269 python_ptr getArrayTypeObject(std::string const & name, PyTypeObject * def = 0)
    00270 {
    00271     ArrayTypeMap *types = getArrayTypeMap();
    00272     if(!types)
    00273         // dummy registry -> handle like empty registry
    00274         return python_ptr((PyObject *)def);
    00275
    00276     python_ptr res;
    00277     ArrayTypeMap::iterator i = types->find(name);
    00278     if(i != types->end())
    00279         res = i->second.first;
    00280     else
    00281         res = python_ptr((PyObject *)def);
    00282 //    std::cerr << "Requested " << name << ", got " << ((PyTypeObject *)res.get())->tp_name << "\n";
    00283     return res;
    00284 }
    00285
    00286 // there are two cases for the return:
    00287 // * if a typecheck function was registered, it is returned
    00288 // * a null pointer is returned if nothing was registered for either key, or if
    00289 //   a type was registered without typecheck function
    00290 inline python_ptr
    00291 getArrayTypecheckFunction(std::string const & keyFull, std::string const & key)
    00292 {
    00293     python_ptr res;
    00294     ArrayTypeMap *types = getArrayTypeMap();
    00295     if(types)
    00296     {
    00297         ArrayTypeMap::iterator i = types->find(keyFull);
    00298         if(i == types->end())
    00299             i = types->find(key);
    00300         if(i != types->end())
    00301             res = i->second.second;
    00302     }
    00303     return res;
    00304 }
    00305
    00306 inline bool
    00307 performCustomizedArrayTypecheck(PyObject * obj, std::string const & keyFull, std::string const & key)
    00308 {
    00309     if(obj == 0 || !PyArray_Check(obj))
    00310         return false;
    00311     python_ptr typecheck = getArrayTypecheckFunction(keyFull, key);
    00312     if(typecheck == 0)
    00313         return true; // no custom test registered
    00314     python_ptr args(PyTuple_Pack(1, obj), python_ptr::keep_count);
    00315     pythonToCppException(args);
    00316     python_ptr res(PyObject_Call(typecheck.get(), args.get(), 0), python_ptr::keep_count);
    00317     pythonToCppException(res);
    00318     vigra_precondition(PyBool_Check(res),
    00319            "NumpyArray conversion: registered typecheck function did not return a boolean.");
    00320     return (void*)res.get() == (void*)Py_True;
    00321 }
    00322
    00323 inline
    00324 python_ptr constructNumpyArrayImpl(
    00325     PyTypeObject * type,
    00326     ArrayVector<npy_intp> const & shape, npy_intp *strides,
    00327     NPY_TYPES typeCode, bool init)
    00328 {
    00329     python_ptr array;
    00330
    00331     if(strides == 0)
    00332     {
    00333         array = python_ptr(PyArray_New(type, shape.size(), (npy_intp *)shape.begin(), typeCode, 0, 0, 0, 1 /* Fortran order */, 0),
    00334                            python_ptr::keep_count);
    00335     }
    00336     else
    00337     {
    00338         int N = shape.size();
    00339         ArrayVector<npy_intp> pshape(N);
    00340         for(int k=0; k<N; ++k)
    00341             pshape[strides[k]] = shape[k];
    00342
    00343         array = python_ptr(PyArray_New(type, N, pshape.begin(), typeCode, 0, 0, 0, 1 /* Fortran order */, 0),
    00344                            python_ptr::keep_count);
    00345         pythonToCppException(array);
    00346
    00347         PyArray_Dims permute = { strides, N };
    00348         array = python_ptr(PyArray_Transpose((PyArrayObject*)array.get(), &permute), python_ptr::keep_count);
    00349     }
    00350     pythonToCppException(array);
    00351
    00352     if(init)
    00353         PyArray_FILLWBYTE((PyArrayObject *)array.get(), 0);
    00354
    00355     return array;
    00356 }
    00357
    00358 // strideOrdering will be ignored unless order == "A"
    00359 // TODO: this function should receive some refactoring in order to make
    00360 //       the rules clear from the code rather than from comments
    00361 inline python_ptr
    00362 constructNumpyArrayImpl(PyTypeObject * type, ArrayVector<npy_intp> const & shape,
    00363                        unsigned int spatialDimensions, unsigned int channels,
    00364                        NPY_TYPES typeCode, std::string order, bool init,
    00365                        ArrayVector<npy_intp> strideOrdering = ArrayVector<npy_intp>())
    00366 {
    00367     // shape must have at least length spatialDimensions, but can also have a channel dimension
    00368     vigra_precondition(shape.size() == spatialDimensions || shape.size() == spatialDimensions + 1,
    00369            "constructNumpyArray(type, shape, ...): shape has wrong length.");
    00370
    00371     // if strideOrdering is given, it must have at least length spatialDimensions,
    00372     // but can also have a channel dimension
    00373     vigra_precondition(strideOrdering.size() == 0 || strideOrdering.size() == spatialDimensions ||
    00374                        strideOrdering.size() == spatialDimensions + 1,
    00375            "constructNumpyArray(type, ..., strideOrdering): strideOrdering has wrong length.");
    00376
    00377     if(channels == 0) // if the requested number of channels is not given ...
    00378     {
    00379         // ... deduce it
    00380         if(shape.size() == spatialDimensions)
    00381             channels = 1;
    00382         else
    00383             channels = shape.back();
    00384     }
    00385     else
    00386     {
    00387         // otherwise, if the shape object also contains a channel dimension, they must be consistent
    00388         if(shape.size() > spatialDimensions)
    00389             vigra_precondition(channels == (unsigned int)shape[spatialDimensions],
    00390                    "constructNumpyArray(type, ...): shape contradicts requested number of channels.");
    00391     }
    00392
    00393     // if we have only one channel, no explicit channel dimension should be in the shape
    00394     unsigned int shapeSize = channels == 1
    00395                                   ? spatialDimensions
    00396                                   : spatialDimensions + 1;
    00397
    00398     // create the shape object with optional channel dimension
    00399     ArrayVector<npy_intp> pshape(shapeSize);
    00400     std::copy(shape.begin(), shape.begin()+std::min(shape.size(), pshape.size()), pshape.begin());
    00401     if(shapeSize > spatialDimensions)
    00402         pshape[spatialDimensions] = channels;
    00403
    00404     // order "A" means "preserve order" when an array is copied, and
    00405     // defaults to "V" when a new array is created without explicit strideOrdering
    00406     //
    00407     if(order == "A")
    00408     {
    00409         if(strideOrdering.size() == 0)
    00410         {
    00411             order = "V";
    00412         }
    00413         else if(strideOrdering.size() > shapeSize)
    00414         {
    00415             // make sure that strideOrdering length matches shape length
    00416             ArrayVector<npy_intp> pstride(strideOrdering.begin(), strideOrdering.begin()+shapeSize);
    00417
    00418             // adjust the ordering when the channel dimension has been dropped because channel == 1
    00419             if(strideOrdering[shapeSize] == 0)
    00420                 for(unsigned int k=0; k<shapeSize; ++k)
    00421                     pstride[k] -= 1;
    00422             pstride.swap(strideOrdering);
    00423         }
    00424         else if(strideOrdering.size() < shapeSize)
    00425         {
    00426             // make sure that strideOrdering length matches shape length
    00427             ArrayVector<npy_intp> pstride(shapeSize);
    00428
    00429             // adjust the ordering when the channel dimension has been dropped because channel == 1
    00430             for(unsigned int k=0; k<shapeSize-1; ++k)
    00431                 pstride[k] = strideOrdering[k] + 1;
    00432             pstride[shapeSize-1] = 0;
    00433             pstride.swap(strideOrdering);
    00434         }
    00435     }
    00436
    00437     // create the appropriate strideOrdering objects for the other memory orders
    00438     // (when strideOrdering already contained data, it is ignored because order != "A")
    00439     if(order == "C")
    00440     {
    00441         strideOrdering.resize(shapeSize);
    00442         for(unsigned int k=0; k<shapeSize; ++k)
    00443             strideOrdering[k] = shapeSize-1-k;
    00444     }
    00445     else if(order == "F" || (order == "V" && channels == 1))
    00446     {
    00447         strideOrdering.resize(shapeSize);
    00448         for(unsigned int k=0; k<shapeSize; ++k)
    00449             strideOrdering[k] = k;
    00450     }
    00451     else if(order == "V")
    00452     {
    00453         strideOrdering.resize(shapeSize);
    00454         for(unsigned int k=0; k<shapeSize-1; ++k)
    00455             strideOrdering[k] = k+1;
    00456         strideOrdering[shapeSize-1] = 0;
    00457     }
    00458
    00459     return constructNumpyArrayImpl(type, pshape, strideOrdering.begin(), typeCode, init);
    00460 }
    00461
    00462 template <class TINY_VECTOR>
    00463 inline
    00464 python_ptr constructNumpyArrayFromData(
    00465     std::string const & typeKeyFull,
    00466     std::string const & typeKey,
    00467     TINY_VECTOR const & shape, npy_intp *strides,
    00468     NPY_TYPES typeCode, void *data)
    00469 {
    00470     ArrayVector<npy_intp> pyShape(shape.begin(), shape.end());
    00471
    00472     python_ptr type = detail::getArrayTypeObject(typeKeyFull);
    00473     if(type == 0)
    00474         type = detail::getArrayTypeObject(typeKey, &PyArray_Type);
    00475
    00476     python_ptr array(PyArray_New((PyTypeObject *)type.ptr(), shape.size(), pyShape.begin(), typeCode, strides, data, 0, NPY_WRITEABLE, 0),
    00477                      python_ptr::keep_count);
    00478     pythonToCppException(array);
    00479
    00480     return array;
    00481 }
    00482
    00483
    00484 } // namespace detail
    00485
    00486 /********************************************************/
    00487 /*                                                      */
    00488 /*               NumpyArrayValuetypeTraits              */
    00489 /*                                                      */
    00490 /********************************************************/
    00491
    00492 template<class ValueType>
    00493 struct ERROR_NumpyArrayValuetypeTraits_not_specialized_for_ { };
    00494
    00495 template<class ValueType>
    00496 struct NumpyArrayValuetypeTraits
    00497 {
    00498     static bool isValuetypeCompatible(PyArrayObject const * obj)
    00499     {
    00500         return ERROR_NumpyArrayValuetypeTraits_not_specialized_for_<ValueType>();
    00501     }
    00502
    00503     static ERROR_NumpyArrayValuetypeTraits_not_specialized_for_<ValueType> typeCode;
    00504
    00505     static std::string typeName()
    00506     {
    00507         return std::string("ERROR: NumpyArrayValuetypeTraits not specialized for this case");
    00508     }
    00509
    00510     static std::string typeNameImpex()
    00511     {
    00512         return std::string("ERROR: NumpyArrayValuetypeTraits not specialized for this case");
    00513     }
    00514
    00515     static PyObject * typeObject()
    00516     {
    00517         return (PyObject *)0;
    00518     }
    00519 };
    00520
    00521 template<class ValueType>
    00522 ERROR_NumpyArrayValuetypeTraits_not_specialized_for_<ValueType> NumpyArrayValuetypeTraits<ValueType>::typeCode;
    00523
    00524 #define VIGRA_NUMPY_VALUETYPE_TRAITS(type, typeID, numpyTypeName, impexTypeName) \
    00525 template <> \
    00526 struct NumpyArrayValuetypeTraits<type > \
    00527 { \
    00528     static bool isValuetypeCompatible(PyArrayObject const * obj) /* obj must not be NULL */ \
    00529     { \
    00530         return PyArray_EquivTypenums(typeID, PyArray_DESCR((PyObject *)obj)->type_num) && \
    00531                PyArray_ITEMSIZE((PyObject *)obj) == sizeof(type); \
    00532     } \
    00533     \
    00534     static NPY_TYPES const typeCode = typeID; \
    00535     \
    00536     static std::string typeName() \
    00537     { \
    00538         return #numpyTypeName; \
    00539     } \
    00540     \
    00541     static std::string typeNameImpex() \
    00542     { \
    00543         return impexTypeName; \
    00544     } \
    00545     \
    00546     static PyObject * typeObject() \
    00547     { \
    00548         return PyArray_TypeObjectFromType(typeID); \
    00549     } \
    00550 };
    00551
    00552 VIGRA_NUMPY_VALUETYPE_TRAITS(bool,           NPY_BOOL, bool, "UINT8")
    00553 VIGRA_NUMPY_VALUETYPE_TRAITS(signed char,    NPY_INT8, int8, "INT16")
    00554 VIGRA_NUMPY_VALUETYPE_TRAITS(unsigned char,  NPY_UINT8, uint8, "UINT8")
    00555 VIGRA_NUMPY_VALUETYPE_TRAITS(short,          NPY_INT16, int16, "INT16")
    00556 VIGRA_NUMPY_VALUETYPE_TRAITS(unsigned short, NPY_UINT16, uint16, "UINT16")
    00557
    00558 #if VIGRA_BITSOF_LONG == 32
    00559 VIGRA_NUMPY_VALUETYPE_TRAITS(long,           NPY_INT32, int32, "INT32")
    00560 VIGRA_NUMPY_VALUETYPE_TRAITS(unsigned long,  NPY_UINT32, uint32, "UINT32")
    00561 #elif VIGRA_BITSOF_LONG == 64
    00562 VIGRA_NUMPY_VALUETYPE_TRAITS(long,           NPY_INT64, int64, "DOUBLE")
    00563 VIGRA_NUMPY_VALUETYPE_TRAITS(unsigned long,  NPY_UINT64, uint64, "DOUBLE")
    00564 #endif
    00565
    00566 #if VIGRA_BITSOF_INT == 32
    00567 VIGRA_NUMPY_VALUETYPE_TRAITS(int,            NPY_INT32, int32, "INT32")
    00568 VIGRA_NUMPY_VALUETYPE_TRAITS(unsigned int,   NPY_UINT32, uint32, "UINT32")
    00569 #elif VIGRA_BITSOF_INT == 64
    00570 VIGRA_NUMPY_VALUETYPE_TRAITS(int,            NPY_INT64, int64, "DOUBLE")
    00571 VIGRA_NUMPY_VALUETYPE_TRAITS(unsigned int,   NPY_UINT64, uint64, "DOUBLE")
    00572 #endif
    00573
    00574 #ifdef PY_LONG_LONG
    00575 # if VIGRA_BITSOF_LONG_LONG == 32
    00576 VIGRA_NUMPY_VALUETYPE_TRAITS(long long,            NPY_INT32, int32, "INT32")
    00577 VIGRA_NUMPY_VALUETYPE_TRAITS(unsigned long long,   NPY_UINT32, uint32, "UINT32")
    00578 # elif VIGRA_BITSOF_LONG_LONG == 64
    00579 VIGRA_NUMPY_VALUETYPE_TRAITS(long long,          NPY_INT64, int64, "DOUBLE")
    00580 VIGRA_NUMPY_VALUETYPE_TRAITS(unsigned long long, NPY_UINT64, uint64, "DOUBLE")
    00581 # endif
    00582 #endif
    00583
    00584 VIGRA_NUMPY_VALUETYPE_TRAITS(npy_float32, NPY_FLOAT32, float32, "FLOAT")
    00585 VIGRA_NUMPY_VALUETYPE_TRAITS(npy_float64, NPY_FLOAT64, float64, "DOUBLE")
    00586 #if NPY_SIZEOF_LONGDOUBLE != NPY_SIZEOF_DOUBLE
    00587 VIGRA_NUMPY_VALUETYPE_TRAITS(npy_longdouble, NPY_LONGDOUBLE, longdouble, "")
    00588 #endif
    00589 VIGRA_NUMPY_VALUETYPE_TRAITS(npy_cfloat, NPY_CFLOAT, complex64, "")
    00590 VIGRA_NUMPY_VALUETYPE_TRAITS(std::complex<npy_float>, NPY_CFLOAT, complex64, "")
    00591 VIGRA_NUMPY_VALUETYPE_TRAITS(npy_cdouble, NPY_CDOUBLE, complex128, "")
    00592 VIGRA_NUMPY_VALUETYPE_TRAITS(std::complex<npy_double>, NPY_CDOUBLE, complex128, "")
    00593 VIGRA_NUMPY_VALUETYPE_TRAITS(npy_clongdouble, NPY_CLONGDOUBLE, clongdouble, "")
    00594 #if NPY_SIZEOF_LONGDOUBLE != NPY_SIZEOF_DOUBLE
    00595 VIGRA_NUMPY_VALUETYPE_TRAITS(std::complex<npy_longdouble>, NPY_CLONGDOUBLE, clongdouble, "")
    00596 #endif
    00597
    00598 #undef VIGRA_NUMPY_VALUETYPE_TRAITS
    00599
    00600 /********************************************************/
    00601 /*                                                      */
    00602 /*                  NumpyArrayTraits                    */
    00603 /*                                                      */
    00604 /********************************************************/
    00605
    00606 template <class U, int N>
    00607 bool stridesAreAscending(TinyVector<U, N> const & strides)
    00608 {
    00609     for(int k=1; k<N; ++k)
    00610         if(strides[k] < strides[k-1])
    00611             return false;
    00612     return true;
    00613 }
    00614
    00615 template<unsigned int N, class T, class Stride>
    00616 struct NumpyArrayTraits;
    00617
    00618 template<unsigned int N, class T>
    00619 struct NumpyArrayTraits<N, T, StridedArrayTag>
    00620 {
    00621     typedef T dtype;
    00622     typedef T value_type;
    00623     typedef NumpyArrayValuetypeTraits<T> ValuetypeTraits;
    00624     static NPY_TYPES const typeCode = ValuetypeTraits::typeCode;
    00625
    00626     enum { spatialDimensions = N, channels = 1 };
    00627
    00628     static bool isArray(PyObject * obj)
    00629     {
    00630         return obj && PyArray_Check(obj);
    00631     }
    00632
    00633     static bool isClassCompatible(PyObject * obj)
    00634     {
    00635         return detail::performCustomizedArrayTypecheck(obj, typeKeyFull(), typeKey());
    00636     }
    00637
    00638     static bool isValuetypeCompatible(PyArrayObject * obj)  /* obj must not be NULL */
    00639     {
    00640         return ValuetypeTraits::isValuetypeCompatible(obj);
    00641     }
    00642
    00643     static bool isShapeCompatible(PyArrayObject * obj) /* obj must not be NULL */
    00644     {
    00645         return PyArray_NDIM((PyObject *)obj) == N-1 ||
    00646                PyArray_NDIM((PyObject *)obj) == N ||
    00647                (PyArray_NDIM((PyObject *)obj) == N+1 && PyArray_DIM((PyObject *)obj, N) == 1);
    00648     }
    00649
    00650     static bool isPropertyCompatible(PyArrayObject * obj) /* obj must not be NULL */
    00651     {
    00652         return ValuetypeTraits::isValuetypeCompatible(obj) &&
    00653                isShapeCompatible(obj);
    00654     }
    00655
    00656     template <class U>
    00657     static python_ptr constructor(TinyVector<U, N> const & shape,
    00658                                   T *data, TinyVector<U, N> const & stride)
    00659     {
    00660         TinyVector<npy_intp, N> npyStride(stride * sizeof(T));
    00661         return detail::constructNumpyArrayFromData(typeKeyFull(), typeKey(), shape, npyStride.begin(), ValuetypeTraits::typeCode, data);
    00662     }
    00663
    00664     static std::string typeKey()
    00665     {
    00666         static std::string key = std::string("NumpyArray<") + asString(N) + ", *>";
    00667         return key;
    00668     }
    00669
    00670     static std::string typeKeyFull()
    00671     {
    00672         static std::string key = std::string("NumpyArray<") + asString(N) + ", " +
    00673                                  ValuetypeTraits::typeName() + ", StridedArrayTag>";
    00674         return key;
    00675     }
    00676 };
    00677
    00678 /********************************************************/
    00679
    00680 template<unsigned int N, class T>
    00681 struct NumpyArrayTraits<N, T, UnstridedArrayTag>
    00682 : public NumpyArrayTraits<N, T, StridedArrayTag>
    00683 {
    00684     typedef NumpyArrayTraits<N, T, StridedArrayTag> BaseType;
    00685     typedef typename BaseType::ValuetypeTraits ValuetypeTraits;
    00686
    00687     static bool isShapeCompatible(PyArrayObject * obj) /* obj must not be NULL */
    00688     {
    00689         return BaseType::isShapeCompatible(obj) &&
    00690                PyArray_STRIDES((PyObject *)obj)[0] == PyArray_ITEMSIZE((PyObject *)obj);
    00691     }
    00692
    00693     static bool isPropertyCompatible(PyArrayObject * obj) /* obj must not be NULL */
    00694     {
    00695         return BaseType::isValuetypeCompatible(obj) &&
    00696                isShapeCompatible(obj);
    00697     }
    00698
    00699     template <class U>
    00700     static python_ptr constructor(TinyVector<U, N> const & shape,
    00701                                   T *data, TinyVector<U, N> const & stride)
    00702     {
    00703         TinyVector<npy_intp, N> npyStride(stride * sizeof(T));
    00704         return detail::constructNumpyArrayFromData(typeKeyFull(), BaseType::typeKey(), shape, npyStride.begin(), ValuetypeTraits::typeCode, data);
    00705     }
    00706
    00707     static std::string typeKeyFull()
    00708     {
    00709         static std::string key = std::string("NumpyArray<") + asString(N) + ", " +
    00710                                  ValuetypeTraits::typeName() + ", UnstridedArrayTag>";
    00711         return key;
    00712     }
    00713 };
    00714
    00715 /********************************************************/
    00716
    00717 template<unsigned int N, class T>
    00718 struct NumpyArrayTraits<N, Singleband<T>, StridedArrayTag>
    00719 : public NumpyArrayTraits<N, T, StridedArrayTag>
    00720 {
    00721     typedef NumpyArrayTraits<N, T, StridedArrayTag> BaseType;
    00722     typedef typename BaseType::ValuetypeTraits ValuetypeTraits;
    00723
    00724     static bool isClassCompatible(PyObject * obj)
    00725     {
    00726         return detail::performCustomizedArrayTypecheck(obj, typeKeyFull(), typeKey());
    00727     }
    00728
    00729     template <class U>
    00730     static python_ptr constructor(TinyVector<U, N> const & shape,
    00731                                   T *data, TinyVector<U, N> const & stride)
    00732     {
    00733         TinyVector<npy_intp, N> npyStride(stride * sizeof(T));
    00734         return detail::constructNumpyArrayFromData(typeKeyFull(), typeKey(), shape, npyStride.begin(), ValuetypeTraits::typeCode, data);
    00735     }
    00736
    00737     static std::string typeKey()
    00738     {
    00739         static std::string key = std::string("NumpyArray<") + asString(N) + ", Singleband<*> >";
    00740         return key;
    00741     }
    00742
    00743     static std::string typeKeyFull()
    00744     {
    00745         static std::string key = std::string("NumpyArray<") + asString(N) + ", Singleband<" +
    00746                                  ValuetypeTraits::typeName() + ">, StridedArrayTag>";
    00747         return key;
    00748     }
    00749 };
    00750
    00751 /********************************************************/
    00752
    00753 template<unsigned int N, class T>
    00754 struct NumpyArrayTraits<N, Singleband<T>, UnstridedArrayTag>
    00755 : public NumpyArrayTraits<N, Singleband<T>, StridedArrayTag>
    00756 {
    00757     typedef NumpyArrayTraits<N, T, UnstridedArrayTag> UnstridedTraits;
    00758     typedef NumpyArrayTraits<N, Singleband<T>, StridedArrayTag> BaseType;
    00759     typedef typename BaseType::ValuetypeTraits ValuetypeTraits;
    00760
    00761     static bool isShapeCompatible(PyArrayObject * obj) /* obj must not be NULL */
    00762     {
    00763         return UnstridedTraits::isShapeCompatible(obj);
    00764     }
    00765
    00766     static bool isPropertyCompatible(PyArrayObject * obj) /* obj must not be NULL */
    00767     {
    00768         return UnstridedTraits::isPropertyCompatible(obj);
    00769     }
    00770
    00771     template <class U>
    00772     static python_ptr constructor(TinyVector<U, N> const & shape,
    00773                                   T *data, TinyVector<U, N> const & stride)
    00774     {
    00775         TinyVector<npy_intp, N> npyStride(stride * sizeof(T));
    00776         return detail::constructNumpyArrayFromData(typeKeyFull(), BaseType::typeKey(), shape, npyStride.begin(), ValuetypeTraits::typeCode, data);
    00777     }
    00778
    00779     static std::string typeKeyFull()
    00780     {
    00781         static std::string key = std::string("NumpyArray<") + asString(N) + ", Singleband<" +
    00782                                  ValuetypeTraits::typeName() + ">, UnstridedArrayTag>";
    00783         return key;
    00784     }
    00785 };
    00786
    00787 /********************************************************/
    00788
    00789 template<unsigned int N, class T>
    00790 struct NumpyArrayTraits<N, Multiband<T>, StridedArrayTag>
    00791 : public NumpyArrayTraits<N, T, StridedArrayTag>
    00792 {
    00793     typedef NumpyArrayTraits<N, T, StridedArrayTag> BaseType;
    00794     typedef typename BaseType::ValuetypeTraits ValuetypeTraits;
    00795
    00796     enum { spatialDimensions = N-1, channels = 0 };
    00797
    00798     static bool isClassCompatible(PyObject * obj)
    00799     {
    00800         return detail::performCustomizedArrayTypecheck(obj, typeKeyFull(), typeKey());
    00801     }
    00802
    00803     static bool isShapeCompatible(PyArrayObject * obj) /* obj must not be NULL */
    00804     {
    00805         return PyArray_NDIM(obj) == N || PyArray_NDIM(obj) == N-1;
    00806     }
    00807
    00808     static bool isPropertyCompatible(PyArrayObject * obj) /* obj must not be NULL */
    00809     {
    00810         return ValuetypeTraits::isValuetypeCompatible(obj) &&
    00811                isShapeCompatible(obj);
    00812     }
    00813
    00814     template <class U>
    00815     static python_ptr constructor(TinyVector<U, N> const & shape,
    00816                                   T *data, TinyVector<U, N> const & stride)
    00817     {
    00818         TinyVector<npy_intp, N> npyStride(stride * sizeof(T));
    00819         return detail::constructNumpyArrayFromData(typeKeyFull(), typeKey(), shape, npyStride.begin(), ValuetypeTraits::typeCode, data);
    00820     }
    00821
    00822     static std::string typeKey()
    00823     {
    00824         static std::string key = std::string("NumpyArray<") + asString(N) + ", Multiband<*> >";
    00825         return key;
    00826     }
    00827
    00828     static std::string typeKeyFull()
    00829     {
    00830         static std::string key = std::string("NumpyArray<") + asString(N) + ", Multiband<" +
    00831                                  ValuetypeTraits::typeName() + ">, StridedArrayTag>";
    00832         return key;
    00833     }
    00834 };
    00835
    00836 /********************************************************/
    00837
    00838 template<unsigned int N, class T>
    00839 struct NumpyArrayTraits<N, Multiband<T>, UnstridedArrayTag>
    00840 : public NumpyArrayTraits<N, Multiband<T>, StridedArrayTag>
    00841 {
    00842     typedef NumpyArrayTraits<N, Multiband<T>, StridedArrayTag> BaseType;
    00843     typedef typename BaseType::ValuetypeTraits ValuetypeTraits;
    00844
    00845     static bool isShapeCompatible(PyArrayObject * obj) /* obj must not be NULL */
    00846     {
    00847         return BaseType::isShapeCompatible(obj) &&
    00848                PyArray_STRIDES((PyObject *)obj)[0] == PyArray_ITEMSIZE((PyObject *)obj);
    00849     }
    00850
    00851     static bool isPropertyCompatible(PyArrayObject * obj) /* obj must not be NULL */
    00852     {
    00853         return BaseType::isValuetypeCompatible(obj) &&
    00854                isShapeCompatible(obj);
    00855     }
    00856
    00857     template <class U>
    00858     static python_ptr constructor(TinyVector<U, N> const & shape,
    00859                                   T *data, TinyVector<U, N> const & stride)
    00860     {
    00861         TinyVector<npy_intp, N> npyStride(stride * sizeof(T));
    00862         return detail::constructNumpyArrayFromData(typeKeyFull(), BaseType::typeKey(), shape, npyStride.begin(), ValuetypeTraits::typeCode, data);
    00863     }
    00864
    00865     static std::string typeKeyFull()
    00866     {
    00867         static std::string key = std::string("NumpyArray<") + asString(N) + ", Multiband<" +
    00868                                  ValuetypeTraits::typeName() + ">, UnstridedArrayTag>";
    00869         return key;
    00870     }
    00871 };
    00872
    00873 /********************************************************/
    00874
    00875 template<unsigned int N, int M, class T>
    00876 struct NumpyArrayTraits<N, TinyVector<T, M>, StridedArrayTag>
    00877 {
    00878     typedef T dtype;
    00879     typedef TinyVector<T, M> value_type;
    00880     typedef NumpyArrayValuetypeTraits<T> ValuetypeTraits;
    00881     static NPY_TYPES const typeCode = ValuetypeTraits::typeCode;
    00882
    00883     enum { spatialDimensions = N, channels = M };
    00884
    00885     static bool isArray(PyObject * obj)
    00886     {
    00887         return obj && PyArray_Check(obj);
    00888     }
    00889
    00890     static bool isClassCompatible(PyObject * obj)
    00891     {
    00892         return detail::performCustomizedArrayTypecheck(obj, typeKeyFull(), typeKey());
    00893     }
    00894
    00895     static bool isValuetypeCompatible(PyArrayObject * obj)  /* obj must not be NULL */
    00896     {
    00897         return ValuetypeTraits::isValuetypeCompatible(obj);
    00898     }
    00899
    00900     static bool isShapeCompatible(PyArrayObject * obj) /* obj must not be NULL */
    00901     {
    00902         return PyArray_NDIM((PyObject *)obj) == N+1 &&
    00903                PyArray_DIM((PyObject *)obj, N) == M &&
    00904                PyArray_STRIDES((PyObject *)obj)[N] == PyArray_ITEMSIZE((PyObject *)obj);
    00905     }
    00906
    00907     static bool isPropertyCompatible(PyArrayObject * obj) /* obj must not be NULL */
    00908     {
    00909         return ValuetypeTraits::isValuetypeCompatible(obj) &&
    00910                isShapeCompatible(obj);
    00911     }
    00912
    00913     template <class U>
    00914     static python_ptr constructor(TinyVector<U, N> const & shape,
    00915                                   T *data, TinyVector<U, N> const & stride)
    00916     {
    00917         TinyVector<npy_intp, N+1> npyShape;
    00918         std::copy(shape.begin(), shape.end(), npyShape.begin());
    00919         npyShape[N] = M;
    00920
    00921         TinyVector<npy_intp, N+1> npyStride;
    00922         std::transform(
    00923             stride.begin(), stride.end(), npyStride.begin(),
    00924             std::bind2nd(std::multiplies<npy_intp>(), sizeof(value_type)));
    00925         npyStride[N] = sizeof(T);
    00926
    00927         return detail::constructNumpyArrayFromData(
    00928             typeKeyFull(), typeKey(), npyShape,
    00929             npyStride.begin(), ValuetypeTraits::typeCode, data);
    00930     }
    00931
    00932     static std::string typeKey()
    00933     {
    00934         static std::string key = std::string("NumpyArray<") + asString(N) + ", TinyVector<*, " + asString(M) + "> >";
    00935         return key;
    00936     }
    00937
    00938     static std::string typeKeyFull()
    00939     {
    00940         static std::string key = std::string("NumpyArray<") + asString(N) +
    00941                       ", TinyVector<" + ValuetypeTraits::typeName() + ", " + asString(M) + ">, StridedArrayTag>";
    00942         return key;
    00943     }
    00944 };
    00945
    00946 /********************************************************/
    00947
    00948 template<unsigned int N, int M, class T>
    00949 struct NumpyArrayTraits<N, TinyVector<T, M>, UnstridedArrayTag>
    00950 : public NumpyArrayTraits<N, TinyVector<T, M>, StridedArrayTag>
    00951 {
    00952     typedef NumpyArrayTraits<N, TinyVector<T, M>, StridedArrayTag> BaseType;
    00953     typedef typename BaseType::value_type value_type;
    00954     typedef typename BaseType::ValuetypeTraits ValuetypeTraits;
    00955
    00956     static bool isShapeCompatible(PyArrayObject * obj) /* obj must not be NULL */
    00957     {
    00958         return BaseType::isShapeCompatible(obj) &&
    00959                PyArray_STRIDES((PyObject *)obj)[0] == sizeof(TinyVector<T, M>);
    00960     }
    00961
    00962     static bool isPropertyCompatible(PyArrayObject * obj) /* obj must not be NULL */
    00963     {
    00964         return BaseType::isValuetypeCompatible(obj) &&
    00965                isShapeCompatible(obj);
    00966     }
    00967
    00968     template <class U>
    00969     static python_ptr constructor(TinyVector<U, N> const & shape,
    00970                                   T *data, TinyVector<U, N> const & stride)
    00971     {
    00972         TinyVector<npy_intp, N+1> npyShape;
    00973         std::copy(shape.begin(), shape.end(), npyShape.begin());
    00974         npyShape[N] = M;
    00975
    00976         TinyVector<npy_intp, N+1> npyStride;
    00977         std::transform(
    00978             stride.begin(), stride.end(), npyStride.begin(),
    00979             std::bind2nd(std::multiplies<npy_intp>(), sizeof(value_type)));
    00980         npyStride[N] = sizeof(T);
    00981
    00982         return detail::constructNumpyArrayFromData(
    00983             typeKeyFull(), BaseType::typeKey(), npyShape,
    00984             npyStride.begin(), ValuetypeTraits::typeCode, data);
    00985     }
    00986
    00987     static std::string typeKeyFull()
    00988     {
    00989         static std::string key = std::string("NumpyArray<") + asString(N) +
    00990                       ", TinyVector<" + ValuetypeTraits::typeName() + ", " + asString(M) + ">, UnstridedArrayTag>";
    00991         return key;
    00992     }
    00993 };
    00994
    00995 /********************************************************/
    00996
    00997 template<unsigned int N, class T>
    00998 struct NumpyArrayTraits<N, RGBValue<T>, StridedArrayTag>
    00999 : public NumpyArrayTraits<N, TinyVector<T, 3>, StridedArrayTag>
    01000 {
    01001     typedef T dtype;
    01002     typedef RGBValue<T> value_type;
    01003     typedef NumpyArrayValuetypeTraits<T> ValuetypeTraits;
    01004
    01005     static bool isClassCompatible(PyObject * obj)
    01006     {
    01007         return detail::performCustomizedArrayTypecheck(obj, typeKeyFull(), typeKey());
    01008     }
    01009
    01010     template <class U>
    01011     static python_ptr constructor(TinyVector<U, N> const & shape,
    01012                                   T *data, TinyVector<U, N> const & stride)
    01013     {
    01014         TinyVector<npy_intp, N+1> npyShape;
    01015         std::copy(shape.begin(), shape.end(), npyShape.begin());
    01016         npyShape[N] = 3;
    01017
    01018         TinyVector<npy_intp, N+1> npyStride;
    01019         std::transform(
    01020             stride.begin(), stride.end(), npyStride.begin(),
    01021             std::bind2nd(std::multiplies<npy_intp>(), sizeof(value_type)));
    01022         npyStride[N] = sizeof(T);
    01023
    01024         return detail::constructNumpyArrayFromData(
    01025             typeKeyFull(), typeKey(), npyShape,
    01026             npyStride.begin(), ValuetypeTraits::typeCode, data);
    01027     }
    01028
    01029     static std::string typeKey()
    01030     {
    01031         static std::string key = std::string("NumpyArray<") + asString(N) + ", RGBValue<*> >";
    01032         return key;
    01033     }
    01034
    01035     static std::string typeKeyFull()
    01036     {
    01037         static std::string key = std::string("NumpyArray<") + asString(N) +
    01038                       ", RGBValue<" + ValuetypeTraits::typeName() + ">, StridedArrayTag>";
    01039         return key;
    01040     }
    01041 };
    01042
    01043 /********************************************************/
    01044
    01045 template<unsigned int N, class T>
    01046 struct NumpyArrayTraits<N, RGBValue<T>, UnstridedArrayTag>
    01047 : public NumpyArrayTraits<N, RGBValue<T>, StridedArrayTag>
    01048 {
    01049     typedef NumpyArrayTraits<N, TinyVector<T, 3>, UnstridedArrayTag> UnstridedTraits;
    01050     typedef NumpyArrayTraits<N, RGBValue<T>, StridedArrayTag> BaseType;
    01051     typedef typename BaseType::value_type value_type;
    01052     typedef typename BaseType::ValuetypeTraits ValuetypeTraits;
    01053
    01054     static bool isShapeCompatible(PyArrayObject * obj) /* obj must not be NULL */
    01055     {
    01056         return UnstridedTraits::isShapeCompatible(obj);
    01057     }
    01058
    01059     static bool isPropertyCompatible(PyArrayObject * obj) /* obj must not be NULL */
    01060     {
    01061         return UnstridedTraits::isPropertyCompatible(obj);
    01062     }
    01063
    01064     template <class U>
    01065     static python_ptr constructor(TinyVector<U, N> const & shape,
    01066                                   T *data, TinyVector<U, N> const & stride)
    01067     {
    01068         TinyVector<npy_intp, N+1> npyShape;
    01069         std::copy(shape.begin(), shape.end(), npyShape.begin());
    01070         npyShape[N] = 3;
    01071
    01072         TinyVector<npy_intp, N+1> npyStride;
    01073         std::transform(
    01074             stride.begin(), stride.end(), npyStride.begin(),
    01075             std::bind2nd(std::multiplies<npy_intp>(), sizeof(value_type)));
    01076         npyStride[N] = sizeof(T);
    01077
    01078         return detail::constructNumpyArrayFromData(
    01079             typeKeyFull(), BaseType::typeKey(), npyShape,
    01080             npyStride.begin(), ValuetypeTraits::typeCode, data);
    01081     }
    01082
    01083     static std::string typeKeyFull()
    01084     {
    01085         static std::string key = std::string("NumpyArray<") + asString(N) +
    01086                       ", RGBValue<" + ValuetypeTraits::typeName() + ">, UnstridedArrayTag>";
    01087         return key;
    01088     }
    01089 };
    01090
    01091 /********************************************************/
    01092 /*                                                      */
    01093 /*                    NumpyAnyArray                     */
    01094 /*                                                      */
    01095 /********************************************************/
    01096
    01097 /** Wrapper class for a Python array.
    01098
    01099     This class stores a reference-counted pointer to an Python numpy array object,
    01100     i.e. an object where <tt>PyArray_Check(object)</tt> returns true (in Python, the
    01101     object is then a subclass of <tt>numpy.ndarray</tt>). This class is mainly used
    01102     as a smart pointer to these arrays, but some basic access and conversion functions
    01103     are also provided.
    01104
    01105     <b>\#include</b> <<a href="numpy__array_8hxx-source.html">vigra/numpy_array.hxx</a>><br>
    01106     Namespace: vigra
    01107 */
    01108 class NumpyAnyArray
    01109 {
    01110   protected:
    01111     python_ptr pyArray_;
    01112
    01113     // We want to apply broadcasting to the channel dimension.
    01114     // Since only leading dimensions can be added during numpy
    01115     // broadcasting, we permute the array accordingly.
    01116     NumpyAnyArray permuteChannelsToFront() const
    01117     {
    01118         MultiArrayIndex M = ndim();
    01119         ArrayVector<npy_intp> permutation(M);
    01120         for(int k=0; k<M; ++k)
    01121             permutation[k] = M-1-k;
    01122         // explicit cast to int is neede here to avoid gcc c++0x compilation
    01123         // error: narrowing conversion of ‘M’ from ‘vigra::MultiArrayIndex’
    01124         //        to ‘int’ inside { }
    01125         // int overflow should not occur here because PyArray_NDIM returns
    01126         // an integer which is converted to long in NumpyAnyArray::ndim()
    01127         PyArray_Dims permute = { permutation.begin(), (int) M };
    01128         python_ptr array(PyArray_Transpose(pyArray(), &permute), python_ptr::keep_count);
    01129         pythonToCppException(array);
    01130         return NumpyAnyArray(array.ptr());
    01131     }
    01132
    01133   public:
    01134
    01135         /// difference type
    01136     typedef ArrayVector<npy_intp> difference_type;
    01137
    01138         /**
    01139          Construct from a Python object. If \a obj is NULL, or is not a subclass
    01140          of numpy.ndarray, the resulting NumpyAnyArray will have no data (i.e.
    01141          hasData() returns false). Otherwise, it creates a new reference to the array
    01142          \a obj, unless \a createCopy is true, where a new array is created by calling
    01143          the C-equivalent of obj->copy().
    01144          */
    01145     explicit NumpyAnyArray(PyObject * obj = 0, bool createCopy = false, PyTypeObject * type = 0)
    01146     {
    01147         if(obj == 0)
    01148             return;
    01149         vigra_precondition(type == 0 || PyType_IsSubtype(type, &PyArray_Type),
    01150              "NumpyAnyArray(obj, createCopy, type): type must be numpy.ndarray or a subclass thereof.");
    01151         if(createCopy)
    01152             makeCopy(obj, type);
    01153         else
    01154             vigra_precondition(makeReference(obj, type), "NumpyAnyArray(obj): obj isn't a numpy array.");
    01155     }
    01156
    01157         /**
    01158          Copy constructor. By default, it creates a new reference to the array
    01159          \a other. When \a createCopy is true, a new array is created by calling
    01160          the C-equivalent of other.copy().
    01161          */
    01162     NumpyAnyArray(NumpyAnyArray const & other, bool createCopy = false, PyTypeObject * type = 0)
    01163     {
    01164         if(!other.hasData())
    01165             return;
    01166         vigra_precondition(type == 0 || PyType_IsSubtype(type, &PyArray_Type),
    01167              "NumpyAnyArray(obj, createCopy, type): type must be numpy.ndarray or a subclass thereof.");
    01168         if(createCopy)
    01169             makeCopy(other.pyObject(), type);
    01170         else
    01171             makeReference(other.pyObject(), type);
    01172     }
    01173
    01174     // auto-generated destructor is ok
    01175
    01176         /**
    01177          * Assignment operator. If this is already a view with data
    01178          * (i.e. hasData() is true) and the shapes match, the RHS
    01179          * array contents are copied via the C-equivalent of
    01180          * 'self[...] = other[...]'. If the shapes don't matched,
    01181          * broadcasting is tried on the trailing (i.e. channel)
    01182          * dimension.
    01183          * If the LHS is an empty view, assignment is identical to
    01184          * makeReference(other.pyObject()).
    01185          */
    01186     NumpyAnyArray & operator=(NumpyAnyArray const & other)
    01187     {
    01188         if(hasData())
    01189         {
    01190             vigra_precondition(other.hasData(),
    01191                 "NumpyArray::operator=(): Cannot assign from empty array.");
    01192             if(PyArray_CopyInto(permuteChannelsToFront().pyArray(), other.permuteChannelsToFront().pyArray()) == -1)
    01193                 pythonToCppException(0);
    01194         }
    01195         else
    01196         {
    01197             pyArray_ = other.pyArray_;
    01198         }
    01199         return *this;
    01200     }
    01201
    01202         /**
    01203          Returns the number of dimensions of this array, or 0 if
    01204          hasData() is false.
    01205          */
    01206     MultiArrayIndex ndim() const
    01207     {
    01208         if(hasData())
    01209             return PyArray_NDIM(pyObject());
    01210         return 0;
    01211     }
    01212
    01213         /**
    01214          Returns the number of spatial dimensions of this array, or 0 if
    01215          hasData() is false. If the enclosed Python array does not define
    01216          the attribute spatialDimensions, ndim() is returned.
    01217          */
    01218     MultiArrayIndex spatialDimensions() const
    01219     {
    01220         if(!hasData())
    01221             return 0;
    01222         MultiArrayIndex s = detail::spatialDimensions(pyObject());
    01223         if(s == -1)
    01224             s = ndim();
    01225         return s;
    01226     }
    01227
    01228         /**
    01229          Returns the shape of this array. The size of
    01230          the returned shape equals ndim().
    01231          */
    01232     difference_type shape() const
    01233     {
    01234         if(hasData())
    01235             return difference_type(PyArray_DIMS(pyObject()), PyArray_DIMS(pyObject()) + ndim());
    01236         return difference_type();
    01237     }
    01238
    01239         /** Compute the ordering of the strides of this array.
    01240             The result is describes the current permutation of the axes relative
    01241             to an ascending stride order.
    01242         */
    01243     difference_type strideOrdering() const
    01244     {
    01245         if(!hasData())
    01246             return difference_type();
    01247         MultiArrayIndex N = ndim();
    01248         difference_type stride(PyArray_STRIDES(pyObject()), PyArray_STRIDES(pyObject()) + N),
    01249                         permutation(N);
    01250         for(MultiArrayIndex k=0; k<N; ++k)
    01251             permutation[k] = k;
    01252         for(MultiArrayIndex k=0; k<N-1; ++k)
    01253         {
    01254             MultiArrayIndex smallest = k;
    01255             for(MultiArrayIndex j=k+1; j<N; ++j)
    01256             {
    01257                 if(stride[j] < stride[smallest])
    01258                     smallest = j;
    01259             }
    01260             if(smallest != k)
    01261             {
    01262                 std::swap(stride[k], stride[smallest]);
    01263                 std::swap(permutation[k], permutation[smallest]);
    01264             }
    01265         }
    01266         difference_type ordering(N);
    01267         for(MultiArrayIndex k=0; k<N; ++k)
    01268             ordering[permutation[k]] = k;
    01269         return ordering;
    01270     }
    01271
    01272         /**
    01273          Returns the value type of the elements in this array, or -1
    01274          when hasData() is false.
    01275          */
    01276     int dtype() const
    01277     {
    01278         if(hasData())
    01279             return PyArray_DESCR(pyObject())->type_num;
    01280         return -1;
    01281     }
    01282
    01283         /**
    01284          * Return a borrowed reference to the internal PyArrayObject.
    01285          */
    01286     PyArrayObject * pyArray() const
    01287     {
    01288         return (PyArrayObject *)pyArray_.get();
    01289     }
    01290
    01291         /**
    01292          * Return a borrowed reference to the internal PyArrayObject
    01293          * (see pyArray()), cast to PyObject for your convenience.
    01294          */
    01295     PyObject * pyObject() const
    01296     {
    01297         return pyArray_.get();
    01298     }
    01299
    01300         /**
    01301            Reset the NumpyAnyArray to the given object. If \a obj is a numpy array object,
    01302            a new reference to that array is created, and the function returns
    01303            true. Otherwise, it returns false and the NumpyAnyArray remains unchanged.
    01304            If \a type is given, the new reference will be a view with that type, provided
    01305            that \a type is a numpy ndarray or a subclass thereof. Otherwise, an
    01306            exception is thrown.
    01307          */
    01308     bool makeReference(PyObject * obj, PyTypeObject * type = 0)
    01309     {
    01310         if(obj == 0 || !PyArray_Check(obj))
    01311             return false;
    01312         if(type != 0)
    01313         {
    01314             vigra_precondition(PyType_IsSubtype(type, &PyArray_Type) != 0,
    01315                 "NumpyAnyArray::makeReference(obj, type): type must be numpy.ndarray or a subclass thereof.");
    01316             obj = PyArray_View((PyArrayObject*)obj, 0, type);
    01317             pythonToCppException(obj);
    01318         }
    01319         pyArray_.reset(obj);
    01320         return true;
    01321     }
    01322
    01323         /**
    01324            Create a copy of the given array object. If \a obj is a numpy array object,
    01325            a copy is created via the C-equivalent of 'obj->copy()'. If
    01326            this call fails, or obj was not an array, an exception is thrown
    01327            and the NumpyAnyArray remains unchanged.
    01328          */
    01329     void makeCopy(PyObject * obj, PyTypeObject * type = 0)
    01330     {
    01331         vigra_precondition(obj && PyArray_Check(obj),
    01332              "NumpyAnyArray::makeCopy(obj): obj is not an array.");
    01333         vigra_precondition(type == 0 || PyType_IsSubtype(type, &PyArray_Type),
    01334              "NumpyAnyArray::makeCopy(obj, type): type must be numpy.ndarray or a subclass thereof.");
    01335         python_ptr array(PyArray_NewCopy((PyArrayObject*)obj, NPY_ANYORDER), python_ptr::keep_count);
    01336         pythonToCppException(array);
    01337         makeReference(array, type);
    01338     }
    01339
    01340          /**
    01341            Check whether this NumpyAnyArray actually points to a Python array.
    01342          */
    01343     bool hasData() const
    01344     {
    01345         return pyArray_ != 0;
    01346     }
    01347 };
    01348
    01349 /********************************************************/
    01350 /*                                                      */
    01351 /*                     NumpyArray                       */
    01352 /*                                                      */
    01353 /********************************************************/
    01354
    01355 /** Provide the MultiArrayView interface for a Python array.
    01356
    01357     This class inherits from both \ref vigra::MultiArrayView and \ref vigra::NumpyAnyArray
    01358     in order to support easy and save application of VIGRA functions to Python arrays.
    01359
    01360     <b>\#include</b> <<a href="numpy__array_8hxx-source.html">vigra/numpy_array.hxx</a>><br>
    01361     Namespace: vigra
    01362 */
    01363 template <unsigned int N, class T, class Stride = StridedArrayTag>
    01364 class NumpyArray
    01365 : public MultiArrayView<N, typename NumpyArrayTraits<N, T, Stride>::value_type, Stride>,
    01366   public NumpyAnyArray
    01367 {
    01368   public:
    01369     typedef NumpyArrayTraits<N, T, Stride> ArrayTraits;
    01370     typedef typename ArrayTraits::dtype dtype;
    01371     typedef T pseudo_value_type;
    01372
    01373     static NPY_TYPES const typeCode = ArrayTraits::typeCode;
    01374
    01375         /** the view type associated with this array.
    01376          */
    01377     typedef MultiArrayView<N, typename ArrayTraits::value_type, Stride> view_type;
    01378
    01379     enum { actual_dimension = view_type::actual_dimension };
    01380
    01381         /** the array's value type
    01382          */
    01383     typedef typename view_type::value_type value_type;
    01384
    01385         /** pointer type
    01386          */
    01387     typedef typename view_type::pointer pointer;
    01388
    01389         /** const pointer type
    01390          */
    01391     typedef typename view_type::const_pointer const_pointer;
    01392
    01393         /** reference type (result of operator[])
    01394          */
    01395     typedef typename view_type::reference reference;
    01396
    01397         /** const reference type (result of operator[] const)
    01398          */
    01399     typedef typename view_type::const_reference const_reference;
    01400
    01401         /** size type
    01402          */
    01403     typedef typename view_type::size_type size_type;
    01404
    01405         /** difference type (used for multi-dimensional offsets and indices)
    01406          */
    01407     typedef typename view_type::difference_type difference_type;
    01408
    01409         /** difference and index type for a single dimension
    01410          */
    01411     typedef typename view_type::difference_type_1 difference_type_1;
    01412
    01413         /** traverser type
    01414          */
    01415     typedef typename view_type::traverser traverser;
    01416
    01417         /** traverser type to const data
    01418          */
    01419     typedef typename view_type::const_traverser const_traverser;
    01420
    01421         /** sequential (random access) iterator type
    01422          */
    01423     typedef value_type * iterator;
    01424
    01425         /** sequential (random access) const iterator type
    01426          */
    01427     typedef value_type * const_iterator;
    01428
    01429     using view_type::shape;   // resolve ambiguity of multiple inheritance
    01430     using view_type::hasData; // resolve ambiguity of multiple inheritance
    01431     using view_type::strideOrdering; // resolve ambiguity of multiple inheritance
    01432
    01433   protected:
    01434
    01435     // this function assumes that pyArray_ has already been set, and compatibility been checked
    01436     void setupArrayView();
    01437
    01438     static python_ptr getArrayTypeObject()
    01439     {
    01440         python_ptr type = detail::getArrayTypeObject(ArrayTraits::typeKeyFull());
    01441         if(type == 0)
    01442             type = detail::getArrayTypeObject(ArrayTraits::typeKey(), &PyArray_Type);
    01443         return type;
    01444     }
    01445
    01446     static python_ptr init(difference_type const & shape, bool init = true)
    01447     {
    01448         ArrayVector<npy_intp> pshape(shape.begin(), shape.end());
    01449         return detail::constructNumpyArrayImpl((PyTypeObject *)getArrayTypeObject().ptr(), pshape,
    01450                        ArrayTraits::spatialDimensions, ArrayTraits::channels,
    01451                        typeCode, "V", init);
    01452     }
    01453
    01454     static python_ptr init(difference_type const & shape, difference_type const & strideOrdering, bool init = true)
    01455     {
    01456         ArrayVector<npy_intp> pshape(shape.begin(), shape.end()),
    01457                               pstrideOrdering(strideOrdering.begin(), strideOrdering.end());
    01458         return detail::constructNumpyArrayImpl((PyTypeObject *)getArrayTypeObject().ptr(), pshape,
    01459                        ArrayTraits::spatialDimensions, ArrayTraits::channels,
    01460                        typeCode, "A", init, pstrideOrdering);
    01461     }
    01462
    01463   public:
    01464
    01465     using view_type::init;
    01466
    01467         /**
    01468          * Construct from a given PyObject pointer. When the given
    01469          * python object is NULL, the internal python array will be
    01470          * NULL and hasData() will return false.
    01471          *
    01472          * Otherwise, the function attempts to create a
    01473          * new reference to the given Python object, unless
    01474          * copying is forced by setting \a createCopy to true.
    01475          * If either of this fails, the function throws an exception.
    01476          * This will not happen if isStrictlyCompatible(obj) (in case
    01477          * of creating a new reference) or isCopyCompatible(obj)
    01478          * (in case of copying) have returned true beforehand.
    01479          */
    01480     explicit NumpyArray(PyObject *obj = 0, bool createCopy = false)
    01481     {
    01482         if(obj == 0)
    01483             return;
    01484         if(createCopy)
    01485             makeCopy(obj);
    01486         else
    01487             vigra_precondition(makeReference(obj),
    01488                   "NumpyArray(obj): Cannot construct from incompatible array.");
    01489     }
    01490
    01491        /**
    01492          * Copy constructor; does not copy the memory, but creates a
    01493          * new reference to the same underlying python object, unless
    01494          * a copy is forced by setting \a createCopy to true.
    01495          * (If the source object has no data, this one will have
    01496          * no data, too.)
    01497          */
    01498     NumpyArray(const NumpyArray &other, bool createCopy = false) :
    01499             MultiArrayView<N, typename NumpyArrayTraits<N, T, Stride>::value_type, Stride>(other),
    01500             NumpyAnyArray(other, createCopy)
    01501     {
    01502         if(!other.hasData())
    01503             return;
    01504         if(createCopy)
    01505             makeCopy(other.pyObject());
    01506         else
    01507             makeReferenceUnchecked(other.pyObject());
    01508     }
    01509
    01510        /**
    01511          * Allocate new memory and copy data from a MultiArrayView.
    01512          */
    01513     explicit NumpyArray(const view_type &other)
    01514     {
    01515         if(!other.hasData())
    01516             return;
    01517         vigra_postcondition(makeReference(init(other.shape(), false)),
    01518                   "NumpyArray(view_type): Python constructor did not produce a compatible array.");
    01519         static_cast<view_type &>(*this) = other;
    01520     }
    01521
    01522         /**
    01523          * Construct a new array object, allocating an internal python
    01524          * ndarray of the given shape (in fortran order), initialized
    01525          * with zeros.
    01526          *
    01527          * An exception is thrown when construction fails.
    01528          */
    01529     explicit NumpyArray(difference_type const & shape)
    01530     {
    01531         vigra_postcondition(makeReference(init(shape)),
    01532                      "NumpyArray(shape): Python constructor did not produce a compatible array.");
    01533     }
    01534
    01535         /**
    01536          * Construct a new array object, allocating an internal python
    01537          * ndarray of the given shape and given stride ordering, initialized
    01538          * with zeros.
    01539          *
    01540          * An exception is thrown when construction fails.
    01541          */
    01542     NumpyArray(difference_type const & shape, difference_type const & strideOrdering)
    01543     {
    01544         vigra_postcondition(makeReference(init(shape, strideOrdering)),
    01545                      "NumpyArray(shape): Python constructor did not produce a compatible array.");
    01546     }
    01547
    01548         /**
    01549          * Constructor from NumpyAnyArray.
    01550          * Equivalent to NumpyArray(other.pyObject())
    01551          */
    01552     NumpyArray(const NumpyAnyArray &other, bool createCopy = false)
    01553     {
    01554         if(!other.hasData())
    01555             return;
    01556         if(createCopy)
    01557             makeCopy(other.pyObject());
    01558         else
    01559             vigra_precondition(makeReference(other.pyObject()), //, false),
    01560                    "NumpyArray(NumpyAnyArray): Cannot construct from incompatible or empty array.");
    01561     }
    01562
    01563         /**
    01564          * Assignment operator. If this is already a view with data
    01565          * (i.e. hasData() is true) and the shapes match, the RHS
    01566          * array contents are copied.  If this is an empty view,
    01567          * assignment is identical to makeReferenceUnchecked(other.pyObject()).
    01568          * See MultiArrayView::operator= for further information on
    01569          * semantics.
    01570          */
    01571     NumpyArray &operator=(const NumpyArray &other)
    01572     {
    01573         if(hasData())
    01574             view_type::operator=(other);
    01575         else
    01576             makeReferenceUnchecked(other.pyObject());
    01577         return *this;
    01578     }
    01579
    01580         /**
    01581          * Assignment operator. If this is already a view with data
    01582          * (i.e. hasData() is true) and the shapes match, the RHS
    01583          * array contents are copied.
    01584          * If this is an empty view, assignment is identical to
    01585          * makeReference(other.pyObject()).
    01586          * Otherwise, an exception is thrown.
    01587          */
    01588     NumpyArray &operator=(const NumpyAnyArray &other)
    01589     {
    01590         if(hasData())
    01591         {
    01592             NumpyAnyArray::operator=(other);
    01593         }
    01594         else if(isStrictlyCompatible(other.pyObject()))
    01595         {
    01596             makeReferenceUnchecked(other.pyObject());
    01597         }
    01598         else
    01599         {
    01600             vigra_precondition(false,
    01601                 "NumpyArray::operator=(): Cannot assign from incompatible array.");
    01602         }
    01603         return *this;
    01604     }
    01605
    01606         /**
    01607          * Test whether a given python object is a numpy array that can be
    01608          * converted (copied) into an array compatible to this NumpyArray type.
    01609          * This means that the array's shape conforms to the requirements of
    01610          * makeCopy().
    01611          */
    01612     static bool isCopyCompatible(PyObject *obj)
    01613     {
    01614         return ArrayTraits::isArray(obj) &&
    01615                ArrayTraits::isShapeCompatible((PyArrayObject *)obj);
    01616     }
    01617
    01618         /**
    01619          * Test whether a given python object is a numpy array with a
    01620          * compatible dtype and the correct shape and strides, so that it
    01621          * can be referenced as a view by this NumpyArray type (i.e.
    01622          * it conforms to the requirements of makeReference()).
    01623          */
    01624     static bool isReferenceCompatible(PyObject *obj)
    01625     {
    01626         return ArrayTraits::isArray(obj) &&
    01627                ArrayTraits::isPropertyCompatible((PyArrayObject *)obj);
    01628     }
    01629
    01630         /**
    01631          * Like isReferenceCompatible(obj), but also executes a customized type compatibility
    01632          * check when such a check has been registered for this class via
    01633          * registerPythonArrayType().
    01634          *
    01635          * This facilitates proper overload resolution between
    01636          * NumpyArray<3, Multiband<T> > (a multiband image) and NumpyArray<3, Singleband<T> > (a scalar volume).
    01637          */
    01638     static bool isStrictlyCompatible(PyObject *obj)
    01639     {
    01640 #if VIGRA_CONVERTER_DEBUG
    01641         std::cerr << "class " << typeid(NumpyArray).name() << " got " << obj->ob_type->tp_name << "\n";
    01642         bool isClassCompatible=ArrayTraits::isClassCompatible(obj);
    01643         bool isPropertyCompatible((PyArrayObject *)obj);
    01644         std::cerr<<"isClassCompatible: "<<isClassCompatible<<std::endl;
    01645         std::cerr<<"isPropertyCompatible: "<<isPropertyCompatible<<std::endl;
    01646 #endif
    01647         return ArrayTraits::isClassCompatible(obj) &&
    01648                ArrayTraits::isPropertyCompatible((PyArrayObject *)obj);
    01649     }
    01650
    01651         /**
    01652          * Create a vector representing the standard stride ordering of a NumpyArray.
    01653          * That is, we get a vector representing the range [0,...,N-1], which
    01654          * denotes the stride ordering for Fortran order.
    01655          */
    01656     static difference_type standardStrideOrdering()
    01657     {
    01658         difference_type strideOrdering;
    01659         for(unsigned int k=0; k<N; ++k)
    01660             strideOrdering[k] = k;
    01661         return strideOrdering;
    01662     }
    01663
    01664         /**
    01665          * Set up a view to the given object without checking compatibility.
    01666          * This function must not be used unless isReferenceCompatible(obj) returned
    01667          * true on the given object (otherwise, a crash is likely).
    01668          */
    01669     void makeReferenceUnchecked(PyObject *obj)
    01670     {
    01671         NumpyAnyArray::makeReference(obj);
    01672         setupArrayView();
    01673     }
    01674
    01675         /**
    01676          * Try to set up a view referencing the given PyObject.
    01677          * Returns false if the python object is not a compatible
    01678          * numpy array (see isReferenceCompatible() or
    01679          * isStrictlyCompatible(), according to the parameter \a
    01680          * strict).
    01681          */
    01682     bool makeReference(PyObject *obj, bool strict = true)
    01683     {
    01684         if(strict)
    01685         {
    01686             if(!isStrictlyCompatible(obj))
    01687                 return false;
    01688         }
    01689         else
    01690         {
    01691             if(!isReferenceCompatible(obj))
    01692                 return false;
    01693         }
    01694         makeReferenceUnchecked(obj);
    01695         return true;
    01696     }
    01697
    01698         /**
    01699          * Try to set up a view referencing the same data as the given
    01700          * NumpyAnyArray.  This overloaded variant simply calls
    01701          * makeReference() on array.pyObject().
    01702          */
    01703     bool makeReference(const NumpyAnyArray &array, bool strict = true)
    01704     {
    01705         return makeReference(array.pyObject(), strict);
    01706     }
    01707
    01708         /**
    01709          * Set up an unsafe reference to the given MultiArrayView.
    01710          * ATTENTION: This creates a numpy.ndarray that points to the
    01711          * same data, but does not own it, so it must be ensured by
    01712          * other means that the memory does not get freed before the
    01713          * end of the ndarray's lifetime!  (One elegant way would be
    01714          * to set the 'base' attribute of the resulting ndarray to a
    01715          * python object which directly or indirectly holds the memory
    01716          * of the given MultiArrayView.)
    01717          */
    01718     void makeReference(const view_type &multiArrayView)
    01719     {
    01720         vigra_precondition(!hasData(), "makeReference(): cannot replace existing view with given buffer");
    01721
    01722         // construct an ndarray that points to our data (taking strides into account):
    01723         python_ptr array(ArrayTraits::constructor(multiArrayView.shape(), multiArrayView.data(), multiArrayView.stride()));
    01724
    01725         view_type::operator=(multiArrayView);
    01726         pyArray_ = array;
    01727     }
    01728
    01729         /**
    01730          Try to create a copy of the given PyObject.
    01731          Raises an exception when obj is not a compatible array
    01732          (see isCopyCompatible() or isStrictlyCompatible(), according to the
    01733          parameter \a strict) or the Python constructor call failed.
    01734          */
    01735     void makeCopy(PyObject *obj, bool strict = false)
    01736     {
    01737         vigra_precondition(strict ? isStrictlyCompatible(obj) : isCopyCompatible(obj),
    01738                      "NumpyArray::makeCopy(obj): Cannot copy an incompatible array.");
    01739
    01740         int M = PyArray_NDIM(obj);
    01741         TinyVector<npy_intp, N> shape;
    01742         std::copy(PyArray_DIMS(obj), PyArray_DIMS(obj)+M, shape.begin());
    01743         if(M == N-1)
    01744             shape[M] = 1;
    01745         vigra_postcondition(makeReference(init(shape, false)),
    01746                      "NumpyArray::makeCopy(obj): Copy created an incompatible array.");
    01747         NumpyAnyArray::operator=(NumpyAnyArray(obj));
    01748 //        if(PyArray_CopyInto(pyArray(), (PyArrayObject*)obj) == -1)
    01749 //            pythonToCppException(0);
    01750     }
    01751
    01752         /**
    01753             Allocate new memory with the given shape and initialize with zeros.<br>
    01754             If a stride ordering is given, the resulting array will have this stride
    01755             ordering, when it is compatible with the array's memory layout (unstrided
    01756             arrays only permit the standard ascending stride ordering).
    01757
    01758             <em>Note:</em> this operation invalidates dependent objects
    01759             (MultiArrayViews and iterators)
    01760          */
    01761     void reshape(difference_type const & shape, difference_type const & strideOrdering = standardStrideOrdering())
    01762     {
    01763         vigra_postcondition(makeReference(init(shape, strideOrdering)),
    01764                      "NumpyArray(shape): Python constructor did not produce a compatible array.");
    01765     }
    01766
    01767         /**
    01768             When this array has no data, allocate new memory with the given \a shape and
    01769             initialize with zeros. Otherwise, check if the new shape matches the old shape
    01770             and throw a precondition exception with the given \a message if not.
    01771          */
    01772     void reshapeIfEmpty(difference_type const & shape, std::string message = "")
    01773     {
    01774         reshapeIfEmpty(shape, standardStrideOrdering(), message);
    01775     }
    01776
    01777         /**
    01778             When this array has no data, allocate new memory with the given \a shape and
    01779             initialize with zeros. Otherwise, check if the new shape matches the old shape
    01780             and throw a precondition exception with the given \a message if not. If strict
    01781             is true, the given stride ordering must also match that of the existing data.
    01782          */
    01783     void reshapeIfEmpty(difference_type const & shape, difference_type const & strideOrdering,
    01784                         std::string message = "", bool strict = false)
    01785     {
    01786         if(hasData())
    01787         {
    01788             if(strict)
    01789             {
    01790                 if(message == "")
    01791                     message = "NumpyArray::reshapeIfEmpty(shape): array was not empty, and shape or stride ordering did not match.";
    01792                 vigra_precondition(shape == this->shape() && strideOrdering == this->strideOrdering(), message.c_str());
    01793             }
    01794             else
    01795             {
    01796                 if(message == "")
    01797                     message = "NumpyArray::reshapeIfEmpty(shape): array was not empty, and shape did not match.";
    01798                 vigra_precondition(shape == this->shape(), message.c_str());
    01799             }
    01800         }
    01801         else
    01802         {
    01803             reshape(shape, strideOrdering);
    01804         }
    01805     }
    01806 };
    01807
    01808     // this function assumes that pyArray_ has already been set, and compatibility been checked
    01809 template <unsigned int N, class T, class Stride>
    01810 void NumpyArray<N, T, Stride>::setupArrayView()
    01811 {
    01812     if(NumpyAnyArray::hasData())
    01813     {
    01814         unsigned int dimension = std::min<unsigned int>(actual_dimension, pyArray()->nd);
    01815         std::copy(pyArray()->dimensions, pyArray()->dimensions + dimension, this->m_shape.begin());
    01816         std::copy(pyArray()->strides, pyArray()->strides + dimension, this->m_stride.begin());
    01817         if(pyArray()->nd < actual_dimension)
    01818         {
    01819             this->m_shape[dimension] = 1;
    01820             this->m_stride[dimension] = sizeof(value_type);
    01821         }
    01822         this->m_stride /= sizeof(value_type);
    01823         this->m_ptr = reinterpret_cast<pointer>(pyArray()->data);
    01824     }
    01825     else
    01826     {
    01827         this->m_ptr = 0;
    01828     }
    01829 }
    01830
    01831
    01832 typedef NumpyArray<2, float >  NumpyFArray2;
    01833 typedef NumpyArray<3, float >  NumpyFArray3;
    01834 typedef NumpyArray<4, float >  NumpyFArray4;
    01835 typedef NumpyArray<2, Singleband<float> >  NumpyFImage;
    01836 typedef NumpyArray<3, Singleband<float> >  NumpyFVolume;
    01837 typedef NumpyArray<2, RGBValue<float> >  NumpyFRGBImage;
    01838 typedef NumpyArray<3, RGBValue<float> >  NumpyFRGBVolume;
    01839 typedef NumpyArray<3, Multiband<float> >  NumpyFMultibandImage;
    01840 typedef NumpyArray<4, Multiband<float> >  NumpyFMultibandVolume;
    01841
    01842 inline void import_vigranumpy()
    01843 {
    01844     if(_import_array() < 0)
    01845         pythonToCppException(0);
    01846     python_ptr module(PyImport_ImportModule("vigra.vigranumpycore"), python_ptr::keep_count);
    01847     pythonToCppException(module);
    01848 }
    01849
    01850 /********************************************************/
    01851 /*                                                      */
    01852 /*   NumpyArray Multiband Argument Object Factories     */
    01853 /*                                                      */
    01854 /********************************************************/
    01855
    01856 template <class PixelType, class Stride>
    01857 inline triple<ConstStridedImageIterator<PixelType>,
    01858               ConstStridedImageIterator<PixelType>,
    01859               MultibandVectorAccessor<PixelType> >
    01860 srcImageRange(NumpyArray<3, Multiband<PixelType>, Stride> const & img)
    01861 {
    01862     ConstStridedImageIterator<PixelType>
    01863         ul(img.data(), 1, img.stride(0), img.stride(1));
    01864     return triple<ConstStridedImageIterator<PixelType>,
    01865                   ConstStridedImageIterator<PixelType>,
    01866                   MultibandVectorAccessor<PixelType> >
    01867         (ul, ul + Size2D(img.shape(0), img.shape(1)), MultibandVectorAccessor<PixelType>(img.shape(2), img.stride(2)));
    01868 }
    01869
    01870 template <class PixelType, class Stride>
    01871 inline pair< ConstStridedImageIterator<PixelType>,
    01872              MultibandVectorAccessor<PixelType> >
    01873 srcImage(NumpyArray<3, Multiband<PixelType>, Stride> const & img)
    01874 {
    01875     ConstStridedImageIterator<PixelType>
    01876         ul(img.data(), 1, img.stride(0), img.stride(1));
    01877     return pair<ConstStridedImageIterator<PixelType>, MultibandVectorAccessor<PixelType> >
    01878         (ul, MultibandVectorAccessor<PixelType>(img.shape(2), img.stride(2)));
    01879 }
    01880
    01881 template <class PixelType, class Stride>
    01882 inline triple< StridedImageIterator<PixelType>,
    01883                StridedImageIterator<PixelType>,
    01884                MultibandVectorAccessor<PixelType> >
    01885 destImageRange(NumpyArray<3, Multiband<PixelType>, Stride> & img)
    01886 {
    01887     StridedImageIterator<PixelType>
    01888         ul(img.data(), 1, img.stride(0), img.stride(1));
    01889     typedef typename AccessorTraits<PixelType>::default_accessor Accessor;
    01890     return triple<StridedImageIterator<PixelType>,
    01891                   StridedImageIterator<PixelType>,
    01892                   MultibandVectorAccessor<PixelType> >
    01893         (ul, ul + Size2D(img.shape(0), img.shape(1)),
    01894         MultibandVectorAccessor<PixelType>(img.shape(2), img.stride(2)));
    01895 }
    01896
    01897 template <class PixelType, class Stride>
    01898 inline pair< StridedImageIterator<PixelType>,
    01899              MultibandVectorAccessor<PixelType> >
    01900 destImage(NumpyArray<3, Multiband<PixelType>, Stride> & img)
    01901 {
    01902     StridedImageIterator<PixelType>
    01903         ul(img.data(), 1, img.stride(0), img.stride(1));
    01904     return pair<StridedImageIterator<PixelType>, MultibandVectorAccessor<PixelType> >
    01905         (ul, MultibandVectorAccessor<PixelType>(img.shape(2), img.stride(2)));
    01906 }
    01907
    01908 template <class PixelType, class Stride>
    01909 inline pair< ConstStridedImageIterator<PixelType>,
    01910              MultibandVectorAccessor<PixelType> >
    01911 maskImage(NumpyArray<3, Multiband<PixelType>, Stride> const & img)
    01912 {
    01913     ConstStridedImageIterator<PixelType>
    01914         ul(img.data(), 1, img.stride(0), img.stride(1));
    01915     typedef typename AccessorTraits<PixelType>::default_accessor Accessor;
    01916     return pair<ConstStridedImageIterator<PixelType>, MultibandVectorAccessor<PixelType> >
    01917         (ul, MultibandVectorAccessor<PixelType>(img.shape(2), img.stride(2)));
    01918 }
    01919
    01920 } // namespace vigra
    01921
    01922 #endif // VIGRA_NUMPY_ARRAY_HXX

© Ullrich Köthe (ullrich.koethe@iwr.uni-heidelberg.de) Heidelberg Collaboratory for Image Processing, University of Heidelberg, Germany | html generated using doxygen and Python vigra 1.7.1 (3 Dec 2010)