Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2023-10-25 10:01:54 

0001 //
0002 // -*- C++ -*-
0003 // Package:    PFTracking
0004 // Class:      PFTrackTransformer
0005 //
0006 // Original Author:  Michele Pioppi
0007 // Other Author: Daniele Benedetti
0008 
0009 #include "RecoParticleFlow/PFTracking/interface/PFTrackTransformer.h"
0010 
0011 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0012 #include "FWCore/Framework/interface/ESHandle.h"
0013 
0014 #include "DataFormats/ParticleFlowReco/interface/PFCluster.h"
0015 #include "DataFormats/ParticleFlowReco/interface/PFRecTrackFwd.h"
0016 #include "DataFormats/TrackReco/interface/Track.h"
0017 
0018 #include "TrackingTools/PatternTools/interface/Trajectory.h"
0019 #include "CommonTools/BaseParticlePropagator/interface/BaseParticlePropagator.h"
0020 
0021 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
0022 // Add by Daniele
0023 #include "TrackingTools/GsfTools/interface/MultiGaussianStateTransform.h"
0024 #include "TrackingTools/GsfTools/interface/MultiGaussianState1D.h"
0025 #include "TrackingTools/GsfTools/interface/GaussianSumUtilities1D.h"
0026 #include "RecoParticleFlow/PFTracking/interface/PFGsfHelper.h"
0027 
0028 using namespace std;
0029 using namespace reco;
0030 using namespace edm;
0031 
0032 PFTrackTransformer::PFTrackTransformer(const math::XYZVector& B) : B_(B) {
0033   LogInfo("PFTrackTransformer") << "PFTrackTransformer built";
0034 
0035   onlyprop_ = false;
0036 }
0037 
0038 PFTrackTransformer::~PFTrackTransformer() {}
0039 
0040 bool PFTrackTransformer::addPoints(reco::PFRecTrack& pftrack,
0041                                    const reco::Track& track,
0042                                    const Trajectory& traj,
0043                                    bool msgwarning) const {
0044   LogDebug("PFTrackTransformer") << "Trajectory propagation started";
0045   using namespace reco;
0046   using namespace std;
0047 
0048   float PT = track.pt();
0049   float pfmass = (pftrack.algoType() == reco::PFRecTrack::KF_ELCAND) ? 0.0005 : 0.139;
0050   float pfenergy = sqrt((pfmass * pfmass) + (track.p() * track.p()));
0051   // closest approach
0052   BaseParticlePropagator theParticle = BaseParticlePropagator(
0053       RawParticle(XYZTLorentzVector(track.px(), track.py(), track.pz(), pfenergy),
0054                   XYZTLorentzVector(track.vertex().x(), track.vertex().y(), track.vertex().z(), 0.),
0055                   track.charge()),
0056       0.,
0057       0.,
0058       B_.z());
0059 
0060   float pfoutenergy = sqrt((pfmass * pfmass) + track.outerMomentum().Mag2());
0061   BaseParticlePropagator theOutParticle;
0062   if (track.outerOk())
0063     theOutParticle = BaseParticlePropagator(
0064         RawParticle(
0065             XYZTLorentzVector(
0066                 track.outerMomentum().x(), track.outerMomentum().y(), track.outerMomentum().z(), pfoutenergy),
0067             XYZTLorentzVector(track.outerPosition().x(), track.outerPosition().y(), track.outerPosition().z(), 0.),
0068             track.charge()),
0069         0.,
0070         0.,
0071         B_.z());
0072   else
0073     theOutParticle = theParticle;  // if outer state is not available, use the information at the closest approach
0074 
0075   math::XYZTLorentzVector momClosest = math::XYZTLorentzVector(track.px(), track.py(), track.pz(), track.p());
0076   const math::XYZPoint& posClosest = track.vertex();
0077 
0078   pftrack.addPoint(PFTrajectoryPoint(-1, PFTrajectoryPoint::ClosestApproach, posClosest, momClosest));
0079 
0080   //BEAMPIPE
0081   theParticle.setPropagationConditions(
0082       pfGeometry_.outerRadius(PFGeometry::BeamPipe), pfGeometry_.outerZ(PFGeometry::BeamPipe), false);
0083   theParticle.propagate();
0084   if (theParticle.getSuccess() != 0)
0085     pftrack.addPoint(PFTrajectoryPoint(-1,
0086                                        PFTrajectoryPoint::BeamPipeOrEndVertex,
0087                                        math::XYZPoint(theParticle.particle().vertex()),
0088                                        math::XYZTLorentzVector(theParticle.particle().momentum())));
0089   else {
0090     PFTrajectoryPoint dummyMaxSh;
0091     pftrack.addPoint(dummyMaxSh);
0092   }
0093 
0094   //trajectory points
0095 
0096   if (!onlyprop_) {
0097     bool direction = (traj.direction() == alongMomentum);
0098     vector<TrajectoryMeasurement> measurements = traj.measurements();
0099     int iTrajFirst = (direction) ? 0 : measurements.size() - 1;
0100     int increment = (direction) ? +1 : -1;
0101     int iTrajLast = (direction) ? int(measurements.size()) : -1;
0102 
0103     for (int iTraj = iTrajFirst; iTraj != iTrajLast; iTraj += increment) {
0104       GlobalPoint v = measurements[iTraj].updatedState().globalPosition();
0105       GlobalVector p = measurements[iTraj].updatedState().globalMomentum();
0106       unsigned int iid = measurements[iTraj].recHit()->det()->geographicalId().rawId();
0107       pftrack.addPoint(PFTrajectoryPoint(
0108           iid, -1, math::XYZPoint(v.x(), v.y(), v.z()), math::XYZTLorentzVector(p.x(), p.y(), p.z(), p.mag())));
0109     }
0110   }
0111 
0112   // ES1
0113   bool isBelowPS = false;
0114   theOutParticle.propagateToPreshowerLayer1(false);
0115   if (theOutParticle.getSuccess() != 0)
0116     pftrack.addPoint(PFTrajectoryPoint(-1,
0117                                        PFTrajectoryPoint::PS1,
0118                                        math::XYZPoint(theOutParticle.particle().vertex()),
0119                                        math::XYZTLorentzVector(theOutParticle.particle().momentum())));
0120   else {
0121     PFTrajectoryPoint dummyPS1;
0122     pftrack.addPoint(dummyPS1);
0123   }
0124 
0125   // ES1
0126   theOutParticle.propagateToPreshowerLayer2(false);
0127   if (theOutParticle.getSuccess() != 0) {
0128     pftrack.addPoint(PFTrajectoryPoint(-1,
0129                                        PFTrajectoryPoint::PS2,
0130                                        math::XYZPoint(theOutParticle.particle().vertex()),
0131                                        math::XYZTLorentzVector(theOutParticle.particle().momentum())));
0132     isBelowPS = true;
0133   } else {
0134     PFTrajectoryPoint dummyPS2;
0135     pftrack.addPoint(dummyPS2);
0136   }
0137 
0138   //ECAL entrance
0139   theOutParticle.propagateToEcalEntrance(false);
0140   if (theOutParticle.getSuccess() != 0) {
0141     pftrack.addPoint(PFTrajectoryPoint(-1,
0142                                        PFTrajectoryPoint::ECALEntrance,
0143                                        math::XYZPoint(theOutParticle.particle().vertex()),
0144                                        math::XYZTLorentzVector(theOutParticle.particle().momentum())));
0145     double ecalShowerDepth = PFCluster::getDepthCorrection(theOutParticle.particle().momentum().E(), isBelowPS, false);
0146 
0147     math::XYZPoint meanShower =
0148         math::XYZPoint(theOutParticle.particle().vertex()) +
0149         math::XYZTLorentzVector(theOutParticle.particle().momentum()).Vect().Unit() * ecalShowerDepth;
0150 
0151     pftrack.addPoint(PFTrajectoryPoint(-1,
0152                                        PFTrajectoryPoint::ECALShowerMax,
0153                                        meanShower,
0154                                        math::XYZTLorentzVector(theOutParticle.particle().momentum())));
0155   } else {
0156     if (PT > 5. && theOutParticle.particle().cos2ThetaV() < cos2ThetaV_Endcap_HiEnd_ && msgwarning)
0157       // Check consistent boundary as in CommonTools/BaseParticlePropagator/src/BaseParticlePropagator.cc
0158       // out of the HB/HE acceptance
0159       // eta = 3.0 -> cos^2(theta) = 0.99014 - cos**2(theta) is faster to determine than eta
0160       LogWarning("PFTrackTransformer") << "KF TRACK " << pftrack << " PROPAGATION TO THE ECAL HAS FAILED\n"
0161                                        << "theOutParticle.particle() pt,eta: " << theOutParticle.particle().pt() << " "
0162                                        << theOutParticle.particle().eta();
0163     PFTrajectoryPoint dummyECAL;
0164     pftrack.addPoint(dummyECAL);
0165     PFTrajectoryPoint dummyMaxSh;
0166     pftrack.addPoint(dummyMaxSh);
0167   }
0168 
0169   //HCAL entrance
0170   theOutParticle.propagateToHcalEntrance(false);
0171   if (theOutParticle.getSuccess() != 0)
0172     pftrack.addPoint(PFTrajectoryPoint(-1,
0173                                        PFTrajectoryPoint::HCALEntrance,
0174                                        math::XYZPoint(theOutParticle.particle().vertex()),
0175                                        math::XYZTLorentzVector(theOutParticle.particle().momentum())));
0176   else {
0177     if (PT > 5. && theOutParticle.particle().cos2ThetaV() < cos2ThetaV_Endcap_HiEnd_ && msgwarning)
0178       LogWarning("PFTrackTransformer") << "KF TRACK " << pftrack << " PROPAGATION TO THE HCAL ENTRANCE HAS FAILED";
0179     PFTrajectoryPoint dummyHCALentrance;
0180     pftrack.addPoint(dummyHCALentrance);
0181   }
0182 
0183   //HCAL exit
0184   // theOutParticle.setMagneticField(0); //Show we propagate as straight line inside HCAL ?
0185   theOutParticle.propagateToHcalExit(false);
0186   if (theOutParticle.getSuccess() != 0)
0187     pftrack.addPoint(PFTrajectoryPoint(-1,
0188                                        PFTrajectoryPoint::HCALExit,
0189                                        math::XYZPoint(theOutParticle.particle().vertex()),
0190                                        math::XYZTLorentzVector(theOutParticle.particle().momentum())));
0191   else {
0192     if (PT > 5. && theOutParticle.particle().cos2ThetaV() < cos2ThetaV_Endcap_HiEnd_ && msgwarning)
0193       LogWarning("PFTrackTransformer") << "KF TRACK " << pftrack << " PROPAGATION TO THE HCAL EXIT HAS FAILED";
0194     PFTrajectoryPoint dummyHCALexit;
0195     pftrack.addPoint(dummyHCALexit);
0196   }
0197 
0198   //HO layer0
0199   // if (abs(theOutParticle.particle().vertex().z())<550) {
0200   if (PT > 3.0) {  //Same value is used in PFBlockAlgo::link( case PFBlockLink::TRACKandHO:
0201     theOutParticle.setMagneticField(0);
0202     theOutParticle.propagateToHOLayer(false);
0203     if (theOutParticle.getSuccess() != 0) {
0204       pftrack.addPoint(PFTrajectoryPoint(-1,
0205                                          PFTrajectoryPoint::HOLayer,
0206                                          math::XYZPoint(theOutParticle.particle().vertex()),
0207                                          math::XYZTLorentzVector(theOutParticle.particle().momentum())));
0208     } else {
0209       if (PT > 5. && abs(theOutParticle.particle().Z()) < 700.25 && msgwarning)
0210         LogWarning("PFTrackTransformer") << "KF TRACK " << pftrack << " PROPAGATION TO THE HO HAS FAILED";
0211       PFTrajectoryPoint dummyHOLayer;
0212       pftrack.addPoint(dummyHOLayer);
0213     }
0214   } else {
0215     PFTrajectoryPoint dummyHOLayer;
0216     pftrack.addPoint(dummyHOLayer);
0217   }
0218 
0219   //VFcal(HF) entrance
0220   theOutParticle.setMagneticField(0);  // assume B=0 works ok for extrapolation to HF
0221   theOutParticle.propagateToVFcalEntrance(false);
0222   if (theOutParticle.getSuccess() != 0) {
0223     pftrack.addPoint(PFTrajectoryPoint(-1,
0224                                        PFTrajectoryPoint::VFcalEntrance,
0225                                        math::XYZPoint(theOutParticle.particle().vertex()),
0226                                        math::XYZTLorentzVector(theOutParticle.particle().momentum())));
0227   } else {
0228     PFTrajectoryPoint dummyVFcalentrance;
0229     pftrack.addPoint(dummyVFcalentrance);
0230   }
0231 
0232   return true;
0233 }
0234 bool PFTrackTransformer::addPointsAndBrems(reco::GsfPFRecTrack& pftrack,
0235                                            const reco::Track& track,
0236                                            const Trajectory& traj,
0237                                            const bool& GetMode) const {
0238   float PT = track.pt();
0239   // Trajectory for each trajectory point
0240 
0241   bool direction = (traj.direction() == alongMomentum);
0242   vector<TrajectoryMeasurement> measurements = traj.measurements();
0243   int iTrajFirst = (direction) ? 0 : measurements.size() - 1;
0244   int increment = (direction) ? +1 : -1;
0245   int iTrajLast = (direction) ? int(measurements.size()) : -1;
0246 
0247   unsigned int iTrajPos = 0;
0248   for (int iTraj = iTrajFirst; iTraj != iTrajLast; iTraj += increment) {
0249     GlobalPoint v = measurements[iTraj].updatedState().globalPosition();
0250     PFGsfHelper* PFGsf = new PFGsfHelper(measurements[iTraj]);
0251     //if (PFGsf->isValid()){
0252     bool ComputeMODE = GetMode;
0253     GlobalVector p = PFGsf->computeP(ComputeMODE);
0254     double DP = PFGsf->fittedDP();
0255     double SigmaDP = PFGsf->sigmafittedDP();
0256     unsigned int iid = measurements[iTraj].recHit()->det()->geographicalId().rawId();
0257     delete PFGsf;
0258 
0259     // --------------------------   Fill GSF Track -------------------------------------
0260 
0261     //    float pfmass= (pftrack.algoType()==reco::PFRecTrack::KF_ELCAND) ? 0.0005 : 0.139;
0262     float ptot = sqrt((p.x() * p.x()) + (p.y() * p.y()) + (p.z() * p.z()));
0263     float pfenergy = ptot;
0264 
0265     if (iTraj == iTrajFirst) {
0266       math::XYZTLorentzVector momClosest = math::XYZTLorentzVector(p.x(), p.y(), p.z(), ptot);
0267       const math::XYZPoint& posClosest = track.vertex();
0268       pftrack.addPoint(PFTrajectoryPoint(-1, PFTrajectoryPoint::ClosestApproach, posClosest, momClosest));
0269 
0270       BaseParticlePropagator theInnerParticle = BaseParticlePropagator(
0271           RawParticle(XYZTLorentzVector(p.x(), p.y(), p.z(), pfenergy),
0272                       XYZTLorentzVector(track.vertex().x(), track.vertex().y(), track.vertex().z(), 0.),
0273                       track.charge()),  //DANIELE Same thing v.x(),v.y(),v.()?
0274           0.,
0275           0.,
0276           B_.z());
0277 
0278       //BEAMPIPE
0279       theInnerParticle.setPropagationConditions(
0280           pfGeometry_.outerRadius(PFGeometry::BeamPipe), pfGeometry_.outerZ(PFGeometry::BeamPipe), false);
0281       theInnerParticle.propagate();
0282       if (theInnerParticle.getSuccess() != 0)
0283         pftrack.addPoint(PFTrajectoryPoint(-1,
0284                                            PFTrajectoryPoint::BeamPipeOrEndVertex,
0285                                            math::XYZPoint(theInnerParticle.particle().vertex()),
0286                                            math::XYZTLorentzVector(theInnerParticle.particle().momentum())));
0287       else {
0288         PFTrajectoryPoint dummyMaxSh;
0289         pftrack.addPoint(dummyMaxSh);
0290       }
0291 
0292       // First Point for the trajectory == Vertex ??
0293       pftrack.addPoint(PFTrajectoryPoint(
0294           iid, -1, math::XYZPoint(v.x(), v.y(), v.z()), math::XYZTLorentzVector(p.x(), p.y(), p.z(), p.mag())));
0295     }
0296     if (iTraj != iTrajFirst && iTraj != (abs(iTrajLast) - 1)) {
0297       pftrack.addPoint(PFTrajectoryPoint(
0298           iid, -1, math::XYZPoint(v.x(), v.y(), v.z()), math::XYZTLorentzVector(p.x(), p.y(), p.z(), p.mag())));
0299     }
0300     if (iTraj == (abs(iTrajLast) - 1)) {
0301       // Last Trajectory Meas
0302       pftrack.addPoint(PFTrajectoryPoint(
0303           iid, -1, math::XYZPoint(v.x(), v.y(), v.z()), math::XYZTLorentzVector(p.x(), p.y(), p.z(), p.mag())));
0304 
0305       BaseParticlePropagator theOutParticle =
0306           BaseParticlePropagator(RawParticle(XYZTLorentzVector(p.x(), p.y(), p.z(), pfenergy),
0307                                              XYZTLorentzVector(v.x(), v.y(), v.z(), 0.),
0308                                              track.charge()),
0309                                  0.,
0310                                  0.,
0311                                  B_.z());
0312 
0313       // ES1
0314       bool isBelowPS = false;
0315       theOutParticle.propagateToPreshowerLayer1(false);
0316       if (theOutParticle.getSuccess() != 0)
0317         pftrack.addPoint(PFTrajectoryPoint(-1,
0318                                            PFTrajectoryPoint::PS1,
0319                                            math::XYZPoint(theOutParticle.particle().vertex()),
0320                                            math::XYZTLorentzVector(theOutParticle.particle().momentum())));
0321       else {
0322         PFTrajectoryPoint dummyPS1;
0323         pftrack.addPoint(dummyPS1);
0324       }
0325 
0326       // ES2
0327       theOutParticle.propagateToPreshowerLayer2(false);
0328       if (theOutParticle.getSuccess() != 0) {
0329         pftrack.addPoint(PFTrajectoryPoint(-1,
0330                                            PFTrajectoryPoint::PS2,
0331                                            math::XYZPoint(theOutParticle.particle().vertex()),
0332                                            math::XYZTLorentzVector(theOutParticle.particle().momentum())));
0333         isBelowPS = true;
0334       } else {
0335         PFTrajectoryPoint dummyPS2;
0336         pftrack.addPoint(dummyPS2);
0337       }
0338 
0339       theOutParticle.propagateToEcalEntrance(false);
0340 
0341       if (theOutParticle.getSuccess() != 0) {
0342         pftrack.addPoint(PFTrajectoryPoint(-1,
0343                                            PFTrajectoryPoint::ECALEntrance,
0344                                            math::XYZPoint(theOutParticle.particle().vertex()),
0345                                            math::XYZTLorentzVector(theOutParticle.particle().momentum())));
0346         double ecalShowerDepth =
0347             PFCluster::getDepthCorrection(theOutParticle.particle().momentum().E(), isBelowPS, false);
0348 
0349         math::XYZPoint meanShower =
0350             math::XYZPoint(theOutParticle.particle().vertex()) +
0351             math::XYZTLorentzVector(theOutParticle.particle().momentum()).Vect().Unit() * ecalShowerDepth;
0352 
0353         pftrack.addPoint(PFTrajectoryPoint(-1,
0354                                            PFTrajectoryPoint::ECALShowerMax,
0355                                            meanShower,
0356                                            math::XYZTLorentzVector(theOutParticle.particle().momentum())));
0357       } else {
0358         if (PT > 5. && theOutParticle.particle().cos2ThetaV() < cos2ThetaV_Endcap_HiEnd_)
0359           LogWarning("PFTrackTransformer") << "GSF TRACK " << pftrack << " PROPAGATION TO THE ECAL HAS FAILED\n"
0360                                            << "theOutParticle.particle() pt,eta: " << theOutParticle.particle().pt()
0361                                            << " " << theOutParticle.particle().eta();
0362         PFTrajectoryPoint dummyECAL;
0363         pftrack.addPoint(dummyECAL);
0364         PFTrajectoryPoint dummyMaxSh;
0365         pftrack.addPoint(dummyMaxSh);
0366       }
0367 
0368       //HCAL entrance
0369       theOutParticle.propagateToHcalEntrance(false);
0370       if (theOutParticle.getSuccess() != 0)
0371         pftrack.addPoint(PFTrajectoryPoint(-1,
0372                                            PFTrajectoryPoint::HCALEntrance,
0373                                            math::XYZPoint(theOutParticle.particle().vertex()),
0374                                            math::XYZTLorentzVector(theOutParticle.particle().momentum())));
0375       else {
0376         if (PT > 5. && theOutParticle.particle().cos2ThetaV() < cos2ThetaV_Endcap_HiEnd_)
0377           LogWarning("PFTrackTransformer") << "GSF TRACK " << pftrack << " PROPAGATION TO THE HCAL ENTRANCE HAS FAILED";
0378         PFTrajectoryPoint dummyHCALentrance;
0379         pftrack.addPoint(dummyHCALentrance);
0380       }
0381       //HCAL exit
0382       theOutParticle.propagateToHcalExit(false);
0383       if (theOutParticle.getSuccess() != 0)
0384         pftrack.addPoint(PFTrajectoryPoint(-1,
0385                                            PFTrajectoryPoint::HCALExit,
0386                                            math::XYZPoint(theOutParticle.particle().vertex()),
0387                                            math::XYZTLorentzVector(theOutParticle.particle().momentum())));
0388       else {
0389         if (PT > 5. && theOutParticle.particle().cos2ThetaV() < cos2ThetaV_Endcap_HiEnd_)
0390           LogWarning("PFTrackTransformer") << "GSF TRACK " << pftrack << " PROPAGATION TO THE HCAL EXIT HAS FAILED";
0391         PFTrajectoryPoint dummyHCALexit;
0392         pftrack.addPoint(dummyHCALexit);
0393       }
0394 
0395       //HO Layer0
0396       if (abs(theOutParticle.particle().vertex().z()) < 550) {
0397         theOutParticle.setMagneticField(0);
0398         theOutParticle.propagateToHOLayer(false);
0399         if (theOutParticle.getSuccess() != 0)
0400           pftrack.addPoint(PFTrajectoryPoint(-1,
0401                                              PFTrajectoryPoint::HOLayer,
0402                                              math::XYZPoint(theOutParticle.particle().vertex()),
0403                                              math::XYZTLorentzVector(theOutParticle.particle().momentum())));
0404         else {
0405           if (PT > 5. && abs(theOutParticle.particle().Z()) < 700.25)
0406             LogWarning("PFTrackTransformer") << "GSF TRACK " << pftrack << " PROPAGATION TO THE HO HAS FAILED";
0407           PFTrajectoryPoint dummyHOLayer;
0408           pftrack.addPoint(dummyHOLayer);
0409         }
0410       }
0411     }
0412 
0413     // --------------------------   END GSF Track -------------------------------------
0414 
0415     // --------------------------   Fill Brem "Track" ---------------------------------
0416     // Fill the brem for each traj point
0417 
0418     //check that the vertex of the brem is in the tracker volume
0419     if ((v.perp() > 110) || (fabs(v.z()) > 280))
0420       continue;
0421     unsigned int iTrajPoint = iTrajPos + 2;
0422     if (iid % 2 == 1)
0423       iTrajPoint = 99;
0424 
0425     PFBrem brem(DP, SigmaDP, iTrajPoint);
0426 
0427     GlobalVector p_gamma = p * (fabs(DP) / p.mag());  // Direction from the electron (tangent), DP without any sign!;
0428     float e_gamma = fabs(DP);                         // DP = pout-pin so could be negative
0429     constexpr int gamma_charge = 0;
0430     BaseParticlePropagator theBremParticle =
0431         BaseParticlePropagator(RawParticle(XYZTLorentzVector(p_gamma.x(), p_gamma.y(), p_gamma.z(), e_gamma),
0432                                            XYZTLorentzVector(v.x(), v.y(), v.z(), 0.),
0433                                            gamma_charge),
0434                                0.,
0435                                0.,
0436                                B_.z());
0437 
0438     // add TrajectoryPoint for Brem, PS, ECAL, ECALShowMax, HCAL
0439     // Brem Entrance PS Layer1
0440 
0441     PFTrajectoryPoint dummyClosest;  // Added just to have the right number order in PFTrack.cc
0442     brem.addPoint(dummyClosest);
0443 
0444     PFTrajectoryPoint dummyBeamPipe;  // Added just to have the right number order in PFTrack.cc
0445     brem.addPoint(dummyBeamPipe);
0446 
0447     // ES1
0448     bool isBelowPS = false;
0449     theBremParticle.propagateToPreshowerLayer1(false);
0450     if (theBremParticle.getSuccess() != 0)
0451       brem.addPoint(PFTrajectoryPoint(-1,
0452                                       PFTrajectoryPoint::PS1,
0453                                       math::XYZPoint(theBremParticle.particle().vertex()),
0454                                       math::XYZTLorentzVector(theBremParticle.particle().momentum())));
0455     else {
0456       PFTrajectoryPoint dummyPS1;
0457       brem.addPoint(dummyPS1);
0458     }
0459 
0460     // ES2
0461     // Brem Entrance PS Layer 2
0462     theBremParticle.propagateToPreshowerLayer2(false);
0463     if (theBremParticle.getSuccess() != 0) {
0464       brem.addPoint(PFTrajectoryPoint(-1,
0465                                       PFTrajectoryPoint::PS2,
0466                                       math::XYZPoint(theBremParticle.particle().vertex()),
0467                                       math::XYZTLorentzVector(theBremParticle.particle().momentum())));
0468       isBelowPS = true;
0469     } else {
0470       PFTrajectoryPoint dummyPS2;
0471       brem.addPoint(dummyPS2);
0472     }
0473 
0474     //ECAL entrance
0475     theBremParticle.propagateToEcalEntrance(false);
0476     if (theBremParticle.getSuccess() != 0) {
0477       brem.addPoint(PFTrajectoryPoint(-1,
0478                                       PFTrajectoryPoint::ECALEntrance,
0479                                       math::XYZPoint(theBremParticle.particle().vertex()),
0480                                       math::XYZTLorentzVector(theBremParticle.particle().momentum())));
0481       double ecalShowerDepth =
0482           PFCluster::getDepthCorrection(theBremParticle.particle().momentum().E(), isBelowPS, false);
0483 
0484       math::XYZPoint meanShower =
0485           math::XYZPoint(theBremParticle.particle().vertex()) +
0486           math::XYZTLorentzVector(theBremParticle.particle().momentum()).Vect().Unit() * ecalShowerDepth;
0487 
0488       brem.addPoint(PFTrajectoryPoint(-1,
0489                                       PFTrajectoryPoint::ECALShowerMax,
0490                                       meanShower,
0491                                       math::XYZTLorentzVector(theBremParticle.particle().momentum())));
0492     } else {
0493       if ((DP > 5.) && ((DP / SigmaDP) > 3) && theBremParticle.particle().cos2ThetaV() < cos2ThetaV_Endcap_HiEnd_)
0494         LogWarning("PFTrackTransformer") << "BREM " << brem << " PROPAGATION TO THE ECAL HAS FAILED\n"
0495                                          << "theBremParticle.particle() pt,eta,cos2ThetaV: "
0496                                          << theBremParticle.particle().pt() << " " << theBremParticle.particle().eta()
0497                                          << " " << theBremParticle.particle().cos2ThetaV();
0498       PFTrajectoryPoint dummyECAL;
0499       brem.addPoint(dummyECAL);
0500       PFTrajectoryPoint dummyMaxSh;
0501       brem.addPoint(dummyMaxSh);
0502     }
0503 
0504     //HCAL entrance
0505     theBremParticle.propagateToHcalEntrance(false);
0506     if (theBremParticle.getSuccess() != 0)
0507       brem.addPoint(PFTrajectoryPoint(-1,
0508                                       PFTrajectoryPoint::HCALEntrance,
0509                                       math::XYZPoint(theBremParticle.particle().vertex()),
0510                                       math::XYZTLorentzVector(theBremParticle.particle().momentum())));
0511     else {
0512       if ((DP > 5.) && ((DP / SigmaDP) > 3) && theBremParticle.particle().cos2ThetaV() < cos2ThetaV_Endcap_HiEnd_)
0513         LogWarning("PFTrackTransformer") << "BREM " << brem << " PROPAGATION TO THE HCAL ENTRANCE HAS FAILED";
0514       PFTrajectoryPoint dummyHCALentrance;
0515       brem.addPoint(dummyHCALentrance);
0516     }
0517 
0518     //HCAL exit
0519     theBremParticle.propagateToHcalExit(false);
0520     if (theBremParticle.getSuccess() != 0)
0521       brem.addPoint(PFTrajectoryPoint(-1,
0522                                       PFTrajectoryPoint::HCALExit,
0523                                       math::XYZPoint(theBremParticle.particle().vertex()),
0524                                       math::XYZTLorentzVector(theBremParticle.particle().momentum())));
0525     else {
0526       if ((DP > 5.) && ((DP / SigmaDP) > 3) && theBremParticle.particle().cos2ThetaV() < cos2ThetaV_Endcap_HiEnd_)
0527         LogWarning("PFTrackTransformer") << "BREM " << brem << " PROPAGATION TO THE HCAL EXIT HAS FAILED";
0528       PFTrajectoryPoint dummyHCALexit;
0529       brem.addPoint(dummyHCALexit);
0530     }
0531 
0532     //HO Layer0
0533     if (abs(theBremParticle.particle().vertex().z()) < 550.0) {
0534       theBremParticle.setMagneticField(0);
0535       theBremParticle.propagateToHOLayer(false);
0536       if (theBremParticle.getSuccess() != 0)
0537         brem.addPoint(PFTrajectoryPoint(-1,
0538                                         PFTrajectoryPoint::HOLayer,
0539                                         math::XYZPoint(theBremParticle.particle().vertex()),
0540                                         math::XYZTLorentzVector(theBremParticle.particle().momentum())));
0541       else {
0542         if ((DP > 5.) && ((DP / SigmaDP) > 3) && abs(theBremParticle.particle().Z()) < 700.25)
0543           LogWarning("PFTrackTransformer") << "BREM " << brem << " PROPAGATION TO THE H0 HAS FAILED";
0544         PFTrajectoryPoint dummyHOLayer;
0545         brem.addPoint(dummyHOLayer);
0546       }
0547     }
0548     pftrack.addBrem(brem);
0549     iTrajPos++;
0550   }
0551   return true;
0552 }
0553 
0554 bool PFTrackTransformer::addPointsAndBrems(reco::GsfPFRecTrack& pftrack,
0555                                            const reco::GsfTrack& track,
0556                                            const MultiTrajectoryStateTransform& mtjstate) const {
0557   //  float PT= track.pt();
0558   unsigned int iTrajPos = 0;
0559   unsigned int iid = 0;  // not anymore saved
0560 
0561   // *****************************   INNER State *************************************
0562   TrajectoryStateOnSurface inTSOS = mtjstate.innerStateOnSurface((track));
0563   TrajectoryStateOnSurface outTSOS = mtjstate.outerStateOnSurface((track));
0564 
0565   if (!inTSOS.isValid() || !outTSOS.isValid()) {
0566     if (!inTSOS.isValid())
0567       LogWarning("PFTrackTransformer") << " INNER TSOS NOT VALID ";
0568     if (!outTSOS.isValid())
0569       LogWarning("PFTrackTransformer") << " OUTER TSOS NOT VALID ";
0570     return false;
0571   }
0572 
0573   GlobalVector InMom;
0574   GlobalPoint InPos;
0575   if (inTSOS.isValid()) {
0576     multiTrajectoryStateMode::momentumFromModeCartesian(inTSOS, InMom);
0577     multiTrajectoryStateMode::positionFromModeCartesian(inTSOS, InPos);
0578   } else {
0579     InMom = GlobalVector(track.pxMode(), track.pyMode(), track.pzMode());
0580     InPos = GlobalPoint(0., 0., 0.);
0581   }
0582 
0583   //  float pfmass= (pftrack.algoType()==reco::PFRecTrack::KF_ELCAND) ? 0.0005 : 0.139;
0584   float ptot = sqrt((InMom.x() * InMom.x()) + (InMom.y() * InMom.y()) + (InMom.z() * InMom.z()));
0585   float pfenergy = ptot;
0586 
0587   math::XYZTLorentzVector momClosest = math::XYZTLorentzVector(InMom.x(), InMom.y(), InMom.z(), ptot);
0588   const math::XYZPoint& posClosest = track.vertex();
0589   pftrack.addPoint(PFTrajectoryPoint(-1, PFTrajectoryPoint::ClosestApproach, posClosest, momClosest));
0590 
0591   BaseParticlePropagator theInnerParticle = BaseParticlePropagator(
0592       RawParticle(XYZTLorentzVector(InMom.x(), InMom.y(), InMom.z(), pfenergy),
0593                   XYZTLorentzVector(track.vertex().x(), track.vertex().y(), track.vertex().z(), 0.),
0594                   track.charge()),  //DANIELE Same thing v.x(),v.y(),v.()?
0595       0.,
0596       0.,
0597       B_.z());
0598   //BEAMPIPE
0599   theInnerParticle.setPropagationConditions(
0600       pfGeometry_.outerRadius(PFGeometry::BeamPipe), pfGeometry_.outerZ(PFGeometry::BeamPipe), false);
0601   theInnerParticle.propagate();
0602   if (theInnerParticle.getSuccess() != 0)
0603     pftrack.addPoint(PFTrajectoryPoint(-1,
0604                                        PFTrajectoryPoint::BeamPipeOrEndVertex,
0605                                        math::XYZPoint(theInnerParticle.particle().vertex()),
0606                                        math::XYZTLorentzVector(theInnerParticle.particle().momentum())));
0607   else {
0608     PFTrajectoryPoint dummyBeam;
0609     pftrack.addPoint(dummyBeam);
0610   }
0611 
0612   // first tjpoint
0613   pftrack.addPoint(PFTrajectoryPoint(iid,
0614                                      -1,
0615                                      math::XYZPoint(InPos.x(), InPos.y(), InPos.z()),
0616                                      math::XYZTLorentzVector(InMom.x(), InMom.y(), InMom.z(), InMom.mag())));
0617 
0618   //######### Photon at INNER State ##########
0619 
0620   unsigned int iTrajPoint = iTrajPos + 2;
0621   double dp_tang = ptot;
0622   double sdp_tang = track.ptModeError() * (track.pMode() / track.ptMode());
0623   PFBrem brem(dp_tang, sdp_tang, iTrajPoint);
0624   constexpr int gamma_charge = 0;
0625   BaseParticlePropagator theBremParticle =
0626       BaseParticlePropagator(RawParticle(XYZTLorentzVector(InMom.x(), InMom.y(), InMom.z(), dp_tang),
0627                                          XYZTLorentzVector(InPos.x(), InPos.y(), InPos.z(), 0.),
0628                                          gamma_charge),
0629                              0.,
0630                              0.,
0631                              B_.z());
0632   // add TrajectoryPoint for Brem, PS, ECAL, ECALShowMax, HCAL
0633   // Brem Entrance PS Layer1
0634   PFTrajectoryPoint dummyClosest;  // Added just to have the right number order in PFTrack.cc
0635   brem.addPoint(dummyClosest);
0636 
0637   PFTrajectoryPoint dummyBeamPipe;  // Added just to have the right number order in PFTrack.cc
0638   brem.addPoint(dummyBeamPipe);
0639 
0640   bool isBelowPS = false;
0641   theBremParticle.propagateToPreshowerLayer1(false);
0642   if (theBremParticle.getSuccess() != 0)
0643     brem.addPoint(PFTrajectoryPoint(-1,
0644                                     PFTrajectoryPoint::PS1,
0645                                     math::XYZPoint(theBremParticle.particle().vertex()),
0646                                     math::XYZTLorentzVector(theBremParticle.particle().momentum())));
0647   else {
0648     PFTrajectoryPoint dummyPS1;
0649     brem.addPoint(dummyPS1);
0650   }
0651 
0652   // Brem Entrance PS Layer 2
0653 
0654   theBremParticle.propagateToPreshowerLayer2(false);
0655   if (theBremParticle.getSuccess() != 0) {
0656     brem.addPoint(PFTrajectoryPoint(-1,
0657                                     PFTrajectoryPoint::PS2,
0658                                     math::XYZPoint(theBremParticle.particle().vertex()),
0659                                     math::XYZTLorentzVector(theBremParticle.particle().momentum())));
0660     isBelowPS = true;
0661   } else {
0662     PFTrajectoryPoint dummyPS2;
0663     brem.addPoint(dummyPS2);
0664   }
0665 
0666   //ECAL entrance
0667   theBremParticle.propagateToEcalEntrance(false);
0668   if (theBremParticle.getSuccess() != 0) {
0669     brem.addPoint(PFTrajectoryPoint(-1,
0670                                     PFTrajectoryPoint::ECALEntrance,
0671                                     math::XYZPoint(theBremParticle.particle().vertex()),
0672                                     math::XYZTLorentzVector(theBremParticle.particle().momentum())));
0673 
0674     //  for the first brem give a low default DP of 100 MeV.
0675     double EDepthCorr = 0.01;
0676     double ecalShowerDepth = PFCluster::getDepthCorrection(EDepthCorr, isBelowPS, false);
0677 
0678     math::XYZPoint meanShower =
0679         math::XYZPoint(theBremParticle.particle().vertex()) +
0680         math::XYZTLorentzVector(theBremParticle.particle().momentum()).Vect().Unit() * ecalShowerDepth;
0681 
0682     brem.addPoint(PFTrajectoryPoint(-1,
0683                                     PFTrajectoryPoint::ECALShowerMax,
0684                                     meanShower,
0685                                     math::XYZTLorentzVector(theBremParticle.particle().momentum())));
0686   } else {
0687     if ((dp_tang > 5.) && ((dp_tang / sdp_tang) > 3) &&
0688         theBremParticle.particle().cos2ThetaV() < cos2ThetaV_Endcap_HiEnd_)
0689       LogWarning("PFTrackTransformer") << "BREM " << brem << " PROPAGATION TO THE ECAL HAS FAILED\n"
0690                                        << "theBremParticle.particle() pt,eta,cos2ThetaV: "
0691                                        << theBremParticle.particle().pt() << " " << theBremParticle.particle().eta()
0692                                        << " " << theBremParticle.particle().cos2ThetaV();
0693     PFTrajectoryPoint dummyECAL;
0694     brem.addPoint(dummyECAL);
0695     PFTrajectoryPoint dummyMaxSh;
0696     brem.addPoint(dummyMaxSh);
0697   }
0698 
0699   //HCAL entrance
0700   theBremParticle.propagateToHcalEntrance(false);
0701   if (theBremParticle.getSuccess() != 0)
0702     brem.addPoint(PFTrajectoryPoint(-1,
0703                                     PFTrajectoryPoint::HCALEntrance,
0704                                     math::XYZPoint(theBremParticle.particle().vertex()),
0705                                     math::XYZTLorentzVector(theBremParticle.particle().momentum())));
0706   else {
0707     if ((dp_tang > 5.) && ((dp_tang / sdp_tang) > 3) &&
0708         theBremParticle.particle().cos2ThetaV() < cos2ThetaV_Endcap_HiEnd_)
0709       LogWarning("PFTrackTransformer") << "BREM " << brem << " PROPAGATION TO THE HCAL ENTRANCE HAS FAILED";
0710     PFTrajectoryPoint dummyHCALentrance;
0711     brem.addPoint(dummyHCALentrance);
0712   }
0713 
0714   //HCAL exit
0715   theBremParticle.propagateToHcalExit(false);
0716   if (theBremParticle.getSuccess() != 0)
0717     brem.addPoint(PFTrajectoryPoint(-1,
0718                                     PFTrajectoryPoint::HCALExit,
0719                                     math::XYZPoint(theBremParticle.particle().vertex()),
0720                                     math::XYZTLorentzVector(theBremParticle.particle().momentum())));
0721   else {
0722     if ((dp_tang > 5.) && ((dp_tang / sdp_tang) > 3) &&
0723         theBremParticle.particle().cos2ThetaV() < cos2ThetaV_Endcap_HiEnd_)
0724       LogWarning("PFTrackTransformer") << "BREM " << brem << " PROPAGATION TO THE HCAL EXIT HAS FAILED";
0725     PFTrajectoryPoint dummyHCALexit;
0726     brem.addPoint(dummyHCALexit);
0727   }
0728 
0729   //HO Layer0
0730   if (abs(theBremParticle.particle().vertex().z()) < 550) {
0731     theBremParticle.setMagneticField(0);
0732     theBremParticle.propagateToHOLayer(false);
0733     if (theBremParticle.getSuccess() != 0)
0734       brem.addPoint(PFTrajectoryPoint(-1,
0735                                       PFTrajectoryPoint::HOLayer,
0736                                       math::XYZPoint(theBremParticle.particle().vertex()),
0737                                       math::XYZTLorentzVector(theBremParticle.particle().momentum())));
0738     else {
0739       if ((dp_tang > 5.) && ((dp_tang / sdp_tang) > 3) && abs(theBremParticle.particle().Z()) < 700.25)
0740         LogWarning("PFTrackTransformer") << "BREM " << brem << " PROPAGATION TO THE H0 HAS FAILED";
0741       PFTrajectoryPoint dummyHOLayer;
0742       brem.addPoint(dummyHOLayer);
0743     }
0744   }
0745 
0746   pftrack.addBrem(brem);
0747   iTrajPos++;
0748 
0749   // *****************************   INTERMIDIATE State *************************************
0750   //From the new Wolfgang code
0751 
0752   // To think if the cout should be removed.
0753   if (track.gsfExtra()->tangentsSize() == 0)
0754     LogError("PFTrackTransformer")
0755         << "BE CAREFUL: Gsf Tangents not stored in the event. You need to re-reco the particle-flow with "
0756            "RecoToDisplay_cfg.py and not RecoToDisplay_NoTracking_cfg.py ";
0757 
0758   vector<GsfTangent> gsftang = track.gsfExtra()->tangents();
0759   for (unsigned int iTang = 0; iTang < track.gsfExtra()->tangentsSize(); iTang++) {
0760     dp_tang = gsftang[iTang].deltaP().value();
0761     sdp_tang = gsftang[iTang].deltaP().error();
0762 
0763     //check that the vertex of the brem is in the tracker volume
0764     if ((sqrt(gsftang[iTang].position().x() * gsftang[iTang].position().x() +
0765               gsftang[iTang].position().y() * gsftang[iTang].position().y()) > 110) ||
0766         (fabs(gsftang[iTang].position().z()) > 280))
0767       continue;
0768 
0769     iTrajPoint = iTrajPos + 2;
0770     PFBrem brem(dp_tang, sdp_tang, iTrajPoint);
0771 
0772     GlobalVector p_tang =
0773         GlobalVector(gsftang[iTang].momentum().x(), gsftang[iTang].momentum().y(), gsftang[iTang].momentum().z());
0774 
0775     // ###### track tj points
0776     pftrack.addPoint(PFTrajectoryPoint(
0777         iid,
0778         -1,
0779         math::XYZPoint(gsftang[iTang].position().x(), gsftang[iTang].position().y(), gsftang[iTang].position().z()),
0780         math::XYZTLorentzVector(p_tang.x(), p_tang.y(), p_tang.z(), p_tang.mag())));
0781 
0782     //rescale
0783     GlobalVector p_gamma = p_tang * (fabs(dp_tang) / p_tang.mag());
0784 
0785     // GlobalVector
0786 
0787     double e_gamma = fabs(dp_tang);  // DP = pout-pin so could be negative
0788     theBremParticle = BaseParticlePropagator(
0789         RawParticle(
0790             XYZTLorentzVector(p_gamma.x(), p_gamma.y(), p_gamma.z(), e_gamma),
0791             XYZTLorentzVector(
0792                 gsftang[iTang].position().x(), gsftang[iTang].position().y(), gsftang[iTang].position().z(), 0.),
0793             gamma_charge),
0794         0.,
0795         0.,
0796         B_.z());
0797 
0798     PFTrajectoryPoint dummyClosest;  // Added just to have the right number order in PFTrack.cc
0799     brem.addPoint(dummyClosest);
0800 
0801     PFTrajectoryPoint dummyBeamPipe;  // Added just to have the right number order in PFTrack.cc
0802     brem.addPoint(dummyBeamPipe);
0803 
0804     isBelowPS = false;
0805     theBremParticle.propagateToPreshowerLayer1(false);
0806     if (theBremParticle.getSuccess() != 0)
0807       brem.addPoint(PFTrajectoryPoint(-1,
0808                                       PFTrajectoryPoint::PS1,
0809                                       math::XYZPoint(theBremParticle.particle().vertex()),
0810                                       math::XYZTLorentzVector(theBremParticle.particle().momentum())));
0811     else {
0812       PFTrajectoryPoint dummyPS1;
0813       brem.addPoint(dummyPS1);
0814     }
0815 
0816     // Brem Entrance PS Layer 2
0817 
0818     theBremParticle.propagateToPreshowerLayer2(false);
0819     if (theBremParticle.getSuccess() != 0) {
0820       brem.addPoint(PFTrajectoryPoint(-1,
0821                                       PFTrajectoryPoint::PS2,
0822                                       math::XYZPoint(theBremParticle.particle().vertex()),
0823                                       math::XYZTLorentzVector(theBremParticle.particle().momentum())));
0824       isBelowPS = true;
0825     } else {
0826       PFTrajectoryPoint dummyPS2;
0827       brem.addPoint(dummyPS2);
0828     }
0829 
0830     //ECAL entrance
0831     theBremParticle.propagateToEcalEntrance(false);
0832     if (theBremParticle.getSuccess() != 0) {
0833       brem.addPoint(PFTrajectoryPoint(-1,
0834                                       PFTrajectoryPoint::ECALEntrance,
0835                                       math::XYZPoint(theBremParticle.particle().vertex()),
0836                                       math::XYZTLorentzVector(theBremParticle.particle().momentum())));
0837 
0838       double ecalShowerDepth =
0839           PFCluster::getDepthCorrection(theBremParticle.particle().momentum().E(), isBelowPS, false);
0840 
0841       math::XYZPoint meanShower =
0842           math::XYZPoint(theBremParticle.particle().vertex()) +
0843           math::XYZTLorentzVector(theBremParticle.particle().momentum()).Vect().Unit() * ecalShowerDepth;
0844 
0845       brem.addPoint(PFTrajectoryPoint(-1,
0846                                       PFTrajectoryPoint::ECALShowerMax,
0847                                       meanShower,
0848                                       math::XYZTLorentzVector(theBremParticle.particle().momentum())));
0849     } else {
0850       if ((dp_tang > 5.) && ((dp_tang / sdp_tang) > 3) &&
0851           theBremParticle.particle().cos2ThetaV() < cos2ThetaV_Endcap_HiEnd_)
0852         LogWarning("PFTrackTransformer") << "BREM " << brem << " PROPAGATION TO THE ECAL HAS FAILED\n"
0853                                          << "theBremParticle.particle() pt,eta,cos2ThetaV: "
0854                                          << theBremParticle.particle().pt() << " " << theBremParticle.particle().eta()
0855                                          << " " << theBremParticle.particle().cos2ThetaV();
0856       PFTrajectoryPoint dummyECAL;
0857       brem.addPoint(dummyECAL);
0858       PFTrajectoryPoint dummyMaxSh;
0859       brem.addPoint(dummyMaxSh);
0860     }
0861 
0862     //HCAL entrance
0863     theBremParticle.propagateToHcalEntrance(false);
0864     if (theBremParticle.getSuccess() != 0)
0865       brem.addPoint(PFTrajectoryPoint(-1,
0866                                       PFTrajectoryPoint::HCALEntrance,
0867                                       math::XYZPoint(theBremParticle.particle().vertex()),
0868                                       math::XYZTLorentzVector(theBremParticle.particle().momentum())));
0869     else {
0870       if ((dp_tang > 5.) && ((dp_tang / sdp_tang) > 3) &&
0871           theBremParticle.particle().cos2ThetaV() < cos2ThetaV_Endcap_HiEnd_)
0872         LogWarning("PFTrackTransformer") << "BREM " << brem << " PROPAGATION TO THE HCAL ENTRANCE HAS FAILED";
0873       PFTrajectoryPoint dummyHCALentrance;
0874       brem.addPoint(dummyHCALentrance);
0875     }
0876 
0877     //HCAL exit
0878     theBremParticle.propagateToHcalExit(false);
0879     if (theBremParticle.getSuccess() != 0)
0880       brem.addPoint(PFTrajectoryPoint(-1,
0881                                       PFTrajectoryPoint::HCALExit,
0882                                       math::XYZPoint(theBremParticle.particle().vertex()),
0883                                       math::XYZTLorentzVector(theBremParticle.particle().momentum())));
0884     else {
0885       if ((dp_tang > 5.) && ((dp_tang / sdp_tang) > 3) &&
0886           theBremParticle.particle().cos2ThetaV() < cos2ThetaV_Endcap_HiEnd_)
0887         LogWarning("PFTrackTransformer") << "BREM " << brem << " PROPAGATION TO THE HCAL EXIT HAS FAILED";
0888       PFTrajectoryPoint dummyHCALexit;
0889       brem.addPoint(dummyHCALexit);
0890     }
0891 
0892     //HO Layer0
0893     if (abs(theBremParticle.particle().vertex().z()) < 550) {
0894       theBremParticle.setMagneticField(0);
0895       theBremParticle.propagateToHOLayer(false);
0896       if (theBremParticle.getSuccess() != 0)
0897         brem.addPoint(PFTrajectoryPoint(-1,
0898                                         PFTrajectoryPoint::HOLayer,
0899                                         math::XYZPoint(theBremParticle.particle().vertex()),
0900                                         math::XYZTLorentzVector(theBremParticle.particle().momentum())));
0901       else {
0902         if ((dp_tang > 5.) && ((dp_tang / sdp_tang) > 3) && abs(theBremParticle.particle().Z()) < 700.25)
0903           LogWarning("PFTrackTransformer") << "BREM " << brem << " PROPAGATION TO THE H0 HAS FAILED";
0904         PFTrajectoryPoint dummyHOLayer;
0905         brem.addPoint(dummyHOLayer);
0906       }
0907     }
0908 
0909     pftrack.addBrem(brem);
0910     iTrajPos++;
0911   }
0912 
0913   // *****************************   OUTER State *************************************
0914 
0915   if (outTSOS.isValid()) {
0916     GlobalVector OutMom;
0917     GlobalPoint OutPos;
0918 
0919     // DANIELE ?????  if the out is not valid maybe take the last tangent?
0920     // From Wolfgang. It should be always valid
0921 
0922     multiTrajectoryStateMode::momentumFromModeCartesian(outTSOS, OutMom);
0923     multiTrajectoryStateMode::positionFromModeCartesian(outTSOS, OutPos);
0924 
0925     // last tjpoint
0926     pftrack.addPoint(PFTrajectoryPoint(iid,
0927                                        -1,
0928                                        math::XYZPoint(OutPos.x(), OutPos.y(), OutPos.z()),
0929                                        math::XYZTLorentzVector(OutMom.x(), OutMom.y(), OutMom.z(), OutMom.mag())));
0930 
0931     float ptot_out = sqrt((OutMom.x() * OutMom.x()) + (OutMom.y() * OutMom.y()) + (OutMom.z() * OutMom.z()));
0932     float pTtot_out = sqrt((OutMom.x() * OutMom.x()) + (OutMom.y() * OutMom.y()));
0933     float pfenergy_out = ptot_out;
0934     BaseParticlePropagator theOutParticle =
0935         BaseParticlePropagator(RawParticle(XYZTLorentzVector(OutMom.x(), OutMom.y(), OutMom.z(), pfenergy_out),
0936                                            XYZTLorentzVector(OutPos.x(), OutPos.y(), OutPos.z(), 0.),
0937                                            track.charge()),
0938                                0.,
0939                                0.,
0940                                B_.z());
0941     isBelowPS = false;
0942     theOutParticle.propagateToPreshowerLayer1(false);
0943     if (theOutParticle.getSuccess() != 0)
0944       pftrack.addPoint(PFTrajectoryPoint(-1,
0945                                          PFTrajectoryPoint::PS1,
0946                                          math::XYZPoint(theOutParticle.particle().vertex()),
0947                                          math::XYZTLorentzVector(theOutParticle.particle().momentum())));
0948     else {
0949       PFTrajectoryPoint dummyPS1;
0950       pftrack.addPoint(dummyPS1);
0951     }
0952 
0953     theOutParticle.propagateToPreshowerLayer2(false);
0954     if (theOutParticle.getSuccess() != 0) {
0955       pftrack.addPoint(PFTrajectoryPoint(-1,
0956                                          PFTrajectoryPoint::PS2,
0957                                          math::XYZPoint(theOutParticle.particle().vertex()),
0958                                          math::XYZTLorentzVector(theOutParticle.particle().momentum())));
0959       isBelowPS = true;
0960     } else {
0961       PFTrajectoryPoint dummyPS2;
0962       pftrack.addPoint(dummyPS2);
0963     }
0964 
0965     //ECAL entrance
0966     theOutParticle.propagateToEcalEntrance(false);
0967     if (theOutParticle.getSuccess() != 0) {
0968       pftrack.addPoint(PFTrajectoryPoint(-1,
0969                                          PFTrajectoryPoint::ECALEntrance,
0970                                          math::XYZPoint(theOutParticle.particle().vertex()),
0971                                          math::XYZTLorentzVector(theOutParticle.particle().momentum())));
0972       double EDepthCorr = 0.01;
0973       double ecalShowerDepth = PFCluster::getDepthCorrection(EDepthCorr, isBelowPS, false);
0974 
0975       math::XYZPoint meanShower =
0976           math::XYZPoint(theOutParticle.particle().vertex()) +
0977           math::XYZTLorentzVector(theOutParticle.particle().momentum()).Vect().Unit() * ecalShowerDepth;
0978 
0979       pftrack.addPoint(PFTrajectoryPoint(-1,
0980                                          PFTrajectoryPoint::ECALShowerMax,
0981                                          meanShower,
0982                                          math::XYZTLorentzVector(theOutParticle.particle().momentum())));
0983     } else {
0984       if (pTtot_out > 5. && theOutParticle.particle().cos2ThetaV() < cos2ThetaV_Endcap_HiEnd_)
0985         LogWarning("PFTrackTransformer") << "GSF TRACK " << pftrack << " PROPAGATION TO THE ECAL HAS FAILED\n"
0986                                          << "theOutParticle.particle() pt,eta: " << theOutParticle.particle().pt()
0987                                          << " " << theOutParticle.particle().eta();
0988 
0989       PFTrajectoryPoint dummyECAL;
0990       pftrack.addPoint(dummyECAL);
0991       PFTrajectoryPoint dummyMaxSh;
0992       pftrack.addPoint(dummyMaxSh);
0993     }
0994 
0995     //HCAL entrance
0996     theOutParticle.propagateToHcalEntrance(false);
0997     if (theOutParticle.getSuccess() != 0)
0998       pftrack.addPoint(PFTrajectoryPoint(-1,
0999                                          PFTrajectoryPoint::HCALEntrance,
1000                                          math::XYZPoint(theOutParticle.particle().vertex()),
1001                                          math::XYZTLorentzVector(theOutParticle.particle().momentum())));
1002     else {
1003       if (pTtot_out > 5. && theOutParticle.particle().cos2ThetaV() < cos2ThetaV_Endcap_HiEnd_)
1004         LogWarning("PFTrackTransformer") << "GSF TRACK " << pftrack << " PROPAGATION TO THE HCAL ENTRANCE HAS FAILED";
1005       PFTrajectoryPoint dummyHCALentrance;
1006       pftrack.addPoint(dummyHCALentrance);
1007     }
1008     //HCAL exit
1009     theOutParticle.propagateToHcalExit(false);
1010     if (theOutParticle.getSuccess() != 0)
1011       pftrack.addPoint(PFTrajectoryPoint(-1,
1012                                          PFTrajectoryPoint::HCALExit,
1013                                          math::XYZPoint(theOutParticle.particle().vertex()),
1014                                          math::XYZTLorentzVector(theOutParticle.particle().momentum())));
1015     else {
1016       if (pTtot_out > 5. && theOutParticle.particle().cos2ThetaV() < cos2ThetaV_Endcap_HiEnd_)
1017         LogWarning("PFTrackTransformer") << "GSF TRACK " << pftrack << " PROPAGATION TO THE HCAL EXIT HAS FAILED";
1018       PFTrajectoryPoint dummyHCALexit;
1019       pftrack.addPoint(dummyHCALexit);
1020     }
1021 
1022     //HO Layer0
1023     if (abs(theOutParticle.particle().vertex().z()) < 550) {
1024       theOutParticle.setMagneticField(0);
1025       theOutParticle.propagateToHOLayer(false);
1026       if (theOutParticle.getSuccess() != 0)
1027         pftrack.addPoint(PFTrajectoryPoint(-1,
1028                                            PFTrajectoryPoint::HOLayer,
1029                                            math::XYZPoint(theOutParticle.particle().vertex()),
1030                                            math::XYZTLorentzVector(theOutParticle.particle().momentum())));
1031       else {
1032         if (pTtot_out > 5. && abs(theOutParticle.particle().Z()) < 700.25)
1033           LogWarning("PFTrackTransformer") << "GSF TRACK " << pftrack << " PROPAGATION TO THE HO HAS FAILED";
1034         PFTrajectoryPoint dummyHOLayer;
1035         pftrack.addPoint(dummyHOLayer);
1036       }
1037     }
1038     //######## Photon at the OUTER State ##########
1039 
1040     dp_tang = OutMom.mag();
1041     // for the moment same inner error just for semplicity
1042     sdp_tang = track.ptModeError() * (track.pMode() / track.ptMode());
1043     iTrajPoint = iTrajPos + 2;
1044     PFBrem brem(dp_tang, sdp_tang, iTrajPoint);
1045 
1046     theBremParticle = BaseParticlePropagator(RawParticle(XYZTLorentzVector(OutMom.x(), OutMom.y(), OutMom.z(), dp_tang),
1047                                                          XYZTLorentzVector(OutPos.x(), OutPos.y(), OutPos.z(), 0.),
1048                                                          gamma_charge),
1049                                              0.,
1050                                              0.,
1051                                              B_.z());
1052 
1053     PFTrajectoryPoint dummyClosest;  // Added just to have the right number order in PFTrack.cc
1054     brem.addPoint(dummyClosest);
1055 
1056     PFTrajectoryPoint dummyBeamPipe;  // Added just to have the right number order in PFTrack.cc
1057     brem.addPoint(dummyBeamPipe);
1058 
1059     isBelowPS = false;
1060     theBremParticle.propagateToPreshowerLayer1(false);
1061     if (theBremParticle.getSuccess() != 0)
1062       brem.addPoint(PFTrajectoryPoint(-1,
1063                                       PFTrajectoryPoint::PS1,
1064                                       math::XYZPoint(theBremParticle.particle().vertex()),
1065                                       math::XYZTLorentzVector(theBremParticle.particle().momentum())));
1066     else {
1067       PFTrajectoryPoint dummyPS1;
1068       brem.addPoint(dummyPS1);
1069     }
1070 
1071     // Brem Entrance PS Layer 2
1072 
1073     theBremParticle.propagateToPreshowerLayer2(false);
1074     if (theBremParticle.getSuccess() != 0) {
1075       brem.addPoint(PFTrajectoryPoint(-1,
1076                                       PFTrajectoryPoint::PS2,
1077                                       math::XYZPoint(theBremParticle.particle().vertex()),
1078                                       math::XYZTLorentzVector(theBremParticle.particle().momentum())));
1079       isBelowPS = true;
1080     } else {
1081       PFTrajectoryPoint dummyPS2;
1082       brem.addPoint(dummyPS2);
1083     }
1084 
1085     //ECAL entrance
1086     theBremParticle.propagateToEcalEntrance(false);
1087     if (theBremParticle.getSuccess() != 0) {
1088       brem.addPoint(PFTrajectoryPoint(-1,
1089                                       PFTrajectoryPoint::ECALEntrance,
1090                                       math::XYZPoint(theBremParticle.particle().vertex()),
1091                                       math::XYZTLorentzVector(theBremParticle.particle().momentum())));
1092       double ecalShowerDepth =
1093           PFCluster::getDepthCorrection(theBremParticle.particle().momentum().E(), isBelowPS, false);
1094 
1095       math::XYZPoint meanShower =
1096           math::XYZPoint(theBremParticle.particle().vertex()) +
1097           math::XYZTLorentzVector(theBremParticle.particle().momentum()).Vect().Unit() * ecalShowerDepth;
1098 
1099       brem.addPoint(PFTrajectoryPoint(-1,
1100                                       PFTrajectoryPoint::ECALShowerMax,
1101                                       meanShower,
1102                                       math::XYZTLorentzVector(theBremParticle.particle().momentum())));
1103     } else {
1104       if ((dp_tang > 5.) && ((dp_tang / sdp_tang) > 3) &&
1105           theBremParticle.particle().cos2ThetaV() < cos2ThetaV_Endcap_HiEnd_)
1106         LogWarning("PFTrackTransformer") << "BREM " << brem << " PROPAGATION TO THE ECAL HAS FAILED\n"
1107                                          << "theBremParticle.particle() pt,eta,cos2ThetaV: "
1108                                          << theBremParticle.particle().pt() << " " << theBremParticle.particle().eta()
1109                                          << " " << theBremParticle.particle().cos2ThetaV();
1110       PFTrajectoryPoint dummyECAL;
1111       brem.addPoint(dummyECAL);
1112       PFTrajectoryPoint dummyMaxSh;
1113       brem.addPoint(dummyMaxSh);
1114     }
1115 
1116     //HCAL entrance
1117     theBremParticle.propagateToHcalEntrance(false);
1118     if (theBremParticle.getSuccess() != 0)
1119       brem.addPoint(PFTrajectoryPoint(-1,
1120                                       PFTrajectoryPoint::HCALEntrance,
1121                                       math::XYZPoint(theBremParticle.particle().vertex()),
1122                                       math::XYZTLorentzVector(theBremParticle.particle().momentum())));
1123     else {
1124       if ((dp_tang > 5.) && ((dp_tang / sdp_tang) > 3) &&
1125           theBremParticle.particle().cos2ThetaV() < cos2ThetaV_Endcap_HiEnd_)
1126         LogWarning("PFTrackTransformer") << "BREM " << brem << " PROPAGATION TO THE HCAL ENTRANCE HAS FAILED";
1127       PFTrajectoryPoint dummyHCALentrance;
1128       brem.addPoint(dummyHCALentrance);
1129     }
1130     //HCAL exit
1131     theBremParticle.propagateToHcalExit(false);
1132     if (theBremParticle.getSuccess() != 0)
1133       brem.addPoint(PFTrajectoryPoint(-1,
1134                                       PFTrajectoryPoint::HCALExit,
1135                                       math::XYZPoint(theBremParticle.particle().vertex()),
1136                                       math::XYZTLorentzVector(theBremParticle.particle().momentum())));
1137     else {
1138       if ((dp_tang > 5.) && ((dp_tang / sdp_tang) > 3) &&
1139           theBremParticle.particle().cos2ThetaV() < cos2ThetaV_Endcap_HiEnd_)
1140         LogWarning("PFTrackTransformer") << "BREM " << brem << " PROPAGATION TO THE HCAL EXIT HAS FAILED";
1141       PFTrajectoryPoint dummyHCALexit;
1142       brem.addPoint(dummyHCALexit);
1143     }
1144 
1145     //HO Layer0
1146     if (abs(theBremParticle.particle().vertex().z()) < 550) {
1147       theBremParticle.setMagneticField(0);
1148       theBremParticle.propagateToHOLayer(false);
1149       if (theBremParticle.getSuccess() != 0)
1150         brem.addPoint(PFTrajectoryPoint(-1,
1151                                         PFTrajectoryPoint::HOLayer,
1152                                         math::XYZPoint(theBremParticle.particle().vertex()),
1153                                         math::XYZTLorentzVector(theBremParticle.particle().momentum())));
1154       else {
1155         if ((dp_tang > 5.) && ((dp_tang / sdp_tang) > 3) && abs(theBremParticle.particle().Z()) < 700.25)
1156           LogWarning("PFTrackTransformer") << "BREM " << brem << " PROPAGATION TO THE H0 HAS FAILED";
1157         PFTrajectoryPoint dummyHOLayer;
1158         brem.addPoint(dummyHOLayer);
1159       }
1160     }
1161     pftrack.addBrem(brem);
1162     iTrajPos++;
1163   }
1164 
1165   return true;
1166 }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!