CMS 3D CMS Logo

FFTJetDijetFilter.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // Package: FFTJetProducers
4 // Class: FFTJetDijetFilter
5 //
13 //
14 // Original Author: Igor Volobouev
15 // Created: Thu Jul 18 19:19:40 CDT 2012
16 //
17 //
18 #include <cmath>
19 #include <cassert>
20 #include <iostream>
21 
22 #include "fftjet/EquidistantSequence.hh"
23 #include "fftjet/ProximityClusteringTree.hh"
24 #include "fftjet/SparseClusteringTree.hh"
25 #include "fftjet/PeakEtaPhiDistance.hh"
26 #include "fftjet/peakEtLifetime.hh"
27 
28 // framework include files
35 
39 
43 
44 #define init_param(type, varname) varname (ps.getParameter< type >( #varname ))
45 
46 using namespace fftjetcms;
47 
48 //
49 // class declaration
50 //
52 {
53 public:
54  typedef fftjet::ProximityClusteringTree<fftjet::Peak,long> ClusteringTree;
55  typedef fftjet::SparseClusteringTree<fftjet::Peak,long> SparseTree;
56 
57  explicit FFTJetDijetFilter(const edm::ParameterSet&);
58  ~FFTJetDijetFilter() override;
59 
60 private:
62 
63  FFTJetDijetFilter() = delete;
64  FFTJetDijetFilter(const FFTJetDijetFilter&) = delete;
65  FFTJetDijetFilter& operator=(const FFTJetDijetFilter&) = delete;
66 
67  void beginJob() override;
68  bool filter(edm::Event& iEvent, const edm::EventSetup& iSetup) override;
69  void endJob() override;
70 
71  template<class Ptr>
72  inline void checkConfig(const Ptr& ptr, const char* message) const
73  {
74  if (ptr.get() == nullptr)
75  throw cms::Exception("FFTJetBadConfig") << message << std::endl;
76  }
77 
78  inline double peakPt(const fftjet::Peak& peak) const
79  {
80  const double s = peak.scale();
81  return ptConversionFactor*s*s*peak.magnitude();
82  }
83 
84  // Module parameters
87  double fixedScale;
90  double minDeltaPhi;
92  double minPt0;
93  double minPt1;
94  double maxPeakEta;
96 
97  // Distance calculator for the clustering tree
98  std::unique_ptr<fftjet::AbsDistanceCalculator<fftjet::Peak> > distanceCalc;
99 
100  // Scales used
101  std::unique_ptr<std::vector<double> > iniScales;
102 
103  // The clustering trees
104  ClusteringTree* clusteringTree;
105  SparseTree* sparseTree;
106 
107  // Space for peaks
108  std::vector<fftjet::Peak> peaks;
109 
110  // Space for sparse tree nodes
111  std::vector<unsigned> nodes;
112 
113  // pass/fail decision counters
114  unsigned long nPassed;
115  unsigned long nRejected;
116 
118 };
119 
120 
121 //
122 // constructors and destructor
123 //
125  : init_param(edm::InputTag, treeLabel),
127  init_param(double, fixedScale),
129  init_param(double, min1to0PtRatio),
130  init_param(double, minDeltaPhi),
132  init_param(double, minPt0),
133  init_param(double, minPt1),
134  init_param(double, maxPeakEta),
136  clusteringTree(nullptr),
137  sparseTree(nullptr)
138 {
139  // Parse the set of scales
141  ps.getParameter<edm::ParameterSet>("InitialScales"));
142  checkConfig(iniScales, "invalid set of scales");
143  std::sort(iniScales->begin(), iniScales->end(), std::greater<double>());
144 
145  // Parse the distance calculator
147  ps.getParameter<edm::ParameterSet>("TreeDistanceCalculator"));
148  distanceCalc = fftjet_DistanceCalculator_parser(TreeDistanceCalculator);
149  checkConfig(distanceCalc, "invalid tree distance calculator");
150 
151  // Create the clustering tree
153  sparseTree = new SparseTree();
154 
155  treeToken = consumes<StoredTree>(treeLabel);
156 }
157 
158 
160 {
161  delete sparseTree;
162  delete clusteringTree;
163 }
164 
165 
167 {
168  nPassed = 0;
169  nRejected = 0;
170 }
171 
172 
174 {
175 // std::cout << "In FTJetDijetFilter::endJob: nPassed = " << nPassed
176 // << ", nRejected = " << nRejected << std::endl;
177 }
178 
179 
180 // ------------ method called to produce the data ------------
182  edm::Event& iEvent, const edm::EventSetup& iSetup)
183 {
185  iEvent.getByToken(treeToken, input);
186 
187  // Convert the stored tree into a normal FFTJet clustering tree
188  // and extract the set of peaks at the requested scale
189  const double eventScale = insertCompleteEvent ? completeEventScale : 0.0;
190  if (input->isSparse())
191  {
193  eventScale, sparseTree);
194  sparseTree->sortNodes();
195  fftjet::updateSplitMergeTimes(*sparseTree, sparseTree->minScale(),
196  sparseTree->maxScale());
197  const unsigned usedLevel = sparseTree->getLevel(fixedScale);
198  sparseTree->getLevelNodes(usedLevel, &nodes);
199  const unsigned numNodes = nodes.size();
200  peaks.clear();
201  peaks.reserve(numNodes);
202  for (unsigned i=0; i<numNodes; ++i)
203  peaks.push_back(sparseTree->uncheckedNode(nodes[i]).getCluster());
204  }
205  else
206  {
207  densePeakTreeFromStorable(*input, iniScales.get(),
208  eventScale, clusteringTree);
209  const unsigned usedLevel = clusteringTree->getLevel(fixedScale);
210  double actualScale = 0.0;
211  long dummyInfo;
212  clusteringTree->getLevelData(usedLevel,&actualScale,&peaks,&dummyInfo);
213  }
214 
215  // Get out if we don't have two clusters
216  const unsigned nClusters = peaks.size();
217  if (nClusters < 2)
218  {
219  ++nRejected;
220  return false;
221  }
222  std::sort(peaks.begin(), peaks.end(), std::greater<fftjet::Peak>());
223 
224  // Calculate all quantities needed to make the pass/fail decision
225  const double pt0 = peakPt(peaks[0]);
226  const double pt1 = peakPt(peaks[1]);
227  const double dphi = reco::deltaPhi(peaks[0].phi(), peaks[1].phi());
228  const double ptratio = pt1/pt0;
229  double thirdJetFraction = 0.0;
230  if (nClusters > 2)
231  thirdJetFraction = peakPt(peaks[2])/(pt0 + pt1);
232 
233  // Calculate the cut
234  const bool pass = pt0 > minPt0 &&
235  pt1 > minPt1 &&
236  ptratio > min1to0PtRatio &&
237  std::abs(dphi) > minDeltaPhi &&
238  thirdJetFraction < maxThirdJetFraction &&
239  std::abs(peaks[0].eta()) < maxPeakEta &&
240  std::abs(peaks[1].eta()) < maxPeakEta;
241  if (pass)
242  ++nPassed;
243  else
244  ++nRejected;
245  return pass;
246 }
247 
248 //define this as a plug-in
constexpr double deltaPhi(double phi1, double phi2)
Definition: deltaPhi.h:22
T getParameter(std::string const &) const
std::vector< fftjet::Peak > peaks
#define init_param(type, varname)
bool isSparse() const
Definition: PattRecoTree.h:28
Class for storing FFTJet sparse clustering trees.
Definition: PattRecoTree.h:20
bool filter(edm::Event &iEvent, const edm::EventSetup &iSetup) override
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:517
~FFTJetDijetFilter() override
std::unique_ptr< std::vector< double > > iniScales
#define nullptr
void beginJob() override
fftjet::ProximityClusteringTree< fftjet::Peak, long > ClusteringTree
void beginJob()
Definition: Breakpoints.cc:14
static std::string const input
Definition: EdmProvDump.cc:48
fftjet::SparseClusteringTree< fftjet::Peak, long > SparseTree
int iEvent
Definition: GenABIO.cc:224
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
std::vector< unsigned > nodes
edm::EDGetTokenT< StoredTree > treeToken
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
std::unique_ptr< fftjet::AbsDistanceCalculator< fftjet::Peak > > fftjet_DistanceCalculator_parser(const edm::ParameterSet &ps)
ClusteringTree * clusteringTree
FFTJetDijetFilter()=delete
reco::PattRecoTree< float, reco::PattRecoPeak< float > > StoredTree
void sparsePeakTreeFromStorable(const reco::PattRecoTree< Real, reco::PattRecoPeak< Real > > &in, const std::vector< double > *scaleSetIfNotAdaptive, double completeEventScale, fftjet::SparseClusteringTree< fftjet::Peak, long > *out)
unsigned long nRejected
void densePeakTreeFromStorable(const reco::PattRecoTree< Real, reco::PattRecoPeak< Real > > &in, const std::vector< double > *scaleSetIfNotAdaptive, double completeEventScale, fftjet::AbsClusteringTree< fftjet::Peak, long > *out)
std::unique_ptr< fftjet::AbsDistanceCalculator< fftjet::Peak > > distanceCalc
HLT enums.
void checkConfig(const Ptr &ptr, const char *message) const
edm::InputTag treeLabel
std::unique_ptr< std::vector< double > > fftjet_ScaleSet_parser(const edm::ParameterSet &ps)
void endJob() override
double peakPt(const fftjet::Peak &peak) const