47 const std::string& bgcurvefile,
48 const std::string& bgsigmafile,
const std::string& pdg,
bool ismakePlots)
51 std::vector<TH1F*> hdedx_bg, hchi_bg, hionzsigma_bg;
52 std::map<int, std::vector<TH1F*>> hchi_inj, hchicos_allbg, hchicos_1by3bg, hchicos_2by3bg, hchicos_3by3bg;
57 defineHisto(hionzsigma_bg,
"ionzsigma",
"bg", pdg.data());
58 for (
int i = 0; i < 2; ++i) {
60 defineHisto(hchi_inj[i],
"chi", Form(
"inj_%d", i), pdg.data());
62 std::string charge =
"pos";
63 if (i == 1) charge =
"neg";
65 defineHisto(hchicos_allbg[i],
"chi", Form(
"%s_allbg_cos", charge.data()), pdg.data());
66 defineHisto(hchicos_1by3bg[i],
"chi", Form(
"%s_1b3bg_cos", charge.data()), pdg.data());
67 defineHisto(hchicos_2by3bg[i],
"chi", Form(
"%s_2b3bg_cos", charge.data()), pdg.data());
68 defineHisto(hchicos_3by3bg[i],
"chi", Form(
"%s_3b3bg_cos", charge.data()), pdg.data());
101 hadron->SetBranchAddress(
"dedxnosat", &dedxnosat);
102 hadron->SetBranchAddress(
"p", &p);
103 hadron->SetBranchAddress(
"costh", &costh);
104 hadron->SetBranchAddress(
"timereso", &timereso);
105 hadron->SetBranchAddress(
"nhits", &nhit);
106 hadron->SetBranchAddress(
"injtime", &injtime);
107 hadron->SetBranchAddress(
"injring", &isher);
108 hadron->SetBranchAddress(
"charge", &charge);
111 if (mass == 0.0) B2FATAL(
"Mass of particle " << pdg.data() <<
" is zero");
124 unsigned int entries = hadron->GetEntries();
126 for (
unsigned int index = 0; index < entries; ++index) {
128 hadron->GetEvent(index);
130 int chg = (charge < 0) ? 1 : 0;
131 double bg = fabs(p) / mass;
134 if (nhit < 0 || nhit > 100)
continue;
135 if (injtime <= 0 || isher < 0)
continue;
137 if (dedxnosat <= 0)
continue;
138 if (costh != costh)
continue;
139 if (bg < m_bgMin || bg >
m_bgMax)
continue;
141 if (fabs(p) > 7.0)
continue;
143 if (pdg ==
"proton")
if ((dedxnosat - 0.45) * abs(p) * abs(p) <
m_cut)
continue;
147 double dedx_new = had.
D2I(costh, had.
I2D(costh, 1.0) * dedxnosat);
149 hdedx_bg[bgBin]->Fill(dedx_new);
152 double dedx_cur = mbg.
getMean(bg);
153 double dedx_res = sbg.
getSigma(dedx_new, nhit, costh, timereso);
156 B2INFO(
"RESOLUTION IS ZERO!!!");
161 double chi_new = (dedx_new - dedx_cur) / dedx_res;
163 hchi_bg[bgBin]->Fill(chi_new);
168 B2INFO(
"RESOLUTION costh * nhit * timereso IS ZERO!!!");
172 hionzsigma_bg[bgBin]->Fill((dedx_new - dedx_cur) / ionz_res);
180 int icos = (int)((costh + 1) / cosstep);
181 hchicos_allbg[chg][icos]->Fill(chi_new);
183 if (bgBin <=
int(
m_bgBins / 3)) hchicos_1by3bg[chg][icos]->Fill(chi_new);
184 else if (bgBin <=
int(2 *
m_bgBins / 3)) hchicos_2by3bg[chg][icos]->Fill(chi_new);
185 else hchicos_3by3bg[chg][icos]->Fill(chi_new);
189 if (isher > 0.5) wr = 1;
190 int injBin = (int)((injtime -
m_injMin) / tstep);
191 hchi_inj[wr][injBin]->Fill(chi_new);
204 setPars(outfile, pdg, hdedx_bg, hchi_bg, hionzsigma_bg, hchi_inj);
209 plotDist(hdedx_bg, Form(
"fits_dedx_inbg_%s_%s", suffix.data(), pdg.data()),
m_bgBins);
210 plotDist(hchi_bg, Form(
"fits_chi_inbg_%s_%s", suffix.data(), pdg.data()),
m_bgBins);
211 plotDist(hionzsigma_bg, Form(
"fits_ionzreso_inbg_%s_%s", suffix.data(), pdg.data()),
m_bgBins);
212 plotDist(hchi_inj, Form(
"fits_chi_ininj_%s_%s", suffix.data(), pdg.data()),
m_injBins);
213 printCanvasCos(hchicos_allbg, hchicos_1by3bg, hchicos_2by3bg, hchicos_3by3bg, pdg, suffix);
221 for (
int i = 0; i < 2; ++i) {
234 int nbdEdx = 200, nbchi = 200;
235 double dedxMax = 20.0;
238 nbdEdx = 400, dedxMax = 4.0;
239 }
else if (pdg ==
"kaon") {
240 nbdEdx = 500, dedxMax = 5.0;
241 }
else if (pdg ==
"proton") {
242 nbdEdx = 1500, dedxMax = 30.0;
244 }
else if (pdg ==
"muon") {
245 nbdEdx = 300, dedxMax = 3.0;
246 }
else if (pdg ==
"electron") {
247 nbdEdx = 200, dedxMax = 2.0;
258 double step = (max - min) / bins;
259 if (stype ==
"nhit") step = (max - min + 1) / bins;
261 for (
int j = 0; j < bins; ++j) {
263 double start = min + j * step;
264 double end = start + step;
265 if (stype ==
"nhit") end = int((start + step) * 0.99999);
268 std::string histname = Form(
"%s_%s_%s_%d", pdg.data(), svar.data(), stype.data(), j);
269 std::string title = Form(
"%s_%s_%s (%.02f, %.02f)", pdg.data(), svar.data(), stype.data(), start, end);
272 htemp.push_back(
new TH1F(histname.data(), title.data(), nbdEdx, 0, dedxMax));
273 else if (svar ==
"chi")
274 htemp.push_back(
new TH1F(histname.data(), title.data(), nbchi, -10.0, 10.0));
276 htemp.push_back(
new TH1F(histname.data(), title.data(), 300, -3.0, 3.0));
278 htemp[j]->GetXaxis()->SetTitle(Form(
"%s", svar.data()));
279 htemp[j]->GetYaxis()->SetTitle(
"Entries");
517 std::map<
int, std::vector<TH1F*>>& hchicos_1by3bg, std::map<
int, std::vector<TH1F*>>& hchicos_2by3bg,
518 std::map<
int, std::vector<TH1F*>>& hchicos_3by3bg,
const std::string& pdg,
const std::string& suffix)
521 std::vector<std::vector<double>> chicos(2, std::vector<double>(
m_cosBins));
522 std::vector<std::vector<double>> sigmacos(2, std::vector<double>(
m_cosBins));
523 std::vector<std::vector<double>> chicoserr(2, std::vector<double>(
m_cosBins));
524 std::vector<std::vector<double>> sigmacoserr(2, std::vector<double>(
m_cosBins));
526 std::vector<std::vector<double>> chicos_1b3bg(2, std::vector<double>(
m_cosBins));
527 std::vector<std::vector<double>> sigmacos_1b3bg(2, std::vector<double>(
m_cosBins));
528 std::vector<std::vector<double>> chicos_1b3bgerr(2, std::vector<double>(
m_cosBins));
529 std::vector<std::vector<double>> sigmacos_1b3bgerr(2, std::vector<double>(
m_cosBins));
531 std::vector<std::vector<double>> chicos_2b3bg(2, std::vector<double>(
m_cosBins));
532 std::vector<std::vector<double>> sigmacos_2b3bg(2, std::vector<double>(
m_cosBins));
533 std::vector<std::vector<double>> chicos_2b3bgerr(2, std::vector<double>(
m_cosBins));
534 std::vector<std::vector<double>> sigmacos_2b3bgerr(2, std::vector<double>(
m_cosBins));
536 std::vector<std::vector<double>> chicos2(2, std::vector<double>(
m_cosBins));
537 std::vector<std::vector<double>> sigmacos2(2, std::vector<double>(
m_cosBins));
538 std::vector<std::vector<double>> chicos2err(2, std::vector<double>(
m_cosBins));
539 std::vector<std::vector<double>> sigmacos2err(2, std::vector<double>(
m_cosBins));
541 for (
int c = 0; c < 2; ++c) {
543 if (hchicos_allbg[c][i]->Integral() > 100) {
544 fit(hchicos_allbg[c][i], pdg.data());
545 chicos[c][i] = hchicos_allbg[c][i]->GetFunction(
"gaus")->GetParameter(1);
546 chicoserr[c][i] = hchicos_allbg[c][i]->GetFunction(
"gaus")->GetParError(1);
547 sigmacos[c][i] = hchicos_allbg[c][i]->GetFunction(
"gaus")->GetParameter(2);
548 sigmacoserr[c][i] = hchicos_allbg[c][i]->GetFunction(
"gaus")->GetParError(2);
551 if (hchicos_1by3bg[c][i]->Integral() > 100) {
552 fit(hchicos_1by3bg[c][i], pdg.data());
553 chicos_1b3bg[c][i] = hchicos_1by3bg[c][i]->GetFunction(
"gaus")->GetParameter(1);
554 chicos_1b3bgerr[c][i] = hchicos_1by3bg[c][i]->GetFunction(
"gaus")->GetParError(1);
555 sigmacos_1b3bg[c][i] = hchicos_1by3bg[c][i]->GetFunction(
"gaus")->GetParameter(2);
556 sigmacos_1b3bgerr[c][i] = hchicos_1by3bg[c][i]->GetFunction(
"gaus")->GetParError(2);
559 if (hchicos_2by3bg[c][i]->Integral() > 100) {
560 fit(hchicos_2by3bg[c][i], pdg.data());
561 chicos_2b3bg[c][i] = hchicos_2by3bg[c][i]->GetFunction(
"gaus")->GetParameter(1);
562 chicos_2b3bgerr[c][i] = hchicos_2by3bg[c][i]->GetFunction(
"gaus")->GetParError(1);
563 sigmacos_2b3bg[c][i] = hchicos_2by3bg[c][i]->GetFunction(
"gaus")->GetParameter(2);
564 sigmacos_2b3bgerr[c][i] = hchicos_2by3bg[c][i]->GetFunction(
"gaus")->GetParError(2);
567 if (hchicos_3by3bg[c][i]->Integral() > 100) {
568 fit(hchicos_3by3bg[c][i], pdg.data());
569 chicos2[c][i] = hchicos_3by3bg[c][i]->GetFunction(
"gaus")->GetParameter(1);
570 chicos2err[c][i] = hchicos_3by3bg[c][i]->GetFunction(
"gaus")->GetParError(1);
571 sigmacos2[c][i] = hchicos_3by3bg[c][i]->GetFunction(
"gaus")->GetParameter(2);
572 sigmacos2err[c][i] = hchicos_3by3bg[c][i]->GetFunction(
"gaus")->GetParError(2);
577 plotDist(hchicos_allbg, Form(
"fits_chi_vscos_allbg_%s_%s", pdg.data(), suffix.data()),
m_cosBins);
578 plotDist(hchicos_1by3bg, Form(
"fits_chi_vscos_1by3bg_%s_%s", pdg.data(), suffix.data()),
m_cosBins);
579 plotDist(hchicos_2by3bg, Form(
"fits_chi_vscos_2by3bg_%s_%s", pdg.data(), suffix.data()),
m_cosBins);
580 plotDist(hchicos_3by3bg, Form(
"fits_chi_vscos_3by3bg_%s_%s", pdg.data(), suffix.data()),
m_cosBins);
585 cosArray[i] = -1.0 + (i * cosstep + cosstep / 2.0);
586 cosArrayErr[i] = 0.0;
589 TGraphErrors grchicos(
m_cosBins, cosArray.data(), chicos[0].data(), cosArrayErr.data(), chicoserr[0].data());
590 TGraphErrors grchicos_1b3bg(
m_cosBins, cosArray.data(), chicos_1b3bg[0].data(), cosArrayErr.data(), chicos_1b3bgerr[0].data());
591 TGraphErrors grchicos_2b3bg(
m_cosBins, cosArray.data(), chicos_2b3bg[0].data(), cosArrayErr.data(), chicos_2b3bgerr[0].data());
592 TGraphErrors grchicos2(
m_cosBins, cosArray.data(), chicos2[0].data(), cosArrayErr.data(), chicos2err[0].data());
594 TGraphErrors grchicosn(
m_cosBins, cosArray.data(), chicos[1].data(), cosArrayErr.data(), chicoserr[1].data());
595 TGraphErrors grchicos_1b3bgn(
m_cosBins, cosArray.data(), chicos_1b3bg[1].data(), cosArrayErr.data(), chicos_1b3bgerr[1].data());
596 TGraphErrors grchicos_2b3bgn(
m_cosBins, cosArray.data(), chicos_2b3bg[1].data(), cosArrayErr.data(), chicos_2b3bgerr[1].data());
597 TGraphErrors grchicos2n(
m_cosBins, cosArray.data(), chicos2[1].data(), cosArrayErr.data(), chicos2err[1].data());
599 TGraphErrors grsigmacos(
m_cosBins, cosArray.data(), sigmacos[0].data(), cosArrayErr.data(), sigmacoserr[0].data());
600 TGraphErrors grsigmacos_1b3bg(
m_cosBins, cosArray.data(), sigmacos_1b3bg[0].data(), cosArrayErr.data(),
601 sigmacos_1b3bgerr[0].data());
602 TGraphErrors grsigmacos_2b3bg(
m_cosBins, cosArray.data(), sigmacos_2b3bg[0].data(), cosArrayErr.data(),
603 sigmacos_2b3bgerr[0].data());
604 TGraphErrors grsigmacos2(
m_cosBins, cosArray.data(), sigmacos2[0].data(), cosArrayErr.data(), sigmacos2err[0].data());
606 TGraphErrors grsigmacosn(
m_cosBins, cosArray.data(), sigmacos[1].data(), cosArrayErr.data(), sigmacoserr[1].data());
607 TGraphErrors grsigmacos_1b3bgn(
m_cosBins, cosArray.data(), sigmacos_1b3bg[1].data(), cosArrayErr.data(),
608 sigmacos_1b3bgerr[1].data());
609 TGraphErrors grsigmacos_2b3bgn(
m_cosBins, cosArray.data(), sigmacos_2b3bg[1].data(), cosArrayErr.data(),
610 sigmacos_2b3bgerr[1].data());
611 TGraphErrors grsigmacos2n(
m_cosBins, cosArray.data(), sigmacos2[1].data(), cosArrayErr.data(), sigmacos2err[1].data());
613 TLine line0(-1, 0, 1, 0);
614 line0.SetLineStyle(kDashed);
615 line0.SetLineColor(kRed);
617 TLine line1(-1, 1, 1, 1);
618 line1.SetLineStyle(kDashed);
619 line1.SetLineColor(kRed);
623 if (pdg ==
"pion") {ptype +=
"#pi"; mass =
Const::pion.getMass();}
624 else if (pdg ==
"kaon") {ptype +=
"K"; mass =
Const::kaon.getMass();}
625 else if (pdg ==
"proton") { ptype +=
"p"; mass =
Const::proton.getMass();}
626 else if (pdg ==
"muon") {ptype +=
"#mu"; mass =
Const::muon.getMass();}
627 else if (pdg ==
"electron") {ptype +=
"e"; mass =
Const::electron.getMass();}
628 if (mass == 0.0) B2FATAL(
"Mass of particle " << pdg.data() <<
" is zero");
630 int bglow = 0, bghigh = 0;
633 TLegend lchi(0.7, 0.75, 0.8, 0.85);
634 lchi.SetBorderSize(0);
635 lchi.SetFillColor(kYellow);
636 lchi.AddEntry(&grchicos, ptype +
"^{+}",
"pc");
637 lchi.AddEntry(&grchicosn, ptype +
"^{-}",
"pc");
639 TCanvas cchi(Form(
"cchi_%s", suffix.data()),
"cchi", 700, 600);
645 grchicosn.Draw(
"P,same");
651 FormatGraph(grchicos_1b3bg, 0, Form(
"first 1/3 bg bins, p =(%0.02f, %0.02f)",
m_bgMin * mass, bghigh * bgstep * mass));
653 grchicos_1b3bg.Draw(
"AP");
654 grchicos_1b3bgn.Draw(
"P,same");
661 FormatGraph(grchicos_2b3bg, 0, Form(
"second 1/3 bg bins, p =(%0.02f, %0.02f)", bglow * bgstep * mass, bghigh * bgstep * mass));
663 grchicos_2b3bg.Draw(
"AP");
664 grchicos_2b3bgn.Draw(
"P,same");
670 FormatGraph(grchicos2, 0, Form(
"third 1/3 bg bins, p =(%0.02f, %0.02f)", bglow * bgstep * mass,
m_bgMax * mass));
672 grchicos2.Draw(
"AP");
673 grchicos2n.Draw(
"P,same");
676 cchi.SaveAs(Form(
"plots/HadronCal/Monitoring/mean_chi_vscos_inbg_%s_%s.pdf", pdg.data(), suffix.data()));
681 grsigmacos.Draw(
"AP");
682 grsigmacosn.Draw(
"P,same");
688 FormatGraph(grsigmacos_1b3bg, 2, Form(
"first 1/3 bg bins, p =(%0.02f, %0.02f)",
m_bgMin * mass, bghigh * bgstep * mass));
690 grsigmacos_1b3bg.Draw(
"AP");
691 grsigmacos_1b3bgn.Draw(
"P,same");
698 FormatGraph(grsigmacos_2b3bg, 2, Form(
"second 1/3 bg bins, p =(%0.02f, %0.02f)", bglow * bgstep * mass, bghigh * bgstep * mass));
700 grsigmacos_2b3bg.Draw(
"AP");
701 grsigmacos_2b3bgn.Draw(
"P,same");
707 FormatGraph(grsigmacos2, 2, Form(
"third 1/3 bg bins, p =(%0.02f, %0.02f)", bglow * bgstep * mass,
m_bgMax * mass));
709 grsigmacos2.Draw(
"AP");
710 grsigmacos2n.Draw(
"P,same");
713 cchi.SaveAs(Form(
"plots/HadronCal/Monitoring/sigma_chi_vscos_inbg_%s_%s.pdf", pdg.data(), suffix.data()));