PCL
FFT2D.h
Go to the documentation of this file.
1 // ____ ______ __
2 // / __ \ / ____// /
3 // / /_/ // / / /
4 // / ____// /___ / /___ PixInsight Class Library
5 // /_/ \____//_____/ PCL 2.7.0
6 // ----------------------------------------------------------------------------
7 // pcl/FFT2D.h - Released 2024-06-18T15:48:54Z
8 // ----------------------------------------------------------------------------
9 // This file is part of the PixInsight Class Library (PCL).
10 // PCL is a multiplatform C++ framework for development of PixInsight modules.
11 //
12 // Copyright (c) 2003-2024 Pleiades Astrophoto S.L. All Rights Reserved.
13 //
14 // Redistribution and use in both source and binary forms, with or without
15 // modification, is permitted provided that the following conditions are met:
16 //
17 // 1. All redistributions of source code must retain the above copyright
18 // notice, this list of conditions and the following disclaimer.
19 //
20 // 2. All redistributions in binary form must reproduce the above copyright
21 // notice, this list of conditions and the following disclaimer in the
22 // documentation and/or other materials provided with the distribution.
23 //
24 // 3. Neither the names "PixInsight" and "Pleiades Astrophoto", nor the names
25 // of their contributors, may be used to endorse or promote products derived
26 // from this software without specific prior written permission. For written
27 // permission, please contact info@pixinsight.com.
28 //
29 // 4. All products derived from this software, in any form whatsoever, must
30 // reproduce the following acknowledgment in the end-user documentation
31 // and/or other materials provided with the product:
32 //
33 // "This product is based on software from the PixInsight project, developed
34 // by Pleiades Astrophoto and its contributors (https://pixinsight.com/)."
35 //
36 // Alternatively, if that is where third-party acknowledgments normally
37 // appear, this acknowledgment must be reproduced in the product itself.
38 //
39 // THIS SOFTWARE IS PROVIDED BY PLEIADES ASTROPHOTO AND ITS CONTRIBUTORS
40 // "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
41 // TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
42 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL PLEIADES ASTROPHOTO OR ITS
43 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
44 // EXEMPLARY OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, BUSINESS
45 // INTERRUPTION; PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; AND LOSS OF USE,
46 // DATA OR PROFITS) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
47 // CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
48 // ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
49 // POSSIBILITY OF SUCH DAMAGE.
50 // ----------------------------------------------------------------------------
51 
52 #ifndef __PCL_FFT2D_h
53 #define __PCL_FFT2D_h
54 
56 
57 /*
58  * The FFT routines used in this PCL version are based on the KISS FFT library
59  * by Mark Borgerding: http://sourceforge.net/projects/kissfft/
60  *
61  * KISS FFT LICENSE INFORMATION
62  * ============================
63  *
64  * Copyright (c) 2003-2009 Mark Borgerding
65  *
66  * All rights reserved.
67  *
68  * Redistribution and use in source and binary forms, with or without
69  * modification, are permitted provided that the following conditions are met:
70  *
71  * * Redistributions of source code must retain the above copyright notice,
72  * this list of conditions and the following disclaimer.
73  * * Redistributions in binary form must reproduce the above copyright notice,
74  * this list of conditions and the following disclaimer in the documentation
75  * and/or other materials provided with the distribution.
76  * * Neither the author nor the names of any contributors may be used to
77  * endorse or promote products derived from this software without specific
78  * prior written permission.
79  *
80  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
81  * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
82  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
83  * ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE
84  * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
85  * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
86  * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
87  * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
88  * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
89  * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
90  * POSSIBILITY OF SUCH DAMAGE.
91  */
92 
93 #include <pcl/Defs.h>
94 
95 #include <pcl/FFT1D.h>
96 #include <pcl/Matrix.h>
97 #include <pcl/ParallelProcess.h>
98 #include <pcl/StatusMonitor.h>
99 
100 namespace pcl
101 {
102 
103 // ----------------------------------------------------------------------------
104 
105 class FFT2DBase
106 {
107 protected:
108 
109  static void Transform( int, int, fcomplex*, const fcomplex*, int, StatusMonitor*, bool, int );
110  static void Transform( int, int, dcomplex*, const dcomplex*, int, StatusMonitor*, bool, int );
111  static void Transform( int, int, fcomplex*, const float*, StatusMonitor*, bool, int );
112  static void Transform( int, int, dcomplex*, const double*, StatusMonitor*, bool, int );
113  static void Transform( int, int, float*, const fcomplex*, StatusMonitor*, bool, int );
114  static void Transform( int, int, double*, const dcomplex*, StatusMonitor*, bool, int );
115 };
116 
117 // ----------------------------------------------------------------------------
118 
128 template <typename T>
129 class PCL_CLASS AbstractFFT2D : public FFT2DBase, public ParallelProcess
130 {
131 public:
132 
136  using scalar = T;
137 
142 
147 
152 
157 
162 
167 
172  AbstractFFT2D( int rows, int cols )
173  : m_rows( rows )
174  , m_cols( cols )
175  {
176  }
177 
185  AbstractFFT2D( int rows, int cols, StatusMonitor& monitor )
186  : m_rows( rows )
187  , m_cols( cols )
188  , m_monitor( &monitor )
189  {
190  }
191 
196  ~AbstractFFT2D() override
197  {
198  Release();
199  }
200 
206  int Rows() const
207  {
208  return m_rows;
209  }
210 
216  int Cols() const
217  {
218  return m_cols;
219  }
220 
226  int NumberOfElements() const
227  {
228  return m_rows*m_cols;
229  }
230 
242  transform DFT() const
243  {
244  return m_dft;
245  }
246 
259  {
260  return m_dft;
261  }
262 
267  virtual void Release()
268  {
269  m_dft = transform();
270  }
271 
272 private:
273 
274  // Transform dimensions
275  int m_rows;
276  int m_cols;
277 
278 protected:
279 
280  transform m_dft; // DFT of complex or real data
281  StatusMonitor* m_monitor = nullptr;
282 };
283 
284 // ----------------------------------------------------------------------------
285 
286 #define m_dft this->m_dft
287 #define m_rows this->Rows()
288 #define m_cols this->Cols()
289 #define m_monitor this->m_monitor
290 #define m_parallel this->m_parallel
291 #define m_maxProcessors this->m_maxProcessors
292 
293 // ----------------------------------------------------------------------------
294 
307 template <typename T>
308 class PCL_CLASS GenericFFT2D : public AbstractFFT2D<T>
309 {
310 public:
311 
316 
320  using scalar = typename base::scalar;
321 
325  using complex = typename base::complex;
326 
330  using vector = typename base::vector;
331 
336 
340  using matrix = typename base::matrix;
341 
346 
350  using transform = typename base::transform;
351 
367  GenericFFT2D( int rows, int cols )
368  : base( rows, cols )
369  {
370  }
371 
381  GenericFFT2D( int rows, int cols, StatusMonitor& status )
382  : base( rows, cols, status )
383  {
384  }
385 
389  ~GenericFFT2D() override
390  {
391  }
392 
405  {
406  if ( m_dft.IsEmpty() )
407  m_dft = transform( m_rows, m_cols );
408  this->Transform( m_rows, m_cols, *m_dft, x, PCL_FFT_FORWARD, m_monitor, m_parallel, m_maxProcessors );
409  return *this;
410  }
411 
427  {
428  if ( m_dft.IsEmpty() )
429  throw Error( "Invalid out-of-place inverse FFT: No FFT has been performed." );
430  this->Transform( m_rows, m_cols, y, *m_dft, PCL_FFT_BACKWARD, m_monitor, m_parallel, m_maxProcessors );
431  return *this;
432  }
433 
444  {
445  PCL_PRECONDITION( x.Rows() == m_rows )
446  PCL_PRECONDITION( x.Cols() == m_cols )
447  if ( x.Rows() != m_rows || x.Cols() != m_cols )
448  throw Error( "Invalid FFT input matrix dimensions." );
449  return operator <<( *x );
450  }
451 
465  {
466  PCL_PRECONDITION( y.Rows() == m_rows )
467  PCL_PRECONDITION( y.Cols() == m_cols )
468  if ( y.Rows() != m_rows || y.Cols() != m_cols )
469  throw Error( "Invalid FFT output matrix dimensions." );
470  return operator >>( *y );
471  }
472 
505  GenericFFT2D& operator ()( complex* y, const complex* x, int dir = PCL_FFT_FORWARD ) const
506  {
507  this->Transform( m_rows, m_cols, y, x, dir, m_monitor, m_parallel, m_maxProcessors );
508  return const_cast<GenericFFT2D&>( *this );
509  }
510 
518  static int OptimizedLength( int n )
519  {
520  return GenericFFT<T>::OptimizedLength( n );
521  }
522 };
523 
524 // ----------------------------------------------------------------------------
525 
540 template <typename T>
541 class PCL_CLASS GenericRealFFT2D : public AbstractFFT2D<T>
542 {
543 public:
544 
549 
553  using scalar = typename base::scalar;
554 
558  using complex = typename base::complex;
559 
563  using vector = typename base::vector;
564 
569 
573  using matrix = typename base::matrix;
574 
579 
583  using transform = typename base::transform;
584 
604  GenericRealFFT2D( int rows, int cols )
605  : base( rows, cols )
606  {
607  }
608 
618  GenericRealFFT2D( int rows, int cols, StatusMonitor& status )
619  : base( rows, cols, status )
620  {
621  }
622 
626  ~GenericRealFFT2D() override
627  {
628  }
629 
642  {
643  if ( m_dft.IsEmpty() )
644  m_dft = transform( m_rows, m_cols/2 + 1 );
645  this->Transform( m_rows, m_cols, *m_dft, x, m_monitor, m_parallel, m_maxProcessors );
646  return *this;
647  }
648 
664  {
665  if ( m_dft.IsEmpty() )
666  throw Error( "Invalid out-of-place inverse FFT: No FFT has been performed." );
667  this->Transform( m_rows, m_cols, y, *m_dft, m_monitor, m_parallel, m_maxProcessors );
668  return *this;
669  }
670 
681  {
682  PCL_PRECONDITION( x.Rows() == m_rows )
683  PCL_PRECONDITION( x.Cols() == m_cols )
684  if ( x.Rows() != m_rows || x.Cols() != m_cols )
685  throw Error( "Invalid FFT input matrix dimensions." );
686  return operator <<( *x );
687  }
688 
702  {
703  PCL_PRECONDITION( y.Rows() == m_rows )
704  PCL_PRECONDITION( y.Cols() == m_cols )
705  if ( y.Rows() != m_rows || y.Cols() != m_cols )
706  throw Error( "Invalid FFT output matrix dimensions." );
707  return operator >>( *y );
708  }
709 
735  GenericRealFFT2D& operator()( complex* y, const scalar* x ) const
736  {
737  this->Transform( m_rows, m_cols, y, x, m_monitor, m_parallel, m_maxProcessors );
738  return const_cast<GenericRealFFT2D&>( *this );
739  }
740 
766  GenericRealFFT2D& operator()( scalar* y, const complex* x ) const
767  {
768  this->Transform( m_rows, m_cols, y, x, m_monitor, m_parallel, m_maxProcessors );
769  return const_cast<GenericRealFFT2D&>( *this );
770  }
771 
779  static int OptimizedLength( int n )
780  {
781  return GenericFFT<T>::OptimizedLength( n );
782  }
783 };
784 
785 // ----------------------------------------------------------------------------
786 
787 #undef m_dft
788 #undef m_rows
789 #undef m_cols
790 #undef m_monitor
791 #undef m_parallel
792 #undef m_maxProcessors
793 
794 // ----------------------------------------------------------------------------
795 
800 #ifndef __PCL_NO_FFT2D_INSTANTIATE
801 
809 using FFFT2D = GenericFFT2D<float>;
810 
818 using DFFT2D = GenericFFT2D<double>;
819 
827 using FRealFFT2D = GenericRealFFT2D<float>;
828 
836 using DRealFFT2D = GenericRealFFT2D<double>;
837 
846 using FFT2D = FFFT2D;
847 
856 using RealFFT2D = FRealFFT2D;
857 
858 #endif // __PCL_NO_FFT2D_INSTANTIATE
859 
860 // ----------------------------------------------------------------------------
861 
862 } // pcl
863 
864 #endif // __PCL_FFT2D_h
865 
866 // ----------------------------------------------------------------------------
867 // EOF pcl/FFT2D.h - Released 2024-06-18T15:48:54Z
Abstract base class of all two-dimensional fast Fourier transform classes.
Definition: FFT2D.h:130
int Cols() const
Definition: FFT2D.h:216
int NumberOfElements() const
Definition: FFT2D.h:226
int Rows() const
Definition: FFT2D.h:206
transform & DFT()
Definition: FFT2D.h:258
AbstractFFT2D(int rows, int cols, StatusMonitor &monitor)
Definition: FFT2D.h:185
virtual void Release()
Definition: FFT2D.h:267
transform DFT() const
Definition: FFT2D.h:242
~AbstractFFT2D() override
Definition: FFT2D.h:196
AbstractFFT2D(int rows, int cols)
Definition: FFT2D.h:172
Generic complex number.
Definition: Complex.h:84
Fast Fourier transform of 64-bit floating point complex data.
Fast Fourier transform of 64-bit floating point real data.
A simple exception with an associated error message.
Definition: Exception.h:239
Fast Fourier transform of 32-bit floating point complex data.
Fast Fourier transform of 32-bit floating point complex data.
Fast Fourier transform of 32-bit floating point real data.
Generic two-dimensional fast Fourier transform of complex data.
Definition: FFT2D.h:309
~GenericFFT2D() override
Definition: FFT2D.h:389
GenericFFT2D(int rows, int cols, StatusMonitor &status)
Definition: FFT2D.h:381
GenericFFT2D(int rows, int cols)
Definition: FFT2D.h:367
static int OptimizedLength(int n)
Definition: FFT2D.h:518
static int OptimizedLength(int n)
Definition: FFT1D.h:520
Generic dynamic matrix of arbitrary dimensions.
Definition: Matrix.h:123
int Cols() const noexcept
Definition: Matrix.h:969
int Rows() const noexcept
Definition: Matrix.h:960
Generic two-dimensional fast Fourier transform of real data.
Definition: FFT2D.h:542
GenericRealFFT2D & operator()(complex *y, const scalar *x) const
Definition: FFT2D.h:735
~GenericRealFFT2D() override
Definition: FFT2D.h:626
GenericRealFFT2D(int rows, int cols)
Definition: FFT2D.h:604
static int OptimizedLength(int n)
Definition: FFT2D.h:779
GenericRealFFT2D(int rows, int cols, StatusMonitor &status)
Definition: FFT2D.h:618
GenericRealFFT2D & operator()(scalar *y, const complex *x) const
Definition: FFT2D.h:766
A process using multiple concurrent execution threads.
Fast Fourier transform of 32-bit floating point real data.
An asynchronous status monitoring system.
A complex number whose components are 64-bit floating point real numbers.
A complex number whose components are 32-bit floating point real numbers.
Array< T, A > & operator<<(Array< T, A > &x, const V &v)
Definition: Array.h:2295
Console & operator>>(Console &o, char &c)
Definition: Console.h:803
#define PCL_FFT_BACKWARD
Indicates an inverse Fourier transform.
Definition: FFT1D.h:150
#define PCL_FFT_FORWARD
Indicates a Fourier transform.
Definition: FFT1D.h:143
PCL root namespace.
Definition: AbstractImage.h:77