Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-04-22 03:16:49

0001 // -*- C++ -*-
0002 //
0003 // Package:    FastSimulation/CTPPSFastTrackingProducer
0004 // Class:      CTPPSFastTrackingProducer
0005 //
0006 /**\class CTPPSFastTrackingProducer CTPPSFastTrackingProducer.cc FastSimulation/CTPPSFastTrackingProducer/plugins/CTPPSFastTrackingProducer.cc
0007 
0008 Description: [one line class summary]
0009 
0010 Implementation:
0011 [Notes on implementation]
0012 */
0013 //
0014 // Original Author:  Sandro Fonseca De Souza
0015 //         Created:  Thu, 29 Sep 2016 16:13:41 GMT
0016 //
0017 //
0018 
0019 #include <memory>
0020 
0021 #include "DataFormats/Common/interface/Handle.h"
0022 #include "DataFormats/Math/interface/Point3D.h"
0023 #include "DataFormats/Math/interface/Vector.h"
0024 #include "DataFormats/Math/interface/Vector3D.h"
0025 #include "FWCore/Framework/interface/ConsumesCollector.h"
0026 #include "FWCore/Framework/interface/Event.h"
0027 #include "FWCore/Framework/interface/MakerMacros.h"
0028 #include "FWCore/Framework/interface/stream/EDProducer.h"
0029 #include "FWCore/ParameterSet/interface/FileInPath.h"
0030 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0031 #include "FWCore/Utilities/interface/isFinite.h"
0032 #include "FWCore/Utilities/interface/StreamID.h"
0033 #include "FastSimDataFormats/CTPPSFastSim/interface/CTPPSFastRecHit.h"
0034 #include "FastSimDataFormats/CTPPSFastSim/interface/CTPPSFastRecHitContainer.h"
0035 #include "FastSimDataFormats/CTPPSFastSim/interface/CTPPSFastTrack.h"
0036 #include "FastSimDataFormats/CTPPSFastSim/interface/CTPPSFastTrackContainer.h"
0037 #include "FastSimulation/CTPPSFastGeometry/interface/CTPPSToFDetector.h"
0038 #include "FastSimulation/CTPPSFastGeometry/interface/CTPPSTrkDetector.h"
0039 #include "Utilities/PPS/interface/PPSUnitConversion.h"
0040 #include "Utilities/PPS/interface/PPSUtilities.h"
0041 
0042 #include "TLorentzVector.h"
0043 
0044 // hector includes
0045 #include "H_Parameters.h"
0046 #include "H_BeamLine.h"
0047 #include "H_RecRPObject.h"
0048 #include "H_BeamParticle.h"
0049 
0050 class CTPPSFastTrackingProducer : public edm::stream::EDProducer<> {
0051 public:
0052   explicit CTPPSFastTrackingProducer(const edm::ParameterSet&);
0053 
0054 private:
0055   void produce(edm::Event&, const edm::EventSetup&) override;
0056 
0057   // ----------member data ---------------------------
0058 
0059   typedef std::vector<CTPPSFastRecHit> CTPPSFastRecHitContainer;
0060   edm::EDGetTokenT<CTPPSFastRecHitContainer> _recHitToken;
0061   void ReadRecHits(edm::Handle<CTPPSFastRecHitContainer>&);
0062   void FastReco(int Direction, H_RecRPObject* station);
0063   void Reconstruction();
0064   void ReconstructArm(
0065       H_RecRPObject* pps_station, double x1, double y1, double x2, double y2, double& tx, double& ty, double& eloss);
0066   void MatchCellId(int cellId, std::vector<int> vrecCellId, std::vector<double> vrecTof, bool& match, double& recTof);
0067   bool SearchTrack(int,
0068                    int,
0069                    int Direction,
0070                    double& xi,
0071                    double& t,
0072                    double& partP,
0073                    double& pt,
0074                    double& thx,
0075                    double& thy,
0076                    double& x0,
0077                    double& y0,
0078                    double& xt,
0079                    double& yt,
0080                    double& X1d,
0081                    double& Y1d,
0082                    double& X2d,
0083                    double& Y2d);
0084   void TrackerStationClear();
0085   void TrackerStationStarting();
0086   void ProjectToToF(const double x1, const double y1, const double x2, const double y2, double& xt, double& yt) {
0087     xt = ((fz_timing - fz_tracker2) * (x2 - x1) / (fz_tracker2 - fz_tracker1)) + x2;
0088     yt = ((fz_timing - fz_tracker2) * (y2 - y1) / (fz_tracker2 - fz_tracker1)) + y2;
0089   };
0090   // Hector objects
0091   bool SetBeamLine();
0092 
0093   std::unique_ptr<H_BeamLine> m_beamlineCTPPS1;
0094   std::unique_ptr<H_BeamLine> m_beamlineCTPPS2;
0095   std::unique_ptr<H_RecRPObject> pps_stationF;
0096   std::unique_ptr<H_RecRPObject> pps_stationB;
0097 
0098   std::string beam1filename;
0099   std::string beam2filename;
0100 
0101   // Defaults
0102   double lengthctpps;
0103   bool m_verbosity;
0104   double fBeamEnergy;
0105   double fBeamMomentum;
0106   bool fCrossAngleCorr;
0107   double fCrossingAngleBeam1;
0108   double fCrossingAngleBeam2;
0109   ////////////////////////////////////////////////
0110   std::unique_ptr<CTPPSTrkStation> TrkStation_F;  // auxiliary object with the tracker geometry
0111   std::unique_ptr<CTPPSTrkStation> TrkStation_B;
0112   std::unique_ptr<CTPPSTrkDetector> det1F;
0113   std::unique_ptr<CTPPSTrkDetector> det1B;
0114   std::unique_ptr<CTPPSTrkDetector> det2F;
0115   std::unique_ptr<CTPPSTrkDetector> det2B;
0116   std::unique_ptr<CTPPSToFDetector> detToF_F;
0117   std::unique_ptr<CTPPSToFDetector> detToF_B;
0118 
0119   std::vector<CTPPSFastTrack> theCTPPSFastTrack;
0120 
0121   CTPPSFastTrack track;
0122 
0123   std::vector<int> recCellId_F, recCellId_B;
0124   std::vector<double> recTof_F, recTof_B;
0125 
0126   double fz_tracker1, fz_tracker2, fz_timing;
0127   double fTrackerWidth, fTrackerHeight, fTrackerInsertion, fBeamXRMS_Trk1, fBeamXRMS_Trk2, fTrk1XOffset, fTrk2XOffset;
0128   std::vector<double> fToFCellWidth;
0129   double fToFCellHeight, fToFPitchX, fToFPitchY;
0130   int fToFNCellX, fToFNCellY;
0131   double fToFInsertion, fBeamXRMS_ToF, fToFXOffset, fTimeSigma, fImpParcut;
0132 };
0133 //////////////////////
0134 // constructors and destructor
0135 //
0136 CTPPSFastTrackingProducer::CTPPSFastTrackingProducer(const edm::ParameterSet& iConfig)
0137     : m_verbosity(false), fBeamMomentum(0.), fCrossAngleCorr(false), fCrossingAngleBeam1(0.), fCrossingAngleBeam2(0.) {
0138   //register your products
0139   produces<edm::CTPPSFastTrackContainer>("CTPPSFastTrack");
0140   using namespace edm;
0141   _recHitToken = consumes<CTPPSFastRecHitContainer>(iConfig.getParameter<edm::InputTag>("recHitTag"));
0142   m_verbosity = iConfig.getParameter<bool>("Verbosity");
0143   // User definitons
0144 
0145   // Read beam parameters needed for Hector reconstruction
0146   lengthctpps = iConfig.getParameter<double>("BeamLineLengthCTPPS");
0147   beam1filename = iConfig.getParameter<string>("Beam1");
0148   beam2filename = iConfig.getParameter<string>("Beam2");
0149   fBeamEnergy = iConfig.getParameter<double>("BeamEnergy");  // beam energy in GeV
0150   fBeamMomentum = sqrt(fBeamEnergy * fBeamEnergy - PPSTools::ProtonMassSQ);
0151   fCrossingAngleBeam1 = iConfig.getParameter<double>("CrossingAngleBeam1");
0152   fCrossingAngleBeam2 = iConfig.getParameter<double>("CrossingAngleBeam2");
0153 
0154   if (fCrossingAngleBeam1 != 0 || fCrossingAngleBeam2 != 0)
0155     fCrossAngleCorr = true;
0156   //Read detectors positions and parameters
0157 
0158   fz_tracker1 = iConfig.getParameter<double>("Z_Tracker1");
0159   fz_tracker2 = iConfig.getParameter<double>("Z_Tracker2");
0160   fz_timing = iConfig.getParameter<double>("Z_Timing");
0161   //
0162   fTrackerWidth = iConfig.getParameter<double>("TrackerWidth");
0163   fTrackerHeight = iConfig.getParameter<double>("TrackerHeight");
0164   fTrackerInsertion = iConfig.getParameter<double>("TrackerInsertion");
0165   fBeamXRMS_Trk1 = iConfig.getParameter<double>("BeamXRMS_Trk1");
0166   fBeamXRMS_Trk2 = iConfig.getParameter<double>("BeamXRMS_Trk2");
0167   fTrk1XOffset = iConfig.getParameter<double>("Trk1XOffset");
0168   fTrk2XOffset = iConfig.getParameter<double>("Trk2XOffset");
0169   fToFCellWidth = iConfig.getUntrackedParameter<std::vector<double> >("ToFCellWidth");
0170   fToFCellHeight = iConfig.getParameter<double>("ToFCellHeight");
0171   fToFPitchX = iConfig.getParameter<double>("ToFPitchX");
0172   fToFPitchY = iConfig.getParameter<double>("ToFPitchY");
0173   fToFNCellX = iConfig.getParameter<int>("ToFNCellX");
0174   fToFNCellY = iConfig.getParameter<int>("ToFNCellY");
0175   fToFInsertion = iConfig.getParameter<double>("ToFInsertion");
0176   fBeamXRMS_ToF = iConfig.getParameter<double>("BeamXRMS_ToF");
0177   fToFXOffset = iConfig.getParameter<double>("ToFXOffset");
0178   fTimeSigma = iConfig.getParameter<double>("TimeSigma");
0179   fImpParcut = iConfig.getParameter<double>("ImpParcut");
0180 
0181   if (!SetBeamLine()) {
0182     if (m_verbosity)
0183       LogDebug("CTPPSFastTrackingProducer") << "CTPPSFastTrackingProducer: WARNING: lengthctpps=  " << lengthctpps;
0184     return;
0185   }
0186 
0187   // Create a particle to get the beam energy from the beam file
0188   // Take care: the z inside the station is in meters
0189   //
0190   //Tracker Detector Description
0191   det1F = std::make_unique<CTPPSTrkDetector>(
0192       fTrackerWidth, fTrackerHeight, fTrackerInsertion * fBeamXRMS_Trk1 + fTrk1XOffset);
0193   det2F = std::make_unique<CTPPSTrkDetector>(
0194       fTrackerWidth, fTrackerHeight, fTrackerInsertion * fBeamXRMS_Trk2 + fTrk2XOffset);
0195   det1B = std::make_unique<CTPPSTrkDetector>(
0196       fTrackerWidth, fTrackerHeight, fTrackerInsertion * fBeamXRMS_Trk1 + fTrk1XOffset);
0197   det2B = std::make_unique<CTPPSTrkDetector>(
0198       fTrackerWidth, fTrackerHeight, fTrackerInsertion * fBeamXRMS_Trk2 + fTrk2XOffset);
0199 
0200   //Timing Detector Description
0201   std::vector<double> vToFCellWidth;
0202   vToFCellWidth.reserve(8);
0203   for (int i = 0; i < 8; i++) {
0204     vToFCellWidth.push_back(fToFCellWidth[i]);
0205   }
0206   double pos_tof = fToFInsertion * fBeamXRMS_ToF + fToFXOffset;
0207   detToF_F = std::make_unique<CTPPSToFDetector>(
0208       fToFNCellX, fToFNCellY, vToFCellWidth, fToFCellHeight, fToFPitchX, fToFPitchY, pos_tof, fTimeSigma);
0209   detToF_B = std::make_unique<CTPPSToFDetector>(
0210       fToFNCellX, fToFNCellY, vToFCellWidth, fToFCellHeight, fToFPitchX, fToFPitchY, pos_tof, fTimeSigma);
0211   //
0212 }
0213 
0214 // ------------ method called to produce the data  ------------
0215 void CTPPSFastTrackingProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0216   using namespace edm;
0217   TrackerStationStarting();
0218   Handle<CTPPSFastRecHitContainer> recHits;
0219   iEvent.getByToken(_recHitToken, recHits);
0220   recCellId_F.clear();
0221   recCellId_B.clear();
0222   recTof_F.clear();
0223   recTof_B.clear();
0224   ReadRecHits(recHits);
0225   Reconstruction();
0226   TrackerStationClear();
0227 
0228   std::unique_ptr<CTPPSFastTrackContainer> output_tracks(new CTPPSFastTrackContainer);
0229   int n = 0;
0230   for (std::vector<CTPPSFastTrack>::const_iterator i = theCTPPSFastTrack.begin(); i != theCTPPSFastTrack.end(); i++) {
0231     output_tracks->push_back(*i);
0232     n += 1;
0233   }
0234 
0235   iEvent.put(std::move(output_tracks), "CTPPSFastTrack");
0236 }  //end
0237 
0238 /////////////////////////
0239 void CTPPSFastTrackingProducer::TrackerStationClear() {
0240   TrkStation_F->first.clear();
0241   TrkStation_F->second.clear();
0242   TrkStation_B->first.clear();
0243   TrkStation_B->second.clear();
0244 }
0245 /////////////////////////
0246 void CTPPSFastTrackingProducer::TrackerStationStarting() {
0247   det1F->clear();
0248   det1B->clear();
0249   det2F->clear();
0250   det2B->clear();
0251   detToF_F->clear();
0252   detToF_B->clear();
0253 }
0254 
0255 ////////////////////////////
0256 void CTPPSFastTrackingProducer::ReadRecHits(edm::Handle<CTPPSFastRecHitContainer>& recHits) {
0257   // DetId codification for PSimHit taken from CTPPSPixel- It will be replaced by CTPPSDetId
0258   // 2014314496 -> Tracker1 zPositive
0259   // 2014838784 -> Tracker2 zPositive
0260   // 2046820352 -> Timing   zPositive
0261   // 2031091712 -> Tracker1 zNegative
0262   // 2031616000 -> Tracker2 zNegative
0263   // 2063597568 -> Timing   zNegative
0264 
0265   for (unsigned int irecHits = 0; irecHits < recHits->size(); ++irecHits) {
0266     const CTPPSFastRecHit* recHitDet = &(*recHits)[irecHits];
0267     unsigned int detlayerId = recHitDet->detUnitId();
0268     double x = recHitDet->entryPoint().x();
0269     double y = recHitDet->entryPoint().y();
0270     double z = recHitDet->entryPoint().z();
0271     float tof = recHitDet->tof();
0272     if (detlayerId == 2014314496)
0273       det1F->AddHit(detlayerId, x, y, z);
0274     else if (detlayerId == 2014838784)
0275       det2F->AddHit(detlayerId, x, y, z);
0276     else if (detlayerId == 2031091712)
0277       det1B->AddHit(detlayerId, x, y, z);
0278     else if (detlayerId == 2031616000)
0279       det2B->AddHit(detlayerId, x, y, z);
0280     else if (detlayerId == 2046820352) {
0281       detToF_F->AddHit(x, y, tof);
0282       recCellId_F.push_back(detToF_F->findCellId(x, y));
0283       recTof_F.push_back(tof);
0284     } else if (detlayerId == 2063597568) {
0285       detToF_B->AddHit(x, y, tof);
0286       recCellId_B.push_back(detToF_B->findCellId(x, y));
0287       recTof_B.push_back(tof);
0288     }
0289 
0290   }  //LOOP TRK
0291   //creating Stations
0292   TrkStation_F = std::make_unique<CTPPSTrkStation>(*det1F, *det2F);
0293   TrkStation_B = std::make_unique<CTPPSTrkStation>(*det1B, *det2B);
0294 }  // end function
0295 
0296 void CTPPSFastTrackingProducer::Reconstruction() {
0297   theCTPPSFastTrack.clear();
0298   int Direction;
0299   Direction = 1;  //cms positive Z / forward
0300   FastReco(Direction, &*pps_stationF);
0301   Direction = -1;  //cms negative Z / backward
0302   FastReco(Direction, &*pps_stationB);
0303 }  //end Reconstruction
0304 
0305 bool CTPPSFastTrackingProducer::SearchTrack(int i,
0306                                             int j,
0307                                             int Direction,
0308                                             double& xi,
0309                                             double& t,
0310                                             double& partP,
0311                                             double& pt,
0312                                             double& thx,
0313                                             double& thy,
0314                                             double& x0,
0315                                             double& y0,
0316                                             double& xt,
0317                                             double& yt,
0318                                             double& X1d,
0319                                             double& Y1d,
0320                                             double& X2d,
0321                                             double& Y2d) {
0322   // Given 1 hit in Tracker1 and 1 hit in Tracker2 try to make a track with Hector
0323   xi = 0;
0324   t = 0;
0325   partP = 0;
0326   pt = 0;
0327   x0 = 0.;
0328   y0 = 0.;
0329   xt = 0.;
0330   yt = 0.;
0331   X1d = 0.;
0332   Y1d = 0.;
0333   X2d = 0.;
0334   Y2d = 0.;
0335   CTPPSTrkDetector* det1 = nullptr;
0336   CTPPSTrkDetector* det2 = nullptr;
0337   H_RecRPObject* station = nullptr;
0338   // Separate in forward and backward stations according to direction
0339   if (Direction > 0) {
0340     det1 = &(TrkStation_F->first);
0341     det2 = &(TrkStation_F->second);
0342     station = &*pps_stationF;
0343   } else {
0344     det1 = &(TrkStation_B->first);
0345     det2 = &(TrkStation_B->second);
0346     station = &*pps_stationB;
0347   }
0348   if (det1->ppsNHits_ <= i || det2->ppsNHits_ <= j)
0349     return false;
0350   //
0351   double x1 = det1->ppsX_.at(i);
0352   double y1 = det1->ppsY_.at(i);
0353   double x2 = det2->ppsX_.at(j);
0354   double y2 = det2->ppsY_.at(j);
0355   double eloss;
0356 
0357   //thx and thy are returned in microrad
0358   ReconstructArm(
0359       station, Direction * x1, y1, Direction * x2, y2, thx, thy, eloss);  // Pass the hits in the LHC ref. frame
0360   thx *= -Direction;                                                      //  invert to the CMS ref frame
0361 
0362   // Protect for unphysical results
0363   if (edm::isNotFinite(eloss) || edm::isNotFinite(thx) || edm::isNotFinite(thy))
0364     return false;
0365   //
0366 
0367   if (m_verbosity)
0368     LogDebug("CTPPSFastTrackingProducer::SearchTrack:") << "thx " << thx << " thy " << thy << " eloss " << eloss;
0369 
0370   // Get the start point of the reconstructed track near the origin made by Hector in the CMS ref. frame
0371   x0 = -Direction * station->getX0() * um_to_cm;
0372   y0 = station->getY0() * um_to_cm;
0373   double ImpPar = sqrt(x0 * x0 + y0 * y0);
0374   if (ImpPar > fImpParcut)
0375     return false;
0376   if (eloss < 0. || eloss > fBeamEnergy)
0377     return false;
0378   //
0379   // Calculate the reconstructed track parameters
0380   double theta = sqrt(thx * thx + thy * thy) * urad;
0381   xi = eloss / fBeamEnergy;
0382   double energy = fBeamEnergy * (1. - xi);
0383   partP = sqrt(energy * energy - PPSTools::ProtonMassSQ);
0384   t = -2. * (PPSTools::ProtonMassSQ - fBeamEnergy * energy + fBeamMomentum * partP * cos(theta));
0385   pt = partP * theta;
0386   if (xi < 0. || xi > 1. || t < 0. || t > 10. || pt <= 0.) {
0387     xi = 0.;
0388     t = 0.;
0389     partP = 0.;
0390     pt = 0.;
0391     x0 = 0.;
0392     y0 = 0.;
0393     return false;  // unphysical values
0394   }
0395   //Try to include the timing detector in the track
0396   ProjectToToF(x1, y1, x2, y2, xt, yt);  // the projections is done in the CMS ref frame
0397   X1d = x1;
0398   Y1d = y1;
0399   X2d = x2;
0400   Y2d = y2;
0401   return true;
0402 }  //end  SearchTrack
0403 
0404 void CTPPSFastTrackingProducer::ReconstructArm(
0405     H_RecRPObject* pps_station, double x1, double y1, double x2, double y2, double& tx, double& ty, double& eloss) {
0406   tx = 0.;
0407   ty = 0.;
0408   eloss = 0.;
0409   if (!pps_station)
0410     return;
0411   x1 *= mm_to_um;
0412   x2 *= mm_to_um;
0413   y1 *= mm_to_um;
0414   y2 *= mm_to_um;
0415   pps_station->setPositions(x1, y1, x2, y2);
0416   double energy = pps_station->getE(AM);  // dummy call needed to calculate some Hector internal parameter
0417   if (edm::isNotFinite(energy))
0418     return;
0419   tx = pps_station->getTXIP();  // change orientation to CMS
0420   ty = pps_station->getTYIP();
0421   eloss = pps_station->getE();
0422 }
0423 
0424 void CTPPSFastTrackingProducer::MatchCellId(
0425     int cellId, std::vector<int> vrecCellId, std::vector<double> vrecTof, bool& match, double& recTof) {
0426   for (unsigned int i = 0; i < vrecCellId.size(); i++) {
0427     if (cellId == vrecCellId.at(i)) {
0428       match = true;
0429       recTof = vrecTof.at(i);
0430       continue;
0431     }
0432   }
0433 }
0434 
0435 void CTPPSFastTrackingProducer::FastReco(int Direction, H_RecRPObject* station) {
0436   double theta = 0.;
0437   double xi, t, partP, pt, phi, x0, y0, thx, thy, xt, yt, X1d, Y1d, X2d, Y2d;
0438   CTPPSTrkDetector* Trk1 = nullptr;
0439   CTPPSTrkDetector* Trk2 = nullptr;
0440   double pos_tof = fToFInsertion * fBeamXRMS_ToF + fToFXOffset;
0441   int cellId = 0;
0442   std::vector<double> vToFCellWidth;
0443   vToFCellWidth.reserve(8);
0444   for (int i = 0; i < 8; i++) {
0445     vToFCellWidth.push_back(fToFCellWidth[i]);
0446   }
0447   CTPPSToFDetector* ToF = new CTPPSToFDetector(
0448       fToFNCellX, fToFNCellY, vToFCellWidth, fToFCellHeight, fToFPitchX, fToFPitchY, pos_tof, fTimeSigma);
0449   if (Direction > 0) {
0450     Trk1 = &(TrkStation_F->first);
0451     Trk2 = &(TrkStation_F->second);
0452   } else {
0453     Trk1 = &(TrkStation_B->first);
0454     Trk2 = &(TrkStation_B->second);
0455   }
0456   // Make a track from EVERY pair of hits combining Tracker1 and Tracker2.
0457   // The tracks may not be independent as 1 hit may belong to more than 1 track.
0458   for (int i = 0; i < (int)Trk1->ppsNHits_; i++) {
0459     for (int j = 0; j < (int)Trk2->ppsNHits_; j++) {
0460       if (SearchTrack(i, j, Direction, xi, t, partP, pt, thx, thy, x0, y0, xt, yt, X1d, Y1d, X2d, Y2d)) {
0461         // Check if the hitted timing cell matches the reconstructed track
0462         cellId = ToF->findCellId(xt, yt);
0463         double recTof = 0.;
0464         bool matchCellId = false;
0465         if (Direction > 0) {
0466           theta = sqrt(thx * thx + thy * thy) * urad;
0467           MatchCellId(cellId, recCellId_F, recTof_F, matchCellId, recTof);
0468         } else if (Direction < 0) {
0469           theta = TMath::Pi() - sqrt(thx * thx + thy * thy) * urad;
0470           MatchCellId(cellId, recCellId_B, recTof_B, matchCellId, recTof);
0471         }
0472         phi = atan2(thy, thx);  // at this point, thx is already in the cms ref. frame
0473 
0474         double px = partP * sin(theta) * cos(phi);
0475         double py = partP * sin(theta) * sin(phi);
0476         double pz = partP * cos(theta);
0477         double e = sqrt(partP * partP + PPSTools::ProtonMassSQ);
0478         TLorentzVector p(px, py, pz, e);
0479         // Invert the Lorentz boost made to take into account the crossing angle during simulation
0480         if (fCrossAngleCorr) {
0481           PPSTools::LorentzBoost(p, "MC", {fCrossingAngleBeam1, fCrossingAngleBeam2, fBeamMomentum, fBeamEnergy});
0482         }
0483         //Getting the Xi and t (squared four momentum transferred) of the reconstructed track
0484         PPSTools::Get_t_and_xi(const_cast<TLorentzVector*>(&p), t, xi, {fBeamMomentum, fBeamEnergy});
0485         double pxx = p.Px();
0486         double pyy = p.Py();
0487         double pzz = p.Pz();
0488         math::XYZVector momentum(pxx, pyy, pzz);
0489         math::XYZPoint vertex(x0, y0, 0);
0490 
0491         track.setp(momentum);
0492         track.setvertex(vertex);
0493         track.sett(t);
0494         track.setxi(xi);
0495         track.setx1(X1d);
0496         track.sety1(Y1d);
0497         track.setx2(X2d);
0498         track.sety2(Y2d);
0499         if (matchCellId) {
0500           track.setcellid(cellId);
0501           track.settof(recTof);
0502         } else {
0503           track.setcellid(0);
0504           track.settof(0.);
0505         }
0506         theCTPPSFastTrack.push_back(track);
0507       }
0508     }
0509   }
0510 }  //end FastReco
0511 
0512 bool CTPPSFastTrackingProducer::SetBeamLine() {
0513   edm::FileInPath b1(beam1filename.c_str());
0514   edm::FileInPath b2(beam2filename.c_str());
0515   if (lengthctpps <= 0)
0516     return false;
0517   m_beamlineCTPPS1 = std::make_unique<H_BeamLine>(-1, lengthctpps + 0.1);  // (direction, length)
0518   m_beamlineCTPPS1->fill(b2.fullPath(), 1, "IP5");
0519   m_beamlineCTPPS2 = std::make_unique<H_BeamLine>(1, lengthctpps + 0.1);  //
0520   m_beamlineCTPPS2->fill(b1.fullPath(), 1, "IP5");
0521   m_beamlineCTPPS1->offsetElements(120, 0.097);
0522   m_beamlineCTPPS2->offsetElements(120, -0.097);
0523   pps_stationF = std::make_unique<H_RecRPObject>(fz_tracker1, fz_tracker2, *m_beamlineCTPPS1);
0524   pps_stationB = std::make_unique<H_RecRPObject>(fz_tracker1, fz_tracker2, *m_beamlineCTPPS2);
0525   return true;
0526 }
0527 //define this as a plug-in
0528 DEFINE_FWK_MODULE(CTPPSFastTrackingProducer);