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

details vigra/splineimageview.hxx VIGRA

00001 /************************************************************************/
00002 /*                                                                      */
00003 /*               Copyright 1998-2004 by Ullrich Koethe                  */
00004 /*       Cognitive Systems Group, University of Hamburg, Germany        */
00005 /*                                                                      */
00006 /*    This file is part of the VIGRA computer vision library.           */
00007 /*    ( Version 1.3.3, Aug 18 2005 )                                    */
00008 /*    You may use, modify, and distribute this software according       */
00009 /*    to the terms stated in the LICENSE file included in               */
00010 /*    the VIGRA distribution.                                           */
00011 /*                                                                      */
00012 /*    The VIGRA Website is                                              */
00013 /*        http://kogs-www.informatik.uni-hamburg.de/~koethe/vigra/      */
00014 /*    Please direct questions, bug reports, and contributions to        */
00015 /*        koethe@informatik.uni-hamburg.de                              */
00016 /*                                                                      */
00017 /*  THIS SOFTWARE IS PROVIDED AS IS AND WITHOUT ANY EXPRESS OR          */
00018 /*  IMPLIED WARRANTIES, INCLUDING, WITHOUT LIMITATION, THE IMPLIED      */
00019 /*  WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. */
00020 /*                                                                      */
00021 /************************************************************************/
00022 
00023 #ifndef VIGRA_SPLINEIMAGEVIEW_HXX
00024 #define VIGRA_SPLINEIMAGEVIEW_HXX
00025 
00026 #include "vigra/mathutil.hxx"
00027 #include "vigra/recursiveconvolution.hxx"
00028 #include "vigra/splines.hxx"
00029 #include "vigra/array_vector.hxx"
00030 #include "vigra/basicimage.hxx"
00031 #include "vigra/copyimage.hxx"
00032 
00033 namespace vigra {
00034 
00035 /********************************************************/
00036 /*                                                      */
00037 /*                    SplineImageView                   */
00038 /*                                                      */
00039 /********************************************************/
00040 /** \brief Create a continuous view onto a discrete image using splines.
00041 
00042     This class is very useful if image values or derivatives at arbitrary
00043     real-valued coordinates are needed. Access at such coordinates is implemented by 
00044     interpolating the given discrete image values with a spline of the 
00045     specified <tt>ORDER</TT>. Continuous derivatives are available up to 
00046     degree <tt>ORDER-1</TT>. If the requested coordinates are near the image border, 
00047     reflective boundary conditions are applied. In principle, this class can also be used 
00048     for image resizing, but here the functions from the <tt>resize...</tt> family are 
00049     more efficient. 
00050     
00051 */
00052 template <int ORDER, class VALUETYPE>
00053 class SplineImageView
00054 {
00055     typedef typename NumericTraits<VALUETYPE>::RealPromote InternalValue;
00056 
00057   public:
00058   
00059         /** The view's value type (return type of access and derivative functions).
00060         */
00061     typedef VALUETYPE value_type;
00062     
00063         /** The view's size type.
00064         */
00065     typedef Size2D size_type;
00066 
00067         /** The view's difference type.
00068         */
00069     typedef TinyVector<double, 2> difference_type;
00070 
00071         /** The order of the spline used.
00072         */
00073     enum StaticOrder { order = ORDER };
00074     
00075         /** The type of the internal image holding the spline coefficients.
00076         */
00077     typedef BasicImage<InternalValue> InternalImage;
00078   
00079   private:
00080     typedef typename InternalImage::traverser InternalTraverser;
00081     typedef typename InternalTraverser::row_iterator InternalRowIterator;
00082     typedef typename InternalTraverser::column_iterator InternalColumnIterator;
00083     typedef BSpline<ORDER, double> Spline;
00084     
00085     enum { ksize_ = ORDER + 1, kcenter_ = ORDER / 2 };
00086  
00087   public:   
00088         /** Construct SplineImageView for the given image.
00089         
00090             If <tt>skipPrefiltering = true</tt> (default: <tt>false</tt>), the recursive
00091             prefilter of the cardinal spline function is not applied, resulting
00092             in an approximating (smoothing) rather than interpolating spline. This is
00093             especially useful if customized prefilters are to be applied.
00094         */
00095     template <class SrcIterator, class SrcAccessor>
00096     SplineImageView(SrcIterator is, SrcIterator iend, SrcAccessor sa, bool skipPrefiltering = false)
00097     : w_(iend.x - is.x), h_(iend.y - is.y), w1_(w_-1), h1_(h_-1),
00098       x0_(kcenter_), x1_(w_ - kcenter_ - 2), y0_(kcenter_), y1_(h_ - kcenter_ - 2),
00099       image_(w_, h_),
00100       x_(-1.0), y_(-1.0),
00101       u_(-1.0), v_(-1.0)
00102     {
00103         copyImage(srcIterRange(is, iend, sa), destImage(image_));
00104         if(!skipPrefiltering)
00105             init();
00106     }
00107     
00108         /** Construct SplineImageView for the given image.
00109         
00110             If <tt>skipPrefiltering = true</tt> (default: <tt>false</tt>), the recursive
00111             prefilter of the cardinal spline function is not applied, resulting
00112             in an approximating (smoothing) rather than interpolating spline. This is
00113             especially useful if customized prefilters are to be applied.
00114         */
00115     template <class SrcIterator, class SrcAccessor>
00116     SplineImageView(triple<SrcIterator, SrcIterator, SrcAccessor> s, bool skipPrefiltering = false)
00117     : w_(s.second.x - s.first.x), h_(s.second.y - s.first.y), w1_(w_-1), h1_(h_-1),
00118       x0_(kcenter_), x1_(w_ - kcenter_ - 2), y0_(kcenter_), y1_(h_ - kcenter_ - 2),
00119       image_(w_, h_),
00120       x_(-1.0), y_(-1.0),
00121       u_(-1.0), v_(-1.0)
00122     {
00123         copyImage(srcIterRange(s.first, s.second, s.third), destImage(image_));
00124         if(!skipPrefiltering)
00125             init();
00126     }
00127     
00128         /** Access interpolated function at real-valued coordinate <tt>(x, y)</tt>.
00129         */
00130     value_type operator()(double x, double y) const;
00131     
00132         /** Access derivative of order <tt>(dx, dy)</tt> at real-valued coordinate <tt>(x, y)</tt>.
00133         */
00134     value_type operator()(double x, double y, unsigned int dx, unsigned int dy) const;
00135     
00136         /** Access 1st derivative in x-direction at real-valued coordinate <tt>(x, y)</tt>.
00137             Equivalent to <tt>splineView(x, y, 1, 0)</tt>.
00138         */
00139     value_type dx(double x, double y) const
00140         { return operator()(x, y, 1, 0); }
00141     
00142         /** Access 1st derivative in y-direction at real-valued coordinate <tt>(x, y)</tt>.
00143             Equivalent to <tt>splineView(x, y, 0, 1)</tt>.
00144         */
00145     value_type dy(double x, double y) const
00146         { return operator()(x, y, 0, 1); }
00147     
00148         /** Access 2nd derivative in x-direction at real-valued coordinate <tt>(x, y)</tt>.
00149             Equivalent to <tt>splineView(x, y, 2, 0)</tt>.
00150         */
00151     value_type dxx(double x, double y) const
00152         { return operator()(x, y, 2, 0); }
00153     
00154         /** Access mixed 2nd derivative at real-valued coordinate <tt>(x, y)</tt>.
00155             Equivalent to <tt>splineView(x, y, 1, 1)</tt>.
00156         */
00157     value_type dxy(double x, double y) const
00158         { return operator()(x, y, 1, 1); }
00159     
00160         /** Access 2nd derivative in y-direction at real-valued coordinate <tt>(x, y)</tt>.
00161             Equivalent to <tt>splineView(x, y, 0, 2)</tt>.
00162         */
00163     value_type dyy(double x, double y) const
00164         { return operator()(x, y, 0, 2); }
00165     
00166         /** Access 3rd derivative in x-direction at real-valued coordinate <tt>(x, y)</tt>.
00167             Equivalent to <tt>splineView(x, y, 3, 0)</tt>.
00168         */
00169     value_type dx3(double x, double y) const
00170         { return operator()(x, y, 3, 0); }
00171     
00172         /** Access 3rd derivative in y-direction at real-valued coordinate <tt>(x, y)</tt>.
00173             Equivalent to <tt>splineView(x, y, 0, 3)</tt>.
00174         */
00175     value_type dy3(double x, double y) const
00176         { return operator()(x, y, 0, 3); }
00177     
00178         /** Access mixed 3rd derivative dxxy at real-valued coordinate <tt>(x, y)</tt>.
00179             Equivalent to <tt>splineView(x, y, 2, 1)</tt>.
00180         */
00181     value_type dxxy(double x, double y) const
00182         { return operator()(x, y, 2, 1); }
00183     
00184         /** Access mixed 3rd derivative dxyy at real-valued coordinate <tt>(x, y)</tt>.
00185             Equivalent to <tt>splineView(x, y, 1, 2)</tt>.
00186         */
00187     value_type dxyy(double x, double y) const
00188         { return operator()(x, y, 1, 2); }
00189         
00190         /** Access interpolated function at real-valued coordinate <tt>d</tt>.
00191             Equivalent to <tt>splineView(d[0], d[1])</tt>.
00192         */
00193     value_type operator()(difference_type const & d) const
00194         { return operator()(d[0], d[1]); }
00195         
00196         /** Access derivative of order <tt>(dx, dy)</tt> at real-valued coordinate <tt>d</tt>.
00197             Equivalent to <tt>splineView(d[0], d[1], dx, dy)</tt>.
00198         */
00199     value_type operator()(difference_type const & d, unsigned int dx, unsigned int dy) const
00200         { return operator()(d[0], d[1], dx, dy); }
00201         
00202         /** Access 1st derivative in x-direction at real-valued coordinate <tt>d</tt>.
00203             Equivalent to <tt>splineView.dx(d[0], d[1])</tt>.
00204         */
00205     value_type dx(difference_type const & d) const
00206         { return dx(d[0], d[1]); }
00207         
00208         /** Access 1st derivative in y-direction at real-valued coordinate <tt>d</tt>.
00209             Equivalent to <tt>splineView.dy(d[0], d[1])</tt>.
00210         */
00211     value_type dy(difference_type const & d) const
00212         { return dy(d[0], d[1]); }
00213         
00214         /** Access 2nd derivative in x-direction at real-valued coordinate <tt>d</tt>.
00215             Equivalent to <tt>splineView.dxx(d[0], d[1])</tt>.
00216         */
00217     value_type dxx(difference_type const & d) const
00218         { return dxx(d[0], d[1]); }
00219         
00220         /** Access mixed 2nd derivative at real-valued coordinate <tt>d</tt>.
00221             Equivalent to <tt>splineView.dxy(d[0], d[1])</tt>.
00222         */
00223     value_type dxy(difference_type const & d) const
00224         { return dxy(d[0], d[1]); }
00225         
00226         /** Access 2nd derivative in y-direction at real-valued coordinate <tt>d</tt>.
00227             Equivalent to <tt>splineView.dyy(d[0], d[1])</tt>.
00228         */
00229     value_type dyy(difference_type const & d) const
00230         { return dyy(d[0], d[1]); }
00231         
00232         /** Access 3rd derivative in x-direction at real-valued coordinate <tt>d</tt>.
00233             Equivalent to <tt>splineView.dx3(d[0], d[1])</tt>.
00234         */
00235     value_type dx3(difference_type const & d) const
00236         { return dx3(d[0], d[1]); }
00237         
00238         /** Access 3rd derivative in y-direction at real-valued coordinate <tt>d</tt>.
00239             Equivalent to <tt>splineView.dy3(d[0], d[1])</tt>.
00240         */
00241     value_type dy3(difference_type const & d) const
00242         { return dy3(d[0], d[1]); }
00243         
00244         /** Access mixed 3rd derivative dxxy at real-valued coordinate <tt>d</tt>.
00245             Equivalent to <tt>splineView.dxxy(d[0], d[1])</tt>.
00246         */
00247     value_type dxxy(difference_type const & d) const
00248         { return dxxy(d[0], d[1]); }
00249         
00250         /** Access mixed 3rd derivative dxyy at real-valued coordinate <tt>d</tt>.
00251             Equivalent to <tt>splineView.dxyy(d[0], d[1])</tt>.
00252         */
00253     value_type dxyy(difference_type const & d) const
00254         { return dxyy(d[0], d[1]); }
00255         
00256         /** Access gradient squared magnitude at real-valued coordinate <tt>(x, y)</tt>.
00257         */
00258     value_type g2(double x, double y) const;
00259         
00260         /** Access 1st derivative in x-direction of gradient squared magnitude 
00261             at real-valued coordinate <tt>(x, y)</tt>.
00262         */
00263     value_type g2x(double x, double y) const;
00264         
00265         /** Access 1st derivative in y-direction of gradient squared magnitude 
00266             at real-valued coordinate <tt>(x, y)</tt>.
00267         */
00268     value_type g2y(double x, double y) const;
00269         
00270         /** Access 2nd derivative in x-direction of gradient squared magnitude 
00271             at real-valued coordinate <tt>(x, y)</tt>.
00272         */
00273     value_type g2xx(double x, double y) const;
00274         
00275         /** Access mixed 2nd derivative of gradient squared magnitude 
00276             at real-valued coordinate <tt>(x, y)</tt>.
00277         */
00278     value_type g2xy(double x, double y) const;
00279         
00280         /** Access 2nd derivative in y-direction of gradient squared magnitude 
00281             at real-valued coordinate <tt>(x, y)</tt>.
00282         */
00283     value_type g2yy(double x, double y) const;
00284         
00285         /** Access gradient squared magnitude at real-valued coordinate <tt>d</tt>.
00286         */
00287     value_type g2(difference_type const & d) const
00288         { return g2(d[0], d[1]); }
00289         
00290         /** Access 1st derivative in x-direction of gradient squared magnitude 
00291             at real-valued coordinate <tt>d</tt>.
00292         */
00293     value_type g2x(difference_type const & d) const
00294         { return g2x(d[0], d[1]); }
00295         
00296         /** Access 1st derivative in y-direction of gradient squared magnitude 
00297             at real-valued coordinate <tt>d</tt>.
00298         */
00299     value_type g2y(difference_type const & d) const
00300         { return g2y(d[0], d[1]); }
00301         
00302         /** Access 2nd derivative in x-direction of gradient squared magnitude 
00303             at real-valued coordinate <tt>d</tt>.
00304         */
00305     value_type g2xx(difference_type const & d) const
00306         { return g2xx(d[0], d[1]); }
00307         
00308         /** Access mixed 2nd derivative of gradient squared magnitude 
00309             at real-valued coordinate <tt>d</tt>.
00310         */
00311     value_type g2xy(difference_type const & d) const
00312         { return g2xy(d[0], d[1]); }
00313         
00314         /** Access 2nd derivative in y-direction of gradient squared magnitude 
00315             at real-valued coordinate <tt>d</tt>.
00316         */
00317     value_type g2yy(difference_type const & d) const
00318         { return g2yy(d[0], d[1]); }
00319     
00320         /** The width of the image.
00321             <tt>0 &lt;= x &lt;= width()-1</tt> is required for all access functions.
00322         */
00323     unsigned int width() const
00324         { return w_; }
00325     
00326         /** The height of the image.
00327             <tt>0 &lt;= y &lt;= height()-1</tt> is required for all access functions.
00328         */
00329     unsigned int height() const
00330         { return h_; }
00331     
00332         /** The size of the image.
00333             <tt>0 &lt;= x &lt;= size().x-1</tt> and <tt>0 &lt;= y &lt;= size().y-1</tt> 
00334             are required for all access functions.
00335         */
00336     size_type size() const
00337         { return size_type(w_, h_); }
00338         
00339         /** The internal image holding the spline coefficients.
00340         */
00341     InternalImage const & image() const
00342     {
00343         return image_;
00344     }
00345     
00346         /** Get the array of polynomial coefficients for the facet containing 
00347             the point <tt>(x, y)</tt>. The array <tt>res</tt> will be resized to
00348             dimension <tt>(ORDER+1)x(ORDER+1)</tt>. From these coefficients, the
00349             value of the interpolated function can be calculated by the following
00350             algorithm
00351             
00352             \code
00353             SplineImageView<ORDER, float> view(...);
00354             double x = ..., y = ...;
00355             double dx, dy;
00356             
00357             // calculate the local facet coordinates of x and y
00358             if(ORDER % 2)
00359             {
00360                 // odd order => facet coordinates between 0 and 1
00361                 dx = x - floor(x);
00362                 dy = y - floor(y);
00363             }
00364             else
00365             {
00366                 // even order => facet coordinates between -0.5 and 0.5
00367                 dx = x - floor(x + 0.5);
00368                 dy = y - floor(y + 0.5);
00369             }
00370             
00371             BasicImage<float> coefficients;
00372             view.coefficientArray(x, y, coefficients);
00373             
00374             float f_x_y = 0.0;
00375             for(int ny = 0; ny < ORDER + 1; ++ny)
00376                 for(int nx = 0; nx < ORDER + 1; ++nx)
00377                     f_x_y += pow(dx, nx) * pow(dy, ny) * coefficients(nx, ny);
00378                     
00379             assert(abs(f_x_y - view(x, y)) < 1e-6);
00380             \endcode
00381         */
00382     template <class Array>
00383     void coefficientArray(double x, double y, Array & res) const;
00384     
00385         /** Check if x is in the valid range.
00386             Equivalent to <tt>0 &lt;= x &lt;= width()-1</tt>.
00387         */
00388     bool isInsideX(double x) const
00389     {
00390         return x >= 0.0 && x <= width()-1.0;
00391     }
00392         
00393         /** Check if y is in the valid range.
00394             Equivalent to <tt>0 &lt;= y &lt;= height()-1</tt>.
00395         */
00396     bool isInsideY(double y) const
00397     {
00398         return y >= 0.0 && y <= height()-1.0;
00399     }
00400         
00401         /** Check if x and y are in the valid range.
00402             Equivalent to <tt>0 &lt;= x &lt;= width()-1</tt> and <tt>0 &lt;= y &lt;= height()-1</tt>.
00403         */
00404     bool isInside(double x, double y) const
00405     {
00406         return isInsideX(x) && isInsideY(y);
00407     }
00408     
00409         /** Check whether the points <tt>(x0, y0)</tt> and <tt>(x1, y1)</tt> are in
00410             the same spline facet. For odd order splines, facets span the range
00411             <tt>(floor(x), floor(x)+1) x (floor(y), floor(y)+1)</tt> (i.e. we have 
00412             integer facet borders), whereas even order splines have facet between
00413             half integer values 
00414             <tt>(floor(x)-0.5, floor(x)+0.5) x (floor(y)-0.5, floor(y)+0.5)</tt>.
00415         */
00416     bool sameFacet(double x0, double y0, double x1, double y1) const
00417     {
00418          x0 = VIGRA_CSTD::floor((ORDER % 2) ? x0 : x0 + 0.5);
00419          y0 = VIGRA_CSTD::floor((ORDER % 2) ? y0 : y0 + 0.5);
00420          x1 = VIGRA_CSTD::floor((ORDER % 2) ? x1 : x1 + 0.5);
00421          y1 = VIGRA_CSTD::floor((ORDER % 2) ? y1 : y1 + 0.5);
00422          return x0 == x1 && y0 == y1;
00423     }
00424         
00425   protected:
00426   
00427     void init();
00428     void calculateIndices(double x, double y) const;
00429     void coefficients(double t, double * const & c) const;
00430     void derivCoefficients(double t, unsigned int d, double * const & c) const;
00431     value_type convolve() const;
00432   
00433     unsigned int w_, h_;
00434     int w1_, h1_;
00435     double x0_, x1_, y0_, y1_;
00436     InternalImage image_;
00437     Spline k_;
00438     mutable double x_, y_, u_, v_, kx_[ksize_], ky_[ksize_];
00439     mutable int ix_[ksize_], iy_[ksize_];
00440 };
00441 
00442 template <int ORDER, class VALUETYPE>
00443 void SplineImageView<ORDER, VALUETYPE>::init()
00444 {
00445     ArrayVector<double> const & b = k_.prefilterCoefficients();
00446     
00447     for(unsigned int i=0; i<b.size(); ++i)
00448     {
00449         recursiveFilterX(srcImageRange(image_), destImage(image_), b[i], BORDER_TREATMENT_REFLECT);
00450         recursiveFilterY(srcImageRange(image_), destImage(image_), b[i], BORDER_TREATMENT_REFLECT);
00451     }
00452 }
00453 
00454 namespace detail
00455 {
00456 
00457 template <int i>
00458 struct SplineImageViewUnrollLoop1
00459 {
00460     template <class Array>
00461     static void exec(int c0, Array c)
00462     {
00463         SplineImageViewUnrollLoop1<i-1>::exec(c0, c);
00464         c[i] = c0 + i;
00465     }
00466 };
00467 
00468 template <>
00469 struct SplineImageViewUnrollLoop1<0>
00470 {
00471     template <class Array>
00472     static void exec(int c0, Array c)
00473     {
00474         c[0] = c0;
00475     }
00476 };
00477 
00478 template <int i>
00479 struct SplineImageViewUnrollLoop2
00480 {
00481     template <class Array1, class Image, class Array2>
00482     static typename Image::value_type
00483     exec(Array1 k, Image const & img, Array2 x, int y)
00484     {
00485         return k[i] * img(x[i], y) + SplineImageViewUnrollLoop2<i-1>::exec(k, img, x, y);
00486     }
00487 };
00488 
00489 template <>
00490 struct SplineImageViewUnrollLoop2<0>
00491 {
00492     template <class Array1, class Image, class Array2>
00493     static typename Image::value_type
00494     exec(Array1 k, Image const & img, Array2 x, int y)
00495     {
00496         return k[0] * img(x[0], y);
00497     }
00498 };
00499 
00500 } // namespace detail
00501 
00502 template <int ORDER, class VALUETYPE>
00503 void 
00504 SplineImageView<ORDER, VALUETYPE>::calculateIndices(double x, double y) const
00505 {
00506     if(x == x_ && y == y_)
00507         return;   // still in cache
00508     
00509     if(x > x0_ && x < x1_ && y > y0_ && y < y1_)
00510     {
00511         detail::SplineImageViewUnrollLoop1<ORDER>::exec(
00512                                 (ORDER % 2) ? int(x - kcenter_) : int(x + 0.5 - kcenter_), ix_);
00513         detail::SplineImageViewUnrollLoop1<ORDER>::exec(
00514                                 (ORDER % 2) ? int(y - kcenter_) : int(y + 0.5 - kcenter_), iy_);
00515     }
00516     else
00517     {
00518         vigra_precondition(isInside(x, y),
00519              "SplineImageView<ORDER, VALUETYPE>::calculateIndices(): index out of bounds.");
00520 
00521         ix_[kcenter_] = (ORDER % 2) ?
00522                              (int)x :
00523                              (int)(x + 0.5);
00524         iy_[kcenter_] = (ORDER % 2) ?
00525                              (int)y :
00526                              (int)(y + 0.5);
00527         
00528         for(int i=0; i<kcenter_; ++i)
00529         {
00530             ix_[i] = vigra::abs(ix_[kcenter_] - (kcenter_ - i));
00531             iy_[i] = vigra::abs(iy_[kcenter_] - (kcenter_ - i));
00532         }
00533         for(int i=kcenter_+1; i<ksize_; ++i)
00534         {
00535             ix_[i] = w1_ - vigra::abs(w1_ - ix_[kcenter_] - (i - kcenter_));
00536             iy_[i] = h1_ - vigra::abs(h1_ - iy_[kcenter_] - (i - kcenter_));
00537         }
00538     }
00539     x_ = x;
00540     y_ = y;
00541     u_ = x - ix_[kcenter_];
00542     v_ = y - iy_[kcenter_];
00543 }
00544 
00545 template <int ORDER, class VALUETYPE>
00546 void SplineImageView<ORDER, VALUETYPE>::coefficients(double t, double * const & c) const
00547 {
00548     t += kcenter_;
00549     for(int i = 0; i<ksize_; ++i)
00550         c[i] = k_(t-i);
00551 }
00552 
00553 template <int ORDER, class VALUETYPE>
00554 void SplineImageView<ORDER, VALUETYPE>::derivCoefficients(double t, 
00555                                                unsigned int d, double * const & c) const
00556 {
00557     t += kcenter_;
00558     for(int i = 0; i<ksize_; ++i)
00559         c[i] = k_(t-i, d);
00560 }
00561 
00562 
00563 template <int ORDER, class VALUETYPE>
00564 VALUETYPE SplineImageView<ORDER, VALUETYPE>::convolve() const
00565 {
00566     InternalValue sum;
00567     sum = ky_[0]*detail::SplineImageViewUnrollLoop2<ORDER>::exec(kx_, image_, ix_, iy_[0]);
00568     
00569     for(int j=1; j<ksize_; ++j)
00570     {
00571         sum += ky_[j]*detail::SplineImageViewUnrollLoop2<ORDER>::exec(kx_, image_, ix_, iy_[j]);
00572     }
00573     return NumericTraits<VALUETYPE>::fromRealPromote(sum);
00574 }
00575 
00576 template <int ORDER, class VALUETYPE>
00577 template <class Array>
00578 void 
00579 SplineImageView<ORDER, VALUETYPE>::coefficientArray(double x, double y, Array & res) const
00580 {
00581     typename Spline::WeightMatrix & weights = Spline::weights();
00582     InternalValue tmp[ksize_][ksize_]; 
00583     
00584     calculateIndices(x, y);
00585     for(int j=0; j<ksize_; ++j)
00586     {
00587         for(int i=0; i<ksize_; ++i)
00588         {
00589             tmp[i][j] = 0.0;
00590             for(int k=0; k<ksize_; ++k)
00591             {
00592                 tmp[i][j] += weights[i][k]*image_(ix_[k], iy_[j]);
00593             }
00594        }       
00595     }
00596     res.resize(ksize_, ksize_);
00597     for(int j=0; j<ksize_; ++j)
00598     {
00599         for(int i=0; i<ksize_; ++i)
00600         {
00601             res(i,j) = 0.0;
00602             for(int k=0; k<ksize_; ++k)
00603             {
00604                 res(i,j) += weights[j][k]*tmp[i][k];
00605             }
00606         }       
00607     }
00608 }
00609 
00610 template <int ORDER, class VALUETYPE>
00611 VALUETYPE SplineImageView<ORDER, VALUETYPE>::operator()(double x, double y) const
00612 {
00613     calculateIndices(x, y);
00614     coefficients(u_, kx_);
00615     coefficients(v_, ky_);
00616     return convolve();
00617 }
00618 
00619 template <int ORDER, class VALUETYPE>
00620 VALUETYPE SplineImageView<ORDER, VALUETYPE>::operator()(double x, double y,
00621                                                  unsigned int dx, unsigned int dy) const
00622 {
00623     calculateIndices(x, y);
00624     derivCoefficients(u_, dx, kx_);
00625     derivCoefficients(v_, dy, ky_);
00626     return convolve();
00627 }
00628 
00629 template <int ORDER, class VALUETYPE>
00630 VALUETYPE SplineImageView<ORDER, VALUETYPE>::g2(double x, double y) const
00631 {
00632     return sq(dx(x,y)) + sq(dy(x,y));
00633 }
00634 
00635 template <int ORDER, class VALUETYPE>
00636 VALUETYPE SplineImageView<ORDER, VALUETYPE>::g2x(double x, double y) const
00637 {
00638     return 2.0*(dx(x,y) * dxx(x,y) + dy(x,y) * dxy(x,y));
00639 }
00640 
00641 template <int ORDER, class VALUETYPE>
00642 VALUETYPE SplineImageView<ORDER, VALUETYPE>::g2y(double x, double y) const
00643 {
00644     return 2.0*(dx(x,y) * dxy(x,y) + dy(x,y) * dyy(x,y));
00645 }
00646 
00647 template <int ORDER, class VALUETYPE>
00648 VALUETYPE SplineImageView<ORDER, VALUETYPE>::g2xx(double x, double y) const
00649 {
00650     return 2.0*(sq(dxx(x,y)) + dx(x,y) * dx3(x,y) + sq(dxy(x,y)) + dy(x,y) * dxxy(x,y));
00651 }
00652 
00653 template <int ORDER, class VALUETYPE>
00654 VALUETYPE SplineImageView<ORDER, VALUETYPE>::g2yy(double x, double y) const
00655 {
00656     return 2.0*(sq(dxy(x,y)) + dx(x,y) * dxyy(x,y) + sq(dyy(x,y)) + dy(x,y) * dy3(x,y));
00657 }
00658 
00659 template <int ORDER, class VALUETYPE>
00660 VALUETYPE SplineImageView<ORDER, VALUETYPE>::g2xy(double x, double y) const
00661 {
00662     return 2.0*(dx(x,y) * dxxy(x,y) + dy(x,y) * dxyy(x,y) + dxy(x,y) * (dxx(x,y) + dyy(x,y)));
00663 }
00664 
00665 } // namespace vigra
00666 
00667 
00668 #endif /* VIGRA_SPLINEIMAGEVIEW_HXX */

© Ullrich Köthe (koethe@informatik.uni-hamburg.de)
Cognitive Systems Group, University of Hamburg, Germany

html generated using doxygen and Python
VIGRA 1.3.3 (18 Aug 2005)