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 > > 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 > >
 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
 
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)
 
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 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 45 of file DTSegment4DQuality.h.

Constructor & Destructor Documentation

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

Constructor.

Definition at line 57 of file DTSegment4DQuality.cc.

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

57  {
58  // Get the debug parameter for verbose output
59  debug_ = pset.getUntrackedParameter<bool>("debug");
61 
62  // the name of the simhit collection
63  simHitLabel_ = pset.getUntrackedParameter<InputTag>("simHitLabel");
64  simHitToken_ = consumes<PSimHitContainer>(pset.getUntrackedParameter<InputTag>("simHitLabel"));
65  // the name of the 2D rec hit collection
66  segment4DLabel_ = pset.getUntrackedParameter<InputTag>("segment4DLabel");
67  segment4DToken_ = consumes<DTRecSegment4DCollection>(pset.getUntrackedParameter<InputTag>("segment4DLabel"));
68 
69  // sigma resolution on position
70  sigmaResX_ = pset.getParameter<double>("sigmaResX");
71  sigmaResY_ = pset.getParameter<double>("sigmaResY");
72  // sigma resolution on angle
73  sigmaResAlpha_ = pset.getParameter<double>("sigmaResAlpha");
74  sigmaResBeta_ = pset.getParameter<double>("sigmaResBeta");
75  doall_ = pset.getUntrackedParameter<bool>("doall", false);
76  local_ = pset.getUntrackedParameter<bool>("local", false);
77 }
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 79 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.

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

107  {
108  const float epsilon = 5e-5; // numerical accuracy on angles [rad}
109 
110  // Get the DT Geometry
111  ESHandle<DTGeometry> dtGeom;
112  setup.get<MuonGeometryRecord>().get(dtGeom);
113 
114  // Get the SimHit collection from the event
116  event.getByToken(simHitToken_, simHits); // FIXME: second string to be removed
117 
118  // Map simHits by chamber
119  map<DTChamberId, PSimHitContainer > simHitsPerCh;
120  for (const auto & simHit : *simHits) {
121 
122  // Consider only muon simhits; the others are not considered elsewhere in this class!
123  if (abs(simHit.particleType())== 13) {
124  // Create the id of the chamber (the simHits in the DT known their wireId)
125  DTChamberId chamberId = (((DTWireId(simHit.detUnitId())).layerId()).superlayerId()).chamberId();
126  // Fill the map
127  simHitsPerCh[chamberId].push_back(simHit);
128  }
129  }
130 
131  // Get the 4D rechits from the event
133  event.getByToken(segment4DToken_, segment4Ds);
134 
135  if (!segment4Ds.isValid()) {
136  if (debug_) {
137  cout << "[DTSegment4DQuality]**Warning: no 4D Segments with label: " << segment4DLabel_ << " in this event, skipping!" << endl;
138 }
139  return;
140  }
141 
142  // Loop over all chambers containing a (muon) simhit
143  for (auto & simHitsInChamber : simHitsPerCh) {
144 
145  DTChamberId chamberId = simHitsInChamber.first;
146  int station = chamberId.station();
147  if (station == 4 && !(local_) ) { continue; // use DTSegment2DSLPhiQuality to analyze MB4 performaces in DQM
148 }
149  int wheel = chamberId.wheel();
150 
151  //------------------------- simHits ---------------------------//
152  // Get simHits of this chamber
153  const PSimHitContainer& simHits = simHitsInChamber.second;
154 
155  // Map simhits per wire
156  auto const& simHitsPerWire = DTHitQualityUtils::mapSimHitsPerWire(simHits);
157  auto const& muSimHitPerWire = DTHitQualityUtils::mapMuSimHitsPerWire(simHitsPerWire);
158  int nMuSimHit = muSimHitPerWire.size();
159  if (nMuSimHit <2) { // Skip chamber with less than 2 cells with mu hits
160  continue;
161  }
162  if (debug_) {
163  cout << "=== Chamber " << chamberId << " has " << nMuSimHit << " SimHits" << endl;
164 }
165 
166  // Find outer and inner mu SimHit to build a segment
167  pair<const PSimHit*, const PSimHit*> inAndOutSimHit = DTHitQualityUtils::findMuSimSegment(muSimHitPerWire);
168 
169  // Consider only sim segments crossing at least 2 SLs
170  if ((DTWireId(inAndOutSimHit.first->detUnitId())).superlayer() == (DTWireId(inAndOutSimHit.second->detUnitId())).superLayer()) { continue;
171 }
172 
173  // Find direction and position of the sim Segment in Chamber RF
174  pair<LocalVector, LocalPoint> dirAndPosSimSegm = DTHitQualityUtils::findMuSimSegmentDirAndPos(inAndOutSimHit,
175  chamberId,&(*dtGeom));
176 
177  LocalVector simSegmLocalDir = dirAndPosSimSegm.first;
178  LocalPoint simSegmLocalPos = dirAndPosSimSegm.second;
179  const DTChamber* chamber = dtGeom->chamber(chamberId);
180  GlobalPoint simSegmGlobalPos = chamber->toGlobal(simSegmLocalPos);
181  GlobalVector simSegmGlobalDir = chamber->toGlobal(simSegmLocalDir);
182 
183  // phi and theta angle of simulated segment in Chamber RF
184  float alphaSimSeg = DTHitQualityUtils::findSegmentAlphaAndBeta(simSegmLocalDir).first;
185  float betaSimSeg = DTHitQualityUtils::findSegmentAlphaAndBeta(simSegmLocalDir).second;
186  // x, y position of simulated segment in Chamber RF
187  float xSimSeg = simSegmLocalPos.x();
188  float ySimSeg = simSegmLocalPos.y();
189  // Position (in eta, phi coordinates) in lobal RF
190  float etaSimSeg = simSegmGlobalPos.eta();
191  float phiSimSeg = simSegmGlobalPos.phi();
192 
193  double count_seg = 0;
194 
195  if (debug_) {
196  cout << " Simulated segment: local direction " << simSegmLocalDir << endl
197  << " local position " << simSegmLocalPos << endl
198  << " alpha " << alphaSimSeg << endl
199  << " beta " << betaSimSeg << endl;
200 }
201 
202  //---------------------------- recHits --------------------------//
203  // Get the range of rechit for the corresponding chamberId
204  bool recHitFound = false;
205  DTRecSegment4DCollection::range range = segment4Ds->get(chamberId);
206  int nsegm = distance(range.first, range.second);
207  if (debug_) {
208  cout << " Chamber: " << chamberId << " has " << nsegm
209  << " 4D segments" << endl;
210 }
211 
212  if (nsegm!= 0) {
213  // Find the best RecHit: look for the 4D RecHit with the phi angle closest
214  // to that of segment made of SimHits.
215  // RecHits must have delta alpha and delta position within 5 sigma of
216  // the residual distribution (we are looking for residuals of segments
217  // usefull to the track fit) for efficency purpose
218  const DTRecSegment4D* bestRecHit = nullptr;
219  double deltaAlpha = 99999;
220  double deltaBeta = 99999;
221 
222  // Loop over the recHits of this chamberId
223  for (DTRecSegment4DCollection::const_iterator segment4D = range.first;
224  segment4D!= range.second;
225  ++segment4D) {
226 
227  // Consider only segments with both projections
228  if (station!= 4 && (*segment4D).dimension() != 4) {
229  continue;
230  }
231  // Segment Local Direction and position (in Chamber RF)
232  LocalVector recSegDirection = (*segment4D).localDirection();
233  LocalPoint recSegPosition = (*segment4D).localPosition();
234 
235  pair<double, double> ab = DTHitQualityUtils::findSegmentAlphaAndBeta(recSegDirection);
236  float recSegAlpha = ab.first;
237  float recSegBeta = ab.second;
238 
239  if (debug_) {
240  cout << &(*segment4D)
241  << " RecSegment direction: " << recSegDirection << endl
242  << " position : " << (*segment4D).localPosition() << endl
243  << " alpha : " << recSegAlpha << endl
244  << " beta : " << recSegBeta << endl
245  << " nhits : " << (*segment4D).phiSegment()->recHits().size() << " " << (((*segment4D).zSegment()!= nullptr)?(*segment4D).zSegment()->recHits().size():0)
246  << endl;
247 }
248 
249 
250  float dAlphaRecSim = fabs(recSegAlpha - alphaSimSeg);
251  float dBetaRecSim = fabs(recSegBeta - betaSimSeg);
252 
253 
254  if ((fabs(recSegPosition.x()-simSegmLocalPos.x())<4) // require rec and sim segments to be ~in the same cell in x
255  && ((fabs(recSegPosition.y()-simSegmLocalPos.y())<4)||(*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 only the phi view.
300  float alphaBestRHitRZ = 0; // angle measured in the RZ SL, in its own frame
301  float alphaSimSegRZ = betaSimSeg;
302  if (zedRecSeg) {
303  bestRecHitLocalPosRZ = zedRecSeg->localPosition();
304  bestRecHitLocalDirRZ = zedRecSeg->localDirection();
305  // Errors on x and y
306  bestRecHitLocalPosErrRZ = zedRecSeg->localPositionError();
307  bestRecHitLocalDirErrRZ = zedRecSeg->localDirectionError();
308 
309  // angle measured in the RZ SL, in its own frame
310  alphaBestRHitRZ = DTHitQualityUtils::findSegmentAlphaAndBeta(bestRecHitLocalDirRZ).first;
311 
312  // Get SimSeg position and Direction in rZ SL frame
313  const DTSuperLayer* sl = dtGeom->superLayer(zedRecSeg->superLayerId());
314  LocalPoint simSegLocalPosRZTmp = sl->toLocal(simSegmGlobalPos);
315  LocalVector simSegLocalDirRZ = sl->toLocal(simSegmGlobalDir);
316  simSegLocalPosRZ =
317  simSegLocalPosRZTmp + simSegLocalDirRZ*(-simSegLocalPosRZTmp.z()/(cos(simSegLocalDirRZ.theta())));
318  alphaSimSegRZ = DTHitQualityUtils::findSegmentAlphaAndBeta(simSegLocalDirRZ).first;
319 
320  if (debug_) {
321  cout <<
322  "RZ SL: recPos " << bestRecHitLocalPosRZ <<
323  "recDir " << bestRecHitLocalDirRZ <<
324  "recAlpha " << alphaBestRHitRZ << endl <<
325  "RZ SL: simPos " << simSegLocalPosRZ <<
326  "simDir " << simSegLocalDirRZ <<
327  "simAlpha " << 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 together
352  if (mirrorMinusWheels && wheel<0) {
353  alphaSimSeg*=-1.;
354  alphaBestRHit*=-1.;
355  // Note: local X (xSimSeg, bestRecHitLocalPos.x() would have to be mirrored as well;
356  // but at the moment there is no plot of dependency 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 = sqrt(DTHitQualityUtils::sigmaAngle(betaBestRHit, bestRecHitLocalDirErr.yy())); // FIXME this misses the contribution from uncertainty in extrapolation!
372  float sigmaAlphaBestRhitRZ = sqrt(DTHitQualityUtils::sigmaAngle(alphaBestRHitRZ, bestRecHitLocalDirErrRZ.xx()));
373 
374  histo->fill(alphaSimSeg,
375  alphaBestRHit,
376  betaSimSeg,
377  betaBestRHit,
378  xSimSeg,
379  bestRecHitLocalPos.x(),
380  ySimSeg,
381  bestRecHitLocalPos.y(),
382  etaSimSeg,
383  phiSimSeg,
384  bestRecHitLocalPosRZ.x(),
385  simSegLocalPosRZ.x(),
386  alphaBestRHitRZ,
387  alphaSimSegRZ,
388  sigmaAlphaBestRhit,
389  sigmaBetaBestRhit,
390  sqrt(bestRecHitLocalPosErr.xx()),
391  sqrt(bestRecHitLocalPosErr.yy()),
392  sigmaAlphaBestRhitRZ,
393  sqrt(bestRecHitLocalPosErrRZ.xx()),
394  nHitPhi, nHitTheta, t0phi, t0theta);
395 
396  histograms.h4DHit->fill(alphaSimSeg,
397  alphaBestRHit,
398  betaSimSeg,
399  betaBestRHit,
400  xSimSeg,
401  bestRecHitLocalPos.x(),
402  ySimSeg,
403  bestRecHitLocalPos.y(),
404  etaSimSeg,
405  phiSimSeg,
406  bestRecHitLocalPosRZ.x(),
407  simSegLocalPosRZ.x(),
408  alphaBestRHitRZ,
409  alphaSimSegRZ,
410  sigmaAlphaBestRhit,
411  sigmaBetaBestRhit,
412  sqrt(bestRecHitLocalPosErr.xx()),
413  sqrt(bestRecHitLocalPosErr.yy()),
414  sigmaAlphaBestRhitRZ,
415  sqrt(bestRecHitLocalPosErrRZ.xx()),
416  nHitPhi, nHitTheta, t0phi, t0theta);
417 
418  if (local_) { histograms.h4DHitWS[abs(wheel)][station-1]->fill(alphaSimSeg,
419  alphaBestRHit,
420  betaSimSeg,
421  betaBestRHit,
422  xSimSeg,
423  bestRecHitLocalPos.x(),
424  ySimSeg,
425  bestRecHitLocalPos.y(),
426  etaSimSeg,
427  phiSimSeg,
428  bestRecHitLocalPosRZ.x(),
429  simSegLocalPosRZ.x(),
430  alphaBestRHitRZ,
431  alphaSimSegRZ,
432  sigmaAlphaBestRhit,
433  sigmaBetaBestRhit,
434  sqrt(bestRecHitLocalPosErr.xx()),
435  sqrt(bestRecHitLocalPosErr.yy()),
436  sigmaAlphaBestRhitRZ,
437  sqrt(bestRecHitLocalPosErrRZ.xx()),
438  nHitPhi, nHitTheta, t0phi, t0theta);
439 }
440 
441  } // end of if (bestRecHit)
442 
443  } // end of if (nsegm!= 0)
444 
445  // Fill Efficiency plot
446  if (doall_) {
447  HEff4DHit *heff = nullptr;
448 
449  if (wheel == 0) {
450  heff = histograms.hEff_W0.get();
451  } else if (abs(wheel) == 1) {
452  heff = histograms.hEff_W1.get();
453  } else if (abs(wheel) == 2) {
454  heff = histograms.hEff_W2.get();
455 }
456  heff->fill(etaSimSeg, phiSimSeg, xSimSeg, ySimSeg, alphaSimSeg, betaSimSeg, recHitFound, count_seg);
457  histograms.hEff_All->fill(etaSimSeg, phiSimSeg, xSimSeg, ySimSeg, alphaSimSeg, betaSimSeg, recHitFound, count_seg);
458  if (local_) {
459  histograms.hEffWS[abs(wheel)][station-1]->fill(etaSimSeg, phiSimSeg, xSimSeg, ySimSeg, alphaSimSeg, betaSimSeg, recHitFound, count_seg);
460 }
461  }
462  } // End of loop over chambers
463 }
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
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:99
std::map< DTWireId, const PSimHit * > mapMuSimHitsPerWire(const std::map< DTWireId, edm::PSimHitContainer > &simHitWireMap)
Create a map between the Mu SimHits and corresponding MuBarWireId ;.
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:1
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
std::map< DTWireId, edm::PSimHitContainer > mapSimHitsPerWire(const edm::PSimHitContainer &simhits)
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:728
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:699
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
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:545
std::pair< LocalVector, LocalPoint > findMuSimSegmentDirAndPos(const std::pair< const PSimHit *, const PSimHit * > &inAndOutSimHit, const DetId detId, const DTGeometry *muonGeom)
Find direction and position of a segment (in local RF) from outer and inner mu SimHit in the RF of ob...
T eta() const
Definition: PV3DBase.h:76
std::pair< const PSimHit *, const PSimHit * > findMuSimSegment(const std::map< DTWireId, const PSimHit * > &mapWireAndMuSimHit)
Find Innermost and outermost SimHit from Mu in a SL (they identify a simulated segment) ...
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:104

Member Data Documentation

bool DTSegment4DQuality::debug_
private

Definition at line 76 of file DTSegment4DQuality.h.

bool DTSegment4DQuality::doall_
private

Definition at line 72 of file DTSegment4DQuality.h.

bool DTSegment4DQuality::local_
private

Definition at line 73 of file DTSegment4DQuality.h.

edm::InputTag DTSegment4DQuality::segment4DLabel_
private

Definition at line 60 of file DTSegment4DQuality.h.

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

Definition at line 62 of file DTSegment4DQuality.h.

double DTSegment4DQuality::sigmaResAlpha_
private

Definition at line 69 of file DTSegment4DQuality.h.

double DTSegment4DQuality::sigmaResBeta_
private

Definition at line 70 of file DTSegment4DQuality.h.

double DTSegment4DQuality::sigmaResX_
private

Definition at line 65 of file DTSegment4DQuality.h.

double DTSegment4DQuality::sigmaResY_
private

Definition at line 66 of file DTSegment4DQuality.h.

edm::InputTag DTSegment4DQuality::simHitLabel_
private

Definition at line 59 of file DTSegment4DQuality.h.

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

Definition at line 61 of file DTSegment4DQuality.h.