CMSSW Reference Manual

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