FairRoot/PandaRoot
Functions
ana_check_psi.C File Reference

Go to the source code of this file.

Functions

int ana_check_psi (int nevts=0)
 
void removeCombinatoric (TCandList &psi_list)
 

Function Documentation

int ana_check_psi ( int  nevts = 0)

Definition at line 5 of file ana_check_psi.C.

References ctime, Double_t, RhoFitterBase::Fit(), Rho4CFitter::FitConserveMasses(), RhoFitterBase::GetChi2(), PndMCTrack::GetMotherID(), PndMCTrack::GetPdgCode(), h_chi2_4c, h_chi2b_4c, i, inFile, inSimFile, Mass, mc_array, n_events, nc, out, p1, p2, printf(), rtime, timer, tree, and TString.

6 {
7  TString OutFile = "output.root";
8 
9  gStyle->SetOptFit(1011);
10 
11  TStopwatch timer;
12  timer.Start();
13 
14  gROOT->Macro("$VMCWORKDIR/gconfig/rootlogon.C");
15 
16  TString inPidFile = "evt_pid.root";
17  TString inSimFile = "evt_points.root";
18 
19  TFile *inFile = TFile::Open(inSimFile,"READ");
20  TTree *tree=(TTree *) inFile->Get("pndsim") ;
21  tree->AddFriend("pndsim",inPidFile);
22 
23  TClonesArray* mc_array=new TClonesArray("PndMCTrack");
24  tree->SetBranchAddress("MCTrack",&mc_array);
25 
26  FairMCEventHeader* evthead;
27  tree->SetBranchAddress("MCEventHeader.", &evthead);
28 
29  TFile *out = TFile::Open(OutFile,"RECREATE");
30 
31  // the PndEventReader takes care about file/event handling
32  PndEventReader evr(inPidFile);
33 
34  TH1F *h_psi_nocut=new TH1F("h_psi_nocut","m(psi), (no cuts);E, GeV",100,3.2,4.4);
35  TH1F *h_mD1_nocuts=new TH1F("h_mD1_nocuts","D+: m(#pi+ #pi+ K-) (no cuts)",100,1.5,2.2);
36  TH1F *h_mD2_nocuts=new TH1F("h_mD2_nocuts","D-: m(#pi- #pi- K+) (no cuts)",100,1.5,2.2);
37  TH1F *h_mD1_pid=new TH1F("h_mD1_pid","D+: m(#pi+ #pi+ K-) (pid)",100,1.5,2.2);
38  TH1F *h_mD2_pid=new TH1F("h_mD2_pid","D-: m(#pi- #pi- K+) (pid)",100,1.5,2.2);
39 
40  TH1F *h_psi_pid=new TH1F("h_psi_pid","m(#psi), (MC PID);E, GeV",100,3.2,4.4);
41  TH1F *h_psi_pid_vtx=new TH1F("h_psi_pid_vtx","m(#psi), (MC PID);E, GeV",100,3.2,4.4);
42 
43  TH1F *h_psi_4c=new TH1F("h_psi_4c","m(#psi), 4C-fit",100,3.2,4.4);
44  TH1F *h_mD_4c=new TH1F("h_mD_4c","D: m(K #pi #pi) (4C-fit)",100,1.5,2.2);
45 
46  TH1F *h_psi_vtx=new TH1F("h_psi_vtx","m(#psi), Vertex fit",100,3.2,4.4);
47  TH1F *h_D1_vtx_bef=new TH1F("h_D1_vtx_bef","#phi: m(K #pi #pi) (Vertex fit)",100,1.5,2.2);
48  TH1F *h_D2_vtx_bef=new TH1F("h_D2_vtx_bef","#phi: m(K #pi #pi) (Vertex fit)",100,1.5,2.2);
49  TH1F *h_D1_vtx=new TH1F("h_D1_vtx","#phi: m(K #pi #pi) (Vertex fit)",100,1.5,2.2);
50  TH1F *h_D2_vtx=new TH1F("h_D2_vtx","#phi: m(K #pi #pi) (Vertex fit)",100,1.5,2.2);
51 
52  TH1F *h_psi_Dmass=new TH1F("h_psi_Dmass","m(#psi), (cut on D mass);E, GeV",100,3.2,4.4);
53  TH1F *h_mD_final=new TH1F("h_mD_final","D: m(K #pi #pi)",100,1.5,2.2);
54 
55  TH1F *nc=new TH1F("nc","n charged",20,0,20);
56 
57  TH1F *h_chi2_4c=new TH1F("h_chi2_4c","#chi^{2} 4C-fit;#chi^{2}/N_{df}",100,0,100);
58  TH1F *h_chi2b_4c=new TH1F("h_chi2b_4c","#chi^{2} 4C-fit;#chi^{2}/N_{df}",100,0,100);
59  TH1F *h_chi2_vtx_D1=new TH1F("h_chi2_vtx_D1","#chi^{2} vertex;#chi^{2}/N_{df}",100,0,100);
60  TH1F *h_chi2_vtx_D2=new TH1F("h_chi2_vtx_D2","#chi^{2} vertex;#chi^{2}/N_{df}",100,0,100);
61  TH2F *hvpos_D1 = new TH2F("hvpos_D1","(x,y) projection of fitted decay vertex",100,-5,5,100,-5,5);
62  TH1F *hvzpos_D1 = new TH1F("hvzpos_D1","z position of fitted decay vertex",100,-10,10);
63  TH1F *hvtxresX_D1 = new TH1F("hvtxresX_D1","X resolution of fitted decay vertex",100,-0.1,0.1);
64  TH1F *hvtxresY_D1 = new TH1F("hvtxresY_D1","Y resolution of fitted decay vertex",100,-0.1,0.1);
65  TH1F *hvtxresZ_D1 = new TH1F("hvtxresZ_D1","Z resolution of fitted decay vertex",100,-0.1,0.1);
66  TH2F *hvpos_D2 = new TH2F("hvpos_D2","(x,y) projection of fitted decay vertex",100,-5,5,100,-5,5);
67  TH1F *hvzpos_D2 = new TH1F("hvzpos_D2","z position of fitted decay vertex",100,-10,10);
68  TH1F *hvtxresX_D2 = new TH1F("hvtxresX_D2","X resolution of fitted decay vertex",100,-0.1,0.1);
69  TH1F *hvtxresY_D2 = new TH1F("hvtxresY_D2","Y resolution of fitted decay vertex",100,-0.1,0.1);
70  TH1F *hvtxresZ_D2 = new TH1F("hvtxresZ_D2","Z resolution of fitted decay vertex",100,-0.1,0.1);
71 
72  TPidMassSelector *psiMassSel=new TPidMassSelector("psi",1.87,0.07);
73 
74  TPidPlusSelector *pplusSel=new TPidPlusSelector("pplus");
75  TPidMinusSelector *pminusSel=new TPidMinusSelector("pminus");
76 
77  // the candidates lists we need
78  TCandList p1, p2, k1, k2, D1, D1_pid, D2, D2_pid, psi, psi_pid, psi_nocut;
79 
80  int n_reco=0;
81  // Number of events in file and number of reconstructed eta_c to store in root file
82  TH1F *n_events=new TH1F("n_events","total number of events",1,0,1);
83  TH1F *n_psi=new TH1F("n_psi","number of reconstructed eta_c",1,0,1);
84 
85  TLorentzVector ini(0,0,6.578800,7.583333);
86 
87  if (nevts==0) nevts=evr.GetEntries();
88 
89  n_events->SetBinContent(1,nevts);
90  // cout << "nevts " << nevts << "\n";
91  int i=0,j=0, k=0, l=0;
92 
93  // *************
94  // this is the loop through the events ... as simple as this...
95  // ****************
96  while (evr.GetEvent() && i++<nevts)
97  {
98  //.cout << "evt: " << i << endl;
99  if ((i%100)==0) cout<<"evt " << i << "\n";
100  //if (!((i+1)%100)) cout<<"evt " << i << "\n";
101  evr.FillList(p1,"Charged");
102  evr.FillList(p2,"Charged");
103  evr.FillList(k1,"Charged");
104  evr.FillList(k2,"Charged");
105 
106  p1.Select(pplusSel);
107  p2.Select(pminusSel);
108  k1.Select(pplusSel);
109  k2.Select(pminusSel);
110 
111  int nchrg=p1.GetLength()+p2.GetLength();
112  nc->Fill(nchrg);
113 
114  for (j=0;j<p1.GetLength();++j) {
115  p1[j].SetMass(TRho::Instance()->GetPDG()->GetParticle(211)->Mass());
116  }
117  for (j=0;j<p2.GetLength();++j) {
118  p2[j].SetMass(TRho::Instance()->GetPDG()->GetParticle(211)->Mass());
119  }
120  for (j=0;j<k1.GetLength();++j) {
121  k1[j].SetMass(TRho::Instance()->GetPDG()->GetParticle(321)->Mass());
122  }
123  for (j=0;j<k2.GetLength();++j) {
124  k2[j].SetMass(TRho::Instance()->GetPDG()->GetParticle(321)->Mass());
125  }
126 
127 
128  D1.Combine(p1,p1, k2);
129  D2.Combine(p2,p2, k1);
130 
131  for (j=0;j<D1.GetLength();++j) h_mD1_nocuts->Fill(D1[j].M());
132  for (j=0;j<D2.GetLength();++j) h_mD2_nocuts->Fill(D2[j].M());
133 
134  D1.Select(psiMassSel);
135  D2.Select(psiMassSel);
136  psi_nocut.Combine(D1,D2);
137 
138  for (l=0;l<psi_nocut.GetLength();++l) {
139  h_psi_nocut->Fill(psi_nocut[l].M());
140  }
141 
142  tree->GetEntry(i-1);
143  TVector3 mcVertex, mcD1Vertex, mcD2vertex;
144 
145  if (((PndMCTrack*)mc_array->At(0))!=0) mcD1Vertex = ((PndMCTrack*)mc_array->At(0))->GetStartVertex();
146  else
147  cout << "Not found k1!" << endl;
148  if (((PndMCTrack*)mc_array->At(3))!=0) mcD2Vertex = ((PndMCTrack*)mc_array->At(3))->GetStartVertex();
149  else
150  cout << "Not found k2!" << endl;
151 
152  evthead->GetVertex(mcVertex);
153 
154  // MC PID
155  // Leave only kaons in particle lists
156  int n_removed=0;
157  int ii=0;
158  for (l=0;l<p1.GetLength();++l) {
159  ii=l-n_removed;
160  if (p1[ii].GetMicroCandidate().GetMcIndex()>-1){
161  PndMCTrack *mcTrack = (PndMCTrack*)mc_array->At(p1[ii].GetMicroCandidate().GetMcIndex());
162  if (mcTrack!=0)
163  {
164  if ((mcTrack->GetPdgCode()!=211) || (mcTrack->GetMotherID()!=-1))
165  {
166  p1.Remove(p1[ii]);
167  n_removed++;
168  }
169  }
170  else
171  {
172  std::cout<<"stt h: " << p1[ii].GetMicroCandidate().GetSttHits() << std::endl;
173  std::cout<<"Kaon list 1, element "<<l<<" has no assosiated mcTRack"<<std::endl;
174  }
175  }
176  }
177 
178  n_removed=0;
179  ii=0;
180  for (l=0;l<p2.GetLength();++l) {
181  ii=l-n_removed;
182  if (p2[ii].GetMicroCandidate().GetMcIndex()>-1){
183  PndMCTrack *mcTrack = (PndMCTrack*)mc_array->At(p2[ii].GetMicroCandidate().GetMcIndex());
184  if (mcTrack!=0)
185  {
186  if ((mcTrack->GetPdgCode()!=-211) || (mcTrack->GetMotherID()!=-1))
187  {
188  p2.Remove(p2[ii]);
189  n_removed++;
190  }
191  }
192  else
193  {
194  std::cout<<"stt h: " << p2[ii].GetMicroCandidate().GetSttHits() << std::endl;
195  std::cout<<"Kaon list 2, element "<<l<<" has no assosiated mcTRack"<<std::endl;
196  }
197  }
198  }
199 
200  n_removed=0;
201  ii=0;
202  for (l=0;l<k1.GetLength();++l) {
203  ii=l-n_removed;
204  if (k1[ii].GetMicroCandidate().GetMcIndex()>-1){
205  PndMCTrack *mcTrack = (PndMCTrack*)mc_array->At(k1[ii].GetMicroCandidate().GetMcIndex());
206  if (mcTrack!=0)
207  {
208  if ((mcTrack->GetPdgCode()!=321) || (mcTrack->GetMotherID()!=-1))
209  {
210  k1.Remove(k1[ii]);
211  n_removed++;
212  }
213  }
214  else
215  {
216  std::cout<<"stt h: " << k1[ii].GetMicroCandidate().GetSttHits() << std::endl;
217  std::cout<<"Kaon list 2, element "<<l<<" has no assosiated mcTRack"<<std::endl;
218  }
219  }
220  }
221 
222  n_removed=0;
223  ii=0;
224  for (l=0;l<k2.GetLength();++l) {
225  ii=l-n_removed;
226  if (k2[ii].GetMicroCandidate().GetMcIndex()>-1){
227  PndMCTrack *mcTrack = (PndMCTrack*)mc_array->At(k2[ii].GetMicroCandidate().GetMcIndex());
228  if (mcTrack!=0)
229  {
230  if ((mcTrack->GetPdgCode()!=-321) || (mcTrack->GetMotherID()!=-1))
231  {
232  k2.Remove(k2[ii]);
233  n_removed++;
234  }
235  }
236  else
237  {
238  std::cout<<"stt h: " << k2[ii].GetMicroCandidate().GetSttHits() << std::endl;
239  std::cout<<"Kaon list 2, element "<<l<<" has no assosiated mcTRack"<<std::endl;
240  }
241  }
242  }
243 
244  D1_pid.Combine(p1,p1, k2);
245  D2_pid.Combine(p2,p2, k1);
246 
247  for (j=0;j<D1_pid.GetLength();++j) h_mD1_pid->Fill(D1_pid[j].M());
248  for (j=0;j<D2_pid.GetLength();++j) h_mD2_pid->Fill(D2_pid[j].M());
249  D1_pid.Select(psiMassSel);
250  D2_pid.Select(psiMassSel);
251  psi.Combine(D1_pid,D2_pid);
252 
253  for (l=0;l<psi.GetLength();++l) {
254  h_psi_pid->Fill(psi[l].M());
255  }
256 
257  cout << " ";// << endl;
258  if (use4cfit)
259  {
261  int best_i=0;
262  double best_chi2=100000;
263  TCandidate *ccfit;// = new TCandidate();
264  double m_D1, m_D2;
265  for (l=0;l<psi.GetLength();++l) {
266  Rho4CFitter fitter(psi[l],ini);
267  fitter.FitConserveMasses();
268  double chi2=fitter.GetChi2();
269  if (chi2<best_chi2)
270  {
271  best_chi2=chi2;
272  best_i = l;
273  ccfit = (TCandidate*)fitter.FittedCand(psi[l]);
274  TCandidate *D1best = fitter.FittedCand(*(psi[l].Daughter(0)));
275  TCandidate *D2best = fitter.FittedCand(*(psi[l].Daughter(1)));
276  TCandidate *p1best = fitter.FittedCand(*(psi[l].Daughter(0)->Daughter(0)));
277  TCandidate *p2best = fitter.FittedCand(*(psi[l].Daughter(0)->Daughter(1)));
278  TCandidate *k1best = fitter.FittedCand(*(psi[l].Daughter(0)->Daughter(2)));
279  TCandidate *p3best = fitter.FittedCand(*(psi[l].Daughter(1)->Daughter(0)));
280  TCandidate *p4best = fitter.FittedCand(*(psi[l].Daughter(1)->Daughter(1)));
281  TCandidate *k2best = fitter.FittedCand(*(psi[l].Daughter(1)->Daughter(2)));
282  TLorentzVector tlvk1 = p1best->P4();
283  TLorentzVector tlvk2 = p2best->P4();
284  TLorentzVector tlvk3 = k1best->P4();
285  TLorentzVector tlvk4 = p3best->P4();
286  TLorentzVector tlvk5 = p4best->P4();
287  TLorentzVector tlvk6 = k2best->P4();
288  m_D1= (tlvk1+tlvk2+tlvk3).M();
289  m_D2= (tlvk4+tlvk5+tlvk6).M();
290  }
291  h_chi2_4c->Fill(chi2/9); // Ndf=3N-3=9
292  }
293 
294  if (psi.GetLength()!=0) cout << "evt: " << i << "\tbest: " << best_chi2 << "\tlen " << psi.GetLength() << endl;
295 
296  if(/*(best_chi2<270)&&*/(psi.GetLength()!=0))
297  {
298  h_chi2b_4c->Fill(best_chi2/9); // Ndf=3N-3=9
299  h_psi_4c->Fill(ccfit->M());
300  h_mD_4c->Fill(m_D1);
301  h_mD_4c->Fill(m_D2);
302  h_mD_final->Fill(m_D1);
303  h_mD_final->Fill(m_D2);
304  if (((m_D1>1.87-0.07)&&(m_D1<1.87+0.07))&&((m_D2>1.87-0.07)&&(m_D2<1.87+0.07)))
305  {
306  h_psi_Dmass->Fill(psi[best_i].M());
307  if ((psi[best_i].M()>3.7)&&(psi[best_i].M()<3.9))
308  n_reco++;
309  }
310  }
311  }
312  else // use vertex fit
313  {
314  TCandList psi_vtx, D1_vtx, D2_vtx;
315  TCandidate *ck1, *ck2, *cp1, *cp2, *cp3, *cp4, *D1tmp, *D2tmp, *D1_tmp, *D2_tmp, *psi_tmp;
316 
317  //Combine 4 kaons directly to candidates
318  for (j=0;j<psi.GetLength();++j)
319  {
320  D1tmp=psi[j].Daughter(0);
321  D2tmp=psi[j].Daughter(1);
322 
323  cp1=D1tmp->Daughter(0);
324  cp2=D1tmp->Daughter(1);
325  ck1=D1tmp->Daughter(2);
326 
327  cp3=D2tmp->Daughter(0);
328  cp4=D2tmp->Daughter(1);
329  ck2=D2tmp->Daughter(2);
330  D1_tmp=cp1->Combine(*cp2,*ck1);
331  D2_tmp=cp3->Combine(*cp4,*ck2);
332  //psi_tmp=cp1->Combine(*cp2,*ck1,*cp3, *cp4, *ck2);
333  D1_vtx.Add(*D1_tmp);
334  D2_vtx.Add(*D2_tmp);
335  //psi_vtx.Add(*psi_tmp);
336  }
337 
338  for (l=0;l<psi.GetLength();++l) {
339  // h_psi_pid_vtx->Fill(psi_vtx[l].M());
340  }
341 
342  {
344  int best_i=0;
345  double best_chi2=1000;
346  TCandidate *D1fit_best=0;
347  TCandidate *k1fit_best, *p1fit_best, *p2fit_best, *k1nofit, *p1nofit, *p2nofit;
348  TVector3 bestPos;
349  Float_t D1vtx_mass, pre_m_D1;
350  for (j=0;j<D1_vtx.GetLength();++j)
351  {
352  RhoKinVtxFitter vtxfitter(D1_vtx[j]); // instantiate a vertex fitter
353  vtxfitter.Fit(); // do the vertex fit
354 
355  TCandidate *D1fit=vtxfitter.FittedCand(D1_vtx[j]); // request the fitted Psi candidate
356  TVector3 D1Vtx=D1fit->Pos(); // and the decay vertex position
357  double chi2_vtx_D1=vtxfitter.GlobalChi2();
358  h_chi2_vtx_D1->Fill(chi2_vtx_D1/5); // Number degree of freedom 2N-3=5
359  // plot mass and vtx x,y projection after fit
360  hvpos_D1->Fill(D1Vtx.X(),D1Vtx.Y());
361  hvzpos_D1->Fill(D1Vtx.Z());
362  if(chi2_vtx_D1<best_chi2)
363  {
364  best_chi2=chi2_vtx_D1;
365  best_i=l;
366  D1fit_best=D1fit;
367  p1fit_best=vtxfitter.FittedCand(*(D1fit_best->Daughter(0)));
368  p2fit_best=vtxfitter.FittedCand(*(D1fit_best->Daughter(1)));
369  k1fit_best=vtxfitter.FittedCand(*(D1fit_best->Daughter(2)));
370  p1nofit=D1_vtx[j].Daughter(0);
371  p2nofit=D1_vtx[j].Daughter(1);
372  k1nofit=D1_vtx[j].Daughter(2);
373  TLorentzVector tlvp1 = p1nofit->P4();
374  TLorentzVector tlvp2 = p2nofit->P4();
375  TLorentzVector tlvk1 = k1nofit->P4();
376  pre_m_D1= (tlvp1+tlvp2+tlvk1).M();
377 
378  D1vtx_mass = D1fit_best->M();
379  bestPos = D1fit->Pos();
380  }
381  }
382  if((psi.GetLength()!=0))
383  {
384  h_D1_vtx->Fill(D1vtx_mass);
385  h_D1_vtx_bef->Fill(pre_m_D1);
386 
387  hvtxresX_D1->Fill(mcD1Vertex.X()-bestPos.X());
388  hvtxresY_D1->Fill(mcD1Vertex.Y()-bestPos.Y());
389  hvtxresZ_D1->Fill(mcD1Vertex.Z()-bestPos.Z());
390 
391  }
392  }
393 
394 
395  {
397  int best_i=0;
398  double best_chi2=1000;
399  TCandidate *D2fit_best=0;
400  TCandidate *k1fit_best, *p1fit_best, *p2fit_best, *k1nofit, *p1nofit, *p2nofit;
401  TVector3 bestPos;
402  Float_t D2vtx_mass, pre_m_D2;
403  for (j=0;j<D2_vtx.GetLength();++j)
404  {
405  RhoKinVtxFitter vtxfitter(D2_vtx[j]); // instantiate a vertex fitter
406  vtxfitter.Fit(); // do the vertex fit
407 
408  TCandidate *D2fit=vtxfitter.FittedCand(D2_vtx[j]); // request the fitted Psi candidate
409  TVector3 D2Vtx=D2fit->Pos(); // and the decay vertex position
410  double chi2_vtx_D2=vtxfitter.GlobalChi2();
411  h_chi2_vtx_D2->Fill(chi2_vtx_D2/5); // Number degree of freedom 2N-3=5
412  // plot mass and vtx x,y projection after fit
413  hvpos_D2->Fill(D2Vtx.X(),D2Vtx.Y());
414  hvzpos_D2->Fill(D2Vtx.Z());
415  if(chi2_vtx_D2<best_chi2)
416  {
417  best_chi2=chi2_vtx_D2;
418  best_i=l;
419  D2fit_best=D2fit;
420  p1fit_best=vtxfitter.FittedCand(*(D2fit_best->Daughter(0)));
421  p2fit_best=vtxfitter.FittedCand(*(D2fit_best->Daughter(1)));
422  k1fit_best=vtxfitter.FittedCand(*(D2fit_best->Daughter(2)));
423  p1nofit=D2_vtx[j].Daughter(0);
424  p2nofit=D2_vtx[j].Daughter(1);
425  k1nofit=D2_vtx[j].Daughter(2);
426  TLorentzVector tlvp1 = p1nofit->P4();
427  TLorentzVector tlvp2 = p2nofit->P4();
428  TLorentzVector tlvk1 = k1nofit->P4();
429  pre_m_D2= (tlvp1+tlvp2+tlvk1).M();
430 
431  D2vtx_mass = D2fit_best->M();
432  bestPos = D2fit->Pos();
433  }
434  }
435  if((psi.GetLength()!=0))
436  {
437  h_D2_vtx->Fill(D2vtx_mass);
438  h_D2_vtx_bef->Fill(pre_m_D2);
439  hvtxresX_D2->Fill(mcD2Vertex.X()-bestPos.X());
440  hvtxresY_D2->Fill(mcD2Vertex.Y()-bestPos.Y());
441  hvtxresZ_D2->Fill(mcD2Vertex.Z()-bestPos.Z());
442 
443  }
444  }
445  }
446 
447  }
448  std::cout<<"Number of reconstructed eta_c = "<<n_reco<<std::endl;
449  n_psi->SetBinContent(1,n_reco);
450 
451 
452  out->cd();
453 
454  h_mD1_nocuts->Write();
455  h_mD1_pid->Write();
456  h_mD2_nocuts->Write();
457  h_mD2_pid->Write();
458  h_psi_nocut->Write();
459  h_psi_pid->Write();
460 
461 
462  n_psi->Write();
463  n_events->Write();
464 
465  h_psi_Dmass->Write();
466  h_psi_vtx->Write();
467  h_D1_vtx_bef->Write();
468  h_D2_vtx_bef->Write();
469  h_D1_vtx->Write();
470  h_D2_vtx->Write();
471 
472  h_psi_4c->Write();
473 
474  //h_mD_vtx->Write();
475  h_mD_4c->Write();
476  h_mD_final->Write();
477 
478  nc->Write();
479 
480  h_chi2_4c->Write();
481  h_chi2b_4c->Write();
482  h_chi2_vtx_D1->Write();
483  hvzpos_D1->Write();
484  hvpos_D1->Write();
485  hvtxresX_D1->Write();
486  hvtxresY_D1->Write();
487  hvtxresZ_D1->Write();
488  hvzpos_D2->Write();
489  hvpos_D2->Write();
490  hvtxresX_D2->Write();
491  hvtxresY_D2->Write();
492  hvtxresZ_D2->Write();
493 
494 
495  out->Save();
496 
497  timer.Stop();
498  Double_t rtime = timer.RealTime();
499  Double_t ctime = timer.CpuTime();
500  printf("RealTime=%f seconds, CpuTime=%f seconds\n",rtime,ctime);
501 
502  return 0;
503 }
printf
printf("RealTime=%f seconds, CpuTime=%f seconds\n", rtime, ctime)
i
Int_t i
Definition: run_full.C:25
tree
TTree * tree
Definition: plot_dirc.C:12
nc
TH1F * nc
Definition: plot_eta_c.C:38
PndMCTrack::GetPdgCode
Int_t GetPdgCode() const
Definition: PndMCTrack.h:73
mc_array
TClonesArray * mc_array
Definition: anaLmdCluster.C:25
inFile
TString inFile
Definition: hit_dirc.C:8
h_chi2_4c
TH1F * h_chi2_4c
Definition: plot_eta_c.C:41
inSimFile
TString inSimFile
Definition: run_reco_llbar_evtgen.C:8
Double_t
Double_t
Definition: SimCompleteLinkDef.h:6
h_chi2b_4c
TH1F * h_chi2b_4c
Definition: plot_eta_c_stt.C:31
timer
TStopwatch timer
Definition: hit_dirc.C:51
out
TFile * out
Definition: reco_muo.C:20
ctime
Double_t ctime
Definition: hit_dirc.C:114
n_events
TVectorD * n_events
Definition: plot_eta_c.C:14
p2
TPad * p2
Definition: hist-t7.C:117
p1
TPad * p1
Definition: hist-t7.C:116
PndMCTrack::GetMotherID
Int_t GetMotherID() const
Definition: PndMCTrack.h:74
rtime
Double_t rtime
Definition: hit_dirc.C:113
void removeCombinatoric ( TCandList &  psi_list)

Definition at line 505 of file ana_check_psi.C.

References i.

506 {
507  TCandidate *phi1, *phi2, *k1, *k2, *k3, *k4;
508  bool isOverlap;
509 
510  int len1=psi_list.GetLength();
511  int ii=0;
512 
513  for (int i=0; i<len1; ++i)
514  {
515  phi1=psi_list[ii].Daughter(0);
516  phi2=psi_list[ii].Daughter(1);
517  k1=phi1->Daughter(0);
518  k2=phi1->Daughter(1);
519  k3=phi2->Daughter(0);
520  k4=phi2->Daughter(1);
521 
522  if ((k1->IsCloneOf(*k3))||(k2->IsCloneOf(*k4)))
523  {
524  psi_list.Remove(psi_list[ii]);
525  std::cout<<"Overlap found, i="<<i<<std::endl;
526  }
527  else
528  ii++;
529  }
530 
531 }
i
Int_t i
Definition: run_full.C:25
Generated on Wed Apr 3 2019 10:03:10 for FairRoot/PandaRoot by   doxygen 1.8.5