FairRoot/PandaRoot
Public Types | Public Member Functions | Protected Attributes | Private Member Functions | Private Attributes | List of all members
PndHypFullAna Class Reference

#include <PndHypFullAna.h>

Inheritance diagram for PndHypFullAna:

Public Types

typedef std::map< Int_t, Float_t > mapper
 

Public Member Functions

 PndHypFullAna ()
 
 ~PndHypFullAna ()
 
virtual InitStatus Init ()
 
virtual void Exec (Option_t *opt)
 
void Finish (TString cat)
 
void SetEnergySpectra (int event, int cluster)
 
void SetTotESpectra (int clus)
 
Int_t GetIonCharge (Int_t Z, Int_t &mass, Int_t &str)
 

Protected Attributes

int evcount
 
RhoNeutralParticleSelectorneutralSel
 
RhoPlusParticleSelectorplusSel
 
RhoMinusParticleSelectorminusSel
 
RhoMassParticleSelectorphiMSel
 
RhoMassParticleSelectorpi0MSel
 
RhoMassParticleSelectordsMSel
 
RhoMassParticleSelectorLambMSel
 
RhoSimpleKaonSelectorkSel
 
RhoSimplePionSelectorpiSel
 
RhoSimpleProtonSelectorpSel
 
TH2F * hvtx2 [10]
 
TH1F * spectra [10]
 
TH1F * ds0mass
 
TH1F * ximass
 
TH1F * Lamb
 
TH1F * ppi2mass
 
TH1F * ppi2
 
TH1F * e
 
TH2F * pid
 
TH2F * pidh
 
TH1F * nmult
 

Private Member Functions

virtual void SetParContainers ()
 
 ClassDef (PndHypFullAna, 1)
 

Private Attributes

TClonesArray * fChargedArray
 
TClonesArray * fMcTr
 
TClonesArray * fMicroArray
 
TClonesArray * fMcCands
 
TClonesArray * fMc
 
TClonesArray * fGe
 

Detailed Description

Definition at line 27 of file PndHypFullAna.h.

Member Typedef Documentation

typedef std::map<Int_t, Float_t> PndHypFullAna::mapper

Definition at line 31 of file PndHypFullAna.h.

Constructor & Destructor Documentation

PndHypFullAna::PndHypFullAna ( )

Default constructor

Definition at line 57 of file PndHypFullAna.cxx.

57  :
58  FairTask("Panda HypFullAna Task") {
59 }
PndHypFullAna::~PndHypFullAna ( )

Destructor

Definition at line 63 of file PndHypFullAna.cxx.

63 { }

Member Function Documentation

PndHypFullAna::ClassDef ( PndHypFullAna  ,
 
)
private
void PndHypFullAna::Exec ( Option_t *  opt)
virtual

Virtual method Exec

hvtx2[1]->Fill(pi3v.Mag(),pp3v.Mag());

Definition at line 178 of file PndHypFullAna.cxx.

References RhoCandList::Add(), RhoCandList::Cleanup(), RhoCandList::Combine(), evcount, fMc, fMcCands, fMcTr, fMicroArray, RhoCandList::Get(), PndHypPoint::GetEventID(), RhoCandList::GetLength(), PndMCTrack::GetMotherID(), PndPidCandidate::GetMvdHits(), PndMCTrack::GetPdgCode(), RhoCandidate::GetRecoCandidate(), PndMCTrack::GetStartVertex(), hit, hvtx2, RhoFactory::Instance(), jj, Lamb, LambMSel, RhoCandidate::Mass(), minusSel, MotherId, Motherpdg, RhoCandListIterator::Next(), npi, p, RhoCandidate::P4(), pi, piSel, plusSel, RhoCandidate::Pos(), ppi2, ppi2mass, pSel, RhoFactory::Reset(), RhoCandList::Select(), SetEnergySpectra(), SetTotESpectra(), t2, t3, and ximass.

178  {
179 
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;
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 }
Double_t Mass() const
void Add(const RhoCandidate *c)
Definition: RhoCandList.h:49
Double_t p
Definition: anasim.C:58
Int_t Motherpdg
RhoMinusParticleSelector * minusSel
Definition: PndHypFullAna.h:67
Int_t t2
Definition: hist-t7.C:106
void Cleanup()
Definition: RhoCandList.cxx:62
TClonesArray * fMc
Definition: PndHypFullAna.h:99
RhoPlusParticleSelector * plusSel
Definition: PndHypFullAna.h:66
Int_t GetPdgCode() const
Definition: PndMCTrack.h:73
Int_t GetLength() const
Definition: RhoCandList.h:46
static void Reset()
Definition: RhoFactory.cxx:28
RhoSimpleProtonSelector * pSel
Definition: PndHypFullAna.h:77
PndPidCandidate * GetRecoCandidate() const
Definition: RhoCandidate.h:376
#define pi
Definition: createSTT.C:60
RhoMassParticleSelector * LambMSel
Definition: PndHypFullAna.h:74
void SetEnergySpectra(int event, int cluster)
TH2F * hvtx2[10]
Definition: PndHypFullAna.h:80
void Combine(RhoCandList &l1, RhoCandList &l2)
TClonesArray * fMcTr
Definition: PndHypFullAna.h:96
TClonesArray * fMicroArray
Definition: PndHypFullAna.h:97
void Select(RhoParticleSelectorBase *pidmgr)
TLorentzVector P4() const
Definition: RhoCandidate.h:195
Int_t t3
Definition: hist-t7.C:106
static RhoFactory * Instance()
Definition: RhoFactory.cxx:34
Int_t GetMvdHits() const
Int_t GetEventID() const
Definition: PndHypPoint.h:90
void SetTotESpectra(int clus)
PndSdsMCPoint * hit
Definition: anasim.C:70
TVector3 GetStartVertex() const
Definition: PndMCTrack.h:76
Int_t GetMotherID() const
Definition: PndMCTrack.h:74
Int_t MotherId
RhoSimplePionSelector * piSel
Definition: PndHypFullAna.h:76
RhoCandidate * Get(Int_t)
Definition: RhoCandList.cxx:94
int npi
Definition: toy_core.C:124
TClonesArray * fMcCands
Definition: PndHypFullAna.h:98
void PndHypFullAna::Finish ( TString  cat)

Definition at line 822 of file PndHypFullAna.cxx.

References e, file, hvtx2, i, Lamb, pid, pidh, ppi2, ppi2mass, spectra, and ximass.

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 }
TH1F * spectra[10]
Definition: PndHypFullAna.h:81
Int_t i
Definition: run_full.C:25
TFile * file
TH2F * hvtx2[10]
Definition: PndHypFullAna.h:80
Int_t PndHypFullAna::GetIonCharge ( Int_t  Z,
Int_t &  mass,
Int_t &  str 
)

Definition at line 791 of file PndHypFullAna.cxx.

References L, and Z.

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 }
double Z
Definition: anaLmdDigi.C:68
InitStatus PndHypFullAna::Init ( )
virtual

Virtual method Init

Definition at line 69 of file PndHypFullAna.cxx.

References e, evcount, fChargedArray, fGe, fMc, fMcCands, fMcTr, fMicroArray, hvtx2, Lamb, LambMSel, minusSel, pid, pidh, piSel, plusSel, ppi2, ppi2mass, pSel, RhoParticleSelectorBase::SetCriterion(), spectra, and ximass.

69  {
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;
145 
146  // **** mass selectors for the resonances/composites
147  //
148 
150  piSel->SetCriterion("veryLoose");
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 }
TH1F * spectra[10]
Definition: PndHypFullAna.h:81
RhoMinusParticleSelector * minusSel
Definition: PndHypFullAna.h:67
TClonesArray * fMc
Definition: PndHypFullAna.h:99
RhoPlusParticleSelector * plusSel
Definition: PndHypFullAna.h:66
TClonesArray * fChargedArray
Definition: PndHypFullAna.h:95
RhoSimpleProtonSelector * pSel
Definition: PndHypFullAna.h:77
RhoMassParticleSelector * LambMSel
Definition: PndHypFullAna.h:74
TH2F * hvtx2[10]
Definition: PndHypFullAna.h:80
TClonesArray * fMcTr
Definition: PndHypFullAna.h:96
TClonesArray * fMicroArray
Definition: PndHypFullAna.h:97
virtual void SetCriterion(const char *crit)
TClonesArray * fGe
RhoSimplePionSelector * piSel
Definition: PndHypFullAna.h:76
TClonesArray * fMcCands
Definition: PndHypFullAna.h:98
void PndHypFullAna::SetEnergySpectra ( int  event,
int  cluster 
)

Definition at line 742 of file PndHypFullAna.cxx.

References fGe, PndHypGePoint::GetEnergyLoss(), PndHypGePoint::GetEventID(), and spectra.

Referenced by Exec().

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 }
TH1F * spectra[10]
Definition: PndHypFullAna.h:81
Int_t GetEventID() const
Definition: PndHypGePoint.h:52
TClonesArray * fGe
Double_t GetEnergyLoss() const
Definition: PndHypGePoint.h:62
void PndHypFullAna::SetParContainers ( )
privatevirtual

Geo file to use Get parameter containers

Definition at line 163 of file PndHypFullAna.cxx.

References run.

163  {
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 }
Int_t run
Definition: autocutx.C:47
void PndHypFullAna::SetTotESpectra ( int  clus)

Definition at line 770 of file PndHypFullAna.cxx.

References fGe, PndHypGePoint::GetEnergyLoss(), and spectra.

Referenced by Exec().

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 }
TH1F * spectra[10]
Definition: PndHypFullAna.h:81
TClonesArray * fGe
Double_t GetEnergyLoss() const
Definition: PndHypGePoint.h:62

Member Data Documentation

TH1F* PndHypFullAna::ds0mass
protected

Definition at line 83 of file PndHypFullAna.h.

RhoMassParticleSelector* PndHypFullAna::dsMSel
protected

Definition at line 73 of file PndHypFullAna.h.

TH1F* PndHypFullAna::e
protected

Definition at line 86 of file PndHypFullAna.h.

Referenced by Finish(), and Init().

int PndHypFullAna::evcount
protected

Definition at line 58 of file PndHypFullAna.h.

Referenced by Exec(), and Init().

TClonesArray* PndHypFullAna::fChargedArray
private

Input array of TpcLheTrack

Definition at line 95 of file PndHypFullAna.h.

Referenced by Init().

TClonesArray* PndHypFullAna::fGe
private

Definition at line 100 of file PndHypFullAna.h.

Referenced by Init(), SetEnergySpectra(), and SetTotESpectra().

TClonesArray* PndHypFullAna::fMc
private

Definition at line 99 of file PndHypFullAna.h.

Referenced by Exec(), and Init().

TClonesArray* PndHypFullAna::fMcCands
private

Definition at line 98 of file PndHypFullAna.h.

Referenced by Exec(), and Init().

TClonesArray* PndHypFullAna::fMcTr
private

Definition at line 96 of file PndHypFullAna.h.

Referenced by Exec(), and Init().

TClonesArray* PndHypFullAna::fMicroArray
private

Definition at line 97 of file PndHypFullAna.h.

Referenced by Exec(), and Init().

TH2F* PndHypFullAna::hvtx2[10]
protected

book all the histograms

Definition at line 80 of file PndHypFullAna.h.

Referenced by Exec(), Finish(), and Init().

RhoSimpleKaonSelector* PndHypFullAna::kSel
protected

Definition at line 75 of file PndHypFullAna.h.

TH1F* PndHypFullAna::Lamb
protected

Definition at line 84 of file PndHypFullAna.h.

Referenced by Exec(), Finish(), and Init().

RhoMassParticleSelector* PndHypFullAna::LambMSel
protected

Definition at line 74 of file PndHypFullAna.h.

Referenced by Exec(), and Init().

RhoMinusParticleSelector* PndHypFullAna::minusSel
protected

Definition at line 67 of file PndHypFullAna.h.

Referenced by Exec(), and Init().

RhoNeutralParticleSelector* PndHypFullAna::neutralSel
protected

Definition at line 65 of file PndHypFullAna.h.

TH1F* PndHypFullAna::nmult
protected

Definition at line 89 of file PndHypFullAna.h.

RhoMassParticleSelector* PndHypFullAna::phiMSel
protected

Definition at line 71 of file PndHypFullAna.h.

RhoMassParticleSelector* PndHypFullAna::pi0MSel
protected

Definition at line 72 of file PndHypFullAna.h.

TH2F* PndHypFullAna::pid
protected

Definition at line 87 of file PndHypFullAna.h.

Referenced by Finish(), and Init().

TH2F* PndHypFullAna::pidh
protected

Definition at line 87 of file PndHypFullAna.h.

Referenced by Finish(), and Init().

RhoSimplePionSelector* PndHypFullAna::piSel
protected

Definition at line 76 of file PndHypFullAna.h.

Referenced by Exec(), and Init().

RhoPlusParticleSelector* PndHypFullAna::plusSel
protected

Definition at line 66 of file PndHypFullAna.h.

Referenced by Exec(), and Init().

TH1F* PndHypFullAna::ppi2
protected

Definition at line 86 of file PndHypFullAna.h.

Referenced by Exec(), Finish(), and Init().

TH1F* PndHypFullAna::ppi2mass
protected

Definition at line 85 of file PndHypFullAna.h.

Referenced by Exec(), Finish(), and Init().

RhoSimpleProtonSelector* PndHypFullAna::pSel
protected

Definition at line 77 of file PndHypFullAna.h.

Referenced by Exec(), and Init().

TH1F* PndHypFullAna::spectra[10]
protected

Definition at line 81 of file PndHypFullAna.h.

Referenced by Finish(), Init(), SetEnergySpectra(), and SetTotESpectra().

TH1F* PndHypFullAna::ximass
protected

Definition at line 84 of file PndHypFullAna.h.

Referenced by Exec(), Finish(), and Init().


The documentation for this class was generated from the following files: