// @(#)root/fft:$Id$
// Author: Anna Kreshuk   07/4/2006

/*************************************************************************
 * Copyright (C) 1995-2006, Rene Brun and Fons Rademakers.               *
 * All rights reserved.                                                  *
 *                                                                       *
 * For the licensing terms see $ROOTSYS/LICENSE.                         *
 * For the list of contributors see $ROOTSYS/README/CREDITS.             *
 *************************************************************************/

//////////////////////////////////////////////////////////////////////////
//
// TFFTComplexReal
//
// One of the interface classes to the FFTW package, can be used directly
// or via the TVirtualFFT class. Only the basic interface of FFTW is implemented.
//
// Computes the inverse of the real-to-complex transforms (class TFFTRealComplex)
// taking complex input (storing the non-redundant half of a logically Hermitian array)
// to real output (see FFTW manual for more details)
//
// How to use it:
// 1) Create an instance of TFFTComplexReal - this will allocate input and output
//    arrays (unless an in-place transform is specified)
// 2) Run the Init() function with the desired flags and settings
// 3) Set the data (via SetPoints(), SetPoint() or SetPointComplex() functions)
// 4) Run the Transform() function
// 5) Get the output (via GetPoints(), GetPoint() or GetPointReal() functions)
// 6) Repeat steps 3)-5) as needed
//
// For a transform of the same size, but with different flags, rerun the Init()
// function and continue with steps 3)-5)
// NOTE: 1) running Init() function will overwrite the input array! Don't set any data
//          before running the Init() function
//       2) FFTW computes unnormalized transform, so doing a transform followed by
//          its inverse will lead to the original array scaled by the transform size
//
//       3) In Complex to Real transform the input array is destroyed. It cannot then
//          be retrieved when using the Get's methods.
//
//////////////////////////////////////////////////////////////////////////

#include "TFFTComplexReal.h"
#include "fftw3.h"
#include "TComplex.h"


ClassImp(TFFTComplexReal)

//_____________________________________________________________________________
TFFTComplexReal::TFFTComplexReal()
{
//default

   fIn   = 0;
   fOut  = 0;
   fPlan = 0;
   fN    = 0;
   fTotalSize = 0;
   fFlags = 0;
   fNdim = 0;
}

//_____________________________________________________________________________
TFFTComplexReal::TFFTComplexReal(Int_t n, Bool_t inPlace)
{
//For 1d transforms
//Allocates memory for the input array, and, if inPlace = kFALSE, for the output array

   if (!inPlace){
      fIn = fftw_malloc(sizeof(fftw_complex)*(n/2+1));
      fOut = fftw_malloc(sizeof(Double_t)*n);

   } else {
      fIn = fftw_malloc(sizeof(fftw_complex)*(n/2+1));
      fOut = 0;
   }
   fNdim = 1;
   fN = new Int_t[1];
   fN[0] = n;
   fPlan = 0;
   fTotalSize = n;
   fFlags = 0;
}

//_____________________________________________________________________________
TFFTComplexReal::TFFTComplexReal(Int_t ndim, Int_t *n, Bool_t inPlace)
{
//For ndim-dimensional transforms
//Second argurment contains sizes of the transform in each dimension

   fNdim = ndim;
   fTotalSize = 1;
   fN = new Int_t[fNdim];
   for (Int_t i=0; i<fNdim; i++){
      fN[i] = n[i];
      fTotalSize*=n[i];
   }
   Int_t sizein = Int_t(Double_t(fTotalSize)*(n[ndim-1]/2+1)/n[ndim-1]);
   if (!inPlace){
      fIn = fftw_malloc(sizeof(fftw_complex)*sizein);
      fOut = fftw_malloc(sizeof(Double_t)*fTotalSize);
   } else {
      fIn = fftw_malloc(sizeof(fftw_complex)*sizein);
      fOut = 0;
   }
   fPlan = 0;
   fFlags = 0;
}


//_____________________________________________________________________________
TFFTComplexReal::~TFFTComplexReal()
{
//Destroys the data arrays and the plan. However, some plan information stays around
//until the root session is over, and is reused if other plans of the same size are
//created

   fftw_destroy_plan((fftw_plan)fPlan);
   fPlan = 0;
   fftw_free((fftw_complex*)fIn);
   if (fOut)
      fftw_free(fOut);
   fIn = 0;
   fOut = 0;
   if (fN)
      delete [] fN;
   fN = 0;
}

//_____________________________________________________________________________
void TFFTComplexReal::Init( Option_t *flags, Int_t /*sign*/,const Int_t* /*kind*/)
{
//Creates the fftw-plan
//
//NOTE:  input and output arrays are overwritten during initialisation,
//       so don't set any points, before running this function!!!!!
//
//Arguments sign and kind are dummy and not need to be specified
//Possible flag_options:
//"ES" (from "estimate") - no time in preparing the transform, but probably sub-optimal
//   performanc
//"M" (from "measure") - some time spend in finding the optimal way to do the transform
//"P" (from "patient") - more time spend in finding the optimal way to do the transform
//"EX" (from "exhaustive") - the most optimal way is found
//This option should be chosen depending on how many transforms of the same size and
//type are going to be done. Planning is only done once, for the first transform of this
//size and type.

   fFlags = flags;

   if (fPlan)
      fftw_destroy_plan((fftw_plan)fPlan);
   fPlan = 0;

   if (fOut)
      fPlan = (void*)fftw_plan_dft_c2r(fNdim, fN,(fftw_complex*)fIn,(Double_t*)fOut, MapFlag(flags));
   else
      fPlan = (void*)fftw_plan_dft_c2r(fNdim, fN, (fftw_complex*)fIn, (Double_t*)fIn, MapFlag(flags));
}

//_____________________________________________________________________________
void TFFTComplexReal::Transform()
{
//Computes the transform, specified in Init() function

   if (fPlan)
      fftw_execute((fftw_plan)fPlan);
   else {
      Error("Transform", "transform was not initialized");
      return;
   }
}

//_____________________________________________________________________________
void TFFTComplexReal::GetPoints(Double_t *data, Bool_t fromInput) const
{
//Fills the argument array with the computed transform
// Works only for output (input array is destroyed in a C2R transform)

   if (fromInput){
      Error("GetPoints", "Input array has been destroyed");
      return;
   }
   const Double_t * array =  (const Double_t*) ( (fOut) ?  fOut :  fIn );
   std::copy(array,array+fTotalSize, data);
}

//_____________________________________________________________________________
Double_t TFFTComplexReal::GetPointReal(Int_t ipoint, Bool_t fromInput) const
{
//Returns the point #ipoint
// Works only for output (input array is destroyed in a C2R transform)

   if (fromInput){
      Error("GetPointReal", "Input array has been destroyed");
      return 0;
   }
   const Double_t * array =  (const Double_t*) ( (fOut) ?  fOut :  fIn );
   return array[ipoint];
}

//_____________________________________________________________________________
Double_t TFFTComplexReal::GetPointReal(const Int_t *ipoint, Bool_t fromInput) const
{
//For multidimensional transforms. Returns the point #ipoint
// Works only for output (input array is destroyed in a C2R transform)

   if (fromInput){
      Error("GetPointReal", "Input array has been destroyed");
      return 0;
   }
   Int_t ireal = ipoint[0];
   for (Int_t i=0; i<fNdim-1; i++)
      ireal=fN[i+1]*ireal + ipoint[i+1];

   const Double_t * array =  (const Double_t*) ( (fOut) ?  fOut :  fIn );
   return array[ireal];
}


//_____________________________________________________________________________
void TFFTComplexReal::GetPointComplex(Int_t ipoint, Double_t &re, Double_t &im, Bool_t fromInput) const
{
// Works only for output (input array is destroyed in a C2R transform)

   if (fromInput){
      Error("GetPointComplex", "Input array has been destroyed");
      return;
   }
   const Double_t * array =  (const Double_t*) ( (fOut) ?  fOut :  fIn );
   re = array[ipoint];
   im = 0;
}

//_____________________________________________________________________________
void TFFTComplexReal::GetPointComplex(const Int_t *ipoint, Double_t &re, Double_t &im, Bool_t fromInput) const
{
//For multidimensional transforms. Returns the point #ipoint
// Works only for output (input array is destroyed in a C2R transform)

   if (fromInput){
      Error("GetPointComplex", "Input array has been destroyed");
      return;
   }
   const Double_t * array =  (const Double_t*) ( (fOut) ?  fOut :  fIn );

   Int_t ireal = ipoint[0];
   for (Int_t i=0; i<fNdim-1; i++)
      ireal=fN[i+1]*ireal + ipoint[i+1];

   re = array[ireal];
   im = 0;
}
//_____________________________________________________________________________
Double_t* TFFTComplexReal::GetPointsReal(Bool_t fromInput) const
{
//Returns the array of computed transform
// Works only for output (input array is destroyed in a C2R transform)

   // we have 2 different cases
   // fromInput = false; fOut = !NULL (transformed is not in place) : return fOut
   // fromInput = false; fOut = NULL (transformed is in place) : return fIn

   if (fromInput) {
      Error("GetPointsReal","Input array was destroyed");
      return 0;
   }
   return  (Double_t*) ( (fOut) ?  fOut :  fIn );
}

//_____________________________________________________________________________
void TFFTComplexReal::GetPointsComplex(Double_t *re, Double_t *im, Bool_t fromInput) const
{
//Fills the argument array with the computed transform
// Works only for output (input array is destroyed in a C2R transform)

   if (fromInput){
      Error("GetPointsComplex", "Input array has been destroyed");
      return;
   }
   const Double_t * array =  (const Double_t*) ( (fOut) ?  fOut :  fIn );
   for (Int_t i=0; i<fTotalSize; i++){
      re[i] = array[i];
      im[i] = 0;
   }
}

//_____________________________________________________________________________
void TFFTComplexReal::GetPointsComplex(Double_t *data, Bool_t fromInput) const
{
//Fills the argument array with the computed transform.
// Works only for output (input array is destroyed in a C2R transform)

   if (fromInput){
      Error("GetPointsComplex", "Input array has been destroyed");
      return;
   }
   const Double_t * array =  (const Double_t*) ( (fOut) ?  fOut :  fIn );
   for (Int_t i=0; i<fTotalSize; i+=2){
      data[i] = array[i/2];
      data[i+1]=0;
   }
}

//_____________________________________________________________________________
void TFFTComplexReal::SetPoint(Int_t ipoint, Double_t re, Double_t im)
{
//since the input must be complex-Hermitian, if the ipoint > n/2, the according
//point before n/2 is set to (re, -im)

   if (ipoint <= fN[0]/2){
      ((fftw_complex*)fIn)[ipoint][0] = re;
      ((fftw_complex*)fIn)[ipoint][1] = im;
   } else {
      ((fftw_complex*)fOut)[2*(fN[0]/2)-ipoint][0] = re;
      ((fftw_complex*)fOut)[2*(fN[0]/2)-ipoint][1] = -im;
   }
}

//_____________________________________________________________________________
void TFFTComplexReal::SetPoint(const Int_t *ipoint, Double_t re, Double_t im)
{
//Set the point #ipoint. Since the input is Hermitian, only the first (roughly)half of
//the points have to be set.

   Int_t ireal = ipoint[0];
   for (Int_t i=0; i<fNdim-2; i++){
      ireal=fN[i+1]*ireal + ipoint[i+1];
   }
   ireal = (fN[fNdim-1]/2+1)*ireal+ipoint[fNdim-1];
   Int_t realN = Int_t(Double_t(fTotalSize)*(fN[fNdim-1]/2+1)/fN[fNdim-1]);

   if (ireal > realN){
      Error("SetPoint", "Illegal index value");
      return;
   }
   ((fftw_complex*)fIn)[ireal][0] = re;
   ((fftw_complex*)fIn)[ireal][1] = im;
}

//_____________________________________________________________________________
void TFFTComplexReal::SetPointComplex(Int_t ipoint, TComplex &c)
{
//since the input must be complex-Hermitian, if the ipoint > n/2, the according
//point before n/2 is set to (re, -im)

   if (ipoint <= fN[0]/2){
      ((fftw_complex*)fIn)[ipoint][0] = c.Re();
      ((fftw_complex*)fIn)[ipoint][1] = c.Im();
   } else {
      ((fftw_complex*)fIn)[2*(fN[0]/2)-ipoint][0] = c.Re();
      ((fftw_complex*)fIn)[2*(fN[0]/2)-ipoint][1] = -c.Im();
   }
}

//_____________________________________________________________________________
void TFFTComplexReal::SetPoints(const Double_t *data)
{
//set all points. the values are copied. points should be ordered as follows:
//[re_0, im_0, re_1, im_1, ..., re_n, im_n)

   Int_t sizein = Int_t(Double_t(fTotalSize)*(fN[fNdim-1]/2+1)/fN[fNdim-1]);

   for (Int_t i=0; i<2*(sizein); i+=2){
      ((fftw_complex*)fIn)[i/2][0]=data[i];
      ((fftw_complex*)fIn)[i/2][1]=data[i+1];
   }
}

//_____________________________________________________________________________
void TFFTComplexReal::SetPointsComplex(const Double_t *re, const Double_t *im)
{
//Set all points. The values are copied.

   Int_t sizein = Int_t(Double_t(fTotalSize)*(fN[fNdim-1]/2+1)/fN[fNdim-1]);
   for (Int_t i=0; i<sizein; i++){
      ((fftw_complex*)fIn)[i][0]=re[i];
      ((fftw_complex*)fIn)[i][1]=im[i];
   }
}

//_____________________________________________________________________________
UInt_t TFFTComplexReal::MapFlag(Option_t *flag)
{
//allowed options:
//"ES" - FFTW_ESTIMATE
//"M" - FFTW_MEASURE
//"P" - FFTW_PATIENT
//"EX" - FFTW_EXHAUSTIVE

   TString opt = flag;
   opt.ToUpper();
   if (opt.Contains("ES"))
      return FFTW_ESTIMATE;
   if (opt.Contains("M"))
      return FFTW_MEASURE;
   if (opt.Contains("P"))
      return FFTW_PATIENT;
   if (opt.Contains("EX"))
      return FFTW_EXHAUSTIVE;
   return FFTW_ESTIMATE;
}
 TFFTComplexReal.cxx:1
 TFFTComplexReal.cxx:2
 TFFTComplexReal.cxx:3
 TFFTComplexReal.cxx:4
 TFFTComplexReal.cxx:5
 TFFTComplexReal.cxx:6
 TFFTComplexReal.cxx:7
 TFFTComplexReal.cxx:8
 TFFTComplexReal.cxx:9
 TFFTComplexReal.cxx:10
 TFFTComplexReal.cxx:11
 TFFTComplexReal.cxx:12
 TFFTComplexReal.cxx:13
 TFFTComplexReal.cxx:14
 TFFTComplexReal.cxx:15
 TFFTComplexReal.cxx:16
 TFFTComplexReal.cxx:17
 TFFTComplexReal.cxx:18
 TFFTComplexReal.cxx:19
 TFFTComplexReal.cxx:20
 TFFTComplexReal.cxx:21
 TFFTComplexReal.cxx:22
 TFFTComplexReal.cxx:23
 TFFTComplexReal.cxx:24
 TFFTComplexReal.cxx:25
 TFFTComplexReal.cxx:26
 TFFTComplexReal.cxx:27
 TFFTComplexReal.cxx:28
 TFFTComplexReal.cxx:29
 TFFTComplexReal.cxx:30
 TFFTComplexReal.cxx:31
 TFFTComplexReal.cxx:32
 TFFTComplexReal.cxx:33
 TFFTComplexReal.cxx:34
 TFFTComplexReal.cxx:35
 TFFTComplexReal.cxx:36
 TFFTComplexReal.cxx:37
 TFFTComplexReal.cxx:38
 TFFTComplexReal.cxx:39
 TFFTComplexReal.cxx:40
 TFFTComplexReal.cxx:41
 TFFTComplexReal.cxx:42
 TFFTComplexReal.cxx:43
 TFFTComplexReal.cxx:44
 TFFTComplexReal.cxx:45
 TFFTComplexReal.cxx:46
 TFFTComplexReal.cxx:47
 TFFTComplexReal.cxx:48
 TFFTComplexReal.cxx:49
 TFFTComplexReal.cxx:50
 TFFTComplexReal.cxx:51
 TFFTComplexReal.cxx:52
 TFFTComplexReal.cxx:53
 TFFTComplexReal.cxx:54
 TFFTComplexReal.cxx:55
 TFFTComplexReal.cxx:56
 TFFTComplexReal.cxx:57
 TFFTComplexReal.cxx:58
 TFFTComplexReal.cxx:59
 TFFTComplexReal.cxx:60
 TFFTComplexReal.cxx:61
 TFFTComplexReal.cxx:62
 TFFTComplexReal.cxx:63
 TFFTComplexReal.cxx:64
 TFFTComplexReal.cxx:65
 TFFTComplexReal.cxx:66
 TFFTComplexReal.cxx:67
 TFFTComplexReal.cxx:68
 TFFTComplexReal.cxx:69
 TFFTComplexReal.cxx:70
 TFFTComplexReal.cxx:71
 TFFTComplexReal.cxx:72
 TFFTComplexReal.cxx:73
 TFFTComplexReal.cxx:74
 TFFTComplexReal.cxx:75
 TFFTComplexReal.cxx:76
 TFFTComplexReal.cxx:77
 TFFTComplexReal.cxx:78
 TFFTComplexReal.cxx:79
 TFFTComplexReal.cxx:80
 TFFTComplexReal.cxx:81
 TFFTComplexReal.cxx:82
 TFFTComplexReal.cxx:83
 TFFTComplexReal.cxx:84
 TFFTComplexReal.cxx:85
 TFFTComplexReal.cxx:86
 TFFTComplexReal.cxx:87
 TFFTComplexReal.cxx:88
 TFFTComplexReal.cxx:89
 TFFTComplexReal.cxx:90
 TFFTComplexReal.cxx:91
 TFFTComplexReal.cxx:92
 TFFTComplexReal.cxx:93
 TFFTComplexReal.cxx:94
 TFFTComplexReal.cxx:95
 TFFTComplexReal.cxx:96
 TFFTComplexReal.cxx:97
 TFFTComplexReal.cxx:98
 TFFTComplexReal.cxx:99
 TFFTComplexReal.cxx:100
 TFFTComplexReal.cxx:101
 TFFTComplexReal.cxx:102
 TFFTComplexReal.cxx:103
 TFFTComplexReal.cxx:104
 TFFTComplexReal.cxx:105
 TFFTComplexReal.cxx:106
 TFFTComplexReal.cxx:107
 TFFTComplexReal.cxx:108
 TFFTComplexReal.cxx:109
 TFFTComplexReal.cxx:110
 TFFTComplexReal.cxx:111
 TFFTComplexReal.cxx:112
 TFFTComplexReal.cxx:113
 TFFTComplexReal.cxx:114
 TFFTComplexReal.cxx:115
 TFFTComplexReal.cxx:116
 TFFTComplexReal.cxx:117
 TFFTComplexReal.cxx:118
 TFFTComplexReal.cxx:119
 TFFTComplexReal.cxx:120
 TFFTComplexReal.cxx:121
 TFFTComplexReal.cxx:122
 TFFTComplexReal.cxx:123
 TFFTComplexReal.cxx:124
 TFFTComplexReal.cxx:125
 TFFTComplexReal.cxx:126
 TFFTComplexReal.cxx:127
 TFFTComplexReal.cxx:128
 TFFTComplexReal.cxx:129
 TFFTComplexReal.cxx:130
 TFFTComplexReal.cxx:131
 TFFTComplexReal.cxx:132
 TFFTComplexReal.cxx:133
 TFFTComplexReal.cxx:134
 TFFTComplexReal.cxx:135
 TFFTComplexReal.cxx:136
 TFFTComplexReal.cxx:137
 TFFTComplexReal.cxx:138
 TFFTComplexReal.cxx:139
 TFFTComplexReal.cxx:140
 TFFTComplexReal.cxx:141
 TFFTComplexReal.cxx:142
 TFFTComplexReal.cxx:143
 TFFTComplexReal.cxx:144
 TFFTComplexReal.cxx:145
 TFFTComplexReal.cxx:146
 TFFTComplexReal.cxx:147
 TFFTComplexReal.cxx:148
 TFFTComplexReal.cxx:149
 TFFTComplexReal.cxx:150
 TFFTComplexReal.cxx:151
 TFFTComplexReal.cxx:152
 TFFTComplexReal.cxx:153
 TFFTComplexReal.cxx:154
 TFFTComplexReal.cxx:155
 TFFTComplexReal.cxx:156
 TFFTComplexReal.cxx:157
 TFFTComplexReal.cxx:158
 TFFTComplexReal.cxx:159
 TFFTComplexReal.cxx:160
 TFFTComplexReal.cxx:161
 TFFTComplexReal.cxx:162
 TFFTComplexReal.cxx:163
 TFFTComplexReal.cxx:164
 TFFTComplexReal.cxx:165
 TFFTComplexReal.cxx:166
 TFFTComplexReal.cxx:167
 TFFTComplexReal.cxx:168
 TFFTComplexReal.cxx:169
 TFFTComplexReal.cxx:170
 TFFTComplexReal.cxx:171
 TFFTComplexReal.cxx:172
 TFFTComplexReal.cxx:173
 TFFTComplexReal.cxx:174
 TFFTComplexReal.cxx:175
 TFFTComplexReal.cxx:176
 TFFTComplexReal.cxx:177
 TFFTComplexReal.cxx:178
 TFFTComplexReal.cxx:179
 TFFTComplexReal.cxx:180
 TFFTComplexReal.cxx:181
 TFFTComplexReal.cxx:182
 TFFTComplexReal.cxx:183
 TFFTComplexReal.cxx:184
 TFFTComplexReal.cxx:185
 TFFTComplexReal.cxx:186
 TFFTComplexReal.cxx:187
 TFFTComplexReal.cxx:188
 TFFTComplexReal.cxx:189
 TFFTComplexReal.cxx:190
 TFFTComplexReal.cxx:191
 TFFTComplexReal.cxx:192
 TFFTComplexReal.cxx:193
 TFFTComplexReal.cxx:194
 TFFTComplexReal.cxx:195
 TFFTComplexReal.cxx:196
 TFFTComplexReal.cxx:197
 TFFTComplexReal.cxx:198
 TFFTComplexReal.cxx:199
 TFFTComplexReal.cxx:200
 TFFTComplexReal.cxx:201
 TFFTComplexReal.cxx:202
 TFFTComplexReal.cxx:203
 TFFTComplexReal.cxx:204
 TFFTComplexReal.cxx:205
 TFFTComplexReal.cxx:206
 TFFTComplexReal.cxx:207
 TFFTComplexReal.cxx:208
 TFFTComplexReal.cxx:209
 TFFTComplexReal.cxx:210
 TFFTComplexReal.cxx:211
 TFFTComplexReal.cxx:212
 TFFTComplexReal.cxx:213
 TFFTComplexReal.cxx:214
 TFFTComplexReal.cxx:215
 TFFTComplexReal.cxx:216
 TFFTComplexReal.cxx:217
 TFFTComplexReal.cxx:218
 TFFTComplexReal.cxx:219
 TFFTComplexReal.cxx:220
 TFFTComplexReal.cxx:221
 TFFTComplexReal.cxx:222
 TFFTComplexReal.cxx:223
 TFFTComplexReal.cxx:224
 TFFTComplexReal.cxx:225
 TFFTComplexReal.cxx:226
 TFFTComplexReal.cxx:227
 TFFTComplexReal.cxx:228
 TFFTComplexReal.cxx:229
 TFFTComplexReal.cxx:230
 TFFTComplexReal.cxx:231
 TFFTComplexReal.cxx:232
 TFFTComplexReal.cxx:233
 TFFTComplexReal.cxx:234
 TFFTComplexReal.cxx:235
 TFFTComplexReal.cxx:236
 TFFTComplexReal.cxx:237
 TFFTComplexReal.cxx:238
 TFFTComplexReal.cxx:239
 TFFTComplexReal.cxx:240
 TFFTComplexReal.cxx:241
 TFFTComplexReal.cxx:242
 TFFTComplexReal.cxx:243
 TFFTComplexReal.cxx:244
 TFFTComplexReal.cxx:245
 TFFTComplexReal.cxx:246
 TFFTComplexReal.cxx:247
 TFFTComplexReal.cxx:248
 TFFTComplexReal.cxx:249
 TFFTComplexReal.cxx:250
 TFFTComplexReal.cxx:251
 TFFTComplexReal.cxx:252
 TFFTComplexReal.cxx:253
 TFFTComplexReal.cxx:254
 TFFTComplexReal.cxx:255
 TFFTComplexReal.cxx:256
 TFFTComplexReal.cxx:257
 TFFTComplexReal.cxx:258
 TFFTComplexReal.cxx:259
 TFFTComplexReal.cxx:260
 TFFTComplexReal.cxx:261
 TFFTComplexReal.cxx:262
 TFFTComplexReal.cxx:263
 TFFTComplexReal.cxx:264
 TFFTComplexReal.cxx:265
 TFFTComplexReal.cxx:266
 TFFTComplexReal.cxx:267
 TFFTComplexReal.cxx:268
 TFFTComplexReal.cxx:269
 TFFTComplexReal.cxx:270
 TFFTComplexReal.cxx:271
 TFFTComplexReal.cxx:272
 TFFTComplexReal.cxx:273
 TFFTComplexReal.cxx:274
 TFFTComplexReal.cxx:275
 TFFTComplexReal.cxx:276
 TFFTComplexReal.cxx:277
 TFFTComplexReal.cxx:278
 TFFTComplexReal.cxx:279
 TFFTComplexReal.cxx:280
 TFFTComplexReal.cxx:281
 TFFTComplexReal.cxx:282
 TFFTComplexReal.cxx:283
 TFFTComplexReal.cxx:284
 TFFTComplexReal.cxx:285
 TFFTComplexReal.cxx:286
 TFFTComplexReal.cxx:287
 TFFTComplexReal.cxx:288
 TFFTComplexReal.cxx:289
 TFFTComplexReal.cxx:290
 TFFTComplexReal.cxx:291
 TFFTComplexReal.cxx:292
 TFFTComplexReal.cxx:293
 TFFTComplexReal.cxx:294
 TFFTComplexReal.cxx:295
 TFFTComplexReal.cxx:296
 TFFTComplexReal.cxx:297
 TFFTComplexReal.cxx:298
 TFFTComplexReal.cxx:299
 TFFTComplexReal.cxx:300
 TFFTComplexReal.cxx:301
 TFFTComplexReal.cxx:302
 TFFTComplexReal.cxx:303
 TFFTComplexReal.cxx:304
 TFFTComplexReal.cxx:305
 TFFTComplexReal.cxx:306
 TFFTComplexReal.cxx:307
 TFFTComplexReal.cxx:308
 TFFTComplexReal.cxx:309
 TFFTComplexReal.cxx:310
 TFFTComplexReal.cxx:311
 TFFTComplexReal.cxx:312
 TFFTComplexReal.cxx:313
 TFFTComplexReal.cxx:314
 TFFTComplexReal.cxx:315
 TFFTComplexReal.cxx:316
 TFFTComplexReal.cxx:317
 TFFTComplexReal.cxx:318
 TFFTComplexReal.cxx:319
 TFFTComplexReal.cxx:320
 TFFTComplexReal.cxx:321
 TFFTComplexReal.cxx:322
 TFFTComplexReal.cxx:323
 TFFTComplexReal.cxx:324
 TFFTComplexReal.cxx:325
 TFFTComplexReal.cxx:326
 TFFTComplexReal.cxx:327
 TFFTComplexReal.cxx:328
 TFFTComplexReal.cxx:329
 TFFTComplexReal.cxx:330
 TFFTComplexReal.cxx:331
 TFFTComplexReal.cxx:332
 TFFTComplexReal.cxx:333
 TFFTComplexReal.cxx:334
 TFFTComplexReal.cxx:335
 TFFTComplexReal.cxx:336
 TFFTComplexReal.cxx:337
 TFFTComplexReal.cxx:338
 TFFTComplexReal.cxx:339
 TFFTComplexReal.cxx:340
 TFFTComplexReal.cxx:341
 TFFTComplexReal.cxx:342
 TFFTComplexReal.cxx:343
 TFFTComplexReal.cxx:344
 TFFTComplexReal.cxx:345
 TFFTComplexReal.cxx:346
 TFFTComplexReal.cxx:347
 TFFTComplexReal.cxx:348
 TFFTComplexReal.cxx:349
 TFFTComplexReal.cxx:350
 TFFTComplexReal.cxx:351
 TFFTComplexReal.cxx:352
 TFFTComplexReal.cxx:353
 TFFTComplexReal.cxx:354
 TFFTComplexReal.cxx:355
 TFFTComplexReal.cxx:356
 TFFTComplexReal.cxx:357
 TFFTComplexReal.cxx:358
 TFFTComplexReal.cxx:359
 TFFTComplexReal.cxx:360
 TFFTComplexReal.cxx:361
 TFFTComplexReal.cxx:362
 TFFTComplexReal.cxx:363
 TFFTComplexReal.cxx:364
 TFFTComplexReal.cxx:365
 TFFTComplexReal.cxx:366
 TFFTComplexReal.cxx:367
 TFFTComplexReal.cxx:368
 TFFTComplexReal.cxx:369
 TFFTComplexReal.cxx:370
 TFFTComplexReal.cxx:371
 TFFTComplexReal.cxx:372
 TFFTComplexReal.cxx:373
 TFFTComplexReal.cxx:374
 TFFTComplexReal.cxx:375
 TFFTComplexReal.cxx:376
 TFFTComplexReal.cxx:377
 TFFTComplexReal.cxx:378
 TFFTComplexReal.cxx:379
 TFFTComplexReal.cxx:380
 TFFTComplexReal.cxx:381
 TFFTComplexReal.cxx:382
 TFFTComplexReal.cxx:383
 TFFTComplexReal.cxx:384
 TFFTComplexReal.cxx:385
 TFFTComplexReal.cxx:386
 TFFTComplexReal.cxx:387
 TFFTComplexReal.cxx:388
 TFFTComplexReal.cxx:389
 TFFTComplexReal.cxx:390
 TFFTComplexReal.cxx:391
 TFFTComplexReal.cxx:392
 TFFTComplexReal.cxx:393
 TFFTComplexReal.cxx:394
 TFFTComplexReal.cxx:395
 TFFTComplexReal.cxx:396
 TFFTComplexReal.cxx:397
 TFFTComplexReal.cxx:398
 TFFTComplexReal.cxx:399
 TFFTComplexReal.cxx:400
 TFFTComplexReal.cxx:401
 TFFTComplexReal.cxx:402
 TFFTComplexReal.cxx:403
 TFFTComplexReal.cxx:404