CMSSW Reference Manual

00001 // -*- C++ -*-
00002 // Package:    SiStripMonitorCluster
00003 // Class:      SiStripMonitorCluster
00006 // Original Author:  Dorian Kcira
00007 //         Created:  Wed Feb  1 16:42:34 CET 2006
00008 // $Id: SiStripMonitorCluster.cc,v 1.49 2008/11/25 17:15:09 dutta Exp $
00009 #include <vector>
00010 #include <numeric>
00011 #include <fstream>
00012 #include <math.h>
00013 #include "TNamed.h"
00014 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00015 #include "FWCore/ServiceRegistry/interface/Service.h"
00016 
00017 #include "CalibTracker/Records/interface/SiStripDetCablingRcd.h"
00018 #include "CalibFormats/SiStripObjects/interface/SiStripDetCabling.h"
00019 #include "CondFormats/DataRecord/interface/SiStripNoisesRcd.h"
00020 #include "CondFormats/SiStripObjects/interface/SiStripNoises.h"
00021 #include "CalibTracker/Records/interface/SiStripGainRcd.h"
00022 #include "CalibFormats/SiStripObjects/interface/SiStripGain.h"
00023 #include "CalibTracker/Records/interface/SiStripQualityRcd.h"
00024 #include "CalibFormats/SiStripObjects/interface/SiStripQuality.h"
00025 
00026 #include "DataFormats/Common/interface/DetSetVector.h"
00027 #include "DataFormats/Common/interface/DetSetVectorNew.h"
00028 #include "DataFormats/SiStripCluster/interface/SiStripCluster.h"
00029 #include "DQM/SiStripCommon/interface/SiStripFolderOrganizer.h"
00030 #include "DQM/SiStripCommon/interface/SiStripHistoId.h"
00031 #include "DQM/SiStripMonitorCluster/interface/SiStripMonitorCluster.h"
00032 #include "DQMServices/Core/interface/DQMStore.h"
00033 #include "DQMServices/Core/interface/MonitorElement.h"
00034 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
00035 #include "DataFormats/SiStripDetId/interface/SiStripSubStructure.h"
00036 #include "DataFormats/SiStripDetId/interface/TECDetId.h"
00037 #include "DataFormats/SiStripDetId/interface/TIBDetId.h"
00038 #include "DataFormats/SiStripDetId/interface/TIDDetId.h"
00039 #include "DataFormats/SiStripDetId/interface/TOBDetId.h"
00040 
00041 #include "TMath.h"
00042 #include <iostream>
00043 
00044 //--------------------------------------------------------------------------------------------
00045 SiStripMonitorCluster::SiStripMonitorCluster(const edm::ParameterSet& iConfig) : dqmStore_(edm::Service<DQMStore>().operator->()), conf_(iConfig), show_mechanical_structure_view(true), show_readout_view(false), show_control_view(false), select_all_detectors(false), reset_each_run(false), m_cacheID_(0)
00046 {
00047 
00048   firstEvent = -1;
00049   eventNb = 0;
00050 
00051 
00052   //get on/off option for every cluster from cfi
00053   edm::ParameterSet ParametersnClusters =  conf_.getParameter<edm::ParameterSet>("TH1nClusters");
00054   layerswitchncluson = ParametersnClusters.getParameter<bool>("layerswitchon");
00055   moduleswitchncluson = ParametersnClusters.getParameter<bool>("moduleswitchon");
00056   
00057   edm::ParameterSet ParametersClusterCharge =  conf_.getParameter<edm::ParameterSet>("TH1ClusterCharge");
00058   layerswitchcluschargeon = ParametersClusterCharge.getParameter<bool>("layerswitchon");
00059   moduleswitchcluschargeon = ParametersClusterCharge.getParameter<bool>("moduleswitchon");
00060   
00061   edm::ParameterSet ParametersClusterStoN =  conf_.getParameter<edm::ParameterSet>("TH1ClusterStoN");
00062   layerswitchclusstonon = ParametersClusterStoN.getParameter<bool>("layerswitchon");
00063   moduleswitchclusstonon = ParametersClusterStoN.getParameter<bool>("moduleswitchon");
00064   
00065   edm::ParameterSet ParametersClusterPos =  conf_.getParameter<edm::ParameterSet>("TH1ClusterPos");
00066   layerswitchclusposon = ParametersClusterPos.getParameter<bool>("layerswitchon");
00067   moduleswitchclusposon = ParametersClusterPos.getParameter<bool>("moduleswitchon");
00068   
00069   edm::ParameterSet ParametersClusterNoise =  conf_.getParameter<edm::ParameterSet>("TH1ClusterNoise");
00070   layerswitchclusnoiseon = ParametersClusterNoise.getParameter<bool>("layerswitchon");
00071   moduleswitchclusnoiseon = ParametersClusterNoise.getParameter<bool>("moduleswitchon");
00072   
00073   edm::ParameterSet ParametersClusterWidth =  conf_.getParameter<edm::ParameterSet>("TH1ClusterWidth");
00074   layerswitchcluswidthon = ParametersClusterWidth.getParameter<bool>("layerswitchon");
00075   moduleswitchcluswidthon = ParametersClusterWidth.getParameter<bool>("moduleswitchon");
00076   
00077   edm::ParameterSet ParametersModuleLocalOccupancy =  conf_.getParameter<edm::ParameterSet>("TH1ModuleLocalOccupancy");
00078   layerswitchlocaloccupancy = ParametersModuleLocalOccupancy.getParameter<bool>("layerswitchon");
00079   moduleswitchlocaloccupancy = ParametersModuleLocalOccupancy.getParameter<bool>("moduleswitchon");
00080 
00081   edm::ParameterSet ParametersNrOfClusterizedStrips =  conf_.getParameter<edm::ParameterSet>("TH1NrOfClusterizedStrips");
00082   layerswitchnrclusterizedstrip = ParametersNrOfClusterizedStrips.getParameter<bool>("layerswitchon");
00083   moduleswitchnrclusterizedstrip = ParametersNrOfClusterizedStrips.getParameter<bool>("moduleswitchon");
00084 
00085   edm::ParameterSet ParametersClusterProf = conf_.getParameter<edm::ParameterSet>("TProfNumberOfCluster");
00086   layerswitchnumclusterprofon = ParametersClusterProf.getParameter<bool>("layerswitchon");
00087 
00088   edm::ParameterSet ParametersClusterWidthProf = conf_.getParameter<edm::ParameterSet>("TProfClusterWidth");
00089   layerswitchclusterwidthprofon = ParametersClusterWidthProf.getParameter<bool>("layerswitchon");
00090 
00091   edm::ParameterSet ParametersTotDigiProf = conf_.getParameter<edm::ParameterSet>("TProfTotalNumberOfClusters");
00092   subdetswitchtotclusterprofon = ParametersTotDigiProf.getParameter<bool>("subdetswitchon");
00093 
00094   edm::ParameterSet ParametersDetsOn =  conf_.getParameter<edm::ParameterSet>("detectorson");
00095   tibon = ParametersDetsOn.getParameter<bool>("tibon");
00096   tidon = ParametersDetsOn.getParameter<bool>("tidon");
00097   tobon = ParametersDetsOn.getParameter<bool>("tobon");
00098   tecon = ParametersDetsOn.getParameter<bool>("tecon");
00099 
00100   createTrendMEs = conf_.getParameter<bool>("CreateTrendMEs");
00101   Mod_On_ = conf_.getParameter<bool>("Mod_On");
00102 } 
00103 
00104 
00105 SiStripMonitorCluster::~SiStripMonitorCluster() { }
00106 
00107 //--------------------------------------------------------------------------------------------
00108 void SiStripMonitorCluster::beginRun(const edm::Run& run, const edm::EventSetup& es){
00109 
00110   if (show_mechanical_structure_view) {
00111     unsigned long long cacheID = es.get<SiStripDetCablingRcd>().cacheIdentifier();
00112     if (m_cacheID_ != cacheID) {
00113       m_cacheID_ = cacheID;       
00114      edm::LogInfo("SiStripMonitorCluster") <<"SiStripMonitorCluster::beginRun: " 
00115                                             << " Creating MEs for new Cabling ";     
00116 
00117       createMEs(es);
00118     } 
00119   } else if (reset_each_run) {
00120     edm::LogInfo("SiStripMonitorCluster") <<"SiStripMonitorCluster::beginRun: " 
00121                                           << " Resetting MEs ";        
00122     for (std::map<uint32_t, ModMEs >::const_iterator idet = ModuleMEMap.begin() ; idet!=ModuleMEMap.end() ; idet++) {
00123       ResetModuleMEs(idet->first);
00124     }
00125   }
00126 }
00127 
00128 //--------------------------------------------------------------------------------------------
00129 void SiStripMonitorCluster::createMEs(const edm::EventSetup& es){
00130   if ( show_mechanical_structure_view ){
00131     // take from eventSetup the SiStripDetCabling object - here will use SiStripDetControl later on
00132     es.get<SiStripDetCablingRcd>().get(SiStripDetCabling_);
00133     
00134     // get list of active detectors from SiStripDetCabling 
00135     std::vector<uint32_t> activeDets;
00136     SiStripDetCabling_->addActiveDetectorsRawIds(activeDets);
00137     
00138     std::vector<uint32_t> SelectedDetIds;
00139     SiStripSubStructure substructure;
00140     if(select_all_detectors){
00141       // select all detectors if appropriate flag is set,  for example for the mtcc
00142       SelectedDetIds = activeDets;
00143     }else{
00144       // use SiStripSubStructure for selecting certain regions
00145       if(tibon) substructure.getTIBDetectors(activeDets, SelectedDetIds, 0, 0, 0, 0); // this adds rawDetIds to SelectedDetIds
00146       if(tobon) substructure.getTOBDetectors(activeDets, SelectedDetIds, 0, 0, 0);    // this adds rawDetIds to SelectedDetIds
00147       if(tidon) substructure.getTIDDetectors(activeDets, SelectedDetIds, 0, 0, 0, 0); // this adds rawDetIds to SelectedDetIds
00148       if(tecon) substructure.getTECDetectors(activeDets, SelectedDetIds, 0, 0, 0, 0, 0, 0); // this adds rawDetIds to SelectedDetIds
00149       
00150     }
00151     SiStripFolderOrganizer folder_organizer;
00152     folder_organizer.setSiStripFolder();
00153 
00154     // loop over detectors and book MEs
00155     edm::LogInfo("SiStripTkDQM|SiStripMonitorCluster")<<"nr. of SelectedDetIds:  "<<SelectedDetIds.size();
00156     for(std::vector<uint32_t>::iterator detid_iterator = SelectedDetIds.begin(); detid_iterator!=SelectedDetIds.end(); detid_iterator++){
00157       uint32_t detid = (*detid_iterator);
00158       // remove any eventual zero elements - there should be none, but just in case
00159       if(detid == 0) {
00160         SelectedDetIds.erase(detid_iterator);
00161         continue;
00162       }
00163 
00164       if (Mod_On_) {
00165         ModMEs mod_single;
00166         // set appropriate folder using SiStripFolderOrganizer
00167         folder_organizer.setDetectorFolder(detid); // pass the detid to this method
00168         if (reset_each_run) ResetModuleMEs(detid);
00169         createModuleMEs(mod_single, detid);
00170         // append to ModuleMEMap
00171         ModuleMEMap.insert( std::make_pair(detid, mod_single));
00172       }
00173 
00174       // Created Layer Level MEs if they are not created already
00175       std::pair<std::string,int32_t> det_layer_pair = folder_organizer.GetSubDetAndLayer(detid);
00176       if (DetectedLayers.find(det_layer_pair) == DetectedLayers.end()){
00177         DetectedLayers[det_layer_pair]=true;
00178 
00179         int32_t lnumber = det_layer_pair.second;
00180         std::vector<uint32_t> layerDetIds;
00181         if (det_layer_pair.first == "TIB") {
00182           substructure.getTIBDetectors(SelectedDetIds,layerDetIds,lnumber,0,0,0);
00183           if (subdetswitchtotclusterprofon && SubDetMEsMap.find("TIB") == SubDetMEsMap.end()) createSubDetMEs("TIB");
00184         } else if (det_layer_pair.first == "TOB") {
00185           substructure.getTOBDetectors(SelectedDetIds,layerDetIds,lnumber,0,0);
00186           if (subdetswitchtotclusterprofon && SubDetMEsMap.find("TOB") == SubDetMEsMap.end()) createSubDetMEs("TOB");
00187         } else if (det_layer_pair.first == "TID" && lnumber > 0) {
00188           substructure.getTIDDetectors(SelectedDetIds,layerDetIds,2,abs(lnumber),0,0);
00189           if (subdetswitchtotclusterprofon && SubDetMEsMap.find("TID") == SubDetMEsMap.end()) createSubDetMEs("TID");
00190         } else if (det_layer_pair.first == "TID" && lnumber < 0) {
00191           substructure.getTIDDetectors(SelectedDetIds,layerDetIds,1,abs(lnumber),0,0);
00192           if (subdetswitchtotclusterprofon && SubDetMEsMap.find("TID") == SubDetMEsMap.end()) createSubDetMEs("TID");
00193         } else if (det_layer_pair.first == "TEC" && lnumber > 0) {
00194           substructure.getTECDetectors(SelectedDetIds,layerDetIds,2,abs(lnumber),0,0,0,0);
00195           if (subdetswitchtotclusterprofon && SubDetMEsMap.find("TEC") == SubDetMEsMap.end()) createSubDetMEs("TEC");
00196         } else if (det_layer_pair.first == "TEC" && lnumber < 0) {
00197           substructure.getTECDetectors(SelectedDetIds,layerDetIds,1,abs(lnumber),0,0,0,0);
00198           if (subdetswitchtotclusterprofon && SubDetMEsMap.find("TEC") == SubDetMEsMap.end()) createSubDetMEs("TEC");
00199         }
00200 
00201         std::string label;
00202         getLayerLabel(detid, label);
00203         LayerDetMap[label] = layerDetIds;
00204 
00205         // book Layer plots
00206         folder_organizer.setLayerFolder(detid,det_layer_pair.second);
00207         createLayerMEs(label, layerDetIds.size());
00208       }
00209 
00210     }//end of loop over detectors
00211 
00212   }//end of if
00213 
00214 }//end of method
00215 
00216 
00217 
00218 //--------------------------------------------------------------------------------------------
00219 void SiStripMonitorCluster::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
00220 {
00221 
00222   using namespace edm;
00223 
00224   runNb   = iEvent.id().run();
00225   //   eventNb = iEvent.id().event();
00226   eventNb++;
00227   float iOrbit = iEvent.orbitNumber()/11223.0;
00228 
00229   //std::cout << " run " << iEvent.id().run() << runNb << " event " << iEvent.id().event() << eventNb << std::endl;
00230 
00231   edm::ESHandle<SiStripNoises> noiseHandle;
00232   iSetup.get<SiStripNoisesRcd>().get(noiseHandle);
00233 
00234   edm::ESHandle<SiStripGain> gainHandle;
00235   iSetup.get<SiStripGainRcd>().get(gainHandle);
00236 
00237   std::string quality_label  = conf_.getParameter<std::string>("StripQualityLabel");
00238   edm::ESHandle<SiStripQuality> qualityHandle;
00239   iSetup.get<SiStripQualityRcd>().get(quality_label,qualityHandle);
00240 
00241   iSetup.get<SiStripDetCablingRcd>().get(SiStripDetCabling_);
00242 
00243   // retrieve producer name of input StripClusterCollection
00244   std::string clusterProducer = conf_.getParameter<std::string>("ClusterProducer");
00245   // get collection of DetSetVector of clusters from Event
00246   edm::Handle< edmNew::DetSetVector<SiStripCluster> > cluster_detsetvektor;
00247   iEvent.getByLabel(clusterProducer, cluster_detsetvektor);
00248 
00249   //if (!cluster_detsetvektor.isValid()) std::cout<<" collection not valid"<<std::endl;
00250   if (!cluster_detsetvektor.isValid()) return;
00251  
00252   int nTotClusterTIB = 0;
00253   int nTotClusterTOB = 0;
00254   int nTotClusterTEC = 0;
00255   int nTotClusterTID = 0;
00256 
00257   bool found_layer_me = false;
00258   for (std::map<std::string, std::vector< uint32_t > >::const_iterator iterLayer = LayerDetMap.begin();
00259        iterLayer != LayerDetMap.end(); iterLayer++) {
00260     
00261     std::string layer_label = iterLayer->first;
00262     
00263     std::vector< uint32_t > layer_dets = iterLayer->second;
00264     int ncluster_layer = 0;
00265     std::map<std::string, LayerMEs>::iterator iLayerME = LayerMEMap.find(layer_label);
00266     
00267     //get Layer MEs 
00268     LayerMEs layer_single;
00269     if(iLayerME != LayerMEMap.end()) {
00270        layer_single = iLayerME->second; 
00271        found_layer_me = true;
00272      } 
00273 
00274     bool found_module_me = false;
00275     uint16_t iDet = 0;
00276     // loop over all modules in the layer
00277     for (std::vector< uint32_t >::const_iterator iterDets = layer_dets.begin() ; 
00278          iterDets != layer_dets.end() ; iterDets++) {
00279       iDet++;
00280       // detid and type of ME
00281       uint32_t detid = (*iterDets);
00282       
00283       // DetId and corresponding set of MEs
00284       ModMEs mod_single;
00285       if (Mod_On_) {
00286         std::map<uint32_t, ModMEs >::iterator imodME = ModuleMEMap.find(detid);
00287         if (imodME != ModuleMEMap.end()) {
00288           mod_single = imodME->second;
00289           found_module_me = true;
00290         } 
00291       } else found_module_me = false;
00292 
00293       edmNew::DetSetVector<SiStripCluster>::const_iterator isearch = cluster_detsetvektor->find(detid); // search  clusters of detid
00294     
00295       if(isearch==cluster_detsetvektor->end()){
00296         if(found_module_me && moduleswitchncluson && (mod_single.NumberOfClusters)){
00297           (mod_single.NumberOfClusters)->Fill(0.); // no clusters for this detector module,fill histogram with 0
00298         }
00299         if (found_layer_me && layerswitchnumclusterprofon) layer_single.LayerNumberOfClusterProfile->Fill(iDet, 0.0);
00300         continue; // no clusters for this detid => jump to next step of loop
00301       }
00302       
00303       //cluster_detset is a structure, cluster_detset.data is a std::vector<SiStripCluster>, cluster_detset.id is uint32_t
00304       edmNew::DetSet<SiStripCluster> cluster_detset = (*cluster_detsetvektor)[detid]; // the statement above makes sure there exists an element with 'detid'
00305       
00306       if(moduleswitchncluson && found_module_me && (mod_single.NumberOfClusters != NULL)){ // nr. of clusters per module
00307         (mod_single.NumberOfClusters)->Fill(static_cast<float>(cluster_detset.size()));
00308       }
00309       if (found_layer_me && layerswitchnumclusterprofon) 
00310         layer_single.LayerNumberOfClusterProfile->Fill(iDet, static_cast<float>(cluster_detset.size()));
00311       ncluster_layer +=  cluster_detset.size();
00312       
00313       short total_clusterized_strips = 0;
00314       //
00315       
00316       SiStripNoises::Range detNoiseRange = noiseHandle->getRange(detid);
00317       SiStripApvGain::Range detGainRange = gainHandle->getRange(detid); 
00318       SiStripQuality::Range qualityRange = qualityHandle->getRange(detid);
00319       
00320       for(edmNew::DetSet<SiStripCluster>::const_iterator clusterIter = cluster_detset.begin(); clusterIter!= cluster_detset.end(); clusterIter++){
00321         const  edm::ParameterSet ps = conf_.getParameter<edm::ParameterSet>("ClusterConditions");
00322         const std::vector<uint8_t>& ampls = clusterIter->amplitudes();
00323         // cluster position
00324         float cluster_position = clusterIter->barycenter();
00325         // width defined as nr. of strips that belong to cluster
00326         short cluster_width    = ampls.size(); 
00327         // add nr of strips of this cluster to total nr. of clusterized strips
00328         total_clusterized_strips = total_clusterized_strips + cluster_width; 
00329         
00330         // cluster signal and noise from the amplitudes
00331         float cluster_signal = 0.0;
00332         float cluster_noise  = 0.0;
00333         int nrnonzeroamplitudes = 0;
00334         float noise2 = 0.0;
00335         float noise  = 0.0;
00336         for(uint iamp=0; iamp<ampls.size(); iamp++){
00337           if(ampls[iamp]>0){ // nonzero amplitude
00338             cluster_signal += ampls[iamp];
00339             try{
00340               if(!qualityHandle->IsStripBad(qualityRange, clusterIter->firstStrip()+iamp)){
00341                 noise = noiseHandle->getNoise(clusterIter->firstStrip()+iamp,detNoiseRange)/gainHandle->getStripGain(clusterIter->firstStrip()+iamp, detGainRange);
00342               }
00343             }
00344             catch(cms::Exception& e){
00345               edm::LogError("SiStripTkDQM|SiStripMonitorCluster|DB")<<" cms::Exception:  detid="<<detid<<" firstStrip="<<clusterIter->firstStrip()<<" iamp="<<iamp<<e.what();
00346             }
00347             noise2 += noise*noise;
00348             nrnonzeroamplitudes++;
00349           }
00350         } // End loop over cluster amplitude
00351         
00352         if (nrnonzeroamplitudes > 0) cluster_noise = sqrt(noise2/nrnonzeroamplitudes);
00353         
00354         if( ps.getParameter<bool>("On") &&
00355             (cluster_signal/cluster_noise < ps.getParameter<double>("minStoN") ||
00356              cluster_signal/cluster_noise > ps.getParameter<double>("maxStoN") ||
00357              cluster_width < ps.getParameter<double>("minWidth") ||
00358              cluster_width  > ps.getParameter<double>("maxWidth")) ) continue;  
00359         
00360         ClusterProperties cluster_properties;
00361         cluster_properties.charge    = cluster_signal;
00362         cluster_properties.position  = cluster_position;
00363         cluster_properties.width     = cluster_width;
00364         cluster_properties.noise     = cluster_noise;
00365         
00366         // Fill Module Level MEs
00367         if (found_module_me) fillModuleMEs(mod_single, cluster_properties);
00368         
00369         // Fill Layer Level MEs
00370         if (found_layer_me) {
00371           fillLayerMEs(layer_single, cluster_properties);
00372           if (layerswitchclusterwidthprofon) 
00373             layer_single.LayerClusterWidthProfile->Fill(iDet, cluster_width);
00374         }
00375       } // end loop over clusters
00376       
00377       short total_nr_strips = SiStripDetCabling_->nApvPairs(detid) * 2 * 128; // get correct # of avp pairs
00378       float local_occupancy = static_cast<float>(total_clusterized_strips)/static_cast<float>(total_nr_strips);
00379       if (found_module_me) {
00380         if(moduleswitchnrclusterizedstrip && mod_single.NrOfClusterizedStrips ){ // nr of clusterized strips
00381           mod_single.NrOfClusterizedStrips->Fill(static_cast<float>(total_clusterized_strips));
00382         }
00383         
00384         if(moduleswitchlocaloccupancy && mod_single.ModuleLocalOccupancy ){ // Occupancy
00385           mod_single.ModuleLocalOccupancy->Fill(local_occupancy);
00386         }
00387       }
00388       if (layerswitchlocaloccupancy && found_layer_me && layer_single.LayerLocalOccupancy) {
00389         fillME(layer_single.LayerLocalOccupancy,local_occupancy);
00390         if (createTrendMEs) fillTrend(layer_single.LayerLocalOccupancyTrend,local_occupancy);
00391       }
00392     }
00393     if (layer_label.find("TIB") != std::string::npos)      nTotClusterTIB += ncluster_layer;
00394     else if (layer_label.find("TOB") != std::string::npos) nTotClusterTOB += ncluster_layer;
00395     else if (layer_label.find("TEC") != std::string::npos) nTotClusterTEC += ncluster_layer;        
00396     else if (layer_label.find("TID") != std::string::npos) nTotClusterTID += ncluster_layer;        
00397   }
00398   if (subdetswitchtotclusterprofon) {
00399     for (std::map<std::string, MonitorElement*>::iterator it = SubDetMEsMap.begin();
00400          it != SubDetMEsMap.end(); it++) {
00401       MonitorElement* me = it->second;
00402       if (!me) continue;
00403       if (it->first == "TIB") me->Fill(iOrbit,nTotClusterTIB);
00404       else if (it->first == "TOB") me->Fill(iOrbit,nTotClusterTOB);
00405       else if (it->first == "TID") me->Fill(iOrbit,nTotClusterTID);
00406       else if (it->first == "TEC") me->Fill(iOrbit,nTotClusterTEC);      
00407     }
00408   }    
00409 }
00410 //
00411 // -- EndJob
00412 //
00413 void SiStripMonitorCluster::endJob(void){
00414   bool outputMEsInRootFile = conf_.getParameter<bool>("OutputMEsInRootFile");
00415   std::string outputFileName = conf_.getParameter<std::string>("OutputFileName");
00416  
00417   // save histos in a file
00418   if(outputMEsInRootFile) dqmStore_->save(outputFileName);
00419 }
00420 //
00421 // -- Reset MEs
00422 //------------------------------------------------------------------------------
00423 void SiStripMonitorCluster::ResetModuleMEs(uint32_t idet){
00424   std::map<uint32_t, ModMEs >::iterator pos = ModuleMEMap.find(idet);
00425   ModMEs mod_me = pos->second;
00426 
00427   if (moduleswitchncluson)            mod_me.NumberOfClusters->Reset();
00428   if (moduleswitchclusposon)          mod_me.ClusterPosition->Reset();
00429   if (moduleswitchclusstonVsposon)    mod_me.ClusterSignalOverNoiseVsPos->Reset();
00430   if (moduleswitchcluswidthon)        mod_me.ClusterWidth->Reset();
00431   if (moduleswitchcluschargeon)       mod_me.ClusterCharge->Reset();
00432   if (moduleswitchclusnoiseon)        mod_me.ClusterNoise->Reset();
00433   if (moduleswitchclusstonon)         mod_me.ClusterSignalOverNoise->Reset();
00434   if (moduleswitchlocaloccupancy)     mod_me.ModuleLocalOccupancy->Reset();
00435   if (moduleswitchnrclusterizedstrip) mod_me.NrOfClusterizedStrips->Reset(); 
00436 }
00437 //
00438 // -- Create Module Level MEs
00439 //
00440 void SiStripMonitorCluster::createModuleMEs(ModMEs& mod_single, uint32_t detid) {
00441 
00442   // use SistripHistoId for producing histogram id (and title)
00443   SiStripHistoId hidmanager;
00444   std::string hid;
00445   
00446   //nr. of clusters per module
00447   if(moduleswitchncluson) {
00448     hid = hidmanager.createHistoId("NumberOfClusters","det",detid);
00449     mod_single.NumberOfClusters = bookME1D("TH1nClusters", hid.c_str());
00450     dqmStore_->tag(mod_single.NumberOfClusters, detid);
00451     mod_single.NumberOfClusters->setAxisTitle("number of clusters in one detector module");
00452     mod_single.NumberOfClusters->getTH1()->StatOverflows(kTRUE);  // over/underflows in Mean calculation
00453   }
00454   
00455   //ClusterPosition
00456   if(moduleswitchclusposon) {
00457     short total_nr_strips = SiStripDetCabling_->nApvPairs(detid) * 2 * 128; // get correct # of avp pairs    
00458     hid = hidmanager.createHistoId("ClusterPosition","det",detid);
00459     mod_single.ClusterPosition = dqmStore_->book1D(hid, hid, total_nr_strips, 0.5, total_nr_strips+0.5);
00460     dqmStore_->tag(mod_single.ClusterPosition, detid);
00461     mod_single.ClusterPosition->setAxisTitle("cluster position [strip number +0.5]");
00462   }
00463   
00464   //ClusterWidth
00465   if(moduleswitchcluswidthon) {
00466     hid = hidmanager.createHistoId("ClusterWidth","det",detid);
00467     mod_single.ClusterWidth = bookME1D("TH1ClusterWidth", hid.c_str());
00468     dqmStore_->tag(mod_single.ClusterWidth, detid);
00469     mod_single.ClusterWidth->setAxisTitle("cluster width [nr strips]");
00470   }
00471   
00472   //ClusterCharge
00473   if(moduleswitchcluschargeon) {
00474     hid = hidmanager.createHistoId("ClusterCharge","det",detid);
00475     mod_single.ClusterCharge = bookME1D("TH1ClusterCharge", hid.c_str());
00476     dqmStore_->tag(mod_single.ClusterCharge, detid);
00477     mod_single.ClusterCharge->setAxisTitle("cluster charge [ADC]");
00478   }
00479   
00480   //ClusterNoise
00481   if(moduleswitchclusnoiseon) {
00482     hid = hidmanager.createHistoId("ClusterNoise","det",detid);
00483     mod_single.ClusterNoise = bookME1D("TH1ClusterNoise", hid.c_str());
00484     dqmStore_->tag(mod_single.ClusterNoise, detid);
00485     mod_single.ClusterNoise->setAxisTitle("cluster noise");
00486   }
00487   
00488   //ClusterSignalOverNoise
00489   if(moduleswitchclusstonon) {
00490     hid = hidmanager.createHistoId("ClusterSignalOverNoise","det",detid);
00491     mod_single.ClusterSignalOverNoise = bookME1D("TH1ClusterStoN", hid.c_str());
00492     dqmStore_->tag(mod_single.ClusterSignalOverNoise, detid);
00493     mod_single.ClusterSignalOverNoise->setAxisTitle("ratio of signal to noise for each cluster");
00494   }
00495 
00496   //ClusterSignalOverNoiseVsPos
00497   if(moduleswitchclusstonVsposon) {
00498     hid = hidmanager.createHistoId("ClusterSignalOverNoiseVsPos","det",detid);
00499     Parameters =  conf_.getParameter<edm::ParameterSet>("TH1ClusterStoNVsPos");
00500     mod_single.ClusterSignalOverNoiseVsPos= dqmStore_->bookProfile(hid.c_str(),hid.c_str(),
00501                                                                    Parameters.getParameter<int32_t>("Nbinx"),
00502                                                                    Parameters.getParameter<double>("xmin"),
00503                                                                    Parameters.getParameter<double>("xmax"),
00504                                                                    Parameters.getParameter<int32_t>("Nbiny"),
00505                                                                    Parameters.getParameter<double>("ymin"),
00506                                                                    Parameters.getParameter<double>("ymax")
00507                                                                    );
00508     dqmStore_->tag(mod_single.ClusterSignalOverNoiseVsPos, detid);
00509     mod_single.ClusterSignalOverNoiseVsPos->setAxisTitle("pos");
00510   }
00511   
00512   //ModuleLocalOccupancy
00513   if (moduleswitchlocaloccupancy) {
00514     hid = hidmanager.createHistoId("ClusterLocalOccupancy","det",detid);
00515     mod_single.ModuleLocalOccupancy = bookME1D("TH1ModuleLocalOccupancy", hid.c_str());
00516     dqmStore_->tag(mod_single.ModuleLocalOccupancy, detid);
00517     mod_single.ModuleLocalOccupancy->setAxisTitle("module local occupancy [% of clusterized strips]");
00518   }
00519   
00520   //NrOfClusterizedStrips
00521   if (moduleswitchnrclusterizedstrip) {
00522     hid = hidmanager.createHistoId("NrOfClusterizedStrips","det",detid);
00523     mod_single.NrOfClusterizedStrips = bookME1D("TH1NrOfClusterizedStrips", hid.c_str());
00524     dqmStore_->tag(mod_single.NrOfClusterizedStrips, detid);
00525     mod_single.NrOfClusterizedStrips->setAxisTitle("number of clusterized strips");
00526   }
00527 }  
00528 //
00529 // -- Create Module Level MEs
00530 //  
00531 void SiStripMonitorCluster::createLayerMEs(std::string label, int ndets) {
00532 
00533   std::map<std::string, LayerMEs>::iterator iLayerME  = LayerMEMap.find(label);
00534   if(iLayerME==LayerMEMap.end()){
00535     SiStripHistoId hidmanager;
00536     LayerMEs layerMEs; 
00537 
00538 
00539     //Cluster Width
00540     if(layerswitchcluswidthon) {
00541       layerMEs.LayerClusterWidth=bookME1D("TH1ClusterWidth", hidmanager.createHistoLayer("Summary_ClusterWidth","layer",label,"").c_str()); 
00542       if (createTrendMEs) layerMEs.LayerClusterWidthTrend=bookMETrend("TH1ClusterWidth", hidmanager.createHistoLayer("Trend_ClusterWidth","layer",label,"").c_str()); 
00543     }
00544 
00545     //Cluster Noise
00546     if(layerswitchclusnoiseon) {
00547       layerMEs.LayerClusterNoise=bookME1D("TH1ClusterNoise", hidmanager.createHistoLayer("Summary_ClusterNoise","layer",label,"").c_str()); 
00548       if (createTrendMEs) layerMEs.LayerClusterNoiseTrend=bookMETrend("TH1ClusterNoise", hidmanager.createHistoLayer("Trend_ClusterNoise","layer",label,"").c_str()); 
00549     }
00550 
00551     //Cluster Charge
00552     if(layerswitchcluschargeon) {
00553       layerMEs.LayerClusterCharge=bookME1D("TH1ClusterCharge", hidmanager.createHistoLayer("Summary_ClusterCharge","layer",label,"").c_str());
00554        if (createTrendMEs) layerMEs.LayerClusterChargeTrend=bookMETrend("TH1ClusterCharge", hidmanager.createHistoLayer("Trend_ClusterCharge","layer",label,"").c_str());
00555     }
00556 
00557     //Cluster StoN
00558     if(layerswitchclusstonon) {
00559       layerMEs.LayerClusterStoN=bookME1D("TH1ClusterStoN", hidmanager.createHistoLayer("Summary_ClusterSignalOverNoise","layer",label,"").c_str());
00560       if (createTrendMEs) layerMEs.LayerClusterStoNTrend=bookMETrend("TH1ClusterStoN", hidmanager.createHistoLayer("Trend_ClusterSignalOverNoise","layer",label,"").c_str());
00561     }
00562 
00563     //Cluster Occupancy
00564     if(layerswitchlocaloccupancy) {
00565       layerMEs.LayerLocalOccupancy=bookME1D("TH1ModuleLocalOccupancy", hidmanager.createHistoLayer("Summary_ClusterLocalOccupancy","layer",label,"").c_str());  
00566       if (createTrendMEs) layerMEs.LayerLocalOccupancyTrend=bookME1D("TH1ModuleLocalOccupancy", hidmanager.createHistoLayer("Trend_ClusterLocalOccupancy","layer",label,"").c_str());  
00567       
00568     }
00569 
00570     // # of Cluster Profile 
00571     if(layerswitchnumclusterprofon) {
00572       std::string hid = hidmanager.createHistoLayer("NumberOfClusterProfile","layer",label,"");
00573       layerMEs.LayerNumberOfClusterProfile = dqmStore_->bookProfile(hid, hid, ndets, 0.5, ndets+0.5,21, -0.5, 20.5);      
00574     }
00575 
00576     // Cluster Width Profile 
00577     if(layerswitchclusterwidthprofon) {
00578       std::string hid = hidmanager.createHistoLayer("ClusterWidthProfile","layer",label,"");      
00579       layerMEs.LayerClusterWidthProfile = dqmStore_->bookProfile(hid, hid, ndets, 0.5, ndets+0.5, 20, -0.5, 19.5);      
00580     }
00581 
00582     LayerMEMap[label]=layerMEs;
00583   }
00584 
00585 }
00586 //
00587 // -- Create SubDetector MEs
00588 //
00589 void SiStripMonitorCluster::createSubDetMEs(std::string label) {
00590   edm::ParameterSet Parameters =  conf_.getParameter<edm::ParameterSet>("TProfTotalNumberOfClusters");
00591   std::string HistoName;
00592   dqmStore_->cd("SiStrip/MechanicalView/"+label);
00593   HistoName = "TotalNumberOfClusterProfile__" + label;
00594   MonitorElement* me = dqmStore_->bookProfile(HistoName,HistoName,
00595                                               Parameters.getParameter<int32_t>("Nbins"),
00596                                               Parameters.getParameter<double>("xmin"),
00597                                               Parameters.getParameter<double>("xmax"),
00598                                               100, //that parameter should not be there !?
00599                                               Parameters.getParameter<double>("ymin"),
00600                                               Parameters.getParameter<double>("ymax"),
00601                                               "" );
00602   if (!me) return;
00603   me->setAxisTitle("Event Time (Seconds)",1);
00604   if (me->kind() == MonitorElement::DQM_KIND_TPROFILE) me->getTH1()->SetBit(TH1::kCanRebin);
00605   SubDetMEsMap.insert( std::make_pair(label, me));
00606 }
00607 
00608 //
00609 // -- Fill Module Level Histograms
00610 //
00611 void SiStripMonitorCluster::fillModuleMEs(ModMEs& mod_mes, ClusterProperties& cluster) {
00612   
00613   if(moduleswitchclusposon && (mod_mes.ClusterPosition)) // position of cluster
00614     (mod_mes.ClusterPosition)->Fill(cluster.position);
00615 
00616   if(moduleswitchcluswidthon && (mod_mes.ClusterWidth)) // width of cluster
00617     (mod_mes.ClusterWidth)->Fill(static_cast<float>(cluster.width));
00618  
00619   if(moduleswitchclusstonon && (mod_mes.ClusterSignalOverNoise)) {// SignalToNoise
00620     if (cluster.noise > 0) 
00621       (mod_mes.ClusterSignalOverNoise)->Fill(cluster.charge/cluster.noise);
00622   }
00623 
00624   if(moduleswitchclusstonVsposon && (mod_mes.ClusterSignalOverNoiseVsPos)) {// SignalToNoise
00625     if (cluster.noise > 0) 
00626       (mod_mes.ClusterSignalOverNoiseVsPos)->Fill(cluster.position,cluster.charge/cluster.noise);
00627   }
00628 
00629   if(moduleswitchclusnoiseon && (mod_mes.ClusterNoise))  // Noise
00630     (mod_mes.ClusterNoise)->Fill(cluster.noise);
00631 
00632   if(moduleswitchcluschargeon && (mod_mes.ClusterCharge)) // charge of cluster
00633     (mod_mes.ClusterCharge)->Fill(cluster.charge);
00634   
00635 } 
00636 //
00637 // -- Fill Layer Level MEs
00638 //
00639 void SiStripMonitorCluster::fillLayerMEs(LayerMEs& layerMEs, ClusterProperties& cluster) { 
00640   if(layerswitchclusstonon) {
00641     fillME(layerMEs.LayerClusterStoN  ,cluster.charge/cluster.noise);
00642     if (createTrendMEs) fillTrend(layerMEs.LayerClusterStoNTrend,cluster.charge/cluster.noise);
00643   }
00644   
00645   if(layerswitchcluschargeon) {
00646     fillME(layerMEs.LayerClusterCharge,cluster.charge);
00647     if (createTrendMEs) fillTrend(layerMEs.LayerClusterChargeTrend,cluster.charge);
00648   }
00649   
00650   if(layerswitchclusnoiseon) {
00651     fillME(layerMEs.LayerClusterNoise ,cluster.noise);
00652     if (createTrendMEs) fillTrend(layerMEs.LayerClusterNoiseTrend,cluster.noise);
00653   }
00654   
00655   if(layerswitchcluswidthon) {
00656     fillME(layerMEs.LayerClusterWidth ,cluster.width);
00657     if (createTrendMEs) fillTrend(layerMEs.LayerClusterWidthTrend,cluster.width);
00658   }
00659   
00660 }
00661 //
00662 // -- Fill Trend
00663 //
00664 void SiStripMonitorCluster::fillTrend(MonitorElement* me ,float value)
00665 {
00666   if(!me) return;
00667   //check the origin and check options
00668   int option = conf_.getParameter<edm::ParameterSet>("Trending").getParameter<int32_t>("UpdateMode");
00669   if(firstEvent==-1) firstEvent = eventNb;
00670   int CurrentStep = atoi(me->getAxisTitle(1).c_str()+8);
00671   int firstEventUsed = firstEvent;
00672   int presentOverflow = (int)me->getBinEntries(me->getNbinsX()+1);
00673   if(option==2) firstEventUsed += CurrentStep * int(me->getBinEntries(me->getNbinsX()+1));
00674   else if(option==3) firstEventUsed += CurrentStep * int(me->getBinEntries(me->getNbinsX()+1)) * me->getNbinsX();
00675   //fill
00676   me->Fill((eventNb-firstEventUsed)/CurrentStep,value);
00677 
00678   if(eventNb-firstEvent<1) return;
00679   // check if we reached the end
00680   if(presentOverflow == me->getBinEntries(me->getNbinsX()+1)) return;
00681   switch(option) {
00682   case 1:
00683     {
00684       // mode 1: rebin and change X scale
00685       int NbinsX = me->getNbinsX();
00686       float entries = 0.;
00687       float content = 0.;
00688       float error = 0.;
00689       int bin = 1;
00690       int totEntries = int(me->getEntries());
00691       for(;bin<=NbinsX/2;++bin) {
00692         content = (me->getBinContent(2*bin-1) + me->getBinContent(2*bin))/2.; 
00693         error   = pow((me->getBinError(2*bin-1)*me->getBinError(2*bin-1)) + (me->getBinError(2*bin)*me->getBinError(2*bin)),0.5)/2.; 
00694         entries = me->getBinEntries(2*bin-1) + me->getBinEntries(2*bin);
00695         me->setBinContent(bin,content*entries);
00696         me->setBinError(bin,error);
00697         me->setBinEntries(bin,entries);
00698       }
00699       for(;bin<=NbinsX+1;++bin) {
00700         me->setBinContent(bin,0);
00701         me->setBinError(bin,0);
00702         me->setBinEntries(bin,0); 
00703       }
00704       me->setEntries(totEntries);
00705       char buffer[256];
00706       sprintf(buffer,"EventId/%d",CurrentStep*2);
00707       me->setAxisTitle(std::string(buffer),1);
00708       break;
00709     }
00710   case 2:
00711     {
00712       // mode 2: slide
00713       int bin=1;
00714       int NbinsX = me->getNbinsX();
00715       for(;bin<=NbinsX;++bin) {
00716         me->setBinContent(bin,me->getBinContent(bin+1)*me->getBinEntries(bin+1));
00717         me->setBinError(bin,me->getBinError(bin+1));
00718         me->setBinEntries(bin,me->getBinEntries(bin+1));
00719       }
00720       break;
00721     }
00722   case 3:
00723     {
00724       // mode 3: reset
00725       int NbinsX = me->getNbinsX();
00726       for(int bin=0;bin<=NbinsX;++bin) {
00727         me->setBinContent(bin,0);
00728         me->setBinError(bin,0);
00729         me->setBinEntries(bin,0); 
00730       }
00731       break;
00732     }
00733   }
00734 }
00735 MonitorElement* SiStripMonitorCluster::bookMETrend(const char* ParameterSetLabel, const char* HistoName)
00736 {
00737   Parameters =  conf_.getParameter<edm::ParameterSet>(ParameterSetLabel);
00738   edm::ParameterSet ParametersTrend =  conf_.getParameter<edm::ParameterSet>("Trending");
00739   MonitorElement* me = dqmStore_->bookProfile(HistoName,HistoName,
00740                                               ParametersTrend.getParameter<int32_t>("Nbins"),
00741                                               //                                              0,
00742                                               ParametersTrend.getParameter<double>("xmin"),
00743                                               ParametersTrend.getParameter<double>("xmax"),
00744                                               //                                              ParametersTrend.getParameter<int32_t>("Nbins"),
00745                                               100, //that parameter should not be there !?
00746                                               ParametersTrend.getParameter<double>("ymin"),
00747                                               ParametersTrend.getParameter<double>("ymax"),
00748                                               "" );
00749   if(!me) return me;
00750   char buffer[256];
00751   sprintf(buffer,"EventId/%d",ParametersTrend.getParameter<int32_t>("Steps"));
00752   me->setAxisTitle(std::string(buffer),1);
00753   return me;
00754 }
00755 
00756 //------------------------------------------------------------------------------------------
00757 MonitorElement* SiStripMonitorCluster::bookME1D(const char* ParameterSetLabel, const char* HistoName)
00758 {
00759   Parameters =  conf_.getParameter<edm::ParameterSet>(ParameterSetLabel);
00760   return dqmStore_->book1D(HistoName,HistoName,
00761                            Parameters.getParameter<int32_t>("Nbinx"),
00762                            Parameters.getParameter<double>("xmin"),
00763                            Parameters.getParameter<double>("xmax")
00764                            );
00765 }
00766 //-------------------------------------------------------------------------------------------
00767 void SiStripMonitorCluster::getLayerLabel(uint32_t detid, std::string& label) {
00768   StripSubdetector subdet(detid);
00769   std::ostringstream label_str;
00770 
00771   if(subdet.subdetId() == StripSubdetector::TIB ){
00772     // ---------------------------  TIB  --------------------------- //
00773     TIBDetId tib1 = TIBDetId(detid);
00774     label_str << "TIB__layer__" << tib1.layer();
00775   }else if(subdet.subdetId() == StripSubdetector::TID){
00776     // ---------------------------  TID  --------------------------- //
00777     TIDDetId tid1 = TIDDetId(detid);
00778     label_str << "TID__side__" << tid1.side() << "__wheel__" << tid1.wheel();
00779   }else if(subdet.subdetId() == StripSubdetector::TOB){
00780     // ---------------------------  TOB  --------------------------- //
00781     TOBDetId tob1 = TOBDetId(detid);
00782     label_str << "TOB__layer__" << tob1.layer();
00783   }else if(subdet.subdetId() == StripSubdetector::TEC) {
00784     // ---------------------------  TEC  --------------------------- //
00785     TECDetId tec1 = TECDetId(detid);
00786     label_str << "TEC__side__"<<tec1.side() << "__wheel__" << tec1.wheel();
00787   }else{
00788     // ---------------------------  ???  --------------------------- //
00789     edm::LogError("SiStripTkDQM|WrongInput")<<"no such subdetector type :"<<subdet.subdetId()<<" no folder set!"<<std::endl;
00790     label_str << "";
00791   }
00792   label = label_str.str();
00793 }
00794 
00795