diff --git a/Generators/ForeseeGenerator/python/Validate.py b/Generators/ForeseeGenerator/python/Validate.py index 8f27b36ab8820a66b19133cf8bb526092cd53ad8..0112555e0d384f7d37ee8e2b93cc7bb074775fe2 100644 --- a/Generators/ForeseeGenerator/python/Validate.py +++ b/Generators/ForeseeGenerator/python/Validate.py @@ -81,7 +81,7 @@ class EvgenValidation(EvgenAnalysisAlg): def initialize(self): # All daughters - self.hists.add("PIDs", 60, -30, 30) + self.hists.add("PIDs", 600, -300, 300) # Daughter i tbins, pbins = self.binning() @@ -93,7 +93,7 @@ class EvgenValidation(EvgenAnalysisAlg): self.hists.add(f"Theta_d{i}", 20, 0, 0.001) self.hists.add(f"Phi_d{i}", 16, -3.2, 3.2) self.hists.add(f"ThetaVsP_d{i}", arrayX = tbins, arrayY = pbins) - self.hists.add(f"Mass_d{i}", 200, 0, 0.01) + self.hists.add(f"Mass_d{i}", 5000, 0, 1) # Mother self.hists.add("E_M", 1000, 0, 10000) @@ -106,11 +106,17 @@ class EvgenValidation(EvgenAnalysisAlg): self.hists.add("ThetaVsP_M", arrayX = tbins, arrayY = pbins) # Vertex - self.hists.add("Vtx_X", 50, -100, 100) - self.hists.add("Vtx_Y", 50, -100, 100) - self.hists.add("Vtx_Z", 500, -5000, 0) - self.hists.add("Vtx_R", 20, 0, 200) - self.hists.add("Vtx_XY", 50, -100, 100, 50, -100, 100) + self.hists.add("Vtx_X_LLP", 50, -100, 100) + self.hists.add("Vtx_Y_LLP", 50, -100, 100) + self.hists.add("Vtx_Z_LLP", 500, -5000, 0) + self.hists.add("Vtx_R_LLP", 20, 0, 200) + self.hists.add("Vtx_XY_LLP", 50, -100, 100, 50, -100, 100) + + self.hists.add("Vtx_X_All", 50, -100, 100) + self.hists.add("Vtx_Y_All", 50, -100, 100) + self.hists.add("Vtx_Z_All", 500, -5000, 0) + self.hists.add("Vtx_R_All", 20, 0, 200) + self.hists.add("Vtx_XY_All", 50, -100, 100, 50, -100, 100) return StatusCode.Success @@ -136,12 +142,12 @@ class EvgenValidation(EvgenAnalysisAlg): self.hists["PIDs"].Fill(p.pdg_id(), self.weight) return - def fillVertex(self, v): - self.hists["Vtx_X"].Fill(v.x(), self.weight) - self.hists["Vtx_Y"].Fill(v.y(), self.weight) - self.hists["Vtx_Z"].Fill(v.z(), self.weight) - self.hists["Vtx_XY"].Fill(v.x(), v.y(), self.weight) - self.hists["Vtx_R"].Fill(sqrt(v.x()**2 + v.y()**2), self.weight) + def fillVertex(self, label, v): + self.hists[f"Vtx_X_{label}"].Fill(v.x(), self.weight) + self.hists[f"Vtx_Y_{label}"].Fill(v.y(), self.weight) + self.hists[f"Vtx_Z_{label}"].Fill(v.z(), self.weight) + self.hists[f"Vtx_XY_{label}"].Fill(v.x(), v.y(), self.weight) + self.hists[f"Vtx_R_{label}"].Fill(sqrt(v.x()**2 + v.y()**2), self.weight) return @@ -151,26 +157,38 @@ class EvgenValidation(EvgenAnalysisAlg): # Loop over all particles in events (assuming mother not stored) momenta = [] + vertices = [] mother = HepMC.FourVector(0,0,0,0) llp_vtx = None for i, p in enumerate(evt.particles): #p.print() self.fillDaughter(p) + momenta.append(p.momentum()) mother += p.momentum() - if i == 0: + + if i == 0 and p.production_vertex(): #p.production_vertex().print() llp_vtx = p.production_vertex().point3d() + if p.production_vertex(): + vertices.append(p.production_vertex().point3d()) + # Fill daughter plots for i in range(self.ndaughters): + if i >= len(momenta): continue self.fillKin(f"d{i}", momenta[i]) # Fill mother plots self.fillKin("M", mother, mass = True) - # Fill vertex plots - self.fillVertex(llp_vtx) + # Fill all vertices + for v in vertices: + self.fillVertex("All", v) + + # Fill LLP vertex plots + if llp_vtx: + self.fillVertex("LLP", llp_vtx) return StatusCode.Success diff --git a/Generators/ForeseeGenerator/share/plot_validation.py b/Generators/ForeseeGenerator/share/plot_validation.py index 5f702f6953dc891b1878719b36daff4714aee36b..49436e95fe348844da228d7eff38059139a37b40 100644 --- a/Generators/ForeseeGenerator/share/plot_validation.py +++ b/Generators/ForeseeGenerator/share/plot_validation.py @@ -90,7 +90,7 @@ if __name__ == "__main__": for i in range(args.ndaughters): config = [Hist(f"P_d{i}", logy = True, xtitle = "p^{0} [GeV]", ndiv = 5, r = 5), Hist(f"Theta_d{i}", xtitle = "#theta [rad]", ndiv = -4), - Hist(f"Mass_d{i}", xtitle = "m^{0} [GeV]", xlo = 0, xhi = 0.001, ndiv = 4), + Hist(f"Mass_d{i}", xtitle = "m^{0} [GeV]", xlo = 0, xhi = 0.1, ndiv = 4), Hist(f"Pt_d{i}", logy = True, xtitle = "p_{T}^{0} [GeV]", ndiv = 10, r = 5), Hist(f"Phi_d{i}", xtitle = "#phi [rad]"), Hist(f"ThetaVsP_d{i}", ytitle = "p^{0} [GeV]", xtitle = "#theta [rad]", logx = True, logy = True, d = "colz") @@ -98,9 +98,9 @@ if __name__ == "__main__": plotn(f, args, config, 3, 2, f"daug{i}") - config = [Hist("P_M", logy = True, xtitle = "p^{0} [GeV]", ndiv = 5, xlo = 0, xhi = 2000), + config = [Hist("P_M", logy = True, xtitle = "p^{0} [GeV]", ndiv = 5, xlo = 0, xhi = 10000), Hist("Theta_M", xtitle = "#theta [rad]", ndiv = -4, r = 10), - Hist("Mass_M", xtitle = "m^{0} [GeV]", xlo = 0, xhi = 0.2, ndiv = 5), + Hist("Mass_M", xtitle = "m^{0} [GeV]", xlo = 0, xhi = 1., ndiv = 5), Hist("Pt_M", logy = True, xtitle = "p_{T}^{0} [GeV]", ndiv = 10, r = 50), Hist("Phi_M", xtitle = "#phi [rad]", r = 2), Hist("ThetaVsP_M", ytitle = "p^{0} [GeV]", xtitle = "#theta [rad]", logx = True, logy = True, d = "colz") @@ -108,7 +108,7 @@ if __name__ == "__main__": plotn(f, args, config, 3, 2, "mother") - plotn(f, args, Hist("PIDs", xtitle="PDG Id"), 1, 1, "pid") + plotn(f, args, Hist("PIDs", xtitle="PDG Id", xlo=-30, xhi=30), 1, 1, "pid") # config = [Hist("ThetaVsP_M", ytitle = "p^{0} [GeV]", xtitle = "#theta [rad]", logx = True, logy = True, d = "colz"), @@ -118,11 +118,34 @@ if __name__ == "__main__": # # plotn(f, args, config, 2, 2, "twod") - config = [Hist("Vtx_X", xtitle = "x [mm]", r = 5), - Hist("Vtx_Y", xtitle = "y [mm]", r = 5), - Hist("Vtx_Z", xtitle = "z [mm]", r = 5, ndiv = 5), - Hist("Vtx_XY", xtitle = "x [mm]", ytitle = "y [mm]", d = "colz", r = (5,5)), - Hist("Vtx_R", xtitle = "r [mm]", r = 5, ndiv = 5) + config = [Hist("Vtx_X_LLP", xtitle = "x [mm]", r = 5), + Hist("Vtx_Y_LLP", xtitle = "y [mm]", r = 5), + Hist("Vtx_Z_LLP", xtitle = "z [mm]", r = 5, ndiv = 5), + Hist("Vtx_XY_LLP", xtitle = "x [mm]", ytitle = "y [mm]", d = "colz", r = (5,5)), + Hist("Vtx_R_LLP", xtitle = "r [mm]", r = 5, ndiv = 5) ] - plotn(f, args, config, 3, 2, "vtx") + plotn(f, args, config, 3, 2, "vtx_llp") + + + config = [Hist("Vtx_X_All", xtitle = "x [mm]", r = 5), + Hist("Vtx_Y_All", xtitle = "y [mm]", r = 5), + Hist("Vtx_Z_All", xtitle = "z [mm]", r = 5, ndiv = 5), + Hist("Vtx_XY_All", xtitle = "x [mm]", ytitle = "y [mm]", d = "colz", r = (5,5)), + Hist("Vtx_R_All", xtitle = "r [mm]", r = 5, ndiv = 5) + ] + + plotn(f, args, config, 3, 2, "vtx_all") + + +# config = [Hist("Vtx_X", xtitle = "x [mm]", r = 5), +# Hist("Vtx_Y", xtitle = "y [mm]", r = 5), +# Hist("Vtx_Z", xtitle = "z [mm]", r = 5, ndiv = 5), +# Hist("Vtx_XY", xtitle = "x [mm]", ytitle = "y [mm]", d = "colz", r = (5,5)), +# Hist("Vtx_R", xtitle = "r [mm]", r = 5, ndiv = 5) +# ] +# +# plotn(f, args, config, 3, 2, "vtx") + + +