CMS 3D CMS Logo

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

#include <CalibTracker/SiStripLorentzAnglePCLMonitor/plugins/SiStripLorentzAnglePCLMonitor.cc>

Inheritance diagram for SiStripLorentzAnglePCLMonitor:
DQMEDAnalyzer edm::stream::EDProducer< edm::GlobalCache< DQMEDAnalyzerGlobalCache >, edm::EndRunProducer, edm::EndLuminosityBlockProducer, edm::Accumulator >

Classes

struct  OnTrackCluster
 

Public Member Functions

 SiStripLorentzAnglePCLMonitor (const edm::ParameterSet &)
 
 ~SiStripLorentzAnglePCLMonitor () override=default
 
- Public Member Functions inherited from DQMEDAnalyzer
void accumulate (edm::Event const &event, edm::EventSetup const &setup) final
 
void beginLuminosityBlock (edm::LuminosityBlock const &lumi, edm::EventSetup const &setup) final
 
void beginRun (edm::Run const &run, edm::EventSetup const &setup) final
 
void beginStream (edm::StreamID id) final
 
 DQMEDAnalyzer ()
 
void endLuminosityBlock (edm::LuminosityBlock const &lumi, edm::EventSetup const &setup) final
 
void endRun (edm::Run const &run, edm::EventSetup const &setup) final
 
virtual bool getCanSaveByLumi ()
 
- Public Member Functions inherited from edm::stream::EDProducer< edm::GlobalCache< DQMEDAnalyzerGlobalCache >, edm::EndRunProducer, edm::EndLuminosityBlockProducer, edm::Accumulator >
 EDProducer ()=default
 
 EDProducer (const EDProducer &)=delete
 
bool hasAbilityToProduceInBeginLumis () const final
 
bool hasAbilityToProduceInBeginProcessBlocks () const final
 
bool hasAbilityToProduceInBeginRuns () const final
 
bool hasAbilityToProduceInEndLumis () const final
 
bool hasAbilityToProduceInEndProcessBlocks () const final
 
bool hasAbilityToProduceInEndRuns () const final
 
const EDProduceroperator= (const EDProducer &)=delete
 

Static Public Member Functions

static void fillDescriptions (edm::ConfigurationDescriptions &descriptions)
 
- Static Public Member Functions inherited from DQMEDAnalyzer
static void globalEndJob (DQMEDAnalyzerGlobalCache const *)
 
static void globalEndLuminosityBlockProduce (edm::LuminosityBlock &lumi, edm::EventSetup const &setup, LuminosityBlockContext const *context)
 
static void globalEndRunProduce (edm::Run &run, edm::EventSetup const &setup, RunContext const *context)
 
static std::unique_ptr< DQMEDAnalyzerGlobalCacheinitializeGlobalCache (edm::ParameterSet const &)
 

Private Member Functions

void analyze (const edm::Event &, const edm::EventSetup &) override
 
void bookHistograms (DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &) override
 
void dqmBeginRun (edm::Run const &, edm::EventSetup const &) override
 
std::string moduleLocationType (const uint32_t &mod, const TrackerTopology *tTopo)
 

Private Attributes

const std::string folder_
 
SiStripLorentzAngleCalibrationHistograms iHists_
 
const edm::EDGetTokenT< TrajTrackAssociationCollectionm_association_token
 
SiStripClusterInfo m_clusterInfo
 
SiStripHashedDetId m_hash
 
const edm::ESGetToken< SiStripLatency, SiStripLatencyRcdm_latencyTokenBR
 
const edm::ESGetToken< SiStripLorentzAngle, SiStripLorentzAngleDepRcdm_lorentzAngleTokenBR
 
const edm::ESGetToken< MagneticField, IdealMagneticFieldRecordm_magFieldTokenBR
 
const edm::ESGetToken< TrackerGeometry, TrackerDigiGeometryRecordm_tkGeomToken
 
const edm::ESGetToken< TrackerGeometry, TrackerDigiGeometryRecordm_tkGeomTokenBR
 
const edm::ESGetToken< TrackerTopology, TrackerTopologyRcdm_topoEsTokenBR
 
const edm::EDGetTokenT< edm::View< reco::Track > > m_tracks_token
 
bool mismatchedBField_
 
bool mismatchedLatency_
 
const bool saveHistosMods_
 

Static Private Attributes

static constexpr float teslaToInverseGeV_ = 2.99792458e-3f
 

Additional Inherited Members

- Public Types inherited from DQMEDAnalyzer
typedef dqm::reco::DQMStore DQMStore
 
typedef dqm::reco::MonitorElement MonitorElement
 
- Public Types inherited from edm::stream::EDProducer< edm::GlobalCache< DQMEDAnalyzerGlobalCache >, edm::EndRunProducer, edm::EndLuminosityBlockProducer, edm::Accumulator >
using CacheTypes = CacheContexts< T... >
 
using GlobalCache = typename CacheTypes::GlobalCache
 
using HasAbility = AbilityChecker< T... >
 
using InputProcessBlockCache = typename CacheTypes::InputProcessBlockCache
 
using LuminosityBlockCache = typename CacheTypes::LuminosityBlockCache
 
using LuminosityBlockContext = LuminosityBlockContextT< LuminosityBlockCache, RunCache, GlobalCache >
 
using LuminosityBlockSummaryCache = typename CacheTypes::LuminosityBlockSummaryCache
 
using RunCache = typename CacheTypes::RunCache
 
using RunContext = RunContextT< RunCache, GlobalCache >
 
using RunSummaryCache = typename CacheTypes::RunSummaryCache
 
- Protected Member Functions inherited from DQMEDAnalyzer
uint64_t meId () const
 
- Protected Attributes inherited from DQMEDAnalyzer
edm::EDPutTokenT< DQMTokenlumiToken_
 
edm::EDPutTokenT< DQMTokenrunToken_
 
unsigned int streamId_
 

Detailed Description

Description: class to book and fill histograms necessary for the online monitoring of the SiStripLorentzAngle

Implementation: Largely taken from https://github.com/robervalwalsh/tracker-la/blob/master/SiStripLAMonitor.cc

Definition at line 60 of file SiStripLorentzAnglePCLMonitor.cc.

Constructor & Destructor Documentation

◆ SiStripLorentzAnglePCLMonitor()

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

Definition at line 111 of file SiStripLorentzAnglePCLMonitor.cc.

112  : m_clusterInfo(consumesCollector()),
113  mismatchedBField_{false},
114  mismatchedLatency_{false},
115  folder_(iConfig.getParameter<std::string>("folder")),
116  saveHistosMods_(iConfig.getParameter<bool>("saveHistoMods")),
117  m_tracks_token(consumes<edm::View<reco::Track>>(iConfig.getParameter<edm::InputTag>("Tracks"))),
118  m_association_token(consumes<TrajTrackAssociationCollection>(iConfig.getParameter<edm::InputTag>("Tracks"))),
119  m_tkGeomToken{esConsumes<>()},
120  m_latencyTokenBR{esConsumes<edm::Transition::BeginRun>()},
121  m_topoEsTokenBR{esConsumes<edm::Transition::BeginRun>()},
122  m_tkGeomTokenBR{esConsumes<edm::Transition::BeginRun>()},
123  m_magFieldTokenBR{esConsumes<edm::Transition::BeginRun>()},
124  m_lorentzAngleTokenBR{esConsumes<edm::Transition::BeginRun>()} {}
T getParameter(std::string const &) const
Definition: ParameterSet.h:307
const edm::ESGetToken< SiStripLorentzAngle, SiStripLorentzAngleDepRcd > m_lorentzAngleTokenBR
const edm::ESGetToken< TrackerGeometry, TrackerDigiGeometryRecord > m_tkGeomTokenBR
const edm::EDGetTokenT< TrajTrackAssociationCollection > m_association_token
const edm::ESGetToken< TrackerGeometry, TrackerDigiGeometryRecord > m_tkGeomToken
const edm::ESGetToken< SiStripLatency, SiStripLatencyRcd > m_latencyTokenBR
const edm::ESGetToken< TrackerTopology, TrackerTopologyRcd > m_topoEsTokenBR
const edm::EDGetTokenT< edm::View< reco::Track > > m_tracks_token
const edm::ESGetToken< MagneticField, IdealMagneticFieldRecord > m_magFieldTokenBR

◆ ~SiStripLorentzAnglePCLMonitor()

SiStripLorentzAnglePCLMonitor::~SiStripLorentzAnglePCLMonitor ( )
overridedefault

Member Function Documentation

◆ analyze()

void SiStripLorentzAnglePCLMonitor::analyze ( const edm::Event iEvent,
const edm::EventSetup iSetup 
)
overrideprivatevirtual

Reimplemented from DQMEDAnalyzer.

Definition at line 222 of file SiStripLorentzAnglePCLMonitor.cc.

References trackingPlots::assoc, bsc_activity_cfg::clusters, HcalObjRepresent::Fill(), SiStripLorentzAngleCalibrationHistograms::h1_, SiStripLorentzAngleCalibrationHistograms::h2_, SiStripLorentzAngleCalibrationHistograms::h2_ct_var2_m_, SiStripLorentzAngleCalibrationHistograms::h2_ct_var3_m_, SiStripLorentzAngleCalibrationHistograms::h2_ct_w_m_, SiStripLorentzAngleCalibrationHistograms::h2_t_var2_m_, SiStripLorentzAngleCalibrationHistograms::h2_t_var3_m_, SiStripLorentzAngleCalibrationHistograms::h2_t_w_m_, ecalpyutils::hashedIndex(), SiStripHashedDetId::hashedIndex(), iEvent, iHists_, LogDebug, m_association_token, m_clusterInfo, m_hash, m_tracks_token, muonTagProbeFilters_cff::matched, mismatchedBField_, mismatchedLatency_, mod(), SiStripLorentzAngleCalibrationHistograms::moduleLocationType_, EcalCalibMonitorClient_cfi::moduleName, SiStripLorentzAngleCalibrationHistograms::orientation_, saveHistosMods_, SiStripClusterInfo::setCluster(), Validation_hcalonly_cfi::sign, AlCaHLTBitMon_QueryRunRegistry::string, theta(), HLT_2023v12_cff::track, pwdgSkimBPark_cfi::tracks, and SiStripClusterInfo::variance().

222  {
223  using namespace edm;
224 
225  // return immediately if the field is not consistent!
226  if (mismatchedBField_)
227  return;
228 
229  if (mismatchedLatency_)
230  return;
231 
233  iEvent.getByToken(m_tracks_token, tracks);
234  edm::Handle<TrajTrackAssociationCollection> trajTrackAssociations;
235  iEvent.getByToken(m_association_token, trajTrackAssociations);
236 
237  LogDebug(moduleDescription().moduleName()) << "I AM IN EVENT" << iEvent.id() << std::endl;
238 
239  std::vector<OnTrackCluster> clusters{};
240 
241  // first collect all the clusters
242  for (const auto& assoc : *trajTrackAssociations) {
243  const auto traj = assoc.key.get();
244  const auto track = assoc.val.get();
245 
246  iHists_.h1_["track_pt"]->Fill(track->pt());
247  iHists_.h1_["track_eta"]->Fill(track->eta());
248  iHists_.h1_["track_phi"]->Fill(track->phi());
249  iHists_.h1_["track_validhits"]->Fill(track->numberOfValidHits());
250 
251  const auto normChi2 = track->ndof() > 0 ? track->chi2() / track->ndof() : -1.;
252  iHists_.h1_["track_chi2ndof"]->Fill(normChi2);
253  iHists_.h2_["track_chi2xhits"]->Fill(normChi2, track->numberOfValidHits());
254  iHists_.h2_["track_ptxhits"]->Fill(track->pt(), track->numberOfValidHits());
255  iHists_.h2_["track_etaxhits"]->Fill(track->eta(), track->numberOfValidHits());
256  iHists_.h2_["track_ptxchi2"]->Fill(track->pt(), normChi2);
257  iHists_.h2_["track_ptxeta"]->Fill(track->pt(), track->eta());
258  iHists_.h2_["track_etaxchi2"]->Fill(track->eta(), normChi2);
259 
260  edm::LogInfo("SiStripLorentzAnglePCLMonitor")
261  << " track pT()" << track->pt() << " track eta()" << track->eta() << std::endl;
262 
263  for (const auto& meas : traj->measurements()) {
264  const auto& trajState = meas.updatedState();
265  if (!trajState.isValid())
266  continue;
267 
268  // there can be 2 (stereo module), 1 (no stereo module), or 0 (no strip hit) clusters per measurement
269  const auto trechit = meas.recHit()->hit();
270  const auto simple1d = dynamic_cast<const SiStripRecHit1D*>(trechit);
271  const auto simple = dynamic_cast<const SiStripRecHit2D*>(trechit);
272  const auto matched = dynamic_cast<const SiStripMatchedRecHit2D*>(trechit);
273  if (matched) {
274  clusters.emplace_back(matched->monoId(), &matched->monoCluster(), traj, track, meas);
275  clusters.emplace_back(matched->stereoId(), &matched->stereoCluster(), traj, track, meas);
276  } else if (simple) {
277  clusters.emplace_back(simple->geographicalId().rawId(), simple->cluster().get(), traj, track, meas);
278  } else if (simple1d) {
279  clusters.emplace_back(simple1d->geographicalId().rawId(), simple1d->cluster().get(), traj, track, meas);
280  }
281  }
282  }
283 
284  for (const auto clus : clusters) {
285  uint32_t c_nstrips = clus.cluster->amplitudes().size();
286  m_clusterInfo.setCluster(*clus.cluster, clus.det);
287  float c_variance = m_clusterInfo.variance();
288  const auto& trajState = clus.measurement.updatedState();
289  const auto trackDir = trajState.localDirection();
290  float c_localdirx = trackDir.x();
291  float c_localdiry = trackDir.y();
292  float c_localdirz = trackDir.z();
293  const auto hit = clus.measurement.recHit()->hit();
294 
295  // not yet needed (might be used for Backplane correction later on
296  /*
297  const auto& tkGeom = iSetup.getData(m_tkGeomToken);
298  const auto stripDet = dynamic_cast<const StripGeomDetUnit*>(tkGeom.idToDet(hit->geographicalId()));
299  float c_barycenter = stripDet->specificTopology().localPosition(clus.cluster->barycenter()).x();
300  float c_localx = stripDet->toLocal(trajState.globalPosition()).x();
301  float c_rhlocalx = hit->localPosition().x();
302  float c_rhlocalxerr = hit->localPositionError().xx();
303  */
304 
305  const uint32_t mod = hit->geographicalId().rawId();
306 
307  std::string locationtype = iHists_.moduleLocationType_[mod];
308  if (locationtype.empty())
309  return;
310 
311  const auto& hashedIndex = m_hash.hashedIndex(mod);
312 
313  TVector3 localdir(c_localdirx, c_localdiry, c_localdirz);
314  int sign = iHists_.orientation_[mod];
315  float tantheta = TMath::Tan(localdir.Theta());
316  float cosphi = TMath::Cos(localdir.Phi());
317  float theta = localdir.Theta();
318 
319  iHists_.h1_[Form("%s_nstrips", locationtype.c_str())]->Fill(c_nstrips);
320  iHists_.h1_[Form("%s_tanthetatrk", locationtype.c_str())]->Fill(sign * tantheta);
321  iHists_.h1_[Form("%s_cosphitrk", locationtype.c_str())]->Fill(cosphi);
322 
323  // nstrips
324  iHists_.h2_[Form("%s_tanthcosphtrk_nstrip", locationtype.c_str())]->Fill(sign * cosphi * tantheta, c_nstrips);
325  iHists_.h2_[Form("%s_thetatrk_nstrip", locationtype.c_str())]->Fill(sign * theta * cosphi, c_nstrips);
326 
327  // variance for width == 2
328  if (c_nstrips == 2) {
329  iHists_.h1_[Form("%s_variance_w2", locationtype.c_str())]->Fill(c_variance);
330  iHists_.h2_[Form("%s_tanthcosphtrk_var2", locationtype.c_str())]->Fill(sign * cosphi * tantheta, c_variance);
331  iHists_.h2_[Form("%s_thcosphtrk_var2", locationtype.c_str())]->Fill(sign * cosphi * theta, c_variance);
332 
333  // not in PCL
334  if (saveHistosMods_) {
335  iHists_.h2_ct_var2_m_[hashedIndex]->Fill(sign * cosphi * tantheta, c_variance);
336  iHists_.h2_t_var2_m_[hashedIndex]->Fill(sign * cosphi * theta, c_variance);
337  }
338  }
339  // variance for width == 3
340  if (c_nstrips == 3) {
341  iHists_.h1_[Form("%s_variance_w3", locationtype.c_str())]->Fill(c_variance);
342  iHists_.h2_[Form("%s_tanthcosphtrk_var3", locationtype.c_str())]->Fill(sign * cosphi * tantheta, c_variance);
343  iHists_.h2_[Form("%s_thcosphtrk_var3", locationtype.c_str())]->Fill(sign * cosphi * theta, c_variance);
344 
345  // not in PCL
346  if (saveHistosMods_) {
347  iHists_.h2_ct_var3_m_[hashedIndex]->Fill(sign * cosphi * tantheta, c_variance);
348  iHists_.h2_t_var3_m_[hashedIndex]->Fill(sign * cosphi * theta, c_variance);
349  }
350  }
351  // not in PCL
352  if (saveHistosMods_) {
353  iHists_.h2_ct_w_m_[hashedIndex]->Fill(sign * cosphi * tantheta, c_nstrips);
354  iHists_.h2_t_w_m_[hashedIndex]->Fill(sign * cosphi * theta, c_nstrips);
355  }
356  }
357 }
void setCluster(const SiStripCluster &cluster, int detId)
std::vector< dqm::reco::MonitorElement * > h2_ct_var3_m_
std::vector< dqm::reco::MonitorElement * > h2_t_w_m_
std::vector< dqm::reco::MonitorElement * > h2_t_var3_m_
uint32_t hashedIndex(uint32_t det_id)
std::vector< dqm::reco::MonitorElement * > h2_ct_w_m_
SiStripLorentzAngleCalibrationHistograms iHists_
int hashedIndex(int ieta, int iphi)
Definition: EcalPyUtils.cc:36
int iEvent
Definition: GenABIO.cc:224
const edm::EDGetTokenT< TrajTrackAssociationCollection > m_association_token
std::map< std::string, dqm::reco::MonitorElement * > h2_
void Fill(HcalDetId &id, double val, std::vector< TH2F > &depth)
std::vector< dqm::reco::MonitorElement * > h2_t_var2_m_
Log< level::Info, false > LogInfo
std::vector< dqm::reco::MonitorElement * > h2_ct_var2_m_
std::map< std::string, dqm::reco::MonitorElement * > h1_
const edm::EDGetTokenT< edm::View< reco::Track > > m_tracks_token
HLT enums.
T mod(const T &a, const T &b)
Definition: ecalDccMap.h:4
Geom::Theta< T > theta() const
#define LogDebug(id)

◆ bookHistograms()

void SiStripLorentzAnglePCLMonitor::bookHistograms ( DQMStore::IBooker ibook,
edm::Run const &  run,
edm::EventSetup const &  iSetup 
)
overrideprivatevirtual

Implements DQMEDAnalyzer.

Definition at line 359 of file SiStripLorentzAnglePCLMonitor.cc.

References SiStripLorentzAngleCalibrationHistograms::apvmode_, SiStripHashedDetId::begin(), SiStripLorentzAngleCalibrationHistograms::bfield_, dqm::implementation::IBooker::book1D(), dqm::implementation::IBooker::book2D(), SiStripHashedDetId::end(), folder_, dqm-mbProfile::format, SiStripLorentzAngleCalibrationHistograms::h1_, SiStripLorentzAngleCalibrationHistograms::h2_, SiStripLorentzAngleCalibrationHistograms::h2_ct_var2_m_, SiStripLorentzAngleCalibrationHistograms::h2_ct_var3_m_, SiStripLorentzAngleCalibrationHistograms::h2_ct_w_m_, SiStripLorentzAngleCalibrationHistograms::h2_t_var2_m_, SiStripLorentzAngleCalibrationHistograms::h2_t_var3_m_, SiStripLorentzAngleCalibrationHistograms::h2_t_w_m_, iHists_, MainPageGenerator::l, hgcalTBTopologyTester_cfi::layers, m_hash, mod(), SiStripLorentzAngleCalibrationHistograms::modtypes_, SiStripLorentzAngleCalibrationHistograms::moduleLocationType_, EcalCalibMonitorClient_cfi::moduleName, SiStripLorentzAngleCalibrationHistograms::nlayers_, saveHistosMods_, dqm::implementation::NavigatorBase::setCurrentFolder(), AlCaHLTBitMon_QueryRunRegistry::string, submitPVValidationJobs::t, and SiStripCondTypes::titles.

361  {
362  std::string bvalue = (iHists_.bfield_ == "3.8") ? "B-ON" : "B-OFF";
363  std::string folderToBook = fmt::format("{}/{}_{}", folder_, bvalue, iHists_.apvmode_);
364 
365  ibook.setCurrentFolder(folderToBook);
366  edm::LogPrint(moduleDescription().moduleName()) << "booking in " << folderToBook;
367 
368  // prepare track histograms
369  // clang-format off
370  iHists_.h1_["track_pt"] = ibook.book1D("track_pt", "track p_{T};track p_{T} [GeV];n. tracks", 2000, 0, 1000);
371  iHists_.h1_["track_eta"] = ibook.book1D("track_eta", "track #eta;track #eta;n. tracks", 100, -4, 4);
372  iHists_.h1_["track_phi"] = ibook.book1D("track_phi", "track #phi;track #phi;n. tracks", 80, -3.2, 3.2);
373  iHists_.h1_["track_validhits"] =
374  ibook.book1D("track_validhits", "track n. valid hits;track n. valid hits;n. tracks", 50, 0, 50);
375  iHists_.h1_["track_chi2ndof"] =
376  ibook.book1D("track_chi2ndof", "track #chi^{2}/ndf;track #chi^{2}/ndf;n. tracks", 100, 0, 5);
377  iHists_.h2_["track_chi2xhits"] =
378  ibook.book2D("track_chi2xhits_2d",
379  "track track n. hits vs track #chi^{2}/ndf;track #chi^{2};track n. valid hits;tracks",
380  100, 0, 5, 50, 0, 50);
381  iHists_.h2_["track_ptxhits"] = ibook.book2D(
382  "track_ptxhits_2d", "track n. hits vs p_{T};track p_{T} [GeV];track n. valid hits;tracks", 200, 0, 100, 50, 0, 50);
383  iHists_.h2_["track_etaxhits"] = ibook.book2D(
384  "track_etaxhits_2d", "track n. hits vs track #eta;track #eta;track n. valid hits;tracks", 60, -3, 3, 50, 0, 50);
385  iHists_.h2_["track_ptxchi2"] =
386  ibook.book2D("track_ptxchi2_2d",
387  "track #chi^{2}/ndf vs track p_{T};track p_{T} [GeV]; track #chi^{2}/ndf;tracks",
388  200, 0, 100, 100, 0, 5);
389  iHists_.h2_["track_ptxeta"] = ibook.book2D(
390  "track_ptxeta_2d", "track #eta vs track p_{T};track p_{T} [GeV];track #eta;tracks", 200, 0, 100, 60, -3, 3);
391  iHists_.h2_["track_etaxchi2"] = ibook.book2D(
392  "track_etaxchi2_2d", "track #chi^{2}/ndf vs track #eta;track #eta;track #chi^{2};tracks", 60, -3, 3, 100, 0, 5);
393  // clang-format on
394 
395  // fill in the module types
396  iHists_.nlayers_["TIB"] = 4;
397  iHists_.nlayers_["TOB"] = 6;
398  iHists_.modtypes_.push_back("s");
399  iHists_.modtypes_.push_back("a");
400 
401  // prepare type histograms
402  for (auto& layers : iHists_.nlayers_) {
403  std::string subdet = layers.first;
404  for (int l = 1; l <= layers.second; ++l) {
405  ibook.setCurrentFolder(folderToBook + Form("/%s/L%d", subdet.c_str(), l));
406  for (auto& t : iHists_.modtypes_) {
407  // do not fill stereo where there aren't
408  if (l > 2 && t == "s")
409  continue;
410  std::string locType = Form("%s_L%d%s", subdet.c_str(), l, t.c_str());
411 
412  // clang-format off
413  const char* titles = Form("n.strips in %s;n.strips;n. clusters", locType.c_str());
414  iHists_.h1_[Form("%s_nstrips", locType.c_str())] = ibook.book1D(Form("%s_nstrips", locType.c_str()), titles, 20, 0, 20);
415 
416  titles = Form("tan(#theta_{trk}) in %s;tan(#theta_{trk});n. clusters", locType.c_str());
417  iHists_.h1_[Form("%s_tanthetatrk", locType.c_str())] = ibook.book1D(Form("%s_tanthetatrk", locType.c_str()), titles, 300, -1.5, 1.5);
418 
419  titles = Form("cos(#phi_{trk}) in %s;cos(#phi_{trk});n. clusters", locType.c_str());
420  iHists_.h1_[Form("%s_cosphitrk", locType.c_str())] = ibook.book1D(Form("%s_cosphitrk", locType.c_str()), titles, 40, -1, 1);
421 
422  titles = Form("Cluster variance (w=2) in %s;cluster variance (w=2);n. clusters", locType.c_str());
423  iHists_.h1_[Form("%s_variance_w2", locType.c_str())] = ibook.book1D(Form("%s_variance_w2", locType.c_str()), titles, 100, 0, 1);
424 
425  titles = Form("Cluster variance (w=3) in %s;cluster variance (w=3);n. clusters", locType.c_str());
426  iHists_.h1_[Form("%s_variance_w3", locType.c_str())] = ibook.book1D(Form("%s_variance_w3", locType.c_str()), titles, 100, 0, 1);
427 
428  titles = Form("n. strips in %s vs tan(#theta_{trk})cos(#phi_{trk});tan(#theta_{trk})cos(#phi_{trk});n. strips;n. clusters", locType.c_str());
429  iHists_.h2_[Form("%s_tanthcosphtrk_nstrip", locType.c_str())] = ibook.book2D(Form("%s_tanthcosphtrk_nstrip", locType.c_str()), titles, 360, -0.9, 0.9, 20, 0, 20);
430 
431  titles = Form("n. strips in %s vs #theta_{trk};#theta_{trk} [rad];n. strips;n. clusters", locType.c_str());
432  iHists_.h2_[Form("%s_thetatrk_nstrip", locType.c_str())] = ibook.book2D(Form("%s_thetatrk_nstrip", locType.c_str()), titles, 360, -0.9, 0.9, 20, 0, 20);
433 
434  titles = Form("cluster variance (w=2) in %s vs tan(#theta_{trk})cos(#phi_{trk});tan(#theta_{trk})cos(#phi_{trk});cluster variance (w=2);n. clusters", locType.c_str());
435  iHists_.h2_[Form("%s_tanthcosphtrk_var2", locType.c_str())] = ibook.book2D(Form("%s_tanthcosphtrk_var2", locType.c_str()), titles, 360, -0.9, 0.9, 50, 0, 1);
436 
437  titles = Form("cluster variance (w=3) in %s vs tan(#theta_{trk})cos(#phi_{trk});tan(#theta_{trk})cos(#phi_{trk});cluster variance (w=3);n. clusters", locType.c_str());
438  iHists_.h2_[Form("%s_tanthcosphtrk_var3", locType.c_str())] = ibook.book2D(Form("%s_tanthcosphtrk_var3", locType.c_str()), titles, 360, -0.9, 0.9, 50, 0, 1);
439 
440  titles = Form("cluster variance (w=2) in %s vs #theta_{trk}cos(#phi_{trk});#theta_{trk}cos(#phi_{trk});cluster variance (w=2);n. clusters", locType.c_str());
441  iHists_.h2_[Form("%s_thcosphtrk_var2", locType.c_str())] = ibook.book2D(Form("%s_thcosphtrk_var2", locType.c_str()), titles, 360, -0.9, 0.9, 50, 0, 1);
442 
443  titles = Form("cluster variance (w=3) in %s vs #theta_{trk}cos(#phi_{trk});#theta_{trk}cos(#phi_{trk});cluster variance (w=3);n. clusters", locType.c_str());
444  iHists_.h2_[Form("%s_thcosphtrk_var3", locType.c_str())] = ibook.book2D(Form("%s_thcosphtrk_var3", locType.c_str()), titles, 360, -0.9, 0.9, 50, 0, 1);
445  // clang-format on
446  }
447  }
448  }
449 
450  // prepare module histograms
451  if (saveHistosMods_) {
452  ibook.setCurrentFolder(folderToBook + "/modules");
453  for (const auto& [mod, locationType] : iHists_.moduleLocationType_) {
454  // histograms for each module
455  iHists_.h1_[Form("%s_%d_nstrips", locationType.c_str(), mod)] =
456  ibook.book1D(Form("%s_%d_nstrips", locationType.c_str(), mod), "", 10, 0, 10);
457  iHists_.h1_[Form("%s_%d_tanthetatrk", locationType.c_str(), mod)] =
458  ibook.book1D(Form("%s_%d_tanthetatrk", locationType.c_str(), mod), "", 40, -1., 1.);
459  iHists_.h1_[Form("%s_%d_cosphitrk", locationType.c_str(), mod)] =
460  ibook.book1D(Form("%s_%d_cosphitrk", locationType.c_str(), mod), "", 40, -1, 1);
461  iHists_.h1_[Form("%s_%d_variance_w2", locationType.c_str(), mod)] =
462  ibook.book1D(Form("%s_%d_variance_w2", locationType.c_str(), mod), "", 20, 0, 1);
463  iHists_.h1_[Form("%s_%d_variance_w3", locationType.c_str(), mod)] =
464  ibook.book1D(Form("%s_%d_variance_w3", locationType.c_str(), mod), "", 20, 0, 1);
465  }
466 
467  int counter{0};
469  for (; iter != m_hash.end(); ++iter) {
470  const auto& locationType = iHists_.moduleLocationType_[(*iter)];
471  iHists_.h2_ct_w_m_.push_back(
472  ibook.book2D(Form("ct_w_m_%s_%d", locationType.c_str(), *iter), "", 90, -0.9, 0.9, 10, 0, 10));
473  iHists_.h2_t_w_m_.push_back(
474  ibook.book2D(Form("t_w_m_%s_%d", locationType.c_str(), *iter), "", 90, -0.9, 0.9, 10, 0, 10));
475  iHists_.h2_ct_var2_m_.push_back(
476  ibook.book2D(Form("ct_var2_m_%s_%d", locationType.c_str(), *iter), "", 90, -0.9, 0.9, 20, 0, 1));
477  iHists_.h2_ct_var3_m_.push_back(
478  ibook.book2D(Form("ct_var3_m_%s_%d", locationType.c_str(), *iter), "", 90, -0.9, 0.9, 20, 0, 1));
479  iHists_.h2_t_var2_m_.push_back(
480  ibook.book2D(Form("t_var2_m_%s_%d", locationType.c_str(), *iter), "", 90, -0.9, 0.9, 20, 0, 1));
481  iHists_.h2_t_var3_m_.push_back(
482  ibook.book2D(Form("t_var3_m_%s_%d", locationType.c_str(), *iter), "", 90, -0.9, 0.9, 20, 0, 1));
483  counter++;
484  }
485  edm::LogPrint(moduleDescription().moduleName())
486  << __PRETTY_FUNCTION__ << " Booked " << counter << " module level histograms!";
487  } // if saveHistoMods
488 }
static const std::array< std::string, 5 > titles
std::vector< uint32_t >::const_iterator const_iterator
virtual void setCurrentFolder(std::string const &fullpath)
Definition: DQMStore.cc:36
std::vector< dqm::reco::MonitorElement * > h2_ct_var3_m_
std::vector< dqm::reco::MonitorElement * > h2_t_w_m_
std::vector< dqm::reco::MonitorElement * > h2_t_var3_m_
std::vector< dqm::reco::MonitorElement * > h2_ct_w_m_
SiStripLorentzAngleCalibrationHistograms iHists_
std::map< std::string, dqm::reco::MonitorElement * > h2_
std::vector< dqm::reco::MonitorElement * > h2_t_var2_m_
Log< level::Warning, true > LogPrint
std::vector< dqm::reco::MonitorElement * > h2_ct_var2_m_
std::map< std::string, dqm::reco::MonitorElement * > h1_
MonitorElement * book2D(TString const &name, TString const &title, int nchX, double lowX, double highX, int nchY, double lowY, double highY, FUNC onbooking=NOOP())
Definition: DQMStore.h:212
const_iterator end() const
MonitorElement * book1D(TString const &name, TString const &title, int const nchX, double const lowX, double const highX, FUNC onbooking=NOOP())
Definition: DQMStore.h:98
T mod(const T &a, const T &b)
Definition: ecalDccMap.h:4
const_iterator begin() const

◆ dqmBeginRun()

void SiStripLorentzAnglePCLMonitor::dqmBeginRun ( edm::Run const &  run,
edm::EventSetup const &  iSetup 
)
overrideprivatevirtual

Reimplemented from DQMEDAnalyzer.

Definition at line 129 of file SiStripLorentzAnglePCLMonitor.cc.

References SiStripLorentzAngleCalibrationHistograms::apvmode_, siStripLACalibration::apvModeAsString(), SiStripLorentzAngleCalibrationHistograms::bfield_, shallow::drift(), siStripLACalibration::fieldAsString(), edm::EventSetup::getData(), SiStripLorentzAngleCalibrationHistograms::h2_ct_var2_m_, SiStripLorentzAngleCalibrationHistograms::h2_ct_var3_m_, SiStripLorentzAngleCalibrationHistograms::h2_ct_w_m_, SiStripLorentzAngleCalibrationHistograms::h2_t_var2_m_, SiStripLorentzAngleCalibrationHistograms::h2_t_var3_m_, SiStripLorentzAngleCalibrationHistograms::h2_t_w_m_, iHists_, SiStripLorentzAngleCalibrationHistograms::la_db_, m_hash, m_latencyTokenBR, m_lorentzAngleTokenBR, m_magFieldTokenBR, m_tkGeomTokenBR, m_topoEsTokenBR, mismatchedBField_, mismatchedLatency_, moduleLocationType(), SiStripLorentzAngleCalibrationHistograms::moduleLocationType_, SiStripLorentzAngleCalibrationHistograms::orientation_, GloballyPositioned< T >::position(), jetUpdater_cfi::sort, GeomDet::surface(), teslaToInverseGeV_, GeomDet::toGlobal(), GloballyPositioned< T >::toLocal(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

129  {
130  const auto& tkGeom = iSetup.getData(m_tkGeomTokenBR);
131  const auto& magField = iSetup.getData(m_magFieldTokenBR);
132  const auto& lorentzAngle = iSetup.getData(m_lorentzAngleTokenBR);
133  const TrackerTopology* tTopo = &iSetup.getData(m_topoEsTokenBR);
134 
135  // fast cachecd access
136  const auto& theMagField = 1.f / (magField.inverseBzAtOriginInGeV() * teslaToInverseGeV_);
137 
138  if (iHists_.bfield_.empty()) {
140  } else {
142  mismatchedBField_ = true;
143  }
144  }
145 
146  const SiStripLatency* apvlat = &iSetup.getData(m_latencyTokenBR);
147  if (iHists_.apvmode_.empty()) {
149  } else {
151  mismatchedLatency_ = true;
152  }
153  }
154 
155  std::vector<uint32_t> c_rawid;
156  std::vector<float> c_globalZofunitlocalY, c_localB, c_BdotY, c_driftx, c_drifty, c_driftz, c_lorentzAngle;
157 
158  auto dets = tkGeom.detsTIB();
159  //dets.insert(dets.end(), tkGeom.detsTID().begin(), tkGeom.detsTID().end()); // no LA in endcaps
160  dets.insert(dets.end(), tkGeom.detsTOB().begin(), tkGeom.detsTOB().end());
161  //dets.insert(dets.end(), tkGeom.detsTEC().begin(), tkGeom.detsTEC().end()); // no LA in endcaps
162 
163  for (auto det : dets) {
164  auto detid = det->geographicalId().rawId();
165  const StripGeomDetUnit* stripDet = dynamic_cast<const StripGeomDetUnit*>(tkGeom.idToDet(det->geographicalId()));
166  if (stripDet) {
167  c_rawid.push_back(detid);
168  c_globalZofunitlocalY.push_back(stripDet->toGlobal(LocalVector(0, 1, 0)).z());
169  iHists_.orientation_[detid] = (stripDet->toGlobal(LocalVector(0, 1, 0)).z() < 0 ? -1 : 1);
170  const auto locB = magField.inTesla(stripDet->surface().position());
171  c_localB.push_back(locB.mag());
172  c_BdotY.push_back(stripDet->surface().toLocal(locB).y());
173  const auto drift = shallow::drift(stripDet, magField, lorentzAngle);
174  c_driftx.push_back(drift.x());
175  c_drifty.push_back(drift.y());
176  c_driftz.push_back(drift.z());
177  c_lorentzAngle.push_back(lorentzAngle.getLorentzAngle(detid));
178  iHists_.la_db_[detid] = lorentzAngle.getLorentzAngle(detid);
179  iHists_.moduleLocationType_[detid] = this->moduleLocationType(detid, tTopo);
180  }
181  }
182 
183  // Sorted DetId list gives max performance, anything else is worse
184  std::sort(c_rawid.begin(), c_rawid.end());
185 
186  // initialized the hash map
187  m_hash = SiStripHashedDetId(c_rawid);
188 
189  //reserve the size of the vector
190  iHists_.h2_ct_w_m_.reserve(c_rawid.size());
191  iHists_.h2_ct_var2_m_.reserve(c_rawid.size());
192  iHists_.h2_ct_var3_m_.reserve(c_rawid.size());
193 
194  iHists_.h2_t_w_m_.reserve(c_rawid.size());
195  iHists_.h2_t_var2_m_.reserve(c_rawid.size());
196  iHists_.h2_t_var3_m_.reserve(c_rawid.size());
197 }
const edm::ESGetToken< SiStripLorentzAngle, SiStripLorentzAngleDepRcd > m_lorentzAngleTokenBR
Local3DVector LocalVector
Definition: LocalVector.h:12
std::string moduleLocationType(const uint32_t &mod, const TrackerTopology *tTopo)
std::vector< dqm::reco::MonitorElement * > h2_ct_var3_m_
std::vector< dqm::reco::MonitorElement * > h2_t_w_m_
T z() const
Definition: PV3DBase.h:61
std::vector< dqm::reco::MonitorElement * > h2_t_var3_m_
LocalVector drift(const StripGeomDetUnit *, const MagneticField &, const SiStripLorentzAngle &)
Definition: ShallowTools.cc:36
const edm::ESGetToken< TrackerGeometry, TrackerDigiGeometryRecord > m_tkGeomTokenBR
LocalPoint toLocal(const GlobalPoint &gp) const
Provides dense hash map in place of DetId.
std::vector< dqm::reco::MonitorElement * > h2_ct_w_m_
SiStripLorentzAngleCalibrationHistograms iHists_
T x() const
Definition: PV3DBase.h:59
T y() const
Definition: PV3DBase.h:60
const std::string fieldAsString(const float &inputField)
std::vector< dqm::reco::MonitorElement * > h2_t_var2_m_
const std::string apvModeAsString(const SiStripLatency *latency)
const edm::ESGetToken< SiStripLatency, SiStripLatencyRcd > m_latencyTokenBR
GlobalPoint toGlobal(const Local2DPoint &lp) const
Conversion to the global R.F. from the R.F. of the GeomDet.
Definition: GeomDet.h:49
std::vector< dqm::reco::MonitorElement * > h2_ct_var2_m_
const Plane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:37
const PositionType & position() const
const edm::ESGetToken< TrackerTopology, TrackerTopologyRcd > m_topoEsTokenBR
const edm::ESGetToken< MagneticField, IdealMagneticFieldRecord > m_magFieldTokenBR

◆ fillDescriptions()

void SiStripLorentzAnglePCLMonitor::fillDescriptions ( edm::ConfigurationDescriptions descriptions)
static

Definition at line 491 of file SiStripLorentzAnglePCLMonitor.cc.

References edm::ConfigurationDescriptions::addWithDefaultLabel(), submitPVResolutionJobs::desc, ProducerED_cfi::InputTag, and AlCaHLTBitMon_QueryRunRegistry::string.

491  {
493  desc.add<std::string>("folder", "AlCaReco/SiStripLorentzAngle")->setComment("DQM folder to write into");
494  desc.add<bool>("saveHistoMods", false)->setComment("save module level hisotgrams. Warning! takes a lot of space!");
495  desc.add<edm::InputTag>("Tracks", edm::InputTag("SiStripCalCosmics"))->setComment("input track collection");
496  descriptions.addWithDefaultLabel(desc);
497 }
void addWithDefaultLabel(ParameterSetDescription const &psetDescription)

◆ moduleLocationType()

std::string SiStripLorentzAnglePCLMonitor::moduleLocationType ( const uint32_t &  mod,
const TrackerTopology tTopo 
)
private

Definition at line 199 of file SiStripLorentzAnglePCLMonitor.cc.

References nano_mu_digi_cff::layer, TrackerTopology::layer(), mod(), SiStripDetId::stereo(), AlCaHLTBitMon_QueryRunRegistry::string, SiStripDetId::subDetector(), SiStripDetId::TIB, and SiStripDetId::TOB.

Referenced by dqmBeginRun().

199  {
200  const SiStripDetId detid(mod);
201  std::string subdet = "";
202  unsigned int layer = 0;
203  if (detid.subDetector() == SiStripDetId::TIB) {
204  subdet = "TIB";
205  layer = tTopo->layer(mod);
206  }
207 
208  if (detid.subDetector() == SiStripDetId::TOB) {
209  subdet = "TOB";
210  layer = tTopo->layer(mod);
211  }
212 
213  std::string type = (detid.stereo() ? "s" : "a");
214  std::string d_l_t = Form("%s_L%d%s", subdet.c_str(), layer, type.c_str());
215 
216  if (layer == 0)
217  return subdet;
218  return d_l_t;
219 }
unsigned int layer(const DetId &id) const
static constexpr auto TOB
Definition: SiStripDetId.h:39
Detector identifier class for the strip tracker.
Definition: SiStripDetId.h:18
static constexpr auto TIB
Definition: SiStripDetId.h:37
T mod(const T &a, const T &b)
Definition: ecalDccMap.h:4

Member Data Documentation

◆ folder_

const std::string SiStripLorentzAnglePCLMonitor::folder_
private

Definition at line 84 of file SiStripLorentzAnglePCLMonitor.cc.

Referenced by bookHistograms().

◆ iHists_

SiStripLorentzAngleCalibrationHistograms SiStripLorentzAnglePCLMonitor::iHists_
private

Definition at line 76 of file SiStripLorentzAnglePCLMonitor.cc.

Referenced by analyze(), bookHistograms(), and dqmBeginRun().

◆ m_association_token

const edm::EDGetTokenT<TrajTrackAssociationCollection> SiStripLorentzAnglePCLMonitor::m_association_token
private

Definition at line 87 of file SiStripLorentzAnglePCLMonitor.cc.

Referenced by analyze().

◆ m_clusterInfo

SiStripClusterInfo SiStripLorentzAnglePCLMonitor::m_clusterInfo
private

Definition at line 75 of file SiStripLorentzAnglePCLMonitor.cc.

Referenced by analyze().

◆ m_hash

SiStripHashedDetId SiStripLorentzAnglePCLMonitor::m_hash
private

Definition at line 77 of file SiStripLorentzAnglePCLMonitor.cc.

Referenced by analyze(), bookHistograms(), and dqmBeginRun().

◆ m_latencyTokenBR

const edm::ESGetToken<SiStripLatency, SiStripLatencyRcd> SiStripLorentzAnglePCLMonitor::m_latencyTokenBR
private

Definition at line 90 of file SiStripLorentzAnglePCLMonitor.cc.

Referenced by dqmBeginRun().

◆ m_lorentzAngleTokenBR

const edm::ESGetToken<SiStripLorentzAngle, SiStripLorentzAngleDepRcd> SiStripLorentzAnglePCLMonitor::m_lorentzAngleTokenBR
private

Definition at line 94 of file SiStripLorentzAnglePCLMonitor.cc.

Referenced by dqmBeginRun().

◆ m_magFieldTokenBR

const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> SiStripLorentzAnglePCLMonitor::m_magFieldTokenBR
private

Definition at line 93 of file SiStripLorentzAnglePCLMonitor.cc.

Referenced by dqmBeginRun().

◆ m_tkGeomToken

const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> SiStripLorentzAnglePCLMonitor::m_tkGeomToken
private

Definition at line 88 of file SiStripLorentzAnglePCLMonitor.cc.

◆ m_tkGeomTokenBR

const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> SiStripLorentzAnglePCLMonitor::m_tkGeomTokenBR
private

Definition at line 92 of file SiStripLorentzAnglePCLMonitor.cc.

Referenced by dqmBeginRun().

◆ m_topoEsTokenBR

const edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> SiStripLorentzAnglePCLMonitor::m_topoEsTokenBR
private

Definition at line 91 of file SiStripLorentzAnglePCLMonitor.cc.

Referenced by dqmBeginRun().

◆ m_tracks_token

const edm::EDGetTokenT<edm::View<reco::Track> > SiStripLorentzAnglePCLMonitor::m_tracks_token
private

Definition at line 86 of file SiStripLorentzAnglePCLMonitor.cc.

Referenced by analyze().

◆ mismatchedBField_

bool SiStripLorentzAnglePCLMonitor::mismatchedBField_
private

Definition at line 82 of file SiStripLorentzAnglePCLMonitor.cc.

Referenced by analyze(), and dqmBeginRun().

◆ mismatchedLatency_

bool SiStripLorentzAnglePCLMonitor::mismatchedLatency_
private

Definition at line 83 of file SiStripLorentzAnglePCLMonitor.cc.

Referenced by analyze(), and dqmBeginRun().

◆ saveHistosMods_

const bool SiStripLorentzAnglePCLMonitor::saveHistosMods_
private

Definition at line 85 of file SiStripLorentzAnglePCLMonitor.cc.

Referenced by analyze(), and bookHistograms().

◆ teslaToInverseGeV_

constexpr float SiStripLorentzAnglePCLMonitor::teslaToInverseGeV_ = 2.99792458e-3f
staticprivate

Definition at line 80 of file SiStripLorentzAnglePCLMonitor.cc.

Referenced by dqmBeginRun().