ROOT logo
// @(#)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.             *
 *************************************************************************/

//////////////////////////////////////////////////////////////////////////
//                                                                      
// TFFTRealComplex                                                       
//                                                                      
// 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 a real input/complex output discrete Fourier transform in 1 or more
// dimensions. However, only out-of-place transforms are now supported for transforms
// in more than 1 dimension. For detailed information about the computed transforms,
// please refer to the FFTW manual
//
// How to use it:
// 1) Create an instance of TFFTRealComplex - 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 (see function
//    comments for possible kind parameters)
// 3) Set the data (via SetPoints()or SetPoint() functions)
// 4) Run the Transform() function
// 5) Get the output (via GetPoints() or GetPoint() 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
// 
//
//////////////////////////////////////////////////////////////////////////

#include "TFFTRealComplex.h"
#include "fftw3.h"
#include "TComplex.h"


ClassImp(TFFTRealComplex)

//_____________________________________________________________________________
TFFTRealComplex::TFFTRealComplex()
{
//default

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

//_____________________________________________________________________________
TFFTRealComplex::TFFTRealComplex(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(Double_t)*n);
      fOut = fftw_malloc(sizeof(fftw_complex)*(n/2+1));
   } else {
      fIn = fftw_malloc(sizeof(Double_t)*(2*(n/2+1)));
      fOut = 0;
   }
   fN = new Int_t[1];
   fN[0] = n;
   fTotalSize = n;
   fNdim = 1;
   fPlan = 0;
   fFlags = 0;
}

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

   if (ndim>1 && inPlace==kTRUE){
      Error("TFFTRealComplex", "multidimensional in-place r2c transforms are not implemented yet");
      return;
   }
   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 sizeout = Int_t(Double_t(fTotalSize)*(n[ndim-1]/2+1)/n[ndim-1]);
   if (!inPlace){
      fIn = fftw_malloc(sizeof(Double_t)*fTotalSize);
      fOut = fftw_malloc(sizeof(fftw_complex)*sizeout);
   } else {
      fIn = fftw_malloc(sizeof(Double_t)*(2*sizeout));
      fOut = 0;
   }
   fPlan = 0;
   fFlags = 0;
}

//_____________________________________________________________________________
TFFTRealComplex::~TFFTRealComplex()
{
//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(fIn);
   fIn = 0;
   if (fOut)
      fftw_free((fftw_complex*)fOut);
   fOut = 0;
   if (fN)
      delete [] fN;
   fN = 0;
}

//_____________________________________________________________________________
void TFFTRealComplex::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_r2c(fNdim, fN, (Double_t*)fIn, (fftw_complex*)fOut,MapFlag(flags));
   else
      fPlan = (void*)fftw_plan_dft_r2c(fNdim, fN, (Double_t*)fIn, (fftw_complex*)fIn, MapFlag(flags));
}

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


   if (fPlan){
      fftw_execute((fftw_plan)fPlan);
   }
   else {
      Error("Transform", "transform hasn't been initialised");
      return;
   }
}

//_____________________________________________________________________________
void TFFTRealComplex::GetPoints(Double_t *data, Bool_t fromInput) const
{
//Fills the array data with the computed transform.
//Only (roughly) a half of the transform is copied (exactly the output of FFTW),
//the rest being Hermitian symmetric with the first half

   if (fromInput){
      for (Int_t i=0; i<fTotalSize; i++)
         data[i] = ((Double_t*)fIn)[i];
   } else {
      Int_t realN = 2*Int_t(Double_t(fTotalSize)*(fN[fNdim-1]/2+1)/fN[fNdim-1]);
      if (fOut){
         for (Int_t i=0; i<realN; i+=2){
            data[i] = ((fftw_complex*)fOut)[i/2][0];
            data[i+1] = ((fftw_complex*)fOut)[i/2][1];
         }
      }
      else {
         for (Int_t i=0; i<realN; i++)
            data[i] = ((Double_t*)fIn)[i];
      }
   }
}

//_____________________________________________________________________________
Double_t TFFTRealComplex::GetPointReal(Int_t ipoint, Bool_t fromInput) const
{
//Returns the real part of the point #ipoint from the output or the point #ipoint
//from the input

   if (fOut && !fromInput){
      Warning("GetPointReal", "Output is complex. Only real part returned");
      return ((fftw_complex*)fOut)[ipoint][0];
   }
   else
      return ((Double_t*)fIn)[ipoint];
}

//_____________________________________________________________________________
Double_t TFFTRealComplex::GetPointReal(const Int_t *ipoint, Bool_t fromInput) const
{
//Returns the real part of the point #ipoint from the output or the point #ipoint
//from the input

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

    if (fOut && !fromInput){
      Warning("GetPointReal", "Output is complex. Only real part returned");
      return ((fftw_complex*)fOut)[ireal][0];
   }
   else
      return ((Double_t*)fIn)[ireal];
}


//_____________________________________________________________________________
void TFFTRealComplex::GetPointComplex(Int_t ipoint, Double_t &re, Double_t &im, Bool_t fromInput) const
{
//Returns the point #ipoint.
//For 1d, if ipoint > fN/2+1 (the point is in the Hermitian symmetric part), it is still
//returned. For >1d, only the first (roughly)half of points can be returned
//For 2d, see function GetPointComplex(Int_t *ipoint,...)

   if (fromInput){
      re = ((Double_t*)fIn)[ipoint];
   } else {
      if (fNdim==1){
         if (fOut){
            if (ipoint<fN[0]/2+1){
               re = ((fftw_complex*)fOut)[ipoint][0];
               im = ((fftw_complex*)fOut)[ipoint][1];
            } else {
               re = ((fftw_complex*)fOut)[fN[0]-ipoint][0];
               im = -((fftw_complex*)fOut)[fN[0]-ipoint][1];
            }
         } else {
            if (ipoint<fN[0]/2+1){
               re = ((Double_t*)fIn)[2*ipoint];
               im = ((Double_t*)fIn)[2*ipoint+1];
            } else {
               re = ((Double_t*)fIn)[2*(fN[0]-ipoint)];
               im = ((Double_t*)fIn)[2*(fN[0]-ipoint)+1];
            }
         }
      }
      else {
         Int_t realN = 2*Int_t(Double_t(fTotalSize)*(fN[fNdim-1]/2+1)/fN[fNdim-1]);
         if (ipoint>realN/2){
            Error("GetPointComplex", "Illegal index value");
            return;
         }
         if (fOut){
            re = ((fftw_complex*)fOut)[ipoint][0];
            im = ((fftw_complex*)fOut)[ipoint][1];
         } else {
            re = ((Double_t*)fIn)[2*ipoint];
            im = ((Double_t*)fIn)[2*ipoint+1];
         }
      }
   }
}
//_____________________________________________________________________________
void TFFTRealComplex::GetPointComplex(const Int_t *ipoint, Double_t &re, Double_t &im, Bool_t fromInput) const
{
//For multidimensional transforms. Returns the point #ipoint.
//In case of transforms of more than 2 dimensions,
//only points from the first (roughly)half are returned, the rest being Hermitian symmetric


   Int_t ireal = ipoint[0];
   for (Int_t i=0; i<fNdim-2; i++)
      ireal=fN[i+1]*ireal + ipoint[i+1];
   //special treatment of the last dimension
   ireal = (fN[fNdim-1]/2+1)*ireal + ipoint[fNdim-1];

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


//_____________________________________________________________________________
Double_t* TFFTRealComplex::GetPointsReal(Bool_t fromInput) const
{
//Returns the input array// 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.

   if (!fromInput){
      Error("GetPointsReal", "Output array is complex");
      return 0;
   }
   return (Double_t*)fIn;
}

//_____________________________________________________________________________
void TFFTRealComplex::GetPointsComplex(Double_t *re, Double_t *im, Bool_t fromInput) const
{
//Fills the argument arrays with the real and imaginary parts of the computed transform.
//Only (roughly) a half of the transform is copied, the rest being Hermitian
//symmetric with the first half

   Int_t nreal;
   if (fOut && !fromInput){
      nreal = Int_t(Double_t(fTotalSize)*(fN[fNdim-1]/2+1)/fN[fNdim-1]);
      for (Int_t i=0; i<nreal; i++){
         re[i] = ((fftw_complex*)fOut)[i][0];
         im[i] = ((fftw_complex*)fOut)[i][1];
      }
   } else if (fromInput) {
      for (Int_t i=0; i<fTotalSize; i++){
         re[i] = ((Double_t *)fIn)[i];
         im[i] = 0;
      }
   }
   else {
      nreal = 2*Int_t(Double_t(fTotalSize)*(fN[fNdim-1]/2+1)/fN[fNdim-1]);
      for (Int_t i=0; i<nreal; i+=2){
         re[i/2] = ((Double_t*)fIn)[i];
         im[i/2] = ((Double_t*)fIn)[i+1];
      }
   }
}

//_____________________________________________________________________________
void TFFTRealComplex::GetPointsComplex(Double_t *data, Bool_t fromInput) const
{
//Fills the argument arrays with the real and imaginary parts of the computed transform.
//Only (roughly) a half of the transform is copied, the rest being Hermitian
//symmetric with the first half

   Int_t nreal;

   if (fOut && !fromInput){
      nreal = Int_t(Double_t(fTotalSize)*(fN[fNdim-1]/2+1)/fN[fNdim-1]);

      for (Int_t i=0; i<nreal; i+=2){
         data[i] = ((fftw_complex*)fOut)[i/2][0];
         data[i+1] = ((fftw_complex*)fOut)[i/2][1];
      }
   } else if (fromInput){
      for (Int_t i=0; i<fTotalSize; i+=2){
         data[i] = ((Double_t*)fIn)[i/2];
         data[i+1] = 0;
      }
   }
   else {

      nreal = 2*Int_t(Double_t(fTotalSize)*(fN[fNdim-1]/2+1)/fN[fNdim-1]);
      for (Int_t i=0; i<nreal; i++)
         data[i] = ((Double_t*)fIn)[i];
   }
}

//_____________________________________________________________________________
void TFFTRealComplex::SetPoint(Int_t ipoint, Double_t re, Double_t /*im*/)
{
//Set the point #ipoint

   ((Double_t *)fIn)[ipoint] = re;
}

//_____________________________________________________________________________
void TFFTRealComplex::SetPoint(const Int_t *ipoint, Double_t re, Double_t /*im*/)
{
//For multidimensional transforms. Set the point #ipoint

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

//_____________________________________________________________________________
void TFFTRealComplex::SetPoints(const Double_t *data)
{
//Set all input points

   for (Int_t i=0; i<fTotalSize; i++){
      ((Double_t*)fIn)[i]=data[i];
   }
}

//_____________________________________________________________________________
void TFFTRealComplex::SetPointComplex(Int_t ipoint, TComplex &c)
{
//Sets the point #ipoint (only the real part of the argument is taken)

   ((Double_t *)fIn)[ipoint]=c.Re();
}

//_____________________________________________________________________________
void TFFTRealComplex::SetPointsComplex(const Double_t *re, const Double_t* /*im*/)
{
//Set all points. Only the real array is used

   for (Int_t i=0; i<fTotalSize; i++)
      ((Double_t *)fIn)[i] = re[i];
}

//_____________________________________________________________________________
UInt_t TFFTRealComplex::MapFlag(Option_t *flag)
{
//allowed options:
//"ES"
//"M"
//"P"
//"EX"


   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;
}
 TFFTRealComplex.cxx:1
 TFFTRealComplex.cxx:2
 TFFTRealComplex.cxx:3
 TFFTRealComplex.cxx:4
 TFFTRealComplex.cxx:5
 TFFTRealComplex.cxx:6
 TFFTRealComplex.cxx:7
 TFFTRealComplex.cxx:8
 TFFTRealComplex.cxx:9
 TFFTRealComplex.cxx:10
 TFFTRealComplex.cxx:11
 TFFTRealComplex.cxx:12
 TFFTRealComplex.cxx:13
 TFFTRealComplex.cxx:14
 TFFTRealComplex.cxx:15
 TFFTRealComplex.cxx:16
 TFFTRealComplex.cxx:17
 TFFTRealComplex.cxx:18
 TFFTRealComplex.cxx:19
 TFFTRealComplex.cxx:20
 TFFTRealComplex.cxx:21
 TFFTRealComplex.cxx:22
 TFFTRealComplex.cxx:23
 TFFTRealComplex.cxx:24
 TFFTRealComplex.cxx:25
 TFFTRealComplex.cxx:26
 TFFTRealComplex.cxx:27
 TFFTRealComplex.cxx:28
 TFFTRealComplex.cxx:29
 TFFTRealComplex.cxx:30
 TFFTRealComplex.cxx:31
 TFFTRealComplex.cxx:32
 TFFTRealComplex.cxx:33
 TFFTRealComplex.cxx:34
 TFFTRealComplex.cxx:35
 TFFTRealComplex.cxx:36
 TFFTRealComplex.cxx:37
 TFFTRealComplex.cxx:38
 TFFTRealComplex.cxx:39
 TFFTRealComplex.cxx:40
 TFFTRealComplex.cxx:41
 TFFTRealComplex.cxx:42
 TFFTRealComplex.cxx:43
 TFFTRealComplex.cxx:44
 TFFTRealComplex.cxx:45
 TFFTRealComplex.cxx:46
 TFFTRealComplex.cxx:47
 TFFTRealComplex.cxx:48
 TFFTRealComplex.cxx:49
 TFFTRealComplex.cxx:50
 TFFTRealComplex.cxx:51
 TFFTRealComplex.cxx:52
 TFFTRealComplex.cxx:53
 TFFTRealComplex.cxx:54
 TFFTRealComplex.cxx:55
 TFFTRealComplex.cxx:56
 TFFTRealComplex.cxx:57
 TFFTRealComplex.cxx:58
 TFFTRealComplex.cxx:59
 TFFTRealComplex.cxx:60
 TFFTRealComplex.cxx:61
 TFFTRealComplex.cxx:62
 TFFTRealComplex.cxx:63
 TFFTRealComplex.cxx:64
 TFFTRealComplex.cxx:65
 TFFTRealComplex.cxx:66
 TFFTRealComplex.cxx:67
 TFFTRealComplex.cxx:68
 TFFTRealComplex.cxx:69
 TFFTRealComplex.cxx:70
 TFFTRealComplex.cxx:71
 TFFTRealComplex.cxx:72
 TFFTRealComplex.cxx:73
 TFFTRealComplex.cxx:74
 TFFTRealComplex.cxx:75
 TFFTRealComplex.cxx:76
 TFFTRealComplex.cxx:77
 TFFTRealComplex.cxx:78
 TFFTRealComplex.cxx:79
 TFFTRealComplex.cxx:80
 TFFTRealComplex.cxx:81
 TFFTRealComplex.cxx:82
 TFFTRealComplex.cxx:83
 TFFTRealComplex.cxx:84
 TFFTRealComplex.cxx:85
 TFFTRealComplex.cxx:86
 TFFTRealComplex.cxx:87
 TFFTRealComplex.cxx:88
 TFFTRealComplex.cxx:89
 TFFTRealComplex.cxx:90
 TFFTRealComplex.cxx:91
 TFFTRealComplex.cxx:92
 TFFTRealComplex.cxx:93
 TFFTRealComplex.cxx:94
 TFFTRealComplex.cxx:95
 TFFTRealComplex.cxx:96
 TFFTRealComplex.cxx:97
 TFFTRealComplex.cxx:98
 TFFTRealComplex.cxx:99
 TFFTRealComplex.cxx:100
 TFFTRealComplex.cxx:101
 TFFTRealComplex.cxx:102
 TFFTRealComplex.cxx:103
 TFFTRealComplex.cxx:104
 TFFTRealComplex.cxx:105
 TFFTRealComplex.cxx:106
 TFFTRealComplex.cxx:107
 TFFTRealComplex.cxx:108
 TFFTRealComplex.cxx:109
 TFFTRealComplex.cxx:110
 TFFTRealComplex.cxx:111
 TFFTRealComplex.cxx:112
 TFFTRealComplex.cxx:113
 TFFTRealComplex.cxx:114
 TFFTRealComplex.cxx:115
 TFFTRealComplex.cxx:116
 TFFTRealComplex.cxx:117
 TFFTRealComplex.cxx:118
 TFFTRealComplex.cxx:119
 TFFTRealComplex.cxx:120
 TFFTRealComplex.cxx:121
 TFFTRealComplex.cxx:122
 TFFTRealComplex.cxx:123
 TFFTRealComplex.cxx:124
 TFFTRealComplex.cxx:125
 TFFTRealComplex.cxx:126
 TFFTRealComplex.cxx:127
 TFFTRealComplex.cxx:128
 TFFTRealComplex.cxx:129
 TFFTRealComplex.cxx:130
 TFFTRealComplex.cxx:131
 TFFTRealComplex.cxx:132
 TFFTRealComplex.cxx:133
 TFFTRealComplex.cxx:134
 TFFTRealComplex.cxx:135
 TFFTRealComplex.cxx:136
 TFFTRealComplex.cxx:137
 TFFTRealComplex.cxx:138
 TFFTRealComplex.cxx:139
 TFFTRealComplex.cxx:140
 TFFTRealComplex.cxx:141
 TFFTRealComplex.cxx:142
 TFFTRealComplex.cxx:143
 TFFTRealComplex.cxx:144
 TFFTRealComplex.cxx:145
 TFFTRealComplex.cxx:146
 TFFTRealComplex.cxx:147
 TFFTRealComplex.cxx:148
 TFFTRealComplex.cxx:149
 TFFTRealComplex.cxx:150
 TFFTRealComplex.cxx:151
 TFFTRealComplex.cxx:152
 TFFTRealComplex.cxx:153
 TFFTRealComplex.cxx:154
 TFFTRealComplex.cxx:155
 TFFTRealComplex.cxx:156
 TFFTRealComplex.cxx:157
 TFFTRealComplex.cxx:158
 TFFTRealComplex.cxx:159
 TFFTRealComplex.cxx:160
 TFFTRealComplex.cxx:161
 TFFTRealComplex.cxx:162
 TFFTRealComplex.cxx:163
 TFFTRealComplex.cxx:164
 TFFTRealComplex.cxx:165
 TFFTRealComplex.cxx:166
 TFFTRealComplex.cxx:167
 TFFTRealComplex.cxx:168
 TFFTRealComplex.cxx:169
 TFFTRealComplex.cxx:170
 TFFTRealComplex.cxx:171
 TFFTRealComplex.cxx:172
 TFFTRealComplex.cxx:173
 TFFTRealComplex.cxx:174
 TFFTRealComplex.cxx:175
 TFFTRealComplex.cxx:176
 TFFTRealComplex.cxx:177
 TFFTRealComplex.cxx:178
 TFFTRealComplex.cxx:179
 TFFTRealComplex.cxx:180
 TFFTRealComplex.cxx:181
 TFFTRealComplex.cxx:182
 TFFTRealComplex.cxx:183
 TFFTRealComplex.cxx:184
 TFFTRealComplex.cxx:185
 TFFTRealComplex.cxx:186
 TFFTRealComplex.cxx:187
 TFFTRealComplex.cxx:188
 TFFTRealComplex.cxx:189
 TFFTRealComplex.cxx:190
 TFFTRealComplex.cxx:191
 TFFTRealComplex.cxx:192
 TFFTRealComplex.cxx:193
 TFFTRealComplex.cxx:194
 TFFTRealComplex.cxx:195
 TFFTRealComplex.cxx:196
 TFFTRealComplex.cxx:197
 TFFTRealComplex.cxx:198
 TFFTRealComplex.cxx:199
 TFFTRealComplex.cxx:200
 TFFTRealComplex.cxx:201
 TFFTRealComplex.cxx:202
 TFFTRealComplex.cxx:203
 TFFTRealComplex.cxx:204
 TFFTRealComplex.cxx:205
 TFFTRealComplex.cxx:206
 TFFTRealComplex.cxx:207
 TFFTRealComplex.cxx:208
 TFFTRealComplex.cxx:209
 TFFTRealComplex.cxx:210
 TFFTRealComplex.cxx:211
 TFFTRealComplex.cxx:212
 TFFTRealComplex.cxx:213
 TFFTRealComplex.cxx:214
 TFFTRealComplex.cxx:215
 TFFTRealComplex.cxx:216
 TFFTRealComplex.cxx:217
 TFFTRealComplex.cxx:218
 TFFTRealComplex.cxx:219
 TFFTRealComplex.cxx:220
 TFFTRealComplex.cxx:221
 TFFTRealComplex.cxx:222
 TFFTRealComplex.cxx:223
 TFFTRealComplex.cxx:224
 TFFTRealComplex.cxx:225
 TFFTRealComplex.cxx:226
 TFFTRealComplex.cxx:227
 TFFTRealComplex.cxx:228
 TFFTRealComplex.cxx:229
 TFFTRealComplex.cxx:230
 TFFTRealComplex.cxx:231
 TFFTRealComplex.cxx:232
 TFFTRealComplex.cxx:233
 TFFTRealComplex.cxx:234
 TFFTRealComplex.cxx:235
 TFFTRealComplex.cxx:236
 TFFTRealComplex.cxx:237
 TFFTRealComplex.cxx:238
 TFFTRealComplex.cxx:239
 TFFTRealComplex.cxx:240
 TFFTRealComplex.cxx:241
 TFFTRealComplex.cxx:242
 TFFTRealComplex.cxx:243
 TFFTRealComplex.cxx:244
 TFFTRealComplex.cxx:245
 TFFTRealComplex.cxx:246
 TFFTRealComplex.cxx:247
 TFFTRealComplex.cxx:248
 TFFTRealComplex.cxx:249
 TFFTRealComplex.cxx:250
 TFFTRealComplex.cxx:251
 TFFTRealComplex.cxx:252
 TFFTRealComplex.cxx:253
 TFFTRealComplex.cxx:254
 TFFTRealComplex.cxx:255
 TFFTRealComplex.cxx:256
 TFFTRealComplex.cxx:257
 TFFTRealComplex.cxx:258
 TFFTRealComplex.cxx:259
 TFFTRealComplex.cxx:260
 TFFTRealComplex.cxx:261
 TFFTRealComplex.cxx:262
 TFFTRealComplex.cxx:263
 TFFTRealComplex.cxx:264
 TFFTRealComplex.cxx:265
 TFFTRealComplex.cxx:266
 TFFTRealComplex.cxx:267
 TFFTRealComplex.cxx:268
 TFFTRealComplex.cxx:269
 TFFTRealComplex.cxx:270
 TFFTRealComplex.cxx:271
 TFFTRealComplex.cxx:272
 TFFTRealComplex.cxx:273
 TFFTRealComplex.cxx:274
 TFFTRealComplex.cxx:275
 TFFTRealComplex.cxx:276
 TFFTRealComplex.cxx:277
 TFFTRealComplex.cxx:278
 TFFTRealComplex.cxx:279
 TFFTRealComplex.cxx:280
 TFFTRealComplex.cxx:281
 TFFTRealComplex.cxx:282
 TFFTRealComplex.cxx:283
 TFFTRealComplex.cxx:284
 TFFTRealComplex.cxx:285
 TFFTRealComplex.cxx:286
 TFFTRealComplex.cxx:287
 TFFTRealComplex.cxx:288
 TFFTRealComplex.cxx:289
 TFFTRealComplex.cxx:290
 TFFTRealComplex.cxx:291
 TFFTRealComplex.cxx:292
 TFFTRealComplex.cxx:293
 TFFTRealComplex.cxx:294
 TFFTRealComplex.cxx:295
 TFFTRealComplex.cxx:296
 TFFTRealComplex.cxx:297
 TFFTRealComplex.cxx:298
 TFFTRealComplex.cxx:299
 TFFTRealComplex.cxx:300
 TFFTRealComplex.cxx:301
 TFFTRealComplex.cxx:302
 TFFTRealComplex.cxx:303
 TFFTRealComplex.cxx:304
 TFFTRealComplex.cxx:305
 TFFTRealComplex.cxx:306
 TFFTRealComplex.cxx:307
 TFFTRealComplex.cxx:308
 TFFTRealComplex.cxx:309
 TFFTRealComplex.cxx:310
 TFFTRealComplex.cxx:311
 TFFTRealComplex.cxx:312
 TFFTRealComplex.cxx:313
 TFFTRealComplex.cxx:314
 TFFTRealComplex.cxx:315
 TFFTRealComplex.cxx:316
 TFFTRealComplex.cxx:317
 TFFTRealComplex.cxx:318
 TFFTRealComplex.cxx:319
 TFFTRealComplex.cxx:320
 TFFTRealComplex.cxx:321
 TFFTRealComplex.cxx:322
 TFFTRealComplex.cxx:323
 TFFTRealComplex.cxx:324
 TFFTRealComplex.cxx:325
 TFFTRealComplex.cxx:326
 TFFTRealComplex.cxx:327
 TFFTRealComplex.cxx:328
 TFFTRealComplex.cxx:329
 TFFTRealComplex.cxx:330
 TFFTRealComplex.cxx:331
 TFFTRealComplex.cxx:332
 TFFTRealComplex.cxx:333
 TFFTRealComplex.cxx:334
 TFFTRealComplex.cxx:335
 TFFTRealComplex.cxx:336
 TFFTRealComplex.cxx:337
 TFFTRealComplex.cxx:338
 TFFTRealComplex.cxx:339
 TFFTRealComplex.cxx:340
 TFFTRealComplex.cxx:341
 TFFTRealComplex.cxx:342
 TFFTRealComplex.cxx:343
 TFFTRealComplex.cxx:344
 TFFTRealComplex.cxx:345
 TFFTRealComplex.cxx:346
 TFFTRealComplex.cxx:347
 TFFTRealComplex.cxx:348
 TFFTRealComplex.cxx:349
 TFFTRealComplex.cxx:350
 TFFTRealComplex.cxx:351
 TFFTRealComplex.cxx:352
 TFFTRealComplex.cxx:353
 TFFTRealComplex.cxx:354
 TFFTRealComplex.cxx:355
 TFFTRealComplex.cxx:356
 TFFTRealComplex.cxx:357
 TFFTRealComplex.cxx:358
 TFFTRealComplex.cxx:359
 TFFTRealComplex.cxx:360
 TFFTRealComplex.cxx:361
 TFFTRealComplex.cxx:362
 TFFTRealComplex.cxx:363
 TFFTRealComplex.cxx:364
 TFFTRealComplex.cxx:365
 TFFTRealComplex.cxx:366
 TFFTRealComplex.cxx:367
 TFFTRealComplex.cxx:368
 TFFTRealComplex.cxx:369
 TFFTRealComplex.cxx:370
 TFFTRealComplex.cxx:371
 TFFTRealComplex.cxx:372
 TFFTRealComplex.cxx:373
 TFFTRealComplex.cxx:374
 TFFTRealComplex.cxx:375
 TFFTRealComplex.cxx:376
 TFFTRealComplex.cxx:377
 TFFTRealComplex.cxx:378
 TFFTRealComplex.cxx:379
 TFFTRealComplex.cxx:380
 TFFTRealComplex.cxx:381
 TFFTRealComplex.cxx:382
 TFFTRealComplex.cxx:383
 TFFTRealComplex.cxx:384
 TFFTRealComplex.cxx:385
 TFFTRealComplex.cxx:386
 TFFTRealComplex.cxx:387
 TFFTRealComplex.cxx:388
 TFFTRealComplex.cxx:389
 TFFTRealComplex.cxx:390
 TFFTRealComplex.cxx:391
 TFFTRealComplex.cxx:392
 TFFTRealComplex.cxx:393
 TFFTRealComplex.cxx:394
 TFFTRealComplex.cxx:395
 TFFTRealComplex.cxx:396
 TFFTRealComplex.cxx:397
 TFFTRealComplex.cxx:398
 TFFTRealComplex.cxx:399
 TFFTRealComplex.cxx:400
 TFFTRealComplex.cxx:401
 TFFTRealComplex.cxx:402
 TFFTRealComplex.cxx:403
 TFFTRealComplex.cxx:404
 TFFTRealComplex.cxx:405
 TFFTRealComplex.cxx:406
 TFFTRealComplex.cxx:407
 TFFTRealComplex.cxx:408
 TFFTRealComplex.cxx:409
 TFFTRealComplex.cxx:410
 TFFTRealComplex.cxx:411
 TFFTRealComplex.cxx:412
 TFFTRealComplex.cxx:413
 TFFTRealComplex.cxx:414
 TFFTRealComplex.cxx:415
 TFFTRealComplex.cxx:416
 TFFTRealComplex.cxx:417
 TFFTRealComplex.cxx:418
 TFFTRealComplex.cxx:419
 TFFTRealComplex.cxx:420
 TFFTRealComplex.cxx:421
 TFFTRealComplex.cxx:422
 TFFTRealComplex.cxx:423
 TFFTRealComplex.cxx:424
 TFFTRealComplex.cxx:425
 TFFTRealComplex.cxx:426
 TFFTRealComplex.cxx:427
 TFFTRealComplex.cxx:428
 TFFTRealComplex.cxx:429
 TFFTRealComplex.cxx:430
 TFFTRealComplex.cxx:431
 TFFTRealComplex.cxx:432
 TFFTRealComplex.cxx:433
 TFFTRealComplex.cxx:434
 TFFTRealComplex.cxx:435
 TFFTRealComplex.cxx:436
 TFFTRealComplex.cxx:437
 TFFTRealComplex.cxx:438
 TFFTRealComplex.cxx:439
 TFFTRealComplex.cxx:440
 TFFTRealComplex.cxx:441
 TFFTRealComplex.cxx:442
 TFFTRealComplex.cxx:443
 TFFTRealComplex.cxx:444
 TFFTRealComplex.cxx:445
 TFFTRealComplex.cxx:446
 TFFTRealComplex.cxx:447
 TFFTRealComplex.cxx:448
 TFFTRealComplex.cxx:449
 TFFTRealComplex.cxx:450
 TFFTRealComplex.cxx:451
 TFFTRealComplex.cxx:452
 TFFTRealComplex.cxx:453
 TFFTRealComplex.cxx:454
 TFFTRealComplex.cxx:455
 TFFTRealComplex.cxx:456
 TFFTRealComplex.cxx:457
 TFFTRealComplex.cxx:458
 TFFTRealComplex.cxx:459
 TFFTRealComplex.cxx:460
 TFFTRealComplex.cxx:461
 TFFTRealComplex.cxx:462
 TFFTRealComplex.cxx:463
 TFFTRealComplex.cxx:464
 TFFTRealComplex.cxx:465
 TFFTRealComplex.cxx:466
 TFFTRealComplex.cxx:467
 TFFTRealComplex.cxx:468
 TFFTRealComplex.cxx:469
 TFFTRealComplex.cxx:470
 TFFTRealComplex.cxx:471
 TFFTRealComplex.cxx:472
 TFFTRealComplex.cxx:473
 TFFTRealComplex.cxx:474
 TFFTRealComplex.cxx:475
 TFFTRealComplex.cxx:476
 TFFTRealComplex.cxx:477
 TFFTRealComplex.cxx:478
 TFFTRealComplex.cxx:479
 TFFTRealComplex.cxx:480
 TFFTRealComplex.cxx:481
 TFFTRealComplex.cxx:482
 TFFTRealComplex.cxx:483
 TFFTRealComplex.cxx:484
 TFFTRealComplex.cxx:485
 TFFTRealComplex.cxx:486
 TFFTRealComplex.cxx:487
 TFFTRealComplex.cxx:488
 TFFTRealComplex.cxx:489
 TFFTRealComplex.cxx:490
 TFFTRealComplex.cxx:491
 TFFTRealComplex.cxx:492
 TFFTRealComplex.cxx:493
 TFFTRealComplex.cxx:494
 TFFTRealComplex.cxx:495
 TFFTRealComplex.cxx:496
 TFFTRealComplex.cxx:497
 TFFTRealComplex.cxx:498
 TFFTRealComplex.cxx:499
 TFFTRealComplex.cxx:500
 TFFTRealComplex.cxx:501
 TFFTRealComplex.cxx:502
 TFFTRealComplex.cxx:503