28 #include <AliAnalysisTask.h>
29 #include <AliAnalysisManager.h>
30 #include <AliCentrality.h>
31 #include <AliVVertex.h>
32 #include <AliESDEvent.h>
33 #include <AliAODEvent.h>
34 #include <AliAODTrack.h>
41 #include "AliMultSelection.h"
64 for(
Int_t i(0); i < 10; i++) {
85 for(
Int_t i(0); i < 10; i++) {
90 DefineInput(0, TChain::Class());
91 DefineOutput(1, TList::Class());
95 DefineOutput(2, TList::Class());
96 DefineOutput(3, TList::Class());
116 #ifdef ALIANALYSISTASKLOCALRHO_DEBUG_FLAG_0
117 printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
124 if(!(InputEvent()->FindListObject(
fLocalRho->GetName()))) {
127 AliFatal(Form(
"%s: Container with same name %s already present. Aborting", GetName(),
fLocalRho->GetName()));
135 AliFatal(Form(
"%s: Couldn't find container for scaled rho. Aborting !", GetName()));
138 if(!
GetJetContainer()) AliFatal(Form(
"%s: Couldn't get jet container. Aborting !", GetName()));
145 #ifdef ALIANALYSISTASKLOCALRHO_DEBUG_FLAG_0
146 printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
153 SetModulationFit(
new TF1(
"fit_kV2",
"[0]*([1]+[2]*[3]*TMath::Cos([2]*(x-[4])))", 0, TMath::TwoPi()));
160 SetModulationFit(
new TF1(
"fit_kV3",
"[0]*([1]+[2]*[3]*TMath::Cos([2]*(x-[4])))", 0, TMath::TwoPi()));
167 SetModulationFit(
new TF1(
"fit_kCombined",
"[0]*([1]+[2]*([3]*TMath::Cos([2]*(x-[4]))+[7]*TMath::Cos([5]*(x-[6]))))", 0, TMath::TwoPi()));
188 #ifdef ALIANALYSISTASKLOCALRHO_DEBUG_FLAG_0
189 printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
191 fHistSwap =
new TH1F(
"fHistSwap",
"fHistSwap", 20, 0, TMath::TwoPi());
193 Int_t c[] = {0, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100};
205 fHistPsi2[i] =
BookTH1F(
"fHistPsi2",
"#Psi_{2}", 100, -.5*TMath::Pi(), .5*TMath::Pi(), i);
206 fHistPsi3[i] =
BookTH1F(
"fHistPsi3",
"#Psi_{3}", 100, -1.*TMath::Pi()/3., TMath::Pi()/3., i);
258 #ifdef ALIANALYSISTASKLOCALRHO_DEBUG_FLAG_0
259 printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
264 name = Form(
"%s_%i", name, c);
267 title += Form(
";%s;[counts]", x);
268 TH1F* histogram =
new TH1F(name, title.Data(), bins, min, max);
279 #ifdef ALIANALYSISTASKLOCALRHO_DEBUG_FLAG_0
280 printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
285 name = Form(
"%s_%i", name, c);
288 title += Form(
";%s;%s", x, y);
289 TH2F* histogram =
new TH2F(name, title.Data(), binsx, minx, maxx, binsy, miny, maxy);
299 #ifdef ALIANALYSISTASKLOCALRHO_DEBUG_FLAG_0
300 printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
308 AliMultSelection *MultSelection = 0x0;
309 MultSelection = (AliMultSelection*)InputEvent()->FindListObject(
"MultSelection");
312 AliWarning(
"AliMultSelection object not found!");
314 lPercentile = MultSelection->GetMultiplicityPercentile(
"V0M");
335 psi2 = tpc[0]; psi3 = tpc[1];
342 psi2 = vzero[0][0]; psi3 = vzero[0][1];
349 psi2 = vzero[1][0]; psi3 = vzero[1][1];
357 psi2 = vzeroComb[0]; psi3 = vzeroComb[1];
452 #ifdef ALIANALYSISTASKLOCALRHO_DEBUG_FLAG_0
453 printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
462 Double_t a(0), b(0),
c(0), d(0), e(0), f(0), g(0), h(0);
463 vzero[0][0] = InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 8, 2, a, b);
464 vzero[1][0] = InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 9, 2,
c, d);
465 vzero[0][1] = InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 8, 3, e, f);
466 vzero[1][1] = InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 9, 3, g, h);
470 if(fDebug > 0) printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
471 Double_t qxa2(0), qya2(0), qxc2(0), qyc2(0);
472 Double_t qxa3(0), qya3(0), qxc3(0), qyc3(0);
473 for(
Int_t iVZERO(0); iVZERO < 64; iVZERO++) {
474 Double_t phi(TMath::PiOver4()*(.5+iVZERO%8)), weight(InputEvent()->GetVZEROEqMultiplicity(iVZERO));
477 qxa2 += weight*TMath::Cos(2.*phi);
478 qya2 += weight*TMath::Sin(2.*phi);
479 qxa3 += weight*TMath::Cos(3.*phi);
480 qya3 += weight*TMath::Sin(3.*phi);
483 qxc2 += weight*TMath::Cos(2.*phi);
484 qyc2 += weight*TMath::Sin(2.*phi);
485 qxc3 += weight*TMath::Cos(3.*phi);
486 qyc3 += weight*TMath::Sin(3.*phi);
489 vzero[0][0] = .5*TMath::ATan2(qya2, qxa2);
490 vzero[1][0] = .5*TMath::ATan2(qyc2, qxc2);
491 vzero[0][1] = (1./3.)*TMath::ATan2(qya3, qxa3);
492 vzero[1][1] = (1./3.)*TMath::ATan2(qyc3, qxc3);
501 #ifdef ALIANALYSISTASKLOCALRHO_DEBUG_FLAG_0
502 printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
511 if(leadingJet) leadingJet->
Eta();
514 for(
Int_t iTPC(0); iTPC < iTracks; iTPC++) {
515 AliVTrack* track =
static_cast<AliVTrack*
>(
fTracks->At(iTPC));
519 qx2+= TMath::Cos(2.*track->Phi());
520 qy2+= TMath::Sin(2.*track->Phi());
521 qx3+= TMath::Cos(3.*track->Phi());
522 qy3+= TMath::Sin(3.*track->Phi());
525 tpc[0] = .5*TMath::ATan2(qy2, qx2);
526 tpc[1] = (1./3.)*TMath::ATan2(qy3, qx3);
533 #ifdef ALIANALYSISTASKLOCALRHO_DEBUG_FLAG_0
534 printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
537 comb[0] = InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 10, 2, a, b);
538 comb[1] = InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 10, 3,
c, d);
545 #ifdef ALIANALYSISTASKLOCALRHO_DEBUG_FLAG_0
546 printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
548 Double_t reQ(0), imQ(0), modQ(0), M11(0), M(0);
550 QCnQnk(harm, 1, reQ, imQ);
551 modQ = reQ*reQ+imQ*imQ;
553 return (M11 > 0) ? ((modQ -
QCnS(1,2))/M11) : -999;
555 QCnQnk(harm, 0, reQ, imQ);
556 modQ = reQ*reQ+imQ*imQ;
558 return (M > 1) ? (modQ - M)/(M*(M-1)) : -999;
565 #ifdef ALIANALYSISTASKLOCALRHO_DEBUG_FLAG_0
566 printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
568 Double_t reQn1(0), imQn1(0), reQ2n2(0), imQ2n2(0), reQn3(0), imQn3(0), M1111(0), M(0);
569 Double_t a(0), b(0),
c(0), d(0), e(0), f(0), g(0);
571 QCnQnk(harm, 1, reQn1, imQn1);
572 QCnQnk(harm*2, 2, reQ2n2, imQ2n2);
573 QCnQnk(harm, 3, reQn3, imQn3);
575 a = (reQn1*reQn1+imQn1*imQn1)*(reQn1*reQn1+imQn1*imQn1);
576 b = reQ2n2*reQ2n2 + imQ2n2*imQ2n2;
577 c = -2.*(reQ2n2*reQn1*reQn1-reQ2n2*imQn1*imQn1+2.*imQ2n2*reQn1*imQn1);
578 d = 8.*(reQn3*reQn1+imQn3*imQn1);
579 e = -4.*
QCnS(1,2)*(reQn1*reQn1+imQn1*imQn1);
583 return (M1111 > 0) ? (a+b+
c+d+e+f+g)/M1111 : -999;
585 Double_t reQn(0), imQn(0), reQ2n(0), imQ2n(0);
586 QCnQnk(harm, 0, reQn, imQn);
587 QCnQnk(harm*2, 0, reQ2n, imQ2n);
590 if(M < 4)
return -999;
591 a = (reQn*reQn+imQn*imQn)*(reQn*reQn+imQn*imQn);
592 b = reQ2n*reQ2n + imQ2n*imQ2n;
593 c = -2.*(reQ2n*reQn*reQn-reQ2n*imQn*imQn+2.*imQ2n*reQn*imQn);
594 e = -4.*(M-2)*(reQn*reQn+imQn*imQn);
596 return (a+b+
c+e+f)/(M*(M-1)*(M-2)*(M-3));
603 #ifdef ALIANALYSISTASKLOCALRHO_DEBUG_FLAG_0
604 printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
609 for(
Int_t iTPC(0); iTPC < iTracks; iTPC++) {
610 AliVTrack* track =
static_cast<AliVTrack*
>(
fTracks->At(iTPC));
614 reQ += TMath::Power(track->Pt(), k) * TMath::Cos(((
double)n)*track->Phi());
615 imQ += TMath::Power(track->Pt(), k) * TMath::Sin(((
double)n)*track->Phi());
623 #ifdef ALIANALYSISTASKLOCALRHO_DEBUG_FLAG_0
624 printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
626 if(!
fTracks || i <= 0 || j <= 0)
return -999;
629 for(
Int_t iTPC(0); iTPC < iTracks; iTPC++) {
630 AliVTrack* track =
static_cast<AliVTrack*
>(
fTracks->At(iTPC));
632 Sij+=TMath::Power(track->Pt(), j);
634 return TMath::Power(Sij, i);
641 #ifdef ALIANALYSISTASKLOCALRHO_DEBUG_FLAG_0
642 printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
651 #ifdef ALIANALYSISTASKLOCALRHO_DEBUG_FLAG_0
652 printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
661 #ifdef ALIANALYSISTASKLOCALRHO_DEBUG_FLAG_0
662 printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
672 #ifdef ALIANALYSISTASKLOCALRHO_DEBUG_FLAG_0
673 printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
691 if( c2 < 0 ) c2 = -1.*TMath::Sqrt(-1.*c2);
692 if( c3 < 0 ) c3 = -1.*TMath::Sqrt(-1.*c3);
706 default :
return kFALSE;
729 #ifdef ALIANALYSISTASKLOCALRHO_DEBUG_FLAG_0
730 printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
756 }
else if (!
QCnRecovery(psi2, psi3))
return kFALSE;
786 }
else if (!
QCnRecovery(psi2, psi3))
return kFALSE;
812 case kTPC : detector+=
"TPC";
814 case kVZEROA : detector+=
"VZEROA";
816 case kVZEROC : detector+=
"VZEROC";
826 if(iTracks <= 0 || fLocalRho->GetVal() <= 0 )
return kFALSE;
830 excludeInEta = leadingJet->
Eta();
831 excludeInPhi = leadingJet->
Phi();
832 excludeInPt = leadingJet->
Pt();
840 _tempSwap = TH1F(
"_tempSwap",
"_tempSwap", TMath::CeilNint(TMath::Sqrt(
fNAcceptedTracks)), 0, TMath::TwoPi());
846 Double_t totalpts(0.), totalptsquares(0.), totalns(0.);
847 for(
Int_t i(0); i < iTracks; i++) {
848 AliVTrack* track =
static_cast<AliVTrack*
>(
fTracks->At(i));
852 _tempSwap.Fill(track->Phi(), track->Pt());
854 totalpts += track->Pt();
855 totalptsquares += track->Pt()*track->Pt();
857 _tempSwapN.Fill(track->Phi());
860 else _tempSwap.Fill(track->Phi());
867 if(totalns < 1)
return kFALSE;
868 for(
Int_t l = 0; l < _tempSwap.GetNbinsX(); l++) {
869 if(_tempSwapN.GetBinContent(l+1) == 0) {
870 _tempSwap.SetBinContent(l+1,0);
871 _tempSwap.SetBinError(l+1,0);
874 Double_t vartimesnsq = totalptsquares*totalns - totalpts*totalpts;
875 Double_t variance = vartimesnsq/(totalns*(totalns-1.));
876 Double_t SDOMSq = variance / _tempSwapN.GetBinContent(l+1);
877 Double_t SDOMSqOverMeanSq = SDOMSq * _tempSwapN.GetBinContent(l+1) * _tempSwapN.GetBinContent(l+1) / (_tempSwapN.GetBinContent(l+1) * _tempSwapN.GetBinContent(l+1));
878 Double_t poissonfrac = 1./_tempSwapN.GetBinContent(l+1);
879 Double_t vartotalfrac = SDOMSqOverMeanSq + poissonfrac;
880 Double_t vartotal = vartotalfrac * _tempSwap.GetBinContent(l+1) * _tempSwap.GetBinContent(l+1);
881 if(vartotal > 0.0001) _tempSwap.SetBinError(l+1,TMath::Sqrt(vartotal));
883 _tempSwap.SetBinContent(l+1,0);
884 _tempSwap.SetBinError(l+1,0);
907 Double_t cos2(0), sin2(0), cos3(0), sin3(0), sumPt(0);
908 for(
Int_t i(0); i < iTracks; i++) {
909 AliVTrack* track =
static_cast<AliVTrack*
>(
fTracks->At(i));
911 sumPt += track->Pt();
912 cos2 += track->Pt()*TMath::Cos(2*
PhaseShift(track->Phi()-psi2));
913 sin2 += track->Pt()*TMath::Sin(2*
PhaseShift(track->Phi()-psi2));
914 cos3 += track->Pt()*TMath::Cos(3*
PhaseShift(track->Phi()-psi3));
915 sin3 += track->Pt()*TMath::Sin(3*
PhaseShift(track->Phi()-psi3));
936 static Int_t didacticCounterBest(0);
937 TProfile* didacticProfile = (TProfile*)_tempSwap.Clone(Form(
"Fit_%i_1-CDF_%.3f_cen_%i_%s", didacticCounterBest, CDF,
fInCentralitySelection, detector.Data()));
939 didacticProfile->GetListOfFunctions()->Add(didactifFit);
941 didacticCounterBest++;
942 TH2F* didacticSurface =
BookTH2F(Form(
"surface_%s", didacticProfile->GetName()),
"#phi",
"#eta", 50, 0, TMath::TwoPi(), 50, -1, 1, -1, kFALSE);
943 for(
Int_t i(0); i < iTracks; i++) {
944 AliVTrack* track =
static_cast<AliVTrack*
>(
fTracks->At(i));
946 if(
fUsePtWeight) didacticSurface->Fill(track->Phi(), track->Eta(), track->Pt());
947 else didacticSurface->Fill(track->Phi(), track->Eta());
951 TF2 *f2 =
new TF2(Form(
"%s_LJ", didacticSurface->GetName()),
"[0]*TMath::Gaus(x,[1],[2])*TMath::Gaus(y,[3],[4])", 0, TMath::TwoPi(), -1, 1);
952 f2->SetParameters(excludeInPt/3.,excludeInPhi,.1,excludeInEta,.1);
953 didacticSurface->GetListOfFunctions()->Add(f2);
963 static Int_t didacticCounterWorst(0);
965 TProfile* didacticProfile = (TProfile*)_tempSwap.Clone(Form(
"Fit_%i_1-CDF_%.3f_cen_%i_%s", didacticCounterWorst, CDF,
fInCentralitySelection, detector.Data() ));
967 didacticProfile->GetListOfFunctions()->Add(didactifFit);
969 didacticCounterWorst++;
991 #ifdef ALIANALYSISTASKLOCALRHO_DEBUG_FLAG_0
992 printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
1040 #ifdef ALIANALYSISTASKLOCALRHO_DEBUG_FLAG_0
1041 printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
1050 #ifdef ALIANALYSISTASKLOCALRHO_DEBUG_FLAG_0
1051 printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
1060 #ifdef ALIANALYSISTASKLOCALRHO_DEBUG_FLAG_0
1061 printf(
"__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
Double_t QCnS(Int_t i, Int_t j)
void CalculateEventPlaneCombinedVZERO(Double_t *comb) const
detectorType fDetectorType
Float_t fPercentageOfFits
virtual Bool_t Run()
Run function. This is the core function of the analysis and contains the user code. Therefore users have to implement this function.
void CalculateEventPlaneTPC(Double_t *tpc)
virtual void Terminate(Option_t *option)
virtual void UserCreateOutputObjects()
Bool_t PassesCuts(AliVTrack *track) const
AliJetContainer * GetJetContainer(Int_t i=0) const
Bool_t QCnRecovery(Double_t psi2, Double_t psi3)
Bool_t fAbsVnHarmonics
status of rho vs centrality
Double_t ChiSquareCDF(Int_t ndf, Double_t x) const
Double_t GetJetEtaMax() const
Bool_t fRebinSwapHistoOnTheFly
TH1F * fHistAnalysisSummary
swap histogram
Float_t fExcludeLeadingJetsFromFit
virtual ~AliAnalysisTaskLocalRho()
TH1F * fHistPsi2[10]
v3 cumulant
Double_t GetJetRadius(Int_t i=0) const
void CalculateEventPlaneVZERO(Double_t vzero[2][2]) const
TString fLocalRhoName
name for local rho
Bool_t CorrectRho(Double_t psi2, Double_t psi3)
TArrayI * fCentralityClasses
Int_t fNAcceptedTracksQCn
number of accepted tracks
AliEmcalJet * GetLeadingJet(const char *opt="")
void SetJetEtaLimits(Float_t min, Float_t max, Int_t c=0)
AliRhoParameter * fRho
! event rho
TList * fOutputListBad
output list for local analysis
Double_t CalculateQC4(Int_t harm)
void QCnQnk(Int_t n, Int_t k, Double_t &reQ, Double_t &imQ)
void ExecOnce()
Perform steps needed to initialize the analysis.
Double_t fCent
!event centrality
AliLocalRhoParameter * fLocalRho
! local event rho
Bool_t fUseV0EventPlaneFromHeader
fitModulationType fFitModulationType
centrality bin, only for QA plots
TClonesArray * fJets
! jets
void SetJetPhiLimits(Float_t min, Float_t max, Int_t c=0)
TList * fOutputListGood
output list
TH2F * fHistRhoStatusCent
cdf value of chisquare p
TProfile * fProfV3
v2 cumulant
TString fFitModulationOptions
Bool_t InitializeAnalysis()
Float_t GetJetRadius() const
void FillEventPlaneHistograms(Double_t psi2, Double_t psi3) const
TH2F * BookTH2F(const char *name, const char *x, const char *y, Int_t binsx, Double_t minx, Double_t maxx, Int_t binsy, Double_t miny, Double_t maxy, Int_t c=-1, Bool_t append=kTRUE)
TClonesArray * fTracks
!tracks
TProfile * fProfV3Cumulant
extracted v3
void FillAnalysisSummaryHistogram() const
void SetModulationFit(TF1 *fit)
Base task in the EMCAL jet framework.
Represent a jet reconstructed using the EMCal jet framework.
Bool_t fUsePtWeightErrorPropagation
TH1F * fHistSwap
output list for local analysis
Int_t fInCentralitySelection
accepted tracks for QCn
Int_t GetRunNumber(TString)
TProfile * fProfV2Cumulant
extracted v2
Double_t PhaseShift(Double_t x) const
TH1F * BookTH1F(const char *name, const char *x, Int_t bins, Double_t min, Double_t max, Int_t c=-1, Bool_t append=kTRUE)
Bool_t fAttachToEvent
is the analysis initialized?
AliAnalysisTaskLocalRho()
Double_t CalculateQC2(Int_t harm)
Bool_t fNoEventWeightsForQC
TH1F * fHistPsi3[10]
psi 2