CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
TreeUtility.cc
Go to the documentation of this file.
1 #include "RecoParticleFlow/PFClusterTools/interface/TreeUtility.h"
2 #include "TBranch.h"
3 #include "TTree.h"
4 #include "DataFormats/ParticleFlowReco/interface/Calibratable.h"
5 #include <cmath>
6 #include <fstream>
7 #include <iostream>
8 #include <TF1.h>
9 
10 using namespace pftools;
11 TreeUtility::TreeUtility() {
12 }
13 
14 TreeUtility::~TreeUtility() {
15 }
16 
17 double deltaR(double eta1, double eta2, double phi1, double phi2) {
18  return sqrt(pow(eta1 - eta2, 2) + pow(phi1 - phi2, 2));
19 }
20 
21 unsigned TreeUtility::getCalibratablesFromRootFile(TChain& tree,
22  std::vector<Calibratable>& toBeFilled) {
23 
24  // f.cd("extraction");
25  // TTree* tree = (TTree*) f.Get("extraction/Extraction");
26  // if (tree == 0) {
27  // PFToolsException me("Couldn't open tree!");
28  // throw me;
29  // }
30  // std::cout << "Successfully opened file. Getting branches..."<< std::endl;
31  CalibratablePtr calib_ptr(new Calibratable());
32  //TBranch* calibBr = tree.GetBranch("Calibratable");
33  //spwBr->SetAddress(&spw);
34  tree.SetBranchAddress("Calibratable", &calib_ptr);
35 
36  std::cout << "Looping over tree's "<< tree.GetEntries() << " entries...\n";
37  for (unsigned entries(0); entries < tree.GetEntries(); entries++) {
38  tree.GetEntry(entries);
39  Calibratable c(*calib_ptr);
40  if (c.cands_num_ == 1 && (c.cluster_ecal_.size() + c.cluster_hcal_.size()) > 0)
41  toBeFilled.push_back(c);
42  }
43  std::cout << "Done." << std::endl;
44  return tree.GetEntries();
45 
46 }
47 
48 void TreeUtility::dumpCaloDataToCSV(TChain& tree, std::string csvFilename, double range, bool gaus) {
49 
50  CalibratablePtr calib_ptr(new Calibratable());
51 
52  tree.SetBranchAddress("Calibratable", &calib_ptr);
53  std::ofstream csvFile;
54  csvFile.open(csvFilename.c_str());
55 
56  std::cout << "Looping over tree's "<< tree.GetEntries() << " entries...\n";
57  unsigned writes(0);
58  TFile freq("freq.root", "recreate");
59  TH1F frequencies("f", "f", 50, 0, range);
60  TF1 g("g", "gaus(0)");
61  g.FixParameter(1, range/2.0);
62  g.FixParameter(0, 1),
63  g.FixParameter(2, range/4.0);
64 
65  for (unsigned entries(0); entries < tree.GetEntries(); entries++) {
66  tree.GetEntry(entries);
67  Calibratable c(*calib_ptr);
68  bool veto(false);
69 
70  //Check vetos as usual
71  if (c.cands_num_ > 1)
72  veto = true;
73  if(c.cluster_ecal_.size() == 0 && c.cluster_hcal_.size() == 0)
74  veto = true;
75  if(!veto) {
76  if(frequencies.GetBinContent(static_cast<int>(floor(c.sim_energyEvent_)) + 1) < (3000) * g.Eval(c.sim_energyEvent_) ) {
77  frequencies.Fill(static_cast<int>(floor(c.sim_energyEvent_)));
78  c.recompute();
79  csvFile << c.sim_energyEvent_ << "\t";
80  /*
81  csvFile << c.sim_energyEcal_ << "\t";
82  csvFile << c.sim_energyHcal_ << "\t";
83 
84 
85  csvFile << c.cluster_energyEcal_/range << "\t";
86  csvFile << c.cluster_energyHcal_/range << "\t";
87 
88  CaloWindow newEcalWindow(c.cluster_meanEcal_.eta_, c.cluster_meanEcal_.phi_, 5, 0.01, 3);
89  const std::vector<CalibratableElement>& ecal = c.cluster_ecal_;
90  std::vector<CalibratableElement>::const_iterator cit = ecal.begin();
91  for(; cit != ecal.end(); ++cit) {
92  const CalibratableElement& hit = *cit;
93  bool added = newEcalWindow.addHit(hit.eta_, hit.phi_, hit.energy_);
94  if(!added)
95  veto = true;
96  }
97  */
98 
99  csvFile << fabs(c.sim_eta_/2) << "\n";
100 
101  ++writes;
102  }
103  }
104 
105 
106  }
107  frequencies.Print("frequencies.eps");
108  frequencies.Write();
109  freq.Close();
110  std::cout << "Closing file " << csvFilename << " with " << writes << " entries.\n" << std::endl;
111 
112  csvFile.close();
113 
114 }
115 
116 unsigned TreeUtility::getParticleDepositsDirectly(TChain& sourceChain,
117  std::vector<ParticleDepositPtr>& toBeFilled, CalibrationTarget target,
118  DetectorElementPtr offset, DetectorElementPtr ecal,
119  DetectorElementPtr hcal, bool includeOffset) {
120 
121  CalibratablePtr calib_ptr(new Calibratable());
122  sourceChain.SetBranchAddress("Calibratable", &calib_ptr);
123  std::cout << __PRETTY_FUNCTION__ << std::endl;
124  std::cout << "WARNING: Using fabs() for eta value assignments!\n";
125  std::cout << "Cutting on > 1 PFCandidate.\n";
126  std::cout << "Looping over tree's "<< sourceChain.GetEntries() << " entries...\n";
127  //neither of these two are supported yet
128  if (target == UNDEFINED || target == PFELEMENT)
129  return 0;
130  unsigned count(0);
131 
132  for (unsigned entries(0); entries < sourceChain.GetEntries(); entries++) {
133  sourceChain.GetEntry(entries);
134  Calibratable c(*calib_ptr);
135 
136  ParticleDepositPtr pd(new ParticleDeposit());
137  bool veto(false);
138  if (c.sim_isMC_) {
139  pd->setTruthEnergy(c.sim_energyEvent_);
140  pd->setEta(fabs(c.sim_eta_));
141  pd->setPhi(c.sim_phi_);
142  //TODO:: sort this out
143  if (c.sim_energyEvent_== 0)
144  veto = true;
145  }
146  if (c.tb_isTB_) {
147  pd->setTruthEnergy(c.sim_energyEvent_);
148  pd->setEta(c.tb_eta_);
149  pd->setPhi(c.tb_phi_);
150  veto = false;
151  }
152 
153  if (c.cands_num_ > 1)
154  veto = true;
155 
156  std::cout << "WARNING: HARD CUT ON 100 GeV SIM PARTICLES!\n";
157  if(c.sim_energyEvent_ > 100)
158  veto = true;
159 
160  if (target == CLUSTER) {
161  if (c.cluster_ecal_.size() == 0 && c.cluster_hcal_.size() ==0)
162  veto = true;
163  // if (c.cluster_numEcal_ > 1|| c.cluster_numHcal_ > 1)
164  // veto = true;
165  //TODO: using fabs for eta! WARNING!!!
166  Deposition decal(ecal, fabs(c.cluster_meanEcal_.eta_),
167  c.cluster_meanEcal_.phi_, c.cluster_energyEcal_, 0);
168  Deposition dhcal(hcal, fabs(c.cluster_meanHcal_.eta_),
169  c.cluster_meanHcal_.phi_, c.cluster_energyHcal_, 0);
170  Deposition doffset(offset, fabs(c.cluster_meanEcal_.eta_),
171  c.cluster_meanEcal_.phi_, 0.001, 0);
172 
173  pd->addTruthDeposition(decal);
174  pd->addRecDeposition(decal);
175 
176  pd->addTruthDeposition(dhcal);
177  pd->addRecDeposition(dhcal);
178 
179  if (includeOffset) {
180  pd->addTruthDeposition(doffset);
181  pd->addRecDeposition(doffset);
182  }
183 
184  }
185 
186  else if (target == PFCANDIDATE) {
187  // if(c.cands_num_ != 1)
188  // veto = true;
189  Deposition decal(ecal, c.cand_eta_, c.cand_phi_,
190  c.cand_energyEcal_, 0);
191  Deposition dhcal(hcal, c.cand_eta_, c.cand_phi_,
192  c.cand_energyHcal_, 0);
193  Deposition doffset(offset, c.cand_eta_, c.cand_phi_, 1.0, 0);
194 
195  pd->addTruthDeposition(decal);
196  pd->addTruthDeposition(dhcal);
197  pd->addRecDeposition(decal);
198  pd->addRecDeposition(dhcal);
199 
200  if (includeOffset) {
201  pd->addTruthDeposition(doffset);
202  pd->addRecDeposition(doffset);
203  }
204  }
205 
206  else if (target == RECHIT) {
207  if (c.rechits_ecal_.size() == 0&& c.rechits_hcal_.size() == 0)
208  veto = true;
209  Deposition decal(ecal, c.rechits_meanEcal_.eta_,
210  c.rechits_meanEcal_.phi_, c.rechits_meanEcal_.energy_
211  * c.rechits_ecal_.size(), 0);
212  Deposition dhcal(hcal, c.rechits_meanHcal_.eta_,
213  c.rechits_meanHcal_.phi_, c.rechits_meanHcal_.energy_
214  * c.rechits_hcal_.size(), 0);
215  Deposition doffset(offset, c.rechits_meanEcal_.eta_,
216  c.rechits_meanEcal_.phi_, 1.0, 0);
217 
218  pd->addTruthDeposition(decal);
219  pd->addTruthDeposition(dhcal);
220  pd->addRecDeposition(decal);
221  pd->addRecDeposition(dhcal);
222 
223  if (includeOffset) {
224  pd->addTruthDeposition(doffset);
225  pd->addRecDeposition(doffset);
226  }
227 
228  }
229  if (!veto)
230  toBeFilled.push_back(pd);
231 
232  ++count;
233  }
234 
235  return toBeFilled.size();
236 }
237 
238 unsigned TreeUtility::convertCalibratablesToParticleDeposits(
239  const std::vector<Calibratable>& input,
240  std::vector<ParticleDepositPtr>& toBeFilled, CalibrationTarget target,
241  DetectorElementPtr offset, DetectorElementPtr ecal,
242  DetectorElementPtr hcal, bool includeOffset) {
243 
244  std::cout << __PRETTY_FUNCTION__ << std::endl;
245  std::cout << "WARNING: Using fabs() for eta value assignments!\n";
246  std::cout << "Input Calibratable has size "<< input.size() << "\n";
247  std::cout << "Cutting on > 1 PFCandidate.\n";
248 
249  //neither of these two are supported yet
250  if (target == UNDEFINED || target == PFELEMENT)
251  return 0;
252  unsigned count(0);
253  for (std::vector<Calibratable>::const_iterator cit = input.begin(); cit
254  != input.end(); ++cit) {
255  Calibratable c = *cit;
256  ParticleDepositPtr pd(new ParticleDeposit());
257  bool veto(false);
258  if (c.sim_isMC_) {
259  pd->setTruthEnergy(c.sim_energyEvent_);
260  pd->setEta(fabs(c.sim_eta_));
261  pd->setPhi(c.sim_phi_);
262  //TODO:: sort this out
263  if (c.sim_energyEvent_== 0)
264  veto = true;
265  }
266  if (c.tb_isTB_) {
267  pd->setTruthEnergy(c.sim_energyEvent_);
268  pd->setEta(c.tb_eta_);
269  pd->setPhi(c.tb_phi_);
270  veto = false;
271  }
272 
273  if (c.cands_num_ > 1)
274  veto = true;
275 
276  if (target == CLUSTER) {
277  if (c.cluster_ecal_.size() == 0&& c.cluster_hcal_.size() ==0)
278  veto = true;
279  // if (c.cluster_numEcal_ > 1|| c.cluster_numHcal_ > 1)
280  // veto = true;
281  //TODO: using fabs for eta! WARNING!!!
282  Deposition decal(ecal, fabs(c.cluster_meanEcal_.eta_),
283  c.cluster_meanEcal_.phi_, c.cluster_energyEcal_, 0);
284  Deposition dhcal(hcal, fabs(c.cluster_meanHcal_.eta_),
285  c.cluster_meanHcal_.phi_, c.cluster_energyHcal_, 0);
286  Deposition doffset(offset, fabs(c.cluster_meanEcal_.eta_),
287  c.cluster_meanEcal_.phi_, 0.001, 0);
288 
289  pd->addTruthDeposition(decal);
290  pd->addRecDeposition(decal);
291 
292  pd->addTruthDeposition(dhcal);
293  pd->addRecDeposition(dhcal);
294 
295  if (includeOffset) {
296  pd->addTruthDeposition(doffset);
297  pd->addRecDeposition(doffset);
298  }
299 
300  }
301 
302  else if (target == PFCANDIDATE) {
303  // if(c.cands_num_ != 1)
304  // veto = true;
305  Deposition decal(ecal, c.cand_eta_, c.cand_phi_,
306  c.cand_energyEcal_, 0);
307  Deposition dhcal(hcal, c.cand_eta_, c.cand_phi_,
308  c.cand_energyHcal_, 0);
309  Deposition doffset(offset, c.cand_eta_, c.cand_phi_, 1.0, 0);
310 
311  pd->addTruthDeposition(decal);
312  pd->addTruthDeposition(dhcal);
313  pd->addRecDeposition(decal);
314  pd->addRecDeposition(dhcal);
315 
316  if (includeOffset) {
317  pd->addTruthDeposition(doffset);
318  pd->addRecDeposition(doffset);
319  }
320  }
321 
322  else if (target == RECHIT) {
323  if (c.rechits_ecal_.size() == 0&& c.rechits_hcal_.size() == 0)
324  veto = true;
325  Deposition decal(ecal, c.rechits_meanEcal_.eta_,
326  c.rechits_meanEcal_.phi_, c.rechits_meanEcal_.energy_
327  * c.rechits_ecal_.size(), 0);
328  Deposition dhcal(hcal, c.rechits_meanHcal_.eta_,
329  c.rechits_meanHcal_.phi_, c.rechits_meanHcal_.energy_
330  * c.rechits_hcal_.size(), 0);
331  Deposition doffset(offset, c.rechits_meanEcal_.eta_,
332  c.rechits_meanEcal_.phi_, 1.0, 0);
333 
334  pd->addTruthDeposition(decal);
335  pd->addTruthDeposition(dhcal);
336  pd->addRecDeposition(decal);
337  pd->addRecDeposition(dhcal);
338 
339  if (includeOffset) {
340  pd->addTruthDeposition(doffset);
341  pd->addRecDeposition(doffset);
342  }
343 
344  }
345  if (!veto)
346  toBeFilled.push_back(pd);
347 
348  ++count;
349  }
350 
351  return toBeFilled.size();
352 
353 }
354 
pftools::Calibratable::rechits_meanHcal_
CalibratableElement rechits_meanHcal_
Definition: Calibratable.h:193
pftools::Calibratable
Wraps essential single particle calibration data ready for export to a Root file. ...
Definition: Calibratable.h:122
pftools::TreeUtility::convertCalibratablesToParticleDeposits
unsigned convertCalibratablesToParticleDeposits(const std::vector< Calibratable > &input, std::vector< ParticleDepositPtr > &toBeFilled, CalibrationTarget target, DetectorElementPtr offset, DetectorElementPtr ecal, DetectorElementPtr hcal, bool includeOffset=false)
Definition: TreeUtility.cc:238
pftools::Calibratable::cluster_meanEcal_
CalibratableElement cluster_meanEcal_
Definition: Calibratable.h:186
pftools::CalibratablePtr
boost::shared_ptr< Calibratable > CalibratablePtr
Definition: TreeUtility.h:26
pftools::TreeUtility::dumpCaloDataToCSV
void dumpCaloDataToCSV(TChain &chain, std::string csvFilename, double range, bool gaus=false)
Definition: TreeUtility.cc:48
pftools::Calibratable::rechits_ecal_
std::vector< CalibratableElement > rechits_ecal_
Definition: Calibratable.h:191
pftools::TreeUtility::~TreeUtility
virtual ~TreeUtility()
Definition: TreeUtility.cc:14
g
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e g
Definition: Activities.doc:4
AlCaHLTBitMon_QueryRunRegistry.string
string string
Definition: AlCaHLTBitMon_QueryRunRegistry.py:255
pftools::TreeUtility::getParticleDepositsDirectly
unsigned getParticleDepositsDirectly(TChain &sourceChain, std::vector< ParticleDepositPtr > &toBeFilled, CalibrationTarget target, DetectorElementPtr offset, DetectorElementPtr ecal, DetectorElementPtr hcal, bool includeOffset=false)
Definition: TreeUtility.cc:116
input
static std::string const input
Definition: EdmProvDump.cc:44
pftools::ParticleDepositPtr
boost::shared_ptr< ParticleDeposit > ParticleDepositPtr
Definition: ParticleDeposit.h:103
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:48
pftools::TreeUtility::getCalibratablesFromRootFile
unsigned getCalibratablesFromRootFile(TChain &tree, std::vector< Calibratable > &toBeFilled)
Definition: TreeUtility.cc:21
pftools::Calibratable::rechits_meanEcal_
CalibratableElement rechits_meanEcal_
Definition: Calibratable.h:193
pftools::Calibratable::rechits_hcal_
std::vector< CalibratableElement > rechits_hcal_
Definition: Calibratable.h:191
evf::evtn::offset
unsigned int offset(bool)
Definition: GlobalEventNumber.cc:60
pftools::DetectorElementPtr
boost::shared_ptr< DetectorElement > DetectorElementPtr
Definition: DetectorElement.h:83
pftools::Calibratable::cluster_ecal_
std::vector< CalibratableElement > cluster_ecal_
Definition: Calibratable.h:184
deltaR
double deltaR(double eta1, double eta2, double phi1, double phi2)
Definition: TreeUtility.cc:17
pftools::Deposition
This class holds an arbitrary energy deposition, specified in terms of angular position, energy, depth (optional) and detector element type.
Definition: Deposition.h:20
patCandidatesForDimuonsSequences_cff.ecal
tuple ecal
Definition: patCandidatesForDimuonsSequences_cff.py:70
pftools::Calibratable::recompute
virtual void recompute()
Definition: Calibratable.cc:45
pftools::Calibratable::cluster_hcal_
std::vector< CalibratableElement > cluster_hcal_
Definition: Calibratable.h:184
pftools::Calibratable::cluster_meanHcal_
CalibratableElement cluster_meanHcal_
Definition: Calibratable.h:186
gather_cfg.cout
tuple cout
Definition: gather_cfg.py:121
patCandidatesForDimuonsSequences_cff.hcal
tuple hcal
Definition: patCandidatesForDimuonsSequences_cff.py:71
funct::pow
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
Generated for CMSSW Reference Manual by doxygen 1.8.5