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

details vigra/distancetransform.hxx VIGRA

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.4.0, Dec 21 2005 )                                    */
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 /*        koethe@informatik.uni-hamburg.de          or                  */
00012 /*        vigra@kogs1.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  
00039 #ifndef VIGRA_DISTANCETRANSFORM_HXX
00040 #define VIGRA_DISTANCETRANSFORM_HXX
00041 
00042 #include <cmath>
00043 #include "vigra/stdimage.hxx"
00044 
00045 namespace vigra {
00046 
00047 /*
00048  * functors to determine the distance norm 
00049  * these functors assume that dx and dy are positive
00050  * (this is OK for use in internalDistanceTransform())
00051  */
00052  
00053 // chessboard metric
00054 struct InternalDistanceTransformLInifinityNormFunctor
00055 {
00056     float operator()(float dx, float dy) const
00057     {
00058         return (dx < dy) ? dy : dx;
00059     }
00060 };
00061 
00062 // Manhattan metric
00063 struct InternalDistanceTransformL1NormFunctor
00064 {
00065     float operator()(float dx, float dy) const
00066     {
00067         return dx + dy;
00068     }
00069 };
00070 
00071 // Euclidean metric
00072 struct InternalDistanceTransformL2NormFunctor
00073 {
00074     float operator()(float dx, float dy) const
00075     {
00076         return VIGRA_CSTD::sqrt(dx*dx + dy*dy);
00077     }
00078 };
00079 
00080 
00081 template <class SrcImageIterator, class SrcAccessor,
00082                    class DestImageIterator, class DestAccessor,
00083                    class ValueType, class NormFunctor>
00084 void
00085 internalDistanceTransform(SrcImageIterator src_upperleft, 
00086                 SrcImageIterator src_lowerright, SrcAccessor sa,
00087                 DestImageIterator dest_upperleft, DestAccessor da,
00088                 ValueType background, NormFunctor norm)
00089 {
00090     int w = src_lowerright.x - src_upperleft.x;  
00091     int h = src_lowerright.y - src_upperleft.y;  
00092     
00093     FImage xdist(w,h), ydist(w,h);
00094     
00095     xdist = w;    // init x and
00096     ydist = h;    // y distances with 'large' values
00097 
00098     SrcImageIterator sy = src_upperleft;
00099     DestImageIterator ry = dest_upperleft;
00100     FImage::Iterator xdy = xdist.upperLeft();
00101     FImage::Iterator ydy = ydist.upperLeft();
00102     SrcImageIterator sx = sy;
00103     DestImageIterator rx = ry;
00104     FImage::Iterator xdx = xdy;
00105     FImage::Iterator ydx = ydy;
00106     
00107     static const Diff2D left(-1, 0);
00108     static const Diff2D right(1, 0);
00109     static const Diff2D top(0, -1);
00110     static const Diff2D bottom(0, 1);
00111         
00112     int x,y;
00113     if(sa(sx) != background)    // first pixel
00114     {
00115         *xdx = 0.0;
00116         *ydx = 0.0;
00117         da.set(0.0, rx);
00118     }
00119     else
00120     {
00121         da.set(norm(*xdx, *ydx), rx);
00122     }
00123     
00124     
00125     for(x=1, ++xdx.x, ++ydx.x, ++sx.x, ++rx.x; 
00126         x<w; 
00127         ++x, ++xdx.x, ++ydx.x, ++sx.x, ++rx.x)   // first row left to right
00128     {
00129         if(sa(sx) != background)
00130         {
00131             *xdx = 0.0;
00132             *ydx = 0.0;
00133             da.set(0.0, rx);
00134         }
00135         else
00136         {
00137             *xdx  = xdx[left] + 1.0;   // propagate x and
00138             *ydx  = ydx[left];         // y components of distance from left pixel
00139             da.set(norm(*xdx, *ydx), rx); // calculate distance from x and y components
00140         }
00141     }
00142     for(x=w-2, xdx.x -= 2, ydx.x -= 2, sx.x -= 2, rx.x -= 2; 
00143         x>=0; 
00144         --x, --xdx.x, --ydx.x, --sx.x, --rx.x)   // first row right to left
00145     {
00146         float d = norm(xdx[right] + 1.0, ydx[right]);
00147         
00148         if(da(rx) < d) continue;
00149         
00150         *xdx = xdx[right] + 1.0;
00151         *ydx = ydx[right];
00152         da.set(d, rx);
00153     }
00154     for(y=1, ++xdy.y, ++ydy.y, ++sy.y, ++ry.y; 
00155         y<h;
00156         ++y, ++xdy.y, ++ydy.y, ++sy.y, ++ry.y)   // top to bottom
00157     {
00158         sx = sy;
00159         rx = ry;
00160         xdx = xdy;
00161         ydx = ydy;
00162         
00163         if(sa(sx) != background)    // first pixel of current row
00164         {
00165             *xdx = 0.0;
00166             *ydx = 0.0;
00167             da.set(0.0, rx);
00168         }
00169         else
00170         {
00171             *xdx = xdx[top];
00172             *ydx = ydx[top] + 1.0;
00173             da.set(norm(*xdx, *ydx), rx);
00174         }
00175         
00176         for(x=1, ++xdx.x, ++ydx.x, ++sx.x, ++rx.x; 
00177             x<w; 
00178             ++x, ++xdx.x, ++ydx.x, ++sx.x, ++rx.x)  // current row left to right
00179         {
00180             if(sa(sx) != background)
00181             {
00182                 *xdx = 0.0;
00183                 *ydx = 0.0;
00184                 da.set(0.0, rx);
00185             }
00186             else
00187             {
00188                 float d1 = norm(xdx[left] + 1.0, ydx[left]);
00189                 float d2 = norm(xdx[top], ydx[top] + 1.0);
00190                 
00191                 if(d1 < d2)
00192                 {
00193                     *xdx = xdx[left] + 1.0;
00194                     *ydx = ydx[left];
00195                     da.set(d1, rx);
00196                 }
00197                 else
00198                 {
00199                     *xdx = xdx[top];
00200                     *ydx = ydx[top] + 1.0;
00201                     da.set(d2, rx);
00202                 }
00203             }
00204         }
00205         for(x=w-2, xdx.x -= 2, ydx.x -= 2, sx.x -= 2, rx.x -= 2; 
00206             x>=0; 
00207             --x, --xdx.x, --ydx.x, --sx.x, --rx.x)  // current row right to left
00208         {
00209             float d1 = norm(xdx[right] + 1.0, ydx[right]);
00210             
00211             if(da(rx) < d1) continue;
00212             
00213             *xdx = xdx[right] + 1.0;
00214             *ydx = ydx[right];
00215             da.set(d1, rx);
00216         }
00217     }
00218     for(y=h-2, xdy.y -= 2, ydy.y -= 2, sy.y -= 2, ry.y -= 2; 
00219         y>=0;
00220         --y, --xdy.y, --ydy.y, --sy.y, --ry.y)    // bottom to top
00221     {
00222         sx = sy;
00223         rx = ry;
00224         xdx = xdy;
00225         ydx = ydy;
00226         
00227         float d = norm(xdx[bottom], ydx[bottom] + 1.0);
00228         if(d < da(rx))    // first pixel of current row
00229         { 
00230             *xdx = xdx[bottom];
00231             *ydx = ydx[bottom] + 1.0;
00232             da.set(d, rx);
00233         }
00234             
00235         for(x=1, ++xdx.x, ++ydx.x, ++sx.x, ++rx.x; 
00236             x<w;
00237             ++x, ++xdx.x, ++ydx.x, ++sx.x, ++rx.x)  // current row left to right
00238         {
00239             float d1 = norm(xdx[left] + 1.0, ydx[left]);
00240             float d2 = norm(xdx[bottom], ydx[bottom] + 1.0);
00241             
00242             if(d1 < d2)
00243             {
00244                 if(da(rx) < d1) continue;
00245                 *xdx = xdx[left] + 1.0;
00246                 *ydx = ydx[left];
00247                 da.set(d1, rx);
00248             }
00249             else
00250             {
00251                 if(da(rx) < d2) continue;
00252                 *xdx = xdx[bottom];
00253                 *ydx = ydx[bottom] + 1.0;
00254                 da.set(d2, rx);
00255             }
00256         }
00257         for(x=w-2, xdx.x -= 2, ydx.x -= 2, sx.x -= 2, rx.x -= 2; 
00258             x>=0; 
00259             --x, --xdx.x, --ydx.x, --sx.x, --rx.x)  // current row right to left
00260         {
00261             float d1 = norm(xdx[right] + 1.0, ydx[right]);
00262 
00263             if(da(rx) < d1) continue;
00264             *xdx = xdx[right] + 1.0;
00265             *ydx = ydx[right];
00266             da.set(d1, rx);
00267         }
00268     }
00269 }
00270 
00271 /********************************************************/
00272 /*                                                      */
00273 /*                 distanceTransform                    */
00274 /*                                                      */
00275 /********************************************************/
00276 
00277 /** \addtogroup DistanceTransform Distance Transform
00278     Perform a distance transform using either the Euclidean, Manhattan, 
00279     or chessboard metrics.
00280 */
00281 //@{
00282 
00283 /** For all background pixels, calculate the distance to 
00284     the nearest object or contour. The label of the pixels to be considered 
00285     background in the source image is passed in the parameter 'background'.
00286     Source pixels with other labels will be considered objects. In the 
00287     destination image, all pixels corresponding to background will be assigned 
00288     the their distance value, all pixels corresponding to objects will be
00289     assigned 0.
00290     
00291     The parameter 'norm' gives the distance norm to be used:
00292     
00293     <ul>
00294 
00295     <li> norm == 0: use chessboard distance (L-infinity norm)
00296     <li> norm == 1: use Manhattan distance (L1 norm)
00297     <li> norm == 2: use Euclidean distance (L2 norm)
00298     
00299     </ul>
00300     
00301     If you use the L2 norm, the destination pixels must be real valued to give
00302     correct results.
00303     
00304     <b> Declarations:</b>
00305     
00306     pass arguments explicitly:
00307     \code
00308     namespace vigra {
00309         template <class SrcImageIterator, class SrcAccessor,
00310                            class DestImageIterator, class DestAccessor,
00311                            class ValueType>
00312         void distanceTransform(SrcImageIterator src_upperleft, 
00313                         SrcImageIterator src_lowerright, SrcAccessor sa,
00314                         DestImageIterator dest_upperleft, DestAccessor da,
00315                         ValueType background, int norm)
00316     }
00317     \endcode
00318     
00319     
00320     use argument objects in conjunction with \ref ArgumentObjectFactories:
00321     \code
00322     namespace vigra {
00323         template <class SrcImageIterator, class SrcAccessor,
00324                            class DestImageIterator, class DestAccessor,
00325                            class ValueType>
00326         void distanceTransform(
00327             triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src,
00328             pair<DestImageIterator, DestAccessor> dest,
00329             ValueType background, int norm)
00330     }
00331     \endcode
00332     
00333     <b> Usage:</b>
00334     
00335     <b>\#include</b> "<a href="distancetransform_8hxx-source.html">vigra/distancetransform.hxx</a>"<br>
00336     Namespace: vigra
00337     
00338     
00339     \code
00340     
00341     vigra::BImage src(w,h), edges(w,h);
00342     vigra::FImage distance(w, h);
00343 
00344     // empty edge image
00345     edges = 0;
00346     ...
00347 
00348     // detect edges in src image (edges will be marked 1, background 0)
00349     vigra::differenceOfExponentialEdgeImage(srcImageRange(src), destImage(edges), 
00350                                      0.8, 4.0);
00351      
00352     // find distance of all pixels from nearest edge
00353     vigra::distanceTransform(srcImageRange(edges), destImage(distance),
00354                              0,                   2);
00355     //                       ^ background label   ^ norm (Euclidean)
00356     \endcode
00357 
00358     <b> Required Interface:</b>
00359     
00360     \code
00361     
00362     SrcImageIterator src_upperleft, src_lowerright;
00363     DestImageIterator dest_upperleft;
00364     
00365     SrcAccessor sa;
00366     DestAccessor da;
00367     
00368     ValueType background;
00369     float distance;
00370     
00371     sa(src_upperleft) != background;
00372     da(dest_upperleft) < distance;
00373     da.set(distance, dest_upperleft);
00374  
00375     \endcode
00376 */
00377 template <class SrcImageIterator, class SrcAccessor,
00378                    class DestImageIterator, class DestAccessor,
00379                    class ValueType>
00380 inline void
00381 distanceTransform(SrcImageIterator src_upperleft, 
00382                 SrcImageIterator src_lowerright, SrcAccessor sa,
00383                 DestImageIterator dest_upperleft, DestAccessor da,
00384                 ValueType background, int norm)
00385 {
00386     if(norm == 1)
00387     {
00388         internalDistanceTransform(src_upperleft, src_lowerright, sa,
00389                                   dest_upperleft, da, background,
00390                                   InternalDistanceTransformL1NormFunctor());
00391     }
00392     else if(norm == 2)
00393     {
00394         internalDistanceTransform(src_upperleft, src_lowerright, sa,
00395                                   dest_upperleft, da, background,
00396                                   InternalDistanceTransformL2NormFunctor());
00397     }
00398     else
00399     {
00400         internalDistanceTransform(src_upperleft, src_lowerright, sa,
00401                                   dest_upperleft, da, background,
00402                                   InternalDistanceTransformLInifinityNormFunctor());
00403     }
00404 }
00405 
00406 template <class SrcImageIterator, class SrcAccessor,
00407                    class DestImageIterator, class DestAccessor,
00408                    class ValueType>
00409 inline void
00410 distanceTransform(
00411     triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src,
00412     pair<DestImageIterator, DestAccessor> dest,
00413     ValueType background, int norm)
00414 {
00415     distanceTransform(src.first, src.second, src.third,
00416                       dest.first, dest.second, background, norm);
00417 }
00418 
00419 //@}
00420 
00421 } // namespace vigra
00422 
00423 #endif // VIGRA_DISTANCETRANSFORM_HXX
00424 

© Ullrich Köthe (koethe@informatik.uni-hamburg.de)
Cognitive Systems Group, University of Hamburg, Germany

html generated using doxygen and Python
VIGRA 1.4.0 (21 Dec 2005)