36 #ifndef VIGRA_FFTW3_HXX
37 #define VIGRA_FFTW3_HXX
41 #include "stdimage.hxx"
42 #include "copyimage.hxx"
43 #include "transformimage.hxx"
44 #include "combineimages.hxx"
45 #include "numerictraits.hxx"
46 #include "imagecontainer.hxx"
51 typedef double fftw_real;
131 data_[0] = o.data_[0];
132 data_[1] = o.data_[1];
156 data_[0] = o.data_[0];
157 data_[1] = o.data_[1];
209 {
return data_[0]*data_[0]+data_[1]*data_[1]; }
240 {
return data_ + 2; }
246 {
return data_ + 2; }
355 struct NumericTraits<fftw_complex>
357 typedef fftw_complex Type;
358 typedef fftw_complex Promote;
359 typedef fftw_complex RealPromote;
360 typedef fftw_complex ComplexPromote;
361 typedef fftw_real ValueType;
363 typedef VigraFalseType isIntegral;
364 typedef VigraFalseType isScalar;
365 typedef NumericTraits<fftw_real>::isSigned isSigned;
366 typedef VigraFalseType isOrdered;
367 typedef VigraTrueType isComplex;
369 static FFTWComplex zero() {
return FFTWComplex(0.0, 0.0); }
370 static FFTWComplex one() {
return FFTWComplex(1.0, 0.0); }
371 static FFTWComplex nonZero() {
return one(); }
373 static const Promote & toPromote(
const Type & v) {
return v; }
374 static const RealPromote & toRealPromote(
const Type & v) {
return v; }
375 static const Type & fromPromote(
const Promote & v) {
return v; }
376 static const Type & fromRealPromote(
const RealPromote & v) {
return v; }
380 struct NumericTraits<FFTWComplex>
382 typedef FFTWComplex Type;
383 typedef FFTWComplex Promote;
384 typedef FFTWComplex RealPromote;
385 typedef FFTWComplex ComplexPromote;
386 typedef fftw_real ValueType;
388 typedef VigraFalseType isIntegral;
389 typedef VigraFalseType isScalar;
390 typedef NumericTraits<fftw_real>::isSigned isSigned;
391 typedef VigraFalseType isOrdered;
392 typedef VigraTrueType isComplex;
394 static FFTWComplex zero() {
return FFTWComplex(0.0, 0.0); }
395 static FFTWComplex one() {
return FFTWComplex(1.0, 0.0); }
396 static FFTWComplex nonZero() {
return one(); }
398 static const Promote & toPromote(
const Type & v) {
return v; }
399 static const RealPromote & toRealPromote(
const Type & v) {
return v; }
400 static const Type & fromPromote(
const Promote & v) {
return v; }
401 static const Type & fromRealPromote(
const RealPromote & v) {
return v; }
405 struct NormTraits<fftw_complex>
407 typedef fftw_complex Type;
408 typedef fftw_real SquaredNormType;
409 typedef fftw_real NormType;
413 struct NormTraits<FFTWComplex>
415 typedef FFTWComplex Type;
416 typedef fftw_real SquaredNormType;
417 typedef fftw_real NormType;
421 struct PromoteTraits<fftw_complex, fftw_complex>
423 typedef fftw_complex Promote;
427 struct PromoteTraits<fftw_complex, double>
429 typedef fftw_complex Promote;
433 struct PromoteTraits<double, fftw_complex>
435 typedef fftw_complex Promote;
439 struct PromoteTraits<FFTWComplex, FFTWComplex>
441 typedef FFTWComplex Promote;
445 struct PromoteTraits<FFTWComplex, double>
447 typedef FFTWComplex Promote;
451 struct PromoteTraits<double, FFTWComplex>
453 typedef FFTWComplex Promote;
457 struct CanSkipInitialization<std::complex<T> >
459 typedef typename CanSkipInitialization<T>::type type;
460 static const bool value = type::asBool;
464 struct CanSkipInitialization<FFTWComplex>
466 typedef CanSkipInitialization<fftw_real>::type type;
467 static const bool value = type::asBool;
493 return a.re() == b.re() && a.im() == b.im();
498 return a.re() != b.re() || a.im() != b.im();
518 a.im() = a.re()*b.im()+a.im()*b.re();
527 a.im() = (b.re()*a.im()-a.re()*b.im())/sm;
651 typedef Iterator iterator;
654 typedef iterator::iterator_category iterator_category;
655 typedef iterator::value_type value_type;
656 typedef iterator::reference reference;
657 typedef iterator::index_reference index_reference;
658 typedef iterator::pointer pointer;
659 typedef iterator::difference_type difference_type;
660 typedef iterator::row_iterator row_iterator;
661 typedef iterator::column_iterator column_iterator;
664 typedef VigraTrueType hasConstantStrides;
668 struct IteratorTraits<
669 ConstBasicImageIterator<FFTWComplex, FFTWComplex **> >
671 typedef ConstBasicImageIterator<FFTWComplex, FFTWComplex **> Iterator;
672 typedef Iterator iterator;
673 typedef BasicImageIterator<FFTWComplex, FFTWComplex **> mutable_iterator;
674 typedef ConstBasicImageIterator<FFTWComplex, FFTWComplex **> const_iterator;
675 typedef iterator::iterator_category iterator_category;
676 typedef iterator::value_type value_type;
677 typedef iterator::reference reference;
678 typedef iterator::index_reference index_reference;
679 typedef iterator::pointer pointer;
680 typedef iterator::difference_type difference_type;
681 typedef iterator::row_iterator row_iterator;
682 typedef iterator::column_iterator column_iterator;
683 typedef VectorAccessor<FFTWComplex> default_accessor;
684 typedef VectorAccessor<FFTWComplex> DefaultAccessor;
685 typedef VigraTrueType hasConstantStrides;
728 template <
class ITERATOR>
734 template <
class ITERATOR,
class DIFFERENCE>
740 template <
class ITERATOR>
746 template <
class ITERATOR,
class DIFFERENCE>
747 void set(
value_type const & v, ITERATOR
const & i, DIFFERENCE d)
const {
752 template <
class ITERATOR>
758 template <
class ITERATOR,
class DIFFERENCE>
759 void set(
FFTWComplex const & v, ITERATOR
const & i, DIFFERENCE d)
const {
777 template <
class ITERATOR>
783 template <
class ITERATOR,
class DIFFERENCE>
789 template <
class ITERATOR>
795 template <
class ITERATOR,
class DIFFERENCE>
796 void set(
value_type const & v, ITERATOR
const & i, DIFFERENCE d)
const {
801 template <
class ITERATOR>
807 template <
class ITERATOR,
class DIFFERENCE>
808 void set(
FFTWComplex const & v, ITERATOR
const & i, DIFFERENCE d)
const {
829 template <
class ITERATOR>
838 template <
class ITERATOR,
class DIFFERENCE>
839 void set(
value_type const & v, ITERATOR
const & i, DIFFERENCE d)
const {
858 template <
class ITERATOR>
860 return (*i).magnitude();
864 template <
class ITERATOR,
class DIFFERENCE>
866 return (i[d]).magnitude();
883 template <
class ITERATOR>
889 template <
class ITERATOR,
class DIFFERENCE>
891 return (i[d]).phase();
1050 template <
class SrcImageIterator,
class SrcAccessor,
1051 class DestImageIterator,
class DestAccessor>
1053 SrcImageIterator src_lowerright, SrcAccessor sa,
1054 DestImageIterator dest_upperleft, DestAccessor da)
1056 int w = int(src_lowerright.x - src_upperleft.x);
1057 int h = int(src_lowerright.y - src_upperleft.y);
1065 src_upperleft + Diff2D(w2, h2), sa),
1066 destIter (dest_upperleft + Diff2D(w1, h1), da));
1069 copyImage(srcIterRange(src_upperleft + Diff2D(w2, h2),
1070 src_lowerright, sa),
1071 destIter (dest_upperleft, da));
1074 copyImage(srcIterRange(src_upperleft + Diff2D(w2, 0),
1075 src_upperleft + Diff2D(w, h2), sa),
1076 destIter (dest_upperleft + Diff2D(0, h1), da));
1079 copyImage(srcIterRange(src_upperleft + Diff2D(0, h2),
1080 src_upperleft + Diff2D(w2, h), sa),
1081 destIter (dest_upperleft + Diff2D(w1, 0), da));
1084 template <
class SrcImageIterator,
class SrcAccessor,
1085 class DestImageIterator,
class DestAccessor>
1087 triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src,
1088 pair<DestImageIterator, DestAccessor> dest)
1091 dest.first, dest.second);
1133 template <
class SrcImageIterator,
class SrcAccessor,
1134 class DestImageIterator,
class DestAccessor>
1136 SrcImageIterator src_lowerright, SrcAccessor sa,
1137 DestImageIterator dest_upperleft, DestAccessor da)
1139 int w = int(src_lowerright.x - src_upperleft.x);
1140 int h = int(src_lowerright.y - src_upperleft.y);
1148 src_upperleft + Diff2D(w2, h2), sa),
1149 destIter (dest_upperleft + Diff2D(w1, h1), da));
1152 copyImage(srcIterRange(src_upperleft + Diff2D(w2, h2),
1153 src_lowerright, sa),
1154 destIter (dest_upperleft, da));
1157 copyImage(srcIterRange(src_upperleft + Diff2D(w2, 0),
1158 src_upperleft + Diff2D(w, h2), sa),
1159 destIter (dest_upperleft + Diff2D(0, h1), da));
1162 copyImage(srcIterRange(src_upperleft + Diff2D(0, h2),
1163 src_upperleft + Diff2D(w2, h), sa),
1164 destIter (dest_upperleft + Diff2D(w1, 0), da));
1167 template <
class SrcImageIterator,
class SrcAccessor,
1168 class DestImageIterator,
class DestAccessor>
1170 triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src,
1171 pair<DestImageIterator, DestAccessor> dest)
1174 dest.first, dest.second);
1177 template <
class DestImageIterator,
class DestAccessor>
1178 void fftShift(DestImageIterator upperleft,
1179 DestImageIterator lowerright, DestAccessor da)
1181 int w = int(lowerright.x - upperleft.x);
1182 int h = int(lowerright.y - upperleft.y);
1189 swapImageData(destIterRange(upperleft,
1190 upperleft + Diff2D(w2, h2), da),
1191 destIter (upperleft + Diff2D(w1, h1), da));
1194 swapImageData(destIterRange(upperleft + Diff2D(w2, 0),
1195 upperleft + Diff2D(w, h2), da),
1196 destIter (upperleft + Diff2D(0, h1), da));
1199 template <
class DestImageIterator,
class DestAccessor>
1200 inline void fftShift(
1201 triple<DestImageIterator, DestImageIterator, DestAccessor> dest)
1203 fftShift(dest.first, dest.second, dest.third);
1215 int w = int(slr.x - sul.x);
1216 int h = int(slr.y - sul.y);
1220 fftw_complex * srcPtr = (fftw_complex *)(&*sul);
1221 fftw_complex * destPtr = (fftw_complex *)(&*dul);
1224 if (h > 1 && &(*(sul + Diff2D(w, 0))) != &(*(sul + Diff2D(0, 1))))
1226 sworkImage.resize(w, h);
1227 copyImage(srcIterRange(sul, slr, src), destImage(sworkImage));
1228 srcPtr = (fftw_complex *)(&(*sworkImage.upperLeft()));
1230 if (h > 1 && &(*(dul + Diff2D(w, 0))) != &(*(dul + Diff2D(0, 1))))
1232 dworkImage.resize(w, h);
1233 destPtr = (fftw_complex *)(&(*dworkImage.upperLeft()));
1236 fftw_plan plan = fftw_plan_dft_2d(h, w, srcPtr, destPtr, sign, FFTW_ESTIMATE );
1238 fftw_destroy_plan(plan);
1240 if (h > 1 && &(*(dul + Diff2D(w, 0))) != &(*(dul + Diff2D(0, 1))))
1242 copyImage(srcImageRange(dworkImage), destIter(dul, dest));
1316 detail::fourierTransformImpl(sul, slr, src, dul, dest, FFTW_FORWARD);
1319 template <
class SrcImageIterator,
class SrcAccessor>
1321 SrcImageIterator srcLowerRight, SrcAccessor sa,
1325 int w= srcLowerRight.x - srcUpperLeft.x;
1326 int h= srcLowerRight.y - srcUpperLeft.y;
1329 copyImage(srcIterRange(srcUpperLeft, srcLowerRight, sa),
1330 destImage(workImage, FFTWWriteRealAccessor()));
1334 fourierTransform(cworkImage.upperLeft(), cworkImage.lowerRight(), cworkImage.accessor(),
1338 template <
class SrcImageIterator,
class SrcAccessor>
1340 void fourierTransform(triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src,
1341 pair<FFTWComplexImage::traverser, FFTWComplexImage::Accessor> dest)
1343 fourierTransform(src.first, src.second, src.third, dest.first, dest.second);
1355 detail::fourierTransformImpl(sul, slr, src, dul, dest, FFTW_BACKWARD);
1358 template <
class DestImageIterator,
class DestAccessor>
1361 DestImageIterator dul, DestAccessor dest)
1363 int w = slr.x - sul.x;
1364 int h = slr.y - sul.y;
1368 copyImage(srcImageRange(workImage), destIter(dul, dest));
1372 template <
class DestImageIterator,
class DestAccessor>
1376 pair<DestImageIterator, DestAccessor> dest)
1463 template <
class SrcImageIterator,
class SrcAccessor,
1464 class FilterImageIterator,
class FilterAccessor,
1465 class DestImageIterator,
class DestAccessor>
1467 SrcImageIterator srcLowerRight, SrcAccessor sa,
1468 FilterImageIterator filterUpperLeft, FilterAccessor fa,
1469 DestImageIterator destUpperLeft, DestAccessor da)
1472 int w = int(srcLowerRight.x - srcUpperLeft.x);
1473 int h = int(srcLowerRight.y - srcUpperLeft.y);
1476 copyImage(srcIterRange(srcUpperLeft, srcLowerRight, sa),
1477 destImage(workImage, FFTWWriteRealAccessor()));
1481 applyFourierFilterImpl(cworkImage.upperLeft(), cworkImage.lowerRight(), cworkImage.accessor(),
1482 filterUpperLeft, fa,
1486 template <
class FilterImageIterator,
class FilterAccessor,
1487 class DestImageIterator,
class DestAccessor>
1493 FilterImageIterator filterUpperLeft, FilterAccessor fa,
1494 DestImageIterator destUpperLeft, DestAccessor da)
1496 int w = srcLowerRight.x - srcUpperLeft.x;
1497 int h = srcLowerRight.y - srcUpperLeft.y;
1500 if (&(*(srcUpperLeft + Diff2D(w, 0))) == &(*(srcUpperLeft + Diff2D(0, 1))))
1501 applyFourierFilterImpl(srcUpperLeft, srcLowerRight, sa,
1502 filterUpperLeft, fa,
1507 copyImage(srcIterRange(srcUpperLeft, srcLowerRight, sa),
1508 destImage(workImage));
1511 applyFourierFilterImpl(cworkImage.upperLeft(), cworkImage.lowerRight(), cworkImage.accessor(),
1512 filterUpperLeft, fa,
1517 template <
class SrcImageIterator,
class SrcAccessor,
1518 class FilterImageIterator,
class FilterAccessor,
1519 class DestImageIterator,
class DestAccessor>
1522 pair<FilterImageIterator, FilterAccessor> filter,
1523 pair<DestImageIterator, DestAccessor> dest)
1526 filter.first, filter.second,
1527 dest.first, dest.second);
1530 template <
class FilterImageIterator,
class FilterAccessor,
1531 class DestImageIterator,
class DestAccessor>
1532 void applyFourierFilterImpl(
1536 FilterImageIterator filterUpperLeft, FilterAccessor fa,
1537 DestImageIterator destUpperLeft, DestAccessor da)
1539 int w = int(srcLowerRight.x - srcUpperLeft.x);
1540 int h = int(srcLowerRight.y - srcUpperLeft.y);
1545 fftw_plan forwardPlan=
1546 fftw_plan_dft_2d(h, w, (fftw_complex *)&(*srcUpperLeft),
1547 (fftw_complex *)complexResultImg.begin(),
1548 FFTW_FORWARD, FFTW_ESTIMATE );
1549 fftw_execute(forwardPlan);
1550 fftw_destroy_plan(forwardPlan);
1553 combineTwoImages(srcImageRange(complexResultImg), srcIter(filterUpperLeft, fa),
1554 destImage(complexResultImg), std::multiplies<FFTWComplex>());
1557 fftw_plan backwardPlan=
1558 fftw_plan_dft_2d(h, w, (fftw_complex *)complexResultImg.begin(),
1559 (fftw_complex *)complexResultImg.begin(),
1560 FFTW_BACKWARD, FFTW_ESTIMATE);
1561 fftw_execute(backwardPlan);
1562 fftw_destroy_plan(backwardPlan);
1565 NumericTraits<typename DestAccessor::value_type>::isScalar
1569 applyFourierFilterImplNormalization(complexResultImg, destUpperLeft, da,
1573 template <
class DestImageIterator,
class DestAccessor>
1575 DestImageIterator destUpperLeft,
1579 double normFactor= 1.0/(srcImage.width() * srcImage.height());
1581 for(
int y=0; y<srcImage.height(); y++, destUpperLeft.y++)
1583 DestImageIterator dIt= destUpperLeft;
1584 for(
int x= 0; x< srcImage.width(); x++, dIt.x++)
1586 da.setComponent(srcImage(x, y).re()*normFactor, dIt, 0);
1587 da.setComponent(srcImage(x, y).im()*normFactor, dIt, 1);
1593 void applyFourierFilterImplNormalization(
FFTWComplexImage const & srcImage,
1598 transformImage(srcImageRange(srcImage), destIter(destUpperLeft, da),
1599 linearIntensityTransform<FFTWComplex>(1.0/(srcImage.width() * srcImage.height())));
1602 template <
class DestImageIterator,
class DestAccessor>
1603 void applyFourierFilterImplNormalization(
FFTWComplexImage const & srcImage,
1604 DestImageIterator destUpperLeft,
1608 double normFactor= 1.0/(srcImage.width() * srcImage.height());
1610 for(
int y=0; y<srcImage.height(); y++, destUpperLeft.y++)
1612 DestImageIterator dIt= destUpperLeft;
1613 for(
int x= 0; x< srcImage.width(); x++, dIt.x++)
1614 da.set(srcImage(x, y).re()*normFactor, dIt);
1679 template <
class SrcImageIterator,
class SrcAccessor,
1680 class FilterType,
class DestImage>
1683 const ImageArray<FilterType> &filters,
1684 ImageArray<DestImage> &results)
1690 template <
class SrcImageIterator,
class SrcAccessor,
1691 class FilterType,
class DestImage>
1693 SrcImageIterator srcLowerRight, SrcAccessor sa,
1694 const ImageArray<FilterType> &filters,
1695 ImageArray<DestImage> &results)
1697 int w = int(srcLowerRight.x - srcUpperLeft.x);
1698 int h = int(srcLowerRight.y - srcUpperLeft.y);
1701 copyImage(srcIterRange(srcUpperLeft, srcLowerRight, sa),
1702 destImage(workImage, FFTWWriteRealAccessor()));
1705 applyFourierFilterFamilyImpl(cworkImage.upperLeft(), cworkImage.lowerRight(), cworkImage.accessor(),
1709 template <
class FilterType,
class DestImage>
1715 const ImageArray<FilterType> &filters,
1716 ImageArray<DestImage> &results)
1718 int w= srcLowerRight.x - srcUpperLeft.x;
1721 if (&(*(srcUpperLeft + Diff2D(w, 0))) == &(*(srcUpperLeft + Diff2D(0, 1))))
1722 applyFourierFilterFamilyImpl(srcUpperLeft, srcLowerRight, sa,
1726 int h = srcLowerRight.y - srcUpperLeft.y;
1728 copyImage(srcIterRange(srcUpperLeft, srcLowerRight, sa),
1729 destImage(workImage));
1732 applyFourierFilterFamilyImpl(cworkImage.upperLeft(), cworkImage.lowerRight(), cworkImage.accessor(),
1737 template <
class FilterType,
class DestImage>
1738 void applyFourierFilterFamilyImpl(
1742 const ImageArray<FilterType> &filters,
1743 ImageArray<DestImage> &results)
1749 vigra_precondition((srcLowerRight - srcUpperLeft) == filters.imageSize(),
1750 "applyFourierFilterFamily called with src image size != filters.imageSize()!");
1753 results.resize(filters.size());
1754 results.resizeImages(filters.imageSize());
1757 int w = int(srcLowerRight.x - srcUpperLeft.x);
1758 int h = int(srcLowerRight.y - srcUpperLeft.y);
1763 fftw_plan forwardPlan=
1764 fftw_plan_dft_2d(h, w, (fftw_complex *)&(*srcUpperLeft),
1765 (fftw_complex *)freqImage.begin(),
1766 FFTW_FORWARD, FFTW_ESTIMATE );
1767 fftw_execute(forwardPlan);
1768 fftw_destroy_plan(forwardPlan);
1770 fftw_plan backwardPlan=
1771 fftw_plan_dft_2d(h, w, (fftw_complex *)result.begin(),
1772 (fftw_complex *)result.begin(),
1773 FFTW_BACKWARD, FFTW_ESTIMATE );
1775 NumericTraits<typename DestImage::Accessor::value_type>::isScalar
1779 for (
unsigned int i= 0; i < filters.size(); i++)
1782 destImage(result), std::multiplies<FFTWComplex>());
1785 fftw_execute(backwardPlan);
1788 applyFourierFilterImplNormalization(result,
1789 results[i].upperLeft(), results[i].accessor(),
1792 fftw_destroy_plan(backwardPlan);
1914 template <
class SrcTraverser,
class SrcAccessor,
1915 class DestTraverser,
class DestAccessor>
1917 fourierTransformRealEE(triple<SrcTraverser, SrcTraverser, SrcAccessor> src,
1918 pair<DestTraverser, DestAccessor> dest, fftw_real
norm)
1920 fourierTransformRealEE(src.first, src.second, src.third,
1921 dest.first, dest.second, norm);
1924 template <
class SrcTraverser,
class SrcAccessor,
1925 class DestTraverser,
class DestAccessor>
1927 fourierTransformRealEE(SrcTraverser sul, SrcTraverser slr, SrcAccessor src,
1928 DestTraverser dul, DestAccessor dest, fftw_real
norm)
1930 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest,
1931 norm, FFTW_REDFT00, FFTW_REDFT00);
1934 template <
class DestTraverser,
class DestAccessor>
1936 fourierTransformRealEE(
1940 DestTraverser dul, DestAccessor dest, fftw_real
norm)
1942 int w = slr.x - sul.x;
1945 if (&(*(sul + Diff2D(w, 0))) == &(*(sul + Diff2D(0, 1))))
1946 fourierTransformRealImpl(sul, slr, dul, dest,
1947 norm, FFTW_REDFT00, FFTW_REDFT00);
1949 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest,
1950 norm, FFTW_REDFT00, FFTW_REDFT00);
1955 template <
class SrcTraverser,
class SrcAccessor,
1956 class DestTraverser,
class DestAccessor>
1958 fourierTransformRealOE(triple<SrcTraverser, SrcTraverser, SrcAccessor> src,
1959 pair<DestTraverser, DestAccessor> dest, fftw_real norm)
1961 fourierTransformRealOE(src.first, src.second, src.third,
1962 dest.first, dest.second, norm);
1965 template <
class SrcTraverser,
class SrcAccessor,
1966 class DestTraverser,
class DestAccessor>
1968 fourierTransformRealOE(SrcTraverser sul, SrcTraverser slr, SrcAccessor src,
1969 DestTraverser dul, DestAccessor dest, fftw_real norm)
1971 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest,
1972 norm, FFTW_RODFT00, FFTW_REDFT00);
1975 template <
class DestTraverser,
class DestAccessor>
1977 fourierTransformRealOE(
1981 DestTraverser dul, DestAccessor dest, fftw_real norm)
1983 int w = slr.x - sul.x;
1986 if (&(*(sul + Diff2D(w, 0))) == &(*(sul + Diff2D(0, 1))))
1987 fourierTransformRealImpl(sul, slr, dul, dest,
1988 norm, FFTW_RODFT00, FFTW_REDFT00);
1990 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest,
1991 norm, FFTW_RODFT00, FFTW_REDFT00);
1996 template <
class SrcTraverser,
class SrcAccessor,
1997 class DestTraverser,
class DestAccessor>
1999 fourierTransformRealEO(triple<SrcTraverser, SrcTraverser, SrcAccessor> src,
2000 pair<DestTraverser, DestAccessor> dest, fftw_real norm)
2002 fourierTransformRealEO(src.first, src.second, src.third,
2003 dest.first, dest.second, norm);
2006 template <
class SrcTraverser,
class SrcAccessor,
2007 class DestTraverser,
class DestAccessor>
2009 fourierTransformRealEO(SrcTraverser sul, SrcTraverser slr, SrcAccessor src,
2010 DestTraverser dul, DestAccessor dest, fftw_real norm)
2012 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest,
2013 norm, FFTW_REDFT00, FFTW_RODFT00);
2016 template <
class DestTraverser,
class DestAccessor>
2018 fourierTransformRealEO(
2022 DestTraverser dul, DestAccessor dest, fftw_real norm)
2024 int w = slr.x - sul.x;
2027 if (&(*(sul + Diff2D(w, 0))) == &(*(sul + Diff2D(0, 1))))
2028 fourierTransformRealImpl(sul, slr, dul, dest,
2029 norm, FFTW_REDFT00, FFTW_RODFT00);
2031 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest,
2032 norm, FFTW_REDFT00, FFTW_RODFT00);
2037 template <
class SrcTraverser,
class SrcAccessor,
2038 class DestTraverser,
class DestAccessor>
2040 fourierTransformRealOO(triple<SrcTraverser, SrcTraverser, SrcAccessor> src,
2041 pair<DestTraverser, DestAccessor> dest, fftw_real norm)
2043 fourierTransformRealOO(src.first, src.second, src.third,
2044 dest.first, dest.second, norm);
2047 template <
class SrcTraverser,
class SrcAccessor,
2048 class DestTraverser,
class DestAccessor>
2050 fourierTransformRealOO(SrcTraverser sul, SrcTraverser slr, SrcAccessor src,
2051 DestTraverser dul, DestAccessor dest, fftw_real norm)
2053 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest,
2054 norm, FFTW_RODFT00, FFTW_RODFT00);
2057 template <
class DestTraverser,
class DestAccessor>
2059 fourierTransformRealOO(
2063 DestTraverser dul, DestAccessor dest, fftw_real norm)
2065 int w = slr.x - sul.x;
2068 if (&(*(sul + Diff2D(w, 0))) == &(*(sul + Diff2D(0, 1))))
2069 fourierTransformRealImpl(sul, slr, dul, dest,
2070 norm, FFTW_RODFT00, FFTW_RODFT00);
2072 fourierTransformRealWorkImageImpl(sul, slr, src, dul, dest,
2073 norm, FFTW_RODFT00, FFTW_RODFT00);
2078 template <
class SrcTraverser,
class SrcAccessor,
2079 class DestTraverser,
class DestAccessor>
2081 fourierTransformRealWorkImageImpl(SrcTraverser sul, SrcTraverser slr, SrcAccessor src,
2082 DestTraverser dul, DestAccessor dest,
2083 fftw_real norm, fftw_r2r_kind kindx, fftw_r2r_kind kindy)
2086 copyImage(srcIterRange(sul, slr, src), destImage(workImage));
2088 fourierTransformRealImpl(cworkImage.upperLeft(), cworkImage.lowerRight(),
2089 dul, dest,
norm, kindx, kindy);
2093 template <
class DestTraverser,
class DestAccessor>
2095 fourierTransformRealImpl(
2098 DestTraverser dul, DestAccessor dest,
2099 fftw_real norm, fftw_r2r_kind kindx, fftw_r2r_kind kindy)
2101 int w = slr.x - sul.x;
2102 int h = slr.y - sul.y;
2103 BasicImage<fftw_real> res(w, h);
2105 fftw_plan plan = fftw_plan_r2r_2d(h, w,
2106 (fftw_real *)&(*sul), (fftw_real *)res.begin(),
2107 kindy, kindx, FFTW_ESTIMATE);
2109 fftw_destroy_plan(plan);
2113 std::bind1st(std::multiplies<fftw_real>(), 1.0 / norm));
2115 copyImage(srcImageRange(res), destIter(dul, dest));
2123 #endif // VIGRA_FFTW3_HXX