36 #ifndef VIGRA_FFTW3_HXX 37 #define VIGRA_FFTW3_HXX 42 #include "stdimage.hxx" 43 #include "copyimage.hxx" 44 #include "transformimage.hxx" 45 #include "combineimages.hxx" 46 #include "numerictraits.hxx" 47 #include "imagecontainer.hxx" 52 typedef double fftw_real;
58 struct FFTWReal<fftw_complex>
64 struct FFTWReal<fftwf_complex>
70 struct FFTWReal<fftwl_complex>
72 typedef long double type;
76 struct FFTWReal2Complex;
79 struct FFTWReal2Complex<double>
81 typedef fftw_complex type;
82 typedef fftw_plan plan_type;
86 struct FFTWReal2Complex<float>
88 typedef fftwf_complex type;
89 typedef fftwf_plan plan_type;
93 struct FFTWReal2Complex<long double>
95 typedef fftwl_complex type;
96 typedef fftwl_plan plan_type;
130 template <
class Real =
double>
169 FFTWComplex(value_type
const & re = 0.0, value_type
const & im = 0.0)
179 data_[0] = o.data_[0];
180 data_[1] = o.data_[1];
188 data_[0] = (Real)o.real();
189 data_[1] = (Real)o.imag();
196 data_[0] = (Real)o[0];
197 data_[1] = (Real)o[1];
204 data_[0] = (Real)o[0];
205 data_[1] = (Real)o[1];
212 data_[0] = (Real)o[0];
213 data_[1] = (Real)o[1];
221 data_[0] = (Real)o.real();
222 data_[1] = (Real)o.imag();
230 data_[0] = (Real)o[0];
231 data_[1] = (Real)o[1];
238 data_[0] = o.data_[0];
239 data_[1] = o.data_[1];
248 data_[0] = (Real)o.real();
249 data_[1] = (Real)o.imag();
257 data_[0] = (Real)o[0];
258 data_[1] = (Real)o[1];
266 data_[0] = (Real)o[0];
267 data_[1] = (Real)o[1];
275 data_[0] = (Real)o[0];
276 data_[1] = (Real)o[1];
312 data_[0] = (Real)o[0];
313 data_[1] = (Real)o[1];
322 data_[0] = (Real)o.real();
323 data_[1] = (Real)o.imag();
330 const_reference re()
const 336 const_reference
real()
const 342 const_reference im()
const 348 const_reference
imag()
const 359 {
return data_[0]*data_[0]+data_[1]*data_[1]; }
390 {
return data_ + 2; }
392 const_iterator begin()
const 395 const_iterator end()
const 396 {
return data_ + 2; }
508 struct NumericTraits<fftw_complex>
510 typedef fftw_complex Type;
511 typedef fftw_complex Promote;
512 typedef fftw_complex RealPromote;
513 typedef fftw_complex ComplexPromote;
514 typedef fftw_real ValueType;
516 typedef VigraFalseType isIntegral;
517 typedef VigraFalseType isScalar;
518 typedef NumericTraits<fftw_real>::isSigned isSigned;
519 typedef VigraFalseType isOrdered;
520 typedef VigraTrueType isComplex;
526 static const Promote & toPromote(
const Type & v) {
return v; }
527 static const RealPromote & toRealPromote(
const Type & v) {
return v; }
528 static const Type & fromPromote(
const Promote & v) {
return v; }
529 static const Type & fromRealPromote(
const RealPromote & v) {
return v; }
533 struct NumericTraits<FFTWComplex<Real> >
539 typedef typename Type::value_type ValueType;
541 typedef VigraFalseType isIntegral;
542 typedef VigraFalseType isScalar;
543 typedef typename NumericTraits<ValueType>::isSigned isSigned;
544 typedef VigraFalseType isOrdered;
545 typedef VigraTrueType isComplex;
547 static Type zero() {
return Type(0.0, 0.0); }
548 static Type one() {
return Type(1.0, 0.0); }
549 static Type nonZero() {
return one(); }
551 static const Promote & toPromote(
const Type & v) {
return v; }
552 static const RealPromote & toRealPromote(
const Type & v) {
return v; }
553 static const Type & fromPromote(
const Promote & v) {
return v; }
554 static const Type & fromRealPromote(
const RealPromote & v) {
return v; }
558 struct NormTraits<fftw_complex>
560 typedef fftw_complex Type;
566 struct NormTraits<FFTWComplex<Real> >
570 typedef typename Type::NormType
NormType;
574 struct PromoteTraits<fftw_complex, fftw_complex>
576 typedef fftw_complex Promote;
580 struct PromoteTraits<fftw_complex, double>
582 typedef fftw_complex Promote;
586 struct PromoteTraits<double, fftw_complex>
588 typedef fftw_complex Promote;
591 template <
class Real>
592 struct PromoteTraits<FFTWComplex<Real>, FFTWComplex<Real> >
597 template <
class Real>
598 struct PromoteTraits<FFTWComplex<Real>, double>
603 template <
class Real>
604 struct PromoteTraits<double, FFTWComplex<Real> >
610 struct CanSkipInitialization<std::complex<T> >
612 typedef typename CanSkipInitialization<T>::type type;
613 static const bool value = type::asBool;
617 struct CanSkipInitialization<FFTWComplex<Real> >
619 typedef typename CanSkipInitialization<Real>::type type;
620 static const bool value = type::asBool;
623 namespace multi_math {
626 struct MultiMathOperand;
628 template <
class Real>
629 struct MultiMathOperand<FFTWComplex<Real> >
631 typedef MultiMathOperand<FFTWComplex<Real> > AllowOverload;
634 static const int ndim = 0;
640 template <
class SHAPE>
641 bool checkShape(SHAPE
const &)
const 646 template <
class SHAPE>
652 void inc(
unsigned int )
const 655 void reset(
unsigned int )
const 672 typedef std::size_t size_type;
673 typedef std::ptrdiff_t difference_type;
675 typedef const Ty *const_pointer;
676 typedef Ty& reference;
677 typedef const Ty& const_reference;
678 typedef Ty value_type;
680 pointer address(reference val)
const 683 const_pointer address(const_reference val)
const 686 template<
class Other>
689 typedef FFTWAllocator<Other> other;
692 FFTWAllocator()
throw()
695 template<
class Other>
696 FFTWAllocator(
const FFTWAllocator<Other>& right)
throw()
699 template<
class Other>
700 FFTWAllocator& operator=(
const FFTWAllocator<Other>& right)
705 pointer allocate(size_type count,
void * = 0)
707 return (pointer)fftw_malloc(count *
sizeof(Ty));
710 void deallocate(pointer ptr, size_type count)
715 void construct(pointer ptr,
const Ty& val)
721 void destroy(pointer ptr)
726 size_type max_size()
const throw()
728 return NumericTraits<std::ptrdiff_t>::max() /
sizeof(Ty);
741 typedef size_t size_type;
742 typedef std::ptrdiff_t difference_type;
743 typedef value_type *pointer;
744 typedef const value_type *const_pointer;
745 typedef value_type& reference;
746 typedef const value_type& const_reference;
748 pointer address(reference val)
const 751 const_pointer address(const_reference val)
const 754 template<
class Other>
757 typedef allocator<Other> other;
763 template<
class Other>
764 allocator(
const allocator<Other>& right)
throw()
767 template<
class Other>
768 allocator& operator=(
const allocator<Other>& right)
773 pointer allocate(size_type count,
void * = 0)
775 return (pointer)fftw_malloc(count *
sizeof(value_type));
778 void deallocate(pointer ptr, size_type )
783 void construct(pointer ptr,
const value_type& val)
785 new(ptr) value_type(val);
789 void destroy(pointer ptr)
794 size_type max_size()
const throw()
796 return vigra::NumericTraits<std::ptrdiff_t>::max() /
sizeof(value_type);
826 return a.re() == b.re() && a.im() == b.im();
831 return a.re() == b && a.im() == 0.0;
836 return a == b.re() && 0.0 == b.im();
842 return a.re() != b.re() || a.im() != b.im();
848 return a.re() != b || a.im() != 0.0;
854 return a != b.re() || 0.0 != b.im();
877 a.im() = a.re()*b.im()+a.im()*b.re();
887 a.im() = (b.re()*a.im()-a.re()*b.im())/sm;
1049 #define VIGRA_DEFINE_FFTW_COMPLEX_FUNCTION(fct) \ 1050 template <class R> \ 1051 inline FFTWComplex<R> fct(const FFTWComplex<R> &a) \ 1053 return std::fct(reinterpret_cast<std::complex<R> const &>(a)); \ 1056 VIGRA_DEFINE_FFTW_COMPLEX_FUNCTION(
cos)
1057 VIGRA_DEFINE_FFTW_COMPLEX_FUNCTION(cosh)
1058 VIGRA_DEFINE_FFTW_COMPLEX_FUNCTION(
exp)
1059 VIGRA_DEFINE_FFTW_COMPLEX_FUNCTION(
log)
1060 VIGRA_DEFINE_FFTW_COMPLEX_FUNCTION(
log10)
1061 VIGRA_DEFINE_FFTW_COMPLEX_FUNCTION(
sin)
1062 VIGRA_DEFINE_FFTW_COMPLEX_FUNCTION(sinh)
1063 VIGRA_DEFINE_FFTW_COMPLEX_FUNCTION(
sqrt)
1064 VIGRA_DEFINE_FFTW_COMPLEX_FUNCTION(
tan)
1065 VIGRA_DEFINE_FFTW_COMPLEX_FUNCTION(tanh)
1067 #undef VIGRA_DEFINE_FFTW_COMPLEX_FUNCTION 1072 return std::pow(
reinterpret_cast<std::complex<R>
const &
>(a), e);
1078 return std::pow(
reinterpret_cast<std::complex<R>
const &
>(a), e);
1084 return std::pow(
reinterpret_cast<std::complex<R>
const &
>(a),
1085 reinterpret_cast<std::complex<R>
const &
>(e));
1091 return std::pow(a,
reinterpret_cast<std::complex<R>
const &
>(e));
1100 template <
class Real>
1101 ostream & operator<<(ostream & s, vigra::FFTWComplex<Real>
const & v)
1103 s << std::complex<Real>(v.re(), v.im());
1146 typedef Iterator iterator;
1149 typedef typename iterator::iterator_category iterator_category;
1150 typedef typename iterator::value_type value_type;
1151 typedef typename iterator::reference reference;
1152 typedef typename iterator::index_reference index_reference;
1153 typedef typename iterator::pointer pointer;
1154 typedef typename iterator::difference_type difference_type;
1155 typedef typename iterator::row_iterator row_iterator;
1156 typedef typename iterator::column_iterator column_iterator;
1159 typedef VigraTrueType hasConstantStrides;
1164 ConstBasicImageIterator<FFTWComplex<R>, FFTWComplex<R> **> >
1168 typedef Iterator iterator;
1171 typedef typename iterator::iterator_category iterator_category;
1172 typedef typename iterator::value_type value_type;
1173 typedef typename iterator::reference reference;
1174 typedef typename iterator::index_reference index_reference;
1175 typedef typename iterator::pointer pointer;
1176 typedef typename iterator::difference_type difference_type;
1177 typedef typename iterator::row_iterator row_iterator;
1178 typedef typename iterator::column_iterator column_iterator;
1181 typedef VigraTrueType hasConstantStrides;
1216 template <
class Real =
double>
1225 template <
class ITERATOR>
1231 template <
class ITERATOR,
class DIFFERENCE>
1237 template <
class ITERATOR>
1238 void set(value_type
const & v, ITERATOR
const & i)
const {
1243 template <
class ITERATOR,
class DIFFERENCE>
1244 void set(value_type
const & v, ITERATOR
const & i, DIFFERENCE d)
const {
1249 template <
class R,
class ITERATOR>
1255 template <
class R,
class ITERATOR,
class DIFFERENCE>
1267 template <
class Real =
double>
1275 template <
class ITERATOR>
1281 template <
class ITERATOR,
class DIFFERENCE>
1287 template <
class ITERATOR>
1288 void set(value_type
const & v, ITERATOR
const & i)
const {
1293 template <
class ITERATOR,
class DIFFERENCE>
1294 void set(value_type
const & v, ITERATOR
const & i, DIFFERENCE d)
const {
1299 template <
class R,
class ITERATOR>
1305 template <
class R,
class ITERATOR,
class DIFFERENCE>
1318 template <
class Real =
double>
1329 template <
class ITERATOR>
1330 void set(value_type
const & v, ITERATOR
const & i)
const {
1338 template <
class ITERATOR,
class DIFFERENCE>
1339 void set(value_type
const & v, ITERATOR
const & i, DIFFERENCE d)
const {
1350 template <
class Real =
double>
1358 template <
class ITERATOR>
1360 return (*i).squaredMagnitude();
1364 template <
class ITERATOR,
class DIFFERENCE>
1366 return (i[d]).squaredMagnitude();
1376 template <
class Real =
double>
1384 template <
class ITERATOR>
1386 return (*i).magnitude();
1390 template <
class ITERATOR,
class DIFFERENCE>
1392 return (i[d]).magnitude();
1402 template <
class Real =
double>
1410 template <
class ITERATOR>
1412 return std::log((*i).magnitude() + 1);
1416 template <
class ITERATOR,
class DIFFERENCE>
1418 return std::log((i[d]).magnitude() + 1);
1428 template <
class Real =
double>
1436 template <
class ITERATOR>
1438 return (*i).phase();
1442 template <
class ITERATOR,
class DIFFERENCE>
1444 return (i[d]).phase();
1620 template <
class SrcImageIterator,
class SrcAccessor,
1621 class DestImageIterator,
class DestAccessor>
1623 SrcImageIterator src_lowerright, SrcAccessor sa,
1624 DestImageIterator dest_upperleft, DestAccessor da)
1626 int w = int(src_lowerright.x - src_upperleft.x);
1627 int h = int(src_lowerright.y - src_upperleft.y);
1635 src_upperleft +
Diff2D(w2, h2), sa),
1636 destIter (dest_upperleft +
Diff2D(w1, h1), da));
1640 src_lowerright, sa),
1641 destIter (dest_upperleft, da));
1645 src_upperleft +
Diff2D(w, h2), sa),
1646 destIter (dest_upperleft +
Diff2D(0, h1), da));
1650 src_upperleft +
Diff2D(w2, h), sa),
1651 destIter (dest_upperleft +
Diff2D(w1, 0), da));
1654 template <
class SrcImageIterator,
class SrcAccessor,
1655 class DestImageIterator,
class DestAccessor>
1657 triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src,
1658 pair<DestImageIterator, DestAccessor> dest)
1661 dest.first, dest.second);
1711 template <
class SrcImageIterator,
class SrcAccessor,
1712 class DestImageIterator,
class DestAccessor>
1714 SrcImageIterator src_lowerright, SrcAccessor sa,
1715 DestImageIterator dest_upperleft, DestAccessor da)
1717 int w = int(src_lowerright.x - src_upperleft.x);
1718 int h = int(src_lowerright.y - src_upperleft.y);
1726 src_upperleft +
Diff2D(w2, h2), sa),
1727 destIter (dest_upperleft +
Diff2D(w1, h1), da));
1731 src_lowerright, sa),
1732 destIter (dest_upperleft, da));
1736 src_upperleft +
Diff2D(w, h2), sa),
1737 destIter (dest_upperleft +
Diff2D(0, h1), da));
1741 src_upperleft +
Diff2D(w2, h), sa),
1742 destIter (dest_upperleft +
Diff2D(w1, 0), da));
1745 template <
class SrcImageIterator,
class SrcAccessor,
1746 class DestImageIterator,
class DestAccessor>
1748 triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src,
1749 pair<DestImageIterator, DestAccessor> dest)
1752 dest.first, dest.second);
1763 int w = int(slr.x - sul.x);
1764 int h = int(slr.y - sul.y);
1766 FFTWComplexImage sworkImage, dworkImage;
1768 fftw_complex * srcPtr = (fftw_complex *)(&*sul);
1769 fftw_complex * destPtr = (fftw_complex *)(&*dul);
1772 if (h > 1 && &(*(sul +
Diff2D(w, 0))) != &(*(sul +
Diff2D(0, 1))))
1775 copyImage(srcIterRange(sul, slr, src), destImage(sworkImage));
1776 srcPtr = (fftw_complex *)(&(*sworkImage.
upperLeft()));
1778 if (h > 1 && &(*(dul +
Diff2D(w, 0))) != &(*(dul +
Diff2D(0, 1))))
1781 destPtr = (fftw_complex *)(&(*dworkImage.
upperLeft()));
1784 fftw_plan plan = fftw_plan_dft_2d(h, w, srcPtr, destPtr, sign, FFTW_ESTIMATE );
1786 fftw_destroy_plan(plan);
1788 if (h > 1 && &(*(dul +
Diff2D(w, 0))) != &(*(dul +
Diff2D(0, 1))))
1790 copyImage(srcImageRange(dworkImage), destIter(dul, dest));
1949 detail::fourierTransformImpl(sul, slr, src, dul, dest, FFTW_FORWARD);
1952 template <
class SrcImageIterator,
class SrcAccessor>
1954 SrcImageIterator srcLowerRight, SrcAccessor sa,
1958 int w= srcLowerRight.x - srcUpperLeft.x;
1959 int h= srcLowerRight.y - srcUpperLeft.y;
1961 FFTWComplexImage workImage(w, h);
1962 copyImage(srcIterRange(srcUpperLeft, srcLowerRight, sa),
1966 FFTWComplexImage
const & cworkImage = workImage;
1971 template <
class SrcImageIterator,
class SrcAccessor>
1973 void fourierTransform(triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src,
1974 pair<FFTWComplexImage::traverser, FFTWComplexImage::Accessor> dest)
1976 fourierTransform(src.first, src.second, src.third, dest.first, dest.second);
1990 detail::fourierTransformImpl(sul, slr, src, dul, dest, FFTW_BACKWARD);
1993 template <
class DestImageIterator,
class DestAccessor>
1996 DestImageIterator dul, DestAccessor dest)
1998 int w = slr.x - sul.x;
1999 int h = slr.y - sul.y;
2001 FFTWComplexImage workImage(w, h);
2003 copyImage(srcImageRange(workImage), destIter(dul, dest));
2007 template <
class DestImageIterator,
class DestAccessor>
2011 pair<DestImageIterator, DestAccessor> dest)
2114 template <
class SrcImageIterator,
class SrcAccessor,
2115 class FilterImageIterator,
class FilterAccessor,
2116 class DestImageIterator,
class DestAccessor>
2118 SrcImageIterator srcLowerRight, SrcAccessor sa,
2119 FilterImageIterator filterUpperLeft, FilterAccessor fa,
2120 DestImageIterator destUpperLeft, DestAccessor da)
2123 int w = int(srcLowerRight.x - srcUpperLeft.x);
2124 int h = int(srcLowerRight.y - srcUpperLeft.y);
2126 FFTWComplexImage workImage(w, h);
2127 copyImage(srcIterRange(srcUpperLeft, srcLowerRight, sa),
2131 FFTWComplexImage
const & cworkImage = workImage;
2133 filterUpperLeft, fa,
2137 template <
class FilterImageIterator,
class FilterAccessor,
2138 class DestImageIterator,
class DestAccessor>
2144 FilterImageIterator filterUpperLeft, FilterAccessor fa,
2145 DestImageIterator destUpperLeft, DestAccessor da)
2147 int w = srcLowerRight.x - srcUpperLeft.x;
2148 int h = srcLowerRight.y - srcUpperLeft.y;
2151 if (&(*(srcUpperLeft +
Diff2D(w, 0))) == &(*(srcUpperLeft +
Diff2D(0, 1))))
2152 applyFourierFilterImpl(srcUpperLeft, srcLowerRight, sa,
2153 filterUpperLeft, fa,
2157 FFTWComplexImage workImage(w, h);
2158 copyImage(srcIterRange(srcUpperLeft, srcLowerRight, sa),
2159 destImage(workImage));
2161 FFTWComplexImage
const & cworkImage = workImage;
2163 filterUpperLeft, fa,
2168 template <
class SrcImageIterator,
class SrcAccessor,
2169 class FilterImageIterator,
class FilterAccessor,
2170 class DestImageIterator,
class DestAccessor>
2173 pair<FilterImageIterator, FilterAccessor> filter,
2174 pair<DestImageIterator, DestAccessor> dest)
2177 filter.first, filter.second,
2178 dest.first, dest.second);
2181 template <
class FilterImageIterator,
class FilterAccessor,
2182 class DestImageIterator,
class DestAccessor>
2183 void applyFourierFilterImpl(
2187 FilterImageIterator filterUpperLeft, FilterAccessor fa,
2188 DestImageIterator destUpperLeft, DestAccessor da)
2190 int w = int(srcLowerRight.x - srcUpperLeft.x);
2191 int h = int(srcLowerRight.y - srcUpperLeft.y);
2193 FFTWComplexImage complexResultImg(srcLowerRight - srcUpperLeft);
2196 fftw_plan forwardPlan=
2197 fftw_plan_dft_2d(h, w, (fftw_complex *)&(*srcUpperLeft),
2198 (fftw_complex *)complexResultImg.
begin(),
2199 FFTW_FORWARD, FFTW_ESTIMATE );
2200 fftw_execute(forwardPlan);
2201 fftw_destroy_plan(forwardPlan);
2204 combineTwoImages(srcImageRange(complexResultImg), srcIter(filterUpperLeft, fa),
2205 destImage(complexResultImg), std::multiplies<
FFTWComplex<> >());
2208 fftw_plan backwardPlan=
2209 fftw_plan_dft_2d(h, w, (fftw_complex *)complexResultImg.
begin(),
2210 (fftw_complex *)complexResultImg.
begin(),
2211 FFTW_BACKWARD, FFTW_ESTIMATE);
2212 fftw_execute(backwardPlan);
2213 fftw_destroy_plan(backwardPlan);
2216 NumericTraits<typename DestAccessor::value_type>::isScalar
2220 applyFourierFilterImplNormalization(complexResultImg, destUpperLeft, da,
2224 template <
class DestImageIterator,
class DestAccessor>
2225 void applyFourierFilterImplNormalization(FFTWComplexImage
const &srcImage,
2226 DestImageIterator destUpperLeft,
2230 double normFactor= 1.0/(srcImage.
width() * srcImage.
height());
2232 for(
int y=0; y<srcImage.
height(); y++, destUpperLeft.y++)
2234 DestImageIterator dIt= destUpperLeft;
2235 for(
int x= 0; x< srcImage.
width(); x++, dIt.x++)
2237 da.setComponent(srcImage(x, y).re()*normFactor, dIt, 0);
2238 da.setComponent(srcImage(x, y).im()*normFactor, dIt, 1);
2244 void applyFourierFilterImplNormalization(FFTWComplexImage
const & srcImage,
2249 transformImage(srcImageRange(srcImage), destIter(destUpperLeft, da),
2253 template <
class DestImageIterator,
class DestAccessor>
2254 void applyFourierFilterImplNormalization(FFTWComplexImage
const & srcImage,
2255 DestImageIterator destUpperLeft,
2259 double normFactor= 1.0/(srcImage.
width() * srcImage.
height());
2261 for(
int y=0; y<srcImage.
height(); y++, destUpperLeft.y++)
2263 DestImageIterator dIt= destUpperLeft;
2264 for(
int x= 0; x< srcImage.
width(); x++, dIt.x++)
2265 da.set(srcImage(x, y).re()*normFactor, dIt);
2341 template <
class SrcImageIterator,
class SrcAccessor,
2342 class FilterType,
class DestImage>
2352 template <
class SrcImageIterator,
class SrcAccessor,
2353 class FilterType,
class DestImage>
2355 SrcImageIterator srcLowerRight, SrcAccessor sa,
2359 int w = int(srcLowerRight.x - srcUpperLeft.x);
2360 int h = int(srcLowerRight.y - srcUpperLeft.y);
2362 FFTWComplexImage workImage(w, h);
2363 copyImage(srcIterRange(srcUpperLeft, srcLowerRight, sa),
2366 FFTWComplexImage
const & cworkImage = workImage;
2371 template <
class FilterType,
class DestImage>
2380 int w= srcLowerRight.x - srcUpperLeft.x;
2383 if (&(*(srcUpperLeft +
Diff2D(w, 0))) == &(*(srcUpperLeft +
Diff2D(0, 1))))
2384 applyFourierFilterFamilyImpl(srcUpperLeft, srcLowerRight, sa,
2388 int h = srcLowerRight.y - srcUpperLeft.y;
2389 FFTWComplexImage workImage(w, h);
2390 copyImage(srcIterRange(srcUpperLeft, srcLowerRight, sa),
2391 destImage(workImage));
2393 FFTWComplexImage
const & cworkImage = workImage;
2399 template <
class FilterType,
class DestImage>
2400 void applyFourierFilterFamilyImpl(
2411 vigra_precondition((srcLowerRight - srcUpperLeft) == filters.
imageSize(),
2412 "applyFourierFilterFamily called with src image size != filters.imageSize()!");
2419 int w = int(srcLowerRight.x - srcUpperLeft.x);
2420 int h = int(srcLowerRight.y - srcUpperLeft.y);
2422 FFTWComplexImage freqImage(w, h);
2423 FFTWComplexImage result(w, h);
2425 fftw_plan forwardPlan=
2426 fftw_plan_dft_2d(h, w, (fftw_complex *)&(*srcUpperLeft),
2427 (fftw_complex *)freqImage.
begin(),
2428 FFTW_FORWARD, FFTW_ESTIMATE );
2429 fftw_execute(forwardPlan);
2430 fftw_destroy_plan(forwardPlan);
2432 fftw_plan backwardPlan=
2433 fftw_plan_dft_2d(h, w, (fftw_complex *)result.
begin(),
2434 (fftw_complex *)result.
begin(),
2435 FFTW_BACKWARD, FFTW_ESTIMATE );
2437 NumericTraits<typename DestImage::Accessor::value_type>::isScalar
2441 for (
unsigned int i= 0; i < filters.
size(); i++)
2447 fftw_execute(backwardPlan);
2450 applyFourierFilterImplNormalization(result,
2451 results[i].upperLeft(), results[i].accessor(),
2454 fftw_destroy_plan(backwardPlan);
2589 template <
class SrcTraverser,
class SrcAccessor,
2590 class DestTraverser,
class DestAccessor>
2592 fourierTransformRealEE(triple<SrcTraverser, SrcTraverser, SrcAccessor> src,
2593 pair<DestTraverser, DestAccessor> dest, fftw_real
norm)
2595 fourierTransformRealEE(src.first, src.second, src.third,
2596 dest.first, dest.second, norm);
2599 template <
class SrcTraverser,
class SrcAccessor,
2600 class DestTraverser,
class DestAccessor>
2602 fourierTransformRealEE(SrcTraverser sul, SrcTraverser slr, SrcAccessor src,
2603 DestTraverser dul, DestAccessor dest, fftw_real norm)
2605 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest,
2606 norm, FFTW_REDFT00, FFTW_REDFT00);
2609 template <
class DestTraverser,
class DestAccessor>
2611 fourierTransformRealEE(
2615 DestTraverser dul, DestAccessor dest, fftw_real norm)
2617 int w = slr.x - sul.x;
2620 if (&(*(sul +
Diff2D(w, 0))) == &(*(sul +
Diff2D(0, 1))))
2621 fourierTransformRealImpl(sul, slr, dul, dest,
2622 norm, FFTW_REDFT00, FFTW_REDFT00);
2624 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest,
2625 norm, FFTW_REDFT00, FFTW_REDFT00);
2630 template <
class SrcTraverser,
class SrcAccessor,
2631 class DestTraverser,
class DestAccessor>
2633 fourierTransformRealOE(triple<SrcTraverser, SrcTraverser, SrcAccessor> src,
2634 pair<DestTraverser, DestAccessor> dest, fftw_real norm)
2636 fourierTransformRealOE(src.first, src.second, src.third,
2637 dest.first, dest.second, norm);
2640 template <
class SrcTraverser,
class SrcAccessor,
2641 class DestTraverser,
class DestAccessor>
2643 fourierTransformRealOE(SrcTraverser sul, SrcTraverser slr, SrcAccessor src,
2644 DestTraverser dul, DestAccessor dest, fftw_real norm)
2646 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest,
2647 norm, FFTW_RODFT00, FFTW_REDFT00);
2650 template <
class DestTraverser,
class DestAccessor>
2652 fourierTransformRealOE(
2656 DestTraverser dul, DestAccessor dest, fftw_real norm)
2658 int w = slr.x - sul.x;
2661 if (&(*(sul +
Diff2D(w, 0))) == &(*(sul +
Diff2D(0, 1))))
2662 fourierTransformRealImpl(sul, slr, dul, dest,
2663 norm, FFTW_RODFT00, FFTW_REDFT00);
2665 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest,
2666 norm, FFTW_RODFT00, FFTW_REDFT00);
2671 template <
class SrcTraverser,
class SrcAccessor,
2672 class DestTraverser,
class DestAccessor>
2674 fourierTransformRealEO(triple<SrcTraverser, SrcTraverser, SrcAccessor> src,
2675 pair<DestTraverser, DestAccessor> dest, fftw_real norm)
2677 fourierTransformRealEO(src.first, src.second, src.third,
2678 dest.first, dest.second, norm);
2681 template <
class SrcTraverser,
class SrcAccessor,
2682 class DestTraverser,
class DestAccessor>
2684 fourierTransformRealEO(SrcTraverser sul, SrcTraverser slr, SrcAccessor src,
2685 DestTraverser dul, DestAccessor dest, fftw_real norm)
2687 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest,
2688 norm, FFTW_REDFT00, FFTW_RODFT00);
2691 template <
class DestTraverser,
class DestAccessor>
2693 fourierTransformRealEO(
2697 DestTraverser dul, DestAccessor dest, fftw_real norm)
2699 int w = slr.x - sul.x;
2702 if (&(*(sul +
Diff2D(w, 0))) == &(*(sul +
Diff2D(0, 1))))
2703 fourierTransformRealImpl(sul, slr, dul, dest,
2704 norm, FFTW_REDFT00, FFTW_RODFT00);
2706 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest,
2707 norm, FFTW_REDFT00, FFTW_RODFT00);
2712 template <
class SrcTraverser,
class SrcAccessor,
2713 class DestTraverser,
class DestAccessor>
2715 fourierTransformRealOO(triple<SrcTraverser, SrcTraverser, SrcAccessor> src,
2716 pair<DestTraverser, DestAccessor> dest, fftw_real norm)
2718 fourierTransformRealOO(src.first, src.second, src.third,
2719 dest.first, dest.second, norm);
2722 template <
class SrcTraverser,
class SrcAccessor,
2723 class DestTraverser,
class DestAccessor>
2725 fourierTransformRealOO(SrcTraverser sul, SrcTraverser slr, SrcAccessor src,
2726 DestTraverser dul, DestAccessor dest, fftw_real norm)
2728 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest,
2729 norm, FFTW_RODFT00, FFTW_RODFT00);
2732 template <
class DestTraverser,
class DestAccessor>
2734 fourierTransformRealOO(
2738 DestTraverser dul, DestAccessor dest, fftw_real norm)
2740 int w = slr.x - sul.x;
2743 if (&(*(sul +
Diff2D(w, 0))) == &(*(sul +
Diff2D(0, 1))))
2744 fourierTransformRealImpl(sul, slr, dul, dest,
2745 norm, FFTW_RODFT00, FFTW_RODFT00);
2747 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest,
2748 norm, FFTW_RODFT00, FFTW_RODFT00);
2753 template <
class SrcTraverser,
class SrcAccessor,
2754 class DestTraverser,
class DestAccessor>
2756 fourierTransformRealWorkImageImpl(SrcTraverser sul, SrcTraverser slr, SrcAccessor src,
2757 DestTraverser dul, DestAccessor dest,
2758 fftw_real norm, fftw_r2r_kind kindx, fftw_r2r_kind kindy)
2760 FFTWRealImage workImage(slr - sul);
2761 copyImage(srcIterRange(sul, slr, src), destImage(workImage));
2762 FFTWRealImage
const & cworkImage = workImage;
2764 dul, dest,
norm, kindx, kindy);
2768 template <
class DestTraverser,
class DestAccessor>
2770 fourierTransformRealImpl(
2773 DestTraverser dul, DestAccessor dest,
2774 fftw_real norm, fftw_r2r_kind kindx, fftw_r2r_kind kindy)
2776 int w = slr.x - sul.x;
2777 int h = slr.y - sul.y;
2780 fftw_plan plan = fftw_plan_r2r_2d(h, w,
2781 (fftw_real *)&(*sul), (fftw_real *)res.
begin(),
2782 kindy, kindx, FFTW_ESTIMATE);
2784 fftw_destroy_plan(plan);
2788 std::bind1st(std::multiplies<fftw_real>(), 1.0 / norm));
2790 copyImage(srcImageRange(res), destIter(dul, dest));
2798 #endif // VIGRA_FFTW3_HXX FFTWComplex(fftwf_complex const &o)
Definition: fftw3.hxx:202
Size2D imageSize() const
Definition: imagecontainer.hxx:411
iterator begin()
Definition: basicimage.hxx:963
Definition: fftw3.hxx:1217
FixedPoint16< 2, OverflowHandling > atan2(FixedPoint16< IntBits, OverflowHandling > y, FixedPoint16< IntBits, OverflowHandling > x)
Arctangent. Accuracy better than 1/3 degree (9 significant bits).
Definition: fixedpoint.hxx:1654
BasicImage< fftw_real > FFTWRealImage
Definition: fftw3.hxx:1132
Fundamental class template for arrays of equal-sized images.
Definition: imagecontainer.hxx:72
Real value_type
Definition: fftw3.hxx:140
value_type operator()(ITERATOR const &i, DIFFERENCE d) const
Read phase at offset from iterator position.
Definition: fftw3.hxx:1443
Export associated information for each image iterator.
Definition: iteratortraits.hxx:109
FFTWComplex< R > conj(const FFTWComplex< R > &a)
complex conjugate
Definition: fftw3.hxx:1030
R imag(const FFTWComplex< R > &a)
imaginary part
Definition: fftw3.hxx:1023
FFTWComplex(value_type const &re=0.0, value_type const &im=0.0)
Definition: fftw3.hxx:169
value_type SquaredNormType
Definition: fftw3.hxx:164
Accessor accessor()
Definition: basicimage.hxx:1064
R arg(const FFTWComplex< R > &a)
phase
Definition: fftw3.hxx:1009
FFTWReal2Complex< Real >::type complex_type
Definition: fftw3.hxx:136
Definition: fftw3.hxx:1351
SquaredNormType squaredMagnitude() const
Definition: fftw3.hxx:358
Diff2D operator-(Diff2D const &a, Diff2D const &b)
Definition: diff2d.hxx:711
Real value_type
The accessor's value type.
Definition: fftw3.hxx:1381
linalg::TemporaryMatrix< T > sin(MultiArrayView< 2, T, C > const &v)
Real value_type
The accessor's value type.
Definition: fftw3.hxx:1324
Definition: fftw3.hxx:1319
value_type phase() const
Definition: fftw3.hxx:368
linalg::TemporaryMatrix< T > exp(MultiArrayView< 2, T, C > const &v)
FFTWComplex & operator=(fftw_complex const &o)
Definition: fftw3.hxx:255
FFTWComplex & operator=(TinyVector< T, 2 > const &o)
Definition: fftw3.hxx:310
value_type const * const_iterator
Definition: fftw3.hxx:156
Real value_type
The accessor's value type.
Definition: fftw3.hxx:1355
Definition: array_vector.hxx:954
value_type operator()(ITERATOR const &i, DIFFERENCE d) const
Read natural log of magnitude at offset from iterator position.
Definition: fftw3.hxx:1417
FFTWComplex< R >::SquaredNormType squaredNorm(const FFTWComplex< R > &a)
squared norm (= squared magnitude)
Definition: fftw3.hxx:1044
FFTWComplex & operator=(fftwl_complex const &o)
Definition: fftw3.hxx:273
int size() const
Definition: fftw3.hxx:383
R real(const FFTWComplex< R > &a)
real part
Definition: fftw3.hxx:1016
Diff2D operator+(Diff2D const &a, Diff2D const &b)
Definition: diff2d.hxx:739
Two dimensional difference vector.
Definition: diff2d.hxx:185
Real value_type
The accessor's value type.
Definition: fftw3.hxx:1222
FFTWComplex & operator=(std::complex< T > const &o)
Definition: fftw3.hxx:320
Real value_type
The accessor's value type.
Definition: fftw3.hxx:1272
FFTWComplex & operator=(FFTWComplex const &o)
Definition: fftw3.hxx:236
Definition: accessor.hxx:43
std::ptrdiff_t height() const
Definition: basicimage.hxx:845
value_type operator()(ITERATOR const &i) const
Read squared magnitude at iterator position.
Definition: fftw3.hxx:1359
FFTWComplex< R > & operator-=(FFTWComplex< R > &a, const FFTWComplex< R > &b)
subtract-assignment
Definition: fftw3.hxx:867
value_type operator()(ITERATOR const &i, DIFFERENCE d) const
Read squared magnitude at offset from iterator position.
Definition: fftw3.hxx:1365
Definition: fftw3.hxx:1429
doxygen_overloaded_function(template<... > void separableConvolveBlockwise) template< unsigned int N
Separated convolution on ChunkedArrays.
value_type operator()(ITERATOR const &i) const
Read natural log of magnitude at iterator position.
Definition: fftw3.hxx:1411
Accessor for items that are STL compatible vectors.
Definition: accessor.hxx:771
FFTWComplex< R >::NormType norm(const FFTWComplex< R > &a)
norm (= magnitude)
Definition: fftw3.hxx:1037
FFTWComplex & operator=(float o)
Definition: fftw3.hxx:291
Definition: basicimage.hxx:262
value_type operator()(ITERATOR const &i, DIFFERENCE d) const
Read real part at offset from iterator position.
Definition: fftw3.hxx:1232
FFTWComplex< R > & operator+=(FFTWComplex< R > &a, const FFTWComplex< R > &b)
add-assignment
Definition: fftw3.hxx:859
linalg::TemporaryMatrix< T > log10(MultiArrayView< 2, T, C > const &v)
FFTWComplex(FFTWComplex const &o)
Definition: fftw3.hxx:177
reference operator[](int i)
Definition: fftw3.hxx:373
value_type operator()(ITERATOR const &i) const
Read imaginary part at iterator position.
Definition: fftw3.hxx:1276
value_type const & const_reference
Definition: fftw3.hxx:148
bool operator!=(FFTWComplex< R > const &a, const FFTWComplex< R > &b)
not equal
Definition: fftw3.hxx:841
value_type operator()(ITERATOR const &i, DIFFERENCE d) const
Read imaginary part at offset from iterator position.
Definition: fftw3.hxx:1282
FFTWComplex(FFTWComplex< U > const &o)
Definition: fftw3.hxx:186
value_type operator()(ITERATOR const &i) const
Read phase at iterator position.
Definition: fftw3.hxx:1437
FFTWComplex & operator=(long double o)
Definition: fftw3.hxx:300
bool operator==(FFTWComplex< R > const &a, const FFTWComplex< R > &b)
equal
Definition: fftw3.hxx:825
FFTWComplex & operator=(FFTWComplex< U > const &o)
Definition: fftw3.hxx:246
void combineTwoImages(...)
Combine two source images into destination image.
FFTWComplex(fftw_complex const &o)
Definition: fftw3.hxx:194
size_type size() const
Definition: imagecontainer.hxx:228
void copyImage(...)
Copy source image into destination image.
FFTWComplex< R > & operator*=(FFTWComplex< R > &a, const FFTWComplex< R > &b)
multiply-assignment
Definition: fftw3.hxx:875
FFTWComplex(fftwl_complex const &o)
Definition: fftw3.hxx:210
value_type operator()(ITERATOR const &i) const
Read real part at iterator position.
Definition: fftw3.hxx:1226
Class for fixed size vectors.This class contains an array of size SIZE of the specified VALUETYPE...
Definition: accessor.hxx:940
value_type operator()(ITERATOR const &i) const
Read magnitude at iterator position.
Definition: fftw3.hxx:1385
value_type NormType
Definition: fftw3.hxx:160
Definition: fftw3.hxx:1268
traverser lowerRight()
Definition: basicimage.hxx:934
Definition: fftw3.hxx:1377
linalg::TemporaryMatrix< T > log(MultiArrayView< 2, T, C > const &v)
Fundamental class template for images.
Definition: basicimage.hxx:473
Real value_type
The accessor's value type.
Definition: fftw3.hxx:1407
FFTWComplex operator-() const
Definition: fftw3.hxx:353
Definition: basicimage.hxx:294
FFTWComplex< R >::NormType abs(const FFTWComplex< R > &a)
absolute value (= magnitude)
Definition: fftw3.hxx:1002
void resize(size_type newSize)
Definition: imagecontainer.hxx:310
void resize(std::ptrdiff_t width, std::ptrdiff_t height)
Definition: basicimage.hxx:776
linalg::TemporaryMatrix< T > tan(MultiArrayView< 2, T, C > const &v)
value_type * iterator
Definition: fftw3.hxx:152
const_reference operator[](int i) const
Definition: fftw3.hxx:378
FFTWComplex< R > & operator/=(FFTWComplex< R > &a, const FFTWComplex< R > &b)
divide-assignment
Definition: fftw3.hxx:884
NormType magnitude() const
Definition: fftw3.hxx:363
Encapsulate access to the values an iterator points to.
Definition: accessor.hxx:133
Definition: fftw3.hxx:1403
linalg::TemporaryMatrix< T > cos(MultiArrayView< 2, T, C > const &v)
std::ptrdiff_t width() const
Definition: basicimage.hxx:838
T sign(T t)
The sign function.
Definition: mathutil.hxx:574
virtual void resizeImages(const Diff2D &newSize)
Definition: imagecontainer.hxx:418
Real value_type
The accessor's value type.
Definition: fftw3.hxx:1433
value_type operator()(ITERATOR const &i, DIFFERENCE d) const
Read magnitude at offset from iterator position.
Definition: fftw3.hxx:1391
Wrapper class for the FFTW complex types 'fftw_complex'.
Definition: fftw3.hxx:131
SquareRootTraits< FixedPoint< IntBits, FracBits > >::SquareRootResult sqrt(FixedPoint< IntBits, FracBits > v)
square root.
Definition: fixedpoint.hxx:616
value_type & reference
Definition: fftw3.hxx:144
FFTWComplex(std::complex< T > const &o)
Definition: fftw3.hxx:219
BasicImage< FFTWComplex<> > FFTWComplexImage
Definition: fftw3.hxx:1192
FFTWComplex(TinyVector< T, 2 > const &o)
Definition: fftw3.hxx:228
FFTWComplex & operator=(double o)
Definition: fftw3.hxx:282
traverser upperLeft()
Definition: basicimage.hxx:923
FFTWComplex & operator=(fftwf_complex const &o)
Definition: fftw3.hxx:264