CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
SiStripClusters2ApproxClusters.cc
Go to the documentation of this file.
31 
32 #include <vector>
33 #include <memory>
34 
36 public:
38  void produce(edm::Event&, const edm::EventSetup&) override;
39 
40  static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
41 
42 private:
45 
46  unsigned int maxNSat;
47  static constexpr double subclusterWindow_ = .7;
48  static constexpr double seedCutMIPs_ = .35;
49  static constexpr double seedCutSN_ = 7.;
50  static constexpr double subclusterCutMIPs_ = .45;
51  static constexpr double subclusterCutSN_ = 12.;
52 
55 
57 
60 
63 
66 };
67 
69  inputClusters = conf.getParameter<edm::InputTag>("inputClusters");
70  maxNSat = conf.getParameter<unsigned int>("maxSaturatedStrips");
71 
72  clusterToken = consumes<edmNew::DetSetVector<SiStripCluster> >(inputClusters);
73 
74  beamSpot_ = conf.getParameter<edm::InputTag>("beamSpot");
75  beamSpotToken_ = consumes<reco::BeamSpot>(beamSpot_);
76 
78 
81 
82  csfLabel_ = conf.getParameter<std::string>("clusterShapeHitFilterLabel");
84 
86 
87  produces<edmNew::DetSetVector<SiStripApproximateCluster> >();
88 }
89 
91  auto result = std::make_unique<edmNew::DetSetVector<SiStripApproximateCluster> >();
92  const auto& clusterCollection = event.get(clusterToken);
93 
94  auto const beamSpotHandle = event.getHandle(beamSpotToken_);
95  auto const& bs = beamSpotHandle.isValid() ? *beamSpotHandle : reco::BeamSpot();
96  if (not beamSpotHandle.isValid()) {
97  edm::LogError("SiStripClusters2ApproxClusters")
98  << "didn't find a valid beamspot with label \"" << beamSpot_.encode() << "\" -> using (0,0,0)";
99  }
100 
101  const auto& tkGeom = &iSetup.getData(tkGeomToken_);
102  const auto& theFilter = &iSetup.getData(csfToken_);
103  const auto& theNoise_ = &iSetup.getData(stripNoiseToken_);
104 
105  for (const auto& detClusters : clusterCollection) {
107 
108  unsigned int detId = detClusters.id();
109  const GeomDet* det = tkGeom->idToDet(detId);
110  double nApvs = detInfo_.getNumberOfApvsAndStripLength(detId).first;
111  double stripLength = detInfo_.getNumberOfApvsAndStripLength(detId).second;
112  double barycenter_ypos = 0.5 * stripLength;
113 
114  const StripGeomDetUnit* stripDet = dynamic_cast<const StripGeomDetUnit*>(det);
115  float mip = 3.9 / (sistrip::MeVperADCStrip / stripDet->surface().bounds().thickness());
116 
117  for (const auto& cluster : detClusters) {
118  const LocalPoint& lp = LocalPoint(((cluster.barycenter() * 10 / (sistrip::STRIPS_PER_APV * nApvs)) -
119  ((stripDet->surface().bounds().width()) * 0.5f)),
120  barycenter_ypos - (0.5f * stripLength),
121  0.);
122  const GlobalPoint& gpos = det->surface().toGlobal(lp);
123  GlobalPoint beamspot(bs.position().x(), bs.position().y(), bs.position().z());
124  const GlobalVector& gdir = gpos - beamspot;
125  const LocalVector& ldir = det->toLocal(gdir);
126 
127  int hitStrips;
128  float hitPredPos;
129  theFilter->getSizes(detId, cluster, lp, ldir, hitStrips, hitPredPos);
130 
131  bool peakFilter = false;
132  SlidingPeakFinder pf(std::max<int>(2, std::ceil(std::abs(hitPredPos) + subclusterWindow_)));
133  float mipnorm = mip / std::abs(ldir.z());
134  PeakFinderTest test(mipnorm,
135  detId,
136  cluster.firstStrip(),
137  theNoise_,
138  seedCutMIPs_,
139  seedCutSN_,
142  peakFilter = pf.apply(cluster.amplitudes(), test);
143 
144  ff.push_back(SiStripApproximateCluster(cluster, maxNSat, hitPredPos, peakFilter));
145  }
146  }
147 
148  event.put(std::move(result));
149 }
150 
153  desc.add<edm::InputTag>("inputClusters", edm::InputTag("siStripClusters"));
154  desc.add<unsigned int>("maxSaturatedStrips", 3);
155  desc.add<std::string>("clusterShapeHitFilterLabel", "ClusterShapeHitFilter"); // add CSF label
156  desc.add<edm::InputTag>("beamSpot", edm::InputTag("offlineBeamSpot")); // add BeamSpot tag
157  descriptions.add("SiStripClusters2ApproxClusters", desc);
158 }
159 
constexpr int32_t ceil(float num)
ESGetTokenH3DDVariant esConsumes(std::string const &Record, edm::ConsumesCollector &)
Definition: DeDxTools.cc:283
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
Point3DBase< Scalar, LocalTag > LocalPoint
Definition: Definitions.h:30
T const & getData(const ESGetToken< T, R > &iToken) const noexcept(false)
Definition: EventSetup.h:119
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
std::string fullPath() const
Definition: FileInPath.cc:161
LocalPoint toLocal(const GlobalPoint &gp) const
Conversion to the R.F. of the GeomDet.
Definition: GeomDet.h:58
T z() const
Definition: PV3DBase.h:61
std::string encode() const
Definition: InputTag.cc:159
SiStripClusters2ApproxClusters(const edm::ParameterSet &conf)
edm::ESGetToken< TrackerGeometry, TrackerDigiGeometryRecord > tkGeomToken_
edm::ESHandle< SiStripNoises > theNoise_
Log< level::Error, false > LogError
edm::ESGetToken< ClusterShapeHitFilter, CkfComponentsRecord > csfToken_
virtual float thickness() const =0
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
double f[11][100]
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
SiStripDetInfo read(std::string filePath)
void produce(edm::Event &, const edm::EventSetup &) override
GlobalPoint toGlobal(const Point2DBase< Scalar, LocalTag > lp) const
Definition: Surface.h:79
const Plane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:37
const std::pair< unsigned short, double > getNumberOfApvsAndStripLength(uint32_t detId) const
Constants and enumerated types for FED/FEC systems.
void add(std::string const &label, ParameterSetDescription const &psetDescription)
edm::ESGetToken< SiStripNoises, SiStripNoisesRcd > stripNoiseToken_
static const uint16_t STRIPS_PER_APV
static constexpr float MeVperADCStrip
edm::EDGetTokenT< reco::BeamSpot > beamSpotToken_
static constexpr char const *const kDefaultFile
virtual float width() const =0
edm::EDGetTokenT< edmNew::DetSetVector< SiStripCluster > > clusterToken
def move(src, dest)
Definition: eostools.py:511
Definition: event.py:1
const Bounds & bounds() const
Definition: Surface.h:87