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

#include "SignalShaping.h"

Public Member Functions

 SignalShaping ()
 
virtual ~SignalShaping ()
 
const std::vector< double > & Response () const
 
const std::vector< double > & Response_save () const
 
const std::vector< TComplex > & ConvKernel () const
 
const std::vector< TComplex > & Filter () const
 
const std::vector< TComplex > & DeconvKernel () const
 
template<class T >
void Convolute (std::vector< T > &func) const
 
template<class T >
void Deconvolute (std::vector< T > &func) const
 
void Reset ()
 
void save_response ()
 
void set_normflag (bool flag)
 
void AddResponseFunction (const std::vector< double > &resp, bool ResetResponse=false)
 
void ShiftResponseTime (double ticks)
 
void SetPeakResponseTime (double tick)
 
void AddFilterFunction (const std::vector< TComplex > &filt)
 
void SetDeconvKernelPolarity (int pol)
 
void LockResponse () const
 
void CalculateDeconvKernel () const
 

Private Attributes

bool fResponseLocked
 
bool fFilterLocked
 
std::vector< double > fResponse
 
std::vector< double > fResponse_save
 
std::vector< TComplex > fConvKernel
 
std::vector< TComplex > fFilter
 
std::vector< TComplex > fDeconvKernel
 
int fDeconvKernelPolarity
 
bool fNorm
 

Detailed Description

Definition at line 69 of file SignalShaping.h.

Constructor & Destructor Documentation

util::SignalShaping::SignalShaping ( )

Definition at line 19 of file SignalShaping.cxx.

20  : fResponseLocked(false)
21  , fFilterLocked (false)
22  , fNorm (true)
23 {}
util::SignalShaping::~SignalShaping ( )
virtual

Definition at line 29 of file SignalShaping.cxx.

30 {}

Member Function Documentation

void util::SignalShaping::AddFilterFunction ( const std::vector< TComplex > &  filt)

Definition at line 168 of file SignalShaping.cxx.

References fFilter, fFilterLocked, util::LArFFT::FFTSize(), filt, min, and n.

Referenced by set_normflag().

169 {
170  // Make sure configuration is not locked.
171 
172  if(fFilterLocked)
173  throw cet::exception("SignalShaping") << "Configuration locked.\n";
174 
175  // Get FFT service.
176 
178 
179  // If this is the first filter function, just copy the filter function.
180  // Otherwise, update the overall filter function.
181 
182  if(fFilter.size() == 0) {
183  fFilter = filt;
184  fFilter.resize(fft->FFTSize() / 2 + 1);
185  }
186  else {
187  unsigned int n = std::min(fFilter.size(), filt.size());
188  for(unsigned int i=0; i<n; ++i)
189  fFilter[i] *= filt[i];
190  for(unsigned int i=n; i<fFilter.size(); ++i)
191  fFilter[i] = 0.;
192  }
193 }
std::vector< TComplex > fFilter
TString filt[ntarg]
Definition: Style.C:28
int FFTSize() const
Definition: LArFFT.h:69
Int_t min
Definition: plot.C:26
Char_t n[5]
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
void util::SignalShaping::AddResponseFunction ( const std::vector< double > &  resp,
bool  ResetResponse = false 
)

Definition at line 62 of file SignalShaping.cxx.

References util::LArFFT::DoFFT(), util::LArFFT::DoInvFFT(), fConvKernel, util::LArFFT::FFTSize(), fResponse, and fResponseLocked.

Referenced by set_normflag().

63 {
64  // Make sure configuration is not locked.
65 
66  if(fResponseLocked)
67  throw cet::exception("SignalShaping") << "Configuration locked.\n";
68 
69  // Get FFT service.
70 
72  int nticks = fft->FFTSize();
73 
74  // Copy new response function into fResponse attribute, and pad or
75  // truncate to correct size.
76 
77  fResponse = resp;
78  fResponse.resize(nticks, 0.);
79 
80  // Is this the first response function?
81 
82  if ( fConvKernel.size() == 0 || ResetResponse ) {
83 
84  // This is the first response function.
85  // Just calculate the fourier transform.
86 
87  fConvKernel.resize(nticks/2 + 1);
89  }
90  else {
91 
92  // Not the first response function.
93  // Calculate the fourier transform of new response function.
94 
95  std::vector<TComplex> kern(nticks/2 + 1);
96  fft->DoFFT(fResponse, kern);
97 
98  // Update overall convolution kernel.
99 
100  if (kern.size() != fConvKernel.size()) {
101  throw cet::exception("SignalShaping") << __func__ << ": inconsistent kernel size, "
102  << kern.size() << " vs. " << fConvKernel.size();
103  }
104  for(unsigned int i=0; i<kern.size(); ++i)
105  fConvKernel[i] *= kern[i];
106 
107  // Recalculate overall response function.
108 
110  }
111 }
std::vector< TComplex > fConvKernel
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
int FFTSize() const
Definition: LArFFT.h:69
std::vector< double > fResponse
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
void util::SignalShaping::CalculateDeconvKernel ( ) const

Definition at line 251 of file SignalShaping.cxx.

References util::LArFFT::DoInvFFT(), fConvKernel, fDeconvKernel, fDeconvKernelPolarity, fFilter, fFilterLocked, util::LArFFT::FFTSize(), fNorm, fResponse, LockResponse(), and n.

Referenced by Deconvolute(), and set_normflag().

252 {
253  // Make sure configuration is not locked.
254 
255  if(fFilterLocked)
256  throw cet::exception("SignalShaping") << "Configuration locked.\n";
257 
258  // Lock response configuration.
259 
260  LockResponse();
261 
262  // Get FFT service.
263 
265 
266  // Make sure filter function has been configured.
267 
268  if(fFilter.size() == 0)
269  throw cet::exception("SignalShaping")
270  << "Filter function has not been configured.\n";
271 
272  // Make sure filter function has the correct size.
273  // (Should always be the case if we get here.)
274 
275  unsigned int n = fft->FFTSize();
276  if (2 * (fFilter.size() - 1) != n)
277  if (fFilter.size() != fConvKernel.size()) {
278  throw cet::exception("SignalShaping") << __func__ << ": inconsistent size, "
279  << fFilter.size() << " vs. " << fConvKernel.size() << "\n";
280  }
281 
282  // Calculate deconvolution kernel as the ratio of the
283  // filter function and the convolution kernel.
284 
286  for(unsigned int i=0; i<fDeconvKernel.size(); ++i) {
287  if(std::abs(fConvKernel[i].Re()) <= 0.0001 && std::abs(fConvKernel[i].Im()) <= 0.0001) {
288  fDeconvKernel[i] = 0.;
289  }
290  else {
291  fDeconvKernel[i] /= fConvKernel[i];
292  }
293  }
294 
295  // Normalize the deconvolution kernel.
296 
297  // Calculate the unnormalized deconvoluted response
298  // (inverse FFT of filter function).
299 
300  std::vector<double> deconv(n, 0.);
301  fft->DoInvFFT(const_cast<std::vector<TComplex>&>(fFilter), deconv);
302 
303  if (fNorm){
304  // Find the peak value of the response
305  // Should normally be at zero, but don't assume that.
306  // Use DeconvKernelPolairty to find what to normalize to
307  double peak_response = 0;
308  if ( fDeconvKernelPolarity == -1 )
309  peak_response = 4096;
310  for(unsigned int i = 0; i < fResponse.size(); ++i) {
311  if( (fResponse[i] > peak_response)
312  and (fDeconvKernelPolarity == 1))
313  peak_response = fResponse[i];
314  else if ( (fResponse[i] < peak_response)
315  and ( fDeconvKernelPolarity == -1) )
316  peak_response = fResponse[i];
317  }
318  if ( fDeconvKernelPolarity == -1 )
319  peak_response *= -1;
320  if (peak_response <= 0.) {
321  throw cet::exception("SignalShaping") << __func__
322  << ": peak should always be positive (got " << peak_response << ")\n";
323  }
324 
325  // Find the peak value of the deconvoluted response
326  // Should normally be at zero, but don't assume that.
327 
328  double peak_deconv = 0.;
329  for(unsigned int i = 0; i < deconv.size(); ++i) {
330  if(deconv[i] > peak_deconv)
331  peak_deconv = deconv[i];
332  }
333  if (peak_deconv <= 0.) {
334  throw cet::exception("SignalShaping") << __func__
335  << ": deconvolution peak should always be positive (got " << peak_deconv << ")\n";
336  }
337 
338  // Multiply the deconvolution kernel by a factor such that
339  // (Peak of response) = (Peak of deconvoluted response).
340 
341  double ratio = peak_response / peak_deconv;
342  for(unsigned int i = 0; i < fDeconvKernel.size(); ++i)
343  fDeconvKernel[i] *= ratio;
344  }
345  // Set the lock flag.
346 
347  fFilterLocked = true;
348 }
std::vector< TComplex > fConvKernel
void LockResponse() const
std::vector< TComplex > fFilter
void DoInvFFT(std::vector< TComplex > &input, std::vector< T > &output)
Definition: LArFFT.h:119
int FFTSize() const
Definition: LArFFT.h:69
std::vector< TComplex > fDeconvKernel
std::vector< double > fResponse
Char_t n[5]
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
const std::vector<TComplex>& util::SignalShaping::ConvKernel ( ) const
inline

Definition at line 79 of file SignalShaping.h.

References fConvKernel.

79 {return fConvKernel;}
std::vector< TComplex > fConvKernel
template<class T >
void util::SignalShaping::Convolute ( std::vector< T > &  func) const
inline

Definition at line 167 of file SignalShaping.h.

References util::LArFFT::Convolute(), fConvKernel, util::LArFFT::FFTSize(), fResponseLocked, LockResponse(), and n.

Referenced by DeconvKernel().

168 {
169  // Make sure response configuration is locked.
170 
171  if(!fResponseLocked)
172  LockResponse();
173 
174  // Get FFT service.
175 
177 
178  // Make sure that time series has the correct size.
179 
180  int n = func.size();
181  if(n != fft->FFTSize())
182  throw cet::exception("SignalShaping") << "Bad time series size = " << n << "\n";
183 
184  // Do convolution.
185 
186  fft->Convolute(func, const_cast<std::vector<TComplex>&>(fConvKernel));
187 }
std::vector< TComplex > fConvKernel
void LockResponse() const
int FFTSize() const
Definition: LArFFT.h:69
void Convolute(std::vector< T > &input, std::vector< T > &respFunc)
Definition: LArFFT.h:172
Char_t n[5]
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
const std::vector<TComplex>& util::SignalShaping::DeconvKernel ( ) const
inline

Definition at line 81 of file SignalShaping.h.

References Convolute(), Deconvolute(), fDeconvKernel, and Reset().

81 {return fDeconvKernel;}
std::vector< TComplex > fDeconvKernel
template<class T >
void util::SignalShaping::Deconvolute ( std::vector< T > &  func) const
inline

Definition at line 191 of file SignalShaping.h.

References CalculateDeconvKernel(), util::LArFFT::Convolute(), fDeconvKernel, fFilterLocked, util::LArFFT::FFTSize(), and n.

Referenced by DeconvKernel().

192 {
193  // Make sure deconvolution kernel is configured.
194 
195  if(!fFilterLocked)
197 
198  // Get FFT service.
199 
201 
202  // Make sure that time series has the correct size.
203 
204  int n = func.size();
205  if(n != fft->FFTSize())
206  throw cet::exception("SignalShaping") << "Bad time series size = " << n << "\n";
207 
208  // Do convolution.
209 
210  fft->Convolute(func, fDeconvKernel);
211 }
int FFTSize() const
Definition: LArFFT.h:69
void CalculateDeconvKernel() const
void Convolute(std::vector< T > &input, std::vector< T > &respFunc)
Definition: LArFFT.h:172
std::vector< TComplex > fDeconvKernel
Char_t n[5]
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
const std::vector<TComplex>& util::SignalShaping::Filter ( ) const
inline

Definition at line 80 of file SignalShaping.h.

References fFilter.

80 {return fFilter;}
std::vector< TComplex > fFilter
void util::SignalShaping::LockResponse ( ) const

Definition at line 214 of file SignalShaping.cxx.

References fConvKernel, util::LArFFT::FFTSize(), fResponse, fResponseLocked, and n.

Referenced by CalculateDeconvKernel(), Convolute(), and set_normflag().

215 {
216  // Do nothing if the response is already locked.
217 
218  if(!fResponseLocked) {
219 
220  // Get FFT service.
221 
223 
224  // Make sure response has been configured.
225 
226  if(fResponse.size() == 0)
227  throw cet::exception("SignalShaping")
228  << "Response has not been configured.\n";
229 
230  // Make sure response and convolution kernel have the correct
231  // size (should always be the case if we get here).
232 
233  unsigned int n = fft->FFTSize();
234  if (fResponse.size() != n)
235  throw cet::exception("SignalShaping") << __func__ << ": inconsistent kernel size, "
236  << fResponse.size() << " vs. " << n << "\n";
237  if (2 * (fConvKernel.size() - 1) != n)
238  throw cet::exception("SignalShaping") << __func__ << ": unexpected FFT size, "
239  << n << " vs. expected " << (2 * (fConvKernel.size() - 1)) << "\n";
240 
241  // Set the lock flag.
242 
243  fResponseLocked = true;
244  }
245 }
std::vector< TComplex > fConvKernel
int FFTSize() const
Definition: LArFFT.h:69
std::vector< double > fResponse
Char_t n[5]
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
void util::SignalShaping::Reset ( )

Definition at line 47 of file SignalShaping.cxx.

References fConvKernel, fDeconvKernel, fDeconvKernelPolarity, fFilter, fFilterLocked, fResponse, and fResponseLocked.

Referenced by DeconvKernel().

48 {
49  fResponseLocked = false;
50  fFilterLocked = false;
51  fResponse.clear();
52  fConvKernel.clear();
53  fFilter.clear();
54  fDeconvKernel.clear();
55  //Set deconvolution polarity to + as default
57 }
std::vector< TComplex > fConvKernel
std::vector< TComplex > fFilter
std::vector< TComplex > fDeconvKernel
std::vector< double > fResponse
const std::vector<double>& util::SignalShaping::Response ( ) const
inline

Definition at line 77 of file SignalShaping.h.

References fResponse.

77 {return fResponse;}
std::vector< double > fResponse
const std::vector<double>& util::SignalShaping::Response_save ( ) const
inline

Definition at line 78 of file SignalShaping.h.

References fResponse_save.

78 {return fResponse_save;}
std::vector< double > fResponse_save
void util::SignalShaping::save_response ( )
inline

Definition at line 100 of file SignalShaping.h.

References fResponse, and fResponse_save.

std::vector< double > fResponse
std::vector< double > fResponse_save
void util::SignalShaping::set_normflag ( bool  flag)
inline
void util::SignalShaping::SetDeconvKernelPolarity ( int  pol)

Definition at line 197 of file SignalShaping.cxx.

References fDeconvKernelPolarity.

Referenced by set_normflag().

198 {
199 
200  if ( (pol != 1) and (pol != -1) ) {
201  throw cet::exception("SignalShaping") << __func__
202  << ": DeconvKernelPolarity should be +1 or -1 (got " << pol << "). Setting to +1\n";
204  }
205 
206  else
207  fDeconvKernelPolarity = pol;
208 
209 }
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
void util::SignalShaping::SetPeakResponseTime ( double  tick)

Definition at line 140 of file SignalShaping.cxx.

References util::LArFFT::FFTSize(), fResponse, fResponseLocked, util::LArFFT::PeakCorrelation(), and ShiftResponseTime().

Referenced by set_normflag().

141 {
142  // Make sure configuration is not locked.
143 
144  if(fResponseLocked)
145  throw cet::exception("SignalShaping") << "Configuration locked.\n";
146 
147  // Get FFT service.
148 
150 
151  // Construct a delta-function response centered at tick zero.
152 
153  std::vector<double> delta(fft->FFTSize(), 0.);
154  delta[0] = 1.;
155 
156  // Figure out peak of current overall response.
157 
158  double peak = fft->PeakCorrelation(delta, fResponse);
159 
160  // Shift peak response to desired tick.
161 
162  ShiftResponseTime(tick - peak);
163 }
T PeakCorrelation(std::vector< T > &shape1, std::vector< T > &shape2)
Definition: LArFFT.h:271
int FFTSize() const
Definition: LArFFT.h:69
void ShiftResponseTime(double ticks)
std::vector< double > fResponse
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
void util::SignalShaping::ShiftResponseTime ( double  ticks)

Definition at line 117 of file SignalShaping.cxx.

References util::LArFFT::DoInvFFT(), fConvKernel, fResponse, fResponseLocked, and util::LArFFT::ShiftData().

Referenced by set_normflag(), and SetPeakResponseTime().

118 {
119  // Make sure configuration is not locked.
120 
121  if(fResponseLocked)
122  throw cet::exception("SignalShaping") << "Configuration locked.\n";
123 
124  // Get FFT service.
125 
127 
128  // Update convolution kernel.
129 
130  fft->ShiftData(fConvKernel, ticks);
131 
132  // Recalculate overall response functiion.
133 
135 }
void ShiftData(std::vector< TComplex > &input, double shift)
std::vector< TComplex > fConvKernel
void DoInvFFT(std::vector< TComplex > &input, std::vector< T > &output)
Definition: LArFFT.h:119
std::vector< double > fResponse
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33

Member Data Documentation

std::vector<TComplex> util::SignalShaping::fConvKernel
private
std::vector<TComplex> util::SignalShaping::fDeconvKernel
mutableprivate

Definition at line 152 of file SignalShaping.h.

Referenced by CalculateDeconvKernel(), DeconvKernel(), Deconvolute(), and Reset().

int util::SignalShaping::fDeconvKernelPolarity
private

Definition at line 157 of file SignalShaping.h.

Referenced by CalculateDeconvKernel(), Reset(), and SetDeconvKernelPolarity().

std::vector<TComplex> util::SignalShaping::fFilter
private

Definition at line 149 of file SignalShaping.h.

Referenced by AddFilterFunction(), CalculateDeconvKernel(), Filter(), and Reset().

bool util::SignalShaping::fFilterLocked
mutableprivate

Definition at line 139 of file SignalShaping.h.

Referenced by AddFilterFunction(), CalculateDeconvKernel(), Deconvolute(), and Reset().

bool util::SignalShaping::fNorm
private

Definition at line 160 of file SignalShaping.h.

Referenced by CalculateDeconvKernel(), and set_normflag().

std::vector<double> util::SignalShaping::fResponse
private
std::vector<double> util::SignalShaping::fResponse_save
private

Definition at line 143 of file SignalShaping.h.

Referenced by Response_save(), and save_response().

bool util::SignalShaping::fResponseLocked
mutableprivate

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