CMS 3D CMS Logo

List of all members | Classes | Public Types | Public Member Functions | Private Member Functions | Private Attributes | Friends
ClusterShapeHitFilter Class Reference

#include <ClusterShapeHitFilter.h>

Classes

struct  PixelData
 
struct  StripData
 

Public Types

typedef TrajectoryFilter::Record Record
 

Public Member Functions

 ClusterShapeHitFilter (const TrackerGeometry *theTracker_, const TrackerTopology *theTkTopol_, const MagneticField *theMagneticField_, const SiPixelLorentzAngle *theSiPixelLorentzAngle_, const SiStripLorentzAngle *theSiStripLorentzAngle_, const std::string &pixelShapeFile_, const std::string &pixelShapeFileL1_)
 
bool getSizes (const SiPixelRecHit &recHit, const LocalVector &ldir, const SiPixelClusterShapeCache &clusterShapeCache, int &part, ClusterData::ArrayType &meas, std::pair< float, float > &predr, PixelData const *pd=0) const
 
bool getSizes (DetId detId, const SiStripCluster &cluster, const LocalPoint &lpos, const LocalVector &ldir, int &meas, float &pred) const
 
bool getSizes (const SiStripRecHit2D &recHit, const LocalPoint &lpos, const LocalVector &ldir, int &meas, float &pred) const
 
bool isCompatible (const SiPixelRecHit &recHit, const LocalVector &ldir, const SiPixelClusterShapeCache &clusterShapeCache, PixelData const *pd=0) const
 
bool isCompatible (const SiPixelRecHit &recHit, const GlobalVector &gdir, const SiPixelClusterShapeCache &clusterShapeCache, PixelData const *pd=0) const
 
bool isCompatible (DetId detId, const SiStripCluster &cluster, const LocalPoint &lpos, const LocalVector &ldir) const
 
bool isCompatible (DetId detId, const SiStripCluster &cluster, const LocalVector &ldir) const
 
bool isCompatible (DetId detId, const SiStripCluster &cluster, const GlobalPoint &gpos, const GlobalVector &gdir) const
 
bool isCompatible (DetId detId, const SiStripCluster &cluster, const GlobalVector &gdir) const
 
bool isCompatible (const SiStripRecHit2D &recHit, const LocalPoint &lpos, const LocalVector &ldir) const
 
bool isCompatible (const SiStripRecHit2D &recHit, const LocalVector &ldir) const
 
bool isCompatible (const SiStripRecHit2D &recHit, const GlobalPoint &gpos, const GlobalVector &gdir) const
 
bool isCompatible (const SiStripRecHit2D &recHit, const GlobalVector &gdir) const
 
void setChargeCuts (bool cutOnPixelCharge, float minGoodPixelCharge, bool cutOnStripCharge, float minGoodStripCharge)
 
void setShapeCuts (bool cutOnPixelShape, bool cutOnStripShape)
 
 ~ClusterShapeHitFilter ()
 

Private Member Functions

bool checkClusterCharge (DetId detId, const SiStripCluster &cluster, const LocalVector &ldir) const
 
bool checkClusterCharge (DetId detId, const SiPixelCluster &cluster, const LocalVector &ldir) const
 
 ClusterShapeHitFilter (std::string const &f1, std::string const &f2)
 
void fillPixelData ()
 
void fillStripData ()
 
std::pair< float, float > getCotangent (const PixelGeomDetUnit *pixelDet) const
 
float getCotangent (const ClusterShapeHitFilter::StripData &sd, const LocalPoint &p) const
 
std::pair< float, float > getDrift (const PixelGeomDetUnit *pixelDet) const
 
float getDrift (const StripGeomDetUnit *stripDet) const
 
const PixelDatagetpd (const SiPixelRecHit &recHit, PixelData const *pd=0) const
 
const StripDatagetsd (DetId id) const
 
bool isNormalOriented (const GeomDetUnit *geomDet) const
 
void loadPixelLimits (std::string const &file, PixelLimits *plim)
 
void loadStripLimits ()
 

Private Attributes

bool cutOnPixelCharge_
 
bool cutOnPixelShape_
 
bool cutOnStripCharge_
 
bool cutOnStripShape_
 
float minGoodPixelCharge_
 
float minGoodStripCharge_
 
std::unordered_map< unsigned int, PixelDatapixelData
 
PixelLimits pixelLimits [PixelKeys::N+1]
 
PixelLimits pixelLimitsL1 [PixelKeys::N+1]
 
std::unordered_map< unsigned int, StripDatastripData
 
StripLimits stripLimits [StripKeys::N+1]
 
float theAngle [6]
 
const MagneticFieldtheMagneticField
 
const SiPixelLorentzAngletheSiPixelLorentzAngle
 
const SiStripLorentzAngletheSiStripLorentzAngle
 
const TrackerTopologytheTkTopol
 
const TrackerGeometrytheTracker
 

Friends

int test::ClusterShapeHitFilterTest::test ()
 

Detailed Description

Definition at line 153 of file ClusterShapeHitFilter.h.

Member Typedef Documentation

Definition at line 177 of file ClusterShapeHitFilter.h.

Constructor & Destructor Documentation

ClusterShapeHitFilter::ClusterShapeHitFilter ( const TrackerGeometry theTracker_,
const TrackerTopology theTkTopol_,
const MagneticField theMagneticField_,
const SiPixelLorentzAngle theSiPixelLorentzAngle_,
const SiStripLorentzAngle theSiStripLorentzAngle_,
const std::string &  pixelShapeFile_,
const std::string &  pixelShapeFileL1_ 
)

Definition at line 37 of file ClusterShapeHitFilter.cc.

44  : theTracker(theTracker_),
45  theTkTopol(theTkTopol_),
46  theMagneticField(theMagneticField_),
47  theSiPixelLorentzAngle(theSiPixelLorentzAngle_),
48  theSiStripLorentzAngle(theSiStripLorentzAngle_)
49 {
50  // Load pixel limits
51  loadPixelLimits(pixelShapeFile_,pixelLimits);
52  loadPixelLimits(pixelShapeFileL1_,pixelLimitsL1);
53  fillPixelData();
54 
55  // Load strip limits
57  fillStripData();
60 }
const MagneticField * theMagneticField
const SiStripLorentzAngle * theSiStripLorentzAngle
PixelLimits pixelLimitsL1[PixelKeys::N+1]
PixelLimits pixelLimits[PixelKeys::N+1]
const TrackerGeometry * theTracker
void loadPixelLimits(std::string const &file, PixelLimits *plim)
const SiPixelLorentzAngle * theSiPixelLorentzAngle
const TrackerTopology * theTkTopol
ClusterShapeHitFilter::~ClusterShapeHitFilter ( )

Definition at line 63 of file ClusterShapeHitFilter.cc.

64 {
65 }
ClusterShapeHitFilter::ClusterShapeHitFilter ( std::string const &  f1,
std::string const &  f2 
)
inlineprivate

Definition at line 261 of file ClusterShapeHitFilter.h.

261  {
264  loadStripLimits();
265  }
PixelLimits pixelLimitsL1[PixelKeys::N+1]
PixelLimits pixelLimits[PixelKeys::N+1]
void loadPixelLimits(std::string const &file, PixelLimits *plim)

Member Function Documentation

bool ClusterShapeHitFilter::checkClusterCharge ( DetId  detId,
const SiStripCluster cluster,
const LocalVector ldir 
) const
private

Definition at line 452 of file ClusterShapeHitFilter.cc.

References siStripClusterTools::chargePerCM().

453 {
454  return siStripClusterTools::chargePerCM(detId, cluster, ldir) > minGoodStripCharge_;
455 }
float chargePerCM(DetId detid, Iter a, Iter b)
bool ClusterShapeHitFilter::checkClusterCharge ( DetId  detId,
const SiPixelCluster cluster,
const LocalVector ldir 
) const
private

Definition at line 458 of file ClusterShapeHitFilter.cc.

References siStripClusterTools::chargePerCM().

459 {
460  return siStripClusterTools::chargePerCM(detId, cluster, ldir) > minGoodPixelCharge_;
461 }
float chargePerCM(DetId detid, Iter a, Iter b)
void ClusterShapeHitFilter::fillPixelData ( )
private

Definition at line 137 of file ClusterShapeHitFilter.cc.

References ClusterShapeHitFilter::PixelData::cotangent, ClusterShapeHitFilter::PixelData::det, ClusterShapeHitFilter::PixelData::drift, GeomDet::geographicalId(), ClusterShapeHitFilter::PixelData::layer, and ClusterShapeHitFilter::PixelData::part.

137  {
138 
139  //barrel
140  for (auto det : theTracker->detsPXB()) {
141  // better not to fail..
142  const PixelGeomDetUnit * pixelDet =
143  dynamic_cast<const PixelGeomDetUnit*>(det);
144  assert(pixelDet);
145  PixelData & pd = pixelData[pixelDet->geographicalId()];
146  pd.det = pixelDet;
147  pd.part=0;
148  pd.layer = theTkTopol->pxbLayer(pixelDet->geographicalId());
149  pd.cotangent=getCotangent(pixelDet);
150  pd.drift=getDrift(pixelDet);
151  }
152 
153  //endcap
154  for (auto det : theTracker->detsPXF()) {
155  // better not to fail..
156  const PixelGeomDetUnit * pixelDet =
157  dynamic_cast<const PixelGeomDetUnit*>(det);
158  assert(pixelDet);
159  PixelData & pd = pixelData[pixelDet->geographicalId()];
160  pd.det = pixelDet;
161  pd.part=1;
162  pd.layer=0;
163  pd.cotangent=getCotangent(pixelDet);
164  pd.drift=getDrift(pixelDet);
165  }
166 
167 }
std::unordered_map< unsigned int, PixelData > pixelData
const TrackerGeometry * theTracker
std::pair< float, float > getCotangent(const PixelGeomDetUnit *pixelDet) const
const DetContainer & detsPXB() const
DetId geographicalId() const
The label of this GeomDet.
Definition: GeomDet.h:79
unsigned int pxbLayer(const DetId &id) const
std::pair< float, float > getDrift(const PixelGeomDetUnit *pixelDet) const
const TrackerTopology * theTkTopol
const DetContainer & detsPXF() const
void ClusterShapeHitFilter::fillStripData ( )
private

Definition at line 170 of file ClusterShapeHitFilter.cc.

References Surface::bounds(), GeomDet::geographicalId(), getCotangent(), mps_fire::i, GeomDet::index(), createfilelist::int, GeomDetType::isTrackerStrip(), PFRecoTauDiscriminationByIsolation_cfi::offset, AlCaHLTBitMon_ParallelJobs::p, GeomDet::specificSurface(), GeomDetEnumerators::tkDetEnum, StripGeomDetUnit::topology(), and StripGeomDetUnit::type().

170  {
171  // copied from StripCPE (FIXME maybe we should move all this in LocalReco)
172  auto const & geom_ = *theTracker;
173  auto const & dus = geom_.detUnits();
174  auto offset=dus.size();
175  for(unsigned int i=1;i<7;++i) {
176  if(geom_.offsetDU(GeomDetEnumerators::tkDetEnum[i]) != dus.size() &&
177  dus[geom_.offsetDU(GeomDetEnumerators::tkDetEnum[i])]->type().isTrackerStrip()) {
178  if(geom_.offsetDU(GeomDetEnumerators::tkDetEnum[i]) < offset) offset = geom_.offsetDU(GeomDetEnumerators::tkDetEnum[i]);
179  }
180  }
181 
182  for (auto i=offset; i!=dus.size();++i) {
183  const StripGeomDetUnit * stripdet=(const StripGeomDetUnit*)(dus[i]);
184  assert(stripdet->index()==int(i));
185  assert(stripdet->type().isTrackerStrip()); // not pixel
186  auto const & bounds = stripdet->specificSurface().bounds();
187  auto detid = stripdet->geographicalId();
188  auto & p = stripData[detid];
189  p.det = stripdet;
190  p.topology=(StripTopology*)(&stripdet->topology());
191  p.drift = getDrift(stripdet);
192  p.thickness = bounds.thickness();
193  p.nstrips = p.topology->nstrips();
194  }
195 }
std::unordered_map< unsigned int, StripData > stripData
const DetContainer & detUnits() const override
Returm a vector of all GeomDet.
const Bounds & bounds() const
Definition: Surface.h:120
bool isTrackerStrip() const
Definition: GeomDetType.cc:24
const TrackerGeometry * theTracker
DetId geographicalId() const
The label of this GeomDet.
Definition: GeomDet.h:79
int index() const
Definition: GeomDet.h:99
SubDetector tkDetEnum[8]
const GeomDetType & type() const override
const Topology & topology() const override
Returns a reference to the strip proxy topology.
std::pair< float, float > getDrift(const PixelGeomDetUnit *pixelDet) const
const Plane & specificSurface() const
Same as surface(), kept for backward compatibility.
Definition: GeomDet.h:45
pair< float, float > ClusterShapeHitFilter::getCotangent ( const PixelGeomDetUnit pixelDet) const
private

Definition at line 199 of file ClusterShapeHitFilter.cc.

References Surface::bounds(), PixelTopology::pitch(), PixelGeomDetUnit::specificTopology(), GeomDet::surface(), and Bounds::thickness().

Referenced by fillStripData().

200 {
201  pair<float,float> cotangent;
202 
203  cotangent.first = pixelDet->surface().bounds().thickness() /
204  pixelDet->specificTopology().pitch().first;
205  cotangent.second = pixelDet->surface().bounds().thickness() /
206  pixelDet->specificTopology().pitch().second;
207 
208  return cotangent;
209 }
virtual std::pair< float, float > pitch() const =0
const Bounds & bounds() const
Definition: Surface.h:120
const Plane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:42
virtual float thickness() const =0
virtual const PixelTopology & specificTopology() const
Returns a reference to the pixel proxy topology.
float ClusterShapeHitFilter::getCotangent ( const ClusterShapeHitFilter::StripData sd,
const LocalPoint p 
) const
private

Definition at line 213 of file ClusterShapeHitFilter.cc.

References getDrift(), StripTopology::localPitch(), ClusterShapeHitFilter::StripData::thickness, and ClusterShapeHitFilter::StripData::topology.

214 {
215  // FIXME may be problematic in case of RadialStriptolopgy
216  return sd.thickness/sd.topology->localPitch(pos);
217 }
virtual float localPitch(const LocalPoint &) const =0
pair< float, float > ClusterShapeHitFilter::getDrift ( const PixelGeomDetUnit pixelDet) const
private

Definition at line 221 of file ClusterShapeHitFilter.cc.

References dir, GeomDet::geographicalId(), GloballyPositioned< T >::position(), DetId::rawId(), GeomDet::surface(), toLocal(), PV3DBase< T, PVType, FrameType >::x(), and PV3DBase< T, PVType, FrameType >::y().

Referenced by getCotangent().

222 {
223  LocalVector lBfield =
224  (pixelDet->surface()).toLocal(
226  pixelDet->surface().position()));
227 
228  double theTanLorentzAnglePerTesla =
230  pixelDet->geographicalId().rawId());
231 
232  pair<float,float> dir;
233  dir.first = - theTanLorentzAnglePerTesla * lBfield.y();
234  dir.second = theTanLorentzAnglePerTesla * lBfield.x();
235 
236  return dir;
237 }
const MagneticField * theMagneticField
constexpr uint32_t rawId() const
get the raw id
Definition: DetId.h:50
T y() const
Definition: PV3DBase.h:63
LocalVector toLocal(const reco::Track::Vector &v, const Surface &s)
const Plane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:42
DetId geographicalId() const
The label of this GeomDet.
Definition: GeomDet.h:79
const SiPixelLorentzAngle * theSiPixelLorentzAngle
virtual GlobalVector inTesla(const GlobalPoint &gp) const =0
Field value ad specified global point, in Tesla.
float getLorentzAngle(const uint32_t &) const
dbl *** dir
Definition: mlp_gen.cc:35
T x() const
Definition: PV3DBase.h:62
const PositionType & position() const
float ClusterShapeHitFilter::getDrift ( const StripGeomDetUnit stripDet) const
private

Definition at line 240 of file ClusterShapeHitFilter.cc.

References dir, GeomDet::geographicalId(), isNormalOriented(), GloballyPositioned< T >::position(), DetId::rawId(), GeomDet::surface(), toLocal(), and PV3DBase< T, PVType, FrameType >::y().

242 {
243  LocalVector lBfield =
244  (stripDet->surface()).toLocal(
246  stripDet->surface().position()));
247 
248  double theTanLorentzAnglePerTesla =
250  stripDet->geographicalId().rawId());
251 
252  float dir = theTanLorentzAnglePerTesla * lBfield.y();
253 
254  return dir;
255 }
const MagneticField * theMagneticField
const SiStripLorentzAngle * theSiStripLorentzAngle
constexpr uint32_t rawId() const
get the raw id
Definition: DetId.h:50
T y() const
Definition: PV3DBase.h:63
LocalVector toLocal(const reco::Track::Vector &v, const Surface &s)
const Plane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:42
float getLorentzAngle(const uint32_t &) const
DetId geographicalId() const
The label of this GeomDet.
Definition: GeomDet.h:79
virtual GlobalVector inTesla(const GlobalPoint &gp) const =0
Field value ad specified global point, in Tesla.
dbl *** dir
Definition: mlp_gen.cc:35
const PositionType & position() const
const PixelData& ClusterShapeHitFilter::getpd ( const SiPixelRecHit recHit,
PixelData const *  pd = 0 
) const
inlineprivate

Definition at line 267 of file ClusterShapeHitFilter.h.

References TrackingRecHit::geographicalId(), and AlCaHLTBitMon_ParallelJobs::p.

267  {
268  if (pd) return *pd;
269  // Get detector
270  DetId id = recHit.geographicalId();
271  auto p = pixelData.find(id);
272  return (*p).second;
273  }
std::unordered_map< unsigned int, PixelData > pixelData
Definition: DetId.h:18
DetId geographicalId() const
const StripData& ClusterShapeHitFilter::getsd ( DetId  id) const
inlineprivate

Definition at line 275 of file ClusterShapeHitFilter.h.

References FrontierConditions_GlobalTag_cff::file, AlCaHLTBitMon_ParallelJobs::p, sd, and AlCaHLTBitMon_QueryRunRegistry::string.

275 {return stripData.find(id)->second;}
std::unordered_map< unsigned int, StripData > stripData
bool ClusterShapeHitFilter::getSizes ( const SiPixelRecHit recHit,
const LocalVector ldir,
const SiPixelClusterShapeCache clusterShapeCache,
int &  part,
ClusterData::ArrayType meas,
std::pair< float, float > &  predr,
PixelData const *  pd = 0 
) const

Definition at line 279 of file ClusterShapeHitFilter.cc.

References edm::VecArray< T, N >::capacity(), edm::VecArray< T, N >::clear(), SiPixelRecHit::cluster(), ClusterShapeHitFilter::PixelData::cotangent, data, ClusterShapeHitFilter::PixelData::det, SoftLeptonByDistance_cfi::distance, shallow::drift(), ClusterShapeHitFilter::PixelData::drift, SiPixelClusterShapeCache::get(), isCompatible(), SiPixelClusterShapeData::isComplete(), SiPixelClusterShapeData::isStraight(), ClusterShapeHitFilter::PixelData::part, edm::VecArray< T, N >::push_back_unchecked(), alignCSCRings::s, SiPixelClusterShapeData::size(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by isCompatible(), isNormalOriented(), and PixelClusterShapeExtractor::processRec().

283 {
284  // Get detector
285  const PixelData & pd = getpd(recHit,ipd);
286 
287  // Get shape information
288  const SiPixelClusterShapeData& data = clusterShapeCache.get(recHit.cluster(), pd.det);
289  bool usable = (data.isStraight() && data.isComplete());
290 
291  // Usable?
292  //if(usable)
293  {
294  part = pd.part;
295 
296  // Predicted size
297  pred.first = ldir.x() / ldir.z();
298  pred.second = ldir.y() / ldir.z();
299 
300  SiPixelClusterShapeData::Range sizeRange = data.size();
301  if(sizeRange.first->second < 0)
302  pred.second = - pred.second;
303 
304  meas.clear();
305  assert(meas.capacity() >= std::distance(sizeRange.first, sizeRange.second));
306  for(auto s=sizeRange.first; s != sizeRange.second; ++s) {
307  meas.push_back_unchecked(*s);
308  }
309  if(sizeRange.first->second < 0) {
310  for(auto& s: meas)
311  s.second = -s.second;
312  }
313 
314  // Take out drift
315  std::pair<float,float> const & drift = pd.drift;
316  pred.first += drift.first;
317  pred.second += drift.second;
318 
319  // Apply cotangent
320  std::pair<float,float> const & cotangent = pd.cotangent;
321  pred.first *= cotangent.first;
322  pred.second *= cotangent.second;
323  }
324 
325  // Usable?
326  return usable;
327 }
const PixelData & getpd(const SiPixelRecHit &recHit, PixelData const *pd=0) const
LocalVector drift(const StripGeomDetUnit *, const MagneticField &, const SiStripLorentzAngle &)
Definition: ShallowTools.cc:37
T y() const
Definition: PV3DBase.h:63
T z() const
Definition: PV3DBase.h:64
std::pair< const_iterator, const_iterator > Range
void clear()
Definition: VecArray.h:73
part
Definition: HCALResponse.h:20
ClusterRef cluster() const
Definition: SiPixelRecHit.h:49
char data[epos_bytes_allocation]
Definition: EPOS_Wrapper.h:82
void push_back_unchecked(const T &value)
Definition: VecArray.h:83
SiPixelClusterShapeData get(const ClusterRef &cluster, const PixelGeomDetUnit *pixDet) const
T x() const
Definition: PV3DBase.h:62
static constexpr size_type capacity() noexcept
Definition: VecArray.h:71
bool ClusterShapeHitFilter::getSizes ( DetId  detId,
const SiStripCluster cluster,
const LocalPoint lpos,
const LocalVector ldir,
int &  meas,
float &  pred 
) const

Definition at line 377 of file ClusterShapeHitFilter.cc.

References SiStripCluster::amplitudes(), shallow::drift(), SiStripCluster::firstStrip(), isCompatible(), AlCaHLTBitMon_ParallelJobs::p, PV3DBase< T, PVType, FrameType >::x(), and PV3DBase< T, PVType, FrameType >::z().

379 {
380  // Get detector
381  auto const & p=getsd(id);
382 
383  // Measured width
384  meas = cluster.amplitudes().size();
385 
386  // Usable?
387  int fs = cluster.firstStrip();
388  int ns = p.nstrips;
389  // bool usable = (fs > 1 && fs + meas - 1 < ns);
390  bool usable = (fs >= 1) & (fs + meas <= ns+1);
391 
392  // Usable?
393  //if(usable)
394  {
395  // Predicted width
396  pred = ldir.x() / ldir.z();
397 
398  // Take out drift
399  float drift = p.drift;;
400  pred += drift;
401 
402  // Apply cotangent
403  pred *= getCotangent(p,lpos);
404  }
405 
406  return usable;
407 }
LocalVector drift(const StripGeomDetUnit *, const MagneticField &, const SiStripLorentzAngle &)
Definition: ShallowTools.cc:37
uint16_t firstStrip() const
const StripData & getsd(DetId id) const
std::pair< float, float > getCotangent(const PixelGeomDetUnit *pixelDet) const
T z() const
Definition: PV3DBase.h:64
T x() const
Definition: PV3DBase.h:62
const std::vector< uint8_t > & amplitudes() const
bool ClusterShapeHitFilter::getSizes ( const SiStripRecHit2D recHit,
const LocalPoint lpos,
const LocalVector ldir,
int &  meas,
float &  pred 
) const
inline

Definition at line 215 of file ClusterShapeHitFilter.h.

References TrackingRecHit::geographicalId(), and TrackerSingleRecHit::stripCluster().

216  {
217  return getSizes(recHit.geographicalId(), recHit.stripCluster(), lpos, ldir, meas, pred);
218  }
SiStripCluster const & stripCluster() const
bool getSizes(const SiPixelRecHit &recHit, const LocalVector &ldir, const SiPixelClusterShapeCache &clusterShapeCache, int &part, ClusterData::ArrayType &meas, std::pair< float, float > &predr, PixelData const *pd=0) const
DetId geographicalId() const
bool ClusterShapeHitFilter::isCompatible ( const SiPixelRecHit recHit,
const LocalVector ldir,
const SiPixelClusterShapeCache clusterShapeCache,
PixelData const *  pd = 0 
) const

Definition at line 330 of file ClusterShapeHitFilter.cc.

References SiPixelRecHit::cluster(), TrackingRecHit::geographicalId(), PixelKeys::isValid(), crabWrapper::key, ClusterShapeHitFilter::PixelData::layer, and funct::m.

Referenced by PixelClusterShapeSeedComparitor::compatibleHit(), getSizes(), MultiHitGeneratorFromChi2::hitSets(), isCompatible(), ClusterShapeTrackFilter::operator()(), and ClusterShapeTrajectoryFilter::toBeContinued().

333 {
334  // Get detector
335  if (cutOnPixelCharge_ && (!checkClusterCharge(recHit.geographicalId(), *(recHit.cluster()), ldir))) return false;
336  if (!cutOnPixelShape_) return true;
337 
338  const PixelData & pd = getpd(recHit,ipd);
339 
340  int part;
342  pair<float,float> pred;
343 
344  PixelLimits const * pl = pd.layer==1 ? pixelLimitsL1 : pixelLimits;
345  if(getSizes(recHit, ldir, clusterShapeCache, part,meas, pred,&pd))
346  {
347  for(const auto& m: meas)
348  {
349  PixelKeys key(part, m.first, m.second);
350  if (!key.isValid()) return true; // FIXME original logic
351  if (pl[key].isInside(pred)) return true;
352  }
353  // none of the choices worked
354  return false;
355  }
356  // not usable
357  return true;
358 }
const PixelData & getpd(const SiPixelRecHit &recHit, PixelData const *pd=0) const
bool checkClusterCharge(DetId detId, const SiStripCluster &cluster, const LocalVector &ldir) const
PixelLimits pixelLimitsL1[PixelKeys::N+1]
PixelLimits pixelLimits[PixelKeys::N+1]
part
Definition: HCALResponse.h:20
ClusterRef cluster() const
Definition: SiPixelRecHit.h:49
bool getSizes(const SiPixelRecHit &recHit, const LocalVector &ldir, const SiPixelClusterShapeCache &clusterShapeCache, int &part, ClusterData::ArrayType &meas, std::pair< float, float > &predr, PixelData const *pd=0) const
DetId geographicalId() const
bool ClusterShapeHitFilter::isCompatible ( const SiPixelRecHit recHit,
const GlobalVector gdir,
const SiPixelClusterShapeCache clusterShapeCache,
PixelData const *  pd = 0 
) const

Definition at line 361 of file ClusterShapeHitFilter.cc.

References ClusterShapeHitFilter::PixelData::det, getSizes(), and GeomDet::toLocal().

364 {
365  // Get detector
366  const PixelData & pd = getpd(recHit,ipd);
367 
368  LocalVector ldir =pd.det->toLocal(gdir);
369 
370  return isCompatible(recHit, ldir, clusterShapeCache, &pd);
371 }
const PixelData & getpd(const SiPixelRecHit &recHit, PixelData const *pd=0) const
bool isCompatible(const SiPixelRecHit &recHit, const LocalVector &ldir, const SiPixelClusterShapeCache &clusterShapeCache, PixelData const *pd=0) const
bool ClusterShapeHitFilter::isCompatible ( DetId  detId,
const SiStripCluster cluster,
const LocalPoint lpos,
const LocalVector ldir 
) const

Definition at line 412 of file ClusterShapeHitFilter.cc.

References isCompatible(), StripKeys::isValid(), and crabWrapper::key.

413 {
414  int meas;
415  float pred;
416 
417  if (cutOnStripCharge_ && (!checkClusterCharge(detId, cluster, ldir))) return false;
418  if (!cutOnStripShape_) return true;
419 
420  if(getSizes(detId, cluster, lpos, ldir, meas, pred))
421  {
422  StripKeys key(meas);
423  if (key.isValid())
424  return stripLimits[key].isInside(pred);
425  }
426 
427  // Not usable or no limits
428  return true;
429 }
bool checkClusterCharge(DetId detId, const SiStripCluster &cluster, const LocalVector &ldir) const
StripLimits stripLimits[StripKeys::N+1]
bool getSizes(const SiPixelRecHit &recHit, const LocalVector &ldir, const SiPixelClusterShapeCache &clusterShapeCache, int &part, ClusterData::ArrayType &meas, std::pair< float, float > &predr, PixelData const *pd=0) const
bool ClusterShapeHitFilter::isCompatible ( DetId  detId,
const SiStripCluster cluster,
const LocalVector ldir 
) const
inline

Definition at line 223 of file ClusterShapeHitFilter.h.

225  {
226  return isCompatible(detId, cluster, LocalPoint(0,0,0), ldir);
227  }
Point3DBase< Scalar, LocalTag > LocalPoint
Definition: Definitions.h:32
bool isCompatible(const SiPixelRecHit &recHit, const LocalVector &ldir, const SiPixelClusterShapeCache &clusterShapeCache, PixelData const *pd=0) const
bool ClusterShapeHitFilter::isCompatible ( DetId  detId,
const SiStripCluster cluster,
const GlobalPoint gpos,
const GlobalVector gdir 
) const

Definition at line 434 of file ClusterShapeHitFilter.cc.

References isCompatible(), GeomDet::toLocal(), and PV3DBase< T, PVType, FrameType >::z().

435 {
436  const GeomDet *det = getsd(detId).det;
437  LocalVector ldir = det->toLocal(gdir);
438  LocalPoint lpos = det->toLocal(gpos);
439  // now here we do the transformation
440  lpos -= ldir * lpos.z()/ldir.z();
441  return isCompatible(detId, cluster, lpos, ldir);
442 }
LocalPoint toLocal(const GlobalPoint &gp) const
Conversion to the R.F. of the GeomDet.
Definition: GeomDet.h:69
const StripData & getsd(DetId id) const
T z() const
Definition: PV3DBase.h:64
bool isCompatible(const SiPixelRecHit &recHit, const LocalVector &ldir, const SiPixelClusterShapeCache &clusterShapeCache, PixelData const *pd=0) const
bool ClusterShapeHitFilter::isCompatible ( DetId  detId,
const SiStripCluster cluster,
const GlobalVector gdir 
) const

Definition at line 444 of file ClusterShapeHitFilter.cc.

445 {
446  return isCompatible(detId, cluster, getsd(detId).det->toLocal(gdir));
447 }
const StripData & getsd(DetId id) const
bool isCompatible(const SiPixelRecHit &recHit, const LocalVector &ldir, const SiPixelClusterShapeCache &clusterShapeCache, PixelData const *pd=0) const
bool ClusterShapeHitFilter::isCompatible ( const SiStripRecHit2D recHit,
const LocalPoint lpos,
const LocalVector ldir 
) const
inline

Definition at line 238 of file ClusterShapeHitFilter.h.

References TrackingRecHit::geographicalId(), and TrackerSingleRecHit::stripCluster().

240  {
241  return isCompatible(recHit.geographicalId(), recHit.stripCluster(), lpos, ldir);
242  }
bool isCompatible(const SiPixelRecHit &recHit, const LocalVector &ldir, const SiPixelClusterShapeCache &clusterShapeCache, PixelData const *pd=0) const
SiStripCluster const & stripCluster() const
DetId geographicalId() const
bool ClusterShapeHitFilter::isCompatible ( const SiStripRecHit2D recHit,
const LocalVector ldir 
) const
inline

Definition at line 243 of file ClusterShapeHitFilter.h.

References TrackingRecHit::geographicalId(), and TrackerSingleRecHit::stripCluster().

244  {
245  return isCompatible(recHit.geographicalId(), recHit.stripCluster(), ldir);
246  }
bool isCompatible(const SiPixelRecHit &recHit, const LocalVector &ldir, const SiPixelClusterShapeCache &clusterShapeCache, PixelData const *pd=0) const
SiStripCluster const & stripCluster() const
DetId geographicalId() const
bool ClusterShapeHitFilter::isCompatible ( const SiStripRecHit2D recHit,
const GlobalPoint gpos,
const GlobalVector gdir 
) const
inline

Definition at line 247 of file ClusterShapeHitFilter.h.

References TrackingRecHit::geographicalId(), and TrackerSingleRecHit::stripCluster().

249  {
250  return isCompatible(recHit.geographicalId(), recHit.stripCluster(), gpos, gdir);
251  }
bool isCompatible(const SiPixelRecHit &recHit, const LocalVector &ldir, const SiPixelClusterShapeCache &clusterShapeCache, PixelData const *pd=0) const
SiStripCluster const & stripCluster() const
DetId geographicalId() const
bool ClusterShapeHitFilter::isCompatible ( const SiStripRecHit2D recHit,
const GlobalVector gdir 
) const
inline

Definition at line 252 of file ClusterShapeHitFilter.h.

References TrackingRecHit::geographicalId(), and TrackerSingleRecHit::stripCluster().

253  {
254  return isCompatible(recHit.geographicalId(), recHit.stripCluster(), gdir);
255  }
bool isCompatible(const SiPixelRecHit &recHit, const LocalVector &ldir, const SiPixelClusterShapeCache &clusterShapeCache, PixelData const *pd=0) const
SiStripCluster const & stripCluster() const
DetId geographicalId() const
bool ClusterShapeHitFilter::isNormalOriented ( const GeomDetUnit geomDet) const
private

Definition at line 259 of file ClusterShapeHitFilter.cc.

References getSizes(), GeomDetType::isBarrel(), PV3DBase< T, PVType, FrameType >::perp(), makeMuonMisalignmentScenario::rot, GeomDet::toGlobal(), GeomDet::type(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by getDrift().

260 {
261  if(geomDet->type().isBarrel())
262  { // barrel
263  float perp0 = geomDet->toGlobal( Local3DPoint(0.,0.,0.) ).perp();
264  float perp1 = geomDet->toGlobal( Local3DPoint(0.,0.,1.) ).perp();
265  return (perp1 > perp0);
266  }
267  else
268  { // endcap
269  float rot = geomDet->toGlobal( LocalVector (0.,0.,1.) ).z();
270  float pos = geomDet->toGlobal( Local3DPoint(0.,0.,0.) ).z();
271  return (rot * pos > 0);
272  }
273 }
bool isBarrel() const
Definition: GeomDetType.cc:13
T perp() const
Definition: PV3DBase.h:72
virtual const GeomDetType & type() const
Definition: GeomDet.cc:85
GlobalPoint toGlobal(const Local2DPoint &lp) const
Conversion to the global R.F. from the R.F. of the GeomDet.
Definition: GeomDet.h:54
T z() const
Definition: PV3DBase.h:64
Point3DBase< float, LocalTag > Local3DPoint
Definition: LocalPoint.h:9
void ClusterShapeHitFilter::loadPixelLimits ( std::string const &  file,
PixelLimits plim 
)
private

Definition at line 68 of file ClusterShapeHitFilter.cc.

References b, edmIntegrityCheck::d, PixelLimits::data, PVValHelper::dx, PVValHelper::dy, f, gen::k, crabWrapper::key, and LogTrace.

69 {
70  edm::FileInPath fileInPath(file.c_str());
71  ifstream inFile(fileInPath.fullPath().c_str());
72 
73 
74  while(inFile.eof() == false)
75  {
76  int part,dx,dy;
77 
78  inFile >> part; // 0or 1
79  inFile >> dx; // 0 to 10
80  inFile >> dy; // 0 to 15 ...
81 
82  const PixelKeys key(part,dx,dy);
83  auto & pl = plim[key];
84 
85  for(int b = 0; b<2 ; b++) // branch
86  for(int d = 0; d<2 ; d++) // direction
87  for(int k = 0; k<2 ; k++) // lower and upper
88  inFile >> pl.data[b][d][k];
89 
90 
91  double f;
92  int d;
93 
94  inFile >> f; // density
95  inFile >> d; // points
96  inFile >> f; // density
97  inFile >> d; // points
98  }
99 
100  inFile.close();
101 
102  LogTrace("MinBiasTracking|ClusterShapeHitFilter")
103  << " [ClusterShapeHitFilter] pixel-cluster-shape filter loaded";
104  }
double f[11][100]
#define LogTrace(id)
int k[5][pyjets_maxn]
part
Definition: HCALResponse.h:20
double b
Definition: hdecay.h:120
float data[2][2][2]
void ClusterShapeHitFilter::loadStripLimits ( )
private

Definition at line 107 of file ClusterShapeHitFilter.cc.

References b, PVValHelper::dx, edm::FileInPath::fullPath(), gen::k, crabWrapper::key, and LogTrace.

108 {
109  // Load strip
111  fileInPath("RecoPixelVertexing/PixelLowPtUtilities/data/stripShape.par");
112  ifstream inFile(fileInPath.fullPath().c_str());
113 
114 
115  while(inFile.eof() == false)
116  {
117  int dx;
118  inFile >> dx;
119 
120  StripKeys key(dx);
121  auto & sl = stripLimits[key];
122 
123  for(int b = 0; b<2 ; b++) // branch
124  for(int k = 0; k<2 ; k++) // lower and upper
125  inFile >> sl.data[b][k];
126 
127  }
128 
129  inFile.close();
130 
131  LogTrace("MinBiasTracking|ClusterShapeHitFilter")
132  << " [ClusterShapeHitFilter] strip-cluster-width filter loaded";
133 }
StripLimits stripLimits[StripKeys::N+1]
#define LogTrace(id)
int k[5][pyjets_maxn]
double b
Definition: hdecay.h:120
void ClusterShapeHitFilter::setChargeCuts ( bool  cutOnPixelCharge,
float  minGoodPixelCharge,
bool  cutOnStripCharge,
float  minGoodStripCharge 
)
inline
void ClusterShapeHitFilter::setShapeCuts ( bool  cutOnPixelShape,
bool  cutOnStripShape 
)
inline

Definition at line 189 of file ClusterShapeHitFilter.h.

189  {
190  cutOnPixelShape_ = cutOnPixelShape; cutOnStripShape_ = cutOnStripShape;}

Friends And Related Function Documentation

Member Data Documentation

bool ClusterShapeHitFilter::cutOnPixelCharge_
private

Definition at line 307 of file ClusterShapeHitFilter.h.

bool ClusterShapeHitFilter::cutOnPixelShape_
private

Definition at line 309 of file ClusterShapeHitFilter.h.

bool ClusterShapeHitFilter::cutOnStripCharge_
private

Definition at line 307 of file ClusterShapeHitFilter.h.

bool ClusterShapeHitFilter::cutOnStripShape_
private

Definition at line 309 of file ClusterShapeHitFilter.h.

float ClusterShapeHitFilter::minGoodPixelCharge_
private

Definition at line 308 of file ClusterShapeHitFilter.h.

float ClusterShapeHitFilter::minGoodStripCharge_
private

Definition at line 308 of file ClusterShapeHitFilter.h.

std::unordered_map<unsigned int, PixelData> ClusterShapeHitFilter::pixelData
private

Definition at line 298 of file ClusterShapeHitFilter.h.

PixelLimits ClusterShapeHitFilter::pixelLimits[PixelKeys::N+1]
private

Definition at line 301 of file ClusterShapeHitFilter.h.

PixelLimits ClusterShapeHitFilter::pixelLimitsL1[PixelKeys::N+1]
private

Definition at line 302 of file ClusterShapeHitFilter.h.

std::unordered_map<unsigned int, StripData> ClusterShapeHitFilter::stripData
private

Definition at line 299 of file ClusterShapeHitFilter.h.

StripLimits ClusterShapeHitFilter::stripLimits[StripKeys::N+1]
private

Definition at line 304 of file ClusterShapeHitFilter.h.

float ClusterShapeHitFilter::theAngle[6]
private

Definition at line 306 of file ClusterShapeHitFilter.h.

const MagneticField* ClusterShapeHitFilter::theMagneticField
private

Definition at line 293 of file ClusterShapeHitFilter.h.

const SiPixelLorentzAngle* ClusterShapeHitFilter::theSiPixelLorentzAngle
private

Definition at line 295 of file ClusterShapeHitFilter.h.

const SiStripLorentzAngle* ClusterShapeHitFilter::theSiStripLorentzAngle
private

Definition at line 296 of file ClusterShapeHitFilter.h.

const TrackerTopology* ClusterShapeHitFilter::theTkTopol
private

Definition at line 292 of file ClusterShapeHitFilter.h.

const TrackerGeometry* ClusterShapeHitFilter::theTracker
private

Definition at line 291 of file ClusterShapeHitFilter.h.