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)
188 ,fh1nRecBckgJetsCuts(0)
190 ,fh2PtRecVsGenPrim(0)
194 ,fhJetPtRefMultEta5(0)
195 ,fhJetPtRefMultEta8(0)
196 ,fhJetPtMultPercent(0)
197 ,fQATrackHistosRecEffGen(0)
198 ,fQATrackHistosRecEffRec(0)
199 ,fQATrackHistosSecRecNS(0)
200 ,fQATrackHistosSecRecS(0)
201 ,fQATrackHistosSecRecSsc(0)
202 ,fFFHistosRecEffRec(0)
203 ,fFFHistosSecRecNS(0)
205 ,fFFHistosSecRecSsc(0)
212 ,fh1FractionPtEmbedded(0)
214 ,fh2DeltaPtVsJetPtEmbedded(0)
215 ,fh2DeltaPtVsRecJetPtEmbedded(0)
216 ,fh1DeltaREmbedded(0)
217 ,fQABckgHisto0RecCuts(0)
219 ,fQABckgHisto1RecCuts(0)
221 ,fQABckgHisto2RecCuts(0)
223 ,fQABckgHisto3RecCuts(0)
225 ,fQABckgHisto4RecCuts(0)
227 ,fFFBckgHisto0RecCuts(0)
229 ,fFFBckgHisto1RecCuts(0)
231 ,fFFBckgHisto2RecCuts(0)
233 ,fFFBckgHisto3RecCuts(0)
235 ,fFFBckgHisto4RecCuts(0)
237 ,fFFBckgHisto0RecEffRec(0)
238 ,fFFBckgHisto0SecRecNS(0)
239 ,fFFBckgHisto0SecRecS(0)
240 ,fFFBckgHisto0SecRecSsc(0)
242 ,fProNtracksLeadingJet(0)
244 ,fProNtracksLeadingJetGen(0)
245 ,fProDelR80pcPtGen(0)
246 ,fProNtracksLeadingJetBgrPerp2(0)
247 ,fProNtracksLeadingJetRecPrim(0)
248 ,fProDelR80pcPtRecPrim(0)
249 ,fProNtracksLeadingJetRecSecNS(0)
250 ,fProNtracksLeadingJetRecSecS(0)
251 ,fProNtracksLeadingJetRecSecSsc(0)
255 ,fOnlyLeadingJets(kFALSE)
261 ,fNumInclusivePIDtasks(0)
263 ,fNameInclusivePIDtask(0x0)
264 ,fNameJetPIDtask(0x0)
265 ,fInclusivePIDtask(0x0)
267 ,fUseInclusivePIDtask(kFALSE)
268 ,fUseJetPIDtask(kFALSE)
271 // default constructor
278 for(Int_t ii=0; ii<5; ii++){
279 fProDelRPtSum[ii] = 0;
280 fProDelRPtSumGen[ii] = 0;
281 fProDelRPtSumBgrPerp2[ii] = 0;
282 fProDelRPtSumRecPrim[ii] = 0;
283 fProDelRPtSumRecSecNS[ii] = 0;
284 fProDelRPtSumRecSecS[ii] = 0;
285 fProDelRPtSumRecSecSsc[ii] = 0;
288 for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
289 fIDFFHistosRecCuts[i] = 0x0;
290 fIDFFHistosGen[i] = 0x0;
292 fhDCA_XY_prim_MCID[i] = 0x0;
293 fhDCA_Z_prim_MCID[i] = 0x0;
295 fhDCA_XY_sec_MCID[i] = 0x0;
296 fhDCA_Z_sec_MCID[i] = 0x0;
300 //_______________________________________________________________________________________________
301 AliAnalysisTaskIDFragmentationFunction::AliAnalysisTaskIDFragmentationFunction(const char *name)
302 : AliAnalysisTaskSE(name)
308 ,fBranchRecJets("jets")
309 ,fBranchRecBckgClusters("")
311 ,fBranchEmbeddedJets("")
315 ,fUseAODInputJets(kTRUE)
317 ,fUsePhysicsSelection(kTRUE)
318 ,fEvtSelectionMask(0)
327 ,fUseExtraTracksBgr(0)
328 ,fCutFractionPtEmbedded(0)
329 ,fUseEmbeddedJetAxis(0)
330 ,fUseEmbeddedJetPt(0)
349 ,fTracksRecCutsEfficiency(0)
351 ,fTracksAODMCCharged(0)
352 ,fTracksAODMCChargedSecNS(0)
353 ,fTracksAODMCChargedSecS(0)
354 ,fTracksRecQualityCuts(0)
363 ,fQATrackHistosRecCuts(0)
364 ,fQATrackHistosGen(0)
366 ,fQAJetHistosRecCuts(0)
367 ,fQAJetHistosRecCutsLeading(0)
369 ,fQAJetHistosGenLeading(0)
370 ,fQAJetHistosRecEffLeading(0)
372 ,fFFHistosRecCutsInc(0)
373 ,fFFHistosRecLeadingTrack(0)
376 ,fFFHistosGenLeadingTrack(0)
377 ,fQATrackHighPtThreshold(0)
411 ,fh1VertexNContributors(0)
424 ,fh1nRecBckgJetsCuts(0)
426 ,fh2PtRecVsGenPrim(0)
430 ,fhJetPtRefMultEta5(0)
431 ,fhJetPtRefMultEta8(0)
432 ,fhJetPtMultPercent(0)
433 ,fQATrackHistosRecEffGen(0)
434 ,fQATrackHistosRecEffRec(0)
435 ,fQATrackHistosSecRecNS(0)
436 ,fQATrackHistosSecRecS(0)
437 ,fQATrackHistosSecRecSsc(0)
438 ,fFFHistosRecEffRec(0)
439 ,fFFHistosSecRecNS(0)
441 ,fFFHistosSecRecSsc(0)
448 ,fh1FractionPtEmbedded(0)
450 ,fh2DeltaPtVsJetPtEmbedded(0)
451 ,fh2DeltaPtVsRecJetPtEmbedded(0)
452 ,fh1DeltaREmbedded(0)
453 ,fQABckgHisto0RecCuts(0)
455 ,fQABckgHisto1RecCuts(0)
457 ,fQABckgHisto2RecCuts(0)
459 ,fQABckgHisto3RecCuts(0)
461 ,fQABckgHisto4RecCuts(0)
463 ,fFFBckgHisto0RecCuts(0)
465 ,fFFBckgHisto1RecCuts(0)
467 ,fFFBckgHisto2RecCuts(0)
469 ,fFFBckgHisto3RecCuts(0)
471 ,fFFBckgHisto4RecCuts(0)
473 ,fFFBckgHisto0RecEffRec(0)
474 ,fFFBckgHisto0SecRecNS(0)
475 ,fFFBckgHisto0SecRecS(0)
476 ,fFFBckgHisto0SecRecSsc(0)
478 ,fProNtracksLeadingJet(0)
480 ,fProNtracksLeadingJetGen(0)
481 ,fProDelR80pcPtGen(0)
482 ,fProNtracksLeadingJetBgrPerp2(0)
483 ,fProNtracksLeadingJetRecPrim(0)
484 ,fProDelR80pcPtRecPrim(0)
485 ,fProNtracksLeadingJetRecSecNS(0)
486 ,fProNtracksLeadingJetRecSecS(0)
487 ,fProNtracksLeadingJetRecSecSsc(0)
489 ,fOnlyLeadingJets(kFALSE)
493 ,fNumInclusivePIDtasks(0)
495 ,fNameInclusivePIDtask(0x0)
496 ,fNameJetPIDtask(0x0)
497 ,fInclusivePIDtask(0x0)
499 ,fUseInclusivePIDtask(kFALSE)
500 ,fUseJetPIDtask(kFALSE)
510 for(Int_t ii=0; ii<5; ii++){
511 fProDelRPtSum[ii] = 0;
512 fProDelRPtSumGen[ii] = 0;
513 fProDelRPtSumBgrPerp2[ii] = 0;
514 fProDelRPtSumRecPrim[ii] = 0;
515 fProDelRPtSumRecSecNS[ii] = 0;
516 fProDelRPtSumRecSecS[ii] = 0;
517 fProDelRPtSumRecSecSsc[ii] = 0;
520 for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
521 fIDFFHistosRecCuts[i] = 0x0;
522 fIDFFHistosGen[i] = 0x0;
524 fhDCA_XY_prim_MCID[i] = 0x0;
525 fhDCA_Z_prim_MCID[i] = 0x0;
527 fhDCA_XY_sec_MCID[i] = 0x0;
528 fhDCA_Z_sec_MCID[i] = 0x0;
531 DefineOutput(1,TList::Class());
534 //__________________________________________________________________________________________________________________________
535 AliAnalysisTaskIDFragmentationFunction::AliAnalysisTaskIDFragmentationFunction(const AliAnalysisTaskIDFragmentationFunction ©)
536 : AliAnalysisTaskSE()
539 ,fAODJets(copy.fAODJets)
540 ,fAODExtension(copy.fAODExtension)
541 ,fNonStdFile(copy.fNonStdFile)
542 ,fBranchRecJets(copy.fBranchRecJets)
543 ,fBranchRecBckgClusters(copy.fBranchRecBckgClusters)
544 ,fBranchGenJets(copy.fBranchGenJets)
545 ,fBranchEmbeddedJets(copy.fBranchEmbeddedJets)
546 ,fTrackTypeGen(copy.fTrackTypeGen)
547 ,fJetTypeGen(copy.fJetTypeGen)
548 ,fJetTypeRecEff(copy.fJetTypeRecEff)
549 ,fUseAODInputJets(copy.fUseAODInputJets)
550 ,fFilterMask(copy.fFilterMask)
551 ,fUsePhysicsSelection(copy.fUsePhysicsSelection)
552 ,fEvtSelectionMask(copy.fEvtSelectionMask)
553 ,fEventClass(copy.fEventClass)
554 ,fMaxVertexZ(copy.fMaxVertexZ)
555 ,fTrackPtCut(copy.fTrackPtCut)
556 ,fTrackEtaMin(copy.fTrackEtaMin)
557 ,fTrackEtaMax(copy.fTrackEtaMax)
558 ,fTrackPhiMin(copy.fTrackPhiMin)
559 ,fTrackPhiMax(copy.fTrackPhiMax)
560 ,fUseExtraTracks(copy.fUseExtraTracks)
561 ,fUseExtraTracksBgr(copy.fUseExtraTracksBgr)
562 ,fCutFractionPtEmbedded(copy.fCutFractionPtEmbedded)
563 ,fUseEmbeddedJetAxis(copy.fUseEmbeddedJetAxis)
564 ,fUseEmbeddedJetPt(copy.fUseEmbeddedJetPt)
565 ,fJetPtCut(copy.fJetPtCut)
566 ,fJetEtaMin(copy.fJetEtaMin)
567 ,fJetEtaMax(copy.fJetEtaMax)
568 ,fJetPhiMin(copy.fJetPhiMin)
569 ,fJetPhiMax(copy.fJetPhiMax)
570 ,fFFRadius(copy.fFFRadius)
571 ,fFFMinLTrackPt(copy.fFFMinLTrackPt)
572 ,fFFMaxTrackPt(copy.fFFMaxTrackPt)
573 ,fFFMinnTracks(copy.fFFMinnTracks)
574 ,fFFBckgRadius(copy.fFFBckgRadius)
575 ,fBckgMode(copy.fBckgMode)
576 ,fQAMode(copy.fQAMode)
577 ,fFFMode(copy.fFFMode)
578 ,fIDFFMode(copy.fIDFFMode)
579 ,fEffMode(copy.fEffMode)
580 ,fJSMode(copy.fJSMode)
581 ,fAvgTrials(copy.fAvgTrials)
582 ,fTracksRecCuts(copy.fTracksRecCuts)
583 ,fTracksRecCutsEfficiency(copy.fTracksRecCutsEfficiency)
584 ,fTracksGen(copy.fTracksGen)
585 ,fTracksAODMCCharged(copy.fTracksAODMCCharged)
586 ,fTracksAODMCChargedSecNS(copy.fTracksAODMCChargedSecNS)
587 ,fTracksAODMCChargedSecS(copy.fTracksAODMCChargedSecS)
588 ,fTracksRecQualityCuts(copy.fTracksRecQualityCuts)
589 ,fJetsRec(copy.fJetsRec)
590 ,fJetsRecCuts(copy.fJetsRecCuts)
591 ,fJetsGen(copy.fJetsGen)
592 ,fJetsRecEff(copy.fJetsRecEff)
593 ,fJetsEmbedded(copy.fJetsEmbedded)
594 ,fBckgJetsRec(copy.fBckgJetsRec)
595 ,fBckgJetsRecCuts(copy.fBckgJetsRecCuts)
596 ,fBckgJetsGen(copy.fBckgJetsGen)
597 ,fQATrackHistosRecCuts(copy.fQATrackHistosRecCuts)
598 ,fQATrackHistosGen(copy.fQATrackHistosGen)
599 ,fQAJetHistosRec(copy.fQAJetHistosRec)
600 ,fQAJetHistosRecCuts(copy.fQAJetHistosRecCuts)
601 ,fQAJetHistosRecCutsLeading(copy.fQAJetHistosRecCutsLeading)
602 ,fQAJetHistosGen(copy.fQAJetHistosGen)
603 ,fQAJetHistosGenLeading(copy.fQAJetHistosGenLeading)
604 ,fQAJetHistosRecEffLeading(copy.fQAJetHistosRecEffLeading)
605 ,fFFHistosRecCuts(copy.fFFHistosRecCuts)
606 ,fFFHistosRecCutsInc(copy.fFFHistosRecCutsInc)
607 ,fFFHistosRecLeadingTrack(copy.fFFHistosRecLeadingTrack)
608 ,fFFHistosGen(copy.fFFHistosGen)
609 ,fFFHistosGenInc(copy.fFFHistosGenInc)
610 ,fFFHistosGenLeadingTrack(copy.fFFHistosGenLeadingTrack)
611 ,fQATrackHighPtThreshold(copy.fQATrackHighPtThreshold)
612 ,fFFNBinsJetPt(copy.fFFNBinsJetPt)
613 ,fFFJetPtMin(copy.fFFJetPtMin)
614 ,fFFJetPtMax(copy.fFFJetPtMax)
615 ,fFFNBinsPt(copy.fFFNBinsPt)
616 ,fFFPtMin(copy.fFFPtMin)
617 ,fFFPtMax(copy.fFFPtMax)
618 ,fFFNBinsXi(copy.fFFNBinsXi)
619 ,fFFXiMin(copy.fFFXiMin)
620 ,fFFXiMax(copy.fFFXiMax)
621 ,fFFNBinsZ(copy.fFFNBinsZ)
622 ,fFFZMin(copy.fFFZMin)
623 ,fFFZMax(copy.fFFZMax)
624 ,fQAJetNBinsPt(copy.fQAJetNBinsPt)
625 ,fQAJetPtMin(copy.fQAJetPtMin)
626 ,fQAJetPtMax(copy.fQAJetPtMax)
627 ,fQAJetNBinsEta(copy.fQAJetNBinsEta)
628 ,fQAJetEtaMin(copy.fQAJetEtaMin)
629 ,fQAJetEtaMax(copy.fQAJetEtaMax)
630 ,fQAJetNBinsPhi(copy.fQAJetNBinsPhi)
631 ,fQAJetPhiMin(copy.fQAJetPhiMin)
632 ,fQAJetPhiMax(copy.fQAJetPhiMax)
633 ,fQATrackNBinsPt(copy.fQATrackNBinsPt)
634 ,fQATrackPtMin(copy.fQATrackPtMin)
635 ,fQATrackPtMax(copy.fQATrackPtMax)
636 ,fQATrackNBinsEta(copy.fQATrackNBinsEta)
637 ,fQATrackEtaMin(copy.fQATrackEtaMin)
638 ,fQATrackEtaMax(copy.fQATrackEtaMax)
639 ,fQATrackNBinsPhi(copy.fQATrackNBinsPhi)
640 ,fQATrackPhiMin(copy.fQATrackPhiMin)
641 ,fQATrackPhiMax(copy.fQATrackPhiMax)
642 ,fCommonHistList(copy.fCommonHistList)
643 ,fh1EvtSelection(copy.fh1EvtSelection)
644 ,fh1VtxSelection(copy.fh1VtxSelection)
645 ,fh1VertexNContributors(copy.fh1VertexNContributors)
646 ,fh1VertexZ(copy.fh1VertexZ)
647 ,fh1EvtMult(copy.fh1EvtMult)
648 ,fh1EvtCent(copy.fh1EvtCent)
649 ,fh1Xsec(copy.fh1Xsec)
650 ,fh1Trials(copy.fh1Trials)
651 ,fh1PtHard(copy.fh1PtHard)
652 ,fh1PtHardTrials(copy.fh1PtHardTrials)
653 ,fh1EvtsPtHardCut(copy.fh1EvtsPtHardCut)
654 ,fh1nRecJetsCuts(copy.fh1nRecJetsCuts)
655 ,fh1nGenJets(copy.fh1nGenJets)
656 ,fh1nRecEffJets(copy.fh1nRecEffJets)
657 ,fh1nEmbeddedJets(copy.fh1nEmbeddedJets)
658 ,fh1nRecBckgJetsCuts(copy.fh1nRecBckgJetsCuts)
659 ,fh1nGenBckgJets(copy.fh1nGenBckgJets)
660 ,fh2PtRecVsGenPrim(copy.fh2PtRecVsGenPrim)
661 ,fh2PtRecVsGenSec(copy.fh2PtRecVsGenSec)
662 ,fhDCA_XY(copy.fhDCA_XY)
663 ,fhDCA_Z(copy.fhDCA_Z)
664 ,fhJetPtRefMultEta5(copy.fhJetPtRefMultEta5)
665 ,fhJetPtRefMultEta8(copy.fhJetPtRefMultEta8)
666 ,fhJetPtMultPercent(copy.fhJetPtMultPercent)
667 ,fQATrackHistosRecEffGen(copy.fQATrackHistosRecEffGen)
668 ,fQATrackHistosRecEffRec(copy.fQATrackHistosRecEffRec)
669 ,fQATrackHistosSecRecNS(copy.fQATrackHistosSecRecNS)
670 ,fQATrackHistosSecRecS(copy.fQATrackHistosSecRecS)
671 ,fQATrackHistosSecRecSsc(copy.fQATrackHistosSecRecSsc)
672 ,fFFHistosRecEffRec(copy.fFFHistosRecEffRec)
673 ,fFFHistosSecRecNS(copy.fFFHistosSecRecNS)
674 ,fFFHistosSecRecS(copy.fFFHistosSecRecS)
675 ,fFFHistosSecRecSsc(copy.fFFHistosSecRecSsc)
677 ,fh1BckgMult0(copy.fh1BckgMult0)
678 ,fh1BckgMult1(copy.fh1BckgMult1)
679 ,fh1BckgMult2(copy.fh1BckgMult2)
680 ,fh1BckgMult3(copy.fh1BckgMult3)
681 ,fh1BckgMult4(copy.fh1BckgMult4)
682 ,fh1FractionPtEmbedded(copy.fh1FractionPtEmbedded)
683 ,fh1IndexEmbedded(copy.fh1IndexEmbedded)
684 ,fh2DeltaPtVsJetPtEmbedded(copy.fh2DeltaPtVsJetPtEmbedded)
685 ,fh2DeltaPtVsRecJetPtEmbedded(copy.fh2DeltaPtVsRecJetPtEmbedded)
686 ,fh1DeltaREmbedded(copy.fh1DeltaREmbedded)
687 ,fQABckgHisto0RecCuts(copy.fQABckgHisto0RecCuts)
688 ,fQABckgHisto0Gen(copy.fQABckgHisto0Gen)
689 ,fQABckgHisto1RecCuts(copy.fQABckgHisto1RecCuts)
690 ,fQABckgHisto1Gen(copy.fQABckgHisto1Gen)
691 ,fQABckgHisto2RecCuts(copy.fQABckgHisto2RecCuts)
692 ,fQABckgHisto2Gen(copy.fQABckgHisto2Gen)
693 ,fQABckgHisto3RecCuts(copy.fQABckgHisto3RecCuts)
694 ,fQABckgHisto3Gen(copy.fQABckgHisto3Gen)
695 ,fQABckgHisto4RecCuts(copy.fQABckgHisto4RecCuts)
696 ,fQABckgHisto4Gen(copy.fQABckgHisto4Gen)
697 ,fFFBckgHisto0RecCuts(copy.fFFBckgHisto0RecCuts)
698 ,fFFBckgHisto0Gen(copy.fFFBckgHisto0Gen)
699 ,fFFBckgHisto1RecCuts(copy.fFFBckgHisto1RecCuts)
700 ,fFFBckgHisto1Gen(copy.fFFBckgHisto1Gen)
701 ,fFFBckgHisto2RecCuts(copy.fFFBckgHisto2RecCuts)
702 ,fFFBckgHisto2Gen(copy.fFFBckgHisto2Gen)
703 ,fFFBckgHisto3RecCuts(copy.fFFBckgHisto3RecCuts)
704 ,fFFBckgHisto3Gen(copy.fFFBckgHisto3Gen)
705 ,fFFBckgHisto4RecCuts(copy.fFFBckgHisto4RecCuts)
706 ,fFFBckgHisto4Gen(copy.fFFBckgHisto4Gen)
707 ,fFFBckgHisto0RecEffRec(copy.fFFBckgHisto0RecEffRec)
708 ,fFFBckgHisto0SecRecNS(copy.fFFBckgHisto0SecRecNS)
709 ,fFFBckgHisto0SecRecS(copy.fFFBckgHisto0SecRecS)
710 ,fFFBckgHisto0SecRecSsc(copy.fFFBckgHisto0SecRecSsc)
712 ,fProNtracksLeadingJet(copy.fProNtracksLeadingJet)
713 ,fProDelR80pcPt(copy.fProDelR80pcPt)
714 ,fProNtracksLeadingJetGen(copy.fProNtracksLeadingJetGen)
715 ,fProDelR80pcPtGen(copy.fProDelR80pcPtGen)
716 ,fProNtracksLeadingJetBgrPerp2(copy.fProNtracksLeadingJetBgrPerp2)
717 ,fProNtracksLeadingJetRecPrim(copy.fProNtracksLeadingJetRecPrim)
718 ,fProDelR80pcPtRecPrim(copy.fProDelR80pcPtRecPrim)
719 ,fProNtracksLeadingJetRecSecNS(copy.fProNtracksLeadingJetRecSecNS)
720 ,fProNtracksLeadingJetRecSecS(copy.fProNtracksLeadingJetRecSecS)
721 ,fProNtracksLeadingJetRecSecSsc(copy.fProNtracksLeadingJetRecSecSsc)
722 ,fRandom(copy.fRandom)
723 ,fOnlyLeadingJets(copy.fOnlyLeadingJets)
724 ,fMCPtHardCut(copy.fMCPtHardCut)
725 ,fAnaUtils(copy.fAnaUtils)
727 ,fNumInclusivePIDtasks(copy.fNumInclusivePIDtasks)
728 ,fNumJetPIDtasks(copy.fNumJetPIDtasks)
729 ,fNameInclusivePIDtask(0x0)
730 ,fNameJetPIDtask(0x0)
731 ,fInclusivePIDtask(0x0)
733 ,fUseInclusivePIDtask(copy.fUseInclusivePIDtask)
734 ,fUseJetPIDtask(copy.fUseJetPIDtask)
738 fBckgType[0] = copy.fBckgType[0];
739 fBckgType[1] = copy.fBckgType[1];
740 fBckgType[2] = copy.fBckgType[2];
741 fBckgType[3] = copy.fBckgType[3];
742 fBckgType[4] = copy.fBckgType[4];
745 for(Int_t ii=0; ii<5; ii++){
746 fProDelRPtSum[ii] = copy.fProDelRPtSum[ii];
747 fProDelRPtSumGen[ii] = copy.fProDelRPtSumGen[ii];
748 fProDelRPtSumBgrPerp2[ii] = copy.fProDelRPtSumBgrPerp2[ii];
749 fProDelRPtSumRecPrim[ii] = copy.fProDelRPtSumRecPrim[ii];
750 fProDelRPtSumRecSecNS[ii] = copy.fProDelRPtSumRecSecNS[ii];
751 fProDelRPtSumRecSecS[ii] = copy.fProDelRPtSumRecSecS[ii];
752 fProDelRPtSumRecSecSsc[ii] = copy.fProDelRPtSumRecSecSsc[ii];
755 for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
756 fIDFFHistosRecCuts[i] = 0x0;
757 if (copy.fIDFFHistosRecCuts[i])
758 fIDFFHistosRecCuts[i] = copy.fIDFFHistosRecCuts[i];
760 fIDFFHistosGen[i] = 0x0;
761 if (copy.fIDFFHistosGen[i])
762 fIDFFHistosGen[i] = copy.fIDFFHistosGen[i];
765 fhDCA_XY_prim_MCID[i] = 0x0;
766 if (copy.fhDCA_XY_prim_MCID[i])
767 fhDCA_XY_prim_MCID[i] = copy.fhDCA_XY_prim_MCID[i];
769 fhDCA_Z_prim_MCID[i] = 0x0;
770 if (copy.fhDCA_Z_prim_MCID[i])
771 fhDCA_Z_prim_MCID[i] = copy.fhDCA_Z_prim_MCID[i];
773 fhDCA_XY_sec_MCID[i] = 0x0;
774 if (copy.fhDCA_XY_sec_MCID[i])
775 fhDCA_XY_sec_MCID[i] = copy.fhDCA_XY_sec_MCID[i];
777 fhDCA_Z_sec_MCID[i] = 0x0;
778 if (copy.fhDCA_Z_sec_MCID[i])
779 fhDCA_Z_sec_MCID[i] = copy.fhDCA_Z_sec_MCID[i];
782 if (fNumInclusivePIDtasks > 0) {
783 fNameInclusivePIDtask = new TString[fNumInclusivePIDtasks];
784 fInclusivePIDtask = new AliAnalysisTaskPID*[fNumInclusivePIDtasks];
786 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
787 fNameInclusivePIDtask[i] = "";
788 fInclusivePIDtask[i] = 0x0;
790 if (copy.fNameInclusivePIDtask[i])
791 fNameInclusivePIDtask[i] = copy.fNameInclusivePIDtask[i];
793 if (copy.fInclusivePIDtask[i])
794 fInclusivePIDtask[i] = copy.fInclusivePIDtask[i];
798 if (fNumJetPIDtasks > 0) {
799 fNameJetPIDtask = new TString[fNumJetPIDtasks];
800 fJetPIDtask = new AliAnalysisTaskPID*[fNumJetPIDtasks];
802 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
803 fNameJetPIDtask[i] = "";
804 fJetPIDtask[i] = 0x0;
806 if (copy.fNameJetPIDtask[i])
807 fNameJetPIDtask[i] = copy.fNameJetPIDtask[i];
809 if (copy.fJetPIDtask[i])
810 fJetPIDtask[i] = copy.fJetPIDtask[i];
815 // _________________________________________________________________________________________________________________________________
816 AliAnalysisTaskIDFragmentationFunction& AliAnalysisTaskIDFragmentationFunction::operator=(const AliAnalysisTaskIDFragmentationFunction& o)
822 AliAnalysisTaskSE::operator=(o);
825 fAODJets = o.fAODJets;
826 fAODExtension = o.fAODExtension;
827 fNonStdFile = o.fNonStdFile;
828 fBranchRecJets = o.fBranchRecJets;
829 fBranchRecBckgClusters = o.fBranchRecBckgClusters;
830 fBranchGenJets = o.fBranchGenJets;
831 fBranchEmbeddedJets = o.fBranchEmbeddedJets;
832 fTrackTypeGen = o.fTrackTypeGen;
833 fJetTypeGen = o.fJetTypeGen;
834 fJetTypeRecEff = o.fJetTypeRecEff;
835 fUseAODInputJets = o.fUseAODInputJets;
836 fFilterMask = o.fFilterMask;
837 fUsePhysicsSelection = o.fUsePhysicsSelection;
838 fEvtSelectionMask = o.fEvtSelectionMask;
839 fEventClass = o.fEventClass;
840 fMaxVertexZ = o.fMaxVertexZ;
841 fTrackPtCut = o.fTrackPtCut;
842 fTrackEtaMin = o.fTrackEtaMin;
843 fTrackEtaMax = o.fTrackEtaMax;
844 fTrackPhiMin = o.fTrackPhiMin;
845 fTrackPhiMax = o.fTrackPhiMax;
846 fUseExtraTracks = o.fUseExtraTracks;
847 fUseExtraTracksBgr = o.fUseExtraTracksBgr;
848 fCutFractionPtEmbedded = o.fCutFractionPtEmbedded;
849 fUseEmbeddedJetAxis = o.fUseEmbeddedJetAxis;
850 fUseEmbeddedJetPt = o.fUseEmbeddedJetPt;
851 fJetPtCut = o.fJetPtCut;
852 fJetEtaMin = o.fJetEtaMin;
853 fJetEtaMax = o.fJetEtaMax;
854 fJetPhiMin = o.fJetPhiMin;
855 fJetPhiMax = o.fJetPhiMin;
856 fFFRadius = o.fFFRadius;
857 fFFMinLTrackPt = o.fFFMinLTrackPt;
858 fFFMaxTrackPt = o.fFFMaxTrackPt;
859 fFFMinnTracks = o.fFFMinnTracks;
860 fFFBckgRadius = o.fFFBckgRadius;
861 fBckgMode = o.fBckgMode;
864 fIDFFMode = o.fIDFFMode;
865 fEffMode = o.fEffMode;
867 fBckgType[0] = o.fBckgType[0];
868 fBckgType[1] = o.fBckgType[1];
869 fBckgType[2] = o.fBckgType[2];
870 fBckgType[3] = o.fBckgType[3];
871 fBckgType[4] = o.fBckgType[4];
872 fAvgTrials = o.fAvgTrials;
873 fTracksRecCuts = o.fTracksRecCuts;
874 fTracksRecCutsEfficiency = o.fTracksRecCutsEfficiency;
875 fTracksGen = o.fTracksGen;
876 fTracksAODMCCharged = o.fTracksAODMCCharged;
877 fTracksAODMCChargedSecNS = o.fTracksAODMCChargedSecNS;
878 fTracksAODMCChargedSecS = o.fTracksAODMCChargedSecS;
879 fTracksRecQualityCuts = o.fTracksRecQualityCuts;
880 fJetsRec = o.fJetsRec;
881 fJetsRecCuts = o.fJetsRecCuts;
882 fJetsGen = o.fJetsGen;
883 fJetsRecEff = o.fJetsRecEff;
884 fJetsEmbedded = o.fJetsEmbedded;
885 fBckgJetsRec = o.fBckgJetsRec;
886 fBckgJetsRecCuts = o.fBckgJetsRecCuts;
887 fBckgJetsGen = o.fBckgJetsGen;
888 fQATrackHistosRecCuts = o.fQATrackHistosRecCuts;
889 fQATrackHistosGen = o.fQATrackHistosGen;
890 fQAJetHistosRec = o.fQAJetHistosRec;
891 fQAJetHistosRecCuts = o.fQAJetHistosRecCuts;
892 fQAJetHistosRecCutsLeading = o.fQAJetHistosRecCutsLeading;
893 fQAJetHistosGen = o.fQAJetHistosGen;
894 fQAJetHistosGenLeading = o.fQAJetHistosGenLeading;
895 fQAJetHistosRecEffLeading = o.fQAJetHistosRecEffLeading;
896 fFFHistosRecCuts = o.fFFHistosRecCuts;
897 fFFHistosRecCutsInc = o.fFFHistosRecCutsInc;
898 fFFHistosRecLeadingTrack = o.fFFHistosRecLeadingTrack;
899 fFFHistosGen = o.fFFHistosGen;
900 fFFHistosGenInc = o.fFFHistosGenInc;
901 fFFHistosGenLeadingTrack = o.fFFHistosGenLeadingTrack;
902 fQATrackHighPtThreshold = o.fQATrackHighPtThreshold;
903 fFFNBinsJetPt = o.fFFNBinsJetPt;
904 fFFJetPtMin = o.fFFJetPtMin;
905 fFFJetPtMax = o.fFFJetPtMax;
906 fFFNBinsPt = o.fFFNBinsPt;
907 fFFPtMin = o.fFFPtMin;
908 fFFPtMax = o.fFFPtMax;
909 fFFNBinsXi = o.fFFNBinsXi;
910 fFFXiMin = o.fFFXiMin;
911 fFFXiMax = o.fFFXiMax;
912 fFFNBinsZ = o.fFFNBinsZ;
915 fQAJetNBinsPt = o.fQAJetNBinsPt;
916 fQAJetPtMin = o.fQAJetPtMin;
917 fQAJetPtMax = o.fQAJetPtMax;
918 fQAJetNBinsEta = o.fQAJetNBinsEta;
919 fQAJetEtaMin = o.fQAJetEtaMin;
920 fQAJetEtaMax = o.fQAJetEtaMax;
921 fQAJetNBinsPhi = o.fQAJetNBinsPhi;
922 fQAJetPhiMin = o.fQAJetPhiMin;
923 fQAJetPhiMax = o.fQAJetPhiMax;
924 fQATrackNBinsPt = o.fQATrackNBinsPt;
925 fQATrackPtMin = o.fQATrackPtMin;
926 fQATrackPtMax = o.fQATrackPtMax;
927 fQATrackNBinsEta = o.fQATrackNBinsEta;
928 fQATrackEtaMin = o.fQATrackEtaMin;
929 fQATrackEtaMax = o.fQATrackEtaMax;
930 fQATrackNBinsPhi = o.fQATrackNBinsPhi;
931 fQATrackPhiMin = o.fQATrackPhiMin;
932 fQATrackPhiMax = o.fQATrackPhiMax;
933 fCommonHistList = o.fCommonHistList;
934 fh1EvtSelection = o.fh1EvtSelection;
935 fh1VtxSelection = o.fh1VtxSelection;
936 fh1VertexNContributors = o.fh1VertexNContributors;
937 fh1VertexZ = o.fh1VertexZ;
938 fh1EvtMult = o.fh1EvtMult;
939 fh1EvtCent = o.fh1EvtCent;
941 fh1Trials = o.fh1Trials;
942 fh1PtHard = o.fh1PtHard;
943 fh1PtHardTrials = o.fh1PtHardTrials;
944 fh1EvtsPtHardCut = o.fh1EvtsPtHardCut;
945 fh1nRecJetsCuts = o.fh1nRecJetsCuts;
946 fh1nGenJets = o.fh1nGenJets;
947 fh1nRecEffJets = o.fh1nRecEffJets;
948 fh1nEmbeddedJets = o.fh1nEmbeddedJets;
949 fh2PtRecVsGenPrim = o.fh2PtRecVsGenPrim;
950 fh2PtRecVsGenSec = o.fh2PtRecVsGenSec;
951 fQATrackHistosRecEffGen = o.fQATrackHistosRecEffGen;
952 fQATrackHistosRecEffRec = o.fQATrackHistosRecEffRec;
953 fQATrackHistosSecRecNS = o.fQATrackHistosSecRecNS;
954 fQATrackHistosSecRecS = o.fQATrackHistosSecRecS;
955 fQATrackHistosSecRecSsc = o.fQATrackHistosSecRecSsc;
956 fFFHistosRecEffRec = o.fFFHistosRecEffRec;
957 fFFHistosSecRecNS = o.fFFHistosSecRecNS;
958 fFFHistosSecRecS = o.fFFHistosSecRecS;
959 fFFHistosSecRecSsc = o.fFFHistosSecRecSsc;
961 fh1BckgMult0 = o.fh1BckgMult0;
962 fh1BckgMult1 = o.fh1BckgMult1;
963 fh1BckgMult2 = o.fh1BckgMult2;
964 fh1BckgMult3 = o.fh1BckgMult3;
965 fh1BckgMult4 = o.fh1BckgMult4;
966 fh1FractionPtEmbedded = o.fh1FractionPtEmbedded;
967 fh1IndexEmbedded = o.fh1IndexEmbedded;
968 fh2DeltaPtVsJetPtEmbedded = o.fh2DeltaPtVsJetPtEmbedded;
969 fh2DeltaPtVsRecJetPtEmbedded = o.fh2DeltaPtVsRecJetPtEmbedded;
970 fh1DeltaREmbedded = o.fh1DeltaREmbedded;
971 fhDCA_XY = o.fhDCA_XY;
973 fhJetPtRefMultEta5 = o.fhJetPtRefMultEta5;
974 fhJetPtRefMultEta8 = o.fhJetPtRefMultEta8;
975 fhJetPtMultPercent = o.fhJetPtMultPercent;
976 fQABckgHisto0RecCuts = o.fQABckgHisto0RecCuts;
977 fQABckgHisto0Gen = o.fQABckgHisto0Gen;
978 fQABckgHisto1RecCuts = o.fQABckgHisto1RecCuts;
979 fQABckgHisto1Gen = o.fQABckgHisto1Gen;
980 fQABckgHisto2RecCuts = o.fQABckgHisto2RecCuts;
981 fQABckgHisto2Gen = o.fQABckgHisto2Gen;
982 fQABckgHisto3RecCuts = o.fQABckgHisto3RecCuts;
983 fQABckgHisto3Gen = o.fQABckgHisto3Gen;
984 fQABckgHisto4RecCuts = o.fQABckgHisto4RecCuts;
985 fQABckgHisto4Gen = o.fQABckgHisto4Gen;
986 fFFBckgHisto0RecCuts = o.fFFBckgHisto0RecCuts;
987 fFFBckgHisto0Gen = o.fFFBckgHisto0Gen;
988 fFFBckgHisto1RecCuts = o.fFFBckgHisto1RecCuts;
989 fFFBckgHisto1Gen = o.fFFBckgHisto1Gen;
990 fFFBckgHisto2RecCuts = o.fFFBckgHisto2RecCuts;
991 fFFBckgHisto2Gen = o.fFFBckgHisto2Gen;
992 fFFBckgHisto3RecCuts = o.fFFBckgHisto3RecCuts;
993 fFFBckgHisto3Gen = o.fFFBckgHisto3Gen;
994 fFFBckgHisto4RecCuts = o.fFFBckgHisto4RecCuts;
995 fFFBckgHisto4Gen = o.fFFBckgHisto4Gen;
996 fFFBckgHisto0RecEffRec = o.fFFBckgHisto0RecEffRec;
997 fFFBckgHisto0SecRecNS = o.fFFBckgHisto0SecRecNS;
998 fFFBckgHisto0SecRecS = o.fFFBckgHisto0SecRecS;
999 fFFBckgHisto0SecRecSsc = o.fFFBckgHisto0SecRecSsc;
1000 fProNtracksLeadingJet = o.fProNtracksLeadingJet;
1001 fProDelR80pcPt = o.fProDelR80pcPt;
1002 fProNtracksLeadingJetGen = o.fProNtracksLeadingJetGen;
1003 fProDelR80pcPtGen = o.fProDelR80pcPtGen;
1004 fProNtracksLeadingJetBgrPerp2 = o.fProNtracksLeadingJetBgrPerp2;
1005 fProNtracksLeadingJetRecPrim = o.fProNtracksLeadingJetRecPrim;
1006 fProDelR80pcPtRecPrim = o.fProDelR80pcPtRecPrim;
1007 fProNtracksLeadingJetRecSecNS = o.fProNtracksLeadingJetRecSecNS;
1008 fProNtracksLeadingJetRecSecS = o.fProNtracksLeadingJetRecSecS;
1009 fProNtracksLeadingJetRecSecSsc = o.fProNtracksLeadingJetRecSecSsc;
1010 fRandom = o.fRandom;
1011 fOnlyLeadingJets = o.fOnlyLeadingJets;
1012 fMCPtHardCut = o.fMCPtHardCut;
1013 fAnaUtils = o.fAnaUtils;
1016 fUseInclusivePIDtask = o.fUseInclusivePIDtask;
1017 fUseJetPIDtask = o.fUseJetPIDtask;
1021 for(Int_t ii=0; ii<5; ii++){
1022 fProDelRPtSum[ii] = o.fProDelRPtSum[ii];
1023 fProDelRPtSumGen[ii] = o.fProDelRPtSumGen[ii];
1024 fProDelRPtSumBgrPerp2[ii] = o.fProDelRPtSumBgrPerp2[ii];
1025 fProDelRPtSumRecPrim[ii] = o.fProDelRPtSumRecPrim[ii];
1026 fProDelRPtSumRecSecNS[ii] = o.fProDelRPtSumRecSecNS[ii];
1027 fProDelRPtSumRecSecS[ii] = o.fProDelRPtSumRecSecS[ii];
1028 fProDelRPtSumRecSecSsc[ii] = o.fProDelRPtSumRecSecSsc[ii];
1031 for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
1032 fIDFFHistosRecCuts[i] = 0x0;
1033 if (o.fIDFFHistosRecCuts[i])
1034 fIDFFHistosRecCuts[i] = o.fIDFFHistosRecCuts[i];
1036 fIDFFHistosGen[i] = 0x0;
1037 if (o.fIDFFHistosGen[i])
1038 fIDFFHistosGen[i] = o.fIDFFHistosGen[i];
1040 fhDCA_XY_prim_MCID[i] = 0x0;
1041 if (o.fhDCA_XY_prim_MCID[i])
1042 fhDCA_XY_prim_MCID[i] = o.fhDCA_XY_prim_MCID[i];
1044 fhDCA_Z_prim_MCID[i] = 0x0;
1045 if (o.fhDCA_Z_prim_MCID[i])
1046 fhDCA_Z_prim_MCID[i] = o.fhDCA_Z_prim_MCID[i];
1048 fhDCA_XY_sec_MCID[i] = 0x0;
1049 if (o.fhDCA_XY_sec_MCID[i])
1050 fhDCA_XY_sec_MCID[i] = o.fhDCA_XY_sec_MCID[i];
1052 fhDCA_Z_sec_MCID[i] = 0x0;
1053 if (o.fhDCA_Z_sec_MCID[i])
1054 fhDCA_Z_sec_MCID[i] = o.fhDCA_Z_sec_MCID[i];
1059 if (fNumJetPIDtasks > 0) {
1060 delete [] fNameJetPIDtask;
1061 fNameJetPIDtask = 0x0;
1063 delete [] fJetPIDtask;
1067 fNumJetPIDtasks = o.fNumJetPIDtasks;
1070 if (fNumInclusivePIDtasks > 0) {
1071 delete [] fNameInclusivePIDtask;
1072 fNameInclusivePIDtask = 0x0;
1074 delete [] fInclusivePIDtask;
1075 fInclusivePIDtask = 0x0;
1078 fNumInclusivePIDtasks = o.fNumInclusivePIDtasks;
1081 if (fNumInclusivePIDtasks > 0) {
1082 fNameInclusivePIDtask = new TString[fNumInclusivePIDtasks];
1083 fInclusivePIDtask = new AliAnalysisTaskPID*[fNumInclusivePIDtasks];
1085 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
1086 fNameInclusivePIDtask[i] = "";
1087 fInclusivePIDtask[i] = 0x0;
1089 if (o.fNameInclusivePIDtask[i])
1090 fNameInclusivePIDtask[i] = o.fNameInclusivePIDtask[i];
1092 if (o.fInclusivePIDtask[i])
1093 fInclusivePIDtask[i] = o.fInclusivePIDtask[i];
1097 if (fNumJetPIDtasks > 0) {
1098 fNameJetPIDtask = new TString[fNumJetPIDtasks];
1099 fJetPIDtask = new AliAnalysisTaskPID*[fNumJetPIDtasks];
1101 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
1102 fNameJetPIDtask[i] = "";
1103 fJetPIDtask[i] = 0x0;
1105 if (o.fNameJetPIDtask[i])
1106 fNameJetPIDtask[i] = o.fNameJetPIDtask[i];
1108 if (o.fJetPIDtask[i])
1109 fJetPIDtask[i] = o.fJetPIDtask[i];
1116 //___________________________________________________________________________
1117 AliAnalysisTaskIDFragmentationFunction::~AliAnalysisTaskIDFragmentationFunction()
1121 if(fTracksRecCuts) delete fTracksRecCuts;
1122 if(fTracksRecCutsEfficiency) delete fTracksRecCutsEfficiency;
1123 if(fTracksGen) delete fTracksGen;
1124 if(fTracksAODMCCharged) delete fTracksAODMCCharged;
1125 if(fTracksAODMCChargedSecNS) delete fTracksAODMCChargedSecNS;
1126 if(fTracksAODMCChargedSecS) delete fTracksAODMCChargedSecS;
1127 if(fTracksRecQualityCuts) delete fTracksRecQualityCuts;
1128 if(fJetsRec) delete fJetsRec;
1129 if(fJetsRecCuts) delete fJetsRecCuts;
1130 if(fJetsGen) delete fJetsGen;
1131 if(fJetsRecEff) delete fJetsRecEff;
1132 if(fJetsEmbedded) delete fJetsEmbedded;
1135 (fBckgType[0]==kBckgClusters || fBckgType[1]==kBckgClusters || fBckgType[2]==kBckgClusters || fBckgType[3]==kBckgClusters || fBckgType[4]==kBckgClusters ||
1136 fBckgType[0]==kBckgClustersOutLeading || fBckgType[1]==kBckgClustersOutLeading || fBckgType[2]==kBckgClustersOutLeading ||
1137 fBckgType[3]==kBckgClustersOutLeading || fBckgType[4]==kBckgClustersOutLeading)){
1139 if(fBckgJetsRec) delete fBckgJetsRec;
1140 if(fBckgJetsRecCuts) delete fBckgJetsRecCuts;
1141 if(fBckgJetsGen) delete fBckgJetsGen;
1143 if(fRandom) delete fRandom;
1145 delete [] fNameInclusivePIDtask;
1146 fNameInclusivePIDtask = 0x0;
1148 delete [] fInclusivePIDtask;
1149 fInclusivePIDtask = 0x0;
1151 delete [] fNameJetPIDtask;
1152 fNameJetPIDtask = 0x0;
1154 delete [] fJetPIDtask;
1161 //______________________________________________________________________________________________________
1162 AliAnalysisTaskIDFragmentationFunction::AliFragFuncHistos::AliFragFuncHistos(const char* name,
1163 Int_t nJetPt, Float_t jetPtMin, Float_t jetPtMax,
1164 Int_t nPt, Float_t ptMin, Float_t ptMax,
1165 Int_t nXi, Float_t xiMin, Float_t xiMax,
1166 Int_t nZ , Float_t zMin , Float_t zMax)
1168 ,fNBinsJetPt(nJetPt)
1169 ,fJetPtMin(jetPtMin)
1170 ,fJetPtMax(jetPtMax)
1186 // default constructor
1190 //___________________________________________________________________________
1191 AliAnalysisTaskIDFragmentationFunction::AliFragFuncHistos::AliFragFuncHistos(const AliFragFuncHistos& copy)
1193 ,fNBinsJetPt(copy.fNBinsJetPt)
1194 ,fJetPtMin(copy.fJetPtMin)
1195 ,fJetPtMax(copy.fJetPtMax)
1196 ,fNBinsPt(copy.fNBinsPt)
1197 ,fPtMin(copy.fPtMin)
1198 ,fPtMax(copy.fPtMax)
1199 ,fNBinsXi(copy.fNBinsXi)
1200 ,fXiMin(copy.fXiMin)
1201 ,fXiMax(copy.fXiMax)
1202 ,fNBinsZ(copy.fNBinsZ)
1205 ,fh2TrackPt(copy.fh2TrackPt)
1208 ,fh1JetPt(copy.fh1JetPt)
1209 ,fNameFF(copy.fNameFF)
1214 //_______________________________________________________________________________________________________________________________________________________________
1215 AliAnalysisTaskIDFragmentationFunction::AliFragFuncHistos& AliAnalysisTaskIDFragmentationFunction::AliFragFuncHistos::operator=(const AliAnalysisTaskIDFragmentationFunction::AliFragFuncHistos& o)
1220 TObject::operator=(o);
1221 fNBinsJetPt = o.fNBinsJetPt;
1222 fJetPtMin = o.fJetPtMin;
1223 fJetPtMax = o.fJetPtMax;
1224 fNBinsPt = o.fNBinsPt;
1227 fNBinsXi = o.fNBinsXi;
1230 fNBinsZ = o.fNBinsZ;
1233 fh2TrackPt = o.fh2TrackPt;
1236 fh1JetPt = o.fh1JetPt;
1237 fNameFF = o.fNameFF;
1243 //_________________________________________________________
1244 AliAnalysisTaskIDFragmentationFunction::AliFragFuncHistos::~AliFragFuncHistos()
1248 if(fh1JetPt) delete fh1JetPt;
1249 if(fh2TrackPt) delete fh2TrackPt;
1250 if(fh2Xi) delete fh2Xi;
1251 if(fh2Z) delete fh2Z;
1254 //_________________________________________________________________
1255 void AliAnalysisTaskIDFragmentationFunction::AliFragFuncHistos::DefineHistos()
1259 fh1JetPt = new TH1F(Form("fh1FFJetPt%s", fNameFF.Data()),"",fNBinsJetPt,fJetPtMin,fJetPtMax);
1260 fh2TrackPt = new TH2F(Form("fh2FFTrackPt%s",fNameFF.Data()),"",fNBinsJetPt, fJetPtMin, fJetPtMax,fNBinsPt, fPtMin, fPtMax);
1261 fh2Z = new TH2F(Form("fh2FFZ%s",fNameFF.Data()),"",fNBinsJetPt, fJetPtMin, fJetPtMax, fNBinsZ, fZMin, fZMax);
1262 fh2Xi = new TH2F(Form("fh2FFXi%s",fNameFF.Data()),"",fNBinsJetPt, fJetPtMin, fJetPtMax, fNBinsXi, fXiMin, fXiMax);
1264 AliAnalysisTaskIDFragmentationFunction::SetProperties(fh1JetPt, "p_{T} [GeV/c]", "entries");
1265 AliAnalysisTaskIDFragmentationFunction::SetProperties(fh2TrackPt,"jet p_{T} [GeV/c]","p_{T} [GeV/c]","entries");
1266 AliAnalysisTaskIDFragmentationFunction::SetProperties(fh2Xi,"jet p_{T} [GeV/c]","#xi", "entries");
1267 AliAnalysisTaskIDFragmentationFunction::SetProperties(fh2Z,"jet p_{T} [GeV/c]","z","entries");
1270 //_______________________________________________________________________________________________________________
1271 void AliAnalysisTaskIDFragmentationFunction::AliFragFuncHistos::FillFF(Float_t trackPt, Float_t jetPt, Bool_t incrementJetPt, Float_t norm,
1272 Bool_t scaleStrangeness, Float_t scaleFacStrangeness)
1276 if(incrementJetPt && norm) fh1JetPt->Fill(jetPt,1/norm);
1277 else if(incrementJetPt) fh1JetPt->Fill(jetPt);
1279 // Added for proper normalization of FF background estimation
1280 // when zero track are found in the background region
1281 if((int)trackPt==-1) return;
1283 if(norm)fh2TrackPt->Fill(jetPt,trackPt,1/norm);
1284 else if(scaleStrangeness) fh2TrackPt->Fill(jetPt,trackPt,scaleFacStrangeness);
1285 else fh2TrackPt->Fill(jetPt,trackPt);
1287 Double_t z = -1., xi = -1.;
1288 AliAnalysisTaskPID::GetJetTrackObservables(trackPt, jetPt, z, xi);
1292 fh2Xi->Fill(jetPt,xi,1/norm);
1293 fh2Z->Fill(jetPt,z,1/norm);
1295 else if(scaleStrangeness){
1296 fh2Xi->Fill(jetPt,xi,scaleFacStrangeness);
1297 fh2Z->Fill(jetPt,z,scaleFacStrangeness);
1300 fh2Xi->Fill(jetPt,xi);
1301 fh2Z->Fill(jetPt,z);
1305 //_________________________________________________________________________________
1306 void AliAnalysisTaskIDFragmentationFunction::AliFragFuncHistos::AddToOutput(TList* list) const
1308 // add histos to list
1310 list->Add(fh1JetPt);
1312 list->Add(fh2TrackPt);
1317 //_________________________________________________________________________________________________________
1318 AliAnalysisTaskIDFragmentationFunction::AliFragFuncQAJetHistos::AliFragFuncQAJetHistos(const char* name,
1319 Int_t nPt, Float_t ptMin, Float_t ptMax,
1320 Int_t nEta, Float_t etaMin, Float_t etaMax,
1321 Int_t nPhi, Float_t phiMin, Float_t phiMax)
1336 // default constructor
1339 //____________________________________________________________________________________
1340 AliAnalysisTaskIDFragmentationFunction::AliFragFuncQAJetHistos::AliFragFuncQAJetHistos(const AliFragFuncQAJetHistos& copy)
1342 ,fNBinsPt(copy.fNBinsPt)
1343 ,fPtMin(copy.fPtMin)
1344 ,fPtMax(copy.fPtMax)
1345 ,fNBinsEta(copy.fNBinsEta)
1346 ,fEtaMin(copy.fEtaMin)
1347 ,fEtaMax(copy.fEtaMax)
1348 ,fNBinsPhi(copy.fNBinsPhi)
1349 ,fPhiMin(copy.fPhiMin)
1350 ,fPhiMax(copy.fPhiMax)
1351 ,fh2EtaPhi(copy.fh2EtaPhi)
1353 ,fNameQAJ(copy.fNameQAJ)
1358 //________________________________________________________________________________________________________________________________________________________________________
1359 AliAnalysisTaskIDFragmentationFunction::AliFragFuncQAJetHistos& AliAnalysisTaskIDFragmentationFunction::AliFragFuncQAJetHistos::operator=(const AliAnalysisTaskIDFragmentationFunction::AliFragFuncQAJetHistos& o)
1364 TObject::operator=(o);
1365 fNBinsPt = o.fNBinsPt;
1368 fNBinsEta = o.fNBinsEta;
1369 fEtaMin = o.fEtaMin;
1370 fEtaMax = o.fEtaMax;
1371 fNBinsPhi = o.fNBinsPhi;
1372 fPhiMin = o.fPhiMin;
1373 fPhiMax = o.fPhiMax;
1374 fh2EtaPhi = o.fh2EtaPhi;
1376 fNameQAJ = o.fNameQAJ;
1382 //______________________________________________________________
1383 AliAnalysisTaskIDFragmentationFunction::AliFragFuncQAJetHistos::~AliFragFuncQAJetHistos()
1387 if(fh2EtaPhi) delete fh2EtaPhi;
1388 if(fh1Pt) delete fh1Pt;
1391 //____________________________________________________________________
1392 void AliAnalysisTaskIDFragmentationFunction::AliFragFuncQAJetHistos::DefineHistos()
1394 // book jet QA histos
1396 fh2EtaPhi = new TH2F(Form("fh2JetQAEtaPhi%s", fNameQAJ.Data()), Form("%s: #eta - #phi distribution", fNameQAJ.Data()), fNBinsEta, fEtaMin, fEtaMax, fNBinsPhi, fPhiMin, fPhiMax);
1397 fh1Pt = new TH1F(Form("fh1JetQAPt%s", fNameQAJ.Data()), Form("%s: p_{T} distribution", fNameQAJ.Data()), fNBinsPt, fPtMin, fPtMax);
1399 AliAnalysisTaskIDFragmentationFunction::SetProperties(fh2EtaPhi, "#eta", "#phi");
1400 AliAnalysisTaskIDFragmentationFunction::SetProperties(fh1Pt, "p_{T} [GeV/c]", "entries");
1403 //____________________________________________________________________________________________________
1404 void AliAnalysisTaskIDFragmentationFunction::AliFragFuncQAJetHistos::FillJetQA(Float_t eta, Float_t phi, Float_t pt)
1406 // fill jet QA histos
1408 fh2EtaPhi->Fill( eta, phi);
1412 //____________________________________________________________________________________
1413 void AliAnalysisTaskIDFragmentationFunction::AliFragFuncQAJetHistos::AddToOutput(TList* list) const
1415 // add histos to list
1417 list->Add(fh2EtaPhi);
1421 //___________________________________________________________________________________________________________
1422 AliAnalysisTaskIDFragmentationFunction::AliFragFuncQATrackHistos::AliFragFuncQATrackHistos(const char* name,
1423 Int_t nPt, Float_t ptMin, Float_t ptMax,
1424 Int_t nEta, Float_t etaMin, Float_t etaMax,
1425 Int_t nPhi, Float_t phiMin, Float_t phiMax,
1437 ,fHighPtThreshold(ptThresh)
1444 // default constructor
1447 //__________________________________________________________________________________________
1448 AliAnalysisTaskIDFragmentationFunction::AliFragFuncQATrackHistos::AliFragFuncQATrackHistos(const AliFragFuncQATrackHistos& copy)
1450 ,fNBinsPt(copy.fNBinsPt)
1451 ,fPtMin(copy.fPtMin)
1452 ,fPtMax(copy.fPtMax)
1453 ,fNBinsEta(copy.fNBinsEta)
1454 ,fEtaMin(copy.fEtaMin)
1455 ,fEtaMax(copy.fEtaMax)
1456 ,fNBinsPhi(copy.fNBinsPhi)
1457 ,fPhiMin(copy.fPhiMin)
1458 ,fPhiMax(copy.fPhiMax)
1459 ,fHighPtThreshold(copy.fHighPtThreshold)
1460 ,fh2EtaPhi(copy.fh2EtaPhi)
1462 ,fh2HighPtEtaPhi(copy.fh2HighPtEtaPhi)
1463 ,fh2PhiPt(copy.fh2PhiPt)
1464 ,fNameQAT(copy.fNameQAT)
1469 // _____________________________________________________________________________________________________________________________________________________________________________
1470 AliAnalysisTaskIDFragmentationFunction::AliFragFuncQATrackHistos& AliAnalysisTaskIDFragmentationFunction::AliFragFuncQATrackHistos::operator=(const AliAnalysisTaskIDFragmentationFunction::AliFragFuncQATrackHistos& o)
1475 TObject::operator=(o);
1476 fNBinsPt = o.fNBinsPt;
1479 fNBinsEta = o.fNBinsEta;
1480 fEtaMin = o.fEtaMin;
1481 fEtaMax = o.fEtaMax;
1482 fNBinsPhi = o.fNBinsPhi;
1483 fPhiMin = o.fPhiMin;
1484 fPhiMax = o.fPhiMax;
1485 fHighPtThreshold = o.fHighPtThreshold;
1486 fh2EtaPhi = o.fh2EtaPhi;
1488 fh2HighPtEtaPhi = o.fh2HighPtEtaPhi;
1489 fh2PhiPt = o.fh2PhiPt;
1490 fNameQAT = o.fNameQAT;
1496 //___________________________________________________________________
1497 AliAnalysisTaskIDFragmentationFunction::AliFragFuncQATrackHistos::~AliFragFuncQATrackHistos()
1501 if(fh2EtaPhi) delete fh2EtaPhi;
1502 if(fh2HighPtEtaPhi) delete fh2HighPtEtaPhi;
1503 if(fh1Pt) delete fh1Pt;
1504 if(fh2PhiPt) delete fh2PhiPt;
1507 //______________________________________________________________________
1508 void AliAnalysisTaskIDFragmentationFunction::AliFragFuncQATrackHistos::DefineHistos()
1510 // book track QA histos
1512 fh2EtaPhi = new TH2F(Form("fh2TrackQAEtaPhi%s", fNameQAT.Data()), Form("%s: #eta - #phi distribution", fNameQAT.Data()), fNBinsEta, fEtaMin, fEtaMax, fNBinsPhi, fPhiMin, fPhiMax);
1513 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);
1514 fh1Pt = new TH1F(Form("fh1TrackQAPt%s", fNameQAT.Data()), Form("%s: p_{T} distribution", fNameQAT.Data()), fNBinsPt, fPtMin, fPtMax);
1515 fh2PhiPt = new TH2F(Form("fh2TrackQAPhiPt%s", fNameQAT.Data()), Form("%s: #phi - p_{T} distribution", fNameQAT.Data()), fNBinsPhi, fPhiMin, fPhiMax, fNBinsPt, fPtMin, fPtMax);
1517 AliAnalysisTaskIDFragmentationFunction::SetProperties(fh2EtaPhi, "#eta", "#phi");
1518 AliAnalysisTaskIDFragmentationFunction::SetProperties(fh2HighPtEtaPhi, "#eta", "#phi");
1519 AliAnalysisTaskIDFragmentationFunction::SetProperties(fh1Pt, "p_{T} [GeV/c]", "entries");
1520 AliAnalysisTaskIDFragmentationFunction::SetProperties(fh2PhiPt, "#phi", "p_{T} [GeV/c]");
1523 //________________________________________________________________________________________________________
1524 void AliAnalysisTaskIDFragmentationFunction::AliFragFuncQATrackHistos::FillTrackQA(Float_t eta, Float_t phi, Float_t pt, Bool_t weightPt, Float_t norm,
1525 Bool_t scaleStrangeness, Float_t scaleFacStrangeness)
1527 // fill track QA histos
1528 Float_t weight = 1.;
1529 if(weightPt) weight = pt;
1530 fh2EtaPhi->Fill( eta, phi, weight);
1531 if(scaleStrangeness) fh2EtaPhi->Fill( eta, phi, scaleFacStrangeness);
1532 if(pt > fHighPtThreshold) fh2HighPtEtaPhi->Fill( eta, phi, weight);
1533 if(pt > fHighPtThreshold && scaleStrangeness) fh2HighPtEtaPhi->Fill( eta, phi, weight);
1534 if(norm) fh1Pt->Fill( pt, 1/norm );
1535 else if(scaleStrangeness) fh1Pt->Fill(pt,scaleFacStrangeness);
1536 else fh1Pt->Fill( pt );
1538 if(scaleFacStrangeness) fh2PhiPt->Fill(phi, pt, scaleFacStrangeness);
1539 else fh2PhiPt->Fill(phi, pt);
1542 //______________________________________________________________________________________
1543 void AliAnalysisTaskIDFragmentationFunction::AliFragFuncQATrackHistos::AddToOutput(TList* list) const
1545 // add histos to list
1547 list->Add(fh2EtaPhi);
1548 list->Add(fh2HighPtEtaPhi);
1550 list->Add(fh2PhiPt);
1553 //_________________________________________________________________________________
1554 Bool_t AliAnalysisTaskIDFragmentationFunction::Notify()
1557 // Implemented Notify() to read the cross sections
1558 // and number of trials from pyxsec.root
1559 // (taken from AliAnalysisTaskJetSpectrum2)
1561 TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
1562 Float_t xsection = 0;
1563 Float_t ftrials = 1;
1567 TFile *curfile = tree->GetCurrentFile();
1569 Error("Notify","No current file");
1573 AliAnalysisHelperJetTasks::PythiaInfoFromFile(curfile->GetName(),xsection,ftrials);
1575 if (fUseInclusivePIDtask) {
1576 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++)
1577 fInclusivePIDtask[i]->FillXsec(xsection);
1580 if (fUseJetPIDtask) {
1581 for (Int_t i = 0; i < fNumJetPIDtasks; i++)
1582 fJetPIDtask[i]->FillXsec(xsection);
1585 if(!fh1Xsec||!fh1Trials){
1586 Printf("%s:%d No Histogram fh1Xsec",(char*)__FILE__,__LINE__);
1590 fh1Xsec->Fill("<#sigma>",xsection);
1591 // construct a poor man average trials
1592 Float_t nEntries = (Float_t)tree->GetTree()->GetEntries();
1593 if(ftrials>=nEntries && nEntries>0.)fAvgTrials = ftrials/nEntries;
1596 // Set seed for backg study
1598 fRandom = new TRandom3();
1599 fRandom->SetSeed(0);
1604 //__________________________________________________________________
1605 void AliAnalysisTaskIDFragmentationFunction::UserCreateOutputObjects()
1607 // create output objects
1609 if(fDebug > 1) Printf("AliAnalysisTaskIDFragmentationFunction::UserCreateOutputObjects()");
1611 // create list of tracks and jets
1613 fTracksRecCuts = new TList();
1614 fTracksRecCuts->SetOwner(kFALSE);
1616 fTracksRecCutsEfficiency = new TList();
1617 fTracksRecCutsEfficiency->SetOwner(kFALSE);
1619 fTracksGen = new TList();
1620 fTracksGen->SetOwner(kFALSE);
1622 fTracksAODMCCharged = new TList();
1623 fTracksAODMCCharged->SetOwner(kFALSE);
1625 fTracksAODMCChargedSecNS = new TList();
1626 fTracksAODMCChargedSecNS->SetOwner(kFALSE);
1628 fTracksAODMCChargedSecS = new TList();
1629 fTracksAODMCChargedSecS->SetOwner(kFALSE);
1631 fTracksRecQualityCuts = new TList();
1632 fTracksRecQualityCuts->SetOwner(kFALSE);
1634 fJetsRec = new TList();
1635 fJetsRec->SetOwner(kFALSE);
1637 fJetsRecCuts = new TList();
1638 fJetsRecCuts->SetOwner(kFALSE);
1640 fJetsGen = new TList();
1641 fJetsGen->SetOwner(kFALSE);
1643 fJetsRecEff = new TList();
1644 fJetsRecEff->SetOwner(kFALSE);
1646 fJetsEmbedded = new TList();
1647 fJetsEmbedded->SetOwner(kFALSE);
1651 (fBckgType[0]==kBckgClusters || fBckgType[1]==kBckgClusters || fBckgType[2]==kBckgClusters || fBckgType[3]==kBckgClusters || fBckgType[4]==kBckgClusters ||
1652 fBckgType[0]==kBckgClustersOutLeading || fBckgType[1]==kBckgClustersOutLeading || fBckgType[2]==kBckgClustersOutLeading ||
1653 fBckgType[3]==kBckgClustersOutLeading || fBckgType[4]==kBckgClustersOutLeading)){
1655 fBckgJetsRec = new TList();
1656 fBckgJetsRec->SetOwner(kFALSE);
1658 fBckgJetsRecCuts = new TList();
1659 fBckgJetsRecCuts->SetOwner(kFALSE);
1661 fBckgJetsGen = new TList();
1662 fBckgJetsGen->SetOwner(kFALSE);
1666 // Create histograms / output container
1670 fCommonHistList = new TList();
1671 fCommonHistList->SetOwner(kTRUE);
1673 Bool_t oldStatus = TH1::AddDirectoryStatus();
1674 TH1::AddDirectory(kFALSE);
1678 fh1EvtSelection = new TH1F("fh1EvtSelection", "Event Selection", 6, -0.5, 5.5);
1679 fh1EvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
1680 fh1EvtSelection->GetXaxis()->SetBinLabel(2,"event selection: rejected");
1681 fh1EvtSelection->GetXaxis()->SetBinLabel(3,"event class: rejected");
1682 fh1EvtSelection->GetXaxis()->SetBinLabel(4,"vertex Ncontr: rejected");
1683 fh1EvtSelection->GetXaxis()->SetBinLabel(5,"vertex z: rejected");
1684 fh1EvtSelection->GetXaxis()->SetBinLabel(6,"vertex type: rejected");
1686 fh1VtxSelection = new TH1F("fh1VtxSelection", "Vertex Selection", 10, -1, 9);
1687 fh1VtxSelection->GetXaxis()->SetBinLabel(1,"Undef");
1688 fh1VtxSelection->GetXaxis()->SetBinLabel(2,"Primary");
1689 fh1VtxSelection->GetXaxis()->SetBinLabel(3,"Kink");
1690 fh1VtxSelection->GetXaxis()->SetBinLabel(4,"V0");
1691 fh1VtxSelection->GetXaxis()->SetBinLabel(5,"Cascade");
1692 fh1VtxSelection->GetXaxis()->SetBinLabel(6,"Multi");
1693 fh1VtxSelection->GetXaxis()->SetBinLabel(7,"SPD");
1694 fh1VtxSelection->GetXaxis()->SetBinLabel(8,"PileUpSPD");
1695 fh1VtxSelection->GetXaxis()->SetBinLabel(9,"PileUpTracks");
1696 fh1VtxSelection->GetXaxis()->SetBinLabel(10,"TPC");
1698 fh1VertexNContributors = new TH1F("fh1VertexNContributors", "Vertex N contributors", 2500,-.5, 2499.5);
1699 fh1VertexZ = new TH1F("fh1VertexZ", "Vertex z distribution", 30, -15., 15.);
1700 fh1EvtMult = new TH1F("fh1EvtMult","Event multiplicity, track pT cut > 150 MeV/c, |#eta| < 0.9",120,0.,12000.);
1701 fh1EvtCent = new TH1F("fh1EvtCent","centrality",100,0.,100.);
1703 fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
1704 fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
1705 fh1Trials = new TH1F("fh1Trials","trials from pyxsec.root",1,0,1);
1706 fh1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
1707 fh1PtHard = new TH1F("fh1PtHard","PYTHIA Pt hard;p_{T,hard}",350,-.5,349.5);
1708 fh1PtHardTrials = new TH1F("fh1PtHardTrials","PYTHIA Pt hard weight with trials;p_{T,hard}",350,-.5,349.5);
1710 fh1EvtsPtHardCut = new TH1F("fh1EvtsPtHardCut", "#events before and after MC #it{p}_{T,hard} cut;;Events",2,0,2);
1711 fh1EvtsPtHardCut->GetXaxis()->SetBinLabel(1, "All");
1712 fh1EvtsPtHardCut->GetXaxis()->SetBinLabel(2, "#it{p}_{T,hard}");
1715 fh1nRecJetsCuts = new TH1F("fh1nRecJetsCuts","reconstructed jets per event",10,-0.5,9.5);
1716 fh1nGenJets = new TH1F("fh1nGenJets","generated jets per event",10,-0.5,9.5);
1717 fh1nRecEffJets = new TH1F("fh1nRecEffJets","reconstruction effiency: jets per event",10,-0.5,9.5);
1718 fh1nEmbeddedJets = new TH1F("fh1nEmbeddedJets","embedded jets per event",10,-0.5,9.5);
1720 fh2PtRecVsGenPrim = new TH2F("fh2PtRecVsGenPrim","rec vs gen pt",fQATrackNBinsPt,fQATrackPtMin,fQATrackPtMax,fQATrackNBinsPt,fQATrackPtMin,fQATrackPtMax);
1721 fh2PtRecVsGenSec = new TH2F("fh2PtRecVsGenSec","rec vs gen pt",fQATrackNBinsPt,fQATrackPtMin,fQATrackPtMax,fQATrackNBinsPt,fQATrackPtMin,fQATrackPtMax);
1724 if(fBranchEmbeddedJets.Length()){
1725 fh1FractionPtEmbedded = new TH1F("fh1FractionPtEmbedded","",200,0,2);
1726 fh1IndexEmbedded = new TH1F("fh1IndexEmbedded","",11,-1,10);
1727 fh2DeltaPtVsJetPtEmbedded = new TH2F("fh2DeltaPtVsJetPtEmbedded","",250,0,250,200,-100,100);
1728 fh2DeltaPtVsRecJetPtEmbedded = new TH2F("fh2DeltaPtVsRecJetPtEmbedded","",250,0,250,200,-100,100);
1729 fh1DeltaREmbedded = new TH1F("fh1DeltaREmbedded","",50,0,0.5);
1730 fh1nEmbeddedJets = new TH1F("fh1nEmbeddedJets","embedded jets per event",10,-0.5,9.5);
1735 if(fQAMode&1){ // track QA
1736 fQATrackHistosRecCuts = new AliFragFuncQATrackHistos("RecCuts", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1737 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1738 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1739 fQATrackHighPtThreshold);
1740 fQATrackHistosGen = new AliFragFuncQATrackHistos("Gen", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1741 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1742 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1743 fQATrackHighPtThreshold);
1746 if(fQAMode&2){ // jet QA
1747 fQAJetHistosRec = new AliFragFuncQAJetHistos("Rec", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax,
1748 fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,
1749 fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1750 fQAJetHistosRecCuts = new AliFragFuncQAJetHistos("RecCuts", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax,
1751 fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,
1752 fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1753 fQAJetHistosRecCutsLeading = new AliFragFuncQAJetHistos("RecCutsLeading", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax,
1754 fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,
1755 fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1756 fQAJetHistosGen = new AliFragFuncQAJetHistos("Gen", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax,
1757 fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,
1758 fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1759 fQAJetHistosGenLeading = new AliFragFuncQAJetHistos("GenLeading", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax,
1760 fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,
1761 fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1762 if(fEffMode) fQAJetHistosRecEffLeading = new AliFragFuncQAJetHistos("RecEffLeading", fQAJetNBinsPt, fQAJetPtMin, fQAJetPtMax,
1763 fQAJetNBinsEta, fQAJetEtaMin, fQAJetEtaMax,fQAJetNBinsPhi, fQAJetPhiMin, fQAJetPhiMax);
1767 if(fFFMode || fIDFFMode){
1769 fFFHistosRecCuts = new AliFragFuncHistos("RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1770 fFFNBinsPt, fFFPtMin, fFFPtMax,
1771 fFFNBinsXi, fFFXiMin, fFFXiMax,
1772 fFFNBinsZ , fFFZMin , fFFZMax );
1775 fFFHistosRecCutsInc = new AliFragFuncHistos("RecCutsInc", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1776 fFFNBinsPt, fFFPtMin, fFFPtMax,
1777 fFFNBinsXi, fFFXiMin, fFFXiMax,
1778 fFFNBinsZ , fFFZMin , fFFZMax );
1781 fFFHistosRecLeadingTrack = new AliFragFuncHistos("RecLeadingTrack", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1782 fFFNBinsPt, fFFPtMin, fFFPtMax,
1783 fFFNBinsXi, fFFXiMin, fFFXiMax,
1784 fFFNBinsZ , fFFZMin , fFFZMax );
1786 fFFHistosGen = new AliFragFuncHistos("Gen", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1787 fFFNBinsPt, fFFPtMin, fFFPtMax,
1788 fFFNBinsXi, fFFXiMin, fFFXiMax,
1789 fFFNBinsZ , fFFZMin , fFFZMax);
1791 fFFHistosGenInc = new AliFragFuncHistos("GenInc", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1792 fFFNBinsPt, fFFPtMin, fFFPtMax,
1793 fFFNBinsXi, fFFXiMin, fFFXiMax,
1794 fFFNBinsZ , fFFZMin , fFFZMax);
1796 fFFHistosGenLeadingTrack = new AliFragFuncHistos("GenLeadingTrack", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1797 fFFNBinsPt, fFFPtMin, fFFPtMax,
1798 fFFNBinsXi, fFFXiMin, fFFXiMax,
1799 fFFNBinsZ , fFFZMin , fFFZMax);
1802 for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
1803 fIDFFHistosRecCuts[i] = new AliFragFuncHistos(Form("RecCuts_%s", AliPID::ParticleName(i)), fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1804 fFFNBinsPt, fFFPtMin, fFFPtMax,
1805 fFFNBinsXi, fFFXiMin, fFFXiMax,
1806 fFFNBinsZ , fFFZMin , fFFZMax );
1807 fIDFFHistosGen[i] = new AliFragFuncHistos(Form("Gen_%s", AliPID::ParticleName(i)), fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1808 fFFNBinsPt, fFFPtMin, fFFPtMax,
1809 fFFNBinsXi, fFFXiMin, fFFXiMax,
1810 fFFNBinsZ , fFFZMin , fFFZMax );
1819 fQATrackHistosRecEffGen = new AliFragFuncQATrackHistos("RecEffGen", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1820 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1821 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1822 fQATrackHighPtThreshold);
1824 fQATrackHistosRecEffRec = new AliFragFuncQATrackHistos("RecEffRec", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1825 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1826 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1827 fQATrackHighPtThreshold);
1829 fQATrackHistosSecRecNS = new AliFragFuncQATrackHistos("SecRecNS", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1830 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1831 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1832 fQATrackHighPtThreshold);
1834 fQATrackHistosSecRecS = new AliFragFuncQATrackHistos("SecRecS", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1835 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1836 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1837 fQATrackHighPtThreshold);
1839 fQATrackHistosSecRecSsc = new AliFragFuncQATrackHistos("SecRecSsc", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1840 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1841 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1842 fQATrackHighPtThreshold);
1846 fFFHistosRecEffRec = new AliFragFuncHistos("RecEffRec", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1847 fFFNBinsPt, fFFPtMin, fFFPtMax,
1848 fFFNBinsXi, fFFXiMin, fFFXiMax,
1849 fFFNBinsZ , fFFZMin , fFFZMax);
1851 fFFHistosSecRecNS = new AliFragFuncHistos("SecRecNS", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1852 fFFNBinsPt, fFFPtMin, fFFPtMax,
1853 fFFNBinsXi, fFFXiMin, fFFXiMax,
1854 fFFNBinsZ , fFFZMin , fFFZMax);
1856 fFFHistosSecRecS = new AliFragFuncHistos("SecRecS", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1857 fFFNBinsPt, fFFPtMin, fFFPtMax,
1858 fFFNBinsXi, fFFXiMin, fFFXiMax,
1859 fFFNBinsZ , fFFZMin , fFFZMax);
1861 fFFHistosSecRecSsc = new AliFragFuncHistos("SecRecSsc", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1862 fFFNBinsPt, fFFPtMin, fFFPtMax,
1863 fFFNBinsXi, fFFXiMin, fFFXiMax,
1864 fFFNBinsZ , fFFZMin , fFFZMax);
1867 } // end: efficiency
1871 if(fBckgType[0]==kBckgNone){
1872 AliError("no bgr method selected !");
1876 for(Int_t i=0; i<5; i++){
1877 if(fBckgType[i]==kBckgPerp) title[i]="Perp";
1878 else if(fBckgType[i]==kBckgPerp2) title[i]="Perp2";
1879 else if(fBckgType[i]==kBckgPerp2Area) title[i]="Perp2Area";
1880 else if(fBckgType[i]==kBckgPerpWindow) title[i]="PerpW";
1881 else if(fBckgType[i]==kBckgASide) title[i]="ASide";
1882 else if(fBckgType[i]==kBckgASideWindow) title[i]="ASideW";
1883 else if(fBckgType[i]==kBckgOutLJ) title[i]="OutLeadingJet";
1884 else if(fBckgType[i]==kBckgOut2J) title[i]="Out2Jets";
1885 else if(fBckgType[i]==kBckgOut3J) title[i]="Out3Jets";
1886 else if(fBckgType[i]==kBckgOutAJ) title[i]="AllJets";
1887 else if(fBckgType[i]==kBckgOutLJStat) title[i]="OutLeadingJetStat";
1888 else if(fBckgType[i]==kBckgOut2JStat) title[i]="Out2JetsStat";
1889 else if(fBckgType[i]==kBckgOut3JStat) title[i]="Out3JetsStat";
1890 else if(fBckgType[i]==kBckgOutAJStat) title[i]="AllJetsStat";
1891 else if(fBckgType[i]==kBckgClustersOutLeading) title[i]="OutClusters";
1892 else if(fBckgType[i]==kBckgClusters) title[i]="MedianClusters";
1893 else if(fBckgType[i]==kBckgNone) title[i]="";
1894 else printf("Please chose background method number %d!",i);
1898 if(fBckgType[0]==kBckgClusters || fBckgType[1]==kBckgClusters || fBckgType[2]==kBckgClusters || fBckgType[3]==kBckgClusters || fBckgType[4]==kBckgClusters ||
1899 fBckgType[0]==kBckgClustersOutLeading || fBckgType[1]==kBckgClustersOutLeading || fBckgType[2]==kBckgClustersOutLeading ||
1900 fBckgType[3]==kBckgClustersOutLeading || fBckgType[4]==kBckgClustersOutLeading){
1902 fh1nRecBckgJetsCuts = new TH1F("fh1nRecBckgJetsCuts","reconstructed background jets per event",10,-0.5,9.5);
1903 fh1nGenBckgJets = new TH1F("fh1nGenBckgJets","generated background jets per event",10,-0.5,9.5);
1907 fh1BckgMult0 = new TH1F("fh1BckgMult0","bckg mult "+title[0],500,0,500);
1908 if(fBckgType[1] != kBckgNone) fh1BckgMult1 = new TH1F("fh1BckgMult1","bckg mult "+title[1],500,0,500);
1909 if(fBckgType[2] != kBckgNone) fh1BckgMult2 = new TH1F("fh1BckgMult2","bckg mult "+title[2],500,0,500);
1910 if(fBckgType[3] != kBckgNone) fh1BckgMult3 = new TH1F("fh1BckgMult3","bckg mult "+title[3],500,0,500);
1911 if(fBckgType[4] != kBckgNone) fh1BckgMult4 = new TH1F("fh1BckgMult4","bckg mult "+title[4],500,0,500);
1915 fQABckgHisto0RecCuts = new AliFragFuncQATrackHistos("Bckg"+title[0]+"RecCuts", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1916 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1917 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1918 fQATrackHighPtThreshold);
1919 fQABckgHisto0Gen = new AliFragFuncQATrackHistos("Bckg"+title[0]+"Gen", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1920 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1921 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1922 fQATrackHighPtThreshold);
1924 if(fBckgType[1] != kBckgNone){
1925 fQABckgHisto1RecCuts = new AliFragFuncQATrackHistos("Bckg"+title[1]+"RecCuts", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1926 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1927 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1928 fQATrackHighPtThreshold);
1929 fQABckgHisto1Gen = new AliFragFuncQATrackHistos("Bckg"+title[1]+"Gen", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1930 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1931 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1932 fQATrackHighPtThreshold);
1934 if(fBckgType[2] != kBckgNone){
1935 fQABckgHisto2RecCuts = new AliFragFuncQATrackHistos("Bckg"+title[2]+"RecCuts", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1936 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1937 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1938 fQATrackHighPtThreshold);
1939 fQABckgHisto2Gen = new AliFragFuncQATrackHistos("Bckg"+title[2]+"Gen", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1940 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1941 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1942 fQATrackHighPtThreshold);
1944 if(fBckgType[3] != kBckgNone){
1945 fQABckgHisto3RecCuts = new AliFragFuncQATrackHistos("Bckg"+title[3]+"RecCuts", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1946 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1947 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1948 fQATrackHighPtThreshold);
1949 fQABckgHisto3Gen = new AliFragFuncQATrackHistos("Bckg"+title[3]+"Gen", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1950 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1951 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1952 fQATrackHighPtThreshold);
1954 if(fBckgType[4] != kBckgNone){
1955 fQABckgHisto4RecCuts = new AliFragFuncQATrackHistos("Bckg"+title[4]+"RecCuts", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1956 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1957 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1958 fQATrackHighPtThreshold);
1959 fQABckgHisto4Gen = new AliFragFuncQATrackHistos("Bckg"+title[4]+"Gen", fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax,
1960 fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
1961 fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax,
1962 fQATrackHighPtThreshold);
1964 } // end: background QA
1967 fFFBckgHisto0RecCuts = new AliFragFuncHistos("Bckg"+title[0]+"RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1968 fFFNBinsPt, fFFPtMin, fFFPtMax,
1969 fFFNBinsXi, fFFXiMin, fFFXiMax,
1970 fFFNBinsZ , fFFZMin , fFFZMax);
1972 fFFBckgHisto0Gen = new AliFragFuncHistos("Bckg"+title[0]+"Gen", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1973 fFFNBinsPt, fFFPtMin, fFFPtMax,
1974 fFFNBinsXi, fFFXiMin, fFFXiMax,
1975 fFFNBinsZ , fFFZMin , fFFZMax);
1977 if(fBckgType[1] != kBckgNone){
1978 fFFBckgHisto1RecCuts = new AliFragFuncHistos("Bckg"+title[1]+"RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1979 fFFNBinsPt, fFFPtMin, fFFPtMax,
1980 fFFNBinsXi, fFFXiMin, fFFXiMax,
1981 fFFNBinsZ , fFFZMin , fFFZMax);
1982 fFFBckgHisto1Gen = new AliFragFuncHistos("Bckg"+title[1]+"Gen", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1983 fFFNBinsPt, fFFPtMin, fFFPtMax,
1984 fFFNBinsXi, fFFXiMin, fFFXiMax,
1985 fFFNBinsZ , fFFZMin , fFFZMax);
1987 if(fBckgType[2] != kBckgNone){
1988 fFFBckgHisto2RecCuts = new AliFragFuncHistos("Bckg"+title[2]+"RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1989 fFFNBinsPt, fFFPtMin, fFFPtMax,
1990 fFFNBinsXi, fFFXiMin, fFFXiMax,
1991 fFFNBinsZ , fFFZMin , fFFZMax);
1993 fFFBckgHisto2Gen = new AliFragFuncHistos("Bckg"+title[2]+"Gen", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
1994 fFFNBinsPt, fFFPtMin, fFFPtMax,
1995 fFFNBinsXi, fFFXiMin, fFFXiMax,
1996 fFFNBinsZ , fFFZMin , fFFZMax);
1998 if(fBckgType[3] != kBckgNone){
1999 fFFBckgHisto3RecCuts = new AliFragFuncHistos("Bckg"+title[3]+"RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
2000 fFFNBinsPt, fFFPtMin, fFFPtMax,
2001 fFFNBinsXi, fFFXiMin, fFFXiMax,
2002 fFFNBinsZ , fFFZMin , fFFZMax);
2004 fFFBckgHisto3Gen = new AliFragFuncHistos("Bckg"+title[3]+"Gen", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
2005 fFFNBinsPt, fFFPtMin, fFFPtMax,
2006 fFFNBinsXi, fFFXiMin, fFFXiMax,
2007 fFFNBinsZ , fFFZMin , fFFZMax);
2009 if(fBckgType[4] != kBckgNone){
2010 fFFBckgHisto4RecCuts = new AliFragFuncHistos("Bckg"+title[4]+"RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
2011 fFFNBinsPt, fFFPtMin, fFFPtMax,
2012 fFFNBinsXi, fFFXiMin, fFFXiMax,
2013 fFFNBinsZ , fFFZMin , fFFZMax);
2015 fFFBckgHisto4Gen = new AliFragFuncHistos("Bckg"+title[4]+"Gen", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
2016 fFFNBinsPt, fFFPtMin, fFFPtMax,
2017 fFFNBinsXi, fFFXiMin, fFFXiMax,
2018 fFFNBinsZ , fFFZMin , fFFZMax);
2021 fFFBckgHisto0RecEffRec = new AliFragFuncHistos("Bckg"+title[0]+"RecEffRec", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
2022 fFFNBinsPt, fFFPtMin, fFFPtMax,
2023 fFFNBinsXi, fFFXiMin, fFFXiMax,
2024 fFFNBinsZ , fFFZMin , fFFZMax);
2026 fFFBckgHisto0SecRecNS = new AliFragFuncHistos("Bckg"+title[0]+"SecRecNS", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
2027 fFFNBinsPt, fFFPtMin, fFFPtMax,
2028 fFFNBinsXi, fFFXiMin, fFFXiMax,
2029 fFFNBinsZ , fFFZMin , fFFZMax);
2031 fFFBckgHisto0SecRecS = new AliFragFuncHistos("Bckg"+title[0]+"SecRecS", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
2032 fFFNBinsPt, fFFPtMin, fFFPtMax,
2033 fFFNBinsXi, fFFXiMin, fFFXiMax,
2034 fFFNBinsZ , fFFZMin , fFFZMax);
2036 fFFBckgHisto0SecRecSsc = new AliFragFuncHistos("Bckg"+title[0]+"SecRecSsc", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax,
2037 fFFNBinsPt, fFFPtMin, fFFPtMax,
2038 fFFNBinsXi, fFFXiMin, fFFXiMax,
2039 fFFNBinsZ , fFFZMin , fFFZMax);
2042 } // end: background FF
2045 } // end: background
2048 // ____________ define histograms ____________________
2051 if(fQAMode&1){ // track QA
2052 fQATrackHistosRecCuts->DefineHistos();
2053 fQATrackHistosGen->DefineHistos();
2056 if(fQAMode&2){ // jet QA
2057 fQAJetHistosRec->DefineHistos();
2058 fQAJetHistosRecCuts->DefineHistos();
2059 fQAJetHistosRecCutsLeading->DefineHistos();
2060 fQAJetHistosGen->DefineHistos();
2061 fQAJetHistosGenLeading->DefineHistos();
2062 if(fEffMode) fQAJetHistosRecEffLeading->DefineHistos();
2066 if(fFFMode || fIDFFMode){
2067 fFFHistosRecCuts->DefineHistos();
2068 fFFHistosRecCutsInc->DefineHistos();
2069 fFFHistosRecLeadingTrack->DefineHistos();
2070 fFFHistosGen->DefineHistos();
2071 fFFHistosGenInc->DefineHistos();
2072 fFFHistosGenLeadingTrack->DefineHistos();
2075 for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
2076 fIDFFHistosRecCuts[i]->DefineHistos();
2077 fIDFFHistosGen[i]->DefineHistos();
2084 fQATrackHistosRecEffGen->DefineHistos();
2085 fQATrackHistosRecEffRec->DefineHistos();
2086 fQATrackHistosSecRecNS->DefineHistos();
2087 fQATrackHistosSecRecS->DefineHistos();
2088 fQATrackHistosSecRecSsc->DefineHistos();
2091 fFFHistosRecEffRec->DefineHistos();
2092 fFFHistosSecRecNS->DefineHistos();
2093 fFFHistosSecRecS->DefineHistos();
2094 fFFHistosSecRecSsc->DefineHistos();
2096 } // end: efficiency
2101 fFFBckgHisto0RecCuts->DefineHistos();
2102 fFFBckgHisto0Gen->DefineHistos();
2103 if(fBckgType[1] != kBckgNone) fFFBckgHisto1RecCuts->DefineHistos();
2104 if(fBckgType[1] != kBckgNone) fFFBckgHisto1Gen->DefineHistos();
2105 if(fBckgType[2] != kBckgNone) fFFBckgHisto2RecCuts->DefineHistos();
2106 if(fBckgType[2] != kBckgNone) fFFBckgHisto2Gen->DefineHistos();
2107 if(fBckgType[3] != kBckgNone) fFFBckgHisto3RecCuts->DefineHistos();
2108 if(fBckgType[3] != kBckgNone) fFFBckgHisto3Gen->DefineHistos();
2109 if(fBckgType[4] != kBckgNone) fFFBckgHisto4RecCuts->DefineHistos();
2110 if(fBckgType[4] != kBckgNone) fFFBckgHisto4Gen->DefineHistos();
2113 fFFBckgHisto0RecEffRec->DefineHistos();
2114 fFFBckgHisto0SecRecNS->DefineHistos();
2115 fFFBckgHisto0SecRecS->DefineHistos();
2116 fFFBckgHisto0SecRecSsc->DefineHistos();
2121 fQABckgHisto0RecCuts->DefineHistos();
2122 fQABckgHisto0Gen->DefineHistos();
2123 if(fBckgType[1] != kBckgNone) fQABckgHisto1RecCuts->DefineHistos();
2124 if(fBckgType[1] != kBckgNone) fQABckgHisto1Gen->DefineHistos();
2125 if(fBckgType[2] != kBckgNone) fQABckgHisto2RecCuts->DefineHistos();
2126 if(fBckgType[2] != kBckgNone) fQABckgHisto2Gen->DefineHistos();
2127 if(fBckgType[3] != kBckgNone) fQABckgHisto3RecCuts->DefineHistos();
2128 if(fBckgType[3] != kBckgNone) fQABckgHisto3Gen->DefineHistos();
2129 if(fBckgType[4] != kBckgNone) fQABckgHisto4RecCuts->DefineHistos();
2130 if(fBckgType[4] != kBckgNone) fQABckgHisto4Gen->DefineHistos();
2132 } // end: background
2135 Bool_t genJets = (fJetTypeGen != kJetsUndef) ? kTRUE : kFALSE;
2136 Bool_t genTracks = (fTrackTypeGen != kTrackUndef) ? kTRUE : kFALSE;
2137 Bool_t recJetsEff = (fJetTypeRecEff != kJetsUndef) ? kTRUE : kFALSE;
2139 fCommonHistList->Add(fh1EvtSelection);
2140 fCommonHistList->Add(fh1VtxSelection);
2141 fCommonHistList->Add(fh1EvtMult);
2142 fCommonHistList->Add(fh1EvtCent);
2143 fCommonHistList->Add(fh1VertexNContributors);
2144 fCommonHistList->Add(fh1VertexZ);
2145 fCommonHistList->Add(fh1nRecJetsCuts);
2146 fCommonHistList->Add(fh1Xsec);
2147 fCommonHistList->Add(fh1Trials);
2148 fCommonHistList->Add(fh1PtHard);
2149 fCommonHistList->Add(fh1PtHardTrials);
2150 fCommonHistList->Add(fh1EvtsPtHardCut);
2152 if(genJets) fCommonHistList->Add(fh1nGenJets);
2156 fFFHistosRecCuts->AddToOutput(fCommonHistList);
2157 fFFHistosRecCutsInc->AddToOutput(fCommonHistList);
2158 fFFHistosRecLeadingTrack->AddToOutput(fCommonHistList);
2160 if(genJets && genTracks){
2161 fFFHistosGen->AddToOutput(fCommonHistList);
2162 fFFHistosGenInc->AddToOutput(fCommonHistList);
2163 fFFHistosGenLeadingTrack->AddToOutput(fCommonHistList);
2167 for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
2168 if(genJets && genTracks)
2169 fIDFFHistosGen[i]->AddToOutput(fCommonHistList);
2170 fIDFFHistosRecCuts[i]->AddToOutput(fCommonHistList);
2178 fFFBckgHisto0RecCuts->AddToOutput(fCommonHistList);
2179 if(fBckgType[1] != kBckgNone) fFFBckgHisto1RecCuts->AddToOutput(fCommonHistList);
2180 if(fBckgType[2] != kBckgNone) fFFBckgHisto2RecCuts->AddToOutput(fCommonHistList);
2181 if(fBckgType[3] != kBckgNone) fFFBckgHisto3RecCuts->AddToOutput(fCommonHistList);
2182 if(fBckgType[4] != kBckgNone) fFFBckgHisto4RecCuts->AddToOutput(fCommonHistList);
2184 if(genJets && genTracks){
2185 fFFBckgHisto0Gen->AddToOutput(fCommonHistList);
2186 if(fBckgType[1] != kBckgNone) fFFBckgHisto1Gen->AddToOutput(fCommonHistList);
2187 if(fBckgType[2] != kBckgNone) fFFBckgHisto2Gen->AddToOutput(fCommonHistList);
2188 if(fBckgType[3] != kBckgNone) fFFBckgHisto3Gen->AddToOutput(fCommonHistList);
2189 if(fBckgType[4] != kBckgNone) fFFBckgHisto4Gen->AddToOutput(fCommonHistList);
2193 fFFBckgHisto0RecEffRec->AddToOutput(fCommonHistList);
2194 fFFBckgHisto0SecRecNS->AddToOutput(fCommonHistList);
2195 fFFBckgHisto0SecRecS->AddToOutput(fCommonHistList);
2196 fFFBckgHisto0SecRecSsc->AddToOutput(fCommonHistList);
2201 fQABckgHisto0RecCuts->AddToOutput(fCommonHistList);
2202 if(fBckgType[1] != kBckgNone) fQABckgHisto1RecCuts->AddToOutput(fCommonHistList);
2203 if(fBckgType[2] != kBckgNone) fQABckgHisto2RecCuts->AddToOutput(fCommonHistList);
2204 if(fBckgType[3] != kBckgNone) fQABckgHisto3RecCuts->AddToOutput(fCommonHistList);
2205 if(fBckgType[4] != kBckgNone) fQABckgHisto4RecCuts->AddToOutput(fCommonHistList);
2206 if(genJets && genTracks){
2207 fQABckgHisto0Gen->AddToOutput(fCommonHistList);
2208 if(fBckgType[1] != kBckgNone) fQABckgHisto1Gen->AddToOutput(fCommonHistList);
2209 if(fBckgType[2] != kBckgNone) fQABckgHisto2Gen->AddToOutput(fCommonHistList);
2210 if(fBckgType[3] != kBckgNone) fQABckgHisto3Gen->AddToOutput(fCommonHistList);
2211 if(fBckgType[4] != kBckgNone) fQABckgHisto4Gen->AddToOutput(fCommonHistList);
2215 if(fh1BckgMult0) fCommonHistList->Add(fh1BckgMult0);
2216 if(fBckgType[1] != kBckgNone) fCommonHistList->Add(fh1BckgMult1);
2217 if(fBckgType[2] != kBckgNone) fCommonHistList->Add(fh1BckgMult2);
2218 if(fBckgType[3] != kBckgNone) fCommonHistList->Add(fh1BckgMult3);
2219 if(fBckgType[4] != kBckgNone) fCommonHistList->Add(fh1BckgMult4);
2223 if(fBranchEmbeddedJets.Length()){
2224 fCommonHistList->Add(fh1FractionPtEmbedded);
2225 fCommonHistList->Add(fh1IndexEmbedded);
2226 fCommonHistList->Add(fh2DeltaPtVsJetPtEmbedded);
2227 fCommonHistList->Add(fh2DeltaPtVsRecJetPtEmbedded);
2228 fCommonHistList->Add(fh1DeltaREmbedded);
2229 fCommonHistList->Add(fh1nEmbeddedJets);
2235 if(fQAMode&1){ // track QA
2236 fQATrackHistosRecCuts->AddToOutput(fCommonHistList);
2237 if(genTracks) fQATrackHistosGen->AddToOutput(fCommonHistList);
2240 if(fQAMode&2){ // jet QA
2241 fQAJetHistosRec->AddToOutput(fCommonHistList);
2242 fQAJetHistosRecCuts->AddToOutput(fCommonHistList);
2243 fQAJetHistosRecCutsLeading->AddToOutput(fCommonHistList);
2244 if(recJetsEff && fEffMode) fQAJetHistosRecEffLeading->AddToOutput(fCommonHistList);
2246 fQAJetHistosGen->AddToOutput(fCommonHistList);
2247 fQAJetHistosGenLeading->AddToOutput(fCommonHistList);
2253 (fBckgType[0]==kBckgClusters || fBckgType[1]==kBckgClusters || fBckgType[2]==kBckgClusters || fBckgType[3]==kBckgClusters || fBckgType[4]==kBckgClusters ||
2254 fBckgType[0]==kBckgClustersOutLeading || fBckgType[1]==kBckgClustersOutLeading || fBckgType[2]==kBckgClustersOutLeading ||
2255 fBckgType[3]==kBckgClustersOutLeading || fBckgType[4]==kBckgClustersOutLeading)) {
2256 fCommonHistList->Add(fh1nRecBckgJetsCuts);
2257 if(genJets) fCommonHistList->Add(fh1nGenBckgJets);
2261 if(fEffMode && recJetsEff && genTracks){
2263 fQATrackHistosRecEffGen->AddToOutput(fCommonHistList);
2264 fQATrackHistosRecEffRec->AddToOutput(fCommonHistList);
2265 fQATrackHistosSecRecNS->AddToOutput(fCommonHistList);
2266 fQATrackHistosSecRecS->AddToOutput(fCommonHistList);
2267 fQATrackHistosSecRecSsc->AddToOutput(fCommonHistList);
2270 fFFHistosRecEffRec->AddToOutput(fCommonHistList);
2271 fFFHistosSecRecNS->AddToOutput(fCommonHistList);
2272 fFFHistosSecRecS->AddToOutput(fCommonHistList);
2273 fFFHistosSecRecSsc->AddToOutput(fCommonHistList);
2275 fCommonHistList->Add(fh1nRecEffJets);
2276 fCommonHistList->Add(fh2PtRecVsGenPrim);
2277 fCommonHistList->Add(fh2PtRecVsGenSec);
2283 fProNtracksLeadingJet = new TProfile("AvgNoOfTracksLeadingJet","AvgNoOfTracksLeadingJet",100,0,250,0,50);
2284 fProDelR80pcPt = new TProfile("AvgdelR80pcPt","AvgdelR80pcPt",100,0,250,0,1);
2286 if(genJets && genTracks){
2287 fProNtracksLeadingJetGen = new TProfile("AvgNoOfTracksLeadingJetGen","AvgNoOfTracksLeadingJetGen",100,0,250,0,50);
2288 fProDelR80pcPtGen = new TProfile("AvgdelR80pcPtGen","AvgdelR80pcPtGen",100,0,250,0,1);
2292 fProNtracksLeadingJetBgrPerp2 = new TProfile("AvgNoOfTracksLeadingJetBgrPerp2","AvgNoOfTracksLeadingJetBgrPerp2",100,0,250,0,50);
2295 fProNtracksLeadingJetRecPrim = new TProfile("AvgNoOfTracksLeadingJetRecPrim","AvgNoOfTracksLeadingJetRecPrim",100,0,250,0,50);
2296 fProDelR80pcPtRecPrim = new TProfile("AvgdelR80pcPtRecPrim","AvgdelR80pcPtRecPrim",100,0,250,0,1);
2297 fProNtracksLeadingJetRecSecNS = new TProfile("AvgNoOfTracksLeadingJetRecSecNS","AvgNoOfTracksLeadingJetRecSecNS",100,0,250,0,50);
2298 fProNtracksLeadingJetRecSecS = new TProfile("AvgNoOfTracksLeadingJetRecSecS","AvgNoOfTracksLeadingJetRecSecS",100,0,250,0,50);
2299 fProNtracksLeadingJetRecSecSsc = new TProfile("AvgNoOfTracksLeadingJetRecSecSsc","AvgNoOfTracksLeadingJetRecSecSsc",100,0,250,0,50);
2303 for(Int_t ii=0; ii<5; ii++){
2304 if(ii==0)strTitJS = "_JetPt20to30";
2305 if(ii==1)strTitJS = "_JetPt30to40";
2306 if(ii==2)strTitJS = "_JetPt40to60";
2307 if(ii==3)strTitJS = "_JetPt60to80";
2308 if(ii==4)strTitJS = "_JetPt80to100";
2310 fProDelRPtSum[ii] = new TProfile(Form("AvgPtSumDelR%s",strTitJS.Data()),Form("AvgPtSumDelR%s",strTitJS.Data()),50,0,1,0,250);
2311 if(genJets && genTracks)
2312 fProDelRPtSumGen[ii] = new TProfile(Form("AvgPtSumDelRGen%s",strTitJS.Data()),Form("AvgPtSumDelRGen%s",strTitJS.Data()),50,0,1,0,250);
2314 fProDelRPtSumBgrPerp2[ii] = new TProfile(Form("AvgPtSumDelRBgrPerp2%s",strTitJS.Data()),Form("AvgPtSumDelRBgrPerp2%s",strTitJS.Data()),50,0,1,0,250);
2316 fProDelRPtSumRecPrim[ii] = new TProfile(Form("AvgPtSumDelRRecPrim%s",strTitJS.Data()),Form("AvgPtSumDelRRecPrim%s",strTitJS.Data()),50,0,1,0,250);
2317 fProDelRPtSumRecSecNS[ii] = new TProfile(Form("AvgPtSumDelRRecSecNS%s",strTitJS.Data()),Form("AvgPtSumDelRRecSecNS%s",strTitJS.Data()),50,0,1,0,250);
2318 fProDelRPtSumRecSecS[ii] = new TProfile(Form("AvgPtSumDelRRecSecS%s",strTitJS.Data()),Form("AvgPtSumDelRRecSecS%s",strTitJS.Data()),50,0,1,0,250);
2319 fProDelRPtSumRecSecSsc[ii] = new TProfile(Form("AvgPtSumDelRRecSecSsc%s",strTitJS.Data()),Form("AvgPtSumDelRRecSecSsc%s",strTitJS.Data()),50,0,1,0,250);
2323 fCommonHistList->Add(fProNtracksLeadingJet);
2324 fCommonHistList->Add(fProDelR80pcPt);
2325 for(int ii=0; ii<5; ii++) fCommonHistList->Add(fProDelRPtSum[ii]);
2327 if(genJets && genTracks){
2328 fCommonHistList->Add(fProNtracksLeadingJetGen);
2329 fCommonHistList->Add(fProDelR80pcPtGen);
2330 for(Int_t ii=0; ii<5; ii++) fCommonHistList->Add(fProDelRPtSumGen[ii]);
2334 fCommonHistList->Add(fProNtracksLeadingJetBgrPerp2);
2335 for(Int_t ii=0; ii<5; ii++) fCommonHistList->Add(fProDelRPtSumBgrPerp2[ii]);
2339 fCommonHistList->Add(fProNtracksLeadingJetRecPrim);
2340 fCommonHistList->Add(fProDelR80pcPtRecPrim);
2341 for(Int_t ii=0; ii<5; ii++) fCommonHistList->Add(fProDelRPtSumRecPrim[ii]);
2343 fCommonHistList->Add(fProNtracksLeadingJetRecSecNS);
2344 for(Int_t ii=0; ii<5; ii++) fCommonHistList->Add(fProDelRPtSumRecSecNS[ii]);
2346 fCommonHistList->Add(fProNtracksLeadingJetRecSecS);
2347 for(Int_t ii=0; ii<5; ii++) fCommonHistList->Add(fProDelRPtSumRecSecS[ii]);
2349 fCommonHistList->Add(fProNtracksLeadingJetRecSecSsc);
2350 for(Int_t ii=0; ii<5; ii++) fCommonHistList->Add(fProDelRPtSumRecSecSsc[ii]);
2354 // Default analysis utils
2355 fAnaUtils = new AliAnalysisUtils();
2357 // Not used yet, but to be save, forward vertex z cut to analysis utils object
2358 fAnaUtils->SetMaxVtxZ(fMaxVertexZ);
2360 // Load PID framework if desired
2361 if(fDebug > 1) Printf("AliAnalysisTaskIDFragmentationFunction::UserCreateOutputObjects() -> Loading PID framework");
2363 fUseJetPIDtask = fIDFFMode || fFFMode;
2364 fUseInclusivePIDtask = fQAMode && (fQAMode&1);
2366 if (fUseJetPIDtask || fUseInclusivePIDtask) {
2367 TObjArray* tasks = AliAnalysisManager::GetAnalysisManager()->GetTasks();
2369 Printf("ERROR loading PID tasks: Failed to retrieve tasks from analysis manager!\n");
2371 fUseJetPIDtask = kFALSE;
2372 fUseInclusivePIDtask = kFALSE;
2375 if (fUseJetPIDtask) {
2376 delete [] fJetPIDtask;
2379 if (fNumJetPIDtasks > 0) {
2380 fJetPIDtask = new AliAnalysisTaskPID*[fNumJetPIDtasks];
2382 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
2383 fJetPIDtask[i] = (AliAnalysisTaskPID*)tasks->FindObject(fNameJetPIDtask[i].Data());
2385 if (!fJetPIDtask[i]) {
2386 Printf("ERROR: Failed to load jet pid task!\n");
2387 fUseJetPIDtask = kFALSE;
2392 Printf("WARNING: zero jet pid tasks!\n");
2393 fUseJetPIDtask = kFALSE;
2397 if (fUseInclusivePIDtask) {
2398 delete [] fInclusivePIDtask;
2399 fInclusivePIDtask = 0x0;
2401 if (fNumInclusivePIDtasks > 0) {
2402 fInclusivePIDtask = new AliAnalysisTaskPID*[fNumInclusivePIDtasks];
2404 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
2405 fInclusivePIDtask[i] = (AliAnalysisTaskPID*)tasks->FindObject(fNameInclusivePIDtask[i].Data());
2407 if (!fInclusivePIDtask[i]) {
2408 Printf("ERROR: Failed to load inclusive pid task!\n");
2409 fUseInclusivePIDtask = kFALSE;
2414 Printf("WARNING: zero inclusive pid tasks!\n");
2415 fUseInclusivePIDtask = kFALSE;
2420 const Int_t nRefMultBins = 100;
2421 const Double_t refMultUp = 100.;
2422 const Double_t refMultDown = 0.;
2424 const Int_t nJetPtBins = 20;
2425 const Double_t jetPtUp = 100.;
2426 const Double_t jetPtDown = 0.;
2428 const Int_t nCentBins = 12;
2429 const Double_t binsCentpp[nCentBins+1] = { 0, 0.01, 0.1, 1, 5, 10, 15, 20, 30, 40, 50, 70, 100};
2431 fhJetPtRefMultEta5 = new TH2F("fhJetPtRefMultEta5",
2432 "Correlation between jet energy and event multiplicity (|#eta| < 0.5);Ref. mult. (|#eta| < 0.5);#it{p}_{T, jet}^{ch} (GeV/#it{c})",
2433 nRefMultBins, refMultDown, refMultUp, nJetPtBins, jetPtDown, jetPtUp);
2434 fhJetPtRefMultEta8 = new TH2F("fhJetPtRefMultEta8",
2435 "Correlation between jet energy and event multiplicity (|#eta| < 0.8);Ref. mult. (|#eta| < 0.8);#it{p}_{T, jet}^{ch} (GeV/#it{c})",
2436 nRefMultBins, refMultDown, refMultUp, nJetPtBins, jetPtDown, jetPtUp);
2437 fhJetPtMultPercent = new TH2F("fhJetPtMultPercent",
2438 "Correlation between jet energy and event multiplicity percentile (V0M);Multiplicity Percentile (V0M);#it{p}_{T, jet}^{ch} (GeV/#it{c})",
2439 nCentBins, binsCentpp, nJetPtBins, jetPtDown, jetPtUp);
2440 fCommonHistList->Add(fhJetPtRefMultEta5);
2441 fCommonHistList->Add(fhJetPtRefMultEta8);
2442 fCommonHistList->Add(fhJetPtMultPercent);
2444 if (fUseJetPIDtask) {
2445 const Int_t nPtBins = 68;
2446 Double_t binsPt[nPtBins+1] = {0. , 0.05, 0.1, 0.15, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45,
2447 0.5, 0.55, 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95,
2448 1.0, 1.1 , 1.2, 1.3 , 1.4, 1.5 , 1.6, 1.7 , 1.8, 1.9 ,
2449 2.0, 2.2 , 2.4, 2.6 , 2.8, 3.0 , 3.2, 3.4 , 3.6, 3.8 ,
2450 4.0, 4.5 , 5.0, 5.5 , 6.0, 6.5 , 7.0, 8.0 , 9.0, 10.0,
2451 11.0, 12.0, 13.0, 14.0, 15.0, 16.0, 18.0, 20.0, 22.0, 24.0,
2452 26.0, 28.0, 30.0, 32.0, 34.0, 36.0, 40.0, 45.0, 50.0 };
2454 const Int_t DCAbins = 320;
2455 const Double_t DCA_Z_max = 3.5;
2456 const Double_t DCA_XY_max = 2.5;
2458 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);
2459 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);
2460 fCommonHistList->Add(fhDCA_XY);
2461 fCommonHistList->Add(fhDCA_Z);
2463 for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
2464 fhDCA_XY_prim_MCID[i] = new TH2F(Form("fhDCA_XY_prim_MCID_%s", AliPID::ParticleShortName(i)),
2465 Form("Rec. %s (prim.);#it{p}_{T} (GeV/#it{c});DCA_{XY}", AliPID::ParticleLatexName(i)),
2466 nPtBins, binsPt, DCAbins, -DCA_XY_max, DCA_XY_max);
2467 fhDCA_Z_prim_MCID[i] = new TH2F(Form("fhDCA_Z_prim_MCID_%s", AliPID::ParticleShortName(i)),
2468 Form("Rec. %s (prim.);#it{p}_{T} (GeV/#it{c});DCA_{Z}", AliPID::ParticleLatexName(i)),
2469 nPtBins, binsPt, DCAbins, -DCA_Z_max, DCA_Z_max);
2470 fCommonHistList->Add(fhDCA_XY_prim_MCID[i]);
2471 fCommonHistList->Add(fhDCA_Z_prim_MCID[i]);
2473 fhDCA_XY_sec_MCID[i] = new TH2F(Form("fhDCA_XY_sec_MCID_%s", AliPID::ParticleShortName(i)),
2474 Form("Rec. %s (sec.);#it{p}_{T} (GeV/#it{c});DCA_{XY}", AliPID::ParticleLatexName(i)),
2475 nPtBins, binsPt, DCAbins, -DCA_XY_max, DCA_XY_max);
2476 fhDCA_Z_sec_MCID[i] = new TH2F(Form("fhDCA_Z_sec_MCID_%s", AliPID::ParticleShortName(i)),
2477 Form("Rec. %s (sec.);#it{p}_{T} (GeV/#it{c});DCA_{Z}", AliPID::ParticleLatexName(i)),
2478 nPtBins, binsPt, DCAbins, -DCA_Z_max, DCA_Z_max);
2479 fCommonHistList->Add(fhDCA_XY_sec_MCID[i]);
2480 fCommonHistList->Add(fhDCA_Z_sec_MCID[i]);
2485 // =========== Switch on Sumw2 for all histos ===========
2486 for (Int_t i=0; i<fCommonHistList->GetEntries(); ++i){
2487 TH1 *h1 = dynamic_cast<TH1*>(fCommonHistList->At(i));
2488 if (h1) h1->Sumw2();
2490 THnSparse *hnSparse = dynamic_cast<THnSparse*>(fCommonHistList->At(i));
2491 if(hnSparse) hnSparse->Sumw2();
2495 TH1::AddDirectory(oldStatus);
2497 if(fDebug > 2) Printf("AliAnalysisTaskIDFragmentationFunction::UserCreateOutputObjects() -> Posting Output");
2499 PostData(1, fCommonHistList);
2501 if(fDebug > 2) Printf("AliAnalysisTaskIDFragmentationFunction::UserCreateOutputObjects() -> Done");
2504 //_______________________________________________
2505 void AliAnalysisTaskIDFragmentationFunction::Init()
2508 if(fDebug > 1) Printf("AliAnalysisTaskIDFragmentationFunction::Init()");
2512 //_____________________________________________________________
2513 void AliAnalysisTaskIDFragmentationFunction::UserExec(Option_t *)
2516 // Called for each event
2518 if(fDebug > 1) Printf("AliAnalysisTaskIDFragmentationFunction::UserExec()");
2521 if(fDebug > 1) Printf("Analysis event #%5d", (Int_t) fEntry);
2523 fMCEvent = MCEvent();
2525 if(fDebug>3) Printf("%s:%d MCEvent not found in the input", (char*)__FILE__,__LINE__);
2528 // Extract pThard and nTrials in case of MC.
2530 Double_t ptHard = 0.;
2531 Double_t nTrials = 1; // trials for MC trigger weight for real data
2532 Bool_t pythiaGenHeaderFound = kFALSE;
2535 AliGenEventHeader* genHeader = fMCEvent->GenEventHeader();
2538 AliGenPythiaEventHeader* pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
2539 AliGenHijingEventHeader* hijingGenHeader = 0x0;
2541 if(pythiaGenHeader){
2542 if(fDebug>3) Printf("%s:%d pythiaGenHeader found", (char*)__FILE__,__LINE__);
2543 pythiaGenHeaderFound = kTRUE;
2544 nTrials = pythiaGenHeader->Trials();
2545 ptHard = pythiaGenHeader->GetPtHard();
2546 } else { // no pythia, hijing?
2548 if(fDebug>3) Printf("%s:%d no pythiaGenHeader found", (char*)__FILE__,__LINE__);
2550 hijingGenHeader = dynamic_cast<AliGenHijingEventHeader*>(genHeader);
2551 if(!hijingGenHeader){
2552 Printf("%s:%d no pythiaGenHeader or hjingGenHeader found", (char*)__FILE__,__LINE__);
2554 if(fDebug>3) Printf("%s:%d hijingGenHeader found", (char*)__FILE__,__LINE__);
2558 //fh1Trials->Fill("#sum{ntrials}",fAvgTrials);
2563 // Cut on pThard if fMCEvent and pThard >= 0 and fill histo with #evt before and after the cut
2564 if (pythiaGenHeaderFound) {
2566 fh1EvtsPtHardCut->Fill(0.);
2568 if (fUseJetPIDtask) {
2569 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
2570 fJetPIDtask[i]->FillCutHisto(0., AliAnalysisTaskPID::kMCPtHardCut);
2574 if (fUseInclusivePIDtask) {
2575 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
2576 fInclusivePIDtask[i]->FillCutHisto(0., AliAnalysisTaskPID::kMCPtHardCut);
2582 if (fMCPtHardCut >= 0. && ptHard >= fMCPtHardCut) {
2583 if (fDebug>3) Printf("%s:%d skipping event with pThard %f (>= %f)", (char*)__FILE__,__LINE__, ptHard, fMCPtHardCut);
2584 PostData(1, fCommonHistList);
2589 fh1EvtsPtHardCut->Fill(1.);
2591 if (fUseJetPIDtask) {
2592 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
2593 fJetPIDtask[i]->FillCutHisto(1., AliAnalysisTaskPID::kMCPtHardCut);
2597 if (fUseInclusivePIDtask) {
2598 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
2599 fInclusivePIDtask[i]->FillCutHisto(1., AliAnalysisTaskPID::kMCPtHardCut);
2604 // Trigger selection
2605 AliInputEventHandler* inputHandler = (AliInputEventHandler*)
2606 ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
2608 if(!(inputHandler->IsEventSelected() & fEvtSelectionMask)){
2609 fh1EvtSelection->Fill(1.);
2610 if (fDebug > 1 ) Printf(" Trigger Selection: event REJECTED ... ");
2611 PostData(1, fCommonHistList);
2615 fESD = dynamic_cast<AliESDEvent*>(InputEvent());
2617 if(fDebug>3) Printf("%s:%d ESDEvent not found in the input", (char*)__FILE__,__LINE__);
2620 // get AOD event from input/ouput
2621 TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
2622 if( handler && handler->InheritsFrom("AliAODInputHandler") ) {
2623 fAOD = ((AliAODInputHandler*)handler)->GetEvent();
2624 if(fUseAODInputJets) fAODJets = fAOD;
2625 if (fDebug > 1) Printf("%s:%d AOD event from input", (char*)__FILE__,__LINE__);
2628 handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
2629 if( handler && handler->InheritsFrom("AliAODHandler") ) {
2630 fAOD = ((AliAODHandler*)handler)->GetAOD();
2632 if (fDebug > 1) Printf("%s:%d AOD event from output", (char*)__FILE__,__LINE__);
2636 if(!fAODJets && !fUseAODInputJets){ // case we have AOD in input & output and want jets from output
2637 TObject* outHandler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
2638 if( outHandler && outHandler->InheritsFrom("AliAODHandler") ) {
2639 fAODJets = ((AliAODHandler*)outHandler)->GetAOD();
2640 if (fDebug > 1) Printf("%s:%d jets from output AOD", (char*)__FILE__,__LINE__);
2644 if(fNonStdFile.Length()!=0){
2645 // case we have an AOD extension - fetch the jets from the extended output
2647 AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
2648 fAODExtension = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);
2650 if(fDebug>1)Printf("AODExtension not found for %s",fNonStdFile.Data());
2655 Printf("%s:%d AODEvent not found", (char*)__FILE__,__LINE__);
2659 Printf("%s:%d AODEvent with jet branch not found", (char*)__FILE__,__LINE__);
2664 // event selection **************************************************
2665 // *** event class ***
2666 AliVEvent* evtForCentDetermination = handler->InheritsFrom("AliAODInputHandler") ? fAOD : InputEvent();
2668 Double_t centPercent = -1;
2671 if(handler->InheritsFrom("AliAODInputHandler")){
2672 // since it is not supported by the helper task define own classes
2673 centPercent = ((AliAODHeader*)fAOD->GetHeader())->GetCentrality();
2675 if(centPercent>10) cl = 2;
2676 if(centPercent>30) cl = 3;
2677 if(centPercent>50) cl = 4;
2680 cl = AliAnalysisHelperJetTasks::EventClass();
2681 if(fESD) centPercent = fESD->GetCentrality()->GetCentralityPercentile("V0M"); // retrieve value 'by hand'
2684 if(cl!=fEventClass){
2685 // event not in selected event class, reject event
2686 if (fDebug > 1) Printf("%s:%d event not in selected event class: event REJECTED ...",(char*)__FILE__,__LINE__);
2687 fh1EvtSelection->Fill(2.);
2688 PostData(1, fCommonHistList);
2693 // Retrieve reference multiplicities in |eta|<0.8 and <0.5
2694 const Int_t refMult5 = ((AliAODHeader*)fAOD->GetHeader())->GetRefMultiplicityComb05();
2695 const Int_t refMult8 = ((AliAODHeader*)fAOD->GetHeader())->GetRefMultiplicityComb08();
2696 const Double_t centPercentPP = fAnaUtils->GetMultiplicityPercentile(fAOD, "V0M");
2699 // Count events with trigger selection, note: Set centrality percentile fix to -1 for pp for PID framework
2700 if (fUseJetPIDtask) {
2701 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
2702 fJetPIDtask[i]->IncrementEventCounter(fIsPP ? -1. : centPercent, AliAnalysisTaskPID::kTriggerSel);
2706 if (fUseInclusivePIDtask) {
2707 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
2708 fInclusivePIDtask[i]->IncrementEventCounter(fIsPP ? -1. : centPercent, AliAnalysisTaskPID::kTriggerSel);
2712 // *** vertex cut ***
2713 AliAODVertex* primVtx = fAOD->GetPrimaryVertex();
2715 Printf("%s:%d Primary vertex not found", (char*)__FILE__,__LINE__);
2719 Int_t nTracksPrim = primVtx->GetNContributors();
2720 fh1VertexNContributors->Fill(nTracksPrim);
2723 if (fDebug > 1) Printf("%s:%d primary vertex selection: %d", (char*)__FILE__,__LINE__,nTracksPrim);
2724 if(nTracksPrim <= 0) {
2725 if (fDebug > 1) Printf("%s:%d primary vertex selection: event REJECTED...",(char*)__FILE__,__LINE__);
2726 fh1EvtSelection->Fill(3.);
2727 PostData(1, fCommonHistList);
2731 TString primVtxName(primVtx->GetName());
2733 if(primVtxName.CompareTo("TPCVertex",TString::kIgnoreCase) == 1){
2734 if (fDebug > 1) Printf("%s:%d primary vertex selection: TPC vertex, event REJECTED...",(char*)__FILE__,__LINE__);
2735 fh1EvtSelection->Fill(5.);
2736 PostData(1, fCommonHistList);
2740 // Count events with trigger selection and vtx cut, note: Set centrality percentile fix to -1 for pp for PID framework
2741 if (fUseJetPIDtask) {
2742 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
2743 fJetPIDtask[i]->IncrementEventCounter(fIsPP ? -1. : centPercent, AliAnalysisTaskPID::kTriggerSelAndVtxCut);
2747 if (fUseInclusivePIDtask) {
2748 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
2749 fInclusivePIDtask[i]->IncrementEventCounter(fIsPP ? -1. : centPercent, AliAnalysisTaskPID::kTriggerSelAndVtxCut);
2753 fh1VertexZ->Fill(primVtx->GetZ());
2755 if(TMath::Abs(primVtx->GetZ())>fMaxVertexZ){
2756 if (fDebug > 1) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ());
2757 fh1EvtSelection->Fill(4.);
2758 PostData(1, fCommonHistList);
2762 // Count events with trigger selection, vtx cut and z vtx cut, note: Set centrality percentile fix to -1 for pp for PID framework
2763 if (fUseJetPIDtask) {
2764 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
2765 fJetPIDtask[i]->IncrementEventCounter(fIsPP ? -1. : centPercent, AliAnalysisTaskPID::kTriggerSelAndVtxCutAndZvtxCutNoPileUpRejection);
2769 if (fUseInclusivePIDtask) {
2770 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
2771 fInclusivePIDtask[i]->IncrementEventCounter(fIsPP ? -1. : centPercent, AliAnalysisTaskPID::kTriggerSelAndVtxCutAndZvtxCutNoPileUpRejection);
2776 // Store for each task, whether this task would tag this event as pile-up or not
2777 const Int_t arrSizeJet = TMath::Max(1, fNumJetPIDtasks);
2778 const Int_t arrSizeInclusive = TMath::Max(1, fNumInclusivePIDtasks);
2779 Bool_t isPileUpJetPIDtask[arrSizeJet];
2780 Bool_t isPileUpInclusivePIDtask[arrSizeInclusive];
2782 for (Int_t i = 0; i < arrSizeJet; i++)
2783 isPileUpJetPIDtask[i] = kFALSE;
2785 for (Int_t i = 0; i < arrSizeInclusive; i++)
2786 isPileUpInclusivePIDtask[i] = kFALSE;
2788 // Check whether there is at least one task that does not reject the event (saves processing time in the following code)
2789 Bool_t isPileUpForAllJetPIDTasks = kTRUE;
2790 Bool_t isPileUpForAllInclusivePIDTasks = kTRUE;
2792 // Count events with trigger selection, vtx cut, z vtx cut and after pile-up rejection (if enabled in that task)
2793 // Note: Set centrality percentile fix to -1 for pp for PID framework
2794 if (fUseJetPIDtask) {
2795 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
2796 isPileUpJetPIDtask[i] = fJetPIDtask[i]->GetIsPileUp(fAOD, fJetPIDtask[i]->GetPileUpRejectionType());
2797 if (!isPileUpJetPIDtask[i])
2798 fJetPIDtask[i]->IncrementEventCounter(fIsPP ? -1. : centPercent, AliAnalysisTaskPID::kTriggerSelAndVtxCutAndZvtxCut);
2800 isPileUpForAllJetPIDTasks = isPileUpForAllJetPIDTasks && isPileUpJetPIDtask[i];
2804 if (fUseInclusivePIDtask) {
2805 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
2806 isPileUpInclusivePIDtask[i] = fInclusivePIDtask[i]->GetIsPileUp(fAOD, fInclusivePIDtask[i]->GetPileUpRejectionType());
2807 if (!isPileUpInclusivePIDtask[i])
2808 fInclusivePIDtask[i]->IncrementEventCounter(fIsPP ? -1. : centPercent, AliAnalysisTaskPID::kTriggerSelAndVtxCutAndZvtxCut);
2810 isPileUpForAllInclusivePIDTasks = isPileUpForAllInclusivePIDTasks && isPileUpInclusivePIDtask[i];
2815 if (fDebug > 1) Printf("%s:%d event ACCEPTED ...",(char*)__FILE__,__LINE__);
2816 fh1EvtSelection->Fill(0.);
2817 fh1VtxSelection->Fill(primVtx->GetType());
2818 fh1EvtCent->Fill(centPercent);
2820 // Set centrality percentile fix to -1 for pp to be used for the PID framework
2826 // Call ConfigureTaskForCurrentEvent of PID tasks to ensure that everything is set up properly for the current event
2827 // (e.g. run/period dependence of max eta variation map)
2828 if (fUseInclusivePIDtask) {
2829 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++)
2830 if (!isPileUpInclusivePIDtask[i])
2831 fInclusivePIDtask[i]->ConfigureTaskForCurrentEvent(fAOD);
2834 if (fUseJetPIDtask) {
2835 for (Int_t i = 0; i < fNumJetPIDtasks; i++)
2836 if (!isPileUpJetPIDtask[i])
2837 fJetPIDtask[i]->ConfigureTaskForCurrentEvent(fAOD);
2842 //___ fill MC information __________________________________________________________________
2844 fh1Trials->Fill("#sum{ntrials}",fAvgTrials);
2846 if (fUseInclusivePIDtask) {
2847 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
2848 if (!isPileUpInclusivePIDtask[i])
2849 fInclusivePIDtask[i]->FillPythiaTrials(fAvgTrials);
2853 if (fUseJetPIDtask) {
2854 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
2855 if (!isPileUpJetPIDtask[i])
2856 fJetPIDtask[i]->FillPythiaTrials(fAvgTrials);
2860 if (pythiaGenHeaderFound) {
2861 fh1PtHard->Fill(ptHard);
2862 fh1PtHardTrials->Fill(ptHard,nTrials);
2865 //___ fetch jets __________________________________________________________________________
2867 Int_t nJ = GetListOfJets(fJetsRec, kJetsRec);
2869 if(nJ>=0) nRecJets = fJetsRec->GetEntries();
2870 if(fDebug>2)Printf("%s:%d Selected Rec jets: %d %d",(char*)__FILE__,__LINE__,nJ,nRecJets);
2871 if(nJ != nRecJets) Printf("%s:%d Mismatch Selected Rec Jets: %d %d",(char*)__FILE__,__LINE__,nJ,nRecJets);
2873 Int_t nJCuts = GetListOfJets(fJetsRecCuts, kJetsRecAcceptance);
2874 Int_t nRecJetsCuts = 0;
2875 if(nJCuts>=0) nRecJetsCuts = fJetsRecCuts->GetEntries();
2876 if(fDebug>2)Printf("%s:%d Selected Rec jets after cuts: %d %d",(char*)__FILE__,__LINE__,nJCuts,nRecJetsCuts);
2877 if(nRecJetsCuts != nJCuts) Printf("%s:%d Mismatch selected Rec jets after cuts: %d %d",(char*)__FILE__,__LINE__,nJCuts,nRecJetsCuts);
2878 fh1nRecJetsCuts->Fill(nRecJetsCuts);
2880 if(fJetTypeGen==kJetsKine || fJetTypeGen == kJetsKineAcceptance) fJetsGen->SetOwner(kTRUE); // kine aod jets allocated on heap, delete them with TList::Clear()
2882 Int_t nJGen = GetListOfJets(fJetsGen, fJetTypeGen);
2884 if(nJGen>=0) nGenJets = fJetsGen->GetEntries();
2885 if(fDebug>2)Printf("%s:%d Selected Gen jets: %d %d",(char*)__FILE__,__LINE__,nJGen,nGenJets);
2887 if(nJGen != nGenJets) Printf("%s:%d Mismatch selected Gen jets: %d %d",(char*)__FILE__,__LINE__,nJGen,nGenJets);
2888 fh1nGenJets->Fill(nGenJets);
2891 if(fJetTypeRecEff==kJetsKine || fJetTypeRecEff == kJetsKineAcceptance) fJetsRecEff->SetOwner(kTRUE); // kine aod jets allocated on heap, delete them with TList::Clear()
2892 Int_t nJRecEff = GetListOfJets(fJetsRecEff, fJetTypeRecEff);
2893 Int_t nRecEffJets = 0;
2894 if(nJRecEff>=0) nRecEffJets = fJetsRecEff->GetEntries();
2895 if(fDebug>2)Printf("%s:%d Selected RecEff jets: %d %d",(char*)__FILE__,__LINE__,nJRecEff,nRecEffJets);
2896 if(nJRecEff != nRecEffJets) Printf("%s:%d Mismatch selected RecEff jets: %d %d",(char*)__FILE__,__LINE__,nJRecEff,nRecEffJets);
2897 fh1nRecEffJets->Fill(nRecEffJets);
2900 Int_t nEmbeddedJets = 0;
2901 TArrayI iEmbeddedMatchIndex;
2902 TArrayF fEmbeddedPtFraction;
2905 if(fBranchEmbeddedJets.Length()){
2906 Int_t nJEmbedded = GetListOfJets(fJetsEmbedded, kJetsEmbedded);
2907 if(nJEmbedded>=0) nEmbeddedJets = fJetsEmbedded->GetEntries();
2908 if(fDebug>2)Printf("%s:%d Selected Embedded jets: %d %d",(char*)__FILE__,__LINE__,nJEmbedded,nEmbeddedJets);
2909 if(nJEmbedded != nEmbeddedJets) Printf("%s:%d Mismatch Selected Embedded Jets: %d %d",(char*)__FILE__,__LINE__,nJEmbedded,nEmbeddedJets);
2910 fh1nEmbeddedJets->Fill(nEmbeddedJets);
2912 Float_t maxDist = 0.3;
2914 iEmbeddedMatchIndex.Set(nEmbeddedJets);
2915 fEmbeddedPtFraction.Set(nEmbeddedJets);
2917 iEmbeddedMatchIndex.Reset(-1);
2918 fEmbeddedPtFraction.Reset(0);
2920 AliAnalysisHelperJetTasks::GetJetMatching(fJetsEmbedded, nEmbeddedJets,
2921 fJetsRecCuts, nRecJetsCuts,
2922 iEmbeddedMatchIndex, fEmbeddedPtFraction,
2927 //____ fetch background clusters ___________________________________________________
2929 (fBckgType[0]==kBckgClusters || fBckgType[1]==kBckgClusters || fBckgType[2]==kBckgClusters || fBckgType[3]==kBckgClusters || fBckgType[4]==kBckgClusters ||
2930 fBckgType[0]==kBckgClustersOutLeading || fBckgType[1]==kBckgClustersOutLeading || fBckgType[2]==kBckgClustersOutLeading ||
2931 fBckgType[3]==kBckgClustersOutLeading || fBckgType[4]==kBckgClustersOutLeading)){
2933 Int_t nBJ = GetListOfBckgJets(fBckgJetsRec, kJetsRec);
2934 Int_t nRecBckgJets = 0;
2935 if(nBJ>=0) nRecBckgJets = fBckgJetsRec->GetEntries();
2936 if(fDebug>2)Printf("%s:%d Selected Rec background jets: %d %d",(char*)__FILE__,__LINE__,nBJ,nRecBckgJets);
2937 if(nBJ != nRecBckgJets) Printf("%s:%d Mismatch Selected Rec background jets: %d %d",(char*)__FILE__,__LINE__,nBJ,nRecBckgJets);
2939 Int_t nBJCuts = GetListOfBckgJets(fBckgJetsRecCuts, kJetsRecAcceptance);
2940 Int_t nRecBckgJetsCuts = 0;
2941 if(nBJCuts>=0) nRecBckgJetsCuts = fBckgJetsRecCuts->GetEntries();
2942 if(fDebug>2)Printf("%s:%d Selected Rec background jets after cuts: %d %d",(char*)__FILE__,__LINE__,nJCuts,nRecJetsCuts);
2943 if(nRecBckgJetsCuts != nBJCuts) Printf("%s:%d Mismatch selected Rec background jets after cuts: %d %d",(char*)__FILE__,__LINE__,nBJCuts,nRecBckgJetsCuts);
2944 fh1nRecBckgJetsCuts->Fill(nRecBckgJetsCuts);
2946 if(0){ // protection OB - not yet implemented
2947 if(fJetTypeGen==kJetsKine || fJetTypeGen == kJetsKineAcceptance) fBckgJetsGen->SetOwner(kTRUE); // kine aod jets allocated on heap, delete them with TList::Clear()
2948 Int_t nBJGen = GetListOfBckgJets(fBckgJetsGen, fJetTypeGen);
2949 Int_t nGenBckgJets = 0;
2950 if(nBJGen>=0) nGenBckgJets = fBckgJetsGen->GetEntries();
2951 if(fDebug>2)Printf("%s:%d Selected Gen background jets: %d %d",(char*)__FILE__,__LINE__,nBJGen,nGenBckgJets);
2952 if(nBJGen != nGenBckgJets) Printf("%s:%d Mismatch selected Gen background jets: %d %d",(char*)__FILE__,__LINE__,nBJGen,nGenBckgJets);
2953 fh1nGenBckgJets->Fill(nGenBckgJets);
2958 //____ fetch particles __________________________________________________________
2961 if(fUseExtraTracks == 1) nTCuts = GetListOfTracks(fTracksRecCuts, kTrackAODExtraCuts);
2962 else if(fUseExtraTracks == -1) nTCuts = GetListOfTracks(fTracksRecCuts, kTrackAODExtraonlyCuts);
2963 else nTCuts = GetListOfTracks(fTracksRecCuts, kTrackAODCuts);
2965 Int_t nRecPartCuts = 0;
2966 if(nTCuts>=0) nRecPartCuts = fTracksRecCuts->GetEntries();
2967 if(fDebug>2)Printf("%s:%d Selected Rec tracks after cuts: %d %d",(char*)__FILE__,__LINE__,nTCuts,nRecPartCuts);
2968 if(nRecPartCuts != nTCuts) Printf("%s:%d Mismatch selected Rec tracks after cuts: %d %d",
2969 (char*)__FILE__,__LINE__,nTCuts,nRecPartCuts);
2970 fh1EvtMult->Fill(nRecPartCuts);
2973 Int_t nTGen = GetListOfTracks(fTracksGen,fTrackTypeGen);
2975 if(nTGen>=0) nGenPart = fTracksGen->GetEntries();
2976 if(fDebug>2)Printf("%s:%d Selected Gen tracks: %d %d",(char*)__FILE__,__LINE__,nTGen,nGenPart);
2977 if(nGenPart != nTGen) Printf("%s:%d Mismatch selected Gen tracks: %d %d",(char*)__FILE__,__LINE__,nTGen,nGenPart);
2979 // Just cut on filterBit, but take the full acceptance
2980 Int_t nTCutsEfficiency = GetListOfTracks(fTracksRecCutsEfficiency, kTrackAODQualityCuts);
2981 Int_t nRecPartCutsEfficiency = 0;
2982 if(nTCutsEfficiency>=0) nRecPartCutsEfficiency = fTracksRecCutsEfficiency->GetEntries();
2983 if(fDebug>2)Printf("%s:%d Selected Rec tracks efficiency after cuts: %d %d",(char*)__FILE__,__LINE__,
2984 nTCutsEfficiency,nRecPartCutsEfficiency);
2985 if(nRecPartCutsEfficiency != nTCutsEfficiency) Printf("%s:%d Mismatch selected Rec tracks after cuts: %d %d",
2986 (char*)__FILE__,__LINE__,nTCutsEfficiency,nRecPartCutsEfficiency);
2988 AliPIDResponse* pidResponse = 0x0;
2989 Bool_t tuneOnDataTPC = kFALSE;
2990 if (fUseJetPIDtask || fUseInclusivePIDtask) {
2991 if (!inputHandler) {
2992 AliFatal("Input handler needed");
2996 // PID response object
2997 pidResponse = inputHandler->GetPIDResponse();
2999 AliFatal("PIDResponse object was not created");
3003 tuneOnDataTPC = pidResponse->IsTunedOnData() &&
3004 ((pidResponse->GetTunedOnDataMask() & AliPIDResponse::kDetTPC) == AliPIDResponse::kDetTPC);
3009 if(fDebug>2)Printf("%s:%d Starting processing...",(char*)__FILE__,__LINE__);
3011 //____ analysis, fill histos ___________________________________________________
3016 TClonesArray *tca = fUseInclusivePIDtask ? dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName())) : 0x0;
3018 // Fill efficiency for generated primaries and also fill histos for generated yields (primaries + all)
3019 // Efficiency, inclusive - particle level
3020 if (fUseInclusivePIDtask && tca && !isPileUpForAllInclusivePIDTasks) {
3021 for (Int_t it = 0; it < tca->GetEntriesFast(); it++) {
3022 AliAODMCParticle* part = dynamic_cast<AliAODMCParticle*>(tca->At(it));
3026 // Define clean MC sample with corresponding particle level track cuts:
3027 // - MC-track must be in desired eta range
3028 // - MC-track must be physical primary
3029 // - Species must be one of those in question
3031 if (part->Eta() > fTrackEtaMax || part->Eta() < fTrackEtaMin)
3034 Int_t mcID = AliAnalysisTaskPID::PDGtoMCID(part->GetPdgCode());
3036 // Following lines are not needed - just keep other species (like casecades) - will end up in overflow bin
3037 // and only affect the efficiencies for all (i.e. not identified) what is desired!
3038 //if (mcID == AliPID::kUnknown)
3041 if (!part->IsPhysicalPrimary())
3044 Int_t iMother = part->GetMother();
3046 continue; // Not a physical primary
3049 Double_t pT = part->Pt();
3051 // AliAODMCParticle->Charge() calls TParticlePDG->Charge(), which returns the charge in units of e0 / 3
3052 Double_t chargeMC = part->Charge() / 3.;
3054 if (TMath::Abs(chargeMC) < 0.01)
3055 continue; // Reject neutral particles (only relevant, if mcID is not used)
3057 Double_t valuesGenYield[AliAnalysisTaskPID::kGenYieldNumAxes] = { static_cast<Double_t>(mcID), pT, centPercent, -1, -1, -1, -1 };
3059 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
3060 if (!isPileUpInclusivePIDtask[i] && fInclusivePIDtask[i]->IsInAcceptedEtaRange(TMath::Abs(part->Eta()))) {
3061 valuesGenYield[fInclusivePIDtask[i]->GetIndexOfChargeAxisGenYield()] = chargeMC;
3062 fInclusivePIDtask[i]->FillGeneratedYield(valuesGenYield);
3066 Double_t valuesEff[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), pT, part->Eta(), chargeMC,
3067 centPercent, -1, -1, -1 };// no jet pT etc since inclusive spectrum
3068 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
3069 if (!isPileUpInclusivePIDtask[i])
3070 fInclusivePIDtask[i]->FillEfficiencyContainer(valuesEff, AliAnalysisTaskPID::kStepGenWithGenCuts);
3075 if (fUseInclusivePIDtask && !isPileUpForAllInclusivePIDTasks) {
3076 //Efficiency, inclusive - detector level
3077 for(Int_t it=0; it<nRecPartCutsEfficiency; ++it){
3078 // fill inclusive tracks XXX, they have the same track cuts!
3079 AliAODTrack * inclusiveaod = dynamic_cast<AliAODTrack*>(fTracksRecCutsEfficiency->At(it));
3081 Double_t dEdxTPC = tuneOnDataTPC ? pidResponse->GetTPCsignalTunedOnData(inclusiveaod)
3082 : inclusiveaod->GetTPCsignal();
3087 Bool_t survivedTPCCutMIGeo = AliAnalysisTaskPID::TPCCutMIGeo(inclusiveaod, InputEvent());
3088 Bool_t survivedTPCnclCut = AliAnalysisTaskPID::TPCnclCut(inclusiveaod);
3090 Int_t label = TMath::Abs(inclusiveaod->GetLabel());
3092 // find MC track in our list, if available
3093 AliAODMCParticle* gentrack = tca ? dynamic_cast<AliAODMCParticle*>(tca->At(label)) : 0x0;
3097 pdg = gentrack->GetPdgCode();
3099 // For efficiency: Reconstructed track has survived all cuts on the detector level (no cut on eta acceptance)
3100 // and has an associated MC track
3101 // -> Check whether associated MC track belongs to the clean MC sample defined above,
3102 // i.e. survives the particle level track cuts
3104 Int_t mcID = AliAnalysisTaskPID::PDGtoMCID(pdg);
3106 // Following lines are not needed - just keep other species (like casecades) - will end up in overflow bin
3107 // and only affect the efficiencies for all (i.e. not identified) what is desired!
3108 //if (mcID == AliPID::kUnknown)
3111 // Fill efficiency for reconstructed primaries
3112 if (!gentrack->IsPhysicalPrimary())
3115 * Int_t iMother = gentrack->GetMother();
3117 * continue; // Not a physical primary
3120 if (gentrack->Eta() > fTrackEtaMax || gentrack->Eta() < fTrackEtaMin)
3123 // AliAODMCParticle->Charge() calls TParticlePDG->Charge(), which returns the charge in units of e0 / 3
3124 Double_t value[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), gentrack->Pt(), gentrack->Eta(), gentrack->Charge() / 3.,
3126 -1, -1, -1 };// no jet pT etc since inclusive spectrum
3127 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
3128 if (!isPileUpInclusivePIDtask[i] &&
3129 ((!fInclusivePIDtask[i]->GetUseTPCCutMIGeo() && !fInclusivePIDtask[i]->GetUseTPCnclCut()) ||
3130 (survivedTPCCutMIGeo && fInclusivePIDtask[i]->GetUseTPCCutMIGeo()) ||
3131 (survivedTPCnclCut && fInclusivePIDtask[i]->GetUseTPCnclCut())))
3132 fInclusivePIDtask[i]->FillEfficiencyContainer(value, AliAnalysisTaskPID::kStepRecWithGenCuts);
3135 Double_t valueMeas[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), inclusiveaod->Pt(), inclusiveaod->Eta(),
3136 static_cast<Double_t>(inclusiveaod->Charge()), centPercent,
3137 -1, -1, -1 };// no jet pT etc since inclusive spectrum
3138 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
3139 if (!isPileUpInclusivePIDtask[i] &&
3140 ((!fInclusivePIDtask[i]->GetUseTPCCutMIGeo() && !fInclusivePIDtask[i]->GetUseTPCnclCut()) ||
3141 (survivedTPCCutMIGeo && fInclusivePIDtask[i]->GetUseTPCCutMIGeo()) ||
3142 (survivedTPCnclCut && fInclusivePIDtask[i]->GetUseTPCnclCut())))
3143 fInclusivePIDtask[i]->FillEfficiencyContainer(valueMeas, AliAnalysisTaskPID::kStepRecWithGenCutsMeasuredObs);
3151 for(Int_t it=0; it<nRecPartCuts; ++it){
3152 AliVParticle *part = dynamic_cast<AliVParticle*>(fTracksRecCuts->At(it));
3153 if(part)fQATrackHistosRecCuts->FillTrackQA( part->Eta(), TVector2::Phi_0_2pi(part->Phi()), part->Pt() );
3155 // fill inclusive tracks XXX, they have the same track cuts!
3156 AliAODTrack * inclusiveaod = dynamic_cast<AliAODTrack*>(fTracksRecCuts->At(it));
3158 if (fUseInclusivePIDtask && !isPileUpForAllInclusivePIDTasks) {
3159 Double_t dEdxTPC = tuneOnDataTPC ? pidResponse->GetTPCsignalTunedOnData(inclusiveaod)
3160 : inclusiveaod->GetTPCsignal();
3165 Bool_t survivedTPCCutMIGeo = AliAnalysisTaskPID::TPCCutMIGeo(inclusiveaod, InputEvent());
3166 Bool_t survivedTPCnclCut = AliAnalysisTaskPID::TPCnclCut(inclusiveaod);
3168 Int_t label = TMath::Abs(inclusiveaod->GetLabel());
3170 // find MC track in our list, if available
3171 AliAODMCParticle* gentrack = tca ? dynamic_cast<AliAODMCParticle*>(tca->At(label)) : 0x0;
3175 pdg = gentrack->GetPdgCode();
3177 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
3178 if (!isPileUpInclusivePIDtask[i] &&
3179 ((!fInclusivePIDtask[i]->GetUseTPCCutMIGeo() && !fInclusivePIDtask[i]->GetUseTPCnclCut()) ||
3180 (survivedTPCCutMIGeo && fInclusivePIDtask[i]->GetUseTPCCutMIGeo()) ||
3181 (survivedTPCnclCut && fInclusivePIDtask[i]->GetUseTPCnclCut()))) {
3182 if (fInclusivePIDtask[i]->IsInAcceptedEtaRange(TMath::Abs(inclusiveaod->Eta())))
3183 fInclusivePIDtask[i]->ProcessTrack(inclusiveaod, pdg, centPercent, -1); // no jet pT since inclusive spectrum
3188 Int_t mcID = AliAnalysisTaskPID::PDGtoMCID(pdg);
3189 Double_t valueRecAllCuts[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), inclusiveaod->Pt(), inclusiveaod->Eta(),
3190 static_cast<Double_t>(inclusiveaod->Charge()), centPercent,
3192 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
3193 if (!isPileUpInclusivePIDtask[i] &&
3194 ((!fInclusivePIDtask[i]->GetUseTPCCutMIGeo() && !fInclusivePIDtask[i]->GetUseTPCnclCut()) ||
3195 (survivedTPCCutMIGeo && fInclusivePIDtask[i]->GetUseTPCCutMIGeo()) ||
3196 (survivedTPCnclCut && fInclusivePIDtask[i]->GetUseTPCnclCut())))
3197 fInclusivePIDtask[i]->FillEfficiencyContainer(valueRecAllCuts, AliAnalysisTaskPID::kStepRecWithRecCutsMeasuredObs);
3200 Double_t weight = IsSecondaryWithStrangeMotherMC(gentrack) ? GetMCStrangenessFactorCMS(gentrack) : 1.0;
3201 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
3202 if (!isPileUpInclusivePIDtask[i] &&
3203 ((!fInclusivePIDtask[i]->GetUseTPCCutMIGeo() && !fInclusivePIDtask[i]->GetUseTPCnclCut()) ||
3204 (survivedTPCCutMIGeo && fInclusivePIDtask[i]->GetUseTPCCutMIGeo()) ||
3205 (survivedTPCnclCut && fInclusivePIDtask[i]->GetUseTPCnclCut())))
3206 fInclusivePIDtask[i]->FillEfficiencyContainer(valueRecAllCuts,
3207 AliAnalysisTaskPID::kStepRecWithRecCutsMeasuredObsStrangenessScaled,
3211 if (gentrack->IsPhysicalPrimary()) {
3212 // AliAODMCParticle->Charge() calls TParticlePDG->Charge(), which returns the charge in units of e0 / 3
3213 Double_t valueGenAllCuts[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), gentrack->Pt(), gentrack->Eta(),
3214 gentrack->Charge() / 3., centPercent, -1, -1,
3217 Double_t valuePtResolution[AliAnalysisTaskPID::kPtResNumAxes] = { -1, gentrack->Pt(), inclusiveaod->Pt(),
3218 gentrack->Charge() / 3., centPercent };
3220 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
3221 if (!isPileUpInclusivePIDtask[i] &&
3222 ((!fInclusivePIDtask[i]->GetUseTPCCutMIGeo() && !fInclusivePIDtask[i]->GetUseTPCnclCut()) ||
3223 (survivedTPCCutMIGeo && fInclusivePIDtask[i]->GetUseTPCCutMIGeo()) ||
3224 (survivedTPCnclCut && fInclusivePIDtask[i]->GetUseTPCnclCut()))) {
3225 fInclusivePIDtask[i]->FillEfficiencyContainer(valueRecAllCuts,
3226 AliAnalysisTaskPID::kStepRecWithRecCutsMeasuredObsPrimaries);
3227 fInclusivePIDtask[i]->FillEfficiencyContainer(valueGenAllCuts,
3228 AliAnalysisTaskPID::kStepRecWithRecCutsPrimaries);
3230 fInclusivePIDtask[i]->FillPtResolution(mcID, valuePtResolution);
3238 for(Int_t it=0; it<nGenPart; ++it){
3239 AliVParticle *part = dynamic_cast<AliVParticle*>(fTracksGen->At(it));
3240 if(part)fQATrackHistosGen->FillTrackQA( part->Eta(), TVector2::Phi_0_2pi(part->Phi()), part->Pt());
3246 for(Int_t ij=0; ij<nRecJets; ++ij){
3247 AliAODJet* jet = dynamic_cast<AliAODJet*>(fJetsRec->At(ij));
3248 if(jet)fQAJetHistosRec->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt());
3253 if(fQAMode || fFFMode){
3255 for(Int_t ij=0; ij<nGenJets; ++ij){
3256 AliAODJet* jet = dynamic_cast<AliAODJet*>(fJetsGen->At(ij));
3260 if(fQAMode&2) fQAJetHistosGen->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt());
3262 if(fQAMode&2 && (ij==0)) fQAJetHistosGenLeading->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt() );
3264 if((ij==0) || !fOnlyLeadingJets){ // leading jets or all jets
3265 TList* jettracklist = new TList();
3266 Double_t sumPt = 0.;
3267 Bool_t isBadJet = kFALSE;
3269 if(GetFFRadius()<=0){
3270 GetJetTracksTrackrefs(jettracklist, jet, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet);
3272 GetJetTracksPointing(fTracksGen, jettracklist, jet, GetFFRadius(), sumPt, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet);
3275 if(GetFFMinNTracks()>0 && jettracklist->GetSize()<=GetFFMinNTracks()) isBadJet = kTRUE;
3276 if(isBadJet) continue;
3278 for(Int_t it=0; it<jettracklist->GetSize(); ++it){
3280 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(jettracklist->At(it));
3281 if(!trackVP)continue;
3282 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
3284 Float_t jetPt = jet->Pt();
3285 Float_t trackPt = trackV->Pt();
3287 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
3289 if(fFFMode && (ij==0)) fFFHistosGen->FillFF( trackPt, jetPt, incrementJetPt );
3290 if(fFFMode) fFFHistosGenInc->FillFF( trackPt, jetPt, incrementJetPt );
3291 if(it==0){ // leading track
3292 if(fFFMode) fFFHistosGenLeadingTrack->FillFF( trackPt, jetPt, kTRUE );
3295 if (fUseJetPIDtask && incrementJetPt) {
3296 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
3297 if (!isPileUpJetPIDtask[i])
3298 fJetPIDtask[i]->FillGenJets(fJetPIDtask[i]->GetCentralityPercentile(evtForCentDetermination), jetPt);
3303 Int_t mcID = AliAnalysisTaskPID::PDGtoMCID(trackVP->PdgCode());
3304 if (mcID != AliPID::kUnknown) {
3305 // WARNING: The number of jets for the different species does not make sense -> One has to take
3306 // the number of jets for ALL particles. Thus, just do not fill the num of jet histos in order
3307 // not to get confused
3308 fIDFFHistosGen[mcID]->FillFF(trackPt, jetPt, kFALSE);
3311 Int_t pidWeightedSpecies = fJetPIDtask->GetRandomParticleTypeAccordingToParticleFractions(trackPt, jetPt, centPercent, kTRUE);
3312 if (pidWeightedSpecies < 0 || pidWeightedSpecies >= AliPID::kSPECIES) {
3313 Printf("Failed to determine particle ID for track in jet (gen) -> ID FF histos not filled with this track!");
3314 Printf("Track details: trackPt %f, jetPt %f, centrality %f!", trackPt, jetPt, centPercent);
3317 // WARNING: The number of jets for the different species does not make sense -> One has to take
3318 // the number of jets for ALL particles. Thus, just do not fill the num of jet histos in order
3319 // not to get confused
3320 fIDFFHistosGen[pidWeightedSpecies]->FillFF(trackPt, jetPt, kFALSE);
3327 // Efficiency, jets - particle level
3328 if (fUseJetPIDtask && !isPileUpForAllJetPIDTasks) {
3329 AliAODMCParticle* part = dynamic_cast<AliAODMCParticle*>(jettracklist->At(it));
3331 AliError("expected ref track not found ");
3334 // Fill efficiency for generated primaries and also fill histos for generated yields (primaries + all)
3335 if (part->Eta() > fTrackEtaMax || part->Eta() < fTrackEtaMin)
3338 Int_t mcID = AliAnalysisTaskPID::PDGtoMCID(part->GetPdgCode());
3340 // Following lines are not needed - just keep other species (like casecades) - will end up in overflow bin
3341 // and only affect the efficiencies for all (i.e. not identified) what is desired!
3342 //if (mcID == AliPID::kUnknown)
3345 if (!part->IsPhysicalPrimary())
3348 // Int_t iMother = part->GetMother();
3349 // if (iMother >= 0)
3350 // continue; // Not a physical primary
3353 Double_t z = -1., xi = -1.;
3354 AliAnalysisTaskPID::GetJetTrackObservables(trackPt, jetPt, z, xi);
3356 // AliAODMCParticle->Charge() calls TParticlePDG->Charge(), which returns the charge in units of e0 / 3
3357 Double_t chargeMC = part->Charge() / 3.;
3359 if (TMath::Abs(chargeMC) < 0.01)
3360 continue; // Reject neutral particles (only relevant, if mcID is not used)
3362 Double_t valuesGenYield[AliAnalysisTaskPID::kGenYieldNumAxes] = { static_cast<Double_t>(mcID), trackPt, centPercent, jetPt, z, xi, chargeMC };
3364 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
3365 if (!isPileUpJetPIDtask[i] && fJetPIDtask[i]->IsInAcceptedEtaRange(TMath::Abs(part->Eta()))) {
3366 valuesGenYield[fJetPIDtask[i]->GetIndexOfChargeAxisGenYield()] = chargeMC;
3367 fJetPIDtask[i]->FillGeneratedYield(valuesGenYield);
3372 Double_t valuesEff[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), trackPt, part->Eta(), chargeMC,
3373 centPercent, jetPt, z, xi };
3374 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
3375 if (!isPileUpJetPIDtask[i])
3376 fJetPIDtask[i]->FillEfficiencyContainer(valuesEff, AliAnalysisTaskPID::kStepGenWithGenCuts);
3382 if(fBckgType[0]!=kBckgNone)
3383 FillBckgHistos(fBckgType[0], fTracksGen, fJetsGen, jet,
3384 fFFBckgHisto0Gen, fQABckgHisto0Gen);
3385 if(fBckgType[1]!=kBckgNone)
3386 FillBckgHistos(fBckgType[1], fTracksGen, fJetsGen, jet,
3387 fFFBckgHisto1Gen, fQABckgHisto1Gen);
3388 if(fBckgType[2]!=kBckgNone)
3389 FillBckgHistos(fBckgType[2], fTracksGen, fJetsGen, jet,
3390 fFFBckgHisto2Gen, fQABckgHisto2Gen);
3391 if(fBckgType[3]!=kBckgNone)
3392 FillBckgHistos(fBckgType[3], fTracksGen, fJetsGen, jet,
3393 fFFBckgHisto3Gen, fQABckgHisto3Gen);
3394 if(fBckgType[4]!=kBckgNone)
3395 FillBckgHistos(fBckgType[4], fTracksGen, fJetsGen, jet,
3396 fFFBckgHisto4Gen, fQABckgHisto4Gen);
3397 } // end if(fBckgMode)
3400 if(fJSMode) FillJetShape(jet, jettracklist, fProNtracksLeadingJetGen, fProDelRPtSumGen, fProDelR80pcPtGen);
3402 delete jettracklist;
3406 for(Int_t ij=0; ij<nRecJetsCuts; ++ij){
3408 AliAODJet* jet = (AliAODJet*)(fJetsRecCuts->At(ij));
3409 if(fQAMode&2) fQAJetHistosRecCuts->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt());
3410 if(fQAMode&2 && (ij==0)) fQAJetHistosRecCutsLeading->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt() );
3412 if((ij==0) || !fOnlyLeadingJets){ // leading jets or all jets
3414 fhJetPtRefMultEta5->Fill(refMult5, jet->Pt());
3415 fhJetPtRefMultEta8->Fill(refMult8, jet->Pt());
3416 fhJetPtMultPercent->Fill(centPercentPP, jet->Pt());
3418 Double_t ptFractionEmbedded = 0;
3419 AliAODJet* embeddedJet = 0;
3421 if(fBranchEmbeddedJets.Length()){ // find embedded jet
3423 Int_t indexEmbedded = -1;
3424 for(Int_t i=0; i<nEmbeddedJets; i++){
3425 if(iEmbeddedMatchIndex[i] == ij){
3427 ptFractionEmbedded = fEmbeddedPtFraction[i];
3431 fh1IndexEmbedded->Fill(indexEmbedded);
3432 fh1FractionPtEmbedded->Fill(ptFractionEmbedded);
3434 if(indexEmbedded>-1){
3436 embeddedJet = dynamic_cast<AliAODJet*>(fJetsEmbedded->At(indexEmbedded));
3437 if(!embeddedJet) continue;
3439 Double_t deltaPt = jet->Pt() - embeddedJet->Pt();
3440 Double_t deltaR = jet->DeltaR((AliVParticle*) (embeddedJet));
3442 fh2DeltaPtVsJetPtEmbedded->Fill(embeddedJet->Pt(),deltaPt);
3443 fh2DeltaPtVsRecJetPtEmbedded->Fill(jet->Pt(),deltaPt);
3444 fh1DeltaREmbedded->Fill(deltaR);
3448 // get tracks in jet
3449 TList* jettracklist = new TList();
3450 Double_t sumPt = 0.;
3451 Bool_t isBadJet = kFALSE;
3453 if(GetFFRadius()<=0){
3454 GetJetTracksTrackrefs(jettracklist, jet, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet);
3456 if(fUseEmbeddedJetAxis){
3457 if(embeddedJet) GetJetTracksPointing(fTracksRecCuts, jettracklist, embeddedJet, GetFFRadius(), sumPt, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet);
3459 else GetJetTracksPointing(fTracksRecCuts, jettracklist, jet, GetFFRadius(), sumPt, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJet);
3462 if(GetFFMinNTracks()>0 && jettracklist->GetSize()<=GetFFMinNTracks()) isBadJet = kTRUE;
3464 if(isBadJet) continue;
3466 if(ptFractionEmbedded>=fCutFractionPtEmbedded){ // if no embedding: ptFraction = cutFraction = 0
3468 TClonesArray *tca = fUseJetPIDtask ? dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName())) : 0x0;
3470 for(Int_t it=0; it<jettracklist->GetSize(); ++it){
3471 AliAODTrack * aodtrack = dynamic_cast<AliAODTrack*>(jettracklist->At(it));
3472 if(!aodtrack) continue;
3474 Double_t pT = aodtrack->Pt();
3475 Float_t jetPt = jet->Pt();
3476 if(fUseEmbeddedJetPt){
3477 if(embeddedJet) jetPt = embeddedJet->Pt();
3481 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
3483 if(fFFMode && (ij==0)) fFFHistosRecCuts->FillFF(pT, jetPt, incrementJetPt);
3484 if(fFFMode) fFFHistosRecCutsInc->FillFF(pT, jetPt, incrementJetPt);
3486 if(it==0){ // leading track
3487 if(fFFMode) fFFHistosRecLeadingTrack->FillFF(pT, jetPt, kTRUE);
3490 if (fUseJetPIDtask && incrementJetPt) {
3491 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
3492 if (!isPileUpJetPIDtask[i])
3493 fJetPIDtask[i]->FillRecJets(fJetPIDtask[i]->GetCentralityPercentile(evtForCentDetermination), jetPt);
3497 if (fUseJetPIDtask && (!isPileUpForAllJetPIDTasks || fIDFFMode)) {
3498 Double_t dEdxTPC = tuneOnDataTPC ? pidResponse->GetTPCsignalTunedOnData(aodtrack)
3499 : aodtrack->GetTPCsignal();
3504 Bool_t survivedTPCCutMIGeo = AliAnalysisTaskPID::TPCCutMIGeo(aodtrack, InputEvent());
3505 Bool_t survivedTPCnclCut = AliAnalysisTaskPID::TPCnclCut(aodtrack);
3507 Int_t label = TMath::Abs(aodtrack->GetLabel());
3509 // Find MC track in our list, if available
3510 AliAODMCParticle* gentrack = tca ? dynamic_cast<AliAODMCParticle*>(tca->At(label)) : 0x0;
3512 Int_t mcID = AliPID::kUnknown;
3515 pdg = gentrack->GetPdgCode();
3517 // Secondaries, jets
3518 mcID = AliAnalysisTaskPID::PDGtoMCID(pdg);
3520 Double_t z = -1., xi = -1.;
3521 AliAnalysisTaskPID::GetJetTrackObservables(pT, jetPt, z, xi);
3523 Double_t valueRecAllCuts[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), pT, aodtrack->Eta(), static_cast<Double_t>(aodtrack->Charge()),
3524 centPercent, jetPt, z, xi };
3525 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
3526 if (!isPileUpJetPIDtask[i] &&
3527 ((!fJetPIDtask[i]->GetUseTPCCutMIGeo() && !fJetPIDtask[i]->GetUseTPCnclCut()) ||
3528 (survivedTPCCutMIGeo && fJetPIDtask[i]->GetUseTPCCutMIGeo()) ||
3529 (survivedTPCnclCut && fJetPIDtask[i]->GetUseTPCnclCut())))
3530 fJetPIDtask[i]->FillEfficiencyContainer(valueRecAllCuts, AliAnalysisTaskPID::kStepRecWithRecCutsMeasuredObs);
3533 Double_t weight = IsSecondaryWithStrangeMotherMC(gentrack) ? GetMCStrangenessFactorCMS(gentrack) : 1.0;
3534 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
3535 if (!isPileUpJetPIDtask[i] &&
3536 ((!fJetPIDtask[i]->GetUseTPCCutMIGeo() && !fJetPIDtask[i]->GetUseTPCnclCut()) ||
3537 (survivedTPCCutMIGeo && fJetPIDtask[i]->GetUseTPCCutMIGeo()) ||
3538 (survivedTPCnclCut && fJetPIDtask[i]->GetUseTPCnclCut())))
3539 fJetPIDtask[i]->FillEfficiencyContainer(valueRecAllCuts,
3540 AliAnalysisTaskPID::kStepRecWithRecCutsMeasuredObsStrangenessScaled,
3544 if (gentrack->IsPhysicalPrimary()) {
3545 // AliAODMCParticle->Charge() calls TParticlePDG->Charge(), which returns the charge in units of e0 / 3
3546 Double_t genPt = gentrack->Pt();
3547 Double_t genZ = -1., genXi = -1.;
3548 AliAnalysisTaskPID::GetJetTrackObservables(genPt, jetPt, genZ, genXi);
3549 Double_t valueGenAllCuts[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), genPt, gentrack->Eta(),
3550 gentrack->Charge() / 3., centPercent, jetPt, genZ,
3553 Double_t valuePtResolution[AliAnalysisTaskPID::kPtResNumAxes] = { jetPt, genPt, pT, gentrack->Charge() / 3., centPercent };
3555 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
3556 if (!isPileUpJetPIDtask[i] &&
3557 ((!fJetPIDtask[i]->GetUseTPCCutMIGeo() && !fJetPIDtask[i]->GetUseTPCnclCut()) ||
3558 (survivedTPCCutMIGeo && fJetPIDtask[i]->GetUseTPCCutMIGeo()) ||
3559 (survivedTPCnclCut && fJetPIDtask[i]->GetUseTPCnclCut()))) {
3560 fJetPIDtask[i]->FillEfficiencyContainer(valueRecAllCuts,
3561 AliAnalysisTaskPID::kStepRecWithRecCutsMeasuredObsPrimaries);
3562 fJetPIDtask[i]->FillEfficiencyContainer(valueGenAllCuts,
3563 AliAnalysisTaskPID::kStepRecWithRecCutsPrimaries);
3565 fJetPIDtask[i]->FillPtResolution(mcID, valuePtResolution);
3571 Bool_t filledDCA = kFALSE;
3573 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
3574 if (!isPileUpJetPIDtask[i] &&
3575 ((!fJetPIDtask[i]->GetUseTPCCutMIGeo() && !fJetPIDtask[i]->GetUseTPCnclCut()) ||
3576 (survivedTPCCutMIGeo && fJetPIDtask[i]->GetUseTPCCutMIGeo()) ||
3577 (survivedTPCnclCut && fJetPIDtask[i]->GetUseTPCnclCut()))) {
3578 if (fJetPIDtask[i]->IsInAcceptedEtaRange(TMath::Abs(aodtrack->Eta()))) {
3579 fJetPIDtask[i]->ProcessTrack(aodtrack, pdg, centPercent, jetPt);
3581 // Fill DCA histo (once) if at least one of the PID tasks accecpts the track
3585 Double_t dca[2] = {0., 0.}; // 0: xy; 1: z
3586 if (aodtrack->IsGlobalConstrained()) {
3587 dca[0] = aodtrack->DCA();
3588 dca[1] = aodtrack->ZAtDCA();
3591 Double_t v[3] = {0, };
3592 Double_t pos[3] = {0, };
3594 aodtrack->GetXYZ(pos);
3595 dca[0] = pos[0] - v[0];
3596 dca[1] = pos[1] - v[1];
3599 // "Unidentified" for data and MC
3600 fhDCA_XY->Fill(pT, dca[0]);
3601 fhDCA_Z->Fill(pT, dca[1]);
3603 // "Identified" for MC
3604 if (gentrack && mcID != AliPID::kUnknown) {
3606 if (gentrack->IsPhysicalPrimary()) {
3607 fhDCA_XY_prim_MCID[mcID]->Fill(pT, dca[0]);
3608 fhDCA_Z_prim_MCID[mcID]->Fill(pT, dca[1]);
3611 fhDCA_XY_sec_MCID[mcID]->Fill(pT, dca[0]);
3612 fhDCA_Z_sec_MCID[mcID]->Fill(pT, dca[1]);
3620 if (fIDFFMode && ((!fJetPIDtask[0]->GetUseTPCCutMIGeo() && !fJetPIDtask[0]->GetUseTPCnclCut()) ||
3621 (survivedTPCCutMIGeo && fJetPIDtask[0]->GetUseTPCCutMIGeo()) ||
3622 (survivedTPCnclCut && fJetPIDtask[0]->GetUseTPCnclCut()))) {
3623 // NOTE: Just take particle fraction from first task (should anyway be the same for all tasks)
3624 Int_t pidWeightedSpecies = fJetPIDtask[0]->GetRandomParticleTypeAccordingToParticleFractions(pT, jetPt,
3625 centPercent, kTRUE);
3626 if (pidWeightedSpecies < 0 || pidWeightedSpecies >= AliPID::kSPECIES) {
3627 Printf("Failed to determine particle ID for track in jet (recCuts) -> ID FF histos not filled with this track!");
3628 Printf("Track details: trackPt %f, jetPt %f, centrality %f!", pT, jetPt, centPercent);
3631 // WARNING: The number of jets for the different species does not make sense -> One has to take
3632 // the number of jets for ALL particles. Thus, just do not fill the num of jet histos in order
3633 // not to get confused
3634 fIDFFHistosRecCuts[pidWeightedSpecies]->FillFF(aodtrack->Pt(), jetPt, kFALSE);
3638 // Efficiency, jets - detector level
3640 // Following lines are not needed - just keep other species (like casecades) - will end up in overflow bin
3641 // and only affect the efficiencies for all (i.e. not identified) what is desired!
3642 //if (mcID == AliPID::kUnknown)
3645 // Fill efficiency for reconstructed primaries
3646 if (!gentrack->IsPhysicalPrimary())
3649 Int_t iMother = gentrack->GetMother();
3651 continue; // Not a physical primary
3654 if (gentrack->Eta() > fTrackEtaMax || gentrack->Eta() < fTrackEtaMin)
3657 Double_t genZ = -1., genXi = -1.;
3658 Double_t genPt = gentrack->Pt();
3659 AliAnalysisTaskPID::GetJetTrackObservables(genPt, jetPt, genZ, genXi);
3661 Double_t measZ = -1., measXi = -1.;
3662 Double_t measPt = aodtrack->Pt();
3663 AliAnalysisTaskPID::GetJetTrackObservables(measPt, jetPt, measZ, measXi);
3666 // AliAODMCParticle->Charge() calls TParticlePDG->Charge(), which returns the charge in units of e0 / 3
3667 Double_t value[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), genPt, gentrack->Eta(), gentrack->Charge() / 3.,
3668 centPercent, jetPt, genZ, genXi };
3669 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
3670 if (!isPileUpJetPIDtask[i] &&
3671 ((!fJetPIDtask[i]->GetUseTPCCutMIGeo() && !fJetPIDtask[i]->GetUseTPCnclCut()) ||
3672 (survivedTPCCutMIGeo && fJetPIDtask[i]->GetUseTPCCutMIGeo()) ||
3673 (survivedTPCnclCut && fJetPIDtask[i]->GetUseTPCnclCut())))
3674 fJetPIDtask[i]->FillEfficiencyContainer(value, AliAnalysisTaskPID::kStepRecWithGenCuts);
3677 Double_t valueMeas[AliAnalysisTaskPID::kEffNumAxes] = { static_cast<Double_t>(mcID), measPt, aodtrack->Eta(), static_cast<Double_t>(aodtrack->Charge()),
3678 centPercent, jetPt, measZ, measXi };
3679 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
3680 if (!isPileUpJetPIDtask[i] &&
3681 ((!fJetPIDtask[i]->GetUseTPCCutMIGeo() && !fJetPIDtask[i]->GetUseTPCnclCut()) ||
3682 (survivedTPCCutMIGeo && fJetPIDtask[i]->GetUseTPCCutMIGeo()) ||
3683 (survivedTPCnclCut && fJetPIDtask[i]->GetUseTPCnclCut())))
3684 fJetPIDtask[i]->FillEfficiencyContainer(valueMeas, AliAnalysisTaskPID::kStepRecWithGenCutsMeasuredObs);
3691 if(fBckgMode && (ij==0)){
3692 if(fBckgType[0]!=kBckgNone)
3693 FillBckgHistos(fBckgType[0], fTracksRecCuts, fJetsRecCuts, jet,
3694 fFFBckgHisto0RecCuts,fQABckgHisto0RecCuts, fh1BckgMult0);
3695 if(fBckgType[1]!=kBckgNone)
3696 FillBckgHistos(fBckgType[1], fTracksRecCuts, fJetsRecCuts, jet,
3697 fFFBckgHisto1RecCuts,fQABckgHisto1RecCuts, fh1BckgMult1);
3698 if(fBckgType[2]!=kBckgNone)
3699 FillBckgHistos(fBckgType[2], fTracksRecCuts, fJetsRecCuts, jet,
3700 fFFBckgHisto2RecCuts,fQABckgHisto2RecCuts, fh1BckgMult1);
3701 if(fBckgType[3]!=kBckgNone)
3702 FillBckgHistos(fBckgType[3], fTracksRecCuts, fJetsRecCuts, jet,
3703 fFFBckgHisto3RecCuts,fQABckgHisto3RecCuts, fh1BckgMult2);
3704 if(fBckgType[4]!=kBckgNone)
3705 FillBckgHistos(fBckgType[4], fTracksRecCuts, fJetsRecCuts, jet,
3706 fFFBckgHisto4RecCuts,fQABckgHisto4RecCuts, fh1BckgMult3);
3707 } // end if(fBckgMode)
3710 if(fJSMode && (ij==0)) FillJetShape(jet, jettracklist, fProNtracksLeadingJet, fProDelRPtSum, fProDelR80pcPt);
3712 delete jettracklist;
3714 } // end: cut embedded ratio
3715 } // end: leading jet or all jets
3716 } // end: rec. jets after cuts
3717 } // end: QA, FF and intra-jet
3720 // ____ efficiency _______________________________
3722 if(fEffMode && (fJetTypeRecEff != kJetsUndef)){
3724 // arrays holding for each generated particle the reconstructed AOD track index & isPrimary flag, are initialized in AssociateGenRec(...) function
3728 // array holding for each reconstructed AOD track generated particle index, initialized in AssociateGenRec(...) function
3731 // ... and another set for secondaries from strange/non strange mothers (secondary MC tracks are stored in different lists)
3732 TArrayI indexAODTrSecNS;
3734 TArrayI indexMCTrSecNS;
3736 TArrayI indexAODTrSecS;
3738 TArrayI indexMCTrSecS;
3740 Int_t nTracksAODMCCharged = GetListOfTracks(fTracksAODMCCharged, kTrackAODMCCharged);
3741 if(fDebug>2)Printf("%s:%d selected AODMC tracks: %d ",(char*)__FILE__,__LINE__,nTracksAODMCCharged);
3743 Int_t nTracksAODMCChargedSecNS = GetListOfTracks(fTracksAODMCChargedSecNS, kTrackAODMCChargedSecNS);
3744 if(fDebug>2)Printf("%s:%d selected AODMC secondary tracks NS: %d ",(char*)__FILE__,__LINE__,nTracksAODMCChargedSecNS);
3746 Int_t nTracksAODMCChargedSecS = GetListOfTracks(fTracksAODMCChargedSecS, kTrackAODMCChargedSecS);
3747 if(fDebug>2)Printf("%s:%d selected AODMC secondary tracks S: %d ",(char*)__FILE__,__LINE__,nTracksAODMCChargedSecS);
3749 Int_t nTracksRecQualityCuts = GetListOfTracks(fTracksRecQualityCuts, kTrackAODQualityCuts);
3750 if(fDebug>2)Printf("%s:%d selected rec tracks quality after cuts, full acceptance/pt : %d ",(char*)__FILE__,__LINE__,nTracksRecQualityCuts);
3752 // associate gen and rec tracks, store indices in TArrays
3753 AssociateGenRec(fTracksAODMCCharged,fTracksRecQualityCuts,indexAODTr,indexMCTr,isGenPrim,fh2PtRecVsGenPrim);
3754 AssociateGenRec(fTracksAODMCChargedSecNS,fTracksRecQualityCuts,indexAODTrSecNS,indexMCTrSecNS,isGenSecNS,fh2PtRecVsGenSec);
3755 AssociateGenRec(fTracksAODMCChargedSecS,fTracksRecQualityCuts,indexAODTrSecS,indexMCTrSecS,isGenSecS,fh2PtRecVsGenSec);
3758 if(fQAMode&1) FillSingleTrackHistosRecGen(fQATrackHistosRecEffGen,fQATrackHistosRecEffRec,fTracksAODMCCharged,indexAODTr,isGenPrim);
3761 if(fQAMode&1) FillSingleTrackHistosRecGen(0x0,fQATrackHistosSecRecNS,fTracksAODMCChargedSecNS,indexAODTrSecNS,isGenSecNS);
3762 if(fQAMode&1) FillSingleTrackHistosRecGen(0x0,fQATrackHistosSecRecS,fTracksAODMCChargedSecS,indexAODTrSecS,isGenSecS);
3763 if(fQAMode&1) FillSingleTrackHistosRecGen(0x0,fQATrackHistosSecRecSsc,fTracksAODMCChargedSecS,indexAODTrSecS,isGenSecS,kTRUE);
3767 Double_t sumPtGenLeadingJetRecEff = 0;
3768 Double_t sumPtGenLeadingJetSec = 0;
3769 Double_t sumPtRecLeadingJetRecEff = 0;
3771 for(Int_t ij=0; ij<nRecEffJets; ++ij){ // jet loop
3773 AliAODJet* jet = (AliAODJet*)(fJetsRecEff->At(ij));
3775 Bool_t isBadJetGenPrim = kFALSE;
3776 Bool_t isBadJetGenSec = kFALSE;
3777 Bool_t isBadJetRec = kFALSE;
3780 if((ij==0) || !fOnlyLeadingJets){ // leading jets or all jets
3782 // for efficiency: gen tracks from pointing with gen/rec jet
3783 TList* jettracklistGenPrim = new TList();
3785 // if radius<0 -> trackRefs: collect gen tracks in wide radius + fill FF recEff rec histos with tracks contained in track refs
3786 // note : FF recEff gen histos will be somewhat useless in this approach
3788 if(GetFFRadius() >0)
3789 GetJetTracksPointing(fTracksAODMCCharged, jettracklistGenPrim, jet, GetFFRadius(), sumPtGenLeadingJetRecEff, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetGenPrim);
3791 GetJetTracksPointing(fTracksAODMCCharged, jettracklistGenPrim, jet, TMath::Abs(GetFFRadius())+0.2, sumPtGenLeadingJetRecEff, GetFFMinLTrackPt(), GetFFMaxTrackPt(), isBadJetGenPrim);
3793 TList* jettracklistGenSecNS = new TList();
3794 if(GetFFRadius() >0)
3795 GetJetTracksPointing(fTracksAODMCChargedSecNS, jettracklistGenSecNS, jet, GetFFRadius(), sumPtGenLeadingJetSec, GetFFMinLTrackPt() , GetFFMaxTrackPt(), isBadJetGenSec);
3797 GetJetTracksPointing(fTracksAODMCChargedSecNS, jettracklistGenSecNS, jet, TMath::Abs(GetFFRadius())+0.2, sumPtGenLeadingJetSec, GetFFMinLTrackPt() , GetFFMaxTrackPt(), isBadJetGenSec);
3799 TList* jettracklistGenSecS = new TList();
3800 if(GetFFRadius() >0)
3801 GetJetTracksPointing(fTracksAODMCChargedSecS, jettracklistGenSecS, jet, GetFFRadius(), sumPtGenLeadingJetSec, GetFFMinLTrackPt() , GetFFMaxTrackPt(), isBadJetGenSec);
3803 GetJetTracksPointing(fTracksAODMCChargedSecS, jettracklistGenSecS, jet, TMath::Abs(GetFFRadius())+0.2, sumPtGenLeadingJetSec, GetFFMinLTrackPt() , GetFFMaxTrackPt(), isBadJetGenSec);
3806 // bin efficiency in jet pt bins using rec tracks
3807 TList* jettracklistRec = new TList();
3808 if(GetFFRadius() >0) GetJetTracksPointing(fTracksRecCuts,jettracklistRec, jet, GetFFRadius(), sumPtRecLeadingJetRecEff, GetFFMinLTrackPt() , GetFFMaxTrackPt(), isBadJetRec);
3809 else GetJetTracksTrackrefs(jettracklistRec, jet, GetFFMinLTrackPt() , GetFFMaxTrackPt(), isBadJetRec);
3812 Double_t jetEta = jet->Eta();
3813 Double_t jetPhi = TVector2::Phi_0_2pi(jet->Phi());
3815 if(GetFFMinNTracks()>0 && jettracklistGenPrim->GetSize()<=GetFFMinNTracks()) isBadJetGenPrim = kTRUE;
3816 if(GetFFMinNTracks()>0 && jettracklistGenSecNS->GetSize()<=GetFFMinNTracks()) isBadJetGenSec = kTRUE;
3817 if(GetFFMinNTracks()>0 && jettracklistRec->GetSize()<=GetFFMinNTracks()) isBadJetRec = kTRUE;
3819 if(isBadJetRec) continue;
3821 if(fQAMode&2) fQAJetHistosRecEffLeading->FillJetQA( jetEta, jetPhi, sumPtGenLeadingJetRecEff );
3825 if(GetFFRadius()>0) FillJetTrackHistosRec(fFFHistosRecEffRec,jet,
3826 jettracklistGenPrim,fTracksAODMCCharged,fTracksRecQualityCuts,indexAODTr,isGenPrim,
3827 0,kFALSE,fJSMode,fProNtracksLeadingJetRecPrim,fProDelRPtSumRecPrim,fProDelR80pcPtRecPrim);
3829 else FillJetTrackHistosRec(fFFHistosRecEffRec,jet,
3830 jettracklistGenPrim,fTracksAODMCCharged,fTracksRecQualityCuts,indexAODTr,isGenPrim,
3831 jettracklistRec,kFALSE,fJSMode,fProNtracksLeadingJetRecPrim,fProDelRPtSumRecPrim,fProDelR80pcPtRecPrim);
3834 // secondaries: use jet pt from primaries
3835 if(GetFFRadius()>0) FillJetTrackHistosRec(fFFHistosSecRecNS,jet,
3836 jettracklistGenSecNS,fTracksAODMCChargedSecNS,fTracksRecQualityCuts, indexAODTrSecNS,isGenSecNS,
3837 0,kFALSE,fJSMode,fProNtracksLeadingJetRecSecNS,fProDelRPtSumRecSecNS);
3839 else FillJetTrackHistosRec(fFFHistosSecRecNS,jet,
3840 jettracklistGenSecNS,fTracksAODMCChargedSecNS,fTracksRecQualityCuts,indexAODTrSecNS,isGenSecNS,
3841 jettracklistRec,kFALSE,fJSMode,fProNtracksLeadingJetRecSecNS,fProDelRPtSumRecSecNS);
3843 if(GetFFRadius()>0) FillJetTrackHistosRec(fFFHistosSecRecS,jet,
3844 jettracklistGenSecS,fTracksAODMCChargedSecS,fTracksRecQualityCuts,indexAODTrSecS,isGenSecS,
3845 0,kFALSE,fJSMode,fProNtracksLeadingJetRecSecS,fProDelRPtSumRecSecS);
3847 else FillJetTrackHistosRec(fFFHistosSecRecS,jet,
3848 jettracklistGenSecS,fTracksAODMCChargedSecS,fTracksRecQualityCuts,indexAODTrSecS,isGenSecS,
3849 jettracklistRec,kFALSE,fJSMode,fProNtracksLeadingJetRecSecS,fProDelRPtSumRecSecS);
3851 if(GetFFRadius()>0) FillJetTrackHistosRec(fFFHistosSecRecSsc,jet,
3852 jettracklistGenSecS,fTracksAODMCChargedSecS,fTracksRecQualityCuts,indexAODTrSecS,isGenSecS,
3853 0,kTRUE,fJSMode,fProNtracksLeadingJetRecSecSsc,fProDelRPtSumRecSecSsc);
3855 else FillJetTrackHistosRec(fFFHistosSecRecSsc,jet,
3856 jettracklistGenSecS,fTracksAODMCChargedSecS,fTracksRecQualityCuts,indexAODTrSecS,isGenSecS,
3857 jettracklistRec,kTRUE,fJSMode,fProNtracksLeadingJetRecSecSsc,fProDelRPtSumRecSecSsc);
3860 delete jettracklistGenPrim;
3861 delete jettracklistGenSecNS;
3862 delete jettracklistGenSecS;
3863 delete jettracklistRec;
3866 if(fBckgMode && fFFMode){
3868 TList* perpjettracklistGen = new TList();
3869 TList* perpjettracklistGen1 = new TList();
3870 TList* perpjettracklistGen2 = new TList();
3872 //Double_t sumPtGenPerp = 0.;
3873 Double_t sumPtGenPerp1 = 0.;
3874 Double_t sumPtGenPerp2 = 0.;
3875 GetTracksTiltedwrpJetAxis(TMath::Pi()/2.,fTracksAODMCCharged, perpjettracklistGen1, jet, TMath::Abs(GetFFRadius()) , sumPtGenPerp1);
3876 GetTracksTiltedwrpJetAxis(-1*TMath::Pi()/2.,fTracksAODMCCharged, perpjettracklistGen2, jet, TMath::Abs(GetFFRadius()) ,
3879 perpjettracklistGen->AddAll(perpjettracklistGen1);
3880 perpjettracklistGen->AddAll(perpjettracklistGen2);
3881 //sumPtGenPerp = 0.5*(sumPtGenPerp1+sumPtGenPerp2);
3883 TList* perpjettracklistGenSecNS = new TList();
3884 TList* perpjettracklistGenSecNS1 = new TList();
3885 TList* perpjettracklistGenSecNS2 = new TList();
3887 //Double_t sumPtGenPerpNS = 0.;
3888 Double_t sumPtGenPerpNS1 = 0.;
3889 Double_t sumPtGenPerpNS2 = 0.;
3890 GetTracksTiltedwrpJetAxis(TMath::Pi()/2.,fTracksAODMCChargedSecNS, perpjettracklistGenSecNS1, jet, TMath::Abs(GetFFRadius()) ,
3892 GetTracksTiltedwrpJetAxis(-1*TMath::Pi()/2.,fTracksAODMCChargedSecNS, perpjettracklistGenSecNS2, jet, TMath::Abs(GetFFRadius()) ,
3895 perpjettracklistGenSecNS->AddAll(perpjettracklistGenSecNS1);
3896 perpjettracklistGenSecNS->AddAll(perpjettracklistGenSecNS2);
3897 //sumPtGenPerpNS = 0.5*(sumPtGenPerpNS1+sumPtGenPerpNS2);
3900 TList* perpjettracklistGenSecS = new TList();
3901 TList* perpjettracklistGenSecS1 = new TList();
3902 TList* perpjettracklistGenSecS2 = new TList();
3904 //Double_t sumPtGenPerpS = 0.;
3905 Double_t sumPtGenPerpS1 = 0.;
3906 Double_t sumPtGenPerpS2 = 0.;
3907 GetTracksTiltedwrpJetAxis(TMath::Pi()/2.,fTracksAODMCChargedSecS, perpjettracklistGenSecS1, jet, TMath::Abs(GetFFRadius()) ,
3909 GetTracksTiltedwrpJetAxis(-1*TMath::Pi()/2.,fTracksAODMCChargedSecS, perpjettracklistGenSecS2, jet, TMath::Abs(GetFFRadius()) ,
3912 perpjettracklistGenSecS->AddAll(perpjettracklistGenSecS1);
3913 perpjettracklistGenSecS->AddAll(perpjettracklistGenSecS2);
3914 //sumPtGenPerpS = 0.5*(sumPtGenPerpS1+sumPtGenPerpS2);
3917 if(perpjettracklistGen->GetSize() != perpjettracklistGen1->GetSize() + perpjettracklistGen2->GetSize()){
3918 cout<<" ERROR: perpjettracklistGen size "<<perpjettracklistGen->GetSize()<<" perp1 "<<perpjettracklistGen1->GetSize()
3919 <<" perp2 "<<perpjettracklistGen2->GetSize()<<endl;
3923 if(perpjettracklistGenSecNS->GetSize() != perpjettracklistGenSecNS1->GetSize() + perpjettracklistGenSecNS2->GetSize()){
3924 cout<<" ERROR: perpjettracklistGenSecNS size "<<perpjettracklistGenSecNS->GetSize()<<" perp1 "<<perpjettracklistGenSecNS1->GetSize()
3925 <<" perp2 "<<perpjettracklistGenSecNS2->GetSize()<<endl;
3929 if(perpjettracklistGenSecS->GetSize() != perpjettracklistGenSecS1->GetSize() + perpjettracklistGenSecS2->GetSize()){
3930 cout<<" ERROR: perpjettracklistGenSecS size "<<perpjettracklistGenSecS->GetSize()<<" perp1 "<<perpjettracklistGenSecS1->GetSize()
3931 <<" perp2 "<<perpjettracklistGenSecS2->GetSize()<<endl;
3936 FillJetTrackHistosRec(fFFBckgHisto0RecEffRec,jet,
3937 perpjettracklistGen,fTracksAODMCCharged,fTracksRecQualityCuts,indexAODTr,isGenPrim);
3939 FillJetTrackHistosRec(fFFBckgHisto0SecRecNS,jet,
3940 perpjettracklistGenSecNS,fTracksAODMCChargedSecNS,fTracksRecQualityCuts,indexAODTrSecNS,isGenSecNS);
3942 FillJetTrackHistosRec(fFFBckgHisto0SecRecS,jet,
3943 perpjettracklistGenSecS,fTracksAODMCChargedSecS,fTracksRecQualityCuts,indexAODTrSecS,isGenSecS);
3945 FillJetTrackHistosRec(fFFBckgHisto0SecRecSsc,jet,
3946 perpjettracklistGenSecS,fTracksAODMCChargedSecS,fTracksRecQualityCuts,indexAODTrSecS,isGenSecS,0,kTRUE);
3949 delete perpjettracklistGen;
3950 delete perpjettracklistGen1;
3951 delete perpjettracklistGen2;
3953 delete perpjettracklistGenSecNS;
3954 delete perpjettracklistGenSecNS1;
3955 delete perpjettracklistGenSecNS2;
3957 delete perpjettracklistGenSecS;
3958 delete perpjettracklistGenSecS1;
3959 delete perpjettracklistGenSecS2;
3966 //___________________
3968 fTracksRecCuts->Clear();
3969 fTracksRecCutsEfficiency->Clear();
3970 fTracksGen->Clear();
3971 fTracksAODMCCharged->Clear();
3972 fTracksAODMCChargedSecNS->Clear();
3973 fTracksAODMCChargedSecS->Clear();
3974 fTracksRecQualityCuts->Clear();
3977 fJetsRecCuts->Clear();
3979 fJetsRecEff->Clear();
3980 fJetsEmbedded->Clear();
3984 (fBckgType[0]==kBckgClusters || fBckgType[1]==kBckgClusters || fBckgType[2]==kBckgClusters || fBckgType[3]==kBckgClusters || fBckgType[4]==kBckgClusters ||
3985 fBckgType[0]==kBckgClustersOutLeading || fBckgType[1]==kBckgClustersOutLeading || fBckgType[2]==kBckgClustersOutLeading ||
3986 fBckgType[3]==kBckgClustersOutLeading || fBckgType[4]==kBckgClustersOutLeading)){
3988 fBckgJetsRec->Clear();
3989 fBckgJetsRecCuts->Clear();
3990 fBckgJetsGen->Clear();
3995 PostData(1, fCommonHistList);
3997 if (fUseJetPIDtask) {
3998 for (Int_t i = 0; i < fNumJetPIDtasks; i++) {
3999 fJetPIDtask[i]->PostOutputData();
4003 if (fUseInclusivePIDtask) {
4004 for (Int_t i = 0; i < fNumInclusivePIDtasks; i++) {
4005 fInclusivePIDtask[i]->PostOutputData();
4010 //______________________________________________________________
4011 void AliAnalysisTaskIDFragmentationFunction::Terminate(Option_t *)
4015 if(fDebug > 1) printf("AliAnalysisTaskIDFragmentationFunction::Terminate() \n");
4018 //_________________________________________________________________________________
4019 Int_t AliAnalysisTaskIDFragmentationFunction::GetListOfTracks(TList *list, Int_t type)
4021 // fill list of tracks selected according to type
4023 if(fDebug > 2) Printf("%s:%d Selecting tracks with %d", (char*)__FILE__,__LINE__,type);
4026 if(fDebug>1) Printf("%s:%d no input list", (char*)__FILE__,__LINE__);
4030 if(!fAOD) return -1;
4032 if(!fAOD->GetTracks()) return 0;
4034 if(type==kTrackUndef) return 0;
4038 if(type==kTrackAODExtraCuts || type==kTrackAODExtraonlyCuts || type==kTrackAODExtra || type==kTrackAODExtraonly){
4040 TClonesArray *aodExtraTracks = dynamic_cast<TClonesArray*>(fAOD->FindListObject("aodExtraTracks"));
4041 if(!aodExtraTracks)return iCount;
4042 for(int it =0; it<aodExtraTracks->GetEntries(); it++) {
4043 AliVParticle *track = dynamic_cast<AliVParticle*> ((*aodExtraTracks)[it]);
4044 if (!track) continue;
4046 AliAODTrack *tr = dynamic_cast<AliAODTrack*> (track);
4049 if(type==kTrackAODExtraCuts || type==kTrackAODExtraonlyCuts){
4051 if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask))) continue;
4053 if(tr->Eta() < fTrackEtaMin || tr->Eta() > fTrackEtaMax) continue;
4054 if(tr->Phi() < fTrackPhiMin || tr->Phi() > fTrackPhiMax) continue;
4055 if(tr->Pt() < fTrackPtCut) continue;
4063 if(type==kTrackAODCuts || type==kTrackAODQualityCuts || type==kTrackAOD || type==kTrackAODExtraCuts || type==kTrackAODExtra){
4065 // all rec. tracks, esd filter mask, eta range
4067 for(Int_t it=0; it<fAOD->GetNumberOfTracks(); ++it){
4068 AliAODTrack *tr = (AliAODTrack*)fAOD->GetTrack(it);
4070 if(type == kTrackAODCuts || type==kTrackAODQualityCuts || type==kTrackAODExtraCuts){
4072 if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask))) continue;
4073 if(type == kTrackAODCuts){
4074 if(tr->Eta() < fTrackEtaMin || tr->Eta() > fTrackEtaMax) continue;
4075 if(tr->Phi() < fTrackPhiMin || tr->Phi() > fTrackPhiMax) continue;
4076 if(tr->Pt() < fTrackPtCut) continue;
4083 else if (type==kTrackKineAll || type==kTrackKineCharged || type==kTrackKineChargedAcceptance){
4084 // kine particles, all or rather charged
4085 if(!fMCEvent) return iCount;
4087 for(Int_t it=0; it<fMCEvent->GetNumberOfTracks(); ++it){
4088 AliMCParticle* part = (AliMCParticle*) fMCEvent->GetTrack(it);
4090 if(type == kTrackKineCharged || type == kTrackKineChargedAcceptance){
4091 if(part->Charge()==0) continue;
4093 if(type == kTrackKineChargedAcceptance &&
4094 ( part->Eta() < fTrackEtaMin
4095 || part->Eta() > fTrackEtaMax
4096 || part->Phi() < fTrackPhiMin
4097 || part->Phi() > fTrackPhiMax
4098 || part->Pt() < fTrackPtCut)) continue;
4105 else if (type==kTrackAODMCCharged || type==kTrackAODMCAll || type==kTrackAODMCChargedAcceptance || type==kTrackAODMCChargedSecNS || type==kTrackAODMCChargedSecS) {
4106 // MC particles (from AOD), physical primaries, all or rather charged or rather charged within acceptance
4107 if(!fAOD) return -1;
4109 TClonesArray *tca = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
4110 if(!tca)return iCount;
4112 for(int it=0; it<tca->GetEntriesFast(); ++it){
4113 AliAODMCParticle *part = dynamic_cast<AliAODMCParticle*>(tca->At(it));
4115 if(type != kTrackAODMCChargedSecNS && type != kTrackAODMCChargedSecS && !part->IsPhysicalPrimary())continue;
4116 if((type == kTrackAODMCChargedSecNS || type == kTrackAODMCChargedSecS) && part->IsPhysicalPrimary())continue;
4118 if (type==kTrackAODMCCharged || type==kTrackAODMCChargedAcceptance || type==kTrackAODMCChargedSecNS || type==kTrackAODMCChargedSecS){
4119 if(part->Charge()==0) continue;
4121 if(type==kTrackAODMCChargedSecNS || type==kTrackAODMCChargedSecS){
4122 Bool_t isFromStrange = kFALSE;
4123 Int_t iMother = part->GetMother();
4125 AliAODMCParticle *partM = dynamic_cast<AliAODMCParticle*>(tca->At(iMother));
4126 if(!partM) continue;
4128 Int_t codeM = TMath::Abs(partM->GetPdgCode());
4129 Int_t mfl = Int_t (codeM/ TMath::Power(10, Int_t(TMath::Log10(codeM))));
4131 if (mfl == 3 && codeM != 3) isFromStrange = kTRUE;
4132 if(codeM == 130) isFromStrange = kTRUE; // K0 long
4133 if(part->IsSecondaryFromMaterial()) isFromStrange = kFALSE; // strange resonances from hadronic showers ?
4136 // cout<<" mfl "<<mfl<<" codeM "<<partM->GetPdgCode()<<" code this track "<<part->GetPdgCode()<<endl;
4137 // cout<<" index this track "<<it<<" index daughter 0 "<<partM->GetDaughter(0)<<" 1 "<<partM->GetDaughter(1)<<endl;
4140 if(type==kTrackAODMCChargedSecNS && isFromStrange) continue;
4141 if(type==kTrackAODMCChargedSecS && !isFromStrange) continue;
4145 if(type==kTrackAODMCChargedAcceptance &&
4146 ( part->Eta() > fTrackEtaMax
4147 || part->Eta() < fTrackEtaMin
4148 || part->Phi() > fTrackPhiMax
4149 || part->Phi() < fTrackPhiMin
4150 || part->Pt() < fTrackPtCut)) continue;
4162 // _______________________________________________________________________________
4163 Int_t AliAnalysisTaskIDFragmentationFunction::GetListOfJets(TList *list, Int_t type)
4165 // fill list of jets selected according to type
4168 if(fDebug>1) Printf("%s:%d no input list", (char*)__FILE__,__LINE__);
4172 if(type == kJetsRec || type == kJetsRecAcceptance){ // reconstructed jets
4174 if(fBranchRecJets.Length()==0){
4175 Printf("%s:%d no rec jet branch specified", (char*)__FILE__,__LINE__);
4176 if(fDebug>1)fAOD->Print();
4180 TClonesArray *aodRecJets = 0;
4181 if(fBranchRecJets.Length()) aodRecJets = dynamic_cast<TClonesArray*>(fAODJets->FindListObject(fBranchRecJets.Data()));
4182 if(!aodRecJets) aodRecJets = dynamic_cast<TClonesArray*>(fAODJets->GetList()->FindObject(fBranchRecJets.Data()));
4183 if(fAODExtension&&!aodRecJets) aodRecJets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fBranchRecJets.Data()));
4186 if(fBranchRecJets.Length()) Printf("%s:%d no reconstructed jet array with name %s in AOD", (char*)__FILE__,__LINE__,fBranchRecJets.Data());
4187 if(fDebug>1)fAOD->Print();
4191 // Reorder jet pt and fill new temporary AliAODJet objects
4194 for(Int_t ij=0; ij<aodRecJets->GetEntries(); ++ij){
4196 AliAODJet *tmp = dynamic_cast<AliAODJet*>(aodRecJets->At(ij));
4199 if( tmp->Pt() < fJetPtCut ) continue;
4200 if( type == kJetsRecAcceptance &&
4201 ( tmp->Eta() < fJetEtaMin
4202 || tmp->Eta() > fJetEtaMax
4203 || tmp->Phi() < fJetPhiMin
4204 || tmp->Phi() > fJetPhiMax )) continue;
4215 else if(type == kJetsKine || type == kJetsKineAcceptance){
4221 if(fDebug>1) Printf("%s:%d no mcEvent",(char*)__FILE__,__LINE__);
4225 AliGenEventHeader* genHeader = fMCEvent->GenEventHeader();
4226 AliGenPythiaEventHeader* pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
4227 AliGenHijingEventHeader* hijingGenHeader = 0x0;
4229 if(!pythiaGenHeader){
4230 hijingGenHeader = dynamic_cast<AliGenHijingEventHeader*>(genHeader);
4232 if(!hijingGenHeader){
4233 Printf("%s:%d no pythiaGenHeader or hijingGenHeader found", (char*)__FILE__,__LINE__);
4236 TLorentzVector mom[4];
4238 hijingGenHeader->GetJets(mom[0], mom[1], mom[2], mom[3]);
4240 for(Int_t i=0; i<2; ++i){
4241 if(!mom[i].Pt()) continue;
4242 jet[i] = new AliAODJet(mom[i]);
4244 if( type == kJetsKineAcceptance &&
4245 ( jet[i]->Eta() < fJetEtaMin
4246 || jet[i]->Eta() > fJetEtaMax
4247 || jet[i]->Phi() < fJetPhiMin
4248 || jet[i]->Phi() > fJetPhiMax )) continue;
4258 // fetch the pythia generated jets
4259 for(int ip=0; ip<pythiaGenHeader->NTriggerJets(); ++ip){
4262 AliAODJet *jet = new AliAODJet();
4263 pythiaGenHeader->TriggerJet(ip, p);
4264 jet->SetPxPyPzE(p[0], p[1], p[2], p[3]);
4266 if( type == kJetsKineAcceptance &&
4267 ( jet->Eta() < fJetEtaMin
4268 || jet->Eta() > fJetEtaMax
4269 || jet->Phi() < fJetPhiMin
4270 || jet->Phi() > fJetPhiMax )) continue;
4278 else if(type == kJetsGen || type == kJetsGenAcceptance ){
4280 if(fBranchGenJets.Length()==0){
4281 if(fDebug>1) Printf("%s:%d no gen jet branch specified", (char*)__FILE__,__LINE__);
4285 TClonesArray *aodGenJets = 0;
4286 if(fBranchGenJets.Length()) aodGenJets = dynamic_cast<TClonesArray*>(fAODJets->FindListObject(fBranchGenJets.Data()));
4287 if(!aodGenJets) aodGenJets = dynamic_cast<TClonesArray*>(fAODJets->GetList()->FindObject(fBranchGenJets.Data()));
4288 if(fAODExtension&&!aodGenJets) aodGenJets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fBranchGenJets.Data()));
4292 if(fBranchGenJets.Length()) Printf("%s:%d Generated jet branch %s not found",(char*)__FILE__,__LINE__,fBranchGenJets.Data());
4294 if(fDebug>1)fAOD->Print();
4300 for(Int_t ig=0; ig<aodGenJets->GetEntries(); ++ig){
4302 AliAODJet *tmp = dynamic_cast<AliAODJet*>(aodGenJets->At(ig));
4305 if( tmp->Pt() < fJetPtCut ) continue;
4306 if( type == kJetsGenAcceptance &&
4307 ( tmp->Eta() < fJetEtaMin
4308 || tmp->Eta() > fJetEtaMax
4309 || tmp->Phi() < fJetPhiMin
4310 || tmp->Phi() > fJetPhiMax )) continue;
4318 else if(type == kJetsEmbedded){ // embedded jets
4320 if(fBranchEmbeddedJets.Length()==0){
4321 Printf("%s:%d no embedded jet branch specified", (char*)__FILE__,__LINE__);
4322 if(fDebug>1)fAOD->Print();
4326 TClonesArray *aodEmbeddedJets = 0;
4327 if(fBranchEmbeddedJets.Length()) aodEmbeddedJets = dynamic_cast<TClonesArray*>(fAODJets->FindListObject(fBranchEmbeddedJets.Data()));
4328 if(!aodEmbeddedJets) aodEmbeddedJets = dynamic_cast<TClonesArray*>(fAODJets->GetList()->FindObject(fBranchEmbeddedJets.Data()));
4329 if(fAODExtension&&!aodEmbeddedJets) aodEmbeddedJets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fBranchEmbeddedJets.Data()));
4331 if(!aodEmbeddedJets){
4332 if(fBranchEmbeddedJets.Length()) Printf("%s:%d no reconstructed jet array with name %s in AOD", (char*)__FILE__,__LINE__,fBranchEmbeddedJets.Data());
4333 if(fDebug>1)fAOD->Print();
4337 // Reorder jet pt and fill new temporary AliAODJet objects
4338 Int_t nEmbeddedJets = 0;
4340 for(Int_t ij=0; ij<aodEmbeddedJets->GetEntries(); ++ij){
4342 AliAODJet *tmp = dynamic_cast<AliAODJet*>(aodEmbeddedJets->At(ij));
4345 if( tmp->Pt() < fJetPtCut ) continue;
4346 if( tmp->Eta() < fJetEtaMin
4347 || tmp->Eta() > fJetEtaMax
4348 || tmp->Phi() < fJetPhiMin
4349 || tmp->Phi() > fJetPhiMax ) continue;
4357 return nEmbeddedJets;
4360 if(fDebug>0)Printf("%s:%d no such type %d",(char*)__FILE__,__LINE__,type);
4365 // ___________________________________________________________________________________
4366 Int_t AliAnalysisTaskIDFragmentationFunction::GetListOfBckgJets(TList *list, Int_t type)
4368 // fill list of bgr clusters selected according to type
4370 if(type == kJetsRec || type == kJetsRecAcceptance){ // reconstructed jets
4372 if(fBranchRecBckgClusters.Length()==0){
4373 Printf("%s:%d no rec jet branch specified", (char*)__FILE__,__LINE__);
4374 if(fDebug>1)fAOD->Print();
4378 TClonesArray *aodRecJets = 0;
4379 if(fBranchRecBckgClusters.Length()) aodRecJets = dynamic_cast<TClonesArray*>(fAODJets->FindListObject(fBranchRecBckgClusters.Data()));
4380 if(!aodRecJets) aodRecJets = dynamic_cast<TClonesArray*>(fAODJets->GetList()->FindObject(fBranchRecBckgClusters.Data()));
4381 if(fAODExtension&&!aodRecJets) aodRecJets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fBranchRecBckgClusters.Data()));
4384 if(fBranchRecBckgClusters.Length()) Printf("%s:%d no reconstructed jet array with name %s in AOD", (char*)__FILE__,__LINE__,fBranchRecBckgClusters.Data());
4385 if(fDebug>1)fAOD->Print();
4389 // Reorder jet pt and fill new temporary AliAODJet objects
4392 for(Int_t ij=0; ij<aodRecJets->GetEntries(); ++ij){
4394 AliAODJet *tmp = dynamic_cast<AliAODJet*>(aodRecJets->At(ij));
4397 // if( tmp->Pt() < fJetPtCut ) continue; // no pt cut on bckg clusters !
4398 if( type == kJetsRecAcceptance &&
4399 ( tmp->Eta() < fJetEtaMin
4400 || tmp->Eta() > fJetEtaMax
4401 || tmp->Phi() < fJetPhiMin
4402 || tmp->Phi() > fJetPhiMax )) continue;
4416 // MC clusters still Under construction
4422 // _________________________________________________________________________________________________________
4423 void AliAnalysisTaskIDFragmentationFunction::SetProperties(THnSparse* h, Int_t dim, const char** labels)
4425 // Set properties of THnSparse
4427 for(Int_t i=0; i<dim; i++){
4428 h->GetAxis(i)->SetTitle(labels[i]);
4429 h->GetAxis(i)->SetTitleColor(1);
4433 // __________________________________________________________________________________________
4434 void AliAnalysisTaskIDFragmentationFunction::SetProperties(TH1* h,const char* x, const char* y)
4436 //Set properties of histos (x and y title)
4440 h->GetXaxis()->SetTitleColor(1);
4441 h->GetYaxis()->SetTitleColor(1);
4444 // _________________________________________________________________________________________________________
4445 void AliAnalysisTaskIDFragmentationFunction::SetProperties(TH1* h,const char* x, const char* y, const char* z)
4447 //Set properties of histos (x,y and z title)
4452 h->GetXaxis()->SetTitleColor(1);
4453 h->GetYaxis()->SetTitleColor(1);
4454 h->GetZaxis()->SetTitleColor(1);
4457 // ________________________________________________________________________________________________________________________________________________________
4458 void AliAnalysisTaskIDFragmentationFunction::GetJetTracksPointing(TList* inputlist, TList* outputlist, const AliAODJet* jet,
4459 Double_t radius, Double_t& sumPt, Double_t minPtL, Double_t maxPt, Bool_t& isBadPt)
4461 // fill list of tracks in cone around jet axis
4464 Bool_t isBadMaxPt = kFALSE;
4465 Bool_t isBadMinPt = kTRUE;
4468 jet->PxPyPz(jetMom);
4469 TVector3 jet3mom(jetMom);
4471 for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++){
4473 AliVParticle* track = dynamic_cast<AliVParticle*>(inputlist->At(itrack));
4475 Double_t trackMom[3];
4476 track->PxPyPz(trackMom);
4477 TVector3 track3mom(trackMom);
4479 Double_t dR = jet3mom.DeltaR(track3mom);
4483 outputlist->Add(track);
4485 sumPt += track->Pt();
4487 if(maxPt>0 && track->Pt()>maxPt) isBadMaxPt = kTRUE;
4488 if(minPtL>0 && track->Pt()>minPtL) isBadMinPt = kFALSE;
4493 if(minPtL>0 && isBadMinPt) isBadPt = kTRUE;
4494 if(maxPt>0 && isBadMaxPt) isBadPt = kTRUE;
4499 // _________________________________________________________________________________________________________________________________________________________________
4500 void AliAnalysisTaskIDFragmentationFunction::GetJetTracksTrackrefs(TList* list, const AliAODJet* jet, Double_t minPtL, Double_t maxPt, Bool_t& isBadPt)
4502 // list of jet tracks from trackrefs
4504 Int_t nTracks = jet->GetRefTracks()->GetEntriesFast();
4506 Bool_t isBadMaxPt = kFALSE;
4507 Bool_t isBadMinPt = kTRUE;
4509 for(Int_t itrack=0; itrack<nTracks; itrack++) {
4511 AliVParticle* track = dynamic_cast<AliVParticle*>(jet->GetRefTracks()->At(itrack));
4513 AliError("expected ref track not found ");
4517 if(track->Pt() < fTrackPtCut) continue; // track refs may contain low pt cut (bug in AliFastJetInput)
4518 if(maxPt>0 && track->Pt()>maxPt) isBadMaxPt = kTRUE;
4519 if(minPtL>0 && track->Pt()>minPtL) isBadMinPt = kFALSE;
4525 if(minPtL>0 && isBadMinPt) isBadPt = kTRUE;
4526 if(maxPt>0 && isBadMaxPt) isBadPt = kTRUE;
4531 // _ ________________________________________________________________________________________________________________________________
4532 void AliAnalysisTaskIDFragmentationFunction::AssociateGenRec(TList* tracksAODMCCharged,TList* tracksRec, TArrayI& indexAODTr,TArrayI& indexMCTr,
4533 TArrayS& isRefGen,TH2F* fh2PtRecVsGen)
4535 // associate generated and reconstructed tracks, fill TArrays of list indices
4537 Int_t nTracksRec = tracksRec->GetSize();
4538 Int_t nTracksGen = tracksAODMCCharged->GetSize();
4539 TClonesArray *tca = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
4542 if(!nTracksGen) return;
4546 indexAODTr.Set(nTracksGen);
4547 indexMCTr.Set(nTracksRec);
4548 isRefGen.Set(nTracksGen);
4550 indexAODTr.Reset(-1);
4551 indexMCTr.Reset(-1);
4554 // loop over reconstructed tracks, get generated track
4556 for(Int_t iRec=0; iRec<nTracksRec; iRec++){
4558 AliAODTrack* rectrack = dynamic_cast<AliAODTrack*>(tracksRec->At(iRec));
4559 if(!rectrack)continue;
4560 Int_t label = TMath::Abs(rectrack->GetLabel());
4562 // find MC track in our list
4563 AliAODMCParticle* gentrack = dynamic_cast<AliAODMCParticle*> (tca->At(label));
4565 Int_t listIndex = -1;
4566 if(gentrack) listIndex = tracksAODMCCharged->IndexOf(gentrack);
4570 indexAODTr[listIndex] = iRec;
4571 indexMCTr[iRec] = listIndex;
4576 // define reference sample of primaries/secondaries (for reconstruction efficiency / contamination)
4578 for(Int_t iGen=0; iGen<nTracksGen; iGen++){
4580 AliAODMCParticle* gentrack = dynamic_cast<AliAODMCParticle*> (tracksAODMCCharged->At(iGen));
4581 if(!gentrack)continue;
4582 Int_t pdg = gentrack->GetPdgCode();
4584 // 211 - pi, 2212 - proton, 321 - Kaon, 11 - electron, 13 - muon
4585 if(TMath::Abs(pdg) == 211 || TMath::Abs(pdg) == 2212 || TMath::Abs(pdg) == 321 ||
4586 TMath::Abs(pdg) == 11 || TMath::Abs(pdg) == 13){
4588 isRefGen[iGen] = kTRUE;
4590 Int_t iRec = indexAODTr[iGen]; // can be -1 if no good reconstructed track
4593 Float_t genPt = gentrack->Pt();
4594 AliAODTrack* vt = dynamic_cast<AliAODTrack*>(tracksRec->At(iRec));
4596 Float_t recPt = vt->Pt();
4597 fh2PtRecVsGen->Fill(genPt,recPt);
4604 // _____________________________________________________________________________________________________________________________________________
4605 void AliAnalysisTaskIDFragmentationFunction::FillSingleTrackHistosRecGen(AliFragFuncQATrackHistos* trackQAGen, AliFragFuncQATrackHistos* trackQARec, TList* tracksGen,
4606 const TArrayI& indexAODTr, const TArrayS& isRefGen, Bool_t scaleStrangeness){
4608 // fill QA for single track reconstruction efficiency
4610 Int_t nTracksGen = tracksGen->GetSize();
4612 if(!nTracksGen) return;
4614 for(Int_t iGen=0; iGen<nTracksGen; iGen++){
4616 if(isRefGen[iGen] != 1) continue; // select primaries
4618 AliAODMCParticle* gentrack = dynamic_cast<AliAODMCParticle*> (tracksGen->At(iGen));
4619 if(!gentrack) continue;
4620 Double_t ptGen = gentrack->Pt();
4621 Double_t etaGen = gentrack->Eta();
4622 Double_t phiGen = TVector2::Phi_0_2pi(gentrack->Phi());
4624 // apply same acc & pt cuts as for FF GetMCStrangenessFactorCMS
4626 if(etaGen < fTrackEtaMin || etaGen > fTrackEtaMax) continue;
4627 if(phiGen < fTrackPhiMin || phiGen > fTrackPhiMax) continue;
4628 if(ptGen < fTrackPtCut) continue;
4630 if(trackQAGen) trackQAGen->FillTrackQA(etaGen, phiGen, ptGen);
4632 Int_t iRec = indexAODTr[iGen]; // can be -1 if no good reconstructed track
4634 if(iRec>=0 && trackQARec){
4635 if(scaleStrangeness){
4636 //Double_t weight = GetMCStrangenessFactor(ptGen);
4637 Double_t weight = GetMCStrangenessFactorCMS(gentrack);
4638 trackQARec->FillTrackQA(etaGen, phiGen, ptGen, kFALSE, 0, kTRUE, weight);
4640 else trackQARec->FillTrackQA(etaGen, phiGen, ptGen);
4645 // ______________________________________________________________________________________________________________________________________________________
4647 void AliAnalysisTaskIDFragmentationFunction::FillJetTrackHistosRec(AliFragFuncHistos* ffhistRec, AliAODJet* jet,
4648 TList* jetTrackList, const TList* tracksGen, const TList* tracksRec, const TArrayI& indexAODTr,
4649 const TArrayS& isRefGen, TList* jetTrackListTR, Bool_t scaleStrangeness,
4650 Bool_t fillJS, TProfile* hProNtracksLeadingJet, TProfile** hProDelRPtSum, TProfile* hProDelR80pcPt)
4652 // fill objects for jet track reconstruction efficiency or secondaries contamination
4653 // arguments histGen/histRec can be of different type: AliFragFuncHistos*, AliFragFuncIntraJetHistos*, ...
4654 // jetTrackListTR pointer: track refs if not NULL
4657 // ensure proper normalization, even for secondaries
4658 Double_t jetPtRec = jet->Pt();
4659 ffhistRec->FillFF(-1, jetPtRec, kTRUE);
4661 Int_t nTracksJet = jetTrackList->GetSize(); // list with AODMC tracks
4662 if(nTracksJet == 0) return;
4664 TList* listRecTracks = new TList();
4665 listRecTracks->Clear();
4667 for(Int_t iTr=0; iTr<nTracksJet; iTr++){ // jet tracks loop
4669 AliAODMCParticle* gentrack = dynamic_cast<AliAODMCParticle*> (jetTrackList->At(iTr));
4670 if(!gentrack)continue;
4671 // find jet track in gen tracks list
4672 Int_t iGen = tracksGen->IndexOf(gentrack);
4675 if(fDebug>0) Printf("%s:%d gen jet track not found ",(char*)__FILE__,__LINE__);
4679 if(isRefGen[iGen] != 1) continue; // select primaries
4681 Double_t ptGen = gentrack->Pt();
4682 Double_t etaGen = gentrack->Eta();
4683 Double_t phiGen = TVector2::Phi_0_2pi(gentrack->Phi());
4685 // gen level acc & pt cuts - skip in case of track refs
4686 if(!jetTrackListTR && (etaGen < fTrackEtaMin || etaGen > fTrackEtaMax)) continue;
4687 if(!jetTrackListTR && (phiGen < fTrackPhiMin || phiGen > fTrackPhiMax)) continue;
4688 if(!jetTrackListTR && ptGen < fTrackPtCut) continue;
4691 Double_t ptRec = -1;
4693 Int_t iRec = indexAODTr[iGen]; // can be -1 if no good reconstructed track
4694 Bool_t isRec = (iRec>=0) ? kTRUE : kFALSE;
4696 Bool_t isJetTrack = kFALSE;
4697 if(!jetTrackListTR) isJetTrack = kTRUE; // skip trackRefs check for tracks in ideal cone
4701 AliAODTrack* rectrack = dynamic_cast<AliAODTrack*> (tracksRec->At(iRec));
4702 if(!rectrack) continue;
4704 ptRec = rectrack->Pt();
4707 Int_t iRecTR = jetTrackListTR->IndexOf(rectrack);
4708 if(iRecTR >=0 ) isJetTrack = kTRUE; // rec tracks assigned to jet
4713 Double_t trackPt = ptRec;
4714 Bool_t incrementJetPt = kFALSE;
4716 if(scaleStrangeness){
4717 //Double_t weight = GetMCStrangenessFactor(ptGen);
4718 Double_t weight = GetMCStrangenessFactorCMS(gentrack);
4720 ffhistRec->FillFF( trackPt, jetPtRec, incrementJetPt, 0, kTRUE, weight );
4723 ffhistRec->FillFF( trackPt, jetPtRec, incrementJetPt );
4726 listRecTracks->Add(rectrack);
4733 if(fillJS) FillJetShape(jet,listRecTracks,hProNtracksLeadingJet, hProDelRPtSum, hProDelR80pcPt,0,0,scaleStrangeness);
4735 delete listRecTracks;
4739 // _____________________________________________________________________________________________________________________________________________________________________
4740 void AliAnalysisTaskIDFragmentationFunction::GetTracksTiltedwrpJetAxis(Float_t alpha, TList* inputlist, TList* outputlist, const AliAODJet* jet, Double_t radius,Double_t& sumPt)
4742 // List of tracks in cone perpendicular to the jet azimuthal direction
4745 jet->PxPyPz(jetMom);
4747 TVector3 jet3mom(jetMom);
4748 // Rotate phi and keep eta unchanged
4749 Double_t etaTilted = jet3mom.Eta();
4750 Double_t phiTilted = TVector2::Phi_0_2pi(jet3mom.Phi()) + alpha;
4751 if(phiTilted > 2*TMath::Pi()) phiTilted = phiTilted - 2*TMath::Pi();
4753 for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++){
4756 if( fUseExtraTracksBgr != 1){
4757 if(AliAODTrack* trackAOD = dynamic_cast<AliAODTrack*> (inputlist->At(itrack))){
4758 if(fUseExtraTracksBgr == 0 && ((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
4759 if(fUseExtraTracksBgr == -1 && !((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
4763 AliVParticle* track = dynamic_cast<AliVParticle*>(inputlist->At(itrack));
4765 Double_t trackMom[3];
4766 track->PxPyPz(trackMom);
4767 TVector3 track3mom(trackMom);
4769 Double_t deta = track3mom.Eta() - etaTilted;
4770 Double_t dphi = TMath::Abs(track3mom.Phi() - phiTilted);
4771 if (dphi > TMath::Pi()) dphi = 2. * TMath::Pi() - dphi;
4772 Double_t dR = TMath::Sqrt(deta * deta + dphi * dphi);
4776 outputlist->Add(track);
4777 sumPt += track->Pt();
4783 // ________________________________________________________________________________________________________________________________________________________
4784 void AliAnalysisTaskIDFragmentationFunction::GetTracksTiltedwrpJetAxisWindow(Float_t alpha, TList* inputlist, TList* outputlist, const AliAODJet* jet, Double_t radius,Double_t& sumPt,Double_t &normFactor)
4786 // List of tracks in cone perpendicular to the jet azimuthal direction
4789 jet->PxPyPz(jetMom);
4791 TVector3 jet3mom(jetMom);
4792 // Rotate phi and keep eta unchanged
4793 Double_t etaTilted = jet3mom.Eta();
4794 Double_t phiTilted = TVector2::Phi_0_2pi(jet3mom.Phi()) + alpha;
4795 if(phiTilted > 2*TMath::Pi()) phiTilted = phiTilted - 2*TMath::Pi();
4797 for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++)
4801 if( fUseExtraTracksBgr != 1){
4802 if(AliAODTrack* trackAOD = dynamic_cast<AliAODTrack*> (inputlist->At(itrack))){
4803 if(fUseExtraTracksBgr == 0 && ((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
4804 if(fUseExtraTracksBgr == -1 && !((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
4808 AliVParticle* track = dynamic_cast<AliVParticle*>(inputlist->At(itrack));
4810 Float_t trackEta = track->Eta();
4811 Float_t trackPhi = track->Phi();
4813 if( ( phiTilted-radius >= 0 ) && ( phiTilted+radius <= 2*TMath::Pi()))
4815 if((trackPhi<=phiTilted+radius) &&
4816 (trackPhi>=phiTilted-radius) &&
4817 (trackEta<=fTrackEtaMax)&&(trackEta>=fTrackEtaMin)) // 0.9 and - 0.9
4819 outputlist->Add(track);
4820 sumPt += track->Pt();
4823 else if( phiTilted-radius < 0 )
4825 if((( trackPhi < phiTilted+radius ) ||
4826 ( trackPhi > 2*TMath::Pi()-(radius-phiTilted) )) &&
4827 (( trackEta <= fTrackEtaMax ) && ( trackEta >= fTrackEtaMin )))
4829 outputlist->Add(track);
4830 sumPt += track->Pt();
4833 else if( phiTilted+radius > 2*TMath::Pi() )
4835 if((( trackPhi > phiTilted-radius ) ||
4836 ( trackPhi < phiTilted+radius-2*TMath::Pi() )) &&
4837 (( trackEta <= fTrackEtaMax ) && ( trackEta >= fTrackEtaMin )))
4839 outputlist->Add(track);
4840 sumPt += track->Pt();
4845 // Jet area - Temporarily added should be modified with the proper jet area value
4846 Float_t areaJet = CalcJetArea(etaTilted,radius);
4847 Float_t areaTilted = 2*radius*(fTrackEtaMax-fTrackEtaMin);
4849 normFactor = (Float_t) 1. / (areaJet / areaTilted);
4854 // ________________________________________________________________________________________________________________________________________________________
4855 void AliAnalysisTaskIDFragmentationFunction::GetTracksOutOfNJets(Int_t nCases, TList* inputlist, TList* outputlist, const TList* jetlist,
4858 // List of tracks outside cone around N jet axis
4859 // Particles taken randomly
4862 // Int_t nj = jetlist->GetSize();
4863 Float_t rc = TMath::Abs(GetFFRadius());
4864 Float_t rcl = GetFFBckgRadius();
4866 // Estimate jet and background areas
4867 Float_t* areaJet = new Float_t[nCases];
4868 memset(areaJet, 0, sizeof(Float_t) * nCases);
4869 Float_t* areaJetLarge = new Float_t[nCases];
4870 memset(areaJetLarge, 0, sizeof(Float_t) * nCases);
4871 Float_t areaFull = (fTrackEtaMax-fTrackEtaMin)*(fTrackPhiMax-fTrackPhiMin);
4872 Float_t areaOut = areaFull;
4874 //estimate jets and background areas
4877 TList* templist = new TList();
4878 TClonesArray *vect3Jet = new TClonesArray("TVector3",nCases);
4880 for(Int_t ij=0; ij<nCases; ++ij)
4882 // Get jet information
4883 AliAODJet* jet = dynamic_cast<AliAODJet*>(jetlist->At(ij));
4886 jet3mom.SetPtEtaPhi(jet->Pt(),jet->Eta(),jet->Phi());
4887 new((*vect3Jet)[ijet]) TVector3((TVector3)jet3mom);
4888 Float_t etaJet = (Float_t)((TVector3*) vect3Jet->At(ij))->Eta();
4891 areaJet[ij] = CalcJetArea(etaJet,rc);
4893 // Area jet larger angle
4894 areaJetLarge[ij] = CalcJetArea(etaJet,rcl);
4897 areaOut = areaOut - areaJetLarge[ij];
4901 // List of all tracks outside jet areas
4902 for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++){
4905 if( fUseExtraTracksBgr != 1){
4906 if(AliAODTrack* trackAOD = dynamic_cast<AliAODTrack*> (inputlist->At(itrack))){
4907 if(fUseExtraTracksBgr == 0 && ((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
4908 if(fUseExtraTracksBgr == -1 && !((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
4912 AliVParticle* track = dynamic_cast<AliVParticle*>(inputlist->At(itrack));
4915 Double_t trackMom[3];
4916 track->PxPyPz(trackMom);
4917 TVector3 track3mom(trackMom);
4919 Double_t *dR = new Double_t[nCases];
4920 for(Int_t ij=0; ij<nCases; ij++)
4921 dR[ij] = (Double_t)((TVector3*) vect3Jet->At(ij))->DeltaR(track3mom);
4923 if((nCases==1 && (dR[0]>rcl)) ||
4924 (nCases==2 && (dR[0]>rcl && dR[1]>rcl)) ||
4925 (nCases==3 && (dR[0]>rcl && dR[1]>rcl && dR[2]>rcl)))
4927 templist->Add(track);
4933 // Take tracks randomly
4934 Int_t nScaled = (Int_t) (nOut * areaJet[0] / areaOut + 0.5);
4935 TArrayI* ar = new TArrayI(nOut);
4937 for(Int_t init=0; init<nOut; init++)
4940 Int_t *randIndex = new Int_t[nScaled];
4941 for(Int_t init2=0; init2<nScaled; init2++)
4942 randIndex[init2] = -1;
4944 // Select nScaled different random numbers in nOut
4945 for(Int_t i=0; i<nScaled; i++)
4947 Int_t* tmpArr = new Int_t[nOut-i];
4948 Int_t temp = fRandom->Integer(nOut-i);
4949 for(Int_t ind = 0; ind< ar->GetSize()-1; ind++)
4951 if(ind<temp) tmpArr[ind] = (*ar)[ind];
4952 else tmpArr[ind] = (*ar)[ind+1];
4954 randIndex[i] = (*ar)[temp];
4956 ar->Set(nOut-i-1,tmpArr);
4962 for(Int_t ipart=0; ipart<nScaled; ipart++)
4964 AliVParticle* track = (AliVParticle*)(templist->At(randIndex[ipart]));
4965 outputlist->Add(track);
4966 sumPt += track->Pt();
4973 delete [] areaJetLarge;
4976 delete [] randIndex;
4980 // ________________________________________________________________________________________________________________________________________________________
4981 void AliAnalysisTaskIDFragmentationFunction::GetTracksOutOfNJetsStat(Int_t nCases, TList* inputlist, TList* outputlist,
4982 const TList* jetlist, Double_t& sumPt, Double_t &normFactor)
4984 // List of tracks outside cone around N jet axis
4985 // All particles taken + final scaling factor
4988 Float_t rc = TMath::Abs(GetFFRadius());
4989 Float_t rcl = GetFFBckgRadius();
4991 // Estimate jet and background areas
4992 Float_t* areaJet = new Float_t[nCases];
4993 memset(areaJet, 0, sizeof(Float_t) * nCases);
4994 Float_t* areaJetLarge = new Float_t[nCases];
4995 memset(areaJetLarge, 0, sizeof(Float_t) * nCases);
4996 Float_t areaFull = (fTrackEtaMax-fTrackEtaMin)*(fTrackPhiMax-fTrackPhiMin);
4997 Float_t areaOut = areaFull;
4999 //estimate jets and background areas
5002 TClonesArray *vect3Jet = new TClonesArray("TVector3",nCases);
5004 for(Int_t ij=0; ij<nCases; ++ij)
5006 // Get jet information
5007 AliAODJet* jet = dynamic_cast<AliAODJet*>(jetlist->At(ij));
5010 jet3mom.SetPtEtaPhi(jet->Pt(),jet->Eta(),jet->Phi());
5011 new((*vect3Jet)[ijet]) TVector3((TVector3)jet3mom);
5012 Float_t etaJet = (Float_t)((TVector3*) vect3Jet->At(ij))->Eta();
5015 areaJet[ij] = CalcJetArea(etaJet,rc);
5017 // Area jet larger angle
5018 areaJetLarge[ij] = CalcJetArea(etaJet,rcl);
5020 // Outside jets area
5021 areaOut = areaOut - areaJetLarge[ij];
5025 for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++){
5028 if( fUseExtraTracksBgr != 1){
5029 if(AliAODTrack* trackAOD = dynamic_cast<AliAODTrack*> (inputlist->At(itrack))){
5030 if(fUseExtraTracksBgr == 0 && ((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
5031 if(fUseExtraTracksBgr == -1 && !((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
5035 AliVParticle* track = dynamic_cast<AliVParticle*>(inputlist->At(itrack));
5037 Double_t trackMom[3];
5038 track->PxPyPz(trackMom);
5039 TVector3 track3mom(trackMom);
5041 Double_t *dR = new Double_t[nCases];
5042 for(Int_t ij=0; ij<nCases; ij++)
5043 dR[ij] = (Double_t)((TVector3*) vect3Jet->At(ij))->DeltaR(track3mom);
5046 (nCases==1 && (dR[0]>rcl)) ||
5047 (nCases==2 && (dR[0]>rcl && dR[1]>rcl)) ||
5048 (nCases==3 && (dR[0]>rcl && dR[1]>rcl && dR[2]>rcl)))
5050 outputlist->Add(track);
5051 sumPt += track->Pt();
5057 if(nCases==0) areaJet[0] = TMath::Pi()*rc*rc;
5058 normFactor = (Float_t) 1./(areaJet[0] / areaOut);
5063 delete [] areaJetLarge;
5068 // ______________________________________________________________________________________________________________________________________________________
5069 Float_t AliAnalysisTaskIDFragmentationFunction::CalcJetArea(Float_t etaJet, Float_t rc) const
5071 // calculate area of jet with eta etaJet and radius rc
5073 Float_t detamax = etaJet + rc;
5074 Float_t detamin = etaJet - rc;
5075 Float_t accmax = 0.0; Float_t accmin = 0.0;
5076 if(detamax > fTrackEtaMax){ // sector outside etamax
5077 Float_t h = fTrackEtaMax - etaJet;
5078 accmax = rc*rc*TMath::ACos(h/rc) - h*TMath::Sqrt(rc*rc - h*h);
5080 if(detamin < fTrackEtaMin){ // sector outside etamin
5081 Float_t h = fTrackEtaMax + etaJet;
5082 accmin = rc*rc*TMath::ACos(h/rc) - h*TMath::Sqrt(rc*rc - h*h);
5084 Float_t areaJet = rc*rc*TMath::Pi() - accmax - accmin;
5090 // ___________________________________________________________________________________________________________________________
5091 void AliAnalysisTaskIDFragmentationFunction::GetClusterTracksOutOf1Jet(AliAODJet* jet, TList* outputlist, Double_t &normFactor)
5093 // fill tracks from bckgCluster branch in list,
5094 // for all clusters outside jet cone
5095 // sum up total area of clusters
5097 Double_t rc = GetFFRadius();
5098 Double_t rcl = GetFFBckgRadius();
5100 Double_t areaTotal = 0;
5101 Double_t sumPtTotal = 0;
5103 for(Int_t ij=0; ij<fBckgJetsRec->GetEntries(); ++ij){
5105 AliAODJet* bgrCluster = (AliAODJet*)(fBckgJetsRec->At(ij)); // not 'recCuts': use all clusters in full eta range
5107 Double_t dR = jet->DeltaR(bgrCluster);
5109 if(dR<rcl) continue;
5111 Double_t clusterPt = bgrCluster->Pt();
5112 Double_t area = bgrCluster->EffectiveAreaCharged();
5114 sumPtTotal += clusterPt;
5116 Int_t nTracksJet = bgrCluster->GetRefTracks()->GetEntries();
5118 for(Int_t it = 0; it<nTracksJet; it++){
5120 // embedded tracks - note: using ref tracks here, fBranchRecBckgClusters has to be consistent
5121 if( fUseExtraTracksBgr != 1){
5122 if(AliAODTrack* trackAOD = dynamic_cast<AliAODTrack*> (bgrCluster->GetTrack(it))){
5123 if(fUseExtraTracksBgr == 0 && ((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
5124 if(fUseExtraTracksBgr == -1 && !((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
5128 AliVParticle* track = dynamic_cast<AliVParticle*>(bgrCluster->GetTrack(it));
5129 if(!track) continue;
5131 Float_t trackPt = track->Pt();
5132 Float_t trackEta = track->Eta();
5133 Float_t trackPhi = TVector2::Phi_0_2pi(track->Phi());
5135 if(trackEta < fTrackEtaMin || trackEta > fTrackEtaMax) continue;
5136 if(trackPhi < fTrackPhiMin || trackPhi > fTrackPhiMax) continue;
5137 if(trackPt < fTrackPtCut) continue;
5139 outputlist->Add(track);
5143 Double_t areaJet = TMath::Pi()*rc*rc;
5144 if(areaTotal) normFactor = (Float_t) 1./(areaJet / areaTotal);
5149 // _______________________________________________________________________________________________________________________
5150 void AliAnalysisTaskIDFragmentationFunction::GetClusterTracksMedian(TList* outputlist, Double_t &normFactor)
5152 // fill tracks from bckgCluster branch,
5153 // using cluster with median density (odd number of clusters)
5154 // or picking randomly one of the two closest to median (even number)
5158 const Int_t nBckgClusters = fBckgJetsRec->GetEntries(); // not 'recCuts': use all clusters in full eta range
5160 if(nBckgClusters<3) return; // need at least 3 clusters (skipping 2 highest)
5162 Double_t* bgrDensity = new Double_t[nBckgClusters];
5163 Int_t* indices = new Int_t[nBckgClusters];
5165 for(Int_t ij=0; ij<nBckgClusters; ++ij){
5167 AliAODJet* bgrCluster = (AliAODJet*)(fBckgJetsRec->At(ij));
5168 Double_t clusterPt = bgrCluster->Pt();
5169 Double_t area = bgrCluster->EffectiveAreaCharged();
5171 Double_t density = 0;
5172 if(area>0) density = clusterPt/area;
5174 bgrDensity[ij] = density;
5178 TMath::Sort(nBckgClusters, bgrDensity, indices);
5180 // get median cluster
5182 AliAODJet* medianCluster = 0;
5183 //Double_t medianDensity = 0;
5185 if(TMath::Odd(nBckgClusters)){
5187 Int_t medianIndex = indices[(Int_t) (0.5*(nBckgClusters-1))];
5188 medianCluster = (AliAODJet*)(fBckgJetsRec->At(medianIndex));
5190 //Double_t clusterPt = medianCluster->Pt();
5191 //Double_t area = medianCluster->EffectiveAreaCharged();
5193 //if(area>0) medianDensity = clusterPt/area;
5197 Int_t medianIndex1 = indices[(Int_t) (0.5*nBckgClusters-1)];
5198 Int_t medianIndex2 = indices[(Int_t) (0.5*nBckgClusters)];
5200 AliAODJet* medianCluster1 = (AliAODJet*)(fBckgJetsRec->At(medianIndex1));
5201 AliAODJet* medianCluster2 = (AliAODJet*)(fBckgJetsRec->At(medianIndex2));
5203 //Double_t density1 = 0;
5204 //Double_t clusterPt1 = medianCluster1->Pt();
5205 //Double_t area1 = medianCluster1->EffectiveAreaCharged();
5206 //if(area1>0) density1 = clusterPt1/area1;
5208 //Double_t density2 = 0;
5209 //Double_t clusterPt2 = medianCluster2->Pt();
5210 //Double_t area2 = medianCluster2->EffectiveAreaCharged();
5211 //if(area2>0) density2 = clusterPt2/area2;
5213 //medianDensity = 0.5*(density1+density2);
5215 medianCluster = ( (fRandom->Rndm()>0.5) ? medianCluster1 : medianCluster2 ); // select one randomly to avoid adding areas
5218 Int_t nTracksJet = medianCluster->GetRefTracks()->GetEntries();
5220 for(Int_t it = 0; it<nTracksJet; it++){
5222 // embedded tracks - note: using ref tracks here, fBranchRecBckgClusters has to be consistent
5223 if( fUseExtraTracksBgr != 1){
5224 if(AliAODTrack* trackAOD = dynamic_cast<AliAODTrack*> (medianCluster->GetTrack(it))){
5225 if(fUseExtraTracksBgr == 0 && ((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
5226 if(fUseExtraTracksBgr == -1 && !((trackAOD->GetFlags() & AliESDtrack::kEmbedded)>0)) continue;
5230 AliVParticle* track = dynamic_cast<AliVParticle*>(medianCluster->GetTrack(it));
5231 if(!track) continue;
5233 Float_t trackPt = track->Pt();
5234 Float_t trackEta = track->Eta();
5235 Float_t trackPhi = TVector2::Phi_0_2pi(track->Phi());
5237 if(trackEta < fTrackEtaMin || trackEta > fTrackEtaMax) continue;
5238 if(trackPhi < fTrackPhiMin || trackPhi > fTrackPhiMax) continue;
5239 if(trackPt < fTrackPtCut) continue;
5241 outputlist->Add(track);
5244 Double_t areaMedian = medianCluster->EffectiveAreaCharged();
5245 Double_t areaJet = TMath::Pi()*GetFFRadius()*GetFFRadius();
5247 if(areaMedian) normFactor = (Float_t) 1./(areaJet / areaMedian);
5251 delete[] bgrDensity;
5255 // ______________________________________________________________________________________________________________________________________________________
5256 void AliAnalysisTaskIDFragmentationFunction::FillBckgHistos(Int_t type, TList* inputtracklist, TList* inputjetlist, AliAODJet* jet,
5257 AliFragFuncHistos* ffbckghistocuts, AliFragFuncQATrackHistos* qabckghistocuts, TH1F* fh1Mult){
5259 // List of tracks outside jets for background study
5260 TList* tracklistout2jets = new TList();
5261 TList* tracklistout3jets = new TList();
5262 TList* tracklistout2jetsStat = new TList();
5263 TList* tracklistout3jetsStat = new TList();
5264 Double_t sumPtOut2Jets = 0.;
5265 Double_t sumPtOut3Jets = 0.;
5266 Double_t sumPtOut2JetsStat = 0.;
5267 Double_t sumPtOut3JetsStat = 0.;
5268 Double_t normFactor2Jets = 0.;
5269 Double_t normFactor3Jets = 0.;
5271 Int_t nRecJetsCuts = inputjetlist->GetEntries();
5273 if(nRecJetsCuts>1) {
5274 GetTracksOutOfNJets(2,inputtracklist, tracklistout2jets, inputjetlist, sumPtOut2Jets);
5275 GetTracksOutOfNJetsStat(2,inputtracklist, tracklistout2jetsStat, inputjetlist,sumPtOut2JetsStat, normFactor2Jets);
5278 if(nRecJetsCuts>2) {
5279 GetTracksOutOfNJets(3,inputtracklist, tracklistout3jets, inputjetlist, sumPtOut3Jets);
5280 GetTracksOutOfNJetsStat(3,inputtracklist, tracklistout3jetsStat, inputjetlist, sumPtOut3JetsStat, normFactor3Jets);
5283 if(type==kBckgOutLJ || type==kBckgOutAJ)
5285 TList* tracklistoutleading = new TList();
5286 Double_t sumPtOutLeading = 0.;
5287 GetTracksOutOfNJets(1,inputtracklist, tracklistoutleading, inputjetlist, sumPtOutLeading);
5288 if(type==kBckgOutLJ && fh1Mult) fh1Mult->Fill(tracklistoutleading->GetSize());
5290 for(Int_t it=0; it<tracklistoutleading->GetSize(); ++it){
5292 AliVParticle* trackVP = (AliVParticle*)(tracklistoutleading->At(it));
5293 if(!trackVP) continue;
5294 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5296 Float_t jetPt = jet->Pt();
5297 Float_t trackPt = trackV->Pt();
5299 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5301 if(type==kBckgOutLJ)
5303 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt);
5305 // Fill track QA for background
5306 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt);
5309 // All cases included
5310 if(nRecJetsCuts==1 && type==kBckgOutAJ)
5312 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt );
5316 // Increment jet pt with one entry in case #tracks outside jets = 0
5317 if(tracklistoutleading->GetSize()==0) {
5318 Float_t jetPt = jet->Pt();
5319 Bool_t incrementJetPt = kTRUE;
5320 if(type==kBckgOutLJ)
5322 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
5324 // All cases included
5325 if(nRecJetsCuts==1 && type==kBckgOutAJ)
5327 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
5330 delete tracklistoutleading;
5332 if(type==kBckgOutLJStat || type==kBckgOutAJStat)
5334 TList* tracklistoutleadingStat = new TList();
5335 Double_t sumPtOutLeadingStat = 0.;
5336 Double_t normFactorLeading = 0.;
5338 GetTracksOutOfNJetsStat(1,inputtracklist, tracklistoutleadingStat, inputjetlist, sumPtOutLeadingStat, normFactorLeading);
5339 if(type==kBckgOutLJStat && fh1Mult) fh1Mult->Fill(tracklistoutleadingStat->GetSize());
5341 for(Int_t it=0; it<tracklistoutleadingStat->GetSize(); ++it){
5343 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistoutleadingStat->At(it));
5344 if(!trackVP) continue;
5345 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5347 Float_t jetPt = jet->Pt();
5348 Float_t trackPt = trackV->Pt();
5349 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5352 if(type==kBckgOutLJStat)
5354 if(fFFMode)ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactorLeading);
5356 // Fill track QA for background
5357 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt); // OB added bgr QA
5360 // All cases included
5361 if(nRecJetsCuts==1 && type==kBckgOutAJStat)
5363 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactorLeading);
5364 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt ); // OB added bgr QA
5369 // Increment jet pt with one entry in case #tracks outside jets = 0
5370 if(tracklistoutleadingStat->GetSize()==0) {
5371 Float_t jetPt = jet->Pt();
5372 Bool_t incrementJetPt = kTRUE;
5373 if(type==kBckgOutLJStat)
5375 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactorLeading);
5377 // All cases included
5378 if(nRecJetsCuts==1 && type==kBckgOutLJStat)
5380 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactorLeading);
5384 delete tracklistoutleadingStat;
5387 if(type==kBckgPerp || type==kBckgPerp2 || type==kBckgPerp2Area)
5389 Double_t sumPtPerp1 = 0.;
5390 Double_t sumPtPerp2 = 0.;
5391 TList* tracklistperp = new TList();
5392 TList* tracklistperp1 = new TList();
5393 TList* tracklistperp2 = new TList();
5396 if(type == kBckgPerp2) norm = 2; // in FillFF() scaleFac = 1/norm = 0.5 - account for double area;
5397 if(type == kBckgPerp2Area) norm = 2*TMath::Pi()*TMath::Abs(GetFFRadius())*TMath::Abs(GetFFRadius()) / jet->EffectiveAreaCharged(); // in FillFF() scaleFac = 1/norm;
5399 GetTracksTiltedwrpJetAxis(TMath::Pi()/2., inputtracklist,tracklistperp1,jet,TMath::Abs(GetFFRadius()),sumPtPerp1);
5400 if(type==kBckgPerp2 || type==kBckgPerp2Area) GetTracksTiltedwrpJetAxis(-1*TMath::Pi()/2., inputtracklist,tracklistperp2,jet,TMath::Abs(GetFFRadius()),sumPtPerp2);
5402 tracklistperp->AddAll(tracklistperp1);
5403 tracklistperp->AddAll(tracklistperp2);
5405 if(tracklistperp->GetSize() != tracklistperp1->GetSize() + tracklistperp2->GetSize()){
5406 cout<<" ERROR: tracklistperp size "<<tracklistperp->GetSize()<<" perp1 "<<tracklistperp1->GetSize()<<" perp2 "<<tracklistperp2->GetSize()<<endl;
5410 if(fh1Mult) fh1Mult->Fill(tracklistperp->GetSize());
5412 for(Int_t it=0; it<tracklistperp->GetSize(); ++it){
5414 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistperp->At(it));
5415 if(!trackVP)continue;
5416 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5418 Float_t jetPt = jet->Pt();
5419 Float_t trackPt = trackV->Pt();
5421 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5423 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, norm );
5425 // Fill track QA for background
5426 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt);
5430 // Increment jet pt with one entry in case #tracks outside jets = 0
5431 if(tracklistperp->GetSize()==0) {
5432 Float_t jetPt = jet->Pt();
5433 Bool_t incrementJetPt = kTRUE;
5434 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
5439 // fill for tracklistperp1/2 separately, divide norm by 2
5440 if(type==kBckgPerp){
5441 FillJetShape(jet, tracklistperp, fProNtracksLeadingJetBgrPerp2, fProDelRPtSumBgrPerp2, 0, TMath::Pi()/2., 0., kFALSE);
5443 if(type==kBckgPerp2){
5444 FillJetShape(jet, tracklistperp1, fProNtracksLeadingJetBgrPerp2, fProDelRPtSumBgrPerp2, 0, TMath::Pi()/2., 0., kFALSE);
5445 FillJetShape(jet, tracklistperp2, fProNtracksLeadingJetBgrPerp2, fProDelRPtSumBgrPerp2, 0, -1*TMath::Pi()/2., 0., kFALSE);
5447 if(type==kBckgPerp2Area){ // divide norm by 2: listperp1/2 filled separately
5448 FillJetShape(jet, tracklistperp1, fProNtracksLeadingJetBgrPerp2, fProDelRPtSumBgrPerp2, 0, TMath::Pi()/2., 0.5*norm, kFALSE);
5449 FillJetShape(jet, tracklistperp2, fProNtracksLeadingJetBgrPerp2, fProDelRPtSumBgrPerp2, 0, -1*TMath::Pi()/2., 0.5*norm, kFALSE);
5453 delete tracklistperp;
5454 delete tracklistperp1;
5455 delete tracklistperp2;
5458 if(type==kBckgASide)
5460 Double_t sumPtASide = 0.;
5461 TList* tracklistaside = new TList();
5462 GetTracksTiltedwrpJetAxis(TMath::Pi(),inputtracklist,tracklistaside,jet,TMath::Abs(GetFFRadius()),sumPtASide);
5463 if(fh1Mult) fh1Mult->Fill(tracklistaside->GetSize());
5465 for(Int_t it=0; it<tracklistaside->GetSize(); ++it){
5467 AliVParticle* trackVP = (AliVParticle*)(tracklistaside->At(it));
5468 if(!trackVP) continue;
5469 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5471 Float_t jetPt = jet->Pt();
5472 Float_t trackPt = trackV->Pt();
5474 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5476 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt );
5478 // Fill track QA for background
5479 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt);
5483 if(tracklistaside->GetSize()==0) {
5484 Float_t jetPt = jet->Pt();
5485 Bool_t incrementJetPt = kTRUE;
5486 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
5489 delete tracklistaside;
5492 if(type==kBckgASideWindow)
5494 Double_t normFactorASide = 0.;
5495 Double_t sumPtASideW = 0.;
5496 TList* tracklistasidew = new TList();
5497 GetTracksTiltedwrpJetAxisWindow(TMath::Pi(),inputtracklist,tracklistasidew,jet,TMath::Abs(GetFFRadius()),sumPtASideW,normFactorASide);
5498 if(fh1Mult) fh1Mult->Fill(tracklistasidew->GetSize());
5500 for(Int_t it=0; it<tracklistasidew->GetSize(); ++it){
5502 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistasidew->At(it));
5503 if(!trackVP) continue;
5504 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5506 Float_t jetPt = jet->Pt();
5507 Float_t trackPt = trackV->Pt();
5508 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5510 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactorASide);
5512 // Fill track QA for background
5513 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt, kFALSE, normFactorASide);
5517 if(tracklistasidew->GetSize()==0) {
5518 Float_t jetPt = jet->Pt();
5519 Bool_t incrementJetPt = kTRUE;
5520 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactorASide);
5523 delete tracklistasidew;
5526 if(type==kBckgPerpWindow)
5528 Double_t normFactorPerp = 0.;
5529 Double_t sumPtPerpW = 0.;
5530 TList* tracklistperpw = new TList();
5531 GetTracksTiltedwrpJetAxisWindow(TMath::Pi()/2.,inputtracklist,tracklistperpw,jet,TMath::Abs(GetFFRadius()),sumPtPerpW,normFactorPerp);
5532 if(fh1Mult) fh1Mult->Fill(tracklistperpw->GetSize());
5534 for(Int_t it=0; it<tracklistperpw->GetSize(); ++it){
5536 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistperpw->At(it));
5537 if(!trackVP) continue;
5538 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5540 Float_t jetPt = jet->Pt();
5541 Float_t trackPt = trackV->Pt();
5542 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5544 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactorPerp);
5546 // Fill track QA for background
5547 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt, kFALSE, normFactorPerp);
5551 if(tracklistperpw->GetSize()==0) {
5552 Float_t jetPt = jet->Pt();
5553 Bool_t incrementJetPt = kTRUE;
5554 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactorPerp);
5557 delete tracklistperpw;
5561 if(type==kBckgOut2J || type==kBckgOutAJ)
5563 if(type==kBckgOut2J && fh1Mult) fh1Mult->Fill(tracklistout2jets->GetSize());
5564 for(Int_t it=0; it<tracklistout2jets->GetSize(); ++it){
5566 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistout2jets->At(it));
5567 if(!trackVP) continue;
5568 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5570 Float_t jetPt = jet->Pt();
5571 Float_t trackPt = trackV->Pt();
5573 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5575 if(type==kBckgOut2J)
5577 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt );
5578 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt);
5581 // All cases included
5582 if(nRecJetsCuts==2 && type==kBckgOutAJ)
5584 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt );
5589 // Increment jet pt with one entry in case #tracks outside jets = 0
5590 if(tracklistout2jets->GetSize()==0) {
5591 Float_t jetPt = jet->Pt();
5592 Bool_t incrementJetPt = kTRUE;
5593 if(type==kBckgOut2J)
5595 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
5597 // All cases included
5598 if(nRecJetsCuts==2 && type==kBckgOutAJ)
5600 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
5605 if(type==kBckgOut2JStat || type==kBckgOutAJStat)
5607 for(Int_t it=0; it<tracklistout2jetsStat->GetSize(); ++it){
5609 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistout2jetsStat->At(it));
5610 if(!trackVP) continue;
5611 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5613 Float_t jetPt = jet->Pt();
5614 Float_t trackPt = trackV->Pt();
5615 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5617 if(type==kBckgOut2JStat)
5619 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactor2Jets);
5621 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt ); // OB added bgr QA
5624 // All cases included
5625 if(nRecJetsCuts==2 && type==kBckgOutAJStat)
5627 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactor2Jets);
5629 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt ); // OB added bgr QA
5633 // Increment jet pt with one entry in case #tracks outside jets = 0
5634 if(tracklistout2jetsStat->GetSize()==0) {
5635 Float_t jetPt = jet->Pt();
5636 Bool_t incrementJetPt = kTRUE;
5637 if(type==kBckgOut2JStat)
5639 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactor2Jets);
5641 // All cases included
5642 if(nRecJetsCuts==2 && type==kBckgOutAJStat)
5644 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactor2Jets);
5650 if(type==kBckgOut3J || type==kBckgOutAJ)
5652 if(type==kBckgOut3J && fh1Mult) fh1Mult->Fill(tracklistout3jets->GetSize());
5654 for(Int_t it=0; it<tracklistout3jets->GetSize(); ++it){
5656 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistout3jets->At(it));
5657 if(!trackVP) continue;
5658 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5660 Float_t jetPt = jet->Pt();
5661 Float_t trackPt = trackV->Pt();
5663 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5665 if(type==kBckgOut3J)
5667 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt );
5669 qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt);
5672 // All cases included
5673 if(nRecJetsCuts==3 && type==kBckgOutAJ)
5675 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt );
5680 // Increment jet pt with one entry in case #tracks outside jets = 0
5681 if(tracklistout3jets->GetSize()==0) {
5682 Float_t jetPt = jet->Pt();
5683 Bool_t incrementJetPt = kTRUE;
5684 if(type==kBckgOut3J)
5686 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
5688 // All cases included
5689 if(nRecJetsCuts==3 && type==kBckgOutAJ)
5691 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt );
5696 if(type==kBckgOut3JStat || type==kBckgOutAJStat)
5698 for(Int_t it=0; it<tracklistout3jetsStat->GetSize(); ++it){
5700 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistout3jetsStat->At(it));
5701 if(!trackVP) continue;
5702 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5704 Float_t jetPt = jet->Pt();
5705 Float_t trackPt = trackV->Pt();
5706 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5708 if(type==kBckgOut3JStat)
5710 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactor3Jets);
5712 //if(fQAMode&1) qabckghistocuts->FillTrackQA( trackEta, TVector2::Phi_0_2pi(trackPhi), trackPt);
5715 // All cases included
5716 if(nRecJetsCuts==3 && type==kBckgOutAJStat)
5718 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactor3Jets );
5720 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt );
5725 // Increment jet pt with one entry in case #tracks outside jets = 0
5726 if(tracklistout3jetsStat->GetSize()==0) {
5727 Float_t jetPt = jet->Pt();
5728 Bool_t incrementJetPt = kTRUE;
5729 if(type==kBckgOut3JStat)
5731 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactor3Jets);
5733 // All cases included
5734 if(nRecJetsCuts==3 && type==kBckgOutAJStat)
5736 if(fFFMode) ffbckghistocuts->FillFF( -1., jetPt, incrementJetPt, normFactor3Jets);
5742 if(type==kBckgClustersOutLeading){ // clusters bgr: all tracks in clusters out of leading jet
5744 TList* tracklistClustersOutLeading = new TList();
5745 Double_t normFactorClusters = 0;
5746 Float_t jetPt = jet->Pt();
5748 GetClusterTracksOutOf1Jet(jet, tracklistClustersOutLeading, normFactorClusters);
5749 if(fh1Mult) fh1Mult->Fill(tracklistClustersOutLeading->GetSize());
5751 for(Int_t it=0; it<tracklistClustersOutLeading->GetSize(); ++it){
5753 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistClustersOutLeading->At(it));
5754 if(!trackVP) continue;
5755 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5757 Float_t trackPt = trackVP->Pt();
5759 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5761 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactorClusters );
5762 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt );
5767 delete tracklistClustersOutLeading;
5771 if(type == kBckgClusters){ // clusters bgr: all tracks in 'median cluster'
5773 TList* tracklistClustersMedian = new TList();
5774 Double_t normFactorClusters = 0;
5775 Float_t jetPt = jet->Pt();
5777 GetClusterTracksMedian(tracklistClustersMedian, normFactorClusters);
5778 if(fh1Mult) fh1Mult->Fill(tracklistClustersMedian->GetSize());
5780 for(Int_t it=0; it<tracklistClustersMedian->GetSize(); ++it){
5782 AliVParticle* trackVP = dynamic_cast<AliVParticle*>(tracklistClustersMedian->At(it));
5783 if(!trackVP) continue;
5784 TLorentzVector* trackV = new TLorentzVector(trackVP->Px(),trackVP->Py(),trackVP->Pz(),trackVP->P());
5786 Float_t trackPt = trackVP->Pt();
5788 Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
5790 if(fFFMode) ffbckghistocuts->FillFF( trackPt, jetPt, incrementJetPt, normFactorClusters );
5791 if(fQAMode&1) qabckghistocuts->FillTrackQA( trackV->Eta(), TVector2::Phi_0_2pi(trackV->Phi()), trackPt );
5796 delete tracklistClustersMedian;
5799 delete tracklistout2jets;
5800 delete tracklistout3jets;
5801 delete tracklistout2jetsStat;
5802 delete tracklistout3jetsStat;
5805 //_____________________________________________________________________________________
5806 Double_t AliAnalysisTaskIDFragmentationFunction::GetMCStrangenessFactor(Double_t pt) const
5808 // factor strangeness data/MC as function of pt from UE analysis (Sara Vallero)
5812 if(0.150<pt && pt<0.200) alpha = 3.639;
5813 if(0.200<pt && pt<0.250) alpha = 2.097;
5814 if(0.250<pt && pt<0.300) alpha = 1.930;
5815 if(0.300<pt && pt<0.350) alpha = 1.932;
5816 if(0.350<pt && pt<0.400) alpha = 1.943;
5817 if(0.400<pt && pt<0.450) alpha = 1.993;
5818 if(0.450<pt && pt<0.500) alpha = 1.989;
5819 if(0.500<pt && pt<0.600) alpha = 1.963;
5820 if(0.600<pt && pt<0.700) alpha = 1.917;
5821 if(0.700<pt && pt<0.800) alpha = 1.861;
5822 if(0.800<pt && pt<0.900) alpha = 1.820;
5823 if(0.900<pt && pt<1.000) alpha = 1.741;
5824 if(1.000<pt && pt<1.500) alpha = 0.878;
5829 //__________________________________________________________________________________________________
5830 Double_t AliAnalysisTaskIDFragmentationFunction::GetMCStrangenessFactorCMS(AliAODMCParticle* daughter) const
5832 // strangeness ratio MC/data as function of mother pt from CMS data in |eta|<2.0
5834 TClonesArray *tca = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
5837 AliAODMCParticle* currentMother = daughter;
5838 AliAODMCParticle* currentDaughter = daughter;
5841 // find first primary mother K0s, Lambda or Xi
5844 Int_t daughterPDG = currentDaughter->GetPdgCode();
5846 Int_t motherLabel = currentDaughter->GetMother();
5847 if(motherLabel >= tca->GetEntriesFast()){ // protection
5848 currentMother = currentDaughter;
5852 currentMother = (AliAODMCParticle*) tca->At(motherLabel);
5855 currentMother = currentDaughter;
5859 Int_t motherPDG = currentMother->GetPdgCode();
5861 // phys. primary found ?
5862 if(currentMother->IsPhysicalPrimary()) break;
5864 if(TMath::Abs(daughterPDG) == 321){ // K+/K- e.g. from phi (ref data not feeddown corrected)
5865 currentMother = currentDaughter; break;
5867 if(TMath::Abs(motherPDG) == 310 ){ // K0s e.g. from phi (ref data not feeddown corrected)
5870 if(TMath::Abs(motherPDG) == 3212 && TMath::Abs(daughterPDG) == 3122){ // mother Sigma0, daughter Lambda (this case not included in feeddown corr.)
5871 currentMother = currentDaughter; break;
5874 currentDaughter = currentMother;
5878 Int_t motherPDG = currentMother->GetPdgCode();
5879 Double_t motherGenPt = currentMother->Pt();
5881 return AliAnalysisTaskPID::GetMCStrangenessFactorCMS(motherPDG, motherGenPt);
5884 // _________________________________________________________________________________
5885 void AliAnalysisTaskIDFragmentationFunction::FillJetShape(AliAODJet* jet, TList* list,
5886 TProfile* hProNtracksLeadingJet, TProfile** hProDelRPtSum, TProfile* hProDelR80pcPt,
5887 Double_t dPhiUE, Double_t normUE, Bool_t scaleStrangeness)
5889 // Fill jet shape histos
5891 const Int_t kNbinsR = 50;
5892 const Float_t kBinWidthR = 0.02;
5894 Int_t nJetTracks = list->GetEntries();
5896 Float_t pTsumA[kNbinsR] = {0.0};
5898 Float_t *delRA = new Float_t[nJetTracks];
5899 Float_t *trackPtA = new Float_t[nJetTracks];
5900 Int_t *index = new Int_t[nJetTracks];
5902 for(Int_t i=0; i<nJetTracks; i++){
5909 jet->PxPyPz(jetMom);
5910 TVector3 jet3mom(jetMom);
5912 if(TMath::Abs(dPhiUE)>0){
5913 Double_t phiTilted = jet3mom.Phi();
5914 phiTilted += dPhiUE;
5915 phiTilted = TVector2::Phi_0_2pi(phiTilted);
5916 jet3mom.SetPhi(phiTilted);
5919 Double_t jetPt = jet->Pt();
5920 Double_t sumWeights = 0;
5922 for (Int_t j =0; j<nJetTracks; j++){
5924 AliVParticle* track = dynamic_cast<AliVParticle*>(list->At(j));
5927 Double_t trackMom[3];
5928 track->PxPyPz(trackMom);
5929 TVector3 track3mom(trackMom);
5931 Double_t dR = jet3mom.DeltaR(track3mom);
5934 trackPtA[j] = track->Pt();
5936 Double_t weight = GetMCStrangenessFactor(track->Pt()); // more correctly should be gen pt
5937 sumWeights += weight;
5939 for(Int_t ibin=1; ibin<=kNbinsR; ibin++){
5940 Float_t xlow = kBinWidthR*(ibin-1);
5941 Float_t xup = kBinWidthR*ibin;
5942 if(xlow <= dR && dR < xup){
5944 if(scaleStrangeness) pTsumA[ibin-1] += track->Pt()*weight;
5945 else pTsumA[ibin-1] += track->Pt();
5953 for(Int_t ibin=0; ibin<kNbinsR; ibin++){
5954 Float_t fR = kBinWidthR*(ibin+0.5);
5956 for(Int_t k=0; k<5; k++){
5957 if(k==0){jetPtMin=20.0;jetPtMax=30.0;}
5958 if(k==1){jetPtMin=30.0;jetPtMax=40.0;}
5959 if(k==2){jetPtMin=40.0;jetPtMax=60.0;}
5960 if(k==3){jetPtMin=60.0;jetPtMax=80.0;}
5961 if(k==4){jetPtMin=80.0;jetPtMax=100.0;}
5962 if(jetPt>jetPtMin && jetPt<jetPtMax){
5964 hProDelRPtSum[k]->Fill(fR,pTsumA[ibin]);
5970 if(scaleStrangeness) hProNtracksLeadingJet->Fill(jetPt,sumWeights);
5971 else hProNtracksLeadingJet->Fill(jetPt,nJetTracks);
5973 if(normUE) hProNtracksLeadingJet->Fill(jetPt,nJetTracks/normUE);
5978 Float_t delRPtSum80pc = 0;
5980 TMath::Sort(nJetTracks,delRA,index,0);
5982 for(Int_t ii=0; ii<nJetTracks; ii++){
5984 if(scaleStrangeness){
5985 Double_t weight = GetMCStrangenessFactor(trackPtA[index[ii]]); // more correctly should be gen pt
5986 pTsum += weight*trackPtA[index[ii]];
5988 else pTsum += trackPtA[index[ii]];
5991 if(pTsum/jetPt >= 0.8000){
5992 delRPtSum80pc = delRA[index[ii]];
5996 hProDelR80pcPt->Fill(jetPt,delRPtSum80pc);
6005 // _________________________________________________________________________________
6006 Bool_t AliAnalysisTaskIDFragmentationFunction::IsSecondaryWithStrangeMotherMC(AliAODMCParticle* part)
6008 // Check whether particle is a secondary with strange mother, i.e. returns kTRUE if a strange mother is found
6009 // and the particle is NOT a physical primary. In all other cases kFALSE is returned
6011 TClonesArray *tca = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
6015 if (part->IsPhysicalPrimary())
6018 Int_t iMother = part->GetMother();
6023 AliAODMCParticle* partM = dynamic_cast<AliAODMCParticle*>(tca->At(iMother));
6027 Int_t codeM = TMath::Abs(partM->GetPdgCode());
6028 Int_t mfl = Int_t(codeM / TMath::Power(10, Int_t(TMath::Log10(codeM))));
6029 if (mfl == 3 && codeM != 3) // codeM = 3 is for s quark