]>
Commit | Line | Data |
---|---|---|
3c40321c | 1 | /* $Id: AliAnalysisTaskPi0V2.cxx 55404 2012-03-29 10:10:19Z fca $ */ |
2 | ||
3c40321c | 3 | #include "AliAnalysisTaskPi0V2.h" |
4 | ||
cd0e5e33 | 5 | #include <Riostream.h> |
6 | #include <TCanvas.h> | |
7 | #include <TChain.h> | |
8 | #include <TClonesArray.h> | |
9 | #include <TH1F.h> | |
10 | #include <TH2F.h> | |
11 | #include <TH3F.h> | |
12 | #include <THnSparse.h> | |
13 | #include <TList.h> | |
14 | #include <TProfile.h> | |
15 | #include <TString.h> | |
16 | #include <TTree.h> | |
aea8063e | 17 | #include <TRandom3.h> |
3c40321c | 18 | |
04b116e8 | 19 | #include "AliAODEvent.h" |
cd0e5e33 | 20 | #include "AliAnalysisManager.h" |
21 | #include "AliAnalysisTaskSE.h" | |
cd0e5e33 | 22 | #include "AliEMCALGeometry.h" |
2fbb6418 | 23 | #include "AliEPFlattener.h" |
cd0e5e33 | 24 | #include "AliESDEvent.h" |
25 | #include "AliESDInputHandler.h" | |
26 | #include "AliESDtrackCuts.h" | |
27 | #include "AliEventplane.h" | |
28 | #include "AliMCEvent.h" | |
79ad78fd | 29 | #include "AliOADBContainer.h" |
cd0e5e33 | 30 | #include "AliStack.h" |
31 | #include "AliVCluster.h" | |
3c40321c | 32 | |
9d6804f4 | 33 | using std::cout; |
34 | using std::endl; | |
35 | ||
3c40321c | 36 | ClassImp(AliAnalysisTaskPi0V2) |
37 | ||
38 | //________________________________________________________________________ | |
cd0e5e33 | 39 | AliAnalysisTaskPi0V2::AliAnalysisTaskPi0V2(const char *name) : |
40 | AliAnalysisTaskSE(name), | |
41 | fOutput(0), | |
42 | fESD(0),fAOD(0), | |
43 | fTracksName("PicoTrack"), fV1ClusName("CaloCluster"), fV2ClusName("CaloCluster"), | |
44 | fTrigClass("CVLN_|CSEMI_|CCENT|CVHN"), | |
45 | fTracks(0), fV1Clus(0), fV2Clus(0), | |
46 | fRunNumber(-999),fInterRunNumber(-999), | |
47 | fVtxCut(15.), | |
3b7c7423 | 48 | fNcellCut(2.), fECut(1.), fEtaCut(0.65), fM02Cut(0.5),fDrCut(0.025), fPi0AsyCut(0), isV1Clus(1), isPhosCali(0), isCentFlat(0), isFullHist(0), |
cd0e5e33 | 49 | fCentrality(99.), |
50 | fEPTPC(-999.), | |
51 | fEPTPCreso(0.), | |
52 | fEPV0(-999.), fEPV0A(-999.), fEPV0C(-999.), fEPV0Ar(-999.), fEPV0Cr(-999.), fEPV0r(-999.), | |
53 | fEPV0AR4(-999.), fEPV0AR5(-999.), fEPV0AR6(-999.), fEPV0AR7(-999.), fEPV0CR0(-999.), fEPV0CR1(-999.), fEPV0CR2(-999.), fEPV0CR3(-999.), | |
900db6d2 | 54 | hEvtCount(0), hCent(0), |
cd0e5e33 | 55 | h2DcosV0A(0), h2DsinV0A(0), h2DcosV0C(0), h2DsinV0C(0), h2DcosTPC(0), h2DsinTPC(0), |
56 | hEPTPC(0), hresoTPC(0), | |
57 | hEPV0(0), hEPV0A(0), hEPV0C(0), hEPV0Ar(0), hEPV0Cr(0), hEPV0r(0), hEPV0AR4(0), hEPV0AR7(0), hEPV0CR0(0), hEPV0CR3(0), | |
58 | hEPTPCCor(0), hEPV0ACor(0), hEPV0CCor(0), | |
59 | hdifV0Ar_V0Cr(0), hdifV0A_V0CR0(0), hdifV0A_V0CR3(0), hdifV0ACR0_V0CR3(0), hdifV0C_V0AR4(0), hdifV0C_V0AR7(0), hdifV0AR4_V0AR7(0), | |
60 | hdifV0A_V0C(0), hdifV0A_TPC(0), hdifV0C_TPC(0), hdifV0C_V0A(0), | |
61 | hM02vsPtA(0), hM02vsPtB(0), hClusDxDZA(0), hClusDxDZB(0), | |
91d49d74 | 62 | hdifEMC_EPV0A(0), hdifEMC_EPV0C(0), |
63 | hdifful_EPV0A(0), hdifful_EPV0C(0), | |
64 | hdifout_EPV0A(0), hdifout_EPV0C(0), | |
65 | hCv2EMC_EPV0A(0), hCv2EMC_EPV0C(0), hCv2ful_EPV0A(0), hCv2ful_EPV0C(0), hCv2out_EPV0A(0), hCv2out_EPV0C(0), | |
66 | hclusDif_EPV0A(0), hclusDif_EPV0C(0), hclusv2_EPV0A(0), hclusv2_EPV0C(0), | |
cd0e5e33 | 67 | fEPcalibFileName("$ALICE_ROOT/OADB/PHOS/PHOSflat.root"), fTPCFlat(0x0), fV0AFlat(0x0), fV0CFlat(0x0), |
68 | fClusterPbV0(0), fClusterPbV0A(0), fClusterPbV0C(0), fClusterPbTPC(0), | |
d118f9e0 | 69 | fHEPV0A(0x0), fHEPV0C(0x0), fHEPTPC(0x0), |
70 | fHEPV0AM2(0x0), fHEPV0CM2(0x0), fHEPTPCM2(0x0) | |
3c40321c | 71 | { |
cd0e5e33 | 72 | // Dummy constructor ALWAYS needed for I/O. |
73 | DefineInput(0, TChain::Class()); | |
74 | DefineOutput(1, TList::Class()); // for output list | |
3c40321c | 75 | } |
76 | ||
77 | //________________________________________________________________________ | |
cd0e5e33 | 78 | AliAnalysisTaskPi0V2::AliAnalysisTaskPi0V2() : |
79 | AliAnalysisTaskSE("default_name"), | |
80 | fOutput(0), | |
81 | fESD(0),fAOD(0), | |
82 | fTracksName("PicoTrack"), fV1ClusName("CaloCluster"), fV2ClusName("CaloCluster"), | |
83 | fTrigClass("CVLN_|CSEMI_|CCENT|CVHN"), | |
84 | fTracks(0), fV1Clus(0), fV2Clus(0), | |
85 | fRunNumber(-999),fInterRunNumber(-999), | |
86 | fVtxCut(15.), | |
3b7c7423 | 87 | fNcellCut(2.), fECut(1.), fEtaCut(0.65), fM02Cut(0.5), fDrCut(0.025), fPi0AsyCut(0), isV1Clus(1),isPhosCali(0),isCentFlat(0), isFullHist(0), |
cd0e5e33 | 88 | fCentrality(99.), |
89 | fEPTPC(-999.), | |
90 | fEPTPCreso(0.), | |
91 | fEPV0(-999.), fEPV0A(-999.), fEPV0C(-999.), fEPV0Ar(-999.), fEPV0Cr(-999.), fEPV0r(-999.), | |
92 | fEPV0AR4(-999.), fEPV0AR5(-999.), fEPV0AR6(-999.), fEPV0AR7(-999.), fEPV0CR0(-999.), fEPV0CR1(-999.), fEPV0CR2(-999.), fEPV0CR3(-999.), | |
900db6d2 | 93 | hEvtCount(0), hCent(0), |
cd0e5e33 | 94 | h2DcosV0A(0), h2DsinV0A(0), h2DcosV0C(0), h2DsinV0C(0), h2DcosTPC(0), h2DsinTPC(0), |
95 | hEPTPC(0), hresoTPC(0), | |
96 | hEPV0(0), hEPV0A(0), hEPV0C(0), hEPV0Ar(0), hEPV0Cr(0), hEPV0r(0), hEPV0AR4(0), hEPV0AR7(0), hEPV0CR0(0), hEPV0CR3(0), | |
97 | hEPTPCCor(0), hEPV0ACor(0), hEPV0CCor(0), | |
98 | hdifV0Ar_V0Cr(0), hdifV0A_V0CR0(0), hdifV0A_V0CR3(0), hdifV0ACR0_V0CR3(0), hdifV0C_V0AR4(0), hdifV0C_V0AR7(0), hdifV0AR4_V0AR7(0), | |
99 | hdifV0A_V0C(0), hdifV0A_TPC(0), hdifV0C_TPC(0), hdifV0C_V0A(0), | |
100 | hM02vsPtA(0), hM02vsPtB(0), hClusDxDZA(0), hClusDxDZB(0), | |
91d49d74 | 101 | hdifEMC_EPV0A(0), hdifEMC_EPV0C(0), |
102 | hdifful_EPV0A(0), hdifful_EPV0C(0), | |
103 | hdifout_EPV0A(0), hdifout_EPV0C(0), | |
104 | hCv2EMC_EPV0A(0), hCv2EMC_EPV0C(0), hCv2ful_EPV0A(0), hCv2ful_EPV0C(0), hCv2out_EPV0A(0), hCv2out_EPV0C(0), | |
105 | hclusDif_EPV0A(0), hclusDif_EPV0C(0), hclusv2_EPV0A(0), hclusv2_EPV0C(0), | |
cd0e5e33 | 106 | fEPcalibFileName("$ALICE_ROOT/OADB/PHOS/PHOSflat.root"), fTPCFlat(0x0), fV0AFlat(0x0), fV0CFlat(0x0), |
107 | fClusterPbV0(0), fClusterPbV0A(0), fClusterPbV0C(0), fClusterPbTPC(0), | |
d118f9e0 | 108 | fHEPV0A(0x0), fHEPV0C(0x0), fHEPTPC(0x0), |
109 | fHEPV0AM2(0x0), fHEPV0CM2(0x0), fHEPTPCM2(0x0) | |
3c40321c | 110 | { |
cd0e5e33 | 111 | // Constructor |
112 | // Define input and output slots here (never in the dummy constructor) | |
113 | // Input slot #0 works with a TChain - it is connected to the default input container | |
114 | // Output slot #1 writes into a TH1 container | |
115 | DefineInput(0, TChain::Class()); | |
116 | DefineOutput(1, TList::Class()); // for output list | |
3c40321c | 117 | } |
118 | ||
119 | //________________________________________________________________________ | |
120 | AliAnalysisTaskPi0V2::~AliAnalysisTaskPi0V2() | |
121 | { | |
cd0e5e33 | 122 | // Destructor. Clean-up the output list, but not the histograms that are put inside |
123 | // (the list is owner and will clean-up these histograms). Protect in PROOF case. | |
124 | if (fTPCFlat) | |
125 | delete fTPCFlat; | |
126 | fTPCFlat=0x0; | |
127 | if (fV0AFlat) | |
128 | delete fV0AFlat; | |
129 | fV0AFlat=0x0; | |
130 | if (fV0CFlat) | |
131 | delete fV0CFlat; | |
132 | fV0CFlat=0x0; | |
133 | delete fOutput; | |
3c40321c | 134 | } |
cd0e5e33 | 135 | |
3c40321c | 136 | //_____________________________________________________________________ |
137 | Double_t AliAnalysisTaskPi0V2::GetMaxCellEnergy(const AliVCluster *cluster, Short_t &id) const | |
138 | { | |
139 | // Get maximum energy of attached cell. | |
140 | ||
141 | id = -1; | |
142 | ||
143 | AliVCaloCells *cells = 0; | |
cd0e5e33 | 144 | if (fESD) { |
3c40321c | 145 | cells = fESD->GetEMCALCells(); |
cd0e5e33 | 146 | } else { |
79ad78fd | 147 | cells = fAOD->GetEMCALCells(); |
148 | } | |
3c40321c | 149 | if (!cells) |
150 | return 0; | |
151 | ||
152 | Double_t maxe = 0; | |
153 | const Int_t ncells = cluster->GetNCells(); | |
154 | for (Int_t i=0; i<ncells; i++) { | |
155 | Double_t e = cells->GetCellAmplitude(TMath::Abs(cluster->GetCellAbsId(i))); | |
156 | if (e>maxe) { | |
157 | maxe = e; | |
158 | id = cluster->GetCellAbsId(i); | |
159 | } | |
160 | } | |
161 | return maxe; | |
162 | } | |
cd0e5e33 | 163 | |
3c40321c | 164 | //_____________________________________________________________________ |
165 | Double_t AliAnalysisTaskPi0V2::GetCrossEnergy(const AliVCluster *cluster, Short_t &idmax) const | |
166 | { | |
167 | // Calculate the energy of cross cells around the leading cell. | |
168 | ||
79ad78fd | 169 | AliVCaloCells *cells; |
cd0e5e33 | 170 | if (fESD) { |
3c40321c | 171 | cells = fESD->GetEMCALCells(); |
cd0e5e33 | 172 | } else { |
79ad78fd | 173 | cells = fAOD->GetEMCALCells(); |
174 | } | |
3c40321c | 175 | if (!cells) |
176 | return 0; | |
177 | ||
178 | AliEMCALGeometry *geom = AliEMCALGeometry::GetInstance(); | |
179 | if (!geom) | |
180 | return 0; | |
181 | ||
182 | Int_t iSupMod = -1; | |
183 | Int_t iTower = -1; | |
184 | Int_t iIphi = -1; | |
185 | Int_t iIeta = -1; | |
186 | Int_t iphi = -1; | |
187 | Int_t ieta = -1; | |
188 | Int_t iphis = -1; | |
189 | Int_t ietas = -1; | |
190 | ||
965c985f | 191 | Double_t crossEnergy = 0.; |
3c40321c | 192 | |
193 | geom->GetCellIndex(idmax,iSupMod,iTower,iIphi,iIeta); | |
194 | geom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi, iIeta,iphis,ietas); | |
195 | ||
196 | Int_t ncells = cluster->GetNCells(); | |
197 | for (Int_t i=0; i<ncells; i++) { | |
198 | Int_t cellAbsId = cluster->GetCellAbsId(i); | |
199 | geom->GetCellIndex(cellAbsId,iSupMod,iTower,iIphi,iIeta); | |
200 | geom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi, iIeta,iphi,ieta); | |
201 | Int_t aphidiff = TMath::Abs(iphi-iphis); | |
202 | if (aphidiff>1) | |
203 | continue; | |
204 | Int_t aetadiff = TMath::Abs(ieta-ietas); | |
205 | if (aetadiff>1) | |
206 | continue; | |
207 | if ( (aphidiff==1 && aetadiff==0) || | |
208 | (aphidiff==0 && aetadiff==1) ) { | |
209 | crossEnergy += cells->GetCellAmplitude(cellAbsId); | |
210 | } | |
211 | } | |
212 | ||
213 | return crossEnergy; | |
214 | } | |
cd0e5e33 | 215 | |
3c40321c | 216 | //_____________________________________________________________________ |
217 | Bool_t AliAnalysisTaskPi0V2::IsWithinFiducialVolume(Short_t id) const | |
218 | { | |
219 | // Check if cell is within given fiducial volume. | |
220 | ||
221 | Double_t fNFiducial = 1; | |
222 | ||
223 | Int_t iSupMod = -1; | |
224 | Int_t iTower = -1; | |
225 | Int_t iIphi = -1; | |
226 | Int_t iIeta = -1; | |
227 | Int_t iphi = -1; | |
228 | Int_t ieta = -1; | |
229 | ||
230 | Bool_t okrow = kFALSE; | |
231 | Bool_t okcol = kFALSE; | |
232 | ||
04b116e8 | 233 | AliEMCALGeometry *geom = AliEMCALGeometry::GetInstance(); |
3c40321c | 234 | if (!geom) |
235 | return kFALSE; | |
236 | ||
237 | Int_t cellAbsId = id; | |
238 | geom->GetCellIndex(cellAbsId,iSupMod,iTower,iIphi,iIeta); | |
239 | geom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi, iIeta,iphi,ieta); | |
240 | ||
241 | // Check rows/phi | |
242 | if (iSupMod < 10) { | |
243 | if (iphi >= fNFiducial && iphi < 24-fNFiducial) | |
244 | okrow = kTRUE; | |
245 | } else { | |
246 | if (iphi >= fNFiducial && iphi < 12-fNFiducial) | |
247 | okrow = kTRUE; | |
248 | } | |
249 | // Check columns/eta | |
250 | Bool_t noEMCALBorderAtEta0 = kTRUE; | |
251 | if (!noEMCALBorderAtEta0) { | |
252 | if (ieta > fNFiducial && ieta < 48-fNFiducial) | |
253 | okcol = kTRUE; | |
254 | } else { | |
255 | if (iSupMod%2==0) { | |
256 | if (ieta >= fNFiducial) | |
257 | okcol = kTRUE; | |
258 | } else { | |
259 | if (ieta < 48-fNFiducial) | |
260 | okcol = kTRUE; | |
261 | } | |
262 | } | |
263 | if (okrow && okcol) | |
264 | return kTRUE; | |
265 | ||
266 | return kFALSE; | |
267 | } | |
cd0e5e33 | 268 | |
3c40321c | 269 | //______________________________________________________________________ |
ceba2d0c | 270 | Bool_t AliAnalysisTaskPi0V2::IsGoodCluster(const AliVCluster *c) const |
3c40321c | 271 | { |
cd0e5e33 | 272 | if (!c) |
3c40321c | 273 | return kFALSE; |
274 | ||
670ffa5c | 275 | if(c->GetNCells() < fNcellCut) |
3c40321c | 276 | return kFALSE; |
277 | ||
670ffa5c | 278 | if(c->E() < fECut) |
3c40321c | 279 | return kFALSE; |
3c40321c | 280 | Short_t id = -1; |
281 | Double_t maxE = GetMaxCellEnergy(c, id); | |
965c985f | 282 | if((1. - double(GetCrossEnergy(c,id))/maxE) > 0.97) |
3c40321c | 283 | return kFALSE; |
284 | ||
965c985f | 285 | |
3c40321c | 286 | Float_t pos1[3]; |
287 | c->GetPosition(pos1); | |
288 | TVector3 clsPos(pos1); | |
289 | Double_t eta = clsPos.Eta(); | |
290 | ||
cd0e5e33 | 291 | if (TMath::Abs(eta) > fEtaCut) |
3c40321c | 292 | return kFALSE; |
293 | ||
294 | if (!IsWithinFiducialVolume(id)) | |
295 | return kFALSE; | |
296 | ||
670ffa5c | 297 | if(c->GetM02() >fM02Cut) |
3c40321c | 298 | return kFALSE; |
299 | ||
04b116e8 | 300 | |
3c40321c | 301 | return kTRUE; |
70d53162 | 302 | |
cc57d293 | 303 | } |
304 | //________________________________________________________________________________________________ | |
ceba2d0c | 305 | Bool_t AliAnalysisTaskPi0V2::IsGoodClusterV1(const AliVCluster *c) const |
cc57d293 | 306 | { |
cd0e5e33 | 307 | if (!c) |
cc57d293 | 308 | return kFALSE; |
309 | ||
cd0e5e33 | 310 | if (c->GetNCells() < fNcellCut) |
cc57d293 | 311 | return kFALSE; |
312 | ||
cd0e5e33 | 313 | if (c->E() < fECut) |
cc57d293 | 314 | return kFALSE; |
315 | ||
316 | Short_t id = -1; | |
317 | Double_t maxE = GetMaxCellEnergy(c, id); | |
9ff4e38b | 318 | if((1. - double(GetCrossEnergy(c,id))/maxE) > 0.97) |
cc57d293 | 319 | return kFALSE; |
320 | ||
321 | ||
322 | Float_t pos1[3]; | |
323 | c->GetPosition(pos1); | |
324 | TVector3 clsPos(pos1); | |
325 | Double_t eta = clsPos.Eta(); | |
326 | ||
cd0e5e33 | 327 | if (TMath::Abs(eta) > fEtaCut) |
cc57d293 | 328 | return kFALSE; |
329 | ||
330 | if (!IsWithinFiducialVolume(id)) | |
331 | return kFALSE; | |
332 | ||
cd0e5e33 | 333 | if (c->GetM02() <fM02Cut) |
cc57d293 | 334 | return kFALSE; |
335 | ||
336 | Double_t dr = TMath::Sqrt(c->GetTrackDx()*c->GetTrackDx() + c->GetTrackDz()*c->GetTrackDz()); | |
c45a0734 | 337 | if(dr<fDrCut) |
cc57d293 | 338 | return kFALSE; |
339 | ||
340 | return kTRUE; | |
04b116e8 | 341 | } |
cd0e5e33 | 342 | |
3c40321c | 343 | //_____________________________________________________________________ |
344 | Bool_t AliAnalysisTaskPi0V2::IsGoodPion(const TLorentzVector &p1, const TLorentzVector &p2) const | |
345 | { | |
346 | // Is good pion? | |
347 | ||
670ffa5c | 348 | if(fPi0AsyCut){ |
349 | Double_t asym = TMath::Abs(p1.E()-p2.E())/(p1.E()+p2.E()); | |
350 | if (asym>0.7) | |
351 | return kFALSE; | |
352 | } | |
3c40321c | 353 | TLorentzVector pion; |
354 | pion = p1 + p2; | |
355 | Double_t eta = pion.Eta(); | |
670ffa5c | 356 | if(TMath::Abs(eta) > fEtaCut) |
3c40321c | 357 | return kFALSE; |
358 | ||
359 | return kTRUE; | |
360 | } | |
cd0e5e33 | 361 | |
3c40321c | 362 | //_______________________________________________________________________ |
79ad78fd | 363 | void AliAnalysisTaskPi0V2::FillPion(const TLorentzVector& p1, const TLorentzVector& p2, Double_t EPV0A, Double_t EPV0C, Double_t EPTPC) |
3c40321c | 364 | { |
365 | // Fill histogram. | |
366 | ||
367 | if (!IsGoodPion(p1,p2)) | |
368 | return; | |
369 | TLorentzVector pion; | |
370 | pion = p1 + p2; | |
371 | ||
372 | Double_t mass = pion.M(); | |
373 | Double_t pt = pion.Pt(); | |
374 | Double_t phi = pion.Phi(); | |
375 | ||
3c40321c | 376 | Double_t dphiV0A = phi-EPV0A; |
377 | Double_t dphiV0C = phi-EPV0C; | |
378 | Double_t dphiTPC = phi-EPTPC; | |
379 | ||
1accedbd | 380 | Double_t cos2phiV0A = TMath::Cos(2.*(dphiV0A)); |
381 | Double_t cos2phiV0C = TMath::Cos(2.*(dphiV0C)); | |
382 | Double_t cos2phiTPC = TMath::Cos(2.*(dphiTPC)); | |
383 | ||
3b7c7423 | 384 | while(dphiV0A<0.) dphiV0A+=TMath::Pi(); while(dphiV0A>TMath::Pi()) dphiV0A-=TMath::Pi(); |
385 | while(dphiV0C<0.) dphiV0C+=TMath::Pi(); while(dphiV0C>TMath::Pi()) dphiV0C-=TMath::Pi(); | |
386 | while(dphiTPC<0.) dphiTPC+=TMath::Pi(); while(dphiTPC>TMath::Pi()) dphiTPC-=TMath::Pi(); | |
387 | ||
3c40321c | 388 | |
d118f9e0 | 389 | Double_t xV0A[4]; // Match ndims in fH V0A EP for method 1 |
3c40321c | 390 | xV0A[0] = mass; |
391 | xV0A[1] = pt; | |
392 | xV0A[2] = fCentrality; | |
393 | xV0A[3] = dphiV0A; | |
3c40321c | 394 | fHEPV0A->Fill(xV0A); |
395 | ||
d118f9e0 | 396 | |
397 | Double_t xV0AM2[4]; // Match ndims in fH V0A EP for method 2 | |
398 | xV0AM2[0] = mass; | |
399 | xV0AM2[1] = pt; | |
400 | xV0AM2[2] = fCentrality; | |
401 | xV0AM2[3] = cos2phiV0A; | |
402 | fHEPV0AM2->Fill(xV0AM2); | |
403 | ||
404 | ||
405 | Double_t xV0C[4]; // Match ndims in fH V0C EP for method 1 | |
3c40321c | 406 | xV0C[0] = mass; |
407 | xV0C[1] = pt; | |
408 | xV0C[2] = fCentrality; | |
409 | xV0C[3] = dphiV0C; | |
3c40321c | 410 | fHEPV0C->Fill(xV0C); |
411 | ||
d118f9e0 | 412 | Double_t xV0CM2[4]; // Match ndims in fH V0C EP for method 2 |
413 | xV0CM2[0] = mass; | |
414 | xV0CM2[1] = pt; | |
415 | xV0CM2[2] = fCentrality; | |
416 | xV0CM2[3] = cos2phiV0C; | |
417 | fHEPV0CM2->Fill(xV0CM2); | |
418 | ||
419 | ||
af1defae | 420 | if (fEPTPC!=-999.){ |
d118f9e0 | 421 | Double_t xTPC[4]; // Match ndims in fH TPC EP for method 1 |
af1defae | 422 | xTPC[0] = mass; |
423 | xTPC[1] = pt; | |
424 | xTPC[2] = fCentrality; | |
91d49d74 | 425 | xTPC[3] = dphiTPC; |
af1defae | 426 | fHEPTPC->Fill(xTPC); |
d118f9e0 | 427 | |
428 | Double_t xTPCM2[4]; // Match ndims in fH TPC EP | |
429 | xTPCM2[0] = mass; | |
430 | xTPCM2[1] = pt; | |
431 | xTPCM2[2] = fCentrality; | |
432 | xTPCM2[3] = cos2phiTPC; | |
433 | fHEPTPCM2->Fill(xTPCM2); | |
434 | ||
af1defae | 435 | } |
7bb82cc2 | 436 | } |
cd0e5e33 | 437 | |
e50db689 | 438 | //________________________________________________________________________________________________________________________________ |
79ad78fd | 439 | void AliAnalysisTaskPi0V2::FillCluster(const TLorentzVector& p1, Double_t EPV0A, Double_t EPV0C, Double_t EPTPC, AliVCluster *c) |
7bb82cc2 | 440 | { |
cd0e5e33 | 441 | // Cluster(photon) v2 method |
442 | ||
8e5ac2d1 | 443 | Double_t Et = p1.Et(); |
7bb82cc2 | 444 | Double_t Phi = p1.Phi(); |
9ff4e38b | 445 | Double_t M02 = c->GetM02(); |
7bb82cc2 | 446 | |
3b7c7423 | 447 | Double_t difClusV0A = Phi-EPV0A; |
448 | Double_t difClusV0C = Phi-EPV0C; | |
449 | Double_t difClusTPC = Phi-EPTPC; | |
450 | while(difClusV0A<0.) difClusV0A+=TMath::Pi(); while(difClusV0A>TMath::Pi()) difClusV0A-=TMath::Pi(); | |
451 | while(difClusV0C<0.) difClusV0C+=TMath::Pi(); while(difClusV0C>TMath::Pi()) difClusV0C-=TMath::Pi(); | |
452 | while(difClusTPC<0.) difClusTPC+=TMath::Pi(); while(difClusTPC>TMath::Pi()) difClusTPC-=TMath::Pi(); | |
7bb82cc2 | 453 | |
900db6d2 | 454 | Double_t DataV0A[4]; |
8e5ac2d1 | 455 | DataV0A[0] = Et; |
456 | DataV0A[1] = M02; | |
457 | DataV0A[2] = fCentrality; | |
79ad78fd | 458 | DataV0A[3] = difClusV0A; |
8e5ac2d1 | 459 | fClusterPbV0A->Fill(DataV0A); |
460 | ||
900db6d2 | 461 | Double_t DataV0C[4]; |
8e5ac2d1 | 462 | DataV0C[0] = Et; |
463 | DataV0C[1] = M02; | |
464 | DataV0C[2] = fCentrality; | |
79ad78fd | 465 | DataV0C[3] = difClusV0C; |
8e5ac2d1 | 466 | fClusterPbV0C->Fill(DataV0C); |
467 | ||
900db6d2 | 468 | Double_t DataTPC[4]; |
8e5ac2d1 | 469 | DataTPC[0] = Et; |
470 | DataTPC[1] = M02; | |
471 | DataTPC[2] = fCentrality; | |
79ad78fd | 472 | DataTPC[3] = difClusTPC; |
8e5ac2d1 | 473 | fClusterPbTPC->Fill(DataTPC); |
04b116e8 | 474 | } |
cd0e5e33 | 475 | |
3c40321c | 476 | //_________________________________________________________________________________________________ |
ceba2d0c | 477 | void AliAnalysisTaskPi0V2::GetMom(TLorentzVector& p, const AliVCluster *c, Double_t *vertex) |
3c40321c | 478 | { |
479 | // Calculate momentum. | |
480 | Float_t posMom[3]; | |
481 | c->GetPosition(posMom); | |
482 | TVector3 clsPos2(posMom); | |
483 | ||
484 | Double_t e = c->E(); | |
485 | Double_t r = clsPos2.Perp(); | |
486 | Double_t eta = clsPos2.Eta(); | |
487 | Double_t phi = clsPos2.Phi(); | |
488 | ||
489 | TVector3 pos; | |
490 | pos.SetPtEtaPhi(r,eta,phi); | |
491 | ||
492 | if (vertex) { //calculate direction relative to vertex | |
493 | pos -= vertex; | |
494 | } | |
495 | ||
496 | Double_t rad = pos.Mag(); | |
497 | p.SetPxPyPzE(e*pos.x()/rad, e*pos.y()/rad, e*pos.z()/rad, e); | |
70d53162 | 498 | |
04b116e8 | 499 | } |
cd0e5e33 | 500 | |
3c40321c | 501 | //________________________________________________________________________ |
502 | void AliAnalysisTaskPi0V2::UserCreateOutputObjects() | |
503 | { | |
cd0e5e33 | 504 | // Create histograms |
505 | // Called once (on the worker node) | |
3c40321c | 506 | |
cd0e5e33 | 507 | fOutput = new TList(); |
508 | fOutput->SetOwner(); // IMPORTANT! | |
509 | ||
510 | hEvtCount = new TH1F("hEvtCount", " Event Plane", 9, 0.5, 9.5); | |
511 | hEvtCount->GetXaxis()->SetBinLabel(1,"All"); | |
512 | hEvtCount->GetXaxis()->SetBinLabel(2,"Evt"); | |
513 | hEvtCount->GetXaxis()->SetBinLabel(3,"Trg Class"); | |
514 | hEvtCount->GetXaxis()->SetBinLabel(4,"Vtx"); | |
515 | hEvtCount->GetXaxis()->SetBinLabel(5,"Cent"); | |
516 | hEvtCount->GetXaxis()->SetBinLabel(6,"EPtask"); | |
517 | hEvtCount->GetXaxis()->SetBinLabel(7,"ClusterTask"); | |
518 | hEvtCount->GetXaxis()->SetBinLabel(8,"Pass"); | |
519 | fOutput->Add(hEvtCount); | |
900db6d2 | 520 | |
521 | hCent = new TH1F("hCent", "centrality dist. before App. flat cut", 100, 0., 100.); | |
522 | fOutput->Add(hCent); | |
523 | ||
cd0e5e33 | 524 | hEPTPC = new TH2F("hEPTPC", "EPTPC vs cent", 100, 0., 100., 100, 0., TMath::Pi()); |
525 | hresoTPC = new TH2F("hresoTPC", "TPc reso vs cent", 100, 0., 100., 100, 0., 1.); | |
cd0e5e33 | 526 | hEPV0A = new TH2F("hEPV0A", "EPV0A vs cent", 100, 0., 100., 100, 0., TMath::Pi()); |
527 | hEPV0C = new TH2F("hEPV0C", "EPV0C vs cent", 100, 0., 100., 100, 0., TMath::Pi()); | |
cd0e5e33 | 528 | fOutput->Add(hEPTPC); |
529 | fOutput->Add(hresoTPC); | |
cd0e5e33 | 530 | fOutput->Add(hEPV0A); |
531 | fOutput->Add(hEPV0C); | |
3b7c7423 | 532 | |
533 | if(isFullHist){ | |
534 | hEPV0 = new TH2F("hEPV0", "EPV0 vs cent", 100, 0., 100., 100, 0., TMath::Pi()); | |
535 | hEPV0Ar = new TH2F("hEPV0Ar", "EPV0Ar vs cent", 100, 0., 100., 100, 0., TMath::Pi()); | |
536 | hEPV0Cr = new TH2F("hEPV0Cr", "EPV0Cr vs cent", 100, 0., 100., 100, 0., TMath::Pi()); | |
537 | hEPV0r = new TH2F("hEPV0r", "EPV0r vs cent", 100, 0., 100., 100, 0., TMath::Pi()); | |
538 | hEPV0AR4 = new TH2F("hEPV0AR4", "EPV0AR4 vs cent", 100, 0., 100., 100, 0., TMath::Pi()); | |
539 | hEPV0AR7 = new TH2F("hEPV0AR7", "EPV0AR7 vs cent", 100, 0., 100., 100, 0., TMath::Pi()); | |
540 | hEPV0CR0 = new TH2F("hEPV0CR0", "EPV0CR0 vs cent", 100, 0., 100., 100, 0., TMath::Pi()); | |
541 | hEPV0CR3 = new TH2F("hEPV0CR3", "EPV0CR3 vs cent", 100, 0., 100., 100, 0., TMath::Pi()); | |
542 | fOutput->Add(hEPV0); | |
543 | fOutput->Add(hEPV0Ar); | |
544 | fOutput->Add(hEPV0Cr); | |
545 | fOutput->Add(hEPV0r); | |
546 | fOutput->Add(hEPV0AR4); | |
547 | fOutput->Add(hEPV0AR7); | |
548 | fOutput->Add(hEPV0CR0); | |
549 | fOutput->Add(hEPV0CR3); | |
550 | } | |
cd0e5e33 | 551 | |
552 | hEPTPCCor = new TH2F("hEPTPCCor", "EPTPC vs cent after PHOS Correct", 100, 0., 100., 100, 0., TMath::Pi()); | |
553 | hEPV0ACor = new TH2F("hEPV0ACor", "EPV0A vs cent after PHOS Correct", 100, 0., 100., 100, 0., TMath::Pi()); | |
554 | hEPV0CCor = new TH2F("hEPV0CCor", "EPV0C vs cent after PHOS Correct", 100, 0., 100., 100, 0., TMath::Pi()); | |
555 | fOutput->Add(hEPTPCCor); | |
556 | fOutput->Add(hEPV0ACor); | |
557 | fOutput->Add(hEPV0CCor); | |
558 | ||
cd0e5e33 | 559 | hdifV0A_V0CR0 = new TH2F("hdifV0A_V0CR0", "EP A-R0 ", 100, 0., 100., 100, -1., 1.); |
560 | hdifV0A_V0CR3 = new TH2F("hdifV0A_V0CR3", "EP A-R3 ", 100, 0., 100., 100, -1., 1.); | |
561 | hdifV0ACR0_V0CR3 = new TH2F("hdifV0ACR0_V0CR3", "EP R0-R3 ", 100, 0., 100., 100, -1., 1.); | |
562 | hdifV0C_V0AR4 = new TH2F("hdifV0C_V0AR4", "EP C-R4 ", 100, 0., 100., 100, -1., 1.); | |
563 | hdifV0C_V0AR7 = new TH2F("hdifV0C_V0AR7", "EP C-R7 ", 100, 0., 100., 100, -1., 1.); | |
564 | hdifV0AR4_V0AR7 = new TH2F("hdifV0AR4_V0AR7", "EP R4-R7 ", 100, 0., 100., 100, -1., 1.); | |
cd0e5e33 | 565 | fOutput->Add(hdifV0A_V0CR0); |
566 | fOutput->Add(hdifV0A_V0CR3); | |
567 | fOutput->Add(hdifV0ACR0_V0CR3); | |
568 | fOutput->Add(hdifV0C_V0AR4); | |
569 | fOutput->Add(hdifV0C_V0AR7); | |
570 | fOutput->Add(hdifV0AR4_V0AR7); | |
571 | ||
3b7c7423 | 572 | if(isFullHist){ |
573 | hdifV0Ar_V0Cr = new TH2F("hdifV0Ar_V0Cr", "EP Ar-Cr ", 100, 0., 100., 100, -1., 1.); | |
574 | fOutput->Add(hdifV0Ar_V0Cr); | |
575 | ||
576 | hdifV0A_V0C = new TH2F("hdifV0A_V0C", "EP A-C ", 100, 0., 100., 100, -1., 1.); | |
577 | hdifV0A_TPC = new TH2F("hdifV0A_TPC", "EP A-TPC", 100, 0., 100., 100, -1., 1.); | |
578 | hdifV0C_TPC = new TH2F("hdifV0C_TPC", "EP C-TPC", 100, 0., 100., 100, -1., 1.); | |
579 | hdifV0C_V0A = new TH2F("hdifV0C_V0A", "EP C-A ", 100, 0., 100., 100, -1., 1.); | |
580 | fOutput->Add(hdifV0A_V0C); | |
581 | fOutput->Add(hdifV0A_TPC); | |
582 | fOutput->Add(hdifV0C_TPC); | |
583 | fOutput->Add(hdifV0C_V0A); | |
584 | } | |
cd0e5e33 | 585 | |
91d49d74 | 586 | hdifEMC_EPV0A = new TH3F("hdifEMC_EPV0A", "dif phi in EMC with EPV0A", 100, 0., 100., 100, 0., TMath::Pi(), 15, 0., 15.); |
587 | hdifEMC_EPV0C = new TH3F("hdifEMC_EPV0C", "dif phi in EMC with EPV0C", 100, 0., 100., 100, 0., TMath::Pi(), 15, 0., 15.); | |
cd0e5e33 | 588 | fOutput->Add(hdifEMC_EPV0A); |
589 | fOutput->Add(hdifEMC_EPV0C); | |
590 | ||
91d49d74 | 591 | hdifful_EPV0A = new TH3F("hdifful_EPV0A", "dif phi in full with EPV0A", 100, 0., 100., 100, 0., TMath::Pi(), 15, 0., 15.); |
592 | hdifful_EPV0C = new TH3F("hdifful_EPV0C", "dif phi in full with EPV0C", 100, 0., 100., 100, 0., TMath::Pi(), 15, 0., 15.); | |
cd0e5e33 | 593 | fOutput->Add(hdifful_EPV0A); |
594 | fOutput->Add(hdifful_EPV0C); | |
595 | ||
91d49d74 | 596 | hdifout_EPV0A = new TH3F("hdifout_EPV0A", "dif phi NOT in EMC with EPV0A", 100, 0., 100., 100, 0., TMath::Pi(), 15, 0., 15.); |
597 | hdifout_EPV0C = new TH3F("hdifout_EPV0C", "dif phi NOT in EMC with EPV0C", 100, 0., 100., 100, 0., TMath::Pi(), 15, 0., 15.); | |
cd0e5e33 | 598 | fOutput->Add(hdifout_EPV0A); |
599 | fOutput->Add(hdifout_EPV0C); | |
600 | ||
91d49d74 | 601 | hCv2EMC_EPV0A = new TH3F("hCv2EMC_EPV0A", " raw v2 of charged trc in EMC with V0A", 100, 0, 100, 50, -1., 1., 15, 0., 15.); |
602 | hCv2EMC_EPV0C = new TH3F("hCv2EMC_EPV0C", " raw v2 of charged trc in EMC with V0C", 100, 0, 100, 50, -1., 1., 15, 0., 15.); | |
603 | fOutput->Add(hCv2EMC_EPV0A); | |
604 | fOutput->Add(hCv2EMC_EPV0C); | |
605 | ||
606 | hCv2ful_EPV0A = new TH3F("hCv2ful_EPV0A", " raw v2 of charged trc in ful with V0A", 100, 0, 100, 50, -1., 1., 15, 0., 15.); | |
607 | hCv2ful_EPV0C = new TH3F("hCv2ful_EPV0C", " raw v2 of charged trc in ful with V0C", 100, 0, 100, 50, -1., 1., 15, 0., 15.); | |
608 | fOutput->Add(hCv2ful_EPV0A); | |
609 | fOutput->Add(hCv2ful_EPV0C); | |
610 | ||
611 | hCv2out_EPV0A = new TH3F("hCv2out_EPV0A", " raw v2 of charged trc out with V0A", 100, 0, 100, 50, -1., 1., 15, 0., 15.); | |
612 | hCv2out_EPV0C = new TH3F("hCv2out_EPV0C", " raw v2 of charged trc out with V0A", 100, 0, 100, 50, -1., 1., 15, 0., 15.); | |
613 | fOutput->Add(hCv2out_EPV0A); | |
614 | fOutput->Add(hCv2out_EPV0C); | |
615 | ||
616 | hclusDif_EPV0A = new TH3F("hclusDif_EPV0A", "dif phi of clus with EP V0A", 100, 0., 100., 100, 0., TMath::Pi(), 15, 0., 15.); | |
617 | hclusDif_EPV0C = new TH3F("hclusDif_EPV0C", "dif phi of clus with EP V0C", 100, 0., 100., 100, 0., TMath::Pi(), 15, 0., 15.); | |
618 | fOutput->Add(hclusDif_EPV0A); | |
619 | fOutput->Add(hclusDif_EPV0C); | |
620 | ||
621 | hclusv2_EPV0A = new TH3F("hclusv2_EPV0A", " raw v2 of clus in ful with V0A", 100, 0, 100, 50, -1., 1., 15, 0., 15.); | |
622 | hclusv2_EPV0C = new TH3F("hclusv2_EPV0C", " raw v2 of clus in ful with V0C", 100, 0, 100, 50, -1., 1., 15, 0., 15.); | |
623 | fOutput->Add(hclusv2_EPV0A); | |
624 | fOutput->Add(hclusv2_EPV0C); | |
625 | ||
cd0e5e33 | 626 | if (isV1Clus) { |
900db6d2 | 627 | // Et M02 spdcent DeltaPhi |
3b7c7423 | 628 | Int_t bins[4] = { 40, 350, 60, 100 }; // binning |
900db6d2 | 629 | Double_t min[4] = { 0.0, 0.0, 0, 0.0 }; // min x |
3b7c7423 | 630 | Double_t max[4] = { 40.0, 3.5, 60, TMath::Pi()}; // max x |
8e5ac2d1 | 631 | |
79ad78fd | 632 | fClusterPbV0A = new THnSparseF("fClusterPbV0A","",5,bins,min,max); |
cd0e5e33 | 633 | fClusterPbV0A->GetAxis(0)->SetTitle("Transverse Energy [GeV]"); |
634 | fClusterPbV0A->GetAxis(1)->SetTitle("M02"); | |
635 | fClusterPbV0A->GetAxis(2)->SetTitle("V0M Centrality"); | |
636 | fClusterPbV0A->GetAxis(3)->SetTitle("Delta(#phi) [rad]"); | |
8e5ac2d1 | 637 | fOutput->Add(fClusterPbV0A); |
638 | ||
79ad78fd | 639 | fClusterPbV0C = new THnSparseF("fClusterPbV0C","",5,bins,min,max); |
cd0e5e33 | 640 | fClusterPbV0C->GetAxis(0)->SetTitle("Transverse Energy [GeV]"); |
641 | fClusterPbV0C->GetAxis(1)->SetTitle("M02"); | |
642 | fClusterPbV0C->GetAxis(2)->SetTitle("V0M Centrality"); | |
643 | fClusterPbV0C->GetAxis(3)->SetTitle("Delta(#phi) [rad]"); | |
8e5ac2d1 | 644 | fOutput->Add(fClusterPbV0C); |
645 | ||
79ad78fd | 646 | fClusterPbTPC = new THnSparseF("fClusterPbTPC","",5,bins,min,max); |
cd0e5e33 | 647 | fClusterPbTPC->GetAxis(0)->SetTitle("Transverse Energy [GeV]"); |
648 | fClusterPbTPC->GetAxis(1)->SetTitle("M02"); | |
649 | fClusterPbTPC->GetAxis(2)->SetTitle("V0M Centrality"); | |
650 | fClusterPbTPC->GetAxis(3)->SetTitle("Delta(#phi) [rad]"); | |
8e5ac2d1 | 651 | fOutput->Add(fClusterPbTPC); |
cd0e5e33 | 652 | } |
79ad78fd | 653 | |
3b7c7423 | 654 | if(isFullHist){ |
655 | h2DcosV0C = new TProfile("h2DcosV0C", "cos(Phi) V0r vs Run NUmber", 200, 0., 200.); | |
656 | h2DsinV0C = new TProfile("h2DsinV0C", "sin(Phi) V0r vs Run NUmber", 200, 0., 200.); | |
657 | h2DcosTPC = new TProfile("h2DcosTPC", "cos(Phi) V0r vs Run NUmber", 200, 0., 200.); | |
658 | h2DsinTPC = new TProfile("h2DsinTPC", "sin(Phi) V0r vs Run NUmber", 200, 0., 200.); | |
659 | fOutput->Add(h2DcosV0C); | |
660 | fOutput->Add(h2DsinV0C); | |
661 | fOutput->Add(h2DcosTPC); | |
662 | fOutput->Add(h2DsinTPC); | |
663 | } | |
664 | ||
cd0e5e33 | 665 | h2DcosV0A = new TProfile("h2DcosV0A", "cos(Phi) V0r vs Run NUmber", 200, 0., 200.); |
666 | h2DsinV0A = new TProfile("h2DsinV0A", "sin(Phi) V0r vs Run NUmber", 200, 0., 200.); | |
cd0e5e33 | 667 | fOutput->Add(h2DcosV0A); |
668 | fOutput->Add(h2DsinV0A); | |
cd0e5e33 | 669 | |
670 | if (isV1Clus) { | |
671 | hM02vsPtA = new TH2F("hM02vsPtA", "M02 vs Et before cut", 5000, 0, 50, 400, 0, 4.); | |
672 | hM02vsPtB = new TH2F("hM02vsPtB", "M02 vs Et before cut", 5000, 0, 50, 400, 0, 4.); | |
673 | fOutput->Add(hM02vsPtA); | |
674 | fOutput->Add(hM02vsPtB); | |
675 | } | |
3b7c7423 | 676 | if(isFullHist){ |
677 | hClusDxDZA = new TH2F("hClusDxDZA", "clus Dx vs Dz", 1000, -1., 1., 1000, -1., 1); | |
678 | hClusDxDZB = new TH2F("hClusDxDZB", "clus Dx vs Dz", 1000, -1., 1., 1000, -1., 1); | |
679 | fOutput->Add(hClusDxDZA); | |
680 | fOutput->Add(hClusDxDZB); | |
681 | } | |
682 | ||
cd0e5e33 | 683 | if (!isV1Clus) { |
d118f9e0 | 684 | const Int_t ndims = 4; |
3b7c7423 | 685 | Int_t nMgg=500, nPt=40, nCent=20, nDeltaPhi=315, ncos2phi=200; |
d118f9e0 | 686 | Int_t binsv1[ndims] = {nMgg, nPt, nCent, nDeltaPhi}; |
91d49d74 | 687 | Double_t xmin[ndims] = { 0, 0., 0, 0. }; |
688 | Double_t xmax[ndims] = { 0.5, 20., 100, 3.15 }; | |
8e5ac2d1 | 689 | fHEPV0A = new THnSparseF("fHEPV0A", "Flow histogram EPV0A", ndims, binsv1, xmin, xmax); |
690 | fHEPV0C = new THnSparseF("fHEPV0C", "Flow histogram EPV0C", ndims, binsv1, xmin, xmax); | |
691 | fHEPTPC = new THnSparseF("fHEPTPC", "Flow histogram EPTPC", ndims, binsv1, xmin, xmax); | |
cd0e5e33 | 692 | fHEPV0A->GetAxis(0)->SetTitle("m_{#gamma#gamma} "); |
693 | fHEPV0A->GetAxis(1)->SetTitle("p_{T}[GeV]"); | |
694 | fHEPV0A->GetAxis(2)->SetTitle("centrality"); | |
695 | fHEPV0A->GetAxis(3)->SetTitle("#delta #phi"); | |
696 | fHEPV0C->GetAxis(0)->SetTitle("m_{#gamma#gamma} "); | |
697 | fHEPV0C->GetAxis(1)->SetTitle("p_{T}[GeV]"); | |
698 | fHEPV0C->GetAxis(2)->SetTitle("centrality"); | |
699 | fHEPV0C->GetAxis(3)->SetTitle("#delta #phi"); | |
700 | fHEPTPC->GetAxis(0)->SetTitle("m_{#gamma#gamma} "); | |
701 | fHEPTPC->GetAxis(1)->SetTitle("p_{T}[GeV]"); | |
702 | fHEPTPC->GetAxis(2)->SetTitle("centrality"); | |
703 | fHEPTPC->GetAxis(3)->SetTitle("#delta #phi"); | |
04b116e8 | 704 | fOutput->Add(fHEPV0A); |
705 | fOutput->Add(fHEPV0C); | |
706 | fOutput->Add(fHEPTPC); | |
d118f9e0 | 707 | |
708 | Int_t binsv2[ndims] = {nMgg, nPt, nCent, ncos2phi}; | |
709 | Double_t xmin2[ndims] = { 0, 0., 0, -1.}; | |
710 | Double_t xmax2[ndims] = { 0.5, 20., 100, 1.}; | |
711 | fHEPV0AM2 = new THnSparseF("fHEPV0AM2", "Flow histogram EPV0A M2", ndims, binsv2, xmin2, xmax2); | |
712 | fHEPV0CM2 = new THnSparseF("fHEPV0CM2", "Flow histogram EPV0C M2", ndims, binsv2, xmin2, xmax2); | |
713 | fHEPTPCM2 = new THnSparseF("fHEPTPCM2", "Flow histogram EPTPC M2", ndims, binsv2, xmin2, xmax2); | |
714 | fHEPV0AM2->GetAxis(0)->SetTitle("m_{#gamma#gamma} "); | |
715 | fHEPV0AM2->GetAxis(1)->SetTitle("p_{T}[GeV]"); | |
716 | fHEPV0AM2->GetAxis(2)->SetTitle("centrality"); | |
717 | fHEPV0AM2->GetAxis(3)->SetTitle("cos(2*#delta #phi)"); | |
718 | fHEPV0CM2->GetAxis(0)->SetTitle("m_{#gamma#gamma} "); | |
719 | fHEPV0CM2->GetAxis(1)->SetTitle("p_{T}[GeV]"); | |
720 | fHEPV0CM2->GetAxis(2)->SetTitle("centrality"); | |
721 | fHEPV0CM2->GetAxis(3)->SetTitle("cos(2*#delta #phi)"); | |
722 | fHEPTPCM2->GetAxis(0)->SetTitle("m_{#gamma#gamma} "); | |
723 | fHEPTPCM2->GetAxis(1)->SetTitle("p_{T}[GeV]"); | |
724 | fHEPTPCM2->GetAxis(2)->SetTitle("centrality"); | |
725 | fHEPTPCM2->GetAxis(3)->SetTitle("cos(2*#delta #phi)"); | |
726 | fOutput->Add(fHEPV0AM2); | |
727 | fOutput->Add(fHEPV0CM2); | |
728 | fOutput->Add(fHEPTPCM2); | |
729 | ||
79ad78fd | 730 | } |
cd0e5e33 | 731 | PostData(1, fOutput); // Post data for ALL output slots >0 here, to get at least an empty histogram |
3c40321c | 732 | } |
733 | ||
734 | //________________________________________________________________________ | |
735 | void AliAnalysisTaskPi0V2::UserExec(Option_t *) | |
736 | { | |
79ad78fd | 737 | // Main loop |
738 | // Called for each event | |
739 | ||
740 | hEvtCount->Fill(1); | |
741 | // Create pointer to reconstructed event | |
cd0e5e33 | 742 | |
79ad78fd | 743 | AliVEvent *event = InputEvent(); |
cd0e5e33 | 744 | if (!event) { |
745 | AliError("Could not retrieve event"); | |
746 | return; | |
747 | } | |
748 | ||
749 | // create pointer to event | |
750 | TString type = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()->GetDataType(); | |
af1defae | 751 | if (type=="ESD") { |
79ad78fd | 752 | fESD = dynamic_cast<AliESDEvent*>(event); |
753 | if (!fESD) { | |
cd0e5e33 | 754 | AliError("Cannot get the ESD event"); |
755 | return; | |
79ad78fd | 756 | } |
af1defae | 757 | } else if (type=="AOD") { |
79ad78fd | 758 | fAOD = dynamic_cast<AliAODEvent*>(event); |
759 | if (!fAOD) { | |
cd0e5e33 | 760 | AliError("Cannot get the AOD event"); |
761 | return; | |
79ad78fd | 762 | } |
cd0e5e33 | 763 | } else { |
764 | AliError("Cannot happen"); | |
765 | return; | |
8071d5b2 | 766 | } |
ea7921ea | 767 | |
4fbbf89d | 768 | hEvtCount->Fill(2); |
cd0e5e33 | 769 | if (!fTrigClass.IsNull()) { |
7bb82cc2 | 770 | TString fired; |
cd0e5e33 | 771 | if (fESD) { |
772 | fired = fESD->GetFiredTriggerClasses(); | |
773 | } else { | |
774 | fired = fAOD->GetFiredTriggerClasses(); | |
79ad78fd | 775 | } |
7bb82cc2 | 776 | if (!fired.Contains("-B-")) |
777 | return; | |
778 | TObjArray *arr = fTrigClass.Tokenize("|"); | |
779 | if (!arr) | |
780 | return; | |
781 | Bool_t match = 0; | |
782 | for (Int_t i=0;i<arr->GetEntriesFast();++i) { | |
783 | TObject *obj = arr->At(i); | |
784 | if (!obj) | |
785 | continue; | |
786 | if (fired.Contains(obj->GetName())) { | |
787 | match = 1; | |
788 | break; | |
789 | } | |
790 | } | |
791 | delete arr; | |
cd0e5e33 | 792 | if (!match) |
793 | return; | |
7bb82cc2 | 794 | } |
cd0e5e33 | 795 | hEvtCount->Fill(3); |
79ad78fd | 796 | |
a20c44e6 | 797 | if (fRunNumber != event->GetRunNumber()) { |
798 | fRunNumber = event->GetRunNumber(); | |
799 | SetFlatteningData(); | |
79ad78fd | 800 | } |
cd0e5e33 | 801 | |
af1defae | 802 | fInterRunNumber = ConvertToInternalRunNumber(fRunNumber); |
803 | ||
cd0e5e33 | 804 | const AliVVertex* fvertex; |
af1defae | 805 | fvertex = event->GetPrimaryVertex(); |
965c985f | 806 | |
cd0e5e33 | 807 | if (TMath::Abs(fvertex->GetZ())>fVtxCut) |
808 | return; | |
809 | Double_t vertex[3] = {fvertex->GetX(), fvertex->GetY(), fvertex->GetZ()}; | |
810 | ||
811 | hEvtCount->Fill(4); | |
ef7e23cf | 812 | |
cd0e5e33 | 813 | fCentrality = event->GetCentrality()->GetCentralityPercentile("CL1"); //spd vertex |
900db6d2 | 814 | hCent->Fill(fCentrality); |
815 | if(isCentFlat){ | |
816 | Bool_t bIsNot = kFALSE; | |
817 | if (fCentrality<=10){ //0-10% | |
818 | TRandom3 *rndm = new TRandom3(0); | |
819 | Double_t Nrndm = rndm->Uniform(0.,1.); | |
820 | if(fCentrality<=1){ | |
3b7c7423 | 821 | if(Nrndm > 0.77308) bIsNot = kTRUE; |
900db6d2 | 822 | } else if(1<fCentrality && fCentrality<=2) { |
3b7c7423 | 823 | if(Nrndm > 0.75863) bIsNot = kTRUE; |
900db6d2 | 824 | } else if (2<fCentrality && fCentrality<=3){ |
3b7c7423 | 825 | if(Nrndm > 0.76365) bIsNot = kTRUE; |
900db6d2 | 826 | } else if (3<fCentrality && fCentrality<=4){ |
3b7c7423 | 827 | if(Nrndm > 0.76763) bIsNot = kTRUE; |
900db6d2 | 828 | } else if (4<fCentrality && fCentrality<=5){ |
3b7c7423 | 829 | if(Nrndm > 0.76251) bIsNot = kTRUE; |
900db6d2 | 830 | } else if (5<fCentrality && fCentrality<=6){ |
3b7c7423 | 831 | if(Nrndm > 0.79069) bIsNot = kTRUE; |
900db6d2 | 832 | } else if (6<fCentrality && fCentrality<=7){ |
3b7c7423 | 833 | if(Nrndm > 0.77669) bIsNot = kTRUE; |
900db6d2 | 834 | } else if (7<fCentrality && fCentrality<=8){ |
3b7c7423 | 835 | if(Nrndm > 0.78537) bIsNot = kTRUE; |
900db6d2 | 836 | } else if (8<fCentrality && fCentrality<=9){ |
3b7c7423 | 837 | if(Nrndm > 0.82727) bIsNot = kTRUE; |
900db6d2 | 838 | } else if (9<fCentrality && fCentrality<=10){ |
839 | if(Nrndm > 1) bIsNot = kTRUE; | |
840 | } | |
841 | delete rndm; rndm = 0; | |
842 | if(bIsNot) | |
843 | return; | |
3b7c7423 | 844 | } |
845 | } | |
846 | if (10<fCentrality && fCentrality<=50){ //10-50% | |
900db6d2 | 847 | TString centfired; |
848 | if (fESD) { | |
849 | centfired = fESD->GetFiredTriggerClasses(); | |
850 | } else { | |
851 | centfired = fAOD->GetFiredTriggerClasses(); | |
852 | } | |
2669a7d3 | 853 | if(!centfired.Contains("CVLN_B2-B-NOPF-ALLNOTRD") && !centfired.Contains("CVLN_R1-B-NOPF-ALLNOTRD") && !centfired.Contains("CSEMI_R1-B-NOPF-ALLNOTRD")) |
900db6d2 | 854 | return; |
855 | } | |
cd0e5e33 | 856 | hEvtCount->Fill(5); |
857 | ||
858 | AliEventplane *ep = event->GetEventplane(); | |
859 | if (ep) { | |
af1defae | 860 | if (ep->GetEventplane("Q") != -1) |
861 | fEPTPC = ep->GetEventplane("Q"); | |
cd0e5e33 | 862 | else |
863 | fEPTPC = -999.; | |
af1defae | 864 | if (ep->GetEventplane("Q") != -1) |
865 | fEPTPCreso = TMath::Cos(2.*(ep->GetQsubRes())); | |
cd0e5e33 | 866 | else |
867 | fEPTPCreso = -1; | |
af1defae | 868 | |
869 | fEPV0 = ep->GetEventplane("V0", event); | |
870 | fEPV0A = ep->GetEventplane("V0A", event); | |
871 | fEPV0C = ep->GetEventplane("V0C", event); | |
cd0e5e33 | 872 | Double_t qx=0, qy=0; |
873 | Double_t qxr=0, qyr=0; | |
af1defae | 874 | fEPV0Ar = ep->CalculateVZEROEventPlane(event, 4, 5, 2, qxr, qyr); |
875 | fEPV0Cr = ep->CalculateVZEROEventPlane(event, 2, 3, 2, qx, qy); | |
cd0e5e33 | 876 | qxr += qx; |
877 | qyr += qy; | |
878 | fEPV0r = TMath::ATan2(qyr,qxr)/2.; | |
af1defae | 879 | fEPV0AR4 = ep->CalculateVZEROEventPlane(event, 4, 2, qx, qy); |
880 | fEPV0AR5 = ep->CalculateVZEROEventPlane(event, 5, 2, qx, qy); | |
881 | fEPV0AR6 = ep->CalculateVZEROEventPlane(event, 6, 2, qx, qy); | |
882 | fEPV0AR7 = ep->CalculateVZEROEventPlane(event, 7, 2, qx, qy); | |
883 | fEPV0CR0 = ep->CalculateVZEROEventPlane(event, 0, 2, qx, qy); | |
884 | fEPV0CR1 = ep->CalculateVZEROEventPlane(event, 1, 2, qx, qy); | |
885 | fEPV0CR2 = ep->CalculateVZEROEventPlane(event, 2, 2, qx, qy); | |
886 | fEPV0CR3 = ep->CalculateVZEROEventPlane(event, 3, 2, qx, qy); | |
887 | } | |
ef7e23cf | 888 | |
cd0e5e33 | 889 | FillEPQA(); //Fill the EP QA |
890 | ||
891 | hEvtCount->Fill(6); | |
892 | ||
3b7c7423 | 893 | while(fEPV0<0.) fEPV0+=TMath::Pi(); while(fEPV0>TMath::Pi()) fEPV0-=TMath::Pi(); |
894 | while(fEPV0r<0.) fEPV0r+=TMath::Pi(); while(fEPV0r>TMath::Pi()) fEPV0r-=TMath::Pi(); | |
895 | while(fEPV0A<0.) fEPV0A+=TMath::Pi(); while(fEPV0A>TMath::Pi()) fEPV0A-=TMath::Pi(); | |
896 | while(fEPV0C<0.) fEPV0C+=TMath::Pi(); while(fEPV0C>TMath::Pi()) fEPV0C-=TMath::Pi(); | |
897 | while(fEPV0Ar<0.) fEPV0Ar+=TMath::Pi(); while(fEPV0Ar>TMath::Pi()) fEPV0Ar-=TMath::Pi(); | |
898 | while(fEPV0Cr<0.) fEPV0Cr+=TMath::Pi(); while(fEPV0Cr>TMath::Pi()) fEPV0Cr-=TMath::Pi(); | |
899 | ||
900 | while(fEPV0AR4<0.) fEPV0AR4+=TMath::Pi(); while(fEPV0AR4>TMath::Pi()) fEPV0AR4-=TMath::Pi(); | |
901 | while(fEPV0AR7<0.) fEPV0AR7+=TMath::Pi(); while(fEPV0AR7>TMath::Pi()) fEPV0AR7-=TMath::Pi(); | |
902 | while(fEPV0CR0<0.) fEPV0CR0+=TMath::Pi(); while(fEPV0CR0>TMath::Pi()) fEPV0CR0-=TMath::Pi(); | |
903 | while(fEPV0CR3<0.) fEPV0CR3+=TMath::Pi(); while(fEPV0CR3>TMath::Pi()) fEPV0CR3-=TMath::Pi(); | |
904 | ||
905 | while(fEPTPC<0.) fEPTPC+=TMath::Pi(); while(fEPTPC>TMath::Pi()) fEPTPC-=TMath::Pi(); | |
af1defae | 906 | if (fEPTPC != -999. ) |
cd0e5e33 | 907 | hEPTPC->Fill(fCentrality, fEPTPC); |
908 | if (fEPTPCreso!=-1) | |
909 | hresoTPC->Fill(fCentrality, fEPTPCreso); | |
3b7c7423 | 910 | if(isFullHist){ |
911 | hEPV0->Fill(fCentrality, fEPV0); | |
912 | hEPV0Ar->Fill(fCentrality, fEPV0Ar); | |
913 | hEPV0Cr->Fill(fCentrality, fEPV0Cr); | |
914 | hEPV0r->Fill(fCentrality, fEPV0r); | |
915 | hEPV0AR4->Fill(fCentrality, fEPV0AR4); | |
916 | hEPV0AR7->Fill(fCentrality, fEPV0AR7); | |
917 | hEPV0CR0->Fill(fCentrality, fEPV0CR0); | |
918 | hEPV0CR3->Fill(fCentrality, fEPV0CR3); | |
919 | } | |
cd0e5e33 | 920 | hEPV0A->Fill(fCentrality, fEPV0A); |
921 | hEPV0C->Fill(fCentrality, fEPV0C); | |
cd0e5e33 | 922 | |
923 | if (isPhosCali) { | |
924 | // PHOS Flattening | |
925 | fEPV0A = ApplyFlatteningV0A(fEPV0A, fCentrality); //V0A after Phos flatten | |
926 | fEPV0C = ApplyFlatteningV0C(fEPV0C, fCentrality); //V0C after Phos flatten | |
af1defae | 927 | if(fEPTPC != -999.) |
928 | fEPTPC = ApplyFlattening(fEPTPC, fCentrality); //TPC after Phos flatten | |
3b7c7423 | 929 | while(fEPV0A <0.) fEPV0A+=TMath::Pi(); while(fEPV0A >TMath::Pi()) fEPV0A-=TMath::Pi(); |
00a93d18 | 930 | while(fEPV0C <0.) fEPV0C+=TMath::Pi(); while(fEPV0C >TMath::Pi()) fEPV0C-=TMath::Pi(); |
931 | while(fEPTPC <0.) fEPTPC+=TMath::Pi(); while(fEPTPC >TMath::Pi()) fEPTPC-=TMath::Pi(); | |
cd0e5e33 | 932 | } |
933 | ||
934 | if (!isPhosCali) { | |
3b7c7423 | 935 | Double_t EPV0ACor = ApplyFlattening(fEPTPC, fCentrality); |
936 | Double_t EPV0CCor = ApplyFlattening(fEPTPC, fCentrality); | |
937 | Double_t EPTPCCor = ApplyFlattening(fEPTPC, fCentrality); | |
00a93d18 | 938 | while(EPV0ACor <0.) EPV0ACor+=TMath::Pi(); while(EPV0ACor >TMath::Pi()) EPV0ACor-=TMath::Pi(); |
939 | while(EPV0CCor <0.) EPV0CCor+=TMath::Pi(); while(EPV0CCor >TMath::Pi()) EPV0CCor-=TMath::Pi(); | |
940 | while(EPTPCCor <0.) EPTPCCor+=TMath::Pi(); while(EPTPCCor >TMath::Pi()) EPTPCCor-=TMath::Pi(); | |
3b7c7423 | 941 | |
af1defae | 942 | if(fEPTPC != -999.) |
3b7c7423 | 943 | hEPTPCCor->Fill(fCentrality, EPTPCCor); |
944 | hEPV0ACor->Fill(fCentrality, EPV0ACor); | |
945 | hEPV0CCor->Fill(fCentrality, EPV0CCor); | |
cd0e5e33 | 946 | } else { |
af1defae | 947 | if(fEPTPC != -999.) |
948 | hEPTPCCor->Fill(fCentrality, fEPTPC); | |
cd0e5e33 | 949 | hEPV0ACor->Fill(fCentrality, fEPV0A); |
950 | hEPV0CCor->Fill(fCentrality, fEPV0C); | |
79ad78fd | 951 | } |
04b116e8 | 952 | |
cd0e5e33 | 953 | hdifV0A_V0CR0->Fill(fCentrality, TMath::Cos(2.*(fEPV0A - fEPV0CR0))); |
954 | hdifV0A_V0CR3->Fill(fCentrality, TMath::Cos(2.*(fEPV0A - fEPV0CR3))); | |
955 | hdifV0ACR0_V0CR3->Fill(fCentrality, TMath::Cos(2*(fEPV0CR0 - fEPV0CR3))); | |
956 | hdifV0C_V0AR4->Fill(fCentrality, TMath::Cos(2*(fEPV0C - fEPV0AR4))); | |
957 | hdifV0C_V0AR7->Fill(fCentrality, TMath::Cos(2*(fEPV0C - fEPV0AR7))); | |
958 | hdifV0AR4_V0AR7->Fill(fCentrality, TMath::Cos(2*(fEPV0AR4 - fEPV0AR7))); | |
3c40321c | 959 | |
3b7c7423 | 960 | if(isFullHist){ |
961 | hdifV0Ar_V0Cr->Fill(fCentrality, TMath::Cos(2.*(fEPV0Ar - fEPV0Cr))); | |
962 | hdifV0A_V0C->Fill(fCentrality, TMath::Cos(2*(fEPV0A - fEPV0C))); | |
963 | if (fEPTPC!=-999.){ | |
964 | hdifV0A_TPC->Fill(fCentrality, TMath::Cos(2*(fEPV0A - fEPTPC))); | |
965 | hdifV0C_TPC->Fill(fCentrality, TMath::Cos(2*(fEPV0C - fEPTPC))); | |
966 | } | |
967 | hdifV0C_V0A->Fill(fCentrality, TMath::Cos(2*(fEPV0C - fEPV0A))); | |
cd0e5e33 | 968 | } |
cd0e5e33 | 969 | // Cluster loop for reconstructed event |
970 | ||
971 | //================ for v2 clusterize analysis============================================== | |
972 | if (!isV1Clus) { | |
13e6ff28 | 973 | if (!fV2ClusName.IsNull() && !fV2Clus) { |
974 | fV2Clus = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fV2ClusName)); | |
975 | if (!fV2Clus) { | |
c0ca428c | 976 | AliError(Form("%s: Could not retrieve v2 cluster name %s!", GetName(), fV2ClusName.Data())); |
13e6ff28 | 977 | return; |
978 | } | |
979 | } | |
980 | Int_t nCluster = fV2Clus->GetEntries(); | |
cd0e5e33 | 981 | for (Int_t i=0; i<nCluster; ++i) { |
13e6ff28 | 982 | AliVCluster *c1 = static_cast<AliVCluster*>(fV2Clus->At(i)); |
cd0e5e33 | 983 | if (!c1) |
984 | continue; | |
3b7c7423 | 985 | if(isFullHist) hClusDxDZA->Fill(c1->GetTrackDz(), c1->GetTrackDx()); |
cd0e5e33 | 986 | if (!c1->IsEMCAL()) |
987 | continue; | |
988 | if (!IsGoodCluster(c1)) | |
989 | continue; | |
3b7c7423 | 990 | if(isFullHist) hClusDxDZB->Fill(c1->GetTrackDz(), c1->GetTrackDx()); |
13e6ff28 | 991 | TLorentzVector p1; |
992 | GetMom(p1, c1, vertex); | |
91d49d74 | 993 | Double_t cluPhi = p1.Phi(); |
994 | Double_t cluPt = p1.Pt(); | |
3b7c7423 | 995 | Double_t difclusV0A = cluPhi-fEPV0A; |
00a93d18 | 996 | if(difclusV0A<0.) difclusV0A+=TMath::Pi(); if(difclusV0A>TMath::Pi()) difclusV0A-=TMath::Pi(); |
3b7c7423 | 997 | Double_t difclusV0C = cluPhi-fEPV0C; |
00a93d18 | 998 | if(difclusV0C<0.) difclusV0C+=TMath::Pi(); if(difclusV0C>TMath::Pi()) difclusV0C-=TMath::Pi(); |
91d49d74 | 999 | hclusDif_EPV0A->Fill(fCentrality, difclusV0A, cluPt); |
1000 | hclusDif_EPV0C->Fill(fCentrality, difclusV0C, cluPt); | |
1001 | hclusv2_EPV0A->Fill(fCentrality, TMath::Cos(2.*difclusV0A), cluPt); | |
1002 | hclusv2_EPV0C->Fill(fCentrality, TMath::Cos(2.*difclusV0C), cluPt); | |
cd0e5e33 | 1003 | for (Int_t j=i+1; j<nCluster; ++j) { |
1004 | AliVCluster *c2 = static_cast<AliVCluster*>(fV2Clus->At(j)); | |
1005 | if (!c2) | |
1006 | continue; | |
1007 | if (!c2->IsEMCAL()) | |
1008 | continue; | |
1009 | if (!IsGoodCluster(c2)) | |
1010 | continue; | |
13e6ff28 | 1011 | TLorentzVector p2; |
1012 | GetMom(p2, c2, vertex); | |
79ad78fd | 1013 | FillPion(p1, p2, fEPV0A, fEPV0C, fEPTPC); |
13e6ff28 | 1014 | } |
1015 | } | |
ceba2d0c | 1016 | } |
cd0e5e33 | 1017 | |
13e6ff28 | 1018 | //================ for v1 clusterize analysis============================================== |
cd0e5e33 | 1019 | if (isV1Clus) { |
ceba2d0c | 1020 | if (!fV2ClusName.IsNull() && !fV1Clus) { |
1021 | fV1Clus = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fV1ClusName)); | |
1022 | if (!fV1Clus) { | |
c0ca428c | 1023 | AliError(Form("%s: Could not retrieve v1 cluster name %s!", GetName(), fV1ClusName.Data())); |
ceba2d0c | 1024 | return; |
1025 | } | |
1026 | } | |
1027 | Int_t nClusterV1 = fV1Clus->GetEntries(); | |
cd0e5e33 | 1028 | for (Int_t i=0; i<nClusterV1; ++i) { |
79ad78fd | 1029 | AliVCluster *c3 = dynamic_cast<AliVCluster*>(fV1Clus->At(i)); |
cd0e5e33 | 1030 | if (!c3) |
1031 | continue; | |
1032 | if (!c3->IsEMCAL()) | |
1033 | continue; | |
9ff4e38b | 1034 | Double_t M02c3 = c3->GetM02(); |
1035 | Double_t Dxc3 = c3->GetTrackDx(); | |
1036 | Double_t Dzc3 = c3->GetTrackDz(); | |
1037 | ||
3b7c7423 | 1038 | if(isFullHist) hClusDxDZA->Fill(Dzc3, Dxc3); |
13e6ff28 | 1039 | Float_t clsPosEt[3] = {0,0,0}; |
1040 | c3->GetPosition(clsPosEt); | |
1041 | TVector3 clsVec(clsPosEt); | |
1042 | Double_t Et = c3->E()*TMath::Sin(clsVec.Theta()); | |
9ff4e38b | 1043 | hM02vsPtA->Fill(Et, M02c3); |
cd0e5e33 | 1044 | if (!IsGoodClusterV1(c3)) |
1045 | continue; | |
9ff4e38b | 1046 | hM02vsPtB->Fill(Et, M02c3); |
3b7c7423 | 1047 | if(isFullHist) hClusDxDZB->Fill(Dzc3, Dxc3); |
8f4922cb | 1048 | TLorentzVector p3; |
1049 | GetMom(p3, c3, vertex); | |
79ad78fd | 1050 | FillCluster(p3, fEPV0A, fEPV0C, fEPTPC, c3); |
8f4922cb | 1051 | } |
1052 | } | |
1053 | ||
cd0e5e33 | 1054 | hEvtCount->Fill(7); |
8f4922cb | 1055 | |
cd0e5e33 | 1056 | if (!fTracksName.IsNull() && !fTracks) { |
1057 | fTracks = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTracksName)); | |
1058 | if (!fTracks) { | |
1059 | AliError(Form("%s: Could not retrieve tracks %s!", GetName(), fTracksName.Data())); | |
1060 | return; | |
1061 | } | |
50ebbe79 | 1062 | } |
e5ba16b4 | 1063 | |
cd0e5e33 | 1064 | Int_t ntracks = fTracks->GetEntries(); |
1065 | for (Int_t i=0; i<ntracks; ++i){ | |
1066 | AliVTrack *track = static_cast<AliVTrack*>(fTracks->At(i)); | |
1067 | if (!track) | |
1068 | continue; | |
1069 | Double_t tPhi = track->Phi(); | |
1070 | Double_t tPt = track->Pt(); | |
1071 | Double_t Eta = track->Eta(); | |
1072 | ||
3b7c7423 | 1073 | Double_t difTrackV0 = tPhi-fEPV0; |
1074 | while(difTrackV0 <0.) difTrackV0+=TMath::Pi(); while(difTrackV0 >TMath::Pi()) difTrackV0-=TMath::Pi(); | |
1075 | Double_t difTrackV0A = tPhi-fEPV0A; | |
1076 | while(difTrackV0A <0.) difTrackV0A+=TMath::Pi(); while(difTrackV0A >TMath::Pi()) difTrackV0A-=TMath::Pi(); | |
1077 | Double_t difTrackV0C = tPhi-fEPV0C; | |
1078 | while(difTrackV0C <0.) difTrackV0C+=TMath::Pi(); while(difTrackV0C >TMath::Pi()) difTrackV0C-=TMath::Pi(); | |
1079 | Double_t difTrackTPC = tPhi-fEPTPC; | |
1080 | while(difTrackTPC <0.) difTrackTPC+=TMath::Pi(); while(difTrackTPC >TMath::Pi()) difTrackTPC-=TMath::Pi(); | |
cd0e5e33 | 1081 | if (tPhi*TMath::RadToDeg()>80. && tPhi*TMath::RadToDeg()<180. && Eta <0.7 && Eta >(-0.7)){ |
cd0e5e33 | 1082 | hdifEMC_EPV0A->Fill(fCentrality, difTrackV0A, tPt); |
1083 | hdifEMC_EPV0C->Fill(fCentrality, difTrackV0C, tPt); | |
91d49d74 | 1084 | hCv2EMC_EPV0A->Fill(fCentrality, TMath::Cos(2.*difTrackV0A), tPt); |
1085 | hCv2EMC_EPV0C->Fill(fCentrality, TMath::Cos(2.*difTrackV0C), tPt); | |
cd0e5e33 | 1086 | } else { |
cd0e5e33 | 1087 | hdifout_EPV0A->Fill(fCentrality, difTrackV0A, tPt); |
1088 | hdifout_EPV0C->Fill(fCentrality, difTrackV0C, tPt); | |
91d49d74 | 1089 | hCv2out_EPV0A->Fill(fCentrality, TMath::Cos(2.*difTrackV0A), tPt); |
1090 | hCv2out_EPV0C->Fill(fCentrality, TMath::Cos(2.*difTrackV0C), tPt); | |
cd0e5e33 | 1091 | } |
cd0e5e33 | 1092 | hdifful_EPV0A->Fill(fCentrality, difTrackV0A, tPt); |
1093 | hdifful_EPV0C->Fill(fCentrality, difTrackV0C, tPt); | |
91d49d74 | 1094 | hCv2ful_EPV0A->Fill(fCentrality, TMath::Cos(2.*difTrackV0A), tPt); |
1095 | hCv2ful_EPV0C->Fill(fCentrality, TMath::Cos(2.*difTrackV0C), tPt); | |
cd0e5e33 | 1096 | } |
1097 | hEvtCount->Fill(8); | |
1098 | ||
1099 | // NEW HISTO should be filled before this point, as PostData puts the | |
1100 | // information for this iteration of the UserExec in the container | |
1101 | PostData(1, fOutput); | |
3c40321c | 1102 | } |
cd0e5e33 | 1103 | |
7bb82cc2 | 1104 | //____________________________________________________________________ |
1105 | Int_t AliAnalysisTaskPi0V2::ConvertToInternalRunNumber(Int_t n) | |
1106 | { | |
cd0e5e33 | 1107 | switch(n) { |
1108 | case 170593 : return 179; | |
1109 | case 170572 : return 178; | |
1110 | case 170556 : return 177; | |
1111 | case 170552 : return 176; | |
1112 | case 170546 : return 175; | |
1113 | case 170390 : return 174; | |
1114 | case 170389 : return 173; | |
1115 | case 170388 : return 172; | |
1116 | case 170387 : return 171; | |
1117 | case 170315 : return 170; | |
1118 | case 170313 : return 169; | |
1119 | case 170312 : return 168; | |
1120 | case 170311 : return 167; | |
1121 | case 170309 : return 166; | |
1122 | case 170308 : return 165; | |
1123 | case 170306 : return 164; | |
1124 | case 170270 : return 163; | |
1125 | case 170269 : return 162; | |
1126 | case 170268 : return 161; | |
1127 | case 170267 : return 160; | |
1128 | case 170264 : return 159; | |
1129 | case 170230 : return 158; | |
1130 | case 170228 : return 157; | |
1131 | case 170208 : return 156; | |
1132 | case 170207 : return 155; | |
1133 | case 170205 : return 154; | |
1134 | case 170204 : return 153; | |
1135 | case 170203 : return 152; | |
1136 | case 170195 : return 151; | |
1137 | case 170193 : return 150; | |
1138 | case 170163 : return 149; | |
1139 | case 170162 : return 148; | |
1140 | case 170159 : return 147; | |
1141 | case 170155 : return 146; | |
1142 | case 170152 : return 145; | |
1143 | case 170091 : return 144; | |
1144 | case 170089 : return 143; | |
1145 | case 170088 : return 142; | |
1146 | case 170085 : return 141; | |
1147 | case 170084 : return 140; | |
1148 | case 170083 : return 139; | |
1149 | case 170081 : return 138; | |
1150 | case 170040 : return 137; | |
1151 | case 170038 : return 136; | |
1152 | case 170036 : return 135; | |
1153 | case 170027 : return 134; | |
1154 | case 169981 : return 133; | |
1155 | case 169975 : return 132; | |
1156 | case 169969 : return 131; | |
1157 | case 169965 : return 130; | |
1158 | case 169961 : return 129; | |
1159 | case 169956 : return 128; | |
1160 | case 169926 : return 127; | |
1161 | case 169924 : return 126; | |
1162 | case 169923 : return 125; | |
1163 | case 169922 : return 124; | |
1164 | case 169919 : return 123; | |
1165 | case 169918 : return 122; | |
1166 | case 169914 : return 121; | |
1167 | case 169859 : return 120; | |
1168 | case 169858 : return 119; | |
1169 | case 169855 : return 118; | |
1170 | case 169846 : return 117; | |
1171 | case 169838 : return 116; | |
1172 | case 169837 : return 115; | |
1173 | case 169835 : return 114; | |
1174 | case 169683 : return 113; | |
1175 | case 169628 : return 112; | |
1176 | case 169591 : return 111; | |
1177 | case 169590 : return 110; | |
1178 | case 169588 : return 109; | |
1179 | case 169587 : return 108; | |
1180 | case 169586 : return 107; | |
1181 | case 169584 : return 106; | |
1182 | case 169557 : return 105; | |
1183 | case 169555 : return 104; | |
1184 | case 169554 : return 103; | |
1185 | case 169553 : return 102; | |
1186 | case 169550 : return 101; | |
1187 | case 169515 : return 100; | |
1188 | case 169512 : return 99; | |
1189 | case 169506 : return 98; | |
1190 | case 169504 : return 97; | |
1191 | case 169498 : return 96; | |
1192 | case 169475 : return 95; | |
1193 | case 169420 : return 94; | |
1194 | case 169419 : return 93; | |
1195 | case 169418 : return 92; | |
1196 | case 169417 : return 91; | |
1197 | case 169415 : return 90; | |
1198 | case 169411 : return 89; | |
1199 | case 169238 : return 88; | |
1200 | case 169236 : return 87; | |
1201 | case 169167 : return 86; | |
1202 | case 169160 : return 85; | |
1203 | case 169156 : return 84; | |
1204 | case 169148 : return 83; | |
1205 | case 169145 : return 82; | |
1206 | case 169144 : return 81; | |
1207 | case 169143 : return 80; | |
1208 | case 169138 : return 79; | |
1209 | case 169099 : return 78; | |
1210 | case 169094 : return 77; | |
1211 | case 169091 : return 76; | |
1212 | case 169045 : return 75; | |
1213 | case 169044 : return 74; | |
1214 | case 169040 : return 73; | |
1215 | case 169035 : return 72; | |
1216 | case 168992 : return 71; | |
1217 | case 168988 : return 70; | |
1218 | case 168984 : return 69; | |
1219 | case 168826 : return 68; | |
1220 | case 168777 : return 67; | |
1221 | case 168514 : return 66; | |
1222 | case 168512 : return 65; | |
1223 | case 168511 : return 64; | |
1224 | case 168467 : return 63; | |
1225 | case 168464 : return 62; | |
1226 | case 168461 : return 61; | |
1227 | case 168460 : return 60; | |
1228 | case 168458 : return 59; | |
1229 | case 168362 : return 58; | |
1230 | case 168361 : return 57; | |
1231 | case 168356 : return 56; | |
1232 | case 168342 : return 55; | |
1233 | case 168341 : return 54; | |
1234 | case 168325 : return 53; | |
1235 | case 168322 : return 52; | |
1236 | case 168318 : return 51; | |
1237 | case 168311 : return 50; | |
1238 | case 168310 : return 49; | |
1239 | case 168213 : return 48; | |
1240 | case 168212 : return 47; | |
1241 | case 168208 : return 46; | |
1242 | case 168207 : return 45; | |
1243 | case 168206 : return 44; | |
1244 | case 168205 : return 43; | |
1245 | case 168204 : return 42; | |
1246 | case 168203 : return 41; | |
1247 | case 168181 : return 40; | |
1248 | case 168177 : return 39; | |
1249 | case 168175 : return 38; | |
1250 | case 168173 : return 37; | |
1251 | case 168172 : return 36; | |
1252 | case 168171 : return 35; | |
1253 | case 168115 : return 34; | |
1254 | case 168108 : return 33; | |
1255 | case 168107 : return 32; | |
1256 | case 168105 : return 31; | |
1257 | case 168104 : return 30; | |
1258 | case 168103 : return 29; | |
1259 | case 168076 : return 28; | |
1260 | case 168069 : return 27; | |
1261 | case 168068 : return 26; | |
1262 | case 168066 : return 25; | |
1263 | case 167988 : return 24; | |
1264 | case 167987 : return 23; | |
1265 | case 167986 : return 22; | |
1266 | case 167985 : return 21; | |
1267 | case 167921 : return 20; | |
1268 | case 167920 : return 19; | |
1269 | case 167915 : return 18; | |
1270 | case 167909 : return 17; | |
1271 | case 167903 : return 16; | |
1272 | case 167902 : return 15; | |
1273 | case 167818 : return 14; | |
1274 | case 167814 : return 13; | |
1275 | case 167813 : return 12; | |
1276 | case 167808 : return 11; | |
1277 | case 167807 : return 10; | |
1278 | case 167806 : return 9; | |
1279 | case 167713 : return 8; | |
1280 | case 167712 : return 7; | |
1281 | case 167711 : return 6; | |
1282 | case 167706 : return 5; | |
1283 | case 167693 : return 4; | |
1284 | case 166532 : return 3; | |
1285 | case 166530 : return 2; | |
1286 | case 166529 : return 1; | |
7bb82cc2 | 1287 | |
1288 | default : return 199; | |
cd0e5e33 | 1289 | } |
7bb82cc2 | 1290 | } |
cd0e5e33 | 1291 | |
7bb82cc2 | 1292 | //_______________________________________________________________________ |
1293 | void AliAnalysisTaskPi0V2::FillEPQA() | |
1294 | { | |
79ad78fd | 1295 | h2DcosV0A->Fill(fInterRunNumber, TMath::Cos(fEPV0A)); |
1296 | h2DsinV0A->Fill(fInterRunNumber, TMath::Sin(fEPV0A)); | |
3b7c7423 | 1297 | if(isFullHist){ |
1298 | h2DcosV0C->Fill(fInterRunNumber, TMath::Cos(fEPV0C)); | |
1299 | h2DsinV0C->Fill(fInterRunNumber, TMath::Sin(fEPV0C)); | |
1300 | if (fEPTPC!=-999.){ | |
1301 | h2DcosTPC->Fill(fInterRunNumber, TMath::Cos(fEPTPC)); | |
1302 | h2DsinTPC->Fill(fInterRunNumber, TMath::Sin(fEPTPC)); | |
1303 | } | |
af1defae | 1304 | } |
79ad78fd | 1305 | } |
cd0e5e33 | 1306 | |
79ad78fd | 1307 | //_________________________________________________________________________________ |
cd0e5e33 | 1308 | void AliAnalysisTaskPi0V2::SetFlatteningData() |
1309 | { | |
79ad78fd | 1310 | //Read objects with flattening parameters |
1311 | AliOADBContainer flatContainer("phosFlat"); | |
1312 | flatContainer.InitFromFile(fEPcalibFileName.Data(),"phosFlat"); | |
1313 | TObjArray *maps = (TObjArray*)flatContainer.GetObject(fRunNumber,"phosFlat"); | |
cd0e5e33 | 1314 | if (!maps) { |
1315 | AliError(Form("Can not read Flattening for run %d. \n From file >%s<\n",fRunNumber,fEPcalibFileName.Data())) ; | |
1316 | } else { | |
79ad78fd | 1317 | AliInfo(Form("Setting PHOS flattening with name %s \n",maps->GetName())) ; |
2fbb6418 | 1318 | AliEPFlattener * h = (AliEPFlattener*)maps->At(0) ; |
79ad78fd | 1319 | if(fTPCFlat) delete fTPCFlat ; |
2fbb6418 | 1320 | fTPCFlat = new AliEPFlattener(); |
79ad78fd | 1321 | fTPCFlat = h ; |
2fbb6418 | 1322 | h = (AliEPFlattener*)maps->At(1); |
79ad78fd | 1323 | if(fV0AFlat) delete fV0AFlat ; |
2fbb6418 | 1324 | fV0AFlat = new AliEPFlattener(); |
79ad78fd | 1325 | fV0AFlat = h ; |
2fbb6418 | 1326 | h = (AliEPFlattener*)maps->At(2); |
79ad78fd | 1327 | if(fV0CFlat) delete fV0CFlat ; |
2fbb6418 | 1328 | fV0CFlat = new AliEPFlattener(); |
cd0e5e33 | 1329 | fV0CFlat = h; |
79ad78fd | 1330 | } |
79ad78fd | 1331 | } |
79ad78fd | 1332 | |
cd0e5e33 | 1333 | //____________________________________________________________________________ |
1334 | Double_t AliAnalysisTaskPi0V2::ApplyFlattening(Double_t phi, Double_t c) | |
1335 | { | |
79ad78fd | 1336 | if(fTPCFlat) |
1337 | return fTPCFlat->MakeFlat(phi,c); | |
cd0e5e33 | 1338 | return phi; |
79ad78fd | 1339 | } |
cd0e5e33 | 1340 | |
79ad78fd | 1341 | //____________________________________________________________________________ |
cd0e5e33 | 1342 | Double_t AliAnalysisTaskPi0V2::ApplyFlatteningV0A(Double_t phi, Double_t c) |
1343 | { | |
79ad78fd | 1344 | if(fV0AFlat) |
1345 | return fV0AFlat->MakeFlat(phi,c); | |
cd0e5e33 | 1346 | return phi; |
79ad78fd | 1347 | } |
cd0e5e33 | 1348 | |
79ad78fd | 1349 | //____________________________________________________________________________ |
1350 | Double_t AliAnalysisTaskPi0V2::ApplyFlatteningV0C(Double_t phi, Double_t c){ | |
1351 | ||
1352 | if(fV0CFlat) | |
1353 | return fV0CFlat->MakeFlat(phi,c); | |
cd0e5e33 | 1354 | return phi; |
7bb82cc2 | 1355 | } |
cd0e5e33 | 1356 | |
3c40321c | 1357 | //________________________________________________________________________ |
1358 | void AliAnalysisTaskPi0V2::Terminate(Option_t *) | |
1359 | { | |
3c40321c | 1360 | } |