| 1d08b6e8 |
1 | /************************************************************************** |
| 2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * |
| 3 | * * |
| 4 | * Author: The ALICE Off-line Project. * |
| 5 | * Contributors are mentioned in the code where appropriate. * |
| 6 | * * |
| 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 | **************************************************************************/ |
| 15 | |
| 16 | // Class for heavy-flavour electron v2 with EMCal triggered events |
| 17 | // Author: Denise Godoy |
| 18 | |
| 19 | |
| 20 | #include "TChain.h" |
| 21 | #include "TTree.h" |
| 22 | #include "TH2F.h" |
| 23 | #include "TMath.h" |
| 24 | #include "TCanvas.h" |
| 25 | #include "THnSparse.h" |
| 26 | #include "TLorentzVector.h" |
| 27 | #include "TString.h" |
| 28 | #include "TFile.h" |
| 29 | |
| 30 | #include "AliAnalysisTask.h" |
| 31 | #include "AliAnalysisManager.h" |
| 32 | |
| 33 | #include "AliESDEvent.h" |
| 34 | #include "AliESDHandler.h" |
| 35 | #include "AliAODEvent.h" |
| 36 | #include "AliAODHandler.h" |
| 37 | |
| 38 | #include "AliAnalysisTaskElecV2.h" |
| 39 | #include "TGeoGlobalMagField.h" |
| 40 | #include "AliLog.h" |
| 41 | #include "AliAnalysisTaskSE.h" |
| 42 | #include "TRefArray.h" |
| 43 | #include "TVector.h" |
| 44 | #include "AliESDInputHandler.h" |
| 45 | #include "AliESDpid.h" |
| 46 | #include "AliESDtrackCuts.h" |
| 47 | #include "AliPhysicsSelection.h" |
| 48 | #include "AliESDCaloCluster.h" |
| 49 | #include "AliAODCaloCluster.h" |
| 50 | #include "AliEMCALRecoUtils.h" |
| 51 | #include "AliEMCALGeometry.h" |
| 52 | #include "AliGeomManager.h" |
| 53 | #include "stdio.h" |
| 54 | #include "TGeoManager.h" |
| 55 | #include "iostream" |
| 56 | #include "fstream" |
| 57 | |
| 58 | #include "AliEMCALTrack.h" |
| 59 | #include "AliMagF.h" |
| 60 | |
| 61 | #include "AliKFParticle.h" |
| 62 | #include "AliKFVertex.h" |
| 63 | |
| 3a94d6fe |
64 | #include "AliMCEventHandler.h" |
| 65 | #include "AliMCEvent.h" |
| 66 | #include "AliMCParticle.h" |
| 67 | #include "AliStack.h" |
| 68 | |
| 1d08b6e8 |
69 | #include "AliPID.h" |
| 70 | #include "AliPIDResponse.h" |
| 71 | #include "AliHFEcontainer.h" |
| 72 | #include "AliHFEcuts.h" |
| 73 | #include "AliHFEpid.h" |
| 74 | #include "AliHFEpidBase.h" |
| 75 | #include "AliHFEpidQAmanager.h" |
| 76 | #include "AliHFEtools.h" |
| 77 | #include "AliCFContainer.h" |
| 78 | #include "AliCFManager.h" |
| 79 | |
| 80 | #include "AliEventplane.h" |
| 81 | #include "AliCentrality.h" |
| 82 | |
| 83 | ClassImp(AliAnalysisTaskElecV2) |
| 84 | //________________________________________________________________________ |
| 85 | AliAnalysisTaskElecV2::AliAnalysisTaskElecV2(const char *name) |
| 86 | : AliAnalysisTaskSE(name) |
| 87 | ,fESD(0) |
| 3a94d6fe |
88 | ,fMC(0) |
| 1d08b6e8 |
89 | ,fOutputList(0) |
| 90 | ,fTrackCuts(0) |
| 91 | ,fCuts(0) |
| 92 | ,fIdentifiedAsOutInz(kFALSE) |
| 93 | ,fPassTheEventCut(kFALSE) |
| 94 | ,fRejectKinkMother(kFALSE) |
| 3a94d6fe |
95 | ,fIsMC(kFALSE) |
| 1d08b6e8 |
96 | ,fVz(0.0) |
| 97 | ,fCFM(0) |
| 98 | ,fPID(0) |
| 99 | ,fPIDqa(0) |
| 100 | ,fOpeningAngleCut(0.1) |
| 101 | ,fInvmassCut(0.01) |
| 102 | ,fNoEvents(0) |
| 103 | ,fTrkpt(0) |
| 104 | ,fTrkEovPBef(0) |
| 105 | ,fTrkEovPAft(0) |
| 106 | ,fdEdxBef(0) |
| 107 | ,fdEdxAft(0) |
| 108 | ,fInvmassLS(0) |
| 109 | ,fInvmassULS(0) |
| 110 | ,fOpeningAngleLS(0) |
| 111 | ,fOpeningAngleULS(0) |
| 112 | ,fPhotoElecPt(0) |
| 113 | ,fSemiInclElecPt(0) |
| 3a94d6fe |
114 | ,fMCphotoElecPt(0) |
| 1d08b6e8 |
115 | ,fTrackPtBefTrkCuts(0) |
| 116 | ,fTrackPtAftTrkCuts(0) |
| 117 | ,fTPCnsigma(0) |
| 118 | ,fCent(0) |
| 44fdf7ac |
119 | ,fevPlaneV0A(0) |
| 120 | ,fevPlaneV0C(0) |
| 121 | ,fevPlaneTPC(0) |
| 1d08b6e8 |
122 | ,fTPCsubEPres(0) |
| 123 | ,fEPres(0) |
| 124 | ,fCorr(0) |
| 125 | ,feTPCV2(0) |
| 126 | ,feV2(0) |
| 127 | ,fphoteV2(0) |
| 128 | ,fChargPartV2(0) |
| 129 | { |
| 130 | //Named constructor |
| 131 | |
| 132 | fPID = new AliHFEpid("hfePid"); |
| 133 | fTrackCuts = new AliESDtrackCuts(); |
| 134 | |
| 135 | // Define input and output slots here |
| 136 | // Input slot #0 works with a TChain |
| 137 | DefineInput(0, TChain::Class()); |
| 138 | // Output slot #0 id reserved by the base class for AOD |
| 139 | // Output slot #1 writes into a TH1 container |
| 140 | // DefineOutput(1, TH1I::Class()); |
| 141 | DefineOutput(1, TList::Class()); |
| 142 | // DefineOutput(3, TTree::Class()); |
| 143 | } |
| 144 | |
| 145 | //________________________________________________________________________ |
| 146 | AliAnalysisTaskElecV2::AliAnalysisTaskElecV2() |
| 147 | : AliAnalysisTaskSE("DefaultAnalysis_AliAnalysisElecHadCorrel") |
| 148 | ,fESD(0) |
| 3a94d6fe |
149 | ,fMC(0) |
| 1d08b6e8 |
150 | ,fOutputList(0) |
| 151 | ,fTrackCuts(0) |
| 152 | ,fCuts(0) |
| 153 | ,fIdentifiedAsOutInz(kFALSE) |
| 154 | ,fPassTheEventCut(kFALSE) |
| 155 | ,fRejectKinkMother(kFALSE) |
| 3a94d6fe |
156 | ,fIsMC(kFALSE) |
| 1d08b6e8 |
157 | ,fVz(0.0) |
| 158 | ,fCFM(0) |
| 159 | ,fPID(0) |
| 160 | ,fPIDqa(0) |
| 161 | ,fOpeningAngleCut(0.1) |
| 162 | ,fInvmassCut(0.01) |
| 163 | ,fNoEvents(0) |
| 164 | ,fTrkpt(0) |
| 165 | ,fTrkEovPBef(0) |
| 166 | ,fTrkEovPAft(0) |
| 167 | ,fdEdxBef(0) |
| 168 | ,fdEdxAft(0) |
| 169 | ,fInvmassLS(0) |
| 170 | ,fInvmassULS(0) |
| 171 | ,fOpeningAngleLS(0) |
| 172 | ,fOpeningAngleULS(0) |
| 173 | ,fPhotoElecPt(0) |
| 174 | ,fSemiInclElecPt(0) |
| 3a94d6fe |
175 | ,fMCphotoElecPt(0) |
| 1d08b6e8 |
176 | ,fTrackPtBefTrkCuts(0) |
| 177 | ,fTrackPtAftTrkCuts(0) |
| 178 | ,fTPCnsigma(0) |
| 179 | ,fCent(0) |
| 44fdf7ac |
180 | ,fevPlaneV0A(0) |
| 181 | ,fevPlaneV0C(0) |
| 182 | ,fevPlaneTPC(0) |
| 1d08b6e8 |
183 | ,fTPCsubEPres(0) |
| 184 | ,fEPres(0) |
| 185 | ,fCorr(0) |
| 186 | ,feTPCV2(0) |
| 187 | ,feV2(0) |
| 188 | ,fphoteV2(0) |
| 189 | ,fChargPartV2(0) |
| 190 | { |
| 191 | //Default constructor |
| 192 | fPID = new AliHFEpid("hfePid"); |
| 193 | |
| 194 | fTrackCuts = new AliESDtrackCuts(); |
| 195 | |
| 196 | // Constructor |
| 197 | // Define input and output slots here |
| 198 | // Input slot #0 works with a TChain |
| 199 | DefineInput(0, TChain::Class()); |
| 200 | // Output slot #0 id reserved by the base class for AOD |
| 201 | // Output slot #1 writes into a TH1 container |
| 202 | // DefineOutput(1, TH1I::Class()); |
| 203 | DefineOutput(1, TList::Class()); |
| 204 | //DefineOutput(3, TTree::Class()); |
| 205 | } |
| 206 | //_________________________________________ |
| 207 | |
| 208 | AliAnalysisTaskElecV2::~AliAnalysisTaskElecV2() |
| 209 | { |
| 210 | //Destructor |
| 211 | |
| 212 | delete fOutputList; |
| 213 | delete fPID; |
| 214 | delete fCFM; |
| 215 | delete fPIDqa; |
| 216 | delete fTrackCuts; |
| 217 | } |
| 218 | //_________________________________________ |
| 219 | |
| 220 | void AliAnalysisTaskElecV2::UserExec(Option_t*) |
| 221 | { |
| 222 | //Main loop |
| 223 | //Called for each event |
| 224 | |
| 225 | // create pointer to event |
| 226 | fESD = dynamic_cast<AliESDEvent*>(InputEvent()); |
| 227 | if (!fESD) { |
| 228 | printf("ERROR: fESD not available\n"); |
| 229 | return; |
| 230 | } |
| 231 | |
| 232 | if(!fCuts){ |
| 233 | AliError("HFE cuts not available"); |
| 234 | return; |
| 235 | } |
| 236 | |
| 237 | if(!fPID->IsInitialized()){ |
| 238 | // Initialize PID with the given run number |
| 239 | AliWarning("PID not initialised, get from Run no"); |
| 240 | fPID->InitializePID(fESD->GetRunNumber()); |
| 241 | } |
| 242 | |
| 3a94d6fe |
243 | if(fIsMC)fMC = MCEvent(); |
| 244 | AliStack* stack = NULL; |
| 245 | if(fIsMC && fMC) stack = fMC->Stack(); |
| 1d08b6e8 |
246 | |
| 247 | Int_t fNOtrks = fESD->GetNumberOfTracks(); |
| 248 | const AliESDVertex *pVtx = fESD->GetPrimaryVertex(); |
| 249 | |
| 250 | Double_t pVtxZ = -999; |
| 251 | pVtxZ = pVtx->GetZ(); |
| 252 | |
| 253 | if(TMath::Abs(pVtxZ)>10) return; |
| 254 | fNoEvents->Fill(0); |
| 255 | |
| 256 | if(fNOtrks<2) return; |
| 257 | |
| 258 | AliPIDResponse *pidResponse = fInputHandler->GetPIDResponse(); |
| 259 | if(!pidResponse){ |
| 260 | AliDebug(1, "Using default PID Response"); |
| 261 | pidResponse = AliHFEtools::GetDefaultPID(kFALSE, fInputEvent->IsA() == AliAODEvent::Class()); |
| 262 | } |
| 263 | |
| 264 | fPID->SetPIDResponse(pidResponse); |
| 265 | |
| 266 | fCFM->SetRecEventInfo(fESD); |
| 267 | |
| 268 | Float_t cent = -1.; |
| 269 | AliCentrality *centrality = fESD->GetCentrality(); |
| 270 | cent = centrality->GetCentralityPercentile("V0M"); |
| 271 | fCent->Fill(cent); |
| 272 | |
| 273 | if(cent>90.) return; |
| 274 | |
| 275 | //Event planes |
| 276 | |
| 277 | Double_t evPlaneV0A = TVector2::Phi_0_2pi(fESD->GetEventplane()->GetEventplane("V0A",fESD,2)); |
| 278 | if(evPlaneV0A > TMath::Pi()) evPlaneV0A = evPlaneV0A - TMath::Pi(); |
| 5945f829 |
279 | fevPlaneV0A->Fill(evPlaneV0A,cent); |
| 1d08b6e8 |
280 | |
| 281 | Double_t evPlaneV0C = TVector2::Phi_0_2pi(fESD->GetEventplane()->GetEventplane("V0C",fESD,2)); |
| 282 | if(evPlaneV0C > TMath::Pi()) evPlaneV0C = evPlaneV0C - TMath::Pi(); |
| 5945f829 |
283 | fevPlaneV0C->Fill(evPlaneV0C,cent); |
| 1d08b6e8 |
284 | |
| 285 | AliEventplane* esdTPCep = fESD->GetEventplane(); |
| 44fdf7ac |
286 | TVector2 *standardQ = 0x0; |
| 1d08b6e8 |
287 | Double_t qx = -999., qy = -999.; |
| 44fdf7ac |
288 | standardQ = esdTPCep->GetQVector(); |
| ee89e744 |
289 | if(!standardQ)return; |
| 290 | |
| 291 | qx = standardQ->X(); |
| 292 | qy = standardQ->Y(); |
| 293 | |
| 1d08b6e8 |
294 | TVector2 qVectorfortrack; |
| 295 | qVectorfortrack.Set(qx,qy); |
| 296 | Float_t evPlaneTPC = TVector2::Phi_0_2pi(qVectorfortrack.Phi())/2.; |
| 5945f829 |
297 | fevPlaneTPC->Fill(evPlaneTPC,cent); |
| 44fdf7ac |
298 | |
| 1d08b6e8 |
299 | TVector2 *qsub1a = esdTPCep->GetQsub1(); |
| 300 | TVector2 *qsub2a = esdTPCep->GetQsub2(); |
| 301 | Double_t evPlaneResTPC = -999.; |
| 302 | if(qsub1a && qsub2a) |
| 303 | { |
| 304 | evPlaneResTPC = TMath::Cos(2.*TVector2::Phi_0_2pi(qsub1a->Phi()/2.- qsub2a->Phi()/2.)); |
| 305 | } |
| 306 | |
| 307 | fTPCsubEPres->Fill(evPlaneResTPC,cent); |
| 308 | |
| 309 | Double_t evPlaneRes[4]={GetCos2DeltaPhi(evPlaneV0A,evPlaneV0C),GetCos2DeltaPhi(evPlaneV0A,evPlaneTPC),GetCos2DeltaPhi(evPlaneV0C,evPlaneTPC),cent}; |
| 310 | fEPres->Fill(evPlaneRes); |
| 311 | |
| 312 | // Track loop |
| 313 | for (Int_t iTracks = 0; iTracks < fESD->GetNumberOfTracks(); iTracks++) { |
| 314 | AliESDtrack* track = fESD->GetTrack(iTracks); |
| 315 | if (!track) { |
| 316 | printf("ERROR: Could not receive track %d\n", iTracks); |
| 317 | continue; |
| 318 | } |
| 319 | |
| 320 | if(TMath::Abs(track->Eta())>0.7) continue; |
| 321 | |
| 322 | fTrackPtBefTrkCuts->Fill(track->Pt()); |
| b1a1875f |
323 | |
| 1d08b6e8 |
324 | if(!ProcessCutStep(AliHFEcuts::kStepRecKineITSTPC, track)) continue; |
| 325 | |
| 1d08b6e8 |
326 | if(fRejectKinkMother) { // Quick and dirty fix to reject both kink mothers and daughters |
| 327 | if(track->GetKinkIndex(0) != 0) continue; |
| 328 | } |
| 329 | |
| 1d08b6e8 |
330 | if(!ProcessCutStep(AliHFEcuts::kStepRecPrim, track)) continue; |
| 331 | |
| 1d08b6e8 |
332 | if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsITS, track)) continue; |
| 333 | |
| 1d08b6e8 |
334 | if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTPC, track)) continue; |
| 335 | |
| 336 | fTrackPtAftTrkCuts->Fill(track->Pt()); |
| 337 | |
| 06d74a7b |
338 | Double_t clsE = -999., p = -999., EovP=-999., pt = -999., dEdx=-999., fTPCnSigma=0, phi=-999., wclsE = -999., wEovP = -999.; |
| 1d08b6e8 |
339 | |
| 1d08b6e8 |
340 | |
| 341 | pt = track->Pt(); |
| 44fdf7ac |
342 | if(pt<2) continue; |
| 1d08b6e8 |
343 | fTrkpt->Fill(pt); |
| 3a94d6fe |
344 | |
| 1d08b6e8 |
345 | Int_t clsId = track->GetEMCALcluster(); |
| 346 | if (clsId>0){ |
| 347 | AliESDCaloCluster *cluster = fESD->GetCaloCluster(clsId); |
| 348 | if(cluster && cluster->IsEMCAL()){ |
| 349 | clsE = cluster->E(); |
| 06d74a7b |
350 | // m20 = cluster->M20(); |
| 351 | // m02 = cluster->M02(); |
| 1d08b6e8 |
352 | } |
| 353 | } |
| b1a1875f |
354 | |
| 355 | |
| 356 | |
| 1d08b6e8 |
357 | p = track->P(); |
| 358 | phi = track->Phi(); |
| 359 | dEdx = track->GetTPCsignal(); |
| 360 | EovP = clsE/p; |
| b1a1875f |
361 | wEovP = wclsE/p; |
| 1d08b6e8 |
362 | fTPCnSigma = fPID->GetPIDResponse() ? fPID->GetPIDResponse()->NumberOfSigmasTPC(track, AliPID::kElectron) : 1000; |
| 363 | fdEdxBef->Fill(p,dEdx); |
| 364 | fTPCnsigma->Fill(p,fTPCnSigma); |
| 365 | |
| 06d74a7b |
366 | Bool_t fFlagPhotonicElec = kFALSE; |
| 367 | Bool_t fFlagPhotonicElecBCG = kFALSE; |
| 5945f829 |
368 | |
| 06d74a7b |
369 | SelectPhotonicElectron(iTracks,track, fFlagPhotonicElec, fFlagPhotonicElecBCG); |
| 370 | |
| 371 | Double_t corr[10]={phi,fTPCnSigma,cent,pt,EovP,GetDeltaPhi(phi,evPlaneTPC),GetDeltaPhi(phi,evPlaneV0A),GetDeltaPhi(phi,evPlaneV0C),fFlagPhotonicElec,fFlagPhotonicElecBCG}; |
| 1d08b6e8 |
372 | fCorr->Fill(corr); |
| 06d74a7b |
373 | |
| 374 | if(fIsMC && fMC && stack){ |
| 375 | Int_t label = track->GetLabel(); |
| 376 | if(label>0){ |
| 377 | TParticle *particle = stack->Particle(label); |
| 378 | if(particle){ |
| 379 | Int_t partPDG = particle->GetPdgCode(); |
| 380 | Double_t partPt = particle->Pt(); |
| 381 | Int_t IsElec = 0; |
| 382 | if (TMath::Abs(partPDG)==11) IsElec = 1; |
| 383 | |
| 384 | Int_t idMother = particle->GetFirstMother(); |
| 385 | if (idMother>0){ |
| 386 | TParticle *mother = stack->Particle(idMother); |
| 387 | Int_t motherPDG = mother->GetPdgCode(); |
| 388 | |
| 389 | Double_t mc[8]={EovP,fTPCnSigma,partPt,fFlagPhotonicElec,fFlagPhotonicElecBCG,IsElec,cent,pt}; |
| 390 | |
| 391 | if (motherPDG==22 || motherPDG==111 || motherPDG==221) fMCphotoElecPt->Fill(mc);// gamma, pi0, eta |
| 392 | } |
| 393 | } |
| 394 | } |
| 395 | } |
| 396 | |
| 397 | |
| 1d08b6e8 |
398 | |
| 399 | if(fTPCnSigma >= 1.5 && fTPCnSigma <= 3)fTrkEovPBef->Fill(pt,EovP); |
| 400 | Int_t pidpassed = 1; |
| 401 | |
| 402 | //--- track accepted |
| 403 | AliHFEpidObject hfetrack; |
| 404 | hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis); |
| 405 | hfetrack.SetRecTrack(track); |
| 406 | hfetrack.SetPbPb(); |
| 407 | if(!fPID->IsSelected(&hfetrack, NULL, "", fPIDqa)) pidpassed = 0; |
| 408 | |
| 44fdf7ac |
409 | Double_t corrV2[7]={phi,cent,pt,EovP,GetCos2DeltaPhi(phi,evPlaneTPC),GetCos2DeltaPhi(phi,evPlaneV0A),GetCos2DeltaPhi(phi,evPlaneV0C)}; |
| 1d08b6e8 |
410 | fChargPartV2->Fill(corrV2); |
| b3bd01b3 |
411 | |
| 412 | if(fTPCnSigma >= -0.5){ |
| 413 | Double_t qX = standardQ->X() - esdTPCep->GetQContributionX(track); |
| 414 | Double_t qY = standardQ->Y() - esdTPCep->GetQContributionY(track); |
| 415 | TVector2 newQVectorfortrack; |
| 416 | newQVectorfortrack.Set(qX,qY); |
| 417 | Double_t corrV2TPC = -999.; |
| 418 | corrV2TPC = TVector2::Phi_0_2pi(newQVectorfortrack.Phi())/2; |
| 419 | |
| 420 | Double_t correctedV2[5]={cent,pt,GetCos2DeltaPhi(phi,corrV2TPC),GetCos2DeltaPhi(phi,evPlaneV0A),GetCos2DeltaPhi(phi,evPlaneV0C)}; |
| 421 | |
| 422 | feTPCV2->Fill(correctedV2); |
| 423 | } |
| 424 | |
| 425 | if(pidpassed==0) continue; |
| 426 | |
| 1d08b6e8 |
427 | Double_t qX = standardQ->X() - esdTPCep->GetQContributionX(track); |
| 428 | Double_t qY = standardQ->Y() - esdTPCep->GetQContributionY(track); |
| 429 | TVector2 newQVectorfortrack; |
| 430 | newQVectorfortrack.Set(qX,qY); |
| 431 | Double_t corrV2TPC = -999.; |
| 432 | corrV2TPC = TVector2::Phi_0_2pi(newQVectorfortrack.Phi())/2; |
| 433 | |
| 434 | Double_t correctedV2[5]={cent,pt,GetCos2DeltaPhi(phi,corrV2TPC),GetCos2DeltaPhi(phi,evPlaneV0A),GetCos2DeltaPhi(phi,evPlaneV0C)}; |
| 435 | |
| 1d08b6e8 |
436 | feV2->Fill(correctedV2); |
| 437 | |
| 438 | fTrkEovPAft->Fill(pt,EovP); |
| 439 | fdEdxAft->Fill(p,dEdx); |
| 440 | |
| 3a94d6fe |
441 | if(fFlagPhotonicElec){ |
| 442 | fphoteV2->Fill(correctedV2); |
| 443 | fPhotoElecPt->Fill(pt); |
| 444 | } |
| 445 | |
| 446 | if(!fFlagPhotonicElec) fSemiInclElecPt->Fill(pt); |
| 447 | |
| 1d08b6e8 |
448 | } |
| 449 | PostData(1, fOutputList); |
| 450 | } |
| 451 | //_________________________________________ |
| 452 | void AliAnalysisTaskElecV2::UserCreateOutputObjects() |
| 453 | { |
| 3a94d6fe |
454 | //--- Check MC |
| 455 | if(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()){ |
| 456 | fIsMC = kTRUE; |
| 457 | printf("+++++ MC Data available"); |
| 458 | } |
| 1d08b6e8 |
459 | //--------Initialize PID |
| 3a94d6fe |
460 | fPID->SetHasMCData(fIsMC); |
| 461 | |
| 1d08b6e8 |
462 | if(!fPID->GetNumberOfPIDdetectors()) |
| 463 | { |
| 464 | fPID->AddDetector("TPC", 0); |
| 465 | fPID->AddDetector("EMCAL", 1); |
| 466 | } |
| 467 | |
| 468 | fPID->SortDetectors(); |
| 469 | fPIDqa = new AliHFEpidQAmanager(); |
| 470 | fPIDqa->Initialize(fPID); |
| 471 | |
| 472 | //--------Initialize correction Framework and Cuts |
| 473 | fCFM = new AliCFManager; |
| 474 | const Int_t kNcutSteps = AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kNcutStepsRecTrack + AliHFEcuts::kNcutStepsDETrack; |
| 475 | fCFM->SetNStepParticle(kNcutSteps); |
| 476 | for(Int_t istep = 0; istep < kNcutSteps; istep++) |
| 477 | fCFM->SetParticleCutsList(istep, NULL); |
| 478 | |
| 479 | if(!fCuts){ |
| 480 | AliWarning("Cuts not available. Default cuts will be used"); |
| 481 | fCuts = new AliHFEcuts; |
| 482 | fCuts->CreateStandardCuts(); |
| 483 | } |
| 484 | fCuts->Initialize(fCFM); |
| 485 | |
| 486 | //---------Output Tlist |
| 487 | fOutputList = new TList(); |
| 488 | fOutputList->SetOwner(); |
| 489 | fOutputList->Add(fPIDqa->MakeList("PIDQA")); |
| 490 | |
| 491 | fNoEvents = new TH1F("fNoEvents","",1,0,1) ; |
| 492 | fOutputList->Add(fNoEvents); |
| 493 | |
| 494 | fTrkpt = new TH1F("fTrkpt","track pt",100,0,50); |
| 495 | fOutputList->Add(fTrkpt); |
| 496 | |
| 497 | fTrackPtBefTrkCuts = new TH1F("fTrackPtBefTrkCuts","track pt before track cuts",100,0,50); |
| 498 | fOutputList->Add(fTrackPtBefTrkCuts); |
| 499 | |
| 500 | fTrackPtAftTrkCuts = new TH1F("fTrackPtAftTrkCuts","track pt after track cuts",100,0,50); |
| 501 | fOutputList->Add(fTrackPtAftTrkCuts); |
| 502 | |
| 503 | fTPCnsigma = new TH2F("fTPCnsigma", "TPC - n sigma",100,0,50,200,-10,10); |
| 504 | fOutputList->Add(fTPCnsigma); |
| 505 | |
| 506 | fTrkEovPBef = new TH2F("fTrkEovPBef","track E/p before HFE pid",100,0,50,100,0,2); |
| 507 | fOutputList->Add(fTrkEovPBef); |
| 508 | |
| 509 | fTrkEovPAft = new TH2F("fTrkEovPAft","track E/p after HFE pid",100,0,50,100,0,2); |
| 510 | fOutputList->Add(fTrkEovPAft); |
| 511 | |
| 512 | fdEdxBef = new TH2F("fdEdxBef","track dEdx vs p before HFE pid",100,0,50,150,0,150); |
| 513 | fOutputList->Add(fdEdxBef); |
| 514 | |
| 515 | fdEdxAft = new TH2F("fdEdxAft","track dEdx vs p after HFE pid",100,0,50,150,0,150); |
| 516 | fOutputList->Add(fdEdxAft); |
| 517 | |
| 518 | fInvmassLS = new TH1F("fInvmassLS", "Inv mass of LS (e,e); mass(GeV/c^2); counts;", 500,0,0.5); |
| 519 | fOutputList->Add(fInvmassLS); |
| 520 | |
| 521 | fInvmassULS = new TH1F("fInvmassULS", "Inv mass of ULS (e,e); mass(GeV/c^2); counts;", 500,0,0.5); |
| 522 | fOutputList->Add(fInvmassULS); |
| 523 | |
| 524 | fOpeningAngleLS = new TH1F("fOpeningAngleLS","Opening angle for LS pairs",100,0,1); |
| 525 | fOutputList->Add(fOpeningAngleLS); |
| 526 | |
| 527 | fOpeningAngleULS = new TH1F("fOpeningAngleULS","Opening angle for ULS pairs",100,0,1); |
| 528 | fOutputList->Add(fOpeningAngleULS); |
| 529 | |
| 530 | fPhotoElecPt = new TH1F("fPhotoElecPt", "photonic electron pt",100,0,50); |
| 531 | fOutputList->Add(fPhotoElecPt); |
| 532 | |
| 533 | fSemiInclElecPt = new TH1F("fSemiInclElecPt", "Semi-inclusive electron pt",100,0,50); |
| 534 | fOutputList->Add(fSemiInclElecPt); |
| 535 | |
| 536 | fCent = new TH1F("fCent","Centrality",100,0,100) ; |
| 537 | fOutputList->Add(fCent); |
| 538 | |
| 5945f829 |
539 | fevPlaneV0A = new TH2F("fevPlaneV0A","V0A EP",100,0,TMath::Pi(),90,0,90); |
| 44fdf7ac |
540 | fOutputList->Add(fevPlaneV0A); |
| 541 | |
| 5945f829 |
542 | fevPlaneV0C = new TH2F("fevPlaneV0C","V0C EP",100,0,TMath::Pi(),90,0,90); |
| 44fdf7ac |
543 | fOutputList->Add(fevPlaneV0C); |
| 544 | |
| 5945f829 |
545 | fevPlaneTPC = new TH2F("fevPlaneTPC","TPC EP",100,0,TMath::Pi(),90,0,90); |
| 44fdf7ac |
546 | fOutputList->Add(fevPlaneTPC); |
| 547 | |
| 1d08b6e8 |
548 | fTPCsubEPres = new TH2F("fTPCsubEPres","TPC subevent plane resolution",100,-1,1,90,0,90); |
| 549 | fOutputList->Add(fTPCsubEPres); |
| 550 | |
| 551 | Int_t binsv1[4]={100,100,100,90}; // V0A-V0C, V0A-TPC, V0C-TPC, cent |
| 552 | Double_t xminv1[4]={-1,-1,-1,0}; |
| 553 | Double_t xmaxv1[4]={1,1,1,90}; |
| 554 | fEPres = new THnSparseD ("fEPres","EP resolution",4,binsv1,xminv1,xmaxv1); |
| 555 | fOutputList->Add(fEPres); |
| 556 | |
| 06d74a7b |
557 | //phi,fTPCnSigma,cent,pt,EovP,GetDeltaPhi(phi,evPlaneTPC),GetDeltaPhi(phi,evPlaneV0A),GetDeltaPhi(phi,evPlaneV0C),fFlagPhotonicElec,fFlagPhotonicElecBCG,m20,m02 |
| 558 | Int_t binsv2[10]={100,100,90,100,100,100,100,100,3,3}; |
| 559 | Double_t xminv2[10]={0,-3.5,0,0,0,0,0,0,-1,-1}; |
| 560 | Double_t xmaxv2[10]={2*TMath::Pi(),3.5,90,50,3,TMath::Pi(),TMath::Pi(),TMath::Pi(),2,2}; |
| 561 | fCorr = new THnSparseD ("fCorr","Correlations",10,binsv2,xminv2,xmaxv2); |
| 1d08b6e8 |
562 | fOutputList->Add(fCorr); |
| 563 | |
| 564 | Int_t binsv3[5]={90,100,100,100,100}; // cent, pt, TPCcos2DeltaPhi, V0Acos2DeltaPhi, V0Ccos2DeltaPhi |
| 565 | Double_t xminv3[5]={0,0,-1,-1,-1}; |
| 566 | Double_t xmaxv3[5]={90,50,1,1,1}; |
| 567 | feV2 = new THnSparseD ("feV2","inclusive electron v2",5,binsv3,xminv3,xmaxv3); |
| 568 | fOutputList->Add(feV2); |
| 569 | |
| 570 | Int_t binsv4[5]={90,100,100,100,100}; // cent, pt, TPCdeltaPhi, V0AdeltaPhi, V0CdeltaPhi |
| 571 | Double_t xminv4[5]={0,0,-1,-1,-1}; |
| 572 | Double_t xmaxv4[5]={90,50,1,1,1}; |
| 573 | fphoteV2 = new THnSparseD ("fphoteV2","photonic electron v2",5,binsv4,xminv4,xmaxv4); |
| 574 | fOutputList->Add(fphoteV2); |
| 575 | |
| 44fdf7ac |
576 | Int_t binsv5[7]={100,90,100,100,100,100,100}; // phi, cent, pt, EovP, TPCdeltaPhi, V0AdeltaPhi, V0CdeltaPhi |
| 577 | Double_t xminv5[7]={0,0,0,0,-1,-1,-1}; |
| 578 | Double_t xmaxv5[7]={2*TMath::Pi(),90,50,3,1,1,1}; |
| 579 | fChargPartV2 = new THnSparseD ("fChargPartV2","Charged particle v2",7,binsv5,xminv5,xmaxv5); |
| 1d08b6e8 |
580 | fOutputList->Add(fChargPartV2); |
| 581 | |
| 582 | Int_t binsv6[5]={90,100,100,100,100}; // cent, pt, TPCdeltaPhi, V0AdeltaPhi, V0CdeltaPhi |
| 583 | Double_t xminv6[5]={0,0,-1,-1,-1}; |
| 584 | Double_t xmaxv6[5]={90,50,1,1,1}; |
| 585 | feTPCV2 = new THnSparseD ("feTPCV2","inclusive electron v2 (TPC)",5,binsv6,xminv6,xmaxv6); |
| 586 | fOutputList->Add(feTPCV2); |
| 587 | |
| 06d74a7b |
588 | //EovP,fTPCnSigma,partPt,fFlagPhotonicElec,fFlagPhotonicElecBCG,IsElec,cent,pt}; |
| 589 | Int_t binsv7[8]={100,100,100,3,3,3,90,100}; |
| 590 | Double_t xminv7[8]={0,-3.5,0,-1,-1,-1,0,0}; |
| 591 | Double_t xmaxv7[8]={3,3.5,50,2,2,2,90,50}; |
| 592 | fMCphotoElecPt = new THnSparseD ("fMCphotoElecPt", "pt distribution (MC)",8,binsv7,xminv7,xmaxv7); |
| 3a94d6fe |
593 | fOutputList->Add(fMCphotoElecPt); |
| 594 | |
| 1d08b6e8 |
595 | PostData(1,fOutputList); |
| 596 | } |
| 597 | |
| 598 | //________________________________________________________________________ |
| 599 | void AliAnalysisTaskElecV2::Terminate(Option_t *) |
| 600 | { |
| 601 | // Info("Terminate"); |
| 602 | AliAnalysisTaskSE::Terminate(); |
| 603 | } |
| 604 | |
| 605 | //________________________________________________________________________ |
| 606 | Bool_t AliAnalysisTaskElecV2::ProcessCutStep(Int_t cutStep, AliVParticle *track) |
| 607 | { |
| 608 | // Check single track cuts for a given cut step |
| 609 | const Int_t kMCOffset = AliHFEcuts::kNcutStepsMCTrack; |
| 610 | if(!fCFM->CheckParticleCuts(cutStep + kMCOffset, track)) return kFALSE; |
| 611 | return kTRUE; |
| 612 | } |
| 613 | //_________________________________________ |
| 06d74a7b |
614 | void AliAnalysisTaskElecV2::SelectPhotonicElectron(Int_t iTracks,AliESDtrack *track,Bool_t &fFlagPhotonicElec, Bool_t &fFlagPhotonicElecBCG) |
| 1d08b6e8 |
615 | { |
| 616 | //Identify non-heavy flavour electrons using Invariant mass method |
| 617 | |
| 618 | fTrackCuts->SetAcceptKinkDaughters(kFALSE); |
| 619 | fTrackCuts->SetRequireTPCRefit(kTRUE); |
| 5945f829 |
620 | fTrackCuts->SetRequireITSRefit(kTRUE); |
| 1d08b6e8 |
621 | fTrackCuts->SetEtaRange(-0.7,0.7); |
| 622 | fTrackCuts->SetRequireSigmaToVertex(kTRUE); |
| 623 | fTrackCuts->SetMaxChi2PerClusterTPC(3.5); |
| 624 | fTrackCuts->SetMinNClustersTPC(100); |
| 625 | |
| 626 | const AliESDVertex *pVtx = fESD->GetPrimaryVertex(); |
| 627 | |
| 628 | Bool_t flagPhotonicElec = kFALSE; |
| 06d74a7b |
629 | Bool_t flagPhotonicElecBCG = kFALSE; |
| 1d08b6e8 |
630 | |
| b1a1875f |
631 | for(Int_t jTracks = 0; jTracks<fESD->GetNumberOfTracks(); jTracks++){ |
| 632 | |
| 5945f829 |
633 | if(jTracks==iTracks) continue; |
| b1a1875f |
634 | |
| 1d08b6e8 |
635 | AliESDtrack* trackAsso = fESD->GetTrack(jTracks); |
| 636 | if (!trackAsso) { |
| 637 | printf("ERROR: Could not receive track %d\n", jTracks); |
| 638 | continue; |
| 639 | } |
| 640 | |
| 5945f829 |
641 | Double_t dEdxAsso = -999., ptAsso=-999.; |
| 06d74a7b |
642 | Bool_t fFlagLS=kFALSE, fFlagULS=kFALSE; |
| 643 | Double_t openingAngle = -999., mass=999., width = -999; |
| 1d08b6e8 |
644 | |
| 645 | dEdxAsso = trackAsso->GetTPCsignal(); |
| 646 | ptAsso = trackAsso->Pt(); |
| 647 | Int_t chargeAsso = trackAsso->Charge(); |
| 648 | Int_t charge = track->Charge(); |
| 649 | |
| 06d74a7b |
650 | if(ptAsso <0.5) continue; |
| 1d08b6e8 |
651 | if(!fTrackCuts->AcceptTrack(trackAsso)) continue; |
| 06d74a7b |
652 | if(dEdxAsso <65 || dEdxAsso>100) continue; |
| 1d08b6e8 |
653 | |
| 654 | Int_t fPDGe1 = 11; Int_t fPDGe2 = 11; |
| 655 | if(charge>0) fPDGe1 = -11; |
| 656 | if(chargeAsso>0) fPDGe2 = -11; |
| 657 | |
| 658 | if(charge == chargeAsso) fFlagLS = kTRUE; |
| 659 | if(charge != chargeAsso) fFlagULS = kTRUE; |
| 660 | |
| 661 | AliKFParticle ge1(*track, fPDGe1); |
| 662 | AliKFParticle ge2(*trackAsso, fPDGe2); |
| 663 | AliKFParticle recg(ge1, ge2); |
| 664 | |
| 665 | if(recg.GetNDF()<1) continue; |
| 666 | Double_t chi2recg = recg.GetChi2()/recg.GetNDF(); |
| 667 | if(TMath::Sqrt(TMath::Abs(chi2recg))>3.) continue; |
| 668 | |
| 669 | AliKFVertex primV(*pVtx); |
| 670 | primV += recg; |
| 671 | recg.SetProductionVertex(primV); |
| 672 | |
| 673 | recg.SetMassConstraint(0,0.0001); |
| 674 | |
| 675 | openingAngle = ge1.GetAngle(ge2); |
| 676 | if(fFlagLS) fOpeningAngleLS->Fill(openingAngle); |
| 677 | if(fFlagULS) fOpeningAngleULS->Fill(openingAngle); |
| 678 | |
| 5945f829 |
679 | if(openingAngle > fOpeningAngleCut) continue; |
| 06d74a7b |
680 | |
| 681 | recg.GetMass(mass,width); |
| 5945f829 |
682 | |
| 1d08b6e8 |
683 | if(fFlagLS) fInvmassLS->Fill(mass); |
| 684 | if(fFlagULS) fInvmassULS->Fill(mass); |
| 685 | |
| 5945f829 |
686 | if(mass<fInvmassCut && fFlagULS && !flagPhotonicElec) flagPhotonicElec = kTRUE; |
| 06d74a7b |
687 | if(mass<fInvmassCut && fFlagLS && !flagPhotonicElecBCG) flagPhotonicElecBCG = kTRUE; |
| 1d08b6e8 |
688 | |
| 689 | } |
| 690 | fFlagPhotonicElec = flagPhotonicElec; |
| 06d74a7b |
691 | fFlagPhotonicElecBCG = flagPhotonicElecBCG; |
| 1d08b6e8 |
692 | |
| 693 | } |
| 694 | //_________________________________________ |
| 695 | Double_t AliAnalysisTaskElecV2::GetCos2DeltaPhi(Double_t phiA,Double_t phiB) const |
| 696 | { |
| 697 | //Get cos[2(phi-psi_EP)] or cos[2(psi_subEP1 - psi_subEP2)] |
| 698 | Double_t dPhi = TVector2::Phi_0_2pi(phiA - phiB); |
| 699 | if(dPhi > TMath::Pi()) dPhi = dPhi - TMath::Pi(); |
| 700 | Double_t cos2DeltaPhi = TMath::Cos(2*dPhi); |
| 701 | |
| 702 | return cos2DeltaPhi; |
| 703 | } |
| 704 | |
| 705 | //_________________________________________ |
| 706 | Double_t AliAnalysisTaskElecV2::GetDeltaPhi(Double_t phiA,Double_t phiB) const |
| 707 | { |
| 708 | //Get phi-psi_EP |
| 709 | Double_t dPhi = TVector2::Phi_0_2pi(phiA - phiB); |
| 710 | if(dPhi > TMath::Pi()) dPhi = dPhi - TMath::Pi(); |
| 711 | |
| 712 | return dPhi; |
| 713 | } |
| b1a1875f |
714 | //_________________________________________ |
| 715 | Double_t AliAnalysisTaskElecV2::GetclusterE(Int_t iTrack, Double_t clsPhi, Double_t clsEta) const |
| 716 | { |
| 717 | //Return E |
| 718 | for (Int_t jTracks = 0; jTracks < fESD->GetNumberOfTracks(); jTracks++){ |
| 719 | |
| 720 | if(jTracks==iTrack) continue; |
| 721 | |
| 722 | AliESDtrack* wtrack = fESD->GetTrack(jTracks); |
| 723 | if (!wtrack) continue; |
| 724 | |
| 725 | Double_t wclsPhi=-999., wclsEta=-999., dPhi=-999., dEta=-999., dR=-999., wclsE=-999.; |
| 726 | |
| 727 | Int_t wclsId = wtrack->GetEMCALcluster(); |
| 728 | if (wclsId>0){ |
| 729 | AliESDCaloCluster *wcluster = fESD->GetCaloCluster(wclsId); |
| 730 | if(wcluster && wcluster->IsEMCAL()){ |
| 731 | Float_t wclusterPosition[3]={0,0,0}; |
| 732 | wcluster->GetPosition(wclusterPosition); |
| 733 | TVector3 clsPosVec(wclusterPosition[0],wclusterPosition[1],wclusterPosition[2]); |
| 734 | wclsPhi = clsPosVec.Phi(); |
| 735 | wclsEta = clsPosVec.Eta(); |
| 736 | |
| 737 | dPhi = TMath::Abs(wclsPhi - clsPhi); |
| 738 | dEta = TMath::Abs(wclsEta - clsEta); |
| 739 | dR = TMath::Sqrt(dPhi*dPhi+dEta*dEta); |
| 740 | if(dR>0.15){ |
| 741 | wclsE = wcluster->E(); |
| 742 | return wclsE; |
| 743 | } |
| 744 | } |
| 745 | } |
| 746 | } |
| 747 | return -999.; |
| 748 | } |
git://git.uio.no / u / mrichter / AliRoot.git / blame