Clustering1D/interface/Cluster1DMerger.h

0001 #ifndef _Cluster1DMerger_H_
0002 #define _Cluster1DMerger_H_
0003
0004 #include "CommonTools/Clustering1D/interface/Cluster1D.h"
0005 #include "CommonTools/Clustering1D/interface/WeightEstimator.h"
0006 #include <cmath>
0007
0008 /**
0009  *  The class that should always be used to merge
0010  *  two Cluster1D into a single Cluster1D.
0011  */
0012
0013 template <class T>
0014 class Cluster1DMerger {
0015 public:
0016   Cluster1DMerger(const WeightEstimator<T>&);
0017   ~Cluster1DMerger();
0018   Cluster1DMerger(const Cluster1DMerger&);
0019   Cluster1D<T> operator()(const Cluster1D<T>& first, const Cluster1D<T>& second) const;
0020
0021 private:
0022   WeightEstimator<T>* theEstimator;
0023 };
0024
0025 /*
0026  *                                implementation
0027  */
0028
0029 template <class T>
0030 Cluster1DMerger<T>::Cluster1DMerger(const WeightEstimator<T>& est) : theEstimator(est.clone()) {}
0031
0032 template <class T>
0033 Cluster1DMerger<T>::~Cluster1DMerger() {
0034   delete theEstimator;
0035 }
0036
0037 template <class T>
0038 Cluster1DMerger<T>::Cluster1DMerger(const Cluster1DMerger& other) : theEstimator(other.theEstimator->clone()) {}
0039
0040 template <class T>
0041 Cluster1D<T> Cluster1DMerger<T>::operator()(const Cluster1D<T>& first, const Cluster1D<T>& second) const {
0042   std::vector<const T*> tracks = first.tracks();
0043   std::vector<const T*> sectracks = second.tracks();
0044   for (typename std::vector<const T*>::const_iterator i = sectracks.begin(); i != sectracks.end(); ++i) {
0045     tracks.push_back(*i);
0046   };
0047   double V1 = first.position().error() * first.position().error();
0048   double V2 = second.position().error() * second.position().error();
0049   double C1 = first.weight() / V1;
0050   double C2 = second.weight() / V2;
0051
0052   double newpos = (first.position().value() * C1 + second.position().value() * C2) / (C1 + C2);
0053
0054   double newerr = sqrt(C1 * C1 * V1 + C2 * C2 * V2) / (C1 + C2);
0055   double newWeight = theEstimator->weight(tracks);
0056
0057   Measurement1D newmeas(newpos, newerr);
0058   return Cluster1D<T>(newmeas, tracks, newWeight);
0059 }
0060
0061 #endif