-
Notifications
You must be signed in to change notification settings - Fork 1
/
dp_dalitz_rdf.C
executable file
·99 lines (64 loc) · 4.07 KB
/
dp_dalitz_rdf.C
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
{//begin
#include <TMath.h>
#include <math.h>
#include "scripts/fit1MeV_GaussianPlusCBWithExp_redo.C"
using namespace std;
using namespace ROOT;
EnableImplicitMT();
//create rdataframe for dp magdown data
RDataFrame dpdf("D2KKpi/DecayTree", "/share/lazy/D2KKpi/dp2kkpi_magdown.root");
const double phi_pm = 12;
const double phiupperbound = (1019.455+phi_pm)*(1019.455+phi_pm)/1000000;
const double philowerbound = (1019.455-phi_pm)*(1019.455-phi_pm)/1000000;
auto inv_m_func = [](double px1, double py1, double pz1, double pe1, double px2, double py2, double pz2, double pe2) {return TLorentzVector(TLorentzVector(px1, py1, pz1, pe1)+TLorentzVector(px2, py2, pz2, pe2)).Mag2()/1000000 ;};
auto prob_func = [](double prob1, double prob2) {return TMath::Log(prob1) - TMath::Log(prob2) ;};
auto cut_ipchi2 = [](double x) {return x < 5 ;};
auto cut_fdchi2 = [](double x) {return x > 175 ;};
auto cut_endvertexchi2 = [](double x) {return x < 3 ;};
auto cut_phi = [phiupperbound, philowerbound](double x) {return x > philowerbound && x < phiupperbound ;};
auto cut_prob_5 = [] (double x) {return x>5 ;};
auto cut_prob_0 = [] (double x) {return x>5 ;};
auto dp_og = dpdf.Define("mkpkm", inv_m_func, {"Kplus_PX", "Kplus_PY", "Kplus_PZ", "Kplus_PE", "Kminus_PX", "Kminus_PY", "Kminus_PZ", "Kminus_PE"})
.Define("mkmpi", inv_m_func, {"Piplus_PX", "Piplus_PY", "Piplus_PZ", "Piplus_PE", "Kminus_PX", "Kminus_PY", "Kminus_PZ", "Kminus_PE"});
auto dp_cut = dp_og.Filter(cut_ipchi2, {"Dplus_IPCHI2_OWNPV"})
.Filter(cut_fdchi2, {"Dplus_FDCHI2_OWNPV"})
.Filter(cut_endvertexchi2, {"Dplus_ENDVERTEX_CHI2"})
.Define("kminuslog", prob_func, {"Kminus_MC15TuneV1_ProbNNk", "Kminus_MC15TuneV1_ProbNNpi"})
.Define("kpluslog", prob_func, {"Kplus_MC15TuneV1_ProbNNk", "Kplus_MC15TuneV1_ProbNNpi"})
.Define("pipluslog", prob_func, {"Piplus_MC15TuneV1_ProbNNpi", "Piplus_MC15TuneV1_ProbNNk"})
.Filter(cut_prob_5, {"kminuslog"})
.Filter(cut_prob_5, {"kpluslog"})
.Filter(cut_prob_0, {"pipluslog"});
auto dp_phicut = dp_cut.Filter(cut_phi, {"mkpkm"});
gStyle->SetPalette(55, 0);
auto dalitz_og = dp_og.Fill<double, double>(TH2D("dalitz_og", "D^{+} #rightarrow K^{+}K^{-}#pi^{+} Dalitz Plot (no cuts)",200,0.8,3.3,200,0.1,2.3), {"mkpkm", "mkmpi"});
dalitz_og->SetStats(0);
dalitz_og->GetXaxis()->SetTitle("m^{2}(K^{+}K^{-}) [GeV/c^{2}]^{2}");
dalitz_og->GetYaxis()->SetTitle("m^{2}(K^{-}#pi^{+}) [GeV/c^{2}]^{2}");
auto dalitz_og_ROOT = dalitz_og->DrawCopy();//converts Rhist into ROOThist
auto dalitz_og_can = new TCanvas("dalitz_og_can", "dalitz_og_can", 1000, 800);
dalitz_og_can->SetRightMargin(0.15);
dalitz_og_can->cd();
dalitz_og_ROOT->Draw("colz");
dalitz_og_can->SaveAs("image/dp_dalitz_nocuts.png");
auto dalitz_cut = dp_cut.Fill<double, double>(TH2D("dalitz_cut", "D^{+} #rightarrow K^{+}K^{-}#pi^{+} Dalitz Plot (cut not on phi)",200,0.8,3.3,200,0.1,2.3), {"mkpkm", "mkmpi"});
dalitz_cut->SetStats(0);
dalitz_cut->GetXaxis()->SetTitle("m^{2}(K^{+}K^{-}) [GeV/c^{2}]^{2}");
dalitz_cut->GetYaxis()->SetTitle("m^{2}(K^{-}#pi^{+}) [GeV/c^{2}]^{2}");
auto dalitz_cut_ROOT = dalitz_cut->DrawCopy();//converts Rhist into ROOThist
auto dalitz_cut_can = new TCanvas("dalitz_cut_can", "dalitz_cut_can", 1000, 800);
dalitz_cut_can->SetRightMargin(0.15);
dalitz_cut_can->cd();
dalitz_cut_ROOT->Draw("colz");
dalitz_cut_can->SaveAs("image/dp_dalitz_somecuts.png");
auto dalitz_phicut = dp_phicut.Fill<double, double>(TH2D("dalitz_phicut", "D^{+} #rightarrow K^{+}K^{-}#pi^{+} Dalitz Plot (cut on phi)",200,0.8,3.3,200,0.1,2.3), {"mkpkm", "mkmpi"});
dalitz_phicut->SetStats(0);
dalitz_phicut->GetXaxis()->SetTitle("m^{2}(K^{+}K^{-}) [GeV/c^{2}]^{2}");
dalitz_phicut->GetYaxis()->SetTitle("m^{2}(K^{-}#pi^{+}) [GeV/c^{2}]^{2}");
auto dalitz_phicut_ROOT = dalitz_phicut->DrawCopy();//converts Rhist into ROOThist
auto dalitz_phicut_can = new TCanvas("dalitz_phicut_can", "dalitz_phicut_can", 1000, 800);
dalitz_phicut_can->SetRightMargin(0.15);
dalitz_phicut_can->cd();
dalitz_phicut_ROOT->Draw("colz");
dalitz_phicut_can->SaveAs("image/dp_dalitz_phicut.png");
}//end