FairRoot/PandaRoot
Functions
master/ana_complete.C File Reference

Go to the source code of this file.

Functions

int SelectTruePid (PndAnalysis *ana, RhoCandList &l)
 
void printCand (RhoCandidate *c)
 
void countDoubles (RhoCandList &l, int &n1, int &n2, int &n3)
 
int ana_complete (int nevts=0, TString prefix="evtcomplete")
 

Function Documentation

int ana_complete ( int  nevts = 0,
TString  prefix = "evtcomplete" 
)

Definition at line 56 of file master/ana_complete.C.

References RhoKinFitter::Add4MomConstraint(), PndMasterRunAna::AddFriend(), RhoKinFitter::AddMassConstraint(), RhoCandList::Combine(), countDoubles(), PndAnalysis::FillList(), PndMasterRunAna::Finish(), RhoKinFitter::Fit(), RhoFitterBase::Fit(), fRun, RhoFitterBase::GetChi2(), PndAnalysis::GetEntries(), PndAnalysis::GetEvent(), RhoCandList::GetLength(), RhoFitterBase::GetProb(), hvpos, i, RhoCandidate::M(), PndAnalysis::McTruthMatch(), muminus, muplus, out, output, P, piminus, piplus, PndAnalysis::PndAnalysis(), RhoCandidate::Pos(), RhoCandList::Select(), SelectTruePid(), PndMasterRunAna::SetInput(), PndMasterRunAna::SetOutput(), PndMasterRunAna::SetParamAsciiFile(), RhoCandList::SetType(), PndMasterRunAna::Setup(), and TString.

57 {
58  //-----User Settings:------------------------------------------------------
59  TString parAsciiFile = "all.par";
60  TString input = "psi2s_Jpsi2pi_Jpsi_mumu.dec";
61  TString output = "ana";
62  TString friend1 = "reco_single";
63  TString friend2 = "pid_single";
64  TString friend3 = "";
65  TString friend4 = "";
66 
67  // ----- Initial Settings --------------------------------------------
68  PndMasterRunAna *fRun= new PndMasterRunAna();
69  fRun->SetInput(input);
70  fRun->SetOutput(output);
71  fRun->AddFriend(friend1);
72  fRun->AddFriend(friend2);
73  fRun->AddFriend(friend3);
74  fRun->AddFriend(friend4);
75  fRun->SetParamAsciiFile(parAsciiFile);
76  fRun->Setup(prefix);
77 
78  // *** some variables
79  int i=0,j=0, k=0, l=0;
80  gStyle->SetOptFit(1011);
81 
82  fRun->Init();
83 
84  // *** create an output file for all histograms
85  TFile *out = TFile::Open(prefix+"_output_ana.root","RECREATE");
86 
87  // *** create some histograms
88  TH1F *hmomtrk = new TH1F("hmomtrk","track momentum (all)",200,0,5);
89  TH1F *hthttrk = new TH1F("hthttrk","track theta (all)",200,0,3.1415);
90 
91  TH1F *hjpsim_all = new TH1F("hjpsim_all","J/#psi mass (all)",200,0,4.5);
92  TH1F *hpsim_all = new TH1F("hpsim_all","#psi(2S) mass (all)",200,0,5);
93 
94  TH1F *hjpsim_lpid = new TH1F("hjpsim_lpid","J/#psi mass (loose pid)",200,0,4.5);
95  TH1F *hpsim_lpid = new TH1F("hpsim_lpid","#psi(2S) mass (loose pid)",200,0,5);
96 
97  TH1F *hjpsim_tpid = new TH1F("hjpsim_tpid","J/#psi mass (tight pid)",200,0,4.5);
98  TH1F *hpsim_tpid = new TH1F("hpsim_tpid","#psi(2S) mass (tight pid)",200,0,5);
99 
100  TH1F *hjpsim_trpid = new TH1F("hjpsim_trpid","J/#psi mass (true pid)",200,0,4.5);
101  TH1F *hpsim_trpid = new TH1F("hpsim_trpid","#psi(2S) mass (true pid)",200,0,5);
102 
103 
104  TH1F *hjpsim_ftm = new TH1F("hjpsim_ftm","J/#psi mass (full truth match)",200,0,4.5);
105  TH1F *hpsim_ftm = new TH1F("hpsim_ftm","#psi(2S) mass (full truth match)",200,0,5);
106 
107  TH1F *hjpsim_nm = new TH1F("hjpsim_nm","J/#psi mass (no truth match)",200,0,4.5);
108  TH1F *hpsim_nm = new TH1F("hpsim_nm","#psi(2S) mass (no truth match)",200,0,5);
109 
110  TH1F *hjpsim_diff = new TH1F("hjpsim_diff","J/#psi mass diff to truth",100,-2,2);
111  TH1F *hpsim_diff = new TH1F("hpsim_diff","#psi(2S) mass diff to truth",100,-2,2);
112 
113 
114  TH1F *hjpsim_vf = new TH1F("hjpsim_vf","J/#psi mass (vertex fit)",200,0,4.5);
115  TH1F *hjpsim_4cf = new TH1F("hjpsim_4cf","J/#psi mass (4C fit)",200,0,4.5);
116  TH1F *hjpsim_mcf = new TH1F("hjpsim_mcf","J/#psi mass (mass constraint fit)",200,0,4.5);
117 
118  TH1F *hjpsi_chi2_vf = new TH1F("hjpsi_chi2_vf", "J/#psi: #chi^{2} vertex fit",100,0,10);
119  TH1F *hpsi_chi2_4c = new TH1F("hpsi_chi2_4c", "#psi(2S): #chi^{2} 4C fit",100,0,250);
120  TH1F *hjpsi_chi2_mf = new TH1F("hjpsi_chi2_mf", "J/#psi: #chi^{2} mass fit",100,0,10);
121 
122  TH1F *hjpsi_prob_vf = new TH1F("hjpsi_prob_vf", "J/#psi: Prob vertex fit",100,0,1);
123  TH1F *hpsi_prob_4c = new TH1F("hpsi_prob_4c", "#psi(2S): Prob 4C fit",100,0,1);
124  TH1F *hjpsi_prob_mf = new TH1F("hjpsi_prob_mf", "J/#psi: Prob mass fit",100,0,1);
125 
126  TH2F *hvpos = new TH2F("hvpos","(x,y) projection of fitted decay vertex",100,-2,2,100,-2,2);
127 
128  //
129  // Now the analysis stuff comes...
130  //
131 
132 
133  // *** the data reader object
134  PndAnalysis* theAnalysis = new PndAnalysis();
135  if (nevts==0) nevts= theAnalysis->GetEntries();
136 
137  // *** RhoCandLists for the analysis
138  RhoCandList chrg, muplus, muminus, piplus, piminus, jpsi, psi2s;
139 
140  // *** Mass selector for the jpsi cands
141  double m0_jpsi = TDatabasePDG::Instance()->GetParticle("J/psi")->Mass(); // Get nominal PDG mass of the J/psi
142  RhoMassParticleSelector *jpsiMassSel=new RhoMassParticleSelector("jpsi",m0_jpsi,1.0);
143 
144  // *** the lorentz vector of the initial psi(2S)
145  TLorentzVector ini(0, 0, 6.231552, 7.240065);
146 
147  // ***
148  // the event loop
149  // ***
150 
151  int cntdbltrk=0, cntdblmc=0, cntdblboth=0, cnttrk=0, cnt_dbl_jpsi=0, cnt_dbl_psip=0;
152 
153  while (theAnalysis->GetEvent() && i++<nevts)
154  {
155  if ((i%100)==0) cout<<"evt " << i << endl;
156 
157  // *** Select with no PID info ('All'); type and mass are set
158  theAnalysis->FillList(chrg, "Charged");
159  theAnalysis->FillList(muplus, "MuonAllPlus");
160  theAnalysis->FillList(muminus, "MuonAllMinus");
161  theAnalysis->FillList(piplus, "PionAllPlus");
162  theAnalysis->FillList(piminus, "PionAllMinus");
163 
164  // *** momentum and theta histograms
165  for (j=0;j<muplus.GetLength();++j)
166  {
167  hmomtrk->Fill(muplus[j]->P());
168  hthttrk->Fill(muplus[j]->P4().Theta());
169  }
170  for (j=0;j<muminus.GetLength();++j)
171  {
172  hmomtrk->Fill(muminus[j]->P());
173  hthttrk->Fill(muminus[j]->P4().Theta());
174  }
175 
176  cnttrk += chrg.GetLength();
177 
178  int n1, n2, n3;
179 
180  countDoubles(chrg,n1,n2,n3);
181  cntdbltrk += n1;
182  cntdblmc += n2;
183  cntdblboth += n3;
184 
185  // *** combinatorics for J/psi -> mu+ mu-
186  jpsi.Combine(muplus, muminus);
187 
188 
189  // ***
190  // *** do the TRUTH MATCH for jpsi
191  // ***
192  jpsi.SetType(443);
193 
194  int nm = 0;
195  for (j=0;j<jpsi.GetLength();++j)
196  {
197  hjpsim_all->Fill( jpsi[j]->M() );
198 
199  if (theAnalysis->McTruthMatch(jpsi[j]))
200  {
201  nm++;
202  hjpsim_ftm->Fill( jpsi[j]->M() );
203  hjpsim_diff->Fill( jpsi[j]->GetMcTruth()->M() - jpsi[j]->M() );
204  }
205  else
206  hjpsim_nm->Fill( jpsi[j]->M() );
207  }
208 
209  if (nm>1) cnt_dbl_jpsi++;
210  // ***
211  // *** do VERTEX FIT (J/psi)
212  // ***
213  for (j=0;j<jpsi.GetLength();++j)
214  {
215  RhoKinVtxFitter vtxfitter(jpsi[j]); // instantiate a vertex fitter
216  vtxfitter.Fit();
217 
218  double chi2_vtx = vtxfitter.GetChi2(); // access chi2 of fit
219  double prob_vtx = vtxfitter.GetProb(); // access probability of fit
220  hjpsi_chi2_vf->Fill(chi2_vtx);
221  hjpsi_prob_vf->Fill(prob_vtx);
222 
223  if ( prob_vtx > 0.01 ) // when good enough, fill some histos
224  {
225  RhoCandidate *jfit = jpsi[j]->GetFit(); // access the fitted cand
226  TVector3 jVtx=jfit->Pos(); // and the decay vertex position
227 
228  hjpsim_vf->Fill(jfit->M());
229  hvpos->Fill(jVtx.X(),jVtx.Y());
230  }
231  }
232 
233  // *** some rough mass selection
234  jpsi.Select(jpsiMassSel);
235 
236  // *** combinatorics for psi(2S) -> J/psi pi+ pi-
237  psi2s.Combine(jpsi, piplus, piminus);
238 
239 
240  // ***
241  // *** do the TRUTH MATCH for psi(2S)
242  // ***
243  psi2s.SetType(88888);
244 
245  nm = 0;
246  for (j=0;j<psi2s.GetLength();++j)
247  {
248  hpsim_all->Fill( psi2s[j]->M() );
249 
250  if (theAnalysis->McTruthMatch(psi2s[j]))
251  {
252  nm++;
253  hpsim_ftm->Fill( psi2s[j]->M() );
254  hpsim_diff->Fill( psi2s[j]->GetMcTruth()->M() - psi2s[j]->M() );
255  }
256  else
257  hpsim_nm->Fill( psi2s[j]->M() );
258  }
259  if (nm>1) cnt_dbl_psip++;
260 
261 
262  // ***
263  // *** do 4C FIT (initial psi(2S) system)
264  // ***
265  for (j=0;j<psi2s.GetLength();++j)
266  {
267  RhoKinFitter fitter(psi2s[j]); // instantiate the kin fitter in psi(2S)
268  fitter.Add4MomConstraint(ini); // set 4 constraint
269  fitter.Fit(); // do fit
270 
271  double chi2_4c = fitter.GetChi2(); // get chi2 of fit
272  double prob_4c = fitter.GetProb(); // access probability of fit
273  hpsi_chi2_4c->Fill(chi2_4c);
274  hpsi_prob_4c->Fill(prob_4c);
275 
276  if ( prob_4c > 0.01 ) // when good enough, fill some histo
277  {
278  RhoCandidate *jfit = psi2s[j]->Daughter(0)->GetFit(); // get fitted J/psi
279 
280  hjpsim_4cf->Fill(jfit->M());
281  }
282  }
283 
284 
285  // ***
286  // *** do MASS CONSTRAINT FIT (J/psi)
287  // ***
288  for (j=0;j<jpsi.GetLength();++j)
289  {
290  RhoKinFitter mfitter(jpsi[j]); // instantiate the RhoKinFitter in psi(2S)
291  mfitter.AddMassConstraint(m0_jpsi); // add the mass constraint
292  mfitter.Fit(); // do fit
293 
294  double chi2_m = mfitter.GetChi2(); // get chi2 of fit
295  double prob_m = mfitter.GetProb(); // access probability of fit
296  hjpsi_chi2_mf->Fill(chi2_m);
297  hjpsi_prob_mf->Fill(prob_m);
298 
299  if ( prob_m > 0.01 ) // when good enough, fill some histo
300  {
301  RhoCandidate *jfit = jpsi[j]->GetFit(); // access the fitted cand
302  hjpsim_mcf->Fill(jfit->M());
303  }
304  }
305 // cout<<"event "<<i<<" DONE"<<endl;
306 
307 
308  // ***
309  // *** TRUE PID combinatorics
310  // ***
311 
312  // *** do MC truth match for PID type
313  SelectTruePid(theAnalysis, muplus);
314  SelectTruePid(theAnalysis, muminus);
315  SelectTruePid(theAnalysis, piplus);
316  SelectTruePid(theAnalysis, piminus);
317 
318  // *** all combinatorics again with true PID
319  jpsi.Combine(muplus, muminus);
320  for (j=0;j<jpsi.GetLength();++j) hjpsim_trpid->Fill( jpsi[j]->M() );
321  jpsi.Select(jpsiMassSel);
322 
323  psi2s.Combine(jpsi, piplus, piminus);
324  for (j=0;j<psi2s.GetLength();++j) hpsim_trpid->Fill( psi2s[j]->M() );
325 
326 
327  // ***
328  // *** LOOSE PID combinatorics
329  // ***
330 
331  // *** and again with PidAlgoMvd;PidAlgoStt;PidAlgoDrc and loose selection
332  theAnalysis->FillList(muplus, "MuonLoosePlus", "PidAlgoMvd;PidAlgoStt;PidAlgoDrc");
333  theAnalysis->FillList(muminus, "MuonLooseMinus", "PidAlgoMvd;PidAlgoStt;PidAlgoDrc");
334  theAnalysis->FillList(piplus, "PionLoosePlus", "PidAlgoMvd;PidAlgoStt;PidAlgoDrc");
335  theAnalysis->FillList(piminus, "PionLooseMinus", "PidAlgoMvd;PidAlgoStt;PidAlgoDrc");
336 
337  jpsi.Combine(muplus, muminus);
338  for (j=0;j<jpsi.GetLength();++j) hjpsim_lpid->Fill( jpsi[j]->M() );
339  jpsi.Select(jpsiMassSel);
340 
341  psi2s.Combine(jpsi, piplus, piminus);
342  for (j=0;j<psi2s.GetLength();++j) hpsim_lpid->Fill( psi2s[j]->M() );
343 
344 
345  // ***
346  // *** TIGHT PID combinatorics
347  // ***
348 
349  // *** and again with PidAlgoMvd;PidAlgoStt and tight selection
350  theAnalysis->FillList(muplus, "MuonTightPlus", "PidAlgoMdtHardCuts");
351  theAnalysis->FillList(muminus, "MuonTightMinus", "PidAlgoMdtHardCuts");
352  theAnalysis->FillList(piplus, "PionLoosePlus", "PidAlgoMvd;PidAlgoStt;PidAlgoDrc");
353  theAnalysis->FillList(piminus, "PionLooseMinus", "PidAlgoMvd;PidAlgoStt;PidAlgoDrc");
354 
355  jpsi.Combine(muplus, muminus);
356  for (j=0;j<jpsi.GetLength();++j) hjpsim_tpid->Fill( jpsi[j]->M() );
357  jpsi.Select(jpsiMassSel);
358 
359  psi2s.Combine(jpsi, piplus, piminus);
360  for (j=0;j<psi2s.GetLength();++j) hpsim_tpid->Fill( psi2s[j]->M() );
361 
362  }
363 
364  // *** write out all the histos
365  out->cd();
366 
367  hmomtrk->Write();
368  hthttrk->Write();
369 
370  hjpsim_all->Write();
371  hpsim_all->Write();
372  hjpsim_lpid->Write();
373  hpsim_lpid->Write();
374  hjpsim_tpid->Write();
375  hpsim_tpid->Write();
376  hjpsim_trpid->Write();
377  hpsim_trpid->Write();
378 
379  hjpsim_ftm->Write();
380  hpsim_ftm->Write();
381  hjpsim_nm->Write();
382  hpsim_nm->Write();
383 
384  hpsim_diff->Write();
385  hjpsim_diff->Write();
386 
387  hjpsim_vf->Write();
388  hjpsim_4cf->Write();
389  hjpsim_mcf->Write();
390 
391  hjpsi_chi2_vf->Write();
392  hpsi_chi2_4c->Write();
393  hjpsi_chi2_mf->Write();
394 
395  hjpsi_prob_vf->Write();
396  hpsi_prob_4c->Write();
397  hjpsi_prob_mf->Write();
398 
399  hvpos->Write();
400 
401  out->Save();
402 
403  fRun->Finish();
404  return 0;
405 
406 }
PndAnalysis::GetEntries
Int_t GetEntries()
Definition: PndAnalysis.cxx:471
SelectTruePid
int SelectTruePid(PndAnalysis *ana, RhoCandList &l)
Definition: master/ana_complete.C:8
PndMasterRunAna
Class for the master reconstruction chain.
Definition: PndMasterRunAna.h:21
i
Int_t i
Definition: run_full.C:25
PndMasterRunAna::Setup
Bool_t Setup(TString outprefix="")
Initial setup.
Definition: PndMasterRunAna.cxx:48
RhoCandidate::Pos
TVector3 Pos() const
Definition: RhoCandidate.h:186
RhoCandList::GetLength
Int_t GetLength() const
Definition: RhoCandList.h:46
PndAnalysis::FillList
Bool_t FillList(RhoCandList &l, TString listkey="All", TString pidTcaNames="", int trackHypothesis=-1)
Definition: PndAnalysis.cxx:319
piplus
CandList piplus
Definition: softtrigger_kin5.C:144
PndMasterRunAna::AddFriend
void AddFriend(TString par)
Setter of friend root files.
Definition: PndMasterRunAna.h:123
output
FairParRootFileIo * output
Definition: sim_emc_apd.C:120
fRun
FairRunAna * fRun
Definition: hit_dirc.C:58
muplus
CandList muplus
Definition: softtrigger_kin5.C:144
PndMasterRunAna::SetInput
void SetInput(TString par)
Input of the macro.
Definition: PndMasterRunAna.h:103
RhoCandList::Combine
void Combine(RhoCandList &l1, RhoCandList &l2)
Definition: RhoCandList.cxx:289
RhoCandList::Select
void Select(RhoParticleSelectorBase *pidmgr)
Definition: RhoCandList.cxx:1359
RhoCandList::SetType
void SetType(const TParticlePDG *pdt, int start=0)
Definition: RhoCandList.cxx:264
out
TFile * out
Definition: reco_muo.C:20
RhoCandidate::M
Double_t M() const
Definition: RhoCandidate.cxx:1636
countDoubles
void countDoubles(RhoCandList &l, int &n1, int &n2, int &n3)
Definition: master/ana_complete.C:31
P
GeV c P
Definition: ReadMesDecCorrRes_PiT.C:42
PndMasterRunAna::Finish
void Finish()
Final diagnostics.
Definition: PndMasterRunAna.cxx:205
piminus
CandList piminus
Definition: softtrigger_kin5.C:144
PndAnalysis::McTruthMatch
Bool_t McTruthMatch(RhoCandidate *cand, Int_t level=2, bool verbose=false)
Definition: PndAnalysis.cxx:932
PndAnalysis::GetEvent
Int_t GetEvent(Int_t n=-1)
Definition: PndAnalysis.cxx:257
PndMasterRunAna::SetOutput
void SetOutput(TString par)
Tag of the output file of the macro.
Definition: PndMasterRunAna.h:108
muminus
CandList muminus
Definition: softtrigger_kin5.C:144
PndMasterRunAna::SetParamAsciiFile
void SetParamAsciiFile(TString par)
Setter of the parameter ascii file.
Definition: PndMasterRunAna.h:118
hvpos
TH2F * hvpos
Definition: plot_eta_c_stt.C:36
void countDoubles ( RhoCandList &  l,
int &  n1,
int &  n2,
int &  n3 
)

Definition at line 31 of file master/ana_complete.C.

References d, fabs(), RhoCandList::GetLength(), and i.

Referenced by ana_complete().

32 {
33  int n_smc = 0;
34  int n_strk = 0;
35  int n_both = 0;
36  double d = 0.00001;
37 
38  for (int i=0;i<l.GetLength()-1;++i)
39  {
40  for (int j=i+1;j<l.GetLength();++j)
41  {
42  TLorentzVector dl = l[i]->P4() - l[j]->P4();
43 
44  bool chkmc = (l[i]->GetMcTruth()==l[j]->GetMcTruth());
45  bool chktrk = (fabs(dl.X())<d) && (fabs(dl.Y())<d) && (fabs(dl.Z())<d) && (fabs(dl.E())<d);
46  if (chkmc) n_smc++;
47  if (chktrk) n_strk++;
48  if (chktrk && chkmc) n_both++;
49  }
50  }
51  n1 = n_strk;
52  n2 = n_smc;
53  n3 = n_both;
54 }
d
TObjArray * d
Definition: Pnd_Hc_etaee7G.C:78
i
Int_t i
Definition: run_full.C:25
RhoCandList::GetLength
Int_t GetLength() const
Definition: RhoCandList.h:46
fabs
friend F32vec4 fabs(const F32vec4 &a)
Definition: P4_F32vec4.h:47
void printCand ( RhoCandidate *  c)

Definition at line 24 of file master/ana_complete.C.

References RhoCandidate::P4(), and RhoCandidate::PdgCode().

25 {
26  TLorentzVector lv=c->P4();
27 
28  cout <<c->PdgCode()<<" ("<<lv.X()<<"/"<<lv.Y()<<"/"<<lv.Z()<<"/"<<lv.E()<<")"<<endl;
29 }
RhoCandidate::P4
TLorentzVector P4() const
Definition: RhoCandidate.h:195
int SelectTruePid ( PndAnalysis *  ana,
RhoCandList &  l 
)

Definition at line 8 of file master/ana_complete.C.

References RhoCandList::GetLength(), PndAnalysis::McTruthMatch(), and RhoCandList::Remove().

Referenced by ana_complete().

9 {
10  int removed = 0;
11 
12  for (int ii=l.GetLength()-1;ii>=0;--ii)
13  {
14  if ( !(ana->McTruthMatch(l[ii])) )
15  {
16  l.Remove(l[ii]);
17  removed++;
18  }
19  }
20 
21  return removed;
22 }
RhoCandList::GetLength
Int_t GetLength() const
Definition: RhoCandList.h:46
RhoCandList::Remove
Int_t Remove(RhoCandidate *)
Definition: RhoCandList.cxx:131
PndAnalysis::McTruthMatch
Bool_t McTruthMatch(RhoCandidate *cand, Int_t level=2, bool verbose=false)
Definition: PndAnalysis.cxx:932
Generated on Sun Aug 18 2019 03:02:20 for FairRoot/PandaRoot by   doxygen 1.8.8