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 > edm::global::EDAnalyzer< edm::RunCache< dtsegment4d::Histograms >, Args... > edm::global::EDAnalyzerBase edm::EDConsumerBase

Public Member Functions

 DTSegment4DQuality (const edm::ParameterSet &pset)
 Constructor. More...
 
- Public Member Functions inherited from edm::global::EDAnalyzer< edm::RunCache< dtsegment4d::Histograms >, Args... >
 EDAnalyzer ()=default
 
bool wantsGlobalLuminosityBlocks () const final
 
bool wantsGlobalRuns () const final
 
bool wantsStreamLuminosityBlocks () const final
 
bool wantsStreamRuns () const final
 
- Public Member Functions inherited from edm::global::EDAnalyzerBase
void callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)
 
 EDAnalyzerBase ()
 
ModuleDescription const & moduleDescription () const
 
 ~EDAnalyzerBase () 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
 
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::vector< ModuleDescription const * > &modules, 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
 
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::ConcurrentBooker &, 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::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 edm::global::EDAnalyzerBase
typedef EDAnalyzerBase ModuleType
 
- Public Types inherited from edm::EDConsumerBase
typedef ProductLabels Labels
 
- Static Public Member Functions inherited from edm::global::EDAnalyzerBase
static const std::string & baseType ()
 
static void fillDescriptions (ConfigurationDescriptions &descriptions)
 
static void prevalidate (ConfigurationDescriptions &descriptions)
 
- Protected Member Functions inherited from edm::EDConsumerBase
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > consumes (edm::InputTag const &tag)
 
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<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)
 

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 ( const edm::ParameterSet pset)

Constructor.

Definition at line 56 of file DTSegment4DQuality.cc.

References DTHitQualityUtils::debug, edm::ParameterSet::getParameter(), and edm::ParameterSet::getUntrackedParameter().

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

Member Function Documentation

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

Book the DQM plots.

Implements DQMGlobalEDAnalyzer< dtsegment4d::Histograms >.

Definition at line 78 of file DTSegment4DQuality.cc.

References dtsegment4d::Histograms::h4DHit, dtsegment4d::Histograms::h4DHit_W0, dtsegment4d::Histograms::h4DHit_W1, dtsegment4d::Histograms::h4DHit_W2, dtsegment4d::Histograms::h4DHitWS, dtsegment4d::Histograms::hEff_All, dtsegment4d::Histograms::hEff_W0, dtsegment4d::Histograms::hEff_W1, dtsegment4d::Histograms::hEff_W2, dtsegment4d::Histograms::hEffWS, dataset::name, alignCSCRings::s, and w.

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

Perform the real analysis.

Definition at line 109 of file DTSegment4DQuality.cc.

References funct::abs(), relativeConstraints::chamber, DTGeometry::chamber(), funct::cos(), gather_cfg::cout, DEFINE_FWK_MODULE, SoftLeptonByDistance_cfi::distance, MillePedeFileConverter_cfg::e, geometryDiff::epsilon, PV3DBase< T, PVType, FrameType >::eta(), HRes4DHit::fill(), HEff4DHit::fill(), DTHitQualityUtils::findMuSimSegment(), DTHitQualityUtils::findMuSimSegmentDirAndPos(), DTHitQualityUtils::findSegmentAlphaAndBeta(), edm::EventSetup::get(), dtsegment4d::Histograms::h4DHit, dtsegment4d::Histograms::h4DHit_W0, dtsegment4d::Histograms::h4DHit_W1, dtsegment4d::Histograms::h4DHit_W2, dtsegment4d::Histograms::h4DHitWS, dtsegment4d::Histograms::hEff_All, dtsegment4d::Histograms::hEff_W0, dtsegment4d::Histograms::hEff_W1, dtsegment4d::Histograms::hEff_W2, dtsegment4d::Histograms::hEffWS, trackerHits::histo, edm::HandleBase::isValid(), DTRecSegment4D::localDirection(), DTRecSegment2D::localDirection(), DTRecSegment4D::localDirectionError(), DTRecSegment2D::localDirectionError(), DTRecSegment4D::localPosition(), DTRecSegment2D::localPosition(), DTRecSegment4D::localPositionError(), DTRecSegment2D::localPositionError(), DTHitQualityUtils::mapMuSimHitsPerWire(), DTHitQualityUtils::mapSimHitsPerWire(), PV3DBase< T, PVType, FrameType >::phi(), DTRecSegment4D::phiSegment(), DTRecSegment2D::recHits(), DTHitQualityUtils::sigmaAngle(), rpcPointValidation_cfi::simHit, trackerHits::simHits, mathSSE::sqrt(), DTChamberId::station(), relativeConstraints::station, DTGeometry::superLayer(), DTSLRecSegment2D::superLayerId(), DTRecSegment2D::t0(), PV3DBase< T, PVType, FrameType >::theta(), GeomDet::toGlobal(), 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().

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

Member Data Documentation

bool DTSegment4DQuality::debug_
private

Definition at line 80 of file DTSegment4DQuality.h.

bool DTSegment4DQuality::doall_
private

Definition at line 76 of file DTSegment4DQuality.h.

bool DTSegment4DQuality::local_
private

Definition at line 77 of file DTSegment4DQuality.h.

edm::InputTag DTSegment4DQuality::segment4DLabel_
private

Definition at line 64 of file DTSegment4DQuality.h.

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

Definition at line 66 of file DTSegment4DQuality.h.

double DTSegment4DQuality::sigmaResAlpha_
private

Definition at line 73 of file DTSegment4DQuality.h.

double DTSegment4DQuality::sigmaResBeta_
private

Definition at line 74 of file DTSegment4DQuality.h.

double DTSegment4DQuality::sigmaResX_
private

Definition at line 69 of file DTSegment4DQuality.h.

double DTSegment4DQuality::sigmaResY_
private

Definition at line 70 of file DTSegment4DQuality.h.

edm::InputTag DTSegment4DQuality::simHitLabel_
private

Definition at line 63 of file DTSegment4DQuality.h.

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

Definition at line 65 of file DTSegment4DQuality.h.