#include <OpenMS/FILTERING/SMOOTHING/SavitzkyGolayFilter.h>
This class represents a Savitzky-Golay lowpass-filter. The idea of the Savitzky-Golay filter is to find filtercoefficients that preserve higher moments, which means to approximate the underlying function within the moving window by a polynomial of higher order (typically quadratic or quartic). Therefore we least-squares fit for each data point a polynomial to all points in the window and set
to be the value of that polynomial at position
. This method is superior to adjacent averaging because it tends to preserve features of the data such as peak height and width, which are usually 'washed out' by adjacent averaging.
Because of the linearity of the problem, we can reduce the work computing by fitting in advance, for fictious data consisiting of all zeros except for a singe 1 and then do the fits on the real data just by taking linear combinations. There are a particular sets of filter coefficients which accomplish the process of polynomial least-squares fit inside a moving window. To get the symmetric coefficient-matrix
with
The first (last) rows of
we need to smooth the first (last)
data points of the signal. So we use for the smoothing of the first data point the data point itself and the next
future points. For the second point we take the first datapoint, the data point itself and
of rightward data points... . We compute the Matrix
by solving the underlying least-squares problems with the singular value decomposition. Here we demonstrate the computation of the first row of a coefficient-matrix
for a Savitzky-Golay Filter of order=3 and frameSize=5: The design-matrix for the least-squares fit of a linear combination of 3 basis functions to 5 data points is:
To smooth the first data point we have to create a design-matrix with . Now we have to solve the over-determined set of
linear equations
where represents the fictious data. Therefore we solve the normal equations of the least-squares problem
Now, it is possible to get
with . Because we only need one row of the inverse matrix, it is possible to use LU decomposition with only a single backsubstitution. The vector
represents the wanted coefficients. Note that the solution of a least-squares problem directly from the normal equations is faster than the singular value decomposition but rather susceptible to roundoff error!
Public Member Functions | |
SavitzkyGolayFilter () | |
Constructor. | |
virtual | ~SavitzkyGolayFilter () |
Destructor. | |
template<typename InputPeakIterator, typename OutputPeakContainer> | |
void | filter (InputPeakIterator first, InputPeakIterator last, OutputPeakContainer &smoothed_data_container) |
Applies the convolution with the filter coefficients to an given iterator range. | |
template<typename InputPeakContainer, typename OutputPeakContainer> | |
void | filter (const InputPeakContainer &input_peak_container, OutputPeakContainer &baseline_filtered_container) |
Convolutes the filter coefficients and the input raw data. | |
template<typename InputSpectrumIterator, typename OutputPeakType> | |
void | filterExperiment (InputSpectrumIterator first, InputSpectrumIterator last, MSExperiment< OutputPeakType > &ms_exp_filtered) |
Filters every MSSpectrum in a given iterator range. | |
template<typename InputPeakType, typename OutputPeakType> | |
void | filterExperiment (const MSExperiment< InputPeakType > &ms_exp_raw, MSExperiment< OutputPeakType > &ms_exp_filtered) |
Filters a MSExperiment. | |
Protected Member Functions | |
virtual void | updateMembers_ () |
This method is used to update extra member variables at the end of the setParam() method. | |
Protected Attributes | |
UInt | frame_size_ |
UInt of the filter kernel (number of pre-tabulated coefficients). | |
UInt | order_ |
The order of the smoothing polynomial. |
Constructor.
virtual ~SavitzkyGolayFilter | ( | ) | [inline, virtual] |
Destructor.
void filter | ( | InputPeakIterator | first, | |
InputPeakIterator | last, | |||
OutputPeakContainer & | smoothed_data_container | |||
) | [inline] |
Applies the convolution with the filter coefficients to an given iterator range.
Convolutes the filter and the raw data in the iterator intervall [first,last) and writes the resulting data to the smoothed_data_container.
The resulting peaks in the smoothed_data_container (e.g. of type MSSpectrum<DRawDataPoint<1> >) can be of type DRawDataPoint<1> or any other class derived from DRawDataPoint.
Reimplemented from SmoothFilter.
void filter | ( | const InputPeakContainer & | input_peak_container, | |
OutputPeakContainer & | baseline_filtered_container | |||
) | [inline] |
Convolutes the filter coefficients and the input raw data.
Convolutes the filter and the raw data in the input_peak_container and writes the resulting data to the smoothed_data_container.
The resulting peaks in the smoothed_data_container (e.g. of type MSSpectrum<DRawDataPoint<1> >) can be of type DRawDataPoint<1> or any other class derived from DRawDataPoint.
Reimplemented from SmoothFilter.
void filterExperiment | ( | InputSpectrumIterator | first, | |
InputSpectrumIterator | last, | |||
MSExperiment< OutputPeakType > & | ms_exp_filtered | |||
) | [inline] |
Filters every MSSpectrum in a given iterator range.
Filters the data successive in every scan in the intervall [first,last). The filtered data are stored in a MSExperiment.
You have to copy the ExperimentalSettings of the raw data on your own.
Reimplemented from SmoothFilter.
void filterExperiment | ( | const MSExperiment< InputPeakType > & | ms_exp_raw, | |
MSExperiment< OutputPeakType > & | ms_exp_filtered | |||
) | [inline] |
Filters a MSExperiment.
Filters the data every scan in the MSExperiment. The filtered data are stored in a MSExperiment.
Reimplemented from SmoothFilter.
virtual void updateMembers_ | ( | ) | [protected, virtual] |
This method is used to update extra member variables at the end of the setParam() method.
Also call it at the end of the derived classes' copy constructor and assignment operator.
The default implementation is empty.
Reimplemented from DefaultParamHandler.
UInt frame_size_ [protected] |
UInt of the filter kernel (number of pre-tabulated coefficients).
Generated Tue Apr 1 15:36:44 2008 -- using doxygen 1.5.4 | OpenMS / TOPP 1.1 |