//////////////////////////////////////////////////////////////////////////////// // // // AliFemtoCorrFctnDEtaDPhiCorrections - A correlation function that analyzes // // two particle correlations with respect to the azimuthal angle (phi) // // and pseudorapidity (eta) difference // // // // Authors: Adam Kisiel Adam.Kisiel@cern.ch // // // //////////////////////////////////////////////////////////////////////////////// #include "AliFemtoCorrFctnDEtaDPhiCorrections.h" #include "AliFemtoModelHiddenInfo.h" //#include "AliFemtoHisto.hh" #include #include #include "THn.h" #ifdef __ROOT__ ClassImp(AliFemtoCorrFctnDEtaDPhiCorrections) #endif #define PIH 1.57079632679489656 #define PIT 6.28318530717958623 //____________________________ AliFemtoCorrFctnDEtaDPhiCorrections::AliFemtoCorrFctnDEtaDPhiCorrections(char* title, const int& aPhiBins=20, const int& aEtaBins=20): AliFemtoCorrFctn(), fDPhiDEtaNumerator(0), fDPhiDEtaDenominator(0), fDPhiNumerator(0), fDPhiDenominator(0), fDCosNumerator(0), fDCosDenominator(0), fDoFullAnalysis(kFALSE), fPhi(0), fEta(0), fPtSumDist(0), fYtYtNumerator(0), fYtYtDenominator(0), fIfCorrection(0), fPtCorrectionsNum(0), fPtCorrectionsDen(0), fEtaCorrectionsNum(0), fEtaCorrectionsDen(0), fCorrFactorTab(0), fpTab(0), fPartType(kNoCorrection), fphiL(0), fphiT(0), ifileCorrTab(0), fdoPtCorr(0), fdoEtaCorr(0), fdoPhiCorr(0), fdoZVertCorr(0), fpartType1(0), fpartType2(0), fhntReco1(0), fhntReco2(0), fh1Reco1(0), fh1Reco2(0), fh2Reco1(0), fh2Reco2(0), fh3Reco1(0), fh3Reco2(0), fhCont1(0), fhCont2(0), fCorr1D(kFALSE) { fphiL = (-(int)(aPhiBins/4)+0.5)*2.*TMath::Pi()/aPhiBins; fphiT = 2*TMath::Pi()+(-(int)(aPhiBins/4)+0.5)*2.*TMath::Pi()/aPhiBins; // set up numerator char tTitNumD[101] = "NumDPhiDEta"; strncat(tTitNumD,title, 100); fDPhiDEtaNumerator = new TH2D(tTitNumD,title,aPhiBins,fphiL,fphiT,aEtaBins,-2.0,2.0); // set up denominator char tTitDenD[101] = "DenDPhiDEta"; strncat(tTitDenD,title, 100); fDPhiDEtaDenominator = new TH2D(tTitDenD,title,aPhiBins,fphiL,fphiT,aEtaBins,-2.0,2.0); // set up numerator char tTitNum[101] = "PtSumDist"; strncat(tTitNum,title, 100); fPtSumDist = new TH1D(tTitNum,title,200,0,10); fDPhiDEtaNumerator->Sumw2(); fDPhiDEtaDenominator->Sumw2(); fPtSumDist->Sumw2(); // to enable error bar calculation... } //____________________________ AliFemtoCorrFctnDEtaDPhiCorrections::AliFemtoCorrFctnDEtaDPhiCorrections(const AliFemtoCorrFctnDEtaDPhiCorrections& aCorrFctn) : AliFemtoCorrFctn(), fDPhiDEtaNumerator(0), fDPhiDEtaDenominator(0), fDPhiNumerator(0), fDPhiDenominator(0), fDCosNumerator(0), fDCosDenominator(0), fDoFullAnalysis(kFALSE), fPhi(0), fEta(0), fPtSumDist(0), fYtYtNumerator(0), fYtYtDenominator(0), fIfCorrection(0), fPtCorrectionsNum(0), fPtCorrectionsDen(0), fEtaCorrectionsNum(0), fEtaCorrectionsDen(0), fCorrFactorTab(0), fpTab(0), fPartType(kNoCorrection), fphiL(0), fphiT(0), ifileCorrTab(0), fdoPtCorr(0), fdoEtaCorr(0), fdoPhiCorr(0), fdoZVertCorr(0), fpartType1(0), fpartType2(0), fhntReco1(0), fhntReco2(0), fh1Reco1(0), fh1Reco2(0), fh2Reco1(0), fh2Reco2(0), fh3Reco1(0), fh3Reco2(0), fhCont1(0), fhCont2(0), fCorr1D(kFALSE) { // copy constructor if (aCorrFctn.fDPhiDEtaNumerator) fDPhiDEtaNumerator = new TH2D(*aCorrFctn.fDPhiDEtaNumerator); else fDPhiDEtaNumerator = 0; if (aCorrFctn.fDPhiDEtaDenominator) fDPhiDEtaDenominator = new TH2D(*aCorrFctn.fDPhiDEtaDenominator); else fDPhiDEtaDenominator = 0; if (aCorrFctn.fDPhiNumerator) fDPhiNumerator = new TH1D(*aCorrFctn.fDPhiNumerator); else fDPhiNumerator = 0; if (aCorrFctn.fDPhiDenominator) fDPhiDenominator = new TH1D(*aCorrFctn.fDPhiDenominator); else fDPhiDenominator = 0; if (aCorrFctn.fDCosNumerator) fDCosNumerator = new TH1D(*aCorrFctn.fDCosNumerator); else fDCosNumerator = 0; if (aCorrFctn.fDCosDenominator) fDCosDenominator = new TH1D(*aCorrFctn.fDCosDenominator); else fDCosDenominator = 0; if (aCorrFctn.fPhi) fPhi = new TH1D(*aCorrFctn.fPhi); else fPhi = 0; if (aCorrFctn.fEta) fEta = new TH1D(*aCorrFctn.fEta); else fEta = 0; if (aCorrFctn.fPtSumDist) fPtSumDist = new TH1D(*aCorrFctn.fPtSumDist); else fPtSumDist = 0; if (aCorrFctn.fYtYtNumerator) fYtYtNumerator = new TH2D(*aCorrFctn.fDPhiDEtaDenominator); else fYtYtNumerator = 0; if (aCorrFctn.fYtYtDenominator) fYtYtDenominator = new TH2D(*aCorrFctn.fDPhiDEtaDenominator); else fYtYtDenominator = 0; fphiL = aCorrFctn.fphiL; fphiT = aCorrFctn.fphiT; fPartType = aCorrFctn.fPartType; // if (aCorrFctn.fPtCorrectionsNum) // fPtCorrectionsNum = new THnSparseF(*aCorrFctn.fPtCorrectionsNum); // else // fPtCorrectionsNum = 0; // if (aCorrFctn.fPtCorrectionsDen) // fPtCorrectionsDen = new THnSparseF(*aCorrFctn.fPtCorrectionsDen); // else // fPtCorrectionsDen = 0; } //____________________________ AliFemtoCorrFctnDEtaDPhiCorrections::~AliFemtoCorrFctnDEtaDPhiCorrections(){ // destructor delete fDPhiDEtaNumerator; delete fDPhiDEtaDenominator; delete fPtSumDist; if (fDoFullAnalysis) { delete fDPhiNumerator; delete fDPhiDenominator; delete fDCosNumerator; delete fDCosDenominator; delete fYtYtNumerator; delete fYtYtDenominator; delete fPhi; delete fEta; delete fPtCorrectionsNum; delete fPtCorrectionsDen; delete fEtaCorrectionsNum; delete fEtaCorrectionsDen; } //corrctions if(fhntReco1){ delete fhntReco1; delete fhntReco2; delete fhCont1; delete fhCont2; } if(fh1Reco1) delete fh1Reco1; if(fh1Reco2) delete fh1Reco2; if(fh2Reco1) delete fh2Reco1; if(fh2Reco2) delete fh2Reco2; if(fh3Reco1) delete fh3Reco1; if(fh3Reco2) delete fh3Reco2; } //_________________________ AliFemtoCorrFctnDEtaDPhiCorrections& AliFemtoCorrFctnDEtaDPhiCorrections::operator=(const AliFemtoCorrFctnDEtaDPhiCorrections& aCorrFctn) { // assignment operator if (this == &aCorrFctn) return *this; if (aCorrFctn.fDPhiDEtaNumerator) fDPhiDEtaNumerator = new TH2D(*aCorrFctn.fDPhiDEtaNumerator); else fDPhiDEtaNumerator = 0; if (aCorrFctn.fDPhiDEtaDenominator) fDPhiDEtaDenominator = new TH2D(*aCorrFctn.fDPhiDEtaDenominator); else fDPhiDEtaDenominator = 0; if (aCorrFctn.fDPhiNumerator) fDPhiNumerator = new TH1D(*aCorrFctn.fDPhiNumerator); else fDPhiNumerator = 0; if (aCorrFctn.fDPhiDenominator) fDPhiDenominator = new TH1D(*aCorrFctn.fDPhiDenominator); else fDPhiDenominator = 0; if (aCorrFctn.fDCosNumerator) fDCosNumerator = new TH1D(*aCorrFctn.fDCosNumerator); else fDCosNumerator = 0; if (aCorrFctn.fDCosDenominator) fDCosDenominator = new TH1D(*aCorrFctn.fDCosDenominator); else fDCosDenominator = 0; if (aCorrFctn.fPhi) fPhi = new TH1D(*aCorrFctn.fPhi); else fPhi = 0; if (aCorrFctn.fEta) fEta = new TH1D(*aCorrFctn.fEta); else fEta = 0; if (aCorrFctn.fPtSumDist) fPtSumDist = new TH1D(*aCorrFctn.fPtSumDist); else fPtSumDist = 0; if (aCorrFctn.fYtYtNumerator) fYtYtNumerator = new TH2D(*aCorrFctn.fDPhiDEtaDenominator); else fYtYtNumerator = 0; if (aCorrFctn.fYtYtDenominator) fYtYtDenominator = new TH2D(*aCorrFctn.fDPhiDEtaDenominator); else fYtYtDenominator = 0; fphiL = aCorrFctn.fphiL; fphiT = aCorrFctn.fphiT; fPartType = aCorrFctn.fPartType; // if (aCorrFctn.fPtCorrectionsNum) // fPtCorrectionsNum = new THnSparseF(*aCorrFctn.fPtCorrectionsNum); // else // fPtCorrectionsNum = 0; // if (aCorrFctn.fPtCorrectionsDen) // fPtCorrectionsDen = new THnSparseF(*aCorrFctn.fPtCorrectionsDen); // else // fPtCorrectionsDen = 0; return *this; } //_________________________ void AliFemtoCorrFctnDEtaDPhiCorrections::Finish(){ // here is where we should normalize, fit, etc... // we should NOT Draw() the histos (as I had done it below), // since we want to insulate ourselves from root at this level // of the code. Do it instead at root command line with browser. // mShareNumerator->Draw(); //mShareDenominator->Draw(); //mRatio->Draw(); } //____________________________ AliFemtoString AliFemtoCorrFctnDEtaDPhiCorrections::Report(){ // create report string stemp = "TPC Ncls Correlation Function Report:\n"; char ctemp[100]; snprintf(ctemp , 100, "Number of entries in numerator:\t%E\n",fDPhiDEtaNumerator->GetEntries()); stemp += ctemp; snprintf(ctemp , 100, "Number of entries in denominator:\t%E\n",fDPhiDEtaDenominator->GetEntries()); stemp += ctemp; // stemp += mCoulombWeight->Report(); AliFemtoString returnThis = stemp; return returnThis; } //____________________________ void AliFemtoCorrFctnDEtaDPhiCorrections::AddRealPair( AliFemtoPair* pair){ // add real (effect) pair if (fPairCut) if (!fPairCut->Pass(pair)) return; /*double phi1 = pair->Track1()->Track()->P().Phi(); double phi2 = pair->Track2()->Track()->P().Phi(); double eta1 = pair->Track1()->Track()->P().PseudoRapidity(); double eta2 = pair->Track2()->Track()->P().PseudoRapidity();*/ double phi1 = pair->Track1()->FourMomentum().Phi(); double phi2 = pair->Track2()->FourMomentum().Phi(); double eta1 = pair->Track1()->FourMomentum().PseudoRapidity(); double eta2 = pair->Track2()->FourMomentum().PseudoRapidity(); double dphi = phi1 - phi2; while (dphifphiT) dphi-=PIT; double deta = eta1 - eta2; double px1 = pair->Track1()->Track()->P().x(); double py1 = pair->Track1()->Track()->P().y(); //double pz1 = pair->Track1()->Track()->P().z(); double px2 = pair->Track2()->Track()->P().x(); double py2 = pair->Track2()->Track()->P().y(); //double pz2 = pair->Track2()->Track()->P().z(); double pt1 = TMath::Hypot(px1, py1); double pt2 = TMath::Hypot(px2, py2); double vert1[3]; pair->Track1()->Track()->GetPrimaryVertex(vert1); double vert2[3]; pair->Track2()->Track()->GetPrimaryVertex(vert2); double corrweight=0; //double corrweightpT1=1; double corrweightpT2=1; //if (fIfCorrection) corrweight = CalculateCorrectionWeight(pt1, pt2); if (fIfCorrection) { corrweight = CalculateCorrectionWeight(pt1, pt2, eta1, eta2, phi1, phi2, vert1[2], vert2[2]); } else if(fCorr1D) { corrweight = CalculateCorrectionWeight(pt1, pt2); //corrweightpT1 = CalculateCorrectionWeight(pt1); //corrweightpT2 = CalculateCorrectionWeight(pt2); } fPtSumDist->Fill(pt1+pt2,corrweight); /* double ptmin = pt1>pt2 ? pt2 : pt1; double cosphi = (px1*px2 + py1*py2 + pz1*pz2)/ sqrt((px1*px1 + py1*py1 + pz1*pz1)*(px2*px2 + py2*py2 + pz2*pz2)); */ if (fIfCorrection || fCorr1D) { fDPhiDEtaNumerator->Fill(dphi, deta, corrweight); } else{ fDPhiDEtaNumerator->Fill(dphi, deta); } if (fDoFullAnalysis) { //fDPhiPtNumerator->Fill(dphi, ptmin); //fDCosPtNumerator->Fill(cosphi, ptmin); fDPhiNumerator->Fill(dphi); //fDCosNumerator->Fill(cosphi); double PionMass = 0.13956995; double yt1 = TMath::Log(sqrt(1+(pt1/PionMass)*(pt1/PionMass))+(pt1/PionMass)); double yt2 = TMath::Log(sqrt(1+(pt2/PionMass)*(pt2/PionMass))+(pt2/PionMass)); fYtYtNumerator->Fill(yt1,yt2); fPhi->Fill(phi1); fEta->Fill(eta1); } } //____________________________ void AliFemtoCorrFctnDEtaDPhiCorrections::AddMixedPair( AliFemtoPair* pair){ // add mixed (background) pair if (fPairCut) if (!fPairCut->Pass(pair)) return; /*double phi1 = pair->Track1()->Track()->P().Phi(); double phi2 = pair->Track2()->Track()->P().Phi(); double eta1 = pair->Track1()->Track()->P().PseudoRapidity(); double eta2 = pair->Track2()->Track()->P().PseudoRapidity();*/ double phi1 = pair->Track1()->FourMomentum().Phi(); double phi2 = pair->Track2()->FourMomentum().Phi(); double eta1 = pair->Track1()->FourMomentum().PseudoRapidity(); double eta2 = pair->Track2()->FourMomentum().PseudoRapidity(); double dphi = phi1 - phi2; while (dphifphiT) dphi-=PIT; double deta = eta1 - eta2; double px1 = pair->Track1()->Track()->P().x(); double py1 = pair->Track1()->Track()->P().y(); //double pz1 = pair->Track1()->Track()->P().z(); double px2 = pair->Track2()->Track()->P().x(); double py2 = pair->Track2()->Track()->P().y(); //double pz2 = pair->Track2()->Track()->P().z(); double pt1 = TMath::Hypot(px1, py1); double pt2 = TMath::Hypot(px2, py2); // double ptmin = pt1>pt2 ? pt2 : pt1; // double cosphi = (px1*px2 + py1*py2 + pz1*pz2)/ // sqrt((px1*px1 + py1*py1 + pz1*pz1)*(px2*px2 + py2*py2 + pz2*pz2)); double vert1[3]; pair->Track1()->Track()->GetPrimaryVertex(vert1); double vert2[3]; pair->Track2()->Track()->GetPrimaryVertex(vert2); double corrweight=-999; //if (fIfCorrection) corrweight = CalculateCorrectionWeight(pt1, pt2); if (fIfCorrection) { corrweight = CalculateCorrectionWeight(pt1, pt2, eta1, eta2, phi1, phi2, vert1[2], vert2[2]); } else if(fCorr1D) { corrweight = CalculateCorrectionWeight(pt1, pt2); } if(fIfCorrection || fCorr1D) fDPhiDEtaDenominator->Fill(dphi, deta, corrweight); else fDPhiDEtaDenominator->Fill(dphi, deta); if (fDoFullAnalysis) { //fDPhiPtDenominator->Fill(dphi, ptmin); //fDCosPtDenominator->Fill(cosphi, ptmin); fDPhiDenominator->Fill(dphi); double PionMass = 0.13956995; double yt1 = TMath::Log(sqrt(1+(pt1/PionMass)*(pt1/PionMass))+(pt1/PionMass)); double yt2 = TMath::Log(sqrt(1+(pt2/PionMass)*(pt2/PionMass))+(pt2/PionMass)); fYtYtDenominator->Fill(yt1,yt2); } } void AliFemtoCorrFctnDEtaDPhiCorrections::WriteHistos() { // Write out result histograms fDPhiDEtaNumerator->Write(); fDPhiDEtaDenominator->Write(); fPtSumDist->Write(); /*fDPhiNumerator->Write(); fDPhiDenominator->Write(); fDCosNumerator->Write(); fDCosDenominator->Write(); if (fDoFullAnalysis) { fDPhiPtNumerator->Write(); fDPhiPtDenominator->Write(); fDCosPtNumerator->Write(); fDCosPtDenominator->Write(); }*/ // fPhi->Write(); // fEta->Write(); } TList* AliFemtoCorrFctnDEtaDPhiCorrections::GetOutputList() { // Prepare the list of objects to be written to the output TList *tOutputList = new TList(); tOutputList->Add(fDPhiDEtaNumerator); tOutputList->Add(fDPhiDEtaDenominator); tOutputList->Add(fPtSumDist); if (fDoFullAnalysis) { // tOutputList->Add(fDPhiPtNumerator); // tOutputList->Add(fDPhiPtDenominator); //tOutputList->Add(fDCosPtNumerator); //tOutputList->Add(fDCosPtDenominator); tOutputList->Add(fDPhiNumerator); tOutputList->Add(fDPhiDenominator); tOutputList->Add(fDCosNumerator); tOutputList->Add(fDCosDenominator); tOutputList->Add(fYtYtNumerator); tOutputList->Add(fYtYtDenominator); tOutputList->Add(fPhi); tOutputList->Add(fEta); } return tOutputList; } void AliFemtoCorrFctnDEtaDPhiCorrections::SetDoFullAnalysis(Bool_t do2d) { fDoFullAnalysis = do2d; if(fDoFullAnalysis) { int aPhiBins = fDPhiDEtaNumerator->GetNbinsX(); int aEtaBins = fDPhiDEtaNumerator->GetNbinsY(); const char *title = fDPhiDEtaNumerator->GetTitle(); // set up numerator char tTitNumDPhi[101] = "NumDPhi"; strncat(tTitNumDPhi,title, 100); fDPhiNumerator = new TH1D(tTitNumDPhi,title,aPhiBins*2,fphiL, fphiT); // set up denominator char tTitDenDPhi[101] = "DenDPhi"; strncat(tTitDenDPhi,title, 100); fDPhiDenominator = new TH1D(tTitDenDPhi,title,aPhiBins*2,fphiL, fphiT); // set up numerator char tTitNumDCos[101] = "NumDCos"; strncat(tTitNumDCos,title, 100); fDCosNumerator = new TH1D(tTitNumDCos,title,aPhiBins*2,-1.0,1.0); // set up denominator char tTitDenDCos[101] = "DenDCos"; strncat(tTitDenDCos,title, 100); fDCosDenominator = new TH1D(tTitDenDCos,title,aPhiBins*2,-1.0,1.0); // set up numerator char tTitYtNum[101] = "NumYtYt"; strncat(tTitYtNum,title, 100); fYtYtNumerator = new TH2D(tTitYtNum,title,aPhiBins,1,5,aEtaBins,1,5); // set up denominator char tTitYtYtDen[101] = "DenYtYt"; strncat(tTitYtYtDen,title, 100); fYtYtDenominator = new TH2D(tTitYtYtDen,title,aPhiBins,1,5,aEtaBins,1,5); char tTitPhi[101] = "Phi"; strncat(tTitPhi,title, 100); fPhi = new TH1D(tTitPhi,title,90,-TMath::Pi(),TMath::Pi()); char tTitEta[101] = "Eta"; strncat(tTitEta,title, 100); fEta = new TH1D(tTitEta,title,90,-1.2,1.2); char tTitPtCorrectionsNum[101] = "NumpT1pT2EtaPhi"; strncat(tTitPtCorrectionsNum,title, 100); char tTitPtCorrectionsDen[101] = "DenpT1pT2EtaPhi"; strncat(tTitPtCorrectionsDen,title, 100); Int_t nbins[4] = {20,20,aPhiBins,aEtaBins}; Double_t xmin[4] = {0,0,-0.5*TMath::Pi(),-2.0}; Double_t xmax[4] = {4,4,1.5*TMath::Pi(),2.0}; fPtCorrectionsNum = new THnSparseF(tTitPtCorrectionsNum,title,4,nbins,xmin,xmax); fPtCorrectionsDen = new THnSparseF(tTitPtCorrectionsDen,title,4,nbins,xmin,xmax); char tTitEtaCorrectionsNum[101] = "NumEta1Eta2EtaPhi"; strncat(tTitEtaCorrectionsNum,title, 100); char tTitEtaCorrectionsDen[101] = "DenEta1Eta2EtaPhi"; strncat(tTitEtaCorrectionsDen,title, 100); Double_t xmineta[4] = {-1,1,-0.5*TMath::Pi(),-2.0}; Double_t xmaxeta[4] = {-1,1,1.5*TMath::Pi(),2.0}; fEtaCorrectionsNum = new THnSparseF(tTitEtaCorrectionsNum,title,4,nbins,xmineta,xmaxeta); fEtaCorrectionsDen = new THnSparseF(tTitEtaCorrectionsDen,title,4,nbins,xmineta,xmaxeta); // THnSparse(const char* name, const char* title, Int_t dim, // const Int_t* nbins, const Double_t* xmin, const Double_t* xmax, // Int_t chunksize); // to enable error bar calculation... fDPhiNumerator->Sumw2(); fDPhiDenominator->Sumw2(); fDCosNumerator->Sumw2(); fDCosDenominator->Sumw2(); fYtYtNumerator->Sumw2(); fPhi->Sumw2(); fEta->Sumw2(); fYtYtDenominator->Sumw2(); fPtCorrectionsNum->Sumw2(); fPtCorrectionsDen->Sumw2(); } } void AliFemtoCorrFctnDEtaDPhiCorrections::LoadCorrectionTabFromROOTFile(const char *file, ParticleType partType1, ParticleType partType2, bool doPtCorr, bool doEtaCorr, bool doPhiCorr, bool doZVertCorr) { fIfCorrection = kTRUE; ifileCorrTab = TFile::Open(file); fdoPtCorr = doPtCorr; fdoEtaCorr = doEtaCorr; fdoPhiCorr = doPhiCorr; fdoZVertCorr = doZVertCorr; fpartType1 = partType1; fpartType2 = partType2; char* type1 = new char[12]; char* type2 = new char[12]; if(fpartType1==kPion) strcpy(type1,"Pion"); else if(fpartType1==kKaon) strcpy(type1,"Kaon"); else if (fpartType1==kProton)strcpy(type1,"Proton"); else if (fpartType1==kAll) strcpy(type1,"All"); else if(fpartType1==kPionMinus) strcpy(type1,"PionMinus"); else if(fpartType1==kKaonMinus) strcpy(type1,"KaonMinus"); else if (fpartType1==kProtonMinus)strcpy(type1,"ProtonMinus"); else strcpy(type1,""); if(fpartType2==kPion) strcpy(type2,"Pion"); else if(fpartType2==kKaon) strcpy(type2,"Kaon"); else if (fpartType2==kProton) strcpy(type2,"Proton"); else if (fpartType2==kAll) strcpy(type2,"All"); else if(fpartType2==kPionMinus) strcpy(type1,"PionMinus"); else if(fpartType2==kKaonMinus) strcpy(type1,"KaonMinus"); else if (fpartType2==kProtonMinus)strcpy(type1,"ProtonMinus"); else strcpy(type1,""); fhntReco1 = (THnT*)(ifileCorrTab->Get(Form("fCorrectionMapData%s",type1)))->Clone(); fhntReco2 = (THnT*)(ifileCorrTab->Get(Form("fCorrectionMapData%s",type2)))->Clone(); fhCont1 = (TH1D*)(ifileCorrTab->Get(Form("SecondariesContamination%s",type1)))->Clone(); fhCont2 = (TH1D*)(ifileCorrTab->Get(Form("SecondariesContamination%s",type2)))->Clone(); delete type1; delete type2; double fhntReco1_nbins = fhntReco1->GetNbins(); double fhntReco2_nbins = fhntReco2->GetNbins(); int boolSum = fdoPtCorr+fdoEtaCorr+fdoPhiCorr+fdoZVertCorr; /*if(boolSum == 0) { return 1; }*/ if(boolSum == 1) { if(fdoPtCorr == 1) { fh1Reco1 = (TH1F*)(fhntReco1->Projection(0))->Clone(); fh1Reco2 = (TH1F*)(fhntReco2->Projection(0))->Clone(); fh1Reco1->Scale(1./fhntReco1_nbins*fh1Reco1->GetNbinsX()); fh1Reco2->Scale(1./fhntReco2_nbins*fh1Reco2->GetNbinsX()); } else if(fdoEtaCorr == 1) { fh1Reco1 = (TH1F*)(fhntReco1->Projection(1))->Clone(); fh1Reco2 = (TH1F*)(fhntReco2->Projection(1))->Clone(); fh1Reco1->Scale(1./fhntReco1_nbins*fh1Reco1->GetNbinsX()); fh1Reco2->Scale(1./fhntReco2_nbins*fh1Reco2->GetNbinsX()); } else if(fdoPhiCorr == 1) { fh1Reco1 = (TH1F*)(fhntReco1->Projection(2))->Clone(); fh1Reco2 = (TH1F*)(fhntReco2->Projection(2))->Clone(); fh1Reco1->Scale(1./fhntReco1_nbins*fh1Reco1->GetNbinsX()); fh1Reco2->Scale(1./fhntReco2_nbins*fh1Reco2->GetNbinsX()); } else if(fdoZVertCorr == 1) { fh1Reco1 = (TH1F*)(fhntReco1->Projection(3))->Clone(); fh1Reco2 = (TH1F*)(fhntReco2->Projection(3))->Clone(); fh1Reco1->Scale(1./fhntReco1_nbins*fh1Reco1->GetNbinsX()); fh1Reco2->Scale(1./fhntReco2_nbins*fh1Reco2->GetNbinsX()); } } else if(boolSum == 2) { if(fdoPtCorr == 1 && fdoEtaCorr == 1) { fh2Reco1 = (TH2F*)(fhntReco1->Projection(1,0))->Clone(); fh2Reco2 = (TH2F*)(fhntReco2->Projection(1,0))->Clone(); fh2Reco1->Scale(1./fhntReco1_nbins*fh2Reco1->GetNbinsX()*fh2Reco1->GetNbinsY()); fh2Reco2->Scale(1./fhntReco2_nbins*fh2Reco2->GetNbinsX()*fh2Reco2->GetNbinsY()); } else if(fdoPtCorr == 1 && fdoPhiCorr == 1) { fh2Reco1 = (TH2F*)(fhntReco1->Projection(2,0))->Clone(); fh2Reco2 = (TH2F*)(fhntReco2->Projection(2,0))->Clone(); fh2Reco1->Scale(1./fhntReco1_nbins*fh2Reco1->GetNbinsX()*fh2Reco1->GetNbinsY()); fh2Reco2->Scale(1./fhntReco2_nbins*fh2Reco2->GetNbinsX()*fh2Reco2->GetNbinsY()); } else if(fdoPtCorr == 1 && fdoZVertCorr == 1) { fh2Reco1 = (TH2F*)(fhntReco1->Projection(3,0))->Clone(); fh2Reco2 = (TH2F*)(fhntReco2->Projection(3,0))->Clone(); fh2Reco1->Scale(1./fhntReco1_nbins*fh2Reco1->GetNbinsX()*fh2Reco1->GetNbinsY()); fh2Reco2->Scale(1./fhntReco2_nbins*fh2Reco2->GetNbinsX()*fh2Reco2->GetNbinsY()); } else if(fdoEtaCorr == 1 && fdoPhiCorr == 1) { fh2Reco1 = (TH2F*)(fhntReco1->Projection(2,1))->Clone(); fh2Reco2 = (TH2F*)(fhntReco2->Projection(2,1))->Clone(); fh2Reco1->Scale(1./fhntReco1_nbins*fh2Reco1->GetNbinsX()*fh2Reco1->GetNbinsY()); fh2Reco2->Scale(1./fhntReco2_nbins*fh2Reco2->GetNbinsX()*fh2Reco2->GetNbinsY()); } else if(fdoEtaCorr == 1 && fdoZVertCorr == 1) { fh2Reco1 = (TH2F*)(fhntReco1->Projection(3,1))->Clone(); fh2Reco2 = (TH2F*)(fhntReco2->Projection(3,1))->Clone(); fh2Reco1->Scale(1./fhntReco1_nbins*fh2Reco1->GetNbinsX()*fh2Reco1->GetNbinsY()); fh2Reco2->Scale(1./fhntReco2_nbins*fh2Reco2->GetNbinsX()*fh2Reco2->GetNbinsY()); } else if(fdoPhiCorr == 1 && fdoZVertCorr == 1) { fh2Reco1 = (TH2F*)(fhntReco1->Projection(3,2))->Clone(); fh2Reco2 = (TH2F*)(fhntReco2->Projection(3,2))->Clone(); fh2Reco1->Scale(1./fhntReco1_nbins*fh2Reco1->GetNbinsX()*fh2Reco1->GetNbinsY()); fh2Reco2->Scale(1./fhntReco2_nbins*fh2Reco2->GetNbinsX()*fh2Reco2->GetNbinsY()); } } else if(boolSum == 3) { if(fdoPtCorr == 1 && fdoEtaCorr == 1 && fdoPhiCorr == 1) { fh3Reco1 = (TH3F*)(fhntReco1->Projection(0,1,2))->Clone(); fh3Reco2 = (TH3F*)(fhntReco2->Projection(0,1,2))->Clone(); fh3Reco1->Scale(1./fhntReco1_nbins*fh3Reco1->GetNbinsX()*fh3Reco1->GetNbinsY()*fh3Reco1->GetNbinsZ()); fh3Reco2->Scale(1./fhntReco2_nbins*fh3Reco2->GetNbinsX()*fh3Reco2->GetNbinsY()*fh3Reco2->GetNbinsZ()); } else if(fdoPtCorr == 1 && fdoEtaCorr == 1 && fdoZVertCorr == 1) { fh3Reco1 = (TH3F*)(fhntReco1->Projection(0,1,3))->Clone(); fh3Reco2 = (TH3F*)(fhntReco2->Projection(0,1,3))->Clone(); fh3Reco1->Scale(1./fhntReco1_nbins*fh3Reco1->GetNbinsX()*fh3Reco1->GetNbinsY()*fh3Reco1->GetNbinsZ()); fh3Reco2->Scale(1./fhntReco2_nbins*fh3Reco2->GetNbinsX()*fh3Reco2->GetNbinsY()*fh3Reco2->GetNbinsZ()); } else if(fdoPtCorr == 1 && fdoPhiCorr == 1 && fdoZVertCorr == 1) { fh3Reco1 = (TH3F*)(fhntReco1->Projection(0,2,3))->Clone(); fh3Reco2 = (TH3F*)(fhntReco2->Projection(0,2,3))->Clone(); fh3Reco1->Scale(1./fhntReco1_nbins*fh3Reco1->GetNbinsX()*fh3Reco1->GetNbinsY()*fh3Reco1->GetNbinsZ()); fh3Reco2->Scale(1./fhntReco2_nbins*fh3Reco2->GetNbinsX()*fh3Reco2->GetNbinsY()*fh3Reco2->GetNbinsZ()); } else if(fdoEtaCorr == 1 && fdoPhiCorr == 1 && fdoZVertCorr == 1) { fh3Reco1 = (TH3F*)(fhntReco1->Projection(1,2,3))->Clone(); fh3Reco2 = (TH3F*)(fhntReco2->Projection(1,2,3))->Clone(); fh3Reco1->Scale(1./fhntReco1_nbins*fh3Reco1->GetNbinsX()*fh3Reco1->GetNbinsY()*fh3Reco1->GetNbinsZ()); fh3Reco2->Scale(1./fhntReco2_nbins*fh3Reco2->GetNbinsX()*fh3Reco2->GetNbinsY()*fh3Reco2->GetNbinsZ()); } } /*else if(boolSum == 4) { }*/ ifileCorrTab->Close(); } void AliFemtoCorrFctnDEtaDPhiCorrections::LoadCorrectionTabFromROOTFile1D(const char *file, ParticleType partType1, ParticleType partType2) { fCorr1D = kTRUE; ifileCorrTab = TFile::Open(file); fpartType1 = partType1; fpartType2 = partType2; char type1[12]; char type2[12]; if(fpartType1==kPion) strcpy(type1,"Pion"); else if(fpartType1==kKaon) strcpy(type1,"Kaon"); else if (fpartType1==kProton)strcpy(type1,"Proton"); else if (fpartType1==kAll) strcpy(type1,"All"); else if(fpartType1==kPionMinus) strcpy(type1,"PionMinus"); else if (fpartType1==kKaonMinus) strcpy(type1,"KaonMinus"); else if (fpartType1==kProtonMinus)strcpy(type1,"ProtonMinus"); else strcpy(type1,""); if(fpartType2==kPion) strcpy(type2,"Pion"); else if(fpartType2==kKaon) strcpy(type2,"Kaon"); else if (fpartType2==kProton) strcpy(type2,"Proton"); else if (fpartType2==kAll) strcpy(type2,"All"); else if(fpartType2==kPionMinus) strcpy(type1,"PionMinus"); else if(fpartType2==kKaonMinus) strcpy(type1,"KaonMinus"); else if (fpartType2==kProtonMinus)strcpy(type1,"ProtonMinus"); else strcpy(type1,""); fhCont1 = (TH1D*)(ifileCorrTab->Get(Form("CorrectionFactorPtEffandCont%s",type1)));//->Clone(); fhCont2 = (TH1D*)(ifileCorrTab->Get(Form("CorrectionFactorPtEffandCont%s",type2)));//->Clone(); ifileCorrTab->Close(); } void AliFemtoCorrFctnDEtaDPhiCorrections::SetCorrectionTab(ParticleType partType) { double pttab[] = {0, 0.025, 0.05, 0.075, 0.1, 0.125, 0.15, 0.175, 0.2, 0.225, 0.25, 0.275, 0.3, 0.325, 0.35, 0.375, 0.4, 0.425, 0.45, 0.475, 0.5, 0.525, 0.55, 0.575, 0.6, 0.625, 0.65, 0.675, 0.7, 0.725, 0.75, 0.775, 0.8, 0.825, 0.85, 0.875, 0.9, 0.925, 0.95, 0.975, 1, 1.025, 1.05, 1.075, 1.1, 1.125, 1.15, 1.175, 1.2, 1.225, 1.25, 1.275, 1.3, 1.325, 1.35, 1.375, 1.4, 1.425, 1.45, 1.475, 1.5, 1.525, 1.55, 1.575, 1.6, 1.625, 1.65, 1.675, 1.7, 1.725, 1.75, 1.775, 1.8, 1.825, 1.85, 1.875, 1.9, 1.925, 1.95, 1.975, 2, 2.025, 2.05, 2.075, 2.1, 2.125, 2.15, 2.175, 2.2, 2.225, 2.25, 2.275, 2.3, 2.325, 2.35, 2.375, 2.4, 2.425, 2.45, 2.475, 2.5, 2.525, 2.55, 2.575, 2.6, 2.625, 2.65, 2.675, 2.7, 2.725, 2.75, 2.775, 2.8, 2.825, 2.85, 2.875, 2.9, 2.925, 2.95, 2.975, 3, 3.025, 3.05, 3.075, 3.1, 3.125, 3.15, 3.175, 3.2, 3.225, 3.25, 3.275, 3.3, 3.325, 3.35, 3.375, 3.4, 3.425, 3.45, 3.475, 3.5, 3.525, 3.55, 3.575, 3.6, 3.625, 3.65, 3.675, 3.7, 3.725, 3.75, 3.775, 3.8, 3.825, 3.85, 3.875, 3.9, 3.925, 3.95, 3.975, 4, 4.025, 4.05, 4.075, 4.1, 4.125, 4.15, 4.175, 4.2, 4.225, 4.25, 4.275, 4.3, 4.325, 4.35, 4.375, 4.4, 4.425, 4.45, 4.475, 4.5, 4.525, 4.55, 4.575, 4.6, 4.625, 4.65, 4.675, 4.7, 4.725, 4.75}; double pioncorrtab[] = {0, 0, 0, 0, 0, 0, 0, 1.40089, 1.40089, 1.29482, 1.29482, 1.25595, 1.22529, 1.22529, 1.23099, 1.32027, 1.32027, 1.44774, 1.44774, 1.74645, 1.8619, 1.8619, 1.82089, 1.78506, 1.78506, 1.75918, 1.75918, 1.74951, 1.74614, 1.74614, 1.74006, 1.73229, 1.73229, 1.72844, 1.72844, 1.72306, 1.71906, 1.71906, 1.71375, 1.71301, 1.71301, 1.70381, 1.70381, 1.69975, 1.69242, 1.69242, 1.69013, 1.67698, 1.67698, 1.6772, 1.6772, 1.67118, 1.66607, 1.66607, 1.66131, 1.67228, 1.67228, 1.66834, 1.66834, 1.66031, 1.6588, 1.6588, 1.6555, 1.64923, 1.64923, 1.6467, 1.6467, 1.63894, 1.63682, 1.63682, 1.6297, 1.62904, 1.62904, 1.63007, 1.63007, 1.62832, 1.62557, 1.62557, 1.62687, 1.62928, 1.62928, 1.62767, 1.62767, 1.62767, 1.62767, 1.62767, 1.62767, 1.63415, 1.63415, 1.63415, 1.63415, 1.63415, 1.63415, 1.63415, 1.64141, 1.64141, 1.64141, 1.64141, 1.64141, 1.64141, 1.65191, 1.65191, 1.65191, 1.65191, 1.65191, 1.65191, 1.65191, 1.66838, 1.66838, 1.66838, 1.66838, 1.66838, 1.66838, 1.6839, 1.6839, 1.6839, 1.6839, 1.6839, 1.6839, 1.69601, 1.69601, 1.69601, 1.69601, 1.69601, 1.69601, 1.69601, 1.70062, 1.70062, 1.70062, 1.70062, 1.70062, 1.70062, 1.68668, 1.68668, 1.68668, 1.68668, 1.68668, 1.68668, 1.68668, 1.68182, 1.68182, 1.68182, 1.68182, 1.68182, 1.68182, 1.681, 1.681, 1.681, 1.681, 1.681, 1.681, 1.67749, 1.67749, 1.67749, 1.67749, 1.67749, 1.67749, 1.67749, 1.66558, 1.67223, 1.67223, 1.67223, 1.67223, 1.67223, 1.67223, 1.67223, 1.67223, 1.67223, 1.67223, 1.67223, 1.67223, 1.67223, 1.67223, 1.67223, 1.66872, 1.66872, 1.66872, 1.66872, 1.66872, 1.66872, 1.66872, 1.66872, 1.66872, 1.66872, 1.66872, 1.66872, 1.66872, 1.66872, 1.66872, 1.66872, 1.64419}; double protoncorrtab[] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 357.585, 357.585, 8.66944, 2.10995, 2.10995, 1.50443, 1.73168, 1.73168, 2.3605, 2.3605, 4.7726, 4.40359, 4.40359, 3.0307, 2.49649, 2.49649, 2.2231, 2.2231, 2.11247, 2.05862, 2.05862, 2.00703, 1.9623, 1.9623, 1.93393, 1.93393, 1.9101, 1.89334, 1.89334, 1.87734, 1.86342, 1.86342, 1.85075, 1.85075, 1.83985, 1.83684, 1.83684, 1.82915, 1.81832, 1.81832, 1.81215, 1.81215, 1.7998, 1.79524, 1.79524, 1.78568, 1.79989, 1.79989, 1.7973, 1.7973, 1.79591, 1.78468, 1.78468, 1.78037, 1.77394, 1.77394, 1.77198, 1.77198, 1.76736, 1.76875, 1.76875, 1.76221, 1.75729, 1.75729, 1.75397, 1.75397, 1.75229, 1.74918, 1.74918, 1.75064, 1.75643, 1.75643, 1.75765, 1.75765, 1.75765, 1.75765, 1.75765, 1.75765, 1.76345, 1.76345, 1.76345, 1.76345, 1.76345, 1.76345, 1.76345, 1.76901, 1.76901, 1.76901, 1.76901, 1.76901, 1.76901, 1.78291, 1.78291, 1.78291, 1.78291, 1.78291, 1.78291, 1.78291, 1.80009, 1.80009, 1.80009, 1.80009, 1.80009, 1.80009, 1.81064, 1.81064, 1.81064, 1.81064, 1.81064, 1.81064, 1.81765, 1.81765, 1.81765, 1.81765, 1.81765, 1.81765, 1.81765, 1.79549, 1.79549, 1.79549, 1.79549, 1.79549, 1.79549, 1.80455, 1.80455, 1.80455, 1.80455, 1.80455, 1.80455, 1.80455, 1.78912, 1.78912, 1.78912, 1.78912, 1.78912, 1.78912, 1.78501, 1.78501, 1.78501, 1.78501, 1.78501, 1.78501, 1.79512, 1.79512, 1.79512, 1.79512, 1.79512, 1.79512, 1.79512, 1.77138, 1.784, 1.784, 1.784, 1.784, 1.784, 1.784, 1.784, 1.784, 1.784, 1.784, 1.784, 1.784, 1.784, 1.784, 1.784, 1.75152, 1.75152, 1.75152, 1.75152, 1.75152, 1.75152, 1.75152, 1.75152, 1.75152, 1.75152, 1.75152, 1.75152, 1.75152, 1.75152, 1.75152}; double kaoncorrtab[] = {0, 0, 0, 0, 0, 0, 0, 8.43268, 8.43268, 3.30657, 3.30657, 2.5102, 2.16256, 2.16256, 2.03757, 2.27166, 2.27166, 2.70432, 2.70432, 4.06234, 4.69199, 4.69199, 4.13074, 3.75139, 3.75139, 3.48381, 3.48381, 3.29762, 3.15261, 3.15261, 3.03022, 2.91874, 2.91874, 2.82421, 2.82421, 2.7388, 2.65961, 2.65961, 2.58426, 2.5174, 2.5174, 2.45378, 2.45378, 2.39687, 2.34699, 2.34699, 2.30247, 2.25299, 2.25299, 2.22443, 2.22443, 2.18303, 2.16012, 2.16012, 2.13083, 2.12806, 2.12806, 2.11376, 2.11376, 2.09566, 2.07526, 2.07526, 2.05378, 2.03252, 2.03252, 2.02466, 2.02466, 2.00531, 1.98945, 1.98945, 1.97877, 1.97226, 1.97226, 1.95475, 1.95475, 1.94838, 1.9314, 1.9314, 1.92571, 1.96346, 1.96346, 1.92849, 1.92849, 1.92849, 1.92849, 1.92849, 1.92849, 1.90949, 1.90949, 1.90949, 1.90949, 1.90949, 1.90949, 1.90949, 1.88743, 1.88743, 1.88743, 1.88743, 1.88743, 1.88743, 1.87486, 1.87486, 1.87486, 1.87486, 1.87486, 1.87486, 1.87486, 1.87785, 1.87785, 1.87785, 1.87785, 1.87785, 1.87785, 1.8757, 1.8757, 1.8757, 1.8757, 1.8757, 1.8757, 1.87948, 1.87948, 1.87948, 1.87948, 1.87948, 1.87948, 1.87948, 1.86148, 1.86148, 1.86148, 1.86148, 1.86148, 1.86148, 1.84329, 1.84329, 1.84329, 1.84329, 1.84329, 1.84329, 1.84329, 1.83105, 1.83105, 1.83105, 1.83105, 1.83105, 1.83105, 1.81955, 1.81955, 1.81955, 1.81955, 1.81955, 1.81955, 1.79944, 1.79944, 1.79944, 1.79944, 1.79944, 1.79944, 1.79944, 1.79345, 1.80077, 1.80077, 1.80077, 1.80077, 1.80077, 1.80077, 1.80077, 1.80077, 1.80077, 1.80077, 1.80077, 1.80077, 1.80077, 1.80077, 1.80077, 1.78333, 1.78333, 1.78333, 1.78333, 1.78333, 1.78333, 1.78333, 1.78333, 1.78333, 1.78333, 1.78333, 1.78333, 1.78333, 1.78333, 1.78333, 1.78333, 1.74958}; double allcorrtab[] = {0, 0, 0, 0, 0, 0, 0, 1.46883, 1.46883, 1.3528, 1.3528, 1.30939, 1.26936, 1.26936, 1.23645, 1.21359, 1.21359, 1.19759, 1.19759, 1.18565, 1.17772, 1.17772, 1.17203, 1.16739, 1.16739, 1.16398, 1.16398, 1.16201, 1.16065, 1.16065, 1.16012, 1.16009, 1.16009, 1.16044, 1.16044, 1.16104, 1.16139, 1.16139, 1.16134, 1.16278, 1.16278, 1.1631, 1.1631, 1.16227, 1.16152, 1.16152, 1.16066, 1.15984, 1.15984, 1.15932, 1.15932, 1.15912, 1.15818, 1.15818, 1.15877, 1.16754, 1.16754, 1.17075, 1.17075, 1.17047, 1.16995, 1.16995, 1.16885, 1.16845, 1.16845, 1.16824, 1.16824, 1.16771, 1.16704, 1.16704, 1.16681, 1.16723, 1.16723, 1.16819, 1.16819, 1.16811, 1.16974, 1.16974, 1.17217, 1.16759, 1.16759, 1.17376, 1.17376, 1.17376, 1.17376, 1.17376, 1.17376, 1.18247, 1.18247, 1.18247, 1.18247, 1.18247, 1.18247, 1.18247, 1.18916, 1.18916, 1.18916, 1.18916, 1.18916, 1.18916, 1.19649, 1.19649, 1.19649, 1.19649, 1.19649, 1.19649, 1.19649, 1.20315, 1.20315, 1.20315, 1.20315, 1.20315, 1.20315, 1.20984, 1.20984, 1.20984, 1.20984, 1.20984, 1.20984, 1.21236, 1.21236, 1.21236, 1.21236, 1.21236, 1.21236, 1.21236, 1.21272, 1.21272, 1.21272, 1.21272, 1.21272, 1.21272, 1.21416, 1.21416, 1.21416, 1.21416, 1.21416, 1.21416, 1.21416, 1.21308, 1.21308, 1.21308, 1.21308, 1.21308, 1.21308, 1.21332, 1.21332, 1.21332, 1.21332, 1.21332, 1.21332, 1.21204, 1.21204, 1.21204, 1.21204, 1.21204, 1.21204, 1.21204, 1.21006, 1.21141, 1.21141, 1.21141, 1.21141, 1.21141, 1.21141, 1.21141, 1.21141, 1.21141, 1.21141, 1.21141, 1.21141, 1.21141, 1.21141, 1.21141, 1.2092, 1.2092, 1.2092, 1.2092, 1.2092, 1.2092, 1.2092, 1.2092, 1.2092, 1.2092, 1.2092, 1.2092, 1.2092, 1.2092, 1.2092, 1.2092, 1.2066}; fpTab = new double[190]; for(int i=0;i<190;i++) fpTab[i]=pttab[i]; if(partType==kPion || partType==kPionMinus) { fCorrFactorTab = new double[190]; for(int i=0;i<190;i++) fCorrFactorTab[i] = pioncorrtab[i]; } else if(partType==kKaon || partType==kKaonMinus) { fCorrFactorTab = new double[190]; for(int i=0;i<190;i++) fCorrFactorTab[i] = kaoncorrtab[i]; } else if(partType==kProton||partType==kProtonMinus) { fCorrFactorTab = new double[190]; for(int i=0;i<190;i++) fCorrFactorTab[i] = protoncorrtab[i]; } else if(partType==kAll) { fCorrFactorTab = new double[190]; for(int i=0;i<190;i++) fCorrFactorTab[i] = allcorrtab[i]; } } double AliFemtoCorrFctnDEtaDPhiCorrections::CalculateCorrectionWeight(double pT1, double pT2) { double w1=0., w2=0.; if(pT1 > fhCont1->GetXaxis()->GetXmin() && pT1 < fhCont1->GetXaxis()->GetXmax() && pT2 > fhCont2->GetXaxis()->GetXmin() && pT2 < fhCont2->GetXaxis()->GetXmax()) { w1 = fhCont1->GetBinContent(fhCont1->FindFixBin(pT1)); w2 = fhCont2->GetBinContent(fhCont2->FindFixBin(pT2)); return w1*w2; } else return 0; } double AliFemtoCorrFctnDEtaDPhiCorrections::CalculateCorrectionWeight(double pT1) { double w1=0.; if(pT1 > fhCont1->GetXaxis()->GetXmin() && pT1 < fhCont1->GetXaxis()->GetXmax()) { w1 = fhCont1->GetBinContent(fhCont1->FindFixBin(pT1)); return w1; } else return 0; } double AliFemtoCorrFctnDEtaDPhiCorrections::CalculateCorrectionWeight(double pT1, double pT2, double eta1, double eta2, double phi1, double phi2, double zvert1, double zvert2) { double w1=0., w2=0.; double eps1=0., eps2=0; double cont1=0., cont2=0; //w=(1-cont)/eps phi1 += TMath::Pi(); phi2 += TMath::Pi(); if(pT1 > fhCont1->GetXaxis()->GetXmin() && pT1 < fhCont1->GetXaxis()->GetXmax() && pT2 > fhCont2->GetXaxis()->GetXmin() && pT2 < fhCont2->GetXaxis()->GetXmax()) { cont1 = fhCont1->GetBinContent(fhCont1->FindFixBin(pT1)); cont2 = fhCont1->GetBinContent(fhCont2->FindFixBin(pT2)); } else return 0; int boolSum = fdoPtCorr+fdoEtaCorr+fdoPhiCorr+fdoZVertCorr; if(boolSum == 0) { return 1; } else if(boolSum == 1) { if(fdoPtCorr == 1) { if(pT1 > fh1Reco1->GetXaxis()->GetXmin() && pT1 < fh1Reco1->GetXaxis()->GetXmax() && pT2 > fh1Reco2->GetXaxis()->GetXmin() && pT2 < fh1Reco2->GetXaxis()->GetXmax()) { eps1 = fh1Reco1->GetBinContent(fh1Reco1->FindFixBin(pT1)); eps2 = fh1Reco2->GetBinContent(fh1Reco2->FindFixBin(pT2)); w1 = (1-cont1)/eps1; w2 = (1-cont2)/eps2; return w1*w2; } else return 0; } else if(fdoEtaCorr == 1) { if(eta1 > fh1Reco1->GetXaxis()->GetXmin() && eta1 < fh1Reco1->GetXaxis()->GetXmax() && eta2 > fh1Reco2->GetXaxis()->GetXmin() && eta2 < fh1Reco2->GetXaxis()->GetXmax()) { eps1 = fh1Reco1->GetBinContent(fh1Reco1->FindFixBin(eta1)); eps2 = fh1Reco2->GetBinContent(fh1Reco2->FindFixBin(eta2)); w1 = (1-cont1)/eps1; w2 = (1-cont2)/eps2; return w1*w2; } else return 0; } else if(fdoPhiCorr == 1) { if(phi1 > fh1Reco1->GetXaxis()->GetXmin() && phi1 < fh1Reco1->GetXaxis()->GetXmax() && phi2 > fh1Reco2->GetXaxis()->GetXmin() && phi2 < fh1Reco2->GetXaxis()->GetXmax()) { eps1 = fh1Reco1->GetBinContent(fh1Reco1->FindFixBin(phi1)); eps2 = fh1Reco2->GetBinContent(fh1Reco2->FindFixBin(phi2)); w1 = (1-cont1)/eps1; w2 = (1-cont2)/eps2; return w1*w2; } else return 0; } else if(fdoZVertCorr == 1) { if(zvert1 > fh1Reco1->GetXaxis()->GetXmin() && zvert1 < fh1Reco1->GetXaxis()->GetXmax() && zvert2 > fh1Reco2->GetXaxis()->GetXmin() && zvert2 < fh1Reco2->GetXaxis()->GetXmax()) { eps1 = fh1Reco1->GetBinContent(fh1Reco1->FindFixBin(zvert1)); eps2 = fh1Reco2->GetBinContent(fh1Reco2->FindFixBin(zvert2)); w1 = (1-cont1)/eps1; w2 = (1-cont2)/eps2; return w1*w2; } else return 0; } } else if(boolSum == 2) { if(fdoPtCorr == 1 && fdoEtaCorr == 1) { if(pT1 > fh2Reco1->GetXaxis()->GetXmin() && pT1 < fh2Reco1->GetXaxis()->GetXmax() && pT2 > fh2Reco2->GetXaxis()->GetXmin() && pT2 < fh2Reco2->GetXaxis()->GetXmax() && eta1 > fh2Reco1->GetYaxis()->GetXmin() && eta1 < fh2Reco1->GetYaxis()->GetXmax() && eta2 > fh2Reco2->GetYaxis()->GetXmin() && eta2 < fh2Reco2->GetYaxis()->GetXmax()) { eps1 = fh2Reco1->GetBinContent(fh2Reco1->GetXaxis()->FindFixBin(pT1),fh2Reco1->GetYaxis()->FindFixBin(eta1)); eps2 = fh2Reco2->GetBinContent(fh2Reco2->GetXaxis()->FindFixBin(pT2),fh2Reco2->GetYaxis()->FindFixBin(eta2)); w1 = (1-cont1)/eps1; w2 = (1-cont2)/eps2; return w1*w2; } else return 0; } if(fdoPtCorr == 1 && fdoPhiCorr == 1) { if(pT1 > fh2Reco1->GetXaxis()->GetXmin() && pT1 < fh2Reco1->GetXaxis()->GetXmax() && pT2 > fh2Reco2->GetXaxis()->GetXmin() && pT2 < fh2Reco2->GetXaxis()->GetXmax() && phi1 > fh2Reco1->GetYaxis()->GetXmin() && phi1 < fh2Reco1->GetYaxis()->GetXmax() && phi2 > fh2Reco2->GetYaxis()->GetXmin() && phi2 < fh2Reco2->GetYaxis()->GetXmax()) { eps1 = fh2Reco1->GetBinContent(fh2Reco1->GetXaxis()->FindFixBin(pT1),fh2Reco1->GetYaxis()->FindFixBin(phi1)); eps2 = fh2Reco2->GetBinContent(fh2Reco2->GetXaxis()->FindFixBin(pT2),fh2Reco2->GetYaxis()->FindFixBin(phi2)); w1 = (1-cont1)/eps1; w2 = (1-cont2)/eps2; return w1*w2; } else return 0; } else if(fdoPtCorr == 1 && fdoZVertCorr == 1) { if(pT1 > fh2Reco1->GetXaxis()->GetXmin() && pT1 < fh2Reco1->GetXaxis()->GetXmax() && pT2 > fh2Reco2->GetXaxis()->GetXmin() && pT2 < fh2Reco2->GetXaxis()->GetXmax() && zvert1 > fh2Reco1->GetYaxis()->GetXmin() && zvert1 < fh2Reco1->GetYaxis()->GetXmax() && zvert2 > fh2Reco2->GetYaxis()->GetXmin() && zvert2 < fh2Reco2->GetYaxis()->GetXmax()) { eps1 = fh2Reco1->GetBinContent(fh2Reco1->GetXaxis()->FindFixBin(pT1),fh2Reco1->GetYaxis()->FindFixBin(zvert1)); eps2 = fh2Reco2->GetBinContent(fh2Reco2->GetXaxis()->FindFixBin(pT2),fh2Reco2->GetYaxis()->FindFixBin(zvert2)); w1 = (1-cont1)/eps1; w2 = (1-cont2)/eps2; return w1*w2; } else return 0; } else if(fdoEtaCorr == 1 && fdoPhiCorr == 1) { if(eta1 > fh2Reco1->GetXaxis()->GetXmin() && eta1 < fh2Reco1->GetXaxis()->GetXmax() && eta2 > fh2Reco2->GetXaxis()->GetXmin() && eta2 < fh2Reco2->GetXaxis()->GetXmax() && phi1 > fh2Reco1->GetYaxis()->GetXmin() && phi1 < fh2Reco1->GetYaxis()->GetXmax() && phi2 > fh2Reco2->GetYaxis()->GetXmin() && phi2 < fh2Reco2->GetYaxis()->GetXmax()) { eps1 = fh2Reco1->GetBinContent(fh2Reco1->GetXaxis()->FindFixBin(eta1),fh2Reco1->GetYaxis()->FindFixBin(phi1)); eps2 = fh2Reco2->GetBinContent(fh2Reco2->GetXaxis()->FindFixBin(eta2),fh2Reco2->GetYaxis()->FindFixBin(phi2)); w1 = (1-cont1)/eps1; w2 = (1-cont2)/eps2; return w1*w2; } else return 0; } else if(fdoEtaCorr == 1 && fdoZVertCorr == 1) { if(eta1 > fh2Reco1->GetXaxis()->GetXmin() && eta1 < fh2Reco1->GetXaxis()->GetXmax() && eta2 > fh2Reco2->GetXaxis()->GetXmin() && eta2 < fh2Reco2->GetXaxis()->GetXmax() && zvert1 > fh2Reco1->GetYaxis()->GetXmin() && zvert1 < fh2Reco1->GetYaxis()->GetXmax() && zvert2 > fh2Reco2->GetXaxis()->GetXmin() && zvert2 < fh2Reco2->GetXaxis()->GetXmax()) { eps1 = fh2Reco1->GetBinContent(fh2Reco1->GetXaxis()->FindFixBin(eta1),fh2Reco1->GetYaxis()->FindFixBin(zvert1)); eps1 = fh2Reco2->GetBinContent(fh2Reco2->GetXaxis()->FindFixBin(eta2),fh2Reco2->GetYaxis()->FindFixBin(zvert2)); w1 = (1-cont1)/eps1; w2 = (1-cont2)/eps2; return w1*w2; } else return 0; } else if(fdoPhiCorr == 1 && fdoZVertCorr == 1) { if(phi1 > fh2Reco1->GetXaxis()->GetXmin() && phi1 < fh2Reco1->GetXaxis()->GetXmax() && phi2 > fh2Reco2->GetXaxis()->GetXmin() && phi2 < fh2Reco2->GetXaxis()->GetXmax() && zvert1 > fh2Reco1->GetYaxis()->GetXmin() && zvert1 < fh2Reco1->GetYaxis()->GetXmax() && zvert2 > fh2Reco2->GetYaxis()->GetXmin() && zvert2 < fh2Reco2->GetYaxis()->GetXmax()) { eps1 = fh2Reco1->GetBinContent(fh2Reco1->GetXaxis()->FindFixBin(phi1),fh2Reco1->GetYaxis()->FindFixBin(zvert1)); eps2 = fh2Reco2->GetBinContent(fh2Reco2->GetXaxis()->FindFixBin(phi2),fh2Reco2->GetYaxis()->FindFixBin(zvert2)); w1 = (1-cont1)/eps1; w2 = (1-cont2)/eps2; return w1*w2; } else return 0; } } else if(boolSum == 3) { if(fdoPtCorr == 1 && fdoEtaCorr == 1 && fdoPhiCorr == 1) { if(pT1 >fh3Reco1->GetXaxis()->GetXmin() && pT1 GetXaxis()->GetXmax() && pT2 > fh3Reco2->GetXaxis()->GetXmin() && pT2 GetXaxis()->GetXmax() && eta1 > fh3Reco1->GetYaxis()->GetXmin() && eta1 GetYaxis()->GetXmax() && eta2 > fh3Reco2->GetYaxis()->GetXmin() && eta2 GetYaxis()->GetXmax() && phi1 > fh3Reco1->GetZaxis()->GetXmin() && phi1 < fh3Reco1->GetZaxis()->GetXmax() && phi2 > fh3Reco2->GetZaxis()->GetXmin() && phi2 < fh3Reco2->GetZaxis()->GetXmax()) { eps1 = fh3Reco1->GetBinContent(fh3Reco1->GetXaxis()->FindFixBin(pT1),fh3Reco1->GetYaxis()->FindFixBin(eta1),fh3Reco1->GetZaxis()->FindFixBin(phi1)); eps2 = fh3Reco2->GetBinContent(fh3Reco2->GetXaxis()->FindFixBin(pT2),fh3Reco2->GetYaxis()->FindFixBin(eta2),fh3Reco2->GetZaxis()->FindFixBin(phi2)); w1 = (1-cont1)/eps1; w2 = (1-cont2)/eps2; return w1*w2; } else return 0; } else if(fdoPtCorr == 1 && fdoEtaCorr == 1 && fdoZVertCorr == 1) { if(pT1 >fh3Reco1->GetXaxis()->GetXmin() && pT1 GetXaxis()->GetXmax() && pT2 > fh3Reco2->GetXaxis()->GetXmin() && pT2 GetXaxis()->GetXmax() && eta1 > fh3Reco1->GetYaxis()->GetXmin() && eta1 GetYaxis()->GetXmax() && eta2 > fh3Reco2->GetYaxis()->GetXmin() && eta2 GetYaxis()->GetXmax() && zvert1 > fh3Reco1->GetZaxis()->GetXmin() && zvert1 < fh3Reco1->GetZaxis()->GetXmax() && zvert2 > fh3Reco2->GetZaxis()->GetXmin() && zvert2 < fh3Reco2->GetZaxis()->GetXmax()) { eps1 = fh3Reco1->GetBinContent(fh3Reco1->GetXaxis()->FindFixBin(pT1),fh3Reco1->GetYaxis()->FindFixBin(eta1),fh3Reco1->GetZaxis()->FindFixBin(zvert1)); eps2 = fh3Reco2->GetBinContent(fh3Reco2->GetXaxis()->FindFixBin(pT2),fh3Reco2->GetYaxis()->FindFixBin(eta2),fh3Reco2->GetZaxis()->FindFixBin(zvert2)); w1 = (1-cont1)/eps1; w2 = (1-cont2)/eps2; return w1*w2; } else return 0; } else if(fdoPtCorr == 1 && fdoPhiCorr == 1 && fdoZVertCorr == 1) { if(pT1 >fh3Reco1->GetXaxis()->GetXmin() && pT1 GetXaxis()->GetXmax() && pT2 > fh3Reco2->GetXaxis()->GetXmin() && pT2 GetXaxis()->GetXmax() && phi1 > fh3Reco1->GetYaxis()->GetXmin() && phi1 GetYaxis()->GetXmax() && phi2 > fh3Reco2->GetYaxis()->GetXmin() && phi2 GetYaxis()->GetXmax() && zvert1 > fh3Reco1->GetZaxis()->GetXmin() && zvert1 < fh3Reco1->GetZaxis()->GetXmax() && zvert2 > fh3Reco2->GetZaxis()->GetXmin() && zvert2 < fh3Reco2->GetZaxis()->GetXmax()) { eps1 = fh3Reco1->GetBinContent(fh3Reco1->GetXaxis()->FindFixBin(pT1),fh3Reco1->GetYaxis()->FindFixBin(phi1),fh3Reco1->GetZaxis()->FindFixBin(zvert1)); eps2 = fh3Reco2->GetBinContent(fh3Reco2->GetXaxis()->FindFixBin(pT2),fh3Reco2->GetYaxis()->FindFixBin(phi2),fh3Reco2->GetZaxis()->FindFixBin(zvert2)); w1 = (1-cont1)/eps1; w2 = (1-cont2)/eps2; return w1*w2; } else return 0; } else if(fdoEtaCorr == 1 && fdoPhiCorr == 1 && fdoZVertCorr == 1) { if(eta1 >fh3Reco1->GetXaxis()->GetXmin() && eta1 GetXaxis()->GetXmax() && eta2 > fh3Reco2->GetXaxis()->GetXmin() && eta2 GetXaxis()->GetXmax() && phi1 > fh3Reco1->GetYaxis()->GetXmin() && phi1 GetYaxis()->GetXmax() && phi2 > fh3Reco2->GetYaxis()->GetXmin() && phi2 GetYaxis()->GetXmax() && zvert1 > fh3Reco1->GetZaxis()->GetXmin() && zvert1 < fh3Reco1->GetZaxis()->GetXmax() && zvert2 > fh3Reco2->GetZaxis()->GetXmin() && zvert2 < fh3Reco2->GetZaxis()->GetXmax()) { eps1 = fh3Reco1->GetBinContent(fh3Reco1->GetXaxis()->FindFixBin(eta1),fh3Reco1->GetYaxis()->FindFixBin(phi1),fh3Reco1->GetZaxis()->FindFixBin(zvert1)); eps2 = fh3Reco2->GetBinContent(fh3Reco2->GetXaxis()->FindFixBin(eta2),fh3Reco2->GetYaxis()->FindFixBin(phi2),fh3Reco2->GetZaxis()->FindFixBin(zvert2)); w1 = (1-cont1)/eps1; w2 = (1-cont2)/eps2; return w1*w2; } else return 0; } } else if(boolSum == 4) { if(pT1 > fhntReco1->GetAxis(0)->GetXmin() && pT1 < fhntReco1->GetAxis(0)->GetXmax() && pT2 > fhntReco2->GetAxis(0)->GetXmin() && pT2 < fhntReco2->GetAxis(0)->GetXmax() && eta1 > fhntReco1->GetAxis(1)->GetXmin() && eta1 GetAxis(1)->GetXmax() && eta2 > fhntReco2->GetAxis(1)->GetXmin() && eta2 < fhntReco2->GetAxis(1)->GetXmax() && phi1 > fhntReco1->GetAxis(2)->GetXmin() && phi2 < fhntReco2->GetAxis(2)->GetXmax() && phi2 > fhntReco2->GetAxis(2)->GetXmin() && phi2 < fhntReco2->GetAxis(2)->GetXmax() && zvert1 > fhntReco1->GetAxis(3)->GetXmin() && zvert1 < fhntReco1->GetAxis(3)->GetXmax() && zvert2 > fhntReco2->GetAxis(3)->GetXmin() && zvert2 < fhntReco2->GetAxis(3)->GetXmax()) { int tab1[] = {fhntReco1->GetAxis(0)->FindFixBin(pT1),fhntReco1->GetAxis(1)->FindFixBin(eta1),fhntReco1->GetAxis(2)->FindFixBin(phi1),fhntReco1->GetAxis(3)->FindFixBin(zvert1)}; int tab2[] = {fhntReco2->GetAxis(0)->FindFixBin(pT2),fhntReco2->GetAxis(1)->FindFixBin(eta2),fhntReco2->GetAxis(2)->FindFixBin(phi2),fhntReco2->GetAxis(3)->FindFixBin(zvert2)}; eps1 = fhntReco1->GetBinContent(tab1); eps2 = fhntReco2->GetBinContent(tab2); w1 = (1-cont1)/eps1; w2 = (1-cont2)/eps2; return w1*w2; } else return 0; } return 0; }