Perform the real analysis.
82 map<DTChamberId, PSimHitContainer> simHitsPerCh;
83 for (
const auto &
simHit : *simHits) {
87 simHitsPerCh[chamberId].push_back(
simHit);
96 cout <<
"[DTSegment2DSLPhiQuality]**Warning: no 4D Segments with label: " <<
segment4DLabel_ 97 <<
" in this event, skipping!" << endl;
104 for (chamberId = segment4Ds->id_begin(); chamberId != segment4Ds->id_end(); ++chamberId) {
112 int nMuSimHit = muSimHitPerWire.size();
113 if (nMuSimHit == 0 || nMuSimHit == 1) {
114 if (
debug_ && nMuSimHit == 1) {
115 cout <<
"[DTSegment2DSLPhiQuality] Only " << nMuSimHit <<
" mu SimHit in this chamber, skipping!" << endl;
120 cout <<
"=== Chamber " << (*chamberId) <<
" has " << nMuSimHit <<
" SimHits" << endl;
127 pair<LocalVector, LocalPoint> dirAndPosSimSegm =
130 LocalVector simSegmLocalDir = dirAndPosSimSegm.first;
131 LocalPoint simSegmLocalPos = dirAndPosSimSegm.second;
137 float posSimSeg = simSegmLocalPos.
x();
139 float etaSimSeg = simSegmGlobalPos.
eta();
140 float phiSimSeg = simSegmGlobalPos.
phi();
143 cout <<
" Simulated segment: local direction " << simSegmLocalDir << endl
144 <<
" local position " << simSegmLocalPos << endl
145 <<
" angle " << angleSimSeg << endl;
150 bool recHitFound =
false;
152 int nsegm =
distance(range.first, range.second);
154 cout <<
" Chamber: " << *chamberId <<
" has " << nsegm <<
" 4D segments" << endl;
164 bool bestRecHitFound =
false;
165 double deltaAlpha = 99999;
170 if ((*segment4D).dimension() != 4) {
172 cout <<
"[DTSegment2DSLPhiQuality]***Error: This is not 4D " 181 if ((*phiSegment2D).dimension() != 2) {
183 cout <<
"[DTSegment2DQuality]***Error: This is not 2D segment!!!" << endl;
189 LocalVector recSegDirection = (*phiSegment2D).localDirection();
193 cout <<
" RecSegment direction: " << recSegDirection << endl
194 <<
" position : " << (*phiSegment2D).localPosition() << endl
195 <<
" alpha : " << recSegAlpha << endl;
198 if (fabs(recSegAlpha - angleSimSeg) < deltaAlpha) {
199 deltaAlpha = fabs(recSegAlpha - angleSimSeg);
200 bestRecHit = &(*phiSegment2D);
201 bestRecHitFound =
true;
205 if (bestRecHitFound) {
215 fabs(bestRecHitLocalPos.
x() - posSimSeg) < 5 *
sigmaResPos_) {
223 bestRecHitLocalPos.
x(),
226 sqrt(bestRecHitLocalPosErr.
xx()),
227 sqrt(bestRecHitLocalDirErr.
xx()));
233 histograms.h2DHitEff_SuperPhi->fill(etaSimSeg, phiSimSeg, posSimSeg, angleSimSeg, recHitFound);
edm::EDGetTokenT< edm::PSimHitContainer > simHitToken_
std::pair< double, double > findSegmentAlphaAndBeta(const LocalVector &direction)
LocalError localPositionError() const override
local position error in SL frame
std::map< DTWireId, edm::PSimHitContainer > mapSimHitsPerWire(const edm::PSimHitContainer &simhits)
std::pair< const_iterator, const_iterator > range
iterator range
const DTChamber * chamber(const DTChamberId &id) const
Return a DTChamber given its id.
GlobalPoint toGlobal(const Local2DPoint &lp) const
Conversion to the global R.F. from the R.F. of the GeomDet.
Geom::Phi< T > phi() const
def setup(process, global_tag, zero_tesla=False)
C::const_iterator const_iterator
constant access iterator type
edm::EDGetTokenT< DTRecSegment4DCollection > segment4DToken_
std::map< DTWireId, const PSimHit * > mapMuSimHitsPerWire(const std::map< DTWireId, edm::PSimHitContainer > &simHitWireMap)
Create a map between the Mu SimHits and corresponding MuBarWireId ;.
std::pair< const PSimHit *, const PSimHit * > findMuSimSegment(const std::map< DTWireId, const PSimHit * > &mapWireAndMuSimHit)
LocalPoint localPosition() const override
local position in SL frame
edm::InputTag segment4DLabel_
std::pair< LocalVector, LocalPoint > findMuSimSegmentDirAndPos(const std::pair< const PSimHit *, const PSimHit * > &inAndOutSimHit, const DetId detId, const DTGeometry *muonGeom)
LocalVector localDirection() const override
the local direction in SL frame
std::vector< PSimHit > PSimHitContainer
LocalError localDirectionError() const override
the local direction error (xx,xy,yy) in SL frame: only xx is not 0.