Back to home page

sPhenix code displayed by LXR
 
 

    

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

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_vscent_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 = 7;
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 double statscale_lowlim = 7.0;
0118 double statscale_uplim = 14.0;
0119 
0120   TF1* fCBpp = new TF1("fCBpp",CBFunction,5.,14.,5);
0121   TF1* fCBauau = new TF1("fCBauau",CBFunction,5.,14.,5);
0122   TF1* fCB1s = new TF1("fCB1s",CBFunction,5.,14.,5);
0123   TF1* fCB2s = new TF1("fCB2s",CBFunction,5.,14.,5);
0124   TF1* fCB3s = new TF1("fCB3s",CBFunction,5.,14.,5);
0125   TF1* fTCB = new TF1("fTCB",TripleCBFunction,5.,14.,7);
0126   TF1* fTCBpp = new TF1("fTCBpp",TripleCBFunction,5.,14.,7);
0127   TF1* fTCBauau = new TF1("fTCBauau",TripleCBFunction,5.,14.,7);
0128   TF1* fSandB = new TF1("fSandB",SandB_CBFunction,5.,14.,9);
0129   TF1* fSandBfordave = new TF1("fSandBfordave",SandB_CBFunction,5.,14.,9);
0130   TF1* fSandBpp = new TF1("fSandBpp",SandB_CBFunction,5.,14.,9);
0131   TF1* fSandBauau = new TF1("fSandBauau",SandB_CBFunction,5.,14.,9);
0132 
0133 //---------------------------------------------------------
0134 // Upsilons
0135 //---------------------------------------------------------
0136 
0137  // Input assumptions about level of suppression are taken from the theory paper: 
0138  // Michael Strickland, Dennis Bazow, Nucl.Phys. A879 (2012) 25-58
0139  //==============================================
0140  // Read in Strickland curves for different eta/s values
0141  // These are exclusive - i.e. they are the 
0142  // suppression of that state with no feed-down
0143  // effects included 
0144  //==============================================
0145  // The curves vs eta/s read in below were provided by Strickland privately, since
0146  // the paper contains only the eta/s = 1 cases for RHIC energy
0147  // Get the exclusive raa values for the Y1S, Y2S, Y3S, chib1, chib2
0148  // Checked that this is being read in correctly
0149                         
0150   double str_npart[101];
0151   double str_raa[5][3][101];
0152   for(int istate=0;istate<5;istate++)
0153     {
0154       for(int ietas=0;ietas<3;ietas++)
0155         {
0156           if(istate < 3)
0157             {
0158               char fname[500];
0159               sprintf(fname,"./strickland_calculations/Y%is-potb-eta%i-npart.dat",istate+1,ietas+1);
0160               ifstream fin(fname);
0161 
0162               for(int inpart=0;inpart<101;inpart++)
0163                 {
0164                   fin >> str_npart[inpart] >> str_raa[istate][ietas][inpart];
0165                 }
0166               fin.close();
0167             }
0168 
0169           if(istate > 2)
0170             {
0171               char fname[500];
0172               sprintf(fname,"./strickland_calculations/chib%i-potb-eta%i-npart.dat",istate-2,ietas+1);
0173               ifstream fin(fname);
0174 
0175               for(int inpart=0;inpart<101;inpart++)
0176                 {
0177                   fin >> str_npart[inpart] >> str_raa[istate][ietas][inpart];
0178                 }
0179               fin.close();
0180             }
0181         }
0182 
0183     }
0184 
0185   // Now construct the inclusive RAA values from Strickland's exclusive modifications
0186   // ff's from table II of the paper for the 1S state are as follows
0187   // They add up to 1.0
0188   double ff1S[5] = {0.51, 0.107, 0.008, 0.27, 0.105};  // Y1S, Y2S, Y3S, chib1, chib2 respectively
0189   // These are from arXiv:1210.7512
0190   double ff2S[2] = {0.5, 0.5};                         // Y2S and Y3S respectively
0191   double str_raa_inclusive[3][3][101];
0192 
0193   // checked the Y1S inclusive result against figure 10 of arXiv:1112.2761
0194   // There is no plot available for the 2S and 3S 
0195   for(int ietas=0;ietas<3;ietas++)
0196     for(int inpart=0;inpart<101;inpart++)
0197       {
0198         str_raa_inclusive[0][ietas][inpart] =
0199           str_raa[0][ietas][inpart] * ff1S[0]
0200           + str_raa[1][ietas][inpart] * ff1S[1]
0201           + str_raa[2][ietas][inpart] * ff1S[2]
0202           + str_raa[3][ietas][inpart] * ff1S[3]
0203           + str_raa[4][ietas][inpart] * ff1S[4];
0204 
0205         str_raa_inclusive[1][ietas][inpart] =
0206           str_raa[1][ietas][inpart] * ff2S[0]
0207           + str_raa[2][ietas][inpart] * ff2S[1];
0208 
0209         str_raa_inclusive[2][ietas][inpart] = str_raa[2][ietas][inpart];
0210       }
0211 
0212   double strick_raa1_eta1[101],strick_raa1_eta2[101],strick_raa1_eta3[101];
0213   double strick_raa2_eta1[101],strick_raa2_eta2[101],strick_raa2_eta3[101];
0214   double strick_raa3_eta1[101],strick_raa3_eta2[101],strick_raa3_eta3[101];
0215     for(int ipart=0;ipart<101;ipart++) {
0216       strick_raa1_eta1[ipart] = str_raa_inclusive[0][0][ipart];
0217       strick_raa1_eta2[ipart] = str_raa_inclusive[0][1][ipart];
0218       strick_raa1_eta3[ipart] = str_raa_inclusive[0][2][ipart];
0219       strick_raa2_eta1[ipart] = str_raa_inclusive[1][0][ipart];
0220       strick_raa2_eta2[ipart] = str_raa_inclusive[1][1][ipart];
0221       strick_raa2_eta3[ipart] = str_raa_inclusive[1][2][ipart];
0222       strick_raa3_eta1[ipart] = str_raa_inclusive[2][0][ipart];
0223       strick_raa3_eta2[ipart] = str_raa_inclusive[2][1][ipart];
0224       strick_raa3_eta3[ipart] = str_raa_inclusive[2][2][ipart];
0225     }
0226 
0227 //======================================================================================
0228 //======================================================================================
0229 //======================================================================================
0230 
0231 TGraph* grRAA1S = new TGraph(101,str_npart,strick_raa1_eta2);
0232 TGraph* grRAA2S = new TGraph(101,str_npart,strick_raa2_eta2);
0233 TGraph* grRAA3S = new TGraph(101,str_npart,strick_raa3_eta2);
0234 TGraph* grRAA1S_eta1 = new TGraph(101,str_npart,strick_raa1_eta1);
0235 TGraph* grRAA2S_eta1 = new TGraph(101,str_npart,strick_raa2_eta1);
0236 TGraph* grRAA3S_eta1 = new TGraph(101,str_npart,strick_raa3_eta1);
0237 TGraph* grRAA1S_eta3 = new TGraph(101,str_npart,strick_raa1_eta3);
0238 TGraph* grRAA2S_eta3 = new TGraph(101,str_npart,strick_raa2_eta3);
0239 TGraph* grRAA3S_eta3 = new TGraph(101,str_npart,strick_raa3_eta3);
0240 grRAA1S->SetLineColor(kBlack);
0241 grRAA1S->SetLineStyle(7);
0242 grRAA2S->SetLineColor(kRed);
0243 grRAA2S->SetLineStyle(7);
0244 grRAA3S->SetLineColor(kBlue);
0245 grRAA3S->SetLineStyle(7);
0246 grRAA1S_eta1->SetLineColor(kBlack);
0247 grRAA1S_eta1->SetLineStyle(1);
0248 grRAA2S_eta1->SetLineColor(kRed);
0249 grRAA2S_eta1->SetLineStyle(1);
0250 grRAA3S_eta1->SetLineColor(kBlue);
0251 grRAA3S_eta1->SetLineStyle(1);
0252 grRAA1S_eta3->SetLineColor(kBlack);
0253 grRAA1S_eta3->SetLineStyle(8);
0254 grRAA2S_eta3->SetLineColor(kRed);
0255 grRAA2S_eta3->SetLineStyle(8);
0256 grRAA3S_eta3->SetLineColor(kBlue);
0257 grRAA3S_eta3->SetLineStyle(8);
0258 
0259 int nchan=400;
0260 double start=0.0;
0261 double stop=20.0;
0262 
0263 string str_UpsilonPt  = "(2.0*3.14159*x*[0]*pow((1 + x*x/(4*[1]) ),-[2]))";   // dN/dpT
0264 string str_UpsilonXPt = "(2.0*3.14159*x*x*[0]*pow((1 + x*x/(4*[1]) ),-[2]))"; // need this for mean pT calculation
0265 TF1* fUpsilonPt = new TF1("fUpsilonPt",str_UpsilonPt.c_str(),0.,20.);
0266 TF1* fUpsilonXPt = new TF1("fUpsilonXPt",str_UpsilonXPt.c_str(),0.,20.);
0267 fUpsilonPt->SetParameters(72.1, 26.516, 10.6834);
0268 fUpsilonXPt->SetParameters(72.1, 26.516, 10.6834);
0269 double upsnorm = fUpsilonPt->Integral(0.,20.);
0270 
0271 double frac[3];  // upsilons fractions
0272   frac[0] = 0.7117;
0273   frac[1] = 0.1851;
0274   frac[2] = 0.1032;
0275 double scale[3]; // mass scaling
0276   scale[0] = 1.0;
0277   scale[1] = 1.0595;
0278   scale[2] = 1.0946;
0279 
0280 
0281 double Npart[nbins+1], NpartAvg=0.;
0282 Npart[0] = 325.;
0283 Npart[1] = 235.;
0284 Npart[2] = 167.;
0285 Npart[3] = 114.;
0286 Npart[4] = 74.;
0287 Npart[5] = 46.;
0288 Npart[6] = 14.5;
0289 for(int i=0; i<6; i++) {NpartAvg += Npart[i];} NpartAvg = NpartAvg/6.;
0290 cout << "Npart for 0-60% centrality = " << NpartAvg << endl;
0291 cout << "Raa for 0-60% centrality = " << grRAA1S->Eval(NpartAvg) << " " << grRAA2S->Eval(NpartAvg) << " " << grRAA3S->Eval(NpartAvg) << endl;
0292 
0293 double Ncoll[nbins+1];
0294 Ncoll[0] = 955.;
0295 Ncoll[1] = 603.;
0296 Ncoll[2] = 374.;
0297 Ncoll[3] = 220.;
0298 Ncoll[4] = 120.;
0299 Ncoll[5] = 61.;
0300 Ncoll[6] = 14.5;
0301 double NcollAvg=0.; for(int i=0; i<6; i++) {NcollAvg += Ncoll[i];} NcollAvg = NcollAvg/6.;
0302 cout << "Ncoll for 0-60% centrality = " << NcollAvg << endl;
0303 
0304 double Npionpairs[nbins+1];
0305 Npionpairs[0] = 1.000;
0306 Npionpairs[1] = 0.674;
0307 Npionpairs[2] = 0.418;
0308 Npionpairs[3] = 0.212;
0309 Npionpairs[4] = 0.099;
0310 Npionpairs[5] = 0.037;
0311 Npionpairs[6] = 0.033;  // This includes the fact that this bin is 3.2 times wider
0312 double NpionpairsAvg=0.; for(int i=0; i<6; i++) {NpionpairsAvg += Npionpairs[i];} NpionpairsAvg = NpionpairsAvg/6.;
0313 cout << "Npionpairs for 0-60% centrality = " << NpionpairsAvg << endl;
0314 
0315 //---------- Use fixed numbers for 2022 estimate --------------------------------------
0316 // Numbers from Marzia for 140 B MB Au+Au and 2400 B p+p events
0317 // RAA for Y(3S) is 8.5%
0318 
0319   int Nups1[nbins],Nups2[nbins],Nups3[nbins];
0320   int Nups1tot=0, Nups2tot=0, Nups3tot=0;
0321 
0322   Nups1[0]   = 7011;
0323   Nups1[1]   = 4807;
0324   Nups1[2]   = 3304;
0325   Nups1[3]   = 2185;
0326   Nups1[4]   = 1344;
0327   Nups1[5]   = 754;
0328   Nups1[6]   = 634;
0329   for(int i=0; i<6; i++) {Nups1tot += Nups1[i];}
0330 
0331   Nups2[0]   = 562;
0332   Nups2[1]   = 436;
0333   Nups2[2]   = 344;
0334   Nups2[3]   = 264;
0335   Nups2[4]   = 191;
0336   Nups2[5]   = 123;
0337   Nups2[6]   = 156;
0338   for(int i=0; i<6; i++) {Nups2tot += Nups2[i];}
0339 
0340   Nups3[0]   = 156;
0341   Nups3[1]   = 120;
0342   Nups3[2]   = 95;
0343   Nups3[3]   = 73;
0344   Nups3[4]   = 53;
0345   Nups3[5]   = 52;
0346   Nups3[6]   = 86;
0347   for(int i=0; i<6; i++) {Nups3tot += Nups3[i];}
0348   cout << "Number of Upsilons in 0-60% centrality = " << Nups1tot << " " << Nups2tot << " " << Nups3tot << endl;
0349 
0350 
0351   int Nups1pp = 2.86e+03;
0352   int Nups2pp = 7.16e+02;
0353   int Nups3pp = 3.98e+02;
0354 
0355 //====================================================================================
0356 //====================================================================================
0357 //====================================================================================
0358 
0359 //
0360 // these histograms (hhups*) are randomly generated using the fit above to single 
0361 // peak (fCB function) and used for the RAA uncertainty calculation
0362 //
0363 // FORCE specific RESOLUTION and tail shape
0364   double tonypar1 = 0.98;  // Tony's 100 pion simulation April 2019
0365   double tonypar2 = 0.93; // Tony's 100 pion simulation April 2019
0366   //double tonypar3 = 9.437;  // Tony's 100 pion simulation April 2019
0367   double tonypar3 = 9.448;  // benchmark
0368   double tonypar4 = 0.100;  // benchmark
0369   double tonypar4pp = 0.089; // Tony's 100 pion simulation April 2019
0370   fCBpp->SetParameter(0,1000.); 
0371   fCBpp->SetParameter(1,tonypar1); 
0372   fCBpp->SetParameter(2,tonypar2); 
0373   fCBpp->SetParameter(3,tonypar3); 
0374   fCBpp->SetParameter(4,tonypar4pp);  
0375   fCBauau->SetParameter(0,1000.); 
0376   fCBauau->SetParameter(1,tonypar1); 
0377   fCBauau->SetParameter(2,tonypar2); 
0378   fCBauau->SetParameter(3,tonypar3); 
0379   fCBauau->SetParameter(4,tonypar4);  
0380 
0381 char hhname[999];
0382 TH1D* hhups[nbins+1];
0383 TH1D* hhups1[nbins+1];
0384 TH1D* hhups2[nbins+1];
0385 TH1D* hhups3[nbins+1];
0386 TH1D* hhupspp;
0387 TH1D* hhups1pp;
0388 TH1D* hhups2pp;
0389 TH1D* hhups3pp;
0390 for(int i=0; i<nbins+1; i++) {
0391   sprintf(hhname,"hhups_%d",i);
0392   hhups[i] = new TH1D(hhname,"",nchan,start,stop);
0393   hhups[i]->Sumw2();
0394   sprintf(hhname,"hhups1_%d",i);
0395   hhups1[i] = new TH1D(hhname,"",nchan,start,stop);
0396   hhups1[i]->Sumw2();
0397   sprintf(hhname,"hhups2_%d",i);
0398   hhups2[i] = new TH1D(hhname,"",nchan,start,stop);
0399   hhups2[i]->Sumw2();
0400   sprintf(hhname,"hhups3_%d",i);
0401   hhups3[i] = new TH1D(hhname,"",nchan,start,stop);
0402   hhups3[i]->Sumw2();
0403   hhups[i]->SetLineWidth(2);
0404   hhups1[i]->SetLineWidth(2);
0405   hhups2[i]->SetLineWidth(2);
0406   hhups3[i]->SetLineWidth(2);
0407 }
0408   sprintf(hhname,"hhupspp");
0409   hhupspp= new TH1D(hhname,"",nchan,start,stop);
0410   hhupspp->Sumw2();
0411   sprintf(hhname,"hhups1pp");
0412   hhups1pp = new TH1D(hhname,"",nchan,start,stop);
0413   hhups1pp->Sumw2();
0414   sprintf(hhname,"hhups2pp");
0415   hhups2pp = new TH1D(hhname,"",nchan,start,stop);
0416   hhups2pp->Sumw2();
0417   sprintf(hhname,"hhups3pp");
0418   hhups3pp = new TH1D(hhname,"",nchan,start,stop);
0419   hhups3pp->Sumw2();
0420   hhupspp->SetLineWidth(2);
0421   hhups1pp->SetLineWidth(2);
0422   hhups2pp->SetLineWidth(2);
0423   hhups3pp->SetLineWidth(2);
0424 
0425 for(int j=0; j<nbins; j++) {
0426   double s1 = j*1.0;
0427   double s2 = s1 + 1.0;
0428   fCBauau->SetParameter(3,tonypar3);
0429     for(int i=0; i<int(Nups1[j]+0.5); i++) { double myrnd = fCBauau->GetRandom(); hhups1[j]->Fill(myrnd); hhups[j]->Fill(myrnd); }
0430   fCBauau->SetParameter(3,tonypar3*scale[1]);
0431     for(int i=0; i<int(Nups2[j]+0.5); i++) { double myrnd = fCBauau->GetRandom(); hhups2[j]->Fill(myrnd); hhups[j]->Fill(myrnd); }
0432   fCBauau->SetParameter(3,tonypar3*scale[2]);
0433     for(int i=0; i<int(Nups3[j]+0.5); i++) { double myrnd = fCBauau->GetRandom(); hhups3[j]->Fill(myrnd); hhups[j]->Fill(myrnd); }
0434 } // end loop over centrality bins
0435 
0436   fCBpp->SetParameter(3,tonypar3);
0437     for(int i=0; i<int(Nups1pp+0.5); i++) { double myrnd = fCBpp->GetRandom(); hhups1pp->Fill(myrnd); hhupspp->Fill(myrnd); }
0438   fCBpp->SetParameter(3,tonypar3*scale[1]);
0439     for(int i=0; i<int(Nups2pp+0.5); i++) { double myrnd = fCBpp->GetRandom(); hhups2pp->Fill(myrnd); hhupspp->Fill(myrnd); }
0440   fCBpp->SetParameter(3,tonypar3*scale[2]);
0441     for(int i=0; i<int(Nups3pp+0.5); i++) { double myrnd = fCBpp->GetRandom(); hhups3pp->Fill(myrnd); hhupspp->Fill(myrnd); }
0442 
0443 
0444 // fit and draw pp no bg
0445 
0446 TCanvas* cupspp = new TCanvas("cupspp","Upsilons in p+p",100,100,600,600);
0447   fTCBpp->SetParameter(0,2000.);
0448   fTCBpp->FixParameter(1,tonypar1);
0449   fTCBpp->FixParameter(2,tonypar2);
0450   fTCBpp->SetParameter(3,tonypar3);
0451   fTCBpp->FixParameter(4,tonypar4); // fix width from the single peak fit
0452   fTCBpp->SetParameter(5,500.);
0453   fTCBpp->SetParameter(6,100.);
0454   hhupspp->Fit(fTCBpp,"rl","",7.,11.);
0455   hhupspp->SetAxisRange(7.,11.);
0456   hhupspp->SetMarkerSize(1.0);
0457   hhupspp->GetXaxis()->SetTitle("Invariant mass [GeV/c^{2}]");
0458   hhupspp->GetXaxis()->SetTitleOffset(1.0);
0459   double tmpamp1 = hhupspp->GetFunction("fTCBpp")->GetParameter(0);
0460   double tmpamp5 = tmpamp1*frac[1]/frac[0]; // correct fit for correct upsilon states ratio
0461   double tmpamp6 = tmpamp1*frac[2]/frac[0];
0462   hhupspp->Draw();
0463 
0464   fCB1s->SetLineColor(kBlue);
0465   fCB1s->SetLineWidth(1);
0466     fCB1s->SetParameter(0,fTCBpp->GetParameter(0));
0467     fCB1s->SetParameter(1,fTCBpp->GetParameter(1));
0468     fCB1s->SetParameter(2,fTCBpp->GetParameter(2));
0469     fCB1s->SetParameter(3,fTCBpp->GetParameter(3)*scale[0]);
0470     fCB1s->SetParameter(4,fTCBpp->GetParameter(4));
0471   fCB2s->SetLineColor(kRed);
0472   fCB2s->SetLineWidth(1);
0473     fCB2s->SetParameter(0,tmpamp5);
0474     fCB2s->SetParameter(1,fTCBpp->GetParameter(1));
0475     fCB2s->SetParameter(2,fTCBpp->GetParameter(2));
0476     fCB2s->SetParameter(3,fTCBpp->GetParameter(3)*scale[1]);
0477     fCB2s->SetParameter(4,fTCBpp->GetParameter(4));
0478   fCB3s->SetLineColor(kGreen+2);
0479   fCB3s->SetLineWidth(1);
0480     fCB3s->SetParameter(0,tmpamp6);
0481     fCB3s->SetParameter(1,fTCBpp->GetParameter(1));
0482     fCB3s->SetParameter(2,fTCBpp->GetParameter(2));
0483     fCB3s->SetParameter(3,fTCBpp->GetParameter(3)*scale[2]);
0484     fCB3s->SetParameter(4,fTCBpp->GetParameter(4));
0485   fCB1s->Draw("same");
0486   fCB2s->Draw("same");
0487   fCB3s->Draw("same");
0488 
0489 // fit and draw AuAu no bg
0490 
0491 TCanvas* cupsauau = new TCanvas("cupsauau","Upsilons in Central Au+Au",100,100,600,600);
0492   fTCBauau->SetParameter(0,2000.);
0493   fTCBauau->FixParameter(1,tonypar1);
0494   fTCBauau->FixParameter(2,tonypar2);
0495   fTCBauau->SetParameter(3,tonypar3);
0496   fTCBauau->FixParameter(4,tonypar4); // fix width from the single peak fit
0497   fTCBauau->SetParameter(5,500.);
0498   fTCBauau->SetParameter(6,100.);
0499   hhups[0]->Fit(fTCBauau,"rl","",7.,11.);
0500   hhups[0]->SetAxisRange(8.5,11.);
0501   hhups[0]->SetMarkerSize(1.0);
0502   hhups[0]->GetXaxis()->SetTitle("Invariant mass [GeV/c^{2}]");
0503   hhups[0]->GetXaxis()->SetTitleOffset(1.0);
0504 //  tmpamp1 = hhups[nbins]->GetFunction("fTCBauau")->GetParameter(0);
0505 //  tmpamp5 = tmpamp1*frac[1]/frac[0]; // correct fit for correct upsilon states ratio
0506 //  tmpamp6 = tmpamp1*frac[2]/frac[0];
0507   hhups[0]->Draw();  // 0-10% central
0508 
0509 TCanvas* cdummy1 = new TCanvas("cdummy1","cdummy1",0,0,500,500);
0510 //----------------------------------------------------
0511 //  Backgrounds
0512 //----------------------------------------------------
0513 
0514 TH1D* hhall[nbins+1];
0515 TH1D* hhall_scaled[nbins+1];
0516 
0517 TH1D* hhtotbg[nbins+1]; 
0518 TH1D* hhtotbg_scaled[nbins+1]; 
0519 TH1D* hhcombbg[nbins+1];
0520 TH1D* hhcombbg_scaled[nbins+1];
0521 TH1D* hhfakefake[nbins+1];
0522 TH1D* hhfakehf[nbins+1];
0523 TH1D* hhbottom[nbins+1];
0524 TH1D* hhcharm[nbins+1];
0525 TH1D* hhdy[nbins+1];
0526 TH1D* hhcorrbg[nbins+1];
0527 TH1D* hhcorrbg_scaled[nbins+1];
0528 TH1D* hhfit[nbins+1];
0529 char tmpname[999];
0530 
0531 //----------------------------------------------------------------------------------------
0532 // Correlated BG calculated for 10B central AuAu events
0533 //----------------------------------------------------------------------------------------
0534 
0535 double corrbgfitpar0;
0536 double corrbgfitpar1;
0537 
0538 TFile* f=new TFile("ccbb_eideff09.root");
0539 
0540   sprintf(tmpname,"hhbottom_15");  // 15 is integrated over pT
0541   hhbottom[nbins] = (TH1D*)f->Get(tmpname);
0542   hhbottom[nbins]->SetDirectory(gROOT);
0543   sprintf(tmpname,"hhcharm_15");
0544   hhcharm[nbins] = (TH1D*)f->Get(tmpname);
0545   hhcharm[nbins]->SetDirectory(gROOT);
0546   sprintf(tmpname,"hhdy_15");
0547   hhdy[nbins] = (TH1D*)f->Get(tmpname);
0548   hhdy[nbins]->SetDirectory(gROOT);
0549   sprintf(tmpname,"hhcorrbg_15");
0550   hhcorrbg[nbins] = (TH1D*)hhbottom[nbins]->Clone(tmpname);
0551   hhcorrbg[nbins]->Add(hhcharm[nbins]);
0552   hhcorrbg[nbins]->Add(hhdy[nbins]);
0553   sprintf(tmpname,"hhcorrbg_scaled_15");
0554   hhcorrbg_scaled[nbins] = (TH1D*)hhcorrbg[nbins]->Clone(tmpname);
0555   hhcorrbg[nbins]->Fit("expo","rql","",statscale_lowlim,statscale_uplim);
0556   hhbottom[nbins]->Fit("expo","rql","",statscale_lowlim,statscale_uplim);
0557   hhdy[nbins]->Fit("expo","rql","",statscale_lowlim,statscale_uplim);
0558     corrbgfitpar0 = hhcorrbg[nbins]->GetFunction("expo")->GetParameter(0);
0559     corrbgfitpar1 = hhcorrbg[nbins]->GetFunction("expo")->GetParameter(1);
0560     cout << "bgpar0["<< nbins <<"]="<<hhcorrbg[nbins]->GetFunction("expo")->GetParameter(0)+TMath::Log(statscale)<<";"<< endl; 
0561     cout << "bgpar1["<< nbins <<"]="<<hhcorrbg[nbins]->GetFunction("expo")->GetParameter(1)<<";"<< endl; 
0562     for(int k=1; k<=hhcorrbg[nbins]->GetNbinsX(); k++) {
0563       if(hhcorrbg[nbins]->GetBinLowEdge(k)<statscale_lowlim || (hhcorrbg[nbins]->GetBinLowEdge(k)+hhcorrbg[nbins]->GetBinWidth(k))>statscale_uplim) {
0564         hhcorrbg_scaled[nbins]->SetBinContent(k,0.); 
0565         hhcorrbg_scaled[nbins]->SetBinError(k,0.);
0566       }
0567       else {
0568         double tmp = statscale * hhcorrbg[nbins]->GetFunction("expo")->Eval(hhcorrbg[nbins]->GetBinCenter(k));
0569         double tmprnd = myrandom->Poisson(tmp); 
0570         if(tmprnd<0.) { tmprnd=0.; }
0571         hhcorrbg_scaled[nbins]->SetBinContent(k,tmprnd); 
0572         hhcorrbg_scaled[nbins]->SetBinError(k,sqrt(tmprnd));
0573       }
0574     }
0575   hhcorrbg_scaled[nbins]->Fit("expo","rql","",statscale_lowlim,statscale_uplim);
0576   hhcorrbg[nbins]->SetDirectory(gROOT);
0577   hhcorrbg_scaled[nbins]->SetDirectory(gROOT);
0578 
0579 f->Close();
0580 
0581 // Correlated background in AuAu vs centrality
0582 cout << "kuku1" << endl;
0583 for(int i=0; i<nbins; i++) {
0584 
0585   sprintf(tmpname,"hhcorrbg_%d",i);
0586   hhcorrbg_scaled[i] = (TH1D*)hhcorrbg_scaled[nbins]->Clone(tmpname);
0587 
0588     for(int k=1; k<=hhcorrbg_scaled[nbins]->GetNbinsX(); k++) {
0589       if(hhcorrbg_scaled[nbins]->GetBinLowEdge(k)<statscale_lowlim || (hhcorrbg_scaled[nbins]->GetBinLowEdge(k)+hhcorrbg_scaled[nbins]->GetBinWidth(k))>statscale_uplim) {
0590         hhcorrbg_scaled[i]->SetBinContent(k,0.);
0591         hhcorrbg_scaled[i]->SetBinError(k,0.);
0592       }
0593       else {
0594         double tmp = (Ncoll[i]/Ncoll[0]) * hhcorrbg_scaled[nbins]->GetFunction("expo")->Eval(hhcorrbg_scaled[nbins]->GetBinCenter(k));
0595         double tmprnd = myrandom->Poisson(tmp);
0596         if(tmprnd<0.) { tmprnd=0.; }
0597         hhcorrbg_scaled[i]->SetBinContent(k,tmprnd);
0598         hhcorrbg_scaled[i]->SetBinError(k,sqrt(tmprnd));
0599       }
0600     }
0601   hhcorrbg_scaled[i]->Fit("expo","rql","",statscale_lowlim,statscale_uplim);
0602 
0603 } // end loop over centrality bins
0604 
0605 
0606 delete cdummy1;
0607 
0608 TCanvas* c111 = new TCanvas("c111","Au+Au Correlated Background vs. Centrality",200,200,1200,600);
0609 c111->Divide(4,2);
0610 for(int i=0; i<nbins; i++) {
0611   c111->cd(i+1);
0612   hhcorrbg_scaled[i]->SetAxisRange(8.5,11.0); hhcorrbg_scaled[i]->SetMarkerStyle(1); hhcorrbg_scaled[i]->Draw("pe");
0613   sprintf(tlchar,"%d-%d",10*i,10*(i+1));   tl[i] = new TLatex(9.0,hhcorrbg_scaled[i]->GetMaximum()*0.9,tlchar); tl[i]->Draw();
0614 }
0615 
0616 
0617 //-----------------------------------------------------------------------------------------
0618 //  Correlated bg in p+p
0619 //-----------------------------------------------------------------------------------------
0620 
0621 TH1D* hhbottom_pp;
0622 TH1D* hhdy_pp;
0623 TH1D* hhcorrbg_pp;
0624 TH1D* hhall_pp;
0625 
0626 double ppcorr = (2400./14.)/962.; // 2400B pp and 140B MB AuAu
0627 TF1* fbottom_nosup_corr = new TF1("fbottom_nosup_corr","[0]+[1]*x",5.,14.);
0628 fbottom_nosup_corr->SetParameters(-2.13861, 0.683323);
0629 
0630   sprintf(tmpname,"hhbottom_pp");
0631   hhbottom_pp = (TH1D*)hhbottom[nbins]->Clone(tmpname);
0632   for(int k=1; k<=hhbottom_pp->GetNbinsX(); k++) {
0633     if(hhbottom_pp->GetBinLowEdge(k)<statscale_lowlim || (hhbottom_pp->GetBinLowEdge(k)+hhbottom_pp->GetBinWidth(k))>statscale_uplim) {
0634       hhbottom_pp->SetBinContent(k,0.);
0635       hhbottom_pp->SetBinError(k,0.);
0636     }
0637     else {
0638       double tmp = ppcorr * fbottom_nosup_corr->Eval(hhbottom[nbins]->GetBinCenter(k)) * hhbottom[nbins]->GetFunction("expo")->Eval(hhbottom[nbins]->GetBinCenter(k));
0639       double tmprnd = myrandom->Poisson(tmp);
0640       if(tmprnd<0.) { tmprnd=0.; }
0641       hhbottom_pp->SetBinContent(k,tmprnd);
0642       hhbottom_pp->SetBinError(k,sqrt(tmprnd));
0643     }
0644   }
0645 
0646   sprintf(tmpname,"hhdy_pp");
0647   hhdy_pp = (TH1D*)hhdy[nbins]->Clone(tmpname);
0648   for(int k=1; k<=hhdy_pp->GetNbinsX(); k++) {
0649     if(hhdy_pp->GetBinLowEdge(k)<statscale_lowlim || (hhdy_pp->GetBinLowEdge(k)+hhdy_pp->GetBinWidth(k))>statscale_uplim) {
0650       hhdy_pp->SetBinContent(k,0.);
0651       hhdy_pp->SetBinError(k,0.);
0652     }
0653     else {
0654       double tmp = ppcorr * hhdy[nbins]->GetFunction("expo")->Eval(hhdy[nbins]->GetBinCenter(k));
0655       double tmprnd = myrandom->Poisson(tmp);
0656       if(tmprnd<0.) { tmprnd=0.; }
0657       hhdy_pp->SetBinContent(k,tmprnd);
0658       hhdy_pp->SetBinError(k,sqrt(tmprnd));
0659     }
0660   }
0661 
0662   sprintf(tmpname,"hhcorrbg_pp");
0663   hhcorrbg_pp = (TH1D*)hhbottom_pp->Clone(tmpname);
0664   hhcorrbg_pp->Add(hhdy_pp);
0665     hhcorrbg_pp->SetMarkerColor(kBlack);
0666     hhcorrbg_pp->SetLineColor(kBlack);
0667     hhbottom_pp->SetLineColor(kBlue);
0668     hhdy_pp->SetLineColor(kGreen+2);
0669   sprintf(tmpname,"hhall_pp");
0670   hhall_pp = (TH1D*)hhcorrbg_pp->Clone(tmpname);
0671   hhall_pp->Add(hhupspp);
0672     hhall_pp->SetLineColor(kMagenta);
0673     hhall_pp->SetMarkerColor(kMagenta);
0674 
0675 
0676 TCanvas* cbginpp = new TCanvas("cbginpp","corr bg in pp",10,10,700,700);
0677 
0678   hhcorrbg_pp->Fit("expo","rql","",statscale_lowlim,statscale_uplim);
0679   hhcorrbg_pp->GetFunction("expo")->SetLineColor(kBlack);
0680   hhbottom_pp->Fit("expo","rql","",statscale_lowlim,statscale_uplim);
0681   hhbottom_pp->GetFunction("expo")->SetLineColor(kBlue);
0682   hhdy_pp->Fit("expo","rql","",statscale_lowlim,statscale_uplim);
0683   hhdy_pp->GetFunction("expo")->SetLineColor(kGreen+2);
0684 
0685   hhall_pp->SetAxisRange(7.,12.);
0686   hhcorrbg_pp->Draw("pehist");
0687   hhbottom_pp->Draw("histsame");
0688   hhdy_pp->Draw("histsame");
0689 
0690 
0691 TCanvas* cpp = new TCanvas("cpp","corr bg + sig in pp",100,100,700,700);
0692   hhall_pp->SetAxisRange(7.,12.);
0693   hhall_pp->Draw("pehist");
0694   hhcorrbg_pp->Draw("pesame");
0695   hhbottom_pp->Draw("same");
0696   hhdy_pp->Draw("same");
0697 
0698 //-----------------------------------------------------------------------------------------
0699 // Combinatorial BG calculated for 10B central AuAu events
0700 //-----------------------------------------------------------------------------------------
0701 
0702 TCanvas* cdummy = new TCanvas("cdummy","cdummy",0,0,500,500);
0703 
0704 f = new TFile("fakee_eideff09.root");
0705   sprintf(tmpname,"hhfakefake_15"); // 15 is integrated over pT
0706   hhfakefake[nbins] = (TH1D*)f->Get(tmpname);
0707   hhfakefake[nbins]->SetDirectory(gROOT);
0708 f->Close();
0709 
0710 f = new TFile("crossterms_eideff09.root");
0711   sprintf(tmpname,"hhfakehf_15");
0712   hhfakehf[nbins] = (TH1D*)f->Get(tmpname);
0713   hhfakehf[nbins]->SetDirectory(gROOT);
0714 f->Close();
0715 
0716 TF1* fbg = new TF1("fbg","exp([0]+[1]*x)+exp([2]+[3]*x)",8.,11.);
0717 fbg->SetParameters(10., -1.0, 4., -0.1);
0718 fbg->SetParLimits(1.,-999.,0.);
0719 fbg->SetParLimits(3.,-999.,0.);
0720 
0721   sprintf(tmpname,"hhcombbg_15");
0722   hhcombbg[nbins] = (TH1D*)hhfakefake[nbins]->Clone(tmpname);
0723   hhcombbg[nbins]->Add(hhfakehf[nbins]);
0724   sprintf(tmpname,"hhcombbg_scaled_15");
0725   hhcombbg_scaled[nbins] = (TH1D*)hhcombbg[nbins]->Clone(tmpname);
0726   fbg->SetParameters(10., -1.0, 4., -0.1); 
0727   hhcombbg[nbins]->Fit(fbg,"qrl","",statscale_lowlim,statscale_uplim);
0728 
0729     for(int k=1; k<=hhcombbg[nbins]->GetNbinsX(); k++) {
0730       if(hhcombbg[nbins]->GetBinLowEdge(k)<statscale_lowlim || (hhcombbg[nbins]->GetBinLowEdge(k)+hhcombbg[nbins]->GetBinWidth(k))>statscale_uplim) {
0731         hhcombbg_scaled[nbins]->SetBinContent(k,0.);
0732         hhcombbg_scaled[nbins]->SetBinError(k,0.);
0733       }
0734       else {
0735         double tmp = statscale * hhcombbg[nbins]->GetFunction("fbg")->Eval(hhcombbg[nbins]->GetBinCenter(k));
0736         double tmprnd = myrandom->Poisson(tmp);
0737         if(tmprnd<0.) { tmprnd=0.; }
0738         hhcombbg_scaled[nbins]->SetBinContent(k,tmprnd);
0739         hhcombbg_scaled[nbins]->SetBinError(k,sqrt(tmprnd));
0740       }
0741     }
0742   hhcombbg_scaled[nbins]->Fit(fbg,"qrl","",statscale_lowlim,statscale_uplim);
0743 
0744 delete cdummy;
0745 
0746 TCanvas* C1 = new TCanvas("C1","Combinatorial BG Central Au+Au",100,100,600,600);
0747 C1->SetLogy();
0748   hhfakefake[nbins]->SetAxisRange(7.0,14.0);
0749   hhfakefake[nbins]->SetMinimum(0.1);
0750   hhfakefake[nbins]->SetMaximum(5000.);
0751   hhfakefake[nbins]->SetLineColor(kGreen+2);
0752   hhfakefake[nbins]->SetLineWidth(2);
0753   hhfakefake[nbins]->GetXaxis()->SetTitle("Transverse momentum [GeV/c]");
0754   hhfakefake[nbins]->GetXaxis()->SetTitleOffset(1.0);
0755   hhfakefake[nbins]->GetXaxis()->SetTitleColor(1);
0756   hhfakefake[nbins]->GetXaxis()->SetTitleSize(0.040);
0757   hhfakefake[nbins]->GetXaxis()->SetLabelSize(0.040);
0758   hhfakefake[nbins]->GetYaxis()->SetTitle("Combinatorial background");
0759   hhfakefake[nbins]->GetYaxis()->SetTitleOffset(1.3);
0760   hhfakefake[nbins]->GetYaxis()->SetTitleSize(0.040);
0761   hhfakefake[nbins]->GetYaxis()->SetLabelSize(0.040);
0762   hhfakefake[nbins]->Draw("e");
0763 
0764   hhfakehf[nbins]->SetLineColor(kOrange+4);
0765   hhfakehf[nbins]->SetLineWidth(2);
0766   hhfakehf[nbins]->Draw("esame");
0767 
0768   hhcombbg[nbins]->SetLineColor(kBlack);
0769   hhcombbg[nbins]->SetLineWidth(2);
0770   hhcombbg[nbins]->Draw("esame");
0771 
0772 TCanvas* C1sc = new TCanvas("C1sc","SCALED Combinatorial BG Central Au+Au",100,100,600,600);
0773 C1sc->SetLogy(); 
0774   hhcombbg_scaled[nbins]->SetAxisRange(7.,14.);
0775   hhcombbg_scaled[nbins]->Draw("esame");
0776 
0777 // Scaled Combinatorial background vs centrality
0778 
0779 for(int i=0; i<nbins; i++) {
0780 
0781   sprintf(tmpname,"hhcombbg_%d",i);
0782   hhcombbg_scaled[i] = (TH1D*)hhcombbg_scaled[nbins]->Clone(tmpname);
0783 
0784     for(int k=1; k<=hhcombbg_scaled[nbins]->GetNbinsX(); k++) {
0785       if(hhcombbg_scaled[nbins]->GetBinLowEdge(k)<statscale_lowlim || (hhcombbg_scaled[nbins]->GetBinLowEdge(k)+hhcombbg_scaled[nbins]->GetBinWidth(k))>statscale_uplim) {
0786         hhcombbg_scaled[i]->SetBinContent(k,0.);
0787         hhcombbg_scaled[i]->SetBinError(k,0.);
0788       }
0789       else {
0790         double tmp = Npionpairs[i] * hhcombbg_scaled[nbins]->GetFunction("fbg")->Eval(hhcombbg_scaled[nbins]->GetBinCenter(k));
0791         double tmprnd = myrandom->Poisson(tmp);
0792         if(tmprnd<0.) { tmprnd=0.; }
0793         hhcombbg_scaled[i]->SetBinContent(k,tmprnd);
0794         hhcombbg_scaled[i]->SetBinError(k,sqrt(tmprnd));
0795       }
0796     }
0797   hhcombbg_scaled[i]->Fit(fbg,"qrl","",statscale_lowlim,statscale_uplim);
0798 
0799 } // end over centrality bins
0800 
0801 TCanvas* c_comb_scaled = new TCanvas("c_comb_scaled","Combinatorial Background vs. Centrality",200,100,1200,600);
0802 c_comb_scaled->Divide(4,2);
0803 for(int i=0; i<nbins; i++) {
0804   c_comb_scaled->cd(i+1);
0805   hhcombbg_scaled[i]->SetAxisRange(8.5,11.0); hhcombbg_scaled[i]->SetMarkerStyle(1); hhcombbg_scaled[i]->Draw("pe");
0806   sprintf(tlchar,"%d-%d",10*i,10*(i+1));   tl[i] = new TLatex(9.0,hhcombbg_scaled[i]->GetMaximum()*0.9,tlchar); tl[i]->Draw();
0807 }
0808 
0809 
0810 //---------------------------------------------
0811 //  Au+Au Signal + all BG
0812 //---------------------------------------------
0813 
0814 for(int i=0; i<nbins; i++) {
0815   sprintf(tmpname,"hhtotbg_scaled_%d",i);
0816   hhtotbg_scaled[i] = (TH1D*)hhcombbg_scaled[i]->Clone(tmpname);
0817   hhtotbg_scaled[i]->Add(hhcorrbg_scaled[i]);
0818 }
0819 
0820 for(int i=0; i<nbins; i++) {
0821   sprintf(tmpname,"hhall_scaled_%d",i);
0822   hhall_scaled[i] = (TH1D*)hhtotbg_scaled[i]->Clone(tmpname);
0823   hhall_scaled[i]->Add(hhups[i]);
0824 }
0825 
0826 TCanvas* c000 = new TCanvas("c000","Au+Au Signal + All Background vs. Centrality",200,200,1200,600);
0827 c000->Divide(4,2);
0828 for(int i=0; i<nbins; i++) {
0829   c000->cd(i+1);
0830   hhall_scaled[i]->SetAxisRange(8.5,11.0); hhall_scaled[i]->SetMarkerStyle(1); hhall_scaled[i]->Draw("pehist");
0831   sprintf(tlchar,"%d-%d",10*i,10*(i+1));   tl[i] = new TLatex(9.0,hhall_scaled[i]->GetMaximum()*0.9,tlchar); tl[i]->Draw();
0832 }
0833 
0834 /*
0835 // Fit Signal + Correlated BG
0836 
0837 // hhfit[] is signal + correlated bg
0838 for(int i=0; i<nbins; i++) {
0839   sprintf(tmpname,"hhfit_%d",i);
0840   hhfit[i] = (TH1D*)hhall_scaled[i]->Clone(tmpname);
0841   for(int j=1; j<=hhall_scaled[i]->GetNbinsX(); j++) {
0842     hhfit[i]->SetBinContent(j,hhfit[i]->GetBinContent(j) - hhcombbg_scaled[i]->GetFunction("fbg")->Eval(hhcombbg_scaled[i]->GetBinCenter(j)));
0843     hhfit[i]->SetBinError(j,sqrt(pow(hhfit[i]->GetBinError(j),2)+hhcombbg_scaled[i]->GetFunction("fbg")->Eval(hhcombbg_scaled[i]->GetBinCenter(j))));
0844   }
0845 }
0846 
0847 //---------------------------------------------------------------------
0848 // plot signal + correlated bg for all pT
0849 //---------------------------------------------------------------------
0850 
0851   double tmppar0 = corrbgfitpar0+TMath::Log(statscale);
0852   double tmppar1 = corrbgfitpar1;
0853 cout << "###### " << tmppar0 << " " << tmppar1 << endl;
0854 
0855 double ppauauscale = fTCBauau->GetParameter(0)/fTCBpp->GetParameter(0)*0.9783;
0856 fSandBpp->SetParameter(0,fTCBpp->GetParameter(0)*ppauauscale);
0857 fSandBpp->SetParameter(1,fTCBpp->GetParameter(1));
0858 fSandBpp->SetParameter(2,fTCBpp->GetParameter(2));
0859 fSandBpp->SetParameter(3,fTCBpp->GetParameter(3));
0860 fSandBpp->SetParameter(4,fTCBpp->GetParameter(4));
0861 fSandBpp->SetParameter(5,fTCBpp->GetParameter(5)*ppauauscale);
0862 fSandBpp->SetParameter(6,fTCBpp->GetParameter(6)*ppauauscale);
0863 fSandBpp->SetParameter(7,tmppar0);
0864 fSandBpp->SetParameter(8,tmppar1);
0865 fSandBpp->SetLineColor(kBlue);
0866 fSandBpp->SetLineStyle(2);
0867 
0868 fSandBauau->SetParameter(0,fTCBauau->GetParameter(0));
0869 fSandBauau->SetParameter(1,fTCBauau->GetParameter(1));
0870 fSandBauau->SetParameter(2,fTCBauau->GetParameter(2));
0871 fSandBauau->SetParameter(3,fTCBauau->GetParameter(3));
0872 fSandBauau->SetParameter(4,fTCBauau->GetParameter(4));
0873 fSandBauau->SetParameter(5,fTCBauau->GetParameter(5));
0874 fSandBauau->SetParameter(6,fTCBauau->GetParameter(6));
0875 fSandBauau->SetParameter(7,tmppar0);
0876 fSandBauau->SetParameter(8,tmppar1);
0877 fSandBauau->SetLineColor(kRed);
0878 
0879 TCanvas* cfitall = new TCanvas("cfitall","Fit Central",270,270,600,600);
0880   hhfit[0]->SetAxisRange(7.0,14.);
0881   hhfit[0]->GetXaxis()->CenterTitle();
0882   hhfit[0]->GetXaxis()->SetTitle("Mass(e^{+}e^{-}) [GeV/c^2]");
0883   hhfit[0]->GetXaxis()->SetTitleOffset(1.1);
0884   hhfit[0]->GetXaxis()->SetLabelSize(0.045);
0885   hhfit[0]->GetYaxis()->CenterTitle();
0886   hhfit[0]->GetYaxis()->SetLabelSize(0.045);
0887   hhfit[0]->GetYaxis()->SetTitle("Events / (50 MeV/c^{2})");
0888   hhfit[0]->GetYaxis()->SetTitleOffset(1.5);
0889   hhfit[0]->Draw("pehist");
0890 //  fSandBpp->Draw("same");
0891 //  fSandBauau->Draw("same");
0892     //hhcorrbg_scaled[0]->SetLineColor(kRed);
0893     //hhcorrbg_scaled[0]->Scale(0.82);
0894     //hhcorrbg_scaled[0]->Scale(0.70);
0895     //hhcorrbg_scaled[0]->Draw("histsame");
0896 
0897   TF1* fmycorrbg = new TF1("fmycorrbg","exp([0]+[1]*x)",7.,14.);
0898   fmycorrbg->SetParameters(tmppar0,tmppar1);
0899   fmycorrbg->SetLineStyle(2);
0900   fmycorrbg->SetLineColor(kRed);
0901   fmycorrbg->Draw("same");
0902 
0903 //double myheight = fTCBauau->GetParameter(0);
0904 //TLatex* ld1 = new TLatex(10.1,myheight,"sPHENIX Simulation");
0905 //ld1->SetTextSize(0.035);
0906 //ld1->Draw();
0907 //TLatex* ld2 = new TLatex(10.1,myheight-100.,"0-10% Au+Au #sqrt{s} = 200 GeV");
0908 //ld2->SetTextSize(0.035);
0909 //ld2->Draw();
0910 
0911 TCanvas* cfitall2 = new TCanvas("cfitall2","FIT all pT",270,270,600,600);
0912 TH1D* hhfit_tmp = (TH1D*)hhfit[0]->Clone("hhfit_tmp");
0913 hhfit_tmp->SetAxisRange(8.0,11.);
0914 hhfit_tmp->Draw("pehist");
0915 //fSandBauau->Draw("same");
0916 fmycorrbg->Draw("same");
0917 */
0918 
0919 //----------------------------------------------------------------------
0920 
0921 // plot signal + both backgrounds for central Au+Au
0922 
0923 TCanvas* callpt = new TCanvas("callpt","Signal + All BG Central Au+Au",300,300,600,600);
0924 
0925   hhall_scaled[0]->GetXaxis()->SetTitle("Invariant mass GeV/c");
0926   hhall_scaled[0]->SetLineColor(kBlack);
0927   hhall_scaled[0]->SetMarkerColor(kBlack);
0928   hhall_scaled[0]->SetMarkerStyle(20);
0929   hhall_scaled[0]->SetAxisRange(8.0,10.8);
0930   hhall_scaled[0]->Draw("pehist");
0931     hhcombbg_scaled[0]->SetLineColor(kBlue);
0932     hhcombbg_scaled[0]->Draw("histsame");
0933       hhcorrbg_scaled[0]->SetLineColor(kRed);
0934       hhcorrbg_scaled[0]->Draw("histsame");
0935 
0936 
0937 
0938 //----------------------------------------------------------------
0939 //  Calculate RAA
0940 //----------------------------------------------------------------
0941 
0942 double u1start = 9.25;
0943 double u1stop  = 9.65;
0944 double u2start = 9.80;
0945 double u2stop  = 10.20;
0946 double u3start = 10.20;
0947 double u3stop  = 10.55;
0948 cout << "kuku2" << endl;
0949   double raa1[nbins+1],raa2[nbins+1],raa3[nbins+1],erraa1[nbins+1],erraa2[nbins+1],erraa3[nbins+1];
0950   for(int i=0; i<nbins; i++) { 
0951     //raa1[i] = raa1s[i];
0952     //raa2[i] = raa2s[i];
0953     //raa3[i] = raa3s[i];
0954     raa1[i] = grRAA1S->Eval(Npart[i]);
0955     raa2[i] = grRAA2S->Eval(Npart[i]);
0956     if(i<5) { raa3[i] = grRAA2S->Eval(Npart[i])/2.; }
0957       else { raa3[i] = (grRAA2S->Eval(Npart[i])+grRAA3S->Eval(Npart[i]))/2.; }
0958   }
0959   int fbin1 = hhall_scaled[0]->FindBin(u1start + 0.001);
0960   int lbin1 = hhall_scaled[0]->FindBin(u1stop  - 0.001);
0961   int fbin2 = hhall_scaled[0]->FindBin(u2start + 0.001);
0962   int lbin2 = hhall_scaled[0]->FindBin(u2stop  - 0.001);
0963   int fbin3 = hhall_scaled[0]->FindBin(u3start + 0.001);
0964   int lbin3 = hhall_scaled[0]->FindBin(u3stop  - 0.001);
0965   cout << "Y(1S) bin range: " << fbin1 << " - " << lbin1 << endl;
0966   cout << "Y(1S) inv. mass range: " << u1start << " - " << u1stop << endl;
0967   cout << "Y(2S) bin range: " << fbin2 << " - " << lbin2 << endl;
0968   cout << "Y(2S) inv. mass range: " << u2start << " - " << u2stop << endl;
0969   cout << "Y(3S) bin range: " << fbin3 << " - " << lbin3 << endl;
0970   cout << "Y(3S) inv. mass range: " << u3start << " - " << u3stop << endl;
0971 
0972   double sum1[99]   = {0.};
0973   double truesum1[99]   = {0.};
0974   double ersum1[99] = {0.};
0975   double sumpp1   = 0.;
0976   double ersumpp1 = 0.;
0977   double sum2[99]   = {0.};
0978   double truesum2[99]   = {0.};
0979   double ersum2[99] = {0.};
0980   double sumpp2   = 0.;
0981   double ersumpp2 = 0.;
0982   double sum3[99]   = {0.};
0983   double truesum3[99]   = {0.};
0984   double ersum3[99] = {0.};
0985   double sumpp3   = 0.;
0986   double ersumpp3 = 0.;
0987 
0988   double sumsum1[99]   = {0.};
0989   double sumsum2[99]   = {0.};
0990   double sumsum3[99]   = {0.};
0991   double sumsum1pp   = 0.;
0992   double sumsum2pp   = 0.;
0993   double sumsum3pp   = 0.;
0994 
0995   for(int j=fbin1; j<=lbin1; j++) {
0996     sumpp1   += hhups1pp->GetBinContent(j);
0997     ersumpp1 += hhupspp->GetBinError(j)*hhupspp->GetBinError(j);
0998   }
0999   for(int j=fbin2; j<=lbin2; j++) {
1000     sumpp2   += hhups2pp->GetBinContent(j);
1001     ersumpp2 += hhupspp->GetBinError(j)*hhupspp->GetBinError(j);
1002   }
1003   for(int j=fbin3; j<=lbin3; j++) {
1004     sumpp3   += hhups3pp->GetBinContent(j);
1005     ersumpp3 += hhupspp->GetBinError(j)*hhupspp->GetBinError(j);
1006   }
1007 
1008   for(int i=0; i<nbins; i++) {
1009 
1010     //double s1 = double(i); double s2 = double(i+1); if(i==nbins) {s1 = 0.;}
1011 
1012     for(int j=fbin1; j<=lbin1; j++) {
1013       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)));
1014       truesum1[i] += hhups1[i]->GetBinContent(j);
1015       ersum1[i] += hhall_scaled[i]->GetBinError(j)*hhall_scaled[i]->GetBinError(j);
1016     }
1017       sumsum1[i]     = truesum1[i];                   // direct count in mass range
1018       sumsum1pp   = sumpp1;       
1019       //sumsum1[i]     = hhups1[i]->GetEntries();       // total number of upsilons in pT bin (rounded up)
1020       //sumsum1pp   = hhups1pp->GetEntries();
1021       //sumsum1[i]     = Nups1*fUpsilonPt->Integral(s1,s2)/upsnorm;  // total number of upsilons in pT bin
1022       //sumsum1pp   = Nups1pp*fUpsilonPt->Integral(s1,s2)/upsnorm;
1023 
1024         if(sumsum1[i]>0. && sumsum1pp>0.) {
1025           erraa1[i] = raa1[i]*sqrt(ersum1[i]/sumsum1[i]/sumsum1[i] + ersumpp1/sumsum1pp/sumsum1pp);
1026         } else {raa1[i]=-1.0; erraa1[i] = 999.; }
1027 
1028     for(int j=fbin2; j<=lbin2; j++) {
1029       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)));
1030       truesum2[i] += hhups2[i]->GetBinContent(j);
1031       ersum2[i] += hhall_scaled[i]->GetBinError(j)*hhall_scaled[i]->GetBinError(j);
1032     }
1033       sumsum2[i]     = truesum2[i];                   // direct count in mass range
1034       sumsum2pp   = sumpp2;       
1035       //sumsum2[i]     = hhups2[i]->GetEntries();       // total number of upsilons in pT bin (rounded up)
1036       //sumsum2pp   = hhups2pp->GetEntries();
1037       //sumsum2[i]     = Nups2*fUpsilonPt->Integral(s1,s2)/upsnorm;  // total number of upsilons in pT bin
1038       //sumsum2pp   = Nups2pp*fUpsilonPt->Integral(s1,s2)/upsnorm;
1039       
1040         if(sumsum2[i]>0. && sumsum2pp>0.) {
1041           erraa2[i] = raa2[i]*sqrt(ersum2[i]/sumsum2[i]/sumsum2[i] + ersumpp2/sumsum2pp/sumsum2pp);
1042         } else {raa2[i]=-1.0; erraa2[i] = 999.; }
1043 
1044     for(int j=fbin3; j<=lbin3; j++) {
1045       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)));
1046       truesum3[i] += hhups3[i]->GetBinContent(j);
1047       ersum3[i] += hhall_scaled[i]->GetBinError(j)*hhall_scaled[i]->GetBinError(j);
1048     }
1049       sumsum3[i]     = truesum3[i];                   // direct count in mass range
1050       sumsum3pp   = sumpp3;       
1051       //sumsum3[i]     = hhups3[i]->GetEntries();       // total number of upsilons in pT bin (rounded up)
1052       //sumsum3pp   = hhups3pp->GetEntries();
1053       //sumsum3[i]     = Nups3*fUpsilonPt->Integral(s1,s2)/upsnorm;   // total number of upsilons in pT bin
1054       //sumsum3pp   = Nups3pp*fUpsilonPt->Integral(s1,s2)/upsnorm;
1055       
1056         if(truesum3[i]>0. && sumpp3>0.) {
1057           erraa3[i] = raa3[i]*sqrt(ersum3[i]/sumsum3[i]/sumsum3[i] + ersumpp3/sumsum3pp/sumsum3pp);
1058         } else {raa3[i]=-1.0; erraa3[i] = 999.; }
1059   
1060   } // end loop over centrality
1061 
1062   erraa3[3] = erraa3[3]*1.2;
1063 
1064   for(int i=0; i<nbins; i++) {
1065     cout << "Npart, Raa = " << Npart[i] << " " << raa1[i] << " " << raa2[i] << " " << raa3[i] << endl;
1066   }
1067 
1068 cout << "====== Y(1S):" << endl;
1069   for(int i=0; i<nbins; i++) {
1070     cout << "   " << i << " " << sumsum1[i] << "(" << Nups1[i] << ")" << " +- " << sqrt(ersum1[i]) 
1071          << " \t\t pp: " << sumsum1pp << " +- " << sqrt(ersumpp1) << endl;
1072   }
1073 cout << "====== Y(2S):" << endl;
1074   for(int i=0; i<nbins; i++) {
1075     cout << "   " << i << " " << sumsum2[i] << "(" << Nups2[i] << ")" << " +- " << sqrt(ersum2[i]) 
1076          << " \t\t pp: " << sumsum2pp << " +- " << sqrt(ersumpp2) << endl;
1077   }
1078 cout << "====== Y(3S):" << endl;
1079   for(int i=0; i<nbins; i++) {
1080     cout << "   " << i << " " << sumsum3[i] << "(" << Nups3[i] << ")" << " +- " << sqrt(ersum3[i]) 
1081          << " \t\t pp: " << sumsum3pp << " +- " << sqrt(ersumpp3) << endl;
1082   }
1083 
1084 
1085 /*
1086 double raa2_rebin[9],raapt_rebin2[9],erraa2_rebin[9];
1087 double raa3_rebin[9],raapt_rebin3[9],erraa3_rebin[9];
1088 double sum2_rebin[9],ersum2_rebin[9],sum2pp_rebin[9],ersumpp2_rebin[9];
1089 double sum3_rebin[9],ersum3_rebin[9],sum3pp_rebin[9],ersumpp3_rebin[9];
1090 
1091 // Rebin Y(2S) by 2
1092 raapt_rebin2[0] = 1.;
1093 raapt_rebin2[1] = 3.;
1094 raapt_rebin2[2] = 5.;
1095 raapt_rebin2[3] = 7.;
1096 raa2_rebin[0] = grRAA2S->Eval(raapt_rebin2[0]);
1097 raa2_rebin[1] = grRAA2S->Eval(raapt_rebin2[1]);
1098 raa2_rebin[2] = grRAA2S->Eval(raapt_rebin2[2]);
1099 raa2_rebin[3] = grRAA2S->Eval(raapt_rebin2[3]);
1100 sum2_rebin[0] = truesum2[0]+truesum2[1];
1101 sum2_rebin[1] = truesum2[2]+truesum2[3];
1102 sum2_rebin[2] = truesum2[4]+truesum2[5];
1103 sum2_rebin[3] = truesum2[6]+truesum2[7]+truesum2[8]+truesum2[9];
1104 ersum2_rebin[0] = ersum2[0]+ersum2[1];
1105 ersum2_rebin[1] = ersum2[2]+ersum2[3];
1106 ersum2_rebin[2] = ersum2[4]+ersum2[5];
1107 ersum2_rebin[3] = ersum2[6]+ersum2[7]+ersum2[8]+ersum2[9];
1108 sum2pp_rebin[0] = sumpp2[0]+sumpp2[1];
1109 sum2pp_rebin[1] = sumpp2[2]+sumpp2[3];
1110 sum2pp_rebin[2] = sumpp2[4]+sumpp2[5];
1111 sum2pp_rebin[3] = sumpp2[6]+sumpp2[7]+sumpp2[8]+sumpp2[9];
1112 ersumpp2_rebin[0] = ersumpp2[0]+ersumpp2[1];
1113 ersumpp2_rebin[1] = ersumpp2[2]+ersumpp2[3];
1114 ersumpp2_rebin[2] = ersumpp2[4]+ersumpp2[5];
1115 ersumpp2_rebin[3] = ersumpp2[6]+ersumpp2[7]+ersumpp2[8]+ersumpp2[9];
1116   erraa2_rebin[0] = raa2[0]*sqrt(ersum2_rebin[0]/sum2_rebin[0]/sum2_rebin[0] + ersumpp2_rebin[0]/sum2pp_rebin[0]/sum2pp_rebin[0]);
1117   erraa2_rebin[1] = raa2[1]*sqrt(ersum2_rebin[1]/sum2_rebin[1]/sum2_rebin[1] + ersumpp2_rebin[1]/sum2pp_rebin[1]/sum2pp_rebin[1]);
1118   erraa2_rebin[2] = raa2[2]*sqrt(ersum2_rebin[2]/sum2_rebin[2]/sum2_rebin[2] + ersumpp2_rebin[2]/sum2pp_rebin[2]/sum2pp_rebin[2]);
1119   erraa2_rebin[3] = raa2[3]*sqrt(ersum2_rebin[3]/sum2_rebin[3]/sum2_rebin[3] + ersumpp2_rebin[3]/sum2pp_rebin[3]/sum2pp_rebin[3]);
1120 // Rebin Y(3S) by 4
1121 raapt_rebin3[0] = 2.;
1122 raapt_rebin3[1] = 6.;
1123 raa3_rebin[0] = grRAA3S->Eval(raapt_rebin3[0]);
1124 raa3_rebin[1] = grRAA3S->Eval(raapt_rebin3[1]);
1125 sum3_rebin[0] = truesum3[0]+truesum3[1]+truesum3[2]+truesum3[3];
1126 sum3_rebin[1] = truesum3[4]+truesum3[5]+truesum3[6]+truesum3[7]+truesum3[8]+truesum3[9];
1127 ersum3_rebin[0] = ersum3[0]+ersum3[1]+ersum3[2]+ersum3[3];
1128 ersum3_rebin[1] = ersum3[4]+ersum3[5]+ersum3[6]+ersum3[7]+ersum3[8]+ersum3[9];
1129 sum3pp_rebin[0] = sumpp3[0]+sumpp3[1]+sumpp3[3]+sumpp3[3];
1130 sum3pp_rebin[1] = sumpp3[4]+sumpp3[5]+sumpp3[6]+sumpp3[7]+sumpp3[8]+sumpp3[9];
1131 ersumpp3_rebin[0] = ersumpp3[0]+ersumpp3[1]+ersumpp3[2]+ersumpp3[3];
1132 ersumpp3_rebin[1] = ersumpp3[4]+ersumpp3[5]+ersumpp3[6]+ersumpp3[7]+ersumpp3[8]+ersumpp3[9];
1133   erraa3_rebin[0] = raa3[0]*sqrt(ersum3_rebin[0]/sum3_rebin[0]/sum3_rebin[0] + ersumpp3_rebin[0]/sum3pp_rebin[0]/sum3pp_rebin[0]);
1134   erraa3_rebin[1] = raa3[1]*sqrt(ersum3_rebin[1]/sum3_rebin[1]/sum3_rebin[1] + ersumpp3_rebin[1]/sum3pp_rebin[1]/sum3pp_rebin[1]);
1135 // Rebin Y(3S) by 8
1136 */
1137 /*
1138 raapt_rebin3[0] = 5.;
1139 raa3_rebin[0] = grRAA3S->Eval(raapt_rebin3[0]);
1140 sum3_rebin[0] = truesum3[0]+truesum3[1]+truesum3[2]+truesum3[3]+
1141                 truesum3[4]+truesum3[5]+truesum3[6]+truesum3[7]+truesum3[8]+truesum3[9];
1142 ersum3_rebin[0] = ersum3[0]+ersum3[1]+ersum3[2]+ersum3[3]+
1143                   ersum3[4]+ersum3[5]+ersum3[6]+ersum3[7]+ersum3[8]+ersum3[9];
1144 sum3pp_rebin[0] = sumpp3[0]+sumpp3[1]+sumpp3[3]+sumpp3[3]+
1145                   sumpp3[4]+sumpp3[5]+sumpp3[6]+sumpp3[7]+sumpp3[8]+sumpp3[9];
1146 ersumpp3_rebin[0] = ersumpp3[0]+ersumpp3[1]+ersumpp3[2]+ersumpp3[3]+
1147                     ersumpp3[4]+ersumpp3[5]+ersumpp3[6]+ersumpp3[7]+ersumpp3[8]+ersumpp3[9];
1148 erraa3_rebin[0] = raa3[0]*sqrt(ersum3_rebin[0]/sum3_rebin[0]/sum3_rebin[0] + ersumpp3_rebin[0]/sum3pp_rebin[0]/sum3pp_rebin[0]);
1149 */
1150 
1151 //-------------------------------------------------
1152 //  Plot RAA
1153 //-------------------------------------------------
1154 
1155 int npts1 = nbins;
1156 int npts2 = nbins;
1157 int npts3 = nbins;
1158 int npts2_rebin = 4;
1159 int npts3_rebin = 2;
1160 
1161 TCanvas* craa = new TCanvas("craa","R_{AA}",120,120,800,600);
1162 TH2F* hh2 = new TH2F("hh2"," ",10,0.,400.,10,0.,1.1);
1163 hh2->GetXaxis()->SetTitle("N_{part}");
1164 hh2->GetXaxis()->SetTitleOffset(0.9);
1165 hh2->GetXaxis()->SetTitleColor(1);
1166 hh2->GetXaxis()->SetTitleSize(0.050);
1167 hh2->GetXaxis()->SetLabelSize(0.040);
1168 hh2->GetYaxis()->SetTitle("R_{AA}");
1169 hh2->GetYaxis()->SetTitleOffset(0.7);
1170 hh2->GetYaxis()->SetTitleSize(0.050);
1171 hh2->GetYaxis()->SetLabelSize(0.040);
1172 hh2->Draw();
1173 
1174 double xx1[nbins+1]; for(int i=0; i<nbins; i++) {xx1[i] = Npart[i];}
1175 double xx2[nbins+1]; for(int i=0; i<nbins; i++) {xx2[i] = Npart[i] - 1.;}
1176 double xx3[nbins+1]; for(int i=0; i<nbins; i++) {xx3[i] = Npart[i] + 1.;}
1177 //double xx3_rebin[nbins+1]; for(int i=0; i<npts3_rebin; i++) {xx3_rebin[i] = raapt_rebin3[i];}
1178 
1179 TGraphErrors* gr1 = new TGraphErrors(npts1,xx1,raa1,0,erraa1);
1180 gr1->SetMarkerStyle(20);
1181 gr1->SetMarkerColor(kBlack);
1182 gr1->SetLineColor(kBlack);
1183 gr1->SetLineWidth(2);
1184 gr1->SetMarkerSize(1.5);
1185 gr1->SetName("gr1");
1186 gr1->Draw("p");
1187 
1188 TGraphErrors* gr2 = new TGraphErrors(npts2,xx2,raa2,0,erraa2);
1189 gr2->SetMarkerStyle(20);
1190 gr2->SetMarkerColor(kRed);
1191 gr2->SetLineColor(kRed);
1192 gr2->SetLineWidth(2);
1193 gr2->SetMarkerSize(1.5);
1194 gr2->SetName("gr2");
1195 gr2->Draw("p");
1196 
1197 //TGraphErrors* gr2_rebin = new TGraphErrors(npts2_rebin,xx2_rebin,raa2_rebin,0,erraa2_rebin);
1198 //gr2_rebin->SetMarkerStyle(20);
1199 //gr2_rebin->SetMarkerColor(kRed);
1200 //gr2_rebin->SetLineColor(kRed);
1201 //gr2_rebin->SetLineWidth(2);
1202 //gr2_rebin->SetMarkerSize(1.5);
1203 //gr2_rebin->SetName("gr2");
1204 //gr2_rebin->Draw("p");
1205 
1206 //--- 3S state -------------------
1207 // double dummy[9]; for(int i=0; i<9; i++) {dummy[i]=-99.;}
1208  TGraphErrors* gr3 = new TGraphErrors(nbins,xx3,raa3,0,erraa3);
1209  gr3->SetMarkerStyle(20);
1210  gr3->SetMarkerColor(kBlue);
1211  gr3->SetLineColor(kBlue);
1212  gr3->SetLineWidth(2);
1213  gr3->SetMarkerSize(1.5);
1214  gr3->SetName("gr3");
1215 // gr3->Draw("p");
1216 /*
1217  TGraphErrors* gr3_rebin = new TGraphErrors(npts3_rebin,xx3_rebin,raa3_rebin,0,erraa3_rebin);
1218  gr3_rebin->SetMarkerStyle(20);
1219  gr3_rebin->SetMarkerColor(kBlue);
1220  gr3_rebin->SetLineColor(kBlue);
1221  gr3_rebin->SetLineWidth(2);
1222  gr3_rebin->SetMarkerSize(1.5);
1223  gr3_rebin->SetName("gr3");
1224  gr3_rebin->Draw("p");
1225 */
1226 /*
1227  TGraphErrors* gr3_rebin = new TGraphErrors(1,xx3_rebin,raa3_rebin,0,erraa3_rebin);
1228  gr3_rebin->SetMarkerStyle(20);
1229  gr3_rebin->SetMarkerColor(kBlue);
1230  gr3_rebin->SetLineColor(kBlue);
1231  gr3_rebin->SetLineWidth(2);
1232  gr3_rebin->SetMarkerSize(1.5);
1233  gr3_rebin->SetName("gr3");
1234  gr3_rebin->Draw("p");
1235 */
1236 //TArrow* aa[9];
1237 //TLine* ll[9];
1238 //erraa3[3] = erraa3[3]*1.5;
1239 //erraa3[4] = erraa3[4]*0.6;
1240 //erraa3[4] = erraa3[4]*1.5;
1241 //erraa3[6] = erraa3[6]*1.2;
1242 /*
1243 for(int i=0; i<npts3; i++) {
1244   aa[i] = new TArrow(xx3[i],1.64*erraa3[i],xx3[i],-0.02,0.02);
1245   aa[i]->SetLineColor(kBlue);
1246   aa[i]->SetLineWidth(2);
1247   aa[i]->Draw();
1248   ll[i] = new TLine(xx3[i]-0.15,1.64*erraa3[i],xx3[i]+0.15,1.64*erraa3[i]);
1249   ll[i]->SetLineColor(kBlue);
1250   ll[i]->SetLineWidth(2);
1251   ll[i]->Draw();
1252 }
1253 */
1254 //--- end 3S state --------------
1255 
1256 //TLegend *leg = new TLegend(0.70,0.70,0.88,0.88);
1257 TLegend *leg = new TLegend(0.73,0.76,0.89,0.88);
1258 leg->SetBorderSize(0);
1259 leg->SetFillColor(10);
1260 leg->SetFillStyle(1001);
1261 TLegendEntry *entry1=leg->AddEntry("gr1","Y(1S)","p");
1262 TLegendEntry *entry2=leg->AddEntry("gr2","Y(2S)","p");
1263 //TLegendEntry *entry3=leg->AddEntry("gr3","Y(3S)","l");
1264 //TLegendEntry *entry3=leg->AddEntry("gr3","Y(3S)","p");
1265 leg->Draw();
1266 
1267 TLatex* l1 = new TLatex(155.,1.02,"#font[72]{sPHENIX} Projection"); l1->SetTextFont(42); l1->Draw();
1268 //TLatex* l11 = new TLatex(20.,0.90,"0-10% cent. Au+Au, Years 1-3"); l11->SetTextFont(42); l11->Draw();
1269 TLatex* l2 = new TLatex(155.,0.93,"21 nb^{-1} rec. Au+Au"); l2->SetTextFont(42); l2->Draw();
1270 TLatex* l3 = new TLatex(155.,0.84,"62 pb^{-1} samp. #it{p+p}"); l3->SetTextFont(42); l3->Draw();
1271 
1272 TLine* lll = new TLine(0.6,0.64,1.3,0.64);
1273 lll->SetLineColor(kBlue);
1274 lll->SetLineWidth(2);
1275 //lll->Draw();
1276 
1277 grRAA1S->Draw("l");
1278 grRAA1S_eta1->Draw("l");
1279 grRAA1S_eta3->Draw("l");
1280 grRAA2S->Draw("l");
1281 grRAA2S_eta1->Draw("l");
1282 grRAA2S_eta3->Draw("l");
1283 //grRAA3S->Draw("l");
1284 //grRAA3S_eta1->Draw("l");
1285 //grRAA3S_eta3->Draw("l");
1286 
1287 //==================================================================================
1288 
1289 /*
1290 fout = new TFile("test.root","recreate");
1291 for(int i=0; i<nbins+1; i++) {
1292   sprintf(tmpname,"hhfit_%d",i);
1293   hhfit[i]->Write();
1294   hhups1pp[i]->Write();
1295   hhups2pp[i]->Write();
1296   hhups3pp[i]->Write();
1297   hhupspp[i]->Write();
1298   hhups1[i]->Write();
1299   hhups2[i]->Write();
1300   hhups3[i]->Write();
1301   hhups[i]->Write();
1302   hhall_pp[i]->Write();
1303   hhcorrbg_pp[i]->Write();
1304   hhcorrbg_scaled[i]->Write();
1305 }
1306 fout->Write();
1307 fout->Close();
1308 */
1309 
1310 }
1311
[ 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!