FairRoot/PandaRoot
PndHypFullAna.cxx
Go to the documentation of this file.
1 /******************************************************
2 
3 Analysis Task created by A.Sanchez
4 Read MicroCandidates and do analysis
5 for Hypernuclei.
6 *******************************************************/
7 
8 #include "TClonesArray.h"
9 
10 #include "FairRootManager.h"
11 #include "FairRunAna.h"
12 #include "FairRuntimeDb.h"
13 //#include "PndTpcLheGFTrack.h"
14 
15 //#include "PndTpcLheGFTrack.h"
16 #include "PndHypHit.h"
17 #include "PndHypPoint.h"
18 #include "../../hypGe/PndHypGePoint.h"
19 #include "PndStack.h"
20 #include "PndMCTrack.h"
21 
22 #include "TVector3.h"
23 #include "TH1F.h"
24 #include "TH2F.h"
25 
26 #include "FairRun.h"
27 #include "FairRuntimeDb.h"
28 #include "PndHypFullAna.h"
29 #include <string>
30 #include <iostream>
31 
32 #include "GFTrack.h"
33 #include "LSLTrackRep.h"
34 
35 //RHO stuff
36 #include "RhoCandidate.h"
37 #include "PndPidCandidate.h"
38 #include "PndPidCandidate.h"
39 #include "RhoCandList.h"
40 #include "RhoCandListIterator.h"
41 #include "RhoFactory.h"
42 #include "RhoSelector/RhoMassParticleSelector.h"
43 #include "RhoSelector/RhoPlusParticleSelector.h"
44 #include "RhoSelector/RhoMinusParticleSelector.h"
45 #include "RhoSelector/RhoSimpleElectronSelector.h"
46 #include "RhoSelector/RhoSimpleKaonSelector.h"
47 #include "RhoSelector/RhoSimpleMuonSelector.h"
48 #include "RhoSelector/RhoSimplePionSelector.h"
49 #include "RhoSelector/RhoSimpleProtonSelector.h"
50 
51 
52 using std::cout;
53 using std::endl;
54 
55 
56 // ----- Default constructor -------------------------------------------
57 PndHypFullAna::PndHypFullAna() :
58  FairTask("Panda HypFullAna Task") {
59 }
60 // -------------------------------------------------------------------------
61 
62 // ----- Destructor ----------------------------------------------------
63 PndHypFullAna::~PndHypFullAna () { }
64 // -------------------------------------------------------------------------
65 
66 
67 
68 // ----- Public method Init --------------------------------------------
69 InitStatus PndHypFullAna::Init() {
70 
71  //cout << " Inside the Init function****" << endl;
72 
73  //FairDetector::Initialize();
74  //FairRun* sim = FairRun::Instance();
75  //FairRuntimeDb* rtdb=sim->GetRuntimeDb();
76 
77  // Get RootManager
78  FairRootManager* ioman = FairRootManager::Instance();
79  if ( ! ioman ) {
80  cout << "-E- PndEmcHitProducer::Init: "
81  << "RootManager not instantiated!" << endl;
82  return kFATAL;
83  }
84 
85  // Get input array
86 
87  fMcTr = (TClonesArray*) ioman->GetObject("MCTrack");
88  fMcCands = (TClonesArray*) ioman->GetObject("HypHit");
89  fMc = (TClonesArray*) ioman->GetObject("HypPoint");
90  fGe = (TClonesArray*) ioman->GetObject("HypGePoint");
91 
92  fChargedArray = (TClonesArray*) ioman->GetObject("PndChargedCandidates");
93  fMicroArray = (TClonesArray*) ioman->GetObject("PndPidCandidates");
94 
95  if ( !fChargedArray && !fMicroArray) {
96  cout << "-W- PndHypFullAna ::Init: "
97  << "No PndChargedCandidates && PndNeutralCandidates array!" << endl;
98  return kERROR;
99  }
100 
101  // Create and register output array
102  cout << "-I- PndHypFullAna : Intialization successfull" << endl;
103 
104  ppi2mass = new TH1F("ppimass","p pi cands",200,0.05,0.3);
105  ppi2 = new TH1F("ppion","p pion cands",200,0.05,1.7);
106  e = new TH1F("evenet","p event",200,10000,50000);
107 
108  pid = new TH2F("pid","cands",200,0,16,200,0,16);
109  pidh = new TH2F("pidh","candsh",200,0,16,200,0,16);
110  ximass = new TH1F("ximass","xi cands",200,0.05,2);
111  Lamb = new TH1F("lamb mass","lamb cands",200,0.05,2);
112 
113  //itTop = new TH1F("p","cands",10000000,1010000000,1020000000);
114  spectra[0] = new TH1F("spectra01","cluster 1",200,0.0,0.01);
115  spectra[1] = new TH1F("spectra02","cluster 2",200,0.0,0.01);
116  spectra[2] = new TH1F("spectra03","cluster 3",200,0.0,0.01);
117  spectra[3] = new TH1F("spectra04","cluster 4",200,0.0,0.01);
118  spectra[4] = new TH1F("spectra05","cluster 5",200,0.0,0.01);
119  spectra[5] = new TH1F("spectra06","cluster 6",200,0.0,0.01);
120  spectra[6] = new TH1F("spectra07","cluster 7",200,0.0,0.01);
121  spectra[7] = new TH1F("spectra08","cluster 8",200,0.0,0.01);
122  spectra[8] = new TH1F("spectra09","cluster 9",200,0.0,0.01);
123  spectra[9] = new TH1F("spectra10","cluster 9",200,0.0,0.01);
124 
125 
126  //hvtx2[0]=new TH2F("hvtx201","vertex positions (x,z)",200,0.03,1.,200,0.03,1.);
127  hvtx2[0]=new TH2F("hvtx201","vertex positions (x,z)",100,0.03,0.15,100,0.03,0.15);
128  hvtx2[1]=new TH2F("hvtx202","vertex positions (x,z)",100,0.03,0.15,100,0.03,0.15);
129  hvtx2[2]=new TH2F("hvtx203","vertex positions (x,z)",100,0.03,0.15,100,0.03,0.15);
130  hvtx2[3]=new TH2F("hvtx204","vertex positions (x,z)",100,0.03,0.15,100,0.03,0.15);
131  hvtx2[4]=new TH2F("hvtx205","vertex positions (x,z)",100,0.03,0.15,100,0.03,0.15);
132  hvtx2[5]=new TH2F("hvtx206","vertex positions (x,z)",100,0.03,0.15,100,0.03,0.15);
133  hvtx2[6]=new TH2F("hvtx207","vertex positions (x,z)",100,0.03,0.15,100,0.03,0.15);
134  hvtx2[7]=new TH2F("hvtx208","vertex positions (x,z)",100,0.03,0.15,100,0.03,0.15);
135  hvtx2[8]=new TH2F("hvtx209","vertex positions (x,z)",100,0.03,0.15,100,0.03,0.15);
136  hvtx2[9]=new TH2F("hvtx210","vertex positions (x,z)",100,0.03,0.15,100,0.03,0.15);
137 
138 
139  // **** create and configure the selectors/filters we'd like to use later
140  //
141  //chargedSel = new RhoChargedParticleSelector;
142  //neutralSel = new RhoNeutralParticleSelector;
143  plusSel = new RhoPlusParticleSelector;
144  minusSel = new RhoMinusParticleSelector;
145 
146  // **** mass selectors for the resonances/composites
147  //
148 
149  piSel = new RhoSimplePionSelector();
150  piSel->SetCriterion("veryLoose");
151  pSel = new RhoSimpleProtonSelector();
152  pSel->SetCriterion("veryLoose");
153 
154  LambMSel = new RhoMassParticleSelector("LambSelector" , 1.115 , 0.04);
155 
156 
157  evcount=0;
158 
159  return kSUCCESS;
160 
161 }
162 
163 void PndHypFullAna::SetParContainers() {
164 
165  // Get run and runtime database
166  FairRunAna* run = FairRunAna::Instance();
167  if ( ! run ) Fatal("SetParContainers", "No analysis run");
168 
169  //FairRuntimeDb* db = run->GetRuntimeDb();
170  //if ( ! db ) Fatal("SetParContainers", "No runtime database");
171 
172 
173 }
174 
175 // -------------------------------------------------------------------------
176 
177 // ----- Public method Exec --------------------------------------------
178 void PndHypFullAna::Exec(Option_t*) {
179 
180  RhoFactory::Instance()->Reset();
181 
182  if (!(++evcount%100)) cout <<"evt "<<evcount<<endl;
183  //cout <<"evt "<<evcount<<endl;++evcount;
184 
185 
186  // **** create all the particle lists we'll need for rebuilding the decay tree
187  //
188  RhoCandList neutralCands,chargedCands, plusCands,minusCands;
189 
190  RhoCandList kpCands,kmCands,piCands,ppiCands;
191 
192  RhoCandList xiCands,nonOvCands,dsCands,ds0Cands,ppCands;
193  std::map<Int_t,Float_t > mapp;
194 
195  //RhoCandidate *tc;
196 
197 
198  // **** loop over all Candidates and add them to the list allCands
199  //
200  chargedCands.Cleanup();
201  //neutralCands.Cleanup();
202 
203  for (Int_t i1=0; i1<fMicroArray->GetEntriesFast(); i1++){
204  PndPidCandidate *mic = (PndPidCandidate *)fMicroArray->At(i1);
205  RhoCandidate tc(*mic,i1);
206  TLorentzVector l=tc.P4();
207  TVector3 p=tc.Pos();
208  //cout<<" micro to tcaaand "<<tc.GetCharge()<<endl;
209 
210  chargedCands.Add(&tc);
211  }
212 
213 
214 
215  //cout <<"c:"<<chargedCands.GetLength()<<" n:"<<neutralCands.GetLength()<<endl;
216 
217 
218  // **** select all the basic lists
219  //
220 
221 
222  plusCands.Select(chargedCands ,plusSel);
223  minusCands.Select(chargedCands ,minusSel);
224 
225  // **** pid selection
226  //
227 
228  piCands.Select(minusCands ,piSel);
229  ppCands.Select(plusCands ,pSel);
230 
231  ppiCands.Combine(minusCands,plusCands);
232  //ppiCands.Combine(pplusCands,piminusCands);
233  ppiCands.Select(LambMSel);
234 
235  for (int la=0;la<ppiCands.GetLength();la++)
236  {
237  Lamb->Fill((ppiCands.Get(la)->P4()).M());
238 
239  }
240 
241  xiCands.Combine(ppiCands,minusCands);
242 
243  RhoCandidate *t4;
244  RhoCandListIterator itX(xiCands);
245  while (t4=itX.Next())
246  {
247  ximass->Fill(t4->Mass());
248  }
249 
250  /* for (int c=0;c<chargedCands.GetLength();c++)
251  {
252  bool notinlist = true;
253  for (int d=0;d<xiCands.GetLength();d++)
254  {
255  if(chargedCands[c].Overlaps(xiCands.Get(d)))notinlist =false;
256  if(notinlist)nonOvCands.Append(chargedCands[c]);
257 
258  }
259  }
260  */
261 
262  //cout <<"charg non Overlap:"<<nonOvCands.GetLength()<<endl;
263 
264  //cout <<"pi-:"<<piCands.GetLength()<<endl;
265 
266 
267  // **** now start combining all composits; inbetween plot masses
268  // before using the mass selectors
269  //
270 
271  RhoCandidate *t2;
272 
273  RhoCandListIterator iterP(piCands);
274  while (t2=iterP.Next())
275  {
276  //ppimass->Fill(l.P());//tc->Mass());
277 
278  TLorentzVector l=t2->P4();
279  // ppimass->Fill(l.M());//tc->Mass());
280  TVector3 pim = l.Vect();
281 
282  //cout <<"fitted = "<<pim.Mag()<<" charge "<<t2->GetCharge()<<endl;
283  ppi2mass->Fill(pim.Mag());
284 
285 
286  }
287 
288  int ii,jj;
289 
290  int npi=piCands.GetLength();
291 
292  //int np=pplusCands.GetLength();
293  PndHypHit *hit;
294  PndHypPoint *po;
295  int Motherpdg,MotherId;
296 
297  //Calculate total energy SPectra
298  SetTotESpectra(9);
299 
300  for (int k=0;k<npi;k++)
301  {
302  RhoCandidate* pion=piCands.Get(k);
303  PndPidCandidate *pionpid =(PndPidCandidate *) pion->GetRecoCandidate();
304 
305  PndHypHit* hp=(PndHypHit*)fMcCands->At((pionpid->GetMvdHits())-1);
306  PndHypPoint* pop=(PndHypPoint*)fMc->At(hp->GetRefIndex());
307  if(pop==0)continue;
308 
309  PndMCTrack* moc=(PndMCTrack*)fMcTr->At(pop->GetTrackID());
310  if(moc==0)continue;
311  if( moc->GetPdgCode()==-211)dsCands.Add(pion);
312  //cout<<" number real "<<pop->GetEventID()<<" pion number real "<<dsCands.GetLength()<<endl;
313 
314  }
315  //cout<<" dsCands.GetLength() "<<dsCands.GetLength()<<endl;
316 
317  int dsi=dsCands.GetLength();
318  RhoCandidate *t3;
319  RhoCandListIterator itP(dsCands);
320  while (t3=itP.Next())
321  {
322  //ppimass->Fill(l.P());//tc->Mass());
323 
324  TLorentzVector v4=t3->P4();
325  // ppimass->Fill(l.M());//tc->Mass());
326  TVector3 v3 = v4.Vect();
327 
328  //cout <<"fitted = "<<pim.Mag()<<" charge "<<t2->GetCharge()<<endl;
329  ppi2->Fill(v3.Mag());
330 
331 
332  }
333 
334 /* for (int k=0;k<npi;k++)
335  {
336  RhoCandidate pion=piCands[k];
337 
338  TLorentzVector v4=pion.P4();
339  TVector3 v3 = v4.Vect();
340  //VAbsMicroCandidate cm;
341  //cm = pi.GetMicroCandidate();
342  PndHypHit* hp=(PndHypHit*)fMcCands->At((pion.GetMicroCandidate().GetMvdHits())-1);
343  PndHypPoint* pop=(PndHypPoint*)fMc->At(hp->GetRefIndex());
344  if(pop==0)continue;
345 
346  PndMCTrack* moc=(PndMCTrack*)fMcTr->At(pop->GetTrackID());
347  if(moc==0)continue;
348  MotherId= moc->GetMotherID();
349  if (MotherId==-1)Motherpdg = moc->GetPdgCode();
350  else {
351  PndMCTrack *mother =(PndMCTrack*)fMcTr->At(MotherId);
352  Motherpdg = mother->GetPdgCode();
353  }
354  //****cut on PCA to primary vertex has to be added ****
355  TVector3 vx=moc->GetStartVertex();
356 
357  //cout<<" event "<<pop->GetEventID()<<" "<<endl;
358  //****PID has to be added ****
359  //if(moc->GetPdgCode()!=-211)continue;
360  //if(vx.x()==0&&vx.y()==0&&vx.z()==-76.5)continue;
361  if((Motherpdg>1020000000||Motherpdg>1010000000)&&moc->GetPdgCode()==-211)ppi2->Fill(v3.Mag());
362 
363  if((Motherpdg>1020000000||Motherpdg>1010000000)&&v3.Mag()>0.5)
364  {
365  Int_t Z,A,L;
366 
367  cout<<" mala part "<<moc->GetPdgCode()<<" "<<Motherpdg<<endl;
368  cout<<" event "<<pop->GetEventID()<<" "<<endl;
369  Z= GetIonCharge(Motherpdg,A,L);
370  e->Fill(pop->GetEventID());
371 
372  if(L==1)pid->Fill(Z,A);
373  if(L==2)pidh->Fill(Z,A);
374  }
375 
376  }
377 */
378 
379  if(dsi==1){
380  //if(npi==1){
381  //RhoCandidate pi=piCands[0];
382  RhoCandidate* pi=dsCands.Get(0);
383  //RhoCandidate pp=piminusCands[jj];
384  TLorentzVector vpim=pi->P4();
385  TVector3 pi3v = vpim.Vect();
386  //cout<<ii<<" "<<pi3v.Mag()<<endl;
387 
388  if(pi3v.Mag()<0.2)hvtx2[0]->Fill(pi3v.Mag(),0.);
389 
390  }
391 
392 
393  if(dsi>1)
394  {
395  for (ii=0;ii<dsi-1;ii++)
396  {
397  //RhoCandidate pi=piCands[ii];
398  RhoCandidate *pi=dsCands.Get(ii);
399  //cout<<" charge pion "<<pi.GetCharge()<<endl;
400  PndPidCandidate* pipid = (PndPidCandidate*)pi->GetRecoCandidate();
401  TLorentzVector vpim=pi->P4();
402  TVector3 pi3v = vpim.Vect();
403 
404  for (jj=ii+1;jj<dsi;jj++)
405  {
406  //RhoCandidate pp=piCands[jj];
407  RhoCandidate *pp=dsCands.Get(jj);
408 
409  TLorentzVector vpp=pp->P4();
410  TVector3 pp3v = vpp.Vect();
411  //VAbsMicroCandidate cm;
412  //cm = pi.GetMicroCandidate();
413 
414  hit=(PndHypHit*)fMcCands->At((pipid->GetMvdHits())-1);
415  po=(PndHypPoint*)fMc->At(hit->GetRefIndex());
416  if(po==0)continue;
417 
418  PndMCTrack* mc=(PndMCTrack*)fMcTr->At(po->GetTrackID());
419  if(mc==0)continue;
420  MotherId= mc->GetMotherID();
421  if (MotherId==-1)Motherpdg = mc->GetPdgCode();
422  else {
423  PndMCTrack *mother =(PndMCTrack*)fMcTr->At(MotherId);
424  Motherpdg = mother->GetPdgCode();
425  }
426  //****cut on PCA to primary vertex has to be added ****
427  TVector3 vertex=mc->GetStartVertex();
428 
429  //if(mc->GetPdgCode()==3112) cout<<" Motherpdg mala "<<Motherpdg<<endl;
430  // cout<<" event "<<po->GetEventID()<<" "<<endl;
431  //****PID has to be added ****
432  //if(mc->GetPdgCode()!=-211)continue;
433 
434  //if(Motherpdg==-211||Motherpdg==310||
435  //if(Motherpdg==3312)continue;
436 
437  //if(vertex.x()==0&&vertex.y()==0&&vertex.z()==-76.5)continue;
438  cout<<" Motherpdg "<<Motherpdg<<" "<<mc->GetPdgCode()<<endl;
439  //if(pi3v.Mag()<0.2&&pp3v.Mag()<0.2)continue;
440 
441  // cout<<" pion selected for clu0 "<<mc->GetPdgCode()<<" mother "<<Motherpdg<<endl;
442  // cout<<" hits on hyp track "<<pi.GetMicroCandidate().GetMvdHits()<<" ev "<<po->GetEventID()<<endl;
443  //cout<<ii<<" "<<pi3v.Mag()<<" "<<jj<<" "<<pp3v.Mag()<<endl;
444 
445  //if(pi3v.Mag()>pp3v.Mag())cout<<" pi3v >pp3v "<<endl;
446  //if(pi3v.Mag()<pp3v.Mag())cout<<" pi3v <pp3v "<<endl;
447 
448  //if( pp.Overlaps(pim)) continue;
449 
450  //TMatrixD pcov=pp.Cov7();
451  //TMatrixD picov=pim.Cov7();
452  if(pi3v.Mag()>pp3v.Mag())
453  {
454 
455  hvtx2[0]->Fill(pi3v.Mag(),pp3v.Mag());
456 
457  if((pi3v.Mag()>0.12&&pi3v.Mag()<0.14)&&(pp3v.Mag()>0.065&&pp3v.Mag()<0.08)) {
458  //if(Motherpdg==1020040110||Motherpdg==1010050110){
459 
460  //cout<<" Be11LL "<<Motherpdg<<" "<<mc->GetPdgCode()<<endl;
461  //cout<<" pi3v >pp3v "<<po->GetEventID()<<endl;
462  hvtx2[1]->Fill(pi3v.Mag(),pp3v.Mag());
463  //if((pi3v.Mag()>0.12&&pi3v.Mag()<0.14)&&(pp3v.Mag()>0.065&&pp3v.Mag()<0.08)) {
464 
465  if(mapp[po->GetEventID()]==0) {
466  mapp[po->GetEventID()]=Motherpdg;
467  //cout<<" Motherpdg "<<Motherpdg<<endl;
468 
469  //
470  SetEnergySpectra(po->GetEventID(),0);
471  }
472 
473  }
474  if(Motherpdg==1020040110||Motherpdg==1010050110){
476  }
477 
478  if(Motherpdg==1020030090||(Motherpdg==1010040090 && MotherId!=-1)){
479 
480  //cout<<" pi3v >pp3v "<<po->GetEventID()<<endl;
481  //hvtx2[2]->Fill(pi3v.Mag(),pp3v.Mag());
482  //if((pi3v.Mag()>0.12&&pi3v.Mag()<0.14)&&(pp3v.Mag()>0.065&&pp3v.Mag()<0.08)) {
483  }
484 
485 
486  //cout<<" Motherpdg Li9LL "<<Motherpdg<<endl;
487 
488  //
489  if((pi3v.Mag()>0.112&&pi3v.Mag()<0.126)&&(pp3v.Mag()>0.09&&pp3v.Mag()<0.103)){
490  hvtx2[2]->Fill(pi3v.Mag(),pp3v.Mag());
491  if(mapp[po->GetEventID()]==0) {
492  mapp[po->GetEventID()]=Motherpdg;
493  SetEnergySpectra(po->GetEventID(),1);//cout<<" Motherpdg Li9LL cl1 "<<Motherpdg<<endl;
494  }
495  }
496 
497 
498  if((pi3v.Mag()>0.128&&pi3v.Mag()<0.147)&&(pp3v.Mag()>0.0898&&pp3v.Mag()<0.109)){
499  hvtx2[9]->Fill(pi3v.Mag(),pp3v.Mag());
500 
501  //cout<<" Motherpdg Li9LL cl2 "<<Motherpdg<<" "<<mc->GetPdgCode()<<endl;
502  if(mapp[po->GetEventID()]==0) {
503  mapp[po->GetEventID()]=Motherpdg;
504  SetEnergySpectra(po->GetEventID(),4);
505  //cout<<" Motherpdg Li9LL cl2 "<<Motherpdg<<endl;
506  }
507 
508  }
509 
510  //}
511 
512 
513  if(Motherpdg==1020040100||Motherpdg==1010050100){
514 
515  //cout<<" pi3v >pp3v "<<po->GetEventID()<<endl;
516  //hvtx2[3]->Fill(pi3v.Mag(),pp3v.Mag());
517  }
518 
519 
520  if((pi3v.Mag()>0.097&&pi3v.Mag()<0.106)&&(pp3v.Mag()>0.094&&pp3v.Mag()<0.103)) {
521  // if(Motherpdg==1020040100||Motherpdg==1010050100){
522 
523  //cout<<" pi3v >pp3v "<<po->GetEventID()<<endl;
524  hvtx2[3]->Fill(pi3v.Mag(),pp3v.Mag());
525  //if((pi3v.Mag()>0.12&&pi3v.Mag()<0.14)&&(pp3v.Mag()>0.065&&pp3v.Mag()<0.08)) {
526  //cout<<" Be10LL "<<Motherpdg<<" "<<mc->GetPdgCode()<<endl;
527 
528  if(mapp[po->GetEventID()]==0) {
529  mapp[po->GetEventID()]=Motherpdg;
530  // cout<<" Motherpdg "<<Motherpdg<<endl;
531 
532  //
533  SetEnergySpectra(po->GetEventID(),2);
534  }
535  }
536 
537 
538  if(Motherpdg==1020040120||Motherpdg==1010050120){
539 
540  //cout<<" pi3v >pp3v "<<po->GetEventID()<<endl;
541 
542  //hvtx2[4]->Fill(pi3v.Mag(),pp3v.Mag());
543  }
544 
545 
546  if((pi3v.Mag()>0.128&&pi3v.Mag()<0.147)&&(pp3v.Mag()>0.110&&pp3v.Mag()<0.124)) {
547  //if(Motherpdg==1020040120||Motherpdg==1010050120){
548 
549  //cout<<" pi3v >pp3v "<<po->GetEventID()<<endl;
550 
551  hvtx2[4]->Fill(pi3v.Mag(),pp3v.Mag());
552  // cout<<" Middle "<<Motherpdg<<" "<<mc->GetPdgCode()<<endl;
553  if(mapp[po->GetEventID()]==0) {
554  mapp[po->GetEventID()]=Motherpdg;
555  //cout<<" Mot be12LL "<<Motherpdg<<endl;
556 
557  //
558  SetEnergySpectra(po->GetEventID(),3);
559  }
560  }
561 
562  if(Motherpdg==1020020060||Motherpdg==1010030060){
563 
564  //cout<<" pi3v >pp3v "<<po->GetEventID()<<endl;
565 
566  //hvtx2[5]->Fill(pi3v.Mag(),pp3v.Mag());
567  }
568 
569 
570  if((pi3v.Mag()>0.128&&pi3v.Mag()<0.147)&&(pp3v.Mag()>0.124&&pp3v.Mag()<0.143)) {
571  //if(Motherpdg==1020020060||Motherpdg==1010030060){
572 
573  //cout<<" pi3v >pp3v "<<po->GetEventID()<<endl;
574 
575  hvtx2[5]->Fill(pi3v.Mag(),pp3v.Mag());
576 
577  //cout<<" Top "<<Motherpdg<<" "<<mc->GetPdgCode()<<endl;
578  if(mapp[po->GetEventID()]==0) {
579  mapp[po->GetEventID()]=Motherpdg;
580  // cout<<" Mot He6LL "<<Motherpdg<<endl;
581 
582  //
583  SetEnergySpectra(po->GetEventID(),5);
584  }
585  }
586 
587 
588  /*
589  if((pi3v.Mag()>0.095&&pi3v.Mag()<0.11)&&(pp3v.Mag()>0.09&&pp3v.Mag()<0.103)) SetEnergySpectra(po->GetEventID(),1);
590  if((pi3v.Mag()>0.112&&pi3v.Mag()<0.126)&&(pp3v.Mag()>0.09&&pp3v.Mag()<0.103)) SetEnergySpectra(po->GetEventID(),2);
591  if((pi3v.Mag()>0.126&&pi3v.Mag()<0.14)&&(pp3v.Mag()>0.09&&pp3v.Mag()<0.103)) SetEnergySpectra(po->GetEventID(),3);*/
592  }
593 
594  if(pi3v.Mag()<pp3v.Mag()) {
595  //cout<<" pi3v <pp3v "<<po->GetEventID()<<endl;
596 
597  hvtx2[0]->Fill(pp3v.Mag(),pi3v.Mag());
598 
599 
600  if(Motherpdg==1020040110||Motherpdg==1010050110){
601 
602  //hvtx2[1]->Fill(pp3v.Mag(),pi3v.Mag());
603  }
604 
605 
606  if((pp3v.Mag()>0.12&&pp3v.Mag()<0.14)&&(pi3v.Mag()>0.065&&pi3v.Mag()<0.08)) {
607  //if(Motherpdg==1020040110||Motherpdg==1010050110){
608 
609  hvtx2[1]->Fill(pp3v.Mag(),pi3v.Mag());
610  //
611  if(mapp[po->GetEventID()]==0) {
612  mapp[po->GetEventID()]=Motherpdg;
613  //cout<<" Motherpdg "<<Motherpdg<<endl;
614  //
615  SetEnergySpectra(po->GetEventID(),0);
616  }
617 
618  }
619 
620  if(Motherpdg==1020030090||(Motherpdg==1010040090&& MotherId!=-1)){
621 
622  //hvtx2[2]->Fill(pp3v.Mag(),pi3v.Mag());
623  }
624 
625  if((pp3v.Mag()>0.112&&pp3v.Mag()<0.126)&&(pi3v.Mag()>0.09&&pi3v.Mag()<0.103))
626  {
627  //if(Motherpdg==1020030090||(Motherpdg==1010040090&& MotherId!=-1)){
628 
629  hvtx2[2]->Fill(pp3v.Mag(),pi3v.Mag());
630  //if((pp3v.Mag()>0.12&&pp3v.Mag()<0.14)&&(pi3v.Mag()>0.065&&pi3v.Mag()<0.08)) {
631  if(mapp[po->GetEventID()]==0) {
632  mapp[po->GetEventID()]=Motherpdg;
633  //cout<<" Motherpdg "<<Motherpdg<<endl;
634  SetEnergySpectra(po->GetEventID(),1);//cout<<" Motherpdg Li9LL cl1 "<<Motherpdg<<endl;
635  }
636  }
637 
638 
639 
640  if((pp3v.Mag()>0.128&&pp3v.Mag()<0.147)&&(pi3v.Mag()>0.0898&&pi3v.Mag()<0.109))
641  {
642  hvtx2[9]->Fill(pp3v.Mag(),pi3v.Mag());
643  //cout<<" Motherpdg Li9LL cl2 "<<Motherpdg<<" "<<mc->GetPdgCode()<<endl;
644  if(mapp[po->GetEventID()]==0) {
645  mapp[po->GetEventID()]=Motherpdg;
646  SetEnergySpectra(po->GetEventID(),4);//
647  }
648 
649  }
650 
651  //}
652 
653  if(Motherpdg==1020040100||Motherpdg==1010050100){
654 
655  //hvtx2[3]->Fill(pp3v.Mag(),pi3v.Mag());
656  }
657 
658 
659  if((pp3v.Mag()>0.097&&pp3v.Mag()<0.106)&&(pi3v.Mag()>0.094&&pi3v.Mag()<0.103)) {
660  // if(Motherpdg==1020040100||Motherpdg==1010050100){
661 
662  hvtx2[3]->Fill(pp3v.Mag(),pi3v.Mag());
663  //
664  if(mapp[po->GetEventID()]==0) {
665  mapp[po->GetEventID()]=Motherpdg;
666  //cout<<" Motherpdg "<<Motherpdg<<endl;
667  //
668  SetEnergySpectra(po->GetEventID(),2);
669  }
670  }
671 
672  if(Motherpdg==1020040120||Motherpdg==1010050120){
673 
674  //hvtx2[4]->Fill(pp3v.Mag(),pi3v.Mag());
675  }
676 
677 
678  if((pp3v.Mag()>0.128&&pp3v.Mag()<0.147)&&(pi3v.Mag()>0.110&&pi3v.Mag()<0.124)) {
679  //if(Motherpdg==1020040120||Motherpdg==1010050120){
680 
681  hvtx2[4]->Fill(pp3v.Mag(),pi3v.Mag());
682  //cout<<" Middle "<<Motherpdg<<" "<<mc->GetPdgCode()<<endl;
683  //
684  if(mapp[po->GetEventID()]==0) {
685  mapp[po->GetEventID()]=Motherpdg;
686  //cout<<" Motherpdg "<<Motherpdg<<endl;
687  //
688  SetEnergySpectra(po->GetEventID(),3);
689  }
690  }
691 
692 
693  if(Motherpdg==1020020060||Motherpdg==1010030060){
694 
695  //cout<<" pi3v >pp3v "<<po->GetEventID()<<endl;
696 
697  //hvtx2[5]->Fill(pp3v.Mag(),pi3v.Mag());
698  }
699 
700 
701  if((pp3v.Mag()>0.128&&pp3v.Mag()<0.147)&&(pi3v.Mag()>0.124&&pi3v.Mag()<0.143)) {
702  // if(Motherpdg==1020020060||Motherpdg==1010030060){
703 
704  //cout<<" pi3v >pp3v "<<po->GetEventID()<<endl;
705 
706  hvtx2[5]->Fill(pp3v.Mag(),pi3v.Mag());
707  //cout<<" top "<<Motherpdg<<" "<<mc->GetPdgCode()<<endl;
708  //SetEnergySpectra(po->GetEventID(),5);
709  if(mapp[po->GetEventID()]==0) {
710  mapp[po->GetEventID()]=Motherpdg;
711  // // cout<<" Mot He6LL "<<Motherpdg<<endl;
712 
713  SetEnergySpectra(po->GetEventID(),5);
714  }
715  }
716 
717  /*
718  if((pp3v.Mag()>0.095&&pp3v.Mag()<0.11)&&(pi3v.Mag()>0.09&&pi3v.Mag()<0.103)) SetEnergySpectra(po->GetEventID(),1);
719  if((pp3v.Mag()>0.112&&pp3v.Mag()<0.126)&&(pi3v.Mag()>0.09&&pi3v.Mag()<0.103)) SetEnergySpectra(po->GetEventID(),2);
720  if((pp3v.Mag()>0.126&&pp3v.Mag()<0.14)&&(pi3v.Mag()>0.09&&pi3v.Mag()<0.103)) SetEnergySpectra(po->GetEventID(),3);*/
721  }
722 
723  //ppimass->Fill(compo->Mass());
724 
725  }
726  }
727  }
728 
729  // cout<<" pdg mother "<<mapp.size()<<endl;
730 
731  // for( std::map<Int_t, TH1F >::const_iterator ci=mapp.begin();
732  // ci != mapp.end();ci++){
733  // TH1F *h = itTop = ci->second;
734  // cout<<" pdg mother "<<ci->first<<" length "<<itTop.size()<<endl;
735 
736  // }
737 
738 
739 }
740 // -------------------------------------------------------------------------
741 
742 void PndHypFullAna::SetEnergySpectra(int event,int cluster)
743 {
744  Float_t E;
745  E=0;
746  if(fGe->GetEntriesFast()!=0)
747  {
748  //cout<<" ge entries "<<fGe->GetEntriesFast()<<endl;
749  PndHypGePoint *evCheck;
750  evCheck=(PndHypGePoint*)fGe->At(0);
751 
752  if(event == evCheck->GetEventID()){
753  cout<<" event "<<event<<" is correlated "<<evCheck->GetEventID()<<endl;
754 
755  for(int entries=0;entries<fGe->GetEntriesFast();entries++)
756  {
757 
758  PndHypGePoint *ge;
759  ge=(PndHypGePoint*)fGe->At(entries);
760 
761  E += ge->GetEnergyLoss();
762  }
763  }
764  }
765 
766 
767  if(E>0.&&E<0.010)spectra[cluster]->Fill(E);
768 }
769 
770 void PndHypFullAna::SetTotESpectra(int clus)
771 {
772  Float_t E;
773  E=0;
774  if(fGe->GetEntriesFast()!=0)
775  {
776  //cout<<" ge entries total "<<fGe->GetEntriesFast()<<endl;
777  for(int entries=0;entries<fGe->GetEntriesFast();entries++)
778  {
779 
780  PndHypGePoint *ge;
781  ge=(PndHypGePoint*)fGe->At(entries);
782 
783  E += ge->GetEnergyLoss();
784  }
785  }
786 
787 
788 
789  if(E>0.&&E<0.010)spectra[clus]->Fill(E);
790 }
791 Int_t PndHypFullAna::GetIonCharge(Int_t ion,Int_t &mass,Int_t &str)
792 {
793  Int_t A,Z,L;
794 
795  if(ion>1000000000&&(ion<1010000000))
796  { ion -= 1000000000;
797  Z = ion/10000;
798  ion -= 10000*Z;
799  A = ion/10;
800  //cout<<" ion charge "<<Z<<endl;
801  mass = A;
802  str =0;
803 
804  return Z;
805 
806  }
807 if((ion>1010000000||ion>1020000000))
808  {
809  ion -= 1000000000;
810  L = ion/10000000;
811  ion -= 10000000*L;
812  Z = ion/10000;
813  ion -= 10000*Z;
814  A = ion/10;
815  //cout<<L<<" hypernuclei charge "<<Z<<endl;
816  mass = A;
817  str =L;
818  return Z;
819 
820  }
821 }
822 void PndHypFullAna::Finish(TString cat)
823 {
824  TFile* file = FairRootManager::Instance()->GetOutFile();
825  file->cd();
826  //file->mkdir(cat.Data());//"HypHitAnaF");
827  //file->cd(cat.Data());//"HypHitAnaF");
828 
829  ppi2mass->Write();
830  delete ppi2mass;
831  ppi2mass=NULL;
832 
833  ximass->Write();
834  delete ximass;
835  ximass=NULL;
836 
837  Lamb->Write();
838  delete Lamb;
839  Lamb=NULL;
840 
841  ppi2->Write();
842  delete ppi2;
843  ppi2=NULL;
844 
845  e->Write();
846  delete e;
847  e=NULL;
848 
849  pid->Write();
850  delete pid;
851  pid=NULL;
852  pidh->Write();
853  delete pidh;
854  pidh=NULL;
855 
856  for(int i=0;i<10;i++){
857  hvtx2[i]->Write();
858  delete hvtx2[i];
859  hvtx2[i] =NULL;
860 
861  spectra[i]->Write();
862  delete spectra[i];
863  spectra[i]=NULL;
864  }
865 
866  /*
867  hvtx2[0]->Write();
868  hvtx2[1]->Write();
869  hvtx2[2]->Write();
870  hvtx2[3]->Write();
871  hvtx2[4]->Write();
872  hvtx2[5]->Write(); hvtx2[6]->Write(); hvtx2[7]->Write();hvtx2[8]->Write();
873  // hvtx2[9]->Write();
874 
875  ppi2mass->Write();
876  spectra[0]->Write();
877  spectra[1]->Write();
878  spectra[2]->Write();
879  spectra[3]->Write();
880  spectra[4]->Write();
881  spectra[5]->Write();spectra[6]->Write();spectra[7]->Write();
882  spectra[8]->Write(); //spectra[9]->Write();
883  */
884 
885 
886 
887 
888 
889 }
890 
891 ClassImp(PndHypFullAna)
RhoCandidate::Mass
Double_t Mass() const
Definition: RhoCandidate.cxx:981
PndHypFullAna::spectra
TH1F * spectra[10]
Definition: PndHypFullAna.h:81
RhoCandList::Add
void Add(const RhoCandidate *c)
Definition: RhoCandList.h:49
Motherpdg
Int_t Motherpdg
Definition: GammaSpectraAnalysis.C:39
PndHypFullAna::minusSel
RhoMinusParticleSelector * minusSel
Definition: PndHypFullAna.h:67
run
Int_t run
Definition: autocutx.C:47
RhoCandList::Cleanup
void Cleanup()
Definition: RhoCandList.cxx:62
PndHypFullAna::fMc
TClonesArray * fMc
Definition: PndHypFullAna.h:99
i
Int_t i
Definition: run_full.C:25
file
TFile * file
Definition: drc/old/draw_geom.C:12
PndHypFullAna::plusSel
RhoPlusParticleSelector * plusSel
Definition: PndHypFullAna.h:66
PndMCTrack::GetPdgCode
Int_t GetPdgCode() const
Definition: PndMCTrack.h:73
RhoCandidate::Pos
TVector3 Pos() const
Definition: RhoCandidate.h:186
PndHypFullAna::fChargedArray
TClonesArray * fChargedArray
Definition: PndHypFullAna.h:95
RhoCandList::GetLength
Int_t GetLength() const
Definition: RhoCandList.h:46
RhoFactory::Reset
static void Reset()
Definition: RhoFactory.cxx:28
PndHypFullAna::pSel
RhoSimpleProtonSelector * pSel
Definition: PndHypFullAna.h:77
RhoCandidate::GetRecoCandidate
PndPidCandidate * GetRecoCandidate() const
Definition: RhoCandidate.h:376
pi
#define pi
Definition: createSTT.C:60
PndHypFullAna::SetParContainers
virtual void SetParContainers()
Definition: PndHypFullAna.cxx:163
jj
Int_t jj
Definition: createRootGeoFileBwEnd_FwEndCapGeo.C:434
p
Double_t p
Definition: anasim.C:58
PndHypFullAna::LambMSel
RhoMassParticleSelector * LambMSel
Definition: PndHypFullAna.h:74
PndHypFullAna::GetIonCharge
Int_t GetIonCharge(Int_t Z, Int_t &mass, Int_t &str)
Definition: PndHypFullAna.cxx:791
PndHypFullAna::SetEnergySpectra
void SetEnergySpectra(int event, int cluster)
Definition: PndHypFullAna.cxx:742
PndHypFullAna::hvtx2
TH2F * hvtx2[10]
Definition: PndHypFullAna.h:80
RhoCandList::Combine
void Combine(RhoCandList &l1, RhoCandList &l2)
Definition: RhoCandList.cxx:289
PndHypFullAna::fMcTr
TClonesArray * fMcTr
Definition: PndHypFullAna.h:96
PndHypFullAna::fMicroArray
TClonesArray * fMicroArray
Definition: PndHypFullAna.h:97
RhoCandList::Select
void Select(RhoParticleSelectorBase *pidmgr)
Definition: RhoCandList.cxx:1359
RhoCandidate::P4
TLorentzVector P4() const
Definition: RhoCandidate.h:195
t3
Int_t t3
Definition: hist-t7.C:106
RhoFactory::Instance
static RhoFactory * Instance()
Definition: RhoFactory.cxx:34
PndPidCandidate::GetMvdHits
Int_t GetMvdHits() const
Definition: PndPidCandidate.h:101
PndHypGePoint::GetEventID
Int_t GetEventID() const
Definition: PndHypGePoint.h:52
PndHypFullAna::Exec
virtual void Exec(Option_t *opt)
Definition: PndHypFullAna.cxx:178
RhoParticleSelectorBase::SetCriterion
virtual void SetCriterion(const char *crit)
Definition: RhoParticleSelectorBase.cxx:61
hit
int hit(Int_t nEvents=0, TString inFile="sim.root", TString parFile="par.root", TString inDigi="digi.root", TString outFile="hit.root", Int_t timeBased=0)
Definition: hit.C:1
t2
Int_t t2
Definition: hist-t7.C:106
PndHypPoint::GetEventID
Int_t GetEventID() const
Definition: PndHypPoint.h:90
PndHypFullAna::SetTotESpectra
void SetTotESpectra(int clus)
Definition: PndHypFullAna.cxx:770
PndHypFullAna::Finish
void Finish(TString cat)
Definition: PndHypFullAna.cxx:822
ClassImp
ClassImp(PndAnaContFact)
Z
double Z
Definition: anaLmdDigi.C:68
PndHypFullAna::fGe
TClonesArray * fGe
Definition: PndHypFullAna.h:100
program_options
PndMCTrack::GetStartVertex
TVector3 GetStartVertex() const
Definition: PndMCTrack.h:76
PndMCTrack::GetMotherID
Int_t GetMotherID() const
Definition: PndMCTrack.h:74
MotherId
Int_t MotherId
Definition: GammaSpectraAnalysis.C:39
PndHypFullAna::piSel
RhoSimplePionSelector * piSel
Definition: PndHypFullAna.h:76
PndHypGePoint::GetEnergyLoss
Double_t GetEnergyLoss() const
Definition: PndHypGePoint.h:62
RhoCandList::Get
RhoCandidate * Get(Int_t)
Definition: RhoCandList.cxx:94
npi
int npi
Definition: toy_core.C:124
PndHypFullAna::fMcCands
TClonesArray * fMcCands
Definition: PndHypFullAna.h:98
PndHypFullAna::Init
virtual InitStatus Init()
Definition: PndHypFullAna.cxx:69
Generated on Sun Aug 18 2019 03:01:31 for FairRoot/PandaRoot by   doxygen 1.8.8