37 #ifndef VIGRA_MULTI_ARRAY_HXX
38 #define VIGRA_MULTI_ARRAY_HXX
42 #include "accessor.hxx"
43 #include "tinyvector.hxx"
44 #include "rgbvalue.hxx"
45 #include "basicimageview.hxx"
46 #include "imageiterator.hxx"
47 #include "numerictraits.hxx"
48 #include "multi_iterator.hxx"
49 #include "metaprogramming.hxx"
50 #include "mathutil.hxx"
53 #ifdef VIGRA_CHECK_BOUNDS
54 #define VIGRA_ASSERT_INSIDE(diff) \
55 vigra_precondition(this->isInside(diff), "Index out of bounds")
57 #define VIGRA_ASSERT_INSIDE(diff)
75 template <
unsigned int N>
76 inline TinyVector <MultiArrayIndex, N>
77 defaultStride(
const TinyVector <MultiArrayIndex, N> &shape)
79 TinyVector <MultiArrayIndex, N> ret;
81 for (
int i = 1; i < (int)N; ++i)
82 ret [i] = ret [i-1] * shape [i-1];
99 struct ScanOrderToOffset
104 const TinyVector <MultiArrayIndex, N> & stride)
106 return stride[N-K] * (d % shape[N-K]) +
107 ScanOrderToOffset<K-1>::exec(d / shape[N-K], shape, stride);
112 struct ScanOrderToOffset<1>
117 const TinyVector <MultiArrayIndex, N> & stride)
119 return stride[N-1] * d;
124 struct ScanOrderToCoordinate
129 TinyVector <MultiArrayIndex, N> & result)
131 result[N-K] = (d % shape[N-K]);
132 ScanOrderToCoordinate<K-1>::exec(d / shape[N-K], shape, result);
137 struct ScanOrderToCoordinate<1>
142 TinyVector <MultiArrayIndex, N> & result)
149 struct CoordinateToScanOrder
153 exec(
const TinyVector <MultiArrayIndex, N> &shape,
154 const TinyVector <MultiArrayIndex, N> & coordinate)
156 return coordinate[N-K] + shape[N-K] * CoordinateToScanOrder<K-1>::exec(shape, coordinate);
161 struct CoordinateToScanOrder<1>
165 exec(
const TinyVector <MultiArrayIndex, N> & ,
166 const TinyVector <MultiArrayIndex, N> & coordinate)
168 return coordinate[N-1];
174 struct CoordinatesToOffest
178 exec(
const TinyVector <MultiArrayIndex, N> & stride,
MultiArrayIndex x)
180 return stride[0] * x;
186 return stride[0] * x + stride[1] * y;
191 struct CoordinatesToOffest<UnstridedArrayTag>
203 return x + stride[1] * y;
222 template <
class Str
ideTag,
unsigned int N>
225 typedef StrideTag type;
228 template <
class Str
ideTag>
229 struct MaybeStrided <StrideTag, 0>
231 typedef StridedArrayTag type;
249 struct MultiIteratorChooser
253 template <
unsigned int N,
class T,
class REFERENCE,
class POINTER>
274 struct MultiIteratorChooser <StridedArrayTag>
276 template <
unsigned int N,
class T,
class REFERENCE,
class POINTER>
279 typedef StridedMultiIterator <N, T, REFERENCE, POINTER> type;
297 struct MultiIteratorChooser <UnstridedArrayTag>
299 template <
unsigned int N,
class T,
class REFERENCE,
class POINTER>
302 typedef MultiIterator <N, T, REFERENCE, POINTER> type;
312 template <
class DestIterator,
class Shape,
class T>
314 initMultiArrayData(DestIterator d, Shape
const & shape, T
const & init, MetaInt<0>)
316 DestIterator dend = d + shape[0];
323 template <
class DestIterator,
class Shape,
class T,
int N>
325 initMultiArrayData(DestIterator d, Shape
const & shape, T
const & init, MetaInt<N>)
327 DestIterator dend = d + shape[N];
330 initMultiArrayData(d.begin(), shape, init, MetaInt<N-1>());
334 #define VIGRA_COPY_MULTI_ARRAY_DATA(name, op) \
335 template <class SrcIterator, class Shape, class DestIterator> \
337 name##MultiArrayData(SrcIterator s, Shape const & shape, DestIterator d, MetaInt<0>) \
339 SrcIterator send = s + shape[0]; \
340 for(; s < send; ++s, ++d) \
342 *d op detail::RequiresExplicitCast<typename DestIterator::value_type>::cast(*s); \
346 template <class SrcIterator, class Shape, class DestIterator, int N> \
348 name##MultiArrayData(SrcIterator s, Shape const & shape, DestIterator d, MetaInt<N>) \
350 SrcIterator send = s + shape[N]; \
351 for(; s < send; ++s, ++d) \
353 name##MultiArrayData(s.begin(), shape, d.begin(), MetaInt<N-1>()); \
357 template <class DestIterator, class Shape, class T> \
359 name##ScalarMultiArrayData(DestIterator d, Shape const & shape, T const & init, MetaInt<0>) \
361 DestIterator dend = d + shape[0]; \
362 for(; d < dend; ++d) \
364 *d op detail::RequiresExplicitCast<typename DestIterator::value_type>::cast(init); \
368 template <class DestIterator, class Shape, class T, int N> \
370 name##ScalarMultiArrayData(DestIterator d, Shape const & shape, T const & init, MetaInt<N>) \
372 DestIterator dend = d + shape[N]; \
373 for(; d < dend; ++d) \
375 name##ScalarMultiArrayData(d.begin(), shape, init, MetaInt<N-1>()); \
379 VIGRA_COPY_MULTI_ARRAY_DATA(copy, =)
380 VIGRA_COPY_MULTI_ARRAY_DATA(copyAdd, +=)
381 VIGRA_COPY_MULTI_ARRAY_DATA(copySub, -=)
382 VIGRA_COPY_MULTI_ARRAY_DATA(copyMul, *=)
383 VIGRA_COPY_MULTI_ARRAY_DATA(copyDiv, /=)
385 #undef VIGRA_COPY_MULTI_ARRAY_DATA
387 template <
class SrcIterator,
class Shape,
class T,
class ALLOC>
389 uninitializedCopyMultiArrayData(SrcIterator s, Shape
const & shape, T * & d, ALLOC & a, MetaInt<0>)
391 SrcIterator send = s + shape[0];
392 for(; s < send; ++s, ++d)
394 a.construct(d, static_cast<T const &>(*s));
398 template <
class SrcIterator,
class Shape,
class T,
class ALLOC,
int N>
400 uninitializedCopyMultiArrayData(SrcIterator s, Shape
const & shape, T * & d, ALLOC & a, MetaInt<N>)
402 SrcIterator send = s + shape[N];
405 uninitializedCopyMultiArrayData(s.begin(), shape, d, a, MetaInt<N-1>());
409 template <
class SrcIterator,
class Shape,
class T,
class Functor>
411 reduceOverMultiArray(SrcIterator s, Shape
const & shape, T & result, Functor
const & f, MetaInt<0>)
413 SrcIterator send = s + shape[0];
420 template <
class SrcIterator,
class Shape,
class T,
class Functor,
int N>
422 reduceOverMultiArray(SrcIterator s, Shape
const & shape, T & result, Functor
const & f, MetaInt<N>)
424 SrcIterator send = s + shape[N];
427 reduceOverMultiArray(s.begin(), shape, result, f, MetaInt<N-1>());
431 struct MaxNormReduceFunctor
433 template <
class T,
class U>
434 void operator()(T & result, U
const & u)
const
442 struct L1NormReduceFunctor
444 template <
class T,
class U>
445 void operator()(T & result, U
const & u)
const
451 struct SquaredL2NormReduceFunctor
453 template <
class T,
class U>
454 void operator()(T & result, U
const & u)
const
461 struct WeightedL2NormReduceFunctor
465 WeightedL2NormReduceFunctor(T s)
470 void operator()(T & result, U
const & u)
const
476 struct SumReduceFunctor
478 template <
class T,
class U>
479 void operator()(T & result, U
const & u)
const
485 struct ProdReduceFunctor
487 template <
class T,
class U>
488 void operator()(T & result, U
const & u)
const
494 struct MinmaxReduceFunctor
496 template <
class T,
class U>
497 void operator()(T & result, U
const & u)
const
501 if(result.second < u)
506 struct MeanVarianceReduceFunctor
508 template <
class T,
class U>
509 void operator()(T & result, U
const & u)
const
512 typename T::second_type t1 = u - result.second;
513 typename T::second_type t2 = t1 / result.first;
515 result.third += (result.first-1.0)*t1*t2;
519 struct AllTrueReduceFunctor
521 template <
class T,
class U>
522 void operator()(T & result, U
const & u)
const
524 result = result && (u != NumericTraits<U>::zero());
528 struct AnyTrueReduceFunctor
530 template <
class T,
class U>
531 void operator()(T & result, U
const & u)
const
533 result = result || (u != NumericTraits<U>::zero());
537 template <
class SrcIterator,
class Shape,
class DestIterator>
539 equalityOfMultiArrays(SrcIterator s, Shape
const & shape, DestIterator d, MetaInt<0>)
541 SrcIterator send = s + shape[0];
542 for(; s < send; ++s, ++d)
550 template <
class SrcIterator,
class Shape,
class DestIterator,
int N>
552 equalityOfMultiArrays(SrcIterator s, Shape
const & shape, DestIterator d, MetaInt<N>)
554 SrcIterator send = s + shape[N];
555 for(; s < send; ++s, ++d)
557 if(!equalityOfMultiArrays(s.begin(), shape, d.begin(), MetaInt<N-1>()))
564 template <
class SrcIterator,
class Shape,
class DestIterator>
566 swapDataImpl(SrcIterator s, Shape
const & shape, DestIterator d, MetaInt<0>)
568 SrcIterator send = s + shape[0];
569 for(; s < send; ++s, ++d)
573 template <
class SrcIterator,
class Shape,
class DestIterator,
int N>
575 swapDataImpl(SrcIterator s, Shape
const & shape, DestIterator d, MetaInt<N>)
577 SrcIterator send = s + shape[N];
578 for(; s < send; ++s, ++d)
579 swapDataImpl(s.begin(), shape, d.begin(), MetaInt<N-1>());
592 template <
unsigned int N,
class T,
class C = Unstr
idedArrayTag>
593 class MultiArrayView;
594 template <
unsigned int N,
class T,
class A = std::allocator<T> >
597 namespace multi_math {
600 struct MultiMathOperand;
604 template <
unsigned int N,
class T,
class C,
class E>
605 void assign(MultiArrayView<N, T, C>, MultiMathOperand<E>
const &);
607 template <
unsigned int N,
class T,
class C,
class E>
608 void plusAssign(MultiArrayView<N, T, C>, MultiMathOperand<E>
const &);
610 template <
unsigned int N,
class T,
class C,
class E>
611 void minusAssign(MultiArrayView<N, T, C>, MultiMathOperand<E>
const &);
613 template <
unsigned int N,
class T,
class C,
class E>
614 void multiplyAssign(MultiArrayView<N, T, C>, MultiMathOperand<E>
const &);
616 template <
unsigned int N,
class T,
class C,
class E>
617 void divideAssign(MultiArrayView<N, T, C>, MultiMathOperand<E>
const &);
619 template <
unsigned int N,
class T,
class A,
class E>
620 void assignOrResize(MultiArray<N, T, A> &, MultiMathOperand<E>
const &);
622 template <
unsigned int N,
class T,
class A,
class E>
623 void plusAssignOrResize(MultiArray<N, T, A> &, MultiMathOperand<E>
const &);
625 template <
unsigned int N,
class T,
class A,
class E>
626 void minusAssignOrResize(MultiArray<N, T, A> &, MultiMathOperand<E>
const &);
628 template <
unsigned int N,
class T,
class A,
class E>
629 void multiplyAssignOrResize(MultiArray<N, T, A> &, MultiMathOperand<E>
const &);
631 template <
unsigned int N,
class T,
class A,
class E>
632 void divideAssignOrResize(MultiArray<N, T, A> &, MultiMathOperand<E>
const &);
638 template <
class T>
class FindSum;
640 struct UnsuitableTypeForExpandElements {};
643 struct ExpandElementResult
645 typedef UnsuitableTypeForExpandElements type;
649 struct ExpandElementResult<std::complex<T> >
659 struct ExpandElementResult<FFTWComplex<T> >
665 template <
class T,
int SIZE>
666 struct ExpandElementResult<TinyVector<T, SIZE> >
669 enum { size = SIZE };
672 template <
class T,
unsigned int R,
unsigned int G,
unsigned int B>
673 struct ExpandElementResult<RGBValue<T, R, G, B> >
679 #define VIGRA_DEFINE_EXPAND_ELEMENT_RESULT(TYPE) \
681 struct ExpandElementResult<TYPE> \
687 VIGRA_DEFINE_EXPAND_ELEMENT_RESULT(
bool)
688 VIGRA_DEFINE_EXPAND_ELEMENT_RESULT(
char)
689 VIGRA_DEFINE_EXPAND_ELEMENT_RESULT(
signed char)
690 VIGRA_DEFINE_EXPAND_ELEMENT_RESULT(
signed short)
691 VIGRA_DEFINE_EXPAND_ELEMENT_RESULT(
signed int)
692 VIGRA_DEFINE_EXPAND_ELEMENT_RESULT(
signed long)
693 VIGRA_DEFINE_EXPAND_ELEMENT_RESULT(
signed long long)
694 VIGRA_DEFINE_EXPAND_ELEMENT_RESULT(
unsigned char)
695 VIGRA_DEFINE_EXPAND_ELEMENT_RESULT(
unsigned short)
696 VIGRA_DEFINE_EXPAND_ELEMENT_RESULT(
unsigned int)
697 VIGRA_DEFINE_EXPAND_ELEMENT_RESULT(
unsigned long)
698 VIGRA_DEFINE_EXPAND_ELEMENT_RESULT(
unsigned long long)
699 VIGRA_DEFINE_EXPAND_ELEMENT_RESULT(
float)
700 VIGRA_DEFINE_EXPAND_ELEMENT_RESULT(
double)
701 VIGRA_DEFINE_EXPAND_ELEMENT_RESULT(
long double)
703 #undef VIGRA_DEFINE_EXPAND_ELEMENT_RESULT
712 template <
unsigned int N,
class T,
class C>
713 struct NormTraits<MultiArrayView<N, T, C> >
715 typedef MultiArrayView<N, T, C> Type;
716 typedef typename NormTraits<T>::SquaredNormType SquaredNormType;
717 typedef typename SquareRootTraits<SquaredNormType>::SquareRootResult NormType;
720 template <
unsigned int N,
class T,
class A>
721 struct NormTraits<MultiArray<N, T, A> >
722 :
public NormTraits<MultiArrayView<N, T, UnstridedArrayTag> >
724 typedef NormTraits<MultiArrayView<N, T, UnstridedArrayTag> > BaseType;
725 typedef MultiArray<N, T, A> Type;
726 typedef typename BaseType::SquaredNormType SquaredNormType;
727 typedef typename BaseType::NormType NormType;
765 template <
unsigned int N,
class T,
class Str
ideTag>
821 typedef typename vigra::detail::MultiIteratorChooser <
822 StrideTag>::template Traverser <actual_dimension, T, T &, T *>::type
traverser;
826 typedef typename vigra::detail::MultiIteratorChooser <
827 StrideTag>::template Traverser <actual_dimension, T, T const &, T const *>::type
const_traverser;
839 typedef typename difference_type::value_type diff_zero_t;
854 template <
class U,
class CN>
857 template <
class U,
class CN>
863 template <
class U,
class CN>
871 return m_stride[0] <= 1;
874 bool checkInnerStride(StridedArrayTag)
885 : m_shape (diff_zero_t(0)), m_stride (diff_zero_t(0)), m_ptr (0)
892 m_stride (detail::defaultStride <
MultiArrayView<N,T>::actual_dimension> (shape)),
902 const difference_type &
stride,
908 vigra_precondition(checkInnerStride(StrideTag()),
909 "MultiArrayView<..., UnstridedArrayTag>::MultiArrayView(): First dimension of given array is not unstrided.");
933 template<
class U,
class C1>
936 vigra_precondition(this->
shape() == rhs.
shape(),
937 "MultiArrayView::operator=() size mismatch.");
945 template<
class U,
class C1>
951 template<
class U,
class C1>
957 template<
class U,
class C1>
963 template<
class U,
class C1>
970 detail::copyAddScalarMultiArrayData(
traverser_begin(),
shape(), rhs, MetaInt<actual_dimension-1>());
978 detail::copySubScalarMultiArrayData(
traverser_begin(),
shape(), rhs, MetaInt<actual_dimension-1>());
986 detail::copyMulScalarMultiArrayData(
traverser_begin(),
shape(), rhs, MetaInt<actual_dimension-1>());
994 detail::copyDivScalarMultiArrayData(
traverser_begin(),
shape(), rhs, MetaInt<actual_dimension-1>());
1001 template<
class Expression>
1004 multi_math::detail::assign(*
this, rhs);
1011 template<
class Expression>
1014 multi_math::detail::plusAssign(*
this, rhs);
1021 template<
class Expression>
1024 multi_math::detail::minusAssign(*
this, rhs);
1031 template<
class Expression>
1034 multi_math::detail::multiplyAssign(*
this, rhs);
1041 template<
class Expression>
1044 multi_math::detail::divideAssign(*
this, rhs);
1052 VIGRA_ASSERT_INSIDE(d);
1053 return m_ptr [
dot (d, m_stride)];
1060 VIGRA_ASSERT_INSIDE(d);
1061 return m_ptr [
dot (d, m_stride)];
1066 template <
unsigned int M>
1085 return m_ptr [detail::ScanOrderToOffset<actual_dimension>::exec(d, m_shape, m_stride)];
1101 return m_ptr [detail::ScanOrderToOffset<actual_dimension>::exec(d, m_shape, m_stride)];
1108 difference_type result;
1109 detail::ScanOrderToCoordinate<actual_dimension>::exec(d, m_shape, result);
1117 return detail::CoordinateToScanOrder<actual_dimension>::exec(m_shape, d);
1125 return m_ptr [detail::CoordinatesToOffest<StrideTag>::exec(m_stride, x)];
1133 return m_ptr [detail::CoordinatesToOffest<StrideTag>::exec(m_stride, x, y)];
1138 reference
operator() (difference_type_1 x, difference_type_1 y, difference_type_1 z)
1141 return m_ptr [m_stride[0]*x + m_stride[1]*y + m_stride[2]*z];
1147 difference_type_1 z, difference_type_1 u)
1150 return m_ptr [m_stride[0]*x + m_stride[1]*y + m_stride[2]*z + m_stride[3]*u];
1155 reference
operator() (difference_type_1 x, difference_type_1 y, difference_type_1 z,
1156 difference_type_1 u, difference_type_1 v)
1159 return m_ptr [m_stride[0]*x + m_stride[1]*y + m_stride[2]*z + m_stride[3]*u + m_stride[4]*v];
1167 return m_ptr [detail::CoordinatesToOffest<StrideTag>::exec(m_stride, x)];
1172 const_reference
operator() (difference_type_1 x, difference_type_1 y)
const
1175 return m_ptr [detail::CoordinatesToOffest<StrideTag>::exec(m_stride, x, y)];
1180 const_reference
operator() (difference_type_1 x, difference_type_1 y, difference_type_1 z)
const
1183 return m_ptr [m_stride[0]*x + m_stride[1]*y + m_stride[2]*z];
1188 const_reference
operator() (difference_type_1 x, difference_type_1 y,
1189 difference_type_1 z, difference_type_1 u)
const
1192 return m_ptr [m_stride[0]*x + m_stride[1]*y + m_stride[2]*z + m_stride[3]*u];
1197 const_reference
operator() (difference_type_1 x, difference_type_1 y, difference_type_1 z,
1198 difference_type_1 u, difference_type_1 v)
const
1201 return m_ptr [m_stride[0]*x + m_stride[1]*y + m_stride[2]*z + m_stride[3]*u + m_stride[4]*v];
1209 detail::copyScalarMultiArrayData(
traverser_begin(),
shape(), init, MetaInt<actual_dimension-1>());
1220 this->copyImpl(rhs);
1225 template <
class U,
class CN>
1228 this->copyImpl(rhs);
1245 template <
class T2,
class C2>
1260 for(
unsigned int k = dimension+1; k < N; ++k)
1279 template <
unsigned int M>
1296 template <
unsigned int M>
1317 template <
unsigned int M>
1318 MultiArrayView <N-1, T,
typename vigra::detail::MaybeStrided<StrideTag, M>::type >
1319 bind (difference_type_1 d)
const;
1368 bindAt (difference_type_1 m, difference_type_1 d)
const;
1389 vigra_precondition(0 <= i && i < ExpandElementResult<T>::size,
1390 "MultiArrayView::bindElementChannel(i): 'i' out of range.");
1461 const difference_type &q)
const
1463 const difference_type_1 offset =
dot (m_stride, p);
1476 difference_type shape =
m_shape;
1477 for (
unsigned int i = 0; i < actual_dimension; ++i)
1500 difference_type
shape(m_shape.begin(), difference_type::ReverseCopy),
1501 stride(m_stride.begin(), difference_type::ReverseCopy);
1556 difference_type_1 ret = m_shape[0];
1557 for(
int i = 1; i < actual_dimension; ++i)
1579 difference_type_1
size (difference_type_1 n)
const
1587 difference_type_1
shape (difference_type_1 n)
const
1603 return m_stride [n];
1608 template <
class U,
class C1>
1619 template <
class U,
class C1>
1629 for(
int d=0; d<actual_dimension; ++d)
1630 if(p[d] < 0 || p[d] >=
shape(d))
1643 detail::AllTrueReduceFunctor(),
1644 MetaInt<actual_dimension-1>());
1656 detail::AnyTrueReduceFunctor(),
1657 MetaInt<actual_dimension-1>());
1668 detail::MinmaxReduceFunctor(),
1669 MetaInt<actual_dimension-1>());
1670 *minimum = res.first;
1671 *maximum = res.second;
1679 typedef typename NumericTraits<U>::RealPromote R;
1680 triple<R, R, R> res(0.0, 0.0, 0.0);
1683 detail::MeanVarianceReduceFunctor(),
1684 MetaInt<actual_dimension-1>());
1686 *variance = res.third / res.first;
1701 U res = NumericTraits<U>::zero();
1704 detail::SumReduceFunctor(),
1705 MetaInt<actual_dimension-1>());
1734 template <
class U,
class S>
1738 destMultiArrayRange(sums),
1754 U res = NumericTraits<U>::one();
1757 detail::ProdReduceFunctor(),
1758 MetaInt<actual_dimension-1>());
1764 typename NormTraits<MultiArrayView>::SquaredNormType
1767 typedef typename NormTraits<MultiArrayView>::SquaredNormType SquaredNormType;
1768 SquaredNormType res = NumericTraits<SquaredNormType>::zero();
1771 detail::SquaredL2NormReduceFunctor(),
1772 MetaInt<actual_dimension-1>());
1789 typename NormTraits<MultiArrayView>::NormType
1790 norm(
int type = 2,
bool useSquaredNorm =
true)
const;
1814 return iterator(m_ptr, m_shape, m_stride);
1830 return begin().getEndIterator();
1838 return begin().getEndIterator();
1846 traverser ret (m_ptr, m_stride.begin (), m_shape.begin ());
1855 const_traverser ret (m_ptr, m_stride.begin (), m_shape.begin ());
1865 traverser ret (m_ptr, m_stride.begin (), m_shape.begin ());
1866 ret += m_shape [actual_dimension-1];
1876 const_traverser ret (m_ptr, m_stride.begin (), m_shape.begin ());
1877 ret += m_shape [actual_dimension-1];
1887 template <
unsigned int N,
class T,
class Str
ideTag>
1888 MultiArrayView<N, T, StrideTag> &
1893 vigra_precondition(this->
shape() == rhs.shape() ||
m_ptr == 0,
1894 "MultiArrayView::operator=(MultiArrayView const &) size mismatch.");
1902 this->copyImpl(rhs);
1906 template <
unsigned int N,
class T,
class Str
ideTag>
1911 vigra_precondition (shape () == rhs.shape (),
1912 "MultiArrayView::arraysOverlap(): shape mismatch.");
1916 rhs_first_element = rhs.data(),
1917 rhs_last_element = rhs_first_element +
dot(rhs.shape() -
difference_type(1), rhs.stride());
1918 return !(last_element < rhs_first_element || rhs_last_element < first_element);
1921 template <
unsigned int N,
class T,
class Str
ideTag>
1922 template <
class U,
class CN>
1924 MultiArrayView <N, T, StrideTag>::copyImpl(
const MultiArrayView <N, U, CN>& rhs)
1926 if(!arraysOverlap(rhs))
1929 detail::copyMultiArrayData(rhs.traverser_begin(),
shape(),
traverser_begin(), MetaInt<actual_dimension-1>());
1935 MultiArray<N, T> tmp(rhs);
1936 detail::copyMultiArrayData(tmp.traverser_begin(),
shape(),
traverser_begin(), MetaInt<actual_dimension-1>());
1940 #define VIGRA_MULTI_ARRAY_COMPUTED_ASSIGNMENT(name, op) \
1941 template <unsigned int N, class T, class StrideTag> \
1942 template<class U, class C1> \
1943 MultiArrayView<N, T, StrideTag> & \
1944 MultiArrayView <N, T, StrideTag>::operator op(MultiArrayView<N, U, C1> const & rhs) \
1946 vigra_precondition(this->shape() == rhs.shape(), "MultiArrayView::operator" #op "() size mismatch."); \
1947 if(!arraysOverlap(rhs)) \
1949 detail::name##MultiArrayData(rhs.traverser_begin(), shape(), traverser_begin(), MetaInt<actual_dimension-1>()); \
1953 MultiArray<N, T> tmp(rhs); \
1954 detail::name##MultiArrayData(tmp.traverser_begin(), shape(), traverser_begin(), MetaInt<actual_dimension-1>()); \
1959 VIGRA_MULTI_ARRAY_COMPUTED_ASSIGNMENT(copyAdd, +=)
1960 VIGRA_MULTI_ARRAY_COMPUTED_ASSIGNMENT(copySub, -=)
1961 VIGRA_MULTI_ARRAY_COMPUTED_ASSIGNMENT(copyMul, *=)
1962 VIGRA_MULTI_ARRAY_COMPUTED_ASSIGNMENT(copyDiv, /=)
1964 #undef VIGRA_MULTI_ARRAY_COMPUTED_ASSIGNMENT
1966 template <
unsigned int N,
class T,
class Str
ideTag>
1967 template <
class U,
class CN>
1969 MultiArrayView <N, T, StrideTag>::swapDataImpl(MultiArrayView <N, U, CN> rhs)
1971 vigra_precondition (shape () == rhs.shape (),
1972 "MultiArrayView::swapData(): shape mismatch.");
1977 typename MultiArrayView <N, U, CN>::const_pointer
1978 rhs_first_element = rhs.
data(),
1979 rhs_last_element = rhs_first_element +
dot(rhs.shape() -
difference_type(1), rhs.stride());
1980 if(last_element < rhs_first_element || rhs_last_element < first_element)
1983 detail::swapDataImpl(
traverser_begin(),
shape(), rhs.traverser_begin(), MetaInt<actual_dimension-1>());
1989 MultiArray<N, T> tmp(*
this);
1995 template <
unsigned int N,
class T,
class Str
ideTag>
1996 MultiArrayView <N, T, StridedArrayTag>
2000 for (
unsigned int i = 0; i < actual_dimension; ++i)
2007 "MultiArrayView::permuteDimensions(): every dimension must occur exactly once.");
2011 template <
unsigned int N,
class T,
class Str
ideTag>
2016 for(
int k=0; k<(int)N; ++k)
2018 for(
int k=0; k<(int)N-1; ++k)
2021 for(
int j=k+1; j<(int)N; ++j)
2029 std::swap(permutation[k], permutation[smallest]);
2033 for(
unsigned int k=0; k<N; ++k)
2034 ordering[permutation[k]] = k;
2038 template <
unsigned int N,
class T,
class Str
ideTag>
2044 permutation[ordering[k]] = k;
2048 template <
unsigned int N,
class T,
class Str
ideTag>
2054 permutation[ordering[N-1-k]] = k;
2058 template <
unsigned int N,
class T,
class Str
ideTag>
2059 template <
unsigned int M>
2066 static const int NNew = (N-M == 0) ? 1 : N-M;
2070 inner_shape [0] = 1;
2071 inner_stride [0] = 0;
2078 return MultiArrayView <N-M, T, StrideTag> (inner_shape, inner_stride, ptr);
2081 template <
unsigned int N,
class T,
class Str
ideTag>
2082 template <
unsigned int M>
2089 static const int NNew = (N-M == 0) ? 1 : N-M;
2093 outer_shape [0] = 1;
2094 outer_stride [0] = 0;
2102 (outer_shape, outer_stride, ptr);
2105 template <
unsigned int N,
class T,
class Str
ideTag>
2106 template <
unsigned int M>
2107 MultiArrayView <N-1, T,
typename detail::MaybeStrided<StrideTag, M>::type >
2110 static const int NNew = (N-1 == 0) ? 1 : N-1;
2122 shape.
begin () + M);
2125 stride.
begin () + M);
2127 return MultiArrayView <N-1, T,
typename detail::MaybeStrided<StrideTag, M>::type>
2131 template <
unsigned int N,
class T,
class Str
ideTag>
2135 static const int NNew = (N-1 == 0) ? 1 : N-1;
2139 inner_shape [0] = 1;
2140 inner_stride [0] = 0;
2147 return MultiArrayView <N-1, T, StrideTag> (inner_shape, inner_stride,
2151 template <
unsigned int N,
class T,
class Str
ideTag>
2155 static const int NNew = (N-1 == 0) ? 1 : N-1;
2159 outer_shape [0] = 1;
2160 outer_stride [0] = 0;
2171 template <
unsigned int N,
class T,
class Str
ideTag>
2175 vigra_precondition (
2176 n < static_cast <int> (N),
2177 "MultiArrayView <N, T, StrideTag>::bindAt(): dimension out of range.");
2178 static const int NNew = (N-1 == 0) ? 1 : N-1;
2190 shape.
begin () + n);
2193 stride.
begin () + n);
2200 template <
unsigned int N,
class T,
class Str
ideTag>
2204 vigra_precondition(0 <= d && d <= static_cast <difference_type_1> (N),
2205 "MultiArrayView<N, ...>::expandElements(d): 0 <= 'd' <= N required.");
2207 int elementSize = ExpandElementResult<T>::size;
2209 for(
int k=0; k<d; ++k)
2212 newStrides[k] =
m_stride[k]*elementSize;
2215 newShape[d] = elementSize;
2218 for(
int k=d; k<N; ++k)
2221 newStrides[k+1] =
m_stride[k]*elementSize;
2224 typedef typename ExpandElementResult<T>::type U;
2226 newShape, newStrides,
reinterpret_cast<U*
>(
m_ptr));
2229 template <
unsigned int N,
class T,
class Str
ideTag>
2233 vigra_precondition (
2234 0 <= i && i <= static_cast <difference_type_1> (N),
2235 "MultiArrayView <N, T, StrideTag>::insertSingletonDimension(): index out of range.");
2247 template <
unsigned int N,
class T,
class Str
ideTag>
2248 typename NormTraits<MultiArrayView <N, T, StrideTag> >::NormType
2251 typedef typename NormTraits<MultiArrayView>::NormType NormType;
2257 NormType res = NumericTraits<NormType>::zero();
2260 detail::MaxNormReduceFunctor(),
2261 MetaInt<actual_dimension-1>());
2266 NormType res = NumericTraits<NormType>::zero();
2269 detail::L1NormReduceFunctor(),
2270 MetaInt<actual_dimension-1>());
2281 NormType normMax = NumericTraits<NormType>::zero();
2284 detail::MaxNormReduceFunctor(),
2285 MetaInt<actual_dimension-1>());
2286 if(normMax == NumericTraits<NormType>::zero())
2288 NormType res = NumericTraits<NormType>::zero();
2291 detail::WeightedL2NormReduceFunctor<NormType>(1.0/normMax),
2292 MetaInt<actual_dimension-1>());
2293 return sqrt(res)*normMax;
2297 vigra_precondition(
false,
"MultiArrayView::norm(): Unknown norm type.");
2298 return NumericTraits<NormType>::zero();
2309 template <
unsigned int N,
class T,
class Str
ideTag>
2310 inline typename NormTraits<MultiArrayView <N, T, StrideTag> >::SquaredNormType
2316 template <
unsigned int N,
class T,
class Str
ideTag>
2317 inline typename NormTraits<MultiArrayView <N, T, StrideTag> >::NormType
2318 norm(MultiArrayView <N, T, StrideTag>
const & a)
2352 template <
unsigned int N,
class T,
class A >
2405 typedef typename vigra::detail::MultiIteratorChooser <
2411 typedef typename vigra::detail::MultiIteratorChooser <
2425 typedef typename difference_type::value_type diff_zero_t;
2445 template <
class U,
class Str
ideTag>
2452 template <
class U,
class Str
ideTag>
2473 allocator_type
const & alloc = allocator_type());
2478 allocator_type
const & alloc = allocator_type());
2483 allocator_type
const & alloc = allocator_type());
2489 m_alloc (rhs.m_alloc)
2496 template<
class Expression>
2497 MultiArray (multi_math::MultiMathOperand<Expression>
const & rhs,
2503 multi_math::detail::assignOrResize(*
this, rhs);
2508 template <
class U,
class Str
ideTag>
2510 allocator_type
const & alloc = allocator_type());
2520 this->copyOrReshape(rhs);
2529 template <
class U,
class Str
ideTag>
2532 this->copyOrReshape(rhs);
2539 template <
class U,
class Str
ideTag>
2549 template <
class U,
class Str
ideTag>
2559 template <
class U,
class Str
ideTag>
2569 template <
class U,
class Str
ideTag>
2610 template<
class Expression>
2613 multi_math::detail::assignOrResize(*
this, rhs);
2620 template<
class Expression>
2623 multi_math::detail::plusAssignOrResize(*
this, rhs);
2630 template<
class Expression>
2633 multi_math::detail::minusAssignOrResize(*
this, rhs);
2640 template<
class Expression>
2643 multi_math::detail::multiplyAssignOrResize(*
this, rhs);
2650 template<
class Expression>
2653 multi_math::detail::divideAssignOrResize(*
this, rhs);
2680 reshape (shape, T());
2701 return this->
data();
2715 return this->
data();
2733 template <
unsigned int N,
class T,
class A>
2749 template <
unsigned int N,
class T,
class A>
2765 template <
unsigned int N,
class T,
class A>
2781 template <
unsigned int N,
class T,
class A>
2782 template <
class U,
class Str
ideTag>
2786 detail::defaultStride <MultiArrayView<N,T>::actual_dimension>(rhs.shape()),
2790 allocate (this->
m_ptr, rhs);
2793 template <
unsigned int N,
class T,
class A>
2794 template <
class U,
class Str
ideTag>
2798 if (this->
shape() == rhs.shape())
2807 template <
unsigned int N,
class T,
class A>
2815 else if(new_shape == this->
shape())
2817 this->
init(initial);
2821 difference_type new_stride = detail::defaultStride <MultiArrayView<N,T>::actual_dimension> (new_shape);
2824 allocate (new_ptr, new_size, initial);
2826 this->
m_ptr = new_ptr;
2833 template <
unsigned int N,
class T,
class A>
2839 std::swap(this->
m_shape, other.m_shape);
2840 std::swap(this->
m_stride, other.m_stride);
2841 std::swap(this->
m_ptr, other.m_ptr);
2842 std::swap(this->m_alloc, other.m_alloc);
2845 template <
unsigned int N,
class T,
class A>
2849 ptr = m_alloc.
allocate ((
typename A::size_type)s);
2852 for (i = 0; i < s; ++i)
2853 m_alloc.construct (ptr + i, init);
2857 m_alloc.destroy (ptr + j);
2858 m_alloc.deallocate (ptr, (
typename A::size_type)s);
2863 template <
unsigned int N,
class T,
class A>
2868 ptr = m_alloc.
allocate ((
typename A::size_type)s);
2871 for (i = 0; i < s; ++i, ++
init)
2872 m_alloc.construct (ptr + i, *init);
2876 m_alloc.destroy (ptr + j);
2877 m_alloc.deallocate (ptr, (
typename A::size_type)s);
2882 template <
unsigned int N,
class T,
class A>
2883 template <
class U,
class Str
ideTag>
2887 ptr = m_alloc.allocate ((
typename A::size_type)s);
2890 detail::uninitializedCopyMultiArrayData(init.traverser_begin(), init.shape(),
2891 p, m_alloc, MetaInt<actual_dimension-1>());
2894 for (
pointer pp = ptr; pp < p; ++pp)
2895 m_alloc.destroy (pp);
2896 m_alloc.deallocate (ptr, (
typename A::size_type)s);
2901 template <
unsigned int N,
class T,
class A>
2907 m_alloc.destroy (ptr + i);
2908 m_alloc.
deallocate (ptr, (
typename A::size_type)s);
2918 template <
unsigned int N,
class T,
class Str
ideTag>
2919 inline triple<typename MultiArrayView<N,T,StrideTag>::const_traverser,
2924 return triple<typename MultiArrayView<N,T,StrideTag>::const_traverser,
2932 template <
unsigned int N,
class T,
class Str
ideTag,
class Accessor>
2933 inline triple<typename MultiArrayView<N,T,StrideTag>::const_traverser,
2934 typename MultiArrayView<N,T,StrideTag>::difference_type,
2936 srcMultiArrayRange( MultiArrayView<N,T,StrideTag>
const & array, Accessor a )
2938 return triple<typename MultiArrayView<N,T,StrideTag>::const_traverser,
2939 typename MultiArrayView<N,T,StrideTag>::difference_type,
2941 ( array.traverser_begin(),
2946 template <
unsigned int N,
class T,
class Str
ideTag>
2947 inline pair<typename MultiArrayView<N,T,StrideTag>::const_traverser,
2948 typename AccessorTraits<T>::default_const_accessor >
2949 srcMultiArray( MultiArrayView<N,T,StrideTag>
const & array )
2951 return pair<typename MultiArrayView<N,T,StrideTag>::const_traverser,
2952 typename AccessorTraits<T>::default_const_accessor >
2953 ( array.traverser_begin(),
2954 typename AccessorTraits<T>::default_const_accessor() );
2957 template <
unsigned int N,
class T,
class Str
ideTag,
class Accessor>
2958 inline pair<typename MultiArrayView<N,T,StrideTag>::const_traverser,
2960 srcMultiArray( MultiArrayView<N,T,StrideTag>
const & array, Accessor a )
2962 return pair<typename MultiArrayView<N,T,StrideTag>::const_traverser,
2964 ( array.traverser_begin(), a );
2967 template <
unsigned int N,
class T,
class Str
ideTag>
2968 inline triple<typename MultiArrayView<N,T,StrideTag>::traverser,
2969 typename MultiArrayView<N,T,StrideTag>::difference_type,
2970 typename AccessorTraits<T>::default_accessor >
2971 destMultiArrayRange( MultiArrayView<N,T,StrideTag> & array )
2973 return triple<typename MultiArrayView<N,T,StrideTag>::traverser,
2974 typename MultiArrayView<N,T,StrideTag>::difference_type,
2975 typename AccessorTraits<T>::default_accessor >
2976 ( array.traverser_begin(),
2978 typename AccessorTraits<T>::default_accessor() );
2981 template <
unsigned int N,
class T,
class Str
ideTag,
class Accessor>
2982 inline triple<typename MultiArrayView<N,T,StrideTag>::traverser,
2983 typename MultiArrayView<N,T,StrideTag>::difference_type,
2985 destMultiArrayRange( MultiArrayView<N,T,StrideTag> & array, Accessor a )
2987 return triple<typename MultiArrayView<N,T,StrideTag>::traverser,
2988 typename MultiArrayView<N,T,StrideTag>::difference_type,
2990 ( array.traverser_begin(),
2995 template <
unsigned int N,
class T,
class Str
ideTag>
2996 inline pair<typename MultiArrayView<N,T,StrideTag>::traverser,
2997 typename AccessorTraits<T>::default_accessor >
2998 destMultiArray( MultiArrayView<N,T,StrideTag> & array )
3000 return pair<typename MultiArrayView<N,T,StrideTag>::traverser,
3001 typename AccessorTraits<T>::default_accessor >
3002 ( array.traverser_begin(),
3003 typename AccessorTraits<T>::default_accessor() );
3006 template <
unsigned int N,
class T,
class Str
ideTag,
class Accessor>
3007 inline pair<typename MultiArrayView<N,T,StrideTag>::traverser,
3009 destMultiArray( MultiArrayView<N,T,StrideTag> & array, Accessor a )
3011 return pair<typename MultiArrayView<N,T,StrideTag>::traverser,
3013 ( array.traverser_begin(), a );
3018 template <
class PixelType,
class Accessor>
3019 inline triple<ConstStridedImageIterator<PixelType>,
3020 ConstStridedImageIterator<PixelType>, Accessor>
3021 srcImageRange(
const MultiArrayView<2, PixelType, StridedArrayTag> & img, Accessor a)
3023 ConstStridedImageIterator<PixelType>
3024 ul(img.data(), 1, img.stride(0), img.stride(1));
3025 return triple<ConstStridedImageIterator<PixelType>,
3026 ConstStridedImageIterator<PixelType>,
3028 ul, ul + Size2D(img.shape(0), img.shape(1)), a);
3031 template <
class PixelType,
class Accessor>
3032 inline pair<ConstStridedImageIterator<PixelType>, Accessor>
3033 srcImage(
const MultiArrayView<2, PixelType, StridedArrayTag> & img, Accessor a)
3035 ConstStridedImageIterator<PixelType>
3036 ul(img.data(), 1, img.stride(0), img.stride(1));
3037 return pair<ConstStridedImageIterator<PixelType>, Accessor>
3041 template <
class PixelType,
class Accessor>
3042 inline triple<StridedImageIterator<PixelType>,
3043 StridedImageIterator<PixelType>, Accessor>
3044 destImageRange(MultiArrayView<2, PixelType, StridedArrayTag> & img, Accessor a)
3046 StridedImageIterator<PixelType>
3047 ul(img.data(), 1, img.stride(0), img.stride(1));
3048 return triple<StridedImageIterator<PixelType>,
3049 StridedImageIterator<PixelType>,
3051 ul, ul + Size2D(img.shape(0), img.shape(1)), a);
3054 template <
class PixelType,
class Accessor>
3055 inline pair<StridedImageIterator<PixelType>, Accessor>
3056 destImage(MultiArrayView<2, PixelType, StridedArrayTag> & img, Accessor a)
3058 StridedImageIterator<PixelType>
3059 ul(img.data(), 1, img.stride(0), img.stride(1));
3060 return pair<StridedImageIterator<PixelType>, Accessor>
3064 template <
class PixelType,
class Accessor>
3065 inline pair<StridedImageIterator<PixelType>, Accessor>
3066 maskImage(MultiArrayView<2, PixelType, StridedArrayTag> & img, Accessor a)
3068 StridedImageIterator<PixelType>
3069 ul(img.data(), 1, img.stride(0), img.stride(1));
3070 return pair<StridedImageIterator<PixelType>, Accessor>
3076 template <
class PixelType>
3077 inline triple<ConstStridedImageIterator<PixelType>,
3078 ConstStridedImageIterator<PixelType>,
3079 typename AccessorTraits<PixelType>::default_const_accessor>
3080 srcImageRange(MultiArrayView<2, PixelType, StridedArrayTag>
const & img)
3082 ConstStridedImageIterator<PixelType>
3083 ul(img.data(), 1, img.stride(0), img.stride(1));
3084 typedef typename AccessorTraits<PixelType>::default_const_accessor Accessor;
3085 return triple<ConstStridedImageIterator<PixelType>,
3086 ConstStridedImageIterator<PixelType>,
3088 (ul, ul + Size2D(img.shape(0), img.shape(1)), Accessor());
3091 template <
class PixelType>
3092 inline triple<ConstImageIterator<PixelType>,
3093 ConstImageIterator<PixelType>,
3094 typename AccessorTraits<PixelType>::default_const_accessor>
3095 srcImageRange(MultiArrayView<2, PixelType, UnstridedArrayTag>
const & img)
3097 ConstImageIterator<PixelType>
3098 ul(img.data(), img.stride(1));
3099 typedef typename AccessorTraits<PixelType>::default_const_accessor Accessor;
3100 return triple<ConstImageIterator<PixelType>,
3101 ConstImageIterator<PixelType>,
3103 (ul, ul + Size2D(img.shape(0), img.shape(1)), Accessor());
3106 template <
class PixelType>
3107 inline pair< ConstStridedImageIterator<PixelType>,
3108 typename AccessorTraits<PixelType>::default_const_accessor>
3109 srcImage(MultiArrayView<2, PixelType, StridedArrayTag>
const & img)
3111 ConstStridedImageIterator<PixelType>
3112 ul(img.data(), 1, img.stride(0), img.stride(1));
3113 typedef typename AccessorTraits<PixelType>::default_const_accessor Accessor;
3114 return pair<ConstStridedImageIterator<PixelType>,
3119 template <
class PixelType>
3120 inline pair< ConstImageIterator<PixelType>,
3121 typename AccessorTraits<PixelType>::default_const_accessor>
3122 srcImage(MultiArrayView<2, PixelType, UnstridedArrayTag>
const & img)
3124 ConstImageIterator<PixelType>
3125 ul(img.data(), img.stride(1));
3126 typedef typename AccessorTraits<PixelType>::default_const_accessor Accessor;
3127 return pair<ConstImageIterator<PixelType>,
3132 template <
class PixelType>
3133 inline triple< StridedImageIterator<PixelType>,
3134 StridedImageIterator<PixelType>,
3135 typename AccessorTraits<PixelType>::default_accessor>
3136 destImageRange(MultiArrayView<2, PixelType, StridedArrayTag> & img)
3138 StridedImageIterator<PixelType>
3139 ul(img.data(), 1, img.stride(0), img.stride(1));
3140 typedef typename AccessorTraits<PixelType>::default_accessor Accessor;
3141 return triple<StridedImageIterator<PixelType>,
3142 StridedImageIterator<PixelType>,
3144 (ul, ul + Size2D(img.shape(0), img.shape(1)), Accessor());
3147 template <
class PixelType>
3148 inline triple< ImageIterator<PixelType>,
3149 ImageIterator<PixelType>,
3150 typename AccessorTraits<PixelType>::default_accessor>
3151 destImageRange(MultiArrayView<2, PixelType, UnstridedArrayTag> & img)
3153 ImageIterator<PixelType>
3154 ul(img.data(), img.stride(1));
3155 typedef typename AccessorTraits<PixelType>::default_accessor Accessor;
3156 return triple<ImageIterator<PixelType>,
3157 ImageIterator<PixelType>,
3159 (ul, ul + Size2D(img.shape(0), img.shape(1)), Accessor());
3162 template <
class PixelType>
3163 inline pair< StridedImageIterator<PixelType>,
3164 typename AccessorTraits<PixelType>::default_accessor>
3165 destImage(MultiArrayView<2, PixelType, StridedArrayTag> & img)
3167 StridedImageIterator<PixelType>
3168 ul(img.data(), 1, img.stride(0), img.stride(1));
3169 typedef typename AccessorTraits<PixelType>::default_accessor Accessor;
3170 return pair<StridedImageIterator<PixelType>, Accessor>
3174 template <
class PixelType>
3175 inline pair< ImageIterator<PixelType>,
3176 typename AccessorTraits<PixelType>::default_accessor>
3177 destImage(MultiArrayView<2, PixelType, UnstridedArrayTag> & img)
3179 ImageIterator<PixelType> ul(img.data(), img.stride(1));
3180 typedef typename AccessorTraits<PixelType>::default_accessor Accessor;
3181 return pair<ImageIterator<PixelType>, Accessor>(ul, Accessor());
3184 template <
class PixelType>
3185 inline pair< ConstStridedImageIterator<PixelType>,
3186 typename AccessorTraits<PixelType>::default_accessor>
3187 maskImage(MultiArrayView<2, PixelType, StridedArrayTag>
const & img)
3189 ConstStridedImageIterator<PixelType>
3190 ul(img.data(), 1, img.stride(0), img.stride(1));
3191 typedef typename AccessorTraits<PixelType>::default_accessor Accessor;
3192 return pair<ConstStridedImageIterator<PixelType>, Accessor>
3196 template <
class PixelType>
3197 inline pair< ConstImageIterator<PixelType>,
3198 typename AccessorTraits<PixelType>::default_accessor>
3199 maskImage(MultiArrayView<2, PixelType, UnstridedArrayTag>
const & img)
3201 ConstImageIterator<PixelType>
3202 ul(img.data(), img.stride(1));
3203 typedef typename AccessorTraits<PixelType>::default_accessor Accessor;
3204 return pair<ConstImageIterator<PixelType>, Accessor>
3257 BasicImageView <RGBValue<T> >
3260 vigra_precondition (
3261 array.
shape (0) == 3,
"makeRGBImageView(): array.shape(0) must be 3.");
3270 #undef VIGRA_ASSERT_INSIDE
3271 #endif // VIGRA_MULTI_ARRAY_HXX