CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
LowPtClusterShapeSeedComparitor.cc
Go to the documentation of this file.
1 #include <cmath>
2 
18 
19 namespace {
21 
22  inline float sqr(float x) { return x * x; }
23 
24  /*****************************************************************************/
25  inline float areaParallelogram(const Vector2D& a, const Vector2D& b) { return a.x() * b.y() - a.y() * b.x(); }
26 
27  /*****************************************************************************/
28 
29  inline bool getGlobalDirs(GlobalPoint const* g, GlobalVector* globalDirs) {
30  // Determine circle
31  CircleFromThreePoints circle(g[0], g[1], g[2]);
32 
33  float curvature = circle.curvature();
34  if (0.f == curvature) {
35  LogDebug("LowPtClusterShapeSeedComparitor")
36  << "the curvature is null:"
37  << "\n point1: " << g[0] << "\n point2: " << g[1] << "\n point3: " << g[2];
38  return false;
39  }
40 
41  // Get 2d points
42  Vector2D p[3];
43  Vector2D c = circle.center();
44  for (int i = 0; i != 3; i++)
45  p[i] = g[i].basicVector().xy() - c;
46 
47  float area = std::abs(areaParallelogram(p[1] - p[0], p[1]));
48 
49  float a12 = std::asin(std::min(area * curvature * curvature, 1.f));
50 
51  float slope = (g[1].z() - g[0].z()) / a12;
52 
53  // Calculate globalDirs
54 
55  float cotTheta = slope * curvature;
56  float sinTheta = 1.f / std::sqrt(1.f + sqr(cotTheta));
57  float cosTheta = cotTheta * sinTheta;
58 
59  if (areaParallelogram(p[0], p[1]) < 0)
60  sinTheta = -sinTheta;
61 
62  for (int i = 0; i != 3; i++) {
63  Vector2D vl = p[i] * (curvature * sinTheta);
64  globalDirs[i] = GlobalVector(-vl.y(), vl.x(), cosTheta);
65  }
66  return true;
67  }
68 
69  /*****************************************************************************/
70 
71  inline void getGlobalPos(const SeedingHitSet& hits, GlobalPoint* globalPoss) {
72  for (unsigned int i = 0; i != hits.size(); ++i)
73  globalPoss[i] = hits[i]->globalPosition();
74  }
75 
76 } // namespace
77 
78 /*****************************************************************************/
81  : thePixelClusterShapeCacheToken(
82  iC.consumes<SiPixelClusterShapeCache>(ps.getParameter<edm::InputTag>("clusterShapeCacheSrc"))),
83  theShapeFilterLabel_(ps.getParameter<std::string>("clusterShapeHitFilter")),
84  clusterShapeHitFilterESToken_(iC.esConsumes(edm::ESInputTag("", theShapeFilterLabel_))),
85  trackerTopologyESToken_(iC.esConsumes()) {}
86 
87 /*****************************************************************************/
92 }
93 
95 //(const reco::Track* track, const vector<const TrackingRecHit *> & recHits) const
96 {
97  assert(hits.size() == 3);
98 
99  if (clusterShapeHitFilter_ == nullptr)
100  throw cms::Exception("LogicError") << "LowPtClusterShapeSeedComparitor: init(EventSetup) method was not called";
101 
102  // Get global positions
103  GlobalPoint globalPoss[3];
104  getGlobalPos(hits, globalPoss);
105 
106  // Get global directions
107  GlobalVector globalDirs[3];
108 
109  bool ok = getGlobalDirs(globalPoss, globalDirs);
110 
111  // Check whether shape of pixel cluster is compatible
112  // with local track direction
113 
114  if (!ok) {
115  LogDebug("LowPtClusterShapeSeedComparitor") << "curvarture 0:"
116  << "\nnHits: " << hits.size() << " will say the seed is good anyway.";
117  return true;
118  }
119 
120  for (int i = 0; i < 3; i++) {
121  const SiPixelRecHit* pixelRecHit = dynamic_cast<const SiPixelRecHit*>(hits[i]->hit());
122 
123  if (!pixelRecHit) {
124  edm::LogError("LowPtClusterShapeSeedComparitor") << "this is not a pixel cluster";
125  ok = false;
126  break;
127  }
128 
129  if (!pixelRecHit->isValid()) {
130  ok = false;
131  break;
132  }
133 
134  LogDebug("LowPtClusterShapeSeedComparitor") << "about to compute compatibility."
135  << "hit ptr: " << pixelRecHit << "global direction:" << globalDirs[i];
136 
137  if (!clusterShapeHitFilter_->isCompatible(*pixelRecHit, globalDirs[i], *thePixelClusterShapeCache)) {
138  LogTrace("LowPtClusterShapeSeedComparitor")
139  << " clusShape is not compatible" << HitInfo::getInfo(*hits[i]->hit(), trackerTopology_);
140 
141  ok = false;
142  break;
143  }
144  }
145 
146  return ok;
147 }
ESGetTokenH3DDVariant esConsumes(std::string const &Record, edm::ConsumesCollector &)
Definition: DeDxTools.cc:283
T const & getData(const ESGetToken< T, R > &iToken) const noexcept(false)
Definition: EventSetup.h:119
bool compatible(const SeedingHitSet &hits) const override
static const double slope[3]
Log< level::Error, false > LogError
assert(be >=bs)
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e g
Definition: Activities.doc:4
T curvature(T InversePt, const MagneticField &field)
#define LogTrace(id)
Basic2DVector< double >::MathVector Vector2D
Definition: LinkByRecHit.cc:7
bool isValid() const
T sqrt(T t)
Definition: SSEVec.h:19
const edm::ESGetToken< ClusterShapeHitFilter, CkfComponentsRecord > clusterShapeHitFilterESToken_
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
LowPtClusterShapeSeedComparitor(const edm::ParameterSet &ps, edm::ConsumesCollector &iC)
double f[11][100]
edm::Handle< SiPixelClusterShapeCache > thePixelClusterShapeCache
edm::EDGetTokenT< SiPixelClusterShapeCache > thePixelClusterShapeCacheToken
void init(const edm::Event &e, const edm::EventSetup &es) override
double b
Definition: hdecay.h:120
HLT enums.
double a
Definition: hdecay.h:121
Square< F >::type sqr(const F &f)
Definition: Square.h:14
float x
const ClusterShapeHitFilter * clusterShapeHitFilter_
static std::string getInfo(const DetId &id, const TrackerTopology *tTopo)
Definition: HitInfo.cc:19
const edm::ESGetToken< TrackerTopology, TrackerTopologyRcd > trackerTopologyESToken_
Global3DVector GlobalVector
Definition: GlobalVector.h:10
Our base class.
Definition: SiPixelRecHit.h:23
#define LogDebug(id)