From 2b7a391bbb269dce1df7087a71d4c1917d1ba2ea Mon Sep 17 00:00:00 2001 From: Andreas Suter Date: Mon, 29 May 2017 17:54:57 +0200 Subject: [PATCH] back-integrate DKS Fourier after DKS upgrade and modularization. Not yet tested other than compilation. --- src/classes/PFourier.cpp | 46 +--------------------------------------- src/include/PFourier.h | 11 ++-------- 2 files changed, 3 insertions(+), 54 deletions(-) diff --git a/src/classes/PFourier.cpp b/src/classes/PFourier.cpp index 9c1cd54f..7c6808f4 100644 --- a/src/classes/PFourier.cpp +++ b/src/classes/PFourier.cpp @@ -8,7 +8,7 @@ ***************************************************************************/ /*************************************************************************** - * Copyright (C) 2007-2016 by Andreas Suter * + * Copyright (C) 2007-2017 by Andreas Suter * * andreas.suter@psi.ch * * * * This program is free software; you can redistribute it and/or modify * @@ -298,12 +298,10 @@ PFourier::PFourier(TH1F *data, Int_t unitTag, Double_t startTime, Double_t endTi fIn = 0; fOut = 0; -/*//as #ifdef HAVE_DKS fInDKS = 0; fOutDKS = 0; #endif -*/ SetUseFFTW(useFFTW); @@ -365,13 +363,11 @@ PFourier::PFourier(TH1F *data, Int_t unitTag, Double_t startTime, Double_t endTi fIn = (fftw_complex *)fftw_malloc(sizeof(fftw_complex)*fNoOfBins); fOut = (fftw_complex *)fftw_malloc(sizeof(fftw_complex)*fNoOfBins); } else { // try DKS -/*//as #ifdef HAVE_DKS fInDKS = new double[fNoOfBins]; unsigned int size=fNoOfBins/2+1; fOutDKS = new complex[size]; #endif -*/ } // check if memory allocation has been successful @@ -381,17 +377,12 @@ PFourier::PFourier(TH1F *data, Int_t unitTag, Double_t startTime, Double_t endTi return; } } else { // try DKS -/*//as #ifdef HAVE_DKS if ((fInDKS == 0) || (fOutDKS == 0)) { fValid = false; return; } -#else - fValid = false; - return; #endif -*/ fValid = false; return; } @@ -405,7 +396,6 @@ PFourier::PFourier(TH1F *data, Int_t unitTag, Double_t startTime, Double_t endTi fValid = false; } } else { // try DKS -/*//as #ifdef HAVE_DKS // init DKSBase fDks.setAPI("Cuda", 4); @@ -430,10 +420,7 @@ PFourier::PFourier(TH1F *data, Int_t unitTag, Double_t startTime, Double_t endTi } if ((fReal_ptr==0) || (fComp_ptr==0)) fValid = false; -#else - fValid = false; #endif -*/ fValid = false; } } @@ -454,7 +441,6 @@ PFourier::~PFourier() if (fOut) fftw_free(fOut); } else { -/*//as #ifdef HAVE_DKS // free accelerator memory fDks.freeMemory(fReal_ptr, (int)fNoOfBins); @@ -465,7 +451,6 @@ PFourier::~PFourier() if (fOut) delete [] fOutDKS; #endif -*/ } } @@ -488,7 +473,6 @@ void PFourier::Transform(UInt_t apodizationTag) if (fUseFFTW) { fftw_execute(fFFTwPlan); } else { -/*//as #ifdef HAVE_DKS int dimsize[3] = {(int)fNoOfBins, 1, 1}; int status=0, size=fNoOfBins/2+1; @@ -498,11 +482,7 @@ void PFourier::Transform(UInt_t apodizationTag) status = fDks.callR2CFFT(fReal_ptr, fComp_ptr, 1, dimsize); // read data from accelerator status = fDks.readData< complex >(fComp_ptr, fOutDKS, size); -#else - fValid = false; - return; #endif -*/ fValid = false; return; } @@ -527,7 +507,6 @@ void PFourier::Transform(UInt_t apodizationTag) fOut[i][1] = im; } } else { // try DKS -/*//as UInt_t size=fNoOfBins/2+1; for (UInt_t i=0; i(re, im); #endif } -*/ } } @@ -609,15 +587,11 @@ TH1F* PFourier::GetRealFourier(const Double_t scale) } } else { for (UInt_t i=0; iSetBinContent(i+1, scale*fOutDKS[i].real()); #else realFourier->SetBinContent(i+1, PMUSR_UNDEFINED); #endif -*/ - realFourier->SetBinContent(i+1, PMUSR_UNDEFINED); //as needs to be removed - realFourier->SetBinError(i+1, 0.0); } } return realFourier; @@ -761,15 +735,11 @@ TH1F* PFourier::GetImaginaryFourier(const Double_t scale) } } else { for (UInt_t i=0; iSetBinContent(i+1, scale*fOutDKS[i].imag()); #else imaginaryFourier->SetBinContent(i+1, PMUSR_UNDEFINED); #endif -*/ - imaginaryFourier->SetBinContent(i+1, PMUSR_UNDEFINED); //as needs to be removed - imaginaryFourier->SetBinError(i+1, 0.0); } } return imaginaryFourier; @@ -812,15 +782,11 @@ TH1F* PFourier::GetPowerFourier(const Double_t scale) } } else { for (UInt_t i=0; iSetBinContent(i+1, scale*sqrt(fOutDKS[i].real()*fOutDKS[i].real()+fOutDKS[i].imag()*fOutDKS[i].imag())); #else pwrFourier->SetBinContent(i+1, PMUSR_UNDEFINED); #endif -*/ - pwrFourier->SetBinContent(i+1, PMUSR_UNDEFINED); //as needs to be removed - pwrFourier->SetBinError(i+1, 0.0); } } return pwrFourier; @@ -863,7 +829,6 @@ TH1F* PFourier::GetPhaseFourier(const Double_t scale) re = fOut[i][0]; im = fOut[i][1]; } else { -/*//as #ifdef HAVE_DKS re = fOutDKS[i].real(); im = fOutDKS[i].imag(); @@ -871,9 +836,6 @@ TH1F* PFourier::GetPhaseFourier(const Double_t scale) re = 1.0; im = 0.0; #endif -*/ - re = 1.0; //as needs to be removed - im = 0.0; //as needs to be removed } // calculate the phase if (re == 0) { @@ -946,18 +908,14 @@ void PFourier::PrepareFFTwInputData(UInt_t apodizationTag) } } else { for (UInt_t i=0; iGetBinContent(i+t0bin) - mean; #endif -*/ } for (UInt_t i=fNoOfData-t0bin; i -using namespace std; -#include "fftw3.h" -/*//as #ifndef HAVE_DKS #include using namespace std; @@ -47,8 +42,8 @@ using namespace std; #include using namespace std; #include "DKSBase.h" +#include "DKSFFT.h" #endif -*/ #include "Minuit2/FCNBase.h" @@ -159,15 +154,13 @@ class PFourier fftw_complex *fIn; ///< real part of the Fourier transform fftw_complex *fOut; ///< imaginary part of the Fourier transform -/*//as #ifdef HAVE_DKS double *fInDKS; ///< real part of the Fourier transform complex *fOutDKS; ///< imaginary part of the Fourier transform - DKSBase fDks; ///< Dynamic Kernel Scheduler + DKSFFT fDks; ///< Dynamic Kernel Scheduler void *fReal_ptr; ///< real part of the Fourier on accelartor void *fComp_ptr; ///< imaginary part of the Fourier on the acclerator #endif -*/ virtual void PrepareFFTwInputData(UInt_t apodizationTag); virtual void ApodizeData(Int_t apodizationTag);