updates from salvatore
[u/mrichter/AliRoot.git] / PWGGA / EMCALJetTasks / AliAnalysisTaskScale.cxx
1 // $Id$
2 //
3 // Scale task.
4 //
5 // Author: R.Reed, M.Connors
6
7 #include "AliAnalysisTaskScale.h"
8
9 #include <TClonesArray.h>
10 #include <TF1.h>
11 #include <TH1F.h>
12 #include <TH2F.h>
13 #include <TList.h>
14 #include <TLorentzVector.h>
15
16 #include "AliAnalysisManager.h"
17 #include "AliCentrality.h"
18 #include "AliEMCALGeometry.h"
19 #include "AliLog.h"
20 #include "AliVCluster.h"
21 #include "AliVEvent.h"
22 #include "AliVTrack.h"
23
24 ClassImp(AliAnalysisTaskScale)
25
26 //________________________________________________________________________
27 AliAnalysisTaskScale::AliAnalysisTaskScale() : 
28   AliAnalysisTaskEmcal("AliAnalysisTaskScale", kTRUE), 
29   fScaleFunction(0),
30   fGeom(0),
31   fHistCentrality(0), 
32   fHistPtTPCvsCent(0), 
33   fHistPtEMCALvsCent(0), 
34   fHistEtvsCent(0),  
35   fHistScalevsCent(0),  
36   fHistDeltaScalevsCent(0), 
37   fHistPtTPCvsNtrack(0), 
38   fHistPtEMCALvsNtrack(0), 
39   fHistEtvsNtrack(0),  
40   fHistScalevsNtrack(0),  
41   fHistDeltaScalevsNtrack(0),
42   fHistTrackPtvsCent(0), 
43   fHistClusterPtvsCent(0),
44   fHistTrackEtaPhi(0), 
45   fHistClusterEtaPhi(0)
46 {
47   // Default constructor.
48 }
49
50 //________________________________________________________________________
51 AliAnalysisTaskScale::AliAnalysisTaskScale(const char *name) :
52   AliAnalysisTaskEmcal(name, kTRUE), 
53   fScaleFunction(0),
54   fGeom(0),
55   fHistCentrality(0), 
56   fHistPtTPCvsCent(0), 
57   fHistPtEMCALvsCent(0), 
58   fHistEtvsCent(0),  
59   fHistScalevsCent(0),  
60   fHistDeltaScalevsCent(0), 
61   fHistPtTPCvsNtrack(0), 
62   fHistPtEMCALvsNtrack(0), 
63   fHistEtvsNtrack(0),  
64   fHistScalevsNtrack(0),  
65   fHistDeltaScalevsNtrack(0),
66   fHistTrackPtvsCent(0), 
67   fHistClusterPtvsCent(0),
68   fHistTrackEtaPhi(0), 
69   fHistClusterEtaPhi(0)
70 {
71   // Constructor.
72
73 }
74
75 //________________________________________________________________________
76 void AliAnalysisTaskScale::UserCreateOutputObjects()
77 {
78   // Create my user objects.
79
80   OpenFile(1);
81   fOutput = new TList();
82   fOutput->SetOwner();
83
84   fHistCentrality         = new TH1F("Centrality","Centrality",              101, -1, 100);
85   fHistPtTPCvsCent        = new TH2F("PtTPCvsCent","rho vs cent",            101, -1, 100, 500,   0, 1000);
86   fHistPtEMCALvsCent      = new TH2F("PtEMCALvsCent","rho vs cent",          101, -1, 100, 500,   0, 1000);
87   fHistEtvsCent           = new TH2F("EtvsCent","rho vs cent",               101, -1, 100, 500,   0, 1000);
88   fHistScalevsCent        = new TH2F("ScalevsCent","rho vs cent",            101, -1, 100, 400,   0, 4);
89   fHistDeltaScalevsCent   = new TH2F("DeltaScalevsCent","rho vs cent",       101, -1, 100, 400,  -2, 2);
90   fHistPtTPCvsNtrack      = new TH2F("PtTPCvsNtrack","rho vs cent",          500,  0, 2500, 500,  0, 1000);
91   fHistPtEMCALvsNtrack    = new TH2F("PtEMCALvsNtrack","rho vs cent",        500,  0, 2500, 500,  0, 1000);
92   fHistEtvsNtrack         = new TH2F("EtvsNtrack","rho vs cent",             500,  0, 2500, 500,  0, 1000);
93   fHistScalevsNtrack      = new TH2F("ScalevsNtrack","rho vs cent",          500,  0, 2500, 400,  0, 4);
94   fHistDeltaScalevsNtrack = new TH2F("DeltaScalevsNtrack","rho vs cent",     500,  0, 2500, 400, -2, 2);
95   fHistTrackPtvsCent      = new TH2F("TrackPtvsCent","Track pt vs cent",     101, -1, 100,  500,  0, 100);
96   fHistClusterPtvsCent    = new TH2F("ClusterPtvsCent","Cluster pt vs cent", 101, -1, 100,  500,  0, 100);
97   fHistTrackEtaPhi        = new TH2F("TrackEtaPhi","Track eta phi",          100, -1.0, 1.0, 64,  0, 6.4);
98   fHistClusterEtaPhi      = new TH2F("ClusterEtaPhi","Cluster eta phi",      100, -1.0, 1.0, 64, -3.2, 3.2);
99
100   fOutput->Add(fHistCentrality);
101   fOutput->Add(fHistPtTPCvsCent);
102   fOutput->Add(fHistPtEMCALvsCent);
103   fOutput->Add(fHistEtvsCent);
104   fOutput->Add(fHistScalevsCent);
105   fOutput->Add(fHistDeltaScalevsCent);
106   fOutput->Add(fHistPtTPCvsNtrack);
107   fOutput->Add(fHistPtEMCALvsNtrack);
108   fOutput->Add(fHistEtvsNtrack);
109   fOutput->Add(fHistScalevsNtrack);
110   fOutput->Add(fHistDeltaScalevsNtrack);
111   fOutput->Add(fHistTrackPtvsCent);
112   fOutput->Add(fHistClusterPtvsCent);
113   fOutput->Add(fHistTrackEtaPhi);
114   fOutput->Add(fHistClusterEtaPhi);
115
116   PostData(1, fOutput);
117 }
118
119 //________________________________________________________________________
120 Double_t AliAnalysisTaskScale::GetScaleFactor(Double_t cent)
121 {
122   // Get scale function.
123
124   Double_t scale = -1;
125   if (fScaleFunction)
126     scale = fScaleFunction->Eval(cent);
127   return scale;
128 }
129
130 //________________________________________________________________________
131 void AliAnalysisTaskScale::ExecOnce() 
132 {
133   // Init the analysis.
134
135   fGeom = AliEMCALGeometry::GetInstance();
136
137   if (!fGeom) {
138     AliFatal("Can not create geometry");
139     return;
140   }
141
142   AliAnalysisTaskEmcal::ExecOnce();
143 }
144
145 //________________________________________________________________________
146 Bool_t AliAnalysisTaskScale::FillHistograms() 
147 {
148   // Execute on each event.
149
150   const Double_t EmcalMinEta = fGeom->GetArm1EtaMin();
151   const Double_t EmcalMaxEta = fGeom->GetArm1EtaMax();
152   const Double_t EmcalMinPhi = fGeom->GetArm1PhiMin() * TMath::DegToRad();
153   const Double_t EmcalMaxPhi = fGeom->GetArm1PhiMax() * TMath::DegToRad();
154
155   const Double_t TpcMinPhi   = 0;
156   const Double_t TpcMaxPhi   = 2 * TMath::Pi();
157
158   const Double_t TpcArea     = (TpcMaxPhi - TpcMinPhi) * (EmcalMinEta - EmcalMaxEta);
159   const Double_t EmcalArea   = (EmcalMaxPhi - EmcalMinPhi) * (EmcalMinEta - EmcalMaxEta);
160
161   Double_t ptTPC   = 0;
162   Double_t ptEMCAL = 0;
163
164   const Int_t Ntracks = fTracks->GetEntries();
165   for (Int_t iTracks = 0; iTracks < Ntracks; ++iTracks) {
166     AliVTrack *track = static_cast<AliVTrack*>(fTracks->At(iTracks));
167
168     if (!track)
169       continue;
170
171     if (!AcceptTrack(track))
172       continue;
173
174     if (TMath::Abs(track->Eta()) > 0.7)   // only accept tracks in the EMCal eta range
175       continue;
176
177     fHistTrackPtvsCent->Fill(fCent,track->Pt());
178     fHistTrackEtaPhi->Fill(track->Eta(),track->Phi());
179
180     ptTPC += track->Pt();
181     if ((track->Phi() > EmcalMaxPhi) || (track->Phi() < EmcalMinPhi))
182       continue;
183
184     ptEMCAL += track->Pt();
185   }
186
187   if (ptTPC == 0) 
188     return kFALSE;
189   
190   Double_t Et = 0;
191   const Int_t Nclus = fCaloClusters->GetEntries();
192   for (Int_t iClus = 0; iClus < Nclus; ++iClus) {
193     AliVCluster *c = static_cast<AliVCluster*>(fCaloClusters->At(iClus));
194     if (!c)
195       continue;
196
197     if (!AcceptCluster(c))
198       continue;
199
200     TLorentzVector nPart;
201     c->GetMomentum(nPart, fVertex);
202
203     fHistClusterPtvsCent->Fill(fCent, nPart.Pt());
204     fHistClusterEtaPhi->Fill(nPart.Eta(), nPart.Phi());
205
206     Et += nPart.Pt();
207   }
208
209   const Double_t scalecalc = ((Et + ptEMCAL) / EmcalArea) * (TpcArea / ptTPC);
210   const Double_t scale     = GetScaleFactor(fCent);
211
212   fHistCentrality->Fill(fCent);
213   fHistPtTPCvsCent->Fill(fCent, ptTPC);
214   fHistPtEMCALvsCent->Fill(fCent, ptEMCAL);
215   fHistEtvsCent->Fill(fCent, Et);
216   fHistScalevsCent->Fill(fCent, scalecalc);
217   fHistDeltaScalevsCent->Fill(fCent, scalecalc - scale);
218   fHistPtTPCvsNtrack->Fill(Ntracks, ptTPC);
219   fHistPtEMCALvsNtrack->Fill(Ntracks, ptEMCAL);
220   fHistEtvsNtrack->Fill(Ntracks, Et);
221   fHistScalevsNtrack->Fill(Ntracks, scalecalc);
222   fHistDeltaScalevsNtrack->Fill(Ntracks, scalecalc - scale);
223
224   return kTRUE;
225 }      
226
227 //________________________________________________________________________
228 void AliAnalysisTaskScale::Terminate(Option_t *) 
229 {
230   // Called once at the end of the analysis.
231 }