CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Chi2ChargeMeasurementEstimatorESProducer.cc
Go to the documentation of this file.
1 
12 
14 
15 #include<limits>
16 
18 
20 
21 
22 namespace {
23 
25 public:
26 
32  explicit Chi2ChargeMeasurementEstimator(double maxChi2, double nSigma,
33  float minGoodPixelCharge, float minGoodStripCharge,
34  float pTChargeCutThreshold) :
35  Chi2MeasurementEstimator( maxChi2, nSigma),
36  minGoodPixelCharge_(minGoodPixelCharge),
37  minGoodStripCharge_(minGoodStripCharge) {
38  if (pTChargeCutThreshold>=0.) pTChargeCutThreshold2_=pTChargeCutThreshold*pTChargeCutThreshold;
39  else pTChargeCutThreshold2_=std::numeric_limits<float>::max();
40  }
41 
42 
43  bool preFilter(const TrajectoryStateOnSurface& ts,
44  const MeasurementEstimator::OpaquePayload & opay) const override;
45 
46 
47  virtual Chi2ChargeMeasurementEstimator* clone() const {
48  return new Chi2ChargeMeasurementEstimator(*this);
49  }
50 private:
51 
52  float minGoodPixelCharge_;
53  float minGoodStripCharge_;
54  float pTChargeCutThreshold2_;
55 
56  bool checkClusterCharge(DetId id, SiStripCluster const & cluster, const TrajectoryStateOnSurface& ts) const {
57  return siStripClusterTools::chargePerCM(id, cluster, ts.localParameters() ) > minGoodStripCharge_;
58 
59  }
60 
61 
62 };
63 
65  const MeasurementEstimator::OpaquePayload & opay) const {
66 
67  // what we got?
68  if (opay.tag != ClusterFilterPayload::myTag) return true; // not mine...
69 
70  auto const & clf = reinterpret_cast<ClusterFilterPayload const &>(opay);
71 
72  if (ts.globalMomentum().perp2()>pTChargeCutThreshold2_) return true;
73 
74  DetId detid = clf.detId;
75  uint32_t subdet = detid.subdetId();
76 
77  if (subdet>2) {
78  return checkClusterCharge(detid, *clf.cluster[0],ts) && ( nullptr==clf.cluster[1] || checkClusterCharge(detid, *clf.cluster[1],ts) ) ;
79 
80  }
81 
82  /* pixel charge not implemented as not used...
83  auto const & thit = static_cast<const SiPixelRecHit &>(hit);
84  thit.cluster()->charge() ...
85 
86  */
87 
88  return true;
89 }
90 
91 }
92 
93 
94 
98 #include <boost/shared_ptr.hpp>
99 
100 namespace {
101 
102 
103 class Chi2ChargeMeasurementEstimatorESProducer: public edm::ESProducer{
104  public:
105  Chi2ChargeMeasurementEstimatorESProducer(const edm::ParameterSet & p);
106  virtual ~Chi2ChargeMeasurementEstimatorESProducer();
107  boost::shared_ptr<Chi2MeasurementEstimatorBase> produce(const TrackingComponentsRecord &);
108  private:
109  boost::shared_ptr<Chi2MeasurementEstimatorBase> _estimator;
110 
111  double maxChi2_;
112  double nSigma_;
113  float minGoodPixelCharge_;
114  float minGoodStripCharge_;
115  float pTChargeCutThreshold_;
116 
117 };
118 
119 Chi2ChargeMeasurementEstimatorESProducer::Chi2ChargeMeasurementEstimatorESProducer(const edm::ParameterSet & pset)
120 {
121  std::string const & myname = pset.getParameter<std::string>("ComponentName");
122  setWhatProduced(this,myname);
123 
124  maxChi2_ = pset.getParameter<double>("MaxChi2");
125  nSigma_ = pset.getParameter<double>("nSigma");
126  minGoodPixelCharge_ = 0;
127  minGoodStripCharge_ = clusterChargeCut(pset);
128  pTChargeCutThreshold_= pset.getParameter<double>("pTChargeCutThreshold");
129 
130 
131 }
132 
133 Chi2ChargeMeasurementEstimatorESProducer::~Chi2ChargeMeasurementEstimatorESProducer() {}
134 
135 boost::shared_ptr<Chi2MeasurementEstimatorBase>
136 Chi2ChargeMeasurementEstimatorESProducer::produce(const TrackingComponentsRecord & iRecord){
137 
138  _estimator = boost::shared_ptr<Chi2MeasurementEstimatorBase>(
139  new Chi2ChargeMeasurementEstimator(maxChi2_,nSigma_,
140  minGoodPixelCharge_, minGoodStripCharge_, pTChargeCutThreshold_));
141  return _estimator;
142 }
143 
144 }
145 
146 
148 DEFINE_FWK_EVENTSETUP_MODULE(Chi2ChargeMeasurementEstimatorESProducer);
149 
T getParameter(std::string const &) const
#define GCC11_FINAL
virtual Chi2MeasurementEstimator * clone() const
virtual bool preFilter(const TrajectoryStateOnSurface &, OpaquePayload const &) const
const LocalTrajectoryParameters & localParameters() const
float chargePerCM(DetId detid, Iter a, Iter b)
float clusterChargeCut(const edm::ParameterSet &conf, const char *name="clusterChargeCut")
T perp2() const
Definition: PV3DBase.h:71
tuple preFilter
Definition: gather_cfg.py:50
int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:37
Definition: DetId.h:18
#define DEFINE_FWK_EVENTSETUP_MODULE(type)
Definition: ModuleFactory.h:60
GlobalVector globalMomentum() const
static constexpr int myTag