Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Draw_HFJetTruth_InvMass.C
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file Draw_HFJetTruth_InvMass.C
1 // $Id: $
2 
11 #include <TFile.h>
12 #include <TGraphAsymmErrors.h>
13 #include <TGraphErrors.h>
14 #include <TLatex.h>
15 #include <TLine.h>
16 #include <TString.h>
17 #include <TTree.h>
18 #include <cassert>
19 #include <cmath>
20 #include "SaveCanvas.C"
21 #include "sPhenixStyle.C"
22 
23 TFile *_file0 = NULL;
24 TTree *T = NULL;
25 
26 void Draw_HFJetTruth_InvMass(const TString infile =
27  // "/sphenix/user/jinhuang/HF-jet/event_gen/200pp_pythia8_CTEQ6L_20GeV/200pp_pythia8_CTEQ6L_20GeV_ALL.cfg_eneg_DSTReader.root",
28  // double int_lumi = 210715 / 5.533e-05 / 1e9, const double dy = 0.6 * 2)
29 
30  "/sphenix/user/jinhuang/HF-jet/event_gen/200pp_pythia8_CTEQ6L_7GeV/200pp_pythia8_CTEQ6L_7GeV_ALL.cfg_eneg_DSTReader.root",
31  double int_lumi = 891093 / 3.332e-02 / 1e9, const double dy = 0.6 * 2)
32 
33 //"/sphenix/user/jinhuang/HF-jet/event_gen/200pp_pythia8_CTEQ6L/200pp_pythia8_CTEQ6L_111ALL.cfg_eneg_DSTReader.root",
34 //double int_lumi = 789908/4.631e-03 / 1e9, const double dy = 0.6*2)
35 
36 //Draw_HFJetTruth_InvMass(const TString infile =
37 // "../macros3/G4sPHENIXCells.root_DSTReader.root",
38 // double int_lumi = 789908/4.631e-03 / 1e9, const double dy = 0.6*2)
39 //Draw_HFJetTruth_InvMass(const TString infile =
40 // "../macros2/G4sPHENIXCells.root_DSTReader.root",
41 // const double int_lumi = 1842152 / 5.801e-03 / 1e9,
42 // const double dy = 0.6 * 2)
43 {
44  SetsPhenixStyle();
45  gStyle->SetOptStat(0);
46  gStyle->SetOptFit(1111);
47  TVirtualFitter::SetDefaultFitter("Minuit2");
48 
49  gSystem->Load("libg4eval.so");
50  gSystem->Load("libg4dst.so");
51 
52  if (!_file0)
53  {
54  TString chian_str = infile;
55  chian_str.ReplaceAll("ALL", "*");
56 
57  TChain *t = new TChain("T");
58  const int n = t->Add(chian_str);
59 
60  int_lumi *= n;
61  // cout << "Loaded " << n << " root files with " << chian_str << endl;
62  cout << "int_lumi = " << int_lumi << " pb^-1 from " << n << " root files with " << chian_str << endl;
63  assert(n > 0);
64 
65  T = t;
66 
67  _file0 = new TFile;
68  _file0->SetName(infile);
69  }
70 
71  // step1: get the cross section
72  // DrawCrossSection(int_lumi, dy);
73 
74  // step2: ploting
75 
76  // Draw_HFJetTruth_InvMass_DrawCrossSection_PR(infile);
77  // CrossSection2RAA_Proposal(infile);
78  // CrossSection2RAA(infile);
79  const double acceptance1 = 2. * (1.1 - .4);
80  // CrossSection2RAA(infile, false, acceptance1);
81  CrossSection2RAA(infile, false, acceptance1, true);
82 
83  // const double acceptance2 = 2.* (0.85 - .4);
84  // CrossSection2RAA(infile, false, acceptance2);
85 
86  // CrossSection2v2(infile, false, .7, acceptance);
87  // CrossSection2v2(infile, false, .4, acceptance);
88 }
89 
90 void DrawCrossSection(double int_lumi, const double dy)
91 {
92  TCut basicDijetEventQA("(AntiKt_Truth_r04[n_AntiKt_Truth_r04-1].get_pt())>15 && (AntiKt_Truth_r04[n_AntiKt_Truth_r04-2].get_pt())>10 && abs(AntiKt_Truth_r04[n_AntiKt_Truth_r04-1].get_eta())<.6 && abs(AntiKt_Truth_r04[n_AntiKt_Truth_r04-2].get_eta())<.6 && cos(AntiKt_Truth_r04[n_AntiKt_Truth_r04-1].get_phi() - AntiKt_Truth_r04[n_AntiKt_Truth_r04-2].get_phi())<-.5");
93 
94  cout << "Build event selection of " << (const char *) basicDijetEventQA << endl;
95 
96  T->Draw(">>EventList", basicDijetEventQA);
97  TEventList *elist = gDirectory->GetObjectChecked("EventList", "TEventList");
98  assert(elist);
99  cout << elist->GetN() << " / " << T->GetEntriesFast() << " events selected" << endl;
100  T->SetEventList(elist);
101 
102  T->SetAlias("DiJetInvMass", "sqrt( (AntiKt_Truth_r04[n_AntiKt_Truth_r04-1].get_e() + AntiKt_Truth_r04[n_AntiKt_Truth_r04-2].get_e())**2 - (AntiKt_Truth_r04[n_AntiKt_Truth_r04-1].get_px() + AntiKt_Truth_r04[n_AntiKt_Truth_r04-2].get_px())**2 - (AntiKt_Truth_r04[n_AntiKt_Truth_r04-1].get_py() + AntiKt_Truth_r04[n_AntiKt_Truth_r04-2].get_py())**2 - (AntiKt_Truth_r04[n_AntiKt_Truth_r04-1].get_pz() + AntiKt_Truth_r04[n_AntiKt_Truth_r04-2].get_pz())**2 )");
103 
104  TH1 *hall = new TH1F("hall", ";m_{12} [GeV]", 200, 0, 200);
105  TH1 *h_b = new TH1F("h_b", ";m_{12} [GeV]", 200, 0, 200);
106  TH1 *h_bq = new TH1F("h_bq", ";m_{12} [GeV]", 200, 0, 200);
107  TH1 *h_bh = new TH1F("h_bh", ";m_{12} [GeV]", 200, 0, 200);
108  TH1 *h_bh5 = new TH1F("h_bh5", ";m_{12} [GeV]", 200, 0, 200);
109  TH1 *h_ch5 = new TH1F("h_ch5", ";m_{12} [GeV]", 200, 0, 200);
110  TH1 *h_c = new TH1F("h_c", ";m_{12} [GeV]", 200, 0, 200);
111 
112  cout << "DiJetInvMass>>hall" << endl;
113  T->Draw("DiJetInvMass>>hall",
114  "",
115  "goff");
116 
117  cout << "DiJetInvMass>>h_b" << endl;
118  T->Draw("DiJetInvMass>>h_b",
119  "abs(AntiKt_Truth_r04[n_AntiKt_Truth_r04-1].get_property(1000))==5 && abs(AntiKt_Truth_r04[n_AntiKt_Truth_r04-2].get_property(1000))==5",
120  "goff");
121  // T->Draw(
122  // "AntiKt_Truth_r04.get_pt()* AntiKt_Truth_r04.get_property(1001) >>h_bq",
123  // "AntiKt_Truth_r04.get_pt()>15 && abs(AntiKt_Truth_r04.get_eta())<0.6 && abs(AntiKt_Truth_r04.get_property(1000))==5",
124  // "goff");
125 
126  // T->Draw("AntiKt_Truth_r04.get_pt()>>h_bh",
127  // "AntiKt_Truth_r04.get_pt()>15 && abs(AntiKt_Truth_r04.get_eta())<0.6 && AntiKt_Truth_r04.get_property(1010)==5",
128  // "goff");
129 
130  // T->Draw("AntiKt_Truth_r04.get_pt()>>h_bh5",
131  // "AntiKt_Truth_r04.get_pt()>15 && abs(AntiKt_Truth_r04.get_eta())<0.6 && abs(AntiKt_Truth_r04.get_property(1010))==5 && AntiKt_Truth_r04.get_property(1011) * AntiKt_Truth_r04.get_pt() > 5",
132  // "goff");
133  //
134  // T->Draw("AntiKt_Truth_r04.get_pt()>>h_c",
135  // "AntiKt_Truth_r04.get_pt()>15 && abs(AntiKt_Truth_r04.get_eta())<0.6 && abs(AntiKt_Truth_r04.get_property(1000))==4",
136  // "goff");
137  // T->Draw("AntiKt_Truth_r04.get_pt()>>h_ch5",
138  // "AntiKt_Truth_r04.get_pt()>15 && abs(AntiKt_Truth_r04.get_eta())<0.6 && abs(AntiKt_Truth_r04.get_property(1010))==4 && AntiKt_Truth_r04.get_property(1011) * AntiKt_Truth_r04.get_pt() > 5",
139  // "goff");
140 
141  Convert2CrossSection(hall, int_lumi, dy);
142  Convert2CrossSection(h_b, int_lumi, dy);
143  // Convert2CrossSection(h_bq, int_lumi, dy);
144  // Convert2CrossSection(h_bh5, int_lumi, dy);
145  // Convert2CrossSection(h_ch5, int_lumi, dy);
146  // Convert2CrossSection(h_c, int_lumi, dy);
147 
148  h_b->SetLineColor(kBlue);
149  h_b->SetMarkerColor(kBlue);
150 
151  h_bq->SetLineColor(kBlue + 3);
152  h_bq->SetMarkerColor(kBlue + 3);
153 
154  h_bh5->SetLineColor(kBlue + 3);
155  h_bh5->SetMarkerColor(kBlue + 3);
156  h_bh5->SetMarkerStyle(kStar);
157 
158  h_c->SetLineColor(kRed);
159  h_c->SetMarkerColor(kRed);
160 
161  h_ch5->SetLineColor(kRed + 3);
162  h_ch5->SetMarkerColor(kRed + 3);
163  h_ch5->SetMarkerStyle(kStar);
164 
165  // TGraphAsymmErrors *gr_fonll_b = GetFONLL_B();
166  // gr_fonll_b->SetFillColor(kBlue - 7);
167  // gr_fonll_b->SetFillStyle(3002);
168  // TGraphAsymmErrors *gr_fonll_c = GetFONLL_C();
169  // gr_fonll_c->SetFillColor(kRed - 7);
170  // gr_fonll_c->SetFillStyle(3003);
171  // TGraph *gr_phenix = GetPHENIX_jet();
172  // gr_phenix->SetLineColor(kGray + 2);
173  // gr_phenix->SetMarkerColor(kBlack);
174  // gr_phenix->SetMarkerStyle(kOpenCross);
175 
176  TCanvas *c1 = new TCanvas("Draw_HFJetTruth_InvMass_DrawCrossSection", "Draw_HFJetTruth_InvMass_DrawCrossSection", 1000, 860);
177  // c1->Divide(2, 2);
178  int idx = 1;
179  TPad *p;
180 
181  p = (TPad *) c1->cd(idx++);
182  c1->Update();
183  p->SetGridx(0);
184  p->SetGridy(0);
185  p->SetLogy();
186 
187  hall->Draw();
188 
189  h_b->Draw("same");
190  // h_bh5->Draw("same");
191  // h_bq->Draw("same");
192  // h_c->Draw("same");
193  // h_ch5->Draw("same");
194 
195  // hall->GetXaxis()->SetRangeUser(30, 160);
196  hall->GetYaxis()->SetTitle(
197  "d^{3}#sigma/(d#eta_{1}d#eta_{2}dm_{12}) [pb/(GeV)]");
198 
199  TLegend *leg = new TLegend(0.45, 0.6, 0.95, 0.93);
200  leg->SetFillColor(kWhite);
201  leg->SetFillStyle(1001);
202  leg->SetHeader("#splitline{#it{#bf{sPHENIX}} fast simulation, #it{p}+#it{p} #sqrt{s} = 200 GeV}{|#eta_{1,2}|<0.6, |#Delta#phi_{1,2}|>2#pi/3, p_{T,1}>15 GeV/c, p_{T,2}>10 GeV/c}");
203  leg->AddEntry(hall, "Inclusive jet, PYTHIA8, Truth, anti-k_{t}, R=0.4",
204  "lpe");
205  // leg->AddEntry(h_c, "c-quark jet, Pythia8, Truth, anti-k_{t}, R=0.4", "lpe");
206  leg->AddEntry(h_b, "b-quark jet, PYTHIA8, Truth, anti-k_{t}, R=0.4", "lpe");
207  // leg->AddEntry(h_bq, "Leading b-quark in jet, Pythia8, anti-k_{t}, R=0.4", "lpe");
208  // leg->AddEntry(h_bh5,
209  // "b-hadron jet, Pythia8, Truth, anti-k_{t}, R=0.4, p_{T, b-hadron}>5 GeV/c",
210  // "lpe");
211  // leg->AddEntry(gr_phenix, "PHENIX inclusive jet, PRL 116, 122301 (2016)", "ple");
212  // leg->AddEntry(gr_fonll_c, "c-quark, FONLL v1.3.2, CTEQ6.6, |y|<0.6", "f");
213  // leg->AddEntry(gr_fonll_b, "b-quark, FONLL v1.3.2, CTEQ6.6, |y|<0.6", "f");
214  leg->Draw();
215 
216  SaveCanvas(c1, TString(_file0->GetName()) + TString(c1->GetName()), kTRUE);
217 }
218 
219 void Draw_HFJetTruth_InvMass_DrawCrossSection_PR(const TString infile)
220 {
221  TFile *f = TFile::Open(infile + "Draw_HFJetTruth_InvMass_DrawCrossSection.root");
222  assert(f);
223 
224  TH1F *hall = (TH1F *) f->GetObjectChecked("hall", "TH1F");
225  assert(hall);
226  TH1F *h_b = (TH1F *) f->GetObjectChecked("h_b", "TH1F");
227  assert(h_b);
228 
229  hall->SetMarkerColor(kBlack);
230  hall->SetLineColor(kBlack);
231  hall->SetMarkerStyle(kFullCircle);
232 
233  h_b->SetMarkerColor(kBlue + 2);
234  h_b->SetLineColor(kBlue + 2);
235  h_b->SetMarkerStyle(kFullCircle);
236 
237  TGraph *gr_star = GetSTAR2019();
238  gr_star->SetLineColor(kGray + 2);
239  gr_star->SetLineWidth(3);
240  gr_star->SetMarkerColor(kGray + 2);
241  gr_star->SetMarkerSize(2);
242  gr_star->SetMarkerStyle(kFullSquare);
243 
244  TCanvas *c1 = new TCanvas("Draw_HFJetTruth_InvMass_DrawCrossSection_PR", "Draw_HFJetTruth_InvMass_DrawCrossSection_PR", 1000, 860);
245  // c1->Divide(2, 2);
246  int idx = 1;
247  TPad *p;
248 
249  p = (TPad *) c1->cd(idx++);
250  c1->Update();
251  p->SetLogy();
252 
253  TH1 *hframe = p->DrawFrame(35, 1e-2, 140, 1e6);
254  hframe->SetTitle(";m_{12} [GeV/c^{2}];d^{3}#sigma/(d#eta_{1}d#eta_{2}dm_{12}) [pb/(GeV/c^{2})]");
255 
256 
257  gr_star->Draw("pe");
258  hall->Draw("same");
259  h_b->Draw("same");
260 
261  TLegend *leg = new TLegend(0.2, 0.7, 0.95, 0.92);
262  leg->SetFillColor(kWhite);
263  leg->SetFillStyle(1001);
264  // leg->SetHeader("#splitline{#it{#bf{sPHENIX }} Simulation}{#it{p}+#it{p}, #sqrt{s} = 200 GeV, |#eta|<0.6}");
265  leg->SetHeader("#splitline{#it{#bf{sPHENIX}} fast simulation, #it{p}+#it{p} #sqrt{s} = 200 GeV}{|#eta_{1,2}|<0.6, |#Delta#phi_{1,2}|>2#pi/3, p_{T,1}>15 GeV/c, p_{T,2}>10 GeV/c}");
266  leg->AddEntry("", "",
267  "");
268  leg->AddEntry(hall, "Inclusive jet, PYTHIA8 + CTEQ6L, anti-k_{T} R=0.4",
269  "lpe");
270  leg->AddEntry(gr_star, "STAR, PRD95 (2017) no.7, 071103, R=0.6", "ple");
271  leg->AddEntry(h_b, "#it{b}-quark jet, PYTHIA8 + CTEQ6L, anti-k_{T} R=0.4", "lpe");
272 
273  leg->Draw();
274 
275  SaveCanvas(c1, TString(_file0->GetName()) + TString(c1->GetName()), kTRUE);
276 }
277 
278 void CrossSection2RAA(const TString infile, const bool use_AA_jet_trigger = true, const double dy = .7 * 2, const bool b3yr = false)
279 {
280  TFile *f = TFile::Open(infile + "Draw_HFJetTruth_InvMass_DrawCrossSection.root");
281  assert(f);
282 
283  TH1F *hall = (TH1F *) f->GetObjectChecked("hall", "TH1F");
284  assert(hall);
285  TH1F *h_b = (TH1F *) f->GetObjectChecked("h_b", "TH1F");
286  assert(h_b);
287 
288  TString s_suffix(use_AA_jet_trigger ? "_AAJetTriggered" : "");
289  s_suffix += b3yr ? "_3yr" : "";
290  s_suffix += Form("_deta%.2f", dy / 2);
291 
292  const double b_jet_RAA = 0.4;
293  const double target_RAA_ratio = 4;;
294 
295  const double pp_eff = 0.6;
296  const double pp_purity = 0.4;
297  const double AuAu_eff = 0.4;
298  const double AuAu_purity = 0.4;
299 
301  // 5-year lumi in [sPH-TRG-000]
303 
304  const double pp_lumi = b3yr ? 62 : 200; // pb^-1 [sPH-TRG-000], rounded up from 197 pb^-1
305  const double pp_inelastic_crosssec = 42e-3 / 1e-12; // 42 mb in pb [sPH-TRG-000]
306 
307  const double AuAu_MB_Evt = use_AA_jet_trigger ? (b3yr ? 320e9 : 550e9) : (b3yr ? 142e9 : 240e9); // [sPH-TRG-000], depending on whether jet trigger applied in AA collisions
308  const double pAu_MB_Evt = 600e9; // [sPH-TRG-000]
309 
310  const double AuAu_Ncoll_C0_10 = 960.2; // [DOI:?10.1103/PhysRevC.87.034911?]
311  const double AuAu_Ncoll_C0_20 = 770.6; // [DOI:?10.1103/PhysRevC.91.064904?]
312  const double AuAu_Ncoll_C0_100 = 250; // pb^-1 [sPH-TRG-000]
313  const double pAu_Ncoll_C0_100 = 4.7; // pb^-1 [sPH-TRG-000]
314 
315  const double AuAu_eq_lumi_C0_10 = AuAu_MB_Evt * 0.1 * AuAu_Ncoll_C0_10 / pp_inelastic_crosssec; //
316  const double AuAu_eq_lumi_C0_20 = AuAu_MB_Evt * 0.2 * AuAu_Ncoll_C0_20 / pp_inelastic_crosssec; //
317  const double AuAu_eq_lumi_C0_100 = AuAu_MB_Evt * 1 * AuAu_Ncoll_C0_100 / pp_inelastic_crosssec; //
318 
319  const double pAu_eq_lumi_C0_100 = pAu_MB_Evt * 1 * pAu_Ncoll_C0_100 / pp_inelastic_crosssec; //
320 
321  cout << "CrossSection2RAA integrated luminosity assumptions in pb^-1: " << endl;
322  cout << "\t"
323  << "pp_lumi = " << pp_lumi << endl;
324  cout << "\t"
325  << "AuAu_eq_lumi_C0_10 = " << AuAu_eq_lumi_C0_10 << endl;
326  cout << "\t"
327  << "AuAu_eq_lumi_C0_20 = " << AuAu_eq_lumi_C0_20 << endl;
328  cout << "\t"
329  << "AuAu_eq_lumi_C0_100 = " << AuAu_eq_lumi_C0_100 << endl;
330  cout << "\t"
331  << "pAu_eq_lumi_C0_100 = " << pAu_eq_lumi_C0_100 << endl;
332 
333  TCanvas *c1 = new TCanvas("Draw_HFJetTruth_InvMass_CrossSection2RAA_Ratio" + s_suffix, "Draw_HFJetTruth_InvMass_CrossSection2RAA_Ratio" + s_suffix, 1100, 800);
334  c1->Divide(1, 1);
335  int idx = 1;
336  TPad *p;
337 
338  p = (TPad *) c1->cd(idx++);
339  c1->Update();
340  // p->SetGridx(0);
341  // p->SetGridy(0);
342  // p->SetLogy();
343 
344  // in sPH-HF-2017-001, dijet purity efficiency was conservatively assumed to be the product of that of two single b-jets,
345  // i.e. without benifit of the enhanace of 2nd b-jet abundance after the 1st jet in the event is tagged as a b-jet
346  TH1 *g_pp = CrossSection2RelUncert(h_b, 1., dy, pp_lumi * pp_eff * pp_purity * pp_eff * pp_purity);
347  TH1 *g_AA = CrossSection2RelUncert(h_b, b_jet_RAA, dy, AuAu_eq_lumi_C0_10 * AuAu_eff * AuAu_purity * AuAu_eff * AuAu_purity);
348 
349  g_pp->SetLineColor(kRed);
350  g_AA->SetLineColor(kBlue);
351 
352  g_pp->Draw();
353  g_AA->Draw("same");
354  SaveCanvas(c1, infile + "_" + TString(c1->GetName()), kTRUE);
355 
356  TCanvas *c1 = new TCanvas("Draw_HFJetTruth_InvMass_CrossSection2RAA" + s_suffix, "Draw_HFJetTruth_InvMass_CrossSection2RAA" + s_suffix, 1100, 800);
357  c1->Divide(1, 1);
358  int idx = 1;
359  TPad *p;
360 
361  p = (TPad *) c1->cd(idx++);
362  c1->Update();
363 
364  p->DrawFrame(30, 0, 75, 1.3)->SetTitle(";Di-jet invariant mass [GeV/#it{c}^{2}];#it{R}^{bb}_{AA}");
365 
366  TGraphErrors *ge_RAA = GetRAA(g_pp, g_AA);
367 
368  ge_RAA->SetLineWidth(4);
369  ge_RAA->SetMarkerStyle(kFullCircle);
370  ge_RAA->SetMarkerSize(3);
371 
372  ge_RAA->Draw("pe");
373  ge_RAA->Print();
374 
375  TLegend *leg = new TLegend(.0, .7, .85, .93);
376  leg->SetFillStyle(0);
377  leg->AddEntry("", "#it{#bf{sPHENIX}} Projection, #it{b}-jet, 0-10% Au+Au, Year 1-3", "");
378  leg->AddEntry("", Form("|#eta_{1,2}|<%.1f, |#Delta#phi_{1,2}|>2#pi/3, p_{T,1}>15 GeV/c, p_{T,2}>10 GeV/c", dy / 2), "");
379  leg->AddEntry("", Form("#it{p}+#it{p}: %.0f pb^{-1} samp., %.0f%% Eff., %.0f%% Pur.", pp_lumi, pp_eff * 100, pp_purity * 100), "");
380  leg->AddEntry("", Form("Au+Au: %.0fnb^{-1} rec., %.0f%% Eff., %.0f%% Pur.", '%', AuAu_MB_Evt/6.8252 / 1e9, AuAu_eff * 100, AuAu_purity * 100), "");
381  leg->Draw();
382 
383 // TLegend *leg2 = new TLegend(.17, .70, .88, .77);
384 // leg2->AddEntry(ge_RAA, "#it{b}-jet #it{R}_{AA}, Au+Au 0-10%C, #sqrt{s_{NN}}=200 GeV", "pl");
385 // leg2->Draw();
386 
387  SaveCanvas(c1, infile + "_" + TString(c1->GetName()), kTRUE);
388 
389  TCanvas *c1 = new TCanvas("Draw_HFJetTruth_InvMass_CrossSection2RAA_Theory" + s_suffix, "Draw_HFJetTruth_InvMass_CrossSection2RAA_Theory" + s_suffix, 1100, 800);
390  c1->Divide(1, 1);
391  int idx = 1;
392  TPad *p;
393 
394  p = (TPad *) c1->cd(idx++);
395  c1->Update();
396 
397  p->DrawFrame(35, 0, 75, 1.4)->SetTitle(";Di-jet invariant mass [GeV];#it{R}^{bb}_{AA}");
398 
399  //
400  TGraph *g18 = GetRAAKang2019(1, 1.8);
401  TGraph *g20 = GetRAAKang2019(1, 2.0);
402  TGraph *g22 = GetRAAKang2019(1, 2.2);
403  //
404  g18->SetLineColor(kBlue + 3);
405  g20->SetLineColor(kGreen + 3);
406  g22->SetLineColor(kRed + 3);
407  g18->SetLineWidth(3);
408  g20->SetLineWidth(3);
409  g22->SetLineWidth(3);
410  //
411  g18->Draw("l");
412  g20->Draw("l");
413  g22->Draw("l");
414  //
415  ge_RAA->DrawClone("pe");
416  leg->DrawClone();
417 // leg2->DrawClone();
418  //
419  TLegend *leg1 = new TLegend(.19, .6, .92, .7);
420 // leg1->AddEntry("", "Kang et al, PRD #bf{99}, 034006 (2019), m=m_{b}, rad.+col.", "");
421  leg1->SetHeader("Kang et al, PRD #bf{99}, 034006 (2019), m=m_{b}, rad.+col.");
422  leg1->Draw();
423  TLegend *leg11 = new TLegend(.35, .5, .55, .62);
424  leg11->AddEntry(g18, "g_{med} = 1.8", "l");
425  leg11->AddEntry(g20, "g_{med} = 2.0", "l");
426  leg11->AddEntry(g22, "g_{med} = 2.2", "l");
427  leg11->Draw();
428 
429  SaveCanvas(c1, infile + "_" + TString(c1->GetName()), kTRUE);
430 
431 
432  TCanvas *c1 = new TCanvas("Draw_HFJetTruth_InvMass_CrossSection2RAARatio_Theory" + s_suffix, "Draw_HFJetTruth_InvMass_CrossSection2RAARatio_Theory" + s_suffix, 1100, 800);
433  c1->Divide(1, 1);
434  int idx = 1;
435  TPad *p;
436 
437  p = (TPad *) c1->cd(idx++);
438  c1->Update();
439 
440  p->DrawFrame(35, 0, 75, 15)->SetTitle(";Di-jet invariant mass [GeV];#it{R}^{bb}_{AA}/#it{R}^{jj}_{AA}");
441  TLine * l = new TLine(35, 1, 75, 1);
442 // l->SetLineStyle(kDashed);
443  l->Draw();
444 
445  //
446  TGraph *g18 = GetRAARatioKang2019(1, 1.8);
447  TGraph *g20 = GetRAARatioKang2019(1, 2.0);
448  TGraph *g22 = GetRAARatioKang2019(1, 2.2);
449  //
450  g18->SetLineColor(kBlue + 3);
451  g20->SetLineColor(kGreen + 3);
452  g22->SetLineColor(kRed + 3);
453  g18->SetLineWidth(3);
454  g20->SetLineWidth(3);
455  g22->SetLineWidth(3);
456  //
457  g18->Draw("l");
458  g20->Draw("l");
459  g22->Draw("l");
460  //
461  TGraphErrors * ge_RAARatio = (TGraphErrors *) ge_RAA->DrawClone("pe");
462  assert(ge_RAARatio);
463  for (int bin = 0; bin<ge_RAARatio->GetN();++bin)
464  {
465  ge_RAARatio->GetY()[bin] *= target_RAA_ratio/b_jet_RAA;
466  ge_RAARatio->GetEY()[bin] *= target_RAA_ratio/b_jet_RAA;
467  }
468 
469  leg->DrawClone();
470 // leg2->DrawClone();
471  //
472  TLegend *leg1 = new TLegend(.2, .6, .92, .70);
473 // leg1->AddEntry("", "Kang et al, PRD #bf{99}, 034006 (2019), m=m_{b}, rad.+col.", "");
474  leg1->SetHeader("Kang et al, PRD #bf{99}, 034006 (2019), m=m_{b}, rad.+col.");
475  leg1->Draw();
476  TLegend *leg11 = new TLegend(.35, .5, .55, .62);
477  leg11->AddEntry(g18, "g_{med} = 1.8", "l");
478  leg11->AddEntry(g20, "g_{med} = 2.0", "l");
479  leg11->AddEntry(g22, "g_{med} = 2.2", "l");
480  leg11->Draw();
481 
482  SaveCanvas(c1, infile + "_" + TString(c1->GetName()), kTRUE);
483 
484 
485 }
486 
487 TGraphErrors *GetRAA(TH1 *h_pp, TH1 *h_AA)
488 {
489  int n_bin = 0;
490 
491  double xs[1000] = {0};
492  double ys[1000] = {0};
493  double eys[1000] = {0};
494 
495  assert(h_pp);
496  assert(h_AA);
497  assert(h_pp->GetNbinsX() == h_AA->GetNbinsX());
498 
499  for (int i = 1; i <= h_pp->GetNbinsX(); ++i)
500  {
501  if (h_pp->GetBinError(i) > 0 && h_pp->GetBinError(i) < 1 && h_AA->GetBinError(i) > 0 && h_AA->GetBinError(i) < 1)
502  {
503  xs[n_bin] = h_pp->GetXaxis()->GetBinCenter(i);
504 
505  ys[n_bin] = h_AA->GetBinContent(i) / h_pp->GetBinContent(i);
506 
507  eys[n_bin] = ys[n_bin] * sqrt(pow(h_AA->GetBinError(i) / h_AA->GetBinContent(i), 2) + pow(h_pp->GetBinError(i) / h_pp->GetBinContent(i), 2));
508 
509  // final uncertainty cut
510  if (eys[n_bin] < .4)
511  n_bin += 1;
512  }
513  }
514 
515  TGraphErrors *ge = new TGraphErrors(n_bin, xs, ys, NULL, eys);
516  ge->SetName(TString("RAA_") + h_AA->GetName());
517 
518  return ge;
519 }
520 
521 TH1 *CrossSection2RelUncert(const TH1F *h_cross,
522  const double suppression,
523  const double deta,
524  const double pp_quiv_int_lum)
525 {
526  assert(h_cross);
527  TH1 *
528  h_ratio = (TH1 *)
529  h_cross->Clone(TString(h_cross->GetName()) + Form("_copy%d", rand()));
530 
531  //convert to count per bin
532  h_ratio->Scale(deta * deta * h_ratio->GetXaxis()->GetBinWidth(0) * pp_quiv_int_lum * suppression);
533 
534  // define the x-binning
535  h_ratio->Rebin(5);
536 
537  for (int i = 1; i <= h_ratio->GetNbinsX(); ++i)
538  {
539  const double yield = h_ratio->GetBinContent(i);
540 
541  if (yield > 0)
542  {
543  h_ratio->SetBinContent(i, suppression);
544 
545  h_ratio->SetBinError(i, suppression / sqrt(yield));
546  }
547  else
548  {
549  h_ratio->SetBinContent(i, 0);
550 
551  h_ratio->SetBinError(i, 0);
552  }
553  }
554 
555  h_ratio->GetYaxis()->SetTitle("Relative Cross Section and Uncertainty");
556 
557  return h_ratio;
558 }
559 
560 void Convert2CrossSection(TH1 *h, const double int_lumi, const double dy)
561 {
562  h->Sumw2();
563 
564  for (int i = 1; i <= h->GetXaxis()->GetNbins(); ++i)
565  {
566  const double m1 = h->GetXaxis()->GetBinLowEdge(i);
567  const double m2 = h->GetXaxis()->GetBinUpEdge(i);
568 
569  // const double dpT2 = pT2*pT2 - pT1*pT1;
570  const double dm = m2 - m1;
571 
572  // const double count2sigma = 1./dpT2/dy/int_lumi;
573  const double count2sigma = 1. / dm / int_lumi / dy / dy;
574 
575  h->SetBinContent(i, h->GetBinContent(i) * count2sigma);
576  h->SetBinError(i, h->GetBinError(i) * count2sigma);
577  }
578 
579  // h->GetYaxis()->SetTitle("#sigma/(dp_{T}^{2} dy) (pb/(GeV/c)^{2})");
580 
581  h->SetMarkerStyle(kFullCircle);
582 }
583 
584 TGraphAsymmErrors *GetSTAR2019(const TString hepdata = "/sphenix/user/jinhuang/HF-jet/event_gen/HEPData-ins1493842-v1-Table_4.root")
585 {
587  // Measurement of the cross section and longitudinal double-spin asymmetry for di-jet production in polarized pp collisions at s = 200 GeV
588  //
589  // The STAR collaboration
590  // Adamczyk, L. , Adkins, J.K. , Agakishiev, G. , Aggarwal, M.M. , Ahammed, Z. , Alekseev, I. , Anderson, D.M. , Aoyama, R. , Aparin, A. , Arkhipkin, D.
591  // No Journal Information, 2016
592  // https://doi.org/10.17182/hepdata.77208
593 
594  TFile *_file0 = TFile::Open(hepdata);
595  _file0->cd("Table 4");
596  TGraphAsymmErrors *gae = (TGraphAsymmErrors *)
597  gDirectory->GetObjectChecked("Graph1D_y1", "TGraphAsymmErrors");
598  assert(gae);
599 
600  for (int bin = 0; bin < gae->GetN(); ++bin)
601  {
602  (gae->GetY())[bin] *= 1e6; // ub -> pb
603  gae->SetPointEYhigh(bin, gae->GetErrorYhigh(bin) * 1e6); // ub -> pb
604  gae->SetPointEYlow(bin, gae->GetErrorYlow(bin) * 1e6); // ub -> pb
605  }
606  gae->SetName("gaeSTARDiJetMass");
607 
608  return gae;
609 }
610 
611 TGraph *GetRAAKang2019(const int mass = 1, const double g = 1.8)
612 {
613  // m12 2mb,g=1.8 1mb,g=1.8 2mb,g=2.0 1mb,g=2.0 2mb,g=2.2 1mb,g=2.2
614  // 17 0.959318 0.813631 0.937975 0.718446 0.904169 0.583889
615  // 19 0.924265 0.734466 0.887045 0.618448 0.831249 0.469604
616  // 21 0.889933 0.676922 0.838959 0.550413 0.765768 0.397484
617  // 23 0.852022 0.624044 0.78657 0.489034 0.696008 0.335986
618  // 25 0.813713 0.572249 0.736696 0.435365 0.634859 0.287682
619  // 27 0.829553 0.588906 0.752964 0.444663 0.647904 0.287249
620  // 29 0.787303 0.554276 0.702613 0.412568 0.592341 0.26188
621  // 31 0.767304 0.519823 0.674582 0.374897 0.555677 0.22774
622  // 33 0.774209 0.531997 0.68111 0.383097 0.559728 0.2299
623  // 35 0.733919 0.503672 0.633725 0.359385 0.509993 0.213934
624  // 37 0.763963 0.528822 0.665819 0.373186 0.537606 0.214043
625  // 39 0.741506 0.564112 0.645253 0.416396 0.52631 0.253036
626  // 41 0.74574 0.497485 0.638774 0.336391 0.500389 0.179804
627  // 43 0.725435 0.496459 0.617035 0.336753 0.479825 0.179778
628  // 45 0.714991 0.472252 0.599281 0.309914 0.453767 0.15839
629  // 47 0.680055 0.443568 0.559106 0.287112 0.413458 0.145234
630  // 49 0.679489 0.450226 0.557737 0.28911 0.410027 0.142904
631  // 51 0.652643 0.420073 0.525378 0.262171 0.375742 0.124363
632  // 53 0.648852 0.421612 0.520426 0.260902 0.368284 0.121749
633  // 55 0.659723 0.424185 0.522925 0.255923 0.36165 0.114938
634  // 57 0.608096 0.394726 0.474039 0.23987 0.324512 0.109253
635  // 59 0.599733 0.385379 0.462397 0.228792 0.309607 0.100507
636  // 61 0.602495 0.398467 0.466165 0.238347 0.312724 0.103957
637  // 63 0.58698 0.369703 0.441268 0.20886 0.281604 0.0844847
638  // 65 0.58804 0.39896 0.449849 0.238734 0.296312 0.102839
639  // 67 0.568857 0.369444 0.423411 0.208794 0.264683 0.0833302
640  // 69 0.569308 0.374558 0.421737 0.210743 0.259612 0.0834087
641  // 71 0.570937 0.379757 0.418615 0.213217 0.253412 0.0841172
642  // 73 0.571426 0.368965 0.404248 0.201111 0.234995 0.0773808
643  // 75 0.444747 0.285571 0.308362 0.153737 0.17623 0.0575338
644  // 77 0.518565 0.353508 0.368773 0.198289 0.215024 0.0772062
645  // 79 0.52314 0.354647 0.36598 0.194516 0.207028 0.0732502
646  // 81 0.462486 0.306947 0.313368 0.163098 0.169892 0.0593907
647  // 83 0.494343 0.338822 0.337503 0.182332 0.182904 0.0665247
648  // 85 0.477053 0.321068 0.315372 0.167256 0.163585 0.0589566
649  // 87 0.425877 0.287827 0.276778 0.14819 0.140167 0.0513555
650  // 89 0.449382 0.310267 0.291343 0.16052 0.146033 0.0555645
651  // 91 0.426898 0.294991 0.270933 0.150014 0.131644 0.0508301
652  // 93 0.58802 0.451299 0.406374 0.248587 0.210967 0.0874347
653  // 95 0.367013 0.245387 0.216773 0.11698 0.0963793 0.0374309
654  // 97 0.363587 0.251356 0.213642 0.120293 0.0932188 0.0381195
655  // 99 0.349846 0.240864 0.197647 0.111828 0.0824488 0.0345796
656  // 101 0.294855 0.20189 0.158814 0.0907104 0.0634189 0.0274697
657  // 103 0.264038 0.18211 0.137505 0.0801353 0.0531663 0.0239006
658 
659  const int n_Data = 44;
660 
661  const double m12[] = {
662  17,
663  19,
664  21,
665  23,
666  25,
667  27,
668  29,
669  31,
670  33,
671  35,
672  37,
673  39,
674  41,
675  43,
676  45,
677  47,
678  49,
679  51,
680  53,
681  55,
682  57,
683  59,
684  61,
685  63,
686  65,
687  67,
688  69,
689  71,
690  73,
691  75,
692  77,
693  79,
694  81,
695  83,
696  85,
697  87,
698  89,
699  91,
700  93,
701  95,
702  97,
703  99,
704  101,
705  103};
706 
707  const double RAA_2mb_g18[] = {
708  0.959318,
709  0.924265,
710  0.889933,
711  0.852022,
712  0.813713,
713  0.829553,
714  0.787303,
715  0.767304,
716  0.774209,
717  0.733919,
718  0.763963,
719  0.741506,
720  0.74574,
721  0.725435,
722  0.714991,
723  0.680055,
724  0.679489,
725  0.652643,
726  0.648852,
727  0.659723,
728  0.608096,
729  0.599733,
730  0.602495,
731  0.58698,
732  0.58804,
733  0.568857,
734  0.569308,
735  0.570937,
736  0.571426,
737  0.444747,
738  0.518565,
739  0.52314,
740  0.462486,
741  0.494343,
742  0.477053,
743  0.425877,
744  0.449382,
745  0.426898,
746  0.58802,
747  0.367013,
748  0.363587,
749  0.349846,
750  0.294855,
751  0.264038};
752 
753  const double RAA_1mb_g_18[] = {
754  0.813631,
755  0.734466,
756  0.676922,
757  0.624044,
758  0.572249,
759  0.588906,
760  0.554276,
761  0.519823,
762  0.531997,
763  0.503672,
764  0.528822,
765  0.564112,
766  0.497485,
767  0.496459,
768  0.472252,
769  0.443568,
770  0.450226,
771  0.420073,
772  0.421612,
773  0.424185,
774  0.394726,
775  0.385379,
776  0.398467,
777  0.369703,
778  0.39896,
779  0.369444,
780  0.374558,
781  0.379757,
782  0.368965,
783  0.285571,
784  0.353508,
785  0.354647,
786  0.306947,
787  0.338822,
788  0.321068,
789  0.287827,
790  0.310267,
791  0.294991,
792  0.451299,
793  0.245387,
794  0.251356,
795  0.240864,
796  0.20189,
797  0.18211};
798 
799  const double RAA_2mb_g_20[] = {
800  0.937975,
801  0.887045,
802  0.838959,
803  0.78657,
804  0.736696,
805  0.752964,
806  0.702613,
807  0.674582,
808  0.68111,
809  0.633725,
810  0.665819,
811  0.645253,
812  0.638774,
813  0.617035,
814  0.599281,
815  0.559106,
816  0.557737,
817  0.525378,
818  0.520426,
819  0.522925,
820  0.474039,
821  0.462397,
822  0.466165,
823  0.441268,
824  0.449849,
825  0.423411,
826  0.421737,
827  0.418615,
828  0.404248,
829  0.308362,
830  0.368773,
831  0.36598,
832  0.313368,
833  0.337503,
834  0.315372,
835  0.276778,
836  0.291343,
837  0.270933,
838  0.406374,
839  0.216773,
840  0.213642,
841  0.197647,
842  0.158814,
843  0.137505};
844 
845  const double RAA_1mb_g_20[] = {
846  0.718446,
847  0.618448,
848  0.550413,
849  0.489034,
850  0.435365,
851  0.444663,
852  0.412568,
853  0.374897,
854  0.383097,
855  0.359385,
856  0.373186,
857  0.416396,
858  0.336391,
859  0.336753,
860  0.309914,
861  0.287112,
862  0.28911,
863  0.262171,
864  0.260902,
865  0.255923,
866  0.23987,
867  0.228792,
868  0.238347,
869  0.20886,
870  0.238734,
871  0.208794,
872  0.210743,
873  0.213217,
874  0.201111,
875  0.153737,
876  0.198289,
877  0.194516,
878  0.163098,
879  0.182332,
880  0.167256,
881  0.14819,
882  0.16052,
883  0.150014,
884  0.248587,
885  0.11698,
886  0.120293,
887  0.111828,
888  0.0907104,
889  0.0801353};
890 
891  const double RAA_2mb_g_22[] = {
892  0.904169,
893  0.831249,
894  0.765768,
895  0.696008,
896  0.634859,
897  0.647904,
898  0.592341,
899  0.555677,
900  0.559728,
901  0.509993,
902  0.537606,
903  0.52631,
904  0.500389,
905  0.479825,
906  0.453767,
907  0.413458,
908  0.410027,
909  0.375742,
910  0.368284,
911  0.36165,
912  0.324512,
913  0.309607,
914  0.312724,
915  0.281604,
916  0.296312,
917  0.264683,
918  0.259612,
919  0.253412,
920  0.234995,
921  0.17623,
922  0.215024,
923  0.207028,
924  0.169892,
925  0.182904,
926  0.163585,
927  0.140167,
928  0.146033,
929  0.131644,
930  0.210967,
931  0.0963793,
932  0.0932188,
933  0.0824488,
934  0.0634189,
935  0.0531663};
936 
937  const double RAA_1mb_g_22[] = {
938  0.583889,
939  0.469604,
940  0.397484,
941  0.335986,
942  0.287682,
943  0.287249,
944  0.26188,
945  0.22774,
946  0.2299,
947  0.213934,
948  0.214043,
949  0.253036,
950  0.179804,
951  0.179778,
952  0.15839,
953  0.145234,
954  0.142904,
955  0.124363,
956  0.121749,
957  0.114938,
958  0.109253,
959  0.100507,
960  0.103957,
961  0.0844847,
962  0.102839,
963  0.0833302,
964  0.0834087,
965  0.0841172,
966  0.0773808,
967  0.0575338,
968  0.0772062,
969  0.0732502,
970  0.0593907,
971  0.0665247,
972  0.0589566,
973  0.0513555,
974  0.0555645,
975  0.0508301,
976  0.0874347,
977  0.0374309,
978  0.0381195,
979  0.0345796,
980  0.0274697,
981  0.0239006};
982 
983  TGraph *gRAAKang2019(NULL);
984 
985  if (mass == 1)
986  {
987  if (g == 1.8)
988  {
989  gRAAKang2019 = new TGraph(n_Data, m12, RAA_1mb_g_18);
990  }
991  else if (g == 2)
992  {
993  gRAAKang2019 = new TGraph(n_Data, m12, RAA_1mb_g_20);
994  }
995  else if (g == 2.2)
996  {
997  gRAAKang2019 = new TGraph(n_Data, m12, RAA_1mb_g_22);
998  }
999  }
1000  else if (mass == 2)
1001  {
1002  if (g == 1.8)
1003  {
1004  gRAAKang2019 = new TGraph(n_Data, m12, RAA_2mb_g_18);
1005  }
1006  else if (g == 2)
1007  {
1008  gRAAKang2019 = new TGraph(n_Data, m12, RAA_2mb_g_20);
1009  }
1010  else if (g == 2.2)
1011  {
1012  gRAAKang2019 = new TGraph(n_Data, m12, RAA_2mb_g_22);
1013  }
1014  }
1015 
1016  cout <<"mass = "<<mass<<", g="<<g<<endl;
1017  assert(gRAAKang2019);
1018  gRAAKang2019->Print();
1019 
1020  return gRAAKang2019;
1021 }
1022 
1023 TGraph *GetRAARatioKang2019(const int mass = 1, const double g = 1.8)
1024 {
1025  // m12 g=1.8 g=2.0 g=2.2
1026  // 17 4.10724 7.40351 16.6273
1027  // 19 4.44519 7.83383 18.0237
1028  // 21 4.57977 7.8564 17.5948
1029  // 23 4.29797 7.07926 15.2511
1030  // 25 3.58437 5.58706 11.445
1031  // 27 3.81735 5.93123 12.1819
1032  // 29 2.88061 4.20894 8.22601
1033  // 31 3.0618 4.4508 8.61545
1034  // 33 2.94732 4.17463 7.85741
1035  // 35 2.92894 4.20026 8.0066
1036  // 37 2.70991 3.76532 6.85238
1037  // 39 2.91857 4.25505 8.23559
1038  // 41 2.39766 3.14728 5.31672
1039  // 43 2.41189 3.19676 5.40647
1040  // 45 2.0471 2.56789 4.11909
1041  // 47 2.04743 2.57154 4.09681
1042  // 49 1.97205 2.43492 3.7728
1043  // 51 1.67187 1.96796 2.89214
1044  // 53 1.75599 2.0696 2.98919
1045  // 55 1.69526 1.94817 2.71746
1046  // 57 1.02735 1.09074 1.46593
1047  // 59 1.56319 1.76551 2.39532
1048  // 61 1.60436 1.82217 2.45074
1049  // 63 1.39696 1.47696 1.82356
1050  // 65 1.45008 1.62795 2.14264
1051  // 67 1.36533 1.43573 1.73491
1052  // 69 1.40584 1.4939 1.80942
1053  // 71 1.38037 1.45466 1.75179
1054  // 73 1.27568 1.28553 1.49002
1055  // 75 1.03127 1.04286 1.19033
1056  // 77 1.20987 1.25797 1.47642
1057  // 79 1.27784 1.31767 1.51129
1058  // 81 1.06928 1.06112 1.1719
1059  // 83 1.22018 1.2414 1.386
1060  // 85 1.15586 1.13797 1.23096
1061  // 87 1.00112 0.966856 1.02192
1062  // 89 1.07953 1.05446 1.1226
1063  // 91 1.08107 1.06196 1.1313
1064  // 93 1.6204 1.70796 1.87276
1065  // 95 0.905063 0.833268 0.839565
1066  // 97 0.94488 0.884202 0.892716
1067  // 99 0.956238 0.882098 0.880081
1068  // 101 0.855383 0.777236 0.768626
1069  // 103 0.832522 0.754264 0.742384
1070 
1071  const int n_Data = 44;
1072 
1073  const double m12[] = {
1074  17,
1075  19,
1076  21,
1077  23,
1078  25,
1079  27,
1080  29,
1081  31,
1082  33,
1083  35,
1084  37,
1085  39,
1086  41,
1087  43,
1088  45,
1089  47,
1090  49,
1091  51,
1092  53,
1093  55,
1094  57,
1095  59,
1096  61,
1097  63,
1098  65,
1099  67,
1100  69,
1101  71,
1102  73,
1103  75,
1104  77,
1105  79,
1106  81,
1107  83,
1108  85,
1109  87,
1110  89,
1111  91,
1112  93,
1113  95,
1114  97,
1115  99,
1116  101,
1117  103};
1118 
1119  const double RAA_1mb_g_18[] = {
1120  4.10724,
1121  4.44519,
1122  4.57977,
1123  4.29797,
1124  3.58437,
1125  3.81735,
1126  2.88061,
1127  3.0618,
1128  2.94732,
1129  2.92894,
1130  2.70991,
1131  2.91857,
1132  2.39766,
1133  2.41189,
1134  2.0471,
1135  2.04743,
1136  1.97205,
1137  1.67187,
1138  1.75599,
1139  1.69526,
1140  1.02735,
1141  1.56319,
1142  1.60436,
1143  1.39696,
1144  1.45008,
1145  1.36533,
1146  1.40584,
1147  1.38037,
1148  1.27568,
1149  1.03127,
1150  1.20987,
1151  1.27784,
1152  1.06928,
1153  1.22018,
1154  1.15586,
1155  1.00112,
1156  1.07953,
1157  1.08107,
1158  1.6204,
1159  0.905063,
1160  0.94488,
1161  0.956238,
1162  0.855383,
1163  0.832522};
1164 
1165  const double RAA_1mb_g_20[] = {
1166  7.40351,
1167  7.83383,
1168  7.8564,
1169  7.07926,
1170  5.58706,
1171  5.93123,
1172  4.20894,
1173  4.4508,
1174  4.17463,
1175  4.20026,
1176  3.76532,
1177  4.25505,
1178  3.14728,
1179  3.19676,
1180  2.56789,
1181  2.57154,
1182  2.43492,
1183  1.96796,
1184  2.0696,
1185  1.94817,
1186  1.09074,
1187  1.76551,
1188  1.82217,
1189  1.47696,
1190  1.62795,
1191  1.43573,
1192  1.4939,
1193  1.45466,
1194  1.28553,
1195  1.04286,
1196  1.25797,
1197  1.31767,
1198  1.06112,
1199  1.2414,
1200  1.13797,
1201  .966856,
1202  1.05446,
1203  1.06196,
1204  1.70796,
1205  .833268,
1206  .884202,
1207  .882098,
1208  0.777236,
1209  0.754264};
1210 
1211  const double RAA_1mb_g_22[] = {
1212  16.6273,
1213  18.0237,
1214  17.5948,
1215  15.2511,
1216  11.445,
1217  12.1819,
1218  8.22601,
1219  8.61545,
1220  7.85741,
1221  8.0066,
1222  6.85238,
1223  8.23559,
1224  5.31672,
1225  5.40647,
1226  4.11909,
1227  4.09681,
1228  3.7728,
1229  2.89214,
1230  2.98919,
1231  2.71746,
1232  1.46593,
1233  2.39532,
1234  2.45074,
1235  1.82356,
1236  2.14264,
1237  1.73491,
1238  1.80942,
1239  1.75179,
1240  1.49002,
1241  1.19033,
1242  1.47642,
1243  1.51129,
1244  1.1719,
1245  1.386,
1246  1.23096,
1247  1.02192,
1248  1.1226,
1249  1.1313,
1250  1.87276,
1251  0.839565,
1252  0.892716,
1253  0.880081,
1254  0.768626,
1255  0.742384};
1256 
1257  TGraph *gRAARatioKang2019(NULL);
1258 
1259  if (mass == 1)
1260  {
1261  if (g == 1.8)
1262  {
1263  gRAARatioKang2019 = new TGraph(n_Data, m12, RAA_1mb_g_18);
1264  }
1265  else if (g == 2)
1266  {
1267  gRAARatioKang2019 = new TGraph(n_Data, m12, RAA_1mb_g_20);
1268  }
1269  else if (g == 2.2)
1270  {
1271  gRAARatioKang2019 = new TGraph(n_Data, m12, RAA_1mb_g_22);
1272  }
1273  }
1274 
1275  assert(gRAARatioKang2019);
1276 
1277  return gRAARatioKang2019;
1278 }