DetLayers/src/MuRingForwardLayer.cc

0001 /** \file
0002  *
0003  *  \author N. Amapane - CERN
0004  */
0005
0006 #include "RecoMuon/DetLayers/interface/MuRingForwardLayer.h"
0007 #include "RecoMuon/DetLayers/interface/MuDetRing.h"
0008 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0009 #include "DataFormats/GeometrySurface/interface/SimpleDiskBounds.h"
0010 #include "TrackingTools/GeomPropagators/interface/Propagator.h"
0011 #include "TrackingTools/DetLayers/interface/MeasurementEstimator.h"
0012 #include "TrackingTools/DetLayers/interface/RBorderFinder.h"
0013 #include "TrackingTools/DetLayers/interface/GeneralBinFinderInR.h"
0014
0015 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0016
0017 #include <algorithm>
0018 #include <iostream>
0019 #include <vector>
0020
0021 using namespace std;
0022
0023 MuRingForwardLayer::MuRingForwardLayer(const vector<const ForwardDetRing*>& rings)
0024     : RingedForwardLayer(false),
0025       theRings(rings),
0026       theComponents(theRings.begin(), theRings.end()),
0027       theBinFinder(nullptr),
0028       isOverlapping(false) {
0029   const std::string metname = "Muon|RecoMuon|RecoMuonDetLayers|MuRingForwardLayer";
0030
0031   // Initial values for R and Z bounds
0032   float theRmin = rings.front()->basicComponents().front()->position().perp();
0033   float theRmax = theRmin;
0034   float theZmin = rings.front()->position().z();
0035   float theZmax = theZmin;
0036
0037   // Cache chamber pointers (the basic components_)
0038   // and find extension in R and Z
0039   for (vector<const ForwardDetRing*>::const_iterator it = rings.begin(); it != rings.end(); it++) {
0040     vector<const GeomDet*> tmp2 = (*it)->basicComponents();
0041     theBasicComps.insert(theBasicComps.end(), tmp2.begin(), tmp2.end());
0042
0043     theRmin = min(theRmin, (*it)->specificSurface().innerRadius());
0044     theRmax = max(theRmax, (*it)->specificSurface().outerRadius());
0045     float halfThick = (*it)->surface().bounds().thickness() / 2.;
0046     float zCenter = (*it)->surface().position().z();
0047     theZmin = min(theZmin, zCenter - halfThick);
0048     theZmax = max(theZmax, zCenter + halfThick);
0049   }
0050
0051   RBorderFinder bf(theRings);
0052   isOverlapping = bf.isROverlapping();
0053   theBinFinder = new GeneralBinFinderInR<double>(bf);
0054
0055   // Build surface
0056
0057   float zPos = (theZmax + theZmin) / 2.;
0058   PositionType pos(0., 0., zPos);
0059   RotationType rot;
0060
0061   setSurface(new BoundDisk(pos, rot, new SimpleDiskBounds(theRmin, theRmax, theZmin - zPos, theZmax - zPos)));
0062
0063   LogTrace(metname) << "Constructing MuRingForwardLayer: " << basicComponents().size() << " Dets " << theRings.size()
0064                     << " Rings "
0065                     << " Z: " << specificSurface().position().z() << " R1: " << specificSurface().innerRadius()
0066                     << " R2: " << specificSurface().outerRadius() << " Per.: " << bf.isRPeriodic()
0067                     << " Overl.: " << bf.isROverlapping();
0068 }
0069
0070 MuRingForwardLayer::~MuRingForwardLayer() {
0071   delete theBinFinder;
0072   for (vector<const ForwardDetRing*>::iterator i = theRings.begin(); i < theRings.end(); i++) {
0073     delete *i;
0074   }
0075 }
0076
0077 vector<GeometricSearchDet::DetWithState> MuRingForwardLayer::compatibleDets(
0078     const TrajectoryStateOnSurface& startingState, const Propagator& prop, const MeasurementEstimator& est) const {
0079   const std::string metname = "Muon|RecoMuon|RecoMuonDetLayers|MuRingForwardLayer";
0080   vector<DetWithState> result;
0081
0082   LogTrace(metname) << "MuRingForwardLayer::compatibleDets,"
0083                     << " R1 " << specificSurface().innerRadius() << " R2: " << specificSurface().outerRadius()
0084                     << " FTS at R: " << startingState.globalPosition().perp();
0085
0086   pair<bool, TrajectoryStateOnSurface> compat = compatible(startingState, prop, est);
0087
0088   if (!compat.first) {
0089     LogTrace(metname) << "     MuRingForwardLayer::compatibleDets: not compatible"
0090                       << " (should not have been selected!)";
0091     return result;
0092   }
0093
0094   TrajectoryStateOnSurface& tsos = compat.second;
0095
0096   int closest = theBinFinder->binIndex(tsos.globalPosition().perp());
0097   const ForwardDetRing* closestRing = theRings[closest];
0098
0099   // Check the closest ring
0100
0101   LogTrace(metname) << "     MuRingForwardLayer::fastCompatibleDets, closestRing: " << closest << " R1 "
0102                     << closestRing->specificSurface().innerRadius()
0103                     << " R2: " << closestRing->specificSurface().outerRadius()
0104                     << " FTS R: " << tsos.globalPosition().perp();
0105   if (tsos.hasError()) {
0106     LogTrace(metname) << " sR: " << sqrt(tsos.localError().positionError().yy())
0107                       << " sX: " << sqrt(tsos.localError().positionError().xx());
0108   }
0109   LogTrace(metname) << endl;
0110
0111   result = closestRing->compatibleDets(tsos, prop, est);
0112
0113   int nclosest = result.size();
0114   int nnextdet = 0;  // MDEBUG counters
0115
0116   //FIXME: if closest is not compatible next cannot be either?
0117
0118   // Use state on layer surface. Note that local coordinates and errors
0119   // are the same on the layer and on all rings surfaces, since
0120   // all BoundDisks are centered in 0,0 and have the same rotation.
0121   // CAVEAT: if the rings are not at the same Z, the local position and error
0122   // will be "Z-projected" to the rings. This is a fairly good approximation.
0123   // However in this case additional propagation will be done when calling
0124   // compatibleDets.
0125   GlobalPoint startPos = tsos.globalPosition();
0126   LocalPoint nextPos(surface().toLocal(startPos));
0127
0128   for (unsigned int idet = closest + 1; idet < theRings.size(); idet++) {
0129     bool inside = false;
0130     if (tsos.hasError()) {
0131       inside = theRings[idet]->specificSurface().bounds().inside(nextPos, tsos.localError().positionError());
0132     } else {
0133       inside = theRings[idet]->specificSurface().bounds().inside(nextPos);
0134     }
0135     if (inside) {
0136       LogTrace(metname) << "     MuRingForwardLayer::fastCompatibleDets:NextRing" << idet << " R1 "
0137                         << theRings[idet]->specificSurface().innerRadius()
0138                         << " R2: " << theRings[idet]->specificSurface().outerRadius() << " FTS R " << nextPos.perp();
0139       nnextdet++;
0140       vector<DetWithState> nextRodDets = theRings[idet]->compatibleDets(tsos, prop, est);
0141       if (!nextRodDets.empty()) {
0142         result.insert(result.end(), nextRodDets.begin(), nextRodDets.end());
0143       } else {
0144         break;
0145       }
0146     }
0147   }
0148
0149   for (int idet = closest - 1; idet >= 0; idet--) {
0150     bool inside = false;
0151     if (tsos.hasError()) {
0152       inside = theRings[idet]->specificSurface().bounds().inside(nextPos, tsos.localError().positionError());
0153     } else {
0154       inside = theRings[idet]->specificSurface().bounds().inside(nextPos);
0155     }
0156     if (inside) {
0157       LogTrace(metname) << "     MuRingForwardLayer::fastCompatibleDets:PreviousRing:" << idet << " R1 "
0158                         << theRings[idet]->specificSurface().innerRadius()
0159                         << " R2: " << theRings[idet]->specificSurface().outerRadius() << " FTS R " << nextPos.perp();
0160       nnextdet++;
0161       vector<DetWithState> nextRodDets = theRings[idet]->compatibleDets(tsos, prop, est);
0162       if (!nextRodDets.empty()) {
0163         result.insert(result.end(), nextRodDets.begin(), nextRodDets.end());
0164       } else {
0165         break;
0166       }
0167     }
0168   }
0169
0170   LogTrace(metname) << "     MuRingForwardLayer::fastCompatibleDets: found: " << result.size()
0171                     << " on closest: " << nclosest << " # checked rings: " << 1 + nnextdet;
0172
0173   return result;
0174 }
0175
0176 vector<DetGroup> MuRingForwardLayer::groupedCompatibleDets(const TrajectoryStateOnSurface& startingState,
0177                                                            const Propagator& prop,
0178                                                            const MeasurementEstimator& est) const {
0179   // FIXME should return only 1 group
0180   cout << "dummy implementation of MuRingForwardLayer::groupedCompatibleDets()" << endl;
0181   return vector<DetGroup>();
0182 }
0183
0184 GeomDetEnumerators::SubDetector MuRingForwardLayer::subDetector() const { return theBasicComps.front()->subDetector(); }
0185
0186 const vector<const GeometricSearchDet*>& MuRingForwardLayer::components() const { return theComponents; }