1 // $Id: AliAnalysisTaskEMCALPi0V2ShSh.cxx$
3 #include "AliAnalysisTaskEMCALPi0V2ShSh.h"
6 //#include <Riostream.h>
15 #include <THnSparse.h>
18 #include <TVirtualFFT.h>
20 //AliRoot include files
21 #include "AliAnalysisTaskSE.h"
22 #include "AliRunLoader.h"
23 #include "AliAnalysisManager.h"
24 #include "AliAnalysisTask.h"
26 #include "AliEMCALGeometry.h"
27 #include "AliESDEvent.h"
28 #include "AliESDVertex.h"
29 #include "AliESDCaloCells.h"
30 #include "AliESDCaloCluster.h"
31 #include "AliESDEvent.h"
32 #include "AliESDHeader.h"
33 #include "AliESDInputHandler.h"
34 #include "AliESDtrack.h"
35 #include "AliKFParticle.h"
36 #include "AliAODEvent.h"
37 #include "AliVCluster.h"
38 #include "AliCentrality.h"
39 #include "AliEventplane.h"
40 #include "AliOADBContainer.h"
44 ClassImp(AliAnalysisTaskEMCALPi0V2ShSh)
46 //________________________________________________________________________
47 AliAnalysisTaskEMCALPi0V2ShSh::AliAnalysisTaskEMCALPi0V2ShSh() :
56 fGeoName("EMCAL_COMPLETEV1"),
83 fHistEPTPCResolution(0),
96 fHistdifV0C0r_V0C3r(0),
99 fHistdifV0A4r_V0A7r(0),
100 fHistdifV0Ar_V0Cr(0),
109 fHistClusterPhiEta(0),
111 fHistClusterEtM02(0),
112 fHistClusterdphiV0(0),
122 // Default constructor.
123 for(Int_t i = 0; i < 12; i++) fGeomMatrix[i] = 0;
126 //________________________________________________________________________
127 AliAnalysisTaskEMCALPi0V2ShSh::AliAnalysisTaskEMCALPi0V2ShSh(const char *name) :
128 AliAnalysisTaskSE(name),
136 fGeoName("EMCAL_COMPLETEV1"),
142 fEPTPCResolution(0.),
163 fHistEPTPCResolution(0),
174 fHistdifV0A_V0C0r(0),
175 fHistdifV0A_V0C3r(0),
176 fHistdifV0C0r_V0C3r(0),
177 fHistdifV0C_V0A4r(0),
178 fHistdifV0C_V0A7r(0),
179 fHistdifV0A4r_V0A7r(0),
180 fHistdifV0Ar_V0Cr(0),
189 fHistClusterPhiEta(0),
191 fHistClusterEtM02(0),
192 fHistClusterdphiV0(0),
203 // Define input and output slots here
204 // Input slot #0 works with a TChain
205 DefineInput(0, TChain::Class());
206 // Output slot #0 id reserved by the base class for AOD
207 // Output slot #1 writes into a TH1 container
208 DefineOutput(1, TList::Class());
211 //________________________________________________________________________
212 void AliAnalysisTaskEMCALPi0V2ShSh::UserCreateOutputObjects()
214 // Create histograms, called once.
216 fESDClusters = new TObjArray();
217 fAODClusters = new TObjArray();
218 fGeom = AliEMCALGeometry::GetInstance(fGeoName.Data());
219 fOADBContainer = new AliOADBContainer("AliEMCALgeo");
220 fOADBContainer->InitFromFile(Form("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root"),"AliEMCALgeo");
222 fOutputList = new TList();
223 fOutputList->SetOwner();// Container cleans up all histos (avoids leaks in merging)
225 fHistTrackPt = new TH1F("fHistTrackPt","Track Transverse Momentum Distribution of Pb+Pb",100,0.0,30.0);
226 fHistTrackPt->GetYaxis()->SetTitle("Entries"); fHistTrackPt->GetXaxis()->SetTitle("P_{t} [GeV/c]");
227 fOutputList->Add(fHistTrackPt);
229 fHistTrackEta = new TH1F("fHistTrackEta","Track Pseudorapidity Distribution of Pb+Pb",100,-1.0,1.0);
230 fHistTrackEta->GetYaxis()->SetTitle("Entries"); fHistTrackEta->GetXaxis()->SetTitle("#eta");
231 fOutputList->Add(fHistTrackEta);
233 fHistTrackPhi = new TH1F("fHistTrackPhi","Track #phi Distribution of Pb+Pb",100,0.0,6.29);
234 fHistTrackPhi->GetYaxis()->SetTitle("Entries"); fHistTrackPhi->GetXaxis()->SetTitle("#phi [rad]");
235 fOutputList->Add(fHistTrackPhi);
237 fHistTrackPhiEta = new TH2F("fHistTrackPhiEta","Track Pseudorapidity vs #phi of Pb+Pb",100,-1.0,1.0,100,0.0,6.29);
238 fHistTrackPhiEta->GetXaxis()->SetTitle("#eta"); fHistTrackPhiEta->GetYaxis()->SetTitle("#phi [rad]");
239 fOutputList->Add(fHistTrackPhiEta);
241 fHistClusterEta = new TH1F("fHistClusterEta","Cluster Pseudorapidity Distribution of Pb+Pb",100,-1.0,1.0);
242 fHistClusterEta->GetYaxis()->SetTitle("Entries"); fHistClusterEta->GetXaxis()->SetTitle("#eta");
243 fOutputList->Add(fHistClusterEta);
245 fHistClusterPhi = new TH1F("fHistClusterPhi","Cluster #phi Distribution of Pb+Pb",100,0.0,6.29);
246 fHistClusterPhi->GetYaxis()->SetTitle("Entries"); fHistClusterPhi->GetXaxis()->SetTitle("#phi [rad]");
247 fOutputList->Add(fHistClusterPhi);
249 fHistClusterPhiEta = new TH2F("fHistClusterPhiEta","Cluster Pseudorapidity vs #phi of Pb+Pb",100,-1.0,1.0,100,0.0,6.29);
250 fHistClusterPhiEta->GetXaxis()->SetTitle("#eta"); fHistClusterPhiEta->GetYaxis()->SetTitle("#phi [rad]");
251 fOutputList->Add(fHistClusterPhiEta);
253 fHistClusterM02 = new TH1F("fHistClusterM02","Cluster M02 Distribution of Pb+Pb",100,0.0,3.0);
254 fHistClusterM02->GetYaxis()->SetTitle("Entries"); fHistClusterM02->GetXaxis()->SetTitle("M02");
255 fOutputList->Add(fHistClusterM02);
257 fHistClusterE = new TH1F("fHistClusterE","Cluster Energy Distribution of Pb+Pb",100,0.0,20.0);
258 fHistClusterE->GetYaxis()->SetTitle("Entries"); fHistClusterE->GetXaxis()->SetTitle("Energy [GeV]");
259 fOutputList->Add(fHistClusterE);
261 fHistClusterEt = new TH1F("fHistClusterEt","Cluster Transverse Energy Distribution of Pb+Pb",100,0.0,20.0);
262 fHistClusterEt->GetYaxis()->SetTitle("Entries"); fHistClusterEt->GetXaxis()->SetTitle("Transverse Energy [GeV]");
263 fOutputList->Add(fHistClusterEt);
265 fHistClusterEM02 = new TH2F("fHistClusterEM02","Cluster Energy vs M02 of Pb+Pb",100,0.0,20.0,100,0.0,3.0);
266 fHistClusterEM02->GetYaxis()->SetTitle("M02"); fHistClusterEM02->GetXaxis()->SetTitle("Energy [GeV]");
267 fOutputList->Add(fHistClusterEM02);
269 fHistClusterEtM02 = new TH2F("fHistClusterEtM02","Cluster Transverse Energy vs M02 of Pb+Pb",100,0.0,20.0,100,0.0,3.0);
270 fHistClusterEtM02->GetYaxis()->SetTitle("M02"); fHistClusterEtM02->GetXaxis()->SetTitle("Transverse Energy [GeV]");
271 fOutputList->Add(fHistClusterEtM02);
273 fHistClusterN = new TH1F("fHistClusterN","Cluster N Distribution of Pb+Pb",30,0.0,30.0);
274 fHistClusterN->GetYaxis()->SetTitle("Entries"); fHistClusterN->GetXaxis()->SetTitle("N");
275 fOutputList->Add(fHistClusterN);
277 fHistClusterEN = new TH2F("fHistClusterEN","N vs Cluster Energy of Pb+Pb",100,0.0,20.0,30,0.0,30.0);
278 fHistClusterEN->GetYaxis()->SetTitle("N"); fHistClusterEN->GetXaxis()->SetTitle("Energy [GeV]");
279 fOutputList->Add(fHistClusterEN);
281 fHistClusterEtN = new TH2F("fHistClusterEtN","N vs Cluster Transverse Energy of Pb+Pb",100,0.0,20.0,30,0.0,30.0);
282 fHistClusterEtN->GetYaxis()->SetTitle("N"); fHistClusterEtN->GetXaxis()->SetTitle("Transverse Energy [GeV]");
283 fOutputList->Add(fHistClusterEtN);
285 fHistClusterdphiV0 = new TH1D("fHistClusterdphiV0","Cluster dphiV0 Distribution of Pb+Pb",100,0.0,TMath::Pi());
286 fHistClusterdphiV0->GetYaxis()->SetTitle("Entries"); fHistClusterdphiV0->GetXaxis()->SetTitle("dphiV0 [rad]");
287 fOutputList->Add(fHistClusterdphiV0);
289 fHistEPV0 = new TH2F("fHistEPV0","V0 Event Plane Angle vs Centrality of Pb+Pb",100,0,100,100,0.0,TMath::Pi());
290 fHistEPV0->GetYaxis()->SetTitle("V0 Event Plane Angle [rad]"); fHistEPV0->GetXaxis()->SetTitle("V0M Centrality");
291 fOutputList->Add(fHistEPV0);
293 fHistEPV0A = new TH2F("fHistEPV0A","V0A Event Plane Angle vs Centrality of Pb+Pb",100,0,100,100,0.0,TMath::Pi());
294 fHistEPV0A->GetYaxis()->SetTitle("V0A Event Plane Angle [rad]"); fHistEPV0A->GetXaxis()->SetTitle("V0M Centrality");
295 fOutputList->Add(fHistEPV0A);
297 fHistEPV0C = new TH2F("fHistEPV0","V0C Event Plane Angle vs Centrality of Pb+Pb",100,0,100,100,0.0,TMath::Pi());
298 fHistEPV0C->GetYaxis()->SetTitle("V0C Event Plane Angle [rad]"); fHistEPV0C->GetXaxis()->SetTitle("V0M Centrality");
299 fOutputList->Add(fHistEPV0C);
301 fHistEPV0r = new TH2F("fHistEPV0r","V0r Event Plane Angle vs Centrality of Pb+Pb",100,0,100,100,0.0,TMath::Pi());
302 fHistEPV0r->GetYaxis()->SetTitle("V0r Event Plane Angle [rad]"); fHistEPV0r->GetXaxis()->SetTitle("V0M Centrality");
303 fOutputList->Add(fHistEPV0r);
305 fHistEPV0Ar = new TH2F("fHistEPV0","V0Ar Event Plane Angle vs Centrality of Pb+Pb",100,0,100,100,0.0,TMath::Pi());
306 fHistEPV0Ar->GetYaxis()->SetTitle("V0Ar Event Plane Angle [rad]"); fHistEPV0Ar->GetXaxis()->SetTitle("V0M Centrality");
307 fOutputList->Add(fHistEPV0Ar);
309 fHistEPV0Cr = new TH2F("fHistEPV0Cr","V0Cr Event Plane Angle vs Centrality of Pb+Pb",100,0,100,100,0.0,TMath::Pi());
310 fHistEPV0Cr->GetYaxis()->SetTitle("V0Cr Event Plane Angle [rad]"); fHistEPV0Cr->GetXaxis()->SetTitle("V0M Centrality");
311 fOutputList->Add(fHistEPV0Cr);
313 fHistEPV0A4r = new TH2F("fHistEPV0Cr","V0A4r Event Plane Angle vs Centrality of Pb+Pb",100,0,100,100,0.0,TMath::Pi());
314 fHistEPV0A4r->GetYaxis()->SetTitle("V0A4r Event Plane Angle [rad]"); fHistEPV0A4r->GetXaxis()->SetTitle("V0M Centrality");
315 fOutputList->Add(fHistEPV0A4r);
317 fHistEPV0A7r = new TH2F("fHistEPV0Cr","V0A7r Event Plane Angle vs Centrality of Pb+Pb",100,0,100,100,0.0,TMath::Pi());
318 fHistEPV0A7r->GetYaxis()->SetTitle("V0A7r Event Plane Angle [rad]"); fHistEPV0A7r->GetXaxis()->SetTitle("V0M Centrality");
319 fOutputList->Add(fHistEPV0A7r);
321 fHistEPV0C0r = new TH2F("fHistEPV0Cr","V0C0r Event Plane Angle vs Centrality of Pb+Pb",100,0,100,100,0.0,TMath::Pi());
322 fHistEPV0C0r->GetYaxis()->SetTitle("V0C0r Event Plane Angle [rad]"); fHistEPV0C0r->GetXaxis()->SetTitle("V0M Centrality");
323 fOutputList->Add(fHistEPV0C0r);
325 fHistEPV0C3r = new TH2F("fHistEPV0Cr","V0C3r Event Plane Angle vs Centrality of Pb+Pb",100,0,100,100,0.0,TMath::Pi());
326 fHistEPV0C3r->GetYaxis()->SetTitle("V0C3r Event Plane Angle [rad]"); fHistEPV0C3r->GetXaxis()->SetTitle("V0M Centrality");
327 fOutputList->Add(fHistEPV0C3r);
329 fHistdifV0A_V0C0r = new TH2F("fHistdifV0A_V0C0r","(V0A - V0C0r) vs Centrality of Pb+Pb",100,0,100,100,-1.0,1.0);
330 fHistdifV0A_V0C0r->GetYaxis()->SetTitle("Cos[2*(V0A - V0C0r)]"); fHistdifV0A_V0C0r->GetXaxis()->SetTitle("V0M Centrality");
331 fOutputList->Add(fHistdifV0A_V0C0r);
333 fHistdifV0A_V0C3r = new TH2F("fHistdifV0A_V0C3r","(V0A - V0C3r) vs Centrality of Pb+Pb",100,0,100,100,-1.0,1.0);
334 fHistdifV0A_V0C3r->GetYaxis()->SetTitle("Cos[2*(V0A - V0C3r)]"); fHistdifV0A_V0C3r->GetXaxis()->SetTitle("V0M Centrality");
335 fOutputList->Add(fHistdifV0A_V0C3r);
337 fHistdifV0C0r_V0C3r = new TH2F("fHistdifV0C0r_V0C3r","(V0C0r - V0C3r) vs Centrality of Pb+Pb",100,0,100,100,-1.0,1.0);
338 fHistdifV0C0r_V0C3r->GetYaxis()->SetTitle("Cos[2*(V0C0r - V0C3r)]"); fHistdifV0C0r_V0C3r->GetXaxis()->SetTitle("V0M Centrality");
339 fOutputList->Add(fHistdifV0C0r_V0C3r);
341 fHistdifV0C_V0A4r = new TH2F("fHistdifV0C_V0A4r","(V0C - V0A4r) vs Centrality of Pb+Pb",100,0,100,100,-1.0,1.0);
342 fHistdifV0C_V0A4r->GetYaxis()->SetTitle("Cos[2*(V0C - V0C4r)]"); fHistdifV0C_V0A4r->GetXaxis()->SetTitle("V0M Centrality");
343 fOutputList->Add(fHistdifV0C_V0A4r);
345 fHistdifV0C_V0A7r = new TH2F("fHistdifV0C_V0A3r","(V0C - V0A7r) vs Centrality of Pb+Pb",100,0,100,100,-1.0,1.0);
346 fHistdifV0C_V0A7r->GetYaxis()->SetTitle("Cos[2*(V0C - V0A7r)]"); fHistdifV0C_V0A7r->GetXaxis()->SetTitle("V0M Centrality");
347 fOutputList->Add(fHistdifV0C_V0A7r);
349 fHistdifV0A4r_V0A7r = new TH2F("fHistdifV0A4r_V0A7r","(V0A4r - V0A7r) vs Centrality of Pb+Pb",100,0,100,100,-1.0,1.0);
350 fHistdifV0A4r_V0A7r->GetYaxis()->SetTitle("Cos[2*(V0A4r - V0A7r)]"); fHistdifV0A4r_V0A7r->GetXaxis()->SetTitle("V0M Centrality");
351 fOutputList->Add(fHistdifV0A4r_V0A7r);
353 fHistdifV0Ar_V0Cr = new TH2F("fHistdifV0Ar_V0Cr","(V0Ar - V0Cr) vs Centrality of Pb+Pb",100,0,100,100,-1.0,1.0);
354 fHistdifV0Ar_V0Cr->GetYaxis()->SetTitle("Cos[2*(V0Ar - V0Cr)]"); fHistdifV0Ar_V0Cr->GetXaxis()->SetTitle("V0M Centrality");
355 fOutputList->Add(fHistdifV0Ar_V0Cr);
357 fHistAllcentV0 = new TH1F("fHistAllcentV0","V0 Event Plane Angle of Pb+Pb",100,0.0,TMath::Pi());
358 fHistAllcentV0->GetXaxis()->SetTitle("V0 Event Plane Angle [rad]"); fHistAllcentV0->GetYaxis()->SetTitle("");
359 fOutputList->Add(fHistAllcentV0);
361 fHistAllcentV0r = new TH1F("fHistAllcentV0r","V0r Event Plane Angle of Pb+Pb",100,0.0,TMath::Pi());
362 fHistAllcentV0r->GetXaxis()->SetTitle("V0r Event Plane Angle [rad]"); fHistAllcentV0r->GetYaxis()->SetTitle("");
363 fOutputList->Add(fHistAllcentV0r);
365 fHistAllcentV0A = new TH1F("fHistAllcentV0A","V0A Event Plane Angle of Pb+Pb",100,0.0,TMath::Pi());
366 fHistAllcentV0A->GetXaxis()->SetTitle("V0A Event Plane Angle [rad]"); fHistAllcentV0A->GetYaxis()->SetTitle("");
367 fOutputList->Add(fHistAllcentV0A);
369 fHistAllcentV0C = new TH1F("fHistAllcentV0C","V0C Event Plane Angle of Pb+Pb",100,0.0,TMath::Pi());
370 fHistAllcentV0C->GetXaxis()->SetTitle("V0C Event Plane Angle [rad]"); fHistAllcentV0C->GetYaxis()->SetTitle("");
371 fOutputList->Add(fHistAllcentV0C);
373 fHistAllcentTPC = new TH1F("fHistAllcentTPC","TPC Event Plane Angle of Pb+Pb",100,0.0,TMath::Pi());
374 fHistAllcentTPC->GetXaxis()->SetTitle("TPC Event Plane Angle [rad]"); fHistAllcentTPC->GetYaxis()->SetTitle("");
375 fOutputList->Add(fHistAllcentTPC);
377 fHistEPTPC = new TH2F("fHistEPTPC","TPC Event Plane Angle vs Centrality of Pb+Pb",100,0,100,100,0.0,TMath::Pi());
378 fHistEPTPC->GetYaxis()->SetTitle("Event Plane Angle [rad]"); fHistEPTPC->GetXaxis()->SetTitle("V0M Centrality");
379 fOutputList->Add(fHistEPTPC);
381 fHistEPTPCResolution = new TH2F("fHistEPTPCResolution","TPC Resolution vs Centrality of Pb+Pb",100,0,100,100,0.0,1.0);
382 fHistEPTPCResolution->GetYaxis()->SetTitle("TPC Resolution"); fHistEPTPCResolution->GetXaxis()->SetTitle("V0M Centrality");
383 fOutputList->Add(fHistEPTPCResolution);
386 // Et M02 V0Mcent DeltaPhi Cos[2*DeltaPhi]
387 Int_t bins[5] = { 500, 250, 100, 100, 100 }; // binning
388 Double_t min[5] = { 0.0, 0.0, 0, 0.0, -1.0}; // min x
389 Double_t max[5] = { 50.0, 2.5, 100, TMath::Pi(), 1.0}; // max x
391 fClusterPbV0 = new THnSparseF("fClusterPbV0","",5,bins,min,max);
392 fClusterPbV0->GetAxis(0)->SetTitle("Transverse Energy [GeV]"); fClusterPbV0->GetAxis(1)->SetTitle("M02"); fClusterPbV0->GetAxis(2)->SetTitle("V0M Centrality");
393 fClusterPbV0->GetAxis(3)->SetTitle("Delta(#phi) [rad]"); fClusterPbV0->GetAxis(4)->SetTitle("Cos[2*Delta(#phi)]");
394 fOutputList->Add(fClusterPbV0);
396 fClusterPbV0A = new THnSparseF("fClusterPbV0A","",5,bins,min,max);
397 fClusterPbV0A->GetAxis(0)->SetTitle("Transverse Energy [GeV]"); fClusterPbV0A->GetAxis(1)->SetTitle("M02"); fClusterPbV0A->GetAxis(2)->SetTitle("V0M Centrality");
398 fClusterPbV0A->GetAxis(3)->SetTitle("Delta(#phi) [rad]"); fClusterPbV0A->GetAxis(4)->SetTitle("Cos[2*Delta(#phi)]");
399 fOutputList->Add(fClusterPbV0A);
401 fClusterPbV0C = new THnSparseF("fClusterPbV0C","",5,bins,min,max);
402 fClusterPbV0C->GetAxis(0)->SetTitle("Transverse Energy [GeV]"); fClusterPbV0C->GetAxis(1)->SetTitle("M02"); fClusterPbV0C->GetAxis(2)->SetTitle("V0M Centrality");
403 fClusterPbV0C->GetAxis(3)->SetTitle("Delta(#phi) [rad]"); fClusterPbV0C->GetAxis(4)->SetTitle("Cos[2*Delta(#phi)]");
404 fOutputList->Add(fClusterPbV0C);
406 fClusterPbTPC = new THnSparseF("fClusterPbTPC","",5,bins,min,max);
407 fClusterPbTPC->GetAxis(0)->SetTitle("Transverse Energy [GeV]"); fClusterPbTPC->GetAxis(1)->SetTitle("M02"); fClusterPbTPC->GetAxis(2)->SetTitle("V0M Centrality");
408 fClusterPbTPC->GetAxis(3)->SetTitle("Delta(#phi) [rad]"); fClusterPbTPC->GetAxis(4)->SetTitle("Cos[2*Delta(#phi)]");
409 fOutputList->Add(fClusterPbTPC);
412 PostData(1, fOutputList);
415 //________________________________________________________________________
416 void AliAnalysisTaskEMCALPi0V2ShSh::UserExec(Option_t *)
418 // Main loop, called for each event.
424 // Create pointer to reconstructed event
425 AliVEvent *event = InputEvent();
426 if (!event) { Printf("ERROR: Could not retrieve event\n"); return; }
428 fESD = dynamic_cast<AliESDEvent*>(event);
430 fAOD = dynamic_cast<AliAODEvent*>(event);
432 printf("ERROR: Could not retrieve the event\n");
437 //Bool_t isSelected =0;
438 //isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & (AliVEvent::kSemiCentral));
439 //if(!isSelected) { return; }
441 if(event->GetCentrality()){
442 fCentralityV0M = event->GetCentrality()->GetCentralityPercentile("V0M");
444 printf("ERROR: Could not retrieve the Centrality\n");
448 fEventPlane = event->GetEventplane();
452 printf("ERROR: Could not retrieve the Centrality\n");
455 Int_t runnumber = InputEvent()->GetRunNumber() ;
457 TObjArray *matEMCAL=(TObjArray*)fOADBContainer->GetObject(runnumber,"EmcalMatrices");
458 for(Int_t mod=0; mod < (fGeom->GetEMCGeometry())->GetNumberOfSuperModules(); mod++){
459 if(fGeoName=="EMCAL_FIRSTYEARV1" && mod>3)
462 fGeom->SetMisalMatrix(fESD->GetEMCALMatrix(mod), mod);
464 // if(fVEvent->GetEMCALMatrix(mod))
465 fGeomMatrix[mod] = (TGeoHMatrix*) matEMCAL->At(mod);
466 fGeom->SetMisalMatrix(fGeomMatrix[mod] , mod);
470 TList *l = fESD->GetList();
472 for(int nk=0;nk<l->GetEntries();nk++){
473 TObject *obj = (TObject*)l->At(nk);
474 TString oname = obj->GetName();
475 //if(oname.Contains("lus"))
476 printf("Object %d has a clus array named %s +++++++++\n",nk,oname.Data());
479 fESDClusters = dynamic_cast<TClonesArray*>(l->FindObject("CaloClusters"));
480 fESDCells = fESD->GetEMCALCells();
482 printf("ESD cluster mult= %d\n",fESDClusters->GetEntriesFast());
483 if(fESDClusters->GetEntriesFast()<1){
484 printf("ERROR: There are no EMCAL clusters in this event\n");
489 fAODClusters = dynamic_cast<TClonesArray*>(fAOD->GetCaloClusters());
490 fAODCells = fAOD->GetEMCALCells();
492 printf("AOD cluster mult= %d\n",fAODClusters->GetEntriesFast());
493 if(fAODClusters->GetEntriesFast()<1){
494 printf("ERROR: There are no EMCAL clusters in this event\n");
502 PostData(1, fOutputList);
505 //________________________________________________________________________
506 void AliAnalysisTaskEMCALPi0V2ShSh::VZEROEventPlane()
507 {// Calculate the V0 Event Plane
509 if (fEventPlane->GetQVector()) {
510 fEPTPC = TVector2::Phi_0_2pi(fEventPlane->GetQVector()->Phi())/2.0; //if(fEPTPC>TMath::Pi()) {fEPTPC-=TMath::Pi();}
511 } else { fEPTPC = -999.; }
513 if (fEventPlane->GetQsub1()&&fEventPlane->GetQsub2()) {
514 fEPTPCResolution = TMath::Cos(2.0*(fEventPlane->GetQsub1()->Phi()/2.0-fEventPlane->GetQsub2()->Phi()/2.0)); }
515 else { fEPTPCResolution = -1; }
517 fEPV0 = TVector2::Phi_0_2pi(fEventPlane->GetEventplane("V0", fESD)); if(fEPV0>TMath::Pi()) {fEPV0-=TMath::Pi();}
518 fEPV0A = TVector2::Phi_0_2pi(fEventPlane->GetEventplane("V0A", fESD)); if(fEPV0A>TMath::Pi()) {fEPV0A-=TMath::Pi();}
519 fEPV0C = TVector2::Phi_0_2pi(fEventPlane->GetEventplane("V0C", fESD)); if(fEPV0C>TMath::Pi()) {fEPV0C-=TMath::Pi();}
521 Double_t qx=0, qy=0, qxr=0, qyr=0;
522 fEPV0Ar = TVector2::Phi_0_2pi(fEventPlane->CalculateVZEROEventPlane(fESD, 4, 5, 2, qxr, qyr)); if(fEPV0Ar>TMath::Pi()) {fEPV0Ar-=TMath::Pi();}
523 fEPV0Cr = TVector2::Phi_0_2pi(fEventPlane->CalculateVZEROEventPlane(fESD, 2, 3, 2, qx, qy)); if(fEPV0Cr>TMath::Pi()) {fEPV0Cr-=TMath::Pi();}
524 qxr += qx; qyr += qy;
525 fEPV0r = TVector2::Phi_0_2pi(TMath::ATan2(qyr,qxr))/2.0;
527 fEPV0A4r = TVector2::Phi_0_2pi(fEventPlane->CalculateVZEROEventPlane(fESD, 4, 2, qx, qy)); if(fEPV0A4r>TMath::Pi()) {fEPV0A4r-=TMath::Pi();}
528 fEPV0A5r = TVector2::Phi_0_2pi(fEventPlane->CalculateVZEROEventPlane(fESD, 5, 2, qx, qy)); if(fEPV0A5r>TMath::Pi()) {fEPV0A5r-=TMath::Pi();}
529 fEPV0A6r = TVector2::Phi_0_2pi(fEventPlane->CalculateVZEROEventPlane(fESD, 6, 2, qx, qy)); if(fEPV0A6r>TMath::Pi()) {fEPV0A6r-=TMath::Pi();}
530 fEPV0A7r = TVector2::Phi_0_2pi(fEventPlane->CalculateVZEROEventPlane(fESD, 7, 2, qx, qy)); if(fEPV0A7r>TMath::Pi()) {fEPV0A7r-=TMath::Pi();}
531 fEPV0C0r = TVector2::Phi_0_2pi(fEventPlane->CalculateVZEROEventPlane(fESD, 0, 2, qx, qy)); if(fEPV0C0r>TMath::Pi()) {fEPV0C0r-=TMath::Pi();}
532 fEPV0C1r = TVector2::Phi_0_2pi(fEventPlane->CalculateVZEROEventPlane(fESD, 1, 2, qx, qy)); if(fEPV0C1r>TMath::Pi()) {fEPV0C1r-=TMath::Pi();}
533 fEPV0C2r = TVector2::Phi_0_2pi(fEventPlane->CalculateVZEROEventPlane(fESD, 2, 2, qx, qy)); if(fEPV0C2r>TMath::Pi()) {fEPV0C2r-=TMath::Pi();}
534 fEPV0C3r = TVector2::Phi_0_2pi(fEventPlane->CalculateVZEROEventPlane(fESD, 3, 2, qx, qy)); if(fEPV0C3r>TMath::Pi()) {fEPV0C3r-=TMath::Pi();}
536 fHistEPTPC->Fill(fCentralityV0M, fEPTPC);
537 if(fEPTPCResolution!=-1) { fHistEPTPCResolution->Fill(fCentralityV0M, fEPTPCResolution); }
538 fHistEPV0->Fill(fCentralityV0M, fEPV0);
539 fHistEPV0A->Fill(fCentralityV0M, fEPV0A);
540 fHistEPV0C->Fill(fCentralityV0M, fEPV0C);
541 fHistEPV0Ar->Fill(fCentralityV0M, fEPV0Ar);
542 fHistEPV0Cr->Fill(fCentralityV0M, fEPV0Cr);
543 fHistEPV0r->Fill(fCentralityV0M, fEPV0r);
544 fHistEPV0A4r->Fill(fCentralityV0M, fEPV0A4r);
545 fHistEPV0A7r->Fill(fCentralityV0M, fEPV0A7r);
546 fHistEPV0C0r->Fill(fCentralityV0M, fEPV0C0r);
547 fHistEPV0C3r->Fill(fCentralityV0M, fEPV0C3r);
549 fHistAllcentV0->Fill(fEPV0);
550 fHistAllcentV0r->Fill(fEPV0r);
551 fHistAllcentV0A->Fill(fEPV0A);
552 fHistAllcentV0C->Fill(fEPV0C);
553 fHistAllcentTPC->Fill(fEPTPC);
555 fHistdifV0A_V0C0r->Fill(fCentralityV0M, TMath::Cos(2.0*(fEPV0A - fEPV0C0r)));
556 fHistdifV0A_V0C3r->Fill(fCentralityV0M, TMath::Cos(2.0*(fEPV0A - fEPV0C3r)));
557 fHistdifV0C0r_V0C3r->Fill(fCentralityV0M, TMath::Cos(2.0*(fEPV0C0r - fEPV0C3r)));
558 fHistdifV0C_V0A4r->Fill(fCentralityV0M, TMath::Cos(2.0*(fEPV0C - fEPV0A4r)));
559 fHistdifV0C_V0A7r->Fill(fCentralityV0M, TMath::Cos(2.0*(fEPV0C - fEPV0A7r)));
560 fHistdifV0A4r_V0A7r->Fill(fCentralityV0M, TMath::Cos(2.0*(fEPV0A4r - fEPV0A7r)));
561 fHistdifV0Ar_V0Cr->Fill(fCentralityV0M, TMath::Cos(2.0*(fEPV0Ar - fEPV0Cr)));
565 //________________________________________________________________________
566 void AliAnalysisTaskEMCALPi0V2ShSh::FillClusterHists()
567 {// Fill cluster histograms.
571 TObjArray *clusters = fESDClusters;
573 clusters = fAODClusters;
577 const Int_t nclus = clusters->GetEntries();
582 for(Int_t iclus=0;iclus<nclus;iclus++){
583 AliVCluster *clus = (AliVCluster *) clusters->At(iclus); // retrieve cluster from esd
586 if(!clus->IsEMCAL()){
587 printf("ERROR: Current Cluster is not an EMCAL Cluster\n");
590 clus->GetPosition(pos);
591 TVector3 vpos(pos[0],pos[1],pos[2]);
592 Double_t Transverse_Energy = ((clus->E())/ (TMath::CosH(vpos.Eta())));
593 Double_t dphiV0 = TVector2::Phi_0_2pi(vpos.Phi()-fEPV0r); if(dphiV0>TMath::Pi()) {dphiV0-=TMath::Pi();}
594 Double_t dphiV0A = TVector2::Phi_0_2pi(vpos.Phi()-fEPV0A); if(dphiV0A>TMath::Pi()) {dphiV0A-=TMath::Pi();}
595 Double_t dphiV0C = TVector2::Phi_0_2pi(vpos.Phi()-fEPV0C); if(dphiV0C>TMath::Pi()) {dphiV0C-=TMath::Pi();}
596 Double_t dphiTPC = TVector2::Phi_0_2pi(vpos.Phi()-fEPTPC); if(dphiTPC>TMath::Pi()) {dphiTPC-=TMath::Pi();}
599 DataV0[0] = Transverse_Energy;
600 DataV0[1] = clus->GetM02();
601 DataV0[2] = fCentralityV0M;
603 DataV0[4] = TMath::Cos(2.0*(dphiV0));
604 fClusterPbV0->Fill(DataV0);
607 DataV0A[0] = Transverse_Energy;
608 DataV0A[1] = clus->GetM02();
609 DataV0A[2] = fCentralityV0M;
610 DataV0A[3] = dphiV0A;
611 DataV0A[4] = TMath::Cos(2.0*(dphiV0A));
612 fClusterPbV0A->Fill(DataV0A);
615 DataV0C[0] = Transverse_Energy;
616 DataV0C[1] = clus->GetM02();
617 DataV0C[2] = fCentralityV0M;
618 DataV0C[3] = dphiV0C;
619 DataV0C[4] = TMath::Cos(2.0*(dphiV0C));
620 fClusterPbV0C->Fill(DataV0C);
623 DataTPC[0] = Transverse_Energy;
624 DataTPC[1] = clus->GetM02();
625 DataTPC[2] = fCentralityV0M;
626 DataTPC[3] = dphiTPC;
627 DataTPC[4] = TMath::Cos(2.0*(dphiTPC));
628 fClusterPbTPC->Fill(DataTPC);
630 fHistClusterE->Fill(clus->E());
631 fHistClusterEt->Fill(Transverse_Energy);
632 fHistClusterM02->Fill(clus->GetM02());
633 fHistClusterN->Fill(clus->GetNCells());
634 fHistClusterPhi->Fill(vpos.Phi());
635 fHistClusterEta->Fill(vpos.Eta());
636 fHistClusterPhiEta->Fill(vpos.Eta(),vpos.Phi());
637 fHistClusterEN->Fill(clus->E(),clus->GetNCells());
638 fHistClusterEM02->Fill(clus->E(),clus->GetM02());
639 fHistClusterEtN->Fill(Transverse_Energy,clus->GetNCells());
640 fHistClusterEtM02->Fill(Transverse_Energy,clus->GetM02());
641 fHistClusterdphiV0->Fill(dphiV0);
646 //________________________________________________________________________
647 void AliAnalysisTaskEMCALPi0V2ShSh::FillTrackHists()
648 {// Fill track histograms.
651 for (Int_t iTracks = 0; iTracks < fESD->GetNumberOfTracks(); iTracks++) {
652 AliESDtrack* track = (AliESDtrack*)fESD->GetTrack(iTracks);
654 printf("ERROR: Could not retreive esdtrack\n");
657 fHistTrackPt->Fill(track->Pt());
658 fHistTrackEta->Fill(track->Eta());
659 fHistTrackPhi->Fill(track->Phi());
660 fHistTrackPhiEta->Fill(track->Eta(),track->Phi());
665 //________________________________________________________________________
666 void AliAnalysisTaskEMCALPi0V2ShSh::Terminate(Option_t *)
668 // Draw result to screen, or perform fitting, normalizations
669 // Called once at the end of the query
671 /* fOutputList = dynamic_cast<TList*> (GetOutputData(1));
672 if(!fOutputList) { Printf("ERROR: could not retrieve TList fOutputList"); return; }
674 fHistTrackPt = dynamic_cast<TH1F*> (fOutputList->FindObject("fHistTrackPt"));
675 if (!fHistTrackPt) { Printf("ERROR: could not retrieve fHistTrackPt"); return;}
676 fHistTrackEta = dynamic_cast<TH1F*> (fOutputList->FindObject("fHistTrackEta"));
677 if (!fHistTrackEta) { Printf("ERROR: could not retrieve fHistTrackEta"); return;}
678 fHistTrackPhi = dynamic_cast<TH1F*> (fOutputList->FindObject("fHistTrackPhi"));
679 if (!fHistTrackPhi) { Printf("ERROR: could not retrieve fHistTrackPhi"); return;}
680 fHistClusterEta = dynamic_cast<TH1F*> (fOutputList->FindObject("fHistClusterEta"));
681 if (!fHistClusterEta) { Printf("ERROR: could not retrieve fHistClusterEta"); return;}
682 fHistClusterPhi = dynamic_cast<TH1F*> (fOutputList->FindObject("fHistClusterPhi"));
683 if (!fHistClusterPhi) { Printf("ERROR: could not retrieve fHistClusterPhi"); return;}
684 fHistClusterE = dynamic_cast<TH1F*> (fOutputList->FindObject("fHistClusterE"));
685 if (!fHistClusterE) { Printf("ERROR: could not retrieve fHistClusterE"); return;}
686 fHistClusterN = dynamic_cast<TH1F*> (fOutputList->FindObject("fHistClusterN"));
687 if (!fHistClusterN) { Printf("ERROR: could not retrieve fHistClusterN"); return;}
688 fHistClusterM02 = dynamic_cast<TH1F*> (fOutputList->FindObject("fHistClusterM02"));
689 if (!fHistClusterM02) { Printf("ERROR: could not retrieve fHistClusterM02"); return;}
691 // Get the physics selection histograms with the selection statistics
692 //AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
693 //AliESDInputHandler *inputH = dynamic_cast<AliESDInputHandler*>(mgr->GetInputEventHandler());
694 //TH2F *histStat = (TH2F*)inputH->GetStatistics();
697 TCanvas *canvas1 = new TCanvas("canvas1","Track P_{T} & #eta & #phi",10,10,1020,510);
698 canvas1->Divide(3,1);
699 canvas1->cd(1)->SetLogy();
700 fHistTrackPt->DrawCopy("");
702 fHistTrackEta->DrawCopy("");
704 fHistTrackPhi->DrawCopy("");
706 TCanvas *canvas2 = new TCanvas("canvas2","Cluster #eta & #phi",10,10,1020,510);
707 canvas2->Divide(2,1);
709 fHistClusterEta->DrawCopy("");
711 fHistClusterPhi->DrawCopy("");
713 TCanvas *canvas3 = new TCanvas("canvas3","Cluster E & N & M02",10,10,1020,510);
714 canvas3->Divide(3,1);
715 canvas3->cd(1)->SetLogy();
716 fHistClusterE->DrawCopy("");
717 canvas3->cd(2)->SetLogy();
718 fHistClusterN->DrawCopy("");
719 canvas3->cd(3)->SetLogy();
720 fHistClusterM02->DrawCopy("");
722 TCanvas *canvas4 = new TCanvas("canvas4","Track #phi vs #eta & Cluster #phi vs #eta",10,10,1020,510);
723 canvas4->Divide(2,1);
724 canvas4->cd(1)->SetLogz();
725 fHistTrackPhiEta->DrawCopy("colorz");
726 canvas4->cd(2)->SetLogz();
727 fHistClusterPhiEta->DrawCopy("colorz");
729 TCanvas *canvas5 = new TCanvas("canvas5","Cluster E vs N & E vs M02",10,10,1020,510);
730 canvas5->Divide(2,1);
731 canvas5->cd(1)->SetLogz();
732 fHistClusterEN->DrawCopy("colorz");
733 canvas5->cd(2)->SetLogz();
734 fHistClusterEM02->DrawCopy("colorz");
736 TCanvas *canvas6 = new TCanvas("canvas6","Cluster Et & Et vs N & Et vs M02",10,10,1020,510);
737 canvas6->cd(1)->SetLogy();
738 fHistClusterEt->DrawCopy("");
739 canvas6->cd(2)->SetLogz();
740 fHistClusterEtN->DrawCopy("colorz");
741 canvas6->cd(3)->SetLogz();
742 fHistClusterEtM02->DrawCopy("colorz");
744 canvas1->SaveAs("lhc11h_2_Track_PT_Eta_Phi.jpg");
745 canvas2->SaveAs("lhc11h_2_Cluster_Eta_Phi.jpg");
746 canvas3->SaveAs("lhc11h_2_Cluster_E_N_M02.jpg");
747 canvas4->SaveAs("lhc11h_2_Track_PhivsEta_Cluster_PhivsEta.jpg");
748 canvas5->SaveAs("lhc11h_2_Cluster_EvsN_EvsM02.jpg");
749 canvas6->SaveAs("lhc11h_2_Cluster_Et_EtvsN_EtvsM02.jpg"); */