diff --git a/RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterInfo.h b/RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterInfo.h new file mode 100644 index 0000000000000..8007b8379ed8e --- /dev/null +++ b/RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterInfo.h @@ -0,0 +1,68 @@ +#ifndef SISTRIPCLUSTERIZER_SISTRIPCLUSTERINFO_H +#define SISTRIPCLUSTERIZER_SISTRIPCLUSTERINFO_H + +#include "FWCore/Framework/interface/Frameworkfwd.h" +#include "FWCore/Framework/interface/ESHandle.h" +#include "DataFormats/SiStripCluster/interface/SiStripCluster.h" +#include + +class SiStripNoises; +class SiStripGain; +class SiStripQuality; + + +class SiStripClusterInfo { + + public: + + SiStripClusterInfo(const SiStripCluster& cluster, + const edm::EventSetup& es, + std::string qualityLabel=""); + + const SiStripCluster * cluster() const {return cluster_ptr;} + + uint32_t detId() const {return cluster()->geographicalId();} + uint16_t width() const {return cluster()->amplitudes().size();} + uint16_t firstStrip() const {return cluster()->firstStrip();} + float baryStrip() const {return cluster()->barycenter();} + uint16_t maxStrip() const {return firstStrip() + maxIndex();} + float variance() const; + + const std::vector& stripCharges() const {return cluster()->amplitudes();} + std::vector stripGains() const; + std::vector stripNoises() const; + std::vector stripNoisesRescaledByGain() const; + std::vector stripQualitiesBad() const; + + uint16_t charge() const {return accumulate( stripCharges().begin(), stripCharges().end(), uint16_t(0));} + uint8_t maxCharge() const {return * max_element(stripCharges().begin(), stripCharges().end());} + uint16_t maxIndex() const {return max_element(stripCharges().begin(), stripCharges().end()) - stripCharges().begin();} + std::pair chargeLR() const; + + float noise() const { return calculate_noise(stripNoises());} + float noiseRescaledByGain() const { return calculate_noise(stripNoisesRescaledByGain());} + + float signalOverNoise() const { return charge()/noiseRescaledByGain(); } + + bool IsAnythingBad() const; + bool IsApvBad() const; + bool IsFiberBad() const; + bool IsModuleBad() const; + bool IsModuleUsable() const; + + std::vector reclusterize(const edm::ParameterSet&) const; + + private: + + float calculate_noise(const std::vector&) const; + + const SiStripCluster* cluster_ptr; + const edm::EventSetup& es; + edm::ESHandle noiseHandle; + edm::ESHandle gainHandle; + edm::ESHandle qualityHandle; + std::string qualityLabel; + +}; + +#endif diff --git a/RecoLocalTracker/SiStripClusterizer/src/SiStripClusterInfo.cc b/RecoLocalTracker/SiStripClusterizer/src/SiStripClusterInfo.cc new file mode 100644 index 0000000000000..44f171d010147 --- /dev/null +++ b/RecoLocalTracker/SiStripClusterizer/src/SiStripClusterInfo.cc @@ -0,0 +1,163 @@ +#include "RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterInfo.h" + +#include "CondFormats/DataRecord/interface/SiStripNoisesRcd.h" +#include "CalibTracker/Records/interface/SiStripGainRcd.h" +#include "CalibTracker/Records/interface/SiStripQualityRcd.h" +#include "FWCore/Framework/interface/EventSetup.h" +#include "RecoLocalTracker/SiStripClusterizer/interface/StripClusterizerAlgorithmFactory.h" +#include "RecoLocalTracker/SiStripClusterizer/interface/StripClusterizerAlgorithm.h" +#include + +SiStripClusterInfo::SiStripClusterInfo(const SiStripCluster& cluster, + const edm::EventSetup& setup, + std::string quality) + : cluster_ptr(&cluster), + es(setup), + qualityLabel(quality) { + es.get().get(noiseHandle); + es.get().get(gainHandle); + es.get().get(qualityLabel,qualityHandle); +} + +std::pair SiStripClusterInfo:: +chargeLR() const { + std::vector::const_iterator + begin( stripCharges().begin() ), + end( stripCharges().end() ), + max; max = max_element(begin,end); + return std::make_pair( accumulate(begin, max, uint16_t(0) ), + accumulate(max+1, end, uint16_t(0) ) ); +} + + +float SiStripClusterInfo:: +variance() const { + float q(0), x1(0), x2(0); + for(std::vector::const_iterator + begin(stripCharges().begin()), end(stripCharges().end()), it(begin); + it!=end; ++it) { + unsigned i = it-begin; + q += (*it); + x1 += (*it) * (i+0.5); + x2 += (*it) * (i*i+i+1./3); + } + return (x2 - x1*x1/q ) / q; +} + +std::vector SiStripClusterInfo:: +stripNoisesRescaledByGain() const { + std::vector noises = stripNoises(); + std::vector gains = stripGains(); + transform(noises.begin(), noises.end(), gains.begin(), + noises.begin(), + std::divides()); + return noises; +} + +std::vector SiStripClusterInfo:: +stripNoises() const { + SiStripNoises::Range detNoiseRange = noiseHandle->getRange(cluster()->geographicalId()); + + std::vector noises; + for(size_t i=0; i < width(); i++){ + noises.push_back( noiseHandle->getNoise( firstStrip()+i, detNoiseRange) ); + } + return noises; +} + +std::vector SiStripClusterInfo:: +stripGains() const { + SiStripApvGain::Range detGainRange = gainHandle->getRange(cluster()->geographicalId()); + + std::vector gains; + for(size_t i=0; i< width(); i++){ + gains.push_back( gainHandle->getStripGain( firstStrip()+i, detGainRange) ); + } + return gains; +} + +std::vector SiStripClusterInfo:: +stripQualitiesBad() const { + std::vector isBad; + for(int i=0; i< width(); i++) { + isBad.push_back( qualityHandle->IsStripBad( cluster()->geographicalId(), + firstStrip()+i) ); + } + return isBad; +} + +float SiStripClusterInfo:: +calculate_noise(const std::vector& noise) const { + float noiseSumInQuadrature = 0; + int numberStripsOverThreshold = 0; + for(int i=0;i stripBad = stripQualitiesBad(); + return + IsApvBad() || + IsFiberBad() || + IsModuleBad() || + accumulate(stripBad.begin(), stripBad.end(), + false, + std::logical_or()); +} + +bool SiStripClusterInfo:: +IsApvBad() const { + return + qualityHandle->IsApvBad( cluster()->geographicalId(), firstStrip()/128 ) || + qualityHandle->IsApvBad( cluster()->geographicalId(), (firstStrip()+width())/128 ) ; +} + +bool SiStripClusterInfo:: +IsFiberBad() const { + return + qualityHandle->IsFiberBad( cluster()->geographicalId(), firstStrip()/256 ) || + qualityHandle->IsFiberBad( cluster()->geographicalId(), (firstStrip()+width())/256 ) ; +} + +bool SiStripClusterInfo:: +IsModuleBad() const { + return qualityHandle->IsModuleBad( cluster()->geographicalId() ); +} + +bool SiStripClusterInfo:: +IsModuleUsable() const { + return qualityHandle->IsModuleUsable( cluster()->geographicalId() ); +} + +std::vector SiStripClusterInfo:: +reclusterize(const edm::ParameterSet& conf) const { + + std::vector clusters; + + std::vector charges = stripCharges(); + std::vector gains = stripGains(); + for(unsigned i=0; i < charges.size(); i++) + charges[i] = (charges[i] < 254) + ? static_cast(charges[i] * gains[i]) + : charges[i]; + + std::auto_ptr + algorithm = StripClusterizerAlgorithmFactory::create(conf); + algorithm->initialize(es); + + if( algorithm->stripByStripBegin( detId() )) { + for(unsigned i = 0; istripByStripAdd( firstStrip()+i, charges[i], clusters ); + algorithm->stripByStripEnd( clusters ); + } + + return clusters; +} +