|
|
|
#include <RecoMuon/L3TrackFinder/interface/MuonRoadTrajectoryBuilder.h>
Definition at line 55 of file MuonRoadTrajectoryBuilder.h.
typedef std::list< trajectory > MuonRoadTrajectoryBuilder::TrajectoryCollection [private] |
Definition at line 103 of file MuonRoadTrajectoryBuilder.h.
typedef flippingPair<TrajectoryCollection> MuonRoadTrajectoryBuilder::TrajectoryCollectionFPair [private] |
Definition at line 104 of file MuonRoadTrajectoryBuilder.h.
| MuonRoadTrajectoryBuilder::MuonRoadTrajectoryBuilder | ( | const edm::ParameterSet & | par, | |
| const MeasurementTracker * | mt, | |||
| const MagneticField * | f, | |||
| const Propagator * | p | |||
| ) |
constructor from PSet and things from record
Definition at line 40 of file MuonRoadTrajectoryBuilder.cc.
References Chi2MeasurementEstimatorESProducer_cfi::Chi2MeasurementEstimator, edm::ParameterSet::getParameter(), KFTrajectorySmootherESProducer_cfi::KFTrajectorySmoother, funct::sqrt(), theBranchonfirstlayer, theCarriedIPatfirstlayer, theCarriedIPatfirstlayerModule, theCategory, theDynamicMaxNumberOfHitPerModule, theField, theHitEstimator, theMaxTrajectories, theMaxTrajectoriesThreshold, theMeasurementTracker, theMinNumberOfHitOnCandidate, theNumberOfHitPerModuleDefault, theNumberOfHitPerModuleThreshold, theOutputAllTraj, thePropagator, theRoadEstimator, theSmoother, and theUpdator.
00043 : 00044 theRoadEstimator(0),theHitEstimator(0),theUpdator(0),theSmoother(0) 00045 { 00046 00047 theMeasurementTracker = mt; 00048 theField = f; 00049 thePropagator = p; 00050 00051 theCategory = "Muon|RecoMuon|MuonRoadTrajectoryBuilder"; 00052 00053 //get parameters from ParameterSet 00054 00055 //propagator name to get from muon service 00056 // thePropagatorName = iConfig.getParameter<std::string>("propagatorName"); 00057 00058 //chi2 estimator 00059 double chi2R=iConfig.getParameter<double>("maxChi2Road"); 00060 theRoadEstimator = new Chi2MeasurementEstimator(chi2R,sqrt(chi2R)); 00061 double chi2H=iConfig.getParameter<double>("maxChi2Hit"); 00062 theHitEstimator = new Chi2MeasurementEstimator(chi2H,sqrt(chi2H)); 00063 00064 //trajectory updator 00065 theUpdator= new KFUpdator(); 00066 00067 //limit the total number of possible trajectories taken into account for a single seed 00068 theMaxTrajectories = iConfig.getParameter<uint>("maxTrajectories"); 00069 00070 //limit the type of module considered to gather rechits 00071 theDynamicMaxNumberOfHitPerModule = iConfig.getParameter<bool>("dynamicMaxNumberOfHitPerModule"); 00072 theNumberOfHitPerModuleDefault = iConfig.getParameter<uint>("numberOfHitPerModule"); 00073 theMaxTrajectoriesThreshold = iConfig.getParameter<std::vector<uint> >("maxTrajectoriesThreshold"); 00074 theNumberOfHitPerModuleThreshold = iConfig.getParameter<std::vector<uint> >("numberOfHitPerModuleThreshold"); 00075 00076 //could be configurable, but not 00077 theBranchonfirstlayer=true; 00078 theCarriedIPatfirstlayer=false; 00079 theCarriedIPatfirstlayerModule =true; 00080 00081 //output track candidate selection 00082 theMinNumberOfHitOnCandidate = iConfig.getParameter<uint>("minNumberOfHitOnCandidate"); 00083 00084 //single or multiple trajectories per muon 00085 theOutputAllTraj = iConfig.getParameter<bool>("outputAllTraj"); 00086 00087 if (!theSmoother) 00088 theSmoother = new KFTrajectorySmoother(thePropagator,theUpdator,theRoadEstimator); 00089 }
| MuonRoadTrajectoryBuilder::~MuonRoadTrajectoryBuilder | ( | ) |
Definition at line 92 of file MuonRoadTrajectoryBuilder.cc.
References theCategory, theHitEstimator, theRoadEstimator, theSmoother, and theUpdator.
00093 { 00094 edm::LogInfo(theCategory)<<"cleaning the object"; 00095 if (theRoadEstimator) delete theRoadEstimator; 00096 if (theHitEstimator) delete theHitEstimator; 00097 if (theUpdator) delete theUpdator; 00098 if (theSmoother) delete theSmoother; 00099 }
| void MuonRoadTrajectoryBuilder::checkDuplicate | ( | TrajectoryCollection & | collection | ) | const [private] |
Definition at line 688 of file MuonRoadTrajectoryBuilder.cc.
References h1, h2, LogDebug, and theCategory.
Referenced by GatherHits().
00689 { 00690 LogDebug(theCategory)<<"checking for duplicates here. list size: "<<collection.size(); 00691 00692 //loop over the trajectory list 00693 TrajectoryCollection::iterator traj1 = collection.begin(); 00694 while(traj1 != collection.end()) 00695 { 00696 LogDebug(theCategory)<<""; 00697 00698 //reloop over the trajectory list from traj1 00699 TrajectoryCollection::iterator traj2 = traj1; 00700 traj2++;//advance one more 00701 bool traj1_removed=false; 00702 while( traj2 != collection.end()) 00703 { 00704 if (traj2 == traj1 ) continue; //skip itself of course 00705 00706 //need to start from the back of the list of measurment 00707 std::list <TrajectoryMeasurement >::reverse_iterator h1 = traj1->measurements.rbegin(); 00708 std::list <TrajectoryMeasurement >::reverse_iterator h2 = traj2->measurements.rbegin(); 00709 00710 bool break_different = false; 00711 while (h1 != traj1->measurements.rend() && h2!=traj2->measurements.rend()) 00712 { 00713 TransientTrackingRecHit::ConstRecHitPointer hit1 = h1->recHit(); 00714 TransientTrackingRecHit::ConstRecHitPointer hit2 = h2->recHit(); 00715 00716 LogDebug(theCategory)<<"(hit1->geographicalId().rawId()) "<<(hit1->geographicalId().rawId())<<"(hit1->globalPosition()) "<<(hit1->globalPosition()) 00717 <<"(hit2->geographicalId().rawId()) "<<(hit2->geographicalId().rawId())<<"(hit2->globalPosition()) "<<(hit2->globalPosition()); 00718 00719 if (hit1 == hit2){/*those are common hits, everything's alright so far*/ h1++;h2++; continue;} 00720 else{break_different =true; 00721 00722 LogDebug(theCategory)<<"list of hits are different"; 00723 00724 break;} 00725 } 00726 if (!break_different) 00727 //meaning one of the list has been exhausted 00728 //one list if the subset of the other 00729 { 00730 LogDebug(theCategory)<<"list of hits are identical/subset. remove one of them."; 00731 //there is a common subset to the two list of rechits. 00732 //remove the one with the fewer hits (iterator that reached the end first) 00733 //in case they are exactly identical (both iterator reached the end at the same time), traj2 will be removed by default 00734 if (h1 != traj1->measurements.rend()) 00735 { 00736 LogDebug(theCategory)<<"I removed traj2"; 00737 //traj2 has been exhausted first. remove it and place the iterator on next item 00738 traj2=collection.erase(traj2); 00739 } 00740 else 00741 { 00742 LogDebug(theCategory)<<"I removed traj1. and decrement traj2"; 00743 //traj1 has been exhausted first. remove it 00744 traj1=collection.erase(traj1); 00745 //and set the iterator traj1 so that next++ will set it to the correct place in the list 00746 traj1_removed=true; 00747 break; // break the traj2 loop, advance to next traj1 item 00748 } 00749 } 00750 else 00751 {traj2++; } 00752 }//loop on traj2 00753 if (!traj1_removed) 00754 {//increment only if you have remove the item at traj1 00755 traj1++;} 00756 }//loop on traj1 00757 }
| bool MuonRoadTrajectoryBuilder::checkStep | ( | TrajectoryCollection & | collection | ) | const [private] |
Definition at line 665 of file MuonRoadTrajectoryBuilder.cc.
References theCategory, theDynamicMaxNumberOfHitPerModule, theMaxTrajectories, theMaxTrajectoriesThreshold, theNumberOfHitPerModule, theNumberOfHitPerModuleThreshold, and trajectoryOrder().
Referenced by makeTrajectories_0().
00666 { 00667 //dynamic cut on the max number of rechit allowed on a single module 00668 if (theDynamicMaxNumberOfHitPerModule) { 00669 for (uint vit = 0; vit!= theMaxTrajectoriesThreshold.size() ; vit++){ 00670 if (collection.size() >theMaxTrajectoriesThreshold[vit]){ 00671 theNumberOfHitPerModule= theNumberOfHitPerModuleThreshold[vit];} 00672 else 00673 break;}} 00674 00675 //reduce the number of possible trajectories if too many 00676 if ( collection.size() > theMaxTrajectories) { 00677 //order with most hits or best chi2 00678 collection.sort(trajectoryOrder); 00679 uint prevSize=collection.size(); 00680 collection.resize(theMaxTrajectories); 00681 edm::LogInfo(theCategory) 00682 <<" too many possible trajectories ("<<prevSize 00683 <<"), reduce the possibilities to "<<theMaxTrajectories<<" bests."; 00684 } 00685 return true; 00686 }
| void MuonRoadTrajectoryBuilder::cleanTrajectory | ( | Trajectory & | traj | ) | const [private] |
Definition at line 394 of file MuonRoadTrajectoryBuilder.cc.
References GeometricSearchTracker::detLayer(), Trajectory::direction(), MeasurementTracker::geometricSearchTracker(), LogDebug, Trajectory::measurements(), Trajectory::seed(), theCategory, and theMeasurementTracker.
Referenced by smooth().
00394 { 00395 //remove the overlapping recHits since the smoother will chock on it 00396 Trajectory::DataContainer meas =traj.measurements(); 00397 00398 const DetLayer* lastDetLayer =0; 00399 Trajectory::DataContainer::iterator mit = meas.begin(); 00400 while (mit!=meas.end()){ 00401 { 00402 const DetLayer * detLayer = theMeasurementTracker->geometricSearchTracker()->detLayer(mit->recHit()->det()->geographicalId()); 00403 LogDebug(theCategory)<<"(mit->recHit()->det()->geographicalId().rawId()) "<<(mit->recHit()->det()->geographicalId().rawId()) 00404 <<"(detLayer) "<<(detLayer)<<"(lastDetLayer) "<<(lastDetLayer); 00405 00406 if (detLayer==lastDetLayer) 00407 { 00408 mit=meas.erase(mit); //serve as mit++ too 00409 } 00410 else {mit++;} 00411 lastDetLayer=detLayer; 00412 } 00413 Trajectory newTraj(traj.seed(),traj.direction()); 00414 for (Trajectory::DataContainer::iterator mit=meas.begin();mit!=meas.end();++mit) 00415 {newTraj.push(*mit);} 00416 traj=newTraj; 00417 } 00418 }
| int MuonRoadTrajectoryBuilder::GatherHits | ( | const TrajectoryStateOnSurface & | step, | |
| const DetLayer * | thislayer, | |||
| TrajectoryCollectionFPair & | Trajectories | |||
| ) | const [private] |
Definition at line 459 of file MuonRoadTrajectoryBuilder.cc.
References checkDuplicate(), MuonRoadTrajectoryBuilder::trajectory::chi2, GeometricSearchDet::compatibleDets(), MuonRoadTrajectoryBuilder::trajectory::duplicate, dws, Chi2MeasurementEstimator::estimate(), MuonRoadTrajectoryBuilder::flippingPair< A >::flip(), Trajectory::foundHits(), TrajectoryStateOnSurface::freeState(), TrajectoryStateOnSurface::globalPosition(), MuonRoadTrajectoryBuilder::flippingPair< A >::head(), MeasurementTracker::idToDet(), TrajectoryStateOnSurface::isValid(), MuonRoadTrajectoryBuilder::trajectory::lastmissed, LogDebug, PV3DBase< T, PVType, FrameType >::mag(), MuonRoadTrajectoryBuilder::trajectory::measurements, MuonRoadTrajectoryBuilder::trajectory::missed, MuonRoadTrajectoryBuilder::trajectory::missedinarow, PV3DBase< T, PVType, FrameType >::perp(), Propagator::propagate(), Trajectory::push(), DetId::rawId(), MeasurementDet::recHits(), MuonRoadTrajectoryBuilder::flippingPair< A >::tail(), theBranchonfirstlayer, theCarriedIPatfirstlayer, theCarriedIPatfirstlayerModule, theCategory, theFirstlayer, theHitEstimator, theMeasurementTracker, theNumberOfHitPerModule, thePropagator, theRoadEstimator, theUpdator, MuonRoadTrajectoryBuilder::trajectory::traj, MuonRoadTrajectoryBuilder::trajectory::TSOS, and TrajectoryStateUpdator::update().
Referenced by makeTrajectories_0().
00460 { 00461 TrajectoryStateOnSurface restep; 00462 bool atleastoneadded=false; 00463 int Nhits=0; 00464 00465 //find compatible modules 00466 std::vector<DetLayer::DetWithState> compatible =thislayer->compatibleDets( step, *thePropagator , *theRoadEstimator); 00467 00468 //loop over compatible modules 00469 for (std::vector< DetLayer::DetWithState > ::iterator dws = compatible.begin();dws != compatible.end();dws++) 00470 { 00471 const DetId presentdetid = dws->first->geographicalId();//get the det Id 00472 restep = dws->second; 00473 00474 00475 LogDebug(theCategory)<<((dws->first->components ().size()!=0) ? /*stereo layer*/"double sided layer":/*single sided*/"single sided layer") 00476 <<"(presentdetid.rawId()) "<<(presentdetid.rawId()); 00477 00478 //get the rechits on this module 00479 TransientTrackingRecHit::ConstRecHitContainer thoseHits = theMeasurementTracker->idToDet(presentdetid)->recHits(restep); 00480 int additionalHits =0; 00481 00482 LogDebug(theCategory)<<"(thoseHits.size()) "<<(thoseHits.size()); 00483 00484 if (thoseHits.size()>theNumberOfHitPerModule) 00485 { edm::LogInfo(theCategory)<<"more than "<<theNumberOfHitPerModule 00486 <<" compatible hits ("<<thoseHits.size()<<")on module "<<presentdetid.rawId() 00487 <<", skip it";continue; } 00488 00489 //loop over the rechit on the module 00490 for (TransientTrackingRecHit::ConstRecHitContainer::const_iterator iTThit = thoseHits.begin();iTThit != thoseHits.end(); iTThit++) 00491 { 00492 if (!(*iTThit)->isValid())/*skip*/{edm::LogInfo(theCategory)<<"rec hit not valid on module "<<presentdetid.rawId() <<". I skip it"; continue;} 00493 00494 LogDebug(theCategory)<<((dynamic_cast<const SiStripMatchedRecHit2D *>((*iTThit)->hit())!=0) ? /*matched rechit*/ "matched rechit" : " r-phi rechit"); 00495 00496 //get the surface of the module Id 00497 const BoundPlane & surface = (*iTThit)->det()->surface(); 00498 00499 //estimate the consistency 00500 MeasurementEstimator::HitReturnType est_road = theRoadEstimator->estimate(restep,**iTThit); 00501 00502 LogDebug(theCategory)<<"(restep.globalPosition().perp()) "<<(restep.globalPosition().perp())<<"(restep.globalPosition().mag()) "<<(restep.globalPosition().mag()) 00503 <<"road step parameters at module: \n"<<restep 00504 << "(est_road.first) "<<(est_road.first)<<"(est_road.second) "<<(est_road.second); 00505 00506 //check consistency 00507 if (est_road.first) 00508 { //hit and propagation are consistent : gather the hit 00509 Nhits++; 00510 additionalHits++; 00511 atleastoneadded=true; 00512 00513 LogDebug(theCategory)<<"hit is consistent with road"<<"(presentdetid.rawId()) "<<(presentdetid.rawId()) 00514 <<"loop over previous trajectories\n" 00515 <<"(Trajectories.tail().size()) "<<(Trajectories.tail().size())<<"(Trajectories.head().size()) "<<(Trajectories.head().size()); 00516 00517 //update the list of trajectory that we have with the consistent rechit. 00518 //loop over the existing list of trajectories and add the hit if necessary 00519 for ( TrajectoryCollection::iterator traj =Trajectories.head().begin();traj != Trajectories.head().end(); traj++) 00520 { 00521 //what is the previous state for this trajectory 00522 const TrajectoryStateOnSurface & previousState = traj->TSOS; 00523 if (!previousState.isValid()) 00524 {edm::LogError(theCategory)<<"previous free state is invalid at trajectory update, this is WRONG"; continue;} 00525 const FreeTrajectoryState & previousFreestate = *previousState.freeState(); 00526 00527 //propagate it to the current surface 00528 TrajectoryStateOnSurface predictedState = thePropagator->propagate(previousFreestate,surface); 00529 if (!predictedState.isValid())/*skip*/{edm::LogError(theCategory) 00530 <<"predicted state is not valid at trajectory update, rechit surface cannot be reached by the previous updated state.";continue;} 00531 00532 MeasurementEstimator::HitReturnType est ; 00533 est= theHitEstimator->estimate(predictedState,**iTThit); 00534 00535 LogDebug(theCategory)<<"(est.first) "<<(est.first)<<"(est.second) "<<(est.second); 00536 00537 //is the current hit consistent with the trajectory 00538 if (est.first ) 00539 { 00540 //update the trajectory state with the rechit 00541 const TrajectoryStateOnSurface & updatedState = theUpdator->update(predictedState,**iTThit); 00542 if (!updatedState.isValid())/*skip*/{edm::LogError(theCategory)<<"updated state is not valid, this is really wrong";continue;} 00543 00544 LogDebug(theCategory)<<"updating a previous state with a rechit"; 00545 00546 //add a combined trajectory to the new list of trajectory, starting from the existing trajecotry 00547 Trajectories.tail().push_back(*traj); //from existing 00548 trajectory & combined = (*Trajectories.tail().rbegin()); 00549 combined.duplicate =false; //this is important 00550 //increment the chi2 00551 //combined.chi2 += est.second; 00552 combined.chi2 += est_road.second;//better to add the road chi2 too unbias the chi2 sum towards first hits 00553 //some history about the trajectory 00554 combined.lastmissed=false; 00555 combined.missedinarow=0; 00556 //add a new hits to the measurements 00557 combined.measurements.push_back(TrajectoryMeasurement(updatedState,*iTThit)); 00558 TrajectoryMeasurement & trajMeasurement = (*combined.measurements.rbegin()); 00559 //assigne updated state 00560 combined.TSOS = updatedState; 00561 //add trajectory measurement 00562 combined.traj.push(trajMeasurement,est.second); 00563 00564 LogDebug(theCategory)<<"(combined.traj.foundHits()) "<<(combined.traj.foundHits()) 00565 <<"muon measurement on previous module: \n"<<traj->TSOS 00566 <<"muon measurement before update: \n"<<predictedState 00567 <<"muon measurement after update: \n"<<updatedState; 00568 00569 }//hit consistent with trajectory 00570 else 00571 { 00572 LogDebug(theCategory)<<"hit failed chi2 test for trajectories update\n"<<"(traj->duplicate) "<<(traj->duplicate); 00573 00574 if (!traj->duplicate){ 00575 Trajectories.tail().push_back(*traj); 00576 trajectory & copied = (*Trajectories.tail().rbegin()); 00577 copied.missed++; 00578 copied.lastmissed=true; 00579 copied.missedinarow++; 00580 traj->duplicate =true; //set this trajectory to already have been copied into next step 00581 }//not a duplicate trajectory 00582 }//hit not consistent with trajectory 00583 }//existing trajectories loop 00584 }//hit is consistent with muon road 00585 }//rechits on module loop 00586 00587 00588 //discard previous list of trajectories 00589 //if something as been done of course 00590 if (Trajectories.tail().size()!=0) 00591 { 00592 //if this is not the first "layer" of detector, set updated trajectories as new seed trajectories 00593 //will branch on every single rechit uncountered for first "layer" 00594 if (!theFirstlayer || theBranchonfirstlayer) 00595 { 00596 LogDebug(theCategory)<<"swapping trajectory list index"; 00597 00598 //always carry on the <IP> alone state in the next list of trajectories to avoid bias from the first rechits 00599 if (theFirstlayer && theCarriedIPatfirstlayerModule) 00600 { 00601 LogDebug(theCategory)<<"push front <IP> to next list of trajectories"; 00602 00603 Trajectories.tail().push_front(*Trajectories.head().begin()); //[0] is <IP> always 00604 trajectory & pushed = *Trajectories.tail().begin(); 00605 pushed.missed+=additionalHits; 00606 pushed.lastmissed = ( additionalHits!=0); 00607 pushed.missedinarow++; 00608 } 00609 //FIXME, there is a candidate leak at this point 00610 //if the IP state is carried on at first layer, without update, then it will be duplicated. this is very unlikely though 00611 00612 Trajectories.head().clear(); 00613 00614 //swap the lists 00615 Trajectories.flip(); 00616 } 00617 }//discard previous list of trajectories 00618 }//module loop 00619 00620 00621 00622 //do some cleaning of the list of trajectories 00623 if(theFirstlayer && atleastoneadded ) 00624 { 00625 theFirstlayer =false; //we are not in the first layer if something has been added 00626 if (!theBranchonfirstlayer) 00627 { 00628 LogDebug(theCategory)<<"swapping trajectory list index (end of first layer)"; 00629 00630 //and for consistency, you have to swap index here, because it could ahve not been done in the module loop above 00631 //always carry on the <IP> alone state in the next list of trajectories to avoid bias from the first rechits 00632 if (theCarriedIPatfirstlayer) 00633 {Trajectories.tail().push_front(*Trajectories.head().begin());} //[0] is <IP> always at this stage 00634 //FIXME, there is a candidate leak at this point 00635 //if the IP state is carried on at first layer, without update, then it will be duplicated. this is very unlikely though 00636 00637 Trajectories.head().clear(); 00638 //swap the switch 00639 Trajectories.flip(); 00640 } 00641 else 00642 { 00643 //actually remove the <IP> from first postion of the next source of trajectories 00644 //since the swaping has already been done. pop up from theTrajectorysource, not !theTrajectorysource 00645 //only if it has been done at the first layer module though 00646 if (!theCarriedIPatfirstlayer && theCarriedIPatfirstlayerModule){ 00647 00648 LogDebug(theCategory)<<"pop up <IP> from trajectories"; 00649 00650 Trajectories.head().pop_front(); } 00651 } 00652 00653 00654 //check an remove trajectories that are subset of the other in next source 00655 //after the first layer only 00656 if (Trajectories.head().size()>=2) 00657 {checkDuplicate(Trajectories.head());} 00658 00659 } //do some cleaning of the list of trajectories 00660 return Nhits; 00661 00662 }
| void MuonRoadTrajectoryBuilder::makeTrajectories | ( | const TrajectorySeed & | seed, | |
| std::vector< Trajectory > & | result, | |||
| int | version = 0 | |||
| ) | const [private] |
Definition at line 143 of file MuonRoadTrajectoryBuilder.cc.
References makeTrajectories_0(), and makeTrajectories_1().
Referenced by trajectories().
00144 { if (version==0) { makeTrajectories_0(seed,result);} 00145 else if (version==1) { makeTrajectories_1(seed,result);}}
| void MuonRoadTrajectoryBuilder::makeTrajectories_0 | ( | const TrajectorySeed & | seed, | |
| std::vector< Trajectory > & | result | |||
| ) | const [private] |
Definition at line 150 of file MuonRoadTrajectoryBuilder.cc.
References alongMomentum, checkStep(), PTrajectoryStateOnDet::detId(), detId, TrajectoryStateOnSurface::freeState(), GatherHits(), MeasurementTracker::geometricSearchTracker(), MeasurementTracker::geomTracker(), TrajectoryStateOnSurface::globalMomentum(), TrajectoryStateOnSurface::globalPosition(), MuonRoadTrajectoryBuilder::flippingPair< A >::head(), TrackingGeometry::idToDet(), SimpleDiskBounds::innerRadius(), TrajectoryStateOnSurface::isValid(), LogDebug, MuonRoadTrajectoryBuilder::trajectory::measurements, GeometricSearchTracker::negTecLayers(), GeometricSearchTracker::negTidLayers(), SimpleDiskBounds::outerRadius(), PV3DBase< T, PVType, FrameType >::perp(), GeometricSearchTracker::posTecLayers(), GeometricSearchTracker::posTidLayers(), Propagator::propagate(), r, DetId::rawId(), smooth(), TrajectorySeed::startingState(), DetId::subdetId(), GeomDet::surface(), GeomDetEnumerators::TEC, theCategory, theField, theFirstlayer, theMeasurementTracker, theMinNumberOfHitOnCandidate, theNumberOfHitPerModule, theNumberOfHitPerModuleDefault, theOutputAllTraj, thePropagator, theTransformer, GeomDetEnumerators::TIB, GeometricSearchTracker::tibLayers(), GeomDetEnumerators::TID, GeomDetEnumerators::TOB, GeometricSearchTracker::tobLayers(), MuonRoadTrajectoryBuilder::trajectory::traj, trajectoryOrder(), TrajectoryStateTransform::transientState(), MuonRoadTrajectoryBuilder::trajectory::TSOS, PV3DBase< T, PVType, FrameType >::z(), and z.
Referenced by makeTrajectories().
00151 { 00152 Trajectory basicTrajectory(seed,alongMomentum); 00153 //add the muon system measurement to the basicTrajectory 00154 //...no yet implemented... 00155 00156 //get the initial state 00157 PTrajectoryStateOnDet PTstart=seed.startingState(); 00158 DetId detId(PTstart.detId()); 00159 const BoundPlane * surface =(&theMeasurementTracker->geomTracker()->idToDet(detId)->surface()); 00160 //start from this point 00161 TrajectoryStateOnSurface TSOS = theTransformer.transientState(PTstart,surface,theField); 00162 if (!TSOS.isValid())/*/abort*/{ edm::LogError(theCategory)<<"TSOS from PTSOD is not valid.";return ;} 00163 00164 LogDebug(theCategory) <<"(detId.rawId()) "<<(detId.rawId())<<"(detId.subdetId()) "<<(detId.subdetId()) 00165 <<"(TSOS.globalPosition()) "<<(TSOS.globalPosition()) 00166 <<"(TSOS.globalMomentum()) "<<(TSOS.globalMomentum()); 00167 00168 //initialization 00169 //---------------------- 00170 // WARNING. this is mutable 00171 // you should never remove this line 00172 theFirstlayer=true; 00173 theNumberOfHitPerModule = theNumberOfHitPerModuleDefault; 00174 //---------------------- 00175 int Nhits=0; 00176 GlobalPoint position; 00177 double z=0,r=0; 00178 //flipping pair of list of trajectories 00179 TrajectoryCollectionFPair Trajectories; 00180 //initialize the head() of pair of list of trajectories 00181 Trajectories.head().push_back(trajectory()); 00182 trajectory & initial = *Trajectories.head().rbegin(); 00183 initial.TSOS=TSOS; 00184 initial.traj=basicTrajectory; 00185 00186 00187 00188 //get the DetLayer in the tracker 00189 std::vector<BarrelDetLayer*> tiblc = theMeasurementTracker->geometricSearchTracker()->tibLayers(); 00190 std::vector<BarrelDetLayer*> toblc = theMeasurementTracker->geometricSearchTracker()->tobLayers(); 00191 std::vector<ForwardDetLayer*> tidlc[2]; 00192 tidlc[0]=theMeasurementTracker->geometricSearchTracker()->negTidLayers(); 00193 tidlc[1]=theMeasurementTracker->geometricSearchTracker()->posTidLayers(); 00194 std::vector<ForwardDetLayer*> teclc[2]; 00195 teclc[0]=theMeasurementTracker->geometricSearchTracker()->negTecLayers(); 00196 teclc[1]=theMeasurementTracker->geometricSearchTracker()->posTecLayers(); 00197 const int Ntib=tiblc.size(); 00198 const int Ntob=toblc.size(); 00199 const int Ntid=tidlc[0].size(); 00200 00201 00202 const int Ntec=teclc[0].size(); 00203 00204 LogDebug(theCategory)<<"(Ntib) "<<(Ntib)<<"(Ntob) "<<(Ntob)<<"(Ntid) "<<(Ntid)<<"(Ntec) "<<(Ntec); 00205 00206 const SimpleDiskBounds * sdbounds=0; 00207 00208 position = TSOS.globalPosition(); 00209 z = position.z(); 00210 r = position.perp(); 00211 00212 //select which part we are in 00213 enum PART { fault , PXB, PXF, TIB , TID , TOB , TEC}; 00214 PART whichpart = fault; 00215 switch(detId.subdetId()){ 00216 case 1: {whichpart=PXB;break;} 00217 case 2: {whichpart=PXF;break;} 00218 00219 case 3: {whichpart=TIB;break;} 00220 case 4: {whichpart=TID;break;} 00221 case 5: {whichpart=TOB;break;} 00222 case 6: {whichpart=TEC;break;}} 00223 00224 bool inapart=true; 00225 bool firstRound =true; 00226 int indexinpart=0; 00227 //loop while in a valid part of the tracker 00228 while(inapart){ 00229 00230 //check on the trajectory list and stuff 00231 if (!checkStep(Trajectories.head())) break; 00232 //................ 00233 00234 switch (whichpart){ 00235 00236 case fault: /*abort*/ {edm::LogError(theCategory)<<"something's wrong with the seed";return;} 00237 case PXB: /*abort*/ {edm::LogError(theCategory)<<"PXB no yet implemented";return;} 00238 case PXF: /*abort*/ {edm::LogError(theCategory)<<"PXF no yet implemented";return;} 00239 00240 //-------------- into TIB---------------- 00241 case TIB:{ 00242 if (indexinpart==Ntib){/*you have reach the last layer of the TIB.let's go to the TOB*/whichpart = TOB; indexinpart=-1;break; } 00243 00244 LogDebug(theCategory)<<"within TIB "<<indexinpart+1; 00245 00246 //propagate to corresponding surface 00247 if (!firstRound) TSOS = thePropagator->propagate(*TSOS.freeState(),tiblc[indexinpart]->surface()); 00248 if (!TSOS.isValid()) {;break;} //go to the next one 00249 00250 z=TSOS.globalPosition().z(); 00251 00252 //have we reached a boundary 00253 if (fabs(z) > fabs(tidlc[(z>0)][0]->surface().position().z())) 00254 {/*z bigger than the TID min z: go to TID*/ 00255 LogDebug(theCategory)<<"|z| ("<<z<<") bigger than the TID min z("<<tidlc[(z>0)][0]->surface().position().z()<<"): go to TID"; 00256 whichpart = TID; indexinpart=-1; break;} 00257 else {/*gather hits in the corresponding TIB layer*/ Nhits+=GatherHits(TSOS,tiblc[indexinpart],Trajectories);} 00258 break;} 00259 00260 //-------------- into TID---------------- 00261 case TID: { 00262 if (indexinpart==Ntid){/*you have reach the last layer of the TID. let's go to the TEC */ whichpart = TEC; indexinpart=-1; break;} 00263 00264 LogDebug(theCategory)<<"within TID "<<indexinpart+1; 00265 00266 //propagate to corresponding surface 00267 if (!firstRound) TSOS = thePropagator->propagate(*TSOS.freeState(),tidlc[(z>0)][indexinpart]->surface()); 00268 if (!TSOS.isValid()){break;}//go to the next one 00269 00270 position = TSOS.globalPosition(); 00271 z = position.z(); 00272 r = position.perp(); 00273 00274 sdbounds = dynamic_cast<const SimpleDiskBounds *>(&tidlc[(z>0)][indexinpart]->surface().bounds()); 00275 if (!sdbounds)/*abort*/{edm::LogError(theCategory)<<" detlayer bounds are not SimpleDiskBounds in tid geometry";return;} 00276 00277 //have we reached a boundary 00278 if (r < sdbounds->innerRadius()) 00279 {/*radius smaller than the TID disk inner radius: next disk please*/ 00280 LogDebug(theCategory)<<"radius ("<<r<<") smaller than the TID disk inner radius ("<<sdbounds->innerRadius()<<"): next disk please"; 00281 break;} 00282 else if (r >sdbounds->outerRadius()) 00283 {/*radius bigger than the TID disk outer radius: go to TOB*/ 00284 LogDebug(theCategory)<<"radius ("<<r<<") bigger than the TID disk outer radius("<<sdbounds->outerRadius()<<"): go to TOB"; 00285 whichpart = TOB; indexinpart=-1;break;} 00286 else {/*gather hits in the corresponding TIB layer*/ 00287 LogDebug(theCategory)<<"collecting hits"; 00288 Nhits+=GatherHits(TSOS,tidlc[(z>0)][indexinpart],Trajectories);} 00289 break;} 00290 00291 //-------------- into TOB---------------- 00292 case TOB: { 00293 if (indexinpart==Ntob){/*you have reach the last layer of the TOB. this is an end*/ inapart=false;break;} 00294 00295 LogDebug(theCategory)<<"within TOB "<<indexinpart+1; 00296 00297 //propagate to corresponding surface 00298 if (!firstRound) TSOS = thePropagator->propagate(*TSOS.freeState(),toblc[indexinpart]->surface()); 00299 if (!TSOS.isValid()) {break;} //go to the next one 00300 00301 z = TSOS.globalPosition().z(); 00302 00303 //have we reached a boundary 00304 if (fabs(z) > fabs(teclc[(z>0)][0]->surface().position().z())) 00305 {/*z bigger than the TOB layer max z: go to TEC*/ 00306 LogDebug(theCategory)<<"|z| ("<<z<<") bigger than the TOB layer max z ("<< teclc[(z>0)][0]->surface().position().z()<<"): go to TEC"; 00307 whichpart = TEC; indexinpart=-1;break;} 00308 else {/*gather hits in the corresponding TOB layer*/Nhits+=GatherHits(TSOS,toblc[indexinpart],Trajectories);} 00309 break;} 00310 00311 //-------------- into TEC---------------- 00312 case TEC: { 00313 if (indexinpart==Ntec){/*you have reach the last layer of the TEC. let's end here*/inapart=false;break;} 00314 00315 LogDebug(theCategory)<<"within TEC "<<indexinpart+1; 00316 00317 //propagate to corresponding TEC disk 00318 if (!firstRound) TSOS = thePropagator->propagate(*TSOS.freeState(),teclc[(z>0)][indexinpart]->surface()); 00319 if (!TSOS.isValid()) {break;} //go to the next one 00320 00321 position = TSOS.globalPosition(); 00322 z = position.z(); 00323 r = position.perp(); 00324 00325 sdbounds = dynamic_cast<const SimpleDiskBounds *>(&teclc[(z>0)][indexinpart]->surface().bounds()); 00326 if (!sdbounds)/*abort*/ {edm::LogError(theCategory)<<" detlayer bounds are not SimpleDiskBounds in tec geometry";return;} 00327 00328 //have we reached a boundary ? 00329 if (r < sdbounds->innerRadius()) 00330 {/*radius smaller than the TEC disk inner radius: next disk please*/ 00331 LogDebug(theCategory)<<"radius ("<<r<<") smaller than the TEC disk inner radius ("<<sdbounds->innerRadius()<<"): next disk please"; 00332 break;} 00333 else if (r > sdbounds->outerRadius()) 00334 {/*radius bigger than the TEC disk outer radius: I can stop here*/ 00335 LogDebug(theCategory)<<"radius ("<<r<<") bigger than the TEC disk outer radius ("<<sdbounds->outerRadius()<<"): I can stop here"; 00336 inapart=false;break;} 00337 else {/*gather hits in the corresponding TEC layer*/Nhits+=GatherHits(TSOS,teclc[(z>0)][indexinpart],Trajectories);} 00338 00339 break;} 00340 00341 }//switch 00342 indexinpart++; 00343 firstRound=false; 00344 }//while inapart 00345 00346 00347 LogDebug(theCategory)<<"propagating through layers is done"; 00348 00349 //--------------------------------------- SWIMMING DONE -------------------------------------- 00350 //the list of trajectory has been build. let's find out which one is the best to use. 00351 00352 edm::LogInfo(theCategory)<<"found: " 00353 <<Nhits<<" hits in the road for this seed \n" 00354 << Trajectories.head().size()<<" possible trajectory(ies)"; 00355 00356 if (Trajectories.head().size() == 0) /*abort*/{edm::LogError(theCategory)<<" no possible trajectory found"; return;} 00357 00358 //order the list to put the best in front 00359 Trajectories.head().sort(trajectoryOrder); 00360 00361 //------------------------------------OUTPUT FILL--------------------------------------------- 00362 if (theOutputAllTraj) { 00363 //output all the possible trajectories found, if they have enough hits on them 00364 //the best is still the first 00365 for (TrajectoryCollection::iterator tit = Trajectories.head().begin(); tit!=Trajectories.head().end();tit++) 00366 {if (tit->measurements.size()>= theMinNumberOfHitOnCandidate){result.push_back(smooth(tit->traj));}} 00367 } 00368 else{ 00369 //output only the best one 00370 trajectory & best = (*Trajectories.head().begin()); 00371 if (best.measurements.size()< theMinNumberOfHitOnCandidate)/*abort*/{edm::LogError(theCategory)<<"best trajectory does not have enough ("<<best.measurements.size()<<") hits on it (<"<<theMinNumberOfHitOnCandidate<<")"; return;} 00372 //output only the best trajectory 00373 result.push_back(smooth(best.traj));} 00374 00375 edm::LogInfo(theCategory)<<result.size()<<" trajectory(ies) output"; 00376 }//makeTrajectories_0
| void MuonRoadTrajectoryBuilder::makeTrajectories_1 | ( | const TrajectorySeed & | seed, | |
| std::vector< Trajectory > & | result | |||
| ) | const [private] |
Definition at line 762 of file MuonRoadTrajectoryBuilder.cc.
References alongMomentum.
Referenced by makeTrajectories().
00763 { 00764 Trajectory basicTrajectory(seed,alongMomentum); 00765 //add the muon system measurement to the basicTrajectory 00766 //... 00767 00768 // //build the trajectories 00769 // std::vector<Trajectory> unsmoothed = theCkfbuilder->trajectories(seed); 00770 00771 // //smoothed them 00772 // if (theOutputAllTraj) { 00773 // for (std::vector<Trajectory>::iterator tit = unsmoothed.begin(); tit!=unsmoothed.end();tit++) 00774 // {result.push_back(smooth(*tit));} 00775 // } 00776 // else{ 00777 // //output only the first one 00778 // result.push_back(smooth(unsmoothed.front()));} 00779 00780 }
| void MuonRoadTrajectoryBuilder::setEvent | ( | const edm::Event & | iEvent | ) | const [virtual] |
Implements TrajectoryBuilder.
Definition at line 101 of file MuonRoadTrajectoryBuilder.cc.
References theMeasurementTracker, and MeasurementTracker::update().
00101 { 00102 theMeasurementTracker->update(iEvent); 00103 }
| Trajectory MuonRoadTrajectoryBuilder::smooth | ( | Trajectory & | traj | ) | const [private] |
Definition at line 440 of file MuonRoadTrajectoryBuilder.cc.
References cleanTrajectory(), sortTrajectoryMeasurements(), theCategory, theSmoother, and KFTrajectorySmoother::trajectories().
Referenced by makeTrajectories_0().
00440 { 00441 00442 //need to order the list of measurements on the trajectory first 00443 sortTrajectoryMeasurements(traj); 00444 00445 std::vector<Trajectory> ret=theSmoother->trajectories(traj); 00446 00447 if (ret.empty()){ 00448 edm::LogError(theCategory)<<"smoother returns an empty vector of trajectories: try cleaning first and resmooth"; 00449 cleanTrajectory(traj); 00450 ret=theSmoother->trajectories(traj); 00451 } 00452 00453 if (ret.empty()){edm::LogError(theCategory)<<"smoother returns an empty vector of trajectories."; 00454 return traj;} 00455 else{return (ret.front());} 00456 }
| void MuonRoadTrajectoryBuilder::trajectories | ( | const TrajectorySeed & | seed, | |
| TrajectoryContainer & | ret | |||
| ) | const |
process the seed, in a faster manner.
Definition at line 119 of file MuonRoadTrajectoryBuilder.cc.
References LogDebug, makeTrajectories(), and theCategory.
00120 { 00121 LogDebug(theCategory)<<"makeTrajectories START"; 00122 //process the seed through the tracker 00123 makeTrajectories(seed,ret); 00124 }
| std::vector< Trajectory > MuonRoadTrajectoryBuilder::trajectories | ( | const TrajectorySeed & | seed | ) | const [virtual] |
Implements TrajectoryBuilder.
Definition at line 127 of file MuonRoadTrajectoryBuilder.cc.
References LogDebug, makeTrajectories(), HLT_VtxMuL3::result, and theCategory.
00128 { 00129 LogDebug(theCategory)<<"makeTrajectories START"; 00130 00131 //default output 00132 std::vector<Trajectory> result; 00133 00134 //process the seed through the tracker 00135 makeTrajectories(seed,result); 00136 00137 //output the result of regional tracking 00138 return result; 00139 }
Definition at line 142 of file MuonRoadTrajectoryBuilder.h.
Referenced by GatherHits(), and MuonRoadTrajectoryBuilder().
Definition at line 143 of file MuonRoadTrajectoryBuilder.h.
Referenced by GatherHits(), and MuonRoadTrajectoryBuilder().
Definition at line 144 of file MuonRoadTrajectoryBuilder.h.
Referenced by GatherHits(), and MuonRoadTrajectoryBuilder().
std::string MuonRoadTrajectoryBuilder::theCategory [private] |
Info/Debug category "Muon|RecoMuon|MuonRoadTrajectoryBuilder".
Definition at line 75 of file MuonRoadTrajectoryBuilder.h.
Referenced by checkDuplicate(), checkStep(), cleanTrajectory(), GatherHits(), makeTrajectories_0(), MuonRoadTrajectoryBuilder(), smooth(), trajectories(), and ~MuonRoadTrajectoryBuilder().
Definition at line 135 of file MuonRoadTrajectoryBuilder.h.
Referenced by checkStep(), and MuonRoadTrajectoryBuilder().
const MagneticField* MuonRoadTrajectoryBuilder::theField [private] |
Definition at line 171 of file MuonRoadTrajectoryBuilder.h.
Referenced by makeTrajectories_0(), and MuonRoadTrajectoryBuilder().
bool MuonRoadTrajectoryBuilder::theFirstlayer [mutable, private] |
Definition at line 126 of file MuonRoadTrajectoryBuilder.h.
Referenced by GatherHits(), and makeTrajectories_0().
Definition at line 161 of file MuonRoadTrajectoryBuilder.h.
Referenced by GatherHits(), MuonRoadTrajectoryBuilder(), and ~MuonRoadTrajectoryBuilder().
uint MuonRoadTrajectoryBuilder::theMaxTrajectories [private] |
Definition at line 129 of file MuonRoadTrajectoryBuilder.h.
Referenced by checkStep(), and MuonRoadTrajectoryBuilder().
std::vector<uint> MuonRoadTrajectoryBuilder::theMaxTrajectoriesThreshold [private] |
Definition at line 138 of file MuonRoadTrajectoryBuilder.h.
Referenced by checkStep(), and MuonRoadTrajectoryBuilder().
const MeasurementTracker* MuonRoadTrajectoryBuilder::theMeasurementTracker [private] |
Definition at line 155 of file MuonRoadTrajectoryBuilder.h.
Referenced by cleanTrajectory(), GatherHits(), makeTrajectories_0(), MuonRoadTrajectoryBuilder(), and setEvent().
unsigned int MuonRoadTrajectoryBuilder::theMinNumberOfHitOnCandidate [private] |
Definition at line 148 of file MuonRoadTrajectoryBuilder.h.
Referenced by makeTrajectories_0(), and MuonRoadTrajectoryBuilder().
uint MuonRoadTrajectoryBuilder::theNumberOfHitPerModule [mutable, private] |
Definition at line 137 of file MuonRoadTrajectoryBuilder.h.
Referenced by checkStep(), GatherHits(), and makeTrajectories_0().
uint MuonRoadTrajectoryBuilder::theNumberOfHitPerModuleDefault [private] |
Definition at line 136 of file MuonRoadTrajectoryBuilder.h.
Referenced by makeTrajectories_0(), and MuonRoadTrajectoryBuilder().
std::vector<uint> MuonRoadTrajectoryBuilder::theNumberOfHitPerModuleThreshold [private] |
Definition at line 139 of file MuonRoadTrajectoryBuilder.h.
Referenced by checkStep(), and MuonRoadTrajectoryBuilder().
Definition at line 150 of file MuonRoadTrajectoryBuilder.h.
Referenced by makeTrajectories_0(), and MuonRoadTrajectoryBuilder().
const Propagator* MuonRoadTrajectoryBuilder::thePropagator [private] |
Definition at line 174 of file MuonRoadTrajectoryBuilder.h.
Referenced by GatherHits(), makeTrajectories_0(), and MuonRoadTrajectoryBuilder().
Definition at line 159 of file MuonRoadTrajectoryBuilder.h.
Referenced by GatherHits(), MuonRoadTrajectoryBuilder(), and ~MuonRoadTrajectoryBuilder().
Definition at line 165 of file MuonRoadTrajectoryBuilder.h.
Referenced by MuonRoadTrajectoryBuilder(), smooth(), and ~MuonRoadTrajectoryBuilder().
Definition at line 163 of file MuonRoadTrajectoryBuilder.h.
Referenced by GatherHits(), MuonRoadTrajectoryBuilder(), and ~MuonRoadTrajectoryBuilder().
1.5.4