36 #ifndef VIGRA_SPLINEIMAGEVIEW_HXX
37 #define VIGRA_SPLINEIMAGEVIEW_HXX
39 #include "mathutil.hxx"
40 #include "recursiveconvolution.hxx"
41 #include "splines.hxx"
42 #include "array_vector.hxx"
43 #include "basicimage.hxx"
44 #include "copyimage.hxx"
45 #include "tinyvector.hxx"
46 #include "fixedpoint.hxx"
47 #include "multi_array.hxx"
98 template <
int ORDER,
class VALUETYPE>
101 typedef typename NumericTraits<VALUETYPE>::RealPromote InternalValue;
127 typedef typename InternalTraverser::row_iterator InternalRowIterator;
131 enum { ksize_ = ORDER + 1, kcenter_ = ORDER / 2 };
141 template <
class SrcIterator,
class SrcAccessor>
142 SplineImageView(SrcIterator is, SrcIterator iend, SrcAccessor sa,
bool skipPrefiltering =
false)
143 : w_(iend.x - is.x), h_(iend.y - is.y), w1_(w_-1), h1_(h_-1),
144 x0_(kcenter_), x1_(w_ - kcenter_ - 2), y0_(kcenter_), y1_(h_ - kcenter_ - 2),
149 copyImage(srcIterRange(is, iend, sa), destImage(image_));
150 if(!skipPrefiltering)
161 template <
class SrcIterator,
class SrcAccessor>
162 SplineImageView(triple<SrcIterator, SrcIterator, SrcAccessor> s,
bool skipPrefiltering =
false)
163 : w_(s.second.x - s.first.x), h_(s.second.y - s.first.y), w1_(w_-1), h1_(h_-1),
164 x0_(kcenter_), x1_(w_ - kcenter_ - 2), y0_(kcenter_), y1_(h_ - kcenter_ - 2),
169 copyImage(srcIterRange(s.first, s.second, s.third), destImage(image_));
170 if(!skipPrefiltering)
258 {
return dx(d[0], d[1]); }
264 {
return dy(d[0], d[1]); }
270 {
return dxx(d[0], d[1]); }
276 {
return dxy(d[0], d[1]); }
282 {
return dyy(d[0], d[1]); }
288 {
return dx3(d[0], d[1]); }
294 {
return dy3(d[0], d[1]); }
300 {
return dxxy(d[0], d[1]); }
306 {
return dxyy(d[0], d[1]); }
340 {
return g2(d[0], d[1]); }
346 {
return g2x(d[0], d[1]); }
352 {
return g2y(d[0], d[1]); }
358 {
return g2xx(d[0], d[1]); }
364 {
return g2xy(d[0], d[1]); }
370 {
return g2yy(d[0], d[1]); }
440 template <
class Array>
448 return x >= 0.0 && x <=
width()-1.0;
456 return y >= 0.0 && y <=
height()-1.0;
474 return x < w1_ + x1_ && x > -x1_ && y < h1_ + y1_ && y > -y1_;
484 bool sameFacet(
double x0,
double y0,
double x1,
double y1)
const
490 return x0 == x1 && y0 == y1;
496 void calculateIndices(
double x,
double y)
const;
497 void coefficients(
double t,
double *
const & c)
const;
498 void derivCoefficients(
double t,
unsigned int d,
double *
const & c)
const;
503 double x0_, x1_, y0_, y1_;
506 mutable double x_, y_, u_, v_, kx_[ksize_], ky_[ksize_];
507 mutable int ix_[ksize_], iy_[ksize_];
510 template <
int ORDER,
class VALUETYPE>
511 void SplineImageView<ORDER, VALUETYPE>::init()
513 ArrayVector<double>
const & b = k_.prefilterCoefficients();
515 for(
unsigned int i=0; i<b.size(); ++i)
517 recursiveFilterX(srcImageRange(image_), destImage(image_), b[i], BORDER_TREATMENT_REFLECT);
518 recursiveFilterY(srcImageRange(image_), destImage(image_), b[i], BORDER_TREATMENT_REFLECT);
526 struct SplineImageViewUnrollLoop1
528 template <
class Array>
529 static void exec(
int c0, Array c)
531 SplineImageViewUnrollLoop1<i-1>::exec(c0, c);
537 struct SplineImageViewUnrollLoop1<0>
539 template <
class Array>
540 static void exec(
int c0, Array c)
546 template <
int i,
class ValueType>
547 struct SplineImageViewUnrollLoop2
549 template <
class Array1,
class RowIterator,
class Array2>
551 exec(Array1 k, RowIterator r, Array2 x)
553 return ValueType(k[i] * r[x[i]]) + SplineImageViewUnrollLoop2<i-1, ValueType>::exec(k, r, x);
557 template <
class ValueType>
558 struct SplineImageViewUnrollLoop2<0, ValueType>
560 template <
class Array1,
class RowIterator,
class Array2>
562 exec(Array1 k, RowIterator r, Array2 x)
564 return ValueType(k[0] * r[x[0]]);
570 template <
int ORDER,
class VALUETYPE>
572 SplineImageView<ORDER, VALUETYPE>::calculateIndices(
double x,
double y)
const
574 if(x == x_ && y == y_)
577 if(x > x0_ && x < x1_ && y > y0_ && y < y1_)
579 detail::SplineImageViewUnrollLoop1<ORDER>::exec(
580 (ORDER % 2) ?
int(x - kcenter_) :
int(x + 0.5 - kcenter_), ix_);
581 detail::SplineImageViewUnrollLoop1<ORDER>::exec(
582 (ORDER % 2) ?
int(y - kcenter_) :
int(y + 0.5 - kcenter_), iy_);
584 u_ = x - ix_[kcenter_];
585 v_ = y - iy_[kcenter_];
589 vigra_precondition(isValid(x,y),
590 "SplineImageView::calculateIndices(): coordinates out of range.");
592 int xCenter = (ORDER % 2) ?
594 (int)VIGRA_CSTD::
floor(x + 0.5);
595 int yCenter = (ORDER % 2) ?
597 (int)VIGRA_CSTD::
floor(y + 0.5);
601 for(
int i = 0; i < ksize_; ++i)
602 ix_[i] = w1_ -
vigra::abs(w1_ - xCenter - (i - kcenter_));
606 for(
int i = 0; i < ksize_; ++i)
607 ix_[i] =
vigra::abs(xCenter - (kcenter_ - i));
611 for(
int i = 0; i < ksize_; ++i)
612 iy_[i] = h1_ -
vigra::abs(h1_ - yCenter - (i - kcenter_));
616 for(
int i = 0; i < ksize_; ++i)
617 iy_[i] =
vigra::abs(yCenter - (kcenter_ - i));
626 template <
int ORDER,
class VALUETYPE>
627 void SplineImageView<ORDER, VALUETYPE>::coefficients(
double t,
double *
const & c)
const
630 for(
int i = 0; i<ksize_; ++i)
634 template <
int ORDER,
class VALUETYPE>
635 void SplineImageView<ORDER, VALUETYPE>::derivCoefficients(
double t,
636 unsigned int d,
double *
const & c)
const
639 for(
int i = 0; i<ksize_; ++i)
643 template <
int ORDER,
class VALUETYPE>
644 VALUETYPE SplineImageView<ORDER, VALUETYPE>::convolve()
const
646 typedef typename NumericTraits<VALUETYPE>::RealPromote RealPromote;
649 ky_[0]*detail::SplineImageViewUnrollLoop2<ORDER, RealPromote>::exec(kx_, image_.rowBegin(iy_[0]), ix_));
651 for(
int j=1; j<ksize_; ++j)
654 ky_[j]*detail::SplineImageViewUnrollLoop2<ORDER, RealPromote>::exec(kx_, image_.rowBegin(iy_[j]), ix_));
656 return detail::RequiresExplicitCast<VALUETYPE>::cast(sum);
659 template <
int ORDER,
class VALUETYPE>
660 template <
class Array>
664 typename Spline::WeightMatrix & weights = Spline::weights();
665 double tmp[ksize_][ksize_];
667 calculateIndices(x, y);
668 for(
int j=0; j<ksize_; ++j)
670 for(
int i=0; i<ksize_; ++i)
673 for(
int k=0; k<ksize_; ++k)
675 tmp[i][j] += weights[i][k]*image_(ix_[k], iy_[j]);
679 res.resize(ksize_, ksize_);
680 for(
int j=0; j<ksize_; ++j)
682 for(
int i=0; i<ksize_; ++i)
685 for(
int k=0; k<ksize_; ++k)
687 res(i,j) += weights[j][k]*tmp[i][k];
693 template <
int ORDER,
class VALUETYPE>
696 calculateIndices(x, y);
697 coefficients(u_, kx_);
698 coefficients(v_, ky_);
702 template <
int ORDER,
class VALUETYPE>
704 unsigned int dx,
unsigned int dy)
const
706 calculateIndices(x, y);
707 derivCoefficients(u_, dx, kx_);
708 derivCoefficients(v_, dy, ky_);
712 template <
int ORDER,
class VALUETYPE>
715 return sq(dx(x,y)) +
sq(dy(x,y));
718 template <
int ORDER,
class VALUETYPE>
721 return VALUETYPE(2.0)*(dx(x,y) * dxx(x,y) + dy(x,y) * dxy(x,y));
724 template <
int ORDER,
class VALUETYPE>
727 return VALUETYPE(2.0)*(dx(x,y) * dxy(x,y) + dy(x,y) * dyy(x,y));
730 template <
int ORDER,
class VALUETYPE>
733 return VALUETYPE(2.0)*(
sq(dxx(x,y)) + dx(x,y) * dx3(x,y) +
sq(dxy(x,y)) + dy(x,y) * dxxy(x,y));
736 template <
int ORDER,
class VALUETYPE>
739 return VALUETYPE(2.0)*(
sq(dxy(x,y)) + dx(x,y) * dxyy(x,y) +
sq(dyy(x,y)) + dy(x,y) * dy3(x,y));
742 template <
int ORDER,
class VALUETYPE>
745 return VALUETYPE(2.0)*(dx(x,y) * dxxy(x,y) + dy(x,y) * dxyy(x,y) + dxy(x,y) * (dxx(x,y) + dyy(x,y)));
753 template <
class VALUETYPE,
class INTERNAL_INDEXER>
754 class SplineImageView0Base
756 typedef typename INTERNAL_INDEXER::value_type InternalValue;
758 typedef VALUETYPE value_type;
761 enum StaticOrder { order = 0 };
765 SplineImageView0Base(
unsigned int w,
unsigned int h)
769 SplineImageView0Base(
int w,
int h, INTERNAL_INDEXER i)
770 : w_(w), h_(h), internalIndexer_(i)
773 template <
unsigned IntBits1,
unsigned FractionalBits1,
774 unsigned IntBits2,
unsigned FractionalBits2>
775 value_type unchecked(FixedPoint<IntBits1, FractionalBits1> x,
776 FixedPoint<IntBits2, FractionalBits2> y)
const
781 template <
unsigned IntBits1,
unsigned FractionalBits1,
782 unsigned IntBits2,
unsigned FractionalBits2>
783 value_type unchecked(FixedPoint<IntBits1, FractionalBits1> x,
784 FixedPoint<IntBits2, FractionalBits2> y,
785 unsigned int dx,
unsigned int dy)
const
787 if((dx != 0) || (dy != 0))
788 return NumericTraits<VALUETYPE>::zero();
789 return unchecked(x, y);
792 value_type unchecked(
double x,
double y)
const
794 return internalIndexer_((
int)(x + 0.5), (
int)(y + 0.5));
797 value_type unchecked(
double x,
double y,
unsigned int dx,
unsigned int dy)
const
799 if((dx != 0) || (dy != 0))
800 return NumericTraits<VALUETYPE>::zero();
801 return unchecked(x, y);
804 value_type operator()(
double x,
double y)
const
809 ix = (int)(-x + 0.5);
810 vigra_precondition(ix <= (
int)w_ - 1,
811 "SplineImageView::operator(): coordinates out of range.");
819 vigra_precondition(ix >= 0,
820 "SplineImageView::operator(): coordinates out of range.");
825 iy = (int)(-y + 0.5);
826 vigra_precondition(iy <= (
int)h_ - 1,
827 "SplineImageView::operator(): coordinates out of range.");
835 vigra_precondition(iy >= 0,
836 "SplineImageView::operator(): coordinates out of range.");
839 return internalIndexer_(ix, iy);
842 value_type operator()(
double x,
double y,
unsigned int dx,
unsigned int dy)
const
844 if((dx != 0) || (dy != 0))
845 return NumericTraits<VALUETYPE>::zero();
846 return operator()(x, y);
849 value_type dx(
double x,
double y)
const
850 {
return NumericTraits<VALUETYPE>::zero(); }
852 value_type dy(
double x,
double y)
const
853 {
return NumericTraits<VALUETYPE>::zero(); }
855 value_type dxx(
double x,
double y)
const
856 {
return NumericTraits<VALUETYPE>::zero(); }
858 value_type dxy(
double x,
double y)
const
859 {
return NumericTraits<VALUETYPE>::zero(); }
861 value_type dyy(
double x,
double y)
const
862 {
return NumericTraits<VALUETYPE>::zero(); }
864 value_type dx3(
double x,
double y)
const
865 {
return NumericTraits<VALUETYPE>::zero(); }
867 value_type dy3(
double x,
double y)
const
868 {
return NumericTraits<VALUETYPE>::zero(); }
870 value_type dxxy(
double x,
double y)
const
871 {
return NumericTraits<VALUETYPE>::zero(); }
873 value_type dxyy(
double x,
double y)
const
874 {
return NumericTraits<VALUETYPE>::zero(); }
876 value_type operator()(difference_type
const & d)
const
877 {
return operator()(d[0], d[1]); }
879 value_type operator()(difference_type
const & d,
unsigned int dx,
unsigned int dy)
const
880 {
return operator()(d[0], d[1], dx, dy); }
882 value_type dx(difference_type
const & d)
const
883 {
return NumericTraits<VALUETYPE>::zero(); }
885 value_type dy(difference_type
const & d)
const
886 {
return NumericTraits<VALUETYPE>::zero(); }
888 value_type dxx(difference_type
const & d)
const
889 {
return NumericTraits<VALUETYPE>::zero(); }
891 value_type dxy(difference_type
const & d)
const
892 {
return NumericTraits<VALUETYPE>::zero(); }
894 value_type dyy(difference_type
const & d)
const
895 {
return NumericTraits<VALUETYPE>::zero(); }
897 value_type dx3(difference_type
const & d)
const
898 {
return NumericTraits<VALUETYPE>::zero(); }
900 value_type dy3(difference_type
const & d)
const
901 {
return NumericTraits<VALUETYPE>::zero(); }
903 value_type dxxy(difference_type
const & d)
const
904 {
return NumericTraits<VALUETYPE>::zero(); }
906 value_type dxyy(difference_type
const & d)
const
907 {
return NumericTraits<VALUETYPE>::zero(); }
909 value_type g2(
double x,
double y)
const
910 {
return NumericTraits<VALUETYPE>::zero(); }
912 value_type g2x(
double x,
double y)
const
913 {
return NumericTraits<VALUETYPE>::zero(); }
915 value_type g2y(
double x,
double y)
const
916 {
return NumericTraits<VALUETYPE>::zero(); }
918 value_type g2xx(
double x,
double y)
const
919 {
return NumericTraits<VALUETYPE>::zero(); }
921 value_type g2xy(
double x,
double y)
const
922 {
return NumericTraits<VALUETYPE>::zero(); }
924 value_type g2yy(
double x,
double y)
const
925 {
return NumericTraits<VALUETYPE>::zero(); }
927 value_type g2(difference_type
const & d)
const
928 {
return NumericTraits<VALUETYPE>::zero(); }
930 value_type g2x(difference_type
const & d)
const
931 {
return NumericTraits<VALUETYPE>::zero(); }
933 value_type g2y(difference_type
const & d)
const
934 {
return NumericTraits<VALUETYPE>::zero(); }
936 value_type g2xx(difference_type
const & d)
const
937 {
return NumericTraits<VALUETYPE>::zero(); }
939 value_type g2xy(difference_type
const & d)
const
940 {
return NumericTraits<VALUETYPE>::zero(); }
942 value_type g2yy(difference_type
const & d)
const
943 {
return NumericTraits<VALUETYPE>::zero(); }
945 unsigned int width()
const
948 unsigned int height()
const
951 size_type size()
const
952 {
return size_type(w_, h_); }
954 TinyVector<unsigned int, 2> shape()
const
955 {
return TinyVector<unsigned int, 2>(w_, h_); }
957 template <
class Array>
958 void coefficientArray(
double x,
double y, Array & res)
const
961 res(0, 0) = operator()(x,y);
964 bool isInsideX(
double x)
const
966 return x >= 0.0 && x <= width() - 1.0;
969 bool isInsideY(
double y)
const
971 return y >= 0.0 && y <= height() - 1.0;
974 bool isInside(
double x,
double y)
const
976 return isInsideX(x) && isInsideY(y);
979 bool isValid(
double x,
double y)
const
981 return x < 2.0*w_-2.0 && x > 1.0-w_ && y < 2.0*h_-2.0 && y > 1.0-h_;
984 bool sameFacet(
double x0,
double y0,
double x1,
double y1)
const
990 return x0 == x1 && y0 == y1;
995 INTERNAL_INDEXER internalIndexer_;
1007 template <class VALUETYPE, class INTERNAL_TRAVERSER = typename BasicImage<VALUETYPE>::const_traverser>
1009 :
public SplineImageView0Base<VALUETYPE, INTERNAL_TRAVERSER>
1011 typedef SplineImageView0Base<VALUETYPE, INTERNAL_TRAVERSER> Base;
1013 typedef typename Base::value_type value_type;
1016 enum StaticOrder { order = Base::order };
1031 : Base(iend.x - is.x, iend.y - is.y, is)
1034 SplineImageView0(triple<InternalTraverser, InternalTraverser, InternalAccessor> s)
1035 : Base(s.second.x - s.first.x, s.second.y - s.first.y, s.first)
1039 : Base(iend.x - is.x, iend.y - is.y, is)
1042 SplineImageView0(triple<InternalConstTraverser, InternalConstTraverser, InternalConstAccessor> s)
1043 : Base(s.second.x - s.first.x, s.second.y - s.first.y, s.first)
1046 template<
class T,
class SU>
1051 for(
unsigned int y=0; y<this->height(); ++y)
1052 for(
unsigned int x=0; x<this->width(); ++x)
1053 image_(x,y) = detail::RequiresExplicitCast<VALUETYPE>::cast(i(x,y));
1054 this->internalIndexer_ = image_.
upperLeft();
1057 template <
class SrcIterator,
class SrcAccessor>
1059 : Base(iend.x - is.x, iend.y - is.y),
1062 copyImage(srcIterRange(is, iend, sa), destImage(image_));
1063 this->internalIndexer_ = image_.
upperLeft();
1066 template <
class SrcIterator,
class SrcAccessor>
1068 : Base(s.second.x - s.first.x, s.second.y - s.first.y),
1069 image_(s.second - s.first)
1072 this->internalIndexer_ = image_.
upperLeft();
1082 template <
class VALUETYPE,
class Str
idedOrUnstr
ided>
1084 :
public SplineImageView0Base<VALUETYPE, MultiArrayView<2, VALUETYPE, StridedOrUnstrided> >
1086 typedef SplineImageView0Base<VALUETYPE, MultiArrayView<2, VALUETYPE, StridedOrUnstrided> > Base;
1088 typedef typename Base::value_type value_type;
1089 typedef typename Base::size_type size_type;
1090 typedef typename Base::difference_type difference_type;
1091 enum StaticOrder { order = Base::order };
1092 typedef BasicImage<VALUETYPE> InternalImage;
1095 typedef MultiArrayView<2, VALUETYPE, StridedOrUnstrided> InternalIndexer;
1102 SplineImageView0(InternalIndexer
const & i)
1103 : Base(i.shape(0), i.shape(1), i)
1106 template<
class T,
class SU>
1107 SplineImageView0(MultiArrayView<2, T, SU>
const & i)
1108 : Base(i.shape(0), i.shape(1)),
1109 image_(i.shape(0), i.shape(1))
1111 for(
unsigned int y=0; y<this->height(); ++y)
1112 for(
unsigned int x=0; x<this->width(); ++x)
1113 image_(x,y) = detail::RequiresExplicitCast<VALUETYPE>::cast(i(x,y));
1114 this->internalIndexer_ = InternalIndexer(
typename InternalIndexer::difference_type(this->width(), this->height()),
1118 template <
class SrcIterator,
class SrcAccessor>
1119 SplineImageView0(SrcIterator is, SrcIterator iend, SrcAccessor sa)
1120 : Base(iend.x - is.x, iend.y - is.y),
1123 copyImage(srcIterRange(is, iend, sa), destImage(image_));
1124 this->internalIndexer_ = InternalIndexer(
typename InternalIndexer::difference_type(this->width(), this->height()),
1128 template <
class SrcIterator,
class SrcAccessor>
1129 SplineImageView0(triple<SrcIterator, SrcIterator, SrcAccessor> s)
1130 : Base(s.second.x - s.first.x, s.second.y - s.first.y),
1131 image_(s.second - s.first)
1134 this->internalIndexer_ = InternalIndexer(
typename InternalIndexer::difference_type(this->width(), this->height()),
1138 InternalImage
const & image()
const
1142 InternalImage image_;
1145 template <
class VALUETYPE>
1146 class SplineImageView<0, VALUETYPE>
1147 :
public SplineImageView0<VALUETYPE>
1149 typedef SplineImageView0<VALUETYPE> Base;
1151 typedef typename Base::value_type
value_type;
1152 typedef typename Base::size_type
size_type;
1158 typedef typename Base::InternalTraverser InternalTraverser;
1159 typedef typename Base::InternalAccessor InternalAccessor;
1160 typedef typename Base::InternalConstTraverser InternalConstTraverser;
1161 typedef typename Base::InternalConstAccessor InternalConstAccessor;
1168 SplineImageView(InternalTraverser is, InternalTraverser iend, InternalAccessor sa,
bool =
false)
1169 : Base(is, iend, sa)
1172 SplineImageView(triple<InternalTraverser, InternalTraverser, InternalAccessor> s,
bool =
false)
1176 SplineImageView(InternalConstTraverser is, InternalConstTraverser iend, InternalConstAccessor sa,
bool =
false)
1177 : Base(is, iend, sa)
1180 SplineImageView(triple<InternalConstTraverser, InternalConstTraverser, InternalConstAccessor> s,
bool =
false)
1184 template <
class SrcIterator,
class SrcAccessor>
1185 SplineImageView(SrcIterator is, SrcIterator iend, SrcAccessor sa,
bool =
false)
1186 : Base(is, iend, sa)
1188 copyImage(srcIterRange(is, iend, sa), destImage(this->image_));
1191 template <
class SrcIterator,
class SrcAccessor>
1192 SplineImageView(triple<SrcIterator, SrcIterator, SrcAccessor> s,
bool =
false)
1204 template <
class VALUETYPE,
class INTERNAL_INDEXER>
1205 class SplineImageView1Base
1207 typedef typename INTERNAL_INDEXER::value_type InternalValue;
1209 typedef VALUETYPE value_type;
1210 typedef Size2D size_type;
1211 typedef TinyVector<double, 2> difference_type;
1212 enum StaticOrder { order = 1 };
1216 SplineImageView1Base(
unsigned int w,
unsigned int h)
1220 SplineImageView1Base(
int w,
int h, INTERNAL_INDEXER i)
1221 : w_(w), h_(h), internalIndexer_(i)
1224 template <
unsigned IntBits1,
unsigned FractionalBits1,
1225 unsigned IntBits2,
unsigned FractionalBits2>
1226 value_type unchecked(FixedPoint<IntBits1, FractionalBits1> x,
1227 FixedPoint<IntBits2, FractionalBits2> y)
const
1230 FixedPoint<0, FractionalBits1> tx =
frac(x - FixedPoint<IntBits1, FractionalBits1>(ix));
1231 FixedPoint<0, FractionalBits1> dtx =
dual_frac(tx);
1232 if(ix == (
int)w_ - 1)
1235 tx.value = FixedPoint<0, FractionalBits1>::ONE;
1239 FixedPoint<0, FractionalBits2> ty =
frac(y - FixedPoint<IntBits2, FractionalBits2>(iy));
1240 FixedPoint<0, FractionalBits2> dty =
dual_frac(ty);
1241 if(iy == (
int)h_ - 1)
1244 ty.value = FixedPoint<0, FractionalBits2>::ONE;
1248 dty*(dtx*fixedPoint(internalIndexer_(ix,iy)) +
1249 tx*fixedPoint(internalIndexer_(ix+1,iy))) +
1250 ty *(dtx*fixedPoint(internalIndexer_(ix,iy+1)) +
1251 tx*fixedPoint(internalIndexer_(ix+1,iy+1))));
1254 template <
unsigned IntBits1,
unsigned FractionalBits1,
1255 unsigned IntBits2,
unsigned FractionalBits2>
1256 value_type unchecked(FixedPoint<IntBits1, FractionalBits1> x,
1257 FixedPoint<IntBits2, FractionalBits2> y,
1258 unsigned int dx,
unsigned int dy)
const
1261 FixedPoint<0, FractionalBits1> tx =
frac(x - FixedPoint<IntBits1, FractionalBits1>(ix));
1262 FixedPoint<0, FractionalBits1> dtx =
dual_frac(tx);
1263 if(ix == (
int)w_ - 1)
1266 tx.value = FixedPoint<0, FractionalBits1>::ONE;
1270 FixedPoint<0, FractionalBits2> ty =
frac(y - FixedPoint<IntBits2, FractionalBits2>(iy));
1271 FixedPoint<0, FractionalBits2> dty =
dual_frac(ty);
1272 if(iy == (
int)h_ - 1)
1275 ty.value = FixedPoint<0, FractionalBits2>::ONE;
1285 dty*(dtx*fixedPoint(internalIndexer_(ix,iy)) +
1286 tx*fixedPoint(internalIndexer_(ix+1,iy))) +
1287 ty *(dtx*fixedPoint(internalIndexer_(ix,iy+1)) +
1288 tx*fixedPoint(internalIndexer_(ix+1,iy+1))));
1291 (dtx*fixedPoint(internalIndexer_(ix,iy+1)) + tx*fixedPoint(internalIndexer_(ix+1,iy+1))) -
1292 (dtx*fixedPoint(internalIndexer_(ix,iy)) + tx*fixedPoint(internalIndexer_(ix+1,iy))));
1294 return NumericTraits<VALUETYPE>::zero();
1301 dty*(fixedPoint(internalIndexer_(ix+1,iy)) - fixedPoint(internalIndexer_(ix,iy))) +
1302 ty *(fixedPoint(internalIndexer_(ix+1,iy+1)) - fixedPoint(internalIndexer_(ix,iy+1))));
1304 return detail::RequiresExplicitCast<value_type>::cast(
1305 (internalIndexer_(ix+1,iy+1) - internalIndexer_(ix,iy+1)) -
1306 (internalIndexer_(ix+1,iy) - internalIndexer_(ix,iy)));
1308 return NumericTraits<VALUETYPE>::zero();
1311 return NumericTraits<VALUETYPE>::zero();
1315 value_type unchecked(
double x,
double y)
const
1318 if(ix == (
int)w_ - 1)
1322 if(iy == (
int)h_ - 1)
1325 return NumericTraits<value_type>::fromRealPromote(
1326 (1.0-ty)*((1.0-tx)*internalIndexer_(ix,iy) + tx*internalIndexer_(ix+1,iy)) +
1327 ty *((1.0-tx)*internalIndexer_(ix,iy+1) + tx*internalIndexer_(ix+1,iy+1)));
1330 value_type unchecked(
double x,
double y,
unsigned int dx,
unsigned int dy)
const
1333 if(ix == (
int)w_ - 1)
1337 if(iy == (
int)h_ - 1)
1346 return detail::RequiresExplicitCast<value_type>::cast(
1347 (1.0-ty)*((1.0-tx)*internalIndexer_(ix,iy) + tx*internalIndexer_(ix+1,iy)) +
1348 ty *((1.0-tx)*internalIndexer_(ix,iy+1) + tx*internalIndexer_(ix+1,iy+1)));
1350 return detail::RequiresExplicitCast<value_type>::cast(
1351 ((1.0-tx)*internalIndexer_(ix,iy+1) + tx*internalIndexer_(ix+1,iy+1)) -
1352 ((1.0-tx)*internalIndexer_(ix,iy) + tx*internalIndexer_(ix+1,iy)));
1354 return NumericTraits<VALUETYPE>::zero();
1360 return detail::RequiresExplicitCast<value_type>::cast(
1361 (1.0-ty)*(internalIndexer_(ix+1,iy) - internalIndexer_(ix,iy)) +
1362 ty *(internalIndexer_(ix+1,iy+1) - internalIndexer_(ix,iy+1)));
1364 return detail::RequiresExplicitCast<value_type>::cast(
1365 (internalIndexer_(ix+1,iy+1) - internalIndexer_(ix,iy+1)) -
1366 (internalIndexer_(ix+1,iy) - internalIndexer_(ix,iy)));
1368 return NumericTraits<VALUETYPE>::zero();
1371 return NumericTraits<VALUETYPE>::zero();
1375 value_type operator()(
double x,
double y)
const
1377 return operator()(x, y, 0, 0);
1380 value_type operator()(
double x,
double y,
unsigned int dx,
unsigned int dy)
const
1382 value_type
mul = NumericTraits<value_type>::one();
1386 vigra_precondition(x <= w_ - 1.0,
1387 "SplineImageView::operator(): coordinates out of range.");
1391 else if(x > w_ - 1.0)
1394 vigra_precondition(x >= 0.0,
1395 "SplineImageView::operator(): coordinates out of range.");
1402 vigra_precondition(y <= h_ - 1.0,
1403 "SplineImageView::operator(): coordinates out of range.");
1407 else if(y > h_ - 1.0)
1410 vigra_precondition(y >= 0.0,
1411 "SplineImageView::operator(): coordinates out of range.");
1415 return mul*unchecked(x, y, dx, dy);
1418 value_type dx(
double x,
double y)
const
1419 {
return operator()(x, y, 1, 0); }
1421 value_type dy(
double x,
double y)
const
1422 {
return operator()(x, y, 0, 1); }
1424 value_type dxx(
double x,
double y)
const
1425 {
return NumericTraits<VALUETYPE>::zero(); }
1427 value_type dxy(
double x,
double y)
const
1428 {
return operator()(x, y, 1, 1); }
1430 value_type dyy(
double x,
double y)
const
1431 {
return NumericTraits<VALUETYPE>::zero(); }
1433 value_type dx3(
double x,
double y)
const
1434 {
return NumericTraits<VALUETYPE>::zero(); }
1436 value_type dy3(
double x,
double y)
const
1437 {
return NumericTraits<VALUETYPE>::zero(); }
1439 value_type dxxy(
double x,
double y)
const
1440 {
return NumericTraits<VALUETYPE>::zero(); }
1442 value_type dxyy(
double x,
double y)
const
1443 {
return NumericTraits<VALUETYPE>::zero(); }
1445 value_type operator()(difference_type
const & d)
const
1446 {
return operator()(d[0], d[1]); }
1448 value_type operator()(difference_type
const & d,
unsigned int dx,
unsigned int dy)
const
1449 {
return operator()(d[0], d[1], dx, dy); }
1451 value_type dx(difference_type
const & d)
const
1452 {
return operator()(d[0], d[1], 1, 0); }
1454 value_type dy(difference_type
const & d)
const
1455 {
return operator()(d[0], d[1], 0, 1); }
1457 value_type dxx(difference_type
const & d)
const
1458 {
return NumericTraits<VALUETYPE>::zero(); }
1460 value_type dxy(difference_type
const & d)
const
1461 {
return operator()(d[0], d[1], 1, 1); }
1463 value_type dyy(difference_type
const & d)
const
1464 {
return NumericTraits<VALUETYPE>::zero(); }
1466 value_type dx3(difference_type
const & d)
const
1467 {
return NumericTraits<VALUETYPE>::zero(); }
1469 value_type dy3(difference_type
const & d)
const
1470 {
return NumericTraits<VALUETYPE>::zero(); }
1472 value_type dxxy(difference_type
const & d)
const
1473 {
return NumericTraits<VALUETYPE>::zero(); }
1475 value_type dxyy(difference_type
const & d)
const
1476 {
return NumericTraits<VALUETYPE>::zero(); }
1478 value_type g2(
double x,
double y)
const
1479 {
return sq(dx(x,y)) +
sq(dy(x,y)); }
1481 value_type g2x(
double x,
double y)
const
1482 {
return NumericTraits<VALUETYPE>::zero(); }
1484 value_type g2y(
double x,
double y)
const
1485 {
return NumericTraits<VALUETYPE>::zero(); }
1487 value_type g2xx(
double x,
double y)
const
1488 {
return NumericTraits<VALUETYPE>::zero(); }
1490 value_type g2xy(
double x,
double y)
const
1491 {
return NumericTraits<VALUETYPE>::zero(); }
1493 value_type g2yy(
double x,
double y)
const
1494 {
return NumericTraits<VALUETYPE>::zero(); }
1496 value_type g2(difference_type
const & d)
const
1497 {
return g2(d[0], d[1]); }
1499 value_type g2x(difference_type
const & d)
const
1500 {
return NumericTraits<VALUETYPE>::zero(); }
1502 value_type g2y(difference_type
const & d)
const
1503 {
return NumericTraits<VALUETYPE>::zero(); }
1505 value_type g2xx(difference_type
const & d)
const
1506 {
return NumericTraits<VALUETYPE>::zero(); }
1508 value_type g2xy(difference_type
const & d)
const
1509 {
return NumericTraits<VALUETYPE>::zero(); }
1511 value_type g2yy(difference_type
const & d)
const
1512 {
return NumericTraits<VALUETYPE>::zero(); }
1514 unsigned int width()
const
1517 unsigned int height()
const
1520 size_type size()
const
1521 {
return size_type(w_, h_); }
1523 TinyVector<unsigned int, 2> shape()
const
1524 {
return TinyVector<unsigned int, 2>(w_, h_); }
1526 template <
class Array>
1527 void coefficientArray(
double x,
double y, Array & res)
const;
1529 void calculateIndices(
double x,
double y,
int & ix,
int & iy,
int & ix1,
int & iy1)
const;
1531 bool isInsideX(
double x)
const
1533 return x >= 0.0 && x <= width() - 1.0;
1536 bool isInsideY(
double y)
const
1538 return y >= 0.0 && y <= height() - 1.0;
1541 bool isInside(
double x,
double y)
const
1543 return isInsideX(x) && isInsideY(y);
1546 bool isValid(
double x,
double y)
const
1548 return x < 2.0*w_-2.0 && x > 1.0-w_ && y < 2.0*h_-2.0 && y > 1.0-h_;
1551 bool sameFacet(
double x0,
double y0,
double x1,
double y1)
const
1557 return x0 == x1 && y0 == y1;
1561 unsigned int w_, h_;
1562 INTERNAL_INDEXER internalIndexer_;
1565 template <
class VALUETYPE,
class INTERNAL_INDEXER>
1566 template <
class Array>
1567 void SplineImageView1Base<VALUETYPE, INTERNAL_INDEXER>::coefficientArray(
double x,
double y, Array & res)
const
1569 int ix, iy, ix1, iy1;
1570 calculateIndices(x, y, ix, iy, ix1, iy1);
1572 res(0,0) = internalIndexer_(ix,iy);
1573 res(1,0) = internalIndexer_(ix1,iy) - internalIndexer_(ix,iy);
1574 res(0,1) = internalIndexer_(ix,iy1) - internalIndexer_(ix,iy);
1575 res(1,1) = internalIndexer_(ix,iy) - internalIndexer_(ix1,iy) -
1576 internalIndexer_(ix,iy1) + internalIndexer_(ix1,iy1);
1579 template <
class VALUETYPE,
class INTERNAL_INDEXER>
1580 void SplineImageView1Base<VALUETYPE, INTERNAL_INDEXER>::calculateIndices(
double x,
double y,
int & ix,
int & iy,
int & ix1,
int & iy1)
const
1585 vigra_precondition(x <= w_ - 1.0,
1586 "SplineImageView::calculateIndices(): coordinates out of range.");
1590 else if(x >= w_ - 1.0)
1593 vigra_precondition(x > 0.0,
1594 "SplineImageView::calculateIndices(): coordinates out of range.");
1606 vigra_precondition(y <= h_ - 1.0,
1607 "SplineImageView::calculateIndices(): coordinates out of range.");
1611 else if(y >= h_ - 1.0)
1614 vigra_precondition(y > 0.0,
1615 "SplineImageView::calculateIndices(): coordinates out of range.");
1641 template <class VALUETYPE, class INTERNAL_TRAVERSER = typename BasicImage<VALUETYPE>::const_traverser>
1643 :
public SplineImageView1Base<VALUETYPE, INTERNAL_TRAVERSER>
1645 typedef SplineImageView1Base<VALUETYPE, INTERNAL_TRAVERSER> Base;
1647 typedef typename Base::value_type value_type;
1650 enum StaticOrder { order = Base::order };
1665 : Base(iend.x - is.x, iend.y - is.y, is)
1668 SplineImageView1(triple<InternalTraverser, InternalTraverser, InternalAccessor> s)
1669 : Base(s.second.x - s.first.x, s.second.y - s.first.y, s.first)
1673 : Base(iend.x - is.x, iend.y - is.y, is)
1676 SplineImageView1(triple<InternalConstTraverser, InternalConstTraverser, InternalConstAccessor> s)
1677 : Base(s.second.x - s.first.x, s.second.y - s.first.y, s.first)
1680 template<
class T,
class SU>
1685 for(
unsigned int y=0; y<this->height(); ++y)
1686 for(
unsigned int x=0; x<this->width(); ++x)
1687 image_(x,y) = detail::RequiresExplicitCast<VALUETYPE>::cast(i(x,y));
1688 this->internalIndexer_ = image_.
upperLeft();
1691 template <
class SrcIterator,
class SrcAccessor>
1693 : Base(iend.x - is.x, iend.y - is.y),
1696 copyImage(srcIterRange(is, iend, sa), destImage(image_));
1697 this->internalIndexer_ = image_.
upperLeft();
1700 template <
class SrcIterator,
class SrcAccessor>
1702 : Base(s.second.x - s.first.x, s.second.y - s.first.y),
1703 image_(s.second - s.first)
1706 this->internalIndexer_ = image_.
upperLeft();
1716 template <
class VALUETYPE,
class Str
idedOrUnstr
ided>
1718 :
public SplineImageView1Base<VALUETYPE, MultiArrayView<2, VALUETYPE, StridedOrUnstrided> >
1720 typedef SplineImageView1Base<VALUETYPE, MultiArrayView<2, VALUETYPE, StridedOrUnstrided> > Base;
1722 typedef typename Base::value_type value_type;
1723 typedef typename Base::size_type size_type;
1724 typedef typename Base::difference_type difference_type;
1725 enum StaticOrder { order = Base::order };
1726 typedef BasicImage<VALUETYPE> InternalImage;
1729 typedef MultiArrayView<2, VALUETYPE, StridedOrUnstrided> InternalIndexer;
1736 SplineImageView1(InternalIndexer
const & i)
1737 : Base(i.shape(0), i.shape(1), i)
1740 template<
class T,
class SU>
1741 SplineImageView1(MultiArrayView<2, T, SU>
const & i)
1742 : Base(i.shape(0), i.shape(1)),
1743 image_(i.shape(0), i.shape(1))
1745 for(
unsigned int y=0; y<this->height(); ++y)
1746 for(
unsigned int x=0; x<this->width(); ++x)
1747 image_(x,y) = detail::RequiresExplicitCast<VALUETYPE>::cast(i(x,y));
1748 this->internalIndexer_ = InternalIndexer(
typename InternalIndexer::difference_type(this->width(), this->height()),
1752 template <
class SrcIterator,
class SrcAccessor>
1753 SplineImageView1(SrcIterator is, SrcIterator iend, SrcAccessor sa)
1754 : Base(iend.x - is.x, iend.y - is.y),
1757 copyImage(srcIterRange(is, iend, sa), destImage(image_));
1758 this->internalIndexer_ = InternalIndexer(
typename InternalIndexer::difference_type(this->width(), this->height()),
1762 template <
class SrcIterator,
class SrcAccessor>
1763 SplineImageView1(triple<SrcIterator, SrcIterator, SrcAccessor> s)
1764 : Base(s.second.x - s.first.x, s.second.y - s.first.y),
1765 image_(s.second - s.first)
1768 this->internalIndexer_ = InternalIndexer(
typename InternalIndexer::difference_type(this->width(), this->height()),
1772 InternalImage
const & image()
const
1776 InternalImage image_;
1779 template <
class VALUETYPE>
1780 class SplineImageView<1, VALUETYPE>
1781 :
public SplineImageView1<VALUETYPE>
1783 typedef SplineImageView1<VALUETYPE> Base;
1785 typedef typename Base::value_type
value_type;
1786 typedef typename Base::size_type
size_type;
1792 typedef typename Base::InternalTraverser InternalTraverser;
1793 typedef typename Base::InternalAccessor InternalAccessor;
1794 typedef typename Base::InternalConstTraverser InternalConstTraverser;
1795 typedef typename Base::InternalConstAccessor InternalConstAccessor;
1802 SplineImageView(InternalTraverser is, InternalTraverser iend, InternalAccessor sa,
bool =
false)
1803 : Base(is, iend, sa)
1806 SplineImageView(triple<InternalTraverser, InternalTraverser, InternalAccessor> s,
bool =
false)
1810 SplineImageView(InternalConstTraverser is, InternalConstTraverser iend, InternalConstAccessor sa,
bool =
false)
1811 : Base(is, iend, sa)
1814 SplineImageView(triple<InternalConstTraverser, InternalConstTraverser, InternalConstAccessor> s,
bool =
false)
1818 template <
class SrcIterator,
class SrcAccessor>
1819 SplineImageView(SrcIterator is, SrcIterator iend, SrcAccessor sa,
bool =
false)
1820 : Base(is, iend, sa)
1822 copyImage(srcIterRange(is, iend, sa), destImage(this->image_));
1825 template <
class SrcIterator,
class SrcAccessor>
1826 SplineImageView(triple<SrcIterator, SrcIterator, SrcAccessor> s,
bool =
false)