Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:28:45

0001 #include "RecoTracker/SpecialSeedGenerators/interface/SeedGeneratorForCosmics.h"
0002 #include "RecoTracker/TkHitPairs/interface/CosmicLayerPairs.h"
0003 #include "RecoTracker/PixelSeeding/interface/CosmicLayerTriplets.h"
0004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0005 #include "FWCore/Utilities/interface/isFinite.h"
0006 #include "TrackingTools/Records/interface/TransientRecHitRecord.h"
0007 #include "RecoTracker/TkSeedGenerator/interface/FastHelix.h"
0008 #include "RecoTracker/Record/interface/TrackerRecoGeometryRecord.h"
0009 void SeedGeneratorForCosmics::init(const SiStripRecHit2DCollection& collstereo,
0010                                    const SiStripRecHit2DCollection& collrphi,
0011                                    const SiStripMatchedRecHit2DCollection& collmatched,
0012                                    const edm::EventSetup& iSetup) {
0013   magfield = iSetup.getHandle(theMagfieldToken);
0014   tracker = iSetup.getHandle(theTrackerToken);
0015   thePropagatorAl = new PropagatorWithMaterial(alongMomentum, 0.1057, &(*magfield));
0016   thePropagatorOp = new PropagatorWithMaterial(oppositeToMomentum, 0.1057, &(*magfield));
0017   theUpdator = new KFUpdator();
0018 
0019   LogDebug("CosmicSeedFinder") << " Hits built with  " << hitsforseeds << " hits";
0020 
0021   GeometricSearchTracker const& track = iSetup.getData(theSearchTrackerToken);
0022   TrackerTopology const& ttopo = iSetup.getData(theTTopoToken);
0023 
0024   CosmicLayerPairs cosmiclayers(geometry, collrphi, collmatched, track, ttopo);
0025 
0026   thePairGenerator = new CosmicHitPairGenerator(cosmiclayers, *tracker);
0027   HitPairs.clear();
0028   if ((hitsforseeds == "pairs") || (hitsforseeds == "pairsandtriplets")) {
0029     thePairGenerator->hitPairs(region, HitPairs);
0030   }
0031 
0032   CosmicLayerTriplets cosmiclayers2(geometry, collrphi, track, ttopo);
0033   theTripletGenerator = new CosmicHitTripletGenerator(cosmiclayers2, *tracker);
0034   HitTriplets.clear();
0035   if ((hitsforseeds == "triplets") || (hitsforseeds == "pairsandtriplets")) {
0036     theTripletGenerator->hitTriplets(region, HitTriplets);
0037   }
0038 }
0039 
0040 SeedGeneratorForCosmics::SeedGeneratorForCosmics(edm::ParameterSet const& conf, edm::ConsumesCollector iCC)
0041     : maxSeeds_(conf.getParameter<int32_t>("maxSeeds")),
0042       theMagfieldToken(iCC.esConsumes()),
0043       theTrackerToken(iCC.esConsumes()),
0044       theSearchTrackerToken(iCC.esConsumes()),
0045       theTTopoToken(iCC.esConsumes()) {
0046   float ptmin = conf.getParameter<double>("ptMin");
0047   float originradius = conf.getParameter<double>("originRadius");
0048   float halflength = conf.getParameter<double>("originHalfLength");
0049   float originz = conf.getParameter<double>("originZPosition");
0050   seedpt = conf.getParameter<double>("SeedPt");
0051 
0052   geometry = conf.getUntrackedParameter<std::string>("GeometricStructure", "STANDARD");
0053   region = GlobalTrackingRegion(ptmin, originradius, halflength, originz);
0054   hitsforseeds = conf.getUntrackedParameter<std::string>("HitsForSeeds", "pairs");
0055   edm::LogInfo("SeedGeneratorForCosmics")
0056       << " PtMin of track is " << ptmin << " The Radius of the cylinder for seeds is " << originradius << "cm"
0057       << " The set Seed Momentum" << seedpt;
0058 
0059   //***top-bottom
0060   positiveYOnly = conf.getParameter<bool>("PositiveYOnly");
0061   negativeYOnly = conf.getParameter<bool>("NegativeYOnly");
0062   //***
0063 }
0064 
0065 void SeedGeneratorForCosmics::run(const SiStripRecHit2DCollection& collstereo,
0066                                   const SiStripRecHit2DCollection& collrphi,
0067                                   const SiStripMatchedRecHit2DCollection& collmatched,
0068                                   const edm::EventSetup& c,
0069                                   TrajectorySeedCollection& output) {
0070   init(collstereo, collrphi, collmatched, c);
0071   seeds(output, region);
0072   delete thePairGenerator;
0073   delete theTripletGenerator;
0074   delete thePropagatorAl;
0075   delete thePropagatorOp;
0076   delete theUpdator;
0077 }
0078 bool SeedGeneratorForCosmics::seeds(TrajectorySeedCollection& output, const TrackingRegion& region) {
0079   LogDebug("CosmicSeedFinder") << "Number of triplets " << HitTriplets.size();
0080   LogDebug("CosmicSeedFinder") << "Number of pairs " << HitPairs.size();
0081 
0082   for (unsigned int it = 0; it < HitTriplets.size(); it++) {
0083     //const TrackingRecHit *hit = &(
0084     //    const TrackingRecHit* hit = it->hits();
0085 
0086     //   GlobalPoint inner = tracker->idToDet(HitTriplets[it].inner().RecHit()->
0087     //                   geographicalId())->surface().
0088     //       toGlobal(HitTriplets[it].inner().RecHit()->localPosition());
0089     //const TrackingRecHit  *innerhit =  &(*HitTriplets[it].inner());
0090     //const TrackingRecHit  *middlehit =  &(*HitTriplets[it].middle());
0091 
0092     GlobalPoint inner = tracker->idToDet((*(HitTriplets[it].inner())).geographicalId())
0093                             ->surface()
0094                             .toGlobal((*(HitTriplets[it].inner())).localPosition());
0095 
0096     GlobalPoint middle = tracker->idToDet((*(HitTriplets[it].middle())).geographicalId())
0097                              ->surface()
0098                              .toGlobal((*(HitTriplets[it].middle())).localPosition());
0099 
0100     GlobalPoint outer = tracker->idToDet((*(HitTriplets[it].outer())).geographicalId())
0101                             ->surface()
0102                             .toGlobal((*(HitTriplets[it].outer())).localPosition());
0103 
0104     SeedingHitSet::ConstRecHitPointer outrhit = HitTriplets[it].outer();
0105     //***top-bottom
0106     SeedingHitSet::ConstRecHitPointer innrhit = HitTriplets[it].inner();
0107     if (positiveYOnly && (outrhit->globalPosition().y() < 0 || innrhit->globalPosition().y() < 0 ||
0108                           outrhit->globalPosition().y() < innrhit->globalPosition().y()))
0109       continue;
0110     if (negativeYOnly && (outrhit->globalPosition().y() > 0 || innrhit->globalPosition().y() > 0 ||
0111                           outrhit->globalPosition().y() > innrhit->globalPosition().y()))
0112       continue;
0113     //***
0114 
0115     edm::OwnVector<TrackingRecHit> hits;
0116     hits.push_back(HitTriplets[it].outer()->hit()->clone());
0117     FastHelix helix(inner, middle, outer, magfield->nominalValue(), &(*magfield));
0118     GlobalVector gv = helix.stateAtVertex().momentum();
0119     float ch = helix.stateAtVertex().charge();
0120     float mom2 = gv.x() * gv.x() + gv.y() * gv.y() + gv.z() * gv.z();
0121     if (mom2 > 1e06f || edm::isNotFinite(mom2))
0122       continue;  // ChangedByDaniele
0123 
0124     if (gv.y() > 0) {
0125       gv = -1. * gv;
0126       ch = -1. * ch;
0127     }
0128 
0129     GlobalTrajectoryParameters Gtp(outer, gv, int(ch), &(*magfield));
0130     FreeTrajectoryState CosmicSeed(Gtp, CurvilinearTrajectoryError(AlgebraicSymMatrix55(AlgebraicMatrixID())));
0131     if ((outer.y() - inner.y()) > 0) {
0132       const TSOS outerState = thePropagatorAl->propagate(
0133           CosmicSeed, tracker->idToDet((*(HitTriplets[it].outer())).geographicalId())->surface());
0134       if (outerState.isValid()) {
0135         LogDebug("CosmicSeedFinder") << "outerState " << outerState;
0136         const TSOS outerUpdated = theUpdator->update(outerState, *outrhit);
0137         if (outerUpdated.isValid()) {
0138           LogDebug("CosmicSeedFinder") << "outerUpdated " << outerUpdated;
0139 
0140           output.push_back(TrajectorySeed(trajectoryStateTransform::persistentState(
0141                                               outerUpdated, (*(HitTriplets[it].outer())).geographicalId().rawId()),
0142                                           hits,
0143                                           alongMomentum));
0144 
0145           if ((maxSeeds_ > 0) && (output.size() > size_t(maxSeeds_))) {
0146             edm::LogError("TooManySeeds") << "Found too many seeds, bailing out.\n";
0147             output.clear();
0148             return false;
0149           }
0150         }
0151       }
0152     } else {
0153       const TSOS outerState = thePropagatorOp->propagate(
0154           CosmicSeed, tracker->idToDet((*(HitTriplets[it].outer())).geographicalId())->surface());
0155       if (outerState.isValid()) {
0156         LogDebug("CosmicSeedFinder") << "outerState " << outerState;
0157         const TSOS outerUpdated = theUpdator->update(outerState, *outrhit);
0158         if (outerUpdated.isValid()) {
0159           LogDebug("CosmicSeedFinder") << "outerUpdated " << outerUpdated;
0160 
0161           output.push_back(TrajectorySeed(trajectoryStateTransform::persistentState(
0162                                               outerUpdated, (*(HitTriplets[it].outer())).geographicalId().rawId()),
0163                                           hits,
0164                                           oppositeToMomentum));
0165 
0166           if ((maxSeeds_ > 0) && (output.size() > size_t(maxSeeds_))) {
0167             edm::LogError("TooManySeeds") << "Found too many seeds, bailing out.\n";
0168             output.clear();
0169             return false;
0170           }
0171         }
0172       }
0173     }
0174   }
0175 
0176   for (unsigned int is = 0; is < HitPairs.size(); is++) {
0177     GlobalPoint inner = tracker->idToDet((*(HitPairs[is].inner())).geographicalId())
0178                             ->surface()
0179                             .toGlobal((*(HitPairs[is].inner())).localPosition());
0180     GlobalPoint outer = tracker->idToDet((*(HitPairs[is].outer())).geographicalId())
0181                             ->surface()
0182                             .toGlobal((*(HitPairs[is].outer())).localPosition());
0183 
0184     LogDebug("CosmicSeedFinder") << "inner point of the seed " << inner << " outer point of the seed " << outer;
0185     SeedingHitSet::ConstRecHitPointer outrhit = HitPairs[is].outer();
0186     //***top-bottom
0187     SeedingHitSet::ConstRecHitPointer innrhit = HitPairs[is].inner();
0188     if (positiveYOnly && (outrhit->globalPosition().y() < 0 || innrhit->globalPosition().y() < 0 ||
0189                           outrhit->globalPosition().y() < innrhit->globalPosition().y()))
0190       continue;
0191     if (negativeYOnly && (outrhit->globalPosition().y() > 0 || innrhit->globalPosition().y() > 0 ||
0192                           outrhit->globalPosition().y() > innrhit->globalPosition().y()))
0193       continue;
0194     //***
0195 
0196     edm::OwnVector<TrackingRecHit> hits;
0197     hits.push_back(HitPairs[is].outer()->hit()->clone());
0198     //    hits.push_back(HitPairs[is].inner()->clone());
0199 
0200     for (int i = 0; i < 2; i++) {
0201       //FIRST STATE IS CALCULATED CONSIDERING THAT THE CHARGE CAN BE POSITIVE OR NEGATIVE
0202       int predsign = (2 * i) - 1;
0203       if ((outer.y() - inner.y()) > 0) {
0204         GlobalTrajectoryParameters Gtp(
0205             outer, (inner - outer) * (seedpt / (inner - outer).mag()), predsign, &(*magfield));
0206 
0207         FreeTrajectoryState CosmicSeed(Gtp, CurvilinearTrajectoryError(AlgebraicSymMatrix55(AlgebraicMatrixID())));
0208 
0209         LogDebug("CosmicSeedFinder") << " FirstTSOS " << CosmicSeed;
0210         //First propagation
0211         const TSOS outerState = thePropagatorAl->propagate(
0212             CosmicSeed, tracker->idToDet((*(HitPairs[is].outer())).geographicalId())->surface());
0213         if (outerState.isValid()) {
0214           LogDebug("CosmicSeedFinder") << "outerState " << outerState;
0215           const TSOS outerUpdated = theUpdator->update(outerState, *outrhit);
0216           if (outerUpdated.isValid()) {
0217             LogDebug("CosmicSeedFinder") << "outerUpdated " << outerUpdated;
0218 
0219             PTrajectoryStateOnDet const& PTraj = trajectoryStateTransform::persistentState(
0220                 outerUpdated, (*(HitPairs[is].outer())).geographicalId().rawId());
0221 
0222             output.push_back(TrajectorySeed(PTraj, hits, alongMomentum));
0223 
0224             if ((maxSeeds_ > 0) && (output.size() > size_t(maxSeeds_))) {
0225               edm::LogError("TooManySeeds") << "Found too many seeds, bailing out.\n";
0226               output.clear();
0227               return false;
0228             }
0229 
0230           } else
0231             edm::LogWarning("CosmicSeedFinder") << " SeedForCosmics first update failed ";
0232         } else
0233           edm::LogWarning("CosmicSeedFinder") << " SeedForCosmics first propagation failed ";
0234 
0235       } else {
0236         GlobalTrajectoryParameters Gtp(
0237             outer, (outer - inner) * (seedpt / (outer - inner).mag()), predsign, &(*magfield));
0238         FreeTrajectoryState CosmicSeed(Gtp, CurvilinearTrajectoryError(AlgebraicSymMatrix55(AlgebraicMatrixID())));
0239         LogDebug("CosmicSeedFinder") << " FirstTSOS " << CosmicSeed;
0240         //First propagation
0241         const TSOS outerState = thePropagatorAl->propagate(
0242             CosmicSeed, tracker->idToDet((*(HitPairs[is].outer())).geographicalId())->surface());
0243         if (outerState.isValid()) {
0244           LogDebug("CosmicSeedFinder") << "outerState " << outerState;
0245           const TSOS outerUpdated = theUpdator->update(outerState, *outrhit);
0246           if (outerUpdated.isValid()) {
0247             LogDebug("CosmicSeedFinder") << "outerUpdated " << outerUpdated;
0248 
0249             PTrajectoryStateOnDet const& PTraj = trajectoryStateTransform::persistentState(
0250                 outerUpdated, (*(HitPairs[is].outer())).geographicalId().rawId());
0251 
0252             output.push_back(TrajectorySeed(PTraj, hits, oppositeToMomentum));
0253 
0254             if ((maxSeeds_ > 0) && (output.size() > size_t(maxSeeds_))) {
0255               edm::LogError("TooManySeeds") << "Found too many seeds, bailing out.\n";
0256               output.clear();
0257               return false;
0258             }
0259 
0260           } else
0261             edm::LogWarning("CosmicSeedFinder") << " SeedForCosmics first update failed ";
0262         } else
0263           edm::LogWarning("CosmicSeedFinder") << " SeedForCosmics first propagation failed ";
0264       }
0265     }
0266   }
0267   return true;
0268 }