LArSoft  v06_85_00
Liquid Argon Software toolkit - http://larsoft.org/
util::LArFFT Class Reference

#include "LArFFT.h"

Public Member Functions

 LArFFT (fhicl::ParameterSet const &pset, art::ActivityRegistry &reg)
 
 ~LArFFT ()
 
template<class T >
void DoFFT (std::vector< T > &input, std::vector< TComplex > &output)
 
template<class T >
void DoInvFFT (std::vector< TComplex > &input, std::vector< T > &output)
 
template<class T >
void Deconvolute (std::vector< T > &input, std::vector< T > &respFunc)
 
template<class T >
void Deconvolute (std::vector< T > &input, std::vector< TComplex > &kern)
 
template<class T >
void Convolute (std::vector< T > &input, std::vector< T > &respFunc)
 
template<class T >
void Convolute (std::vector< T > &input, std::vector< TComplex > &kern)
 
template<class T >
void Correlate (std::vector< T > &input, std::vector< T > &respFunc)
 
template<class T >
void Correlate (std::vector< T > &input, std::vector< TComplex > &kern)
 
template<class T >
void AlignedSum (std::vector< T > &input, std::vector< T > &output, bool add=true)
 
void ShiftData (std::vector< TComplex > &input, double shift)
 
template<class T >
void ShiftData (std::vector< T > &input, double shift)
 
template<class T >
PeakCorrelation (std::vector< T > &shape1, std::vector< T > &shape2)
 
int FFTSize () const
 
std::string FFTOptions () const
 
int FFTFitBins () const
 
void ReinitializeFFT (int, std::string, int)
 

Private Member Functions

void InitializeFFT ()
 

Private Attributes

int fSize
 
int fFreqSize
 
std::string fOption
 
int fFitBins
 
TF1 * fPeakFit
 
TH1D * fConvHist
 
std::vector< TComplex > fCompTemp
 
std::vector< TComplex > fKern
 
TFFTRealComplex * fFFT
 object to do FFT More...
 
TFFTComplexReal * fInverseFFT
 object to do Inverse FF More...
 

Detailed Description

Definition at line 26 of file LArFFT.h.

Constructor & Destructor Documentation

util::LArFFT::LArFFT ( fhicl::ParameterSet const &  pset,
art::ActivityRegistry reg 
)

Definition at line 25 of file LArFFT_service.cc.

References fSize, and InitializeFFT().

26  : fSize (pset.get< int > ("FFTSize", 0))
27  , fOption (pset.get< std::string >("FFTOption"))
28  , fFitBins (pset.get< int >("FitBins"))
29 {
30 
31  // Default to the readout window size if the user didn't input
32  // a specific size
33  if(fSize <= 0)
34  fSize = art::ServiceHandle<detinfo::DetectorPropertiesService>()->provider()->ReadOutWindowSize();
35 
36  InitializeFFT();
37 }
int fSize
Definition: LArFFT.h:77
std::string fOption
Definition: LArFFT.h:79
int fFitBins
Definition: LArFFT.h:80
void InitializeFFT()
util::LArFFT::~LArFFT ( )

Definition at line 66 of file LArFFT_service.cc.

References fConvHist, fFFT, fInverseFFT, and fPeakFit.

67 {
68  delete fFFT;
69  delete fInverseFFT;
70  delete fPeakFit;
71  delete fConvHist;
72 }
TFFTRealComplex * fFFT
object to do FFT
Definition: LArFFT.h:86
TF1 * fPeakFit
Definition: LArFFT.h:81
TFFTComplexReal * fInverseFFT
object to do Inverse FF
Definition: LArFFT.h:87
TH1D * fConvHist
Definition: LArFFT.h:82

Member Function Documentation

template<class T >
void util::LArFFT::AlignedSum ( std::vector< T > &  input,
std::vector< T > &  output,
bool  add = true 
)
inline

Definition at line 243 of file LArFFT.h.

References fSize, PeakCorrelation(), and ShiftData().

Referenced by caldata::CalWireAna::analyze(), detsim::SimWire::ConvoluteResponseFunctions(), and detsim::SimWireT962::ConvoluteResponseFunctions().

246 {
247  double shift = PeakCorrelation(shape1,shape2);
248 
249  ShiftData(shape1,shift);
250 
251  if(add)for(int i = 0; i < fSize; i++) shape1[i]+=shape2[i];
252 
253  return;
254 }
void ShiftData(std::vector< TComplex > &input, double shift)
int fSize
Definition: LArFFT.h:77
T PeakCorrelation(std::vector< T > &shape1, std::vector< T > &shape2)
Definition: LArFFT.h:271
template<class T >
void util::LArFFT::Convolute ( std::vector< T > &  input,
std::vector< T > &  respFunc 
)
inline

Definition at line 172 of file LArFFT.h.

References DoFFT(), DoInvFFT(), fCompTemp, fFreqSize, and fKern.

Referenced by util::SignalShaping::Convolute(), util::SignalShaping::Deconvolute(), caldata::CalWire::produce(), caldata::CalWireT962::produce(), detsim::SimWire::produce(), and detsim::SimWireT962::produce().

174 {
175  DoFFT(shape1, fKern);
176  DoFFT(shape2, fCompTemp);
177 
178  for(int i = 0; i < fFreqSize; i++)
179  fCompTemp[i]*=fKern[i];
180 
181  DoInvFFT(fCompTemp, shape1);
182 
183  return;
184 }
std::vector< TComplex > fKern
Definition: LArFFT.h:84
void DoFFT(std::vector< T > &input, std::vector< TComplex > &output)
Definition: LArFFT.h:97
void DoInvFFT(std::vector< TComplex > &input, std::vector< T > &output)
Definition: LArFFT.h:119
std::vector< TComplex > fCompTemp
Definition: LArFFT.h:83
int fFreqSize
Definition: LArFFT.h:78
template<class T >
void util::LArFFT::Convolute ( std::vector< T > &  input,
std::vector< TComplex > &  kern 
)
inline

Definition at line 191 of file LArFFT.h.

References DoFFT(), DoInvFFT(), fCompTemp, and fFreqSize.

193 {
194  DoFFT(input, fCompTemp);
195 
196  for(int i = 0; i < fFreqSize; i++)
197  fCompTemp[i]*=kern[i];
198 
199  DoInvFFT(fCompTemp, input);
200 
201  return;
202 }
void DoFFT(std::vector< T > &input, std::vector< TComplex > &output)
Definition: LArFFT.h:97
void DoInvFFT(std::vector< TComplex > &input, std::vector< T > &output)
Definition: LArFFT.h:119
std::vector< TComplex > fCompTemp
Definition: LArFFT.h:83
int fFreqSize
Definition: LArFFT.h:78
template<class T >
void util::LArFFT::Correlate ( std::vector< T > &  input,
std::vector< T > &  respFunc 
)
inline

Definition at line 206 of file LArFFT.h.

References DoFFT(), DoInvFFT(), fCompTemp, fFreqSize, and fKern.

Referenced by PeakCorrelation().

208 {
209  DoFFT(shape1, fKern);
210  DoFFT(shape2, fCompTemp);
211 
212  for(int i = 0; i < fFreqSize; i++)
213  fCompTemp[i]*=TComplex::Conjugate(fKern[i]);
214 
215  DoInvFFT(fCompTemp, shape1);
216 
217  return;
218 }
std::vector< TComplex > fKern
Definition: LArFFT.h:84
void DoFFT(std::vector< T > &input, std::vector< TComplex > &output)
Definition: LArFFT.h:97
void DoInvFFT(std::vector< TComplex > &input, std::vector< T > &output)
Definition: LArFFT.h:119
std::vector< TComplex > fCompTemp
Definition: LArFFT.h:83
int fFreqSize
Definition: LArFFT.h:78
template<class T >
void util::LArFFT::Correlate ( std::vector< T > &  input,
std::vector< TComplex > &  kern 
)
inline

Definition at line 225 of file LArFFT.h.

References DoFFT(), DoInvFFT(), fCompTemp, and fFreqSize.

227 {
228  DoFFT(input, fCompTemp);
229 
230  for(int i = 0; i < fFreqSize; i++)
231  fCompTemp[i]*=TComplex::Conjugate(kern[i]);
232 
233  DoInvFFT(fCompTemp, input);
234 
235  return;
236 }
void DoFFT(std::vector< T > &input, std::vector< TComplex > &output)
Definition: LArFFT.h:97
void DoInvFFT(std::vector< TComplex > &input, std::vector< T > &output)
Definition: LArFFT.h:119
std::vector< TComplex > fCompTemp
Definition: LArFFT.h:83
int fFreqSize
Definition: LArFFT.h:78
template<class T >
void util::LArFFT::Deconvolute ( std::vector< T > &  input,
std::vector< T > &  respFunc 
)
inline

Definition at line 137 of file LArFFT.h.

References DoFFT(), DoInvFFT(), fCompTemp, fFreqSize, and fKern.

139 {
140  DoFFT(respFunction, fKern);
141  DoFFT(input, fCompTemp);
142 
143  for(int i = 0; i < fFreqSize; i++)
144  fCompTemp[i]/=fKern[i];
145 
146  DoInvFFT(fCompTemp, input);
147 
148  return;
149 }
std::vector< TComplex > fKern
Definition: LArFFT.h:84
void DoFFT(std::vector< T > &input, std::vector< TComplex > &output)
Definition: LArFFT.h:97
void DoInvFFT(std::vector< TComplex > &input, std::vector< T > &output)
Definition: LArFFT.h:119
std::vector< TComplex > fCompTemp
Definition: LArFFT.h:83
int fFreqSize
Definition: LArFFT.h:78
template<class T >
void util::LArFFT::Deconvolute ( std::vector< T > &  input,
std::vector< TComplex > &  kern 
)
inline

Definition at line 156 of file LArFFT.h.

References DoFFT(), DoInvFFT(), fCompTemp, and fFreqSize.

158 {
159  DoFFT(input, fCompTemp);
160 
161  for(int i = 0; i < fFreqSize; i++)
162  fCompTemp[i]/=kern[i];
163 
164  DoInvFFT(fCompTemp, input);
165 
166  return;
167 }
void DoFFT(std::vector< T > &input, std::vector< TComplex > &output)
Definition: LArFFT.h:97
void DoInvFFT(std::vector< TComplex > &input, std::vector< T > &output)
Definition: LArFFT.h:119
std::vector< TComplex > fCompTemp
Definition: LArFFT.h:83
int fFreqSize
Definition: LArFFT.h:78
template<class T >
void util::LArFFT::DoFFT ( std::vector< T > &  input,
std::vector< TComplex > &  output 
)
inline

Definition at line 97 of file LArFFT.h.

References fFFT, and fFreqSize.

Referenced by util::SignalShaping::AddResponseFunction(), caldata::CalWireAna::analyze(), Convolute(), detsim::SimWire::ConvoluteResponseFunctions(), detsim::SimWireT962::ConvoluteResponseFunctions(), Correlate(), Deconvolute(), and ShiftData().

99 {
100  double real = 0.; //real value holder
101  double imaginary = 0.; //imaginary value hold
102 
103  // set the points
104  for(size_t p = 0; p < input.size(); ++p)
105  fFFT->SetPoint(p, input[p]);
106 
107  fFFT->Transform();
108 
109  for(int i = 0; i < fFreqSize; ++i){
110  fFFT->GetPointComplex(i, real, imaginary);
111  output[i]=TComplex(real, imaginary);
112  }
113 
114  return;
115 }
TFFTRealComplex * fFFT
object to do FFT
Definition: LArFFT.h:86
int fFreqSize
Definition: LArFFT.h:78
template<class T >
void util::LArFFT::DoInvFFT ( std::vector< TComplex > &  input,
std::vector< T > &  output 
)
inline

Definition at line 119 of file LArFFT.h.

References fFreqSize, fInverseFFT, and fSize.

Referenced by util::SignalShaping::AddResponseFunction(), util::SignalShaping::CalculateDeconvKernel(), Convolute(), Correlate(), Deconvolute(), detsim::SimWire::GenNoise(), detsim::SimWireT962::GenNoise(), ShiftData(), and util::SignalShaping::ShiftResponseTime().

121 {
122  for(int i = 0; i < fFreqSize; ++i)
123  fInverseFFT->SetPointComplex(i, input[i]);
124 
125  fInverseFFT->Transform();
126  double factor = 1.0/(double) fSize;
127 
128  for(int i = 0; i < fSize; ++i)
129  output[i] = factor*fInverseFFT->GetPointReal(i,false);
130 
131  return;
132 }
int fSize
Definition: LArFFT.h:77
TFFTComplexReal * fInverseFFT
object to do Inverse FF
Definition: LArFFT.h:87
int fFreqSize
Definition: LArFFT.h:78
int util::LArFFT::FFTFitBins ( ) const
inline

Definition at line 71 of file LArFFT.h.

References fFitBins, and ReinitializeFFT().

71 { return fFitBins; }
int fFitBins
Definition: LArFFT.h:80
std::string util::LArFFT::FFTOptions ( ) const
inline

Definition at line 70 of file LArFFT.h.

References fOption.

70 { return fOption; }
std::string fOption
Definition: LArFFT.h:79
void util::LArFFT::InitializeFFT ( )
private

Definition at line 40 of file LArFFT_service.cc.

References fCompTemp, fConvHist, fFFT, fFitBins, fFreqSize, fInverseFFT, fKern, fOption, fPeakFit, and fSize.

Referenced by LArFFT(), and ReinitializeFFT().

40  {
41 
42  int i;
43  for(i = 1; i < fSize; i *= 2){ }
44  // mf::LogInfo("LArFFt") << "Requested size: " << fSize << " FFT size: " << i ;
45 
46  fSize=i;
47  fFreqSize = fSize/2+1;
48 
49  // allocate and setup Transform objects
50  fFFT = new TFFTRealComplex(fSize, false);
51  fInverseFFT = new TFFTComplexReal(fSize, false);
52 
53  int dummy[1] = {0};
54  // appears to be dummy argument from root page
55  fFFT->Init(fOption.c_str(),-1,dummy);
56  fInverseFFT->Init(fOption.c_str(),1,dummy);
57 
58  fPeakFit = new TF1("fPeakFit","gaus"); //allocate function used for peak fitting
59  fConvHist = new TH1D("fConvHist","Convolution Peak Data",fFitBins,0,fFitBins); //allocate histogram for peak fitting
60  //allocate other data vectors
61  fCompTemp.resize(fFreqSize);
62  fKern.resize(fFreqSize);
63 }
int fSize
Definition: LArFFT.h:77
std::vector< TComplex > fKern
Definition: LArFFT.h:84
TFFTRealComplex * fFFT
object to do FFT
Definition: LArFFT.h:86
TF1 * fPeakFit
Definition: LArFFT.h:81
std::vector< TComplex > fCompTemp
Definition: LArFFT.h:83
std::string fOption
Definition: LArFFT.h:79
TFFTComplexReal * fInverseFFT
object to do Inverse FF
Definition: LArFFT.h:87
TH1D * fConvHist
Definition: LArFFT.h:82
int fFitBins
Definition: LArFFT.h:80
int fFreqSize
Definition: LArFFT.h:78
template<class T >
T util::LArFFT::PeakCorrelation ( std::vector< T > &  shape1,
std::vector< T > &  shape2 
)
inline

Definition at line 271 of file LArFFT.h.

References Correlate(), DECLARE_ART_SERVICE, fConvHist, fFitBins, fPeakFit, and fSize.

Referenced by AlignedSum(), and util::SignalShaping::SetPeakResponseTime().

273 {
274  fConvHist->Reset("ICE");
275  std::vector<T> holder = shape1;
276  Correlate(holder,shape2);
277 
278  int maxT = max_element(holder.begin(), holder.end())-holder.begin();
279  float startT = maxT-fFitBins/2;
280  int offset = 0;
281 
282  for(int i = 0; i < fFitBins; i++) {
283  if(startT+i < 0) offset=fSize;
284  else if(startT+i > fSize) offset=-fSize;
285  else offset = 0;
286  fConvHist->Fill(i,holder[i+startT+offset]);
287  }
288 
289  fPeakFit->SetParameters(fConvHist->GetMaximum(),fFitBins/2,fFitBins/2);
290  fConvHist->Fit(fPeakFit,"QWNR","",0,fFitBins);
291  return fPeakFit->GetParameter(1)+startT;
292 }
int fSize
Definition: LArFFT.h:77
void Correlate(std::vector< T > &input, std::vector< T > &respFunc)
Definition: LArFFT.h:206
TF1 * fPeakFit
Definition: LArFFT.h:81
TH1D * fConvHist
Definition: LArFFT.h:82
int fFitBins
Definition: LArFFT.h:80
void util::LArFFT::ReinitializeFFT ( int  size,
std::string  option,
int  fitbins 
)

Definition at line 75 of file LArFFT_service.cc.

References fConvHist, fFFT, fFitBins, fInverseFFT, fOption, fPeakFit, fSize, and InitializeFFT().

Referenced by FFTFitBins().

75  {
76 
77  //delete these, which will be remade
78  delete fFFT;
79  delete fInverseFFT;
80  delete fPeakFit;
81  delete fConvHist;
82 
83  //set members
84  fSize = size;
85  fOption = option;
86  fFitBins = fitbins;
87 
88  //now initialize
89  InitializeFFT();
90 }
int fSize
Definition: LArFFT.h:77
TFFTRealComplex * fFFT
object to do FFT
Definition: LArFFT.h:86
TF1 * fPeakFit
Definition: LArFFT.h:81
std::string fOption
Definition: LArFFT.h:79
TFFTComplexReal * fInverseFFT
object to do Inverse FF
Definition: LArFFT.h:87
TH1D * fConvHist
Definition: LArFFT.h:82
int fFitBins
Definition: LArFFT.h:80
void InitializeFFT()
void util::LArFFT::ShiftData ( std::vector< TComplex > &  input,
double  shift 
)

Definition at line 107 of file LArFFT_service.cc.

References DEFINE_ART_SERVICE, fFreqSize, and fSize.

Referenced by AlignedSum(), ShiftData(), and util::SignalShaping::ShiftResponseTime().

109 {
110  double factor = -2.0*TMath::Pi()*shift/(double)fSize;
111 
112  for(int i = 0; i < fFreqSize; i++)
113  input[i] *= TComplex::Exp(TComplex(0,factor*(double)i));
114 
115  return;
116 }
int fSize
Definition: LArFFT.h:77
int fFreqSize
Definition: LArFFT.h:78
template<class T >
void util::LArFFT::ShiftData ( std::vector< T > &  input,
double  shift 
)
inline

Definition at line 258 of file LArFFT.h.

References DoFFT(), DoInvFFT(), fCompTemp, and ShiftData().

260 {
261  DoFFT(input,fCompTemp);
262  ShiftData(fCompTemp,shift);
263  DoInvFFT(fCompTemp,input);
264 
265  return;
266 }
void ShiftData(std::vector< TComplex > &input, double shift)
void DoFFT(std::vector< T > &input, std::vector< TComplex > &output)
Definition: LArFFT.h:97
void DoInvFFT(std::vector< TComplex > &input, std::vector< T > &output)
Definition: LArFFT.h:119
std::vector< TComplex > fCompTemp
Definition: LArFFT.h:83

Member Data Documentation

std::vector<TComplex> util::LArFFT::fCompTemp
private

Definition at line 83 of file LArFFT.h.

Referenced by Convolute(), Correlate(), Deconvolute(), InitializeFFT(), and ShiftData().

TH1D* util::LArFFT::fConvHist
private

Definition at line 82 of file LArFFT.h.

Referenced by InitializeFFT(), PeakCorrelation(), ReinitializeFFT(), and ~LArFFT().

TFFTRealComplex* util::LArFFT::fFFT
private

object to do FFT

Definition at line 86 of file LArFFT.h.

Referenced by DoFFT(), InitializeFFT(), ReinitializeFFT(), and ~LArFFT().

int util::LArFFT::fFitBins
private

Definition at line 80 of file LArFFT.h.

Referenced by FFTFitBins(), InitializeFFT(), PeakCorrelation(), and ReinitializeFFT().

int util::LArFFT::fFreqSize
private

Definition at line 78 of file LArFFT.h.

Referenced by Convolute(), Correlate(), Deconvolute(), DoFFT(), DoInvFFT(), InitializeFFT(), and ShiftData().

TFFTComplexReal* util::LArFFT::fInverseFFT
private

object to do Inverse FF

Definition at line 87 of file LArFFT.h.

Referenced by DoInvFFT(), InitializeFFT(), ReinitializeFFT(), and ~LArFFT().

std::vector<TComplex> util::LArFFT::fKern
private

Definition at line 84 of file LArFFT.h.

Referenced by Convolute(), Correlate(), Deconvolute(), and InitializeFFT().

std::string util::LArFFT::fOption
private

Definition at line 79 of file LArFFT.h.

Referenced by FFTOptions(), InitializeFFT(), and ReinitializeFFT().

TF1* util::LArFFT::fPeakFit
private

Definition at line 81 of file LArFFT.h.

Referenced by InitializeFFT(), PeakCorrelation(), ReinitializeFFT(), and ~LArFFT().

int util::LArFFT::fSize
private

The documentation for this class was generated from the following files: