VertexTools/src/FsmwModeFinder3d.cc

0001 #include "CommonTools/Clustering1D/interface/FsmwClusterizer1D.h"
0002 #include "RecoVertex/VertexTools/interface/FsmwModeFinder3d.h"
0003
0004 #include <cmath>
0005 #include <cassert>
0006
0007 /** Half sample mode in 3d, as a functional class.
0008  */
0009
0010 FsmwModeFinder3d::FsmwModeFinder3d(float fraction, float weightExp, float cutoff, int no_w_a)
0011     : theFraction(fraction), theWeightExponent(weightExp), theCutoff(cutoff), theNoWeightsAbove(no_w_a) {
0012   assert(theFraction > 0. && theFraction < 1.);
0013 }
0014
0015 GlobalPoint FsmwModeFinder3d::operator()(const std::vector<PointAndDistance>& values) const {
0016   typedef Cluster1D<void> SimpleCluster;
0017   std::vector<SimpleCluster> vx, vy, vz;
0018   vx.reserve(values.size() - 1);
0019   vy.reserve(values.size() - 1);
0020   vz.reserve(values.size() - 1);
0021   std::vector<const void*> emptyvec;
0022
0023   for (std::vector<PointAndDistance>::const_iterator i = values.begin(); i != values.end(); ++i) {
0024     float weight = 1.;
0025     if (static_cast<int>(values.size()) < theNoWeightsAbove) {
0026       // compute weights if we have fewer than theNoWeightsAbove
0027       // data points
0028       weight = pow(theCutoff + 10000 * i->second, theWeightExponent);
0029     };
0030
0031     SimpleCluster tmp_x(Measurement1D(i->first.x(), 1.0), emptyvec, weight);
0032     SimpleCluster tmp_y(Measurement1D(i->first.y(), 1.0), emptyvec, weight);
0033     SimpleCluster tmp_z(Measurement1D(i->first.z(), 1.0), emptyvec, weight);
0034     vx.push_back(tmp_x);
0035     vy.push_back(tmp_y);
0036     vz.push_back(tmp_z);
0037   };
0038
0039   FsmwClusterizer1D<void> algo(theFraction);
0040   std::vector<SimpleCluster> cresx = algo(vx).first;
0041   std::vector<SimpleCluster> cresy = algo(vy).first;
0042   std::vector<SimpleCluster> cresz = algo(vz).first;
0043   assert(!cresx.empty() && !cresy.empty() && !cresz.empty());
0044
0045   GlobalPoint ret(
0046       cresx.begin()->position().value(), cresy.begin()->position().value(), cresz.begin()->position().value());
0047   return ret;
0048 }
0049
0050 FsmwModeFinder3d* FsmwModeFinder3d::clone() const { return new FsmwModeFinder3d(*this); }