CTPPSFastTrackingProducer/plugins/CTPPSFastTrackingProducer.cc

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   for (std::vector<CTPPSFastTrack>::const_iterator i = theCTPPSFastTrack.begin(); i != theCTPPSFastTrack.end(); i++) {
0230     output_tracks->push_back(*i);
0231   }
0232
0233   iEvent.put(std::move(output_tracks), "CTPPSFastTrack");
0234 }  //end
0235
0236 /////////////////////////
0237 void CTPPSFastTrackingProducer::TrackerStationClear() {
0238   TrkStation_F->first.clear();
0239   TrkStation_F->second.clear();
0240   TrkStation_B->first.clear();
0241   TrkStation_B->second.clear();
0242 }
0243 /////////////////////////
0244 void CTPPSFastTrackingProducer::TrackerStationStarting() {
0245   det1F->clear();
0246   det1B->clear();
0247   det2F->clear();
0248   det2B->clear();
0249   detToF_F->clear();
0250   detToF_B->clear();
0251 }
0252
0253 ////////////////////////////
0254 void CTPPSFastTrackingProducer::ReadRecHits(edm::Handle<CTPPSFastRecHitContainer>& recHits) {
0255   // DetId codification for PSimHit taken from CTPPSPixel- It will be replaced by CTPPSDetId
0256   // 2014314496 -> Tracker1 zPositive
0257   // 2014838784 -> Tracker2 zPositive
0258   // 2046820352 -> Timing   zPositive
0259   // 2031091712 -> Tracker1 zNegative
0260   // 2031616000 -> Tracker2 zNegative
0261   // 2063597568 -> Timing   zNegative
0262
0263   for (unsigned int irecHits = 0; irecHits < recHits->size(); ++irecHits) {
0264     const CTPPSFastRecHit* recHitDet = &(*recHits)[irecHits];
0265     unsigned int detlayerId = recHitDet->detUnitId();
0266     double x = recHitDet->entryPoint().x();
0267     double y = recHitDet->entryPoint().y();
0268     double z = recHitDet->entryPoint().z();
0269     float tof = recHitDet->tof();
0270     if (detlayerId == 2014314496)
0271       det1F->AddHit(detlayerId, x, y, z);
0272     else if (detlayerId == 2014838784)
0273       det2F->AddHit(detlayerId, x, y, z);
0274     else if (detlayerId == 2031091712)
0275       det1B->AddHit(detlayerId, x, y, z);
0276     else if (detlayerId == 2031616000)
0277       det2B->AddHit(detlayerId, x, y, z);
0278     else if (detlayerId == 2046820352) {
0279       detToF_F->AddHit(x, y, tof);
0280       recCellId_F.push_back(detToF_F->findCellId(x, y));
0281       recTof_F.push_back(tof);
0282     } else if (detlayerId == 2063597568) {
0283       detToF_B->AddHit(x, y, tof);
0284       recCellId_B.push_back(detToF_B->findCellId(x, y));
0285       recTof_B.push_back(tof);
0286     }
0287
0288   }  //LOOP TRK
0289   //creating Stations
0290   TrkStation_F = std::make_unique<CTPPSTrkStation>(*det1F, *det2F);
0291   TrkStation_B = std::make_unique<CTPPSTrkStation>(*det1B, *det2B);
0292 }  // end function
0293
0294 void CTPPSFastTrackingProducer::Reconstruction() {
0295   theCTPPSFastTrack.clear();
0296   int Direction;
0297   Direction = 1;  //cms positive Z / forward
0298   FastReco(Direction, &*pps_stationF);
0299   Direction = -1;  //cms negative Z / backward
0300   FastReco(Direction, &*pps_stationB);
0301 }  //end Reconstruction
0302
0303 bool CTPPSFastTrackingProducer::SearchTrack(int i,
0304                                             int j,
0305                                             int Direction,
0306                                             double& xi,
0307                                             double& t,
0308                                             double& partP,
0309                                             double& pt,
0310                                             double& thx,
0311                                             double& thy,
0312                                             double& x0,
0313                                             double& y0,
0314                                             double& xt,
0315                                             double& yt,
0316                                             double& X1d,
0317                                             double& Y1d,
0318                                             double& X2d,
0319                                             double& Y2d) {
0320   // Given 1 hit in Tracker1 and 1 hit in Tracker2 try to make a track with Hector
0321   xi = 0;
0322   t = 0;
0323   partP = 0;
0324   pt = 0;
0325   x0 = 0.;
0326   y0 = 0.;
0327   xt = 0.;
0328   yt = 0.;
0329   X1d = 0.;
0330   Y1d = 0.;
0331   X2d = 0.;
0332   Y2d = 0.;
0333   CTPPSTrkDetector* det1 = nullptr;
0334   CTPPSTrkDetector* det2 = nullptr;
0335   H_RecRPObject* station = nullptr;
0336   // Separate in forward and backward stations according to direction
0337   if (Direction > 0) {
0338     det1 = &(TrkStation_F->first);
0339     det2 = &(TrkStation_F->second);
0340     station = &*pps_stationF;
0341   } else {
0342     det1 = &(TrkStation_B->first);
0343     det2 = &(TrkStation_B->second);
0344     station = &*pps_stationB;
0345   }
0346   if (det1->ppsNHits_ <= i || det2->ppsNHits_ <= j)
0347     return false;
0348   //
0349   double x1 = det1->ppsX_.at(i);
0350   double y1 = det1->ppsY_.at(i);
0351   double x2 = det2->ppsX_.at(j);
0352   double y2 = det2->ppsY_.at(j);
0353   double eloss;
0354
0355   //thx and thy are returned in microrad
0356   ReconstructArm(
0357       station, Direction * x1, y1, Direction * x2, y2, thx, thy, eloss);  // Pass the hits in the LHC ref. frame
0358   thx *= -Direction;                                                      //  invert to the CMS ref frame
0359
0360   // Protect for unphysical results
0361   if (edm::isNotFinite(eloss) || edm::isNotFinite(thx) || edm::isNotFinite(thy))
0362     return false;
0363   //
0364
0365   if (m_verbosity)
0366     LogDebug("CTPPSFastTrackingProducer::SearchTrack:") << "thx " << thx << " thy " << thy << " eloss " << eloss;
0367
0368   // Get the start point of the reconstructed track near the origin made by Hector in the CMS ref. frame
0369   x0 = -Direction * station->getX0() * um_to_cm;
0370   y0 = station->getY0() * um_to_cm;
0371   double ImpPar = sqrt(x0 * x0 + y0 * y0);
0372   if (ImpPar > fImpParcut)
0373     return false;
0374   if (eloss < 0. || eloss > fBeamEnergy)
0375     return false;
0376   //
0377   // Calculate the reconstructed track parameters
0378   double theta = sqrt(thx * thx + thy * thy) * urad;
0379   xi = eloss / fBeamEnergy;
0380   double energy = fBeamEnergy * (1. - xi);
0381   partP = sqrt(energy * energy - PPSTools::ProtonMassSQ);
0382   t = -2. * (PPSTools::ProtonMassSQ - fBeamEnergy * energy + fBeamMomentum * partP * cos(theta));
0383   pt = partP * theta;
0384   if (xi < 0. || xi > 1. || t < 0. || t > 10. || pt <= 0.) {
0385     xi = 0.;
0386     t = 0.;
0387     partP = 0.;
0388     pt = 0.;
0389     x0 = 0.;
0390     y0 = 0.;
0391     return false;  // unphysical values
0392   }
0393   //Try to include the timing detector in the track
0394   ProjectToToF(x1, y1, x2, y2, xt, yt);  // the projections is done in the CMS ref frame
0395   X1d = x1;
0396   Y1d = y1;
0397   X2d = x2;
0398   Y2d = y2;
0399   return true;
0400 }  //end  SearchTrack
0401
0402 void CTPPSFastTrackingProducer::ReconstructArm(
0403     H_RecRPObject* pps_station, double x1, double y1, double x2, double y2, double& tx, double& ty, double& eloss) {
0404   tx = 0.;
0405   ty = 0.;
0406   eloss = 0.;
0407   if (!pps_station)
0408     return;
0409   x1 *= mm_to_um;
0410   x2 *= mm_to_um;
0411   y1 *= mm_to_um;
0412   y2 *= mm_to_um;
0413   pps_station->setPositions(x1, y1, x2, y2);
0414   double energy = pps_station->getE(AM);  // dummy call needed to calculate some Hector internal parameter
0415   if (edm::isNotFinite(energy))
0416     return;
0417   tx = pps_station->getTXIP();  // change orientation to CMS
0418   ty = pps_station->getTYIP();
0419   eloss = pps_station->getE();
0420 }
0421
0422 void CTPPSFastTrackingProducer::MatchCellId(
0423     int cellId, std::vector<int> vrecCellId, std::vector<double> vrecTof, bool& match, double& recTof) {
0424   for (unsigned int i = 0; i < vrecCellId.size(); i++) {
0425     if (cellId == vrecCellId.at(i)) {
0426       match = true;
0427       recTof = vrecTof.at(i);
0428       continue;
0429     }
0430   }
0431 }
0432
0433 void CTPPSFastTrackingProducer::FastReco(int Direction, H_RecRPObject* station) {
0434   double theta = 0.;
0435   double xi, t, partP, pt, phi, x0, y0, thx, thy, xt, yt, X1d, Y1d, X2d, Y2d;
0436   CTPPSTrkDetector* Trk1 = nullptr;
0437   CTPPSTrkDetector* Trk2 = nullptr;
0438   double pos_tof = fToFInsertion * fBeamXRMS_ToF + fToFXOffset;
0439   int cellId = 0;
0440   std::vector<double> vToFCellWidth;
0441   vToFCellWidth.reserve(8);
0442   for (int i = 0; i < 8; i++) {
0443     vToFCellWidth.push_back(fToFCellWidth[i]);
0444   }
0445   CTPPSToFDetector* ToF = new CTPPSToFDetector(
0446       fToFNCellX, fToFNCellY, vToFCellWidth, fToFCellHeight, fToFPitchX, fToFPitchY, pos_tof, fTimeSigma);
0447   if (Direction > 0) {
0448     Trk1 = &(TrkStation_F->first);
0449     Trk2 = &(TrkStation_F->second);
0450   } else {
0451     Trk1 = &(TrkStation_B->first);
0452     Trk2 = &(TrkStation_B->second);
0453   }
0454   // Make a track from EVERY pair of hits combining Tracker1 and Tracker2.
0455   // The tracks may not be independent as 1 hit may belong to more than 1 track.
0456   for (int i = 0; i < (int)Trk1->ppsNHits_; i++) {
0457     for (int j = 0; j < (int)Trk2->ppsNHits_; j++) {
0458       if (SearchTrack(i, j, Direction, xi, t, partP, pt, thx, thy, x0, y0, xt, yt, X1d, Y1d, X2d, Y2d)) {
0459         // Check if the hitted timing cell matches the reconstructed track
0460         cellId = ToF->findCellId(xt, yt);
0461         double recTof = 0.;
0462         bool matchCellId = false;
0463         if (Direction > 0) {
0464           theta = sqrt(thx * thx + thy * thy) * urad;
0465           MatchCellId(cellId, recCellId_F, recTof_F, matchCellId, recTof);
0466         } else if (Direction < 0) {
0467           theta = TMath::Pi() - sqrt(thx * thx + thy * thy) * urad;
0468           MatchCellId(cellId, recCellId_B, recTof_B, matchCellId, recTof);
0469         }
0470         phi = atan2(thy, thx);  // at this point, thx is already in the cms ref. frame
0471
0472         double px = partP * sin(theta) * cos(phi);
0473         double py = partP * sin(theta) * sin(phi);
0474         double pz = partP * cos(theta);
0475         double e = sqrt(partP * partP + PPSTools::ProtonMassSQ);
0476         TLorentzVector p(px, py, pz, e);
0477         // Invert the Lorentz boost made to take into account the crossing angle during simulation
0478         if (fCrossAngleCorr) {
0479           PPSTools::LorentzBoost(p, "MC", {fCrossingAngleBeam1, fCrossingAngleBeam2, fBeamMomentum, fBeamEnergy});
0480         }
0481         //Getting the Xi and t (squared four momentum transferred) of the reconstructed track
0482         PPSTools::Get_t_and_xi(const_cast<TLorentzVector*>(&p), t, xi, {fBeamMomentum, fBeamEnergy});
0483         double pxx = p.Px();
0484         double pyy = p.Py();
0485         double pzz = p.Pz();
0486         math::XYZVector momentum(pxx, pyy, pzz);
0487         math::XYZPoint vertex(x0, y0, 0);
0488
0489         track.setp(momentum);
0490         track.setvertex(vertex);
0491         track.sett(t);
0492         track.setxi(xi);
0493         track.setx1(X1d);
0494         track.sety1(Y1d);
0495         track.setx2(X2d);
0496         track.sety2(Y2d);
0497         if (matchCellId) {
0498           track.setcellid(cellId);
0499           track.settof(recTof);
0500         } else {
0501           track.setcellid(0);
0502           track.settof(0.);
0503         }
0504         theCTPPSFastTrack.push_back(track);
0505       }
0506     }
0507   }
0508 }  //end FastReco
0509
0510 bool CTPPSFastTrackingProducer::SetBeamLine() {
0511   edm::FileInPath b1(beam1filename.c_str());
0512   edm::FileInPath b2(beam2filename.c_str());
0513   if (lengthctpps <= 0)
0514     return false;
0515   m_beamlineCTPPS1 = std::make_unique<H_BeamLine>(-1, lengthctpps + 0.1);  // (direction, length)
0516   m_beamlineCTPPS1->fill(b2.fullPath(), 1, "IP5");
0517   m_beamlineCTPPS2 = std::make_unique<H_BeamLine>(1, lengthctpps + 0.1);  //
0518   m_beamlineCTPPS2->fill(b1.fullPath(), 1, "IP5");
0519   m_beamlineCTPPS1->offsetElements(120, 0.097);
0520   m_beamlineCTPPS2->offsetElements(120, -0.097);
0521   pps_stationF = std::make_unique<H_RecRPObject>(fz_tracker1, fz_tracker2, *m_beamlineCTPPS1);
0522   pps_stationB = std::make_unique<H_RecRPObject>(fz_tracker1, fz_tracker2, *m_beamlineCTPPS2);
0523   return true;
0524 }
0525 //define this as a plug-in
0526 DEFINE_FWK_MODULE(CTPPSFastTrackingProducer);