Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TFFTComplex.cxx
Go to the documentation of this file.
1// @(#)root/fft:$Id$
2// Author: Anna Kreshuk 07/4/2006
3
4/*************************************************************************
5 * Copyright (C) 1995-2006, Rene Brun and Fons Rademakers. *
6 * All rights reserved. *
7 * *
8 * For the licensing terms see $ROOTSYS/LICENSE. *
9 * For the list of contributors see $ROOTSYS/README/CREDITS. *
10 *************************************************************************/
11
12////////////////////////////////////////////////////////////////////////////////
13/// \class TFFTComplex
14///
15/// One of the interface classes to the FFTW package, can be used directly
16/// or via the TVirtualFFT class. Only the basic interface of FFTW is implemented.
17/// Computes complex input/output discrete Fourier transforms (DFT)
18/// in one or more dimensions. For the detailed information on the computed
19/// transforms please refer to the FFTW manual, chapter "What FFTW really computes".
20///
21/// How to use it:
22///
23/// 1. Create an instance of TFFTComplex - this will allocate input and output
24/// arrays (unless an in-place transform is specified)
25/// 2. Run the Init() function with the desired flags and settings
26/// 3. Set the data (via SetPoints(), SetPoint() or SetPointComplex() functions)
27/// 4. Run the Transform() function
28/// 5. Get the output (via GetPoints(), GetPoint() or GetPointComplex() functions)
29/// 6. Repeat steps 3)-5) as needed
30///
31/// For a transform of the same size, but with different flags or sign, rerun the Init()
32/// function and continue with steps 3)-5)
33///
34/// NOTE:
35/// 1. running Init() function will overwrite the input array! Don't set any data
36/// before running the Init() function
37/// 2. FFTW computes unnormalized transform, so doing a transform followed by
38/// its inverse will lead to the original array scaled by the transform size
39///
40////////////////////////////////////////////////////////////////////////////////
41
42#include "TFFTComplex.h"
43#include "fftw3.h"
44#include "TComplex.h"
45
46
47
48////////////////////////////////////////////////////////////////////////////////
49///default
50
52{
53 fIn = nullptr;
54 fOut = nullptr;
55 fPlan = nullptr;
56 fN = nullptr;
57 fNdim = 0;
58 fTotalSize = 0;
59 fSign = 1;
60}
61
62////////////////////////////////////////////////////////////////////////////////
63///For 1d transforms
64///Allocates memory for the input array, and, if inPlace = kFALSE, for the output array
65
67{
68 fIn = fftw_malloc(sizeof(fftw_complex) *n);
69 if (!inPlace)
70 fOut = fftw_malloc(sizeof(fftw_complex) * n);
71 else
72 fOut = nullptr;
73 fN = new Int_t[1];
74 fN[0] = n;
75 fTotalSize = n;
76 fNdim = 1;
77 fSign = 1;
78 fPlan = nullptr;
79}
80
81////////////////////////////////////////////////////////////////////////////////
82///For multidim. transforms
83///Allocates memory for the input array, and, if inPlace = kFALSE, for the output array
84
86{
87 fNdim = ndim;
88 fTotalSize = 1;
89 fN = new Int_t[fNdim];
90 for (Int_t i=0; i<fNdim; i++){
91 fN[i] = n[i];
92 fTotalSize*=n[i];
93 }
95 if (!inPlace)
97 else
98 fOut = nullptr;
99 fSign = 1;
100 fPlan = nullptr;
101}
102
103////////////////////////////////////////////////////////////////////////////////
104///Destroys the data arrays and the plan. However, some plan information stays around
105///until the root session is over, and is reused if other plans of the same size are
106///created
107
109{
111 fPlan = nullptr;
113 if (fOut)
115 if (fN)
116 delete [] fN;
117}
118
119////////////////////////////////////////////////////////////////////////////////
120///Creates the fftw-plan
121///
122///NOTE: input and output arrays are overwritten during initialisation,
123/// so don't set any points, before running this function!!!!!
124///
125///2nd parameter: +1
126///
127///Argument kind is dummy and doesn't need to be specified
128///
129///Possible flag_options:
130/// - "ES" (from "estimate") - no time in preparing the transform, but probably sub-optimal
131/// performance
132/// - "M" (from "measure") - some time spend in finding the optimal way to do the transform
133/// - "P" (from "patient") - more time spend in finding the optimal way to do the transform
134/// - "EX" (from "exhaustive") - the most optimal way is found
135///This option should be chosen depending on how many transforms of the same size and
136///type are going to be done. Planning is only done once, for the first transform of this
137///size and type.
138
139void TFFTComplex::Init( Option_t *flags, Int_t sign,const Int_t* /*kind*/)
140{
141 fSign = sign;
142 fFlags = flags;
143
144 if (fPlan)
146 fPlan = nullptr;
147
148 if (fOut)
150 else
152}
153
154////////////////////////////////////////////////////////////////////////////////
155///Computes the transform, specified in Init() function
156
158{
159 if (fPlan)
161 else {
162 Error("Transform", "transform not initialised");
163 return;
164 }
165}
166
167////////////////////////////////////////////////////////////////////////////////
168///Copies the output(or input) into the argument array
169
171{
172 if (!fromInput){
173 for (Int_t i=0; i<2*fTotalSize; i+=2){
174 data[i] = ((fftw_complex*)fOut)[i/2][0];
175 data[i+1] = ((fftw_complex*)fOut)[i/2][1];
176 }
177 } else {
178 for (Int_t i=0; i<2*fTotalSize; i+=2){
179 data[i] = ((fftw_complex*)fIn)[i/2][0];
180 data[i+1] = ((fftw_complex*)fIn)[i/2][1];
181 }
182 }
183}
184
185////////////////////////////////////////////////////////////////////////////////
186///returns real and imaginary parts of the point `#ipoint`
187
189{
190 if (fOut && !fromInput){
191 re = ((fftw_complex*)fOut)[ipoint][0];
192 im = ((fftw_complex*)fOut)[ipoint][1];
193 } else {
194 re = ((fftw_complex*)fIn)[ipoint][0];
195 im = ((fftw_complex*)fIn)[ipoint][1];
196 }
197}
198
199////////////////////////////////////////////////////////////////////////////////
200///For multidimensional transforms. Returns real and imaginary parts of the point `#ipoint`
201
203{
204 Int_t ireal = ipoint[0];
205 for (Int_t i=0; i<fNdim-1; i++)
206 ireal=fN[i+1]*ireal + ipoint[i+1];
207
208 if (fOut && !fromInput){
209 re = ((fftw_complex*)fOut)[ireal][0];
210 im = ((fftw_complex*)fOut)[ireal][1];
211 } else {
212 re = ((fftw_complex*)fIn)[ireal][0];
213 im = ((fftw_complex*)fIn)[ireal][1];
214 }
215}
216
217////////////////////////////////////////////////////////////////////////////////
218///Copies real and imaginary parts of the output (input) into the argument arrays
219
221{
222 if (fOut && !fromInput){
223 for (Int_t i=0; i<fTotalSize; i++){
224 re[i] = ((fftw_complex*)fOut)[i][0];
225 im[i] = ((fftw_complex*)fOut)[i][1];
226 }
227 } else {
228 for (Int_t i=0; i<fTotalSize; i++){
229 re[i] = ((fftw_complex*)fIn)[i][0];
230 im[i] = ((fftw_complex*)fIn)[i][1];
231 }
232 }
233}
234
235////////////////////////////////////////////////////////////////////////////////
236///Copies the output(input) into the argument array
237
239{
240 if (fOut && !fromInput){
241 for (Int_t i=0; i<fTotalSize; i+=2){
242 data[i] = ((fftw_complex*)fOut)[i/2][0];
243 data[i+1] = ((fftw_complex*)fOut)[i/2][1];
244 }
245 } else {
246 for (Int_t i=0; i<fTotalSize; i+=2){
247 data[i] = ((fftw_complex*)fIn)[i/2][0];
248 data[i+1] = ((fftw_complex*)fIn)[i/2][1];
249 }
250 }
251}
252
253////////////////////////////////////////////////////////////////////////////////
254///sets real and imaginary parts of point # ipoint
255
261
262////////////////////////////////////////////////////////////////////////////////
263///For multidim. transforms. Sets real and imaginary parts of point # ipoint
264
266{
267 Int_t ireal = ipoint[0];
268 for (Int_t i=0; i<fNdim-1; i++)
269 ireal=fN[i+1]*ireal + ipoint[i+1];
270
271 ((fftw_complex*)fIn)[ireal][0]=re;
272 ((fftw_complex*)fIn)[ireal][1]=im;
273}
274
275////////////////////////////////////////////////////////////////////////////////
276
278{
279 ((fftw_complex*)fIn)[ipoint][0] = c.Re();
280 ((fftw_complex*)fIn)[ipoint][1] = c.Im();
281}
282
283////////////////////////////////////////////////////////////////////////////////
284///set all points. the values are copied. points should be ordered as follows:
285///[re_0, im_0, re_1, im_1, ..., re_n, im_n)
286
288{
289 for (Int_t i=0; i<2*fTotalSize-1; i+=2){
290 ((fftw_complex*)fIn)[i/2][0]=data[i];
291 ((fftw_complex*)fIn)[i/2][1]=data[i+1];
292 }
293}
294
295////////////////////////////////////////////////////////////////////////////////
296///set all points. the values are copied
297
299{
300 if (!fIn){
301 Error("SetPointsComplex", "Size is not set yet");
302 return;
303 }
304 for (Int_t i=0; i<fTotalSize; i++){
305 ((fftw_complex*)fIn)[i][0]=re_data[i];
306 ((fftw_complex*)fIn)[i][1]=im_data[i];
307 }
308}
309
310////////////////////////////////////////////////////////////////////////////////
311///allowed options:
312/// - "ES" - FFTW_ESTIMATE
313/// - "M" - FFTW_MEASURE
314/// - "P" - FFTW_PATIENT
315/// - "EX" - FFTW_EXHAUSTIVE
316
318{
319 TString opt = flag;
320 opt.ToUpper();
321 if (opt.Contains("ES"))
322 return FFTW_ESTIMATE;
323 if (opt.Contains("M"))
324 return FFTW_MEASURE;
325 if (opt.Contains("P"))
326 return FFTW_PATIENT;
327 if (opt.Contains("EX"))
328 return FFTW_EXHAUSTIVE;
329 return FFTW_ESTIMATE;
330}
#define c(i)
Definition RSha256.hxx:101
const char Option_t
Option string (const char)
Definition RtypesCore.h:80
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
Int_t fTotalSize
Definition TFFTComplex.h:26
void * fPlan
Definition TFFTComplex.h:24
UInt_t MapFlag(Option_t *flag)
allowed options:
void SetPoints(const Double_t *data) override
set all points.
void SetPointComplex(Int_t ipoint, TComplex &c) override
void GetPointsComplex(Double_t *re, Double_t *im, Bool_t fromInput=kFALSE) const override
Copies real and imaginary parts of the output (input) into the argument arrays.
void GetPoints(Double_t *data, Bool_t fromInput=kFALSE) const override
Copies the output(or input) into the argument array.
void SetPointsComplex(const Double_t *re, const Double_t *im) override
set all points. the values are copied
void Transform() override
Computes the transform, specified in Init() function.
TString fFlags
Definition TFFTComplex.h:29
TFFTComplex()
default
void Init(Option_t *flags, Int_t sign, const Int_t *) override
Creates the fftw-plan.
void * fOut
Definition TFFTComplex.h:23
Int_t * fN
Definition TFFTComplex.h:27
~TFFTComplex() override
Destroys the data arrays and the plan.
void GetPointComplex(Int_t ipoint, Double_t &re, Double_t &im, Bool_t fromInput=kFALSE) const override
returns real and imaginary parts of the point #ipoint
void SetPoint(Int_t ipoint, Double_t re, Double_t im=0) override
sets real and imaginary parts of point # ipoint
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition TObject.cxx:1071
Basic string class.
Definition TString.h:138
void ToUpper()
Change string to upper case.
Definition TString.cxx:1202
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition TString.h:640
const Int_t n
Definition legend1.C:16