CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_6_1_2_SLHC4_patch1/src/RecoTracker/RoadMapRecord/src/Roads.cc

Go to the documentation of this file.
00001 //
00002 // Package:         RecoTracker/RoadMapRecord
00003 // Class:           Roads
00004 // 
00005 // Description:     The Roads object holds the RoadSeeds
00006 //                  and the RoadSets of all Roads through 
00007 //                  the detector. A RoadSeed consists
00008 //                  of the inner and outer SeedRing,
00009 //                  a RoadSet consists of all Rings in
00010 //                  in the Road.
00011 //
00012 // Original Author: Oliver Gutsche, gutsche@fnal.gov
00013 // Created:         Thu Jan 12 21:00:00 UTC 2006
00014 //
00015 // $Author: gutsche $
00016 // $Date: 2007/02/05 19:22:46 $
00017 // $Revision: 1.6 $
00018 //
00019 
00020 #include <iostream>
00021 #include <sstream>
00022 
00023 #include "RecoTracker/RoadMapRecord/interface/Roads.h"
00024 
00025 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00026 
00027 #include "DataFormats/DetId/interface/DetId.h"
00028 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
00029 #include "DataFormats/SiPixelDetId/interface/PixelSubdetector.h"
00030 
00031 #include "TrackingTools/RoadSearchHitAccess/interface/RoadSearchDetIdHelper.h"
00032 
00033 Roads::Roads() {
00034 
00035 }
00036 
00037 Roads::Roads(std::string ascii_filename, const Rings *rings) : rings_(rings) {
00038 
00039   readInFromAsciiFile(ascii_filename);
00040 
00041 }
00042 
00043 Roads::~Roads() { 
00044 
00045 }
00046 
00047 void Roads::readInFromAsciiFile(std::string ascii_filename) {
00048 
00049   // input file
00050   std::ifstream input(ascii_filename.c_str());
00051 
00052   // variable declaration
00053   unsigned int counter         = 0;
00054   std::istringstream stream;
00055   std::string line;
00056   unsigned int nroads          = 0;
00057   unsigned int nrings          = 0;
00058   unsigned int nlayers         = 0;
00059   unsigned int index           = 0;
00060 
00061   // read in number of roads
00062   std::getline(input,line);
00063   while (std::isspace(line[0]) || (line[0] == 35) ) {
00064     std::getline(input,line);
00065   }
00066   stream.str(line);
00067   stream.clear();
00068   stream >> nroads;
00069 
00070   for (unsigned int road = 0; 
00071        road < nroads;
00072        ++road ) {
00073     // read in number of inner seed rings
00074     std::getline(input,line);
00075     while (std::isspace(line[0]) || (line[0] == 35) ) {
00076       std::getline(input,line);
00077     }
00078     std::vector<const Ring*> innerSeedRings;
00079     stream.str(line);
00080     stream.clear();
00081     stream >> nrings;
00082     for ( unsigned int i = 0;
00083           i < nrings;
00084           ++i ) {
00085 
00086       // read in ring indices for inner seed rings
00087       std::getline(input,line);
00088       while (std::isspace(line[0]) || (line[0] == 35) ) {
00089         std::getline(input,line);
00090       }
00091       stream.str(line);
00092       stream.clear();
00093       stream >> index;
00094       innerSeedRings.push_back(rings_->getRing(index));
00095     }
00096 
00097     // read in number of outer seed rings
00098     std::getline(input,line);
00099     while (std::isspace(line[0]) || (line[0] == 35) ) {
00100       std::getline(input,line);
00101     }
00102     std::vector<const Ring*> outerSeedRings;
00103     stream.str(line);
00104     stream.clear();
00105     stream >> nrings;
00106     for ( unsigned int i = 0;
00107           i < nrings;
00108           ++i ) {
00109 
00110       // read in ring indices for outer seed rings
00111       std::getline(input,line);
00112       while (std::isspace(line[0]) || (line[0] == 35) ) {
00113         std::getline(input,line);
00114       }
00115       stream.str(line);
00116       stream.clear();
00117       stream >> index;
00118       outerSeedRings.push_back(rings_->getRing(index));
00119     }
00120 
00121     // RoadSeed
00122     RoadSeed seed(innerSeedRings,outerSeedRings);
00123 
00124     // RoadSet
00125     RoadSet set;
00126 
00127     // number of layers in road set
00128     std::getline(input,line);
00129     while (std::isspace(line[0]) || (line[0] == 35) ) {
00130       std::getline(input,line);
00131     }
00132     stream.str(line);
00133     stream.clear();
00134     stream >> nlayers;
00135 
00136     for ( unsigned int i = 0;
00137           i < nlayers;
00138           ++i ) {
00139 
00140       std::vector<const Ring*> layer;
00141 
00142       // number of rings in layer
00143       std::getline(input,line);
00144       while (std::isspace(line[0]) || (line[0] == 35) ) {
00145         std::getline(input,line);
00146       }
00147       stream.str(line);
00148       stream.clear();
00149       stream >> nrings;
00150       for ( unsigned int j = 0; j < nrings; ++j ) {
00151         std::getline(input,line);
00152         while (std::isspace(line[0]) || (line[0] == 35) ) {
00153           std::getline(input,line);
00154         }
00155         stream.str(line);
00156         stream.clear();
00157         stream >> index;
00158         layer.push_back(rings_->getRing(index));
00159       }
00160       set.push_back(layer);
00161     }
00162       
00163 
00164     // add seed and set to map
00165     roadMap_.insert(make_pair(seed,set));
00166     ++counter;
00167   }
00168   
00169   edm::LogInfo("RoadSearch") << "Read in: " << counter << " RoadSets from file: " << ascii_filename;
00170 
00171 }
00172 
00173 void Roads::dump(std::string ascii_filename) const {
00174 
00175   std::ofstream stream(ascii_filename.c_str());
00176   
00177   dumpHeader(stream);
00178 
00179   stream << "### Road information ###" << std::endl;
00180   stream << roadMap_.size() << std::endl;
00181 
00182   unsigned int counter = 0;
00183 
00184   for ( const_iterator roaditerator = roadMap_.begin(); 
00185         roaditerator != roadMap_.end(); 
00186         ++roaditerator ) {
00187 
00188     ++counter;
00189 
00190     stream << "### RoadMap Entry " << counter << " ###" << std::endl;
00191 
00192     RoadSeed seed = (*roaditerator).first;
00193     RoadSet  set  = (*roaditerator).second;
00194 
00195     stream << "### RoadSeed First Ring ###" << std::endl;
00196     stream << seed.first.size() << std::endl;
00197     for (std::vector<const Ring*>::const_iterator ring = seed.first.begin();
00198          ring != seed.first.end();
00199          ++ring ) {
00200       stream << (*ring)->getindex() << std::endl;
00201     }
00202     stream << "### RoadSeed Second Ring ###" << std::endl;
00203     stream << seed.second.size() << std::endl;
00204     for (std::vector<const Ring*>::const_iterator ring = seed.second.begin();
00205          ring != seed.second.end();
00206          ++ring ) {
00207       stream << (*ring)->getindex() << std::endl;
00208     }
00209 
00210     stream << "### RoadSet ###" << std::endl;
00211     stream << set.size() << std::endl;
00212     for ( RoadSet::const_iterator layer = set.begin(); layer != set.end(); ++layer ) {
00213       stream << "### Layer ###" << std::endl;
00214       stream << layer->size() << std::endl;
00215       for ( std::vector<const Ring*>::const_iterator ring = layer->begin();
00216             ring != layer->end();
00217             ++ring ) {
00218         stream << (*ring)->getindex() << std::endl;
00219       }
00220     }
00221   }
00222 }
00223 
00224 void Roads::dumpHeader(std::ofstream &stream) const {
00225 
00226   stream << "#" << std::endl;
00227   stream << "# Roads for the RoadSearch tracking algorithm" << std::endl;
00228   stream << "# Ascii Dump" << std::endl;
00229   stream << "# " << std::endl;
00230   stream << "# Content:" << std::endl;
00231   stream << "# " << std::endl;
00232   stream << "# a dump of the RoadMap structure:" << std::endl;
00233   stream << "#" << std::endl;
00234   stream << "# Road Information: <number of roads>" << std::endl;
00235   stream << "# Ring: index, rmin, rmax, zmin, zmax, std::vector<DetId>: Ring of DetUnits in phi taken from ring service" << std::endl;
00236   stream << "# RoadSeed: std::pair<std::vector<const Ring*>,std::vector<const Ring*> >: inner and outer Ring Seed for the Road" << std::endl;
00237   stream << "# RoadSet : std::vector<std::vectro<const Ring*> >: all Rings belonging to a road structured in layers" << std::endl;
00238   stream << "# RoadMap: std::multimap<RoadSeed,RoadSet>: main container for the Roads" << std::endl;
00239   stream << "# " << std::endl;
00240   stream << "# Ascii-Format:" << std::endl;
00241   stream << "# " << std::endl;
00242   stream << "# Road Information:" << std::endl;
00243   stream << "#       <number of roads>" << std::endl;
00244   stream << "#" << std::endl;
00245   stream << "# RoadMap for each road:" << std::endl;
00246   stream << "#" << std::endl;
00247   stream << "#       ### RoadMap Entry ###" << std::endl;
00248   stream << "#       ### RoadSeed First Ring ###" << std::endl;
00249   stream << "#       <number of inner seed rings>" << std::endl;
00250   stream << "#       <index>" << std::endl;
00251   stream << "#       <index>" << std::endl;
00252   stream << "#       ..." << std::endl;
00253   stream << "#       ### RoadSeed Second Ring ###" << std::endl;
00254   stream << "#       <number of outer seed rings>" << std::endl;
00255   stream << "#       <index>" << std::endl;
00256   stream << "#       <index>" << std::endl;
00257   stream << "#       ..." << std::endl;
00258   stream << "#       ### RoadSet ###" << std::endl;
00259   stream << "#       <number of Layers in RoadSet>" << std::endl;
00260   stream << "#       ### Layer ###" << std::endl;
00261   stream << "#       <number of rings in layer>" << std::endl;
00262   stream << "#       <index>" << std::endl;
00263   stream << "#       <index>" << std::endl;
00264   stream << "#        ..." << std::endl;
00265   stream << "#       ### Layer ###" << std::endl;
00266   stream << "#        ..." << std::endl;
00267   stream << "#" << std::endl;
00268   stream << "#" << std::endl;
00269   
00270 }
00271 
00272 const Roads::RoadSeed* Roads::getRoadSeed(DetId InnerSeedRing, 
00273                                           DetId OuterSeedRing, 
00274                                           double InnerSeedRingPhi,
00275                                           double OuterSeedRingPhi,
00276                                           double dphi_scalefactor) const {
00277 
00278   // loop over seed Ring pairs
00279 
00280   // determine ringtype for inner seed ring detid
00281   Ring::type innerSeedRingType = getRingType(InnerSeedRing);
00282   Ring::type outerSeedRingType = getRingType(OuterSeedRing);
00283 
00284   for ( const_iterator road = roadMap_.begin(); road != roadMap_.end(); ++road ) {
00285     for ( std::vector<const Ring*>::const_iterator innerRing = road->first.first.begin();
00286           innerRing != road->first.first.end();
00287           ++innerRing ) {
00288       if ( (*innerRing)->getType() == innerSeedRingType ) {
00289         for ( std::vector<const Ring*>::const_iterator outerRing = road->first.second.begin();
00290               outerRing != road->first.second.end();
00291               ++outerRing ) {
00292           if ( (*outerRing)->getType() == outerSeedRingType ) {
00293             if ( (*innerRing)->containsDetId(InnerSeedRing,InnerSeedRingPhi,dphi_scalefactor) &&
00294                  (*outerRing)->containsDetId(OuterSeedRing,OuterSeedRingPhi,dphi_scalefactor) ) {
00295               return &(road->first);
00296             }
00297           }
00298         }
00299       }
00300     }
00301   }
00302       
00303   edm::LogError("RoadSearch") << "RoadSeed could not be found for inner SeedRing type: " << innerSeedRingType << " DetId: " << InnerSeedRing.rawId() 
00304                               << " at " << InnerSeedRingPhi
00305                               << " and outer SeedRing type : " << outerSeedRingType << " DetID: " << OuterSeedRing.rawId() 
00306                               << " at " << OuterSeedRingPhi;
00307   return 0;
00308 }
00309 
00310 const Roads::RoadSeed* Roads::getRoadSeed(std::vector<DetId> seedRingDetIds,
00311                                           std::vector<double> seedRingHitsPhi,
00312                                           double dphi_scalefactor) const {
00313   //
00314   // loop over roads and return first road which contains all seedRingDetIds
00315   //
00316 
00317   for ( const_iterator road = roadMap_.begin(); road != roadMap_.end(); ++road ) {
00318     unsigned int found = 0;
00319     for ( unsigned int detIdCounter = 0;
00320           detIdCounter < seedRingDetIds.size();
00321           ++detIdCounter ) {
00322       DetId      id   = RoadSearchDetIdHelper::ReturnRPhiId(seedRingDetIds[detIdCounter]);
00323       double     phi  = seedRingHitsPhi[detIdCounter];
00324       Ring::type type = getRingType(id);
00325 
00326       bool foundInInnerRing = false;
00327       for ( std::vector<const Ring*>::const_iterator innerRing = road->first.first.begin();
00328             innerRing != road->first.first.end();
00329             ++innerRing ) {
00330         if ( (*innerRing)->getType() == type ) {
00331           if ( (*innerRing)->containsDetId(id,phi,dphi_scalefactor) ) {
00332             ++found;
00333             foundInInnerRing = true;
00334           }
00335         }
00336       }
00337 
00338       if ( !foundInInnerRing ) {
00339         for ( std::vector<const Ring*>::const_iterator outerRing = road->first.second.begin();
00340               outerRing != road->first.second.end();
00341               ++outerRing ) {
00342           if ( (*outerRing)->getType() == type ) {
00343             if ( (*outerRing)->containsDetId(id,phi,dphi_scalefactor) ) {
00344               ++found;
00345             }
00346           }
00347         }
00348       }
00349 
00350       if ( found == seedRingDetIds.size() ) {
00351               return &(road->first);
00352       }
00353     }
00354   }
00355 
00356   std::ostringstream ost;
00357   
00358   ost << "RoadSeed could not be found for following hits:\n";
00359   for ( unsigned int detIdCounter = 0;
00360         detIdCounter < seedRingDetIds.size();
00361         ++detIdCounter ) {
00362     ost << "Hit DetId: " << seedRingDetIds[detIdCounter].rawId() << " phi: " << seedRingHitsPhi[detIdCounter] << "\n";
00363   }
00364   
00365   edm::LogError("RoadSearch") << ost.str();
00366   
00367   return 0;
00368 }
00369 
00370 const Roads::type Roads::getRoadType(const RoadSeed *const seed) const {
00371   //
00372   // check if one of the outer rings is in TOB, then mark as RPhi
00373   // problematic for transition region
00374   bool TOBRing = false;
00375   for ( std::vector<const Ring*>::const_iterator ring = seed->second.begin();
00376         ring != seed->second.end();
00377         ++ring) {
00378     if ( (*ring)->getType() == Ring::TOBRing) {
00379       TOBRing = true;
00380     }
00381   }
00382   if ( TOBRing ) {
00383     return Roads::RPhi;
00384   } else {
00385     return Roads::ZPhi;
00386   }
00387 }
00388 
00389 
00390 const Ring::type Roads::getRingType(DetId id) const {
00391 
00392   Ring::type type = Ring::Unspecified;
00393 
00394   if ( (unsigned int)id.subdetId() == StripSubdetector::TIB ) {
00395     type = Ring::TIBRing;
00396   } else if ( (unsigned int)id.subdetId() == StripSubdetector::TOB ) {
00397     type = Ring::TOBRing;
00398   } else if ( (unsigned int)id.subdetId() == StripSubdetector::TID ) {
00399     type = Ring::TIDRing;
00400   } else if ( (unsigned int)id.subdetId() == StripSubdetector::TEC ) {
00401     type = Ring::TECRing;
00402   } else if ( (unsigned int)id.subdetId() == PixelSubdetector::PixelBarrel ) {
00403     type = Ring::PXBRing;
00404   } else if ( (unsigned int)id.subdetId() == PixelSubdetector::PixelEndcap ) {
00405     type = Ring::PXFRing;
00406   }
00407 
00408   return type;
00409 
00410 }