CMS 3D CMS Logo

/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_5_3_13_patch3/src/RecoTracker/MeasurementDet/plugins/OnDemandMeasurementTracker.cc

Go to the documentation of this file.
00001 #include "OnDemandMeasurementTracker.h"
00002 
00003 #include "FWCore/Framework/interface/ESHandle.h"
00004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00005 
00006 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
00007 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
00008 #include "Geometry/CommonDetUnit/interface/GeomDetUnit.h"
00009 #include "Geometry/TrackerGeometryBuilder/interface/GluedGeomDet.h"
00010 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
00011 
00012 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
00013 
00014 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
00015 #include "DataFormats/SiPixelDetId/interface/PixelSubdetector.h"
00016 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
00017 #include "DataFormats/SiStripCluster/interface/SiStripClusterCollection.h"
00018 #include "DataFormats/Common/interface/ContainerMask.h"
00019 
00020 #include "TrackingTools/MeasurementDet/interface/MeasurementDetException.h"
00021 
00022 #include "RecoLocalTracker/ClusterParameterEstimator/interface/PixelClusterParameterEstimator.h"
00023 #include "RecoLocalTracker/Records/interface/TrackerCPERecord.h"
00024 #include "RecoLocalTracker/SiStripRecHitConverter/interface/SiStripRecHitMatcher.h"
00025 #include "RecoLocalTracker/SiStripRecHitConverter/interface/StripCPE.h"  
00026 
00027 #include "RecoTracker/TkDetLayers/interface/GeometricSearchTracker.h"
00028 #include "RecoTracker/Record/interface/TrackerRecoGeometryRecord.h"
00029 #include "TkStripMeasurementDet.h"
00030 #include "TkPixelMeasurementDet.h"
00031 #include "TkGluedMeasurementDet.h"
00032 
00033 #include "CondFormats/SiStripObjects/interface/SiStripNoises.h"
00034 #include "CondFormats/DataRecord/interface/SiStripNoisesRcd.h"
00035 
00036 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
00037 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
00038 #include "CalibFormats/SiStripObjects/interface/SiStripDetCabling.h"
00039 #include "CalibTracker/Records/interface/SiStripDetCablingRcd.h"
00040 #include "CalibTracker/Records/interface/SiStripRegionCablingRcd.h"
00041 
00042 #include <iostream>
00043 #include <typeinfo>
00044 #include <map>
00045 
00046 #include <DataFormats/GeometrySurface/interface/BoundPlane.h>
00047 #include "DataFormats/Math/interface/deltaR.h"
00048 
00049 #include "FWCore/ServiceRegistry/interface/Service.h"
00050 #include "FWCore/Services/interface/UpdaterService.h"
00051 
00052 using namespace std;
00053 
00054 OnDemandMeasurementTracker::OnDemandMeasurementTracker(const edm::ParameterSet&              conf,
00055                                                        const PixelClusterParameterEstimator* pixelCPE,
00056                                                        const StripClusterParameterEstimator* stripCPE,
00057                                                        const SiStripRecHitMatcher*  hitMatcher,
00058                                                        const TrackerGeometry*  trackerGeom,
00059                                                        const GeometricSearchTracker* geometricSearchTracker,
00060                                                        const SiStripQuality *stripQuality,
00061                                                        int   stripQualityFlags,
00062                                                        int   stripQualityDebugFlags,
00063                                                        const SiPixelQuality *pixelQuality,
00064                                                        const SiPixelFedCabling *pixelCabling,
00065                                                        int   pixelQualityFlags,
00066                                                        int   pixelQualityDebugFlags,
00067                                                        const SiStripRegionCabling * stripRegionCabling,
00068                                                        bool isRegional):
00069   MeasurementTrackerImpl(conf,pixelCPE,stripCPE,hitMatcher,trackerGeom,geometricSearchTracker,
00070         stripQuality,stripQualityFlags,stripQualityDebugFlags,
00071         pixelQuality,pixelCabling,pixelQualityFlags,pixelQualityDebugFlags,
00072         isRegional)
00073   , category_("OnDemandMeasurementTracker")
00074   , StayPacked_(true)
00075   , StripOnDemand_(true)
00076   , PixelOnDemand_(false)
00077   , theStripRegionCabling(stripRegionCabling)
00078 {
00079   //  the constructor does construct the regular MeasurementTracker
00080   //  then a smart copy of the DetMap is made into DetODMap: this could be avoided with modification to MeasurementDet interface
00081   //  the elementIndex to be defined in the refgetter is mapped to the detId
00082   //  flags are set to initialize the DetODMap
00083   
00084   std::map<SiStripRegionCabling::ElementIndex, std::vector< DetODContainer::const_iterator> > local_mapping;
00085 
00086   for (DetContainer::iterator it=theDetMap.begin(); it!= theDetMap.end();++it)
00087     {
00088       DetODContainer::iterator inserted = theDetODMap.insert(make_pair(it->first,DetODStatus(const_cast<MeasurementDet*>(it->second)))).first;
00089 
00090 
00091       GeomDet::SubDetector subdet = it->second->geomDet().subDetector();
00092       if (subdet == GeomDetEnumerators::PixelBarrel  || subdet == GeomDetEnumerators::PixelEndcap ){
00093         //special flag treatement for pixels
00094         //one can never be updated and not defined. except if we don't need to care about them: pixels
00095         inserted->second.defined=false;
00096         inserted->second.updated=true;
00097       }//pixel module
00098       else if (subdet == GeomDetEnumerators::TIB || subdet == GeomDetEnumerators::TOB ||
00099           subdet == GeomDetEnumerators::TID || subdet == GeomDetEnumerators::TEC )
00100         {
00101           //set flag to false
00102           inserted->second.defined=false;
00103           inserted->second.updated=false;
00104 
00105           //what will be the element index in the refgetter
00106           GlobalPoint center = it->second->geomDet().position();
00107           double eta = center.eta();
00108           double phi = center.phi();
00109           uint32_t id = it->first;
00110           SiStripRegionCabling::ElementIndex eIndex = theStripRegionCabling->elementIndex(SiStripRegionCabling::Position(eta,phi),
00111                                                                                           SiStripRegionCabling::subdetFromDetId(id),
00112                                                                                           SiStripRegionCabling::layerFromDetId(id));
00113           LogDebug(category_)<<"region selected (from "<<id<<" center) is:\n"
00114                              <<"position: "<<center
00115                              <<"\n center absolute index: "<<theStripRegionCabling->region(theStripRegionCabling->positionIndex(SiStripRegionCabling::Position(eta,phi)))
00116                              <<"\n center position index: "<<theStripRegionCabling->positionIndex(SiStripRegionCabling::Position(eta,phi)).first<<
00117             " "<<theStripRegionCabling->positionIndex(SiStripRegionCabling::Position(eta,phi)).second
00118                              <<"\n center postion: "<<theStripRegionCabling->position(theStripRegionCabling->positionIndex(SiStripRegionCabling::Position(eta,phi))).first<<
00119             " "<<theStripRegionCabling->position(theStripRegionCabling->positionIndex(SiStripRegionCabling::Position(eta,phi))).second
00120                              <<"\n eta: "<<eta
00121                              <<"\n phi: "<<phi
00122                              <<"\n subedet: "<<SiStripRegionCabling::subdetFromDetId(id)
00123                              <<" layer: "<<SiStripRegionCabling::layerFromDetId(id);
00124 
00125           //      register those in a map
00126           //to be able to know what are the detid in a given elementIndex
00127           local_mapping[eIndex].push_back(inserted);
00128         }//strip module
00129       else{
00130         //abort
00131         edm::LogError(category_)<<"not a tracker geomdet in constructor: "<<it->first;
00132         throw MeasurementDetException("OnDemandMeasurementTracker dealing with a non tracker GeomDet.");
00133       }//abort
00134     }//loop over DetMap
00135   if (theInactiveStripDetectorLabels.size()!=0)
00136     theRawInactiveStripDetIds.reserve(200);
00137 
00138   //move into a vector
00139   region_mapping.reserve(local_mapping.size());
00140   for( auto eIt= local_mapping.begin();
00141        eIt!=local_mapping.end();++eIt)
00142     region_mapping.push_back(std::make_pair((*eIt).first,(*eIt).second));
00143 }
00144 
00145 
00146 void OnDemandMeasurementTracker::define( const edm::Handle< LazyGetter> & aLazyGetterH,
00147                                          std::auto_ptr< RefGetter > & aGetter ) const
00148 {
00149   //  define is supposed to be call by an EDProducer module, which wil put the RefGetter in the event
00150   //  so that reference can be made to it.
00151   //  the lazy getter is retrieved by the calling module and passed along with the event
00152   //  the map is cleared, except for pixel
00153   //  then the known elementIndex are defined to the RefGetter. no unpacking is done at this time
00154   //  the defined region range is registered in the DetODMap for further use.
00155 
00156   //clear all defined tags to start from scratch (except for pixel)
00157   for (DetODContainer::iterator it=theDetODMap.begin(); it!= theDetODMap.end();++it)
00158     {
00159       if (it->second.updated && !it->second.defined) continue; //special treatement for pixels
00160       it->second.defined= false; 
00161       it->second.updated = false;
00162     }
00163 
00164   // nedeed??
00165   theStDets.setLazyGetter(aLazyGetterH);
00166 
00167 
00168   //define all the elementindex in the refgetter
00169   for(auto eIt= region_mapping.begin();
00170        eIt!=region_mapping.end();++eIt){
00171     std::pair<unsigned int, unsigned int> region_range; 
00172     
00173     //before update of the refgetter
00174     region_range.first = aGetter->size();
00175     //update the refegetter with the elementindex
00176     theStripRegionCabling->updateSiStripRefGetter<SiStripCluster> (*aGetter, aLazyGetterH, eIt->first);
00177     //after update of the refgetter
00178     region_range.second = aGetter->size();
00179 
00180     LogDebug(category_)<<"between index: "<<region_range.first<<" "<<region_range.second
00181                        <<"\n"<<dumpRegion(region_range,*aGetter,StayPacked_);
00182     
00183     //now assign to each measurement det for that element index
00184     for (auto dIt=eIt->second.begin();
00185          dIt!=eIt->second.end();++dIt){
00186       DetODStatus & elem = const_cast<DetODStatus &>((*dIt)->second);
00187       elem.region_range = region_range;
00188       elem.defined=true;
00189       LogDebug(category_)<<"detId: "<<(*dIt)->first<<" in region range: "<<region_range.first<<" "<<region_range.second;
00190     }//loop over MeasurementDet attached to that elementIndex
00191   }//loop over know elementindex
00192 }
00193 
00194 void OnDemandMeasurementTracker::updateStrips( const edm::Event& event) const 
00195 {
00196   bool oncePerEvent= edm::Service<UpdaterService>()->checkOnce("OnDemandMeasurementTracker::updateStrips::"+name_);
00197   bool failedToGet = false;
00198   if (!oncePerEvent)
00199     failedToGet = theRefGetterH.failedToGet() || theLazyGetterH.failedToGet();
00200 
00201   if (oncePerEvent || failedToGet)
00202     {
00203       LogDebug(category_)<<"Updating siStrip on event: "<< (unsigned int) event.id().run() <<" : "<<(unsigned int) event.id().event();
00204       
00205       //get the ref getter back from the event
00206       std::string stripClusterProducer = pset_.getParameter<std::string>("stripClusterProducer");
00207       event.getByLabel(stripClusterProducer,theRefGetterH);
00208       
00209       std::string stripLazyGetter = pset_.getParameter<std::string>("stripLazyGetterProducer");
00210       event.getByLabel(stripLazyGetter,theLazyGetterH);
00211 
00212       theStDets.setLazyGetter(theLazyGetterH);
00213 
00214 
00215       //get the skip clusters
00216       if (selfUpdateSkipClusters_){
00217         theSkipClusterRefs=true;
00218         event.getByLabel(pset_.getParameter<edm::InputTag>("skipClusters"),theStripClusterMask);
00219         theStripClusterMask->copyMaskTo(theStDets.clusterToSkip());
00220       } else {
00221         theStDets.clusterToSkip().clear();
00222       }
00223 
00224       //get the detid that are inactive
00225       theRawInactiveStripDetIds.clear();
00226       getInactiveStrips(event,theRawInactiveStripDetIds);
00227     }
00228 }
00229 
00230 void OnDemandMeasurementTracker::update( const edm::Event& event) const
00231 {
00232   //  update is supposed to be called by any module that is useing the MeasurementTracker
00233   //  after checking the the event has not yet been seen
00234   //  update the pixel using MeasurementTracekr specific function
00235   //  retreive the RefGetter from the event: the very one that has been pass to define(...) and put into the event
00236 
00237   if (!PixelOnDemand_) {
00238     LogDebug(category_)<<"pixel are not OnDemand. updating them a la MeasurmentTracker.";
00239     MeasurementTrackerImpl::updatePixels(event);}
00240   else{
00241     edm::LogError(category_)<<"trying to update siPixel as on-demand. Not Implemented yet.";
00242   }
00243 
00244   if (!StripOnDemand_) {
00245     LogDebug(category_)<<"strip are not OnDemand. updating them a la MeasurmentTracker.";
00246     MeasurementTrackerImpl::updateStrips(event);}
00247   else{
00248     LogDebug(category_)<<"strip are OnDemand. updating them a la OnDemandMeasurmentTracker."; 
00249     updateStrips(event);
00250   }
00251 }
00252 
00253 #include <sstream>
00254 
00255 std::string OnDemandMeasurementTracker::dumpCluster(const std::vector<SiStripCluster> ::const_iterator & begin,const  std::vector<SiStripCluster> ::const_iterator & end)const
00256 {
00257   //  dumpCluster is a printout of all the clusters between the iterator. returns a string
00258   std::string tab="      ";
00259   std::stringstream ss;
00260   std::vector<SiStripCluster> ::const_iterator it = begin;
00261   unsigned int i=0;
00262   for (;it!=end;++it){
00263     ss<<tab<<i++<<") center: "<<it->barycenter()<<",id: "<<it->geographicalId()<<" with: "<<it->amplitudes().size()<<" strips\n"<<tab<<tab<<"{";
00264     for (unsigned int is=0;is!=it->amplitudes().size();++is){
00265       ss<<it->amplitudes()[is]<<" ";
00266     }ss<<"}\n";
00267   }
00268   return ss.str();
00269 }
00270 
00271 std::string OnDemandMeasurementTracker::dumpRegion(std::pair<unsigned int,unsigned int> indexes,
00272                                              const RefGetter & theGetter,
00273                                              bool stayPacked)const
00274 {
00275   //  dumpRegion is a printout of all the clusters in a region defined on the RefGetter. returns a string
00276   std::stringstream ss;
00277   ss<<"cluster between: "<<indexes.first<<" and: "<<indexes.second<<"\n";
00278   for (unsigned int iRegion = indexes.first; iRegion != indexes.second; ++iRegion){    
00279     uint32_t reg = SiStripRegionCabling::region((theGetter)[iRegion].region());
00280     SiStripRegionCabling::Position pos = theStripRegionCabling->position(reg);
00281     SiStripRegionCabling::PositionIndex posI = theStripRegionCabling->positionIndex(reg);
00282     
00283     ss<<"Clusters for region:["<<iRegion<<"]"
00284       <<"\n element index: "<<(theGetter)[iRegion].region()
00285       <<"\n region absolute index: "<<reg
00286       <<"\n region position index: "<<posI.first<<" "<<posI.second
00287       <<"\n region position: "<<pos.first<<" "<<pos.second
00288       <<"\n"<< (stayPacked? " hidden to avoid unpacking." : dumpCluster((theGetter)[iRegion].begin(),(theGetter)[iRegion].end()));
00289   }
00290   return ss.str();
00291 }
00292 
00293 void OnDemandMeasurementTracker::assign(const TkStripMeasurementDet * csmdet,
00294                                   DetODContainer::iterator * alreadyFound) const {
00295   //  assign is using the handle to the refgetter and the region index range to update the MeasurementDet with their clusters
00296   
00297   TkStripMeasurementDet * smdet = const_cast<TkStripMeasurementDet *>(csmdet);
00298   DetId id = smdet->rawId();
00299   
00300   LogDebug(category_)<<"assigning: "<<id.rawId();
00301 
00302   // what is the iterator. do not look again if already found and provided with
00303   DetODContainer::iterator elementInMap;
00304   if (alreadyFound){ elementInMap=*alreadyFound;}
00305   else{ elementInMap = theDetODMap.find(id);}
00306   
00307   if  ( elementInMap != theDetODMap.end()){
00308     //flag it as updated
00309     elementInMap->second.updated = true;
00310 
00311     if (!theRawInactiveStripDetIds.empty() && std::binary_search(theRawInactiveStripDetIds.begin(), theRawInactiveStripDetIds.end(), id)) {
00312       smdet->setActiveThisEvent(false); 
00313       return;
00314     }
00315 
00316     //retrieve the region range index for this module
00317     std::pair<unsigned int,unsigned int> & indexes =elementInMap->second.region_range;
00318 
00319     //this printout will trigger unpacking. no problem. it is done on the next regular line (find(id.rawId())
00320     LogDebug(category_)<<"between index: "<<indexes.first<<" and: "<<indexes.second
00321                        <<"\nretrieved for module: "<<id.rawId()
00322                        <<"\n"<<dumpRegion(indexes,*theRefGetterH);
00323     
00324     //look for iterator range in the regions defined for that module
00325     for (unsigned int iRegion = indexes.first; iRegion != indexes.second; ++iRegion){
00326       RefGetter::record_pair range = (*theRefGetterH)[iRegion].find(id.rawId());
00327       if (range.first!=range.second){
00328         //      found something not empty
00329         //update the measurementDet
00330         smdet->update(range.first, range.second);
00331         LogDebug(category_)<<"Valid clusters for: "<<id.rawId()
00332                            <<"\nnumber of regions defined here: "<< indexes.second-indexes.first
00333                            <<"\n"<<dumpCluster(range.first,range.second);
00334         /* since theStripsToSkip is a "static" pointer of the MT, no need to set it at all time.
00335           if (selfUpdateSkipClusters_){
00336           //assign skip clusters
00337           smdet->setClusterToSkip(&theStripsToSkip);
00338         }
00339         */
00340         //and you are done
00341         return;}
00342     }//loop over regions, between indexes
00343 
00344     //if reached. no cluster are found. set the TkStripMeasurementDet to be empty
00345     smdet->setEmpty();
00346 
00347   }//found in the map
00348   else{
00349     //throw excpetion
00350     edm::LogError(category_)<<"failed to find the MeasurementDet for: "<<id.rawId();
00351     throw MeasurementDetException("failed to find the MeasurementDet for: <see message logger>");
00352   }
00353 }
00354 
00355 
00356 
00357 const MeasurementDet* 
00358 OnDemandMeasurementTracker::idToDet(const DetId& id) const
00359 {
00360   //  overloaded from MeasurementTracker
00361   //  find the detid. if not found throw exception
00362   //  if already updated: always for pixel or strip already queried. return it
00363   DetODContainer::iterator it = theDetODMap.find(id);
00364   if ( it != theDetODMap.end()) {
00365     
00366     //it has already been queried. and so already updated: nothing to be done here (valid for pixels too)
00367     if (it->second.updated){
00368       LogDebug(category_)<<"found id: "<<id.rawId()<<" as aleardy updated."; 
00369       return it->second.mdet;
00370     }
00371     
00372     if (StripOnDemand_){
00373       //check glued or single
00374       if (it->second.glued){
00375         //glued det
00376         LogDebug(category_)<<"updating glued id: "<<id.rawId();
00377         //cast to specific type
00378         TkGluedMeasurementDet*  theConcreteDet = static_cast<TkGluedMeasurementDet*>(it->second.mdet);
00379                 
00380         //      update the two components
00381         //update the mono
00382         assign(theConcreteDet->monoDet());
00383         //update the stereo
00384         assign(theConcreteDet->stereoDet());
00385               
00386         //flag the glued det as updated (components are flagged in assign)
00387         it->second.updated=true;
00388       }//glued det
00389       else{
00390         //single layer 
00391         LogDebug(category_)<<"updating singel id: "<<id.rawId();
00392         //cast to specific type
00393         TkStripMeasurementDet*  theConcreteDet = static_cast<TkStripMeasurementDet*>(it->second.mdet);
00394 
00395         //update the single layer
00396         assign(theConcreteDet,&it);
00397       }//single det
00398     }
00399     //eventually return it
00400     return it->second.mdet;
00401   }//found in the DetODMap
00402   else{
00403     //throw excpetion
00404     edm::LogError(category_)<<"failed to find the MeasurementDet for: "<<id.rawId();
00405     throw MeasurementDetException("failed to find the MeasurementDet for: <see message logger>");
00406   }
00407   return 0;
00408 }
00409 
00410