[ 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 res.get() == 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 inline
00463 python_ptr constructNumpyArrayFromData(
00464     PyTypeObject * type,
00465     ArrayVector<npy_intp> const & shape, npy_intp *strides,
00466     NPY_TYPES typeCode, void *data)
00467 {
00468     python_ptr array(PyArray_New(type, shape.size(), (npy_intp *)shape.begin(), typeCode, strides, data, 0, NPY_WRITEABLE, 0),
00469                      python_ptr::keep_count);
00470     pythonToCppException(array);
00471 
00472     return array;
00473 }
00474 
00475 
00476 } // namespace detail
00477 
00478 /********************************************************/
00479 /*                                                      */
00480 /*               NumpyArrayValuetypeTraits              */
00481 /*                                                      */
00482 /********************************************************/
00483 
00484 template<class ValueType>
00485 struct ERROR_NumpyArrayValuetypeTraits_not_specialized_for_ { };
00486 
00487 template<class ValueType>
00488 struct NumpyArrayValuetypeTraits
00489 {
00490     static bool isValuetypeCompatible(PyArrayObject const * obj)
00491     {
00492         return ERROR_NumpyArrayValuetypeTraits_not_specialized_for_<ValueType>();
00493     }
00494 
00495     static ERROR_NumpyArrayValuetypeTraits_not_specialized_for_<ValueType> typeCode;
00496 
00497     static std::string typeName()
00498     {
00499         return std::string("ERROR: NumpyArrayValuetypeTraits not specialized for this case");
00500     }
00501 
00502     static std::string typeNameImpex()
00503     {
00504         return std::string("ERROR: NumpyArrayValuetypeTraits not specialized for this case");
00505     }
00506 
00507     static PyObject * typeObject()
00508     {
00509         return (PyObject *)0;
00510     }
00511 };
00512 
00513 template<class ValueType>
00514 ERROR_NumpyArrayValuetypeTraits_not_specialized_for_<ValueType> NumpyArrayValuetypeTraits<ValueType>::typeCode;
00515 
00516 #define VIGRA_NUMPY_VALUETYPE_TRAITS(type, typeID, numpyTypeName, impexTypeName) \
00517 template <> \
00518 struct NumpyArrayValuetypeTraits<type > \
00519 { \
00520     static bool isValuetypeCompatible(PyArrayObject const * obj) /* obj must not be NULL */ \
00521     { \
00522         return PyArray_EquivTypenums(typeID, PyArray_DESCR((PyObject *)obj)->type_num) && \
00523                PyArray_ITEMSIZE((PyObject *)obj) == sizeof(type); \
00524     } \
00525     \
00526     static NPY_TYPES const typeCode = typeID; \
00527     \
00528     static std::string typeName() \
00529     { \
00530         return #numpyTypeName; \
00531     } \
00532     \
00533     static std::string typeNameImpex() \
00534     { \
00535         return impexTypeName; \
00536     } \
00537     \
00538     static PyObject * typeObject() \
00539     { \
00540         return PyArray_TypeObjectFromType(typeID); \
00541     } \
00542 };
00543 
00544 VIGRA_NUMPY_VALUETYPE_TRAITS(bool,           NPY_BOOL, bool, "UINT8")
00545 VIGRA_NUMPY_VALUETYPE_TRAITS(signed char,    NPY_INT8, int8, "INT16")
00546 VIGRA_NUMPY_VALUETYPE_TRAITS(unsigned char,  NPY_UINT8, uint8, "UINT8")
00547 VIGRA_NUMPY_VALUETYPE_TRAITS(short,          NPY_INT16, int16, "INT16")
00548 VIGRA_NUMPY_VALUETYPE_TRAITS(unsigned short, NPY_UINT16, uint16, "UINT16")
00549 
00550 #if VIGRA_BITSOF_LONG == 32
00551 VIGRA_NUMPY_VALUETYPE_TRAITS(long,           NPY_INT32, int32, "INT32")
00552 VIGRA_NUMPY_VALUETYPE_TRAITS(unsigned long,  NPY_UINT32, uint32, "UINT32")
00553 #elif VIGRA_BITSOF_LONG == 64
00554 VIGRA_NUMPY_VALUETYPE_TRAITS(long,           NPY_INT64, int64, "DOUBLE")
00555 VIGRA_NUMPY_VALUETYPE_TRAITS(unsigned long,  NPY_UINT64, uint64, "DOUBLE")
00556 #endif
00557 
00558 #if VIGRA_BITSOF_INT == 32
00559 VIGRA_NUMPY_VALUETYPE_TRAITS(int,            NPY_INT32, int32, "INT32")
00560 VIGRA_NUMPY_VALUETYPE_TRAITS(unsigned int,   NPY_UINT32, uint32, "UINT32")
00561 #elif VIGRA_BITSOF_INT == 64
00562 VIGRA_NUMPY_VALUETYPE_TRAITS(int,            NPY_INT64, int64, "DOUBLE")
00563 VIGRA_NUMPY_VALUETYPE_TRAITS(unsigned int,   NPY_UINT64, uint64, "DOUBLE")
00564 #endif
00565 
00566 #ifdef PY_LONG_LONG
00567 # if VIGRA_BITSOF_LONG_LONG == 32
00568 VIGRA_NUMPY_VALUETYPE_TRAITS(long long,            NPY_INT32, int32, "INT32")
00569 VIGRA_NUMPY_VALUETYPE_TRAITS(unsigned long long,   NPY_UINT32, uint32, "UINT32")
00570 # elif VIGRA_BITSOF_LONG_LONG == 64
00571 VIGRA_NUMPY_VALUETYPE_TRAITS(long long,          NPY_INT64, int64, "DOUBLE")
00572 VIGRA_NUMPY_VALUETYPE_TRAITS(unsigned long long, NPY_UINT64, uint64, "DOUBLE")
00573 # endif
00574 #endif
00575 
00576 VIGRA_NUMPY_VALUETYPE_TRAITS(npy_float32, NPY_FLOAT32, float32, "FLOAT")
00577 VIGRA_NUMPY_VALUETYPE_TRAITS(npy_float64, NPY_FLOAT64, float64, "DOUBLE")
00578 #if NPY_SIZEOF_LONGDOUBLE != NPY_SIZEOF_DOUBLE
00579 VIGRA_NUMPY_VALUETYPE_TRAITS(npy_longdouble, NPY_LONGDOUBLE, longdouble, "")
00580 #endif
00581 VIGRA_NUMPY_VALUETYPE_TRAITS(npy_cfloat, NPY_CFLOAT, complex64, "")
00582 VIGRA_NUMPY_VALUETYPE_TRAITS(std::complex<npy_float>, NPY_CFLOAT, complex64, "")
00583 VIGRA_NUMPY_VALUETYPE_TRAITS(npy_cdouble, NPY_CDOUBLE, complex128, "")
00584 VIGRA_NUMPY_VALUETYPE_TRAITS(std::complex<npy_double>, NPY_CDOUBLE, complex128, "")
00585 VIGRA_NUMPY_VALUETYPE_TRAITS(npy_clongdouble, NPY_CLONGDOUBLE, clongdouble, "")
00586 #if NPY_SIZEOF_LONGDOUBLE != NPY_SIZEOF_DOUBLE
00587 VIGRA_NUMPY_VALUETYPE_TRAITS(std::complex<npy_longdouble>, NPY_CLONGDOUBLE, clongdouble, "")
00588 #endif
00589 
00590 #undef VIGRA_NUMPY_VALUETYPE_TRAITS
00591 
00592 /********************************************************/
00593 /*                                                      */
00594 /*                  NumpyArrayTraits                    */
00595 /*                                                      */
00596 /********************************************************/
00597 
00598 template <class U, int N>
00599 bool stridesAreAscending(TinyVector<U, N> const & strides)
00600 {
00601     for(int k=1; k<N; ++k)
00602         if(strides[k] < strides[k-1])
00603             return false;
00604     return true;
00605 }
00606 
00607 template<unsigned int N, class T, class Stride>
00608 struct NumpyArrayTraits;
00609 
00610 template<unsigned int N, class T>
00611 struct NumpyArrayTraits<N, T, StridedArrayTag>
00612 {
00613     typedef T dtype;
00614     typedef T value_type;
00615     typedef NumpyArrayValuetypeTraits<T> ValuetypeTraits;
00616     static NPY_TYPES const typeCode = ValuetypeTraits::typeCode;
00617     
00618     enum { spatialDimensions = N, channels = 1 };
00619 
00620     static bool isArray(PyObject * obj)
00621     {
00622         return obj && PyArray_Check(obj);
00623     }
00624 
00625     static bool isClassCompatible(PyObject * obj)
00626     {
00627         return detail::performCustomizedArrayTypecheck(obj, typeKeyFull(), typeKey());
00628     }
00629 
00630     static bool isValuetypeCompatible(PyArrayObject * obj)  /* obj must not be NULL */
00631     {
00632         return ValuetypeTraits::isValuetypeCompatible(obj);
00633     }
00634 
00635     static bool isShapeCompatible(PyArrayObject * obj) /* obj must not be NULL */
00636     {
00637         return PyArray_NDIM((PyObject *)obj) == N-1 ||
00638                PyArray_NDIM((PyObject *)obj) == N ||
00639                (PyArray_NDIM((PyObject *)obj) == N+1 && PyArray_DIM((PyObject *)obj, N) == 1);
00640     }
00641 
00642     static bool isPropertyCompatible(PyArrayObject * obj) /* obj must not be NULL */
00643     {
00644         return ValuetypeTraits::isValuetypeCompatible(obj) &&
00645                isShapeCompatible(obj);
00646     }
00647 
00648     template <class U>
00649     static python_ptr constructor(TinyVector<U, N> const & shape,
00650                                   T *data, TinyVector<U, N> const & stride)
00651     {
00652         TinyVector<npy_intp, N> npyStride(stride * sizeof(T));
00653         return detail::constructNumpyArrayFromData(typeKeyFull(), typeKey(), shape, ValuetypeTraits::typeCode, data, npyStride.begin());
00654     }
00655 
00656     static std::string typeKey()
00657     {
00658         static std::string key = std::string("NumpyArray<") + asString(N) + ", *>";
00659         return key;
00660     }
00661 
00662     static std::string typeKeyFull()
00663     {
00664         static std::string key = std::string("NumpyArray<") + asString(N) + ", " +
00665                                  ValuetypeTraits::typeName() + ", StridedArrayTag>";
00666         return key;
00667     }
00668 };
00669 
00670 /********************************************************/
00671 
00672 template<unsigned int N, class T>
00673 struct NumpyArrayTraits<N, T, UnstridedArrayTag>
00674 : public NumpyArrayTraits<N, T, StridedArrayTag>
00675 {
00676     typedef NumpyArrayTraits<N, T, StridedArrayTag> BaseType;
00677     typedef typename BaseType::ValuetypeTraits ValuetypeTraits;
00678 
00679     static bool isShapeCompatible(PyArrayObject * obj) /* obj must not be NULL */
00680     {
00681         return BaseType::isShapeCompatible(obj) &&
00682                PyArray_STRIDES((PyObject *)obj)[0] == PyArray_ITEMSIZE((PyObject *)obj);
00683     }
00684 
00685     static bool isPropertyCompatible(PyArrayObject * obj) /* obj must not be NULL */
00686     {
00687         return BaseType::isValuetypeCompatible(obj) &&
00688                isShapeCompatible(obj);
00689     }
00690 
00691     template <class U>
00692     static python_ptr constructor(TinyVector<U, N> const & shape,
00693                                   T *data, TinyVector<U, N> const & stride)
00694     {
00695         TinyVector<npy_intp, N> npyStride(stride * sizeof(T));
00696         return detail::constructNumpyArrayFromData(typeKeyFull(), BaseType::typeKey(), shape, ValuetypeTraits::typeCode, data, npyStride.begin());
00697     }
00698 
00699     static std::string typeKeyFull()
00700     {
00701         static std::string key = std::string("NumpyArray<") + asString(N) + ", " +
00702                                  ValuetypeTraits::typeName() + ", UnstridedArrayTag>";
00703         return key;
00704     }
00705 };
00706 
00707 /********************************************************/
00708 
00709 template<unsigned int N, class T>
00710 struct NumpyArrayTraits<N, Singleband<T>, StridedArrayTag>
00711 : public NumpyArrayTraits<N, T, StridedArrayTag>
00712 {
00713     typedef NumpyArrayTraits<N, T, StridedArrayTag> BaseType;
00714     typedef typename BaseType::ValuetypeTraits ValuetypeTraits;
00715 
00716     static bool isClassCompatible(PyObject * obj)
00717     {
00718         return detail::performCustomizedArrayTypecheck(obj, typeKeyFull(), typeKey());
00719     }
00720 
00721     template <class U>
00722     static python_ptr constructor(TinyVector<U, N> const & shape,
00723                                   T *data, TinyVector<U, N> const & stride)
00724     {
00725         TinyVector<npy_intp, N> npyStride(stride * sizeof(T));
00726         return detail::constructNumpyArrayFromData(typeKeyFull(), typeKey(), shape, ValuetypeTraits::typeCode, data, npyStride.begin());
00727     }
00728 
00729     static std::string typeKey()
00730     {
00731         static std::string key = std::string("NumpyArray<") + asString(N) + ", Singleband<*> >";
00732         return key;
00733     }
00734 
00735     static std::string typeKeyFull()
00736     {
00737         static std::string key = std::string("NumpyArray<") + asString(N) + ", Singleband<" +
00738                                  ValuetypeTraits::typeName() + ">, StridedArrayTag>";
00739         return key;
00740     }
00741 };
00742 
00743 /********************************************************/
00744 
00745 template<unsigned int N, class T>
00746 struct NumpyArrayTraits<N, Singleband<T>, UnstridedArrayTag>
00747 : public NumpyArrayTraits<N, Singleband<T>, StridedArrayTag>
00748 {
00749     typedef NumpyArrayTraits<N, T, UnstridedArrayTag> UnstridedTraits;
00750     typedef NumpyArrayTraits<N, Singleband<T>, StridedArrayTag> BaseType;
00751     typedef typename BaseType::ValuetypeTraits ValuetypeTraits;
00752 
00753     static bool isShapeCompatible(PyArrayObject * obj) /* obj must not be NULL */
00754     {
00755         return UnstridedTraits::isShapeCompatible(obj);
00756     }
00757 
00758     static bool isPropertyCompatible(PyArrayObject * obj) /* obj must not be NULL */
00759     {
00760         return UnstridedTraits::isPropertyCompatible(obj);
00761     }
00762 
00763     template <class U>
00764     static python_ptr constructor(TinyVector<U, N> const & shape,
00765                                   T *data, TinyVector<U, N> const & stride)
00766     {
00767         TinyVector<npy_intp, N> npyStride(stride * sizeof(T));
00768         return detail::constructNumpyArrayFromData(typeKeyFull(), BaseType::typeKey(), shape, ValuetypeTraits::typeCode, data, npyStride.begin());
00769     }
00770 
00771     static std::string typeKeyFull()
00772     {
00773         static std::string key = std::string("NumpyArray<") + asString(N) + ", Singleband<" +
00774                                  ValuetypeTraits::typeName() + ">, UnstridedArrayTag>";
00775         return key;
00776     }
00777 };
00778 
00779 /********************************************************/
00780 
00781 template<unsigned int N, class T>
00782 struct NumpyArrayTraits<N, Multiband<T>, StridedArrayTag>
00783 : public NumpyArrayTraits<N, T, StridedArrayTag>
00784 {
00785     typedef NumpyArrayTraits<N, T, StridedArrayTag> BaseType;
00786     typedef typename BaseType::ValuetypeTraits ValuetypeTraits;
00787 
00788     enum { spatialDimensions = N-1, channels = 0 };
00789 
00790     static bool isClassCompatible(PyObject * obj)
00791     {
00792         return detail::performCustomizedArrayTypecheck(obj, typeKeyFull(), typeKey());
00793     }
00794 
00795     static bool isShapeCompatible(PyArrayObject * obj) /* obj must not be NULL */
00796     {
00797         return PyArray_NDIM(obj) == N || PyArray_NDIM(obj) == N-1;
00798     }
00799 
00800     static bool isPropertyCompatible(PyArrayObject * obj) /* obj must not be NULL */
00801     {
00802         return ValuetypeTraits::isValuetypeCompatible(obj) &&
00803                isShapeCompatible(obj);
00804     }
00805 
00806     template <class U>
00807     static python_ptr constructor(TinyVector<U, N> const & shape,
00808                                   T *data, TinyVector<U, N> const & stride)
00809     {
00810         TinyVector<npy_intp, N> npyStride(stride * sizeof(T));
00811         return detail::constructNumpyArrayFromData(typeKeyFull(), typeKey(), shape, ValuetypeTraits::typeCode, data, npyStride.begin());
00812     }
00813 
00814     static std::string typeKey()
00815     {
00816         static std::string key = std::string("NumpyArray<") + asString(N) + ", Multiband<*> >";
00817         return key;
00818     }
00819 
00820     static std::string typeKeyFull()
00821     {
00822         static std::string key = std::string("NumpyArray<") + asString(N) + ", Multiband<" +
00823                                  ValuetypeTraits::typeName() + ">, StridedArrayTag>";
00824         return key;
00825     }
00826 };
00827 
00828 /********************************************************/
00829 
00830 template<unsigned int N, class T>
00831 struct NumpyArrayTraits<N, Multiband<T>, UnstridedArrayTag>
00832 : public NumpyArrayTraits<N, Multiband<T>, StridedArrayTag>
00833 {
00834     typedef NumpyArrayTraits<N, Multiband<T>, StridedArrayTag> BaseType;
00835     typedef typename BaseType::ValuetypeTraits ValuetypeTraits;
00836 
00837     static bool isShapeCompatible(PyArrayObject * obj) /* obj must not be NULL */
00838     {
00839         return BaseType::isShapeCompatible(obj) &&
00840                PyArray_STRIDES((PyObject *)obj)[0] == PyArray_ITEMSIZE((PyObject *)obj);
00841     }
00842 
00843     static bool isPropertyCompatible(PyArrayObject * obj) /* obj must not be NULL */
00844     {
00845         return BaseType::isValuetypeCompatible(obj) &&
00846                isShapeCompatible(obj);
00847     }
00848 
00849     template <class U>
00850     static python_ptr constructor(TinyVector<U, N> const & shape,
00851                                   T *data, TinyVector<U, N> const & stride)
00852     {
00853         TinyVector<npy_intp, N> npyStride(stride * sizeof(T));
00854         return detail::constructNumpyArrayFromData(typeKeyFull(), BaseType::typeKey(), shape, ValuetypeTraits::typeCode, data, npyStride.begin());
00855     }
00856 
00857     static std::string typeKeyFull()
00858     {
00859         static std::string key = std::string("NumpyArray<") + asString(N) + ", Multiband<" +
00860                                  ValuetypeTraits::typeName() + ">, UnstridedArrayTag>";
00861         return key;
00862     }
00863 };
00864 
00865 /********************************************************/
00866 
00867 template<unsigned int N, int M, class T>
00868 struct NumpyArrayTraits<N, TinyVector<T, M>, StridedArrayTag>
00869 {
00870     typedef T dtype;
00871     typedef TinyVector<T, M> value_type;
00872     typedef NumpyArrayValuetypeTraits<T> ValuetypeTraits;
00873     static NPY_TYPES const typeCode = ValuetypeTraits::typeCode;
00874 
00875     enum { spatialDimensions = N, channels = M };
00876 
00877     static bool isArray(PyObject * obj)
00878     {
00879         return obj && PyArray_Check(obj);
00880     }
00881 
00882     static bool isClassCompatible(PyObject * obj)
00883     {
00884         return detail::performCustomizedArrayTypecheck(obj, typeKeyFull(), typeKey());
00885     }
00886 
00887     static bool isValuetypeCompatible(PyArrayObject * obj)  /* obj must not be NULL */
00888     {
00889         return ValuetypeTraits::isValuetypeCompatible(obj);
00890     }
00891 
00892     static bool isShapeCompatible(PyArrayObject * obj) /* obj must not be NULL */
00893     {
00894         return PyArray_NDIM((PyObject *)obj) == N+1 &&
00895                PyArray_DIM((PyObject *)obj, N) == M &&
00896                PyArray_STRIDES((PyObject *)obj)[N] == PyArray_ITEMSIZE((PyObject *)obj);
00897     }
00898 
00899     static bool isPropertyCompatible(PyArrayObject * obj) /* obj must not be NULL */
00900     {
00901         return ValuetypeTraits::isValuetypeCompatible(obj) &&
00902                isShapeCompatible(obj);
00903     }
00904 
00905     template <class U>
00906     static python_ptr constructor(TinyVector<U, N> const & shape,
00907                                   T *data, TinyVector<U, N> const & stride)
00908     {
00909         TinyVector<npy_intp, N+1> npyShape;
00910         std::copy(shape.begin(), shape.end(), npyShape.begin());
00911         npyShape[N] = M;
00912 
00913         TinyVector<npy_intp, N+1> npyStride;
00914         std::transform(
00915             stride.begin(), stride.end(), npyStride.begin(),
00916             std::bind2nd(std::multiplies<npy_intp>(), sizeof(value_type)));
00917         npyStride[N] = sizeof(T);
00918 
00919         return detail::constructNumpyArrayFromData(
00920             typeKeyFull(), typeKey(), npyShape,
00921             ValuetypeTraits::typeCode, data, npyStride.begin());
00922     }
00923 
00924     static std::string typeKey()
00925     {
00926         static std::string key = std::string("NumpyArray<") + asString(N) + ", TinyVector<*, " + asString(M) + "> >";
00927         return key;
00928     }
00929 
00930     static std::string typeKeyFull()
00931     {
00932         static std::string key = std::string("NumpyArray<") + asString(N) +
00933                       ", TinyVector<" + ValuetypeTraits::typeName() + ", " + asString(M) + ">, StridedArrayTag>";
00934         return key;
00935     }
00936 };
00937 
00938 /********************************************************/
00939 
00940 template<unsigned int N, int M, class T>
00941 struct NumpyArrayTraits<N, TinyVector<T, M>, UnstridedArrayTag>
00942 : public NumpyArrayTraits<N, TinyVector<T, M>, StridedArrayTag>
00943 {
00944     typedef NumpyArrayTraits<N, TinyVector<T, M>, StridedArrayTag> BaseType;
00945     typedef typename BaseType::value_type value_type;
00946     typedef typename BaseType::ValuetypeTraits ValuetypeTraits;
00947 
00948     static bool isShapeCompatible(PyArrayObject * obj) /* obj must not be NULL */
00949     {
00950         return BaseType::isShapeCompatible(obj) &&
00951                PyArray_STRIDES((PyObject *)obj)[0] == sizeof(TinyVector<T, M>);
00952     }
00953 
00954     static bool isPropertyCompatible(PyArrayObject * obj) /* obj must not be NULL */
00955     {
00956         return BaseType::isValuetypeCompatible(obj) &&
00957                isShapeCompatible(obj);
00958     }
00959 
00960     template <class U>
00961     static python_ptr constructor(TinyVector<U, N> const & shape,
00962                                   T *data, TinyVector<U, N> const & stride)
00963     {
00964         TinyVector<npy_intp, N+1> npyShape;
00965         std::copy(shape.begin(), shape.end(), npyShape.begin());
00966         npyShape[N] = M;
00967 
00968         TinyVector<npy_intp, N+1> npyStride;
00969         std::transform(
00970             stride.begin(), stride.end(), npyStride.begin(),
00971             std::bind2nd(std::multiplies<npy_intp>(), sizeof(value_type)));
00972         npyStride[N] = sizeof(T);
00973 
00974         return detail::constructNumpyArrayFromData(
00975             typeKeyFull(), BaseType::typeKey(), npyShape,
00976             ValuetypeTraits::typeCode, data, npyStride.begin());
00977     }
00978 
00979     static std::string typeKeyFull()
00980     {
00981         static std::string key = std::string("NumpyArray<") + asString(N) +
00982                       ", TinyVector<" + ValuetypeTraits::typeName() + ", " + asString(M) + ">, UnstridedArrayTag>";
00983         return key;
00984     }
00985 };
00986 
00987 /********************************************************/
00988 
00989 template<unsigned int N, class T>
00990 struct NumpyArrayTraits<N, RGBValue<T>, StridedArrayTag>
00991 : public NumpyArrayTraits<N, TinyVector<T, 3>, StridedArrayTag>
00992 {
00993     typedef T dtype;
00994     typedef RGBValue<T> value_type;
00995     typedef NumpyArrayValuetypeTraits<T> ValuetypeTraits;
00996 
00997     static bool isClassCompatible(PyObject * obj)
00998     {
00999         return detail::performCustomizedArrayTypecheck(obj, typeKeyFull(), typeKey());
01000     }
01001 
01002     template <class U>
01003     static python_ptr constructor(TinyVector<U, N> const & shape,
01004                                   T *data, TinyVector<U, N> const & stride)
01005     {
01006         TinyVector<npy_intp, N+1> npyShape;
01007         std::copy(shape.begin(), shape.end(), npyShape.begin());
01008         npyShape[N] = 3;
01009 
01010         TinyVector<npy_intp, N+1> npyStride;
01011         std::transform(
01012             stride.begin(), stride.end(), npyStride.begin(),
01013             std::bind2nd(std::multiplies<npy_intp>(), sizeof(value_type)));
01014         npyStride[N] = sizeof(T);
01015 
01016         return detail::constructNumpyArrayFromData(
01017             typeKeyFull(), typeKey(), npyShape,
01018             ValuetypeTraits::typeCode, data, npyStride.begin());
01019     }
01020 
01021     static std::string typeKey()
01022     {
01023         static std::string key = std::string("NumpyArray<") + asString(N) + ", RGBValue<*> >";
01024         return key;
01025     }
01026 
01027     static std::string typeKeyFull()
01028     {
01029         static std::string key = std::string("NumpyArray<") + asString(N) +
01030                       ", RGBValue<" + ValuetypeTraits::typeName() + ">, StridedArrayTag>";
01031         return key;
01032     }
01033 };
01034 
01035 /********************************************************/
01036 
01037 template<unsigned int N, class T>
01038 struct NumpyArrayTraits<N, RGBValue<T>, UnstridedArrayTag>
01039 : public NumpyArrayTraits<N, RGBValue<T>, StridedArrayTag>
01040 {
01041     typedef NumpyArrayTraits<N, TinyVector<T, 3>, UnstridedArrayTag> UnstridedTraits;
01042     typedef NumpyArrayTraits<N, RGBValue<T>, StridedArrayTag> BaseType;
01043     typedef typename BaseType::value_type value_type;
01044     typedef typename BaseType::ValuetypeTraits ValuetypeTraits;
01045 
01046     static bool isShapeCompatible(PyArrayObject * obj) /* obj must not be NULL */
01047     {
01048         return UnstridedTraits::isShapeCompatible(obj);
01049     }
01050 
01051     static bool isPropertyCompatible(PyArrayObject * obj) /* obj must not be NULL */
01052     {
01053         return UnstridedTraits::isPropertyCompatible(obj);
01054     }
01055 
01056     template <class U>
01057     static python_ptr constructor(TinyVector<U, N> const & shape,
01058                                   T *data, TinyVector<U, N> const & stride)
01059     {
01060         TinyVector<npy_intp, N+1> npyShape;
01061         std::copy(shape.begin(), shape.end(), npyShape.begin());
01062         npyShape[N] = 3;
01063 
01064         TinyVector<npy_intp, N+1> npyStride;
01065         std::transform(
01066             stride.begin(), stride.end(), npyStride.begin(),
01067             std::bind2nd(std::multiplies<npy_intp>(), sizeof(value_type)));
01068         npyStride[N] = sizeof(T);
01069 
01070         return detail::constructNumpyArrayFromData(
01071             typeKeyFull(), BaseType::typeKey(), npyShape,
01072             ValuetypeTraits::typeCode, data, npyStride.begin());
01073     }
01074 
01075     static std::string typeKeyFull()
01076     {
01077         static std::string key = std::string("NumpyArray<") + asString(N) +
01078                       ", RGBValue<" + ValuetypeTraits::typeName() + ">, UnstridedArrayTag>";
01079         return key;
01080     }
01081 };
01082 
01083 /********************************************************/
01084 /*                                                      */
01085 /*                    NumpyAnyArray                     */
01086 /*                                                      */
01087 /********************************************************/
01088 
01089 /** Wrapper class for a Python array.
01090     
01091     This class stores a reference-counted pointer to an Python numpy array object, 
01092     i.e. an object where <tt>PyArray_Check(object)</tt> returns true (in Python, the
01093     object is then a subclass of <tt>numpy.ndarray</tt>). This class is mainly used 
01094     as a smart pointer to these arrays, but some basic access and conversion functions
01095     are also provided.
01096 
01097     <b>\#include</b> <<a href="numpy__array_8hxx-source.html">vigra/numpy_array.hxx</a>><br>
01098     Namespace: vigra
01099 */
01100 class NumpyAnyArray
01101 {
01102   protected:
01103     python_ptr pyArray_;
01104     
01105     // We want to apply broadcasting to the channel dimension.
01106     // Since only leading dimensions can be added during numpy
01107     // broadcasting, we permute the array accordingly.
01108     NumpyAnyArray permuteChannelsToFront() const
01109     {
01110         MultiArrayIndex M = ndim();
01111         ArrayVector<npy_intp> permutation(M);
01112         for(int k=0; k<M; ++k)
01113             permutation[k] = M-1-k;
01114         // explicit cast to int is neede here to avoid gcc c++0x compilation
01115         // error: narrowing conversion of ‘M’ from ‘vigra::MultiArrayIndex’
01116         //        to ‘int’ inside { }
01117         // int overflow should not occur here because PyArray_NDIM returns
01118         // an integer which is converted to long in NumpyAnyArray::ndim()
01119         PyArray_Dims permute = { permutation.begin(), (int) M };
01120         python_ptr array(PyArray_Transpose(pyArray(), &permute), python_ptr::keep_count);
01121         pythonToCppException(array);
01122         return NumpyAnyArray(array.ptr());
01123     }
01124 
01125   public:
01126 
01127         /// difference type
01128     typedef ArrayVector<npy_intp> difference_type;
01129 
01130         /**
01131          Construct from a Python object. If \a obj is NULL, or is not a subclass
01132          of numpy.ndarray, the resulting NumpyAnyArray will have no data (i.e.
01133          hasData() returns false). Otherwise, it creates a new reference to the array
01134          \a obj, unless \a createCopy is true, where a new array is created by calling
01135          the C-equivalent of obj->copy().
01136          */
01137     explicit NumpyAnyArray(PyObject * obj = 0, bool createCopy = false, PyTypeObject * type = 0)
01138     {
01139         if(obj == 0)
01140             return;
01141         vigra_precondition(type == 0 || PyType_IsSubtype(type, &PyArray_Type),
01142              "NumpyAnyArray(obj, createCopy, type): type must be numpy.ndarray or a subclass thereof.");
01143         if(createCopy)
01144             makeCopy(obj, type);
01145         else
01146             vigra_precondition(makeReference(obj, type), "NumpyAnyArray(obj): obj isn't a numpy array.");
01147     }
01148 
01149         /**
01150          Copy constructor. By default, it creates a new reference to the array
01151          \a other. When \a createCopy is true, a new array is created by calling
01152          the C-equivalent of other.copy().
01153          */
01154     NumpyAnyArray(NumpyAnyArray const & other, bool createCopy = false, PyTypeObject * type = 0)
01155     {
01156         if(!other.hasData())
01157             return;
01158         vigra_precondition(type == 0 || PyType_IsSubtype(type, &PyArray_Type),
01159              "NumpyAnyArray(obj, createCopy, type): type must be numpy.ndarray or a subclass thereof.");
01160         if(createCopy)
01161             makeCopy(other.pyObject(), type);
01162         else
01163             makeReference(other.pyObject(), type);
01164     }
01165 
01166     // auto-generated destructor is ok
01167 
01168         /**
01169          * Assignment operator. If this is already a view with data
01170          * (i.e. hasData() is true) and the shapes match, the RHS
01171          * array contents are copied via the C-equivalent of 
01172          * 'self[...] = other[...]'. If the shapes don't matched,
01173          * broadcasting is tried on the trailing (i.e. channel)
01174          * dimension.
01175          * If the LHS is an empty view, assignment is identical to
01176          * makeReference(other.pyObject()).
01177          */
01178     NumpyAnyArray & operator=(NumpyAnyArray const & other)
01179     {
01180         if(hasData())
01181         {
01182             vigra_precondition(other.hasData(),
01183                 "NumpyArray::operator=(): Cannot assign from empty array.");
01184             if(PyArray_CopyInto(permuteChannelsToFront().pyArray(), other.permuteChannelsToFront().pyArray()) == -1)
01185                 pythonToCppException(0);
01186         }
01187         else
01188         {
01189             pyArray_ = other.pyArray_;
01190         }
01191         return *this;
01192     }
01193 
01194         /**
01195          Returns the number of dimensions of this array, or 0 if
01196          hasData() is false.
01197          */
01198     MultiArrayIndex ndim() const
01199     {
01200         if(hasData())
01201             return PyArray_NDIM(pyObject());
01202         return 0;
01203     }
01204 
01205         /**
01206          Returns the number of spatial dimensions of this array, or 0 if
01207          hasData() is false. If the enclosed Python array does not define
01208          the attribute spatialDimensions, ndim() is returned.
01209          */
01210     MultiArrayIndex spatialDimensions() const
01211     {
01212         if(!hasData())
01213             return 0;
01214         MultiArrayIndex s = detail::spatialDimensions(pyObject());
01215         if(s == -1)
01216             s = ndim();
01217         return s;
01218     }
01219 
01220         /**
01221          Returns the shape of this array. The size of
01222          the returned shape equals ndim().
01223          */
01224     difference_type shape() const
01225     {
01226         if(hasData())
01227             return difference_type(PyArray_DIMS(pyObject()), PyArray_DIMS(pyObject()) + ndim());
01228         return difference_type();
01229     }
01230 
01231         /** Compute the ordering of the strides of this array.
01232             The result is describes the current permutation of the axes relative
01233             to an ascending stride order.
01234         */
01235     difference_type strideOrdering() const
01236     {
01237         if(!hasData())
01238             return difference_type();
01239         MultiArrayIndex N = ndim();
01240         difference_type stride(PyArray_STRIDES(pyObject()), PyArray_STRIDES(pyObject()) + N),
01241                         permutation(N);
01242         for(MultiArrayIndex k=0; k<N; ++k)
01243             permutation[k] = k;
01244         for(MultiArrayIndex k=0; k<N-1; ++k)
01245         {
01246             MultiArrayIndex smallest = k;
01247             for(MultiArrayIndex j=k+1; j<N; ++j)
01248             {
01249                 if(stride[j] < stride[smallest])
01250                     smallest = j;
01251             }
01252             if(smallest != k)
01253             {
01254                 std::swap(stride[k], stride[smallest]);
01255                 std::swap(permutation[k], permutation[smallest]);
01256             }
01257         }
01258         difference_type ordering(N);
01259         for(MultiArrayIndex k=0; k<N; ++k)
01260             ordering[permutation[k]] = k;
01261         return ordering;
01262     }
01263 
01264         /**
01265          Returns the value type of the elements in this array, or -1
01266          when hasData() is false.
01267          */
01268     int dtype() const
01269     {
01270         if(hasData())
01271             return PyArray_DESCR(pyObject())->type_num;
01272         return -1;
01273     }
01274 
01275         /**
01276          * Return a borrowed reference to the internal PyArrayObject.
01277          */
01278     PyArrayObject * pyArray() const
01279     {
01280         return (PyArrayObject *)pyArray_.get();
01281     }
01282 
01283         /**
01284          * Return a borrowed reference to the internal PyArrayObject
01285          * (see pyArray()), cast to PyObject for your convenience.
01286          */
01287     PyObject * pyObject() const
01288     {
01289         return pyArray_.get();
01290     }
01291 
01292         /**
01293            Reset the NumpyAnyArray to the given object. If \a obj is a numpy array object,
01294            a new reference to that array is created, and the function returns
01295            true. Otherwise, it returns false and the NumpyAnyArray remains unchanged.
01296            If \a type is given, the new reference will be a view with that type, provided
01297            that \a type is a numpy ndarray or a subclass thereof. Otherwise, an
01298            exception is thrown.
01299          */
01300     bool makeReference(PyObject * obj, PyTypeObject * type = 0)
01301     {
01302         if(obj == 0 || !PyArray_Check(obj))
01303             return false;
01304         if(type != 0)
01305         {
01306             vigra_precondition(PyType_IsSubtype(type, &PyArray_Type) != 0,
01307                 "NumpyAnyArray::makeReference(obj, type): type must be numpy.ndarray or a subclass thereof.");
01308             obj = PyArray_View((PyArrayObject*)obj, 0, type);
01309             pythonToCppException(obj);
01310         }
01311         pyArray_.reset(obj);
01312         return true;
01313     }
01314 
01315         /**
01316            Create a copy of the given array object. If \a obj is a numpy array object,
01317            a copy is created via the C-equivalent of 'obj->copy()'. If
01318            this call fails, or obj was not an array, an exception is thrown
01319            and the NumpyAnyArray remains unchanged.
01320          */
01321     void makeCopy(PyObject * obj, PyTypeObject * type = 0)
01322     {
01323         vigra_precondition(obj && PyArray_Check(obj),
01324              "NumpyAnyArray::makeCopy(obj): obj is not an array.");
01325         vigra_precondition(type == 0 || PyType_IsSubtype(type, &PyArray_Type),
01326              "NumpyAnyArray::makeCopy(obj, type): type must be numpy.ndarray or a subclass thereof.");
01327         python_ptr array(PyArray_NewCopy((PyArrayObject*)obj, NPY_ANYORDER), python_ptr::keep_count);
01328         pythonToCppException(array);
01329         makeReference(array, type);
01330     }
01331 
01332          /**
01333            Check whether this NumpyAnyArray actually points to a Python array.
01334          */
01335     bool hasData() const
01336     {
01337         return pyArray_ != 0;
01338     }
01339 };
01340 
01341 /********************************************************/
01342 /*                                                      */
01343 /*                     NumpyArray                       */
01344 /*                                                      */
01345 /********************************************************/
01346 
01347 /** Provide the MultiArrayView interface for a Python array.
01348     
01349     This class inherits from both \ref vigra::MultiArrayView and \ref vigra::NumpyAnyArray 
01350     in order to support easy and save application of VIGRA functions to Python arrays.
01351 
01352     <b>\#include</b> <<a href="numpy__array_8hxx-source.html">vigra/numpy_array.hxx</a>><br>
01353     Namespace: vigra
01354 */
01355 template <unsigned int N, class T, class Stride = StridedArrayTag>
01356 class NumpyArray
01357 : public MultiArrayView<N, typename NumpyArrayTraits<N, T, Stride>::value_type, Stride>,
01358   public NumpyAnyArray
01359 {
01360   public:
01361     typedef NumpyArrayTraits<N, T, Stride> ArrayTraits;
01362     typedef typename ArrayTraits::dtype dtype;
01363     typedef T pseudo_value_type;
01364     
01365     static NPY_TYPES const typeCode = ArrayTraits::typeCode;
01366 
01367         /** the view type associated with this array.
01368          */
01369     typedef MultiArrayView<N, typename ArrayTraits::value_type, Stride> view_type;
01370 
01371     enum { actual_dimension = view_type::actual_dimension };
01372 
01373         /** the array's value type
01374          */
01375     typedef typename view_type::value_type value_type;
01376 
01377         /** pointer type
01378          */
01379     typedef typename view_type::pointer pointer;
01380 
01381         /** const pointer type
01382          */
01383     typedef typename view_type::const_pointer const_pointer;
01384 
01385         /** reference type (result of operator[])
01386          */
01387     typedef typename view_type::reference reference;
01388 
01389         /** const reference type (result of operator[] const)
01390          */
01391     typedef typename view_type::const_reference const_reference;
01392 
01393         /** size type
01394          */
01395     typedef typename view_type::size_type size_type;
01396 
01397         /** difference type (used for multi-dimensional offsets and indices)
01398          */
01399     typedef typename view_type::difference_type difference_type;
01400 
01401         /** difference and index type for a single dimension
01402          */
01403     typedef typename view_type::difference_type_1 difference_type_1;
01404 
01405         /** traverser type
01406          */
01407     typedef typename view_type::traverser traverser;
01408 
01409         /** traverser type to const data
01410          */
01411     typedef typename view_type::const_traverser const_traverser;
01412 
01413         /** sequential (random access) iterator type
01414          */
01415     typedef value_type * iterator;
01416 
01417         /** sequential (random access) const iterator type
01418          */
01419     typedef value_type * const_iterator;
01420 
01421     using view_type::shape;   // resolve ambiguity of multiple inheritance
01422     using view_type::hasData; // resolve ambiguity of multiple inheritance
01423     using view_type::strideOrdering; // resolve ambiguity of multiple inheritance
01424 
01425   protected:
01426 
01427     // this function assumes that pyArray_ has already been set, and compatibility been checked
01428     void setupArrayView();
01429     
01430     static python_ptr getArrayTypeObject()
01431     {
01432         python_ptr type = detail::getArrayTypeObject(ArrayTraits::typeKeyFull());
01433         if(type == 0)
01434             type = detail::getArrayTypeObject(ArrayTraits::typeKey(), &PyArray_Type);
01435         return type;
01436     }
01437        
01438     static python_ptr init(difference_type const & shape, bool init = true)
01439     {
01440         ArrayVector<npy_intp> pshape(shape.begin(), shape.end());
01441         return detail::constructNumpyArrayImpl((PyTypeObject *)getArrayTypeObject().ptr(), pshape, 
01442                        ArrayTraits::spatialDimensions, ArrayTraits::channels,
01443                        typeCode, "V", init);
01444     }
01445 
01446     static python_ptr init(difference_type const & shape, difference_type const & strideOrdering, bool init = true)
01447     {
01448         ArrayVector<npy_intp> pshape(shape.begin(), shape.end()), 
01449                               pstrideOrdering(strideOrdering.begin(), strideOrdering.end());
01450         return detail::constructNumpyArrayImpl((PyTypeObject *)getArrayTypeObject().ptr(), pshape, 
01451                        ArrayTraits::spatialDimensions, ArrayTraits::channels,
01452                        typeCode, "A", init, pstrideOrdering);
01453     }
01454 
01455   public:
01456   
01457     using view_type::init;
01458     
01459         /**
01460          * Construct from a given PyObject pointer. When the given
01461          * python object is NULL, the internal python array will be
01462          * NULL and hasData() will return false.
01463          *
01464          * Otherwise, the function attempts to create a
01465          * new reference to the given Python object, unless
01466          * copying is forced by setting \a createCopy to true.
01467          * If either of this fails, the function throws an exception.
01468          * This will not happen if isStrictlyCompatible(obj) (in case
01469          * of creating a new reference) or isCopyCompatible(obj)
01470          * (in case of copying) have returned true beforehand.
01471          */
01472     explicit NumpyArray(PyObject *obj = 0, bool createCopy = false)
01473     {
01474         if(obj == 0)
01475             return;
01476         if(createCopy)
01477             makeCopy(obj);
01478         else
01479             vigra_precondition(makeReference(obj),
01480                   "NumpyArray(obj): Cannot construct from incompatible array.");
01481     }
01482 
01483        /**
01484          * Copy constructor; does not copy the memory, but creates a
01485          * new reference to the same underlying python object, unless
01486          * a copy is forced by setting \a createCopy to true.
01487          * (If the source object has no data, this one will have
01488          * no data, too.)
01489          */
01490     NumpyArray(const NumpyArray &other, bool createCopy = false) :
01491             MultiArrayView<N, typename NumpyArrayTraits<N, T, Stride>::value_type, Stride>(other),
01492             NumpyAnyArray(other, createCopy)
01493     {
01494         if(!other.hasData())
01495             return;
01496         if(createCopy)
01497             makeCopy(other.pyObject());
01498         else
01499             makeReferenceUnchecked(other.pyObject());
01500     }
01501 
01502        /**
01503          * Allocate new memory and copy data from a MultiArrayView.
01504          */
01505     explicit NumpyArray(const view_type &other)
01506     {
01507         if(!other.hasData())
01508             return;
01509         vigra_postcondition(makeReference(init(other.shape(), false)),
01510                   "NumpyArray(view_type): Python constructor did not produce a compatible array.");
01511         static_cast<view_type &>(*this) = other;
01512     }
01513 
01514         /**
01515          * Construct a new array object, allocating an internal python
01516          * ndarray of the given shape (in fortran order), initialized
01517          * with zeros.
01518          *
01519          * An exception is thrown when construction fails.
01520          */
01521     explicit NumpyArray(difference_type const & shape)
01522     {
01523         vigra_postcondition(makeReference(init(shape)),
01524                      "NumpyArray(shape): Python constructor did not produce a compatible array.");
01525     }
01526 
01527         /**
01528          * Construct a new array object, allocating an internal python
01529          * ndarray of the given shape and given stride ordering, initialized
01530          * with zeros.
01531          *
01532          * An exception is thrown when construction fails.
01533          */
01534     NumpyArray(difference_type const & shape, difference_type const & strideOrdering)
01535     {
01536         vigra_postcondition(makeReference(init(shape, strideOrdering)),
01537                      "NumpyArray(shape): Python constructor did not produce a compatible array.");
01538     }
01539 
01540         /**
01541          * Constructor from NumpyAnyArray.
01542          * Equivalent to NumpyArray(other.pyObject())
01543          */
01544     NumpyArray(const NumpyAnyArray &other, bool createCopy = false)
01545     {
01546         if(!other.hasData())
01547             return;
01548         if(createCopy)
01549             makeCopy(other.pyObject());
01550         else
01551             vigra_precondition(makeReference(other.pyObject()), //, false),
01552                    "NumpyArray(NumpyAnyArray): Cannot construct from incompatible or empty array.");
01553     }
01554 
01555         /**
01556          * Assignment operator. If this is already a view with data
01557          * (i.e. hasData() is true) and the shapes match, the RHS
01558          * array contents are copied.  If this is an empty view,
01559          * assignment is identical to makeReferenceUnchecked(other.pyObject()).
01560          * See MultiArrayView::operator= for further information on
01561          * semantics.
01562          */
01563     NumpyArray &operator=(const NumpyArray &other)
01564     {
01565         if(hasData())
01566             view_type::operator=(other);
01567         else
01568             makeReferenceUnchecked(other.pyObject());
01569         return *this;
01570     }
01571 
01572         /**
01573          * Assignment operator. If this is already a view with data
01574          * (i.e. hasData() is true) and the shapes match, the RHS
01575          * array contents are copied.
01576          * If this is an empty view, assignment is identical to
01577          * makeReference(other.pyObject()).
01578          * Otherwise, an exception is thrown.
01579          */
01580     NumpyArray &operator=(const NumpyAnyArray &other)
01581     {
01582         if(hasData())
01583         {
01584             NumpyAnyArray::operator=(other);
01585         }
01586         else if(isStrictlyCompatible(other.pyObject()))
01587         {
01588             makeReferenceUnchecked(other.pyObject());
01589         }
01590         else
01591         {
01592             vigra_precondition(false,
01593                 "NumpyArray::operator=(): Cannot assign from incompatible array.");
01594         }
01595         return *this;
01596     }
01597 
01598         /**
01599          * Test whether a given python object is a numpy array that can be
01600          * converted (copied) into an array compatible to this NumpyArray type.
01601          * This means that the array's shape conforms to the requirements of
01602          * makeCopy().
01603          */
01604     static bool isCopyCompatible(PyObject *obj)
01605     {
01606         return ArrayTraits::isArray(obj) &&
01607                ArrayTraits::isShapeCompatible((PyArrayObject *)obj);
01608     }
01609 
01610         /**
01611          * Test whether a given python object is a numpy array with a
01612          * compatible dtype and the correct shape and strides, so that it
01613          * can be referenced as a view by this NumpyArray type (i.e.
01614          * it conforms to the requirements of makeReference()).
01615          */
01616     static bool isReferenceCompatible(PyObject *obj)
01617     {
01618         return ArrayTraits::isArray(obj) &&
01619                ArrayTraits::isPropertyCompatible((PyArrayObject *)obj);
01620     }
01621 
01622         /**
01623          * Like isReferenceCompatible(obj), but also executes a customized type compatibility
01624          * check when such a check has been registered for this class via
01625          * registerPythonArrayType().
01626          *
01627          * This facilitates proper overload resolution between
01628          * NumpyArray<3, Multiband<T> > (a multiband image) and NumpyArray<3, Singleband<T> > (a scalar volume).
01629          */
01630     static bool isStrictlyCompatible(PyObject *obj)
01631     {
01632 #if VIGRA_CONVERTER_DEBUG
01633         std::cerr << "class " << typeid(NumpyArray).name() << " got " << obj->ob_type->tp_name << "\n";
01634         bool isClassCompatible=ArrayTraits::isClassCompatible(obj);
01635         bool isPropertyCompatible((PyArrayObject *)obj);
01636         std::cerr<<"isClassCompatible: "<<isClassCompatible<<std::endl;
01637         std::cerr<<"isPropertyCompatible: "<<isPropertyCompatible<<std::endl;
01638 #endif
01639         return ArrayTraits::isClassCompatible(obj) &&
01640                ArrayTraits::isPropertyCompatible((PyArrayObject *)obj);
01641     }
01642 
01643         /**
01644          * Create a vector representing the standard stride ordering of a NumpyArray.
01645          * That is, we get a vector representing the range [0,...,N-1], which
01646          * denotes the stride ordering for Fortran order.
01647          */
01648     static difference_type standardStrideOrdering()
01649     {
01650         difference_type strideOrdering;
01651         for(unsigned int k=0; k<N; ++k)
01652             strideOrdering[k] = k;
01653         return strideOrdering;
01654     }
01655 
01656         /**
01657          * Set up a view to the given object without checking compatibility.
01658          * This function must not be used unless isReferenceCompatible(obj) returned
01659          * true on the given object (otherwise, a crash is likely).
01660          */
01661     void makeReferenceUnchecked(PyObject *obj)
01662     {
01663         NumpyAnyArray::makeReference(obj);
01664         setupArrayView();
01665     }
01666 
01667         /**
01668          * Try to set up a view referencing the given PyObject.
01669          * Returns false if the python object is not a compatible
01670          * numpy array (see isReferenceCompatible() or
01671          * isStrictlyCompatible(), according to the parameter \a
01672          * strict).
01673          */
01674     bool makeReference(PyObject *obj, bool strict = true)
01675     {
01676         if(strict)
01677         {
01678             if(!isStrictlyCompatible(obj))
01679                 return false;
01680         }
01681         else
01682         {
01683             if(!isReferenceCompatible(obj))
01684                 return false;
01685         }
01686         makeReferenceUnchecked(obj);
01687         return true;
01688     }
01689 
01690         /**
01691          * Try to set up a view referencing the same data as the given
01692          * NumpyAnyArray.  This overloaded variant simply calls
01693          * makeReference() on array.pyObject().
01694          */
01695     bool makeReference(const NumpyAnyArray &array, bool strict = true)
01696     {
01697         return makeReference(array.pyObject(), strict);
01698     }
01699 
01700         /**
01701          * Set up an unsafe reference to the given MultiArrayView.
01702          * ATTENTION: This creates a numpy.ndarray that points to the
01703          * same data, but does not own it, so it must be ensured by
01704          * other means that the memory does not get freed before the
01705          * end of the ndarray's lifetime!  (One elegant way would be
01706          * to set the 'base' attribute of the resulting ndarray to a
01707          * python object which directly or indirectly holds the memory
01708          * of the given MultiArrayView.)
01709          */
01710     void makeReference(const view_type &multiArrayView)
01711     {
01712         vigra_precondition(!hasData(), "makeReference(): cannot replace existing view with given buffer");
01713 
01714         // construct an ndarray that points to our data (taking strides into account):
01715         python_ptr array(ArrayTraits::constructor(multiArrayView.shape(), multiArrayView.data(), multiArrayView.stride()));
01716 
01717         view_type::operator=(multiArrayView);
01718         pyArray_ = array;
01719     }
01720 
01721         /**
01722          Try to create a copy of the given PyObject.
01723          Raises an exception when obj is not a compatible array
01724          (see isCopyCompatible() or isStrictlyCompatible(), according to the
01725          parameter \a strict) or the Python constructor call failed.
01726          */
01727     void makeCopy(PyObject *obj, bool strict = false)
01728     {
01729         vigra_precondition(strict ? isStrictlyCompatible(obj) : isCopyCompatible(obj),
01730                      "NumpyArray::makeCopy(obj): Cannot copy an incompatible array.");
01731 
01732         int M = PyArray_NDIM(obj);
01733         TinyVector<npy_intp, N> shape;
01734         std::copy(PyArray_DIMS(obj), PyArray_DIMS(obj)+M, shape.begin());
01735         if(M == N-1)
01736             shape[M] = 1;
01737         vigra_postcondition(makeReference(init(shape, false)),
01738                      "NumpyArray::makeCopy(obj): Copy created an incompatible array.");
01739         NumpyAnyArray::operator=(NumpyAnyArray(obj));
01740 //        if(PyArray_CopyInto(pyArray(), (PyArrayObject*)obj) == -1)
01741 //            pythonToCppException(0);
01742     }
01743 
01744         /**
01745             Allocate new memory with the given shape and initialize with zeros.<br>
01746             If a stride ordering is given, the resulting array will have this stride
01747             ordering, when it is compatible with the array's memory layout (unstrided
01748             arrays only permit the standard ascending stride ordering).
01749 
01750             <em>Note:</em> this operation invalidates dependent objects
01751             (MultiArrayViews and iterators)
01752          */
01753     void reshape(difference_type const & shape, difference_type const & strideOrdering = standardStrideOrdering())
01754     {
01755         vigra_postcondition(makeReference(init(shape, strideOrdering)),
01756                      "NumpyArray(shape): Python constructor did not produce a compatible array.");
01757     }
01758 
01759         /**
01760             When this array has no data, allocate new memory with the given \a shape and
01761             initialize with zeros. Otherwise, check if the new shape matches the old shape
01762             and throw a precondition exception with the given \a message if not.
01763          */
01764     void reshapeIfEmpty(difference_type const & shape, std::string message = "")
01765     {
01766         reshapeIfEmpty(shape, standardStrideOrdering(), message);
01767     }
01768 
01769         /**
01770             When this array has no data, allocate new memory with the given \a shape and
01771             initialize with zeros. Otherwise, check if the new shape matches the old shape
01772             and throw a precondition exception with the given \a message if not. If strict
01773             is true, the given stride ordering must also match that of the existing data.
01774          */
01775     void reshapeIfEmpty(difference_type const & shape, difference_type const & strideOrdering,
01776                         std::string message = "", bool strict = false)
01777     {
01778         if(hasData())
01779         {
01780             if(strict)
01781             {
01782                 if(message == "")
01783                     message = "NumpyArray::reshapeIfEmpty(shape): array was not empty, and shape or stride ordering did not match.";
01784                 vigra_precondition(shape == this->shape() && strideOrdering == this->strideOrdering(), message.c_str());
01785             }
01786             else
01787             {
01788                 if(message == "")
01789                     message = "NumpyArray::reshapeIfEmpty(shape): array was not empty, and shape did not match.";
01790                 vigra_precondition(shape == this->shape(), message.c_str());
01791             }
01792         }
01793         else
01794         {
01795             reshape(shape, strideOrdering);
01796         }
01797     }
01798 };
01799 
01800     // this function assumes that pyArray_ has already been set, and compatibility been checked
01801 template <unsigned int N, class T, class Stride>
01802 void NumpyArray<N, T, Stride>::setupArrayView()
01803 {
01804     if(NumpyAnyArray::hasData())
01805     {
01806         unsigned int dimension = std::min<unsigned int>(actual_dimension, pyArray()->nd);
01807         std::copy(pyArray()->dimensions, pyArray()->dimensions + dimension, this->m_shape.begin());
01808         std::copy(pyArray()->strides, pyArray()->strides + dimension, this->m_stride.begin());
01809         if(pyArray()->nd < actual_dimension)
01810         {
01811             this->m_shape[dimension] = 1;
01812             this->m_stride[dimension] = sizeof(value_type);
01813         }
01814         this->m_stride /= sizeof(value_type);
01815         this->m_ptr = reinterpret_cast<pointer>(pyArray()->data);
01816     }
01817     else
01818     {
01819         this->m_ptr = 0;
01820     }
01821 }
01822 
01823 
01824 typedef NumpyArray<2, float >  NumpyFArray2;
01825 typedef NumpyArray<3, float >  NumpyFArray3;
01826 typedef NumpyArray<4, float >  NumpyFArray4;
01827 typedef NumpyArray<2, Singleband<float> >  NumpyFImage;
01828 typedef NumpyArray<3, Singleband<float> >  NumpyFVolume;
01829 typedef NumpyArray<2, RGBValue<float> >  NumpyFRGBImage;
01830 typedef NumpyArray<3, RGBValue<float> >  NumpyFRGBVolume;
01831 typedef NumpyArray<3, Multiband<float> >  NumpyFMultibandImage;
01832 typedef NumpyArray<4, Multiband<float> >  NumpyFMultibandVolume;
01833 
01834 inline void import_vigranumpy()
01835 {
01836     if(_import_array() < 0)
01837         pythonToCppException(0);
01838     python_ptr module(PyImport_ImportModule("vigra.vigranumpycore"), python_ptr::keep_count);
01839     pythonToCppException(module);
01840 }
01841 
01842 /********************************************************/
01843 /*                                                      */
01844 /*   NumpyArray Multiband Argument Object Factories     */
01845 /*                                                      */
01846 /********************************************************/
01847 
01848 template <class PixelType, class Stride>
01849 inline triple<ConstStridedImageIterator<PixelType>,
01850               ConstStridedImageIterator<PixelType>,
01851               MultibandVectorAccessor<PixelType> >
01852 srcImageRange(NumpyArray<3, Multiband<PixelType>, Stride> const & img)
01853 {
01854     ConstStridedImageIterator<PixelType>
01855         ul(img.data(), 1, img.stride(0), img.stride(1));
01856     return triple<ConstStridedImageIterator<PixelType>,
01857                   ConstStridedImageIterator<PixelType>,
01858                   MultibandVectorAccessor<PixelType> >
01859         (ul, ul + Size2D(img.shape(0), img.shape(1)), MultibandVectorAccessor<PixelType>(img.shape(2), img.stride(2)));
01860 }
01861 
01862 template <class PixelType, class Stride>
01863 inline pair< ConstStridedImageIterator<PixelType>,
01864              MultibandVectorAccessor<PixelType> >
01865 srcImage(NumpyArray<3, Multiband<PixelType>, Stride> const & img)
01866 {
01867     ConstStridedImageIterator<PixelType>
01868         ul(img.data(), 1, img.stride(0), img.stride(1));
01869     return pair<ConstStridedImageIterator<PixelType>, MultibandVectorAccessor<PixelType> >
01870         (ul, MultibandVectorAccessor<PixelType>(img.shape(2), img.stride(2)));
01871 }
01872 
01873 template <class PixelType, class Stride>
01874 inline triple< StridedImageIterator<PixelType>,
01875                StridedImageIterator<PixelType>,
01876                MultibandVectorAccessor<PixelType> >
01877 destImageRange(NumpyArray<3, Multiband<PixelType>, Stride> & img)
01878 {
01879     StridedImageIterator<PixelType>
01880         ul(img.data(), 1, img.stride(0), img.stride(1));
01881     typedef typename AccessorTraits<PixelType>::default_accessor Accessor;
01882     return triple<StridedImageIterator<PixelType>,
01883                   StridedImageIterator<PixelType>,
01884                   MultibandVectorAccessor<PixelType> >
01885         (ul, ul + Size2D(img.shape(0), img.shape(1)),
01886         MultibandVectorAccessor<PixelType>(img.shape(2), img.stride(2)));
01887 }
01888 
01889 template <class PixelType, class Stride>
01890 inline pair< StridedImageIterator<PixelType>,
01891              MultibandVectorAccessor<PixelType> >
01892 destImage(NumpyArray<3, Multiband<PixelType>, Stride> & img)
01893 {
01894     StridedImageIterator<PixelType>
01895         ul(img.data(), 1, img.stride(0), img.stride(1));
01896     return pair<StridedImageIterator<PixelType>, MultibandVectorAccessor<PixelType> >
01897         (ul, MultibandVectorAccessor<PixelType>(img.shape(2), img.stride(2)));
01898 }
01899 
01900 template <class PixelType, class Stride>
01901 inline pair< ConstStridedImageIterator<PixelType>,
01902              MultibandVectorAccessor<PixelType> >
01903 maskImage(NumpyArray<3, Multiband<PixelType>, Stride> const & img)
01904 {
01905     ConstStridedImageIterator<PixelType>
01906         ul(img.data(), 1, img.stride(0), img.stride(1));
01907     typedef typename AccessorTraits<PixelType>::default_accessor Accessor;
01908     return pair<ConstStridedImageIterator<PixelType>, MultibandVectorAccessor<PixelType> >
01909         (ul, MultibandVectorAccessor<PixelType>(img.shape(2), img.stride(2)));
01910 }
01911 
01912 } // namespace vigra
01913 
01914 #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.0 (Thu Aug 25 2011)