#include <TotemTransport.h>

Inheritance diagram for TotemTransport:

Public Member Functions
void	process (const HepMC::GenEvent ev, const edm::EventSetup &es, CLHEP::HepRandomEngine engine) override

	TotemTransport (const edm::ParameterSet &ps, bool verbosity)

	TotemTransport ()

	~TotemTransport () override

Public Member Functions inherited from ProtonTransport
void	addPartToHepMC (HepMC::GenEvent *)

void	ApplyBeamCorrection (HepMC::GenParticle *p)

void	ApplyBeamCorrection (TLorentzVector &p)

void	clear ()

std::vector< LHCTransportLink > &	getCorrespondenceMap ()

	ProtonTransport ()

virtual	~ProtonTransport ()

Private Member Functions
LHCOpticsApproximator *	ReadParameterization (const std::string &, const std::string &)

bool	transportProton (const HepMC::GenParticle *)

Private Attributes
LHCOpticsApproximator *	m_aprox_ip_150_l = 0

LHCOpticsApproximator *	m_aprox_ip_150_r = 0

double	m_beampipe_aperture_radius

double	m_fEtacut

double	m_fMomentumMin

std::string	m_model_ip_150_l_name

double	m_model_ip_150_l_zmax

double	m_model_ip_150_l_zmin

std::string	m_model_ip_150_r_name

double	m_model_ip_150_r_zmax

double	m_model_ip_150_r_zmin

std::string	m_model_root_file_l

std::string	m_model_root_file_r

edm::FileInPath	m_opticsFileBeam1_

edm::FileInPath	m_opticsFileBeam2_

edm::ParameterSet	m_parameters

edm::ESHandle< ParticleDataTable >	m_pdt

edm::EDGetTokenT< edm::HepMCProduct >	m_protonsToken_

bool	m_verbosity

double	m_Zin_

double	m_Zout_

Additional Inherited Members
Protected Types inherited from ProtonTransport
enum	TransportMode { TransportMode::HECTOR, TransportMode::TOTEM }

Protected Attributes inherited from ProtonTransport
bool	bApplyZShift

CLHEP::HepRandomEngine *	engine

double	fBeamEnergy

double	fBeamMomentum

double	fBeamXatIP

double	fBeamYatIP

double	fCrossingAngle_45

double	fCrossingAngle_56

double	fPPSRegionStart_45

double	fPPSRegionStart_56

double	fVtxMeanX

double	fVtxMeanY

double	fVtxMeanZ

std::map< unsigned int, TLorentzVector >	m_beamPart

std::vector< LHCTransportLink >	m_CorrespondenceMap

double	m_sig_E

double	m_sigmaSTX

double	m_sigmaSTY

double	m_sigmaSX

double	m_sigmaSY

bool	m_verbosity

std::map< unsigned int, double >	m_xAtTrPoint

std::map< unsigned int, double >	m_yAtTrPoint

TransportMode	MODE

int	NEvent

Detailed Description

Definition at line 24 of file TotemTransport.h.

Constructor & Destructor Documentation

TotemTransport::TotemTransport	(	const edm::ParameterSet &	ps,
		bool	verbosity
	)

Definition at line 15 of file TotemTransport.cc.

                                                                              :ProtonTransport(),
         m_parameters(iConfig.getParameter<edm::ParameterSet>("BeamProtTransportSetup")),
         m_verbosity(iConfig.getParameter<bool>("Verbosity")),
         m_model_root_file_r(m_parameters.getParameter<std::string>("ModelRootFile_R")),
         m_model_root_file_l(m_parameters.getParameter<std::string>("ModelRootFile_L")),
         m_model_ip_150_r_name(m_parameters.getParameter<std::string>("Model_IP_150_R_Name")),
         m_model_ip_150_l_name(m_parameters.getParameter<std::string>("Model_IP_150_L_Name")),
         m_model_ip_150_r_zmin(m_parameters.getParameter<double>("Model_IP_150_R_Zmin")),
         m_model_ip_150_r_zmax(m_parameters.getParameter<double>("Model_IP_150_R_Zmax")),
         m_model_ip_150_l_zmin(m_parameters.getParameter<double>("Model_IP_150_L_Zmin")),
         m_model_ip_150_l_zmax(m_parameters.getParameter<double>("Model_IP_150_L_Zmax")), 
         m_beampipe_aperture_radius(m_parameters.getParameter<double>("BeampipeApertureRadius"))
 {
         fBeamEnergy= m_parameters.getParameter<double>("sqrtS");
         m_sigmaSTX = m_parameters.getParameter<double>("beamDivergenceX");
         m_sigmaSTY = m_parameters.getParameter<double>("beamDivergenceY");
         m_sig_E    = m_parameters.getParameter<double>("beamEnergyDispersion");
         fCrossingAngle_45 = m_parameters.getParameter<double>("halfCrossingAngleSector45");
         fCrossingAngle_56 = m_parameters.getParameter<double>("halfCrossingAngleSector56");
         fVtxMeanX       = iConfig.getParameter<double>("VtxMeanX");
         fVtxMeanY       = iConfig.getParameter<double>("VtxMeanY");
         fVtxMeanZ       = iConfig.getParameter<double>("VtxMeanZ");
         fBeamXatIP      = m_parameters.getUntrackedParameter<double>("BeamXatIP",fVtxMeanX);
         fBeamYatIP      = m_parameters.getUntrackedParameter<double>("BeamYatIP",fVtxMeanY);
         bApplyZShift    = m_parameters.getParameter<bool>("ApplyZShift");
 
         MODE = TransportMode::TOTEM;
         edm::LogInfo("ProtonTransport")
              <<"=============================================================================\n"
              <<"             Bulding LHC Proton transporter based on TOTEM model\n"
              <<"=============================================================================\n";
 
         fBeamMomentum        = sqrt(fBeamEnergy*fBeamEnergy - ProtonMassSQ);
 
         fPPSRegionStart_56=m_model_ip_150_r_zmax;
         fPPSRegionStart_45=m_model_ip_150_l_zmax;
         
         m_aprox_ip_150_r = ReadParameterization(m_model_ip_150_r_name,m_model_root_file_r);
         m_aprox_ip_150_l = ReadParameterization(m_model_ip_150_l_name,m_model_root_file_l);
 
         if (m_aprox_ip_150_r == nullptr || m_aprox_ip_150_l == nullptr) {
                 edm::LogError("ProtonTransport") << "Parameterisation "
                         << m_model_ip_150_r_name << " or " << m_model_ip_150_l_name << " missing in file. Cannot proceed. ";
                 exit(1);
         }
         edm::LogInfo("TotemRPProtonTransportSetup") <<
                 "Parameterizations read from file, pointers:" << m_aprox_ip_150_r << " " << m_aprox_ip_150_l << " ";
 }

TotemTransport::TotemTransport ( )

Definition at line 10 of file TotemTransport.cc.

References ProtonTransport::MODE, and ProtonTransport::TOTEM.

10 :ProtonTransport(){MODE=TransportMode::TOTEM;};

ProtonTransport::TransportMode::TOTEM

ProtonTransport::MODE

TransportMode MODE

Definition: ProtonTransport.h:29

ProtonTransport::ProtonTransport

ProtonTransport()

Definition: ProtonTransport.cc:7

TotemTransport::~TotemTransport ( )

override

Definition at line 11 of file TotemTransport.cc.

References ProtonTransport::clear().

 {
       this->clear();
 }

Member Function Documentation

void TotemTransport::process	(	const HepMC::GenEvent *	ev,
		const edm::EventSetup &	es,
		CLHEP::HepRandomEngine *	engine
	)

overridevirtual

Adds the stable protons from the event ev to a beamline

Implements ProtonTransport.

Definition at line 63 of file TotemTransport.cc.

References ProtonTransport::addPartToHepMC(), ProtonTransport::engine, GenParticle::GenParticle, mps_splice::line, ProtonTransport::m_beamPart, and transportProton().

 {
         engine=_engine;
 
         for (HepMC::GenEvent::particle_const_iterator eventParticle =evt->particles_begin(); eventParticle != evt->particles_end(); ++eventParticle ) {
          if (!((*eventParticle)->status() == 1 && (*eventParticle)->pdg_id()==2212 )) continue;
          unsigned int line = (*eventParticle)->barcode();
          HepMC::GenParticle * gpart = (*eventParticle);
          if ( gpart->pdg_id()!=2212 ) continue; // only transport stable protons
          if ( gpart->status()!=1 /*&& gpart->status()<83 */) continue;
          if ( m_beamPart.find(line) != m_beamPart.end() ) continue;
 
 
          transportProton(gpart);
  
      }
      addPartToHepMC(const_cast<HepMC::GenEvent*>(evt));
 }

LHCOpticsApproximator * TotemTransport::ReadParameterization	(	const std::string &	m_model_name,
		const std::string &	rootfile
	)

private

Definition at line 151 of file TotemTransport.cc.

References f, and MillePedeFileConverter_cfg::fileName.

Referenced by TotemTransport().

 {
     edm::FileInPath fileName(rootfile.c_str());
     TFile *f = TFile::Open(fileName.fullPath().c_str(), "read");
     if (!f) {
         edm::LogError("TotemRPProtonTransportSetup") << "File " << fileName << " not found. Exiting.";
         return nullptr;
     }
     edm::LogInfo("TotemRPProtonTransportSetup") << "Root file opened, pointer:" << f;
 
     // read parametrization
     LHCOpticsApproximator* aprox = (LHCOpticsApproximator *) f->Get(m_model_name.c_str());
     f->Close();
     return aprox;
 }

bool TotemTransport::transportProton ( const HepMC::GenParticle * in_trk )

private

Definition at line 81 of file TotemTransport.cc.

References ProtonTransport::ApplyBeamCorrection(), MillePedeFileConverter_cfg::e, ProtonTransport::fBeamXatIP, ProtonTransport::fBeamYatIP, ProtonTransport::fCrossingAngle_45, ProtonTransport::fCrossingAngle_56, ProtonTransport::fVtxMeanX, ProtonTransport::fVtxMeanY, ProtonTransport::fVtxMeanZ, mps_splice::line, LogDebug, m_aprox_ip_150_l, m_aprox_ip_150_r, ProtonTransport::m_beamPart, m_beampipe_aperture_radius, m_model_ip_150_l_zmax, m_model_ip_150_l_zmin, m_model_ip_150_r_zmax, m_model_ip_150_r_zmin, m_verbosity, ProtonTransport::m_xAtTrPoint, ProtonTransport::m_yAtTrPoint, m_Zin_, m_Zout_, ProtonMassSQ, mathSSE::sqrt(), funct::tan(), LHCOpticsApproximator::Transport_m_GeV(), and urad.

Referenced by process().

 {
      edm::LogInfo("ProtonTransport")<<"Starting proton transport using TOTEM method\n";
      ApplyBeamCorrection(const_cast<HepMC::GenParticle*>(in_trk));
 
      const HepMC::GenVertex* in_pos = in_trk->production_vertex();
      const HepMC::FourVector in_mom = in_trk->momentum();
 //
 // ATTENTION: HepMC uses mm, vertex config of CMS uses cm and SimTransport uses mm
 //
      double in_position[3] = {(in_pos->position().x()-fVtxMeanX*cm) / meter+fBeamXatIP*mm/meter,
                               (in_pos->position().y()-fVtxMeanY*cm) / meter+fBeamYatIP*mm/meter,
                               (in_pos->position().z()-fVtxMeanZ*cm) / meter};  // move to z=0 if configured below
 
 // (bApplyZShift) -- The TOTEM parameterization requires the shift to z=0
      double fCrossingAngle = (in_mom.z()>0)?fCrossingAngle_45:-fCrossingAngle_56;
      in_position[0] = in_position[0]+(tan((long double)fCrossingAngle*urad)-((long double)in_mom.x())/((long double)in_mom.z()))*in_position[2];
      in_position[1] = in_position[1]-((long double)in_mom.y())/((long double)in_mom.z())*in_position[2];
      in_position[2] = 0.;
 //
      double in_momentum[3] = {in_mom.x(), in_mom.y() , in_mom.z()};
      double out_position[3];
      double out_momentum[3];
      edm::LogInfo("ProtonTransport") << "before transport ->" <<
              " position: " << in_position[0] << ", " << in_position[1] << ", " << in_position[2] <<
              " momentum: " << in_momentum[0] << ", " << in_momentum[1] << ", " << in_momentum[2];
 
      LHCOpticsApproximator* approximator_= nullptr;
      if (in_mom.z()>0) {approximator_ = m_aprox_ip_150_l; m_Zin_ = m_model_ip_150_l_zmin; m_Zout_ = m_model_ip_150_l_zmax;}
      else              {approximator_ = m_aprox_ip_150_r; m_Zin_ = m_model_ip_150_r_zmin; m_Zout_ = m_model_ip_150_r_zmax;}
 
      bool invert_beam_coord_system=true; // it doesn't matter the option here, it is hard coded as TRUE inside LHCOpticsApproximator!
 
      bool tracked = approximator_->Transport_m_GeV(in_position, in_momentum, out_position, out_momentum, invert_beam_coord_system, m_Zout_ - m_Zin_);
 
      if (!tracked) return false;
 
      edm::LogInfo("ProtonTransport") << "after transport -> " <<
              "position: " << out_position[0] << ", " << out_position[1] << ", " << out_position[2] <<
              "momentum: " << out_momentum[0] << ", " << out_momentum[1] << ", " << out_momentum[2];
 
      if (out_position[0] * out_position[0] + out_position[1] * out_position[1] >
                      m_beampipe_aperture_radius * m_beampipe_aperture_radius) {
              edm::LogInfo("ProtonTransport") << "Proton ouside beampipe";
              edm::LogInfo("ProtonTransport") << "===== END Transport " << "====================";
              return false;
      }
 
     TVector3 out_pos(out_position[0] * meter, out_position[1] * meter, out_position[2] * meter);
     TVector3 out_mom(out_momentum[0], out_momentum[1], out_momentum[2]);
     edm::LogInfo("TotemRPProtonTransportModel") << "output -> " << "position: ";out_pos.Print();
     edm::LogInfo("TotemRPProtonTransportModel") << " momentum: ";out_mom.Print();
      double px = -out_momentum[0];
      double py = out_momentum[1];  // this need to be checked again, since it seems an invertion is occuring in  the prop.
      double pz = out_momentum[2];
      double e = sqrt(px*px+py*py+pz*pz+ProtonMassSQ);
      TLorentzVector p_out(px,py,pz,e);
      double x1_ctpps = -out_position[0]*meter; // Totem parameterization uses meter, one need it in millimeter
      double y1_ctpps = -out_position[1]*meter;
 
      unsigned int line = in_trk->barcode();
 
      if(m_verbosity) LogDebug("ProtonTransportEventProcessing") <<
              "ProtonTransport:filterPPS: barcode = " << line << " x=  "<< x1_ctpps <<" y= " << y1_ctpps;
 
      m_beamPart[line]    = p_out;
      m_xAtTrPoint[line]  = x1_ctpps;
      m_yAtTrPoint[line]  = y1_ctpps;
      return true;
 }