30 #include "THnSparse.h" 38 #include "AliAODInputHandler.h" 39 #include "AliAODHandler.h" 40 #include "AliESDEvent.h" 41 #include "AliAODMCParticle.h" 43 #include "AliAODJetEventBackground.h" 44 #include "AliGenPythiaEventHeader.h" 45 #include "AliGenHijingEventHeader.h" 46 #include "AliInputEventHandler.h" 50 #include "AliAnalysisManager.h" 51 #include "AliAnalysisTaskSE.h" 52 #include "AliAnalysisUtils.h" 53 #include "AliVParticle.h" 54 #include "AliVEvent.h" 57 #include "AliPIDResponse.h" 75 ,fCentralityEstimator("V0M")
76 ,fNameTrackContainer("Tracks")
77 ,fNameTrackContainerEfficiency("")
78 ,fNameMCParticleContainer("")
79 ,fNameJetContainer("")
80 ,fNameMCParticleJetContainer("")
81 ,fUseAODInputJets(kTRUE)
82 ,fUsePhysicsSelection(kTRUE)
92 ,fStudyTransversalJetStructure(kFALSE)
96 ,fh1VertexNContributors(0)
109 ,fhJetPtRefMultEta5(0)
110 ,fhJetPtRefMultEta8(0)
111 ,fhJetPtMultPercent(0)
115 ,fOnlyLeadingJets(kFALSE)
121 ,fNumInclusivePIDtasks(0)
123 ,fNumJetUEPIDtasks(0)
124 ,fNameInclusivePIDtask(0x0)
125 ,fNameJetPIDtask(0x0)
126 ,fNameJetUEPIDtask(0x0)
127 ,fInclusivePIDtask(0x0)
130 ,fUseInclusivePIDtask(kFALSE)
131 ,fUseJetPIDtask(kFALSE)
132 ,fUseJetUEPIDtask(kFALSE)
141 for (
Int_t i = 0; i < AliPID::kSPECIES; i++) {
143 fhDCA_XY_prim_MCID[i] = 0x0;
144 fhDCA_Z_prim_MCID[i] = 0x0;
146 fhDCA_XY_sec_MCID[i] = 0x0;
147 fhDCA_Z_sec_MCID[i] = 0x0;
160 ,fCentralityEstimator(
"V0M")
161 ,fNameTrackContainer(
"Tracks")
162 ,fNameTrackContainerEfficiency(
"")
163 ,fNameMCParticleContainer(
"")
164 ,fNameJetContainer(
"")
165 ,fNameMCParticleJetContainer(
"")
166 ,fUseAODInputJets(kTRUE)
167 ,fUsePhysicsSelection(kTRUE)
168 ,fEvtSelectionMask(0)
177 ,fStudyTransversalJetStructure(kFALSE)
181 ,fh1VertexNContributors(0)
194 ,fhJetPtRefMultEta5(0)
195 ,fhJetPtRefMultEta8(0)
196 ,fhJetPtMultPercent(0)
198 ,fOnlyLeadingJets(kFALSE)
202 ,fNumInclusivePIDtasks(0)
204 ,fNumJetUEPIDtasks(0)
205 ,fNameInclusivePIDtask(0x0)
206 ,fNameJetPIDtask(0x0)
207 ,fNameJetUEPIDtask(0x0)
208 ,fInclusivePIDtask(0x0)
211 ,fUseInclusivePIDtask(kFALSE)
212 ,fUseJetPIDtask(kFALSE)
213 ,fUseJetUEPIDtask(kFALSE)
222 for (
Int_t i = 0; i < AliPID::kSPECIES; i++) {
232 DefineOutput(1,TList::Class());
272 TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
278 TFile *curfile = tree->GetCurrentFile();
280 Error(
"Notify",
"No current file");
302 Printf(
"%s:%d No Histogram fh1Xsec",(
char*)__FILE__,__LINE__);
306 fh1Xsec->Fill(
"<#sigma>",xsection);
309 if(ftrials>=nEntries && nEntries>0.)
fAvgTrials = ftrials/nEntries;
325 if(fDebug > 1) Printf(
"AliAnalysisTaskIDFragmentationFunction::UserCreateOutputObjects()");
335 Bool_t oldStatus = TH1::AddDirectoryStatus();
336 TH1::AddDirectory(kFALSE);
340 fh1EvtSelection =
new TH1F(
"fh1EvtSelection",
"Event Selection", 6, -0.5, 5.5);
342 fh1EvtSelection->GetXaxis()->SetBinLabel(2,
"event selection: rejected");
348 fh1VtxSelection =
new TH1F(
"fh1VtxSelection",
"Vertex Selection", 10, -1, 9);
361 fh1VertexZ =
new TH1F(
"fh1VertexZ",
"Vertex z distribution", 30, -15., 15.);
362 fh1EvtMult =
new TH1F(
"fh1EvtMult",
"Event multiplicity, track pT cut > 150 MeV/c, |#eta| < 0.9",120,0.,12000.);
363 fh1EvtCent =
new TH1F(
"fh1EvtCent",
"centrality",100,0.,100.);
365 fh1Xsec =
new TProfile(
"fh1Xsec",
"xsec from pyxsec.root",1,0,1);
366 fh1Xsec->GetXaxis()->SetBinLabel(1,
"<#sigma>");
367 fh1Trials =
new TH1F(
"fh1Trials",
"trials from pyxsec.root",1,0,1);
368 fh1Trials->GetXaxis()->SetBinLabel(1,
"#sum{ntrials}");
369 fh1PtHard =
new TH1F(
"fh1PtHard",
"PYTHIA Pt hard;p_{T,hard}",350,-.5,349.5);
370 fh1PtHardTrials =
new TH1F(
"fh1PtHardTrials",
"PYTHIA Pt hard weight with trials;p_{T,hard}",350,-.5,349.5);
372 fh1EvtsPtHardCut =
new TH1F(
"fh1EvtsPtHardCut",
"#events before and after MC #it{p}_{T,hard} cut;;Events",2,0,2);
377 fh1nRecJetsCuts =
new TH1F(
"fh1nRecJetsCuts",
"reconstructed jets per event",10,-0.5,9.5);
378 fh1nGenJets =
new TH1F(
"fh1nGenJets",
"generated jets per event",10,-0.5,9.5);
405 if(fDebug > 1) Printf(
"AliAnalysisTaskIDFragmentationFunction::UserCreateOutputObjects() -> Loading PID framework");
408 TObjArray* tasks = AliAnalysisManager::GetAnalysisManager()->GetTasks();
410 Printf(
"ERROR loading PID tasks: Failed to retrieve tasks from analysis manager!\n");
428 Printf(
"ERROR: Failed to load inclusive pid task!\n");
434 Printf(
"WARNING: zero inclusive pid tasks!\n");
450 Printf(
"ERROR: Failed to load jet pid task!\n");
456 Printf(
"WARNING: zero jet pid tasks!\n");
472 Printf(
"ERROR: Failed to load jet underlying event pid task!\n");
478 Printf(
"WARNING: zero jet underlying event pid tasks!\n");
485 const Int_t nRefMultBins = 100;
489 const Int_t nJetPtBins = 20;
494 const Double_t binsCentpp[nCentBins+1] = { 0, 0.01, 0.1, 1, 5, 10, 15, 20, 30, 40, 50, 70, 100};
497 "Correlation between jet energy and event multiplicity (|#eta| < 0.5);Ref. mult. (|#eta| < 0.5);#it{p}_{T, jet}^{ch} (GeV/#it{c})",
498 nRefMultBins, refMultDown, refMultUp, nJetPtBins, jetPtDown, jetPtUp);
500 "Correlation between jet energy and event multiplicity (|#eta| < 0.8);Ref. mult. (|#eta| < 0.8);#it{p}_{T, jet}^{ch} (GeV/#it{c})",
501 nRefMultBins, refMultDown, refMultUp, nJetPtBins, jetPtDown, jetPtUp);
503 "Correlation between jet energy and event multiplicity percentile (V0M);Multiplicity Percentile (V0M);#it{p}_{T, jet}^{ch} (GeV/#it{c})",
504 nCentBins, binsCentpp, nJetPtBins, jetPtDown, jetPtUp);
511 Double_t binsPt[nPtBins+1] = {0. , 0.05, 0.1, 0.15, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45,
512 0.5, 0.55, 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95,
513 1.0, 1.1 , 1.2, 1.3 , 1.4, 1.5 , 1.6, 1.7 , 1.8, 1.9 ,
514 2.0, 2.2 , 2.4, 2.6 , 2.8, 3.0 , 3.2, 3.4 , 3.6, 3.8 ,
515 4.0, 4.5 , 5.0, 5.5 , 6.0, 6.5 , 7.0, 8.0 , 9.0, 10.0,
516 11.0, 12.0, 13.0, 14.0, 15.0, 16.0, 18.0, 20.0, 22.0, 24.0,
517 26.0, 28.0, 30.0, 32.0, 34.0, 36.0, 40.0, 45.0, 50.0 };
519 const Int_t DCAbins = 320;
523 fhDCA_XY =
new TH2F(
"fhDCA_XY",
"All rec. tracks;#it{p}_{T} (GeV/#it{c});DCA_{XY}", nPtBins, binsPt, DCAbins, -DCA_XY_max, DCA_XY_max);
524 fhDCA_Z =
new TH2F(
"fhDCA_Z",
"All rec. tracks;#it{p}_{T} (GeV/#it{c});DCA_{Z}", nPtBins, binsPt, DCAbins, -DCA_Z_max, DCA_Z_max);
528 for (
Int_t i = 0; i < AliPID::kSPECIES; i++) {
530 Form(
"Rec. %s (prim.);#it{p}_{T} (GeV/#it{c});DCA_{XY}", AliPID::ParticleLatexName(i)),
531 nPtBins, binsPt, DCAbins, -DCA_XY_max, DCA_XY_max);
533 Form(
"Rec. %s (prim.);#it{p}_{T} (GeV/#it{c});DCA_{Z}", AliPID::ParticleLatexName(i)),
534 nPtBins, binsPt, DCAbins, -DCA_Z_max, DCA_Z_max);
539 Form(
"Rec. %s (sec.);#it{p}_{T} (GeV/#it{c});DCA_{XY}", AliPID::ParticleLatexName(i)),
540 nPtBins, binsPt, DCAbins, -DCA_XY_max, DCA_XY_max);
542 Form(
"Rec. %s (sec.);#it{p}_{T} (GeV/#it{c});DCA_{Z}", AliPID::ParticleLatexName(i)),
543 nPtBins, binsPt, DCAbins, -DCA_Z_max, DCA_Z_max);
555 THnSparse *hnSparse =
dynamic_cast<THnSparse*
>(
fCommonHistList->At(i));
556 if(hnSparse) hnSparse->Sumw2();
560 TH1::AddDirectory(oldStatus);
562 if(fDebug > 2) Printf(
"AliAnalysisTaskIDFragmentationFunction::UserCreateOutputObjects() -> Posting Output");
566 if(fDebug > 2) Printf(
"AliAnalysisTaskIDFragmentationFunction::UserCreateOutputObjects() -> Done");
573 if(fDebug > 1) Printf(
"AliAnalysisTaskIDFragmentationFunction::Init()");
581 if(fDebug > 1) Printf(
"AliAnalysisTaskIDFragmentationFunction::FillHistograms()");
584 if(fDebug > 1) Printf(
"Analysis event #%5d", (
Int_t) fEntry);
586 fMCEvent = MCEvent();
588 if(fDebug>3) Printf(
"%s:%d MCEvent not found in the input", (
char*)__FILE__,__LINE__);
596 Bool_t pythiaGenHeaderFound = kFALSE;
599 AliGenEventHeader* genHeader = fMCEvent->GenEventHeader();
602 AliGenPythiaEventHeader* pythiaGenHeader =
dynamic_cast<AliGenPythiaEventHeader*
>(genHeader);
603 AliGenHijingEventHeader* hijingGenHeader = 0x0;
606 if(fDebug>3) Printf(
"%s:%d pythiaGenHeader found", (
char*)__FILE__,__LINE__);
607 pythiaGenHeaderFound = kTRUE;
608 nTrials = pythiaGenHeader->Trials();
609 ptHard = pythiaGenHeader->GetPtHard();
612 if(fDebug>3) Printf(
"%s:%d no pythiaGenHeader found", (
char*)__FILE__,__LINE__);
614 hijingGenHeader =
dynamic_cast<AliGenHijingEventHeader*
>(genHeader);
615 if(!hijingGenHeader){
616 Printf(
"%s:%d no pythiaGenHeader or hjingGenHeader found", (
char*)__FILE__,__LINE__);
618 if(fDebug>3) Printf(
"%s:%d hijingGenHeader found", (
char*)__FILE__,__LINE__);
652 if (fDebug>3) Printf(
"%s:%d skipping event with pThard %f (>= %f)", (
char*)__FILE__,__LINE__, ptHard,
fMCPtHardCut);
681 AliInputEventHandler* inputHandler = (AliInputEventHandler*)
682 ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
686 if (fDebug > 1 ) Printf(
" Trigger Selection: event REJECTED ... ");
693 if(fDebug>3) Printf(
"%s:%d ESDEvent not found in the input", (
char*)__FILE__,__LINE__);
697 TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
698 if( handler && handler->InheritsFrom(
"AliAODInputHandler") ) {
699 fAOD = ((AliAODInputHandler*)handler)->GetEvent();
701 if (fDebug > 1) Printf(
"%s:%d AOD event from input", (
char*)__FILE__,__LINE__);
704 handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
705 if( handler && handler->InheritsFrom(
"AliAODHandler") ) {
706 fAOD = ((AliAODHandler*)handler)->GetAOD();
708 if (fDebug > 1) Printf(
"%s:%d AOD event from output", (
char*)__FILE__,__LINE__);
713 TObject* outHandler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
714 if( outHandler && outHandler->InheritsFrom(
"AliAODHandler") ) {
715 fAODJets = ((AliAODHandler*)outHandler)->GetAOD();
716 if (fDebug > 1) Printf(
"%s:%d jets from output AOD", (
char*)__FILE__,__LINE__);
723 AliAODHandler *aodH =
dynamic_cast<AliAODHandler*
>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
726 if(fDebug>1)Printf(
"AODExtension not found for %s",
fNonStdFile.Data());
731 Printf(
"%s:%d AODEvent not found", (
char*)__FILE__,__LINE__);
735 Printf(
"%s:%d AODEvent with jet branch not found", (
char*)__FILE__,__LINE__);
742 AliVEvent* evtForCentDetermination = handler->InheritsFrom(
"AliAODInputHandler") ?
fAOD : InputEvent();
748 if(handler->InheritsFrom(
"AliAODInputHandler")){
750 centPercent = ((AliAODHeader*)
fAOD->GetHeader())->GetCentrality();
752 if(centPercent>10) cl = 2;
753 if(centPercent>30) cl = 3;
754 if(centPercent>50) cl = 4;
763 if (fDebug > 1) Printf(
"%s:%d event not in selected event class: event REJECTED ...",(
char*)__FILE__,__LINE__);
774 if (!
fIsPP) centPercent = evtForCentDetermination->GetCentrality()->GetCentralityPercentile(
fCentralityEstimator.Data());
777 AliPIDResponse *fPIDResponse = inputHandler->GetPIDResponse();
779 AliError(
"PIDResponse object was not created");
782 fPIDResponse->SetCurrentCentrality(centPercent);
785 const Int_t refMult5 = ((AliAODHeader*)
fAOD->GetHeader())->GetRefMultiplicityComb05();
786 const Int_t refMult8 = ((AliAODHeader*)
fAOD->GetHeader())->GetRefMultiplicityComb08();
810 AliAODVertex* primVtx =
fAOD->GetPrimaryVertex();
812 Printf(
"%s:%d Primary vertex not found", (
char*)__FILE__,__LINE__);
816 Int_t nTracksPrim = primVtx->GetNContributors();
820 if (fDebug > 1) Printf(
"%s:%d primary vertex selection: %d", (
char*)__FILE__,__LINE__,nTracksPrim);
821 if(nTracksPrim <= 0) {
822 if (fDebug > 1) Printf(
"%s:%d primary vertex selection: event REJECTED...",(
char*)__FILE__,__LINE__);
828 TString primVtxName(primVtx->GetName());
830 if(primVtxName.CompareTo(
"TPCVertex",TString::kIgnoreCase) == 1){
831 if (fDebug > 1) Printf(
"%s:%d primary vertex selection: TPC vertex, event REJECTED...",(
char*)__FILE__,__LINE__);
860 if (fDebug > 1) Printf(
"%s:%d primary vertex z = %f: event REJECTED...",(
char*)__FILE__,__LINE__,primVtx->GetZ());
890 Bool_t isPileUpInclusivePIDtask[arrSizeInclusive];
891 Bool_t isPileUpJetPIDtask[arrSizeJet];
892 Bool_t isPileUpJetUEPIDtask[arrSizeJetUE];
894 for (
Int_t i = 0; i < arrSizeInclusive; i++)
895 isPileUpInclusivePIDtask[i] = kFALSE;
897 for (
Int_t i = 0; i < arrSizeJet; i++)
898 isPileUpJetPIDtask[i] = kFALSE;
900 for (
Int_t i = 0; i < arrSizeJetUE; i++)
901 isPileUpJetUEPIDtask[i] = kFALSE;
904 Bool_t isPileUpForAllInclusivePIDTasks = kTRUE;
905 Bool_t isPileUpForAllJetPIDTasks = kTRUE;
906 Bool_t isPileUpForAllJetUEPIDTasks = kTRUE;
913 if (!isPileUpInclusivePIDtask[i])
916 isPileUpForAllInclusivePIDTasks = isPileUpForAllInclusivePIDTasks && isPileUpInclusivePIDtask[i];
923 if (!isPileUpJetPIDtask[i])
926 isPileUpForAllJetPIDTasks = isPileUpForAllJetPIDTasks && isPileUpJetPIDtask[i];
933 if (!isPileUpJetUEPIDtask[i])
936 isPileUpForAllJetUEPIDTasks = isPileUpForAllJetUEPIDTasks && isPileUpJetUEPIDtask[i];
942 if (fDebug > 1) Printf(
"%s:%d event ACCEPTED ...",(
char*)__FILE__,__LINE__);
955 if (!isPileUpInclusivePIDtask[i])
961 if (!isPileUpJetPIDtask[i])
967 if (!isPileUpJetUEPIDtask[i])
978 if (!isPileUpInclusivePIDtask[i])
985 if (!isPileUpJetPIDtask[i])
992 if (!isPileUpJetUEPIDtask[i])
1002 AliPIDResponse* pidResponse = 0x0;
1003 Bool_t tuneOnDataTPC = kFALSE;
1005 if (!inputHandler) {
1006 AliFatal(
"Input handler needed");
1011 pidResponse = inputHandler->GetPIDResponse();
1013 AliFatal(
"PIDResponse object was not created");
1017 tuneOnDataTPC = pidResponse->IsTunedOnData() &&
1018 ((pidResponse->GetTunedOnDataMask() & AliPIDResponse::kDetTPC) == AliPIDResponse::kDetTPC);
1023 if(fDebug>2)Printf(
"%s:%d Starting processing...",(
char*)__FILE__,__LINE__);
1031 if (trackContainer) {
1038 if(fDebug>2)Printf(
"%s:%d Starting Inclusive Efficiency...",(
char*)__FILE__,__LINE__);
1050 if (
part->Eta() > trackEtaMax ||
part->Eta() < trackEtaMin)
1060 if (!
part->IsPhysicalPrimary())
1073 if (TMath::Abs(chargeMC) < 0.01)
1077 -1, -1, -1, -1, -1, -1 };
1080 if (!isPileUpInclusivePIDtask[i] &&
fInclusivePIDtask[i]->IsInAcceptedEtaRange(TMath::Abs(
part->Eta()))) {
1082 fInclusivePIDtask[i]->FillGeneratedYield(valuesGenYield);
1088 centPercent, -1, -1, -1, -1, -1 };
1090 if (!isPileUpInclusivePIDtask[i])
1099 for(
auto inclusiveaod : trackContainerEfficiency->
accepted()) {
1103 Double_t dEdxTPC = tuneOnDataTPC ? pidResponse->GetTPCsignalTunedOnData(inclusiveaod)
1104 : inclusiveaod->GetTPCsignal();
1112 Int_t label = TMath::Abs(inclusiveaod->GetLabel());
1119 pdg = gentrack->GetPdgCode();
1134 if (!gentrack->IsPhysicalPrimary())
1142 if (gentrack->Eta() > trackEtaMax || gentrack->Eta() < trackEtaMin)
1148 gentrack->Charge() / 3., centPercent, -1, -1, -1, -1, -1 };
1150 if (!isPileUpInclusivePIDtask[i] &&
1158 static_cast<Double_t>(inclusiveaod->Charge()), centPercent,
1159 -1, -1, -1, -1, -1 };
1161 if (!isPileUpInclusivePIDtask[i] &&
1172 if(fDebug>2)Printf(
"%s:%d Process inclusive tracks...",(
char*)__FILE__,__LINE__);
1178 Double_t dEdxTPC = tuneOnDataTPC ? pidResponse->GetTPCsignalTunedOnData(
part)
1179 :
part->GetTPCsignal();
1187 Int_t label = TMath::Abs(
part->GetLabel());
1190 AliAODMCParticle* gentrack = mcParticleContainer ? mcParticleContainer->
GetMCParticleWithLabel(label) : 0x0;
1194 pdg = gentrack->GetPdgCode();
1197 if (!isPileUpInclusivePIDtask[i] &&
1211 -1, -1, -1, -1, -1 };
1213 if (!isPileUpInclusivePIDtask[i] &&
1223 if (!isPileUpInclusivePIDtask[i] &&
1232 if (gentrack->IsPhysicalPrimary()) {
1236 gentrack->Charge() / 3., centPercent, -1, -1, -1, -1, -1 };
1239 gentrack->Charge() / 3., centPercent };
1242 if (!isPileUpInclusivePIDtask[i] &&
1259 if(fDebug>2)Printf(
"%s:%d Process Jets - efficiency, particle level..",(
char*)__FILE__,__LINE__);
1262 mcJetContainer->SortArray();
1264 if (
fUseJetPIDtask && mcJetContainer && !isPileUpForAllJetPIDTasks) {
1267 if (!isPileUpJetPIDtask[i]) {
1269 for(
auto jet : mcJetContainer->
accepted()) {
1276 jetPt = jetPt - mcJetContainer->
GetRhoVal() * jet->Area();
1280 TClonesArray *arrayWithTracks = mcParticleContainer ? mcParticleContainer->GetArray() : 0x0;
1282 for(
Int_t it=0; it<jet->GetNumberOfTracks(); ++it) {
1284 AliVParticle* trackVP =
dynamic_cast<AliVParticle*
>(jet->TrackAt(it,arrayWithTracks));
1295 if(fDebug>2)Printf(
"%s:%d Process jets...",(
char*)__FILE__,__LINE__);
1298 jetContainer->SortArray();
1300 if (
fUseJetPIDtask && jetContainer && !isPileUpForAllJetPIDTasks) {
1302 for (
auto jet : jetContainer->
accepted()) {
1308 jetPt = jetPt - jetContainer->
GetRhoVal() * jet->Area();
1311 if (!isPileUpJetPIDtask[i])
1319 TClonesArray *arrayWithTracks = trackContainer->GetArray();
1321 for(
Int_t j=0; j<jet->GetNumberOfTracks(); ++j) {
1323 AliVTrack *track =
dynamic_cast<AliVTrack*
>(jet->TrackAt(j,arrayWithTracks));
1324 Bool_t trackRejectedByTask[arrSizeJet];
1325 Bool_t trackRejectedByAllTasks = kTRUE;
1329 Double_t dEdxTPC = tuneOnDataTPC ? pidResponse->GetTPCsignalTunedOnData(track) : track->GetTPCsignal();
1332 trackRejectedByAllTasks = trackRejectedByAllTasks && trackRejectedByTask[i];
1335 if (!trackRejectedByAllTasks) {
1336 AnalyseJetTrack(track, jet, 0x0, trackRejectedByTask, centPercent, mcParticleContainer);
1337 FillDCA(track, mcParticleContainer);
1347 if(fDebug>2)Printf(
"%s:%d Process Underlying event...",(
char*)__FILE__,__LINE__);
1350 if (!isPileUpJetUEPIDtask[i]) {
1354 TList *jetUElist = 0x0;
1355 TList* mcJetUElist = 0x0;
1357 if (mcJetContainer) {
1359 mcJetContainer->SortArray();
1361 if (method.Contains(
"RC",TString::kIgnoreCase) || method.Contains(
"Random",TString::kIgnoreCase)) {
1364 else if (method.Contains(
"PC",TString::kIgnoreCase) || method.Contains(
"Perpendicular",TString::kIgnoreCase)) {
1372 jetContainer->SortArray();
1374 if (method.Contains(
"RC",TString::kIgnoreCase) || method.Contains(
"Random",TString::kIgnoreCase)) {
1378 else if (method.Contains(
"PC",TString::kIgnoreCase) || method.Contains(
"Perpendicular",TString::kIgnoreCase)) {
1384 for (
Int_t i=0;i<jetUElist->GetEntries();++i) {
1394 if (!jetUEtracklist) {
1395 cout <<
"No track list for an underlying event jet. Check Code!" << endl;
1399 jetUEtracklist->SetOwner(kFALSE);
1401 for (
Int_t j=0;j<jetUEtracklist->GetEntries();++j) {
1402 AliVTrack* UEtrack = (AliVTrack*)(jetUEtracklist->At(j));
1407 Double_t dEdxTPC = tuneOnDataTPC ? pidResponse->GetTPCsignalTunedOnData(UEtrack) : UEtrack->GetTPCsignal();
1410 AnalyseJetTrack(UEtrack, jet, task, 0x0, centPercent, mcParticleContainer);
1412 UEPtDensity += UEtrack->Pt();
1415 delete jetUEtracklist;
1416 jetUEtracklist = 0x0;
1424 for (
Int_t i=0;i<mcJetUElist->GetEntries();++i) {
1429 if (!mcJetUEtracklist) {
1430 cout <<
"No track list for MC Tracks in Underlying event. Check Code!" << endl;
1433 mcJetUEtracklist->SetOwner(kFALSE);
1434 for (
Int_t j=0;j<mcJetUEtracklist->GetEntries();++j) {
1435 AliVParticle* mcUEParticle =
dynamic_cast<AliVParticle*
>(mcJetUEtracklist->At(j));
1447 if(fDebug > 2) Printf(
"%s:%d Processing done, posting output data now...",(
char*)__FILE__,__LINE__);
1470 if(fDebug > 2) Printf(
"%s:%d Done",(
char*)__FILE__,__LINE__);
1485 if(fDebug > 1) printf(
"AliAnalysisTaskIDFragmentationFunction::Terminate() \n");
1493 for(
Int_t i=0; i<dim; i++){
1494 h->GetAxis(i)->SetTitle(labels[i]);
1495 h->GetAxis(i)->SetTitleColor(1);
1506 h->GetXaxis()->SetTitleColor(1);
1507 h->GetYaxis()->SetTitleColor(1);
1518 h->GetXaxis()->SetTitleColor(1);
1519 h->GetYaxis()->SetTitleColor(1);
1520 h->GetZaxis()->SetTitleColor(1);
1526 Bool_t useInternalFlag)
const 1537 TVector3 v_track(track->Px(), track->Py(), track->Pz());
1538 TVector3 v_jet(jet->
Px(), jet->
Py(), jet->
Pz());
1540 return v_track.DeltaR(v_jet);
1551 if(0.150<pt && pt<0.200) alpha = 3.639;
1552 if(0.200<pt && pt<0.250) alpha = 2.097;
1553 if(0.250<pt && pt<0.300) alpha = 1.930;
1554 if(0.300<pt && pt<0.350) alpha = 1.932;
1555 if(0.350<pt && pt<0.400) alpha = 1.943;
1556 if(0.400<pt && pt<0.450) alpha = 1.993;
1557 if(0.450<pt && pt<0.500) alpha = 1.989;
1558 if(0.500<pt && pt<0.600) alpha = 1.963;
1559 if(0.600<pt && pt<0.700) alpha = 1.917;
1560 if(0.700<pt && pt<0.800) alpha = 1.861;
1561 if(0.800<pt && pt<0.900) alpha = 1.820;
1562 if(0.900<pt && pt<1.000) alpha = 1.741;
1563 if(1.000<pt && pt<1.500) alpha = 0.878;
1574 if(!mcParticleContainer)
return 1;
1576 AliAODMCParticle* currentMother = daughter;
1577 AliAODMCParticle* currentDaughter = daughter;
1583 Int_t daughterPDG = currentDaughter->GetPdgCode();
1585 Int_t motherLabel = currentDaughter->GetMother();
1590 currentMother = currentDaughter;
1594 Int_t motherPDG = currentMother->GetPdgCode();
1597 if(currentMother->IsPhysicalPrimary())
break;
1599 if(TMath::Abs(daughterPDG) == 321){
1600 currentMother = currentDaughter;
break;
1602 if(TMath::Abs(motherPDG) == 310 ){
1605 if(TMath::Abs(motherPDG) == 3212 && TMath::Abs(daughterPDG) == 3122){
1606 currentMother = currentDaughter;
break;
1609 currentDaughter = currentMother;
1613 Int_t motherPDG = currentMother->GetPdgCode();
1614 Double_t motherGenPt = currentMother->Pt();
1622 Bool_t useInternalFlag)
const 1633 TVector3 v_track(track->Px(), track->Py(), track->Pz());
1634 TVector3 v_jet(jet->
Px(), jet->
Py(), jet->
Pz());
1637 return v_track.Perp(v_jet);
1647 if (!mcParticleContainer || !part)
1650 if (part->IsPhysicalPrimary())
1653 Int_t iMother = part->GetMother();
1658 AliAODMCParticle* partM =
dynamic_cast<AliAODMCParticle*
>(mcParticleContainer->
GetMCParticleWithLabel(iMother));
1662 Int_t codeM = TMath::Abs(partM->GetPdgCode());
1663 Int_t mfl =
Int_t(codeM / TMath::Power(10,
Int_t(TMath::Log10(codeM))));
1664 if (mfl == 3 && codeM != 3)
1672 jetUElist->SetOwner(kTRUE);
1675 for (
auto jet : jetContainer->
accepted()) {
1676 jetRC =
GetRandomCone(jet, maxEtaTrack - coneRadius, 2.0 * coneRadius);
1678 jetUElist->Add(jetRC);
1689 jetUElist->SetOwner(kTRUE);
1691 Double_t perpAngle = TMath::Pi()/2.0;
1692 for (
auto jet : jetContainer->
accepted()) {
1696 jetUElist->Add(jetPC);
1700 jetUElist->Add(jetPC);
1713 TLorentzVector vecRdCone;
1723 dEta = dEtaConeMax * (2 *
fRandom->Rndm() - 1.0);
1724 dPhi = TMath::TwoPi() *
fRandom->Rndm();
1760 if(!part)
return kFALSE;
1761 if(!dDistance)
return kFALSE;
1765 for(
auto jet : jetContainer->
accepted()){
1767 if (fDebug>2) std::cout <<
"AliAnalysisTaskJetChem::IsRCJCOverlap jet pointer invalid!" << std::endl;
1778 if (!part1 || !part2)
1781 TVector3 vecMom2(part2->Px(),part2->Py(),part2->Pz());
1782 TVector3 vecMom1(part1->Px(),part1->Py(),part1->Pz());
1783 Double_t dR = vecMom2.DeltaR(vecMom1);
1791 jetTrackList->SetOwner(kFALSE);
1793 return jetTrackList;
1795 if (!particleContainer)
1798 for (
auto track : particleContainer->
accepted()) {
1804 jetTrackList->Add(track);
1807 return jetTrackList;
1813 if (!jet || !trackVP || !task)
1816 if (!mcJetContainer)
1819 if (!mcJetContainer) {
1820 std::cout <<
"Monte-Carlo Jet Container not found." << std::endl;
1826 Float_t trackPt = trackVP->Pt();
1829 AliAODMCParticle*
part =
dynamic_cast<AliAODMCParticle*
>(trackVP);
1831 AliError(
"expected ref track not found ");
1839 if (mcParticleContainer) {
1845 if (part->Eta() > trackEtaMax || part->Eta() < trackEtaMin)
1855 if (!part->IsPhysicalPrimary())
1869 Double_t chargeMC = part->Charge() / 3.;
1871 if (TMath::Abs(chargeMC) < 0.01)
1875 chargeMC, distance, jT };
1884 centPercent, jetPt, z, xi, distance, jT };
1890 if (!track || (!task && !trackRejectedByTask)) {
1891 std::cout <<
"Cannot analyse track" << std::endl;
1892 std::cout <<
"Track: " << track << std::endl;
1893 std::cout <<
"Task: " << task << std::endl;
1894 std::cout <<
"Rejection Array: " << trackRejectedByTask << std::endl;
1906 Int_t label = TMath::Abs(track->GetLabel());
1909 Int_t mcID = AliPID::kUnknown;
1916 if (!mcParticleContainer)
1919 AliAODMCParticle* gentrack = mcParticleContainer ? mcParticleContainer->
GetMCParticleWithLabel(label) : 0x0;
1923 pdg = gentrack->GetPdgCode();
1930 if (!trackRejectedByTask[i])
1942 static_cast<Double_t>(track->Charge()),
1943 centPercent, jetPt, z, xi, distance, jT };
1955 if (!trackRejectedByTask[i])
1962 if (gentrack->IsPhysicalPrimary()) {
1971 gentrack->Charge() / 3., centPercent, jetPt, genZ,
1972 genXi, genDistance, genJt };
1983 if (!trackRejectedByTask[i])
2001 if (mcParticleContainer) {
2006 if (gentrack->Eta() <= trackEtaMax || gentrack->Eta() >= trackEtaMin) {
2014 Double_t measZ = -1., measXi = -1.;
2024 static_cast<Double_t>(track->Charge()),
2025 centPercent, jetPt, measZ, measXi, measDistance, measJt };
2033 if (!trackRejectedByTask[i])
2051 Int_t label = TMath::Abs(track->GetLabel());
2053 if (!mcParticleContainer)
2056 Int_t mcID = AliPID::kUnknown;
2058 AliAODMCParticle* gentrack = mcParticleContainer ? mcParticleContainer->
GetMCParticleWithLabel(label) : 0x0;
2061 Int_t pdg = gentrack->GetPdgCode();
2069 AliAODVertex* primVtx =
fAOD->GetPrimaryVertex();
2071 std::cout <<
"Primary Vertex not found, do not fill DCA" << std::endl;
2076 dca[0] = pos[0] - v[0];
2077 dca[1] = pos[1] - v[1];
2081 fhDCA_Z->Fill(trackPt, dca[1]);
2084 if (gentrack && mcID != AliPID::kUnknown) {
2086 if (gentrack->IsPhysicalPrimary()) {
virtual TString GetNameJetContainer() const
TString fNonStdFile
where we take the jets from can be input or output AOD
TH1F * fh1VertexZ
NContributors to prim vertex.
TH1F * fh1EvtCent
number of reconstructed tracks after cuts
Double_t GetRhoVal() const
void NormalizeJetArea(Double_t jetParameter)
static void SetProperties(TH1 *h, const char *x, const char *y)
TH1F * fh1PtHardTrials
pt hard of the event
TH1F * fh1nGenJets
number of jets from reconstructed tracks per event
virtual TString GetNameMCParticleContainer() const
Bool_t FillEfficiencyContainer(const Double_t *values, EffSteps step, Double_t weight=1.0)
AliAODEvent * fAOD
ESD event.
TH1F * fh1VtxSelection
event cuts
AliJetContainer * GetJetContainer(Int_t i=0) const
TH1F * fh1PtHard
sum of trials
virtual void AnalyseJetTrack(AliVTrack *track, AliEmcalJet *jet, AliAnalysisTaskPID *task, const Bool_t *trackRejectedByTask, Double_t centPercent, AliMCParticleContainer *mcParticleContainer=0x0)
static Bool_t TPCCutMIGeo(const AliVTrack *track, const AliVEvent *evt, TTreeStream *streamer=0x0)
Bool_t FillGenJets(Double_t centralityPercentile, Double_t jetPt, Double_t norm=-1.)
TString * fNameJetPIDtask
TH2F * fhDCA_XY_sec_MCID[AliPID::kSPECIES]
DCA Z for all rec. prim. particles sorted by MC-ID.
const AliTrackIterableContainer accepted() const
Double_t GetEtaAbsCutUp() const
Container with name, TClonesArray and cuts for particles.
AliAnalysisTaskPID ** fJetUEPIDtask
Pointer to tasks for jet PID spectra.
virtual TString GetNameMCParticleJetContainer() const
Int_t fNumInclusivePIDtasks
Object to use analysis utils like pile-up rejection.
TH1F * fh1EvtsPtHardCut
pt hard of the event
static Double_t GetMCStrangenessFactorCMS(Int_t motherPDG, Double_t motherGenPt)
virtual Bool_t GetUseTPCCutMIGeo() const
TH1F * fh1VertexNContributors
type of accepted vertices
AliAnalysisTaskPID ** fInclusivePIDtask
virtual Bool_t GetUseTPCnclCut() const
Bool_t IsInAcceptedEtaRange(Double_t etaAbs) const
virtual Bool_t GetIsPileUp(AliVEvent *event, PileUpRejectionType pileUpRejection=kPileUpRejectionClass) const
Double_t GetDistanceJetTrack(const AliEmcalJet *jet, const AliVParticle *track, Bool_t useInternalFlag=kTRUE) const
Double_t GetParticleEtaMin() const
virtual void UserCreateOutputObjects()
Bool_t IsSecondaryWithStrangeMotherMC(AliAODMCParticle *part, AliMCParticleContainer *mcParticleContainer)
TString part
use mixed event to constrain combinatorial background
Container for particles within the EMCAL framework.
virtual TList * GetUEJetsWithRandomConeMethod(AliJetContainer *jetContainer, Double_t coneRadius, Double_t maxEtaTrack)
virtual void Terminate(Option_t *)
AliAnalysisUtils * fAnaUtils
TProfile * fh1Xsec
centrality percentile
TString * fNameJetUEPIDtask
Bool_t FillPythiaTrials(Double_t avgTrials)
Bool_t FillGeneratedYield(const Double_t *values, Double_t weight=1.0)
Double_t GetMCStrangenessFactorCMS(AliAODMCParticle *daughter, AliMCParticleContainer *mcParticleContainer) const
virtual Bool_t IsRCJCOverlap(const AliVParticle *part, Double_t dDistance) const
static void GetJetTrackObservables(Double_t trackPt, Double_t jetPt, Double_t &z, Double_t &xi, Bool_t storeXi=kTRUE)
static Int_t EventClass(Bool_t bSet=kFALSE, Int_t iNew=0)
TH2F * fhDCA_Z_sec_MCID[AliPID::kSPECIES]
DCA XY for all rec. sec. particles sorted by MC-ID.
TH2F * fhDCA_Z
DCA XY for all rec. particles.
AliAODEvent * fAODJets
AOD event.
Bool_t GetStoreCharge() const
Bool_t ProcessTrack(const AliVTrack *track, Int_t particlePDGcode, Double_t centralityPercentile, Double_t jetPt, Bool_t isMBSelected=kFALSE, Bool_t isMultSelected=kTRUE, Bool_t storeXi=kTRUE, Double_t radialDistanceToJet=-1, Double_t jT=-1)
void FillUEDensity(Double_t cent, Double_t UEpt)
Bool_t fUseInclusivePIDtask
Pointer to tasks for jet UE PID spectra.
TRandom3 * fRandom
DCA Z for all rec. sec. particles sorted by MC-ID.
Int_t GetIndexOfChargeAxisGenYield() const
virtual void ConfigureTaskForCurrentEvent(AliVEvent *event)
AliAnalysisTaskIDFragmentationFunction()
Float_t GetFFRadius() const
TH1F * fh1EvtMult
prim vertex z distribution
TH1F * fh1Trials
pythia cross section and trials
static Bool_t TPCnclCut(const AliVTrack *track)
Double_t GetPerpendicularMomentumTrackJet(const AliEmcalJet *jet, const AliVParticle *track, Bool_t useInternalFlag=kTRUE) const
AliAnalysisTaskPID ** fJetPIDtask
Pointer to tasks for inclusive PID spectra.
virtual void FillDCA(AliVTrack *track, AliMCParticleContainer *mcParticleContainer)
Bool_t FillHistograms()
Function filling histograms.
TH2F * fhJetPtRefMultEta8
Jet pT vs. reference multiplicity (|eta|<0.5)
Bool_t fStudyTransversalJetStructure
AliRhoParameter * GetRhoParameter()
TH2F * fhJetPtRefMultEta5
DCA Z for all rec. particles.
TString * fNameInclusivePIDtask
virtual void PerformJetMonteCarloAnalysisGeneratedYield(AliEmcalJet *jet, AliVParticle *trackVP, AliAnalysisTaskPID *task, Double_t centPercent, AliJetContainer *mcJetContainer=0x0)
virtual AliEmcalJet * GetPerpendicularCone(AliEmcalJet *processedJet, Double_t perpAngle) const
TString fCentralityEstimator
TH2F * fhJetPtMultPercent
Jet pT vs. reference multiplicity (|eta|<0.8)
AliMCParticleContainer * GetMCParticleContainer(Int_t i=0) const
TH2F * fhDCA_XY_prim_MCID[AliPID::kSPECIES]
Jet pT vs. multiplicity percentile (usually V0M)
virtual TList * GetUEJetsWithPerpendicularConeMethod(AliJetContainer *jetContainer)
TH1F * fh1nRecJetsCuts
Number events before and after the cut on MC pT hard.
Double_t GetParticleEtaMax() const
AliTrackContainer * GetTrackContainer(Int_t i=0) const
AliAODExtension * fAODExtension
AOD event with jet branch (case we have AOD both in input and output)
Base task in the EMCAL jet framework.
Bool_t IsParticleInCone(const AliVParticle *part1, const AliVParticle *part2, Double_t dRMax) const
Represent a jet reconstructed using the EMCal jet framework.
static Int_t PDGtoMCID(Int_t pdg)
static Bool_t PythiaInfoFromFile(const char *currFile, Float_t &fXsec, Float_t &fTrials)
void FillJetArea(Double_t cent, Double_t area)
virtual AliEmcalJet * GetRandomCone(AliEmcalJet *processedJet, Double_t dEtaConeMax, Double_t dDistance) const
const AliParticleIterableContainer accepted() const
Double_t GetMCStrangenessFactor(Double_t pt) const
virtual TString GetNameTrackContainer() const
virtual TList * GetTracksInCone(const AliEmcalJet *jet, AliParticleContainer *particleContainer=0x0) const
const AliJetIterableContainer accepted() const
Bool_t fOnlyLeadingJets
TRandom3 for background estimation.
TH2F * fhDCA_XY
number of jets from generated tracks per event
TH2F * fhDCA_Z_prim_MCID[AliPID::kSPECIES]
DCA XY for all rec. prim. particles sorted by MC-ID.
Container for MC-true particles within the EMCAL framework.
virtual AliAODMCParticle * GetMCParticleWithLabel(Int_t lab) const
Bool_t FillRecJets(Double_t centralityPercentile, Double_t jetPt, Double_t norm=-1.)
const AliMCParticleIterableContainer accepted() const
Container for jet within the EMCAL jet framework.
Bool_t FillCutHisto(Double_t value, CutHistoType type)
virtual TString GetNameTrackContainerEfficiency() const
Bool_t IncrementEventCounter(Double_t centralityPercentile, EventCounterType type)
TList * OpenFile(const char *fname)
virtual ~AliAnalysisTaskIDFragmentationFunction()
Bool_t FillPtResolution(Int_t mcID, const Double_t *values)