MFFM FFTw Wrapper
real2DFFT.H
1/* Copyright 2001,2002 Matt Flax <flatmax@ieee.org>
2 This file is part of the MFFM FFTw Wrapper library.
3
4 MFFM MFFM FFTw Wrapper library is free software; you can
5 redistribute it and/or modify
6 it under the terms of the GNU General Public License as published by
7 the Free Software Foundation; either version 2 of the License, or
8 (at your option) any later version.
9
10 MFFM FFTw Wrapper library is distributed in the hope that it will be useful,
11 but WITHOUT ANY WARRANTY; without even the implied warranty of
12 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 GNU General Public License for more details.
14
15 You have received a copy of the GNU General Public License
16 along with the MFFM FFTw Wrapper library
17*/
18#ifndef REAL2DFFT_H_
19#define REAL2DFFT_H_
20
21#include <string.h>
22
23#include <fftw3.h>
24#ifndef fftw_real
25#define fftw_real double
26#endif
27#define c_re(c) ((c)[0])
28#define c_im(c) ((c)[1])
29
30#include <iomanip>
31using namespace std;
32
33#define PLANTYPE FFTW_ESTIMATE
34
35/// class real2DFFTData controls and manipulates real 2D fft data
37 /// x=row y=column
38 int x, y;
39 /// The total memory used by this class
40 fftw_real *mem;
41 /// Free the memory
42 void memDeInit(void);
43public:
44 /// The input data and power spectrum
45 fftw_real *in, *power;
46 /// The output data
47 fftw_complex *out;
48 /// Arrays which sum across rows (x) and columns (y)
49 fftw_real *xSum, *ySum;
50 /// A sum across the input time signal
51 fftw_real *timeXSum;
52 /// Power spectral sums across rows (x) and columns (y)
53 fftw_real *realXSum, *imagXSum;
54
55 /// The total power in the power spectrum, the maximum and minimum powers too
56 double totalPower, maxPower, minPower;
57 /// The minimum/maximum row (x) and column (y) sums
58 double xSumMin, xSumMax, ySumMin, ySumMax;
59 /// Row (x) and Column (y) max sum indexes
60 int maxXSumIndex, maxYSumIndex;
61
62 /// Constructor with all memory to be allocated internally
63 real2DFFTData(int r, int c);
64 /// Deconstructor
66
67 /// The row count
68 int getXSize(){return x;}
69 /// The column count
70 int getYSize(){return y;}
71 /// The half row count
72 int getXHalfSize(){ if (!(x%2)) return x/2; else return x/2+1;}
73 /// The half column count
74 int getYHalfSize(){ if (!(y%2)) return y/2; else return y/2+1;}
75
76 /// Scales the output down by the number of elements
77 void reScale(void);
78 /// This function computes the power spectrum and updates the totalPower, maxPower and minPower
79 void compPowerSpec(); // Find the power spectrum
80 /// Finds 10*log10(power spectrum) and updates the totalPower, maxPower and minPower
81 void compLogPowerSpec(); // Find the log power spectrum
82
83 /// Updates timeXSum
85 /// Updates realXSum and imagXSum
87 /// Finds the power Spectrum averages and
88 /// updates the xSumMin, xSumMax, ySumMin, ySumMax, xSum, ySum, maxXSumIndex, maxYSumIndex
90 /// Finds the y-sum between columns start and stop
91 void findYSum(int start, int stop);
92 /// Finds the y-max for the ySum array, updates ySumMin, ySumMax, maxYSumIndex
93 void findYMax(void);
94
95 /// Zeros the in array
96 void clearInput(void){memset(in, 0, x*2*(y/2+1)*sizeof(fftw_real));}
97 /// Zeros the out awway
98 void clearOutput(void){memset(out, 0, x*(y/2+1)*sizeof(fftw_complex));}
99};
100
101///class real2DFFT controls fftw plans and executes fwd/inv transforms
103 /// The forward and inverse plans
104 fftw_plan fwdPlan, invPlan;
105protected:
106 /// The pointer to the relevant data
108public:
109 /// fft init ... data pointed to by 'd'
111 /// fft deconstructor
113
114 /// Forward transform the data (in to out)
115 void fwdTransform(); // Forward 2D fft
116 /// Inverse transform the data (out to in)
117 void invTransform(); // Inverse 2D fft
118};
119/** \example real2DFFTExample.cc
120 * This is an example of how to use the class.
121 */
122#endif // REAL2DFFT_H_
class real2DFFTData controls and manipulates real 2D fft data
Definition real2DFFT.H:36
int getXSize()
The row count.
Definition real2DFFT.H:68
fftw_real * xSum
Arrays which sum across rows (x) and columns (y)
Definition real2DFFT.H:49
double totalPower
The total power in the power spectrum, the maximum and minimum powers too.
Definition real2DFFT.H:56
fftw_real * realXSum
Power spectral sums across rows (x) and columns (y)
Definition real2DFFT.H:53
int getYSize()
The column count.
Definition real2DFFT.H:70
void reScale(void)
Scales the output down by the number of elements.
void clearOutput(void)
Zeros the out awway.
Definition real2DFFT.H:98
void compPowerSpec()
This function computes the power spectrum and updates the totalPower, maxPower and minPower.
void findYSum(int start, int stop)
Finds the y-sum between columns start and stop.
void timeSpecAverage()
Updates timeXSum.
~real2DFFTData()
Deconstructor.
double xSumMin
The minimum/maximum row (x) and column (y) sums.
Definition real2DFFT.H:58
void powerSpecAverage()
int maxXSumIndex
Row (x) and Column (y) max sum indexes.
Definition real2DFFT.H:60
void clearInput(void)
Zeros the in array.
Definition real2DFFT.H:96
int getYHalfSize()
The half column count.
Definition real2DFFT.H:74
int getXHalfSize()
The half row count.
Definition real2DFFT.H:72
void compLogPowerSpec()
Finds 10*log10(power spectrum) and updates the totalPower, maxPower and minPower.
real2DFFTData(int r, int c)
Constructor with all memory to be allocated internally.
fftw_complex * out
The output data.
Definition real2DFFT.H:47
void findYMax(void)
Finds the y-max for the ySum array, updates ySumMin, ySumMax, maxYSumIndex.
fftw_real * in
The input data and power spectrum.
Definition real2DFFT.H:45
void complexSpecAverage()
Updates realXSum and imagXSum.
fftw_real * timeXSum
A sum across the input time signal.
Definition real2DFFT.H:51
class real2DFFT controls fftw plans and executes fwd/inv transforms
Definition real2DFFT.H:102
void fwdTransform()
Forward transform the data (in to out)
real2DFFT(real2DFFTData *d)
fft init ... data pointed to by 'd'
void invTransform()
Inverse transform the data (out to in)
real2DFFTData * data
The pointer to the relevant data.
Definition real2DFFT.H:107
~real2DFFT()
fft deconstructor