CMS 3D CMS Logo

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

#include <DTSegment4DQuality.h>

Inheritance diagram for DTSegment4DQuality:
DQMGlobalEDAnalyzer< dtsegment4d::Histograms > DQMGlobalEDAnalyzerBase< dtsegment4d::Histograms, Args... > edm::global::EDProducer< edm::RunCache< dtsegment4d::Histograms >, edm::EndRunProducer, edm::Accumulator, Args... > edm::global::EDProducerBase edm::ProducerBase edm::EDConsumerBase edm::ProductRegistryHelper

Public Member Functions

 DTSegment4DQuality (const edm::ParameterSet &pset)
 Constructor. More...
 
- Public Member Functions inherited from DQMGlobalEDAnalyzer< dtsegment4d::Histograms >
virtual void dqmEndRun (edm::Run const &, edm::EventSetup const &, dtsegment4d::Histograms const &) const
 
void globalEndRunProduce (edm::Run &run, edm::EventSetup const &setup) const final
 
- Public Member Functions inherited from DQMGlobalEDAnalyzerBase< dtsegment4d::Histograms, Args... >
void accumulate (edm::StreamID id, edm::Event const &event, edm::EventSetup const &setup) const final
 
virtual void dqmBeginRun (edm::Run const &, edm::EventSetup const &, dtsegment4d::Histograms &) const
 
 DQMGlobalEDAnalyzerBase ()
 
std::shared_ptr< dtsegment4d::HistogramsglobalBeginRun (edm::Run const &run, edm::EventSetup const &setup) const final
 
void globalEndRun (edm::Run const &, edm::EventSetup const &) const final
 
- Public Member Functions inherited from edm::global::EDProducer< edm::RunCache< dtsegment4d::Histograms >, edm::EndRunProducer, edm::Accumulator, Args... >
 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
 
EDProduceroperator= (const EDProducer &)=delete
 
bool wantsGlobalLuminosityBlocks () const final
 
bool wantsGlobalRuns () const final
 
bool wantsInputProcessBlocks () const final
 
bool wantsProcessBlocks () const final
 
bool wantsStreamLuminosityBlocks () const final
 
bool wantsStreamRuns () const final
 
- Public Member Functions inherited from edm::global::EDProducerBase
 EDProducerBase ()
 
ModuleDescription const & moduleDescription () const
 
 ~EDProducerBase () override
 
- Public Member Functions inherited from edm::ProducerBase
void callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)
 
std::vector< edm::ProductResolverIndex > const & indiciesForPutProducts (BranchType iBranchType) const
 
 ProducerBase ()
 
std::vector< edm::ProductResolverIndex > const & putTokenIndexToProductResolverIndex () const
 
std::vector< bool > const & recordProvenanceList () const
 
void registerProducts (ProducerBase *, ProductRegistry *, ModuleDescription const &)
 
std::function< void(BranchDescription const &)> registrationCallback () const
 used by the fwk to register list of products More...
 
void resolvePutIndicies (BranchType iBranchType, ModuleToResolverIndicies const &iIndicies, std::string const &moduleLabel)
 
TypeLabelList const & typeLabelList () const
 used by the fwk to register the list of products of this module More...
 
 ~ProducerBase () noexcept(false) override
 
- Public Member Functions inherited from edm::EDConsumerBase
std::vector< ConsumesInfoconsumesInfo () const
 
void convertCurrentProcessAlias (std::string const &processName)
 Convert "@currentProcess" in InputTag process names to the actual current process name. More...
 
 EDConsumerBase ()
 
 EDConsumerBase (EDConsumerBase const &)=delete
 
 EDConsumerBase (EDConsumerBase &&)=default
 
ESProxyIndex const * esGetTokenIndices (edm::Transition iTrans) const
 
std::vector< ESProxyIndex > const & esGetTokenIndicesVector (edm::Transition iTrans) const
 
std::vector< ESRecordIndex > const & esGetTokenRecordIndicesVector (edm::Transition iTrans) const
 
ProductResolverIndexAndSkipBit indexFrom (EDGetToken, BranchType, TypeID const &) const
 
void itemsMayGet (BranchType, std::vector< ProductResolverIndexAndSkipBit > &) const
 
void itemsToGet (BranchType, std::vector< ProductResolverIndexAndSkipBit > &) const
 
std::vector< ProductResolverIndexAndSkipBit > const & itemsToGetFrom (BranchType iType) const
 
void labelsForToken (EDGetToken iToken, Labels &oLabels) const
 
void modulesWhoseProductsAreConsumed (std::array< std::vector< ModuleDescription const *> *, NumBranchTypes > &modulesAll, std::vector< ModuleProcessName > &modulesInPreviousProcesses, ProductRegistry const &preg, std::map< std::string, ModuleDescription const *> const &labelsToDesc, std::string const &processName) const
 
EDConsumerBase const & operator= (EDConsumerBase const &)=delete
 
EDConsumerBaseoperator= (EDConsumerBase &&)=default
 
bool registeredToConsume (ProductResolverIndex, bool, BranchType) const
 
bool registeredToConsumeMany (TypeID const &, BranchType) const
 
void selectInputProcessBlocks (ProductRegistry const &productRegistry, ProcessBlockHelperBase const &processBlockHelperBase)
 
ProductResolverIndexAndSkipBit uncheckedIndexFrom (EDGetToken) const
 
void updateLookup (BranchType iBranchType, ProductResolverIndexHelper const &, bool iPrefetchMayGet)
 
void updateLookup (eventsetup::ESRecordsToProxyIndices const &)
 
virtual ~EDConsumerBase () noexcept(false)
 

Private Member Functions

void bookHistograms (DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &, dtsegment4d::Histograms &) const override
 Book the DQM plots. More...
 
void dqmAnalyze (edm::Event const &, edm::EventSetup const &, dtsegment4d::Histograms const &) const override
 Perform the real analysis. More...
 

Private Attributes

bool debug_
 
bool doall_
 
bool local_
 
edm::ESGetToken< DTGeometry, MuonGeometryRecordmuonGeomToken_
 
edm::InputTag segment4DLabel_
 
edm::EDGetTokenT< DTRecSegment4DCollectionsegment4DToken_
 
double sigmaResAlpha_
 
double sigmaResBeta_
 
double sigmaResX_
 
double sigmaResY_
 
edm::InputTag simHitLabel_
 
edm::EDGetTokenT< edm::PSimHitContainersimHitToken_
 

Additional Inherited Members

- Public Types inherited from DQMGlobalEDAnalyzerBase< dtsegment4d::Histograms, Args... >
typedef dqm::reco::DQMStore DQMStore
 
typedef dqm::reco::MonitorElement MonitorElement
 
- Public Types inherited from edm::global::EDProducerBase
typedef EDProducerBase ModuleType
 
- Public Types inherited from edm::ProducerBase
template<typename T >
using BranchAliasSetterT = ProductRegistryHelper::BranchAliasSetterT< T >
 
using ModuleToResolverIndicies = std::unordered_multimap< std::string, std::tuple< edm::TypeID const *, const char *, edm::ProductResolverIndex > >
 
typedef ProductRegistryHelper::TypeLabelList TypeLabelList
 
- Public Types inherited from edm::EDConsumerBase
typedef ProductLabels Labels
 
- Static Public Member Functions inherited from edm::global::EDProducerBase
static const std::string & baseType ()
 
static void fillDescriptions (ConfigurationDescriptions &descriptions)
 
static void prevalidate (ConfigurationDescriptions &descriptions)
 
- Protected Member Functions inherited from DQMGlobalEDAnalyzerBase< dtsegment4d::Histograms, Args... >
uint64_t meId (edm::Run const &run) const
 
- Protected Member Functions inherited from edm::ProducerBase
template<Transition Tr = Transition::Event>
auto produces (std::string instanceName) noexcept
 declare what type of product will make and with which optional label More...
 
template<Transition B>
BranchAliasSetter produces (const TypeID &id, std::string instanceName=std::string(), bool recordProvenance=true)
 
template<BranchType B>
BranchAliasSetter produces (const TypeID &id, std::string instanceName=std::string(), bool recordProvenance=true)
 
BranchAliasSetter produces (const TypeID &id, std::string instanceName=std::string(), bool recordProvenance=true)
 
template<typename ProductType , Transition B>
BranchAliasSetterT< ProductType > produces (std::string instanceName)
 
template<class ProductType >
BranchAliasSetterT< ProductType > produces ()
 
template<typename ProductType , BranchType B>
BranchAliasSetterT< ProductType > produces (std::string instanceName)
 
template<typename ProductType , BranchType B>
BranchAliasSetterT< ProductType > produces ()
 
template<class ProductType >
BranchAliasSetterT< ProductType > produces (std::string instanceName)
 
template<typename ProductType , Transition B>
BranchAliasSetterT< ProductType > produces ()
 
template<Transition Tr = Transition::Event>
auto produces () noexcept
 
ProducesCollector producesCollector ()
 
- Protected Member Functions inherited from edm::EDConsumerBase
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > consumes (edm::InputTag const &tag)
 
template<BranchType B = InEvent>
EDConsumerBaseAdaptor< Bconsumes (edm::InputTag tag) noexcept
 
EDGetToken consumes (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken consumes (TypeToGet const &id, edm::InputTag const &tag)
 
ConsumesCollector consumesCollector ()
 Use a ConsumesCollector to gather consumes information from helper functions. More...
 
template<typename ProductType , BranchType B = InEvent>
void consumesMany ()
 
void consumesMany (const TypeToGet &id)
 
template<BranchType B>
void consumesMany (const TypeToGet &id)
 
template<typename ESProduct , typename ESRecord , Transition Tr = Transition::Event>
auto esConsumes ()
 
template<typename ESProduct , typename ESRecord , Transition Tr = Transition::Event>
auto esConsumes (ESInputTag const &tag)
 
template<Transition Tr = Transition::Event>
constexpr auto esConsumes ()
 
template<Transition Tr = Transition::Event>
auto esConsumes (ESInputTag tag)
 
template<Transition Tr = Transition::Event>
ESGetTokenGeneric esConsumes (eventsetup::EventSetupRecordKey const &iRecord, eventsetup::DataKey const &iKey)
 Used with EventSetupRecord::doGet. More...
 
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > mayConsume (edm::InputTag const &tag)
 
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 
void resetItemsToGetFrom (BranchType iType)
 
- Protected Attributes inherited from DQMGlobalEDAnalyzerBase< dtsegment4d::Histograms, Args... >
DQMStoredqmstore_
 
edm::EDPutTokenT< DQMTokenrunToken_
 

Detailed Description

Basic analyzer class which accesses 4D DTSegments and plots resolution comparing reconstructed and simulated quantities

Only true 4D segments are considered. Station 4 segments are not looked at. FIXME: Add flag to consider also segments with only phi view? Possible bias?

Residual/pull plots are filled for the reco segment with alpha closest to the simulated muon direction (defined from muon simhits in the chamber).

Efficiencies are defined as reconstructed 4D segments with alpha, beta, x, y, within 5 sigma relative to the sim muon, with sigmas specified in the config. Note that loss of even only one of the two views is considered as inefficiency!

Author
S. Bolognesi and G. Cerminara - INFN Torino

Definition at line 46 of file DTSegment4DQuality.h.

Constructor & Destructor Documentation

◆ DTSegment4DQuality()

DTSegment4DQuality::DTSegment4DQuality ( const edm::ParameterSet pset)

Constructor.

Definition at line 55 of file DTSegment4DQuality.cc.

References DTHitQualityUtils::debug, debug_, doall_, local_, muonDTDigis_cfi::pset, segment4DLabel_, segment4DToken_, sigmaResAlpha_, sigmaResBeta_, sigmaResX_, sigmaResY_, simHitLabel_, and simHitToken_.

56  // Get the debug parameter for verbose output
57  debug_ = pset.getUntrackedParameter<bool>("debug");
59 
60  // the name of the simhit collection
61  simHitLabel_ = pset.getUntrackedParameter<InputTag>("simHitLabel");
62  simHitToken_ = consumes<PSimHitContainer>(pset.getUntrackedParameter<InputTag>("simHitLabel"));
63  // the name of the 2D rec hit collection
64  segment4DLabel_ = pset.getUntrackedParameter<InputTag>("segment4DLabel");
65  segment4DToken_ = consumes<DTRecSegment4DCollection>(pset.getUntrackedParameter<InputTag>("segment4DLabel"));
66 
67  // sigma resolution on position
68  sigmaResX_ = pset.getParameter<double>("sigmaResX");
69  sigmaResY_ = pset.getParameter<double>("sigmaResY");
70  // sigma resolution on angle
71  sigmaResAlpha_ = pset.getParameter<double>("sigmaResAlpha");
72  sigmaResBeta_ = pset.getParameter<double>("sigmaResBeta");
73  doall_ = pset.getUntrackedParameter<bool>("doall", false);
74  local_ = pset.getUntrackedParameter<bool>("local", false);
75 }
edm::ESGetToken< DTGeometry, MuonGeometryRecord > muonGeomToken_
edm::InputTag segment4DLabel_
edm::InputTag simHitLabel_
std::atomic< bool > debug
edm::EDGetTokenT< DTRecSegment4DCollection > segment4DToken_
edm::EDGetTokenT< edm::PSimHitContainer > simHitToken_

Member Function Documentation

◆ bookHistograms()

void DTSegment4DQuality::bookHistograms ( DQMStore::IBooker booker,
edm::Run const &  run,
edm::EventSetup const &  setup,
dtsegment4d::Histograms histograms 
) const
overrideprivatevirtual

Book the DQM plots.

Implements DQMGlobalEDAnalyzerBase< dtsegment4d::Histograms, Args... >.

Definition at line 77 of file DTSegment4DQuality.cc.

References doall_, local_, Skims_PA_cff::name, alignCSCRings::s, and w().

80  {
81  histograms.h4DHit = std::make_unique<HRes4DHit>("All", booker, doall_, local_);
82  histograms.h4DHit_W0 = std::make_unique<HRes4DHit>("W0", booker, doall_, local_);
83  histograms.h4DHit_W1 = std::make_unique<HRes4DHit>("W1", booker, doall_, local_);
84  histograms.h4DHit_W2 = std::make_unique<HRes4DHit>("W2", booker, doall_, local_);
85 
86  if (doall_) {
87  histograms.hEff_All = std::make_unique<HEff4DHit>("All", booker);
88  histograms.hEff_W0 = std::make_unique<HEff4DHit>("W0", booker);
89  histograms.hEff_W1 = std::make_unique<HEff4DHit>("W1", booker);
90  histograms.hEff_W2 = std::make_unique<HEff4DHit>("W2", booker);
91  }
92 
93  if (local_) {
94  // Plots with finer granularity, not to be included in DQM
95  TString name = "W";
96  for (long w = 0; w <= 2; ++w) {
97  for (long s = 1; s <= 4; ++s) {
98  // FIXME station 4 is not filled
99  TString nameWS = (name + w + "_St" + s);
100  histograms.h4DHitWS[w][s - 1] = std::make_unique<HRes4DHit>(nameWS.Data(), booker, doall_, local_);
101  histograms.hEffWS[w][s - 1] = std::make_unique<HEff4DHit>(nameWS.Data(), booker);
102  }
103  }
104  }
105 };
T w() const

◆ dqmAnalyze()

void DTSegment4DQuality::dqmAnalyze ( edm::Event const &  event,
edm::EventSetup const &  setup,
dtsegment4d::Histograms const &  histograms 
) const
overrideprivatevirtual

Perform the real analysis.

Implements DQMGlobalEDAnalyzerBase< dtsegment4d::Histograms, Args... >.

Definition at line 108 of file DTSegment4DQuality.cc.

References funct::abs(), relativeConstraints::chamber, DTGeometry::chamber(), funct::cos(), gather_cfg::cout, debug_, HLT_2023v12_cff::distance, doall_, MillePedeFileConverter_cfg::e, geometryDiff::epsilon, PV3DBase< T, PVType, FrameType >::eta(), HEff4DHit::fill(), DTHitQualityUtils::findMuSimSegment(), DTHitQualityUtils::findMuSimSegmentDirAndPos(), DTHitQualityUtils::findSegmentAlphaAndBeta(), timingPdfMaker::histo, edm::HandleBase::isValid(), local_, DTRecSegment4D::localDirection(), DTRecSegment2D::localDirection(), DTRecSegment4D::localDirectionError(), DTRecSegment2D::localDirectionError(), DTRecSegment4D::localPosition(), DTRecSegment2D::localPosition(), DTRecSegment4D::localPositionError(), DTRecSegment2D::localPositionError(), DTHitQualityUtils::mapMuSimHitsPerWire(), DTHitQualityUtils::mapSimHitsPerWire(), muonGeomToken_, PV3DBase< T, PVType, FrameType >::phi(), DTRecSegment4D::phiSegment(), FastTimerService_cff::range, DTRecSegment2D::recHits(), segment4DLabel_, segment4DToken_, singleTopDQM_cfi::setup, DTHitQualityUtils::sigmaAngle(), rpcPointValidation_cfi::simHit, FastTrackerRecHitCombiner_cfi::simHits, simHitToken_, mathSSE::sqrt(), DTChamberId::station(), relativeConstraints::station, DTGeometry::superLayer(), DTSLRecSegment2D::superLayerId(), DTRecSegment2D::t0(), PV3DBase< T, PVType, FrameType >::theta(), GeomDet::toLocal(), DTChamberId::wheel(), makeMuonMisalignmentScenario::wheel, PV3DBase< T, PVType, FrameType >::x(), LocalError::xx(), PV3DBase< T, PVType, FrameType >::y(), LocalError::yy(), PV3DBase< T, PVType, FrameType >::z(), and DTRecSegment4D::zSegment().

110  {
111  const float epsilon = 5e-5; // numerical accuracy on angles [rad}
112 
113  // Get the DT Geometry
114  const DTGeometry &dtGeom = setup.getData(muonGeomToken_);
115 
116  // Get the SimHit collection from the event
118  event.getByToken(simHitToken_, simHits); // FIXME: second string to be removed
119 
120  // Map simHits by chamber
121  map<DTChamberId, PSimHitContainer> simHitsPerCh;
122  for (const auto &simHit : *simHits) {
123  // Consider only muon simhits; the others are not considered elsewhere in
124  // this class!
125  if (abs(simHit.particleType()) == 13) {
126  // Create the id of the chamber (the simHits in the DT known their wireId)
127  DTChamberId chamberId = (((DTWireId(simHit.detUnitId())).layerId()).superlayerId()).chamberId();
128  // Fill the map
129  simHitsPerCh[chamberId].push_back(simHit);
130  }
131  }
132 
133  // Get the 4D rechits from the event
135  event.getByToken(segment4DToken_, segment4Ds);
136 
137  if (!segment4Ds.isValid()) {
138  if (debug_) {
139  cout << "[DTSegment4DQuality]**Warning: no 4D Segments with label: " << segment4DLabel_
140  << " in this event, skipping!" << endl;
141  }
142  return;
143  }
144 
145  // Loop over all chambers containing a (muon) simhit
146  for (auto &simHitsInChamber : simHitsPerCh) {
147  DTChamberId chamberId = simHitsInChamber.first;
148  int station = chamberId.station();
149  if (station == 4 && !(local_)) {
150  continue; // use DTSegment2DSLPhiQuality to analyze MB4 performaces in DQM
151  }
152  int wheel = chamberId.wheel();
153 
154  //------------------------- simHits ---------------------------//
155  // Get simHits of this chamber
156  const PSimHitContainer &simHits = simHitsInChamber.second;
157 
158  // Map simhits per wire
159  auto const &simHitsPerWire = DTHitQualityUtils::mapSimHitsPerWire(simHits);
160  auto const &muSimHitPerWire = DTHitQualityUtils::mapMuSimHitsPerWire(simHitsPerWire);
161  int nMuSimHit = muSimHitPerWire.size();
162  if (nMuSimHit < 2) { // Skip chamber with less than 2 cells with mu hits
163  continue;
164  }
165  if (debug_) {
166  cout << "=== Chamber " << chamberId << " has " << nMuSimHit << " SimHits" << endl;
167  }
168 
169  // Find outer and inner mu SimHit to build a segment
170  pair<const PSimHit *, const PSimHit *> inAndOutSimHit = DTHitQualityUtils::findMuSimSegment(muSimHitPerWire);
171 
172  // Consider only sim segments crossing at least 2 SLs
173  if ((DTWireId(inAndOutSimHit.first->detUnitId())).superlayer() ==
174  (DTWireId(inAndOutSimHit.second->detUnitId())).superLayer()) {
175  continue;
176  }
177 
178  // Find direction and position of the sim Segment in Chamber RF
179  pair<LocalVector, LocalPoint> dirAndPosSimSegm =
180  DTHitQualityUtils::findMuSimSegmentDirAndPos(inAndOutSimHit, chamberId, dtGeom);
181 
182  LocalVector simSegmLocalDir = dirAndPosSimSegm.first;
183  LocalPoint simSegmLocalPos = dirAndPosSimSegm.second;
184  const DTChamber *chamber = dtGeom.chamber(chamberId);
185  GlobalPoint simSegmGlobalPos = chamber->toGlobal(simSegmLocalPos);
186  GlobalVector simSegmGlobalDir = chamber->toGlobal(simSegmLocalDir);
187 
188  // phi and theta angle of simulated segment in Chamber RF
189  float alphaSimSeg = DTHitQualityUtils::findSegmentAlphaAndBeta(simSegmLocalDir).first;
190  float betaSimSeg = DTHitQualityUtils::findSegmentAlphaAndBeta(simSegmLocalDir).second;
191  // x, y position of simulated segment in Chamber RF
192  float xSimSeg = simSegmLocalPos.x();
193  float ySimSeg = simSegmLocalPos.y();
194  // Position (in eta, phi coordinates) in lobal RF
195  float etaSimSeg = simSegmGlobalPos.eta();
196  float phiSimSeg = simSegmGlobalPos.phi();
197 
198  double count_seg = 0;
199 
200  if (debug_) {
201  cout << " Simulated segment: local direction " << simSegmLocalDir << endl
202  << " local position " << simSegmLocalPos << endl
203  << " alpha " << alphaSimSeg << endl
204  << " beta " << betaSimSeg << endl;
205  }
206 
207  //---------------------------- recHits --------------------------//
208  // Get the range of rechit for the corresponding chamberId
209  bool recHitFound = false;
210  DTRecSegment4DCollection::range range = segment4Ds->get(chamberId);
211  int nsegm = distance(range.first, range.second);
212  if (debug_) {
213  cout << " Chamber: " << chamberId << " has " << nsegm << " 4D segments" << endl;
214  }
215 
216  if (nsegm != 0) {
217  // Find the best RecHit: look for the 4D RecHit with the phi angle closest
218  // to that of segment made of SimHits.
219  // RecHits must have delta alpha and delta position within 5 sigma of
220  // the residual distribution (we are looking for residuals of segments
221  // usefull to the track fit) for efficency purpose
222  const DTRecSegment4D *bestRecHit = nullptr;
223  double deltaAlpha = 99999;
224  double deltaBeta = 99999;
225 
226  // Loop over the recHits of this chamberId
227  for (DTRecSegment4DCollection::const_iterator segment4D = range.first; segment4D != range.second; ++segment4D) {
228  // Consider only segments with both projections
229  if (station != 4 && (*segment4D).dimension() != 4) {
230  continue;
231  }
232  // Segment Local Direction and position (in Chamber RF)
233  LocalVector recSegDirection = (*segment4D).localDirection();
234  LocalPoint recSegPosition = (*segment4D).localPosition();
235 
236  pair<double, double> ab = DTHitQualityUtils::findSegmentAlphaAndBeta(recSegDirection);
237  float recSegAlpha = ab.first;
238  float recSegBeta = ab.second;
239 
240  if (debug_) {
241  cout << &(*segment4D) << " RecSegment direction: " << recSegDirection << endl
242  << " position : " << (*segment4D).localPosition() << endl
243  << " alpha : " << recSegAlpha << endl
244  << " beta : " << recSegBeta << endl
245  << " nhits : " << (*segment4D).phiSegment()->recHits().size() << " "
246  << (((*segment4D).zSegment() != nullptr) ? (*segment4D).zSegment()->recHits().size() : 0) << endl;
247  }
248 
249  float dAlphaRecSim = fabs(recSegAlpha - alphaSimSeg);
250  float dBetaRecSim = fabs(recSegBeta - betaSimSeg);
251 
252  if ((fabs(recSegPosition.x() - simSegmLocalPos.x()) <
253  4) // require rec and sim segments to be ~in the same cell in x
254  && ((fabs(recSegPosition.y() - simSegmLocalPos.y()) < 4) ||
255  (*segment4D).dimension() < 4)) { // ~in the same cell in y, if segment has the theta view
256  ++count_seg;
257 
258  if (fabs(dAlphaRecSim - deltaAlpha) < epsilon) { // Numerically equivalent alphas, choose based on beta
259  if (dBetaRecSim < deltaBeta) {
260  deltaAlpha = dAlphaRecSim;
261  deltaBeta = dBetaRecSim;
262  bestRecHit = &(*segment4D);
263  }
264 
265  } else if (dAlphaRecSim < deltaAlpha) {
266  deltaAlpha = dAlphaRecSim;
267  deltaBeta = dBetaRecSim;
268  bestRecHit = &(*segment4D);
269  }
270  }
271 
272  } // End of Loop over all 4D RecHits
273 
274  if (bestRecHit) {
275  if (debug_) {
276  cout << endl << "Chosen: " << bestRecHit << endl;
277  }
278  // Best rechit direction and position in Chamber RF
279  LocalPoint bestRecHitLocalPos = bestRecHit->localPosition();
280  LocalVector bestRecHitLocalDir = bestRecHit->localDirection();
281  // Errors on x and y
282  LocalError bestRecHitLocalPosErr = bestRecHit->localPositionError();
283  LocalError bestRecHitLocalDirErr = bestRecHit->localDirectionError();
284 
285  pair<double, double> ab = DTHitQualityUtils::findSegmentAlphaAndBeta(bestRecHitLocalDir);
286  float alphaBestRHit = ab.first;
287  float betaBestRHit = ab.second;
288  // Errors on alpha and beta
289 
290  // Get position and direction using the rx projection (so in SL
291  // reference frame). Note that x (and y) are swapped wrt to Chamber
292  // frame
293  // if (bestRecHit->hasZed()) //
294  const DTSLRecSegment2D *zedRecSeg = bestRecHit->zSegment();
295  LocalPoint bestRecHitLocalPosRZ;
296  LocalVector bestRecHitLocalDirRZ;
297  LocalError bestRecHitLocalPosErrRZ;
298  LocalError bestRecHitLocalDirErrRZ;
299  LocalPoint simSegLocalPosRZ; // FIXME: this is not set for segments with
300  // only the phi view.
301  float alphaBestRHitRZ = 0; // angle measured in the RZ SL, in its own frame
302  float alphaSimSegRZ = betaSimSeg;
303  if (zedRecSeg) {
304  bestRecHitLocalPosRZ = zedRecSeg->localPosition();
305  bestRecHitLocalDirRZ = zedRecSeg->localDirection();
306  // Errors on x and y
307  bestRecHitLocalPosErrRZ = zedRecSeg->localPositionError();
308  bestRecHitLocalDirErrRZ = zedRecSeg->localDirectionError();
309 
310  // angle measured in the RZ SL, in its own frame
311  alphaBestRHitRZ = DTHitQualityUtils::findSegmentAlphaAndBeta(bestRecHitLocalDirRZ).first;
312 
313  // Get SimSeg position and Direction in rZ SL frame
314  const DTSuperLayer *sl = dtGeom.superLayer(zedRecSeg->superLayerId());
315  LocalPoint simSegLocalPosRZTmp = sl->toLocal(simSegmGlobalPos);
316  LocalVector simSegLocalDirRZ = sl->toLocal(simSegmGlobalDir);
317  simSegLocalPosRZ =
318  simSegLocalPosRZTmp + simSegLocalDirRZ * (-simSegLocalPosRZTmp.z() / (cos(simSegLocalDirRZ.theta())));
319  alphaSimSegRZ = DTHitQualityUtils::findSegmentAlphaAndBeta(simSegLocalDirRZ).first;
320 
321  if (debug_) {
322  cout << "RZ SL: recPos " << bestRecHitLocalPosRZ << "recDir " << bestRecHitLocalDirRZ << "recAlpha "
323  << alphaBestRHitRZ << endl
324  << "RZ SL: simPos " << simSegLocalPosRZ << "simDir " << simSegLocalDirRZ << "simAlpha "
325  << alphaSimSegRZ << endl;
326  }
327  }
328 
329  // get nhits and t0
330  const DTChamberRecSegment2D *phiSeg = bestRecHit->phiSegment();
331 
332  float t0phi = -999;
333  float t0theta = -999;
334  int nHitPhi = 0;
335  int nHitTheta = 0;
336 
337  if (phiSeg) {
338  t0phi = phiSeg->t0();
339  nHitPhi = phiSeg->recHits().size();
340  }
341 
342  if (zedRecSeg) {
343  t0theta = zedRecSeg->t0();
344  nHitTheta = zedRecSeg->recHits().size();
345  }
346 
347  recHitFound = true;
348 
349  // Mirror alpha in phi SLs so that + and - wheels can be plotted
350  // together
351  if (mirrorMinusWheels && wheel < 0) {
352  alphaSimSeg *= -1.;
353  alphaBestRHit *= -1.;
354  // Note: local X (xSimSeg, bestRecHitLocalPos.x() would have to be
355  // mirrored as well; but at the moment there is no plot of dependency
356  // vs X, except for efficiency.
357  }
358 
359  // Fill Residual histos
360  HRes4DHit *histo = nullptr;
361 
362  if (wheel == 0) {
363  histo = histograms.h4DHit_W0.get();
364  } else if (abs(wheel) == 1) {
365  histo = histograms.h4DHit_W1.get();
366  } else if (abs(wheel) == 2) {
367  histo = histograms.h4DHit_W2.get();
368  }
369 
370  float sigmaAlphaBestRhit = sqrt(DTHitQualityUtils::sigmaAngle(alphaBestRHit, bestRecHitLocalDirErr.xx()));
371  float sigmaBetaBestRhit =
372  sqrt(DTHitQualityUtils::sigmaAngle(betaBestRHit,
373  bestRecHitLocalDirErr.yy())); // FIXME this misses the contribution
374  // from uncertainty in extrapolation!
375  float sigmaAlphaBestRhitRZ = sqrt(DTHitQualityUtils::sigmaAngle(alphaBestRHitRZ, bestRecHitLocalDirErrRZ.xx()));
376 
377  histo->fill(alphaSimSeg,
378  alphaBestRHit,
379  betaSimSeg,
380  betaBestRHit,
381  xSimSeg,
382  bestRecHitLocalPos.x(),
383  ySimSeg,
384  bestRecHitLocalPos.y(),
385  etaSimSeg,
386  phiSimSeg,
387  bestRecHitLocalPosRZ.x(),
388  simSegLocalPosRZ.x(),
389  alphaBestRHitRZ,
390  alphaSimSegRZ,
391  sigmaAlphaBestRhit,
392  sigmaBetaBestRhit,
393  sqrt(bestRecHitLocalPosErr.xx()),
394  sqrt(bestRecHitLocalPosErr.yy()),
395  sigmaAlphaBestRhitRZ,
396  sqrt(bestRecHitLocalPosErrRZ.xx()),
397  nHitPhi,
398  nHitTheta,
399  t0phi,
400  t0theta);
401 
402  histograms.h4DHit->fill(alphaSimSeg,
403  alphaBestRHit,
404  betaSimSeg,
405  betaBestRHit,
406  xSimSeg,
407  bestRecHitLocalPos.x(),
408  ySimSeg,
409  bestRecHitLocalPos.y(),
410  etaSimSeg,
411  phiSimSeg,
412  bestRecHitLocalPosRZ.x(),
413  simSegLocalPosRZ.x(),
414  alphaBestRHitRZ,
415  alphaSimSegRZ,
416  sigmaAlphaBestRhit,
417  sigmaBetaBestRhit,
418  sqrt(bestRecHitLocalPosErr.xx()),
419  sqrt(bestRecHitLocalPosErr.yy()),
420  sigmaAlphaBestRhitRZ,
421  sqrt(bestRecHitLocalPosErrRZ.xx()),
422  nHitPhi,
423  nHitTheta,
424  t0phi,
425  t0theta);
426 
427  if (local_) {
428  histograms.h4DHitWS[abs(wheel)][station - 1]->fill(alphaSimSeg,
429  alphaBestRHit,
430  betaSimSeg,
431  betaBestRHit,
432  xSimSeg,
433  bestRecHitLocalPos.x(),
434  ySimSeg,
435  bestRecHitLocalPos.y(),
436  etaSimSeg,
437  phiSimSeg,
438  bestRecHitLocalPosRZ.x(),
439  simSegLocalPosRZ.x(),
440  alphaBestRHitRZ,
441  alphaSimSegRZ,
442  sigmaAlphaBestRhit,
443  sigmaBetaBestRhit,
444  sqrt(bestRecHitLocalPosErr.xx()),
445  sqrt(bestRecHitLocalPosErr.yy()),
446  sigmaAlphaBestRhitRZ,
447  sqrt(bestRecHitLocalPosErrRZ.xx()),
448  nHitPhi,
449  nHitTheta,
450  t0phi,
451  t0theta);
452  }
453 
454  } // end of if (bestRecHit)
455 
456  } // end of if (nsegm!= 0)
457 
458  // Fill Efficiency plot
459  if (doall_) {
460  HEff4DHit *heff = nullptr;
461 
462  if (wheel == 0) {
463  heff = histograms.hEff_W0.get();
464  } else if (abs(wheel) == 1) {
465  heff = histograms.hEff_W1.get();
466  } else if (abs(wheel) == 2) {
467  heff = histograms.hEff_W2.get();
468  }
469  heff->fill(etaSimSeg, phiSimSeg, xSimSeg, ySimSeg, alphaSimSeg, betaSimSeg, recHitFound, count_seg);
470  histograms.hEff_All->fill(
471  etaSimSeg, phiSimSeg, xSimSeg, ySimSeg, alphaSimSeg, betaSimSeg, recHitFound, count_seg);
472  if (local_) {
473  histograms.hEffWS[abs(wheel)][station - 1]->fill(
474  etaSimSeg, phiSimSeg, xSimSeg, ySimSeg, alphaSimSeg, betaSimSeg, recHitFound, count_seg);
475  }
476  }
477  } // End of loop over chambers
478 }
int station() const
Return the station number.
Definition: DTChamberId.h:42
std::pair< double, double > findSegmentAlphaAndBeta(const LocalVector &direction)
std::map< DTWireId, edm::PSimHitContainer > mapSimHitsPerWire(const edm::PSimHitContainer &simhits)
std::pair< const_iterator, const_iterator > range
iterator range
Definition: RangeMap.h:50
LocalPoint localPosition() const override
local position in SL frame
LocalVector localDirection() const override
the local direction in SL frame
LocalPoint toLocal(const GlobalPoint &gp) const
Conversion to the R.F. of the GeomDet.
Definition: GeomDet.h:58
double t0() const
Get the segment t0 (if recomputed, 0 is returned otherwise)
edm::ESGetToken< DTGeometry, MuonGeometryRecord > muonGeomToken_
T z() const
Definition: PV3DBase.h:61
Geom::Phi< T > phi() const
Definition: PV3DBase.h:66
LocalError localDirectionError() const override
Local direction error in the Chamber frame.
T eta() const
Definition: PV3DBase.h:73
std::pair< LocalVector, LocalPoint > findMuSimSegmentDirAndPos(const std::pair< const PSimHit *, const PSimHit *> &inAndOutSimHit, const DetId detId, const DTGeometry &muonGeom)
LocalVector localDirection() const override
Local direction in Chamber frame.
edm::InputTag segment4DLabel_
LocalPoint localPosition() const override
Local position in Chamber frame.
const DTSLRecSegment2D * zSegment() const
The Z segment: 0 if not zed projection available.
edm::EDGetTokenT< DTRecSegment4DCollection > segment4DToken_
C::const_iterator const_iterator
constant access iterator type
Definition: RangeMap.h:43
float yy() const
Definition: LocalError.h:24
T x() const
Definition: PV3DBase.h:59
T y() const
Definition: PV3DBase.h:60
const DTSuperLayer * superLayer(const DTSuperLayerId &id) const
Return a DTSuperLayer given its id.
Definition: DTGeometry.cc:92
T sqrt(T t)
Definition: SSEVec.h:19
void fill(float etaSimSegm, float phiSimSegm, float xSimSegm, float ySimSegm, float alphaSimSegm, float betaSimSegm, bool fillRecHit, int nSeg)
Definition: Histograms.h:975
LocalError localDirectionError() const override
the local direction error (xx,xy,yy) in SL frame: only xx is not 0.
std::map< DTWireId, const PSimHit * > mapMuSimHitsPerWire(const std::map< DTWireId, edm::PSimHitContainer > &simHitWireMap)
Create a map between the Mu SimHits and corresponding MuBarWireId ;.
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
A set of histograms for efficiency 4D RecHits (producer)
Definition: Histograms.h:940
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
LocalError localPositionError() const override
local position error in SL frame
std::vector< const TrackingRecHit * > recHits() const override
Access to component RecHits (if any)
DTSuperLayerId superLayerId() const
The id of the superlayer on which reside the segment.
const DTChamberRecSegment2D * phiSegment() const
The superPhi segment: 0 if no phi projection available.
std::pair< const PSimHit *, const PSimHit * > findMuSimSegment(const std::map< DTWireId, const PSimHit *> &mapWireAndMuSimHit)
bool isValid() const
Definition: HandleBase.h:70
int wheel() const
Return the wheel number.
Definition: DTChamberId.h:39
LocalError localPositionError() const override
Local position error in Chamber frame.
std::vector< PSimHit > PSimHitContainer
double sigmaAngle(double Angle, double sigma2TanAngle)
edm::EDGetTokenT< edm::PSimHitContainer > simHitToken_
const DTChamber * chamber(const DTChamberId &id) const
Return a DTChamber given its id.
Definition: DTGeometry.cc:90
float xx() const
Definition: LocalError.h:22
Geom::Theta< T > theta() const
Definition: PV3DBase.h:72

Member Data Documentation

◆ debug_

bool DTSegment4DQuality::debug_
private

Definition at line 83 of file DTSegment4DQuality.h.

Referenced by dqmAnalyze(), and DTSegment4DQuality().

◆ doall_

bool DTSegment4DQuality::doall_
private

Definition at line 79 of file DTSegment4DQuality.h.

Referenced by bookHistograms(), dqmAnalyze(), and DTSegment4DQuality().

◆ local_

bool DTSegment4DQuality::local_
private

Definition at line 80 of file DTSegment4DQuality.h.

Referenced by bookHistograms(), dqmAnalyze(), and DTSegment4DQuality().

◆ muonGeomToken_

edm::ESGetToken<DTGeometry, MuonGeometryRecord> DTSegment4DQuality::muonGeomToken_
private

Definition at line 69 of file DTSegment4DQuality.h.

Referenced by dqmAnalyze().

◆ segment4DLabel_

edm::InputTag DTSegment4DQuality::segment4DLabel_
private

Definition at line 64 of file DTSegment4DQuality.h.

Referenced by dqmAnalyze(), and DTSegment4DQuality().

◆ segment4DToken_

edm::EDGetTokenT<DTRecSegment4DCollection> DTSegment4DQuality::segment4DToken_
private

Definition at line 66 of file DTSegment4DQuality.h.

Referenced by dqmAnalyze(), and DTSegment4DQuality().

◆ sigmaResAlpha_

double DTSegment4DQuality::sigmaResAlpha_
private

Definition at line 76 of file DTSegment4DQuality.h.

Referenced by DTSegment4DQuality().

◆ sigmaResBeta_

double DTSegment4DQuality::sigmaResBeta_
private

Definition at line 77 of file DTSegment4DQuality.h.

Referenced by DTSegment4DQuality().

◆ sigmaResX_

double DTSegment4DQuality::sigmaResX_
private

Definition at line 72 of file DTSegment4DQuality.h.

Referenced by DTSegment4DQuality().

◆ sigmaResY_

double DTSegment4DQuality::sigmaResY_
private

Definition at line 73 of file DTSegment4DQuality.h.

Referenced by DTSegment4DQuality().

◆ simHitLabel_

edm::InputTag DTSegment4DQuality::simHitLabel_
private

Definition at line 63 of file DTSegment4DQuality.h.

Referenced by DTSegment4DQuality().

◆ simHitToken_

edm::EDGetTokenT<edm::PSimHitContainer> DTSegment4DQuality::simHitToken_
private

Definition at line 65 of file DTSegment4DQuality.h.

Referenced by dqmAnalyze(), and DTSegment4DQuality().