[ VIGRA Homepage | Class Index | Function Index | File Index | Main Page ]
![]() |
vigra/gaborfilter.hxx | ![]() |
---|
00001 /************************************************************************/ 00002 /* */ 00003 /* Copyright 2002-2004 by Ullrich Koethe and 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.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_GABORFILTER_HXX 00040 #define VIGRA_GABORFILTER_HXX 00041 00042 #include "vigra/imagecontainer.hxx" 00043 #include "vigra/config.hxx" 00044 #include "vigra/stdimage.hxx" 00045 #include "vigra/copyimage.hxx" 00046 #include "vigra/transformimage.hxx" 00047 #include "vigra/combineimages.hxx" 00048 #include "vigra/utilities.hxx" 00049 00050 #include <functional> 00051 #include <vector> 00052 #include <cmath> 00053 00054 namespace vigra { 00055 00056 /** \addtogroup GaborFilter Gabor Filter 00057 Functions to create or apply gabor filter (latter based on FFTW). 00058 */ 00059 //@{ 00060 00061 /********************************************************/ 00062 /* */ 00063 /* createGaborFilter */ 00064 /* */ 00065 /********************************************************/ 00066 00067 /** \brief Create a gabor filter in frequency space. 00068 00069 The orientation is given in radians, the other parameters are the 00070 center frequency (for example 0.375 or smaller) and the two 00071 angular and radial sigmas of the gabor filter. (See \ref 00072 angularGaborSigma() for an explanation of possible values.) 00073 00074 The energy of the filter is explicitely normalized to 1.0. 00075 00076 <b> Declarations:</b> 00077 00078 pass arguments explicitly: 00079 \code 00080 namespace vigra { 00081 template <class DestImageIterator, class DestAccessor> 00082 void createGaborFilter(DestImageIterator destUpperLeft, 00083 DestImageIterator destLowerRight, 00084 DestAccessor da, 00085 double orientation, double centerFrequency, 00086 double angularSigma, double radialSigma) 00087 } 00088 \endcode 00089 00090 use argument objects in conjunction with \ref ArgumentObjectFactories: 00091 \code 00092 namespace vigra { 00093 template <class DestImageIterator, class DestAccessor> 00094 inline 00095 void createGaborFilter(triple<DestImageIterator, 00096 DestImageIterator, 00097 DestAccessor> dest, 00098 double orientation, double centerFrequency, 00099 double angularSigma, double radialSigma) 00100 } 00101 \endcode 00102 00103 <b> Usage:</b> 00104 00105 <b>\#include</b> "<a href="gaborfilter_8hxx-source.html">vigra/gaborfilter.hxx</a>"<br> 00106 00107 Namespace: vigra 00108 00109 \code 00110 vigra::FImage gabor(w,h); 00111 00112 vigra::createGaborFilter(destImageRange(gabor), orient, freq, 00113 angularGaborSigma(directionCount, freq) 00114 radialGaborSigma(freq)); 00115 \endcode 00116 */ 00117 00118 template <class DestImageIterator, class DestAccessor> 00119 void createGaborFilter(DestImageIterator destUpperLeft, 00120 DestImageIterator destLowerRight, DestAccessor da, 00121 double orientation, double centerFrequency, 00122 double angularSigma, double radialSigma) 00123 { 00124 int w= destLowerRight.x - destUpperLeft.x; 00125 int h= destLowerRight.y - destUpperLeft.y; 00126 00127 double squaredSum = 0.0; 00128 double cosTheta= VIGRA_CSTD::cos(orientation); 00129 double sinTheta= VIGRA_CSTD::sin(orientation); 00130 00131 double radialSigma2 = radialSigma*radialSigma; 00132 double angularSigma2 = angularSigma*angularSigma; 00133 00134 double wscale = w % 1 ? 00135 1.0f / (w-1) : 00136 1.0f / w; 00137 double hscale = h % 1 ? 00138 1.0f / (h-1) : 00139 1.0f / h; 00140 00141 int dcX= (w+1)/2, dcY= (h+1)/2; 00142 00143 double u, v; 00144 for ( int y=0; y<h; y++, destUpperLeft.y++ ) 00145 { 00146 typename DestImageIterator::row_iterator dix = destUpperLeft.rowIterator(); 00147 00148 v = hscale * ((h - (y - dcY))%h - dcY); 00149 for ( int x=0; x<w; x++, dix++ ) 00150 { 00151 u= wscale*((x - dcX + w)%w - dcX); 00152 00153 double uu = cosTheta*u + sinTheta*v - centerFrequency; 00154 double vv = -sinTheta*u + cosTheta*v; 00155 double gabor; 00156 00157 gabor = VIGRA_CSTD::exp(-0.5*(uu*uu / radialSigma2 + vv*vv / angularSigma2)); 00158 squaredSum += gabor * gabor; 00159 da.set( gabor, dix ); 00160 } 00161 } 00162 destUpperLeft.y -= h; 00163 00164 // clear out DC value and remove it from the squared sum 00165 double dcValue = da(destUpperLeft); 00166 squaredSum -= dcValue * dcValue; 00167 da.set( 0.0, destUpperLeft ); 00168 00169 // normalize energy to one 00170 double factor = VIGRA_CSTD::sqrt(squaredSum); 00171 for ( int y=0; y<h; y++, destUpperLeft.y++ ) 00172 { 00173 typename DestImageIterator::row_iterator dix = destUpperLeft.rowIterator(); 00174 00175 for ( int x=0; x<w; x++, dix++ ) 00176 { 00177 da.set( da(dix) / factor, dix ); 00178 } 00179 } 00180 } 00181 00182 template <class DestImageIterator, class DestAccessor> 00183 inline 00184 void createGaborFilter(triple<DestImageIterator, DestImageIterator, 00185 DestAccessor> dest, 00186 double orientation, double centerFrequency, 00187 double angularSigma, double radialSigma) 00188 { 00189 createGaborFilter(dest.first, dest.second, dest.third, 00190 orientation, centerFrequency, 00191 angularSigma, radialSigma); 00192 } 00193 00194 /********************************************************/ 00195 /* */ 00196 /* radialGaborSigma */ 00197 /* */ 00198 /********************************************************/ 00199 00200 /** \brief Calculate sensible radial sigma for given parameters. 00201 00202 For a brief introduction what is meant with "sensible" sigmas, see 00203 \ref angularGaborSigma(). 00204 00205 <b> Declaration:</b> 00206 00207 \code 00208 namespace vigra { 00209 inline 00210 double radialGaborSigma(double centerFrequency) 00211 } 00212 \endcode 00213 */ 00214 00215 inline double radialGaborSigma(double centerFrequency) 00216 { 00217 static double sfactor = 3.0 * VIGRA_CSTD::sqrt(VIGRA_CSTD::log(4.0)); 00218 return centerFrequency / sfactor; 00219 } 00220 00221 /********************************************************/ 00222 /* */ 00223 /* angularGaborSigma */ 00224 /* */ 00225 /********************************************************/ 00226 00227 /** \brief Calculate sensible angular sigma for given parameters. 00228 00229 "Sensible" means: If you use a range of gabor filters for feature 00230 detection, you are interested in minimal redundance. This is hard 00231 to define but one possible try is to arrange the filters in 00232 frequency space, so that the half-peak-magnitude ellipses touch 00233 each other. 00234 00235 To do so, you must know the number of directions (first parameter 00236 for the angular sigma function) and the center frequency of the 00237 filter you want to calculate the sigmas for. 00238 00239 The exact formulas are: 00240 \code 00241 sigma_radial= 1/sqrt(ln(4)) * centerFrequency/3 00242 \endcode 00243 00244 \code 00245 sigma_angular= 1/sqrt(ln(4)) * tan(pi/(directions*2)) 00246 * sqrt(8/9) * centerFrequency 00247 \endcode 00248 00249 <b> Declaration:</b> 00250 00251 \code 00252 namespace vigra { 00253 inline 00254 double angularGaborSigma(int directionCount, double centerFrequency) 00255 } 00256 \endcode 00257 */ 00258 00259 inline double angularGaborSigma(int directionCount, double centerFrequency) 00260 { 00261 return VIGRA_CSTD::tan(M_PI/directionCount/2.0) * centerFrequency 00262 * VIGRA_CSTD::sqrt(8.0 / (9 * VIGRA_CSTD::log(4.0))); 00263 } 00264 00265 /********************************************************/ 00266 /* */ 00267 /* GaborFilterFamily */ 00268 /* */ 00269 /********************************************************/ 00270 00271 /** \brief Family of gabor filters of different scale and direction. 00272 00273 A GaborFilterFamily can be used to quickly create a whole family 00274 of gabor filters in frequency space. Especially useful in 00275 conjunction with \ref applyFourierFilterFamily, since it's derived 00276 from \ref ImageArray. 00277 00278 The filter parameters are chosen to make the center frequencies 00279 decrease in octaves with increasing scale indices, and to make the 00280 half-peak-magnitude ellipses touch each other to somewhat reduce 00281 redundancy in the filter answers. This is done by using \ref 00282 angularGaborSigma() and \ref radialGaborSigma(), you'll find more 00283 information there. 00284 00285 The template parameter ImageType should be a scalar image type suitable for filling in 00286 00287 <b>\#include</b> "<a href="gaborfilter_8hxx-source.html">vigra/gaborfilter.hxx</a>" 00288 00289 Namespace: vigra 00290 */ 00291 template <class ImageType, 00292 class Alloc = typename ImageType::allocator_type::template rebind<ImageType>::other > 00293 class GaborFilterFamily 00294 : public ImageArray<ImageType, Alloc> 00295 { 00296 typedef ImageArray<ImageType, Alloc> ParentClass; 00297 int scaleCount_, directionCount_; 00298 double maxCenterFrequency_; 00299 00300 protected: 00301 void initFilters() 00302 { 00303 for(int direction= 0; direction<directionCount_; direction++) 00304 for(int scale= 0; scale<scaleCount_; scale++) 00305 { 00306 double angle = direction * M_PI / directionCount(); 00307 double centerFrequency = 00308 maxCenterFrequency_ / VIGRA_CSTD::pow(2.0, (double)scale); 00309 createGaborFilter(destImageRange(this->images_[filterIndex(direction, scale)]), 00310 angle, centerFrequency, 00311 angularGaborSigma(directionCount(), centerFrequency), 00312 radialGaborSigma(centerFrequency)); 00313 } 00314 } 00315 00316 public: 00317 enum { stdFilterSize= 128, stdDirectionCount= 6, stdScaleCount= 4 }; 00318 00319 /** Constructs a family of gabor filters in frequency 00320 space. The filters will be calculated on construction, so 00321 it makes sense to provide good parameters right now 00322 although they can be changed later, too. If you leave them 00323 out, the defaults are a \ref directionCount of 6, a \ref 00324 scaleCount of 4 and a \ref maxCenterFrequency of 00325 3/8(=0.375). 00326 */ 00327 GaborFilterFamily(const Diff2D & size, 00328 int directionCount = stdDirectionCount, int scaleCount = stdScaleCount, 00329 double maxCenterFrequency = 3.0/8.0, 00330 Alloc const & alloc = Alloc()) 00331 : ParentClass(directionCount*scaleCount, size, alloc), 00332 scaleCount_(scaleCount), 00333 directionCount_(directionCount), 00334 maxCenterFrequency_(maxCenterFrequency) 00335 { 00336 initFilters(); 00337 } 00338 00339 /** Convenience variant of the above constructor taking width 00340 and height separately. Also, this one serves as default 00341 constructor constructing 128x128 pixel filters. 00342 */ 00343 GaborFilterFamily(int width= stdFilterSize, int height= -1, 00344 int directionCount = stdDirectionCount, int scaleCount = stdScaleCount, 00345 double maxCenterFrequency = 3.0/8.0, 00346 Alloc const & alloc = Alloc()) 00347 : ParentClass(directionCount*scaleCount, 00348 Size2D(width, height > 0 ? height : width), alloc), 00349 scaleCount_(scaleCount), 00350 directionCount_(directionCount), 00351 maxCenterFrequency_(maxCenterFrequency) 00352 { 00353 initFilters(); 00354 } 00355 00356 /** Return the index of the filter with the given direction and 00357 scale in this ImageArray. direction must in the range 00358 0..directionCount()-1 and scale in the range 00359 0..rangeCount()-1. This is useful for example if you used 00360 \ref applyFourierFilterFamily() and got a resulting 00361 ImageArray which still has the same order of images, but no 00362 \ref getFilter() method anymore. 00363 */ 00364 int filterIndex(int direction, int scale) const 00365 { 00366 return scale*directionCount()+direction; 00367 } 00368 00369 /** Return the filter with the given direction and 00370 scale. direction must in the range 0..directionCount()-1 00371 and scale in the range 0..rangeCount()-1. 00372 <tt>filters.getFilter(direction, scale)</tt> is the same as 00373 <tt>filters[filterIndex(direction, scale)]</tt>. 00374 */ 00375 ImageType const & getFilter(int direction, int scale) const 00376 { 00377 return this->images_[filterIndex(direction, scale)]; 00378 } 00379 00380 /** Resize all filters (causing their recalculation). 00381 */ 00382 virtual void resizeImages(const Diff2D &newSize) 00383 { 00384 ParentClass::resizeImages(newSize); 00385 initFilters(); 00386 } 00387 00388 /** Query the number of filter scales available. 00389 */ 00390 int scaleCount() const 00391 { return scaleCount_; } 00392 00393 /** Query the number of filter directions available. 00394 */ 00395 int directionCount() const 00396 { return directionCount_; } 00397 00398 /** Change the number of directions / scales. This causes the 00399 recalculation of all filters. 00400 */ 00401 void setDirectionScaleCounts(int directionCount, int scaleCount) 00402 { 00403 this->resize(directionCount * scaleCount); 00404 scaleCount_ = scaleCount; 00405 directionCount_ = directionCount; 00406 initFilters(); 00407 } 00408 00409 /** Return the center frequency of the filter(s) with 00410 scale==0. Filters with scale>0 will have a center frequency 00411 reduced in octaves: 00412 <tt>centerFrequency= maxCenterFrequency / 2.0^scale</tt> 00413 */ 00414 double maxCenterFrequency() 00415 { return maxCenterFrequency_; } 00416 00417 /** Change the center frequency of the filter(s) with 00418 scale==0. See \ref maxCenterFrequency(). 00419 */ 00420 void setMaxCenterFrequency(double maxCenterFrequency) 00421 { 00422 maxCenterFrequency_ = maxCenterFrequency; 00423 initFilters(); 00424 } 00425 }; 00426 00427 //@} 00428 00429 } // namespace vigra 00430 00431 #endif // VIGRA_GABORFILTER_HXX
© Ullrich Köthe (koethe@informatik.uni-hamburg.de) |
html generated using doxygen and Python
|