[ VIGRA Homepage | Class Index | Function Index | File Index | Main Page ]
![]() |
vigra/fftw3.hxx | ![]() |
---|
00001 /************************************************************************/ 00002 /* */ 00003 /* Copyright 1998-2004 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.3, Aug 18 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 #ifndef VIGRA_FFTW3_HXX 00024 #define VIGRA_FFTW3_HXX 00025 00026 #include <cmath> 00027 #include <functional> 00028 #include "vigra/stdimage.hxx" 00029 #include "vigra/copyimage.hxx" 00030 #include "vigra/transformimage.hxx" 00031 #include "vigra/combineimages.hxx" 00032 #include "vigra/numerictraits.hxx" 00033 #include "vigra/imagecontainer.hxx" 00034 #include <fftw3.h> 00035 00036 namespace vigra { 00037 00038 typedef double fftw_real; 00039 00040 /********************************************************/ 00041 /* */ 00042 /* FFTWComplex */ 00043 /* */ 00044 /********************************************************/ 00045 00046 /** \brief Wrapper class for the FFTW type '<TT>fftw_complex</TT>'. 00047 00048 This class provides constructors and other member functions 00049 for the C struct '<TT>fftw_complex</TT>'. This struct is the basic 00050 pixel type of the <a href="http://www.fftw.org/">FFTW Fast Fourier Transform</a> 00051 library. It inherits the data members '<TT>re</TT>' and '<TT>im</TT>' 00052 that denote the real and imaginary part of the number. In addition it 00053 defines transformations to polar coordinates, 00054 as well as \ref FFTWComplexOperators "arithmetic operators" and 00055 \ref FFTWComplexAccessors "accessors". 00056 00057 FFTWComplex implements the concepts \ref AlgebraicField and 00058 \ref DivisionAlgebra. The standard image types <tt>FFTWRealImage</tt> 00059 and <tt>FFTWComplexImage</tt> are defined. 00060 00061 <b>See also:</b> 00062 <ul> 00063 <li> \ref FFTWComplexTraits 00064 <li> \ref FFTWComplexOperators 00065 <li> \ref FFTWComplexAccessors 00066 </ul> 00067 00068 <b>\#include</b> "<a href="fftw3_8hxx-source.html">vigra/fftw3.hxx</a>" (for FFTW 3) or<br> 00069 <b>\#include</b> "<a href="fftw_8hxx-source.html">vigra/fftw.hxx</a>" (for deprecated FFTW 2)<br> 00070 Namespace: vigra 00071 */ 00072 class FFTWComplex 00073 { 00074 fftw_complex data_; 00075 00076 public: 00077 /** The complex' component type, as defined in '<TT>fftw3.h</TT>' 00078 */ 00079 typedef fftw_real value_type; 00080 00081 /** reference type (result of operator[]) 00082 */ 00083 typedef fftw_real & reference; 00084 00085 /** const reference type (result of operator[] const) 00086 */ 00087 typedef fftw_real const & const_reference; 00088 00089 /** iterator type (result of begin() ) 00090 */ 00091 typedef fftw_real * iterator; 00092 00093 /** const iterator type (result of begin() const) 00094 */ 00095 typedef fftw_real const * const_iterator; 00096 00097 /** The norm type (result of magnitde()) 00098 */ 00099 typedef fftw_real NormType; 00100 00101 /** The squared norm type (result of squaredMagnitde()) 00102 */ 00103 typedef fftw_real SquaredNormType; 00104 00105 /** Construct from real and imaginary part. 00106 Default: 0. 00107 */ 00108 FFTWComplex(value_type const & re = 0.0, value_type const & im = 0.0) 00109 { 00110 data_[0] = re; 00111 data_[1] = im; 00112 } 00113 00114 /** Copy constructor. 00115 */ 00116 FFTWComplex(FFTWComplex const & o) 00117 { 00118 data_[0] = o.data_[0]; 00119 data_[1] = o.data_[1]; 00120 } 00121 00122 /** Construct from plain <TT>fftw_complex</TT>. 00123 */ 00124 FFTWComplex(fftw_complex const & o) 00125 { 00126 data_[0] = o[0]; 00127 data_[1] = o[1]; 00128 } 00129 00130 /** Construct from TinyVector. 00131 */ 00132 template <class T> 00133 FFTWComplex(TinyVector<T, 2> const & o) 00134 { 00135 data_[0] = o[0]; 00136 data_[1] = o[1]; 00137 } 00138 00139 /** Assignment. 00140 */ 00141 FFTWComplex& operator=(FFTWComplex const & o) 00142 { 00143 data_[0] = o.data_[0]; 00144 data_[1] = o.data_[1]; 00145 return *this; 00146 } 00147 00148 /** Assignment. 00149 */ 00150 FFTWComplex& operator=(fftw_complex const & o) 00151 { 00152 data_[0] = o[0]; 00153 data_[1] = o[1]; 00154 return *this; 00155 } 00156 00157 /** Assignment. 00158 */ 00159 FFTWComplex& operator=(fftw_real const & o) 00160 { 00161 data_[0] = o; 00162 data_[1] = 0.0; 00163 return *this; 00164 } 00165 00166 /** Assignment. 00167 */ 00168 template <class T> 00169 FFTWComplex& operator=(TinyVector<T, 2> const & o) 00170 { 00171 data_[0] = o[0]; 00172 data_[1] = o[1]; 00173 return *this; 00174 } 00175 00176 reference re() 00177 { return data_[0]; } 00178 00179 const_reference re() const 00180 { return data_[0]; } 00181 00182 reference im() 00183 { return data_[1]; } 00184 00185 const_reference im() const 00186 { return data_[1]; } 00187 00188 /** Unary negation. 00189 */ 00190 FFTWComplex operator-() const 00191 { return FFTWComplex(-data_[0], -data_[1]); } 00192 00193 /** Squared magnitude x*conj(x) 00194 */ 00195 SquaredNormType squaredMagnitude() const 00196 { return data_[0]*data_[0]+data_[1]*data_[1]; } 00197 00198 /** Magnitude (length of radius vector). 00199 */ 00200 NormType magnitude() const 00201 { return VIGRA_CSTD::sqrt(squaredMagnitude()); } 00202 00203 /** Phase angle. 00204 */ 00205 value_type phase() const 00206 { return VIGRA_CSTD::atan2(data_[1], data_[0]); } 00207 00208 /** Access components as if number were a vector. 00209 */ 00210 reference operator[](int i) 00211 { return data_[i]; } 00212 00213 /** Read components as if number were a vector. 00214 */ 00215 const_reference operator[](int i) const 00216 { return data_[i]; } 00217 00218 /** Length of complex number (always 2). 00219 */ 00220 int size() const 00221 { return 2; } 00222 00223 iterator begin() 00224 { return data_; } 00225 00226 iterator end() 00227 { return data_ + 2; } 00228 00229 const_iterator begin() const 00230 { return data_; } 00231 00232 const_iterator end() const 00233 { return data_ + 2; } 00234 }; 00235 00236 /********************************************************/ 00237 /* */ 00238 /* FFTWComplexTraits */ 00239 /* */ 00240 /********************************************************/ 00241 00242 /** \page FFTWComplexTraits Numeric and Promote Traits of FFTWComplex 00243 00244 The numeric and promote traits for fftw_complex and FFTWComplex follow 00245 the general specifications for \ref NumericPromotionTraits and 00246 \ref AlgebraicField. They are explicitly specialized for the types 00247 involved: 00248 00249 \code 00250 00251 template<> 00252 struct NumericTraits<fftw_complex> 00253 { 00254 typedef fftw_complex Promote; 00255 typedef fftw_complex RealPromote; 00256 typedef fftw_complex ComplexPromote; 00257 typedef fftw_real ValueType; 00258 00259 typedef VigraFalseType isIntegral; 00260 typedef VigraFalseType isScalar; 00261 typedef VigraFalseType isOrdered; 00262 typedef VigraTrueType isComplex; 00263 00264 // etc. 00265 }; 00266 00267 template<> 00268 struct NumericTraits<FFTWComplex> 00269 { 00270 typedef FFTWComplex Promote; 00271 typedef FFTWComplex RealPromote; 00272 typedef FFTWComplex ComplexPromote; 00273 typedef fftw_real ValueType; 00274 00275 typedef VigraFalseType isIntegral; 00276 typedef VigraFalseType isScalar; 00277 typedef VigraFalseType isOrdered; 00278 typedef VigraTrueType isComplex; 00279 00280 // etc. 00281 }; 00282 00283 template<> 00284 struct NormTraits<fftw_complex> 00285 { 00286 typedef fftw_complex Type; 00287 typedef fftw_real SquaredNormType; 00288 typedef fftw_real NormType; 00289 }; 00290 00291 template<> 00292 struct NormTraits<FFTWComplex> 00293 { 00294 typedef FFTWComplex Type; 00295 typedef fftw_real SquaredNormType; 00296 typedef fftw_real NormType; 00297 }; 00298 00299 template <> 00300 struct PromoteTraits<fftw_complex, fftw_complex> 00301 { 00302 typedef fftw_complex Promote; 00303 }; 00304 00305 template <> 00306 struct PromoteTraits<fftw_complex, double> 00307 { 00308 typedef fftw_complex Promote; 00309 }; 00310 00311 template <> 00312 struct PromoteTraits<double, fftw_complex> 00313 { 00314 typedef fftw_complex Promote; 00315 }; 00316 00317 template <> 00318 struct PromoteTraits<FFTWComplex, FFTWComplex> 00319 { 00320 typedef FFTWComplex Promote; 00321 }; 00322 00323 template <> 00324 struct PromoteTraits<FFTWComplex, double> 00325 { 00326 typedef FFTWComplex Promote; 00327 }; 00328 00329 template <> 00330 struct PromoteTraits<double, FFTWComplex> 00331 { 00332 typedef FFTWComplex Promote; 00333 }; 00334 \endcode 00335 00336 <b>\#include</b> "<a href="fftw3_8hxx-source.html">vigra/fftw3.hxx</a>" (for FFTW 3) or<br> 00337 <b>\#include</b> "<a href="fftw_8hxx-source.html">vigra/fftw.hxx</a>" (for deprecated FFTW 2)<br> 00338 Namespace: vigra 00339 00340 */ 00341 template<> 00342 struct NumericTraits<fftw_complex> 00343 { 00344 typedef fftw_complex Type; 00345 typedef fftw_complex Promote; 00346 typedef fftw_complex RealPromote; 00347 typedef fftw_complex ComplexPromote; 00348 typedef fftw_real ValueType; 00349 00350 typedef VigraFalseType isIntegral; 00351 typedef VigraFalseType isScalar; 00352 typedef VigraFalseType isOrdered; 00353 typedef VigraTrueType isComplex; 00354 00355 static FFTWComplex zero() { return FFTWComplex(0.0, 0.0); } 00356 static FFTWComplex one() { return FFTWComplex(1.0, 0.0); } 00357 static FFTWComplex nonZero() { return one(); } 00358 00359 static const Promote & toPromote(const Type & v) { return v; } 00360 static const RealPromote & toRealPromote(const Type & v) { return v; } 00361 static const Type & fromPromote(const Promote & v) { return v; } 00362 static const Type & fromRealPromote(const RealPromote & v) { return v; } 00363 }; 00364 00365 template<> 00366 struct NumericTraits<FFTWComplex> 00367 { 00368 typedef FFTWComplex Type; 00369 typedef FFTWComplex Promote; 00370 typedef FFTWComplex RealPromote; 00371 typedef FFTWComplex ComplexPromote; 00372 typedef fftw_real ValueType; 00373 00374 typedef VigraFalseType isIntegral; 00375 typedef VigraFalseType isScalar; 00376 typedef VigraFalseType isOrdered; 00377 typedef VigraTrueType isComplex; 00378 00379 static FFTWComplex zero() { return FFTWComplex(0.0, 0.0); } 00380 static FFTWComplex one() { return FFTWComplex(1.0, 0.0); } 00381 static FFTWComplex nonZero() { return one(); } 00382 00383 static const Promote & toPromote(const Type & v) { return v; } 00384 static const RealPromote & toRealPromote(const Type & v) { return v; } 00385 static const Type & fromPromote(const Promote & v) { return v; } 00386 static const Type & fromRealPromote(const RealPromote & v) { return v; } 00387 }; 00388 00389 template<> 00390 struct NormTraits<fftw_complex> 00391 { 00392 typedef fftw_complex Type; 00393 typedef fftw_real SquaredNormType; 00394 typedef fftw_real NormType; 00395 }; 00396 00397 template<> 00398 struct NormTraits<FFTWComplex> 00399 { 00400 typedef FFTWComplex Type; 00401 typedef fftw_real SquaredNormType; 00402 typedef fftw_real NormType; 00403 }; 00404 00405 template <> 00406 struct PromoteTraits<fftw_complex, fftw_complex> 00407 { 00408 typedef fftw_complex Promote; 00409 }; 00410 00411 template <> 00412 struct PromoteTraits<fftw_complex, double> 00413 { 00414 typedef fftw_complex Promote; 00415 }; 00416 00417 template <> 00418 struct PromoteTraits<double, fftw_complex> 00419 { 00420 typedef fftw_complex Promote; 00421 }; 00422 00423 template <> 00424 struct PromoteTraits<FFTWComplex, FFTWComplex> 00425 { 00426 typedef FFTWComplex Promote; 00427 }; 00428 00429 template <> 00430 struct PromoteTraits<FFTWComplex, double> 00431 { 00432 typedef FFTWComplex Promote; 00433 }; 00434 00435 template <> 00436 struct PromoteTraits<double, FFTWComplex> 00437 { 00438 typedef FFTWComplex Promote; 00439 }; 00440 00441 00442 /********************************************************/ 00443 /* */ 00444 /* FFTWComplex Operations */ 00445 /* */ 00446 /********************************************************/ 00447 00448 /** \addtogroup FFTWComplexOperators Functions for FFTWComplex 00449 00450 <b>\#include</b> "<a href="fftw3_8hxx-source.html">vigra/fftw3.hxx</a>" (for FFTW 3) or<br> 00451 <b>\#include</b> "<a href="fftw_8hxx-source.html">vigra/fftw.hxx</a>" (for deprecated FFTW 2)<br> 00452 00453 These functions fulfill the requirements of an Algebraic Field. 00454 Return types are determined according to \ref FFTWComplexTraits. 00455 00456 Namespace: vigra 00457 <p> 00458 00459 */ 00460 //@{ 00461 /// equal 00462 inline bool operator ==(FFTWComplex const &a, const FFTWComplex &b) { 00463 return a.re() == b.re() && a.im() == b.im(); 00464 } 00465 00466 /// not equal 00467 inline bool operator !=(FFTWComplex const &a, const FFTWComplex &b) { 00468 return a.re() != b.re() || a.im() != b.im(); 00469 } 00470 00471 /// add-assignment 00472 inline FFTWComplex &operator +=(FFTWComplex &a, const FFTWComplex &b) { 00473 a.re() += b.re(); 00474 a.im() += b.im(); 00475 return a; 00476 } 00477 00478 /// subtract-assignment 00479 inline FFTWComplex &operator -=(FFTWComplex &a, const FFTWComplex &b) { 00480 a.re() -= b.re(); 00481 a.im() -= b.im(); 00482 return a; 00483 } 00484 00485 /// multiply-assignment 00486 inline FFTWComplex &operator *=(FFTWComplex &a, const FFTWComplex &b) { 00487 FFTWComplex::value_type t = a.re()*b.re()-a.im()*b.im(); 00488 a.im() = a.re()*b.im()+a.im()*b.re(); 00489 a.re() = t; 00490 return a; 00491 } 00492 00493 /// divide-assignment 00494 inline FFTWComplex &operator /=(FFTWComplex &a, const FFTWComplex &b) { 00495 FFTWComplex::value_type sm = b.squaredMagnitude(); 00496 FFTWComplex::value_type t = (a.re()*b.re()+a.im()*b.im())/sm; 00497 a.im() = (b.re()*a.im()-a.re()*b.im())/sm; 00498 a.re() = t; 00499 return a; 00500 } 00501 00502 /// multiply-assignment with scalar double 00503 inline FFTWComplex &operator *=(FFTWComplex &a, const double &b) { 00504 a.re() *= b; 00505 a.im() *= b; 00506 return a; 00507 } 00508 00509 /// divide-assignment with scalar double 00510 inline FFTWComplex &operator /=(FFTWComplex &a, const double &b) { 00511 a.re() /= b; 00512 a.im() /= b; 00513 return a; 00514 } 00515 00516 /// addition 00517 inline FFTWComplex operator +(FFTWComplex a, const FFTWComplex &b) { 00518 a += b; 00519 return a; 00520 } 00521 00522 /// subtraction 00523 inline FFTWComplex operator -(FFTWComplex a, const FFTWComplex &b) { 00524 a -= b; 00525 return a; 00526 } 00527 00528 /// multiplication 00529 inline FFTWComplex operator *(FFTWComplex a, const FFTWComplex &b) { 00530 a *= b; 00531 return a; 00532 } 00533 00534 /// right multiplication with scalar double 00535 inline FFTWComplex operator *(FFTWComplex a, const double &b) { 00536 a *= b; 00537 return a; 00538 } 00539 00540 /// left multiplication with scalar double 00541 inline FFTWComplex operator *(const double &a, FFTWComplex b) { 00542 b *= a; 00543 return b; 00544 } 00545 00546 /// division 00547 inline FFTWComplex operator /(FFTWComplex a, const FFTWComplex &b) { 00548 a /= b; 00549 return a; 00550 } 00551 00552 /// right division with scalar double 00553 inline FFTWComplex operator /(FFTWComplex a, const double &b) { 00554 a /= b; 00555 return a; 00556 } 00557 00558 using VIGRA_CSTD::abs; 00559 00560 /// absolute value (= magnitude) 00561 inline FFTWComplex::value_type abs(const FFTWComplex &a) 00562 { 00563 return a.magnitude(); 00564 } 00565 00566 /// complex conjugate 00567 inline FFTWComplex conj(const FFTWComplex &a) 00568 { 00569 return FFTWComplex(a.re(), -a.im()); 00570 } 00571 00572 /// norm (= magnitude) 00573 inline FFTWComplex::NormType norm(const FFTWComplex &a) 00574 { 00575 return a.magnitude(); 00576 } 00577 00578 /// squared norm (= squared magnitude) 00579 inline FFTWComplex::SquaredNormType squaredNorm(const FFTWComplex &a) 00580 { 00581 return a.squaredMagnitude(); 00582 } 00583 00584 //@} 00585 00586 00587 /** \addtogroup StandardImageTypes 00588 */ 00589 //@{ 00590 00591 /********************************************************/ 00592 /* */ 00593 /* FFTWRealImage */ 00594 /* */ 00595 /********************************************************/ 00596 00597 /** Float (<tt>fftw_real</tt>) image. 00598 00599 The type <tt>fftw_real</tt> is defined as <tt>double</tt> (in FFTW 2 it used to be 00600 either <tt>float</tt> or <tt>double</tt>, as specified during compilation of FFTW). 00601 FFTWRealImage uses \ref vigra::BasicImageIterator and \ref vigra::StandardAccessor and 00602 their const counterparts to access the data. 00603 00604 <b>\#include</b> "<a href="fftw3_8hxx-source.html">vigra/fftw3.hxx</a>" (for FFTW 3) or<br> 00605 <b>\#include</b> "<a href="fftw_8hxx-source.html">vigra/fftw.hxx</a>" (for deprecated FFTW 2)<br> 00606 Namespace: vigra 00607 */ 00608 typedef BasicImage<fftw_real> FFTWRealImage; 00609 00610 /********************************************************/ 00611 /* */ 00612 /* FFTWComplexImage */ 00613 /* */ 00614 /********************************************************/ 00615 00616 template<> 00617 struct IteratorTraits< 00618 BasicImageIterator<FFTWComplex, FFTWComplex **> > 00619 { 00620 typedef BasicImageIterator<FFTWComplex, FFTWComplex **> Iterator; 00621 typedef Iterator iterator; 00622 typedef iterator::iterator_category iterator_category; 00623 typedef iterator::value_type value_type; 00624 typedef iterator::reference reference; 00625 typedef iterator::index_reference index_reference; 00626 typedef iterator::pointer pointer; 00627 typedef iterator::difference_type difference_type; 00628 typedef iterator::row_iterator row_iterator; 00629 typedef iterator::column_iterator column_iterator; 00630 typedef VectorAccessor<FFTWComplex> default_accessor; 00631 typedef VectorAccessor<FFTWComplex> DefaultAccessor; 00632 typedef VigraTrueType hasConstantStrides; 00633 }; 00634 00635 template<> 00636 struct IteratorTraits< 00637 ConstBasicImageIterator<FFTWComplex, FFTWComplex **> > 00638 { 00639 typedef ConstBasicImageIterator<FFTWComplex, FFTWComplex **> Iterator; 00640 typedef Iterator iterator; 00641 typedef iterator::iterator_category iterator_category; 00642 typedef iterator::value_type value_type; 00643 typedef iterator::reference reference; 00644 typedef iterator::index_reference index_reference; 00645 typedef iterator::pointer pointer; 00646 typedef iterator::difference_type difference_type; 00647 typedef iterator::row_iterator row_iterator; 00648 typedef iterator::column_iterator column_iterator; 00649 typedef VectorAccessor<FFTWComplex> default_accessor; 00650 typedef VectorAccessor<FFTWComplex> DefaultAccessor; 00651 typedef VigraTrueType hasConstantStrides; 00652 }; 00653 00654 /** Complex (FFTWComplex) image. 00655 It uses \ref vigra::BasicImageIterator and \ref vigra::StandardAccessor and 00656 their const counterparts to access the data. 00657 00658 <b>\#include</b> "<a href="fftw3_8hxx-source.html">vigra/fftw3.hxx</a>" (for FFTW 3) or<br> 00659 <b>\#include</b> "<a href="fftw_8hxx-source.html">vigra/fftw.hxx</a>" (for deprecated FFTW 2)<br> 00660 Namespace: vigra 00661 */ 00662 typedef BasicImage<FFTWComplex> FFTWComplexImage; 00663 00664 //@} 00665 00666 /********************************************************/ 00667 /* */ 00668 /* FFTWComplex-Accessors */ 00669 /* */ 00670 /********************************************************/ 00671 00672 /** \addtogroup DataAccessors 00673 */ 00674 //@{ 00675 /** \defgroup FFTWComplexAccessors Accessors for FFTWComplex 00676 00677 Encapsulate access to pixels of type FFTWComplex 00678 */ 00679 //@{ 00680 /** Encapsulate access to the the real part of a complex number. 00681 00682 <b>\#include</b> "<a href="fftw3_8hxx-source.html">vigra/fftw3.hxx</a>" (for FFTW 3) or<br> 00683 <b>\#include</b> "<a href="fftw_8hxx-source.html">vigra/fftw.hxx</a>" (for deprecated FFTW 2)<br> 00684 Namespace: vigra 00685 */ 00686 class FFTWRealAccessor 00687 { 00688 public: 00689 00690 /// The accessor's value type. 00691 typedef fftw_real value_type; 00692 00693 /// Read real part at iterator position. 00694 template <class ITERATOR> 00695 value_type operator()(ITERATOR const & i) const { 00696 return (*i).re(); 00697 } 00698 00699 /// Read real part at offset from iterator position. 00700 template <class ITERATOR, class DIFFERENCE> 00701 value_type operator()(ITERATOR const & i, DIFFERENCE d) const { 00702 return i[d].re(); 00703 } 00704 00705 /// Write real part at iterator position. 00706 template <class ITERATOR> 00707 void set(value_type const & v, ITERATOR const & i) const { 00708 (*i).re()= v; 00709 } 00710 00711 /// Write real part at offset from iterator position. 00712 template <class ITERATOR, class DIFFERENCE> 00713 void set(value_type const & v, ITERATOR const & i, DIFFERENCE d) const { 00714 i[d].re()= v; 00715 } 00716 }; 00717 00718 /** Encapsulate access to the the imaginary part of a complex number. 00719 00720 <b>\#include</b> "<a href="fftw3_8hxx-source.html">vigra/fftw3.hxx</a>" (for FFTW 3) or<br> 00721 <b>\#include</b> "<a href="fftw_8hxx-source.html">vigra/fftw.hxx</a>" (for deprecated FFTW 2)<br> 00722 Namespace: vigra 00723 */ 00724 class FFTWImaginaryAccessor 00725 { 00726 public: 00727 /// The accessor's value type. 00728 typedef fftw_real value_type; 00729 00730 /// Read imaginary part at iterator position. 00731 template <class ITERATOR> 00732 value_type operator()(ITERATOR const & i) const { 00733 return (*i).im(); 00734 } 00735 00736 /// Read imaginary part at offset from iterator position. 00737 template <class ITERATOR, class DIFFERENCE> 00738 value_type operator()(ITERATOR const & i, DIFFERENCE d) const { 00739 return i[d].im(); 00740 } 00741 00742 /// Write imaginary part at iterator position. 00743 template <class ITERATOR> 00744 void set(value_type const & v, ITERATOR const & i) const { 00745 (*i).im()= v; 00746 } 00747 00748 /// Write imaginary part at offset from iterator position. 00749 template <class ITERATOR, class DIFFERENCE> 00750 void set(value_type const & v, ITERATOR const & i, DIFFERENCE d) const { 00751 i[d].im()= v; 00752 } 00753 }; 00754 00755 /** Write a real number into a complex one. The imaginary part is set 00756 to 0. 00757 00758 <b>\#include</b> "<a href="fftw3_8hxx-source.html">vigra/fftw3.hxx</a>" (for FFTW 3) or<br> 00759 <b>\#include</b> "<a href="fftw_8hxx-source.html">vigra/fftw.hxx</a>" (for deprecated FFTW 2)<br> 00760 Namespace: vigra 00761 */ 00762 class FFTWWriteRealAccessor: public FFTWRealAccessor 00763 { 00764 public: 00765 /// The accessor's value type. 00766 typedef fftw_real value_type; 00767 00768 /** Write real number at iterator position. Set imaginary part 00769 to 0. 00770 */ 00771 template <class ITERATOR> 00772 void set(value_type const & v, ITERATOR const & i) const { 00773 (*i).re()= v; 00774 (*i).im()= 0; 00775 } 00776 00777 /** Write real number at offset from iterator position. Set imaginary part 00778 to 0. 00779 */ 00780 template <class ITERATOR, class DIFFERENCE> 00781 void set(value_type const & v, ITERATOR const & i, DIFFERENCE d) const { 00782 i[d].re()= v; 00783 i[d].im()= 0; 00784 } 00785 }; 00786 00787 /** Calculate magnitude of complex number on the fly. 00788 00789 <b>\#include</b> "<a href="fftw3_8hxx-source.html">vigra/fftw3.hxx</a>" (for FFTW 3) or<br> 00790 <b>\#include</b> "<a href="fftw_8hxx-source.html">vigra/fftw.hxx</a>" (for deprecated FFTW 2)<br> 00791 Namespace: vigra 00792 */ 00793 class FFTWMagnitudeAccessor 00794 { 00795 public: 00796 /// The accessor's value type. 00797 typedef fftw_real value_type; 00798 00799 /// Read magnitude at iterator position. 00800 template <class ITERATOR> 00801 value_type operator()(ITERATOR const & i) const { 00802 return (*i).magnitude(); 00803 } 00804 00805 /// Read magnitude at offset from iterator position. 00806 template <class ITERATOR, class DIFFERENCE> 00807 value_type operator()(ITERATOR const & i, DIFFERENCE d) const { 00808 return (i[d]).magnitude(); 00809 } 00810 }; 00811 00812 /** Calculate phase of complex number on the fly. 00813 00814 <b>\#include</b> "<a href="fftw3_8hxx-source.html">vigra/fftw3.hxx</a>" (for FFTW 3) or<br> 00815 <b>\#include</b> "<a href="fftw_8hxx-source.html">vigra/fftw.hxx</a>" (for deprecated FFTW 2)<br> 00816 Namespace: vigra 00817 */ 00818 class FFTWPhaseAccessor 00819 { 00820 public: 00821 /// The accessor's value type. 00822 typedef fftw_real value_type; 00823 00824 /// Read phase at iterator position. 00825 template <class ITERATOR> 00826 value_type operator()(ITERATOR const & i) const { 00827 return (*i).phase(); 00828 } 00829 00830 /// Read phase at offset from iterator position. 00831 template <class ITERATOR, class DIFFERENCE> 00832 value_type operator()(ITERATOR const & i, DIFFERENCE d) const { 00833 return (i[d]).phase(); 00834 } 00835 }; 00836 00837 //@} 00838 //@} 00839 00840 /********************************************************/ 00841 /* */ 00842 /* Fourier Transform */ 00843 /* */ 00844 /********************************************************/ 00845 00846 /** \addtogroup FourierTransform Fast Fourier Transform 00847 00848 This documentation describes the VIGRA interface to FFTW 3. There is also an 00849 \link FourierTransformFFTW2 interface to the older version FFTW 2\endlink 00850 00851 VIGRA uses the <a href="http://www.fftw.org/">FFTW Fast Fourier 00852 Transform</a> package to perform Fourier transformations. VIGRA 00853 provides a wrapper for FFTW's complex number type (FFTWComplex), 00854 but FFTW's functions are used verbatim. If the image is stored as 00855 a FFTWComplexImage, the simplest call to an FFT function is like this: 00856 00857 \code 00858 vigra::FFTWComplexImage spatial(width,height), fourier(width,height); 00859 ... // fill image with data 00860 00861 // create a plan with estimated performance optimization 00862 fftw_plan forwardPlan = fftw_plan_dft_2d(height, width, 00863 (fftw_complex *)spatial.begin(), (fftw_complex *)fourier.begin(), 00864 FFTW_FORWARD, FFTW_ESTIMATE ); 00865 // calculate FFT (this can be repeated as often as needed, 00866 // with fresh data written into the source array) 00867 fftw_execute(forwardPlan); 00868 00869 // release the plan memory 00870 fftw_destroy_plan(forwardPlan); 00871 00872 // likewise for the inverse transform 00873 fftw_plan backwardPlan = fftw_plan_dft_2d(height, width, 00874 (fftw_complex *)fourier.begin(), (fftw_complex *)spatial.begin(), 00875 FFTW_BACKWARD, FFTW_ESTIMATE); 00876 fftw_execute(backwardPlan); 00877 fftw_destroy_plan(backwardPlan); 00878 00879 // do not forget to normalize the result according to the image size 00880 transformImage(srcImageRange(spatial), destImage(spatial), 00881 std::bind1st(std::multiplies<FFTWComplex>(), 1.0 / width / height)); 00882 \endcode 00883 00884 Note that in the creation of a plan, the height must be given 00885 first. Note also that <TT>spatial.begin()</TT> may only be passed 00886 to <TT>fftw_plan_dft_2d</TT> if the transform shall be applied to the 00887 entire image. When you want to restrict operation to an ROI, you 00888 can create a copy of the ROI in an image of appropriate size, or 00889 you may use the Guru interface to FFTW. 00890 00891 More information on using FFTW can be found <a href="http://www.fftw.org/doc/">here</a>. 00892 00893 FFTW produces fourier images that have the DC component (the 00894 origin of the Fourier space) in the upper left corner. Often, one 00895 wants the origin in the center of the image, so that frequencies 00896 always increase towards the border of the image. This can be 00897 achieved by calling \ref moveDCToCenter(). The inverse 00898 transformation is done by \ref moveDCToUpperLeft(). 00899 00900 <b>\#include</b> "<a href="fftw3_8hxx-source.html">vigra/fftw3.hxx</a>"<br> 00901 Namespace: vigra 00902 */ 00903 00904 /** \addtogroup FourierTransform 00905 */ 00906 //@{ 00907 00908 /********************************************************/ 00909 /* */ 00910 /* moveDCToCenter */ 00911 /* */ 00912 /********************************************************/ 00913 00914 /** \brief Rearrange the quadrants of a Fourier image so that the origin is 00915 in the image center. 00916 00917 FFTW produces fourier images where the DC component (origin of 00918 fourier space) is located in the upper left corner of the 00919 image. The quadrants are placed like this (using a 4x4 image for 00920 example): 00921 00922 \code 00923 DC 4 3 3 00924 4 4 3 3 00925 1 1 2 2 00926 1 1 2 2 00927 \endcode 00928 00929 After applying the function, the quadrants are at their usual places: 00930 00931 \code 00932 2 2 1 1 00933 2 2 1 1 00934 3 3 DC 4 00935 3 3 4 4 00936 \endcode 00937 00938 This transformation can be reversed by \ref moveDCToUpperLeft(). 00939 Note that the transformation must not be executed in place - input 00940 and output images must be different. 00941 00942 <b> Declarations:</b> 00943 00944 pass arguments explicitly: 00945 \code 00946 namespace vigra { 00947 template <class SrcImageIterator, class SrcAccessor, 00948 class DestImageIterator, class DestAccessor> 00949 void moveDCToCenter(SrcImageIterator src_upperleft, 00950 SrcImageIterator src_lowerright, SrcAccessor sa, 00951 DestImageIterator dest_upperleft, DestAccessor da); 00952 } 00953 \endcode 00954 00955 00956 use argument objects in conjunction with \ref ArgumentObjectFactories: 00957 \code 00958 namespace vigra { 00959 template <class SrcImageIterator, class SrcAccessor, 00960 class DestImageIterator, class DestAccessor> 00961 void moveDCToCenter( 00962 triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src, 00963 pair<DestImageIterator, DestAccessor> dest); 00964 } 00965 \endcode 00966 00967 <b> Usage:</b> 00968 00969 <b>\#include</b> "<a href="fftw3_8hxx-source.html">vigra/fftw3.hxx</a>"<br> 00970 Namespace: vigra 00971 00972 \code 00973 vigra::FFTWComplexImage spatial(width,height), fourier(width,height); 00974 ... // fill image with data 00975 00976 // create a plan with estimated performance optimization 00977 fftw_plan forwardPlan = fftw_plan_dft_2d(height, width, 00978 (fftw_complex *)spatial.begin(), (fftw_complex *)fourier.begin(), 00979 FFTW_FORWARD, FFTW_ESTIMATE ); 00980 // calculate FFT 00981 fftw_execute(forwardPlan); 00982 00983 vigra::FFTWComplexImage rearrangedFourier(width, height); 00984 moveDCToCenter(srcImageRange(fourier), destImage(rearrangedFourier)); 00985 00986 // delete the plan 00987 fftw_destroy_plan(forwardPlan); 00988 \endcode 00989 */ 00990 template <class SrcImageIterator, class SrcAccessor, 00991 class DestImageIterator, class DestAccessor> 00992 void moveDCToCenter(SrcImageIterator src_upperleft, 00993 SrcImageIterator src_lowerright, SrcAccessor sa, 00994 DestImageIterator dest_upperleft, DestAccessor da) 00995 { 00996 int w= src_lowerright.x - src_upperleft.x; 00997 int h= src_lowerright.y - src_upperleft.y; 00998 int w1 = w/2; 00999 int h1 = h/2; 01000 int w2 = (w+1)/2; 01001 int h2 = (h+1)/2; 01002 01003 // 2. Quadrant zum 4. 01004 copyImage(srcIterRange(src_upperleft, 01005 src_upperleft + Diff2D(w2, h2), sa), 01006 destIter (dest_upperleft + Diff2D(w1, h1), da)); 01007 01008 // 4. Quadrant zum 2. 01009 copyImage(srcIterRange(src_upperleft + Diff2D(w2, h2), 01010 src_lowerright, sa), 01011 destIter (dest_upperleft, da)); 01012 01013 // 1. Quadrant zum 3. 01014 copyImage(srcIterRange(src_upperleft + Diff2D(w2, 0), 01015 src_upperleft + Diff2D(w, h2), sa), 01016 destIter (dest_upperleft + Diff2D(0, h1), da)); 01017 01018 // 3. Quadrant zum 1. 01019 copyImage(srcIterRange(src_upperleft + Diff2D(0, h2), 01020 src_upperleft + Diff2D(w2, h), sa), 01021 destIter (dest_upperleft + Diff2D(w1, 0), da)); 01022 } 01023 01024 template <class SrcImageIterator, class SrcAccessor, 01025 class DestImageIterator, class DestAccessor> 01026 inline void moveDCToCenter( 01027 triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src, 01028 pair<DestImageIterator, DestAccessor> dest) 01029 { 01030 moveDCToCenter(src.first, src.second, src.third, 01031 dest.first, dest.second); 01032 } 01033 01034 /********************************************************/ 01035 /* */ 01036 /* moveDCToUpperLeft */ 01037 /* */ 01038 /********************************************************/ 01039 01040 /** \brief Rearrange the quadrants of a Fourier image so that the origin is 01041 in the image's upper left. 01042 01043 This function is the inversion of \ref moveDCToCenter(). See there 01044 for declarations and a usage example. 01045 01046 <b> Declarations:</b> 01047 01048 pass arguments explicitly: 01049 \code 01050 namespace vigra { 01051 template <class SrcImageIterator, class SrcAccessor, 01052 class DestImageIterator, class DestAccessor> 01053 void moveDCToUpperLeft(SrcImageIterator src_upperleft, 01054 SrcImageIterator src_lowerright, SrcAccessor sa, 01055 DestImageIterator dest_upperleft, DestAccessor da); 01056 } 01057 \endcode 01058 01059 01060 use argument objects in conjunction with \ref ArgumentObjectFactories: 01061 \code 01062 namespace vigra { 01063 template <class SrcImageIterator, class SrcAccessor, 01064 class DestImageIterator, class DestAccessor> 01065 void moveDCToUpperLeft( 01066 triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src, 01067 pair<DestImageIterator, DestAccessor> dest); 01068 } 01069 \endcode 01070 */ 01071 template <class SrcImageIterator, class SrcAccessor, 01072 class DestImageIterator, class DestAccessor> 01073 void moveDCToUpperLeft(SrcImageIterator src_upperleft, 01074 SrcImageIterator src_lowerright, SrcAccessor sa, 01075 DestImageIterator dest_upperleft, DestAccessor da) 01076 { 01077 int w= src_lowerright.x - src_upperleft.x; 01078 int h= src_lowerright.y - src_upperleft.y; 01079 int w2 = w/2; 01080 int h2 = h/2; 01081 int w1 = (w+1)/2; 01082 int h1 = (h+1)/2; 01083 01084 // 2. Quadrant zum 4. 01085 copyImage(srcIterRange(src_upperleft, 01086 src_upperleft + Diff2D(w2, h2), sa), 01087 destIter (dest_upperleft + Diff2D(w1, h1), da)); 01088 01089 // 4. Quadrant zum 2. 01090 copyImage(srcIterRange(src_upperleft + Diff2D(w2, h2), 01091 src_lowerright, sa), 01092 destIter (dest_upperleft, da)); 01093 01094 // 1. Quadrant zum 3. 01095 copyImage(srcIterRange(src_upperleft + Diff2D(w2, 0), 01096 src_upperleft + Diff2D(w, h2), sa), 01097 destIter (dest_upperleft + Diff2D(0, h1), da)); 01098 01099 // 3. Quadrant zum 1. 01100 copyImage(srcIterRange(src_upperleft + Diff2D(0, h2), 01101 src_upperleft + Diff2D(w2, h), sa), 01102 destIter (dest_upperleft + Diff2D(w1, 0), da)); 01103 } 01104 01105 template <class SrcImageIterator, class SrcAccessor, 01106 class DestImageIterator, class DestAccessor> 01107 inline void moveDCToUpperLeft( 01108 triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src, 01109 pair<DestImageIterator, DestAccessor> dest) 01110 { 01111 moveDCToUpperLeft(src.first, src.second, src.third, 01112 dest.first, dest.second); 01113 } 01114 01115 /********************************************************/ 01116 /* */ 01117 /* applyFourierFilter */ 01118 /* */ 01119 /********************************************************/ 01120 01121 /** \brief Apply a filter (defined in the frequency domain) to an image. 01122 01123 After transferring the image into the frequency domain, it is 01124 multiplied pixel-wise with the filter and transformed back. The 01125 result is always put into the given FFTWComplexImage 01126 <TT>destImg</TT> which must have the right size. Finally, the 01127 result will be normalized to compensate for the two FFTs. 01128 01129 The input and filter images can be scalar images (more precisely, 01130 <TT>SrcAccessor::value_type</TT> must be scalar) or 01131 FFTWComplexImages (in this case, one conversion is saved). 01132 01133 The DC entry of the filter must be in the upper left, which is the 01134 position where FFTW expects it (see \ref moveDCToUpperLeft()). 01135 01136 You can optionally pass two fftwnd_plans as last parameters, the 01137 forward plan and the in-place backward plan. Both must have been 01138 created for the right image size (see sample code). 01139 01140 <b> Declarations:</b> 01141 01142 pass arguments explicitly: 01143 \code 01144 namespace vigra { 01145 template <class SrcImageIterator, class SrcAccessor, 01146 class FilterImageIterator, class FilterAccessor> 01147 void applyFourierFilter(SrcImageIterator srcUpperLeft, 01148 SrcImageIterator srcLowerRight, SrcAccessor sa, 01149 FilterImageIterator filterUpperLeft, FilterAccessor fa, 01150 FFTWComplexImage & destImg) 01151 01152 template <class SrcImageIterator, class SrcAccessor, 01153 class FilterImageIterator, class FilterAccessor> 01154 void applyFourierFilter(SrcImageIterator srcUpperLeft, 01155 SrcImageIterator srcLowerRight, SrcAccessor sa, 01156 FilterImageIterator filterUpperLeft, FilterAccessor fa, 01157 FFTWComplexImage & destImg, 01158 const fftwnd_plan & forwardPlan, const fftwnd_plan & backwardPlan) 01159 } 01160 \endcode 01161 01162 use argument objects in conjunction with \ref ArgumentObjectFactories: 01163 \code 01164 namespace vigra { 01165 template <class SrcImageIterator, class SrcAccessor, 01166 class FilterImageIterator, class FilterAccessor> 01167 inline 01168 void applyFourierFilter(triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src, 01169 pair<FilterImageIterator, FilterAccessor> filter, 01170 FFTWComplexImage &destImg) 01171 01172 template <class SrcImageIterator, class SrcAccessor, 01173 class FilterImageIterator, class FilterAccessor> 01174 inline 01175 void applyFourierFilter(triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src, 01176 pair<FilterImageIterator, FilterAccessor> filter, 01177 FFTWComplexImage &destImg, 01178 const fftwnd_plan &forwardPlan, const fftwnd_plan &backwardPlan) 01179 } 01180 \endcode 01181 01182 <b> Usage:</b> 01183 01184 <b>\#include</b> "<a href="fftw3_8hxx-source.html">vigra/fftw3.hxx</a>"<br> 01185 Namespace: vigra 01186 01187 \code 01188 // create a Gaussian filter in Fourier space 01189 vigra::FImage gaussFilter(w, h), filter(w, h); 01190 for(int y=0; y<h; ++y) 01191 for(int x=0; x<w; ++x) 01192 { 01193 xx = float(x - w / 2) / w; 01194 yy = float(y - h / 2) / h; 01195 01196 gaussFilter(x,y) = std::exp(-(xx*xx + yy*yy) / 2.0 * scale); 01197 } 01198 01199 // applyFourierFilter() expects the filter's DC in the upper left 01200 moveDCToUpperLeft(srcImageRange(gaussFilter), destImage(filter)); 01201 01202 vigra::FFTWComplexImage result(w, h); 01203 01204 vigra::applyFourierFilter(srcImageRange(image), srcImage(filter), result); 01205 \endcode 01206 01207 For inspection of the result, \ref FFTWMagnitudeAccessor might be 01208 useful. If you want to apply the same filter repeatedly, it may be more 01209 efficient to use the FFTW functions directly with FFTW plans optimized 01210 for good performance. 01211 */ 01212 template <class SrcImageIterator, class SrcAccessor, 01213 class FilterImageIterator, class FilterAccessor, 01214 class DestImageIterator, class DestAccessor> 01215 inline 01216 void applyFourierFilter(triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src, 01217 pair<FilterImageIterator, FilterAccessor> filter, 01218 pair<DestImageIterator, DestAccessor> dest) 01219 { 01220 applyFourierFilter(src.first, src.second, src.third, 01221 filter.first, filter.second, 01222 dest.first, dest.second); 01223 } 01224 01225 template <class SrcImageIterator, class SrcAccessor, 01226 class FilterImageIterator, class FilterAccessor, 01227 class DestImageIterator, class DestAccessor> 01228 void applyFourierFilter(SrcImageIterator srcUpperLeft, 01229 SrcImageIterator srcLowerRight, SrcAccessor sa, 01230 FilterImageIterator filterUpperLeft, FilterAccessor fa, 01231 DestImageIterator destUpperLeft, DestAccessor da) 01232 { 01233 // copy real input images into a complex one... 01234 int w= srcLowerRight.x - srcUpperLeft.x; 01235 int h= srcLowerRight.y - srcUpperLeft.y; 01236 01237 FFTWComplexImage workImage(w, h); 01238 copyImage(srcIterRange(srcUpperLeft, srcLowerRight, sa), 01239 destImage(workImage, FFTWWriteRealAccessor())); 01240 01241 // ...and call the impl 01242 FFTWComplexImage const & cworkImage = workImage; 01243 applyFourierFilterImpl(cworkImage.upperLeft(), cworkImage.lowerRight(), cworkImage.accessor(), 01244 filterUpperLeft, fa, 01245 destUpperLeft, da); 01246 } 01247 01248 template <class FilterImageIterator, class FilterAccessor, 01249 class DestImageIterator, class DestAccessor> 01250 inline 01251 void applyFourierFilter( 01252 FFTWComplexImage::const_traverser srcUpperLeft, 01253 FFTWComplexImage::const_traverser srcLowerRight, 01254 FFTWComplexImage::ConstAccessor sa, 01255 FilterImageIterator filterUpperLeft, FilterAccessor fa, 01256 DestImageIterator destUpperLeft, DestAccessor da) 01257 { 01258 int w = srcLowerRight.x - srcUpperLeft.x; 01259 int h = srcLowerRight.y - srcUpperLeft.y; 01260 01261 // test for right memory layout (fftw expects a 2*width*height floats array) 01262 if (&(*(srcUpperLeft + Diff2D(w, 0))) == &(*(srcUpperLeft + Diff2D(0, 1)))) 01263 applyFourierFilterImpl(srcUpperLeft, srcLowerRight, sa, 01264 filterUpperLeft, fa, 01265 destUpperLeft, da); 01266 else 01267 { 01268 FFTWComplexImage workImage(w, h); 01269 copyImage(srcIterRange(srcUpperLeft, srcLowerRight, sa), 01270 destImage(workImage)); 01271 01272 FFTWComplexImage const & cworkImage = workImage; 01273 applyFourierFilterImpl(cworkImage.upperLeft(), cworkImage.lowerRight(), cworkImage.accessor(), 01274 filterUpperLeft, fa, 01275 destUpperLeft, da); 01276 } 01277 } 01278 01279 template <class FilterImageIterator, class FilterAccessor, 01280 class DestImageIterator, class DestAccessor> 01281 void applyFourierFilterImpl( 01282 FFTWComplexImage::const_traverser srcUpperLeft, 01283 FFTWComplexImage::const_traverser srcLowerRight, 01284 FFTWComplexImage::ConstAccessor sa, 01285 FilterImageIterator filterUpperLeft, FilterAccessor fa, 01286 DestImageIterator destUpperLeft, DestAccessor da) 01287 { 01288 int w = srcLowerRight.x - srcUpperLeft.x; 01289 int h = srcLowerRight.y - srcUpperLeft.y; 01290 01291 FFTWComplexImage complexResultImg(srcLowerRight - srcUpperLeft); 01292 01293 // FFT from srcImage to complexResultImg 01294 fftw_plan forwardPlan= 01295 fftw_plan_dft_2d(h, w, (fftw_complex *)&(*srcUpperLeft), 01296 (fftw_complex *)complexResultImg.begin(), 01297 FFTW_FORWARD, FFTW_ESTIMATE ); 01298 fftw_execute(forwardPlan); 01299 fftw_destroy_plan(forwardPlan); 01300 01301 // convolve in freq. domain (in complexResultImg) 01302 combineTwoImages(srcImageRange(complexResultImg), srcIter(filterUpperLeft, fa), 01303 destImage(complexResultImg), std::multiplies<FFTWComplex>()); 01304 01305 // FFT back into spatial domain (inplace in complexResultImg) 01306 fftw_plan backwardPlan= 01307 fftw_plan_dft_2d(h, w, (fftw_complex *)complexResultImg.begin(), 01308 (fftw_complex *)complexResultImg.begin(), 01309 FFTW_BACKWARD, FFTW_ESTIMATE); 01310 fftw_execute(backwardPlan); 01311 fftw_destroy_plan(backwardPlan); 01312 01313 typedef typename 01314 NumericTraits<typename DestAccessor::value_type>::isScalar 01315 isScalarResult; 01316 01317 // normalization (after FFTs), maybe stripping imaginary part 01318 applyFourierFilterImplNormalization(complexResultImg, destUpperLeft, da, 01319 isScalarResult()); 01320 } 01321 01322 template <class DestImageIterator, class DestAccessor> 01323 void applyFourierFilterImplNormalization(FFTWComplexImage const &srcImage, 01324 DestImageIterator destUpperLeft, 01325 DestAccessor da, 01326 VigraFalseType) 01327 { 01328 double normFactor= 1.0/(srcImage.width() * srcImage.height()); 01329 01330 for(int y=0; y<srcImage.height(); y++, destUpperLeft.y++) 01331 { 01332 DestImageIterator dIt= destUpperLeft; 01333 for(int x= 0; x< srcImage.width(); x++, dIt.x++) 01334 { 01335 da.setComponent(srcImage(x, y).re()*normFactor, dIt, 0); 01336 da.setComponent(srcImage(x, y).im()*normFactor, dIt, 1); 01337 } 01338 } 01339 } 01340 01341 inline 01342 void applyFourierFilterImplNormalization(FFTWComplexImage const & srcImage, 01343 FFTWComplexImage::traverser destUpperLeft, 01344 FFTWComplexImage::Accessor da, 01345 VigraFalseType) 01346 { 01347 transformImage(srcImageRange(srcImage), destIter(destUpperLeft, da), 01348 linearIntensityTransform<FFTWComplex>(1.0/(srcImage.width() * srcImage.height()))); 01349 } 01350 01351 template <class DestImageIterator, class DestAccessor> 01352 void applyFourierFilterImplNormalization(FFTWComplexImage const & srcImage, 01353 DestImageIterator destUpperLeft, 01354 DestAccessor da, 01355 VigraTrueType) 01356 { 01357 double normFactor= 1.0/(srcImage.width() * srcImage.height()); 01358 01359 for(int y=0; y<srcImage.height(); y++, destUpperLeft.y++) 01360 { 01361 DestImageIterator dIt= destUpperLeft; 01362 for(int x= 0; x< srcImage.width(); x++, dIt.x++) 01363 da.set(srcImage(x, y).re()*normFactor, dIt); 01364 } 01365 } 01366 01367 /**********************************************************/ 01368 /* */ 01369 /* applyFourierFilterFamily */ 01370 /* */ 01371 /**********************************************************/ 01372 01373 /** \brief Apply an array of filters (defined in the frequency domain) to an image. 01374 01375 This provides the same functionality as \ref applyFourierFilter(), 01376 but applying several filters at once allows to avoid 01377 repeated Fourier transforms of the source image. 01378 01379 Filters and result images must be stored in \ref vigra::ImageArray data 01380 structures. In contrast to \ref applyFourierFilter(), this function adjusts 01381 the size of the result images and the the length of the array. 01382 01383 <b> Declarations:</b> 01384 01385 pass arguments explicitly: 01386 \code 01387 namespace vigra { 01388 template <class SrcImageIterator, class SrcAccessor, class FilterType> 01389 void applyFourierFilterFamily(SrcImageIterator srcUpperLeft, 01390 SrcImageIterator srcLowerRight, SrcAccessor sa, 01391 const ImageArray<FilterType> &filters, 01392 ImageArray<FFTWComplexImage> &results) 01393 } 01394 \endcode 01395 01396 use argument objects in conjunction with \ref ArgumentObjectFactories: 01397 \code 01398 namespace vigra { 01399 template <class SrcImageIterator, class SrcAccessor, class FilterType> 01400 inline 01401 void applyFourierFilterFamily(triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src, 01402 const ImageArray<FilterType> &filters, 01403 ImageArray<FFTWComplexImage> &results) 01404 } 01405 \endcode 01406 01407 <b> Usage:</b> 01408 01409 <b>\#include</b> "<a href="fftw3_8hxx-source.html">vigra/fftw3.hxx</a>"<br> 01410 Namespace: vigra 01411 01412 \code 01413 // assuming the presence of a real-valued image named "image" to 01414 // be filtered in this example 01415 01416 vigra::ImageArray<vigra::FImage> filters(16, image.size()); 01417 01418 for (int i=0; i<filters.size(); i++) 01419 // create some meaningful filters here 01420 createMyFilterOfScale(i, destImage(filters[i])); 01421 01422 vigra::ImageArray<vigra::FFTWComplexImage> results(); 01423 01424 vigra::applyFourierFilterFamily(srcImageRange(image), filters, results); 01425 \endcode 01426 */ 01427 template <class SrcImageIterator, class SrcAccessor, 01428 class FilterType, class DestImage> 01429 inline 01430 void applyFourierFilterFamily(triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src, 01431 const ImageArray<FilterType> &filters, 01432 ImageArray<DestImage> &results) 01433 { 01434 applyFourierFilterFamily(src.first, src.second, src.third, 01435 filters, results); 01436 } 01437 01438 template <class SrcImageIterator, class SrcAccessor, 01439 class FilterType, class DestImage> 01440 void applyFourierFilterFamily(SrcImageIterator srcUpperLeft, 01441 SrcImageIterator srcLowerRight, SrcAccessor sa, 01442 const ImageArray<FilterType> &filters, 01443 ImageArray<DestImage> &results) 01444 { 01445 int w= srcLowerRight.x - srcUpperLeft.x; 01446 int h= srcLowerRight.y - srcUpperLeft.y; 01447 01448 FFTWComplexImage workImage(w, h); 01449 copyImage(srcIterRange(srcUpperLeft, srcLowerRight, sa), 01450 destImage(workImage, FFTWWriteRealAccessor())); 01451 01452 FFTWComplexImage const & cworkImage = workImage; 01453 applyFourierFilterFamilyImpl(cworkImage.upperLeft(), cworkImage.lowerRight(), cworkImage.accessor(), 01454 filters, results); 01455 } 01456 01457 template <class FilterType, class DestImage> 01458 inline 01459 void applyFourierFilterFamily( 01460 FFTWComplexImage::const_traverser srcUpperLeft, 01461 FFTWComplexImage::const_traverser srcLowerRight, 01462 FFTWComplexImage::ConstAccessor sa, 01463 const ImageArray<FilterType> &filters, 01464 ImageArray<DestImage> &results) 01465 { 01466 int w= srcLowerRight.x - srcUpperLeft.x; 01467 01468 // test for right memory layout (fftw expects a 2*width*height floats array) 01469 if (&(*(srcUpperLeft + Diff2D(w, 0))) == &(*(srcUpperLeft + Diff2D(0, 1)))) 01470 applyFourierFilterFamilyImpl(srcUpperLeft, srcLowerRight, sa, 01471 filters, results); 01472 else 01473 { 01474 int h = srcLowerRight.y - srcUpperLeft.y; 01475 FFTWComplexImage workImage(w, h); 01476 copyImage(srcIterRange(srcUpperLeft, srcLowerRight, sa), 01477 destImage(workImage)); 01478 01479 FFTWComplexImage const & cworkImage = workImage; 01480 applyFourierFilterFamilyImpl(cworkImage.upperLeft(), cworkImage.lowerRight(), cworkImage.accessor(), 01481 filters, results); 01482 } 01483 } 01484 01485 template <class FilterType, class DestImage> 01486 void applyFourierFilterFamilyImpl( 01487 FFTWComplexImage::const_traverser srcUpperLeft, 01488 FFTWComplexImage::const_traverser srcLowerRight, 01489 FFTWComplexImage::ConstAccessor sa, 01490 const ImageArray<FilterType> &filters, 01491 ImageArray<DestImage> &results) 01492 { 01493 // make sure the filter images have the right dimensions 01494 vigra_precondition((srcLowerRight - srcUpperLeft) == filters.imageSize(), 01495 "applyFourierFilterFamily called with src image size != filters.imageSize()!"); 01496 01497 // make sure the result image array has the right dimensions 01498 results.resize(filters.size()); 01499 results.resizeImages(filters.imageSize()); 01500 01501 // FFT from srcImage to freqImage 01502 int w= srcLowerRight.x - srcUpperLeft.x; 01503 int h= srcLowerRight.y - srcUpperLeft.y; 01504 01505 FFTWComplexImage freqImage(w, h); 01506 FFTWComplexImage result(w, h); 01507 01508 fftw_plan forwardPlan= 01509 fftw_plan_dft_2d(h, w, (fftw_complex *)&(*srcUpperLeft), 01510 (fftw_complex *)freqImage.begin(), 01511 FFTW_FORWARD, FFTW_ESTIMATE ); 01512 fftw_execute(forwardPlan); 01513 fftw_destroy_plan(forwardPlan); 01514 01515 fftw_plan backwardPlan= 01516 fftw_plan_dft_2d(h, w, (fftw_complex *)result.begin(), 01517 (fftw_complex *)result.begin(), 01518 FFTW_BACKWARD, FFTW_ESTIMATE ); 01519 typedef typename 01520 NumericTraits<typename DestImage::Accessor::value_type>::isScalar 01521 isScalarResult; 01522 01523 // convolve with filters in freq. domain 01524 for (unsigned int i= 0; i < filters.size(); i++) 01525 { 01526 combineTwoImages(srcImageRange(freqImage), srcImage(filters[i]), 01527 destImage(result), std::multiplies<FFTWComplex>()); 01528 01529 // FFT back into spatial domain (inplace in destImage) 01530 fftw_execute(backwardPlan); 01531 01532 // normalization (after FFTs), maybe stripping imaginary part 01533 applyFourierFilterImplNormalization(result, 01534 results[i].upperLeft(), results[i].accessor(), 01535 isScalarResult()); 01536 } 01537 fftw_destroy_plan(backwardPlan); 01538 } 01539 01540 /********************************************************/ 01541 /* */ 01542 /* fourierTransformReal */ 01543 /* */ 01544 /********************************************************/ 01545 01546 /** \brief Real Fourier transforms for even and odd boundary conditions 01547 (aka. cosine and sine transforms). 01548 01549 01550 If the image is real and has even symmetry, its Fourier transform 01551 is also real and has even symmetry. The Fourier transform of a real image with odd 01552 symmetry is imaginary and has odd symmetry. In either case, only about a quarter 01553 of the pixels need to be stored because the rest can be calculated from the symmetry 01554 properties. This is especially useful, if the original image is implicitly assumed 01555 to have reflective or anti-reflective boundary conditions. Then the "negative" 01556 pixel locations are defined as 01557 01558 \code 01559 even (reflective boundary conditions): f[-x] = f[x] (x = 1,...,N-1) 01560 odd (anti-reflective boundary conditions): f[-1] = 0 01561 f[-x] = -f[x-2] (x = 2,...,N-1) 01562 \endcode 01563 01564 end similar at the other boundary (see the FFTW documentation for details). 01565 This has the advantage that more efficient Fourier transforms that use only 01566 real numbers can be implemented. These are also known as cosine and sine transforms 01567 respectively. 01568 01569 If you use the odd transform it is important to note that in the Fourier domain, 01570 the DC component is always zero and is therefore dropped from the data structure. 01571 This means that index 0 in an odd symmetric Fourier domain image refers to 01572 the <i>first</i> harmonic. This is especially important if an image is first 01573 cosine transformed (even symmetry), then in the Fourier domain multiplied 01574 with an odd symmetric filter (e.g. a first derivative) and finally transformed 01575 back to the spatial domain with a sine transform (odd symmetric). For this to work 01576 properly the image must be shifted left or up by one pixel (depending on whether 01577 the x- or y-axis is odd symmetric) before the inverse transform can be applied. 01578 (see example below). 01579 01580 The real Fourier transform functions are named <tt>fourierTransformReal??</tt> 01581 where the questions marks stand for either <tt>E</tt> or <tt>O</tt> indicating 01582 whether the x- and y-axis is to be transformed using even or odd symmetry. 01583 The same functions can be used for both the forward and inverse transforms, 01584 only the normalization changes. For signal processing, the following 01585 normalization factors are most appropriate: 01586 01587 \code 01588 forward inverse 01589 ------------------------------------------------------------ 01590 X even, Y even 1.0 4.0 * (w-1) * (h-1) 01591 X even, Y odd -1.0 -4.0 * (w-1) * (h+1) 01592 X odd, Y even -1.0 -4.0 * (w+1) * (h-1) 01593 X odd, Y odd 1.0 4.0 * (w+1) * (h+1) 01594 \endcode 01595 01596 where <tt>w</tt> and <tt>h</tt> denote the image width and height. 01597 01598 <b> Declarations:</b> 01599 01600 pass arguments explicitly: 01601 \code 01602 namespace vigra { 01603 template <class SrcTraverser, class SrcAccessor, 01604 class DestTraverser, class DestAccessor> 01605 void 01606 fourierTransformRealEE(SrcTraverser sul, SrcTraverser slr, SrcAccessor src, 01607 DestTraverser dul, DestAccessor dest, fftw_real norm); 01608 01609 fourierTransformRealEO, fourierTransformRealOE, fourierTransformRealOO likewise 01610 } 01611 \endcode 01612 01613 01614 use argument objects in conjunction with \ref ArgumentObjectFactories: 01615 \code 01616 namespace vigra { 01617 template <class SrcTraverser, class SrcAccessor, 01618 class DestTraverser, class DestAccessor> 01619 void 01620 fourierTransformRealEE(triple<SrcTraverser, SrcTraverser, SrcAccessor> src, 01621 pair<DestTraverser, DestAccessor> dest, fftw_real norm); 01622 01623 fourierTransformRealEO, fourierTransformRealOE, fourierTransformRealOO likewise 01624 } 01625 \endcode 01626 01627 <b> Usage:</b> 01628 01629 <b>\#include</b> "<a href="fftw3_8hxx-source.html">vigra/fftw3.hxx</a>"<br> 01630 Namespace: vigra 01631 01632 \code 01633 vigra::FImage spatial(width,height), fourier(width,height); 01634 ... // fill image with data 01635 01636 // forward cosine transform == reflective boundary conditions 01637 fourierTransformRealEE(srcImageRange(spatial), destImage(fourier), (fftw_real)1.0); 01638 01639 // multiply with a first derivative of Gaussian in x-direction 01640 for(int y = 0; y < height; ++y) 01641 { 01642 for(int x = 1; x < width; ++x) 01643 { 01644 double dx = x * M_PI / (width - 1); 01645 double dy = y * M_PI / (height - 1); 01646 fourier(x-1, y) = fourier(x, y) * dx * exp(-(dx*dx + dy*dy) * scale*scale / 2.0); 01647 } 01648 fourier(width-1, y) = 0.0; 01649 } 01650 01651 // inverse transform -- odd symmetry in x-direction, even in y, 01652 // due to symmetry of the filter 01653 fourierTransformRealOE(srcImageRange(fourier), destImage(spatial), 01654 (fftw_real)-4.0 * (width+1) * (height-1)); 01655 \endcode 01656 */ 01657 template <class SrcTraverser, class SrcAccessor, 01658 class DestTraverser, class DestAccessor> 01659 inline void 01660 fourierTransformRealEE(triple<SrcTraverser, SrcTraverser, SrcAccessor> src, 01661 pair<DestTraverser, DestAccessor> dest, fftw_real norm) 01662 { 01663 fourierTransformRealEE(src.first, src.second, src.third, 01664 dest.first, dest.second, norm); 01665 } 01666 01667 template <class SrcTraverser, class SrcAccessor, 01668 class DestTraverser, class DestAccessor> 01669 inline void 01670 fourierTransformRealEE(SrcTraverser sul, SrcTraverser slr, SrcAccessor src, 01671 DestTraverser dul, DestAccessor dest, fftw_real norm) 01672 { 01673 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest, 01674 norm, FFTW_REDFT00, FFTW_REDFT00); 01675 } 01676 01677 template <class DestTraverser, class DestAccessor> 01678 inline void 01679 fourierTransformRealEE( 01680 FFTWRealImage::const_traverser sul, 01681 FFTWRealImage::const_traverser slr, 01682 FFTWRealImage::Accessor src, 01683 DestTraverser dul, DestAccessor dest, fftw_real norm) 01684 { 01685 int w = slr.x - sul.x; 01686 01687 // test for right memory layout (fftw expects a width*height fftw_real array) 01688 if (&(*(sul + Diff2D(w, 0))) == &(*(sul + Diff2D(0, 1)))) 01689 fourierTransformRealImpl(sul, slr, dul, dest, 01690 norm, FFTW_REDFT00, FFTW_REDFT00); 01691 else 01692 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest, 01693 norm, FFTW_REDFT00, FFTW_REDFT00); 01694 } 01695 01696 /********************************************************************/ 01697 01698 template <class SrcTraverser, class SrcAccessor, 01699 class DestTraverser, class DestAccessor> 01700 inline void 01701 fourierTransformRealOE(triple<SrcTraverser, SrcTraverser, SrcAccessor> src, 01702 pair<DestTraverser, DestAccessor> dest, fftw_real norm) 01703 { 01704 fourierTransformRealOE(src.first, src.second, src.third, 01705 dest.first, dest.second, norm); 01706 } 01707 01708 template <class SrcTraverser, class SrcAccessor, 01709 class DestTraverser, class DestAccessor> 01710 inline void 01711 fourierTransformRealOE(SrcTraverser sul, SrcTraverser slr, SrcAccessor src, 01712 DestTraverser dul, DestAccessor dest, fftw_real norm) 01713 { 01714 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest, 01715 norm, FFTW_RODFT00, FFTW_REDFT00); 01716 } 01717 01718 template <class DestTraverser, class DestAccessor> 01719 inline void 01720 fourierTransformRealOE( 01721 FFTWRealImage::const_traverser sul, 01722 FFTWRealImage::const_traverser slr, 01723 FFTWRealImage::Accessor src, 01724 DestTraverser dul, DestAccessor dest, fftw_real norm) 01725 { 01726 int w = slr.x - sul.x; 01727 01728 // test for right memory layout (fftw expects a width*height fftw_real array) 01729 if (&(*(sul + Diff2D(w, 0))) == &(*(sul + Diff2D(0, 1)))) 01730 fourierTransformRealImpl(sul, slr, dul, dest, 01731 norm, FFTW_RODFT00, FFTW_REDFT00); 01732 else 01733 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest, 01734 norm, FFTW_RODFT00, FFTW_REDFT00); 01735 } 01736 01737 /********************************************************************/ 01738 01739 template <class SrcTraverser, class SrcAccessor, 01740 class DestTraverser, class DestAccessor> 01741 inline void 01742 fourierTransformRealEO(triple<SrcTraverser, SrcTraverser, SrcAccessor> src, 01743 pair<DestTraverser, DestAccessor> dest, fftw_real norm) 01744 { 01745 fourierTransformRealEO(src.first, src.second, src.third, 01746 dest.first, dest.second, norm); 01747 } 01748 01749 template <class SrcTraverser, class SrcAccessor, 01750 class DestTraverser, class DestAccessor> 01751 inline void 01752 fourierTransformRealEO(SrcTraverser sul, SrcTraverser slr, SrcAccessor src, 01753 DestTraverser dul, DestAccessor dest, fftw_real norm) 01754 { 01755 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest, 01756 norm, FFTW_REDFT00, FFTW_RODFT00); 01757 } 01758 01759 template <class DestTraverser, class DestAccessor> 01760 inline void 01761 fourierTransformRealEO( 01762 FFTWRealImage::const_traverser sul, 01763 FFTWRealImage::const_traverser slr, 01764 FFTWRealImage::Accessor src, 01765 DestTraverser dul, DestAccessor dest, fftw_real norm) 01766 { 01767 int w = slr.x - sul.x; 01768 01769 // test for right memory layout (fftw expects a width*height fftw_real array) 01770 if (&(*(sul + Diff2D(w, 0))) == &(*(sul + Diff2D(0, 1)))) 01771 fourierTransformRealImpl(sul, slr, dul, dest, 01772 norm, FFTW_REDFT00, FFTW_RODFT00); 01773 else 01774 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest, 01775 norm, FFTW_REDFT00, FFTW_RODFT00); 01776 } 01777 01778 /********************************************************************/ 01779 01780 template <class SrcTraverser, class SrcAccessor, 01781 class DestTraverser, class DestAccessor> 01782 inline void 01783 fourierTransformRealOO(triple<SrcTraverser, SrcTraverser, SrcAccessor> src, 01784 pair<DestTraverser, DestAccessor> dest, fftw_real norm) 01785 { 01786 fourierTransformRealOO(src.first, src.second, src.third, 01787 dest.first, dest.second, norm); 01788 } 01789 01790 template <class SrcTraverser, class SrcAccessor, 01791 class DestTraverser, class DestAccessor> 01792 inline void 01793 fourierTransformRealOO(SrcTraverser sul, SrcTraverser slr, SrcAccessor src, 01794 DestTraverser dul, DestAccessor dest, fftw_real norm) 01795 { 01796 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest, 01797 norm, FFTW_RODFT00, FFTW_RODFT00); 01798 } 01799 01800 template <class DestTraverser, class DestAccessor> 01801 inline void 01802 fourierTransformRealOO( 01803 FFTWRealImage::const_traverser sul, 01804 FFTWRealImage::const_traverser slr, 01805 FFTWRealImage::Accessor src, 01806 DestTraverser dul, DestAccessor dest, fftw_real norm) 01807 { 01808 int w = slr.x - sul.x; 01809 01810 // test for right memory layout (fftw expects a width*height fftw_real array) 01811 if (&(*(sul + Diff2D(w, 0))) == &(*(sul + Diff2D(0, 1)))) 01812 fourierTransformRealImpl(sul, slr, dul, dest, 01813 norm, FFTW_RODFT00, FFTW_RODFT00); 01814 else 01815 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest, 01816 norm, FFTW_RODFT00, FFTW_RODFT00); 01817 } 01818 01819 /*******************************************************************/ 01820 01821 template <class SrcTraverser, class SrcAccessor, 01822 class DestTraverser, class DestAccessor> 01823 void 01824 fourierTransformRealWorkImageImpl(SrcTraverser sul, SrcTraverser slr, SrcAccessor src, 01825 DestTraverser dul, DestAccessor dest, 01826 fftw_real norm, fftw_r2r_kind kindx, fftw_r2r_kind kindy) 01827 { 01828 FFTWRealImage workImage(slr - sul); 01829 copyImage(srcIterRange(sul, slr, src), destImage(workImage)); 01830 FFTWRealImage const & cworkImage = workImage; 01831 fourierTransformRealImpl(cworkImage.upperLeft(), cworkImage.lowerRight(), 01832 dul, dest, norm, kindx, kindy); 01833 } 01834 01835 01836 template <class DestTraverser, class DestAccessor> 01837 void 01838 fourierTransformRealImpl( 01839 FFTWRealImage::const_traverser sul, 01840 FFTWRealImage::const_traverser slr, 01841 DestTraverser dul, DestAccessor dest, 01842 fftw_real norm, fftw_r2r_kind kindx, fftw_r2r_kind kindy) 01843 { 01844 int w = slr.x - sul.x; 01845 int h = slr.y - sul.y; 01846 BasicImage<fftw_real> res(w, h); 01847 01848 fftw_plan plan = fftw_plan_r2r_2d(h, w, 01849 (fftw_real *)&(*sul), (fftw_real *)res.begin(), 01850 kindy, kindx, FFTW_ESTIMATE); 01851 fftw_execute(plan); 01852 fftw_destroy_plan(plan); 01853 01854 if(norm != 1.0) 01855 transformImage(srcImageRange(res), destIter(dul, dest), 01856 std::bind1st(std::multiplies<fftw_real>(), 1.0 / norm)); 01857 else 01858 copyImage(srcImageRange(res), destIter(dul, dest)); 01859 } 01860 01861 01862 //@} 01863 01864 } // namespace vigra 01865 01866 #endif // VIGRA_FFTW3_HXX
© Ullrich Köthe (koethe@informatik.uni-hamburg.de) |
html generated using doxygen and Python
|