CMS 3D CMS Logo

List of all members | Public Member Functions | Static Public Member Functions | Private Member Functions | Private Attributes
BtlDigiHitsValidation Class Reference

#include <Validation/MtdValidation/plugins/BtlDigiHitsValidation.cc>

Inheritance diagram for BtlDigiHitsValidation:
one::DQMEDAnalyzer< T > one::dqmimplementation::DQMBaseClass< T... >

Public Member Functions

 BtlDigiHitsValidation (const edm::ParameterSet &)
 
 ~BtlDigiHitsValidation () override
 
- Public Member Functions inherited from one::DQMEDAnalyzer< T >
 DQMEDAnalyzer ()=default
 
 DQMEDAnalyzer (DQMEDAnalyzer< T... > const &)=delete
 
 DQMEDAnalyzer (DQMEDAnalyzer< T... > &&)=delete
 
 ~DQMEDAnalyzer () override=default
 

Static Public Member Functions

static void fillDescriptions (edm::ConfigurationDescriptions &descriptions)
 

Private Member Functions

void analyze (const edm::Event &, const edm::EventSetup &) override
 
void bookHistograms (DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &) override
 

Private Attributes

edm::EDGetTokenT< BTLDigiCollectionbtlDigiHitsToken_
 
const std::string folder_
 
MonitorElementmeHitCharge_ [2]
 
MonitorElementmeHitEta_ [2]
 
MonitorElementmeHitPhi_ [2]
 
MonitorElementmeHitQvsEta_ [2]
 
MonitorElementmeHitQvsPhi_ [2]
 
MonitorElementmeHitQvsZ_ [2]
 
MonitorElementmeHitTime_ [2]
 
MonitorElementmeHitTvsEta_ [2]
 
MonitorElementmeHitTvsPhi_ [2]
 
MonitorElementmeHitTvsQ_ [2]
 
MonitorElementmeHitTvsZ_ [2]
 
MonitorElementmeHitX_ [2]
 
MonitorElementmeHitY_ [2]
 
MonitorElementmeHitZ_ [2]
 
MonitorElementmeNhits_ [2]
 
MonitorElementmeOccupancy_ [2]
 

Detailed Description

Description: BTL DIGI hits validation

Implementation: [Notes on implementation]

Definition at line 36 of file BtlDigiHitsValidation.cc.

Constructor & Destructor Documentation

BtlDigiHitsValidation::BtlDigiHitsValidation ( const edm::ParameterSet iConfig)
explicit

Definition at line 79 of file BtlDigiHitsValidation.cc.

References btlDigiHitsToken_, and edm::ParameterSet::getParameter().

80  : folder_(iConfig.getParameter<std::string>("folder")) {
81  btlDigiHitsToken_ = consumes<BTLDigiCollection>(iConfig.getParameter<edm::InputTag>("inputTag"));
82 }
T getParameter(std::string const &) const
edm::EDGetTokenT< BTLDigiCollection > btlDigiHitsToken_
BtlDigiHitsValidation::~BtlDigiHitsValidation ( )
override

Definition at line 84 of file BtlDigiHitsValidation.cc.

84 {}

Member Function Documentation

void BtlDigiHitsValidation::analyze ( const edm::Event iEvent,
const edm::EventSetup iSetup 
)
overrideprivate

Definition at line 87 of file BtlDigiHitsValidation.cc.

References ecalMGPA::adc(), btlDigiHitsToken_, BTLDetId::column(), TauDecayModes::dec, Exception, MonitorElement::Fill(), BTLDetId::geographicalId(), relativeConstraints::geom, edm::EventSetup::get(), edm::Event::getHandle(), MTDTopology::getMTDTopologyMode(), MTDGeometry::idToDet(), edm::makeValid(), meHitCharge_, meHitEta_, meHitPhi_, meHitQvsEta_, meHitQvsPhi_, meHitQvsZ_, meHitTime_, meHitTvsEta_, meHitTvsPhi_, meHitTvsQ_, meHitTvsZ_, meHitX_, meHitY_, meHitZ_, meNhits_, meOccupancy_, RectangularMTDTopology::nrows(), RectangularMTDTopology::pixelToModuleLocalPoint(), edm::ESHandle< T >::product(), DetId::rawId(), BTLDetId::row(), ProxyMTDTopology::specificTopology(), GeomDet::toGlobal(), GeomDet::topology(), and ecaldqm::topology().

87  {
88  using namespace edm;
89 
90  edm::ESHandle<MTDGeometry> geometryHandle;
91  iSetup.get<MTDDigiGeometryRecord>().get(geometryHandle);
92  const MTDGeometry* geom = geometryHandle.product();
93 
94  edm::ESHandle<MTDTopology> topologyHandle;
95  iSetup.get<MTDTopologyRcd>().get(topologyHandle);
96  const MTDTopology* topology = topologyHandle.product();
97 
98  auto btlDigiHitsHandle = makeValid(iEvent.getHandle(btlDigiHitsToken_));
99 
100  // --- Loop over the BLT DIGI hits
101 
102  unsigned int n_digi_btl[2] = {0, 0};
103 
104  for (const auto& dataFrame : *btlDigiHitsHandle) {
105  BTLDetId detId = dataFrame.id();
106  DetId geoId = detId.geographicalId(static_cast<BTLDetId::CrysLayout>(topology->getMTDTopologyMode()));
107  const MTDGeomDet* thedet = geom->idToDet(geoId);
108  if (thedet == nullptr)
109  throw cms::Exception("BtlDigiHitsValidation") << "GeographicalID: " << std::hex << geoId.rawId() << " ("
110  << detId.rawId() << ") is invalid!" << std::dec << std::endl;
111  const ProxyMTDTopology& topoproxy = static_cast<const ProxyMTDTopology&>(thedet->topology());
112  const RectangularMTDTopology& topo = static_cast<const RectangularMTDTopology&>(topoproxy.specificTopology());
113 
114  Local3DPoint local_point(0., 0., 0.);
115  local_point = topo.pixelToModuleLocalPoint(local_point, detId.row(topo.nrows()), detId.column(topo.nrows()));
116  const auto& global_point = thedet->toGlobal(local_point);
117 
118  const auto& sample_L = dataFrame.sample(0);
119  const auto& sample_R = dataFrame.sample(1);
120 
121  uint32_t adc[2] = {sample_L.data(), sample_R.data()};
122  uint32_t tdc[2] = {sample_L.toa(), sample_R.toa()};
123 
124  for (int iside = 0; iside < 2; ++iside) {
125  if (adc[iside] == 0)
126  continue;
127 
128  meHitCharge_[iside]->Fill(adc[iside]);
129  meHitTime_[iside]->Fill(tdc[iside]);
130 
131  meOccupancy_[iside]->Fill(global_point.z(), global_point.phi());
132 
133  meHitX_[iside]->Fill(global_point.x());
134  meHitY_[iside]->Fill(global_point.y());
135  meHitZ_[iside]->Fill(global_point.z());
136  meHitPhi_[iside]->Fill(global_point.phi());
137  meHitEta_[iside]->Fill(global_point.eta());
138 
139  meHitTvsQ_[iside]->Fill(adc[iside], tdc[iside]);
140  meHitQvsPhi_[iside]->Fill(global_point.phi(), adc[iside]);
141  meHitQvsEta_[iside]->Fill(global_point.eta(), adc[iside]);
142  meHitQvsZ_[iside]->Fill(global_point.z(), adc[iside]);
143  meHitTvsPhi_[iside]->Fill(global_point.phi(), tdc[iside]);
144  meHitTvsEta_[iside]->Fill(global_point.eta(), tdc[iside]);
145  meHitTvsZ_[iside]->Fill(global_point.z(), tdc[iside]);
146 
147  n_digi_btl[iside]++;
148 
149  } // iside loop
150 
151  } // dataFrame loop
152 
153  meNhits_[0]->Fill(n_digi_btl[0]);
154  meNhits_[1]->Fill(n_digi_btl[1]);
155 }
MonitorElement * meHitZ_[2]
CaloTopology const * topology(0)
virtual const Topology & topology() const
Definition: GeomDet.cc:81
GlobalPoint toGlobal(const Local2DPoint &lp) const
Conversion to the global R.F. from the R.F. of the GeomDet.
Definition: GeomDet.h:54
int getMTDTopologyMode() const
Definition: MTDTopology.h:74
MonitorElement * meHitQvsEta_[2]
constexpr uint32_t rawId() const
get the raw id
Definition: DetId.h:50
MonitorElement * meHitTvsZ_[2]
MonitorElement * meHitY_[2]
MonitorElement * meHitX_[2]
edm::EDGetTokenT< BTLDigiCollection > btlDigiHitsToken_
int nrows() const override
void Fill(long long x)
Handle< PROD > getHandle(EDGetTokenT< PROD > token) const
Definition: Event.h:539
const MTDGeomDet * idToDet(DetId) const override
Definition: MTDGeometry.cc:184
BTLDetId geographicalId(CrysLayout lay) const
Definition: BTLDetId.cc:162
MonitorElement * meHitTvsEta_[2]
virtual const PixelTopology & specificTopology() const
constexpr int adc(sample_type sample)
get the ADC sample (12 bits)
LocalPoint pixelToModuleLocalPoint(const LocalPoint &plp, int row, int col) const
MonitorElement * meNhits_[2]
MonitorElement * meOccupancy_[2]
Definition: DetId.h:18
MonitorElement * meHitTvsQ_[2]
MonitorElement * meHitQvsZ_[2]
MonitorElement * meHitEta_[2]
MonitorElement * meHitPhi_[2]
HLT enums.
T get() const
Definition: EventSetup.h:71
Detector identifier class for the Barrel Timing Layer. The crystal count must start from 0...
Definition: BTLDetId.h:18
MonitorElement * meHitTime_[2]
MonitorElement * meHitTvsPhi_[2]
int column(unsigned nrows=16) const
Definition: BTLDetId.h:110
MonitorElement * meHitCharge_[2]
auto makeValid(const U &iOtherHandleType) noexcept(false)
Definition: ValidHandle.h:59
T const * product() const
Definition: ESHandle.h:86
MonitorElement * meHitQvsPhi_[2]
int row(unsigned nrows=16) const
Definition: BTLDetId.h:105
void BtlDigiHitsValidation::bookHistograms ( DQMStore::IBooker ibook,
edm::Run const &  run,
edm::EventSetup const &  iSetup 
)
overrideprivate

Definition at line 158 of file BtlDigiHitsValidation.cc.

References DQMStore::IBooker::book1D(), DQMStore::IBooker::book2D(), DQMStore::IBooker::bookProfile(), folder_, meHitCharge_, meHitEta_, meHitPhi_, meHitQvsEta_, meHitQvsPhi_, meHitQvsZ_, meHitTime_, meHitTvsEta_, meHitTvsPhi_, meHitTvsQ_, meHitTvsZ_, meHitX_, meHitY_, meHitZ_, meNhits_, meOccupancy_, and DQMStore::IBooker::setCurrentFolder().

160  {
161  ibook.setCurrentFolder(folder_);
162 
163  // --- histograms booking
164 
165  meNhits_[0] = ibook.book1D("BtlNhitsL", "Number of BTL DIGI hits (L);N_{DIGI}", 100, 0., 5000.);
166  meNhits_[1] = ibook.book1D("BtlNhitsR", "Number of BTL DIGI hits (R);N_{DIGI}", 100, 0., 5000.);
167 
168  meHitCharge_[0] = ibook.book1D("BtlHitChargeL", "BTL DIGI hits charge (L);Q_{DIGI} [ADC counts]", 100, 0., 1024.);
169  meHitCharge_[1] = ibook.book1D("BtlHitChargeR", "BTL DIGI hits charge (R);Q_{DIGI} [ADC counts]", 100, 0., 1024.);
170  meHitTime_[0] = ibook.book1D("BtlHitTimeL", "BTL DIGI hits ToA (L);ToA_{DIGI} [TDC counts]", 100, 0., 1024.);
171  meHitTime_[1] = ibook.book1D("BtlHitTimeR", "BTL DIGI hits ToA (R);ToA_{DIGI} [TDC counts]", 100, 0., 1024.);
172 
173  meOccupancy_[0] = ibook.book2D("BtlOccupancyL",
174  "BTL DIGI hits occupancy (L);Z_{DIGI} [cm]; #phi_{DIGI} [rad]",
175  65,
176  -260.,
177  260.,
178  126,
179  -3.15,
180  3.15);
181  meOccupancy_[1] = ibook.book2D("BtlOccupancyR",
182  "BTL DIGI hits occupancy (R);Z_{DIGI} [cm]; #phi_{DIGI} [rad]",
183  65,
184  -260.,
185  260.,
186  126,
187  -3.15,
188  3.15);
189 
190  meHitX_[0] = ibook.book1D("BtlHitXL", "BTL DIGI hits X (L);X_{DIGI} [cm]", 60, -120., 120.);
191  meHitX_[1] = ibook.book1D("BtlHitXR", "BTL DIGI hits X (R);X_{DIGI} [cm]", 60, -120., 120.);
192  meHitY_[0] = ibook.book1D("BtlHitYL", "BTL DIGI hits Y (L);Y_{DIGI} [cm]", 60, -120., 120.);
193  meHitY_[1] = ibook.book1D("BtlHitYR", "BTL DIGI hits Y (R);Y_{DIGI} [cm]", 60, -120., 120.);
194  meHitZ_[0] = ibook.book1D("BtlHitZL", "BTL DIGI hits Z (L);Z_{DIGI} [cm]", 100, -260., 260.);
195  meHitZ_[1] = ibook.book1D("BtlHitZR", "BTL DIGI hits Z (R);Z_{DIGI} [cm]", 100, -260., 260.);
196  meHitPhi_[0] = ibook.book1D("BtlHitPhiL", "BTL DIGI hits #phi (L);#phi_{DIGI} [rad]", 126, -3.15, 3.15);
197  meHitPhi_[1] = ibook.book1D("BtlHitPhiR", "BTL DIGI hits #phi (R);#phi_{DIGI} [rad]", 126, -3.15, 3.15);
198  meHitEta_[0] = ibook.book1D("BtlHitEtaL", "BTL DIGI hits #eta (L);#eta_{DIGI}", 100, -1.55, 1.55);
199  meHitEta_[1] = ibook.book1D("BtlHitEtaR", "BTL DIGI hits #eta (R);#eta_{DIGI}", 100, -1.55, 1.55);
200 
201  meHitTvsQ_[0] = ibook.bookProfile("BtlHitTvsQL",
202  "BTL DIGI ToA vs charge (L);Q_{DIGI} [ADC counts];ToA_{DIGI} [TDC counts]",
203  50,
204  0.,
205  1024.,
206  0.,
207  1024.);
208  meHitTvsQ_[1] = ibook.bookProfile("BtlHitTvsQR",
209  "BTL DIGI ToA vs charge (R);Q_{DIGI} [ADC counts];ToA_{DIGI} [TDC counts]",
210  50,
211  0.,
212  1024.,
213  0.,
214  1024.);
215  meHitQvsPhi_[0] = ibook.bookProfile("BtlHitQvsPhiL",
216  "BTL DIGI charge vs #phi (L);#phi_{DIGI} [rad];Q_{DIGI} [ADC counts]",
217  50,
218  -3.15,
219  3.15,
220  0.,
221  1024.);
222  meHitQvsPhi_[1] = ibook.bookProfile("BtlHitQvsPhiR",
223  "BTL DIGI charge vs #phi (R);#phi_{DIGI} [rad];Q_{DIGI} [ADC counts]",
224  50,
225  -3.15,
226  3.15,
227  0.,
228  1024.);
229  meHitQvsEta_[0] = ibook.bookProfile(
230  "BtlHitQvsEtaL", "BTL DIGI charge vs #eta (L);#eta_{DIGI};Q_{DIGI} [ADC counts]", 50, -1.55, 1.55, 0., 1024.);
231  meHitQvsEta_[1] = ibook.bookProfile(
232  "BtlHitQvsEtaR", "BTL DIGI charge vs #eta (R);#eta_{DIGI};Q_{DIGI} [ADC counts]", 50, -1.55, 1.55, 0., 1024.);
233  meHitQvsZ_[0] = ibook.bookProfile(
234  "BtlHitQvsZL", "BTL DIGI charge vs Z (L);Z_{DIGI} [cm];Q_{DIGI} [ADC counts]", 50, -260., 260., 0., 1024.);
235  meHitQvsZ_[1] = ibook.bookProfile(
236  "BtlHitQvsZR", "BTL DIGI charge vs Z (R);Z_{DIGI} [cm];Q_{DIGI} [ADC counts]", 50, -260., 260., 0., 1024.);
237  meHitTvsPhi_[0] = ibook.bookProfile(
238  "BtlHitTvsPhiL", "BTL DIGI ToA vs #phi (L);#phi_{DIGI} [rad];ToA_{DIGI} [TDC counts]", 50, -3.15, 3.15, 0., 1024.);
239  meHitTvsPhi_[1] = ibook.bookProfile(
240  "BtlHitTvsPhiR", "BTL DIGI ToA vs #phi (R);#phi_{DIGI} [rad];ToA_{DIGI} [TDC counts]", 50, -3.15, 3.15, 0., 1024.);
241  meHitTvsEta_[0] = ibook.bookProfile(
242  "BtlHitTvsEtaL", "BTL DIGI ToA vs #eta (L);#eta_{DIGI};ToA_{DIGI} [TDC counts]", 50, -1.55, 1.55, 0., 1024.);
243  meHitTvsEta_[1] = ibook.bookProfile(
244  "BtlHitTvsEtaR", "BTL DIGI ToA vs #eta (R);#eta_{DIGI};ToA_{DIGI} [TDC counts]", 50, -1.55, 1.55, 0., 1024.);
245  meHitTvsZ_[0] = ibook.bookProfile(
246  "BtlHitTvsZL", "BTL DIGI ToA vs Z (L);Z_{DIGI} [cm];ToA_{DIGI} [TDC counts]", 50, -260., 260., 0., 1024.);
247  meHitTvsZ_[1] = ibook.bookProfile(
248  "BtlHitTvsZR", "BTL DIGI ToA vs Z (R);Z_{DIGI} [cm];ToA_{DIGI} [TDC counts]", 50, -260., 260., 0., 1024.);
249 }
MonitorElement * meHitZ_[2]
MonitorElement * bookProfile(Args &&...args)
Definition: DQMStore.h:113
MonitorElement * meHitQvsEta_[2]
MonitorElement * meHitTvsZ_[2]
MonitorElement * meHitY_[2]
MonitorElement * meHitX_[2]
void setCurrentFolder(std::string const &fullpath)
Definition: DQMStore.cc:268
MonitorElement * meHitTvsEta_[2]
MonitorElement * book1D(Args &&...args)
Definition: DQMStore.h:106
MonitorElement * meNhits_[2]
MonitorElement * meOccupancy_[2]
MonitorElement * book2D(Args &&...args)
Definition: DQMStore.h:109
MonitorElement * meHitTvsQ_[2]
MonitorElement * meHitQvsZ_[2]
MonitorElement * meHitEta_[2]
MonitorElement * meHitPhi_[2]
MonitorElement * meHitTime_[2]
MonitorElement * meHitTvsPhi_[2]
MonitorElement * meHitCharge_[2]
MonitorElement * meHitQvsPhi_[2]
void BtlDigiHitsValidation::fillDescriptions ( edm::ConfigurationDescriptions descriptions)
static

Definition at line 252 of file BtlDigiHitsValidation.cc.

References edm::ConfigurationDescriptions::add(), edm::ParameterSetDescription::add(), DEFINE_FWK_MODULE, and AlCaHLTBitMon_QueryRunRegistry::string.

252  {
254 
255  desc.add<std::string>("folder", "MTD/BTL/DigiHits");
256  desc.add<edm::InputTag>("inputTag", edm::InputTag("mix", "FTLBarrel"));
257 
258  descriptions.add("btlDigiHitsDefault", desc);
259 }
ParameterDescriptionBase * add(U const &iLabel, T const &value)
void add(std::string const &label, ParameterSetDescription const &psetDescription)

Member Data Documentation

edm::EDGetTokenT<BTLDigiCollection> BtlDigiHitsValidation::btlDigiHitsToken_
private

Definition at line 52 of file BtlDigiHitsValidation.cc.

Referenced by analyze(), and BtlDigiHitsValidation().

const std::string BtlDigiHitsValidation::folder_
private

Definition at line 50 of file BtlDigiHitsValidation.cc.

Referenced by bookHistograms().

MonitorElement* BtlDigiHitsValidation::meHitCharge_[2]
private

Definition at line 58 of file BtlDigiHitsValidation.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* BtlDigiHitsValidation::meHitEta_[2]
private

Definition at line 67 of file BtlDigiHitsValidation.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* BtlDigiHitsValidation::meHitPhi_[2]
private

Definition at line 66 of file BtlDigiHitsValidation.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* BtlDigiHitsValidation::meHitQvsEta_[2]
private

Definition at line 71 of file BtlDigiHitsValidation.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* BtlDigiHitsValidation::meHitQvsPhi_[2]
private

Definition at line 70 of file BtlDigiHitsValidation.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* BtlDigiHitsValidation::meHitQvsZ_[2]
private

Definition at line 72 of file BtlDigiHitsValidation.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* BtlDigiHitsValidation::meHitTime_[2]
private

Definition at line 59 of file BtlDigiHitsValidation.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* BtlDigiHitsValidation::meHitTvsEta_[2]
private

Definition at line 74 of file BtlDigiHitsValidation.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* BtlDigiHitsValidation::meHitTvsPhi_[2]
private

Definition at line 73 of file BtlDigiHitsValidation.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* BtlDigiHitsValidation::meHitTvsQ_[2]
private

Definition at line 69 of file BtlDigiHitsValidation.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* BtlDigiHitsValidation::meHitTvsZ_[2]
private

Definition at line 75 of file BtlDigiHitsValidation.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* BtlDigiHitsValidation::meHitX_[2]
private

Definition at line 63 of file BtlDigiHitsValidation.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* BtlDigiHitsValidation::meHitY_[2]
private

Definition at line 64 of file BtlDigiHitsValidation.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* BtlDigiHitsValidation::meHitZ_[2]
private

Definition at line 65 of file BtlDigiHitsValidation.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* BtlDigiHitsValidation::meNhits_[2]
private

Definition at line 56 of file BtlDigiHitsValidation.cc.

Referenced by analyze(), and bookHistograms().

MonitorElement* BtlDigiHitsValidation::meOccupancy_[2]
private

Definition at line 61 of file BtlDigiHitsValidation.cc.

Referenced by analyze(), and bookHistograms().