#include <SignalShaper.h>

Public Member Functions
	SignalShaper (int fftsize, std::string fftopt)

virtual	~SignalShaper ()

const std::vector< double > &	Response () const

const std::vector< double > &	Response_save () const

const std::vector < std::complex< double > > &	ConvKernel () const

const std::vector < std::complex< double > > &	Filter () const

const std::vector < std::complex< double > > &	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< std::complex< double >> &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< std::complex < double > >	fConvKernel

std::vector< std::complex < double > >	fFilter

std::vector< std::complex < double > >	fDeconvKernel

int	fDeconvKernelPolarity

bool	fNorm

int	fFFTSize

const void *	fPlan

const void *	rPlan

std::unique_ptr < util::LArFFTWPlan >	fFFTPlan

std::unique_ptr< util::LArFFTW >	fFFT

Detailed Description

Definition at line 13 of file SignalShaper.h.

Constructor & Destructor Documentation

util::SignalShaper::SignalShaper	(	int	fftsize,
		std::string	fftopt
	)

Definition at line 9 of file SignalShaper.cxx.

   : fResponseLocked(false)
   , fFilterLocked  (false)
   , fNorm (true)
   , fFFTSize (fftsize)
   , fFFTPlan(new util::LArFFTWPlan(fftsize,fftopt))
   , fFFT(new util::LArFFTW(fftsize,fFFTPlan->fPlan,fFFTPlan->rPlan,0))
 {}

util::SignalShaper::~SignalShaper ( )

virtual

Definition at line 22 of file SignalShaper.cxx.

23 {}

Member Function Documentation

void util::SignalShaper::AddFilterFunction ( const std::vector< std::complex< double >> & filt )

Definition at line 138 of file SignalShaper.cxx.

 {
   // Make sure configuration is not locked.
 
   if(fFilterLocked)
     throw cet::exception("SignalShaper") << "Configuration locked.\n";
 
   // If this is the first filter function, just copy the filter function.
   // Otherwise, update the overall filter function.
 
   if(fFilter.size() == 0) {
     fFilter = filt;
     fFilter.resize(fFFTSize / 2 + 1);
   }
   else {
     unsigned int n = std::min(fFilter.size(), filt.size());
     for(unsigned int i=0; i<n; ++i)
       fFilter[i] *= filt[i];
     for(unsigned int i=n; i<fFilter.size(); ++i)
       fFilter[i] = 0.;
   }
 }

void util::SignalShaper::AddResponseFunction	(	const std::vector< double > &	resp,
		bool	ResetResponse = `false`
	)

Definition at line 43 of file SignalShaper.cxx.

 {
   // Make sure configuration is not locked.
 
   if(fResponseLocked)
     throw cet::exception("SignalShaper") << "Configuration locked.\n";
 
   int nticks = fFFTSize;
 
   // Copy new response function into fResponse attribute, and pad or
   // truncate to correct size.
 
   fResponse = resp;
   fResponse.resize(nticks, 0.);
 
   // Is this the first response function?
 
   if ( fConvKernel.size() == 0 || ResetResponse ) {
 
     // This is the first response function.
     // Just calculate the fourier transform.
 
     fConvKernel.resize(nticks/2 + 1);
     fFFT->DoFFT(fResponse, fConvKernel);
   }
   else {
 
     // Not the first response function.
     // Calculate the fourier transform of new response function.
 
     std::vector<std::complex<double>> kern(nticks/2 + 1);
     fFFT->DoFFT(fResponse, kern);
 
     // Update overall convolution kernel.
 
     if (kern.size() != fConvKernel.size()) {
       throw cet::exception("SignalShaper") << __func__ << ": inconsistent kernel size, "
         << kern.size() << " vs. " << fConvKernel.size();
     }
     for(unsigned int i=0; i<kern.size(); ++i)
       fConvKernel[i] *= kern[i];
 
     // Recalculate overall response function.
 
     fFFT->DoInvFFT(fConvKernel, fResponse);
   }
 }

void util::SignalShaper::CalculateDeconvKernel ( ) const

Definition at line 213 of file SignalShaper.cxx.

 {
   // Make sure configuration is not locked.
 
   if(fFilterLocked)
     throw cet::exception("SignalShaper") << "Configuration locked.\n";
 
   // Lock response configuration.
 
   LockResponse();
 
   // Make sure filter function has been configured.
 
   if(fFilter.size() == 0)
     throw cet::exception("SignalShaper")
       << "Filter function has not been configured.\n";
 
   // Make sure filter function has the correct size.
   // (Should always be the case if we get here.)
 
   unsigned int n = fFFTSize;
   if (2 * (fFilter.size() - 1) != n)
   if (fFilter.size() != fConvKernel.size()) {
     throw cet::exception("SignalShaper") << __func__ << ": inconsistent size, "
       << fFilter.size() << " vs. " << fConvKernel.size() << "\n";
   }
   
   // Calculate deconvolution kernel as the ratio of the 
   // filter function and the convolution kernel.
 
   fDeconvKernel = fFilter;
   for(unsigned int i=0; i<fDeconvKernel.size(); ++i) {
     if(std::abs(fConvKernel[i].real()) <= 0.0001 && std::abs(fConvKernel[i].imag()) <= 0.0001) {
       fDeconvKernel[i] = 0.; 
     }
     else {
       fDeconvKernel[i] /= fConvKernel[i]; 
     }
   }
 
   // Normalize the deconvolution kernel.
 
   // Calculate the unnormalized deconvoluted response
   // (inverse FFT of filter function).
 
   std::vector<double> deconv(n, 0.);
   fFFT->DoInvFFT(const_cast<std::vector<std::complex<double>>&>(fFilter), deconv);
   //fFFT->DoInvFFT(fFilter, deconv);
 
   if (fNorm){
     // Find the peak value of the response
     // Should normally be at zero, but don't assume that.
     // Use DeconvKernelPolairty to find what to normalize to
     double peak_response = 0;
     if ( fDeconvKernelPolarity == -1 )
       peak_response = 4096;
     for(unsigned int i = 0; i < fResponse.size(); ++i) {
       if( (fResponse[i] > peak_response) 
           and (fDeconvKernelPolarity == 1))
         peak_response = fResponse[i];
       else if ( (fResponse[i] < peak_response)
                 and ( fDeconvKernelPolarity == -1) )
         peak_response = fResponse[i];
     }
     if ( fDeconvKernelPolarity == -1 )
       peak_response *= -1;
     if (peak_response <= 0.) {
       throw cet::exception("SignalShaper") << __func__
                                             << ": peak should always be positive (got " << peak_response << ")\n";
     }
     
     // Find the peak value of the deconvoluted response
     // Should normally be at zero, but don't assume that.
     
     double peak_deconv = 0.;
     for(unsigned int i = 0; i < deconv.size(); ++i) {
       if(deconv[i] > peak_deconv)
         peak_deconv = deconv[i];
     }
     if (peak_deconv <= 0.) {
       throw cet::exception("SignalShaper") << __func__
                                             << ": deconvolution peak should always be positive (got " << peak_deconv << ")\n";
     }
     
     // Multiply the deconvolution kernel by a factor such that
     // (Peak of response) = (Peak of deconvoluted response).
     
     double ratio = peak_response / peak_deconv;
     for(unsigned int i = 0; i < fDeconvKernel.size(); ++i)
       fDeconvKernel[i] *= ratio;
   }
 
   // Set the lock flag.
 
   fFilterLocked = true;
 }

const std::vector<std::complex<double> >& util::SignalShaper::ConvKernel ( ) const

inline

Definition at line 23 of file SignalShaper.h.

23 {return fConvKernel;}

util::SignalShaper::fConvKernel

std::vector< std::complex< double > > fConvKernel

Definition: SignalShaper.h:84

template<class T >

void util::SignalShaper::Convolute ( std::vector< T > & func ) const

const std::vector<std::complex<double> >& util::SignalShaper::DeconvKernel ( ) const

inline

Definition at line 25 of file SignalShaper.h.

25 {return fDeconvKernel;}

util::SignalShaper::fDeconvKernel

std::vector< std::complex< double > > fDeconvKernel

Definition: SignalShaper.h:90

template<class T >

void util::SignalShaper::Deconvolute ( std::vector< T > & func ) const

const std::vector<std::complex<double> >& util::SignalShaper::Filter ( ) const

inline

Definition at line 24 of file SignalShaper.h.

24 {return fFilter;}

util::SignalShaper::fFilter

std::vector< std::complex< double > > fFilter

Definition: SignalShaper.h:87

void util::SignalShaper::LockResponse ( ) const

Definition at line 180 of file SignalShaper.cxx.

 {
   // Do nothing if the response is already locked.
 
   if(!fResponseLocked) {
 
     // Make sure response has been configured.
 
     if(fResponse.size() == 0)
       throw cet::exception("SignalShaper")
         << "Response has not been configured.\n";
 
     // Make sure response and convolution kernel have the correct
     // size (should always be the case if we get here).
 
     unsigned int n = fFFTSize;
     if (fResponse.size() != n)
       throw cet::exception("SignalShaper") << __func__ << ": inconsistent kernel size, "
         << fResponse.size() << " vs. " << n << "\n";
     if (2 * (fConvKernel.size() - 1) != n)
       throw cet::exception("SignalShaper") << __func__ << ": unexpected FFT size, "
         << n << " vs. expected " << (2 * (fConvKernel.size() - 1)) << "\n";
 
     // Set the lock flag.
 
     fResponseLocked = true;
   }
 }

void util::SignalShaper::Reset ( )

Definition at line 28 of file SignalShaper.cxx.

 {
   fResponseLocked = false;
   fFilterLocked = false;
   fResponse.clear();
   fConvKernel.clear();
   fFilter.clear();
   fDeconvKernel.clear();
   //Set deconvolution polarity to + as default
   fDeconvKernelPolarity = +1;
 }

const std::vector<double>& util::SignalShaper::Response ( ) const

inline

Definition at line 21 of file SignalShaper.h.

21 {return fResponse;}

util::SignalShaper::fResponse

std::vector< double > fResponse

Definition: SignalShaper.h:80

const std::vector<double>& util::SignalShaper::Response_save ( ) const

inline

Definition at line 22 of file SignalShaper.h.

22 {return fResponse_save;}

util::SignalShaper::fResponse_save

std::vector< double > fResponse_save

Definition: SignalShaper.h:81

void util::SignalShaper::save_response ( )

inline

Definition at line 41 of file SignalShaper.h.

41 { fResponse_save.clear(); fResponse_save=fResponse;}

util::SignalShaper::fResponse_save

std::vector< double > fResponse_save

Definition: SignalShaper.h:81

util::SignalShaper::fResponse

std::vector< double > fResponse

Definition: SignalShaper.h:80

void util::SignalShaper::set_normflag ( bool flag )

inline

Definition at line 42 of file SignalShaper.h.

42 {fNorm = flag;}

flag

then echo unknown compiler flag

Definition: build-dunepdsprce.sh:66

util::SignalShaper::fNorm

bool fNorm

Definition: SignalShaper.h:98

void util::SignalShaper::SetDeconvKernelPolarity ( int pol )

Definition at line 163 of file SignalShaper.cxx.

 {
 
   if ( (pol != 1) and (pol != -1) ) {
     throw cet::exception("SignalShaper") << __func__
       << ": DeconvKernelPolarity should be +1 or -1 (got " << pol << "). Setting to +1\n";
     fDeconvKernelPolarity = +1;
   }
 
   else
     fDeconvKernelPolarity = pol;
 
 }

void util::SignalShaper::SetPeakResponseTime ( double tick )

Definition at line 114 of file SignalShaper.cxx.

 {
   // Make sure configuration is not locked.
 
   if(fResponseLocked)
     throw cet::exception("SignalShaper") << "Configuration locked.\n";
 
   // Construct a delta-function response centered at tick zero.
 
   std::vector<double> delta(fFFTSize, 0.);
   delta[0] = 1.;
 
   // Figure out peak of current overall response.
 
   double peak = fFFT->PeakCorrelation(delta, fResponse);
 
   // Shift peak response to desired tick.
 
   ShiftResponseTime(tick - peak);
 }

void util::SignalShaper::ShiftResponseTime ( double ticks )

Definition at line 95 of file SignalShaper.cxx.

 {
   // Make sure configuration is not locked.
 
   if(fResponseLocked)
     throw cet::exception("SignalShaper") << "Configuration locked.\n";
 
   // Update convolution kernel.
 
   fFFT->ShiftData(fConvKernel, ticks);
 
   // Recalculate overall response functiion.
 
   fFFT->DoInvFFT(fConvKernel, fResponse);
 }

Member Data Documentation

std::vector<std::complex<double> > util::SignalShaper::fConvKernel

private

Definition at line 84 of file SignalShaper.h.

std::vector<std::complex<double> > util::SignalShaper::fDeconvKernel

mutableprivate

Definition at line 90 of file SignalShaper.h.

int util::SignalShaper::fDeconvKernelPolarity

private

Definition at line 95 of file SignalShaper.h.

std::unique_ptr<util::LArFFTW> util::SignalShaper::fFFT

private

Definition at line 104 of file SignalShaper.h.

std::unique_ptr<util::LArFFTWPlan> util::SignalShaper::fFFTPlan

private

Definition at line 103 of file SignalShaper.h.

int util::SignalShaper::fFFTSize

private

Definition at line 100 of file SignalShaper.h.

std::vector<std::complex<double> > util::SignalShaper::fFilter

private

Definition at line 87 of file SignalShaper.h.

bool util::SignalShaper::fFilterLocked

mutableprivate

Definition at line 77 of file SignalShaper.h.

bool util::SignalShaper::fNorm

private

Definition at line 98 of file SignalShaper.h.

const void* util::SignalShaper::fPlan

private

Definition at line 101 of file SignalShaper.h.

std::vector<double> util::SignalShaper::fResponse

private

Definition at line 80 of file SignalShaper.h.

std::vector<double> util::SignalShaper::fResponse_save

private

Definition at line 81 of file SignalShaper.h.

bool util::SignalShaper::fResponseLocked

mutableprivate

Definition at line 76 of file SignalShaper.h.

const void* util::SignalShaper::rPlan

private

Definition at line 102 of file SignalShaper.h.

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

Public Member Functions

Private Attributes

Detailed Description

Constructor & Destructor Documentation

Member Function Documentation

Member Data Documentation