00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028 #ifndef OPENMS_FILTERING_BASELINE_TOPHATFILTER_H
00029 #define OPENMS_FILTERING_BASELINE_TOPHATFILTER_H
00030
00031 #include <OpenMS/FILTERING/BASELINE/MorphFilter.h>
00032 #include <OpenMS/KERNEL/MSExperiment.h>
00033 #include <OpenMS/MATH/MISC/MathFunctions.h>
00034
00035 #include <algorithm>
00036
00037 namespace OpenMS
00038 {
00039
00057 class TopHatFilter
00058 : public MorphFilter
00059 {
00060 public:
00061 typedef MorphFilter BaseClass;
00062 using BaseClass::struc_size_;
00063
00065 inline TopHatFilter()
00066 : MorphFilter()
00067 {
00068 }
00069
00071 virtual ~TopHatFilter()
00072 {
00073 }
00074
00087 template <typename InputPeakIterator, typename OutputPeakContainer >
00088 void filter(InputPeakIterator first, InputPeakIterator last, OutputPeakContainer& baseline_filtered_container)
00089 {
00090 typedef typename InputPeakIterator::value_type PeakType;
00091
00092
00093
00094 DoubleReal spacing= ((last-1)->getMZ() - first->getMZ()) / (distance(first,last)-1);
00095 int struc_elem_number_of_points = (int) ceil(struc_size_ / spacing );
00096
00097
00098 if (!Math::isOdd(struc_elem_number_of_points))
00099 {
00100 struc_elem_number_of_points += 1;
00101 }
00102
00103
00104
00105
00106
00107
00108 std::vector<PeakType> erosion_result;
00109
00110 this->template erosion
00111 (first, last, erosion_result, struc_elem_number_of_points);
00112
00113 this->template dilatation
00114 (erosion_result.begin(),erosion_result.end(), baseline_filtered_container, struc_elem_number_of_points);
00115
00116 this->template minusIntensities_
00117 (first,last,baseline_filtered_container);
00118 }
00119
00120
00133 template <typename InputPeakContainer, typename OutputPeakContainer >
00134 void filter(const InputPeakContainer& input_peak_container, OutputPeakContainer& baseline_filtered_container)
00135 {
00136
00137 static_cast<SpectrumSettings&>(baseline_filtered_container) = input_peak_container;
00138
00139 filter(input_peak_container.begin(), input_peak_container.end(), baseline_filtered_container);
00140 }
00141
00142
00154 template <typename InputSpectrumIterator, typename OutputPeakType, typename OutputAllocType >
00155 void filterExperiment(InputSpectrumIterator first, InputSpectrumIterator last, MSExperiment<OutputPeakType, OutputAllocType>& ms_exp_filtered)
00156 {
00157 UInt n = distance(first,last);
00158 ms_exp_filtered.reserve(n);
00159 startProgress(0,n,"filtering baseline of data");
00160
00161 for (UInt i = 0; i < n; ++i)
00162 {
00163 InputSpectrumIterator input_it(first+i);
00164
00165 if ( struc_size_ < fabs((input_it->end()-1)->getMZ()- input_it->begin()->getMZ()))
00166 {
00167 MSSpectrum< OutputPeakType > spectrum;
00168
00169
00170 filter(*input_it,spectrum);
00171 setProgress(i);
00172
00173
00174 if (spectrum.size() > 0)
00175 {
00176
00177 static_cast<SpectrumSettings&>(spectrum) = *input_it;
00178 spectrum.setType(SpectrumSettings::RAWDATA);
00179
00180
00181 spectrum.getPrecursorPeak() = input_it->getPrecursorPeak();
00182 spectrum.setRT(input_it->getRT());
00183 spectrum.setMSLevel(input_it->getMSLevel());
00184 spectrum.getName() = input_it->getName();
00185
00186 ms_exp_filtered.push_back(spectrum);
00187 }
00188 }
00189 }
00190 endProgress();
00191 }
00192
00202 template <typename InputPeakType, typename InputAllocType, typename OutputPeakType, typename OutputAllocType >
00203 void filterExperiment(const MSExperiment< InputPeakType, InputAllocType>& ms_exp_raw, MSExperiment<OutputPeakType, OutputAllocType>& ms_exp_filtered)
00204 {
00205
00206 static_cast<ExperimentalSettings&>(ms_exp_filtered) = ms_exp_raw;
00207
00208 filterExperiment(ms_exp_raw.begin(), ms_exp_raw.end(), ms_exp_filtered);
00209 }
00210 };
00211
00212 }
00213 #endif