CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
BaseProtonTransport.h
Go to the documentation of this file.
1 #ifndef BASEPROTONTRANSPORT
2 #define BASEPROTONTRANSPORT
5 #include "HepMC/GenEvent.h"
9 #include <vector>
10 #include <map>
11 #include <string>
12 #include "TLorentzVector.h"
13 
14 namespace CLHEP {
15  class HepRandomEngine;
16 }
18 public:
19  BaseProtonTransport(const edm::ParameterSet& iConfig);
20  virtual ~BaseProtonTransport();
21 
22  std::vector<LHCTransportLink>& getCorrespondenceMap() { return m_CorrespondenceMap; }
23  virtual void process(const HepMC::GenEvent* ev, const edm::EventSetup& es, CLHEP::HepRandomEngine* engine) = 0;
24 
25  void clear();
26 
28 
30  void ApplyBeamCorrection(TLorentzVector& p);
31 
32  double beamEnergy() { return beamEnergy_; };
33  double beamMomentum() { return beamMomentum_; };
34 
35 protected:
36  void setBeamFileNames(const std::string& nam1, const std::string& nam2) {
37  beam1Filename_ = nam1;
38  beam2Filename_ = nam2;
39  };
40 
41  void setBeamParameters(double stx, double sty, double sx, double sy, double se) {
42  m_sigmaSTX = stx;
43  m_sigmaSTY = sty;
44  m_sigmaSX = sx;
45  m_sigmaSY = sy;
46  m_sig_E = se;
47  };
48 
49  void setCrossingAngles(double cx45, double cx56, double cy45, double cy56) {
50  fCrossingAngleX_45_ = cx45;
51  fCrossingAngleX_56_ = cx56;
52  fCrossingAngleY_45_ = cy45;
53  fCrossingAngleY_56_ = cy56;
54  };
55 
58 
59  CLHEP::HepRandomEngine* engine_{nullptr};
60  std::vector<LHCTransportLink> m_CorrespondenceMap;
61  std::map<unsigned int, TLorentzVector> m_beamPart;
62  std::map<unsigned int, double> m_xAtTrPoint;
63  std::map<unsigned int, double> m_yAtTrPoint;
64 
65  bool verbosity_;
69 
72 
75  double fCrossingAngleX_45_{0.0};
76  double fCrossingAngleX_56_{0.0};
77  double fCrossingAngleY_45_{0.0};
78  double fCrossingAngleY_56_{0.0};
79 
80  double etaCut_;
81  double momentumCut_;
82 
83 private:
84  double beamMomentum_;
85  double beamEnergy_;
86 
87  double m_sigmaSTX{0.0};
88  double m_sigmaSTY{0.0};
89  double m_sigmaSX{0.0};
90  double m_sigmaSY{0.0};
91  double m_sig_E{0.0};
92 };
93 #endif
virtual void process(const HepMC::GenEvent *ev, const edm::EventSetup &es, CLHEP::HepRandomEngine *engine)=0
BaseProtonTransport(const edm::ParameterSet &iConfig)
std::vector< LHCTransportLink > & getCorrespondenceMap()
void setBeamParameters(double stx, double sty, double sx, double sy, double se)
CLHEP::HepRandomEngine * engine_
void setCrossingAngles(double cx45, double cx56, double cy45, double cy56)
void setBeamFileNames(const std::string &nam1, const std::string &nam2)
std::map< unsigned int, double > m_xAtTrPoint
void ApplyBeamCorrection(HepMC::GenParticle *p)
std::map< unsigned int, TLorentzVector > m_beamPart
std::map< unsigned int, double > m_yAtTrPoint
std::vector< LHCTransportLink > m_CorrespondenceMap
void addPartToHepMC(const HepMC::GenEvent *, HepMC::GenEvent *)