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