CMSSW Reference Manual

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