455 string stepname =
fHMean.fAnaStepNames[
static_cast<int>(step)];
456 for (
const string cat : {
"",
"_trueEl"}) {
459 fHMean.CreateHByClone(
fHMean.H1Clone(
"hfsc_MinvBg" + cat, step),
fHMean.GetName(
"hSignal" + cat, step));
460 TH1D* sigred =
fHMean.CreateHByClone(
fHFastSim.H1Clone(
"hfsc_MinvBg" + cat, step),
461 fHMean.GetName(
"hSignalReduced" + cat, step));
462 TH1D* reBg =
fHFastSim.H1Clone(
"hfsc_MinvBg" + cat, step);
465 TH1D* cbc =
fHMean.H1Clone(
"hCbBg" + cat, step);
467 double bW = sig->GetBinWidth(1);
468 int ChangeBin = sig->FindBin(0.4);
469 for (
int iB = 1; iB <= sig->GetNbinsX(); iB++) {
470 double vMcBg = mcBg->GetBinContent(iB);
471 double vReBg = reBg->GetBinContent(iB);
472 double vCbBg = cbc->GetBinContent(iB);
473 double vCoc = coc->GetBinContent(iB);
474 double vSig = (iB < ChangeBin) ? vMcBg + vCoc - vCbBg : vReBg + vCoc - vCbBg;
475 sig->SetBinContent(iB, vSig);
476 sigred->SetBinContent(iB, vSig);
492 sigred->Add(omega, -1.);
493 sigred->Add(omegadalitz, -1.);
494 sigred->Add(phi, -1.);
495 sigred->Add(pi0, -1.);
496 sigred->Add(eta, -1.);
514 LOG(error) <<
"LmvmDrawAll::CalculateSignal: Check bwVar-tag assignment!";
517 TH1D* bre =
VaryBinWidth(
"mean",
"hfsc_MinvBg" + cat, step, bwv);
520 TH1D* errSig = (TH1D*) sigV->Clone();
521 TH1D* errBre = (TH1D*) sigV->Clone();
522 TH1D* errBcb = (TH1D*) sigV->Clone();
528 TH1D* h_s_pp =
VaryBinWidth(
"evmix",
fHMean.GetName(
"hmx_MinvPP" + cat +
"_sameEv", step), bwv);
529 TH1D* h_s_mm =
VaryBinWidth(
"evmix",
fHMean.GetName(
"hmx_MinvMM" + cat +
"_sameEv", step), bwv);
530 TH1D* h_m_pp =
VaryBinWidth(
"evmix",
fHMean.GetName(
"hmx_MinvPP" + cat +
"_mixedEv", step), bwv);
531 TH1D* h_m_mm =
VaryBinWidth(
"evmix",
fHMean.GetName(
"hmx_MinvMM" + cat +
"_mixedEv", step), bwv);
533 for (
int iB = 1; iB <= bre->GetNbinsX(); iB++) {
534 double bW2 = bre->GetBinWidth(iB);
547 double D_bcb = (bre->GetBinCenter(iB) <=
fCbChange) ? D_bcb1 : D_bcb2;
549 double D_dNdM = std::sqrt(D_bre * D_bre + D_bcb * D_bcb);
551 sigV->SetBinError(iB, D_dNdM);
552 sigredV->SetBinError(iB, D_dNdM);
554 errSig->SetBinContent(iB, D_dNdM);
555 errBre->SetBinContent(iB, D_bre);
556 errBcb->SetBinContent(iB, D_bcb);
558 fHMean.CreateHByClone(sigV,
fHMean.GetName(
"hSignal" + cat +
"_" + tag, step));
559 fHMean.CreateHByClone(sigredV,
fHMean.GetName(
"hSignalReduced" + cat +
"_" + tag, step));
561 if (cat !=
"")
continue;
563 sigV->GetYaxis()->SetRangeUser(1e-9, 2e-1);
564 string cNameSig =
"Signal/Standard/" + tag +
"_" + stepname;
565 string cNameErr =
"Signal/Standard/Error/" + tag +
"_" + stepname;
566 fHMean.fHM.CreateCanvas(cNameSig, cNameSig, 800, 800);
569 fHMean.fHM.CreateCanvas(cNameErr, cNameErr, 800, 800);
570 errSig->GetYaxis()->SetRangeUser(1e-10, 5e-2);
571 DrawH1({errSig, errBre, errBcb}, {
"S = B_{FS} - B_{CB}",
"B_{FS}",
"CB_{calc}"},
kLinear,
kLog,
true, 0.65, 0.7,
729 vector<double> paramsA = {5e-5, 1e-4, 5e-4, 1e-3, 2e-3, 5e-3, 0.01, 0.02, 0.05, 0.1, 0.2, 0.5, 1.};
730 vector<double> paramsB = {0.01, 0.02, 0.05, 0.1, 0.2, 0.5, 1., 2., 5., 10., 20., 50., 100.};
732 string h2descr =
"hTempVsInitParams; Parameter #it{a}; Parameter #it{b}; T [MeV]";
733 TH2D* h2 =
new TH2D(
"hTempVsInitParams", h2descr.c_str(), paramsA.size(), paramsA[0], paramsA[paramsA.size() - 1],
734 paramsB.size(), paramsB[0], paramsB[paramsB.size() - 1]);
735 for (
size_t iX = 0; iX < paramsA.size(); iX++) {
737 h2->GetXaxis()->SetBinLabel(iX + 1, startstr.c_str());
739 for (
size_t iY = 0; iY < paramsB.size(); iY++) {
741 h2->GetYaxis()->SetBinLabel(iY + 1, endstr.c_str());
744 for (
size_t iA = 0; iA < paramsA.size(); iA++) {
745 for (
size_t iB = 0; iB < paramsB.size(); iB++) {
749 TF1* TFitCoc =
new TF1(
"TFitCoc",
"[0]*TMath::Exp(-[1]*x)", 1.1, 2.5);
750 TFitCoc->SetParameters(paramsA[iA], paramsB[iB]);
751 string cNameInit =
"Temperature/CheckoutInitParams/Init_" +
Cbm::NumberToString(paramsA[iA], 2) +
"_"
753 fHMean.fHM.CreateCanvas(cNameInit, cNameInit, 1600, 1100);
755 h->Fit(
"TFitCoc",
"QR");
756 double b = TFitCoc->GetParameter(1);
758 DrawTextOnPad(
"Fit: #it{a#upointe^{-bx}}", 0.6, 0.76, 0.8, 0.81);
761 TF1* funcInit =
new TF1(
"TFuncTestInitParams",
"[0]*TMath::Power(x*[1], 3./2.)*TMath::Exp(-x/[1])", 1.1, 2.5);
762 funcInit->SetParameters(0.02, 0.15);
763 string cName =
"Temperature/CheckoutInitParams/Temp_" +
Cbm::NumberToString(paramsA[iA], 2) +
"_"
765 fHMean.fHM.CreateCanvas(cName, cName, 800, 800);
766 int nBins =
h->GetNbinsX();
767 double x[nBins],
y[nBins], ex[nBins], ey[nBins];
768 for (
int iBin = 1; iBin <= nBins; iBin++) {
769 double bW =
h->GetBinWidth(iBin);
770 double x1 =
h->GetBinLowEdge(iBin);
771 x[iBin - 1] = x1 + (TMath::Log(b * bW) - TMath::Log(1 - TMath::Exp(-b * bW))) / b;
772 y[iBin - 1] =
h->GetBinContent(iBin);
774 ey[iBin - 1] =
h->GetBinError(iBin);
776 auto graph =
new TGraphErrors(nBins,
x,
y, ex, ey);
777 graph->SetTitle(
"Temperature; M_{ee} [GeV/c^{2}]; dN/dM [Gev/c^{2}]^{-1}");
778 graph->GetXaxis()->SetRangeUser(0., 2.5);
779 graph->GetYaxis()->SetRangeUser(1e-10, 1e-1);
780 graph->SetMarkerStyle(21);
783 graph->Fit(
"TFuncTestInitParams",
"QR");
784 double T = 1000. * funcInit->GetParameter(1);
785 double TErr = T * funcInit->GetParError(1) / funcInit->GetParameter(1);
790 h2->SetBinContent(iA + 1, iB + 1, T);
793 h2->GetZaxis()->SetRangeUser(163.5, 164.5);
794 string cNameAll =
"Temperature/CheckoutInitParams/All";
795 fHMean.fHM.CreateCanvas(cNameAll, cNameAll, 1200, 800);
797 DrawTextOnPad(
"Temperature in Dependance on Choice Start Params of Fit", 0.0, 0.9, 0.85, 0.99);
859 vector<double> params1 = {0.0001, 0.001, 0.005, 0.01, 0.02, 0.04, 0.1, 0.15, 0.2, 0.4, 1., 2., 5.};
860 vector<double> params2 = {0.001, 0.01, 0.05, 0.1, 0.15, 0.2, 0.25, 0.4, 0.6, 1., 1.5, 2., 5., 10.};
862 string h2descr =
"hTempVsParams; Parameter 1; Parameter 2; T [MeV]";
863 TH2D* h2 =
new TH2D(
"hTempVsParams", h2descr.c_str(), params1.size(), params1[0], params1[params1.size() - 1],
864 params2.size(), params2[0], params2[params2.size() - 1]);
865 for (
size_t iX = 0; iX < params1.size(); iX++) {
867 h2->GetXaxis()->SetBinLabel(iX + 1, startstr.c_str());
869 for (
size_t iY = 0; iY < params2.size(); iY++) {
871 h2->GetYaxis()->SetBinLabel(iY + 1, endstr.c_str());
874 for (
size_t iP1 = 0; iP1 < params1.size(); iP1++) {
875 for (
size_t iP2 = 0; iP2 < params2.size(); iP2++) {
876 TF1* func =
new TF1(
"TFuncTestParams",
"[0]*TMath::Power(x*[1], 3./2.)*TMath::Exp(-x/[1])", 1.1, 2.5);
877 func->SetParameters(params1[iP1], params2[iP2]);
880 fHMean.fHM.CreateCanvas(cName, cName, 800, 800);
881 int nBins =
h->GetNbinsX();
882 double x[nBins],
y[nBins], ex[nBins], ey[nBins];
883 for (
int iB = 1; iB <= nBins; iB++) {
884 double bW =
h->GetBinWidth(iB);
885 double x1 =
h->GetBinLowEdge(iB);
886 x[iB - 1] = x1 + (TMath::Log(b * bW) - TMath::Log(1 - TMath::Exp(-b * bW))) / b;
887 y[iB - 1] =
h->GetBinContent(iB);
889 ey[iB - 1] =
h->GetBinError(iB);
891 auto graph =
new TGraphErrors(nBins,
x,
y, ex, ey);
892 graph->SetTitle(
"Temperature; M_{ee} [GeV/c^{2}]; dN/dM [Gev/c^{2}]^{-1}");
893 graph->GetXaxis()->SetRangeUser(0., 2.5);
894 graph->GetYaxis()->SetRangeUser(1e-10, 1e-1);
895 graph->SetMarkerStyle(21);
898 graph->Fit(
"TFuncTestParams",
"QR");
899 double T = 1000. * func->GetParameter(1);
900 double TErr = T * func->GetParError(1) / func->GetParameter(1);
904 0.12, 0.12, 0.5, 0.22);
905 h2->SetBinContent(iP1 + 1, iP2 + 1, T);
908 h2->GetZaxis()->SetRangeUser(163.5, 164.5);
909 string cNameAll =
"Temperature/CheckoutParams/All";
910 fHMean.fHM.CreateCanvas(cNameAll, cNameAll, 1200, 800);
912 DrawTextOnPad(
"Temperature in Dependance on Start Params of Fit", 0.0, 0.9, 0.85, 0.99);
1089 string cName =
"Temperature/";
1100 fHMean.fHM.CreateCanvas(cName +
"Functions", cName +
"Functions", 950, 800);
1101 TF1* cocfit =
new TF1(
"cocfit",
"0.000177 * TMath::Exp(-5.25 * x)", 1.1, 2.5);
1102 TF1* tempfit =
new TF1(
"tempfit",
"0.005255 * TMath::Power(0.1641 * x, 3./2.) * TMath::Exp(-x/0.1641)", 1.1, 2.5);
1103 cocfit->SetTitle(
"TempFuncs;M_{ee} [GeV/c^{2}];dN/dM [GeV/c^{2}]^{-1}");
1104 tempfit->SetLineColor(4);
1106 tempfit->Draw(
"same");
1107 gPad->SetLogy(
true);
1108 gPad->RangeAxis(0., 1e-10, 2.5, 1e-2);
1109 TLegend* legend =
new TLegend(0.55, 0.75, 0.88, 0.88);
1110 legend->AddEntry(cocfit,
"a #upoint e^{-bx}",
"l");
1111 legend->AddEntry(tempfit,
"p_{1} #upoint (T #upoint x)^{3/2} #upoint e^{-x/T}",
"l");
1117 vector<pair<double, double>> varbin =
fBwVarSig3;
1121 TCanvas* c =
fHMean.fHM.CreateCanvas(cName +
"SignalVsCocktail", cName +
"SignalVsCocktail", 2100, 800);
1147 fHMean.fHM.CreateCanvas(cName +
"Cocktail/CheckCorrections", cName +
"Cocktail/CheckCorrections", 3800, 1600);
1152 DrawTextOnPad(
"Cocktail (#it{X} uncorrected)", 0.1, 0.9, 0.7, 0.99);
1155 DrawTextOnPad(
"Cocktail (#it{X} uncorrected)", 0.1, 0.9, 0.7, 0.99);
1158 DrawTextOnPad(
"Cocktail (#it{X} uncorrected)", 0.1, 0.9, 0.7, 0.99);
1161 DrawTextOnPad(
"Cocktail (#it{X} corrected)", 0.11, 0.9, 0.69, 0.99);
1164 DrawTextOnPad(
"Cocktail (#it{X} corrected)", 0.11, 0.9, 0.69, 0.99);
1167 DrawTextOnPad(
"Cocktail (#it{X} corrected)", 0.11, 0.9, 0.69, 0.99);
1176 coc->GetXaxis()->SetRangeUser(1e-10, 1e-1);
1177 coc2->GetXaxis()->SetRangeUser(1e-10, 1e-1);
1178 coc5->GetXaxis()->SetRangeUser(1e-10, 1e-1);
1179 coc20->GetXaxis()->SetRangeUser(1e-10, 1e-1);
1180 for (
int iB = 1; iB <= coc->GetNbinsX(); iB++) {
1181 double m = coc->GetBinCenter(iB);
1182 double y0 = coc->GetBinContent(iB);
1183 double y2 = y0 * (1 + (2 - 1) * (m - 0.) / (2.5 - 0.));
1184 double y5 = y0 * (1 + (5 - 1) * (m - 0.) / (2.5 - 0.));
1185 double y20 = y0 * (1 + (20 - 1) * (m - 0.) / (2.5 - 0.));
1186 coc2->SetBinContent(iB, y2);
1187 coc5->SetBinContent(iB, y5);
1188 coc20->SetBinContent(iB, y20);
1190 string cNameFull = cName +
"Cocktail/" +
"TempVsSlope";
1191 TCanvas* c =
fHMean.fHM.CreateCanvas(cNameFull, cNameFull, 2200, 1600);
1213 TH1D* hFull =
fHMean.H1Clone(
"hSignal_bwvSig3_elid");
1214 TH1D* hRed =
fHMean.H1Clone(
"hSignalReduced_bwvSig3_elid");
1215 TH1D* hRed2 =
fHMean.H1Clone(
"hSignalReduced_bwvSig3_elid");
1216 hFull->GetYaxis()->SetRangeUser(1e-10, 5e-2);
1217 hRed->GetYaxis()->SetRangeUser(1e-10, 5e-2);
1218 hRed2->GetYaxis()->SetRangeUser(1e-10, 5e-2);
1219 TCanvas* c =
fHMean.fHM.CreateCanvas(cName +
"Signal/Signal", cName +
"Signal/Signal", 3300, 800);
1234 vector<pair<double, double>> bwv0 = {make_pair(1.1, 0.1), make_pair(1.5, 0.4), make_pair(2.5, 1.0)};
1235 vector<pair<double, double>> bwv1 = {make_pair(1.1, 0.1), make_pair(1.5, 0.2), make_pair(2.5, 1.0)};
1236 vector<pair<double, double>> bwv2 = {make_pair(1.1, 0.1), make_pair(1.5, 0.4), make_pair(2.5, 0.5)};
1237 vector<pair<double, double>> bwv3 = {make_pair(1.1, 0.1), make_pair(1.5, 0.2), make_pair(2.5, 0.5)};
1238 vector<pair<double, double>> bwv4 = {make_pair(1.1, 0.1), make_pair(1.5, 0.2), make_pair(2.0, 0.25),
1239 make_pair(2.5, 0.5)};
1240 vector<pair<double, double>> bwv5 = {make_pair(1.1, 0.1), make_pair(1.5, 0.2), make_pair(2.0, 0.25),
1241 make_pair(2.5, 0.25)};
1242 vector<pair<double, double>> bwv6 = {make_pair(1.1, 0.1), make_pair(1.5, 0.2), make_pair(2.0, 0.25),
1243 make_pair(2.5, 0.1)};
1245 fHMean.fHM.CreateCanvas(
"Temperature/Signal_bwVariations",
"Temperature/Signal_bwVariations", 1900, 800);
1248 for (vector<pair<double, double>> bwv : {bwv1, bwv2, bwv3, bwv4, bwv5, bwv6}) {
1258 TH1D*
h =
fHMean.H1Clone(
"hSignal_bwvSig_elid");
1259 TCanvas* c =
fHMean.fHM.CreateCanvas(cName +
"Signal/BinCorrection", cName +
"Signal/BinCorrection", 2100, 800);
1263 DrawTextOnPad(
"Signal (uncorrected in #it{X})", 0.1, 0.9, 0.64, 0.99);
1266 DrawTextOnPad(
"Signal (corrected in #it{X})", 0.11, 0.9, 0.62, 0.99);
1273 TH1D* phiElid2 = (TH1D*) phiElid->Clone();
1274 phiElid2->GetXaxis()->SetRangeUser(1.0, 1.2);
1275 TCanvas* c =
fHMean.fHM.CreateCanvas(cName +
"Phi", cName +
"Phi", 1600, 800);
1278 DrawH1({phiElid, phiAcc}, {
"El-ID",
"Acc"},
kLinear,
kLog,
true, 0.75, 0.78, 0.9, 0.9,
"hist");
1280 DrawH1({phiElid2}, {
"El-ID zoom"},
kLinear,
kLog,
true, 0.67, 0.82, 0.9, 0.9,
"hist");
1286 string cName =
"ChargeSymmetry/";
1291 vector<pair<double, double>> varbin = {make_pair(2.5, 0.05)};
1292 TCanvas* c =
fHMean.fHM.CreateCanvas(cName +
"master", cName +
"master", 1600, 2400);
1295 for (
const string type : {
"hmx",
"hfsc"}) {
1296 for (
const string ev : {
"sameEv",
"mixedEv"}) {
1297 if (type ==
"hmx" && ev ==
"mixedEv")
1299 else if (type ==
"hmx" && ev ==
"sameEv")
1301 string mode = (type ==
"hmx") ?
"evmix" :
"fastsim";
1302 string hName = type +
"_" +
"Minv";
1306 int nBins = pm->GetNbinsX();
1307 TH1D* hr_ppmm = (TH1D*) pp->Clone();
1308 TH1D* hr_pppm = (TH1D*) pp->Clone();
1309 TH1D* hr_mmpm = (TH1D*) mm->Clone();
1310 hr_ppmm->Divide(mm);
1311 hr_pppm->Divide(pm);
1312 hr_mmpm->Divide(pm);
1313 int startbin = pm->FindBin(0.7);
1314 double Ipm = pm->Integral(startbin, nBins);
1315 double Ipp = pp->Integral(startbin, nBins);
1316 double Imm = mm->Integral(startbin, nBins);
1317 double r_ppmm = Ipp / Imm;
1318 double r_pppm = Ipp / Ipm;
1319 double r_mmpm = Imm / Ipm;
1320 double k = Ipm / (2 * std::sqrt(Ipp * Imm));
1326 hr_ppmm->GetYaxis()->SetTitle(
"Ratio");
1327 hr_ppmm->GetYaxis()->SetRangeUser(0., 1.5);
1328 string texttype = (type ==
"hmx") ?
"Simulation" :
"FastSim";
1329 string textmode = (ev ==
"sameEv") ?
"(same Events)" :
"(mixed Events)";
1330 DrawH1({hr_ppmm, hr_pppm, hr_mmpm}, {
"++/--",
"++/+-",
"--/+-"},
kLinear,
kLinear,
true, 0.2, 0.6, 0.38, 0.85,
1333 DrawTextOnPad(texttype +
" " + textmode, 0.1, 0.9, 0.9, 0.99);
1339 TH1D* hr_cbc_bre = (TH1D*) cbc->Clone();
1340 TH1D* hr_cbre_bre = (TH1D*) cbre->Clone();
1341 hr_cbc_bre->Divide(bre);
1342 hr_cbre_bre->Divide(bre);
1343 int nBins = bre->GetNbinsX();
1344 int startbin = bre->FindBin(0.7);
1345 double I_bre = bre->Integral(startbin, nBins);
1346 double I_cbc = cbc->Integral(startbin, nBins);
1347 double I_cbre = cbre->Integral(startbin, nBins);
1348 double r_cbc_bre = I_cbc / I_bre;
1349 double r_cbre_bre = I_cbre / I_bre;
1353 DrawH1({hr_cbc_bre, hr_cbre_bre}, {
"CB_{Sim} / B_{FS}",
"CB_{FS} / B_{FS}"},
kLinear,
kLinear,
true, 0.6, 0.7, 0.92,
1355 DrawTextOnPad(
"CB_{Sim} / B_{FS} = " + sr_cbc_bre +
", CB_{FS} / B_{FS} = " + sr_cbre_bre, 0.08, 0.12, 0.89, 0.25);
1359 for (
const string cat : {
""}) {
1363 TH1D* rat1 = (TH1D*) cbc->Clone();
1365 TH1D* rat2 = (TH1D*) cbre->Clone();
1367 rat1->GetYaxis()->SetTitle(
"Ratio");
1368 rat2->GetYaxis()->SetTitle(
"Ratio");
1369 rat1->GetYaxis()->SetRangeUser(0.95, 1.05);
1370 rat2->GetYaxis()->SetRangeUser(0.95, 1.05);
1371 string cNameFull = (cat ==
"") ? cName +
"ReCb_all" : cName +
"ReCb_trueEl";
1372 string label = (cat ==
"") ?
"(all)" :
"(electrons)";
1373 TCanvas* c =
fHMean.fHM.CreateCanvas(cNameFull, cNameFull, 1600, 1600);
1379 if (cat ==
"")
DrawH1({rat2}, {
"CB_{FS} / B_{FS}"},
kLinear,
kLinear,
true, 0.7, 0.8, 0.95, 0.9,
"hist");
1382 for (
const string ev : {
"sameEv",
"mixedEv"}) {
1383 vector<pair<double, double>> varbin =
fBwVarBg4;
1384 if (ev ==
"mixedEv") varbin =
fBwVarBg3;
1385 TH1D* pm =
VaryBinWidth(
"evmix",
"hmx_MinvPM" + cat +
"_" + ev +
"_elid", varbin);
1386 TH1D* pp =
VaryBinWidth(
"evmix",
"hmx_MinvPP" + cat +
"_" + ev +
"_elid", varbin);
1387 TH1D* mm =
VaryBinWidth(
"evmix",
"hmx_MinvMM" + cat +
"_" + ev +
"_elid", varbin);
1388 TH1D* ppmm = (TH1D*) pp->Clone();
1389 TH1D* pppm = (TH1D*) pp->Clone();
1390 TH1D* mmpm = (TH1D*) mm->Clone();
1394 ppmm->GetYaxis()->SetTitle(
"Ratio");
1395 ppmm->GetYaxis()->SetRangeUser(0., 1.6);
1397 DrawH1({ppmm, pppm, mmpm}, {
"N--/N++",
"N++/N+-",
"N--/N+-"},
kLinear,
kLinear,
true, 0.77, 0.77, 0.95, 0.95,
1399 string title = (ev ==
"sameEv") ?
"Same Events" :
"Mixed Events";
1406 string hName =
"hfsc_Minv";
1408 TH1D* pmSame =
fHFastSim.H1Clone(hName +
"PM_sameEv", step);
1409 TH1D* ppSame =
fHFastSim.H1Clone(hName +
"PP_sameEv", step);
1410 TH1D* mmSame =
fHFastSim.H1Clone(hName +
"MM_sameEv", step);
1411 TH1D* pmMixed =
fHFastSim.H1Clone(hName +
"PM_mixedEv", step);
1412 TH1D* ppMixed =
fHFastSim.H1Clone(hName +
"PP_mixedEv", step);
1413 TH1D* mmMixed =
fHFastSim.H1Clone(hName +
"MM_mixedEv", step);
1414 TH1D* bre =
fHMean.H1Clone(
"hfsc_MinvBg", step);
1415 TH1D* cbre =
fHMean.H1Clone(
"hReCb");
1416 TH1D* k =
fHMean.H1Clone(
"hReCbK");
1417 TH1D* rat = (TH1D*) cbre->Clone();
1419 rat->GetYaxis()->SetTitle(
"Ratio");
1420 rat->GetYaxis()->SetRangeUser(0.95, 1.05);
1421 string cNameFull = cName +
"FastSim";
1422 TCanvas* c =
fHMean.fHM.CreateCanvas(cNameFull, cNameFull, 1600, 1600);
1425 DrawH1({pmSame, ppSame, mmSame}, {
"e+e-",
"e+e+",
"e-e-"},
kLinear,
kLog,
true, 0.7, 0.7, 0.91, 0.9,
"hist");
1426 DrawTextOnPad(
"Pair Yields (Same Events)", 0.2, 0.9, 0.8, 0.99);
1428 DrawH1({pmMixed, ppMixed, mmMixed}, {
"e+e-",
"e+e+",
"e-e-"},
kLinear,
kLog,
true, 0.7, 0.7, 0.91, 0.9,
"hist");
1429 DrawTextOnPad(
"Pair Yields (Mixed Events)", 0.19, 0.9, 0.81, 0.99);
1751 string cName =
"CheckVaryBinWidth/";
1752 vector<pair<double, double>> binValues = {make_pair(0.8, 0.05), make_pair(1.5, 0.1), make_pair(2.5, 0.5)};
1763 inmed->Scale(1. / inmed->GetBinWidth(1));
1766 bg->GetYaxis()->SetRangeUser(1e-7, 2e-2);
1767 TCanvas* c =
fHMean.fHM.CreateCanvas(cName +
"single", cName +
"single", 1600, 800);
1770 DrawH1({bg, bg2}, {
"CB_{calc}",
"CB_{calc}2"},
kLinear,
kLog,
true, 0.8, 0.8, 0.95, 0.91,
"pe");
1772 DrawH1({inmed, inmed2}, {
"inmed",
"inmed2"},
kLinear,
kLog,
true, 0.8, 0.8, 0.95, 0.91,
"pe");
1777 TCanvas* c1 =
fHMean.fHM.CreateCanvas(cName +
"minv", cName +
"minv", 1600, 800);
1785 TH1D* sbgO =
static_cast<TH1D*
>(bgO->Clone());
1786 sbgO->Add(cocktailO);
1787 TH1D* cbO =
nullptr;
1790 cbO =
fHMean.H1Clone(
"hCbBg", step);
1794 double binWidth = bgO->GetBinWidth(1);
1795 vector<TH1D*> sHists(
fHMean.fNofSignals);
1798 sHist->Scale(1. /
H(signal)->H1(
"hEventNumber")->GetEntries());
1800 sHist->Scale(1. / binWidth);
1801 sHists[
static_cast<int>(signal)] = sHist;
1804 vector<LmvmDrawMinvData> drawData;
1806 drawData.emplace_back(sbgO, kBlack, kBlack, 1, -1,
"Cocktail+B_{MC}");
1807 drawData.emplace_back(bgO, kGray, kBlack, 1, -1,
"B_{MC}");
1810 drawData.emplace_back(pi0O, kGreen - 3, kGreen + 3, 2, -1,
"#pi^{0} #rightarrow #gammae^{+}e^{-}");
1811 drawData.emplace_back(etaO, kRed - 4, kRed + 2, 2, -1,
"#eta #rightarrow #gammae^{+}e^{-}");
1812 drawData.emplace_back(sHists[
static_cast<int>(
ELmvmSignal::OmegaD)], kCyan + 2, kCyan + 4, 2, -1,
1813 "#omega #rightarrow #pi^{0}e^{+}e^{-}");
1814 drawData.emplace_back(gammaO, -1, kYellow, 4, -1,
"#gamma #rightarrow e^{+}e^{-}");
1815 drawData.emplace_back(sHists[
static_cast<int>(
ELmvmSignal::Omega)], kOrange + 7, kOrange + 4, 2, -1,
1816 "#omega #rightarrow e^{+}e^{-}");
1817 drawData.emplace_back(sHists[
static_cast<int>(
ELmvmSignal::Phi)], kAzure + 2, kAzure + 3, 2, -1,
1818 "#phi #rightarrow e^{+}e^{-}");
1819 drawData.emplace_back(sHists[
static_cast<int>(
ELmvmSignal::Qgp)], kOrange - 2, kOrange - 3, 4, 3112,
1821 drawData.emplace_back(sHists[
static_cast<int>(
ELmvmSignal::Inmed)], kMagenta - 3, kMagenta - 2, 4, 3018,
1823 drawData.emplace_back(cocktailO, -1, kRed + 2, 4, -1,
"Cocktail");
1825 drawData.emplace_back(cbO, -1, kCyan + 1, 4, -1,
"CB_{calc}");
1830 TLegend* leg =
new TLegend(0.72, 0.55, 0.97, 0.99);
1831 for (
size_t i = 0; i < drawData.size(); i++) {
1832 const auto& d = drawData[i];
1833 d.fH->GetYaxis()->SetLabelSize(0.04);
1835 d.fH->GetYaxis()->SetRangeUser(5e-10, 5e-2);
1837 d.fH->GetYaxis()->SetRangeUser(5e-10, 1e2);
1838 if (d.fFillColor != -1) d.fH->SetFillColor(d.fFillColor);
1839 if (d.fFillStyle != -1) d.fH->SetFillStyle(d.fFillStyle);
1840 leg->AddEntry(d.fH, d.fLegend.c_str(),
"f");
1841 DrawH1(d.fH,
kLinear,
kLog, (h0 ==
nullptr) ?
"hist" :
"hist,same", d.fLineColor, d.fLineWidth, 0);
1842 if (h0 ==
nullptr) h0 = d.fH;
1844 fHMean.DrawAnaStepOnPad(step);
1846 leg->SetFillColor(kWhite);
1848 gPad->SetLogy(
true);
1852 fHMean.DrawAnaStepOnPad(step);
1857 TH1D* bre =
fHFastSim.H1Clone(
"hfsc_MinvBg_elid");
1858 TH1D* cbc =
fHMean.H1Clone(
"hCbBg_elid");
1860 TH1D* sig = (TH1D*) bre->Clone();
1863 TH1D* sig2 = (TH1D*) sig->Clone();
1864 TH1D* sig5 = (TH1D*) sig->Clone();
1865 TH1D* sig10 = (TH1D*) sig->Clone();
1866 TH1D* sig25 = (TH1D*) sig->Clone();
1871 sig2->Scale(1. / 2.);
1872 sig5->Scale(1. / 5.);
1873 sig10->Scale(1. / 10.);
1874 sig25->Scale(1. / 25.);
1875 fHMean.fHM.CreateCanvas(cName +
"Signal", cName +
"Signal", 800, 800);
1876 DrawH1({coc, sig2, sig5, sig10, sig25}, {
"Cocktail",
"Signal2",
"Signal5",
"Signal10",
"Signal25"},
kLinear,
kLog,
1877 true, 0.7, 0.8, 0.95, 0.91,
"hist");
1882 vector<double> v_pm, v_pp, v_mm, v_pm_el, v_pp_el, v_mm_el;
1883 for (
const string comb : {
"PM",
"PP",
"MM"}) {
1884 for (
int iMinv = 0; iMinv <= 24; iMinv++) {
1886 TH2D*
h =
fHMean.H2Clone(hName);
1887 double I =
h->Integral(1,
h->GetNbinsX(), 1,
h->GetNbinsX());
1888 double ITrue =
h->Integral(1, 4, 1, 4);
1891 v_pm_el.push_back(ITrue);
1893 else if (comb ==
"PP") {
1895 v_pp_el.push_back(ITrue);
1897 else if (comb ==
"MM") {
1899 v_mm_el.push_back(ITrue);
1904 TH1D* hSymm_ppmm =
new TH1D(
"hSymm_ppmm-vb",
"hSymm_ppmm-vb; M_{ee} [GeV/c^{2}]; Ratio", 25, 0., 2.5);
1905 TH1D* hSymm_pppm =
new TH1D(
"hSymm_pppm-vb",
"hSymm_pppm-vb; M_{ee} [GeV/c^{2}]; Ratio", 25, 0., 2.5);
1906 TH1D* hSymm_mmpm =
new TH1D(
"hSymm_mmpm-vb",
"hSymm_mmpm-vb; M_{ee} [GeV/c^{2}]; Ratio", 25, 0., 2.5);
1907 for (
int iB = 1; iB <= hSymm_ppmm->GetNbinsX(); iB++) {
1908 double ppmm = v_pp[iB - 1] / v_mm[iB - 1];
1909 double pppm = v_pp[iB - 1] / v_pm[iB - 1];
1910 double mmpm = v_mm[iB - 1] / v_pm[iB - 1];
1911 hSymm_ppmm->SetBinContent(iB, ppmm);
1912 hSymm_pppm->SetBinContent(iB, pppm);
1913 hSymm_mmpm->SetBinContent(iB, mmpm);
1915 TH1D* hSymm_ppmm5 = (TH1D*) hSymm_ppmm->Clone();
1916 TH1D* hSymm_pppm5 = (TH1D*) hSymm_pppm->Clone();
1917 TH1D* hSymm_mmpm5 = (TH1D*) hSymm_mmpm->Clone();
1918 hSymm_ppmm5->Rebin(5);
1919 hSymm_pppm5->Rebin(5);
1920 hSymm_mmpm5->Rebin(5);
1921 hSymm_ppmm5->Scale(1. / 5.);
1922 hSymm_pppm5->Scale(1. / 5.);
1923 hSymm_mmpm5->Scale(1. / 5.);
1924 vector<pair<double, double>> varbin5 = {make_pair(2.5, 0.5)};
1928 TH1D* hSymm_ppmmVBg3 =
VaryBinWidth(hSymm_ppmm, varbin5);
1929 TH1D* hSymm_pppmVBg3 =
VaryBinWidth(hSymm_pppm, varbin5);
1930 TH1D* hSymm_mmpmVBg3 =
VaryBinWidth(hSymm_mmpm, varbin5);
1939 hSymm_ppmm->GetYaxis()->SetRangeUser(0., ymax);
1940 hSymm_pppm->GetYaxis()->SetRangeUser(0., ymax);
1941 hSymm_mmpm->GetYaxis()->SetRangeUser(0., ymax);
1942 hSymm_ppmm5->GetYaxis()->SetRangeUser(0., ymax);
1943 hSymm_pppm5->GetYaxis()->SetRangeUser(0., ymax);
1944 hSymm_mmpm5->GetYaxis()->SetRangeUser(0., ymax);
1945 hSymm_ppmmVReg->GetYaxis()->SetRangeUser(0., ymax);
1946 hSymm_pppmVReg->GetYaxis()->SetRangeUser(0., ymax);
1947 hSymm_mmpmVReg->GetYaxis()->SetRangeUser(0., ymax);
1948 hSymm_ppmmVBg3->GetYaxis()->SetRangeUser(0., ymax);
1949 hSymm_pppmVBg3->GetYaxis()->SetRangeUser(0., ymax);
1950 hSymm_mmpmVBg3->GetYaxis()->SetRangeUser(0., ymax);
1951 hSymm_ppmmVBg4->GetYaxis()->SetRangeUser(0., ymax);
1952 hSymm_pppmVBg4->GetYaxis()->SetRangeUser(0., ymax);
1953 hSymm_mmpmVBg4->GetYaxis()->SetRangeUser(0., ymax);
1954 hSymm_ppmmVBg5->GetYaxis()->SetRangeUser(0., ymax);
1955 hSymm_pppmVBg5->GetYaxis()->SetRangeUser(0., ymax);
1956 hSymm_mmpmVBg5->GetYaxis()->SetRangeUser(0., ymax);
1958 string cNameSymm = cName +
"ChargeSymmetries";
1963 TCanvas* c =
fHMean.fHM.CreateCanvas(cNameSymm, cNameSymm, 2400, 1600);
1966 DrawH1({hSymm_ppmm, hSymm_pppm, hSymm_mmpm}, {
"++/--",
"++/+-",
"--/+-"},
kLinear,
kLinear,
true, x1, y1, x2, y2,
1970 DrawH1({hSymm_ppmm5, hSymm_pppm5, hSymm_mmpm5}, {
"++/--",
"++/+-",
"--/+-"},
kLinear,
kLinear,
true, x1, y1, x2, y2,
1974 DrawH1({hSymm_ppmmVReg, hSymm_pppmVReg, hSymm_mmpmVReg}, {
"++/--",
"++/+-",
"--/+-"},
kLinear,
kLinear,
true, x1,
1975 y1, x2, y2,
"hist");
1978 DrawH1({hSymm_ppmmVBg3, hSymm_pppmVBg3, hSymm_mmpmVBg3}, {
"++/--",
"++/+-",
"--/+-"},
kLinear,
kLinear,
true, x1,
1979 y1, x2, y2,
"hist");
1982 DrawH1({hSymm_ppmmVBg4, hSymm_pppmVBg4, hSymm_mmpmVBg4}, {
"++/--",
"++/+-",
"--/+-"},
kLinear,
kLinear,
true, x1,
1983 y1, x2, y2,
"hist");
1986 DrawH1({hSymm_ppmmVBg5, hSymm_pppmVBg5, hSymm_mmpmVBg5}, {
"++/--",
"++/+-",
"--/+-"},
kLinear,
kLinear,
true, x1,
1987 y1, x2, y2,
"hist");
2148 int nBinsMult =
fHFastSim.H1(
"hfsp_mult_urqmd_gammaEl_plus_elid")->GetNbinsX();
2149 vector<string> xLabel;
2150 for (
int iB = 0; iB < nBinsMult; iB++) {
2153 TH1D* hRatP = (TH1D*) plus->Clone();
2154 TH1D* hRatM = (TH1D*) minus->Clone();
2155 hRatP->Divide(plusrandom);
2156 hRatM->Divide(minusrandom);
2157 plus->GetYaxis()->SetRangeUser(1e-8, 1);
2158 plusrandom->GetYaxis()->SetRangeUser(1e-7, 1);
2160 hRatP->GetYaxis()->SetTitle(
"Ratio true/random");
2161 for (
size_t iL = 0; iL < xLabel.size(); iL++) {
2162 plus->GetXaxis()->SetBinLabel(iL + 1, xLabel[iL].c_str());
2163 hRatP->GetXaxis()->SetBinLabel(iL + 1, xLabel[iL].c_str());
2165 double IP = plus->Integral(1, nBinsMult);
2166 double IM = minus->Integral(1, nBinsMult);
2167 double IPrndm = plusrandom->Integral(1, nBinsMult);
2168 double IMrndm = minusrandom->Integral(1, nBinsMult);
2170 minusrandom->Scale(0.1);
2171 TCanvas* c =
fHMean.fHM.CreateCanvas(cName, cName, 1700, 800);
2174 DrawH1({plus, plusrandom, minus, minusrandom},
2175 {
"e^{+}_{true}",
"e^{+}_{random}",
"e^{-}_{true} #upoint 0.1",
"e^{-}_{random} #upoint 0.1"},
kLinear,
kLog,
2176 true, 0.63, 0.62, 0.88, 0.87,
"hist");
2183 DrawH1({hRatP, hRatM}, {
"Ratio e+",
"Ratio e-"},
kLinear,
kLinear,
true, 0.67, 0.74, 0.88, 0.87,
"hist");
2270 string cName =
"FastSim/";
2273 string nEvText =
"NofFastSimEvents = " + nEvStr + (
fNofFastSimEvents >= 1e9 ?
"#times10^{9}" :
"#times10^{6}");
2275 int nBinsMult =
fHFastSim.H1(
"hfsc_mult_urqmd_urqmd_plus_elid")->GetNbinsX();
2276 vector<string> xLabel;
2277 for (
int iB = 0; iB < nBinsMult; iB++) {
2283 string stepname =
fHMean.fAnaStepNames[
static_cast<int>(step)];
2284 string stepLatex =
fHMean.fAnaStepLatex[
static_cast<int>(step)];
2288 TH1D* hP =
fHFastSim.H1Clone(
"hfsc_mult_urqmd_urqmd_plus", step);
2289 TH1D* hM =
fHFastSim.H1Clone(
"hfsc_mult_urqmd_urqmd_minus", step);
2290 TH1D* hPrndm =
fHFastSim.H1Clone(
"hfsc_multRndm_urqmd_urqmd_plus", step);
2291 TH1D* hMrndm =
fHFastSim.H1Clone(
"hfsc_multRndm_urqmd_urqmd_minus", step);
2292 TH1D* hRatP = (TH1D*) hP->Clone();
2293 TH1D* hRatM = (TH1D*) hM->Clone();
2294 hRatP->Divide(hPrndm);
2295 hRatM->Divide(hMrndm);
2296 hP->GetYaxis()->SetRangeUser(1e-8, 1);
2297 hP->GetYaxis()->SetTitle(
"Yield/Event");
2298 hPrndm->GetYaxis()->SetRangeUser(1e-8, 1);
2299 hRatP->GetYaxis()->SetRangeUser(0.5, 1.3);
2300 hRatP->GetYaxis()->SetTitle(
"Ratio true/random");
2301 for (
size_t iL = 0; iL < xLabel.size(); iL++) {
2302 hP->GetXaxis()->SetBinLabel(iL + 1, xLabel[iL].c_str());
2303 hPrndm->GetXaxis()->SetBinLabel(iL + 1, xLabel[iL].c_str());
2304 hRatP->GetXaxis()->SetBinLabel(iL + 1, xLabel[iL].c_str());
2308 double IPrndm = hPrndm->Integral(1, nBinsMult);
2309 double IMrndm = hMrndm->Integral(1, nBinsMult);
2312 for (
int iB = 2; iB <= hP->GetNbinsX(); iB++) {
2316 if (step ==
ELmvmAnaStep::ElId) LOG(info) <<
"DrawFastSim: nP = " << nP <<
", nM = " << nM;
2317 string cNameFull = cName +
"Multiplicities/urqmd_" + stepname;
2318 TCanvas* c =
fHMean.fHM.CreateCanvas(cNameFull, cNameFull, 3200, 800);
2321 DrawH1({hP, hM}, {
"e^{+}",
"e^{-}"},
kLinear,
kLog,
true, 0.8, 0.68, 0.91, 0.9,
"hist");
2327 DrawH1({hPrndm, hMrndm}, {
"e^{+}",
"e^{-}"},
kLinear,
kLog,
true, 0.8, 0.68, 0.91, 0.9,
"hist");
2333 DrawH1({hRatP, hRatM}, {
"e^{+}",
"e^{-}"},
kLinear,
kLinear,
true, 0.8, 0.68, 0.91, 0.9,
"hist");
2334 DrawTextOnPad(
"True vs. Random Multiplicity", 0.1, 0.9, 0.9, 0.99);
2341 string sigName =
fHMean.fSignalNames[
static_cast<int>(signal)];
2342 TH1D* hP =
fHFastSim.H1Clone(
"hfsc_multCorr_pluto_" + sigName +
"_plus", step);
2343 TH1D* hM =
fHFastSim.H1Clone(
"hfsc_multCorr_pluto_" + sigName +
"_minus", step);
2344 TH1D* hPrndm =
fHFastSim.H1Clone(
"hfsc_multRndm_pluto_" + sigName +
"_plus", step);
2345 TH1D* hMrndm =
fHFastSim.H1Clone(
"hfsc_multRndm_pluto_" + sigName +
"_minus", step);
2348 TH1D* hRatP = (TH1D*) hP->Clone();
2349 TH1D* hRatM = (TH1D*) hM->Clone();
2350 hRatP->Divide(hPrndm);
2351 hRatM->Divide(hMrndm);
2352 hP->GetYaxis()->SetRangeUser(1e-6, 1);
2353 hPrndm->GetYaxis()->SetRangeUser(1e-6, 1);
2354 hRatP->GetYaxis()->SetTitle(
"Ratio true/random");
2355 for (
size_t iL = 0; iL < xLabel.size(); iL++) {
2356 hP->GetXaxis()->SetBinLabel(iL + 1, xLabel[iL].c_str());
2357 hPrndm->GetXaxis()->SetBinLabel(iL + 1, xLabel[iL].c_str());
2358 hRatP->GetXaxis()->SetBinLabel(iL + 1, xLabel[iL].c_str());
2360 string cNameFull = cName +
"Multiplicities/" + sigName +
"_" + stepname;
2361 TCanvas* c =
fHMean.fHM.CreateCanvas(cNameFull, cNameFull, 3500, 800);
2364 DrawH1({hP, hM}, {
"e^{+}",
"e^{-}"},
kLinear,
kLog,
true, 0.8, 0.68, 0.91, 0.9,
"hist");
2368 DrawH1({hPrndm, hMrndm}, {
"e^{+}",
"e^{-}"},
kLinear,
kLog,
true, 0.8, 0.68, 0.91, 0.9,
"hist");
2372 DrawH1({hRatP, hRatM}, {
"e^{+}",
"e^{-}"},
kLinear,
kLinear,
true, 0.8, 0.68, 0.91, 0.9,
"hist");
2373 DrawTextOnPad(
"True vs. Random Multiplicity", 0.1, 0.9, 0.9, 0.99);
2380 string cNameFull = cName +
"Multiplicities/all_" + stepname;
2381 TCanvas* c =
fHMean.fHM.CreateCanvas(cNameFull, cNameFull, 3500, 1600);
2384 TH1D* hPU =
fHFastSim.H1Clone(
"hfsc_mult_urqmd_urqmd_plus", step);
2385 TH1D* hMU =
fHFastSim.H1Clone(
"hfsc_mult_urqmd_urqmd_minus", step);
2386 for (
size_t iL = 0; iL < xLabel.size(); iL++) {
2387 hPU->GetXaxis()->SetBinLabel(iL + 1, xLabel[iL].c_str());
2390 DrawH1({hPU, hMU}, {
"e^{+}",
"e^{-}"},
kLinear,
kLog,
true, 0.8, 0.68, 0.91, 0.9,
"hist");
2391 DrawTextOnPad(
"True Multiplicity (UrQMD)", 0.17, 0.9, 0.83, 0.99);
2395 string sigName =
fHMean.fSignalNames[
static_cast<int>(signal)];
2396 TH1D* hP =
fHFastSim.H1Clone(
"hfsc_multCorr_pluto_" + sigName +
"_plus", step);
2397 TH1D* hM =
fHFastSim.H1Clone(
"hfsc_multCorr_pluto_" + sigName +
"_minus", step);
2399 for (
size_t iL = 0; iL < xLabel.size(); iL++) {
2400 hP->GetXaxis()->SetBinLabel(iL + 1, xLabel[iL].c_str());
2403 DrawH1({hP, hM}, {
"e^{+}",
"e^{-}"},
kLinear,
kLog,
true, 0.8, 0.68, 0.91, 0.9,
"hist");
2404 DrawTextOnPad(
"True Multiplicity (" + sigName +
")", 0.17, 0.9, 0.83, 0.99);
2411 static ProcInfo_t info;
2412 const float toMB = 1.f / 1024.f;
2413 gSystem->GetProcInfo(&info);
2414 double sizeRes1 = info.fMemResident * toMB;
2415 double sizeVir1 = info.fMemVirtual * toMB;
2417 THnSparseD* hPlusT = (THnSparseD*)
fHFastSim.fHM.GetObject(
"hfsc_mom_urqmd_urqmd_plus_" + stepname)->Clone();
2418 THnSparseD* hMinusT = (THnSparseD*)
fHFastSim.fHM.GetObject(
"hfsc_mom_urqmd_urqmd_minus_" + stepname)->Clone();
2422 THnSparseD* hPlusR = (THnSparseD*)
fHFastSim.fHM.GetObject(
"hfsc_momRndm_urqmd_urqmd_plus_" + stepname)->Clone();
2423 THnSparseD* hMinusR =
2424 (THnSparseD*)
fHFastSim.fHM.GetObject(
"hfsc_momRndm_urqmd_urqmd_minus_" + stepname)->Clone();
2428 gSystem->GetProcInfo(&info);
2429 double sizeRes2 = info.fMemResident * toMB;
2430 double sizeVir2 = info.fMemVirtual * toMB;
2431 double sizeRes = (sizeRes2 - sizeRes1) / 4.;
2432 double sizeVir = (sizeVir2 - sizeVir1) / 4.;
2433 LOG(info) <<
"Size of one THnSparseD histogram: virtual: " << sizeVir <<
" MB, residential: " << sizeRes
2436 TH1D* hPlusXt = (TH1D*) hPlusT->Projection(0,
"");
2437 TH1D* hPlusYt = (TH1D*) hPlusT->Projection(1,
"");
2438 TH1D* hPlusZt = (TH1D*) hPlusT->Projection(2,
"");
2439 TH1D* hPlusXr = (TH1D*) hPlusR->Projection(0,
"");
2440 TH1D* hPlusYr = (TH1D*) hPlusR->Projection(1,
"");
2441 TH1D* hPlusZr = (TH1D*) hPlusR->Projection(2,
"");
2442 TH1D* hMinusXt = (TH1D*) hMinusT->Projection(0,
"");
2443 TH1D* hMinusYt = (TH1D*) hMinusT->Projection(1,
"");
2444 TH1D* hMinusZt = (TH1D*) hMinusT->Projection(2,
"");
2445 TH1D* hMinusXr = (TH1D*) hMinusR->Projection(0,
"");
2446 TH1D* hMinusYr = (TH1D*) hMinusR->Projection(1,
"");
2447 TH1D* hMinusZr = (TH1D*) hMinusR->Projection(2,
"");
2454 hPlusXt->GetYaxis()->SetTitle(
"Yield / Bin");
2455 hPlusYt->GetYaxis()->SetTitle(
"Yield / Bin");
2456 hPlusZt->GetYaxis()->SetTitle(
"Yield / Bin");
2457 hMinusXt->GetYaxis()->SetTitle(
"Yield / Bin");
2458 hMinusYt->GetYaxis()->SetTitle(
"Yield / Bin");
2459 hMinusZt->GetYaxis()->SetTitle(
"Yield / Bin");
2462 hPlusXt->GetYaxis()->SetRangeUser(
min,
max);
2463 hPlusYt->GetYaxis()->SetRangeUser(
min,
max);
2464 hPlusZt->GetYaxis()->SetRangeUser(
min,
max);
2465 hMinusXt->GetYaxis()->SetRangeUser(
min,
max);
2466 hMinusYt->GetYaxis()->SetRangeUser(
min,
max);
2467 hMinusZt->GetYaxis()->SetRangeUser(
min,
max);
2469 string cNametrP = cName +
"Momentum/TrueVsRndm_plus_" + stepname;
2470 TCanvas* cP =
fHMean.fHM.CreateCanvas(cNametrP, cNametrP, 2400, 800);
2473 DrawH1({hPlusXt, hPlusXr}, {
"P_{x} true",
"P_{x} random"},
kLinear,
kLog,
true, 0.75, 0.75, 0.91, 0.9,
"hist");
2474 DrawTextOnPad(
"True vs. random Momentum (e^{+})", 0.17, 0.9, 0.83, 0.99);
2477 DrawH1({hPlusYt, hPlusYr}, {
"P_{y} true",
"P_{y} random"},
kLinear,
kLog,
true, 0.75, 0.75, 0.91, 0.9,
"hist");
2478 DrawTextOnPad(
"True vs. random Momentum (e^{+})", 0.17, 0.9, 0.83, 0.99);
2481 DrawH1({hPlusZt, hPlusZr}, {
"P_{z} true",
"P_{z} random"},
kLinear,
kLog,
true, 0.75, 0.75, 0.91, 0.9,
"hist");
2482 DrawTextOnPad(
"True vs. random Momentum (e^{+})", 0.17, 0.9, 0.83, 0.99);
2485 string cNametrM = cName +
"Momentum/TrueVsRndm_minus_" + stepname;
2486 TCanvas* cM =
fHMean.fHM.CreateCanvas(cNametrM, cNametrM, 2400, 800);
2489 DrawH1({hMinusXt, hMinusXr}, {
"P_{x} true",
"P_{x} random"},
kLinear,
kLog,
true, 0.75, 0.75, 0.91, 0.9,
"hist");
2490 DrawTextOnPad(
"True vs. random Momentum (e^{-})", 0.17, 0.9, 0.83, 0.99);
2493 DrawH1({hMinusYt, hMinusYr}, {
"P_{y} true",
"P_{y} random"},
kLinear,
kLog,
true, 0.75, 0.75, 0.91, 0.9,
"hist");
2494 DrawTextOnPad(
"True vs. random Momentum (e^{-})", 0.17, 0.9, 0.83, 0.99);
2497 DrawH1({hMinusZt, hMinusZr}, {
"P_{z} true",
"P_{z} random"},
kLinear,
kLog,
true, 0.75, 0.75, 0.91, 0.9,
"hist");
2498 DrawTextOnPad(
"True vs. random Momentum (e^{-})", 0.17, 0.9, 0.83, 0.99);
2502 TH1D* hPlusXt2 = (TH1D*) hPlusXt->Clone();
2503 TH1D* hPlusYt2 = (TH1D*) hPlusYt->Clone();
2504 TH1D* hPlusZt2 = (TH1D*) hPlusZt->Clone();
2505 TH1D* hMinusXt2 = (TH1D*) hMinusXt->Clone();
2506 TH1D* hMinusYt2 = (TH1D*) hMinusYt->Clone();
2507 TH1D* hMinusZt2 = (TH1D*) hMinusZt->Clone();
2508 hPlusXt2->GetYaxis()->SetRangeUser(
min,
max);
2509 hPlusYt2->GetYaxis()->SetRangeUser(
min,
max);
2510 hPlusZt2->GetYaxis()->SetRangeUser(
min,
max);
2511 string cNamePm = cName +
"Momentum/PlusVsMinus_" + stepname;
2512 TCanvas* cPm =
fHMean.fHM.CreateCanvas(cNamePm, cNamePm, 2400, 800);
2515 DrawH1({hPlusXt2, hMinusXt2}, {
"e^{+}",
"e^{-}"},
kLinear,
kLog,
true, 0.8, 0.75, 0.91, 0.9,
"hist");
2519 DrawH1({hPlusYt2, hMinusYt2}, {
"e^{+}",
"e^{-}"},
kLinear,
kLog,
true, 0.8, 0.75, 0.91, 0.9,
"hist");
2523 DrawH1({hPlusZt2, hMinusZt2}, {
"e^{+}",
"e^{-}"},
kLinear,
kLog,
true, 0.8, 0.75, 0.91, 0.9,
"hist");
2531 TCanvas* c =
fHMean.fHM.CreateCanvas(cName +
"Background/Background", cName +
"Background/Background", 2400, 800);
2535 string stepname =
fHMean.fAnaStepNames[
static_cast<int>(step)];
2537 bg->GetYaxis()->SetRangeUser(1e-9, 1e-2);
2540 fHMean.DrawAnaStepOnPad(step);
2545 for (
const string cat : {
""}) {
2547 fHMean.fHM.CreateCanvas(cName +
"Background/VsMcBg" + cat, cName +
"Background/VsMcBg" + cat, 1600, 2400);
2551 string stepname =
fHMean.fAnaStepNames[
static_cast<int>(step)];
2554 bfs->GetYaxis()->SetRangeUser(1e-7, 2e-2);
2555 TH1D* r = (TH1D*) bfs->Clone();
2557 r->GetYaxis()->SetRangeUser(0.5, 1.4);
2558 r->GetYaxis()->SetTitle(
"Ratio");
2560 DrawH1({bfs, bmc}, {
"B_{FS}",
"B_{MC}"},
kLinear,
kLog,
true, 0.65, 0.75, 0.95, 0.9,
"p");
2561 fHMean.DrawAnaStepOnPad(step);
2564 fHMean.DrawAnaStepOnPad(step);
2587 fHMean.fHM.CreateCanvas(cName +
"Background/SmoothEffect", cName +
"Background/SmoothEffect", 2400, 1600);
2589 vector<int> canvas = {1, 4, 2, 5, 3, 6};
2591 for (
const string cat : {
""}) {
2592 TH1D* mc =
fHMean.H1Clone(
"hMinv" + cat +
"_bg", step);
2593 TH1D* fs =
fHFastSim.H1Clone(
"hfsc_MinvBg" + cat, step);
2594 TH1D* cb =
fHMean.H1Clone(
"hCbBg" + cat, step);
2595 TH1D* rat = (TH1D*) cb->Clone();
2598 rat2->GetYaxis()->SetRangeUser(0.97, 1.1);
2599 rat2->GetYaxis()->SetTitle(
"CB_{calc} / B_{FS}");
2600 mc->GetXaxis()->SetTitle(
"M_{ee} [GeV/c^{2}]");
2601 mc->GetYaxis()->SetRangeUser(1e-9, 1e-2);
2602 string label = (cat ==
"") ?
"PLUTO + UrQMD" : (cat ==
"_urqmdAll") ?
"UrQMD" :
"UrQMD electrons";
2605 DrawH1({mc, fs, cb}, {
"MC Background",
"Fast Sim. Background",
"CB_{calc}"},
kLinear,
kLog,
true, 0.52, 0.7, 0.91,
2607 DrawH1({mc, fs}, {
"MC Background",
"Fast Sim. Background"},
kLinear,
kLog,
true, 0.52, 0.77, 0.91, 0.9,
"hist");
2609 fHMean.DrawSimDataLabel(0.14, 1.02e-8, 0.033);
2619 for (
const string cat : {
""}) {
2620 string cNameFull = cName +
"Signal/signal" + cat;
2621 TCanvas* c =
fHMean.fHM.CreateCanvas(cNameFull, cNameFull, 2400, 800);
2628 TH1D* coc =
VaryBinWidth(
"cocktail" + cat +
"_" +
fHMean.fAnaStepNames[
static_cast<int>(step)],
"",
2632 sig->GetYaxis()->SetRangeUser(1e-9, 3e-2);
2634 DrawH1({sig, coc}, {
"Signal = B_{FS} + Cocktail - CB_{calc}",
"Cocktail"},
kLinear,
kLog,
true, 0.47, 0.77,
2636 fHMean.DrawAnaStepOnPad(step);
2637 fHMean.DrawSimDataLabel(.08, 3e-9, 0.033);
2660 vector<pair<double, double>> bwvar =
fBwVarSig3;
2661 for (
const string cat : {
""}) {
2662 string cNameFull = cName +
"Signal/Reduced/signal_reduced" + cat;
2663 string legSig =
"B_{FS} + Cocktail - CB_{calc} - #pi^{0} - #eta - #omega - #omega_{D} - #phi ";
2664 TCanvas* c =
fHMean.fHM.CreateCanvas(cNameFull, cNameFull, 3200, 800);
2668 string stepname =
fHMean.fAnaStepNames[
static_cast<int>(step)];
2669 TH1D* sig =
VaryBinWidth(
"fastsim",
"hfsc_MinvBg" + cat, step, bwvar);
2670 TH1D* fs =
VaryBinWidth(
"fastsim",
"hfsc_MinvBg" + cat, step, bwvar);
2680 TH1D* coc =
VaryBinWidth(
"cocktail" + cat +
"_" +
fHMean.fAnaStepNames[
static_cast<int>(step)],
"", bwvar);
2690 sig->Add(omega, -1.);
2691 sig->Add(omegadalitz, -1.);
2693 TH1D* cocRed = (TH1D*) inmed->Clone();
2697 sig->GetYaxis()->SetRangeUser(1e-10, 5e-2);
2698 fs->GetYaxis()->SetRangeUser(1e-10, 5e-2);
2701 DrawH1({fs, cbc, coc, omega, omegadalitz, phi, pi0, eta, gamma},
2702 {
"B_{FS}",
"CB_{calc}",
"Cocktail",
"#omega",
"#omega_{D}",
"#phi",
"#pi^{0}",
"#eta",
"#gamma"},
2706 DrawH1({sig, cocRed}, {legSig,
"#rho_{i.m.} + QGP"},
kLinear,
kLog,
true, 0.38, 0.77, 0.91, 0.9,
"pe");
2707 DrawTextOnPad(
"Reduced Cocktail (" + stepname +
")", 0.1, 0.9, 0.9, 0.99);
2716 fHMean.fHM.CreateCanvas(cName +
"Time", cName +
"Time", 800, 800);
2719 h->Fit(
"gaus",
"Q",
"", 0., 35.);
2720 TF1* func =
h->GetFunction(
"gaus");
2721 double mean = (func !=
nullptr) ? func->GetParameter(
"Mean") : 0.;
2727 TH1D* hRes =
fHFastSim.H1Clone(
"hfs_memory_res");
2728 TH1D* hVir =
fHFastSim.H1Clone(
"hfs_memory_vir");
2729 fHMean.fHM.CreateCanvas(cName +
"Memory", cName +
"Memory", 800, 800);
2730 DrawH1({hRes, hVir}, {
"residential",
"virtual"},
kLinear,
kLog,
true, 0.8, 0.75, 0.91, 0.9,
"hist");
2765 TH1D* pmSame = HM.
H1Clone(hName +
"PM_sameEv", step);
2766 TH1D* ppSame = HM.
H1Clone(hName +
"PP_sameEv", step);
2767 TH1D* mmSame = HM.
H1Clone(hName +
"MM_sameEv", step);
2768 TH1D* pmMixed = HM.
H1Clone(hName +
"PM_mixedEv", step);
2769 TH1D* ppMixed = HM.
H1Clone(hName +
"PP_mixedEv", step);
2770 TH1D* mmMixed = HM.
H1Clone(hName +
"MM_mixedEv", step);
2771 TH1D* geomSame = (TH1D*) pmSame->Clone();
2772 TH1D* geomMixed = (TH1D*) pmSame->Clone();
2773 TH1D* k = (TH1D*) pmMixed->Clone();
2774 TH1D* cb = (TH1D*) pmSame->Clone();
2778 for (
int iB = 1; iB <= geomSame->GetNbinsX(); iB++) {
2779 double cppSame = ppSame->GetBinContent(iB);
2780 double cmmSame = mmSame->GetBinContent(iB);
2781 double cSame = std::sqrt(cppSame * cmmSame);
2782 geomSame->SetBinContent(iB, cSame);
2783 double cppMix = ppMixed->GetBinContent(iB);
2784 double cmmMix = mmMixed->GetBinContent(iB);
2785 double cMix = std::sqrt(cppMix * cmmMix);
2786 geomMixed->SetBinContent(iB, cMix);
2788 k->Divide(geomMixed);
2790 int ChangeBin = pmSame->FindBin(
fCbChange);
2791 int NormRangeMinBin = pmSame->FindBin(rangeStart);
2792 int NormRangeMaxBin = pmSame->FindBin(rangeEnd);
2794 ppSame->Integral(NormRangeMinBin, NormRangeMaxBin) / ppMixed->Integral(NormRangeMinBin, NormRangeMaxBin);
2796 mmSame->Integral(NormRangeMinBin, NormRangeMaxBin) / mmMixed->Integral(NormRangeMinBin, NormRangeMaxBin);
2797 for (
int iB = 1; iB <= cb->GetNbinsX(); iB++) {
2798 double kval = k->GetBinContent(iB);
2799 double pm = pmMixed->GetBinContent(iB);
2800 double cbval = (iB < ChangeBin) ? 2. * kval * geomSame->GetBinContent(iB) : pm * std::sqrt(normPPInt * normMMInt);
2801 cb->SetBinContent(iB, cbval);
2805 double bW = cb->GetBinWidth(1);
2806 for (
int iB = 1; iB <= cb->GetNbinsX(); iB++) {
2808 double s_pp = ppSame->GetBinContent(iB) *
fNofSimEvents * bW;
2809 double s_mm = mmSame->GetBinContent(iB) *
fNofSimEvents * bW;
2810 double m_pm = pmMixed->GetBinContent(iB) *
fNofSimEvents * bW;
2811 double m_pp = pmMixed->GetBinContent(iB) *
fNofSimEvents * bW;
2812 double m_mm = pmMixed->GetBinContent(iB) *
fNofSimEvents * bW;
2815 double d_m_pm = std::sqrt(s_pp * s_mm / (m_pp * m_mm));
2816 double d_m_pp = -0.5 * m_pm * std::pow(m_mm, -1. / 2.) * std::pow(m_pp, -3. / 2.) * std::sqrt(s_pp * s_mm);
2817 double d_m_mm = -0.5 * m_pm * std::pow(m_mm, -3. / 2.) * std::pow(m_pp, -1. / 2.) * std::sqrt(s_pp * s_mm);
2818 double d_s_pp = 0.5 * m_pm * std::sqrt(s_mm) * std::pow(m_pp * m_mm * s_pp, -1. / 2.);
2819 double d_s_mm = 0.5 * m_pm * std::sqrt(s_pp) * std::pow(m_pp * m_mm * s_mm, -1. / 2.);
2822 double f_m_pm = std::pow(d_m_pm * std::sqrt(m_pm), 2);
2823 double f_m_pp = std::pow(d_m_pp * std::sqrt(m_pp), 2);
2824 double f_m_mm = std::pow(d_m_mm * std::sqrt(m_mm), 2);
2825 double f_s_pp = std::pow(d_s_pp * std::sqrt(s_pp), 2);
2826 double f_s_mm = std::pow(d_s_mm * std::sqrt(s_mm), 2);
2829 double errorBc = std::sqrt(f_m_pm + f_m_pp + f_m_mm + f_s_pp + f_s_mm) / (
fNofSimEvents * bW);
2830 double errorBc2 = std::sqrt(normPPInt * normMMInt) * std::sqrt(m_pm) / (
fNofSimEvents * bW);
2831 if (iB < ChangeBin) cb->SetBinError(iB, errorBc);
2832 if (iB >= ChangeBin) cb->SetBinError(iB, errorBc2);
2839 string hName =
"hfsc_Minv";
2841 TH1D* pmSame =
fHFastSim.H1Clone(hName +
"Bg", step);
2842 TH1D* ppSame =
fHFastSim.H1Clone(hName +
"PP_sameEv", step);
2843 TH1D* mmSame =
fHFastSim.H1Clone(hName +
"MM_sameEv", step);
2844 TH1D* pmMixed =
fHFastSim.H1Clone(hName +
"PM_mixedEv", step);
2845 TH1D* ppMixed =
fHFastSim.H1Clone(hName +
"PP_mixedEv", step);
2846 TH1D* mmMixed =
fHFastSim.H1Clone(hName +
"MM_mixedEv", step);
2847 TH1D* geomSame = (TH1D*) pmSame->Clone();
2848 TH1D* geomMixed = (TH1D*) pmSame->Clone();
2849 TH1D* k = (TH1D*) pmMixed->Clone();
2850 TH1D* cb = (TH1D*) pmSame->Clone();
2854 k->GetYaxis()->SetTitle(
"#it{k} Factor");
2855 for (
int iB = 1; iB <= geomSame->GetNbinsX(); iB++) {
2856 double cppSame = ppSame->GetBinContent(iB);
2857 double cmmSame = mmSame->GetBinContent(iB);
2858 double cSame = std::sqrt(cppSame * cmmSame);
2859 geomSame->SetBinContent(iB, cSame);
2860 double cppMix = ppMixed->GetBinContent(iB);
2861 double cmmMix = mmMixed->GetBinContent(iB);
2862 double cMix = std::sqrt(cppMix * cmmMix);
2863 geomMixed->SetBinContent(iB, cMix);
2865 k->Divide(geomMixed);
2867 int ChangeBin = pmSame->FindBin(
fCbChange);
2871 ppSame->Integral(NormRangeMinBin, NormRangeMaxBin) / ppMixed->Integral(NormRangeMinBin, NormRangeMaxBin);
2873 mmSame->Integral(NormRangeMinBin, NormRangeMaxBin) / mmMixed->Integral(NormRangeMinBin, NormRangeMaxBin);
2874 for (
int iB = 1; iB <= cb->GetNbinsX(); iB++) {
2875 double kval = k->GetBinContent(iB);
2876 double geomMeanMixNormed = geomMixed->GetBinContent(iB) * std::sqrt(normPPInt * normMMInt);
2877 double cbval = (iB < ChangeBin) ? 2. * kval * geomSame->GetBinContent(iB) : 2. * kval * geomMeanMixNormed;
2878 cb->SetBinContent(iB, cbval);
2882 double bW = cb->GetBinWidth(1);
2883 for (
int iB = 1; iB <= cb->GetNbinsX(); iB++) {
2885 double s_pp = ppSame->GetBinContent(iB) *
fNofSimEvents * bW;
2886 double s_mm = mmSame->GetBinContent(iB) *
fNofSimEvents * bW;
2887 double m_pm = pmMixed->GetBinContent(iB) *
fNofSimEvents * bW;
2888 double m_pp = pmMixed->GetBinContent(iB) *
fNofSimEvents * bW;
2889 double m_mm = pmMixed->GetBinContent(iB) *
fNofSimEvents * bW;
2892 double d_m_pm = std::sqrt(s_pp * s_mm / (m_pp * m_mm));
2893 double d_m_pp = -0.5 * m_pm * std::pow(m_mm, -1. / 2.) * std::pow(m_pp, -3. / 2.) * std::sqrt(s_pp * s_mm);
2894 double d_m_mm = -0.5 * m_pm * std::pow(m_mm, -3. / 2.) * std::pow(m_pp, -1. / 2.) * std::sqrt(s_pp * s_mm);
2895 double d_s_pp = 0.5 * m_pm * std::sqrt(s_mm) * std::pow(m_pp * m_mm * s_pp, -1. / 2.);
2896 double d_s_mm = 0.5 * m_pm * std::sqrt(s_pp) * std::pow(m_pp * m_mm * s_mm, -1. / 2.);
2899 double f_m_pm = std::pow(d_m_pm * std::sqrt(m_pm), 2);
2900 double f_m_pp = std::pow(d_m_pp * std::sqrt(m_pp), 2);
2901 double f_m_mm = std::pow(d_m_mm * std::sqrt(m_mm), 2);
2902 double f_s_pp = std::pow(d_s_pp * std::sqrt(s_pp), 2);
2903 double f_s_mm = std::pow(d_s_mm * std::sqrt(s_mm), 2);
2906 double errorBc = std::sqrt(f_m_pm + f_m_pp + f_m_mm + f_s_pp + f_s_mm) / (
fNofSimEvents * bW);
2907 double errorBc2 = std::sqrt(normPPInt * normMMInt) * std::sqrt(m_pm) / (
fNofSimEvents * bW);
2908 if (iB < ChangeBin) cb->SetBinError(iB, errorBc);
2909 if (iB >= ChangeBin) cb->SetBinError(iB, errorBc2);
2911 fHMean.fHM.Add(
"hReCbGeom_sameEv", geomSame);
2912 fHMean.fHM.Add(
"hReCbGeom_mixedEv", geomMixed);
2913 fHMean.fHM.Add(
"hReCbK", k);
2914 fHMean.fHM.Add(
"hReCb", cb);
2915 string cName =
"FastSim/CombinatorialBackground/BG";
2916 TCanvas* c =
fHMean.fHM.CreateCanvas(cName, cName, 3200, 800);
2919 DrawH1({pmSame, ppSame, mmSame}, {
"e+e-",
"e+e+",
"e-e-"},
kLinear,
kLog,
true, 0.8, 0.75, 0.91, 0.9,
"hist");
2922 DrawH1({pmMixed, ppMixed, mmMixed}, {
"e+e-",
"e+e+",
"e-e-"},
kLinear,
kLog,
true, 0.8, 0.75, 0.91, 0.9,
"hist");
2929 DrawTextOnPad(
"Combinatorial Background from Fast Simulations", 0.1, 0.9, 0.84, 0.99);
2935 string stepname =
fHMean.fAnaStepNames[
static_cast<int>(step)];
2937 vector<TH1*> HistsCB;
2938 vector<TH1*> HistsRE;
2939 vector<TH1*> HistsSig;
2941 for (
const string cutset : {
"cut_standard_8",
"cut_looser_10",
"cut_stricter_17"}) {
2942 string dir =
"/lustre/cbm/users/criesen/data/studies/RandomEvent/SysErrorSignal/" + cutset +
"/";
2943 TH1D* coc =
nullptr;
2944 TH1D* cocPluto =
nullptr;
2945 TH1D* cocUrqmd =
nullptr;
2950 for (
const string pf : {
"inmed",
"qgp",
"omegadalitz",
"omegaepem",
"phi"}) {
2952 TFile* filePluto =
new TFile((dir + pf +
"/analysis.all.root").c_str());
2954 int nEventsSig = HM.
H1(
"hEventNumber")->GetEntries();
2955 nEvents += nEventsSig;
2957 TH1D* pluto = (TH1D*) HM.
fHM.
GetObject(
"hMinv_signal_" + stepname)->Clone();
2958 TH1D* eta = (TH1D*) HM.
fHM.
GetObject(
"hMinv_eta_" + stepname)->Clone();
2959 TH1D* pi0 = (TH1D*) HM.
fHM.
GetObject(
"hMinv_pi0_" + stepname)->Clone();
2960 bW = pluto->GetBinWidth(1);
2961 pluto->Scale(1. / (nEventsSig * bW));
2962 eta->Scale(1. / bW);
2963 pi0->Scale(1. / bW);
2964 if (cocPluto ==
nullptr)
2967 cocPluto->Add(pluto);
2968 if (cocUrqmd ==
nullptr)
2975 cocUrqmd->Scale(1. / nEvents);
2981 TFile* fileRandom =
new TFile((dir +
"randomevent.all.root").c_str());
2983 double nofFastSimEv = HMFastSim.
H1(
"hre_EventNumber")->GetEntries();
2985 TH1D* bre = HMFastSim.
H1Clone(
"hrec_MinvBg");
2989 TFile* filecb =
new TFile((dir +
"evmix.all.root").c_str());
2992 string hName =
"hmx_Minv";
2993 TH1D* pmSame = HMCB.
H1Clone(hName +
"PM_sameEv", step);
2994 TH1D* ppSame = HMCB.
H1Clone(hName +
"PP_sameEv", step);
2995 TH1D* mmSame = HMCB.
H1Clone(hName +
"MM_sameEv", step);
2996 TH1D* pmMixed = HMCB.
H1Clone(hName +
"PM_mixedEv", step);
2997 TH1D* ppMixed = HMCB.
H1Clone(hName +
"PP_mixedEv", step);
2998 TH1D* mmMixed = HMCB.
H1Clone(hName +
"MM_mixedEv", step);
2999 TH1D* geomSame = (TH1D*) pmSame->Clone();
3000 TH1D* geomMixed = (TH1D*) pmSame->Clone();
3001 TH1D* k = (TH1D*) pmMixed->Clone();
3002 TH1D* cb = (TH1D*) pmSame->Clone();
3007 for (
int iB = 1; iB <= geomSame->GetNbinsX(); iB++) {
3008 double cppSame = ppSame->GetBinContent(iB);
3009 double cmmSame = mmSame->GetBinContent(iB);
3010 double cSame = std::sqrt(cppSame * cmmSame);
3011 geomSame->SetBinContent(iB, cSame);
3012 double cppMix = ppMixed->GetBinContent(iB);
3013 double cmmMix = mmMixed->GetBinContent(iB);
3014 double cMix = std::sqrt(cppMix * cmmMix);
3015 geomMixed->SetBinContent(iB, cMix);
3017 k->Divide(geomMixed);
3019 int ChangeBin = pmSame->FindBin(
fCbChange);
3023 ppSame->Integral(NormRangeMinBin, NormRangeMaxBin) / ppMixed->Integral(NormRangeMinBin, NormRangeMaxBin);
3025 mmSame->Integral(NormRangeMinBin, NormRangeMaxBin) / mmMixed->Integral(NormRangeMinBin, NormRangeMaxBin);
3027 for (
int iB = 1; iB <= cb->GetNbinsX(); iB++) {
3028 double kval = k->GetBinContent(iB);
3029 double pp = ppMixed->GetBinContent(iB);
3030 double mm = mmMixed->GetBinContent(iB);
3031 double cbval = (iB < ChangeBin) ? 2. * kval * geomSame->GetBinContent(iB)
3032 : 2. * kval * std::sqrt(normPPInt * pp * normMMInt * mm);
3033 cb->SetBinContent(iB, cbval);
3037 for (
int iB = 1; iB <= cb->GetNbinsX(); iB++) {
3039 double s_pp = ppSame->GetBinContent(iB) *
fNofSimEvents * bW;
3040 double s_mm = mmSame->GetBinContent(iB) *
fNofSimEvents * bW;
3041 double m_pm = pmMixed->GetBinContent(iB) *
fNofSimEvents * bW;
3042 double m_pp = pmMixed->GetBinContent(iB) *
fNofSimEvents * bW;
3043 double m_mm = pmMixed->GetBinContent(iB) *
fNofSimEvents * bW;
3046 double d_m_pm = std::sqrt(s_pp * s_mm / (m_pp * m_mm));
3047 double d_m_pp = -0.5 * m_pm * std::pow(m_mm, -1. / 2.) * std::pow(m_pp, -3. / 2.) * std::sqrt(s_pp * s_mm);
3048 double d_m_mm = -0.5 * m_pm * std::pow(m_mm, -3. / 2.) * std::pow(m_pp, -1. / 2.) * std::sqrt(s_pp * s_mm);
3049 double d_s_pp = 0.5 * m_pm * std::sqrt(s_mm) * std::pow(m_pp * m_mm * s_pp, -1. / 2.);
3050 double d_s_mm = 0.5 * m_pm * std::sqrt(s_pp) * std::pow(m_pp * m_mm * s_mm, -1. / 2.);
3053 double f_m_pm = std::pow(d_m_pm * std::sqrt(m_pm), 2);
3054 double f_m_pp = std::pow(d_m_pp * std::sqrt(m_pp), 2);
3055 double f_m_mm = std::pow(d_m_mm * std::sqrt(m_mm), 2);
3056 double f_s_pp = std::pow(d_s_pp * std::sqrt(s_pp), 2);
3057 double f_s_mm = std::pow(d_s_mm * std::sqrt(s_mm), 2);
3060 double errorBc = std::sqrt(f_m_pm + f_m_pp + f_m_mm + f_s_pp + f_s_mm) / (
fNofSimEvents * bW);
3062 std::sqrt(normPPInt * normMMInt) * std::sqrt(m_pm) / (
fNofSimEvents * bW);
3063 if (iB < ChangeBin) cb->SetBinError(iB, errorBc);
3064 if (iB >= ChangeBin) cb->SetBinError(iB, errorBc2);
3068 TH1D* sig = (TH1D*) bre->Clone();
3071 for (
int iB = 1; iB <= sig->GetNbinsX(); iB++) {
3072 double D_bre = bre->GetBinError(iB);
3073 double D_cbc = cb->GetBinError(iB);
3074 double D_dNdM = std::sqrt(D_bre * D_bre + D_cbc * D_cbc);
3075 sig->SetBinError(iB, D_dNdM);
3078 TH1D* ratio = (TH1D*) sig->Clone();
3081 ratio->GetYaxis()->SetRangeUser(-2., 2.);
3082 ratio->GetYaxis()->SetTitle(
"Ratio Sig_{RE} / MC Cocktail");
3084 HistsSig.push_back(ratio);
3085 HistsCB.push_back(cb);
3086 HistsRE.push_back(bre);
3089 fHMean.fHM.CreateCanvas(
"FastSim/SysErrors/SigToCock",
"FastSim/SysErrors/SigToCock", 800, 800);
3090 DrawH1(HistsSig, {
"cutset_0",
"cutset_looser",
"cutset_stricter"},
kLinear,
kLinear,
true, 0.7, 0.8, 0.95, 0.91,
3095 TH1D* ratCbLoose = (TH1D*) HistsCB[1];
3096 TH1D* ratCbStrict = (TH1D*) HistsCB[2];
3097 ratCbLoose->Divide(HistsCB[0]);
3098 ratCbStrict->Divide(HistsCB[0]);
3100 TH1D* ratReLoose = (TH1D*) HistsRE[1];
3101 TH1D* ratReStrict = (TH1D*) HistsRE[2];
3102 ratReLoose->Divide(HistsRE[0]);
3103 ratReStrict->Divide(HistsRE[0]);
3110 ratCbLoose->GetYaxis()->SetRangeUser(0.5, 1.5);
3111 ratReLoose->GetYaxis()->SetRangeUser(0.5, 1.5);
3112 ratCbLoose->GetYaxis()->SetTitle(
"Ratio");
3113 ratReLoose->GetYaxis()->SetTitle(
"Ratio");
3115 TCanvas* c =
fHMean.fHM.CreateCanvas(
"FastSim/SysErrors/Backgrounds",
"FastSim/SysErrors/Backgrounds", 1600, 800);
3118 DrawH1({ratCbLoose, ratCbStrict}, {
"loose / standard",
"strict / standard"},
kLinear,
kLinear,
true, 0.7, 0.8, 0.95,
3120 DrawTextOnPad(
"Combinatorial Background", 0.25, 0.9, 0.75, 0.99);
3122 DrawH1({ratReLoose, ratReStrict}, {
"loose / standard",
"strict / standard"},
kLinear,
kLinear,
true, 0.7, 0.8, 0.95,
3124 DrawTextOnPad(
"Fast Simulation Background", 0.25, 0.9, 0.75, 0.99);
3129 TH1D* r0 = (TH1D*) HistsRE[0];
3130 TH1D* rL = (TH1D*) HistsRE[1];
3131 TH1D* rS = (TH1D*) HistsRE[2];
3132 r0->Divide(HistsCB[0]);
3133 rL->Divide(HistsCB[1]);
3134 rS->Divide(HistsCB[2]);
3140 r0->GetYaxis()->SetRangeUser(0.97, 1.03);
3141 r0->GetYaxis()->SetTitle(
"Ratio");
3143 fHMean.fHM.CreateCanvas(
"FastSim/SysErrors/BreOverCbcalc",
"FastSim/SysErrors/BreOverCbcalc", 800, 800);
3144 DrawH1({r0, rL, rS}, {
"cutset_0",
"cutset_looser",
"cutset_stricter"},
kLinear,
kLinear,
true, 0.7, 0.8, 0.95, 0.91,
3154 TH1D*
h =
new TH1D(
"hNofEvents",
"hNofEvents; PLUTO File; Number",
fH.size(), 0.,
fH.size());
3157 for (
size_t iS = 0; iS <
fH.size(); iS++) {
3159 double c =
fH[iS]->H1(
"hEventNumber")->GetEntries();
3162 h->SetBinContent(iB, c);
3163 h->GetXaxis()->SetBinLabel(iB,
fHMean.fSignalNames[iS].c_str());
3165 h->GetYaxis()->SetRangeUser(0.9 *
min, 1.1 *
max);
3166 fHMean.fHM.CreateCanvas(
"Numbers/NofEvents",
"Numbers/NofEvents", 800, 800);
3168 h->SetMarkerSize(1.);
3174 vector<string> xLabel = {
"e^{#pm}",
"#pi^{#pm}",
"#pi^{0}",
"p",
"n",
"#gamma_{ch}",
"o."};
3175 fHMean.fHM.CreateCanvas(
"Numbers/MCTracks",
"Numbers/MCTracks", 800, 800);
3176 TH1D*
h =
fHMean.H1Clone(
"hNofMcTracks");
3178 for (
size_t iL = 0; iL < xLabel.size(); iL++) {
3179 h->GetXaxis()->SetBinLabel(iL + 1, xLabel[iL].c_str());
3181 double nAll =
h->Integral(1,
h->GetNbinsX());
3182 double ratEl = 100 *
h->GetBinContent(1) / nAll;
3183 double ratPiC = 100 *
h->GetBinContent(2) / nAll;
3184 double ratPi0 = 100 *
h->GetBinContent(3) / nAll;
3185 double ratProt = 100 *
h->GetBinContent(4) / nAll;
3186 double ratNeut = 100 *
h->GetBinContent(5) / nAll;
3187 double ratGCh = 100 *
h->GetBinContent(6) / nAll;
3188 double ratioO = 100 *
h->GetBinContent(7) / nAll;
3201 vector<string> xLabel = {
"e_{Pluto}",
"e_{#gamma}",
"e_{o.}",
"#pi^{#pm}",
"p",
"o."};
3202 fHMean.fHM.CreateCanvas(
"Numbers/GlobalTracks",
"Numbers/GlobalTracks", 800, 800);
3203 TH1D*
h =
fHMean.H1Clone(
"hNofGTracks");
3205 for (
size_t iL = 0; iL < xLabel.size(); iL++) {
3206 h->GetXaxis()->SetBinLabel(iL + 1, xLabel[iL].c_str());
3208 double nAll =
h->Integral(1,
h->GetNbinsX());
3209 double ratPlutoEl = 100 *
h->GetBinContent(1) / nAll;
3210 double ratGammaEl = 100 *
h->GetBinContent(2) / nAll;
3211 double ratOtherEl = 100 *
h->GetBinContent(3) / nAll;
3212 double ratPiC = 100 *
h->GetBinContent(4) / nAll;
3213 double ratProton = 100 *
h->GetBinContent(5) / nAll;
3214 double ratioO = 100 *
h->GetBinContent(6) / nAll;
3226 gStyle->SetPaintTextFormat(
"4.2f");
3227 fHMean.fHM.CreateCanvas(
"Numbers/Candidates",
"Numbers/Candidates", 1100, 800);
3228 TH1D*
h =
fHMean.H1Clone(
"hNofBgTracks");
3230 h->SetMarkerSize(1.5);
3231 fHMean.SetAnalysisStepAxis(
h);
3426 fHMean.DrawAllCandsAndSteps(2,
"STS/E-Loss/",
"hELossSts", {
""}, {
""},
false,
false, 1e-12, 1e-2);
3427 fHMean.DrawAllCandsAndSteps(2,
"RICH/ANN/",
"hRichAnn", {
""}, {
""},
false,
false, 1e-13, 1.);
3428 fHMean.DrawAllCandsAndSteps(2,
"TRD/El-Like/",
"hTrdElLike", {
""}, {
""},
false,
false, 1e-13, 1.);
3429 fHMean.DrawAllCandsAndSteps(2,
"TRD/Chi2/",
"hChi2Dets_trd", {
""}, {
""},
false,
false, 1e-13, 1e-2);
3430 fHMean.DrawAllCandsAndSteps(2,
"TOF/M2/",
"hTofM2", {
""}, {
""},
false,
false, 1e-12, 1e-1);
3431 fHMean.DrawAllCandsAndSteps(2,
"TOF/Dist/",
"hTofDist", {
""}, {
""},
false,
false, 1e-12, 1e-1);
3432 fHMean.DrawAllCands(2,
"TOF/Dist/truematch",
"hTofDist_truematch", {
""}, {
""},
false,
false, 1e-12, 1e-1);
3433 fHMean.DrawAllCands(2,
"TOF/Dist/mismatch",
"hTofDist_mismatch", {
""}, {
""},
false,
false, 1e-12, 1e-1);
3441 TH2D* stsELoss =
fHMean.GetRatioElectronsOthers(
"hELossSts", step);
3442 TH2D* richAnn =
fHMean.GetRatioElectronsOthers(
"hRichAnn", step);
3443 TH2D* trdLike =
fHMean.GetRatioElectronsOthers(
"hTrdElLike", step);
3444 TH2D* trdChi2 =
fHMean.GetRatioElectronsOthers(
"hChi2Dets_trd", step);
3445 TH2D* tofM2 =
fHMean.GetRatioElectronsOthers(
"hTofM2", step);
3446 TH2D* tofM2z =
fHMean.GetRatioElectronsOthers(
"hTofM2", step);
3447 TH2D* tofDist =
fHMean.GetRatioElectronsOthers(
"hTofDist", step);
3450 stsELoss->GetZaxis()->SetRangeUser(
min,
max);
3451 richAnn->GetZaxis()->SetRangeUser(
min,
max);
3452 trdLike->GetZaxis()->SetRangeUser(
min,
max);
3453 trdChi2->GetZaxis()->SetRangeUser(
min,
max);
3454 tofM2->GetZaxis()->SetRangeUser(
min,
max);
3455 tofM2z->GetZaxis()->SetRangeUser(
min,
max);
3456 tofM2z->GetXaxis()->SetRangeUser(0., 4.);
3457 tofM2z->GetYaxis()->SetRangeUser(-0.1, 0.1);
3458 string cName =
"ElidCuts/" +
fHMean.fAnaStepNames[
static_cast<int>(step)];
3459 TCanvas* c =
fHMean.fHM.CreateCanvas(cName, cName, 3200, 1600);
3478 DrawTextOnPad(
"TOF #it{m}^{2} (zoom)", 0.12, 0.9, 0.64, 0.99);
3488 string stepname =
fHMean.fAnaStepNames[
static_cast<int>(step)];
3489 string cName =
"TRD/El-Like/Steps/1D_" + stepname;
3490 TCanvas* c =
fHMean.fHM.CreateCanvas(cName, cName, 2400, 1800);
3493 for (
const string& ptcl :
fHMean.fCandNames) {
3494 string hName =
fHMean.GetName(
"hTrdElLike_" + ptcl, step);
3495 string projname = hName +
"_y";
3496 TH2D* h2 =
fHMean.H2Clone(hName);
3497 TH1D* h1 = (TH1D*) h2->ProjectionY(projname.c_str(), 1, h2->GetYaxis()->GetNbins(),
"");
3507 string cName =
"TOF/M2/Steps/" +
fHMean.fAnaStepNames[
static_cast<int>(step)] +
"_zoom";
3508 TCanvas* c =
fHMean.fHM.CreateCanvas(cName, cName, 3200, 2400);
3511 for (
const string& ptcl :
fHMean.fCandNames) {
3512 TH2D*
h =
fHMean.H2Clone(
"hTofM2_" + ptcl, step);
3513 h->GetXaxis()->SetRangeUser(0., 4.);
3514 h->GetYaxis()->SetRangeUser(-0.1, 0.1);
3515 h->GetZaxis()->SetRangeUser(1e-12, 1e-1);
3837 string cName =
"Efficiency/";
3838 string hName =
"hMom_";
3843 TCanvas* c1 =
fHMean.fHM.CreateCanvas(cName +
"hMinv", cName +
"hMinv", 2400, 1600);
3847 TH1D* mc =
H(signal)->
H1Clone(
"hMinv_signal_mc");
3848 TH1D* reco =
H(signal)->
H1Clone(
"hMinv_signal_reco");
3849 TH1D* elid =
H(signal)->
H1Clone(
"hMinv_signal_elid");
3853 TH1D* effMc = (TH1D*) elid->Clone();
3854 TH1D* effReco = (TH1D*) elid->Clone();
3856 effReco->Divide(reco);
3857 effMc->GetYaxis()->SetTitle(
"Efficiency");
3858 int nBinsXMinv = mc->GetXaxis()->GetNbins();
3859 double eMc = 100 * elid->Integral(1, nBinsXMinv) / mc->Integral(1, nBinsXMinv);
3860 double eReco = 100 * elid->Integral(1, nBinsXMinv) / reco->Integral(1, nBinsXMinv);
3864 effMc->GetYaxis()->SetRangeUser(0., 0.3);
3865 effReco->GetYaxis()->SetRangeUser(0., 0.3);
3866 DrawH1({effMc, effReco}, {
"Eff_{MC}: " + eMcString,
"Eff_{Reco}: " + eRecoString},
kLinear,
kLinear,
true, 0.6, 0.7,
3867 0.97, 0.92,
"hist");
3872 TH1D* hRef =
fHMean.H1Clone(hName +
fHMean.fCandNames[0], refstep);
3873 hRef->Add(
fHMean.H1(hName +
fHMean.fCandNames[1], refstep));
3881 gamma->Divide(hRef);
3891 int nBinsX = hRef->GetXaxis()->GetNbins();
3893 double effElid = 100 * elid->Integral(1, nBinsX) / nBinsX;
3894 double effGamma = 100 * gamma->Integral(1, nBinsX) / nBinsX;
3895 double effTt = 100 * tt->Integral(1, nBinsX) / nBinsX;
3896 double effPt = 100 * pt->Integral(1, nBinsX) / nBinsX;
3898 fHMean.fHM.CreateCanvas(cName +
"all_steps", cName +
"all_steps", 900, 900);
3899 elid->GetYaxis()->SetTitle(
"Efficiency");
3900 DrawH1({elid, gamma, tt, pt},
3908 string cNameFull = cName +
"plutos_" +
fHMean.fAnaStepNames[
static_cast<int>(step)];
3909 fHMean.fHM.CreateCanvas(cNameFull, cNameFull, 900, 900);
3911 vector<TH1D*> sHists(
fHMean.fNofSignals);
3913 TH1D* hRef2 =
H(signal)->
H1Clone(hName +
fHMean.fCandNames[0], refstep);
3914 hRef2->Add(
H(signal)->H1Clone(hName +
fHMean.fCandNames[1], refstep));
3915 TH1D* hEff =
H(signal)->
H1Clone(hName +
fHMean.fCandNames[0], step);
3916 hEff->Add(
H(signal)->H1Clone(hName +
fHMean.fCandNames[1], step));
3917 hEff->Divide(hRef2);
3918 hEff->GetYaxis()->SetTitle(
"Efficiency");
3919 sHists[
static_cast<int>(signal)] = hEff;
3922 double effInmed = 100 * sHists[0]->Integral(1, nBinsX) / nBinsX;
3923 double effQgp = 100 * sHists[1]->Integral(1, nBinsX) / nBinsX;
3924 double effOmega = 100 * sHists[2]->Integral(1, nBinsX) / nBinsX;
3925 double effPhi = 100 * sHists[3]->Integral(1, nBinsX) / nBinsX;
3926 double effOmegaD = 100 * sHists[4]->Integral(1, nBinsX) / nBinsX;
3928 vector<LmvmDrawMinvData> data;
3931 data.emplace_back(sHists[
static_cast<int>(
ELmvmSignal::Omega)], -1, kOrange + 7, 2, -1,
3933 data.emplace_back(sHists[
static_cast<int>(
ELmvmSignal::Phi)], -1, kAzure + 2, 2, -1,
3935 data.emplace_back(sHists[
static_cast<int>(
ELmvmSignal::Qgp)], -1, kOrange - 2, 2, -1,
3937 data.emplace_back(sHists[
static_cast<int>(
ELmvmSignal::Inmed)], -1, kMagenta - 3, 2, -1,
3940 TLegend* leg =
new TLegend(0.75, 0.65, .99, 0.99);
3941 for (
size_t i = 0; i < data.size(); i++) {
3942 const auto& d = data[i];
3943 d.fH->GetYaxis()->SetLabelSize(0.05);
3944 if (d.fFillColor != -1) d.fH->SetFillColor(d.fFillColor);
3945 if (d.fFillStyle != -1) d.fH->SetFillStyle(d.fFillStyle);
3946 leg->AddEntry(d.fH, d.fLegend.c_str(),
"f");
3947 DrawH1(d.fH,
kLinear,
kLinear, (i == 0) ?
"hist" :
"hist,same", d.fLineColor, d.fLineWidth, 0);
3949 DrawTextOnPad(
"Efficiency (average of each both electrons)", 0.01, 0.9, 0.8, 0.99);
3950 leg->SetFillColor(kWhite);
3954 vector<TH1D*> sHistsStepP(
fHMean.fNofSignals);
3955 vector<TH1D*> sHistsStepE(
fHMean.fNofSignals);
3956 vector<TH1D*> sHistsRefP(
fHMean.fNofSignals);
3957 vector<TH1D*> sHistsRefE(
fHMean.fNofSignals);
3959 TH1D* hStepP =
H(signal)->
H1Clone(hName +
fHMean.fCandNames[0], step);
3960 TH1D* hStepE =
H(signal)->
H1Clone(hName +
fHMean.fCandNames[1], step);
3961 TH1D* hRefP =
H(signal)->
H1Clone(hName +
fHMean.fCandNames[0], refstep);
3962 TH1D* hRefE =
H(signal)->
H1Clone(hName +
fHMean.fCandNames[1], refstep);
3963 sHistsStepP[
static_cast<int>(signal)] = hStepP;
3964 sHistsStepE[
static_cast<int>(signal)] = hStepE;
3965 sHistsRefP[
static_cast<int>(signal)] = hRefP;
3966 sHistsRefE[
static_cast<int>(signal)] = hRefE;
3970 string cNameFull2 = cName +
"Components_" +
fHMean.fSignalNames[
static_cast<int>(signal)];
3971 TCanvas* c =
fHMean.fHM.CreateCanvas(cNameFull2, cNameFull2, 1600, 800);
3974 DrawH1({sHistsStepP[
static_cast<int>(signal)], sHistsStepE[
static_cast<int>(signal)]}, {
"positron",
"electron"},
3975 kLinear,
kLog,
true, 0.85, 0.7, 0.99, 0.99,
"hist");
3976 fHMean.DrawAnaStepOnPad(step);
3978 DrawH1({sHistsRefP[
static_cast<int>(signal)], sHistsRefE[
static_cast<int>(signal)]}, {
"positron",
"electron"},
3979 kLinear,
kLog,
true, 0.85, 0.7, 0.99, 0.99,
"hist");
3980 fHMean.DrawAnaStepOnPad(refstep);
3985 vector<pair<double, double>> bwvar = {make_pair(10., 1.)};
3986 TH1D* hElRef =
VaryBinWidth(
"mean",
"hMom_plutoEl+_reco", bwvar);
3987 TH1D* hElCounter =
VaryBinWidth(
"mean",
"hMom_plutoEl+_elid", bwvar);
3988 for (
int iP = 1; iP <= 3; iP++) {
3989 string ptcl =
fHMean.fCandNames[iP];
3990 TH1D* href =
VaryBinWidth(
"mean",
"hMom_" + ptcl +
"_reco", bwvar);
3991 TH1D* hcounter =
VaryBinWidth(
"mean",
"hMom_" + ptcl +
"_elid", bwvar);
3993 hElCounter->Add(hcounter);
3995 TH1D* hEff = (TH1D*) hElCounter->Clone();
3996 hEff->Divide(hElRef);
3997 int nX = hElRef->GetNbinsX();
3998 double eff = hEff->Integral(1, nX) / nX;
3999 TCanvas* c =
fHMean.fHM.CreateCanvas(cName +
"efficiency_electrons", cName +
"Efficiency_electrons", 1600, 800);
4002 DrawH1({hElRef, hElCounter}, {
"Ref",
"Counter"},
kLinear,
kLog,
true, 0.85, 0.85, 0.95, 0.95,
"p");
4200 vector<string> pLabel = {
"e^{#pm}_{PLUTO}",
"e^{#pm}_{UrQMD}",
"e_{#gamma}",
"#pi^{#pm}",
"p",
"K^{#pm}",
"o."};
4204 for (
const string suff : {
"",
"+",
"-"}) {
4210 if (suff ==
"+" || suff ==
"-") {
4211 hName =
"hCandPdgVsMom" + suff;
4212 cName =
"Purity/purity" + suff;
4215 hName =
"hCandPdgVsMom_" +
fHMean.fAnaStepNames[
static_cast<int>(step)];
4216 cName =
"Purity/purity_" +
fHMean.fAnaStepNames[
static_cast<int>(step)];
4218 fHMean.fHM.CreateCanvas(cName, cName, 800, 600);
4219 TH2D*
h =
fHMean.H2Clone(hName);
4220 h->GetZaxis()->SetRangeUser(1e-11, 1e-2);
4222 for (
size_t y = 1;
y <= pLabel.size();
y++) {
4223 h->GetYaxis()->SetBinLabel(
y, pLabel[
y - 1].c_str());
4225 double nSigEl =
h->Integral(1,
h->GetXaxis()->GetNbins(), 1, 1);
4226 double nUrqmdEl =
h->Integral(1,
h->GetXaxis()->GetNbins(), 2, 2);
4227 double nGammaEl =
h->Integral(1,
h->GetXaxis()->GetNbins(), 3, 3);
4228 double nEl = nSigEl + nUrqmdEl + nGammaEl;
4229 double nPions =
h->Integral(1,
h->GetXaxis()->GetNbins(), 4, 4);
4230 double nProtons =
h->Integral(1,
h->GetXaxis()->GetNbins(), 5, 5);
4231 double nKaons =
h->Integral(1,
h->GetXaxis()->GetNbins(), 6, 6);
4232 double nOther =
h->Integral(1,
h->GetXaxis()->GetNbins(), 7, 7);
4233 double nAll =
h->Integral(1,
h->GetXaxis()->GetNbins(), 1,
h->GetYaxis()->GetNbins());
4234 double purSig = (nSigEl / nAll) * 100;
4235 double purUrqmdEl = (nUrqmdEl / nAll) * 100;
4236 double purGammaEl = (nGammaEl / nAll) * 100;
4237 double purEl = (nEl / nAll) * 100;
4238 double purPions = (nPions / nAll) * 100;
4239 double purProtons = (nProtons / nAll) * 100;
4240 double purKaons = (nKaons / nAll) * 100;
4241 double purOther = (nOther / nAll) * 100;
4243 0.1, 0.9, 0.45, 0.99);
4244 if (suff ==
"+" || suff ==
"-")
DrawTextOnPad(
"(" + suff +
")", 0.5, 0.9, 0.6, 0.99);
4249 xMin, 0.15, xMax, 0.25);
4252 xMin, 0.27, xMax, 0.35);
4255 xMin, 0.37, xMax, 0.46);
4257 xMin, 0.47, xMax, 0.56);
4260 xMin, 0.58, xMax, 0.67);
4262 xMin, 0.70, xMax, 0.78);
4264 xMin, 0.80, xMax, 0.89);
4272 int nBins =
fHMean.H2(
"hCandPdgVsMom", step)->GetXaxis()->GetNbins();
4273 double xMin =
fHMean.H2(
"hCandPdgVsMom", step)->GetXaxis()->GetXmin();
4274 double xMax =
fHMean.H2(
"hCandPdgVsMom", step)->GetXaxis()->GetXmax();
4275 string hNameEl =
"purityEl_mom_" +
fHMean.fAnaStepNames[
static_cast<int>(step)];
4276 string hNameSig =
"puritySig_mom_" +
fHMean.fAnaStepNames[
static_cast<int>(step)];
4277 TH1D* purityEl =
new TH1D(hNameEl.c_str(),
"purityEl_mom; P [GeV/c]; Purity [%]", nBins, xMin, xMax);
4278 TH1D* puritySig =
new TH1D(hNameSig.c_str(),
"puritySig_mom; P [GeV/c]; Purity [%]", nBins, xMin, xMax);
4279 for (
int i = 1; i <= purityEl->GetNbinsX(); i++) {
4280 double nEl =
fHMean.H2(
"hCandPdgVsMom", step)->GetBinContent(i, 2);
4281 double nSig =
fHMean.H2(
"hCandPdgVsMom", step)->GetBinContent(i, 1);
4283 fHMean.H2(
"hCandPdgVsMom", step)->Integral(i, i, 1,
fHMean.H2(
"hCandPdgVsMom_elid")->GetYaxis()->GetNbins());
4284 double pEl = (nAll != 0) ? 100 * nEl / nAll : 0.;
4285 double pSig = (nAll != 0) ? 100 * nSig / nAll : 0.;
4286 purityEl->SetBinContent(i, pEl);
4287 puritySig->SetBinContent(i, pSig);
4290 purityEl->Scale(1. / 5.);
4291 puritySig->Rebin(5);
4292 puritySig->Scale(1. / 5.);
4293 puritySig->GetYaxis()->SetRangeUser(0., 100.);
4294 string cName =
"Purity/PurityVsMom_" +
fHMean.fAnaStepNames[
static_cast<int>(step)];
4295 fHMean.fHM.CreateCanvas(cName, cName, 800, 800);
4296 DrawH1({puritySig, purityEl}, {
"Signal Electrons",
"All Electrons"},
kLinear,
kLinear,
true, 0.7, 0.8, 0.95, 0.91,
4303 for (
auto step :
fHMean.fAnaSteps) {
4305 string hName =
"hCandPdgInt_" +
fHMean.fAnaStepNames[
static_cast<int>(step)];
4306 TH2D* h2D =
fHMean.H2Clone(
"hCandPdgVsMom_" +
fHMean.fAnaStepNames[
static_cast<int>(step)]);
4308 h2D->ProjectionY(hName.c_str(), 0, h2D->GetXaxis()->GetNbins() + 1,
"");
4309 h1D->GetYaxis()->SetRangeUser(1e-6, 5.);
4310 h1D->GetYaxis()->SetTitle(
"N / Event");
4311 for (
size_t iL = 0; iL < pLabel.size(); iL++) {
4312 h1D->GetXaxis()->SetBinLabel(iL + 1, pLabel[iL].c_str());
4314 string cName =
"Purity/1D/" +
fHMean.fAnaStepNames[
static_cast<int>(step)];
4315 fHMean.fHM.CreateCanvas(cName, cName, 800, 800);
4317 fHMean.DrawAnaStepOnPad(step);
4319 double nAll = h1D->Integral(1, h1D->GetNbinsX(),
"width");
4321 double pElP = 100 * (h1D->GetBinContent(1) / nAll);
4322 double pElU = 100 * (h1D->GetBinContent(2) / nAll);
4323 double pElG = 100 * (h1D->GetBinContent(3) / nAll);
4324 double pPiC = 100 * (h1D->GetBinContent(4) / nAll);
4325 double pProt = 100 * (h1D->GetBinContent(5) / nAll);
4326 double pK = 100 * (h1D->GetBinContent(6) / nAll);
4327 double pOth = 100 * (h1D->GetBinContent(7) / nAll);
4341 vector<string> yLabel = {
"#pi^{+}",
"#pi^{-}",
"p",
"K^{+}",
"K^{-}",
"o."};
4345 for (
const string cat : {
"rec",
"acc",
"chi2prim"}) {
4346 TH2D* rich =
fHMean.H2Clone(
"hPdgVsMom_rich_" + cat);
4347 TH2D* trd =
fHMean.H2Clone(
"hPdgVsMom_trd_" + cat);
4348 TH2D* tof =
fHMean.H2Clone(
"hPdgVsMom_tof_" + cat);
4350 rich->GetZaxis()->SetRangeUser(
min,
max);
4351 trd->GetZaxis()->SetRangeUser(
min,
max);
4352 tof->GetZaxis()->SetRangeUser(
min,
max);
4354 for (
size_t y = 1;
y <= yLabel.size();
y++) {
4355 rich->GetYaxis()->SetBinLabel(
y, yLabel[
y - 1].c_str());
4356 trd->GetYaxis()->SetBinLabel(
y, yLabel[
y - 1].c_str());
4357 tof->GetYaxis()->SetBinLabel(
y, yLabel[
y - 1].c_str());
4360 TCanvas* c =
fHMean.fHM.CreateCanvas(
"Purity/misid_gTracks_" + cat,
"Purity/misid_gTracks_" + cat, 2400, 800);
4559 vector<TH1D*> sHists(
fHMean.fNofSignals);
4561 string sigName =
fHMean.fSignalNames[
static_cast<int>(signal)];
4563 sHist->Scale(1. /
H(signal)->H1(
"hEventNumber")->GetEntries());
4564 sHists[
static_cast<int>(signal)] = sHist;
4576 vector<LmvmDrawMinvData> drawData;
4577 sbg->GetYaxis()->ChangeLabel(0, -1., 10., -1, -1, -1);
4579 drawData.emplace_back(sbg, kBlack, kBlack, 1, -1,
"Cocktail+B_{MC}");
4580 drawData.emplace_back(bg, kGray, kBlack, 1, -1,
"B_{MC}");
4583 drawData.emplace_back(pi0, kGreen - 3, kGreen + 3, 2, -1,
"#pi^{0} #rightarrow #gammae^{+}e^{-}");
4584 drawData.emplace_back(eta, kRed - 4, kRed + 2, 2, -1,
"#eta #rightarrow #gammae^{+}e^{-}");
4585 drawData.emplace_back(sHists[
static_cast<int>(
ELmvmSignal::OmegaD)], kCyan + 2, kCyan + 4, 2, -1,
4586 "#omega #rightarrow #pi^{0}e^{+}e^{-}");
4587 drawData.emplace_back(gamma, -1, kYellow, 4, -1,
"#gamma #rightarrow e^{+}e^{-}");
4588 drawData.emplace_back(sHists[
static_cast<int>(
ELmvmSignal::Omega)], kOrange + 7, kOrange + 4, 2, -1,
4589 "#omega #rightarrow e^{+}e^{-}");
4590 drawData.emplace_back(sHists[
static_cast<int>(
ELmvmSignal::Phi)], kAzure + 2, kAzure + 3, 2, -1,
4591 "#phi #rightarrow e^{+}e^{-}");
4592 drawData.emplace_back(sHists[
static_cast<int>(
ELmvmSignal::Qgp)], kOrange - 2, kOrange - 3, 4, 3112,
"QGP radiation");
4593 drawData.emplace_back(sHists[
static_cast<int>(
ELmvmSignal::Inmed)], kMagenta - 3, kMagenta - 2, 4, 3018,
4595 drawData.emplace_back(cocktail, -1, kRed + 2, 4, -1,
"Cocktail");
4596 if (step >=
ELmvmAnaStep::ElId) drawData.emplace_back(cb, -1, kCyan + 1, 4, -1,
"CB_{calc}");
4599 TLegend* leg =
new TLegend(0.72, 0.55, 0.97, 0.99);
4600 for (
size_t i = 0; i < drawData.size(); i++) {
4601 const auto& d = drawData[i];
4602 d.fH->GetYaxis()->SetLabelSize(0.04);
4604 d.fH->GetYaxis()->SetRangeUser(5e-10, 5e-2);
4606 d.fH->GetYaxis()->SetRangeUser(5e-10, 1e2);
4607 if (d.fFillColor != -1) d.fH->SetFillColor(d.fFillColor);
4608 if (d.fFillStyle != -1) d.fH->SetFillStyle(d.fFillStyle);
4609 leg->AddEntry(d.fH, d.fLegend.c_str(),
"f");
4610 DrawH1(d.fH,
kLinear,
kLinear, (h0 ==
nullptr) ?
"hist" :
"hist,same", d.fLineColor, d.fLineWidth, 0);
4611 if (h0 ==
nullptr) h0 = d.fH;
4614 leg->SetFillColor(kWhite);
4616 gPad->SetLogy(
true);
4803 string cName =
"Background/";
4812 TH1D* rat = (TH1D*) pm2->Clone();
4814 string cNameFull = cName +
"McBg_vs_McPM";
4815 TCanvas* c =
fHMean.fHM.CreateCanvas(cNameFull, cNameFull, 1600, 800);
4818 DrawH1({pm1, pm2}, {
"hMinv_bg + Coc + #gamma",
"hMinvPM"},
kLinear,
kLog,
true, 0.6, 0.8, 0.91, 0.91,
"hist");
4824 for (
const string cat : {
"",
"_trueEl"}) {
4825 string cNameFull = cName +
"Types/Backgrounds" + cat;
4826 TCanvas* c =
fHMean.fHM.CreateCanvas(cNameFull, cNameFull, 2400, 800);
4829 for (
auto step :
fHMean.fAnaStepsFS) {
4830 string stepname =
fHMean.fAnaStepNames[
static_cast<int>(step)];
4835 mc->GetYaxis()->SetRangeUser(1e-7, 1e-2);
4837 DrawH1({mc, cb, re}, {
"B_{MC}",
"CB_{calc}",
"B_{FS}"},
kLinear,
kLog,
true, 0.8, 0.76, 0.91, 0.93,
"pe");
4838 fHMean.DrawAnaStepOnPad(step);
4843 vector<pair<double, double>> varbinZoom = {make_pair(2.5, 0.02)};
4844 for (
const string cat : {
""}) {
4845 string cNameFull = cName +
"Types/Ratios" + cat;
4846 TCanvas* c =
fHMean.fHM.CreateCanvas(cNameFull, cNameFull, 3200, 2400);
4849 for (
auto step :
fHMean.fAnaStepsFS) {
4850 string stepname =
fHMean.fAnaStepNames[
static_cast<int>(step)];
4860 mcZoom->GetXaxis()->SetRangeUser(0., 0.5);
4861 cbZoom->GetXaxis()->SetRangeUser(0., 0.5);
4862 reZoom->GetXaxis()->SetRangeUser(0., 0.5);
4863 TH1D* remc = (TH1D*) re->Clone();
4864 TH1D* cbmc = (TH1D*) cb->Clone();
4865 TH1D* cbre = (TH1D*) cb->Clone();
4887 remc->GetYaxis()->SetTitle(
"Ratio");
4888 cbre->GetYaxis()->SetTitle(
"Ratio");
4889 remc->GetYaxis()->SetRangeUser(0.5, 2.);
4890 cbmc->GetYaxis()->SetRangeUser(0.5, 2.);
4891 cbre->GetYaxis()->SetRangeUser(0.8, 1.2);
4892 mc->GetYaxis()->SetRangeUser(1e-7, 1e-2);
4895 DrawH1({mc, re, cb}, {
"B_{MC}",
"B_{FS}",
"CB_{calc}"},
kLinear,
kLog,
true, 0.78, 0.76, 0.95, 0.93,
"pe");
4896 fHMean.DrawAnaStepOnPad(step);
4899 DrawH1({reZoom, cbZoom}, {
"B_{FS}",
"CB_{calc}"},
kLinear,
kLog,
true, 0.78, 0.82, 0.95, 0.93,
"pe");
4900 fHMean.DrawAnaStepOnPad(step);
4902 DrawH1({cbmc, remc}, {
"CB_{calc} / B_{MC}",
"B_{FS} / B_{MC}"},
kLinear,
kLinear,
true, 0.7, 0.75, 0.95, 0.9,
4904 fHMean.DrawAnaStepOnPad(step);
4907 fHMean.DrawAnaStepOnPad(step);
4919 TH1D* rat = (TH1D*) cbre->Clone();
4921 ksim->GetYaxis()->SetTitle(
"#it{k} Factor");
4922 rat->GetYaxis()->SetTitle(
"Ratio");
4923 string cNameFull = cName +
"Types/CombinatorialBackground";
4924 TCanvas* c =
fHMean.fHM.CreateCanvas(cNameFull, cNameFull, 2400, 800);
4927 DrawH1({cbsim, cbre}, {
"Simulation",
"Fast Simulations"},
kLinear,
kLog,
true, 0.65, 0.8, 0.92, 0.92,
"hist");
4928 DrawTextOnPad(
"Combinatorial Background", 0.2, 0.9, 0.8, 0.99);
4930 DrawH1({ksim, kre}, {
"Simulation",
"Fast Simulations"},
kLinear,
kLinear,
true, 0.65, 0.8, 0.92, 0.92,
"hist");
4933 DrawH1({rat}, {
"CB_{calc, RE} / CB_{calc, Sim}"},
kLinear,
kLinear,
true, 0.65, 0.82, 0.92, 0.92,
"hist");
4934 DrawTextOnPad(
"Ratio Combinatorial Backgrounds", 0.2, 0.9, 0.85, 0.99);
4940 string stepname =
fHMean.fAnaStepNames[
static_cast<int>(step)];
4941 vector<pair<double, double>> varbin = {make_pair(1.6, 0.2), make_pair(2.0, 0.4), make_pair(2.5, 0.5)};
4942 for (
const string ev : {
"sameEv",
"mixedEv"}) {
4943 string cNameFull = cName +
"Types/PairYields_" + ev;
4944 TCanvas* c =
fHMean.fHM.CreateCanvas(cNameFull, cNameFull, 2400, 2400);
4947 for (
const string comb : {
"PM",
"PP",
"MM"}) {
4951 TH1D* remc = (TH1D*) re->Clone();
4952 TH1D* cbmc = (TH1D*) cb->Clone();
4953 TH1D* cbre = (TH1D*) cb->Clone();
4957 remc->GetYaxis()->SetTitle(
"Ratio");
4958 cbre->GetYaxis()->SetTitle(
"Ratio");
4959 mc->GetYaxis()->SetRangeUser(1e-7, 1e-2);
4960 remc->GetYaxis()->SetRangeUser(0.0, 2.);
4961 string combMode = (comb ==
"PP") ?
"e^{+}e^{+}" : (comb ==
"MM") ?
"e^{-}e^{-}" :
"e^{+}e^{-}";
4962 string label = combMode +
" (" + ev +
")";
4964 DrawH1({mc, re, cb}, {
"MC Simulation",
"Fast Simulations",
"Comb. BG Procedure"},
kLinear,
kLog,
true, 0.6, 0.7,
4968 DrawH1({remc, cbmc}, {
"Fast Simulations / MC Simulation",
"Comb. BG Procedure / MC Simulation"},
kLinear,
4969 kLinear,
true, 0.55, 0.75, 0.92, 0.92,
"hist");
4972 DrawH1({cbre}, {
"Comb. BG Procedure / Fast Sim."},
kLinear,
kLinear,
true, 0.5, 0.82, 0.95, 0.9,
"hist");
4978 TH1D* rat = (TH1D*) gmre->Clone();
4980 rat->GetYaxis()->SetTitle(
"Ratio");
4981 string cNameFull2 = cName +
"Types/GeomMean_" + ev;
4982 TCanvas* c2 =
fHMean.fHM.CreateCanvas(cNameFull2, cNameFull2, 1600, 800);
4985 DrawH1({gmre, gmcb}, {
"Fast Simulations",
"Comb. BG Procedure"},
kLinear,
kLog,
true, 0.65, 0.75, 0.95, 0.9,
4987 DrawTextOnPad(
"Geometric Mean (" + ev +
")", 0.25, 0.9, 0.75, 0.99);
4989 DrawH1({rat}, {
"Fast Simulations / Comb. BG Procedure"},
kLinear,
kLinear,
true, 0.65, 0.82, 0.95, 0.9,
"hist");
4997 string stepname =
fHMean.fAnaStepNames[
static_cast<int>(step)];
4998 vector<pair<double, double>> varbin = {make_pair(1.6, 0.2), make_pair(2.0, 0.4), make_pair(2.5, 0.5)};
5003 TH1D* rat1 = (TH1D*) re->Clone();
5004 TH1D* rat2 = (TH1D*) re->Clone();
5008 fHMean.fHM.CreateCanvas(cName +
"Types/PM_wwo-Cocktail", cName +
"Types/PM_wwo-Cocktail", 800, 800);
5009 DrawH1({rat1, rat2}, {
"B_{FS}/B_{MC}",
"(B_{FS} + Coc)/B_{MC}"},
kLinear,
kLog,
true, 0.65, 0.75, 0.95, 0.9,
5017 TH1D* rat = (TH1D*) bg1->Clone();
5019 rat->GetYaxis()->SetTitle(
"Ratio");
5020 TCanvas* c =
fHMean.fHM.CreateCanvas(cName +
"Types/FastSimBgs", cName +
"Types/FastSimBgs", 1600, 800);
5023 DrawH1({bg1, bg2}, {
"hMinv_bg_elid",
"hMinvPM_sameEv"},
kLinear,
kLog,
true, 0.6, 0.7, 0.95, 0.9,
"hist");
5028 for (
auto step :
fHMean.fAnaSteps) {
5030 string stepname =
fHMean.fAnaStepNames[
static_cast<int>(step)];
5032 vector<pair<double, double>> varbin =
fBwVarBg2;
5036 string hName =
"hMinvBgSource2_" +
fHMean.fAnaStepNames[
static_cast<int>(step)] +
"_";
5038 TH1D* pipi =
VaryBinWidth(
"mean", hName +
"pipi", varbin);
5039 TH1D* pi0pi0 =
VaryBinWidth(
"mean", hName +
"pi0pi0", varbin);
5041 TH1D* gpi =
VaryBinWidth(
"mean", hName +
"gpi", varbin);
5042 TH1D* gpi0 =
VaryBinWidth(
"mean", hName +
"gpi0", varbin);
5044 TH1D* pipi0 =
VaryBinWidth(
"mean", hName +
"pipi0", varbin);
5045 TH1D* pio =
VaryBinWidth(
"mean", hName +
"pio", varbin);
5046 TH1D* pi0o =
VaryBinWidth(
"mean", hName +
"pi0o", varbin);
5048 vector<LmvmDrawMinvData> drawData;
5049 drawData.emplace_back(gpi0, kBlack, kBlack, 1, -1,
"#gamma - #pi^{0}");
5050 drawData.emplace_back(pi0pi0, kGray + 1, kGray + 1, 1, -1,
"#pi^{0} - #pi^{0}");
5051 drawData.emplace_back(gg, kCyan + 2, kCyan + 4, 2, -1,
"#gamma - #gamma");
5052 drawData.emplace_back(pipi0, kGreen - 3, kGreen + 3, 2, -1,
"#pi^{#pm}_{mis} - #pi^{0}");
5053 drawData.emplace_back(pi0o, kRed - 4, kRed + 2, 1, -1,
"#pi^{0} - o.");
5054 drawData.emplace_back(gpi, kOrange + 7, kOrange + 4, 2, -1,
"#gamma - #pi^{#pm}_{mis}");
5055 drawData.emplace_back(go, kAzure + 2, kAzure + 3, 2, -1,
"#gamma - o.");
5056 drawData.emplace_back(pipi, kOrange - 2, kOrange - 3, 4, 3112,
"#pi^{#pm}_{mis} - #pi^{#pm}_{mis}");
5057 drawData.emplace_back(pio, kMagenta - 3, kMagenta - 2, 4, 3018,
"#pi^{#pm}_{mis} - o.");
5058 drawData.emplace_back(oo, -1, kRed + 2, 4, -1,
"o. - o.");
5060 string cNameFull = cName +
"Composition/comp_" +
fHMean.fAnaStepNames[
static_cast<int>(step)];
5061 fHMean.fHM.CreateCanvas(cNameFull, cNameFull, 1000, 1000);
5064 TLegend* leg =
new TLegend(0.75, 0.35, 0.92, 0.9);
5065 for (
size_t i = 0; i < drawData.size(); i++) {
5066 const auto& d = drawData[i];
5067 d.fH->GetYaxis()->SetLabelSize(0.05);
5068 d.fH->GetYaxis()->SetRangeUser(2e-7, 3e-3);
5069 if (d.fFillColor != -1) d.fH->SetFillColor(d.fFillColor);
5070 if (d.fFillStyle != -1) d.fH->SetFillStyle(d.fFillStyle);
5071 leg->AddEntry(d.fH, d.fLegend.c_str(),
"f");
5072 DrawH1(d.fH,
kLinear,
kLinear, (h0 ==
nullptr) ?
"hist" :
"hist,same", d.fLineColor, d.fLineWidth, 0);
5073 if (h0 ==
nullptr) h0 = d.fH;
5076 leg->SetFillColor(kWhite);
5078 fHMean.DrawSimDataLabel(.57, 1.5e-3, 0.03);
5079 gPad->SetLogy(
true);
5084 string hName =
"hMinvBgSource2_" +
fHMean.fAnaStepNames[
static_cast<int>(step)] +
"_";
5085 string cNameFull = cName +
"Composition/sum_" +
fHMean.fAnaStepNames[
static_cast<int>(step)];
5087 TH1D* physBg =
fHMean.H1Clone(hName +
"gg");
5088 physBg->Add(
fHMean.H1(hName +
"gpi0"));
5089 physBg->Add(
fHMean.H1(hName +
"pi0pi0"));
5091 TH1D* cPi =
fHMean.H1Clone(hName +
"pipi");
5092 cPi->Add(
fHMean.H1Clone(hName +
"gpi"));
5093 cPi->Add(
fHMean.H1Clone(hName +
"pipi0"));
5094 cPi->Add(
fHMean.H1Clone(hName +
"pio"));
5096 TH1D* rest =
fHMean.H1Clone(hName +
"oo");
5097 rest->Add(
fHMean.H1Clone(hName +
"go"));
5098 rest->Add(
fHMean.H1Clone(hName +
"pi0o"));
5103 physBg2->GetYaxis()->SetRangeUser(1e-7, 1e-2);
5105 fHMean.fHM.CreateCanvas(cNameFull, cNameFull, 1000, 1000);
5106 DrawH1({physBg2, cPi2, rest2}, {
"Phys. BG",
"BG w. misid. #pi^{#pm}",
"Rest"},
kLinear,
kLog,
true, 0.7, 0.8,
5108 fHMean.DrawSimDataLabel(.15, 2.7e-7, 0.03);
5112 TH1D* bgRat =
static_cast<TH1D*
>(physBg2->Clone());
5120 int nBins = bgRat->GetNbinsX();
5121 for (
int i = 1; i <= nBins; i++) {
5122 if (bgRat->GetBinContent(i) == 0)
5123 bgRat->SetBinContent(i, 2);
5125 double r = bgRat->GetBinContent(i) / bg2->GetBinContent(i);
5126 bgRat->SetBinContent(i, r);
5130 bgRat->SetMinimum(0.95);
5131 bgRat->SetMaximum(1.05);
5133 string cNameFull2 = cName +
"Composition/Ratio_" +
fHMean.fAnaStepNames[
static_cast<int>(step)];
5134 fHMean.fHM.CreateCanvas(cNameFull2, cNameFull2, 800, 800);
5181 TLine* lineH =
new TLine(0., 4., 4., 4.);
5182 TLine* lineV =
new TLine(4., 0., 4., 4.);
5183 lineH->SetLineWidth(6.);
5184 lineV->SetLineWidth(6.);
5185 for (
const string type : {
"pRec"}) {
5186 vector<double> v_pm, v_pp, v_mm, v_pm_el, v_pp_el, v_mm_el;
5187 for (
const string comb : {
"PM",
"PP",
"MM"}) {
5188 for (
int iMinv = 0; iMinv <= 24; iMinv++) {
5192 TH2D*
h =
fHMean.H2Clone(hName);
5193 h->GetZaxis()->SetRangeUser(1e-15, 1e-3);
5194 vector<string> label = {
"e_{signal}",
"e_{#pi^{0}}",
"e_{#gamma}",
"e_{o.}",
5195 "#pi^{#pm}_{mis}",
"proton",
"kaon",
"other"};
5196 for (
size_t y = 1;
y <= label.size();
y++) {
5197 h->GetXaxis()->SetBinLabel(
y, label[
y - 1].c_str());
5198 h->GetYaxis()->SetBinLabel(
y, label[
y - 1].c_str());
5200 double I =
h->Integral(1,
h->GetNbinsX(), 1,
h->GetNbinsX());
5201 double ITrue =
h->Integral(1, 4, 1, 4);
5204 v_pm_el.push_back(ITrue);
5206 else if (comb ==
"PP") {
5208 v_pp_el.push_back(ITrue);
5210 else if (comb ==
"MM") {
5212 v_mm_el.push_back(ITrue);
5214 string xTitle =
"PDG +";
5215 string yTitle =
"PDG -";
5216 if (comb ==
"PP") yTitle =
"PDG +";
5217 if (comb ==
"MM") xTitle =
"PDG -";
5218 h->GetXaxis()->SetTitle(xTitle.c_str());
5219 h->GetYaxis()->SetTitle(yTitle.c_str());
5220 string cNameFull = cName +
"PairPdg/" + type +
"/" + comb +
"/" +
Cbm::NumberToString(iMinv, 0);
5221 fHMean.fHM.CreateCanvas(cNameFull, cNameFull, 900, 600);
5230 LOG(warn) <<
"PairYieldSymmetry: Check Calculation of Symmetry Ratios!!!";
5232 TH1D* hSymm_ppmm =
new TH1D((
"hSymm_ppmm_" + type).c_str(),
"hSymm_ppmm; M_{ee} [GeV/c^{2}]; Ratio", 25, 0., 2.5);
5233 TH1D* hSymm_pppm =
new TH1D((
"hSymm_pppm_" + type).c_str(),
"hSymm_pppm; M_{ee} [GeV/c^{2}]; Ratio", 25, 0., 2.5);
5234 TH1D* hSymm_mmpm =
new TH1D((
"hSymm_mmpm_" + type).c_str(),
"hSymm_mmpm; M_{ee} [GeV/c^{2}]; Ratio", 25, 0., 2.5);
5235 TH1D* hSymm_ppmm_el =
5236 new TH1D((
"hSymm_ppmm_" + type +
"_el").c_str(),
"hSymm_ppmm_el; M_{ee} [GeV/c^{2}]; Ratio", 25, 0., 2.5);
5237 TH1D* hSymm_pppm_el =
5238 new TH1D((
"hSymm_pppm_" + type +
"_el").c_str(),
"hSymm_pppm_el; M_{ee} [GeV/c^{2}]; Ratio", 25, 0., 2.5);
5239 TH1D* hSymm_mmpm_el =
5240 new TH1D((
"hSymm_mmpm_" + type +
"_el").c_str(),
"hSymm_mmpm_el; M_{ee} [GeV/c^{2}]; Ratio", 25, 0., 2.5);
5241 for (
int iB = 1; iB <= hSymm_ppmm->GetNbinsX(); iB++) {
5242 double ppmm = v_pp[iB - 1] / v_mm[iB - 1];
5243 double pppm = v_pp[iB - 1] / v_pm[iB - 1];
5244 double mmpm = v_mm[iB - 1] / v_pm[iB - 1];
5245 double ppmmEl = v_pp_el[iB - 1] / v_mm_el[iB - 1];
5246 double pppmEl = v_pp_el[iB - 1] / v_pm_el[iB - 1];
5247 double mmpmEl = v_mm_el[iB - 1] / v_pm_el[iB - 1];
5249 hSymm_ppmm->SetBinContent(iB, ppmm);
5250 hSymm_pppm->SetBinContent(iB, pppm);
5251 hSymm_mmpm->SetBinContent(iB, mmpm);
5252 hSymm_ppmm_el->SetBinContent(iB, ppmmEl);
5253 hSymm_pppm_el->SetBinContent(iB, pppmEl);
5254 hSymm_mmpm_el->SetBinContent(iB, mmpmEl);
5256 vector<pair<double, double>> varbin =
fBwVarBg5;
5260 TH1D* hSymm_ppmm_el2 =
VaryBinWidth(hSymm_ppmm_el, varbin);
5261 TH1D* hSymm_pppm_el2 =
VaryBinWidth(hSymm_pppm_el, varbin);
5262 TH1D* hSymm_mmpm_el2 =
VaryBinWidth(hSymm_mmpm_el, varbin);
5263 hSymm_ppmm2->GetYaxis()->SetRangeUser(0., 1.2);
5264 hSymm_pppm2->GetYaxis()->SetRangeUser(0., 1.2);
5265 hSymm_mmpm2->GetYaxis()->SetRangeUser(0., 1.2);
5266 hSymm_ppmm_el2->GetYaxis()->SetRangeUser(0., 1.2);
5267 hSymm_pppm_el2->GetYaxis()->SetRangeUser(0., 1.2);
5268 hSymm_mmpm_el2->GetYaxis()->SetRangeUser(0., 1.2);
5269 hSymm_ppmm->GetYaxis()->SetRangeUser(0., 1.2);
5270 hSymm_pppm->GetYaxis()->SetRangeUser(0., 1.2);
5271 hSymm_mmpm->GetYaxis()->SetRangeUser(0., 1.2);
5272 hSymm_ppmm_el->GetYaxis()->SetRangeUser(0., 1.2);
5273 hSymm_pppm_el->GetYaxis()->SetRangeUser(0., 1.2);
5274 hSymm_mmpm_el->GetYaxis()->SetRangeUser(0., 1.2);
5275 string cNameSymm = cName +
"PairPdg/ChargeSymmetries_" + type;
5276 TCanvas* c =
fHMean.fHM.CreateCanvas(cNameSymm, cNameSymm, 1600, 800);
5279 DrawH1({hSymm_ppmm2, hSymm_pppm2, hSymm_mmpm2}, {
"++/--",
"++/+-",
"--/+-"},
kLinear,
kLinear,
true, 0.72, 0.7,
5280 0.91, 0.91,
"hist");
5281 DrawTextOnPad(
"Charge Symmetries (all)", 0.1, 0.9, 0.8, 0.99);
5283 DrawH1({hSymm_ppmm_el2, hSymm_pppm_el2, hSymm_mmpm_el2}, {
"++/--",
"++/+-",
"--/+-"},
kLinear,
kLinear,
true,
5284 0.72, 0.7, 0.91, 0.91,
"hist");
5285 DrawTextOnPad(
"Charge Symmetries (electrons)", 0.05, 0.9, 0.85, 0.99);
5287 TCanvas* ct =
fHMean.fHM.CreateCanvas(cNameSymm +
"_test", cNameSymm +
"_test", 1600, 800);
5290 DrawH1({hSymm_ppmm, hSymm_pppm, hSymm_mmpm}, {
"++/--",
"++/+-",
"--/+-"},
kLinear,
kLinear,
true, 0.72, 0.7, 0.91,
5292 DrawTextOnPad(
"Charge Symmetries (all)", 0.1, 0.9, 0.8, 0.99);
5294 DrawH1({hSymm_ppmm_el, hSymm_pppm_el, hSymm_mmpm_el}, {
"++/--",
"++/+-",
"--/+-"},
kLinear,
kLinear,
true, 0.72,
5295 0.7, 0.91, 0.91,
"hist");
5296 DrawTextOnPad(
"Charge Symmetries (electrons)", 0.05, 0.9, 0.85, 0.99);
5306 TH1D* rat1 = (TH1D*) bgU->Clone();
5307 TH1D* rat2 = (TH1D*) bgUEl->Clone();
5310 rat1->GetYaxis()->SetTitle(
"Ratio");
5311 rat2->GetYaxis()->SetTitle(
"Ratio");
5312 bg->GetYaxis()->SetRangeUser(2e-7, 1e-2);
5313 rat1->GetYaxis()->SetRangeUser(0., 1.2);
5314 TCanvas* c =
fHMean.fHM.CreateCanvas(
"Background/Categories",
"Background/Categories", 1600, 800);
5317 DrawH1({bg, bgU, bgUEl}, {
"B_{MC} (all)",
"B_{MC} (UrQMD)",
"B_{MC} (UrQMD electrons)"},
kLinear,
kLog,
true, 0.65,
5318 0.75, 0.91, 0.9,
"p");
5320 DrawH1({rat1, rat2}, {
"B_{MC}(UrQMD) / B_{MC}(all)",
"B_{MC}(UrQMD-El) / B_{MC}(all)"},
kLinear,
kLinear,
true,
5321 0.65, 0.75, 0.91, 0.9,
"p");
5461 string folder =
"CombinatorialBackground/";
5462 TH1D* time =
fHEvMix.H1Clone(
"hmx_time");
5463 fHMean.fHM.CreateCanvas(
"CombinatorialBackground/time",
"CombinatorialBackground/time", 800, 800);
5468 vector<int> canvas = {1, 4, 2, 5, 3, 6};
5469 for (
const string evMode : {
"sameEv",
"mixedEv"}) {
5470 string cName = folder +
"PairYields/CheckWithSim_" + evMode;
5471 TCanvas* c =
fHMean.fHM.CreateCanvas(cName, cName, 2400, 1600);
5474 for (
const string comb : {
"PM",
"PP",
"MM"}) {
5477 TH1D* r = (TH1D*) cb->Clone();
5479 r->GetYaxis()->SetTitle(
"Ratio");
5480 r->GetYaxis()->SetRangeUser(.95, 1.05);
5483 DrawH1({sim, cb}, {
"Simulation",
"CB Procedure"},
kLinear,
kLog,
true, 0.6, 0.7, 0.92, 0.9,
"hist");
5496 for (
const string cat : {
"",
"_trueEl"}) {
5497 for (
const string ev : {
"sameEv",
"mixedEv"}) {
5498 string cName = folder +
"PairYields/AllCombs_" + ev + cat;
5499 TCanvas* c =
fHMean.fHM.CreateCanvas(cName.c_str(), cName.c_str(), 1600, 800);
5503 vector<pair<double, double>> varbin =
fBwVarBg4;
5505 TH1D* pm =
VaryBinWidth(
"evmix",
fHMean.GetName(
"hmx_MinvPM" + cat +
"_" + ev, step), nEvents, varbin);
5506 TH1D* pp =
VaryBinWidth(
"evmix",
fHMean.GetName(
"hmx_MinvPP" + cat +
"_" + ev, step), nEvents, varbin);
5507 TH1D* mm =
VaryBinWidth(
"evmix",
fHMean.GetName(
"hmx_MinvMM" + cat +
"_" + ev, step), nEvents, varbin);
5508 TH1D* ppmm = (TH1D*) pp->Clone();
5510 TH1D* pppm = (TH1D*) pp->Clone();
5511 TH1D* mmpm = (TH1D*) mm->Clone();
5514 ppmm->GetYaxis()->SetTitle(
"Ratio");
5515 ppmm->GetYaxis()->SetRangeUser(0., 1.2);
5517 DrawH1({pm, pp, mm}, {
"e+e-",
"e+e+",
"e-e-"},
kLinear,
kLog,
true, 0.8, 0.7, 0.95, 0.9,
"pe");
5518 fHMean.DrawAnaStepOnPad(step);
5522 DrawH1({ppmm, pppm, mmpm}, {
"++/--",
"++/+-",
"--/+-"},
kLinear,
kLinear,
true, 0.77, 0.77, 0.95, 0.95,
"pe");
5523 fHMean.DrawAnaStepOnPad(step);
5531 string stepstr =
fHMean.fAnaStepNames[
static_cast<int>(step)];
5532 for (
const string cat : {
"",
"_trueEl"}) {
5535 TH1D* ppMix =
VaryBinWidth(
"evmix",
"hmx_MinvPP" + cat +
"_mixedEv_" + stepstr,
5537 TH1D* mmMix =
VaryBinWidth(
"evmix",
"hmx_MinvMM" + cat +
"_mixedEv_" + stepstr,
5542 double nP = ppSame->Integral(minBin, maxBin) / ppMix->Integral(minBin, maxBin);
5543 double nM = mmSame->Integral(minBin, maxBin) / mmMix->Integral(minBin, maxBin);
5548 ppSame->GetYaxis()->SetRangeUser(1e-8, 5e-3);
5549 TH1D* ppSame2 = (TH1D*) ppSame->Clone();
5550 TH1D* mmSame2 = (TH1D*) mmSame->Clone();
5551 TH1D* ppMix2 = (TH1D*) ppMix->Clone();
5552 TH1D* mmMix2 = (TH1D*) mmMix->Clone();
5553 ppSame2->GetXaxis()->SetRangeUser(0.3, 1.2);
5554 mmSame2->GetXaxis()->SetRangeUser(0.3, 1.2);
5555 ppMix2->GetXaxis()->SetRangeUser(0.3, 1.2);
5556 mmMix2->GetXaxis()->SetRangeUser(0.3, 1.2);
5557 ppSame2->GetYaxis()->SetRangeUser(1e-6, 5e-3);
5558 TCanvas* c =
fHMean.fHM.CreateCanvas(folder +
"PairYields/Likesign/likesign" + cat,
5559 folder +
"PairYields/Likesign/likesign" + cat, 1600, 800);
5562 DrawH1({ppSame, ppMix, mmSame, mmMix},
5563 {
"N_{e^{+}e^{+}}^{same}",
"N_{e^{+}e^{+}}^{mixed} #upoint n_{+}",
"N_{e^{-}e^{-}}^{same} #upoint 0.1",
5564 "N_{e^{-}e^{-}}^{mixed} #upoint n_{-} #upoint 0.1"},
5565 kLinear,
kLog,
true, 0.7, 0.75, 0.95, 0.99,
"pe");
5566 fHMean.DrawAnaStepOnPad(step);
5568 DrawH1({ppSame2, ppMix2, mmSame2, mmMix2},
5569 {
"N_{e^{+}e^{+}}^{same}",
"N_{e^{+}e^{+}}^{mixed} #upoint n_{+}",
"N_{e^{-}e^{-}}^{same} #upoint 0.1",
5570 "N_{e^{-}e^{-}}^{mixed} #upoint n_{-} #upoint 0.1"},
5571 kLinear,
kLog,
true, 0.7, 0.75, 0.95, 0.99,
"pe");
5572 fHMean.DrawAnaStepOnPad(step);
5575 0.06, 0.15, 0.55, 0.32);
5583 for (
const string cat : {
"",
"_trueEl"}) {
5584 TCanvas* c =
fHMean.fHM.CreateCanvas(folder +
"PairYields/Likesign/Ratio_" + cat,
5585 folder +
"PairYields/Likesign/Ratio_" + cat, 1800, 1800);
5588 for (
auto step :
fHMean.fAnaSteps) {
5592 TH1D* r = (TH1D*) mm->Clone();
5594 r->GetYaxis()->SetRangeUser(0., 5.);
5595 r->GetYaxis()->SetTitle(
"Ratio e^{-}e^{-}/e^{+}e^{+}");
5598 fHMean.DrawAnaStepOnPad(step);
5605 for (
const string cat : {
"",
"_trueEl"}) {
5606 fHMean.fHM.CreateCanvas(folder +
"GeometricMean/GeomMean" + cat, folder +
"GeometricMean/GeomMean" + cat, 800,
5611 for (
auto step :
fHMean.fAnaSteps) {
5618 double scale = same->Integral(minBin, maxBin) / mixed->Integral(minBin, maxBin);
5619 mixed->Scale(scale);
5620 same->GetYaxis()->SetRangeUser(1e-7, 5e-3);
5621 same->GetXaxis()->SetTitle(
"M_{ee} [GeV/c^{2}]");
5622 DrawH1({same, mixed}, {
"geom. mean (same)",
"geom. mean (mixed)"},
kLinear,
kLog,
true, 0.6, 0.8, 0.99, 0.95,
5624 fHMean.DrawAnaStepOnPad(step);
5632 for (
const string cat : {
"",
"_trueEl"}) {
5634 fHMean.fHM.CreateCanvas(folder +
"/k-Factor/steps" + cat, folder +
"/k-Factor/steps" + cat, 2400, 2400);
5641 k->GetYaxis()->SetRangeUser(0.6, 1.4);
5642 k->GetYaxis()->SetTitle(
"#it{k} Factor");
5645 fHMean.DrawAnaStepOnPad(step);
5658 TH1D* g2 = (TH1D*) g->Clone();
5660 TH1D* k2 = (TH1D*) pm->Clone();
5661 for (
int iB = 1; iB <= pm->GetNbinsX(); iB++) {
5662 double con = pm->GetBinContent(iB) / (2 * std::sqrt(pp->GetBinContent(iB) * mm->GetBinContent(iB)));
5663 k2->SetBinContent(iB, con);
5665 k->GetYaxis()->SetTitle(
"#it{k} Factor");
5666 k->GetYaxis()->SetRangeUser(0.6, 1.2);
5667 TCanvas* c =
fHMean.fHM.CreateCanvas(folder +
"/k-Factor/components", folder +
"/k-Factor/components", 2400, 800);
5670 DrawH1({pm, pp, mm, g, g2}, {
"#LT b^{#pm} #GT",
"#LT b^{++} #GT",
"#LT b^{--} #GT",
"Geom. Mean",
"2 * Geom. Mean"},
5671 kLinear,
kLog,
true, 0.75, 0.65, 0.99, 0.95,
"p");
5672 fHMean.DrawAnaStepOnPad(step);
5675 DrawH1({k, k2}, {
"k",
"k2"},
kLinear,
kLinear,
true, 0.8, 0.8, 0.9, 0.9,
"hist");
5676 fHMean.DrawAnaStepOnPad(step);
5679 TH1D* gmean2 = (TH1D*) pm->Clone();
5680 for (
int iB = 1; iB <= gmean2->GetNbinsX(); iB++) {
5681 double con = std::sqrt(pp->GetBinContent(iB) * mm->GetBinContent(iB));
5682 gmean2->SetBinContent(iB, con);
5684 DrawH1({gmean1, gmean2}, {
"standard",
"calc. new"},
kLinear,
kLog,
true, 0.75, 0.65, 0.99, 0.95,
"p");
5685 DrawTextOnPad(
"Geometric Mean (mixed ev.)", 0.2, 0.9, 0.8, 0.99);
5691 TCanvas* c =
fHMean.fHM.CreateCanvas(folder +
"/k-Factor/categories", folder +
"/k-Factor/categories", 1600, 800);
5694 for (
const string cat : {
"_raw",
"_fit"}) {
5699 kAll->GetYaxis()->SetTitle(
"#it{k} Factor");
5700 kAll->GetYaxis()->SetRangeUser(0.6, 1.2);
5705 fHMean.DrawAnaStepOnPad(step);
5712 for (
const string cat : {
"",
"_mother"}) {
5713 TH2D* pxPM =
fHMean.H2Clone(
"hMinvCombPM_px" + cat);
5714 TH2D* pxPP =
fHMean.H2Clone(
"hMinvCombPP_px" + cat);
5715 TH2D* pxMM =
fHMean.H2Clone(
"hMinvCombMM_px" + cat);
5716 TH2D* pyPM =
fHMean.H2Clone(
"hMinvCombPM_py" + cat);
5717 TH2D* pyPP =
fHMean.H2Clone(
"hMinvCombPP_py" + cat);
5718 TH2D* pyMM =
fHMean.H2Clone(
"hMinvCombMM_py" + cat);
5719 TH2D* kPx = (TH2D*) pxPM->Clone();
5720 TH2D* kPy = (TH2D*) pyPM->Clone();
5721 for (
int iX = 1; iX <= pxPM->GetXaxis()->GetNbins(); iX++) {
5722 for (
int iY = 1; iY <= pxPM->GetYaxis()->GetNbins(); iY++) {
5724 pxPM->GetBinContent(iX, iY) / (2 * std::sqrt(pxPP->GetBinContent(iX, iY) * pxMM->GetBinContent(iX, iY)));
5726 pyPM->GetBinContent(iX, iY) / (2 * std::sqrt(pyPP->GetBinContent(iX, iY) * pyMM->GetBinContent(iX, iY)));
5727 kPx->SetBinContent(iX, iY, kx);
5728 kPy->SetBinContent(iX, iY, ky);
5731 kPx->GetZaxis()->SetRangeUser(0., 2.);
5732 kPy->GetZaxis()->SetRangeUser(0., 2.);
5733 kPx->GetXaxis()->SetTitle(
"P_{X_{1}} [GeV/c]");
5734 kPx->GetYaxis()->SetTitle(
"P_{X_{2}} [GeV/c]");
5735 kPy->GetXaxis()->SetTitle(
"P_{Y_{1}} [GeV/c]");
5736 kPy->GetYaxis()->SetTitle(
"P_{Y_{2}} [GeV/c]");
5737 kPx->GetZaxis()->SetTitle(
"#it{k} Factor");
5738 kPy->GetZaxis()->SetTitle(
"#it{k} Factor");
5740 fHMean.fHM.CreateCanvas(folder +
"k-Factor/dep_Pt" + cat, folder +
"k-Factor/dep_Pt" + cat, 1600, 800);
5748 TH2D* RPM =
fHMean.H2Clone(
"hMinvCombPM_R" + cat);
5749 TH2D* RPP =
fHMean.H2Clone(
"hMinvCombPP_R" + cat);
5750 TH2D* RMM =
fHMean.H2Clone(
"hMinvCombMM_R" + cat);
5751 TH2D* kR = (TH2D*) RPM->Clone();
5752 for (
int iX = 1; iX <= RPM->GetXaxis()->GetNbins(); iX++) {
5753 for (
int iY = 1; iY <= RPM->GetYaxis()->GetNbins(); iY++) {
5755 RPM->GetBinContent(iX, iY) / (2 * std::sqrt(RPP->GetBinContent(iX, iY) * RMM->GetBinContent(iX, iY)));
5756 kR->SetBinContent(iX, iY, k);
5759 kR->GetZaxis()->SetRangeUser(0., 2.);
5760 kR->GetZaxis()->SetTitle(
"#it{k} Factor");
5761 fHMean.fHM.CreateCanvas(folder +
"k-Factor/dep_Vertex" + cat, folder +
"k-Factor/dep_Vertex" + cat, 800, 800);
5766 TH2D* aPM =
fHMean.H2Clone(
"hMinvCombPM_theta-minv");
5767 TH2D* aPP =
fHMean.H2Clone(
"hMinvCombPP_theta-minv");
5768 TH2D* aMM =
fHMean.H2Clone(
"hMinvCombMM_theta-minv");
5769 TH2D* kAngle = (TH2D*) aPM->Clone();
5770 for (
int iA = 1; iA <= aPM->GetXaxis()->GetNbins(); iA++) {
5771 for (
int iM = 1; iM <= aPM->GetYaxis()->GetNbins(); iM++) {
5772 double den = (2 * std::sqrt(aPP->GetBinContent(iA, iM) * aMM->GetBinContent(iA, iM)));
5773 double k = (den == 0) ? 0 : aPM->GetBinContent(iA, iM) / den;
5774 kAngle->SetBinContent(iA, iM, k);
5777 kAngle->GetZaxis()->SetRangeUser(0., 2.);
5778 kAngle->GetZaxis()->SetTitle(
"#it{k} Factor");
5779 fHMean.fHM.CreateCanvas(folder +
"k-Factor/dep_theta-minv", folder +
"k-Factor/dep_theta-minv", 800, 800);
5781 kAngle->SetMarkerSize(2.);
5786 for (
const string cat : {
"",
"_trueEl",
"_urqmdAll",
"_urqmdEl"}) {
5788 fHMean.fHM.CreateCanvas(folder +
"Background/VsMcBg" + cat, folder +
"Background/VsMcBg" + cat, 1600, 2400);
5792 string stepname =
fHMean.fAnaStepNames[
static_cast<int>(step)];
5796 cbc->GetYaxis()->SetRangeUser(1e-7, 2e-2);
5797 TH1D* r = (TH1D*) cbc->Clone();
5799 r->GetYaxis()->SetRangeUser(0.5, 1.4);
5800 r->GetYaxis()->SetTitle(
"Ratio");
5802 DrawH1({cbc, mcb}, {
"CB_{calc}",
"B_{MC}"},
kLinear,
kLog,
true, 0.65, 0.75, 0.95, 0.9,
"p");
5803 fHMean.DrawAnaStepOnPad(step);
5806 fHMean.DrawAnaStepOnPad(step);
5822 pu->GetYaxis()->SetRangeUser(1e-7, 1e-2);
5823 fHMean.fHM.CreateCanvas(folder +
"Background/CbTypes", folder +
"Background/CbTypes", 800, 800);
5824 DrawH1({pu, puEl, u, uEl}, {
"PLUTO + UrQMD",
"PLUTO + UrQMD (electrons)",
"UrQMD",
"UrQMD (electrons)"},
kLinear,
5825 kLog,
true, 0.55, 0.78, 0.99, 0.95,
"pe");
5826 fHMean.DrawAnaStepOnPad(step);
5834 for (
const string type : {
"Mc",
"Npm"}) {
5835 for (
const string cat : {
""}) {
5839 fHMean.fHM.CreateCanvas(folder +
"Signal/VsInput" + type + cat, folder +
"Signal/VsInput_" + type + cat, 800,
5841 for (
auto step :
fHMean.fAnaSteps) {
5849 for (
int iB = 1; iB <= sig->GetNbinsX(); iB++) {
5850 double bW = sig->GetBinWidth(iB);
5853 sig->SetBinError(iB, err);
5855 sig->GetYaxis()->SetRangeUser(1e-10, 5e-2);
5857 string sigLabel =
"Signal";
5859 DrawH1({sig, coc}, {sigLabel,
"MC Cocktail"},
kLinear,
kLog,
true, 0.65, 0.78, 0.91, 0.91,
"pe");
5870 for (
const string cat : {
"",
"_trueEl"}) {
5875 TH1D* coc =
VaryBinWidth(
"cocktail" + cat +
"_" +
fHMean.fAnaStepNames[
static_cast<int>(step)],
"",
5877 TH1D* sbg = (TH1D*) bg->Clone();
5880 TH1D* ratS = (TH1D*) npm->Clone();
5881 TH1D* ratB = (TH1D*) cb->Clone();
5886 fHMean.fHM.CreateCanvas(folder +
"Signal/N+-VsBgCoc" + cat, folder +
"Signal/N+-VsBgCoc" + cat, 1800, 900);
5889 npm->GetYaxis()->SetRangeUser(1e-9, 2e-2);
5890 DrawH1({npm, sbg, cb, bg}, {
"N^{+-}_{same}",
"B_{MC} + Cocktail",
"CB_{calc}",
"B_{MC}"},
kLinear,
kLog,
true,
5891 0.7, 0.7, 0.99, 0.95,
"p");
5892 fHMean.DrawAnaStepOnPad(step);
5894 ratS->GetYaxis()->SetTitle(
"Ratio");
5895 ratS->GetYaxis()->SetRangeUser(0., 2.);
5896 ratB->GetYaxis()->SetRangeUser(0., 2.);
5897 DrawH1({ratS, ratB}, {
"N^{+-}_{same} / B_{MC}+Coc",
"CB_{calc} / B_{MC}"},
kLinear,
kLinear,
true, 0.7, 0.8, 0.99,
5899 fHMean.DrawAnaStepOnPad(step);
5975 string hName =
"hmx_Minv";
5976 string cbName =
"hCb";
5978 TFile* file =
new TFile((
fDataDir +
"Background/EventMix/eventmix.all.root").c_str());
5979 fHEvMix.fHM.ReadFromFile(file);
5985 LOG(info) <<
"fEvMixDepth: " <<
fEvMixDepth <<
" events.";
5988 if (diff <= 0.995 || diff >= 1.005)
5989 LOG(warn) <<
"Number of LmvmEventMix Events differs significantly from number of simulated events!";
5993 if ((h1->GetNbinsX() != h2->GetNbinsX()) || (h1->GetXaxis()->GetBinLowEdge(1) != h2->GetXaxis()->GetBinLowEdge(1))
5994 || (h1->GetXaxis()->GetBinUpEdge(h1->GetNbinsX()) != (h2->GetXaxis()->GetBinUpEdge(h2->GetNbinsX()))))
5995 LOG(error) <<
"LmvmDrawAll::CalculateCombBGHistos: Minv histograms from LmvmTask and LmvmEventMix do not agree!";
6000 for (
const string cat : {
"",
"_trueEl",
"_urqmdAll",
"_urqmdEl"}) {
6001 TH1D* pmSame =
fHEvMix.H1Clone(hName +
"PM" + cat +
"_sameEv", step);
6002 TH1D* ppSame =
fHEvMix.H1Clone(hName +
"PP" + cat +
"_sameEv", step);
6003 TH1D* mmSame =
fHEvMix.H1Clone(hName +
"MM" + cat +
"_sameEv", step);
6004 TH1D* pmMixed =
fHEvMix.H1Clone(hName +
"PM" + cat +
"_mixedEv", step);
6005 TH1D* ppMixed =
fHEvMix.H1Clone(hName +
"PP" + cat +
"_mixedEv", step);
6006 TH1D* mmMixed =
fHEvMix.H1Clone(hName +
"MM" + cat +
"_mixedEv", step);
6009 TH1D* hGeomSame =
fHMean.CreateHByClone(
fHEvMix.H1Clone(
fHMean.GetName(hName +
"MM" + cat +
"_sameEv", step)),
6010 fHMean.GetName(cbName +
"Geom" + cat +
"_sameEv", step));
6011 TH1D* hGeomMixed =
fHMean.CreateHByClone(
fHEvMix.H1Clone(
fHMean.GetName(hName +
"MM" + cat +
"_sameEv", step)),
6012 fHMean.GetName(cbName +
"Geom" + cat +
"_mixedEv", step));
6014 for (
int iB = 1; iB <= hGeomSame->GetNbinsX(); iB++) {
6015 double cppSame = ppSame->GetBinContent(iB);
6016 double cmmSame = mmSame->GetBinContent(iB);
6017 double cSame = std::sqrt(cppSame * cmmSame);
6018 hGeomSame->SetBinContent(iB, cSame);
6019 double cppMix = ppMixed->GetBinContent(iB);
6020 double cmmMix = mmMixed->GetBinContent(iB);
6021 double cMix = std::sqrt(cppMix * cmmMix);
6022 hGeomMixed->SetBinContent(iB, cMix);
6027 TH1D* hKRaw = (TH1D*) pmMixed->Clone();
6028 hKRaw->Divide(hGeomMixed);
6030 fHMean.fHM.Add(
fHMean.GetName(cbName +
"K_raw" + cat, step), hKRaw);
6033 TH1D* hKFit = (TH1D*) hKRaw->Clone();
6034 double kFitMin = 1.0;
6035 int changeBinK = hKFit->FindBin(kFitMin);
6036 hKFit->Fit(
"pol0",
"Q",
"", kFitMin, hKFit->GetMaximum());
6037 TF1* func = hKFit->GetFunction(
"pol0");
6038 double mean = (func !=
nullptr) ? func->GetParameter(0) : 0;
6039 for (
int iB = changeBinK; iB <= hKFit->GetNbinsX(); iB++) {
6040 hKFit->SetBinContent(iB, mean);
6042 fHMean.fHM.Add(
fHMean.GetName(cbName +
"K_fit" + cat, step), hKRaw);
6045 TH1D* hKOne = (TH1D*) hKRaw->Clone();
6046 double kOneStart = 0.7;
6047 int binOneStart = hKRaw->FindBin(kOneStart);
6048 for (
int iB = binOneStart; iB <= hKRaw->GetNbinsX(); iB++) {
6049 hKOne->SetBinContent(iB, 1.);
6051 fHMean.fHM.Add(
fHMean.GetName(cbName +
"K_one" + cat, step), hKRaw);
6054 TH1D* hK = (TH1D*) hKRaw->Clone();
6055 fHMean.fHM.Add(
fHMean.GetName(cbName +
"K" + cat, step), hK);
6062 int ChangeBin = pmSame->FindBin(
fCbChange);
6066 ppSame->Integral(NormRangeMinBin, NormRangeMaxBin) / ppMixed->Integral(NormRangeMinBin, NormRangeMaxBin);
6068 mmSame->Integral(NormRangeMinBin, NormRangeMaxBin) / mmMixed->Integral(NormRangeMinBin, NormRangeMaxBin);
6070 fCbNormFactor[
static_cast<int>(step)] = 1. / ((normPPInt + normMMInt) / 2.);
6073 TH1D* hBc =
fHMean.CreateHByClone(
fHEvMix.H1Clone(
fHMean.GetName(hName +
"MM" + cat +
"_sameEv", step)),
6074 fHMean.GetName(cbName +
"Bg" + cat, step));
6075 for (
int iB = 1; iB <= hBc->GetNbinsX(); iB++) {
6076 double k = hK->GetBinContent(iB);
6077 double geomMeanMixNormed = hGeomMixed->GetBinContent(iB) * std::sqrt(normPPInt * normMMInt);
6078 double cb = (iB < ChangeBin) ? 2. * k * hGeomSame->GetBinContent(iB) : 2. * k * geomMeanMixNormed;
6079 hBc->SetBinContent(iB, cb);
6083 TH1D* hBcSame =
fHMean.CreateHByClone<TH1D>(cbName +
"Geom" + cat +
"_sameEv", cbName +
"BgSame" + cat, step);
6084 hBcSame->Multiply(hK);
6089 TH1D* hSigNpm =
fHEvMix.H1Clone(hName +
"PM" + cat +
"_sameEv", step);
6090 fHMean.fHM.Add(
fHMean.GetName(cbName +
"SigNpm" + cat, step), hSigNpm);
6091 hSigNpm->Add(
fHMean.H1(cbName +
"Bg" + cat, step), -1.);
6094 TH1D* hSigCoc =
fHMean.H1Clone(
"hMinv" + cat +
"_bg", step);
6095 fHMean.fHM.Add(
fHMean.GetName(cbName +
"SigMc" + cat, step), hSigCoc);
6097 hSigCoc->Add(
fHMean.H1(cbName +
"Bg" + cat, step), -1.);
6102 TH1D* h_s_pm =
VaryBinWidth(
"evmix",
fHMean.GetName(hName +
"PM" + cat +
"_sameEv", step), bwv);
6103 TH1D* h_s_pp =
VaryBinWidth(
"evmix",
fHMean.GetName(hName +
"PP" + cat +
"_sameEv", step), bwv);
6104 TH1D* h_s_mm =
VaryBinWidth(
"evmix",
fHMean.GetName(hName +
"MM" + cat +
"_sameEv", step), bwv);
6105 TH1D* h_m_pm =
VaryBinWidth(
"evmix",
fHMean.GetName(hName +
"PM" + cat +
"_mixedEv", step), bwv);
6106 TH1D* h_m_pp =
VaryBinWidth(
"evmix",
fHMean.GetName(hName +
"PP" + cat +
"_mixedEv", step), bwv);
6107 TH1D* h_m_mm =
VaryBinWidth(
"evmix",
fHMean.GetName(hName +
"MM" + cat +
"_mixedEv", step), bwv);
6123 LOG(error) <<
"LmvmDrawAll::CalculateCombBGHistos: Check bwVar-tag assignment!";
6124 TH1D* bg =
VaryBinWidth(
"mean", cbName +
"Bg" + cat, step, bwv);
6125 TH1D* sigMc =
VaryBinWidth(
"mean", cbName +
"SigMc" + cat, step, bwv);
6126 TH1D* sigNpm =
VaryBinWidth(
"mean", cbName +
"SigNpm" + cat, step, bwv);
6129 (TH1D*) bg->Clone();
6130 TH1D* trueerror = (TH1D*) bg->Clone();
6131 approxerror->Reset();
6134 int nofBins = h_s_pm->GetNbinsX();
6135 double bW2 = h_s_pm->GetBinWidth(1);
6136 for (
int iB = 1; iB <= nofBins; iB++) {
6145 double d_m_pm = std::sqrt(N_s_pp * N_s_mm / (N_m_pp * N_m_mm));
6147 -0.5 * N_m_pm * std::pow(N_m_mm, -1. / 2.) * std::pow(N_m_pp, -3. / 2.) * std::sqrt(N_s_pp * N_s_mm);
6149 -0.5 * N_m_pm * std::pow(N_m_mm, -3. / 2.) * std::pow(N_m_pp, -1. / 2.) * std::sqrt(N_s_pp * N_s_mm);
6150 double d_s_pp = 0.5 * N_m_pm * std::sqrt(N_s_mm) * std::pow(N_m_pp * N_m_mm * N_s_pp, -1. / 2.);
6151 double d_s_mm = 0.5 * N_m_pm * std::sqrt(N_s_pp) * std::pow(N_m_pp * N_m_mm * N_s_mm, -1. / 2.);
6154 double f_m_pm = std::pow(d_m_pm * std::sqrt(N_m_pm), 2);
6155 double f_m_pp = std::pow(d_m_pp * std::sqrt(N_m_pp), 2);
6156 double f_m_mm = std::pow(d_m_mm * std::sqrt(N_m_mm), 2);
6157 double f_s_pp = std::pow(d_s_pp * std::sqrt(N_s_pp), 2);
6158 double f_s_mm = std::pow(d_s_mm * std::sqrt(N_s_mm), 2);
6161 double errorBc = std::sqrt(f_m_pm + f_m_pp + f_m_mm + f_s_pp + f_s_mm) / (bW2 *
fNofEvMixEvents);
6164 double errorSig = std::sqrt(N_s_pm + std::pow(errorBc, 2)) / (bW2 *
fNofEvMixEvents);
6165 double errorSig2 = std::sqrt(N_s_pm / std::pow(bW2 *
fNofEvMixEvents, 2)
6168 if (iB < ChangeBin) bg->SetBinError(iB, errorBc);
6169 if (iB >= ChangeBin) bg->SetBinError(iB, errorBc2);
6170 if (iB < ChangeBin) sigMc->SetBinError(iB, errorSig);
6171 if (iB >= ChangeBin) sigMc->SetBinError(iB, errorSig2);
6172 if (iB < ChangeBin) sigNpm->SetBinError(iB, errorSig);
6173 if (iB >= ChangeBin) sigNpm->SetBinError(iB, errorSig2);
6176 double errappr = std::sqrt(N_s_mm + N_s_pp) / (bW2 *
fNofEvMixEvents);
6179 if (iB < ChangeBin) {
6180 trueerror->SetBinContent(iB, errorBc);
6181 approxerror->SetBinContent(iB, errappr);
6183 if (iB >= ChangeBin) {
6184 trueerror->SetBinContent(iB, errorBc2);
6185 approxerror->SetBinContent(iB, errappr2);
6188 fHMean.CreateHByClone(bg,
fHMean.GetName(cbName +
"Bg" + cat +
"_" + tag, step));
6189 fHMean.CreateHByClone(sigMc,
fHMean.GetName(cbName +
"SigMc" + cat +
"_" + tag, step));
6190 fHMean.CreateHByClone(sigNpm,
fHMean.GetName(cbName +
"SigNpm" + cat +
"_" + tag, step));
6192 fHMean.fHM.CreateCanvas(
"CombinatorialBackground/Error",
"CombinatorialBackground/Error", 800, 800);
6193 DrawH1({trueerror, approxerror}, {
"exact error",
"approx. error"},
kLinear,
kLog,
true, 0.65, 0.75, 0.89,
6201 LOG(info) <<
"fCbNormFactor[" <<
fHMean.fAnaStepNames[
static_cast<int>(step)] <<
"]: " << f;