6 #include <TClonesArray.h> 10 #include <THnSparse.h> 13 #include <TLorentzVector.h> 19 #include <AliAnalysisDataSlot.h> 20 #include <AliAnalysisDataContainer.h> 22 #include "TMatrixDSym.h" 23 #include "TMatrixDSymEigen.h" 26 #include "AliVCluster.h" 27 #include "AliVTrack.h" 32 #include "AliMCEvent.h" 33 #include "AliGenPythiaEventHeader.h" 34 #include "AliAODMCHeader.h" 35 #include "AliAODHandler.h" 36 #include "AliAODEvent.h" 37 #include "AliMCEvent.h" 38 #include "AliAnalysisManager.h" 43 #include "AliAODInputHandler.h" 44 #include "AliPIDResponse.h" 45 #include "AliCFContainer.h" 46 #include "AliCFManager.h" 47 #include "AliMultiEventInputHandler.h" 48 #include "AliKFParticle.h" 49 #include "AliKFVertex.h" 50 #include "AliMultSelection.h" 51 #include "AliAnalysisUtils.h" 53 #include "AliAODEvent.h" 76 fMinFractionShared(0),
79 fJetSelection(kInclusive),
85 fangWindowRecoil(0.6),
92 fOneConstSelectOn(kFALSE),
99 fAssITSrefitCut(kTRUE),
126 fnULSmLSpairsPerElectron(0),
128 fnElecOverPartPerJet(0),
131 fnIncSubPhotElecPerJet(0),
133 fnTrueHFElecPerJet(0),
134 fnTruePElecPerJet(0),
150 fTreeObservableTagging(0)
153 for(
Int_t i=0;i<21;i++){
157 for(
Int_t i=0;i<4;i++){
158 fptTrueHFE[i] = NULL;
161 for(
Int_t i=0;i<5;i++){
162 fInvmassLS[i] = NULL;
163 fInvmassULS[i] = NULL;
164 fUlsLsElecPt[i] = NULL;
165 fTotElecPt[i] = NULL;
166 for(
Int_t j=0;j<18;j++){
167 fnTPCTrueParticles[i][j] = NULL;
170 SetMakeGeneralHistograms(kTRUE);
171 DefineOutput(1, TList::Class());
172 DefineOutput(2, TTree::Class());
191 fMinFractionShared(0),
192 fJetShapeType(
kData),
193 fJetShapeSub(kNoSub),
194 fJetSelection(kInclusive),
195 fPtThreshold(-9999.),
200 fangWindowRecoil(0.6),
204 fCentSelectOn(kTRUE),
207 fOneConstSelectOn(kFALSE),
214 fAssITSrefitCut(kTRUE),
241 fnULSmLSpairsPerElectron(0),
243 fnElecOverPartPerJet(0),
246 fnIncSubPhotElecPerJet(0),
248 fnTrueHFElecPerJet(0),
249 fnTruePElecPerJet(0),
265 fTreeObservableTagging(0)
269 for(
Int_t i=0;i<21;i++){
273 for(
Int_t i=0;i<4;i++){
277 for(
Int_t i=0;i<5;i++){
282 for(
Int_t j=0;j<18;j++){
288 DefineOutput(1, TList::Class());
289 DefineOutput(2, TTree::Class());
306 Bool_t oldStatus = TH1::AddDirectoryStatus();
307 TH1::AddDirectory(kFALSE);
309 fh2ResponseUW=
new TH2F(
"fh2ResponseUW",
"fh2ResponseUW", 100, 0, 200, 100, 0, 200);
311 fh2ResponseW=
new TH2F(
"fh2ResponseW",
"fh2ResponseW", 100, 0, 200, 100, 0, 200);
313 fPhiJetCorr6=
new TH2F(
"fPhiJetCorr6",
"fPhiJetCorr6", 50, 0, 2*TMath::Pi(), 50, 0, 2*TMath::Pi());
315 fEtaJetCorr6=
new TH2F(
"fEtaJetCorr6",
"fEtaJetCorr6", 50, -1.5, 1.5, 50, -1.5, 1.5);
318 fPhiJetCorr7=
new TH2F(
"fPhiJetCorr7",
"fPhiJetCorr7", 50, 0, 2*TMath::Pi(), 50, 0, 2*TMath::Pi());
320 fEtaJetCorr7=
new TH2F(
"fEtaJetCorr7",
"fEtaJetCorr7", 50, -1.5, 1.5, 50, -1.5, 1.5);
323 fPtJetCorr=
new TH2F(
"fPtJetCorr",
"fPtJetCorr", 100, 0, 200, 100, 0, 200);
326 fPtJet=
new TH1F(
"fPtJet",
"fPtJet", 100, 0, 200);
329 fPtGenJet=
new TH1F(
"fPtGenJet",
"fPtGenJet", 100, 0, 200);
332 fPhiJet=
new TH2F(
"fPhiJet",
"fPhiJet", 100, 0, 200, 100, 0, TMath::TwoPi());
335 fEtaJet=
new TH2F(
"fEtaJet",
"fEtaJet", 100, 0, 200, 100, -1,1);
338 fhpTjetpT=
new TH2F(
"fhpTjetpT",
"fhpTjetpT", 200, 0, 200, 200, 0, 200);
341 fhPt=
new TH1F(
"fhPt",
"fhPt", 200, 0, 200);
344 fhPhi=
new TH1F(
"fhPhi",
"fhPhi", 100, -TMath::Pi(), TMath::Pi());
347 fNbOfConstvspT=
new TH2F(
"fNbOfConstvspT",
"fNbOfConstvspT", 100, 0, 100, 200, 0, 200);
350 fnTPCnTOFnocut=
new TH2F(
"fnTPCnTOFnocut",
"fnTPCnTOFnocut", 200, -10, 10, 200, -10, 10);
353 fnTPCnocut=
new TH2F(
"fnTPCnocut",
"fnTPCnocut", 50, 0, 5, 200, -10, 10);
356 fnTOFnocut=
new TH2F(
"fnTOFnocut",
"fnTOFnocut", 50, 0, 5, 200, -10, 10);
359 fnTPCcut=
new TH2F(
"fnTPCcut",
"fnTPCcut", 50, 0, 5, 200, -10, 10);
365 fnPartPerJet=
new TH1F(
"fnPartPerJet",
"fnPartPerJet", 500,0,500);
371 fnInclElecPerJet=
new TH1F(
"fnInclElecPerJet",
"fnInclElecPerJet", 100,0,100);
374 fnPhotElecPerJet=
new TH1F(
"fnPhotElecPerJet",
"fnPhotElecPerJet", 100,0,100);
380 fnTrueElecPerJet=
new TH1F(
"fnTrueElecPerJet",
"fnTrueElecPerJet", 100,0,100);
390 double xbins[60] = {0.01,0.1,0.12,0.14,0.16,0.18,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.55,0.6,0.65,0.7,0.75,
391 0.8,0.85,0.9,0.95,1,1.1,1.2,1.3,1.4,1.5,1.6,1.7,1.8,1.9,2,2.2,2.4,2.6,2.8,3,3.2,3.4,3.6,3.8,4,4.5,5,
392 5.5,6,6.5,7,8,9,10,11,12,13,14,15,16,18,20};
394 fPi0Pt =
new TH2F(
"fPi0Pt",
"fPi0Pt",4,0,4,nbin,xbins);
397 fEtaPt =
new TH2F(
"fEtaPt",
"fEtaPt",4,0,4,nbin,xbins);
400 fGenHfePt =
new TH1F(
"fGenHfePt",
"fGenHfePt",nbin,xbins);
403 fGenPePt =
new TH1F(
"fGenPePt",
"fGenPePt",nbin,xbins);
406 Double_t ptRange[19] = {0.5,0.7,0.9,1.1,1.3,1.5,1.7,1.9,2.1,2.3,2.5,2.7,2.9,3.1,3.3,3.5,3.7,3.9,5};
408 fPtP=
new TH2F(
"fPtP",
"fPtP", 18,ptRange,18,ptRange);
411 for(
Int_t i=0;i<5;i++){
413 fInvmassLS[i] =
new TH2F(Form(
"fInvmassLS%d",i), Form(
"fInvmassLS%d",i), 50, 0, 5, 100, 0, 0.5);
416 fInvmassULS[i] =
new TH2F(Form(
"fInvmassULS%d",i), Form(
"fInvmassULS%d",i), 50, 0, 5, 100, 0, 0.5);
419 fUlsLsElecPt[i] =
new TH1F(Form(
"fUlsLsElecPt%d",i), Form(
"fUlsLsElecPt%d",i), 18,ptRange);
422 fTotElecPt[i] =
new TH1F(Form(
"fTotElecPt%d",i), Form(
"fTotElecPt%d",i), 18,ptRange);
426 for(
Int_t j=0;j<18;j++){
427 fnTPCTrueParticles[i][j] =
new TH2F(Form(
"fnTPCTrueParticles%d%d",i,j),Form(
"fnTPCTrueParticles%d%d",i,j), 7,0,7,100,-15,15);
432 for(
Int_t i=0;i<4;i++){
433 fptTrueHFE[i] =
new TH1F(Form(
"fptTrueHFE%d",i), Form(
"fptTrueHFE%d",i), 18,ptRange);
437 fptJetIE=
new TH1F(
"fptJetIE",
"fptJetIE", 100, 0, 200);
440 fptJetPE=
new TH1F(
"fptJetPE",
"fptJetPE", 100, 0, 200);
443 fptJetHFE=
new TH1F(
"fptJetHFE",
"fptJetHFE", 100, 0, 200);
446 fptRecPE=
new TH1F(
"fptRecPE",
"fptRecPE", 100, 0, 50);
449 fptTruePE=
new TH1F(
"fptTruePE",
"fptTruePE", 100, 0, 50);
452 fptWrongPE=
new TH1F(
"fptWrongPE",
"fptWrongPE", 100, 0, 50);
455 fPhiTrack=
new TH2F(
"fPhiTrack",
"fPhiTrack", 100, 0, 50, 100, 0, TMath::TwoPi());
458 fEtaTrack=
new TH2F(
"fEtaTrack",
"fEtaTrack", 100, 0, 50, 100, -1,1);
461 fPhiElec=
new TH2F(
"fPhiElec",
"fPhiElec", 100, 0, 50, 100, 0, TMath::TwoPi());
464 fEtaElec=
new TH2F(
"fEtaElec",
"fEtaElec", 100, 0, 50, 100, -1,1);
474 THnSparse *hn =
dynamic_cast<THnSparse*
>(
fOutput->At(i));
479 TH1::AddDirectory(oldStatus);
486 fShapesVarNames[0] =
"partonCode";
487 fShapesVarNames[1] =
"ptJet";
488 fShapesVarNames[2] =
"ptDJet";
489 fShapesVarNames[3] =
"mJet";
491 fShapesVarNames[4] =
"angularity";
492 fShapesVarNames[5] =
"circularity";
493 fShapesVarNames[6] =
"lesub";
494 fShapesVarNames[7] =
"coronna";
496 fShapesVarNames[8] =
"ptJetMatch";
497 fShapesVarNames[9] =
"ptDJetMatch";
498 fShapesVarNames[10] =
"mJetMatch";
500 fShapesVarNames[11] =
"angularityMatch";
501 fShapesVarNames[12] =
"circularityMatch";
502 fShapesVarNames[13] =
"lesubMatch";
503 fShapesVarNames[14] =
"coronnaMatch";
504 fShapesVarNames[15]=
"weightPythia";
506 fShapesVarNames[16]=
"rhoVal";
507 fShapesVarNames[17]=
"rhoMassVal";
508 fShapesVarNames[18]=
"ptUnsub";
509 fShapesVarNames[19]=
"pElec";
510 fShapesVarNames[20]=
"ptElec";
511 fShapesVarNames[21]=
"nInclElec";
512 fShapesVarNames[22]=
"nPhotElec";
513 fShapesVarNames[23]=
"hasElec";
524 delete [] fShapesVarNames;
534 printf(
"ERROR: fAOD not available\n");
538 fVevent =
dynamic_cast<AliVEvent*
>(InputEvent());
540 printf(
"ERROR: fVEvent not available\n");
548 fMCarray =
dynamic_cast<TClonesArray*
>(
fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
549 fMCheader =
dynamic_cast<AliAODMCHeader*
>(
fAOD->GetList()->FindObject(AliAODMCHeader::StdBranchName()));
556 AliGenEventHeader* gh=(AliGenEventHeader*)lh->At(0);
557 NpureMC = gh->NProduced();
562 Bool_t isPileupfromSPDmulbins = kFALSE, isPileupFromMV = kFALSE;
563 isPileupfromSPDmulbins =
fAOD->IsPileupFromSPDInMultBins();
564 if (isPileupfromSPDmulbins)
return kFALSE;
566 AliAnalysisUtils utils;
567 utils.SetMinPlpContribMV(5);
568 utils.SetMaxPlpChi2MV(5.);
569 utils.SetMinWDistMV(15);
570 utils.SetCheckPlpFromDifferentBCMV(kFALSE);
571 isPileupFromMV = utils.IsPileUpMV(
fAOD);
572 if (isPileupFromMV)
return kFALSE;
576 Bool_t isPrimary = kFALSE, isFromLMdecay = kTRUE, isFromHFdecay=kTRUE;
583 for (
Int_t iParticle = 0; iParticle < nParticles; iParticle++) {
584 AliAODMCParticle* particle = (AliAODMCParticle*)
fMCarray->At(iParticle);
585 int fPDG = particle->GetPdgCode();
586 double pTMC = particle->Pt();
589 if(iParticle>=NpureMC)iEnhance = kTRUE;
604 if (TMath::Abs(etaMC)<0.8 && TMath::Abs(fPDG)==11){
606 if (isFromHFdecay)
fGenHfePt->Fill(pTMC);
607 if (iDecay>0 && iDecay<7)
fGenPePt->Fill(pTMC);
612 if (TMath::Abs(yMC)>1.0)
continue;
615 if(fPDG==111)
fPi0Pt->Fill(iEnhance,pTMC);
616 if(fPDG==221)
fEtaPt->Fill(iEnhance,pTMC);
618 if (!isFromHFdecay && !isFromLMdecay){
619 if(fPDG==111)
fPi0Pt->Fill(iEnhance+2,pTMC);
620 if(fPDG==221)
fEtaPt->Fill(iEnhance+2,pTMC);
640 AliAODTrack *triggerHadron = 0x0;
647 if (triggerHadronLabel==-99999) {
662 TClonesArray *TrackArray = NULL;
663 TClonesArray *TrackArrayMC = NULL;
665 else TrackArray = PartCont->GetArray();
667 else triggerHadron =
static_cast<AliAODTrack*
>(TrackArray->At(triggerHadronLabel));
671 if (!triggerHadron) {
682 fhPt->Fill(triggerHadron->Pt());
688 nVertices =
fAOD->GetNumberOfVertices();
689 Double_t listofmotherkink[nVertices];
690 Int_t nMotherKink = 0;
691 for(
Int_t ivertex=0; ivertex < nVertices; ivertex++) {
692 AliAODVertex *vertex =
fAOD->GetVertex(ivertex);
693 if(!vertex)
continue;
694 if(vertex->GetType()==AliAODVertex::kKink) {
695 AliAODTrack *mother = (AliAODTrack *) vertex->GetParent();
696 if(!mother)
continue;
697 Int_t idmother = mother->GetID();
698 listofmotherkink[nMotherKink] = idmother;
707 Float_t rhoVal=0, rhoMassVal = 0.;
710 jetCont->ResetCurrentID();
715 Printf(
"%s: Could not retrieve rho %s (some histograms will be filled with zero)!", GetName(), jetCont->
GetRhoName().Data());
716 }
else rhoVal = rhoParam->GetVal();
720 Printf(
"%s: Could not retrieve rho_m %s (some histograms will be filled with zero)!", GetName(), jetCont->
GetRhoMassName().Data());
721 }
else rhoMassVal = rhomParam->GetVal();
732 Int_t ifound=0, jfound=0;
733 Int_t ilab=-1, jlab=-1;
763 jetUS = jetContUS->
GetJet(i);
766 if(ifound==1) ilab = i;
769 if(ilab==-1)
continue;
770 jetUS=jetContUS->
GetJet(ilab);
776 Printf(
"jet2 does not exist, returning");
784 Printf(
"jet3 does not exist, returning");
813 jetUS = jetContUS->
GetJet(i);
816 if(ifound==1) ilab = i;
819 if(ilab==-1)
continue;
820 jetUS=jetContUS->
GetJet(ilab);
824 Printf(
"jet2 does not exist, returning");
830 jet3 = jetContPart->
GetJet(j);
834 if(jfound==1) jlab = j;
837 if(jlab==-1)
continue;
838 jet3=jetContPart->
GetJet(jlab);
840 Printf(
"jet3 does not exist, returning");
846 Printf(
"jet3 does not exist, returning");
868 Float_t dRp1 = TMath::Sqrt(jp1 * jp1 + detap1 * detap1);
878 dRp1 = TMath::Sqrt(jp1 * jp1 + detap1 * detap1);
921 Float_t ptMatch=0., ptDMatch=0., massMatch=0., angulMatch=0.,circMatch=0., lesubMatch=0., coronnaMatch=0;
983 Int_t nInclusiveElectrons = 0, nPhotonicElectrons = 0, nTrueElectronsMC, nTrueHFElecMC;
985 Bool_t hasElectrons = kFALSE;
989 GetNumberOfElectrons(jet1, 0,nMotherKink,listofmotherkink,nInclusiveElectrons,nPhotonicElectrons,pElec,ptElec,hasElectrons);
996 AliVParticle *vp1 = 0x0;
998 Int_t elecCounter = 0;
1005 Printf(
"AliVParticle associated to constituent not found");
1009 pdgCode = vp1->PdgCode();
1010 if (TMath::Abs(pdgCode)==11) elecCounter++;
1044 AliVParticle *vp1 = 0x0;
1045 Int_t nIE=0, nPairs=0, nPE=0, sub = -1;
1048 Bool_t hasElecCand = kFALSE;
1056 printf(
"ERROR: pid response not available\n");
1063 Printf(
"AliVParticle associated to constituent not found");
1067 AliVTrack *vtrack =
dynamic_cast<AliVTrack*
>(vp1);
1069 printf(
"ERROR: Could not receive track%d\n", i);
1073 AliAODTrack *track =
dynamic_cast<AliAODTrack*
>(vtrack);
1075 if (!track)
continue;
1078 Bool_t passTrackCut = kFALSE;
1080 if (!passTrackCut)
continue;
1084 Double_t p=-9., pt=-9., fTPCnSigma=-99., fTOFnSigma=-99., phi = -9., eta = -99.;
1093 fTPCnSigma =
fpidResponse->NumberOfSigmasTPC(track, AliPID::kElectron);
1094 fTOFnSigma =
fpidResponse->NumberOfSigmasTOF(track, AliPID::kElectron);
1101 if(TMath::Abs(fTPCnSigma)<3.5){
1102 hasElecCand = kTRUE;
1114 if (nPairs>0) nPE++;
1130 hasElec = hasElecCand;
1136 AliVParticle *vp1 = 0x0;
1137 Int_t nIE=0, nHFE=0, nPE=0, PartId=0, nPairs=0, iDecay = 0;
1138 Double_t p=-9., pt=-9., fTPCnSigma=-99., fTOFnSigma=-99., MCweight = 1.;
1139 Double_t ptRange[19] = {0.5,0.7,0.9,1.1,1.3,1.5,1.7,1.9,2.1,2.3,2.5,2.7,2.9,3.1,3.3,3.5,3.7,3.9,5};
1140 Double_t ptJetRange[6] = {5,15,25,35,45,60};
1142 Bool_t isFromHFdecay=kFALSE;
1143 Bool_t isFromLMdecay=kFALSE;
1144 Bool_t passTrackCut = kFALSE;
1154 Printf(
"AliVParticle associated to constituent not found");
1158 AliVTrack *vtrack =
dynamic_cast<AliVTrack*
>(vp1);
1160 printf(
"ERROR: Could not receive track%d\n", i);
1164 AliAODTrack *track =
dynamic_cast<AliAODTrack*
>(vtrack);
1166 passTrackCut = kFALSE;
1167 isFromHFdecay=kFALSE;
1168 isFromLMdecay=kFALSE;
1173 fTPCnSigma =
fpidResponse->NumberOfSigmasTPC(track, AliPID::kElectron);
1174 fTOFnSigma =
fpidResponse->NumberOfSigmasTOF(track, AliPID::kElectron);
1181 Int_t label = track->GetLabel();
1191 if (TMath::Abs(partPDG)==11 && isFromHFdecay)
fptTrueHFE[0]->Fill(pt);
1195 if (!passTrackCut)
continue;
1197 if (TMath::Abs(partPDG)==11 && isFromHFdecay)
fptTrueHFE[1]->Fill(pt);
1201 if (TMath::Abs(partPDG)==11) PartId = 1;
1202 if (TMath::Abs(partPDG)==13) PartId = 2;
1203 if (TMath::Abs(partPDG)==321) PartId = 3;
1204 if (TMath::Abs(partPDG)==2212) PartId = 4;
1205 if (TMath::Abs(partPDG)==211) PartId = 5;
1207 for (
Int_t l=0;l<5;l++){
1208 for(
Int_t k=0;k<18;k++){
1209 if (jetPt>=ptJetRange[l] && jetPt<ptJetRange[l+1] && p>=ptRange[k] && p<ptRange[k+1])
fnTPCTrueParticles[l][k]->Fill(PartId,fTPCnSigma);
1214 if (TMath::Abs(partPDG)!=11)
continue;
1221 if (isFromHFdecay) nHFE++;
1222 if (isFromLMdecay) nPE++;
1227 if (nPairs>0)
fptRecPE->Fill(pt,MCweight);
1228 if (nPairs>0 && (iDecay<1 || iDecay>6))
fptWrongPE->Fill(pt,MCweight);
1230 if (iDecay>0 && iDecay<7){
1232 for (
Int_t l=0;l<5;l++){
1234 if (jetPt>=ptJetRange[l] && jetPt<ptJetRange[l+1]){
1247 if (nIE==1 && nPE==1)
fptJetPE->Fill(jet->
Pt());
1248 if (nIE==1 && nHFE==1)
fptJetHFE->Fill(jet->
Pt());
1265 ntracks =
fVevent->GetNumberOfTracks();
1267 Int_t nULSpairs = 0;
1271 Double_t ptJetRange[6] = {5,15,25,35,45,60};
1274 for(
Int_t jTracks = 0; jTracks < ntracks; jTracks++){
1276 AliVParticle* Vassotrack =
fVevent->GetTrack(jTracks);
1279 printf(
"ERROR: Could not receive associated track %d\n", jTracks);
1282 AliAODTrack *trackAsso =
dynamic_cast<AliAODTrack*
>(Vassotrack);
1283 if(!trackAsso)
continue;
1285 if((track->GetID())==(trackAsso->GetID()))
continue;
1288 Bool_t passAssoTrackCut = kFALSE;
1290 if (!passAssoTrackCut)
continue;
1293 Double_t p=-9.,pt=-9.,eta =-9.,phi=-9.;
1300 charge = track->Charge();
1304 Double_t pAsso=-9.,ptAsso=-9.,etaAsso =-9.,fITSnSigmaAsso=-9.,fTPCnSigmaAsso=-9.,phiAsso=-9.;
1305 Int_t chargeAsso = 0;
1307 ptAsso = trackAsso->Pt();
1308 pAsso = trackAsso->P();
1309 phiAsso = trackAsso->Phi();
1310 etaAsso = trackAsso->Eta();
1311 chargeAsso = trackAsso->Charge();
1315 fITSnSigmaAsso =
fpidResponse->NumberOfSigmasITS(trackAsso, AliPID::kElectron);
1316 fTPCnSigmaAsso =
fpidResponse->NumberOfSigmasTPC(trackAsso, AliPID::kElectron);
1318 if(TMath::Abs(fTPCnSigmaAsso)>3)
continue;
1321 Bool_t fFlagLS=kFALSE, fFlagULS=kFALSE;
1322 Double_t openingAngle = -999.,
mass=999., width = -999;
1325 if(charge>0) fPDGe1 = -11;
1326 if(chargeAsso>0) fPDGe2 = -11;
1328 if(charge == chargeAsso) fFlagLS = kTRUE;
1329 if(charge != chargeAsso) fFlagULS = kTRUE;
1331 AliKFParticle::SetField(
fVevent->GetMagneticField());
1332 AliKFParticle ge1(*track, fPDGe1);
1333 AliKFParticle ge2(*trackAsso, fPDGe2);
1334 AliKFParticle recg(ge1, ge2);
1336 if(recg.GetNDF()<1)
continue;
1337 Double_t chi2recg = recg.GetChi2()/recg.GetNDF();
1338 if(TMath::Sqrt(TMath::Abs(chi2recg))>3.)
continue;
1340 openingAngle = ge1.GetAngle(ge2);
1343 Int_t MassCorrect=-9;
1344 MassCorrect = recg.GetMass(
mass,width);
1346 for (
Int_t l=0;l<5;l++){
1347 if (jetPt>=ptJetRange[l] && jetPt<ptJetRange[l+1] && fFlagULS)
fInvmassULS[l]->Fill(pt,
mass);
1348 if (jetPt>=ptJetRange[l] && jetPt<ptJetRange[l+1] && fFlagLS)
fInvmassLS[l]->Fill(pt,
mass);
1359 sub = nULSpairs-nLSpairs;
1370 Bool_t isHFdecay = kFALSE;
1373 Int_t idMother = particle->GetMother();
1375 AliAODMCParticle* mother = (AliAODMCParticle*)
fMCarray->At(idMother);
1376 Int_t motherPDG = TMath::Abs(mother->GetPdgCode());
1379 if((motherPDG % 1000) / 100 == 4) isHFdecay = kTRUE;
1380 if(motherPDG / 1000 == 4) isHFdecay = kTRUE;
1383 if((motherPDG % 1000) / 100 == 5) isHFdecay = kTRUE;
1384 if(motherPDG / 1000 == 5) isHFdecay = kTRUE;
1387 if(motherPDG==4 || motherPDG==5 || motherPDG == 211 || motherPDG ==13 || motherPDG ==2112 || motherPDG ==130 || motherPDG ==3122 ||
1388 motherPDG ==310 || motherPDG ==3222 || motherPDG ==2212 || motherPDG ==3112 || motherPDG ==321 ||
1389 motherPDG ==11 || motherPDG ==3212 || motherPDG ==311 || motherPDG ==20213 || motherPDG ==3312 ||
1390 motherPDG ==3334 || motherPDG ==3324 || motherPDG ==3322 || motherPDG ==1000010020 || motherPDG ==15
1391 || motherPDG ==10323 || motherPDG ==2114 || motherPDG ==1000010030 || motherPDG ==2214 || motherPDG ==2224
1392 || motherPDG ==1114 || motherPDG == 2214 || motherPDG == 3114 || motherPDG ==3224 || motherPDG ==3124
1393 || motherPDG ==3314 || motherPDG ==10323 || motherPDG == 3214) isHFdecay = kTRUE;
1403 Bool_t isLMdecay = kFALSE;
1406 Int_t idMother = particle->GetMother();
1408 AliAODMCParticle* mother = (AliAODMCParticle*)
fMCarray->At(idMother);
1409 Int_t motherPDG = TMath::Abs(mother->GetPdgCode());
1411 if(motherPDG == 111 || motherPDG == 221 || motherPDG==223 || motherPDG==333 || motherPDG==331 ||
1412 motherPDG==113 || motherPDG==213 || motherPDG==313 || motherPDG==323) isLMdecay = kTRUE;
1421 Bool_t isprimary = kFALSE;
1425 if (particle->IsPrimary()) isprimary = kTRUE;
1433 double weight = 1.0;
1435 double parPi0_enh[5] = {0.530499,0.732775,0.000997414,3.46894,4.84342};
1436 if (
fMCweight==1) weight = parPi0_enh[0] / TMath::Power((exp(parPi0_enh[1]*mcPi0pT - parPi0_enh[2]*mcPi0pT*mcPi0pT) + (mcPi0pT/parPi0_enh[3])) , parPi0_enh[4]);
1444 double weight = 1.0;
1446 double parEta_enh[5] = {0.78512,-0.606822,0.0326254,3.13959,4.83715};
1447 if (
fMCweight==1) weight = parEta_enh[0] / TMath::Power((exp(parEta_enh[1]*mcEtapT - parEta_enh[2]*mcEtapT*mcEtapT) + (mcEtapT/parEta_enh[3])), parEta_enh[4]);
1457 Int_t partPDG = particle->GetPdgCode();
1459 if (TMath::Abs(partPDG)==11){
1460 Int_t idMother = particle->GetMother();
1463 AliAODMCParticle* mother = (AliAODMCParticle*)
fMCarray->At(idMother);
1464 Int_t motherPDG = mother->GetPdgCode();
1471 if (motherPDG==111 && (isMotherPrimary || (!isMotherFromHF && !isMotherFromLM))){
1476 if (motherPDG==221 && (isMotherPrimary || (!isMotherFromHF && !isMotherFromLM))){
1482 Int_t idSecondMother = mother->GetMother();
1484 if (idSecondMother>0){
1485 AliAODMCParticle* secondMother = (AliAODMCParticle*)
fMCarray->At(idSecondMother);
1486 Int_t secondMotherPDG = secondMother->GetPdgCode();
1487 Double_t secondMotherPt = secondMother->Pt();
1493 if (motherPDG==22 && secondMotherPDG==111 && (isSecondMotherPrimary || (!isSecondMotherFromHF && !isSecondMotherFromLM))){
1498 if (motherPDG==22 && secondMotherPDG==221 && (isSecondMotherPrimary || (!isSecondMotherFromHF && !isSecondMotherFromLM))){
1503 if (motherPDG==111 && secondMotherPDG==221 && (isSecondMotherPrimary || (!isSecondMotherFromHF && !isSecondMotherFromLM))){
1508 Int_t idThirdMother = secondMother->GetMother();
1509 if (idThirdMother>0){
1510 AliAODMCParticle* thirdMother = (AliAODMCParticle*)
fMCarray->At(idThirdMother);
1511 Int_t thirdMotherPDG = thirdMother->GetPdgCode();
1512 Double_t thirdMotherPt = thirdMother->Pt();
1518 if (motherPDG==22 && secondMotherPDG==111 && thirdMotherPDG==221 && (isThirdMotherPrimary || (!isThirdMotherFromHF && !isThirdMotherFromLM))){
1537 AliVParticle *vp1 = 0x0;
1542 Printf(
"AliVParticle associated to constituent not found");
1547 Double_t dr2 = (vp1->Eta()-jet->
Eta())*(vp1->Eta()-jet->
Eta()) + dphi*dphi;
1549 num=num+vp1->Pt()*dr;
1583 TClonesArray *TracksArray = NULL;
1584 TClonesArray *TracksArrayMC = NULL;
1587 else TracksArray = PartCont->GetArray();
1590 if(!PartContMC || !TracksArrayMC)
return -2;
1593 if(!PartCont || !TracksArray)
return -2;
1597 AliAODTrack *
Track = 0x0;
1601 else NTracks = TracksArray->GetEntriesFast();
1603 for(
Int_t i=0; i < NTracks; i++){
1606 if (!Track)
continue;
1607 if(TMath::Abs(Track->Eta())>0.9)
continue;
1609 Double_t dr2 = (Track->Eta()-jet->
Eta())*(Track->Eta()-jet->
Eta()) + dphi*dphi;
1611 if((dr>=0.8) && (dr<1)) sumpt=sumpt+Track->Pt();
1615 if((Track = static_cast<AliAODTrack*>(PartCont->
GetAcceptTrack(i)))){
1616 if (!Track)
continue;
1617 if(TMath::Abs(Track->Eta())>0.9)
continue;
1618 if (Track->Pt()<0.15)
continue;
1620 Double_t dr2 = (Track->Eta()-jet->
Eta())*(Track->Eta()-jet->
Eta()) + dphi*dphi;
1622 if((dr>=0.8) && (dr<1)) sumpt=sumpt+Track->Pt();
1635 return Coronna(jet, jetContNb);
1647 AliVParticle *vp1 = 0x0;
1652 Printf(
"AliVParticle associated to constituent not found");
1656 num=num+vp1->Pt()*vp1->Pt();
1659 return TMath::Sqrt(num)/den;
1669 return PTD(jet, jetContNb);
1690 TVector3 ppJ1(pxjet, pyjet, pzjet);
1691 TVector3 ppJ3(- pxjet* pzjet, - pyjet * pzjet, pxjet * pxjet + pyjet * pyjet);
1693 TVector3 ppJ2(-pyjet, pxjet, 0);
1695 AliVParticle *vp1 = 0x0;
1700 Printf(
"AliVParticle associated to constituent not found");
1704 TVector3 pp(vp1->Px(), vp1->Py(), vp1->Pz());
1707 TVector3 pLong = pp.Dot(ppJ1) / ppJ1.Mag2() * ppJ1;
1708 TVector3 pPerp = pp - pLong;
1711 Float_t ppjX = pPerp.Dot(ppJ2);
1712 Float_t ppjY = pPerp.Dot(ppJ3);
1713 Float_t ppjT = TMath::Sqrt(ppjX * ppjX + ppjY * ppjY);
1714 if(ppjT<=0)
return 0;
1716 mxx += (ppjX * ppjX / ppjT);
1717 myy += (ppjY * ppjY / ppjT);
1718 mxy += (ppjX * ppjY / ppjT);
1723 if(sump2==0)
return 0;
1725 Double_t ele[4] = {mxx / sump2, mxy / sump2, mxy / sump2, myy / sump2};
1726 TMatrixDSym m0(2,ele);
1729 TMatrixDSymEigen m(m0);
1731 TMatrixD evecm = m.GetEigenVectors();
1732 eval = m.GetEigenValues();
1736 if (eval[0] < eval[1]) jev = 1;
1739 evec0 = TMatrixDColumn(evecm, jev);
1742 TVector2 evec(compx, compy);
1744 if(jev==1) circ=2*eval[0];
1745 if(jev==0) circ=2*eval[1];
1773 AliVParticle *vp1 = 0x0;
1774 AliVParticle *vp2 = 0x0;
1775 std::vector<int> ordindex;
1780 if(ordindex.size()<2)
return -1;
1784 Printf(
"AliVParticle associated to Leading constituent not found");
1790 Printf(
"AliVParticle associated to Subleading constituent not found");
1810 return LeSub(jet, jetContNb);
1839 AliVParticle *vp1 = 0x0;
1844 Printf(
"AliVParticle associated to constituent not found");
1852 mxx += ppt*ppt*deta*deta;
1853 myy += ppt*ppt*dphi*dphi;
1854 mxy -= ppt*ppt*deta*TMath::Abs(dphi);
1860 if(sump2==0)
return 0;
1862 Double_t ele[4] = {mxx , mxy , mxy , myy };
1863 TMatrixDSym m0(2,ele);
1866 TMatrixDSymEigen m(m0);
1868 TMatrixD evecm = m.GetEigenVectors();
1869 eval = m.GetEigenValues();
1873 if (eval[0] < eval[1]) jev = 1;
1876 evec0 = TMatrixDColumn(evecm, jev);
1879 TVector2 evec(compx, compy);
1881 if(jev==1) sig=TMath::Sqrt(TMath::Abs(eval[0])/sump2);
1882 if(jev==0) sig=TMath::Sqrt(TMath::Abs(eval[1])/sump2);
1896 return Sigma2(jet, jetContNb);
1916 TClonesArray *TracksArray = NULL;
1917 TClonesArray *TracksArrayMC = NULL;
1920 else TracksArray = PartCont->GetArray();
1923 if(!PartContMC || !TracksArrayMC)
return -99999;
1926 if(!PartCont || !TracksArray)
return -99999;
1930 AliAODTrack *
Track = 0x0;
1935 Int_t triggers[100];
1936 for (
Int_t iTrigger=0; iTrigger<100; iTrigger++) triggers[iTrigger] = 0;
1940 else NTracks = TracksArray->GetEntriesFast();
1942 for(
Int_t i=0; i < NTracks; i++){
1945 if (!Track)
continue;
1946 if(TMath::Abs(Track->Eta())>0.9)
continue;
1947 if (Track->Pt()<0.15)
continue;
1948 if ((Track->Pt() >= minpT) && (Track->Pt()< maxpT)) {
1955 if((Track = static_cast<AliAODTrack*>(PartCont->
GetAcceptTrack(i)))){
1956 if (!Track)
continue;
1957 if(TMath::Abs(Track->Eta())>0.9)
continue;
1958 if (Track->Pt()<0.15)
continue;
1959 if ((Track->Pt() >= minpT) && (Track->Pt()< maxpT)) {
1968 if (iTT == 0)
return -99999;
1969 Int_t nbRn = 0, index = 0 ;
1970 TRandom3* random =
new TRandom3(0);
1971 nbRn = random->Integer(iTT);
1972 index = triggers[nbRn];
1981 if (vphi < -1*TMath::Pi()) vphi += (2*TMath::Pi());
1982 else if (vphi > TMath::Pi()) vphi -= (2*TMath::Pi());
1983 if (mphi < -1*TMath::Pi()) mphi += (2*TMath::Pi());
1984 else if (mphi > TMath::Pi()) mphi -= (2*TMath::Pi());
1985 double dphi = mphi-vphi;
1986 if (dphi < -1*TMath::Pi()) dphi += (2*TMath::Pi());
1987 else if (dphi > TMath::Pi()) dphi -= (2*TMath::Pi());
2008 if(!ietrack->TestFilterMask(AliAODTrack::kTrkGlobalNoDCA))
return kFALSE;
2009 if(TMath::Abs(ietrack->Eta())>0.8)
return kFALSE;
2010 if(ietrack->GetTPCNcls() <
fTPCnCut)
return kFALSE;
2011 if (ietrack->GetITSNcls() <
fITSncut)
return kFALSE;
2012 if(!ietrack->IsOn(AliAODTrack::kITSrefit))
return kFALSE;
2013 if(!ietrack->IsOn(AliAODTrack::kTPCrefit))
return kFALSE;
2014 if(!(ietrack->HasPointOnITSLayer(0) && ietrack->HasPointOnITSLayer(1)))
return kFALSE;
2015 if ((ietrack->Pt()<0.5) || (ietrack->Pt()>4))
return kFALSE;
2017 Bool_t kinkmotherpass = kTRUE;
2018 for(
Int_t kinkmother = 0; kinkmother < nMother; kinkmother++) {
2019 if(ietrack->GetID() == listMother[kinkmother]) {
2020 kinkmotherpass = kFALSE;
2024 if(!kinkmotherpass)
return kFALSE;
2026 Double_t d0z0[2]={-999,-999}, cov[3];
2028 if(ietrack->PropagateToDCA(pVietx,
fVevent->GetMagneticField(), 20., d0z0, cov))
2029 if(TMath::Abs(d0z0[0]) >
fDcaXYcut || TMath::Abs(d0z0[1]) >
fDcaZcut)
return kFALSE;
2039 if(!aetrack->TestFilterMask(AliAODTrack::kTrkTPCOnly))
return kFALSE;
2040 if(aetrack->Pt() <
fAssPtCut)
return kFALSE;
2041 if(TMath::Abs(aetrack->Eta())>0.9)
return kFALSE;
2042 if(aetrack->GetTPCNcls() <
fAssTPCnCut)
return kFALSE;
2043 if (
fAssITSrefitCut && !(aetrack->GetStatus()&AliAODTrack::kITSrefit))
return kFALSE;
2044 if(!(aetrack->GetStatus()&AliAODTrack::kTPCrefit))
return kFALSE;
Float_t Angularity(AliEmcalJet *jet, Int_t jetContNb)
TClonesArray * fMCarray
MC particle.
Double_t GetFirstOrderSubtractedAngularity() const
Float_t GetJetLeSub(AliEmcalJet *jet, Int_t jetContNb)
Bool_t IsFromHFdecay(AliAODMCParticle *particle)
const TString & GetRhoName() const
Double_t GetSecondOrderSubtractedSigma2() const
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.
Bool_t FillHistograms()
Function filling histograms.
Float_t GetPartonEta6() const
AliEmcalJet * ClosestJet() const
AliJetContainer * GetJetContainer(Int_t i=0) const
Float_t GetPartonEta7() const
Double_t GetSecondOrderSubtractedConstituent() const
Float_t GetPythiaEventWeight() const
Container with name, TClonesArray and cuts for particles.
Double_t GetPi0weight(Double_t mcPi0pT) const
TH2F * fnTPCTrueParticles[5][18]
Float_t GetPartonPhi7() const
Bool_t RetrieveEventObjects()
Retrieve common objects from event.
Bool_t PhotElecTrackCuts(const AliVVertex *pVtx, AliAODTrack *aetrack, Int_t nMother, Double_t listMother[])
Float_t Circularity(AliEmcalJet *jet, Int_t jetContNb)
Double_t GetSecondOrderSubtractedpTD() const
Bool_t RetrieveEventObjects()
Double_t GetFirstOrderSubtractedLeSub() const
AliAnalysisTaskEmcalHfeTagging()
Int_t GetNumberOfPairs(AliEmcalJet *jet, AliAODTrack *track, const AliVVertex *pVtx, Int_t nMother, Double_t listMother[])
void UserCreateOutputObjects()
Float_t GetJetCircularity(AliEmcalJet *jet, Int_t jetContNb)
Int_t GetPartonFlag7() const
Container for particles within the EMCAL framework.
Float_t Sigma2(AliEmcalJet *jet, Int_t jetContNb)
Float_t GetPartonPhi6() const
Int_t TrackAt(Int_t idx) const
Double_t GetFirstOrderSubtractedpTD() const
UShort_t GetNumberOfTracks() const
TString kData
Declare data MC or deltaAOD.
Float_t PTD(AliEmcalJet *jet, Int_t jetContNb)
AliParticleContainer * GetParticleContainer(Int_t i=0) const
Get particle container attached to this task.
const TString & GetRhoMassName() const
Double_t GetSecondOrderSubtractedLeSub() const
Int_t GetPartonFlag6() const
Double_t GetEtaweight(Double_t mcEtapT) const
AliParticleContainer * GetParticleContainer() const
Double_t GetFirstOrderSubtractedConstituent() const
JetSelectionType fJetSelection
Bool_t IsPrimary(AliAODMCParticle *particle)
TH2F * fnULSmLSpairsPerElectron
Float_t GetJetMass(AliEmcalJet *jet, Int_t jetContNb)
Double_t GetSecondOrderSubtractedAngularity() const
Bool_t InclElecTrackCuts(const AliVVertex *pVtx, AliAODTrack *ietrack, Int_t nMother, Double_t listMother[])
Double_t RelativePhi(Double_t mphi, Double_t vphi)
virtual AliVParticle * GetAcceptParticle(Int_t i=-1) const
Double_t fCent
!event centrality
AliVEvent * fVevent
AOD object.
Double_t GetSecondOrderSubtracted() const
AliEmcalJet * GetNextAcceptJet()
void Terminate(Option_t *option)
void GetNumberOfElectrons(AliEmcalJet *jet, Int_t jetContNb, Int_t nMother, Double_t listMother[], Int_t &nIncElec, Int_t &nPhotElec, Double_t &pElec, Double_t &ptElec, Bool_t &hasElec)
Bool_t Data(TH1F *h, Double_t *rangefit, Bool_t writefit, Double_t &sgn, Double_t &errsgn, Double_t &bkg, Double_t &errbkg, Double_t &sgnf, Double_t &errsgnf, Double_t &sigmafit, Int_t &status)
Float_t LeSub(AliEmcalJet *jet, Int_t jetContNb)
virtual ~AliAnalysisTaskEmcalHfeTagging()
Double_t GetRhoVal(Int_t i=0) const
void GetWeightAndDecay(AliAODMCParticle *particle, Int_t &decay, Double_t &weight)
Double_t GetFirstOrderSubtractedCircularity() const
Float_t fMinFractionShared
AliAODMCParticle * fMCparticle
stack
AliEmcalList * fOutput
!output list
TH1F * fnIncSubPhotElecPerJet
Float_t GetJetCoronna(AliEmcalJet *jet, Int_t jetContNb)
Double_t GetSecondOrderSubtractedCircularity() const
AliTrackContainer * GetTrackContainer(Int_t i=0) const
Store some informaion about a Pythia eventThis class is used to store some information about a Pythia...
void SetMakeGeneralHistograms(Bool_t g)
Float_t GetJetAngularity(AliEmcalJet *jet, Int_t jetContNb)
Float_t GetSigma2(AliEmcalJet *jet, Int_t jetContNb)
virtual AliVTrack * GetAcceptTrack(Int_t i=-1) const
Base task in the EMCAL jet framework.
Represent a jet reconstructed using the EMCal jet framework.
void GetNumberOfTrueElectrons(AliEmcalJet *jet, Int_t jetContNb, Int_t nMother, Double_t listMother[], Int_t &nTrueElec, Int_t &nTrueHFElec)
const AliEmcalPythiaInfo * GetPythiaInfo() const
AliAODMCHeader * fMCheader
MC array.
Double_t GetFractionSharedPt(const AliEmcalJet *jet, AliParticleContainer *cont2=0x0) const
Declaration of class AliEmcalPythiaInfo.
JetShapeType fJetShapeType
void UserCreateOutputObjects()
Main initialization function on the worker.
Double_t GetFirstOrderSubtractedSigma2() const
Bool_t IsFromLMdecay(AliAODMCParticle *particle)
AliEmcalJetShapeProperties * GetShapeProperties() const
Int_t SelectTrigger(Float_t minpT, Float_t maxpT)
TH1F * fnElecOverPartPerJet
AliPIDResponse * fpidResponse
VEvent.
Float_t GetJetpTD(AliEmcalJet *jet, Int_t jetContNb)
const AliVVertex * spdVtx
Float_t Coronna(AliEmcalJet *jet, Int_t jetContNb)
Float_t GetPartonPt7() const
TTree * fTreeObservableTagging
Double_t GetFirstOrderSubtracted() const
Container for jet within the EMCAL jet framework.
std::vector< int > GetPtSortedTrackConstituentIndexes(TClonesArray *tracks) const
Float_t GetPartonPt6() const
Float_t GetJetNumberOfConstituents(AliEmcalJet *jet, Int_t jetContNb)
TH1F * fnTrueHFElecPerJet
AliEmcalJet * GetJet(Int_t i) const