test/DBTools/DTVDriftAnalyzer.cc

0001
0002 /*
0003  *  See header file for a description of this class.
0004  *
0005  *  \author S. Bolognesi - INFN Torino
0006  */
0007
0008 #include "DTVDriftAnalyzer.h"
0009
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011 #include "FWCore/Framework/interface/EventSetup.h"
0012 #include "CondFormats/DTObjects/interface/DTMtime.h"
0013 #include "CondFormats/DataRecord/interface/DTMtimeRcd.h"
0014 #include "CondFormats/DTObjects/interface/DTRecoConditions.h"
0015 #include "CondFormats/DataRecord/interface/DTRecoConditionsVdriftRcd.h"
0016 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0017 #include <iostream>
0018 #include "TFile.h"
0019 #include "TH1D.h"
0020 #include "TString.h"
0021
0022 using namespace edm;
0023 using namespace std;
0024
0025 DTVDriftAnalyzer::DTVDriftAnalyzer(const ParameterSet& pset)
0026     : readLegacyVDriftDB(pset.getParameter<bool>("readLegacyVDriftDB")) {
0027   // The root file which will contain the histos
0028   string rootFileName = pset.getUntrackedParameter<string>("rootFileName");
0029   theFile = new TFile(rootFileName.c_str(), "RECREATE");
0030   theFile->cd();
0031   mTimeMapToken_ = esConsumes<edm::Transition::BeginRun>();
0032   vDriftMapToken_ = esConsumes<edm::Transition::BeginRun>();
0033 }
0034
0035 DTVDriftAnalyzer::~DTVDriftAnalyzer() { theFile->Close(); }
0036
0037 void DTVDriftAnalyzer::beginRun(const edm::Run& run, const edm::EventSetup& eventSetup) {
0038   if (readLegacyVDriftDB) {
0039     ESHandle<DTMtime> mTime = eventSetup.getHandle(mTimeMapToken_);
0040     mTimeMap = &*mTime;
0041     vDriftMap_ = nullptr;
0042     edm::LogVerbatim("DTVDriftAnalyzer") << "[DTVDriftAnalyzer] MTime version: " << mTime->version() << endl;
0043   } else {
0044     ESHandle<DTRecoConditions> hVdrift = eventSetup.getHandle(vDriftMapToken_);
0045     vDriftMap_ = &*hVdrift;
0046     mTimeMap = nullptr;
0047     // Consistency check: no parametrization is implemented for the time being
0048     int version = vDriftMap_->version();
0049     if (version != 1) {
0050       throw cms::Exception("Configuration") << "only version 1 is presently supported for VDriftDB";
0051     }
0052   }
0053 }
0054
0055 void DTVDriftAnalyzer::endJob() {
0056   // Loop over DB entries
0057
0058   map<uint32_t, pair<float, float>> values;
0059
0060   if (readLegacyVDriftDB) {
0061     for (DTMtime::const_iterator mtime = mTimeMap->begin(); mtime != mTimeMap->end(); ++mtime) {
0062       DTWireId wireId(
0063           (*mtime).first.wheelId, (*mtime).first.stationId, (*mtime).first.sectorId, (*mtime).first.slId, 0, 0);
0064       float vdrift;
0065       float reso;
0066       DetId detId(wireId.rawId());
0067       // vdrift is cm/ns , resolution is cm
0068       mTimeMap->get(detId, vdrift, reso, DTVelocityUnits::cm_per_ns);
0069       values[wireId.rawId()] = make_pair(vdrift, reso);
0070     }
0071   } else {
0072     for (DTRecoConditions::const_iterator vd = vDriftMap_->begin(); vd != vDriftMap_->end(); ++vd) {
0073       DTWireId wireId(vd->first);
0074       float vdrift = vDriftMap_->get(wireId);
0075       values[vd->first] = make_pair(vdrift, 0.f);
0076     }
0077   }
0078
0079   for (map<uint32_t, pair<float, float>>::const_iterator it = values.begin(); it != values.end(); ++it) {
0080     float vdrift = it->second.first;
0081     float reso = it->second.second;
0082     DTWireId wireId(it->first);
0083     // vdrift is cm/ns , resolution is cm
0084     edm::LogVerbatim("DTVDriftAnalyzer") << "Wire: " << wireId << endl
0085                                          << " vdrift (cm/ns): " << vdrift << endl
0086                                          << " reso (cm): " << reso << endl;
0087
0088     //Define an histo for each wheel and each superlayer type
0089     TH1D* hVDriftHisto = theVDriftHistoMap[make_pair(wireId.wheel(), wireId.superlayer())];
0090     TH1D* hResoHisto = theResoHistoMap[make_pair(wireId.wheel(), wireId.superlayer())];
0091     if (hVDriftHisto == 0) {
0092       theFile->cd();
0093       TString name = getHistoName(wireId).c_str();
0094       if (wireId.superlayer() != 2) {
0095         hVDriftHisto = new TH1D(name + "_VDrift", "VDrift calibrated from MT per superlayer", 50, 0, 50);
0096         hResoHisto = new TH1D(name + "_Reso", "Reso calibrated from MT per superlayer", 50, 0, 50);
0097       } else {
0098         hVDriftHisto = new TH1D(name + "_VDrift", "VDrift calibrated from MT per superlayer", 36, 0, 36);
0099         hResoHisto = new TH1D(name + "_Reso", "Reso calibrated from MT per superlayer", 36, 0, 36);
0100       }
0101       theVDriftHistoMap[make_pair(wireId.wheel(), wireId.superlayer())] = hVDriftHisto;
0102       theResoHistoMap[make_pair(wireId.wheel(), wireId.superlayer())] = hResoHisto;
0103     }
0104
0105     //Fill the histos and set the bin label
0106     int binNumber = wireId.sector() + 12 * (wireId.station() - 1);
0107     hVDriftHisto->SetBinContent(binNumber, vdrift);
0108     hResoHisto->SetBinContent(binNumber, reso);
0109     string labelName;
0110     stringstream theStream;
0111     if (wireId.sector() == 1)
0112       theStream << "MB" << wireId.station() << "_Sec" << wireId.sector();
0113     else
0114       theStream << "Sec" << wireId.sector();
0115     theStream >> labelName;
0116     hVDriftHisto->GetXaxis()->SetBinLabel(binNumber, labelName.c_str());
0117     hResoHisto->GetXaxis()->SetBinLabel(binNumber, labelName.c_str());
0118
0119     //Define a distribution for each wheel,station and each superlayer type
0120     vector<int> Wh_St_SL;
0121     Wh_St_SL.push_back(wireId.wheel());
0122     Wh_St_SL.push_back(wireId.station());
0123     Wh_St_SL.push_back(wireId.superlayer());
0124     TH1D* hVDriftDistrib = theVDriftDistribMap[Wh_St_SL];
0125     TH1D* hResoDistrib = theResoDistribMap[Wh_St_SL];
0126     if (hVDriftDistrib == 0) {
0127       theFile->cd();
0128       TString name = getDistribName(wireId).c_str();
0129       hVDriftDistrib = new TH1D(name + "_VDrift", "VDrift calibrated from MT per superlayer", 100, 0.00530, 0.00580);
0130       hResoDistrib = new TH1D(name + "_Reso", "Reso calibrated from MT per superlayer", 300, 0.015, 0.045);
0131       theVDriftDistribMap[Wh_St_SL] = hVDriftDistrib;
0132       theResoDistribMap[Wh_St_SL] = hResoDistrib;
0133     }
0134     //Fill the distributions
0135     hVDriftDistrib->Fill(vdrift);
0136     hResoDistrib->Fill(reso);
0137   }
0138
0139   //Write histos in a .root file
0140   theFile->cd();
0141   for (map<pair<int, int>, TH1D*>::const_iterator lHisto = theVDriftHistoMap.begin(); lHisto != theVDriftHistoMap.end();
0142        ++lHisto) {
0143     (*lHisto).second->GetXaxis()->LabelsOption("v");
0144     (*lHisto).second->Write();
0145   }
0146   for (map<pair<int, int>, TH1D*>::const_iterator lHisto = theResoHistoMap.begin(); lHisto != theResoHistoMap.end();
0147        ++lHisto) {
0148     (*lHisto).second->GetXaxis()->LabelsOption("v");
0149     (*lHisto).second->Write();
0150   }
0151   for (map<vector<int>, TH1D*>::const_iterator lDistrib = theVDriftDistribMap.begin();
0152        lDistrib != theVDriftDistribMap.end();
0153        ++lDistrib) {
0154     (*lDistrib).second->Write();
0155   }
0156   for (map<vector<int>, TH1D*>::const_iterator lDistrib = theResoDistribMap.begin();
0157        lDistrib != theResoDistribMap.end();
0158        ++lDistrib) {
0159     (*lDistrib).second->Write();
0160   }
0161 }
0162
0163 string DTVDriftAnalyzer::getHistoName(const DTWireId& wId) const {
0164   string histoName;
0165   stringstream theStream;
0166   theStream << "Wheel" << wId.wheel() << "_SL" << wId.superlayer();
0167   theStream >> histoName;
0168   return histoName;
0169 }
0170
0171 string DTVDriftAnalyzer::getDistribName(const DTWireId& wId) const {
0172   string histoName;
0173   stringstream theStream;
0174   theStream << "Wheel" << wId.wheel() << "_Station" << wId.station() << "_SL" << wId.superlayer();
0175   theStream >> histoName;
0176   return histoName;
0177 }