npstat is hosted by Hepforge, IPPP Durham
NPStat  5.10.0
npstat::ConvolutionEngine1D Class Reference

#include <ConvolutionEngine1D.hh>

Public Member Functions

 ConvolutionEngine1D (unsigned dataLen, unsigned optimization=FFTW_ESTIMATE)
 
template<typename Real >
void setFilter (const Real *data, unsigned dataLen, unsigned shift=0)
 
void setFilterImage (const std::complex< double > *image, unsigned imLen)
 
template<typename Real1 , typename Real2 >
void convolveWithFilter (const Real1 *in, Real2 *out, unsigned dataLen, unsigned shift=0)
 
template<typename Real1 , typename Real2 >
void deconvolveFilter (const Real1 *in, Real2 *out, unsigned dataLen, unsigned shift=0)
 
template<typename Real >
void depositFilter (unsigned long id, const Real *data, unsigned dataLen, unsigned shift=0)
 
void depositFilterImage (unsigned long id, const std::complex< double > *image, unsigned imLen)
 
template<typename Real1 , typename Real2 >
void convolveWithDeposit (unsigned long id, const Real1 *in, Real2 *out, unsigned dataLen, unsigned shift=0)
 
template<typename Real1 , typename Real2 >
void deconvolveDeposit (unsigned long id, const Real1 *in, Real2 *out, unsigned dataLen, unsigned shift=0)
 
template<typename Real2 >
void convolveFilterAndDeposit (unsigned long id, Real2 *out, unsigned dataLen)
 
template<typename Real2 >
void deconvolveFilterAndDeposit (unsigned long id, Real2 *out, unsigned dataLen)
 
bool discardFilter (unsigned long id)
 
void clearFilterBank ()
 
unsigned dataLength () const
 
unsigned imageLength () const
 
bool isFilterSet () const
 
unsigned filterBankSize () const
 

Detailed Description

Class for performing convolutions by FFT using FFTW library. Internally, transforms are calculated in double precision.

Constructor & Destructor Documentation

◆ ConvolutionEngine1D()

npstat::ConvolutionEngine1D::ConvolutionEngine1D ( unsigned  dataLen,
unsigned  optimization = FFTW_ESTIMATE 
)
explicit

Constructor takes data length and FFTW optimization flag as parameters. Make sure the data length is sufficient to avoid the convolution wrap-around.

Member Function Documentation

◆ clearFilterBank()

void npstat::ConvolutionEngine1D::clearFilterBank ( )

Discard all filters deposited by "depositFilter" method

◆ convolveFilterAndDeposit()

template<typename Real2 >
void npstat::ConvolutionEngine1D::convolveFilterAndDeposit ( unsigned long  id,
Real2 *  out,
unsigned  dataLen 
)

Convolve the filter established by the "setFilter" method with a previously deposited filter

◆ convolveWithDeposit()

template<typename Real1 , typename Real2 >
void npstat::ConvolutionEngine1D::convolveWithDeposit ( unsigned long  id,
const Real1 *  in,
Real2 *  out,
unsigned  dataLen,
unsigned  shift = 0 
)

Convolve provided data with one of the filters previously established by the depositFilter method. The "id" argument is used to identify a particular filter.

The "shift" argument allows you to rotate the input data buffer so that the element with index "shift" becomes new element with index 0 (essentially, the buffer is rotated left by "shift" elements).

◆ convolveWithFilter()

template<typename Real1 , typename Real2 >
void npstat::ConvolutionEngine1D::convolveWithFilter ( const Real1 *  in,
Real2 *  out,
unsigned  dataLen,
unsigned  shift = 0 
)

Convolve provided data with the filter previously established by setFilter.

The "shift" argument allows you to rotate the input data buffer so that the element with index "shift" becomes new element with index 0 (essentially, the buffer is rotated left by "shift" elements).

◆ dataLength()

unsigned npstat::ConvolutionEngine1D::dataLength ( ) const
inline

Expected data length

◆ deconvolveDeposit()

template<typename Real1 , typename Real2 >
void npstat::ConvolutionEngine1D::deconvolveDeposit ( unsigned long  id,
const Real1 *  in,
Real2 *  out,
unsigned  dataLen,
unsigned  shift = 0 
)

Deconvolve the data with a previously deposited filter

◆ deconvolveFilter()

template<typename Real1 , typename Real2 >
void npstat::ConvolutionEngine1D::deconvolveFilter ( const Real1 *  in,
Real2 *  out,
unsigned  dataLen,
unsigned  shift = 0 
)

Deconvolve the input using the previously established filter

◆ deconvolveFilterAndDeposit()

template<typename Real2 >
void npstat::ConvolutionEngine1D::deconvolveFilterAndDeposit ( unsigned long  id,
Real2 *  out,
unsigned  dataLen 
)

Deconvolve the filter established by the "setFilter" method with a previously deposited filter

◆ depositFilter()

template<typename Real >
void npstat::ConvolutionEngine1D::depositFilter ( unsigned long  id,
const Real *  data,
unsigned  dataLen,
unsigned  shift = 0 
)

Provide a filter for subsequent convolutions.

The "shift" argument allows you to rotate the filter data buffer so that the element with index "shift" becomes new element with index 0 (essentially, the buffer is rotated left by "shift" elements).

"id" is the user id for the filter which can be later reused with the method "convolveWithDeposit".

◆ depositFilterImage()

void npstat::ConvolutionEngine1D::depositFilterImage ( unsigned long  id,
const std::complex< double > *  image,
unsigned  imLen 
)

Provide a filter for subsequent convolutions using its Fourier space representation

"id" is the user id for the filter which can be later reused with the method "convolveWithDeposit".

◆ discardFilter()

bool npstat::ConvolutionEngine1D::discardFilter ( unsigned long  id)

Discard a filter previously deposited by "depositFilter" method. This method returns "true" if the filter with the given id was indeed discarded and "false" if such a filter was not found.

◆ imageLength()

unsigned npstat::ConvolutionEngine1D::imageLength ( ) const
inline

Expected image length

◆ setFilter()

template<typename Real >
void npstat::ConvolutionEngine1D::setFilter ( const Real *  data,
unsigned  dataLen,
unsigned  shift = 0 
)

Provide the filter for subsequent convolutions.

The "shift" argument allows you to rotate the filter data buffer so that the element with index "shift" becomes new element with index 0 (essentially, the buffer is rotated left by "shift" elements).

◆ setFilterImage()

void npstat::ConvolutionEngine1D::setFilterImage ( const std::complex< double > *  image,
unsigned  imLen 
)

Provide the filter for subsequent convolutions using its Fourier space representation


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