]>
Commit | Line | Data |
---|---|---|
ea3fd2d5 | 1 | // $Id$ |
2 | ||
ea3fd2d5 | 3 | #include "AliAnalysisTaskEMCALPi0PbPb.h" |
fa443410 | 4 | #include <TAxis.h> |
717fe7de | 5 | #include <TChain.h> |
6 | #include <TClonesArray.h> | |
f5d4ab70 | 7 | #include <TFile.h> |
296ea9b4 | 8 | #include <TGeoGlobalMagField.h> |
717fe7de | 9 | #include <TH1F.h> |
10 | #include <TH2F.h> | |
11 | #include <TList.h> | |
12 | #include <TLorentzVector.h> | |
f5d4ab70 | 13 | #include <TNtuple.h> |
296ea9b4 | 14 | #include <TProfile.h> |
b3ee6797 | 15 | #include <TString.h> |
296ea9b4 | 16 | #include <TVector2.h> |
717fe7de | 17 | #include "AliAODEvent.h" |
18 | #include "AliAODVertex.h" | |
19 | #include "AliAnalysisManager.h" | |
20 | #include "AliAnalysisTaskEMCALClusterizeFast.h" | |
296ea9b4 | 21 | #include "AliCDBManager.h" |
ea3fd2d5 | 22 | #include "AliCentrality.h" |
ea3fd2d5 | 23 | #include "AliEMCALGeoUtils.h" |
0ec74551 | 24 | #include "AliEMCALGeometry.h" |
0fbe8d4f | 25 | #include "AliEMCALRecPoint.h" |
296ea9b4 | 26 | #include "AliEMCALRecoUtils.h" |
3a952328 | 27 | #include "AliESDCaloTrigger.h" |
717fe7de | 28 | #include "AliESDEvent.h" |
5fe1ca23 | 29 | #include "AliESDUtils.h" |
717fe7de | 30 | #include "AliESDVertex.h" |
296ea9b4 | 31 | #include "AliESDtrack.h" |
0ec74551 | 32 | #include "AliESDtrackCuts.h" |
b6c599fe | 33 | #include "AliEventplane.h" |
296ea9b4 | 34 | #include "AliGeomManager.h" |
b3ee6797 | 35 | #include "AliInputEventHandler.h" |
43cfaa06 | 36 | #include "AliLog.h" |
296ea9b4 | 37 | #include "AliMagF.h" |
5fe1ca23 | 38 | #include "AliMultiplicity.h" |
0ec74551 | 39 | #include "AliTrackerBase.h" |
ea3fd2d5 | 40 | |
41 | ClassImp(AliAnalysisTaskEMCALPi0PbPb) | |
717fe7de | 42 | |
ea3fd2d5 | 43 | //________________________________________________________________________ |
44 | AliAnalysisTaskEMCALPi0PbPb::AliAnalysisTaskEMCALPi0PbPb(const char *name) | |
45 | : AliAnalysisTaskSE(name), | |
717fe7de | 46 | fCentVar("V0M"), |
47 | fCentFrom(0), | |
48 | fCentTo(100), | |
76332037 | 49 | fVtxZMin(-10), |
50 | fVtxZMax(+10), | |
43cfaa06 | 51 | fUseQualFlag(1), |
717fe7de | 52 | fClusName(), |
f5d4ab70 | 53 | fDoNtuple(0), |
a49742b5 | 54 | fDoAfterburner(0), |
f224d35b | 55 | fAsymMax(1), |
a49742b5 | 56 | fNminCells(1), |
3c661da5 | 57 | fMinE(0.100), |
f224d35b | 58 | fMinErat(0), |
59 | fMinEcc(0), | |
6bf90832 | 60 | fGeoName("EMCAL_FIRSTYEARV1"), |
b6c599fe | 61 | fMinNClusPerTr(50), |
296ea9b4 | 62 | fIsoDist(0.2), |
b3ee6797 | 63 | fTrClassNames(""), |
0ec74551 | 64 | fTrCuts(0), |
b6c599fe | 65 | fPrimTrCuts(0), |
66 | fDoTrMatGeom(0), | |
3a952328 | 67 | fTrainMode(0), |
68 | fMarkCells(), | |
69 | fMinL0Time(-1), | |
70 | fMaxL0Time(1024), | |
27c2e3d9 | 71 | fIsGeoMatsSet(0), |
d9f26424 | 72 | fNEvs(0), |
d595acbb | 73 | fGeom(0), |
296ea9b4 | 74 | fReco(0), |
717fe7de | 75 | fOutput(0), |
b3ee6797 | 76 | fTrClassNamesArr(0), |
717fe7de | 77 | fEsdEv(0), |
78 | fAodEv(0), | |
43cfaa06 | 79 | fRecPoints(0), |
717fe7de | 80 | fEsdClusters(0), |
81 | fEsdCells(0), | |
82 | fAodClusters(0), | |
286b47a5 | 83 | fAodCells(0), |
fa443410 | 84 | fPtRanges(0), |
296ea9b4 | 85 | fSelTracks(0), |
b6c599fe | 86 | fSelPrimTracks(0), |
3a952328 | 87 | fNAmpInTrigger(0), |
88 | fAmpInTrigger(0), | |
788ca675 | 89 | fNtuple(0), |
90 | fHeader(0), | |
91 | fPrimVert(0), | |
92 | fSpdVert(0), | |
93 | fTpcVert(0), | |
94 | fClusters(0), | |
3a952328 | 95 | fTriggers(0), |
fa443410 | 96 | fHCuts(0x0), |
97 | fHVertexZ(0x0), | |
76332037 | 98 | fHVertexZ2(0x0), |
fa443410 | 99 | fHCent(0x0), |
100 | fHCentQual(0x0), | |
b3ee6797 | 101 | fHTclsBeforeCuts(0x0), |
102 | fHTclsAfterCuts(0x0), | |
d595acbb | 103 | fHColuRow(0x0), |
104 | fHColuRowE(0x0), | |
105 | fHCellMult(0x0), | |
106 | fHCellE(0x0), | |
107 | fHCellH(0x0), | |
6eb6260e | 108 | fHCellM(0x0), |
a49742b5 | 109 | fHCellM2(0x0), |
296ea9b4 | 110 | fHCellFreqNoCut(0x0), |
2e4d8148 | 111 | fHCellFreqCut100M(0x0), |
112 | fHCellFreqCut300M(0x0), | |
113 | fHCellFreqE(0x0), | |
296ea9b4 | 114 | fHCellCheckE(0x0), |
6eb6260e | 115 | fHClustEccentricity(0), |
fa443410 | 116 | fHClustEtaPhi(0x0), |
117 | fHClustEnergyPt(0x0), | |
118 | fHClustEnergySigma(0x0), | |
d595acbb | 119 | fHClustSigmaSigma(0x0), |
6eb6260e | 120 | fHClustNCellEnergyRatio(0x0), |
b6c599fe | 121 | fHMatchDr(0x0), |
122 | fHMatchDz(0x0), | |
123 | fHMatchEp(0x0), | |
fa443410 | 124 | fHPionEtaPhi(0x0), |
125 | fHPionMggPt(0x0), | |
6eb6260e | 126 | fHPionMggAsym(0x0), |
127 | fHPionMggDgg(0x0) | |
ea3fd2d5 | 128 | { |
717fe7de | 129 | // Constructor. |
ea3fd2d5 | 130 | |
d595acbb | 131 | if (!name) |
132 | return; | |
133 | SetName(name); | |
ea3fd2d5 | 134 | DefineInput(0, TChain::Class()); |
ea3fd2d5 | 135 | DefineOutput(1, TList::Class()); |
3a952328 | 136 | fBranchNames="ESD:AliESDRun.,AliESDHeader.,PrimaryVertex.,SPDVertex.,TPCVertex.,EMCALCells.,Tracks EMCALTrigger." |
296ea9b4 | 137 | "AOD:header,vertices,emcalCells,tracks"; |
ea3fd2d5 | 138 | } |
ea3fd2d5 | 139 | |
717fe7de | 140 | //________________________________________________________________________ |
141 | AliAnalysisTaskEMCALPi0PbPb::~AliAnalysisTaskEMCALPi0PbPb() | |
142 | { | |
143 | // Destructor. | |
ea3fd2d5 | 144 | |
b3ee6797 | 145 | if (!AliAnalysisManager::GetAnalysisManager()->IsProofMode()) { |
146 | delete fOutput; fOutput = 0; | |
147 | } | |
fa443410 | 148 | delete fPtRanges; fPtRanges = 0; |
d595acbb | 149 | delete fGeom; fGeom = 0; |
296ea9b4 | 150 | delete fReco; fReco = 0; |
b3ee6797 | 151 | delete fTrClassNamesArr; |
296ea9b4 | 152 | delete fSelTracks; |
b6c599fe | 153 | delete fSelPrimTracks; |
3a952328 | 154 | delete [] fAmpInTrigger; |
d595acbb | 155 | delete [] fHColuRow; |
156 | delete [] fHColuRowE; | |
157 | delete [] fHCellMult; | |
296ea9b4 | 158 | delete [] fHCellFreqNoCut; |
2e4d8148 | 159 | delete [] fHCellFreqCut100M; |
160 | delete [] fHCellFreqCut300M; | |
ea3fd2d5 | 161 | } |
717fe7de | 162 | |
ea3fd2d5 | 163 | //________________________________________________________________________ |
164 | void AliAnalysisTaskEMCALPi0PbPb::UserCreateOutputObjects() | |
165 | { | |
717fe7de | 166 | // Create user objects here. |
ea3fd2d5 | 167 | |
f3582e89 | 168 | cout << "AliAnalysisTaskEMCALPi0PbPb: Input settings" << endl; |
169 | cout << " fCentVar: " << fCentVar << endl; | |
170 | cout << " fCentFrom: " << fCentFrom << endl; | |
171 | cout << " fCentTo: " << fCentTo << endl; | |
172 | cout << " fVtxZMin: " << fVtxZMin << endl; | |
173 | cout << " fVtxZMax: " << fVtxZMax << endl; | |
174 | cout << " fUseQualFlag: " << fUseQualFlag << endl; | |
175 | cout << " fClusName: \"" << fClusName << "\"" << endl; | |
176 | cout << " fDoNtuple: " << fDoNtuple << endl; | |
177 | cout << " fDoAfterburner: " << fDoAfterburner << endl; | |
178 | cout << " fAsymMax: " << fAsymMax << endl; | |
179 | cout << " fNminCells: " << fNminCells << endl; | |
180 | cout << " fMinE: " << fMinE << endl; | |
181 | cout << " fMinErat: " << fMinErat << endl; | |
182 | cout << " fMinEcc: " << fMinEcc << endl; | |
183 | cout << " fGeoName: \"" << fGeoName << "\"" << endl; | |
184 | cout << " fMinNClusPerTr: " << fMinNClusPerTr << endl; | |
185 | cout << " fTrClassNames: \"" << fTrClassNames << "\"" << endl; | |
f3582e89 | 186 | cout << " fTrCuts: " << fTrCuts << endl; |
187 | cout << " fPrimTrCuts: " << fPrimTrCuts << endl; | |
188 | cout << " fDoTrMatGeom: " << fDoTrMatGeom << endl; | |
189 | cout << " fTrainMode: " << fTrainMode << endl; | |
190 | cout << " fMarkCells: " << fMarkCells << endl; | |
191 | cout << " fMinL0Time: " << fMinL0Time << endl; | |
192 | cout << " fMaxL0Time: " << fMaxL0Time << endl; | |
193 | ||
296ea9b4 | 194 | fGeom = new AliEMCALGeoUtils(fGeoName,"EMCAL"); |
195 | fReco = new AliEMCALRecoUtils(); | |
b3ee6797 | 196 | fTrClassNamesArr = fTrClassNames.Tokenize(" "); |
717fe7de | 197 | fOutput = new TList(); |
198 | fOutput->SetOwner(); | |
296ea9b4 | 199 | fSelTracks = new TObjArray; |
b6c599fe | 200 | fSelPrimTracks = new TObjArray; |
ea3fd2d5 | 201 | |
f5d4ab70 | 202 | if (fDoNtuple) { |
203 | TFile *f = OpenFile(1); | |
6bf90832 | 204 | if (f) { |
205 | f->SetCompressionLevel(2); | |
788ca675 | 206 | fNtuple = new TTree(Form("tree%.0fto%.0f",fCentFrom,fCentTo), "StandaloneTree"); |
6bf90832 | 207 | fNtuple->SetDirectory(f); |
3a952328 | 208 | if (fTrainMode) { |
209 | fNtuple->SetAutoFlush(-2*1024*1024); | |
210 | fNtuple->SetAutoSave(0); | |
211 | } else { | |
212 | fNtuple->SetAutoFlush(-32*1024*1024); | |
213 | fNtuple->SetAutoSave(0); | |
214 | } | |
215 | ||
788ca675 | 216 | fHeader = new AliStaHeader; |
217 | fNtuple->Branch("header", &fHeader, 16*1024, 99); | |
218 | fPrimVert = new AliStaVertex; | |
219 | fNtuple->Branch("primv", &fPrimVert, 16*1024, 99); | |
220 | fSpdVert = new AliStaVertex; | |
221 | fNtuple->Branch("spdv", &fSpdVert, 16*1024, 99); | |
222 | fTpcVert = new AliStaVertex; | |
223 | fNtuple->Branch("tpcv", &fTpcVert, 16*1024, 99); | |
224 | if (TClass::GetClass("AliStaCluster")) | |
225 | TClass::GetClass("AliStaCluster")->IgnoreTObjectStreamer(); | |
3a952328 | 226 | fClusters = new TClonesArray("AliStaCluster"); |
788ca675 | 227 | fNtuple->Branch("clusters", &fClusters, 8*16*1024, 99); |
3a952328 | 228 | if (TClass::GetClass("AliStaTrigger")) |
229 | TClass::GetClass("AliStaTrigger")->IgnoreTObjectStreamer(); | |
230 | fTriggers = new TClonesArray("AliStaTrigger"); | |
231 | fNtuple->Branch("l0prim", &fTriggers, 16*1024, 99); | |
6bf90832 | 232 | } |
f5d4ab70 | 233 | } |
234 | ||
002dcebe | 235 | AliEMCALEMCGeometry *emc = fGeom->GetEMCGeometry(); |
b6c599fe | 236 | Double_t phimin = emc->GetArm1PhiMin()*TMath::DegToRad(); |
237 | Double_t phimax = emc->GetArm1PhiMax()*TMath::DegToRad(); | |
002dcebe | 238 | |
d595acbb | 239 | // histograms |
a49742b5 | 240 | TH1::SetDefaultSumw2(kTRUE); |
241 | TH2::SetDefaultSumw2(kTRUE); | |
b3ee6797 | 242 | fHCuts = new TH1F("hEventCuts","",5,0.5,5.5); |
243 | fHCuts->GetXaxis()->SetBinLabel(1,"All"); | |
244 | fHCuts->GetXaxis()->SetBinLabel(2,"PS"); | |
245 | fHCuts->GetXaxis()->SetBinLabel(3,Form("%s: %.0f-%.0f",fCentVar.Data(),fCentFrom,fCentTo)); | |
246 | fHCuts->GetXaxis()->SetBinLabel(4,"QFlag"); | |
247 | fHCuts->GetXaxis()->SetBinLabel(5,Form("zvtx: %.0f-%.0f",fVtxZMin,fVtxZMax)); | |
fa443410 | 248 | fOutput->Add(fHCuts); |
249 | fHVertexZ = new TH1F("hVertexZBeforeCut","",100,-25,25); | |
250 | fHVertexZ->SetXTitle("z [cm]"); | |
251 | fOutput->Add(fHVertexZ); | |
76332037 | 252 | fHVertexZ2 = new TH1F("hVertexZAfterCut","",100,-25,25); |
253 | fHVertexZ2->SetXTitle("z [cm]"); | |
254 | fOutput->Add(fHVertexZ2); | |
f2b8fca6 | 255 | fHCent = new TH1F("hCentBeforeCut","",102,-1,101); |
fa443410 | 256 | fHCent->SetXTitle(fCentVar.Data()); |
76332037 | 257 | fOutput->Add(fHCent); |
f2b8fca6 | 258 | fHCentQual = new TH1F("hCentAfterCut","",102,-1,101); |
fa443410 | 259 | fHCentQual->SetXTitle(fCentVar.Data()); |
260 | fOutput->Add(fHCentQual); | |
2e4d8148 | 261 | fHTclsBeforeCuts = new TH1F("hTclsBeforeCuts","",fTrClassNamesArr->GetEntries(),0.5,0.5+fTrClassNamesArr->GetEntries()); |
262 | fHTclsAfterCuts = new TH1F("hTclsAfterCuts","",fTrClassNamesArr->GetEntries(),0.5,0.5+fTrClassNamesArr->GetEntries()); | |
b3ee6797 | 263 | for (Int_t i = 0; i<fTrClassNamesArr->GetEntries(); ++i) { |
264 | const char *name = fTrClassNamesArr->At(i)->GetName(); | |
265 | fHTclsBeforeCuts->GetXaxis()->SetBinLabel(1+i,name); | |
266 | fHTclsAfterCuts->GetXaxis()->SetBinLabel(1+i,name); | |
267 | } | |
268 | fOutput->Add(fHTclsBeforeCuts); | |
269 | fOutput->Add(fHTclsAfterCuts); | |
90d5b88b | 270 | |
d595acbb | 271 | // histograms for cells |
296ea9b4 | 272 | Int_t nsm = fGeom->GetEMCGeometry()->GetNumberOfSuperModules(); |
273 | fHColuRow = new TH2*[nsm]; | |
274 | fHColuRowE = new TH2*[nsm]; | |
275 | fHCellMult = new TH1*[nsm]; | |
276 | for (Int_t i = 0; i<nsm; ++i) { | |
2e4d8148 | 277 | fHColuRow[i] = new TH2F(Form("hColRow_Mod%d",i),"",48,0,48,24,0.,24); |
d595acbb | 278 | fHColuRow[i]->SetTitle(Form("Module %d: Occupancy", i)); |
279 | fHColuRow[i]->SetXTitle("col (i#eta)"); | |
280 | fHColuRow[i]->SetYTitle("row (i#phi)"); | |
2e4d8148 | 281 | fHColuRowE[i] = new TH2F(Form("hColRowE_Mod%d", i),"",48,0,48,24,0,24); |
90d5b88b | 282 | fHColuRowE[i]->SetTitle(Form("Module %d: Cell energy",i)); |
d595acbb | 283 | fHColuRowE[i]->SetXTitle("col (i#eta)"); |
284 | fHColuRowE[i]->SetYTitle("row (i#phi)"); | |
285 | fHCellMult[i] = new TH1F(Form("hCellMult_Mod%d",i),"",1000,0,1000); | |
286 | fHCellMult[i]->SetTitle(Form("Module %d: Cell multiplicity",i)); | |
90d5b88b | 287 | fHCellMult[i]->SetXTitle("# of cells"); |
d595acbb | 288 | fOutput->Add(fHColuRow[i]); |
289 | fOutput->Add(fHColuRowE[i]); | |
290 | fOutput->Add(fHCellMult[i]); | |
291 | } | |
2e4d8148 | 292 | fHCellE = new TH1F("hCellE","",250,0.,25.); |
d595acbb | 293 | fHCellE->SetXTitle("E_{cell} [GeV]"); |
294 | fOutput->Add(fHCellE); | |
fb5a0b69 | 295 | fHCellH = new TH1F ("hCellHighestE","",250,0.,25.); |
6eb6260e | 296 | fHCellH->SetXTitle("E^{max}_{cell} [GeV]"); |
d595acbb | 297 | fOutput->Add(fHCellH); |
fb5a0b69 | 298 | fHCellM = new TH1F ("hCellMeanEperHitCell","",250,0.,2.5); |
6eb6260e | 299 | fHCellM->SetXTitle("#LT E_{cell}#GT [GeV]"); |
300 | fOutput->Add(fHCellM); | |
fb5a0b69 | 301 | fHCellM2 = new TH1F ("hCellMeanEperAllCells","",250,0.,1); |
f4ec884e | 302 | fHCellM2->SetXTitle("1/N_{cells} #Sigma E_{cell} [GeV]"); |
a49742b5 | 303 | fOutput->Add(fHCellM2); |
90d5b88b | 304 | |
2e4d8148 | 305 | fHCellFreqNoCut = new TH1*[nsm]; |
306 | fHCellFreqCut100M = new TH1*[nsm]; | |
307 | fHCellFreqCut300M = new TH1*[nsm]; | |
308 | fHCellFreqE = new TH1*[nsm]; | |
296ea9b4 | 309 | for (Int_t i = 0; i<nsm; ++i){ |
310 | Double_t lbin = i*24*48-0.5; | |
311 | Double_t hbin = lbin+24*48; | |
2e4d8148 | 312 | fHCellFreqNoCut[i] = new TH1F(Form("hCellFreqNoCut_SM%d",i), |
313 | Form("Frequency SM%d (no cut);id;#",i), 1152, lbin, hbin); | |
314 | fHCellFreqCut100M[i] = new TH1F(Form("hCellFreqCut100M_SM%d",i), | |
315 | Form("Frequency SM%d (>0.1GeV);id;#",i), 1152, lbin, hbin); | |
316 | fHCellFreqCut300M[i] = new TH1F(Form("hCellFreqCut300M_SM%d",i), | |
317 | Form("Frequency SM%d (>0.3GeV);id;#",i), 1152, lbin, hbin); | |
318 | fHCellFreqE[i] = new TH1F(Form("hCellFreqE_SM%d",i), | |
319 | Form("Frequency SM%d (E weighted);id;#",i), 1152, lbin, hbin); | |
296ea9b4 | 320 | fOutput->Add(fHCellFreqNoCut[i]); |
2e4d8148 | 321 | fOutput->Add(fHCellFreqCut100M[i]); |
322 | fOutput->Add(fHCellFreqCut300M[i]); | |
323 | fOutput->Add(fHCellFreqE[i]); | |
296ea9b4 | 324 | } |
3a952328 | 325 | if (!fMarkCells.IsNull()) { |
296ea9b4 | 326 | fHCellCheckE = new TH1*[24*48*nsm]; |
408dc04c | 327 | memset(fHCellCheckE,0,24*48*nsm*sizeof(TH1*)); |
3a952328 | 328 | TObjArray *cells = fMarkCells.Tokenize(" "); |
329 | Int_t n = cells->GetEntries(); | |
330 | Int_t *tcs = new Int_t[n]; | |
331 | for (Int_t i=0;i<n;++i) { | |
332 | TString name(cells->At(i)->GetName()); | |
333 | tcs[i]=name.Atoi(); | |
334 | } | |
335 | for (Int_t i = 0; i<n; ++i) { | |
296ea9b4 | 336 | Int_t c=tcs[i]; |
337 | if (c<24*48*nsm) { | |
3a952328 | 338 | fHCellCheckE[i] = new TH1F(Form("hCellE_id%d",c), Form("Cell %d;E [GeV/c];#",c), 1000, 0, 10); |
296ea9b4 | 339 | fOutput->Add(fHCellCheckE[i]); |
340 | } | |
341 | } | |
3a952328 | 342 | delete cells; |
343 | delete [] tcs; | |
296ea9b4 | 344 | } |
345 | ||
d595acbb | 346 | // histograms for clusters |
3a952328 | 347 | if (!fTrainMode) { |
348 | fHClustEccentricity = new TH1F("hClustEccentricity","",100,-0.1,1.1); | |
349 | fHClustEccentricity->SetXTitle("#epsilon_{C}"); | |
350 | fOutput->Add(fHClustEccentricity); | |
351 | fHClustEtaPhi = new TH2F("hClustEtaPhi","",500,-0.8,0.8,100*nsm,phimin,phimax); | |
352 | fHClustEtaPhi->SetXTitle("#eta"); | |
353 | fHClustEtaPhi->SetYTitle("#varphi"); | |
354 | fOutput->Add(fHClustEtaPhi); | |
355 | fHClustEnergyPt = new TH2F("hClustEnergyPt","",250,0,50,250,0,50); | |
356 | fHClustEnergyPt->SetXTitle("E [GeV]"); | |
357 | fHClustEnergyPt->SetYTitle("p_{T} [GeV/c]"); | |
358 | fOutput->Add(fHClustEnergyPt); | |
359 | fHClustEnergySigma = new TH2F("hClustEnergySigma","",100,0,100,500,0,50); | |
360 | fHClustEnergySigma->SetXTitle("E_{C} * #sigma_{max} [GeV*cm]"); | |
361 | fHClustEnergySigma->SetYTitle("E_{C} [GeV]"); | |
362 | fOutput->Add(fHClustEnergySigma); | |
363 | fHClustSigmaSigma = new TH2F("hClustSigmaSigma","",500,0,50,500,0,50); | |
364 | fHClustSigmaSigma->SetXTitle("#lambda_{0} [cm]"); | |
365 | fHClustSigmaSigma->SetYTitle("#sigma_{max} [cm]"); | |
366 | fOutput->Add(fHClustSigmaSigma); | |
367 | fHClustNCellEnergyRatio = new TH2F("hClustNCellEnergyRatio","",27,-0.5,26.5,101,-0.05,1.05); | |
368 | fHClustNCellEnergyRatio->SetXTitle("N_{cells}"); | |
369 | fHClustNCellEnergyRatio->SetYTitle("E^{max}_{cell}/E_{clus}"); | |
370 | fOutput->Add(fHClustNCellEnergyRatio); | |
371 | } | |
90d5b88b | 372 | |
b6c599fe | 373 | // histograms for track matching |
374 | fHMatchDr = new TH1F("hMatchDrDist",";dR [cm]",500,0,200); | |
375 | fOutput->Add(fHMatchDr); | |
376 | fHMatchDz = new TH1F("hMatchDzDist",";dZ [cm]",500,-100,100); | |
377 | fOutput->Add(fHMatchDz); | |
378 | fHMatchEp = new TH1F("hMatchEpDist",";E/p",100,0,10); | |
379 | fOutput->Add(fHMatchEp); | |
3a952328 | 380 | |
d595acbb | 381 | // histograms for pion candidates |
3a952328 | 382 | if (!fTrainMode) { |
383 | fHPionEtaPhi = new TH2F("hPionEtaPhi","",100,-0.8,0.8,100*nsm,phimin,phimax); | |
384 | fHPionEtaPhi->SetXTitle("#eta_{#gamma#gamma}"); | |
385 | fHPionEtaPhi->SetYTitle("#varphi_{#gamma#gamma}"); | |
386 | fOutput->Add(fHPionEtaPhi); | |
387 | fHPionMggPt = new TH2F("hPionMggPt","",1000,0,2,100,0,20.0); | |
388 | fHPionMggPt->SetXTitle("M_{#gamma#gamma} [GeV/c^{2}]"); | |
389 | fHPionMggPt->SetYTitle("p_{T}^{#gamma#gamma} [GeV/c]"); | |
390 | fOutput->Add(fHPionMggPt); | |
391 | fHPionMggAsym = new TH2F("hPionMggAsym","",1000,0,2,100,0,1); | |
392 | fHPionMggAsym->SetXTitle("M_{#gamma#gamma} [GeV/c^{2}]"); | |
393 | fHPionMggAsym->SetYTitle("Z_{#gamma#gamma} [GeV]"); | |
394 | fOutput->Add(fHPionMggAsym); | |
395 | fHPionMggDgg = new TH2F("hPionMggDgg","",1000,0,2,100,0,10); | |
396 | fHPionMggDgg->SetXTitle("M_{#gamma#gamma} [GeV/c^{2}]"); | |
397 | fHPionMggDgg->SetYTitle("opening angle [grad]"); | |
398 | fOutput->Add(fHPionMggDgg); | |
399 | const Int_t nbins = 20; | |
400 | Double_t xbins[nbins] = {0.5,1,1.5,2,2.5,3,3.5,4,4.5,5,6,7,8,9,10,12.5,15,20,25,50}; | |
401 | fPtRanges = new TAxis(nbins-1,xbins); | |
402 | for (Int_t i = 0; i<=nbins; ++i) { | |
403 | fHPionInvMasses[i] = new TH1F(Form("hPionInvMass%d",i),"",1000,0,2); | |
404 | fHPionInvMasses[i]->SetXTitle("M_{#gamma#gamma} [GeV/c^{2}]"); | |
405 | if (i==0) | |
406 | fHPionInvMasses[i]->SetTitle(Form("0 < p_{T}^{#gamma#gamma} <%.1f",xbins[0])); | |
407 | else if (i==nbins) | |
408 | fHPionInvMasses[i]->SetTitle(Form("p_{T}^{#gamma#gamma} > 50")); | |
409 | else | |
410 | fHPionInvMasses[i]->SetTitle(Form("%.1f < p_{T}^{#gamma#gamma} <%.1f",xbins[i-1],xbins[i])); | |
411 | fOutput->Add(fHPionInvMasses[i]); | |
412 | } | |
fa443410 | 413 | } |
296ea9b4 | 414 | |
717fe7de | 415 | PostData(1, fOutput); |
ea3fd2d5 | 416 | } |
417 | ||
418 | //________________________________________________________________________ | |
419 | void AliAnalysisTaskEMCALPi0PbPb::UserExec(Option_t *) | |
420 | { | |
717fe7de | 421 | // Called for each event. |
422 | ||
423 | if (!InputEvent()) | |
424 | return; | |
425 | ||
43cfaa06 | 426 | AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager(); |
427 | fEsdEv = dynamic_cast<AliESDEvent*>(InputEvent()); | |
44310bce | 428 | UInt_t offtrigger = 0; |
43cfaa06 | 429 | if (fEsdEv) { |
430 | am->LoadBranch("AliESDRun."); | |
431 | am->LoadBranch("AliESDHeader."); | |
9809ed8c | 432 | UInt_t mask1 = fEsdEv->GetESDRun()->GetDetectorsInDAQ(); |
433 | UInt_t mask2 = fEsdEv->GetESDRun()->GetDetectorsInReco(); | |
434 | Bool_t desc1 = (mask1 >> 18) & 0x1; | |
435 | Bool_t desc2 = (mask2 >> 18) & 0x1; | |
436 | if (desc1==0 || desc2==0) { //AliDAQ::OfflineModuleName(180=="EMCAL" | |
437 | AliError(Form("EMCAL not in DAQ/RECO: %u (%u)/%u (%u)", | |
438 | mask1, fEsdEv->GetESDRun()->GetDetectorsInReco(), | |
439 | mask2, fEsdEv->GetESDRun()->GetDetectorsInDAQ())); | |
44310bce | 440 | return; |
441 | } | |
442 | offtrigger = ((AliInputEventHandler*)(am->GetInputEventHandler()))->IsEventSelected(); | |
43cfaa06 | 443 | } else { |
444 | fAodEv = dynamic_cast<AliAODEvent*>(InputEvent()); | |
408dc04c | 445 | if (!fAodEv) { |
446 | AliFatal("Neither ESD nor AOD event found"); | |
447 | return; | |
448 | } | |
43cfaa06 | 449 | am->LoadBranch("header"); |
44310bce | 450 | offtrigger = fAodEv->GetHeader()->GetOfflineTrigger(); |
451 | } | |
452 | if (offtrigger & AliVEvent::kFastOnly) { | |
453 | AliWarning(Form("EMCAL not in fast only partition")); | |
454 | return; | |
43cfaa06 | 455 | } |
b6c599fe | 456 | if (fDoTrMatGeom && !AliGeomManager::GetGeometry()) { // get geometry |
27c2e3d9 | 457 | AliWarning("Accessing geometry from OCDB, this is not very efficient!"); |
458 | AliCDBManager *cdb = AliCDBManager::Instance(); | |
459 | if (!cdb->IsDefaultStorageSet()) | |
460 | cdb->SetDefaultStorage("raw://"); | |
461 | Int_t runno = InputEvent()->GetRunNumber(); | |
462 | if (runno != cdb->GetRun()) | |
463 | cdb->SetRun(runno); | |
464 | AliGeomManager::LoadGeometry(); | |
465 | } | |
466 | ||
467 | if (!AliGeomManager::GetGeometry()&&!fIsGeoMatsSet) { // set misalignment matrices (stored in first event) | |
468 | Int_t nsm = fGeom->GetEMCGeometry()->GetNumberOfSuperModules(); | |
9809ed8c | 469 | for (Int_t i=0; i<nsm; ++i) { |
470 | const TGeoHMatrix *geom = 0; | |
471 | if (fEsdEv) | |
472 | geom = fEsdEv->GetESDRun()->GetEMCALMatrix(i); | |
473 | else | |
474 | geom = fAodEv->GetHeader()->GetEMCALMatrix(i); | |
475 | if (!geom) | |
476 | continue; | |
477 | geom->Print(); | |
478 | fGeom->SetMisalMatrix(geom,i); | |
27c2e3d9 | 479 | } |
480 | fIsGeoMatsSet = kTRUE; | |
481 | } | |
482 | ||
483 | if (!TGeoGlobalMagField::Instance()->GetField()) { // construct field map | |
484 | if (fEsdEv) { | |
485 | const AliESDRun *erun = fEsdEv->GetESDRun(); | |
486 | AliMagF *field = AliMagF::CreateFieldMap(erun->GetCurrentL3(), | |
487 | erun->GetCurrentDip(), | |
488 | AliMagF::kConvLHC, | |
489 | kFALSE, | |
490 | erun->GetBeamEnergy(), | |
491 | erun->GetBeamType()); | |
492 | TGeoGlobalMagField::Instance()->SetField(field); | |
493 | } else { | |
494 | Double_t pol = -1; //polarity | |
495 | Double_t be = -1; //beam energy | |
496 | AliMagF::BeamType_t btype = AliMagF::kBeamTypepp; | |
497 | Int_t runno = fAodEv->GetRunNumber(); | |
498 | if (runno>=136851 && runno<138275) { | |
499 | pol = -1; | |
500 | be = 2760; | |
501 | btype = AliMagF::kBeamTypeAA; | |
502 | } else if (runno>=138275 && runno<=139517) { | |
503 | pol = +1; | |
504 | be = 2760; | |
505 | btype = AliMagF::kBeamTypeAA; | |
506 | } else { | |
507 | AliError(Form("Do not know the bfield parameters for run %d! Using defaults!!!", runno)); | |
508 | } | |
509 | TGeoGlobalMagField::Instance()->SetField(new AliMagF("Maps","Maps", pol, pol, AliMagF::k5kG, btype, be)); | |
510 | } | |
511 | } | |
512 | ||
b3ee6797 | 513 | Int_t cut = 1; |
514 | fHCuts->Fill(cut++); | |
515 | ||
516 | TString trgclasses; | |
517 | AliESDHeader *h = dynamic_cast<AliESDHeader*>(InputEvent()->GetHeader()); | |
518 | if (h) { | |
519 | trgclasses = fEsdEv->GetFiredTriggerClasses(); | |
520 | } else { | |
521 | AliAODHeader *h2 = dynamic_cast<AliAODHeader*>(InputEvent()->GetHeader()); | |
522 | if (h2) | |
523 | trgclasses = h2->GetFiredTriggerClasses(); | |
524 | } | |
525 | for (Int_t i = 0; i<fTrClassNamesArr->GetEntries(); ++i) { | |
526 | const char *name = fTrClassNamesArr->At(i)->GetName(); | |
527 | if (trgclasses.Contains(name)) | |
528 | fHTclsBeforeCuts->Fill(1+i); | |
529 | } | |
530 | ||
531 | UInt_t res = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected(); | |
532 | if (res==0) | |
533 | return; | |
534 | ||
fa443410 | 535 | fHCuts->Fill(cut++); |
286b47a5 | 536 | |
717fe7de | 537 | const AliCentrality *centP = InputEvent()->GetCentrality(); |
538 | Double_t cent = centP->GetCentralityPercentileUnchecked(fCentVar); | |
fa443410 | 539 | fHCent->Fill(cent); |
717fe7de | 540 | if (cent<fCentFrom||cent>fCentTo) |
286b47a5 | 541 | return; |
43cfaa06 | 542 | |
fa443410 | 543 | fHCuts->Fill(cut++); |
286b47a5 | 544 | |
43cfaa06 | 545 | if (fUseQualFlag) { |
546 | if (centP->GetQuality()>0) | |
547 | return; | |
548 | } | |
286b47a5 | 549 | |
fa443410 | 550 | fHCentQual->Fill(cent); |
551 | fHCuts->Fill(cut++); | |
717fe7de | 552 | |
43cfaa06 | 553 | if (fEsdEv) { |
fa443410 | 554 | am->LoadBranch("PrimaryVertex."); |
555 | am->LoadBranch("SPDVertex."); | |
556 | am->LoadBranch("TPCVertex."); | |
43cfaa06 | 557 | } else { |
558 | fAodEv = dynamic_cast<AliAODEvent*>(InputEvent()); | |
559 | am->LoadBranch("vertices"); | |
3c661da5 | 560 | if (!fAodEv) return; |
43cfaa06 | 561 | } |
562 | ||
563 | const AliVVertex *vertex = InputEvent()->GetPrimaryVertex(); | |
717fe7de | 564 | if (!vertex) |
565 | return; | |
566 | ||
fa443410 | 567 | fHVertexZ->Fill(vertex->GetZ()); |
286b47a5 | 568 | |
717fe7de | 569 | if(vertex->GetZ()<fVtxZMin||vertex->GetZ()>fVtxZMax) |
570 | return; | |
286b47a5 | 571 | |
fa443410 | 572 | fHCuts->Fill(cut++); |
76332037 | 573 | fHVertexZ2->Fill(vertex->GetZ()); |
717fe7de | 574 | |
b3ee6797 | 575 | // count number of accepted events |
576 | ++fNEvs; | |
577 | ||
578 | for (Int_t i = 0; i<fTrClassNamesArr->GetEntries(); ++i) { | |
579 | const char *name = fTrClassNamesArr->At(i)->GetName(); | |
580 | if (trgclasses.Contains(name)) | |
581 | fHTclsAfterCuts->Fill(1+i); | |
582 | } | |
583 | ||
43cfaa06 | 584 | fRecPoints = 0; // will be set if fClusName is given and AliAnalysisTaskEMCALClusterizeFast is used |
b6c599fe | 585 | fEsdClusters = 0; // will be set if ESD input used and if fRecPoints are not set or if clusters are attached |
43cfaa06 | 586 | fEsdCells = 0; // will be set if ESD input used |
b6c599fe | 587 | fAodClusters = 0; // will be set if AOD input used and if fRecPoints are not set or if clusters are attached |
43cfaa06 | 588 | // or if fClusName is given and AliAnalysisTaskEMCALClusterizeFast in AOD output mode |
589 | fAodCells = 0; // will be set if AOD input used | |
717fe7de | 590 | |
43cfaa06 | 591 | // deal with special output from AliAnalysisTaskEMCALClusterizeFast first |
592 | Bool_t clusattached = 0; | |
593 | Bool_t recalibrated = 0; | |
594 | if (1 && !fClusName.IsNull()) { | |
595 | AliAnalysisTaskEMCALClusterizeFast *cltask = 0; | |
596 | TObjArray *ts = am->GetTasks(); | |
597 | cltask = dynamic_cast<AliAnalysisTaskEMCALClusterizeFast*>(ts->FindObject(fClusName)); | |
598 | if (cltask && cltask->GetClusters()) { | |
d595acbb | 599 | fRecPoints = const_cast<TObjArray*>(cltask->GetClusters()); |
43cfaa06 | 600 | clusattached = cltask->GetAttachClusters(); |
601 | if (cltask->GetCalibData()!=0) | |
602 | recalibrated = kTRUE; | |
603 | } | |
604 | } | |
605 | if (1 && AODEvent() && !fClusName.IsNull()) { | |
717fe7de | 606 | TList *l = AODEvent()->GetList(); |
607 | TClonesArray *clus = 0; | |
608 | if (l) { | |
609 | clus = dynamic_cast<TClonesArray*>(l->FindObject(fClusName)); | |
610 | fAodClusters = clus; | |
611 | } | |
ea3fd2d5 | 612 | } |
717fe7de | 613 | |
614 | if (fEsdEv) { // ESD input mode | |
43cfaa06 | 615 | if (1 && (!fRecPoints||clusattached)) { |
616 | if (!clusattached) | |
617 | am->LoadBranch("CaloClusters"); | |
717fe7de | 618 | TList *l = fEsdEv->GetList(); |
619 | TClonesArray *clus = 0; | |
620 | if (l) { | |
621 | clus = dynamic_cast<TClonesArray*>(l->FindObject("CaloClusters")); | |
622 | fEsdClusters = clus; | |
ea3fd2d5 | 623 | } |
624 | } | |
43cfaa06 | 625 | if (1) { |
626 | if (!recalibrated) | |
627 | am->LoadBranch("EMCALCells."); | |
717fe7de | 628 | fEsdCells = fEsdEv->GetEMCALCells(); |
629 | } | |
630 | } else if (fAodEv) { // AOD input mode | |
43cfaa06 | 631 | if (1 && (!fAodClusters || clusattached)) { |
632 | if (!clusattached) | |
633 | am->LoadBranch("caloClusters"); | |
717fe7de | 634 | TList *l = fAodEv->GetList(); |
635 | TClonesArray *clus = 0; | |
636 | if (l) { | |
637 | clus = dynamic_cast<TClonesArray*>(l->FindObject("caloClusters")); | |
638 | fAodClusters = clus; | |
ea3fd2d5 | 639 | } |
640 | } | |
43cfaa06 | 641 | if (1) { |
642 | if (!recalibrated) | |
643 | am->LoadBranch("emcalCells"); | |
717fe7de | 644 | fAodCells = fAodEv->GetEMCALCells(); |
645 | } | |
646 | } else { | |
647 | AliFatal("Impossible to not have either pointer to ESD or AOD event"); | |
648 | } | |
ea3fd2d5 | 649 | |
43cfaa06 | 650 | if (1) { |
651 | AliDebug(2,Form("fRecPoints set: %p", fRecPoints)); | |
652 | AliDebug(2,Form("fEsdClusters set: %p", fEsdClusters)); | |
653 | AliDebug(2,Form("fEsdCells set: %p", fEsdCells)); | |
654 | AliDebug(2,Form("fAodClusters set: %p", fAodClusters)); | |
655 | AliDebug(2,Form("fAodCells set: %p", fAodCells)); | |
656 | } | |
657 | ||
a49742b5 | 658 | if (fDoAfterburner) |
659 | ClusterAfterburner(); | |
6eb6260e | 660 | |
3a952328 | 661 | CalcCaloTriggers(); |
b6c599fe | 662 | CalcPrimTracks(); |
3a952328 | 663 | CalcTracks(); |
296ea9b4 | 664 | CalcClusterProps(); |
665 | ||
286b47a5 | 666 | FillCellHists(); |
3a952328 | 667 | if (!fTrainMode) { |
668 | FillClusHists(); | |
669 | FillPionHists(); | |
670 | FillOtherHists(); | |
671 | } | |
788ca675 | 672 | FillNtuple(); |
ea3fd2d5 | 673 | |
3a952328 | 674 | if (fTrainMode) { |
675 | fSelTracks->Clear(); | |
676 | fSelPrimTracks->Clear(); | |
677 | fTriggers->Clear(); | |
678 | } | |
de4cee41 | 679 | |
717fe7de | 680 | PostData(1, fOutput); |
ea3fd2d5 | 681 | } |
682 | ||
683 | //________________________________________________________________________ | |
684 | void AliAnalysisTaskEMCALPi0PbPb::Terminate(Option_t *) | |
685 | { | |
717fe7de | 686 | // Terminate called at the end of analysis. |
f5d4ab70 | 687 | |
688 | if (fNtuple) { | |
689 | TFile *f = OpenFile(1); | |
690 | if (f) | |
691 | fNtuple->Write(); | |
692 | } | |
d9f26424 | 693 | |
b3ee6797 | 694 | AliInfo(Form("%s: Accepted %lld events", GetName(), fNEvs)); |
ea3fd2d5 | 695 | } |
ea3fd2d5 | 696 | |
296ea9b4 | 697 | //________________________________________________________________________ |
3a952328 | 698 | void AliAnalysisTaskEMCALPi0PbPb::CalcCaloTriggers() |
296ea9b4 | 699 | { |
3a952328 | 700 | // Calculate triggers |
296ea9b4 | 701 | |
3a952328 | 702 | fTriggers->Clear(); |
296ea9b4 | 703 | |
3a952328 | 704 | AliVCaloCells *cells = fEsdCells; |
705 | if (!cells) | |
706 | cells = fAodCells; | |
707 | if (!cells) | |
708 | return; | |
709 | ||
710 | Int_t ncells = cells->GetNumberOfCells(); | |
711 | if (ncells<=0) | |
712 | return; | |
713 | ||
714 | if (ncells>fNAmpInTrigger) { | |
715 | delete [] fAmpInTrigger; | |
716 | fAmpInTrigger = new Float_t[ncells]; | |
717 | fNAmpInTrigger = ncells; | |
296ea9b4 | 718 | } |
3a952328 | 719 | for (Int_t i=0;i<ncells;++i) |
720 | fAmpInTrigger[i] = 0; | |
296ea9b4 | 721 | |
3a952328 | 722 | std::map<Short_t,Short_t> map; |
723 | for (Short_t pos=0;pos<ncells;++pos) { | |
724 | Short_t id = cells->GetCellNumber(pos); | |
725 | map[id]=pos; | |
726 | } | |
727 | ||
728 | AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager(); | |
729 | am->LoadBranch("EMCALTrigger."); | |
730 | ||
731 | AliESDCaloTrigger *triggers = fEsdEv->GetCaloTrigger("EMCAL"); | |
732 | if (!triggers) | |
733 | return; | |
734 | if (triggers->GetEntries()<=0) | |
b6c599fe | 735 | return; |
736 | ||
3a952328 | 737 | triggers->Reset(); |
738 | Int_t ntrigs=0; | |
739 | while (triggers->Next()) { | |
740 | Int_t gCol=0, gRow=0, ntimes=0; | |
741 | triggers->GetPosition(gCol,gRow); | |
742 | triggers->GetNL0Times(ntimes); | |
743 | if (ntimes<1) | |
744 | continue; | |
745 | Float_t amp=0; | |
746 | triggers->GetAmplitude(amp); | |
747 | Int_t find = -1; | |
748 | fGeom->GetAbsFastORIndexFromPositionInEMCAL(gCol,gRow,find); | |
749 | if (find<0) | |
750 | continue; | |
751 | Int_t cidx[4] = {-1}; | |
752 | Bool_t ret = fGeom->GetCellIndexFromFastORIndex(find, cidx); | |
753 | if (!ret) | |
754 | continue; | |
755 | Int_t trgtimes[25]; | |
756 | triggers->GetL0Times(trgtimes); | |
757 | Int_t mintime = trgtimes[0]; | |
758 | Int_t maxtime = trgtimes[0]; | |
759 | Bool_t trigInTimeWindow = 0; | |
760 | for (Int_t i=0;i<ntimes;++i) { | |
761 | if (trgtimes[i]<mintime) | |
762 | mintime = trgtimes[i]; | |
763 | if (maxtime<trgtimes[i]) | |
764 | maxtime = trgtimes[i]; | |
f2b8fca6 | 765 | if ((fMinL0Time<=trgtimes[i]) && (fMaxL0Time>=trgtimes[i])) |
3a952328 | 766 | trigInTimeWindow = 1; |
767 | } | |
768 | ||
769 | Double_t tenergy = 0; | |
770 | Double_t tphi=0; | |
771 | Double_t teta=0; | |
772 | for (Int_t i=0;i<3;++i) { | |
773 | Short_t pos = -1; | |
774 | std::map<Short_t,Short_t>::iterator it = map.find(cidx[i]); | |
775 | if (it!=map.end()) | |
776 | pos = it->second; | |
777 | if (pos<0) | |
b6c599fe | 778 | continue; |
3a952328 | 779 | if (trigInTimeWindow) |
5fe1ca23 | 780 | fAmpInTrigger[pos] = amp; |
3a952328 | 781 | Float_t eta=-1, phi=-1; |
782 | fGeom->EtaPhiFromIndex(cidx[i],eta,phi); | |
783 | Double_t en= cells->GetAmplitude(pos); | |
784 | tenergy+=en; | |
785 | teta+=eta*en; | |
786 | tphi+=phi*en; | |
787 | } | |
788 | teta/=tenergy; | |
789 | tphi/=tenergy; | |
790 | AliStaTrigger *trignew = static_cast<AliStaTrigger*>(fTriggers->New(ntrigs++)); | |
791 | trignew->fE = tenergy; | |
792 | trignew->fEta = teta; | |
793 | trignew->fPhi = tphi; | |
794 | trignew->fAmp = amp; | |
795 | trignew->fMinTime = mintime; | |
796 | trignew->fMaxTime = maxtime; | |
797 | } | |
798 | } | |
799 | ||
800 | //________________________________________________________________________ | |
801 | void AliAnalysisTaskEMCALPi0PbPb::CalcClusterProps() | |
802 | { | |
803 | // Calculate cluster properties | |
804 | ||
805 | fClusters->Clear(); | |
806 | ||
807 | AliVCaloCells *cells = fEsdCells; | |
808 | if (!cells) | |
809 | cells = fAodCells; | |
810 | if (!cells) | |
811 | return; | |
812 | ||
813 | TObjArray *clusters = fEsdClusters; | |
814 | if (!clusters) | |
815 | clusters = fAodClusters; | |
816 | if (!clusters) | |
817 | return; | |
818 | ||
819 | Int_t ncells = cells->GetNumberOfCells(); | |
820 | Int_t nclus = clusters->GetEntries(); | |
821 | Int_t ntrks = fSelTracks->GetEntries(); | |
822 | Bool_t btracks[6][ntrks]; | |
823 | memset(btracks,0,sizeof(btracks)); | |
824 | ||
825 | std::map<Short_t,Short_t> map; | |
826 | for (Short_t pos=0;pos<ncells;++pos) { | |
827 | Short_t id = cells->GetCellNumber(pos); | |
828 | map[id]=pos; | |
829 | } | |
830 | ||
831 | for(Int_t i=0, ncl=0; i<nclus; ++i) { | |
832 | AliVCluster *clus = static_cast<AliVCluster*>(clusters->At(i)); | |
833 | if (!clus) | |
834 | continue; | |
835 | if (!clus->IsEMCAL()) | |
836 | continue; | |
837 | if (clus->E()<fMinE) | |
838 | continue; | |
839 | ||
840 | Float_t clsPos[3] = {0}; | |
841 | clus->GetPosition(clsPos); | |
842 | TVector3 clsVec(clsPos); | |
843 | Double_t vertex[3] = {0}; | |
844 | InputEvent()->GetPrimaryVertex()->GetXYZ(vertex); | |
845 | TLorentzVector clusterVec; | |
846 | clus->GetMomentum(clusterVec,vertex); | |
847 | Double_t clsEta = clusterVec.Eta(); | |
848 | ||
849 | AliStaCluster *cl = static_cast<AliStaCluster*>(fClusters->New(ncl++)); | |
850 | cl->fE = clus->E(); | |
851 | cl->fR = clsVec.Perp(); | |
852 | cl->fEta = clsVec.Eta(); | |
853 | cl->fPhi = clsVec.Phi(); | |
854 | cl->fN = clus->GetNCells(); | |
855 | cl->fN1 = GetNCells(clus,0.100); | |
856 | cl->fN3 = GetNCells(clus,0.300); | |
857 | Short_t id = -1; | |
858 | Double_t emax = GetMaxCellEnergy(clus, id); | |
859 | cl->fIdMax = id; | |
860 | cl->fEmax = emax; | |
861 | if (clus->GetDistanceToBadChannel()<10000) | |
862 | cl->fDbc = clus->GetDistanceToBadChannel(); | |
863 | if (!TMath::IsNaN(clus->GetDispersion())) | |
864 | cl->fDisp = clus->GetDispersion(); | |
865 | if (!TMath::IsNaN(clus->GetM20())) | |
0fbe8d4f | 866 | cl->fM20 = clus->GetM20(); |
3a952328 | 867 | if (!TMath::IsNaN(clus->GetM02())) |
0fbe8d4f | 868 | cl->fM02 = clus->GetM02(); |
3a952328 | 869 | Double_t maxAxis = 0, minAxis = 0; |
870 | GetSigma(clus,maxAxis,minAxis); | |
871 | clus->SetTOF(maxAxis); // store sigma in TOF | |
872 | cl->fSig = maxAxis; | |
873 | Double_t clusterEcc = 0; | |
874 | if (maxAxis > 0) | |
875 | clusterEcc = TMath::Sqrt(1.0 - minAxis*minAxis/(maxAxis*maxAxis)); | |
876 | clus->SetChi2(clusterEcc); // store ecc in chi2 | |
877 | cl->fEcc = clusterEcc; | |
878 | cl->fTrIso = GetTrackIsolation(clusterVec.Eta(),clusterVec.Phi(),fIsoDist); | |
879 | cl->fTrIso1 = GetTrackIsolation(clusterVec.Eta(),clusterVec.Phi(),fIsoDist, 1); | |
880 | cl->fTrIso2 = GetTrackIsolation(clusterVec.Eta(),clusterVec.Phi(),fIsoDist, 2); | |
881 | cl->fCeCore = GetCellIsolation(clsVec.Eta(),clsVec.Phi(),0.05); | |
882 | cl->fCeIso = GetCellIsolation(clsVec.Eta(),clsVec.Phi(),fIsoDist); | |
0fbe8d4f | 883 | |
3a952328 | 884 | Double_t trigpen = 0; |
885 | Double_t trignen = 0; | |
886 | for(Int_t j=0; j<cl->fN; ++j) { | |
887 | Short_t cid = TMath::Abs(clus->GetCellAbsId(j)); | |
888 | Short_t pos = -1; | |
889 | std::map<Short_t,Short_t>::iterator it = map.find(cid); | |
890 | if (it!=map.end()) | |
891 | pos = it->second; | |
892 | if (pos<0) | |
b6c599fe | 893 | continue; |
3a952328 | 894 | if (fAmpInTrigger[pos]>0) |
895 | trigpen += cells->GetAmplitude(pos); | |
896 | else if (fAmpInTrigger[pos]<0) | |
897 | trignen += cells->GetAmplitude(pos); | |
b6c599fe | 898 | } |
3a952328 | 899 | if (trigpen>0) { |
f3582e89 | 900 | cl->fIsTrigM = 1; |
901 | cl->fTrigE = trigpen; | |
3a952328 | 902 | } |
903 | if (trignen>0) { | |
f3582e89 | 904 | cl->fIsTrigM = 1; |
3a952328 | 905 | cl->fTrigMaskE = trignen; |
906 | } | |
0fbe8d4f | 907 | cl->fIsShared = IsShared(clus); |
3a952328 | 908 | |
909 | // track matching | |
910 | Double_t mind2 = 1e10; | |
911 | for(Int_t j = 0; j<ntrks; ++j) { | |
912 | AliVTrack *track = static_cast<AliVTrack*>(fSelTracks->At(j)); | |
b6c599fe | 913 | if (!track) |
914 | continue; | |
3a952328 | 915 | |
916 | if (TMath::Abs(clsEta-track->Eta())>0.5) | |
b6c599fe | 917 | continue; |
3a952328 | 918 | |
919 | TVector3 vec(clsPos); | |
920 | Int_t index = (Int_t)(vec.Phi()*TMath::RadToDeg()/20); | |
921 | if (btracks[index-4][j]) { | |
b6c599fe | 922 | continue; |
3a952328 | 923 | } |
b6c599fe | 924 | |
3a952328 | 925 | Float_t tmpR=-1, tmpZ=-1; |
926 | if (!fDoTrMatGeom) { | |
927 | AliExternalTrackParam *tParam = 0; | |
928 | if (fEsdEv) { | |
929 | AliESDtrack *esdTrack = static_cast<AliESDtrack*>(track); | |
930 | tParam = new AliExternalTrackParam(*esdTrack->GetTPCInnerParam()); | |
931 | } else | |
932 | tParam = new AliExternalTrackParam(track); | |
933 | ||
934 | Double_t bfield[3] = {0}; | |
935 | track->GetBxByBz(bfield); | |
936 | Double_t alpha = (index+0.5)*20*TMath::DegToRad(); | |
937 | vec.RotateZ(-alpha); //Rotate the cluster to the local extrapolation coordinate system | |
938 | tParam->Rotate(alpha); //Rotate the track to the same local extrapolation system | |
939 | Bool_t ret = tParam->PropagateToBxByBz(vec.X(), bfield); | |
940 | if (!ret) { | |
941 | btracks[index-4][j]=1; | |
942 | delete tParam; | |
943 | continue; | |
b6c599fe | 944 | } |
3a952328 | 945 | Double_t trkPos[3] = {0}; |
946 | tParam->GetXYZ(trkPos); //Get the extrapolated global position | |
947 | tmpR = TMath::Sqrt( TMath::Power(clsPos[0]-trkPos[0],2)+TMath::Power(clsPos[1]-trkPos[1],2)+TMath::Power(clsPos[2]-trkPos[2],2) ); | |
948 | tmpZ = clsPos[2]-trkPos[2]; | |
949 | delete tParam; | |
950 | } else { | |
951 | if (TMath::Abs(clsEta-track->Eta())>fIsoDist) | |
952 | continue; | |
953 | AliExternalTrackParam tParam(track); | |
954 | if (!fReco->ExtrapolateTrackToCluster(&tParam, clus, tmpR, tmpZ)) | |
955 | continue; | |
b6c599fe | 956 | } |
3a952328 | 957 | |
958 | Double_t d2 = tmpR; | |
959 | if (mind2>d2) { | |
960 | mind2=d2; | |
961 | cl->fTrDz = tmpZ; | |
962 | cl->fTrDr = TMath::Sqrt(tmpR*tmpR-tmpZ*tmpZ); | |
963 | cl->fTrEp = clus->E()/track->P(); | |
f3582e89 | 964 | cl->fIsTrackM = 1; |
3a952328 | 965 | } |
966 | } | |
967 | ||
f3582e89 | 968 | if (cl->fIsTrackM) { |
3a952328 | 969 | fHMatchDr->Fill(cl->fTrDr); |
970 | fHMatchDz->Fill(cl->fTrDz); | |
971 | fHMatchEp->Fill(cl->fTrEp); | |
b6c599fe | 972 | } |
973 | } | |
974 | } | |
975 | ||
976 | //________________________________________________________________________ | |
977 | void AliAnalysisTaskEMCALPi0PbPb::CalcPrimTracks() | |
978 | { | |
979 | // Calculate track properties. | |
980 | ||
981 | fSelPrimTracks->Clear(); | |
982 | ||
983 | AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager(); | |
984 | TClonesArray *tracks = 0; | |
985 | if (fEsdEv) { | |
986 | am->LoadBranch("Tracks"); | |
987 | TList *l = fEsdEv->GetList(); | |
988 | tracks = dynamic_cast<TClonesArray*>(l->FindObject("Tracks")); | |
989 | } else { | |
990 | am->LoadBranch("tracks"); | |
991 | TList *l = fAodEv->GetList(); | |
992 | tracks = dynamic_cast<TClonesArray*>(l->FindObject("tracks")); | |
993 | } | |
994 | ||
296ea9b4 | 995 | if (!tracks) |
996 | return; | |
997 | ||
0ec74551 | 998 | AliEMCALEMCGeometry *emc = fGeom->GetEMCGeometry(); |
b6c599fe | 999 | Double_t phimin = emc->GetArm1PhiMin()*TMath::DegToRad()-fIsoDist*1.25; |
1000 | Double_t phimax = emc->GetArm1PhiMax()*TMath::DegToRad()+fIsoDist*1.25; | |
0ec74551 | 1001 | |
296ea9b4 | 1002 | if (fEsdEv) { |
b6c599fe | 1003 | fSelPrimTracks->SetOwner(kTRUE); |
0ec74551 | 1004 | am->LoadBranch("PrimaryVertex."); |
1005 | const AliESDVertex *vtx = fEsdEv->GetPrimaryVertexSPD(); | |
1006 | am->LoadBranch("SPDVertex."); | |
1007 | const AliESDVertex *vtxSPD = fEsdEv->GetPrimaryVertexSPD(); | |
1008 | am->LoadBranch("Tracks"); | |
1009 | const Int_t Ntracks = tracks->GetEntries(); | |
1010 | for (Int_t iTracks = 0; iTracks < Ntracks; ++iTracks) { | |
1011 | AliESDtrack *track = static_cast<AliESDtrack*>(tracks->At(iTracks)); | |
1012 | if (!track) { | |
1013 | AliWarning(Form("Could not receive track %d\n", iTracks)); | |
1014 | continue; | |
1015 | } | |
1016 | if (fTrCuts && !fTrCuts->IsSelected(track)) | |
296ea9b4 | 1017 | continue; |
1018 | Double_t eta = track->Eta(); | |
1019 | if (eta<-1||eta>1) | |
1020 | continue; | |
0ec74551 | 1021 | if (track->Phi()<phimin||track->Phi()>phimax) |
1022 | continue; | |
2e4d8148 | 1023 | |
0ec74551 | 1024 | AliESDtrack copyt(*track); |
1025 | Double_t bfield[3]; | |
1026 | copyt.GetBxByBz(bfield); | |
1027 | AliExternalTrackParam tParam; | |
1028 | Bool_t relate = copyt.RelateToVertexBxByBz(vtxSPD,bfield,kVeryBig,&tParam); | |
1029 | if (!relate) | |
1030 | continue; | |
1031 | copyt.Set(tParam.GetX(),tParam.GetAlpha(),tParam.GetParameter(),tParam.GetCovariance()); | |
2e4d8148 | 1032 | |
0ec74551 | 1033 | Double_t p[3] = { 0. }; |
1034 | copyt.GetPxPyPz(p); | |
1035 | Double_t pos[3] = { 0. }; | |
1036 | copyt.GetXYZ(pos); | |
1037 | Double_t covTr[21] = { 0. }; | |
1038 | copyt.GetCovarianceXYZPxPyPz(covTr); | |
1039 | Double_t pid[10] = { 0. }; | |
1040 | copyt.GetESDpid(pid); | |
1041 | AliAODTrack *aTrack = new AliAODTrack(copyt.GetID(), | |
1042 | copyt.GetLabel(), | |
1043 | p, | |
1044 | kTRUE, | |
1045 | pos, | |
1046 | kFALSE, | |
1047 | covTr, | |
1048 | (Short_t)copyt.GetSign(), | |
1049 | copyt.GetITSClusterMap(), | |
1050 | pid, | |
1051 | 0,/*fPrimaryVertex,*/ | |
1052 | kTRUE, // check if this is right | |
1053 | vtx->UsesTrack(copyt.GetID())); | |
1054 | aTrack->SetTPCClusterMap(copyt.GetTPCClusterMap()); | |
1055 | aTrack->SetTPCSharedMap (copyt.GetTPCSharedMap()); | |
1056 | Float_t ndf = copyt.GetTPCNcls() + 1 - 5 ; | |
1057 | if(ndf>0) | |
1058 | aTrack->SetChi2perNDF(copyt.GetTPCchi2()/ndf); | |
1059 | else | |
1060 | aTrack->SetChi2perNDF(-1); | |
1061 | aTrack->SetFlags(copyt.GetStatus()); | |
1062 | aTrack->SetTPCPointsF(copyt.GetTPCNclsF()); | |
b6c599fe | 1063 | fSelPrimTracks->Add(aTrack); |
296ea9b4 | 1064 | } |
1065 | } else { | |
1066 | Int_t ntracks = tracks->GetEntries(); | |
1067 | for (Int_t i=0; i<ntracks; ++i) { | |
1068 | AliAODTrack *track = static_cast<AliAODTrack*>(tracks->At(i)); | |
1069 | if (!track) | |
1070 | continue; | |
1071 | Double_t eta = track->Eta(); | |
1072 | if (eta<-1||eta>1) | |
1073 | continue; | |
0ec74551 | 1074 | if (track->Phi()<phimin||track->Phi()>phimax) |
1075 | continue; | |
b6c599fe | 1076 | if(track->GetTPCNcls()<fMinNClusPerTr) |
296ea9b4 | 1077 | continue; |
b6c599fe | 1078 | //todo: Learn how to set/filter AODs for prim/sec tracks |
1079 | fSelPrimTracks->Add(track); | |
296ea9b4 | 1080 | } |
1081 | } | |
1082 | } | |
1083 | ||
1084 | //________________________________________________________________________ | |
3a952328 | 1085 | void AliAnalysisTaskEMCALPi0PbPb::CalcTracks() |
296ea9b4 | 1086 | { |
3a952328 | 1087 | // Calculate track properties (including secondaries). |
b6c599fe | 1088 | |
3a952328 | 1089 | fSelTracks->Clear(); |
296ea9b4 | 1090 | |
3a952328 | 1091 | AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager(); |
1092 | TClonesArray *tracks = 0; | |
1093 | if (fEsdEv) { | |
1094 | am->LoadBranch("Tracks"); | |
1095 | TList *l = fEsdEv->GetList(); | |
1096 | tracks = dynamic_cast<TClonesArray*>(l->FindObject("Tracks")); | |
1097 | } else { | |
1098 | am->LoadBranch("tracks"); | |
1099 | TList *l = fAodEv->GetList(); | |
1100 | tracks = dynamic_cast<TClonesArray*>(l->FindObject("tracks")); | |
1101 | } | |
296ea9b4 | 1102 | |
3a952328 | 1103 | if (!tracks) |
1104 | return; | |
296ea9b4 | 1105 | |
3a952328 | 1106 | if (fEsdEv) { |
1107 | const Int_t Ntracks = tracks->GetEntries(); | |
1108 | for (Int_t iTracks = 0; iTracks < Ntracks; ++iTracks) { | |
1109 | AliESDtrack *track = static_cast<AliESDtrack*>(tracks->At(iTracks)); | |
1110 | if (!track) { | |
1111 | AliWarning(Form("Could not receive track %d\n", iTracks)); | |
1112 | continue; | |
1113 | } | |
1114 | if (fTrCuts && !fTrCuts->IsSelected(track)) | |
1115 | continue; | |
1116 | Double_t eta = track->Eta(); | |
1117 | if (eta<-1||eta>1) | |
1118 | continue; | |
1119 | fSelTracks->Add(track); | |
1120 | } | |
1121 | } else { | |
1122 | Int_t ntracks = tracks->GetEntries(); | |
1123 | for (Int_t i=0; i<ntracks; ++i) { | |
1124 | AliAODTrack *track = static_cast<AliAODTrack*>(tracks->At(i)); | |
296ea9b4 | 1125 | if (!track) |
1126 | continue; | |
3a952328 | 1127 | Double_t eta = track->Eta(); |
1128 | if (eta<-1||eta>1) | |
c4236a58 | 1129 | continue; |
3a952328 | 1130 | if(track->GetTPCNcls()<fMinNClusPerTr) |
b6c599fe | 1131 | continue; |
2e4d8148 | 1132 | |
3a952328 | 1133 | if (0 && (track->Pt()>=0.6) && (track->PxAtDCA()==-999)) { // compute position on EMCAL |
2e4d8148 | 1134 | AliExternalTrackParam tParam(track); |
3a952328 | 1135 | if (AliTrackerBase::PropagateTrackToBxByBz(&tParam, 438, 0.139, 1, kTRUE)) { |
1136 | Double_t trkPos[3]; | |
1137 | tParam.GetXYZ(trkPos); | |
1138 | track->SetPxPyPzAtDCA(trkPos[0],trkPos[1],trkPos[2]); | |
1139 | } | |
b6c599fe | 1140 | } |
3a952328 | 1141 | fSelTracks->Add(track); |
c4236a58 | 1142 | } |
296ea9b4 | 1143 | } |
1144 | } | |
1145 | ||
323834f0 | 1146 | //________________________________________________________________________ |
3a952328 | 1147 | void AliAnalysisTaskEMCALPi0PbPb::ClusterAfterburner() |
323834f0 | 1148 | { |
296ea9b4 | 1149 | // Run custer reconstruction afterburner. |
323834f0 | 1150 | |
1151 | AliVCaloCells *cells = fEsdCells; | |
1152 | if (!cells) | |
1153 | cells = fAodCells; | |
1154 | ||
1155 | if (!cells) | |
1156 | return; | |
1157 | ||
1158 | Int_t ncells = cells->GetNumberOfCells(); | |
1159 | if (ncells<=0) | |
1160 | return; | |
1161 | ||
1162 | Double_t cellMeanE = 0, cellSigE = 0; | |
1163 | for (Int_t i = 0; i<ncells; ++i) { | |
1164 | Double_t cellE = cells->GetAmplitude(i); | |
1165 | cellMeanE += cellE; | |
1166 | cellSigE += cellE*cellE; | |
1167 | } | |
1168 | cellMeanE /= ncells; | |
1169 | cellSigE /= ncells; | |
1170 | cellSigE -= cellMeanE*cellMeanE; | |
1171 | if (cellSigE<0) | |
1172 | cellSigE = 0; | |
1173 | cellSigE = TMath::Sqrt(cellSigE / ncells); | |
1174 | ||
1175 | Double_t subE = cellMeanE - 7*cellSigE; | |
1176 | if (subE<0) | |
1177 | return; | |
1178 | ||
1179 | for (Int_t i = 0; i<ncells; ++i) { | |
1180 | Short_t id=-1; | |
1181 | Double_t amp=0,time=0; | |
1182 | if (!cells->GetCell(i, id, amp, time)) | |
1183 | continue; | |
1184 | amp -= cellMeanE; | |
1185 | if (amp<0.001) | |
1186 | amp = 0; | |
1187 | cells->SetCell(i, id, amp, time); | |
1188 | } | |
1189 | ||
1190 | TObjArray *clusters = fEsdClusters; | |
1191 | if (!clusters) | |
1192 | clusters = fAodClusters; | |
323834f0 | 1193 | if (!clusters) |
1194 | return; | |
1195 | ||
1196 | Int_t nclus = clusters->GetEntries(); | |
1197 | for (Int_t i = 0; i<nclus; ++i) { | |
1198 | AliVCluster *clus = static_cast<AliVCluster*>(clusters->At(i)); | |
1199 | if (!clus->IsEMCAL()) | |
1200 | continue; | |
1201 | Int_t nc = clus->GetNCells(); | |
1202 | Double_t clusE = 0; | |
1203 | UShort_t ids[100] = {0}; | |
1204 | Double_t fra[100] = {0}; | |
1205 | for (Int_t j = 0; j<nc; ++j) { | |
1206 | Short_t id = TMath::Abs(clus->GetCellAbsId(j)); | |
1207 | Double_t cen = cells->GetCellAmplitude(id); | |
1208 | clusE += cen; | |
1209 | if (cen>0) { | |
1210 | ids[nc] = id; | |
1211 | ++nc; | |
1212 | } | |
1213 | } | |
1214 | if (clusE<=0) { | |
1215 | clusters->RemoveAt(i); | |
1216 | continue; | |
1217 | } | |
1218 | ||
1219 | for (Int_t j = 0; j<nc; ++j) { | |
1220 | Short_t id = ids[j]; | |
1221 | Double_t cen = cells->GetCellAmplitude(id); | |
1222 | fra[j] = cen/clusE; | |
1223 | } | |
1224 | clus->SetE(clusE); | |
1225 | AliAODCaloCluster *aodclus = dynamic_cast<AliAODCaloCluster*>(clus); | |
1226 | if (aodclus) { | |
1227 | aodclus->Clear(""); | |
1228 | aodclus->SetNCells(nc); | |
1229 | aodclus->SetCellsAmplitudeFraction(fra); | |
1230 | aodclus->SetCellsAbsId(ids); | |
1231 | } | |
1232 | } | |
1233 | clusters->Compress(); | |
1234 | } | |
1235 | ||
286b47a5 | 1236 | //________________________________________________________________________ |
1237 | void AliAnalysisTaskEMCALPi0PbPb::FillCellHists() | |
1238 | { | |
1239 | // Fill histograms related to cell properties. | |
d595acbb | 1240 | |
90d5b88b | 1241 | AliVCaloCells *cells = fEsdCells; |
1242 | if (!cells) | |
1243 | cells = fAodCells; | |
1244 | ||
296ea9b4 | 1245 | if (!cells) |
1246 | return; | |
90d5b88b | 1247 | |
296ea9b4 | 1248 | Int_t cellModCount[12] = {0}; |
1249 | Double_t cellMaxE = 0; | |
1250 | Double_t cellMeanE = 0; | |
1251 | Int_t ncells = cells->GetNumberOfCells(); | |
1252 | if (ncells==0) | |
1253 | return; | |
1254 | ||
1255 | Int_t nsm = fGeom->GetEMCGeometry()->GetNumberOfSuperModules(); | |
1256 | ||
1257 | for (Int_t i = 0; i<ncells; ++i) { | |
1258 | Int_t absID = TMath::Abs(cells->GetCellNumber(i)); | |
1259 | Double_t cellE = cells->GetAmplitude(i); | |
1260 | fHCellE->Fill(cellE); | |
1261 | if (cellE>cellMaxE) | |
1262 | cellMaxE = cellE; | |
1263 | cellMeanE += cellE; | |
1264 | ||
1265 | Int_t iSM=-1, iTower=-1, nIphi=-1, nIeta=-1; | |
1266 | Bool_t ret = fGeom->GetCellIndex(absID, iSM, iTower, nIphi, nIeta); | |
1267 | if (!ret) { | |
1268 | AliError(Form("Could not get cell index for %d", absID)); | |
1269 | continue; | |
1270 | } | |
1271 | ++cellModCount[iSM]; | |
1272 | Int_t iPhi=-1, iEta=-1; | |
1273 | fGeom->GetCellPhiEtaIndexInSModule(iSM, iTower, nIphi, nIeta, iPhi, iEta); | |
1274 | fHColuRow[iSM]->Fill(iEta,iPhi,1); | |
1275 | fHColuRowE[iSM]->Fill(iEta,iPhi,cellE); | |
1276 | fHCellFreqNoCut[iSM]->Fill(absID); | |
2e4d8148 | 1277 | if (cellE > 0.1) fHCellFreqCut100M[iSM]->Fill(absID); |
1278 | if (cellE > 0.3) fHCellFreqCut300M[iSM]->Fill(absID); | |
296ea9b4 | 1279 | if (fHCellCheckE && fHCellCheckE[absID]) |
1280 | fHCellCheckE[absID]->Fill(cellE); | |
2e4d8148 | 1281 | fHCellFreqE[iSM]->Fill(absID, cellE); |
296ea9b4 | 1282 | } |
1283 | fHCellH->Fill(cellMaxE); | |
1284 | cellMeanE /= ncells; | |
1285 | fHCellM->Fill(cellMeanE); | |
1286 | fHCellM2->Fill(cellMeanE*ncells/24/48/nsm); | |
1287 | for (Int_t i=0; i<nsm; ++i) | |
1288 | fHCellMult[i]->Fill(cellModCount[i]); | |
286b47a5 | 1289 | } |
1290 | ||
1291 | //________________________________________________________________________ | |
1292 | void AliAnalysisTaskEMCALPi0PbPb::FillClusHists() | |
1293 | { | |
90d5b88b | 1294 | // Fill histograms related to cluster properties. |
fa443410 | 1295 | |
90d5b88b | 1296 | TObjArray *clusters = fEsdClusters; |
1297 | if (!clusters) | |
1298 | clusters = fAodClusters; | |
296ea9b4 | 1299 | if (!clusters) |
1300 | return; | |
90d5b88b | 1301 | |
296ea9b4 | 1302 | Double_t vertex[3] = {0}; |
1303 | InputEvent()->GetPrimaryVertex()->GetXYZ(vertex); | |
1304 | ||
1305 | Int_t nclus = clusters->GetEntries(); | |
1306 | for(Int_t i = 0; i<nclus; ++i) { | |
1307 | AliVCluster *clus = static_cast<AliVCluster*>(clusters->At(i)); | |
1308 | if (!clus) | |
1309 | continue; | |
1310 | if (!clus->IsEMCAL()) | |
1311 | continue; | |
90d5b88b | 1312 | TLorentzVector clusterVec; |
296ea9b4 | 1313 | clus->GetMomentum(clusterVec,vertex); |
3a952328 | 1314 | Double_t maxAxis = clus->GetTOF(); //sigma |
1315 | Double_t clusterEcc = clus->Chi2(); //eccentricity | |
296ea9b4 | 1316 | fHClustEccentricity->Fill(clusterEcc); |
1317 | fHClustEtaPhi->Fill(clusterVec.Eta(),clusterVec.Phi()); | |
1318 | fHClustEnergyPt->Fill(clusterVec.E(),clusterVec.Pt()); | |
1319 | fHClustEnergySigma->Fill(clus->E()*maxAxis,clus->E()); | |
1320 | fHClustSigmaSigma->Fill(max(clus->GetM02(),clus->GetM20()),clus->E()*maxAxis); | |
1321 | fHClustNCellEnergyRatio->Fill(clus->GetNCells(),GetMaxCellEnergy(clus)/clus->E()); | |
788ca675 | 1322 | } |
1323 | } | |
b3ee6797 | 1324 | |
788ca675 | 1325 | //________________________________________________________________________ |
1326 | void AliAnalysisTaskEMCALPi0PbPb::FillNtuple() | |
1327 | { | |
1328 | // Fill ntuple. | |
1329 | ||
1330 | if (!fNtuple) | |
1331 | return; | |
1332 | ||
1333 | AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager(); | |
1334 | if (fAodEv) { | |
1335 | fHeader->fRun = fAodEv->GetRunNumber(); | |
1336 | fHeader->fOrbit = fAodEv->GetHeader()->GetOrbitNumber(); | |
1337 | fHeader->fPeriod = fAodEv->GetHeader()->GetPeriodNumber(); | |
1338 | fHeader->fBx = fAodEv->GetHeader()->GetBunchCrossNumber(); | |
1339 | fHeader->fL0 = fAodEv->GetHeader()->GetL0TriggerInputs(); | |
1340 | fHeader->fL1 = fAodEv->GetHeader()->GetL1TriggerInputs(); | |
1341 | fHeader->fL2 = fAodEv->GetHeader()->GetL2TriggerInputs(); | |
1342 | fHeader->fTrClassMask = fAodEv->GetHeader()->GetTriggerMask(); | |
1343 | fHeader->fTrCluster = fAodEv->GetHeader()->GetTriggerCluster(); | |
1344 | fHeader->fOffTriggers = fAodEv->GetHeader()->GetOfflineTrigger(); | |
1345 | fHeader->fFiredTriggers = fAodEv->GetHeader()->GetFiredTriggerClasses(); | |
1346 | } else { | |
1347 | fHeader->fRun = fEsdEv->GetRunNumber(); | |
1348 | fHeader->fOrbit = fEsdEv->GetHeader()->GetOrbitNumber(); | |
1349 | fHeader->fPeriod = fEsdEv->GetESDRun()->GetPeriodNumber(); | |
1350 | fHeader->fBx = fEsdEv->GetHeader()->GetBunchCrossNumber(); | |
1351 | fHeader->fL0 = fEsdEv->GetHeader()->GetL0TriggerInputs(); | |
1352 | fHeader->fL1 = fEsdEv->GetHeader()->GetL1TriggerInputs(); | |
1353 | fHeader->fL2 = fEsdEv->GetHeader()->GetL2TriggerInputs(); | |
1354 | fHeader->fTrClassMask = fEsdEv->GetHeader()->GetTriggerMask(); | |
1355 | fHeader->fTrCluster = fEsdEv->GetHeader()->GetTriggerCluster(); | |
1356 | fHeader->fOffTriggers = ((AliInputEventHandler*)(am->GetInputEventHandler()))->IsEventSelected(); | |
1357 | fHeader->fFiredTriggers = fEsdEv->GetFiredTriggerClasses(); | |
5fe1ca23 | 1358 | Float_t v0CorrR = 0; |
1359 | fHeader->fV0 = AliESDUtils::GetCorrV0(fEsdEv,v0CorrR); | |
1360 | const AliMultiplicity *mult = fEsdEv->GetMultiplicity(); | |
1361 | if (mult) | |
1362 | fHeader->fCl1 = mult->GetNumberOfITSClusters(1); | |
1363 | fHeader->fTr = AliESDtrackCuts::GetReferenceMultiplicity(fEsdEv,1); | |
788ca675 | 1364 | } |
1365 | AliCentrality *cent = InputEvent()->GetCentrality(); | |
1366 | fHeader->fV0Cent = cent->GetCentralityPercentileUnchecked("V0M"); | |
1367 | fHeader->fCl1Cent = cent->GetCentralityPercentileUnchecked("CL1"); | |
1368 | fHeader->fTrCent = cent->GetCentralityPercentileUnchecked("TRK"); | |
1369 | fHeader->fCqual = cent->GetQuality(); | |
b6c599fe | 1370 | |
1371 | AliEventplane *ep = InputEvent()->GetEventplane(); | |
1372 | if (ep) { | |
1373 | if (ep->GetQVector()) | |
1374 | fHeader->fPsi = ep->GetQVector()->Phi()/2. ; | |
f2b8fca6 | 1375 | else |
1376 | fHeader->fPsi = -1; | |
b6c599fe | 1377 | if (ep->GetQsub1()&&ep->GetQsub2()) |
1378 | fHeader->fPsiRes = ep->GetQsub1()->Phi()/2.-ep->GetQsub2()->Phi()/2.; | |
1379 | else | |
1380 | fHeader->fPsiRes = 0; | |
1381 | } | |
1382 | ||
788ca675 | 1383 | Double_t val = 0; |
1384 | TString trgclasses(fHeader->fFiredTriggers); | |
1385 | for (Int_t j = 0; j<fTrClassNamesArr->GetEntries(); ++j) { | |
1386 | const char *name = fTrClassNamesArr->At(j)->GetName(); | |
1387 | if (trgclasses.Contains(name)) | |
1388 | val += TMath::Power(2,j); | |
1389 | } | |
1390 | fHeader->fTcls = (UInt_t)val; | |
1391 | ||
5fe1ca23 | 1392 | fHeader->fNSelTr = fSelTracks->GetEntries(); |
1393 | fHeader->fNSelPrimTr = fSelPrimTracks->GetEntries(); | |
1394 | ||
1395 | fHeader->fNCells = 0; | |
1396 | fHeader->fNCells1 = 0; | |
1397 | fHeader->fNCells2 = 0; | |
1398 | fHeader->fNCells5 = 0; | |
1399 | fHeader->fMaxCellE = 0; | |
1400 | ||
1401 | AliVCaloCells *cells = fEsdCells; | |
1402 | if (!cells) | |
1403 | cells = fAodCells; | |
1404 | ||
1405 | if (cells) { | |
1406 | Int_t ncells = cells->GetNumberOfCells(); | |
1407 | for(Int_t j=0; j<ncells; ++j) { | |
1408 | Double_t cellen = cells->GetAmplitude(j); | |
1409 | if (cellen>1) | |
1410 | ++fHeader->fNCells1; | |
1411 | if (cellen>2) | |
1412 | ++fHeader->fNCells2; | |
1413 | if (cellen>5) | |
1414 | ++fHeader->fNCells5; | |
1415 | if (cellen>fHeader->fMaxCellE) | |
1416 | fHeader->fMaxCellE = cellen; | |
1417 | } | |
1418 | fHeader->fNCells = ncells; | |
1419 | } | |
1420 | ||
1421 | fHeader->fNClus = 0; | |
1422 | fHeader->fNClus1 = 0; | |
1423 | fHeader->fNClus2 = 0; | |
1424 | fHeader->fNClus5 = 0; | |
1425 | fHeader->fMaxClusE = 0; | |
1426 | ||
1427 | TObjArray *clusters = fEsdClusters; | |
1428 | if (!clusters) | |
1429 | clusters = fAodClusters; | |
1430 | ||
1431 | if (clusters) { | |
1432 | Int_t nclus = clusters->GetEntries(); | |
1433 | for(Int_t j=0; j<nclus; ++j) { | |
1434 | AliVCluster *clus = static_cast<AliVCluster*>(clusters->At(j)); | |
1435 | Double_t clusen = clus->E(); | |
1436 | if (clusen>1) | |
1437 | ++fHeader->fNClus1; | |
1438 | if (clusen>2) | |
1439 | ++fHeader->fNClus2; | |
1440 | if (clusen>5) | |
1441 | ++fHeader->fNClus5; | |
1442 | if (clusen>fHeader->fMaxClusE) | |
1443 | fHeader->fMaxClusE = clusen; | |
1444 | } | |
1445 | fHeader->fNClus = nclus; | |
1446 | } | |
1447 | ||
1448 | ||
788ca675 | 1449 | if (fAodEv) { |
1450 | am->LoadBranch("vertices"); | |
1451 | AliAODVertex *pv = fAodEv->GetPrimaryVertex(); | |
1452 | FillVertex(fPrimVert, pv); | |
1453 | AliAODVertex *sv = fAodEv->GetPrimaryVertexSPD(); | |
1454 | FillVertex(fSpdVert, sv); | |
1455 | } else { | |
1456 | am->LoadBranch("PrimaryVertex."); | |
1457 | const AliESDVertex *pv = fEsdEv->GetPrimaryVertexTracks(); | |
1458 | FillVertex(fPrimVert, pv); | |
1459 | am->LoadBranch("SPDVertex."); | |
1460 | const AliESDVertex *sv = fEsdEv->GetPrimaryVertexSPD(); | |
1461 | FillVertex(fSpdVert, sv); | |
1462 | am->LoadBranch("TPCVertex."); | |
1463 | const AliESDVertex *tv = fEsdEv->GetPrimaryVertexTPC(); | |
1464 | FillVertex(fTpcVert, tv); | |
fa443410 | 1465 | } |
788ca675 | 1466 | |
788ca675 | 1467 | fNtuple->Fill(); |
286b47a5 | 1468 | } |
1469 | ||
1470 | //________________________________________________________________________ | |
1471 | void AliAnalysisTaskEMCALPi0PbPb::FillPionHists() | |
1472 | { | |
1473 | // Fill histograms related to pions. | |
286b47a5 | 1474 | |
296ea9b4 | 1475 | Double_t vertex[3] = {0}; |
fa443410 | 1476 | InputEvent()->GetPrimaryVertex()->GetXYZ(vertex); |
1477 | ||
90d5b88b | 1478 | TObjArray *clusters = fEsdClusters; |
1479 | if (!clusters) | |
1480 | clusters = fAodClusters; | |
fa443410 | 1481 | |
90d5b88b | 1482 | if (clusters) { |
1483 | TLorentzVector clusterVec1; | |
1484 | TLorentzVector clusterVec2; | |
1485 | TLorentzVector pionVec; | |
1486 | ||
1487 | Int_t nclus = clusters->GetEntries(); | |
fa443410 | 1488 | for (Int_t i = 0; i<nclus; ++i) { |
90d5b88b | 1489 | AliVCluster *clus1 = static_cast<AliVCluster*>(clusters->At(i)); |
fa443410 | 1490 | if (!clus1) |
1491 | continue; | |
1492 | if (!clus1->IsEMCAL()) | |
1493 | continue; | |
3c661da5 | 1494 | if (clus1->E()<fMinE) |
6eb6260e | 1495 | continue; |
a49742b5 | 1496 | if (clus1->GetNCells()<fNminCells) |
1497 | continue; | |
f224d35b | 1498 | if (GetMaxCellEnergy(clus1)/clus1->E()<fMinErat) |
1499 | continue; | |
3c661da5 | 1500 | if (clus1->Chi2()<fMinEcc) // eccentricity cut |
f224d35b | 1501 | continue; |
fa443410 | 1502 | clus1->GetMomentum(clusterVec1,vertex); |
1503 | for (Int_t j = i+1; j<nclus; ++j) { | |
90d5b88b | 1504 | AliVCluster *clus2 = static_cast<AliVCluster*>(clusters->At(j)); |
fa443410 | 1505 | if (!clus2) |
1506 | continue; | |
1507 | if (!clus2->IsEMCAL()) | |
1508 | continue; | |
3c661da5 | 1509 | if (clus2->E()<fMinE) |
6eb6260e | 1510 | continue; |
a49742b5 | 1511 | if (clus2->GetNCells()<fNminCells) |
1512 | continue; | |
f224d35b | 1513 | if (GetMaxCellEnergy(clus2)/clus2->E()<fMinErat) |
1514 | continue; | |
3c661da5 | 1515 | if (clus2->Chi2()<fMinEcc) // eccentricity cut |
f224d35b | 1516 | continue; |
fa443410 | 1517 | clus2->GetMomentum(clusterVec2,vertex); |
1518 | pionVec = clusterVec1 + clusterVec2; | |
1519 | Double_t pionZgg = TMath::Abs(clusterVec1.E()-clusterVec2.E())/pionVec.E(); | |
6eb6260e | 1520 | Double_t pionOpeningAngle = clusterVec1.Angle(clusterVec2.Vect()); |
d595acbb | 1521 | if (pionZgg < fAsymMax) { |
1522 | fHPionEtaPhi->Fill(pionVec.Eta(),pionVec.Phi()); | |
1523 | fHPionMggPt->Fill(pionVec.M(),pionVec.Pt()); | |
1524 | fHPionMggAsym->Fill(pionVec.M(),pionZgg); | |
6eb6260e | 1525 | fHPionMggDgg->Fill(pionVec.M(),pionOpeningAngle); |
d595acbb | 1526 | Int_t bin = fPtRanges->FindBin(pionVec.Pt()); |
1527 | fHPionInvMasses[bin]->Fill(pionVec.M()); | |
1528 | } | |
fa443410 | 1529 | } |
1530 | } | |
90d5b88b | 1531 | } |
fa443410 | 1532 | } |
d595acbb | 1533 | |
6eb6260e | 1534 | //________________________________________________________________________ |
323834f0 | 1535 | void AliAnalysisTaskEMCALPi0PbPb::FillOtherHists() |
6eb6260e | 1536 | { |
788ca675 | 1537 | // Fill other histograms. |
1538 | } | |
1539 | ||
1540 | //__________________________________________________________________________________________________ | |
1541 | void AliAnalysisTaskEMCALPi0PbPb::FillVertex(AliStaVertex *v, const AliESDVertex *esdv) | |
1542 | { | |
1543 | // Fill vertex from ESD vertex info. | |
1544 | ||
1545 | v->fVx = esdv->GetX(); | |
1546 | v->fVy = esdv->GetY(); | |
1547 | v->fVz = esdv->GetZ(); | |
1548 | v->fVc = esdv->GetNContributors(); | |
1549 | v->fDisp = esdv->GetDispersion(); | |
1550 | v->fZres = esdv->GetZRes(); | |
1551 | v->fChi2 = esdv->GetChi2(); | |
1552 | v->fSt = esdv->GetStatus(); | |
1553 | v->fIs3D = esdv->IsFromVertexer3D(); | |
1554 | v->fIsZ = esdv->IsFromVertexerZ(); | |
1555 | } | |
1556 | ||
1557 | //__________________________________________________________________________________________________ | |
1558 | void AliAnalysisTaskEMCALPi0PbPb::FillVertex(AliStaVertex *v, const AliAODVertex *aodv) | |
1559 | { | |
1560 | // Fill vertex from AOD vertex info. | |
1561 | ||
1562 | v->fVx = aodv->GetX(); | |
1563 | v->fVy = aodv->GetY(); | |
1564 | v->fVz = aodv->GetZ(); | |
1565 | v->fVc = aodv->GetNContributors(); | |
1566 | v->fChi2 = aodv->GetChi2(); | |
6eb6260e | 1567 | } |
1568 | ||
d595acbb | 1569 | //________________________________________________________________________ |
296ea9b4 | 1570 | Double_t AliAnalysisTaskEMCALPi0PbPb::GetCellIsolation(Double_t cEta, Double_t cPhi, Double_t radius) const |
1571 | { | |
1572 | // Compute isolation based on cell content. | |
1573 | ||
1574 | AliVCaloCells *cells = fEsdCells; | |
1575 | if (!cells) | |
1576 | cells = fAodCells; | |
1577 | if (!cells) | |
1578 | return 0; | |
1579 | ||
1580 | Double_t cellIsolation = 0; | |
1581 | Double_t rad2 = radius*radius; | |
1582 | Int_t ncells = cells->GetNumberOfCells(); | |
1583 | for (Int_t i = 0; i<ncells; ++i) { | |
1584 | Int_t absID = TMath::Abs(cells->GetCellNumber(i)); | |
0fbe8d4f | 1585 | Float_t eta=-1, phi=-1; |
296ea9b4 | 1586 | fGeom->EtaPhiFromIndex(absID,eta,phi); |
0fbe8d4f | 1587 | Double_t phidiff = TVector2::Phi_mpi_pi(phi-cPhi); |
296ea9b4 | 1588 | Double_t dist = (eta-cEta)*(eta-cEta)+phidiff*phidiff; |
1589 | if(dist>rad2) | |
1590 | continue; | |
0fbe8d4f | 1591 | Double_t cellE = cells->GetAmplitude(i); |
296ea9b4 | 1592 | cellIsolation += cellE; |
1593 | } | |
1594 | return cellIsolation; | |
1595 | } | |
1596 | ||
1597 | //________________________________________________________________________ | |
0fbe8d4f | 1598 | Double_t AliAnalysisTaskEMCALPi0PbPb::GetCellEnergy(const AliVCluster *cluster) const |
1599 | { | |
1600 | // Get maximum energy of attached cell. | |
1601 | ||
1602 | Double_t ret = 0; | |
1603 | Int_t ncells = cluster->GetNCells(); | |
1604 | if (fEsdCells) { | |
1605 | for (Int_t i=0; i<ncells; i++) { | |
1606 | Double_t e = fEsdCells->GetCellAmplitude(TMath::Abs(cluster->GetCellAbsId(i))); | |
1607 | ret += e; | |
1608 | } | |
1609 | } else { | |
1610 | for (Int_t i=0; i<ncells; i++) { | |
1611 | Double_t e = fAodCells->GetCellAmplitude(TMath::Abs(cluster->GetCellAbsId(i))); | |
1612 | ret += e; | |
1613 | } | |
1614 | } | |
1615 | return ret; | |
1616 | } | |
1617 | ||
1618 | //________________________________________________________________________ | |
1619 | Double_t AliAnalysisTaskEMCALPi0PbPb::GetMaxCellEnergy(const AliVCluster *cluster, Short_t &id) const | |
d595acbb | 1620 | { |
90d5b88b | 1621 | // Get maximum energy of attached cell. |
1622 | ||
788ca675 | 1623 | id = -1; |
90d5b88b | 1624 | Double_t maxe = 0; |
90d5b88b | 1625 | Int_t ncells = cluster->GetNCells(); |
1626 | if (fEsdCells) { | |
1627 | for (Int_t i=0; i<ncells; i++) { | |
27422847 | 1628 | Double_t e = fEsdCells->GetCellAmplitude(TMath::Abs(cluster->GetCellAbsId(i))); |
f0897c18 | 1629 | if (e>maxe) { |
90d5b88b | 1630 | maxe = e; |
788ca675 | 1631 | id = cluster->GetCellAbsId(i); |
f0897c18 | 1632 | } |
90d5b88b | 1633 | } |
1634 | } else { | |
1635 | for (Int_t i=0; i<ncells; i++) { | |
27422847 | 1636 | Double_t e = fAodCells->GetCellAmplitude(TMath::Abs(cluster->GetCellAbsId(i))); |
90d5b88b | 1637 | if (e>maxe) |
1638 | maxe = e; | |
788ca675 | 1639 | id = cluster->GetCellAbsId(i); |
90d5b88b | 1640 | } |
1641 | } | |
6eb6260e | 1642 | return maxe; |
d595acbb | 1643 | } |
1644 | ||
1645 | //________________________________________________________________________ | |
0fbe8d4f | 1646 | void AliAnalysisTaskEMCALPi0PbPb::GetSigma(const AliVCluster *c, Double_t& sigmaMax, Double_t &sigmaMin) const |
d595acbb | 1647 | { |
90d5b88b | 1648 | // Calculate the (E) weighted variance along the longer (eigen) axis. |
1649 | ||
6bf90832 | 1650 | sigmaMax = 0; // cluster variance along its longer axis |
1651 | sigmaMin = 0; // cluster variance along its shorter axis | |
1652 | Double_t Ec = c->E(); // cluster energy | |
296ea9b4 | 1653 | if(Ec<=0) |
1654 | return; | |
6bf90832 | 1655 | Double_t Xc = 0 ; // cluster first moment along X |
1656 | Double_t Yc = 0 ; // cluster first moment along Y | |
1657 | Double_t Sxx = 0 ; // cluster second central moment along X (variance_X^2) | |
1658 | Double_t Sxy = 0 ; // cluster second central moment along Y (variance_Y^2) | |
1659 | Double_t Syy = 0 ; // cluster covariance^2 | |
90d5b88b | 1660 | |
1661 | AliVCaloCells *cells = fEsdCells; | |
1662 | if (!cells) | |
1663 | cells = fAodCells; | |
1664 | ||
6eb6260e | 1665 | if (!cells) |
1666 | return; | |
1667 | ||
6bf90832 | 1668 | Int_t ncells = c->GetNCells(); |
6eb6260e | 1669 | if (ncells==1) |
1670 | return; | |
1671 | ||
6eb6260e | 1672 | for(Int_t j=0; j<ncells; ++j) { |
6bf90832 | 1673 | Int_t id = TMath::Abs(c->GetCellAbsId(j)); |
6eb6260e | 1674 | Double_t cellen = cells->GetCellAmplitude(id); |
3a952328 | 1675 | TVector3 pos; |
1676 | fGeom->GetGlobal(id,pos); | |
6eb6260e | 1677 | Xc += cellen*pos.X(); |
1678 | Yc += cellen*pos.Y(); | |
1679 | Sxx += cellen*pos.X()*pos.X(); | |
1680 | Syy += cellen*pos.Y()*pos.Y(); | |
1681 | Sxy += cellen*pos.X()*pos.Y(); | |
1682 | } | |
1683 | Xc /= Ec; | |
1684 | Yc /= Ec; | |
1685 | Sxx /= Ec; | |
1686 | Syy /= Ec; | |
1687 | Sxy /= Ec; | |
1688 | Sxx -= Xc*Xc; | |
1689 | Syy -= Yc*Yc; | |
1690 | Sxy -= Xc*Yc; | |
f0897c18 | 1691 | Sxx = TMath::Abs(Sxx); |
1692 | Syy = TMath::Abs(Syy); | |
296ea9b4 | 1693 | sigmaMax = (Sxx + Syy + TMath::Sqrt(TMath::Abs((Sxx-Syy)*(Sxx-Syy)+4.0*Sxy*Sxy)))/2.0; |
1694 | sigmaMax = TMath::Sqrt(TMath::Abs(sigmaMax)); | |
1695 | sigmaMin = TMath::Abs(Sxx + Syy - TMath::Sqrt(TMath::Abs((Sxx-Syy)*(Sxx-Syy)+4.0*Sxy*Sxy)))/2.0; | |
1696 | sigmaMin = TMath::Sqrt(TMath::Abs(sigmaMin)); | |
d595acbb | 1697 | } |
90d5b88b | 1698 | |
6bf90832 | 1699 | //________________________________________________________________________ |
0fbe8d4f | 1700 | Int_t AliAnalysisTaskEMCALPi0PbPb::GetNCells(const AliVCluster *c, Double_t emin) const |
6bf90832 | 1701 | { |
1702 | // Calculate number of attached cells above emin. | |
1703 | ||
6bf90832 | 1704 | AliVCaloCells *cells = fEsdCells; |
1705 | if (!cells) | |
1706 | cells = fAodCells; | |
6bf90832 | 1707 | if (!cells) |
296ea9b4 | 1708 | return 0; |
6bf90832 | 1709 | |
296ea9b4 | 1710 | Int_t n = 0; |
6bf90832 | 1711 | Int_t ncells = c->GetNCells(); |
1712 | for(Int_t j=0; j<ncells; ++j) { | |
1713 | Int_t id = TMath::Abs(c->GetCellAbsId(j)); | |
1714 | Double_t cellen = cells->GetCellAmplitude(id); | |
1715 | if (cellen>=emin) | |
1716 | ++n; | |
1717 | } | |
1718 | return n; | |
1719 | } | |
1720 | ||
296ea9b4 | 1721 | //________________________________________________________________________ |
b6c599fe | 1722 | Double_t AliAnalysisTaskEMCALPi0PbPb::GetTrackIsolation(Double_t cEta, Double_t cPhi, Double_t radius, Double_t pt) const |
296ea9b4 | 1723 | { |
1724 | // Compute isolation based on tracks. | |
1725 | ||
1726 | Double_t trkIsolation = 0; | |
1727 | Double_t rad2 = radius*radius; | |
b6c599fe | 1728 | Int_t ntrks = fSelPrimTracks->GetEntries(); |
296ea9b4 | 1729 | for(Int_t j = 0; j<ntrks; ++j) { |
1730 | AliVTrack *track = static_cast<AliVTrack*>(fSelTracks->At(j)); | |
1731 | if (!track) | |
1732 | continue; | |
b6c599fe | 1733 | if (track->Pt()<pt) |
1734 | continue; | |
296ea9b4 | 1735 | Float_t eta = track->Eta(); |
1736 | Float_t phi = track->Phi(); | |
0fbe8d4f | 1737 | Double_t phidiff = TVector2::Phi_mpi_pi(phi-cPhi); |
296ea9b4 | 1738 | Double_t dist = (eta-cEta)*(eta-cEta)+phidiff*phidiff; |
1739 | if(dist>rad2) | |
1740 | continue; | |
1741 | trkIsolation += track->Pt(); | |
1742 | } | |
1743 | return trkIsolation; | |
1744 | } | |
3a952328 | 1745 | |
0fbe8d4f | 1746 | //________________________________________________________________________ |
1747 | Bool_t AliAnalysisTaskEMCALPi0PbPb::IsShared(const AliVCluster *c) const | |
1748 | { | |
1749 | // Returns if cluster shared across super modules. | |
1750 | ||
1751 | AliVCaloCells *cells = fEsdCells; | |
1752 | if (!cells) | |
1753 | cells = fAodCells; | |
1754 | if (!cells) | |
1755 | return 0; | |
1756 | ||
1757 | Int_t n = -1; | |
1758 | Int_t ncells = c->GetNCells(); | |
1759 | for(Int_t j=0; j<ncells; ++j) { | |
1760 | Int_t id = TMath::Abs(c->GetCellAbsId(j)); | |
1761 | Int_t got = id / (24*48); | |
1762 | if (n==-1) { | |
1763 | n = got; | |
1764 | continue; | |
1765 | } | |
1766 | if (got!=n) | |
1767 | return 1; | |
1768 | } | |
1769 | return 0; | |
1770 | } | |
3a952328 | 1771 |