4747a8b3d00b4fb557c031ea42dc41371f08daaa
[u/mrichter/AliRoot.git] / PWG4 / totEt / AliAnalysisHadEtMonteCarlo.cxx
1 //_________________________________________________________________________
2 //  Utility Class for transverse energy studies, charged hadrons
3 //  Base class for MC analysis
4 //  - MC output
5 // implementation file
6 //
7 //Created by Christine Nattrass, Rebecca Scott, Irakli Martashvili
8 //University of Tennessee at Knoxville
9 //_________________________________________________________________________
10 #include "AliAnalysisHadEtMonteCarlo.h"
11 #include "AliAnalysisEtCuts.h"
12
13 #include "AliStack.h"
14 #include "AliMCEvent.h"
15 #include "AliESDEvent.h"
16 #include "AliESDtrackCuts.h"
17 #include "AliESDpid.h"
18 #include <iostream>
19
20 using namespace std;
21
22 ClassImp(AliAnalysisHadEtMonteCarlo);
23
24
25 Int_t AliAnalysisHadEtMonteCarlo::AnalyseEvent(AliVEvent* ev,AliVEvent* ev2)
26 { // analyse MC and real event info
27   FillHisto1D("NEvents",0.5,1);
28
29   AnalyseEvent(ev);
30   AliMCEvent *mcEvent = dynamic_cast<AliMCEvent*>(ev);
31   AliESDEvent *realEvent = dynamic_cast<AliESDEvent*>(ev2);
32   AliStack *stack = mcEvent->Stack();
33
34   //for PID
35   AliESDpid *pID = new AliESDpid();
36   pID->MakePID(realEvent);
37
38   //This code taken from https://twiki.cern.ch/twiki/bin/view/ALICE/SelectionOfPrimaryTracksForPp2009DataAnalysis
39   //Gets good tracks
40   //=============================================
41   // Primary vertex
42   const AliESDVertex *vertex = realEvent->GetPrimaryVertexTracks();
43   if(vertex->GetNContributors()<1) {
44     // SPD vertex
45     vertex = realEvent->GetPrimaryVertexSPD();
46     if(vertex->GetNContributors()<1) {
47       // NO GOOD VERTEX, SKIP EVENT 
48     }
49   }
50   // apply a cut |zVertex| < CUT, if needed
51
52   
53
54   //fEsdtrackCutsITSTPC->SetEtaRange(-0.8,0.8); // normally, |eta|<0.8
55   //=============================================
56
57   //Roughly following $ALICE_ROOT/PWG0/dNdEta/AlidNdEtaCorrectionTask
58
59   //=============================================TPC&&ITS=============================================
60   TString *strTPC = new TString("TPC");
61   TString *strITS = new TString("ITS");
62   TString *strTPCITS = new TString("TPCITS");
63   for(Int_t cutset=0;cutset<3;cutset++){
64     TString *cutName;
65     TObjArray* list;
66     switch(cutset){
67     case 0:
68       cutName = strTPC;
69       list = fEsdtrackCutsTPC->GetAcceptedTracks(realEvent);
70       break;
71     case 1:
72       cutName = strITS;
73       list = fEsdtrackCutsITS->GetAcceptedTracks(realEvent);
74       break;
75     case 2:
76       cutName = strTPCITS;
77       list = fEsdtrackCutsITSTPC->GetAcceptedTracks(realEvent);
78       break;
79     default:
80       cerr<<"Error:  cannot fill histograms!"<<endl;
81       return -1;
82     }
83     Int_t nGoodTracks = list->GetEntries();
84     for (Int_t iTrack = 0; iTrack < nGoodTracks; iTrack++)
85       {
86         AliESDtrack *track = dynamic_cast<AliESDtrack*> (list->At(iTrack));
87         if (!track)
88           {
89             Printf("ERROR: Could not get track %d", iTrack);
90             continue;
91           }
92         else{
93           Float_t nSigmaPion,nSigmaProton,nSigmaKaon,nSigmaElectron;
94           if(cutset!=1){
95             nSigmaPion = TMath::Abs(pID->NumberOfSigmasTPC(track,AliPID::kPion));
96             nSigmaProton = TMath::Abs(pID->NumberOfSigmasTPC(track,AliPID::kProton));
97             nSigmaKaon = TMath::Abs(pID->NumberOfSigmasTPC(track,AliPID::kKaon));
98             nSigmaElectron = TMath::Abs(pID->NumberOfSigmasTPC(track,AliPID::kElectron));
99           }
100           else{
101             nSigmaPion = TMath::Abs(pID->NumberOfSigmasITS(track,AliPID::kPion));
102             nSigmaProton = TMath::Abs(pID->NumberOfSigmasITS(track,AliPID::kProton));
103             nSigmaKaon = TMath::Abs(pID->NumberOfSigmasITS(track,AliPID::kKaon));
104             nSigmaElectron = TMath::Abs(pID->NumberOfSigmasITS(track,AliPID::kElectron));
105           }
106           bool isPion = (nSigmaPion<3.0 && nSigmaProton>2.0 && nSigmaKaon>2.0);
107           bool isElectron = (nSigmaElectron<2.0 && nSigmaPion>4.0 && nSigmaProton>3.0 && nSigmaKaon>3.0);
108           bool isKaon = (nSigmaPion>3.0 && nSigmaProton>2.0 && nSigmaKaon<2.0);
109           bool isProton = (nSigmaPion>3.0 && nSigmaProton<2.0 && nSigmaKaon>2.0);
110
111           //bool IsElectron = false;
112           bool unidentified = (!isProton && !isKaon && !isElectron);
113           Float_t dEdx = track->GetTPCsignal();
114           if(cutset==1) dEdx = track->GetITSsignal();
115
116           FillHisto2D(Form("dEdxAll%s",cutName->Data()),track->P(),dEdx,1.0);
117           //if(cutset==1) cout<<"filling "<<track->P()<<" "<<dEdx<<endl;
118
119           UInt_t label = (UInt_t)TMath::Abs(track->GetLabel());
120           TParticle  *simPart  = stack->Particle(label);
121           if(!simPart) {
122             Printf("no MC particle\n");                 
123             continue;           
124           }
125           else{//analysis
126             if(stack->IsPhysicalPrimary(label)){
127               if (TMath::Abs(simPart->Eta()) < fCuts->GetCommonEtaCut())            {
128
129                 Int_t pdgCode =  simPart->GetPDG(0)->PdgCode();
130                 Int_t mypid = 0;
131                 if(pdgCode==fPiPlusCode) mypid = 1;
132                 if(pdgCode==fProtonCode) mypid = 2;
133                 if(pdgCode==fKPlusCode) mypid = 3;
134                 if(pdgCode==fEPlusCode) mypid = 4;
135                 if(pdgCode==fPiMinusCode) mypid = 1;
136                 if(pdgCode==fAntiProtonCode) mypid = 2;
137                 if(pdgCode==fKMinusCode) mypid = 3;
138                 if(pdgCode==fEMinusCode) mypid = 4;
139                 //cout<<pdgCode->PdgCode()<<" ";
140                 //fPdgDB->GetSimParticle("pi+")->PdgCode();
141                 bool filled = false;      
142                 //============Charged hadrons===================================
143                 //identified...
144                 if(isPion){
145                   if(pdgCode!=fPiPlusCode && pdgCode!=fPiMinusCode){
146                     FillHisto2D(Form("MisidentifiedPIDs%s",cutName->Data()),1,mypid,1);
147                     //if(mypid==0)cerr<<"I was misidentified! I'm not a pion! I am a "<<simPart->GetName()<<endl;
148                   }
149                   float myEt = Et(simPart);
150                   if(track->Charge()>0){ FillHisto2D(Form("EtReconstructed%sIdentifiedPiPlus",cutName->Data()),track->Pt(),track->Eta(),myEt);}
151                   else{ FillHisto2D(Form("EtReconstructed%sIdentifiedPiMinus",cutName->Data()),track->Pt(),track->Eta(),myEt);}
152                   FillHisto2D(Form("dEdxPion%s",cutName->Data()),track->P(),dEdx,1.0);
153                 }
154                 if(isProton){
155                   if(pdgCode!=fProtonCode && pdgCode!=fAntiProtonCode){
156                     FillHisto2D(Form("MisidentifiedPIDs%s",cutName->Data()),2,mypid,1);
157                     // if(mypid==0)cerr<<"I was misidentified!  I'm not a proton! I am a "<<simPart->GetName()<<endl;
158                   }
159                   float myEt = Et(simPart);
160                   if(track->Charge()>0){ FillHisto2D(Form("EtReconstructed%sIdentifiedProton",cutName->Data()),track->Pt(),track->Eta(),myEt);}
161                   else{ FillHisto2D(Form("EtReconstructed%sIdentifiedAntiProton",cutName->Data()),track->Pt(),track->Eta(),myEt);}
162                   FillHisto2D(Form("dEdxProton%s",cutName->Data()),track->P(),dEdx,1.0);
163                 }
164                 if(isKaon){
165                   if(pdgCode!=fKMinusCode && pdgCode!=fKPlusCode){
166                     FillHisto2D(Form("MisidentifiedPIDs%s",cutName->Data()),3,mypid,1);
167                     //if(mypid==0)cerr<<"I was misidentified!  I'm not a kaon! I am a "<<simPart->GetName()<<" p "<<track->P()<<" nSigmaProton "<<nSigmaProton<<" nSigmaPion "<<nSigmaPion<<" nSigmaKaon "<<nSigmaKaon<<endl;
168                   }
169                   float myEt = Et(simPart);
170                   if(track->Charge()>0){ FillHisto2D(Form("EtReconstructed%sIdentifiedKPlus",cutName->Data()),track->Pt(),track->Eta(),myEt);}
171                   else{ FillHisto2D(Form("EtReconstructed%sIdentifiedKMinus",cutName->Data()),track->Pt(),track->Eta(),myEt);}
172                   FillHisto2D(Form("dEdxKaon%s",cutName->Data()),track->P(),dEdx,1.0);
173                 }
174                 if(isElectron){
175                   if(pdgCode!=fEMinusCode && pdgCode!=fEPlusCode){
176                     FillHisto2D(Form("MisidentifiedPIDs%s",cutName->Data()),4,mypid,1);
177                     //cerr<<"I was misidentified!  I'm not an electron! I am a "<<simPart->GetName()<<endl;
178                   }
179                   float myEt = Et(simPart);
180                   if(track->Charge()>0){ FillHisto2D(Form("EtReconstructed%sIdentifiedEPlus",cutName->Data()),track->Pt(),track->Eta(),myEt);}
181                   else{ FillHisto2D(Form("EtReconstructed%sIdentifiedEMinus",cutName->Data()),track->Pt(),track->Eta(),myEt);}
182                   FillHisto2D(Form("dEdxElectron%s",cutName->Data()),track->P(),dEdx,1.0);
183                 }
184                 if(unidentified){
185                   if(pdgCode!=fEMinusCode && pdgCode!=fEPlusCode){
186                     float myEtPi = Et(simPart,fPionMass);
187                     float myEt = Et(simPart);
188                     FillHisto2D(Form("EtReconstructed%sUnidentifiedAssumingPion",cutName->Data()),track->Pt(),track->Eta(),myEtPi);
189                     FillHisto2D(Form("EtReconstructed%sUnidentified",cutName->Data()),track->Pt(),track->Eta(),myEt);
190                     FillHisto2D(Form("EtNReconstructed%sUnidentified",cutName->Data()),track->Pt(),track->Eta(),1.0);
191                   }
192                   FillHisto2D(Form("dEdxUnidentified%s",cutName->Data()),track->P(),dEdx,1.0);
193                   //cout<<"I was not identified.  I am a "<<simPart->GetName()<<" PID "<<pdgCode<<endl;
194                   //track what was not identified successfully
195                   FillHisto1D(Form("UnidentifiedPIDs%s",cutName->Data()),mypid,1);
196                 }
197                 //...simulated
198                 if(pdgCode == fPiPlusCode){
199                   //cout<<"I'm a real primary "<<simPart->GetName()<<"! "<<"my label is "<<simPart->GetFirstMother()<<" track no "<<iTrack<<"/"<<realEvent->GetNumberOfTracks()<<endl;//<<" "<<label<<" "<<pdgCode<<endl;
200                 
201                   float myEt = Et(simPart);
202                   FillHisto2D(Form("EtReconstructed%sPiPlus",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
203                   FillHisto2D(Form("EtReconstructed%sChargedHadron",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
204                   FillHisto2D(Form("EtNReconstructed%sPiPlus",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
205                   FillHisto2D(Form("EtNReconstructed%sChargedHadron",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
206                   FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingPion",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
207                   filled = true;
208                 }
209                 if(pdgCode == fPiMinusCode){
210                   float myEt = Et(simPart);
211                   FillHisto2D(Form("EtReconstructed%sPiMinus",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
212                   FillHisto2D(Form("EtReconstructed%sChargedHadron",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
213                   FillHisto2D(Form("EtNReconstructed%sPiMinus",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
214                   FillHisto2D(Form("EtNReconstructed%sChargedHadron",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
215                   FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingPion",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
216                   filled = true;
217                 }
218                 if(pdgCode == fKPlusCode){
219                   float myEt = Et(simPart);
220                   float myEtPi = Et(simPart,fPionMass);
221                   FillHisto2D(Form("EtReconstructed%sKPlus",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
222                   FillHisto2D(Form("EtReconstructed%sChargedHadron",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
223                   FillHisto2D(Form("EtNReconstructed%sKPlus",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
224                   FillHisto2D(Form("EtNReconstructed%sChargedHadron",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
225                   FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingPion",cutName->Data()),simPart->Pt(),simPart->Eta(),myEtPi);
226                   FillHisto2D(Form("EtReconstructed%sKPlusAssumingPion",cutName->Data()),simPart->Pt(),simPart->Eta(),myEtPi);
227                   filled = true;
228                 }
229                 if(pdgCode == fKMinusCode){
230                   float myEt = Et(simPart);
231                   float myEtPi = Et(simPart,fPionMass);
232                   FillHisto2D(Form("EtReconstructed%sKMinus",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
233                   FillHisto2D(Form("EtReconstructed%sChargedHadron",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
234                   FillHisto2D(Form("EtNReconstructed%sKMinus",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
235                   FillHisto2D(Form("EtNReconstructed%sChargedHadron",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
236                   FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingPion",cutName->Data()),simPart->Pt(),simPart->Eta(),myEtPi);
237                   FillHisto2D(Form("EtReconstructed%sKMinusAssumingPion",cutName->Data()),simPart->Pt(),simPart->Eta(),myEtPi);
238                   filled = true;
239                 }
240                 if(pdgCode == fProtonCode){
241                   float myEt = Et(simPart);
242                   float myEtPi = Et(simPart,fPionMass);
243                   FillHisto2D(Form("EtReconstructed%sProton",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
244                   FillHisto2D(Form("EtReconstructed%sChargedHadron",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
245                   FillHisto2D(Form("EtNReconstructed%sProton",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
246                   FillHisto2D(Form("EtNReconstructed%sChargedHadron",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
247                   FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingPion",cutName->Data()),simPart->Pt(),simPart->Eta(),myEtPi);
248                   FillHisto2D(Form("EtReconstructed%sProtonAssumingPion",cutName->Data()),simPart->Pt(),simPart->Eta(),myEtPi);
249                   filled = true;
250                 }
251                 if(pdgCode == fAntiProtonCode){
252                   float myEt = Et(simPart);
253                   float myEtPi = Et(simPart,fPionMass);
254                   FillHisto2D(Form("EtReconstructed%sAntiProton",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
255                   FillHisto2D(Form("EtReconstructed%sChargedHadron",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
256                   FillHisto2D(Form("EtNReconstructed%sAntiProton",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
257                   FillHisto2D(Form("EtNReconstructed%sChargedHadron",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
258                   FillHisto2D(Form("EtReconstructed%sChargedHadronAssumingPion",cutName->Data()),simPart->Pt(),simPart->Eta(),myEtPi);
259                   FillHisto2D(Form("EtReconstructed%sAntiProtonAssumingPion",cutName->Data()),simPart->Pt(),simPart->Eta(),myEtPi);
260                   filled = true;
261                 }
262                 if(pdgCode == fEPlusCode){
263                   float myEt = Et(simPart);
264                   FillHisto2D(Form("EtReconstructed%sEPlus",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
265                   if(!isElectron || unidentified){
266                     float myEtPi = Et(simPart,fPionMass);
267                     FillHisto2D(Form("EtReconstructed%sMisidentifiedElectrons",cutName->Data()),track->Pt(),track->Eta(),myEtPi);
268                   }
269                   filled = true;
270                 }
271                 if(pdgCode == fEMinusCode){
272                   if(!isElectron || unidentified){
273                     float myEtPi = Et(simPart,fPionMass);
274                     FillHisto2D(Form("EtReconstructed%sMisidentifiedElectrons",cutName->Data()),track->Pt(),track->Eta(),myEtPi);
275                   }
276                   float myEt = Et(simPart);
277                     FillHisto2D(Form("EtReconstructed%sEMinus",cutName->Data()),simPart->Pt(),simPart->Eta(),myEt);
278                     filled = true;
279                 }
280                 //if(!filled){
281                   //TParticlePDG *pc = simPart->GetPDG(0);
282                   //if( strcmp(pc->ParticleClass(),"Baryon")==0 || strcmp(pc->ParticleClass(),"Meson")==0 ){
283                   //cout<<"Did not find a place for "<<simPart->GetName()<<" "<<pdgCode<<" which is a "<<pc->ParticleClass()<<endl;
284                   //}
285                   //}
286               }
287               
288             }
289             else{//not a primary - we're after V0 daughters!
290               //cout<<"I'm a secondary "<<simPart->GetName()<<"!";//<<endl;
291               TParticle *mom = stack->Particle(simPart->GetFirstMother());
292               if(mom){
293                 TParticlePDG *pc = mom->GetPDG(0);
294                 if(pc){
295                   Int_t pdgCode =  mom->GetPDG(0)->PdgCode();
296                   if(pdgCode == fLambdaCode){
297                     float myEt = Et(simPart);
298                     FillHisto2D(Form("EtReconstructed%sLambdaDaughters",cutName->Data()),track->Pt(),track->Eta(),myEt);
299                   }
300                   if(pdgCode == fAntiLambdaCode){
301                     float myEt = Et(simPart);
302                     FillHisto2D(Form("EtReconstructed%sAntiLambdaDaughters",cutName->Data()),track->Pt(),track->Eta(),myEt);
303                   }
304                   if(pdgCode == fK0SCode){
305                     float myEt = Et(simPart);
306                     FillHisto2D(Form("EtReconstructed%sK0SDaughters",cutName->Data()),track->Pt(),track->Eta(),myEt);
307                   }
308                   if(pdgCode == fXiCode){
309                     float myEt = Et(simPart);
310                     FillHisto2D(Form("EtReconstructed%sXiDaughters",cutName->Data()),track->Pt(),track->Eta(),myEt);
311                   }
312                   if(pdgCode == fAntiXiCode){
313                     float myEt = Et(simPart);
314                     FillHisto2D(Form("EtReconstructed%sAntiXiDaughters",cutName->Data()),track->Pt(),track->Eta(),myEt);
315                   }
316                   if(pdgCode == fOmegaCode){
317                     float myEt = Et(simPart);
318                     FillHisto2D(Form("EtReconstructed%sOmegaDaughters",cutName->Data()),track->Pt(),track->Eta(),myEt);
319                   }
320                   if(pdgCode == fXiCode){
321                     float myEt = Et(simPart);
322                     FillHisto2D(Form("EtReconstructed%sAntiOmegaDaughters",cutName->Data()),track->Pt(),track->Eta(),myEt);
323                   }
324
325                   if(mom->GetFirstMother()>0){
326                     TParticle *grandma = stack->Particle(mom->GetFirstMother());
327                     if(grandma){
328                       Int_t pdgCodeMom =  mom->GetPDG(0)->PdgCode();
329                       if(pdgCodeMom==fPiPlusCode || pdgCodeMom==fPiMinusCode || pdgCodeMom==fProtonCode ||pdgCodeMom==fAntiProtonCode || pdgCodeMom==fKPlusCode || pdgCode==fKMinusCode){
330                         //cout<<" my grandmother is "<<grandma->GetName()<<" "<<endl;
331                         Int_t pdgCodeGrandma =  grandma->GetPDG(0)->PdgCode();
332                       
333                         if(pdgCodeGrandma == fXiCode){
334                           float myEt = Et(simPart);
335                           FillHisto2D(Form("EtReconstructed%sXiDaughters",cutName->Data()),track->Pt(),track->Eta(),myEt);
336                         }
337                         if(pdgCodeGrandma == fAntiXiCode){
338                           float myEt = Et(simPart);
339                           FillHisto2D(Form("EtReconstructed%sAntiXiDaughters",cutName->Data()),track->Pt(),track->Eta(),myEt);
340                         }
341                         if(pdgCodeGrandma == fOmegaCode){
342                           float myEt = Et(simPart);
343                           FillHisto2D(Form("EtReconstructed%sOmegaDaughters",cutName->Data()),track->Pt(),track->Eta(),myEt);
344                         }
345                         if(pdgCodeGrandma == fXiCode){
346                           float myEt = Et(simPart);
347                           FillHisto2D(Form("EtReconstructed%sAntiOmegaDaughters",cutName->Data()),track->Pt(),track->Eta(),myEt);
348                         }
349
350                       }
351                     }
352                   }
353                 }
354               }
355             }
356           }
357
358         }
359       }
360     delete list;
361   }
362   //delete AliESDpid;
363   return 1;
364 }
365 Int_t AliAnalysisHadEtMonteCarlo::AnalyseEvent(AliVEvent* ev)
366 { // analyse MC event
367      ResetEventValues();
368      
369     // Get us an mc event
370     AliMCEvent *mcEvent = dynamic_cast<AliMCEvent*>(ev);
371
372     // Let's play with the stack!
373     AliStack *stack = mcEvent->Stack();
374
375     Int_t nPrim = stack->GetNtrack();
376
377     //=================Tracks which may or may not have been reconstructed=================
378
379     for (Int_t iPart = 0; iPart < nPrim; iPart++)
380     {
381
382         TParticle *part = stack->Particle(iPart);
383
384         if (!part)
385           {
386             Printf("ERROR: Could not get particle %d", iPart);
387             continue;
388           }
389
390         TParticlePDG *pc = part->GetPDG(0);
391
392         // Check if it is a primary particle
393         if (stack->IsPhysicalPrimary(iPart)){//primaries
394
395           if (TMath::Abs(part->Eta()) < fCuts->GetCommonEtaCut())           {
396
397             Int_t pdgCode =  part->GetPDG(0)->PdgCode();
398             //cout<<pdgCode->PdgCode()<<" ";
399             //fPdgDB->GetParticle("pi+")->PdgCode();
400             bool filled = false;
401             //============Charged hadrons===================================
402             if(pdgCode == fPiPlusCode){
403               //cout<<"I'm a simulated primary "<<part->GetName()<<"! "<<"my label is "<<part->GetFirstMother()<<" pt "<<part->Pt()<<endl;
404               float myEt = Et(part);
405               FillHisto2D("EtSimulatedPiPlus",part->Pt(),part->Eta(),myEt);
406               FillHisto2D("EtNSimulatedPiPlus",part->Pt(),part->Eta(),1.0);
407               FillHisto2D("EtSimulatedChargedHadron",part->Pt(),part->Eta(),myEt);
408               FillHisto2D("EtNSimulatedChargedHadron",part->Pt(),part->Eta(),1.0);
409               FillHisto2D("EtSimulatedChargedHadronAssumingPion",part->Pt(),part->Eta(),myEt);
410               FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
411               filled = true;
412             }
413             if(pdgCode == fPiMinusCode){
414               float myEt = Et(part);
415               FillHisto2D("EtSimulatedPiMinus",part->Pt(),part->Eta(),myEt);
416               FillHisto2D("EtNSimulatedPiMinus",part->Pt(),part->Eta(),1.0);
417               FillHisto2D("EtSimulatedChargedHadron",part->Pt(),part->Eta(),myEt);
418               FillHisto2D("EtNSimulatedChargedHadron",part->Pt(),part->Eta(),1.0);
419               FillHisto2D("EtSimulatedChargedHadronAssumingPion",part->Pt(),part->Eta(),myEt);
420               FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
421               filled = true;
422             }
423             if(pdgCode == fKPlusCode){
424               float myEt = Et(part);
425               float myEtPi = Et(part,fPionMass);
426               FillHisto2D("EtSimulatedKPlus",part->Pt(),part->Eta(),myEt);
427               FillHisto2D("EtNSimulatedKPlus",part->Pt(),part->Eta(),1.0);
428               FillHisto2D("EtSimulatedChargedHadron",part->Pt(),part->Eta(),myEt);
429               FillHisto2D("EtNSimulatedChargedHadron",part->Pt(),part->Eta(),1.0);
430               FillHisto2D("EtSimulatedChargedHadronAssumingPion",part->Pt(),part->Eta(),myEtPi);
431               FillHisto2D("EtSimulatedKPlusAssumingPion",part->Pt(),part->Eta(),myEtPi);
432               FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
433               filled = true;
434             }
435             if(pdgCode == fKMinusCode){
436               float myEt = Et(part);
437               float myEtPi = Et(part,fPionMass);
438               FillHisto2D("EtSimulatedKMinus",part->Pt(),part->Eta(),myEt);
439               FillHisto2D("EtNSimulatedKMinus",part->Pt(),part->Eta(),1.0);
440               FillHisto2D("EtSimulatedChargedHadron",part->Pt(),part->Eta(),myEt);
441               FillHisto2D("EtNSimulatedChargedHadron",part->Pt(),part->Eta(),1.0);
442               FillHisto2D("EtSimulatedChargedHadronAssumingPion",part->Pt(),part->Eta(),myEtPi);
443               FillHisto2D("EtSimulatedKMinusAssumingPion",part->Pt(),part->Eta(),myEtPi);
444               FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
445               filled = true;
446             }
447             if(pdgCode == fProtonCode){
448               float myEt = Et(part);
449               float myEtPi = Et(part,fPionMass);
450               FillHisto2D("EtSimulatedProton",part->Pt(),part->Eta(),myEt);
451               FillHisto2D("EtNSimulatedProton",part->Pt(),part->Eta(),1.0);
452               FillHisto2D("EtSimulatedChargedHadron",part->Pt(),part->Eta(),myEt);
453               FillHisto2D("EtNSimulatedChargedHadron",part->Pt(),part->Eta(),1.0);
454               FillHisto2D("EtSimulatedChargedHadronAssumingPion",part->Pt(),part->Eta(),myEtPi);
455               FillHisto2D("EtSimulatedProtonAssumingPion",part->Pt(),part->Eta(),myEtPi);
456               FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
457               filled = true;
458             }
459             if(pdgCode == fAntiProtonCode){
460               float myEt = Et(part);
461               float myEtPi = Et(part,fPionMass);
462               FillHisto2D("EtSimulatedAntiProton",part->Pt(),part->Eta(),myEt);
463               FillHisto2D("EtNSimulatedAntiProton",part->Pt(),part->Eta(),1.0);
464               FillHisto2D("EtSimulatedChargedHadron",part->Pt(),part->Eta(),myEt);
465               FillHisto2D("EtNSimulatedChargedHadron",part->Pt(),part->Eta(),1.0);
466               FillHisto2D("EtSimulatedChargedHadronAssumingPion",part->Pt(),part->Eta(),myEtPi);
467               FillHisto2D("EtSimulatedAntiProtonAssumingPion",part->Pt(),part->Eta(),myEtPi);
468               FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
469               filled = true;
470             }
471             //============Other hadrons===================================
472
473             if(pdgCode == fNeutronCode){
474               float myEt = Et(part);
475               FillHisto2D("EtSimulatedNeutron",part->Pt(),part->Eta(),myEt);
476               FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
477               filled = true;
478             }
479             if(pdgCode == fAntiNeutronCode){
480               float myEt = Et(part);
481               FillHisto2D("EtSimulatedAntiNeutron",part->Pt(),part->Eta(),myEt);
482               FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
483               filled = true;
484             }
485             if(pdgCode == fLambdaCode){
486               float myEt = Et(part);
487               //cout<<"I am a simulated lambda! pt "<<part->Pt()<<" eta "<<part->Eta()<<endl;
488               FillHisto2D("EtSimulatedLambda",part->Pt(),part->Eta(),myEt);
489               FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
490               Int_t ndaughters = part->GetNDaughters();
491               for(Int_t idaughter = 0;idaughter<ndaughters;idaughter++){
492                 Int_t daughterindex = part->GetDaughter(idaughter);
493                 if(daughterindex<0 || daughterindex>1e5) continue;
494                 TParticle *daughter = stack->ParticleFromTreeK(daughterindex);
495                 if(daughter){
496                   if(daughter->GetPDG(0)){
497                     Int_t daughtercode = daughter->GetPDG(0)->PdgCode();
498                     if(daughtercode==fPiMinusCode || daughtercode==fProtonCode){
499                       myEt = Et(daughter);
500                       FillHisto2D("EtSimulatedLambdaDaughters",daughter->Pt(),daughter->Eta(),myEt);
501                       //cout<<"Lambda daughter is a "<<daughter->GetName()<<endl;
502                     }
503                   }
504                   else{
505                     //cout<<"Lambda daughter is a "<<daughter->GetName()<<endl;
506                   }
507                 }
508               }
509               filled = true;
510             }
511             if(pdgCode == fAntiLambdaCode){
512               float myEt = Et(part);
513               FillHisto2D("EtSimulatedAntiLambda",part->Pt(),part->Eta(),myEt);
514               FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
515               Int_t ndaughters = part->GetNDaughters();
516               for(Int_t idaughter = 0;idaughter<ndaughters;idaughter++){
517                 Int_t daughterindex = part->GetDaughter(idaughter);
518                 if(daughterindex<0 || daughterindex>1e5) continue;
519                 TParticle *daughter = stack->ParticleFromTreeK(daughterindex);
520                 if(daughter){
521                   if(daughter->GetPDG(0)){
522                     Int_t daughtercode = daughter->GetPDG(0)->PdgCode();
523                     if(daughtercode==fPiPlusCode || daughtercode==fAntiProtonCode){
524                       myEt = Et(daughter);
525                       FillHisto2D("EtSimulatedAntiLambdaDaughters",daughter->Pt(),daughter->Eta(),myEt);
526                       //cout<<"AntiLambda daughter is a "<<daughter->GetName()<<endl;
527                     }
528                   }
529                   else{
530                     //cout<<"AntiLambda daughter is a "<<daughter->GetName()<<endl;
531                   }
532                 }
533               }
534               filled = true;
535             }
536             if(pdgCode == fK0SCode){
537               float myEt = Et(part);
538               FillHisto2D("EtSimulatedK0S",part->Pt(),part->Eta(),myEt);
539               FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
540               Int_t ndaughters = part->GetNDaughters();
541               for(Int_t idaughter = 0;idaughter<ndaughters;idaughter++){
542                 Int_t daughterindex = part->GetDaughter(idaughter);
543                 if(daughterindex<0 || daughterindex>1e5) continue;
544                 TParticle *daughter = stack->ParticleFromTreeK(daughterindex);
545                 if(daughter){
546                   if(daughter->GetPDG(0)){
547
548                     Int_t daughtercode = daughter->GetPDG(0)->PdgCode();
549                     if(daughtercode==fPiMinusCode || daughtercode==fPiPlusCode){
550                       myEt = Et(daughter);
551                       FillHisto2D("EtSimulatedK0SDaughters",daughter->Pt(),daughter->Eta(),myEt);
552                       //cout<<"K0S daughter is a "<<daughter->GetName()<<endl;
553                     }
554                   }
555                   else{
556                     //cout<<"K0S daughter is a "<<daughter->GetName()<<endl;
557                   }
558                 }
559               }
560               filled = true;
561             }
562             if(pdgCode == fK0LCode){
563               float myEt = Et(part);
564               FillHisto2D("EtSimulatedK0L",part->Pt(),part->Eta(),myEt);
565               FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
566               filled = true;
567             }
568             if(pdgCode == fOmegaCode){
569               float myEt = Et(part);
570               FillHisto2D("EtSimulatedOmega",part->Pt(),part->Eta(),myEt);
571               FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
572               Int_t ndaughters = part->GetNDaughters();
573               for(Int_t idaughter = 0;idaughter<ndaughters;idaughter++){
574                 Int_t daughterindex = part->GetDaughter(idaughter);
575                 if(daughterindex<0 || daughterindex>1e5) continue;
576                 TParticle *daughter = stack->Particle(daughterindex);
577                 if(daughter){
578                   if(daughter->GetPDG(0)){
579
580                     Int_t daughtercode = daughter->GetPDG(0)->PdgCode();
581                     if(daughtercode==fPiPlusCode || daughtercode==fProtonCode || daughtercode==fKMinusCode){
582                       myEt = Et(daughter);
583                       FillHisto2D("EtSimulatedOmegaDaughters",daughter->Pt(),daughter->Eta(),myEt);
584                     //cout<<"Omega daughter is a "<<daughter->GetName()<<endl;
585                     }
586                   }
587                   else{
588                     //cout<<"Omega daughter is a "<<daughter->GetName()<<endl;
589                   }
590                 }
591               }
592               filled = true;
593             }
594             if(pdgCode == fAntiOmegaCode){
595               float myEt = Et(part);
596               FillHisto2D("EtSimulatedOmega",part->Pt(),part->Eta(),myEt);
597               FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
598               Int_t ndaughters = part->GetNDaughters();
599               for(Int_t idaughter = 0;idaughter<ndaughters;idaughter++){
600                 Int_t daughterindex = part->GetDaughter(idaughter);
601                 if(daughterindex<0 || daughterindex>1e5) continue;
602                 TParticle *daughter = stack->ParticleFromTreeK(daughterindex);
603                 if(daughter){
604                   if(daughter->GetPDG(0)){
605                     Int_t daughtercode = daughter->GetPDG(0)->PdgCode();
606                     if(daughtercode==fPiMinusCode || daughtercode==fAntiProtonCode || daughtercode==fKPlusCode){
607                       myEt = Et(daughter);
608                       FillHisto2D("EtSimulatedAntiOmegaDaughters",daughter->Pt(),daughter->Eta(),myEt);
609                       //cout<<"AntiOmega daughter is a "<<daughter->GetName()<<endl;
610                     }
611                   }
612                   else{
613                     //cout<<"AntiOmega daughter is a "<<daughter->GetName()<<endl;
614                   }
615                 }
616               }
617               filled = true;
618             }
619             //There are two codes for Sigmas
620             if(pdgCode == fSigmaCode || pdgCode == -3222){
621               float myEt = Et(part);
622               FillHisto2D("EtSimulatedSigma",part->Pt(),part->Eta(),myEt);
623               FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
624               filled = true;
625             }
626             if(pdgCode == fAntiSigmaCode || pdgCode == 3222){
627               float myEt = Et(part);
628               FillHisto2D("EtSimulatedAntiSigma",part->Pt(),part->Eta(),myEt);
629               FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
630               filled = true;
631             }
632             if(pdgCode == fXiCode){
633               float myEt = Et(part);
634               FillHisto2D("EtSimulatedXi",part->Pt(),part->Eta(),myEt);
635               FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
636               Int_t ndaughters = part->GetNDaughters();
637               for(Int_t idaughter = 0;idaughter<ndaughters;idaughter++){
638                 Int_t daughterindex = part->GetDaughter(idaughter);
639                 if(daughterindex<0 || daughterindex>1e5 || daughterindex>1e5) continue;
640                 //cerr<<"Daughter index "<<daughterindex<<" npart "<<nPrim<<endl;
641                 TParticle *daughter = stack->ParticleFromTreeK(daughterindex);
642                 if(daughter){
643                   if(daughter->GetPDG(0)){
644
645                     Int_t daughtercode = daughter->GetPDG(0)->PdgCode();
646                     if(daughtercode==fPiPlusCode || daughtercode==fProtonCode || daughtercode==fPiMinusCode){
647                       myEt = Et(daughter);
648                       FillHisto2D("EtSimulatedXiDaughters",daughter->Pt(),daughter->Eta(),myEt);
649                     //cout<<"Xi daughter is a "<<daughter->GetName()<<endl;
650                     }
651                   }
652                   else{
653                     //cout<<"Xi daughter is a "<<daughter->GetName()<<endl;
654                   }
655                 }
656               }
657               filled = true;
658             }
659             if(pdgCode == fAntiXiCode){
660               float myEt = Et(part);
661               FillHisto2D("EtSimulatedAntiXi",part->Pt(),part->Eta(),myEt);
662               FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
663               Int_t ndaughters = part->GetNDaughters();
664               for(Int_t idaughter = 0;idaughter<ndaughters;idaughter++){
665                 Int_t daughterindex = part->GetDaughter(idaughter);
666                 if(daughterindex<0 || daughterindex>1e5) continue;
667                 TParticle *daughter = stack->ParticleFromTreeK(daughterindex);
668                 if(daughter){
669                   if(daughter->GetPDG(0)){
670                     Int_t daughtercode = daughter->GetPDG(0)->PdgCode();
671                     if(daughtercode==fPiPlusCode || daughtercode==fAntiProtonCode || daughtercode==fPiMinusCode){
672                       myEt = Et(daughter);
673                       FillHisto2D("EtSimulatedAntiXiDaughters",daughter->Pt(),daughter->Eta(),myEt);
674                       //cout<<"AntiXi daughter is a "<<daughter->GetName()<<endl;
675                     }
676                   }
677                   else{
678                     //cout<<"AntiXi daughter is a "<<daughter->GetName()<<endl;
679                   }
680                 }
681               }
682               filled = true;
683             }
684             if(pdgCode == fXi0Code){
685               float myEt = Et(part);
686               FillHisto2D("EtSimulatedXi0",part->Pt(),part->Eta(),myEt);
687               FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
688               filled = true;
689             }
690             if(pdgCode == fAntiXi0Code){
691               float myEt = Et(part);
692               FillHisto2D("EtSimulatedAntiXi0",part->Pt(),part->Eta(),myEt);
693               FillHisto2D("EtSimulatedAllHadron",part->Pt(),part->Eta(),myEt);
694               filled = true;
695             }
696             //============electrons===================================
697
698             if(pdgCode == fEPlusCode){
699               float myEt = Et(part);
700               FillHisto2D("EtSimulatedEPlus",part->Pt(),part->Eta(),myEt);
701               filled = true;
702             }
703             if(pdgCode == fEMinusCode){
704               float myEt = Et(part);
705               FillHisto2D("EtSimulatedEMinus",part->Pt(),part->Eta(),myEt);
706               filled = true;
707             }
708             //============neutrals===================================
709             if(pdgCode == fGammaCode){
710               TParticle *mom = stack->Particle(part->GetFirstMother());
711               Int_t pdgCodeMom =  mom->GetPDG(0)->PdgCode();
712               //cout<<"I am a gamma and my mom is "<<mom->GetName()<<endl;
713               if(pdgCodeMom == fEtaCode){
714                 float myEt = Et(mom);
715                 FillHisto2D("EtSimulatedEta",mom->Pt(),mom->Eta(),myEt);
716                 filled = true;
717               }
718               if(pdgCodeMom == fPi0Code){
719                 float myEt = Et(mom);
720                 FillHisto2D("EtSimulatedPi0",mom->Pt(),mom->Eta(),myEt);
721                 filled = true;
722               }
723               if(pdgCodeMom == fOmega0Code){
724                 float myEt = Et(mom);
725                 FillHisto2D("EtSimulatedOmega0",mom->Pt(),mom->Eta(),myEt);
726                 filled = true;
727               }
728               if(!filled){
729                 float myEt = Et(part);
730                 FillHisto2D("EtSimulatedGamma",part->Pt(),part->Eta(),myEt);
731                 filled = true;
732               }
733             }
734             if(pdgCode == fEtaCode){
735               float myEt = Et(part);
736               FillHisto2D("EtSimulatedEta",part->Pt(),part->Eta(),myEt);
737               filled = true;
738             }
739             if(pdgCode == fPi0Code){
740               float myEt = Et(part);
741               FillHisto2D("EtSimulatedPi0",part->Pt(),part->Eta(),myEt);
742               filled = true;
743             }
744             if(pdgCode == fOmega0Code){
745               float myEt = Et(part);
746               FillHisto2D("EtSimulatedOmega0",part->Pt(),part->Eta(),myEt);
747               filled = true;
748             }
749             if(!filled){
750               //if( strcmp(pc->ParticleClass(),"Baryon")==0 || strcmp(pc->ParticleClass(),"Meson")==0 ){
751                 cout<<"Did not find a place for "<<part->GetName()<<" "<<pdgCode<<" which is a "<<pc->ParticleClass()<<endl;
752                 //}
753             }
754           }
755         }
756     }
757
758
759
760
761 //     fTotNeutralEtAcc = fTotNeutralEt;
762 //     fTotEt = fTotChargedEt + fTotNeutralEt;
763 //     fTotEtAcc = fTotChargedEtAcc + fTotNeutralEtAcc;
764     
765 //     FillHistograms();
766
767     return 1;
768     
769 }
770
771 void AliAnalysisHadEtMonteCarlo::Init()
772 { // Init
773     AliAnalysisHadEt::Init();
774 }
775
776 void AliAnalysisHadEtMonteCarlo::CreateHistograms(){
777   //for simulated Et only (no reconstruction)
778   CreateEtaPtHisto2D(TString("EtSimulatedPiPlus"),TString("Simulated E_{T} from #pi^{+}"));
779   CreateEtaPtHisto2D("EtSimulatedPiMinus","Simulated E_{T} from #pi^{-}");
780   CreateEtaPtHisto2D("EtSimulatedKPlus","Simulated E_{T} from K^{+}");
781   CreateEtaPtHisto2D("EtSimulatedKMinus","Simulated E_{T} from K^{-}");
782   CreateEtaPtHisto2D("EtSimulatedProton","Simulated E_{T} from p");
783   CreateEtaPtHisto2D("EtSimulatedAntiProton","Simulated E_{T} from #bar{p}");
784   CreateEtaPtHisto2D("EtSimulatedChargedHadron","Simulated E_{T} from charged hadrons");
785   CreateEtaPtHisto2D(TString("EtNSimulatedPiPlus"),TString("Number of Simulated #pi^{+}"));
786   CreateEtaPtHisto2D("EtNSimulatedPiMinus","Number of simulated #pi^{-}");
787   CreateEtaPtHisto2D("EtNSimulatedKPlus","Number of simulated K^{+}");
788   CreateEtaPtHisto2D("EtNSimulatedKMinus","Number of simulated K^{-}");
789   CreateEtaPtHisto2D("EtNSimulatedProton","Number of simulated p");
790   CreateEtaPtHisto2D("EtNSimulatedAntiProton","Number of simulated #bar{p}");
791   CreateEtaPtHisto2D("EtNSimulatedChargedHadron","Number of simulated charged hadrons");
792
793   CreateEtaPtHisto2D("EtSimulatedChargedHadronAssumingPion","Simulated E_{T} from charged hadrons assuming they are all pions");
794   CreateEtaPtHisto2D("EtSimulatedKPlusAssumingPion","Simulated E_{T} from K^{+} assuming #pi mass");
795   CreateEtaPtHisto2D("EtSimulatedKMinusAssumingPion","Simulated E_{T} from K^{-} assuming #pi mass");
796   CreateEtaPtHisto2D("EtSimulatedProtonAssumingPion","Simulated E_{T} from p assuming #pi mass");
797   CreateEtaPtHisto2D("EtSimulatedAntiProtonAssumingPion","Simulated E_{T} from #bar{p} assuming #pi mass");
798
799   CreateEtaPtHisto2D("EtSimulatedLambda","Simulated E_{T} from #Lambda");
800   CreateEtaPtHisto2D("EtSimulatedAntiLambda","Simulated E_{T} from #bar{#Lambda}");
801   CreateEtaPtHisto2D("EtSimulatedK0S","Simulated E_{T} from K^{0}_{S}");
802   CreateEtaPtHisto2D("EtSimulatedK0L","Simulated E_{T} from K^{0}_{L}");
803   CreateEtaPtHisto2D("EtSimulatedNeutron","Simulated E_{T} from neutrons");
804   CreateEtaPtHisto2D("EtSimulatedAntiNeutron","Simulated E_{T} from #bar{n}");
805   CreateEtaPtHisto2D("EtSimulatedEPlus","Simulated E_{T} from e^{+}");
806   CreateEtaPtHisto2D("EtSimulatedEMinus","Simulated E_{T} from e^{-}");
807   CreateEtaPtHisto2D("EtSimulatedOmega","Simulated E_{T} from #Omega^{-}");
808   CreateEtaPtHisto2D("EtSimulatedAntiOmega","Simulated E_{T} from #Omega^{+}");
809   CreateEtaPtHisto2D("EtSimulatedXi","Simulated E_{T} from #Xi^{-}");
810   CreateEtaPtHisto2D("EtSimulatedAntiXi","Simulated E_{T} from #Xi^{+}");
811   CreateEtaPtHisto2D("EtSimulatedSigma","Simulated E_{T} from #Xi^{-}");
812   CreateEtaPtHisto2D("EtSimulatedAntiSigma","Simulated E_{T} from #Xi^{+}");
813   CreateEtaPtHisto2D("EtSimulatedXi0","Simulated E_{T} from #Xi^{0}");
814   CreateEtaPtHisto2D("EtSimulatedAntiXi0","Simulated E_{T} from #Xi^{0}");
815   CreateEtaPtHisto2D("EtSimulatedAllHadron","Simulated E_{T} from all hadrons");
816
817
818   CreateEtaPtHisto2D("EtSimulatedLambdaDaughters","Simulated E_{T} from #Lambda Daughters");
819   CreateEtaPtHisto2D("EtSimulatedAntiLambdaDaughters","Simulated E_{T} from #bar{#Lambda} Daughters");
820   CreateEtaPtHisto2D("EtSimulatedK0SDaughters","Simulated E_{T} from K^{0}_{S} Daughters");
821   CreateEtaPtHisto2D("EtSimulatedOmegaDaughters","Simulated E_{T} from #Omega^{-} Daughters");
822   CreateEtaPtHisto2D("EtSimulatedAntiOmegaDaughters","Simulated E_{T} from #Omega^{+} Daughters");
823   CreateEtaPtHisto2D("EtSimulatedXiDaughters","Simulated E_{T} from #Xi^{-} Daughters");
824   CreateEtaPtHisto2D("EtSimulatedAntiXiDaughters","Simulated E_{T} from #Xi^{+} Daughters");
825
826
827   CreateEtaPtHisto2D("EtSimulatedGamma","Simulated E_{T} from #gamma");
828   CreateEtaPtHisto2D("EtSimulatedEta","Simulated E_{T} from #eta");
829   CreateEtaPtHisto2D("EtSimulatedPi0","Simulated E_{T} from #pi^{0}");
830   CreateEtaPtHisto2D("EtSimulatedOmega0","Simulated E_{T} from #omega");
831
832   TString *strTPC = new TString("TPC");
833   TString *strITS = new TString("ITS");
834   TString *strTPCITS = new TString("TPCITS");
835   for(Int_t i=0;i<3;i++){
836     TString *cutName;
837     Float_t maxPtdEdx = 10;
838     Float_t mindEdx = 35;
839     Float_t maxdEdx = 150.0;
840     switch(i){
841     case 0:
842       cutName = strTPC;
843       break;
844     case 1:
845       cutName = strITS;
846       maxPtdEdx = 5;
847       maxdEdx = 500.0;
848       break;
849     case 2:
850       cutName = strTPCITS;
851       break;
852     default:
853       cerr<<"Error:  cannot make histograms!"<<endl;
854       return;
855     }
856
857     CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedPiPlus",cutName->Data()),"Reconstructed E_{T} from identified #pi^{+}");
858     CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedPiMinus",cutName->Data()),"Reconstructed E_{T} from identified #pi^{-}");
859     CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedKPlus",cutName->Data()),"Reconstructed E_{T} from identified K^{+}");
860     CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedEMinus",cutName->Data()),"Reconstructed E_{T} from identified e^{-}");
861     CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedEPlus",cutName->Data()),"Reconstructed E_{T} from identified e^{+}");
862     CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedKMinus",cutName->Data()),"Reconstructed E_{T} from identified K^{-}");
863     CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedProton",cutName->Data()),"Reconstructed E_{T} from identified p");
864     CreateEtaPtHisto2D(Form("EtReconstructed%sIdentifiedAntiProton",cutName->Data()),"Reconstructed E_{T} from identified #bar{p}");
865     CreateEtaPtHisto2D(Form("EtNReconstructed%sUnidentified",cutName->Data()),"Number of Reconstructed unidentified particles");
866     CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentifiedAssumingPion",cutName->Data()),"Reconstructed E_{T} from unidentified particles assuming pion mass");
867     CreateEtaPtHisto2D(Form("EtReconstructed%sUnidentified",cutName->Data()),"Reconstructed E_{T} from unidentified particles using real mass");
868     CreateEtaPtHisto2D(Form("EtReconstructed%sMisidentifiedElectrons",cutName->Data()),"Reconstructed E_{T} from misidentified electrons");
869
870
871     CreateEtaPtHisto2D(Form("EtReconstructed%sPiPlus",cutName->Data()),"Reconstructed E_{T} from #pi^{+}");
872     CreateEtaPtHisto2D(Form("EtReconstructed%sPiMinus",cutName->Data()),"Reconstructed E_{T} from #pi^{-}");
873     CreateEtaPtHisto2D(Form("EtReconstructed%sKPlus",cutName->Data()),"Reconstructed E_{T} from K^{+}");
874     CreateEtaPtHisto2D(Form("EtReconstructed%sKMinus",cutName->Data()),"Reconstructed E_{T} from K^{-}");
875     CreateEtaPtHisto2D(Form("EtReconstructed%sProton",cutName->Data()),"Reconstructed E_{T} from p");
876     CreateEtaPtHisto2D(Form("EtReconstructed%sAntiProton",cutName->Data()),"Reconstructed E_{T} from #bar{p}");
877     CreateEtaPtHisto2D(Form("EtReconstructed%sChargedHadron",cutName->Data()),"Reconstructed E_{T} from charged hadrons");
878     CreateEtaPtHisto2D(Form("EtNReconstructed%sPiPlus",cutName->Data()),"Reconstructed E_{T} from #pi^{+}");
879     CreateEtaPtHisto2D(Form("EtNReconstructed%sPiMinus",cutName->Data()),"Reconstructed E_{T} from #pi^{-}");
880     CreateEtaPtHisto2D(Form("EtNReconstructed%sKPlus",cutName->Data()),"Reconstructed E_{T} from K^{+}");
881     CreateEtaPtHisto2D(Form("EtNReconstructed%sKMinus",cutName->Data()),"Reconstructed E_{T} from K^{-}");
882     CreateEtaPtHisto2D(Form("EtNReconstructed%sProton",cutName->Data()),"Reconstructed E_{T} from p");
883     CreateEtaPtHisto2D(Form("EtNReconstructed%sAntiProton",cutName->Data()),"Reconstructed E_{T} from #bar{p}");
884     CreateEtaPtHisto2D(Form("EtNReconstructed%sChargedHadron",cutName->Data()),"Reconstructed E_{T} from charged hadrons");
885
886     CreateEtaPtHisto2D(Form("EtReconstructed%sChargedHadronAssumingPion",cutName->Data()),"Reconstructed E_{T} from charged hadrons assuming they are all pions");
887     CreateEtaPtHisto2D(Form("EtReconstructed%sKPlusAssumingPion",cutName->Data()),"Reconstructed E_{T} from K^{+} assuming #pi mass");
888     CreateEtaPtHisto2D(Form("EtReconstructed%sKMinusAssumingPion",cutName->Data()),"Reconstructed E_{T} from K^{-} assuming #pi mass");
889     CreateEtaPtHisto2D(Form("EtReconstructed%sProtonAssumingPion",cutName->Data()),"Reconstructed E_{T} from p assuming #pi mass");
890     CreateEtaPtHisto2D(Form("EtReconstructed%sAntiProtonAssumingPion",cutName->Data()),"Reconstructed E_{T} from #bar{p} assuming #pi mass");
891
892     CreateEtaPtHisto2D(Form("EtReconstructed%sEPlus",cutName->Data()),"Reconstructed E_{T} from e^{+}");
893     CreateEtaPtHisto2D(Form("EtReconstructed%sEMinus",cutName->Data()),"Reconstructed E_{T} from e^{-}");
894
895
896
897     CreateEtaPtHisto2D(Form("EtReconstructed%sLambdaDaughters",cutName->Data()),"Reconstructed E_{T} from #Lambda Daughters");
898   CreateEtaPtHisto2D(Form("EtReconstructed%sAntiLambdaDaughters",cutName->Data()),"Reconstructed E_{T} from #bar{#Lambda} Daughters");
899   CreateEtaPtHisto2D(Form("EtReconstructed%sK0SDaughters",cutName->Data()),"Reconstructed E_{T} from K^{0}_{S} Daughters");
900   CreateEtaPtHisto2D(Form("EtReconstructed%sOmegaDaughters",cutName->Data()),"Reconstructed E_{T} from #Omega^{-} Daughters");
901   CreateEtaPtHisto2D(Form("EtReconstructed%sAntiOmegaDaughters",cutName->Data()),"Reconstructed E_{T} from #Omega^{+} Daughters");
902   CreateEtaPtHisto2D(Form("EtReconstructed%sXiDaughters",cutName->Data()),"Reconstructed E_{T} from #Xi^{-} Daughters");
903   CreateEtaPtHisto2D(Form("EtReconstructed%sAntiXiDaughters",cutName->Data()),"Reconstructed E_{T} from #Xi^{+} Daughters");
904
905     CreateIntHisto1D(Form("UnidentifiedPIDs%s",cutName->Data()),"PIDs of unidentified particles", "PID", "Number of particles",9, -4,4);
906     CreateHisto2D(Form("MisidentifiedPIDs%s",cutName->Data()),"PIDs of misidentified particles", "PID real","PID identified",5, -.5,4.5,5, -.5,4.5);
907     CreateHisto2D(Form("dEdxAll%s",cutName->Data()),"dE/dx for all particles","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
908     CreateHisto2D(Form("dEdxPion%s",cutName->Data()),"dE/dx for #pi^{#pm}","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
909     CreateHisto2D(Form("dEdxKaon%s",cutName->Data()),"dE/dx for K^{#pm}","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
910     CreateHisto2D(Form("dEdxProton%s",cutName->Data()),"dE/dx for p(#bar{p})","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
911     CreateHisto2D(Form("dEdxElectron%s",cutName->Data()),"dE/dx for e^{#pm}","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
912     CreateHisto2D(Form("dEdxUnidentified%s",cutName->Data()),"dE/dx for unidentified particles","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
913   }
914
915   CreateIntHisto1D("NEvents","Number of events","number of events","Number of events",1,0,1);
916
917   //CreateHisto1D("MisidentifiedPIDs","PIDs for particles misidentified that are not a #pi,K,p","PID","number of entries",3000,0.5,3000.5);
918 }
919