[ 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     kernels.resize(xperiod);
00766     createResamplingKernels(spline, xmapCoordinate, kernels);
00767 
00768     for(y=0; y < height_new; ++y, ++y_tmp.y, ++dest_iter.y)
00769     {
00770         typename DestIterator::row_iterator r_dest = dest_iter.rowIterator();
00771         typename TmpImageIterator::row_iterator r_tmp = y_tmp.rowIterator();
00772 
00773         if(prefilterCoeffs.size() == 0)
00774         {
00775             if(width_new >= width_old)
00776             {
00777                 resamplingConvolveLine(r_tmp, r_tmp + width_old, tmp.accessor(),
00778                                        r_dest, r_dest + width_new, dest_acc,
00779                                        kernels, xmapCoordinate);
00780             }
00781             else
00782             {
00783                 recursiveSmoothLine(r_tmp, r_tmp + width_old, tmp.accessor(),
00784                                   line_tmp, line.accessor(), (double)width_old/width_new/scale);
00785                 resamplingConvolveLine(line_tmp, line_tmp + width_old, line.accessor(),
00786                                        r_dest, r_dest + width_new, dest_acc,
00787                                        kernels, xmapCoordinate);
00788             }
00789         }
00790         else
00791         {
00792             recursiveFilterLine(r_tmp, r_tmp + width_old, tmp.accessor(),
00793                                 line_tmp, line.accessor(),
00794                                 prefilterCoeffs[0], BORDER_TREATMENT_REFLECT);
00795             for(unsigned int b = 1; b < prefilterCoeffs.size(); ++b)
00796             {
00797                 recursiveFilterLine(line_tmp, line_tmp + width_old, line.accessor(),
00798                                     line_tmp, line.accessor(),
00799                                     prefilterCoeffs[b], BORDER_TREATMENT_REFLECT);
00800             }
00801             if(width_new < width_old)
00802             {
00803                 recursiveSmoothLine(line_tmp, line_tmp + width_old, line.accessor(),
00804                                     line_tmp, line.accessor(), (double)width_old/width_new/scale);
00805             }
00806             resamplingConvolveLine(line_tmp, line_tmp + width_old, line.accessor(),
00807                                    r_dest, r_dest + width_new, dest_acc,
00808                                    kernels, xmapCoordinate);
00809         }
00810     }
00811 }
00812 
00813 template <class SrcIterator, class SrcAccessor,
00814           class DestIterator, class DestAccessor,
00815           class SPLINE>
00816 inline
00817 void
00818 resizeImageSplineInterpolation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
00819                       triple<DestIterator, DestIterator, DestAccessor> dest,
00820                       SPLINE const & spline)
00821 {
00822     resizeImageSplineInterpolation(src.first, src.second, src.third,
00823                                    dest.first, dest.second, dest.third, spline);
00824 }
00825 
00826 template <class SrcIterator, class SrcAccessor,
00827           class DestIterator, class DestAccessor>
00828 void
00829 resizeImageSplineInterpolation(SrcIterator is, SrcIterator iend, SrcAccessor sa,
00830                       DestIterator id, DestIterator idend, DestAccessor da)
00831 {
00832     resizeImageSplineInterpolation(is, iend, sa, id, idend, da, BSpline<3, double>());
00833 }
00834 
00835 template <class SrcIterator, class SrcAccessor,
00836           class DestIterator, class DestAccessor>
00837 inline
00838 void
00839 resizeImageSplineInterpolation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
00840                       triple<DestIterator, DestIterator, DestAccessor> dest)
00841 {
00842     resizeImageSplineInterpolation(src.first, src.second, src.third,
00843                                    dest.first, dest.second, dest.third);
00844 }
00845 
00846 /*****************************************************************/
00847 /*                                                               */
00848 /*              resizeImageCatmullRomInterpolation               */
00849 /*                                                               */
00850 /*****************************************************************/
00851 
00852 /** \brief Resize image using the Catmull/Rom interpolation function.
00853 
00854     The function calls like \ref resizeImageSplineInterpolation() with
00855     \ref vigra::CatmullRomSpline as an interpolation kernel.
00856     The interpolated function has one continuous derivative.
00857     (See \ref resizeImageSplineInterpolation() for more documentation)
00858 
00859     <b> Declarations:</b>
00860 
00861     pass arguments explicitly:
00862     \code
00863     namespace vigra {
00864         template <class SrcIterator, class SrcAccessor,
00865                   class DestIterator, class DestAccessor>
00866         void
00867         resizeImageCatmullRomInterpolation(SrcIterator src_iter, SrcIterator src_iter_end, SrcAccessor src_acc,
00868                               DestIterator dest_iter, DestIterator dest_iter_end, DestAccessor dest_acc);
00869     }
00870     \endcode
00871 
00872 
00873     use argument objects in conjunction with \ref ArgumentObjectFactories :
00874     \code
00875     namespace vigra {
00876         template <class SrcIterator, class SrcAccessor,
00877                   class DestIterator, class DestAccessor>
00878         void
00879         resizeImageCatmullRomInterpolation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
00880                               triple<DestIterator, DestIterator, DestAccessor> dest);
00881     }
00882     \endcode
00883 
00884 
00885     <b>\#include</b> <vigra/resizeimage.hxx><br>
00886     Namespace: vigra
00887 
00888 */
00889 doxygen_overloaded_function(template <...> void resizeImageCatmullRomInterpolation)
00890 
00891 template <class SrcIterator, class SrcAccessor,
00892           class DestIterator, class DestAccessor>
00893 inline void
00894 resizeImageCatmullRomInterpolation(SrcIterator src_iter, SrcIterator src_iter_end, SrcAccessor src_acc,
00895                       DestIterator dest_iter, DestIterator dest_iter_end, DestAccessor dest_acc)
00896 {
00897     resizeImageSplineInterpolation(src_iter, src_iter_end, src_acc, dest_iter, dest_iter_end, dest_acc,
00898                                   CatmullRomSpline<double>());
00899 }
00900 
00901 template <class SrcIterator, class SrcAccessor,
00902           class DestIterator, class DestAccessor>
00903 inline
00904 void
00905 resizeImageCatmullRomInterpolation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
00906                       triple<DestIterator, DestIterator, DestAccessor> dest)
00907 {
00908     resizeImageCatmullRomInterpolation(src.first, src.second, src.third,
00909                                      dest.first, dest.second, dest.third);
00910 }
00911 
00912 #if 0
00913 /*****************************************************************/
00914 /*                                                               */
00915 /*              resizeImageCubicInterpolation                    */
00916 /*                                                               */
00917 /*****************************************************************/
00918 
00919 /** \brief Resize image using the cardinal B-spline interpolation function.
00920 
00921     The function calls like \ref resizeImageSplineInterpolation() with
00922     \ref vigra::BSpline<3, double> and prefiltering as an interpolation kernel.
00923     The interpolated function has two continuous derivatives.
00924     (See \ref resizeImageSplineInterpolation() for more documentation)
00925 
00926     <b>\#include</b> <vigra/resizeimage.hxx><br>
00927     Namespace: vigra
00928 
00929 */
00930 template <class SrcIterator, class SrcAccessor,
00931           class DestIterator, class DestAccessor>
00932 void
00933 resizeImageCubicInterpolation(SrcIterator src_iter, SrcIterator src_iter_end, SrcAccessor src_acc,
00934                       DestIterator dest_iter, DestIterator dest_iter_end, DestAccessor dest_acc)
00935 {
00936     resizeImageSplineInterpolation(src_iter, src_iter_end, src_acc, dest_iter, dest_iter_end, dest_acc,
00937                                   BSpline<3, double>());
00938 }
00939 
00940 template <class SrcIterator, class SrcAccessor,
00941           class DestIterator, class DestAccessor>
00942 inline
00943 void
00944 resizeImageCubicInterpolation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
00945                       triple<DestIterator, DestIterator, DestAccessor> dest)
00946 {
00947     resizeImageCubicInterpolation(src.first, src.second, src.third,
00948                                    dest.first, dest.second, dest.third);
00949 }
00950 #endif
00951 
00952 /*****************************************************************/
00953 /*                                                               */
00954 /*              resizeImageCoscotInterpolation                   */
00955 /*                                                               */
00956 /*****************************************************************/
00957 
00958 /** \brief Resize image using the Coscot interpolation function.
00959 
00960     The function calls \ref resizeImageSplineInterpolation() with
00961     \ref vigra::CoscotFunction as an interpolation kernel.
00962     The interpolated function has one continuous derivative.
00963     (See \ref resizeImageSplineInterpolation() for more documentation)
00964 
00965     <b> Declarations:</b>
00966 
00967     pass arguments explicitly:
00968     \code
00969     namespace vigra {
00970         template <class SrcIterator, class SrcAccessor,
00971                   class DestIterator, class DestAccessor>
00972         void
00973         resizeImageCoscotInterpolation(SrcIterator src_iter, SrcIterator src_iter_end, SrcAccessor src_acc,
00974                               DestIterator dest_iter, DestIterator dest_iter_end, DestAccessor dest_acc);
00975     }
00976     \endcode
00977 
00978 
00979     use argument objects in conjunction with \ref ArgumentObjectFactories :
00980     \code
00981     namespace vigra {
00982         template <class SrcIterator, class SrcAccessor,
00983                   class DestIterator, class DestAccessor>
00984         void
00985         resizeImageCoscotInterpolation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
00986                               triple<DestIterator, DestIterator, DestAccessor> dest);
00987     }
00988     \endcode
00989 
00990 
00991     <b>\#include</b> <vigra/resizeimage.hxx><br>
00992     Namespace: vigra
00993 
00994 */
00995 doxygen_overloaded_function(template <...> void resizeImageCoscotInterpolation)
00996 
00997 template <class SrcIterator, class SrcAccessor,
00998           class DestIterator, class DestAccessor>
00999 void
01000 resizeImageCoscotInterpolation(SrcIterator src_iter, SrcIterator src_iter_end, SrcAccessor src_acc,
01001                       DestIterator dest_iter, DestIterator dest_iter_end, DestAccessor dest_acc)
01002 {
01003     resizeImageSplineInterpolation(src_iter, src_iter_end, src_acc, dest_iter, dest_iter_end, dest_acc,
01004                                    CoscotFunction<double>());
01005 }
01006 
01007 template <class SrcIterator, class SrcAccessor,
01008           class DestIterator, class DestAccessor>
01009 inline
01010 void
01011 resizeImageCoscotInterpolation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
01012                       triple<DestIterator, DestIterator, DestAccessor> dest)
01013 {
01014     resizeImageCoscotInterpolation(src.first, src.second, src.third,
01015                                    dest.first, dest.second, dest.third);
01016 }
01017 
01018 
01019 #if 0 // old version of the spline interpolation algorithm
01020 
01021 /********************************************************/
01022 /*                                                      */
01023 /*           resizeCalculateSplineCoefficients          */
01024 /*         (internally used by resize functions)        */
01025 /*                                                      */
01026 /********************************************************/
01027 
01028 template <class SrcIterator, class SrcAccessor, class VALUETYPE>
01029 void
01030 resizeCalculateSplineCoefficients(SrcIterator i1, SrcIterator iend,
01031                 SrcAccessor a, VALUETYPE * i2)
01032 {
01033     int n = iend - i1;
01034 
01035     if(n <= 0) return;
01036 
01037     VALUETYPE zero = NumericTraits<VALUETYPE>::zero();
01038     VALUETYPE two = 2.0 * NumericTraits<VALUETYPE>::one();
01039     VALUETYPE half = 0.5 * NumericTraits<VALUETYPE>::one();
01040 
01041     *i2 = zero;
01042     if(n == 1) return;
01043 
01044     std::vector<VALUETYPE> vec(n);
01045     typename std::vector<VALUETYPE>::iterator u = vec.begin();
01046 
01047     *u = zero;
01048 
01049     for(++i1, ++i2, ++u, --iend; i1 != iend; ++i1, ++i2, ++u)
01050     {
01051         VALUETYPE p = 0.5 * i2[-1] + two;
01052         *i2 = half / p;
01053         *u = 3.0 *(a(i1,1) - 2.0 * a(i1) + a(i1, -1)) - 0.5 * u[-1] / p;
01054     }
01055 
01056     *i2 = zero;
01057 
01058     for(--i2, --u; u != vec; --u, --i2)
01059     {
01060         *i2 = *i2 * i2[1] + *u;
01061     }
01062 }
01063 
01064 /********************************************************/
01065 /*                                                      */
01066 /*         resizeImageInternalSplineGradient            */
01067 /*                                                      */
01068 /********************************************************/
01069 
01070 template <class SrcIterator, class SrcAccessor,
01071           class DoubleIterator, class TempIterator, class DestIterator>
01072 void
01073 resizeImageInternalSplineGradient(SrcIterator in, SrcIterator inend, SrcAccessor sa,
01074                          DoubleIterator tmp, TempIterator r, DestIterator id)
01075 {
01076     int w = inend - in;
01077 
01078     int x;
01079 
01080     typedef typename SrcAccessor::value_type SRCVT;
01081     typedef typename NumericTraits<SRCVT>::RealPromote TMPTYPE;
01082 
01083     // calculate border derivatives
01084     SrcIterator xs = in;
01085     TMPTYPE p0 = -11.0/6.0 * sa(xs);  ++xs;
01086             p0 += 3.0 * sa(xs);  ++xs;
01087             p0 += -1.5 * sa(xs);  ++xs;
01088             p0 += 1.0/3.0 * sa(xs);
01089 
01090     xs = in + w-1;
01091     TMPTYPE pw = 11.0/6.0 * sa(xs);  --xs;
01092             pw += -3.0 * sa(xs);  --xs;
01093             pw +=  1.5 * sa(xs);  --xs;
01094             pw += -1.0/3.0 * sa(xs);
01095 
01096     xs = in + 2;
01097     SrcIterator xs1 = in;
01098 
01099     for(x=1; x<w-1; ++x, ++xs, ++xs1)
01100     {
01101         r[x] = 3.0 * (sa(xs) - sa(xs1));
01102     }
01103 
01104     r[1] -= p0;
01105     r[w-2] -= pw;
01106 
01107     double q = 0.25;
01108 
01109     id[0] = p0;
01110     id[w-1] = pw;
01111     id[1] = 0.25 * r[1];
01112 
01113     for(x=2; x<w-1; ++x)
01114     {
01115         tmp[x] = q;
01116         q = 1.0 / (4.0 - q);
01117         id[x] = q * (r[x] - id[x-1]);
01118     }
01119 
01120     for(x=w-3; x>=1; --x)
01121     {
01122         id[x] -= tmp[x+1]*id[x+1];
01123     }
01124 }
01125 
01126 /********************************************************/
01127 /*                                                      */
01128 /*         resizeImageInternalSplineInterpolation       */
01129 /*                                                      */
01130 /********************************************************/
01131 
01132 template <class SrcIterator, class SrcAccessor,
01133           class DestIterator, class DestAccessor>
01134 void
01135 resizeImageInternalSplineInterpolation(SrcIterator is, SrcIterator iend, SrcAccessor sa,
01136                       DestIterator id, DestIterator idend, DestAccessor da)
01137 {
01138     int w = iend.x - is.x;
01139     int h = iend.y - is.y;
01140 
01141     int wnew = idend.x - id.x;
01142     int hnew = idend.y - id.y;
01143 
01144     typedef typename SrcAccessor::value_type SRCVT;
01145     typedef typename NumericTraits<SRCVT>::RealPromote TMPTYPE;
01146     typedef typename BasicImage<TMPTYPE>::Iterator TMPITER;
01147     typedef
01148         NumericTraits<typename DestAccessor::value_type> DestTraits;
01149 
01150     BasicImage<TMPTYPE> dx(w,h);
01151     BasicImage<TMPTYPE> dy(w,h);
01152     BasicImage<TMPTYPE> dxy(w,h);
01153     BasicImage<TMPTYPE> W(4,4), W1(4,4);
01154     std::vector<TMPTYPE> R(w > h ? w : h);
01155     std::vector<double> tmp(w > h ? w : h);
01156 
01157     typename BasicImage<TMPTYPE>::Accessor ta;
01158 
01159     SrcIterator in = is;
01160 
01161     TMPITER idx = dx.upperLeft();
01162     TMPITER idy = dy.upperLeft();
01163     TMPITER idxy = dxy.upperLeft();
01164     typename std::vector<TMPTYPE>::iterator r = R.begin();
01165     typename std::vector<double>::iterator it = tmp.begin();
01166 
01167     double ig[] = { 1.0, 0.0, -3.0,  2.0,
01168                     0.0, 1.0, -2.0,  1.0,
01169                     0.0, 0.0,  3.0, -2.0,
01170                     0.0, 0.0, -1.0,  1.0 };
01171 
01172     int x, y, i, j, k;
01173 
01174 
01175     // calculate x derivatives
01176     for(y=0; y<h; ++y, ++in.y, ++idx.y)
01177     {
01178         typename SrcIterator::row_iterator sr = in.rowIterator();
01179         typename TMPITER::row_iterator dr = idx.rowIterator();
01180         resizeImageInternalSplineGradient(sr, sr+w, sa,
01181                                           it, r, dr);
01182     }
01183 
01184     in = is;
01185 
01186     // calculate y derivatives
01187     for(x=0; x<w; ++x, ++in.x, ++idy.x)
01188     {
01189         typename SrcIterator::column_iterator sc = in.columnIterator();
01190         typename TMPITER::column_iterator dc = idy.columnIterator();
01191         resizeImageInternalSplineGradient(sc, sc+h, sa,
01192                                           it, r, dc);
01193     }
01194 
01195     in = is;
01196     idy = dy.upperLeft();
01197 
01198     // calculate mixed derivatives
01199     for(y=0; y<h; ++y, ++idy.y, ++idxy.y)
01200     {
01201         typename TMPITER::row_iterator sr = idy.rowIterator();
01202         typename TMPITER::row_iterator dr = idxy.rowIterator();
01203         resizeImageInternalSplineGradient(sr, sr+w, ta,
01204                                           it, r, dr);
01205     }
01206 
01207     double du = (double)(w-1) / (wnew-1);
01208     double dv = (double)(h-1) / (hnew-1);
01209     double ov = 0.0;
01210     int oy = 0;
01211     int yy = oy;
01212 
01213     DestIterator xxd = id, yyd = id;
01214 
01215     static Diff2D down(0,1), right(1,0), downright(1,1);
01216 
01217     for(y=0; y<h-1; ++y, ++in.y, ov -= 1.0)
01218     {
01219         if(y < h-2 && ov >= 1.0) continue;
01220         int y1 = y+1;
01221         double v = ov;
01222         double ou = 0.0;
01223         int ox = 0;
01224         int xx = ox;
01225 
01226         SrcIterator xs = in;
01227         for(x=0; x<w-1; ++x, ++xs.x, ou -= 1.0)
01228         {
01229             if(x < w-2 && ou >= 1.0) continue;
01230             int x1 = x+1;
01231             double u = ou;
01232 
01233             DestIterator xd = id + Diff2D(ox,oy);
01234             W[0][0] = sa(xs);
01235             W[0][1] = dy(x, y);
01236             W[0][2] = sa(xs, down);
01237             W[0][3] = dy(x, y1);
01238             W[1][0] = dx(x, y);
01239             W[1][1] = dxy(x, y);
01240             W[1][2] = dx(x, y1);
01241             W[1][3] = dxy(x, y1);
01242             W[2][0] = sa(xs, right);
01243             W[2][1] = dy(x1,y);
01244             W[2][2] = sa(xs, downright);
01245             W[2][3] = dy(x1, y1);
01246             W[3][0] = dx(x1, y);
01247             W[3][1] = dxy(x1, y);
01248             W[3][2] = dx(x1, y1);
01249             W[3][3] = dxy(x1, y1);
01250 
01251             for(i=0; i<4; ++i)
01252             {
01253                 for(j=0; j<4; ++j)
01254                 {
01255                     W1[j][i] = ig[j] * W[0][i];
01256                     for(k=1; k<4; ++k)
01257                     {
01258                         W1[j][i] += ig[j+4*k] * W[k][i];
01259                     }
01260                 }
01261             }
01262             for(i=0; i<4; ++i)
01263             {
01264                 for(j=0; j<4; ++j)
01265                 {
01266                     W[j][i] = ig[i] * W1[j][0];
01267                     for(k=1; k<4; ++k)
01268                     {
01269                        W[j][i] += ig[4*k+i] * W1[j][k];
01270                     }
01271                 }
01272             }
01273 
01274             TMPTYPE a1,a2,a3,a4;
01275 
01276             yyd = xd;
01277             for(v=ov, yy=oy; v<1.0; v+=dv, ++yyd.y, ++yy)
01278             {
01279                 a1 = W[0][0] + v * (W[0][1] +
01280                                v * (W[0][2] + v * W[0][3]));
01281                 a2 = W[1][0] + v * (W[1][1] +
01282                                v * (W[1][2] + v * W[1][3]));
01283                 a3 = W[2][0] + v * (W[2][1] +
01284                                v * (W[2][2] + v * W[2][3]));
01285                 a4 = W[3][0] + v * (W[3][1] +
01286                                v * (W[3][2] + v * W[3][3]));
01287 
01288                 xxd = yyd;
01289                 for(u=ou, xx=ox; u<1.0; u+=du, ++xxd.x, ++xx)
01290                 {
01291                     da.set(DestTraits::fromRealPromote(a1 + u * (a2 + u * (a3 + u * a4))), xxd);
01292                 }
01293 
01294                 if(xx == wnew-1)
01295                 {
01296                     da.set(DestTraits::fromRealPromote(a1 + a2 + a3 + a4), xxd);
01297                 }
01298             }
01299 
01300             if(yy == hnew-1)
01301             {
01302                 a1 = W[0][0] + W[0][1] + W[0][2] + W[0][3];
01303                 a2 = W[1][0] + W[1][1] + W[1][2] + W[1][3];
01304                 a3 = W[2][0] + W[2][1] + W[2][2] + W[2][3];
01305                 a4 = W[3][0] + W[3][1] + W[3][2] + W[3][3];
01306 
01307                 DestIterator xxd = yyd;
01308                 for(u=ou, xx=ox; u<1.0; u+=du, ++xxd.x, ++xx)
01309                 {
01310                     da.set(DestTraits::fromRealPromote(a1 + u * (a2 + u * (a3 + u * a4))), xxd);
01311                 }
01312 
01313                 if(xx == wnew-1)
01314                 {
01315                     da.set(DestTraits::fromRealPromote(a1 + a2 + a3 + a4), xxd);
01316                 }
01317             }
01318 
01319             ou = u;
01320             ox = xx;
01321         }
01322         ov = v;
01323         oy = yy;
01324     }
01325 }
01326 
01327 template <class SrcIterator, class SrcAccessor,
01328           class DestIterator, class DestAccessor>
01329 void
01330 resizeImageSplineInterpolation(SrcIterator is, SrcIterator iend, SrcAccessor sa,
01331                       DestIterator id, DestIterator idend, DestAccessor da)
01332 {
01333     int w = iend.x - is.x;
01334     int h = iend.y - is.y;
01335 
01336     int wnew = idend.x - id.x;
01337     int hnew = idend.y - id.y;
01338 
01339     vigra_precondition((w > 3) && (h > 3),
01340                  "resizeImageSplineInterpolation(): "
01341                  "Source image too small.\n");
01342     vigra_precondition((wnew > 1) && (hnew > 1),
01343                  "resizeImageSplineInterpolation(): "
01344                  "Destination image too small.\n");
01345 
01346     double scale = 2.0;
01347 
01348     if(wnew < w || hnew < h)
01349     {
01350         typedef typename SrcAccessor::value_type SRCVT;
01351         typedef typename NumericTraits<SRCVT>::RealPromote TMPTYPE;
01352         typedef typename BasicImage<TMPTYPE>::Iterator TMPITER;
01353 
01354         BasicImage<TMPTYPE> t(w,h);
01355         TMPITER it = t.upperLeft();
01356 
01357         if(wnew < w)
01358         {
01359             recursiveSmoothX(is, iend, sa,
01360                     it, t.accessor(), (double)w/wnew/scale);
01361 
01362             if(hnew < h)
01363             {
01364                recursiveSmoothY(it, t.lowerRight(), t.accessor(),
01365                     it, t.accessor(), (double)h/hnew/scale);
01366             }
01367         }
01368         else
01369         {
01370            recursiveSmoothY(is, iend, sa,
01371                     it, t.accessor(), (double)h/hnew/scale);
01372         }
01373 
01374         resizeImageInternalSplineInterpolation(it, t.lowerRight(), t.accessor(),
01375                                                id, idend, da);
01376     }
01377     else
01378     {
01379         resizeImageInternalSplineInterpolation(is, iend, sa, id, idend, da);
01380     }
01381 }
01382 
01383 template <class SrcIterator, class SrcAccessor,
01384           class DestIterator, class DestAccessor>
01385 inline
01386 void
01387 resizeImageSplineInterpolation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
01388                       triple<DestIterator, DestIterator, DestAccessor> dest)
01389 {
01390     resizeImageSplineInterpolation(src.first, src.second, src.third,
01391                                    dest.first, dest.second, dest.third);
01392 }
01393 #endif  // old algorithm version
01394 
01395 //@}
01396 
01397 } // namespace vigra
01398 
01399 #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.9.0 (Tue Nov 6 2012)