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

| vigra/resampling_convolution.hxx |
--- | --- | ---

    00001 /************************************************************************/
    00002 /*                                                                      */
    00003 /*               Copyright 1998-2004 by Ullrich Koethe                  */
    00004 /*       Cognitive Systems Group, University of Hamburg, Germany        */
    00005 /*                                                                      */
    00006 /*    This file is part of the VIGRA computer vision library.           */
    00007 /*    ( Version 1.3.2, Jan 27 2005 )                                    */
    00008 /*    You may use, modify, and distribute this software according       */
    00009 /*    to the terms stated in the LICENSE file included in               */
    00010 /*    the VIGRA distribution.                                           */
    00011 /*                                                                      */
    00012 /*    The VIGRA Website is                                              */
    00013 /*        http://kogs-www.informatik.uni-hamburg.de/~koethe/vigra/      */
    00014 /*    Please direct questions, bug reports, and contributions to        */
    00015 /*        koethe@informatik.uni-hamburg.de                              */
    00016 /*                                                                      */
    00017 /*  THIS SOFTWARE IS PROVIDED AS IS AND WITHOUT ANY EXPRESS OR          */
    00018 /*  IMPLIED WARRANTIES, INCLUDING, WITHOUT LIMITATION, THE IMPLIED      */
    00019 /*  WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. */
    00020 /*                                                                      */
    00021 /************************************************************************/
    00022
    00023 #ifndef VIGRA_RESAMPLING_CONVOLUTION_HXX
    00024 #define VIGRA_RESAMPLING_CONVOLUTION_HXX
    00025
    00026 #include <cmath>
    00027 #include "vigra/stdimage.hxx"
    00028 #include "vigra/array_vector.hxx"
    00029 #include "vigra/rational.hxx"
    00030 #include "vigra/functortraits.hxx"
    00031
    00032 namespace vigra {
    00033
    00034 namespace resampling_detail
    00035 {
    00036
    00037 struct MapTargetToSourceCoordinate
    00038 {
    00039     MapTargetToSourceCoordinate(Rational<int> const & samplingRatio,
    00040                                 Rational<int> const & offset)
    00041     : a(samplingRatio.denominator()*offset.denominator()),
    00042       b(samplingRatio.numerator()*offset.numerator()),
    00043       c(samplingRatio.numerator()*offset.denominator())
    00044     {}
    00045
    00046 //        the following funcions are more efficient realizations of:
    00047 //             rational_cast<T>(i / samplingRatio + offset);
    00048 //        we need efficiency because this may be called in the inner loop
    00049
    00050     int operator()(int i) const
    00051     {
    00052         return (i * a + b) / c;
    00053     }
    00054
    00055     double toDouble(int i) const
    00056     {
    00057         return double(i * a + b) / c;
    00058     }
    00059
    00060     Rational<int> toRational(int i) const
    00061     {
    00062         return Rational<int>(i * a + b, c);
    00063     }
    00064
    00065     int a, b, c;
    00066 };
    00067
    00068 } // namespace resampling_detail
    00069
    00070 template <>
    00071 class FunctorTraits<resampling_detail::MapTargetToSourceCoordinate>
    00072 : public FunctorTraitsBase<resampling_detail::MapTargetToSourceCoordinate>
    00073 {
    00074   public:
    00075     typedef VigraTrueType isUnaryFunctor;
    00076 };
    00077
    00078 template <class SrcIter, class SrcAcc,
    00079           class DestIter, class DestAcc,
    00080           class KernelArray,
    00081           class Functor>
    00082 void
    00083 resamplingConvolveLine(SrcIter s, SrcIter send, SrcAcc src,
    00084                        DestIter d, DestIter dend, DestAcc dest,
    00085                        KernelArray const & kernels,
    00086                        Functor mapTargetToSourceCoordinate)
    00087 {
    00088     typedef typename
    00089         NumericTraits<typename SrcAcc::value_type>::RealPromote
    00090         TmpType;
    00091     typedef typename KernelArray::value_type Kernel;
    00092     typedef typename Kernel::const_iterator KernelIter;
    00093
    00094     int wo = send - s;
    00095     int wn = dend - d;
    00096     int wo2 = 2*wo - 2;
    00097
    00098     int i;
    00099     typename KernelArray::const_iterator kernel = kernels.begin();
    00100     for(i=0; i<wn; ++i, ++d, ++kernel)
    00101     {
    00102         // use the kernels periodically
    00103         if(kernel == kernels.end())
    00104             kernel = kernels.begin();
    00105
    00106         // calculate current target point into source location
    00107         int is = mapTargetToSourceCoordinate(i);
    00108
    00109         TmpType sum = NumericTraits<TmpType>::zero();
    00110
    00111         int lbound = is - kernel->right(),
    00112             hbound = is - kernel->left();
    00113
    00114         KernelIter k = kernel->center() + kernel->right();
    00115         if(lbound < 0 || hbound >= wo)
    00116         {
    00117             vigra_precondition(-lbound < wo && wo2 - hbound >= 0,
    00118                 "resamplingConvolveLine(): kernel or offset larger than image.");
    00119             for(int m=lbound; m <= hbound; ++m, --k)
    00120             {
    00121                 int mm = (m < 0) ?
    00122                             -m :
    00123                             (m >= wo) ?
    00124                                 wo2 - m :
    00125                                 m;
    00126                 sum += *k * src(s, mm);
    00127             }
    00128         }
    00129         else
    00130         {
    00131             SrcIter ss = s + lbound;
    00132             SrcIter ssend = s + hbound;
    00133
    00134             for(; ss <= ssend; ++ss, --k)
    00135             {
    00136                 sum += *k * src(ss);
    00137             }
    00138         }
    00139
    00140         dest.set(sum, d);
    00141     }
    00142 }
    00143
    00144 template <class Kernel, class MapCoordinate, class KernelArray>
    00145 void
    00146 createResamplingKernels(Kernel const & kernel,
    00147              MapCoordinate const & mapCoordinate, KernelArray & kernels)
    00148 {
    00149     for(unsigned int idest = 0; idest < kernels.size(); ++idest)
    00150     {
    00151         int isrc = mapCoordinate(idest);
    00152         double idsrc = mapCoordinate.toDouble(idest);
    00153         double offset = idsrc - isrc;
    00154         double radius = kernel.radius();
    00155         int left = int(ceil(-radius - offset));
    00156         int right = int(floor(radius - offset));
    00157         kernels[idest].initExplicitly(left, right);
    00158
    00159         double x = left + offset;
    00160         for(int i = left; i <= right; ++i, ++x)
    00161             kernels[idest][i] = kernel(x);
    00162         kernels[idest].normalize(1.0, kernel.derivativeOrder(), offset);
    00163     }
    00164 }
    00165
    00166 /** \addtogroup ResamplingConvolutionFilters Resampling Convolution Filters
    00167
    00168     These functions implement the convolution operation when the source and target images
    00169     have different sizes. This is realized by accessing a continous kernel at the
    00170     appropriate non-integer positions. The technique is, for example, described in
    00171     D. Schumacher: <i>General Filtered Image Rescaling</i>, in: Graphics Gems III,
    00172     Academic Press, 1992.
    00173 */
    00174 //@{
    00175
    00176 /********************************************************/
    00177 /*                                                      */
    00178 /*                  resamplingConvolveX                 */
    00179 /*                                                      */
    00180 /********************************************************/
    00181
    00182 /** \brief Apply a resampling filter in the x-direction.
    00183
    00184     This function implements a convolution operation in x-direction
    00185     (i.e. applies a 1D filter to every row) where the width of the source
    00186     and destination images differ. This is typically used to avoid aliasing if
    00187     the image is scaled down, or to interpolate smoothly if the image is scaled up.
    00188     The target coordinates are transformed into source coordinates by
    00189
    00190     \code
    00191     xsource = (xtarget - offset) / samplingRatio
    00192     \endcode
    00193
    00194     The <tt>samplingRatio</tt> and <tt>offset</tt> must be given as \ref vigra::Rational
    00195     in order to avoid rounding errors in this transformation. It is required that for all
    00196     pixels of the target image, <tt>xsource</tt> remains within the range of the source
    00197     image (i.e. <tt>0 <= xsource <= sourceWidth-1</tt>. Since <tt>xsource</tt> is
    00198     in general not an integer, the <tt>kernel</tt> must be a functor that can be accessed at
    00199     arbitrary (<tt>double</tt>) coordinates. It must also provide a member function <tt>radius()</tt>
    00200     which specifies the support (non-zero interval) of the kernel. VIGRA already
    00201     provides a number of suitable functors, e.g. \ref vigra::Gaussian, \ref vigra::BSpline
    00202     \ref vigra::CatmullRomSpline, and \ref vigra::CoscotFunction. The function
    00203     \ref resizeImageSplineInterpolation() is implemented by means resamplingConvolveX() and
    00204     resamplingConvolveY().
    00205
    00206     <b> Declarations:</b>
    00207
    00208     pass arguments explicitly:
    00209     \code
    00210     namespace vigra {
    00211         template <class SrcIter, class SrcAcc,
    00212                   class DestIter, class DestAcc,
    00213                   class Kernel>
    00214         void
    00215         resamplingConvolveX(SrcIter sul, SrcIter slr, SrcAcc src,
    00216                             DestIter dul, DestIter dlr, DestAcc dest,
    00217                             Kernel const & kernel,
    00218                             Rational<int> const & samplingRatio, Rational<int> const & offset);
    00219     }
    00220     \endcode
    00221
    00222
    00223     use argument objects in conjunction with \ref ArgumentObjectFactories:
    00224     \code
    00225     namespace vigra {
    00226         template <class SrcIter, class SrcAcc,
    00227                   class DestIter, class DestAcc,
    00228                   class Kernel>
    00229         void
    00230         resamplingConvolveX(triple<SrcIter, SrcIter, SrcAcc> src,
    00231                             triple<DestIter, DestIter, DestAcc> dest,
    00232                             Kernel const & kernel,
    00233                             Rational<int> const & samplingRatio, Rational<int> const & offset);
    00234     }
    00235     \endcode
    00236
    00237     <b> Usage:</b>
    00238
    00239     <b>\#include</b> "<a href="resampling_convolution_8hxx-source.html">vigra/resampling_convolution.hxx</a>"
    00240
    00241
    00242     \code
    00243     Rational<int> ratio(2), offset(0);
    00244
    00245     FImage src(w,h),
    00246            dest(rational_cast<int>(ratio*w), h);
    00247
    00248     float sigma = 2.0;
    00249     Gaussian<float> smooth(sigma);
    00250     ...
    00251
    00252     // simpultaneously enlarge and smooth source image
    00253     resamplingConvolveX(srcImageRange(src), destImageRange(dest),
    00254                         smooth, ratio, offset);
    00255     \endcode
    00256
    00257     <b> Required Interface:</b>
    00258
    00259     \code
    00260     Kernel kernel;
    00261     int kernelRadius = kernel.radius();
    00262     double x = ...;  // must be <= radius()
    00263     double value = kernel(x);
    00264     \endcode
    00265 */
    00266 template <class SrcIter, class SrcAcc,
    00267           class DestIter, class DestAcc,
    00268           class Kernel>
    00269 void
    00270 resamplingConvolveX(SrcIter sul, SrcIter slr, SrcAcc src,
    00271                     DestIter dul, DestIter dlr, DestAcc dest,
    00272                     Kernel const & kernel,
    00273                     Rational<int> const & samplingRatio, Rational<int> const & offset)
    00274 {
    00275     int wold = slr.x - sul.x;
    00276     int wnew = dlr.x - dul.x;
    00277
    00278     vigra_precondition(!samplingRatio.is_inf() && samplingRatio > 0,
    00279                 "resamplingConvolveX(): sampling ratio must be > 0 and < infinity");
    00280     vigra_precondition(!offset.is_inf(),
    00281                 "resamplingConvolveX(): offset must be < infinity");
    00282
    00283     int period = lcm(samplingRatio.numerator(), samplingRatio.denominator());
    00284     resampling_detail::MapTargetToSourceCoordinate mapCoordinate(samplingRatio, offset);
    00285
    00286     ArrayVector<Kernel1D<double> > kernels(period);
    00287
    00288     createResamplingKernels(kernel, mapCoordinate, kernels);
    00289
    00290     for(; sul.y < slr.y; ++sul.y, ++dul.y)
    00291     {
    00292         typename SrcIter::row_iterator sr = sul.rowIterator();
    00293         typename DestIter::row_iterator dr = dul.rowIterator();
    00294         resamplingConvolveLine(sr, sr+wold, src, dr, dr+wnew, dest,
    00295                                kernels, mapCoordinate);
    00296     }
    00297 }
    00298
    00299 template <class SrcIter, class SrcAcc,
    00300           class DestIter, class DestAcc,
    00301           class Kernel>
    00302 inline void
    00303 resamplingConvolveX(triple<SrcIter, SrcIter, SrcAcc> src,
    00304                     triple<DestIter, DestIter, DestAcc> dest,
    00305                     Kernel const & kernel,
    00306                     Rational<int> const & samplingRatio, Rational<int> const & offset)
    00307 {
    00308     resamplingConvolveX(src.first, src.second, src.third,
    00309                         dest.first, dest.second, dest.third,
    00310                         kernel, samplingRatio, offset);
    00311 }
    00312
    00313 /********************************************************/
    00314 /*                                                      */
    00315 /*                  resamplingConvolveY                 */
    00316 /*                                                      */
    00317 /********************************************************/
    00318
    00319 /** \brief Apply a resampling filter in the y-direction.
    00320
    00321     This function implements a convolution operation in y-direction
    00322     (i.e. applies a 1D filter to every column) where the height of the source
    00323     and destination images differ. This is typically used to avoid aliasing if
    00324     the image is scaled down, or to interpolate smoothly if the image is scaled up.
    00325     The target coordinates are transformed into source coordinates by
    00326
    00327     \code
    00328     ysource = (ytarget - offset) / samplingRatio
    00329     \endcode
    00330
    00331     The <tt>samplingRatio</tt> and <tt>offset</tt> must be given as \ref vigra::Rational
    00332     in order to avoid rounding errors in this transformation. It is required that for all
    00333     pixels of the target image, <tt>ysource</tt> remains within the range of the source
    00334     image (i.e. <tt>0 <= ysource <= sourceHeight-1</tt>. Since <tt>ysource</tt> is
    00335     in general not an integer, the <tt>kernel</tt> must be a functor that can be accessed at
    00336     arbitrary (<tt>double</tt>) coordinates. It must also provide a member function <tt>radius()</tt>
    00337     which specifies the support (non-zero interval) of the kernel. VIGRA already
    00338     provides a number of suitable functors, e.g. \ref vigra::Gaussian, \ref vigra::BSpline
    00339     \ref vigra::CatmullRomSpline, and \ref vigra::CoscotFunction. The function
    00340     \ref resizeImageSplineInterpolation() is implemented by means resamplingConvolveX() and
    00341     resamplingConvolveY().
    00342
    00343     <b> Declarations:</b>
    00344
    00345     pass arguments explicitly:
    00346     \code
    00347     namespace vigra {
    00348         template <class SrcIter, class SrcAcc,
    00349                   class DestIter, class DestAcc,
    00350                   class Kernel>
    00351         void
    00352         resamplingConvolveY(SrcIter sul, SrcIter slr, SrcAcc src,
    00353                             DestIter dul, DestIter dlr, DestAcc dest,
    00354                             Kernel const & kernel,
    00355                             Rational<int> const & samplingRatio, Rational<int> const & offset);
    00356     }
    00357     \endcode
    00358
    00359
    00360     use argument objects in conjunction with \ref ArgumentObjectFactories:
    00361     \code
    00362     namespace vigra {
    00363         template <class SrcIter, class SrcAcc,
    00364                   class DestIter, class DestAcc,
    00365                   class Kernel>
    00366         void
    00367         resamplingConvolveY(triple<SrcIter, SrcIter, SrcAcc> src,
    00368                             triple<DestIter, DestIter, DestAcc> dest,
    00369                             Kernel const & kernel,
    00370                             Rational<int> const & samplingRatio, Rational<int> const & offset);
    00371     }
    00372     \endcode
    00373
    00374     <b> Usage:</b>
    00375
    00376     <b>\#include</b> "<a href="resampling_convolution_8hxx-source.html">vigra/resampling_convolution.hxx</a>"
    00377
    00378
    00379     \code
    00380     Rational<int> ratio(2), offset(0);
    00381
    00382     FImage src(w,h),
    00383            dest(w, rational_cast<int>(ratio*h));
    00384
    00385     float sigma = 2.0;
    00386     Gaussian<float> smooth(sigma);
    00387     ...
    00388
    00389     // simpultaneously enlarge and smooth source image
    00390     resamplingConvolveY(srcImageRange(src), destImageRange(dest),
    00391                         smooth, ratio, offset);
    00392     \endcode
    00393
    00394     <b> Required Interface:</b>
    00395
    00396     \code
    00397     Kernel kernel;
    00398     int kernelRadius = kernel.radius();
    00399     double y = ...;  // must be <= radius()
    00400     double value = kernel(y);
    00401     \endcode
    00402 */
    00403 template <class SrcIter, class SrcAcc,
    00404           class DestIter, class DestAcc,
    00405           class Kernel>
    00406 void
    00407 resamplingConvolveY(SrcIter sul, SrcIter slr, SrcAcc src,
    00408                     DestIter dul, DestIter dlr, DestAcc dest,
    00409                     Kernel const & kernel,
    00410                     Rational<int> const & samplingRatio, Rational<int> const & offset)
    00411 {
    00412     int hold = slr.y - sul.y;
    00413     int hnew = dlr.y - dul.y;
    00414
    00415     vigra_precondition(!samplingRatio.is_inf() && samplingRatio > 0,
    00416                 "resamplingConvolveY(): sampling ratio must be > 0 and < infinity");
    00417     vigra_precondition(!offset.is_inf(),
    00418                 "resamplingConvolveY(): offset must be < infinity");
    00419
    00420     int period = lcm(samplingRatio.numerator(), samplingRatio.denominator());
    00421
    00422     resampling_detail::MapTargetToSourceCoordinate mapCoordinate(samplingRatio, offset);
    00423
    00424     ArrayVector<Kernel1D<double> > kernels(period);
    00425
    00426     createResamplingKernels(kernel, mapCoordinate, kernels);
    00427
    00428     for(; sul.x < slr.x; ++sul.x, ++dul.x)
    00429     {
    00430         typename SrcIter::column_iterator sc = sul.columnIterator();
    00431         typename DestIter::column_iterator dc = dul.columnIterator();
    00432         resamplingConvolveLine(sc, sc+hold, src, dc, dc+hnew, dest,
    00433                                kernels, mapCoordinate);
    00434     }
    00435 }
    00436
    00437 template <class SrcIter, class SrcAcc,
    00438           class DestIter, class DestAcc,
    00439           class Kernel>
    00440 inline void
    00441 resamplingConvolveY(triple<SrcIter, SrcIter, SrcAcc> src,
    00442                     triple<DestIter, DestIter, DestAcc> dest,
    00443                     Kernel const & kernel,
    00444                     Rational<int> const & samplingRatio, Rational<int> const & offset)
    00445 {
    00446     resamplingConvolveY(src.first, src.second, src.third,
    00447                         dest.first, dest.second, dest.third,
    00448                         kernel, samplingRatio, offset);
    00449 }
    00450
    00451 /********************************************************/
    00452 /*                                                      */
    00453 /*               resamplingConvolveImage                */
    00454 /*                                                      */
    00455 /********************************************************/
    00456
    00457 /** \brief Apply two separable resampling filters successively, the first in x-direction,
    00458            the second in y-direction.
    00459
    00460     This function is a shorthand for the concatenation of a call to
    00461     \link ResamplingConvolutionFilters#resamplingConvolveX resamplingConvolveX\endlink()
    00462     and \link ResamplingConvolutionFilters#resamplingConvolveY resamplingConvolveY\endlink()
    00463     with the given kernels. See there for detailed documentation.
    00464
    00465     <b> Declarations:</b>
    00466
    00467     pass arguments explicitly:
    00468     \code
    00469     namespace vigra {
    00470         template <class SrcIterator, class SrcAccessor,
    00471                   class DestIterator, class DestAccessor,
    00472                   class KernelX, class KernelY>
    00473         void resamplingConvolveImage(SrcIterator sul,SrcIterator slr, SrcAccessor src,
    00474                            DestIterator dul, DestIterator dlr, DestAccessor dest,
    00475                            KernelX const & kx,
    00476                            Rational<int> const & samplingRatioX, Rational<int> const & offsetX,
    00477                            KernelY const & ky,
    00478                            Rational<int> const & samplingRatioY, Rational<int> const & offsetY);
    00479     }
    00480     \endcode
    00481
    00482
    00483     use argument objects in conjunction with \ref ArgumentObjectFactories:
    00484     \code
    00485     namespace vigra {
    00486         template <class SrcIterator, class SrcAccessor,
    00487                   class DestIterator, class DestAccessor,
    00488                   class KernelX, class KernelY>
    00489         void
    00490         resamplingConvolveImage(triple<SrcIterator, SrcIterator, SrcAccessor> src,
    00491                            triple<DestIterator, DestIterator, DestAccessor> dest,
    00492                            KernelX const & kx,
    00493                            Rational<int> const & samplingRatioX, Rational<int> const & offsetX,
    00494                            KernelY const & ky,
    00495                            Rational<int> const & samplingRatioY, Rational<int> const & offsetY);
    00496     }
    00497     \endcode
    00498
    00499     <b> Usage:</b>
    00500
    00501     <b>\#include</b> "<a href="resampling_convolution_8hxx-source.html">vigra/resampling_convolution.hxx</a>"
    00502
    00503
    00504     \code
    00505     Rational<int> xratio(2), yratio(3), offset(0);
    00506
    00507     FImage src(w,h),
    00508            dest(rational_cast<int>(xratio*w), rational_cast<int>(yratio*h));
    00509
    00510     float sigma = 2.0;
    00511     Gaussian<float> smooth(sigma);
    00512     ...
    00513
    00514     // simpultaneously enlarge and smooth source image
    00515     resamplingConvolveImage(srcImageRange(src), destImageRange(dest),
    00516                             smooth, xratio, offset,
    00517                             smooth, yratio, offset);
    00518
    00519     \endcode
    00520
    00521 */
    00522 template <class SrcIterator, class SrcAccessor,
    00523           class DestIterator, class DestAccessor,
    00524           class KernelX, class KernelY>
    00525 void resamplingConvolveImage(SrcIterator sul,SrcIterator slr, SrcAccessor src,
    00526                    DestIterator dul, DestIterator dlr, DestAccessor dest,
    00527                    KernelX const & kx,
    00528                    Rational<int> const & samplingRatioX, Rational<int> const & offsetX,
    00529                    KernelY const & ky,
    00530                    Rational<int> const & samplingRatioY, Rational<int> const & offsetY)
    00531 {
    00532     typedef typename
    00533         NumericTraits<typename SrcAccessor::value_type>::RealPromote
    00534         TmpType;
    00535
    00536     BasicImage<TmpType> tmp(dlr.x - dul.x, slr.y - sul.y);
    00537
    00538     resamplingConvolveX(srcIterRange(sul, slr, src),
    00539                         destImageRange(tmp),
    00540                         kx, samplingRatioX, offsetX);
    00541     resamplingConvolveY(srcImageRange(tmp),
    00542                         destIterRange(dul, dlr, dest),
    00543                         ky, samplingRatioY, offsetY);
    00544 }
    00545
    00546 template <class SrcIterator, class SrcAccessor,
    00547           class DestIterator, class DestAccessor,
    00548           class KernelX, class KernelY>
    00549 inline void
    00550 resamplingConvolveImage(triple<SrcIterator, SrcIterator, SrcAccessor> src,
    00551                    triple<DestIterator, DestIterator, DestAccessor> dest,
    00552                    KernelX const & kx,
    00553                    Rational<int> const & samplingRatioX, Rational<int> const & offsetX,
    00554                    KernelY const & ky,
    00555                    Rational<int> const & samplingRatioY, Rational<int> const & offsetY)
    00556 {
    00557     resamplingConvolveImage(src.first, src.second, src.third,
    00558                             dest.first, dest.second, dest.third,
    00559                             kx, samplingRatioX, offsetX,
    00560                             ky, samplingRatioY, offsetY);
    00561 }
    00562
    00563 } // namespace vigra
    00564
    00565
    00566 #endif /* VIGRA_RESAMPLING_CONVOLUTION_HXX */

© Ullrich Köthe (koethe@informatik.uni-hamburg.de) Cognitive Systems Group, University of Hamburg, Germany | html generated using doxygen and Python VIGRA 1.3.2 (27 Jan 2005)