1 /**************************************************************************
2 * Copyright(c) 1998-2014, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
16 // Analysis task for identified PHOS cluster from pi0 and take korrelation betwen hadron-pi0 angel's.
17 // Authors: Daniil Ponomarenko (Daniil.Ponomarenko@cern.ch)
21 #include <Riostream.h>
22 #include "THashList.h"
23 #include "TObjArray.h"
24 #include "TClonesArray.h"
34 #include "AliAnalysisManager.h"
35 #include "AliAnalysisTaskSE.h"
36 #include "AliPHOSCorrelations.h"
37 #include "AliPHOSGeometry.h"
38 #include "AliESDEvent.h"
39 #include "AliESDCaloCluster.h"
40 #include "AliESDVertex.h"
41 #include "AliESDtrackCuts.h"
42 #include "AliESDtrack.h"
43 #include "AliAODTrack.h"
44 #include "AliVTrack.h"
46 #include "AliTriggerAnalysis.h"
47 #include "AliPIDResponse.h"
48 #include "AliPHOSEsdCluster.h"
49 #include "AliCDBManager.h"
50 #include "AliPHOSCalibData.h"
51 #include "AliCentrality.h"
52 #include "AliAnalysisTaskSE.h"
53 #include "AliEventplane.h"
54 #include "AliOADBContainer.h"
55 #include "AliAODEvent.h"
56 #include "AliAODCaloCells.h"
57 #include "AliAODCaloCluster.h"
58 #include "AliCaloPhoton.h"
59 #include "AliAODVertex.h"
64 ClassImp(AliPHOSCorrelations)
66 //_______________________________________________________________________________
67 AliPHOSCorrelations::AliPHOSCorrelations()
70 fOutputContainer(0x0),
71 fMinClusterEnergy(0.3),
82 fCheckHibridGlobal(kOnlyHibridTracks),
84 fPeriod(kUndefinedPeriod),
85 fInternalTriggerSelection(kNoSelection),
87 fManualV0EPCalc(false),
93 fMassMeanP0(-20.9476),
96 fMassSigmaP1(-0.0001),
102 fInternalRunNumber(0),
104 fV0Cpol(0.),fV0Apol(0.),
105 fEPcalibFileName("$ALICE_ROOT/OADB/PHOS/PHOSflat.root"),
108 fCentralityEstimator("V0M"),
114 fCaloPhotonsPHOS(0x0),
116 fCaloPhotonsPHOSLists(0x0),
119 //Deafult constructor, no memory allocations here
122 //_______________________________________________________________________________
123 AliPHOSCorrelations::AliPHOSCorrelations(const char *name, Period period)
124 :AliAnalysisTaskSE(name),
126 fOutputContainer(0x0),
127 fMinClusterEnergy(0.3),
138 fCheckHibridGlobal(kOnlyHibridTracks),
141 fInternalTriggerSelection(kNoSelection),
143 fManualV0EPCalc(false),
149 fMassMeanP0(-20.9476),
152 fMassSigmaP1(-0.0001),
158 fInternalRunNumber(0),
160 fV0Cpol(0.),fV0Apol(0.),
161 fEPcalibFileName("$ALICE_ROOT/OADB/PHOS/PHOSflat.root"),
164 fCentralityEstimator("V0M"),
170 fCaloPhotonsPHOS(0x0),
172 fCaloPhotonsPHOSLists(0x0),
176 // Output slots #0 write into a TH1 container
177 DefineOutput(1,THashList::Class());
179 const Int_t nPtAssoc=10 ;
180 Double_t ptAssocBins[nPtAssoc]={0.,0.5,1.0,1.5,2.0,3.,5.,7.,10.,16} ;
181 fAssocBins.Set(nPtAssoc,ptAssocBins) ;
184 Double_t edges[nbins+1] = {0., 5., 10., 20., 30., 40., 50., 60., 70., 80.};
185 TArrayD centEdges(nbins+1, edges);
186 Int_t nMixed[nbins] = {4,4,6,10,20,30,50,100,100};
187 //Int_t nMixed[nbins] = {100,100,100,100,100,100,100,100,100};
188 TArrayI centNMixed(nbins, nMixed);
189 SetCentralityBinning(centEdges, centNMixed);
191 fVertex[0]=0; fVertex[1]=0; fVertex[2]=0;
197 //_______________________________________________________________________________
198 AliPHOSCorrelations::~AliPHOSCorrelations()
200 if(fCaloPhotonsPHOS){
201 delete fCaloPhotonsPHOS;
202 fCaloPhotonsPHOS=0x0;
210 if(fCaloPhotonsPHOSLists){
211 fCaloPhotonsPHOSLists->SetOwner() ;
212 delete fCaloPhotonsPHOSLists;
213 fCaloPhotonsPHOSLists=0x0;
217 fTracksTPCLists->SetOwner() ;
218 delete fTracksTPCLists;
219 fTracksTPCLists=0x0 ;
223 delete fESDtrackCuts;
227 if(fOutputContainer){
228 delete fOutputContainer;
229 fOutputContainer=0x0;
232 //_______________________________________________________________________________
233 void AliPHOSCorrelations::UserCreateOutputObjects()
237 const Int_t nRuns=200 ;
238 const Int_t ptMult = 200;
239 const Double_t ptMin = 0.;
240 const Double_t ptMax = 20.;
243 if(fOutputContainer != NULL) { delete fOutputContainer; }
244 fOutputContainer = new THashList();
245 fOutputContainer->SetOwner(kTRUE);
248 fOutputContainer->Add(new TH1F("hTriggerPassedEvents","Event selection passed Cuts", 20, 0., 20.) );
250 fOutputContainer->Add(new TH1F("hTotSelEvents","Event selection", kTotalSelected+3, 0., double(kTotalSelected+3))) ;
252 fOutputContainer->Add(new TH2F("hSelEvents","Event selection", kTotalSelected+1, 0., double(kTotalSelected+1), nRuns,0.,float(nRuns))) ;
253 fOutputContainer->Add(new TH2F("hCentrality","Event centrality", 100,0.,100.,nRuns,0.,float(nRuns))) ;
254 fOutputContainer->Add(new TH2F("phiRPflat","RP distribution with TPC flat", 100, 0., 2.*TMath::Pi(),20,0.,100.)) ;
255 fOutputContainer->Add(new TH2F("massWindow","mean & sigma", 100,0.095,0.185,500,0.,0.05));
256 fOutputContainer->Add(new TH2F("hCluEvsClu","ClusterMult vs E",200,0.,10.,100,0.,100.)) ;
259 // Set hists, with track's and cluster's angle distributions.
260 SetHistPtNumTrigger(ptMult, ptMin, ptMax);
262 SetHistPHOSClusterMap();
263 SetHistMass(ptMult, ptMin, ptMax);
264 SetHistPtAssoc(ptMult, ptMin, ptMax);
266 // Setup photon lists
267 Int_t kapacity = fNVtxZBins * GetNumberOfCentralityBins() * fNEMRPBins;
268 fCaloPhotonsPHOSLists = new TObjArray(kapacity);
269 fCaloPhotonsPHOSLists->SetOwner();
271 fTracksTPCLists = new TObjArray(kapacity);
272 fTracksTPCLists->SetOwner();
274 PostData(1, fOutputContainer);
276 //_______________________________________________________________________________
277 void AliPHOSCorrelations::SetHistPtNumTrigger(Int_t ptMult, Double_t ptMin, Double_t ptMax)
279 TString spid[4]={"all","cpv","disp","both"} ;
280 for(Int_t ipid=0; ipid<4; ipid++)
282 fOutputContainer->Add(new TH2F(Form("nTrigger_%s", spid[ipid].Data()), Form("Num of trigger particle %s", spid[ipid].Data()), ptMult+300, ptMin, ptMax, 10000, 0, 1 ) );
283 TH2F *h = static_cast<TH2F*>(fOutputContainer->Last()) ;
284 h->GetXaxis()->SetTitle("Pt [GEV]");
285 h->GetYaxis()->SetTitle("#varepsilon"); // 1/efficiensy
288 //_______________________________________________________________________________
289 void AliPHOSCorrelations::SetHistEtaPhi()
291 // Set hists, with track's and cluster's angle distributions.
293 Float_t pi = TMath::Pi();
296 fOutputContainer->Add(new TH2F("clu_phieta","Cluster's #phi & #eta distribution", 300, double(-1.8), double(-0.6), 300, double(-0.2), double(0.2) ) );
297 TH2F * h = static_cast<TH2F*>(fOutputContainer->Last()) ;
298 h->GetXaxis()->SetTitle("#phi [rad]");
299 h->GetYaxis()->SetTitle("#eta");
302 fOutputContainer->Add(new TH2F("clusingle_phieta","Cluster's #phi & #eta distribution", 300, double(-1.8), double(-0.6), 300, double(-0.2), double(0.2) ) );
303 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
304 h->GetXaxis()->SetTitle("#phi [rad]");
305 h->GetYaxis()->SetTitle("#eta");
308 fOutputContainer->Add(new TH2F("track_phieta","TPC track's #phi & #eta distribution", 200, double(-pi-0.3), double(pi+0.3), 200, double(-0.9), double(0.9) ) );
309 h = static_cast<TH2F*>(fOutputContainer->FindObject("track_phieta")) ;
310 h->GetXaxis()->SetTitle("#phi [rad]");
311 h->GetYaxis()->SetTitle("#eta");
313 //_______________________________________________________________________________
314 void AliPHOSCorrelations::SetHistMass(Int_t ptMult, Double_t ptMin, Double_t ptMax)
316 // Set other histograms.
317 // cout<<"\nSetting output SetHist_CutDistribution...";
319 Double_t massMin = fMassInvMean-fMassInvSigma;
320 Double_t massMax = fMassInvMean+fMassInvSigma;
322 TString spid[4]={"all","cpv","disp","both"} ;
326 for(Int_t ipid=0; ipid<4; ipid++)
328 // Real ++++++++++++++++++++++++++++++
330 fOutputContainer->Add(new TH2F(Form("%s_mpt",spid[ipid].Data() )," real ", 100, massMin, massMax, ptMult, ptMin, ptMax ) );
331 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
332 h->GetXaxis()->SetTitle("Mass [GeV]");
333 h->GetYaxis()->SetTitle("Pt [GEV]");
335 // MIX +++++++++++++++++++++++++
337 fOutputContainer->Add(new TH2F(Form("mix_%s_mpt",spid[ipid].Data() )," mix ", 100, massMin, massMax, ptMult, ptMin, ptMax ) );
338 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
339 h->GetXaxis()->SetTitle("Mass [GeV]");
340 h->GetYaxis()->SetTitle("Pt [GEV]");
342 // Real ++++++++++++++++++++++++++++++
344 fOutputContainer->Add(new TH2F(Form("%s_mpt_left",spid[ipid].Data() )," real ", 100, 0.05, 0.1, ptMult, ptMin, ptMax ) );
345 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
346 h->GetXaxis()->SetTitle("Mass [GeV]");
347 h->GetYaxis()->SetTitle("Pt [GEV]");
349 fOutputContainer->Add(new TH2F(Form("%s_mpt_right",spid[ipid].Data() )," real ", 100, 0.2, 0.4, ptMult, ptMin, ptMax ) );
350 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
351 h->GetXaxis()->SetTitle("Mass [GeV]");
352 h->GetYaxis()->SetTitle("Pt [GEV]");
354 // MIX +++++++++++++++++++++++++
356 fOutputContainer->Add(new TH2F(Form("mix_%s_mpt_left",spid[ipid].Data() )," mix ", 100, 0.05, 0.1, ptMult, ptMin, ptMax ) );
357 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
358 h->GetXaxis()->SetTitle("Mass [GeV]");
359 h->GetYaxis()->SetTitle("Pt [GEV]");
361 fOutputContainer->Add(new TH2F(Form("mix_%s_mpt_right",spid[ipid].Data() )," mix ", 100, 0.2, 0.4, ptMult, ptMin, ptMax ) );
362 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
363 h->GetXaxis()->SetTitle("Mass [GeV]");
364 h->GetYaxis()->SetTitle("Pt [GEV]");
367 // Calibration PHOS Module Pi0peak {REAL}
368 for(Int_t mod=1; mod<4; mod++){
369 fOutputContainer->Add(new TH2F(Form("both%d_mpt",mod),Form("Both cuts (CPV + Disp) mod[%d]",mod), 100, massMin, massMax, ptMult, ptMin, ptMax ) );
370 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
371 h->GetXaxis()->SetTitle("Mass [GeV]");
372 h->GetYaxis()->SetTitle("Pt [GEV]");
374 // Calibration PHOS Module Pi0peak {MIX}
375 fOutputContainer->Add(new TH2F(Form("mix_both%d_mpt",mod),Form(" Both cuts (CPV + Disp) mod[%d]",mod), 100, massMin, massMax, ptMult, ptMin, ptMax ) );
376 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
377 h->GetXaxis()->SetTitle("Mass [GeV]");
378 h->GetYaxis()->SetTitle("Pt [GEV]");
383 for(Int_t ipid=0; ipid<4; ipid++)
385 // Real ++++++++++++++++++++++++++++++
387 fOutputContainer->Add(new TH2F(Form("%s_mpt_eff",spid[ipid].Data() )," real ", 100, massMin, massMax, ptMult, ptMin, ptMax ) );
388 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
389 h->GetXaxis()->SetTitle("Mass [GeV]");
390 h->GetYaxis()->SetTitle("Pt [GEV]");
392 // MIX +++++++++++++++++++++++++
394 fOutputContainer->Add(new TH2F(Form("mix_%s_mpt_eff",spid[ipid].Data() )," mix ", 100, massMin, massMax, ptMult, ptMin, ptMax ) );
395 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
396 h->GetXaxis()->SetTitle("Mass [GeV]");
397 h->GetYaxis()->SetTitle("Pt [GEV]");
399 // Real ++++++++++++++++++++++++++++++
401 fOutputContainer->Add(new TH2F(Form("%s_mpt_left_eff",spid[ipid].Data() )," real ", 100, 0.05, 0.1, ptMult, ptMin, ptMax ) );
402 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
403 h->GetXaxis()->SetTitle("Mass [GeV]");
404 h->GetYaxis()->SetTitle("Pt [GEV]");
406 fOutputContainer->Add(new TH2F(Form("%s_mpt_right_eff",spid[ipid].Data() )," real ", 100, 0.2, 0.4, ptMult, ptMin, ptMax ) );
407 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
408 h->GetXaxis()->SetTitle("Mass [GeV]");
409 h->GetYaxis()->SetTitle("Pt [GEV]");
411 // MIX +++++++++++++++++++++++++
413 fOutputContainer->Add(new TH2F(Form("mix_%s_mpt_left_eff",spid[ipid].Data() )," mix ", 100, 0.05, 0.1, ptMult, ptMin, ptMax ) );
414 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
415 h->GetXaxis()->SetTitle("Mass [GeV]");
416 h->GetYaxis()->SetTitle("Pt [GEV]");
418 fOutputContainer->Add(new TH2F(Form("mix_%s_mpt_right_eff",spid[ipid].Data() )," mix ", 100, 0.2, 0.4, ptMult, ptMin, ptMax ) );
419 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
420 h->GetXaxis()->SetTitle("Mass [GeV]");
421 h->GetYaxis()->SetTitle("Pt [GEV]");
424 // Calibration PHOS Module Pi0peak {REAL}
425 for(Int_t mod=1; mod<4; mod++){
426 fOutputContainer->Add(new TH2F(Form("both%d_mpt_eff",mod),Form("Both cuts (CPV + Disp) mod[%d]",mod), 100, massMin, massMax, ptMult, ptMin, ptMax ) );
427 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
428 h->GetXaxis()->SetTitle("Mass [GeV]");
429 h->GetYaxis()->SetTitle("Pt [GEV]");
431 // Calibration PHOS Module Pi0peak {MIX}
432 fOutputContainer->Add(new TH2F(Form("mix_both%d_mpt_eff",mod),Form(" Both cuts (CPV + Disp) mod[%d]",mod), 100, massMin, massMax, ptMult, ptMin, ptMax ) );
433 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
434 h->GetXaxis()->SetTitle("Mass [GeV]");
435 h->GetYaxis()->SetTitle("Pt [GEV]");
439 // cout<<" OK!"<<endl;
441 //_______________________________________________________________________________
442 void AliPHOSCorrelations::SetHistPtAssoc(Int_t ptMult, Double_t ptMin, Double_t ptMax)
444 Double_t pi = TMath::Pi();
447 Float_t PhiMin = -0.5*pi;
448 Float_t PhiMax = 1.5*pi;
450 Float_t EtaMin = -1.;
453 TString spid[4]={"all","cpv","disp","both"} ;
454 Int_t PhotonsInMod[6] = {1, 2, 3, 12, 13, 23};
456 for (int i = 0; i<fAssocBins.GetSize()-1; i++){
457 for(Int_t ipid=0; ipid<4; ipid++){
458 fOutputContainer->Add(new TH3F(Form("%s_ptphieta_ptAssoc_%3.1f",spid[ipid].Data(),fAssocBins.At(i+1)),
459 Form("%s_ptphieta_ptAssoc_%3.1f",spid[ipid].Data(),fAssocBins.At(i+1)),
460 ptMult, ptMin, ptMax, PhiMult, PhiMin, PhiMax, EtaMult, EtaMin, EtaMax ) );
461 TH3F * h = static_cast<TH3F*>(fOutputContainer->Last()) ;
462 h->GetXaxis()->SetTitle("Pt_{triger} [GEV]");
463 h->GetYaxis()->SetTitle("#phi [rad]");
464 h->GetZaxis()->SetTitle("#eta");
466 fOutputContainer->Add(new TH3F(Form("mix_%s_ptphieta_ptAssoc_%3.1f",spid[ipid].Data(),fAssocBins.At(i+1)),
467 Form("Mixed %s_ptphieta_ptAssoc_%3.1f",spid[ipid].Data(),fAssocBins.At(i+1)),
468 ptMult, ptMin, ptMax, PhiMult, PhiMin, PhiMax, EtaMult, EtaMin, EtaMax ) );
469 h = static_cast<TH3F*>(fOutputContainer->Last()) ;
470 h->GetXaxis()->SetTitle("Pt_{triger} [GEV]");
471 h->GetYaxis()->SetTitle("#phi [rad]");
472 h->GetZaxis()->SetTitle("#eta");
475 for(Int_t m=0; m<6; m++)
477 fOutputContainer->Add(new TH3F(Form("mix_%s_ptphieta_ptAssoc_%3.1f_mod%i",spid[ipid].Data(),fAssocBins.At(i+1), PhotonsInMod[m]),
478 Form("Mixed %s_ptphieta_ptAssoc_%3.1f_mod%i",spid[ipid].Data(),fAssocBins.At(i+1), PhotonsInMod[m]),
479 ptMult, ptMin, ptMax, PhiMult, PhiMin, PhiMax, EtaMult, EtaMin, EtaMax ) );
480 h = static_cast<TH3F*>(fOutputContainer->Last()) ;
481 h->GetXaxis()->SetTitle("Pt_{triger} [GEV]");
482 h->GetYaxis()->SetTitle("#phi [rad]");
483 h->GetZaxis()->SetTitle("#eta");
486 for(Int_t itpc=1; itpc<3; itpc++)
488 fOutputContainer->Add(new TH3F(Form("mix_%s_ptphieta_ptAssoc_%3.1f_tpc%i",spid[ipid].Data(),fAssocBins.At(i+1), itpc),
489 Form("Mixed %s_ptphieta_ptAssoc_%3.1f_tpc%i",spid[ipid].Data(),fAssocBins.At(i+1), itpc),
490 ptMult, ptMin, ptMax, PhiMult, PhiMin, PhiMax, EtaMult, EtaMin, EtaMax ) );
491 h = static_cast<TH3F*>(fOutputContainer->Last()) ;
492 h->GetXaxis()->SetTitle("Pt_{triger} [GEV]");
493 h->GetYaxis()->SetTitle("#phi [rad]");
494 h->GetZaxis()->SetTitle("#eta");
500 void AliPHOSCorrelations::SetHistPHOSClusterMap()
502 for(int i = 0; i<3; i++)
504 // Cluster X/Z/E distribution.
505 fOutputContainer->Add(new TH3F(Form("QA_cluXZE_mod%i", i+1),Form("PHOS Clusters XZE distribution of module %i", i+1), 100, 0, 100, 100, 0, 100, 100, 0, 10 ) );
506 TH3F *h = static_cast<TH3F*>(fOutputContainer->Last()) ;
507 h->GetXaxis()->SetTitle("X");
508 h->GetYaxis()->SetTitle("Z");
509 h->GetZaxis()->SetTitle("E");
512 //_______________________________________________________________________________
513 void AliPHOSCorrelations::UserExec(Option_t *)
515 // Main loop, called for each event analyze ESD/AOD
516 // Step 0: Event Objects
518 fEvent = InputEvent();
521 AliError("Event could not be retrieved");
522 PostData(1, fOutputContainer);
528 fEventESD = dynamic_cast<AliESDEvent*>(fEvent);
529 fEventAOD = dynamic_cast<AliAODEvent*>(fEvent);
533 // Step 1(done once):
534 if( fRunNumber != fEvent->GetRunNumber() )
536 fRunNumber = fEvent->GetRunNumber();
537 fInternalRunNumber = ConvertToInternalRunNumber(fRunNumber);
543 if( RejectTriggerMaskSelection() )
545 PostData(1, fOutputContainer);
551 // fVertex, fVertexVector, fVtxBin
553 if( RejectEventVertex() )
555 PostData(1, fOutputContainer);
560 // Step 3: Centrality
561 // fCentrality, fCentBin
563 if( RejectEventCentrality() )
565 PostData(1, fOutputContainer);
568 FillHistogram("hCentrality",fCentrality,fInternalRunNumber-0.5) ;
571 // Step 4: Reaction Plane
572 // fHaveTPCRP, fRP, fRPV0A, fRPV0C, fRPBin
574 fEMRPBin = GetRPBin();
577 // Step 5: Event Photons (PHOS Clusters) selectionMakeFlat
578 SelectPhotonClusters();
579 if( ! fCaloPhotonsPHOS->GetEntriesFast() )
580 LogSelection(kHasPHOSClusters, fInternalRunNumber);
583 // Step 6: Event Associated particles (TPC Tracks) selection
584 SelectAccosiatedTracks();
586 if( ! fTracksTPC->GetEntriesFast() )
587 LogSelection(kHasTPCTracks, fInternalRunNumber);
588 LogSelection(kTotalSelected, fInternalRunNumber);
591 // Step 7: Consider pi0 (photon/cluster) pairs.
598 //ConsiderTracksMixME();
601 // Step 9: Make TPC's mask
605 // Step 10: Update lists
612 PostData(1, fOutputContainer);
614 //_______________________________________________________________________________
615 void AliPHOSCorrelations::SetESDTrackCuts()
618 // Create ESD track cut
619 fESDtrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts() ;
620 //fESDtrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010();
621 fESDtrackCuts->SetRequireTPCRefit(kTRUE);
624 //_______________________________________________________________________________
625 Int_t AliPHOSCorrelations::ConvertToInternalRunNumber(Int_t run){
626 if(fPeriod== kLHC11h){
629 case 170593 : return 179 ;
630 case 170572 : return 178 ;
631 case 170556 : return 177 ;
632 case 170552 : return 176 ;
633 case 170546 : return 175 ;
634 case 170390 : return 174 ;
635 case 170389 : return 173 ;
636 case 170388 : return 172 ;
637 case 170387 : return 171 ;
638 case 170315 : return 170 ;
639 case 170313 : return 169 ;
640 case 170312 : return 168 ;
641 case 170311 : return 167 ;
642 case 170309 : return 166 ;
643 case 170308 : return 165 ;
644 case 170306 : return 164 ;
645 case 170270 : return 163 ;
646 case 170269 : return 162 ;
647 case 170268 : return 161 ;
648 case 170267 : return 160 ;
649 case 170264 : return 159 ;
650 case 170230 : return 158 ;
651 case 170228 : return 157 ;
652 case 170208 : return 156 ;
653 case 170207 : return 155 ;
654 case 170205 : return 154 ;
655 case 170204 : return 153 ;
656 case 170203 : return 152 ;
657 case 170195 : return 151 ;
658 case 170193 : return 150 ;
659 case 170163 : return 149 ;
660 case 170162 : return 148 ;
661 case 170159 : return 147 ;
662 case 170155 : return 146 ;
663 case 170152 : return 145 ;
664 case 170091 : return 144 ;
665 case 170089 : return 143 ;
666 case 170088 : return 142 ;
667 case 170085 : return 141 ;
668 case 170084 : return 140 ;
669 case 170083 : return 139 ;
670 case 170081 : return 138 ;
671 case 170040 : return 137 ;
672 case 170038 : return 136 ;
673 case 170036 : return 135 ;
674 case 170027 : return 134 ;
675 case 169981 : return 133 ;
676 case 169975 : return 132 ;
677 case 169969 : return 131 ;
678 case 169965 : return 130 ;
679 case 169961 : return 129 ;
680 case 169956 : return 128 ;
681 case 169926 : return 127 ;
682 case 169924 : return 126 ;
683 case 169923 : return 125 ;
684 case 169922 : return 124 ;
685 case 169919 : return 123 ;
686 case 169918 : return 122 ;
687 case 169914 : return 121 ;
688 case 169859 : return 120 ;
689 case 169858 : return 119 ;
690 case 169855 : return 118 ;
691 case 169846 : return 117 ;
692 case 169838 : return 116 ;
693 case 169837 : return 115 ;
694 case 169835 : return 114 ;
695 case 169683 : return 113 ;
696 case 169628 : return 112 ;
697 case 169591 : return 111 ;
698 case 169590 : return 110 ;
699 case 169588 : return 109 ;
700 case 169587 : return 108 ;
701 case 169586 : return 107 ;
702 case 169584 : return 106 ;
703 case 169557 : return 105 ;
704 case 169555 : return 104 ;
705 case 169554 : return 103 ;
706 case 169553 : return 102 ;
707 case 169550 : return 101 ;
708 case 169515 : return 100 ;
709 case 169512 : return 99 ;
710 case 169506 : return 98 ;
711 case 169504 : return 97 ;
712 case 169498 : return 96 ;
713 case 169475 : return 95 ;
714 case 169420 : return 94 ;
715 case 169419 : return 93 ;
716 case 169418 : return 92 ;
717 case 169417 : return 91 ;
718 case 169415 : return 90 ;
719 case 169411 : return 89 ;
720 case 169238 : return 88 ;
721 case 169236 : return 87 ;
722 case 169167 : return 86 ;
723 case 169160 : return 85 ;
724 case 169156 : return 84 ;
725 case 169148 : return 83 ;
726 case 169145 : return 82 ;
727 case 169144 : return 81 ;
728 case 169143 : return 80 ;
729 case 169138 : return 79 ;
730 case 169099 : return 78 ;
731 case 169094 : return 77 ;
732 case 169091 : return 76 ;
733 case 169045 : return 75 ;
734 case 169044 : return 74 ;
735 case 169040 : return 73 ;
736 case 169035 : return 72 ;
737 case 168992 : return 71 ;
738 case 168988 : return 70 ;
739 case 168984 : return 69 ;
740 case 168826 : return 68 ;
741 case 168777 : return 67 ;
742 case 168514 : return 66 ;
743 case 168512 : return 65 ;
744 case 168511 : return 64 ;
745 case 168467 : return 63 ;
746 case 168464 : return 62 ;
747 case 168461 : return 61 ;
748 case 168460 : return 60 ;
749 case 168458 : return 59 ;
750 case 168362 : return 58 ;
751 case 168361 : return 57 ;
752 case 168356 : return 56 ;
753 case 168342 : return 55 ;
754 case 168341 : return 54 ;
755 case 168325 : return 53 ;
756 case 168322 : return 52 ;
757 case 168318 : return 51 ;
758 case 168311 : return 50 ;
759 case 168310 : return 49 ;
760 case 168213 : return 48 ;
761 case 168212 : return 47 ;
762 case 168208 : return 46 ;
763 case 168207 : return 45 ;
764 case 168206 : return 44 ;
765 case 168205 : return 43 ;
766 case 168204 : return 42 ;
767 case 168203 : return 41 ;
768 case 168181 : return 40 ;
769 case 168177 : return 39 ;
770 case 168175 : return 38 ;
771 case 168173 : return 37 ;
772 case 168172 : return 36 ;
773 case 168171 : return 35 ;
774 case 168115 : return 34 ;
775 case 168108 : return 33 ;
776 case 168107 : return 32 ;
777 case 168105 : return 31 ;
778 case 168104 : return 30 ;
779 case 168103 : return 29 ;
780 case 168076 : return 28 ;
781 case 168069 : return 27 ;
782 case 168068 : return 26 ;
783 case 168066 : return 25 ;
784 case 167988 : return 24 ;
785 case 167987 : return 23 ;
786 case 167986 : return 22 ;
787 case 167985 : return 21 ;
788 case 167921 : return 20 ;
789 case 167920 : return 19 ;
790 case 167915 : return 18 ;
791 case 167909 : return 17 ;
792 case 167903 : return 16 ;
793 case 167902 : return 15 ;
794 case 167818 : return 14 ;
795 case 167814 : return 13 ;
796 case 167813 : return 12 ;
797 case 167808 : return 11 ;
798 case 167807 : return 10 ;
799 case 167806 : return 9 ;
800 case 167713 : return 8 ;
801 case 167712 : return 7 ;
802 case 167711 : return 6 ;
803 case 167706 : return 5 ;
804 case 167693 : return 4 ;
805 case 166532 : return 3 ;
806 case 166530 : return 2 ;
807 case 166529 : return 1 ;
809 default : return 199;
812 if(fPeriod== kLHC10h){
814 case 139517 : return 137;
815 case 139514 : return 136;
816 case 139513 : return 135;
817 case 139511 : return 134;
818 case 139510 : return 133;
819 case 139507 : return 132;
820 case 139505 : return 131;
821 case 139504 : return 130;
822 case 139503 : return 129;
823 case 139470 : return 128;
824 case 139467 : return 127;
825 case 139466 : return 126;
826 case 139465 : return 125;
827 case 139440 : return 124;
828 case 139439 : return 123;
829 case 139438 : return 122;
830 case 139437 : return 121;
831 case 139360 : return 120;
832 case 139329 : return 119;
833 case 139328 : return 118;
834 case 139314 : return 117;
835 case 139311 : return 116;
836 case 139310 : return 115;
837 case 139309 : return 114;
838 case 139308 : return 113;
839 case 139173 : return 112;
840 case 139172 : return 111;
841 case 139110 : return 110;
842 case 139107 : return 109;
843 case 139105 : return 108;
844 case 139104 : return 107;
845 case 139042 : return 106;
846 case 139038 : return 105;
847 case 139037 : return 104;
848 case 139036 : return 103;
849 case 139029 : return 102;
850 case 139028 : return 101;
851 case 138983 : return 100;
852 case 138982 : return 99;
853 case 138980 : return 98;
854 case 138979 : return 97;
855 case 138978 : return 96;
856 case 138977 : return 95;
857 case 138976 : return 94;
858 case 138973 : return 93;
859 case 138972 : return 92;
860 case 138965 : return 91;
861 case 138924 : return 90;
862 case 138872 : return 89;
863 case 138871 : return 88;
864 case 138870 : return 87;
865 case 138837 : return 86;
866 case 138830 : return 85;
867 case 138828 : return 84;
868 case 138826 : return 83;
869 case 138796 : return 82;
870 case 138795 : return 81;
871 case 138742 : return 80;
872 case 138732 : return 79;
873 case 138730 : return 78;
874 case 138666 : return 77;
875 case 138662 : return 76;
876 case 138653 : return 75;
877 case 138652 : return 74;
878 case 138638 : return 73;
879 case 138624 : return 72;
880 case 138621 : return 71;
881 case 138583 : return 70;
882 case 138582 : return 69;
883 case 138579 : return 68;
884 case 138578 : return 67;
885 case 138534 : return 66;
886 case 138469 : return 65;
887 case 138442 : return 64;
888 case 138439 : return 63;
889 case 138438 : return 62;
890 case 138396 : return 61;
891 case 138364 : return 60;
892 case 138359 : return 59;
893 case 138275 : return 58;
894 case 138225 : return 57;
895 case 138201 : return 56;
896 case 138200 : return 55;
897 case 138197 : return 54;
898 case 138192 : return 53;
899 case 138190 : return 52;
900 case 138154 : return 51;
901 case 138153 : return 50;
902 case 138151 : return 49;
903 case 138150 : return 48;
904 case 138126 : return 47;
905 case 138125 : return 46;
906 case 137848 : return 45;
907 case 137847 : return 44;
908 case 137844 : return 43;
909 case 137843 : return 42;
910 case 137752 : return 41;
911 case 137751 : return 40;
912 case 137748 : return 39;
913 case 137724 : return 38;
914 case 137722 : return 37;
915 case 137718 : return 36;
916 case 137704 : return 35;
917 case 137693 : return 34;
918 case 137692 : return 33;
919 case 137691 : return 32;
920 case 137689 : return 31;
921 case 137686 : return 30;
922 case 137685 : return 29;
923 case 137639 : return 28;
924 case 137638 : return 27;
925 case 137608 : return 26;
926 case 137595 : return 25;
927 case 137549 : return 24;
928 case 137546 : return 23;
929 case 137544 : return 22;
930 case 137541 : return 21;
931 case 137539 : return 20;
932 case 137531 : return 19;
933 case 137530 : return 18;
934 case 137443 : return 17;
935 case 137441 : return 16;
936 case 137440 : return 15;
937 case 137439 : return 14;
938 case 137434 : return 13;
939 case 137432 : return 12;
940 case 137431 : return 11;
941 case 137430 : return 10;
942 case 137366 : return 9;
943 case 137243 : return 8;
944 case 137236 : return 7;
945 case 137235 : return 6;
946 case 137232 : return 5;
947 case 137231 : return 4;
948 case 137165 : return 3;
949 case 137162 : return 2;
950 case 137161 : return 1;
951 default : return 199;
954 if( kLHC13 == fPeriod )
958 case 195344 : return 1;
959 case 195346 : return 2;
960 case 195351 : return 3;
961 case 195389 : return 4;
962 case 195390 : return 5;
963 case 195391 : return 6;
964 case 195478 : return 7;
965 case 195479 : return 8;
966 case 195480 : return 9;
967 case 195481 : return 10;
968 case 195482 : return 11;
969 case 195483 : return 12;
970 case 195529 : return 13;
971 case 195531 : return 14;
972 case 195532 : return 15;
973 case 195566 : return 16;
974 case 195567 : return 17;
975 case 195568 : return 18;
976 case 195592 : return 19;
977 case 195593 : return 20;
978 case 195596 : return 21;
979 case 195633 : return 22;
980 case 195635 : return 23;
981 case 195644 : return 24;
982 case 195673 : return 25;
983 case 195675 : return 26;
984 case 195676 : return 27;
985 case 195677 : return 28;
986 case 195681 : return 29;
987 case 195682 : return 30;
988 case 195720 : return 31;
989 case 195721 : return 32;
990 case 195722 : return 33;
991 case 195724 : return 34;
992 case 195725 : return 34;
993 case 195726 : return 35;
994 case 195727 : return 36;
995 case 195760 : return 37;
996 case 195761 : return 38;
997 case 195765 : return 39;
998 case 195767 : return 40;
999 case 195783 : return 41;
1000 case 195787 : return 42;
1001 case 195826 : return 43;
1002 case 195827 : return 44;
1003 case 195829 : return 45;
1004 case 195830 : return 46;
1005 case 195831 : return 47;
1006 case 195867 : return 48;
1007 case 195869 : return 49;
1008 case 195871 : return 50;
1009 case 195872 : return 51;
1010 case 195873 : return 52;
1011 case 195935 : return 53;
1012 case 195949 : return 54;
1013 case 195950 : return 55;
1014 case 195954 : return 56;
1015 case 195955 : return 57;
1016 case 195958 : return 58;
1017 case 195989 : return 59;
1018 case 195994 : return 60;
1019 case 195998 : return 61;
1020 case 196000 : return 62;
1021 case 196006 : return 63;
1022 case 196085 : return 64;
1023 case 196089 : return 65;
1024 case 196090 : return 66;
1025 case 196091 : return 67;
1026 case 196099 : return 68;
1027 case 196105 : return 69;
1028 case 196107 : return 70;
1029 case 196185 : return 71;
1030 case 196187 : return 72;
1031 case 196194 : return 73;
1032 case 196197 : return 74;
1033 case 196199 : return 75;
1034 case 196200 : return 76;
1035 case 196201 : return 77;
1036 case 196203 : return 78;
1037 case 196208 : return 79;
1038 case 196214 : return 80;
1039 case 196308 : return 81;
1040 case 196309 : return 82;
1041 case 196310 : return 83;
1042 case 196311 : return 84;
1043 case 196433 : return 85;
1044 case 196474 : return 86;
1045 case 196475 : return 87;
1046 case 196477 : return 88;
1047 case 196528 : return 89;
1048 case 196533 : return 90;
1049 case 196535 : return 91;
1050 case 196563 : return 92;
1051 case 196564 : return 93;
1052 case 196566 : return 94;
1053 case 196568 : return 95;
1054 case 196601 : return 96;
1055 case 196605 : return 97;
1056 case 196608 : return 98;
1057 case 196646 : return 99;
1058 case 196648 : return 100;
1059 case 196701 : return 101;
1060 case 196702 : return 102;
1061 case 196703 : return 103;
1062 case 196706 : return 104;
1063 case 196714 : return 105;
1064 case 196720 : return 106;
1065 case 196721 : return 107;
1066 case 196722 : return 108;
1067 case 196772 : return 109;
1068 case 196773 : return 110;
1069 case 196774 : return 111;
1070 case 196869 : return 112;
1071 case 196870 : return 113;
1072 case 196874 : return 114;
1073 case 196876 : return 115;
1074 case 196965 : return 116;
1075 case 196967 : return 117;
1076 case 196972 : return 118;
1077 case 196973 : return 119;
1078 case 196974 : return 120;
1079 case 197003 : return 121;
1080 case 197011 : return 122;
1081 case 197012 : return 123;
1082 case 197015 : return 124;
1083 case 197027 : return 125;
1084 case 197031 : return 126;
1085 case 197089 : return 127;
1086 case 197090 : return 128;
1087 case 197091 : return 129;
1088 case 197092 : return 130;
1089 case 197094 : return 131;
1090 case 197098 : return 132;
1091 case 197099 : return 133;
1092 case 197138 : return 134;
1093 case 197139 : return 135;
1094 case 197142 : return 136;
1095 case 197143 : return 137;
1096 case 197144 : return 138;
1097 case 197145 : return 139;
1098 case 197146 : return 140;
1099 case 197147 : return 140;
1100 case 197148 : return 141;
1101 case 197149 : return 142;
1102 case 197150 : return 143;
1103 case 197152 : return 144;
1104 case 197153 : return 145;
1105 case 197184 : return 146;
1106 case 197189 : return 147;
1107 case 197247 : return 148;
1108 case 197248 : return 149;
1109 case 197254 : return 150;
1110 case 197255 : return 151;
1111 case 197256 : return 152;
1112 case 197258 : return 153;
1113 case 197260 : return 154;
1114 case 197296 : return 155;
1115 case 197297 : return 156;
1116 case 197298 : return 157;
1117 case 197299 : return 158;
1118 case 197300 : return 159;
1119 case 197302 : return 160;
1120 case 197341 : return 161;
1121 case 197342 : return 162;
1122 case 197348 : return 163;
1123 case 197349 : return 164;
1124 case 197351 : return 165;
1125 case 197386 : return 166;
1126 case 197387 : return 167;
1127 case 197388 : return 168;
1128 default : return 199;
1131 if((fPeriod == kUndefinedPeriod) && (fDebug >= 1) )
1133 AliWarning("Period not defined");
1138 //_______________________________________________________________________________
1139 Bool_t AliPHOSCorrelations::RejectTriggerMaskSelection()
1141 const Bool_t REJECT = true;
1142 const Bool_t ACCEPT = false;
1144 // No need to check trigger mask if no selection is done
1145 if( kNoSelection == fInternalTriggerSelection )
1148 Bool_t reject = REJECT;
1150 Bool_t isMB = (fEvent->GetTriggerMask() & (ULong64_t(1)<<1));
1151 Bool_t isCentral = (fEvent->GetTriggerMask() & (ULong64_t(1)<<4));
1152 Bool_t isSemiCentral = (fEvent->GetTriggerMask() & (ULong64_t(1)<<7));
1155 if( kCentralInclusive == fInternalTriggerSelection
1156 && isCentral ) reject = ACCEPT; // accept event.
1157 else if( kCentralExclusive == fInternalTriggerSelection
1158 && isCentral && !isSemiCentral && !isMB ) reject = ACCEPT; // accept event.
1160 else if( kSemiCentralInclusive == fInternalTriggerSelection
1161 && isSemiCentral ) reject = ACCEPT; // accept event
1162 else if( kSemiCentralExclusive == fInternalTriggerSelection
1163 && isSemiCentral && !isCentral && !isMB ) reject = ACCEPT; // accept event.
1165 else if( kMBInclusive == fInternalTriggerSelection
1166 && isMB ) reject = ACCEPT; // accept event.
1167 else if( kMBExclusive == fInternalTriggerSelection
1168 && isMB && !isCentral && !isSemiCentral ) reject = ACCEPT; // accept event.
1170 if( REJECT == reject )
1173 LogSelection(kInternalTriggerMaskSelection, fInternalRunNumber);
1177 //_______________________________________________________________________________
1178 void AliPHOSCorrelations::SetVertex()
1180 const AliVVertex *primaryVertex = fEvent->GetPrimaryVertex();
1183 fVertex[0] = primaryVertex->GetX();
1184 fVertex[1] = primaryVertex->GetY();
1185 fVertex[2] = primaryVertex->GetZ();
1189 //AliError("Event has 0x0 Primary Vertex, defaulting to origo");
1194 fVertexVector = TVector3(fVertex);
1196 fVtxBin=0 ;// No support for vtx binning implemented.
1198 //_______________________________________________________________________________
1199 Bool_t AliPHOSCorrelations::RejectEventVertex()
1201 if( ! fEvent->GetPrimaryVertex() )
1202 return true; // reject
1203 LogSelection(kHasVertex, fInternalRunNumber);
1205 if ( TMath::Abs(fVertexVector.z()) > fMaxAbsVertexZ )
1206 return true; // reject
1207 LogSelection(kHasAbsVertex, fInternalRunNumber);
1209 return false; // accept event.
1211 //_______________________________________________________________________________
1212 void AliPHOSCorrelations::SetCentrality()
1214 AliCentrality *centrality = fEvent->GetCentrality();
1216 fCentrality=centrality->GetCentralityPercentile(fCentralityEstimator);
1219 AliError("Event has 0x0 centrality");
1223 //cout<<"fCentrality: "<<fCentrality<<endl;
1224 //FillHistogram("hCentrality",fCentrality,fInternalRunNumber-0.5) ;
1225 fCentBin = GetCentralityBin(fCentrality);
1227 //_______________________________________________________________________________
1228 Bool_t AliPHOSCorrelations::RejectEventCentrality()
1230 if (fCentrality<fCentCutoffDown)
1231 return true; //reject
1232 if(fCentrality>fCentCutoffUp)
1235 return false; // accept event.
1237 //_______________________________________________________________________________
1238 void AliPHOSCorrelations::SetCentralityBinning(const TArrayD& edges, const TArrayI& nMixed){
1239 // Define centrality bins by their edges
1240 for(int i=0; i<edges.GetSize()-1; ++i)
1241 if(edges.At(i) > edges.At(i+1)) AliFatal("edges are not sorted");
1242 if( edges.GetSize() != nMixed.GetSize()+1) AliFatal("edges and nMixed don't have appropriate relative sizes");
1245 fCentNMixed = nMixed;
1247 //_______________________________________________________________________________
1248 Int_t AliPHOSCorrelations::GetCentralityBin(Float_t centralityV0M){
1249 int lastBinUpperIndex = fCentEdges.GetSize() -1;
1250 if( centralityV0M > fCentEdges[lastBinUpperIndex] ) {
1252 AliWarning( Form("centrality (%f) larger then upper edge of last centrality bin (%f)!", centralityV0M, fCentEdges[lastBinUpperIndex]) );
1253 return lastBinUpperIndex-1;
1255 if( centralityV0M < fCentEdges[0] ) {
1257 AliWarning( Form("centrality (%f) smaller then lower edge of first bin (%f)!", centralityV0M, fCentEdges[0]) );
1261 fCentBin = TMath::BinarySearch<Double_t> ( GetNumberOfCentralityBins(), fCentEdges.GetArray(), centralityV0M );
1264 //_______________________________________________________________________________
1265 void AliPHOSCorrelations::SetCentralityBorders (double down , double up ){
1266 if (down < 0. || up > 100 || up<=down)
1267 AliError( Form("Warning. Bad value of centrality borders. Setting as default: fCentCutoffDown=%2.f, fCentCutoffUp=%2.f",fCentCutoffDown,fCentCutoffUp) );
1269 fCentCutoffDown = down;
1271 AliInfo( Form("Centrality border was set as fCentCutoffDown=%2.f, fCentCutoffUp=%2.f",fCentCutoffDown,fCentCutoffUp) );
1275 //_______________________________________________________________________________
1276 void AliPHOSCorrelations::EvalReactionPlane()
1278 // assigns: fHaveTPCRP and fRP
1279 // also does a few histogram fills
1281 AliEventplane *eventPlane = fEvent->GetEventplane();
1282 if( ! eventPlane ) { AliError("Event has no event plane"); return; }
1284 Double_t reactionPlaneQ = eventPlane->GetEventplane("Q");
1286 if(reactionPlaneQ>=999 || reactionPlaneQ < 0.)
1288 //reaction plain was not defined
1289 fHaveTPCRP = kFALSE;
1297 fRP = reactionPlaneQ;
1301 FillHistogram("phiRPflat",fRP,fCentrality) ;
1303 //_______________________________________________________________________________
1304 Int_t AliPHOSCorrelations::GetRPBin()
1307 averageRP = fRP ; // If possible, it is better to have EP bin from TPC
1308 // to have similar events for miximng (including jets etc) (fRPV0A+fRPV0C+fRP) /3.;
1310 fEMRPBin = Int_t(fNEMRPBins*(averageRP)/TMath::Pi());
1312 if( fEMRPBin > (Int_t)fNEMRPBins-1 )
1313 fEMRPBin = fNEMRPBins-1 ;
1315 if(fEMRPBin < 0) fEMRPBin=0;
1319 //_______________________________________________________________________________
1320 void AliPHOSCorrelations::SelectPhotonClusters()
1322 //Selects PHOS clusters
1324 // clear (or create) array for holding events photons/clusters
1325 if(fCaloPhotonsPHOS)
1326 fCaloPhotonsPHOS->Clear();
1329 fCaloPhotonsPHOS = new TClonesArray("AliCaloPhoton",200);
1330 fCaloPhotonsPHOS->SetOwner();
1335 for (Int_t i=0; i<fEvent->GetNumberOfCaloClusters(); i++)
1337 AliVCluster *clu = fEvent->GetCaloCluster(i);
1338 if (!clu->IsPHOS() || clu->E()< fMinClusterEnergy) continue; // reject cluster
1340 Float_t position[3];
1341 clu->GetPosition(position);
1342 TVector3 global(position) ;
1344 fPHOSGeo->GlobalPos2RelId(global,relId) ;
1345 Int_t modPHOS = relId[0] ;
1346 Int_t cellXPHOS = relId[2];
1347 Int_t cellZPHOS = relId[3] ;
1349 Double_t distBC=clu->GetDistanceToBadChannel();
1350 if(distBC<fMinBCDistance)
1353 if(clu->GetNCells() < fMinNCells) continue ;
1354 if(clu->GetM02() < fMinM02) continue ;
1357 Double_t tof = clu->GetTOF();
1358 if(TMath::Abs(tof) > fTOFCut ) continue ;
1360 TLorentzVector lorentzMomentum;
1361 Double_t ecore = clu->GetCoreEnergy();
1362 //Double_t ecore = clu->E();
1364 FillHistogram("hCluEvsClu", clu->E(), clu->GetNCells()) ;
1366 Double_t origo[3] = {0,0,0}; // don't rely on event vertex, assume (0,0,0) ?
1367 //clu->GetMomentum(lorentzMomentum, fVertex);
1368 clu->GetMomentum(lorentzMomentum, origo);
1371 if(inPHOS>=fCaloPhotonsPHOS->GetSize()){
1372 fCaloPhotonsPHOS->Expand(inPHOS+50) ;
1375 AliCaloPhoton * ph =new((*fCaloPhotonsPHOS)[inPHOS]) AliCaloPhoton(lorentzMomentum.X(),lorentzMomentum.Py(),lorentzMomentum.Z(),lorentzMomentum.E());
1377 ph->SetCluster(clu);
1379 /*Float_t cellId=clu->GetCellAbsId(0) ;
1380 Int_t mod = (Int_t)TMath:: Ceil(cellId/(56*64) ) ; */
1381 ph->SetModule(modPHOS) ;
1383 lorentzMomentum*=ecore/lorentzMomentum.E() ;
1385 //ph->SetNCells(clu->GetNCells());
1386 ph->SetMomV2(&lorentzMomentum) ;
1387 ph->SetDispBit(clu->GetDispersion()<2.5) ;
1388 ph->SetCPVBit(clu->GetEmcCpvDistance()>2.) ;
1390 FillHistogram(Form("QA_cluXZE_mod%i", modPHOS), cellXPHOS, cellZPHOS, lorentzMomentum.E() ) ;
1393 //_______________________________________________________________________________
1394 void AliPHOSCorrelations::SelectAccosiatedTracks()
1396 // clear (or create) array for holding events tracks
1398 fTracksTPC->Clear();
1401 fTracksTPC = new TClonesArray("TLorentzVector",12000);
1404 for (Int_t i=0; i<fEvent->GetNumberOfTracks(); i++)
1407 AliVParticle *track = fEvent->GetTrack(i);
1409 if(!SelectESDTrack((AliESDtrack*)track)) continue ;
1412 if(!SelectAODTrack((AliAODTrack*)track)) continue ;
1414 Double_t px = track->Px();
1415 Double_t py = track->Py();
1416 Double_t pz = track->Pz() ;
1417 Double_t e = track->E() ;
1419 if(iTracks>=fTracksTPC->GetSize())
1420 fTracksTPC->Expand(iTracks+50) ;
1422 new((*fTracksTPC)[iTracks]) TLorentzVector(px, py, pz,e);
1426 //_______________________________________________________________________________
1427 void AliPHOSCorrelations::ConsiderPi0s()
1429 // Must consider only PHOS events in real distribution.
1432 const Int_t nPHOS=fCaloPhotonsPHOS->GetEntriesFast() ;
1433 for(Int_t i1=0; i1 < nPHOS-1; i1++)
1435 AliCaloPhoton * ph1=(AliCaloPhoton*)fCaloPhotonsPHOS->At(i1) ;
1436 for (Int_t i2=i1+1; i2<nPHOS; i2++)
1438 AliCaloPhoton * ph2=(AliCaloPhoton*)fCaloPhotonsPHOS->At(i2) ;
1439 TLorentzVector p12 = *ph1 + *ph2;
1441 Double_t phiTrigger=p12.Phi() ;
1442 Double_t etaTrigger=p12.Eta() ;
1444 Double_t m=p12.M() ;
1445 Double_t pt=p12.Pt() ;
1446 Double_t eff = 1./GetEfficiency(pt);
1447 int mod1 = ph1->Module() ;
1448 int mod2 = ph2->Module() ;
1450 FillHistogram("clu_phieta",phiTrigger,etaTrigger);
1451 FillHistogram("clusingle_phieta",ph1->Phi(), ph1->Eta());
1452 FillHistogram("clusingle_phieta",ph2->Phi(), ph2->Eta());
1455 FillHistogram("all_mpt",m, pt);
1456 FillHistogram("all_mpt_left",m, pt);
1457 FillHistogram("all_mpt_right",m, pt);
1459 FillHistogram("all_mpt_eff",m, pt, eff);
1460 FillHistogram("all_mpt_left_eff",m, pt, eff);
1461 FillHistogram("all_mpt_right_eff",m, pt, eff);
1463 if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1465 FillHistogram("cpv_mpt",m, pt);
1466 FillHistogram("cpv_mpt_left",m, pt);
1467 FillHistogram("cpv_mpt_right",m, pt);
1469 FillHistogram("cpv_mpt_eff",m, pt, eff);
1470 FillHistogram("cpv_mpt_left_eff",m, pt, eff);
1471 FillHistogram("cpv_mpt_right_eff",m, pt, eff);
1474 if ( ph1->IsDispOK() && ph2->IsDispOK() )
1476 FillHistogram("disp_mpt",m, pt);
1477 FillHistogram("disp_mpt_left",m, pt);
1478 FillHistogram("disp_mpt_right",m, pt);
1480 FillHistogram("disp_mpt_eff",m, pt, eff);
1481 FillHistogram("disp_mpt_left_eff",m, pt, eff);
1482 FillHistogram("disp_mpt_right_eff",m, pt, eff);
1483 if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1485 FillHistogram("both_mpt",m, pt);
1486 FillHistogram("both_mpt_left",m, pt);
1487 FillHistogram("both_mpt_right",m, pt);
1489 FillHistogram("both_mpt_eff",m, pt, eff);
1490 FillHistogram("both_mpt_left_eff",m, pt, eff);
1491 FillHistogram("both_mpt_right_eff",m, pt, eff);
1492 if(mod1 == mod2) // for each module
1494 FillHistogram(Form("both%d_mpt",mod1),m, pt);
1495 FillHistogram(Form("both%d_mpt_eff",mod1),m, pt, eff);
1500 if(!TestMass(m,pt)) continue;
1502 FillHistogram("nTrigger_all", pt, 1./eff);
1503 if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1504 FillHistogram("nTrigger_cpv", pt, 1./eff);
1505 if ( ph1->IsDispOK() && ph2->IsDispOK() )
1507 FillHistogram("nTrigger_disp", pt, 1./eff);
1508 if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1509 FillHistogram("nTrigger_both", pt, 1./eff);
1512 // Take track's angles and compare with cluster's angles.
1513 for(Int_t i3=0; i3<fTracksTPC->GetEntriesFast(); i3++){
1514 TLorentzVector * track = (TLorentzVector*)fTracksTPC->At(i3);
1516 Double_t phiAssoc = track->Phi();
1517 Double_t etaAssoc = track->Eta();
1518 Double_t ptAssoc = track->Pt();
1520 Double_t dPhi = phiTrigger - phiAssoc;
1521 while (dPhi > 1.5*TMath::Pi()) dPhi-=2*TMath::Pi();
1522 while (dPhi < -.5*TMath::Pi()) dPhi+=2*TMath::Pi();
1524 Double_t dEta = etaTrigger - etaAssoc;
1526 Double_t ptAssocBin=GetAssocBin(ptAssoc) ;
1527 FillHistogram(Form("all_ptphieta_ptAssoc_%3.1f",ptAssocBin), pt, dPhi, dEta, eff);
1528 if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1529 FillHistogram(Form("cpv_ptphieta_ptAssoc_%3.1f",ptAssocBin), pt, dPhi, dEta, eff);
1531 if ( ph1->IsDispOK() && ph2->IsDispOK() ){
1532 FillHistogram(Form("disp_ptphieta_ptAssoc_%3.1f",ptAssocBin), pt, dPhi, dEta, eff);
1533 if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1534 FillHistogram(Form("both_ptphieta_ptAssoc_%3.1f",ptAssocBin), pt, dPhi, dEta, eff);
1542 //_______________________________________________________________________________
1543 void AliPHOSCorrelations::ConsiderPi0sMix()
1545 // We must consider only PHOS events in real distribution.
1546 //UInt_t currentOfflineTriggerMask = GetCollisionCandidates();
1550 TList * arrayList = GetCaloPhotonsPHOSList(fVtxBin, fCentBin, fEMRPBin);
1551 for(Int_t evi=0; evi<arrayList->GetEntries();evi++)
1553 TClonesArray * mixPHOS = static_cast<TClonesArray*>(arrayList->At(evi));
1554 for (Int_t i1=0; i1 < fCaloPhotonsPHOS->GetEntriesFast(); i1++)
1556 AliCaloPhoton * ph1 = (AliCaloPhoton*)fCaloPhotonsPHOS->At(i1) ;
1557 for(Int_t i2=0; i2<mixPHOS->GetEntriesFast(); i2++)
1559 AliCaloPhoton * ph2 = (AliCaloPhoton*)mixPHOS->At(i2) ;
1560 TLorentzVector p12 = *ph1 + *ph2;
1561 Double_t m=p12.M() ;
1562 Double_t pt=p12.Pt() ;
1563 Double_t eff = 1./GetEfficiency(pt);
1565 int mod1 = ph1->Module() ;
1566 int mod2 = ph2->Module() ;
1568 FillHistogram("mix_all_mpt", m, pt);
1569 FillHistogram("mix_all_mpt_left",m, pt);
1570 FillHistogram("mix_all_mpt_right",m, pt);
1572 FillHistogram("mix_all_mpt_eff", m, pt, eff);
1573 FillHistogram("mix_all_mpt_left_eff",m, pt, eff);
1574 FillHistogram("mix_all_mpt_right_eff",m, pt, eff);
1576 if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1578 FillHistogram("mix_cpv_mpt",m, pt);
1579 FillHistogram("mix_cpv_mpt_left",m, pt);
1580 FillHistogram("mix_cpv_mpt_right",m, pt);
1582 FillHistogram("mix_cpv_mpt_eff",m, pt, eff);
1583 FillHistogram("mix_cpv_mpt_left_eff",m, pt, eff);
1584 FillHistogram("mix_cpv_mpt_right_eff",m, pt, eff);
1586 if ( ph1->IsDispOK() && ph2->IsDispOK() )
1588 FillHistogram("mix_disp_mpt",m, pt);
1589 FillHistogram("mix_disp_mpt_left",m, pt);
1590 FillHistogram("mix_disp_mpt_right",m, pt);
1592 FillHistogram("mix_disp_mpt_eff",m, pt, eff);
1593 FillHistogram("mix_disp_mpt_left_eff",m, pt, eff);
1594 FillHistogram("mix_disp_mpt_right_eff",m, pt, eff);
1596 if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1598 FillHistogram("mix_both_mpt",m, pt);
1599 FillHistogram("mix_both_mpt_left",m, pt);
1600 FillHistogram("mix_both_mpt_right",m, pt);
1602 FillHistogram("mix_both_mpt_eff",m, pt, eff);
1603 FillHistogram("mix_both_mpt_left_eff",m, pt, eff);
1604 FillHistogram("mix_both_mpt_right_eff",m, pt, eff);
1606 if (mod1 == mod2) // for each module
1608 FillHistogram(Form("mix_both%d_mpt",mod1),m, pt);
1609 FillHistogram(Form("mix_both%d_mpt_eff",mod1),m, pt, eff);
1618 //_______________________________________________________________________________
1619 void AliPHOSCorrelations::ConsiderTracksMix()
1621 TList * arrayList = GetTracksTPCList(fVtxBin, fCentBin, fEMRPBin);
1622 for (Int_t i1=0; i1 < fCaloPhotonsPHOS->GetEntriesFast(); i1++) {
1623 AliCaloPhoton * ph1=(AliCaloPhoton*)fCaloPhotonsPHOS->At(i1) ;
1624 for (Int_t i2=0; i2<fCaloPhotonsPHOS->GetEntriesFast(); i2++){
1625 AliCaloPhoton * ph2=(AliCaloPhoton*)fCaloPhotonsPHOS->At(i2) ;
1626 TLorentzVector p12 = *ph1 + *ph2;
1627 Double_t phiTrigger=p12.Phi() ;
1628 Double_t etaTrigger=p12.Eta() ;
1630 Double_t m=p12.M() ;
1631 Double_t pt=p12.Pt() ;
1632 Double_t eff = 1./GetEfficiency(pt);
1633 Int_t mod1 = ph1->Module();
1634 Int_t mod2 = ph2->Module();
1637 if(!TestMass(m,pt)) continue;
1639 for(Int_t evi=0; evi<arrayList->GetEntries();evi++){
1640 TClonesArray * mixTracks = static_cast<TClonesArray*>(arrayList->At(evi));
1641 for(Int_t i3=0; i3<mixTracks->GetEntriesFast(); i3++){
1642 TLorentzVector * track = (TLorentzVector*)mixTracks->At(i3);
1644 Double_t phiAssoc = track->Phi();
1645 Double_t etaAssoc = track->Eta();
1646 Double_t ptAssoc = track->Pt();
1648 Double_t ptAssocBin=GetAssocBin(ptAssoc) ;
1650 Double_t dPhi = phiTrigger - phiAssoc;
1651 while (dPhi > 1.5*TMath::Pi()) dPhi-=2*TMath::Pi();
1652 while (dPhi < -.5*TMath::Pi()) dPhi+=2*TMath::Pi();
1654 Double_t dEta = etaTrigger - etaAssoc;
1656 FillHistogram(Form("mix_all_ptphieta_ptAssoc_%3.1f",ptAssocBin), pt, dPhi, dEta, eff);
1657 FillHistogram(Form("mix_all_ptphieta_ptAssoc_%3.1f_mod%i",ptAssocBin, GetModCase(mod1, mod2)), pt, dPhi, dEta, eff);
1658 FillHistogram(Form("mix_all_ptphieta_ptAssoc_%3.1f_tpc%i",ptAssocBin, CheckTriggerEta(etaTrigger)), pt, dPhi, dEta, eff);
1660 if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1662 FillHistogram(Form("mix_cpv_ptphieta_ptAssoc_%3.1f",ptAssocBin), pt, dPhi, dEta, eff);
1663 FillHistogram(Form("mix_cpv_ptphieta_ptAssoc_%3.1f_mod%i",ptAssocBin, GetModCase(mod1, mod2)), pt, dPhi, dEta, eff);
1664 FillHistogram(Form("mix_cpv_ptphieta_ptAssoc_%3.1f_tpc%i",ptAssocBin, CheckTriggerEta(etaTrigger)), pt, dPhi, dEta, eff);
1667 if ( ph1->IsDispOK() && ph2->IsDispOK() ) {
1668 FillHistogram(Form("mix_disp_ptphieta_ptAssoc_%3.1f",ptAssocBin), pt, dPhi, dEta, eff);
1669 FillHistogram(Form("mix_disp_ptphieta_ptAssoc_%3.1f_mod%i",ptAssocBin, GetModCase(mod1, mod2)), pt, dPhi, dEta, eff);
1670 FillHistogram(Form("mix_disp_ptphieta_ptAssoc_%3.1f_tpc%i",ptAssocBin, CheckTriggerEta(etaTrigger)), pt, dPhi, dEta, eff);
1671 if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1673 FillHistogram(Form("mix_both_ptphieta_ptAssoc_%3.1f",ptAssocBin), pt, dPhi, dEta, eff);
1674 FillHistogram(Form("mix_both_ptphieta_ptAssoc_%3.1f_mod%i",ptAssocBin, GetModCase(mod1, mod2)), pt, dPhi, dEta, eff);
1675 FillHistogram(Form("mix_both_ptphieta_ptAssoc_%3.1f_tpc%i",ptAssocBin, CheckTriggerEta(etaTrigger)), pt, dPhi, dEta, eff);
1684 //_______________________________________________________________________________
1685 void AliPHOSCorrelations::ConsiderPi0sME()
1687 //Seek Most Energetic (ME) Pi0 and work whit it.
1688 // Must consider only PHOS events in real distribution.
1691 const Int_t nPHOS=fCaloPhotonsPHOS->GetEntriesFast() ;
1692 for(Int_t i1=0; i1 < nPHOS-1; i1++)
1694 AliCaloPhoton * ph1=(AliCaloPhoton*)fCaloPhotonsPHOS->At(i1) ;
1695 for (Int_t i2=i1+1; i2<nPHOS; i2++)
1697 AliCaloPhoton * ph2=(AliCaloPhoton*)fCaloPhotonsPHOS->At(i2) ;
1698 TLorentzVector p12 = *ph1 + *ph2;
1700 Double_t phiTrigger=p12.Phi() ;
1701 Double_t etaTrigger=p12.Eta() ;
1703 Double_t m=p12.M() ;
1704 Double_t pt=p12.Pt() ;
1705 Double_t eff = 1./GetEfficiency(pt);
1706 int mod1 = ph1->Module() ;
1707 int mod2 = ph2->Module() ;
1709 FillHistogram("clu_phieta",phiTrigger,etaTrigger);
1710 FillHistogram("clusingle_phieta",ph1->Phi(), ph1->Eta());
1711 FillHistogram("clusingle_phieta",ph2->Phi(), ph2->Eta());
1714 FillHistogram("all_mpt",m, pt);
1715 FillHistogram("all_mpt_left",m, pt);
1716 FillHistogram("all_mpt_right",m, pt);
1718 FillHistogram("all_mpt_eff",m, pt, eff);
1719 FillHistogram("all_mpt_left_eff",m, pt, eff);
1720 FillHistogram("all_mpt_right_eff",m, pt, eff);
1722 if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1724 FillHistogram("cpv_mpt",m, pt);
1725 FillHistogram("cpv_mpt_left",m, pt);
1726 FillHistogram("cpv_mpt_right",m, pt);
1728 FillHistogram("cpv_mpt_eff",m, pt, eff);
1729 FillHistogram("cpv_mpt_left_eff",m, pt, eff);
1730 FillHistogram("cpv_mpt_right_eff",m, pt, eff);
1733 if ( ph1->IsDispOK() && ph2->IsDispOK() )
1735 FillHistogram("disp_mpt",m, pt);
1736 FillHistogram("disp_mpt_left",m, pt);
1737 FillHistogram("disp_mpt_right",m, pt);
1739 FillHistogram("disp_mpt_eff",m, pt, eff);
1740 FillHistogram("disp_mpt_left_eff",m, pt, eff);
1741 FillHistogram("disp_mpt_right_eff",m, pt, eff);
1742 if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1744 FillHistogram("both_mpt",m, pt);
1745 FillHistogram("both_mpt_left",m, pt);
1746 FillHistogram("both_mpt_right",m, pt);
1748 FillHistogram("both_mpt_eff",m, pt, eff);
1749 FillHistogram("both_mpt_left_eff",m, pt, eff);
1750 FillHistogram("both_mpt_right_eff",m, pt, eff);
1751 if(mod1 == mod2) // for each module
1753 FillHistogram(Form("both%d_mpt",mod1),m, pt);
1754 FillHistogram(Form("both%d_mpt_eff",mod1),m, pt, eff);
1759 if(!TestMass(m,pt)) continue;
1761 Int_t modCase = GetModCase(mod1, mod2);
1763 TestPi0ME(kPidAll, p12, modCase);
1764 if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1765 TestPi0ME(kPidCPV, p12, modCase);
1766 if ( ph1->IsDispOK() && ph2->IsDispOK() )
1768 TestPi0ME(kPidDisp, p12, modCase);
1769 if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1770 TestPi0ME(kPidBoth, p12, modCase);
1775 TString spid[4]={"all","cpv","disp","both"} ;
1776 for (int ipid = 0; ipid < 4; ipid++)
1778 if (fMEExists[ipid])
1779 FillHistogram(Form("nTrigger_%s", spid[ipid].Data()), fMEPt[ipid], GetEfficiency(fMEPt[ipid]));
1782 // Take track's angles and compare with cluster's angles.
1783 for(Int_t i3=0; i3<fTracksTPC->GetEntriesFast(); i3++){
1784 TLorentzVector * track = (TLorentzVector*)fTracksTPC->At(i3);
1786 Double_t phiAssoc = track->Phi();
1787 Double_t etaAssoc = track->Eta();
1788 Double_t ptAssoc = track->Pt();
1790 Double_t ptAssocBin=GetAssocBin(ptAssoc) ;
1791 Double_t dPhi(0.), dEta(0.);
1793 for (int ipid = 0; ipid < 4; ipid++)
1795 if (fMEExists[ipid])
1797 dPhi = fMEPhi[ipid] - phiAssoc;
1798 while (dPhi > 1.5*TMath::Pi()) dPhi-=2*TMath::Pi();
1799 while (dPhi < -.5*TMath::Pi()) dPhi+=2*TMath::Pi();
1800 dEta = fMEEta[ipid] - etaAssoc;
1801 FillHistogram(Form("%s_ptphieta_ptAssoc_%3.1f", spid[ipid].Data(), ptAssocBin), fMEPt[ipid], dPhi, dEta, 1./GetEfficiency(fMEPt[ipid]) );
1807 //_______________________________________________________________________________
1808 void AliPHOSCorrelations::ConsiderTracksMixME()
1810 TString spid[4]={"all","cpv","disp","both"} ;
1812 TList * arrayList = GetTracksTPCList(fVtxBin, fCentBin, fEMRPBin);
1814 for(Int_t evi=0; evi<arrayList->GetEntries();evi++){
1815 TClonesArray * mixTracks = static_cast<TClonesArray*>(arrayList->At(evi));
1816 for(Int_t i3=0; i3<mixTracks->GetEntriesFast(); i3++){
1817 TLorentzVector * track = (TLorentzVector*)mixTracks->At(i3);
1819 Double_t phiAssoc = track->Phi();
1820 Double_t etaAssoc = track->Eta();
1821 Double_t ptAssoc = track->Pt();
1823 Double_t ptAssocBin=GetAssocBin(ptAssoc) ;
1825 Double_t dPhi(0.), dEta(0.);
1827 for (int ipid = 0; ipid < 4; ipid++)
1829 if (fMEExists[ipid])
1831 dPhi = fMEPhi[ipid] - phiAssoc;
1832 while (dPhi > 1.5*TMath::Pi()) dPhi-=2*TMath::Pi();
1833 while (dPhi < -.5*TMath::Pi()) dPhi+=2*TMath::Pi();
1834 dEta = fMEEta[ipid] - etaAssoc;
1836 FillHistogram(Form("mix_%s_ptphieta_ptAssoc_%3.1f", spid[ipid].Data(), ptAssocBin), fMEPt[ipid], dPhi, dEta, 1./GetEfficiency(fMEPt[ipid]));
1837 FillHistogram(Form("mix_%s_ptphieta_ptAssoc_%3.1f_mod%i", spid[ipid].Data(), ptAssocBin, fMEModCase[ipid]), fMEPt[ipid], dPhi, dEta, 1./GetEfficiency(fMEPt[ipid]));
1838 FillHistogram(Form("mix_%s_ptphieta_ptAssoc_%3.1f_tpc%i", spid[ipid].Data(), ptAssocBin, CheckTriggerEta(fMEEta[ipid])), fMEPt[ipid], dPhi, dEta, 1./GetEfficiency(fMEPt[ipid]));
1845 //_______________________________________________________________________________
1846 TList* AliPHOSCorrelations::GetCaloPhotonsPHOSList(UInt_t vtxBin, UInt_t centBin, UInt_t rpBin){
1848 int offset = vtxBin * GetNumberOfCentralityBins() * fNEMRPBins + centBin * fNEMRPBins + rpBin;
1849 if( fCaloPhotonsPHOSLists->At(offset) ) {
1850 TList* list = dynamic_cast<TList*> (fCaloPhotonsPHOSLists->At(offset));
1853 else{ // no list for this bin has been created, yet
1854 TList* list = new TList();
1855 fCaloPhotonsPHOSLists->AddAt(list, offset);
1859 //_______________________________________________________________________________
1860 TList* AliPHOSCorrelations::GetTracksTPCList(UInt_t vtxBin, UInt_t centBin, UInt_t rpBin){
1862 int offset = vtxBin * GetNumberOfCentralityBins() * fNEMRPBins + centBin * fNEMRPBins + rpBin;
1863 if( fTracksTPCLists->At(offset) ) { // list exists
1864 TList* list = dynamic_cast<TList*> (fTracksTPCLists->At(offset));
1867 else { // no list for this bin has been created, yet
1868 TList* list = new TList();
1869 fTracksTPCLists->AddAt(list, offset);
1873 //_______________________________________________________________________________
1874 Double_t AliPHOSCorrelations::GetAssocBin(Double_t pt) const
1877 for(Int_t i=1; i<fAssocBins.GetSize(); i++){
1878 if(pt>fAssocBins.At(i-1) && pt<fAssocBins.At(i))
1879 return fAssocBins.At(i) ;
1881 return fAssocBins.At(fAssocBins.GetSize()-1) ;
1883 //_______________________________________________________________________________
1884 void AliPHOSCorrelations::FillTrackEtaPhi()
1886 // Distribution TPC's tracks by angles.
1887 for (Int_t i1=0; i1<fTracksTPC->GetEntriesFast(); i1++){
1888 TLorentzVector * track = (TLorentzVector*)fTracksTPC->At(i1);
1889 FillHistogram( "track_phieta", track->Phi(), track->Eta() );
1893 //_______________________________________________________________________________
1894 void AliPHOSCorrelations::UpdatePhotonLists()
1896 //Now we either add current events to stack or remove
1897 //If no photons in current event - no need to add it to mixed
1899 TList * arrayList = GetCaloPhotonsPHOSList(fVtxBin, fCentBin, fEMRPBin);
1901 AliInfo( Form("fCentBin=%d, fCentNMixed[]=%d",fCentBin,fCentNMixed[fCentBin]) );
1902 if(fCaloPhotonsPHOS->GetEntriesFast()>0)
1904 arrayList->AddFirst(fCaloPhotonsPHOS) ;
1905 fCaloPhotonsPHOS=0x0;
1906 if(arrayList->GetEntries() > fCentNMixed[fCentBin])
1907 { // Remove redundant events
1908 TClonesArray * tmp = static_cast<TClonesArray*>(arrayList->Last()) ;
1909 arrayList->RemoveLast() ;
1914 //_______________________________________________________________________________
1915 void AliPHOSCorrelations::UpdateTrackLists()
1917 //Now we either add current events to stack or remove
1918 //If no photons in current event - no need to add it to mixed
1920 TList * arrayList = GetTracksTPCList(fVtxBin, fCentBin, fEMRPBin);
1923 AliInfo( Form("fCentBin=%d, fCentNMixed[]=%d",fCentBin,fCentNMixed[fCentBin]) );
1924 if(fTracksTPC->GetEntriesFast()>0)
1927 arrayList->AddFirst(fTracksTPC) ;
1929 if(arrayList->GetEntries() > fCentNMixed[fCentBin])
1930 { // Remove redundant events
1931 TClonesArray * tmp = static_cast<TClonesArray*>(arrayList->Last()) ;
1932 arrayList->RemoveLast() ;
1937 //_______________________________________________________________________________
1938 Bool_t AliPHOSCorrelations::SelectESDTrack(AliESDtrack * t) const
1939 // Estimate if this track can be used for the RP calculation. If all right - return "TRUE"
1942 if(pt<0.5 || pt>20.) return kFALSE ;
1943 if(fabs( t->Eta() )>0.8) return kFALSE;
1944 if(!fESDtrackCuts->AcceptTrack(t)) return kFALSE ;
1947 //_______________________________________________________________________________
1948 Bool_t AliPHOSCorrelations::SelectAODTrack(AliAODTrack * t) const
1949 // Estimate if this track can be used for the RP calculation. If all right - return "TRUE"
1952 if(pt<0.5 || pt>20.) return kFALSE ;
1953 if(fabs( t->Eta() )>0.8) return kFALSE;
1954 if(fCheckHibridGlobal == kOnlyHibridTracks)
1956 if(!t->IsHybridGlobalConstrainedGlobal())
1960 if (fCheckHibridGlobal == kWithOutHibridTracks)
1962 if(t->IsHybridGlobalConstrainedGlobal())
1968 //_______________________________________________________________________________
1969 void AliPHOSCorrelations::LogProgress(int step)
1970 // Fill "step by step" hist
1972 //FillHistogram("hSelEvents", step+0.5, internalRunNumber-0.5);
1973 FillHistogram("hTotSelEvents", step+0.5);
1975 //_______________________________________________________________________________
1976 void AliPHOSCorrelations::LogSelection(int step, int internalRunNumber)
1978 // the +0.5 is not realy neccisarry, but oh well... -henrik
1979 FillHistogram("hSelEvents", step+0.5, internalRunNumber-0.5);
1980 //FillHistogram("hTotSelEvents", step+0.5);
1982 //_______________________________________________________________________________
1983 Bool_t AliPHOSCorrelations::TestMass(Double_t m, Double_t pt)
1985 //Check if mair in pi0 peak window
1986 //To make pT-dependent
1987 if (!fSigmaWidth) // Default big window
1989 FillHistogram("massWindow", fMassInvMean, fMassInvSigma);
1990 return (fMassInvMean-fMassInvSigma<m && m<fMassInvMean+fMassInvSigma) ;
1992 else // Parametrization
1994 FillHistogram("massWindow", MassMeanFunktion(pt), MassSigmaFunktion(pt)*fSigmaWidth);
1995 /*cout <<"MinMass: " << MassMeanFunktion(pt)-MassSigmaFunktion(pt)*fSigmaWidth
1998 <<" MaxMass "<< MassMeanFunktion(pt)+MassSigmaFunktion(pt)*fSigmaWidth<<endl;*/
1999 return ( MassMeanFunktion(pt)-MassSigmaFunktion(pt)*fSigmaWidth<m && m<MassMeanFunktion(pt)+MassSigmaFunktion(pt)*fSigmaWidth );
2002 //_______________________________________________________________________________
2003 Double_t AliPHOSCorrelations::MassMeanFunktion(Double_t &pt) const
2005 // Parametrization mean of mass window
2006 return ( fMassMeanP1+TMath::Power(1.25,-pt+fMassMeanP0) );
2008 //_______________________________________________________________________________
2009 Double_t AliPHOSCorrelations::MassSigmaFunktion(Double_t &pt) const
2011 // Parametrization sigma of mass window
2012 //TODO:: Kill falling at large pT.
2013 return ( fabs(fMassSigmaP0 + fMassSigmaP1*pt) );
2015 //_____________________________________________________________________________
2016 void AliPHOSCorrelations::FillHistogram(const char * key,Double_t x)const{
2018 TH1 * hist = dynamic_cast<TH1*>(fOutputContainer->FindObject(key)) ;
2022 AliError(Form("can not find histogram (of instance TH1) <%s> ",key)) ;
2024 //_____________________________________________________________________________
2025 void AliPHOSCorrelations::FillHistogram(const char * key,Double_t x,Double_t y)const{
2027 TH1 * th1 = dynamic_cast<TH1*> (fOutputContainer->FindObject(key));
2031 AliError(Form("can not find histogram (of instance TH1) <%s> ",key)) ;
2034 //_____________________________________________________________________________
2035 void AliPHOSCorrelations::FillHistogram(const char * key,Double_t x,Double_t y, Double_t z) const{
2036 //Fills 1D histograms with key
2037 TObject * obj = fOutputContainer->FindObject(key);
2039 TH2 * th2 = dynamic_cast<TH2*> (obj);
2041 th2->Fill(x, y, z) ;
2045 TH3 * th3 = dynamic_cast<TH3*> (obj);
2047 th3->Fill(x, y, z) ;
2051 AliError(Form("can not find histogram (of instance TH2) <%s> ",key)) ;
2053 //_____________________________________________________________________________
2054 void AliPHOSCorrelations::FillHistogram(const char * key,Double_t x,Double_t y, Double_t z, Double_t w) const{
2055 //Fills 1D histograms with key
2056 TObject * obj = fOutputContainer->FindObject(key);
2058 TH3 * th3 = dynamic_cast<TH3*> (obj);
2060 th3->Fill(x, y, z, w) ;
2064 AliError(Form("can not find histogram (of instance TH3) <%s> ",key)) ;
2066 //_____________________________________________________________________________
2067 void AliPHOSCorrelations::SetGeometry()
2069 // Initialize the PHOS geometry
2072 AliOADBContainer geomContainer("phosGeo");
2073 geomContainer.InitFromFile("$ALICE_ROOT/OADB/PHOS/PHOSGeometry.root","PHOSRotationMatrixes");
2074 TObjArray *matrixes = (TObjArray*)geomContainer.GetObject(fRunNumber,"PHOSRotationMatrixes");
2075 fPHOSGeo = AliPHOSGeometry::GetInstance("IHEP") ;
2076 for(Int_t mod=0; mod<5; mod++) {
2077 if(!matrixes->At(mod)) {
2079 AliInfo(Form("No PHOS Matrix for mod:%d, geo=%p\n", mod, fPHOSGeo));
2083 fPHOSGeo->SetMisalMatrix(((TGeoHMatrix*)matrixes->At(mod)),mod) ;
2085 AliInfo(Form("Adding PHOS Matrix for mod:%d, geo=%p\n", mod, fPHOSGeo));
2090 //_____________________________________________________________________________
2091 Double_t AliPHOSCorrelations::GetEfficiency(Double_t x) const {
2092 //Efficiency for Both2core only!
2095 // From 0 to 5 - 11h for different centrality.
2104 Double_t par0[9] = {-798863, 339.714, 6407.1, -457.778, 1283.65, -117.075, -19.3764, 0, 0};
2105 Double_t par1[9] = {-799344, -1852.1, 3326.29, -384.229, 504.046, 562.608, 130.518, 0, 0};
2106 Double_t par2[9] = {-858904, -1923.28, 5350.74, -568.946, 945.497, 419.647, 101.911, 0, 0};
2107 Double_t par3[9] = {-795652, -1495.97, 2926.46, -357.804, 478.961, 551.127, 128.86, 0, 0};
2108 Double_t par4[9] = {-891951, 279626, -123110, -5464.75, 27470.8, 283264, 15355.1, 192762, 44828.6};
2109 Double_t par5[9] = {-1.1094e+06, -986.915, 2127.71, -268.908, 375.594, 380.791, 89.4053, 0, 0};
2110 // Double_t par6[7] = {4.86106e+09, 4.47013e+08, -1.48079e+09, 1.47233e+08, -2.62356e+08, -1.00639e+08, -2.45629e+07, 0, 0};
2111 // Double_t par7[7] = {-1.36243e+06, -26011.1, 135838, -12161.3, 24956.8, 4985.4, 1285.57, 0, 0};
2113 // 8 for pPb13 and 0-100%
2114 Double_t par8[9] = {6.87095e+06, 8.36553e+06, -3.29572e+06, 2.18688e+06, -739490, 521666, 106661, 0, 0};
2117 Double_t* pFitPoint;
2119 if(fPeriod == kLHC11h)
2123 if (fCentrality<=5) pFitPoint = &par0[0];
2124 if (fCentrality>5 && fCentrality<=10) pFitPoint = &par1[0];
2125 if (fCentrality>10 && fCentrality<=20) pFitPoint = &par2[0];
2126 if (fCentrality>20 && fCentrality<=40) pFitPoint = &par3[0];
2127 if (fCentrality>40 && fCentrality<=60) pFitPoint = &par4[0];
2128 if (fCentrality>60) pFitPoint = &par5[0];
2131 for (int i = 0; i < 10; ++i)
2133 pFit[i] = *(pFitPoint+i);
2136 if (fCentrality>40 && fCentrality<=60)
2137 e = TMath::Exp(-(((((1.+(pFit[1]*x))+(pFit[2]*(x*x)))+(pFit[5]*(x*(x*x))))+(pFit[7]*(x*(x*(x*x)))))/((((pFit[3]*x)+(pFit[4]*(x*x)))+(pFit[6]*(x*(x*x))))+(pFit[8]*(x*(x*(x*x))))))) ;
2139 e = TMath::Exp(-((((1.+(pFit[1]*x))+(pFit[2]*(x*x)))+(pFit[5]*(x*(x*x))))/(((pFit[3]*x)+(pFit[4]*(x*x)))+(pFit[6]*(x*(x*x)))))) ;
2142 if( fPeriod == kLHC13 )
2144 pFitPoint = &par8[0];
2146 for( int i = 0; i < 10; i++ )
2148 pFit[i] = *(pFitPoint+i);
2151 e = TMath::Exp(-((((pFit[0]+(pFit[1]*x))+(pFit[2]*(x*x)))+(pFit[5]*(x*(x*x))))/(((1.+(pFit[3]*x))+(pFit[4]*(x*x)))+(pFit[6]*(x*(x*x)))))) ;
2156 AliWarning(Form("No efficiensy choise."));
2162 //_____________________________________________________________________________
2163 Int_t AliPHOSCorrelations::GetModCase(Int_t &mod1, Int_t &mod2) const {
2165 // Return modules pair namber.
2168 if(mod1 == 1) return 1;
2169 if(mod1 == 2) return 2;
2170 if(mod1 == 3) return 3;
2174 if(mod1 == 1 || mod2 == 1)
2175 if(mod1 == 2 || mod2 == 2)
2178 if(mod1 == 1 || mod2 == 1)
2179 if(mod1 == 3 || mod2 == 3)
2181 if(mod1 == 2 || mod2 == 2)
2182 if(mod1 == 3 || mod2 == 3)
2186 AliError(Form("No choise for mod1 = %i, mod2 = %i", mod1, mod2));
2189 //_____________________________________________________________________________
2190 void AliPHOSCorrelations::TestTrigger(){
2191 FillHistogram("hTriggerPassedEvents", 0); // All events
2192 if (fEvent->GetFiredTriggerClasses().Contains("PHI7") )
2193 FillHistogram("hTriggerPassedEvents", 1.); // 13 events
2194 if (fEvent->GetFiredTriggerClasses().Contains("PHS") )
2195 FillHistogram("hTriggerPassedEvents", 2.); // 11h events
2197 if (fEvent->GetFiredTriggerClasses().Contains("PHI7") || fEvent->GetFiredTriggerClasses().Contains("PHS"))
2201 AliInfo( Form("Event passed offline phos trigger test: %s ", fEvent->GetFiredTriggerClasses().Data() ) );
2203 //fPHOSEvent = true;
2205 //_____________________________________________________________________________
2206 Int_t AliPHOSCorrelations::CheckTriggerEta(Double_t eta){
2211 //_____________________________________________________________________________
2212 void AliPHOSCorrelations::TestPi0ME(Int_t ipid, TLorentzVector p12, Int_t modCase)
2214 Double_t phiTrigger=p12.Phi() ;
2215 Double_t etaTrigger=p12.Eta() ;
2216 Double_t pt=p12.Pt() ;
2218 if (pt >= fMEPt[ipid])
2221 fMEPhi[ipid] = phiTrigger;
2222 fMEEta[ipid] = etaTrigger;
2223 fMEModCase[ipid] = modCase;
2224 fMEExists[ipid] = true;
2227 //_____________________________________________________________________________
2228 void AliPHOSCorrelations::ZeroingVariables(){
2229 for (int i = 0; i < 4; ++i)
2231 fMEPhi[i] = fMEEta[i] = fMEPt[i] = -99;
2233 fMEExists[i] = false;