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

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