00001
00002
00003
00004
00005
00006
00007
00008
00009
00010 #include <iostream>
00011 #include <fstream>
00012 #include <string>
00013 #include <cmath>
00014
00015 #include "FWCore/ServiceRegistry/interface/Service.h"
00016 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00017 #include "FWCore/Framework/interface/ESHandle.h"
00018
00019 #include "DQMServices/Core/interface/MonitorElement.h"
00020
00021 #include "DQMServices/Core/interface/DQMStore.h"
00022
00023 #include "DataFormats/EcalRawData/interface/EcalRawDataCollections.h"
00024 #include "DataFormats/EcalDetId/interface/EBDetId.h"
00025 #include "DataFormats/EcalRecHit/interface/EcalUncalibratedRecHit.h"
00026 #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
00027 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerEvmReadoutRecord.h"
00028
00029 #include "RecoLocalCalo/EcalRecAlgos/interface/EcalSeverityLevelAlgo.h"
00030 #include "RecoLocalCalo/EcalRecAlgos/interface/EcalSeverityLevelAlgoRcd.h"
00031
00032 #include "DQM/EcalCommon/interface/Numbers.h"
00033
00034 #include "DQM/EcalBarrelMonitorTasks/interface/EBTimingTask.h"
00035
00036 EBTimingTask::EBTimingTask(const edm::ParameterSet& ps){
00037
00038 init_ = false;
00039
00040 dqmStore_ = edm::Service<DQMStore>().operator->();
00041
00042 prefixME_ = ps.getUntrackedParameter<std::string>("prefixME", "");
00043
00044 enableCleanup_ = ps.getUntrackedParameter<bool>("enableCleanup", false);
00045
00046 mergeRuns_ = ps.getUntrackedParameter<bool>("mergeRuns", false);
00047
00048 energyThreshold_ = ps.getUntrackedParameter<double>("energyTreshold",1.0);
00049
00050 EcalRawDataCollection_ = ps.getParameter<edm::InputTag>("EcalRawDataCollection");
00051 EcalRecHitCollection_ = ps.getParameter<edm::InputTag>("EcalRecHitCollection");
00052
00053 L1GtEvmReadoutRecord_ = ps.getParameter<edm::InputTag>("L1GtEvmReadoutRecord");
00054
00055 useBeamStatus_ = ps.getUntrackedParameter<bool>("useBeamStatus", false);
00056
00057 for (int i = 0; i < 36; i++) {
00058 meTime_[i] = 0;
00059 meTimeMap_[i] = 0;
00060 meTimeAmpli_[i] = 0;
00061 }
00062
00063 meTimeAmpliSummary_ = 0;
00064 meTimeSummary1D_ = 0;
00065 meTimeSummaryMap_ = 0;
00066
00067 stableBeamsDeclared_ = false;
00068
00069 }
00070
00071 EBTimingTask::~EBTimingTask(){
00072
00073 }
00074
00075 void EBTimingTask::beginJob(void){
00076
00077 ievt_ = 0;
00078
00079 if ( dqmStore_ ) {
00080 dqmStore_->setCurrentFolder(prefixME_ + "/EBTimingTask");
00081 dqmStore_->rmdir(prefixME_ + "/EBTimingTask");
00082 }
00083
00084 }
00085
00086 void EBTimingTask::beginRun(const edm::Run& r, const edm::EventSetup& c) {
00087
00088 Numbers::initGeometry(c, false);
00089
00090 if ( ! mergeRuns_ ) this->reset();
00091
00092 stableBeamsDeclared_ = false;
00093
00094 }
00095
00096 void EBTimingTask::endRun(const edm::Run& r, const edm::EventSetup& c) {
00097
00098 }
00099
00100 void EBTimingTask::reset(void) {
00101
00102 for (int i = 0; i < 36; i++) {
00103 if ( meTime_[i] ) meTime_[i]->Reset();
00104 if ( meTimeMap_[i] ) meTimeMap_[i]->Reset();
00105 if ( meTimeAmpli_[i] ) meTimeAmpli_[i]->Reset();
00106 }
00107
00108 if ( meTimeAmpliSummary_ ) meTimeAmpliSummary_->Reset();
00109 if ( meTimeSummary1D_ ) meTimeSummary1D_->Reset();
00110 if ( meTimeSummaryMap_ ) meTimeSummaryMap_->Reset();
00111
00112 }
00113
00114 void EBTimingTask::setup(void){
00115
00116 init_ = true;
00117
00118 std::string name;
00119
00120
00121 const int nbinsE = 25;
00122 const float minlogE = -0.5;
00123 const float maxlogE = 2.;
00124 float binEdgesE[nbinsE + 1];
00125 for(int i = 0; i <= nbinsE; i++)
00126 binEdgesE[i] = std::pow((float)10., minlogE + (maxlogE - minlogE) / nbinsE * i);
00127
00128 const int nbinsT = 200;
00129 const float minT = -50.;
00130 const float maxT = 50.;
00131 float binEdgesT[nbinsT + 1];
00132 for(int i = 0; i <= nbinsT; i++)
00133 binEdgesT[i] = minT + (maxT - minT) / nbinsT * i;
00134
00135 if ( dqmStore_ ) {
00136 dqmStore_->setCurrentFolder(prefixME_ + "/EBTimingTask");
00137
00138 for (int i = 0; i < 36; i++) {
00139 name = "EBTMT timing 1D " + Numbers::sEB(i+1);
00140 meTime_[i] = dqmStore_->book1D(name, name, 50, -25., 25.);
00141 meTime_[i]->setAxisTitle("time (ns)", 1);
00142 dqmStore_->tag(meTime_[i], i+1);
00143
00144 name = "EBTMT timing " + Numbers::sEB(i+1);
00145 meTimeMap_[i] = dqmStore_->bookProfile2D(name, name, 85, 0., 85., 20, 0., 20., -20.+shiftProf2D_, 20.+shiftProf2D_, "s");
00146 meTimeMap_[i]->setAxisTitle("ieta", 1);
00147 meTimeMap_[i]->setAxisTitle("iphi", 2);
00148 meTimeMap_[i]->setAxisTitle("time (ns)", 3);
00149 dqmStore_->tag(meTimeMap_[i], i+1);
00150
00151 name = "EBTMT timing vs amplitude " + Numbers::sEB(i+1);
00152 meTimeAmpli_[i] = dqmStore_->book2D(name, name, nbinsE, binEdgesE, nbinsT, binEdgesT);
00153 meTimeAmpli_[i]->setAxisTitle("energy (GeV)", 1);
00154 meTimeAmpli_[i]->setAxisTitle("time (ns)", 2);
00155 dqmStore_->tag(meTimeAmpli_[i], i+1);
00156 }
00157
00158 name = "EBTMT timing vs amplitude summary";
00159 meTimeAmpliSummary_ = dqmStore_->book2D(name, name, nbinsE, binEdgesE, nbinsT, binEdgesT);
00160 meTimeAmpliSummary_->setAxisTitle("energy (GeV)", 1);
00161 meTimeAmpliSummary_->setAxisTitle("time (ns)", 2);
00162
00163 name = "EBTMT timing 1D summary";
00164 meTimeSummary1D_ = dqmStore_->book1D(name, name, 50, -25., 25.);
00165 meTimeSummary1D_->setAxisTitle("time (ns)", 1);
00166
00167 name = "EBTMT timing map";
00168 meTimeSummaryMap_ = dqmStore_->bookProfile2D(name, name, 72, 0., 360., 34, -85, 85, -20.+shiftProf2D_, 20.+shiftProf2D_, "s");
00169 meTimeSummaryMap_->setAxisTitle("jphi", 1);
00170 meTimeSummaryMap_->setAxisTitle("jeta", 2);
00171 meTimeSummaryMap_->setAxisTitle("time (ns)", 3);
00172
00173 }
00174
00175 }
00176
00177 void EBTimingTask::cleanup(void){
00178
00179 if ( ! init_ ) return;
00180
00181 if ( dqmStore_ ) {
00182 dqmStore_->setCurrentFolder(prefixME_ + "/EBTimingTask");
00183
00184 for ( int i = 0; i < 36; i++ ) {
00185 if ( meTime_[i] ) dqmStore_->removeElement( meTime_[i]->getName() );
00186 meTime_[i] = 0;
00187
00188 if ( meTimeMap_[i] ) dqmStore_->removeElement( meTimeMap_[i]->getName() );
00189 meTimeMap_[i] = 0;
00190
00191 if ( meTimeAmpli_[i] ) dqmStore_->removeElement( meTimeAmpli_[i]->getName() );
00192 meTimeAmpli_[i] = 0;
00193 }
00194
00195 if ( meTimeAmpliSummary_ ) dqmStore_->removeElement( meTimeAmpliSummary_->getName() );
00196 meTimeAmpliSummary_ = 0;
00197
00198 if ( meTimeSummary1D_ ) dqmStore_->removeElement( meTimeSummary1D_->getName() );
00199 meTimeSummary1D_ = 0;
00200
00201 if ( meTimeSummaryMap_ ) dqmStore_->removeElement( meTimeSummaryMap_->getName() );
00202 meTimeSummaryMap_ = 0;
00203
00204 }
00205
00206 init_ = false;
00207
00208 }
00209
00210 void EBTimingTask::endJob(void){
00211
00212 edm::LogInfo("EBTimingTask") << "analyzed " << ievt_ << " events";
00213
00214 if ( enableCleanup_ ) this->cleanup();
00215
00216 }
00217
00218 void EBTimingTask::analyze(const edm::Event& e, const edm::EventSetup& c){
00219
00220 const unsigned STABLE_BEAMS = 11;
00221
00222 bool isData = true;
00223 bool enable = false;
00224 int runType[36];
00225 for (int i=0; i<36; i++) runType[i] = -1;
00226
00227 edm::Handle<EcalRawDataCollection> dcchs;
00228
00229 if ( e.getByLabel(EcalRawDataCollection_, dcchs) ) {
00230
00231 for ( EcalRawDataCollection::const_iterator dcchItr = dcchs->begin(); dcchItr != dcchs->end(); ++dcchItr ) {
00232
00233 if ( Numbers::subDet( *dcchItr ) != EcalBarrel ) continue;
00234
00235 int ism = Numbers::iSM( *dcchItr, EcalBarrel );
00236
00237 runType[ism-1] = dcchItr->getRunType();
00238
00239 if ( dcchItr->getRunType() == EcalDCCHeaderBlock::COSMIC ||
00240 dcchItr->getRunType() == EcalDCCHeaderBlock::MTCC ||
00241 dcchItr->getRunType() == EcalDCCHeaderBlock::COSMICS_GLOBAL ||
00242 dcchItr->getRunType() == EcalDCCHeaderBlock::PHYSICS_GLOBAL ||
00243 dcchItr->getRunType() == EcalDCCHeaderBlock::COSMICS_LOCAL ||
00244 dcchItr->getRunType() == EcalDCCHeaderBlock::PHYSICS_LOCAL ) enable = true;
00245
00246 }
00247
00248 } else {
00249
00250 isData = false; enable = true;
00251 edm::LogWarning("EBTimingTask") << EcalRawDataCollection_ << " not available";
00252
00253 }
00254
00255 if ( ! enable ) return;
00256
00257 if ( ! init_ ) this->setup();
00258
00259 ievt_++;
00260
00261
00262 if( useBeamStatus_ && !stableBeamsDeclared_ ) {
00263 edm::Handle<L1GlobalTriggerEvmReadoutRecord> gtRecord;
00264 if( e.getByLabel(L1GtEvmReadoutRecord_, gtRecord) ) {
00265
00266 unsigned lhcBeamMode = gtRecord->gtfeWord().beamMode();
00267
00268 if( lhcBeamMode == STABLE_BEAMS ){
00269
00270 reset();
00271
00272 stableBeamsDeclared_ = true;
00273
00274 }
00275 }
00276 }
00277
00278 edm::ESHandle<EcalSeverityLevelAlgo> sevlv;
00279 c.get<EcalSeverityLevelAlgoRcd>().get(sevlv);
00280
00281 edm::Handle<EcalRecHitCollection> hits;
00282
00283 if ( e.getByLabel(EcalRecHitCollection_, hits) ) {
00284
00285 int neh = hits->size();
00286 LogDebug("EBTimingTask") << "event " << ievt_ << " hits collection size " << neh;
00287
00288 for ( EcalRecHitCollection::const_iterator hitItr = hits->begin(); hitItr != hits->end(); ++hitItr ) {
00289
00290 EBDetId id = hitItr->id();
00291
00292 int ic = id.ic();
00293 int ie = (ic-1)/20 + 1;
00294 int ip = (ic-1)%20 + 1;
00295
00296 int ism = Numbers::iSM( id );
00297
00298 float xie = ie - 0.5;
00299 float xip = ip - 0.5;
00300
00301 if ( isData ) {
00302
00303 if ( ! ( runType[ism-1] == EcalDCCHeaderBlock::COSMIC ||
00304 runType[ism-1] == EcalDCCHeaderBlock::MTCC ||
00305 runType[ism-1] == EcalDCCHeaderBlock::COSMICS_GLOBAL ||
00306 runType[ism-1] == EcalDCCHeaderBlock::PHYSICS_GLOBAL ||
00307 runType[ism-1] == EcalDCCHeaderBlock::COSMICS_LOCAL ||
00308 runType[ism-1] == EcalDCCHeaderBlock::PHYSICS_LOCAL ) ) continue;
00309
00310 }
00311
00312 MonitorElement* meTime = 0;
00313 MonitorElement* meTimeMap = 0;
00314 MonitorElement* meTimeAmpli = 0;
00315
00316 meTime = meTime_[ism-1];
00317 meTimeMap = meTimeMap_[ism-1];
00318 meTimeAmpli = meTimeAmpli_[ism-1];
00319
00320 float xval = hitItr->energy();
00321 float yval = hitItr->time();
00322
00323 uint32_t flag = hitItr->recoFlag();
00324
00325 uint32_t sev = sevlv->severityLevel(id, *hits);
00326
00327 if ( (flag == EcalRecHit::kGood || flag == EcalRecHit::kOutOfTime) && sev != EcalSeverityLevel::kWeird ) {
00328 if ( meTimeAmpli ) meTimeAmpli->Fill(xval, yval);
00329 if ( meTimeAmpliSummary_ ) meTimeAmpliSummary_->Fill(xval, yval);
00330
00331 if ( xval > energyThreshold_ ) {
00332 if ( meTime ) meTime->Fill(yval);
00333 if ( meTimeMap ) meTimeMap->Fill(xie, xip, yval+shiftProf2D_);
00334 if ( meTimeSummary1D_ ) meTimeSummary1D_->Fill(yval);
00335
00336 float xebeta = id.ieta() - 0.5 * id.zside();
00337 float xebphi = id.iphi() - 0.5;
00338 if ( meTimeSummaryMap_ ) meTimeSummaryMap_->Fill(xebphi, xebeta, yval+shiftProf2D_);
00339 }
00340
00341 }
00342 }
00343
00344 } else {
00345
00346 edm::LogWarning("EBTimingTask") << EcalRecHitCollection_ << " not available";
00347
00348 }
00349
00350 }
00351