50 bool enableDTMeasurement = par.
getParameter<
bool>(
"EnableDTMeasurement");
51 bool enableCSCMeasurement = par.
getParameter<
bool>(
"EnableCSCMeasurement");
52 bool enableRPCMeasurement = par.
getParameter<
bool>(
"EnableRPCMeasurement");
68 enableRPCMeasurement);
94 category_ =
"Muon|RecoMuon|CosmicMuon|CosmicMuonTrajectoryBuilder";
136 vector<Trajectory*> trajL = vector<Trajectory*>();
141 DetId did(ptsd1.detId());
151 vector<const DetLayer*> navLayers;
164 LogTrace(
category_) <<
"found "<<navLayers.size()<<
" compatible DetLayers for the Seed";
166 if (navLayers.empty())
return trajL;
168 vector<DetWithState> detsWithStates;
170 for ( vector<const DetLayer*>::const_iterator layer = navLayers.begin();
171 layer != navLayers.end(); layer++) {
177 LogTrace(
category_) <<
"Number of compatible dets: " << detsWithStates.size() << endl;
179 if ( !detsWithStates.empty() ) {
181 if ( detsWithStates.front().second.isValid() ) {
185 lastTsos = detsWithStates.front().second;
190 detsWithStates.clear();
191 if ( !lastTsos.
isValid() )
return trajL;
212 int DTChamberUsedBack = 0;
213 int CSCChamberUsedBack = 0;
214 int RPCChamberUsedBack = 0;
215 int TotalChamberUsedBack = 0;
217 vector<TrajectoryMeasurement> measL;
221 for ( vector<const DetLayer*>::const_iterator rnxtlayer = navLayers.begin(); rnxtlayer!= navLayers.end(); ++rnxtlayer) {
230 if ( measL.empty() && (fabs(
theService->magneticField()->inTesla(
GlobalPoint(0,0,0)).z()) < 0.01) && (
theService->propagator(
"StraightLinePropagator").isValid() ) ) {
234 if ( measL.empty() )
continue;
236 for (vector<TrajectoryMeasurement>::const_iterator theMeas = measL.begin(); theMeas != measL.end(); ++theMeas) {
241 incrementChamberCounters((*rnxtlayer), DTChamberUsedBack, CSCChamberUsedBack, RPCChamberUsedBack, TotalChamberUsedBack);
242 secondLast = lastTsos;
243 if ( (!theTraj->
empty()) && result.second.isValid() ) {
244 lastTsos = result.second;
246 }
else if ((theMeas)->predictedState().isValid()) lastTsos = (theMeas)->predictedState();
251 while (!theTraj->
empty()) {
255 if (!theTraj->
isValid() || TotalChamberUsedBack < 2 || (DTChamberUsedBack+CSCChamberUsedBack) == 0 || !lastTsos.
isValid()) {
263 DTChamberUsedBack = 0;
264 CSCChamberUsedBack = 0;
265 RPCChamberUsedBack = 0;
266 TotalChamberUsedBack = 0;
276 if ( lastPos.basicVector().dot(momDir.
basicVector()) > 0 ) {
286 std::reverse(navLayers.begin(), navLayers.end());
291 LogTrace(
category_) <<
"Begin backward refitting, with " << navLayers.size() <<
" layers" << endl;
293 for (vector<const DetLayer*>::const_iterator rnxtlayer = navLayers.begin();
294 rnxtlayer!= navLayers.end(); ++rnxtlayer) {
300 if ( measL.empty() ) {
302 for (MuonRecHitContainer::const_iterator ihit = tmpHits.begin();
303 ihit != tmpHits.end(); ++ihit ) {
304 allUnusedHits.push_back(*ihit);
309 for (vector<TrajectoryMeasurement>::const_iterator theMeas = measL.begin(); theMeas != measL.end(); ++theMeas) {
312 if (rnxtlayer != navLayers.begin()) {
314 vector<const DetLayer*>::const_iterator lastlayer = rnxtlayer;
317 if ( (*rnxtlayer)->location() != (*lastlayer)->location() ) {
320 GlobalPoint thisPos = (theMeas)->predictedState().globalPosition();
324 if ( momDir.
mag() > 0.01 ) {
332 (lastTsos.
globalPosition().
z() * (theMeas)->predictedState().globalPosition().z() < 0) ) {
340 pair<bool,TrajectoryStateOnSurface> bkresult
343 if (bkresult.first ) {
345 incrementChamberCounters((*rnxtlayer), DTChamberUsedBack, CSCChamberUsedBack, RPCChamberUsedBack, TotalChamberUsedBack);
350 allUnusedHits.insert(allUnusedHits.end(),tmpUnusedHits.begin(),tmpUnusedHits.end());
352 if ( (!myTraj.
empty()) && bkresult.second.isValid() )
353 lastTsos = bkresult.second;
354 else if ((theMeas)->predictedState().isValid())
355 lastTsos = (theMeas)->predictedState();
360 for ( vector<Trajectory*>::iterator
t = trajL.begin();
t != trajL.end(); ++
t )
delete *
t;
364 if (( !myTraj.
isValid() ) || ( myTraj.
empty() ) || ( (
selfDuplicate(myTraj)) )|| TotalChamberUsedBack < 2 || (DTChamberUsedBack+CSCChamberUsedBack) < 1) {
381 if ( myTraj.
empty() )
return trajL;
386 if ( !smoothed.empty() && smoothed.front().foundHits()> 3 ) {
388 myTraj = smoothed.front();
417 if ( !myTraj.
isValid() )
return trajL;
437 if (dir != propDir ) {
441 if ( myTraj.
empty() )
return trajL;
457 for (MuonRecHitContainer::const_iterator ihit = tmpHits.begin();
458 ihit != tmpHits.end(); ++ihit ) {
459 if ((*ihit)->geographicalId() != meas.
recHit()->geographicalId() ){
460 result.push_back(*ihit);
482 vector<const DetLayer*> navLayers;
509 if (navLayers.empty())
return;
513 int DTChamberUsedBack = 0;
514 int CSCChamberUsedBack = 0;
515 int RPCChamberUsedBack = 0;
516 int TotalChamberUsedBack = 0;
529 vector<TrajectoryMeasurement> measL;
530 for (vector<const DetLayer*>::const_iterator rnxtlayer = navLayers.begin();
531 rnxtlayer!= navLayers.end(); ++rnxtlayer) {
536 if ( measL.empty() )
continue;
538 for (vector<TrajectoryMeasurement>::const_iterator theMeas = measL.begin(); theMeas != measL.end(); ++theMeas) {
540 pair<bool,TrajectoryStateOnSurface> bkresult
542 if (bkresult.first ) {
545 incrementChamberCounters((*rnxtlayer), DTChamberUsedBack, CSCChamberUsedBack, RPCChamberUsedBack, TotalChamberUsedBack);
547 if ( (!traj.
empty()) && bkresult.second.isValid() )
548 lastTsos = bkresult.second;
549 else if ((theMeas)->predictedState().isValid())
550 lastTsos = (theMeas)->predictedState();
576 if (traj.
empty())
return;
602 if ( predTsos.isValid() )
603 LogTrace(
category_)<<
"intermediateState: a intermediate state: pos: "<<predTsos.globalPosition() <<
"mom: " << predTsos.globalMomentum();
615 if ( hits.size() < 2 )
return;
618 vector<bool>
keep(hits.size(),
true);
622 for (MuonRecHitContainer::const_iterator ihit = hits.begin();
623 ihit != hits.end(); ++ihit ) {
624 if ( !
keep[i] ) { i++;
continue; };
626 for (MuonRecHitContainer::const_iterator ihit2 = ihit + 1;
627 ihit2 != hits.end(); ++ihit2 ) {
628 if ( !
keep[j] ) { j++;
continue; }
629 if ((*ihit)->geographicalId() == (*ihit2)->geographicalId() ) {
630 if ( (*ihit)->dimension() > (*ihit2)->dimension() ) {
632 }
else if ( (*ihit)->dimension() < (*ihit2)->dimension() ) {
635 if ( (*ihit)->transientHits().size()>(*ihit2)->transientHits().size() ) {
637 }
else if ( (*ihit)->transientHits().size()<(*ihit2)->transientHits().size() ) {
640 else if ( (*ihit)->degreesOfFreedom() != 0 && (*ihit2)->degreesOfFreedom() != 0) {
641 if (((*ihit)->chi2()/(*ihit)->degreesOfFreedom()) > ((*ihit2)->chi2()/(*ihit)->degreesOfFreedom()))
keep[
i] =
false;
642 else keep[
j] =
false;
652 for (MuonRecHitContainer::const_iterator ihit = hits.begin();
653 ihit != hits.end(); ++ihit ) {
654 if (
keep[i] ) tmp.push_back(*ihit);
673 if (traj.
empty())
return true;
676 for (ConstRecHitContainer::const_iterator ir = hits.begin(); ir != hits.end(); ir++ ) {
677 if ( !(*ir)->isValid() )
continue;
678 for (ConstRecHitContainer::const_iterator ir2 = ir+1; ir2 != hits.end(); ir2++ ) {
679 if ( !(*ir2)->isValid() )
continue;
680 if ( (*ir) == (*ir2) ) result =
true;
700 const std::vector<TrajectoryMeasurement>& meas = traj.
measurements();
702 for (std::vector<TrajectoryMeasurement>::const_reverse_iterator itm = meas.rbegin();
703 itm != meas.rend(); ++itm ) {
721 std::reverse(hits.begin(), hits.end());
728 if ( refittedback.empty() ) {
732 LogTrace(
category_) <<
"flipTrajectory: first "<< refittedback.front().firstMeasurement().updatedState()
733 <<
"\nflipTrajectory: last "<<refittedback.front().lastMeasurement().updatedState();
735 traj = refittedback.front();
750 const std::vector<TrajectoryMeasurement>& meas = traj.
measurements();
752 for (std::vector<TrajectoryMeasurement>::const_iterator itm = meas.begin(); itm != meas.end(); ++itm) {
756 while (!traj.
empty()) {
793 if ( !refitted.empty() ) traj = refitted.front();
796 std::reverse(hits.begin(), hits.end());
799 if ( !refittedback.empty() ) traj = refittedback.front();
815 for (TransientTrackingRecHit::ConstRecHitContainer::const_iterator ir = hits.begin(); ir != hits.end(); ir++ ) {
816 if ( !(*ir)->isValid() ) {
824 <<
"radius " << pos.
perp()
825 <<
" dim " << (*ir)->dimension()
826 <<
" det " << (*ir)->det()->geographicalId().det()
827 <<
" sub det " << (*ir)->det()->subDetector()<<endl;
829 if ((*ir)->det()->geographicalId().det() == 2 && (*ir)->det()->subDetector() == 6) {
834 if ((*rechit).isValid())
LogTrace(
category_)<<
"csc from collection tpeak "<<(*rechit).tpeak();
837 if ((*ir)->det()->geographicalId().det() == 2 && (*ir)->det()->subDetector() == 7) {
842 if ((*rechit).isValid())
LogTrace(
category_)<<
"dt from collection digitime "<<(*rechit).digiTime();
855 std::vector<TrajectoryMeasurement>
861 std::vector<TrajectoryMeasurement>
result;
862 std::vector<TrajectoryMeasurement> measurements;
865 std::vector<TrajectoryMeasurementGroup> measurementGroups =
868 for (std::vector<TrajectoryMeasurementGroup>::const_iterator tmGroupItr = measurementGroups.begin();
869 tmGroupItr != measurementGroups.end(); ++tmGroupItr) {
871 measurements = tmGroupItr->measurements();
875 if (bestMeasurement) result.push_back(*bestMeasurement);
883 if (bestMeasurement) result.push_back(*bestMeasurement);
885 measurements.clear();
899 int& totalChambers) {
914 if ( (dtseg == 0) || (!dtseg->
hasPhi()) )
return 0;
bool empty() const
True if trajectory has no measurements.
void rescaleError(double factor)
void reverseTrajectoryPropagationDirection(Trajectory &) const
reverse the propagation direction of a trajectory
T getParameter(std::string const &) const
const Propagator * propagator() const
double t0(const DTRecSegment4D *deseg) const
std::pair< const_iterator, const_iterator > range
iterator range
const MeasurementEstimator * estimator() const
accasso at the propagator
std::string dumpLayer(const DetLayer *layer) const
TrajectorySeed const & seed() const
Access to the seed used to reconstruct the Trajectory.
std::vector< const DetLayer * > compatibleEndcapLayers(const FreeTrajectoryState &fts, PropagationDirection timeDirection) const
CosmicMuonSmoother * theSmoother
void estimateDirection(Trajectory &) const
check the direction of trajectory by checking eta spread
void makeFirstTime()
reset the theFirstTSOSFlag
TrajectoryStateOnSurface intermediateState(const TrajectoryStateOnSurface &) const
edm::Handle< CSCRecHit2DCollection > cschits_
const DTChamberRecSegment2D * phiSegment() const
The superPhi segment: 0 if no phi projection available.
unsigned long long theCacheId_DG
Global3DPoint GlobalPoint
virtual SubDetector subDetector() const =0
The type of detector (PixelBarrel, PixelEndcap, TIB, TOB, TID, TEC, CSC, DT, RPCBarrel, RPCEndcap)
MuonTrajectoryUpdator * backwardUpdator() const
GlobalPoint globalPosition() const
std::vector< Trajectory * > trajectories(const TrajectorySeed &)
build trajectories from seed
ConstRecHitPointer recHit() const
MuonCandidate::TrajectoryContainer TrajectoryContainer
std::string thePropagatorName
edm::Handle< DTRecHitCollection > dthits_
const Propagator * propagatorOpposite() const
MuonDetLayerMeasurements * theLayerMeasurements
bool theTraversingMuonFlag
DirectMuonNavigation * theNavigation
std::vector< TrajectoryMeasurementGroup > groupedMeasurements(const DetLayer *layer, const TrajectoryStateOnSurface &startingState, const Propagator &prop, const MeasurementEstimator &est, const edm::Event &iEvent)
std::vector< const DetLayer * > compatibleLayers(const FreeTrajectoryState &fts, PropagationDirection timeDirection) const
ConstRecHitContainer recHits(bool splitting=false) const
std::string dumpMuonId(const DetId &id) const
PropagationDirection const & direction() const
C::const_iterator const_iterator
constant access iterator type
DataContainer const & measurements() const
DirectMuonNavigation * navigation() const
void flipTrajectory(Trajectory &) const
flip a trajectory with refit (the momentum direction is opposite)
FreeTrajectoryState * freeState(bool withErrors=true) const
MuonTrajectoryUpdator * updator() const
MuonRecHitContainer recHits(const DetLayer *layer, const edm::Event &iEvent)
returns the rechits which are on the layer
void reverseTrajectory(Trajectory &) const
reverse a trajectory without refit (out the measurements order changed)
MuonBestMeasurementFinder * theBestMeasurementFinder
void setFitDirection(NavigationDirection fitDirection)
set fit direction
TrajectoryMeasurement const & lastMeasurement() const
MuonTrajectoryUpdator * theBKUpdator
TrajectoryStateOnSurface updatedState() const
void incrementChamberCounters(const DetLayer *layer, int &dtChambers, int &cscChambers, int &rpcChambers, int &totalChambers)
bool hasPhi() const
Does it have the Phi projection?
TrajectoryStateOnSurface predictedState() const
How EventSelector::AcceptEvent() decides whether to accept an event for output otherwise it is excluding the probing of A single or multiple positive and the trigger will pass if any such matching triggers are PASS or EXCEPTION[A criterion thatmatches no triggers at all is detected and causes a throw.] A single negative with an expectation of appropriate bit checking in the decision and the trigger will pass if any such matching triggers are FAIL or EXCEPTION A wildcarded negative criterion that matches more than one trigger in the trigger but the state exists so we define the behavior If all triggers are the negative crieriion will lead to accepting the event(this again matches the behavior of"!*"before the partial wildcard feature was incorporated).The per-event"cost"of each negative criterion with multiple relevant triggers is about the same as!*was in the past
MeasurementContainer measurements(const DetLayer *layer, const GeomDet *det, const TrajectoryStateOnSurface &stateOnDet, const MeasurementEstimator &est, const edm::Event &iEvent)
GlobalVector momentum() const
std::vector< DTRecHit1D > specificRecHits() const
Access to specific components.
virtual TrajectoryStateOnSurface propagate(const FreeTrajectoryState &, const Surface &) const
void getDirectionByTime(Trajectory &) const
check the direction of trajectory by checking the timing
std::vector< ConstRecHitPointer > ConstRecHitContainer
GlobalPoint position() const
PropagationDirection checkDirectionByT0(const DTRecSegment4D *, const DTRecSegment4D *) const
TrajectoryMeasurement const & firstMeasurement() const
PTrajectoryStateOnDet const & startingState() const
std::vector< Trajectory > fit(const Trajectory &) const
edm::ParameterSet theNavigationPSet
virtual bool hasGroups() const =0
MuonTrajectoryUpdator * theUpdator
CosmicMuonUtilities * utilities() const
const MuonServiceProxy * theService
virtual ~CosmicMuonTrajectoryBuilder()
Destructor.
T dot(const Basic3DVector &v) const
Scalar product, or "dot" product, with a vector of same type.
void selectHits(MuonTransientTrackingRecHit::MuonRecHitContainer &) const
std::vector< std::vector< double > > tmp
void setEvent(const edm::Event &)
set event
GlobalVector globalMomentum() const
const Propagator * propagatorAlong() const
TrajectoryMeasurement * findBestMeasurement(std::vector< TrajectoryMeasurement > &measC, const Propagator *propagator)
return the Tm with the best chi2: no cut applied.
MuonTransientTrackingRecHit::MuonRecHitContainer MuonRecHitContainer
perl if(1 lt scalar(@::datatypes))
std::vector< TrajectoryMeasurement > findBestMeasurements(const DetLayer *, const TrajectoryStateOnSurface &, const Propagator *, const MeasurementEstimator *)
CosmicMuonTrajectoryBuilder(const edm::ParameterSet &, const MuonServiceProxy *service)
Constructor.
virtual std::pair< bool, TrajectoryStateOnSurface > update(const TrajectoryMeasurement *measurement, Trajectory &trajectory, const Propagator *propagator)
update the Trajectory with the TrajectoryMeasurement
void buildSecondHalf(Trajectory &)
virtual std::vector< Trajectory > trajectories(const Trajectory &) const
bool selfDuplicate(const Trajectory &) const
check if the trajectory iterates the same hit more than once
MuonTransientTrackingRecHit::MuonRecHitContainer unusedHits(const DetLayer *, const TrajectoryMeasurement &) const
double t0() const
Get the segment t0 (if recomputed, 0 is returned otherwise)
void reverseDirection(TrajectoryStateOnSurface &, const MagneticField *) const
void build(const TrajectoryStateOnSurface &, const NavigationDirection &, Trajectory &)
const BasicVectorType & basicVector() const
std::vector< MuonRecHitPointer > MuonRecHitContainer
double chiSquared() const
GlobalVector globalDirection() const
virtual void setEvent(const edm::Event &)
pass the Event to the algo at each event
const double par[8 *NPar][4]