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

SmoothFilter.h (Maintainer: Eva Lange)

Go to the documentation of this file.
00001 // -*- Mode: C++; tab-width: 2; -*-
00002 // vi: set ts=2:
00003 //
00004 // --------------------------------------------------------------------------
00005 //                   OpenMS Mass Spectrometry Framework
00006 // --------------------------------------------------------------------------
00007 //  Copyright (C) 2003-2008 -- Oliver Kohlbacher, Knut Reinert
00008 //
00009 //  This library is free software; you can redistribute it and/or
00010 //  modify it under the terms of the GNU Lesser General Public
00011 //  License as published by the Free Software Foundation; either
00012 //  version 2.1 of the License, or (at your option) any later version.
00013 //
00014 //  This library is distributed in the hope that it will be useful,
00015 //  but WITHOUT ANY WARRANTY; without even the implied warranty of
00016 //  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00017 //  Lesser General Public License for more details.
00018 //
00019 //  You should have received a copy of the GNU Lesser General Public
00020 //  License along with this library; if not, write to the Free Software
00021 //  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
00022 //
00023 // --------------------------------------------------------------------------
00024 // $Maintainer: Eva Lange $
00025 // --------------------------------------------------------------------------
00026 //
00027 
00028 #ifndef OPENMS_FILTERING_SMOOTHING_SMOOTHFILTER_H
00029 #define OPENMS_FILTERING_SMOOTHING_SMOOTHFILTER_H
00030 
00031 #include <OpenMS/KERNEL/MSExperiment.h>
00032 #include <OpenMS/CONCEPT/ProgressLogger.h>
00033 
00034 
00035 namespace OpenMS
00036 {
00040   class SmoothFilter : public ProgressLogger
00041   {
00042     public:
00043 
00045       inline SmoothFilter()
00046       : coeffs_(0)
00047       {
00048       }
00049 
00051       virtual ~SmoothFilter()
00052       {
00053       }
00054 
00067       template <typename InputPeakIterator, typename OutputPeakContainer  >
00068       void filter(InputPeakIterator first, InputPeakIterator last, OutputPeakContainer& smoothed_data_container)
00069       {
00070         smoothed_data_container.resize(distance(first,last));
00071 
00072         // needed for multiply the signal with the filter coefficients
00073         InputPeakIterator it_back;
00074         typename OutputPeakContainer::iterator out_it = smoothed_data_container.begin();
00075 
00076         int m,i,j;
00077         DoubleReal help;
00078 
00079         int frame_size = coeffs_.size();
00080         // compute the transient on
00081         for (i=0; i<frame_size;++i)
00082         {
00083           it_back=first;
00084           help=0;
00085           m=0;
00086 
00087           for (j=i; j>=0; --j)
00088           {
00089             help+=it_back->getIntensity()*coeffs_[m];
00090             --it_back;
00091             ++m;
00092           }
00093 
00094           out_it->setPosition(first->getPosition());
00095           out_it->setIntensity(help);
00096           ++out_it;
00097           ++first;
00098         }
00099 
00100         // compute the steady state output
00101         while (first!=last)
00102         {
00103           it_back=first;
00104           help=0;
00105 
00106           for (j=0; j<frame_size; ++j)
00107           {
00108             help+=it_back->getIntensity()*coeffs_[j];
00109             --it_back;
00110           }
00111 
00112           out_it->setPosition(first->getPosition());
00113           out_it->setIntensity(help);
00114           ++out_it;
00115           ++first;
00116         }
00117       }
00118 
00119 
00132       template <typename InputPeakContainer, typename OutputPeakContainer >
00133       void filter(const InputPeakContainer& input_peak_container, OutputPeakContainer& smoothed_data_container)
00134       {
00135         // copy the experimental settings
00136         static_cast<SpectrumSettings&>(smoothed_data_container) = input_peak_container;
00137         
00138         filter(input_peak_container.begin(), input_peak_container.end(), smoothed_data_container);
00139       }
00140 
00141 
00153       template <typename InputSpectrumIterator, typename OutputPeakType >
00154       void filterExperiment(InputSpectrumIterator first, InputSpectrumIterator last, MSExperiment<OutputPeakType>& ms_exp_filtered)
00155       {
00156         UInt n = distance(first,last);
00157         ms_exp_filtered.reserve(n);
00158         // pick peaks on each scan
00159         for (UInt i = 0; i < n; ++i)
00160         {
00161           MSSpectrum< OutputPeakType > spectrum;
00162           InputSpectrumIterator input_it(first+i);
00163 
00164           // pick the peaks in scan i
00165           filter(*input_it,spectrum);
00166 
00167           // if any peaks are found copy the spectrum settings
00168           if (spectrum.size() > 0)
00169           {
00170             // copy the spectrum settings
00171             static_cast<SpectrumSettings&>(spectrum) = *input_it;
00172             spectrum.setType(SpectrumSettings::RAWDATA);
00173 
00174             // copy the spectrum information
00175             spectrum.getPrecursorPeak() = input_it->getPrecursorPeak();
00176             spectrum.setRT(input_it->getRT());
00177             spectrum.setMSLevel(input_it->getMSLevel());
00178             spectrum.getName() = input_it->getName();
00179 
00180             ms_exp_filtered.push_back(spectrum);
00181           }
00182         }
00183       }
00184 
00185 
00186 
00196       template <typename InputPeakType, typename OutputPeakType >
00197       void filterExperiment(const MSExperiment< InputPeakType >& ms_exp_raw, MSExperiment<OutputPeakType>& ms_exp_filtered)
00198       {
00199         // copy the experimental settings
00200         static_cast<ExperimentalSettings&>(ms_exp_filtered) = ms_exp_raw;
00201 
00202         filterExperiment(ms_exp_raw.begin(), ms_exp_raw.end(), ms_exp_filtered);
00203       }
00204 
00205     protected:
00207       std::vector<DoubleReal> coeffs_;
00208   };
00209 
00210 }// namespace OpenMS
00211 
00212 
00213 #endif

Generated Tue Apr 1 15:36:38 2008 -- using doxygen 1.5.4 OpenMS / TOPP 1.1