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