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

# vigra/seededregiongrowing.hxx

    00001 /************************************************************************/
    00002 /*                                                                      */
    00003 /*         Copyright 1998-2003 by Ullrich Koethe, Hans Meine            */
    00004 /*       Cognitive Systems Group, University of Hamburg, Germany        */
    00005 /*                                                                      */
    00006 /*    This file is part of the VIGRA computer vision library.           */
    00007 /*    ( Version 1.6.0, Aug 13 2008 )                                    */
    00008 /*    The VIGRA Website is                                              */
    00009 /*        http://kogs-www.informatik.uni-hamburg.de/~koethe/vigra/      */
    00010 /*    Please direct questions, bug reports, and contributions to        */
    00011 /*        ullrich.koethe@iwr.uni-heidelberg.de    or                    */
    00012 /*        vigra@informatik.uni-hamburg.de                               */
    00013 /*                                                                      */
    00014 /*    Permission is hereby granted, free of charge, to any person       */
    00015 /*    obtaining a copy of this software and associated documentation    */
    00016 /*    files (the "Software"), to deal in the Software without           */
    00017 /*    restriction, including without limitation the rights to use,      */
    00018 /*    copy, modify, merge, publish, distribute, sublicense, and/or      */
    00019 /*    sell copies of the Software, and to permit persons to whom the    */
    00020 /*    Software is furnished to do so, subject to the following          */
    00021 /*    conditions:                                                       */
    00022 /*                                                                      */
    00023 /*    The above copyright notice and this permission notice shall be    */
    00024 /*    included in all copies or substantial portions of the             */
    00025 /*    Software.                                                         */
    00026 /*                                                                      */
    00027 /*    THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND    */
    00028 /*    EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES   */
    00029 /*    OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND          */
    00030 /*    NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT       */
    00031 /*    HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,      */
    00032 /*    WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING      */
    00033 /*    FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR     */
    00034 /*    OTHER DEALINGS IN THE SOFTWARE.                                   */
    00035 /*                                                                      */
    00036 /************************************************************************/
    00037
    00038 #ifndef VIGRA_SEEDEDREGIONGROWING_HXX
    00039 #define VIGRA_SEEDEDREGIONGROWING_HXX
    00040
    00041 #include <vector>
    00042 #include <stack>
    00043 #include <queue>
    00044 #include "utilities.hxx"
    00045 #include "stdimage.hxx"
    00046 #include "stdimagefunctions.hxx"
    00047
    00048 namespace vigra {
    00049
    00050 namespace detail {
    00051
    00052 template <class COST>
    00053 class SeedRgPixel
    00054 {
    00055 public:
    00056     Point2D location_, nearest_;
    00057     COST cost_;
    00058     int count_;
    00059     int label_;
    00060     int dist_;
    00061
    00062     SeedRgPixel()
    00063     : location_(0,0), nearest_(0,0), cost_(0), count_(0), label_(0)
    00064     {}
    00065
    00066     SeedRgPixel(Point2D const & location, Point2D const & nearest,
    00067                 COST const & cost, int const & count, int const & label)
    00068     : location_(location), nearest_(nearest),
    00069       cost_(cost), count_(count), label_(label)
    00070     {
    00071         int dx = location_.x - nearest_.x;
    00072         int dy = location_.y - nearest_.y;
    00073         dist_ = dx * dx + dy * dy;
    00074     }
    00075
    00076     void set(Point2D const & location, Point2D const & nearest,
    00077              COST const & cost, int const & count, int const & label)
    00078     {
    00079         location_ = location;
    00080         nearest_ = nearest;
    00081         cost_ = cost;
    00082         count_ = count;
    00083         label_ = label;
    00084
    00085         int dx = location_.x - nearest_.x;
    00086         int dy = location_.y - nearest_.y;
    00087         dist_ = dx * dx + dy * dy;
    00088     }
    00089
    00090     struct Compare
    00091     {
    00092         // must implement > since priority_queue looks for largest element
    00093         bool operator()(SeedRgPixel const & l,
    00094                         SeedRgPixel const & r) const
    00095         {
    00096             if(r.cost_ == l.cost_)
    00097             {
    00098                 if(r.dist_ == l.dist_) return r.count_ < l.count_;
    00099
    00100                 return r.dist_ < l.dist_;
    00101             }
    00102
    00103             return r.cost_ < l.cost_;
    00104         }
    00105         bool operator()(SeedRgPixel const * l,
    00106                         SeedRgPixel const * r) const
    00107         {
    00108             if(r->cost_ == l->cost_)
    00109             {
    00110                 if(r->dist_ == l->dist_) return r->count_ < l->count_;
    00111
    00112                 return r->dist_ < l->dist_;
    00113             }
    00114
    00115             return r->cost_ < l->cost_;
    00116         }
    00117     };
    00118
    00119     struct Allocator
    00120     {
    00121         ~Allocator()
    00122         {
    00123             while(!freelist_.empty())
    00124             {
    00125                 delete freelist_.top();
    00126                 freelist_.pop();
    00127             }
    00128         }
    00129
    00130         SeedRgPixel *
    00131         create(Point2D const & location, Point2D const & nearest,
    00132                COST const & cost, int const & count, int const & label)
    00133         {
    00134             if(!freelist_.empty())
    00135             {
    00136                 SeedRgPixel * res = freelist_.top();
    00137                 freelist_.pop();
    00138                 res->set(location, nearest, cost, count, label);
    00139                 return res;
    00140             }
    00141
    00142             return new SeedRgPixel(location, nearest, cost, count, label);
    00143         }
    00144
    00145         void dismiss(SeedRgPixel * p)
    00146         {
    00147             freelist_.push(p);
    00148         }
    00149
    00150         std::stack<SeedRgPixel<COST> *> freelist_;
    00151     };
    00152 };
    00153
    00154 struct UnlabelWatersheds
    00155 {
    00156     int operator()(int label) const
    00157     {
    00158         return label < 0 ? 0 : label;
    00159     }
    00160 };
    00161
    00162 } // namespace detail
    00163
    00164 enum SRGType { KeepContours, CompleteGrow, SRGWatershedLabel = -1 };
    00165
    00166 /** \addtogroup SeededRegionGrowing Region Segmentation Algorithms
    00167     Region growing, watersheds, and voronoi tesselation
    00168 */
    00169 //@{
    00170
    00171 /********************************************************/
    00172 /*                                                      */
    00173 /*                    seededRegionGrowing               */
    00174 /*                                                      */
    00175 /********************************************************/
    00176
    00177 /** \brief Region Segmentation by means of Seeded Region Growing.
    00178
    00179     This algorithm implements seeded region growing as described in
    00180
    00181     R. Adams, L. Bischof: "<em> Seeded Region Growing</em>", IEEE Trans. on Pattern
    00182     Analysis and Maschine Intelligence, vol 16, no 6, 1994, and
    00183
    00184     Ullrich K&ouml;the:
    00185     <em><a href="http://kogs-www.informatik.uni-hamburg.de/~koethe/papers/#primary">Primary Image Segmentation</a></em>,
    00186     in: G. Sagerer, S.
    00187     Posch, F. Kummert (eds.): Mustererkennung 1995, Proc. 17. DAGM-Symposium,
    00188     Springer 1995
    00189
    00190     The seed image is a partly segmented image which contains uniquely
    00191     labeled regions (the seeds) and unlabeled pixels (the candidates, label 0).
    00192     Seed regions can be as large as you wish and as small as one pixel. If
    00193     there are no candidates, the algorithm will simply copy the seed image
    00194     into the output image. Otherwise it will aggregate the candidates into
    00195     the existing regions so that a cost function is minimized. This
    00196     works as follows:
    00197
    00198     <ol>
    00199
    00200     <li> Find all candidate pixels that are 4-adjacent to a seed region.
    00201     Calculate the cost for aggregating each candidate into its adajacent region
    00202     and put the candidates into a priority queue.
    00203
    00204     <li> While( priority queue is not empty)
    00205
    00206         <ol>
    00207
    00208         <li> Take the candidate with least cost from the queue. If it has not
    00209         already been merged, merge it with it's adjacent region.
    00210
    00211         <li> Put all candidates that are 4-adjacent to the pixel just processed
    00212         into the priority queue.
    00213
    00214         </ol>
    00215
    00216     </ol>
    00217
    00218     If <tt>SRGType == CompleteGrow</tt> (the default), this algorithm
    00219     will produce a complete 4-connected tesselation of the image.
    00220     If <tt>SRGType == KeepContours</tt>, a one-pixel-wide border will be left
    00221     between the regions. The border pixels get label 0 (zero).
    00222
    00223     The cost is determined jointly by the source image and the
    00224     region statistics functor. The source image contains feature values for each
    00225     pixel which will be used by the region statistics functor to calculate and
    00226     update statistics for each region and to calculate the cost for each
    00227     candidate. The <TT>RegionStatisticsArray</TT> must be compatible to the
    00228     \ref ArrayOfRegionStatistics functor and contains an <em> array</em> of
    00229     statistics objects for each region. The indices must correspond to the
    00230     labels of the seed regions. The statistics for the initial regions must have
    00231     been calculated prior to calling <TT>seededRegionGrowing()</TT> (for example by
    00232     means of \ref inspectTwoImagesIf()).
    00233
    00234     For each candidate
    00235     <TT>x</TT> that is adjacent to region <TT>i</TT>, the algorithm will call
    00236     <TT>stats[i].cost(as(x))</TT> to get the cost (where <TT>x</TT> is a <TT>SrcImageIterator</TT>
    00237     and <TT>as</TT> is
    00238     the SrcAccessor). When a candidate has been merged with a region, the
    00239     statistics are updated by calling <TT>stats[i].operator()(as(x))</TT>. Since
    00240     the <TT>RegionStatisticsArray</TT> is passed by reference, this will overwrite
    00241     the original statistics.
    00242
    00243     If a candidate could be merged into more than one regions with identical
    00244     cost, the algorithm will favour the nearest region.
    00245
    00246     In some cases, the cost only depends on the feature value of the current
    00247     pixel. Then the update operation will simply be a no-op, and the <TT>cost()</TT>
    00248     function returns its argument. This behavior is implemented by the
    00249     \ref SeedRgDirectValueFunctor. With <tt>SRGType == KeepContours</tt>,
    00250     this is equivalent to the watershed algorithm.
    00251
    00252     <b> Declarations:</b>
    00253
    00254     pass arguments explicitly:
    00255     \code
    00256     namespace vigra {
    00257         template <class SrcImageIterator, class SrcAccessor,
    00258                   class SeedImageIterator, class SeedAccessor,
    00259                   class DestImageIterator, class DestAccessor,
    00260                   class RegionStatisticsArray>
    00261         void seededRegionGrowing(SrcImageIterator srcul,
    00262                                  SrcImageIterator srclr, SrcAccessor as,
    00263                                  SeedImageIterator seedsul, SeedAccessor aseeds,
    00264                                  DestImageIterator destul, DestAccessor ad,
    00265                                  RegionStatisticsArray & stats,
    00266                                  SRGType srgType = CompleteGrow);
    00267     }
    00268     \endcode
    00269
    00270     use argument objects in conjunction with \ref ArgumentObjectFactories :
    00271     \code
    00272     namespace vigra {
    00273         template <class SrcImageIterator, class SrcAccessor,
    00274                   class SeedImageIterator, class SeedAccessor,
    00275                   class DestImageIterator, class DestAccessor,
    00276                   class RegionStatisticsArray>
    00277         void
    00278         seededRegionGrowing(triple<SrcImageIterator, SrcImageIterator, SrcAccessor> img1,
    00279                             pair<SeedImageIterator, SeedAccessor> img3,
    00280                             pair<DestImageIterator, DestAccessor> img4,
    00281                             RegionStatisticsArray & stats,
    00282                             SRGType srgType = CompleteGrow);
    00283     }
    00284     \endcode
    00285
    00286     <b> Usage:</b>
    00287
    00288     <b>\#include</b> <<a href="seededregiongrowing_8hxx-source.html">vigra/seededregiongrowing.hxx</a>><br>
    00289     Namespace: vigra
    00290
    00291     Example: implementation of the voronoi tesselation
    00292
    00293     \code
    00294     vigra::BImage points(w,h);
    00295     vigra::FImage dist(x,y);
    00296
    00297     // empty edge image
    00298     points = 0;
    00299     dist = 0;
    00300
    00301     int max_region_label = 100;
    00302
    00303     // throw in some random points:
    00304     for(int i = 1; i <= max_region_label; ++i)
    00305            points(w * rand() / RAND_MAX , h * rand() / RAND_MAX) = i;
    00306
    00307     // calculate Euclidean distance transform
    00308     vigra::distanceTransform(srcImageRange(points), destImage(dist), 2);
    00309
    00310     // init statistics functor
    00311     vigra::ArrayOfRegionStatistics<vigra::SeedRgDirectValueFunctor<float> >
    00312                                               stats(max_region_label);
    00313
    00314     // find voronoi region of each point
    00315    vigra:: seededRegionGrowing(srcImageRange(dist), srcImage(points),
    00316                                destImage(points), stats);
    00317     \endcode
    00318
    00319     <b> Required Interface:</b>
    00320
    00321     \code
    00322     SrcImageIterator src_upperleft, src_lowerright;
    00323     SeedImageIterator seed_upperleft;
    00324     DestImageIterator dest_upperleft;
    00325
    00326     SrcAccessor src_accessor;
    00327     SeedAccessor seed_accessor;
    00328     DestAccessor dest_accessor;
    00329
    00330     RegionStatisticsArray stats;
    00331
    00332     // calculate costs
    00333     RegionStatisticsArray::value_type::cost_type cost =
    00334         stats[seed_accessor(seed_upperleft)].cost(src_accessor(src_upperleft));
    00335
    00336     // compare costs
    00337     cost < cost;
    00338
    00339     // update statistics
    00340     stats[seed_accessor(seed_upperleft)](src_accessor(src_upperleft));
    00341
    00342     // set result
    00343     dest_accessor.set(seed_accessor(seed_upperleft), dest_upperleft);
    00344     \endcode
    00345
    00346     Further requirements are determined by the <TT>RegionStatisticsArray</TT>.
    00347 */
    00348 doxygen_overloaded_function(template <...> void seededRegionGrowing)
    00349
    00350 template <class SrcImageIterator, class SrcAccessor,
    00351           class SeedImageIterator, class SeedAccessor,
    00352           class DestImageIterator, class DestAccessor,
    00353           class RegionStatisticsArray>
    00354 void seededRegionGrowing(SrcImageIterator srcul,
    00355                          SrcImageIterator srclr, SrcAccessor as,
    00356                          SeedImageIterator seedsul, SeedAccessor aseeds,
    00357                          DestImageIterator destul, DestAccessor ad,
    00358                          RegionStatisticsArray & stats,
    00359                          const SRGType srgType)
    00360 {
    00361     int w = srclr.x - srcul.x;
    00362     int h = srclr.y - srcul.y;
    00363     int count = 0;
    00364
    00365     SrcImageIterator isy = srcul, isx = srcul;  // iterators for the src image
    00366
    00367     typedef typename RegionStatisticsArray::value_type RegionStatistics;
    00368     typedef typename RegionStatistics::cost_type CostType;
    00369     typedef detail::SeedRgPixel<CostType> Pixel;
    00370
    00371     typename Pixel::Allocator allocator;
    00372
    00373     typedef std::priority_queue<Pixel *, std::vector<Pixel *>,
    00374                                 typename Pixel::Compare>  SeedRgPixelHeap;
    00375
    00376     // copy seed image in an image with border
    00377     IImage regions(w+2, h+2);
    00378     IImage::Iterator ir = regions.upperLeft() + Diff2D(1,1);
    00379     IImage::Iterator iry, irx;
    00380
    00381     initImageBorder(destImageRange(regions), 1, SRGWatershedLabel);
    00382     copyImage(seedsul, seedsul+Diff2D(w,h), aseeds, ir, regions.accessor());
    00383
    00384     // allocate and init memory for the results
    00385
    00386     SeedRgPixelHeap pheap;
    00387     int cneighbor;
    00388
    00389     static const Diff2D dist[] = { Diff2D(-1,0), Diff2D(0,-1),
    00390                                    Diff2D(1,0),  Diff2D(0,1) };
    00391
    00392     Point2D pos(0,0);
    00393     for(isy=srcul, iry=ir, pos.y=0; pos.y<h;
    00394         ++pos.y, ++isy.y, ++iry.y)
    00395     {
    00396         for(isx=isy, irx=iry, pos.x=0; pos.x<w;
    00397             ++pos.x, ++isx.x, ++irx.x)
    00398         {
    00399             if(*irx == 0)
    00400             {
    00401                 // find candidate pixels for growing and fill heap
    00402                 for(int i=0; i<4; i++)
    00403                 {
    00404                     cneighbor = irx[dist[i]];
    00405                     if(cneighbor > 0)
    00406                     {
    00407                         CostType cost = stats[cneighbor].cost(as(isx));
    00408
    00409                         Pixel * pixel =
    00410                             allocator.create(pos, pos+dist[i], cost, count++, cneighbor);
    00411                         pheap.push(pixel);
    00412                     }
    00413                 }
    00414             }
    00415         }
    00416     }
    00417
    00418     // perform region growing
    00419     while(pheap.size() != 0)
    00420     {
    00421         Pixel * pixel = pheap.top();
    00422         pheap.pop();
    00423
    00424         Point2D pos = pixel->location_;
    00425         Point2D nearest = pixel->nearest_;
    00426         int lab = pixel->label_;
    00427
    00428         allocator.dismiss(pixel);
    00429
    00430         irx = ir + pos;
    00431         isx = srcul + pos;
    00432
    00433         if(*irx) // already labelled region / watershed?
    00434             continue;
    00435
    00436         if(srgType == KeepContours)
    00437         {
    00438             for(int i=0; i<4; i++)
    00439             {
    00440                 cneighbor = irx[dist[i]];
    00441                 if((cneighbor>0) && (cneighbor != lab))
    00442                 {
    00443                     lab = SRGWatershedLabel;
    00444                     break;
    00445                 }
    00446             }
    00447         }
    00448
    00449         *irx = lab;
    00450
    00451         if((srgType != KeepContours) || (lab > 0))
    00452         {
    00453             // update statistics
    00454             stats[*irx](as(isx));
    00455
    00456             // search neighborhood
    00457             // second pass: find new candidate pixels
    00458             for(int i=0; i<4; i++)
    00459             {
    00460                 if(irx[dist[i]] == 0)
    00461                 {
    00462                     CostType cost = stats[lab].cost(as(isx, dist[i]));
    00463
    00464                     Pixel * new_pixel =
    00465                         allocator.create(pos+dist[i], nearest, cost, count++, lab);
    00466                     pheap.push(new_pixel);
    00467                 }
    00468             }
    00469         }
    00470     }
    00471
    00472     // write result
    00473     transformImage(ir, ir+Point2D(w,h), regions.accessor(), destul, ad,
    00474                    detail::UnlabelWatersheds());
    00475 }
    00476
    00477 template <class SrcImageIterator, class SrcAccessor,
    00478           class SeedImageIterator, class SeedAccessor,
    00479           class DestImageIterator, class DestAccessor,
    00480           class RegionStatisticsArray>
    00481 inline void
    00482 seededRegionGrowing(SrcImageIterator srcul,
    00483                     SrcImageIterator srclr, SrcAccessor as,
    00484                     SeedImageIterator seedsul, SeedAccessor aseeds,
    00485                     DestImageIterator destul, DestAccessor ad,
    00486                     RegionStatisticsArray & stats)
    00487 {
    00488     seededRegionGrowing(srcul, srclr, as,
    00489                         seedsul, aseeds,
    00490                         destul, ad,
    00491                         stats, CompleteGrow);
    00492 }
    00493
    00494 template <class SrcImageIterator, class SrcAccessor,
    00495           class SeedImageIterator, class SeedAccessor,
    00496           class DestImageIterator, class DestAccessor,
    00497           class RegionStatisticsArray>
    00498 inline void
    00499 seededRegionGrowing(triple<SrcImageIterator, SrcImageIterator, SrcAccessor> img1,
    00500                     pair<SeedImageIterator, SeedAccessor> img3,
    00501                     pair<DestImageIterator, DestAccessor> img4,
    00502                     RegionStatisticsArray & stats,
    00503                     SRGType srgType)
    00504 {
    00505     seededRegionGrowing(img1.first, img1.second, img1.third,
    00506                         img3.first, img3.second,
    00507                         img4.first, img4.second,
    00508                         stats, srgType);
    00509 }
    00510
    00511 template <class SrcImageIterator, class SrcAccessor,
    00512           class SeedImageIterator, class SeedAccessor,
    00513           class DestImageIterator, class DestAccessor,
    00514           class RegionStatisticsArray>
    00515 inline void
    00516 seededRegionGrowing(triple<SrcImageIterator, SrcImageIterator, SrcAccessor> img1,
    00517                     pair<SeedImageIterator, SeedAccessor> img3,
    00518                     pair<DestImageIterator, DestAccessor> img4,
    00519                     RegionStatisticsArray & stats)
    00520 {
    00521     seededRegionGrowing(img1.first, img1.second, img1.third,
    00522                         img3.first, img3.second,
    00523                         img4.first, img4.second,
    00524                         stats, CompleteGrow);
    00525 }
    00526
    00527 /********************************************************/
    00528 /*                                                      */
    00529 /*               SeedRgDirectValueFunctor               */
    00530 /*                                                      */
    00531 /********************************************************/
    00532
    00533 /** \brief Statistics functor to be used for seeded region growing.
    00534
    00535     This functor can be used if the cost of a candidate during
    00536     \ref seededRegionGrowing() is equal to the feature value of that
    00537     candidate and does not depend on properties of the region it is going to
    00538     be merged with.
    00539
    00540     <b>\#include</b> <<a href="seededregiongrowing_8hxx-source.html">vigra/seededregiongrowing.hxx</a>><br>
    00541     Namespace: vigra
    00542
    00543
    00544      <b> Required Interface:</b>
    00545
    00546      no requirements
    00547 */
    00548 template <class Value>
    00549 class SeedRgDirectValueFunctor
    00550 {
    00551   public:
    00552         /** the functor's argument type
    00553         */
    00554     typedef Value argument_type;
    00555
    00556         /** the functor's result type (unused, only necessary for
    00557             use of SeedRgDirectValueFunctor in \ref vigra::ArrayOfRegionStatistics
    00558         */
    00559     typedef Value result_type;
    00560
    00561         /** \deprecated use argument_type
    00562         */
    00563     typedef Value value_type;
    00564
    00565         /** the return type of the cost() function
    00566         */
    00567     typedef Value cost_type;
    00568
    00569         /** Do nothing (since we need not update region statistics).
    00570         */
    00571     void operator()(argument_type const &) const {}
    00572
    00573         /** Return argument (since cost is identical to feature value)
    00574         */
    00575     cost_type const & cost(argument_type const & v) const
    00576     {
    00577         return v;
    00578     }
    00579 };
    00580
    00581 //@}
    00582
    00583 } // namespace vigra
    00584
    00585 #endif // VIGRA_SEEDEDREGIONGROWING_HXX
    00586

© 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.6.0 (13 Aug 2008)