CMS 3D CMS Logo

List of all members | Public Types | Public Member Functions | Static Public Member Functions | Private Attributes
SiPixelCompareTracks< T > Class Template Reference

#include <SiPixelCompareTracks.cc>

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

Public Types

using PixelTrackSoA = TracksHost< T >
 
- 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
 

Public Member Functions

void analyze (const edm::Event &iEvent, const edm::EventSetup &iSetup) override
 
template<typename U , typename V >
void analyzeSeparate (U tokenRef, V tokenTar, const edm::Event &iEvent)
 
void bookHistograms (DQMStore::IBooker &ibooker, edm::Run const &iRun, edm::EventSetup const &iSetup) override
 
 SiPixelCompareTracks (const edm::ParameterSet &)
 
 ~SiPixelCompareTracks () 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
 
virtual void dqmBeginRun (edm::Run const &, edm::EventSetup const &)
 
 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 Attributes

const float dr2cut_
 
MonitorElementhCharge_
 
MonitorElementhchi2_
 
MonitorElementhChi2VsEta_
 
MonitorElementhChi2VsPhi_
 
MonitorElementhCurvature_
 
MonitorElementhCurvdiffMatched_
 
MonitorElementhDeltaNLooseAndAboveTracks_
 
MonitorElementhDeltaNLooseAndAboveTracks_matched_
 
MonitorElementhDeltaNTracks_
 
MonitorElementheta_
 
MonitorElementhetadiffMatched_
 
MonitorElementhnHits_
 
MonitorElementhnHitsVsEta_
 
MonitorElementhnHitsVsPhi_
 
MonitorElementhnLayers_
 
MonitorElementhnLayersVsEta_
 
MonitorElementhnLayersVsPhi_
 
MonitorElementhnLooseAndAboveTracks_
 
MonitorElementhnLooseAndAboveTracks_matched_
 
MonitorElementhnTracks_
 
MonitorElementhphi_
 
MonitorElementhphi_z_tkAllRef_
 
MonitorElementhphi_z_tkAllRefMatched_
 
MonitorElementhphidiffMatched_
 
MonitorElementhpt_
 
MonitorElementhpt_eta_tkAllRef_
 
MonitorElementhpt_eta_tkAllRefMatched_
 
MonitorElementhptdiffMatched_
 
MonitorElementhptLogLog_
 
MonitorElementhquality_
 
MonitorElementhtip_
 
MonitorElementhtipdiffMatched_
 
MonitorElementhz_
 
MonitorElementhzdiffMatched_
 
const pixelTrack::Quality minQuality_
 
const edm::EDGetTokenT< PixelTrackSoAtokenSoATrackReference_
 
const edm::EDGetTokenT< PixelTrackSoAtokenSoATrackTarget_
 
const std::string topFolderName_
 
const bool useQualityCut_
 

Additional Inherited Members

- 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

template<typename T>
class SiPixelCompareTracks< T >

Definition at line 73 of file SiPixelCompareTracks.cc.

Member Typedef Documentation

◆ PixelTrackSoA

template<typename T >
using SiPixelCompareTracks< T >::PixelTrackSoA = TracksHost<T>

Definition at line 75 of file SiPixelCompareTracks.cc.

Constructor & Destructor Documentation

◆ SiPixelCompareTracks()

template<typename T >
SiPixelCompareTracks< T >::SiPixelCompareTracks ( const edm::ParameterSet iConfig)
explicit

Definition at line 139 of file SiPixelCompareTracks.cc.

140  : tokenSoATrackReference_(consumes<PixelTrackSoA>(iConfig.getParameter<edm::InputTag>("pixelTrackReferenceSoA"))),
141  tokenSoATrackTarget_(consumes<PixelTrackSoA>(iConfig.getParameter<edm::InputTag>("pixelTrackTargetSoA"))),
142  topFolderName_(iConfig.getParameter<std::string>("topFolderName")),
143  useQualityCut_(iConfig.getParameter<bool>("useQualityCut")),
145  dr2cut_(iConfig.getParameter<double>("deltaR2cut")) {}
T getParameter(std::string const &) const
Definition: ParameterSet.h:307
Quality qualityByName(std::string const &name)
const edm::EDGetTokenT< PixelTrackSoA > tokenSoATrackReference_
const pixelTrack::Quality minQuality_
const std::string topFolderName_
const edm::EDGetTokenT< PixelTrackSoA > tokenSoATrackTarget_

◆ ~SiPixelCompareTracks()

template<typename T >
SiPixelCompareTracks< T >::~SiPixelCompareTracks ( )
overridedefault

Member Function Documentation

◆ analyze()

template<typename T >
void SiPixelCompareTracks< T >::analyze ( const edm::Event iEvent,
const edm::EventSetup iSetup 
)
overridevirtual

Reimplemented from DQMEDAnalyzer.

Definition at line 283 of file SiPixelCompareTracks.cc.

References iEvent.

283  {
284  // The default use case is to use vertices from Alpaka reconstructed on CPU and GPU;
285  // The function is left templated if any other cases need to be added
287 }
const edm::EDGetTokenT< PixelTrackSoA > tokenSoATrackReference_
void analyzeSeparate(U tokenRef, V tokenTar, const edm::Event &iEvent)
int iEvent
Definition: GenABIO.cc:224
const edm::EDGetTokenT< PixelTrackSoA > tokenSoATrackTarget_

◆ analyzeSeparate()

template<typename T >
template<typename U , typename V >
void SiPixelCompareTracks< T >::analyzeSeparate ( tokenRef,
tokenTar,
const edm::Event iEvent 
)

Definition at line 149 of file SiPixelCompareTracks.cc.

References reco::charge(), reco::deltaPhi(), reco::deltaR2(), PV3DBase< T, PVType, FrameType >::eta(), iEvent, ALPAKA_ACCELERATOR_NAMESPACE::vertexFinder::it, SiStripPI::min, nHits, notFound, or, MillePedeFileConverter_cfg::out, reco::phi(), reco::tip(), and reco::zip().

149  {
150  using helper = TracksUtilities<T>;
151 
152  const auto& tsoaHandleRef = iEvent.getHandle(tokenRef);
153  const auto& tsoaHandleTar = iEvent.getHandle(tokenTar);
154 
155  if (not tsoaHandleRef or not tsoaHandleTar) {
156  edm::LogWarning out("SiPixelCompareTracks");
157  if (not tsoaHandleRef) {
158  out << "reference tracks not found; ";
159  }
160  if (not tsoaHandleTar) {
161  out << "target tracks not found; ";
162  }
163  out << "the comparison will not run.";
164  return;
165  }
166 
167  auto const& tsoaRef = *tsoaHandleRef;
168  auto const& tsoaTar = *tsoaHandleTar;
169 
170  auto maxTracksRef = tsoaRef.view().metadata().size(); //this should be same for both?
171  auto maxTracksTar = tsoaTar.view().metadata().size(); //this should be same for both?
172 
173  auto const* qualityRef = tsoaRef.view().quality();
174  auto const* qualityTar = tsoaTar.view().quality();
175 
176  int32_t nTracksRef = 0;
177  int32_t nTracksTar = 0;
178  int32_t nLooseAndAboveTracksRef = 0;
179  int32_t nLooseAndAboveTracksRef_matchedTar = 0;
180  int32_t nLooseAndAboveTracksTar = 0;
181 
182  //Loop over Tar tracks and store the indices of the loose tracks. Whats happens if useQualityCut_ is false?
183  std::vector<int32_t> looseTrkidxTar;
184  for (int32_t jt = 0; jt < maxTracksTar; ++jt) {
185  if (helper::nHits(tsoaTar.view(), jt) == 0)
186  break; // this is a guard
187  if (!(tsoaTar.view()[jt].pt() > 0.))
188  continue;
189  nTracksTar++;
190  if (useQualityCut_ && qualityTar[jt] < minQuality_)
191  continue;
192  nLooseAndAboveTracksTar++;
193  looseTrkidxTar.emplace_back(jt);
194  }
195 
196  //Now loop over Ref tracks//nested loop for loose gPU tracks
197  for (int32_t it = 0; it < maxTracksRef; ++it) {
198  int nHitsRef = helper::nHits(tsoaRef.view(), it);
199 
200  if (nHitsRef == 0)
201  break; // this is a guard
202 
203  float ptRef = tsoaRef.view()[it].pt();
204  float etaRef = tsoaRef.view()[it].eta();
205  float phiRef = reco::phi(tsoaRef.view(), it);
206  float zipRef = reco::zip(tsoaRef.view(), it);
207  float tipRef = reco::tip(tsoaRef.view(), it);
208  auto qRef = reco::charge(tsoaRef.view(), it);
209 
210  if (!(ptRef > 0.))
211  continue;
212  nTracksRef++;
213  if (useQualityCut_ && qualityRef[it] < minQuality_)
214  continue;
215  nLooseAndAboveTracksRef++;
216  //Now loop over loose Tar trk and find the closest in DeltaR//do we need pt cut?
217  const int32_t notFound = -1;
218  int32_t closestTkidx = notFound;
219  float mindr2 = dr2cut_;
220 
221  for (auto gid : looseTrkidxTar) {
222  float etaTar = tsoaTar.view()[gid].eta();
223  float phiTar = reco::phi(tsoaTar.view(), gid);
224  float dr2 = reco::deltaR2(etaRef, phiRef, etaTar, phiTar);
225  if (dr2 > dr2cut_)
226  continue; // this is arbitrary
227  if (mindr2 > dr2) {
228  mindr2 = dr2;
229  closestTkidx = gid;
230  }
231  }
232 
233  hpt_eta_tkAllRef_->Fill(etaRef, ptRef); //all Ref tk
234  hphi_z_tkAllRef_->Fill(phiRef, zipRef);
235  if (closestTkidx == notFound)
236  continue;
237  nLooseAndAboveTracksRef_matchedTar++;
238 
239  hchi2_->Fill(tsoaRef.view()[it].chi2(), tsoaTar.view()[closestTkidx].chi2());
240  hCharge_->Fill(qRef, reco::charge(tsoaTar.view(), closestTkidx));
241  hnHits_->Fill(helper::nHits(tsoaRef.view(), it), helper::nHits(tsoaTar.view(), closestTkidx));
242  hnLayers_->Fill(tsoaRef.view()[it].nLayers(), tsoaTar.view()[closestTkidx].nLayers());
243  hpt_->Fill(ptRef, tsoaTar.view()[closestTkidx].pt());
244  hCurvature_->Fill(qRef / ptRef, reco::charge(tsoaTar.view(), closestTkidx) / tsoaTar.view()[closestTkidx].pt());
245  hptLogLog_->Fill(ptRef, tsoaTar.view()[closestTkidx].pt());
246  heta_->Fill(etaRef, tsoaTar.view()[closestTkidx].eta());
247  hphi_->Fill(phiRef, reco::phi(tsoaTar.view(), closestTkidx));
248  hz_->Fill(zipRef, reco::zip(tsoaTar.view(), closestTkidx));
249  htip_->Fill(tipRef, reco::tip(tsoaTar.view(), closestTkidx));
250  hptdiffMatched_->Fill(ptRef - tsoaTar.view()[closestTkidx].pt());
251  hCurvdiffMatched_->Fill(qRef / ptRef -
252  (reco::charge(tsoaTar.view(), closestTkidx) / tsoaTar.view()[closestTkidx].pt()));
253  hetadiffMatched_->Fill(etaRef - tsoaTar.view()[closestTkidx].eta());
254  hphidiffMatched_->Fill(reco::deltaPhi(phiRef, reco::phi(tsoaTar.view(), closestTkidx)));
255  hzdiffMatched_->Fill(zipRef - reco::zip(tsoaTar.view(), closestTkidx));
256  htipdiffMatched_->Fill(tipRef - reco::tip(tsoaTar.view(), closestTkidx));
257  hpt_eta_tkAllRefMatched_->Fill(etaRef, tsoaRef.view()[it].pt()); //matched to gpu
258  hphi_z_tkAllRefMatched_->Fill(etaRef, zipRef);
259  }
260 
261  // Define a lambda function for filling the histograms
262  auto fillHistogram = [](auto& histogram, auto xValue, auto yValue) { histogram->Fill(xValue, yValue); };
263 
264  // Define a lambda for filling delta histograms
265  auto fillDeltaHistogram = [](auto& histogram, int cpuValue, int gpuValue) {
266  histogram->Fill(std::min(cpuValue, 1000), std::clamp(gpuValue - cpuValue, -100, 100));
267  };
268 
269  // Fill the histograms
270  fillHistogram(hnTracks_, nTracksRef, nTracksTar);
271  fillHistogram(hnLooseAndAboveTracks_, nLooseAndAboveTracksRef, nLooseAndAboveTracksTar);
272  fillHistogram(hnLooseAndAboveTracks_matched_, nLooseAndAboveTracksRef, nLooseAndAboveTracksRef_matchedTar);
273 
274  fillDeltaHistogram(hDeltaNTracks_, nTracksRef, nTracksTar);
275  fillDeltaHistogram(hDeltaNLooseAndAboveTracks_, nLooseAndAboveTracksRef, nLooseAndAboveTracksTar);
276  fillDeltaHistogram(hDeltaNLooseAndAboveTracks_matched_, nLooseAndAboveTracksRef, nLooseAndAboveTracksRef_matchedTar);
277 }
constexpr double deltaPhi(double phi1, double phi2)
Definition: deltaPhi.h:26
MonitorElement * hDeltaNTracks_
MonitorElement * hnLooseAndAboveTracks_
Definition: helper.py:1
MonitorElement * hptdiffMatched_
MonitorElement * hphi_z_tkAllRef_
ALPAKA_FN_HOST_ACC ALPAKA_FN_INLINE constexpr float zip(ConstView const &tracks, int32_t i)
Definition: TracksSoA.h:90
MonitorElement * hphi_z_tkAllRefMatched_
MonitorElement * hCurvdiffMatched_
MonitorElement * hnTracks_
ALPAKA_FN_HOST_ACC ALPAKA_FN_INLINE constexpr float charge(ConstView const &tracks, int32_t i)
Definition: TracksSoA.h:73
MonitorElement * hpt_eta_tkAllRefMatched_
MonitorElement * hnLooseAndAboveTracks_matched_
MonitorElement * hphidiffMatched_
void Fill(long long x)
ALPAKA_FN_HOST_ACC ALPAKA_FN_INLINE constexpr float tip(ConstView const &tracks, int32_t i)
Definition: TracksSoA.h:85
int iEvent
Definition: GenABIO.cc:224
MonitorElement * hpt_eta_tkAllRef_
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e< void, edm::EventID const &, edm::Timestamp const & > We also list in braces which AR_WATCH_USING_METHOD_ is used for those or
Definition: Activities.doc:12
ALPAKA_FN_HOST_ACC ALPAKA_FN_INLINE constexpr float phi(ConstView const &tracks, int32_t i)
Definition: TracksSoA.h:80
MonitorElement * hDeltaNLooseAndAboveTracks_
MonitorElement * hDeltaNLooseAndAboveTracks_matched_
MonitorElement * hzdiffMatched_
MonitorElement * hetadiffMatched_
constexpr auto deltaR2(const T1 &t1, const T2 &t2) -> decltype(t1.eta())
Definition: deltaR.h:16
MonitorElement * hCurvature_
static const GlobalPoint notFound(0, 0, 0)
const pixelTrack::Quality minQuality_
MonitorElement * htipdiffMatched_
TupleMultiplicity< TrackerTraits > const *__restrict__ uint32_t nHits

◆ bookHistograms()

template<typename T >
void SiPixelCompareTracks< T >::bookHistograms ( DQMStore::IBooker ibooker,
edm::Run const &  iRun,
edm::EventSetup const &  iSetup 
)
overridevirtual

Implements DQMEDAnalyzer.

Definition at line 293 of file SiPixelCompareTracks.cc.

References dqm::implementation::IBooker::book1D(), dqm::implementation::IBooker::book2I(), dqm::implementation::NavigatorBase::cd(), ALPAKA_ACCELERATOR_NAMESPACE::brokenline::constexpr(), M_PI, mergeVDriftHistosByStation::name, dqm::implementation::NavigatorBase::setCurrentFolder(), AlCaHLTBitMon_QueryRunRegistry::string, runGCPTkAlMap::title, multiplicitycorr_cfi::xBins, multiplicitycorr_cfi::xMax, photonAnalyzer_cfi::xMin, multiplicitycorr_cfi::yBins, multiplicitycorr_cfi::yMax, and photonAnalyzer_cfi::yMin.

295  {
296  iBook.cd();
297  iBook.setCurrentFolder(topFolderName_);
298 
299  // Define a helper function for booking histograms
300  std::string toRep = "Number of tracks";
301  auto bookTracksTH2I = [&](const std::string& name,
302  const std::string& title,
303  int xBins,
304  double xMin,
305  double xMax,
306  int yBins,
307  double yMin,
308  double yMax) {
309  return iBook.book2I(name, fmt::sprintf(title, toRep), xBins, xMin, xMax, yBins, yMin, yMax);
310  };
311 
312  // Define common parameters for different histogram types
313  constexpr int xBins = 501;
314  constexpr double xMin = -0.5;
315  constexpr double xMax = 1001.5;
316 
317  constexpr int dXBins = 1001;
318  constexpr double dXMin = -0.5;
319  constexpr double dXMax = 1000.5;
320 
321  constexpr int dYBins = 201;
322  constexpr double dYMin = -100.5;
323  constexpr double dYMax = 100.5;
324 
325  // FIXME: all the 2D correlation plots are quite heavy in terms of memory consumption, so a as soon as DQM supports THnSparse
326  // these should be moved to a less resource consuming format
327 
328  // Book histograms using the helper function
329  // clang-format off
330  hnTracks_ = bookTracksTH2I("nTracks", "%s per event; Reference; Target", xBins, xMin, xMax, xBins, xMin, xMax);
331  hnLooseAndAboveTracks_ = bookTracksTH2I("nLooseAndAboveTracks", "%s (quality #geq loose) per event; Reference; Target", xBins, xMin, xMax, xBins, xMin, xMax);
332  hnLooseAndAboveTracks_matched_ = bookTracksTH2I("nLooseAndAboveTracks_matched", "%s (quality #geq loose) per event; Reference; Target", xBins, xMin, xMax, xBins, xMin, xMax);
333 
334  hDeltaNTracks_ = bookTracksTH2I("deltaNTracks", "%s per event; Reference; Target - Reference", dXBins, dXMin, dXMax, dYBins, dYMin, dYMax);
335  hDeltaNLooseAndAboveTracks_ = bookTracksTH2I("deltaNLooseAndAboveTracks", "%s (quality #geq loose) per event; Reference; Target - Reference", dXBins, dXMin, dXMax, dYBins, dYMin, dYMax);
336  hDeltaNLooseAndAboveTracks_matched_ = bookTracksTH2I("deltaNLooseAndAboveTracks_matched", "%s (quality #geq loose) per event; Reference; Target - Reference", dXBins, dXMin, dXMax, dYBins, dYMin, dYMax);
337 
338  toRep = "Number of all RecHits per track (quality #geq loose)";
339  hnHits_ = iBook.book2I("nRecHits", fmt::sprintf("%s;Reference;Target",toRep), 15, -0.5, 14.5, 15, -0.5, 14.5);
340 
341  toRep = "Number of all layers per track (quality #geq loose)";
342  hnLayers_ = iBook.book2I("nLayers", fmt::sprintf("%s;Reference;Target",toRep), 15, -0.5, 14.5, 15, -0.5, 14.5);
343 
344  toRep = "Track (quality #geq loose) #chi^{2}/ndof";
345  hchi2_ = iBook.book2I("nChi2ndof", fmt::sprintf("%s;Reference;Target",toRep), 40, 0., 20., 40, 0., 20.);
346 
347  toRep = "Track (quality #geq loose) charge";
348  hCharge_ = iBook.book2I("charge",fmt::sprintf("%s;Reference;Target",toRep),3, -1.5, 1.5, 3, -1.5, 1.5);
349 
350  hpt_ = iBook.book2I("pt", "Track (quality #geq loose) p_{T} [GeV];Reference;Target", 200, 0., 200., 200, 0., 200.);
351  hCurvature_ = iBook.book2I("curvature", "Track (quality #geq loose) q/p_{T} [GeV^{-1}];Reference;Target", 60,- 3., 3., 60, -3., 3. );
352  hptLogLog_ = make2DIfLog(iBook, true, true, "ptLogLog", "Track (quality #geq loose) p_{T} [GeV];Reference;Target", 200, log10(0.5), log10(200.), 200, log10(0.5), log10(200.));
353  heta_ = iBook.book2I("eta", "Track (quality #geq loose) #eta;Reference;Target", 30, -3., 3., 30, -3., 3.);
354  hphi_ = iBook.book2I("phi", "Track (quality #geq loose) #phi;Reference;Target", 30, -M_PI, M_PI, 30, -M_PI, M_PI);
355  hz_ = iBook.book2I("z", "Track (quality #geq loose) z [cm];Reference;Target", 30, -30., 30., 30, -30., 30.);
356  htip_ = iBook.book2I("tip", "Track (quality #geq loose) TIP [cm];Reference;Target", 100, -0.5, 0.5, 100, -0.5, 0.5);
357 
358  //1D difference plots
359  hptdiffMatched_ = iBook.book1D("ptdiffmatched", " p_{T} diff [GeV] between matched tracks; #Delta p_{T} [GeV]", 61, -30.5, 30.5);
360  hCurvdiffMatched_ = iBook.book1D("curvdiffmatched", "q/p_{T} diff [GeV^{-1}] between matched tracks; #Delta q/p_{T} [GeV^{-1}]", 61, -3.05, 3.05);
361  hetadiffMatched_ = iBook.book1D("etadiffmatched", " #eta diff between matched tracks; #Delta #eta", 161, -0.045 ,0.045);
362  hphidiffMatched_ = iBook.book1D("phidiffmatched", " #phi diff between matched tracks; #Delta #phi", 161, -0.045 ,0.045);
363  hzdiffMatched_ = iBook.book1D("zdiffmatched", " z diff between matched tracks; #Delta z [cm]", 301, -1.55, 1.55);
364  htipdiffMatched_ = iBook.book1D("tipdiffmatched", " TIP diff between matched tracks; #Delta TIP [cm]", 301, -1.55, 1.55);
365  //2D plots for eff
366  hpt_eta_tkAllRef_ = iBook.book2I("ptetatrkAllReference", "Track (quality #geq loose) on Reference; #eta; p_{T} [GeV];", 30, -M_PI, M_PI, 200, 0., 200.);
367  hpt_eta_tkAllRefMatched_ = iBook.book2I("ptetatrkAllReferencematched", "Track (quality #geq loose) on Reference matched to Target track; #eta; p_{T} [GeV];", 30, -M_PI, M_PI, 200, 0., 200.);
368 
369  hphi_z_tkAllRef_ = iBook.book2I("phiztrkAllReference", "Track (quality #geq loose) on Reference; #phi; z [cm];", 30, -M_PI, M_PI, 30, -30., 30.);
370  hphi_z_tkAllRefMatched_ = iBook.book2I("phiztrkAllReferencematched", "Track (quality #geq loose) on Reference; #phi; z [cm];", 30, -M_PI, M_PI, 30, -30., 30.);
371 
372 }
MonitorElement * hDeltaNTracks_
MonitorElement * hnLooseAndAboveTracks_
MonitorElement * hptdiffMatched_
MonitorElement * hphi_z_tkAllRef_
MonitorElement * hphi_z_tkAllRefMatched_
MonitorElement * hCurvdiffMatched_
MonitorElement * hnTracks_
MonitorElement * hpt_eta_tkAllRefMatched_
MonitorElement * hnLooseAndAboveTracks_matched_
MonitorElement * hphidiffMatched_
MonitorElement * hpt_eta_tkAllRef_
MonitorElement * hDeltaNLooseAndAboveTracks_
MonitorElement * hDeltaNLooseAndAboveTracks_matched_
MonitorElement * hzdiffMatched_
#define M_PI
MonitorElement * hetadiffMatched_
MonitorElement * hCurvature_
MonitorElement * htipdiffMatched_
const std::string topFolderName_

◆ fillDescriptions()

template<typename T >
void SiPixelCompareTracks< T >::fillDescriptions ( edm::ConfigurationDescriptions descriptions)
static

Definition at line 375 of file SiPixelCompareTracks.cc.

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

375  {
376  // monitorpixelTrackSoA
378  desc.add<edm::InputTag>("pixelTrackReferenceSoA", edm::InputTag("pixelTracksAlpakaSerial"));
379  desc.add<edm::InputTag>("pixelTrackTargetSoA", edm::InputTag("pixelTracksAlpaka"));
380  desc.add<std::string>("topFolderName", "SiPixelHeterogeneous/PixelTrackCompareDeviceVSHost");
381  desc.add<bool>("useQualityCut", true);
382  desc.add<std::string>("minQuality", "loose");
383  desc.add<double>("deltaR2cut", 0.02 * 0.02)->setComment("deltaR2 cut between track on device and host");
384  descriptions.addWithDefaultLabel(desc);
385 }
void addWithDefaultLabel(ParameterSetDescription const &psetDescription)

Member Data Documentation

◆ dr2cut_

template<typename T >
const float SiPixelCompareTracks< T >::dr2cut_
private

Definition at line 94 of file SiPixelCompareTracks.cc.

◆ hCharge_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hCharge_
private

Definition at line 107 of file SiPixelCompareTracks.cc.

◆ hchi2_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hchi2_
private

Definition at line 108 of file SiPixelCompareTracks.cc.

◆ hChi2VsEta_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hChi2VsEta_
private

Definition at line 110 of file SiPixelCompareTracks.cc.

◆ hChi2VsPhi_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hChi2VsPhi_
private

Definition at line 109 of file SiPixelCompareTracks.cc.

◆ hCurvature_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hCurvature_
private

Definition at line 112 of file SiPixelCompareTracks.cc.

◆ hCurvdiffMatched_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hCurvdiffMatched_
private

Definition at line 121 of file SiPixelCompareTracks.cc.

◆ hDeltaNLooseAndAboveTracks_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hDeltaNLooseAndAboveTracks_
private

Definition at line 99 of file SiPixelCompareTracks.cc.

◆ hDeltaNLooseAndAboveTracks_matched_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hDeltaNLooseAndAboveTracks_matched_
private

Definition at line 100 of file SiPixelCompareTracks.cc.

◆ hDeltaNTracks_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hDeltaNTracks_
private

Definition at line 98 of file SiPixelCompareTracks.cc.

◆ heta_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::heta_
private

Definition at line 114 of file SiPixelCompareTracks.cc.

◆ hetadiffMatched_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hetadiffMatched_
private

Definition at line 122 of file SiPixelCompareTracks.cc.

◆ hnHits_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hnHits_
private

Definition at line 101 of file SiPixelCompareTracks.cc.

◆ hnHitsVsEta_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hnHitsVsEta_
private

Definition at line 103 of file SiPixelCompareTracks.cc.

◆ hnHitsVsPhi_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hnHitsVsPhi_
private

Definition at line 102 of file SiPixelCompareTracks.cc.

◆ hnLayers_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hnLayers_
private

Definition at line 104 of file SiPixelCompareTracks.cc.

◆ hnLayersVsEta_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hnLayersVsEta_
private

Definition at line 106 of file SiPixelCompareTracks.cc.

◆ hnLayersVsPhi_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hnLayersVsPhi_
private

Definition at line 105 of file SiPixelCompareTracks.cc.

◆ hnLooseAndAboveTracks_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hnLooseAndAboveTracks_
private

Definition at line 96 of file SiPixelCompareTracks.cc.

◆ hnLooseAndAboveTracks_matched_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hnLooseAndAboveTracks_matched_
private

Definition at line 97 of file SiPixelCompareTracks.cc.

◆ hnTracks_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hnTracks_
private

Definition at line 95 of file SiPixelCompareTracks.cc.

◆ hphi_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hphi_
private

Definition at line 115 of file SiPixelCompareTracks.cc.

◆ hphi_z_tkAllRef_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hphi_z_tkAllRef_
private

Definition at line 130 of file SiPixelCompareTracks.cc.

◆ hphi_z_tkAllRefMatched_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hphi_z_tkAllRefMatched_
private

Definition at line 131 of file SiPixelCompareTracks.cc.

◆ hphidiffMatched_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hphidiffMatched_
private

Definition at line 123 of file SiPixelCompareTracks.cc.

◆ hpt_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hpt_
private

Definition at line 111 of file SiPixelCompareTracks.cc.

◆ hpt_eta_tkAllRef_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hpt_eta_tkAllRef_
private

Definition at line 128 of file SiPixelCompareTracks.cc.

◆ hpt_eta_tkAllRefMatched_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hpt_eta_tkAllRefMatched_
private

Definition at line 129 of file SiPixelCompareTracks.cc.

◆ hptdiffMatched_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hptdiffMatched_
private

Definition at line 120 of file SiPixelCompareTracks.cc.

◆ hptLogLog_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hptLogLog_
private

Definition at line 113 of file SiPixelCompareTracks.cc.

◆ hquality_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hquality_
private

Definition at line 118 of file SiPixelCompareTracks.cc.

◆ htip_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::htip_
private

Definition at line 117 of file SiPixelCompareTracks.cc.

◆ htipdiffMatched_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::htipdiffMatched_
private

Definition at line 125 of file SiPixelCompareTracks.cc.

◆ hz_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hz_
private

Definition at line 116 of file SiPixelCompareTracks.cc.

◆ hzdiffMatched_

template<typename T >
MonitorElement* SiPixelCompareTracks< T >::hzdiffMatched_
private

Definition at line 124 of file SiPixelCompareTracks.cc.

◆ minQuality_

template<typename T >
const pixelTrack::Quality SiPixelCompareTracks< T >::minQuality_
private

Definition at line 93 of file SiPixelCompareTracks.cc.

◆ tokenSoATrackReference_

template<typename T >
const edm::EDGetTokenT<PixelTrackSoA> SiPixelCompareTracks< T >::tokenSoATrackReference_
private

Definition at line 89 of file SiPixelCompareTracks.cc.

◆ tokenSoATrackTarget_

template<typename T >
const edm::EDGetTokenT<PixelTrackSoA> SiPixelCompareTracks< T >::tokenSoATrackTarget_
private

Definition at line 90 of file SiPixelCompareTracks.cc.

◆ topFolderName_

template<typename T >
const std::string SiPixelCompareTracks< T >::topFolderName_
private

Definition at line 91 of file SiPixelCompareTracks.cc.

◆ useQualityCut_

template<typename T >
const bool SiPixelCompareTracks< T >::useQualityCut_
private

Definition at line 92 of file SiPixelCompareTracks.cc.