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

vigra/stdconvolution.hxx VIGRA

00001 /************************************************************************/
00002 /*                                                                      */
00003 /*               Copyright 1998-2002 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_STDCONVOLUTION_HXX
00038 #define VIGRA_STDCONVOLUTION_HXX
00039 
00040 #include <cmath>
00041 #include "stdimage.hxx"
00042 #include "bordertreatment.hxx"
00043 #include "separableconvolution.hxx"
00044 #include "utilities.hxx"
00045 #include "sized_int.hxx"
00046 
00047 namespace vigra {
00048 
00049 /** \addtogroup StandardConvolution Two-dimensional convolution functions
00050 
00051 Perform 2D non-separable convolution, with and without ROI mask.
00052 
00053 These generic convolution functions implement
00054 the standard 2D convolution operation for images that fit
00055 into the required interface. Arbitrary ROI's are supported
00056 by the mask version of the algorithm.
00057 The functions need a suitable 2D kernel to operate.
00058 */
00059 //@{
00060 
00061 /** \brief Performs a 2 dimensional convolution of the source image using the given
00062     kernel.
00063 
00064     The KernelIterator must point to the center of the kernel, and
00065     the kernel's size is given by its upper left (x and y of distance <= 0) and
00066     lower right (distance >= 0) corners. The image must always be larger than the
00067     kernel. At those positions where the kernel does not completely fit
00068     into the image, the specified \ref BorderTreatmentMode is
00069     applied. You can choice between following BorderTreatmentModes:
00070     <ul>
00071     <li>BORDER_TREATMENT_CLIP</li>
00072     <li>BORDER_TREATMENT_AVOID</li>
00073     <li>BORDER_TREATMENT_WRAP</li>
00074     <li>BORDER_TREATMENT_REFLECT</li>
00075     <li>BORDER_TREATMENT_REPEAT</li>
00076     </ul><br>
00077     The images's pixel type (SrcAccessor::value_type) must be a
00078     linear space over the kernel's value_type (KernelAccessor::value_type),
00079     i.e. addition of source values, multiplication with kernel values,
00080     and NumericTraits must be defined.
00081     The kernel's value_type must be an algebraic field,
00082     i.e. the arithmetic operations (+, -, *, /) and NumericTraits must
00083     be defined.
00084 
00085     <b> Declarations:</b>
00086 
00087     pass arguments explicitly:
00088     \code
00089     namespace vigra {
00090         template <class SrcIterator, class SrcAccessor,
00091                   class DestIterator, class DestAccessor,
00092                   class KernelIterator, class KernelAccessor>
00093         void convolveImage(SrcIterator src_ul, SrcIterator src_lr, SrcAccessor src_acc,
00094                            DestIterator dest_ul, DestAccessor dest_acc,
00095                            KernelIterator ki, KernelAccessor ak,
00096                            Diff2D kul, Diff2D klr, BorderTreatmentMode border);
00097     }
00098     \endcode
00099 
00100 
00101     use argument objects in conjunction with \ref ArgumentObjectFactories :
00102     \code
00103     namespace vigra {
00104         template <class SrcIterator, class SrcAccessor,
00105                   class DestIterator, class DestAccessor,
00106                   class KernelIterator, class KernelAccessor>
00107         void convolveImage(triple<SrcIterator, SrcIterator, SrcAccessor> src,
00108                            pair<DestIterator, DestAccessor> dest,
00109                            tuple5<KernelIterator, KernelAccessor, Diff2D, Diff2D,
00110                            BorderTreatmentMode> kernel);
00111     }
00112     \endcode
00113 
00114     <b> Usage:</b>
00115 
00116     <b>\#include</b> <vigra/stdconvolution.hxx><br>
00117     Namespace: vigra
00118 
00119 
00120     \code
00121     vigra::FImage src(w,h), dest(w,h);
00122     ...
00123 
00124     // define horizontal Sobel filter
00125     vigra::Kernel2D<float> sobel;
00126 
00127     sobel.initExplicitly(Diff2D(-1,-1), Diff2D(1,1)) =  // upper left and lower right
00128                          0.125, 0.0, -0.125,
00129                          0.25,  0.0, -0.25,
00130                          0.125, 0.0, -0.125;
00131     sobel.setBorderTreatment(vigra::BORDER_TREATMENT_REFLECT);
00132 
00133     vigra::convolveImage(srcImageRange(src), destImage(dest), kernel2d(sobel));
00134     \endcode
00135 
00136     <b> Required Interface:</b>
00137 
00138     \code
00139     ImageIterator src_ul, src_lr;
00140     ImageIterator dest_ul;
00141     ImageIterator ik;
00142 
00143     SrcAccessor src_accessor;
00144     DestAccessor dest_accessor;
00145     KernelAccessor kernel_accessor;
00146 
00147     NumericTraits<SrcAccessor::value_type>::RealPromote s = src_accessor(src_ul);
00148 
00149     s = s + s;
00150     s = kernel_accessor(ik) * s;
00151     s -= s;
00152 
00153     dest_accessor.set(
00154     NumericTraits<DestAccessor::value_type>::fromRealPromote(s), dest_ul);
00155 
00156     NumericTraits<KernelAccessor::value_type>::RealPromote k = kernel_accessor(ik);
00157 
00158     k += k;
00159     k -= k;
00160     k = k / k;
00161 
00162     \endcode
00163 
00164     <b> Preconditions:</b>
00165 
00166     \code
00167     kul.x <= 0
00168     kul.y <= 0
00169     klr.x >= 0
00170     klr.y >= 0
00171     src_lr.x - src_ul.x >= klr.x + kul.x + 1
00172     src_lr.y - src_ul.y >= klr.y + kul.y + 1
00173     \endcode
00174 
00175     If border == BORDER_TREATMENT_CLIP: Sum of kernel elements must be
00176     != 0.
00177 
00178 */
00179 template <class SrcIterator, class SrcAccessor,
00180           class DestIterator, class DestAccessor,
00181           class KernelIterator, class KernelAccessor>
00182 void convolveImage(SrcIterator src_ul, SrcIterator src_lr, SrcAccessor src_acc,
00183                    DestIterator dest_ul, DestAccessor dest_acc,
00184                    KernelIterator ki, KernelAccessor ak,
00185                    Diff2D kul, Diff2D klr, BorderTreatmentMode border)
00186 {
00187     vigra_precondition((border == BORDER_TREATMENT_CLIP    ||
00188                         border == BORDER_TREATMENT_AVOID   ||
00189                         border == BORDER_TREATMENT_REFLECT ||
00190                         border == BORDER_TREATMENT_REPEAT  ||
00191                         border == BORDER_TREATMENT_WRAP),
00192                        "convolveImage():\n"
00193                        "  Border treatment must be one of follow treatments:\n"
00194                        "  - BORDER_TREATMENT_CLIP\n"
00195                        "  - BORDER_TREATMENT_AVOID\n"
00196                        "  - BORDER_TREATMENT_REFLECT\n"
00197                        "  - BORDER_TREATMENT_REPEAT\n"
00198                        "  - BORDER_TREATMENT_WRAP\n");
00199 
00200     vigra_precondition(kul.x <= 0 && kul.y <= 0,
00201                        "convolveImage(): coordinates of "
00202                        "kernel's upper left must be <= 0.");
00203     vigra_precondition(klr.x >= 0 && klr.y >= 0,
00204                        "convolveImage(): coordinates of "
00205                        "kernel's lower right must be >= 0.");
00206 
00207     // use traits to determine SumType as to prevent possible overflow
00208     typedef typename
00209         PromoteTraits<typename SrcAccessor::value_type,
00210                       typename KernelAccessor::value_type>::Promote SumType;
00211     typedef typename
00212         NumericTraits<typename KernelAccessor::value_type>::RealPromote KernelSumType;
00213     typedef typename DestAccessor::value_type DestType;
00214 
00215     // calculate width and height of the image
00216     int w = src_lr.x - src_ul.x;
00217     int h = src_lr.y - src_ul.y;
00218 
00219     // calculate width and height of the kernel
00220     int kernel_width  = klr.x - kul.x + 1;
00221     int kernel_height = klr.y - kul.y + 1;
00222 
00223     vigra_precondition(w >= std::max(klr.x, -kul.x) + 1 && h >= std::max(klr.y, -kul.y) + 1,
00224                        "convolveImage(): kernel larger than image.");
00225 
00226     KernelSumType norm = NumericTraits<KernelSumType>::zero();
00227     if(border == BORDER_TREATMENT_CLIP)
00228     {
00229         // calculate the sum of the kernel elements for renormalization
00230         KernelIterator yk  = ki + klr;
00231 
00232         // determine sum within kernel (= norm)
00233         for(int y = 0; y < kernel_height; ++y, --yk.y)
00234         {
00235             KernelIterator xk  = yk;
00236             for(int x = 0; x < kernel_width; ++x, --xk.x)
00237             {
00238                 norm += ak(xk);
00239             }
00240         }
00241         vigra_precondition(norm != NumericTraits<KernelSumType>::zero(),
00242             "convolveImage(): Cannot use BORDER_TREATMENT_CLIP with a DC-free kernel");
00243     }
00244 
00245     DestIterator yd = dest_ul;
00246     SrcIterator ys = src_ul;
00247     SrcIterator send = src_lr;
00248 
00249     // iterate over the interior part
00250     for(int y=0; y<h; ++y, ++ys.y, ++yd.y)
00251     {
00252         // create x iterators
00253         DestIterator xd(yd);
00254         SrcIterator xs(ys);
00255 
00256         for(int x=0; x < w; ++x, ++xs.x, ++xd.x)
00257         {
00258             // init the sum
00259             SumType sum = NumericTraits<SumType>::zero();
00260             KernelIterator ykernel  = ki + klr;
00261             
00262             if(x >= klr.x && y >= klr.y && x < w + kul.x && y < h + kul.y)
00263             {
00264                 // kernel is entirely inside the image
00265                 SrcIterator yys = xs - klr;
00266                 SrcIterator yyend = xs - kul;
00267 
00268                 for(; yys.y <= yyend.y; ++yys.y, --ykernel.y)
00269                 {
00270                     typename SrcIterator::row_iterator xxs = yys.rowIterator();
00271                     typename SrcIterator::row_iterator xxe = xxs + kernel_width;
00272                     typename KernelIterator::row_iterator xkernel= ykernel.rowIterator();
00273 
00274                     for(; xxs < xxe; ++xxs, --xkernel)
00275                     {
00276                         sum += ak(xkernel) * src_acc(xxs);
00277                     }
00278                 }
00279             }
00280             else if(border == BORDER_TREATMENT_REPEAT)
00281             {
00282                 Diff2D diff;
00283                 for(int yk = klr.y; yk >= kul.y; --yk, --ykernel.y)
00284                 {
00285                     diff.y = std::min(std::max(y - yk, 0), h-1);
00286                     typename KernelIterator::row_iterator xkernel  = ykernel.rowIterator();
00287 
00288                     for(int xk = klr.x; xk >= kul.x; --xk, --xkernel)
00289                     {
00290                         diff.x = std::min(std::max(x - xk, 0), w-1);
00291                         sum += ak(xkernel) * src_acc(src_ul, diff);
00292                     }
00293                 }
00294             }
00295             else if(border == BORDER_TREATMENT_REFLECT)
00296             {
00297                 Diff2D diff;
00298                 for(int yk = klr.y; yk >= kul.y; --yk , --ykernel.y)
00299                 {
00300                     diff.y = abs(y - yk);
00301                     if(diff.y >= h)
00302                         diff.y = 2*h - 2 - diff.y;
00303                     typename KernelIterator::row_iterator xkernel  = ykernel.rowIterator();
00304 
00305                     for(int xk = klr.x; xk >= kul.x; --xk, --xkernel)
00306                     {
00307                         diff.x = abs(x - xk);
00308                         if(diff.x >= w)
00309                             diff.x = 2*w - 2 - diff.x;
00310                         sum += ak(xkernel) * src_acc(src_ul, diff);
00311                     }
00312                 }
00313             }
00314             else if(border == BORDER_TREATMENT_WRAP)
00315             {
00316                 Diff2D diff;
00317                 for(int yk = klr.y; yk >= kul.y; --yk, --ykernel.y)
00318                 {
00319                     diff.y = (y - yk + h) % h;
00320                     typename KernelIterator::row_iterator xkernel  = ykernel.rowIterator();
00321 
00322                     for(int xk = klr.x; xk >= kul.x; --xk, --xkernel)
00323                     {
00324                         diff.x = (x - xk + w) % w;
00325                         sum += ak(xkernel) * src_acc(src_ul, diff);
00326                     }
00327                 }
00328             }
00329             else if(border == BORDER_TREATMENT_CLIP)
00330             {
00331                 KernelSumType ksum = NumericTraits<KernelSumType>::zero();
00332                 Diff2D diff;
00333                 for(int yk = klr.y; yk >= kul.y; --yk, --ykernel.y)
00334                 {
00335                     diff.y = y - yk;
00336                     if(diff.y < 0 || diff.y >= h)
00337                         continue;
00338                     typename KernelIterator::row_iterator xkernel  = ykernel.rowIterator();
00339 
00340                     for(int xk = klr.x; xk >= kul.x; --xk, --xkernel)
00341                     {
00342                         diff.x = x - xk;
00343                         if(diff.x < 0 || diff.x >= w)
00344                             continue;
00345                         ksum += ak(xkernel);
00346                         sum += ak(xkernel) * src_acc(src_ul, diff);
00347                     }
00348                 }
00349                 
00350                 sum *= norm / ksum;
00351             }
00352             else if(border == BORDER_TREATMENT_AVOID)
00353             {
00354                 continue;
00355             }
00356 
00357             // store convolution result in destination pixel
00358             dest_acc.set(detail::RequiresExplicitCast<DestType>::cast(sum), xd);
00359         }
00360     }
00361 }
00362 
00363 template <class SrcIterator, class SrcAccessor,
00364           class DestIterator, class DestAccessor,
00365           class KernelIterator, class KernelAccessor>
00366 inline
00367 void convolveImage(
00368                    triple<SrcIterator, SrcIterator, SrcAccessor> src,
00369                    pair<DestIterator, DestAccessor> dest,
00370                    tuple5<KernelIterator, KernelAccessor, Diff2D, Diff2D,
00371                    BorderTreatmentMode> kernel)
00372 {
00373     convolveImage(src.first, src.second, src.third,
00374                   dest.first, dest.second,
00375                   kernel.first, kernel.second, kernel.third,
00376                   kernel.fourth, kernel.fifth);
00377 }
00378 
00379 
00380 /** \brief Performs a 2-dimensional normalized convolution, i.e. convolution with a mask image.
00381 
00382     This functions computes
00383     <a href ="http://homepages.inf.ed.ac.uk/rbf/CVonline/LOCAL_COPIES/PIRODDI1/NormConv/NormConv.html">normalized
00384     convolution</a> as defined in
00385     Knutsson, H. and Westin, C-F.: <i>Normalized and differential convolution:
00386     Methods for Interpolation and Filtering of incomplete and uncertain data</i>.
00387     Proc. of the IEEE Conf. on Computer Vision and Pattern Recognition, 1993, 515-523.
00388 
00389     The mask image must be binary and encodes which pixels of the original image
00390     are valid. It is used as follows:
00391     Only pixel under the mask are used in the calculations. Whenever a part of the
00392     kernel lies outside the mask, it is ignored, and the kernel is renormalized to its
00393     original norm (analogous to the CLIP \ref BorderTreatmentMode). Thus, a useful convolution
00394     result is computed whenever <i>at least one valid pixel is within the current window</i>
00395     Thus, destination pixels not under the mask still receive a value if they are <i>near</i>
00396     the mask. Therefore, this algorithm is useful as an interpolator of sparse input data.
00397     If you are only interested in the destination values under the mask, you can perform
00398     a subsequent \ref copyImageIf().
00399 
00400     The KernelIterator must point to the center of the kernel, and
00401     the kernel's size is given by its upper left (x and y of distance <= 0) and
00402     lower right (distance >= 0) corners. The image must always be larger than the
00403     kernel. At those positions where the kernel does not completely fit
00404     into the image, the specified \ref BorderTreatmentMode is
00405     applied. Only BORDER_TREATMENT_CLIP and BORDER_TREATMENT_AVOID are currently
00406     supported.
00407 
00408     The images's pixel type (SrcAccessor::value_type) must be a
00409     linear space over the kernel's value_type (KernelAccessor::value_type),
00410     i.e. addition of source values, multiplication with kernel values,
00411     and NumericTraits must be defined.
00412     The kernel's value_type must be an algebraic field,
00413     i.e. the arithmetic operations (+, -, *, /) and NumericTraits must
00414     be defined.
00415 
00416     <b> Declarations:</b>
00417 
00418     pass arguments explicitly:
00419     \code
00420     namespace vigra {
00421         template <class SrcIterator, class SrcAccessor,
00422                   class MaskIterator, class MaskAccessor,
00423                   class DestIterator, class DestAccessor,
00424                   class KernelIterator, class KernelAccessor>
00425         void
00426         normalizedConvolveImage(SrcIterator src_ul, SrcIterator src_lr, SrcAccessor src_acc,
00427                                 MaskIterator mul, MaskAccessor am,
00428                                 DestIterator dest_ul, DestAccessor dest_acc,
00429                                 KernelIterator ki, KernelAccessor ak,
00430                                 Diff2D kul, Diff2D klr, BorderTreatmentMode border);
00431     }
00432     \endcode
00433 
00434 
00435     use argument objects in conjunction with \ref ArgumentObjectFactories :
00436     \code
00437     namespace vigra {
00438         template <class SrcIterator, class SrcAccessor,
00439                   class MaskIterator, class MaskAccessor,
00440                   class DestIterator, class DestAccessor,
00441                   class KernelIterator, class KernelAccessor>
00442         void normalizedConvolveImage(triple<SrcIterator, SrcIterator, SrcAccessor> src,
00443                                      pair<MaskIterator, MaskAccessor> mask,
00444                                      pair<DestIterator, DestAccessor> dest,
00445                                      tuple5<KernelIterator, KernelAccessor, Diff2D, Diff2D,
00446                                      BorderTreatmentMode> kernel);
00447     }
00448     \endcode
00449 
00450     <b> Usage:</b>
00451 
00452     <b>\#include</b> <vigra/stdconvolution.hxx><br>
00453     Namespace: vigra
00454 
00455 
00456     \code
00457     vigra::FImage src(w,h), dest(w,h);
00458     vigra::CImage mask(w,h);
00459     ...
00460 
00461     // define 3x3 binomial filter
00462     vigra::Kernel2D<float> binom;
00463 
00464     binom.initExplicitly(Diff2D(-1,-1), Diff2D(1,1)) =   // upper left and lower right
00465                          0.0625, 0.125, 0.0625,
00466                          0.125,  0.25,  0.125,
00467                          0.0625, 0.125, 0.0625;
00468 
00469     vigra::normalizedConvolveImage(srcImageRange(src), maskImage(mask), destImage(dest), kernel2d(binom));
00470     \endcode
00471 
00472     <b> Required Interface:</b>
00473 
00474     \code
00475     ImageIterator src_ul, src_lr;
00476     ImageIterator mul;
00477     ImageIterator dest_ul;
00478     ImageIterator ik;
00479 
00480     SrcAccessor src_accessor;
00481     MaskAccessor mask_accessor;
00482     DestAccessor dest_accessor;
00483     KernelAccessor kernel_accessor;
00484 
00485     NumericTraits<SrcAccessor::value_type>::RealPromote s = src_accessor(src_ul);
00486 
00487     s = s + s;
00488     s = kernel_accessor(ik) * s;
00489     s -= s;
00490 
00491     if(mask_accessor(mul)) ...;
00492 
00493     dest_accessor.set(
00494     NumericTraits<DestAccessor::value_type>::fromRealPromote(s), dest_ul);
00495 
00496     NumericTraits<KernelAccessor::value_type>::RealPromote k = kernel_accessor(ik);
00497 
00498     k += k;
00499     k -= k;
00500     k = k / k;
00501 
00502     \endcode
00503 
00504     <b> Preconditions:</b>
00505 
00506     \code
00507     kul.x <= 0
00508     kul.y <= 0
00509     klr.x >= 0
00510     klr.y >= 0
00511     src_lr.x - src_ul.x >= klr.x + kul.x + 1
00512     src_lr.y - src_ul.y >= klr.y + kul.y + 1
00513     border == BORDER_TREATMENT_CLIP || border == BORDER_TREATMENT_AVOID
00514     \endcode
00515 
00516     Sum of kernel elements must be != 0.
00517 
00518 */
00519 doxygen_overloaded_function(template <...> void normalizedConvolveImage)
00520 
00521 template <class SrcIterator, class SrcAccessor,
00522           class DestIterator, class DestAccessor,
00523           class MaskIterator, class MaskAccessor,
00524           class KernelIterator, class KernelAccessor>
00525 void
00526 normalizedConvolveImage(SrcIterator src_ul, SrcIterator src_lr, SrcAccessor src_acc,
00527                         MaskIterator mul, MaskAccessor am,
00528                         DestIterator dest_ul, DestAccessor dest_acc,
00529                         KernelIterator ki, KernelAccessor ak,
00530                         Diff2D kul, Diff2D klr, BorderTreatmentMode border)
00531 {
00532     vigra_precondition((border == BORDER_TREATMENT_CLIP  ||
00533                         border == BORDER_TREATMENT_AVOID),
00534                        "normalizedConvolveImage(): "
00535                        "Border treatment must be BORDER_TREATMENT_CLIP or BORDER_TREATMENT_AVOID.");
00536 
00537     vigra_precondition(kul.x <= 0 && kul.y <= 0,
00538                        "normalizedConvolveImage(): left borders must be <= 0.");
00539     vigra_precondition(klr.x >= 0 && klr.y >= 0,
00540                        "normalizedConvolveImage(): right borders must be >= 0.");
00541 
00542     // use traits to determine SumType as to prevent possible overflow
00543     typedef typename
00544         NumericTraits<typename SrcAccessor::value_type>::RealPromote SumType;
00545     typedef typename
00546         NumericTraits<typename KernelAccessor::value_type>::RealPromote KSumType;
00547     typedef
00548         NumericTraits<typename DestAccessor::value_type> DestTraits;
00549 
00550     // calculate width and height of the image
00551     int w = src_lr.x - src_ul.x;
00552     int h = src_lr.y - src_ul.y;
00553     int kernel_width = klr.x - kul.x + 1;
00554     int kernel_height = klr.y - kul.y + 1;
00555 
00556     int x,y;
00557     int ystart = (border == BORDER_TREATMENT_AVOID) ?  klr.y : 0;
00558     int yend   = (border == BORDER_TREATMENT_AVOID) ? h+kul.y : h;
00559     int xstart = (border == BORDER_TREATMENT_AVOID) ?  klr.x : 0;
00560     int xend   = (border == BORDER_TREATMENT_AVOID) ? w+kul.x : w;
00561 
00562     // create y iterators
00563     DestIterator yd = dest_ul + Diff2D(xstart, ystart);
00564     SrcIterator ys = src_ul + Diff2D(xstart, ystart);
00565     MaskIterator ym = mul + Diff2D(xstart, ystart);
00566 
00567     KSumType norm = ak(ki);
00568     int xx, yy;
00569     KernelIterator yk  = ki + klr;
00570     for(yy=0; yy<kernel_height; ++yy, --yk.y)
00571     {
00572         KernelIterator xk  = yk;
00573 
00574         for(xx=0; xx<kernel_width; ++xx, --xk.x)
00575         {
00576             norm += ak(xk);
00577         }
00578     }
00579     norm -= ak(ki);
00580 
00581 
00582     for(y=ystart; y < yend; ++y, ++ys.y, ++yd.y, ++ym.y)
00583     {
00584         // create x iterators
00585         DestIterator xd(yd);
00586         SrcIterator xs(ys);
00587         MaskIterator xm(ym);
00588 
00589         for(x=xstart; x < xend; ++x, ++xs.x, ++xd.x, ++xm.x)
00590         {
00591             // how much of the kernel fits into the image ?
00592             int x0, y0, x1, y1;
00593 
00594             y0 = (y<klr.y) ? -y : -klr.y;
00595             y1 = (h-y-1<-kul.y) ? h-y-1 : -kul.y;
00596             x0 = (x<klr.x) ? -x : -klr.x;
00597             x1 = (w-x-1<-kul.x) ? w-x-1 : -kul.x;
00598 
00599             bool first = true;
00600             // init the sum
00601             SumType sum = NumericTraits<SumType>::zero();
00602             KSumType ksum = NumericTraits<KSumType>::zero();
00603 
00604             SrcIterator yys = xs + Diff2D(x0, y0);
00605             MaskIterator yym = xm + Diff2D(x0, y0);
00606             KernelIterator yk  = ki - Diff2D(x0, y0);
00607 
00608             int xx, kernel_width, kernel_height;
00609             kernel_width = x1 - x0 + 1;
00610             kernel_height = y1 - y0 + 1;
00611             for(yy=0; yy<kernel_height; ++yy, ++yys.y, --yk.y, ++yym.y)
00612             {
00613                 typename SrcIterator::row_iterator xxs = yys.rowIterator();
00614                 typename SrcIterator::row_iterator xxend = xxs + kernel_width;
00615                 typename MaskIterator::row_iterator xxm = yym.rowIterator();
00616                 typename KernelIterator::row_iterator xk  = yk.rowIterator();
00617 
00618                 for(xx=0; xxs < xxend; ++xxs, --xk, ++xxm)
00619                 {
00620                     if(!am(xxm)) continue;
00621 
00622                     if(first)
00623                     {
00624                         sum = detail::RequiresExplicitCast<SumType>::cast(ak(xk) * src_acc(xxs));
00625                         ksum = ak(xk);
00626                         first = false;
00627                     }
00628                     else
00629                     {
00630                         sum = detail::RequiresExplicitCast<SumType>::cast(sum + ak(xk) * src_acc(xxs));
00631                         ksum += ak(xk);
00632                     }
00633                 }
00634             }
00635             // store average in destination pixel
00636             if(ksum != NumericTraits<KSumType>::zero())
00637             {
00638                 dest_acc.set(DestTraits::fromRealPromote(
00639                              detail::RequiresExplicitCast<SumType>::cast((norm / ksum) * sum)), xd);
00640             }
00641         }
00642     }
00643 }
00644 
00645 
00646 template <class SrcIterator, class SrcAccessor,
00647           class DestIterator, class DestAccessor,
00648           class MaskIterator, class MaskAccessor,
00649           class KernelIterator, class KernelAccessor>
00650 inline
00651 void normalizedConvolveImage(
00652                            triple<SrcIterator, SrcIterator, SrcAccessor> src,
00653                            pair<MaskIterator, MaskAccessor> mask,
00654                            pair<DestIterator, DestAccessor> dest,
00655                            tuple5<KernelIterator, KernelAccessor, Diff2D, Diff2D,
00656                            BorderTreatmentMode> kernel)
00657 {
00658     normalizedConvolveImage(src.first, src.second, src.third,
00659                             mask.first, mask.second,
00660                             dest.first, dest.second,
00661                             kernel.first, kernel.second, kernel.third,
00662                             kernel.fourth, kernel.fifth);
00663 }
00664 
00665 /** \brief Deprecated name of 2-dimensional normalized convolution, i.e. convolution with a mask image.
00666 
00667     See \ref normalizedConvolveImage() for documentation.
00668 
00669     <b> Declarations:</b>
00670 
00671     pass arguments explicitly:
00672     \code
00673     namespace vigra {
00674         template <class SrcIterator, class SrcAccessor,
00675                   class MaskIterator, class MaskAccessor,
00676                   class DestIterator, class DestAccessor,
00677                   class KernelIterator, class KernelAccessor>
00678         void
00679         convolveImageWithMask(SrcIterator src_ul, SrcIterator src_lr, SrcAccessor src_acc,
00680                               MaskIterator mul, MaskAccessor am,
00681                               DestIterator dest_ul, DestAccessor dest_acc,
00682                               KernelIterator ki, KernelAccessor ak,
00683                               Diff2D kul, Diff2D klr, BorderTreatmentMode border);
00684     }
00685     \endcode
00686 
00687 
00688     use argument objects in conjunction with \ref ArgumentObjectFactories :
00689     \code
00690     namespace vigra {
00691         template <class SrcIterator, class SrcAccessor,
00692                   class MaskIterator, class MaskAccessor,
00693                   class DestIterator, class DestAccessor,
00694                   class KernelIterator, class KernelAccessor>
00695         void convolveImageWithMask(triple<SrcIterator, SrcIterator, SrcAccessor> src,
00696                                    pair<MaskIterator, MaskAccessor> mask,
00697                                    pair<DestIterator, DestAccessor> dest,
00698                                    tuple5<KernelIterator, KernelAccessor, Diff2D, Diff2D,
00699                                    BorderTreatmentMode> kernel);
00700     }
00701     \endcode
00702 */
00703 doxygen_overloaded_function(template <...> void convolveImageWithMask)
00704 
00705 template <class SrcIterator, class SrcAccessor,
00706           class DestIterator, class DestAccessor,
00707           class MaskIterator, class MaskAccessor,
00708           class KernelIterator, class KernelAccessor>
00709 inline void
00710 convolveImageWithMask(SrcIterator src_ul, SrcIterator src_lr, SrcAccessor src_acc,
00711                       MaskIterator mul, MaskAccessor am,
00712                       DestIterator dest_ul, DestAccessor dest_acc,
00713                       KernelIterator ki, KernelAccessor ak,
00714                       Diff2D kul, Diff2D klr, BorderTreatmentMode border)
00715 {
00716     normalizedConvolveImage(src_ul, src_lr, src_acc,
00717                             mul, am,
00718                             dest_ul, dest_acc,
00719                             ki, ak, kul, klr, border);
00720 }
00721 
00722 template <class SrcIterator, class SrcAccessor,
00723           class DestIterator, class DestAccessor,
00724           class MaskIterator, class MaskAccessor,
00725           class KernelIterator, class KernelAccessor>
00726 inline
00727 void convolveImageWithMask(
00728                            triple<SrcIterator, SrcIterator, SrcAccessor> src,
00729                            pair<MaskIterator, MaskAccessor> mask,
00730                            pair<DestIterator, DestAccessor> dest,
00731                            tuple5<KernelIterator, KernelAccessor, Diff2D, Diff2D,
00732                            BorderTreatmentMode> kernel)
00733 {
00734     normalizedConvolveImage(src.first, src.second, src.third,
00735                             mask.first, mask.second,
00736                             dest.first, dest.second,
00737                             kernel.first, kernel.second, kernel.third,
00738                             kernel.fourth, kernel.fifth);
00739 }
00740 
00741 //@}
00742 
00743 /********************************************************/
00744 /*                                                      */
00745 /*                      Kernel2D                        */
00746 /*                                                      */
00747 /********************************************************/
00748 
00749 /** \brief Generic 2 dimensional convolution kernel.
00750 
00751     This kernel may be used for convolution of 2 dimensional signals.
00752 
00753     Convolution functions access the kernel via an ImageIterator
00754     which they get by calling \ref center(). This iterator
00755     points to the center of the kernel. The kernel's size is given by its upperLeft()
00756     (upperLeft().x <= 0, upperLeft().y <= 0)
00757     and lowerRight() (lowerRight().x >= 0, lowerRight().y >= 0) methods.
00758     The desired border treatment mode is returned by borderTreatment().
00759     (Note that the \ref StandardConvolution "2D convolution functions" don't currently
00760     support all modes.)
00761 
00762     The different init functions create a kernel with the specified
00763     properties. The requirements for the kernel's value_type depend
00764     on the init function used. At least NumericTraits must be defined.
00765 
00766     The kernel defines a factory function kernel2d() to create an argument object
00767     (see \ref KernelArgumentObjectFactories).
00768 
00769     <b> Usage:</b>
00770 
00771     <b>\#include</b> <vigra/stdconvolution.hxx><br>
00772     Namespace: vigra
00773 
00774     \code
00775     vigra::FImage src(w,h), dest(w,h);
00776     ...
00777 
00778     // define horizontal Sobel filter
00779     vigra::Kernel2D<float> sobel;
00780 
00781     sobel.initExplicitly(Diff2D(-1,-1), Diff2D(1,1)) =  // upper left and lower right
00782                          0.125, 0.0, -0.125,
00783                          0.25,  0.0, -0.25,
00784                          0.125, 0.0, -0.125;
00785 
00786     vigra::convolveImage(srcImageRange(src), destImage(dest), kernel2d(sobel));
00787     \endcode
00788 
00789     <b> Required Interface:</b>
00790 
00791     \code
00792     value_type v = NumericTraits<value_type>::one();
00793     \endcode
00794 
00795     See also the init functions.
00796 */
00797 template <class ARITHTYPE>
00798 class Kernel2D
00799 {
00800 public:
00801         /** the kernel's value type
00802          */
00803     typedef ARITHTYPE value_type;
00804 
00805         /** 2D random access iterator over the kernel's values
00806          */
00807     typedef typename BasicImage<value_type>::traverser Iterator;
00808 
00809         /** const 2D random access iterator over the kernel's values
00810          */
00811     typedef typename BasicImage<value_type>::const_traverser ConstIterator;
00812 
00813         /** the kernel's accessor
00814          */
00815     typedef typename BasicImage<value_type>::Accessor Accessor;
00816 
00817         /** the kernel's const accessor
00818          */
00819     typedef typename BasicImage<value_type>::ConstAccessor ConstAccessor;
00820 
00821     struct InitProxy
00822     {
00823         typedef typename
00824         BasicImage<value_type>::ScanOrderIterator Iterator;
00825 
00826         InitProxy(Iterator i, int count, value_type & norm)
00827             : iter_(i), base_(i),
00828               count_(count), sum_(count),
00829               norm_(norm)
00830         {}
00831 
00832         ~InitProxy()
00833         {
00834             vigra_precondition(count_ == 1 || count_ == sum_,
00835                                "Kernel2D::initExplicitly(): "
00836                                "Too few init values.");
00837         }
00838 
00839         InitProxy & operator,(value_type const & v)
00840         {
00841             if(count_ == sum_)  norm_ = *iter_;
00842 
00843             --count_;
00844             vigra_precondition(count_ > 0,
00845                                "Kernel2D::initExplicitly(): "
00846                                "Too many init values.");
00847 
00848             norm_ += v;
00849 
00850             ++iter_;
00851             *iter_ = v;
00852 
00853             return *this;
00854         }
00855 
00856         Iterator iter_, base_;
00857         int count_, sum_;
00858         value_type & norm_;
00859     };
00860 
00861     static value_type one() { return NumericTraits<value_type>::one(); }
00862 
00863         /** Default constructor.
00864             Creates a kernel of size 1x1 which would copy the signal
00865             unchanged.
00866         */
00867     Kernel2D()
00868         : kernel_(1, 1, one()),
00869           left_(0, 0),
00870           right_(0, 0),
00871           norm_(one()),
00872           border_treatment_(BORDER_TREATMENT_REFLECT)
00873     {}
00874 
00875         /** Copy constructor.
00876          */
00877     Kernel2D(Kernel2D const & k)
00878         : kernel_(k.kernel_),
00879           left_(k.left_),
00880           right_(k.right_),
00881           norm_(k.norm_),
00882           border_treatment_(k.border_treatment_)
00883     {}
00884 
00885         /** Copy assignment.
00886          */
00887     Kernel2D & operator=(Kernel2D const & k)
00888     {
00889         if(this != &k)
00890         {
00891         kernel_ = k.kernel_;
00892             left_ = k.left_;
00893             right_ = k.right_;
00894             norm_ = k.norm_;
00895         border_treatment_ = k.border_treatment_;
00896         }
00897         return *this;
00898     }
00899 
00900         /** Initialization.
00901             This initializes the kernel with the given constant. The norm becomes
00902             v*width()*height().
00903 
00904             Instead of a single value an initializer list of length width()*height()
00905             can be used like this:
00906 
00907             \code
00908             vigra::Kernel2D<float> binom;
00909 
00910             binom.initExplicitly(Diff2D(-1,-1), Diff2D(1,1)) =
00911             0.0625, 0.125, 0.0625,
00912             0.125,  0.25,  0.125,
00913             0.0625, 0.125, 0.0625;
00914             \endcode
00915 
00916             In this case, the norm will be set to the sum of the init values.
00917             An initializer list of wrong length will result in a run-time error.
00918         */
00919     InitProxy operator=(value_type const & v)
00920     {
00921         int size = (right_.x - left_.x + 1) *
00922                    (right_.y - left_.y + 1);
00923         kernel_ = v;
00924         norm_ = (double)size*v;
00925 
00926         return InitProxy(kernel_.begin(), size, norm_);
00927     }
00928 
00929         /** Destructor.
00930          */
00931     ~Kernel2D()
00932     {}
00933 
00934         /** Init the 2D kernel as the cartesian product of two 1D kernels
00935             of type \ref Kernel1D. The norm becomes the product of the two original
00936             norms.
00937 
00938             <b> Required Interface:</b>
00939 
00940             The kernel's value_type must be a linear algebra.
00941 
00942             \code
00943             vigra::Kernel2D<...>::value_type v;
00944             v = v * v;
00945             \endcode
00946         */
00947     void initSeparable(Kernel1D<value_type> const & kx,
00948                        Kernel1D<value_type> const & ky)
00949     {
00950         left_ = Diff2D(kx.left(), ky.left());
00951         right_ = Diff2D(kx.right(), ky.right());
00952         int w = right_.x - left_.x + 1;
00953         int h = right_.y - left_.y + 1;
00954         kernel_.resize(w, h);
00955 
00956         norm_ = kx.norm() * ky.norm();
00957 
00958         typedef typename Kernel1D<value_type>::const_iterator KIter;
00959         typename Kernel1D<value_type>::Accessor ka;
00960 
00961         KIter kiy = ky.center() + left_.y;
00962         Iterator iy = center() + left_;
00963 
00964         for(int y=left_.y; y<=right_.y; ++y, ++kiy, ++iy.y)
00965         {
00966             KIter kix = kx.center() + left_.x;
00967             Iterator ix = iy;
00968             for(int x=left_.x; x<=right_.x; ++x, ++kix, ++ix.x)
00969             {
00970                 *ix = ka(kix) * ka(kiy);
00971             }
00972         }
00973     }
00974 
00975         /** Init the 2D kernel as the cartesian product of two 1D kernels
00976             given explicitly by iterators and sizes. The norm becomes the
00977             sum of the resulting kernel values.
00978 
00979             <b> Required Interface:</b>
00980 
00981             The kernel's value_type must be a linear algebra.
00982 
00983             \code
00984             vigra::Kernel2D<...>::value_type v;
00985             v = v * v;
00986             v += v;
00987             \endcode
00988 
00989             <b> Preconditions:</b>
00990 
00991             \code
00992             xleft <= 0;
00993             xright >= 0;
00994             yleft <= 0;
00995             yright >= 0;
00996             \endcode
00997         */
00998     template <class KernelIterator>
00999     void initSeparable(KernelIterator kxcenter, int xleft, int xright,
01000                        KernelIterator kycenter, int yleft, int yright)
01001     {
01002         vigra_precondition(xleft <= 0 && yleft <= 0,
01003                            "Kernel2D::initSeparable(): left borders must be <= 0.");
01004         vigra_precondition(xright >= 0 && yright >= 0,
01005                            "Kernel2D::initSeparable(): right borders must be >= 0.");
01006 
01007         left_ = Point2D(xleft, yleft);
01008         right_ = Point2D(xright, yright);
01009 
01010         int w = right_.x - left_.x + 1;
01011         int h = right_.y - left_.y + 1;
01012         kernel_.resize(w, h);
01013 
01014         KernelIterator kiy = kycenter + left_.y;
01015         Iterator iy = center() + left_;
01016 
01017         for(int y=left_.y; y<=right_.y; ++y, ++kiy, ++iy.y)
01018         {
01019             KernelIterator kix = kxcenter + left_.x;
01020             Iterator ix = iy;
01021             for(int x=left_.x; x<=right_.x; ++x, ++kix, ++ix.x)
01022             {
01023                 *ix = *kix * *kiy;
01024             }
01025         }
01026 
01027         typename BasicImage<value_type>::iterator i = kernel_.begin();
01028         typename BasicImage<value_type>::iterator iend = kernel_.end();
01029         norm_ = *i;
01030         ++i;
01031 
01032         for(; i!= iend; ++i)
01033         {
01034             norm_ += *i;
01035         }
01036     }
01037 
01038         /** Init as a 2D Gaussian function with given standard deviation and norm.
01039          */    
01040     void initGaussian(double std_dev, value_type norm)
01041     {
01042         Kernel1D<value_type> gauss;
01043         gauss.initGaussian(std_dev, norm);
01044         initSeparable(gauss, gauss);
01045     }
01046 
01047         /** Init as a 2D Gaussian function with given standard deviation and unit norm.
01048          */
01049     void initGaussian(double std_dev)
01050     {
01051         initGaussian(std_dev, NumericTraits<value_type>::one());
01052     }
01053 
01054         /** Init the 2D kernel as a circular averaging filter. The norm will be
01055             calculated as
01056             <TT>NumericTraits<value_type>::one() / (number of non-zero kernel values)</TT>.
01057             The kernel's value_type must be a linear space.
01058 
01059             <b> Required Interface:</b>
01060 
01061             \code
01062             value_type v = vigra::NumericTraits<value_type>::one();
01063 
01064             double d;
01065             v = d * v;
01066             \endcode
01067 
01068             <b> Precondition:</b>
01069 
01070             \code
01071             radius > 0;
01072             \endcode
01073         */
01074     void initDisk(int radius)
01075     {
01076         vigra_precondition(radius > 0,
01077                            "Kernel2D::initDisk(): radius must be > 0.");
01078 
01079         left_ = Point2D(-radius, -radius);
01080         right_ = Point2D(radius, radius);
01081         int w = right_.x - left_.x + 1;
01082         int h = right_.y - left_.y + 1;
01083         kernel_.resize(w, h);
01084         norm_ = NumericTraits<value_type>::one();
01085 
01086         kernel_ = NumericTraits<value_type>::zero();
01087         double count = 0.0;
01088 
01089         Iterator k = center();
01090         double r2 = (double)radius*radius;
01091 
01092         int i;
01093         for(i=0; i<= radius; ++i)
01094         {
01095             double r = (double) i - 0.5;
01096             int w = (int)(VIGRA_CSTD::sqrt(r2 - r*r) + 0.5);
01097             for(int j=-w; j<=w; ++j)
01098             {
01099                 k(j, i) = NumericTraits<value_type>::one();
01100                 k(j, -i) = NumericTraits<value_type>::one();
01101                 count += (i != 0) ? 2.0 : 1.0;
01102             }
01103         }
01104 
01105         count = 1.0 / count;
01106 
01107         for(int y=-radius; y<=radius; ++y)
01108         {
01109             for(int x=-radius; x<=radius; ++x)
01110             {
01111                 k(x,y) = count * k(x,y);
01112             }
01113         }
01114     }
01115 
01116         /** Init the kernel by an explicit initializer list.
01117             The upper left and lower right corners of the kernel must be passed.
01118             A comma-separated initializer list is given after the assignment operator.
01119             This function is used like this:
01120 
01121             \code
01122             // define horizontal Sobel filter
01123             vigra::Kernel2D<float> sobel;
01124 
01125             sobel.initExplicitly(Diff2D(-1,-1), Diff2D(1,1)) =
01126             0.125, 0.0, -0.125,
01127             0.25,  0.0, -0.25,
01128             0.125, 0.0, -0.125;
01129             \endcode
01130 
01131             The norm is set to the sum of the initializer values. If the wrong number of
01132             values is given, a run-time error results. It is, however, possible to give
01133             just one initializer. This creates an averaging filter with the given constant:
01134 
01135             \code
01136             vigra::Kernel2D<float> average3x3;
01137 
01138             average3x3.initExplicitly(Diff2D(-1,-1), Diff2D(1,1)) = 1.0/9.0;
01139             \endcode
01140 
01141             Here, the norm is set to value*width()*height().
01142 
01143             <b> Preconditions:</b>
01144 
01145             \code
01146             1. upperleft.x <= 0;
01147             2. upperleft.y <= 0;
01148             3. lowerright.x >= 0;
01149             4. lowerright.y >= 0;
01150             5. the number of values in the initializer list
01151             is 1 or equals the size of the kernel.
01152             \endcode
01153         */
01154     Kernel2D & initExplicitly(Diff2D upperleft, Diff2D lowerright)
01155     {
01156         vigra_precondition(upperleft.x <= 0 && upperleft.y <= 0,
01157                            "Kernel2D::initExplicitly(): left borders must be <= 0.");
01158         vigra_precondition(lowerright.x >= 0 && lowerright.y >= 0,
01159                            "Kernel2D::initExplicitly(): right borders must be >= 0.");
01160 
01161         left_ = Point2D(upperleft);
01162         right_ = Point2D(lowerright);
01163 
01164         int w = right_.x - left_.x + 1;
01165         int h = right_.y - left_.y + 1;
01166         kernel_.resize(w, h);
01167 
01168         return *this;
01169     }
01170 
01171         /** Coordinates of the upper left corner of the kernel.
01172          */
01173     Point2D upperLeft() const { return left_; }
01174 
01175         /** Coordinates of the lower right corner of the kernel.
01176          */
01177     Point2D lowerRight() const { return right_; }
01178 
01179         /** Width of the kernel.
01180          */
01181     int width() const { return right_.x - left_.x + 1; }
01182 
01183         /** Height of the kernel.
01184          */
01185     int height() const { return right_.y - left_.y + 1; }
01186 
01187         /** ImageIterator that points to the center of the kernel (coordinate (0,0)).
01188          */
01189     Iterator center() { return kernel_.upperLeft() - left_; }
01190 
01191         /** ImageIterator that points to the center of the kernel (coordinate (0,0)).
01192          */
01193     ConstIterator center() const { return kernel_.upperLeft() - left_; }
01194 
01195         /** Access kernel entry at given position.
01196          */
01197     value_type & operator()(int x, int y)
01198     { return kernel_[Diff2D(x,y) - left_]; }
01199 
01200         /** Read kernel entry at given position.
01201          */
01202     value_type operator()(int x, int y) const
01203     { return kernel_[Diff2D(x,y) - left_]; }
01204 
01205         /** Access kernel entry at given position.
01206          */
01207     value_type & operator[](Diff2D const & d)
01208     { return kernel_[d - left_]; }
01209 
01210         /** Read kernel entry at given position.
01211          */
01212     value_type operator[](Diff2D const & d) const
01213     { return kernel_[d - left_]; }
01214 
01215         /** Norm of the kernel (i.e. sum of its elements).
01216          */
01217     value_type norm() const { return norm_; }
01218 
01219         /** The kernels default accessor.
01220          */
01221     Accessor accessor() { return Accessor(); }
01222 
01223         /** The kernels default const accessor.
01224          */
01225     ConstAccessor accessor() const { return ConstAccessor(); }
01226 
01227         /** Normalize the kernel to the given value. (The norm is the sum of all kernel
01228             elements.) The kernel's value_type must be a division algebra or
01229             algebraic field.
01230 
01231             <b> Required Interface:</b>
01232 
01233             \code
01234             value_type v = vigra::NumericTraits<value_type>::one(); // if norm is not
01235                                                                     // given explicitly
01236 
01237             v += v;
01238             v = v * v;
01239             v = v / v;
01240             \endcode
01241         */
01242     void normalize(value_type norm)
01243     {
01244         typename BasicImage<value_type>::iterator i = kernel_.begin();
01245         typename BasicImage<value_type>::iterator iend = kernel_.end();
01246         typename NumericTraits<value_type>::RealPromote sum = *i;
01247         ++i;
01248 
01249         for(; i!= iend; ++i)
01250         {
01251             sum += *i;
01252         }
01253 
01254         sum = norm / sum;
01255         i = kernel_.begin();
01256         for(; i != iend; ++i)
01257         {
01258             *i = *i * sum;
01259         }
01260 
01261         norm_ = norm;
01262     }
01263 
01264         /** Normalize the kernel to norm 1.
01265          */
01266     void normalize()
01267     {
01268         normalize(one());
01269     }
01270 
01271         /** current border treatment mode
01272          */
01273     BorderTreatmentMode borderTreatment() const
01274     { return border_treatment_; }
01275 
01276         /** Set border treatment mode.
01277             Only <TT>BORDER_TREATMENT_CLIP</TT> and <TT>BORDER_TREATMENT_AVOID</TT> are currently
01278             allowed.
01279         */
01280     void setBorderTreatment( BorderTreatmentMode new_mode)
01281     {
01282         vigra_precondition((new_mode == BORDER_TREATMENT_CLIP    ||
01283                             new_mode == BORDER_TREATMENT_AVOID   ||
01284                             new_mode == BORDER_TREATMENT_REFLECT ||
01285                             new_mode == BORDER_TREATMENT_REPEAT  ||
01286                             new_mode == BORDER_TREATMENT_WRAP),
01287                            "convolveImage():\n"
01288                            "  Border treatment must be one of follow treatments:\n"
01289                            "  - BORDER_TREATMENT_CLIP\n"
01290                            "  - BORDER_TREATMENT_AVOID\n"
01291                            "  - BORDER_TREATMENT_REFLECT\n"
01292                            "  - BORDER_TREATMENT_REPEAT\n"
01293                            "  - BORDER_TREATMENT_WRAP\n");
01294 
01295         border_treatment_ = new_mode;
01296     }
01297 
01298 
01299 private:
01300     BasicImage<value_type> kernel_;
01301     Point2D left_, right_;
01302     value_type norm_;
01303     BorderTreatmentMode border_treatment_;
01304 };
01305 
01306 /**************************************************************/
01307 /*                                                            */
01308 /*         Argument object factories for Kernel2D             */
01309 /*                                                            */
01310 /*     (documentation: see vigra/convolution.hxx)             */
01311 /*                                                            */
01312 /**************************************************************/
01313 
01314 template <class KernelIterator, class KernelAccessor>
01315 inline
01316 tuple5<KernelIterator, KernelAccessor, Diff2D, Diff2D, BorderTreatmentMode>
01317 kernel2d(KernelIterator ik, KernelAccessor ak, Diff2D kul, Diff2D klr,
01318          BorderTreatmentMode border)
01319 
01320 {
01321     return
01322         tuple5<KernelIterator, KernelAccessor, Diff2D, Diff2D, BorderTreatmentMode> (
01323                                                              ik, ak, kul, klr, border);
01324 }
01325 
01326 template <class T>
01327 inline
01328 tuple5<typename Kernel2D<T>::ConstIterator,
01329        typename Kernel2D<T>::ConstAccessor,
01330        Diff2D, Diff2D, BorderTreatmentMode>
01331 kernel2d(Kernel2D<T> const & k)
01332 
01333 {
01334     return
01335         tuple5<typename Kernel2D<T>::ConstIterator,
01336                typename Kernel2D<T>::ConstAccessor,
01337                Diff2D, Diff2D, BorderTreatmentMode>(
01338             k.center(),
01339             k.accessor(),
01340             k.upperLeft(), k.lowerRight(),
01341             k.borderTreatment());
01342 }
01343 
01344 template <class T>
01345 inline
01346 tuple5<typename Kernel2D<T>::ConstIterator,
01347        typename Kernel2D<T>::ConstAccessor,
01348        Diff2D, Diff2D, BorderTreatmentMode>
01349 kernel2d(Kernel2D<T> const & k, BorderTreatmentMode border)
01350 
01351 {
01352     return
01353         tuple5<typename Kernel2D<T>::ConstIterator,
01354                typename Kernel2D<T>::ConstAccessor,
01355                Diff2D, Diff2D, BorderTreatmentMode>(
01356             k.center(),
01357             k.accessor(),
01358             k.upperLeft(), k.lowerRight(),
01359             border);
01360 }
01361 
01362 
01363 } // namespace vigra
01364 
01365 #endif // VIGRA_STDCONVOLUTION_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)