1 // *************************************************************************
3 // * Task for ID Fragmentation Function Analysis *
5 // *************************************************************************
8 /**************************************************************************
9 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
11 * Author: The ALICE Off-line Project. *
12 * Contributors are mentioned in the code where appropriate. *
14 * Permission to use, copy, modify and distribute this software and its *
15 * documentation strictly for non-commercial purposes is hereby granted *
16 * without fee, provided that the above copyright notice appears in all *
17 * copies and that both the copyright notice and this permission notice *
18 * appear in the supporting documentation. The authors make no claims *
19 * about the suitability of this software for any purpose. It is *
20 * provided "as is" without express or implied warranty. *
21 **************************************************************************/
30 #include "THnSparse.h"
37 #include "AliAODInputHandler.h"
38 #include "AliAODHandler.h"
39 #include "AliESDEvent.h"
40 #include "AliAODMCParticle.h"
41 #include "AliAODJet.h"
42 #include "AliAODJetEventBackground.h"
43 #include "AliGenPythiaEventHeader.h"
44 #include "AliGenHijingEventHeader.h"
45 #include "AliInputEventHandler.h"
47 #include "AliAnalysisHelperJetTasks.h"
48 #include "AliAnalysisManager.h"
49 #include "AliAnalysisTaskSE.h"
50 #include "AliAnalysisUtils.h"
51 #include "AliVParticle.h"
52 #include "AliVEvent.h"
54 #include "AliAnalysisTaskPID.h"
55 #include "AliPIDResponse.h"
57 #include "AliAnalysisTaskIDFragmentationFunction.h"
62 ClassImp(AliAnalysisTaskIDFragmentationFunction)
64 //____________________________________________________________________________
65 AliAnalysisTaskIDFragmentationFunction::AliAnalysisTaskIDFragmentationFunction()
72 ,fBranchRecJets("jets")
73 ,fBranchRecBckgClusters("")
75 ,fBranchEmbeddedJets("")
79 ,fUseAODInputJets(kTRUE)
81 ,fUsePhysicsSelection(kTRUE)
91 ,fUseExtraTracksBgr(0)
92 ,fCutFractionPtEmbedded(0)
93 ,fUseEmbeddedJetAxis(0)
113 ,fTracksRecCutsEfficiency(0)
115 ,fTracksAODMCCharged(0)
116 ,fTracksAODMCChargedSecNS(0)
117 ,fTracksAODMCChargedSecS(0)
118 ,fTracksRecQualityCuts(0)
127 ,fQATrackHistosRecCuts(0)
128 ,fQATrackHistosGen(0)
130 ,fQAJetHistosRecCuts(0)
131 ,fQAJetHistosRecCutsLeading(0)
133 ,fQAJetHistosGenLeading(0)
134 ,fQAJetHistosRecEffLeading(0)
136 ,fFFHistosRecCutsInc(0)
137 ,fFFHistosRecLeadingTrack(0)
140 ,fFFHistosGenLeadingTrack(0)
141 ,fQATrackHighPtThreshold(0)
175 ,fh1VertexNContributors(0)
187 ,fh1nRecBckgJetsCuts(0)
189 ,fh2PtRecVsGenPrim(0)
193 ,fhJetPtRefMultEta5(0)
194 ,fhJetPtRefMultEta8(0)
195 ,fhJetPtMultPercent(0)
196 ,fQATrackHistosRecEffGen(0)
197 ,fQATrackHistosRecEffRec(0)
198 ,fQATrackHistosSecRecNS(0)
199 ,fQATrackHistosSecRecS(0)
200 ,fQATrackHistosSecRecSsc(0)
201 ,fFFHistosRecEffRec(0)
202 ,fFFHistosSecRecNS(0)
204 ,fFFHistosSecRecSsc(0)
211 ,fh1FractionPtEmbedded(0)
213 ,fh2DeltaPtVsJetPtEmbedded(0)
214 ,fh2DeltaPtVsRecJetPtEmbedded(0)
215 ,fh1DeltaREmbedded(0)
216 ,fQABckgHisto0RecCuts(0)
218 ,fQABckgHisto1RecCuts(0)
220 ,fQABckgHisto2RecCuts(0)
222 ,fQABckgHisto3RecCuts(0)
224 ,fQABckgHisto4RecCuts(0)
226 ,fFFBckgHisto0RecCuts(0)
228 ,fFFBckgHisto1RecCuts(0)
230 ,fFFBckgHisto2RecCuts(0)
232 ,fFFBckgHisto3RecCuts(0)
234 ,fFFBckgHisto4RecCuts(0)
236 ,fFFBckgHisto0RecEffRec(0)
237 ,fFFBckgHisto0SecRecNS(0)
238 ,fFFBckgHisto0SecRecS(0)
239 ,fFFBckgHisto0SecRecSsc(0)
241 ,fProNtracksLeadingJet(0)
243 ,fProNtracksLeadingJetGen(0)
244 ,fProDelR80pcPtGen(0)
245 ,fProNtracksLeadingJetBgrPerp2(0)
246 ,fProNtracksLeadingJetRecPrim(0)
247 ,fProDelR80pcPtRecPrim(0)
248 ,fProNtracksLeadingJetRecSecNS(0)
249 ,fProNtracksLeadingJetRecSecS(0)
250 ,fProNtracksLeadingJetRecSecSsc(0)
254 ,fOnlyLeadingJets(kFALSE)
259 ,fNumInclusivePIDtasks(0)
261 ,fNameInclusivePIDtask(0x0)
262 ,fNameJetPIDtask(0x0)
263 ,fInclusivePIDtask(0x0)
265 ,fUseInclusivePIDtask(kFALSE)
266 ,fUseJetPIDtask(kFALSE)
269 // default constructor
276 for(Int_t ii=0; ii<5; ii++){
277 fProDelRPtSum[ii] = 0;
278 fProDelRPtSumGen[ii] = 0;
279 fProDelRPtSumBgrPerp2[ii] = 0;
280 fProDelRPtSumRecPrim[ii] = 0;
281 fProDelRPtSumRecSecNS[ii] = 0;
282 fProDelRPtSumRecSecS[ii] = 0;
283 fProDelRPtSumRecSecSsc[ii] = 0;
286 for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
287 fIDFFHistosRecCuts[i] = 0x0;
288 fIDFFHistosGen[i] = 0x0;
290 fhDCA_XY_prim_MCID[i] = 0x0;
291 fhDCA_Z_prim_MCID[i] = 0x0;
293 fhDCA_XY_sec_MCID[i] = 0x0;
294 fhDCA_Z_sec_MCID[i] = 0x0;
298 //_______________________________________________________________________________________________
299 AliAnalysisTaskIDFragmentationFunction::AliAnalysisTaskIDFragmentationFunction(const char *name)
300 : AliAnalysisTaskSE(name)
306 ,fBranchRecJets("jets")
307 ,fBranchRecBckgClusters("")
309 ,fBranchEmbeddedJets("")
313 ,fUseAODInputJets(kTRUE)
315 ,fUsePhysicsSelection(kTRUE)
316 ,fEvtSelectionMask(0)
325 ,fUseExtraTracksBgr(0)
326 ,fCutFractionPtEmbedded(0)
327 ,fUseEmbeddedJetAxis(0)
328 ,fUseEmbeddedJetPt(0)
347 ,fTracksRecCutsEfficiency(0)
349 ,fTracksAODMCCharged(0)
350 ,fTracksAODMCChargedSecNS(0)
351 ,fTracksAODMCChargedSecS(0)
352 ,fTracksRecQualityCuts(0)
361 ,fQATrackHistosRecCuts(0)
362 ,fQATrackHistosGen(0)
364 ,fQAJetHistosRecCuts(0)
365 ,fQAJetHistosRecCutsLeading(0)
367 ,fQAJetHistosGenLeading(0)
368 ,fQAJetHistosRecEffLeading(0)
370 ,fFFHistosRecCutsInc(0)
371 ,fFFHistosRecLeadingTrack(0)
374 ,fFFHistosGenLeadingTrack(0)
375 ,fQATrackHighPtThreshold(0)
409 ,fh1VertexNContributors(0)
421 ,fh1nRecBckgJetsCuts(0)
423 ,fh2PtRecVsGenPrim(0)
427 ,fhJetPtRefMultEta5(0)
428 ,fhJetPtRefMultEta8(0)
429 ,fhJetPtMultPercent(0)
430 ,fQATrackHistosRecEffGen(0)
431 ,fQATrackHistosRecEffRec(0)
432 ,fQATrackHistosSecRecNS(0)
433 ,fQATrackHistosSecRecS(0)
434 ,fQATrackHistosSecRecSsc(0)
435 ,fFFHistosRecEffRec(0)
436 ,fFFHistosSecRecNS(0)
438 ,fFFHistosSecRecSsc(0)
445 ,fh1FractionPtEmbedded(0)
447 ,fh2DeltaPtVsJetPtEmbedded(0)
448 ,fh2DeltaPtVsRecJetPtEmbedded(0)
449 ,fh1DeltaREmbedded(0)
450 ,fQABckgHisto0RecCuts(0)
452 ,fQABckgHisto1RecCuts(0)
454 ,fQABckgHisto2RecCuts(0)
456 ,fQABckgHisto3RecCuts(0)
458 ,fQABckgHisto4RecCuts(0)
460 ,fFFBckgHisto0RecCuts(0)
462 ,fFFBckgHisto1RecCuts(0)
464 ,fFFBckgHisto2RecCuts(0)
466 ,fFFBckgHisto3RecCuts(0)
468 ,fFFBckgHisto4RecCuts(0)
470 ,fFFBckgHisto0RecEffRec(0)
471 ,fFFBckgHisto0SecRecNS(0)
472 ,fFFBckgHisto0SecRecS(0)
473 ,fFFBckgHisto0SecRecSsc(0)
475 ,fProNtracksLeadingJet(0)
477 ,fProNtracksLeadingJetGen(0)
478 ,fProDelR80pcPtGen(0)
479 ,fProNtracksLeadingJetBgrPerp2(0)
480 ,fProNtracksLeadingJetRecPrim(0)
481 ,fProDelR80pcPtRecPrim(0)
482 ,fProNtracksLeadingJetRecSecNS(0)
483 ,fProNtracksLeadingJetRecSecS(0)
484 ,fProNtracksLeadingJetRecSecSsc(0)
486 ,fOnlyLeadingJets(kFALSE)
489 ,fNumInclusivePIDtasks(0)
491 ,fNameInclusivePIDtask(0x0)
492 ,fNameJetPIDtask(0x0)
493 ,fInclusivePIDtask(0x0)
495 ,fUseInclusivePIDtask(kFALSE)
496 ,fUseJetPIDtask(kFALSE)
506 for(Int_t ii=0; ii<5; ii++){
507 fProDelRPtSum[ii] = 0;
508 fProDelRPtSumGen[ii] = 0;
509 fProDelRPtSumBgrPerp2[ii] = 0;
510 fProDelRPtSumRecPrim[ii] = 0;
511 fProDelRPtSumRecSecNS[ii] = 0;
512 fProDelRPtSumRecSecS[ii] = 0;
513 fProDelRPtSumRecSecSsc[ii] = 0;
516 for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
517 fIDFFHistosRecCuts[i] = 0x0;
518 fIDFFHistosGen[i] = 0x0;
520 fhDCA_XY_prim_MCID[i] = 0x0;
521 fhDCA_Z_prim_MCID[i] = 0x0;
523 fhDCA_XY_sec_MCID[i] = 0x0;
524 fhDCA_Z_sec_MCID[i] = 0x0;
527 DefineOutput(1,TList::Class());
530 //__________________________________________________________________________________________________________________________
531 AliAnalysisTaskIDFragmentationFunction::AliAnalysisTaskIDFragmentationFunction(const AliAnalysisTaskIDFragmentationFunction ©)
532 : AliAnalysisTaskSE()
535 ,fAODJets(copy.fAODJets)
536 ,fAODExtension(copy.fAODExtension)
537 ,fNonStdFile(copy.fNonStdFile)
538 ,fBranchRecJets(copy.fBranchRecJets)
539 ,fBranchRecBckgClusters(copy.fBranchRecBckgClusters)
540 ,fBranchGenJets(copy.fBranchGenJets)
541 ,fBranchEmbeddedJets(copy.fBranchEmbeddedJets)
542 ,fTrackTypeGen(copy.fTrackTypeGen)
543 ,fJetTypeGen(copy.fJetTypeGen)
544 ,fJetTypeRecEff(copy.fJetTypeRecEff)
545 ,fUseAODInputJets(copy.fUseAODInputJets)
546 ,fFilterMask(copy.fFilterMask)
547 ,fUsePhysicsSelection(copy.fUsePhysicsSelection)
548 ,fEvtSelectionMask(copy.fEvtSelectionMask)
549 ,fEventClass(copy.fEventClass)
550 ,fMaxVertexZ(copy.fMaxVertexZ)
551 ,fTrackPtCut(copy.fTrackPtCut)
552 ,fTrackEtaMin(copy.fTrackEtaMin)
553 ,fTrackEtaMax(copy.fTrackEtaMax)
554 ,fTrackPhiMin(copy.fTrackPhiMin)
555 ,fTrackPhiMax(copy.fTrackPhiMax)
556 ,fUseExtraTracks(copy.fUseExtraTracks)
557 ,fUseExtraTracksBgr(copy.fUseExtraTracksBgr)
558 ,fCutFractionPtEmbedded(copy.fCutFractionPtEmbedded)
559 ,fUseEmbeddedJetAxis(copy.fUseEmbeddedJetAxis)
560 ,fUseEmbeddedJetPt(copy.fUseEmbeddedJetPt)
561 ,fJetPtCut(copy.fJetPtCut)
562 ,fJetEtaMin(copy.fJetEtaMin)
563 ,fJetEtaMax(copy.fJetEtaMax)
564 ,fJetPhiMin(copy.fJetPhiMin)
565 ,fJetPhiMax(copy.fJetPhiMax)
566 ,fFFRadius(copy.fFFRadius)
567 ,fFFMinLTrackPt(copy.fFFMinLTrackPt)
568 ,fFFMaxTrackPt(copy.fFFMaxTrackPt)
569 ,fFFMinnTracks(copy.fFFMinnTracks)
570 ,fFFBckgRadius(copy.fFFBckgRadius)
571 ,fBckgMode(copy.fBckgMode)
572 ,fQAMode(copy.fQAMode)
573 ,fFFMode(copy.fFFMode)
574 ,fIDFFMode(copy.fIDFFMode)
575 ,fEffMode(copy.fEffMode)
576 ,fJSMode(copy.fJSMode)
577 ,fAvgTrials(copy.fAvgTrials)
578 ,fTracksRecCuts(copy.fTracksRecCuts)
579 ,fTracksRecCutsEfficiency(copy.fTracksRecCutsEfficiency)
580 ,fTracksGen(copy.fTracksGen)
581 ,fTracksAODMCCharged(copy.fTracksAODMCCharged)
582 ,fTracksAODMCChargedSecNS(copy.fTracksAODMCChargedSecNS)
583 ,fTracksAODMCChargedSecS(copy.fTracksAODMCChargedSecS)
584 ,fTracksRecQualityCuts(copy.fTracksRecQualityCuts)
585 ,fJetsRec(copy.fJetsRec)
586 ,fJetsRecCuts(copy.fJetsRecCuts)
587 ,fJetsGen(copy.fJetsGen)
588 ,fJetsRecEff(copy.fJetsRecEff)
589 ,fJetsEmbedded(copy.fJetsEmbedded)
590 ,fBckgJetsRec(copy.fBckgJetsRec)
591 ,fBckgJetsRecCuts(copy.fBckgJetsRecCuts)
592 ,fBckgJetsGen(copy.fBckgJetsGen)
593 ,fQATrackHistosRecCuts(copy.fQATrackHistosRecCuts)
594 ,fQATrackHistosGen(copy.fQATrackHistosGen)
595 ,fQAJetHistosRec(copy.fQAJetHistosRec)
596 ,fQAJetHistosRecCuts(copy.fQAJetHistosRecCuts)
597 ,fQAJetHistosRecCutsLeading(copy.fQAJetHistosRecCutsLeading)
598 ,fQAJetHistosGen(copy.fQAJetHistosGen)
599 ,fQAJetHistosGenLeading(copy.fQAJetHistosGenLeading)
600 ,fQAJetHistosRecEffLeading(copy.fQAJetHistosRecEffLeading)
601 ,fFFHistosRecCuts(copy.fFFHistosRecCuts)
602 ,fFFHistosRecCutsInc(copy.fFFHistosRecCutsInc)
603 ,fFFHistosRecLeadingTrack(copy.fFFHistosRecLeadingTrack)
604 ,fFFHistosGen(copy.fFFHistosGen)
605 ,fFFHistosGenInc(copy.fFFHistosGenInc)
606 ,fFFHistosGenLeadingTrack(copy.fFFHistosGenLeadingTrack)
607 ,fQATrackHighPtThreshold(copy.fQATrackHighPtThreshold)
608 ,fFFNBinsJetPt(copy.fFFNBinsJetPt)
609 ,fFFJetPtMin(copy.fFFJetPtMin)
610 ,fFFJetPtMax(copy.fFFJetPtMax)
611 ,fFFNBinsPt(copy.fFFNBinsPt)
612 ,fFFPtMin(copy.fFFPtMin)
613 ,fFFPtMax(copy.fFFPtMax)
614 ,fFFNBinsXi(copy.fFFNBinsXi)
615 ,fFFXiMin(copy.fFFXiMin)
616 ,fFFXiMax(copy.fFFXiMax)
617 ,fFFNBinsZ(copy.fFFNBinsZ)
618 ,fFFZMin(copy.fFFZMin)
619 ,fFFZMax(copy.fFFZMax)
620 ,fQAJetNBinsPt(copy.fQAJetNBinsPt)
621 ,fQAJetPtMin(copy.fQAJetPtMin)
622 ,fQAJetPtMax(copy.fQAJetPtMax)
623 ,fQAJetNBinsEta(copy.fQAJetNBinsEta)
624 ,fQAJetEtaMin(copy.fQAJetEtaMin)
625 ,fQAJetEtaMax(copy.fQAJetEtaMax)
626 ,fQAJetNBinsPhi(copy.fQAJetNBinsPhi)
627 ,fQAJetPhiMin(copy.fQAJetPhiMin)
628 ,fQAJetPhiMax(copy.fQAJetPhiMax)
629 ,fQATrackNBinsPt(copy.fQATrackNBinsPt)
630 ,fQATrackPtMin(copy.fQATrackPtMin)
631 ,fQATrackPtMax(copy.fQATrackPtMax)
632 ,fQATrackNBinsEta(copy.fQATrackNBinsEta)
633 ,fQATrackEtaMin(copy.fQATrackEtaMin)
634 ,fQATrackEtaMax(copy.fQATrackEtaMax)
635 ,fQATrackNBinsPhi(copy.fQATrackNBinsPhi)
636 ,fQATrackPhiMin(copy.fQATrackPhiMin)
637 ,fQATrackPhiMax(copy.fQATrackPhiMax)
638 ,fCommonHistList(copy.fCommonHistList)
639 ,fh1EvtSelection(copy.fh1EvtSelection)
640 ,fh1VtxSelection(copy.fh1VtxSelection)
641 ,fh1VertexNContributors(copy.fh1VertexNContributors)
642 ,fh1VertexZ(copy.fh1VertexZ)
643 ,fh1EvtMult(copy.fh1EvtMult)
644 ,fh1EvtCent(copy.fh1EvtCent)
645 ,fh1Xsec(copy.fh1Xsec)
646 ,fh1Trials(copy.fh1Trials)
647 ,fh1PtHard(copy.fh1PtHard)
648 ,fh1PtHardTrials(copy.fh1PtHardTrials)
649 ,fh1nRecJetsCuts(copy.fh1nRecJetsCuts)
650 ,fh1nGenJets(copy.fh1nGenJets)
651 ,fh1nRecEffJets(copy.fh1nRecEffJets)
652 ,fh1nEmbeddedJets(copy.fh1nEmbeddedJets)
653 ,fh1nRecBckgJetsCuts(copy.fh1nRecBckgJetsCuts)
654 ,fh1nGenBckgJets(copy.fh1nGenBckgJets)
655 ,fh2PtRecVsGenPrim(copy.fh2PtRecVsGenPrim)
656 ,fh2PtRecVsGenSec(copy.fh2PtRecVsGenSec)
657 ,fhDCA_XY(copy.fhDCA_XY)
658 ,fhDCA_Z(copy.fhDCA_Z)
659 ,fhJetPtRefMultEta5(copy.fhJetPtRefMultEta5)
660 ,fhJetPtRefMultEta8(copy.fhJetPtRefMultEta8)
661 ,fhJetPtMultPercent(copy.fhJetPtMultPercent)
662 ,fQATrackHistosRecEffGen(copy.fQATrackHistosRecEffGen)
663 ,fQATrackHistosRecEffRec(copy.fQATrackHistosRecEffRec)
664 ,fQATrackHistosSecRecNS(copy.fQATrackHistosSecRecNS)
665 ,fQATrackHistosSecRecS(copy.fQATrackHistosSecRecS)
666 ,fQATrackHistosSecRecSsc(copy.fQATrackHistosSecRecSsc)
667 ,fFFHistosRecEffRec(copy.fFFHistosRecEffRec)
668 ,fFFHistosSecRecNS(copy.fFFHistosSecRecNS)
669 ,fFFHistosSecRecS(copy.fFFHistosSecRecS)
670 ,fFFHistosSecRecSsc(copy.fFFHistosSecRecSsc)
672 ,fh1BckgMult0(copy.fh1BckgMult0)
673 ,fh1BckgMult1(copy.fh1BckgMult1)
674 ,fh1BckgMult2(copy.fh1BckgMult2)
675 ,fh1BckgMult3(copy.fh1BckgMult3)
676 ,fh1BckgMult4(copy.fh1BckgMult4)
677 ,fh1FractionPtEmbedded(copy.fh1FractionPtEmbedded)
678 ,fh1IndexEmbedded(copy.fh1IndexEmbedded)
679 ,fh2DeltaPtVsJetPtEmbedded(copy.fh2DeltaPtVsJetPtEmbedded)
680 ,fh2DeltaPtVsRecJetPtEmbedded(copy.fh2DeltaPtVsRecJetPtEmbedded)
681 ,fh1DeltaREmbedded(copy.fh1DeltaREmbedded)
682 ,fQABckgHisto0RecCuts(copy.fQABckgHisto0RecCuts)
683 ,fQABckgHisto0Gen(copy.fQABckgHisto0Gen)
684 ,fQABckgHisto1RecCuts(copy.fQABckgHisto1RecCuts)
685 ,fQABckgHisto1Gen(copy.fQABckgHisto1Gen)
686 ,fQABckgHisto2RecCuts(copy.fQABckgHisto2RecCuts)
687 ,fQABckgHisto2Gen(copy.fQABckgHisto2Gen)
688 ,fQABckgHisto3RecCuts(copy.fQABckgHisto3RecCuts)
689 ,fQABckgHisto3Gen(copy.fQABckgHisto3Gen)
690 ,fQABckgHisto4RecCuts(copy.fQABckgHisto4RecCuts)
691 ,fQABckgHisto4Gen(copy.fQABckgHisto4Gen)
692 ,fFFBckgHisto0RecCuts(copy.fFFBckgHisto0RecCuts)
693 ,fFFBckgHisto0Gen(copy.fFFBckgHisto0Gen)
694 ,fFFBckgHisto1RecCuts(copy.fFFBckgHisto1RecCuts)
695 ,fFFBckgHisto1Gen(copy.fFFBckgHisto1Gen)
696 ,fFFBckgHisto2RecCuts(copy.fFFBckgHisto2RecCuts)
697 ,fFFBckgHisto2Gen(copy.fFFBckgHisto2Gen)
698 ,fFFBckgHisto3RecCuts(copy.fFFBckgHisto3RecCuts)
699 ,fFFBckgHisto3Gen(copy.fFFBckgHisto3Gen)
700 ,fFFBckgHisto4RecCuts(copy.fFFBckgHisto4RecCuts)
701 ,fFFBckgHisto4Gen(copy.fFFBckgHisto4Gen)
702 ,fFFBckgHisto0RecEffRec(copy.fFFBckgHisto0RecEffRec)
703 ,fFFBckgHisto0SecRecNS(copy.fFFBckgHisto0SecRecNS)
704 ,fFFBckgHisto0SecRecS(copy.fFFBckgHisto0SecRecS)
705 ,fFFBckgHisto0SecRecSsc(copy.fFFBckgHisto0SecRecSsc)
707 ,fProNtracksLeadingJet(copy.fProNtracksLeadingJet)
708 ,fProDelR80pcPt(copy.fProDelR80pcPt)
709 ,fProNtracksLeadingJetGen(copy.fProNtracksLeadingJetGen)
710 ,fProDelR80pcPtGen(copy.fProDelR80pcPtGen)
711 ,fProNtracksLeadingJetBgrPerp2(copy.fProNtracksLeadingJetBgrPerp2)
712 ,fProNtracksLeadingJetRecPrim(copy.fProNtracksLeadingJetRecPrim)
713 ,fProDelR80pcPtRecPrim(copy.fProDelR80pcPtRecPrim)
714 ,fProNtracksLeadingJetRecSecNS(copy.fProNtracksLeadingJetRecSecNS)
715 ,fProNtracksLeadingJetRecSecS(copy.fProNtracksLeadingJetRecSecS)
716 ,fProNtracksLeadingJetRecSecSsc(copy.fProNtracksLeadingJetRecSecSsc)
717 ,fRandom(copy.fRandom)
718 ,fOnlyLeadingJets(copy.fOnlyLeadingJets)
719 ,fAnaUtils(copy.fAnaUtils)
721 ,fNumInclusivePIDtasks(copy.fNumInclusivePIDtasks)
722 ,fNumJetPIDtasks(copy.fNumJetPIDtasks)
723 ,fNameInclusivePIDtask(0x0)
724 ,fNameJetPIDtask(0x0)
725 ,fInclusivePIDtask(0x0)
727 ,fUseInclusivePIDtask(copy.fUseInclusivePIDtask)
728 ,fUseJetPIDtask(copy.fUseJetPIDtask)
732 fBckgType[0] = copy.fBckgType[0];
733 fBckgType[1] = copy.fBckgType[1];
734 fBckgType[2] = copy.fBckgType[2];
735 fBckgType[3] = copy.fBckgType[3];
736 fBckgType[4] = copy.fBckgType[4];
739 for(Int_t ii=0; ii<5; ii++){
740 fProDelRPtSum[ii] = copy.fProDelRPtSum[ii];
741 fProDelRPtSumGen[ii] = copy.fProDelRPtSumGen[ii];
742 fProDelRPtSumBgrPerp2[ii] = copy.fProDelRPtSumBgrPerp2[ii];
743 fProDelRPtSumRecPrim[ii] = copy.fProDelRPtSumRecPrim[ii];
744 fProDelRPtSumRecSecNS[ii] = copy.fProDelRPtSumRecSecNS[ii];
745 fProDelRPtSumRecSecS[ii] = copy.fProDelRPtSumRecSecS[ii];
746 fProDelRPtSumRecSecSsc[ii] = copy.fProDelRPtSumRecSecSsc[ii];
749 for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
750 fIDFFHistosRecCuts[i] = 0x0;
751 if (copy.fIDFFHistosRecCuts[i])
752 fIDFFHistosRecCuts[i] = copy.fIDFFHistosRecCuts[i];
754 fIDFFHistosGen[i] = 0x0;
755 if (copy.fIDFFHistosGen[i])
756 fIDFFHistosGen[i] = copy.fIDFFHistosGen[i];
759 fhDCA_XY_prim_MCID[i] = 0x0;
760 if (copy.fhDCA_XY_prim_MCID[i])
761 fhDCA_XY_prim_MCID[i] = copy.fhDCA_XY_prim_MCID[i];
763 fhDCA_Z_prim_MCID[i] = 0x0;
764 if (copy.fhDCA_Z_prim_MCID[i])
765 fhDCA_Z_prim_MCID[i] = copy.fhDCA_Z_prim_MCID[i];
767 fhDCA_XY_sec_MCID[i] = 0x0;
768 if (copy.fhDCA_XY_sec_MCID[i])
769 fhDCA_XY_sec_MCID[i] = copy.fhDCA_XY_sec_MCID[i];
771 fhDCA_Z_sec_MCID[i] = 0x0;
772 if (copy.fhDCA_Z_sec_MCID[i])
773 fhDCA_Z_sec_MCID[i] = copy.fhDCA_Z_sec_MCID[i];
776 if (fNumInclusivePIDtasks > 0) {
777 fNameInclusivePIDtask = new TString[fNumInclusivePIDtasks];
778 fInclusivePIDtask = new AliAnalysisTaskPID*[fNumInclusivePIDtasks];
780 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
781 fNameInclusivePIDtask[i] = "";
782 fInclusivePIDtask[i] = 0x0;
784 if (copy.fNameInclusivePIDtask[i])
785 fNameInclusivePIDtask[i] = copy.fNameInclusivePIDtask[i];
787 if (copy.fInclusivePIDtask[i])
788 fInclusivePIDtask[i] = copy.fInclusivePIDtask[i];
792 if (fNumJetPIDtasks > 0) {
793 fNameJetPIDtask = new TString[fNumJetPIDtasks];
794 fJetPIDtask = new AliAnalysisTaskPID*[fNumJetPIDtasks];
796 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
797 fNameJetPIDtask[i] = "";
798 fJetPIDtask[i] = 0x0;
800 if (copy.fNameJetPIDtask[i])
801 fNameJetPIDtask[i] = copy.fNameJetPIDtask[i];
803 if (copy.fJetPIDtask[i])
804 fJetPIDtask[i] = copy.fJetPIDtask[i];
809 // _________________________________________________________________________________________________________________________________
810 AliAnalysisTaskIDFragmentationFunction& AliAnalysisTaskIDFragmentationFunction::operator=(const AliAnalysisTaskIDFragmentationFunction& o)
816 AliAnalysisTaskSE::operator=(o);
819 fAODJets = o.fAODJets;
820 fAODExtension = o.fAODExtension;
821 fNonStdFile = o.fNonStdFile;
822 fBranchRecJets = o.fBranchRecJets;
823 fBranchRecBckgClusters = o.fBranchRecBckgClusters;
824 fBranchGenJets = o.fBranchGenJets;
825 fBranchEmbeddedJets = o.fBranchEmbeddedJets;
826 fTrackTypeGen = o.fTrackTypeGen;
827 fJetTypeGen = o.fJetTypeGen;
828 fJetTypeRecEff = o.fJetTypeRecEff;
829 fUseAODInputJets = o.fUseAODInputJets;
830 fFilterMask = o.fFilterMask;
831 fUsePhysicsSelection = o.fUsePhysicsSelection;
832 fEvtSelectionMask = o.fEvtSelectionMask;
833 fEventClass = o.fEventClass;
834 fMaxVertexZ = o.fMaxVertexZ;
835 fTrackPtCut = o.fTrackPtCut;
836 fTrackEtaMin = o.fTrackEtaMin;
837 fTrackEtaMax = o.fTrackEtaMax;
838 fTrackPhiMin = o.fTrackPhiMin;
839 fTrackPhiMax = o.fTrackPhiMax;
840 fUseExtraTracks = o.fUseExtraTracks;
841 fUseExtraTracksBgr = o.fUseExtraTracksBgr;
842 fCutFractionPtEmbedded = o.fCutFractionPtEmbedded;
843 fUseEmbeddedJetAxis = o.fUseEmbeddedJetAxis;
844 fUseEmbeddedJetPt = o.fUseEmbeddedJetPt;
845 fJetPtCut = o.fJetPtCut;
846 fJetEtaMin = o.fJetEtaMin;
847 fJetEtaMax = o.fJetEtaMax;
848 fJetPhiMin = o.fJetPhiMin;
849 fJetPhiMax = o.fJetPhiMin;
850 fFFRadius = o.fFFRadius;
851 fFFMinLTrackPt = o.fFFMinLTrackPt;
852 fFFMaxTrackPt = o.fFFMaxTrackPt;
853 fFFMinnTracks = o.fFFMinnTracks;
854 fFFBckgRadius = o.fFFBckgRadius;
855 fBckgMode = o.fBckgMode;
858 fIDFFMode = o.fIDFFMode;
859 fEffMode = o.fEffMode;
861 fBckgType[0] = o.fBckgType[0];
862 fBckgType[1] = o.fBckgType[1];
863 fBckgType[2] = o.fBckgType[2];
864 fBckgType[3] = o.fBckgType[3];
865 fBckgType[4] = o.fBckgType[4];
866 fAvgTrials = o.fAvgTrials;
867 fTracksRecCuts = o.fTracksRecCuts;
868 fTracksRecCutsEfficiency = o.fTracksRecCutsEfficiency;
869 fTracksGen = o.fTracksGen;
870 fTracksAODMCCharged = o.fTracksAODMCCharged;
871 fTracksAODMCChargedSecNS = o.fTracksAODMCChargedSecNS;
872 fTracksAODMCChargedSecS = o.fTracksAODMCChargedSecS;
873 fTracksRecQualityCuts = o.fTracksRecQualityCuts;
874 fJetsRec = o.fJetsRec;
875 fJetsRecCuts = o.fJetsRecCuts;
876 fJetsGen = o.fJetsGen;
877 fJetsRecEff = o.fJetsRecEff;
878 fJetsEmbedded = o.fJetsEmbedded;
879 fBckgJetsRec = o.fBckgJetsRec;
880 fBckgJetsRecCuts = o.fBckgJetsRecCuts;
881 fBckgJetsGen = o.fBckgJetsGen;
882 fQATrackHistosRecCuts = o.fQATrackHistosRecCuts;
883 fQATrackHistosGen = o.fQATrackHistosGen;
884 fQAJetHistosRec = o.fQAJetHistosRec;
885 fQAJetHistosRecCuts = o.fQAJetHistosRecCuts;
886 fQAJetHistosRecCutsLeading = o.fQAJetHistosRecCutsLeading;
887 fQAJetHistosGen = o.fQAJetHistosGen;
888 fQAJetHistosGenLeading = o.fQAJetHistosGenLeading;
889 fQAJetHistosRecEffLeading = o.fQAJetHistosRecEffLeading;
890 fFFHistosRecCuts = o.fFFHistosRecCuts;
891 fFFHistosRecCutsInc = o.fFFHistosRecCutsInc;
892 fFFHistosRecLeadingTrack = o.fFFHistosRecLeadingTrack;
893 fFFHistosGen = o.fFFHistosGen;
894 fFFHistosGenInc = o.fFFHistosGenInc;
895 fFFHistosGenLeadingTrack = o.fFFHistosGenLeadingTrack;
896 fQATrackHighPtThreshold = o.fQATrackHighPtThreshold;
897 fFFNBinsJetPt = o.fFFNBinsJetPt;
898 fFFJetPtMin = o.fFFJetPtMin;
899 fFFJetPtMax = o.fFFJetPtMax;
900 fFFNBinsPt = o.fFFNBinsPt;
901 fFFPtMin = o.fFFPtMin;
902 fFFPtMax = o.fFFPtMax;
903 fFFNBinsXi = o.fFFNBinsXi;
904 fFFXiMin = o.fFFXiMin;
905 fFFXiMax = o.fFFXiMax;
906 fFFNBinsZ = o.fFFNBinsZ;
909 fQAJetNBinsPt = o.fQAJetNBinsPt;
910 fQAJetPtMin = o.fQAJetPtMin;
911 fQAJetPtMax = o.fQAJetPtMax;
912 fQAJetNBinsEta = o.fQAJetNBinsEta;
913 fQAJetEtaMin = o.fQAJetEtaMin;
914 fQAJetEtaMax = o.fQAJetEtaMax;
915 fQAJetNBinsPhi = o.fQAJetNBinsPhi;
916 fQAJetPhiMin = o.fQAJetPhiMin;
917 fQAJetPhiMax = o.fQAJetPhiMax;
918 fQATrackNBinsPt = o.fQATrackNBinsPt;
919 fQATrackPtMin = o.fQATrackPtMin;
920 fQATrackPtMax = o.fQATrackPtMax;
921 fQATrackNBinsEta = o.fQATrackNBinsEta;
922 fQATrackEtaMin = o.fQATrackEtaMin;
923 fQATrackEtaMax = o.fQATrackEtaMax;
924 fQATrackNBinsPhi = o.fQATrackNBinsPhi;
925 fQATrackPhiMin = o.fQATrackPhiMin;
926 fQATrackPhiMax = o.fQATrackPhiMax;
927 fCommonHistList = o.fCommonHistList;
928 fh1EvtSelection = o.fh1EvtSelection;
929 fh1VtxSelection = o.fh1VtxSelection;
930 fh1VertexNContributors = o.fh1VertexNContributors;
931 fh1VertexZ = o.fh1VertexZ;
932 fh1EvtMult = o.fh1EvtMult;
933 fh1EvtCent = o.fh1EvtCent;
935 fh1Trials = o.fh1Trials;
936 fh1PtHard = o.fh1PtHard;
937 fh1PtHardTrials = o.fh1PtHardTrials;
938 fh1nRecJetsCuts = o.fh1nRecJetsCuts;
939 fh1nGenJets = o.fh1nGenJets;
940 fh1nRecEffJets = o.fh1nRecEffJets;
941 fh1nEmbeddedJets = o.fh1nEmbeddedJets;
942 fh2PtRecVsGenPrim = o.fh2PtRecVsGenPrim;
943 fh2PtRecVsGenSec = o.fh2PtRecVsGenSec;
944 fQATrackHistosRecEffGen = o.fQATrackHistosRecEffGen;
945 fQATrackHistosRecEffRec = o.fQATrackHistosRecEffRec;
946 fQATrackHistosSecRecNS = o.fQATrackHistosSecRecNS;
947 fQATrackHistosSecRecS = o.fQATrackHistosSecRecS;
948 fQATrackHistosSecRecSsc = o.fQATrackHistosSecRecSsc;
949 fFFHistosRecEffRec = o.fFFHistosRecEffRec;
950 fFFHistosSecRecNS = o.fFFHistosSecRecNS;
951 fFFHistosSecRecS = o.fFFHistosSecRecS;
952 fFFHistosSecRecSsc = o.fFFHistosSecRecSsc;
954 fh1BckgMult0 = o.fh1BckgMult0;
955 fh1BckgMult1 = o.fh1BckgMult1;
956 fh1BckgMult2 = o.fh1BckgMult2;
957 fh1BckgMult3 = o.fh1BckgMult3;
958 fh1BckgMult4 = o.fh1BckgMult4;
959 fh1FractionPtEmbedded = o.fh1FractionPtEmbedded;
960 fh1IndexEmbedded = o.fh1IndexEmbedded;
961 fh2DeltaPtVsJetPtEmbedded = o.fh2DeltaPtVsJetPtEmbedded;
962 fh2DeltaPtVsRecJetPtEmbedded = o.fh2DeltaPtVsRecJetPtEmbedded;
963 fh1DeltaREmbedded = o.fh1DeltaREmbedded;
964 fhDCA_XY = o.fhDCA_XY;
966 fhJetPtRefMultEta5 = o.fhJetPtRefMultEta5;
967 fhJetPtRefMultEta8 = o.fhJetPtRefMultEta8;
968 fhJetPtMultPercent = o.fhJetPtMultPercent;
969 fQABckgHisto0RecCuts = o.fQABckgHisto0RecCuts;
970 fQABckgHisto0Gen = o.fQABckgHisto0Gen;
971 fQABckgHisto1RecCuts = o.fQABckgHisto1RecCuts;
972 fQABckgHisto1Gen = o.fQABckgHisto1Gen;
973 fQABckgHisto2RecCuts = o.fQABckgHisto2RecCuts;
974 fQABckgHisto2Gen = o.fQABckgHisto2Gen;
975 fQABckgHisto3RecCuts = o.fQABckgHisto3RecCuts;
976 fQABckgHisto3Gen = o.fQABckgHisto3Gen;
977 fQABckgHisto4RecCuts = o.fQABckgHisto4RecCuts;
978 fQABckgHisto4Gen = o.fQABckgHisto4Gen;
979 fFFBckgHisto0RecCuts = o.fFFBckgHisto0RecCuts;
980 fFFBckgHisto0Gen = o.fFFBckgHisto0Gen;
981 fFFBckgHisto1RecCuts = o.fFFBckgHisto1RecCuts;
982 fFFBckgHisto1Gen = o.fFFBckgHisto1Gen;
983 fFFBckgHisto2RecCuts = o.fFFBckgHisto2RecCuts;
984 fFFBckgHisto2Gen = o.fFFBckgHisto2Gen;
985 fFFBckgHisto3RecCuts = o.fFFBckgHisto3RecCuts;
986 fFFBckgHisto3Gen = o.fFFBckgHisto3Gen;
987 fFFBckgHisto4RecCuts = o.fFFBckgHisto4RecCuts;
988 fFFBckgHisto4Gen = o.fFFBckgHisto4Gen;
989 fFFBckgHisto0RecEffRec = o.fFFBckgHisto0RecEffRec;
990 fFFBckgHisto0SecRecNS = o.fFFBckgHisto0SecRecNS;
991 fFFBckgHisto0SecRecS = o.fFFBckgHisto0SecRecS;
992 fFFBckgHisto0SecRecSsc = o.fFFBckgHisto0SecRecSsc;
993 fProNtracksLeadingJet = o.fProNtracksLeadingJet;
994 fProDelR80pcPt = o.fProDelR80pcPt;
995 fProNtracksLeadingJetGen = o.fProNtracksLeadingJetGen;
996 fProDelR80pcPtGen = o.fProDelR80pcPtGen;
997 fProNtracksLeadingJetBgrPerp2 = o.fProNtracksLeadingJetBgrPerp2;
998 fProNtracksLeadingJetRecPrim = o.fProNtracksLeadingJetRecPrim;
999 fProDelR80pcPtRecPrim = o.fProDelR80pcPtRecPrim;
1000 fProNtracksLeadingJetRecSecNS = o.fProNtracksLeadingJetRecSecNS;
1001 fProNtracksLeadingJetRecSecS = o.fProNtracksLeadingJetRecSecS;
1002 fProNtracksLeadingJetRecSecSsc = o.fProNtracksLeadingJetRecSecSsc;
1003 fRandom = o.fRandom;
1004 fOnlyLeadingJets = o.fOnlyLeadingJets;
1005 fAnaUtils = o.fAnaUtils;
1008 fUseInclusivePIDtask = o.fUseInclusivePIDtask;
1009 fUseJetPIDtask = o.fUseJetPIDtask;
1013 for(Int_t ii=0; ii<5; ii++){
1014 fProDelRPtSum[ii] = o.fProDelRPtSum[ii];
1015 fProDelRPtSumGen[ii] = o.fProDelRPtSumGen[ii];
1016 fProDelRPtSumBgrPerp2[ii] = o.fProDelRPtSumBgrPerp2[ii];
1017 fProDelRPtSumRecPrim[ii] = o.fProDelRPtSumRecPrim[ii];
1018 fProDelRPtSumRecSecNS[ii] = o.fProDelRPtSumRecSecNS[ii];
1019 fProDelRPtSumRecSecS[ii] = o.fProDelRPtSumRecSecS[ii];
1020 fProDelRPtSumRecSecSsc[ii] = o.fProDelRPtSumRecSecSsc[ii];
1023 for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
1024 fIDFFHistosRecCuts[i] = 0x0;
1025 if (o.fIDFFHistosRecCuts[i])
1026 fIDFFHistosRecCuts[i] = o.fIDFFHistosRecCuts[i];
1028 fIDFFHistosGen[i] = 0x0;
1029 if (o.fIDFFHistosGen[i])
1030 fIDFFHistosGen[i] = o.fIDFFHistosGen[i];
1032 fhDCA_XY_prim_MCID[i] = 0x0;
1033 if (o.fhDCA_XY_prim_MCID[i])
1034 fhDCA_XY_prim_MCID[i] = o.fhDCA_XY_prim_MCID[i];
1036 fhDCA_Z_prim_MCID[i] = 0x0;
1037 if (o.fhDCA_Z_prim_MCID[i])
1038 fhDCA_Z_prim_MCID[i] = o.fhDCA_Z_prim_MCID[i];
1040 fhDCA_XY_sec_MCID[i] = 0x0;
1041 if (o.fhDCA_XY_sec_MCID[i])
1042 fhDCA_XY_sec_MCID[i] = o.fhDCA_XY_sec_MCID[i];
1044 fhDCA_Z_sec_MCID[i] = 0x0;
1045 if (o.fhDCA_Z_sec_MCID[i])
1046 fhDCA_Z_sec_MCID[i] = o.fhDCA_Z_sec_MCID[i];
1051 if (fNumJetPIDtasks > 0) {
1052 delete [] fNameJetPIDtask;
1053 fNameJetPIDtask = 0x0;
1055 delete [] fJetPIDtask;
1059 fNumJetPIDtasks = o.fNumJetPIDtasks;
1062 if (fNumInclusivePIDtasks > 0) {
1063 delete [] fNameInclusivePIDtask;
1064 fNameInclusivePIDtask = 0x0;
1066 delete [] fInclusivePIDtask;
1067 fInclusivePIDtask = 0x0;
1070 fNumInclusivePIDtasks = o.fNumInclusivePIDtasks;
1073 if (fNumInclusivePIDtasks > 0) {
1074 fNameInclusivePIDtask = new TString[fNumInclusivePIDtasks];
1075 fInclusivePIDtask = new AliAnalysisTaskPID*[fNumInclusivePIDtasks];
1077 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
1078 fNameInclusivePIDtask[i] = "";
1079 fInclusivePIDtask[i] = 0x0;
1081 if (o.fNameInclusivePIDtask[i])
1082 fNameInclusivePIDtask[i] = o.fNameInclusivePIDtask[i];
1084 if (o.fInclusivePIDtask[i])
1085 fInclusivePIDtask[i] = o.fInclusivePIDtask[i];
1089 if (fNumJetPIDtasks > 0) {
1090 fNameJetPIDtask = new TString[fNumJetPIDtasks];
1091 fJetPIDtask = new AliAnalysisTaskPID*[fNumJetPIDtasks];
1093 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
1094 fNameJetPIDtask[i] = "";
1095 fJetPIDtask[i] = 0x0;
1097 if (o.fNameJetPIDtask[i])
1098 fNameJetPIDtask[i] = o.fNameJetPIDtask[i];
1100 if (o.fJetPIDtask[i])
1101 fJetPIDtask[i] = o.fJetPIDtask[i];
1108 //___________________________________________________________________________
1109 AliAnalysisTaskIDFragmentationFunction::~AliAnalysisTaskIDFragmentationFunction()
1113 if(fTracksRecCuts) delete fTracksRecCuts;
1114 if(fTracksRecCutsEfficiency) delete fTracksRecCutsEfficiency;
1115 if(fTracksGen) delete fTracksGen;
1116 if(fTracksAODMCCharged) delete fTracksAODMCCharged;
1117 if(fTracksAODMCChargedSecNS) delete fTracksAODMCChargedSecNS;
1118 if(fTracksAODMCChargedSecS) delete fTracksAODMCChargedSecS;
1119 if(fTracksRecQualityCuts) delete fTracksRecQualityCuts;
1120 if(fJetsRec) delete fJetsRec;
1121 if(fJetsRecCuts) delete fJetsRecCuts;
1122 if(fJetsGen) delete fJetsGen;
1123 if(fJetsRecEff) delete fJetsRecEff;
1124 if(fJetsEmbedded) delete fJetsEmbedded;
1127 (fBckgType[0]==kBckgClusters || fBckgType[1]==kBckgClusters || fBckgType[2]==kBckgClusters || fBckgType[3]==kBckgClusters || fBckgType[4]==kBckgClusters ||
1128 fBckgType[0]==kBckgClustersOutLeading || fBckgType[1]==kBckgClustersOutLeading || fBckgType[2]==kBckgClustersOutLeading ||
1129 fBckgType[3]==kBckgClustersOutLeading || fBckgType[4]==kBckgClustersOutLeading)){
1131 if(fBckgJetsRec) delete fBckgJetsRec;
1132 if(fBckgJetsRecCuts) delete fBckgJetsRecCuts;
1133 if(fBckgJetsGen) delete fBckgJetsGen;
1135 if(fRandom) delete fRandom;
1137 delete [] fNameInclusivePIDtask;
1138 fNameInclusivePIDtask = 0x0;
1140 delete [] fInclusivePIDtask;
1141 fInclusivePIDtask = 0x0;
1143 delete [] fNameJetPIDtask;
1144 fNameJetPIDtask = 0x0;
1146 delete [] fJetPIDtask;
1153 //______________________________________________________________________________________________________
1154 AliAnalysisTaskIDFragmentationFunction::AliFragFuncHistos::AliFragFuncHistos(const char* name,
1155 Int_t nJetPt, Float_t jetPtMin, Float_t jetPtMax,
1156 Int_t nPt, Float_t ptMin, Float_t ptMax,
1157 Int_t nXi, Float_t xiMin, Float_t xiMax,
1158 Int_t nZ , Float_t zMin , Float_t zMax)
1160 ,fNBinsJetPt(nJetPt)
1161 ,fJetPtMin(jetPtMin)
1162 ,fJetPtMax(jetPtMax)
1178 // default constructor
1182 //___________________________________________________________________________
1183 AliAnalysisTaskIDFragmentationFunction::AliFragFuncHistos::AliFragFuncHistos(const AliFragFuncHistos& copy)
1185 ,fNBinsJetPt(copy.fNBinsJetPt)
1186 ,fJetPtMin(copy.fJetPtMin)
1187 ,fJetPtMax(copy.fJetPtMax)
1188 ,fNBinsPt(copy.fNBinsPt)
1189 ,fPtMin(copy.fPtMin)
1190 ,fPtMax(copy.fPtMax)
1191 ,fNBinsXi(copy.fNBinsXi)
1192 ,fXiMin(copy.fXiMin)
1193 ,fXiMax(copy.fXiMax)
1194 ,fNBinsZ(copy.fNBinsZ)
1197 ,fh2TrackPt(copy.fh2TrackPt)
1200 ,fh1JetPt(copy.fh1JetPt)
1201 ,fNameFF(copy.fNameFF)
1206 //_______________________________________________________________________________________________________________________________________________________________
1207 AliAnalysisTaskIDFragmentationFunction::AliFragFuncHistos& AliAnalysisTaskIDFragmentationFunction::AliFragFuncHistos::operator=(const AliAnalysisTaskIDFragmentationFunction::AliFragFuncHistos& o)
1212 TObject::operator=(o);
1213 fNBinsJetPt = o.fNBinsJetPt;
1214 fJetPtMin = o.fJetPtMin;
1215 fJetPtMax = o.fJetPtMax;
1216 fNBinsPt = o.fNBinsPt;
1219 fNBinsXi = o.fNBinsXi;
1222 fNBinsZ = o.fNBinsZ;
1225 fh2TrackPt = o.fh2TrackPt;
1228 fh1JetPt = o.fh1JetPt;
1229 fNameFF = o.fNameFF;
1235 //_________________________________________________________
1236 AliAnalysisTaskIDFragmentationFunction::AliFragFuncHistos::~AliFragFuncHistos()
1240 if(fh1JetPt) delete fh1JetPt;
1241 if(fh2TrackPt) delete fh2TrackPt;
1242 if(fh2Xi) delete fh2Xi;
1243 if(fh2Z) delete fh2Z;
1246 //_________________________________________________________________
1247 void AliAnalysisTaskIDFragmentationFunction::AliFragFuncHistos::DefineHistos()
1251 fh1JetPt = new TH1F(Form("fh1FFJetPt%s", fNameFF.Data()),"",fNBinsJetPt,fJetPtMin,fJetPtMax);
1252 fh2TrackPt = new TH2F(Form("fh2FFTrackPt%s",fNameFF.Data()),"",fNBinsJetPt, fJetPtMin, fJetPtMax,fNBinsPt, fPtMin, fPtMax);
1253 fh2Z = new TH2F(Form("fh2FFZ%s",fNameFF.Data()),"",fNBinsJetPt, fJetPtMin, fJetPtMax, fNBinsZ, fZMin, fZMax);
1254 fh2Xi = new TH2F(Form("fh2FFXi%s",fNameFF.Data()),"",fNBinsJetPt, fJetPtMin, fJetPtMax, fNBinsXi, fXiMin, fXiMax);
1256 AliAnalysisTaskIDFragmentationFunction::SetProperties(fh1JetPt, "p_{T} [GeV/c]", "entries");
1257 AliAnalysisTaskIDFragmentationFunction::SetProperties(fh2TrackPt,"jet p_{T} [GeV/c]","p_{T} [GeV/c]","entries");
1258 AliAnalysisTaskIDFragmentationFunction::SetProperties(fh2Xi,"jet p_{T} [GeV/c]","#xi", "entries");
1259 AliAnalysisTaskIDFragmentationFunction::SetProperties(fh2Z,"jet p_{T} [GeV/c]","z","entries");
1262 //_______________________________________________________________________________________________________________
1263 void AliAnalysisTaskIDFragmentationFunction::AliFragFuncHistos::FillFF(Float_t trackPt, Float_t jetPt, Bool_t incrementJetPt, Float_t norm,
1264 Bool_t scaleStrangeness, Float_t scaleFacStrangeness)
1268 if(incrementJetPt && norm) fh1JetPt->Fill(jetPt,1/norm);
1269 else if(incrementJetPt) fh1JetPt->Fill(jetPt);
1271 // Added for proper normalization of FF background estimation
1272 // when zero track are found in the background region
1273 if((int)trackPt==-1) return;
1275 if(norm)fh2TrackPt->Fill(jetPt,trackPt,1/norm);
1276 else if(scaleStrangeness) fh2TrackPt->Fill(jetPt,trackPt,scaleFacStrangeness);
1277 else fh2TrackPt->Fill(jetPt,trackPt);
1279 Double_t z = -1., xi = -1.;
1280 AliAnalysisTaskPID::GetJetTrackObservables(trackPt, jetPt, z, xi);
1284 fh2Xi->Fill(jetPt,xi,1/norm);
1285 fh2Z->Fill(jetPt,z,1/norm);
1287 else if(scaleStrangeness){
1288 fh2Xi->Fill(jetPt,xi,scaleFacStrangeness);
1289 fh2Z->Fill(jetPt,z,scaleFacStrangeness);
1292 fh2Xi->Fill(jetPt,xi);
1293 fh2Z->Fill(jetPt,z);
1297 //_________________________________________________________________________________
1298 void AliAnalysisTaskIDFragmentationFunction::AliFragFuncHistos::AddToOutput(TList* list) const
1300 // add histos to list
1302 list->Add(fh1JetPt);
1304 list->Add(fh2TrackPt);
1309 //_________________________________________________________________________________________________________
1310 AliAnalysisTaskIDFragmentationFunction::AliFragFuncQAJetHistos::AliFragFuncQAJetHistos(const char* name,
1311 Int_t nPt, Float_t ptMin, Float_t ptMax,
1312 Int_t nEta, Float_t etaMin, Float_t etaMax,
1313 Int_t nPhi, Float_t phiMin, Float_t phiMax)
1328 // default constructor
1331 //____________________________________________________________________________________
1332 AliAnalysisTaskIDFragmentationFunction::AliFragFuncQAJetHistos::AliFragFuncQAJetHistos(const AliFragFuncQAJetHistos& copy)
1334 ,fNBinsPt(copy.fNBinsPt)
1335 ,fPtMin(copy.fPtMin)
1336 ,fPtMax(copy.fPtMax)
1337 ,fNBinsEta(copy.fNBinsEta)
1338 ,fEtaMin(copy.fEtaMin)
1339 ,fEtaMax(copy.fEtaMax)
1340 ,fNBinsPhi(copy.fNBinsPhi)
1341 ,fPhiMin(copy.fPhiMin)
1342 ,fPhiMax(copy.fPhiMax)
1343 ,fh2EtaPhi(copy.fh2EtaPhi)
1345 ,fNameQAJ(copy.fNameQAJ)
1350 //________________________________________________________________________________________________________________________________________________________________________
1351 AliAnalysisTaskIDFragmentationFunction::AliFragFuncQAJetHistos& AliAnalysisTaskIDFragmentationFunction::AliFragFuncQAJetHistos::operator=(const AliAnalysisTaskIDFragmentationFunction::AliFragFuncQAJetHistos& o)
1356 TObject::operator=(o);
1357 fNBinsPt = o.fNBinsPt;
1360 fNBinsEta = o.fNBinsEta;
1361 fEtaMin = o.fEtaMin;
1362 fEtaMax = o.fEtaMax;
1363 fNBinsPhi = o.fNBinsPhi;
1364 fPhiMin = o.fPhiMin;
1365 fPhiMax = o.fPhiMax;
1366 fh2EtaPhi = o.fh2EtaPhi;
1368 fNameQAJ = o.fNameQAJ;
1374 //______________________________________________________________
1375 AliAnalysisTaskIDFragmentationFunction::AliFragFuncQAJetHistos::~AliFragFuncQAJetHistos()
1379 if(fh2EtaPhi) delete fh2EtaPhi;
1380 if(fh1Pt) delete fh1Pt;
1383 //____________________________________________________________________
1384 void AliAnalysisTaskIDFragmentationFunction::AliFragFuncQAJetHistos::DefineHistos()
1386 // book jet QA histos
1388 fh2EtaPhi = new TH2F(Form("fh2JetQAEtaPhi%s", fNameQAJ.Data()), Form("%s: #eta - #phi distribution", fNameQAJ.Data()), fNBinsEta, fEtaMin, fEtaMax, fNBinsPhi, fPhiMin, fPhiMax);
1389 fh1Pt = new TH1F(Form("fh1JetQAPt%s", fNameQAJ.Data()), Form("%s: p_{T} distribution", fNameQAJ.Data()), fNBinsPt, fPtMin, fPtMax);
1391 AliAnalysisTaskIDFragmentationFunction::SetProperties(fh2EtaPhi, "#eta", "#phi");
1392 AliAnalysisTaskIDFragmentationFunction::SetProperties(fh1Pt, "p_{T} [GeV/c]", "entries");
1395 //____________________________________________________________________________________________________
1396 void AliAnalysisTaskIDFragmentationFunction::AliFragFuncQAJetHistos::FillJetQA(Float_t eta, Float_t phi, Float_t pt)
1398 // fill jet QA histos
1400 fh2EtaPhi->Fill( eta, phi);
1404 //____________________________________________________________________________________
1405 void AliAnalysisTaskIDFragmentationFunction::AliFragFuncQAJetHistos::AddToOutput(TList* list) const
1407 // add histos to list
1409 list->Add(fh2EtaPhi);
1413 //___________________________________________________________________________________________________________
1414 AliAnalysisTaskIDFragmentationFunction::AliFragFuncQATrackHistos::AliFragFuncQATrackHistos(const char* name,
1415 Int_t nPt, Float_t ptMin, Float_t ptMax,
1416 Int_t nEta, Float_t etaMin, Float_t etaMax,
1417 Int_t nPhi, Float_t phiMin, Float_t phiMax,
1429 ,fHighPtThreshold(ptThresh)
1436 // default constructor
1439 //__________________________________________________________________________________________
1440 AliAnalysisTaskIDFragmentationFunction::AliFragFuncQATrackHistos::AliFragFuncQATrackHistos(const AliFragFuncQATrackHistos& copy)
1442 ,fNBinsPt(copy.fNBinsPt)
1443 ,fPtMin(copy.fPtMin)
1444 ,fPtMax(copy.fPtMax)
1445 ,fNBinsEta(copy.fNBinsEta)
1446 ,fEtaMin(copy.fEtaMin)
1447 ,fEtaMax(copy.fEtaMax)
1448 ,fNBinsPhi(copy.fNBinsPhi)
1449 ,fPhiMin(copy.fPhiMin)
1450 ,fPhiMax(copy.fPhiMax)
1451 ,fHighPtThreshold(copy.fHighPtThreshold)
1452 ,fh2EtaPhi(copy.fh2EtaPhi)
1454 ,fh2HighPtEtaPhi(copy.fh2HighPtEtaPhi)
1455 ,fh2PhiPt(copy.fh2PhiPt)
1456 ,fNameQAT(copy.fNameQAT)
1461 // _____________________________________________________________________________________________________________________________________________________________________________
1462 AliAnalysisTaskIDFragmentationFunction::AliFragFuncQATrackHistos& AliAnalysisTaskIDFragmentationFunction::AliFragFuncQATrackHistos::operator=(const AliAnalysisTaskIDFragmentationFunction::AliFragFuncQATrackHistos& o)
1467 TObject::operator=(o);
1468 fNBinsPt = o.fNBinsPt;
1471 fNBinsEta = o.fNBinsEta;
1472 fEtaMin = o.fEtaMin;
1473 fEtaMax = o.fEtaMax;
1474 fNBinsPhi = o.fNBinsPhi;
1475 fPhiMin = o.fPhiMin;
1476 fPhiMax = o.fPhiMax;
1477 fHighPtThreshold = o.fHighPtThreshold;
1478 fh2EtaPhi = o.fh2EtaPhi;
1480 fh2HighPtEtaPhi = o.fh2HighPtEtaPhi;
1481 fh2PhiPt = o.fh2PhiPt;
1482 fNameQAT = o.fNameQAT;
1488 //___________________________________________________________________
1489 AliAnalysisTaskIDFragmentationFunction::AliFragFuncQATrackHistos::~AliFragFuncQATrackHistos()
1493 if(fh2EtaPhi) delete fh2EtaPhi;
1494 if(fh2HighPtEtaPhi) delete fh2HighPtEtaPhi;
1495 if(fh1Pt) delete fh1Pt;
1496 if(fh2PhiPt) delete fh2PhiPt;
1499 //______________________________________________________________________
1500 void AliAnalysisTaskIDFragmentationFunction::AliFragFuncQATrackHistos::DefineHistos()
1502 // book track QA histos
1504 fh2EtaPhi = new TH2F(Form("fh2TrackQAEtaPhi%s", fNameQAT.Data()), Form("%s: #eta - #phi distribution", fNameQAT.Data()), fNBinsEta, fEtaMin, fEtaMax, fNBinsPhi, fPhiMin, fPhiMax);
1505 fh2HighPtEtaPhi = new TH2F(Form("fh2TrackQAHighPtEtaPhi%s", fNameQAT.Data()), Form("%s: #eta - #phi distribution for high-p_{T}", fNameQAT.Data()), fNBinsEta, fEtaMin, fEtaMax, fNBinsPhi, fPhiMin, fPhiMax);
1506 fh1Pt = new TH1F(Form("fh1TrackQAPt%s", fNameQAT.Data()), Form("%s: p_{T} distribution", fNameQAT.Data()), fNBinsPt, fPtMin, fPtMax);
1507 fh2PhiPt = new TH2F(Form("fh2TrackQAPhiPt%s", fNameQAT.Data()), Form("%s: #phi - p_{T} distribution", fNameQAT.Data()), fNBinsPhi, fPhiMin, fPhiMax, fNBinsPt, fPtMin, fPtMax);
1509 AliAnalysisTaskIDFragmentationFunction::SetProperties(fh2EtaPhi, "#eta", "#phi");
1510 AliAnalysisTaskIDFragmentationFunction::SetProperties(fh2HighPtEtaPhi, "#eta", "#phi");
1511 AliAnalysisTaskIDFragmentationFunction::SetProperties(fh1Pt, "p_{T} [GeV/c]", "entries");
1512 AliAnalysisTaskIDFragmentationFunction::SetProperties(fh2PhiPt, "#phi", "p_{T} [GeV/c]");
1515 //________________________________________________________________________________________________________
1516 void AliAnalysisTaskIDFragmentationFunction::AliFragFuncQATrackHistos::FillTrackQA(Float_t eta, Float_t phi, Float_t pt, Bool_t weightPt, Float_t norm,
1517 Bool_t scaleStrangeness, Float_t scaleFacStrangeness)
1519 // fill track QA histos
1520 Float_t weight = 1.;
1521 if(weightPt) weight = pt;
1522 fh2EtaPhi->Fill( eta, phi, weight);
1523 if(scaleStrangeness) fh2EtaPhi->Fill( eta, phi, scaleFacStrangeness);
1524 if(pt > fHighPtThreshold) fh2HighPtEtaPhi->Fill( eta, phi, weight);
1525 if(pt > fHighPtThreshold && scaleStrangeness) fh2HighPtEtaPhi->Fill( eta, phi, weight);
1526 if(norm) fh1Pt->Fill( pt, 1/norm );
1527 else if(scaleStrangeness) fh1Pt->Fill(pt,scaleFacStrangeness);
1528 else fh1Pt->Fill( pt );
1530 if(scaleFacStrangeness) fh2PhiPt->Fill(phi, pt, scaleFacStrangeness);
1531 else fh2PhiPt->Fill(phi, pt);
1534 //______________________________________________________________________________________
1535 void AliAnalysisTaskIDFragmentationFunction::AliFragFuncQATrackHistos::AddToOutput(TList* list) const
1537 // add histos to list
1539 list->Add(fh2EtaPhi);
1540 list->Add(fh2HighPtEtaPhi);
1542 list->Add(fh2PhiPt);
1545 //_________________________________________________________________________________
1546 Bool_t AliAnalysisTaskIDFragmentationFunction::Notify()
1549 // Implemented Notify() to read the cross sections
1550 // and number of trials from pyxsec.root
1551 // (taken from AliAnalysisTaskJetSpectrum2)
1553 TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
1554 Float_t xsection = 0;
1555 Float_t ftrials = 1;
1559 TFile *curfile = tree->GetCurrentFile();
1561 Error("Notify","No current file");
1565 AliAnalysisHelperJetTasks::PythiaInfoFromFile(curfile->GetName(),xsection,ftrials);
1567 if (fUseInclusivePIDtask) {
1568 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++)
1569 fInclusivePIDtask[i]->FillXsec(xsection);
1572 if (fUseJetPIDtask) {
1573 for (Int_t i = 0; i < fNumJetPIDtasks; i++)
1574 fJetPIDtask[i]->FillXsec(xsection);
1577 if(!fh1Xsec||!fh1Trials){
1578 Printf("%s:%d No Histogram fh1Xsec",(char*)__FILE__,__LINE__);
1582 fh1Xsec->Fill("<#sigma>",xsection);
1583 // construct a poor man average trials
1584 Float_t nEntries = (Float_t)tree->GetTree()->GetEntries();
1585 if(ftrials>=nEntries && nEntries>0.)fAvgTrials = ftrials/nEntries;
1588 // Set seed for backg study
1590 fRandom = new TRandom3();
1591 fRandom->SetSeed(0);
1596 //__________________________________________________________________
1597 void AliAnalysisTaskIDFragmentationFunction::UserCreateOutputObjects()
1599 // create output objects
1601 if(fDebug > 1) Printf("AliAnalysisTaskIDFragmentationFunction::UserCreateOutputObjects()");
1603 // create list of tracks and jets
1605 fTracksRecCuts = new TList();
1606 fTracksRecCuts->SetOwner(kFALSE);
1608 fTracksRecCutsEfficiency = new TList();
1609 fTracksRecCutsEfficiency->SetOwner(kFALSE);
1611 fTracksGen = new TList();
1612 fTracksGen->SetOwner(kFALSE);
1614 fTracksAODMCCharged = new TList();
1615 fTracksAODMCCharged->SetOwner(kFALSE);
1617 fTracksAODMCChargedSecNS = new TList();
1618 fTracksAODMCChargedSecNS->SetOwner(kFALSE);
1620 fTracksAODMCChargedSecS = new TList();
1621 fTracksAODMCChargedSecS->SetOwner(kFALSE);
1623 fTracksRecQualityCuts = new TList();
1624 fTracksRecQualityCuts->SetOwner(kFALSE);
1626 fJetsRec = new TList();
1627 fJetsRec->SetOwner(kFALSE);
1629 fJetsRecCuts = new TList();
1630 fJetsRecCuts->SetOwner(kFALSE);
1632 fJetsGen = new TList();
1633 fJetsGen->SetOwner(kFALSE);
1635 fJetsRecEff = new TList();
1636 fJetsRecEff->SetOwner(kFALSE);
1638 fJetsEmbedded = new TList();
1639 fJetsEmbedded->SetOwner(kFALSE);
1643 (fBckgType[0]==kBckgClusters || fBckgType[1]==kBckgClusters || fBckgType[2]==kBckgClusters || fBckgType[3]==kBckgClusters || fBckgType[4]==kBckgClusters ||
1644 fBckgType[0]==kBckgClustersOutLeading || fBckgType[1]==kBckgClustersOutLeading || fBckgType[2]==kBckgClustersOutLeading ||
1645 fBckgType[3]==kBckgClustersOutLeading || fBckgType[4]==kBckgClustersOutLeading)){
1647 fBckgJetsRec = new TList();
1648 fBckgJetsRec->SetOwner(kFALSE);
1650 fBckgJetsRecCuts = new TList();
1651 fBckgJetsRecCuts->SetOwner(kFALSE);
1653 fBckgJetsGen = new TList();
1654 fBckgJetsGen->SetOwner(kFALSE);
1658 // Create histograms / output container
1662 fCommonHistList = new TList();
1663 fCommonHistList->SetOwner(kTRUE);
1665 Bool_t oldStatus = TH1::AddDirectoryStatus();
1666 TH1::AddDirectory(kFALSE);
1670 fh1EvtSelection = new TH1F("fh1EvtSelection", "Event Selection", 6, -0.5, 5.5);
1671 fh1EvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
1672 fh1EvtSelection->GetXaxis()->SetBinLabel(2,"event selection: rejected");
1673 fh1EvtSelection->GetXaxis()->SetBinLabel(3,"event class: rejected");
1674 fh1EvtSelection->GetXaxis()->SetBinLabel(4,"vertex Ncontr: rejected");
1675 fh1EvtSelection->GetXaxis()->SetBinLabel(5,"vertex z: rejected");
1676 fh1EvtSelection->GetXaxis()->SetBinLabel(6,"vertex type: rejected");
1678 fh1VtxSelection = new TH1F("fh1VtxSelection", "Vertex Selection", 10, -1, 9);
1679 fh1VtxSelection->GetXaxis()->SetBinLabel(1,"Undef");
1680 fh1VtxSelection->GetXaxis()->SetBinLabel(2,"Primary");
1681 fh1VtxSelection->GetXaxis()->SetBinLabel(3,"Kink");
1682 fh1VtxSelection->GetXaxis()->SetBinLabel(4,"V0");
1683 fh1VtxSelection->GetXaxis()->SetBinLabel(5,"Cascade");
1684 fh1VtxSelection->GetXaxis()->SetBinLabel(6,"Multi");
1685 fh1VtxSelection->GetXaxis()->SetBinLabel(7,"SPD");
1686 fh1VtxSelection->GetXaxis()->SetBinLabel(8,"PileUpSPD");
1687 fh1VtxSelection->GetXaxis()->SetBinLabel(9,"PileUpTracks");
1688 fh1VtxSelection->GetXaxis()->SetBinLabel(10,"TPC");
1690 fh1VertexNContributors = new TH1F("fh1VertexNContributors", "Vertex N contributors", 2500,-.5, 2499.5);
1691 fh1VertexZ = new TH1F("fh1VertexZ", "Vertex z distribution", 30, -15., 15.);
1692 fh1EvtMult = new TH1F("fh1EvtMult","Event multiplicity, track pT cut > 150 MeV/c, |#eta| < 0.9",120,0.,12000.);
1693 fh1EvtCent = new TH1F("fh1EvtCent","centrality",100,0.,100.);
1695 fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
1696 fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
1697 fh1Trials = new TH1F("fh1Trials","trials from pyxsec.root",1,0,1);
1698 fh1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
1699 fh1PtHard = new TH1F("fh1PtHard","PYTHIA Pt hard;p_{T,hard}",350,-.5,349.5);
1700 fh1PtHardTrials = new TH1F("fh1PtHardTrials","PYTHIA Pt hard weight with trials;p_{T,hard}",350,-.5,349.5);
1702 fh1nRecJetsCuts = new TH1F("fh1nRecJetsCuts","reconstructed jets per event",10,-0.5,9.5);
1703 fh1nGenJets = new TH1F("fh1nGenJets","generated jets per event",10,-0.5,9.5);
1704 fh1nRecEffJets = new TH1F("fh1nRecEffJets","reconstruction effiency: jets per event",10,-0.5,9.5);
1705 fh1nEmbeddedJets = new TH1F("fh1nEmbeddedJets","embedded jets per event",10,-0.5,9.5);
1707 fh2PtRecVsGenPrim = new TH2F("fh2PtRecVsGenPrim","rec vs gen pt",fQATrackNBinsPt,fQATrackPtMin,fQATrackPtMax,fQATrackNBinsPt,fQATrackPtMin,fQATrackPtMax);
1708 fh2PtRecVsGenSec = new TH2F("fh2PtRecVsGenSec","rec vs gen pt",fQATrackNBinsPt,fQATrackPtMin,fQATrackPtMax,fQATrackNBinsPt,fQATrackPtMin,fQATrackPtMax);
1711 if(fBranchEmbeddedJets.Length()){
1712 fh1FractionPtEmbedded = new TH1F("fh1FractionPtEmbedded","",200,0,2);
1713 fh1IndexEmbedded = new TH1F("fh1IndexEmbedded","",11,-1,10);
1714 fh2DeltaPtVsJetPtEmbedded = new TH2F("fh2DeltaPtVsJetPtEmbedded","",250,0,250,200,-100,100);
1715 fh2DeltaPtVsRecJetPtEmbedded = new TH2F("fh2DeltaPtVsRecJetPtEmbedded","",250,0,250,200,-100,100);
1716 fh1DeltaREmbedded = new TH1F("fh1DeltaREmbedded","",50,0,0.5);
1717 fh1nEmbeddedJets = new TH1F("fh1nEmbeddedJets","embedded jets per event",10,-0.5,9.5);
1722 if(fQAMode&1){ // track QA
1723 fQATrackHistosRecCuts = new AliFragFuncQATrackHistos("RecCuts", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1724 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1725 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1726 fQATrackHighPtThreshold);
1727 fQATrackHistosGen = new AliFragFuncQATrackHistos("Gen", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1728 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1729 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1730 fQATrackHighPtThreshold);
1733 if(fQAMode&2){ // jet QA
1734 fQAJetHistosRec = new AliFragFuncQAJetHistos("Rec", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax,
1735 fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,
1736 fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1737 fQAJetHistosRecCuts = new AliFragFuncQAJetHistos("RecCuts", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax,
1738 fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,
1739 fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1740 fQAJetHistosRecCutsLeading = new AliFragFuncQAJetHistos("RecCutsLeading", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax,
1741 fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,
1742 fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1743 fQAJetHistosGen = new AliFragFuncQAJetHistos("Gen", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax,
1744 fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,
1745 fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1746 fQAJetHistosGenLeading = new AliFragFuncQAJetHistos("GenLeading", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax,
1747 fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,
1748 fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1749 if(fEffMode) fQAJetHistosRecEffLeading = new AliFragFuncQAJetHistos("RecEffLeading", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax,
1750 fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1754 if(fFFMode || fIDFFMode){
1756 fFFHistosRecCuts = new AliFragFuncHistos("RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1757 fFFNBinsPt, fFFPtMin, fFFPtMax,
1758 fFFNBinsXi, fFFXiMin, fFFXiMax,
1759 fFFNBinsZ , fFFZMin , fFFZMax );
1762 fFFHistosRecCutsInc = new AliFragFuncHistos("RecCutsInc", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1763 fFFNBinsPt, fFFPtMin, fFFPtMax,
1764 fFFNBinsXi, fFFXiMin, fFFXiMax,
1765 fFFNBinsZ , fFFZMin , fFFZMax );
1768 fFFHistosRecLeadingTrack = new AliFragFuncHistos("RecLeadingTrack", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1769 fFFNBinsPt, fFFPtMin, fFFPtMax,
1770 fFFNBinsXi, fFFXiMin, fFFXiMax,
1771 fFFNBinsZ , fFFZMin , fFFZMax );
1773 fFFHistosGen = new AliFragFuncHistos("Gen", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1774 fFFNBinsPt, fFFPtMin, fFFPtMax,
1775 fFFNBinsXi, fFFXiMin, fFFXiMax,
1776 fFFNBinsZ , fFFZMin , fFFZMax);
1778 fFFHistosGenInc = new AliFragFuncHistos("GenInc", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1779 fFFNBinsPt, fFFPtMin, fFFPtMax,
1780 fFFNBinsXi, fFFXiMin, fFFXiMax,
1781 fFFNBinsZ , fFFZMin , fFFZMax);
1783 fFFHistosGenLeadingTrack = new AliFragFuncHistos("GenLeadingTrack", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1784 fFFNBinsPt, fFFPtMin, fFFPtMax,
1785 fFFNBinsXi, fFFXiMin, fFFXiMax,
1786 fFFNBinsZ , fFFZMin , fFFZMax);
1789 for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
1790 fIDFFHistosRecCuts[i] = new AliFragFuncHistos(Form("RecCuts_%s", AliPID::ParticleName(i)), fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1791 fFFNBinsPt, fFFPtMin, fFFPtMax,
1792 fFFNBinsXi, fFFXiMin, fFFXiMax,
1793 fFFNBinsZ , fFFZMin , fFFZMax );
1794 fIDFFHistosGen[i] = new AliFragFuncHistos(Form("Gen_%s", AliPID::ParticleName(i)), fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1795 fFFNBinsPt, fFFPtMin, fFFPtMax,
1796 fFFNBinsXi, fFFXiMin, fFFXiMax,
1797 fFFNBinsZ , fFFZMin , fFFZMax );
1806 fQATrackHistosRecEffGen = new AliFragFuncQATrackHistos("RecEffGen", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1807 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1808 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1809 fQATrackHighPtThreshold);
1811 fQATrackHistosRecEffRec = new AliFragFuncQATrackHistos("RecEffRec", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1812 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1813 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1814 fQATrackHighPtThreshold);
1816 fQATrackHistosSecRecNS = new AliFragFuncQATrackHistos("SecRecNS", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1817 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1818 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1819 fQATrackHighPtThreshold);
1821 fQATrackHistosSecRecS = new AliFragFuncQATrackHistos("SecRecS", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1822 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1823 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1824 fQATrackHighPtThreshold);
1826 fQATrackHistosSecRecSsc = new AliFragFuncQATrackHistos("SecRecSsc", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1827 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1828 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1829 fQATrackHighPtThreshold);
1833 fFFHistosRecEffRec = new AliFragFuncHistos("RecEffRec", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1834 fFFNBinsPt, fFFPtMin, fFFPtMax,
1835 fFFNBinsXi, fFFXiMin, fFFXiMax,
1836 fFFNBinsZ , fFFZMin , fFFZMax);
1838 fFFHistosSecRecNS = new AliFragFuncHistos("SecRecNS", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1839 fFFNBinsPt, fFFPtMin, fFFPtMax,
1840 fFFNBinsXi, fFFXiMin, fFFXiMax,
1841 fFFNBinsZ , fFFZMin , fFFZMax);
1843 fFFHistosSecRecS = new AliFragFuncHistos("SecRecS", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1844 fFFNBinsPt, fFFPtMin, fFFPtMax,
1845 fFFNBinsXi, fFFXiMin, fFFXiMax,
1846 fFFNBinsZ , fFFZMin , fFFZMax);
1848 fFFHistosSecRecSsc = new AliFragFuncHistos("SecRecSsc", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1849 fFFNBinsPt, fFFPtMin, fFFPtMax,
1850 fFFNBinsXi, fFFXiMin, fFFXiMax,
1851 fFFNBinsZ , fFFZMin , fFFZMax);
1854 } // end: efficiency
1858 if(fBckgType[0]==kBckgNone){
1859 AliError("no bgr method selected !");
1863 for(Int_t i=0; i<5; i++){
1864 if(fBckgType[i]==kBckgPerp) title[i]="Perp";
1865 else if(fBckgType[i]==kBckgPerp2) title[i]="Perp2";
1866 else if(fBckgType[i]==kBckgPerp2Area) title[i]="Perp2Area";
1867 else if(fBckgType[i]==kBckgPerpWindow) title[i]="PerpW";
1868 else if(fBckgType[i]==kBckgASide) title[i]="ASide";
1869 else if(fBckgType[i]==kBckgASideWindow) title[i]="ASideW";
1870 else if(fBckgType[i]==kBckgOutLJ) title[i]="OutLeadingJet";
1871 else if(fBckgType[i]==kBckgOut2J) title[i]="Out2Jets";
1872 else if(fBckgType[i]==kBckgOut3J) title[i]="Out3Jets";
1873 else if(fBckgType[i]==kBckgOutAJ) title[i]="AllJets";
1874 else if(fBckgType[i]==kBckgOutLJStat) title[i]="OutLeadingJetStat";
1875 else if(fBckgType[i]==kBckgOut2JStat) title[i]="Out2JetsStat";
1876 else if(fBckgType[i]==kBckgOut3JStat) title[i]="Out3JetsStat";
1877 else if(fBckgType[i]==kBckgOutAJStat) title[i]="AllJetsStat";
1878 else if(fBckgType[i]==kBckgClustersOutLeading) title[i]="OutClusters";
1879 else if(fBckgType[i]==kBckgClusters) title[i]="MedianClusters";
1880 else if(fBckgType[i]==kBckgNone) title[i]="";
1881 else printf("Please chose background method number %d!",i);
1885 if(fBckgType[0]==kBckgClusters || fBckgType[1]==kBckgClusters || fBckgType[2]==kBckgClusters || fBckgType[3]==kBckgClusters || fBckgType[4]==kBckgClusters ||
1886 fBckgType[0]==kBckgClustersOutLeading || fBckgType[1]==kBckgClustersOutLeading || fBckgType[2]==kBckgClustersOutLeading ||
1887 fBckgType[3]==kBckgClustersOutLeading || fBckgType[4]==kBckgClustersOutLeading){
1889 fh1nRecBckgJetsCuts = new TH1F("fh1nRecBckgJetsCuts","reconstructed background jets per event",10,-0.5,9.5);
1890 fh1nGenBckgJets = new TH1F("fh1nGenBckgJets","generated background jets per event",10,-0.5,9.5);
1894 fh1BckgMult0 = new TH1F("fh1BckgMult0","bckg mult "+title[0],500,0,500);
1895 if(fBckgType[1] != kBckgNone) fh1BckgMult1 = new TH1F("fh1BckgMult1","bckg mult "+title[1],500,0,500);
1896 if(fBckgType[2] != kBckgNone) fh1BckgMult2 = new TH1F("fh1BckgMult2","bckg mult "+title[2],500,0,500);
1897 if(fBckgType[3] != kBckgNone) fh1BckgMult3 = new TH1F("fh1BckgMult3","bckg mult "+title[3],500,0,500);
1898 if(fBckgType[4] != kBckgNone) fh1BckgMult4 = new TH1F("fh1BckgMult4","bckg mult "+title[4],500,0,500);
1902 fQABckgHisto0RecCuts = new AliFragFuncQATrackHistos("Bckg"+title[0]+"RecCuts", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1903 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1904 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1905 fQATrackHighPtThreshold);
1906 fQABckgHisto0Gen = new AliFragFuncQATrackHistos("Bckg"+title[0]+"Gen", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1907 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1908 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1909 fQATrackHighPtThreshold);
1911 if(fBckgType[1] != kBckgNone){
1912 fQABckgHisto1RecCuts = new AliFragFuncQATrackHistos("Bckg"+title[1]+"RecCuts", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1913 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1914 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1915 fQATrackHighPtThreshold);
1916 fQABckgHisto1Gen = new AliFragFuncQATrackHistos("Bckg"+title[1]+"Gen", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1917 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1918 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1919 fQATrackHighPtThreshold);
1921 if(fBckgType[2] != kBckgNone){
1922 fQABckgHisto2RecCuts = new AliFragFuncQATrackHistos("Bckg"+title[2]+"RecCuts", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1923 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1924 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1925 fQATrackHighPtThreshold);
1926 fQABckgHisto2Gen = new AliFragFuncQATrackHistos("Bckg"+title[2]+"Gen", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1927 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1928 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1929 fQATrackHighPtThreshold);
1931 if(fBckgType[3] != kBckgNone){
1932 fQABckgHisto3RecCuts = new AliFragFuncQATrackHistos("Bckg"+title[3]+"RecCuts", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1933 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1934 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1935 fQATrackHighPtThreshold);
1936 fQABckgHisto3Gen = new AliFragFuncQATrackHistos("Bckg"+title[3]+"Gen", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1937 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1938 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1939 fQATrackHighPtThreshold);
1941 if(fBckgType[4] != kBckgNone){
1942 fQABckgHisto4RecCuts = new AliFragFuncQATrackHistos("Bckg"+title[4]+"RecCuts", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1943 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1944 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1945 fQATrackHighPtThreshold);
1946 fQABckgHisto4Gen = new AliFragFuncQATrackHistos("Bckg"+title[4]+"Gen", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1947 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1948 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1949 fQATrackHighPtThreshold);
1951 } // end: background QA
1954 fFFBckgHisto0RecCuts = new AliFragFuncHistos("Bckg"+title[0]+"RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1955 fFFNBinsPt, fFFPtMin, fFFPtMax,
1956 fFFNBinsXi, fFFXiMin, fFFXiMax,
1957 fFFNBinsZ , fFFZMin , fFFZMax);
1959 fFFBckgHisto0Gen = new AliFragFuncHistos("Bckg"+title[0]+"Gen", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1960 fFFNBinsPt, fFFPtMin, fFFPtMax,
1961 fFFNBinsXi, fFFXiMin, fFFXiMax,
1962 fFFNBinsZ , fFFZMin , fFFZMax);
1964 if(fBckgType[1] != kBckgNone){
1965 fFFBckgHisto1RecCuts = new AliFragFuncHistos("Bckg"+title[1]+"RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1966 fFFNBinsPt, fFFPtMin, fFFPtMax,
1967 fFFNBinsXi, fFFXiMin, fFFXiMax,
1968 fFFNBinsZ , fFFZMin , fFFZMax);
1969 fFFBckgHisto1Gen = new AliFragFuncHistos("Bckg"+title[1]+"Gen", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1970 fFFNBinsPt, fFFPtMin, fFFPtMax,
1971 fFFNBinsXi, fFFXiMin, fFFXiMax,
1972 fFFNBinsZ , fFFZMin , fFFZMax);
1974 if(fBckgType[2] != kBckgNone){
1975 fFFBckgHisto2RecCuts = new AliFragFuncHistos("Bckg"+title[2]+"RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1976 fFFNBinsPt, fFFPtMin, fFFPtMax,
1977 fFFNBinsXi, fFFXiMin, fFFXiMax,
1978 fFFNBinsZ , fFFZMin , fFFZMax);
1980 fFFBckgHisto2Gen = new AliFragFuncHistos("Bckg"+title[2]+"Gen", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1981 fFFNBinsPt, fFFPtMin, fFFPtMax,
1982 fFFNBinsXi, fFFXiMin, fFFXiMax,
1983 fFFNBinsZ , fFFZMin , fFFZMax);
1985 if(fBckgType[3] != kBckgNone){
1986 fFFBckgHisto3RecCuts = new AliFragFuncHistos("Bckg"+title[3]+"RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1987 fFFNBinsPt, fFFPtMin, fFFPtMax,
1988 fFFNBinsXi, fFFXiMin, fFFXiMax,
1989 fFFNBinsZ , fFFZMin , fFFZMax);
1991 fFFBckgHisto3Gen = new AliFragFuncHistos("Bckg"+title[3]+"Gen", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1992 fFFNBinsPt, fFFPtMin, fFFPtMax,
1993 fFFNBinsXi, fFFXiMin, fFFXiMax,
1994 fFFNBinsZ , fFFZMin , fFFZMax);
1996 if(fBckgType[4] != kBckgNone){
1997 fFFBckgHisto4RecCuts = new AliFragFuncHistos("Bckg"+title[4]+"RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1998 fFFNBinsPt, fFFPtMin, fFFPtMax,
1999 fFFNBinsXi, fFFXiMin, fFFXiMax,
2000 fFFNBinsZ , fFFZMin , fFFZMax);
2002 fFFBckgHisto4Gen = new AliFragFuncHistos("Bckg"+title[4]+"Gen", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
2003 fFFNBinsPt, fFFPtMin, fFFPtMax,
2004 fFFNBinsXi, fFFXiMin, fFFXiMax,
2005 fFFNBinsZ , fFFZMin , fFFZMax);
2008 fFFBckgHisto0RecEffRec = new AliFragFuncHistos("Bckg"+title[0]+"RecEffRec", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
2009 fFFNBinsPt, fFFPtMin, fFFPtMax,
2010 fFFNBinsXi, fFFXiMin, fFFXiMax,
2011 fFFNBinsZ , fFFZMin , fFFZMax);
2013 fFFBckgHisto0SecRecNS = new AliFragFuncHistos("Bckg"+title[0]+"SecRecNS", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
2014 fFFNBinsPt, fFFPtMin, fFFPtMax,
2015 fFFNBinsXi, fFFXiMin, fFFXiMax,
2016 fFFNBinsZ , fFFZMin , fFFZMax);
2018 fFFBckgHisto0SecRecS = new AliFragFuncHistos("Bckg"+title[0]+"SecRecS", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
2019 fFFNBinsPt, fFFPtMin, fFFPtMax,
2020 fFFNBinsXi, fFFXiMin, fFFXiMax,
2021 fFFNBinsZ , fFFZMin , fFFZMax);
2023 fFFBckgHisto0SecRecSsc = new AliFragFuncHistos("Bckg"+title[0]+"SecRecSsc", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
2024 fFFNBinsPt, fFFPtMin, fFFPtMax,
2025 fFFNBinsXi, fFFXiMin, fFFXiMax,
2026 fFFNBinsZ , fFFZMin , fFFZMax);
2029 } // end: background FF
2032 } // end: background
2035 // ____________ define histograms ____________________
2038 if(fQAMode&1){ // track QA
2039 fQATrackHistosRecCuts->DefineHistos();
2040 fQATrackHistosGen->DefineHistos();
2043 if(fQAMode&2){ // jet QA
2044 fQAJetHistosRec->DefineHistos();
2045 fQAJetHistosRecCuts->DefineHistos();
2046 fQAJetHistosRecCutsLeading->DefineHistos();
2047 fQAJetHistosGen->DefineHistos();
2048 fQAJetHistosGenLeading->DefineHistos();
2049 if(fEffMode) fQAJetHistosRecEffLeading->DefineHistos();
2053 if(fFFMode || fIDFFMode){
2054 fFFHistosRecCuts->DefineHistos();
2055 fFFHistosRecCutsInc->DefineHistos();
2056 fFFHistosRecLeadingTrack->DefineHistos();
2057 fFFHistosGen->DefineHistos();
2058 fFFHistosGenInc->DefineHistos();
2059 fFFHistosGenLeadingTrack->DefineHistos();
2062 for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
2063 fIDFFHistosRecCuts[i]->DefineHistos();
2064 fIDFFHistosGen[i]->DefineHistos();
2071 fQATrackHistosRecEffGen->DefineHistos();
2072 fQATrackHistosRecEffRec->DefineHistos();
2073 fQATrackHistosSecRecNS->DefineHistos();
2074 fQATrackHistosSecRecS->DefineHistos();
2075 fQATrackHistosSecRecSsc->DefineHistos();
2078 fFFHistosRecEffRec->DefineHistos();
2079 fFFHistosSecRecNS->DefineHistos();
2080 fFFHistosSecRecS->DefineHistos();
2081 fFFHistosSecRecSsc->DefineHistos();
2083 } // end: efficiency
2088 fFFBckgHisto0RecCuts->DefineHistos();
2089 fFFBckgHisto0Gen->DefineHistos();
2090 if(fBckgType[1] != kBckgNone) fFFBckgHisto1RecCuts->DefineHistos();
2091 if(fBckgType[1] != kBckgNone) fFFBckgHisto1Gen->DefineHistos();
2092 if(fBckgType[2] != kBckgNone) fFFBckgHisto2RecCuts->DefineHistos();
2093 if(fBckgType[2] != kBckgNone) fFFBckgHisto2Gen->DefineHistos();
2094 if(fBckgType[3] != kBckgNone) fFFBckgHisto3RecCuts->DefineHistos();
2095 if(fBckgType[3] != kBckgNone) fFFBckgHisto3Gen->DefineHistos();
2096 if(fBckgType[4] != kBckgNone) fFFBckgHisto4RecCuts->DefineHistos();
2097 if(fBckgType[4] != kBckgNone) fFFBckgHisto4Gen->DefineHistos();
2100 fFFBckgHisto0RecEffRec->DefineHistos();
2101 fFFBckgHisto0SecRecNS->DefineHistos();
2102 fFFBckgHisto0SecRecS->DefineHistos();
2103 fFFBckgHisto0SecRecSsc->DefineHistos();
2108 fQABckgHisto0RecCuts->DefineHistos();
2109 fQABckgHisto0Gen->DefineHistos();
2110 if(fBckgType[1] != kBckgNone) fQABckgHisto1RecCuts->DefineHistos();
2111 if(fBckgType[1] != kBckgNone) fQABckgHisto1Gen->DefineHistos();
2112 if(fBckgType[2] != kBckgNone) fQABckgHisto2RecCuts->DefineHistos();
2113 if(fBckgType[2] != kBckgNone) fQABckgHisto2Gen->DefineHistos();
2114 if(fBckgType[3] != kBckgNone) fQABckgHisto3RecCuts->DefineHistos();
2115 if(fBckgType[3] != kBckgNone) fQABckgHisto3Gen->DefineHistos();
2116 if(fBckgType[4] != kBckgNone) fQABckgHisto4RecCuts->DefineHistos();
2117 if(fBckgType[4] != kBckgNone) fQABckgHisto4Gen->DefineHistos();
2119 } // end: background
2122 Bool_t genJets = (fJetTypeGen != kJetsUndef) ? kTRUE : kFALSE;
2123 Bool_t genTracks = (fTrackTypeGen != kTrackUndef) ? kTRUE : kFALSE;
2124 Bool_t recJetsEff = (fJetTypeRecEff != kJetsUndef) ? kTRUE : kFALSE;
2126 fCommonHistList->Add(fh1EvtSelection);
2127 fCommonHistList->Add(fh1VtxSelection);
2128 fCommonHistList->Add(fh1EvtMult);
2129 fCommonHistList->Add(fh1EvtCent);
2130 fCommonHistList->Add(fh1VertexNContributors);
2131 fCommonHistList->Add(fh1VertexZ);
2132 fCommonHistList->Add(fh1nRecJetsCuts);
2133 fCommonHistList->Add(fh1Xsec);
2134 fCommonHistList->Add(fh1Trials);
2135 fCommonHistList->Add(fh1PtHard);
2136 fCommonHistList->Add(fh1PtHardTrials);
2138 if(genJets) fCommonHistList->Add(fh1nGenJets);
2142 fFFHistosRecCuts->AddToOutput(fCommonHistList);
2143 fFFHistosRecCutsInc->AddToOutput(fCommonHistList);
2144 fFFHistosRecLeadingTrack->AddToOutput(fCommonHistList);
2146 if(genJets && genTracks){
2147 fFFHistosGen->AddToOutput(fCommonHistList);
2148 fFFHistosGenInc->AddToOutput(fCommonHistList);
2149 fFFHistosGenLeadingTrack->AddToOutput(fCommonHistList);
2153 for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
2154 if(genJets && genTracks)
2155 fIDFFHistosGen[i]->AddToOutput(fCommonHistList);
2156 fIDFFHistosRecCuts[i]->AddToOutput(fCommonHistList);
2164 fFFBckgHisto0RecCuts->AddToOutput(fCommonHistList);
2165 if(fBckgType[1] != kBckgNone) fFFBckgHisto1RecCuts->AddToOutput(fCommonHistList);
2166 if(fBckgType[2] != kBckgNone) fFFBckgHisto2RecCuts->AddToOutput(fCommonHistList);
2167 if(fBckgType[3] != kBckgNone) fFFBckgHisto3RecCuts->AddToOutput(fCommonHistList);
2168 if(fBckgType[4] != kBckgNone) fFFBckgHisto4RecCuts->AddToOutput(fCommonHistList);
2170 if(genJets && genTracks){
2171 fFFBckgHisto0Gen->AddToOutput(fCommonHistList);
2172 if(fBckgType[1] != kBckgNone) fFFBckgHisto1Gen->AddToOutput(fCommonHistList);
2173 if(fBckgType[2] != kBckgNone) fFFBckgHisto2Gen->AddToOutput(fCommonHistList);
2174 if(fBckgType[3] != kBckgNone) fFFBckgHisto3Gen->AddToOutput(fCommonHistList);
2175 if(fBckgType[4] != kBckgNone) fFFBckgHisto4Gen->AddToOutput(fCommonHistList);
2179 fFFBckgHisto0RecEffRec->AddToOutput(fCommonHistList);
2180 fFFBckgHisto0SecRecNS->AddToOutput(fCommonHistList);
2181 fFFBckgHisto0SecRecS->AddToOutput(fCommonHistList);
2182 fFFBckgHisto0SecRecSsc->AddToOutput(fCommonHistList);
2187 fQABckgHisto0RecCuts->AddToOutput(fCommonHistList);
2188 if(fBckgType[1] != kBckgNone) fQABckgHisto1RecCuts->AddToOutput(fCommonHistList);
2189 if(fBckgType[2] != kBckgNone) fQABckgHisto2RecCuts->AddToOutput(fCommonHistList);
2190 if(fBckgType[3] != kBckgNone) fQABckgHisto3RecCuts->AddToOutput(fCommonHistList);
2191 if(fBckgType[4] != kBckgNone) fQABckgHisto4RecCuts->AddToOutput(fCommonHistList);
2192 if(genJets && genTracks){
2193 fQABckgHisto0Gen->AddToOutput(fCommonHistList);
2194 if(fBckgType[1] != kBckgNone) fQABckgHisto1Gen->AddToOutput(fCommonHistList);
2195 if(fBckgType[2] != kBckgNone) fQABckgHisto2Gen->AddToOutput(fCommonHistList);
2196 if(fBckgType[3] != kBckgNone) fQABckgHisto3Gen->AddToOutput(fCommonHistList);
2197 if(fBckgType[4] != kBckgNone) fQABckgHisto4Gen->AddToOutput(fCommonHistList);
2201 if(fh1BckgMult0) fCommonHistList->Add(fh1BckgMult0);
2202 if(fBckgType[1] != kBckgNone) fCommonHistList->Add(fh1BckgMult1);
2203 if(fBckgType[2] != kBckgNone) fCommonHistList->Add(fh1BckgMult2);
2204 if(fBckgType[3] != kBckgNone) fCommonHistList->Add(fh1BckgMult3);
2205 if(fBckgType[4] != kBckgNone) fCommonHistList->Add(fh1BckgMult4);
2209 if(fBranchEmbeddedJets.Length()){
2210 fCommonHistList->Add(fh1FractionPtEmbedded);
2211 fCommonHistList->Add(fh1IndexEmbedded);
2212 fCommonHistList->Add(fh2DeltaPtVsJetPtEmbedded);
2213 fCommonHistList->Add(fh2DeltaPtVsRecJetPtEmbedded);
2214 fCommonHistList->Add(fh1DeltaREmbedded);
2215 fCommonHistList->Add(fh1nEmbeddedJets);
2221 if(fQAMode&1){ // track QA
2222 fQATrackHistosRecCuts->AddToOutput(fCommonHistList);
2223 if(genTracks) fQATrackHistosGen->AddToOutput(fCommonHistList);
2226 if(fQAMode&2){ // jet QA
2227 fQAJetHistosRec->AddToOutput(fCommonHistList);
2228 fQAJetHistosRecCuts->AddToOutput(fCommonHistList);
2229 fQAJetHistosRecCutsLeading->AddToOutput(fCommonHistList);
2230 if(recJetsEff && fEffMode) fQAJetHistosRecEffLeading->AddToOutput(fCommonHistList);
2232 fQAJetHistosGen->AddToOutput(fCommonHistList);
2233 fQAJetHistosGenLeading->AddToOutput(fCommonHistList);
2239 (fBckgType[0]==kBckgClusters || fBckgType[1]==kBckgClusters || fBckgType[2]==kBckgClusters || fBckgType[3]==kBckgClusters || fBckgType[4]==kBckgClusters ||
2240 fBckgType[0]==kBckgClustersOutLeading || fBckgType[1]==kBckgClustersOutLeading || fBckgType[2]==kBckgClustersOutLeading ||
2241 fBckgType[3]==kBckgClustersOutLeading || fBckgType[4]==kBckgClustersOutLeading)) {
2242 fCommonHistList->Add(fh1nRecBckgJetsCuts);
2243 if(genJets) fCommonHistList->Add(fh1nGenBckgJets);
2247 if(fEffMode && recJetsEff && genTracks){
2249 fQATrackHistosRecEffGen->AddToOutput(fCommonHistList);
2250 fQATrackHistosRecEffRec->AddToOutput(fCommonHistList);
2251 fQATrackHistosSecRecNS->AddToOutput(fCommonHistList);
2252 fQATrackHistosSecRecS->AddToOutput(fCommonHistList);
2253 fQATrackHistosSecRecSsc->AddToOutput(fCommonHistList);
2256 fFFHistosRecEffRec->AddToOutput(fCommonHistList);
2257 fFFHistosSecRecNS->AddToOutput(fCommonHistList);
2258 fFFHistosSecRecS->AddToOutput(fCommonHistList);
2259 fFFHistosSecRecSsc->AddToOutput(fCommonHistList);
2261 fCommonHistList->Add(fh1nRecEffJets);
2262 fCommonHistList->Add(fh2PtRecVsGenPrim);
2263 fCommonHistList->Add(fh2PtRecVsGenSec);
2269 fProNtracksLeadingJet = new TProfile("AvgNoOfTracksLeadingJet","AvgNoOfTracksLeadingJet",100,0,250,0,50);
2270 fProDelR80pcPt = new TProfile("AvgdelR80pcPt","AvgdelR80pcPt",100,0,250,0,1);
2272 if(genJets && genTracks){
2273 fProNtracksLeadingJetGen = new TProfile("AvgNoOfTracksLeadingJetGen","AvgNoOfTracksLeadingJetGen",100,0,250,0,50);
2274 fProDelR80pcPtGen = new TProfile("AvgdelR80pcPtGen","AvgdelR80pcPtGen",100,0,250,0,1);
2278 fProNtracksLeadingJetBgrPerp2 = new TProfile("AvgNoOfTracksLeadingJetBgrPerp2","AvgNoOfTracksLeadingJetBgrPerp2",100,0,250,0,50);
2281 fProNtracksLeadingJetRecPrim = new TProfile("AvgNoOfTracksLeadingJetRecPrim","AvgNoOfTracksLeadingJetRecPrim",100,0,250,0,50);
2282 fProDelR80pcPtRecPrim = new TProfile("AvgdelR80pcPtRecPrim","AvgdelR80pcPtRecPrim",100,0,250,0,1);
2283 fProNtracksLeadingJetRecSecNS = new TProfile("AvgNoOfTracksLeadingJetRecSecNS","AvgNoOfTracksLeadingJetRecSecNS",100,0,250,0,50);
2284 fProNtracksLeadingJetRecSecS = new TProfile("AvgNoOfTracksLeadingJetRecSecS","AvgNoOfTracksLeadingJetRecSecS",100,0,250,0,50);
2285 fProNtracksLeadingJetRecSecSsc = new TProfile("AvgNoOfTracksLeadingJetRecSecSsc","AvgNoOfTracksLeadingJetRecSecSsc",100,0,250,0,50);
2289 for(Int_t ii=0; ii<5; ii++){
2290 if(ii==0)strTitJS = "_JetPt20to30";
2291 if(ii==1)strTitJS = "_JetPt30to40";
2292 if(ii==2)strTitJS = "_JetPt40to60";
2293 if(ii==3)strTitJS = "_JetPt60to80";
2294 if(ii==4)strTitJS = "_JetPt80to100";
2296 fProDelRPtSum[ii] = new TProfile(Form("AvgPtSumDelR%s",strTitJS.Data()),Form("AvgPtSumDelR%s",strTitJS.Data()),50,0,1,0,250);
2297 if(genJets && genTracks)
2298 fProDelRPtSumGen[ii] = new TProfile(Form("AvgPtSumDelRGen%s",strTitJS.Data()),Form("AvgPtSumDelRGen%s",strTitJS.Data()),50,0,1,0,250);
2300 fProDelRPtSumBgrPerp2[ii] = new TProfile(Form("AvgPtSumDelRBgrPerp2%s",strTitJS.Data()),Form("AvgPtSumDelRBgrPerp2%s",strTitJS.Data()),50,0,1,0,250);
2302 fProDelRPtSumRecPrim[ii] = new TProfile(Form("AvgPtSumDelRRecPrim%s",strTitJS.Data()),Form("AvgPtSumDelRRecPrim%s",strTitJS.Data()),50,0,1,0,250);
2303 fProDelRPtSumRecSecNS[ii] = new TProfile(Form("AvgPtSumDelRRecSecNS%s",strTitJS.Data()),Form("AvgPtSumDelRRecSecNS%s",strTitJS.Data()),50,0,1,0,250);
2304 fProDelRPtSumRecSecS[ii] = new TProfile(Form("AvgPtSumDelRRecSecS%s",strTitJS.Data()),Form("AvgPtSumDelRRecSecS%s",strTitJS.Data()),50,0,1,0,250);
2305 fProDelRPtSumRecSecSsc[ii] = new TProfile(Form("AvgPtSumDelRRecSecSsc%s",strTitJS.Data()),Form("AvgPtSumDelRRecSecSsc%s",strTitJS.Data()),50,0,1,0,250);
2309 fCommonHistList->Add(fProNtracksLeadingJet);
2310 fCommonHistList->Add(fProDelR80pcPt);
2311 for(int ii=0; ii<5; ii++) fCommonHistList->Add(fProDelRPtSum[ii]);
2313 if(genJets && genTracks){
2314 fCommonHistList->Add(fProNtracksLeadingJetGen);
2315 fCommonHistList->Add(fProDelR80pcPtGen);
2316 for(Int_t ii=0; ii<5; ii++) fCommonHistList->Add(fProDelRPtSumGen[ii]);
2320 fCommonHistList->Add(fProNtracksLeadingJetBgrPerp2);
2321 for(Int_t ii=0; ii<5; ii++) fCommonHistList->Add(fProDelRPtSumBgrPerp2[ii]);
2325 fCommonHistList->Add(fProNtracksLeadingJetRecPrim);
2326 fCommonHistList->Add(fProDelR80pcPtRecPrim);
2327 for(Int_t ii=0; ii<5; ii++) fCommonHistList->Add(fProDelRPtSumRecPrim[ii]);
2329 fCommonHistList->Add(fProNtracksLeadingJetRecSecNS);
2330 for(Int_t ii=0; ii<5; ii++) fCommonHistList->Add(fProDelRPtSumRecSecNS[ii]);
2332 fCommonHistList->Add(fProNtracksLeadingJetRecSecS);
2333 for(Int_t ii=0; ii<5; ii++) fCommonHistList->Add(fProDelRPtSumRecSecS[ii]);
2335 fCommonHistList->Add(fProNtracksLeadingJetRecSecSsc);
2336 for(Int_t ii=0; ii<5; ii++) fCommonHistList->Add(fProDelRPtSumRecSecSsc[ii]);
2340 // Default analysis utils
2341 fAnaUtils = new AliAnalysisUtils();
2343 // Not used yet, but to be save, forward vertex z cut to analysis utils object
2344 fAnaUtils->SetMaxVtxZ(fMaxVertexZ);
2346 // Load PID framework if desired
2347 if(fDebug > 1) Printf("AliAnalysisTaskIDFragmentationFunction::UserCreateOutputObjects() -> Loading PID framework");
2349 fUseJetPIDtask = fIDFFMode || fFFMode;
2350 fUseInclusivePIDtask = fQAMode && (fQAMode&1);
2352 if (fUseJetPIDtask || fUseInclusivePIDtask) {
2353 TObjArray* tasks = AliAnalysisManager::GetAnalysisManager()->GetTasks();
2355 Printf("ERROR loading PID tasks: Failed to retrieve tasks from analysis manager!\n");
2357 fUseJetPIDtask = kFALSE;
2358 fUseInclusivePIDtask = kFALSE;
2361 if (fUseJetPIDtask) {
2362 delete [] fJetPIDtask;
2365 if (fNumJetPIDtasks > 0) {
2366 fJetPIDtask = new AliAnalysisTaskPID*[fNumJetPIDtasks];
2368 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
2369 fJetPIDtask[i] = (AliAnalysisTaskPID*)tasks->FindObject(fNameJetPIDtask[i].Data());
2371 if (!fJetPIDtask[i]) {
2372 Printf("ERROR: Failed to load jet pid task!\n");
2373 fUseJetPIDtask = kFALSE;
2378 Printf("WARNING: zero jet pid tasks!\n");
2379 fUseJetPIDtask = kFALSE;
2383 if (fUseInclusivePIDtask) {
2384 delete [] fInclusivePIDtask;
2385 fInclusivePIDtask = 0x0;
2387 if (fNumInclusivePIDtasks > 0) {
2388 fInclusivePIDtask = new AliAnalysisTaskPID*[fNumInclusivePIDtasks];
2390 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
2391 fInclusivePIDtask[i] = (AliAnalysisTaskPID*)tasks->FindObject(fNameInclusivePIDtask[i].Data());
2393 if (!fInclusivePIDtask[i]) {
2394 Printf("ERROR: Failed to load inclusive pid task!\n");
2395 fUseInclusivePIDtask = kFALSE;
2400 Printf("WARNING: zero inclusive pid tasks!\n");
2401 fUseInclusivePIDtask = kFALSE;
2406 const Int_t nRefMultBins = 100;
2407 const Double_t refMultUp = 100.;
2408 const Double_t refMultDown = 0.;
2410 const Int_t nJetPtBins = 20;
2411 const Double_t jetPtUp = 100.;
2412 const Double_t jetPtDown = 0.;
2414 const Int_t nCentBins = 12;
2415 const Double_t binsCentpp[nCentBins+1] = { 0, 0.01, 0.1, 1, 5, 10, 15, 20, 30, 40, 50, 70, 100};
2417 fhJetPtRefMultEta5 = new TH2F("fhJetPtRefMultEta5",
2418 "Correlation between jet energy and event multiplicity (|#eta| < 0.5);Ref. mult. (|#eta| < 0.5);#it{p}_{T, jet}^{ch} (GeV/#it{c})",
2419 nRefMultBins, refMultDown, refMultUp, nJetPtBins, jetPtDown, jetPtUp);
2420 fhJetPtRefMultEta8 = new TH2F("fhJetPtRefMultEta8",
2421 "Correlation between jet energy and event multiplicity (|#eta| < 0.8);Ref. mult. (|#eta| < 0.8);#it{p}_{T, jet}^{ch} (GeV/#it{c})",
2422 nRefMultBins, refMultDown, refMultUp, nJetPtBins, jetPtDown, jetPtUp);
2423 fhJetPtMultPercent = new TH2F("fhJetPtMultPercent",
2424 "Correlation between jet energy and event multiplicity percentile (V0M);Multiplicity Percentile (V0M);#it{p}_{T, jet}^{ch} (GeV/#it{c})",
2425 nCentBins, binsCentpp, nJetPtBins, jetPtDown, jetPtUp);
2426 fCommonHistList->Add(fhJetPtRefMultEta5);
2427 fCommonHistList->Add(fhJetPtRefMultEta8);
2428 fCommonHistList->Add(fhJetPtMultPercent);
2430 if (fUseJetPIDtask) {
2431 const Int_t nPtBins = 68;
2432 Double_t binsPt[nPtBins+1] = {0. , 0.05, 0.1, 0.15, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45,
2433 0.5, 0.55, 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95,
2434 1.0, 1.1 , 1.2, 1.3 , 1.4, 1.5 , 1.6, 1.7 , 1.8, 1.9 ,
2435 2.0, 2.2 , 2.4, 2.6 , 2.8, 3.0 , 3.2, 3.4 , 3.6, 3.8 ,
2436 4.0, 4.5 , 5.0, 5.5 , 6.0, 6.5 , 7.0, 8.0 , 9.0, 10.0,
2437 11.0, 12.0, 13.0, 14.0, 15.0, 16.0, 18.0, 20.0, 22.0, 24.0,
2438 26.0, 28.0, 30.0, 32.0, 34.0, 36.0, 40.0, 45.0, 50.0 };
2440 const Int_t DCAbins = 320;
2441 const Double_t DCA_Z_max = 3.5;
2442 const Double_t DCA_XY_max = 2.5;
2444 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);
2445 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);
2446 fCommonHistList->Add(fhDCA_XY);
2447 fCommonHistList->Add(fhDCA_Z);
2449 for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
2450 fhDCA_XY_prim_MCID[i] = new TH2F(Form("fhDCA_XY_prim_MCID_%s", AliPID::ParticleShortName(i)),
2451 Form("Rec. %s (prim.);#it{p}_{T} (GeV/#it{c});DCA_{XY}", AliPID::ParticleLatexName(i)),
2452 nPtBins, binsPt, DCAbins, -DCA_XY_max, DCA_XY_max);
2453 fhDCA_Z_prim_MCID[i] = new TH2F(Form("fhDCA_Z_prim_MCID_%s", AliPID::ParticleShortName(i)),
2454 Form("Rec. %s (prim.);#it{p}_{T} (GeV/#it{c});DCA_{Z}", AliPID::ParticleLatexName(i)),
2455 nPtBins, binsPt, DCAbins, -DCA_Z_max, DCA_Z_max);
2456 fCommonHistList->Add(fhDCA_XY_prim_MCID[i]);
2457 fCommonHistList->Add(fhDCA_Z_prim_MCID[i]);
2459 fhDCA_XY_sec_MCID[i] = new TH2F(Form("fhDCA_XY_sec_MCID_%s", AliPID::ParticleShortName(i)),
2460 Form("Rec. %s (sec.);#it{p}_{T} (GeV/#it{c});DCA_{XY}", AliPID::ParticleLatexName(i)),
2461 nPtBins, binsPt, DCAbins, -DCA_XY_max, DCA_XY_max);
2462 fhDCA_Z_sec_MCID[i] = new TH2F(Form("fhDCA_Z_sec_MCID_%s", AliPID::ParticleShortName(i)),
2463 Form("Rec. %s (sec.);#it{p}_{T} (GeV/#it{c});DCA_{Z}", AliPID::ParticleLatexName(i)),
2464 nPtBins, binsPt, DCAbins, -DCA_Z_max, DCA_Z_max);
2465 fCommonHistList->Add(fhDCA_XY_sec_MCID[i]);
2466 fCommonHistList->Add(fhDCA_Z_sec_MCID[i]);
2471 // =========== Switch on Sumw2 for all histos ===========
2472 for (Int_t i=0; i<fCommonHistList->GetEntries(); ++i){
2473 TH1 *h1 = dynamic_cast<TH1*>(fCommonHistList->At(i));
2474 if (h1) h1->Sumw2();
2476 THnSparse *hnSparse = dynamic_cast<THnSparse*>(fCommonHistList->At(i));
2477 if(hnSparse) hnSparse->Sumw2();
2481 TH1::AddDirectory(oldStatus);
2483 if(fDebug > 2) Printf("AliAnalysisTaskIDFragmentationFunction::UserCreateOutputObjects() -> Posting Output");
2485 PostData(1, fCommonHistList);
2487 if(fDebug > 2) Printf("AliAnalysisTaskIDFragmentationFunction::UserCreateOutputObjects() -> Done");
2490 //_______________________________________________
2491 void AliAnalysisTaskIDFragmentationFunction::Init()
2494 if(fDebug > 1) Printf("AliAnalysisTaskIDFragmentationFunction::Init()");
2498 //_____________________________________________________________
2499 void AliAnalysisTaskIDFragmentationFunction::UserExec(Option_t *)
2502 // Called for each event
2504 if(fDebug > 1) Printf("AliAnalysisTaskIDFragmentationFunction::UserExec()");
2507 if(fDebug > 1) Printf("Analysis event #%5d", (Int_t) fEntry);
2509 // Trigger selection
2510 AliInputEventHandler* inputHandler = (AliInputEventHandler*)
2511 ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
2513 if(!(inputHandler->IsEventSelected() & fEvtSelectionMask)){
2514 fh1EvtSelection->Fill(1.);
2515 if (fDebug > 1 ) Printf(" Trigger Selection: event REJECTED ... ");
2516 PostData(1, fCommonHistList);
2520 fESD = dynamic_cast<AliESDEvent*>(InputEvent());
2522 if(fDebug>3) Printf("%s:%d ESDEvent not found in the input", (char*)__FILE__,__LINE__);
2525 fMCEvent = MCEvent();
2527 if(fDebug>3) Printf("%s:%d MCEvent not found in the input", (char*)__FILE__,__LINE__);
2530 // get AOD event from input/ouput
2531 TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
2532 if( handler && handler->InheritsFrom("AliAODInputHandler") ) {
2533 fAOD = ((AliAODInputHandler*)handler)->GetEvent();
2534 if(fUseAODInputJets) fAODJets = fAOD;
2535 if (fDebug > 1) Printf("%s:%d AOD event from input", (char*)__FILE__,__LINE__);
2538 handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
2539 if( handler && handler->InheritsFrom("AliAODHandler") ) {
2540 fAOD = ((AliAODHandler*)handler)->GetAOD();
2542 if (fDebug > 1) Printf("%s:%d AOD event from output", (char*)__FILE__,__LINE__);
2546 if(!fAODJets && !fUseAODInputJets){ // case we have AOD in input & output and want jets from output
2547 TObject* outHandler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
2548 if( outHandler && outHandler->InheritsFrom("AliAODHandler") ) {
2549 fAODJets = ((AliAODHandler*)outHandler)->GetAOD();
2550 if (fDebug > 1) Printf("%s:%d jets from output AOD", (char*)__FILE__,__LINE__);
2554 if(fNonStdFile.Length()!=0){
2555 // case we have an AOD extension - fetch the jets from the extended output
2557 AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
2558 fAODExtension = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);
2560 if(fDebug>1)Printf("AODExtension not found for %s",fNonStdFile.Data());
2565 Printf("%s:%d AODEvent not found", (char*)__FILE__,__LINE__);
2569 Printf("%s:%d AODEvent with jet branch not found", (char*)__FILE__,__LINE__);
2574 // event selection **************************************************
2575 // *** event class ***
2576 AliVEvent* evtForCentDetermination = handler->InheritsFrom("AliAODInputHandler") ? fAOD : InputEvent();
2578 Double_t centPercent = -1;
2581 if(handler->InheritsFrom("AliAODInputHandler")){
2582 // since it is not supported by the helper task define own classes
2583 centPercent = fAOD->GetHeader()->GetCentrality();
2585 if(centPercent>10) cl = 2;
2586 if(centPercent>30) cl = 3;
2587 if(centPercent>50) cl = 4;
2590 cl = AliAnalysisHelperJetTasks::EventClass();
2591 if(fESD) centPercent = fESD->GetCentrality()->GetCentralityPercentile("V0M"); // retrieve value 'by hand'
2594 if(cl!=fEventClass){
2595 // event not in selected event class, reject event
2596 if (fDebug > 1) Printf("%s:%d event not in selected event class: event REJECTED ...",(char*)__FILE__,__LINE__);
2597 fh1EvtSelection->Fill(2.);
2598 PostData(1, fCommonHistList);
2603 // Retrieve reference multiplicities in |eta|<0.8 and <0.5
2604 const Int_t refMult5 = fAOD->GetHeader()->GetRefMultiplicityComb05();
2605 const Int_t refMult8 = fAOD->GetHeader()->GetRefMultiplicityComb08();
2606 const Double_t centPercentPP = fAnaUtils->GetMultiplicityPercentile(fAOD, "V0M");
2609 // Count events with trigger selection, note: Set centrality percentile fix to -1 for pp for PID framework
2610 if (fUseJetPIDtask) {
2611 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
2612 fJetPIDtask[i]->IncrementEventCounter(fIsPP ? -1. : centPercent, AliAnalysisTaskPID::kTriggerSel);
2616 if (fUseInclusivePIDtask) {
2617 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
2618 fInclusivePIDtask[i]->IncrementEventCounter(fIsPP ? -1. : centPercent, AliAnalysisTaskPID::kTriggerSel);
2622 // *** vertex cut ***
2623 AliAODVertex* primVtx = fAOD->GetPrimaryVertex();
2625 Printf("%s:%d Primary vertex not found", (char*)__FILE__,__LINE__);
2629 Int_t nTracksPrim = primVtx->GetNContributors();
2630 fh1VertexNContributors->Fill(nTracksPrim);
2633 if (fDebug > 1) Printf("%s:%d primary vertex selection: %d", (char*)__FILE__,__LINE__,nTracksPrim);
2634 if(nTracksPrim <= 0) {
2635 if (fDebug > 1) Printf("%s:%d primary vertex selection: event REJECTED...",(char*)__FILE__,__LINE__);
2636 fh1EvtSelection->Fill(3.);
2637 PostData(1, fCommonHistList);
2641 TString primVtxName(primVtx->GetName());
2643 if(primVtxName.CompareTo("TPCVertex",TString::kIgnoreCase) == 1){
2644 if (fDebug > 1) Printf("%s:%d primary vertex selection: TPC vertex, event REJECTED...",(char*)__FILE__,__LINE__);
2645 fh1EvtSelection->Fill(5.);
2646 PostData(1, fCommonHistList);
2650 // Count events with trigger selection and vtx cut, note: Set centrality percentile fix to -1 for pp for PID framework
2651 if (fUseJetPIDtask) {
2652 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
2653 fJetPIDtask[i]->IncrementEventCounter(fIsPP ? -1. : centPercent, AliAnalysisTaskPID::kTriggerSelAndVtxCut);
2657 if (fUseInclusivePIDtask) {
2658 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
2659 fInclusivePIDtask[i]->IncrementEventCounter(fIsPP ? -1. : centPercent, AliAnalysisTaskPID::kTriggerSelAndVtxCut);
2663 fh1VertexZ->Fill(primVtx->GetZ());
2665 if(TMath::Abs(primVtx->GetZ())>fMaxVertexZ){
2666 if (fDebug > 1) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ());
2667 fh1EvtSelection->Fill(4.);
2668 PostData(1, fCommonHistList);
2672 // Count events with trigger selection, vtx cut and z vtx cut, note: Set centrality percentile fix to -1 for pp for PID framework
2673 if (fUseJetPIDtask) {
2674 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
2675 fJetPIDtask[i]->IncrementEventCounter(fIsPP ? -1. : centPercent, AliAnalysisTaskPID::kTriggerSelAndVtxCutAndZvtxCutNoPileUpRejection);
2679 if (fUseInclusivePIDtask) {
2680 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
2681 fInclusivePIDtask[i]->IncrementEventCounter(fIsPP ? -1. : centPercent, AliAnalysisTaskPID::kTriggerSelAndVtxCutAndZvtxCutNoPileUpRejection);
2686 // Store for each task, whether this task would tag this event as pile-up or not
2687 const Int_t arrSizeJet = TMath::Max(1, fNumJetPIDtasks);
2688 const Int_t arrSizeInclusive = TMath::Max(1, fNumInclusivePIDtasks);
2689 Bool_t isPileUpJetPIDtask[arrSizeJet];
2690 Bool_t isPileUpInclusivePIDtask[arrSizeInclusive];
2692 for (Int_t i = 0; i < arrSizeJet; i++)
2693 isPileUpJetPIDtask[i] = kFALSE;
2695 for (Int_t i = 0; i < arrSizeInclusive; i++)
2696 isPileUpInclusivePIDtask[i] = kFALSE;
2698 // Check whether there is at least one task that does not reject the event (saves processing time in the following code)
2699 Bool_t isPileUpForAllJetPIDTasks = kTRUE;
2700 Bool_t isPileUpForAllInclusivePIDTasks = kTRUE;
2702 // Count events with trigger selection, vtx cut, z vtx cut and after pile-up rejection (if enabled in that task)
2703 // Note: Set centrality percentile fix to -1 for pp for PID framework
2704 if (fUseJetPIDtask) {
2705 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
2706 isPileUpJetPIDtask[i] = fJetPIDtask[i]->GetIsPileUp(fAOD, fJetPIDtask[i]->GetPileUpRejectionType());
2707 if (!isPileUpJetPIDtask[i])
2708 fJetPIDtask[i]->IncrementEventCounter(fIsPP ? -1. : centPercent, AliAnalysisTaskPID::kTriggerSelAndVtxCutAndZvtxCut);
2710 isPileUpForAllJetPIDTasks = isPileUpForAllJetPIDTasks && isPileUpJetPIDtask[i];
2714 if (fUseInclusivePIDtask) {
2715 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
2716 isPileUpInclusivePIDtask[i] = fInclusivePIDtask[i]->GetIsPileUp(fAOD, fInclusivePIDtask[i]->GetPileUpRejectionType());
2717 if (!isPileUpInclusivePIDtask[i])
2718 fInclusivePIDtask[i]->IncrementEventCounter(fIsPP ? -1. : centPercent, AliAnalysisTaskPID::kTriggerSelAndVtxCutAndZvtxCut);
2720 isPileUpForAllInclusivePIDTasks = isPileUpForAllInclusivePIDTasks && isPileUpInclusivePIDtask[i];
2725 if (fDebug > 1) Printf("%s:%d event ACCEPTED ...",(char*)__FILE__,__LINE__);
2726 fh1EvtSelection->Fill(0.);
2727 fh1VtxSelection->Fill(primVtx->GetType());
2728 fh1EvtCent->Fill(centPercent);
2730 // Set centrality percentile fix to -1 for pp to be used for the PID framework
2736 // Call ConfigureTaskForCurrentEvent of PID tasks to ensure that everything is set up properly for the current event
2737 // (e.g. run/period dependence of max eta variation map)
2738 if (fUseInclusivePIDtask) {
2739 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++)
2740 if (!isPileUpInclusivePIDtask[i])
2741 fInclusivePIDtask[i]->ConfigureTaskForCurrentEvent(fAOD);
2744 if (fUseJetPIDtask) {
2745 for (Int_t i = 0; i < fNumJetPIDtasks; i++)
2746 if (!isPileUpJetPIDtask[i])
2747 fJetPIDtask[i]->ConfigureTaskForCurrentEvent(fAOD);
2752 //___ get MC information __________________________________________________________________
2754 fh1Trials->Fill("#sum{ntrials}",fAvgTrials);
2756 if (fUseInclusivePIDtask) {
2757 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
2758 if (!isPileUpInclusivePIDtask[i])
2759 fInclusivePIDtask[i]->FillPythiaTrials(fAvgTrials);
2763 if (fUseJetPIDtask) {
2764 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
2765 if (!isPileUpJetPIDtask[i])
2766 fJetPIDtask[i]->FillPythiaTrials(fAvgTrials);
2770 Double_t ptHard = 0.;
2771 Double_t nTrials = 1; // trials for MC trigger weight for real data
2774 AliGenEventHeader* genHeader = fMCEvent->GenEventHeader();
2778 AliGenPythiaEventHeader* pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
2779 AliGenHijingEventHeader* hijingGenHeader = 0x0;
2781 if(pythiaGenHeader){
2782 if(fDebug>3) Printf("%s:%d pythiaGenHeader found", (char*)__FILE__,__LINE__);
2783 nTrials = pythiaGenHeader->Trials();
2784 ptHard = pythiaGenHeader->GetPtHard();
2786 fh1PtHard->Fill(ptHard);
2787 fh1PtHardTrials->Fill(ptHard,nTrials);
2789 } else { // no pythia, hijing?
2791 if(fDebug>3) Printf("%s:%d no pythiaGenHeader found", (char*)__FILE__,__LINE__);
2793 hijingGenHeader = dynamic_cast<AliGenHijingEventHeader*>(genHeader);
2794 if(!hijingGenHeader){
2795 Printf("%s:%d no pythiaGenHeader or hjingGenHeader found", (char*)__FILE__,__LINE__);
2797 if(fDebug>3) Printf("%s:%d hijingGenHeader found", (char*)__FILE__,__LINE__);
2801 //fh1Trials->Fill("#sum{ntrials}",fAvgTrials);
2805 //___ fetch jets __________________________________________________________________________
2807 Int_t nJ = GetListOfJets(fJetsRec, kJetsRec);
2809 if(nJ>=0) nRecJets = fJetsRec->GetEntries();
2810 if(fDebug>2)Printf("%s:%d Selected Rec jets: %d %d",(char*)__FILE__,__LINE__,nJ,nRecJets);
2811 if(nJ != nRecJets) Printf("%s:%d Mismatch Selected Rec Jets: %d %d",(char*)__FILE__,__LINE__,nJ,nRecJets);
2813 Int_t nJCuts = GetListOfJets(fJetsRecCuts, kJetsRecAcceptance);
2814 Int_t nRecJetsCuts = 0;
2815 if(nJCuts>=0) nRecJetsCuts = fJetsRecCuts->GetEntries();
2816 if(fDebug>2)Printf("%s:%d Selected Rec jets after cuts: %d %d",(char*)__FILE__,__LINE__,nJCuts,nRecJetsCuts);
2817 if(nRecJetsCuts != nJCuts) Printf("%s:%d Mismatch selected Rec jets after cuts: %d %d",(char*)__FILE__,__LINE__,nJCuts,nRecJetsCuts);
2818 fh1nRecJetsCuts->Fill(nRecJetsCuts);
2820 if(fJetTypeGen==kJetsKine || fJetTypeGen == kJetsKineAcceptance) fJetsGen->SetOwner(kTRUE); // kine aod jets allocated on heap, delete them with TList::Clear()
2822 Int_t nJGen = GetListOfJets(fJetsGen, fJetTypeGen);
2824 if(nJGen>=0) nGenJets = fJetsGen->GetEntries();
2825 if(fDebug>2)Printf("%s:%d Selected Gen jets: %d %d",(char*)__FILE__,__LINE__,nJGen,nGenJets);
2827 if(nJGen != nGenJets) Printf("%s:%d Mismatch selected Gen jets: %d %d",(char*)__FILE__,__LINE__,nJGen,nGenJets);
2828 fh1nGenJets->Fill(nGenJets);
2831 if(fJetTypeRecEff==kJetsKine || fJetTypeRecEff == kJetsKineAcceptance) fJetsRecEff->SetOwner(kTRUE); // kine aod jets allocated on heap, delete them with TList::Clear()
2832 Int_t nJRecEff = GetListOfJets(fJetsRecEff, fJetTypeRecEff);
2833 Int_t nRecEffJets = 0;
2834 if(nJRecEff>=0) nRecEffJets = fJetsRecEff->GetEntries();
2835 if(fDebug>2)Printf("%s:%d Selected RecEff jets: %d %d",(char*)__FILE__,__LINE__,nJRecEff,nRecEffJets);
2836 if(nJRecEff != nRecEffJets) Printf("%s:%d Mismatch selected RecEff jets: %d %d",(char*)__FILE__,__LINE__,nJRecEff,nRecEffJets);
2837 fh1nRecEffJets->Fill(nRecEffJets);
2840 Int_t nEmbeddedJets = 0;
2841 TArrayI iEmbeddedMatchIndex;
2842 TArrayF fEmbeddedPtFraction;
2845 if(fBranchEmbeddedJets.Length()){
2846 Int_t nJEmbedded = GetListOfJets(fJetsEmbedded, kJetsEmbedded);
2847 if(nJEmbedded>=0) nEmbeddedJets = fJetsEmbedded->GetEntries();
2848 if(fDebug>2)Printf("%s:%d Selected Embedded jets: %d %d",(char*)__FILE__,__LINE__,nJEmbedded,nEmbeddedJets);
2849 if(nJEmbedded != nEmbeddedJets) Printf("%s:%d Mismatch Selected Embedded Jets: %d %d",(char*)__FILE__,__LINE__,nJEmbedded,nEmbeddedJets);
2850 fh1nEmbeddedJets->Fill(nEmbeddedJets);
2852 Float_t maxDist = 0.3;
2854 iEmbeddedMatchIndex.Set(nEmbeddedJets);
2855 fEmbeddedPtFraction.Set(nEmbeddedJets);
2857 iEmbeddedMatchIndex.Reset(-1);
2858 fEmbeddedPtFraction.Reset(0);
2860 AliAnalysisHelperJetTasks::GetJetMatching(fJetsEmbedded, nEmbeddedJets,
2861 fJetsRecCuts, nRecJetsCuts,
2862 iEmbeddedMatchIndex, fEmbeddedPtFraction,
2867 //____ fetch background clusters ___________________________________________________
2869 (fBckgType[0]==kBckgClusters || fBckgType[1]==kBckgClusters || fBckgType[2]==kBckgClusters || fBckgType[3]==kBckgClusters || fBckgType[4]==kBckgClusters ||
2870 fBckgType[0]==kBckgClustersOutLeading || fBckgType[1]==kBckgClustersOutLeading || fBckgType[2]==kBckgClustersOutLeading ||
2871 fBckgType[3]==kBckgClustersOutLeading || fBckgType[4]==kBckgClustersOutLeading)){
2873 Int_t nBJ = GetListOfBckgJets(fBckgJetsRec, kJetsRec);
2874 Int_t nRecBckgJets = 0;
2875 if(nBJ>=0) nRecBckgJets = fBckgJetsRec->GetEntries();
2876 if(fDebug>2)Printf("%s:%d Selected Rec background jets: %d %d",(char*)__FILE__,__LINE__,nBJ,nRecBckgJets);
2877 if(nBJ != nRecBckgJets) Printf("%s:%d Mismatch Selected Rec background jets: %d %d",(char*)__FILE__,__LINE__,nBJ,nRecBckgJets);
2879 Int_t nBJCuts = GetListOfBckgJets(fBckgJetsRecCuts, kJetsRecAcceptance);
2880 Int_t nRecBckgJetsCuts = 0;
2881 if(nBJCuts>=0) nRecBckgJetsCuts = fBckgJetsRecCuts->GetEntries();
2882 if(fDebug>2)Printf("%s:%d Selected Rec background jets after cuts: %d %d",(char*)__FILE__,__LINE__,nJCuts,nRecJetsCuts);
2883 if(nRecBckgJetsCuts != nBJCuts) Printf("%s:%d Mismatch selected Rec background jets after cuts: %d %d",(char*)__FILE__,__LINE__,nBJCuts,nRecBckgJetsCuts);
2884 fh1nRecBckgJetsCuts->Fill(nRecBckgJetsCuts);
2886 if(0){ // protection OB - not yet implemented
2887 if(fJetTypeGen==kJetsKine || fJetTypeGen == kJetsKineAcceptance) fBckgJetsGen->SetOwner(kTRUE); // kine aod jets allocated on heap, delete them with TList::Clear()
2888 Int_t nBJGen = GetListOfBckgJets(fBckgJetsGen, fJetTypeGen);
2889 Int_t nGenBckgJets = 0;
2890 if(nBJGen>=0) nGenBckgJets = fBckgJetsGen->GetEntries();
2891 if(fDebug>2)Printf("%s:%d Selected Gen background jets: %d %d",(char*)__FILE__,__LINE__,nBJGen,nGenBckgJets);
2892 if(nBJGen != nGenBckgJets) Printf("%s:%d Mismatch selected Gen background jets: %d %d",(char*)__FILE__,__LINE__,nBJGen,nGenBckgJets);
2893 fh1nGenBckgJets->Fill(nGenBckgJets);
2898 //____ fetch particles __________________________________________________________
2901 if(fUseExtraTracks == 1) nTCuts = GetListOfTracks(fTracksRecCuts, kTrackAODExtraCuts);
2902 else if(fUseExtraTracks == -1) nTCuts = GetListOfTracks(fTracksRecCuts, kTrackAODExtraonlyCuts);
2903 else nTCuts = GetListOfTracks(fTracksRecCuts, kTrackAODCuts);
2905 Int_t nRecPartCuts = 0;
2906 if(nTCuts>=0) nRecPartCuts = fTracksRecCuts->GetEntries();
2907 if(fDebug>2)Printf("%s:%d Selected Rec tracks after cuts: %d %d",(char*)__FILE__,__LINE__,nTCuts,nRecPartCuts);
2908 if(nRecPartCuts != nTCuts) Printf("%s:%d Mismatch selected Rec tracks after cuts: %d %d",
2909 (char*)__FILE__,__LINE__,nTCuts,nRecPartCuts);
2910 fh1EvtMult->Fill(nRecPartCuts);
2913 Int_t nTGen = GetListOfTracks(fTracksGen,fTrackTypeGen);
2915 if(nTGen>=0) nGenPart = fTracksGen->GetEntries();
2916 if(fDebug>2)Printf("%s:%d Selected Gen tracks: %d %d",(char*)__FILE__,__LINE__,nTGen,nGenPart);
2917 if(nGenPart != nTGen) Printf("%s:%d Mismatch selected Gen tracks: %d %d",(char*)__FILE__,__LINE__,nTGen,nGenPart);
2919 // Just cut on filterBit, but take the full acceptance
2920 Int_t nTCutsEfficiency = GetListOfTracks(fTracksRecCutsEfficiency, kTrackAODQualityCuts);
2921 Int_t nRecPartCutsEfficiency = 0;
2922 if(nTCutsEfficiency>=0) nRecPartCutsEfficiency = fTracksRecCutsEfficiency->GetEntries();
2923 if(fDebug>2)Printf("%s:%d Selected Rec tracks efficiency after cuts: %d %d",(char*)__FILE__,__LINE__,
2924 nTCutsEfficiency,nRecPartCutsEfficiency);
2925 if(nRecPartCutsEfficiency != nTCutsEfficiency) Printf("%s:%d Mismatch selected Rec tracks after cuts: %d %d",
2926 (char*)__FILE__,__LINE__,nTCutsEfficiency,nRecPartCutsEfficiency);
2928 AliPIDResponse* pidResponse = 0x0;
2929 Bool_t tuneOnDataTPC = kFALSE;
2930 if (fUseJetPIDtask || fUseInclusivePIDtask) {
2931 if (!inputHandler) {
2932 AliFatal("Input handler needed");
2936 // PID response object
2937 pidResponse = inputHandler->GetPIDResponse();
2939 AliFatal("PIDResponse object was not created");
2943 tuneOnDataTPC = pidResponse->IsTunedOnData() &&
2944 ((pidResponse->GetTunedOnDataMask() & AliPIDResponse::kDetTPC) == AliPIDResponse::kDetTPC);
2949 if(fDebug>2)Printf("%s:%d Starting processing...",(char*)__FILE__,__LINE__);
2951 //____ analysis, fill histos ___________________________________________________
2956 TClonesArray *tca = fUseInclusivePIDtask ? dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName())) : 0x0;
2958 // Fill efficiency for generated primaries and also fill histos for generated yields (primaries + all)
2959 // Efficiency, inclusive - particle level
2960 if (fUseInclusivePIDtask && tca && !isPileUpForAllInclusivePIDTasks) {
2961 for (Int_t it = 0; it < tca->GetEntriesFast(); it++) {
2962 AliAODMCParticle* part = dynamic_cast<AliAODMCParticle*>(tca->At(it));
2966 // Define clean MC sample with corresponding particle level track cuts:
2967 // - MC-track must be in desired eta range
2968 // - MC-track must be physical primary
2969 // - Species must be one of those in question
2971 if (part->Eta() > fTrackEtaMax || part->Eta() < fTrackEtaMin)
2974 Int_t mcID = AliAnalysisTaskPID::PDGtoMCID(part->GetPdgCode());
2976 // Following lines are not needed - just keep other species (like casecades) - will end up in overflow bin
2977 // and only affect the efficiencies for all (i.e. not identified) what is desired!
2978 //if (mcID == AliPID::kUnknown)
2981 if (!part->IsPhysicalPrimary())
2984 Int_t iMother = part->GetMother();
2986 continue; // Not a physical primary
2989 Double_t pT = part->Pt();
2991 // AliAODMCParticle->Charge() calls TParticlePDG->Charge(), which returns the charge in units of e0 / 3
2992 Double_t chargeMC = part->Charge() / 3.;
2994 if (TMath::Abs(chargeMC) < 0.01)
2995 continue; // Reject neutral particles (only relevant, if mcID is not used)
2997 Double_t valuesGenYield[AliAnalysisTaskPID::kGenYieldNumAxes] = { static_cast<Double_t>(mcID), pT, centPercent, -1, -1, -1, -1 };
2999 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
3000 if (!isPileUpInclusivePIDtask[i] && fInclusivePIDtask[i]->IsInAcceptedEtaRange(TMath::Abs(part->Eta()))) {
3001 valuesGenYield[fInclusivePIDtask[i]->GetIndexOfChargeAxisGenYield()] = chargeMC;
3002 fInclusivePIDtask[i]->FillGeneratedYield(valuesGenYield);
3006 Double_t valuesEff[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), pT, part->Eta(), chargeMC,
3007 centPercent, -1, -1, -1 };// no jet pT etc since inclusive spectrum
3008 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
3009 if (!isPileUpInclusivePIDtask[i])
3010 fInclusivePIDtask[i]->FillEfficiencyContainer(valuesEff, AliAnalysisTaskPID::kStepGenWithGenCuts);
3015 if (fUseInclusivePIDtask && !isPileUpForAllInclusivePIDTasks) {
3016 //Efficiency, inclusive - detector level
3017 for(Int_t it=0; it<nRecPartCutsEfficiency; ++it){
3018 // fill inclusive tracks XXX, they have the same track cuts!
3019 AliAODTrack * inclusiveaod = dynamic_cast<AliAODTrack*>(fTracksRecCutsEfficiency->At(it));
3021 Double_t dEdxTPC = tuneOnDataTPC ? pidResponse->GetTPCsignalTunedOnData(inclusiveaod)
3022 : inclusiveaod->GetTPCsignal();
3027 Bool_t survivedTPCCutMIGeo = AliAnalysisTaskPID::TPCCutMIGeo(inclusiveaod, InputEvent());
3028 Bool_t survivedTPCnclCut = AliAnalysisTaskPID::TPCnclCut(inclusiveaod);
3030 Int_t label = TMath::Abs(inclusiveaod->GetLabel());
3032 // find MC track in our list, if available
3033 AliAODMCParticle* gentrack = tca ? dynamic_cast<AliAODMCParticle*>(tca->At(label)) : 0x0;
3037 pdg = gentrack->GetPdgCode();
3039 // For efficiency: Reconstructed track has survived all cuts on the detector level (no cut on eta acceptance)
3040 // and has an associated MC track
3041 // -> Check whether associated MC track belongs to the clean MC sample defined above,
3042 // i.e. survives the particle level track cuts
3044 Int_t mcID = AliAnalysisTaskPID::PDGtoMCID(pdg);
3046 // Following lines are not needed - just keep other species (like casecades) - will end up in overflow bin
3047 // and only affect the efficiencies for all (i.e. not identified) what is desired!
3048 //if (mcID == AliPID::kUnknown)
3051 // Fill efficiency for reconstructed primaries
3052 if (!gentrack->IsPhysicalPrimary())
3055 * Int_t iMother = gentrack->GetMother();
3057 * continue; // Not a physical primary
3060 if (gentrack->Eta() > fTrackEtaMax || gentrack->Eta() < fTrackEtaMin)
3063 // AliAODMCParticle->Charge() calls TParticlePDG->Charge(), which returns the charge in units of e0 / 3
3064 Double_t value[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), gentrack->Pt(), gentrack->Eta(), gentrack->Charge() / 3.,
3066 -1, -1, -1 };// no jet pT etc since inclusive spectrum
3067 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
3068 if (!isPileUpInclusivePIDtask[i] &&
3069 ((!fInclusivePIDtask[i]->GetUseTPCCutMIGeo() && !fInclusivePIDtask[i]->GetUseTPCnclCut()) ||
3070 (survivedTPCCutMIGeo && fInclusivePIDtask[i]->GetUseTPCCutMIGeo()) ||
3071 (survivedTPCnclCut && fInclusivePIDtask[i]->GetUseTPCnclCut())))
3072 fInclusivePIDtask[i]->FillEfficiencyContainer(value, AliAnalysisTaskPID::kStepRecWithGenCuts);
3075 Double_t valueMeas[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), inclusiveaod->Pt(), inclusiveaod->Eta(),
3076 static_cast<Double_t>(inclusiveaod->Charge()), centPercent,
3077 -1, -1, -1 };// no jet pT etc since inclusive spectrum
3078 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
3079 if (!isPileUpInclusivePIDtask[i] &&
3080 ((!fInclusivePIDtask[i]->GetUseTPCCutMIGeo() && !fInclusivePIDtask[i]->GetUseTPCnclCut()) ||
3081 (survivedTPCCutMIGeo && fInclusivePIDtask[i]->GetUseTPCCutMIGeo()) ||
3082 (survivedTPCnclCut && fInclusivePIDtask[i]->GetUseTPCnclCut())))
3083 fInclusivePIDtask[i]->FillEfficiencyContainer(valueMeas, AliAnalysisTaskPID::kStepRecWithGenCutsMeasuredObs);
3091 for(Int_t it=0; it<nRecPartCuts; ++it){
3092 AliVParticle *part = dynamic_cast<AliVParticle*>(fTracksRecCuts->At(it));
3093 if(part)fQATrackHistosRecCuts->FillTrackQA( part->Eta(), TVector2::Phi_0_2pi(part->Phi()), part->Pt() );
3095 // fill inclusive tracks XXX, they have the same track cuts!
3096 AliAODTrack * inclusiveaod = dynamic_cast<AliAODTrack*>(fTracksRecCuts->At(it));
3098 if (fUseInclusivePIDtask && !isPileUpForAllInclusivePIDTasks) {
3099 Double_t dEdxTPC = tuneOnDataTPC ? pidResponse->GetTPCsignalTunedOnData(inclusiveaod)
3100 : inclusiveaod->GetTPCsignal();
3105 Bool_t survivedTPCCutMIGeo = AliAnalysisTaskPID::TPCCutMIGeo(inclusiveaod, InputEvent());
3106 Bool_t survivedTPCnclCut = AliAnalysisTaskPID::TPCnclCut(inclusiveaod);
3108 Int_t label = TMath::Abs(inclusiveaod->GetLabel());
3110 // find MC track in our list, if available
3111 AliAODMCParticle* gentrack = tca ? dynamic_cast<AliAODMCParticle*>(tca->At(label)) : 0x0;
3115 pdg = gentrack->GetPdgCode();
3117 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
3118 if (!isPileUpInclusivePIDtask[i] &&
3119 ((!fInclusivePIDtask[i]->GetUseTPCCutMIGeo() && !fInclusivePIDtask[i]->GetUseTPCnclCut()) ||
3120 (survivedTPCCutMIGeo && fInclusivePIDtask[i]->GetUseTPCCutMIGeo()) ||
3121 (survivedTPCnclCut && fInclusivePIDtask[i]->GetUseTPCnclCut()))) {
3122 if (fInclusivePIDtask[i]->IsInAcceptedEtaRange(TMath::Abs(inclusiveaod->Eta())))
3123 fInclusivePIDtask[i]->ProcessTrack(inclusiveaod, pdg, centPercent, -1); // no jet pT since inclusive spectrum
3128 Int_t mcID = AliAnalysisTaskPID::PDGtoMCID(pdg);
3129 Double_t valueRecAllCuts[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), inclusiveaod->Pt(), inclusiveaod->Eta(),
3130 static_cast<Double_t>(inclusiveaod->Charge()), centPercent,
3132 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
3133 if (!isPileUpInclusivePIDtask[i] &&
3134 ((!fInclusivePIDtask[i]->GetUseTPCCutMIGeo() && !fInclusivePIDtask[i]->GetUseTPCnclCut()) ||
3135 (survivedTPCCutMIGeo && fInclusivePIDtask[i]->GetUseTPCCutMIGeo()) ||
3136 (survivedTPCnclCut && fInclusivePIDtask[i]->GetUseTPCnclCut())))
3137 fInclusivePIDtask[i]->FillEfficiencyContainer(valueRecAllCuts, AliAnalysisTaskPID::kStepRecWithRecCutsMeasuredObs);
3140 Double_t weight = IsSecondaryWithStrangeMotherMC(gentrack) ? GetMCStrangenessFactorCMS(gentrack) : 1.0;
3141 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
3142 if (!isPileUpInclusivePIDtask[i] &&
3143 ((!fInclusivePIDtask[i]->GetUseTPCCutMIGeo() && !fInclusivePIDtask[i]->GetUseTPCnclCut()) ||
3144 (survivedTPCCutMIGeo && fInclusivePIDtask[i]->GetUseTPCCutMIGeo()) ||
3145 (survivedTPCnclCut && fInclusivePIDtask[i]->GetUseTPCnclCut())))
3146 fInclusivePIDtask[i]->FillEfficiencyContainer(valueRecAllCuts,
3147 AliAnalysisTaskPID::kStepRecWithRecCutsMeasuredObsStrangenessScaled,
3151 if (gentrack->IsPhysicalPrimary()) {
3152 // AliAODMCParticle->Charge() calls TParticlePDG->Charge(), which returns the charge in units of e0 / 3
3153 Double_t valueGenAllCuts[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), gentrack->Pt(), gentrack->Eta(),
3154 gentrack->Charge() / 3., centPercent, -1, -1,
3157 Double_t valuePtResolution[AliAnalysisTaskPID::kPtResNumAxes] = { -1, gentrack->Pt(), inclusiveaod->Pt(),
3158 gentrack->Charge() / 3., centPercent };
3160 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
3161 if (!isPileUpInclusivePIDtask[i] &&
3162 ((!fInclusivePIDtask[i]->GetUseTPCCutMIGeo() && !fInclusivePIDtask[i]->GetUseTPCnclCut()) ||
3163 (survivedTPCCutMIGeo && fInclusivePIDtask[i]->GetUseTPCCutMIGeo()) ||
3164 (survivedTPCnclCut && fInclusivePIDtask[i]->GetUseTPCnclCut()))) {
3165 fInclusivePIDtask[i]->FillEfficiencyContainer(valueRecAllCuts,
3166 AliAnalysisTaskPID::kStepRecWithRecCutsMeasuredObsPrimaries);
3167 fInclusivePIDtask[i]->FillEfficiencyContainer(valueGenAllCuts,
3168 AliAnalysisTaskPID::kStepRecWithRecCutsPrimaries);
3170 fInclusivePIDtask[i]->FillPtResolution(mcID, valuePtResolution);
3178 for(Int_t it=0; it<nGenPart; ++it){
3179 AliVParticle *part = dynamic_cast<AliVParticle*>(fTracksGen->At(it));
3180 if(part)fQATrackHistosGen->FillTrackQA( part->Eta(), TVector2::Phi_0_2pi(part->Phi()), part->Pt());
3186 for(Int_t ij=0; ij<nRecJets; ++ij){
3187 AliAODJet* jet = dynamic_cast<AliAODJet*>(fJetsRec->At(ij));
3188 if(jet)fQAJetHistosRec->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt());
3193 if(fQAMode || fFFMode){
3195 for(Int_t ij=0; ij<nGenJets; ++ij){
3196 AliAODJet* jet = dynamic_cast<AliAODJet*>(fJetsGen->At(ij));
3200 if(fQAMode&2) fQAJetHistosGen->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt());
3202 if(fQAMode&2 && (ij==0)) fQAJetHistosGenLeading->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt() );
3204 if((ij==0) || !fOnlyLeadingJets){ // leading jets or all jets
3205 TList* jettracklist = new TList();
3206 Double_t sumPt = 0.;
3207 Bool_t isBadJet = kFALSE;
3209 if(GetFFRadius()<=0){
3210 GetJetTracksTrackrefs(jettracklist, jet, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet);
3212 GetJetTracksPointing(fTracksGen, jettracklist, jet, GetFFRadius(), sumPt, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet);
3215 if(GetFFMinNTracks()>0 && jettracklist->GetSize()<=GetFFMinNTracks()) isBadJet = kTRUE;
3216 if(isBadJet) continue;
3218 for(Int_t it=0; it<jettracklist->GetSize(); ++it){
3220 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(jettracklist->At(it));
3221 if(!trackVP)continue;
3222 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
3224 Float_t jetPt = jet->Pt();
3225 Float_t trackPt = trackV->Pt();
3227 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
3229 if(fFFMode && (ij==0)) fFFHistosGen->FillFF( trackPt, jetPt, incrementJetPt );
3230 if(fFFMode) fFFHistosGenInc->FillFF( trackPt, jetPt, incrementJetPt );
3231 if(it==0){ // leading track
3232 if(fFFMode) fFFHistosGenLeadingTrack->FillFF( trackPt, jetPt, kTRUE );
3235 if (fUseJetPIDtask && incrementJetPt) {
3236 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
3237 if (!isPileUpJetPIDtask[i])
3238 fJetPIDtask[i]->FillGenJets(fJetPIDtask[i]->GetCentralityPercentile(evtForCentDetermination), jetPt);
3243 Int_t mcID = AliAnalysisTaskPID::PDGtoMCID(trackVP->PdgCode());
3244 if (mcID != AliPID::kUnknown) {
3245 // WARNING: The number of jets for the different species does not make sense -> One has to take
3246 // the number of jets for ALL particles. Thus, just do not fill the num of jet histos in order
3247 // not to get confused
3248 fIDFFHistosGen[mcID]->FillFF(trackPt, jetPt, kFALSE);
3251 Int_t pidWeightedSpecies = fJetPIDtask->GetRandomParticleTypeAccordingToParticleFractions(trackPt, jetPt, centPercent, kTRUE);
3252 if (pidWeightedSpecies < 0 || pidWeightedSpecies >= AliPID::kSPECIES) {
3253 Printf("Failed to determine particle ID for track in jet (gen) -> ID FF histos not filled with this track!");
3254 Printf("Track details: trackPt %f, jetPt %f, centrality %f!", trackPt, jetPt, centPercent);
3257 // WARNING: The number of jets for the different species does not make sense -> One has to take
3258 // the number of jets for ALL particles. Thus, just do not fill the num of jet histos in order
3259 // not to get confused
3260 fIDFFHistosGen[pidWeightedSpecies]->FillFF(trackPt, jetPt, kFALSE);
3267 // Efficiency, jets - particle level
3268 if (fUseJetPIDtask && !isPileUpForAllJetPIDTasks) {
3269 AliAODMCParticle* part = dynamic_cast<AliAODMCParticle*>(jettracklist->At(it));
3271 AliError("expected ref track not found ");
3274 // Fill efficiency for generated primaries and also fill histos for generated yields (primaries + all)
3275 if (part->Eta() > fTrackEtaMax || part->Eta() < fTrackEtaMin)
3278 Int_t mcID = AliAnalysisTaskPID::PDGtoMCID(part->GetPdgCode());
3280 // Following lines are not needed - just keep other species (like casecades) - will end up in overflow bin
3281 // and only affect the efficiencies for all (i.e. not identified) what is desired!
3282 //if (mcID == AliPID::kUnknown)
3285 if (!part->IsPhysicalPrimary())
3288 // Int_t iMother = part->GetMother();
3289 // if (iMother >= 0)
3290 // continue; // Not a physical primary
3293 Double_t z = -1., xi = -1.;
3294 AliAnalysisTaskPID::GetJetTrackObservables(trackPt, jetPt, z, xi);
3296 // AliAODMCParticle->Charge() calls TParticlePDG->Charge(), which returns the charge in units of e0 / 3
3297 Double_t chargeMC = part->Charge() / 3.;
3299 if (TMath::Abs(chargeMC) < 0.01)
3300 continue; // Reject neutral particles (only relevant, if mcID is not used)
3302 Double_t valuesGenYield[AliAnalysisTaskPID::kGenYieldNumAxes] = { static_cast<Double_t>(mcID), trackPt, centPercent, jetPt, z, xi, chargeMC };
3304 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
3305 if (!isPileUpJetPIDtask[i] && fJetPIDtask[i]->IsInAcceptedEtaRange(TMath::Abs(part->Eta()))) {
3306 valuesGenYield[fJetPIDtask[i]->GetIndexOfChargeAxisGenYield()] = chargeMC;
3307 fJetPIDtask[i]->FillGeneratedYield(valuesGenYield);
3312 Double_t valuesEff[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), trackPt, part->Eta(), chargeMC,
3313 centPercent, jetPt, z, xi };
3314 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
3315 if (!isPileUpJetPIDtask[i])
3316 fJetPIDtask[i]->FillEfficiencyContainer(valuesEff, AliAnalysisTaskPID::kStepGenWithGenCuts);
3322 if(fBckgType[0]!=kBckgNone)
3323 FillBckgHistos(fBckgType[0], fTracksGen, fJetsGen, jet,
3324 fFFBckgHisto0Gen, fQABckgHisto0Gen);
3325 if(fBckgType[1]!=kBckgNone)
3326 FillBckgHistos(fBckgType[1], fTracksGen, fJetsGen, jet,
3327 fFFBckgHisto1Gen, fQABckgHisto1Gen);
3328 if(fBckgType[2]!=kBckgNone)
3329 FillBckgHistos(fBckgType[2], fTracksGen, fJetsGen, jet,
3330 fFFBckgHisto2Gen, fQABckgHisto2Gen);
3331 if(fBckgType[3]!=kBckgNone)
3332 FillBckgHistos(fBckgType[3], fTracksGen, fJetsGen, jet,
3333 fFFBckgHisto3Gen, fQABckgHisto3Gen);
3334 if(fBckgType[4]!=kBckgNone)
3335 FillBckgHistos(fBckgType[4], fTracksGen, fJetsGen, jet,
3336 fFFBckgHisto4Gen, fQABckgHisto4Gen);
3337 } // end if(fBckgMode)
3340 if(fJSMode) FillJetShape(jet, jettracklist, fProNtracksLeadingJetGen, fProDelRPtSumGen, fProDelR80pcPtGen);
3342 delete jettracklist;
3346 for(Int_t ij=0; ij<nRecJetsCuts; ++ij){
3348 AliAODJet* jet = (AliAODJet*)(fJetsRecCuts->At(ij));
3349 if(fQAMode&2) fQAJetHistosRecCuts->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt());
3350 if(fQAMode&2 && (ij==0)) fQAJetHistosRecCutsLeading->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt() );
3352 if((ij==0) || !fOnlyLeadingJets){ // leading jets or all jets
3354 fhJetPtRefMultEta5->Fill(refMult5, jet->Pt());
3355 fhJetPtRefMultEta8->Fill(refMult8, jet->Pt());
3356 fhJetPtMultPercent->Fill(centPercentPP, jet->Pt());
3358 Double_t ptFractionEmbedded = 0;
3359 AliAODJet* embeddedJet = 0;
3361 if(fBranchEmbeddedJets.Length()){ // find embedded jet
3363 Int_t indexEmbedded = -1;
3364 for(Int_t i=0; i<nEmbeddedJets; i++){
3365 if(iEmbeddedMatchIndex[i] == ij){
3367 ptFractionEmbedded = fEmbeddedPtFraction[i];
3371 fh1IndexEmbedded->Fill(indexEmbedded);
3372 fh1FractionPtEmbedded->Fill(ptFractionEmbedded);
3374 if(indexEmbedded>-1){
3376 embeddedJet = dynamic_cast<AliAODJet*>(fJetsEmbedded->At(indexEmbedded));
3377 if(!embeddedJet) continue;
3379 Double_t deltaPt = jet->Pt() - embeddedJet->Pt();
3380 Double_t deltaR = jet->DeltaR((AliVParticle*) (embeddedJet));
3382 fh2DeltaPtVsJetPtEmbedded->Fill(embeddedJet->Pt(),deltaPt);
3383 fh2DeltaPtVsRecJetPtEmbedded->Fill(jet->Pt(),deltaPt);
3384 fh1DeltaREmbedded->Fill(deltaR);
3388 // get tracks in jet
3389 TList* jettracklist = new TList();
3390 Double_t sumPt = 0.;
3391 Bool_t isBadJet = kFALSE;
3393 if(GetFFRadius()<=0){
3394 GetJetTracksTrackrefs(jettracklist, jet, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet);
3396 if(fUseEmbeddedJetAxis){
3397 if(embeddedJet) GetJetTracksPointing(fTracksRecCuts, jettracklist, embeddedJet, GetFFRadius(), sumPt, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet);
3399 else GetJetTracksPointing(fTracksRecCuts, jettracklist, jet, GetFFRadius(), sumPt, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet);
3402 if(GetFFMinNTracks()>0 && jettracklist->GetSize()<=GetFFMinNTracks()) isBadJet = kTRUE;
3404 if(isBadJet) continue;
3406 if(ptFractionEmbedded>=fCutFractionPtEmbedded){ // if no embedding: ptFraction = cutFraction = 0
3408 TClonesArray *tca = fUseJetPIDtask ? dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName())) : 0x0;
3410 for(Int_t it=0; it<jettracklist->GetSize(); ++it){
3411 AliAODTrack * aodtrack = dynamic_cast<AliAODTrack*>(jettracklist->At(it));
3412 if(!aodtrack) continue;
3414 Double_t pT = aodtrack->Pt();
3415 Float_t jetPt = jet->Pt();
3416 if(fUseEmbeddedJetPt){
3417 if(embeddedJet) jetPt = embeddedJet->Pt();
3421 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
3423 if(fFFMode && (ij==0)) fFFHistosRecCuts->FillFF(pT, jetPt, incrementJetPt);
3424 if(fFFMode) fFFHistosRecCutsInc->FillFF(pT, jetPt, incrementJetPt);
3426 if(it==0){ // leading track
3427 if(fFFMode) fFFHistosRecLeadingTrack->FillFF(pT, jetPt, kTRUE);
3430 if (fUseJetPIDtask && incrementJetPt) {
3431 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
3432 if (!isPileUpJetPIDtask[i])
3433 fJetPIDtask[i]->FillRecJets(fJetPIDtask[i]->GetCentralityPercentile(evtForCentDetermination), jetPt);
3437 if (fUseJetPIDtask && (!isPileUpForAllJetPIDTasks || fIDFFMode)) {
3438 Double_t dEdxTPC = tuneOnDataTPC ? pidResponse->GetTPCsignalTunedOnData(aodtrack)
3439 : aodtrack->GetTPCsignal();
3444 Bool_t survivedTPCCutMIGeo = AliAnalysisTaskPID::TPCCutMIGeo(aodtrack, InputEvent());
3445 Bool_t survivedTPCnclCut = AliAnalysisTaskPID::TPCnclCut(aodtrack);
3447 Int_t label = TMath::Abs(aodtrack->GetLabel());
3449 // Find MC track in our list, if available
3450 AliAODMCParticle* gentrack = tca ? dynamic_cast<AliAODMCParticle*>(tca->At(label)) : 0x0;
3452 Int_t mcID = AliPID::kUnknown;
3455 pdg = gentrack->GetPdgCode();
3457 // Secondaries, jets
3458 mcID = AliAnalysisTaskPID::PDGtoMCID(pdg);
3460 Double_t z = -1., xi = -1.;
3461 AliAnalysisTaskPID::GetJetTrackObservables(pT, jetPt, z, xi);
3463 Double_t valueRecAllCuts[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), pT, aodtrack->Eta(), static_cast<Double_t>(aodtrack->Charge()),
3464 centPercent, jetPt, z, xi };
3465 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
3466 if (!isPileUpJetPIDtask[i] &&
3467 ((!fJetPIDtask[i]->GetUseTPCCutMIGeo() && !fJetPIDtask[i]->GetUseTPCnclCut()) ||
3468 (survivedTPCCutMIGeo && fJetPIDtask[i]->GetUseTPCCutMIGeo()) ||
3469 (survivedTPCnclCut && fJetPIDtask[i]->GetUseTPCnclCut())))
3470 fJetPIDtask[i]->FillEfficiencyContainer(valueRecAllCuts, AliAnalysisTaskPID::kStepRecWithRecCutsMeasuredObs);
3473 Double_t weight = IsSecondaryWithStrangeMotherMC(gentrack) ? GetMCStrangenessFactorCMS(gentrack) : 1.0;
3474 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
3475 if (!isPileUpJetPIDtask[i] &&
3476 ((!fJetPIDtask[i]->GetUseTPCCutMIGeo() && !fJetPIDtask[i]->GetUseTPCnclCut()) ||
3477 (survivedTPCCutMIGeo && fJetPIDtask[i]->GetUseTPCCutMIGeo()) ||
3478 (survivedTPCnclCut && fJetPIDtask[i]->GetUseTPCnclCut())))
3479 fJetPIDtask[i]->FillEfficiencyContainer(valueRecAllCuts,
3480 AliAnalysisTaskPID::kStepRecWithRecCutsMeasuredObsStrangenessScaled,
3484 if (gentrack->IsPhysicalPrimary()) {
3485 // AliAODMCParticle->Charge() calls TParticlePDG->Charge(), which returns the charge in units of e0 / 3
3486 Double_t genPt = gentrack->Pt();
3487 Double_t genZ = -1., genXi = -1.;
3488 AliAnalysisTaskPID::GetJetTrackObservables(genPt, jetPt, genZ, genXi);
3489 Double_t valueGenAllCuts[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), genPt, gentrack->Eta(),
3490 gentrack->Charge() / 3., centPercent, jetPt, genZ,
3493 Double_t valuePtResolution[AliAnalysisTaskPID::kPtResNumAxes] = { jetPt, genPt, pT, gentrack->Charge() / 3., centPercent };
3495 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
3496 if (!isPileUpJetPIDtask[i] &&
3497 ((!fJetPIDtask[i]->GetUseTPCCutMIGeo() && !fJetPIDtask[i]->GetUseTPCnclCut()) ||
3498 (survivedTPCCutMIGeo && fJetPIDtask[i]->GetUseTPCCutMIGeo()) ||
3499 (survivedTPCnclCut && fJetPIDtask[i]->GetUseTPCnclCut()))) {
3500 fJetPIDtask[i]->FillEfficiencyContainer(valueRecAllCuts,
3501 AliAnalysisTaskPID::kStepRecWithRecCutsMeasuredObsPrimaries);
3502 fJetPIDtask[i]->FillEfficiencyContainer(valueGenAllCuts,
3503 AliAnalysisTaskPID::kStepRecWithRecCutsPrimaries);
3505 fJetPIDtask[i]->FillPtResolution(mcID, valuePtResolution);
3511 Bool_t filledDCA = kFALSE;
3513 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
3514 if (!isPileUpJetPIDtask[i] &&
3515 ((!fJetPIDtask[i]->GetUseTPCCutMIGeo() && !fJetPIDtask[i]->GetUseTPCnclCut()) ||
3516 (survivedTPCCutMIGeo && fJetPIDtask[i]->GetUseTPCCutMIGeo()) ||
3517 (survivedTPCnclCut && fJetPIDtask[i]->GetUseTPCnclCut()))) {
3518 if (fJetPIDtask[i]->IsInAcceptedEtaRange(TMath::Abs(aodtrack->Eta()))) {
3519 fJetPIDtask[i]->ProcessTrack(aodtrack, pdg, centPercent, jetPt);
3521 // Fill DCA histo (once) if at least one of the PID tasks accecpts the track
3525 Double_t dca[2] = {0., 0.}; // 0: xy; 1: z
3526 if (aodtrack->IsGlobalConstrained()) {
3527 dca[0] = aodtrack->DCA();
3528 dca[1] = aodtrack->ZAtDCA();
3531 Double_t v[3] = {0, };
3532 Double_t pos[3] = {0, };
3534 aodtrack->GetXYZ(pos);
3535 dca[0] = pos[0] - v[0];
3536 dca[1] = pos[1] - v[1];
3539 // "Unidentified" for data and MC
3540 fhDCA_XY->Fill(pT, dca[0]);
3541 fhDCA_Z->Fill(pT, dca[1]);
3543 // "Identified" for MC
3544 if (gentrack && mcID != AliPID::kUnknown) {
3546 if (gentrack->IsPhysicalPrimary()) {
3547 fhDCA_XY_prim_MCID[mcID]->Fill(pT, dca[0]);
3548 fhDCA_Z_prim_MCID[mcID]->Fill(pT, dca[1]);
3551 fhDCA_XY_sec_MCID[mcID]->Fill(pT, dca[0]);
3552 fhDCA_Z_sec_MCID[mcID]->Fill(pT, dca[1]);
3560 if (fIDFFMode && ((!fJetPIDtask[0]->GetUseTPCCutMIGeo() && !fJetPIDtask[0]->GetUseTPCnclCut()) ||
3561 (survivedTPCCutMIGeo && fJetPIDtask[0]->GetUseTPCCutMIGeo()) ||
3562 (survivedTPCnclCut && fJetPIDtask[0]->GetUseTPCnclCut()))) {
3563 // NOTE: Just take particle fraction from first task (should anyway be the same for all tasks)
3564 Int_t pidWeightedSpecies = fJetPIDtask[0]->GetRandomParticleTypeAccordingToParticleFractions(pT, jetPt,
3565 centPercent, kTRUE);
3566 if (pidWeightedSpecies < 0 || pidWeightedSpecies >= AliPID::kSPECIES) {
3567 Printf("Failed to determine particle ID for track in jet (recCuts) -> ID FF histos not filled with this track!");
3568 Printf("Track details: trackPt %f, jetPt %f, centrality %f!", pT, jetPt, centPercent);
3571 // WARNING: The number of jets for the different species does not make sense -> One has to take
3572 // the number of jets for ALL particles. Thus, just do not fill the num of jet histos in order
3573 // not to get confused
3574 fIDFFHistosRecCuts[pidWeightedSpecies]->FillFF(aodtrack->Pt(), jetPt, kFALSE);
3578 // Efficiency, jets - detector level
3580 // Following lines are not needed - just keep other species (like casecades) - will end up in overflow bin
3581 // and only affect the efficiencies for all (i.e. not identified) what is desired!
3582 //if (mcID == AliPID::kUnknown)
3585 // Fill efficiency for reconstructed primaries
3586 if (!gentrack->IsPhysicalPrimary())
3589 Int_t iMother = gentrack->GetMother();
3591 continue; // Not a physical primary
3594 if (gentrack->Eta() > fTrackEtaMax || gentrack->Eta() < fTrackEtaMin)
3597 Double_t genZ = -1., genXi = -1.;
3598 Double_t genPt = gentrack->Pt();
3599 AliAnalysisTaskPID::GetJetTrackObservables(genPt, jetPt, genZ, genXi);
3601 Double_t measZ = -1., measXi = -1.;
3602 Double_t measPt = aodtrack->Pt();
3603 AliAnalysisTaskPID::GetJetTrackObservables(measPt, jetPt, measZ, measXi);
3606 // AliAODMCParticle->Charge() calls TParticlePDG->Charge(), which returns the charge in units of e0 / 3
3607 Double_t value[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), genPt, gentrack->Eta(), gentrack->Charge() / 3.,
3608 centPercent, jetPt, genZ, genXi };
3609 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
3610 if (!isPileUpJetPIDtask[i] &&
3611 ((!fJetPIDtask[i]->GetUseTPCCutMIGeo() && !fJetPIDtask[i]->GetUseTPCnclCut()) ||
3612 (survivedTPCCutMIGeo && fJetPIDtask[i]->GetUseTPCCutMIGeo()) ||
3613 (survivedTPCnclCut && fJetPIDtask[i]->GetUseTPCnclCut())))
3614 fJetPIDtask[i]->FillEfficiencyContainer(value, AliAnalysisTaskPID::kStepRecWithGenCuts);
3617 Double_t valueMeas[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), measPt, aodtrack->Eta(), static_cast<Double_t>(aodtrack->Charge()),
3618 centPercent, jetPt, measZ, measXi };
3619 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
3620 if (!isPileUpJetPIDtask[i] &&
3621 ((!fJetPIDtask[i]->GetUseTPCCutMIGeo() && !fJetPIDtask[i]->GetUseTPCnclCut()) ||
3622 (survivedTPCCutMIGeo && fJetPIDtask[i]->GetUseTPCCutMIGeo()) ||
3623 (survivedTPCnclCut && fJetPIDtask[i]->GetUseTPCnclCut())))
3624 fJetPIDtask[i]->FillEfficiencyContainer(valueMeas, AliAnalysisTaskPID::kStepRecWithGenCutsMeasuredObs);
3631 if(fBckgMode && (ij==0)){
3632 if(fBckgType[0]!=kBckgNone)
3633 FillBckgHistos(fBckgType[0], fTracksRecCuts, fJetsRecCuts, jet,
3634 fFFBckgHisto0RecCuts,fQABckgHisto0RecCuts, fh1BckgMult0);
3635 if(fBckgType[1]!=kBckgNone)
3636 FillBckgHistos(fBckgType[1], fTracksRecCuts, fJetsRecCuts, jet,
3637 fFFBckgHisto1RecCuts,fQABckgHisto1RecCuts, fh1BckgMult1);
3638 if(fBckgType[2]!=kBckgNone)
3639 FillBckgHistos(fBckgType[2], fTracksRecCuts, fJetsRecCuts, jet,
3640 fFFBckgHisto2RecCuts,fQABckgHisto2RecCuts, fh1BckgMult1);
3641 if(fBckgType[3]!=kBckgNone)
3642 FillBckgHistos(fBckgType[3], fTracksRecCuts, fJetsRecCuts, jet,
3643 fFFBckgHisto3RecCuts,fQABckgHisto3RecCuts, fh1BckgMult2);
3644 if(fBckgType[4]!=kBckgNone)
3645 FillBckgHistos(fBckgType[4], fTracksRecCuts, fJetsRecCuts, jet,
3646 fFFBckgHisto4RecCuts,fQABckgHisto4RecCuts, fh1BckgMult3);
3647 } // end if(fBckgMode)
3650 if(fJSMode && (ij==0)) FillJetShape(jet, jettracklist, fProNtracksLeadingJet, fProDelRPtSum, fProDelR80pcPt);
3652 delete jettracklist;
3654 } // end: cut embedded ratio
3655 } // end: leading jet or all jets
3656 } // end: rec. jets after cuts
3657 } // end: QA, FF and intra-jet
3660 // ____ efficiency _______________________________
3662 if(fEffMode && (fJetTypeRecEff != kJetsUndef)){
3664 // arrays holding for each generated particle the reconstructed AOD track index & isPrimary flag, are initialized in AssociateGenRec(...) function
3668 // array holding for each reconstructed AOD track generated particle index, initialized in AssociateGenRec(...) function
3671 // ... and another set for secondaries from strange/non strange mothers (secondary MC tracks are stored in different lists)
3672 TArrayI indexAODTrSecNS;
3674 TArrayI indexMCTrSecNS;
3676 TArrayI indexAODTrSecS;
3678 TArrayI indexMCTrSecS;
3680 Int_t nTracksAODMCCharged = GetListOfTracks(fTracksAODMCCharged, kTrackAODMCCharged);
3681 if(fDebug>2)Printf("%s:%d selected AODMC tracks: %d ",(char*)__FILE__,__LINE__,nTracksAODMCCharged);
3683 Int_t nTracksAODMCChargedSecNS = GetListOfTracks(fTracksAODMCChargedSecNS, kTrackAODMCChargedSecNS);
3684 if(fDebug>2)Printf("%s:%d selected AODMC secondary tracks NS: %d ",(char*)__FILE__,__LINE__,nTracksAODMCChargedSecNS);
3686 Int_t nTracksAODMCChargedSecS = GetListOfTracks(fTracksAODMCChargedSecS, kTrackAODMCChargedSecS);
3687 if(fDebug>2)Printf("%s:%d selected AODMC secondary tracks S: %d ",(char*)__FILE__,__LINE__,nTracksAODMCChargedSecS);
3689 Int_t nTracksRecQualityCuts = GetListOfTracks(fTracksRecQualityCuts, kTrackAODQualityCuts);
3690 if(fDebug>2)Printf("%s:%d selected rec tracks quality after cuts, full acceptance/pt : %d ",(char*)__FILE__,__LINE__,nTracksRecQualityCuts);
3692 // associate gen and rec tracks, store indices in TArrays
3693 AssociateGenRec(fTracksAODMCCharged,fTracksRecQualityCuts,indexAODTr,indexMCTr,isGenPrim,fh2PtRecVsGenPrim);
3694 AssociateGenRec(fTracksAODMCChargedSecNS,fTracksRecQualityCuts,indexAODTrSecNS,indexMCTrSecNS,isGenSecNS,fh2PtRecVsGenSec);
3695 AssociateGenRec(fTracksAODMCChargedSecS,fTracksRecQualityCuts,indexAODTrSecS,indexMCTrSecS,isGenSecS,fh2PtRecVsGenSec);
3698 if(fQAMode&1) FillSingleTrackHistosRecGen(fQATrackHistosRecEffGen,fQATrackHistosRecEffRec,fTracksAODMCCharged,indexAODTr,isGenPrim);
3701 if(fQAMode&1) FillSingleTrackHistosRecGen(0x0,fQATrackHistosSecRecNS,fTracksAODMCChargedSecNS,indexAODTrSecNS,isGenSecNS);
3702 if(fQAMode&1) FillSingleTrackHistosRecGen(0x0,fQATrackHistosSecRecS,fTracksAODMCChargedSecS,indexAODTrSecS,isGenSecS);
3703 if(fQAMode&1) FillSingleTrackHistosRecGen(0x0,fQATrackHistosSecRecSsc,fTracksAODMCChargedSecS,indexAODTrSecS,isGenSecS,kTRUE);
3707 Double_t sumPtGenLeadingJetRecEff = 0;
3708 Double_t sumPtGenLeadingJetSec = 0;
3709 Double_t sumPtRecLeadingJetRecEff = 0;
3711 for(Int_t ij=0; ij<nRecEffJets; ++ij){ // jet loop
3713 AliAODJet* jet = (AliAODJet*)(fJetsRecEff->At(ij));
3715 Bool_t isBadJetGenPrim = kFALSE;
3716 Bool_t isBadJetGenSec = kFALSE;
3717 Bool_t isBadJetRec = kFALSE;
3720 if((ij==0) || !fOnlyLeadingJets){ // leading jets or all jets
3722 // for efficiency: gen tracks from pointing with gen/rec jet
3723 TList* jettracklistGenPrim = new TList();
3725 // if radius<0 -> trackRefs: collect gen tracks in wide radius + fill FF recEff rec histos with tracks contained in track refs
3726 // note : FF recEff gen histos will be somewhat useless in this approach
3728 if(GetFFRadius() >0)
3729 GetJetTracksPointing(fTracksAODMCCharged, jettracklistGenPrim, jet, GetFFRadius(), sumPtGenLeadingJetRecEff, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetGenPrim);
3731 GetJetTracksPointing(fTracksAODMCCharged, jettracklistGenPrim, jet, TMath::Abs(GetFFRadius())+0.2, sumPtGenLeadingJetRecEff, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetGenPrim);
3733 TList* jettracklistGenSecNS = new TList();
3734 if(GetFFRadius() >0)
3735 GetJetTracksPointing(fTracksAODMCChargedSecNS, jettracklistGenSecNS, jet, GetFFRadius(), sumPtGenLeadingJetSec, GetFFMinLTrackPt() , GetFFMaxTrackPt(), isBadJetGenSec);
3737 GetJetTracksPointing(fTracksAODMCChargedSecNS, jettracklistGenSecNS, jet, TMath::Abs(GetFFRadius())+0.2, sumPtGenLeadingJetSec, GetFFMinLTrackPt() , GetFFMaxTrackPt(), isBadJetGenSec);
3739 TList* jettracklistGenSecS = new TList();
3740 if(GetFFRadius() >0)
3741 GetJetTracksPointing(fTracksAODMCChargedSecS, jettracklistGenSecS, jet, GetFFRadius(), sumPtGenLeadingJetSec, GetFFMinLTrackPt() , GetFFMaxTrackPt(), isBadJetGenSec);
3743 GetJetTracksPointing(fTracksAODMCChargedSecS, jettracklistGenSecS, jet, TMath::Abs(GetFFRadius())+0.2, sumPtGenLeadingJetSec, GetFFMinLTrackPt() , GetFFMaxTrackPt(), isBadJetGenSec);
3746 // bin efficiency in jet pt bins using rec tracks
3747 TList* jettracklistRec = new TList();
3748 if(GetFFRadius() >0) GetJetTracksPointing(fTracksRecCuts,jettracklistRec, jet, GetFFRadius(), sumPtRecLeadingJetRecEff, GetFFMinLTrackPt() , GetFFMaxTrackPt(), isBadJetRec);
3749 else GetJetTracksTrackrefs(jettracklistRec, jet, GetFFMinLTrackPt() , GetFFMaxTrackPt(), isBadJetRec);
3752 Double_t jetEta = jet->Eta();
3753 Double_t jetPhi = TVector2::Phi_0_2pi(jet->Phi());
3755 if(GetFFMinNTracks()>0 && jettracklistGenPrim->GetSize()<=GetFFMinNTracks()) isBadJetGenPrim = kTRUE;
3756 if(GetFFMinNTracks()>0 && jettracklistGenSecNS->GetSize()<=GetFFMinNTracks()) isBadJetGenSec = kTRUE;
3757 if(GetFFMinNTracks()>0 && jettracklistRec->GetSize()<=GetFFMinNTracks()) isBadJetRec = kTRUE;
3759 if(isBadJetRec) continue;
3761 if(fQAMode&2) fQAJetHistosRecEffLeading->FillJetQA( jetEta, jetPhi, sumPtGenLeadingJetRecEff );
3765 if(GetFFRadius()>0) FillJetTrackHistosRec(fFFHistosRecEffRec,jet,
3766 jettracklistGenPrim,fTracksAODMCCharged,fTracksRecQualityCuts,indexAODTr,isGenPrim,
3767 0,kFALSE,fJSMode,fProNtracksLeadingJetRecPrim,fProDelRPtSumRecPrim,fProDelR80pcPtRecPrim);
3769 else FillJetTrackHistosRec(fFFHistosRecEffRec,jet,
3770 jettracklistGenPrim,fTracksAODMCCharged,fTracksRecQualityCuts,indexAODTr,isGenPrim,
3771 jettracklistRec,kFALSE,fJSMode,fProNtracksLeadingJetRecPrim,fProDelRPtSumRecPrim,fProDelR80pcPtRecPrim);
3774 // secondaries: use jet pt from primaries
3775 if(GetFFRadius()>0) FillJetTrackHistosRec(fFFHistosSecRecNS,jet,
3776 jettracklistGenSecNS,fTracksAODMCChargedSecNS,fTracksRecQualityCuts, indexAODTrSecNS,isGenSecNS,
3777 0,kFALSE,fJSMode,fProNtracksLeadingJetRecSecNS,fProDelRPtSumRecSecNS);
3779 else FillJetTrackHistosRec(fFFHistosSecRecNS,jet,
3780 jettracklistGenSecNS,fTracksAODMCChargedSecNS,fTracksRecQualityCuts,indexAODTrSecNS,isGenSecNS,
3781 jettracklistRec,kFALSE,fJSMode,fProNtracksLeadingJetRecSecNS,fProDelRPtSumRecSecNS);
3783 if(GetFFRadius()>0) FillJetTrackHistosRec(fFFHistosSecRecS,jet,
3784 jettracklistGenSecS,fTracksAODMCChargedSecS,fTracksRecQualityCuts,indexAODTrSecS,isGenSecS,
3785 0,kFALSE,fJSMode,fProNtracksLeadingJetRecSecS,fProDelRPtSumRecSecS);
3787 else FillJetTrackHistosRec(fFFHistosSecRecS,jet,
3788 jettracklistGenSecS,fTracksAODMCChargedSecS,fTracksRecQualityCuts,indexAODTrSecS,isGenSecS,
3789 jettracklistRec,kFALSE,fJSMode,fProNtracksLeadingJetRecSecS,fProDelRPtSumRecSecS);
3791 if(GetFFRadius()>0) FillJetTrackHistosRec(fFFHistosSecRecSsc,jet,
3792 jettracklistGenSecS,fTracksAODMCChargedSecS,fTracksRecQualityCuts,indexAODTrSecS,isGenSecS,
3793 0,kTRUE,fJSMode,fProNtracksLeadingJetRecSecSsc,fProDelRPtSumRecSecSsc);
3795 else FillJetTrackHistosRec(fFFHistosSecRecSsc,jet,
3796 jettracklistGenSecS,fTracksAODMCChargedSecS,fTracksRecQualityCuts,indexAODTrSecS,isGenSecS,
3797 jettracklistRec,kTRUE,fJSMode,fProNtracksLeadingJetRecSecSsc,fProDelRPtSumRecSecSsc);
3800 delete jettracklistGenPrim;
3801 delete jettracklistGenSecNS;
3802 delete jettracklistGenSecS;
3803 delete jettracklistRec;
3806 if(fBckgMode && fFFMode){
3808 TList* perpjettracklistGen = new TList();
3809 TList* perpjettracklistGen1 = new TList();
3810 TList* perpjettracklistGen2 = new TList();
3812 //Double_t sumPtGenPerp = 0.;
3813 Double_t sumPtGenPerp1 = 0.;
3814 Double_t sumPtGenPerp2 = 0.;
3815 GetTracksTiltedwrpJetAxis(TMath::Pi()/2.,fTracksAODMCCharged, perpjettracklistGen1, jet, TMath::Abs(GetFFRadius()) , sumPtGenPerp1);
3816 GetTracksTiltedwrpJetAxis(-1*TMath::Pi()/2.,fTracksAODMCCharged, perpjettracklistGen2, jet, TMath::Abs(GetFFRadius()) ,
3819 perpjettracklistGen->AddAll(perpjettracklistGen1);
3820 perpjettracklistGen->AddAll(perpjettracklistGen2);
3821 //sumPtGenPerp = 0.5*(sumPtGenPerp1+sumPtGenPerp2);
3823 TList* perpjettracklistGenSecNS = new TList();
3824 TList* perpjettracklistGenSecNS1 = new TList();
3825 TList* perpjettracklistGenSecNS2 = new TList();
3827 //Double_t sumPtGenPerpNS = 0.;
3828 Double_t sumPtGenPerpNS1 = 0.;
3829 Double_t sumPtGenPerpNS2 = 0.;
3830 GetTracksTiltedwrpJetAxis(TMath::Pi()/2.,fTracksAODMCChargedSecNS, perpjettracklistGenSecNS1, jet, TMath::Abs(GetFFRadius()) ,
3832 GetTracksTiltedwrpJetAxis(-1*TMath::Pi()/2.,fTracksAODMCChargedSecNS, perpjettracklistGenSecNS2, jet, TMath::Abs(GetFFRadius()) ,
3835 perpjettracklistGenSecNS->AddAll(perpjettracklistGenSecNS1);
3836 perpjettracklistGenSecNS->AddAll(perpjettracklistGenSecNS2);
3837 //sumPtGenPerpNS = 0.5*(sumPtGenPerpNS1+sumPtGenPerpNS2);
3840 TList* perpjettracklistGenSecS = new TList();
3841 TList* perpjettracklistGenSecS1 = new TList();
3842 TList* perpjettracklistGenSecS2 = new TList();
3844 //Double_t sumPtGenPerpS = 0.;
3845 Double_t sumPtGenPerpS1 = 0.;
3846 Double_t sumPtGenPerpS2 = 0.;
3847 GetTracksTiltedwrpJetAxis(TMath::Pi()/2.,fTracksAODMCChargedSecS, perpjettracklistGenSecS1, jet, TMath::Abs(GetFFRadius()) ,
3849 GetTracksTiltedwrpJetAxis(-1*TMath::Pi()/2.,fTracksAODMCChargedSecS, perpjettracklistGenSecS2, jet, TMath::Abs(GetFFRadius()) ,
3852 perpjettracklistGenSecS->AddAll(perpjettracklistGenSecS1);
3853 perpjettracklistGenSecS->AddAll(perpjettracklistGenSecS2);
3854 //sumPtGenPerpS = 0.5*(sumPtGenPerpS1+sumPtGenPerpS2);
3857 if(perpjettracklistGen->GetSize() != perpjettracklistGen1->GetSize() + perpjettracklistGen2->GetSize()){
3858 cout<<" ERROR: perpjettracklistGen size "<<perpjettracklistGen->GetSize()<<" perp1 "<<perpjettracklistGen1->GetSize()
3859 <<" perp2 "<<perpjettracklistGen2->GetSize()<<endl;
3863 if(perpjettracklistGenSecNS->GetSize() != perpjettracklistGenSecNS1->GetSize() + perpjettracklistGenSecNS2->GetSize()){
3864 cout<<" ERROR: perpjettracklistGenSecNS size "<<perpjettracklistGenSecNS->GetSize()<<" perp1 "<<perpjettracklistGenSecNS1->GetSize()
3865 <<" perp2 "<<perpjettracklistGenSecNS2->GetSize()<<endl;
3869 if(perpjettracklistGenSecS->GetSize() != perpjettracklistGenSecS1->GetSize() + perpjettracklistGenSecS2->GetSize()){
3870 cout<<" ERROR: perpjettracklistGenSecS size "<<perpjettracklistGenSecS->GetSize()<<" perp1 "<<perpjettracklistGenSecS1->GetSize()
3871 <<" perp2 "<<perpjettracklistGenSecS2->GetSize()<<endl;
3876 FillJetTrackHistosRec(fFFBckgHisto0RecEffRec,jet,
3877 perpjettracklistGen,fTracksAODMCCharged,fTracksRecQualityCuts,indexAODTr,isGenPrim);
3879 FillJetTrackHistosRec(fFFBckgHisto0SecRecNS,jet,
3880 perpjettracklistGenSecNS,fTracksAODMCChargedSecNS,fTracksRecQualityCuts,indexAODTrSecNS,isGenSecNS);
3882 FillJetTrackHistosRec(fFFBckgHisto0SecRecS,jet,
3883 perpjettracklistGenSecS,fTracksAODMCChargedSecS,fTracksRecQualityCuts,indexAODTrSecS,isGenSecS);
3885 FillJetTrackHistosRec(fFFBckgHisto0SecRecSsc,jet,
3886 perpjettracklistGenSecS,fTracksAODMCChargedSecS,fTracksRecQualityCuts,indexAODTrSecS,isGenSecS,0,kTRUE);
3889 delete perpjettracklistGen;
3890 delete perpjettracklistGen1;
3891 delete perpjettracklistGen2;
3893 delete perpjettracklistGenSecNS;
3894 delete perpjettracklistGenSecNS1;
3895 delete perpjettracklistGenSecNS2;
3897 delete perpjettracklistGenSecS;
3898 delete perpjettracklistGenSecS1;
3899 delete perpjettracklistGenSecS2;
3906 //___________________
3908 fTracksRecCuts->Clear();
3909 fTracksRecCutsEfficiency->Clear();
3910 fTracksGen->Clear();
3911 fTracksAODMCCharged->Clear();
3912 fTracksAODMCChargedSecNS->Clear();
3913 fTracksAODMCChargedSecS->Clear();
3914 fTracksRecQualityCuts->Clear();
3917 fJetsRecCuts->Clear();
3919 fJetsRecEff->Clear();
3920 fJetsEmbedded->Clear();
3924 (fBckgType[0]==kBckgClusters || fBckgType[1]==kBckgClusters || fBckgType[2]==kBckgClusters || fBckgType[3]==kBckgClusters || fBckgType[4]==kBckgClusters ||
3925 fBckgType[0]==kBckgClustersOutLeading || fBckgType[1]==kBckgClustersOutLeading || fBckgType[2]==kBckgClustersOutLeading ||
3926 fBckgType[3]==kBckgClustersOutLeading || fBckgType[4]==kBckgClustersOutLeading)){
3928 fBckgJetsRec->Clear();
3929 fBckgJetsRecCuts->Clear();
3930 fBckgJetsGen->Clear();
3935 PostData(1, fCommonHistList);
3937 if (fUseJetPIDtask) {
3938 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
3939 fJetPIDtask[i]->PostOutputData();
3943 if (fUseInclusivePIDtask) {
3944 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
3945 fInclusivePIDtask[i]->PostOutputData();
3950 //______________________________________________________________
3951 void AliAnalysisTaskIDFragmentationFunction::Terminate(Option_t *)
3955 if(fDebug > 1) printf("AliAnalysisTaskIDFragmentationFunction::Terminate() \n");
3958 //_________________________________________________________________________________
3959 Int_t AliAnalysisTaskIDFragmentationFunction::GetListOfTracks(TList *list, Int_t type)
3961 // fill list of tracks selected according to type
3963 if(fDebug > 2) Printf("%s:%d Selecting tracks with %d", (char*)__FILE__,__LINE__,type);
3966 if(fDebug>1) Printf("%s:%d no input list", (char*)__FILE__,__LINE__);
3970 if(!fAOD) return -1;
3972 if(!fAOD->GetTracks()) return 0;
3974 if(type==kTrackUndef) return 0;
3978 if(type==kTrackAODExtraCuts || type==kTrackAODExtraonlyCuts || type==kTrackAODExtra || type==kTrackAODExtraonly){
3980 TClonesArray *aodExtraTracks = dynamic_cast<TClonesArray*>(fAOD->FindListObject("aodExtraTracks"));
3981 if(!aodExtraTracks)return iCount;
3982 for(int it =0; it<aodExtraTracks->GetEntries(); it++) {
3983 AliVParticle *track = dynamic_cast<AliVParticle*> ((*aodExtraTracks)[it]);
3984 if (!track) continue;
3986 AliAODTrack *tr = dynamic_cast<AliAODTrack*> (track);
3989 if(type==kTrackAODExtraCuts || type==kTrackAODExtraonlyCuts){
3991 if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask))) continue;
3993 if(tr->Eta() < fTrackEtaMin || tr->Eta() > fTrackEtaMax) continue;
3994 if(tr->Phi() < fTrackPhiMin || tr->Phi() > fTrackPhiMax) continue;
3995 if(tr->Pt() < fTrackPtCut) continue;
4003 if(type==kTrackAODCuts || type==kTrackAODQualityCuts || type==kTrackAOD || type==kTrackAODExtraCuts || type==kTrackAODExtra){
4005 // all rec. tracks, esd filter mask, eta range
4007 for(Int_t it=0; it<fAOD->GetNumberOfTracks(); ++it){
4008 AliAODTrack *tr = fAOD->GetTrack(it);
4010 if(type == kTrackAODCuts || type==kTrackAODQualityCuts || type==kTrackAODExtraCuts){
4012 if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask))) continue;
4013 if(type == kTrackAODCuts){
4014 if(tr->Eta() < fTrackEtaMin || tr->Eta() > fTrackEtaMax) continue;
4015 if(tr->Phi() < fTrackPhiMin || tr->Phi() > fTrackPhiMax) continue;
4016 if(tr->Pt() < fTrackPtCut) continue;
4023 else if (type==kTrackKineAll || type==kTrackKineCharged || type==kTrackKineChargedAcceptance){
4024 // kine particles, all or rather charged
4025 if(!fMCEvent) return iCount;
4027 for(Int_t it=0; it<fMCEvent->GetNumberOfTracks(); ++it){
4028 AliMCParticle* part = (AliMCParticle*) fMCEvent->GetTrack(it);
4030 if(type == kTrackKineCharged || type == kTrackKineChargedAcceptance){
4031 if(part->Charge()==0) continue;
4033 if(type == kTrackKineChargedAcceptance &&
4034 ( part->Eta() < fTrackEtaMin
4035 || part->Eta() > fTrackEtaMax
4036 || part->Phi() < fTrackPhiMin
4037 || part->Phi() > fTrackPhiMax
4038 || part->Pt() < fTrackPtCut)) continue;
4045 else if (type==kTrackAODMCCharged || type==kTrackAODMCAll || type==kTrackAODMCChargedAcceptance || type==kTrackAODMCChargedSecNS || type==kTrackAODMCChargedSecS) {
4046 // MC particles (from AOD), physical primaries, all or rather charged or rather charged within acceptance
4047 if(!fAOD) return -1;
4049 TClonesArray *tca = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
4050 if(!tca)return iCount;
4052 for(int it=0; it<tca->GetEntriesFast(); ++it){
4053 AliAODMCParticle *part = dynamic_cast<AliAODMCParticle*>(tca->At(it));
4055 if(type != kTrackAODMCChargedSecNS && type != kTrackAODMCChargedSecS && !part->IsPhysicalPrimary())continue;
4056 if((type == kTrackAODMCChargedSecNS || type == kTrackAODMCChargedSecS) && part->IsPhysicalPrimary())continue;
4058 if (type==kTrackAODMCCharged || type==kTrackAODMCChargedAcceptance || type==kTrackAODMCChargedSecNS || type==kTrackAODMCChargedSecS){
4059 if(part->Charge()==0) continue;
4061 if(type==kTrackAODMCChargedSecNS || type==kTrackAODMCChargedSecS){
4062 Bool_t isFromStrange = kFALSE;
4063 Int_t iMother = part->GetMother();
4065 AliAODMCParticle *partM = dynamic_cast<AliAODMCParticle*>(tca->At(iMother));
4066 if(!partM) continue;
4068 Int_t codeM = TMath::Abs(partM->GetPdgCode());
4069 Int_t mfl = Int_t (codeM/ TMath::Power(10, Int_t(TMath::Log10(codeM))));
4071 if (mfl == 3 && codeM != 3) isFromStrange = kTRUE;
4072 if(codeM == 130) isFromStrange = kTRUE; // K0 long
4073 if(part->IsSecondaryFromMaterial()) isFromStrange = kFALSE; // strange resonances from hadronic showers ?
4076 // cout<<" mfl "<<mfl<<" codeM "<<partM->GetPdgCode()<<" code this track "<<part->GetPdgCode()<<endl;
4077 // cout<<" index this track "<<it<<" index daughter 0 "<<partM->GetDaughter(0)<<" 1 "<<partM->GetDaughter(1)<<endl;
4080 if(type==kTrackAODMCChargedSecNS && isFromStrange) continue;
4081 if(type==kTrackAODMCChargedSecS && !isFromStrange) continue;
4085 if(type==kTrackAODMCChargedAcceptance &&
4086 ( part->Eta() > fTrackEtaMax
4087 || part->Eta() < fTrackEtaMin
4088 || part->Phi() > fTrackPhiMax
4089 || part->Phi() < fTrackPhiMin
4090 || part->Pt() < fTrackPtCut)) continue;
4102 // _______________________________________________________________________________
4103 Int_t AliAnalysisTaskIDFragmentationFunction::GetListOfJets(TList *list, Int_t type)
4105 // fill list of jets selected according to type
4108 if(fDebug>1) Printf("%s:%d no input list", (char*)__FILE__,__LINE__);
4112 if(type == kJetsRec || type == kJetsRecAcceptance){ // reconstructed jets
4114 if(fBranchRecJets.Length()==0){
4115 Printf("%s:%d no rec jet branch specified", (char*)__FILE__,__LINE__);
4116 if(fDebug>1)fAOD->Print();
4120 TClonesArray *aodRecJets = 0;
4121 if(fBranchRecJets.Length()) aodRecJets = dynamic_cast<TClonesArray*>(fAODJets->FindListObject(fBranchRecJets.Data()));
4122 if(!aodRecJets) aodRecJets = dynamic_cast<TClonesArray*>(fAODJets->GetList()->FindObject(fBranchRecJets.Data()));
4123 if(fAODExtension&&!aodRecJets) aodRecJets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fBranchRecJets.Data()));
4126 if(fBranchRecJets.Length()) Printf("%s:%d no reconstructed jet array with name %s in AOD", (char*)__FILE__,__LINE__,fBranchRecJets.Data());
4127 if(fDebug>1)fAOD->Print();
4131 // Reorder jet pt and fill new temporary AliAODJet objects
4134 for(Int_t ij=0; ij<aodRecJets->GetEntries(); ++ij){
4136 AliAODJet *tmp = dynamic_cast<AliAODJet*>(aodRecJets->At(ij));
4139 if( tmp->Pt() < fJetPtCut ) continue;
4140 if( type == kJetsRecAcceptance &&
4141 ( tmp->Eta() < fJetEtaMin
4142 || tmp->Eta() > fJetEtaMax
4143 || tmp->Phi() < fJetPhiMin
4144 || tmp->Phi() > fJetPhiMax )) continue;
4155 else if(type == kJetsKine || type == kJetsKineAcceptance){
4161 if(fDebug>1) Printf("%s:%d no mcEvent",(char*)__FILE__,__LINE__);
4165 AliGenEventHeader* genHeader = fMCEvent->GenEventHeader();
4166 AliGenPythiaEventHeader* pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
4167 AliGenHijingEventHeader* hijingGenHeader = 0x0;
4169 if(!pythiaGenHeader){
4170 hijingGenHeader = dynamic_cast<AliGenHijingEventHeader*>(genHeader);
4172 if(!hijingGenHeader){
4173 Printf("%s:%d no pythiaGenHeader or hijingGenHeader found", (char*)__FILE__,__LINE__);
4176 TLorentzVector mom[4];
4178 hijingGenHeader->GetJets(mom[0], mom[1], mom[2], mom[3]);
4180 for(Int_t i=0; i<2; ++i){
4181 if(!mom[i].Pt()) continue;
4182 jet[i] = new AliAODJet(mom[i]);
4184 if( type == kJetsKineAcceptance &&
4185 ( jet[i]->Eta() < fJetEtaMin
4186 || jet[i]->Eta() > fJetEtaMax
4187 || jet[i]->Phi() < fJetPhiMin
4188 || jet[i]->Phi() > fJetPhiMax )) continue;
4198 // fetch the pythia generated jets
4199 for(int ip=0; ip<pythiaGenHeader->NTriggerJets(); ++ip){
4202 AliAODJet *jet = new AliAODJet();
4203 pythiaGenHeader->TriggerJet(ip, p);
4204 jet->SetPxPyPzE(p[0], p[1], p[2], p[3]);
4206 if( type == kJetsKineAcceptance &&
4207 ( jet->Eta() < fJetEtaMin
4208 || jet->Eta() > fJetEtaMax
4209 || jet->Phi() < fJetPhiMin
4210 || jet->Phi() > fJetPhiMax )) continue;
4218 else if(type == kJetsGen || type == kJetsGenAcceptance ){
4220 if(fBranchGenJets.Length()==0){
4221 if(fDebug>1) Printf("%s:%d no gen jet branch specified", (char*)__FILE__,__LINE__);
4225 TClonesArray *aodGenJets = 0;
4226 if(fBranchGenJets.Length()) aodGenJets = dynamic_cast<TClonesArray*>(fAODJets->FindListObject(fBranchGenJets.Data()));
4227 if(!aodGenJets) aodGenJets = dynamic_cast<TClonesArray*>(fAODJets->GetList()->FindObject(fBranchGenJets.Data()));
4228 if(fAODExtension&&!aodGenJets) aodGenJets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fBranchGenJets.Data()));
4232 if(fBranchGenJets.Length()) Printf("%s:%d Generated jet branch %s not found",(char*)__FILE__,__LINE__,fBranchGenJets.Data());
4234 if(fDebug>1)fAOD->Print();
4240 for(Int_t ig=0; ig<aodGenJets->GetEntries(); ++ig){
4242 AliAODJet *tmp = dynamic_cast<AliAODJet*>(aodGenJets->At(ig));
4245 if( tmp->Pt() < fJetPtCut ) continue;
4246 if( type == kJetsGenAcceptance &&
4247 ( tmp->Eta() < fJetEtaMin
4248 || tmp->Eta() > fJetEtaMax
4249 || tmp->Phi() < fJetPhiMin
4250 || tmp->Phi() > fJetPhiMax )) continue;
4258 else if(type == kJetsEmbedded){ // embedded jets
4260 if(fBranchEmbeddedJets.Length()==0){
4261 Printf("%s:%d no embedded jet branch specified", (char*)__FILE__,__LINE__);
4262 if(fDebug>1)fAOD->Print();
4266 TClonesArray *aodEmbeddedJets = 0;
4267 if(fBranchEmbeddedJets.Length()) aodEmbeddedJets = dynamic_cast<TClonesArray*>(fAODJets->FindListObject(fBranchEmbeddedJets.Data()));
4268 if(!aodEmbeddedJets) aodEmbeddedJets = dynamic_cast<TClonesArray*>(fAODJets->GetList()->FindObject(fBranchEmbeddedJets.Data()));
4269 if(fAODExtension&&!aodEmbeddedJets) aodEmbeddedJets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fBranchEmbeddedJets.Data()));
4271 if(!aodEmbeddedJets){
4272 if(fBranchEmbeddedJets.Length()) Printf("%s:%d no reconstructed jet array with name %s in AOD", (char*)__FILE__,__LINE__,fBranchEmbeddedJets.Data());
4273 if(fDebug>1)fAOD->Print();
4277 // Reorder jet pt and fill new temporary AliAODJet objects
4278 Int_t nEmbeddedJets = 0;
4280 for(Int_t ij=0; ij<aodEmbeddedJets->GetEntries(); ++ij){
4282 AliAODJet *tmp = dynamic_cast<AliAODJet*>(aodEmbeddedJets->At(ij));
4285 if( tmp->Pt() < fJetPtCut ) continue;
4286 if( tmp->Eta() < fJetEtaMin
4287 || tmp->Eta() > fJetEtaMax
4288 || tmp->Phi() < fJetPhiMin
4289 || tmp->Phi() > fJetPhiMax ) continue;
4297 return nEmbeddedJets;
4300 if(fDebug>0)Printf("%s:%d no such type %d",(char*)__FILE__,__LINE__,type);
4305 // ___________________________________________________________________________________
4306 Int_t AliAnalysisTaskIDFragmentationFunction::GetListOfBckgJets(TList *list, Int_t type)
4308 // fill list of bgr clusters selected according to type
4310 if(type == kJetsRec || type == kJetsRecAcceptance){ // reconstructed jets
4312 if(fBranchRecBckgClusters.Length()==0){
4313 Printf("%s:%d no rec jet branch specified", (char*)__FILE__,__LINE__);
4314 if(fDebug>1)fAOD->Print();
4318 TClonesArray *aodRecJets = 0;
4319 if(fBranchRecBckgClusters.Length()) aodRecJets = dynamic_cast<TClonesArray*>(fAODJets->FindListObject(fBranchRecBckgClusters.Data()));
4320 if(!aodRecJets) aodRecJets = dynamic_cast<TClonesArray*>(fAODJets->GetList()->FindObject(fBranchRecBckgClusters.Data()));
4321 if(fAODExtension&&!aodRecJets) aodRecJets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fBranchRecBckgClusters.Data()));
4324 if(fBranchRecBckgClusters.Length()) Printf("%s:%d no reconstructed jet array with name %s in AOD", (char*)__FILE__,__LINE__,fBranchRecBckgClusters.Data());
4325 if(fDebug>1)fAOD->Print();
4329 // Reorder jet pt and fill new temporary AliAODJet objects
4332 for(Int_t ij=0; ij<aodRecJets->GetEntries(); ++ij){
4334 AliAODJet *tmp = dynamic_cast<AliAODJet*>(aodRecJets->At(ij));
4337 // if( tmp->Pt() < fJetPtCut ) continue; // no pt cut on bckg clusters !
4338 if( type == kJetsRecAcceptance &&
4339 ( tmp->Eta() < fJetEtaMin
4340 || tmp->Eta() > fJetEtaMax
4341 || tmp->Phi() < fJetPhiMin
4342 || tmp->Phi() > fJetPhiMax )) continue;
4356 // MC clusters still Under construction
4362 // _________________________________________________________________________________________________________
4363 void AliAnalysisTaskIDFragmentationFunction::SetProperties(THnSparse* h, Int_t dim, const char** labels)
4365 // Set properties of THnSparse
4367 for(Int_t i=0; i<dim; i++){
4368 h->GetAxis(i)->SetTitle(labels[i]);
4369 h->GetAxis(i)->SetTitleColor(1);
4373 // __________________________________________________________________________________________
4374 void AliAnalysisTaskIDFragmentationFunction::SetProperties(TH1* h,const char* x, const char* y)
4376 //Set properties of histos (x and y title)
4380 h->GetXaxis()->SetTitleColor(1);
4381 h->GetYaxis()->SetTitleColor(1);
4384 // _________________________________________________________________________________________________________
4385 void AliAnalysisTaskIDFragmentationFunction::SetProperties(TH1* h,const char* x, const char* y, const char* z)
4387 //Set properties of histos (x,y and z title)
4392 h->GetXaxis()->SetTitleColor(1);
4393 h->GetYaxis()->SetTitleColor(1);
4394 h->GetZaxis()->SetTitleColor(1);
4397 // ________________________________________________________________________________________________________________________________________________________
4398 void AliAnalysisTaskIDFragmentationFunction::GetJetTracksPointing(TList* inputlist, TList* outputlist, const AliAODJet* jet,
4399 Double_t radius, Double_t& sumPt, Double_t minPtL, Double_t maxPt, Bool_t& isBadPt)
4401 // fill list of tracks in cone around jet axis
4404 Bool_t isBadMaxPt = kFALSE;
4405 Bool_t isBadMinPt = kTRUE;
4408 jet->PxPyPz(jetMom);
4409 TVector3 jet3mom(jetMom);
4411 for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++){
4413 AliVParticle* track = dynamic_cast<AliVParticle*>(inputlist->At(itrack));
4415 Double_t trackMom[3];
4416 track->PxPyPz(trackMom);
4417 TVector3 track3mom(trackMom);
4419 Double_t dR = jet3mom.DeltaR(track3mom);
4423 outputlist->Add(track);
4425 sumPt += track->Pt();
4427 if(maxPt>0 && track->Pt()>maxPt) isBadMaxPt = kTRUE;
4428 if(minPtL>0 && track->Pt()>minPtL) isBadMinPt = kFALSE;
4433 if(minPtL>0 && isBadMinPt) isBadPt = kTRUE;
4434 if(maxPt>0 && isBadMaxPt) isBadPt = kTRUE;
4439 // _________________________________________________________________________________________________________________________________________________________________
4440 void AliAnalysisTaskIDFragmentationFunction::GetJetTracksTrackrefs(TList* list, const AliAODJet* jet, Double_t minPtL, Double_t maxPt, Bool_t& isBadPt)
4442 // list of jet tracks from trackrefs
4444 Int_t nTracks = jet->GetRefTracks()->GetEntriesFast();
4446 Bool_t isBadMaxPt = kFALSE;
4447 Bool_t isBadMinPt = kTRUE;
4449 for(Int_t itrack=0; itrack<nTracks; itrack++) {
4451 AliVParticle* track = dynamic_cast<AliVParticle*>(jet->GetRefTracks()->At(itrack));
4453 AliError("expected ref track not found ");
4457 if(track->Pt() < fTrackPtCut) continue; // track refs may contain low pt cut (bug in AliFastJetInput)
4458 if(maxPt>0 && track->Pt()>maxPt) isBadMaxPt = kTRUE;
4459 if(minPtL>0 && track->Pt()>minPtL) isBadMinPt = kFALSE;
4465 if(minPtL>0 && isBadMinPt) isBadPt = kTRUE;
4466 if(maxPt>0 && isBadMaxPt) isBadPt = kTRUE;
4471 // _ ________________________________________________________________________________________________________________________________
4472 void AliAnalysisTaskIDFragmentationFunction::AssociateGenRec(TList* tracksAODMCCharged,TList* tracksRec, TArrayI& indexAODTr,TArrayI& indexMCTr,
4473 TArrayS& isRefGen,TH2F* fh2PtRecVsGen)
4475 // associate generated and reconstructed tracks, fill TArrays of list indices
4477 Int_t nTracksRec = tracksRec->GetSize();
4478 Int_t nTracksGen = tracksAODMCCharged->GetSize();
4479 TClonesArray *tca = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
4482 if(!nTracksGen) return;
4486 indexAODTr.Set(nTracksGen);
4487 indexMCTr.Set(nTracksRec);
4488 isRefGen.Set(nTracksGen);
4490 indexAODTr.Reset(-1);
4491 indexMCTr.Reset(-1);
4494 // loop over reconstructed tracks, get generated track
4496 for(Int_t iRec=0; iRec<nTracksRec; iRec++){
4498 AliAODTrack* rectrack = dynamic_cast<AliAODTrack*>(tracksRec->At(iRec));
4499 if(!rectrack)continue;
4500 Int_t label = TMath::Abs(rectrack->GetLabel());
4502 // find MC track in our list
4503 AliAODMCParticle* gentrack = dynamic_cast<AliAODMCParticle*> (tca->At(label));
4505 Int_t listIndex = -1;
4506 if(gentrack) listIndex = tracksAODMCCharged->IndexOf(gentrack);
4510 indexAODTr[listIndex] = iRec;
4511 indexMCTr[iRec] = listIndex;
4516 // define reference sample of primaries/secondaries (for reconstruction efficiency / contamination)
4518 for(Int_t iGen=0; iGen<nTracksGen; iGen++){
4520 AliAODMCParticle* gentrack = dynamic_cast<AliAODMCParticle*> (tracksAODMCCharged->At(iGen));
4521 if(!gentrack)continue;
4522 Int_t pdg = gentrack->GetPdgCode();
4524 // 211 - pi, 2212 - proton, 321 - Kaon, 11 - electron, 13 - muon
4525 if(TMath::Abs(pdg) == 211 || TMath::Abs(pdg) == 2212 || TMath::Abs(pdg) == 321 ||
4526 TMath::Abs(pdg) == 11 || TMath::Abs(pdg) == 13){
4528 isRefGen[iGen] = kTRUE;
4530 Int_t iRec = indexAODTr[iGen]; // can be -1 if no good reconstructed track
4533 Float_t genPt = gentrack->Pt();
4534 AliAODTrack* vt = dynamic_cast<AliAODTrack*>(tracksRec->At(iRec));
4536 Float_t recPt = vt->Pt();
4537 fh2PtRecVsGen->Fill(genPt,recPt);
4544 // _____________________________________________________________________________________________________________________________________________
4545 void AliAnalysisTaskIDFragmentationFunction::FillSingleTrackHistosRecGen(AliFragFuncQATrackHistos* trackQAGen, AliFragFuncQATrackHistos* trackQARec, TList* tracksGen,
4546 const TArrayI& indexAODTr, const TArrayS& isRefGen, Bool_t scaleStrangeness){
4548 // fill QA for single track reconstruction efficiency
4550 Int_t nTracksGen = tracksGen->GetSize();
4552 if(!nTracksGen) return;
4554 for(Int_t iGen=0; iGen<nTracksGen; iGen++){
4556 if(isRefGen[iGen] != 1) continue; // select primaries
4558 AliAODMCParticle* gentrack = dynamic_cast<AliAODMCParticle*> (tracksGen->At(iGen));
4559 if(!gentrack) continue;
4560 Double_t ptGen = gentrack->Pt();
4561 Double_t etaGen = gentrack->Eta();
4562 Double_t phiGen = TVector2::Phi_0_2pi(gentrack->Phi());
4564 // apply same acc & pt cuts as for FF GetMCStrangenessFactorCMS
4566 if(etaGen < fTrackEtaMin || etaGen > fTrackEtaMax) continue;
4567 if(phiGen < fTrackPhiMin || phiGen > fTrackPhiMax) continue;
4568 if(ptGen < fTrackPtCut) continue;
4570 if(trackQAGen) trackQAGen->FillTrackQA(etaGen, phiGen, ptGen);
4572 Int_t iRec = indexAODTr[iGen]; // can be -1 if no good reconstructed track
4574 if(iRec>=0 && trackQARec){
4575 if(scaleStrangeness){
4576 //Double_t weight = GetMCStrangenessFactor(ptGen);
4577 Double_t weight = GetMCStrangenessFactorCMS(gentrack);
4578 trackQARec->FillTrackQA(etaGen, phiGen, ptGen, kFALSE, 0, kTRUE, weight);
4580 else trackQARec->FillTrackQA(etaGen, phiGen, ptGen);
4585 // ______________________________________________________________________________________________________________________________________________________
4587 void AliAnalysisTaskIDFragmentationFunction::FillJetTrackHistosRec(AliFragFuncHistos* ffhistRec, AliAODJet* jet,
4588 TList* jetTrackList, const TList* tracksGen, const TList* tracksRec, const TArrayI& indexAODTr,
4589 const TArrayS& isRefGen, TList* jetTrackListTR, Bool_t scaleStrangeness,
4590 Bool_t fillJS, TProfile* hProNtracksLeadingJet, TProfile** hProDelRPtSum, TProfile* hProDelR80pcPt)
4592 // fill objects for jet track reconstruction efficiency or secondaries contamination
4593 // arguments histGen/histRec can be of different type: AliFragFuncHistos*, AliFragFuncIntraJetHistos*, ...
4594 // jetTrackListTR pointer: track refs if not NULL
4597 // ensure proper normalization, even for secondaries
4598 Double_t jetPtRec = jet->Pt();
4599 ffhistRec->FillFF(-1, jetPtRec, kTRUE);
4601 Int_t nTracksJet = jetTrackList->GetSize(); // list with AODMC tracks
4602 if(nTracksJet == 0) return;
4604 TList* listRecTracks = new TList();
4605 listRecTracks->Clear();
4607 for(Int_t iTr=0; iTr<nTracksJet; iTr++){ // jet tracks loop
4609 AliAODMCParticle* gentrack = dynamic_cast<AliAODMCParticle*> (jetTrackList->At(iTr));
4610 if(!gentrack)continue;
4611 // find jet track in gen tracks list
4612 Int_t iGen = tracksGen->IndexOf(gentrack);
4615 if(fDebug>0) Printf("%s:%d gen jet track not found ",(char*)__FILE__,__LINE__);
4619 if(isRefGen[iGen] != 1) continue; // select primaries
4621 Double_t ptGen = gentrack->Pt();
4622 Double_t etaGen = gentrack->Eta();
4623 Double_t phiGen = TVector2::Phi_0_2pi(gentrack->Phi());
4625 // gen level acc & pt cuts - skip in case of track refs
4626 if(!jetTrackListTR && (etaGen < fTrackEtaMin || etaGen > fTrackEtaMax)) continue;
4627 if(!jetTrackListTR && (phiGen < fTrackPhiMin || phiGen > fTrackPhiMax)) continue;
4628 if(!jetTrackListTR && ptGen < fTrackPtCut) continue;
4631 Double_t ptRec = -1;
4633 Int_t iRec = indexAODTr[iGen]; // can be -1 if no good reconstructed track
4634 Bool_t isRec = (iRec>=0) ? kTRUE : kFALSE;
4636 Bool_t isJetTrack = kFALSE;
4637 if(!jetTrackListTR) isJetTrack = kTRUE; // skip trackRefs check for tracks in ideal cone
4641 AliAODTrack* rectrack = dynamic_cast<AliAODTrack*> (tracksRec->At(iRec));
4642 if(!rectrack) continue;
4644 ptRec = rectrack->Pt();
4647 Int_t iRecTR = jetTrackListTR->IndexOf(rectrack);
4648 if(iRecTR >=0 ) isJetTrack = kTRUE; // rec tracks assigned to jet
4653 Double_t trackPt = ptRec;
4654 Bool_t incrementJetPt = kFALSE;
4656 if(scaleStrangeness){
4657 //Double_t weight = GetMCStrangenessFactor(ptGen);
4658 Double_t weight = GetMCStrangenessFactorCMS(gentrack);
4660 ffhistRec->FillFF( trackPt, jetPtRec, incrementJetPt, 0, kTRUE, weight );
4663 ffhistRec->FillFF( trackPt, jetPtRec, incrementJetPt );
4666 listRecTracks->Add(rectrack);
4673 if(fillJS) FillJetShape(jet,listRecTracks,hProNtracksLeadingJet, hProDelRPtSum, hProDelR80pcPt,0,0,scaleStrangeness);
4675 delete listRecTracks;
4679 // _____________________________________________________________________________________________________________________________________________________________________
4680 void AliAnalysisTaskIDFragmentationFunction::GetTracksTiltedwrpJetAxis(Float_t alpha, TList* inputlist, TList* outputlist, const AliAODJet* jet, Double_t radius,Double_t& sumPt)
4682 // List of tracks in cone perpendicular to the jet azimuthal direction
4685 jet->PxPyPz(jetMom);
4687 TVector3 jet3mom(jetMom);
4688 // Rotate phi and keep eta unchanged
4689 Double_t etaTilted = jet3mom.Eta();
4690 Double_t phiTilted = TVector2::Phi_0_2pi(jet3mom.Phi()) + alpha;
4691 if(phiTilted > 2*TMath::Pi()) phiTilted = phiTilted - 2*TMath::Pi();
4693 for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++){
4696 if( fUseExtraTracksBgr != 1){
4697 if(AliAODTrack* trackAOD = dynamic_cast<AliAODTrack*> (inputlist->At(itrack))){
4698 if(fUseExtraTracksBgr == 0 && ((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
4699 if(fUseExtraTracksBgr == -1 && !((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
4703 AliVParticle* track = dynamic_cast<AliVParticle*>(inputlist->At(itrack));
4705 Double_t trackMom[3];
4706 track->PxPyPz(trackMom);
4707 TVector3 track3mom(trackMom);
4709 Double_t deta = track3mom.Eta() - etaTilted;
4710 Double_t dphi = TMath::Abs(track3mom.Phi() - phiTilted);
4711 if (dphi > TMath::Pi()) dphi = 2. * TMath::Pi() - dphi;
4712 Double_t dR = TMath::Sqrt(deta * deta + dphi * dphi);
4716 outputlist->Add(track);
4717 sumPt += track->Pt();
4723 // ________________________________________________________________________________________________________________________________________________________
4724 void AliAnalysisTaskIDFragmentationFunction::GetTracksTiltedwrpJetAxisWindow(Float_t alpha, TList* inputlist, TList* outputlist, const AliAODJet* jet, Double_t radius,Double_t& sumPt,Double_t &normFactor)
4726 // List of tracks in cone perpendicular to the jet azimuthal direction
4729 jet->PxPyPz(jetMom);
4731 TVector3 jet3mom(jetMom);
4732 // Rotate phi and keep eta unchanged
4733 Double_t etaTilted = jet3mom.Eta();
4734 Double_t phiTilted = TVector2::Phi_0_2pi(jet3mom.Phi()) + alpha;
4735 if(phiTilted > 2*TMath::Pi()) phiTilted = phiTilted - 2*TMath::Pi();
4737 for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++)
4741 if( fUseExtraTracksBgr != 1){
4742 if(AliAODTrack* trackAOD = dynamic_cast<AliAODTrack*> (inputlist->At(itrack))){
4743 if(fUseExtraTracksBgr == 0 && ((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
4744 if(fUseExtraTracksBgr == -1 && !((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
4748 AliVParticle* track = dynamic_cast<AliVParticle*>(inputlist->At(itrack));
4750 Float_t trackEta = track->Eta();
4751 Float_t trackPhi = track->Phi();
4753 if( ( phiTilted-radius >= 0 ) && ( phiTilted+radius <= 2*TMath::Pi()))
4755 if((trackPhi<=phiTilted+radius) &&
4756 (trackPhi>=phiTilted-radius) &&
4757 (trackEta<=fTrackEtaMax)&&(trackEta>=fTrackEtaMin)) // 0.9 and - 0.9
4759 outputlist->Add(track);
4760 sumPt += track->Pt();
4763 else if( phiTilted-radius < 0 )
4765 if((( trackPhi < phiTilted+radius ) ||
4766 ( trackPhi > 2*TMath::Pi()-(radius-phiTilted) )) &&
4767 (( trackEta <= fTrackEtaMax ) && ( trackEta >= fTrackEtaMin )))
4769 outputlist->Add(track);
4770 sumPt += track->Pt();
4773 else if( phiTilted+radius > 2*TMath::Pi() )
4775 if((( trackPhi > phiTilted-radius ) ||
4776 ( trackPhi < phiTilted+radius-2*TMath::Pi() )) &&
4777 (( trackEta <= fTrackEtaMax ) && ( trackEta >= fTrackEtaMin )))
4779 outputlist->Add(track);
4780 sumPt += track->Pt();
4785 // Jet area - Temporarily added should be modified with the proper jet area value
4786 Float_t areaJet = CalcJetArea(etaTilted,radius);
4787 Float_t areaTilted = 2*radius*(fTrackEtaMax-fTrackEtaMin);
4789 normFactor = (Float_t) 1. / (areaJet / areaTilted);
4794 // ________________________________________________________________________________________________________________________________________________________
4795 void AliAnalysisTaskIDFragmentationFunction::GetTracksOutOfNJets(Int_t nCases, TList* inputlist, TList* outputlist, const TList* jetlist,
4798 // List of tracks outside cone around N jet axis
4799 // Particles taken randomly
4802 // Int_t nj = jetlist->GetSize();
4803 Float_t rc = TMath::Abs(GetFFRadius());
4804 Float_t rcl = GetFFBckgRadius();
4806 // Estimate jet and background areas
4807 Float_t* areaJet = new Float_t[nCases];
4808 memset(areaJet, 0, sizeof(Float_t) * nCases);
4809 Float_t* areaJetLarge = new Float_t[nCases];
4810 memset(areaJetLarge, 0, sizeof(Float_t) * nCases);
4811 Float_t areaFull = (fTrackEtaMax-fTrackEtaMin)*(fTrackPhiMax-fTrackPhiMin);
4812 Float_t areaOut = areaFull;
4814 //estimate jets and background areas
4817 TList* templist = new TList();
4818 TClonesArray *vect3Jet = new TClonesArray("TVector3",nCases);
4820 for(Int_t ij=0; ij<nCases; ++ij)
4822 // Get jet information
4823 AliAODJet* jet = dynamic_cast<AliAODJet*>(jetlist->At(ij));
4826 jet3mom.SetPtEtaPhi(jet->Pt(),jet->Eta(),jet->Phi());
4827 new((*vect3Jet)[ijet]) TVector3((TVector3)jet3mom);
4828 Float_t etaJet = (Float_t)((TVector3*) vect3Jet->At(ij))->Eta();
4831 areaJet[ij] = CalcJetArea(etaJet,rc);
4833 // Area jet larger angle
4834 areaJetLarge[ij] = CalcJetArea(etaJet,rcl);
4837 areaOut = areaOut - areaJetLarge[ij];
4841 // List of all tracks outside jet areas
4842 for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++){
4845 if( fUseExtraTracksBgr != 1){
4846 if(AliAODTrack* trackAOD = dynamic_cast<AliAODTrack*> (inputlist->At(itrack))){
4847 if(fUseExtraTracksBgr == 0 && ((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
4848 if(fUseExtraTracksBgr == -1 && !((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
4852 AliVParticle* track = dynamic_cast<AliVParticle*>(inputlist->At(itrack));
4855 Double_t trackMom[3];
4856 track->PxPyPz(trackMom);
4857 TVector3 track3mom(trackMom);
4859 Double_t *dR = new Double_t[nCases];
4860 for(Int_t ij=0; ij<nCases; ij++)
4861 dR[ij] = (Double_t)((TVector3*) vect3Jet->At(ij))->DeltaR(track3mom);
4863 if((nCases==1 && (dR[0]>rcl)) ||
4864 (nCases==2 && (dR[0]>rcl && dR[1]>rcl)) ||
4865 (nCases==3 && (dR[0]>rcl && dR[1]>rcl && dR[2]>rcl)))
4867 templist->Add(track);
4873 // Take tracks randomly
4874 Int_t nScaled = (Int_t) (nOut * areaJet[0] / areaOut + 0.5);
4875 TArrayI* ar = new TArrayI(nOut);
4877 for(Int_t init=0; init<nOut; init++)
4880 Int_t *randIndex = new Int_t[nScaled];
4881 for(Int_t init2=0; init2<nScaled; init2++)
4882 randIndex[init2] = -1;
4884 // Select nScaled different random numbers in nOut
4885 for(Int_t i=0; i<nScaled; i++)
4887 Int_t* tmpArr = new Int_t[nOut-i];
4888 Int_t temp = fRandom->Integer(nOut-i);
4889 for(Int_t ind = 0; ind< ar->GetSize()-1; ind++)
4891 if(ind<temp) tmpArr[ind] = (*ar)[ind];
4892 else tmpArr[ind] = (*ar)[ind+1];
4894 randIndex[i] = (*ar)[temp];
4896 ar->Set(nOut-i-1,tmpArr);
4902 for(Int_t ipart=0; ipart<nScaled; ipart++)
4904 AliVParticle* track = (AliVParticle*)(templist->At(randIndex[ipart]));
4905 outputlist->Add(track);
4906 sumPt += track->Pt();
4913 delete [] areaJetLarge;
4916 delete [] randIndex;
4920 // ________________________________________________________________________________________________________________________________________________________
4921 void AliAnalysisTaskIDFragmentationFunction::GetTracksOutOfNJetsStat(Int_t nCases, TList* inputlist, TList* outputlist,
4922 const TList* jetlist, Double_t& sumPt, Double_t &normFactor)
4924 // List of tracks outside cone around N jet axis
4925 // All particles taken + final scaling factor
4928 Float_t rc = TMath::Abs(GetFFRadius());
4929 Float_t rcl = GetFFBckgRadius();
4931 // Estimate jet and background areas
4932 Float_t* areaJet = new Float_t[nCases];
4933 memset(areaJet, 0, sizeof(Float_t) * nCases);
4934 Float_t* areaJetLarge = new Float_t[nCases];
4935 memset(areaJetLarge, 0, sizeof(Float_t) * nCases);
4936 Float_t areaFull = (fTrackEtaMax-fTrackEtaMin)*(fTrackPhiMax-fTrackPhiMin);
4937 Float_t areaOut = areaFull;
4939 //estimate jets and background areas
4942 TClonesArray *vect3Jet = new TClonesArray("TVector3",nCases);
4944 for(Int_t ij=0; ij<nCases; ++ij)
4946 // Get jet information
4947 AliAODJet* jet = dynamic_cast<AliAODJet*>(jetlist->At(ij));
4950 jet3mom.SetPtEtaPhi(jet->Pt(),jet->Eta(),jet->Phi());
4951 new((*vect3Jet)[ijet]) TVector3((TVector3)jet3mom);
4952 Float_t etaJet = (Float_t)((TVector3*) vect3Jet->At(ij))->Eta();
4955 areaJet[ij] = CalcJetArea(etaJet,rc);
4957 // Area jet larger angle
4958 areaJetLarge[ij] = CalcJetArea(etaJet,rcl);
4960 // Outside jets area
4961 areaOut = areaOut - areaJetLarge[ij];
4965 for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++){
4968 if( fUseExtraTracksBgr != 1){
4969 if(AliAODTrack* trackAOD = dynamic_cast<AliAODTrack*> (inputlist->At(itrack))){
4970 if(fUseExtraTracksBgr == 0 && ((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
4971 if(fUseExtraTracksBgr == -1 && !((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
4975 AliVParticle* track = dynamic_cast<AliVParticle*>(inputlist->At(itrack));
4977 Double_t trackMom[3];
4978 track->PxPyPz(trackMom);
4979 TVector3 track3mom(trackMom);
4981 Double_t *dR = new Double_t[nCases];
4982 for(Int_t ij=0; ij<nCases; ij++)
4983 dR[ij] = (Double_t)((TVector3*) vect3Jet->At(ij))->DeltaR(track3mom);
4986 (nCases==1 && (dR[0]>rcl)) ||
4987 (nCases==2 && (dR[0]>rcl && dR[1]>rcl)) ||
4988 (nCases==3 && (dR[0]>rcl && dR[1]>rcl && dR[2]>rcl)))
4990 outputlist->Add(track);
4991 sumPt += track->Pt();
4997 if(nCases==0) areaJet[0] = TMath::Pi()*rc*rc;
4998 normFactor = (Float_t) 1./(areaJet[0] / areaOut);
5003 delete [] areaJetLarge;
5008 // ______________________________________________________________________________________________________________________________________________________
5009 Float_t AliAnalysisTaskIDFragmentationFunction::CalcJetArea(Float_t etaJet, Float_t rc) const
5011 // calculate area of jet with eta etaJet and radius rc
5013 Float_t detamax = etaJet + rc;
5014 Float_t detamin = etaJet - rc;
5015 Float_t accmax = 0.0; Float_t accmin = 0.0;
5016 if(detamax > fTrackEtaMax){ // sector outside etamax
5017 Float_t h = fTrackEtaMax - etaJet;
5018 accmax = rc*rc*TMath::ACos(h/rc) - h*TMath::Sqrt(rc*rc - h*h);
5020 if(detamin < fTrackEtaMin){ // sector outside etamin
5021 Float_t h = fTrackEtaMax + etaJet;
5022 accmin = rc*rc*TMath::ACos(h/rc) - h*TMath::Sqrt(rc*rc - h*h);
5024 Float_t areaJet = rc*rc*TMath::Pi() - accmax - accmin;
5030 // ___________________________________________________________________________________________________________________________
5031 void AliAnalysisTaskIDFragmentationFunction::GetClusterTracksOutOf1Jet(AliAODJet* jet, TList* outputlist, Double_t &normFactor)
5033 // fill tracks from bckgCluster branch in list,
5034 // for all clusters outside jet cone
5035 // sum up total area of clusters
5037 Double_t rc = GetFFRadius();
5038 Double_t rcl = GetFFBckgRadius();
5040 Double_t areaTotal = 0;
5041 Double_t sumPtTotal = 0;
5043 for(Int_t ij=0; ij<fBckgJetsRec->GetEntries(); ++ij){
5045 AliAODJet* bgrCluster = (AliAODJet*)(fBckgJetsRec->At(ij)); // not 'recCuts': use all clusters in full eta range
5047 Double_t dR = jet->DeltaR(bgrCluster);
5049 if(dR<rcl) continue;
5051 Double_t clusterPt = bgrCluster->Pt();
5052 Double_t area = bgrCluster->EffectiveAreaCharged();
5054 sumPtTotal += clusterPt;
5056 Int_t nTracksJet = bgrCluster->GetRefTracks()->GetEntries();
5058 for(Int_t it = 0; it<nTracksJet; it++){
5060 // embedded tracks - note: using ref tracks here, fBranchRecBckgClusters has to be consistent
5061 if( fUseExtraTracksBgr != 1){
5062 if(AliAODTrack* trackAOD = dynamic_cast<AliAODTrack*> (bgrCluster->GetTrack(it))){
5063 if(fUseExtraTracksBgr == 0 && ((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
5064 if(fUseExtraTracksBgr == -1 && !((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
5068 AliVParticle* track = dynamic_cast<AliVParticle*>(bgrCluster->GetTrack(it));
5069 if(!track) continue;
5071 Float_t trackPt = track->Pt();
5072 Float_t trackEta = track->Eta();
5073 Float_t trackPhi = TVector2::Phi_0_2pi(track->Phi());
5075 if(trackEta < fTrackEtaMin || trackEta > fTrackEtaMax) continue;
5076 if(trackPhi < fTrackPhiMin || trackPhi > fTrackPhiMax) continue;
5077 if(trackPt < fTrackPtCut) continue;
5079 outputlist->Add(track);
5083 Double_t areaJet = TMath::Pi()*rc*rc;
5084 if(areaTotal) normFactor = (Float_t) 1./(areaJet / areaTotal);
5089 // _______________________________________________________________________________________________________________________
5090 void AliAnalysisTaskIDFragmentationFunction::GetClusterTracksMedian(TList* outputlist, Double_t &normFactor)
5092 // fill tracks from bckgCluster branch,
5093 // using cluster with median density (odd number of clusters)
5094 // or picking randomly one of the two closest to median (even number)
5098 const Int_t nBckgClusters = fBckgJetsRec->GetEntries(); // not 'recCuts': use all clusters in full eta range
5100 if(nBckgClusters<3) return; // need at least 3 clusters (skipping 2 highest)
5102 Double_t* bgrDensity = new Double_t[nBckgClusters];
5103 Int_t* indices = new Int_t[nBckgClusters];
5105 for(Int_t ij=0; ij<nBckgClusters; ++ij){
5107 AliAODJet* bgrCluster = (AliAODJet*)(fBckgJetsRec->At(ij));
5108 Double_t clusterPt = bgrCluster->Pt();
5109 Double_t area = bgrCluster->EffectiveAreaCharged();
5111 Double_t density = 0;
5112 if(area>0) density = clusterPt/area;
5114 bgrDensity[ij] = density;
5118 TMath::Sort(nBckgClusters, bgrDensity, indices);
5120 // get median cluster
5122 AliAODJet* medianCluster = 0;
5123 //Double_t medianDensity = 0;
5125 if(TMath::Odd(nBckgClusters)){
5127 Int_t medianIndex = indices[(Int_t) (0.5*(nBckgClusters-1))];
5128 medianCluster = (AliAODJet*)(fBckgJetsRec->At(medianIndex));
5130 //Double_t clusterPt = medianCluster->Pt();
5131 //Double_t area = medianCluster->EffectiveAreaCharged();
5133 //if(area>0) medianDensity = clusterPt/area;
5137 Int_t medianIndex1 = indices[(Int_t) (0.5*nBckgClusters-1)];
5138 Int_t medianIndex2 = indices[(Int_t) (0.5*nBckgClusters)];
5140 AliAODJet* medianCluster1 = (AliAODJet*)(fBckgJetsRec->At(medianIndex1));
5141 AliAODJet* medianCluster2 = (AliAODJet*)(fBckgJetsRec->At(medianIndex2));
5143 //Double_t density1 = 0;
5144 //Double_t clusterPt1 = medianCluster1->Pt();
5145 //Double_t area1 = medianCluster1->EffectiveAreaCharged();
5146 //if(area1>0) density1 = clusterPt1/area1;
5148 //Double_t density2 = 0;
5149 //Double_t clusterPt2 = medianCluster2->Pt();
5150 //Double_t area2 = medianCluster2->EffectiveAreaCharged();
5151 //if(area2>0) density2 = clusterPt2/area2;
5153 //medianDensity = 0.5*(density1+density2);
5155 medianCluster = ( (fRandom->Rndm()>0.5) ? medianCluster1 : medianCluster2 ); // select one randomly to avoid adding areas
5158 Int_t nTracksJet = medianCluster->GetRefTracks()->GetEntries();
5160 for(Int_t it = 0; it<nTracksJet; it++){
5162 // embedded tracks - note: using ref tracks here, fBranchRecBckgClusters has to be consistent
5163 if( fUseExtraTracksBgr != 1){
5164 if(AliAODTrack* trackAOD = dynamic_cast<AliAODTrack*> (medianCluster->GetTrack(it))){
5165 if(fUseExtraTracksBgr == 0 && ((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
5166 if(fUseExtraTracksBgr == -1 && !((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
5170 AliVParticle* track = dynamic_cast<AliVParticle*>(medianCluster->GetTrack(it));
5171 if(!track) continue;
5173 Float_t trackPt = track->Pt();
5174 Float_t trackEta = track->Eta();
5175 Float_t trackPhi = TVector2::Phi_0_2pi(track->Phi());
5177 if(trackEta < fTrackEtaMin || trackEta > fTrackEtaMax) continue;
5178 if(trackPhi < fTrackPhiMin || trackPhi > fTrackPhiMax) continue;
5179 if(trackPt < fTrackPtCut) continue;
5181 outputlist->Add(track);
5184 Double_t areaMedian = medianCluster->EffectiveAreaCharged();
5185 Double_t areaJet = TMath::Pi()*GetFFRadius()*GetFFRadius();
5187 if(areaMedian) normFactor = (Float_t) 1./(areaJet / areaMedian);
5191 delete[] bgrDensity;
5195 // ______________________________________________________________________________________________________________________________________________________
5196 void AliAnalysisTaskIDFragmentationFunction::FillBckgHistos(Int_t type, TList* inputtracklist, TList* inputjetlist, AliAODJet* jet,
5197 AliFragFuncHistos* ffbckghistocuts, AliFragFuncQATrackHistos* qabckghistocuts, TH1F* fh1Mult){
5199 // List of tracks outside jets for background study
5200 TList* tracklistout2jets = new TList();
5201 TList* tracklistout3jets = new TList();
5202 TList* tracklistout2jetsStat = new TList();
5203 TList* tracklistout3jetsStat = new TList();
5204 Double_t sumPtOut2Jets = 0.;
5205 Double_t sumPtOut3Jets = 0.;
5206 Double_t sumPtOut2JetsStat = 0.;
5207 Double_t sumPtOut3JetsStat = 0.;
5208 Double_t normFactor2Jets = 0.;
5209 Double_t normFactor3Jets = 0.;
5211 Int_t nRecJetsCuts = inputjetlist->GetEntries();
5213 if(nRecJetsCuts>1) {
5214 GetTracksOutOfNJets(2,inputtracklist, tracklistout2jets, inputjetlist, sumPtOut2Jets);
5215 GetTracksOutOfNJetsStat(2,inputtracklist, tracklistout2jetsStat, inputjetlist,sumPtOut2JetsStat, normFactor2Jets);
5218 if(nRecJetsCuts>2) {
5219 GetTracksOutOfNJets(3,inputtracklist, tracklistout3jets, inputjetlist, sumPtOut3Jets);
5220 GetTracksOutOfNJetsStat(3,inputtracklist, tracklistout3jetsStat, inputjetlist, sumPtOut3JetsStat, normFactor3Jets);
5223 if(type==kBckgOutLJ || type==kBckgOutAJ)
5225 TList* tracklistoutleading = new TList();
5226 Double_t sumPtOutLeading = 0.;
5227 GetTracksOutOfNJets(1,inputtracklist, tracklistoutleading, inputjetlist, sumPtOutLeading);
5228 if(type==kBckgOutLJ && fh1Mult) fh1Mult->Fill(tracklistoutleading->GetSize());
5230 for(Int_t it=0; it<tracklistoutleading->GetSize(); ++it){
5232 AliVParticle* trackVP = (AliVParticle*)(tracklistoutleading->At(it));
5233 if(!trackVP) continue;
5234 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5236 Float_t jetPt = jet->Pt();
5237 Float_t trackPt = trackV->Pt();
5239 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5241 if(type==kBckgOutLJ)
5243 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt);
5245 // Fill track QA for background
5246 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt);
5249 // All cases included
5250 if(nRecJetsCuts==1 && type==kBckgOutAJ)
5252 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt );
5256 // Increment jet pt with one entry in case #tracks outside jets = 0
5257 if(tracklistoutleading->GetSize()==0) {
5258 Float_t jetPt = jet->Pt();
5259 Bool_t incrementJetPt = kTRUE;
5260 if(type==kBckgOutLJ)
5262 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
5264 // All cases included
5265 if(nRecJetsCuts==1 && type==kBckgOutAJ)
5267 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
5270 delete tracklistoutleading;
5272 if(type==kBckgOutLJStat || type==kBckgOutAJStat)
5274 TList* tracklistoutleadingStat = new TList();
5275 Double_t sumPtOutLeadingStat = 0.;
5276 Double_t normFactorLeading = 0.;
5278 GetTracksOutOfNJetsStat(1,inputtracklist, tracklistoutleadingStat, inputjetlist, sumPtOutLeadingStat, normFactorLeading);
5279 if(type==kBckgOutLJStat && fh1Mult) fh1Mult->Fill(tracklistoutleadingStat->GetSize());
5281 for(Int_t it=0; it<tracklistoutleadingStat->GetSize(); ++it){
5283 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistoutleadingStat->At(it));
5284 if(!trackVP) continue;
5285 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5287 Float_t jetPt = jet->Pt();
5288 Float_t trackPt = trackV->Pt();
5289 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5292 if(type==kBckgOutLJStat)
5294 if(fFFMode)ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactorLeading);
5296 // Fill track QA for background
5297 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt); // OB added bgr QA
5300 // All cases included
5301 if(nRecJetsCuts==1 && type==kBckgOutAJStat)
5303 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactorLeading);
5304 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt ); // OB added bgr QA
5309 // Increment jet pt with one entry in case #tracks outside jets = 0
5310 if(tracklistoutleadingStat->GetSize()==0) {
5311 Float_t jetPt = jet->Pt();
5312 Bool_t incrementJetPt = kTRUE;
5313 if(type==kBckgOutLJStat)
5315 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactorLeading);
5317 // All cases included
5318 if(nRecJetsCuts==1 && type==kBckgOutLJStat)
5320 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactorLeading);
5324 delete tracklistoutleadingStat;
5327 if(type==kBckgPerp || type==kBckgPerp2 || type==kBckgPerp2Area)
5329 Double_t sumPtPerp1 = 0.;
5330 Double_t sumPtPerp2 = 0.;
5331 TList* tracklistperp = new TList();
5332 TList* tracklistperp1 = new TList();
5333 TList* tracklistperp2 = new TList();
5336 if(type == kBckgPerp2) norm = 2; // in FillFF() scaleFac = 1/norm = 0.5 - account for double area;
5337 if(type == kBckgPerp2Area) norm = 2*TMath::Pi()*TMath::Abs(GetFFRadius())*TMath::Abs(GetFFRadius()) / jet->EffectiveAreaCharged(); // in FillFF() scaleFac = 1/norm;
5339 GetTracksTiltedwrpJetAxis(TMath::Pi()/2., inputtracklist,tracklistperp1,jet,TMath::Abs(GetFFRadius()),sumPtPerp1);
5340 if(type==kBckgPerp2 || type==kBckgPerp2Area) GetTracksTiltedwrpJetAxis(-1*TMath::Pi()/2., inputtracklist,tracklistperp2,jet,TMath::Abs(GetFFRadius()),sumPtPerp2);
5342 tracklistperp->AddAll(tracklistperp1);
5343 tracklistperp->AddAll(tracklistperp2);
5345 if(tracklistperp->GetSize() != tracklistperp1->GetSize() + tracklistperp2->GetSize()){
5346 cout<<" ERROR: tracklistperp size "<<tracklistperp->GetSize()<<" perp1 "<<tracklistperp1->GetSize()<<" perp2 "<<tracklistperp2->GetSize()<<endl;
5350 if(fh1Mult) fh1Mult->Fill(tracklistperp->GetSize());
5352 for(Int_t it=0; it<tracklistperp->GetSize(); ++it){
5354 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistperp->At(it));
5355 if(!trackVP)continue;
5356 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5358 Float_t jetPt = jet->Pt();
5359 Float_t trackPt = trackV->Pt();
5361 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5363 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, norm );
5365 // Fill track QA for background
5366 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt);
5370 // Increment jet pt with one entry in case #tracks outside jets = 0
5371 if(tracklistperp->GetSize()==0) {
5372 Float_t jetPt = jet->Pt();
5373 Bool_t incrementJetPt = kTRUE;
5374 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
5379 // fill for tracklistperp1/2 separately, divide norm by 2
5380 if(type==kBckgPerp){
5381 FillJetShape(jet, tracklistperp, fProNtracksLeadingJetBgrPerp2, fProDelRPtSumBgrPerp2, 0, TMath::Pi()/2., 0., kFALSE);
5383 if(type==kBckgPerp2){
5384 FillJetShape(jet, tracklistperp1, fProNtracksLeadingJetBgrPerp2, fProDelRPtSumBgrPerp2, 0, TMath::Pi()/2., 0., kFALSE);
5385 FillJetShape(jet, tracklistperp2, fProNtracksLeadingJetBgrPerp2, fProDelRPtSumBgrPerp2, 0, -1*TMath::Pi()/2., 0., kFALSE);
5387 if(type==kBckgPerp2Area){ // divide norm by 2: listperp1/2 filled separately
5388 FillJetShape(jet, tracklistperp1, fProNtracksLeadingJetBgrPerp2, fProDelRPtSumBgrPerp2, 0, TMath::Pi()/2., 0.5*norm, kFALSE);
5389 FillJetShape(jet, tracklistperp2, fProNtracksLeadingJetBgrPerp2, fProDelRPtSumBgrPerp2, 0, -1*TMath::Pi()/2., 0.5*norm, kFALSE);
5393 delete tracklistperp;
5394 delete tracklistperp1;
5395 delete tracklistperp2;
5398 if(type==kBckgASide)
5400 Double_t sumPtASide = 0.;
5401 TList* tracklistaside = new TList();
5402 GetTracksTiltedwrpJetAxis(TMath::Pi(),inputtracklist,tracklistaside,jet,TMath::Abs(GetFFRadius()),sumPtASide);
5403 if(fh1Mult) fh1Mult->Fill(tracklistaside->GetSize());
5405 for(Int_t it=0; it<tracklistaside->GetSize(); ++it){
5407 AliVParticle* trackVP = (AliVParticle*)(tracklistaside->At(it));
5408 if(!trackVP) continue;
5409 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5411 Float_t jetPt = jet->Pt();
5412 Float_t trackPt = trackV->Pt();
5414 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5416 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt );
5418 // Fill track QA for background
5419 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt);
5423 if(tracklistaside->GetSize()==0) {
5424 Float_t jetPt = jet->Pt();
5425 Bool_t incrementJetPt = kTRUE;
5426 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
5429 delete tracklistaside;
5432 if(type==kBckgASideWindow)
5434 Double_t normFactorASide = 0.;
5435 Double_t sumPtASideW = 0.;
5436 TList* tracklistasidew = new TList();
5437 GetTracksTiltedwrpJetAxisWindow(TMath::Pi(),inputtracklist,tracklistasidew,jet,TMath::Abs(GetFFRadius()),sumPtASideW,normFactorASide);
5438 if(fh1Mult) fh1Mult->Fill(tracklistasidew->GetSize());
5440 for(Int_t it=0; it<tracklistasidew->GetSize(); ++it){
5442 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistasidew->At(it));
5443 if(!trackVP) continue;
5444 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5446 Float_t jetPt = jet->Pt();
5447 Float_t trackPt = trackV->Pt();
5448 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5450 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactorASide);
5452 // Fill track QA for background
5453 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt, kFALSE, normFactorASide);
5457 if(tracklistasidew->GetSize()==0) {
5458 Float_t jetPt = jet->Pt();
5459 Bool_t incrementJetPt = kTRUE;
5460 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactorASide);
5463 delete tracklistasidew;
5466 if(type==kBckgPerpWindow)
5468 Double_t normFactorPerp = 0.;
5469 Double_t sumPtPerpW = 0.;
5470 TList* tracklistperpw = new TList();
5471 GetTracksTiltedwrpJetAxisWindow(TMath::Pi()/2.,inputtracklist,tracklistperpw,jet,TMath::Abs(GetFFRadius()),sumPtPerpW,normFactorPerp);
5472 if(fh1Mult) fh1Mult->Fill(tracklistperpw->GetSize());
5474 for(Int_t it=0; it<tracklistperpw->GetSize(); ++it){
5476 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistperpw->At(it));
5477 if(!trackVP) continue;
5478 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5480 Float_t jetPt = jet->Pt();
5481 Float_t trackPt = trackV->Pt();
5482 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5484 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactorPerp);
5486 // Fill track QA for background
5487 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt, kFALSE, normFactorPerp);
5491 if(tracklistperpw->GetSize()==0) {
5492 Float_t jetPt = jet->Pt();
5493 Bool_t incrementJetPt = kTRUE;
5494 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactorPerp);
5497 delete tracklistperpw;
5501 if(type==kBckgOut2J || type==kBckgOutAJ)
5503 if(type==kBckgOut2J && fh1Mult) fh1Mult->Fill(tracklistout2jets->GetSize());
5504 for(Int_t it=0; it<tracklistout2jets->GetSize(); ++it){
5506 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistout2jets->At(it));
5507 if(!trackVP) continue;
5508 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5510 Float_t jetPt = jet->Pt();
5511 Float_t trackPt = trackV->Pt();
5513 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5515 if(type==kBckgOut2J)
5517 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt );
5518 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt);
5521 // All cases included
5522 if(nRecJetsCuts==2 && type==kBckgOutAJ)
5524 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt );
5529 // Increment jet pt with one entry in case #tracks outside jets = 0
5530 if(tracklistout2jets->GetSize()==0) {
5531 Float_t jetPt = jet->Pt();
5532 Bool_t incrementJetPt = kTRUE;
5533 if(type==kBckgOut2J)
5535 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
5537 // All cases included
5538 if(nRecJetsCuts==2 && type==kBckgOutAJ)
5540 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
5545 if(type==kBckgOut2JStat || type==kBckgOutAJStat)
5547 for(Int_t it=0; it<tracklistout2jetsStat->GetSize(); ++it){
5549 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistout2jetsStat->At(it));
5550 if(!trackVP) continue;
5551 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5553 Float_t jetPt = jet->Pt();
5554 Float_t trackPt = trackV->Pt();
5555 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5557 if(type==kBckgOut2JStat)
5559 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactor2Jets);
5561 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt ); // OB added bgr QA
5564 // All cases included
5565 if(nRecJetsCuts==2 && type==kBckgOutAJStat)
5567 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactor2Jets);
5569 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt ); // OB added bgr QA
5573 // Increment jet pt with one entry in case #tracks outside jets = 0
5574 if(tracklistout2jetsStat->GetSize()==0) {
5575 Float_t jetPt = jet->Pt();
5576 Bool_t incrementJetPt = kTRUE;
5577 if(type==kBckgOut2JStat)
5579 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactor2Jets);
5581 // All cases included
5582 if(nRecJetsCuts==2 && type==kBckgOutAJStat)
5584 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactor2Jets);
5590 if(type==kBckgOut3J || type==kBckgOutAJ)
5592 if(type==kBckgOut3J && fh1Mult) fh1Mult->Fill(tracklistout3jets->GetSize());
5594 for(Int_t it=0; it<tracklistout3jets->GetSize(); ++it){
5596 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistout3jets->At(it));
5597 if(!trackVP) continue;
5598 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5600 Float_t jetPt = jet->Pt();
5601 Float_t trackPt = trackV->Pt();
5603 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5605 if(type==kBckgOut3J)
5607 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt );
5609 qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt);
5612 // All cases included
5613 if(nRecJetsCuts==3 && type==kBckgOutAJ)
5615 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt );
5620 // Increment jet pt with one entry in case #tracks outside jets = 0
5621 if(tracklistout3jets->GetSize()==0) {
5622 Float_t jetPt = jet->Pt();
5623 Bool_t incrementJetPt = kTRUE;
5624 if(type==kBckgOut3J)
5626 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
5628 // All cases included
5629 if(nRecJetsCuts==3 && type==kBckgOutAJ)
5631 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
5636 if(type==kBckgOut3JStat || type==kBckgOutAJStat)
5638 for(Int_t it=0; it<tracklistout3jetsStat->GetSize(); ++it){
5640 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistout3jetsStat->At(it));
5641 if(!trackVP) continue;
5642 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5644 Float_t jetPt = jet->Pt();
5645 Float_t trackPt = trackV->Pt();
5646 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5648 if(type==kBckgOut3JStat)
5650 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactor3Jets);
5652 //if(fQAMode&1) qabckghistocuts->FillTrackQA( trackEta, TVector2::Phi_0_2pi(trackPhi), trackPt);
5655 // All cases included
5656 if(nRecJetsCuts==3 && type==kBckgOutAJStat)
5658 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactor3Jets );
5660 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt );
5665 // Increment jet pt with one entry in case #tracks outside jets = 0
5666 if(tracklistout3jetsStat->GetSize()==0) {
5667 Float_t jetPt = jet->Pt();
5668 Bool_t incrementJetPt = kTRUE;
5669 if(type==kBckgOut3JStat)
5671 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactor3Jets);
5673 // All cases included
5674 if(nRecJetsCuts==3 && type==kBckgOutAJStat)
5676 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactor3Jets);
5682 if(type==kBckgClustersOutLeading){ // clusters bgr: all tracks in clusters out of leading jet
5684 TList* tracklistClustersOutLeading = new TList();
5685 Double_t normFactorClusters = 0;
5686 Float_t jetPt = jet->Pt();
5688 GetClusterTracksOutOf1Jet(jet, tracklistClustersOutLeading, normFactorClusters);
5689 if(fh1Mult) fh1Mult->Fill(tracklistClustersOutLeading->GetSize());
5691 for(Int_t it=0; it<tracklistClustersOutLeading->GetSize(); ++it){
5693 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistClustersOutLeading->At(it));
5694 if(!trackVP) continue;
5695 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5697 Float_t trackPt = trackVP->Pt();
5699 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5701 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactorClusters );
5702 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt );
5707 delete tracklistClustersOutLeading;
5711 if(type == kBckgClusters){ // clusters bgr: all tracks in 'median cluster'
5713 TList* tracklistClustersMedian = new TList();
5714 Double_t normFactorClusters = 0;
5715 Float_t jetPt = jet->Pt();
5717 GetClusterTracksMedian(tracklistClustersMedian, normFactorClusters);
5718 if(fh1Mult) fh1Mult->Fill(tracklistClustersMedian->GetSize());
5720 for(Int_t it=0; it<tracklistClustersMedian->GetSize(); ++it){
5722 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistClustersMedian->At(it));
5723 if(!trackVP) continue;
5724 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5726 Float_t trackPt = trackVP->Pt();
5728 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5730 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactorClusters );
5731 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt );
5736 delete tracklistClustersMedian;
5739 delete tracklistout2jets;
5740 delete tracklistout3jets;
5741 delete tracklistout2jetsStat;
5742 delete tracklistout3jetsStat;
5745 //_____________________________________________________________________________________
5746 Double_t AliAnalysisTaskIDFragmentationFunction::GetMCStrangenessFactor(Double_t pt) const
5748 // factor strangeness data/MC as function of pt from UE analysis (Sara Vallero)
5752 if(0.150<pt && pt<0.200) alpha = 3.639;
5753 if(0.200<pt && pt<0.250) alpha = 2.097;
5754 if(0.250<pt && pt<0.300) alpha = 1.930;
5755 if(0.300<pt && pt<0.350) alpha = 1.932;
5756 if(0.350<pt && pt<0.400) alpha = 1.943;
5757 if(0.400<pt && pt<0.450) alpha = 1.993;
5758 if(0.450<pt && pt<0.500) alpha = 1.989;
5759 if(0.500<pt && pt<0.600) alpha = 1.963;
5760 if(0.600<pt && pt<0.700) alpha = 1.917;
5761 if(0.700<pt && pt<0.800) alpha = 1.861;
5762 if(0.800<pt && pt<0.900) alpha = 1.820;
5763 if(0.900<pt && pt<1.000) alpha = 1.741;
5764 if(1.000<pt && pt<1.500) alpha = 0.878;
5769 //__________________________________________________________________________________________________
5770 Double_t AliAnalysisTaskIDFragmentationFunction::GetMCStrangenessFactorCMS(AliAODMCParticle* daughter) const
5772 // strangeness ratio MC/data as function of mother pt from CMS data in |eta|<2.0
5774 TClonesArray *tca = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
5777 AliAODMCParticle* currentMother = daughter;
5778 AliAODMCParticle* currentDaughter = daughter;
5781 // find first primary mother K0s, Lambda or Xi
5784 Int_t daughterPDG = currentDaughter->GetPdgCode();
5786 Int_t motherLabel = currentDaughter->GetMother();
5787 if(motherLabel >= tca->GetEntriesFast()){ // protection
5788 currentMother = currentDaughter;
5792 currentMother = (AliAODMCParticle*) tca->At(motherLabel);
5795 currentMother = currentDaughter;
5799 Int_t motherPDG = currentMother->GetPdgCode();
5801 // phys. primary found ?
5802 if(currentMother->IsPhysicalPrimary()) break;
5804 if(TMath::Abs(daughterPDG) == 321){ // K+/K- e.g. from phi (ref data not feeddown corrected)
5805 currentMother = currentDaughter; break;
5807 if(TMath::Abs(motherPDG) == 310 ){ // K0s e.g. from phi (ref data not feeddown corrected)
5810 if(TMath::Abs(motherPDG) == 3212 && TMath::Abs(daughterPDG) == 3122){ // mother Sigma0, daughter Lambda (this case not included in feeddown corr.)
5811 currentMother = currentDaughter; break;
5814 currentDaughter = currentMother;
5818 Int_t motherPDG = currentMother->GetPdgCode();
5819 Double_t motherGenPt = currentMother->Pt();
5821 return AliAnalysisTaskPID::GetMCStrangenessFactorCMS(motherPDG, motherGenPt);
5824 // _________________________________________________________________________________
5825 void AliAnalysisTaskIDFragmentationFunction::FillJetShape(AliAODJet* jet, TList* list,
5826 TProfile* hProNtracksLeadingJet, TProfile** hProDelRPtSum, TProfile* hProDelR80pcPt,
5827 Double_t dPhiUE, Double_t normUE, Bool_t scaleStrangeness)
5829 // Fill jet shape histos
5831 const Int_t kNbinsR = 50;
5832 const Float_t kBinWidthR = 0.02;
5834 Int_t nJetTracks = list->GetEntries();
5836 Float_t pTsumA[kNbinsR] = {0.0};
5838 Float_t *delRA = new Float_t[nJetTracks];
5839 Float_t *trackPtA = new Float_t[nJetTracks];
5840 Int_t *index = new Int_t[nJetTracks];
5842 for(Int_t i=0; i<nJetTracks; i++){
5849 jet->PxPyPz(jetMom);
5850 TVector3 jet3mom(jetMom);
5852 if(TMath::Abs(dPhiUE)>0){
5853 Double_t phiTilted = jet3mom.Phi();
5854 phiTilted += dPhiUE;
5855 phiTilted = TVector2::Phi_0_2pi(phiTilted);
5856 jet3mom.SetPhi(phiTilted);
5859 Double_t jetPt = jet->Pt();
5860 Double_t sumWeights = 0;
5862 for (Int_t j =0; j<nJetTracks; j++){
5864 AliVParticle* track = dynamic_cast<AliVParticle*>(list->At(j));
5867 Double_t trackMom[3];
5868 track->PxPyPz(trackMom);
5869 TVector3 track3mom(trackMom);
5871 Double_t dR = jet3mom.DeltaR(track3mom);
5874 trackPtA[j] = track->Pt();
5876 Double_t weight = GetMCStrangenessFactor(track->Pt()); // more correctly should be gen pt
5877 sumWeights += weight;
5879 for(Int_t ibin=1; ibin<=kNbinsR; ibin++){
5880 Float_t xlow = kBinWidthR*(ibin-1);
5881 Float_t xup = kBinWidthR*ibin;
5882 if(xlow <= dR && dR < xup){
5884 if(scaleStrangeness) pTsumA[ibin-1] += track->Pt()*weight;
5885 else pTsumA[ibin-1] += track->Pt();
5893 for(Int_t ibin=0; ibin<kNbinsR; ibin++){
5894 Float_t fR = kBinWidthR*(ibin+0.5);
5896 for(Int_t k=0; k<5; k++){
5897 if(k==0){jetPtMin=20.0;jetPtMax=30.0;}
5898 if(k==1){jetPtMin=30.0;jetPtMax=40.0;}
5899 if(k==2){jetPtMin=40.0;jetPtMax=60.0;}
5900 if(k==3){jetPtMin=60.0;jetPtMax=80.0;}
5901 if(k==4){jetPtMin=80.0;jetPtMax=100.0;}
5902 if(jetPt>jetPtMin && jetPt<jetPtMax){
5904 hProDelRPtSum[k]->Fill(fR,pTsumA[ibin]);
5910 if(scaleStrangeness) hProNtracksLeadingJet->Fill(jetPt,sumWeights);
5911 else hProNtracksLeadingJet->Fill(jetPt,nJetTracks);
5913 if(normUE) hProNtracksLeadingJet->Fill(jetPt,nJetTracks/normUE);
5918 Float_t delRPtSum80pc = 0;
5920 TMath::Sort(nJetTracks,delRA,index,0);
5922 for(Int_t ii=0; ii<nJetTracks; ii++){
5924 if(scaleStrangeness){
5925 Double_t weight = GetMCStrangenessFactor(trackPtA[index[ii]]); // more correctly should be gen pt
5926 pTsum += weight*trackPtA[index[ii]];
5928 else pTsum += trackPtA[index[ii]];
5931 if(pTsum/jetPt >= 0.8000){
5932 delRPtSum80pc = delRA[index[ii]];
5936 hProDelR80pcPt->Fill(jetPt,delRPtSum80pc);
5945 // _________________________________________________________________________________
5946 Bool_t AliAnalysisTaskIDFragmentationFunction::IsSecondaryWithStrangeMotherMC(AliAODMCParticle* part)
5948 // Check whether particle is a secondary with strange mother, i.e. returns kTRUE if a strange mother is found
5949 // and the particle is NOT a physical primary. In all other cases kFALSE is returned
5951 TClonesArray *tca = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
5955 if (part->IsPhysicalPrimary())
5958 Int_t iMother = part->GetMother();
5963 AliAODMCParticle* partM = dynamic_cast<AliAODMCParticle*>(tca->At(iMother));
5967 Int_t codeM = TMath::Abs(partM->GetPdgCode());
5968 Int_t mfl = Int_t(codeM / TMath::Power(10, Int_t(TMath::Log10(codeM))));
5969 if (mfl == 3 && codeM != 3) // codeM = 3 is for s quark