51 h_count =
new TH1F(
"h_count",
"", 10, 0., 10.);
52 for (
int i = 0; i < 10; i++) {
53 h_prodvtx[i] =
new TH2F(TString::Format(
"h_prodvtx_%d", i),
"", 200, -400., 400., 200, 0., 400.);
54 h_decavtx[i] =
new TH2F(TString::Format(
"h_decavtx_%d", i),
"", 200, -400., 400., 200, 0., 400.);
55 h_kineticvz[i] =
new TH2F(TString::Format(
"h_kineticvz_%d", i),
"", 200, -400., 400., 1000, 0., 10.);
56 h_kineticvz1[i] =
new TH2F(TString::Format(
"h_kineticvz1_%d", i),
"", 200, -400., 400., 1000, 0., 10.);
57 h_kineticvz2[i] =
new TH2F(TString::Format(
"h_kineticvz2_%d", i),
"", 200, -400., 400., 1000, 0., 10.);
58 h_kineticvz_zoom[i] =
new TH2F(TString::Format(
"h_kineticvz_zoom_%d", i),
"", 200, -400., 400., 1000, 0., 0.01);
59 h_Wkineticvz[i] =
new TH2F(TString::Format(
"h_Wkineticvz_%d", i),
"", 200, -400., 400., 1000, 0., 10.);
60 h_Wkineticvz_zoom[i] =
new TH2F(TString::Format(
"h_Wkineticvz_zoom_%d", i),
"", 200, -400., 400., 1000, 0., 0.01);
61 h_phive[i] =
new TH2F(TString::Format(
"h_phive_%d", i),
"", 360, -180., 180., 1000, 0., 10.);
62 h_phivz[i] =
new TH2F(TString::Format(
"h_phivz_%d", i),
"", 200, -400., 400., 360, -180., 180.);
63 h_rve[i] =
new TH2F(TString::Format(
"h_rve_%d", i),
"", 200, 0., 200., 1000, 0., 10.);
64 h_thetavz[i] =
new TH2F(TString::Format(
"h_thetavz_%d", i),
"", 200, -400., 400., 180, 0., 180.);
75 h_g4_xy =
new TH2F(
"h_g4_xy",
"", 100, -5.99, 5.99, 100, -5.99, 5.99);
76 for (
int i = 0; i < 2; i++) {
77 h_sad_xy[i] =
new TH2F(TString::Format(
"h_sad_xy_%d", i),
"", 100, -5.99, 5.99, 100, -5.99, 5.99);
78 h_sad_sir[i] =
new TH1F(TString::Format(
"h_sad_sir_%d", i),
"", 100, -3.99, 3.99);
79 h_z[i] =
new TH1F(TString::Format(
"h_z_%d", i),
"", 100, -3.99, 3.99);
80 h_sad_sall[i] =
new TH1F(TString::Format(
"h_sad_sall_%d", i),
"", 100, -1499.99, 1499.99);
81 h_sad_sE[i] =
new TH2F(TString::Format(
"h_sad_sE_%d", i),
"", 100, -3.99, 3.99, 1000, 0., 10.);
82 h_sad_sraw[i] =
new TH1F(TString::Format(
"h_sad_sraw_%d", i),
"", 100, -1499.99, 1499.99);
83 h_sad_E[i] =
new TH1F(TString::Format(
"h_sad_E_%d", i),
"", 8000, 0, 7.99);
84 h_sad_s[i] =
new TH1F(TString::Format(
"h_sad_s_%d", i),
"", 400, -3.99, 3.99);
86 h_dpx =
new TH1F(
"h_dpx",
"", 1000, -1., 1.);
87 h_dpy =
new TH1F(
"h_dpy",
"", 1000, -1., 1.);
88 h_dE =
new TH1F(
"h_dE",
"", 1000, -1., 1.);
89 h_px =
new TH1F(
"h_px",
"", 10000, -10., 10.);
90 h_py =
new TH1F(
"h_py",
"", 10000, -10., 10.);
91 h_pz =
new TH1F(
"h_pz",
"", 10000, -10., 10.);
92 h_dx =
new TH1F(
"h_dx",
"", 10000, -400., 400.);
93 h_dy =
new TH1F(
"h_dy",
"", 10000, -400., 400.);
94 h_dz =
new TH1F(
"h_dz",
"", 10000, -400., 400.);
95 h_E =
new TH1F(
"h_E",
"", 8000, 0., 7.99);
96 h_P =
new TH1F(
"h_P",
"", 8000, 0., 7.99);
126 for (
const auto& sadMetaHit : sadMetaHits) {
127 px = sadMetaHit.getpx();
128 py = sadMetaHit.getpy();
129 x = sadMetaHit.getx();
130 y = sadMetaHit.gety();
131 s = sadMetaHit.gets();
134 double sraw = sadMetaHit.getsraw();
135 E = sadMetaHit.getE();
136 rate = sadMetaHit.getrate();
143 if (-400. < s && s < 400.) {
155 for (
const auto& mcParticle : mcParticles) {
156 int PDG = mcParticle.getPDG();
158 float Mass = mcParticle.getMass();
159 float Energy = mcParticle.getEnergy();
160 float Kinetic = Energy - Mass;
162 prodvtx[0] = mcParticle.getProductionVertex().X();
163 prodvtx[1] = mcParticle.getProductionVertex().Y();
164 prodvtx[2] = mcParticle.getProductionVertex().Z();
165 float prodr =
sqrt(prodvtx[0] * prodvtx[0] + prodvtx[1] * prodvtx[1]);
167 decavtx[0] = mcParticle.getDecayVertex().X();
168 decavtx[1] = mcParticle.getDecayVertex().Y();
169 decavtx[2] = mcParticle.getDecayVertex().Z();
170 float decar =
sqrt(decavtx[0] * decavtx[0] + decavtx[1] * decavtx[1]);
172 mom[0] = mcParticle.getMomentum().X();
173 mom[1] = mcParticle.getMomentum().Y();
174 mom[2] = mcParticle.getMomentum().Z();
175 float theta = mcParticle.getMomentum().Theta() * TMath::RadToDeg();
176 float phi = mcParticle.getMomentum().Phi() * TMath::RadToDeg();
177 int partID[9] = {0, 0, 0, 0, 0, 0, 0, 0, 0};
184 else if (PDG == 1000080160) partID[5] = 1;
185 else if (PDG == 1000060120) partID[6] = 1;
186 else if (PDG == 1000020040) partID[7] = 1;
188 for (
int i = 0; i < 8; i++) {
189 if (partID[i] == 1) {
198 h_phivz[i]->Fill(prodvtx[2], phi);
199 h_phive[i]->Fill(phi, Kinetic);
200 h_rve[i]->Fill(prodr, Kinetic);
218 h_phivz[8]->Fill(prodvtx[2], phi);
219 h_phive[8]->Fill(phi, Kinetic);
220 h_rve[8]->Fill(prodr, Kinetic);
236 h_phivz[9]->Fill(prodvtx[2], phi);
238 h_dx->Fill(x - prodvtx[0]);
239 h_dy->Fill(-y - prodvtx[1]);
240 h_dz->Fill(-s - prodvtx[2]);
242 h_z[0]->Fill(prodvtx[2] / 100.);
243 h_z[1]->Fill(prodvtx[2] / 100., rate);
245 h_dpx->Fill(px - mom[0]);
246 h_dpy->Fill(-py - mom[1]);
247 h_dE->Fill(
E - Energy);
249 h_g4_xy->Fill(prodvtx[0], prodvtx[1]);
255 h_P->Fill(
sqrt(mom[0] * mom[0] + mom[1] * mom[1] + mom[2] * mom[2]));