Back to home page

sPhenix code displayed by LXR
 
 

    

File indexing completed on 2025-08-03 08:15:28 

0001 double CBFunction(double *x, double *p)
0002 {
0003   double norm = p[0];
0004   double alpha = p[1];  // tail start
0005   double n = p[2];      // tail shape 
0006   double mu = p[3];     // upsilon mass
0007   double sigma = p[4];  // upsilon width
0008 
0009   double A = pow(n/fabs(alpha),n)*TMath::Exp(-pow(fabs(alpha),2)/2.);
0010   double B = n/fabs(alpha) - fabs(alpha);
0011   double k = (x[0]-mu)/sigma;
0012 
0013   double val;
0014   if( k > -alpha )
0015     val = norm*TMath::Exp(-0.5*pow(k,2));
0016   else
0017     val = norm*A*pow(B-k,-n);
0018 
0019   if( TMath::IsNaN(val) ) val = 0.0;
0020 
0021   return val;
0022 }
0023 
0024 double TripleCBFunction(double *x, double *p)
0025 {
0026   double norm1 = p[0];  // amplitude of Y(1S) peak
0027   double alpha = p[1];  // tail start
0028   double n = p[2];      // tail shape 
0029   double mu1 = p[3];    // upsilon (1S) mass
0030   double sigma = p[4];  // upsilon width
0031   double norm2 = p[5];
0032   double norm3 = p[6];
0033   double mu2 = mu1*1.0595;
0034   double mu3 = mu1*1.0946;
0035 
0036   double A = pow(n/fabs(alpha),n)*TMath::Exp(-pow(fabs(alpha),2)/2.);
0037   double B = n/fabs(alpha) - fabs(alpha);
0038   double k1 = (x[0]-mu1)/sigma;
0039   double k2 = (x[0]-mu2)/sigma;
0040   double k3 = (x[0]-mu3)/sigma;
0041 
0042   double val,val1,val2,val3;
0043 
0044   if( k1 > -alpha ) { val1 = norm1*TMath::Exp(-0.5*pow(k1,2)); }
0045   else { val1 = norm1*A*pow(B-k1,-n); }
0046   if( k2 > -alpha ) { val2 = norm2*TMath::Exp(-0.5*pow(k2,2)); }
0047   else { val2 = norm2*A*pow(B-k2,-n); }
0048   if( k3 > -alpha ) { val3 = norm3*TMath::Exp(-0.5*pow(k3,2)); }
0049   else { val3 = norm3*A*pow(B-k3,-n); }
0050 
0051   val = val1 + val2 + val3;
0052 
0053   if( TMath::IsNaN(val) ) val = 0.0;
0054 
0055   return val;
0056 }
0057 
0058 double SandB_CBFunction(double *x, double *p)
0059 {
0060   double norm1 = p[0];  // amplitude of Y(1S) peak
0061   double alpha = p[1];  // tail start
0062   double n = p[2];      // tail shape 
0063   double mu1 = p[3];     // upsilon (1S) mass
0064   double sigma = p[4];  // upsilon width
0065   double norm2 = p[5];
0066   double norm3 = p[6];
0067   double mu2 = mu1*1.0595;
0068   double mu3 = mu1*1.0946;
0069 
0070   double A = pow(n/fabs(alpha),n)*TMath::Exp(-pow(fabs(alpha),2)/2.);
0071   double B = n/fabs(alpha) - fabs(alpha);
0072   double k1 = (x[0]-mu1)/sigma;
0073   double k2 = (x[0]-mu2)/sigma;
0074   double k3 = (x[0]-mu3)/sigma;
0075 
0076   double val,val1,val2,val3;
0077 
0078   if( k1 > -alpha ) { val1 = norm1*TMath::Exp(-0.5*pow(k1,2)); }
0079   else { val1 = norm1*A*pow(B-k1,-n); }
0080   if( k2 > -alpha ) { val2 = norm2*TMath::Exp(-0.5*pow(k2,2)); }
0081   else { val2 = norm2*A*pow(B-k2,-n); }
0082   if( k3 > -alpha ) { val3 = norm3*TMath::Exp(-0.5*pow(k3,2)); }
0083   else { val3 = norm3*A*pow(B-k3,-n); }
0084 
0085   double bgnorm1 = p[7];
0086   double bgslope1 = p[8];
0087 
0088   double bg = exp(bgnorm1+x[0]*bgslope1);
0089 
0090   val = val1 + val2 + val3 + bg;
0091   if( TMath::IsNaN(val) ) val = 0.0;
0092 
0093   return val;
0094 }
0095 
0096 //----------------------------------------------------------------
0097 //----------------------------------------------------------------
0098 //----------------------------------------------------------------
0099 
0100 void newplotbg_newsupp_2022() {
0101 
0102 //gROOT->LoadMacro("sPHENIXStyle/sPhenixStyle.C");
0103 //SetsPhenixStyle();
0104 
0105 gStyle->SetOptStat(0);
0106 gStyle->SetOptFit(0);
0107 
0108 TRandom* myrandom = new TRandom3();
0109 const int nbins = 15;
0110 double eideff = 0.9;
0111 string times = "*";
0112 TLatex* tl[15];
0113 char tlchar[999];
0114 
0115 double statscale;
0116 statscale = 1.4; // from 10B to 14B
0117 /*
0118 //double statscalepp;
0119 //statscalepp = 0.3265; // from 7350B in the sPHENIX proposal to 2400B
0120 //statscalepp = 0.3947; // from 13.49k in the sPHENIX proposal to 5.32k (Marzia)
0121 //double statscale = 2.4;  // 5-year plan all correladed bg plots are for 10B
0122 //double statscale = 1.43;  // 3-year plan all correladed bg plots are for 10B
0123 //double statscalepp = 1.78; // old wrong???
0124 //double statscalepp = 1.13; // 5-year plan
0125 */
0126 double statscale_lowlim = 7.0;
0127 double statscale_uplim = 14.0;
0128 
0129   TF1* fCBpp = new TF1("fCBpp",CBFunction,5.,14.,5);
0130   TF1* fCBauau = new TF1("fCBauau",CBFunction,5.,14.,5);
0131   TF1* fCB1s = new TF1("fCB1s",CBFunction,5.,14.,5);
0132   TF1* fCB2s = new TF1("fCB2s",CBFunction,5.,14.,5);
0133   TF1* fCB3s = new TF1("fCB3s",CBFunction,5.,14.,5);
0134   TF1* fTCB = new TF1("fTCB",TripleCBFunction,5.,14.,7);
0135   TF1* fTCBpp = new TF1("fTCBpp",TripleCBFunction,5.,14.,7);
0136   TF1* fTCBauau = new TF1("fTCBauau",TripleCBFunction,5.,14.,7);
0137   TF1* fSandB = new TF1("fSandB",SandB_CBFunction,5.,14.,9);
0138   TF1* fSandBfordave = new TF1("fSandBfordave",SandB_CBFunction,5.,14.,9);
0139   TF1* fSandBpp = new TF1("fSandBpp",SandB_CBFunction,5.,14.,9);
0140   TF1* fSandBauau = new TF1("fSandBauau",SandB_CBFunction,5.,14.,9);
0141 
0142 //---------------------------------------------------------
0143 // Upsilons
0144 //---------------------------------------------------------
0145 
0146 double raapt[9],raa1s[9],raa2s[9],raa3s[9];
0147 raapt[0] = 1.5;
0148 raapt[1] = 4.5;
0149 raapt[2] = 7.5;
0150 raapt[3] = 10.5;
0151 raapt[4] = 13.5;
0152 raa1s[0] = 0.535;
0153 raa1s[1] = 0.535;
0154 raa1s[2] = 0.535;
0155 raa1s[3] = 0.535;
0156 raa1s[4] = 0.535;
0157 raa2s[0] = 0.170;
0158 raa2s[1] = 0.170;
0159 raa2s[2] = 0.170;
0160 raa2s[3] = 0.170;
0161 raa2s[4] = 0.170;
0162 raa3s[0] = 0.085;
0163 raa3s[1] = 0.085;
0164 raa3s[2] = 0.085;
0165 raa3s[3] = 0.085;
0166 raa3s[4] = 0.085;
0167 
0168 TGraph* grRAA1S = new TGraph(5,raapt,raa1s);
0169 TGraph* grRAA2S = new TGraph(5,raapt,raa2s);
0170 TGraph* grRAA3S = new TGraph(5,raapt,raa3s);
0171 
0172 int nchan=400;
0173 double start=0.0;
0174 double stop=20.0;
0175 
0176 string str_UpsilonPt  = "(2.0*3.14159*x*[0]*pow((1 + x*x/(4*[1]) ),-[2]))";   // dN/dpT
0177 string str_UpsilonXPt = "(2.0*3.14159*x*x*[0]*pow((1 + x*x/(4*[1]) ),-[2]))"; // need this for mean pT calculation
0178 TF1* fUpsilonPt = new TF1("fUpsilonPt",str_UpsilonPt.c_str(),0.,20.);
0179 TF1* fUpsilonXPt = new TF1("fUpsilonXPt",str_UpsilonXPt.c_str(),0.,20.);
0180 fUpsilonPt->SetParameters(72.1, 26.516, 10.6834);
0181 fUpsilonXPt->SetParameters(72.1, 26.516, 10.6834);
0182 double upsnorm = fUpsilonPt->Integral(0.,20.);
0183 
0184 double frac[3];  // upsilons fractions
0185   frac[0] = 0.7117;
0186   frac[1] = 0.1851;
0187   frac[2] = 0.1032;
0188 double scale[3]; // mass scaling
0189   scale[0] = 1.0;
0190   scale[1] = 1.0595;
0191   scale[2] = 1.0946;
0192 
0193 //---------- Use fixed numbers for 2022 estimate --------------------------------------
0194 // Numbers from Marzia for 140 B MB Au+Au and 2400 B p+p events
0195 // RAA for Y(3S) is 8.5%
0196   int Nups1   = 7011;
0197   int Nups2   = 562;
0198   int Nups3   = 156;
0199   int Nups1pp = 2.86e+03;
0200   int Nups2pp = 7.16e+02;
0201   int Nups3pp = 3.98e+02;
0202 
0203 //-------------------------------------------------------------------------------------
0204 
0205 //====================================================================================
0206 //====================================================================================
0207 //====================================================================================
0208 
0209 //
0210 // these histograms (hhups*) are randomly generated using the fit above to single 
0211 // peak (fCB function) and used for the RAA uncertainty calculation
0212 //
0213 // FORCE specific RESOLUTION and tail shape
0214   double tonypar1 = 0.98;  // Tony's 100 pion simulation April 2019
0215   double tonypar2 = 0.93; // Tony's 100 pion simulation April 2019
0216   //double tonypar3 = 9.437;  // Tony's 100 pion simulation April 2019
0217   double tonypar3 = 9.448;  // benchmark
0218   double tonypar4 = 0.100;  // benchmark
0219   double tonypar4pp = 0.089; // Tony's 100 pion simulation April 2019
0220   fCBpp->SetParameter(0,1000.); 
0221   fCBpp->SetParameter(1,tonypar1); 
0222   fCBpp->SetParameter(2,tonypar2); 
0223   fCBpp->SetParameter(3,tonypar3); 
0224   fCBpp->SetParameter(4,tonypar4pp);  
0225   fCBauau->SetParameter(0,1000.); 
0226   fCBauau->SetParameter(1,tonypar1); 
0227   fCBauau->SetParameter(2,tonypar2); 
0228   fCBauau->SetParameter(3,tonypar3); 
0229   fCBauau->SetParameter(4,tonypar4);  
0230 
0231 
0232 char hhname[999];
0233 TH1D* hhups[nbins+1];
0234 TH1D* hhups1[nbins+1];
0235 TH1D* hhups2[nbins+1];
0236 TH1D* hhups3[nbins+1];
0237 TH1D* hhupspp[nbins+1];
0238 TH1D* hhups1pp[nbins+1];
0239 TH1D* hhups2pp[nbins+1];
0240 TH1D* hhups3pp[nbins+1];
0241 for(int i=0; i<nbins+1; i++) {
0242   sprintf(hhname,"hhups_%d",i);
0243   hhups[i] = new TH1D(hhname,"",nchan,start,stop);
0244   hhups[i]->Sumw2();
0245   sprintf(hhname,"hhups1_%d",i);
0246   hhups1[i] = new TH1D(hhname,"",nchan,start,stop);
0247   hhups1[i]->Sumw2();
0248   sprintf(hhname,"hhups2_%d",i);
0249   hhups2[i] = new TH1D(hhname,"",nchan,start,stop);
0250   hhups2[i]->Sumw2();
0251   sprintf(hhname,"hhups3_%d",i);
0252   hhups3[i] = new TH1D(hhname,"",nchan,start,stop);
0253   hhups3[i]->Sumw2();
0254   hhups[i]->SetLineWidth(2);
0255   hhups1[i]->SetLineWidth(2);
0256   hhups2[i]->SetLineWidth(2);
0257   hhups3[i]->SetLineWidth(2);
0258   sprintf(hhname,"hhupspp_%d",i);
0259   hhupspp[i] = new TH1D(hhname,"",nchan,start,stop);
0260   hhupspp[i]->Sumw2();
0261   sprintf(hhname,"hhups1pp_%d",i);
0262   hhups1pp[i] = new TH1D(hhname,"",nchan,start,stop);
0263   hhups1pp[i]->Sumw2();
0264   sprintf(hhname,"hhups2pp_%d",i);
0265   hhups2pp[i] = new TH1D(hhname,"",nchan,start,stop);
0266   hhups2pp[i]->Sumw2();
0267   sprintf(hhname,"hhups3pp_%d",i);
0268   hhups3pp[i] = new TH1D(hhname,"",nchan,start,stop);
0269   hhups3pp[i]->Sumw2();
0270   hhupspp[i]->SetLineWidth(2);
0271   hhups1pp[i]->SetLineWidth(2);
0272   hhups2pp[i]->SetLineWidth(2);
0273   hhups3pp[i]->SetLineWidth(2);
0274 }
0275 
0276 for(int j=0; j<nbins; j++) {
0277 
0278   double s1 = j*1.0;
0279   double s2 = s1 + 1.0;
0280   //double tmpnups1 = Nups1nosup*fUpsilonPt->Integral(s1,s2)/upsnorm * grRAA1S->Eval((s1+s2)/2.);
0281   //double tmpnups2 = Nups2nosup*fUpsilonPt->Integral(s1,s2)/upsnorm * grRAA2S->Eval((s1+s2)/2.);
0282   //double tmpnups3 = Nups3nosup*fUpsilonPt->Integral(s1,s2)/upsnorm * grRAA3S->Eval((s1+s2)/2.);
0283   double tmpnups1 = Nups1*fUpsilonPt->Integral(s1,s2)/upsnorm;
0284   double tmpnups2 = Nups2*fUpsilonPt->Integral(s1,s2)/upsnorm;
0285   double tmpnups3 = Nups3*fUpsilonPt->Integral(s1,s2)/upsnorm;
0286   fCBauau->SetParameter(3,tonypar3);
0287     for(int i=0; i<int(tmpnups1+0.5); i++) { double myrnd = fCBauau->GetRandom(); hhups1[j]->Fill(myrnd); hhups[j]->Fill(myrnd); }
0288   fCBauau->SetParameter(3,tonypar3*scale[1]);
0289     for(int i=0; i<int(tmpnups2+0.5); i++) { double myrnd = fCBauau->GetRandom(); hhups2[j]->Fill(myrnd); hhups[j]->Fill(myrnd); }
0290   fCBauau->SetParameter(3,tonypar3*scale[2]);
0291     for(int i=0; i<int(tmpnups3+0.5); i++) { double myrnd = fCBauau->GetRandom(); hhups3[j]->Fill(myrnd); hhups[j]->Fill(myrnd); }
0292 
0293   double tmpnups1pp = Nups1pp*fUpsilonPt->Integral(s1,s2)/upsnorm;
0294   double tmpnups2pp = Nups2pp*fUpsilonPt->Integral(s1,s2)/upsnorm;
0295   double tmpnups3pp = Nups3pp*fUpsilonPt->Integral(s1,s2)/upsnorm;
0296   fCBpp->SetParameter(3,tonypar3);
0297     for(int i=0; i<int(tmpnups1pp+0.5); i++) { double myrnd = fCBpp->GetRandom(); hhups1pp[j]->Fill(myrnd); hhupspp[j]->Fill(myrnd); }
0298   fCBpp->SetParameter(3,tonypar3*scale[1]);
0299     for(int i=0; i<int(tmpnups2pp+0.5); i++) { double myrnd = fCBpp->GetRandom(); hhups2pp[j]->Fill(myrnd); hhupspp[j]->Fill(myrnd); }
0300   fCBpp->SetParameter(3,tonypar3*scale[2]);
0301     for(int i=0; i<int(tmpnups3pp+0.5); i++) { double myrnd = fCBpp->GetRandom(); hhups3pp[j]->Fill(myrnd); hhupspp[j]->Fill(myrnd); }
0302 
0303 } // end loop over pT bins
0304 
0305 // all pT
0306 
0307   fCBpp->SetParameter(3,tonypar3);
0308   fCBauau->SetParameter(3,tonypar3);
0309     for(int i=0; i<int(Nups1+0.5); i++) { double myrnd = fCBauau->GetRandom(); hhups1[nbins]->Fill(myrnd); hhups[nbins]->Fill(myrnd); }
0310     for(int i=0; i<int(Nups1pp+0.5); i++) { double myrnd = fCBpp->GetRandom(); hhups1pp[nbins]->Fill(myrnd); hhupspp[nbins]->Fill(myrnd); }
0311   fCBpp->SetParameter(3,tonypar3*scale[1]);
0312   fCBauau->SetParameter(3,tonypar3*scale[1]);
0313     for(int i=0; i<int(Nups2+0.5); i++) { double myrnd = fCBauau->GetRandom(); hhups2[nbins]->Fill(myrnd); hhups[nbins]->Fill(myrnd); }
0314     for(int i=0; i<int(Nups2pp+0.5); i++) { double myrnd = fCBpp->GetRandom(); hhups2pp[nbins]->Fill(myrnd); hhupspp[nbins]->Fill(myrnd); }
0315   fCBpp->SetParameter(3,tonypar3*scale[2]);
0316   fCBauau->SetParameter(3,tonypar3*scale[2]);
0317     for(int i=0; i<int(Nups3+0.5); i++) { double myrnd = fCBauau->GetRandom(); hhups3[nbins]->Fill(myrnd); hhups[nbins]->Fill(myrnd); }
0318     for(int i=0; i<int(Nups3pp+0.5); i++) { double myrnd = fCBpp->GetRandom(); hhups3pp[nbins]->Fill(myrnd); hhupspp[nbins]->Fill(myrnd); }
0319 
0320 // fit and draw pp no bg
0321 
0322 TCanvas* cupspp = new TCanvas("cupspp","Upsilons in p+p",100,100,600,600);
0323   fTCBpp->SetParameter(0,2000.);
0324   fTCBpp->FixParameter(1,tonypar1);
0325   fTCBpp->FixParameter(2,tonypar2);
0326   fTCBpp->SetParameter(3,tonypar3);
0327   fTCBpp->FixParameter(4,tonypar4); // fix width from the single peak fit
0328   fTCBpp->SetParameter(5,500.);
0329   fTCBpp->SetParameter(6,100.);
0330   hhupspp[nbins]->Fit(fTCBpp,"rl","",7.,11.);
0331   hhupspp[nbins]->SetAxisRange(7.,11.);
0332   hhupspp[nbins]->SetMarkerSize(1.0);
0333   hhupspp[nbins]->GetXaxis()->SetTitle("Invariant mass [GeV/c^{2}]");
0334   hhupspp[nbins]->GetXaxis()->SetTitleOffset(1.0);
0335   double tmpamp1 = hhupspp[nbins]->GetFunction("fTCBpp")->GetParameter(0);
0336   double tmpamp5 = tmpamp1*frac[1]/frac[0]; // correct fit for correct upsilon states ratio
0337   double tmpamp6 = tmpamp1*frac[2]/frac[0];
0338   hhupspp[nbins]->Draw();
0339 
0340   fCB1s->SetLineColor(kBlue);
0341   fCB1s->SetLineWidth(1);
0342     fCB1s->SetParameter(0,fTCBpp->GetParameter(0));
0343     fCB1s->SetParameter(1,fTCBpp->GetParameter(1));
0344     fCB1s->SetParameter(2,fTCBpp->GetParameter(2));
0345     fCB1s->SetParameter(3,fTCBpp->GetParameter(3)*scale[0]);
0346     fCB1s->SetParameter(4,fTCBpp->GetParameter(4));
0347   fCB2s->SetLineColor(kRed);
0348   fCB2s->SetLineWidth(1);
0349     fCB2s->SetParameter(0,tmpamp5);
0350     fCB2s->SetParameter(1,fTCBpp->GetParameter(1));
0351     fCB2s->SetParameter(2,fTCBpp->GetParameter(2));
0352     fCB2s->SetParameter(3,fTCBpp->GetParameter(3)*scale[1]);
0353     fCB2s->SetParameter(4,fTCBpp->GetParameter(4));
0354   fCB3s->SetLineColor(kGreen+2);
0355   fCB3s->SetLineWidth(1);
0356     fCB3s->SetParameter(0,tmpamp6);
0357     fCB3s->SetParameter(1,fTCBpp->GetParameter(1));
0358     fCB3s->SetParameter(2,fTCBpp->GetParameter(2));
0359     fCB3s->SetParameter(3,fTCBpp->GetParameter(3)*scale[2]);
0360     fCB3s->SetParameter(4,fTCBpp->GetParameter(4));
0361   fCB1s->Draw("same");
0362   fCB2s->Draw("same");
0363   fCB3s->Draw("same");
0364 
0365 // fit and draw AuAu no bg
0366 
0367 TCanvas* cupsauau = new TCanvas("cupsauau","Upsilons in Au+Au",100,100,600,600);
0368   fTCBauau->SetParameter(0,2000.);
0369   fTCBauau->FixParameter(1,tonypar1);
0370   fTCBauau->FixParameter(2,tonypar2);
0371   fTCBauau->SetParameter(3,tonypar3);
0372   fTCBauau->FixParameter(4,tonypar4); // fix width from the single peak fit
0373   fTCBauau->SetParameter(5,500.);
0374   fTCBauau->SetParameter(6,100.);
0375   hhups[nbins]->Fit(fTCBauau,"rl","",7.,11.);
0376   hhups[nbins]->SetAxisRange(7.,11.);
0377   hhups[nbins]->SetMarkerSize(1.0);
0378   hhups[nbins]->GetXaxis()->SetTitle("Invariant mass [GeV/c^{2}]");
0379   hhups[nbins]->GetXaxis()->SetTitleOffset(1.0);
0380   tmpamp1 = hhups[nbins]->GetFunction("fTCBauau")->GetParameter(0);
0381   tmpamp5 = tmpamp1*frac[1]/frac[0]; // correct fit for correct upsilon states ratio
0382   tmpamp6 = tmpamp1*frac[2]/frac[0];
0383   hhups[nbins]->Draw();
0384 
0385 TCanvas* cupsvspt = new TCanvas("cupsvspt","Upsilons vs. p_{T}",100,100,1200,900);
0386 cupsvspt->Divide(4,3);
0387 for(int i=0; i<12; i++) {
0388 if(i==0)  {cupsvspt->cd(1);}
0389 if(i==1)  {cupsvspt->cd(2);}
0390 if(i==2)  {cupsvspt->cd(3);}
0391 if(i==3)  {cupsvspt->cd(4);}
0392 if(i==4)  {cupsvspt->cd(5);}
0393 if(i==5)  {cupsvspt->cd(6);}
0394 if(i==6)  {cupsvspt->cd(7);}
0395 if(i==7)  {cupsvspt->cd(8);}
0396 if(i==8)  {cupsvspt->cd(9);}
0397 if(i==9)  {cupsvspt->cd(10);}
0398 if(i==10) {cupsvspt->cd(11);}
0399 if(i==11) {cupsvspt->cd(12);}
0400 hhups[i]->SetAxisRange(7.0,11.0); hhups[i]->Draw();
0401   sprintf(tlchar,"%d-%d GeV",i,i+1);   tl[i] = new TLatex(8.0,hhups[i]->GetMaximum()*0.95,tlchar); tl[i]->Draw();
0402 }
0403 
0404 //----------------------------------------------------
0405 //  Backgrounds
0406 //----------------------------------------------------
0407 
0408 TH1D* hhall[nbins+1];
0409 TH1D* hhall_scaled[nbins+1];
0410 
0411 TH1D* hhtotbg[nbins+1]; 
0412 TH1D* hhtotbg_scaled[nbins+1]; 
0413 TH1D* hhcombbg[nbins+1];
0414 TH1D* hhcombbg_scaled[nbins+1];
0415 TH1D* hhfakefake[nbins+1];
0416 TH1D* hhfakehf[nbins+1];
0417 TH1D* hhbottom[nbins+1];
0418 TH1D* hhcharm[nbins+1];
0419 TH1D* hhdy[nbins+1];
0420 TH1D* hhcorrbg[nbins+1];
0421 TH1D* hhcorrbg_scaled[nbins+1];
0422 TH1D* hhfit[nbins+1];
0423 char tmpname[999];
0424 
0425 //----------------------------------------------------------------------------------------
0426 // Correlated BG calculated for 10B central AuAu events
0427 //----------------------------------------------------------------------------------------
0428 
0429 double corrbgfitpar0;
0430 double corrbgfitpar1;
0431 
0432 TFile* f=new TFile("ccbb_eideff09.root");
0433 for(int i=0; i<nbins+1; i++) {
0434   sprintf(tmpname,"hhbottom_%d",i);
0435   hhbottom[i] = (TH1D*)f->Get(tmpname);
0436   hhbottom[i]->SetDirectory(gROOT);
0437   sprintf(tmpname,"hhcharm_%d",i);
0438   hhcharm[i] = (TH1D*)f->Get(tmpname);
0439   hhcharm[i]->SetDirectory(gROOT);
0440   sprintf(tmpname,"hhdy_%d",i);
0441   hhdy[i] = (TH1D*)f->Get(tmpname);
0442   hhdy[i]->SetDirectory(gROOT);
0443   sprintf(tmpname,"hhcorrbg_%d",i);
0444   hhcorrbg[i] = (TH1D*)hhbottom[i]->Clone(tmpname);
0445   hhcorrbg[i]->Add(hhcharm[i]);
0446   hhcorrbg[i]->Add(hhdy[i]);
0447   sprintf(tmpname,"hhcorrbg_scaled_%d",i);
0448   hhcorrbg_scaled[i] = (TH1D*)hhcorrbg[i]->Clone(tmpname);
0449   hhcorrbg[i]->Fit("expo","rql","",statscale_lowlim,statscale_uplim);
0450   hhbottom[i]->Fit("expo","rql","",statscale_lowlim,statscale_uplim);
0451   hhdy[i]->Fit("expo","rql","",statscale_lowlim,statscale_uplim);
0452     if(i==nbins) {
0453        corrbgfitpar0 = hhcorrbg[i]->GetFunction("expo")->GetParameter(0);
0454        corrbgfitpar1 = hhcorrbg[i]->GetFunction("expo")->GetParameter(1);
0455     }
0456     cout << "bgpar0["<< i <<"]="<<hhcorrbg[i]->GetFunction("expo")->GetParameter(0)+TMath::Log(statscale)<<";"<< endl; 
0457     cout << "bgpar1["<< i <<"]="<<hhcorrbg[i]->GetFunction("expo")->GetParameter(1)<<";"<< endl; 
0458     for(int k=1; k<=hhcorrbg[i]->GetNbinsX(); k++) {
0459       if(hhcorrbg[i]->GetBinLowEdge(k)<statscale_lowlim || (hhcorrbg[i]->GetBinLowEdge(k)+hhcorrbg[i]->GetBinWidth(k))>statscale_uplim) {
0460         hhcorrbg_scaled[i]->SetBinContent(k,0.); 
0461         hhcorrbg_scaled[i]->SetBinError(k,0.);
0462       }
0463       else {
0464         double tmp = statscale * hhcorrbg[i]->GetFunction("expo")->Eval(hhcorrbg[i]->GetBinCenter(k));
0465         double tmprnd = myrandom->Poisson(tmp); 
0466         if(tmprnd<0.) { tmprnd=0.; }
0467         hhcorrbg_scaled[i]->SetBinContent(k,tmprnd); 
0468         hhcorrbg_scaled[i]->SetBinError(k,sqrt(tmprnd));
0469       }
0470     }
0471   hhcorrbg_scaled[i]->Fit("expo","rql","",statscale_lowlim,statscale_uplim);
0472   hhcorrbg[i]->SetDirectory(gROOT);
0473   hhcorrbg_scaled[i]->SetDirectory(gROOT);
0474 }
0475 f->Close();
0476 
0477 TCanvas* c2 = new TCanvas("c2","Correlated BG (all p_{T})",100,100,600,600);
0478   
0479 hhbottom[nbins]->SetAxisRange(7.0,14.0);
0480 hhbottom[nbins]->SetLineColor(kBlue);
0481 hhbottom[nbins]->SetLineWidth(2);
0482 
0483 hhbottom[nbins]->GetXaxis()->SetTitle("Invariant mass [GeV/c^{2}]");
0484 hhbottom[nbins]->GetXaxis()->SetTitleOffset(1.0);
0485 hhbottom[nbins]->GetXaxis()->SetTitleColor(1);
0486 hhbottom[nbins]->GetXaxis()->SetTitleSize(0.040);
0487 hhbottom[nbins]->GetXaxis()->SetLabelSize(0.040);
0488 //hhbottom[nbins]->GetYaxis()->SetTitle("Correlated background");
0489 hhbottom[nbins]->GetYaxis()->SetTitleOffset(1.3);
0490 hhbottom[nbins]->GetYaxis()->SetTitleSize(0.040);
0491 hhbottom[nbins]->GetYaxis()->SetLabelSize(0.040);
0492 
0493 hhcorrbg[nbins]->SetAxisRange(7.0,14.0);
0494 hhcorrbg[nbins]->SetMinimum(0.1);
0495 hhcorrbg[nbins]->SetLineColor(kBlack);
0496 hhcorrbg[nbins]->SetLineWidth(2);
0497 hhcorrbg[nbins]->Draw("hist");
0498 
0499 hhbottom[nbins]->SetMarkerColor(kBlue);
0500 hhbottom[nbins]->SetLineColor(kBlue);
0501 hhbottom[nbins]->Draw("same");
0502 
0503 hhdy[nbins]->SetMarkerColor(kGreen+2);
0504 hhdy[nbins]->SetLineColor(kGreen+2);
0505 hhdy[nbins]->SetLineWidth(2);
0506 hhdy[nbins]->Draw("same");
0507 
0508 hhcharm[nbins]->SetMarkerColor(kRed);
0509 hhcharm[nbins]->SetLineColor(kRed);
0510 hhcharm[nbins]->SetLineWidth(2);
0511 hhcharm[nbins]->Draw("same");
0512 
0513 TCanvas* c0 = new TCanvas("c0","Correlated BG vs. p_{T} 10B events",100,100,1200,900);
0514 c0->Divide(4,3);
0515 for(int i=0; i<nbins; i++) {
0516 if(i==0)  {c0->cd(1);}
0517 if(i==1)  {c0->cd(2);}
0518 if(i==2)  {c0->cd(3);}
0519 if(i==3)  {c0->cd(4);}
0520 if(i==4)  {c0->cd(5);}
0521 if(i==5)  {c0->cd(6);}
0522 if(i==6)  {c0->cd(7);}
0523 if(i==7)  {c0->cd(8);}
0524 if(i==8)  {c0->cd(9);}
0525 if(i==9)  {c0->cd(10);}
0526 if(i==10) {c0->cd(11);}
0527 if(i==11) {c0->cd(12);}
0528   hhcorrbg[i]->SetAxisRange(7.,14.);
0529   hhcorrbg[i]->GetFunction("expo")->SetLineColor(kBlack); 
0530   hhbottom[i]->GetFunction("expo")->SetLineColor(kBlue); 
0531   hhdy[i]->GetFunction("expo")->SetLineColor(kGreen+2); 
0532   hhcorrbg[i]->SetLineColor(kBlack); 
0533   hhbottom[i]->SetLineColor(kBlue); 
0534   hhdy[i]->SetLineColor(kGreen+2); 
0535   hhcorrbg[i]->Draw();
0536   hhbottom[i]->Draw("same");
0537   hhdy[i]->Draw("same");
0538   sprintf(tlchar,"%d-%d GeV",i,i+1);   tl[i] = new TLatex(8.0,hhcorrbg[i]->GetMaximum()*0.95,tlchar); tl[i]->Draw();
0539 }
0540 
0541 TCanvas* c0scaled = new TCanvas("c0scaled","SCALED Correlated BG vs. p_{T}",100,100,1200,900);
0542 c0scaled->Divide(4,3);
0543 for(int i=0; i<nbins; i++) {
0544 if(i==0)  {c0scaled->cd(1);}
0545 if(i==1)  {c0scaled->cd(2);}
0546 if(i==2)  {c0scaled->cd(3);}
0547 if(i==3)  {c0scaled->cd(4);}
0548 if(i==4)  {c0scaled->cd(5);}
0549 if(i==5)  {c0scaled->cd(6);}
0550 if(i==6)  {c0scaled->cd(7);}
0551 if(i==7)  {c0scaled->cd(8);}
0552 if(i==8)  {c0scaled->cd(9);}
0553 if(i==9)  {c0scaled->cd(10);}
0554 if(i==10) {c0scaled->cd(11);}
0555 if(i==11) {c0scaled->cd(12);}
0556 hhcorrbg_scaled[i]->SetAxisRange(7.0,14.0); hhcorrbg_scaled[i]->Draw();
0557   sprintf(tlchar,"%d-%d GeV",i,i+1);   tl[i] = new TLatex(8.0,hhcorrbg_scaled[i]->GetMaximum()*0.95,tlchar); tl[i]->Draw();
0558 }
0559 
0560 //-----------------------------------------------------------------------------------------
0561 //  Correlated bg in p+p
0562 //-----------------------------------------------------------------------------------------
0563 
0564 TH1D* hhbottom_pp[nbins+1];
0565 TH1D* hhdy_pp[nbins+1];
0566 TH1D* hhcorrbg_pp[nbins+1];
0567 TH1D* hhall_pp[nbins+1];
0568 
0569 double ppcorr = (2400./14.)/962.; // 2400B pp and 140B MB AuAu
0570 TF1* fbottom_nosup_corr = new TF1("fbottom_nosup_corr","[0]+[1]*x",5.,14.);
0571 fbottom_nosup_corr->SetParameters(-2.13861, 0.683323);
0572 
0573 for(int i=0; i<nbins+1; i++) {
0574 
0575   sprintf(tmpname,"hhbottom_pp_%d",i);
0576   hhbottom_pp[i] = (TH1D*)hhbottom[i]->Clone(tmpname);
0577   for(int k=1; k<=hhbottom_pp[i]->GetNbinsX(); k++) {
0578     if(hhbottom_pp[i]->GetBinLowEdge(k)<statscale_lowlim || (hhbottom_pp[i]->GetBinLowEdge(k)+hhbottom_pp[i]->GetBinWidth(k))>statscale_uplim) {
0579       hhbottom_pp[i]->SetBinContent(k,0.);
0580       hhbottom_pp[i]->SetBinError(k,0.);
0581     }
0582     else {
0583       double tmp = ppcorr * fbottom_nosup_corr->Eval(hhbottom[i]->GetBinCenter(k)) * hhbottom[i]->GetFunction("expo")->Eval(hhbottom[i]->GetBinCenter(k));
0584       double tmprnd = myrandom->Poisson(tmp);
0585       if(tmprnd<0.) { tmprnd=0.; }
0586       hhbottom_pp[i]->SetBinContent(k,tmprnd);
0587       hhbottom_pp[i]->SetBinError(k,sqrt(tmprnd));
0588     }
0589   }
0590 
0591   sprintf(tmpname,"hhdy_pp_%d",i);
0592   hhdy_pp[i] = (TH1D*)hhdy[i]->Clone(tmpname);
0593   for(int k=1; k<=hhdy_pp[i]->GetNbinsX(); k++) {
0594     if(hhdy_pp[i]->GetBinLowEdge(k)<statscale_lowlim || (hhdy_pp[i]->GetBinLowEdge(k)+hhdy_pp[i]->GetBinWidth(k))>statscale_uplim) {
0595       hhdy_pp[i]->SetBinContent(k,0.);
0596       hhdy_pp[i]->SetBinError(k,0.);
0597     }
0598     else {
0599       double tmp = ppcorr * hhdy[i]->GetFunction("expo")->Eval(hhdy[i]->GetBinCenter(k));
0600       double tmprnd = myrandom->Poisson(tmp);
0601       if(tmprnd<0.) { tmprnd=0.; }
0602       hhdy_pp[i]->SetBinContent(k,tmprnd);
0603       hhdy_pp[i]->SetBinError(k,sqrt(tmprnd));
0604     }
0605   }
0606 
0607   sprintf(tmpname,"hhcorrbg_pp_%d",i);
0608   hhcorrbg_pp[i] = (TH1D*)hhbottom_pp[i]->Clone(tmpname);
0609   hhcorrbg_pp[i]->Add(hhdy_pp[i]);
0610     hhcorrbg_pp[i]->SetMarkerColor(kBlack);
0611     hhcorrbg_pp[i]->SetLineColor(kBlack);
0612     hhbottom_pp[i]->SetLineColor(kBlue);
0613     hhdy_pp[i]->SetLineColor(kGreen+2);
0614   sprintf(tmpname,"hhall_pp_%d",i);
0615   hhall_pp[i] = (TH1D*)hhcorrbg_pp[i]->Clone(tmpname);
0616   hhall_pp[i]->Add(hhupspp[i]);
0617     hhall_pp[i]->SetLineColor(kMagenta);
0618     hhall_pp[i]->SetMarkerColor(kMagenta);
0619 
0620   hhcorrbg_pp[i]->Fit("expo","rql","",statscale_lowlim,statscale_uplim);
0621   hhcorrbg_pp[i]->GetFunction("expo")->SetLineColor(kBlack);
0622   hhbottom_pp[i]->Fit("expo","rql","",statscale_lowlim,statscale_uplim);
0623   hhbottom_pp[i]->GetFunction("expo")->SetLineColor(kBlue);
0624   hhdy_pp[i]->Fit("expo","rql","",statscale_lowlim,statscale_uplim);
0625   hhdy_pp[i]->GetFunction("expo")->SetLineColor(kGreen+2);
0626 
0627 } // end loop over pT bins
0628 
0629 TCanvas* cbginpp = new TCanvas("cbginpp","corr bg in pp",10,10,700,700);
0630   hhall_pp[nbins]->SetAxisRange(7.,12.);
0631   hhcorrbg_pp[nbins]->Draw("pehist");
0632   hhbottom_pp[nbins]->Draw("histsame");
0633   hhdy_pp[nbins]->Draw("histsame");
0634 
0635 
0636 TCanvas* cpp = new TCanvas("cpp","corr bg + sig in pp",100,100,700,700);
0637   hhall_pp[nbins]->SetAxisRange(7.,12.);
0638   hhall_pp[nbins]->Draw("pehist");
0639   hhcorrbg_pp[nbins]->Draw("pesame");
0640   hhbottom_pp[nbins]->Draw("same");
0641   hhdy_pp[nbins]->Draw("same");
0642 
0643 TCanvas* cpp_vspt = new TCanvas("cpp_vspt","corr bg + sig vs pt in pp",50,50,1200,900);
0644 cpp_vspt->Divide(4,3);
0645 for(int i=0; i<nbins; i++) {
0646 cpp_vspt->cd(i+1);
0647   hhall_pp[i]->SetAxisRange(7.0,12.0); hhall_pp[i]->Draw("hist"); hhcorrbg_pp[i]->Draw("histsame");
0648   sprintf(tlchar,"%d-%d GeV",i,i+1);   tl[i] = new TLatex(8.0,hhcorrbg_pp[i]->GetMaximum()*0.95,tlchar); tl[i]->Draw();
0649 }
0650 
0651 //-----------------------------------------------------------------------------------------
0652 // Combinatorial BG calculated for 10B central AuAu events
0653 //-----------------------------------------------------------------------------------------
0654 TCanvas* cdummy = new TCanvas("cdummy","cdummy",0,0,500,500);
0655 
0656 f = new TFile("fakee_eideff09.root");
0657 for(int i=0; i<nbins+1; i++) {
0658   sprintf(tmpname,"hhfakefake_%d",i);
0659   hhfakefake[i] = (TH1D*)f->Get(tmpname);
0660   hhfakefake[i]->SetDirectory(gROOT);
0661 }
0662 f->Close();
0663 
0664 f = new TFile("crossterms_eideff09.root");
0665 for(int i=0; i<nbins+1; i++) {
0666   sprintf(tmpname,"hhfakehf_%d",i);
0667   hhfakehf[i] = (TH1D*)f->Get(tmpname);
0668   hhfakehf[i]->SetDirectory(gROOT);
0669 }
0670 f->Close();
0671 
0672 TF1* fbg = new TF1("fbg","exp([0]+[1]*x)+exp([2]+[3]*x)",8.,11.);
0673 fbg->SetParameters(10., -1.0, 4., -0.1);
0674 fbg->SetParLimits(1.,-999.,0.);
0675 fbg->SetParLimits(3.,-999.,0.);
0676 
0677 for(int i=0; i<nbins+1; i++) {
0678   sprintf(tmpname,"hhcombbg_%d",i);
0679   hhcombbg[i] = (TH1D*)hhfakefake[i]->Clone(tmpname);
0680   hhcombbg[i]->Add(hhfakehf[i]);
0681   sprintf(tmpname,"hhcombbg_scaled_%d",i);
0682   hhcombbg_scaled[i] = (TH1D*)hhcombbg[i]->Clone(tmpname);
0683   if(i==nbins) { fbg->SetParameters(10., -1.0, 4., -0.1); }
0684   hhcombbg[i]->Fit(fbg,"qrl","",statscale_lowlim,statscale_uplim);
0685 
0686     for(int k=1; k<=hhcombbg[i]->GetNbinsX(); k++) {
0687       if(hhcombbg[i]->GetBinLowEdge(k)<statscale_lowlim || (hhcombbg[i]->GetBinLowEdge(k)+hhcombbg[i]->GetBinWidth(k))>statscale_uplim) {
0688         hhcombbg_scaled[i]->SetBinContent(k,0.);
0689         hhcombbg_scaled[i]->SetBinError(k,0.);
0690       }
0691       else {
0692         double tmp = statscale * hhcombbg[i]->GetFunction("fbg")->Eval(hhcombbg[i]->GetBinCenter(k));
0693         double tmprnd = myrandom->Poisson(tmp);
0694         if(tmprnd<0.) { tmprnd=0.; }
0695         hhcombbg_scaled[i]->SetBinContent(k,tmprnd);
0696         hhcombbg_scaled[i]->SetBinError(k,sqrt(tmprnd));
0697       }
0698     }
0699   hhcombbg_scaled[i]->Fit(fbg,"qrl","",statscale_lowlim,statscale_uplim);
0700 }
0701 
0702 delete cdummy;
0703 
0704 TCanvas* C1 = new TCanvas("C1","Combinatorial BG (ALL p_{T})",100,100,600,600);
0705 C1->SetLogy();
0706   hhfakefake[nbins]->SetAxisRange(7.0,14.0);
0707   hhfakefake[nbins]->SetMinimum(0.1);
0708   hhfakefake[nbins]->SetMaximum(5000.);
0709   hhfakefake[nbins]->SetLineColor(kGreen+2);
0710   hhfakefake[nbins]->SetLineWidth(2);
0711   hhfakefake[nbins]->GetXaxis()->SetTitle("Transverse momentum [GeV/c]");
0712   hhfakefake[nbins]->GetXaxis()->SetTitleOffset(1.0);
0713   hhfakefake[nbins]->GetXaxis()->SetTitleColor(1);
0714   hhfakefake[nbins]->GetXaxis()->SetTitleSize(0.040);
0715   hhfakefake[nbins]->GetXaxis()->SetLabelSize(0.040);
0716   hhfakefake[nbins]->GetYaxis()->SetTitle("Combinatorial background");
0717   hhfakefake[nbins]->GetYaxis()->SetTitleOffset(1.3);
0718   hhfakefake[nbins]->GetYaxis()->SetTitleSize(0.040);
0719   hhfakefake[nbins]->GetYaxis()->SetLabelSize(0.040);
0720   hhfakefake[nbins]->Draw("e");
0721 
0722   hhfakehf[nbins]->SetLineColor(kOrange+4);
0723   hhfakehf[nbins]->SetLineWidth(2);
0724   hhfakehf[nbins]->Draw("esame");
0725 
0726   hhcombbg[nbins]->SetLineColor(kBlack);
0727   hhcombbg[nbins]->SetLineWidth(2);
0728   hhcombbg[nbins]->Draw("esame");
0729 
0730 TCanvas* C1sc = new TCanvas("C1sc","SCALED Combinatorial BG (ALL p_{T})",100,100,600,600);
0731 C1sc->SetLogy(); 
0732   hhcombbg_scaled[nbins]->SetAxisRange(7.,14.);
0733   hhcombbg_scaled[nbins]->Draw("esame");
0734 
0735 TCanvas* c00 = new TCanvas("c00","Combinatorial BG vs. p_{T}",150,150,1200,900);
0736 c00->Divide(4,3);
0737 
0738 for(int i=0; i<nbins; i++) {
0739   c00->cd(i+1);
0740   hhcombbg[i]->SetAxisRange(7.0,14.0); hhcombbg[i]->Draw();
0741   sprintf(tlchar,"%d-%d GeV",i,i+1);   tl[i] = new TLatex(8.0,hhcombbg[i]->GetMaximum()*0.95,tlchar); tl[i]->Draw();
0742 }
0743 
0744 TCanvas* c00scaled = new TCanvas("c00scaled","SCALED Combinatorial BG vs. p_{T}",150,150,1200,900);
0745 c00scaled->Divide(4,3);
0746 
0747 for(int i=0; i<nbins; i++) {
0748   c00scaled->cd(i+1);
0749   hhcombbg_scaled[i]->SetAxisRange(statscale_lowlim,statscale_uplim); hhcombbg_scaled[i]->Draw();
0750   sprintf(tlchar,"%d-%d GeV",i,i+1);   tl[i] = new TLatex(8.0,hhcombbg_scaled[i]->GetMaximum()*0.95,tlchar); tl[i]->Draw();
0751 }
0752 
0753 //---------------------------------------------
0754 //  Signal + BG
0755 //---------------------------------------------
0756 
0757 for(int i=0; i<nbins+1; i++) {
0758   sprintf(tmpname,"hhtotbg_scaled_%d",i);
0759   hhtotbg_scaled[i] = (TH1D*)hhcombbg_scaled[i]->Clone(tmpname);
0760   hhtotbg_scaled[i]->Add(hhcorrbg_scaled[i]);
0761 }
0762 
0763 for(int i=0; i<nbins+1; i++) {
0764 sprintf(tmpname,"hhall_scaled_%d",i);
0765 hhall_scaled[i] = (TH1D*)hhtotbg_scaled[i]->Clone(tmpname);
0766 hhall_scaled[i]->Add(hhups[i]);
0767 }
0768 
0769 TCanvas* c000 = new TCanvas("c000","Signal + Background vs. p_{T}",200,200,1200,900);
0770 c000->Divide(4,3);
0771 for(int i=0; i<nbins; i++) {
0772   c000->cd(i+1);
0773   hhall_scaled[i]->SetAxisRange(7.0,14.0); hhall_scaled[i]->SetMarkerStyle(1); hhall_scaled[i]->Draw("pehist");
0774   sprintf(tlchar,"%d-%d GeV",i,i+1);   tl[i] = new TLatex(8.0,hhall_scaled[i]->GetMaximum()*0.9,tlchar); tl[i]->Draw();
0775 }
0776 
0777 // Fit Signal + Correlated BG
0778 
0779 // hhfit[] is signal + correlated bg
0780 for(int i=0; i<nbins+1; i++) {
0781   sprintf(tmpname,"hhfit_%d",i);
0782   hhfit[i] = (TH1D*)hhall_scaled[i]->Clone(tmpname);
0783   for(int j=1; j<=hhall_scaled[i]->GetNbinsX(); j++) {
0784     hhfit[i]->SetBinContent(j,hhfit[i]->GetBinContent(j) - hhcombbg_scaled[i]->GetFunction("fbg")->Eval(hhcombbg_scaled[i]->GetBinCenter(j)));
0785     hhfit[i]->SetBinError(j,sqrt(pow(hhfit[i]->GetBinError(j),2)+hhcombbg_scaled[i]->GetFunction("fbg")->Eval(hhcombbg_scaled[i]->GetBinCenter(j))));
0786   }
0787 }
0788 
0789 //---------------------------------------------------------------------
0790 // plot signal + correlated bg for all pT
0791 //---------------------------------------------------------------------
0792 
0793   double tmppar0 = corrbgfitpar0+TMath::Log(statscale);
0794   double tmppar1 = corrbgfitpar1;
0795 //cout << "###### " << tmppar0 << " " << tmppar1 << endl;
0796 
0797 double ppauauscale = fTCBauau->GetParameter(0)/fTCBpp->GetParameter(0)*0.9783;
0798 fSandBpp->SetParameter(0,fTCBpp->GetParameter(0)*ppauauscale);
0799 fSandBpp->SetParameter(1,fTCBpp->GetParameter(1));
0800 fSandBpp->SetParameter(2,fTCBpp->GetParameter(2));
0801 fSandBpp->SetParameter(3,fTCBpp->GetParameter(3));
0802 fSandBpp->SetParameter(4,fTCBpp->GetParameter(4));
0803 fSandBpp->SetParameter(5,fTCBpp->GetParameter(5)*ppauauscale);
0804 fSandBpp->SetParameter(6,fTCBpp->GetParameter(6)*ppauauscale);
0805 fSandBpp->SetParameter(7,tmppar0);
0806 fSandBpp->SetParameter(8,tmppar1);
0807 fSandBpp->SetLineColor(kBlue);
0808 fSandBpp->SetLineStyle(2);
0809 
0810 //cout << "######### " << fTCBauau->GetParameter(0) << " " << fTCBauau->GetParameter(5) << " " << fTCBauau->GetParameter(6) << endl;
0811 fSandBauau->SetParameter(0,fTCBauau->GetParameter(0));
0812 fSandBauau->SetParameter(1,fTCBauau->GetParameter(1));
0813 fSandBauau->SetParameter(2,fTCBauau->GetParameter(2));
0814 fSandBauau->SetParameter(3,fTCBauau->GetParameter(3));
0815 fSandBauau->SetParameter(4,fTCBauau->GetParameter(4));
0816 fSandBauau->SetParameter(5,fTCBauau->GetParameter(5));
0817 fSandBauau->SetParameter(6,fTCBauau->GetParameter(6));
0818 fSandBauau->SetParameter(7,tmppar0);
0819 fSandBauau->SetParameter(8,tmppar1);
0820 fSandBauau->SetLineColor(kRed);
0821 
0822 
0823 TCanvas* cfitall = new TCanvas("cfitall","FIT all pT",270,270,600,600);
0824   hhfit[nbins]->SetAxisRange(7.0,14.);
0825   hhfit[nbins]->GetXaxis()->CenterTitle();
0826   hhfit[nbins]->GetXaxis()->SetTitle("Mass(e^{+}e^{-}) [GeV/c^2]");
0827   hhfit[nbins]->GetXaxis()->SetTitleOffset(1.1);
0828   hhfit[nbins]->GetXaxis()->SetLabelSize(0.045);
0829   hhfit[nbins]->GetYaxis()->CenterTitle();
0830   hhfit[nbins]->GetYaxis()->SetLabelSize(0.045);
0831   hhfit[nbins]->GetYaxis()->SetTitle("Events / (50 MeV/c^{2})");
0832   hhfit[nbins]->GetYaxis()->SetTitleOffset(1.5);
0833   hhfit[nbins]->Draw("pehist");
0834   fSandBpp->Draw("same");
0835   fSandBauau->Draw("same");
0836     //hhcorrbg_scaled[nbins]->SetLineColor(kRed);
0837     //hhcorrbg_scaled[nbins]->Scale(0.82);
0838     //hhcorrbg_scaled[nbins]->Scale(0.70);
0839     //hhcorrbg_scaled[nbins]->Draw("histsame");
0840 
0841   TF1* fmycorrbg = new TF1("fmycorrbg","exp([0]+[1]*x)",7.,14.);
0842   fmycorrbg->SetParameters(tmppar0,tmppar1);
0843   fmycorrbg->SetLineStyle(2);
0844   fmycorrbg->SetLineColor(kRed);
0845   fmycorrbg->Draw("same");
0846 
0847 
0848 double myheight = fTCBauau->GetParameter(0);
0849 TLatex* ld1 = new TLatex(10.1,myheight,"sPHENIX Simulation");
0850 ld1->SetTextSize(0.035);
0851 ld1->Draw();
0852 TLatex* ld2 = new TLatex(10.1,myheight-100.,"0-10% Au+Au #sqrt{s} = 200 GeV");
0853 ld2->SetTextSize(0.035);
0854 ld2->Draw();
0855 
0856 TCanvas* cfitall2 = new TCanvas("cfitall2","FIT all pT",270,270,600,600);
0857 TH1D* hhfit_tmp = (TH1D*)hhfit[nbins]->Clone("hhfit_tmp");
0858 hhfit_tmp->SetAxisRange(8.0,11.);
0859 hhfit_tmp->Draw("pehist");
0860 fSandBauau->Draw("same");
0861 fmycorrbg->Draw("same");
0862 
0863 
0864 //----------------------------------------------------------------------
0865 
0866 // plot signal + both backgrounds for all pT
0867 
0868 TCanvas* callpt = new TCanvas("callpt","Signal+BG (all p_{T})",300,300,600,600);
0869 
0870   hhall_scaled[nbins]->GetXaxis()->SetTitle("Invariant mass GeV/c");
0871   hhall_scaled[nbins]->SetLineColor(kBlack);
0872   hhall_scaled[nbins]->SetMarkerColor(kBlack);
0873   hhall_scaled[nbins]->SetMarkerStyle(20);
0874   hhall_scaled[nbins]->SetAxisRange(8.0,10.8);
0875   hhall_scaled[nbins]->Draw("pehist");
0876     hhcombbg_scaled[nbins]->SetLineColor(kBlue);
0877     hhcombbg_scaled[nbins]->Draw("histsame");
0878       hhcorrbg_scaled[nbins]->SetLineColor(kRed);
0879       hhcorrbg_scaled[nbins]->Draw("histsame");
0880 
0881 //----------------------------------------------------------------
0882 //  Calculate RAA
0883 //----------------------------------------------------------------
0884 
0885 double u1start = 9.25;
0886 double u1stop  = 9.65;
0887 double u2start = 9.80;
0888 double u2stop  = 10.20;
0889 double u3start = 10.20;
0890 double u3stop  = 10.55;
0891 
0892   double raa1[nbins+1],raa2[nbins+1],raa3[nbins+1],erraa1[nbins+1],erraa2[nbins+1],erraa3[nbins+1];
0893   for(int i=0; i<nbins; i++) { 
0894     //raa1[i] = supcor[0]; raa2[i] = supcor[1]; raa3[i] = supcor[2]; 
0895     raa1[i] = grRAA1S->Eval(0.5+i*1.0);
0896     raa2[i] = grRAA2S->Eval(0.5+i*1.0);
0897     raa3[i] = grRAA3S->Eval(0.5+i*1.0);
0898   }
0899   int fbin1 = hhall_scaled[nbins]->FindBin(u1start + 0.001);
0900   int lbin1 = hhall_scaled[nbins]->FindBin(u1stop  - 0.001);
0901   int fbin2 = hhall_scaled[nbins]->FindBin(u2start + 0.001);
0902   int lbin2 = hhall_scaled[nbins]->FindBin(u2stop  - 0.001);
0903   int fbin3 = hhall_scaled[nbins]->FindBin(u3start + 0.001);
0904   int lbin3 = hhall_scaled[nbins]->FindBin(u3stop  - 0.001);
0905   cout << "Y(1S) bin range: " << fbin1 << " - " << lbin1 << endl;
0906   cout << "Y(1S) inv. mass range: " << u1start << " - " << u1stop << endl;
0907   cout << "Y(2S) bin range: " << fbin2 << " - " << lbin2 << endl;
0908   cout << "Y(2S) inv. mass range: " << u2start << " - " << u2stop << endl;
0909   cout << "Y(3S) bin range: " << fbin3 << " - " << lbin3 << endl;
0910   cout << "Y(3S) inv. mass range: " << u3start << " - " << u3stop << endl;
0911 
0912   double sum1[99]   = {0.};
0913   double truesum1[99]   = {0.};
0914   double ersum1[99] = {0.};
0915   double sumpp1[99]   = {0.};
0916   double ersumpp1[99] = {0.};
0917   double sum2[99]   = {0.};
0918   double truesum2[99]   = {0.};
0919   double ersum2[99] = {0.};
0920   double sumpp2[99]   = {0.};
0921   double ersumpp2[99] = {0.};
0922   double sum3[99]   = {0.};
0923   double truesum3[99]   = {0.};
0924   double ersum3[99] = {0.};
0925   double sumpp3[99]   = {0.};
0926   double ersumpp3[99] = {0.};
0927 
0928   double sumsum1[99]   = {0.};
0929   double sumsum2[99]   = {0.};
0930   double sumsum3[99]   = {0.};
0931   double sumsum1pp[99]   = {0.};
0932   double sumsum2pp[99]   = {0.};
0933   double sumsum3pp[99]   = {0.};
0934 
0935   for(int i=0; i<nbins+1; i++) {
0936 
0937     double s1 = double(i); double s2 = double(i+1); if(i==nbins) {s1 = 0.;}
0938 
0939     for(int j=fbin1; j<=lbin1; j++) {
0940       sum1[i]   += (hhall_scaled[i]->GetBinContent(j) - hhcombbg_scaled[i]->GetFunction("fbg")->Eval(hhall_scaled[i]->GetBinCenter(j)) - hhcorrbg_scaled[i]->GetFunction("expo")->Eval(hhall_scaled[i]->GetBinCenter(j)));
0941       truesum1[i] += hhups1[i]->GetBinContent(j);
0942       ersum1[i] += hhall_scaled[i]->GetBinError(j)*hhall_scaled[i]->GetBinError(j);
0943       sumpp1[i]   += hhups1pp[i]->GetBinContent(j);
0944       ersumpp1[i] += hhupspp[i]->GetBinError(j)*hhupspp[i]->GetBinError(j);
0945     }
0946       sumsum1[i]     = truesum1[i];                   // direct count in mass range
0947       sumsum1pp[i]   = sumpp1[i];       
0948       //sumsum1[i]     = hhups1[i]->GetEntries();       // total number of upsilons in pT bin (rounded up)
0949       //sumsum1pp[i]   = hhups1pp[i]->GetEntries();
0950       //sumsum1[i]     = Nups1*fUpsilonPt->Integral(s1,s2)/upsnorm;  // total number of upsilons in pT bin
0951       //sumsum1pp[i]   = Nups1pp*fUpsilonPt->Integral(s1,s2)/upsnorm;
0952 
0953         if(sumsum1[i]>0. && sumsum1pp[i]>0.) {
0954           erraa1[i] = raa1[i]*sqrt(ersum1[i]/sumsum1[i]/sumsum1[i] + ersumpp1[i]/sumsum1pp[i]/sumsum1pp[i]);
0955         } else {raa1[i]=-1.0; erraa1[i] = 999.; }
0956 
0957     for(int j=fbin2; j<=lbin2; j++) {
0958       sum2[i]   += (hhall_scaled[i]->GetBinContent(j) - hhcombbg_scaled[i]->GetFunction("fbg")->Eval(hhall_scaled[i]->GetBinCenter(j)) - hhcorrbg_scaled[i]->GetFunction("expo")->Eval(hhall_scaled[i]->GetBinCenter(j)));
0959       truesum2[i] += hhups2[i]->GetBinContent(j);
0960       ersum2[i] += hhall_scaled[i]->GetBinError(j)*hhall_scaled[i]->GetBinError(j);
0961       sumpp2[i]   += hhups2pp[i]->GetBinContent(j);
0962       ersumpp2[i] += hhupspp[i]->GetBinError(j)*hhupspp[i]->GetBinError(j);
0963     }
0964       sumsum2[i]     = truesum2[i];                   // direct count in mass range
0965       sumsum2pp[i]   = sumpp2[i];       
0966       //sumsum2[i]     = hhups2[i]->GetEntries();       // total number of upsilons in pT bin (rounded up)
0967       //sumsum2pp[i]   = hhups2pp[i]->GetEntries();
0968       //sumsum2[i]     = Nups2*fUpsilonPt->Integral(s1,s2)/upsnorm;  // total number of upsilons in pT bin
0969       //sumsum2pp[i]   = Nups2pp*fUpsilonPt->Integral(s1,s2)/upsnorm;
0970       
0971         if(sumsum2[i]>0. && sumsum2pp[i]>0.) {
0972           erraa2[i] = raa2[i]*sqrt(ersum2[i]/sumsum2[i]/sumsum2[i] + ersumpp2[i]/sumsum2pp[i]/sumsum2pp[i]);
0973         } else {raa2[i]=-1.0; erraa2[i] = 999.; }
0974 
0975     for(int j=fbin3; j<=lbin3; j++) {
0976       sum3[i]   += (hhall_scaled[i]->GetBinContent(j) - hhcombbg_scaled[i]->GetFunction("fbg")->Eval(hhall_scaled[i]->GetBinCenter(j)) - hhcorrbg_scaled[i]->GetFunction("expo")->Eval(hhall_scaled[i]->GetBinCenter(j)));
0977       truesum3[i] += hhups3[i]->GetBinContent(j);
0978       ersum3[i] += hhall_scaled[i]->GetBinError(j)*hhall_scaled[i]->GetBinError(j);
0979       sumpp3[i]   += hhups3pp[i]->GetBinContent(j);
0980       ersumpp3[i] += hhupspp[i]->GetBinError(j)*hhupspp[i]->GetBinError(j);
0981     }
0982       sumsum3[i]     = truesum3[i];                   // direct count in mass range
0983       sumsum3pp[i]   = sumpp3[i];       
0984       //sumsum3[i]     = hhups3[i]->GetEntries();       // total number of upsilons in pT bin (rounded up)
0985       //sumsum3pp[i]   = hhups3pp[i]->GetEntries();
0986       //sumsum3[i]     = Nups3*fUpsilonPt->Integral(s1,s2)/upsnorm;   // total number of upsilons in pT bin
0987       //sumsum3pp[i]   = Nups3pp*fUpsilonPt->Integral(s1,s2)/upsnorm;
0988       
0989         if(truesum3[i]>0. && sumpp3[i]>0.) {
0990           erraa3[i] = raa3[i]*sqrt(ersum3[i]/sumsum3[i]/sumsum3[i] + ersumpp3[i]/sumsum3pp[i]/sumsum3pp[i]);
0991         } else {raa3[i]=-1.0; erraa3[i] = 999.; }
0992   
0993   }
0994 
0995 double raa2_rebin[9],raapt_rebin2[9],erraa2_rebin[9];
0996 double raa3_rebin[9],raapt_rebin3[9],erraa3_rebin[9];
0997 double sum2_rebin[9],ersum2_rebin[9],sum2pp_rebin[9],ersumpp2_rebin[9];
0998 double sum3_rebin[9],ersum3_rebin[9],sum3pp_rebin[9],ersumpp3_rebin[9];
0999 
1000 // Rebin Y(2S) by 2
1001 raapt_rebin2[0] = 1.;
1002 raapt_rebin2[1] = 3.;
1003 raapt_rebin2[2] = 5.;
1004 raapt_rebin2[3] = 7.;
1005 raa2_rebin[0] = grRAA2S->Eval(raapt_rebin2[0]);
1006 raa2_rebin[1] = grRAA2S->Eval(raapt_rebin2[1]);
1007 raa2_rebin[2] = grRAA2S->Eval(raapt_rebin2[2]);
1008 raa2_rebin[3] = grRAA2S->Eval(raapt_rebin2[3]);
1009 sum2_rebin[0] = truesum2[0]+truesum2[1];
1010 sum2_rebin[1] = truesum2[2]+truesum2[3];
1011 sum2_rebin[2] = truesum2[4]+truesum2[5];
1012 sum2_rebin[3] = truesum2[6]+truesum2[7]+truesum2[8]+truesum2[9];
1013 ersum2_rebin[0] = ersum2[0]+ersum2[1];
1014 ersum2_rebin[1] = ersum2[2]+ersum2[3];
1015 ersum2_rebin[2] = ersum2[4]+ersum2[5];
1016 ersum2_rebin[3] = ersum2[6]+ersum2[7]+ersum2[8]+ersum2[9];
1017 sum2pp_rebin[0] = sumpp2[0]+sumpp2[1];
1018 sum2pp_rebin[1] = sumpp2[2]+sumpp2[3];
1019 sum2pp_rebin[2] = sumpp2[4]+sumpp2[5];
1020 sum2pp_rebin[3] = sumpp2[6]+sumpp2[7]+sumpp2[8]+sumpp2[9];
1021 ersumpp2_rebin[0] = ersumpp2[0]+ersumpp2[1];
1022 ersumpp2_rebin[1] = ersumpp2[2]+ersumpp2[3];
1023 ersumpp2_rebin[2] = ersumpp2[4]+ersumpp2[5];
1024 ersumpp2_rebin[3] = ersumpp2[6]+ersumpp2[7]+ersumpp2[8]+ersumpp2[9];
1025   erraa2_rebin[0] = raa2[0]*sqrt(ersum2_rebin[0]/sum2_rebin[0]/sum2_rebin[0] + ersumpp2_rebin[0]/sum2pp_rebin[0]/sum2pp_rebin[0]);
1026   erraa2_rebin[1] = raa2[1]*sqrt(ersum2_rebin[1]/sum2_rebin[1]/sum2_rebin[1] + ersumpp2_rebin[1]/sum2pp_rebin[1]/sum2pp_rebin[1]);
1027   erraa2_rebin[2] = raa2[2]*sqrt(ersum2_rebin[2]/sum2_rebin[2]/sum2_rebin[2] + ersumpp2_rebin[2]/sum2pp_rebin[2]/sum2pp_rebin[2]);
1028   erraa2_rebin[3] = raa2[3]*sqrt(ersum2_rebin[3]/sum2_rebin[3]/sum2_rebin[3] + ersumpp2_rebin[3]/sum2pp_rebin[3]/sum2pp_rebin[3]);
1029 // Rebin Y(3S) by 4
1030 raapt_rebin3[0] = 2.;
1031 raapt_rebin3[1] = 6.;
1032 raa3_rebin[0] = grRAA3S->Eval(raapt_rebin3[0]);
1033 raa3_rebin[1] = grRAA3S->Eval(raapt_rebin3[1]);
1034 sum3_rebin[0] = truesum3[0]+truesum3[1]+truesum3[2]+truesum3[3];
1035 sum3_rebin[1] = truesum3[4]+truesum3[5]+truesum3[6]+truesum3[7]+truesum3[8]+truesum3[9];
1036 ersum3_rebin[0] = ersum3[0]+ersum3[1]+ersum3[2]+ersum3[3];
1037 ersum3_rebin[1] = ersum3[4]+ersum3[5]+ersum3[6]+ersum3[7]+ersum3[8]+ersum3[9];
1038 sum3pp_rebin[0] = sumpp3[0]+sumpp3[1]+sumpp3[3]+sumpp3[3];
1039 sum3pp_rebin[1] = sumpp3[4]+sumpp3[5]+sumpp3[6]+sumpp3[7]+sumpp3[8]+sumpp3[9];
1040 ersumpp3_rebin[0] = ersumpp3[0]+ersumpp3[1]+ersumpp3[2]+ersumpp3[3];
1041 ersumpp3_rebin[1] = ersumpp3[4]+ersumpp3[5]+ersumpp3[6]+ersumpp3[7]+ersumpp3[8]+ersumpp3[9];
1042   erraa3_rebin[0] = raa3[0]*sqrt(ersum3_rebin[0]/sum3_rebin[0]/sum3_rebin[0] + ersumpp3_rebin[0]/sum3pp_rebin[0]/sum3pp_rebin[0]);
1043   erraa3_rebin[1] = raa3[1]*sqrt(ersum3_rebin[1]/sum3_rebin[1]/sum3_rebin[1] + ersumpp3_rebin[1]/sum3pp_rebin[1]/sum3pp_rebin[1]);
1044 // Rebin Y(3S) by 8
1045 /*
1046 raapt_rebin3[0] = 5.;
1047 raa3_rebin[0] = grRAA3S->Eval(raapt_rebin3[0]);
1048 sum3_rebin[0] = truesum3[0]+truesum3[1]+truesum3[2]+truesum3[3]+
1049                 truesum3[4]+truesum3[5]+truesum3[6]+truesum3[7]+truesum3[8]+truesum3[9];
1050 ersum3_rebin[0] = ersum3[0]+ersum3[1]+ersum3[2]+ersum3[3]+
1051                   ersum3[4]+ersum3[5]+ersum3[6]+ersum3[7]+ersum3[8]+ersum3[9];
1052 sum3pp_rebin[0] = sumpp3[0]+sumpp3[1]+sumpp3[3]+sumpp3[3]+
1053                   sumpp3[4]+sumpp3[5]+sumpp3[6]+sumpp3[7]+sumpp3[8]+sumpp3[9];
1054 ersumpp3_rebin[0] = ersumpp3[0]+ersumpp3[1]+ersumpp3[2]+ersumpp3[3]+
1055                     ersumpp3[4]+ersumpp3[5]+ersumpp3[6]+ersumpp3[7]+ersumpp3[8]+ersumpp3[9];
1056 erraa3_rebin[0] = raa3[0]*sqrt(ersum3_rebin[0]/sum3_rebin[0]/sum3_rebin[0] + ersumpp3_rebin[0]/sum3pp_rebin[0]/sum3pp_rebin[0]);
1057 */
1058 
1059 cout << "====== Y(1S):" << endl;
1060   for(int i=0; i<nbins+1; i++) {
1061     double s1 = double(i); double s2 = double(i+1); if(i==nbins) {s1 = 0.;}
1062     cout << "   " << i << " " << truesum1[i] << "(" << Nups1*fUpsilonPt->Integral(s1,s2)/upsnorm << ")" << " +- " 
1063          << sqrt(ersum1[i]) << " \t\t pp: " << sumpp1[i] << " +- " << sqrt(ersumpp1[i]) << endl;
1064   }
1065 cout << "====== Y(2S):" << endl;
1066   for(int i=0; i<nbins+1; i++) {
1067     double s1 = double(i); double s2 = double(i+1); if(i==nbins) {s1 = 0.;}
1068     cout << "   " << i << " " << truesum2[i] << "(" << Nups2*fUpsilonPt->Integral(s1,s2)/upsnorm << ")" << " +- " 
1069          << sqrt(ersum2[i]) << " \t\t pp: " << sumpp2[i] << " +- " << sqrt(ersumpp2[i]) << endl;
1070   }
1071 cout << "====== Y(3S):" << endl;
1072   for(int i=0; i<nbins+1; i++) {
1073 
1074     double s1 = double(i); double s2 = double(i+1); if(i==nbins) {s1 = 0.;}
1075     cout << "   " << i << " " << truesum3[i] << "(" << Nups3*fUpsilonPt->Integral(s1,s2)/upsnorm << ")" << " +- " 
1076          << sqrt(ersum3[i]) << " \t\t pp: " << sumpp3[i] << " +- " << sqrt(ersumpp3[i]) << endl;
1077   }
1078 
1079 //-------------------------------------------------
1080 //  Plot RAA
1081 //-------------------------------------------------
1082 
1083 int npts1 = 10;
1084 int npts2 = 8;
1085 int npts3 = 7;
1086 int npts2_rebin = 4;
1087 int npts3_rebin = 2;
1088 
1089 TCanvas* craa = new TCanvas("craa","R_{AA}",120,120,800,600);
1090 TH2F* hh2 = new TH2F("hh2"," ",10,0.,float(npts1),10,0.,1.1);
1091 hh2->GetXaxis()->SetTitle("Transverse momentum [GeV/c]");
1092 hh2->GetXaxis()->SetTitleOffset(1.0);
1093 hh2->GetXaxis()->SetTitleColor(1);
1094 hh2->GetXaxis()->SetTitleSize(0.040);
1095 hh2->GetXaxis()->SetLabelSize(0.040);
1096 hh2->GetYaxis()->SetTitle("R_{AA}");
1097 hh2->GetYaxis()->SetTitleOffset(0.7);
1098 hh2->GetYaxis()->SetTitleSize(0.050);
1099 hh2->GetYaxis()->SetLabelSize(0.040);
1100 hh2->Draw();
1101 
1102 double xx1[nbins+1]; for(int i=0; i<nbins+1; i++) {xx1[i] = 0.5 + double(i);}
1103 double xx2[nbins+1]; for(int i=0; i<nbins+1; i++) {xx2[i] = 0.5 + double(i);}
1104 double xx2_rebin[nbins+1]; for(int i=0; i<npts2_rebin; i++) {xx2_rebin[i] = raapt_rebin2[i];}
1105 double xx3[nbins+1]; for(int i=0; i<nbins+1; i++) {xx3[i] = 0.5 + double(i);}
1106 double xx3_rebin[nbins+1]; for(int i=0; i<npts3_rebin; i++) {xx3_rebin[i] = raapt_rebin3[i];}
1107 
1108 xx3_rebin[0] = fUpsilonXPt->Integral(0.,4.)/fUpsilonPt->Integral(0.,4.);
1109 xx3_rebin[1] = fUpsilonXPt->Integral(4.,10.)/fUpsilonPt->Integral(4.,10.);
1110 
1111 TGraphErrors* gr1 = new TGraphErrors(npts1-1,xx1,raa1,0,erraa1);
1112 gr1->SetMarkerStyle(20);
1113 gr1->SetMarkerColor(kBlack);
1114 gr1->SetLineColor(kBlack);
1115 gr1->SetLineWidth(2);
1116 gr1->SetMarkerSize(1.5);
1117 gr1->SetName("gr1");
1118 gr1->Draw("p");
1119 
1120 //erraa2[7] = erraa2[7]*0.85;
1121 
1122 TGraphErrors* gr2 = new TGraphErrors(npts2-1,xx2,raa2,0,erraa2);
1123 gr2->SetMarkerStyle(20);
1124 gr2->SetMarkerColor(kRed);
1125 gr2->SetLineColor(kRed);
1126 gr2->SetLineWidth(2);
1127 gr2->SetMarkerSize(1.5);
1128 gr2->SetName("gr2");
1129 gr2->Draw("p");
1130 
1131 TGraphErrors* gr2_rebin = new TGraphErrors(npts2_rebin,xx2_rebin,raa2_rebin,0,erraa2_rebin);
1132 gr2_rebin->SetMarkerStyle(20);
1133 gr2_rebin->SetMarkerColor(kRed);
1134 gr2_rebin->SetLineColor(kRed);
1135 gr2_rebin->SetLineWidth(2);
1136 gr2_rebin->SetMarkerSize(1.5);
1137 gr2_rebin->SetName("gr2");
1138 //gr2_rebin->Draw("p");
1139 
1140 
1141 //--- 3S state -------------------
1142 // double dummy[9]; for(int i=0; i<9; i++) {dummy[i]=-99.;}
1143 // TGraphErrors* gr3 = new TGraphErrors(6,xx3,dummy,0,erraa3);
1144  TGraphErrors* gr3 = new TGraphErrors(6,xx3,raa3,0,erraa3);
1145  gr3->SetMarkerStyle(20);
1146  gr3->SetMarkerColor(kBlue);
1147  gr3->SetLineColor(kBlue);
1148  gr3->SetLineWidth(2);
1149  gr3->SetMarkerSize(1.5);
1150  gr3->SetName("gr3");
1151 // gr3->Draw("p");
1152  TGraphErrors* gr3_rebin = new TGraphErrors(npts3_rebin,xx3_rebin,raa3_rebin,0,erraa3_rebin);
1153  gr3_rebin->SetMarkerStyle(20);
1154  gr3_rebin->SetMarkerColor(kBlue);
1155  gr3_rebin->SetLineColor(kBlue);
1156  gr3_rebin->SetLineWidth(2);
1157  gr3_rebin->SetMarkerSize(1.5);
1158  gr3_rebin->SetName("gr3");
1159  gr3_rebin->Draw("p");
1160 /*
1161  TGraphErrors* gr3_rebin = new TGraphErrors(1,xx3_rebin,raa3_rebin,0,erraa3_rebin);
1162  gr3_rebin->SetMarkerStyle(20);
1163  gr3_rebin->SetMarkerColor(kBlue);
1164  gr3_rebin->SetLineColor(kBlue);
1165  gr3_rebin->SetLineWidth(2);
1166  gr3_rebin->SetMarkerSize(1.5);
1167  gr3_rebin->SetName("gr3");
1168  gr3_rebin->Draw("p");
1169 */
1170 //TArrow* aa[9];
1171 //TLine* ll[9];
1172 //erraa3[3] = erraa3[3]*1.5;
1173 //erraa3[4] = erraa3[4]*0.6;
1174 //erraa3[4] = erraa3[4]*1.5;
1175 //erraa3[6] = erraa3[6]*1.2;
1176 /*
1177 for(int i=0; i<npts3; i++) {
1178   aa[i] = new TArrow(xx3[i],1.64*erraa3[i],xx3[i],-0.02,0.02);
1179   aa[i]->SetLineColor(kBlue);
1180   aa[i]->SetLineWidth(2);
1181   aa[i]->Draw();
1182   ll[i] = new TLine(xx3[i]-0.15,1.64*erraa3[i],xx3[i]+0.15,1.64*erraa3[i]);
1183   ll[i]->SetLineColor(kBlue);
1184   ll[i]->SetLineWidth(2);
1185   ll[i]->Draw();
1186 }
1187 */
1188 //--- end 3S state --------------
1189 
1190 
1191 TLegend *leg = new TLegend(0.65,0.70,0.88,0.88);
1192 leg->SetBorderSize(1);
1193 leg->SetFillColor(10);
1194 leg->SetFillStyle(1001);
1195 TLegendEntry *entry1=leg->AddEntry("gr1","Y(1S)","p");
1196 TLegendEntry *entry2=leg->AddEntry("gr2","Y(2S)","p");
1197 //TLegendEntry *entry3=leg->AddEntry("gr3","Y(3S)","l");
1198 TLegendEntry *entry3=leg->AddEntry("gr3","Y(3S)","p");
1199 leg->Draw();
1200 
1201 TLatex* l1 = new TLatex(0.5,1.0,"#font[72]{sPHENIX} Projection"); l1->SetTextFont(42); l1->Draw();
1202 TLatex* l11 = new TLatex(0.5,0.90,"0-10% cent. Au+Au, Years 1-3"); l11->SetTextFont(42); l11->Draw();
1203 TLatex* l2 = new TLatex(0.5,0.80,"21 nb^{-1} rec. Au+Au"); l2->SetTextFont(42); l2->Draw();
1204 TLatex* l3 = new TLatex(0.5,0.70,"62 pb^{-1} samp. #it{p+p}"); l3->SetTextFont(42); l3->Draw();
1205 
1206 TLine* lll = new TLine(0.6,0.64,1.3,0.64);
1207 lll->SetLineColor(kBlue);
1208 lll->SetLineWidth(2);
1209 //lll->Draw();
1210 
1211 //==================================================================================
1212 
1213 /*
1214 fout = new TFile("test.root","recreate");
1215 for(int i=0; i<nbins+1; i++) {
1216   sprintf(tmpname,"hhfit_%d",i);
1217   hhfit[i]->Write();
1218   hhups1pp[i]->Write();
1219   hhups2pp[i]->Write();
1220   hhups3pp[i]->Write();
1221   hhupspp[i]->Write();
1222   hhups1[i]->Write();
1223   hhups2[i]->Write();
1224   hhups3[i]->Write();
1225   hhups[i]->Write();
1226   hhall_pp[i]->Write();
1227   hhcorrbg_pp[i]->Write();
1228   hhcorrbg_scaled[i]->Write();
1229 }
1230 fout->Write();
1231 fout->Close();
1232 */
1233 
1234 }
1235
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!