CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_4_1_8_patch9/src/RecoMuon/MuonSeedGenerator/src/RPCCosmicSeedrecHitFinder.cc

Go to the documentation of this file.
00001 
00007 #include "RecoMuon/MuonSeedGenerator/src/RPCCosmicSeedrecHitFinder.h"
00008 #include <DataFormats/TrackingRecHit/interface/TrackingRecHit.h>
00009 #include <DataFormats/RPCRecHit/interface/RPCRecHit.h>
00010 #include <FWCore/Framework/interface/ESHandle.h>
00011 #include <Geometry/RPCGeometry/interface/RPCGeometry.h>
00012 #include <Geometry/Records/interface/MuonGeometryRecord.h>
00013 
00014 using namespace std;
00015 using namespace edm;
00016 
00017 MuonTransientTrackingRecHit::MuonRecHitContainer::const_iterator find(MuonTransientTrackingRecHit::MuonRecHitContainer::const_iterator firstIter, MuonTransientTrackingRecHit::MuonRecHitContainer::const_iterator lastIter, const MuonTransientTrackingRecHit::MuonRecHitPointer& recHitRef) {
00018 
00019     MuonTransientTrackingRecHit::MuonRecHitContainer::const_iterator index = lastIter;
00020     for(MuonTransientTrackingRecHit::MuonRecHitContainer::const_iterator iter = firstIter; iter != lastIter; iter++)
00021         if((*iter) == recHitRef)
00022             index = iter;
00023     return index;
00024 }
00025 
00026 RPCCosmicSeedrecHitFinder::RPCCosmicSeedrecHitFinder() {
00027 
00028     // Initiate the member
00029     isLayerset = false;
00030     isConfigured = false;
00031     isInputset = false;
00032     isOutputset = false;
00033     isEdgeset = false;
00034     BxRange = 0;
00035     MaxDeltaPhi = 0;
00036     ClusterSet.clear();
00037     innerBounds.clear();
00038     isLayersmixed = false;
00039     LayersinRPC.clear();
00040     therecHits.clear();
00041     isOuterLayerfilled = false;
00042 }
00043 
00044 RPCCosmicSeedrecHitFinder::~RPCCosmicSeedrecHitFinder() {
00045 
00046 }
00047 
00048 void RPCCosmicSeedrecHitFinder::configure(const edm::ParameterSet& iConfig) {
00049 
00050     // Set the configuration
00051     BxRange = iConfig.getParameter<unsigned int>("BxRange");
00052     MaxDeltaPhi = iConfig.getParameter<double>("MaxDeltaPhi");
00053     ClusterSet = iConfig.getParameter< std::vector<int> >("ClusterSet");
00054 
00055     // Set the signal open
00056     isConfigured = true;
00057 }
00058 
00059 void RPCCosmicSeedrecHitFinder::setInput(MuonRecHitContainer (&recHits)[RPCLayerNumber]) {
00060 
00061     for(unsigned int i = 0; i < RPCLayerNumber; i++) {
00062         AllrecHits[i].clear();
00063         AllrecHits[i] = recHits[i];
00064     }
00065 
00066     // Set the signal open
00067     isInputset = true;
00068 }
00069 
00070 void RPCCosmicSeedrecHitFinder::setEdge(const edm::EventSetup& iSetup) {
00071   
00072     // Get RPCGeometry
00073     edm::ESHandle<RPCGeometry> rpcGeometry;
00074     iSetup.get<MuonGeometryRecord>().get(rpcGeometry);
00075 
00076     // Find all chamber in RB1in and collect their surface
00077     const std::vector<DetId> AllRPCId = rpcGeometry->detIds();
00078     for(std::vector<DetId>::const_iterator it = AllRPCId.begin(); it != AllRPCId.end(); it++) {
00079         RPCDetId RPCId(it->rawId());
00080         int Region = RPCId.region();
00081         int Station = RPCId.station();
00082         int Layer = RPCId.layer();
00083         if(Region == 0 && Station == 1 && Layer == 1) {
00084             const BoundPlane RPCChamberSurface = rpcGeometry->chamber(RPCId)->surface();
00085             innerBounds.push_back(RPCChamberSurface);
00086         }
00087     }
00088 
00089     // Set the signal open
00090     isEdgeset = true;
00091 }
00092 
00093 void RPCCosmicSeedrecHitFinder::unsetEdge() {
00094 
00095     // Clear all surfaces of chambers in RB1in
00096     innerBounds.clear();
00097     isEdgeset = false;
00098 }
00099 
00100 void RPCCosmicSeedrecHitFinder::unsetInput() {
00101 
00102     for(unsigned int i = 0; i < RPCLayerNumber; i++)
00103         AllrecHits[i].clear();
00104     isInputset = false;
00105 }
00106 
00107 void RPCCosmicSeedrecHitFinder::setOutput(RPCSeedFinder *Seed) {
00108 
00109     theSeed = Seed;
00110     // Set the signal open
00111     isOutputset = true;
00112 }
00113 
00114 void RPCCosmicSeedrecHitFinder::setLayers(const std::vector<unsigned int>& Layers) {
00115 
00116     LayersinRPC = Layers;
00117     // Set the signal open
00118     isLayerset = true;
00119 }
00120 
00121 void RPCCosmicSeedrecHitFinder::fillrecHits() {
00122 
00123     if(isLayerset == false || isConfigured == false || isOutputset == false || isInputset == false || isEdgeset == false) {
00124         cout << "Not set the IO or not configured yet" << endl;
00125         return;
00126     } 
00127 
00128     therecHits.clear();
00129 
00130     if(LayersinRPC.size() == 0) {
00131         cout << "Not set with any layers" << endl;
00132         LayersinRPC.clear();
00133         therecHits.clear();
00134         isLayerset = false;
00135     }
00136 
00137     // check the layers, 1=all barrel, 2=all endcap, 3=mix
00138     unsigned int Component = LayerComponent();
00139     if(Component == 3)
00140         isLayersmixed = true;
00141     else
00142         isLayersmixed = false;
00143 
00144     GlobalVector initVector(0, 0, 0);
00145     const MuonRecHitPointer recHitRef;
00146     isOuterLayerfilled = false;
00147     complete(initVector, recHitRef);
00148 
00149     // Unset the signal
00150     LayersinRPC.clear();
00151     therecHits.clear();
00152     isLayerset = false;
00153 }
00154 
00155 int RPCCosmicSeedrecHitFinder::LayerComponent() {
00156 
00157     bool isBarrel = false;
00158     bool isEndcap = false;
00159     for(std::vector<unsigned int>::const_iterator it = LayersinRPC.begin(); it != LayersinRPC.end(); it++) {
00160         if((*it) < BarrelLayerNumber)
00161             isBarrel = true;
00162         if((*it) >= BarrelLayerNumber && (*it) < (BarrelLayerNumber+EachEndcapLayerNumber*2))
00163             isEndcap = true;
00164     }
00165     if(isBarrel == true && isEndcap == true)
00166         return 3;
00167     if(isEndcap == true)
00168         return 2;
00169     if(isBarrel == true)
00170         return 1;
00171     return 0;
00172 }
00173 
00174 bool RPCCosmicSeedrecHitFinder::complete(const GlobalVector& lastSegment, const MuonRecHitPointer& lastrecHitRef) {
00175 
00176     bool isrecHitsfound = false;
00177 
00178     for(unsigned int i = 0; i < RPCLayerNumber; i++)
00179         for(MuonRecHitContainer::const_iterator it = AllrecHits[i].begin(); it != AllrecHits[i].end(); it++) {
00180 
00181             cout << "Finding recHits from " << i << " th layer" << endl;
00182             // information for recHits
00183             GlobalPoint currentPosition = (*it)->globalPosition();
00184             int currentBX;
00185 
00186             // Check validation
00187             if(!(*it)->isValid())
00188                 continue;
00189 
00190             // Check BX range, be sure there is only RPCRecHit in the MuonRecHitContainer when use the dynamic_cast
00191             TrackingRecHit* thisTrackingRecHit = (*it)->hit()->clone();
00192             // Should also delete the RPCRecHit object cast by dynamic_cast<> ?
00193             RPCRecHit* thisRPCRecHit = dynamic_cast<RPCRecHit*>(thisTrackingRecHit);
00194             currentBX = thisRPCRecHit->BunchX();
00195             int ClusterSize = thisRPCRecHit->clusterSize();
00196             delete thisTrackingRecHit;
00197             // Check BX
00198             if((unsigned int)abs(currentBX) > BxRange)
00199                 continue;
00200 
00201             // Check cluster size
00202             bool Clustercheck = false;
00203             if(ClusterSet.size() == 0)
00204                 Clustercheck = true;
00205             for(std::vector<int>::const_iterator CluIter = ClusterSet.begin(); CluIter != ClusterSet.end(); CluIter++)
00206                 if(ClusterSize == (*CluIter))
00207                     Clustercheck = true;
00208             if(Clustercheck != true)
00209                 continue;
00210  
00211             cout << "Candidate recHit's position: " << currentPosition.x() << ", " << currentPosition.y() << ", " << currentPosition.z() << ". BX : " << currentBX << endl;
00212             // Fill 1st recHit from outer layers and rest recHits from all layers
00213             if(isOuterLayerfilled == false) {
00214                 // Pick out the recHit from outer layers and fill it
00215                 if(!isouterLayer(*it))
00216                     continue;
00217 
00218                 // If pass all, add to the seed
00219                 GlobalVector currentSegment = GlobalVector(0, 0, 0); 
00220                 cout << "1st recHit's global position: " << currentPosition.x() << ", " << currentPosition.y() << ", " << currentPosition.z() << ". BX: " << currentBX << endl;
00221                 isrecHitsfound = true;
00222                 therecHits.push_back(*it);
00223                 isOuterLayerfilled = true;
00224                 complete(currentSegment, *it);
00225                 // Remember to pop the recHit before add another one from the same layer!
00226                 therecHits.pop_back();
00227                 isOuterLayerfilled = false;
00228             }
00229             else {
00230                 GlobalPoint lastPosition = lastrecHitRef->globalPosition();
00231                 TrackingRecHit* lastTrackingRecHit = lastrecHitRef->hit()->clone();
00232                 // Should also delete the RPCRecHit object cast by dynamic_cast<> ?
00233                 RPCRecHit* lastRPCRecHit = dynamic_cast<RPCRecHit*>(lastTrackingRecHit);
00234                 int lastBX = lastRPCRecHit->BunchX();
00235                 delete lastTrackingRecHit;
00236 
00237                 // Check the Y coordinate, shoule be lower than current one
00238                 if(currentPosition.y() >= lastPosition.y())
00239                     continue;
00240 
00241                 // Check the BX, should be larger than current one
00242                 if(currentBX < lastBX)
00243                     continue;
00244 
00245                 // If be the 2nd recHit, just fill it
00246                 bool isinsideRegion = isinsideAngleRange(lastSegment, lastPosition, currentPosition);
00247                 cout << "Check isinsideRegion: " << isinsideRegion << endl;
00248                 if(!isinsideRegion)
00249                     continue;
00250 
00251                 // If cross the edge the recHit should belong to another seed
00252                 bool iscrossanyEdge = iscorssEdge(lastrecHitRef, *it);
00253                 cout << "Check iscrossanyEdge: " << iscrossanyEdge << endl;
00254                 if(iscrossanyEdge)
00255                     continue;
00256 
00257                 // If pass all, add to the seed
00258                 unsigned int NumberinSeed = therecHits.size();
00259                 GlobalVector currentSegment = (GlobalVector)(currentPosition - lastPosition);
00260                 cout << (NumberinSeed + 1) << "th recHit's global position: " << currentPosition.x() << ", " << currentPosition.y() << ", " << currentPosition.z() << ". BX: " << currentBX << endl;
00261                 isrecHitsfound = true;
00262                 therecHits.push_back(*it);
00263 
00264                 // if could not find next recHit in the search path, and have enough recHits already, that is the candidate
00265                 bool findNext = complete(currentSegment, *it);
00266                 if(findNext == false && therecHits.size() > 3) {
00267                     for(ConstMuonRecHitContainer::const_iterator iter = therecHits.begin(); iter != therecHits.end(); iter++)
00268                         cout << "Find recHit in seed candidate : " << (*iter)->globalPosition().x() << ", " << (*iter)->globalPosition().y() << ", " << (*iter)->globalPosition().z() << endl;
00269                     checkandfill();
00270                 }
00271 
00272                 // Remember to pop the recHit before add another one from the same layer!
00273                 therecHits.pop_back();
00274             }
00275         }
00276     return isrecHitsfound;
00277 }
00278 
00279 bool RPCCosmicSeedrecHitFinder::isouterLayer(const MuonRecHitPointer& recHitRef) {
00280 
00281     bool isinsideLayers = false;
00282     for(std::vector<unsigned int>::const_iterator it = LayersinRPC.begin(); it != LayersinRPC.end(); it++) {
00283         MuonRecHitContainer::const_iterator index = find(AllrecHits[*it].begin(), AllrecHits[*it].end(), recHitRef);
00284         if(index != AllrecHits[*it].end())
00285             isinsideLayers = true;
00286     }
00287     return isinsideLayers;
00288 }
00289 
00290 bool RPCCosmicSeedrecHitFinder::isinsideAngleRange(const GlobalVector& lastSegment, const GlobalPoint& lastPosition, const GlobalPoint& currentPosition) {
00291 
00292     bool isinsideAngle = true;
00293     GlobalVector SegVec = currentPosition - lastPosition;
00294     if(lastSegment.mag() != 0)
00295         if(fabs((lastSegment.phi()-SegVec.phi()).value()) > MaxDeltaPhi)
00296             isinsideAngle = false;
00297         
00298     return isinsideAngle;
00299 }
00300 
00301 bool RPCCosmicSeedrecHitFinder::iscorssEdge(const MuonRecHitPointer& lastrecHitRef, const MuonRecHitPointer& currentrecHitRef) {
00302 
00303     bool iscorss = false;
00304 
00305     
00306     // Check if 2 recHits corss the inner bounds
00307     GlobalPoint lastPosition = lastrecHitRef->globalPosition();
00308     GlobalPoint currentPosition = currentrecHitRef->globalPosition();
00309     GlobalPoint testPosition((lastPosition.x()+currentPosition.x())/2, (lastPosition.y()+currentPosition.y())/2, (lastPosition.z()+currentPosition.z())/2);
00310     /*
00311     for(std::vector<BoundPlane>::const_iterator it = innerBounds.begin(); it != innerBounds.end(); it++) {
00312         //SurfaceOrientation::Side TestSide0 = it->side(currentPosition, 0);
00313         //SurfaceOrientation::Side TestSide1 = it->side(lastPosition, 0);
00314         SurfaceOrientation::Side TestSide = it->side(testPosition, 0);
00315         //cout << "Side of currentPosition: " << TestSide0 << ", Side of lastPosition: " << TestSide1 << ", Side of middlePosition: " << TestSide << endl;
00316         //if(TestSide != SurfaceOrientation::positiveSide)
00317             //iscorss = true;
00318     }
00319     */
00320 
00321     // Check when mixLayer is not set
00322     if(isLayersmixed == false) {
00323         DetId lastId = lastrecHitRef->geographicalId();
00324         RPCDetId lastRPCId(lastId.rawId());
00325         int lastRegion = lastRPCId.region();
00326         DetId currentId = currentrecHitRef->geographicalId();
00327         RPCDetId currentRPCId(currentId.rawId());
00328         int currentRegion = currentRPCId.region();
00329         // Check if 2 recHits from different regions
00330         if((lastRegion - currentRegion) != 0)
00331             iscorss = true;
00332     }
00333 
00334     return iscorss;
00335 }
00336 
00337 void RPCCosmicSeedrecHitFinder::checkandfill() {
00338 
00339     if(therecHits.size() >= 3) {
00340         theSeed->setrecHits(therecHits); 
00341         theSeed->seed();
00342     }
00343     else
00344         cout << "Layer less than 3, could not fill a RPCSeedFinder" << endl;
00345 }