CbmRoot
CbmKresSelectGoodEvents.cxx
Go to the documentation of this file.
1 
15 
16 #include "CbmMCTrack.h"
17 
18 #include "FairLogger.h"
19 #include "FairRootManager.h"
20 #include "FairRunSim.h"
21 
22 #include <iostream>
23 
24 using namespace std;
25 
27  : FairTask(), fMcTracks(nullptr), fApp(nullptr) {}
28 
30 
32 
33  FairRunSim* sim = FairRunSim::Instance();
34  if (sim) { fApp = FairMCApplication::Instance(); }
35 
36  FairRootManager* ioman = FairRootManager::Instance();
37  if (nullptr == ioman) {
38  Fatal("CbmKresEta::Init", "RootManager not instantised!");
39  }
40 
41  fMcTracks = (TClonesArray*) ioman->GetObject("MCTrack");
42  if (nullptr == fMcTracks) {
43  Fatal("CbmKresSelectGoodEvents::Init", "No MCTrack array!");
44  }
45 
46  return kSUCCESS;
47 }
48 
49 
51  Electrons.clear();
52  Int_t nofMcTracks = fMcTracks->GetEntriesFast();
53  for (int i = 0; i < nofMcTracks; i++) {
54  CbmMCTrack* mctrack = (CbmMCTrack*) fMcTracks->At(i);
55  if (mctrack == nullptr) continue;
56  if (mctrack->GetMotherId() == -1) continue;
57  CbmMCTrack* mcMotherTrack =
58  (CbmMCTrack*) fMcTracks->At(mctrack->GetMotherId());
59  if (mcMotherTrack == nullptr) continue;
60 
61  if (TMath::Abs(mctrack->GetPdgCode()) == 11
62  && mcMotherTrack->GetPdgCode() == 22) {
63  if (mcMotherTrack->GetMotherId() == -1) continue;
64  CbmMCTrack* mcGrTrack =
65  (CbmMCTrack*) fMcTracks->At(mcMotherTrack->GetMotherId());
66  if (mcGrTrack == nullptr) continue;
67  if (mcGrTrack->GetPdgCode() == 221) { Electrons.push_back(mctrack); }
68  }
69  }
70 
71  int EtaConversion = 0;
72  if (Electrons.size() >= 4) {
73  for (size_t i = 0; i < Electrons.size(); i++) {
74  for (size_t j = i + 1; j < Electrons.size(); j++) {
75  for (size_t k = j + 1; k < Electrons.size(); k++) {
76  for (size_t l = k + 1; l < Electrons.size(); l++) {
77 
78  int pdg1 = Electrons.at(i)->GetPdgCode();
79  int pdg2 = Electrons.at(j)->GetPdgCode();
80  int pdg3 = Electrons.at(k)->GetPdgCode();
81  int pdg4 = Electrons.at(l)->GetPdgCode();
82 
83  if (pdg1 + pdg2 + pdg3 + pdg4 != 0) continue;
84  if (TMath::Abs(pdg1) != 11 || TMath::Abs(pdg2) != 11
85  || TMath::Abs(pdg3) != 11 || TMath::Abs(pdg4) != 11)
86  continue;
87 
88  int motherId1 = Electrons.at(i)->GetMotherId();
89  int motherId2 = Electrons.at(j)->GetMotherId();
90  int motherId3 = Electrons.at(k)->GetMotherId();
91  int motherId4 = Electrons.at(l)->GetMotherId();
92 
93  if (motherId1 == -1 || motherId2 == -1 || motherId3 == -1
94  || motherId4 == -1)
95  continue;
96 
97  CbmMCTrack* mother1 = (CbmMCTrack*) fMcTracks->At(motherId1);
98  CbmMCTrack* mother2 = (CbmMCTrack*) fMcTracks->At(motherId2);
99  CbmMCTrack* mother3 = (CbmMCTrack*) fMcTracks->At(motherId3);
100  CbmMCTrack* mother4 = (CbmMCTrack*) fMcTracks->At(motherId4);
101 
102  int mcMotherPdg1 = mother1->GetPdgCode();
103  int mcMotherPdg2 = mother2->GetPdgCode();
104  int mcMotherPdg3 = mother3->GetPdgCode();
105  int mcMotherPdg4 = mother4->GetPdgCode();
106 
107  if (mcMotherPdg1 != 22 || mcMotherPdg2 != 22 || mcMotherPdg3 != 22
108  || mcMotherPdg4 != 22)
109  continue;
110 
111  int grandmotherId1 = mother1->GetMotherId();
112  int grandmotherId2 = mother2->GetMotherId();
113  int grandmotherId3 = mother3->GetMotherId();
114  int grandmotherId4 = mother4->GetMotherId();
115 
116  if (grandmotherId1 == -1) continue;
117  CbmMCTrack* GrTrack = (CbmMCTrack*) fMcTracks->At(grandmotherId1);
118 
119  if (grandmotherId1 == grandmotherId2
120  && grandmotherId1 == grandmotherId3
121  && grandmotherId1 == grandmotherId4
122  && GrTrack->GetPdgCode() == 221) {
123  EtaConversion++;
124  cout << "Decay eta -> gamma gamma -> e+e- e+e- detected!\t\t mc "
125  "mass: "
126  << GrTrack->GetMass() << endl;
127  cout << "motherids: " << motherId1 << "/" << motherId2 << "/"
128  << motherId3 << "/" << motherId4 << endl;
129  cout << "grandmotherid: " << grandmotherId1 << "/"
130  << grandmotherId2 << "/" << grandmotherId3 << "/"
131  << grandmotherId4 << endl;
132  }
133  }
134  }
135  }
136  }
137  }
138 
139  cout << "CbmKresSelectGoodEvents, EtaConversion = " << EtaConversion << endl;
140 
141  // if (fApp && EtaConversion == 0) {
142  // LOG(WARNING) << "No double converted Eta";
143  // fApp->SetSaveCurrentEvent(kFALSE);
144  // }
145 }
146 
147 
CbmMCTrack::GetMotherId
Int_t GetMotherId() const
Definition: CbmMCTrack.h:71
CbmMCTrack::GetMass
Double_t GetMass() const
Mass of the associated particle.
Definition: CbmMCTrack.cxx:114
CbmKresSelectGoodEvents::CbmKresSelectGoodEvents
CbmKresSelectGoodEvents()
Definition: CbmKresSelectGoodEvents.cxx:26
CbmKresSelectGoodEvents::Exec
virtual void Exec(Option_t *)
Definition: CbmKresSelectGoodEvents.cxx:50
CbmKresSelectGoodEvents::~CbmKresSelectGoodEvents
virtual ~CbmKresSelectGoodEvents()
Definition: CbmKresSelectGoodEvents.cxx:29
CbmMCTrack::GetPdgCode
Int_t GetPdgCode() const
Definition: CbmMCTrack.h:70
i
int i
Definition: L1/vectors/P4_F32vec4.h:25
CbmKresSelectGoodEvents::Electrons
std::vector< CbmMCTrack * > Electrons
Definition: CbmKresSelectGoodEvents.h:30
CbmKresSelectGoodEvents::Finish
virtual void Finish()
Definition: CbmKresSelectGoodEvents.cxx:148
CbmKresSelectGoodEvents::fMcTracks
TClonesArray * fMcTracks
Definition: CbmKresSelectGoodEvents.h:28
CbmMCTrack.h
CbmMCTrack
Definition: CbmMCTrack.h:34
CbmKresSelectGoodEvents.h
CbmKresSelectGoodEvents::Init
virtual InitStatus Init()
Definition: CbmKresSelectGoodEvents.cxx:31
CbmKresSelectGoodEvents::fApp
FairMCApplication * fApp
Definition: CbmKresSelectGoodEvents.h:29