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

DelaunayPairFinder.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_ANALYSIS_MAPMATCHING_DELAUNAYPAIRFINDER_H
00029 #define OPENMS_ANALYSIS_MAPMATCHING_DELAUNAYPAIRFINDER_H
00030 
00031 #include <OpenMS/ANALYSIS/MAPMATCHING/BasePairFinder.h>
00032 #include <OpenMS/SYSTEM/StopWatch.h>
00033 #include <OpenMS/ANALYSIS/MAPMATCHING/IndexTuple.h>
00034 
00035 #include <CGAL/Cartesian.h>
00036 #include <CGAL/Point_set_2.h>
00037 
00038 #define V_DelaunayPairFinder(bla) // std::cout << bla << std::endl;
00039 #define V_DelaunayConsenus(bla) // std::cout << bla << std::endl;
00040 
00041 namespace OpenMS
00042 {
00043 
00067   template < typename ConsensusMapT = FeatureMap< Feature >, typename ElementMapT = FeatureMap< > >
00068   class DelaunayPairFinder
00069         : public BasePairFinder<ConsensusMapT>
00070   {
00071   public:
00072 
00076     enum Maps
00077     {
00078       MODEL = 0,
00079       SCENE = 1
00080     };
00081 
00082     typedef BasePairFinder< ConsensusMapT > Base;
00083 
00084     // The base knows it all...
00085     typedef typename Base::QualityType            QualityType;
00086     typedef typename Base::PositionType           PositionType;
00087     typedef typename Base::IntensityType          IntensityType;
00088     typedef typename Base::PointType              PointType;
00089     typedef typename Base::PointMapType           PointMapType;
00090     typedef typename Base::ElementPairType        ElementPairType;
00091 
00092     using Base::param_;
00093     using Base::defaults_;
00094     using Base::element_map_;
00095     using Base::element_pairs_;
00096     using Base::transformation_;
00097 
00099     DelaunayPairFinder()
00100         : Base()
00101     {
00102       //set the name for DefaultParamHandler error messages
00103       Base::setName(getProductName());
00104 
00105       defaults_.setValue("similarity:max_pair_distance:RT",20.0,"Consider element el_1 of a map map_1 and element el_2  el_2 of a different map map_2. Then el_1 is assigned to el_2 if they are nearest neighbors and if all other elements in map_1 have a distance greater max_pair_distance::RT to el_2 as well as all other elements in map_2 have a distance greater max_pair_distance::RT to el_1.",true);
00106       defaults_.setValue("similarity:max_pair_distance:MZ",1.0,"Consider element el_1 of a map map_1 and element el_2  el_2 of a different map map_2. Then el_1 is assigned to el_2 if they are nearest neighbors and if all other elements in map_1 have a distance greater max_pair_distance::MZ to el_2 as well as all other elements in map_2 have a distance greater max_pair_distance::MZ to el_1.",true);
00107       defaults_.setValue("similarity:precision:RT",60.0,"Maximal deviation in RT of two elements in different maps to be considered as a corresponding element.");
00108       defaults_.setValue("similarity:precision:MZ",0.5,"Maximal deviation in m/z of two elements in different maps to be considered as a corresponding element.");
00109       defaults_.setValue("similarity:diff_intercept:RT",1.0,"Factor for RT position used to balance the influence of RT and m/z deviations");
00110       defaults_.setValue("similarity:diff_intercept:MZ",0.01,"Factor for m/z position used to balance the influence of RT and m/z deviations");
00111 
00112       Base::defaultsToParam_();
00113     }
00114 
00116     virtual ~DelaunayPairFinder()
00117     {}
00118 
00120     static BasePairFinder<PointMapType>* create()
00121     {
00122       return new DelaunayPairFinder();
00123     }
00124 
00126     static const String getProductName()
00127     {
00128       return "DelaunayPairFinder";
00129     }
00130 
00133     class Point 
00134     : public CGAL::Point_2< CGAL::Cartesian<double> >
00135     {
00136     public:
00137 
00138       typedef CGAL::Point_2< CGAL::Cartesian<double> > Base;
00139 
00140       inline Point() : Base()
00141       {
00142         element = 0;
00143         key = 0;
00144       }
00145 
00146       inline Point(const Base& cgal_point) : Base(cgal_point)
00147       {
00148         element = 0;
00149         key = 0;
00150       }
00151 
00152       inline Point(Base::RT hx, Base::RT hy, const PointType& f, UInt k=0)
00153           : Base(hx,hy)
00154       {
00155         element = &f;
00156         key = k;
00157       }
00158 
00159       inline Point(Base::RT hx, Base::RT hy) : Base(hx,hy)
00160       {
00161         element = 0;
00162         key = 0;
00163       }
00164 
00165       ~Point()
00166       {}
00167 
00169       Point(const Point& source)
00170           : Point_2(source),
00171           key(source.key)
00172       {
00173         element = source.element;
00174       }
00175 
00177       Point& operator = (const Point& source)
00178       {
00179         if (this==&source)
00180           return *this;
00181 
00182         element = source.element;
00183         key = source.key;
00184         Base::operator=(source);
00185         return *this;
00186       }
00187 
00188       const PointType* element;
00189       UInt key;
00190     };
00191 
00194     class  GeometricTraits 
00195       : public CGAL::Cartesian<double>
00196     {
00197     public:
00198       typedef Point Point_2;
00200       class Construct_center_2
00201       {
00202         typedef Point   Point_2;
00203         typedef CGAL::Cartesian<double>::Circle_2  Circle_2;
00204       public:
00205         typedef Point_2          result_type;
00206         typedef CGAL::Arity_tag< 1 >   Arity;
00207 
00208         Point_2
00209         operator()(const Circle_2& c) const
00210         {
00211           return c.center();
00212         }
00213       };
00214     };
00215 
00216     typedef CGAL::Point_set_2< GeometricTraits, CGAL::Triangulation_data_structure_2< CGAL::Triangulation_vertex_base_2< GeometricTraits > > > Point_set_2;
00217     typedef typename Point_set_2::Vertex_handle Vertex_handle;
00218 
00220     double getDiffIntercept(UInt dim)
00221     {
00222       return diff_intercept_[dim];
00223     }
00224 
00226     void setDiffIntercept(UInt dim, DoubleReal intercept)
00227     {
00228          param_.setValue(String("similarity:diff_intercept:") + RawDataPoint2D::shortDimensionName(dim), intercept);
00229          updateMembers_();   
00230     }
00231 
00233     float getMaxPairDistance(UInt dim)
00234     {
00235       if (dim == RawDataPoint2D::RT) return max_pair_distance_[dim];
00236       else 
00237         {
00238           double diff = diff_intercept_[RawDataPoint2D::RT] / diff_intercept_[RawDataPoint2D::MZ];
00239      
00240           return max_pair_distance_[dim] / diff;      
00241         }
00242     }
00243 
00245     void setMaxPairDistance(UInt dim, Real max_pair_distance)
00246     {
00247       if (dim == RawDataPoint2D::RT) 
00248       {
00249         max_pair_distance_[dim] = max_pair_distance; 
00250       }
00251       else 
00252       {
00253         double diff = diff_intercept_[RawDataPoint2D::RT] / diff_intercept_[RawDataPoint2D::MZ];
00254         max_pair_distance_[dim] = max_pair_distance * diff;     
00255        }
00256       param_.setValue(String("similarity:max_pair_distance:") + RawDataPoint2D::shortDimensionName(dim), max_pair_distance);
00257     }
00258 
00260     float getPrecision(UInt dim)
00261     {
00262       if (dim == RawDataPoint2D::RT) return precision_[dim];
00263       else 
00264         {
00265           double diff = diff_intercept_[RawDataPoint2D::RT] / diff_intercept_[RawDataPoint2D::MZ];
00266      
00267           return precision_[dim] / diff;      
00268         }
00269       
00270     }
00271 
00273     void setPrecision(UInt dim, Real precision)
00274     {
00275       if (dim == RawDataPoint2D::RT) 
00276       {
00277         precision_[dim] = precision; 
00278       }
00279       else 
00280       {
00281         double diff = diff_intercept_[RawDataPoint2D::RT] / diff_intercept_[RawDataPoint2D::MZ];
00282         precision_[dim] = precision * diff;     
00283       }
00284       param_.setValue(String("similarity:precision:") + RawDataPoint2D::shortDimensionName(dim), precision);
00285     }
00286 
00288     void findElementPairs()
00289     {
00290       const PointMapType& reference_map = *(element_map_[MODEL]);
00291       const PointMapType& transformed_map = *(element_map_[SCENE]);
00292 
00293 #define V_findElementPairs(bla) V_DelaunayPairFinder(bla)
00294 
00295       V_findElementPairs("@@@ findElementPairs_()");
00296 
00297       UInt n = reference_map.size();
00298 
00299       // Vector to fill the point set for triangulation
00300       // Penalize a deviation in mz more than in rt: deviation(diff_intercept_[RawDataPoint2D::RT]) ~ deviation(diff_intercept_[RawDataPoint2D::MZ])
00301       std::vector< Point > positions_reference_map;
00302       for (UInt i = 0; i < n; ++i)
00303       {
00304         positions_reference_map.push_back(Point(reference_map[i].getRT(),
00305                                                 reference_map[i].getMZ() / (diff_intercept_[RawDataPoint2D::MZ] / diff_intercept_[RawDataPoint2D::RT]),reference_map[i],i));
00306       }
00307 
00308       // compute the delaunay triangulation
00309       Point_set_2 p_set(positions_reference_map.begin(),positions_reference_map.end());
00310 
00311       V_findElementPairs("Translation rt " << transformation_[RawDataPoint2D::RT].getParam());
00312       V_findElementPairs("Translation mz " << transformation_[RawDataPoint2D::MZ].getParam());
00313 
00314       // Initialize a hash map for the elements of reference_map to avoid that elements of the reference map occur in several element pairs
00315       std::vector< Int > lookup_table(n,-1);
00316       std::vector< std::pair< const PointType*,const PointType*> > all_element_pairs;
00317 
00318       UInt index_act_element_pair = 0;
00319       // take each point in the first data map and search for its neighbours in the second element map (within a given (transformed) range)
00320       for ( UInt fi1 = 0; fi1 < transformed_map.size(); ++fi1 )
00321       {
00322         // compute the transformed iso-rectangle (upper_left,bottom_left,bottom_right,upper_right) for the range query
00323         double rt_pos = transformed_map[fi1].getRT();
00324         double mz_pos = transformed_map[fi1].getMZ();
00325 
00326         V_findElementPairs("Search for two nearest neighbours of " << rt_pos << ' ' << transformed_map[fi1].getMZ() );
00327         transformation_[RawDataPoint2D::RT].apply(rt_pos);
00328         transformation_[RawDataPoint2D::MZ].apply(mz_pos);
00329 
00330         mz_pos /= (diff_intercept_[RawDataPoint2D::MZ] / diff_intercept_[RawDataPoint2D::RT]);
00331         Point transformed_pos(rt_pos,mz_pos);
00332 
00333         V_findElementPairs("Transformed Position is : " << transformed_pos );
00334 
00335         std::vector< Vertex_handle > resulting_range;
00336         p_set.nearest_neighbors(transformed_pos,2,std::back_inserter(resulting_range));
00337 
00338         V_findElementPairs("Neighbouring points : ");
00339         for (typename std::vector< Vertex_handle >::const_iterator it = resulting_range.begin(); it != resulting_range.end(); it++)
00340         {
00341           V_findElementPairs((*it)->point());
00342           V_findElementPairs(*((*it)->point().element));
00343         }
00344 
00345         // if the first neighbour is close enough to act_pos and the second_nearest neighbour lies far enough from the nearest neighbour
00346         Point nearest = resulting_range[0]->point();
00347         Point second_nearest = resulting_range[1]->point();
00348 
00349         if (((fabs(transformed_pos[RawDataPoint2D::RT] - nearest.hx())  < precision_[RawDataPoint2D::RT])
00350              &&  (fabs(transformed_pos[RawDataPoint2D::MZ] - nearest.hy())  < precision_[RawDataPoint2D::MZ]))
00351             && ((fabs(second_nearest.hx() - nearest.hx())  > max_pair_distance_[RawDataPoint2D::RT])
00352                 || (fabs(second_nearest.hy() - nearest.hy())  > max_pair_distance_[RawDataPoint2D::MZ])))
00353         {
00354           all_element_pairs.push_back(std::pair<const PointType*,const PointType*>(nearest.element,&transformed_map[fi1]));
00355 
00356           Int element_key = resulting_range[0]->point().key;
00357           // if the element already part of a ElementPair the value in the lookup_table becomes -2
00358           if ( lookup_table[element_key] > -1)
00359           {
00360             lookup_table[element_key] = -2;
00361           }
00362           // otherwise if the element is until now no part of a element pair,
00363           // set the value in the lookup_table to the index of the pair in the all_element_pairs vector
00364           else
00365           {
00366             if ( lookup_table[element_key] == -1)
00367             {
00368               lookup_table[element_key] = index_act_element_pair;
00369             }
00370           }
00371           ++index_act_element_pair;
00372         }
00373       }
00374 
00375 //       std::ofstream out("pairs.dat",std::ios::out);
00376       for (UInt i = 0; i < n; ++i)
00377       {
00378         Int pair_key = lookup_table[i];
00379         if ( pair_key > -1 )
00380         {
00381           element_pairs_->push_back(ElementPairType(*(all_element_pairs[pair_key].second),*(all_element_pairs[pair_key].first)));
00382         }
00383       }
00384 #undef V_findElementPairs
00385 
00386     } // findElementPairs_
00389     template < typename ResultMapType >
00390     void computeConsensusMap(const PointMapType& first_map, ResultMapType& second_map)
00391     {
00392 #define V_computeConsensusMap(bla) V_DelaunayConsenus(bla)
00393       V_computeConsensusMap("@@@ computeConsensusMap()");
00394 
00395       // Vector to fill the point set for triangulation
00396       std::vector< Point > positions_reference_map;
00397       UInt n = first_map.size();
00398       for (UInt i = 0; i < n; ++i)
00399       {
00400         positions_reference_map.push_back(Point((double)(first_map[i].getRT()),
00401                                                 (double)(first_map[i].getMZ() / (diff_intercept_[RawDataPoint2D::MZ] / diff_intercept_[RawDataPoint2D::RT])),first_map[i],i));
00402       }
00403       StopWatch timer;
00404       V_computeConsensusMap("Start delaunay triangulation for " << positions_reference_map.size() << " elements");
00405       // compute the delaunay triangulation
00406       timer.start();
00407       Point_set_2 p_set(positions_reference_map.begin(),positions_reference_map.end());
00408       timer.stop();
00409       V_computeConsensusMap("End delaunay triangulation after " << timer.getCPUTime() << "s");
00410 
00411       // Initialize a hash map for the elements of reference_map to avoid that elements of the reference map occur in several element pairs
00412       std::vector< Int > lookup_table(n,-1);
00413       std::vector< std::pair< const PointType*, PointType*> > all_element_pairs;
00414 
00415       UInt trans_single = 0;
00416       UInt ref_single = 0;
00417       UInt pairs = 0;
00418       UInt index_act_element_pair = 0;
00419       // take each point in the first data map and search for its neighbours in the second element map (within a given (transformed) range)
00420       for ( UInt fi1 = 0; fi1 < second_map.size(); ++fi1 )
00421       {
00422         // compute the transformed iso-rectangle (upper_left,bottom_left,bottom_right,upper_right) for the range query
00423         double rt_pos = (double)(second_map[fi1].getRT());
00424         double mz_pos = (double)(second_map[fi1].getMZ() / (diff_intercept_[RawDataPoint2D::MZ]/diff_intercept_[RawDataPoint2D::RT]));
00425 
00426         V_computeConsensusMap("Search for two nearest neighbours of " << rt_pos << ' ' << second_map[fi1].getMZ() << ' ' << second_map[fi1].getIntensity() );
00427         Point transformed_pos(rt_pos,mz_pos,second_map[fi1]);
00428 
00429         V_computeConsensusMap("Transformed Position is : " << transformed_pos );
00430         std::vector< Vertex_handle > resulting_range;
00431         p_set.nearest_neighbors(transformed_pos,2,std::back_inserter(resulting_range));
00432 
00433         Point nearest;
00434         Point second_nearest;
00435         if (resulting_range.size() == 1)
00436         {
00437           nearest = resulting_range[0]->point();
00438           if ((fabs(transformed_pos[RawDataPoint2D::RT] - nearest.hx())  < precision_[RawDataPoint2D::RT])
00439               &&  (fabs(transformed_pos[RawDataPoint2D::MZ] - nearest.hy())  < precision_[RawDataPoint2D::MZ]))
00440           {
00441             all_element_pairs.push_back(std::pair<const PointType*,PointType*>(nearest.element,&(second_map[fi1])));
00442           }
00443         }
00444         else
00445           if (resulting_range.size() > 1)
00446           {
00447             nearest = resulting_range[0]->point();
00448             second_nearest = resulting_range[1]->point();
00449             V_computeConsensusMap("Nearest: " << (nearest.element)->getRT() << ' ' << fabs(transformed_pos[RawDataPoint2D::RT] - nearest.hx()) << ' ' 
00450                 << (nearest.element)->getMZ() << ' ' << fabs(transformed_pos[RawDataPoint2D::MZ] - nearest.hy()) << ' '
00451                                               << (nearest.element)->getIntensity());
00452                 
00453             V_computeConsensusMap("Second nearest: " << (second_nearest.element)->getRT() << ' ' << second_nearest.hx() << ' ' 
00454                                                      << (second_nearest.element)->getMZ() << ' ' << second_nearest.hy() << ' '
00455                                                      << (second_nearest.element)->getIntensity());
00456 
00457             if (((fabs(transformed_pos[RawDataPoint2D::RT] - nearest.hx())  < precision_[RawDataPoint2D::RT])
00458                  &&  (fabs(transformed_pos[RawDataPoint2D::MZ] - nearest.hy())  < precision_[RawDataPoint2D::MZ]))
00459                 && ((fabs(second_nearest.hx() - nearest.hx())  > max_pair_distance_[RawDataPoint2D::RT])
00460                     || (fabs(second_nearest.hy() - nearest.hy())  > max_pair_distance_[RawDataPoint2D::MZ])))
00461             {
00462               all_element_pairs.push_back(std::pair<const PointType*,PointType*>(nearest.element,&(second_map[fi1])));
00463               V_computeConsensusMap("Push first: " << *(nearest.element))
00464               V_computeConsensusMap("Push second: " << second_map[fi1])
00465 
00466               Int element_key = resulting_range[0]->point().key;
00467 
00468               // if the element a is already part of a ElementPair (a,b) do:
00469               //    if (the element c closer to a than b to a) and (the distance between c and b is > a given threshold) do:
00470               //    --> push (a,c)
00471               //    else
00472               //    --> the value in the lookup_table becomes -2 because the mapping is not unique
00473               if ( lookup_table[element_key] > -1)
00474               {
00475                 Int pair_key = lookup_table[element_key];
00476                 const PointType& first_map_a = *(all_element_pairs[pair_key].first);
00477                 PointType& second_map_b = *(all_element_pairs[pair_key].second);
00478                 PointType& second_map_c = second_map[fi1];
00479 
00480                 V_computeConsensusMap("The element " << first_map_a.getPosition() << ' ' << first_map_a.getIntensity() << " has two element partners \n");
00481                 V_computeConsensusMap(second_map_b.getPosition() << ' ' << second_map_b.getIntensity() 
00482                     << "  and  " << second_map_c.getPosition() << ' ' << second_map_c.getIntensity()) ;
00483 
00484                 V_computeConsensusMap("Range " << second_map_b.getPositionRange() << "  and  " << second_map_c.getPositionRange());
00485 
00486                 if (second_map_c.getPositionRange().encloses(first_map_a.getPosition())
00487                     && !second_map_b.getPositionRange().encloses(first_map_a.getPosition()))
00488                 {
00489                   lookup_table[element_key] = index_act_element_pair;
00490                   V_computeConsensusMap(second_map_c.getPosition() << " and " << first_map_a.getPosition() << " are a pair");
00491                 }
00492                 else
00493                 {
00494                   // if second_map_b and second_map_c do not enclose first_map_a
00495                   if (!(second_map_b.getPositionRange().encloses(first_map_a.getPosition())
00496                         && !second_map_c.getPositionRange().encloses(first_map_a.getPosition())))
00497                   {
00498                     V_computeConsensusMap(second_map_b.getPosition() << " and " << first_map_a.getPosition() << " are a pair, but check the distance between c and b");
00499                     // check the distance between second_map_b and second_map_c
00500                     if (fabs(second_map_b.getMZ() / (diff_intercept_[RawDataPoint2D::MZ]/diff_intercept_[RawDataPoint2D::RT])
00501                              - second_map_c.getMZ() / (diff_intercept_[RawDataPoint2D::MZ]/diff_intercept_[RawDataPoint2D::RT]))
00502                         > max_pair_distance_[RawDataPoint2D::MZ])
00503                     {
00504                       V_computeConsensusMap("distance ok");
00505                       // and check which one of the elements lies closer to first_map_
00506                       if( sqrt(pow((first_map_a.getRT() - second_map_b.getRT()), 2)
00507                                + pow((first_map_a.getMZ() - second_map_b.getMZ()), 2))
00508                           > sqrt(pow((first_map_a.getRT() - second_map_c.getRT()), 2)
00509                                  + pow((first_map_a.getMZ() - second_map_c.getMZ()), 2)))
00510                       {
00511                         lookup_table[element_key] = index_act_element_pair;
00512                         V_computeConsensusMap("take a and c");
00513                       }
00514                     }
00515                     else
00516                     {
00517                       lookup_table[element_key] = -2;
00518                       ++trans_single;
00519                     }
00520                   }
00521                 }
00522               }
00523               // otherwise if the element is until now no part of a element pair,
00524               // set the value in the lookup_table to the index of the pair in the all_element_pairs vector
00525               else
00526               {
00527                 if ( lookup_table[element_key] == -1)
00528                 {
00529                   lookup_table[element_key] = index_act_element_pair;
00530                 }
00531               }
00532               ++index_act_element_pair;
00533             }
00534             // no corresponding element in reference map
00535             // add a singleton consensus element
00536             else
00537             {
00538               ++trans_single;
00539             }
00540           }
00541       }
00542       V_computeConsensusMap("Insert elements ");
00543       std::vector< const PointType* > single_elements_first_map;
00544       for (UInt i = 0; i < n; ++i)
00545       {
00546         Int pair_key = lookup_table[i];
00547         if ( pair_key > -1 )
00548         {
00549           IndexTuple< ElementMapT > index_tuple(*((all_element_pairs[pair_key].first)->begin()));
00550           V_computeConsensusMap("First: " << *((all_element_pairs[pair_key].first)))
00551           V_computeConsensusMap("Second: " << *((all_element_pairs[pair_key].second)))
00552           (all_element_pairs[pair_key].second)->insert(index_tuple);
00553           ++pairs;
00554         }
00555         // add a singleton consensus element
00556         else
00557         {
00558           single_elements_first_map.push_back(positions_reference_map[i].element);
00559           ++ref_single;
00560         }
00561       }
00562 
00563       UInt length = single_elements_first_map.size();
00564       for (UInt i = 0; i < length; ++i)
00565       {
00566         second_map.push_back(*(single_elements_first_map[i]));
00567       }
00568 
00569       V_computeConsensusMap("SINGLE TRANS: " << trans_single);
00570       V_computeConsensusMap("SINGLE REF: " << ref_single);
00571       V_computeConsensusMap("PAIRS: " << pairs);
00572 
00573 #undef V_computeConsensusMap
00574 
00575     } // computeConsensusMap
00576 
00577 
00578   protected:
00579     virtual void updateMembers_()
00580     {
00581       diff_intercept_[RawDataPoint2D::RT] = (double)param_.getValue("similarity:diff_intercept:RT");
00582       diff_intercept_[RawDataPoint2D::MZ] = (double)param_.getValue("similarity:diff_intercept:MZ");
00583       
00584       double diff = diff_intercept_[RawDataPoint2D::RT] / diff_intercept_[RawDataPoint2D::MZ];
00585       
00586       max_pair_distance_[RawDataPoint2D::RT] = (float)param_.getValue("similarity:max_pair_distance:RT");
00587       max_pair_distance_[RawDataPoint2D::MZ] = (float)param_.getValue("similarity:max_pair_distance:MZ") * diff;
00588 
00589       precision_[RawDataPoint2D::RT] = (float)param_.getValue("similarity:precision:RT");
00590       precision_[RawDataPoint2D::MZ] = (float)param_.getValue("similarity:precision:MZ") * diff;
00591 
00592     }
00593 
00595     double diff_intercept_[2];
00599     float max_pair_distance_[2];
00601     float precision_[2];
00602   }
00603   ; // DelaunayPairFinder
00604 } // namespace OpenMS
00605 
00606 #endif  // OPENMS_ANALYSIS_MAPMATCHING_DELAUNAYPAIRFINDER_H

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