Home  · Classes  · Annotated Classes  · Modules  · Members  · Namespaces  · Related Pages

SavitzkyGolayFilter Class Reference
[Spectrum filters]

#include <OpenMS/FILTERING/SMOOTHING/SavitzkyGolayFilter.h>

Inheritance diagram for SavitzkyGolayFilter:

SmoothFilter DefaultParamHandler ProgressLogger

List of all members.


Detailed Description

Computes the Savitzky-Golay filter coefficients using QR decomposition.

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 $ f_i $ in the window and set $g_i$ to be the value of that polynomial at position $ i $. 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 $ c_n $ which accomplish the process of polynomial least-squares fit inside a moving window. To get the symmetric coefficient-matrix $C \in R^{frameSize \times frameSize}$ with

\[ C= \left[ \begin{array}{cccc} c_{0,0} & c_{0,1} & \cdots & c_{0,frameSize-1} \\ \vdots & & & \vdots \\ c_{frameSize-1,0} & c_{frameSize-1,2} & \ldots & c_{frameSize-1,frameSize-1} \end{array} \right]\]

The first (last) $ \frac{frameSize}{2} $ rows of $ C $ we need to smooth the first (last) $ frameSize $ data points of the signal. So we use for the smoothing of the first data point the data point itself and the next $ frameSize-1 $ future points. For the second point we take the first datapoint, the data point itself and $ frameSize-2 $ of rightward data points... . We compute the Matrix $ C $ 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 $ C $ 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:

\[ A=\left[ \begin{array}{ccc} x_0^0 & x_0^1 & x_0^2 \\ \\ x_1^0 & x_1^1 & x_1^2 \\ \\ x_2^0 & x_2^1 & x_2^2 \\ \\ x_3^0 & x_3^1 & x_3^2 \\ \\ x_4^0 & x_4^1 & x_4^2 \end{array} \right]. \]

To smooth the first data point we have to create a design-matrix with $ x=[0,\ldots, frameSize-1] $. Now we have to solve the over-determined set of $ frameSize $ linear equations

\[ Ac=b \]

where $ b=[1,0,\ldots,0] $ represents the fictious data. Therefore we solve the normal equations of the least-squares problem

\[ A^TAc=A^Tb. \]

Now, it is possible to get

\[ c_n=\sum_{m=0}^8 \{(A^TA)^{-1}\}_{0,m} n^m, \]

with $ 0\le n \le 8$. Because we only need one row of the inverse matrix, it is possible to use LU decomposition with only a single backsubstitution. The vector $c=[c_0,\ldots,c_8] $ 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!

Note:
This filter works only for uniform raw data! A polynom order of 4 is recommended. The bigger the frame size the smoother the signal (the more detail information get lost!). The frame size corresponds to the number of filter coefficients, so the width of the smoothing intervall is given by frame_size*spacing of the raw data.
SavitzkyGolayFilter Parameters are explained on a separate page.

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 & Destructor Documentation

SavitzkyGolayFilter (  ) 

Constructor.

virtual ~SavitzkyGolayFilter (  )  [inline, virtual]

Destructor.


Member Function Documentation

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.

Note:
This method assumes that the InputPeakIterator (e.g. of type MSSpectrum<DRawDataPoint<1> >const_iterator) points to a data point of type DRawDataPoint<1> or any other class derived from DRawDataPoint<1>.

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.

Note:
This method assumes that the elements of the InputPeakContainer (e.g. of type MSSpectrum<DRawDataPoint<1> >) are of type DRawDataPoint<1> or any other class derived from DRawDataPoint<1>.

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.

Note:
The InputSpectrumIterator should point to a MSSpectrum. Elements of the input spectra should be of type DRawDataPoint<1> or any other derived class of DRawDataPoint.

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.

Note:
The InputPeakType as well as the OutputPeakType should be of type DRawDataPoint<1> or any other derived class of DRawDataPoint.

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.


Member Data Documentation

UInt frame_size_ [protected]

UInt of the filter kernel (number of pre-tabulated coefficients).

UInt order_ [protected]

The order of the smoothing polynomial.


The documentation for this class was generated from the following file:
Generated Tue Apr 1 15:36:44 2008 -- using doxygen 1.5.4 OpenMS / TOPP 1.1