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

vigra/resizeimage.hxx

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 /*    The VIGRA Website is                                              */
00008 /*        http://kogs-www.informatik.uni-hamburg.de/~koethe/vigra/      */
00009 /*    Please direct questions, bug reports, and contributions to        */
00010 /*        ullrich.koethe@iwr.uni-heidelberg.de    or                    */
00011 /*        vigra@informatik.uni-hamburg.de                               */
00012 /*                                                                      */
00013 /*    Permission is hereby granted, free of charge, to any person       */
00014 /*    obtaining a copy of this software and associated documentation    */
00015 /*    files (the "Software"), to deal in the Software without           */
00016 /*    restriction, including without limitation the rights to use,      */
00017 /*    copy, modify, merge, publish, distribute, sublicense, and/or      */
00018 /*    sell copies of the Software, and to permit persons to whom the    */
00019 /*    Software is furnished to do so, subject to the following          */
00020 /*    conditions:                                                       */
00021 /*                                                                      */
00022 /*    The above copyright notice and this permission notice shall be    */
00023 /*    included in all copies or substantial portions of the             */
00024 /*    Software.                                                         */
00025 /*                                                                      */
00026 /*    THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND    */
00027 /*    EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES   */
00028 /*    OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND          */
00029 /*    NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT       */
00030 /*    HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,      */
00031 /*    WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING      */
00032 /*    FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR     */
00033 /*    OTHER DEALINGS IN THE SOFTWARE.                                   */
00034 /*                                                                      */
00035 /************************************************************************/
00036 
00037 
00038 #ifndef VIGRA_RESIZEIMAGE_HXX
00039 #define VIGRA_RESIZEIMAGE_HXX
00040 
00041 #include <vector>
00042 #include "utilities.hxx"
00043 #include "numerictraits.hxx"
00044 #include "stdimage.hxx"
00045 #include "recursiveconvolution.hxx"
00046 #include "separableconvolution.hxx"
00047 #include "resampling_convolution.hxx"
00048 #include "splines.hxx"
00049 
00050 namespace vigra {
00051 
00052 /*****************************************************************/
00053 /*                                                               */
00054 /*                         CoscotFunction                        */
00055 /*                                                               */
00056 /*****************************************************************/
00057 
00058 /*! The Coscot interpolation function.
00059 
00060     Implements the Coscot interpolation function proposed by Maria Magnusson Seger
00061     (maria@isy.liu.se) in the context of tomographic reconstruction. It provides a fast
00062     transition between the pass- and stop-bands and minimal ripple outside the transition
00063     region. Both properties are important for this application and can be tuned by the parameters
00064     <i>m</i> and <i>h</i> (with defaults 3 and 0.5). The function is defined by
00065 
00066     \f[ f_{m,h}(x) = \left\{ \begin{array}{ll}
00067                                    \frac{1}{2m}\sin(\pi x)\cot(\pi x / (2 m))(h + (1-h)\cos(\pi x/m)) & |x| \leq m \\
00068                                   0 & \mbox{otherwise}
00069                         \end{array}\right.
00070     \f]
00071 
00072     It can be used as a functor, and as a kernel for
00073     \ref resamplingConvolveImage() to create a differentiable interpolant
00074     of an image.
00075 
00076     <b>\#include</b> <<a href="resizeimage_8hxx-source.html">vigra/resizeimage.hxx</a>><br>
00077     Namespace: vigra
00078 
00079     \ingroup MathFunctions
00080 */
00081 template <class T>
00082 class CoscotFunction
00083 {
00084   public:
00085 
00086         /** the kernel's value type
00087         */
00088     typedef T            value_type;
00089         /** the unary functor's argument type
00090         */
00091     typedef T            argument_type;
00092         /** the splines polynomial order
00093         */
00094     typedef T            result_type;
00095 
00096     CoscotFunction(unsigned int m = 3, double h = 0.5)
00097     : m_(m),
00098       h_(h)
00099     {}
00100 
00101         /** function (functor) call
00102         */
00103     result_type operator()(argument_type x) const
00104     {
00105         return x == 0.0 ?
00106                     1.0
00107                   : abs(x) < m_ ?
00108                         VIGRA_CSTD::sin(M_PI*x) / VIGRA_CSTD::tan(M_PI * x / 2.0 / m_) *
00109                              (h_ + (1.0 - h_) * VIGRA_CSTD::cos(M_PI * x / m_)) / 2.0 / m_
00110                       : 0.0;
00111     }
00112 
00113         /** index operator -- same as operator()
00114         */
00115     value_type operator[](value_type x) const
00116         { return operator()(x); }
00117 
00118         /** Radius of the function's support.
00119             Needed for  \ref resamplingConvolveImage(), equals m.
00120         */
00121     double radius() const
00122         { return m_; }
00123 
00124         /** Derivative order of the function: always 0.
00125         */
00126     unsigned int derivativeOrder() const
00127         { return 0; }
00128 
00129         /** Prefilter coefficients for compatibility with \ref vigra::BSpline.
00130             (array has zero length, since prefiltering is not necessary).
00131         */
00132     ArrayVector<double> const & prefilterCoefficients() const
00133     {
00134         static ArrayVector<double> b;
00135         return b;
00136     }
00137 
00138   protected:
00139 
00140     unsigned int m_;
00141     double h_;
00142 };
00143 
00144 /** \addtogroup GeometricTransformations Geometric Transformations
00145     Zoom up and down by repeating pixels, or using various interpolation schemes.
00146 
00147     See also: \ref resamplingConvolveImage(), \ref resampleImage(), \ref resizeMultiArraySplineInterpolation()
00148 
00149     <b>\#include</b> <<a href="stdimagefunctions_8hxx-source.html">vigra/stdimagefunctions.hxx</a>><br>
00150     <b>or</b><br>
00151     <b>\#include</b> <<a href="resizeimage_8hxx-source.html">vigra/resizeimage.hxx</a>><br>
00152 */
00153 //@{
00154 
00155 /********************************************************/
00156 /*                                                      */
00157 /*               resizeLineNoInterpolation              */
00158 /*                                                      */
00159 /********************************************************/
00160 
00161 template <class SrcIterator, class SrcAccessor,
00162           class DestIterator, class DestAccessor>
00163 void
00164 resizeLineNoInterpolation(SrcIterator i1, SrcIterator iend, SrcAccessor as,
00165                            DestIterator id, DestIterator idend, DestAccessor ad)
00166 {
00167     int wold = iend - i1;
00168     int wnew = idend - id;
00169 
00170     if((wold <= 1) || (wnew <= 1)) return; // oder error ?
00171 
00172     ad.set(as(i1), id);
00173     ++id;
00174 
00175     --iend, --idend;
00176     ad.set(as(iend), idend);
00177 
00178     double dx = (double)(wold - 1) / (wnew - 1);
00179     double x = dx;
00180 
00181     for(; id != idend; ++id, x += dx)
00182     {
00183     if(x >= 1.0)
00184     {
00185         int xx = (int)x;
00186         i1 += xx;
00187         x -= (double)xx;
00188     }
00189 
00190     ad.set(as(i1), id);
00191     }
00192 }
00193 
00194 /********************************************************/
00195 /*                                                      */
00196 /*              resizeImageNoInterpolation              */
00197 /*                                                      */
00198 /********************************************************/
00199 
00200 /** \brief Resize image by repeating the nearest pixel values.
00201 
00202     This algorithm is very fast and does not require any arithmetic on
00203     the pixel types.
00204 
00205     The range of both the input and output images (resp. regions) must
00206     be given. Both images must have a size of at least 2x2 pixels. The
00207     scaling factors are then calculated accordingly. Destination
00208     pixels are directly copied from the appropriate source pixels.
00209 
00210     The function uses accessors.
00211 
00212     <b> Declarations:</b>
00213 
00214     pass arguments explicitly:
00215     \code
00216     namespace vigra {
00217         template <class SrcImageIterator, class SrcAccessor,
00218                   class DestImageIterator, class DestAccessor>
00219         void
00220         resizeImageNoInterpolation(
00221               SrcImageIterator is, SrcImageIterator iend, SrcAccessor sa,
00222               DestImageIterator id, DestImageIterator idend, DestAccessor da)
00223     }
00224     \endcode
00225 
00226 
00227     use argument objects in conjunction with \ref ArgumentObjectFactories :
00228     \code
00229     namespace vigra {
00230         template <class SrcImageIterator, class SrcAccessor,
00231                   class DestImageIterator, class DestAccessor>
00232         void
00233         resizeImageNoInterpolation(
00234               triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src,
00235               triple<DestImageIterator, DestImageIterator, DestAccessor> dest)
00236     }
00237     \endcode
00238 
00239     <b> Usage:</b>
00240 
00241         <b>\#include</b> <<a href="resizeimage_8hxx-source.html">vigra/resizeimage.hxx</a>><br>
00242         Namespace: vigra
00243 
00244     \code
00245     vigra::resizeImageNoInterpolation(
00246                src.upperLeft(), src.lowerRight(), src.accessor(),
00247                dest.upperLeft(), dest.lowerRight(), dest.accessor());
00248 
00249     \endcode
00250 
00251     <b> Required Interface:</b>
00252 
00253     \code
00254     SrcImageIterator src_upperleft, src_lowerright;
00255     DestImageIterator dest_upperleft, src_lowerright;
00256 
00257     SrcAccessor src_accessor;
00258     DestAccessor dest_accessor;
00259 
00260     dest_accessor.set(src_accessor(src_upperleft), dest_upperleft);
00261 
00262     \endcode
00263 
00264     <b> Preconditions:</b>
00265 
00266     \code
00267     src_lowerright.x - src_upperleft.x > 1
00268     src_lowerright.y - src_upperleft.y > 1
00269     dest_lowerright.x - dest_upperleft.x > 1
00270     dest_lowerright.y - dest_upperleft.y > 1
00271     \endcode
00272 
00273 */
00274 doxygen_overloaded_function(template <...> void resizeImageNoInterpolation)
00275 
00276 template <class SrcIterator, class SrcAccessor,
00277           class DestIterator, class DestAccessor>
00278 void
00279 resizeImageNoInterpolation(SrcIterator is, SrcIterator iend, SrcAccessor sa,
00280                       DestIterator id, DestIterator idend, DestAccessor da)
00281 {
00282     int w = iend.x - is.x;
00283     int h = iend.y - is.y;
00284 
00285     int wnew = idend.x - id.x;
00286     int hnew = idend.y - id.y;
00287 
00288     vigra_precondition((w > 1) && (h > 1),
00289                  "resizeImageNoInterpolation(): "
00290                  "Source image to small.\n");
00291     vigra_precondition((wnew > 1) && (hnew > 1),
00292                  "resizeImageNoInterpolation(): "
00293                  "Destination image to small.\n");
00294 
00295     typedef BasicImage<typename SrcAccessor::value_type> TmpImage;
00296     typedef typename TmpImage::traverser TmpImageIterator;
00297 
00298     TmpImage tmp(w, hnew);
00299 
00300     TmpImageIterator yt = tmp.upperLeft();
00301 
00302     for(int x=0; x<w; ++x, ++is.x, ++yt.x)
00303     {
00304         typename SrcIterator::column_iterator c1 = is.columnIterator();
00305         typename TmpImageIterator::column_iterator ct = yt.columnIterator();
00306 
00307         resizeLineNoInterpolation(c1, c1 + h, sa, ct, ct + hnew, tmp.accessor());
00308     }
00309 
00310     yt = tmp.upperLeft();
00311 
00312     for(int y=0; y < hnew; ++y, ++yt.y, ++id.y)
00313     {
00314         typename DestIterator::row_iterator rd = id.rowIterator();
00315         typename TmpImageIterator::row_iterator rt = yt.rowIterator();
00316 
00317         resizeLineNoInterpolation(rt, rt + w, tmp.accessor(), rd, rd + wnew, da);
00318     }
00319 }
00320 
00321 template <class SrcIterator, class SrcAccessor,
00322           class DestIterator, class DestAccessor>
00323 inline
00324 void
00325 resizeImageNoInterpolation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
00326                            triple<DestIterator, DestIterator, DestAccessor> dest)
00327 {
00328     resizeImageNoInterpolation(src.first, src.second, src.third,
00329                                    dest.first, dest.second, dest.third);
00330 }
00331 
00332 /********************************************************/
00333 /*                                                      */
00334 /*             resizeLineLinearInterpolation            */
00335 /*                                                      */
00336 /********************************************************/
00337 
00338 template <class SrcIterator, class SrcAccessor,
00339           class DestIterator, class DestAccessor>
00340 void
00341 resizeLineLinearInterpolation(SrcIterator i1, SrcIterator iend, SrcAccessor as,
00342                            DestIterator id, DestIterator idend, DestAccessor ad)
00343 {
00344     int wold = iend - i1;
00345     int wnew = idend - id;
00346 
00347     if((wold <= 1) || (wnew <= 1)) return; // oder error ?
00348 
00349     typedef
00350         NumericTraits<typename DestAccessor::value_type> DestTraits;
00351 
00352     ad.set(DestTraits::fromRealPromote(as(i1)), id);
00353     ++id;
00354 
00355     --iend, --idend;
00356     ad.set(DestTraits::fromRealPromote(as(iend)), idend);
00357 
00358     double dx = (double)(wold - 1) / (wnew - 1);
00359     double x = dx;
00360 
00361     for(; id != idend; ++id, x += dx)
00362     {
00363         if(x >= 1.0)
00364         {
00365             int xx = (int)x;
00366             i1 += xx;
00367             x -= (double)xx;
00368         }
00369         double x1 = 1.0 - x;
00370 
00371         ad.set(DestTraits::fromRealPromote(x1 * as(i1) + x * as(i1, 1)), id);
00372     }
00373 }
00374 
00375 /********************************************************/
00376 /*                                                      */
00377 /*           resizeImageLinearInterpolation             */
00378 /*                                                      */
00379 /********************************************************/
00380 
00381 /** \brief Resize image using linear interpolation.
00382 
00383     The function uses the standard separable bilinear interpolation algorithm to
00384     obtain a good compromize between quality and speed.
00385 
00386     The range must of both the input and output images (resp. regions)
00387     must be given. Both images must have a size of at
00388     least 2x2. The scaling factors are then calculated
00389     accordingly. If the source image is larger than the destination, it
00390     is smoothed (band limited) using a recursive
00391     exponential filter. The source value_type (SrcAccessor::value_type) must
00392     be a linear space, i.e. it must support addition, multiplication
00393     with a scalar real number and \ref NumericTraits "NumericTraits".
00394     The function uses accessors.
00395 
00396     <b> Declarations:</b>
00397 
00398     pass arguments explicitly:
00399     \code
00400     namespace vigra {
00401         template <class SrcImageIterator, class SrcAccessor,
00402                   class DestImageIterator, class DestAccessor>
00403         void
00404         resizeImageLinearInterpolation(
00405               SrcImageIterator is, SrcImageIterator iend, SrcAccessor sa,
00406               DestImageIterator id, DestImageIterator idend, DestAccessor da)
00407     }
00408     \endcode
00409 
00410 
00411     use argument objects in conjunction with \ref ArgumentObjectFactories :
00412     \code
00413     namespace vigra {
00414         template <class SrcImageIterator, class SrcAccessor,
00415                   class DestImageIterator, class DestAccessor>
00416         void
00417         resizeImageLinearInterpolation(
00418               triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src,
00419               triple<DestImageIterator, DestImageIterator, DestAccessor> dest)
00420     }
00421     \endcode
00422 
00423     <b> Usage:</b>
00424 
00425         <b>\#include</b> <<a href="resizeimage_8hxx-source.html">vigra/resizeimage.hxx</a>><br>
00426         Namespace: vigra
00427 
00428     \code
00429     vigra::resizeImageLinearInterpolation(
00430                src.upperLeft(), src.lowerRight(), src.accessor(),
00431                dest.upperLeft(), dest.lowerRight(), dest.accessor());
00432 
00433     \endcode
00434 
00435     <b> Required Interface:</b>
00436 
00437     \code
00438     SrcImageIterator src_upperleft, src_lowerright;
00439     DestImageIterator dest_upperleft, src_lowerright;
00440 
00441     SrcAccessor src_accessor;
00442     DestAccessor dest_accessor;
00443 
00444     NumericTraits<SrcAccessor::value_type>::RealPromote
00445                              u = src_accessor(src_upperleft),
00446                  v = src_accessor(src_upperleft, 1);
00447     double d;
00448 
00449     u = d * v;
00450     u = u + v;
00451 
00452     dest_accessor.set(
00453         NumericTraits<DestAccessor::value_type>::fromRealPromote(u),
00454     dest_upperleft);
00455 
00456     \endcode
00457 
00458     <b> Preconditions:</b>
00459 
00460     \code
00461     src_lowerright.x - src_upperleft.x > 1
00462     src_lowerright.y - src_upperleft.y > 1
00463     dest_lowerright.x - dest_upperleft.x > 1
00464     dest_lowerright.y - dest_upperleft.y > 1
00465     \endcode
00466 
00467 */
00468 doxygen_overloaded_function(template <...> void resizeImageLinearInterpolation)
00469 
00470 template <class SrcIterator, class SrcAccessor,
00471           class DestIterator, class DestAccessor>
00472 void
00473 resizeImageLinearInterpolation(SrcIterator is, SrcIterator iend, SrcAccessor sa,
00474                       DestIterator id, DestIterator idend, DestAccessor da)
00475 {
00476     int w = iend.x - is.x;
00477     int h = iend.y - is.y;
00478 
00479     int wnew = idend.x - id.x;
00480     int hnew = idend.y - id.y;
00481 
00482     vigra_precondition((w > 1) && (h > 1),
00483                  "resizeImageLinearInterpolation(): "
00484                  "Source image to small.\n");
00485     vigra_precondition((wnew > 1) && (hnew > 1),
00486                  "resizeImageLinearInterpolation(): "
00487                  "Destination image to small.\n");
00488 
00489     double const scale = 2.0;
00490 
00491     typedef typename SrcAccessor::value_type SRCVT;
00492     typedef typename NumericTraits<SRCVT>::RealPromote TMPTYPE;
00493     typedef BasicImage<TMPTYPE> TmpImage;
00494     typedef typename TmpImage::traverser TmpImageIterator;
00495 
00496     BasicImage<TMPTYPE> tmp(w, hnew);
00497     BasicImage<TMPTYPE> line((h > w) ? h : w, 1);
00498 
00499     int x,y;
00500 
00501     typename BasicImage<TMPTYPE>::Iterator yt = tmp.upperLeft();
00502     typename TmpImageIterator::row_iterator lt = line.upperLeft().rowIterator();
00503 
00504     for(x=0; x<w; ++x, ++is.x, ++yt.x)
00505     {
00506         typename SrcIterator::column_iterator c1 = is.columnIterator();
00507         typename TmpImageIterator::column_iterator ct = yt.columnIterator();
00508 
00509         if(hnew < h)
00510         {
00511             recursiveSmoothLine(c1, c1 + h, sa,
00512                  lt, line.accessor(), (double)h/hnew/scale);
00513 
00514             resizeLineLinearInterpolation(lt, lt + h, line.accessor(),
00515                                           ct, ct + hnew, tmp.accessor());
00516         }
00517         else
00518         {
00519             resizeLineLinearInterpolation(c1, c1 + h, sa,
00520                                           ct, ct + hnew, tmp.accessor());
00521         }
00522     }
00523 
00524     yt = tmp.upperLeft();
00525 
00526     for(y=0; y < hnew; ++y, ++yt.y, ++id.y)
00527     {
00528         typename DestIterator::row_iterator rd = id.rowIterator();
00529         typename TmpImageIterator::row_iterator rt = yt.rowIterator();
00530 
00531         if(wnew < w)
00532         {
00533             recursiveSmoothLine(rt, rt + w, tmp.accessor(),
00534                               lt, line.accessor(), (double)w/wnew/scale);
00535 
00536             resizeLineLinearInterpolation(lt, lt + w, line.accessor(),
00537                                           rd, rd + wnew, da);
00538         }
00539         else
00540         {
00541             resizeLineLinearInterpolation(rt, rt + w, tmp.accessor(),
00542                                           rd, rd + wnew, da);
00543         }
00544     }
00545 }
00546 
00547 template <class SrcIterator, class SrcAccessor,
00548           class DestIterator, class DestAccessor>
00549 inline
00550 void
00551 resizeImageLinearInterpolation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
00552                                triple<DestIterator, DestIterator, DestAccessor> dest)
00553 {
00554     resizeImageLinearInterpolation(src.first, src.second, src.third,
00555                                    dest.first, dest.second, dest.third);
00556 }
00557 
00558 /***************************************************************/
00559 /*                                                             */
00560 /*                resizeImageSplineInterpolation               */
00561 /*                                                             */
00562 /***************************************************************/
00563 
00564 /** \brief Resize image using B-spline interpolation.
00565 
00566     The function implements separable spline interpolation algorithm described in
00567 
00568     M. Unser, A. Aldroubi, M. Eden, <i>"B-Spline Signal Processing"</i>
00569     IEEE Transactions on Signal Processing, vol. 41, no. 2, pp. 821-833 (part I),
00570     pp. 834-848 (part II), 1993.
00571 
00572     to obtain optimal interpolation quality and speed. You may pass the funcion
00573     a spline of arbitrary order (e.g. <TT>BSpline<ORDER, double></tt> or
00574     <TT>CatmullRomSpline<double></tt>). The default is a third order spline
00575     which gives a twice continuously differentiable interpolant.
00576     The implementation ensures that image values are interpolated rather
00577     than smoothed by first calling a recursive (sharpening) prefilter as
00578     described in the above paper. Then the actual interpolation is done
00579     using \ref resamplingConvolveLine().
00580 
00581     The range of both the input and output images (resp. regions)
00582     must be given. The input image must have a size of at
00583     least 4x4, the destination of at least 2x2. The scaling factors are then calculated
00584     accordingly. If the source image is larger than the destination, it
00585     is smoothed (band limited) using a recursive
00586     exponential filter. The source value_type (SrcAccessor::value_type) must
00587     be a linear algebra, i.e. it must support addition, subtraction,
00588     and multiplication (+, -, *), multiplication with a scalar
00589     real number and \ref NumericTraits "NumericTraits".
00590     The function uses accessors.
00591 
00592     <b> Declarations:</b>
00593 
00594     pass arguments explicitly:
00595     \code
00596     namespace vigra {
00597         template <class SrcImageIterator, class SrcAccessor,
00598                   class DestImageIterator, class DestAccessor,
00599                   class SPLINE>
00600         void
00601         resizeImageSplineInterpolation(
00602               SrcImageIterator is, SrcImageIterator iend, SrcAccessor sa,
00603               DestImageIterator id, DestImageIterator idend, DestAccessor da,
00604               SPLINE spline = BSpline<3, double>())
00605     }
00606     \endcode
00607 
00608 
00609     use argument objects in conjunction with \ref ArgumentObjectFactories :
00610     \code
00611     namespace vigra {
00612         template <class SrcImageIterator, class SrcAccessor,
00613                   class DestImageIterator, class DestAccessor,
00614                   class SPLINE>
00615         void
00616         resizeImageSplineInterpolation(
00617               triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src,
00618               triple<DestImageIterator, DestImageIterator, DestAccessor> dest,
00619               SPLINE spline = BSpline<3, double>())
00620     }
00621     \endcode
00622 
00623     <b> Usage:</b>
00624 
00625         <b>\#include</b> <<a href="resizeimage_8hxx-source.html">vigra/resizeimage.hxx</a>><br>
00626         Namespace: vigra
00627 
00628     \code
00629     vigra::resizeImageSplineInterpolation(
00630                src.upperLeft(), src.lowerRight(), src.accessor(),
00631                dest.upperLeft(), dest.lowerRight(), dest.accessor());
00632 
00633     \endcode
00634 
00635     <b> Required Interface:</b>
00636 
00637     \code
00638     SrcImageIterator src_upperleft, src_lowerright;
00639     DestImageIterator dest_upperleft, src_lowerright;
00640 
00641     SrcAccessor src_accessor;
00642     DestAccessor dest_accessor;
00643 
00644     NumericTraits<SrcAccessor::value_type>::RealPromote
00645                              u = src_accessor(src_upperleft),
00646                  v = src_accessor(src_upperleft, 1);
00647     double d;
00648 
00649     u = d * v;
00650     u = u + v;
00651     u = u - v;
00652     u = u * v;
00653     u += v;
00654     u -= v;
00655 
00656     dest_accessor.set(
00657         NumericTraits<DestAccessor::value_type>::fromRealPromote(u),
00658     dest_upperleft);
00659 
00660     \endcode
00661 
00662     <b> Preconditions:</b>
00663 
00664     \code
00665     src_lowerright.x - src_upperleft.x > 3
00666     src_lowerright.y - src_upperleft.y > 3
00667     dest_lowerright.x - dest_upperleft.x > 1
00668     dest_lowerright.y - dest_upperleft.y > 1
00669     \endcode
00670 
00671 */
00672 doxygen_overloaded_function(template <...> void resizeImageSplineInterpolation)
00673 
00674 template <class SrcIterator, class SrcAccessor,
00675           class DestIterator, class DestAccessor,
00676           class SPLINE>
00677 void
00678 resizeImageSplineInterpolation(
00679     SrcIterator src_iter, SrcIterator src_iter_end, SrcAccessor src_acc,
00680     DestIterator dest_iter, DestIterator dest_iter_end, DestAccessor dest_acc,
00681     SPLINE const & spline)
00682 {
00683 
00684     int width_old = src_iter_end.x - src_iter.x;
00685     int height_old = src_iter_end.y - src_iter.y;
00686 
00687     int width_new = dest_iter_end.x - dest_iter.x;
00688     int height_new = dest_iter_end.y - dest_iter.y;
00689 
00690     vigra_precondition((width_old > 1) && (height_old > 1),
00691                  "resizeImageSplineInterpolation(): "
00692                  "Source image to small.\n");
00693 
00694     vigra_precondition((width_new > 1) && (height_new > 1),
00695                  "resizeImageSplineInterpolation(): "
00696                  "Destination image to small.\n");
00697 
00698     Rational<int> xratio(width_new - 1, width_old - 1);
00699     Rational<int> yratio(height_new - 1, height_old - 1);
00700     Rational<int> offset(0);
00701     resampling_detail::MapTargetToSourceCoordinate xmapCoordinate(xratio, offset);
00702     resampling_detail::MapTargetToSourceCoordinate ymapCoordinate(yratio, offset);
00703     int xperiod = lcm(xratio.numerator(), xratio.denominator());
00704     int yperiod = lcm(yratio.numerator(), yratio.denominator());
00705 
00706     double const scale = 2.0;
00707 
00708     typedef typename SrcAccessor::value_type SRCVT;
00709     typedef typename NumericTraits<SRCVT>::RealPromote TMPTYPE;
00710     typedef BasicImage<TMPTYPE> TmpImage;
00711     typedef typename TmpImage::traverser TmpImageIterator;
00712 
00713     BasicImage<TMPTYPE> tmp(width_old, height_new);
00714 
00715     BasicImage<TMPTYPE> line((height_old > width_old) ? height_old : width_old, 1);
00716     typename BasicImage<TMPTYPE>::Accessor tmp_acc = tmp.accessor();
00717     ArrayVector<double> const & prefilterCoeffs = spline.prefilterCoefficients();
00718 
00719     int x,y;
00720 
00721     ArrayVector<Kernel1D<double> > kernels(yperiod);
00722     createResamplingKernels(spline, ymapCoordinate, kernels);
00723 
00724     typename BasicImage<TMPTYPE>::Iterator y_tmp = tmp.upperLeft();
00725     typename TmpImageIterator::row_iterator line_tmp = line.upperLeft().rowIterator();
00726 
00727     for(x=0; x<width_old; ++x, ++src_iter.x, ++y_tmp.x)
00728     {
00729 
00730         typename SrcIterator::column_iterator c_src = src_iter.columnIterator();
00731         typename TmpImageIterator::column_iterator c_tmp = y_tmp.columnIterator();
00732 
00733         if(prefilterCoeffs.size() == 0)
00734         {
00735             if(height_new >= height_old)
00736             {
00737                 resamplingConvolveLine(c_src, c_src + height_old, src_acc,
00738                                        c_tmp, c_tmp + height_new, tmp_acc,
00739                                        kernels, ymapCoordinate);
00740             }
00741             else
00742             {
00743                 recursiveSmoothLine(c_src, c_src + height_old, src_acc,
00744                      line_tmp, line.accessor(), (double)height_old/height_new/scale);
00745                 resamplingConvolveLine(line_tmp, line_tmp + height_old, line.accessor(),
00746                                        c_tmp, c_tmp + height_new, tmp_acc,
00747                                        kernels, ymapCoordinate);
00748             }
00749         }
00750         else
00751         {
00752             recursiveFilterLine(c_src, c_src + height_old, src_acc,
00753                                 line_tmp, line.accessor(),
00754                                 prefilterCoeffs[0], BORDER_TREATMENT_REFLECT);
00755             for(unsigned int b = 1; b < prefilterCoeffs.size(); ++b)
00756             {
00757                 recursiveFilterLine(line_tmp, line_tmp + height_old, line.accessor(),
00758                                     line_tmp, line.accessor(),
00759                                     prefilterCoeffs[b], BORDER_TREATMENT_REFLECT);
00760             }
00761             if(height_new < height_old)
00762             {
00763                 recursiveSmoothLine(line_tmp, line_tmp + height_old, line.accessor(),
00764                      line_tmp, line.accessor(), (double)height_old/height_new/scale);
00765             }
00766             resamplingConvolveLine(line_tmp, line_tmp + height_old, line.accessor(),
00767                                    c_tmp, c_tmp + height_new, tmp_acc,
00768                                    kernels, ymapCoordinate);
00769         }
00770     }
00771 
00772     y_tmp = tmp.upperLeft();
00773 
00774     DestIterator dest = dest_iter;
00775 
00776     kernels.resize(xperiod);
00777     createResamplingKernels(spline, xmapCoordinate, kernels);
00778 
00779     for(y=0; y < height_new; ++y, ++y_tmp.y, ++dest_iter.y)
00780     {
00781         typename DestIterator::row_iterator r_dest = dest_iter.rowIterator();
00782         typename TmpImageIterator::row_iterator r_tmp = y_tmp.rowIterator();
00783 
00784         if(prefilterCoeffs.size() == 0)
00785         {
00786             if(width_new >= width_old)
00787             {
00788                 resamplingConvolveLine(r_tmp, r_tmp + width_old, tmp.accessor(),
00789                                        r_dest, r_dest + width_new, dest_acc,
00790                                        kernels, xmapCoordinate);
00791             }
00792             else
00793             {
00794                 recursiveSmoothLine(r_tmp, r_tmp + width_old, tmp.accessor(),
00795                                   line_tmp, line.accessor(), (double)width_old/width_new/scale);
00796                 resamplingConvolveLine(line_tmp, line_tmp + width_old, line.accessor(),
00797                                        r_dest, r_dest + width_new, dest_acc,
00798                                        kernels, xmapCoordinate);
00799             }
00800         }
00801         else
00802         {
00803             recursiveFilterLine(r_tmp, r_tmp + width_old, tmp.accessor(),
00804                                 line_tmp, line.accessor(),
00805                                 prefilterCoeffs[0], BORDER_TREATMENT_REFLECT);
00806             for(unsigned int b = 1; b < prefilterCoeffs.size(); ++b)
00807             {
00808                 recursiveFilterLine(line_tmp, line_tmp + width_old, line.accessor(),
00809                                     line_tmp, line.accessor(),
00810                                     prefilterCoeffs[b], BORDER_TREATMENT_REFLECT);
00811             }
00812             if(width_new < width_old)
00813             {
00814                 recursiveSmoothLine(line_tmp, line_tmp + width_old, line.accessor(),
00815                                     line_tmp, line.accessor(), (double)width_old/width_new/scale);
00816             }
00817             resamplingConvolveLine(line_tmp, line_tmp + width_old, line.accessor(),
00818                                    r_dest, r_dest + width_new, dest_acc,
00819                                    kernels, xmapCoordinate);
00820         }
00821     }
00822 }
00823 
00824 template <class SrcIterator, class SrcAccessor,
00825           class DestIterator, class DestAccessor,
00826           class SPLINE>
00827 inline
00828 void
00829 resizeImageSplineInterpolation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
00830                       triple<DestIterator, DestIterator, DestAccessor> dest,
00831                       SPLINE const & spline)
00832 {
00833     resizeImageSplineInterpolation(src.first, src.second, src.third,
00834                                    dest.first, dest.second, dest.third, spline);
00835 }
00836 
00837 template <class SrcIterator, class SrcAccessor,
00838           class DestIterator, class DestAccessor>
00839 void
00840 resizeImageSplineInterpolation(SrcIterator is, SrcIterator iend, SrcAccessor sa,
00841                       DestIterator id, DestIterator idend, DestAccessor da)
00842 {
00843     resizeImageSplineInterpolation(is, iend, sa, id, idend, da, BSpline<3, double>());
00844 }
00845 
00846 template <class SrcIterator, class SrcAccessor,
00847           class DestIterator, class DestAccessor>
00848 inline
00849 void
00850 resizeImageSplineInterpolation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
00851                       triple<DestIterator, DestIterator, DestAccessor> dest)
00852 {
00853     resizeImageSplineInterpolation(src.first, src.second, src.third,
00854                                    dest.first, dest.second, dest.third);
00855 }
00856 
00857 /*****************************************************************/
00858 /*                                                               */
00859 /*              resizeImageCatmullRomInterpolation               */
00860 /*                                                               */
00861 /*****************************************************************/
00862 
00863 /** \brief Resize image using the Catmull/Rom interpolation function.
00864 
00865     The function calls like \ref resizeImageSplineInterpolation() with
00866     \ref vigra::CatmullRomSpline as an interpolation kernel.
00867     The interpolated function has one continuous derivative.
00868     (See \ref resizeImageSplineInterpolation() for more documentation)
00869 
00870     <b> Declarations:</b>
00871 
00872     pass arguments explicitly:
00873     \code
00874     namespace vigra {
00875         template <class SrcIterator, class SrcAccessor,
00876                   class DestIterator, class DestAccessor>
00877         void
00878         resizeImageCatmullRomInterpolation(SrcIterator src_iter, SrcIterator src_iter_end, SrcAccessor src_acc,
00879                               DestIterator dest_iter, DestIterator dest_iter_end, DestAccessor dest_acc);
00880     }
00881     \endcode
00882 
00883 
00884     use argument objects in conjunction with \ref ArgumentObjectFactories :
00885     \code
00886     namespace vigra {
00887         template <class SrcIterator, class SrcAccessor,
00888                   class DestIterator, class DestAccessor>
00889         void
00890         resizeImageCatmullRomInterpolation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
00891                               triple<DestIterator, DestIterator, DestAccessor> dest);
00892     }
00893     \endcode
00894 
00895 
00896     <b>\#include</b> <<a href="resizeimage_8hxx-source.html">vigra/resizeimage.hxx</a>><br>
00897     Namespace: vigra
00898 
00899 */
00900 doxygen_overloaded_function(template <...> void resizeImageCatmullRomInterpolation)
00901 
00902 template <class SrcIterator, class SrcAccessor,
00903           class DestIterator, class DestAccessor>
00904 inline void
00905 resizeImageCatmullRomInterpolation(SrcIterator src_iter, SrcIterator src_iter_end, SrcAccessor src_acc,
00906                       DestIterator dest_iter, DestIterator dest_iter_end, DestAccessor dest_acc)
00907 {
00908     resizeImageSplineInterpolation(src_iter, src_iter_end, src_acc, dest_iter, dest_iter_end, dest_acc,
00909                                   CatmullRomSpline<double>());
00910 }
00911 
00912 template <class SrcIterator, class SrcAccessor,
00913           class DestIterator, class DestAccessor>
00914 inline
00915 void
00916 resizeImageCatmullRomInterpolation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
00917                       triple<DestIterator, DestIterator, DestAccessor> dest)
00918 {
00919     resizeImageCatmullRomInterpolation(src.first, src.second, src.third,
00920                                      dest.first, dest.second, dest.third);
00921 }
00922 
00923 #if 0
00924 /*****************************************************************/
00925 /*                                                               */
00926 /*              resizeImageCubicInterpolation                    */
00927 /*                                                               */
00928 /*****************************************************************/
00929 
00930 /** \brief Resize image using the cardinal B-spline interpolation function.
00931 
00932     The function calls like \ref resizeImageSplineInterpolation() with
00933     \ref vigra::BSpline<3, double> and prefiltering as an interpolation kernel.
00934     The interpolated function has two continuous derivatives.
00935     (See \ref resizeImageSplineInterpolation() for more documentation)
00936 
00937     <b>\#include</b> <<a href="resizeimage_8hxx-source.html">vigra/resizeimage.hxx</a>><br>
00938     Namespace: vigra
00939 
00940 */
00941 template <class SrcIterator, class SrcAccessor,
00942           class DestIterator, class DestAccessor>
00943 void
00944 resizeImageCubicInterpolation(SrcIterator src_iter, SrcIterator src_iter_end, SrcAccessor src_acc,
00945                       DestIterator dest_iter, DestIterator dest_iter_end, DestAccessor dest_acc)
00946 {
00947     resizeImageSplineInterpolation(src_iter, src_iter_end, src_acc, dest_iter, dest_iter_end, dest_acc,
00948                                   BSpline<3, double>());
00949 }
00950 
00951 template <class SrcIterator, class SrcAccessor,
00952           class DestIterator, class DestAccessor>
00953 inline
00954 void
00955 resizeImageCubicInterpolation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
00956                       triple<DestIterator, DestIterator, DestAccessor> dest)
00957 {
00958     resizeImageCubicInterpolation(src.first, src.second, src.third,
00959                                    dest.first, dest.second, dest.third);
00960 }
00961 #endif
00962 
00963 /*****************************************************************/
00964 /*                                                               */
00965 /*              resizeImageCoscotInterpolation                   */
00966 /*                                                               */
00967 /*****************************************************************/
00968 
00969 /** \brief Resize image using the Coscot interpolation function.
00970 
00971     The function calls \ref resizeImageSplineInterpolation() with
00972     \ref vigra::CoscotFunction as an interpolation kernel.
00973     The interpolated function has one continuous derivative.
00974     (See \ref resizeImageSplineInterpolation() for more documentation)
00975 
00976     <b> Declarations:</b>
00977 
00978     pass arguments explicitly:
00979     \code
00980     namespace vigra {
00981         template <class SrcIterator, class SrcAccessor,
00982                   class DestIterator, class DestAccessor>
00983         void
00984         resizeImageCoscotInterpolation(SrcIterator src_iter, SrcIterator src_iter_end, SrcAccessor src_acc,
00985                               DestIterator dest_iter, DestIterator dest_iter_end, DestAccessor dest_acc);
00986     }
00987     \endcode
00988 
00989 
00990     use argument objects in conjunction with \ref ArgumentObjectFactories :
00991     \code
00992     namespace vigra {
00993         template <class SrcIterator, class SrcAccessor,
00994                   class DestIterator, class DestAccessor>
00995         void
00996         resizeImageCoscotInterpolation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
00997                               triple<DestIterator, DestIterator, DestAccessor> dest);
00998     }
00999     \endcode
01000 
01001 
01002     <b>\#include</b> <<a href="resizeimage_8hxx-source.html">vigra/resizeimage.hxx</a>><br>
01003     Namespace: vigra
01004 
01005 */
01006 doxygen_overloaded_function(template <...> void resizeImageCoscotInterpolation)
01007 
01008 template <class SrcIterator, class SrcAccessor,
01009           class DestIterator, class DestAccessor>
01010 void
01011 resizeImageCoscotInterpolation(SrcIterator src_iter, SrcIterator src_iter_end, SrcAccessor src_acc,
01012                       DestIterator dest_iter, DestIterator dest_iter_end, DestAccessor dest_acc)
01013 {
01014     resizeImageSplineInterpolation(src_iter, src_iter_end, src_acc, dest_iter, dest_iter_end, dest_acc,
01015                                    CoscotFunction<double>());
01016 }
01017 
01018 template <class SrcIterator, class SrcAccessor,
01019           class DestIterator, class DestAccessor>
01020 inline
01021 void
01022 resizeImageCoscotInterpolation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
01023                       triple<DestIterator, DestIterator, DestAccessor> dest)
01024 {
01025     resizeImageCoscotInterpolation(src.first, src.second, src.third,
01026                                    dest.first, dest.second, dest.third);
01027 }
01028 
01029 
01030 #if 0 // old version of the spline interpolation algorithm
01031 
01032 /********************************************************/
01033 /*                                                      */
01034 /*           resizeCalculateSplineCoefficients          */
01035 /*         (internally used by resize functions)        */
01036 /*                                                      */
01037 /********************************************************/
01038 
01039 template <class SrcIterator, class SrcAccessor, class VALUETYPE>
01040 void
01041 resizeCalculateSplineCoefficients(SrcIterator i1, SrcIterator iend,
01042                 SrcAccessor a, VALUETYPE * i2)
01043 {
01044     int n = iend - i1;
01045 
01046     if(n <= 0) return;
01047 
01048     VALUETYPE zero = NumericTraits<VALUETYPE>::zero();
01049     VALUETYPE two = 2.0 * NumericTraits<VALUETYPE>::one();
01050     VALUETYPE half = 0.5 * NumericTraits<VALUETYPE>::one();
01051 
01052     *i2 = zero;
01053     if(n == 1) return;
01054 
01055     std::vector<VALUETYPE> vec(n);
01056     typename std::vector<VALUETYPE>::iterator u = vec.begin();
01057 
01058     *u = zero;
01059 
01060     for(++i1, ++i2, ++u, --iend; i1 != iend; ++i1, ++i2, ++u)
01061     {
01062         VALUETYPE p = 0.5 * i2[-1] + two;
01063         *i2 = half / p;
01064         *u = 3.0 *(a(i1,1) - 2.0 * a(i1) + a(i1, -1)) - 0.5 * u[-1] / p;
01065     }
01066 
01067     *i2 = zero;
01068 
01069     for(--i2, --u; u != vec; --u, --i2)
01070     {
01071         *i2 = *i2 * i2[1] + *u;
01072     }
01073 }
01074 
01075 /********************************************************/
01076 /*                                                      */
01077 /*         resizeImageInternalSplineGradient            */
01078 /*                                                      */
01079 /********************************************************/
01080 
01081 template <class SrcIterator, class SrcAccessor,
01082           class DoubleIterator, class TempIterator, class DestIterator>
01083 void
01084 resizeImageInternalSplineGradient(SrcIterator in, SrcIterator inend, SrcAccessor sa,
01085                          DoubleIterator tmp, TempIterator r, DestIterator id)
01086 {
01087     int w = inend - in;
01088 
01089     int x;
01090 
01091     typedef typename SrcAccessor::value_type SRCVT;
01092     typedef typename NumericTraits<SRCVT>::RealPromote TMPTYPE;
01093 
01094     // calculate border derivatives
01095     SrcIterator xs = in;
01096     TMPTYPE p0 = -11.0/6.0 * sa(xs);  ++xs;
01097             p0 += 3.0 * sa(xs);  ++xs;
01098             p0 += -1.5 * sa(xs);  ++xs;
01099             p0 += 1.0/3.0 * sa(xs);
01100 
01101     xs = in + w-1;
01102     TMPTYPE pw = 11.0/6.0 * sa(xs);  --xs;
01103             pw += -3.0 * sa(xs);  --xs;
01104             pw +=  1.5 * sa(xs);  --xs;
01105             pw += -1.0/3.0 * sa(xs);
01106 
01107     xs = in + 2;
01108     SrcIterator xs1 = in;
01109 
01110     for(x=1; x<w-1; ++x, ++xs, ++xs1)
01111     {
01112         r[x] = 3.0 * (sa(xs) - sa(xs1));
01113     }
01114 
01115     r[1] -= p0;
01116     r[w-2] -= pw;
01117 
01118     double q = 0.25;
01119 
01120     id[0] = p0;
01121     id[w-1] = pw;
01122     id[1] = 0.25 * r[1];
01123 
01124     for(x=2; x<w-1; ++x)
01125     {
01126         tmp[x] = q;
01127         q = 1.0 / (4.0 - q);
01128         id[x] = q * (r[x] - id[x-1]);
01129     }
01130 
01131     for(x=w-3; x>=1; --x)
01132     {
01133         id[x] -= tmp[x+1]*id[x+1];
01134     }
01135 }
01136 
01137 /********************************************************/
01138 /*                                                      */
01139 /*         resizeImageInternalSplineInterpolation       */
01140 /*                                                      */
01141 /********************************************************/
01142 
01143 template <class SrcIterator, class SrcAccessor,
01144           class DestIterator, class DestAccessor>
01145 void
01146 resizeImageInternalSplineInterpolation(SrcIterator is, SrcIterator iend, SrcAccessor sa,
01147                       DestIterator id, DestIterator idend, DestAccessor da)
01148 {
01149     int w = iend.x - is.x;
01150     int h = iend.y - is.y;
01151 
01152     int wnew = idend.x - id.x;
01153     int hnew = idend.y - id.y;
01154 
01155     typedef typename SrcAccessor::value_type SRCVT;
01156     typedef typename NumericTraits<SRCVT>::RealPromote TMPTYPE;
01157     typedef typename BasicImage<TMPTYPE>::Iterator TMPITER;
01158     typedef
01159         NumericTraits<typename DestAccessor::value_type> DestTraits;
01160 
01161     BasicImage<TMPTYPE> dx(w,h);
01162     BasicImage<TMPTYPE> dy(w,h);
01163     BasicImage<TMPTYPE> dxy(w,h);
01164     BasicImage<TMPTYPE> W(4,4), W1(4,4);
01165     std::vector<TMPTYPE> R(w > h ? w : h);
01166     std::vector<double> tmp(w > h ? w : h);
01167 
01168     typename BasicImage<TMPTYPE>::Accessor ta;
01169 
01170     SrcIterator in = is;
01171 
01172     TMPITER idx = dx.upperLeft();
01173     TMPITER idy = dy.upperLeft();
01174     TMPITER idxy = dxy.upperLeft();
01175     typename std::vector<TMPTYPE>::iterator r = R.begin();
01176     typename std::vector<double>::iterator it = tmp.begin();
01177 
01178     double ig[] = { 1.0, 0.0, -3.0,  2.0,
01179                     0.0, 1.0, -2.0,  1.0,
01180                     0.0, 0.0,  3.0, -2.0,
01181                     0.0, 0.0, -1.0,  1.0 };
01182 
01183     int x, y, i, j, k;
01184 
01185 
01186     // calculate x derivatives
01187     for(y=0; y<h; ++y, ++in.y, ++idx.y)
01188     {
01189         typename SrcIterator::row_iterator sr = in.rowIterator();
01190         typename TMPITER::row_iterator dr = idx.rowIterator();
01191         resizeImageInternalSplineGradient(sr, sr+w, sa,
01192                                           it, r, dr);
01193     }
01194 
01195     in = is;
01196 
01197     // calculate y derivatives
01198     for(x=0; x<w; ++x, ++in.x, ++idy.x)
01199     {
01200         typename SrcIterator::column_iterator sc = in.columnIterator();
01201         typename TMPITER::column_iterator dc = idy.columnIterator();
01202         resizeImageInternalSplineGradient(sc, sc+h, sa,
01203                                           it, r, dc);
01204     }
01205 
01206     in = is;
01207     idy = dy.upperLeft();
01208 
01209     // calculate mixed derivatives
01210     for(y=0; y<h; ++y, ++idy.y, ++idxy.y)
01211     {
01212         typename TMPITER::row_iterator sr = idy.rowIterator();
01213         typename TMPITER::row_iterator dr = idxy.rowIterator();
01214         resizeImageInternalSplineGradient(sr, sr+w, ta,
01215                                           it, r, dr);
01216     }
01217 
01218     double du = (double)(w-1) / (wnew-1);
01219     double dv = (double)(h-1) / (hnew-1);
01220     double ov = 0.0;
01221     int oy = 0;
01222     int yy = oy;
01223 
01224     DestIterator xxd = id, yyd = id;
01225 
01226     static Diff2D down(0,1), right(1,0), downright(1,1);
01227 
01228     for(y=0; y<h-1; ++y, ++in.y, ov -= 1.0)
01229     {
01230         if(y < h-2 && ov >= 1.0) continue;
01231         int y1 = y+1;
01232         double v = ov;
01233         double ou = 0.0;
01234         int ox = 0;
01235         int xx = ox;
01236 
01237         SrcIterator xs = in;
01238         for(x=0; x<w-1; ++x, ++xs.x, ou -= 1.0)
01239         {
01240             if(x < w-2 && ou >= 1.0) continue;
01241             int x1 = x+1;
01242             double u = ou;
01243 
01244             DestIterator xd = id + Diff2D(ox,oy);
01245             W[0][0] = sa(xs);
01246             W[0][1] = dy(x, y);
01247             W[0][2] = sa(xs, down);
01248             W[0][3] = dy(x, y1);
01249             W[1][0] = dx(x, y);
01250             W[1][1] = dxy(x, y);
01251             W[1][2] = dx(x, y1);
01252             W[1][3] = dxy(x, y1);
01253             W[2][0] = sa(xs, right);
01254             W[2][1] = dy(x1,y);
01255             W[2][2] = sa(xs, downright);
01256             W[2][3] = dy(x1, y1);
01257             W[3][0] = dx(x1, y);
01258             W[3][1] = dxy(x1, y);
01259             W[3][2] = dx(x1, y1);
01260             W[3][3] = dxy(x1, y1);
01261 
01262             for(i=0; i<4; ++i)
01263             {
01264                 for(j=0; j<4; ++j)
01265                 {
01266                     W1[j][i] = ig[j] * W[0][i];
01267                     for(k=1; k<4; ++k)
01268                     {
01269                         W1[j][i] += ig[j+4*k] * W[k][i];
01270                     }
01271                 }
01272             }
01273             for(i=0; i<4; ++i)
01274             {
01275                 for(j=0; j<4; ++j)
01276                 {
01277                     W[j][i] = ig[i] * W1[j][0];
01278                     for(k=1; k<4; ++k)
01279                     {
01280                        W[j][i] += ig[4*k+i] * W1[j][k];
01281                     }
01282                 }
01283             }
01284 
01285             TMPTYPE a1,a2,a3,a4;
01286 
01287             yyd = xd;
01288             for(v=ov, yy=oy; v<1.0; v+=dv, ++yyd.y, ++yy)
01289             {
01290                 a1 = W[0][0] + v * (W[0][1] +
01291                                v * (W[0][2] + v * W[0][3]));
01292                 a2 = W[1][0] + v * (W[1][1] +
01293                                v * (W[1][2] + v * W[1][3]));
01294                 a3 = W[2][0] + v * (W[2][1] +
01295                                v * (W[2][2] + v * W[2][3]));
01296                 a4 = W[3][0] + v * (W[3][1] +
01297                                v * (W[3][2] + v * W[3][3]));
01298 
01299                 xxd = yyd;
01300                 for(u=ou, xx=ox; u<1.0; u+=du, ++xxd.x, ++xx)
01301                 {
01302                     da.set(DestTraits::fromRealPromote(a1 + u * (a2 + u * (a3 + u * a4))), xxd);
01303                 }
01304 
01305                 if(xx == wnew-1)
01306                 {
01307                     da.set(DestTraits::fromRealPromote(a1 + a2 + a3 + a4), xxd);
01308                 }
01309             }
01310 
01311             if(yy == hnew-1)
01312             {
01313                 a1 = W[0][0] + W[0][1] + W[0][2] + W[0][3];
01314                 a2 = W[1][0] + W[1][1] + W[1][2] + W[1][3];
01315                 a3 = W[2][0] + W[2][1] + W[2][2] + W[2][3];
01316                 a4 = W[3][0] + W[3][1] + W[3][2] + W[3][3];
01317 
01318                 DestIterator xxd = yyd;
01319                 for(u=ou, xx=ox; u<1.0; u+=du, ++xxd.x, ++xx)
01320                 {
01321                     da.set(DestTraits::fromRealPromote(a1 + u * (a2 + u * (a3 + u * a4))), xxd);
01322                 }
01323 
01324                 if(xx == wnew-1)
01325                 {
01326                     da.set(DestTraits::fromRealPromote(a1 + a2 + a3 + a4), xxd);
01327                 }
01328             }
01329 
01330             ou = u;
01331             ox = xx;
01332         }
01333         ov = v;
01334         oy = yy;
01335     }
01336 }
01337 
01338 template <class SrcIterator, class SrcAccessor,
01339           class DestIterator, class DestAccessor>
01340 void
01341 resizeImageSplineInterpolation(SrcIterator is, SrcIterator iend, SrcAccessor sa,
01342                       DestIterator id, DestIterator idend, DestAccessor da)
01343 {
01344     int w = iend.x - is.x;
01345     int h = iend.y - is.y;
01346 
01347     int wnew = idend.x - id.x;
01348     int hnew = idend.y - id.y;
01349 
01350     vigra_precondition((w > 3) && (h > 3),
01351                  "resizeImageSplineInterpolation(): "
01352                  "Source image to small.\n");
01353     vigra_precondition((wnew > 1) && (hnew > 1),
01354                  "resizeImageSplineInterpolation(): "
01355                  "Destination image to small.\n");
01356 
01357     double scale = 2.0;
01358 
01359     if(wnew < w || hnew < h)
01360     {
01361         typedef typename SrcAccessor::value_type SRCVT;
01362         typedef typename NumericTraits<SRCVT>::RealPromote TMPTYPE;
01363         typedef typename BasicImage<TMPTYPE>::Iterator TMPITER;
01364 
01365         BasicImage<TMPTYPE> t(w,h);
01366         TMPITER it = t.upperLeft();
01367 
01368         if(wnew < w)
01369         {
01370             recursiveSmoothX(is, iend, sa,
01371                     it, t.accessor(), (double)w/wnew/scale);
01372 
01373             if(hnew < h)
01374             {
01375                recursiveSmoothY(it, t.lowerRight(), t.accessor(),
01376                     it, t.accessor(), (double)h/hnew/scale);
01377             }
01378         }
01379         else
01380         {
01381            recursiveSmoothY(is, iend, sa,
01382                     it, t.accessor(), (double)h/hnew/scale);
01383         }
01384 
01385         resizeImageInternalSplineInterpolation(it, t.lowerRight(), t.accessor(),
01386                                                id, idend, da);
01387     }
01388     else
01389     {
01390         resizeImageInternalSplineInterpolation(is, iend, sa, id, idend, da);
01391     }
01392 }
01393 
01394 template <class SrcIterator, class SrcAccessor,
01395           class DestIterator, class DestAccessor>
01396 inline
01397 void
01398 resizeImageSplineInterpolation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
01399                       triple<DestIterator, DestIterator, DestAccessor> dest)
01400 {
01401     resizeImageSplineInterpolation(src.first, src.second, src.third,
01402                                    dest.first, dest.second, dest.third);
01403 }
01404 #endif  // old alghorithm version
01405 
01406 //@}
01407 
01408 } // namespace vigra
01409 
01410 #endif // VIGRA_RESIZEIMAGE_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.6.0 (5 Nov 2009)