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

vigra/multi_array.hxx
00001 /************************************************************************/
00002 /*                                                                      */
00003 /*     Copyright 2003-2008 by Gunnar Kedenburg and Ullrich Koethe       */
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 
00037 #ifndef VIGRA_MULTI_ARRAY_HXX
00038 #define VIGRA_MULTI_ARRAY_HXX
00039 
00040 #include <memory>
00041 #include <algorithm>
00042 #include "accessor.hxx"
00043 #include "tinyvector.hxx"
00044 #include "rgbvalue.hxx"
00045 #include "basicimageview.hxx"
00046 #include "imageiterator.hxx"
00047 #include "numerictraits.hxx"
00048 #include "multi_iterator.hxx"
00049 #include "metaprogramming.hxx"
00050 #include "mathutil.hxx"
00051 
00052 // Bounds checking Macro used if VIGRA_CHECK_BOUNDS is defined.
00053 #ifdef VIGRA_CHECK_BOUNDS
00054 #define VIGRA_ASSERT_INSIDE(diff) \
00055   vigra_precondition(this->isInside(diff), "Index out of bounds")
00056 #else
00057 #define VIGRA_ASSERT_INSIDE(diff)
00058 #endif
00059 
00060 namespace vigra
00061 {
00062 
00063 namespace detail
00064 {
00065 /********************************************************/
00066 /*                                                      */
00067 /*                    defaultStride                     */
00068 /*                                                      */
00069 /********************************************************/
00070 
00071 /* generates the stride for a gapless shape.
00072 
00073     Namespace: vigra::detail
00074 */
00075 template <unsigned int N>
00076 inline TinyVector <MultiArrayIndex, N>
00077 defaultStride(const TinyVector <MultiArrayIndex, N> &shape)
00078 {
00079     TinyVector <MultiArrayIndex, N> ret;
00080     ret [0] = 1;
00081     for (int i = 1; i < (int)N; ++i)
00082         ret [i] = ret [i-1] * shape [i-1];
00083     return ret;
00084 }
00085 
00086 /********************************************************/
00087 /*                                                      */
00088 /*                 ScanOrderToOffset                    */
00089 /*                                                      */
00090 /********************************************************/
00091 
00092 /* transforms an index in scan order sense to a pointer offset in a possibly
00093    strided, multi-dimensional array.
00094 
00095     Namespace: vigra::detail
00096 */
00097 
00098 template <int K>
00099 struct ScanOrderToOffset
00100 {
00101     template <int N>
00102     static MultiArrayIndex
00103     exec(MultiArrayIndex d, const TinyVector <MultiArrayIndex, N> &shape,
00104          const TinyVector <MultiArrayIndex, N> & stride)
00105     {
00106         return stride[N-K] * (d % shape[N-K]) +
00107                ScanOrderToOffset<K-1>::exec(d / shape[N-K], shape, stride);
00108     }
00109 };
00110 
00111 template <>
00112 struct ScanOrderToOffset<1>
00113 {
00114     template <int N>
00115     static MultiArrayIndex
00116     exec(MultiArrayIndex d, const TinyVector <MultiArrayIndex, N> & /*shape*/,
00117          const TinyVector <MultiArrayIndex, N> & stride)
00118     {
00119         return stride[N-1] * d;
00120     }
00121 };
00122 
00123 template <int K>
00124 struct ScanOrderToCoordinate
00125 {
00126     template <int N>
00127     static void
00128     exec(MultiArrayIndex d, const TinyVector <MultiArrayIndex, N> &shape,
00129          TinyVector <MultiArrayIndex, N> & result)
00130     {
00131         result[N-K] = (d % shape[N-K]);
00132         ScanOrderToCoordinate<K-1>::exec(d / shape[N-K], shape, result);
00133     }
00134 };
00135 
00136 template <>
00137 struct ScanOrderToCoordinate<1>
00138 {
00139     template <int N>
00140     static void
00141     exec(MultiArrayIndex d, const TinyVector <MultiArrayIndex, N> & /*shape*/,
00142          TinyVector <MultiArrayIndex, N> & result)
00143     {
00144         result[N-1] = d;
00145     }
00146 };
00147 
00148 template <int K>
00149 struct CoordinateToScanOrder
00150 {
00151     template <int N>
00152     static MultiArrayIndex
00153     exec(const TinyVector <MultiArrayIndex, N> &shape,
00154          const TinyVector <MultiArrayIndex, N> & coordinate)
00155     {
00156         return coordinate[N-K] + shape[N-K] * CoordinateToScanOrder<K-1>::exec(shape, coordinate);
00157     }
00158 };
00159 
00160 template <>
00161 struct CoordinateToScanOrder<1>
00162 {
00163     template <int N>
00164     static MultiArrayIndex
00165     exec(const TinyVector <MultiArrayIndex, N> & /*shape*/,
00166          const TinyVector <MultiArrayIndex, N> & coordinate)
00167     {
00168         return coordinate[N-1];
00169     }
00170 };
00171 
00172 
00173 template <class C>
00174 struct CoordinatesToOffest
00175 {
00176     template <int N>
00177     static MultiArrayIndex
00178     exec(const TinyVector <MultiArrayIndex, N> & stride, MultiArrayIndex x)
00179     {
00180         return stride[0] * x;
00181     }
00182     template <int N>
00183     static MultiArrayIndex
00184     exec(const TinyVector <MultiArrayIndex, N> & stride, MultiArrayIndex x, MultiArrayIndex y)
00185     {
00186         return stride[0] * x + stride[1] * y;
00187     }
00188 };
00189 
00190 template <>
00191 struct CoordinatesToOffest<UnstridedArrayTag>
00192 {
00193     template <int N>
00194     static MultiArrayIndex
00195     exec(const TinyVector <MultiArrayIndex, N> & /*stride*/, MultiArrayIndex x)
00196     {
00197         return x;
00198     }
00199     template <int N>
00200     static MultiArrayIndex
00201     exec(const TinyVector <MultiArrayIndex, N> & stride, MultiArrayIndex x, MultiArrayIndex y)
00202     {
00203         return x + stride[1] * y;
00204     }
00205 };
00206 
00207 /********************************************************/
00208 /*                                                      */
00209 /*                     MaybeStrided                     */
00210 /*                                                      */
00211 /********************************************************/
00212 
00213 /* metatag implementing a test for marking MultiArrays that were
00214     indexed at the zero'th dimension as strided, and all others as
00215     unstrided.
00216 
00217 <b>\#include</b>
00218 <<a href="multi__array_8hxx-source.html">vigra/multi_array.hxx</a>>
00219 
00220 Namespace: vigra::detail
00221 */
00222 template <unsigned int N>
00223 struct MaybeStrided
00224 {
00225     typedef UnstridedArrayTag type;
00226 };
00227 
00228 template <>
00229 struct MaybeStrided <0>
00230 {
00231     typedef StridedArrayTag type;
00232 };
00233 
00234 /********************************************************/
00235 /*                                                      */
00236 /*                MultiIteratorChooser                  */
00237 /*                                                      */
00238 /********************************************************/
00239 
00240 /* metatag implementing a test (by pattern matching) for marking
00241     MultiArrays that were indexed at the zero'th dimension as strided.
00242 
00243 <b>\#include</b>
00244 <<a href="multi__array_8hxx-source.html">vigra/multi_array.hxx</a>>
00245 
00246 Namespace: vigra::detail
00247 */
00248 template <class O>
00249 struct MultiIteratorChooser
00250 {
00251     struct Nil {};
00252 
00253     template <unsigned int N, class T, class REFERENCE, class POINTER>
00254     struct Traverser
00255     {
00256         typedef Nil type;
00257     };
00258 };
00259 
00260 /********************************************************/
00261 /*                                                      */
00262 /*       MultiIteratorChooser <StridedArrayTag>         */
00263 /*                                                      */
00264 /********************************************************/
00265 
00266 /* specialization of the MultiIteratorChooser for strided arrays.
00267 
00268 <b>\#include</b>
00269 <<a href="multi__array_8hxx-source.html">vigra/multi_array.hxx</a>>
00270 
00271 Namespace: vigra::detail
00272 */
00273 template <>
00274 struct MultiIteratorChooser <StridedArrayTag>
00275 {
00276     template <unsigned int N, class T, class REFERENCE, class POINTER>
00277     struct Traverser
00278     {
00279         typedef StridedMultiIterator <N, T, REFERENCE, POINTER> type;
00280     };
00281 };
00282 
00283 /********************************************************/
00284 /*                                                      */
00285 /*      MultiIteratorChooser <UnstridedArrayTag>        */
00286 /*                                                      */
00287 /********************************************************/
00288 
00289 /* specialization of the MultiIteratorChooser for unstrided arrays.
00290 
00291 <b>\#include</b>
00292 <<a href="multi__array_8hxx-source.html">vigra/multi_array.hxx</a>>
00293 
00294 Namespace: vigra::detail
00295 */
00296 template <>
00297 struct MultiIteratorChooser <UnstridedArrayTag>
00298 {
00299     template <unsigned int N, class T, class REFERENCE, class POINTER>
00300     struct Traverser
00301     {
00302         typedef MultiIterator <N, T, REFERENCE, POINTER> type;
00303     };
00304 };
00305 
00306 /********************************************************/
00307 /*                                                      */
00308 /*                   helper functions                   */
00309 /*                                                      */
00310 /********************************************************/
00311 
00312 template <class DestIterator, class Shape, class T>
00313 inline void
00314 initMultiArrayData(DestIterator d, Shape const & shape, T const & init, MetaInt<0>)
00315 {
00316     DestIterator dend = d + shape[0];
00317     for(; d < dend; ++d)
00318     {
00319         *d = init;
00320     }
00321 }
00322 
00323 template <class DestIterator, class Shape, class T, int N>
00324 void
00325 initMultiArrayData(DestIterator d, Shape const & shape, T const & init, MetaInt<N>)
00326 {
00327     DestIterator dend = d + shape[N];
00328     for(; d < dend; ++d)
00329     {
00330         initMultiArrayData(d.begin(), shape, init, MetaInt<N-1>());
00331     }
00332 }
00333 
00334 #define VIGRA_COPY_MULTI_ARRAY_DATA(name, op) \
00335 template <class SrcIterator, class Shape, class DestIterator> \
00336 inline void \
00337 name##MultiArrayData(SrcIterator s, Shape const & shape, DestIterator d, MetaInt<0>) \
00338 {     \
00339     SrcIterator send = s + shape[0]; \
00340     for(; s < send; ++s, ++d) \
00341     { \
00342         *d op detail::RequiresExplicitCast<typename DestIterator::value_type>::cast(*s); \
00343     } \
00344 } \
00345  \
00346 template <class SrcIterator, class Shape, class DestIterator, int N> \
00347 void \
00348 name##MultiArrayData(SrcIterator s, Shape const & shape, DestIterator d, MetaInt<N>) \
00349 { \
00350     SrcIterator send = s + shape[N]; \
00351     for(; s < send; ++s, ++d) \
00352     { \
00353         name##MultiArrayData(s.begin(), shape, d.begin(), MetaInt<N-1>()); \
00354     } \
00355 } \
00356 \
00357 template <class DestIterator, class Shape, class T> \
00358 inline void \
00359 name##ScalarMultiArrayData(DestIterator d, Shape const & shape, T const & init, MetaInt<0>) \
00360 {     \
00361     DestIterator dend = d + shape[0]; \
00362     for(; d < dend; ++d) \
00363     { \
00364         *d op detail::RequiresExplicitCast<typename DestIterator::value_type>::cast(init); \
00365     } \
00366 } \
00367  \
00368 template <class DestIterator, class Shape, class T, int N> \
00369 void \
00370 name##ScalarMultiArrayData(DestIterator d, Shape const & shape, T const & init, MetaInt<N>) \
00371 {     \
00372     DestIterator dend = d + shape[N]; \
00373     for(; d < dend; ++d) \
00374     { \
00375         name##ScalarMultiArrayData(d.begin(), shape, init, MetaInt<N-1>()); \
00376     } \
00377 }
00378 
00379 VIGRA_COPY_MULTI_ARRAY_DATA(copy, =)
00380 VIGRA_COPY_MULTI_ARRAY_DATA(copyAdd, +=)
00381 VIGRA_COPY_MULTI_ARRAY_DATA(copySub, -=)
00382 VIGRA_COPY_MULTI_ARRAY_DATA(copyMul, *=)
00383 VIGRA_COPY_MULTI_ARRAY_DATA(copyDiv, /=)
00384 
00385 #undef VIGRA_COPY_MULTI_ARRAY_DATA
00386 
00387 template <class SrcIterator, class Shape, class T, class ALLOC>
00388 inline void
00389 uninitializedCopyMultiArrayData(SrcIterator s, Shape const & shape, T * & d, ALLOC & a, MetaInt<0>)
00390 {
00391     SrcIterator send = s + shape[0];
00392     for(; s < send; ++s, ++d)
00393     {
00394         a.construct(d, static_cast<T const &>(*s));
00395     }
00396 }
00397 
00398 template <class SrcIterator, class Shape, class T, class ALLOC, int N>
00399 void
00400 uninitializedCopyMultiArrayData(SrcIterator s, Shape const & shape, T * & d, ALLOC & a, MetaInt<N>)
00401 {
00402     SrcIterator send = s + shape[N];
00403     for(; s < send; ++s)
00404     {
00405         uninitializedCopyMultiArrayData(s.begin(), shape, d, a, MetaInt<N-1>());
00406     }
00407 }
00408 
00409 template <class SrcIterator, class Shape, class T>
00410 inline void
00411 normMaxOfMultiArray(SrcIterator s, Shape const & shape, T & result, MetaInt<0>)
00412 {
00413     SrcIterator send = s + shape[0];
00414     for(; s < send; ++s)
00415     {
00416         T v = norm(*s);
00417         if(result < v)
00418             result = v;
00419     }
00420 }
00421 
00422 template <class SrcIterator, class Shape, class T, int N>
00423 void
00424 normMaxOfMultiArray(SrcIterator s, Shape const & shape, T & result, MetaInt<N>)
00425 {
00426     SrcIterator send = s + shape[N];
00427     for(; s < send; ++s)
00428     {
00429         normMaxOfMultiArray(s.begin(), shape, result, MetaInt<N-1>());
00430     }
00431 }
00432 
00433 template <class T>
00434 struct MultiArrayL1Functor
00435 {
00436     template <class U>
00437     T operator()(U t) const
00438     { return norm(t); }
00439 };
00440 
00441 template <class T>
00442 struct MultiArrayL2Functor
00443 {
00444     template <class U>
00445     T operator()(U t) const
00446     { return squaredNorm(t); }
00447 };
00448 
00449 template <class T>
00450 struct MultiArrayScaledL2Functor
00451 {
00452     T scale;
00453 
00454     MultiArrayScaledL2Functor(T s)
00455     : scale(s)
00456     {}
00457 
00458     template <class U>
00459     T operator()(U t) const
00460     { return squaredNorm(T(t) / scale); }
00461 };
00462 
00463 template <class SrcIterator, class Shape, class Functor, class T>
00464 inline void
00465 sumOverMultiArray(SrcIterator s, Shape const & shape, Functor f, T & result, MetaInt<0>)
00466 {
00467     SrcIterator send = s + shape[0];
00468     for(; s < send; ++s)
00469     {
00470         result += f(*s);
00471     }
00472 }
00473 
00474 template <class SrcIterator, class Shape, class Functor, class T, int N>
00475 void
00476 sumOverMultiArray(SrcIterator s, Shape const & shape, Functor f, T & result, MetaInt<N>)
00477 {
00478     SrcIterator send = s + shape[N];
00479     for(; s < send; ++s)
00480     {
00481         sumOverMultiArray(s.begin(), shape, f, result, MetaInt<N-1>());
00482     }
00483 }
00484 
00485 template <class SrcIterator, class Shape, class DestIterator>
00486 inline bool
00487 equalityOfMultiArrays(SrcIterator s, Shape const & shape, DestIterator d, MetaInt<0>)
00488 {
00489     SrcIterator send = s + shape[0];
00490     for(; s < send; ++s, ++d)
00491     {
00492         if(!(*s == *d))
00493             return false;
00494     }
00495     return true;
00496 }
00497 
00498 template <class SrcIterator, class Shape, class DestIterator, int N>
00499 bool
00500 equalityOfMultiArrays(SrcIterator s, Shape const & shape, DestIterator d, MetaInt<N>)
00501 {
00502     SrcIterator send = s + shape[N];
00503     for(; s < send; ++s, ++d)
00504     {
00505         if(!equalityOfMultiArrays(s.begin(), shape, d.begin(), MetaInt<N-1>()))
00506             return false;
00507     }
00508     return true;
00509 }
00510 
00511 
00512 template <class SrcIterator, class Shape, class DestIterator>
00513 inline void
00514 swapDataImpl(SrcIterator s, Shape const & shape, DestIterator d, MetaInt<0>)
00515 {
00516     SrcIterator send = s + shape[0];
00517     for(; s < send; ++s, ++d)
00518         std::swap(*s, *d);
00519 }
00520 
00521 template <class SrcIterator, class Shape, class DestIterator, int N>
00522 void
00523 swapDataImpl(SrcIterator s, Shape const & shape, DestIterator d, MetaInt<N>)
00524 {
00525     SrcIterator send = s + shape[N];
00526     for(; s < send; ++s, ++d)
00527         swapDataImpl(s.begin(), shape, d.begin(), MetaInt<N-1>());
00528 }
00529 
00530 } // namespace detail
00531 
00532 /********************************************************/
00533 /*                                                      */
00534 /*                     MultiArrayView                   */
00535 /*                                                      */
00536 /********************************************************/
00537 
00538 // forward declaration
00539 template <unsigned int N, class T, class C = UnstridedArrayTag>
00540 class MultiArrayView;
00541 template <unsigned int N, class T, class A = std::allocator<T> >
00542 class MultiArray;
00543 
00544 /********************************************************/
00545 /*                                                      */
00546 /*                       NormTraits                     */
00547 /*                                                      */
00548 /********************************************************/
00549 
00550 template <unsigned int N, class T, class C>
00551 struct NormTraits<MultiArrayView<N, T, C> >
00552 {
00553     typedef MultiArrayView<N, T, C>                                      Type;
00554     typedef typename NormTraits<T>::SquaredNormType SquaredNormType;
00555     typedef typename SquareRootTraits<SquaredNormType>::SquareRootResult NormType;
00556 };
00557 
00558 template <unsigned int N, class T, class A>
00559 struct NormTraits<MultiArray<N, T, A> >
00560 : public NormTraits<MultiArrayView<N, T, UnstridedArrayTag> >
00561 {
00562     typedef NormTraits<MultiArrayView<N, T, UnstridedArrayTag> > BaseType;
00563     typedef MultiArray<N, T, A>                                  Type;
00564     typedef typename BaseType::SquaredNormType                   SquaredNormType;
00565     typedef typename BaseType::NormType                          NormType;
00566 };
00567 
00568 /** \brief Base class for, and view to, \ref vigra::MultiArray.
00569 
00570 This class implements the interface of both MultiArray and
00571 MultiArrayView.  By default, MultiArrayViews are tagged as
00572 unstrided. If necessary, strided arrays are constructed automatically
00573 by calls to a variant of the bind...() function.
00574 
00575 If you want to apply an algorithm requiring an image to a
00576 <tt>MultiArrayView</tt> of appropriate (2-dimensional) shape, you can
00577 create a \ref vigra::BasicImageView that acts as a wrapper with the
00578 necessary interface -- see \ref MultiArrayToImage.
00579 
00580 The template parameter are as follows
00581 \code
00582     N: the array dimension
00583 
00584     T: the type of the array elements
00585 
00586     C: a tag determining whether the array's inner dimension is strided
00587        or not. An array is unstrided if the array elements occupy consecutive
00588        memory location, strided if there is an offset in between (e.g.
00589        when a view is created that skips every other array element).
00590        The compiler can generate faster code for unstrided arrays.
00591        Possible values: UnstridedArrayTag (default), StridedArrayTag
00592 \endcode
00593 
00594 <b>\#include</b>
00595 <<a href="multi__array_8hxx-source.html">vigra/multi_array.hxx</a>>
00596 
00597 Namespace: vigra
00598 */
00599 template <unsigned int N, class T, class C>
00600 class MultiArrayView
00601 {
00602 public:
00603 
00604         /** the array's actual dimensionality.
00605             This ensures that MultiArrayView can also be used for
00606             scalars (that is, when <tt>N == 0</tt>). Calculated as:<br>
00607             \code
00608             actual_dimension = (N==0) ? 1 : N
00609             \endcode
00610          */
00611     enum ActualDimension { actual_dimension = (N==0) ? 1 : N };
00612 
00613         /** the array's value type
00614          */
00615     typedef T value_type;
00616 
00617         /** reference type (result of operator[])
00618          */
00619     typedef value_type &reference;
00620 
00621         /** const reference type (result of operator[] const)
00622          */
00623     typedef const value_type &const_reference;
00624 
00625         /** pointer type
00626          */
00627     typedef value_type *pointer;
00628 
00629         /** const pointer type
00630          */
00631     typedef const value_type *const_pointer;
00632 
00633         /** difference type (used for multi-dimensional offsets and indices)
00634          */
00635     typedef typename MultiArrayShape<actual_dimension>::type difference_type;
00636 
00637         /** size type
00638          */
00639     typedef difference_type size_type;
00640 
00641         /** difference and index type for a single dimension
00642          */
00643     typedef MultiArrayIndex difference_type_1;
00644 
00645         /** traverser (MultiIterator) type
00646          */
00647     typedef typename vigra::detail::MultiIteratorChooser <
00648         C>::template Traverser <actual_dimension, T, T &, T *>::type traverser;
00649 
00650         /** const traverser (MultiIterator) type
00651          */
00652     typedef typename vigra::detail::MultiIteratorChooser <
00653         C>::template Traverser <actual_dimension, T, T const &, T const *>::type const_traverser;
00654 
00655         /** the view type associated with this array.
00656          */
00657     typedef MultiArrayView <N, T, C> view_type;
00658 
00659         /** the matrix type associated with this array.
00660          */
00661     typedef MultiArray <N, T> matrix_type;
00662 
00663 protected:
00664 
00665     typedef typename difference_type::value_type diff_zero_t;
00666 
00667         /** the shape of the image pointed to is stored here.
00668         */
00669     difference_type m_shape;
00670 
00671         /** the strides (offset of a sample to the next) for every dimension
00672             are stored here.
00673         */
00674     difference_type m_stride;
00675 
00676         /** pointer to the image.
00677          */
00678     pointer m_ptr;
00679 
00680     template <class U, class CN>
00681     void copyImpl(const MultiArrayView <N, U, CN>& rhs);
00682 
00683     template <class U, class CN>
00684     void swapDataImpl(MultiArrayView <N, U, CN> rhs);
00685 
00686     template <class CN>
00687     bool arraysOverlap(const MultiArrayView <N, T, CN>& rhs) const;
00688 
00689     template <class U, class CN>
00690     bool arraysOverlap(const MultiArrayView <N, U, CN>&) const
00691     {
00692         return false;
00693     }
00694 
00695 public:
00696 
00697         /** default constructor: create an invalid view,
00698          * i.e. hasData() returns false and size() is zero.
00699          */
00700     MultiArrayView ()
00701         : m_shape (diff_zero_t(0)), m_stride (diff_zero_t(0)), m_ptr (0)
00702     {}
00703 
00704         /** construct from shape and pointer
00705          */
00706     MultiArrayView (const difference_type &shape, pointer ptr)
00707     : m_shape (shape),
00708       m_stride (detail::defaultStride <MultiArrayView<N,T>::actual_dimension> (shape)),
00709       m_ptr (ptr)
00710     {}
00711 
00712         /** Construct from shape, strides (offset of a sample to the
00713             next) for every dimension, and pointer.  (Note that
00714             strides are not given in bytes, but in offsets of the
00715             respective pointer type.)
00716          */
00717     MultiArrayView (const difference_type &shape,
00718                     const difference_type &stride,
00719                     pointer ptr)
00720     : m_shape (shape),
00721       m_stride (stride),
00722       m_ptr (ptr)
00723     {}
00724 
00725 
00726         /** Assignment. There are 3 cases:
00727 
00728             <ul>
00729             <li> When this <tt>MultiArrayView</tt> does not point to valid data
00730                  (e.g. after default construction), it becomes a copy of \a rhs.
00731             <li> When the shapes of the two arrays match, the array contents are copied.
00732             <li> Otherwise, a <tt>PreconditionViolation</tt> exception is thrown.
00733             </ul>
00734          */
00735     MultiArrayView & operator=(MultiArrayView const & rhs);
00736 
00737         /** Assignment of a differently typed MultiArrayView. Fails with
00738             <tt>PreconditionViolation</tt> exception when the shapes do not match.
00739          */
00740     template<class U, class C1>
00741     MultiArrayView & operator=(MultiArrayView<N, U, C1> const & rhs)
00742     {
00743         vigra_precondition(this->shape() == rhs.shape(),
00744             "MultiArrayView::operator=() size mismatch.");
00745         this->copyImpl(rhs);
00746         return *this;
00747     }
00748 
00749         /** Add-assignment of a compatible MultiArrayView. Fails with
00750             <tt>PreconditionViolation</tt> exception when the shapes do not match.
00751          */
00752     template<class U, class C1>
00753     MultiArrayView & operator+=(MultiArrayView<N, U, C1> const & rhs);
00754 
00755         /** Subtract-assignment of a compatible MultiArrayView. Fails with
00756             <tt>PreconditionViolation</tt> exception when the shapes do not match.
00757          */
00758     template<class U, class C1>
00759     MultiArrayView & operator-=(MultiArrayView<N, U, C1> const & rhs);
00760 
00761         /** Multiply-assignment of a compatible MultiArrayView. Fails with
00762             <tt>PreconditionViolation</tt> exception when the shapes do not match.
00763          */
00764     template<class U, class C1>
00765     MultiArrayView & operator*=(MultiArrayView<N, U, C1> const & rhs);
00766 
00767         /** Divide-assignment of a compatible MultiArrayView. Fails with
00768             <tt>PreconditionViolation</tt> exception when the shapes do not match.
00769          */
00770     template<class U, class C1>
00771     MultiArrayView & operator/=(MultiArrayView<N, U, C1> const & rhs);
00772 
00773         /** Add-assignment of a scalar.
00774          */
00775     MultiArrayView & operator+=(T const & rhs)
00776     {
00777         detail::copyAddScalarMultiArrayData(traverser_begin(), shape(), rhs, MetaInt<actual_dimension-1>());
00778         return *this;
00779     }
00780 
00781         /** Subtract-assignment of a scalar.
00782          */
00783     MultiArrayView & operator-=(T const & rhs)
00784     {
00785         detail::copySubScalarMultiArrayData(traverser_begin(), shape(), rhs, MetaInt<actual_dimension-1>());
00786         return *this;
00787     }
00788 
00789         /** Multiply-assignment of a scalar.
00790          */
00791     MultiArrayView & operator*=(T const & rhs)
00792     {
00793         detail::copyMulScalarMultiArrayData(traverser_begin(), shape(), rhs, MetaInt<actual_dimension-1>());
00794         return *this;
00795     }
00796 
00797         /** Divide-assignment of a scalar.
00798          */
00799     MultiArrayView & operator/=(T const & rhs)
00800     {
00801         detail::copyDivScalarMultiArrayData(traverser_begin(), shape(), rhs, MetaInt<actual_dimension-1>());
00802         return *this;
00803     }
00804 
00805         /** array access.
00806          */
00807     reference operator[] (const difference_type &d)
00808     {
00809         VIGRA_ASSERT_INSIDE(d);
00810         return m_ptr [dot (d, m_stride)];
00811     }
00812 
00813         /** array access.
00814          */
00815     const_reference operator[] (const difference_type &d) const
00816     {
00817         VIGRA_ASSERT_INSIDE(d);
00818         return m_ptr [dot (d, m_stride)];
00819     }
00820 
00821         /** equvalent to bindInner(), when M < N.
00822          */
00823     template <unsigned int M>
00824     MultiArrayView <N-M, T, StridedArrayTag> operator[] (const TinyVector<MultiArrayIndex, M> &d) const
00825     {
00826         return bindInner(d);
00827     }
00828 
00829         /** array access in scan-order sense.
00830             Mostly useful to support standard indexing for 1-dimensional multi-arrays,
00831             but works for any N. Use scanOrderIndexToCoordinate() and
00832             coordinateToScanOrderIndex() for conversion between indices and coordinates.
00833          */
00834     reference operator[](difference_type_1 d)
00835     {
00836         VIGRA_ASSERT_INSIDE(scanOrderIndexToCoordinate(d));
00837         return m_ptr [detail::ScanOrderToOffset<actual_dimension>::exec(d, m_shape, m_stride)];
00838     }
00839 
00840         /** array access in scan-order sense.
00841             Mostly useful to support standard indexing for 1-dimensional multi-arrays,
00842             but works for any N. Use scanOrderIndexToCoordinate() and
00843             coordinateToScanOrderIndex() for conversion between indices and coordinates.
00844          */
00845     const_reference operator[](difference_type_1 d) const
00846     {
00847         VIGRA_ASSERT_INSIDE(scanOrderIndexToCoordinate(d));
00848         return m_ptr [detail::ScanOrderToOffset<actual_dimension>::exec(d, m_shape, m_stride)];
00849     }
00850 
00851         /** convert scan-order index to coordinate.
00852          */
00853     difference_type scanOrderIndexToCoordinate(difference_type_1 d) const
00854     {
00855         difference_type result;
00856         detail::ScanOrderToCoordinate<actual_dimension>::exec(d, m_shape, result);
00857         return result;
00858     }
00859 
00860         /** convert coordinate to scan-order index.
00861          */
00862     difference_type_1 coordinateToScanOrderIndex(const difference_type &d) const
00863     {
00864         return detail::CoordinateToScanOrder<actual_dimension>::exec(m_shape, d);
00865     }
00866 
00867         /** 1D array access. Use only if N == 1.
00868          */
00869     reference operator() (difference_type_1 x)
00870     {
00871         VIGRA_ASSERT_INSIDE(difference_type(x));
00872         return m_ptr [detail::CoordinatesToOffest<C>::exec(m_stride, x)];
00873     }
00874 
00875         /** 2D array access. Use only if N == 2.
00876          */
00877     reference operator() (difference_type_1 x, difference_type_1 y)
00878     {
00879         VIGRA_ASSERT_INSIDE(difference_type(x, y));
00880         return m_ptr [detail::CoordinatesToOffest<C>::exec(m_stride, x, y)];
00881     }
00882 
00883         /** 3D array access. Use only if N == 3.
00884          */
00885     reference operator() (difference_type_1 x, difference_type_1 y, difference_type_1 z)
00886     {
00887         VIGRA_ASSERT_INSIDE(difference_type(x, y, z));
00888         return m_ptr [m_stride[0]*x + m_stride[1]*y + m_stride[2]*z];
00889     }
00890 
00891         /** 4D array access. Use only if N == 4.
00892          */
00893     reference operator() (difference_type_1 x, difference_type_1 y,
00894                           difference_type_1 z, difference_type_1 u)
00895     {
00896         VIGRA_ASSERT_INSIDE(difference_type(x, y, z, u));
00897         return m_ptr [m_stride[0]*x + m_stride[1]*y + m_stride[2]*z + m_stride[3]*u];
00898     }
00899 
00900         /** 5D array access. Use only if N == 5.
00901          */
00902     reference operator() (difference_type_1 x, difference_type_1 y, difference_type_1 z,
00903                           difference_type_1 u, difference_type_1 v)
00904     {
00905         VIGRA_ASSERT_INSIDE(difference_type(x, y,z, u,v));
00906         return m_ptr [m_stride[0]*x + m_stride[1]*y + m_stride[2]*z + m_stride[3]*u + m_stride[4]*v];
00907     }
00908 
00909         /** 1D const array access. Use only if N == 1.
00910          */
00911     const_reference operator() (difference_type_1 x) const
00912     {
00913         VIGRA_ASSERT_INSIDE(difference_type(x));
00914         return m_ptr [detail::CoordinatesToOffest<C>::exec(m_stride, x)];
00915     }
00916 
00917         /** 2D const array access. Use only if N == 2.
00918          */
00919     const_reference operator() (difference_type_1 x, difference_type_1 y) const
00920     {
00921         VIGRA_ASSERT_INSIDE(difference_type(x, y));
00922         return m_ptr [detail::CoordinatesToOffest<C>::exec(m_stride, x, y)];
00923     }
00924 
00925         /** 3D const array access. Use only if N == 3.
00926          */
00927     const_reference operator() (difference_type_1 x, difference_type_1 y, difference_type_1 z) const
00928     {
00929         VIGRA_ASSERT_INSIDE(difference_type(x,y,z));
00930         return m_ptr [m_stride[0]*x + m_stride[1]*y + m_stride[2]*z];
00931     }
00932 
00933         /** 4D const array access. Use only if N == 4.
00934          */
00935     const_reference operator() (difference_type_1 x, difference_type_1 y,
00936                                 difference_type_1 z, difference_type_1 u) const
00937     {
00938         VIGRA_ASSERT_INSIDE(difference_type(x,y,z,u));
00939         return m_ptr [m_stride[0]*x + m_stride[1]*y + m_stride[2]*z + m_stride[3]*u];
00940     }
00941 
00942         /** 5D const array access. Use only if N == 5.
00943          */
00944     const_reference operator() (difference_type_1 x, difference_type_1 y, difference_type_1 z,
00945                                 difference_type_1 u, difference_type_1 v) const
00946     {
00947         VIGRA_ASSERT_INSIDE(difference_type(x,y,z,u,v));
00948         return m_ptr [m_stride[0]*x + m_stride[1]*y + m_stride[2]*z + m_stride[3]*u + m_stride[4]*v];
00949     }
00950 
00951         /** Init with a constant.
00952          */
00953     template <class U>
00954     MultiArrayView & init(const U & init)
00955     {
00956         detail::copyScalarMultiArrayData(traverser_begin(), shape(), init, MetaInt<actual_dimension-1>());
00957         return *this;
00958     }
00959 
00960 
00961         /** Copy the data of the right-hand array (array shapes must match).
00962          */
00963     void copy(const MultiArrayView & rhs)
00964     {
00965         if(this == &rhs)
00966             return;
00967         this->copyImpl(rhs);
00968     }
00969 
00970         /** Copy the data of the right-hand array (array shapes must match).
00971          */
00972     template <class U, class CN>
00973     void copy(const MultiArrayView <N, U, CN>& rhs)
00974     {
00975         this->copyImpl(rhs);
00976     }
00977 
00978         /** swap the data between two MultiArrayView objects.
00979 
00980             The shapes of the two array must match.
00981         */
00982     void swapData(MultiArrayView rhs)
00983     {
00984         if(this != &rhs)
00985             swapDataImpl(rhs);
00986     }
00987 
00988         /** swap the data between two MultiArrayView objects.
00989 
00990             The shapes of the two array must match.
00991         */
00992     template <class T2, class C2>
00993     void swapData(MultiArrayView <N, T2, C2> rhs)
00994     {
00995         swapDataImpl(rhs);
00996     }
00997 
00998         /** bind the M outmost dimensions to certain indices.
00999             this reduces the dimensionality of the image to
01000             max { 1, N-M }.
01001 
01002             <b>Usage:</b>
01003             \code
01004             // create a 3D array of size 40x30x20
01005             typedef MultiArray<3, double>::difference_type Shape;
01006             MultiArray<3, double> array3(Shape(40, 30, 20));
01007 
01008             // get a 1D array by fixing index 1 to 12, and index 2 to 10
01009             MultiArrayView <1, double> array1 = array3.bindOuter(TinyVector<int, 2>(12, 10));
01010             \endcode
01011         */
01012     template <unsigned int M>
01013     MultiArrayView <N-M, T, C> bindOuter (const TinyVector <MultiArrayIndex, M> &d) const;
01014 
01015         /** bind the M innermost dimensions to certain indices.
01016             this reduces the dimensionality of the image to
01017             max { 1, N-M }.
01018 
01019             <b>Usage:</b>
01020             \code
01021             // create a 3D array of size 40x30x20
01022             typedef MultiArray<3, double>::difference_type Shape;
01023             MultiArray<3, double> array3(Shape(40, 30, 20));
01024 
01025             // get a 1D array by fixing index 0 to 12, and index 1 to 10
01026             MultiArrayView <1, double, StridedArrayTag> array1 = array3.bindInner(TinyVector<int, 2>(12, 10));
01027             \endcode
01028         */
01029     template <unsigned int M>
01030     MultiArrayView <N-M, T, StridedArrayTag>
01031     bindInner (const TinyVector <MultiArrayIndex, M> &d) const;
01032 
01033         /** bind dimension M to index d.
01034             this reduces the dimensionality of the image to
01035             max { 1, N-1 }.
01036 
01037             <b>Usage:</b>
01038             \code
01039             // create a 3D array of size 40x30x20
01040             typedef MultiArray<3, double>::difference_type Shape;
01041             MultiArray<3, double> array3(Shape(40, 30, 20));
01042 
01043             // get a 2D array by fixing index 1 to 12
01044             MultiArrayView <2, double> array2 = array3.bind<1>(12);
01045 
01046             // get a 2D array by fixing index 0 to 23
01047             MultiArrayView <2, double, StridedArrayTag> array2a = array3.bind<0>(23);
01048             \endcode
01049          */
01050     template <unsigned int M>
01051     MultiArrayView <N-1, T, typename vigra::detail::MaybeStrided <M>::type >
01052     bind (difference_type_1 d) const;
01053 
01054         /** bind the outmost dimension to a certain index.
01055             this reduces the dimensionality of the image to
01056             max { 1, N-1 }.
01057 
01058             <b>Usage:</b>
01059             \code
01060             // create a 3D array of size 40x30x20
01061             typedef MultiArray<3, double>::difference_type Shape;
01062             MultiArray<3, double> array3(Shape(40, 30, 20));
01063 
01064             // get a 2D array by fixing the outermost index (i.e. index 2) to 12
01065             MultiArrayView <2, double> array2 = array3.bindOuter(12);
01066             \endcode
01067         */
01068     MultiArrayView <N-1, T, C> bindOuter (difference_type_1 d) const;
01069 
01070         /** bind the innermost dimension to a certain index.
01071             this reduces the dimensionality of the image to
01072             max { 1, N-1 }.
01073 
01074             <b>Usage:</b>
01075             \code
01076             // create a 3D array of size 40x30x20
01077             typedef MultiArray<3, double>::difference_type Shape;
01078             MultiArray<3, double> array3(Shape(40, 30, 20));
01079 
01080             // get a 2D array by fixing the innermost index (i.e. index 0) to 23
01081             MultiArrayView <2, double, StridedArrayTag> array2 = array3.bindInner(23);
01082             \endcode
01083         */
01084     MultiArrayView <N-1, T, StridedArrayTag> bindInner (difference_type_1 d) const;
01085 
01086         /** bind dimension m to index d.
01087             this reduces the dimensionality of the image to
01088             max { 1, N-1 }.
01089 
01090             <b>Usage:</b>
01091             \code
01092             // create a 3D array of size 40x30x20
01093             typedef MultiArray<3, double>::difference_type Shape;
01094             MultiArray<3, double> array3(Shape(40, 30, 20));
01095 
01096             // get a 2D array by fixing index 2 to 15
01097             MultiArrayView <2, double, StridedArrayTag> array2 = array3.bindAt(2, 15);
01098             \endcode
01099          */
01100     MultiArrayView <N-1, T, StridedArrayTag>
01101     bindAt (difference_type_1 m, difference_type_1 d) const;
01102 
01103         /** Add a singleton dimension (dimension of legth 1).
01104 
01105             Singleton dimensions don't change the size of the data, but introduce
01106             a new index that can only take the value 0. This is mainly useful for
01107             the 'reduce mode' of transformMultiArray() and combineTwoMultiArrays(),
01108             because these functions require the source and destination arrays to
01109             have the same number of dimensions.
01110 
01111             The range of \a i must be <tt>0 <= i <= N</tt>. The new dimension will become
01112             the i'th index, and the old indices from i upwards will shift one
01113             place to the right.
01114 
01115             <b>Usage:</b>
01116 
01117             Suppose we want have a 2D array and want to create a 1D array that contains
01118             the row average of the first array.
01119             \code
01120             typedef MultiArrayShape<2>::type Shape2;
01121             MultiArray<2, double> original(Shape2(40, 30));
01122 
01123             typedef MultiArrayShape<1>::type Shape1;
01124             MultiArray<1, double> rowAverages(Shape1(30));
01125 
01126             // temporarily add a singleton dimension to the destination array
01127             transformMultiArray(srcMultiArrayRange(original),
01128                                 destMultiArrayRange(rowAverages.insertSingletonDimension(0)),
01129                                 FindAverage<double>());
01130             \endcode
01131          */
01132     MultiArrayView <N+1, T, C>
01133     insertSingletonDimension (difference_type_1 i) const;
01134 
01135         /** create a rectangular subarray that spans between the
01136             points p and q, where p is in the subarray, q not.
01137 
01138             <b>Usage:</b>
01139             \code
01140             // create a 3D array of size 40x30x20
01141             typedef MultiArray<3, double>::difference_type Shape;
01142             MultiArray<3, double> array3(Shape(40, 30, 20));
01143 
01144             // get a subarray set is smaller by one element at all sides
01145             MultiArrayView <3, double> subarray = array3.subarray(Shape(1,1,1), Shape(39, 29, 19));
01146             \endcode
01147         */
01148     MultiArrayView subarray (const difference_type &p,
01149                              const difference_type &q) const
01150     {
01151         const difference_type_1 offset = dot (m_stride, p);
01152         return MultiArrayView (q - p, m_stride, m_ptr + offset);
01153     }
01154 
01155         /** apply an additional striding to the image, thereby reducing
01156             the shape of the array.
01157             for example, multiplying the stride of dimension one by three
01158             turns an appropriately layed out (interleaved) rgb image into
01159             a single band image.
01160         */
01161     MultiArrayView <N, T, StridedArrayTag>
01162     stridearray (const difference_type &s) const
01163     {
01164         difference_type shape = m_shape;
01165         for (unsigned int i = 0; i < actual_dimension; ++i)
01166             shape [i] /= s [i];
01167         return MultiArrayView <N, T, StridedArrayTag>(shape, m_stride * s, m_ptr);
01168     }
01169 
01170         /** Transpose an array. If N==2, this implements the usual matrix transposition.
01171             For N > 2, it reverses the order of the indices.
01172 
01173             <b>Usage:</b><br>
01174             \code
01175             typedef MultiArray<2, double>::difference_type Shape;
01176             MultiArray<2, double> array(10, 20);
01177 
01178             MultiArray<2, double, StridedArrayTag> transposed = array.transpose();
01179 
01180             for(int i=0; i<array.shape(0), ++i)
01181                 for(int j=0; j<array.shape(1); ++j)
01182                     assert(array(i, j) == transposed(j, i));
01183             \endcode
01184         */
01185     MultiArrayView <N, T, StridedArrayTag>
01186     transpose () const
01187     {
01188         difference_type shape(m_shape.begin(), difference_type::ReverseCopy),
01189                         stride(m_stride.begin(), difference_type::ReverseCopy);
01190         return MultiArrayView <N, T, StridedArrayTag>(shape, stride, m_ptr);
01191     }
01192 
01193         /** permute the dimensions of the array.
01194             The function exchanges the meaning of the dimensions without copying the data.
01195             In case of a 2-dimensional array, this is simply array transposition. In higher dimensions,
01196             there are more possibilities.
01197 
01198             <b>Usage:</b><br>
01199             \code
01200             typedef MultiArray<2, double>::difference_type Shape;
01201             MultiArray<2, double> array(10, 20);
01202 
01203             MultiArray<2, double, StridedArrayTag> transposed = array.permuteDimensions(Shape(1,0));
01204 
01205             for(int i=0; i<array.shape(0), ++i)
01206                 for(int j=0; j<array.shape(1); ++j)
01207                     assert(array(i, j) == transposed(j, i));
01208             \endcode
01209         */
01210     MultiArrayView <N, T, StridedArrayTag>
01211     permuteDimensions (const difference_type &s) const;
01212 
01213         /** Permute the dimensions of the array so that the strides are in ascending order.
01214             Determines the appropriate permutation and then calls permuteDimensions().
01215         */
01216     MultiArrayView <N, T, StridedArrayTag>
01217     permuteStridesAscending() const;
01218     
01219         /** Permute the dimensions of the array so that the strides are in descending order.
01220             Determines the appropriate permutation and then calls permuteDimensions().
01221         */
01222     MultiArrayView <N, T, StridedArrayTag>
01223     permuteStridesDescending() const;
01224     
01225         /** Compute the ordering of the strides in this array.
01226             The result is describes the current permutation of the axes relative 
01227             to the standard ascending stride order.
01228         */
01229     difference_type strideOrdering() const
01230     {
01231         return strideOrdering(m_stride);
01232     }
01233     
01234         /** Compute the ordering of the given strides.
01235             The result is describes the current permutation of the axes relative 
01236             to the standard ascending stride order.
01237         */
01238     static difference_type strideOrdering(difference_type strides);
01239 
01240         /** number of the elements in the array.
01241          */
01242     difference_type_1 elementCount () const
01243     {
01244         difference_type_1 ret = m_shape[0];
01245         for(int i = 1; i < actual_dimension; ++i)
01246             ret *= m_shape[i];
01247         return ret;
01248     }
01249 
01250         /** number of the elements in the array.
01251             Same as <tt>elementCount()</tt>. Mostly useful to support the std::vector interface.
01252          */
01253     difference_type_1 size () const
01254     {
01255         return elementCount();
01256     }
01257 
01258         /** return the array's shape.
01259          */
01260     const difference_type & shape () const
01261     {
01262         return m_shape;
01263     }
01264 
01265         /** return the array's size at a certain dimension.
01266          */
01267     difference_type_1 size (difference_type_1 n) const
01268     {
01269         return m_shape [n];
01270     }
01271 
01272         /** return the array's shape at a certain dimension
01273             (same as <tt>size(n)</tt>).
01274          */
01275     difference_type_1 shape (difference_type_1 n) const
01276     {
01277         return m_shape [n];
01278     }
01279 
01280         /** return the array's stride for every dimension.
01281          */
01282     const difference_type & stride () const
01283     {
01284         return m_stride;
01285     }
01286 
01287         /** return the array's stride at a certain dimension.
01288          */
01289     difference_type_1 stride (int n) const
01290     {
01291         return m_stride [n];
01292     }
01293 
01294         /** check whether two arrays are elementwise equal.
01295          */
01296     template <class U, class C1>
01297     bool operator==(MultiArrayView<N, U, C1> const & rhs) const
01298     {
01299         if(this->shape() != rhs.shape())
01300             return false;
01301         return detail::equalityOfMultiArrays(traverser_begin(), shape(), rhs.traverser_begin(), MetaInt<actual_dimension-1>());
01302     }
01303 
01304         /** check whether two arrays are not elementwise equal.
01305             Also true when the two arrays have different shapes.
01306          */
01307     template <class U, class C1>
01308     bool operator!=(MultiArrayView<N, U, C1> const & rhs) const
01309     {
01310         return !operator==(rhs);
01311     }
01312 
01313         /** check whether the given point is in the array range.
01314          */
01315     bool isInside (difference_type const & p) const
01316     {
01317         for(int d=0; d<actual_dimension; ++d)
01318             if(p[d] < 0 || p[d] >= shape(d))
01319                 return false;
01320         return true;
01321     }
01322 
01323         /** Compute the squared Euclidean norm of the array (sum of squares of the array elements).
01324          */
01325     typename NormTraits<MultiArrayView>::SquaredNormType squaredNorm() const
01326     {
01327         typedef typename NormTraits<MultiArrayView>::SquaredNormType SquaredNormType;
01328         SquaredNormType res = NumericTraits<SquaredNormType>::zero();
01329         detail::sumOverMultiArray(traverser_begin(), shape(), detail::MultiArrayL2Functor<SquaredNormType>(),
01330                                   res, MetaInt<actual_dimension-1>());
01331         return res;
01332     }
01333 
01334         /** Compute various norms of the array.
01335             The norm is determined by parameter \a type:
01336 
01337             <ul>
01338             <li> type == 0: maximum norm (L-infinity): maximum of absolute values of the array elements
01339             <li> type == 1: Manhattan norm (L1): sum of absolute values of the array elements
01340             <li> type == 2: Euclidean norm (L2): square root of <tt>squaredNorm()</tt> when \a useSquaredNorm is <tt>true</tt>,<br>
01341                  or direct algorithm that avoids underflow/overflow otherwise.
01342             </ul>
01343 
01344             Parameter \a useSquaredNorm has no effect when \a type != 2. Defaults: compute L2 norm as square root of
01345             <tt>squaredNorm()</tt>.
01346          */
01347     typename NormTraits<MultiArrayView>::NormType norm(int type = 2, bool useSquaredNorm = true) const;
01348 
01349         /** return the pointer to the image data
01350          */
01351     pointer data () const
01352     {
01353         return m_ptr;
01354     }
01355 
01356         /**
01357          * returns true iff this view refers to valid data,
01358          * i.e. data() is not a NULL pointer.  (this is false after
01359          * default construction.)
01360          */
01361     bool hasData () const
01362     {
01363         return m_ptr != 0;
01364     }
01365 
01366         /** returns the N-dimensional MultiIterator pointing
01367             to the first element in every dimension.
01368         */
01369     traverser traverser_begin ()
01370     {
01371         traverser ret (m_ptr, m_stride.begin (), m_shape.begin ());
01372         return ret;
01373     }
01374 
01375         /** returns the N-dimensional MultiIterator pointing
01376             to the const first element in every dimension.
01377         */
01378     const_traverser traverser_begin () const
01379     {
01380         const_traverser ret (m_ptr, m_stride.begin (), m_shape.begin ());
01381         return ret;
01382     }
01383 
01384         /** returns the N-dimensional MultiIterator pointing
01385             beyond the last element in dimension N, and to the
01386             first element in every other dimension.
01387         */
01388     traverser traverser_end ()
01389     {
01390         traverser ret (m_ptr, m_stride.begin (), m_shape.begin ());
01391         ret += m_shape [actual_dimension-1];
01392         return ret;
01393     }
01394 
01395         /** returns the N-dimensional const MultiIterator pointing
01396             beyond the last element in dimension N, and to the
01397             first element in every other dimension.
01398         */
01399     const_traverser traverser_end () const
01400     {
01401         const_traverser ret (m_ptr, m_stride.begin (), m_shape.begin ());
01402         ret += m_shape [actual_dimension-1];
01403         return ret;
01404     }
01405 
01406     view_type view ()
01407     {
01408         return *this;
01409     }
01410 };
01411 
01412 template <unsigned int N, class T, class C>
01413 MultiArrayView<N, T, C> &
01414 MultiArrayView <N, T, C>::operator=(MultiArrayView<N, T, C> const & rhs)
01415 {
01416     if(this == &rhs)
01417         return *this;
01418     vigra_precondition(this->shape() == rhs.shape() || m_ptr == 0,
01419         "MultiArrayView::operator=(MultiArrayView const &) size mismatch - use MultiArrayView::reset().");
01420     if(m_ptr == 0)
01421     {
01422         m_shape  = rhs.m_shape;
01423         m_stride = rhs.m_stride;
01424         m_ptr    = rhs.m_ptr;
01425     }
01426     else
01427         this->copyImpl(rhs);
01428     return *this;
01429 }
01430 
01431 template <unsigned int N, class T, class C>
01432 template <class CN>
01433 bool
01434 MultiArrayView <N, T, C>::arraysOverlap(const MultiArrayView <N, T, CN>& rhs) const
01435 {
01436     vigra_precondition (shape () == rhs.shape (),
01437         "MultiArrayView::arraysOverlap(): shape mismatch.");
01438     const_pointer first_element = this->m_ptr,
01439                   last_element = first_element + dot(this->m_shape - difference_type(1), this->m_stride);
01440     typename MultiArrayView <N, T, CN>::const_pointer
01441            rhs_first_element = rhs.data(),
01442            rhs_last_element = rhs_first_element + dot(rhs.shape() - difference_type(1), rhs.stride());
01443     return !(last_element < rhs_first_element || rhs_last_element < first_element);
01444 }
01445 
01446 template <unsigned int N, class T, class C>
01447 template <class U, class CN>
01448 void
01449 MultiArrayView <N, T, C>::copyImpl(const MultiArrayView <N, U, CN>& rhs)
01450 {
01451     if(!arraysOverlap(rhs))
01452     {
01453         // no overlap -- can copy directly
01454         detail::copyMultiArrayData(rhs.traverser_begin(), shape(), traverser_begin(), MetaInt<actual_dimension-1>());
01455     }
01456     else
01457     {
01458         // overlap: we got different views to the same data -- copy to intermediate memory in order to avoid
01459         // overwriting elements that are still needed on the rhs.
01460         MultiArray<N, T> tmp(rhs);
01461         detail::copyMultiArrayData(tmp.traverser_begin(), shape(), traverser_begin(), MetaInt<actual_dimension-1>());
01462     }
01463 }
01464 
01465 #define VIGRA_MULTI_ARRY_COMPUTED_ASSIGNMENT(name, op) \
01466 template <unsigned int N, class T, class C> \
01467 template<class U, class C1> \
01468 MultiArrayView<N, T, C> &  \
01469 MultiArrayView <N, T, C>::operator op(MultiArrayView<N, U, C1> const & rhs) \
01470 { \
01471     vigra_precondition(this->shape() == rhs.shape(), "MultiArrayView::operator" #op "() size mismatch."); \
01472     if(!arraysOverlap(rhs)) \
01473     { \
01474         detail::name##MultiArrayData(rhs.traverser_begin(), shape(), traverser_begin(), MetaInt<actual_dimension-1>()); \
01475     } \
01476     else \
01477     { \
01478         MultiArray<N, T> tmp(rhs); \
01479         detail::name##MultiArrayData(tmp.traverser_begin(), shape(), traverser_begin(), MetaInt<actual_dimension-1>()); \
01480     } \
01481     return *this; \
01482 }
01483 
01484 VIGRA_MULTI_ARRY_COMPUTED_ASSIGNMENT(copyAdd, +=)
01485 VIGRA_MULTI_ARRY_COMPUTED_ASSIGNMENT(copySub, -=)
01486 VIGRA_MULTI_ARRY_COMPUTED_ASSIGNMENT(copyMul, *=)
01487 VIGRA_MULTI_ARRY_COMPUTED_ASSIGNMENT(copyDiv, /=)
01488 
01489 #undef VIGRA_MULTI_ARRY_COMPUTED_ASSIGNMENT
01490 
01491 template <unsigned int N, class T, class C>
01492 template <class U, class CN>
01493 void
01494 MultiArrayView <N, T, C>::swapDataImpl(MultiArrayView <N, U, CN> rhs)
01495 {
01496     vigra_precondition (shape () == rhs.shape (),
01497         "MultiArrayView::swapData(): shape mismatch.");
01498 
01499     // check for overlap of this and rhs
01500     const_pointer first_element = this->m_ptr,
01501                   last_element = first_element + dot(this->m_shape - difference_type(1), this->m_stride);
01502     typename MultiArrayView <N, U, CN>::const_pointer
01503            rhs_first_element = rhs.data(),
01504            rhs_last_element = rhs_first_element + dot(rhs.shape() - difference_type(1), rhs.stride());
01505     if(last_element < rhs_first_element || rhs_last_element < first_element)
01506     {
01507         // no overlap -- can swap directly
01508         detail::swapDataImpl(traverser_begin(), shape(), rhs.traverser_begin(), MetaInt<actual_dimension-1>());
01509     }
01510     else
01511     {
01512         // overlap: we got different views to the same data -- copy to intermediate memory in order to avoid
01513         // overwriting elements that are still needed.
01514         MultiArray<N, T> tmp(*this);
01515         copy(rhs);
01516         rhs.copy(tmp);
01517     }
01518 }
01519 
01520 template <unsigned int N, class T, class C>
01521 MultiArrayView <N, T, StridedArrayTag>
01522 MultiArrayView <N, T, C>::permuteDimensions (const difference_type &s) const
01523 {
01524     difference_type shape, stride, check((typename difference_type::value_type)0);
01525     for (unsigned int i = 0; i < actual_dimension; ++i)
01526     {
01527         shape[i]  = m_shape[s[i]];
01528         stride[i] = m_stride[s[i]];
01529         ++check[s[i]];
01530     }
01531     vigra_precondition(check == difference_type(1),
01532        "MultiArrayView::permuteDimensions(): every dimension must occur exactly once.");
01533     return MultiArrayView <N, T, StridedArrayTag>(shape, stride, m_ptr);
01534 }
01535 
01536 template <unsigned int N, class T, class C>
01537 typename MultiArrayView <N, T, C>::difference_type 
01538 MultiArrayView <N, T, C>::strideOrdering(difference_type stride)
01539 {
01540     difference_type permutation;
01541     for(unsigned int k=0; k<N; ++k)
01542         permutation[k] = k;
01543     for(unsigned int k=0; k<N-1; ++k)
01544     {
01545         unsigned int smallest = k;
01546         for(unsigned int j=k+1; j<N; ++j)
01547         {
01548             if(stride[j] < stride[smallest])
01549                 smallest = j;
01550         }
01551         if(smallest != k)
01552         {
01553             std::swap(stride[k], stride[smallest]);
01554             std::swap(permutation[k], permutation[smallest]);
01555         }
01556     }
01557     difference_type ordering;
01558     for(unsigned int k=0; k<N; ++k)
01559         ordering[permutation[k]] = k;
01560     return ordering;
01561 }
01562 
01563 template <unsigned int N, class T, class C>
01564 MultiArrayView <N, T, StridedArrayTag>
01565 MultiArrayView <N, T, C>::permuteStridesAscending() const
01566 {
01567     difference_type ordering(strideOrdering(m_stride)), permutation;
01568     for(MultiArrayIndex k=0; k<N; ++k)
01569         permutation[ordering[k]] = k;
01570     return permuteDimensions(permutation);
01571 }
01572 
01573 template <unsigned int N, class T, class C>
01574 MultiArrayView <N, T, StridedArrayTag>
01575 MultiArrayView <N, T, C>::permuteStridesDescending() const
01576 {
01577     difference_type ordering(strideOrdering(m_stride)), permutation;
01578     for(MultiArrayIndex k=0; k<N; ++k)
01579         permutation[ordering[N-1-k]] = k;
01580     return permuteDimensions(permutation);
01581 }
01582 
01583 template <unsigned int N, class T, class C>
01584 template <unsigned int M>
01585 MultiArrayView <N-M, T, C>
01586 MultiArrayView <N, T, C>::bindOuter (const TinyVector <MultiArrayIndex, M> &d) const
01587 {
01588     TinyVector <MultiArrayIndex, M> stride;
01589     stride.init (m_stride.begin () + N-M, m_stride.end ());
01590     pointer ptr = m_ptr + dot (d, stride);
01591     static const int NNew = (N-M == 0) ? 1 : N-M;
01592     TinyVector <MultiArrayIndex, NNew> inner_shape, inner_stride;
01593     if (N-M == 0)
01594     {
01595         inner_shape [0] = 1;
01596         inner_stride [0] = 0;
01597     }
01598     else
01599     {
01600         inner_shape.init (m_shape.begin (), m_shape.end () - M);
01601         inner_stride.init (m_stride.begin (), m_stride.end () - M);
01602     }
01603     return MultiArrayView <N-M, T, C> (inner_shape, inner_stride, ptr);
01604 }
01605 
01606 template <unsigned int N, class T, class C>
01607 template <unsigned int M>
01608 MultiArrayView <N - M, T, StridedArrayTag>
01609 MultiArrayView <N, T, C>::bindInner (const TinyVector <MultiArrayIndex, M> &d) const
01610 {
01611     TinyVector <MultiArrayIndex, M> stride;
01612     stride.init (m_stride.begin (), m_stride.end () - N + M);
01613     pointer ptr = m_ptr + dot (d, stride);
01614     static const int NNew = (N-M == 0) ? 1 : N-M;
01615     TinyVector <MultiArrayIndex, NNew> outer_shape, outer_stride;
01616     if (N-M == 0)
01617     {
01618         outer_shape [0] = 1;
01619         outer_stride [0] = 0;
01620     }
01621     else
01622     {
01623         outer_shape.init (m_shape.begin () + M, m_shape.end ());
01624         outer_stride.init (m_stride.begin () + M, m_stride.end ());
01625     }
01626     return MultiArrayView <N-M, T, StridedArrayTag>
01627         (outer_shape, outer_stride, ptr);
01628 }
01629 
01630 template <unsigned int N, class T, class C>
01631 template <unsigned int M>
01632 MultiArrayView <N-1, T, typename detail::MaybeStrided <M>::type >
01633 MultiArrayView <N, T, C>::bind (difference_type_1 d) const
01634 {
01635     static const int NNew = (N-1 == 0) ? 1 : N-1;
01636     TinyVector <MultiArrayIndex, NNew> shape, stride;
01637     // the remaining dimensions are 0..n-1,n+1..N-1
01638     if (N-1 == 0)
01639     {
01640         shape[0] = 1;
01641         stride[0] = 0;
01642     }
01643     else
01644     {
01645         std::copy (m_shape.begin (), m_shape.begin () + M, shape.begin ());
01646         std::copy (m_shape.begin () + M+1, m_shape.end (),
01647                    shape.begin () + M);
01648         std::copy (m_stride.begin (), m_stride.begin () + M, stride.begin ());
01649         std::copy (m_stride.begin () + M+1, m_stride.end (),
01650                    stride.begin () + M);
01651     }
01652     return MultiArrayView <N-1, T, typename detail::MaybeStrided <M>::type>
01653         (shape, stride, m_ptr + d * m_stride[M]);
01654 }
01655 
01656 template <unsigned int N, class T, class C>
01657 MultiArrayView <N - 1, T, C>
01658 MultiArrayView <N, T, C>::bindOuter (difference_type_1 d) const
01659 {
01660     static const int NNew = (N-1 == 0) ? 1 : N-1;
01661     TinyVector <MultiArrayIndex, NNew> inner_shape, inner_stride;
01662     if (N-1 == 0)
01663     {
01664         inner_shape [0] = 1;
01665         inner_stride [0] = 0;
01666     }
01667     else
01668     {
01669         inner_shape.init (m_shape.begin (), m_shape.end () - 1);
01670         inner_stride.init (m_stride.begin (), m_stride.end () - 1);
01671     }
01672     return MultiArrayView <N-1, T, C> (inner_shape, inner_stride,
01673                                        m_ptr + d * m_stride [N-1]);
01674 }
01675 
01676 template <unsigned int N, class T, class C>
01677 MultiArrayView <N - 1, T, StridedArrayTag>
01678 MultiArrayView <N, T, C>::bindInner (difference_type_1 d) const
01679 {
01680     static const int NNew = (N-1 == 0) ? 1 : N-1;
01681     TinyVector <MultiArrayIndex, NNew> outer_shape, outer_stride;
01682     if (N-1 == 0)
01683     {
01684         outer_shape [0] = 1;
01685         outer_stride [0] = 0;
01686     }
01687     else
01688     {
01689         outer_shape.init (m_shape.begin () + 1, m_shape.end ());
01690         outer_stride.init (m_stride.begin () + 1, m_stride.end ());
01691     }
01692     return MultiArrayView <N-1, T, StridedArrayTag>
01693         (outer_shape, outer_stride, m_ptr + d * m_stride [0]);
01694 }
01695 
01696 template <unsigned int N, class T, class C>
01697 MultiArrayView <N - 1, T, StridedArrayTag>
01698 MultiArrayView <N, T, C>::bindAt (difference_type_1 n, difference_type_1 d) const
01699 {
01700     vigra_precondition (
01701         n < static_cast <int> (N),
01702         "MultiArrayView <N, T, C>::bindAt(): dimension out of range.");
01703     static const int NNew = (N-1 == 0) ? 1 : N-1;
01704     TinyVector <MultiArrayIndex, NNew> shape, stride;
01705     // the remaining dimensions are 0..n-1,n+1..N-1
01706     if (N-1 == 0)
01707     {
01708         shape [0] = 1;
01709         stride [0] = 0;
01710     }
01711     else
01712     {
01713         std::copy (m_shape.begin (), m_shape.begin () + n, shape.begin ());
01714         std::copy (m_shape.begin () + n+1, m_shape.end (),
01715                    shape.begin () + n);
01716         std::copy (m_stride.begin (), m_stride.begin () + n, stride.begin ());
01717         std::copy (m_stride.begin () + n+1, m_stride.end (),
01718                    stride.begin () + n);
01719     }
01720     return MultiArrayView <N-1, T, StridedArrayTag>
01721         (shape, stride, m_ptr + d * m_stride[n]);
01722 }
01723 
01724 template <unsigned int N, class T, class C>
01725 MultiArrayView <N + 1, T, C>
01726 MultiArrayView <N, T, C>::insertSingletonDimension (difference_type_1 i) const
01727 {
01728     vigra_precondition (
01729         0 <= i && i <= static_cast <difference_type_1> (N),
01730         "MultiArrayView <N, T, C>::insertSingletonDimension(): index out of range.");
01731     TinyVector <MultiArrayIndex, N+1> shape, stride;
01732     std::copy (m_shape.begin (), m_shape.begin () + i, shape.begin ());
01733     std::copy (m_shape.begin () + i, m_shape.end (), shape.begin () + i + 1);
01734     std::copy (m_stride.begin (), m_stride.begin () + i, stride.begin ());
01735     std::copy (m_stride.begin () + i, m_stride.end (), stride.begin () + i + 1);
01736     shape[i] = 1;
01737     stride[i] = 1;
01738 
01739     return MultiArrayView <N+1, T, C>(shape, stride, m_ptr);
01740 }
01741 
01742 template <unsigned int N, class T, class C>
01743 typename NormTraits<MultiArrayView <N, T, C> >::NormType
01744 MultiArrayView <N, T, C>::norm(int type, bool useSquaredNorm) const
01745 {
01746     typedef typename NormTraits<MultiArrayView>::NormType NormType;
01747 
01748     switch(type)
01749     {
01750       case 0:
01751       {
01752         NormType res = NumericTraits<NormType>::zero();
01753         detail::normMaxOfMultiArray(traverser_begin(), shape(), res, MetaInt<actual_dimension-1>());
01754         return res;
01755       }
01756       case 1:
01757       {
01758         NormType res = NumericTraits<NormType>::zero();
01759         detail::sumOverMultiArray(traverser_begin(), shape(), detail::MultiArrayL1Functor<NormType>(),
01760                                 res, MetaInt<actual_dimension-1>());
01761         return res;
01762       }
01763       case 2:
01764       {
01765         if(useSquaredNorm)
01766         {
01767             return sqrt((NormType)squaredNorm());
01768         }
01769         else
01770         {
01771             NormType normMax = NumericTraits<NormType>::zero();
01772             detail::normMaxOfMultiArray(traverser_begin(), shape(), normMax, MetaInt<actual_dimension-1>());
01773             if(normMax == NumericTraits<NormType>::zero())
01774                 return normMax;
01775             NormType res  = NumericTraits<NormType>::zero();
01776             detail::sumOverMultiArray(traverser_begin(), shape(), detail::MultiArrayScaledL2Functor<NormType>(normMax),
01777                                     res, MetaInt<actual_dimension-1>());
01778             return sqrt(res)*normMax;
01779         }
01780       }
01781       default:
01782         vigra_precondition(false, "MultiArrayView::norm(): Unknown norm type.");
01783         return NumericTraits<NormType>::zero(); // unreachable
01784     }
01785 }
01786 
01787 
01788 /********************************************************/
01789 /*                                                      */
01790 /*                          norm                        */
01791 /*                                                      */
01792 /********************************************************/
01793 
01794 template <unsigned int N, class T, class C>
01795 inline typename NormTraits<MultiArrayView <N, T, C> >::SquaredNormType
01796 squaredNorm(MultiArrayView <N, T, C> const & a)
01797 {
01798     return a.squaredNorm();
01799 }
01800 
01801 template <unsigned int N, class T, class C>
01802 inline typename NormTraits<MultiArrayView <N, T, C> >::NormType
01803 norm(MultiArrayView <N, T, C> const & a)
01804 {
01805     return a.norm();
01806 }
01807 
01808 /********************************************************/
01809 /*                                                      */
01810 /*                       MultiArray                     */
01811 /*                                                      */
01812 /********************************************************/
01813 
01814 /** \brief Main <TT>MultiArray</TT> class containing the memory
01815     management.
01816 
01817 This class inherits the interface of MultiArrayView, and implements
01818 the memory ownership.
01819 MultiArray's are always unstrided, striding them creates a MultiArrayView.
01820 
01821 
01822 The template parameters are as follows
01823 \code
01824     N: the array dimension
01825 
01826     T: the type of the array elements
01827 
01828     A: the allocator used for internal storage management
01829        (default: std::allocator<T>)
01830 \endcode
01831 
01832 <b>\#include</b>
01833 <<a href="multi__array_8hxx-source.html">vigra/multi_array.hxx</a>>
01834 
01835 Namespace: vigra
01836 */
01837 template <unsigned int N, class T, class A /* default already declared above */>
01838 class MultiArray : public MultiArrayView <N, T>
01839 {
01840 
01841 public:
01842     using MultiArrayView <N, T>::actual_dimension;
01843 
01844         /** the allocator type used to allocate the memory
01845          */
01846     typedef A allocator_type;
01847 
01848         /** the view type associated with this array.
01849          */
01850     typedef MultiArrayView <N, T> view_type;
01851 
01852         /** the matrix type associated with this array.
01853          */
01854     typedef MultiArray <N, T, A> matrix_type;
01855 
01856         /** the array's value type
01857          */
01858     typedef typename view_type::value_type value_type;
01859 
01860         /** pointer type
01861          */
01862     typedef typename view_type::pointer pointer;
01863 
01864         /** const pointer type
01865          */
01866     typedef typename view_type::const_pointer const_pointer;
01867 
01868         /** reference type (result of operator[])
01869          */
01870     typedef typename view_type::reference reference;
01871 
01872         /** const reference type (result of operator[] const)
01873          */
01874     typedef typename view_type::const_reference const_reference;
01875 
01876         /** size type
01877          */
01878     typedef typename view_type::size_type size_type;
01879 
01880         /** difference type (used for multi-dimensional offsets and indices)
01881          */
01882     typedef typename view_type::difference_type difference_type;
01883 
01884         /** difference and index type for a single dimension
01885          */
01886     typedef typename view_type::difference_type_1 difference_type_1;
01887 
01888         /** traverser type
01889          */
01890     typedef typename vigra::detail::MultiIteratorChooser <
01891         UnstridedArrayTag>::template Traverser <N, T, T &, T *>::type
01892     traverser;
01893 
01894         /** traverser type to const data
01895          */
01896     typedef typename vigra::detail::MultiIteratorChooser <
01897         UnstridedArrayTag>::template Traverser <N, T, T const &, T const *>::type
01898     const_traverser;
01899 
01900         /** sequential (random access) iterator type
01901          */
01902     typedef T * iterator;
01903 
01904         /** sequential (random access) const iterator type
01905          */
01906     typedef T * const_iterator;
01907 
01908 protected:
01909 
01910     typedef typename difference_type::value_type diff_zero_t;
01911 
01912         /** the allocator used to allocate the memory
01913          */
01914     allocator_type m_alloc;
01915 
01916         /** allocate memory for s pixels, write its address into the given
01917             pointer and initialize the pixels with init.
01918         */
01919     void allocate (pointer &ptr, difference_type_1 s, const_reference init);
01920 
01921         /** allocate memory for s pixels, write its address into the given
01922             pointer and initialize the linearized pixels to the values of init.
01923         */
01924     template <class U>
01925     void allocate (pointer &ptr, difference_type_1 s, U const * init);
01926 
01927         /** allocate memory, write its address into the given
01928             pointer and initialize it by copying the data from the given MultiArrayView.
01929         */
01930     template <class U, class C>
01931     void allocate (pointer &ptr, MultiArrayView<N, U, C> const & init);
01932 
01933         /** deallocate the memory (of length s) starting at the given address.
01934          */
01935     void deallocate (pointer &ptr, difference_type_1 s);
01936 
01937     template <class U, class C>
01938     void copyOrReshape (const MultiArrayView<N, U, C> &rhs);
01939 public:
01940         /** default constructor
01941          */
01942     MultiArray ()
01943     : MultiArrayView <N, T> (difference_type (diff_zero_t(0)),
01944                              difference_type (diff_zero_t(0)), 0)
01945     {}
01946 
01947         /** construct with given allocator
01948          */
01949     MultiArray (allocator_type const & alloc)
01950     : MultiArrayView <N, T> (difference_type (diff_zero_t(0)),
01951                              difference_type (diff_zero_t(0)), 0),
01952       m_alloc(alloc)
01953     {}
01954 
01955         /** construct with given shape
01956          */
01957     explicit MultiArray (const difference_type &shape,
01958                          allocator_type const & alloc = allocator_type());
01959 
01960         /** construct from shape with an initial value
01961          */
01962     MultiArray (const difference_type &shape, const_reference init,
01963                 allocator_type const & alloc = allocator_type());
01964 
01965         /** construct from shape and copy values from the given array
01966          */
01967     MultiArray (const difference_type &shape, const_pointer init,
01968                          allocator_type const & alloc = allocator_type());
01969 
01970         /** copy constructor
01971          */
01972     MultiArray (const MultiArray &rhs)
01973     : MultiArrayView <N, T> (rhs.m_shape, rhs.m_stride, 0),
01974       m_alloc (rhs.m_alloc)
01975     {
01976         allocate (this->m_ptr, this->elementCount (), rhs.data ());
01977     }
01978 
01979         /** construct by copying from a MultiArrayView
01980          */
01981     template <class U, class C>
01982     MultiArray (const MultiArrayView<N, U, C>  &rhs,
01983                 allocator_type const & alloc = allocator_type());
01984 
01985         /** assignment.<br>
01986             If the size of \a rhs is the same as the left-hand side arrays's old size, only
01987             the data are copied. Otherwise, new storage is allocated, which invalidates all
01988             objects (array views, iterators) depending on the lhs array.
01989          */
01990     MultiArray & operator= (const MultiArray &rhs)
01991     {
01992         if (this != &rhs)
01993             this->copyOrReshape(rhs);
01994         return *this;
01995     }
01996 
01997         /** assignment from arbitrary MultiArrayView.<br>
01998             If the size of \a rhs is the same as the left-hand side arrays's old size, only
01999             the data are copied. Otherwise, new storage is allocated, which invalidates all
02000             objects (array views, iterators) depending on the lhs array.
02001          */
02002     template <class U, class C>
02003     MultiArray &operator= (const MultiArrayView<N, U, C> &rhs)
02004     {
02005         this->copyOrReshape(rhs);
02006         return *this;
02007     }
02008 
02009         /** Add-assignment from arbitrary MultiArrayView. Fails with
02010             <tt>PreconditionViolation</tt> exception when the shapes do not match.
02011          */
02012     template <class U, class C>
02013     MultiArray &operator+= (const MultiArrayView<N, U, C> &rhs)
02014     {
02015         view_type::operator+=(rhs);
02016         return *this;
02017     }
02018 
02019         /** Subtract-assignment from arbitrary MultiArrayView. Fails with
02020             <tt>PreconditionViolation</tt> exception when the shapes do not match.
02021          */
02022     template <class U, class C>
02023     MultiArray &operator-= (const MultiArrayView<N, U, C> &rhs)
02024     {
02025         view_type::operator-=(rhs);
02026         return *this;
02027     }
02028 
02029         /** Multiply-assignment from arbitrary MultiArrayView. Fails with
02030             <tt>PreconditionViolation</tt> exception when the shapes do not match.
02031          */
02032     template <class U, class C>
02033     MultiArray &operator*= (const MultiArrayView<N, U, C> &rhs)
02034     {
02035         view_type::operator*=(rhs);
02036         return *this;
02037     }
02038 
02039         /** Divide-assignment from arbitrary MultiArrayView. Fails with
02040             <tt>PreconditionViolation</tt> exception when the shapes do not match.
02041          */
02042     template <class U, class C>
02043     MultiArray &operator/= (const MultiArrayView<N, U, C> &rhs)
02044     {
02045         view_type::operator/=(rhs);
02046         return *this;
02047     }
02048 
02049         /** Add-assignment of a scalar.
02050          */
02051     MultiArray &operator+= (const T &rhs)
02052     {
02053         view_type::operator+=(rhs);
02054         return *this;
02055     }
02056 
02057         /** Subtract-assignment of a scalar.
02058          */
02059     MultiArray &operator-= (const T &rhs)
02060     {
02061         view_type::operator-=(rhs);
02062         return *this;
02063     }
02064 
02065         /** Multiply-assignment of a scalar.
02066          */
02067     MultiArray &operator*= (const T &rhs)
02068     {
02069         view_type::operator*=(rhs);
02070         return *this;
02071     }
02072 
02073         /** Divide-assignment of a scalar.
02074          */
02075     MultiArray &operator/= (const T &rhs)
02076     {
02077         view_type::operator/=(rhs);
02078         return *this;
02079     }
02080 
02081         /** destructor
02082          */
02083     ~MultiArray ()
02084     {
02085         deallocate (this->m_ptr, this->elementCount ());
02086     }
02087 
02088 
02089         /** init elements with a constant
02090          */
02091     template <class U>
02092     MultiArray & init(const U & init)
02093     {
02094         view_type::init(init);
02095         return *this;
02096     }
02097 
02098         /** Allocate new memory with the given shape and initialize with zeros.<br>
02099             <em>Note:</em> this operation invalidates all dependent objects
02100             (array views and iterators)
02101          */
02102     void reshape (const difference_type &shape)
02103     {
02104         reshape (shape, NumericTraits <T>::zero ());
02105     }
02106 
02107         /** Allocate new memory with the given shape and initialize it
02108             with the given value.<br>
02109             <em>Note:</em> this operation invalidates all dependent objects
02110             (array views and iterators)
02111          */
02112     void reshape (const difference_type &shape, const_reference init);
02113 
02114         /** Swap the contents with another MultiArray. This is fast,
02115             because no data are copied, but only pointers and shapes swapped.
02116             <em>Note:</em> this operation invalidates all dependent objects
02117             (array views and iterators)
02118          */
02119     void swap (MultiArray & other);
02120 
02121         /** sequential iterator pointing to the first array element.
02122          */
02123     iterator begin ()
02124     {
02125         return this->data();
02126     }
02127 
02128         /** sequential iterator pointing beyond the last array element.
02129          */
02130     iterator end ()
02131     {
02132         return this->data() + this->elementCount();
02133     }
02134 
02135         /** sequential const iterator pointing to the first array element.
02136          */
02137     const_iterator begin () const
02138     {
02139         return this->data();
02140     }
02141 
02142         /** sequential const iterator pointing beyond the last array element.
02143          */
02144     const_iterator end () const
02145     {
02146         return this->data() + this->elementCount();
02147     }
02148 
02149         /** get the allocator.
02150          */
02151     allocator_type const & allocator () const
02152     {
02153         return m_alloc;
02154     }
02155 };
02156 
02157 template <unsigned int N, class T, class A>
02158 MultiArray <N, T, A>::MultiArray (const difference_type &shape,
02159                                   allocator_type const & alloc)
02160 : MultiArrayView <N, T> (shape,
02161                          detail::defaultStride <MultiArrayView<N,T>::actual_dimension> (shape),
02162                          0),
02163   m_alloc(alloc)
02164 {
02165     if (N == 0)
02166     {
02167         this->m_shape [0] = 1;
02168         this->m_stride [0] = 0;
02169     }
02170     allocate (this->m_ptr, this->elementCount (), NumericTraits<T>::zero ());
02171 }
02172 
02173 template <unsigned int N, class T, class A>
02174 MultiArray <N, T, A>::MultiArray (const difference_type &shape, const_reference init,
02175                                   allocator_type const & alloc)
02176 : MultiArrayView <N, T> (shape,
02177                          detail::defaultStride <MultiArrayView<N,T>::actual_dimension> (shape),
02178                          0),
02179   m_alloc(alloc)
02180 {
02181     if (N == 0)
02182     {
02183         this->m_shape [0] = 1;
02184         this->m_stride [0] = 0;
02185     }
02186     allocate (this->m_ptr, this->elementCount (), init);
02187 }
02188 
02189 template <unsigned int N, class T, class A>
02190 MultiArray <N, T, A>::MultiArray (const difference_type &shape, const_pointer init,
02191                                   allocator_type const & alloc)
02192 : MultiArrayView <N, T> (shape,
02193                          detail::defaultStride <MultiArrayView<N,T>::actual_dimension> (shape),
02194                          0),
02195   m_alloc(alloc)
02196 {
02197     if (N == 0)
02198     {
02199         this->m_shape [0] = 1;
02200         this->m_stride [0] = 0;
02201     }
02202     allocate (this->m_ptr, this->elementCount (), init);
02203 }
02204 
02205 template <unsigned int N, class T, class A>
02206 template <class U, class C>
02207 MultiArray <N, T, A>::MultiArray(const MultiArrayView<N, U, C>  &rhs,
02208                                  allocator_type const & alloc)
02209 : MultiArrayView <N, T> (rhs.shape(),
02210                          detail::defaultStride <MultiArrayView<N,T>::actual_dimension>(rhs.shape()),
02211                          0),
02212   m_alloc (alloc)
02213 {
02214     allocate (this->m_ptr, rhs);
02215 }
02216 
02217 template <unsigned int N, class T, class A>
02218 template <class U, class C>
02219 void
02220 MultiArray <N, T, A>::copyOrReshape(const MultiArrayView<N, U, C> &rhs)
02221 {
02222     if (this->shape() == rhs.shape())
02223         this->copy(rhs);
02224     else
02225     {
02226         MultiArray t(rhs);
02227         this->swap(t);
02228     }
02229 }
02230 
02231 template <unsigned int N, class T, class A>
02232 void MultiArray <N, T, A>::reshape (const difference_type & new_shape,
02233                                     const_reference initial)
02234 {
02235     if (N== 0)
02236     {
02237         return;
02238     }
02239     else if(new_shape == this->shape())
02240     {
02241         this->init(initial);
02242     }
02243     else
02244     {
02245         difference_type new_stride = detail::defaultStride <MultiArrayView<N,T>::actual_dimension> (new_shape);
02246         difference_type_1 new_size = new_shape [MultiArrayView<N,T>::actual_dimension-1] * new_stride [MultiArrayView<N,T>::actual_dimension-1];
02247         T *new_ptr;
02248         allocate (new_ptr, new_size, initial);
02249         deallocate (this->m_ptr, this->elementCount ());
02250         this->m_ptr = new_ptr;
02251         this->m_shape = new_shape;
02252         this->m_stride = new_stride;
02253     }
02254 }
02255 
02256 
02257 template <unsigned int N, class T, class A>
02258 inline void
02259 MultiArray <N, T, A>::swap (MultiArray <N, T, A> & other)
02260 {
02261     if (this == &other)
02262         return;
02263     std::swap(this->m_shape,  other.m_shape);
02264     std::swap(this->m_stride, other.m_stride);
02265     std::swap(this->m_ptr,    other.m_ptr);
02266     std::swap(this->m_alloc,  other.m_alloc);
02267 }
02268 
02269 template <unsigned int N, class T, class A>
02270 void MultiArray <N, T, A>::allocate (pointer & ptr, difference_type_1 s,
02271                                      const_reference init)
02272 {
02273     ptr = m_alloc.allocate ((typename A::size_type)s);
02274     difference_type_1 i;
02275     try {
02276         for (i = 0; i < s; ++i)
02277             m_alloc.construct (ptr + i, init);
02278     }
02279     catch (...) {
02280         for (difference_type_1 j = 0; j < i; ++j)
02281             m_alloc.destroy (ptr + j);
02282         m_alloc.deallocate (ptr, (typename A::size_type)s);
02283         throw;
02284     }
02285 }
02286 
02287 template <unsigned int N, class T, class A>
02288 template <class U>
02289 void MultiArray <N, T, A>::allocate (pointer & ptr, difference_type_1 s,
02290                                      U const * init)
02291 {
02292     ptr = m_alloc.allocate ((typename A::size_type)s);
02293     difference_type_1 i;
02294     try {
02295         for (i = 0; i < s; ++i, ++init)
02296             m_alloc.construct (ptr + i, *init);
02297     }
02298     catch (...) {
02299         for (difference_type_1 j = 0; j < i; ++j)
02300             m_alloc.destroy (ptr + j);
02301         m_alloc.deallocate (ptr, (typename A::size_type)s);
02302         throw;
02303     }
02304 }
02305 
02306 template <unsigned int N, class T, class A>
02307 template <class U, class C>
02308 void MultiArray <N, T, A>::allocate (pointer & ptr, MultiArrayView<N, U, C> const & init)
02309 {
02310     difference_type_1 s = init.elementCount();
02311     ptr = m_alloc.allocate ((typename A::size_type)s);
02312     pointer p = ptr;
02313     try {
02314         detail::uninitializedCopyMultiArrayData(init.traverser_begin(), init.shape(),
02315                                                 p, m_alloc, MetaInt<actual_dimension-1>());
02316     }
02317     catch (...) {
02318         for (pointer pp = ptr; pp < p; ++pp)
02319             m_alloc.destroy (pp);
02320         m_alloc.deallocate (ptr, (typename A::size_type)s);
02321         throw;
02322     }
02323 }
02324 
02325 template <unsigned int N, class T, class A>
02326 inline void MultiArray <N, T, A>::deallocate (pointer & ptr, difference_type_1 s)
02327 {
02328     if (ptr == 0)
02329         return;
02330     for (difference_type_1 i = 0; i < s; ++i)
02331         m_alloc.destroy (ptr + i);
02332     m_alloc.deallocate (ptr, (typename A::size_type)s);
02333     ptr = 0;
02334 }
02335 
02336 /********************************************************/
02337 /*                                                      */
02338 /*              argument object factories               */
02339 /*                                                      */
02340 /********************************************************/
02341 
02342 template <unsigned int N, class T, class C>
02343 inline triple<typename MultiArrayView<N,T,C>::const_traverser,
02344               typename MultiArrayView<N,T,C>::difference_type,
02345               typename AccessorTraits<T>::default_const_accessor >
02346 srcMultiArrayRange( MultiArrayView<N,T,C> const & array )
02347 {
02348     return triple<typename MultiArrayView<N,T,C>::const_traverser,
02349                   typename MultiArrayView<N,T,C>::difference_type,
02350                   typename AccessorTraits<T>::default_const_accessor >
02351       ( array.traverser_begin(),
02352         array.shape(),
02353         typename AccessorTraits<T>::default_const_accessor() );
02354 }
02355 
02356 template <unsigned int N, class T, class C, class Accessor>
02357 inline triple<typename MultiArrayView<N,T,C>::const_traverser,
02358               typename MultiArrayView<N,T,C>::difference_type,
02359               Accessor >
02360 srcMultiArrayRange( MultiArrayView<N,T,C> const & array, Accessor a )
02361 {
02362     return triple<typename MultiArrayView<N,T,C>::const_traverser,
02363                   typename MultiArrayView<N,T,C>::difference_type,
02364                   Accessor >
02365       ( array.traverser_begin(),
02366         array.shape(),
02367         a);
02368 }
02369 
02370 template <unsigned int N, class T, class C>
02371 inline pair<typename MultiArrayView<N,T,C>::const_traverser,
02372             typename AccessorTraits<T>::default_const_accessor >
02373 srcMultiArray( MultiArrayView<N,T,C> const & array )
02374 {
02375     return pair<typename MultiArrayView<N,T,C>::const_traverser,
02376                 typename AccessorTraits<T>::default_const_accessor >
02377       ( array.traverser_begin(),
02378         typename AccessorTraits<T>::default_const_accessor() );
02379 }
02380 
02381 template <unsigned int N, class T, class C, class Accessor>
02382 inline pair<typename MultiArrayView<N,T,C>::const_traverser,
02383             Accessor >
02384 srcMultiArray( MultiArrayView<N,T,C> const & array, Accessor a )
02385 {
02386     return pair<typename MultiArrayView<N,T,C>::const_traverser,
02387                 Accessor >
02388       ( array.traverser_begin(), a );
02389 }
02390 
02391 template <unsigned int N, class T, class C>
02392 inline triple<typename MultiArrayView<N,T,C>::traverser,
02393               typename MultiArrayView<N,T,C>::difference_type,
02394               typename AccessorTraits<T>::default_accessor >
02395 destMultiArrayRange( MultiArrayView<N,T,C> & array )
02396 {
02397     return triple<typename MultiArrayView<N,T,C>::traverser,
02398                   typename MultiArrayView<N,T,C>::difference_type,
02399                   typename AccessorTraits<T>::default_accessor >
02400       ( array.traverser_begin(),
02401         array.shape(),
02402         typename AccessorTraits<T>::default_accessor() );
02403 }
02404 
02405 template <unsigned int N, class T, class C, class Accessor>
02406 inline triple<typename MultiArrayView<N,T,C>::traverser,
02407               typename MultiArrayView<N,T,C>::difference_type,
02408               Accessor >
02409 destMultiArrayRange( MultiArrayView<N,T,C> & array, Accessor a )
02410 {
02411     return triple<typename MultiArrayView<N,T,C>::traverser,
02412                   typename MultiArrayView<N,T,C>::difference_type,
02413                   Accessor >
02414       ( array.traverser_begin(),
02415         array.shape(),
02416         a );
02417 }
02418 
02419 template <unsigned int N, class T, class C>
02420 inline pair<typename MultiArrayView<N,T,C>::traverser,
02421             typename AccessorTraits<T>::default_accessor >
02422 destMultiArray( MultiArrayView<N,T,C> & array )
02423 {
02424     return pair<typename MultiArrayView<N,T,C>::traverser,
02425                 typename AccessorTraits<T>::default_accessor >
02426         ( array.traverser_begin(),
02427           typename AccessorTraits<T>::default_accessor() );
02428 }
02429 
02430 template <unsigned int N, class T, class C, class Accessor>
02431 inline pair<typename MultiArrayView<N,T,C>::traverser,
02432             Accessor >
02433 destMultiArray( MultiArrayView<N,T,C> & array, Accessor a )
02434 {
02435     return pair<typename MultiArrayView<N,T,C>::traverser,
02436                 Accessor >
02437         ( array.traverser_begin(), a );
02438 }
02439 
02440 /********************************************************************/
02441 
02442 template <class PixelType, class Accessor>
02443 inline triple<ConstStridedImageIterator<PixelType>,
02444               ConstStridedImageIterator<PixelType>, Accessor>
02445 srcImageRange(const MultiArrayView<2, PixelType, StridedArrayTag> & img, Accessor a)
02446 {
02447     ConstStridedImageIterator<PixelType>
02448         ul(img.data(), 1, img.stride(0), img.stride(1));
02449     return triple<ConstStridedImageIterator<PixelType>,
02450                   ConstStridedImageIterator<PixelType>,
02451                   Accessor>(
02452                       ul, ul + Size2D(img.shape(0), img.shape(1)), a);
02453 }
02454 
02455 template <class PixelType, class Accessor>
02456 inline pair<ConstStridedImageIterator<PixelType>, Accessor>
02457 srcImage(const MultiArrayView<2, PixelType, StridedArrayTag> & img, Accessor a)
02458 {
02459     ConstStridedImageIterator<PixelType>
02460         ul(img.data(), 1, img.stride(0), img.stride(1));
02461     return pair<ConstStridedImageIterator<PixelType>, Accessor>
02462         (ul, a);
02463 }
02464 
02465 template <class PixelType, class Accessor>
02466 inline triple<StridedImageIterator<PixelType>,
02467               StridedImageIterator<PixelType>, Accessor>
02468 destImageRange(MultiArrayView<2, PixelType, StridedArrayTag> & img, Accessor a)
02469 {
02470     StridedImageIterator<PixelType>
02471         ul(img.data(), 1, img.stride(0), img.stride(1));
02472     return triple<StridedImageIterator<PixelType>,
02473                   StridedImageIterator<PixelType>,
02474                   Accessor>(
02475                       ul, ul + Size2D(img.shape(0), img.shape(1)), a);
02476 }
02477 
02478 template <class PixelType, class Accessor>
02479 inline pair<StridedImageIterator<PixelType>, Accessor>
02480 destImage(MultiArrayView<2, PixelType, StridedArrayTag> & img, Accessor a)
02481 {
02482     StridedImageIterator<PixelType>
02483         ul(img.data(), 1, img.stride(0), img.stride(1));
02484     return pair<StridedImageIterator<PixelType>, Accessor>
02485         (ul, a);
02486 }
02487 
02488 template <class PixelType, class Accessor>
02489 inline pair<StridedImageIterator<PixelType>, Accessor>
02490 maskImage(MultiArrayView<2, PixelType, StridedArrayTag> & img, Accessor a)
02491 {
02492     StridedImageIterator<PixelType>
02493         ul(img.data(), 1, img.stride(0), img.stride(1));
02494     return pair<StridedImageIterator<PixelType>, Accessor>
02495         (ul, a);
02496 }
02497 
02498 // -------------------------------------------------------------------
02499 
02500 template <class PixelType>
02501 inline triple<ConstStridedImageIterator<PixelType>,
02502               ConstStridedImageIterator<PixelType>,
02503               typename AccessorTraits<PixelType>::default_const_accessor>
02504 srcImageRange(MultiArrayView<2, PixelType, StridedArrayTag> const & img)
02505 {
02506     ConstStridedImageIterator<PixelType>
02507         ul(img.data(), 1, img.stride(0), img.stride(1));
02508     typedef typename AccessorTraits<PixelType>::default_const_accessor Accessor;
02509     return triple<ConstStridedImageIterator<PixelType>,
02510                   ConstStridedImageIterator<PixelType>,
02511                   Accessor>
02512         (ul, ul + Size2D(img.shape(0), img.shape(1)), Accessor());
02513 }
02514 
02515 template <class PixelType>
02516 inline triple<ConstImageIterator<PixelType>,
02517               ConstImageIterator<PixelType>,
02518               typename AccessorTraits<PixelType>::default_const_accessor>
02519 srcImageRange(MultiArrayView<2, PixelType, UnstridedArrayTag> const & img)
02520 {
02521     ConstImageIterator<PixelType>
02522         ul(img.data(), img.stride(1));
02523     typedef typename AccessorTraits<PixelType>::default_const_accessor Accessor;
02524     return triple<ConstImageIterator<PixelType>,
02525                   ConstImageIterator<PixelType>,
02526                   Accessor>
02527         (ul, ul + Size2D(img.shape(0), img.shape(1)), Accessor());
02528 }
02529 
02530 template <class PixelType>
02531 inline pair< ConstStridedImageIterator<PixelType>,
02532              typename AccessorTraits<PixelType>::default_const_accessor>
02533 srcImage(MultiArrayView<2, PixelType, StridedArrayTag> const & img)
02534 {
02535     ConstStridedImageIterator<PixelType>
02536         ul(img.data(), 1, img.stride(0), img.stride(1));
02537     typedef typename AccessorTraits<PixelType>::default_const_accessor Accessor;
02538     return pair<ConstStridedImageIterator<PixelType>,
02539                 Accessor>
02540         (ul, Accessor());
02541 }
02542 
02543 template <class PixelType>
02544 inline pair< ConstImageIterator<PixelType>,
02545              typename AccessorTraits<PixelType>::default_const_accessor>
02546 srcImage(MultiArrayView<2, PixelType, UnstridedArrayTag> const & img)
02547 {
02548     ConstImageIterator<PixelType>
02549         ul(img.data(), img.stride(1));
02550     typedef typename AccessorTraits<PixelType>::default_const_accessor Accessor;
02551     return pair<ConstImageIterator<PixelType>,
02552                 Accessor>
02553         (ul, Accessor());
02554 }
02555 
02556 template <class PixelType>
02557 inline triple< StridedImageIterator<PixelType>,
02558                StridedImageIterator<PixelType>,
02559                typename AccessorTraits<PixelType>::default_accessor>
02560 destImageRange(MultiArrayView<2, PixelType, StridedArrayTag> & img)
02561 {
02562     StridedImageIterator<PixelType>
02563         ul(img.data(), 1, img.stride(0), img.stride(1));
02564     typedef typename AccessorTraits<PixelType>::default_accessor Accessor;
02565     return triple<StridedImageIterator<PixelType>,
02566                   StridedImageIterator<PixelType>,
02567                   Accessor>
02568         (ul, ul + Size2D(img.shape(0), img.shape(1)), Accessor());
02569 }
02570 
02571 template <class PixelType>
02572 inline triple< ImageIterator<PixelType>,
02573                ImageIterator<PixelType>,
02574                typename AccessorTraits<PixelType>::default_accessor>
02575 destImageRange(MultiArrayView<2, PixelType, UnstridedArrayTag> & img)
02576 {
02577     ImageIterator<PixelType>
02578         ul(img.data(), img.stride(1));
02579     typedef typename AccessorTraits<PixelType>::default_accessor Accessor;
02580     return triple<ImageIterator<PixelType>,
02581                   ImageIterator<PixelType>,
02582                   Accessor>
02583         (ul, ul + Size2D(img.shape(0), img.shape(1)), Accessor());
02584 }
02585 
02586 template <class PixelType>
02587 inline pair< StridedImageIterator<PixelType>,
02588              typename AccessorTraits<PixelType>::default_accessor>
02589 destImage(MultiArrayView<2, PixelType, StridedArrayTag> & img)
02590 {
02591     StridedImageIterator<PixelType>
02592         ul(img.data(), 1, img.stride(0), img.stride(1));
02593     typedef typename AccessorTraits<PixelType>::default_accessor Accessor;
02594     return pair<StridedImageIterator<PixelType>, Accessor>
02595         (ul, Accessor());
02596 }
02597 
02598 template <class PixelType>
02599 inline pair< ImageIterator<PixelType>,
02600              typename AccessorTraits<PixelType>::default_accessor>
02601 destImage(MultiArrayView<2, PixelType, UnstridedArrayTag> & img)
02602 {
02603     ImageIterator<PixelType> ul(img.data(), img.stride(1));
02604     typedef typename AccessorTraits<PixelType>::default_accessor Accessor;
02605     return pair<ImageIterator<PixelType>, Accessor>(ul, Accessor());
02606 }
02607 
02608 template <class PixelType>
02609 inline pair< ConstStridedImageIterator<PixelType>,
02610              typename AccessorTraits<PixelType>::default_accessor>
02611 maskImage(MultiArrayView<2, PixelType, StridedArrayTag> const & img)
02612 {
02613     ConstStridedImageIterator<PixelType>
02614         ul(img.data(), 1, img.stride(0), img.stride(1));
02615     typedef typename AccessorTraits<PixelType>::default_accessor Accessor;
02616     return pair<ConstStridedImageIterator<PixelType>, Accessor>
02617         (ul, Accessor());
02618 }
02619 
02620 template <class PixelType>
02621 inline pair< ConstImageIterator<PixelType>,
02622              typename AccessorTraits<PixelType>::default_accessor>
02623 maskImage(MultiArrayView<2, PixelType, UnstridedArrayTag> const & img)
02624 {
02625     ConstImageIterator<PixelType>
02626         ul(img.data(), img.stride(1));
02627     typedef typename AccessorTraits<PixelType>::default_accessor Accessor;
02628     return pair<ConstImageIterator<PixelType>, Accessor>
02629         (ul, Accessor());
02630 }
02631 
02632 /********************************************************/
02633 /*                                                      */
02634 /*                  makeBasicImageView                  */
02635 /*                                                      */
02636 /********************************************************/
02637 
02638 /** \addtogroup MultiArrayToImage Wrap a \ref vigra::MultiArrayView in
02639                                   a \ref vigra::BasicImageView
02640 */
02641 //@{
02642 /** Create a \ref vigra::BasicImageView from an unstrided 2-dimensional
02643     \ref vigra::MultiArrayView.
02644 
02645     The \ref vigra::BasicImageView will have the same <tt>value_type </tt>
02646     as the original \ref vigra::MultiArrayView.
02647 */
02648 template <class T>
02649 BasicImageView <T>
02650 makeBasicImageView (MultiArrayView <2, T, UnstridedArrayTag> const &array)
02651 {
02652     return BasicImageView <T> (array.data (), array.shape (0),
02653                                array.shape (1));
02654 }
02655 
02656 /** Create a \ref vigra::BasicImageView from a 3-dimensional
02657     \ref vigra::MultiArray.
02658 
02659     This wrapper flattens the two innermost dimensions of the array
02660     into single rows of the resulting image.
02661     The \ref vigra::BasicImageView will have the same <tt>value_type </tt>
02662     as the original \ref vigra::MultiArray.
02663 */
02664 template <class T>
02665 BasicImageView <T>
02666 makeBasicImageView (MultiArray <3, T> const &array)
02667 {
02668     return BasicImageView <T> (array.data (),
02669                                array.shape (0)*array.shape (1), array.shape (2));
02670 }
02671 
02672 /** Create a \ref vigra::BasicImageView from a 3-dimensional
02673     \ref vigra::MultiArray.
02674 
02675     This wrapper only works if <tt>T</tt> is a scalar type and the
02676     array's innermost dimension has size 3. It then re-interprets
02677     the data array as a 2-dimensional array with value_type
02678     <tt>RGBValue<T></tt>.
02679 */
02680 template <class T>
02681 BasicImageView <RGBValue<T> >
02682 makeRGBImageView (MultiArray<3, T> const &array)
02683 {
02684     vigra_precondition (
02685         array.shape (0) == 3, "makeRGBImageView(): array.shape(0) must be 3.");
02686     return BasicImageView <RGBValue<T> > (
02687         reinterpret_cast <RGBValue <T> *> (array.data ()),
02688         array.shape (1), array.shape (2));
02689 }
02690 
02691 //@}
02692 
02693 } // namespace vigra
02694 #undef VIGRA_ASSERT_INSIDE
02695 #endif // VIGRA_MULTI_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)