#include <DQM/TrigXMonitor/interface/L1Scalers.h>
Definition at line 14 of file L1Scalers.h.
L1Scalers::L1Scalers | ( | const edm::ParameterSet & | ps | ) |
Constructors.
Definition at line 34 of file L1Scalers.cc.
References LogDebug.
00034 : 00035 dbe_(0), nev_(0), 00036 verbose_(ps.getUntrackedParameter < bool > ("verbose", false)), 00037 l1GtDataSource_( ps.getParameter< edm::InputTag >("l1GtData")), 00038 folderName_( ps.getUntrackedParameter< std::string>("dqmFolder", 00039 std::string("L1T/L1Scalers_EvF"))), 00040 l1scalers_(0), 00041 l1techScalers_(0), 00042 l1Correlations_(0), 00043 bxNum_(0), 00044 l1scalersBx_(0), 00045 l1techScalersBx_(0), 00046 nLumiBlock_(0), 00047 l1AlgoCounter_(0), 00048 l1TtCounter_(0), 00049 rateAlgoCounter_(0), 00050 rateTtCounter_(0), 00051 pixFedSize_(0), 00052 hfEnergy_(0), 00053 fedStart_(ps.getUntrackedParameter<unsigned int>("firstFED", 0)), 00054 fedStop_(ps.getUntrackedParameter<unsigned int>("lastFED", 931)), 00055 fedRawCollection_(ps.getParameter<edm::InputTag>("fedRawData")), 00056 maskedList_(ps.getUntrackedParameter<std::vector<int> >("maskedChannels", maskedList_)), //this is using the ashed index 00057 HcalRecHitCollection_(ps.getParameter<edm::InputTag>("HFRecHitCollection")) 00058 { 00059 LogDebug("Status") << "constructor" ; 00060 }
virtual L1Scalers::~L1Scalers | ( | ) | [inline, virtual] |
void L1Scalers::analyze | ( | const edm::Event & | e, | |
const edm::EventSetup & | c | |||
) | [virtual] |
Implements edm::EDAnalyzer.
Definition at line 118 of file L1Scalers.cc.
References L1GtfeWord::bxNr(), bxNum_, CaloRecHit::energy(), fedRawCollection_, fedStart_, fedStop_, MonitorElement::Fill(), find(), edm::Event::getByLabel(), HcalRecHitCollection_, hfEnergy_, i, HFRecHit::id(), j, l1AlgoCounter_, l1Correlations_, l1GtDataSource_, l1scalers_, l1scalersBx_, l1techScalers_, l1techScalersBx_, l1TtCounter_, edm::InputTag::label(), LogDebug, maskedList_, nev_, pixFedSize_, rateAlgoCounter_, rateTtCounter_, HLT_VtxMuL3::result, t, and tw.
00119 { 00120 nev_++; 00121 LogDebug("Status") << "L1Scalers::analyze event " << nev_ ; 00122 00123 // get Global Trigger decision and the decision word 00124 // these are locally derived 00125 edm::Handle<L1GlobalTriggerReadoutRecord> gtRecord; 00126 bool t = e.getByLabel(l1GtDataSource_,gtRecord); 00127 if ( ! t ) { 00128 LogDebug("Product") << "can't find L1GlobalTriggerReadoutRecord " 00129 << "with label " << l1GtDataSource_.label() ; 00130 } 00131 else { 00132 00133 // DEBUG 00134 //gtRecord->print(std::cout); 00135 // DEBUG 00136 00137 L1GtfeWord gtfeWord = gtRecord->gtfeWord(); 00138 int gtfeBx = gtfeWord.bxNr(); 00139 bxNum_->Fill(gtfeBx); 00140 00141 // First, the default 00142 // vector of bool 00143 DecisionWord gtDecisionWord = gtRecord->decisionWord(); 00144 if ( ! gtDecisionWord.empty() ) { // if board not there this is zero 00145 // loop over dec. bit to get total rate (no overlap) 00146 for ( int i = 0; i < 128; ++i ) { 00147 if ( gtDecisionWord[i] ) { 00148 rateAlgoCounter_++; 00149 l1AlgoCounter_->Fill(rateAlgoCounter_); 00150 break; 00151 } 00152 } 00153 // loop over decision bits 00154 for ( int i = 0; i < 128; ++i ) { 00155 if ( gtDecisionWord[i] ) { 00156 l1scalers_->Fill(i); 00157 l1scalersBx_->Fill(gtfeBx,i); 00158 for ( int j = i + 1; j < 128; ++j ) { 00159 if ( gtDecisionWord[j] ) { 00160 l1Correlations_->Fill(i,j); 00161 l1Correlations_->Fill(j,i); 00162 } 00163 } 00164 } 00165 } 00166 } 00167 // loop over technical triggers 00168 // vector of bool again. 00169 TechnicalTriggerWord tw = gtRecord->technicalTriggerWord(); 00170 if ( ! tw.empty() ) { 00171 // loop over dec. bit to get total rate (no overlap) 00172 for ( int i = 0; i < 64; ++i ) { 00173 if ( tw[i] ) { 00174 rateTtCounter_++; 00175 l1TtCounter_->Fill(rateTtCounter_); 00176 break; 00177 } 00178 } 00179 for ( int i = 0; i < 64; ++i ) { 00180 if ( tw[i] ) { 00181 l1techScalers_->Fill(i); 00182 l1techScalersBx_->Fill(gtfeBx, i); 00183 } 00184 } 00185 } // ! tw.empty 00186 } // getbylabel succeeded 00187 00188 00189 // HACK 00190 // getting very basic uncalRH 00191 edm::Handle<FEDRawDataCollection> theRaw; 00192 bool getFed = e.getByLabel(fedRawCollection_, theRaw); 00193 if ( ! getFed ) { 00194 edm::LogInfo("FEDSizeFilter") << fedRawCollection_ << " not available"; 00195 } 00196 else { // got the fed raw data 00197 unsigned int totalFEDsize = 0 ; 00198 for (unsigned int i=fedStart_; i<=fedStop_; ++i) { 00199 LogDebug("Parameter") << "Examining fed " << i << " with size " 00200 << theRaw->FEDData(i).size() ; 00201 totalFEDsize += theRaw->FEDData(i).size() ; 00202 } 00203 pixFedSize_->Fill(totalFEDsize); 00204 00205 LogDebug("Parameter") << "Total FED size: " << totalFEDsize; 00206 } 00207 00208 // HF - stolen from HLTrigger/special 00209 // getting very basic uncalRH 00210 edm::Handle<HFRecHitCollection> crudeHits; 00211 bool getHF = e.getByLabel(HcalRecHitCollection_, crudeHits); 00212 if ( ! getHF ) { 00213 LogDebug("Status") << HcalRecHitCollection_ << " not available"; 00214 } 00215 else { 00216 00217 LogDebug("Status") << "Filtering, with " << crudeHits->size() 00218 << " recHits to consider" ; 00219 for ( HFRecHitCollection::const_iterator hitItr = crudeHits->begin(); 00220 hitItr != crudeHits->end(); ++hitItr ) { 00221 HFRecHit hit = (*hitItr); 00222 00223 // masking noisy channels 00224 std::vector<int>::iterator result; 00225 result = std::find( maskedList_.begin(), maskedList_.end(), 00226 HcalDetId(hit.id()).hashed_index() ); 00227 if (result != maskedList_.end()) 00228 continue; 00229 hfEnergy_->Fill(hit.energy()); 00230 00231 } 00232 } 00233 00234 // END HACK 00235 00236 return; 00237 00238 }
void L1Scalers::beginJob | ( | const edm::EventSetup & | c | ) | [virtual] |
BeginJob.
Reimplemented from edm::EDAnalyzer.
Definition at line 65 of file L1Scalers.cc.
References DQMStore::book1D(), DQMStore::book2D(), DQMStore::bookInt(), bxNum_, dbe_, folderName_, hfEnergy_, l1AlgoCounter_, l1Correlations_, l1scalers_, l1scalersBx_, l1techScalers_, l1techScalersBx_, l1TtCounter_, LogDebug, nLumiBlock_, pixFedSize_, DQMStore::setCurrentFolder(), and DQMStore::setVerbose().
00066 { 00067 LogDebug("Status") << "L1Scalers::beginJob()..."; 00068 00069 dbe_ = Service<DQMStore>().operator->(); 00070 if (dbe_ ) { 00071 dbe_->setVerbose(0); 00072 dbe_->setCurrentFolder(folderName_); 00073 00074 00075 l1scalers_ = dbe_->book1D("l1AlgoBits", "L1 Algorithm Bits", 00076 128, -0.5, 127.5); 00077 l1scalersBx_ = dbe_->book2D("l1AlgoBits_Vs_Bx", "L1 Algorithm Bits vs " 00078 "Bunch Number", 00079 3600, -0.5, 3599.5, 00080 128, -0.5, 127.5); 00081 l1Correlations_ = dbe_->book2D("l1Correlations", "L1 Algorithm Bits " 00082 "Correlations", 00083 128, -0.5, 127.5, 00084 128, -0.5, 127.5); 00085 l1techScalers_ = dbe_->book1D("l1TechAlgoBits", "L1 Tech. Trigger Bits", 00086 64, -0.5, 63.5); 00087 l1techScalersBx_ = dbe_->book2D("l1TechAlgoBits_Vs_Bx", "L1 Technical " 00088 "Trigger " 00089 "Bits vs Bunch Number", 00090 3600, -0.5, 3599.5, 64, -0.5, 63.5); 00091 bxNum_ = dbe_->book1D("bxNum", "Bunch number from GTFE", 00092 3600, -0.5, 3599.5); 00093 00094 nLumiBlock_ = dbe_->bookInt("nLumiBlock"); 00095 00096 00097 // l1 total rate 00098 00099 l1AlgoCounter_ = dbe_->bookInt("l1AlgoCounter"); 00100 l1TtCounter_ = dbe_->bookInt("l1TtCounter"); 00101 00102 // early triggers 00103 pixFedSize_ = dbe_->book1D("pixFedSize", "Size of Pixel FED data", 00104 200, 0., 20000.); 00105 hfEnergy_ = dbe_->book1D("hfEnergy", "HF energy", 00106 100, 0., 500.); 00107 00108 } 00109 00110 00111 return; 00112 }
void L1Scalers::beginRun | ( | const edm::Run & | run, | |
const edm::EventSetup & | c | |||
) | [virtual] |
void L1Scalers::endLuminosityBlock | ( | const edm::LuminosityBlock & | lumiSeg, | |
const edm::EventSetup & | c | |||
) | [virtual] |
End LumiBlock DQM Client Diagnostic should be performed here.
Reimplemented from edm::EDAnalyzer.
Definition at line 240 of file L1Scalers.cc.
References MonitorElement::Fill(), edm::LuminosityBlock::id(), edm::LuminosityBlockID::luminosityBlock(), and nLumiBlock_.
00242 { 00243 nLumiBlock_->Fill(lumiSeg.id().luminosityBlock()); 00244 00245 }
void L1Scalers::endRun | ( | const edm::Run & | run, | |
const edm::EventSetup & | c | |||
) | [virtual] |
MonitorElement* L1Scalers::bxNum_ [private] |
DQMStore* L1Scalers::dbe_ [private] |
edm::InputTag L1Scalers::fedRawCollection_ [private] |
unsigned int L1Scalers::fedStart_ [private] |
unsigned int L1Scalers::fedStop_ [private] |
std::string L1Scalers::folderName_ [private] |
MonitorElement* L1Scalers::hfEnergy_ [private] |
MonitorElement* L1Scalers::l1AlgoCounter_ [private] |
MonitorElement* L1Scalers::l1Correlations_ [private] |
edm::InputTag L1Scalers::l1GtDataSource_ [private] |
MonitorElement* L1Scalers::l1scalers_ [private] |
MonitorElement* L1Scalers::l1scalersBx_ [private] |
MonitorElement* L1Scalers::l1techScalers_ [private] |
MonitorElement* L1Scalers::l1techScalersBx_ [private] |
MonitorElement* L1Scalers::l1TtCounter_ [private] |
std::vector<int> L1Scalers::maskedList_ [private] |
int L1Scalers::nev_ [private] |
MonitorElement* L1Scalers::nLumiBlock_ [private] |
MonitorElement* L1Scalers::pixFedSize_ [private] |
unsigned int L1Scalers::rateAlgoCounter_ [private] |
unsigned int L1Scalers::rateTtCounter_ [private] |
unsigned int L1Scalers::threshold_ [private] |
Definition at line 70 of file L1Scalers.h.
bool L1Scalers::verbose_ [private] |
Definition at line 49 of file L1Scalers.h.