00001
00002
00003
00004
00005
00006
00007
00008
00009
00010 #include <iostream>
00011 #include <fstream>
00012
00013 #include "FWCore/ServiceRegistry/interface/Service.h"
00014 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00015 #include "FWCore/Framework/interface/ESHandle.h"
00016
00017 #include "DQMServices/Core/interface/MonitorElement.h"
00018
00019 #include "DQMServices/Core/interface/DQMStore.h"
00020
00021 #include "DataFormats/EcalRawData/interface/EcalRawDataCollections.h"
00022 #include "DataFormats/EcalDetId/interface/EEDetId.h"
00023 #include "DataFormats/EcalRecHit/interface/EcalUncalibratedRecHit.h"
00024 #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
00025
00026 #include "CondFormats/EcalObjects/interface/EcalChannelStatus.h"
00027 #include "CondFormats/DataRecord/interface/EcalChannelStatusRcd.h"
00028 #include "RecoLocalCalo/EcalRecAlgos/interface/EcalSeverityLevelAlgo.h"
00029
00030 #include "DQM/EcalCommon/interface/Numbers.h"
00031
00032 #include "DQM/EcalEndcapMonitorTasks/interface/EETimingTask.h"
00033
00034 EETimingTask::EETimingTask(const edm::ParameterSet& ps){
00035
00036 init_ = false;
00037
00038 initCaloGeometry_ = 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 EcalRawDataCollection_ = ps.getParameter<edm::InputTag>("EcalRawDataCollection");
00049 EcalRecHitCollection_ = ps.getParameter<edm::InputTag>("EcalRecHitCollection");
00050
00051 for (int i = 0; i < 18; i++) {
00052 meTime_[i] = 0;
00053 meTimeMap_[i] = 0;
00054 meTimeAmpli_[i] = 0;
00055 }
00056
00057 for (int i = 0; i < 2; i++) {
00058 meTimeAmpliSummary_[i] = 0;
00059 meTimeSummary1D_[i] = 0;
00060 meTimeSummaryMap_[i] = 0;
00061 meTimeSummaryMapProjEta_[i] = 0;
00062 meTimeSummaryMapProjPhi_[i] = 0;
00063 }
00064
00065 meTimeDelta_ = 0;
00066 meTimeDelta2D_ = 0;
00067
00068 }
00069
00070 EETimingTask::~EETimingTask(){
00071
00072 }
00073
00074 void EETimingTask::beginJob(void){
00075
00076 ievt_ = 0;
00077
00078 if ( dqmStore_ ) {
00079 dqmStore_->setCurrentFolder(prefixME_ + "/EETimingTask");
00080 dqmStore_->rmdir(prefixME_ + "/EETimingTask");
00081 }
00082
00083 }
00084
00085 void EETimingTask::beginRun(const edm::Run& r, const edm::EventSetup& c) {
00086
00087 Numbers::initGeometry(c, true);
00088
00089 if( !initCaloGeometry_ ) {
00090 c.get<CaloGeometryRecord>().get(pGeometry_);
00091 initCaloGeometry_ = true;
00092 }
00093
00094 if ( ! mergeRuns_ ) this->reset();
00095
00096 }
00097
00098 void EETimingTask::endRun(const edm::Run& r, const edm::EventSetup& c) {
00099
00100 }
00101
00102 void EETimingTask::reset(void) {
00103
00104 for (int i = 0; i < 18; i++) {
00105 if ( meTime_[i] ) meTime_[i]->Reset();
00106 if ( meTimeMap_[i] ) meTimeMap_[i]->Reset();
00107 if ( meTimeAmpli_[i] ) meTimeAmpli_[i]->Reset();
00108 }
00109
00110 for (int i = 0; i < 2; i++) {
00111 if ( meTimeAmpliSummary_[i] ) meTimeAmpliSummary_[i]->Reset();
00112 if ( meTimeSummary1D_[i] ) meTimeSummary1D_[i]->Reset();
00113 if ( meTimeSummaryMap_[i] ) meTimeSummaryMap_[i]->Reset();
00114 if ( meTimeSummaryMapProjEta_[i] ) meTimeSummaryMapProjEta_[i]->Reset();
00115 if ( meTimeSummaryMapProjPhi_[i] ) meTimeSummaryMapProjPhi_[i]->Reset();
00116 }
00117
00118 if ( meTimeDelta_ ) meTimeDelta_->Reset();
00119 if ( meTimeDelta2D_ ) meTimeDelta2D_->Reset();
00120
00121 }
00122
00123 void EETimingTask::setup(void){
00124
00125 init_ = true;
00126
00127 char histo[200];
00128
00129 if ( dqmStore_ ) {
00130 dqmStore_->setCurrentFolder(prefixME_ + "/EETimingTask");
00131
00132 for (int i = 0; i < 18; i++) {
00133 sprintf(histo, "EETMT timing 1D %s", Numbers::sEE(i+1).c_str());
00134 meTime_[i] = dqmStore_->book1D(histo, histo, 50, -50., 50.);
00135 meTime_[i]->setAxisTitle("time (ns)", 1);
00136 dqmStore_->tag(meTime_[i], i+1);
00137
00138 sprintf(histo, "EETMT timing %s", Numbers::sEE(i+1).c_str());
00139 meTimeMap_[i] = dqmStore_->bookProfile2D(histo, histo, 50, Numbers::ix0EE(i+1)+0., Numbers::ix0EE(i+1)+50., 50, Numbers::iy0EE(i+1)+0., Numbers::iy0EE(i+1)+50., 25., 75., "s");
00140 meTimeMap_[i]->setAxisTitle("ix", 1);
00141 if ( i+1 >= 1 && i+1 <= 9 ) meTimeMap_[i]->setAxisTitle("101-ix", 1);
00142 meTimeMap_[i]->setAxisTitle("iy", 2);
00143 meTimeMap_[i]->setAxisTitle("time (ns)", 3);
00144 dqmStore_->tag(meTimeMap_[i], i+1);
00145
00146 sprintf(histo, "EETMT timing vs amplitude %s", Numbers::sEE(i+1).c_str());
00147 meTimeAmpli_[i] = dqmStore_->book2D(histo, histo, 100, 0., 10., 50, -50., 50.);
00148 meTimeAmpli_[i]->setAxisTitle("energy (GeV)", 1);
00149 meTimeAmpli_[i]->setAxisTitle("time (ns)", 2);
00150 dqmStore_->tag(meTimeAmpli_[i], i+1);
00151 }
00152
00153 sprintf(histo, "EETMT timing vs amplitude summary EE -");
00154 meTimeAmpliSummary_[0] = dqmStore_->book2D(histo, histo, 100, 0., 10., 50, -50., 50.);
00155 meTimeAmpliSummary_[0]->setAxisTitle("energy (GeV)", 1);
00156 meTimeAmpliSummary_[0]->setAxisTitle("time (ns)", 2);
00157
00158 sprintf(histo, "EETMT timing vs amplitude summary EE +");
00159 meTimeAmpliSummary_[1] = dqmStore_->book2D(histo, histo, 100, 0., 10., 50, -50., 50.);
00160 meTimeAmpliSummary_[1]->setAxisTitle("energy (GeV)", 1);
00161 meTimeAmpliSummary_[1]->setAxisTitle("time (ns)", 2);
00162
00163 sprintf(histo, "EETMT timing 1D summary EE -");
00164 meTimeSummary1D_[0] = dqmStore_->book1D(histo, histo, 50, -50., 50.);
00165 meTimeSummary1D_[0]->setAxisTitle("time (ns)", 1);
00166
00167 sprintf(histo, "EETMT timing 1D summary EE +");
00168 meTimeSummary1D_[1] = dqmStore_->book1D(histo, histo, 50, -50., 50.);
00169 meTimeSummary1D_[1]->setAxisTitle("time (ns)", 1);
00170
00171 sprintf(histo, "EETMT timing map EE -");
00172 meTimeSummaryMap_[0] = dqmStore_->bookProfile2D(histo, histo, 20, 0., 100., 20, 0., 100., 50, 25., 75., "s");
00173 meTimeSummaryMap_[0]->setAxisTitle("ix'", 1);
00174 meTimeSummaryMap_[0]->setAxisTitle("101-iy'", 2);
00175 meTimeSummaryMap_[0]->setAxisTitle("time (ns)", 3);
00176
00177 sprintf(histo, "EETMT timing map EE +");
00178 meTimeSummaryMap_[1] = dqmStore_->bookProfile2D(histo, histo, 20, 0., 100., 20, 0., 100., 50, 25., 75., "s");
00179 meTimeSummaryMap_[1]->setAxisTitle("ix'", 1);
00180 meTimeSummaryMap_[1]->setAxisTitle("iy'", 2);
00181 meTimeSummaryMap_[1]->setAxisTitle("time (ns)", 3);
00182
00183 sprintf(histo, "EETMT timing projection eta EE -");
00184 meTimeSummaryMapProjEta_[0] = dqmStore_->bookProfile(histo, histo, 20, -3.0, -1.479, 50, -50., 50., "s");
00185 meTimeSummaryMapProjEta_[0]->setAxisTitle("eta", 1);
00186 meTimeSummaryMapProjEta_[0]->setAxisTitle("time (ns)", 2);
00187
00188 sprintf(histo, "EETMT timing projection eta EE +");
00189 meTimeSummaryMapProjEta_[1] = dqmStore_->bookProfile(histo, histo, 20, 1.479, 3.0, 50, -50., 50., "s");
00190 meTimeSummaryMapProjEta_[1]->setAxisTitle("phi", 1);
00191 meTimeSummaryMapProjEta_[1]->setAxisTitle("time (ns)", 2);
00192
00193 sprintf(histo, "EETMT timing projection phi EE -");
00194 meTimeSummaryMapProjPhi_[0] = dqmStore_->bookProfile(histo, histo, 50, -M_PI, M_PI, 50, -50., 50., "s");
00195 meTimeSummaryMapProjPhi_[0]->setAxisTitle("phi", 1);
00196 meTimeSummaryMapProjPhi_[0]->setAxisTitle("time (ns)", 2);
00197
00198 sprintf(histo, "EETMT timing projection phi EE +");
00199 meTimeSummaryMapProjPhi_[1] = dqmStore_->bookProfile(histo, histo, 50, -M_PI, M_PI, 50, -50., 50., "s");
00200 meTimeSummaryMapProjPhi_[1]->setAxisTitle("phi", 1);
00201 meTimeSummaryMapProjPhi_[1]->setAxisTitle("time (ns)", 2);
00202
00203 sprintf(histo, "EETMT timing EE+ - EE-");
00204 meTimeDelta_ = dqmStore_->book1D(histo, histo, 100, -3., 3.);
00205 meTimeDelta_->setAxisTitle("time (ns)", 1);
00206
00207 sprintf(histo, "EETMT timing EE+ vs EE-");
00208 meTimeDelta2D_ = dqmStore_->book2D(histo, histo, 50, -50., 50., 50, -50., 50.);
00209 meTimeDelta2D_->setAxisTitle("EE+ average time (ns)", 1);
00210 meTimeDelta2D_->setAxisTitle("EE- average time (ns)", 2);
00211
00212 }
00213
00214 }
00215
00216 void EETimingTask::cleanup(void){
00217
00218 if ( ! init_ ) return;
00219
00220 if ( dqmStore_ ) {
00221 dqmStore_->setCurrentFolder(prefixME_ + "/EETimingTask");
00222
00223 for ( int i = 0; i < 18; i++ ) {
00224 if ( meTime_[i] ) dqmStore_->removeElement( meTime_[i]->getName() );
00225 meTime_[i] = 0;
00226
00227 if ( meTimeMap_[i] ) dqmStore_->removeElement( meTimeMap_[i]->getName() );
00228 meTimeMap_[i] = 0;
00229
00230 if ( meTimeAmpli_[i] ) dqmStore_->removeElement( meTimeAmpli_[i]->getName() );
00231 meTimeAmpli_[i] = 0;
00232 }
00233
00234 for (int i = 0; i < 2; i++) {
00235 if ( meTimeAmpliSummary_[i] ) dqmStore_->removeElement( meTimeAmpliSummary_[i]->getName() );
00236 meTimeAmpliSummary_[i] = 0;
00237
00238 if ( meTimeSummary1D_[i] ) dqmStore_->removeElement( meTimeSummary1D_[i]->getName() );
00239 meTimeSummary1D_[i] = 0;
00240
00241 if ( meTimeSummaryMap_[i] ) dqmStore_->removeElement( meTimeSummaryMap_[i]->getName() );
00242 meTimeSummaryMap_[i] = 0;
00243
00244 if ( meTimeSummaryMapProjEta_[i] ) dqmStore_->removeElement( meTimeSummaryMapProjEta_[i]->getName() );
00245 meTimeSummaryMapProjEta_[i] = 0;
00246
00247 if ( meTimeSummaryMapProjPhi_[i] ) dqmStore_->removeElement( meTimeSummaryMapProjPhi_[i]->getName() );
00248 meTimeSummaryMapProjPhi_[i] = 0;
00249 }
00250
00251 if ( meTimeDelta_ ) dqmStore_->removeElement( meTimeDelta_->getName() );
00252 meTimeDelta_ = 0;
00253
00254 if ( meTimeDelta2D_ ) dqmStore_->removeElement( meTimeDelta2D_->getName() );
00255 meTimeDelta2D_ = 0;
00256
00257 }
00258
00259 init_ = false;
00260
00261 }
00262
00263 void EETimingTask::endJob(void){
00264
00265 edm::LogInfo("EETimingTask") << "analyzed " << ievt_ << " events";
00266
00267 if ( enableCleanup_ ) this->cleanup();
00268
00269 }
00270
00271 void EETimingTask::analyze(const edm::Event& e, const edm::EventSetup& c){
00272
00273 bool isData = true;
00274 bool enable = false;
00275 int runType[18];
00276 for (int i=0; i<18; i++) runType[i] = -1;
00277
00278 edm::Handle<EcalRawDataCollection> dcchs;
00279
00280 if ( e.getByLabel(EcalRawDataCollection_, dcchs) ) {
00281
00282 for ( EcalRawDataCollection::const_iterator dcchItr = dcchs->begin(); dcchItr != dcchs->end(); ++dcchItr ) {
00283
00284 if ( Numbers::subDet( *dcchItr ) != EcalEndcap ) continue;
00285
00286 int ism = Numbers::iSM( *dcchItr, EcalEndcap );
00287
00288 runType[ism-1] = dcchItr->getRunType();
00289
00290 if ( dcchItr->getRunType() == EcalDCCHeaderBlock::COSMIC ||
00291 dcchItr->getRunType() == EcalDCCHeaderBlock::MTCC ||
00292 dcchItr->getRunType() == EcalDCCHeaderBlock::COSMICS_GLOBAL ||
00293 dcchItr->getRunType() == EcalDCCHeaderBlock::PHYSICS_GLOBAL ||
00294 dcchItr->getRunType() == EcalDCCHeaderBlock::COSMICS_LOCAL ||
00295 dcchItr->getRunType() == EcalDCCHeaderBlock::PHYSICS_LOCAL ) enable = true;
00296
00297 }
00298
00299 } else {
00300
00301 isData = false; enable = true;
00302 edm::LogWarning("EETimingTask") << EcalRawDataCollection_ << " not available";
00303
00304 }
00305
00306 if ( ! enable ) return;
00307
00308 if ( ! init_ ) this->setup();
00309
00310 ievt_++;
00311
00312
00313 edm::ESHandle<EcalChannelStatus> pChannelStatus;
00314 c.get<EcalChannelStatusRcd>().get(pChannelStatus);
00315 const EcalChannelStatus* chStatus = pChannelStatus.product();
00316
00317 float sumTime_hithr[2] = {0.,0.};
00318 int n_hithr[2] = {0,0};
00319
00320 edm::Handle<EcalRecHitCollection> hits;
00321
00322 if ( e.getByLabel(EcalRecHitCollection_, hits) ) {
00323
00324 int neh = hits->size();
00325 LogDebug("EETimingTask") << "event " << ievt_ << " hits collection size " << neh;
00326
00327 for ( EcalRecHitCollection::const_iterator hitItr = hits->begin(); hitItr != hits->end(); ++hitItr ) {
00328
00329 EEDetId id = hitItr->id();
00330
00331 int ix = id.ix();
00332 int iy = id.iy();
00333 int iz = ( id.positiveZ() ) ? 1 : 0;
00334
00335 int ism = Numbers::iSM( id );
00336
00337 if ( ism >= 1 && ism <= 9 ) ix = 101 - ix;
00338
00339 float xix = ix - 0.5;
00340 float xiy = iy - 0.5;
00341
00342 if ( isData ) {
00343
00344 if ( ! ( runType[ism-1] == EcalDCCHeaderBlock::COSMIC ||
00345 runType[ism-1] == EcalDCCHeaderBlock::MTCC ||
00346 runType[ism-1] == EcalDCCHeaderBlock::COSMICS_GLOBAL ||
00347 runType[ism-1] == EcalDCCHeaderBlock::PHYSICS_GLOBAL ||
00348 runType[ism-1] == EcalDCCHeaderBlock::COSMICS_LOCAL ||
00349 runType[ism-1] == EcalDCCHeaderBlock::PHYSICS_LOCAL ) ) continue;
00350
00351 }
00352
00353 MonitorElement* meTime = 0;
00354 MonitorElement* meTimeMap = 0;
00355 MonitorElement* meTimeAmpli = 0;
00356
00357 meTime = meTime_[ism-1];
00358 meTimeMap = meTimeMap_[ism-1];
00359 meTimeAmpli = meTimeAmpli_[ism-1];
00360
00361 float xval = hitItr->energy();
00362 float yval = hitItr->time();
00363
00364 uint32_t flag = hitItr->recoFlag();
00365 uint32_t sev = EcalSeverityLevelAlgo::severityLevel(id, *hits, *chStatus );
00366
00367 const GlobalPoint& pos = pGeometry_->getGeometry(id)->getPosition();
00368
00369 float theta = pos.theta();
00370 float eta = pos.eta();
00371 float phi = pos.phi();
00372
00373 float et = hitItr->energy() * std::abs(sin(theta));
00374
00375 if ( (flag == EcalRecHit::kGood || flag == EcalRecHit::kOutOfTime) && sev != EcalSeverityLevelAlgo::kWeird ) {
00376 if ( meTimeAmpli ) meTimeAmpli->Fill(xval, yval);
00377 if ( meTimeAmpliSummary_[iz] ) meTimeAmpliSummary_[iz]->Fill(xval, yval);
00378 if ( et > 0.600 ) {
00379 if ( meTimeMap ) meTimeMap->Fill(xix, xiy, yval+50.);
00380 if ( meTime ) meTime->Fill(yval);
00381 if ( meTimeSummary1D_[iz] ) meTimeSummary1D_[iz]->Fill(yval);
00382
00383 if ( meTimeSummaryMap_[iz] ) meTimeSummaryMap_[iz]->Fill(xix, xiy, yval+50.);
00384 if ( meTimeSummaryMapProjEta_[iz] ) meTimeSummaryMapProjEta_[iz]->Fill(eta, yval);
00385 if ( meTimeSummaryMapProjPhi_[iz] ) meTimeSummaryMapProjPhi_[iz]->Fill(phi, yval);
00386
00387 sumTime_hithr[iz] += yval;
00388 n_hithr[iz]++;
00389 }
00390 }
00391 }
00392
00393 if (n_hithr[0] > 0 && n_hithr[1] > 0 ) {
00394 if ( meTimeDelta_ ) meTimeDelta_->Fill( sumTime_hithr[1]/n_hithr[1] - sumTime_hithr[0]/n_hithr[0] );
00395 if ( meTimeDelta2D_ ) meTimeDelta2D_->Fill( sumTime_hithr[1]/n_hithr[1], sumTime_hithr[0]/n_hithr[0] );
00396 }
00397
00398 } else {
00399
00400 edm::LogWarning("EETimingTask") << EcalRecHitCollection_ << " not available";
00401
00402 }
00403
00404 }
00405
1.7.3