X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=PWG4%2FJetTasks%2FAliAnalysisTaskJetSpectrum2.cxx;h=8cd6b4f1d8fa7badeab2674ca2d0710f0c1e0a42;hb=c8eabe2400d8b6b0ae9ab54a7236d90b21bb60c6;hp=67ce106e39e832985665db2dcb9f24b1398df055;hpb=bb3d13aa372b570062cd719290700f50db38045d;p=u%2Fmrichter%2FAliRoot.git diff --git a/PWG4/JetTasks/AliAnalysisTaskJetSpectrum2.cxx b/PWG4/JetTasks/AliAnalysisTaskJetSpectrum2.cxx index 67ce106e39e..8cd6b4f1d8f 100644 --- a/PWG4/JetTasks/AliAnalysisTaskJetSpectrum2.cxx +++ b/PWG4/JetTasks/AliAnalysisTaskJetSpectrum2.cxx @@ -63,25 +63,29 @@ ClassImp(AliAnalysisTaskJetSpectrum2) AliAnalysisTaskJetSpectrum2::AliAnalysisTaskJetSpectrum2(): AliAnalysisTaskSE(), - fJetHeaderRec(0x0), - fJetHeaderGen(0x0), - fAOD(0x0), - fhnCorrelation(0x0), - fBranchRec("jets"), - fBranchGen(""), - fUseAODJetInput(kFALSE), - fUseAODTrackInput(kFALSE), - fUseAODMCInput(kFALSE), - fUseGlobalSelection(kFALSE), - fUseExternalWeightOnly(kFALSE), - fLimitGenJetEta(kFALSE), - fFilterMask(0), - fAnalysisType(0), + fJetHeaderRec(0x0), + fJetHeaderGen(0x0), + fAOD(0x0), + fhnCorrelation(0x0), + fhnCorrelationPhiZRec(0x0), + f1PtScale(0x0), + fBranchRec("jets"), + fBranchGen(""), + fUseAODJetInput(kFALSE), + fUseAODTrackInput(kFALSE), + fUseAODMCInput(kFALSE), + fUseGlobalSelection(kFALSE), + fUseExternalWeightOnly(kFALSE), + fLimitGenJetEta(kFALSE), + fFilterMask(0), + fAnalysisType(0), fTrackTypeRec(kTrackUndef), fTrackTypeGen(kTrackUndef), fAvgTrials(1), fExternalWeight(1), - fRecEtaWindow(0.5), + fRecEtaWindow(0.5), + fMinJetPt(0), + fDeltaPhiWindow(20./180.*TMath::Pi()), fh1Xsec(0x0), fh1Trials(0x0), fh1PtHard(0x0), @@ -92,6 +96,7 @@ AliAnalysisTaskJetSpectrum2::AliAnalysisTaskJetSpectrum2(): AliAnalysisTaskSE(), fh1PtTrackRec(0x0), fh1SumPtTrackRec(0x0), fh1SumPtTrackAreaRec(0x0), + fh1TmpRho(0x0), fh1PtJetsRecIn(0x0), fh1PtJetsLeadingRecIn(0x0), fh1PtTracksRecIn(0x0), @@ -104,13 +109,23 @@ AliAnalysisTaskJetSpectrum2::AliAnalysisTaskJetSpectrum2(): AliAnalysisTaskSE(), fh2TracksLeadingPhiEta(0x0), fh2TracksLeadingPhiPt(0x0), fh2TracksLeadingJetPhiPt(0x0), + fh2JetPtJetPhi(0x0), + fh2TrackPtTrackPhi(0x0), + fh2DijetDeltaPhiPt(0x0), + fh2DijetAsymPt(0x0), + fh2DijetAsymPtCut(0x0), + fh2DijetDeltaPhiDeltaEta(0x0), + fh2DijetPt2vsPt1(0x0), + fh2DijetDifvsSum(0x0), + fh1DijetMinv(0x0), + fh1DijetMinvCut(0x0), fHistList(0x0) { for(int i = 0;i < kMaxStep*2;++i){ fhnJetContainer[i] = 0; } for(int i = 0;i < kMaxJets;++i){ - fh2PhiPt[i] = fh2PhiEta[i] = fh2FragRec[i] = fh2FragLnRec[i] = fh2FragGen[i] = fh2FragLnGen[i] = 0; + fh2PhiPt[i] = fh2PhiEta[i] = fh2RhoPtRec[i] = fh2PsiPtRec[i] = fh2FragRec[i] = fh2PsiPtGen[i] = fh2FragGen[i] = fh2FragLnRec[i] = fh2FragGen[i] = fh2FragLnGen[i] = 0; fh1PtRecIn[i] = fh1PtGenIn[i] = 0; } @@ -122,6 +137,8 @@ AliAnalysisTaskJetSpectrum2::AliAnalysisTaskJetSpectrum2(const char* name): fJetHeaderGen(0x0), fAOD(0x0), fhnCorrelation(0x0), + fhnCorrelationPhiZRec(0x0), + f1PtScale(0x0), fBranchRec("jets"), fBranchGen(""), fUseAODJetInput(kFALSE), @@ -137,6 +154,8 @@ AliAnalysisTaskJetSpectrum2::AliAnalysisTaskJetSpectrum2(const char* name): fAvgTrials(1), fExternalWeight(1), fRecEtaWindow(0.5), + fMinJetPt(0), + fDeltaPhiWindow(20./180.*TMath::Pi()), fh1Xsec(0x0), fh1Trials(0x0), fh1PtHard(0x0), @@ -147,6 +166,7 @@ AliAnalysisTaskJetSpectrum2::AliAnalysisTaskJetSpectrum2(const char* name): fh1PtTrackRec(0x0), fh1SumPtTrackRec(0x0), fh1SumPtTrackAreaRec(0x0), + fh1TmpRho(0x0), fh1PtJetsRecIn(0x0), fh1PtJetsLeadingRecIn(0x0), fh1PtTracksRecIn(0x0), @@ -159,6 +179,16 @@ AliAnalysisTaskJetSpectrum2::AliAnalysisTaskJetSpectrum2(const char* name): fh2TracksLeadingPhiEta(0x0), fh2TracksLeadingPhiPt(0x0), fh2TracksLeadingJetPhiPt(0x0), + fh2JetPtJetPhi(0x0), + fh2TrackPtTrackPhi(0x0), + fh2DijetDeltaPhiPt(0x0), + fh2DijetAsymPt(0x0), + fh2DijetAsymPtCut(0x0), + fh2DijetDeltaPhiDeltaEta(0x0), + fh2DijetPt2vsPt1(0x0), + fh2DijetDifvsSum(0x0), + fh1DijetMinv(0x0), + fh1DijetMinvCut(0x0), fHistList(0x0) { @@ -166,7 +196,7 @@ AliAnalysisTaskJetSpectrum2::AliAnalysisTaskJetSpectrum2(const char* name): fhnJetContainer[i] = 0; } for(int i = 0;i < kMaxJets;++i){ - fh2PhiPt[i] = fh2PhiEta[i] = fh2FragRec[i] = fh2FragLnRec[i] = fh2FragGen[i] = fh2FragLnGen[i] = 0; + fh2PhiPt[i] = fh2PhiEta[i] = fh2RhoPtRec[i] = fh2PsiPtRec[i] = fh2FragRec[i] = fh2PsiPtGen[i] = fh2RhoPtGen[i] = fh2FragGen[i] = fh2FragLnRec[i] = fh2FragGen[i] = fh2FragLnGen[i] = 0; fh1PtRecIn[i] = fh1PtGenIn[i] = 0; } DefineOutput(1, TList::Class()); @@ -230,7 +260,7 @@ void AliAnalysisTaskJetSpectrum2::UserCreateOutputObjects() // // Histogram - const Int_t nBinPt = 100; + const Int_t nBinPt = 240; Double_t binLimitsPt[nBinPt+1]; for(Int_t iPt = 0;iPt <= nBinPt;iPt++){ if(iPt == 0){ @@ -253,6 +283,18 @@ void AliAnalysisTaskJetSpectrum2::UserCreateOutputObjects() } + const Int_t nBinPhi2 = 360; + Double_t binLimitsPhi2[nBinPhi2+1]; + for(Int_t iPhi2 = 0;iPhi2<=nBinPhi2;iPhi2++){ + if(iPhi2==0){ + binLimitsPhi2[iPhi2] = 0.; + } + else{ + binLimitsPhi2[iPhi2] = binLimitsPhi2[iPhi2-1] + 1/(Float_t)nBinPhi2 * TMath::Pi()*2; + } + } + + const Int_t nBinEta = 40; Double_t binLimitsEta[nBinEta+1]; @@ -304,6 +346,10 @@ void AliAnalysisTaskJetSpectrum2::UserCreateOutputObjects() fh2TracksLeadingPhiPt = new TH2F("fh2TracksLeadingPhiPt","leading p_T vs delta phi to leading jet;#Delta#phi;p_{T} (GeV/c)", nBinPhi,binLimitsPhi,nBinPt,binLimitsPt); + fh2JetPtJetPhi = new TH2F("fh2JetPtJetPhi","Reconstructed jet phi vs. pt",nBinPt,binLimitsPt,nBinPhi2,binLimitsPhi2); + fh2TrackPtTrackPhi = new TH2F("fh2TrackPtTrackPhi","Reconstructed track phi vs. pt",nBinPt,binLimitsPt,nBinPhi2,binLimitsPhi2); + + for(int ij = 0;ij < kMaxJets;++ij){ fh1PtRecIn[ij] = new TH1F(Form("fh1PtRecIn_j%d",ij),"rec p_T input ;p_{T,rec}",nBinPt,binLimitsPt); fh1PtGenIn[ij] = new TH1F(Form("fh1PtGenIn_j%d",ij),"found p_T input ;p_{T,gen}",nBinPt,binLimitsPt); @@ -314,6 +360,17 @@ void AliAnalysisTaskJetSpectrum2::UserCreateOutputObjects() fh2PhiEta[ij] = new TH2F(Form("fh2PhiEtaRec_j%d",ij),"delta eta vs delta phi for jets;#Delta#phi;#Delta#eta", nBinPhi,binLimitsPhi,nBinEta,binLimitsEta); + fh2RhoPtRec[ij] = new TH2F(Form("fh2RhoPtRec_j%d",ij),"jet shape rho for jets;r;p_{T,rec};", + 20,0.,1.,nBinPt,binLimitsPt); + fh2PsiPtRec[ij] = new TH2F(Form("fh2PsiPtRec_j%d",ij),"jet shape psi for jets;r;p_{T,rec};", + 20,0.,1.,nBinPt,binLimitsPt); + + fh2RhoPtGen[ij] = new TH2F(Form("fh2RhoPtGen_j%d",ij),"jet shape rho for jets;r;p_{T,gen};", + 20,0.,1.,nBinPt,binLimitsPt); + fh2PsiPtGen[ij] = new TH2F(Form("fh2PsiPtGen_j%d",ij),"jet shape psi for jets;r;p_{T,gen};", + 20,0.,1.,nBinPt,binLimitsPt); + if(!fh1TmpRho)fh1TmpRho = new TH1F("fh1TmpRho","tmp histo for jet shape",20,0.,1); + fh2FragRec[ij] = new TH2F(Form("fh2FragRec_j%d",ij),"Jet Fragmentation;x=p_{T,i}/p_{T,jet};p_{T,jet};1/N_{jet}dN_{ch}/dx", nBinFrag,0.,1.,nBinPt,binLimitsPt); @@ -321,11 +378,24 @@ void AliAnalysisTaskJetSpectrum2::UserCreateOutputObjects() nBinFrag,0.,10.,nBinPt,binLimitsPt); fh2FragGen[ij] = new TH2F(Form("fh2FragGen_j%d",ij),"Jet Fragmentation;x=p_{T,i}/p_{T,jet};p_{T,jet};1/N_{jet}dN_{ch}/dx", - nBinFrag,0.,1.,nBinPt,binLimitsPt); + nBinFrag,0.,1.,nBinPt,binLimitsPt); fh2FragLnGen[ij] = new TH2F(Form("fh2FragLnGen_j%d",ij),"Jet Fragmentation Ln;#xi=ln(p_{T,jet}/p_{T,i});p_{T,jet}(GeV);1/N_{jet}dN_{ch}/d#xi", - nBinFrag,0.,10.,nBinPt,binLimitsPt); + nBinFrag,0.,10.,nBinPt,binLimitsPt); } + // Dijet histograms + + fh2DijetDeltaPhiPt = new TH2F("fh2DeltaPhiPt","Difference in the azimuthal angle;#Delta#phi;p_{T,1};Entries",180,0.,TMath::Pi(),nBinPt,binLimitsPt); + fh2DijetAsymPt = new TH2F("fh2DijetAsym","Pt asymmetry;#Deltap_{T}/(p_{T,1}+p_{T,2});p_{T,1};Entries",50,0.,1.,nBinPt,binLimitsPt); + fh2DijetAsymPtCut = new TH2F("fh2DijetAsymCut","Pt asymmetry after delta phi cut;#Deltap_{T}/(p_{T,1}+p_{T,2});p_{T,1};Entries",50,0.,1.,nBinPt,binLimitsPt); + fh2DijetDeltaPhiDeltaEta = new TH2F("fh2DijetDeltaPhiDeltaEta","Difference in the azimuthal angle;#Delta#phi;Entries",180,0.,TMath::Pi(),20,-2.,2.); + fh2DijetPt2vsPt1 = new TH2F("fh2DijetPt2vsPt1","Pt2 versus Pt1;p_{T,1} (GeV/c);p_{T,2} (GeV/c)",250,0.,250.,250,0.,250.); + fh2DijetDifvsSum = new TH2F("fh2DijetDifvsSum","Pt difference vs Pt sum;p_{T,1}+p_{T,2} (GeV/c);#Deltap_{T} (GeV/c)",400,0.,400.,150,0.,150.); + fh1DijetMinv = new TH1F("fh1DijetMinv","Dijet invariant mass;m_{JJ}",nBinPt,binLimitsPt); + fh1DijetMinvCut = new TH1F("fh1DijetMinvCut","Dijet invariant mass;m_{JJ}",nBinPt,binLimitsPt); + + + const Int_t saveLevel = 3; // large save level more histos if(saveLevel>0){ @@ -355,17 +425,35 @@ void AliAnalysisTaskJetSpectrum2::UserCreateOutputObjects() for(int i = 0;iAdd(fhnJetContainer[i]); for(int ij = 0;ijAdd( fh1PtRecIn[ij]); + if(fBranchGen.Length()>0){ fHistList->Add( fh1PtGenIn[ij]); fHistList->Add( fh2FragGen[ij]); fHistList->Add( fh2FragLnGen[ij]); + fHistList->Add(fh2RhoPtGen[ij]); + fHistList->Add(fh2PsiPtGen[ij]); + fHistList->Add( fh2FragGen[ij]); } fHistList->Add( fh2PhiPt[ij]); fHistList->Add( fh2PhiEta[ij]); + fHistList->Add(fh2RhoPtRec[ij]); + fHistList->Add(fh2PsiPtRec[ij]); fHistList->Add( fh2FragRec[ij]); fHistList->Add( fh2FragLnRec[ij]); } fHistList->Add(fhnCorrelation); + fHistList->Add(fhnCorrelationPhiZRec); + fHistList->Add(fh2JetPtJetPhi); + fHistList->Add(fh2TrackPtTrackPhi); + + fHistList->Add(fh2DijetDeltaPhiPt); + fHistList->Add(fh2DijetAsymPt); + fHistList->Add(fh2DijetAsymPtCut); + fHistList->Add(fh2DijetDeltaPhiDeltaEta); + fHistList->Add(fh2DijetPt2vsPt1); + fHistList->Add(fh2DijetDifvsSum); + fHistList->Add(fh1DijetMinv); + fHistList->Add(fh1DijetMinvCut); } // =========== Switch on Sumw2 for all histos =========== @@ -387,7 +475,6 @@ void AliAnalysisTaskJetSpectrum2::Init() // Initialization // - Printf(">>> AnalysisTaskJetSpectrum2::Init() debug level %d\n",fDebug); if (fDebug > 1) printf("AnalysisTaskJetSpectrum2::Init() \n"); } @@ -500,19 +587,19 @@ void AliAnalysisTaskJetSpectrum2::UserExec(Option_t */*option*/) pythiaGenHeader->TriggerJet(ip,p); pythiaGenJets[iCount].SetPxPyPzE(p[0],p[1],p[2],p[3]); - /* - if(fLimitGenJetEta){ - if(pythiaGenJets[iCount].Eta()>fJetHeaderRec->GetJetEtaMax()|| - pythiaGenJets[iCount].Eta()GetJetEtaMin())continue; - } - */ - if(fBranchGen.Length()==0){ + /* + if(fLimitGenJetEta){ + if(pythiaGenJets[iCount].Eta()>fJetHeaderRec->GetJetEtaMax()|| + pythiaGenJets[iCount].Eta()GetJetEtaMin())continue; + } + */ + if(fMinJetPt>0&&pythiaGenJets[iCount].Pt()Eta()GetJetEtaMin())continue; } */ + if(fMinJetPt>0&&tmp->Pt()At(ir)); Float_t tmpPt = tmp->Pt(); if(tmpPt>ptOld){ - Printf("%s:%d Jets Not Sorted!! %d:%.3E %d%.3E",(char*)__FILE__,__LINE__,ir,tmpPt,ir-1,ptOld); + Printf("%s:%d Jets Not Sorted %s !! %d:%.3E %d:%.3E",(char*)__FILE__,__LINE__,fBranchRec.Data(),ir,tmpPt,ir-1,ptOld); } ptOld = tmpPt; } @@ -694,14 +782,17 @@ void AliAnalysisTaskJetSpectrum2::UserExec(Option_t */*option*/) delete genIter; } - fh1NRecJets->Fill(nRecJets); nRecJets = TMath::Min(nRecJets,kMaxJets); - + + Int_t iCountRec = 0; for(int ir = 0;ir < nRecJets;++ir){ AliAODJet *tmp = dynamic_cast(aodRecJets->At(ir)); if(!tmp)continue; + if(tmp->Pt() 10)Printf("%s:%d",(char*)__FILE__,__LINE__); @@ -729,6 +820,7 @@ void AliAnalysisTaskJetSpectrum2::UserExec(Option_t */*option*/) Double_t container[6]; + Double_t containerPhiZ[6]; // loop over generated jets @@ -750,30 +842,38 @@ void AliAnalysisTaskJetSpectrum2::UserExec(Option_t */*option*/) Int_t ir = iRecIndex[ig]; if(TMath::Abs(etaGen)Reset(); + fhnJetContainer[kStep1]->Fill(&container[3],eventW); fh1PtGenIn[ig]->Fill(ptGen,eventW); // fill the fragmentation function for(int it = 0;itFill(deltaR,part->Pt()/ptGen); + if(deltaRPt()/ptGen; Float_t lnz = -1.*TMath::Log(z); fh2FragGen[ig]->Fill(z,ptGen,eventW); fh2FragLnGen[ig]->Fill(lnz,ptGen,eventW); } + } - if(ir>=0&&irFill(&container[3],eventW); + Float_t rhoSum = 0; + for(int ibx = 1;ibxGetNbinsX();ibx++){ + Float_t r = fh2RhoPtGen[ir]->GetXaxis()->GetBinCenter(ibx); + Float_t rho = fh1TmpRho->GetBinContent(ibx); + rhoSum += rho; + fh2RhoPtGen[ig]->Fill(r,ptGen,rho); + fh2PsiPtGen[ig]->Fill(r,ptGen,rhoSum); } - } - + } if(ir>=0&&irFill(&container[3],eventW); Double_t etaRec = recJets[ir].Eta(); if(TMath::Abs(etaRec)Fill(&container[3],eventW); } }// loop over generated jets - Float_t sumPt = 0; for(int it = 0;itEta())<0.9){ Float_t pt = part->Pt(); fh1PtTrackRec->Fill(pt,eventW); + fh2TrackPtTrackPhi->Fill(pt,part->Phi()); sumPt += pt; } } @@ -795,14 +896,56 @@ void AliAnalysisTaskJetSpectrum2::UserExec(Option_t */*option*/) Double_t ptRec = recJets[ir].Pt(); Double_t phiRec = recJets[ir].Phi(); if(phiRec<0)phiRec+=TMath::Pi()*2.; + + + // do something with dijets... + if(ir==1){ + Double_t etaRec1 = recJets[0].Eta(); + Double_t ptRec1 = recJets[0].Pt(); + Double_t phiRec1 = recJets[0].Phi(); + if(phiRec1<0)phiRec1+=TMath::Pi()*2.; + + + if(TMath::Abs(etaRec1)TMath::Pi())deltaPhi = deltaPhi - 2.*TMath::Pi(); + if(deltaPhi<(-1.*TMath::Pi()))deltaPhi = deltaPhi + 2.*TMath::Pi(); + deltaPhi = TMath::Abs(deltaPhi); + fh2DijetDeltaPhiPt->Fill(deltaPhi,ptRec1); + Float_t asym = (ptRec1-ptRec)/(ptRec1+ptRec); + fh2DijetAsymPt->Fill(asym,ptRec1); + fh2DijetDeltaPhiDeltaEta->Fill(deltaPhi,etaRec1-etaRec); + fh2DijetPt2vsPt1->Fill(ptRec1,ptRec); + fh2DijetDifvsSum->Fill(ptRec1+ptRec,ptRec1-ptRec); + Float_t minv = 2.*(recJets[0].P()*recJets[1].P()- + recJets[0].Px()*recJets[1].Px()- + recJets[0].Py()*recJets[1].Py()- + recJets[0].Pz()*recJets[1].Py()); + minv = TMath::Sqrt(minv); + // with mass == 0; + + fh1DijetMinv->Fill(minv); + if((TMath::Pi()-deltaPhi)Fill(minv); + fh2DijetAsymPtCut->Fill(asym,ptRec1); + } + } + } + + container[0] = ptRec; container[1] = etaRec; container[2] = phiRec; - - if(ptRec>10.&&fDebug>0){ + containerPhiZ[0] = ptRec; + containerPhiZ[1] = phiRec; + if(ptRec>30.&&fDebug>0){ // need to cast to int, otherwise the printf overwrites Printf("Jet found in Event %d with p_T, %E",(int)Entry(),ptRec); - Printf("%s read event, %d",fInputHandler->GetTree()->GetCurrentFile()->GetName(),fInputHandler->GetTree()->GetReadEntry()); + Printf("%s read event, %d",fInputHandler->GetTree()->GetCurrentFile()->GetName(),(int)fInputHandler->GetTree()->GetTree()->GetReadEntry()); + if(fESD)Printf("ESDEvent GetEventNumberInFile(): %d",fESD->GetEventNumberInFile()); // aodH->SetFillAOD(kTRUE); fAOD->GetHeader()->Print(); Printf("TriggerClasses: %s",fAOD->GetFiredTriggerClasses().Data()); @@ -822,10 +965,16 @@ void AliAnalysisTaskJetSpectrum2::UserExec(Option_t */*option*/) fhnJetContainer[kStep0+kMaxStep]->Fill(container,eventW); if (fDebug > 10)Printf("%s:%d",(char*)__FILE__,__LINE__); + + Float_t zLeading = -1; if(TMath::Abs(etaRec)Fill(ptRec,phiRec); fhnJetContainer[kStep1+kMaxStep]->Fill(container,eventW); fh1PtRecIn[ir]->Fill(ptRec,eventW); // fill the fragmentation function + + fh1TmpRho->Reset(); + for(int it = 0;itEta(); @@ -834,20 +983,38 @@ void AliAnalysisTaskJetSpectrum2::UserExec(Option_t */*option*/) if(phi<0)phi+=TMath::Pi()*2.; Float_t dPhi = phi - phiRec; Float_t dEta = eta - etaRec; - if(dPhi<(-1.*TMath::Pi()))phiRec+=TMath::Pi()*2.; + if(dPhi>TMath::Pi())dPhi = dPhi - 2.*TMath::Pi(); + if(dPhi<(-1.*TMath::Pi()))dPhi = dPhi + 2.*TMath::Pi(); fh2PhiPt[ir]->Fill(dPhi,ptRec,eventW); fh2PhiEta[ir]->Fill(dPhi,dEta,eventW); } - if(recJets[ir].DeltaR(part)Fill(deltaR,part->Pt()/ptRec); + + + if(deltaRPt()/ptRec; + if(z>zLeading)zLeading=z; Float_t lnz = -1.*TMath::Log(z); fh2FragRec[ir]->Fill(z,ptRec,eventW); fh2FragLnRec[ir]->Fill(lnz,ptRec,eventW); } } + // fill the jet shapes + Float_t rhoSum = 0; + for(int ibx = 1;ibxGetNbinsX();ibx++){ + Float_t r = fh2RhoPtRec[ir]->GetXaxis()->GetBinCenter(ibx); + Float_t rho = fh1TmpRho->GetBinContent(ibx); + rhoSum += rho; + fh2RhoPtRec[ir]->Fill(r,ptRec,rho); + fh2PsiPtRec[ir]->Fill(r,ptRec,rhoSum); + } } + containerPhiZ[2] = zLeading; + // Fill Correlation Int_t ig = iGenIndex[ir]; if(ig>=0 && igFill(container,eventW); fhnCorrelation->Fill(container); + if(fhnCorrelationPhiZRec)fhnCorrelationPhiZRec->Fill(containerPhiZ); + }// if etarec in window } - //////////////////////////////////////////////////// else{ - + containerPhiZ[3] = 0; + if(fhnCorrelationPhiZRec)fhnCorrelationPhiZRec->Fill(containerPhiZ); } }// loop over reconstructed jets @@ -892,19 +1061,19 @@ void AliAnalysisTaskJetSpectrum2::MakeJetContainer(){ // link it // const Int_t kNvar = 3 ; //number of variables on the grid:pt,eta, phi - const Double_t kPtmin = 5.0, kPtmax = 105.; // we do not want to have empty bins at the beginning... + const Double_t kPtmin = 0.0, kPtmax = 160.; // we do not want to have empty bins at the beginning... const Double_t kEtamin = -3.0, kEtamax = 3.0; const Double_t kPhimin = 0., kPhimax = 2. * TMath::Pi(); + const Double_t kZmin = 0., kZmax = 1; // can we neglect migration in eta and phi? // phi should be no problem since we cover full phi and are phi symmetric // eta migration is more difficult due to needed acceptance correction // in limited eta range - //arrays for the number of bins in each dimension Int_t iBin[kNvar]; - iBin[0] = 100; //bins in pt + iBin[0] = 160; //bins in pt iBin[1] = 1; //bins in eta iBin[2] = 1; // bins in phi @@ -942,14 +1111,42 @@ void AliAnalysisTaskJetSpectrum2::MakeJetContainer(){ } fhnCorrelation->Sumw2(); + // for second correlation histogram + + + const Int_t kNvarPhiZ = 4; + //arrays for the number of bins in each dimension + Int_t iBinPhiZ[kNvarPhiZ]; + iBinPhiZ[0] = 80; //bins in pt + iBinPhiZ[1] = 72; //bins in phi + iBinPhiZ[2] = 20; // bins in Z + iBinPhiZ[3] = 80; //bins in ptgen + + //arrays for lower bounds : + Double_t* binEdgesPhiZ[kNvarPhiZ]; + for(Int_t ivar = 0; ivar < kNvarPhiZ; ivar++) + binEdgesPhiZ[ivar] = new Double_t[iBinPhiZ[ivar] + 1]; + + for(Int_t i=0; i<=iBinPhiZ[0]; i++) binEdgesPhiZ[0][i]=(Double_t)kPtmin + (kPtmax-kPtmin)/(Double_t)iBinPhiZ[0]*(Double_t)i; + for(Int_t i=0; i<=iBinPhiZ[1]; i++) binEdgesPhiZ[1][i]=(Double_t)kPhimin + (kPhimax-kPhimin)/iBinPhiZ[1]*(Double_t)i; + for(Int_t i=0; i<=iBinPhiZ[2]; i++) binEdgesPhiZ[2][i]=(Double_t)kZmin + (kZmax-kZmin)/iBinPhiZ[2]*(Double_t)i; + for(Int_t i=0; i<=iBinPhiZ[3]; i++) binEdgesPhiZ[3][i]=(Double_t)kPtmin + (kPtmax-kPtmin)/(Double_t)iBinPhiZ[3]*(Double_t)i; + + fhnCorrelationPhiZRec = new THnSparseF("fhnCorrelationPhiZRec","THnSparse with correlations",kNvarPhiZ,iBinPhiZ); + for (int k=0; kSetBinEdges(k,binEdgesPhiZ[k]); + } + fhnCorrelationPhiZRec->Sumw2(); + + // Add a histogram for Fake jets - // thnDim[3] = AliPID::kSPECIES; - // fFakeElectrons = new THnSparseF("fakeEkectrons", "Output for Fake Electrons", kNvar + 1, thnDim); - // for(Int_t idim = 0; idim < kNvar; idim++) - // fFakeElectrons->SetBinEdges(idim, binEdges[idim]); + for(Int_t ivar = 0; ivar < kNvar; ivar++) delete [] binEdges[ivar]; + for(Int_t ivar = 0; ivar < kNvarPhiZ; ivar++) + delete [] binEdgesPhiZ[ivar]; + } void AliAnalysisTaskJetSpectrum2::Terminate(Option_t */*option*/)