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

PeakPickerCWT Class Reference
[PeakPicking]

#include <OpenMS/TRANSFORMATIONS/RAW2PEAK/PeakPickerCWT.h>

Inheritance diagram for PeakPickerCWT:

PeakPicker DefaultParamHandler ProgressLogger

List of all members.


Detailed Description

This class implements a peak picking algorithm using wavelet techniques (as described by Lange et al. (2006) Proc. PSB-06).

This peak picking algorithm uses the continuous wavelet transform of a raw data signal to detect mass peaks. Afterwards a given asymmetric peak function is fitted to the raw data and important peak parameters (e.g. fwhm) are extracted. In an optional step these parameters can be optimized using a non-linear opimization method.

PeakPickerCWT Parameters are explained on a separate page.

Public Types

typedef RawDataPoint1D RawDataPointType
 Raw data point type.
typedef DPeakArray
< RawDataPointType
RawDataArrayType
 Raw data container type using for the temporary storage of the input data.
typedef RawDataArrayType::iterator RawDataPointIterator
 Raw data iterator type.
typedef DPosition< 1 > PositionType
 Position type.

Public Member Functions

 PeakPickerCWT ()
 Constructor.
virtual ~PeakPickerCWT ()
 Destructor.
template<typename InputPeakIterator, typename OutputPeakContainer>
void pick (InputPeakIterator first, InputPeakIterator last, OutputPeakContainer &picked_peak_container, int ms_level=1)
 Applies the peak picking algorithm to an given iterator range.
template<typename InputPeakContainer, typename OutputPeakContainer>
void pick (const InputPeakContainer &input_peak_container, OutputPeakContainer &picked_peaks_container, int ms_level=1)
 Applies the peak picking algorithm to a raw data point container.
template<typename InputSpectrumIterator, typename OutputPeakType>
void pickExperiment (InputSpectrumIterator first, InputSpectrumIterator last, MSExperiment< OutputPeakType > &ms_exp_peaks)
 Picks the peaks in a range of MSSpectra.
template<typename InputPeakType, typename OutputPeakType>
void pickExperiment (const MSExperiment< InputPeakType > &ms_exp_raw, MSExperiment< OutputPeakType > &ms_exp_peaks)
 Picks the peaks in a MSExperiment.

Protected Member Functions

void updateMembers_ ()
 This method is used to update extra member variables at the end of the setParam() method.
void init_ ()
 Initializes the members and parses the parameter object.
template<typename OutputPeakType>
void fillPeak_ (const PeakShape &, OutputPeakType &)
 This function fills the members of a picked peak of type OutputPeakType.
void getPeakArea_ (const PeakArea_ &area, double &area_left, double &area_right)
 Computes the peak's left and right area.
PeakShape fitPeakShape_ (const PeakArea_ &area, bool enable_centroid_fit)
 Returns the best fitting peakshape.
double correlate_ (const PeakShape &peak, const PeakArea_ &area, int direction=0) const
 Returns the squared pearson coefficient.
bool getMaxPosition_ (RawDataPointIterator first, RawDataPointIterator last, const ContinuousWaveletTransform &wt, PeakArea_ &area, int distance_from_scan_border, int ms_level, int direction=1)
 Finds the next maximum position in the wavelet transform wt.
bool getPeakEndPoints_ (RawDataPointIterator first, RawDataPointIterator last, PeakArea_ &area, int distance_from_scan_border, int &peak_left_index, int &peak_right_index)
 Determines a peaks's endpoints.
void getPeakCentroid_ (PeakArea_ &area)
 Estimates a peak's centroid position.
double lorentz_ (double height, double lambda, double pos, double x)
 Computes the value of a theroretical lorentz peak at position x.
void initializeWT_ ()
 Computes the threshold for the peak height in the wavelet transform and initializes the wavelet transform.
template<>
void fillPeak_ (const PeakShape &peak_shape, PickedPeak1D &picked_peak)
 Fills the members of a PickedPeak1D given an PeakShape.
Methods needed for separation of overlapping peaks
bool deconvolutePeak_ (PeakShape &shape)
 Separates overlapping peaks.
int getNumberOfPeaks_ (RawDataPointIterator &first, RawDataPointIterator &last, std::vector< double > &peak_values, int direction, DoubleReal resolution, ContinuousWaveletTransformNumIntegration &wt)
 Determines the number of peaks in the given mass range using the cwt.
int determineChargeState_ (std::vector< double > &peak_values)
 Estimate the charge state of the peaks.
void addPeak_ (std::vector< PeakShape > &peaks_DC, PeakArea_ &area, double left_width, double right_width)
 Add a peak.

Protected Attributes

std::vector< PeakShapepeak_shapes_
 Container the determined peak shapes.
ContinuousWaveletTransformNumIntegration wt_
 The continuous wavelet "transformer".
ContinuousWaveletTransformNumIntegration wtDC_
 The continuous wavelet "transformer" for the deconvolution.
UInt radius_
 The search radius for the determination of a peak's maximum position.
float scale_
 The dilation of the wavelet.
float peak_bound_cwt_
 The minimal height which defines a peak in the CWT (MS 1 level).
float peak_bound_ms2_level_cwt_
 The minimal height which defines a peak in the CWT (MS 2 level).
float peak_corr_bound_
 The threshold for correlation.
float noise_level_
 The threshold for the noise level (TODO: Use the information of the signal to noise estimator).
bool optimization_
 Switch for the optimization of peak parameters.
bool deconvolution_
 Switch for the deconvolution of peak parameters.
bool two_d_optimization_
 Switch for the 2D optimization of peak parameters.

Classes

class  PeakArea_
 Class for the internal peak representation. More...


Member Typedef Documentation

typedef RawDataPoint1D RawDataPointType

Raw data point type.

typedef DPeakArray<RawDataPointType > RawDataArrayType

Raw data container type using for the temporary storage of the input data.

typedef RawDataArrayType::iterator RawDataPointIterator

Raw data iterator type.

typedef DPosition<1> PositionType

Position type.


Constructor & Destructor Documentation

PeakPickerCWT (  ) 

Constructor.

virtual ~PeakPickerCWT (  )  [virtual]

Destructor.


Member Function Documentation

void pick ( InputPeakIterator  first,
InputPeakIterator  last,
OutputPeakContainer &  picked_peak_container,
int  ms_level = 1 
) [inline]

Applies the peak picking algorithm to an given iterator range.

Picks the peaks in the given iterator intervall [first,last) and writes the resulting peaks to the picked_peak_container. The ms_level should be one if the spectrum is a normal mass spectrum, or two if it is a tandem mass spectrum.

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

The resulting peaks in the picked_peak_container (e.g. of type MSSpectrum<PickedPeak1D>) can be of type RawDataPoint1D or any other class derived from DRawDataPoint. We recommend to use the PickedPeak1Dbecause it stores important information gained during the peak picking algorithm.

void pick ( const InputPeakContainer &  input_peak_container,
OutputPeakContainer &  picked_peaks_container,
int  ms_level = 1 
) [inline]

Applies the peak picking algorithm to a raw data point container.

Picks the peaks in the input container (e.g. of type MSSpectrum<RawDataPoint1D >) and writes the resulting peaks to the picked_peak_container (e.g. MSSpectrum<PickedPeak1D>).

The ms_level should be one if the spectrum is a normal mass spectrum, or two if it is a tandem mass spectrum.

Note:
This method assumes that the input_peak_container contains data points of type RawDataPoint1D or any other class derived from DRawDataPoint.

The resulting peaks in the picked_peak_container (e.g. of type MSSpectrum<PickedPeak1D>) can be of type RawDataPoint1D or any other class derived from DRawDataPoint. We recommend to use the PickedPeak1Dbecause it stores important information gained during the peak picking algorithm.

void pickExperiment ( InputSpectrumIterator  first,
InputSpectrumIterator  last,
MSExperiment< OutputPeakType > &  ms_exp_peaks 
) [inline]

Picks the peaks in a range of MSSpectra.

Picks the peaks successive in every scan in the intervall [first,last). The detected peaks are stored in a MSExperiment.

Note:
The InputSpectrumIterator should point to a MSSpectrum. Elements of the input spectra should be of type RawDataPoint1D or any other derived class of DRawDataPoint. For the resulting peaks we recommend to use the PickedPeak1Dbecause it stores important information gained during the peak picking algorithm.

You have to copy the ExperimentalSettings of the raw data by your own.

void pickExperiment ( const MSExperiment< InputPeakType > &  ms_exp_raw,
MSExperiment< OutputPeakType > &  ms_exp_peaks 
) [inline]

Picks the peaks in a MSExperiment.

Picks the peaks on every scan in the MSExperiment. The detected peaks are stored in a MSExperiment.

Note:
The input peaks should be of type RawDataPoint1D or any other derived class of DRawDataPoint. For the resulting peaks we recommend to use the PickedPeak1Dbecause it stores important information gained during the peak picking algorithm.

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 PeakPicker.

void init_ (  )  [protected]

Initializes the members and parses the parameter object.

void fillPeak_ ( const PeakShape ,
OutputPeakType &   
) [inline, protected]

This function fills the members of a picked peak of type OutputPeakType.

void getPeakArea_ ( const PeakArea_ area,
double &  area_left,
double &  area_right 
) [protected]

Computes the peak's left and right area.

PeakShape fitPeakShape_ ( const PeakArea_ area,
bool  enable_centroid_fit 
) [protected]

Returns the best fitting peakshape.

double correlate_ ( const PeakShape peak,
const PeakArea_ area,
int  direction = 0 
) const [protected]

Returns the squared pearson coefficient.

Computes the correlation of the peak and the original data given by the peak enpoints area.left and area.right. If the value is near 1, the fitted peakshape and the raw data are expected to be very similar.

bool getMaxPosition_ ( RawDataPointIterator  first,
RawDataPointIterator  last,
const ContinuousWaveletTransform wt,
PeakArea_ area,
int  distance_from_scan_border,
int  ms_level,
int  direction = 1 
) [protected]

Finds the next maximum position in the wavelet transform wt.

If the maximum is greater than peak_bound_cwt we search for the corresponding maximum in the raw data interval [first,last) given a predefined search radius radius. Only peaks with intensities greater than peak_bound_ are relevant. If no peak is detected the method return false. For direction=1, the method runs from first to last given direction=-1 it runs the other way around.

bool getPeakEndPoints_ ( RawDataPointIterator  first,
RawDataPointIterator  last,
PeakArea_ area,
int  distance_from_scan_border,
int &  peak_left_index,
int &  peak_right_index 
) [protected]

Determines a peaks's endpoints.

The algorithm does the following:

(1) starting from x_l', walk left until one of the following happens

(2) analogous procedure to the right of x_r

void getPeakCentroid_ ( PeakArea_ area  )  [protected]

Estimates a peak's centroid position.

Computes the centroid position of the peak using all raw data points which are greater than 60% of the most intensive raw data point.

double lorentz_ ( double  height,
double  lambda,
double  pos,
double  x 
) [protected]

Computes the value of a theroretical lorentz peak at position x.

void initializeWT_ (  )  [protected]

Computes the threshold for the peak height in the wavelet transform and initializes the wavelet transform.

Given the threshold for the peak height a corresponding value peak_bound_cwt can be computed for the continious wavelet transform. Therefore we compute a theoretical lorentzian peakshape with height=peak_bound_ and a width which is similar to the width of the wavelet. Taking the maximum in the wavelet transform of the lorentzian peak we have a peak bound in the wavelet transform.

bool deconvolutePeak_ ( PeakShape shape  )  [protected]

Separates overlapping peaks.

It determines the number of peaks lying underneath the initial peak using the cwt with different scales. Then a nonlinear optimzation procedure is applied to optimize the peak parameters.

int getNumberOfPeaks_ ( RawDataPointIterator first,
RawDataPointIterator last,
std::vector< double > &  peak_values,
int  direction,
DoubleReal  resolution,
ContinuousWaveletTransformNumIntegration wt 
) [protected]

Determines the number of peaks in the given mass range using the cwt.

int determineChargeState_ ( std::vector< double > &  peak_values  )  [protected]

Estimate the charge state of the peaks.

void addPeak_ ( std::vector< PeakShape > &  peaks_DC,
PeakArea_ area,
double  left_width,
double  right_width 
) [protected]

Add a peak.

void fillPeak_ ( const PeakShape peak_shape,
PickedPeak1D picked_peak 
) [inline, protected]

Fills the members of a PickedPeak1D given an PeakShape.


Member Data Documentation

std::vector<PeakShape> peak_shapes_ [protected]

Container the determined peak shapes.

ContinuousWaveletTransformNumIntegration wt_ [protected]

The continuous wavelet "transformer".

ContinuousWaveletTransformNumIntegration wtDC_ [protected]

The continuous wavelet "transformer" for the deconvolution.

UInt radius_ [protected]

The search radius for the determination of a peak's maximum position.

float scale_ [protected]

The dilation of the wavelet.

float peak_bound_cwt_ [protected]

The minimal height which defines a peak in the CWT (MS 1 level).

float peak_bound_ms2_level_cwt_ [protected]

The minimal height which defines a peak in the CWT (MS 2 level).

float peak_corr_bound_ [protected]

The threshold for correlation.

float noise_level_ [protected]

The threshold for the noise level (TODO: Use the information of the signal to noise estimator).

bool optimization_ [protected]

Switch for the optimization of peak parameters.

bool deconvolution_ [protected]

Switch for the deconvolution of peak parameters.

bool two_d_optimization_ [protected]

Switch for the 2D optimization of peak parameters.


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