Hall-D Software  alpha
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
PIDStudySelector.C
Go to the documentation of this file.
1 #define PIDStudySelector_cxx
2 // The class definition in PIDStudySelector.h has been generated automatically
3 // by the ROOT utility TTree::MakeSelector(). This class is derived
4 // from the ROOT class TSelector. For more information on the TSelector
5 // framework see $ROOTSYS/README/README.SELECTOR or the ROOT User Manual.
6 
7 // The following methods are defined in this file:
8 // Begin(): called every time a loop on the tree starts,
9 // a convenient place to create your histograms.
10 // SlaveBegin(): called after Begin(), when on PROOF called only on the
11 // slave servers.
12 // Process(): called for each event, in this function you decide what
13 // to read and fill your histograms.
14 // SlaveTerminate: called at the end of the loop on the tree, when on PROOF
15 // called only on the slave servers.
16 // Terminate(): called at the end of the loop on the tree,
17 // a convenient place to draw/fit your histograms.
18 //
19 // To use this file, try the following session on your Tree T:
20 //
21 // Root > T->Process("PIDStudySelector.C")
22 // Root > T->Process("PIDStudySelector.C","some options")
23 // Root > T->Process("PIDStudySelector.C+")
24 //
25 
26 #include "PIDStudySelector.h"
27 #include <TH2.h>
28 #include <TStyle.h>
29 
30 double Convert_RadiansToDegrees(double locRadians){
31  return locRadians*180.0/3.141592654;
32 }
33 
34 void PIDStudySelector::Begin(TTree * /*tree*/)
35 {
36  // The Begin() function is called at the start of the query.
37  // When running with PROOF Begin() is only called on the client.
38  // The tree argument is deprecated (on PROOF 0 is passed).
39 
40  TString option = GetOption();
41  TH2F *loc2FHist;
42  TH1F *loc1FHist;
43  ostringstream locHistName, locHistTitle;
44  unsigned int loc_i, loc_j;
45 
46  vector<int> locParticleIDVector; //p, pi+, pi-, k+
47  vector<int> locParticleIDVector_Results; //p, pi+, pi-, k+, pbar, k-
48  vector<string> locParticleSymbolVector, locParticleSymbolVector_Results;
49  int locPID;
50  string locParticleSymbol;
51 
52  locParticleIDVector.push_back(14); locParticleIDVector.push_back(8); locParticleIDVector.push_back(9); locParticleIDVector.push_back(11);
53  locParticleIDVector_Results.push_back(14); locParticleIDVector_Results.push_back(8); locParticleIDVector_Results.push_back(9); locParticleIDVector_Results.push_back(11);
54  locParticleIDVector_Results.push_back(15); locParticleIDVector_Results.push_back(12);
55 
56  locParticleSymbolVector.push_back("p"); locParticleSymbolVector.push_back("#pi^{+}"); locParticleSymbolVector.push_back("#pi^{-}"); locParticleSymbolVector.push_back("K^{+}");
57  locParticleSymbolVector_Results.push_back("p"); locParticleSymbolVector_Results.push_back("#pi^{+}"); locParticleSymbolVector_Results.push_back("#pi^{-}"); locParticleSymbolVector_Results.push_back("K^{+}");
58  locParticleSymbolVector_Results.push_back("#bar{p}"); locParticleSymbolVector_Results.push_back("K^{-}");
59 
60  int locNumBinsP = 900;
61  int locNumBinsTheta = 600;
62  int locNumBinsConfidenceLevel = 1000;
63  float locRangeMinP = 0.0, locRangeMaxP = 9.0;
64  float locRangeMinTheta = 0.0, locRangeMaxTheta = 150.0;
65  float locRangeMinConfidenceLevel = 0.0, locRangeMaxConfidenceLevel = 1.0;
66 
67 
68  dOutputFile = new TFile("dh_PIDStudies.root", "RECREATE");
69 
70  for(loc_i = 0; loc_i < locParticleIDVector.size(); loc_i++){
71  locPID = locParticleIDVector[loc_i];
72  locParticleSymbol = locParticleSymbolVector[loc_i];
73 
74  locHistName.str("");
75  locHistName << "dh_GeneratedTracksPVsTheta_" << locPID;
76  locHistTitle.str("");
77  locHistTitle << "Generated " << locParticleSymbol << "'s;#theta (degrees);p (GeV/c)";
78  loc2FHist = new TH2F(locHistName.str().c_str(), locHistTitle.str().c_str(), locNumBinsTheta, locRangeMinTheta, locRangeMaxTheta, locNumBinsP, locRangeMinP, locRangeMaxP);
79 
80  locHistName.str("");
81  locHistName << "dh_NonReconstructedTracksPVsTheta_" << locPID;
82  locHistTitle.str("");
83  locHistTitle << "Non-Reconstructed " << locParticleSymbol << "'s;#theta (degrees);p (GeV/c)";
84  loc2FHist = new TH2F(locHistName.str().c_str(), locHistTitle.str().c_str(), locNumBinsTheta, locRangeMinTheta, locRangeMaxTheta, locNumBinsP, locRangeMinP, locRangeMaxP);
85 
86  locHistName.str("");
87  locHistName << "dh_ReconstructedTracksPVsTheta_" << locPID;
88  locHistTitle.str("");
89  locHistTitle << "Reconstructed " << locParticleSymbol << "'s;#theta (degrees);p (GeV/c)";
90  loc2FHist = new TH2F(locHistName.str().c_str(), locHistTitle.str().c_str(), locNumBinsTheta, locRangeMinTheta, locRangeMaxTheta, locNumBinsP, locRangeMinP, locRangeMaxP);
91 
92  locHistName.str("");
93  locHistName << "dh_IdentifiedTracksPVsTheta_" << locPID;
94  locHistTitle.str("");
95  locHistTitle << "Identified " << locParticleSymbol << "'s;#theta (degrees);p (GeV/c)";
96  loc2FHist = new TH2F(locHistName.str().c_str(), locHistTitle.str().c_str(), locNumBinsTheta, locRangeMinTheta, locRangeMaxTheta, locNumBinsP, locRangeMinP, locRangeMaxP);
97 
98  locHistName.str("");
99  locHistName << "dh_MisIdentifiedTracksPVsTheta_" << locPID;
100  locHistTitle.str("");
101  locHistTitle << "Misidentified " << locParticleSymbol << "'s;#theta (degrees);p (GeV/c)";
102  loc2FHist = new TH2F(locHistName.str().c_str(), locHistTitle.str().c_str(), locNumBinsTheta, locRangeMinTheta, locRangeMaxTheta, locNumBinsP, locRangeMinP, locRangeMaxP);
103 
104  for(loc_j = 0; loc_j < locParticleIDVector_Results.size(); loc_j++){
105  if(loc_j == loc_i)
106  continue;
107  locHistName.str("");
108  locHistName << "dh_MisIdentifiedAsTracksPVsTheta_" << locPID << "_" << locParticleIDVector_Results[loc_j];
109  locHistTitle.str("");
110  locHistTitle << "MisIdentified " << locParticleSymbol << "as " << locParticleSymbolVector_Results[loc_j] << ";#theta (degrees);p (GeV/c)";
111  loc2FHist = new TH2F(locHistName.str().c_str(), locHistTitle.str().c_str(), locNumBinsTheta, locRangeMinTheta, locRangeMaxTheta, locNumBinsP, locRangeMinP, locRangeMaxP);
112  }
113 
114  for(loc_j = 0; loc_j < locParticleIDVector_Results.size(); loc_j++){
115  locHistName.str("");
116  locHistName << "dh_ConfidenceLevel_ProblemArea_Timing_" << locPID << "_" << locParticleIDVector_Results[loc_j];
117  locHistTitle.str("");
118  locHistTitle << locParticleSymbol << " as " << locParticleSymbolVector_Results[loc_j] << ", 1.5 < p (GeV/c) < 3.0, 30#circ < #theta < 120#circ;Timing Confidence Level";
119  loc1FHist = new TH1F(locHistName.str().c_str(), locHistTitle.str().c_str(), locNumBinsConfidenceLevel, locRangeMinConfidenceLevel, locRangeMaxConfidenceLevel);
120 
121  locHistName.str("");
122  locHistName << "dh_ConfidenceLevel_ProblemArea_DCdEdx_" << locPID << "_" << locParticleIDVector_Results[loc_j];
123  locHistTitle.str("");
124  locHistTitle << locParticleSymbol << " as " << locParticleSymbolVector_Results[loc_j] << ", 1.5 < p (GeV/c) < 3.0, 30#circ < #theta < 120#circ;DC dEdx Confidence Level";
125  loc1FHist = new TH1F(locHistName.str().c_str(), locHistTitle.str().c_str(), locNumBinsConfidenceLevel, locRangeMinConfidenceLevel, locRangeMaxConfidenceLevel);
126  }
127 
128  //CL hists
129  locHistName.str("");
130  locHistName << "dh_ConfidenceLevel_Overall_" << locPID;
131  locHistTitle.str("");
132  locHistTitle << "Overall " << locParticleSymbol << " PID;Confidence Level";
133  loc1FHist = new TH1F(locHistName.str().c_str(), locHistTitle.str().c_str(), locNumBinsConfidenceLevel, locRangeMinConfidenceLevel, locRangeMaxConfidenceLevel);
134 
135  locHistName.str("");
136  locHistName << "dh_ConfidenceLevel_Tracking_" << locPID;
137  locHistTitle.str("");
138  locHistTitle << "Tracking " << locParticleSymbol << " PID;Confidence Level";
139  loc1FHist = new TH1F(locHistName.str().c_str(), locHistTitle.str().c_str(), locNumBinsConfidenceLevel, locRangeMinConfidenceLevel, locRangeMaxConfidenceLevel);
140 
141  locHistName.str("");
142  locHistName << "dh_ConfidenceLevel_DCdEdx_" << locPID;
143  locHistTitle.str("");
144  locHistTitle << "DC dEdx " << locParticleSymbol << " PID;Confidence Level";
145  loc1FHist = new TH1F(locHistName.str().c_str(), locHistTitle.str().c_str(), locNumBinsConfidenceLevel, locRangeMinConfidenceLevel, locRangeMaxConfidenceLevel);
146 
147  locHistName.str("");
148  locHistName << "dh_ConfidenceLevel_Timing_" << locPID;
149  locHistTitle.str("");
150  locHistTitle << "Timing " << locParticleSymbol << " PID;Confidence Level";
151  loc1FHist = new TH1F(locHistName.str().c_str(), locHistTitle.str().c_str(), locNumBinsConfidenceLevel, locRangeMinConfidenceLevel, locRangeMaxConfidenceLevel);
152 
153  locHistName.str("");
154  locHistName << "dh_ConfidenceLevel_TOFdEdx_" << locPID;
155  locHistTitle.str("");
156  locHistTitle << "TOFdEdx " << locParticleSymbol << " PID;Confidence Level";
157  loc1FHist = new TH1F(locHistName.str().c_str(), locHistTitle.str().c_str(), locNumBinsConfidenceLevel, locRangeMinConfidenceLevel, locRangeMaxConfidenceLevel);
158 
159  locHistName.str("");
160  locHistName << "dh_ConfidenceLevel_BCALdEdx_" << locPID;
161  locHistTitle.str("");
162  locHistTitle << "BCALdEdx " << locParticleSymbol << " PID;Confidence Level";
163  loc1FHist = new TH1F(locHistName.str().c_str(), locHistTitle.str().c_str(), locNumBinsConfidenceLevel, locRangeMinConfidenceLevel, locRangeMaxConfidenceLevel);
164  }
165 }
166 
167 void PIDStudySelector::SlaveBegin(TTree * /*tree*/)
168 {
169  // The SlaveBegin() function is called after the Begin() function.
170  // When running with PROOF SlaveBegin() is called on each slave server.
171  // The tree argument is deprecated (on PROOF 0 is passed).
172 
173  TString option = GetOption();
174 
175 }
176 
177 Bool_t PIDStudySelector::Process(Long64_t entry)
178 {
179  // The Process() function is called for each entry in the tree (or possibly
180  // keyed object in the case of PROOF) to be processed. The entry argument
181  // specifies which entry in the currently loaded tree is to be processed.
182  // It can be passed to either PIDStudySelector::GetEntry() or TBranch::GetEntry()
183  // to read either all or the required parts of the data. When processing
184  // keyed objects with PROOF, the object is already loaded and is available
185  // via the fObject pointer.
186  //
187  // This function should contain the "body" of the analysis. It can contain
188  // simple or elaborate selection criteria, run algorithms on the data
189  // of the event and typically fill histograms.
190  //
191  // The processing can be stopped by calling Abort().
192  //
193  // Use fStatus to set the return value of TTree::Process().
194  //
195  // The return value is currently not used.
196 
197 
198  GetEntry(entry);
199  if((entry % 1000) == 0)
200  cout << "entry " << (entry + 1) << " of " << fChain->GetEntries() << endl;
201 
202  unsigned int loc_i, loc_j;
203  const MCReconstructionStatus *locMCReconstructionStatus;
204  const ReconstructedHypothesis *locReconstructedHypothesis;
205  float locFOM, locThrownP, locThrownTheta_Degrees;
206  Particle_t locThrownID, locPID;
207  DLorentzVector locThrownFourMomentum, locThrownSpacetimeVertex, locFourMomentum, locSpacetimeVertex;
208  TH2F *loc2FHist;
209  TH1F *loc1FHist;
210  ostringstream locHistName;
211 
212  for(loc_i = 0; loc_i < dMCReconstructionStatusVector.size(); loc_i++){
213  locMCReconstructionStatus = dMCReconstructionStatusVector[loc_i];
214  locThrownID = locMCReconstructionStatus->dThrownID;
215 
216  locThrownFourMomentum = locMCReconstructionStatus->dThrownFourMomentum;
217  locThrownSpacetimeVertex = locMCReconstructionStatus->dThrownSpacetimeVertex;
218 
219  locThrownP = locThrownFourMomentum.P();
220  locThrownTheta_Degrees = Convert_RadiansToDegrees(locThrownFourMomentum.Theta());
221 
222  locHistName.str("");
223  locHistName << "dh_GeneratedTracksPVsTheta_" << int(locThrownID);
224  loc2FHist = (TH2F*)gDirectory->Get(locHistName.str().c_str());
225  loc2FHist->Fill(locThrownTheta_Degrees, locThrownP);
226 
227  if(locMCReconstructionStatus->dReconstructedHypothesisVector.size() == 0){
228  locHistName.str("");
229  locHistName << "dh_NonReconstructedTracksPVsTheta_" << int(locThrownID);
230  loc2FHist = (TH2F*)gDirectory->Get(locHistName.str().c_str());
231  loc2FHist->Fill(locThrownTheta_Degrees, locThrownP);
232  }
233 
234  for(loc_j = 0; loc_j < locMCReconstructionStatus->dReconstructedHypothesisVector.size(); loc_j++){
235  locReconstructedHypothesis = locMCReconstructionStatus->dReconstructedHypothesisVector[loc_j];
236  locPID = locReconstructedHypothesis->dPID;
237  locFourMomentum = locReconstructedHypothesis->dFourMomentum;
238  locSpacetimeVertex = locReconstructedHypothesis->dSpacetimeVertex;
239 
240  if(loc_j == 0){
241  locHistName.str("");
242  locHistName << "dh_ReconstructedTracksPVsTheta_" << int(locThrownID);
243  loc2FHist = (TH2F*)gDirectory->Get(locHistName.str().c_str());
244  loc2FHist->Fill(locThrownTheta_Degrees, locThrownP);
245 
246  if(locThrownID == locPID){
247  locHistName.str("");
248  locHistName << "dh_IdentifiedTracksPVsTheta_" << int(locThrownID);
249  loc2FHist = (TH2F*)gDirectory->Get(locHistName.str().c_str());
250  loc2FHist->Fill(locThrownTheta_Degrees, locThrownP);
251  }else{
252  locHistName.str("");
253  locHistName << "dh_MisIdentifiedTracksPVsTheta_" << int(locThrownID);
254  loc2FHist = (TH2F*)gDirectory->Get(locHistName.str().c_str());
255  loc2FHist->Fill(locThrownTheta_Degrees, locThrownP);
256 
257  locHistName.str("");
258  locHistName << "dh_MisIdentifiedAsTracksPVsTheta_" << int(locThrownID) << "_" << int(locPID);
259  loc2FHist = (TH2F*)gDirectory->Get(locHistName.str().c_str());
260  if(loc2FHist != NULL)
261  loc2FHist->Fill(locThrownTheta_Degrees, locThrownP);
262  }
263  }
264 
265 
266  if(locThrownID == locPID){
267  locFOM = TMath::Prob(locReconstructedHypothesis->dChiSq_Overall, locReconstructedHypothesis->dNDF_Overall);
268  locHistName.str("");
269  locHistName << "dh_ConfidenceLevel_Overall_" << int(locThrownID);
270  loc1FHist = (TH1F*)gDirectory->Get(locHistName.str().c_str());
271  loc1FHist->Fill(locFOM);
272 
273  locFOM = TMath::Prob(locReconstructedHypothesis->dChiSq_Tracking, locReconstructedHypothesis->dNDF_Tracking);
274  locHistName.str("");
275  locHistName << "dh_ConfidenceLevel_Tracking_" << int(locThrownID);
276  loc1FHist = (TH1F*)gDirectory->Get(locHistName.str().c_str());
277  loc1FHist->Fill(locFOM);
278 
279  locFOM = TMath::Prob(locReconstructedHypothesis->dChiSq_DCdEdx, locReconstructedHypothesis->dNDF_DCdEdx);
280  locHistName.str("");
281  locHistName << "dh_ConfidenceLevel_DCdEdx_" << int(locThrownID);
282  loc1FHist = (TH1F*)gDirectory->Get(locHistName.str().c_str());
283  loc1FHist->Fill(locFOM);
284 
285  locFOM = TMath::Prob(locReconstructedHypothesis->dChiSq_Timing, locReconstructedHypothesis->dNDF_Timing);
286  locHistName.str("");
287  locHistName << "dh_ConfidenceLevel_Timing_" << int(locThrownID);
288  loc1FHist = (TH1F*)gDirectory->Get(locHistName.str().c_str());
289  loc1FHist->Fill(locFOM);
290 /*
291  locFOM = TMath::Prob(locReconstructedHypothesis->dChiSq_TOFdEdx, locReconstructedHypothesis->dNDF_TOFdEdx);
292  locHistName.str("");
293  locHistName << "dh_ConfidenceLevel_TOFdEdx_" << int(locThrownID);
294  loc1FHist = (TH1F*)gDirectory->Get(locHistName.str().c_str());
295  loc1FHist->Fill(locFOM);
296 
297  locFOM = TMath::Prob(locReconstructedHypothesis->dChiSq_BCALdEdx, locReconstructedHypothesis->dNDF_BCALdEdx);
298  locHistName.str("");
299  locHistName << "dh_ConfidenceLevel_BCALdEdx_" << int(locThrownID);
300  loc1FHist = (TH1F*)gDirectory->Get(locHistName.str().c_str());
301  loc1FHist->Fill(locFOM);
302 */
303  }
304 
305  //problem area
306  if((int(locThrownID) == 8) && (locThrownP > 1.5) && (locThrownP < 3.0) && (locThrownTheta_Degrees > 30.0) && (locThrownTheta_Degrees < 120.0)){
307  locFOM = TMath::Prob(locReconstructedHypothesis->dChiSq_Timing, locReconstructedHypothesis->dNDF_Timing);
308  locHistName.str("");
309  locHistName << "dh_ConfidenceLevel_ProblemArea_Timing_" << int(locThrownID) << "_" << int(locPID);
310  loc1FHist = (TH1F*)gDirectory->Get(locHistName.str().c_str());
311  if(loc1FHist != NULL)
312  loc1FHist->Fill(locFOM);
313 
314  locFOM = TMath::Prob(locReconstructedHypothesis->dChiSq_DCdEdx, locReconstructedHypothesis->dNDF_DCdEdx);
315  locHistName.str("");
316  locHistName << "dh_ConfidenceLevel_ProblemArea_DCdEdx_" << int(locThrownID) << "_" << int(locPID);
317  loc1FHist = (TH1F*)gDirectory->Get(locHistName.str().c_str());
318  if(loc1FHist != NULL)
319  loc1FHist->Fill(locFOM);
320  }
321 
322  }
323 
324  }
325 
326 
327  //fraction of generated tracks reconstructed of each type as function of p & theta
328 
329 
330  //fraction of generated tracks id'd of each type, as each type, as function of p & theta
331 
332  //fraction of reconstructed tracks id'd of each type, as each type, as function of p & theta
333 
334  //CL of each PID source for each type in bins of p & theta //& % of events failing 1% CL cut
335 
336  return kTRUE;
337 }
338 
340 {
341  // The SlaveTerminate() function is called after all entries or objects
342  // have been processed. When running with PROOF SlaveTerminate() is called
343  // on each slave server.
344 
345 }
346 
348 {
349  // The Terminate() function is the last function to be called during
350  // a query. It always runs on the client, it can be used to present
351  // the results graphically or save the results to file.
352 
353  dOutputFile->Write();
354 }
355 
virtual Int_t GetEntry(Long64_t entry, Int_t getall=0)
virtual void SlaveTerminate()
TLorentzVector DLorentzVector
DLorentzVector dThrownSpacetimeVertex
vector< ReconstructedHypothesis * > dReconstructedHypothesisVector
double Convert_RadiansToDegrees(double locRadians)
vector< MCReconstructionStatus * > dMCReconstructionStatusVector
virtual void Begin(TTree *tree)
virtual void SlaveBegin(TTree *tree)
virtual void Terminate()
virtual Bool_t Process(Long64_t entry)
Particle_t
Definition: particleType.h:12