CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
GapClusterizerInZ.cc
Go to the documentation of this file.
4 
5 using namespace std;
6 
7 namespace {
8 
9  bool recTrackLessZ(const reco::TransientTrack& tk1, const reco::TransientTrack& tk2) {
10  return tk1.stateAtBeamLine().trackStateAtPCA().position().z() <
12  }
13 
14 } // namespace
15 
17  // some defaults to avoid uninitialized variables
18  verbose_ = conf.getUntrackedParameter<bool>("verbose", false);
19  zSep = conf.getParameter<double>("zSeparation");
20  if (verbose_) {
21  std::cout << "TrackClusterizerInZ: algorithm=gap, zSeparation=" << zSep << std::endl;
22  }
23 }
24 
25 float GapClusterizerInZ::zSeparation() const { return zSep; }
26 
27 vector<vector<reco::TransientTrack> > GapClusterizerInZ::clusterize(const vector<reco::TransientTrack>& tracks) const {
28  vector<reco::TransientTrack> tks = tracks; // copy to be sorted
29 
30  vector<vector<reco::TransientTrack> > clusters;
31  if (tks.empty())
32  return clusters;
33 
34  // sort in increasing order of z
35  stable_sort(tks.begin(), tks.end(), recTrackLessZ);
36 
37  // init first cluster
38  vector<reco::TransientTrack>::const_iterator it = tks.begin();
39  vector<reco::TransientTrack> currentCluster;
40  currentCluster.push_back(*it);
41 
42  it++;
43  for (; it != tks.end(); it++) {
44  double zPrev = currentCluster.back().stateAtBeamLine().trackStateAtPCA().position().z();
45  double zCurr = (*it).stateAtBeamLine().trackStateAtPCA().position().z();
46 
47  if (abs(zCurr - zPrev) < zSeparation()) {
48  // close enough ? cluster together
49  currentCluster.push_back(*it);
50  } else {
51  // store current cluster, start new one
52  clusters.push_back(currentCluster);
53  currentCluster.clear();
54  currentCluster.push_back(*it);
55  }
56  }
57 
58  // store last cluster
59  clusters.push_back(currentCluster);
60 
61  return clusters;
62 }
63 
65  desc.add<double>("zSeparation", 1.0);
66  desc.addUntracked<bool>("verbose", false);
67 }
HLT_FULL_cff.zSeparation
zSeparation
Definition: HLT_FULL_cff.py:104297
GapClusterizerInZ::clusterize
std::vector< std::vector< reco::TransientTrack > > clusterize(const std::vector< reco::TransientTrack > &tracks) const override
Definition: GapClusterizerInZ.cc:27
reco::TransientTrack::stateAtBeamLine
TrajectoryStateClosestToBeamLine stateAtBeamLine() const
Definition: TransientTrack.h:119
gather_cfg.cout
cout
Definition: gather_cfg.py:144
edm::ParameterSetDescription
Definition: ParameterSetDescription.h:52
Measurement1D.h
GapClusterizerInZ::GapClusterizerInZ
GapClusterizerInZ(const edm::ParameterSet &conf)
Definition: GapClusterizerInZ.cc:16
edm::ParameterSet::getUntrackedParameter
T getUntrackedParameter(std::string const &, T const &) const
FreeTrajectoryState::position
GlobalPoint position() const
Definition: FreeTrajectoryState.h:67
PV3DBase::z
T z() const
Definition: PV3DBase.h:61
bsc_activity_cfg.clusters
clusters
Definition: bsc_activity_cfg.py:36
edm::ParameterSet
Definition: ParameterSet.h:47
tracks
const uint32_t *__restrict__ const HitContainer *__restrict__ TkSoA *__restrict__ tracks
Definition: CAHitNtupletGeneratorKernelsImpl.h:176
GapClusterizerInZ::zSeparation
float zSeparation() const
Definition: GapClusterizerInZ.cc:25
submitPVResolutionJobs.desc
string desc
Definition: submitPVResolutionJobs.py:251
std
Definition: JetResolutionObject.h:76
reco::TransientTrack
Definition: TransientTrack.h:19
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
VertexException.h
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
GapClusterizerInZ::fillPSetDescription
static void fillPSetDescription(edm::ParameterSetDescription &desc)
Definition: GapClusterizerInZ.cc:64
GapClusterizerInZ.h
TrajectoryStateClosestToBeamLine::trackStateAtPCA
FTS const & trackStateAtPCA() const
Definition: TrajectoryStateClosestToBeamLine.h:32