DetLayers/src/BarrelDetLayer.cc

0001 #include "TrackingTools/DetLayers/interface/BarrelDetLayer.h"
0002 #include "DataFormats/GeometrySurface/interface/Surface.h"
0003 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
0004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0005 #include "DataFormats/GeometrySurface/interface/SimpleCylinderBounds.h"
0006 #include "DataFormats/GeometrySurface/interface/BoundingBox.h"
0007 #include "FWCore/Utilities/interface/Likely.h"
0008
0009 using namespace std;
0010
0011 BarrelDetLayer::~BarrelDetLayer() {}
0012
0013 //--- Extension of the interface
0014 void BarrelDetLayer::setSurface(BoundCylinder* cp) { theCylinder = cp; }
0015
0016 bool BarrelDetLayer::contains(const Local3DPoint& p) const { return surface().bounds().inside(p); }
0017
0018 void BarrelDetLayer::initialize() { setSurface(computeSurface()); }
0019
0020 //--- protected methods
0021 BoundCylinder* BarrelDetLayer::computeSurface() {
0022   vector<const GeomDet*> comps = basicComponents();
0023
0024   // Find extension in Z
0025   float theRmin = comps.front()->position().perp();
0026   float theRmax = theRmin;
0027   float theZmin = comps.front()->position().z();
0028   float theZmax = theZmin;
0029   for (vector<const GeomDet*>::const_iterator deti = comps.begin(); deti != comps.end(); deti++) {
0030     vector<GlobalPoint> corners = BoundingBox().corners(dynamic_cast<const Plane&>((*deti)->surface()));
0031     for (vector<GlobalPoint>::const_iterator ic = corners.begin(); ic != corners.end(); ic++) {
0032       float r = ic->perp();
0033       float z = ic->z();
0034       theRmin = min(theRmin, r);
0035       theRmax = max(theRmax, r);
0036       theZmin = min(theZmin, z);
0037       theZmax = max(theZmax, z);
0038     }
0039     // in addition to the corners we have to check the middle of the
0040     // det +/- thickness/2
0041     // , since the min  radius for some barrel dets is reached there
0042     float rdet = (**deti).position().perp();
0043     float thick = (**deti).surface().bounds().thickness();
0044     theRmin = min(theRmin, rdet - thick / 2.F);
0045     theRmax = max(theRmax, rdet + thick / 2.F);
0046   }
0047
0048   // By default the barrel layers are positioned at the center of the
0049   // global frame, and the axes of their local frame coincide with
0050   // those of the global grame (z along the cylinder axis)
0051   PositionType pos(0., 0., 0.);
0052   RotationType rot;
0053
0054   auto scp = new SimpleCylinderBounds(theRmin, theRmax, theZmin, theZmax);
0055   return new Cylinder(Cylinder::computeRadius(*scp), pos, rot, scp);
0056 }
0057
0058 pair<bool, TrajectoryStateOnSurface> BarrelDetLayer::compatible(const TrajectoryStateOnSurface& ts,
0059                                                                 const Propagator& prop,
0060                                                                 const MeasurementEstimator&) const {
0061   if UNLIKELY (theCylinder == nullptr)
0062     edm::LogError("DetLayers") << "ERROR: BarrelDetLayer::compatible() is used before the layer surface is initialized";
0063   // throw an exception? which one?
0064
0065   TrajectoryStateOnSurface myState = prop.propagate(ts, specificSurface());
0066   if UNLIKELY (!myState.isValid())
0067     return make_pair(false, myState);
0068
0069   // check z assuming symmetric bounds around position().z()
0070   auto z0 = std::abs(myState.globalPosition().z() - specificSurface().position().z());
0071   auto deltaZ = 0.5f * bounds().length();
0072   if (z0 < deltaZ)
0073     return make_pair(true, myState);
0074
0075   // take into account the thickness of the layer
0076   deltaZ += 0.5f * bounds().thickness() * std::abs(myState.globalDirection().z()) / myState.globalDirection().perp();
0077
0078   // take also into account the error on the predicted state
0079   const float nSigma = 3.;
0080   if (myState.hasError())
0081     deltaZ += nSigma * sqrt(myState.cartesianError().position().czz());
0082   //
0083   //  check z again
0084   return make_pair(z0 < deltaZ, myState);
0085 }