Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2025-01-18 03:42:23

0001 #include "DataFormats/TrackCandidate/interface/TrackCandidate.h"
0002 #include "DataFormats/TrackReco/interface/Track.h"
0003 #include "DataFormats/TrackerRecHit2D/interface/OmniClusterRef.h"
0004 #include "DataFormats/TrackerRecHit2D/interface/SiTrackerMultiRecHit.h"
0005 #include "DataFormats/TrackerRecHit2D/interface/TkCloner.h"
0006 #include "DataFormats/TrackingRecHit/interface/InvalidTrackingRecHit.h"
0007 #include "DataFormats/TrackingRecHit/interface/TrackingRecHitFwd.h"
0008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0009 #include "Geometry/CommonDetUnit/interface/TrackingGeometry.h"
0010 #include "MagneticField/Engine/interface/MagneticField.h"
0011 #include "RecoTracker/SiTrackerMRHTools/interface/MultiRecHitCollector.h"
0012 #include "RecoTracker/SiTrackerMRHTools/interface/SiTrackerMultiRecHitUpdator.h"
0013 #include "RecoTracker/TrackProducer/interface/DAFTrackProducerAlgorithm.h"
0014 #include "RecoTracker/TransientTrackingRecHit/interface/TkClonerImpl.h"
0015 #include "RecoTracker/TransientTrackingRecHit/interface/TkTransientTrackingRecHitBuilder.h"
0016 #include "TrackingTools/PatternTools/interface/TSCBLBuilderNoMaterial.h"
0017 #include "TrackingTools/PatternTools/interface/TrajAnnealing.h"
0018 #include "TrackingTools/PatternTools/interface/TransverseImpactPointExtrapolator.h"
0019 #include "TrackingTools/TrackFitters/interface/TrajectoryFitter.h"
0020 #include "TrackingTools/TrackFitters/interface/TrajectoryStateCombiner.h"
0021 #include "TrackingTools/TrackFitters/interface/TrajectoryStateWithArbitraryError.h"
0022 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
0023 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateTransform.h"
0024 #include "TrackingTools/TransientTrackingRecHit/interface/TransientTrackingRecHit.h"
0025 #include "TrackingTools/TransientTrackingRecHit/interface/TransientTrackingRecHitBuilder.h"
0026 
0027 DAFTrackProducerAlgorithm::DAFTrackProducerAlgorithm(const edm::ParameterSet& conf)
0028     : conf_(conf), minHits_(conf.getParameter<int>("MinHits")) {}
0029 
0030 void DAFTrackProducerAlgorithm::fillPSetDescription(edm::ParameterSetDescription& desc) { desc.add<int>("MinHits", 3); }
0031 
0032 void DAFTrackProducerAlgorithm::runWithCandidate(const TrackingGeometry* theG,
0033                                                  const MagneticField* theMF,
0034                                                  const TrajTrackAssociationCollection& TTmap,
0035                                                  const MeasurementTrackerEvent* measTk,
0036                                                  const TrajectoryFitter* theFitter,
0037                                                  const TransientTrackingRecHitBuilder* builder,
0038                                                  const MultiRecHitCollector* measurementCollector,
0039                                                  const SiTrackerMultiRecHitUpdator* updator,
0040                                                  const reco::BeamSpot& bs,
0041                                                  AlgoProductCollection& algoResults,
0042                                                  TrajAnnealingCollection& trajann,
0043                                                  bool TrajAnnSaving_,
0044                                                  AlgoProductCollection& algoResultsBeforeDAF,
0045                                                  AlgoProductCollection& algoResultsAfterDAF) const {
0046   LogDebug("DAFTrackProducerAlgorithm") << "Size of map: " << TTmap.size() << "\n";
0047 
0048   int cont = 0;
0049   int nTracksChanged = 0;
0050 
0051   for (TrajTrackAssociationCollection::const_iterator itTTmap = TTmap.begin(); itTTmap != TTmap.end(); itTTmap++) {
0052     const edm::Ref<std::vector<Trajectory> > BeforeDAFTraj = itTTmap->key;
0053     std::vector<TrajectoryMeasurement> BeforeDAFTrajMeas = BeforeDAFTraj->measurements();
0054     const reco::TrackRef BeforeDAFTrack = itTTmap->val;
0055 
0056     float ndof = 0;
0057     Trajectory CurrentTraj;
0058 
0059     if (BeforeDAFTraj->isValid()) {
0060       LogDebug("DAFTrackProducerAlgorithm") << "The trajectory #" << cont + 1 << " is valid. \n";
0061 
0062       //getting the MultiRecHit collection and the trajectory with a first fit-smooth round
0063       std::pair<TransientTrackingRecHit::RecHitContainer, TrajectoryStateOnSurface> hits =
0064           collectHits(*BeforeDAFTraj, measurementCollector, &*measTk);
0065 
0066       //new initial fit
0067       CurrentTraj = fit(hits, theFitter, *BeforeDAFTraj);
0068 
0069       //starting the annealing program
0070       for (std::vector<double>::const_iterator ian = updator->getAnnealingProgram().begin();
0071            ian != updator->getAnnealingProgram().end();
0072            ian++) {
0073         if (CurrentTraj.isValid()) {
0074           LogDebug("DAFTrackProducerAlgorithm") << "Seed direction is " << CurrentTraj.seed().direction()
0075                                                 << ".Traj direction is " << CurrentTraj.direction() << std::endl;
0076 
0077           //updating MultiRecHits and fit-smooth again
0078           std::pair<TransientTrackingRecHit::RecHitContainer, TrajectoryStateOnSurface> curiterationhits =
0079               updateHits(CurrentTraj, updator, &*measTk, *ian);
0080           if (curiterationhits.first.size() < 3) {
0081             LogDebug("DAFTrackProducerAlgorithm")
0082                 << "Rejecting trajectory with " << curiterationhits.first.size() << " hits" << std::endl;
0083             CurrentTraj = Trajectory();
0084             break;
0085           }
0086 
0087           CurrentTraj = fit(curiterationhits, theFitter, CurrentTraj);
0088 
0089           //saving trajectory for each annealing cycle ...
0090           if (TrajAnnSaving_) {
0091             TrajAnnealing temp(CurrentTraj, *ian);
0092             trajann.push_back(temp);
0093           }
0094 
0095           LogDebug("DAFTrackProducerAlgorithm") << "done annealing value " << (*ian);
0096 
0097         } else
0098           break;
0099       }  //end annealing cycle
0100 
0101       int percOfHitsUnchangedAfterDAF =
0102           (1. * checkHits(*BeforeDAFTraj, CurrentTraj) / (1. * BeforeDAFTrajMeas.size())) * 100.;
0103       LogDebug("DAFTrackProducerAlgorithm")
0104           << "Ended annealing program with " << percOfHitsUnchangedAfterDAF << " unchanged." << std::endl;
0105 
0106       //computing the ndof keeping into account the weights
0107       ndof = calculateNdof(CurrentTraj);
0108 
0109       //checking if the trajectory has the minimum number of valid hits ( weight (>1e-6) )
0110       //in order to remove tracks with too many outliers.
0111       int goodHits = countingGoodHits(CurrentTraj);
0112 
0113       if (goodHits >= minHits_) {
0114         bool ok = buildTrack(CurrentTraj, algoResults, ndof, bs, &BeforeDAFTrack);
0115         // or filtered?
0116         if (ok)
0117           cont++;
0118 
0119         //saving tracks before and after DAF
0120         if ((100. - percOfHitsUnchangedAfterDAF) > 0.) {
0121           bool okBefore = buildTrack(*BeforeDAFTraj, algoResultsBeforeDAF, ndof, bs, &BeforeDAFTrack);
0122           bool okAfter = buildTrack(CurrentTraj, algoResultsAfterDAF, ndof, bs, &BeforeDAFTrack);
0123           if (okBefore && okAfter)
0124             nTracksChanged++;
0125         }
0126       } else {
0127         LogDebug("DAFTrackProducerAlgorithm") << "Rejecting trajectory with " << CurrentTraj.foundHits() << " hits";
0128       }
0129     }  //end run on track collection
0130     else {
0131       LogDebug("DAFTrackProducerAlgorithm") << "Rejecting empty trajectory" << std::endl;
0132     }
0133   }  //end run on track collection
0134 
0135   LogDebug("DAFTrackProducerAlgorithm") << "Number of Tracks found:   " << cont << "\n";
0136   LogDebug("DAFTrackProducerAlgorithm") << "Number of Tracks changed: " << nTracksChanged << "\n";
0137 }
0138 /*------------------------------------------------------------------------------------------------------*/
0139 std::pair<TransientTrackingRecHit::RecHitContainer, TrajectoryStateOnSurface> DAFTrackProducerAlgorithm::collectHits(
0140     const Trajectory vtraj,
0141     const MultiRecHitCollector* measurementCollector,
0142     const MeasurementTrackerEvent* measTk) const {
0143   LogDebug("DAFTrackProducerAlgorithm") << "Calling DAFTrackProducerAlgorithm::collectHits";
0144 
0145   //getting the traj measurements from the MeasurementCollector
0146   TransientTrackingRecHit::RecHitContainer hits;
0147   std::vector<TrajectoryMeasurement> collectedmeas = measurementCollector->recHits(vtraj, measTk);
0148 
0149   //if the MeasurementCollector is empty, make an "empty" pair
0150   //else taking the collected measured hits and building the pair
0151   if (collectedmeas.empty())
0152     return std::make_pair(TransientTrackingRecHit::RecHitContainer(), TrajectoryStateOnSurface());
0153 
0154   for (std::vector<TrajectoryMeasurement>::const_iterator iter = collectedmeas.begin(); iter != collectedmeas.end();
0155        iter++) {
0156     hits.push_back(iter->recHit());
0157   }
0158 
0159   //TrajectoryStateWithArbitraryError() == Creates a TrajectoryState with the same parameters
0160   //     as the input one, but with "infinite" errors, i.e. errors so big that they don't
0161   //     bias a fit starting with this state.
0162   //return std::make_pair(hits,TrajectoryStateWithArbitraryError()(collectedmeas.front().predictedState()));
0163 
0164   // I do not have to rescale the error because it is already rescaled in the fit code
0165   TrajectoryStateOnSurface initialStateFromTrack = collectedmeas.front().predictedState();
0166 
0167   LogDebug("DAFTrackProducerAlgorithm")
0168       << "Pair (hits, TSOS)  with TSOS predicted(collectedmeas.front().predictedState())";
0169   return std::make_pair(hits, initialStateFromTrack);
0170 }
0171 /*------------------------------------------------------------------------------------------------------*/
0172 std::pair<TransientTrackingRecHit::RecHitContainer, TrajectoryStateOnSurface> DAFTrackProducerAlgorithm::updateHits(
0173     const Trajectory vtraj,
0174     const SiTrackerMultiRecHitUpdator* updator,
0175     const MeasurementTrackerEvent* theMTE,
0176     double annealing) const {
0177   LogDebug("DAFTrackProducerAlgorithm") << "Calling DAFTrackProducerAlgorithm::updateHits";
0178   TransientTrackingRecHit::RecHitContainer hits;
0179   std::vector<TrajectoryMeasurement> vmeas = vtraj.measurements();
0180   std::vector<TrajectoryMeasurement>::reverse_iterator imeas;
0181   unsigned int hitcounter = 1;
0182 
0183   //I run inversely on the trajectory obtained and update the state
0184   for (imeas = vmeas.rbegin(); imeas != vmeas.rend(); imeas++, hitcounter++) {
0185     DetId id = imeas->recHit()->geographicalId();
0186     MeasurementDetWithData measDet = theMTE->idToDet(id);
0187 
0188     TrajectoryStateCombiner combiner;
0189     TrajectoryStateOnSurface combtsos;
0190     if (hitcounter == vmeas.size())
0191       combtsos = imeas->predictedState();  //fwd
0192     else if (hitcounter == 1)
0193       combtsos = imeas->backwardPredictedState();  //bwd
0194     else
0195       combtsos = combiner(imeas->backwardPredictedState(), imeas->predictedState());
0196 
0197     PrintHit(&*imeas->recHit(), combtsos);
0198     if (imeas->recHit()->isValid()) {
0199       TransientTrackingRecHit::RecHitPointer updated = updator->update(imeas->recHit(), combtsos, measDet, annealing);
0200       hits.push_back(updated);
0201     } else {
0202       hits.push_back(imeas->recHit());
0203     }
0204   }
0205 
0206   TrajectoryStateOnSurface updatedStateFromTrack = vmeas.back().predictedState();
0207 
0208   //return std::make_pair(hits,TrajectoryStateWithArbitraryError()(vmeas.back().updatedState()));
0209   LogDebug("DAFTrackProducerAlgorithm") << "Pair (hits, TSOS)  with TSOS predicted (vmeas.back().predictedState())";
0210 
0211   return std::make_pair(hits, updatedStateFromTrack);
0212 }
0213 /*------------------------------------------------------------------------------------------------------*/
0214 Trajectory DAFTrackProducerAlgorithm::fit(
0215     const std::pair<TransientTrackingRecHit::RecHitContainer, TrajectoryStateOnSurface>& hits,
0216     const TrajectoryFitter* theFitter,
0217     Trajectory vtraj) const {
0218   //creating a new trajectory starting from the direction of the seed of the input one and the hits
0219   Trajectory newVec = theFitter->fitOne(TrajectorySeed({}, {}, vtraj.seed().direction()), hits.first, hits.second);
0220 
0221   if (newVec.isValid())
0222     return newVec;
0223   else {
0224     LogDebug("DAFTrackProducerAlgorithm") << "Fit no valid.";
0225     return Trajectory();
0226   }
0227 }
0228 /*------------------------------------------------------------------------------------------------------*/
0229 bool DAFTrackProducerAlgorithm::buildTrack(const Trajectory vtraj,
0230                                            AlgoProductCollection& algoResults,
0231                                            float ndof,
0232                                            const reco::BeamSpot& bs,
0233                                            const reco::TrackRef* BeforeDAFTrack) const {
0234   LogDebug("DAFTrackProducerAlgorithm") << " BUILDER " << std::endl;
0235   ;
0236   TrajectoryStateOnSurface innertsos;
0237 
0238   if (vtraj.isValid()) {
0239     std::unique_ptr<Trajectory> theTraj(new Trajectory(vtraj));
0240 
0241     if (vtraj.direction() == alongMomentum) {
0242       //if (theTraj->direction() == oppositeToMomentum) {
0243       innertsos = vtraj.firstMeasurement().updatedState();
0244     } else {
0245       innertsos = vtraj.lastMeasurement().updatedState();
0246     }
0247 
0248     TSCBLBuilderNoMaterial tscblBuilder;
0249     TrajectoryStateClosestToBeamLine tscbl = tscblBuilder(*(innertsos.freeState()), bs);
0250 
0251     if (tscbl.isValid() == false)
0252       return false;
0253 
0254     GlobalPoint v = tscbl.trackStateAtPCA().position();
0255     math::XYZPoint pos(v.x(), v.y(), v.z());
0256     GlobalVector p = tscbl.trackStateAtPCA().momentum();
0257     math::XYZVector mom(p.x(), p.y(), p.z());
0258 
0259     //    LogDebug("TrackProducer") <<v<<p<<std::endl;
0260 
0261     auto algo = (*BeforeDAFTrack)->algo();
0262     std::unique_ptr<reco::Track> theTrack(new reco::Track(vtraj.chiSquared(),
0263                                                           ndof,  //in the DAF the ndof is not-integer
0264                                                           pos,
0265                                                           mom,
0266                                                           tscbl.trackStateAtPCA().charge(),
0267                                                           tscbl.trackStateAtPCA().curvilinearError(),
0268                                                           algo));
0269     theTrack->setQualityMask((*BeforeDAFTrack)->qualityMask());
0270 
0271     AlgoProduct aProduct{theTraj.release(), theTrack.release(), vtraj.direction(), 0};
0272     algoResults.push_back(aProduct);
0273 
0274     return true;
0275   } else {
0276     LogDebug("DAFTrackProducerAlgorithm") << " BUILDER NOT POSSIBLE: traj is not valid" << std::endl;
0277     ;
0278     return false;
0279   }
0280 }
0281 /*------------------------------------------------------------------------------------------------------*/
0282 int DAFTrackProducerAlgorithm::countingGoodHits(const Trajectory traj) const {
0283   int ngoodhits = 0;
0284   std::vector<TrajectoryMeasurement> vtm = traj.measurements();
0285 
0286   for (std::vector<TrajectoryMeasurement>::const_iterator tm = vtm.begin(); tm != vtm.end(); tm++) {
0287     //if the rechit is valid
0288     if (tm->recHit()->isValid()) {
0289       SiTrackerMultiRecHit const& mHit = dynamic_cast<SiTrackerMultiRecHit const&>(*tm->recHit());
0290       std::vector<const TrackingRecHit*> components = mHit.recHits();
0291 
0292       int iComp = 0;
0293 
0294       for (std::vector<const TrackingRecHit*>::const_iterator iter = components.begin(); iter != components.end();
0295            iter++, iComp++) {
0296         //if there is at least one component with weight higher than 1e-6 then the hit is not an outlier
0297         if (mHit.weight(iComp) > 1e-6) {
0298           ngoodhits++;
0299           iComp++;
0300           break;
0301         }
0302       }
0303     }
0304   }
0305 
0306   LogDebug("DAFTrackProducerAlgorithm") << "Original number of valid hits " << traj.foundHits()
0307                                         << " -> hit with good weight (>1e-6) are " << ngoodhits;
0308   return ngoodhits;
0309 }
0310 /*------------------------------------------------------------------------------------------------------*/
0311 
0312 void DAFTrackProducerAlgorithm::filter(const TrajectoryFitter* fitter,
0313                                        std::vector<Trajectory>& input,
0314                                        int minhits,
0315                                        std::vector<Trajectory>& output,
0316                                        const TransientTrackingRecHitBuilder* builder) const {
0317   if (input.empty())
0318     return;
0319 
0320   int ngoodhits = 0;
0321   std::vector<TrajectoryMeasurement> vtm = input[0].measurements();
0322   TransientTrackingRecHit::RecHitContainer hits;
0323 
0324   //count the number of non-outlier and non-invalid hits
0325   for (std::vector<TrajectoryMeasurement>::reverse_iterator tm = vtm.rbegin(); tm != vtm.rend(); tm++) {
0326     //if the rechit is valid
0327     if (tm->recHit()->isValid()) {
0328       SiTrackerMultiRecHit const& mHit = dynamic_cast<SiTrackerMultiRecHit const&>(*tm->recHit());
0329       std::vector<const TrackingRecHit*> components = mHit.recHits();
0330       int iComp = 0;
0331       bool isGood = false;
0332       for (std::vector<const TrackingRecHit*>::const_iterator iter = components.begin(); iter != components.end();
0333            iter++, iComp++) {
0334         //if there is at least one component with weight higher than 1e-6 then the hit is not an outlier
0335         if (mHit.weight(iComp) > 1e-6) {
0336           ngoodhits++;
0337           iComp++;
0338           isGood = true;
0339           break;
0340         }
0341       }
0342       if (isGood) {
0343         TkClonerImpl hc = static_cast<TkTransientTrackingRecHitBuilder const*>(builder)->cloner();
0344         auto tempHit = hc.makeShared(tm->recHit(), tm->updatedState());
0345         hits.push_back(tempHit);
0346       } else
0347         hits.push_back(std::make_shared<InvalidTrackingRecHit>(*tm->recHit()->det(), TrackingRecHit::missing));
0348     } else {
0349       TkClonerImpl hc = static_cast<TkTransientTrackingRecHitBuilder const*>(builder)->cloner();
0350       auto tempHit = hc.makeShared(tm->recHit(), tm->updatedState());
0351       hits.push_back(tempHit);
0352     }
0353   }
0354 
0355   LogDebug("DAFTrackProducerAlgorithm") << "Original number of valid hits " << input[0].foundHits()
0356                                         << "; after filtering " << ngoodhits;
0357   if (ngoodhits > input[0].foundHits())
0358     edm::LogError("DAFTrackProducerAlgorithm")
0359         << "Something wrong: the number of good hits from DAFTrackProducerAlgorithm::filter " << ngoodhits
0360         << " is higher than the original one " << input[0].foundHits();
0361 
0362   if (ngoodhits < minhits)
0363     return;
0364 
0365   TrajectoryStateOnSurface curstartingTSOS = input.front().lastMeasurement().updatedState();
0366   LogDebug("DAFTrackProducerAlgorithm") << "starting tsos for final refitting " << curstartingTSOS;
0367   //curstartingTSOS.rescaleError(100);
0368 
0369   output = fitter->fit(TrajectorySeed({}, {}, input.front().seed().direction()),
0370                        hits,
0371                        TrajectoryStateWithArbitraryError()(curstartingTSOS));
0372 
0373   LogDebug("DAFTrackProducerAlgorithm") << "After filtering " << output.size() << " trajectories";
0374 }
0375 /*------------------------------------------------------------------------------------------------------*/
0376 float DAFTrackProducerAlgorithm::calculateNdof(const Trajectory vtraj) const {
0377   if (!vtraj.isValid())
0378     return 0;
0379   float ndof = 0;
0380   const std::vector<TrajectoryMeasurement>& meas = vtraj.measurements();
0381   for (std::vector<TrajectoryMeasurement>::const_iterator iter = meas.begin(); iter != meas.end(); iter++) {
0382     if (iter->recHit()->isValid()) {
0383       SiTrackerMultiRecHit const& mHit = dynamic_cast<SiTrackerMultiRecHit const&>(*iter->recHit());
0384       std::vector<const TrackingRecHit*> components = mHit.recHits();
0385       int iComp = 0;
0386       for (std::vector<const TrackingRecHit*>::const_iterator iter2 = components.begin(); iter2 != components.end();
0387            iter2++, iComp++) {
0388         if ((*iter2)->isValid())
0389           ndof += ((*iter2)->dimension()) * mHit.weight(iComp);
0390       }
0391     }
0392   }
0393 
0394   return ndof - 5;
0395 }
0396 //------------------------------------------------------------------------------------------------
0397 int DAFTrackProducerAlgorithm::checkHits(Trajectory iInitTraj, const Trajectory iFinalTraj) const {
0398   std::vector<TrajectoryMeasurement> initmeasurements = iInitTraj.measurements();
0399   std::vector<TrajectoryMeasurement> finalmeasurements = iFinalTraj.measurements();
0400   std::vector<TrajectoryMeasurement>::iterator jmeas;
0401   int nSame = 0;
0402   int ihit = 0;
0403 
0404   if (initmeasurements.empty() || finalmeasurements.empty()) {
0405     LogDebug("DAFTrackProducerAlgorithm") << "Initial or Final Trajectory empty.";
0406     return 0;
0407   }
0408 
0409   if (initmeasurements.size() != finalmeasurements.size()) {
0410     LogDebug("DAFTrackProducerAlgorithm")
0411         << "Initial Trajectory size(" << initmeasurements.size() << " hits) "
0412         << "is different to final traj size (" << finalmeasurements.size() << ")! No checkHits possible! ";
0413     return 0;
0414   }
0415 
0416   for (jmeas = initmeasurements.begin(); jmeas != initmeasurements.end(); jmeas++) {
0417     const TrackingRecHit* initHit = jmeas->recHit()->hit();
0418     if (!initHit->isValid() && ihit == 0)
0419       continue;
0420 
0421     if (initHit->isValid()) {
0422       TrajectoryMeasurement imeas = finalmeasurements.at(ihit);
0423       const TrackingRecHit* finalHit = imeas.recHit()->hit();
0424       const TrackingRecHit* MaxWeightHit = nullptr;
0425       float maxweight = 0;
0426 
0427       const SiTrackerMultiRecHit* mrh = dynamic_cast<const SiTrackerMultiRecHit*>(finalHit);
0428       if (mrh) {
0429         std::vector<const TrackingRecHit*> components = mrh->recHits();
0430         std::vector<const TrackingRecHit*>::const_iterator icomp;
0431         int hitcounter = 0;
0432 
0433         for (icomp = components.begin(); icomp != components.end(); icomp++) {
0434           if ((*icomp)->isValid()) {
0435             double weight = mrh->weight(hitcounter);
0436             if (weight > maxweight) {
0437               MaxWeightHit = *icomp;
0438               maxweight = weight;
0439             }
0440           }
0441           hitcounter++;
0442         }
0443         if (!MaxWeightHit)
0444           continue;
0445 
0446         auto myref1 = reinterpret_cast<const BaseTrackerRecHit*>(initHit)->firstClusterRef();
0447         auto myref2 = reinterpret_cast<const BaseTrackerRecHit*>(MaxWeightHit)->firstClusterRef();
0448 
0449         if (myref1 == myref2) {
0450           nSame++;
0451         } else {
0452           LogDebug("DAFTrackProducerAlgorithm") << "diverso hit!" << std::endl;
0453           TrajectoryStateOnSurface dummState;
0454           LogTrace("DAFTrackProducerAlgorithm") << "  This hit was:\n ";
0455           PrintHit(initHit, dummState);
0456           LogTrace("DAFTrackProducerAlgorithm") << "  instead now is:\n ";
0457           PrintHit(MaxWeightHit, dummState);
0458         }
0459       }
0460     } else {
0461       TrajectoryMeasurement imeas = finalmeasurements.at(ihit);
0462       const TrackingRecHit* finalHit = imeas.recHit()->hit();
0463       if (!finalHit->isValid()) {
0464         nSame++;
0465       }
0466     }
0467 
0468     ihit++;
0469   }
0470 
0471   return nSame;
0472 }
0473 
0474 void DAFTrackProducerAlgorithm::PrintHit(const TrackingRecHit* const& hit, TrajectoryStateOnSurface& tsos) const {
0475   if (hit->isValid()) {
0476     LogTrace("DAFTrackProducerAlgorithm")
0477         << "  Valid Hit with DetId " << hit->geographicalId().rawId() << " and dim:" << hit->dimension()
0478         << " local position " << hit->localPosition() << " global position " << hit->globalPosition() << " and r "
0479         << hit->globalPosition().perp();
0480     if (tsos.isValid())
0481       LogTrace("DAFTrackProducerAlgorithm") << "  TSOS combtsos " << tsos.localPosition();
0482 
0483   } else {
0484     LogTrace("DAFTrackProducerAlgorithm") << "  Invalid Hit with DetId " << hit->geographicalId().rawId();
0485   }
0486 }