//
///////////////////////////////////////////////////////////////////////////////
-enum centrality{ kpp, k07half, kpPb0100, k010, k1020, k020, k2040, k2030, k3040, k4050, k3050, k5060, k4060, k6080, k4080, k80100 };
+enum centrality{ kpp, k07half, kpPb0100, k010, k1020, k020, k2040, k2030, k3040, k4050, k3050, k5060, k4060, k6080, k4080, k5080, k80100, kpPb020, kpPb2040, kpPb4060, kpPb60100 };
+enum centestimator{ kV0M, kV0A, kZNA };
enum energy{ k276, k5dot023, k55 };
enum BFDSubtrMethod { kfc, kNb };
enum RaavsEP {kPhiIntegrated, kInPlane, kOutOfPlane};
Bool_t printout = false;
+Double_t ptprintout = 1.5;
Double_t NormPPUnc = 0.035;
+Double_t NormABUnc = 0.037;
Bool_t elossFDQuadSum = true;
//____________________________________________________________
Int_t fdMethod = kNb, Int_t cc=kpp, Int_t Energy=k276,
Double_t MinHypo=1./3., Double_t MaxHypo=3.0, Double_t MaxRb=6.0,
Bool_t isRbHypo=false, Double_t CentralHypo = 1.0,
+ Int_t ccestimator = kV0M,
+ Bool_t isUseTaaForRaa=true, const char *shadRbcFile="", Int_t nSigmaShad=3.0,
Int_t isRaavsEP=kPhiIntegrated, Bool_t isScaledAndExtrapRef=kFALSE)
{
//
// Defining the TAB values for the given centrality class
//
- Double_t Tab = 1., TabSyst = 0.;
+ Double_t Tab = 1., TabSyst = 0., A=207.2, B=207.2;
if ( Energy!=k276 && Energy!=k5dot023) {
printf("\n The Tab values for this cms energy have not yet been implemented, please do it ! \n");
return;
}
// Values from Alberica's twiki:
// https://twiki.cern.ch/twiki/bin/viewauth/ALICE/CentStudies
- if ( cc == k07half ) {
- Tab = 24.81; TabSyst = 0.8037;
- } else if ( cc == k010 ) {
- Tab = 23.48; TabSyst = 0.97;
- } else if ( cc == k1020 ) {
- Tab = 14.4318; TabSyst = 0.5733;
- } else if ( cc == k020 ) {
- Tab = 18.93; TabSyst = 0.74;
- } else if ( cc == k2040 ) {
- Tab = 6.86; TabSyst = 0.28;
- } else if ( cc == k2030 ) {
- Tab = 8.73769; TabSyst = 0.370219;
- } else if ( cc == k3040 ) {
- Tab = 5.02755; TabSyst = 0.22099;
- } else if ( cc == k4050 ) {
- Tab = 2.68327; TabSyst = 0.137073;
- } else if ( cc == k3050 ) {
- Tab = 3.87011; TabSyst = 0.183847;
- } else if ( cc == k4060 ) {
- Tab = 2.00; TabSyst= 0.11;
- } else if ( cc == k4080 ) {
- Tab = 1.20451; TabSyst = 0.071843;
- } else if ( cc == k5060 ) {
- Tab = 1.32884; TabSyst = 0.0929536;
- } else if ( cc == k6080 ) {
- Tab = 0.419; TabSyst = 0.033;
- } else if ( cc == k80100 ){
- Tab = 0.0690; TabSyst = 0.0062;
+ if( ccestimator == kV0M ) {
+ if ( cc == k07half ) {
+ Tab = 24.81; TabSyst = 0.8037;
+ } else if ( cc == k010 ) {
+ Tab = 23.48; TabSyst = 0.97;
+ } else if ( cc == k1020 ) {
+ Tab = 14.4318; TabSyst = 0.5733;
+ } else if ( cc == k020 ) {
+ Tab = 18.93; TabSyst = 0.74;
+ } else if ( cc == k2040 ) {
+ Tab = 6.86; TabSyst = 0.28;
+ } else if ( cc == k2030 ) {
+ Tab = 8.73769; TabSyst = 0.370219;
+ } else if ( cc == k3040 ) {
+ Tab = 5.02755; TabSyst = 0.22099;
+ } else if ( cc == k4050 ) {
+ Tab = 2.68327; TabSyst = 0.137073;
+ } else if ( cc == k3050 ) {
+ Tab = 3.87011; TabSyst = 0.183847;
+ } else if ( cc == k4060 ) {
+ Tab = 2.00; TabSyst= 0.11;
+ } else if ( cc == k4080 ) {
+ Tab = 1.20451; TabSyst = 0.071843;
+ } else if ( cc == k5060 ) {
+ Tab = 1.32884; TabSyst = 0.0929536;
+ } else if ( cc == k6080 ) {
+ Tab = 0.419; TabSyst = 0.033;
+ } else if ( cc == k5080 ) {
+ Tab = 0.719; TabSyst = 0.054;
+ } else if ( cc == k80100 ){
+ Tab = 0.0690; TabSyst = 0.0062;
+ }
}
// pPb Glauber (A. Toia)
// https://twiki.cern.ch/twiki/bin/viewauth/ALICE/PACentStudies#Glauber_Calculations_with_sigma
- else if( cc == kpPb0100 ){
+ if( cc == kpPb0100 ){
Tab = 0.098334; TabSyst = 0.0070679;
+ A=207.2; B=1.;
+ }
+ else if( ccestimator == kV0A ){
+ if ( cc == kpPb020 ) {
+ Tab = 0.183; TabSyst = 0.006245;
+ } else if ( cc == kpPb2040 ) {
+ Tab = 0.134; TabSyst = 0.004899;
+ } else if ( cc == kpPb4060 ) {
+ Tab = 0.092; TabSyst = 0.004796;
+ } else if ( cc == kpPb60100 ) {
+ Tab = 0.041; TabSyst = 0.008832;
+ }
+ }
+ else if( ccestimator == kZNA ){
+ if ( cc == kpPb020 ) {
+ Tab = 0.164; TabSyst = 0.010724;
+ } else if ( cc == kpPb2040 ) {
+ Tab = 0.137; TabSyst = 0.005099;
+ } else if ( cc == kpPb4060 ) {
+ Tab = 0.1011; TabSyst = 0.006;
+ } else if ( cc == kpPb60100 ) {
+ Tab = 0.0459; TabSyst = 0.003162;
+ }
}
//
// Reading the pp file
//
TFile * ppf = new TFile(ppfile,"read");
- TH1D * hSigmaPP;// = (TH1D*)ppf->Get("fhScaledData");
- TGraphAsymmErrors * gSigmaPPSyst;// = (TGraphAsymmErrors*)ppf->Get("gScaledData");
+ TH1D * hSigmaPP;
+ TGraphAsymmErrors * gSigmaPPSyst;
TGraphAsymmErrors * gSigmaPPSystData = (TGraphAsymmErrors*)ppf->Get("gScaledDataSystData");
TGraphAsymmErrors * gSigmaPPSystTheory = (TGraphAsymmErrors*)ppf->Get("gScaledDataSystExtrap");
TGraphAsymmErrors * gSigmaPPSystFeedDown = (TGraphAsymmErrors*)ppf->Get("gScaledDataSystFeedDown");
//
TFile * ABf = new TFile(ABfile,"read");
TH1D *hSigmaAB = (TH1D*)ABf->Get("histoSigmaCorr");
- TH2D *hSigmaABRcb = (TH2D*)ABf->Get("histoSigmaCorrRcb");
- TGraphAsymmErrors * gSigmaABSyst = (TGraphAsymmErrors*)ABf->Get("gSigmaCorr");
+ // TH2D *hSigmaABRcb = (TH2D*)ABf->Get("histoSigmaCorrRcb");
+ // TGraphAsymmErrors * gSigmaABSyst = (TGraphAsymmErrors*)ABf->Get("gSigmaCorr");
TGraphAsymmErrors * gSigmaABSystFeedDown = (TGraphAsymmErrors*)ABf->Get("gSigmaCorrConservative");
TNtuple * nSigmaAB = (TNtuple*)ABf->Get("fnSigma");
//
systematicsAB->SetIsPbPb2010EnergyScan(true);
}
else if ( cc == k4060 || cc == k4050 || cc == k5060 ) systematicsAB->SetCentrality("4060");
- else if ( cc == k6080 ) systematicsAB->SetCentrality("6080");
+ else if ( cc == k6080 || cc == k5080 ) systematicsAB->SetCentrality("6080");
else if ( cc == k4080 ) systematicsAB->SetCentrality("4080");
else if ( cc == k3050 ) {
if (isRaavsEP == kPhiIntegrated) systematicsAB->SetCentrality("4080");
else if (isRaavsEP == kOutOfPlane) systematicsAB->SetCentrality("3050OutOfPlane");
}
//
- else if ( cc == kpPb0100 ){
- systematicsAB->SetCollisionType(0);
- cout <<endl<<" Beware pPb systematics not yet implemented, using pp at 7 TeV !! "<<endl<<endl;
+ else if ( cc == kpPb0100 || cc == kpPb020 || cc == kpPb2040 || cc == kpPb4060 || cc == kpPb60100 ) {
+ systematicsAB->SetCollisionType(2);
+ if(ccestimator==kV0A) {
+ if(cc == kpPb020) systematicsAB->SetCentrality("020V0A");
+ else if(cc == kpPb2040) systematicsAB->SetCentrality("2040V0A");
+ else if(cc == kpPb4060) systematicsAB->SetCentrality("4060V0A");
+ else if(cc == kpPb60100) systematicsAB->SetCentrality("60100V0A");
+ } else if (ccestimator==kZNA) {
+ if(cc == kpPb020) systematicsAB->SetCentrality("020ZNA");
+ else if(cc == kpPb2040) systematicsAB->SetCentrality("2040ZNA");
+ else if(cc == kpPb4060) systematicsAB->SetCentrality("4060ZNA");
+ else if(cc == kpPb60100) systematicsAB->SetCentrality("60100ZNA");
+ } else {
+ if(!(cc == kpPb0100)) {
+ cout <<" Error on the pPb options"<<endl;
+ return;
+ }
+ }
}
else {
cout << " Systematics not yet implemented " << endl;
binwidths[i-1] = binwidth;
}
limits[nbins] = xlow + binwidth;
+
+
+ //
+ // Read the shadowing file if given as input
+ //
+ Double_t centralRbcShad[nbins+1], minRbcShad[nbins+1], maxRbcShad[nbins+1];
+ for(Int_t i=0; i<=nbins; i++) { centralRbcShad[i]=1.0; minRbcShad[i]=6.0; maxRbcShad[i]=0.0; }
+ Bool_t isShadHypothesis = false;
+ if( strcmp(shadRbcFile,"")!=0 ) {
+ isShadHypothesis = true;
+ cout<<endl<<">> Beware, using the shadowing prediction file with an "<<nSigmaShad<<"*sigma <<"<<endl<<endl;
+ TFile *fshad = new TFile(shadRbcFile,"read");
+ if(!fshad){ cout <<" >> Shadowing file not properly opened!!!"<<endl<<endl; return;}
+ // TH1D *hRbcShadCentral = (TH1D*)fshad->Get("hDfromBoverPromptD_Shadowing_central");
+ // TH1D *hRbcShadMin = (TH1D*)fshad->Get("hDfromBoverPromptD_Shadowing_upper");
+ // TH1D *hRbcShadMax = (TH1D*)fshad->Get("hDfromBoverPromptD_Shadowing_lower");
+ TH1D *hRbcShadCentral = (TH1D*)fshad->Get("hDfromBoverDfromc_L0");
+ TH1D *hRbcShadMin = (TH1D*)fshad->Get("hDfromBoverDfromc_L0");
+ TH1D *hRbcShadMax = (TH1D*)fshad->Get("hDfromBoverDfromc_L1");
+ if(!hRbcShadCentral || !hRbcShadMin || !hRbcShadMax) {
+ cout<< endl <<">> Shadowing input histograms are not ok !! "<<endl<<endl;
+ return;
+ }
+ // nSigmaShad
+ // nSigmaShad
+ for(Int_t i=1; i<=nbins; i++) {
+ Double_t xpt = hSigmaAB->GetBinCenter(i);
+ if(xpt>24) xpt = 20;
+ centralRbcShad[i] = hRbcShadCentral->GetBinContent( hRbcShadCentral->FindBin(xpt) );
+ Double_t minValue0 = hRbcShadMin->GetBinContent( hRbcShadMin->FindBin(xpt) );
+ Double_t maxValue0 = hRbcShadMax->GetBinContent( hRbcShadMax->FindBin(xpt) );
+ Double_t arrayEl[3] = {minValue0,maxValue0, centralRbcShad[i]};
+ Double_t minValue = TMath::MinElement(3,arrayEl);
+ Double_t maxValue = TMath::MaxElement(3,arrayEl);
+ cout<<">> Shadowing pt="<<xpt<<" central="<<centralRbcShad[i]<<" min="<<minValue<<" max="<<maxValue<<endl;
+ if(minValue>centralRbcShad[i]){ minValue = centralRbcShad[i]; }
+ if(maxValue<centralRbcShad[i]){ maxValue = centralRbcShad[i]; }
+ minRbcShad[i] = centralRbcShad[i] - nSigmaShad*(centralRbcShad[i] - minValue);
+ maxRbcShad[i] = centralRbcShad[i] + nSigmaShad*(maxValue - centralRbcShad[i]);
+ cout<<">> Shadowing hypothesis pt="<<xpt<<" central="<<centralRbcShad[i]<<" min="<<minRbcShad[i]<<" max="<<maxRbcShad[i]<<endl;
+ }
+ }
+
//
// define the bins correspondence bw histos/files/graphs
//
//
TH2D * hRABvsRcb = new TH2D("hRABvsRcb"," R_{AB}(c) vs Rcb Eloss hypothesis; p_{T} [GeV/c] ; R_{AB}(c) ; Rcb Eloss hypothesis ",nbins,limits,800,0.,4.);
TH2D * hRABvsRb = new TH2D("hRABvsRb"," R_{AB}(c) vs Rb Eloss hypothesis; p_{T} [GeV/c] ; R_{AB}(c) ; Rb Eloss hypothesis ",nbins,limits,800,0.,4.);
- TH2D * hRABBeautyvsRCharm = new TH2D("hRABBeautyvsRCharm"," R_{AB}(c) vs Rb Eloss hypothesis; p_{T} [GeV/c] ; R_{AB}(b) ; R_{AB}(c) ",nbins,limits,800,0.,4.);
+ // TH2D * hRABBeautyvsRCharm = new TH2D("hRABBeautyvsRCharm"," R_{AB}(c) vs Rb Eloss hypothesis; p_{T} [GeV/c] ; R_{AB}(b) ; R_{AB}(c) ",nbins,limits,800,0.,4.);
Int_t nbinsHypo=800;//200;
Double_t *limitsHypo = new Double_t[nbinsHypo+1];
for(Int_t i=1; i<=nbinsHypo+1; i++) limitsHypo[i-1]= i*4./800.;
//
TNtuple *ntupleRAB=0x0 ;
if (fdMethod==kNb) {
- ntupleRAB = new TNtuple("ntupleRAB","ntupleRAB (Nb)","pt:TAB:sigmaPP:sigmaAB:invyieldAB:invyieldABFDHigh:invyieldABFDLow:RABCharm:RABCharmFDHigh:RABCharmFDLow:RABBeauty:fc");
+ ntupleRAB = new TNtuple("ntupleRAB","ntupleRAB (Nb)","pt:TAB:sigmaPP:sigmaAB:invyieldAB:invyieldABFDHigh:invyieldABFDLow:RABCharm:RABCharmFDHigh:RABCharmFDLow:RABBeauty:fc",100000);
} else if (fdMethod==kfc) {
- ntupleRAB = new TNtuple("ntupleRAB","ntupleRAB (fc)","pt:TAB:sigmaPP:sigmaAB:invyieldAB:invyieldABFDHigh:invyieldABFDLow:Rcb:RABCharm:RABCharmFDHigh:RABCharmFDLow:RABBeauty:RABBeautyFDHigh:RABBeautyFDLow:fc");
+ ntupleRAB = new TNtuple("ntupleRAB","ntupleRAB (fc)","pt:TAB:sigmaPP:sigmaAB:invyieldAB:invyieldABFDHigh:invyieldABFDLow:Rcb:RABCharm:RABCharmFDHigh:RABCharmFDLow:RABBeauty:RABBeautyFDHigh:RABBeautyFDLow:fc",100000);
}
if(!ntupleRAB) printf("ERROR: Wrong method option");
TGraphAsymmErrors *gRAB_GlobalSystematics = new TGraphAsymmErrors(nbins+1);
gRAB_GlobalSystematics->SetNameTitle("gRAB_GlobalSystematics","RAB Measurement global (data, FD, Eloss) systematics");
Double_t ElossMax[nbins+1], ElossMin[nbins+1];
- for(Int_t i=0; i<=nbins; i++) { ElossMax[i]=0.; ElossMin[i]=1.; }
+ for(Int_t i=0; i<=nbins; i++) { ElossMax[i]=0.; ElossMin[i]=6.; }
Double_t fcElossMax[nbins+1], fcElossMin[nbins+1];
- for(Int_t i=0; i<=nbins; i++) { fcElossMax[i]=0.; fcElossMin[i]=1.; }
+ for(Int_t i=0; i<=nbins; i++) { fcElossMax[i]=0.; fcElossMin[i]=6.; }
Double_t FDElossMax[nbins+1], FDElossMin[nbins+1];
- for(Int_t i=0; i<=nbins; i++) { FDElossMax[i]=0.; FDElossMin[i]=1.; }
+ for(Int_t i=0; i<=nbins; i++) { FDElossMax[i]=0.; FDElossMin[i]=6.; }
TGraphAsymmErrors *gRAB_Norm = new TGraphAsymmErrors(1);
gRAB_Norm->SetNameTitle("gRAB_Norm","RAB Normalization systematics (pp norm + Tab)");
Double_t normUnc = TMath::Sqrt ( NormPPUnc*NormPPUnc + (TabSyst/Tab)*(TabSyst/Tab) );
+ if(!isUseTaaForRaa) normUnc = TMath::Sqrt ( NormPPUnc*NormPPUnc + NormABUnc*NormABUnc );
gRAB_Norm->SetPoint(1,0.5,1.);
gRAB_Norm->SetPointError(1,0.25,0.25,normUnc,normUnc);
nSigmaAB->GetEntry(ientry);
// cout << " pt="<< pt<<" sigma-AB="<<sigmaAB<<endl;
if ( !(sigmaAB>0.) ) continue;
+ //if(decay==2 && pt<2.) continue;
// Compute RAB and the statistical uncertainty
Int_t hppbin = hSigmaPP->FindBin( pt );
+ Int_t hABbin = hSigmaAB->FindBin( pt );
Double_t sigmapp = hSigmaPP->GetBinContent( hppbin );
// cout << " pt="<< pt<<", sigma-pp="<< sigmapp<<endl;
if (isRaavsEP>0.) sigmapp = 0.5*sigmapp;
if ( !(sigmapp>0.) ) continue;
+
RaaCharm = ( sigmaAB / sigmaABCINT1B ) / ((Tab*1e3) * sigmapp *1e-12 ) ;
+ if(!isUseTaaForRaa) {
+ RaaCharm = ( sigmaAB ) / ( (A*B) * sigmapp ) ;
+ }
+
if (fdMethod==kNb) {
RaaBeauty = Rb ;
}
if (fdMethod==kfc) { ElossHypo = 1. / Rcb; }
else { ElossHypo = 1. / (RaaCharm / RaaBeauty) ; }
if(isRbHypo) ElossHypo = RaaBeauty;
+
+ // If using shadowing hypothesis, change the central hypothesis too
+ if(isShadHypothesis) CentralHypo = centralRbcShad[hABbin];
+
// cout <<" pt "<< pt << " Raa charm " << RaaCharm << " Raa beauty " << RaaBeauty << " eloss hypo "<< ElossHypo<<endl;
//
// Find the bin for the central Eloss hypo
//
if( TMath::Abs( ElossHypo - CentralHypo ) < 0.075 ){
- Int_t hABbin = hSigmaAB->FindBin( pt );
Double_t DeltaIni = TMath::Abs( ElossCentral[ hABbin ] - CentralHypo );
Double_t DeltaV = TMath::Abs( ElossHypo - CentralHypo );
// cout << " pt " << pt << " ECentral " << ElossCentral[ hABbin ] << " Ehypo "<< ElossHypo ;
Double_t sigmapp = hSigmaPP->GetBinContent( hppbin );
if (isRaavsEP>0.) sigmapp = 0.5*sigmapp;
if ( !(sigmapp>0.) ) continue;
+
RaaCharm = ( sigmaAB / sigmaABCINT1B ) / ((Tab*1e3) * sigmapp *1e-12 );
+ if(!isUseTaaForRaa) {
+ RaaCharm = ( sigmaAB ) / ( (A*B) * sigmapp ) ;
+ }
// Flag to know if it is an scaled or extrapolated point of the pp reference
Bool_t isExtrapolatedBin = kFALSE;
RaaCharmFDhigh = ( sigmaABMax / sigmaABCINT1B ) / ((Tab*1e3) * (sigmapp+yPPh) *1e-12 ) ;
RaaCharmFDlow = ( sigmaABMin / sigmaABCINT1B ) / ((Tab*1e3) * (sigmapp-yPPl) *1e-12 ) ;
- if(printout) cout << endl<<" pt "<< pt << " Raa " << RaaCharm <<" high "<< RaaCharmFDhigh << " low "<< RaaCharmFDlow<<endl;
+ if(printout && TMath::Abs(ptprintout-pt)<0.1 ) cout << endl<<" pt "<< pt << " Raa " << RaaCharm <<" high "<< RaaCharmFDhigh << " low "<< RaaCharmFDlow<<endl;
+ if(!isUseTaaForRaa) {
+ RaaCharmFDhigh = ( sigmaABMax ) / ( (A*B)* (sigmapp+yPPh) ) ;
+ RaaCharmFDlow = ( sigmaABMin ) / ( (A*B)* (sigmapp-yPPl) ) ;
+ }
if (fdMethod==kNb) {
hRABvsRb->Fill( pt, RaaCharm, RaaBeauty );
hRABvsRbFDlow->Fill( pt, RaaCharmFDlow, RaaBeautyFDlow );
hRABvsRbFDhigh->Fill( pt, RaaCharmFDhigh, RaaBeautyFDhigh );
- if(printout) cout << " pt "<< pt << " Rb " << RaaBeauty <<" high "<< RaaBeautyFDhigh << " low "<< RaaBeautyFDlow <<endl;
+ if(printout && TMath::Abs(ptprintout-pt)<0.1) cout << " pt "<< pt << " Rb " << RaaBeauty <<" high "<< RaaBeautyFDhigh << " low "<< RaaBeautyFDlow <<endl;
hRABCharmVsRBeautyVsPt->Fill( pt, RaaBeauty, RaaCharm );
Int_t ptbin = hRABvsPt->FindBin( pt );
Int_t hABbin = hMassAB->FindBin( pt );
- if(printout)
- if ( fdMethod==kNb && TMath::Abs(Rb -1.0)< 0.05) {
+ if(isShadHypothesis) CentralHypo = centralRbcShad[hABbin];
+
+ if(printout && TMath::Abs(ptprintout-pt)<0.1)
+ if ( fdMethod==kNb && TMath::Abs(Rb -CentralHypo)< 0.05) {
cout << " pt "<< pt <<", at bin "<<hABbin<<endl;
cout<<" entries "<<entries<<", i="<<ientry<<", pt="<<pt<<", Rb="<<Rb<<", Tab="<<Tab<<", sigmaAB="<<sigmaAB<<", sigmapp="<<sigmapp<<", Raacharm="<<RaaCharm<<", RaaBeauty="<<RaaBeauty<<endl;
cout << " AB basis: mass "<< hMassAB->GetBinContent(hABbin)<<", eff "<< hDirectEffptAB->GetBinContent(hABbin)<<endl;
cout<<" FD low, err low AB "<< (sigmaAB-sigmaABMin)<<" err low PP "<< yPPl<<" Raacharm="<<RaaCharmFDlow<<", RaaBeauty="<<RaaBeautyFDlow<<endl;
cout<<" FD high, err high AB "<< (sigmaABMax-sigmaAB)<<" err high PP "<< yPPh<<" Raacharm="<<RaaCharmFDhigh<<", RaaBeauty="<<RaaBeautyFDhigh<<endl;
}
- if(printout)
- if ( fdMethod==kfc) if(TMath::Abs(Rcb -1.0)< 0.05 ){
+ if(printout && TMath::Abs(ptprintout-pt)<0.1)
+ if ( fdMethod==kfc) if(TMath::Abs(Rcb -CentralHypo)< 0.05 ){
cout << " pt "<< pt <<", at bin "<<hABbin<<endl;
cout<<" entries "<<entries<<", i="<<ientry<<", pt="<<pt<<", Rcb="<<Rcb<<", Tab="<<Tab<<", sigmaAB="<<sigmaAB<<", sigmapp="<<sigmapp<<", Raacharm="<<RaaCharm<<", RaaBeauty="<<RaaBeauty<<endl;
cout << " AB basis: mass "<< hMassAB->GetBinContent(hABbin)<<", eff "<< hDirectEffptAB->GetBinContent(hABbin)<<", fc "<<histofcAB->GetBinContent(hABbin)<< endl;
if(isRbHypo) ElossHypo = RaaBeauty;
hRCharmVsElossHypo[ptbin]->Fill( ElossHypo, RaaCharm );
- // cout <<" pt "<< pt << " Raa charm " << RaaCharm << " Raa beauty " << RaaBeauty << " eloss hypo "<< ElossHypo<<endl;
+ // If using shadowing hypothesis, change the limit hypothesis too
+ if(isShadHypothesis) {
+ MinHypo = minRbcShad[ hABbin ];
+ MaxHypo = maxRbcShad[ hABbin ];
+ }
+
+ // cout <<" pt "<< pt << " Raa charm " << RaaCharm << " Raa beauty " << RaaBeauty << " eloss hypo "<< ElossHypo
+ if(ientry==0) cout<<" pt"<< pt<< " ElossCentral "<< ElossCentral[hABbin] << " min-hypo "<<MinHypo << " max-hypo "<<MaxHypo<<endl;
+
//
// Fill in histos charm (null Eloss hypothesis)
//
hYieldABvsPt->SetBinContent( hRABbin, sigmaAB/sigmaABCINT1B );
hYieldABvsPt->SetBinError( hRABbin, statUncSigmaAB/sigmaABCINT1B );
- if(printout)
+ cout << "pt="<< pt<< " Raa " << RaaCharm << " stat unc. "<<
+ " sigma-pp "<< sigmapp <<" sigma-AB "<< sigmaAB<<endl;
+ if(printout && TMath::Abs(ptprintout-pt)<0.1) {
cout << " Raa " << RaaCharm << " stat unc. "<< stat << " is "<< stat/RaaCharm * 100. <<
"%, stat-pp "<< sigmappStat/sigmapp*100. <<"% stat-AB "<< statUncSigmaAB/sigmaAB*100.<<"%"<<endl;
+ }
Double_t errstatEff = fhStatUncEffcSigmaAB->GetBinError( hRABbin );
fhStatUncEffcSigmaAB_Raa->SetBinError( hRABbin, errstatEff*RaaCharm );
}
- if(printout) {
+ if(printout && TMath::Abs(ptprintout-pt)<0.1) {
cout << " pt : "<< pt<<" Syst-pp-data "<< dataPPUp/sigmapp << "%, ";
if(!isExtrapolatedBin){
if (isRaavsEP>0.) cout <<" extr unc + "<< 0.5*gSigmaPPSystTheory->GetErrorYhigh(istartPPextr)/sigmapp <<" - "<< 0.5*gSigmaPPSystTheory->GetErrorYlow(istartPPextr)/sigmapp <<" %";
FDL = RaaCharmFDhigh; FDH = RaaCharmFDlow;
}
- if(printout) cout<<" Raa "<<RaaCharm<<", Raa-fd-low "<<RaaCharmFDlow <<", Raa-fd-high "<<RaaCharmFDhigh <<endl;
+ if(printout && TMath::Abs(ptprintout-pt)<0.1) cout<<" Raa "<<RaaCharm<<", Raa-fd-low "<<RaaCharmFDlow <<", Raa-fd-high "<<RaaCharmFDhigh <<endl;
maxFdSyst = TMath::Abs(FDH - RaaCharm);
minFdSyst = TMath::Abs(RaaCharm - FDL);
if ( RaaCharm>0 ) {
//
// Filling Eloss scenarii information
//
+ // trick in case not fine enough Rb hypothesis to cope with the min/max range
+ // if( RaaCharm>0 && ( (ElossHypo >= MinHypo && ElossHypo <=MaxHypo) || ElossHypo == ElossCentral[ hABbin ] ) && RaaBeauty<=MaxRb ) {
+ // by default better not use it, to monitor when this happens (could affect results)
if( RaaCharm>0 && ElossHypo >= MinHypo && ElossHypo <=MaxHypo && RaaBeauty<=MaxRb ) {
+
Double_t Ehigh = ElossMax[ hABbin ] ;
Double_t Elow = ElossMin[ hABbin ] ;
if ( RaaCharm > Ehigh ) ElossMax[ hABbin ] = RaaCharm ;
if ( RaaCharm < Elow ) ElossMin[ hABbin ] = RaaCharm ;
- if(printout) {
+ if(printout && TMath::Abs(ptprintout-pt)<0.1) {
cout<<" Hypothesis " << ElossHypo << " sigma-AB "<< sigmaAB <<", Raa "<< RaaCharm <<", Raa Eloss max "<< ElossMax[hABbin] <<" Raa Eloss min "<< ElossMin[hABbin] << " Rb="<< RaaBeauty <<endl;
cout<<" Rb="<< RaaBeauty <<" max "<< RaaBeautyFDhigh <<" min "<< RaaBeautyFDlow <<endl;
}
Double_t RFDlow = RaaCharmFDlow<RaaCharmFDhigh ? RaaCharmFDlow : RaaCharmFDhigh;
if ( RFDhigh > FDEhigh ) FDElossMax[ hABbin ] = RFDhigh ;
if ( RFDlow < FDEmin ) FDElossMin[ hABbin ] = RFDlow ;
- if(printout)
+ if(printout && TMath::Abs(ptprintout-pt)<0.1)
cout<<" Hypothesis " << ElossHypo << " sigma-AB "<< sigmaAB <<", Raa FD-max Eloss max "<< FDElossMax[hABbin] <<" Raa FD-min Eloss min "<< FDElossMin[hABbin] <<endl;
}
RaaPlot->SetBottomMargin(0.1);
RaaPlot->SetTickx();
RaaPlot->SetTicky();
- TH2D *hRaaCanvas = new TH2D("hRaaCanvas"," R_{AB}(c) vs p_{T} (no Eloss hypothesis); p_{t} [GeV/c] ; R_{AA} prompt D",40,0.,40.,100,0.,2.0);
+ TH2D *hRaaCanvas = new TH2D("hRaaCanvas"," R_{AB}(c) vs p_{T} (no Eloss hypothesis); p_{t} [GeV/c] ; R_{AA} prompt D",40,0.,40.,100,0.,3.0);
hRaaCanvas->GetXaxis()->SetTitleSize(0.05);
hRaaCanvas->GetXaxis()->SetTitleOffset(0.9);
hRaaCanvas->GetYaxis()->SetTitleSize(0.05);
legrcb->AddEntry(gRAB_FeedDownSystematics,"Syst. from B feed-down","f");
legrcb->Draw();
TLatex* tc;
- if(decay==1) tc =new TLatex(0.18,0.82,"D^{0}, Pb-Pb #sqrt{s_{NN}}=2.76 TeV");
- else if(decay==2) tc =new TLatex(0.18,0.82,"D^{+}, Pb-Pb #sqrt{s_{NN}}=2.76 TeV");
- else if(decay==3) tc =new TLatex(0.18,0.82,"D^{*+}, Pb-Pb #sqrt{s_{NN}}=2.76 TeV");
- else if(decay==4) tc =new TLatex(0.18,0.82,"D_{s}^{+}, Pb-Pb #sqrt{s_{NN}}=2.76 TeV");
- else tc =new TLatex(0.18,0.82,"any (?) D meson, Pb-Pb #sqrt{s_{NN}}=2.76 TeV");
+ TString system = "Pb-Pb #sqrt{s_{NN}}=2.76 TeV";
+ if( cc==kpPb0100 || cc==kpPb020 || cc==kpPb2040 || cc==kpPb4060 || cc==kpPb60100 ) system = "p-Pb #sqrt{s_{NN}}=5.023 TeV";
+ if(decay==1) tc =new TLatex(0.18,0.82,Form("D^{0}, %s ",system.Data()));
+ else if(decay==2) tc =new TLatex(0.18,0.82,Form("D^{+}, %s ",system.Data()));
+ else if(decay==3) tc =new TLatex(0.18,0.82,Form("D^{*+}, %s ",system.Data()));
+ else if(decay==4) tc =new TLatex(0.18,0.82,Form("D_{s}^{+}, %s ",system.Data()));
+ else tc =new TLatex(0.18,0.82,Form("any (?) D meson, %s ",system.Data()));
tc->SetNDC();
tc->SetTextSize(0.038);
tc->SetTextFont(42);
Double_t pidunc=0.;
err = syst->GetTotalSystErr(pt)*syst->GetTotalSystErr(pt);
+ errUp = err ;
errDown = err ;
- if( cc==k07half || cc==k020 || cc==k010 || cc==k1020 || cc==k2040 ) {
- if(pt<6) pidunc = 0.15;
- else pidunc = 0.05;
- errUp = err + pidunc*pidunc - syst->GetPIDEffErr(pt)*syst->GetPIDEffErr(pt);
- isOk = true;
- }
- else if ( cc==k3050 || cc==k4080 || cc==k4060 || cc==k6080 ){
- if(pt<3.1) pidunc = 0.10;
- else pidunc = 0.05;
- errUp = err + pidunc*pidunc - syst->GetPIDEffErr(pt)*syst->GetPIDEffErr(pt);
- isOk = true;
+ // Apply an asymetric PID uncertainty for 2010 PbPb data only
+ if( syst->GetRunNumber()==10 && syst->GetCollisionType()==1 ) {
+ if( cc==k07half || cc==k020 || cc==k010 || cc==k1020 || cc==k2040 ) {
+ if(pt<6) pidunc = 0.15;
+ else pidunc = 0.05;
+ errUp = err + pidunc*pidunc - syst->GetPIDEffErr(pt)*syst->GetPIDEffErr(pt);
+ isOk = true;
+ }
+ else if ( cc==k3050 || cc==k4080 || cc==k4060 || cc==k6080 ){
+ if(pt<3.1) pidunc = 0.10;
+ else pidunc = 0.05;
+ errUp = err + pidunc*pidunc - syst->GetPIDEffErr(pt)*syst->GetPIDEffErr(pt);
+ isOk = true;
+ }
}
+ else { isOk = true; }
dataSystUp = TMath::Sqrt(errUp);
dataSystDown = TMath::Sqrt(errDown);
- /// << ------------------------- PATCH FOR 2011 LHC11h data ------------------------------
- dataSystUp = dataSystDown;
-
- if ( cc==kpPb0100 ){
- dataSystUp = TMath::Sqrt(err);
- dataSystDown = TMath::Sqrt(err);
- isOk = true;
- }
-
return isOk;
}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff