#include <StripSubClusterShapeTrajectoryFilter.h>

Inheritance diagram for StripSubClusterShapeFilterBase:

Public Member Functions
	StripSubClusterShapeFilterBase (const edm::ParameterSet &iConfig, edm::ConsumesCollector &iC)

virtual	~StripSubClusterShapeFilterBase ()

Protected Member Functions
void	setEventBase (const edm::Event &, const edm::EventSetup &)

bool	testLastHit (const TrackingRecHit *hit, const TrajectoryStateOnSurface &tsos, bool mustProject=false) const

bool	testLastHit (const TrackingRecHit *hit, const GlobalPoint &gpos, const GlobalVector &gdir, bool mustProject=false) const

Protected Attributes
std::string	label_

std::array< std::array< uint8_t, 10 >, 7 >	layerMask_

uint32_t	maxNSat_

float	maxTrimmedSizeDiffNeg_

float	maxTrimmedSizeDiffPos_

float	seedCutMIPs_

float	seedCutSN_

float	subclusterCutMIPs_

float	subclusterCutSN_

float	subclusterWindow_

edm::ESHandle< ClusterShapeHitFilter >	theFilter

edm::ESHandle< SiStripNoises >	theNoise

edm::ESHandle< TrackerTopology >	theTopology

edm::ESHandle< TrackerGeometry >	theTracker

uint8_t	trimMaxADC_

float	trimMaxFracNeigh_

float	trimMaxFracTotal_

Detailed Description

Definition at line 26 of file StripSubClusterShapeTrajectoryFilter.h.

Constructor & Destructor Documentation

StripSubClusterShapeFilterBase::StripSubClusterShapeFilterBase	(	const edm::ParameterSet &	iConfig,
		edm::ConsumesCollector &	iC
	)

Definition at line 125 of file StripSubClusterShapeTrajectoryFilter.cc.

References begin, end, edm::ParameterSet::exists(), edm::ParameterSet::existsAs(), lumiContext::fill, edm::ParameterSet::getParameter(), and mps_fire::i.

                                                            :
    label_(iCfg.getUntrackedParameter<std::string>("label","")),
    maxNSat_(iCfg.getParameter<uint32_t>("maxNSat")),
    trimMaxADC_(iCfg.getParameter<double>("trimMaxADC")),
    trimMaxFracTotal_(iCfg.getParameter<double>("trimMaxFracTotal")),
    trimMaxFracNeigh_(iCfg.getParameter<double>("trimMaxFracNeigh")),
    maxTrimmedSizeDiffPos_(iCfg.getParameter<double>("maxTrimmedSizeDiffPos")),
    maxTrimmedSizeDiffNeg_(iCfg.getParameter<double>("maxTrimmedSizeDiffNeg")),
    subclusterWindow_(iCfg.getParameter<double>("subclusterWindow")),
    seedCutMIPs_(iCfg.getParameter<double>("seedCutMIPs")),
    seedCutSN_(iCfg.getParameter<double>("seedCutSN")),
    subclusterCutMIPs_(iCfg.getParameter<double>("subclusterCutMIPs")),
    subclusterCutSN_(iCfg.getParameter<double>("subclusterCutSN"))
 #ifdef StripSubClusterShapeFilterBase_COUNTERS
    ,called_(0),saturated_(0),test_(0),passTrim_(0),failTooLarge_(0),passSC_(0),failTooNarrow_(0)
 #endif
 {
     if (iCfg.exists("layerMask")) {
         const edm::ParameterSet &iLM = iCfg.getParameter<edm::ParameterSet>("layerMask");
         const char *ndets[4] =  { "TIB", "TID", "TOB", "TEC" };
         const int   idets[4] =  {   3,     4,     5,     6   };
         for (unsigned int i = 0; i < 4; ++i) {
             if (iLM.existsAs<bool>(ndets[i])) {
                 std::fill(layerMask_[idets[i]].begin(), layerMask_[idets[i]].end(), iLM.getParameter<bool>(ndets[i]));
             } else {
                 layerMask_[idets[i]][0] = 2;
                 std::fill(layerMask_[idets[i]].begin()+1, layerMask_[idets[i]].end(), 0);
                 for (uint32_t lay : iLM.getParameter<std::vector<uint32_t>>(ndets[i])) {
                     layerMask_[idets[i]][lay] = 1;
                 }
             }
         }
     } else {
         for (auto & arr : layerMask_) {
             std::fill(arr.begin(), arr.end(), 1);
         }
     }
 }

StripSubClusterShapeFilterBase::~StripSubClusterShapeFilterBase ( )

virtual

Definition at line 165 of file StripSubClusterShapeTrajectoryFilter.cc.

References gather_cfg::cout, and testLastHit().

 {
 #if 0
     std::cout << "StripSubClusterShapeFilterBase " << label_ <<": called        " << called_ << std::endl;
     std::cout << "StripSubClusterShapeFilterBase " << label_ <<": saturated     " << saturated_ << std::endl;
     std::cout << "StripSubClusterShapeFilterBase " << label_ <<": test          " << test_ << std::endl;
     std::cout << "StripSubClusterShapeFilterBase " << label_ <<": failTooNarrow " << failTooNarrow_ << std::endl;
     std::cout << "StripSubClusterShapeFilterBase " << label_ <<": passTrim      " << passTrim_ << std::endl;
     std::cout << "StripSubClusterShapeFilterBase " << label_ <<": passSC        " << passSC_ << std::endl;
     std::cout << "StripSubClusterShapeFilterBase " << label_ <<": failTooLarge  " << failTooLarge_ << std::endl;
 #endif
 }

Member Function Documentation

void StripSubClusterShapeFilterBase::setEventBase	(	const edm::Event &	event,
		const edm::EventSetup &	es
	)

protected

Definition at line 288 of file StripSubClusterShapeTrajectoryFilter.cc.

References edm::EventSetup::get(), and StripSubClusterShapeTrajectoryFilter::toBeContinued().

Referenced by testLastHit().

 {
   // Get tracker geometry
   es.get<TrackerDigiGeometryRecord>().get(theTracker);
 
   es.get<CkfComponentsRecord>().get("ClusterShapeHitFilter",theFilter);
 
   //Retrieve tracker topology from geometry
   es.get<TrackerTopologyRcd>().get(theTopology);
 
   es.get<SiStripNoisesRcd>().get(theNoise);
 
 }

bool StripSubClusterShapeFilterBase::testLastHit	(	const TrackingRecHit *	hit,
		const TrajectoryStateOnSurface &	tsos,
		bool	mustProject = `false`
	)		const

protected

Definition at line 180 of file StripSubClusterShapeTrajectoryFilter.cc.

References TrajectoryStateOnSurface::globalDirection(), and TrajectoryStateOnSurface::globalPosition().

Referenced by StripSubClusterShapeTrajectoryFilter::setEvent(), and ~StripSubClusterShapeFilterBase().

 {
     return testLastHit(hit, tsos.globalPosition(), tsos.globalDirection(), mustProject);
 
 }

bool StripSubClusterShapeFilterBase::testLastHit	(	const TrackingRecHit *	hit,
		const GlobalPoint &	gpos,
		const GlobalVector &	gdir,
		bool	mustProject = `false`
	)		const

protected

Definition at line 186 of file StripSubClusterShapeTrajectoryFilter.cc.

References funct::abs(), SiStripCluster::amplitudes(), begin, Surface::bounds(), f, SiStripCluster::firstStrip(), TrackingRecHit::geographicalId(), mps_fire::i, INC_COUNTER, plotBeamSpotDB::last, dedxEstimators_cff::MeVperADCStrip, SiStripMatchedRecHit2D::monoHit(), gen::n, ProjectedSiStripRecHit2D::originalHit(), packedPFCandidateRefMixer_cfi::pf, DetId::rawId(), setEventBase(), SiStripMatchedRecHit2D::stereoHit(), TrackerSingleRecHit::stripCluster(), DetId::subdetId(), GeomDet::surface(), PFTauMVAInputDiscriminatorTranslator_cfi::test, Bounds::thickness(), GeomDet::toLocal(), and PV3DBase< T, PVType, FrameType >::z().

 {
    const TrackerSingleRecHit *stripHit = nullptr;
    if (typeid(*hit) == typeid(SiStripMatchedRecHit2D)) {
       const SiStripMatchedRecHit2D & mhit = static_cast<const SiStripMatchedRecHit2D &>(*hit);
       SiStripRecHit2D mono = mhit.monoHit();
       SiStripRecHit2D stereo = mhit.stereoHit();
       return testLastHit(&mono, gpos, gdir, true) && testLastHit(&stereo, gpos, gdir, true);
    } else if (typeid(*hit) == typeid(ProjectedSiStripRecHit2D)) {
       const ProjectedSiStripRecHit2D & mhit = static_cast<const ProjectedSiStripRecHit2D &>(*hit);
       const SiStripRecHit2D &orig = mhit.originalHit();
       return testLastHit(&orig,   gpos, gdir, true);
    } else if ((stripHit = dynamic_cast<const TrackerSingleRecHit *>(hit)) != nullptr) {
       DetId detId = hit->geographicalId();
 
       if (layerMask_[detId.subdetId()][0] == 0) {
           return true; // no filtering here
       } else if (layerMask_[detId.subdetId()][0] == 2) {
           unsigned int ilayer = theTopology->layer(detId);          
           if (layerMask_[detId.subdetId()][ilayer] == 0) {
               return true; // no filtering here
           }
       }
 
       const GeomDet *det = theTracker->idToDet(detId);
       LocalVector ldir = det->toLocal(gdir);
       LocalPoint  lpos = det->toLocal(gpos);
       if (mustProject) {
          lpos -= ldir * lpos.z()/ldir.z();
       }
       int hitStrips; float hitPredPos;
       const SiStripCluster &cluster = stripHit->stripCluster();
       bool usable = theFilter->getSizes(detId, cluster, lpos, ldir, hitStrips, hitPredPos);
       if (!usable) return true;
 
       INC_COUNTER(called_)
       const auto  &ampls = cluster.amplitudes();
       
       // pass-through of trivial case 
       if (std::abs(hitPredPos) < 1.5f && hitStrips <= 2) {
         return true;
       }
 
 
       // compute number of consecutive saturated strips 
       // (i.e. with adc count >= 254, see SiStripCluster class for documentation)
       unsigned int thisSat = (ampls[0] >= 254), maxSat = thisSat;
       for (unsigned int i = 1, n = ampls.size(); i < n; ++i) {
           if (ampls[i] >= 254) {
               thisSat++;
           } else if (thisSat > 0) {
               maxSat = std::max<int>(maxSat, thisSat);
               thisSat = 0;
           }
       }
       if (thisSat > 0) {
           maxSat = std::max<int>(maxSat, thisSat);
       }
       if (maxSat >= maxNSat_) {
           INC_COUNTER(saturated_)
           return true;
       }
      
       // trimming
       INC_COUNTER(test_)
       unsigned int hitStripsTrim = ampls.size();
       int sum = std::accumulate(ampls.begin(), ampls.end(), 0);
       uint8_t trimCut = std::min<uint8_t>(trimMaxADC_, std::floor(trimMaxFracTotal_ * sum));
       auto begin = ampls.begin();
       auto last = ampls.end()-1;
       while (hitStripsTrim > 1 && (*begin < std::max<uint8_t>(trimCut, trimMaxFracNeigh_*(*(begin+1)))) ) { hitStripsTrim--; ++begin; }
       while (hitStripsTrim > 1 && (*last  < std::max<uint8_t>(trimCut, trimMaxFracNeigh_*(*(last -1)))) ) { hitStripsTrim--; --last; }
 
       if (hitStripsTrim < std::floor(std::abs(hitPredPos)-maxTrimmedSizeDiffNeg_)) {
           INC_COUNTER(failTooNarrow_)
           return false;
       } else if (hitStripsTrim <= std::ceil(std::abs(hitPredPos)+maxTrimmedSizeDiffPos_)) {
           INC_COUNTER(passTrim_)
           return true;
       } 
 
       const StripGeomDetUnit* stripDetUnit = dynamic_cast<const StripGeomDetUnit *>(det);
       if (det == nullptr) { edm::LogError("Strip not a StripGeomDetUnit?") << " on " << detId.rawId() << "\n"; return true; }
 
       float MeVperADCStrip = 9.5665E-4; // conversion constant from ADC counts to MeV for the SiStrip detector
       float mip = 3.9 / ( MeVperADCStrip/stripDetUnit->surface().bounds().thickness() ); // 3.9 MeV/cm = ionization in silicon 
       float mipnorm = mip/std::abs(ldir.z());
       ::SlidingPeakFinder pf(std::max<int>(2,std::ceil(std::abs(hitPredPos)+subclusterWindow_)));
       ::PeakFinderTest test(mipnorm, detId(), cluster.firstStrip(), &*theNoise, seedCutMIPs_, seedCutSN_, subclusterCutMIPs_, subclusterCutSN_);
       if (pf.apply(cluster.amplitudes(), test)) {
           INC_COUNTER(passSC_)
           return true;
       } else {
           INC_COUNTER(failTooLarge_)
           return false;
       }
 
     }
     return true; 
 }