File indexing completed on 2024-04-06 11:58:30
0001
0002
0003
0004
0005
0006
0007
0008 #include "DTT0Analyzer.h"
0009
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011 #include "FWCore/Framework/interface/EventSetup.h"
0012 #include "Geometry/Records/interface/MuonGeometryRecord.h"
0013 #include "CondFormats/DTObjects/interface/DTT0.h"
0014 #include "CondFormats/DataRecord/interface/DTT0Rcd.h"
0015 #include <iostream>
0016 #include "TFile.h"
0017 #include "TH1D.h"
0018 #include "TString.h"
0019
0020 using namespace edm;
0021 using namespace std;
0022
0023 DTT0Analyzer::DTT0Analyzer(const ParameterSet& pset) {
0024
0025 string rootFileName = pset.getUntrackedParameter<string>("rootFileName");
0026 theFile = new TFile(rootFileName.c_str(), "RECREATE");
0027 theFile->cd();
0028 t0Token_ = esConsumes<edm::Transition::BeginRun>();
0029 dtGeomToken_ = esConsumes<edm::Transition::BeginRun>();
0030 }
0031
0032 DTT0Analyzer::~DTT0Analyzer() { theFile->Close(); }
0033
0034 void DTT0Analyzer::beginRun(const edm::Run&, const edm::EventSetup& eventSetup) {
0035
0036 ESHandle<DTT0> t0 = eventSetup.getHandle(t0Token_);
0037 tZeroMap = &*t0;
0038
0039
0040 dtGeom = eventSetup.getHandle(dtGeomToken_);
0041 }
0042
0043 void DTT0Analyzer::endJob() {
0044
0045 for (DTT0::const_iterator tzero = tZeroMap->begin(); tzero != tZeroMap->end(); ++tzero) {
0046
0047
0048
0049
0050
0051
0052
0053 int channelId = tzero->channelId;
0054 if (channelId == 0)
0055 continue;
0056 DTWireId wireId(channelId);
0057
0058 float t0mean;
0059 float t0rms;
0060
0061 tZeroMap->get(wireId, t0mean, t0rms, DTTimeUnits::counts);
0062 cout << "Wire: " << wireId << endl
0063 << " T0 mean (TDC counts): " << t0mean << " T0_rms (TDC counts): " << t0rms << endl;
0064
0065 DTLayerId layerId = wireId.layerId();
0066 const int nWires = dtGeom->layer(layerId)->specificTopology().channels();
0067
0068
0069 TH1D* hT0Histo = theMeanHistoMap[layerId];
0070 if (hT0Histo == 0) {
0071 theFile->cd();
0072 TString name = getHistoName(layerId).c_str();
0073 hT0Histo = new TH1D(name + "_t0Mean",
0074 "mean T0 from pulses by Channel",
0075 nWires,
0076 dtGeom->layer(layerId)->specificTopology().firstChannel(),
0077 dtGeom->layer(layerId)->specificTopology().firstChannel() + nWires);
0078 theMeanHistoMap[layerId] = hT0Histo;
0079 }
0080
0081 TH1D* hSigmaT0Histo = theSigmaHistoMap[layerId];
0082 if (hSigmaT0Histo == 0) {
0083 theFile->cd();
0084 TString name = getHistoName(layerId).c_str();
0085 hSigmaT0Histo = new TH1D(name + "_t0Sigma",
0086 "sigma T0 from pulses by Channel",
0087 nWires,
0088 dtGeom->layer(layerId)->specificTopology().firstChannel(),
0089 dtGeom->layer(layerId)->specificTopology().firstChannel() + nWires);
0090 theSigmaHistoMap[layerId] = hSigmaT0Histo;
0091 }
0092
0093
0094 int nBin = hT0Histo->GetXaxis()->FindFixBin(wireId.wire());
0095 hT0Histo->SetBinContent(nBin, t0mean);
0096 hSigmaT0Histo->SetBinContent(nBin, t0rms);
0097 }
0098
0099
0100 theFile->cd();
0101 for (map<DTLayerId, TH1D*>::const_iterator lHisto = theMeanHistoMap.begin(); lHisto != theMeanHistoMap.end();
0102 ++lHisto) {
0103 (*lHisto).second->Write();
0104 }
0105
0106 for (map<DTLayerId, TH1D*>::const_iterator lHisto = theSigmaHistoMap.begin(); lHisto != theSigmaHistoMap.end();
0107 ++lHisto) {
0108 (*lHisto).second->Write();
0109 }
0110 }
0111
0112 string DTT0Analyzer::getHistoName(const DTLayerId& lId) const {
0113 string histoName;
0114 stringstream theStream;
0115 theStream << "Ch_" << lId.wheel() << "_" << lId.station() << "_" << lId.sector() << "_SL" << lId.superlayer() << "_L"
0116 << lId.layer();
0117 theStream >> histoName;
0118 return histoName;
0119 }