vigra - vigra/flatmorphology.hxx Source File

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

vigra/flatmorphology.hxx |

    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_FLATMORPHOLOGY_HXX
    00038 #define VIGRA_FLATMORPHOLOGY_HXX
    00039
    00040 #include <cmath>
    00041 #include <vector>
    00042 #include "utilities.hxx"
    00043
    00044 namespace vigra {
    00045
    00046 /** \addtogroup Morphology Basic Morphological Operations
    00047     Perform erosion, dilation, and median with disc structuring functions
    00048
    00049     See also: \ref MultiArrayMorphology Separable morphology with parabola structuring functions in arbitrary dimensions
    00050 */
    00051 //@{
    00052
    00053 /********************************************************/
    00054 /*                                                      */
    00055 /*                  discRankOrderFilter                 */
    00056 /*                                                      */
    00057 /********************************************************/
    00058
    00059 /** \brief Apply rank order filter with disc structuring function to the image.
    00060
    00061     The pixel values of the source image <b> must</b> be in the range
    00062     0...255. Radius must be >= 0. Rank must be in the range 0.0 <= rank
    00063     <= 1.0. The filter acts as a minimum filter if rank = 0.0,
    00064     as a median if rank = 0.5, and as a maximum filter if rank = 1.0.
    00065     Accessor are used to access the pixel data.
    00066
    00067     <b> Declarations:</b>
    00068
    00069     pass arguments explicitly:
    00070     \code
    00071     namespace vigra {
    00072         template <class SrcIterator, class SrcAccessor,
    00073                   class DestIterator, class DestAccessor>
    00074         void
    00075         discRankOrderFilter(SrcIterator upperleft1,
    00076                             SrcIterator lowerright1, SrcAccessor sa,
    00077                             DestIterator upperleft2, DestAccessor da,
    00078                             int radius, float rank)
    00079     }
    00080     \endcode
    00081
    00082
    00083     use argument objects in conjunction with \ref ArgumentObjectFactories :
    00084     \code
    00085     namespace vigra {
    00086         template <class SrcIterator, class SrcAccessor,
    00087                   class DestIterator, class DestAccessor>
    00088         void
    00089         discRankOrderFilter(triple<SrcIterator, SrcIterator, SrcAccessor> src,
    00090                             pair<DestIterator, DestAccessor> dest,
    00091                             int radius, float rank)
    00092     }
    00093     \endcode
    00094
    00095     <b> Usage:</b>
    00096
    00097         <b>\#include</b> <vigra/flatmorphology.hxx><br>
    00098     Namespace: vigra
    00099
    00100     \code
    00101     vigra::CImage src, dest;
    00102
    00103     // do median filtering
    00104     vigra::discRankOrderFilter(srcImageRange(src), destImage(dest), 10, 0.5);
    00105     \endcode
    00106
    00107     <b> Required Interface:</b>
    00108
    00109     \code
    00110     SrcIterator src_upperleft;
    00111     DestIterator dest_upperleft;
    00112     int x, y;
    00113     unsigned char value;
    00114
    00115     SrcAccessor src_accessor;
    00116     DestAccessor dest_accessor;
    00117
    00118     // value_type of accessor must be convertible to unsigned char
    00119     value = src_accessor(src_upperleft, x, y);
    00120
    00121     dest_accessor.set(value, dest_upperleft, x, y);
    00122     \endcode
    00123
    00124     <b> Preconditions:</b>
    00125
    00126     \code
    00127     for all source pixels: 0 <= value <= 255
    00128
    00129     (rank >= 0.0) && (rank <= 1.0)
    00130     radius >= 0
    00131     \endcode
    00132
    00133 */
    00134 doxygen_overloaded_function(template <...> void discRankOrderFilter)
    00135
    00136 template <class SrcIterator, class SrcAccessor,
    00137           class DestIterator, class DestAccessor>
    00138 void
    00139 discRankOrderFilter(SrcIterator upperleft1,
    00140                     SrcIterator lowerright1, SrcAccessor sa,
    00141                     DestIterator upperleft2, DestAccessor da,
    00142                     int radius, float rank)
    00143 {
    00144     vigra_precondition((rank >= 0.0) && (rank <= 1.0),
    00145             "discRankOrderFilter(): Rank must be between 0 and 1"
    00146             " (inclusive).");
    00147
    00148     vigra_precondition(radius >= 0,
    00149             "discRankOrderFilter(): Radius must be >= 0.");
    00150
    00151     int i, x, y, xmax, ymax, xx, yy;
    00152     int rankpos, winsize, leftsum;
    00153
    00154     long hist[256];
    00155
    00156     // prepare structuring function
    00157     std::vector<int> struct_function(radius+1);
    00158     struct_function[0] = radius;
    00159
    00160     double r2 = (double)radius*radius;
    00161     for(i=1; i<=radius; ++i)
    00162     {
    00163         double r = (double) i - 0.5;
    00164         struct_function[i] = (int)(VIGRA_CSTD::sqrt(r2 - r*r) + 0.5);
    00165     }
    00166
    00167     int w = lowerright1.x - upperleft1.x;
    00168     int h = lowerright1.y - upperleft1.y;
    00169
    00170     SrcIterator ys(upperleft1);
    00171     DestIterator yd(upperleft2);
    00172
    00173     for(y=0; y<h; ++y, ++ys.y, ++yd.y)
    00174     {
    00175         SrcIterator xs(ys);
    00176         DestIterator xd(yd);
    00177
    00178         // first column
    00179         int x0 = 0;
    00180         int y0 = y;
    00181         int x1 = w - 1;
    00182         int y1 = h - y - 1;
    00183
    00184         // clear histogram
    00185         for(i=0; i<256; ++i) hist[i] = 0;
    00186         winsize = 0;
    00187
    00188         // init histogram
    00189         ymax = (y1 < radius) ? y1 : radius;
    00190         for(yy=0; yy<=ymax; ++yy)
    00191         {
    00192             xmax = (x1 < struct_function[yy]) ? x1 : struct_function[yy];
    00193             for(xx=0; xx<=xmax; ++xx)
    00194             {
    00195                 hist[sa(xs, Diff2D(xx, yy))]++;
    00196                 winsize++;
    00197             }
    00198         }
    00199
    00200         ymax = (y0 < radius) ? y0 : radius;
    00201         for(yy=1; yy<=ymax; ++yy)
    00202         {
    00203             xmax = (x1 < struct_function[yy]) ? x1 : struct_function[yy];
    00204             for(xx=0; xx<=xmax; ++xx)
    00205             {
    00206                 hist[sa(xs, Diff2D(xx, -yy))]++;
    00207                 winsize++;
    00208             }
    00209         }
    00210
    00211         // find the desired histogram bin
    00212         leftsum = 0;
    00213         if(rank == 0.0)
    00214         {
    00215             for(i=0; i<256; i++)
    00216             {
    00217                 if(hist[i]) break;
    00218             }
    00219             rankpos = i;
    00220         }
    00221         else
    00222         {
    00223             for(i=0; i<256; i++)
    00224             {
    00225                 if((float)(hist[i]+leftsum) / winsize >= rank) break;
    00226                 leftsum += hist[i];
    00227             }
    00228             rankpos = i;
    00229         }
    00230
    00231         da.set(rankpos, xd);
    00232
    00233         ++xs.x;
    00234         ++xd.x;
    00235
    00236         // inner columns
    00237         for(x=1; x<w; ++x, ++xs.x, ++xd.x)
    00238         {
    00239             x0 = x;
    00240             y0 = y;
    00241             x1 = w - x - 1;
    00242             y1 = h - y - 1;
    00243
    00244             // update histogram
    00245             // remove pixels at left border
    00246             yy = (y1 < radius) ? y1 : radius;
    00247             for(; yy>=0; yy--)
    00248             {
    00249                 unsigned char cur;
    00250                 xx = struct_function[yy]+1;
    00251                 if(xx > x0) break;
    00252
    00253                 cur = sa(xs, Diff2D(-xx, yy));
    00254
    00255                 hist[cur]--;
    00256                 if(cur < rankpos) leftsum--;
    00257                 winsize--;
    00258             }
    00259             yy = (y0 < radius) ? y0 : radius;
    00260             for(; yy>=1; yy--)
    00261             {
    00262                 unsigned char cur;
    00263                 xx = struct_function[yy]+1;
    00264                 if(xx > x0) break;
    00265
    00266                 cur = sa(xs, Diff2D(-xx, -yy));
    00267
    00268                 hist[cur]--;
    00269                 if(cur < rankpos) leftsum--;
    00270                 winsize--;
    00271             }
    00272
    00273             // add pixels at right border
    00274             yy = (y1 < radius) ? y1 : radius;
    00275             for(; yy>=0; yy--)
    00276             {
    00277                 unsigned char cur;
    00278                 xx = struct_function[yy];
    00279                 if(xx > x1) break;
    00280
    00281                 cur = sa(xs, Diff2D(xx, yy));
    00282
    00283                 hist[cur]++;
    00284                 if(cur < rankpos) leftsum++;
    00285                 winsize++;
    00286             }
    00287             yy = (y0 < radius) ? y0 : radius;
    00288             for(; yy>=1; yy--)
    00289             {
    00290                 unsigned char cur;
    00291                 xx = struct_function[yy];
    00292                 if(xx > x1) break;
    00293
    00294                 cur = sa(xs, Diff2D(xx, -yy));
    00295
    00296                 hist[cur]++;
    00297                 if(cur < rankpos) leftsum++;
    00298                 winsize++;
    00299             }
    00300
    00301             // find the desired histogram bin
    00302             if(rank == 0.0)
    00303             {
    00304                 if(leftsum == 0)
    00305                 {
    00306                     // search to the right
    00307                     for(i=rankpos; i<256; i++)
    00308                     {
    00309                         if(hist[i]) break;
    00310                     }
    00311                     rankpos = i;
    00312                 }
    00313                 else
    00314                 {
    00315                     // search to the left
    00316                     for(i=rankpos-1; i>=0; i--)
    00317                     {
    00318                         leftsum -= hist[i];
    00319                         if(leftsum == 0) break;
    00320                     }
    00321                     rankpos = i;
    00322                 }
    00323             }
    00324             else  // rank > 0.0
    00325             {
    00326                 if((float)leftsum / winsize < rank)
    00327                 {
    00328                     // search to the right
    00329                     for(i=rankpos; i<256; i++)
    00330                     {
    00331                         if((float)(hist[i]+leftsum) / winsize >= rank) break;
    00332                         leftsum+=hist[i];
    00333                     }
    00334                     rankpos = i;
    00335                 }
    00336                 else
    00337                 {
    00338                     // search to the left
    00339                     for(i=rankpos-1; i>=0; i--)
    00340                     {
    00341                         leftsum-=hist[i];
    00342                         if((float)leftsum / winsize < rank) break;
    00343                     }
    00344                     rankpos = i;
    00345                 }
    00346             }
    00347
    00348             da.set(rankpos, xd);
    00349         }
    00350     }
    00351 }
    00352
    00353 template <class SrcIterator, class SrcAccessor,
    00354           class DestIterator, class DestAccessor>
    00355 void
    00356 discRankOrderFilter(triple<SrcIterator, SrcIterator, SrcAccessor> src,
    00357                     pair<DestIterator, DestAccessor> dest,
    00358                     int radius, float rank)
    00359 {
    00360     discRankOrderFilter(src.first, src.second, src.third,
    00361                         dest.first, dest.second,
    00362                         radius, rank);
    00363 }
    00364
    00365 /********************************************************/
    00366 /*                                                      */
    00367 /*                      discErosion                     */
    00368 /*                                                      */
    00369 /********************************************************/
    00370
    00371 /** \brief Apply erosion (minimum) filter with disc of given radius to image.
    00372
    00373     This is an abbreviation vor the rank order filter with rank = 0.0.
    00374     See \ref discRankOrderFilter() for more information.
    00375
    00376     <b> Declarations:</b>
    00377
    00378     pass arguments explicitly:
    00379     \code
    00380     namespace vigra {
    00381         template <class SrcIterator, class SrcAccessor,
    00382                   class DestIterator, class DestAccessor>
    00383         void
    00384         discErosion(SrcIterator upperleft1,
    00385                     SrcIterator lowerright1, SrcAccessor sa,
    00386                     DestIterator upperleft2, DestAccessor da,
    00387                     int radius)
    00388     }
    00389     \endcode
    00390
    00391
    00392     use argument objects in conjunction with \ref ArgumentObjectFactories :
    00393     \code
    00394     namespace vigra {
    00395         template <class SrcIterator, class SrcAccessor,
    00396                   class DestIterator, class DestAccessor>
    00397         void
    00398         discErosion(triple<SrcIterator, SrcIterator, SrcAccessor> src,
    00399                     pair<DestIterator, DestAccessor> dest,
    00400                     int radius)
    00401     }
    00402     \endcode
    00403
    00404 */
    00405 doxygen_overloaded_function(template <...> void discErosion)
    00406
    00407 template <class SrcIterator, class SrcAccessor,
    00408           class DestIterator, class DestAccessor>
    00409 inline void
    00410 discErosion(SrcIterator upperleft1,
    00411             SrcIterator lowerright1, SrcAccessor sa,
    00412             DestIterator upperleft2, DestAccessor da,
    00413             int radius)
    00414 {
    00415     vigra_precondition(radius >= 0, "discErosion(): Radius must be >= 0.");
    00416
    00417     discRankOrderFilter(upperleft1, lowerright1, sa,
    00418                         upperleft2, da, radius, 0.0);
    00419 }
    00420
    00421 template <class SrcIterator, class SrcAccessor,
    00422           class DestIterator, class DestAccessor>
    00423 void
    00424 discErosion(triple<SrcIterator, SrcIterator, SrcAccessor> src,
    00425             pair<DestIterator, DestAccessor> dest,
    00426             int radius)
    00427 {
    00428     vigra_precondition(radius >= 0, "discErosion(): Radius must be >= 0.");
    00429
    00430     discRankOrderFilter(src.first, src.second, src.third,
    00431                         dest.first, dest.second,
    00432                         radius, 0.0);
    00433 }
    00434
    00435 /********************************************************/
    00436 /*                                                      */
    00437 /*                     discDilation                     */
    00438 /*                                                      */
    00439 /********************************************************/
    00440
    00441 /** \brief Apply dilation (maximum) filter with disc of given radius to image.
    00442
    00443     This is an abbreviation vor the rank order filter with rank = 1.0.
    00444     See \ref discRankOrderFilter() for more information.
    00445
    00446     <b> Declarations:</b>
    00447
    00448     pass arguments explicitly:
    00449     \code
    00450     namespace vigra {
    00451         template <class SrcIterator, class SrcAccessor,
    00452                   class DestIterator, class DestAccessor>
    00453         void
    00454         discDilation(SrcIterator upperleft1,
    00455                     SrcIterator lowerright1, SrcAccessor sa,
    00456                     DestIterator upperleft2, DestAccessor da,
    00457                     int radius)
    00458     }
    00459     \endcode
    00460
    00461
    00462     use argument objects in conjunction with \ref ArgumentObjectFactories :
    00463     \code
    00464     namespace vigra {
    00465         template <class SrcIterator, class SrcAccessor,
    00466                   class DestIterator, class DestAccessor>
    00467         void
    00468         discDilation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
    00469                     pair<DestIterator, DestAccessor> dest,
    00470                     int radius)
    00471     }
    00472     \endcode
    00473
    00474 */
    00475 doxygen_overloaded_function(template <...> void discDilation)
    00476
    00477 template <class SrcIterator, class SrcAccessor,
    00478           class DestIterator, class DestAccessor>
    00479 inline void
    00480 discDilation(SrcIterator upperleft1,
    00481             SrcIterator lowerright1, SrcAccessor sa,
    00482             DestIterator upperleft2, DestAccessor da,
    00483             int radius)
    00484 {
    00485     vigra_precondition(radius >= 0, "discDilation(): Radius must be >= 0.");
    00486
    00487     discRankOrderFilter(upperleft1, lowerright1, sa,
    00488                         upperleft2, da, radius, 1.0);
    00489 }
    00490
    00491 template <class SrcIterator, class SrcAccessor,
    00492           class DestIterator, class DestAccessor>
    00493 void
    00494 discDilation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
    00495             pair<DestIterator, DestAccessor> dest,
    00496             int radius)
    00497 {
    00498     vigra_precondition(radius >= 0, "discDilation(): Radius must be >= 0.");
    00499
    00500     discRankOrderFilter(src.first, src.second, src.third,
    00501                         dest.first, dest.second,
    00502                         radius, 1.0);
    00503 }
    00504
    00505 /********************************************************/
    00506 /*                                                      */
    00507 /*                      discMedian                      */
    00508 /*                                                      */
    00509 /********************************************************/
    00510
    00511 /** \brief Apply median filter with disc of given radius to image.
    00512
    00513     This is an abbreviation vor the rank order filter with rank = 0.5.
    00514     See \ref discRankOrderFilter() for more information.
    00515
    00516     <b> Declarations:</b>
    00517
    00518     pass arguments explicitly:
    00519     \code
    00520     namespace vigra {
    00521         template <class SrcIterator, class SrcAccessor,
    00522                   class DestIterator, class DestAccessor>
    00523         void
    00524         discMedian(SrcIterator upperleft1,
    00525                     SrcIterator lowerright1, SrcAccessor sa,
    00526                     DestIterator upperleft2, DestAccessor da,
    00527                     int radius)
    00528     }
    00529     \endcode
    00530
    00531
    00532     use argument objects in conjunction with \ref ArgumentObjectFactories :
    00533     \code
    00534     namespace vigra {
    00535         template <class SrcIterator, class SrcAccessor,
    00536                   class DestIterator, class DestAccessor>
    00537         void
    00538         discMedian(triple<SrcIterator, SrcIterator, SrcAccessor> src,
    00539                     pair<DestIterator, DestAccessor> dest,
    00540                     int radius)
    00541     }
    00542     \endcode
    00543
    00544 */
    00545 doxygen_overloaded_function(template <...> void discMedian)
    00546
    00547 template <class SrcIterator, class SrcAccessor,
    00548           class DestIterator, class DestAccessor>
    00549 inline void
    00550 discMedian(SrcIterator upperleft1,
    00551             SrcIterator lowerright1, SrcAccessor sa,
    00552             DestIterator upperleft2, DestAccessor da,
    00553             int radius)
    00554 {
    00555     vigra_precondition(radius >= 0, "discMedian(): Radius must be >= 0.");
    00556
    00557     discRankOrderFilter(upperleft1, lowerright1, sa,
    00558                         upperleft2, da, radius, 0.5);
    00559 }
    00560
    00561 template <class SrcIterator, class SrcAccessor,
    00562           class DestIterator, class DestAccessor>
    00563 void
    00564 discMedian(triple<SrcIterator, SrcIterator, SrcAccessor> src,
    00565             pair<DestIterator, DestAccessor> dest,
    00566             int radius)
    00567 {
    00568     vigra_precondition(radius >= 0, "discMedian(): Radius must be >= 0.");
    00569
    00570     discRankOrderFilter(src.first, src.second, src.third,
    00571                         dest.first, dest.second,
    00572                         radius, 0.5);
    00573 }
    00574
    00575 /********************************************************/
    00576 /*                                                      */
    00577 /*            discRankOrderFilterWithMask               */
    00578 /*                                                      */
    00579 /********************************************************/
    00580
    00581 /** \brief Apply rank order filter with disc structuring function to the image
    00582     using a mask.
    00583
    00584     The pixel values of the source image <b> must</b> be in the range
    00585     0...255. Radius must be >= 0. Rank must be in the range 0.0 <= rank
    00586     <= 1.0. The filter acts as a minimum filter if rank = 0.0,
    00587     as a median if rank = 0.5, and as a maximum filter if rank = 1.0.
    00588     Accessor are used to access the pixel data.
    00589
    00590     The mask is only applied to th input image, i.e. the function
    00591     generates an output wherever the current disc contains at least
    00592     one pixel with mask value 'true'. Source pixels with mask value
    00593     'false' are ignored during the calculation of the rank order.
    00594
    00595     <b> Declarations:</b>
    00596
    00597     pass arguments explicitly:
    00598     \code
    00599     namespace vigra {
    00600         template <class SrcIterator, class SrcAccessor,
    00601                   class MaskIterator, class MaskAccessor,
    00602                   class DestIterator, class DestAccessor>
    00603         void
    00604         discRankOrderFilterWithMask(SrcIterator upperleft1,
    00605                                     SrcIterator lowerright1, SrcAccessor sa,
    00606                                     MaskIterator upperleftm, MaskAccessor mask,
    00607                                     DestIterator upperleft2, DestAccessor da,
    00608                                     int radius, float rank)
    00609     }
    00610     \endcode
    00611
    00612
    00613     group arguments (use in conjunction with \ref ArgumentObjectFactories):
    00614     \code
    00615     namespace vigra {
    00616         template <class SrcIterator, class SrcAccessor,
    00617                   class MaskIterator, class MaskAccessor,
    00618                   class DestIterator, class DestAccessor>
    00619         void
    00620         discRankOrderFilterWithMask(triple<SrcIterator, SrcIterator, SrcAccessor> src,
    00621                                     pair<MaskIterator, MaskAccessor> mask,
    00622                                     pair<DestIterator, DestAccessor> dest,
    00623                                     int radius, float rank)
    00624     }
    00625     \endcode
    00626
    00627     <b> Usage:</b>
    00628
    00629     <b>\#include</b> <vigra/flatmorphology.hxx><br>
    00630     Namespace: vigra
    00631
    00632     \code
    00633     vigra::CImage src, dest, mask;
    00634
    00635     // do median filtering
    00636     vigra::discRankOrderFilterWithMask(srcImageRange(src),
    00637                                 maskImage(mask), destImage(dest), 10, 0.5);
    00638     \endcode
    00639
    00640     <b> Required Interface:</b>
    00641
    00642     \code
    00643     SrcIterator src_upperleft;
    00644     DestIterator dest_upperleft;
    00645     MaskIterator mask_upperleft;
    00646     int x, y;
    00647     unsigned char value;
    00648
    00649     SrcAccessor src_accessor;
    00650     DestAccessor dest_accessor;
    00651     MaskAccessor mask_accessor;
    00652
    00653     mask_accessor(mask_upperleft, x, y) // convertible to bool
    00654
    00655     // value_type of accessor must be convertible to unsigned char
    00656     value = src_accessor(src_upperleft, x, y);
    00657
    00658     dest_accessor.set(value, dest_upperleft, x, y);
    00659     \endcode
    00660
    00661     <b> Preconditions:</b>
    00662
    00663     \code
    00664     for all source pixels: 0 <= value <= 255
    00665
    00666     (rank >= 0.0) && (rank <= 1.0)
    00667     radius >= 0
    00668     \endcode
    00669
    00670 */
    00671 doxygen_overloaded_function(template <...> void discRankOrderFilterWithMask)
    00672
    00673 template <class SrcIterator, class SrcAccessor,
    00674           class MaskIterator, class MaskAccessor,
    00675           class DestIterator, class DestAccessor>
    00676 void
    00677 discRankOrderFilterWithMask(SrcIterator upperleft1,
    00678                             SrcIterator lowerright1, SrcAccessor sa,
    00679                             MaskIterator upperleftm, MaskAccessor mask,
    00680                             DestIterator upperleft2, DestAccessor da,
    00681                             int radius, float rank)
    00682 {
    00683     vigra_precondition((rank >= 0.0) && (rank <= 1.0),
    00684                  "discRankOrderFilter(): Rank must be between 0 and 1"
    00685                  " (inclusive).");
    00686
    00687     vigra_precondition(radius >= 0, "discRankOrderFilter(): Radius must be >= 0.");
    00688
    00689     int i, x, y, xmax, ymax, xx, yy;
    00690     int rankpos, winsize, leftsum;
    00691
    00692     long hist[256];
    00693
    00694     // prepare structuring function
    00695     std::vector<int> struct_function(radius+1);
    00696     struct_function[0] = radius;
    00697
    00698     double r2 = (double)radius*radius;
    00699     for(i=1; i<=radius; ++i)
    00700     {
    00701         double r = (double) i - 0.5;
    00702         struct_function[i] = (int)(VIGRA_CSTD::sqrt(r2 - r*r) + 0.5);
    00703     }
    00704
    00705     int w = lowerright1.x - upperleft1.x;
    00706     int h = lowerright1.y - upperleft1.y;
    00707
    00708     SrcIterator ys(upperleft1);
    00709     MaskIterator ym(upperleftm);
    00710     DestIterator yd(upperleft2);
    00711
    00712     for(y=0; y<h; ++y, ++ys.y, ++yd.y, ++ym.y)
    00713     {
    00714         SrcIterator xs(ys);
    00715         MaskIterator xm(ym);
    00716         DestIterator xd(yd);
    00717
    00718         // first column
    00719         int x0 = 0;
    00720         int y0 = y;
    00721         int x1 = w - 1;
    00722         int y1 = h - y - 1;
    00723
    00724         // clear histogram
    00725         for(i=0; i<256; ++i) hist[i] = 0;
    00726         winsize = 0;
    00727         leftsum = 0;
    00728         rankpos = 0;
    00729
    00730         // init histogram
    00731         ymax = (y1 < radius) ? y1 : radius;
    00732         for(yy=0; yy<=ymax; ++yy)
    00733         {
    00734             xmax = (x1 < struct_function[yy]) ? x1 : struct_function[yy];
    00735             for(xx=0; xx<=xmax; ++xx)
    00736             {
    00737                 Diff2D pos(xx, yy);
    00738                 if(mask(xm, pos))
    00739                 {
    00740                     hist[sa(xs, pos)]++;
    00741                     winsize++;
    00742                 }
    00743             }
    00744         }
    00745
    00746         ymax = (y0 < radius) ? y0 : radius;
    00747         for(yy=1; yy<=ymax; ++yy)
    00748         {
    00749             xmax = (x1 < struct_function[yy]) ? x1 : struct_function[yy];
    00750             for(xx=0; xx<=xmax; ++xx)
    00751             {
    00752                 Diff2D pos(xx, -yy);
    00753                 if(mask(xm, pos))
    00754                 {
    00755                     hist[sa(xs, pos)]++;
    00756                     winsize++;
    00757                 }
    00758             }
    00759         }
    00760
    00761         // find the desired histogram bin
    00762         if(winsize)
    00763         {
    00764             if(rank == 0.0)
    00765             {
    00766                 for(i=0; i<256; i++)
    00767                 {
    00768                     if(hist[i]) break;
    00769                 }
    00770                 rankpos = i;
    00771             }
    00772             else
    00773             {
    00774                 for(i=0; i<256; i++)
    00775                 {
    00776                     if((float)(hist[i]+leftsum) / winsize >= rank) break;
    00777                     leftsum += hist[i];
    00778                 }
    00779                 rankpos = i;
    00780             }
    00781
    00782             da.set(rankpos, xd);
    00783         }
    00784
    00785         ++xs.x;
    00786         ++xd.x;
    00787         ++xm.x;
    00788
    00789         // inner columns
    00790         for(x=1; x<w; ++x, ++xs.x, ++xd.x, ++xm.x)
    00791         {
    00792             x0 = x;
    00793             y0 = y;
    00794             x1 = w - x - 1;
    00795             y1 = h - y - 1;
    00796
    00797             // update histogram
    00798             // remove pixels at left border
    00799             yy = (y1 < radius) ? y1 : radius;
    00800             for(; yy>=0; yy--)
    00801             {
    00802                 unsigned char cur;
    00803                 xx = struct_function[yy]+1;
    00804                 if(xx > x0) break;
    00805
    00806                 Diff2D pos(-xx, yy);
    00807                 if(mask(xm, pos))
    00808                 {
    00809                     cur = sa(xs, pos);
    00810
    00811                     hist[cur]--;
    00812                     if(cur < rankpos) leftsum--;
    00813                     winsize--;
    00814                 }
    00815             }
    00816             yy = (y0 < radius) ? y0 : radius;
    00817             for(; yy>=1; yy--)
    00818             {
    00819                 unsigned char cur;
    00820                 xx = struct_function[yy]+1;
    00821                 if(xx > x0) break;
    00822
    00823                 Diff2D pos(-xx, -yy);
    00824                 if(mask(xm, pos))
    00825                 {
    00826                     cur = sa(xs, pos);
    00827
    00828                     hist[cur]--;
    00829                     if(cur < rankpos) leftsum--;
    00830                     winsize--;
    00831                 }
    00832             }
    00833
    00834             // add pixels at right border
    00835             yy = (y1 < radius) ? y1 : radius;
    00836             for(; yy>=0; yy--)
    00837             {
    00838                 unsigned char cur;
    00839                 xx = struct_function[yy];
    00840                 if(xx > x1) break;
    00841
    00842                 Diff2D pos(xx, yy);
    00843                 if(mask(xm, pos))
    00844                 {
    00845                     cur = sa(xs, pos);
    00846
    00847                     hist[cur]++;
    00848                     if(cur < rankpos) leftsum++;
    00849                     winsize++;
    00850                 }
    00851             }
    00852             yy = (y0 < radius) ? y0 : radius;
    00853             for(; yy>=1; yy--)
    00854             {
    00855                 unsigned char cur;
    00856                 xx = struct_function[yy];
    00857                 if(xx > x1) break;
    00858
    00859                 Diff2D pos(xx, -yy);
    00860                 if(mask(xm, pos))
    00861                 {
    00862                     cur = sa(xs, pos);
    00863
    00864                     hist[cur]++;
    00865                     if(cur < rankpos) leftsum++;
    00866                     winsize++;
    00867                 }
    00868             }
    00869
    00870             // find the desired histogram bin
    00871             if(winsize)
    00872             {
    00873                 if(rank == 0.0)
    00874                 {
    00875                     if(leftsum == 0)
    00876                     {
    00877                         // search to the right
    00878                         for(i=rankpos; i<256; i++)
    00879                         {
    00880                             if(hist[i]) break;
    00881                         }
    00882                         rankpos = i;
    00883                     }
    00884                     else
    00885                     {
    00886                         // search to the left
    00887                         for(i=rankpos-1; i>=0; i--)
    00888                         {
    00889                             leftsum -= hist[i];
    00890                             if(leftsum == 0) break;
    00891                         }
    00892                         rankpos = i;
    00893                     }
    00894                 }
    00895                 else  // rank > 0.0
    00896                 {
    00897                     if((float)leftsum / winsize < rank)
    00898                     {
    00899                         // search to the right
    00900                         for(i=rankpos; i<256; i++)
    00901                         {
    00902                             if((float)(hist[i]+leftsum) / winsize >= rank) break;
    00903                             leftsum+=hist[i];
    00904                         }
    00905                         rankpos = i;
    00906                     }
    00907                     else
    00908                     {
    00909                         // search to the left
    00910                         for(i=rankpos-1; i>=0; i--)
    00911                         {
    00912                             leftsum-=hist[i];
    00913                             if((float)leftsum / winsize < rank) break;
    00914                         }
    00915                         rankpos = i;
    00916                     }
    00917                 }
    00918
    00919                 da.set(rankpos, xd);
    00920             }
    00921             else
    00922             {
    00923                 leftsum = 0;
    00924                 rankpos = 0;
    00925             }
    00926         }
    00927     }
    00928 }
    00929
    00930 template <class SrcIterator, class SrcAccessor,
    00931           class MaskIterator, class MaskAccessor,
    00932           class DestIterator, class DestAccessor>
    00933 void
    00934 discRankOrderFilterWithMask(triple<SrcIterator, SrcIterator, SrcAccessor> src,
    00935                             pair<MaskIterator, MaskAccessor> mask,
    00936                             pair<DestIterator, DestAccessor> dest,
    00937                             int radius, float rank)
    00938 {
    00939     discRankOrderFilterWithMask(src.first, src.second, src.third,
    00940                         mask.first, mask.second,
    00941                         dest.first, dest.second,
    00942                         radius, rank);
    00943 }
    00944
    00945 /********************************************************/
    00946 /*                                                      */
    00947 /*                 discErosionWithMask                  */
    00948 /*                                                      */
    00949 /********************************************************/
    00950
    00951 /** \brief Apply erosion (minimum) filter with disc of given radius to image
    00952     using a mask.
    00953
    00954     This is an abbreviation vor the masked rank order filter with
    00955     rank = 0.0. See \ref discRankOrderFilterWithMask() for more information.
    00956
    00957     <b> Declarations:</b>
    00958
    00959     pass arguments explicitly:
    00960     \code
    00961     namespace vigra {
    00962         template <class SrcIterator, class SrcAccessor,
    00963                   class MaskIterator, class MaskAccessor,
    00964                   class DestIterator, class DestAccessor>
    00965         void
    00966         discErosionWithMask(SrcIterator upperleft1,
    00967                             SrcIterator lowerright1, SrcAccessor sa,
    00968                             MaskIterator upperleftm, MaskAccessor mask,
    00969                             DestIterator upperleft2, DestAccessor da,
    00970                             int radius)
    00971     }
    00972     \endcode
    00973
    00974
    00975     group arguments (use in conjunction with \ref ArgumentObjectFactories):
    00976     \code
    00977     namespace vigra {
    00978         template <class SrcIterator, class SrcAccessor,
    00979                   class MaskIterator, class MaskAccessor,
    00980                   class DestIterator, class DestAccessor>
    00981         void
    00982         discErosionWithMask(triple<SrcIterator, SrcIterator, SrcAccessor> src,
    00983                             pair<MaskIterator, MaskAccessor> mask,
    00984                             pair<DestIterator, DestAccessor> dest,
    00985                             int radius)
    00986     }
    00987     \endcode
    00988
    00989 */
    00990 doxygen_overloaded_function(template <...> void discErosionWithMask)
    00991
    00992 template <class SrcIterator, class SrcAccessor,
    00993           class MaskIterator, class MaskAccessor,
    00994           class DestIterator, class DestAccessor>
    00995 inline void
    00996 discErosionWithMask(SrcIterator upperleft1,
    00997                     SrcIterator lowerright1, SrcAccessor sa,
    00998                     MaskIterator upperleftm, MaskAccessor mask,
    00999                     DestIterator upperleft2, DestAccessor da,
    01000                     int radius)
    01001 {
    01002     vigra_precondition(radius >= 0, "discErosionWithMask(): Radius must be >= 0.");
    01003
    01004     discRankOrderFilterWithMask(upperleft1, lowerright1, sa,
    01005                                 upperleftm, mask,
    01006                                 upperleft2, da,
    01007                                 radius, 0.0);
    01008 }
    01009
    01010 template <class SrcIterator, class SrcAccessor,
    01011           class MaskIterator, class MaskAccessor,
    01012           class DestIterator, class DestAccessor>
    01013 inline void
    01014 discErosionWithMask(triple<SrcIterator, SrcIterator, SrcAccessor> src,
    01015                     pair<MaskIterator, MaskAccessor> mask,
    01016                     pair<DestIterator, DestAccessor> dest,
    01017                     int radius)
    01018 {
    01019     vigra_precondition(radius >= 0, "discErosionWithMask(): Radius must be >= 0.");
    01020
    01021     discRankOrderFilterWithMask(src.first, src.second, src.third,
    01022                         mask.first, mask.second,
    01023                         dest.first, dest.second,
    01024                         radius, 0.0);
    01025 }
    01026
    01027 /********************************************************/
    01028 /*                                                      */
    01029 /*                discDilationWithMask                  */
    01030 /*                                                      */
    01031 /********************************************************/
    01032
    01033 /** \brief Apply dilation (maximum) filter with disc of given radius to image
    01034     using a mask.
    01035
    01036     This is an abbreviation vor the masked rank order filter with
    01037     rank = 1.0. See \ref discRankOrderFilterWithMask() for more information.
    01038
    01039     <b> Declarations:</b>
    01040
    01041     pass arguments explicitly:
    01042     \code
    01043     namespace vigra {
    01044         template <class SrcIterator, class SrcAccessor,
    01045                   class MaskIterator, class MaskAccessor,
    01046                   class DestIterator, class DestAccessor>
    01047         void
    01048         discDilationWithMask(SrcIterator upperleft1,
    01049                             SrcIterator lowerright1, SrcAccessor sa,
    01050                             MaskIterator upperleftm, MaskAccessor mask,
    01051                             DestIterator upperleft2, DestAccessor da,
    01052                             int radius)
    01053     }
    01054     \endcode
    01055
    01056
    01057     group arguments (use in conjunction with \ref ArgumentObjectFactories):
    01058     \code
    01059     namespace vigra {
    01060         template <class SrcIterator, class SrcAccessor,
    01061                   class MaskIterator, class MaskAccessor,
    01062                   class DestIterator, class DestAccessor>
    01063         void
    01064         discDilationWithMask(triple<SrcIterator, SrcIterator, SrcAccessor> src,
    01065                             pair<MaskIterator, MaskAccessor> mask,
    01066                             pair<DestIterator, DestAccessor> dest,
    01067                             int radius)
    01068     }
    01069     \endcode
    01070
    01071 */
    01072 doxygen_overloaded_function(template <...> void discDilationWithMask)
    01073
    01074 template <class SrcIterator, class SrcAccessor,
    01075           class MaskIterator, class MaskAccessor,
    01076           class DestIterator, class DestAccessor>
    01077 inline void
    01078 discDilationWithMask(SrcIterator upperleft1,
    01079                     SrcIterator lowerright1, SrcAccessor sa,
    01080                     MaskIterator upperleftm, MaskAccessor mask,
    01081                     DestIterator upperleft2, DestAccessor da,
    01082                     int radius)
    01083 {
    01084     vigra_precondition(radius >= 0, "discDilationWithMask(): Radius must be >= 0.");
    01085
    01086     discRankOrderFilterWithMask(upperleft1, lowerright1, sa,
    01087                                 upperleftm, mask,
    01088                                 upperleft2, da,
    01089                                 radius, 1.0);
    01090 }
    01091
    01092 template <class SrcIterator, class SrcAccessor,
    01093           class MaskIterator, class MaskAccessor,
    01094           class DestIterator, class DestAccessor>
    01095 inline void
    01096 discDilationWithMask(triple<SrcIterator, SrcIterator, SrcAccessor> src,
    01097                     pair<MaskIterator, MaskAccessor> mask,
    01098                     pair<DestIterator, DestAccessor> dest,
    01099                     int radius)
    01100 {
    01101     vigra_precondition(radius >= 0, "discDilationWithMask(): Radius must be >= 0.");
    01102
    01103     discRankOrderFilterWithMask(src.first, src.second, src.third,
    01104                         mask.first, mask.second,
    01105                         dest.first, dest.second,
    01106                         radius, 1.0);
    01107 }
    01108
    01109 /********************************************************/
    01110 /*                                                      */
    01111 /*                 discMedianWithMask                   */
    01112 /*                                                      */
    01113 /********************************************************/
    01114
    01115 /** \brief Apply median filter with disc of given radius to image
    01116     using a mask.
    01117
    01118     This is an abbreviation vor the masked rank order filter with
    01119     rank = 0.5. See \ref discRankOrderFilterWithMask() for more information.
    01120
    01121     <b> Declarations:</b>
    01122
    01123     pass arguments explicitly:
    01124     \code
    01125     namespace vigra {
    01126         template <class SrcIterator, class SrcAccessor,
    01127                   class MaskIterator, class MaskAccessor,
    01128                   class DestIterator, class DestAccessor>
    01129         void
    01130         discMedianWithMask(SrcIterator upperleft1,
    01131                             SrcIterator lowerright1, SrcAccessor sa,
    01132                             MaskIterator upperleftm, MaskAccessor mask,
    01133                             DestIterator upperleft2, DestAccessor da,
    01134                             int radius)
    01135     }
    01136     \endcode
    01137
    01138
    01139     group arguments (use in conjunction with \ref ArgumentObjectFactories):
    01140     \code
    01141     namespace vigra {
    01142         template <class SrcIterator, class SrcAccessor,
    01143                   class MaskIterator, class MaskAccessor,
    01144                   class DestIterator, class DestAccessor>
    01145         void
    01146         discMedianWithMask(triple<SrcIterator, SrcIterator, SrcAccessor> src,
    01147                             pair<MaskIterator, MaskAccessor> mask,
    01148                             pair<DestIterator, DestAccessor> dest,
    01149                             int radius)
    01150     }
    01151     \endcode
    01152
    01153 */
    01154 doxygen_overloaded_function(template <...> void discMedianWithMask)
    01155
    01156 template <class SrcIterator, class SrcAccessor,
    01157           class MaskIterator, class MaskAccessor,
    01158           class DestIterator, class DestAccessor>
    01159 inline void
    01160 discMedianWithMask(SrcIterator upperleft1,
    01161                     SrcIterator lowerright1, SrcAccessor sa,
    01162                     MaskIterator upperleftm, MaskAccessor mask,
    01163                     DestIterator upperleft2, DestAccessor da,
    01164                     int radius)
    01165 {
    01166     vigra_precondition(radius >= 0, "discMedianWithMask(): Radius must be >= 0.");
    01167
    01168     discRankOrderFilterWithMask(upperleft1, lowerright1, sa,
    01169                                 upperleftm, mask,
    01170                                 upperleft2, da,
    01171                                 radius, 0.5);
    01172 }
    01173
    01174 template <class SrcIterator, class SrcAccessor,
    01175           class MaskIterator, class MaskAccessor,
    01176           class DestIterator, class DestAccessor>
    01177 inline void
    01178 discMedianWithMask(triple<SrcIterator, SrcIterator, SrcAccessor> src,
    01179                     pair<MaskIterator, MaskAccessor> mask,
    01180                     pair<DestIterator, DestAccessor> dest,
    01181                     int radius)
    01182 {
    01183     vigra_precondition(radius >= 0, "discMedianWithMask(): Radius must be >= 0.");
    01184
    01185     discRankOrderFilterWithMask(src.first, src.second, src.third,
    01186                         mask.first, mask.second,
    01187                         dest.first, dest.second,
    01188                         radius, 0.5);
    01189 }
    01190
    01191 //@}
    01192
    01193 } // namespace vigra
    01194
    01195 #endif // VIGRA_FLATMORPHOLOGY_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)