SiTrackerMRHTools/src/SimpleDAFHitCollector.cc

0001 #include "RecoTracker/SiTrackerMRHTools/interface/SimpleDAFHitCollector.h"
0002 #include "RecoTracker/SiTrackerMRHTools/interface/SiTrackerMultiRecHitUpdator.h"
0003 #include "TrackingTools/DetLayers/interface/MeasurementEstimator.h"
0004 #include "TrackingTools/PatternTools/interface/Trajectory.h"
0005 #include "DataFormats/TrackingRecHit/interface/InvalidTrackingRecHit.h"
0006 #include "TrackingTools/MeasurementDet/interface/MeasurementDet.h"
0007 #include "RecoTracker/MeasurementDet/interface/MeasurementTrackerEvent.h"
0008 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit1D.h"
0009 #include "TrackingTools/TrackFitters/interface/TrajectoryStateCombiner.h"
0010
0011 #ifdef EDM_ML_DEBUG
0012 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
0013 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0014 #include "Geometry/Records/interface/IdealGeometryRecord.h"
0015 #include "DataFormats/MuonDetId/interface/CSCDetId.h"
0016 #include "DataFormats/MuonDetId/interface/DTWireId.h"
0017 #include "DataFormats/MuonDetId/interface/RPCDetId.h"
0018 #include "DataFormats/MuonDetId/interface/MuonSubdetId.h"
0019 #endif
0020
0021 #include <vector>
0022 #include <map>
0023
0024 #define _debug_SimpleDAFHitCollector_
0025
0026 using namespace std;
0027
0028 vector<TrajectoryMeasurement> SimpleDAFHitCollector::recHits(const Trajectory& traj,
0029                                                              const MeasurementTrackerEvent* theMTE) const {
0030   LogTrace("MultiRecHitCollector") << " Calling SimpleDAFHitCollector::recHits" << std::endl;
0031
0032   //WARNING: At the moment the trajectories has the measurements
0033   //with reversed sorting after the track smoothing
0034   const vector<TrajectoryMeasurement>& meas = traj.measurements();
0035   unsigned int hitcounter = 1;
0036
0037   if (meas.empty())
0038     return vector<TrajectoryMeasurement>();
0039
0040   LogTrace("MultiRecHitCollector") << "  Original measurements are:";
0041   Debug(meas);
0042
0043   //groups hits on a sensor by sensor with same Id of previous TM
0044   //we have to sort the TrajectoryMeasurements in the opposite way in the fitting direction
0045   vector<TrajectoryMeasurement> result;
0046   for (vector<TrajectoryMeasurement>::const_reverse_iterator itrajmeas = meas.rbegin(); itrajmeas < meas.rend();
0047        itrajmeas++, hitcounter++) {
0048     DetId id = itrajmeas->recHit()->geographicalId();
0049     MeasurementDetWithData measDet = theMTE->idToDet(id);
0050     tracking::TempMeasurements tmps;
0051
0052     std::vector<const TrackingRecHit*> hits;
0053     std::vector<std::unique_ptr<const TrackingRecHit>> hitsOwner;
0054
0055     TrajectoryStateOnSurface smoothtsos = itrajmeas->updatedState();
0056     //the error is scaled in order to take more "compatible" hits
0057     //      if( smoothtsos.isValid() ) smoothtsos.rescaleError(10);
0058
0059     TrajectoryStateOnSurface predtsos_fwd = itrajmeas->predictedState();
0060     TrajectoryStateOnSurface predtsos_bwd = itrajmeas->backwardPredictedState();
0061     if (!predtsos_fwd.isValid() || !predtsos_bwd.isValid()) {
0062       LogTrace("MultiRecHitCollector") << "Something wrong! no valid TSOS found in current group ";
0063       continue;
0064     }
0065
0066     TrajectoryStateCombiner combiner;
0067     TrajectoryStateOnSurface combtsos;
0068     if (hitcounter == meas.size())
0069       combtsos = predtsos_fwd;
0070     else if (hitcounter == 1)
0071       combtsos = predtsos_bwd;
0072     else
0073       combtsos = combiner(predtsos_bwd, predtsos_fwd);
0074
0075     //collected hits compatible with the itrajmeas
0076     if (measDet.measurements(smoothtsos, *(getEstimator()), tmps)) {
0077       LogTrace("MultiRecHitCollector") << "  Found " << tmps.size() << " compatible measurements";
0078
0079       for (std::size_t i = 0; i != tmps.size(); ++i) {
0080         DetId idtemps = tmps.hits[i]->geographicalId();
0081
0082         if (idtemps == id && tmps.hits[i]->hit()->isValid()) {
0083           LogTrace("MultiRecHitCollector") << "  This is valid with position " << tmps.hits[i]->hit()->localPosition()
0084                                            << " and error " << tmps.hits[i]->hit()->localPositionError();
0085
0086           TransientTrackingRecHit::RecHitPointer transient = theUpdator->getBuilder()->build(tmps.hits[i]->hit());
0087           TrackingRecHit::ConstRecHitPointer preciseHit = theHitCloner.makeShared(transient, combtsos);
0088           auto righthit = rightdimension(*preciseHit);
0089           hitsOwner.push_back(std::move(righthit));
0090           hits.push_back(hitsOwner.back().get());
0091         }
0092       }
0093
0094       //the error was scaled, now is scaled back (even if this is not tightly necessary)
0095       //        if (smoothtsos.isValid()) smoothtsos.rescaleError(0.1);
0096
0097       //I will keep the Invalid hit, IF this is not the first one
0098       if (hits.empty()) {
0099         LogTrace("MultiRecHitCollector") << " -> but no valid hits found in current group.";
0100
0101         if (result.empty())
0102           continue;
0103
0104         result.push_back(TrajectoryMeasurement(
0105             predtsos_fwd, std::make_shared<InvalidTrackingRecHit>(measDet.mdet().geomDet(), TrackingRecHit::missing)));
0106       } else {
0107         //measurements in groups are sorted with increating chi2
0108         //sort( *hits.begin(), *hits.end(), TrajMeasLessEstim());
0109         if (!itrajmeas->recHit()->isValid())
0110           LogTrace("MultiRecHitCollector")
0111               << "  -> " << hits.size() << " valid hits for this sensor. (IT WAS INVALID!!!)";
0112         else
0113           LogTrace("MultiRecHitCollector") << "  -> " << hits.size() << " valid hits for this sensor.";
0114
0115         //building a MultiRecHit out of each sensor group
0116         result.push_back(TrajectoryMeasurement(predtsos_fwd, theUpdator->buildMultiRecHit(hits, combtsos, measDet)));
0117       }
0118     } else {
0119       LogTrace("MultiRecHitCollector") << "  No measurements found in current group.";
0120       //the error was scaled, now is scaled back (even if this is not tightly necessary)
0121       //          if (smoothtsos.isValid()) smoothtsos.rescaleError(0.1);
0122
0123       if (result.empty())
0124         continue;
0125
0126       result.push_back(TrajectoryMeasurement(
0127           predtsos_fwd, std::make_shared<InvalidTrackingRecHit>(measDet.mdet().geomDet(), TrackingRecHit::missing)));
0128     }
0129   }
0130   LogTrace("MultiRecHitCollector") << " Ending SimpleDAFHitCollector::recHits >> " << result.size();
0131
0132   //LogTrace("MultiRecHitCollector") << "  New measurements are:";
0133   //Debug(result);
0134
0135   //adding a protection against too few hits and invalid hits
0136   //(due to failed propagation on the same surface of the original hits)
0137   if (result.size() > 2) {
0138     int hitcounter = 0;
0139     //check if the vector result has more than 3 valid hits
0140     for (vector<TrajectoryMeasurement>::const_iterator iimeas = result.begin(); iimeas != result.end(); ++iimeas) {
0141       if (iimeas->recHit()->isValid())
0142         hitcounter++;
0143     }
0144
0145     if (hitcounter > 2)
0146       return result;
0147     else
0148       return vector<TrajectoryMeasurement>();
0149   }
0150
0151   else {
0152     return vector<TrajectoryMeasurement>();
0153   }
0154 }
0155
0156 void SimpleDAFHitCollector::Debug(const std::vector<TrajectoryMeasurement> TM) const {
0157 #ifdef EDM_ML_DEBUG
0158   for (vector<TrajectoryMeasurement>::const_iterator itrajmeas = TM.begin(); itrajmeas < TM.end(); itrajmeas++) {
0159     if (itrajmeas->recHit()->isValid()) {
0160       LogTrace("MultiRecHitCollector") << "  Valid Hit with DetId " << itrajmeas->recHit()->geographicalId().rawId()
0161                                        << " and dim:"
0162                                        << itrajmeas->recHit()->dimension()
0163                                        //<< " type " << typeid(itrajmeas->recHit()).name()
0164                                        << " local position " << itrajmeas->recHit()->hit()->localPosition()
0165                                        << " global position " << itrajmeas->recHit()->hit()->globalPosition()
0166                                        << " and r " << itrajmeas->recHit()->hit()->globalPosition().perp();
0167
0168       DetId hitId = itrajmeas->recHit()->geographicalId();
0169
0170       if (hitId.det() == DetId::Tracker) {
0171         switch (hitId.subdetId()) {
0172           case StripSubdetector::TIB:
0173             LogTrace("MultiRecHitCollector") << "  I am TIB " << theTopology->tibLayer(hitId);
0174             break;
0175           case StripSubdetector::TOB:
0176             LogTrace("MultiRecHitCollector") << "  I am TOB " << theTopology->tobLayer(hitId);
0177             break;
0178           case StripSubdetector::TEC:
0179             LogTrace("MultiRecHitCollector") << "  I am TEC " << theTopology->tecWheel(hitId);
0180             break;
0181           case StripSubdetector::TID:
0182             LogTrace("MultiRecHitCollector") << "  I am TID " << theTopology->tidWheel(hitId);
0183             break;
0184           case PixelSubdetector::PixelBarrel:
0185             LogTrace("MultiRecHitCollector") << "  I am PixBar " << theTopology->pxbLayer(hitId);
0186             break;
0187           case PixelSubdetector::PixelEndcap:
0188             LogTrace("MultiRecHitCollector") << "  I am PixFwd " << theTopology->pxfDisk(hitId);
0189             break;
0190           default:
0191             LogTrace("MultiRecHitCollector") << "  UNKNOWN TRACKER HIT TYPE ";
0192             break;
0193         }
0194       } else if (hitId.det() == DetId::Muon) {
0195         if (hitId.subdetId() == MuonSubdetId::DT)
0196           LogTrace("MultiRecHitCollector") << "  I am DT " << DTWireId(hitId);
0197         else if (hitId.subdetId() == MuonSubdetId::CSC)
0198           LogTrace("MultiRecHitCollector") << "  I am CSC " << CSCDetId(hitId);
0199         else if (hitId.subdetId() == MuonSubdetId::RPC)
0200           LogTrace("MultiRecHitCollector") << "  I am RPC " << RPCDetId(hitId);
0201         else
0202           LogTrace("MultiRecHitCollector") << "  UNKNOWN MUON HIT TYPE ";
0203       } else
0204         LogTrace("MultiRecHitCollector") << "  UNKNOWN HIT TYPE ";
0205
0206       LogTrace("MultiRecHitCollector") << "  TSOS predicted_fwd " << itrajmeas->predictedState().localPosition();
0207       LogTrace("MultiRecHitCollector") << "  TSOS predicted_bwd "
0208                                        << itrajmeas->backwardPredictedState().localPosition();
0209       LogTrace("MultiRecHitCollector") << "  TSOS smoothtsos " << itrajmeas->updatedState().localPosition();
0210     } else {
0211       LogTrace("MultiRecHitCollector") << "  Invalid Hit with DetId " << itrajmeas->recHit()->geographicalId().rawId();
0212     }
0213     LogTrace("MultiRecHitCollector") << "\n";
0214   }
0215 #endif
0216 }