35 #ifndef OPENMS_FILTERING_SMOOTHING_SAVITZKYGOLAYFILTER_H
36 #define OPENMS_FILTERING_SMOOTHING_SAVITZKYGOLAYFILTER_H
115 template <
typename PeakType>
133 int mid = (frame_size_ / 2);
141 for (i = 0; i <= mid; ++i)
143 it_forward = (first - i);
146 for (j = 0; j < frame_size_; ++j)
148 help += it_forward->getIntensity() * coeffs_[(i + 1) * frame_size_ - 1 - j];
153 out_it->setPosition(first->getPosition());
154 out_it->setIntensity(std::max(0.0, help));
160 it_help = (last - mid);
161 while (first != it_help)
163 it_forward = (first - mid);
166 for (j = 0; j < frame_size_; ++j)
168 help += it_forward->getIntensity() * coeffs_[mid * frame_size_ + j];
173 out_it->setPosition(first->getPosition());
174 out_it->setIntensity(std::max(0.0, help));
180 for (i = (mid - 1); i >= 0; --i)
182 it_forward = (first - (frame_size_ - i - 1));
185 for (j = 0; j < frame_size_; ++j)
187 help += it_forward->getIntensity() * coeffs_[i * frame_size_ + j];
191 out_it->setPosition(first->getPosition());
192 out_it->setIntensity(std::max(0.0, help));
200 template <
typename PeakType>
207 filter_spectra.push_back(*it);
209 filter(filter_spectra);
210 chromatogram.
clear(
false);
213 chromatogram.push_back(*it);
221 template <
typename PeakType>
226 for (
Size i = 0; i < map.
size(); ++i)
229 setProgress(++progress);
234 setProgress(++progress);
247 virtual void updateMembers_();
252 #endif // OPENMS_FILTERING_SMOOTHING_SAVITZKYGOLAYFILTER_H
Size size() const
Definition: MSExperiment.h:117
The representation of a chromatogram.
Definition: MSChromatogram.h:52
unsigned int UInt
Unsigned integer type.
Definition: Types.h:88
UInt frame_size_
UInt of the filter kernel (number of pre-tabulated coefficients)
Definition: SavitzkyGolayFilter.h:243
void filter(MSSpectrum< PeakType > &spectrum)
Removed the noise from an MSSpectrum containing profile data.
Definition: SavitzkyGolayFilter.h:116
Main OpenMS namespace.
Definition: FeatureDeconvolution.h:47
UInt order_
The order of the smoothing polynomial.
Definition: SavitzkyGolayFilter.h:245
Computes the Savitzky-Golay filter coefficients using QR decomposition.
Definition: SavitzkyGolayFilter.h:101
MSChromatogram< ChromatogramPeakType > & getChromatogram(Size id)
returns a single chromatogram
Definition: MSExperiment.h:796
void filter(MSChromatogram< PeakType > &chromatogram)
Definition: SavitzkyGolayFilter.h:201
In-Memory representation of a mass spectrometry experiment.
Definition: MSExperiment.h:69
void filterExperiment(MSExperiment< PeakType > &map)
Removed the noise from an MSExperiment containing profile data.
Definition: SavitzkyGolayFilter.h:222
Base class for all classes that want to report their progress.
Definition: ProgressLogger.h:55
A base class for all classes handling default parameters.
Definition: DefaultParamHandler.h:92
std::vector< double > coeffs_
Coefficients.
Definition: SavitzkyGolayFilter.h:241
void clear(bool clear_meta_data)
Clears all data and meta data.
Definition: MSChromatogram.h:587
const std::vector< MSChromatogram< ChromatogramPeakType > > & getChromatograms() const
returns the chromatogram list
Definition: MSExperiment.h:788