Reference Manual

Go to the documentation of this file.
00001 /*
00002  * \file EEBeamHodoTask.cc
00003  *
00004  * $Date: 2010/08/08 08:46:09 $
00005  * $Revision: 1.30 $
00006  * \author G. Della Ricca
00007  * \author G. Franzoni
00008  *
00009  */
00010 
00011 #include <iostream>
00012 #include <fstream>
00013 
00014 #include "FWCore/ServiceRegistry/interface/Service.h"
00015 #include "FWCore/MessageLogger/interface/MessageLogger.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/EBDetId.h"
00023 #include "DataFormats/EcalRecHit/interface/EcalUncalibratedRecHit.h"
00024 #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
00025 
00026 #include "TBDataFormats/EcalTBObjects/interface/EcalTBCollections.h"
00027 
00028 #include "DQM/EcalCommon/interface/Numbers.h"
00029 
00030 #include "DQM/EcalEndcapMonitorTasks/interface/EEBeamHodoTask.h"
00031 
00032 EEBeamHodoTask::EEBeamHodoTask(const edm::ParameterSet& ps){
00033 
00034   init_ = false;
00035 
00036   dqmStore_ = edm::Service<DQMStore>().operator->();
00037 
00038   prefixME_ = ps.getUntrackedParameter<std::string>("prefixME", "");
00039 
00040   enableCleanup_ = ps.getUntrackedParameter<bool>("enableCleanup", false);
00041 
00042   mergeRuns_ = ps.getUntrackedParameter<bool>("mergeRuns", false);
00043 
00044   EcalTBEventHeader_ = ps.getParameter<edm::InputTag>("EcalTBEventHeader");
00045   EcalRawDataCollection_ = ps.getParameter<edm::InputTag>("EcalRawDataCollection");
00046   EcalUncalibratedRecHitCollection_ = ps.getParameter<edm::InputTag>("EcalUncalibratedRecHitCollection");
00047   EcalTBTDCRawInfo_ = ps.getParameter<edm::InputTag>("EcalTBTDCRawInfo");
00048   EcalTBHodoscopeRawInfo_ = ps.getParameter<edm::InputTag>("EcalTBHodoscopeRawInfo");
00049   EcalTBTDCRecInfo_ = ps.getParameter<edm::InputTag>("EcalTBTDCRecInfo");
00050   EcalTBHodoscopeRecInfo_ = ps.getParameter<edm::InputTag>("EcalTBHodoscopeRecInfo");
00051 
00052   tableIsMoving_ = false;
00053   cryInBeam_ =0;
00054   previousCryInBeam_ = -99999;
00055 
00056   // to be filled all the time
00057   for (int i=0; i<4; i++) {
00058     meHodoOcc_[i] =0;
00059     meHodoRaw_[i] =0;
00060   }
00061   meTDCRec_      =0;
00062 
00063   // filled only when: the table does not move
00064   meHodoPosRecX_    =0;
00065   meHodoPosRecY_    =0;
00066   meHodoPosRecXY_ =0;
00067   meHodoSloXRec_  =0;
00068   meHodoSloYRec_  =0;
00069   meHodoQuaXRec_ =0;
00070   meHodoQuaYRec_ =0;
00071   meHodoPosXMinusCaloPosXVsCry_   =0;
00072   meHodoPosYMinusCaloPosYVsCry_   =0;
00073   meTDCTimeMinusCaloTimeVsCry_ =0;
00074   meMissingCollections_        =0;
00075 
00076   meEvsXRecProf_     =0;
00077   meEvsYRecProf_     =0;
00078   meEvsXRecHis_     =0;
00079   meEvsYRecHis_     =0;
00080 
00081   //                       and matrix 5x5 available
00082   meCaloVsHodoXPos_ =0;
00083   meCaloVsHodoYPos_ =0;
00084   meCaloVsTDCTime_    =0;
00085 
00086 }
00087 
00088 EEBeamHodoTask::~EEBeamHodoTask(){
00089 
00090 }
00091 
00092 void EEBeamHodoTask::beginJob(void){
00093 
00094   ievt_  = 0;
00095 
00096   LV1_ = 0;
00097   cryInBeamCounter_ =0;
00098   resetNow_                =false;
00099 
00100   if ( dqmStore_ ) {
00101     dqmStore_->setCurrentFolder(prefixME_ + "/EEBeamHodoTask");
00102     dqmStore_->rmdir(prefixME_ + "/EEBeamHodoTask");
00103   }
00104 
00105 }
00106 
00107 void EEBeamHodoTask::beginRun(const edm::Run& r, const edm::EventSetup& c) {
00108 
00109   Numbers::initGeometry(c, false);
00110 
00111   if ( ! mergeRuns_ ) this->reset();
00112 
00113 }
00114 
00115 void EEBeamHodoTask::endRun(const edm::Run& r, const edm::EventSetup& c) {
00116 
00117 }
00118 
00119 void EEBeamHodoTask::reset(void) {
00120 
00121   for (int i=0; i<4; i++) {
00122     if ( meHodoOcc_[i] ) meHodoOcc_[i]->Reset();
00123     if ( meHodoRaw_[i] ) meHodoRaw_[i]->Reset();
00124   }
00125 
00126   if ( meHodoPosRecX_ ) meHodoPosRecX_->Reset();
00127   if ( meHodoPosRecY_ ) meHodoPosRecY_->Reset();
00128   if ( meHodoPosRecXY_ ) meHodoPosRecXY_->Reset();
00129   if ( meHodoSloXRec_ ) meHodoSloXRec_->Reset();
00130   if ( meHodoSloYRec_ ) meHodoSloYRec_->Reset();
00131   if ( meHodoQuaXRec_ ) meHodoQuaXRec_->Reset();
00132   if ( meHodoQuaYRec_ ) meHodoQuaYRec_->Reset();
00133   if ( meTDCRec_ ) meTDCRec_->Reset();
00134   if ( meEvsXRecProf_ ) meEvsXRecProf_->Reset();
00135   if ( meEvsYRecProf_ ) meEvsYRecProf_->Reset();
00136   if ( meEvsXRecHis_ ) meEvsXRecHis_->Reset();
00137   if ( meEvsYRecHis_ ) meEvsYRecHis_->Reset();
00138   if ( meCaloVsHodoXPos_ ) meCaloVsHodoXPos_->Reset();
00139   if ( meCaloVsHodoYPos_ ) meCaloVsHodoYPos_->Reset();
00140   if ( meCaloVsTDCTime_ ) meCaloVsTDCTime_->Reset();
00141   if ( meHodoPosXMinusCaloPosXVsCry_  ) meHodoPosXMinusCaloPosXVsCry_->Reset();
00142   if ( meHodoPosYMinusCaloPosYVsCry_  ) meHodoPosYMinusCaloPosYVsCry_->Reset();
00143   if ( meTDCTimeMinusCaloTimeVsCry_  ) meTDCTimeMinusCaloTimeVsCry_->Reset();
00144   if ( meMissingCollections_  ) meMissingCollections_->Reset();
00145 
00146 }
00147 
00148 void EEBeamHodoTask::setup(void){
00149 
00150   init_ = true;
00151 
00152   smId =1;
00153 
00154   char histo[200];
00155 
00156   if ( dqmStore_ ) {
00157     dqmStore_->setCurrentFolder(prefixME_ + "/EEBeamHodoTask");
00158 
00159     // following ME (type I):
00160     //  *** do not need to be ever reset
00161     //  *** can be filled regardless of the moving/notMoving status of the table
00162 
00163     for (int i=0; i<4; i++) {
00164       sprintf(histo, "EEBHT occup %s %02d", Numbers::sEE(smId).c_str(), i+1);
00165       meHodoOcc_[i] = dqmStore_->book1D(histo, histo, 30, 0., 30.);
00166       meHodoOcc_[i]->setAxisTitle("hits per event", 1);
00167       sprintf(histo, "EEBHT raw %s %02d", Numbers::sEE(smId).c_str(), i+1);
00168       meHodoRaw_[i] = dqmStore_->book1D(histo, histo, 64, 0., 64.);
00169       meHodoRaw_[i]->setAxisTitle("hodo fiber number", 1);
00170     }
00171 
00172     sprintf(histo, "EEBHT PosX rec %s", Numbers::sEE(smId).c_str());
00173     meHodoPosRecX_ = dqmStore_->book1D(histo, histo, 100, -20, 20);
00174     meHodoPosRecX_->setAxisTitle("reconstructed position    (mm)", 1);
00175 
00176     sprintf(histo, "EEBHT PosY rec %s", Numbers::sEE(smId).c_str());
00177     meHodoPosRecY_ = dqmStore_->book1D(histo, histo, 100, -20, 20);
00178     meHodoPosRecY_->setAxisTitle("reconstructed position    (mm)", 1);
00179 
00180     sprintf(histo, "EEBHT PosYX rec %s", Numbers::sEE(smId).c_str());
00181     meHodoPosRecXY_ = dqmStore_->book2D(histo, histo, 100, -20, 20,100, -20, 20);
00182     meHodoPosRecXY_->setAxisTitle("reconstructed X position    (mm)", 1);
00183     meHodoPosRecXY_->setAxisTitle("reconstructed Y position    (mm)", 2);
00184 
00185     sprintf(histo, "EEBHT SloX %s", Numbers::sEE(smId).c_str());
00186     meHodoSloXRec_ = dqmStore_->book1D(histo, histo, 50, -0.005, 0.005);
00187     meHodoSloXRec_->setAxisTitle("reconstructed track slope", 1);
00188 
00189     sprintf(histo, "EEBHT SloY %s", Numbers::sEE(smId).c_str());
00190     meHodoSloYRec_ = dqmStore_->book1D(histo, histo, 50, -0.005, 0.005);
00191     meHodoSloYRec_->setAxisTitle("reconstructed track slope", 1);
00192 
00193     sprintf(histo, "EEBHT QualX %s", Numbers::sEE(smId).c_str());
00194     meHodoQuaXRec_ = dqmStore_->book1D(histo, histo, 50, 0, 5);
00195     meHodoQuaXRec_->setAxisTitle("track fit quality", 1);
00196 
00197     sprintf(histo, "EEBHT QualY %s", Numbers::sEE(smId).c_str());
00198     meHodoQuaYRec_ = dqmStore_->book1D(histo, histo, 50, 0, 5);
00199     meHodoQuaYRec_->setAxisTitle("track fit quality", 1);
00200 
00201     sprintf(histo, "EEBHT TDC rec %s", Numbers::sEE(smId).c_str());
00202     meTDCRec_  = dqmStore_->book1D(histo, histo, 25, 0, 1);
00203     meTDCRec_->setAxisTitle("offset", 1);
00204 
00205     sprintf(histo, "EEBHT Hodo-Calo X vs Cry %s", Numbers::sEE(smId).c_str());
00206     meHodoPosXMinusCaloPosXVsCry_  = dqmStore_->book1D(histo, histo, 50, 0, 50);
00207     meHodoPosXMinusCaloPosXVsCry_->setAxisTitle("scan step number", 1);
00208     meHodoPosXMinusCaloPosXVsCry_->setAxisTitle("PosX_{hodo} - PosX_{calo}    (mm)", 2);
00209 
00210     sprintf(histo, "EEBHT Hodo-Calo Y vs Cry %s", Numbers::sEE(smId).c_str());
00211     meHodoPosYMinusCaloPosYVsCry_  = dqmStore_->book1D(histo, histo, 50, 0, 50);
00212     meHodoPosYMinusCaloPosYVsCry_->setAxisTitle("scan step number", 1);
00213     meHodoPosYMinusCaloPosYVsCry_->setAxisTitle("PosY_{hodo} - PosY_{calo}    (mm)", 2);
00214 
00215     sprintf(histo, "EEBHT TDC-Calo vs Cry %s", Numbers::sEE(smId).c_str());
00216     meTDCTimeMinusCaloTimeVsCry_  = dqmStore_->book1D(histo, histo, 50, 0, 50);
00217     meTDCTimeMinusCaloTimeVsCry_->setAxisTitle("scan step number", 1);
00218     meTDCTimeMinusCaloTimeVsCry_->setAxisTitle("Time_{TDC} - Time_{calo}    (sample)", 2);
00219 
00220     sprintf(histo, "EEBHT Missing Collections %s", Numbers::sEE(smId).c_str());
00221     meMissingCollections_ = dqmStore_->book1D(histo, histo, 7, 0, 7);
00222     meMissingCollections_->setAxisTitle("missing collection", 1);
00223 
00224     // following ME (type II):
00225     //  *** can be filled only when table is **not** moving
00226     //  *** need to be reset once table goes from 'moving'->notMoving
00227 
00228     sprintf(histo, "EEBHT prof E1 vs X %s", Numbers::sEE(smId).c_str());
00229     meEvsXRecProf_    = dqmStore_-> bookProfile(histo, histo, 100, -20, 20, 500, 0, 5000, "s");
00230     meEvsXRecProf_->setAxisTitle("PosX    (mm)", 1);
00231     meEvsXRecProf_->setAxisTitle("E1 (ADC)", 2);
00232 
00233     sprintf(histo, "EEBHT prof E1 vs Y %s", Numbers::sEE(smId).c_str());
00234     meEvsYRecProf_    = dqmStore_-> bookProfile(histo, histo, 100, -20, 20, 500, 0, 5000, "s");
00235     meEvsYRecProf_->setAxisTitle("PosY    (mm)", 1);
00236     meEvsYRecProf_->setAxisTitle("E1 (ADC)", 2);
00237 
00238     sprintf(histo, "EEBHT his E1 vs X %s", Numbers::sEE(smId).c_str());
00239     meEvsXRecHis_    = dqmStore_-> book2D(histo, histo, 100, -20, 20, 500, 0, 5000);
00240     meEvsXRecHis_->setAxisTitle("PosX    (mm)", 1);
00241     meEvsXRecHis_->setAxisTitle("E1 (ADC)", 2);
00242 
00243     sprintf(histo, "EEBHT his E1 vs Y %s", Numbers::sEE(smId).c_str());
00244     meEvsYRecHis_    = dqmStore_-> book2D(histo, histo, 100, -20, 20, 500, 0, 5000);
00245     meEvsYRecHis_->setAxisTitle("PosY    (mm)", 1);
00246     meEvsYRecHis_->setAxisTitle("E1 (ADC)", 2);
00247 
00248     sprintf(histo, "EEBHT PosX Hodo-Calo %s", Numbers::sEE(smId).c_str());
00249     meCaloVsHodoXPos_   = dqmStore_->book1D(histo, histo, 40, -20, 20);
00250     meCaloVsHodoXPos_->setAxisTitle("PosX_{hodo} - PosX_{calo}     (mm)", 1);
00251 
00252     sprintf(histo, "EEBHT PosY Hodo-Calo %s", Numbers::sEE(smId).c_str());
00253     meCaloVsHodoYPos_   = dqmStore_->book1D(histo, histo, 40, -20, 20);
00254     meCaloVsHodoYPos_->setAxisTitle("PosY_{hodo} - PosY_{calo}     (mm)", 1);
00255 
00256     sprintf(histo, "EEBHT TimeMax TDC-Calo %s", Numbers::sEE(smId).c_str());
00257     meCaloVsTDCTime_  = dqmStore_->book1D(histo, histo, 100, -1, 1);//tentative
00258     meCaloVsTDCTime_->setAxisTitle("Time_{TDC} - Time_{calo} (samples)", 1);
00259 
00260   }
00261 
00262 }
00263 
00264 void EEBeamHodoTask::cleanup(void){
00265 
00266   if ( ! init_ ) return;
00267 
00268   if ( dqmStore_ ) {
00269     dqmStore_->setCurrentFolder(prefixME_ + "/EEBeamHodoTask");
00270 
00271     for (int i=0; i<4; i++) {
00272       if ( meHodoOcc_[i] ) dqmStore_->removeElement( meHodoOcc_[i]->getName() );
00273       meHodoOcc_[i] = 0;
00274       if ( meHodoRaw_[i] ) dqmStore_->removeElement( meHodoRaw_[i]->getName() );
00275       meHodoRaw_[i] = 0;
00276     }
00277 
00278     if ( meHodoPosRecX_ ) dqmStore_->removeElement( meHodoPosRecX_->getName() );
00279     meHodoPosRecX_ = 0;
00280     if ( meHodoPosRecY_ ) dqmStore_->removeElement( meHodoPosRecY_->getName() );
00281     meHodoPosRecY_ = 0;
00282     if ( meHodoPosRecXY_ ) dqmStore_->removeElement( meHodoPosRecXY_->getName() );
00283     meHodoPosRecXY_ = 0;
00284     if ( meHodoSloXRec_ ) dqmStore_->removeElement( meHodoSloXRec_->getName() );
00285     meHodoSloXRec_ = 0;
00286     if ( meHodoSloYRec_ ) dqmStore_->removeElement( meHodoSloYRec_->getName() );
00287     meHodoSloYRec_ = 0;
00288     if ( meHodoQuaXRec_ ) dqmStore_->removeElement( meHodoQuaXRec_->getName() );
00289     meHodoQuaXRec_ = 0;
00290     if ( meHodoQuaYRec_ ) dqmStore_->removeElement( meHodoQuaYRec_->getName() );
00291     meHodoQuaYRec_ = 0;
00292     if ( meTDCRec_ ) dqmStore_->removeElement( meTDCRec_->getName() );
00293     meTDCRec_ = 0;
00294     if ( meEvsXRecProf_ ) dqmStore_->removeElement( meEvsXRecProf_->getName() );
00295     meEvsXRecProf_ = 0;
00296     if ( meEvsYRecProf_ ) dqmStore_->removeElement( meEvsYRecProf_->getName() );
00297     meEvsYRecProf_ = 0;
00298     if ( meEvsXRecHis_ ) dqmStore_->removeElement( meEvsXRecHis_->getName() );
00299     meEvsXRecHis_ = 0;
00300     if ( meEvsYRecHis_ ) dqmStore_->removeElement( meEvsYRecHis_->getName() );
00301     meEvsYRecHis_ = 0;
00302     if ( meCaloVsHodoXPos_ ) dqmStore_->removeElement( meCaloVsHodoXPos_->getName() );
00303     meCaloVsHodoXPos_ = 0;
00304     if ( meCaloVsHodoYPos_ ) dqmStore_->removeElement( meCaloVsHodoYPos_->getName() );
00305     meCaloVsHodoYPos_ = 0;
00306     if ( meCaloVsTDCTime_ ) dqmStore_->removeElement( meCaloVsTDCTime_->getName() );
00307     meCaloVsTDCTime_ = 0;
00308     if ( meHodoPosXMinusCaloPosXVsCry_  ) dqmStore_->removeElement( meHodoPosXMinusCaloPosXVsCry_->getName() );
00309     meHodoPosXMinusCaloPosXVsCry_  = 0;
00310     if ( meHodoPosYMinusCaloPosYVsCry_  ) dqmStore_->removeElement( meHodoPosYMinusCaloPosYVsCry_->getName() );
00311     meHodoPosYMinusCaloPosYVsCry_  = 0;
00312     if ( meTDCTimeMinusCaloTimeVsCry_  ) dqmStore_->removeElement( meTDCTimeMinusCaloTimeVsCry_->getName() );
00313     meTDCTimeMinusCaloTimeVsCry_  = 0;
00314     if ( meMissingCollections_  ) dqmStore_->removeElement( meMissingCollections_->getName() );
00315     meMissingCollections_  = 0;
00316 
00317   }
00318 
00319   init_ = false;
00320 
00321 }
00322 
00323 void EEBeamHodoTask::endJob(void){
00324 
00325   edm::LogInfo("EEBeamHodoTask") << "analyzed " << ievt_ << " events";
00326 
00327   if ( enableCleanup_ ) this->cleanup();
00328 
00329 }
00330 
00331 void EEBeamHodoTask::analyze(const edm::Event& e, const edm::EventSetup& c){
00332 
00333   bool enable = false;
00334   edm::Handle<EcalTBEventHeader> pHeader;
00335   const EcalTBEventHeader* Header =0;
00336 
00337   if ( e.getByLabel(EcalTBEventHeader_, pHeader) ) {
00338     Header = pHeader.product(); // get a ptr to the product
00339     if (!Header) {
00340       edm::LogWarning("EEBeamHodoTask") << "Event header not found. Returning. ";
00341       meMissingCollections_-> Fill(0); // bin1: missing CMSSW edm::Event header
00342       return;
00343     }
00344     tableIsMoving_     = Header->tableIsMoving();
00345     cryInBeam_           = Header->crystalInBeam();  //  cryInBeam_         = Header->nominalCrystalInBeam();
00346     if (previousCryInBeam_ == -99999 )
00347       {      previousCryInBeam_ = cryInBeam_ ;    }
00348 
00349     LogDebug("EEBeamHodoTask") << "event: " << ievt_ << " event header found ";
00350     if (tableIsMoving_){
00351       LogDebug("EEBeamHodoTask") << "Table is moving. ";  }
00352     else
00353       {      LogDebug("EEBeamHodoTask") << "Table is not moving. ";    }
00354   } else {
00355     edm::LogWarning("EEBeamHodoTask") << "Event header not found (exception caught). Returning. ";
00356     return;
00357   }
00358 
00359   edm::Handle<EcalRawDataCollection> dcchs;
00360 
00361   if ( e.getByLabel(EcalRawDataCollection_, dcchs) ) {
00362 
00363     for ( EcalRawDataCollection::const_iterator dcchItr = dcchs->begin(); dcchItr != dcchs->end(); ++dcchItr ) {
00364 
00365       if ( Numbers::subDet( *dcchItr ) != EcalEndcap ) continue;
00366 
00367       if ( dcchItr->getRunType() == EcalDCCHeaderBlock::BEAMH4 ||
00368            dcchItr->getRunType() == EcalDCCHeaderBlock::BEAMH2  ) enable = true;
00369     }
00370 
00371   } else {
00372     edm::LogWarning("EcalBeamTask") << EcalRawDataCollection_ << " not available";
00373     meMissingCollections_-> Fill(1); // bin2: missing DCC headers
00374     return;
00375     // see bottom of cc file for compatibility to 2004 data [***]
00376   }
00377 
00378   if ( ! enable ) return;
00379   if ( ! init_ ) this->setup();
00380   ievt_++;
00381 
00382   LV1_ = Header->eventNumber();
00383 
00384   edm::Handle<EcalUncalibratedRecHitCollection> pUncalRH;
00385   const EcalUncalibratedRecHitCollection* uncalRecH =0;
00386 
00387   if ( e.getByLabel(EcalUncalibratedRecHitCollection_, pUncalRH) ) {
00388     uncalRecH = pUncalRH.product(); // get a ptr to the product
00389     int neh = pUncalRH->size();
00390     LogDebug("EEBeamHodoTask") << EcalUncalibratedRecHitCollection_ << " found in event " << ievt_ << "; hits collection size " << neh;
00391   } else {
00392     edm::LogWarning("EEBeamHodoTask") << EcalUncalibratedRecHitCollection_ << " not available";
00393     meMissingCollections_-> Fill(2); // bin3: missing uncalibRecHits
00394     return;
00395   }
00396 
00397   edm::Handle<EcalTBTDCRawInfo> pTDCRaw;
00398   const EcalTBTDCRawInfo* rawTDC=0;
00399 
00400   if ( e.getByLabel(EcalTBTDCRawInfo_, pTDCRaw) ) {
00401     rawTDC = pTDCRaw.product();
00402   } else {
00403     edm::LogError("EcalBeamTask") << "Error! Can't get the product EcalTBTDCRawInfo. Returning.";
00404     meMissingCollections_-> Fill(4); // bin5: missing raw TDC
00405     return;
00406   }
00407 
00408   edm::Handle<EcalTBHodoscopeRawInfo> pHodoRaw;
00409   const EcalTBHodoscopeRawInfo* rawHodo=0;
00410 
00411   if ( e.getByLabel(EcalTBHodoscopeRawInfo_, pHodoRaw) ) {
00412     rawHodo = pHodoRaw.product();
00413     if(rawHodo->planes() ){
00414     LogDebug("EcalBeamTask") << "hodoscopeRaw:  num planes: " <<  rawHodo->planes()
00415                              << " channels in plane 1: "  <<  rawHodo->channels(0);
00416     }
00417   } else {
00418     edm::LogError("EcalBeamTask") << "Error! Can't get the product EcalTBHodoscopeRawInfo. Returning";
00419     meMissingCollections_-> Fill(3); // bin4: missing raw hodo hits collection
00420     return;
00421   }
00422 
00423   if ( !rawTDC ||!rawHodo || !uncalRecH  || !( rawHodo->planes() )) {
00424       edm::LogWarning("EcalBeamTask") << "analyze: missing a needed collection or hodo collection empty. Returning.";
00425       return;
00426   }
00427   LogDebug("EEBeamHodoTask") << " TDC raw, Hodo raw, uncalRecH and DCCheader found.";
00428 
00429 
00430 
00431   // table has come to a stop is identified by new value of cry_in_beam
00432   //   - increase counter of crystals that have been on beam
00433   //   - set flag for resetting
00434   if (cryInBeam_ != previousCryInBeam_ ) {
00435       previousCryInBeam_ = cryInBeam_ ;
00436       cryInBeamCounter_++;
00437       resetNow_ = true;
00438 
00439       // since flag "tableIsMoving==false" is reliable (as we can tell, so far),
00440       // operations due when "table has started moving"
00441       // can be done after the change in crystal in beam
00442 
00443       LogDebug("EcalBeamTask")  << "At event number : " << LV1_ <<  " switching table status: from moving to still. "
00444                                 << " cry in beam is: " << cryInBeam_ << ", step being: " << cryInBeamCounter_ ;
00445 
00446       // fill here plots which keep history of beamed crystals
00447       float HodoPosXMinusCaloPosXVsCry_mean  =0;
00448       float HodoPosXMinusCaloPosXVsCry_rms   =0;
00449       float HodoPosYMinusCaloPosYVsCry_mean  =0;
00450       float HodoPosYMinusCaloPosYVsCry_rms   =0;
00451       float TDCTimeMinusCaloTimeVsCry_mean   =0;
00452       float TDCTimeMinusCaloTimeVsCry_rms    =0;
00453 
00454       // min number of entries chosen assuming:
00455       //                 prescaling = 100 X 2FU
00456       //                 that we want at leas 2k events per crystal
00457 
00458       if (meCaloVsHodoXPos_-> getEntries()  > 10){
00459         HodoPosXMinusCaloPosXVsCry_mean = meCaloVsHodoXPos_-> getMean(1);
00460         HodoPosXMinusCaloPosXVsCry_rms  = meCaloVsHodoXPos_-> getRMS(1);
00461         meHodoPosXMinusCaloPosXVsCry_->  setBinContent( cryInBeamCounter_, HodoPosXMinusCaloPosXVsCry_mean);
00462         meHodoPosXMinusCaloPosXVsCry_->  setBinError(      cryInBeamCounter_, HodoPosXMinusCaloPosXVsCry_rms);
00463         LogDebug("EcalBeamTask")  << "At event number: " << LV1_ << " step: " << cryInBeamCounter_
00464                                   <<  " DeltaPosX is: " << (meCaloVsHodoXPos_-> getMean(1))
00465                                   << " +-" << ( meCaloVsHodoXPos_-> getRMS(1));
00466       }
00467       if (meCaloVsHodoYPos_-> getEntries()  > 10){
00468         HodoPosYMinusCaloPosYVsCry_mean = meCaloVsHodoYPos_-> getMean(1);
00469         HodoPosYMinusCaloPosYVsCry_rms  = meCaloVsHodoYPos_-> getRMS(1);
00470         meHodoPosYMinusCaloPosYVsCry_->  setBinContent( cryInBeamCounter_, HodoPosYMinusCaloPosYVsCry_mean);
00471         meHodoPosYMinusCaloPosYVsCry_->  setBinError(      cryInBeamCounter_, HodoPosYMinusCaloPosYVsCry_rms);
00472         LogDebug("EcalBeamTask")  << "At event number: " << LV1_ << " step: " << cryInBeamCounter_
00473                                   <<  " DeltaPosY is: " << (meCaloVsHodoYPos_-> getMean(1))
00474                                   << " +-" << ( meCaloVsHodoYPos_-> getRMS(1));
00475       }
00476       if (meCaloVsTDCTime_-> getEntries()  > 10){
00477         TDCTimeMinusCaloTimeVsCry_mean     = meCaloVsTDCTime_-> getMean(1);
00478         TDCTimeMinusCaloTimeVsCry_rms      = meCaloVsTDCTime_-> getRMS(1);
00479         meTDCTimeMinusCaloTimeVsCry_->  setBinContent(cryInBeamCounter_, TDCTimeMinusCaloTimeVsCry_mean);
00480         meTDCTimeMinusCaloTimeVsCry_->  setBinError(cryInBeamCounter_, TDCTimeMinusCaloTimeVsCry_rms);
00481         LogDebug("EcalBeamTask")  << "At event number: " << LV1_ << " step: " << cryInBeamCounter_
00482                                   <<  " DeltaT is: " << (meCaloVsTDCTime_-> getMean(1))
00483                                   << " +-" << ( meCaloVsTDCTime_-> getRMS(1));
00484       }
00485 
00486       LogDebug("EcalBeamTask")  << "At event number: " << LV1_ <<  " trace histos filled ( cryInBeamCounter_="
00487                                 << cryInBeamCounter_ << ")";
00488 
00489     }
00490 
00491 
00492 
00493 
00494   // if table has come to rest (from movement), reset concerned ME's
00495   if (resetNow_)
00496     {
00497       meEvsXRecProf_->Reset();
00498       meEvsYRecProf_->Reset();
00499       meEvsXRecHis_->Reset();
00500       meEvsYRecHis_->Reset();
00501       meCaloVsHodoXPos_->Reset();
00502       meCaloVsHodoYPos_->Reset();
00503       meCaloVsTDCTime_->Reset();
00504 
00505       resetNow_ = false;
00506     }
00507 
00508 
00509 
00510 
00511   /**************************************
00512   // handling histos type I:
00513   **************************************/
00514   for (unsigned int planeId=0; planeId <4; planeId++){
00515 
00516     const EcalTBHodoscopePlaneRawHits& planeRaw = rawHodo->getPlaneRawHits(planeId);
00517     LogDebug("EcalBeamTask")  << "\t plane: " << (planeId+1)
00518                               << "\t number of fibers: " << planeRaw.channels()
00519                               << "\t number of hits: " << planeRaw.numberOfFiredHits();
00520     meHodoOcc_[planeId]-> Fill( planeRaw.numberOfFiredHits() );
00521 
00522     for (unsigned int i=0;i<planeRaw.channels();i++)
00523       {
00524         if (planeRaw.isChannelFired(i))
00525           {
00526             edm::LogInfo("EcalBeamTask")<< " channel " << (i+1) << " has fired";
00527             meHodoRaw_[planeId]-> Fill(i+0.5);
00528           }
00529       }
00530   }
00531 
00532 
00533 
00534   edm::Handle<EcalTBTDCRecInfo> pTDC;
00535   const EcalTBTDCRecInfo* recTDC=0;
00536 
00537   if ( e.getByLabel(EcalTBTDCRecInfo_, pTDC) ) {
00538     recTDC = pTDC.product();
00539     LogDebug("EEBeamHodoTask") << " TDC offset is: " << recTDC->offset();
00540   } else {
00541     edm::LogError("EcalBeamTask") << "Error! Can't get the product EcalTBTDCRecInfo. Returning";
00542     meMissingCollections_-> Fill(5); // bin6: missing reconstructed TDC
00543     return;
00544   }
00545 
00546   edm::Handle<EcalTBHodoscopeRecInfo> pHodo;
00547   const EcalTBHodoscopeRecInfo* recHodo=0;
00548 
00549   if ( e.getByLabel(EcalTBHodoscopeRecInfo_, pHodo) ) {
00550     recHodo = pHodo.product();
00551     LogDebug("EcalBeamTask") << "hodoscopeReco:    x: " << recHodo->posX()
00552                              << "\ty: " << recHodo->posY()
00553                              << "\t sx: " << recHodo->slopeX() << "\t qualx: " << recHodo->qualX()
00554                              << "\t sy: " << recHodo->slopeY() << "\t qualy: " << recHodo->qualY();
00555   } else {
00556     edm::LogError("EcalBeamTask") << "Error! Can't get the product EcalTBHodoscopeRecInfo";
00557     meMissingCollections_-> Fill(6); // bin7: missing reconstructed hodoscopes
00558     return;
00559   }
00560 
00561   if ( (!recHodo) || (!recTDC) ) {
00562       edm::LogWarning("EcalBeamTask") << "analyze: missing a needed collection, recHodo or recTDC. Returning.";
00563       return;
00564   }
00565   LogDebug("EEBeamHodoTask") << " Hodo reco and TDC reco found.";
00566 
00567   meTDCRec_->Fill( recTDC->offset());
00568 
00569   meHodoPosRecXY_->Fill( recHodo->posX(), recHodo->posY() );
00570   meHodoPosRecX_->Fill( recHodo->posX());
00571   meHodoPosRecY_->Fill( recHodo->posY() );
00572   meHodoSloXRec_->Fill( recHodo->slopeX());
00573   meHodoSloYRec_->Fill( recHodo->slopeY());
00574   meHodoQuaXRec_->Fill( recHodo->qualX());
00575   meHodoQuaYRec_->Fill( recHodo->qualY());
00576 
00577   /**************************************
00578   // handling histos type II:
00579   **************************************/
00580 
00581   if (tableIsMoving_) {
00582       LogDebug("EcalBeamTask")<< "At event number:" << LV1_ << " table is moving. Not filling concerned monitoring elements. ";
00583       return;
00584   } else {
00585       LogDebug("EcalBeamTask")<< "At event number:" << LV1_ << " table is not moving - thus filling alos monitoring elements requiring so.";
00586   }
00587 
00588   float maxE =0;
00589   EBDetId maxHitId(0);
00590   for (  EEUncalibratedRecHitCollection::const_iterator uncalHitItr = pUncalRH->begin();  uncalHitItr!= pUncalRH->end(); uncalHitItr++ ) {
00591     double e = (*uncalHitItr).amplitude();
00592     if ( e <= 0. ) e = 0.0;
00593     if ( e > maxE )  {
00594       maxE       = e;
00595       maxHitId = (*uncalHitItr).id();
00596     }
00597   }
00598   if (  maxHitId == EBDetId(0) ) {
00599       edm::LogError("EEBeamHodoTask") << "No positive UncalRecHit found in ECAL in event " << ievt_ << " - returning.";
00600       return;
00601   }
00602 
00603   meEvsXRecProf_-> Fill(recHodo->posX(), maxE);
00604   meEvsYRecProf_-> Fill(recHodo->posY(), maxE);
00605   meEvsXRecHis_-> Fill(recHodo->posX(), maxE);
00606   meEvsYRecHis_-> Fill(recHodo->posY(), maxE);
00607 
00608   edm::LogInfo("EcalBeamTask")<< " channel with max is " << maxHitId;
00609 
00610   bool mat5x5 =true;
00611   int ietaMax = maxHitId.ieta();
00612   int iphiMax = (maxHitId.iphi() % 20);
00613   if (ietaMax ==1 || ietaMax ==2 || ietaMax ==84 || ietaMax == 85 ||
00614       iphiMax ==1 || iphiMax ==2 || iphiMax ==19 || iphiMax == 20 )
00615     {mat5x5 =false;}
00616 
00617   if (!mat5x5) return;
00618 
00619   EBDetId Xtals5x5[25];
00620   double   ene5x5[25];
00621   double   e25    =0;
00622   for (unsigned int icry=0;icry<25;icry++)
00623     {
00624       unsigned int row = icry / 5;
00625       unsigned int column= icry %5;
00626       if ( EBDetId::validDetId(maxHitId.ieta()+column-2,maxHitId.iphi()+row-2) ) {
00627         Xtals5x5[icry]=EBDetId(maxHitId.ieta()+column-2,maxHitId.iphi()+row-2,EBDetId::ETAPHIMODE);
00628         double e = ( *pUncalRH->find( Xtals5x5[icry] )  ).amplitude();
00629         if ( e <= 0. ) e = 0.0;
00630         ene5x5[icry] =e;
00631         e25 +=e;
00632       } else {
00633         LogDebug("EcalBeamTask")<< "Cannot construct 5x5 matrix around EBDetId " << maxHitId;
00634         mat5x5 =false;
00635       }
00636     }
00637 
00638   if (!mat5x5) return;
00639   LogDebug("EcalBeamTask")<< "Could construct 5x5 matrix around EBDetId " << maxHitId;
00640 
00641   // im mm
00642   float sideX =24.06;
00643   float sideY =22.02;
00644   float caloX =0;
00645   float caloY =0;
00646   float weight=0;
00647   float sumWeight=0;
00648   // X and Y calculated from log-weighted energy center of mass
00649   for (unsigned int icry=0;icry<25;icry++)
00650     {
00651       weight = log( ene5x5[icry] / e25) + 3.8;
00652       if (weight>0)
00653         {
00654           unsigned int row      = icry / 5;
00655           unsigned int column = icry %5;
00656           caloX +=  (column-2) * sideX * weight;
00657           caloY -=  (row-2) * sideY * weight;
00658           sumWeight += weight;
00659         }
00660     }
00661   caloX /=sumWeight;
00662   caloY /=sumWeight;
00663 
00664   meCaloVsHodoXPos_->Fill( recHodo->posX()-caloX );
00665   meCaloVsHodoYPos_->Fill( recHodo->posY()-caloY );
00666   meCaloVsTDCTime_->Fill( (*pUncalRH->find( maxHitId )  ).jitter()  -  recTDC->offset() - 3);
00667   LogDebug("EcalBeamTask")<< "jiitter from uncalRecHit: " <<  (*pUncalRH->find( maxHitId )  ).jitter();
00668 
00669 }
00670
/data/refman/pasoursint/CMSSW_4_1_8_patch13/src/DQM/EcalEndcapMonitorTasks/src/EEBeamHodoTask.cc
