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

| vigra/flatmorphology.hxx |
--- | --- | ---

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

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