3 #include <CLHEP/Units/SystemOfUnits.h> 4 #include "TLorentzVector.h" 11 m_model_ip_150_r_name_(iConfig.getParameter<
std::
string>(
"Model_IP_150_R_Name")),
12 m_model_ip_150_l_name_(iConfig.getParameter<
std::
string>(
"Model_IP_150_L_Name")),
13 m_beampipe_aperture_radius_(iConfig.getParameter<double>(
"BeampipeApertureRadius")) {
18 double cax45 = iConfig.
getParameter<
double>(
"halfCrossingAngleXSector45");
19 double cax56 = iConfig.
getParameter<
double>(
"halfCrossingAngleXSector56");
21 double stx = iConfig.
getParameter<
double>(
"BeamDivergenceX");
22 double sty = iConfig.
getParameter<
double>(
"BeamDivergenceY");
25 double se = iConfig.
getParameter<
double>(
"BeamEnergyDispersion");
32 <<
"=============================================================================\n" 33 <<
" Bulding LHC Proton transporter based on TOTEM model\n" 34 <<
"=============================================================================\n";
41 <<
" missing in file. Cannot proceed. ";
54 CLHEP::HepRandomEngine* engine) {
60 for (HepMC::GenEvent::particle_const_iterator eventParticle = evt.particles_begin();
61 eventParticle != evt.particles_end();
63 if (!((*eventParticle)->status() == 1 && (*eventParticle)->pdg_id() == 2212))
66 if (!(fabs((*eventParticle)->momentum().eta()) >
etaCut_ && fabs((*eventParticle)->momentum().pz()) >
momentumCut_))
69 unsigned int line = (*eventParticle)->barcode();
71 if (gpart->pdg_id() != 2212)
73 if (gpart->status() != 1)
88 edm::LogVerbatim(
"TotemTransport") <<
"Starting proton transport using TOTEM method\n";
92 const HepMC::GenVertex* in_pos = in_trk->production_vertex();
93 const HepMC::FourVector in_mom = in_trk->momentum();
97 double in_position[3] = {in_pos->position().x(), in_pos->position().y(), in_pos->position().z()};
103 in_position[0] - in_position[2] * (in_mom.x() / in_mom.z() - crossingAngleX *
urad);
104 in_position[1] = in_position[1] - in_position[2] * (in_mom.y() / (in_mom.z()));
106 double in_momentum[3] = {in_mom.x(), in_mom.y(), in_mom.z()};
107 double out_position[3];
108 double out_momentum[3];
110 <<
" position: " << in_position[0] <<
", " << in_position[1] <<
", " 111 << in_position[2] <<
" momentum: " << in_momentum[0] <<
", " << in_momentum[1]
112 <<
", " << in_momentum[2];
117 if (in_mom.z() > 0) {
127 bool invert_beam_coord_system =
131 in_position, in_momentum, out_position, out_momentum, invert_beam_coord_system, zout - zin);
137 <<
"position: " << out_position[0] <<
", " << out_position[1] <<
", " 138 << out_position[2] <<
"momentum: " << out_momentum[0] <<
", " << out_momentum[1]
139 <<
", " << out_momentum[2];
141 if (out_position[0] * out_position[0] + out_position[1] * out_position[1] >
144 <<
"===== END Transport " 145 <<
"====================";
150 TVector3 out_pos(out_position[0] * meter, out_position[1] * meter, out_position[2] * meter);
151 TVector3 out_mom(out_momentum[0], out_momentum[1], out_momentum[2]);
154 LogDebug(
"TotemTransport") <<
"output -> position: ";
156 LogDebug(
"TotemTransport") <<
" momentum: ";
160 double px = -out_momentum[0];
161 double py = out_momentum[1];
165 TLorentzVector p_out(
px,
py, pz,
e);
166 double x1_ctpps = -out_position[0] * meter;
167 double y1_ctpps = out_position[1] * meter;
169 unsigned int line = in_trk->barcode();
172 LogDebug(
"TotemTransport") <<
"TotemTransport:filterPPS: barcode = " <<
line <<
" x= " << x1_ctpps
173 <<
" y= " << y1_ctpps;
184 TFile*
f = TFile::Open(
fileName.fullPath().c_str(),
"read");
Log< level::Info, true > LogVerbatim
T getParameter(std::string const &) const
std::string beam1Filename_
std::string m_model_ip_150_l_name_
std::string beam2Filename_
bool Transport_m_GeV(double in_pos[3], double in_momentum[3], double out_pos[3], double out_momentum[3], bool check_apertures, double z2_z1_dist) const
pos, momentum: x,y,z; pos in m, momentum in GeV/c
double fCrossingAngleX_56_
~TotemTransport() override
Log< level::Error, false > LogError
void setBeamParameters(double stx, double sty, double sx, double sy, double se)
double fCrossingAngleX_45_
double fPPSRegionStart_56_
LHCOpticsApproximator * m_aprox_ip_150_r_
TotemTransport(const edm::ParameterSet &ps)
std::string m_model_ip_150_r_name_
double m_beampipe_aperture_radius_
Class finds the parametrisation of MADX proton transport and transports the protons according to it 5...
LHCOpticsApproximator * m_aprox_ip_150_l_
CLHEP::HepRandomEngine * engine_
void setCrossingAngles(double cx45, double cx56, double cy45, double cy56)
void setBeamFileNames(const std::string &nam1, const std::string &nam2)
bool transportProton(HepMC::GenParticle *)
std::map< unsigned int, double > m_xAtTrPoint
void ApplyBeamCorrection(HepMC::GenParticle *p)
std::map< unsigned int, TLorentzVector > m_beamPart
LHCOpticsApproximator * ReadParameterization(const std::string &, const std::string &)
double fPPSRegionStart_45_
void process(const HepMC::GenEvent *ev, const edm::EventSetup &es, CLHEP::HepRandomEngine *engine) override
static const double ProtonMassSQ
std::map< unsigned int, double > m_yAtTrPoint