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

vigra/splineimageview.hxx VIGRA

00001 /************************************************************************/
00002 /*                                                                      */
00003 /*               Copyright 1998-2004 by 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 #ifndef VIGRA_SPLINEIMAGEVIEW_HXX
00037 #define VIGRA_SPLINEIMAGEVIEW_HXX
00038 
00039 #include "mathutil.hxx"
00040 #include "recursiveconvolution.hxx"
00041 #include "splines.hxx"
00042 #include "array_vector.hxx"
00043 #include "basicimage.hxx"
00044 #include "copyimage.hxx"
00045 #include "tinyvector.hxx"
00046 #include "fixedpoint.hxx"
00047 #include "multi_array.hxx"
00048 
00049 namespace vigra {
00050 
00051 /********************************************************/
00052 /*                                                      */
00053 /*                    SplineImageView                   */
00054 /*                                                      */
00055 /********************************************************/
00056 /** \brief Create a continuous view onto a discrete image using splines.
00057 
00058     This class is very useful if image values or derivatives at arbitrary
00059     real-valued coordinates are needed. Access at such coordinates is implemented by
00060     interpolating the given discrete image values with a spline of the
00061     specified <tt>ORDER</TT>. Continuous derivatives are available up to
00062     degree <tt>ORDER-1</TT>. If the requested coordinates are near the image border,
00063     reflective boundary conditions are applied. In principle, this class can also be used
00064     for image resizing, but here the functions from the <tt>resize...</tt> family are
00065     more efficient, since they exploit the regularity of the sampling grid.
00066 
00067     The <tt>SplineImageView</tt> template is explicitly specialized to make it as efficient as possible.
00068     In particular, unnecessary copying of the image is avoided when the iterators passed
00069     in the constructor originate from a \ref vigra::BasicImage. In addition, these specializations
00070     provide function <tt>unchecked(...)</tt> that do not perform bounds checking. If the original image
00071     is not a variant of \ref vigra::BasicImage, one can customize the internal representation by
00072     using \ref vigra::SplineImageView0 or \ref vigra::SplineImageView1.
00073 
00074     <b>Usage:</b>
00075 
00076     <b>\#include</b> <vigra/splineimageview.hxx><br>
00077     Namespace: vigra
00078 
00079     \code
00080     BImage img(w,h);
00081     ... // fill img
00082 
00083     // construct spline view for quadratic interpolation
00084     SplineImageView<2, double> spi2(srcImageRange(img));
00085 
00086     double x = ..., y = ...;
00087     double v2 = spi2(x, y);
00088 
00089     // construct spline view for linear interpolation
00090     SplineImageView<1, UInt32> spi1(srcImageRange(img));
00091 
00092     UInt32 v1 = spi1(x, y);
00093 
00094     FixedPoint<16, 15> fx(...), fy(...);
00095     UInt32 vf = spi1.unchecked(fx, fy); // caller is sure that (fx, fy) are valid coordinates
00096     \endcode
00097 */
00098 template <int ORDER, class VALUETYPE>
00099 class SplineImageView
00100 {
00101     typedef typename NumericTraits<VALUETYPE>::RealPromote InternalValue;
00102 
00103   public:
00104 
00105         /** The view's value type (return type of access and derivative functions).
00106         */
00107     typedef VALUETYPE value_type;
00108 
00109         /** The view's size type.
00110         */
00111     typedef Size2D size_type;
00112 
00113         /** The view's difference type.
00114         */
00115     typedef TinyVector<double, 2> difference_type;
00116 
00117         /** The order of the spline used.
00118         */
00119     enum StaticOrder { order = ORDER };
00120 
00121         /** The type of the internal image holding the spline coefficients.
00122         */
00123     typedef BasicImage<InternalValue> InternalImage;
00124 
00125   private:
00126     typedef typename InternalImage::traverser InternalTraverser;
00127     typedef typename InternalTraverser::row_iterator InternalRowIterator;
00128     typedef typename InternalTraverser::column_iterator InternalColumnIterator;
00129     typedef BSpline<ORDER, double> Spline;
00130 
00131     enum { ksize_ = ORDER + 1, kcenter_ = ORDER / 2 };
00132 
00133   public:
00134         /** Construct SplineImageView for the given image.
00135 
00136             If <tt>skipPrefiltering = true</tt> (default: <tt>false</tt>), the recursive
00137             prefilter of the cardinal spline function is not applied, resulting
00138             in an approximating (smoothing) rather than interpolating spline. This is
00139             especially useful if customized prefilters are to be applied.
00140         */
00141     template <class SrcIterator, class SrcAccessor>
00142     SplineImageView(SrcIterator is, SrcIterator iend, SrcAccessor sa, bool skipPrefiltering = false)
00143     : w_(iend.x - is.x), h_(iend.y - is.y), w1_(w_-1), h1_(h_-1),
00144       x0_(kcenter_), x1_(w_ - kcenter_ - 2), y0_(kcenter_), y1_(h_ - kcenter_ - 2),
00145       image_(w_, h_),
00146       x_(-1.0), y_(-1.0),
00147       u_(-1.0), v_(-1.0)
00148     {
00149         copyImage(srcIterRange(is, iend, sa), destImage(image_));
00150         if(!skipPrefiltering)
00151             init();
00152     }
00153 
00154         /** Construct SplineImageView for the given image.
00155 
00156             If <tt>skipPrefiltering = true</tt> (default: <tt>false</tt>), the recursive
00157             prefilter of the cardinal spline function is not applied, resulting
00158             in an approximating (smoothing) rather than interpolating spline. This is
00159             especially useful if customized prefilters are to be applied.
00160         */
00161     template <class SrcIterator, class SrcAccessor>
00162     SplineImageView(triple<SrcIterator, SrcIterator, SrcAccessor> s, bool skipPrefiltering = false)
00163     : w_(s.second.x - s.first.x), h_(s.second.y - s.first.y), w1_(w_-1), h1_(h_-1),
00164       x0_(kcenter_), x1_(w_ - kcenter_ - 2), y0_(kcenter_), y1_(h_ - kcenter_ - 2),
00165       image_(w_, h_),
00166       x_(-1.0), y_(-1.0),
00167       u_(-1.0), v_(-1.0)
00168     {
00169         copyImage(srcIterRange(s.first, s.second, s.third), destImage(image_));
00170         if(!skipPrefiltering)
00171             init();
00172     }
00173 
00174         /** Access interpolated function at real-valued coordinate <tt>(x, y)</tt>.
00175             If <tt>(x, y)</tt> is near the image border or outside the image, the value
00176             is calculated with reflective boundary conditions. An exception is thrown if the
00177             coordinate is outside the first reflection.
00178         */
00179     value_type operator()(double x, double y) const;
00180 
00181         /** Access derivative of order <tt>(dx, dy)</tt> at real-valued coordinate <tt>(x, y)</tt>.
00182             If <tt>(x, y)</tt> is near the image border or outside the image, the value
00183             is calculated with reflective boundary conditions. An exception is thrown if the
00184             coordinate is outside the first reflection.
00185         */
00186     value_type operator()(double x, double y, unsigned int dx, unsigned int dy) const;
00187 
00188         /** Access 1st derivative in x-direction at real-valued coordinate <tt>(x, y)</tt>.
00189             Equivalent to <tt>splineView(x, y, 1, 0)</tt>.
00190         */
00191     value_type dx(double x, double y) const
00192         { return operator()(x, y, 1, 0); }
00193 
00194         /** Access 1st derivative in y-direction at real-valued coordinate <tt>(x, y)</tt>.
00195             Equivalent to <tt>splineView(x, y, 0, 1)</tt>.
00196         */
00197     value_type dy(double x, double y) const
00198         { return operator()(x, y, 0, 1); }
00199 
00200         /** Access 2nd derivative in x-direction at real-valued coordinate <tt>(x, y)</tt>.
00201             Equivalent to <tt>splineView(x, y, 2, 0)</tt>.
00202         */
00203     value_type dxx(double x, double y) const
00204         { return operator()(x, y, 2, 0); }
00205 
00206         /** Access mixed 2nd derivative at real-valued coordinate <tt>(x, y)</tt>.
00207             Equivalent to <tt>splineView(x, y, 1, 1)</tt>.
00208         */
00209     value_type dxy(double x, double y) const
00210         { return operator()(x, y, 1, 1); }
00211 
00212         /** Access 2nd derivative in y-direction at real-valued coordinate <tt>(x, y)</tt>.
00213             Equivalent to <tt>splineView(x, y, 0, 2)</tt>.
00214         */
00215     value_type dyy(double x, double y) const
00216         { return operator()(x, y, 0, 2); }
00217 
00218         /** Access 3rd derivative in x-direction at real-valued coordinate <tt>(x, y)</tt>.
00219             Equivalent to <tt>splineView(x, y, 3, 0)</tt>.
00220         */
00221     value_type dx3(double x, double y) const
00222         { return operator()(x, y, 3, 0); }
00223 
00224         /** Access 3rd derivative in y-direction at real-valued coordinate <tt>(x, y)</tt>.
00225             Equivalent to <tt>splineView(x, y, 0, 3)</tt>.
00226         */
00227     value_type dy3(double x, double y) const
00228         { return operator()(x, y, 0, 3); }
00229 
00230         /** Access mixed 3rd derivative dxxy at real-valued coordinate <tt>(x, y)</tt>.
00231             Equivalent to <tt>splineView(x, y, 2, 1)</tt>.
00232         */
00233     value_type dxxy(double x, double y) const
00234         { return operator()(x, y, 2, 1); }
00235 
00236         /** Access mixed 3rd derivative dxyy at real-valued coordinate <tt>(x, y)</tt>.
00237             Equivalent to <tt>splineView(x, y, 1, 2)</tt>.
00238         */
00239     value_type dxyy(double x, double y) const
00240         { return operator()(x, y, 1, 2); }
00241 
00242         /** Access interpolated function at real-valued coordinate <tt>d</tt>.
00243             Equivalent to <tt>splineView(d[0], d[1])</tt>.
00244         */
00245     value_type operator()(difference_type const & d) const
00246         { return operator()(d[0], d[1]); }
00247 
00248         /** Access derivative of order <tt>(dx, dy)</tt> at real-valued coordinate <tt>d</tt>.
00249             Equivalent to <tt>splineView(d[0], d[1], dx, dy)</tt>.
00250         */
00251     value_type operator()(difference_type const & d, unsigned int dx, unsigned int dy) const
00252         { return operator()(d[0], d[1], dx, dy); }
00253 
00254         /** Access 1st derivative in x-direction at real-valued coordinate <tt>d</tt>.
00255             Equivalent to <tt>splineView.dx(d[0], d[1])</tt>.
00256         */
00257     value_type dx(difference_type const & d) const
00258         { return dx(d[0], d[1]); }
00259 
00260         /** Access 1st derivative in y-direction at real-valued coordinate <tt>d</tt>.
00261             Equivalent to <tt>splineView.dy(d[0], d[1])</tt>.
00262         */
00263     value_type dy(difference_type const & d) const
00264         { return dy(d[0], d[1]); }
00265 
00266         /** Access 2nd derivative in x-direction at real-valued coordinate <tt>d</tt>.
00267             Equivalent to <tt>splineView.dxx(d[0], d[1])</tt>.
00268         */
00269     value_type dxx(difference_type const & d) const
00270         { return dxx(d[0], d[1]); }
00271 
00272         /** Access mixed 2nd derivative at real-valued coordinate <tt>d</tt>.
00273             Equivalent to <tt>splineView.dxy(d[0], d[1])</tt>.
00274         */
00275     value_type dxy(difference_type const & d) const
00276         { return dxy(d[0], d[1]); }
00277 
00278         /** Access 2nd derivative in y-direction at real-valued coordinate <tt>d</tt>.
00279             Equivalent to <tt>splineView.dyy(d[0], d[1])</tt>.
00280         */
00281     value_type dyy(difference_type const & d) const
00282         { return dyy(d[0], d[1]); }
00283 
00284         /** Access 3rd derivative in x-direction at real-valued coordinate <tt>d</tt>.
00285             Equivalent to <tt>splineView.dx3(d[0], d[1])</tt>.
00286         */
00287     value_type dx3(difference_type const & d) const
00288         { return dx3(d[0], d[1]); }
00289 
00290         /** Access 3rd derivative in y-direction at real-valued coordinate <tt>d</tt>.
00291             Equivalent to <tt>splineView.dy3(d[0], d[1])</tt>.
00292         */
00293     value_type dy3(difference_type const & d) const
00294         { return dy3(d[0], d[1]); }
00295 
00296         /** Access mixed 3rd derivative dxxy at real-valued coordinate <tt>d</tt>.
00297             Equivalent to <tt>splineView.dxxy(d[0], d[1])</tt>.
00298         */
00299     value_type dxxy(difference_type const & d) const
00300         { return dxxy(d[0], d[1]); }
00301 
00302         /** Access mixed 3rd derivative dxyy at real-valued coordinate <tt>d</tt>.
00303             Equivalent to <tt>splineView.dxyy(d[0], d[1])</tt>.
00304         */
00305     value_type dxyy(difference_type const & d) const
00306         { return dxyy(d[0], d[1]); }
00307 
00308         /** Access gradient squared magnitude at real-valued coordinate <tt>(x, y)</tt>.
00309         */
00310     value_type g2(double x, double y) const;
00311 
00312         /** Access 1st derivative in x-direction of gradient squared magnitude
00313             at real-valued coordinate <tt>(x, y)</tt>.
00314         */
00315     value_type g2x(double x, double y) const;
00316 
00317         /** Access 1st derivative in y-direction of gradient squared magnitude
00318             at real-valued coordinate <tt>(x, y)</tt>.
00319         */
00320     value_type g2y(double x, double y) const;
00321 
00322         /** Access 2nd derivative in x-direction of gradient squared magnitude
00323             at real-valued coordinate <tt>(x, y)</tt>.
00324         */
00325     value_type g2xx(double x, double y) const;
00326 
00327         /** Access mixed 2nd derivative of gradient squared magnitude
00328             at real-valued coordinate <tt>(x, y)</tt>.
00329         */
00330     value_type g2xy(double x, double y) const;
00331 
00332         /** Access 2nd derivative in y-direction of gradient squared magnitude
00333             at real-valued coordinate <tt>(x, y)</tt>.
00334         */
00335     value_type g2yy(double x, double y) const;
00336 
00337         /** Access gradient squared magnitude at real-valued coordinate <tt>d</tt>.
00338         */
00339     value_type g2(difference_type const & d) const
00340         { return g2(d[0], d[1]); }
00341 
00342         /** Access 1st derivative in x-direction of gradient squared magnitude
00343             at real-valued coordinate <tt>d</tt>.
00344         */
00345     value_type g2x(difference_type const & d) const
00346         { return g2x(d[0], d[1]); }
00347 
00348         /** Access 1st derivative in y-direction of gradient squared magnitude
00349             at real-valued coordinate <tt>d</tt>.
00350         */
00351     value_type g2y(difference_type const & d) const
00352         { return g2y(d[0], d[1]); }
00353 
00354         /** Access 2nd derivative in x-direction of gradient squared magnitude
00355             at real-valued coordinate <tt>d</tt>.
00356         */
00357     value_type g2xx(difference_type const & d) const
00358         { return g2xx(d[0], d[1]); }
00359 
00360         /** Access mixed 2nd derivative of gradient squared magnitude
00361             at real-valued coordinate <tt>d</tt>.
00362         */
00363     value_type g2xy(difference_type const & d) const
00364         { return g2xy(d[0], d[1]); }
00365 
00366         /** Access 2nd derivative in y-direction of gradient squared magnitude
00367             at real-valued coordinate <tt>d</tt>.
00368         */
00369     value_type g2yy(difference_type const & d) const
00370         { return g2yy(d[0], d[1]); }
00371 
00372         /** The width of the image.
00373             <tt>0 <= x <= width()-1</tt> is required for all access functions.
00374         */
00375     unsigned int width() const
00376         { return w_; }
00377 
00378         /** The height of the image.
00379             <tt>0 <= y <= height()-1</tt> is required for all access functions.
00380         */
00381     unsigned int height() const
00382         { return h_; }
00383 
00384         /** The size of the image.
00385             <tt>0 <= x <= size().x-1</tt> and <tt>0 <= y <= size().y-1</tt>
00386             are required for all access functions.
00387         */
00388     size_type size() const
00389         { return size_type(w_, h_); }
00390 
00391         /** The shape of the image.
00392             Same as size(), except for the return type.
00393         */
00394     TinyVector<unsigned int, 2> shape() const
00395         { return TinyVector<unsigned int, 2>(w_, h_); }
00396 
00397         /** The internal image holding the spline coefficients.
00398         */
00399     InternalImage const & image() const
00400     {
00401         return image_;
00402     }
00403 
00404         /** Get the array of polynomial coefficients for the facet containing
00405             the point <tt>(x, y)</tt>. The array <tt>res</tt> will be resized to
00406             dimension <tt>(ORDER+1)x(ORDER+1)</tt>. From these coefficients, the
00407             value of the interpolated function can be calculated by the following
00408             algorithm
00409 
00410             \code
00411             SplineImageView<ORDER, float> view(...);
00412             double x = ..., y = ...;
00413             double dx, dy;
00414 
00415             // calculate the local facet coordinates of x and y
00416             if(ORDER % 2)
00417             {
00418                 // odd order => facet coordinates between 0 and 1
00419                 dx = x - floor(x);
00420                 dy = y - floor(y);
00421             }
00422             else
00423             {
00424                 // even order => facet coordinates between -0.5 and 0.5
00425                 dx = x - floor(x + 0.5);
00426                 dy = y - floor(y + 0.5);
00427             }
00428 
00429             BasicImage<float> coefficients;
00430             view.coefficientArray(x, y, coefficients);
00431 
00432             float f_x_y = 0.0;
00433             for(int ny = 0; ny < ORDER + 1; ++ny)
00434                 for(int nx = 0; nx < ORDER + 1; ++nx)
00435                     f_x_y += pow(dx, nx) * pow(dy, ny) * coefficients(nx, ny);
00436 
00437             assert(abs(f_x_y - view(x, y)) < 1e-6);
00438             \endcode
00439         */
00440     template <class Array>
00441     void coefficientArray(double x, double y, Array & res) const;
00442 
00443         /** Check if x is in the original image range.
00444             Equivalent to <tt>0 <= x <= width()-1</tt>.
00445         */
00446     bool isInsideX(double x) const
00447     {
00448         return x >= 0.0 && x <= width()-1.0;
00449     }
00450 
00451         /** Check if y is in the original image range.
00452             Equivalent to <tt>0 <= y <= height()-1</tt>.
00453         */
00454     bool isInsideY(double y) const
00455     {
00456         return y >= 0.0 && y <= height()-1.0;
00457     }
00458 
00459         /** Check if x and y are in the original image range.
00460             Equivalent to <tt>0 <= x <= width()-1</tt> and <tt>0 <= y <= height()-1</tt>.
00461         */
00462     bool isInside(double x, double y) const
00463     {
00464         return isInsideX(x) && isInsideY(y);
00465     }
00466 
00467         /** Check if x and y are in the valid range. Points outside the original image range are computed
00468             by reflective boundary conditions, but only within the first reflection.
00469             Equivalent to <tt>-width() + ORDER/2 + 2 < x < 2*width() - ORDER/2 - 2</tt> and
00470             <tt>-height() + ORDER/2 + 2 < y < 2*height() - ORDER/2 - 2</tt>.
00471         */
00472     bool isValid(double x, double y) const
00473     {
00474         return x < w1_ + x1_ && x > -x1_ && y < h1_ + y1_ && y > -y1_;
00475     }
00476 
00477         /** Check whether the points <tt>(x0, y0)</tt> and <tt>(x1, y1)</tt> are in
00478             the same spline facet. For odd order splines, facets span the range
00479             <tt>(floor(x), floor(x)+1) x (floor(y), floor(y)+1)</tt> (i.e. we have
00480             integer facet borders), whereas even order splines have facet between
00481             half integer values
00482             <tt>(floor(x)-0.5, floor(x)+0.5) x (floor(y)-0.5, floor(y)+0.5)</tt>.
00483         */
00484     bool sameFacet(double x0, double y0, double x1, double y1) const
00485     {
00486          x0 = VIGRA_CSTD::floor((ORDER % 2) ? x0 : x0 + 0.5);
00487          y0 = VIGRA_CSTD::floor((ORDER % 2) ? y0 : y0 + 0.5);
00488          x1 = VIGRA_CSTD::floor((ORDER % 2) ? x1 : x1 + 0.5);
00489          y1 = VIGRA_CSTD::floor((ORDER % 2) ? y1 : y1 + 0.5);
00490          return x0 == x1 && y0 == y1;
00491     }
00492 
00493   protected:
00494 
00495     void init();
00496     void calculateIndices(double x, double y) const;
00497     void coefficients(double t, double * const & c) const;
00498     void derivCoefficients(double t, unsigned int d, double * const & c) const;
00499     value_type convolve() const;
00500 
00501     unsigned int w_, h_;
00502     int w1_, h1_;
00503     double x0_, x1_, y0_, y1_;
00504     InternalImage image_;
00505     Spline k_;
00506     mutable double x_, y_, u_, v_, kx_[ksize_], ky_[ksize_];
00507     mutable int ix_[ksize_], iy_[ksize_];
00508 };
00509 
00510 template <int ORDER, class VALUETYPE>
00511 void SplineImageView<ORDER, VALUETYPE>::init()
00512 {
00513     ArrayVector<double> const & b = k_.prefilterCoefficients();
00514 
00515     for(unsigned int i=0; i<b.size(); ++i)
00516     {
00517         recursiveFilterX(srcImageRange(image_), destImage(image_), b[i], BORDER_TREATMENT_REFLECT);
00518         recursiveFilterY(srcImageRange(image_), destImage(image_), b[i], BORDER_TREATMENT_REFLECT);
00519     }
00520 }
00521 
00522 namespace detail
00523 {
00524 
00525 template <int i>
00526 struct SplineImageViewUnrollLoop1
00527 {
00528     template <class Array>
00529     static void exec(int c0, Array c)
00530     {
00531         SplineImageViewUnrollLoop1<i-1>::exec(c0, c);
00532         c[i] = c0 + i;
00533     }
00534 };
00535 
00536 template <>
00537 struct SplineImageViewUnrollLoop1<0>
00538 {
00539     template <class Array>
00540     static void exec(int c0, Array c)
00541     {
00542         c[0] = c0;
00543     }
00544 };
00545 
00546 template <int i, class ValueType>
00547 struct SplineImageViewUnrollLoop2
00548 {
00549     template <class Array1, class RowIterator, class Array2>
00550     static ValueType
00551     exec(Array1 k, RowIterator r, Array2 x)
00552     {
00553         return ValueType(k[i] * r[x[i]]) + SplineImageViewUnrollLoop2<i-1, ValueType>::exec(k, r, x);
00554     }
00555 };
00556 
00557 template <class ValueType>
00558 struct SplineImageViewUnrollLoop2<0, ValueType>
00559 {
00560     template <class Array1, class RowIterator, class Array2>
00561     static ValueType
00562     exec(Array1 k, RowIterator r, Array2 x)
00563     {
00564         return ValueType(k[0] * r[x[0]]);
00565     }
00566 };
00567 
00568 } // namespace detail
00569 
00570 template <int ORDER, class VALUETYPE>
00571 void
00572 SplineImageView<ORDER, VALUETYPE>::calculateIndices(double x, double y) const
00573 {
00574     if(x == x_ && y == y_)
00575         return;   // still in cache
00576 
00577     if(x > x0_ && x < x1_ && y > y0_ && y < y1_)
00578     {
00579         detail::SplineImageViewUnrollLoop1<ORDER>::exec(
00580                                 (ORDER % 2) ? int(x - kcenter_) : int(x + 0.5 - kcenter_), ix_);
00581         detail::SplineImageViewUnrollLoop1<ORDER>::exec(
00582                                 (ORDER % 2) ? int(y - kcenter_) : int(y + 0.5 - kcenter_), iy_);
00583 
00584         u_ = x - ix_[kcenter_];
00585         v_ = y - iy_[kcenter_];
00586     }
00587     else
00588     {
00589         vigra_precondition(isValid(x,y),
00590                     "SplineImageView::calculateIndices(): coordinates out of range.");
00591 
00592         int xCenter = (ORDER % 2) ?
00593                       (int)VIGRA_CSTD::floor(x) :
00594                       (int)VIGRA_CSTD::floor(x + 0.5);
00595         int yCenter = (ORDER % 2) ?
00596                       (int)VIGRA_CSTD::floor(y) :
00597                       (int)VIGRA_CSTD::floor(y + 0.5);
00598 
00599         if(x >= x1_)
00600         {
00601             for(int i = 0; i < ksize_; ++i)
00602                 ix_[i] = w1_ - vigra::abs(w1_ - xCenter - (i - kcenter_));
00603         }
00604         else
00605         {
00606             for(int i = 0; i < ksize_; ++i)
00607                 ix_[i] = vigra::abs(xCenter - (kcenter_ - i));
00608         }
00609         if(y >= y1_)
00610         {
00611             for(int i = 0; i < ksize_; ++i)
00612                 iy_[i] = h1_ - vigra::abs(h1_ - yCenter - (i - kcenter_));
00613         }
00614         else
00615         {
00616             for(int i = 0; i < ksize_; ++i)
00617                 iy_[i] = vigra::abs(yCenter - (kcenter_ - i));
00618         }
00619         u_ = x - xCenter;
00620         v_ = y - yCenter;
00621     }
00622     x_ = x;
00623     y_ = y;
00624 }
00625 
00626 template <int ORDER, class VALUETYPE>
00627 void SplineImageView<ORDER, VALUETYPE>::coefficients(double t, double * const & c) const
00628 {
00629     t += kcenter_;
00630     for(int i = 0; i<ksize_; ++i)
00631         c[i] = k_(t-i);
00632 }
00633 
00634 template <int ORDER, class VALUETYPE>
00635 void SplineImageView<ORDER, VALUETYPE>::derivCoefficients(double t,
00636                                                unsigned int d, double * const & c) const
00637 {
00638     t += kcenter_;
00639     for(int i = 0; i<ksize_; ++i)
00640         c[i] = k_(t-i, d);
00641 }
00642 
00643 template <int ORDER, class VALUETYPE>
00644 VALUETYPE SplineImageView<ORDER, VALUETYPE>::convolve() const
00645 {
00646     typedef typename NumericTraits<VALUETYPE>::RealPromote RealPromote;
00647     RealPromote sum;
00648     sum = RealPromote(
00649       ky_[0]*detail::SplineImageViewUnrollLoop2<ORDER, RealPromote>::exec(kx_, image_.rowBegin(iy_[0]), ix_));
00650 
00651     for(int j=1; j<ksize_; ++j)
00652     {
00653         sum += RealPromote(
00654           ky_[j]*detail::SplineImageViewUnrollLoop2<ORDER, RealPromote>::exec(kx_, image_.rowBegin(iy_[j]), ix_));
00655     }
00656     return detail::RequiresExplicitCast<VALUETYPE>::cast(sum);
00657 }
00658 
00659 template <int ORDER, class VALUETYPE>
00660 template <class Array>
00661 void
00662 SplineImageView<ORDER, VALUETYPE>::coefficientArray(double x, double y, Array & res) const
00663 {
00664     typename Spline::WeightMatrix & weights = Spline::weights();
00665     double tmp[ksize_][ksize_];
00666 
00667     calculateIndices(x, y);
00668     for(int j=0; j<ksize_; ++j)
00669     {
00670         for(int i=0; i<ksize_; ++i)
00671         {
00672             tmp[i][j] = 0.0;
00673             for(int k=0; k<ksize_; ++k)
00674             {
00675                 tmp[i][j] += weights[i][k]*image_(ix_[k], iy_[j]);
00676             }
00677        }
00678     }
00679     res.resize(ksize_, ksize_);
00680     for(int j=0; j<ksize_; ++j)
00681     {
00682         for(int i=0; i<ksize_; ++i)
00683         {
00684             res(i,j) = 0.0;
00685             for(int k=0; k<ksize_; ++k)
00686             {
00687                 res(i,j) += weights[j][k]*tmp[i][k];
00688             }
00689         }
00690     }
00691 }
00692 
00693 template <int ORDER, class VALUETYPE>
00694 VALUETYPE SplineImageView<ORDER, VALUETYPE>::operator()(double x, double y) const
00695 {
00696     calculateIndices(x, y);
00697     coefficients(u_, kx_);
00698     coefficients(v_, ky_);
00699     return convolve();
00700 }
00701 
00702 template <int ORDER, class VALUETYPE>
00703 VALUETYPE SplineImageView<ORDER, VALUETYPE>::operator()(double x, double y,
00704                                                  unsigned int dx, unsigned int dy) const
00705 {
00706     calculateIndices(x, y);
00707     derivCoefficients(u_, dx, kx_);
00708     derivCoefficients(v_, dy, ky_);
00709     return convolve();
00710 }
00711 
00712 template <int ORDER, class VALUETYPE>
00713 VALUETYPE SplineImageView<ORDER, VALUETYPE>::g2(double x, double y) const
00714 {
00715     return sq(dx(x,y)) + sq(dy(x,y));
00716 }
00717 
00718 template <int ORDER, class VALUETYPE>
00719 VALUETYPE SplineImageView<ORDER, VALUETYPE>::g2x(double x, double y) const
00720 {
00721     return VALUETYPE(2.0)*(dx(x,y) * dxx(x,y) + dy(x,y) * dxy(x,y));
00722 }
00723 
00724 template <int ORDER, class VALUETYPE>
00725 VALUETYPE SplineImageView<ORDER, VALUETYPE>::g2y(double x, double y) const
00726 {
00727     return VALUETYPE(2.0)*(dx(x,y) * dxy(x,y) + dy(x,y) * dyy(x,y));
00728 }
00729 
00730 template <int ORDER, class VALUETYPE>
00731 VALUETYPE SplineImageView<ORDER, VALUETYPE>::g2xx(double x, double y) const
00732 {
00733     return VALUETYPE(2.0)*(sq(dxx(x,y)) + dx(x,y) * dx3(x,y) + sq(dxy(x,y)) + dy(x,y) * dxxy(x,y));
00734 }
00735 
00736 template <int ORDER, class VALUETYPE>
00737 VALUETYPE SplineImageView<ORDER, VALUETYPE>::g2yy(double x, double y) const
00738 {
00739     return VALUETYPE(2.0)*(sq(dxy(x,y)) + dx(x,y) * dxyy(x,y) + sq(dyy(x,y)) + dy(x,y) * dy3(x,y));
00740 }
00741 
00742 template <int ORDER, class VALUETYPE>
00743 VALUETYPE SplineImageView<ORDER, VALUETYPE>::g2xy(double x, double y) const
00744 {
00745     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)));
00746 }
00747 
00748 /********************************************************/
00749 /*                                                      */
00750 /*                    SplineImageView0                  */
00751 /*                                                      */
00752 /********************************************************/
00753 template <class VALUETYPE, class INTERNAL_INDEXER>
00754 class SplineImageView0Base
00755 {
00756     typedef typename INTERNAL_INDEXER::value_type InternalValue;
00757   public:
00758     typedef VALUETYPE value_type;
00759     typedef Size2D size_type;
00760     typedef TinyVector<double, 2> difference_type;
00761     enum StaticOrder { order = 0 };
00762 
00763   public:
00764 
00765     SplineImageView0Base(unsigned int w, unsigned int h)
00766     : w_(w), h_(h)
00767     {}
00768 
00769     SplineImageView0Base(int w, int h, INTERNAL_INDEXER i)
00770     : w_(w), h_(h), internalIndexer_(i)
00771     {}
00772 
00773     template <unsigned IntBits1, unsigned FractionalBits1,
00774               unsigned IntBits2, unsigned FractionalBits2>
00775     value_type unchecked(FixedPoint<IntBits1, FractionalBits1> x,
00776                          FixedPoint<IntBits2, FractionalBits2> y) const
00777     {
00778         return internalIndexer_(round(x), round(y));
00779     }
00780 
00781     template <unsigned IntBits1, unsigned FractionalBits1,
00782               unsigned IntBits2, unsigned FractionalBits2>
00783     value_type unchecked(FixedPoint<IntBits1, FractionalBits1> x,
00784                          FixedPoint<IntBits2, FractionalBits2> y,
00785                          unsigned int dx, unsigned int dy) const
00786     {
00787         if((dx != 0) || (dy != 0))
00788             return NumericTraits<VALUETYPE>::zero();
00789         return unchecked(x, y);
00790     }
00791 
00792     value_type unchecked(double x, double y) const
00793     {
00794         return internalIndexer_((int)(x + 0.5), (int)(y + 0.5));
00795     }
00796 
00797     value_type unchecked(double x, double y, unsigned int dx, unsigned int dy) const
00798     {
00799         if((dx != 0) || (dy != 0))
00800             return NumericTraits<VALUETYPE>::zero();
00801         return unchecked(x, y);
00802     }
00803 
00804     value_type operator()(double x, double y) const
00805     {
00806         int ix, iy;
00807         if(x < 0.0)
00808         {
00809             ix = (int)(-x + 0.5);
00810             vigra_precondition(ix <= (int)w_ - 1,
00811                     "SplineImageView::operator(): coordinates out of range.");
00812         }
00813         else
00814         {
00815             ix = (int)(x + 0.5);
00816             if(ix >= (int)w_)
00817             {
00818                 ix = 2*w_-2-ix;
00819                 vigra_precondition(ix >= 0,
00820                         "SplineImageView::operator(): coordinates out of range.");
00821             }
00822         }
00823         if(y < 0.0)
00824         {
00825             iy = (int)(-y + 0.5);
00826             vigra_precondition(iy <= (int)h_ - 1,
00827                     "SplineImageView::operator(): coordinates out of range.");
00828         }
00829         else
00830         {
00831             iy = (int)(y + 0.5);
00832             if(iy >= (int)h_)
00833             {
00834                 iy = 2*h_-2-iy;
00835                 vigra_precondition(iy >= 0,
00836                         "SplineImageView::operator(): coordinates out of range.");
00837             }
00838         }
00839         return internalIndexer_(ix, iy);
00840     }
00841 
00842     value_type operator()(double x, double y, unsigned int dx, unsigned int dy) const
00843     {
00844         if((dx != 0) || (dy != 0))
00845             return NumericTraits<VALUETYPE>::zero();
00846         return operator()(x, y);
00847     }
00848 
00849     value_type dx(double x, double y) const
00850         { return NumericTraits<VALUETYPE>::zero(); }
00851 
00852     value_type dy(double x, double y) const
00853         { return NumericTraits<VALUETYPE>::zero(); }
00854 
00855     value_type dxx(double x, double y) const
00856         { return NumericTraits<VALUETYPE>::zero(); }
00857 
00858     value_type dxy(double x, double y) const
00859         { return NumericTraits<VALUETYPE>::zero(); }
00860 
00861     value_type dyy(double x, double y) const
00862         { return NumericTraits<VALUETYPE>::zero(); }
00863 
00864     value_type dx3(double x, double y) const
00865         { return NumericTraits<VALUETYPE>::zero(); }
00866 
00867     value_type dy3(double x, double y) const
00868         { return NumericTraits<VALUETYPE>::zero(); }
00869 
00870     value_type dxxy(double x, double y) const
00871         { return NumericTraits<VALUETYPE>::zero(); }
00872 
00873     value_type dxyy(double x, double y) const
00874         { return NumericTraits<VALUETYPE>::zero(); }
00875 
00876     value_type operator()(difference_type const & d) const
00877         { return operator()(d[0], d[1]); }
00878 
00879     value_type operator()(difference_type const & d, unsigned int dx, unsigned int dy) const
00880         { return operator()(d[0], d[1], dx, dy); }
00881 
00882     value_type dx(difference_type const & d) const
00883         { return NumericTraits<VALUETYPE>::zero(); }
00884 
00885     value_type dy(difference_type const & d) const
00886         { return NumericTraits<VALUETYPE>::zero(); }
00887 
00888     value_type dxx(difference_type const & d) const
00889         { return NumericTraits<VALUETYPE>::zero(); }
00890 
00891     value_type dxy(difference_type const & d) const
00892         { return NumericTraits<VALUETYPE>::zero(); }
00893 
00894     value_type dyy(difference_type const & d) const
00895         { return NumericTraits<VALUETYPE>::zero(); }
00896 
00897     value_type dx3(difference_type const & d) const
00898         { return NumericTraits<VALUETYPE>::zero(); }
00899 
00900     value_type dy3(difference_type const & d) const
00901         { return NumericTraits<VALUETYPE>::zero(); }
00902 
00903     value_type dxxy(difference_type const & d) const
00904         { return NumericTraits<VALUETYPE>::zero(); }
00905 
00906     value_type dxyy(difference_type const & d) const
00907         { return NumericTraits<VALUETYPE>::zero(); }
00908 
00909     value_type g2(double x, double y) const
00910         { return NumericTraits<VALUETYPE>::zero(); }
00911 
00912     value_type g2x(double x, double y) const
00913         { return NumericTraits<VALUETYPE>::zero(); }
00914 
00915     value_type g2y(double x, double y) const
00916         { return NumericTraits<VALUETYPE>::zero(); }
00917 
00918     value_type g2xx(double x, double y) const
00919         { return NumericTraits<VALUETYPE>::zero(); }
00920 
00921     value_type g2xy(double x, double y) const
00922         { return NumericTraits<VALUETYPE>::zero(); }
00923 
00924     value_type g2yy(double x, double y) const
00925         { return NumericTraits<VALUETYPE>::zero(); }
00926 
00927     value_type g2(difference_type const & d) const
00928         { return NumericTraits<VALUETYPE>::zero(); }
00929 
00930     value_type g2x(difference_type const & d) const
00931         { return NumericTraits<VALUETYPE>::zero(); }
00932 
00933     value_type g2y(difference_type const & d) const
00934         { return NumericTraits<VALUETYPE>::zero(); }
00935 
00936     value_type g2xx(difference_type const & d) const
00937         { return NumericTraits<VALUETYPE>::zero(); }
00938 
00939     value_type g2xy(difference_type const & d) const
00940         { return NumericTraits<VALUETYPE>::zero(); }
00941 
00942     value_type g2yy(difference_type const & d) const
00943         { return NumericTraits<VALUETYPE>::zero(); }
00944 
00945     unsigned int width() const
00946         { return w_; }
00947 
00948     unsigned int height() const
00949         { return h_; }
00950 
00951     size_type size() const
00952         { return size_type(w_, h_); }
00953 
00954     TinyVector<unsigned int, 2> shape() const
00955         { return TinyVector<unsigned int, 2>(w_, h_); }
00956 
00957     template <class Array>
00958     void coefficientArray(double x, double y, Array & res) const
00959     {
00960         res.resize(1, 1);
00961         res(0, 0) = operator()(x,y);
00962     }
00963 
00964     bool isInsideX(double x) const
00965     {
00966         return x >= 0.0 && x <= width() - 1.0;
00967     }
00968 
00969     bool isInsideY(double y) const
00970     {
00971         return y >= 0.0 && y <= height() - 1.0;
00972     }
00973 
00974     bool isInside(double x, double y) const
00975     {
00976         return isInsideX(x) && isInsideY(y);
00977     }
00978 
00979     bool isValid(double x, double y) const
00980     {
00981         return x < 2.0*w_-2.0 && x > 1.0-w_ && y < 2.0*h_-2.0 && y > 1.0-h_;
00982     }
00983 
00984     bool sameFacet(double x0, double y0, double x1, double y1) const
00985     {
00986          x0 = VIGRA_CSTD::floor(x0 + 0.5);
00987          y0 = VIGRA_CSTD::floor(y0 + 0.5);
00988          x1 = VIGRA_CSTD::floor(x1 + 0.5);
00989          y1 = VIGRA_CSTD::floor(y1 + 0.5);
00990          return x0 == x1 && y0 == y1;
00991     }
00992 
00993   protected:
00994     unsigned int w_, h_;
00995     INTERNAL_INDEXER internalIndexer_;
00996 };
00997 
00998 /** \brief Create an image view for nearest-neighbor interpolation.
00999 
01000     This class behaves like \ref vigra::SplineImageView&lt;0, ...&gt;, but one can pass
01001     an additional template argument that determined the internal representation of the image.
01002     If this is equal to the argument type passed in the constructor, the image is not copied.
01003     By default, this works for \ref vigra::BasicImage, \ref vigra::BasicImageView,
01004     \ref vigra::MultiArray&lt;2, ...&gt;, and \ref vigra::MultiArrayView&lt;2, ...&gt;.
01005 
01006 */
01007 template <class VALUETYPE, class INTERNAL_TRAVERSER = typename BasicImage<VALUETYPE>::const_traverser>
01008 class SplineImageView0
01009 : public SplineImageView0Base<VALUETYPE, INTERNAL_TRAVERSER>
01010 {
01011     typedef SplineImageView0Base<VALUETYPE, INTERNAL_TRAVERSER> Base;
01012   public:
01013     typedef typename Base::value_type value_type;
01014     typedef typename Base::size_type size_type;
01015     typedef typename Base::difference_type difference_type;
01016     enum StaticOrder { order = Base::order };
01017     typedef BasicImage<VALUETYPE> InternalImage;
01018 
01019   protected:
01020     typedef typename IteratorTraits<INTERNAL_TRAVERSER>::mutable_iterator InternalTraverser;
01021     typedef typename IteratorTraits<InternalTraverser>::DefaultAccessor InternalAccessor;
01022     typedef typename IteratorTraits<INTERNAL_TRAVERSER>::const_iterator InternalConstTraverser;
01023     typedef typename IteratorTraits<InternalConstTraverser>::DefaultAccessor InternalConstAccessor;
01024 
01025   public:
01026 
01027         /* when traverser and accessor types passed to the constructor are the same as the corresponding
01028            internal types, we need not copy the image (speed up)
01029         */
01030     SplineImageView0(InternalTraverser is, InternalTraverser iend, InternalAccessor sa)
01031     : Base(iend.x - is.x, iend.y - is.y, is)
01032     {}
01033 
01034     SplineImageView0(triple<InternalTraverser, InternalTraverser, InternalAccessor> s)
01035     : Base(s.second.x - s.first.x, s.second.y - s.first.y, s.first)
01036     {}
01037 
01038     SplineImageView0(InternalConstTraverser is, InternalConstTraverser iend, InternalConstAccessor sa)
01039     : Base(iend.x - is.x, iend.y - is.y, is)
01040     {}
01041 
01042     SplineImageView0(triple<InternalConstTraverser, InternalConstTraverser, InternalConstAccessor> s)
01043     : Base(s.second.x - s.first.x, s.second.y - s.first.y, s.first)
01044     {}
01045 
01046     template<class T, class SU>
01047     SplineImageView0(MultiArrayView<2, T, SU> const & i)
01048     : Base(i.shape(0), i.shape(1)),
01049       image_(i.shape(0), i.shape(1))
01050     {
01051         for(unsigned int y=0; y<this->height(); ++y)
01052             for(unsigned int x=0; x<this->width(); ++x)
01053                 image_(x,y) = detail::RequiresExplicitCast<VALUETYPE>::cast(i(x,y));
01054         this->internalIndexer_ = image_.upperLeft();
01055     }
01056 
01057     template <class SrcIterator, class SrcAccessor>
01058     SplineImageView0(SrcIterator is, SrcIterator iend, SrcAccessor sa)
01059     : Base(iend.x - is.x, iend.y - is.y),
01060       image_(iend - is)
01061     {
01062         copyImage(srcIterRange(is, iend, sa), destImage(image_));
01063         this->internalIndexer_ = image_.upperLeft();
01064     }
01065 
01066     template <class SrcIterator, class SrcAccessor>
01067     SplineImageView0(triple<SrcIterator, SrcIterator, SrcAccessor> s)
01068     : Base(s.second.x - s.first.x, s.second.y - s.first.y),
01069       image_(s.second - s.first)
01070     {
01071         copyImage(s, destImage(image_));
01072         this->internalIndexer_ = image_.upperLeft();
01073     }
01074 
01075     InternalImage const & image() const
01076         { return image_; }
01077 
01078   protected:
01079     InternalImage image_;
01080 };
01081 
01082 template <class VALUETYPE, class StridedOrUnstrided>
01083 class SplineImageView0<VALUETYPE, MultiArrayView<2, VALUETYPE, StridedOrUnstrided> >
01084 : public SplineImageView0Base<VALUETYPE, MultiArrayView<2, VALUETYPE, StridedOrUnstrided> >
01085 {
01086     typedef SplineImageView0Base<VALUETYPE, MultiArrayView<2, VALUETYPE, StridedOrUnstrided> > Base;
01087   public:
01088     typedef typename Base::value_type value_type;
01089     typedef typename Base::size_type size_type;
01090     typedef typename Base::difference_type difference_type;
01091     enum StaticOrder { order = Base::order };
01092     typedef BasicImage<VALUETYPE> InternalImage;
01093 
01094   protected:
01095     typedef MultiArrayView<2, VALUETYPE, StridedOrUnstrided> InternalIndexer;
01096 
01097   public:
01098 
01099         /* when traverser and accessor types passed to the constructor are the same as the corresponding
01100            internal types, we need not copy the image (speed up)
01101         */
01102     SplineImageView0(InternalIndexer const & i)
01103     : Base(i.shape(0), i.shape(1), i)
01104     {}
01105 
01106     template<class T, class SU>
01107     SplineImageView0(MultiArrayView<2, T, SU> const & i)
01108     : Base(i.shape(0), i.shape(1)),
01109       image_(i.shape(0), i.shape(1))
01110     {
01111         for(unsigned int y=0; y<this->height(); ++y)
01112             for(unsigned int x=0; x<this->width(); ++x)
01113                 image_(x,y) = detail::RequiresExplicitCast<VALUETYPE>::cast(i(x,y));
01114         this->internalIndexer_ = InternalIndexer(typename InternalIndexer::difference_type(this->width(), this->height()),
01115                                                  image_.data());
01116     }
01117 
01118     template <class SrcIterator, class SrcAccessor>
01119     SplineImageView0(SrcIterator is, SrcIterator iend, SrcAccessor sa)
01120     : Base(iend.x - is.x, iend.y - is.y),
01121       image_(iend-is)
01122     {
01123         copyImage(srcIterRange(is, iend, sa), destImage(image_));
01124         this->internalIndexer_ = InternalIndexer(typename InternalIndexer::difference_type(this->width(), this->height()),
01125                                                  image_.data());
01126     }
01127 
01128     template <class SrcIterator, class SrcAccessor>
01129     SplineImageView0(triple<SrcIterator, SrcIterator, SrcAccessor> s)
01130     : Base(s.second.x - s.first.x, s.second.y - s.first.y),
01131       image_(s.second - s.first)
01132     {
01133         copyImage(s, destImage(image_));
01134         this->internalIndexer_ = InternalIndexer(typename InternalIndexer::difference_type(this->width(), this->height()),
01135                                                  image_.data());
01136     }
01137 
01138     InternalImage const & image() const
01139         { return image_; }
01140 
01141   protected:
01142     InternalImage image_;
01143 };
01144 
01145 template <class VALUETYPE>
01146 class SplineImageView<0, VALUETYPE>
01147 : public SplineImageView0<VALUETYPE>
01148 {
01149     typedef SplineImageView0<VALUETYPE> Base;
01150   public:
01151     typedef typename Base::value_type value_type;
01152     typedef typename Base::size_type size_type;
01153     typedef typename Base::difference_type difference_type;
01154     enum StaticOrder { order = Base::order };
01155     typedef typename Base::InternalImage InternalImage;
01156 
01157   protected:
01158     typedef typename Base::InternalTraverser InternalTraverser;
01159     typedef typename Base::InternalAccessor InternalAccessor;
01160     typedef typename Base::InternalConstTraverser InternalConstTraverser;
01161     typedef typename Base::InternalConstAccessor InternalConstAccessor;
01162 
01163 public:
01164 
01165         /* when traverser and accessor types passed to the constructor are the same as the corresponding
01166            internal types, we need not copy the image (speed up)
01167         */
01168     SplineImageView(InternalTraverser is, InternalTraverser iend, InternalAccessor sa, bool /* unused */ = false)
01169     : Base(is, iend, sa)
01170     {}
01171 
01172     SplineImageView(triple<InternalTraverser, InternalTraverser, InternalAccessor> s, bool /* unused */ = false)
01173     : Base(s)
01174     {}
01175 
01176     SplineImageView(InternalConstTraverser is, InternalConstTraverser iend, InternalConstAccessor sa, bool /* unused */ = false)
01177     : Base(is, iend, sa)
01178     {}
01179 
01180     SplineImageView(triple<InternalConstTraverser, InternalConstTraverser, InternalConstAccessor> s, bool /* unused */ = false)
01181     : Base(s)
01182     {}
01183 
01184     template <class SrcIterator, class SrcAccessor>
01185     SplineImageView(SrcIterator is, SrcIterator iend, SrcAccessor sa, bool /* unused */ = false)
01186     : Base(is, iend, sa)
01187     {
01188         copyImage(srcIterRange(is, iend, sa), destImage(this->image_));
01189     }
01190 
01191     template <class SrcIterator, class SrcAccessor>
01192     SplineImageView(triple<SrcIterator, SrcIterator, SrcAccessor> s, bool /* unused */ = false)
01193     : Base(s)
01194     {
01195         copyImage(s, destImage(this->image_));
01196     }
01197 };
01198 
01199 /********************************************************/
01200 /*                                                      */
01201 /*                    SplineImageView1                  */
01202 /*                                                      */
01203 /********************************************************/
01204 template <class VALUETYPE, class INTERNAL_INDEXER>
01205 class SplineImageView1Base
01206 {
01207     typedef typename INTERNAL_INDEXER::value_type InternalValue;
01208   public:
01209     typedef VALUETYPE value_type;
01210     typedef Size2D size_type;
01211     typedef TinyVector<double, 2> difference_type;
01212     enum StaticOrder { order = 1 };
01213 
01214   public:
01215 
01216     SplineImageView1Base(unsigned int w, unsigned int h)
01217     : w_(w), h_(h)
01218     {}
01219 
01220     SplineImageView1Base(int w, int h, INTERNAL_INDEXER i)
01221     : w_(w), h_(h), internalIndexer_(i)
01222     {}
01223 
01224     template <unsigned IntBits1, unsigned FractionalBits1,
01225               unsigned IntBits2, unsigned FractionalBits2>
01226     value_type unchecked(FixedPoint<IntBits1, FractionalBits1> x,
01227                          FixedPoint<IntBits2, FractionalBits2> y) const
01228     {
01229         int ix = floor(x);
01230         FixedPoint<0, FractionalBits1> tx = frac(x - FixedPoint<IntBits1, FractionalBits1>(ix));
01231         FixedPoint<0, FractionalBits1> dtx = dual_frac(tx);
01232         if(ix == (int)w_ - 1)
01233         {
01234             --ix;
01235             tx.value = FixedPoint<0, FractionalBits1>::ONE;
01236             dtx.value = 0;
01237         }
01238         int iy = floor(y);
01239         FixedPoint<0, FractionalBits2> ty = frac(y - FixedPoint<IntBits2, FractionalBits2>(iy));
01240         FixedPoint<0, FractionalBits2> dty = dual_frac(ty);
01241         if(iy == (int)h_ - 1)
01242         {
01243             --iy;
01244             ty.value = FixedPoint<0, FractionalBits2>::ONE;
01245             dty.value = 0;
01246         }
01247         return fixed_point_cast<value_type>(
01248                     dty*(dtx*fixedPoint(internalIndexer_(ix,iy)) +
01249                                    tx*fixedPoint(internalIndexer_(ix+1,iy))) +
01250                     ty *(dtx*fixedPoint(internalIndexer_(ix,iy+1)) +
01251                                    tx*fixedPoint(internalIndexer_(ix+1,iy+1))));
01252     }
01253 
01254     template <unsigned IntBits1, unsigned FractionalBits1,
01255               unsigned IntBits2, unsigned FractionalBits2>
01256     value_type unchecked(FixedPoint<IntBits1, FractionalBits1> x,
01257                          FixedPoint<IntBits2, FractionalBits2> y,
01258                          unsigned int dx, unsigned int dy) const
01259     {
01260         int ix = floor(x);
01261         FixedPoint<0, FractionalBits1> tx = frac(x - FixedPoint<IntBits1, FractionalBits1>(ix));
01262         FixedPoint<0, FractionalBits1> dtx = dual_frac(tx);
01263         if(ix == (int)w_ - 1)
01264         {
01265             --ix;
01266             tx.value = FixedPoint<0, FractionalBits1>::ONE;
01267             dtx.value = 0;
01268         }
01269         int iy = floor(y);
01270         FixedPoint<0, FractionalBits2> ty = frac(y - FixedPoint<IntBits2, FractionalBits2>(iy));
01271         FixedPoint<0, FractionalBits2> dty = dual_frac(ty);
01272         if(iy == (int)h_ - 1)
01273         {
01274             --iy;
01275             ty.value = FixedPoint<0, FractionalBits2>::ONE;
01276             dty.value = 0;
01277         }
01278         switch(dx)
01279         {
01280           case 0:
01281               switch(dy)
01282               {
01283                 case 0:
01284                     return fixed_point_cast<value_type>(
01285                                 dty*(dtx*fixedPoint(internalIndexer_(ix,iy)) +
01286                                                tx*fixedPoint(internalIndexer_(ix+1,iy))) +
01287                                 ty *(dtx*fixedPoint(internalIndexer_(ix,iy+1)) +
01288                                                tx*fixedPoint(internalIndexer_(ix+1,iy+1))));
01289                 case 1:
01290                     return fixed_point_cast<value_type>(
01291                            (dtx*fixedPoint(internalIndexer_(ix,iy+1)) + tx*fixedPoint(internalIndexer_(ix+1,iy+1))) -
01292                            (dtx*fixedPoint(internalIndexer_(ix,iy)) + tx*fixedPoint(internalIndexer_(ix+1,iy))));
01293                 default:
01294                     return NumericTraits<VALUETYPE>::zero();
01295               }
01296           case 1:
01297               switch(dy)
01298               {
01299                 case 0:
01300                     return fixed_point_cast<value_type>(
01301                                 dty*(fixedPoint(internalIndexer_(ix+1,iy)) - fixedPoint(internalIndexer_(ix,iy))) +
01302                                 ty *(fixedPoint(internalIndexer_(ix+1,iy+1)) - fixedPoint(internalIndexer_(ix,iy+1))));
01303                 case 1:
01304                     return detail::RequiresExplicitCast<value_type>::cast(
01305                                 (internalIndexer_(ix+1,iy+1) - internalIndexer_(ix,iy+1)) -
01306                                 (internalIndexer_(ix+1,iy) - internalIndexer_(ix,iy)));
01307                 default:
01308                     return NumericTraits<VALUETYPE>::zero();
01309               }
01310           default:
01311               return NumericTraits<VALUETYPE>::zero();
01312         }
01313     }
01314 
01315     value_type unchecked(double x, double y) const
01316     {
01317         int ix = (int)std::floor(x);
01318         if(ix == (int)w_ - 1)
01319             --ix;
01320         double tx = x - ix;
01321         int iy = (int)std::floor(y);
01322         if(iy == (int)h_ - 1)
01323             --iy;
01324         double ty = y - iy;
01325         return NumericTraits<value_type>::fromRealPromote(
01326                    (1.0-ty)*((1.0-tx)*internalIndexer_(ix,iy) + tx*internalIndexer_(ix+1,iy)) +
01327                     ty *((1.0-tx)*internalIndexer_(ix,iy+1) + tx*internalIndexer_(ix+1,iy+1)));
01328     }
01329 
01330     value_type unchecked(double x, double y, unsigned int dx, unsigned int dy) const
01331     {
01332         int ix = (int)std::floor(x);
01333         if(ix == (int)w_ - 1)
01334             --ix;
01335         double tx = x - ix;
01336         int iy = (int)std::floor(y);
01337         if(iy == (int)h_ - 1)
01338             --iy;
01339         double ty = y - iy;
01340         switch(dx)
01341         {
01342           case 0:
01343               switch(dy)
01344               {
01345                 case 0:
01346                     return detail::RequiresExplicitCast<value_type>::cast(
01347                                (1.0-ty)*((1.0-tx)*internalIndexer_(ix,iy) + tx*internalIndexer_(ix+1,iy)) +
01348                                 ty *((1.0-tx)*internalIndexer_(ix,iy+1) + tx*internalIndexer_(ix+1,iy+1)));
01349                 case 1:
01350                     return detail::RequiresExplicitCast<value_type>::cast(
01351                                ((1.0-tx)*internalIndexer_(ix,iy+1) + tx*internalIndexer_(ix+1,iy+1)) -
01352                                ((1.0-tx)*internalIndexer_(ix,iy) + tx*internalIndexer_(ix+1,iy)));
01353                 default:
01354                     return NumericTraits<VALUETYPE>::zero();
01355               }
01356           case 1:
01357               switch(dy)
01358               {
01359                 case 0:
01360                     return detail::RequiresExplicitCast<value_type>::cast(
01361                                (1.0-ty)*(internalIndexer_(ix+1,iy) - internalIndexer_(ix,iy)) +
01362                                 ty *(internalIndexer_(ix+1,iy+1) - internalIndexer_(ix,iy+1)));
01363                 case 1:
01364                     return detail::RequiresExplicitCast<value_type>::cast(
01365                               (internalIndexer_(ix+1,iy+1) - internalIndexer_(ix,iy+1)) -
01366                               (internalIndexer_(ix+1,iy) - internalIndexer_(ix,iy)));
01367                 default:
01368                     return NumericTraits<VALUETYPE>::zero();
01369               }
01370           default:
01371               return NumericTraits<VALUETYPE>::zero();
01372         }
01373     }
01374 
01375     value_type operator()(double x, double y) const
01376     {
01377         return operator()(x, y, 0, 0);
01378     }
01379 
01380     value_type operator()(double x, double y, unsigned int dx, unsigned int dy) const
01381     {
01382         value_type mul = NumericTraits<value_type>::one();
01383         if(x < 0.0)
01384         {
01385             x = -x;
01386             vigra_precondition(x <= w_ - 1.0,
01387                     "SplineImageView::operator(): coordinates out of range.");
01388             if(dx % 2)
01389                 mul = -mul;
01390         }
01391         else if(x > w_ - 1.0)
01392         {
01393             x = 2.0*w_-2.0-x;
01394             vigra_precondition(x >= 0.0,
01395                     "SplineImageView::operator(): coordinates out of range.");
01396             if(dx % 2)
01397                 mul = -mul;
01398         }
01399         if(y < 0.0)
01400         {
01401             y = -y;
01402             vigra_precondition(y <= h_ - 1.0,
01403                     "SplineImageView::operator(): coordinates out of range.");
01404             if(dy % 2)
01405                 mul = -mul;
01406         }
01407         else if(y > h_ - 1.0)
01408         {
01409             y = 2.0*h_-2.0-y;
01410             vigra_precondition(y >= 0.0,
01411                     "SplineImageView::operator(): coordinates out of range.");
01412             if(dy % 2)
01413                 mul = -mul;
01414         }
01415         return mul*unchecked(x, y, dx, dy);
01416     }
01417 
01418     value_type dx(double x, double y) const
01419         { return operator()(x, y, 1, 0); }
01420 
01421     value_type dy(double x, double y) const
01422         { return operator()(x, y, 0, 1); }
01423 
01424     value_type dxx(double x, double y) const
01425         { return NumericTraits<VALUETYPE>::zero(); }
01426 
01427     value_type dxy(double x, double y) const
01428         { return operator()(x, y, 1, 1); }
01429 
01430     value_type dyy(double x, double y) const
01431         { return NumericTraits<VALUETYPE>::zero(); }
01432 
01433     value_type dx3(double x, double y) const
01434         { return NumericTraits<VALUETYPE>::zero(); }
01435 
01436     value_type dy3(double x, double y) const
01437         { return NumericTraits<VALUETYPE>::zero(); }
01438 
01439     value_type dxxy(double x, double y) const
01440         { return NumericTraits<VALUETYPE>::zero(); }
01441 
01442     value_type dxyy(double x, double y) const
01443         { return NumericTraits<VALUETYPE>::zero(); }
01444 
01445     value_type operator()(difference_type const & d) const
01446         { return operator()(d[0], d[1]); }
01447 
01448     value_type operator()(difference_type const & d, unsigned int dx, unsigned int dy) const
01449         { return operator()(d[0], d[1], dx, dy); }
01450 
01451     value_type dx(difference_type const & d) const
01452         { return operator()(d[0], d[1], 1, 0); }
01453 
01454     value_type dy(difference_type const & d) const
01455         { return operator()(d[0], d[1], 0, 1); }
01456 
01457     value_type dxx(difference_type const & d) const
01458         { return NumericTraits<VALUETYPE>::zero(); }
01459 
01460     value_type dxy(difference_type const & d) const
01461         { return operator()(d[0], d[1], 1, 1); }
01462 
01463     value_type dyy(difference_type const & d) const
01464         { return NumericTraits<VALUETYPE>::zero(); }
01465 
01466     value_type dx3(difference_type const & d) const
01467         { return NumericTraits<VALUETYPE>::zero(); }
01468 
01469     value_type dy3(difference_type const & d) const
01470         { return NumericTraits<VALUETYPE>::zero(); }
01471 
01472     value_type dxxy(difference_type const & d) const
01473         { return NumericTraits<VALUETYPE>::zero(); }
01474 
01475     value_type dxyy(difference_type const & d) const
01476         { return NumericTraits<VALUETYPE>::zero(); }
01477 
01478     value_type g2(double x, double y) const
01479         { return sq(dx(x,y)) + sq(dy(x,y)); }
01480 
01481     value_type g2x(double x, double y) const
01482         { return NumericTraits<VALUETYPE>::zero(); }
01483 
01484     value_type g2y(double x, double y) const
01485         { return NumericTraits<VALUETYPE>::zero(); }
01486 
01487     value_type g2xx(double x, double y) const
01488         { return NumericTraits<VALUETYPE>::zero(); }
01489 
01490     value_type g2xy(double x, double y) const
01491         { return NumericTraits<VALUETYPE>::zero(); }
01492 
01493     value_type g2yy(double x, double y) const
01494         { return NumericTraits<VALUETYPE>::zero(); }
01495 
01496     value_type g2(difference_type const & d) const
01497         { return g2(d[0], d[1]); }
01498 
01499     value_type g2x(difference_type const & d) const
01500         { return NumericTraits<VALUETYPE>::zero(); }
01501 
01502     value_type g2y(difference_type const & d) const
01503         { return NumericTraits<VALUETYPE>::zero(); }
01504 
01505     value_type g2xx(difference_type const & d) const
01506         { return NumericTraits<VALUETYPE>::zero(); }
01507 
01508     value_type g2xy(difference_type const & d) const
01509         { return NumericTraits<VALUETYPE>::zero(); }
01510 
01511     value_type g2yy(difference_type const & d) const
01512         { return NumericTraits<VALUETYPE>::zero(); }
01513 
01514     unsigned int width() const
01515         { return w_; }
01516 
01517     unsigned int height() const
01518         { return h_; }
01519 
01520     size_type size() const
01521         { return size_type(w_, h_); }
01522 
01523     TinyVector<unsigned int, 2> shape() const
01524         { return TinyVector<unsigned int, 2>(w_, h_); }
01525 
01526     template <class Array>
01527     void coefficientArray(double x, double y, Array & res) const;
01528 
01529     void calculateIndices(double x, double y, int & ix, int & iy, int & ix1, int & iy1) const;
01530 
01531     bool isInsideX(double x) const
01532     {
01533         return x >= 0.0 && x <= width() - 1.0;
01534     }
01535 
01536     bool isInsideY(double y) const
01537     {
01538         return y >= 0.0 && y <= height() - 1.0;
01539     }
01540 
01541     bool isInside(double x, double y) const
01542     {
01543         return isInsideX(x) && isInsideY(y);
01544     }
01545 
01546     bool isValid(double x, double y) const
01547     {
01548         return x < 2.0*w_-2.0 && x > 1.0-w_ && y < 2.0*h_-2.0 && y > 1.0-h_;
01549     }
01550 
01551     bool sameFacet(double x0, double y0, double x1, double y1) const
01552     {
01553          x0 = VIGRA_CSTD::floor(x0);
01554          y0 = VIGRA_CSTD::floor(y0);
01555          x1 = VIGRA_CSTD::floor(x1);
01556          y1 = VIGRA_CSTD::floor(y1);
01557          return x0 == x1 && y0 == y1;
01558     }
01559 
01560   protected:
01561     unsigned int w_, h_;
01562     INTERNAL_INDEXER internalIndexer_;
01563 };
01564 
01565 template <class VALUETYPE, class INTERNAL_INDEXER>
01566 template <class Array>
01567 void SplineImageView1Base<VALUETYPE, INTERNAL_INDEXER>::coefficientArray(double x, double y, Array & res) const
01568 {
01569     int ix, iy, ix1, iy1;
01570     calculateIndices(x, y, ix, iy, ix1, iy1);
01571     res.resize(2, 2);
01572     res(0,0) = internalIndexer_(ix,iy);
01573     res(1,0) = internalIndexer_(ix1,iy) - internalIndexer_(ix,iy);
01574     res(0,1) = internalIndexer_(ix,iy1) - internalIndexer_(ix,iy);
01575     res(1,1) = internalIndexer_(ix,iy) - internalIndexer_(ix1,iy) -
01576                internalIndexer_(ix,iy1) + internalIndexer_(ix1,iy1);
01577 }
01578 
01579 template <class VALUETYPE, class INTERNAL_INDEXER>
01580 void SplineImageView1Base<VALUETYPE, INTERNAL_INDEXER>::calculateIndices(double x, double y, int & ix, int & iy, int & ix1, int & iy1) const
01581 {
01582     if(x < 0.0)
01583     {
01584         x = -x;
01585         vigra_precondition(x <= w_ - 1.0,
01586                 "SplineImageView::calculateIndices(): coordinates out of range.");
01587         ix = (int)VIGRA_CSTD::ceil(x);
01588         ix1 = ix - 1;
01589     }
01590     else if(x >= w_ - 1.0)
01591     {
01592         x = 2.0*w_-2.0-x;
01593         vigra_precondition(x > 0.0,
01594                 "SplineImageView::calculateIndices(): coordinates out of range.");
01595         ix = (int)VIGRA_CSTD::ceil(x);
01596         ix1 = ix - 1;
01597     }
01598     else
01599     {
01600         ix = (int)VIGRA_CSTD::floor(x);
01601         ix1 = ix + 1;
01602     }
01603     if(y < 0.0)
01604     {
01605         y = -y;
01606         vigra_precondition(y <= h_ - 1.0,
01607                 "SplineImageView::calculateIndices(): coordinates out of range.");
01608         iy = (int)VIGRA_CSTD::ceil(y);
01609         iy1 = iy - 1;
01610     }
01611     else if(y >= h_ - 1.0)
01612     {
01613         y = 2.0*h_-2.0-y;
01614         vigra_precondition(y > 0.0,
01615                 "SplineImageView::calculateIndices(): coordinates out of range.");
01616         iy = (int)VIGRA_CSTD::ceil(y);
01617         iy1 = iy - 1;
01618     }
01619     else
01620     {
01621         iy = (int)VIGRA_CSTD::floor(y);
01622         iy1 = iy + 1;
01623     }
01624 }
01625 
01626 /** \brief Create an image view for bi-linear interpolation.
01627 
01628     This class behaves like \ref vigra::SplineImageView&lt;1, ...&gt;, but one can pass
01629     an additional template argument that determined the internal representation of the image.
01630     If this is equal to the argument type passed in the constructor, the image is not copied.
01631     By default, this works for \ref vigra::BasicImage, \ref vigra::BasicImageView,
01632     \ref vigra::MultiArray&lt;2, ...&gt;, and \ref vigra::MultiArrayView&lt;2, ...&gt;.
01633 
01634     In addition to the function provided by  \ref vigra::SplineImageView, there are functions
01635     <tt>unchecked(x,y)</tt> and <tt>unchecked(x,y, xorder, yorder)</tt> which improve speed by
01636     not applying bounds checking and reflective border treatment (<tt>isInside(x, y)</tt> must
01637     be <tt>true</tt>), but otherwise behave identically to their checked counterparts.
01638     In addition, <tt>x</tt> and <tt>y</tt> can have type \ref vigra::FixedPoint instead of
01639     <tt>double</tt>.
01640 */
01641 template <class VALUETYPE, class INTERNAL_TRAVERSER = typename BasicImage<VALUETYPE>::const_traverser>
01642 class SplineImageView1
01643 : public SplineImageView1Base<VALUETYPE, INTERNAL_TRAVERSER>
01644 {
01645     typedef SplineImageView1Base<VALUETYPE, INTERNAL_TRAVERSER> Base;
01646   public:
01647     typedef typename Base::value_type value_type;
01648     typedef typename Base::size_type size_type;
01649     typedef typename Base::difference_type difference_type;
01650     enum StaticOrder { order = Base::order };
01651     typedef BasicImage<VALUETYPE> InternalImage;
01652 
01653   protected:
01654     typedef typename IteratorTraits<INTERNAL_TRAVERSER>::mutable_iterator InternalTraverser;
01655     typedef typename IteratorTraits<InternalTraverser>::DefaultAccessor InternalAccessor;
01656     typedef typename IteratorTraits<INTERNAL_TRAVERSER>::const_iterator InternalConstTraverser;
01657     typedef typename IteratorTraits<InternalConstTraverser>::DefaultAccessor InternalConstAccessor;
01658 
01659   public:
01660 
01661         /* when traverser and accessor types passed to the constructor are the same as the corresponding
01662            internal types, we need not copy the image (speed up)
01663         */
01664     SplineImageView1(InternalTraverser is, InternalTraverser iend, InternalAccessor sa)
01665     : Base(iend.x - is.x, iend.y - is.y, is)
01666     {}
01667 
01668     SplineImageView1(triple<InternalTraverser, InternalTraverser, InternalAccessor> s)
01669     : Base(s.second.x - s.first.x, s.second.y - s.first.y, s.first)
01670     {}
01671 
01672     SplineImageView1(InternalConstTraverser is, InternalConstTraverser iend, InternalConstAccessor sa)
01673     : Base(iend.x - is.x, iend.y - is.y, is)
01674     {}
01675 
01676     SplineImageView1(triple<InternalConstTraverser, InternalConstTraverser, InternalConstAccessor> s)
01677     : Base(s.second.x - s.first.x, s.second.y - s.first.y, s.first)
01678     {}
01679 
01680     template<class T, class SU>
01681     SplineImageView1(MultiArrayView<2, T, SU> const & i)
01682     : Base(i.shape(0), i.shape(1)),
01683       image_(i.shape(0), i.shape(1))
01684     {
01685         for(unsigned int y=0; y<this->height(); ++y)
01686             for(unsigned int x=0; x<this->width(); ++x)
01687                 image_(x,y) = detail::RequiresExplicitCast<VALUETYPE>::cast(i(x,y));
01688         this->internalIndexer_ = image_.upperLeft();
01689     }
01690 
01691     template <class SrcIterator, class SrcAccessor>
01692     SplineImageView1(SrcIterator is, SrcIterator iend, SrcAccessor sa)
01693     : Base(iend.x - is.x, iend.y - is.y),
01694       image_(iend - is)
01695     {
01696         copyImage(srcIterRange(is, iend, sa), destImage(image_));
01697         this->internalIndexer_ = image_.upperLeft();
01698     }
01699 
01700     template <class SrcIterator, class SrcAccessor>
01701     SplineImageView1(triple<SrcIterator, SrcIterator, SrcAccessor> s)
01702     : Base(s.second.x - s.first.x, s.second.y - s.first.y),
01703       image_(s.second - s.first)
01704     {
01705         copyImage(s, destImage(image_));
01706         this->internalIndexer_ = image_.upperLeft();
01707     }
01708 
01709     InternalImage const & image() const
01710         { return image_; }
01711 
01712   protected:
01713     InternalImage image_;
01714 };
01715 
01716 template <class VALUETYPE, class StridedOrUnstrided>
01717 class SplineImageView1<VALUETYPE, MultiArrayView<2, VALUETYPE, StridedOrUnstrided> >
01718 : public SplineImageView1Base<VALUETYPE, MultiArrayView<2, VALUETYPE, StridedOrUnstrided> >
01719 {
01720     typedef SplineImageView1Base<VALUETYPE, MultiArrayView<2, VALUETYPE, StridedOrUnstrided> > Base;
01721   public:
01722     typedef typename Base::value_type value_type;
01723     typedef typename Base::size_type size_type;
01724     typedef typename Base::difference_type difference_type;
01725     enum StaticOrder { order = Base::order };
01726     typedef BasicImage<VALUETYPE> InternalImage;
01727 
01728   protected:
01729     typedef MultiArrayView<2, VALUETYPE, StridedOrUnstrided> InternalIndexer;
01730 
01731   public:
01732 
01733         /* when traverser and accessor types passed to the constructor are the same as the corresponding
01734            internal types, we need not copy the image (speed up)
01735         */
01736     SplineImageView1(InternalIndexer const & i)
01737     : Base(i.shape(0), i.shape(1), i)
01738     {}
01739 
01740     template<class T, class SU>
01741     SplineImageView1(MultiArrayView<2, T, SU> const & i)
01742     : Base(i.shape(0), i.shape(1)),
01743       image_(i.shape(0), i.shape(1))
01744     {
01745         for(unsigned int y=0; y<this->height(); ++y)
01746             for(unsigned int x=0; x<this->width(); ++x)
01747                 image_(x,y) = detail::RequiresExplicitCast<VALUETYPE>::cast(i(x,y));
01748         this->internalIndexer_ = InternalIndexer(typename InternalIndexer::difference_type(this->width(), this->height()),
01749                                                  image_.data());
01750     }
01751 
01752     template <class SrcIterator, class SrcAccessor>
01753     SplineImageView1(SrcIterator is, SrcIterator iend, SrcAccessor sa)
01754     : Base(iend.x - is.x, iend.y - is.y),
01755       image_(iend-is)
01756     {
01757         copyImage(srcIterRange(is, iend, sa), destImage(image_));
01758         this->internalIndexer_ = InternalIndexer(typename InternalIndexer::difference_type(this->width(), this->height()),
01759                                                  image_.data());
01760     }
01761 
01762     template <class SrcIterator, class SrcAccessor>
01763     SplineImageView1(triple<SrcIterator, SrcIterator, SrcAccessor> s)
01764     : Base(s.second.x - s.first.x, s.second.y - s.first.y),
01765       image_(s.second - s.first)
01766     {
01767         copyImage(s, destImage(image_));
01768         this->internalIndexer_ = InternalIndexer(typename InternalIndexer::difference_type(this->width(), this->height()),
01769                                                  image_.data());
01770     }
01771 
01772     InternalImage const & image() const
01773         { return image_; }
01774 
01775   protected:
01776     InternalImage image_;
01777 };
01778 
01779 template <class VALUETYPE>
01780 class SplineImageView<1, VALUETYPE>
01781 : public SplineImageView1<VALUETYPE>
01782 {
01783     typedef SplineImageView1<VALUETYPE> Base;
01784   public:
01785     typedef typename Base::value_type value_type;
01786     typedef typename Base::size_type size_type;
01787     typedef typename Base::difference_type difference_type;
01788     enum StaticOrder { order = Base::order };
01789     typedef typename Base::InternalImage InternalImage;
01790 
01791   protected:
01792     typedef typename Base::InternalTraverser InternalTraverser;
01793     typedef typename Base::InternalAccessor InternalAccessor;
01794     typedef typename Base::InternalConstTraverser InternalConstTraverser;
01795     typedef typename Base::InternalConstAccessor InternalConstAccessor;
01796 
01797 public:
01798 
01799         /* when traverser and accessor types passed to the constructor are the same as the corresponding
01800            internal types, we need not copy the image (speed up)
01801         */
01802     SplineImageView(InternalTraverser is, InternalTraverser iend, InternalAccessor sa, bool /* unused */ = false)
01803     : Base(is, iend, sa)
01804     {}
01805 
01806     SplineImageView(triple<InternalTraverser, InternalTraverser, InternalAccessor> s, bool /* unused */ = false)
01807     : Base(s)
01808     {}
01809 
01810     SplineImageView(InternalConstTraverser is, InternalConstTraverser iend, InternalConstAccessor sa, bool /* unused */ = false)
01811     : Base(is, iend, sa)
01812     {}
01813 
01814     SplineImageView(triple<InternalConstTraverser, InternalConstTraverser, InternalConstAccessor> s, bool /* unused */ = false)
01815     : Base(s)
01816     {}
01817 
01818     template <class SrcIterator, class SrcAccessor>
01819     SplineImageView(SrcIterator is, SrcIterator iend, SrcAccessor sa, bool /* unused */ = false)
01820     : Base(is, iend, sa)
01821     {
01822         copyImage(srcIterRange(is, iend, sa), destImage(this->image_));
01823     }
01824 
01825     template <class SrcIterator, class SrcAccessor>
01826     SplineImageView(triple<SrcIterator, SrcIterator, SrcAccessor> s, bool /* unused */ = false)
01827     : Base(s)
01828     {
01829         copyImage(s, destImage(this->image_));
01830     }
01831 };
01832 
01833 } // namespace vigra
01834 
01835 
01836 #endif /* VIGRA_SPLINEIMAGEVIEW_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.8.0 (20 Sep 2011)