CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_4_1_8_patch9/src/CalibTracker/SiStripQuality/plugins/SiStripQualityHotStripIdentifier.cc

Go to the documentation of this file.
00001 #include "CalibTracker/SiStripQuality/plugins/SiStripQualityHotStripIdentifier.h"
00002 
00003 #include "CalibFormats/SiStripObjects/interface/SiStripQuality.h"
00004 
00005 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit2D.h"
00006 #include "DataFormats/TrackerRecHit2D/interface/ProjectedSiStripRecHit2D.h"
00007 #include "DataFormats/TrackerRecHit2D/interface/SiStripMatchedRecHit2D.h"
00008 #include "DataFormats/TrackReco/interface/TrackFwd.h"
00009 #include "DataFormats/TrackReco/interface/Track.h"
00010 #include "DataFormats/TrackingRecHit/interface/TrackingRecHit.h"
00011 #include <iostream>
00012 #include <fstream>
00013 #include <sstream>
00014 
00015 //Insert here the include to the algos
00016 #include "CalibTracker/SiStripQuality/interface/SiStripHotStripAlgorithmFromClusterOccupancy.h"
00017 
00018 
00019 
00020 SiStripQualityHotStripIdentifier::SiStripQualityHotStripIdentifier(const edm::ParameterSet& iConfig) : ConditionDBWriter<SiStripBadStrip>(iConfig),
00021   m_cacheID_(0), 
00022   dataLabel_(iConfig.getUntrackedParameter<std::string>("dataLabel","")),
00023   conf_(iConfig), 
00024   fp_(iConfig.getUntrackedParameter<edm::FileInPath>("file",edm::FileInPath("CalibTracker/SiStripCommon/data/SiStripDetInfo.dat"))),
00025   Cluster_src_(iConfig.getParameter<edm::InputTag>( "Cluster_src" )),
00026   Track_src_(iConfig.getUntrackedParameter<edm::InputTag>( "Track_src" )),
00027   tracksCollection_in_EventTree(iConfig.getUntrackedParameter<bool>("RemoveTrackClusters",false))
00028 {
00029   reader = new SiStripDetInfoFileReader(fp_.fullPath());  
00030 
00031   edm::ParameterSet pset=iConfig.getUntrackedParameter< edm::ParameterSet > ("ClusterSelection",edm::ParameterSet());
00032   MinClusterWidth_=pset.getUntrackedParameter<uint32_t>("minWidth",1);
00033   MaxClusterWidth_=pset.getUntrackedParameter<uint32_t>("maxWidth",1000);
00034 
00035   bookHistos();
00036 }
00037 
00038 SiStripQualityHotStripIdentifier::~SiStripQualityHotStripIdentifier(){
00039 }
00040 
00041 SiStripBadStrip* SiStripQualityHotStripIdentifier::getNewObject(){
00042 
00043   SiStripBadStrip* obj=new SiStripBadStrip();
00044   
00045   edm::ParameterSet parameters=conf_.getParameter<edm::ParameterSet>("AlgoParameters");
00046   std::string AlgoName = parameters.getParameter<std::string>("AlgoName");
00047   if (AlgoName=="SiStripHotStripAlgorithmFromClusterOccupancy"){
00048     
00049     edm::LogInfo("SiStripQualityHotStripIdentifier") <<" [SiStripQualityHotStripIdentifier::getNewObject] call to SiStripHotStripAlgorithmFromClusterOccupancy"<<std::endl;
00050 
00051     SiStripHotStripAlgorithmFromClusterOccupancy theIdentifier(conf_);
00052     theIdentifier.setProbabilityThreshold(parameters.getUntrackedParameter<double>("ProbabilityThreshold",1.E-7));
00053     theIdentifier.setMinNumEntries(parameters.getUntrackedParameter<uint32_t>("MinNumEntries",100));
00054     theIdentifier.setMinNumEntriesPerStrip(parameters.getUntrackedParameter<uint32_t>("MinNumEntriesPerStrip",5));
00055 
00056     SiStripQuality* qobj = new SiStripQuality();
00057     theIdentifier.extractBadStrips(qobj,ClusterPositionHistoMap,SiStripQuality_);
00058 
00059     edm::LogInfo("SiStripQualityHotStripIdentifier") <<" [SiStripQualityHotStripIdentifier::getNewObject] copy SiStripObject in SiStripBadStrip"<<std::endl;
00060 
00061     std::stringstream ss;  
00062   
00063     SiStripBadStrip::RegistryIterator rIter=qobj->getRegistryVectorBegin();
00064     SiStripBadStrip::RegistryIterator rIterEnd=qobj->getRegistryVectorEnd();
00065     for(;rIter!=rIterEnd;++rIter){
00066       SiStripBadStrip::Range range(qobj->getDataVectorBegin()+rIter->ibegin,qobj->getDataVectorBegin()+rIter->iend);
00067       if ( ! obj->put(rIter->detid,range) )
00068         edm::LogError("SiStripQualityHotStripIdentifier")<<"[SiStripQualityHotStripIdentifier::getNewObject] detid already exists"<<std::endl;
00069     }
00070     edm::LogInfo("SiStripQualityHotStripIdentifier") <<" [SiStripQualityHotStripIdentifier::getNewObject] " << ss.str() << std::endl;
00071 
00072   } else {
00073     edm::LogError("SiStripQualityHotStripIdentifier") <<" [SiStripQualityHotStripIdentifier::getNewObject] call for a unknow HotStrip identification algoritm"<<std::endl;
00074 
00075     std::vector<uint32_t> a;
00076     SiStripBadStrip::Range range(a.begin(),a.end());
00077     if ( ! obj->put(0xFFFFFFFF,range) )
00078       edm::LogError("SiStripQualityHotStripIdentifier")<<"[SiStripQualityHotStripIdentifier::getNewObject] detid already exists"<<std::endl;
00079   }
00080   
00081   return obj;
00082 }
00083 
00084 void SiStripQualityHotStripIdentifier::algoBeginRun(const edm::Run& run, const edm::EventSetup& iSetup){
00085   resetHistos(); 
00086   unsigned long long cacheID = iSetup.get<SiStripQualityRcd>().cacheIdentifier();
00087 
00088   if (m_cacheID_ == cacheID) 
00089     return;
00090 
00091   m_cacheID_ = cacheID; 
00092 
00093   iSetup.get<SiStripQualityRcd>().get(dataLabel_,SiStripQuality_);
00094 
00095 }
00096 
00097 void SiStripQualityHotStripIdentifier::algoEndJob(){
00098   //Clear map
00099   ClusterPositionHistoMap.clear();
00100 }
00101 
00102 void SiStripQualityHotStripIdentifier::resetHistos(){
00103   edm::LogInfo("SiStripQualityHotStripIdentifier") <<" [SiStripQualityHotStripIdentifier::resetHistos] " << std::endl;
00104   SiStrip::QualityHistosMap::iterator it=ClusterPositionHistoMap.begin();
00105   SiStrip::QualityHistosMap::iterator iEnd=ClusterPositionHistoMap.end();
00106   for (;it!=iEnd;++it){
00107     it->second->Reset();
00108   }
00109 }
00110 
00111 void SiStripQualityHotStripIdentifier::bookHistos(){
00112   edm::LogInfo("SiStripQualityHotStripIdentifier") <<" [SiStripQualityHotStripIdentifier::bookHistos] " << std::endl;
00113   char hname[1024];
00114   std::map<uint32_t, SiStripDetInfoFileReader::DetInfo >::const_iterator it =reader->getAllData().begin();
00115   std::map<uint32_t, SiStripDetInfoFileReader::DetInfo >::const_iterator iEnd =reader->getAllData().end();
00116   for(; it!=iEnd; ++it){    
00117     sprintf(hname,"h_%d",it->first);
00118     SiStrip::QualityHistosMap::iterator ref=ClusterPositionHistoMap.find(it->first);
00119     if (ref==ClusterPositionHistoMap.end()){
00120       ClusterPositionHistoMap[it->first]=boost::shared_ptr<TH1F>(new TH1F(hname,hname,it->second.nApvs*128,-0.5,it->second.nApvs*128-0.5));
00121     }
00122     else 
00123       edm::LogError("SiStripQualityHotStripIdentifier") <<" [SiStripQualityHotStripIdentifier::bookHistos] DetId " << it->first << " already found in map. Ignoring new data"<<std::endl;
00124   }
00125 }
00126 
00127 void SiStripQualityHotStripIdentifier::fillHisto(uint32_t detid,float value){
00128   
00129   SiStrip::QualityHistosMap::iterator ref=ClusterPositionHistoMap.find(detid);
00130   if (ref!=ClusterPositionHistoMap.end())
00131     ref->second->Fill(value);
00132   else
00133     edm::LogError("SiStripQualityHotStripIdentifier") <<" [SiStripQualityHotStripIdentifier::fillHisto] Histogram not found in the list for DetId " << detid << " Ignoring data value "<< value <<std::endl;
00134 }
00135 
00136 
00137 void SiStripQualityHotStripIdentifier::algoAnalyze(const edm::Event& e, const edm::EventSetup& eSetup){
00138   edm::Handle< edm::DetSetVector<SiStripCluster> >  dsv_SiStripCluster;
00139   e.getByLabel( Cluster_src_, dsv_SiStripCluster);    
00140       
00141   edm::Handle<reco::TrackCollection> trackCollection;
00142   if(tracksCollection_in_EventTree){
00143     e.getByLabel(Track_src_, trackCollection);
00144     if(!trackCollection.isValid()){
00145       edm::LogError("SiStripQualityHotStripIdentifier")<<" [SiStripQualityHotStripIdentifier::algoAnalyze] missing trackCollection with label " << Track_src_ <<std::endl;
00146     }
00147   }
00148 
00149   std::set<const void*> vPSiStripCluster;    
00150   //Perform track study
00151   if (tracksCollection_in_EventTree){
00152    
00153     const reco::TrackCollection tC = *(trackCollection.product());
00154     int i=0;
00155     for (reco::TrackCollection::const_iterator track=tC.begin(); track!=tC.end(); track++){
00156       LogTrace("SiStripQualityHotStripIdentifier")
00157         << "Track number "<< i+1 
00158         << "\n\tmomentum: " << track->momentum()
00159         << "\n\tPT: " << track->pt()
00160         << "\n\tvertex: " << track->vertex()
00161         << "\n\timpact parameter: " << track->d0()
00162         << "\n\tcharge: " << track->charge()
00163         << "\n\tnormalizedChi2: " << track->normalizedChi2() 
00164         <<"\n\tFrom EXTRA : "
00165         <<"\n\t\touter PT "<< track->outerPt()<<std::endl;
00166 
00167       //Loop on rechits
00168       for (trackingRecHit_iterator it = track->recHitsBegin();  it != track->recHitsEnd(); ++it){
00169         const TrackingRecHit* recHit = &(**it);
00170         
00171         if (!recHit->isValid()){
00172           LogTrace("SiStripQualityHotStripIdentifier") <<"\t\t Invalid Hit "<<std::endl;
00173           continue;
00174         }
00175         
00176         const SiStripRecHit2D* singleHit=dynamic_cast<const SiStripRecHit2D*>(recHit);
00177         const SiStripMatchedRecHit2D* matchedHit=dynamic_cast<const SiStripMatchedRecHit2D*>(recHit);
00178         const ProjectedSiStripRecHit2D* projectedHit=dynamic_cast<const ProjectedSiStripRecHit2D*>(recHit);
00179 
00180         if(matchedHit){
00181           vPSiStripCluster.insert((void *)&*(matchedHit->monoHit()->cluster())); 
00182           vPSiStripCluster.insert((void *)&*(matchedHit->stereoHit()->cluster()));         
00183         } else if(projectedHit){
00184           vPSiStripCluster.insert((void *)&*(projectedHit->originalHit().cluster())); 
00185         } else if(singleHit){
00186           vPSiStripCluster.insert((void*)&*(singleHit->cluster())); 
00187         }else{
00188           LogTrace("SiStripQualityHotStripIdentifier") << "NULL hit" << std::endl;
00189         }
00190       }
00191     }
00192   }
00193         
00194   std::stringstream ss;   
00195   //Loop on Det Clusters
00196   edm::DetSetVector<SiStripCluster>::const_iterator DSViter=dsv_SiStripCluster->begin();
00197   for (; DSViter!=dsv_SiStripCluster->end();DSViter++){
00198     edm::DetSet<SiStripCluster>::const_iterator ClusIter = DSViter->data.begin();
00199     edm::DetSet<SiStripCluster>::const_iterator ClusIterEnd = DSViter->data.end();
00200     for(; ClusIter!=ClusIterEnd; ++ClusIter) {    
00201       if (MinClusterWidth_<=ClusIter->amplitudes().size() && ClusIter->amplitudes().size()<=MaxClusterWidth_) {
00202         if (std::find(vPSiStripCluster.begin(),vPSiStripCluster.end(),(void*) &*ClusIter) == vPSiStripCluster.end()){
00203           if ( edm::isDebugEnabled() ) 
00204             ss << " adding cluster to histo for detid " << DSViter->id << " with barycenter "<< ClusIter->barycenter() << std::endl;
00205           fillHisto(DSViter->id,ClusIter->barycenter());
00206         }
00207       }
00208     }
00209   }
00210   LogTrace("SiStripQualityHotStripIdentifier") << ss.str();
00211 }