OnDemandMeasurementTracker.cc

Go to the documentation of this file.

#include "RecoTracker/MeasurementDet/interface/OnDemandMeasurementTracker.h"

#include "FWCore/Framework/interface/ESHandle.h"
#include "FWCore/MessageLogger/interface/MessageLogger.h"

#include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
#include "Geometry/CommonDetUnit/interface/GeomDet.h"
#include "Geometry/CommonDetUnit/interface/GeomDetUnit.h"
#include "Geometry/TrackerGeometryBuilder/interface/GluedGeomDet.h"
#include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"

#include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"

#include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
#include "DataFormats/SiPixelDetId/interface/PixelSubdetector.h"
#include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
#include "DataFormats/SiStripCluster/interface/SiStripClusterCollection.h"

#include "TrackingTools/MeasurementDet/interface/MeasurementDetException.h"

#include "RecoLocalTracker/ClusterParameterEstimator/interface/PixelClusterParameterEstimator.h"
#include "RecoLocalTracker/Records/interface/TrackerCPERecord.h"
#include "RecoLocalTracker/SiStripRecHitConverter/interface/SiStripRecHitMatcher.h"
#include "RecoLocalTracker/SiStripRecHitConverter/interface/StripCPE.h"  

#include "RecoTracker/TkDetLayers/interface/GeometricSearchTracker.h"
#include "RecoTracker/Record/interface/TrackerRecoGeometryRecord.h"
#include "RecoTracker/MeasurementDet/interface/TkStripMeasurementDet.h"
#include "RecoTracker/MeasurementDet/interface/TkPixelMeasurementDet.h"
#include "RecoTracker/MeasurementDet/interface/TkGluedMeasurementDet.h"

#include "CondFormats/SiStripObjects/interface/SiStripNoises.h"
#include "CondFormats/DataRecord/interface/SiStripNoisesRcd.h"

#include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
#include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
#include "CalibFormats/SiStripObjects/interface/SiStripDetCabling.h"
#include "CalibTracker/Records/interface/SiStripDetCablingRcd.h"
#include "CalibTracker/Records/interface/SiStripRegionCablingRcd.h"

#include <iostream>
#include <typeinfo>
#include <map>

#include <DataFormats/GeometrySurface/interface/BoundPlane.h>
#include "DataFormats/Math/interface/deltaR.h"

#include "FWCore/ServiceRegistry/interface/Service.h"
#include "FWCore/Services/interface/UpdaterService.h"

using namespace std;

OnDemandMeasurementTracker::OnDemandMeasurementTracker(const edm::ParameterSet&              conf,
                               const PixelClusterParameterEstimator* pixelCPE,
                               const StripClusterParameterEstimator* stripCPE,
                               const SiStripRecHitMatcher*  hitMatcher,
                               const TrackerGeometry*  trackerGeom,
                               const GeometricSearchTracker* geometricSearchTracker,
                               const SiStripQuality *stripQuality,
                                                       int   stripQualityFlags,
                                                       int   stripQualityDebugFlags,
                                                       const SiPixelQuality *pixelQuality,
                                                       const SiPixelFedCabling *pixelCabling,
                                                       int   pixelQualityFlags,
                                                       int   pixelQualityDebugFlags,
                               const SiStripRegionCabling * stripRegionCabling,
                               bool isRegional):
  MeasurementTracker(conf,pixelCPE,stripCPE,hitMatcher,trackerGeom,geometricSearchTracker,
        stripQuality,stripQualityFlags,stripQualityDebugFlags,
        pixelQuality,pixelCabling,pixelQualityFlags,pixelQualityDebugFlags,
        isRegional)
  , category_("OnDemandMeasurementTracker")
  , StayPacked_(true)
  , StripOnDemand_(true)
  , PixelOnDemand_(false)
  , theStripRegionCabling(stripRegionCabling)
{
  //  the constructor does construct the regular MeasurementTracker
  //  then a smart copy of the DetMap is made into DetODMap: this could be avoided with modification to MeasurementDet interface
  //  the elementIndex to be defined in the refgetter is mapped to the detId
  //  flags are set to initialize the DetODMap

  for (DetContainer::iterator it=theDetMap.begin(); it!= theDetMap.end();++it)
    {
      DetODContainer::iterator inserted = theDetODMap.insert(make_pair(it->first,DetODStatus(const_cast<MeasurementDet*>(it->second)))).first;


      GeomDet::SubDetector subdet = it->second->geomDet().subDetector();
      if (subdet == GeomDetEnumerators::PixelBarrel  || subdet == GeomDetEnumerators::PixelEndcap ){
    //special flag treatement for pixels
    //one can never be updated and not defined. except if we don't need to care about them: pixels
    inserted->second.defined=false;
    inserted->second.updated=true;
      }//pixel module
      else if (subdet == GeomDetEnumerators::TIB || subdet == GeomDetEnumerators::TOB ||
      subdet == GeomDetEnumerators::TID || subdet == GeomDetEnumerators::TEC )
    {
      //set flag to false
      inserted->second.defined=false;
      inserted->second.updated=false;

      //what will be the element index in the refgetter
      GlobalPoint center = it->second->geomDet().position();
      double eta = center.eta();
      double phi = center.phi();
      uint32_t id = it->first.rawId();
      SiStripRegionCabling::ElementIndex eIndex = theStripRegionCabling->elementIndex(SiStripRegionCabling::Position(eta,phi),
                                              SiStripRegionCabling::subdetFromDetId(id),
                                              SiStripRegionCabling::layerFromDetId(id));
      LogDebug(category_)<<"region selected (from "<<id<<" center) is:\n"
                 <<"position: "<<center
                 <<"\n center absolute index: "<<theStripRegionCabling->region(theStripRegionCabling->positionIndex(SiStripRegionCabling::Position(eta,phi)))
                 <<"\n center position index: "<<theStripRegionCabling->positionIndex(SiStripRegionCabling::Position(eta,phi)).first<<
        " "<<theStripRegionCabling->positionIndex(SiStripRegionCabling::Position(eta,phi)).second
                 <<"\n center postion: "<<theStripRegionCabling->position(theStripRegionCabling->positionIndex(SiStripRegionCabling::Position(eta,phi))).first<<
        " "<<theStripRegionCabling->position(theStripRegionCabling->positionIndex(SiStripRegionCabling::Position(eta,phi))).second
                 <<"\n eta: "<<eta
                 <<"\n phi: "<<phi
                 <<"\n subedet: "<<SiStripRegionCabling::subdetFromDetId(id)
                 <<" layer: "<<SiStripRegionCabling::layerFromDetId(id);

      //      register those in a map
      //to be able to know what are the detid in a given elementIndex
      region_mapping[eIndex].push_back(inserted);
    }//strip module
      else{
    //abort
    edm::LogError(category_)<<"not a tracker geomdet in constructor: "<<it->first.rawId();
    throw MeasurementDetException("OnDemandMeasurementTracker dealing with a non tracker GeomDet.");
      }//abort
    }//loop over DetMap
}


void OnDemandMeasurementTracker::define( const edm::Handle< LazyGetter> & theLazyGetterH,
                     std::auto_ptr< RefGetter > &  theGetter ) const
{
  //  define is supposed to be call by an EDProducer module, which wil put the RefGetter in the event
  //  so that reference can be made to it.
  //  the lazy getter is retrieved by the calling module and passed along with the event
  //  the map is cleared, except for pixel
  //  then the known elementIndex are defined to the RefGetter. no unpacking is done at this time
  //  the defined region range is registered in the DetODMap for further use.

  //clear all defined tags to start from scratch (except for pixel)
  for (DetODContainer::iterator it=theDetODMap.begin(); it!= theDetODMap.end();++it)
    {
      if (it->second.updated && !it->second.defined) continue; //special treatement for pixels
      it->second.defined= false; 
      it->second.updated = false;
    }

  //define all the elementindex in the refgetter
  for( std::map<SiStripRegionCabling::ElementIndex, std::vector<DetODContainer::iterator> >::iterator eIt= region_mapping.begin();
       eIt!=region_mapping.end();++eIt){
    std::pair<unsigned int, unsigned int> region_range; 
    
    //before update of the refgetter
    region_range.first = theGetter->size();
    //update the refegetter with the elementindex
    theStripRegionCabling->updateSiStripRefGetter<SiStripCluster> (*theGetter, theLazyGetterH, eIt->first);
    //after update of the refgetter
    region_range.second = theGetter->size();

    LogDebug(category_)<<"between index: "<<region_range.first<<" "<<region_range.second
               <<"\n"<<dumpRegion(region_range,*theGetter,StayPacked_);
    
    //now assign to each measurement det for that element index
    for (std::vector<DetODContainer::iterator>::iterator dIt=eIt->second.begin();
     dIt!=eIt->second.end();++dIt){
      (*dIt)->second.region_range = region_range;
      (*dIt)->second.defined=true;
      LogDebug(category_)<<"detId: "<<(*dIt)->first.rawId()<<" in region range: "<<region_range.first<<" "<<region_range.second;
    }//loop over MeasurementDet attached to that elementIndex
  }//loop over know elementindex
}

void OnDemandMeasurementTracker::updateStrips( const edm::Event& event) const 
{
  bool oncePerEvent= edm::Service<UpdaterService>()->checkOnce("OnDemandMeasurementTracker::updateStrips::"+name_);
  bool failedToGet = false;
  if (!oncePerEvent)
    failedToGet = theRefGetterH.failedToGet() || theLazyGetterH.failedToGet();

  if (oncePerEvent || failedToGet)
    {
      LogDebug(category_)<<"Updating siStrip on event: "<< (unsigned int) event.id().run() <<" : "<<(unsigned int) event.id().event();
      
      //get the ref getter back from the event
      std::string stripClusterProducer = pset_.getParameter<std::string>("stripClusterProducer");
      event.getByLabel(stripClusterProducer,theRefGetterH);
      
      std::string stripLazyGetter = pset_.getParameter<std::string>("stripLazyGetterProducer");
      event.getByLabel(stripLazyGetter,theLazyGetterH);
    }
}

void OnDemandMeasurementTracker::update( const edm::Event& event) const
{
  //  update is supposed to be called by any module that is useing the MeasurementTracker
  //  after checking the the event has not yet been seen
  //  update the pixel using MeasurementTracekr specific function
  //  retreive the RefGetter from the event: the very one that has been pass to define(...) and put into the event

  if (!PixelOnDemand_) {
    LogDebug(category_)<<"pixel are not OnDemand. updating them a la MeasurmentTracker.";
    MeasurementTracker::updatePixels(event);}
  else{
    edm::LogError(category_)<<"trying to update siPixel as on-demand. Not Implemented yet.";
  }

  if (!StripOnDemand_) {
    LogDebug(category_)<<"strip are not OnDemand. updating them a la MeasurmentTracker.";
    MeasurementTracker::updateStrips(event);}
  else{
    LogDebug(category_)<<"strip are OnDemand. updating them a la OnDemandMeasurmentTracker."; 
    updateStrips(event);
  }
}

#include <sstream>

std::string OnDemandMeasurementTracker::dumpCluster(const std::vector<SiStripCluster> ::const_iterator & begin,const  std::vector<SiStripCluster> ::const_iterator & end)const
{
  //  dumpCluster is a printout of all the clusters between the iterator. returns a string
  std::string tab="      ";
  std::stringstream ss;
  std::vector<SiStripCluster> ::const_iterator it = begin;
  unsigned int i=0;
  for (;it!=end;++it){
    ss<<tab<<i++<<") center: "<<it->barycenter()<<",id: "<<it->geographicalId()<<" with: "<<it->amplitudes().size()<<" strips\n"<<tab<<tab<<"{";
    for (unsigned int is=0;is!=it->amplitudes().size();++is){
      ss<<it->amplitudes()[is]<<" ";
    }ss<<"}\n";
  }
  return ss.str();
}

std::string OnDemandMeasurementTracker::dumpRegion(std::pair<unsigned int,unsigned int> indexes,
                         const RefGetter & theGetter,
                         bool stayPacked)const
{
  //  dumpRegion is a printout of all the clusters in a region defined on the RefGetter. returns a string
  std::stringstream ss;
  ss<<"cluster between: "<<indexes.first<<" and: "<<indexes.second<<"\n";
  for (unsigned int iRegion = indexes.first; iRegion != indexes.second; ++iRegion){    
    uint32_t reg = SiStripRegionCabling::region((theGetter)[iRegion].region());
    SiStripRegionCabling::Position pos = theStripRegionCabling->position(reg);
    SiStripRegionCabling::PositionIndex posI = theStripRegionCabling->positionIndex(reg);
    
    ss<<"Clusters for region:["<<iRegion<<"]"
      <<"\n element index: "<<(theGetter)[iRegion].region()
      <<"\n region absolute index: "<<reg
      <<"\n region position index: "<<posI.first<<" "<<posI.second
      <<"\n region position: "<<pos.first<<" "<<pos.second
      <<"\n"<< (stayPacked? " hidden to avoid unpacking." : dumpCluster((theGetter)[iRegion].begin(),(theGetter)[iRegion].end()));
  }
  return ss.str();
}

void OnDemandMeasurementTracker::assign(const TkStripMeasurementDet * csmdet,
                  DetODContainer::iterator * alreadyFound)const {
  //  assign is using the handle to the refgetter and the region index range to update the MeasurementDet with their clusters
  
  TkStripMeasurementDet * smdet = const_cast<TkStripMeasurementDet *>(csmdet);
  DetId id = smdet->geomDet().geographicalId();
  
  LogDebug(category_)<<"assigning: "<<id.rawId();

  // what is the iterator. do not look again if already found and provided with
  DetODContainer::iterator elementInMap;
  if (alreadyFound){ elementInMap=*alreadyFound;}
  else{ elementInMap = theDetODMap.find(id);}
  
  if  ( elementInMap != theDetODMap.end()){
    //flag it as updated
    elementInMap->second.updated = true;
    //retrieve the region range index for this module
    std::pair<unsigned int,unsigned int> & indexes =elementInMap->second.region_range;

    //this printout will trigger unpacking. no problem. it is done on the next regular line (find(id.rawId())
    LogDebug(category_)<<"between index: "<<indexes.first<<" and: "<<indexes.second
               <<"\nretrieved for module: "<<id.rawId()
               <<"\n"<<dumpRegion(indexes,*theRefGetterH);
    
    //look for iterator range in the regions defined for that module
    for (unsigned int iRegion = indexes.first; iRegion != indexes.second; ++iRegion){
      RefGetter::record_pair range = (*theRefGetterH)[iRegion].find(id.rawId());
      if (range.first!=range.second){
    //  found something not empty
    //update the measurementDet
    smdet->update(range.first, range.second, theLazyGetterH, id);
    LogDebug(category_)<<"Valid clusters for: "<<id.rawId()
               <<"\nnumber of regions defined here: "<< indexes.second-indexes.first
               <<"\n"<<dumpCluster(range.first,range.second);
    //and you are done
    return;}
    }//loop over regions, between indexes

    //if reached. no cluster are found. set the TkStripMeasurementDet to be empty
    smdet->setEmpty();

  }//found in the map
  else{
    //throw excpetion
    edm::LogError(category_)<<"failed to find the MeasurementDet for: "<<id.rawId();
    throw MeasurementDetException("failed to find the MeasurementDet for: <see message logger>");
  }
}



const MeasurementDet* 
OnDemandMeasurementTracker::idToDet(const DetId& id) const
{
  //  overloaded from MeasurementTracker
  //  find the detid. if not found throw exception
  //  if already updated: always for pixel or strip already queried. return it
  DetODContainer::iterator it = theDetODMap.find(id);
  if ( it != theDetODMap.end()) {
    
    //it has already been queried. and so already updated: nothing to be done here (valid for pixels too)
    if (it->second.updated){LogDebug(category_)<<"found id: "<<id.rawId()<<" as aleardy updated."; return it->second.mdet;}
    
    if (StripOnDemand_){
      //check glued or single
      std::vector< const GeomDet*> comp = it->second.mdet->geomDet().components();
      if (!comp.empty()){
    //glued det
    LogDebug(category_)<<"updating glued id: "<<id.rawId()<<" ("<<comp.size()<<").";
    //cast to specific type
    TkGluedMeasurementDet*  theConcreteDet = dynamic_cast<TkGluedMeasurementDet*>(it->second.mdet);
    //FIXME. do not require the dynamic cast, since we know what we do here
    //              const TkGluedMeasurementDet*  theConcreteDet = (const TkGluedMeasurementDet*) it->second.mdet;
    if (!theConcreteDet){
      edm::LogError(category_)<<"failed to cast to glued measurement det:"<<it->first.rawId();
      throw MeasurementDetException("failed to cast. see log file.");}
        
    //  update the two components
    //update the mono
    assign(theConcreteDet->monoDet());
    //update the stereo
    assign(theConcreteDet->stereoDet());
          
    //flag the glued det as updated (components are flagged in assign)
    it->second.updated=true;
      }//glued det
      else{
    //single layer 
    LogDebug(category_)<<"updating singel id: "<<id.rawId();
    //cast to specific type
    TkStripMeasurementDet*  theConcreteDet = dynamic_cast<TkStripMeasurementDet*>(it->second.mdet);
    //FIXME. do not require the dynamic cast, since we know what we do here
    //              const TkStripMeasurementDet*  theConcreteDet = (const TkStripMeasurementDet*)it->second.mdet;
    if (!theConcreteDet){
      edm::LogError(category_)<<"failed to cast to strip measurement det:"<<it->first.rawId();
      throw MeasurementDetException("failed to cast. see log file.");}

    //update the single layer
    assign(theConcreteDet,&it);
      }//single det
    }
    //eventually return it
    return it->second.mdet;
  }//found in the DetODMap
  else{
    //throw excpetion
    edm::LogError(category_)<<"failed to find the MeasurementDet for: "<<id.rawId();
    throw MeasurementDetException("failed to find the MeasurementDet for: <see message logger>");
  }
  return 0;
}