Go to the documentation of this file.00001 #include "RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterInfo.h"
00002
00003 #include "CondFormats/DataRecord/interface/SiStripNoisesRcd.h"
00004 #include "CalibTracker/Records/interface/SiStripGainRcd.h"
00005 #include "CalibTracker/Records/interface/SiStripQualityRcd.h"
00006 #include "FWCore/Framework/interface/EventSetup.h"
00007 #include "RecoLocalTracker/SiStripClusterizer/interface/StripClusterizerAlgorithmFactory.h"
00008 #include "RecoLocalTracker/SiStripClusterizer/interface/StripClusterizerAlgorithm.h"
00009 #include <cmath>
00010
00011 SiStripClusterInfo::SiStripClusterInfo(const SiStripCluster& cluster,
00012 const edm::EventSetup& setup,
00013 std::string quality)
00014 : cluster_ptr(&cluster),
00015 es(setup),
00016 qualityLabel(quality) {
00017 es.get<SiStripNoisesRcd>().get(noiseHandle);
00018 es.get<SiStripGainRcd>().get(gainHandle);
00019 es.get<SiStripQualityRcd>().get(qualityLabel,qualityHandle);
00020 }
00021
00022 std::pair<uint16_t,uint16_t > SiStripClusterInfo::
00023 chargeLR() const {
00024 std::vector<uint8_t>::const_iterator
00025 begin( stripCharges().begin() ),
00026 end( stripCharges().end() ),
00027 max; max = max_element(begin,end);
00028 return std::make_pair( accumulate(begin, max, uint16_t(0) ),
00029 accumulate(max+1, end, uint16_t(0) ) );
00030 }
00031
00032
00033 float SiStripClusterInfo::
00034 variance() const {
00035 float q(0), x1(0), x2(0);
00036 for(std::vector<uint8_t>::const_iterator
00037 begin(stripCharges().begin()), end(stripCharges().end()), it(begin);
00038 it!=end; ++it) {
00039 unsigned i = it-begin;
00040 q += (*it);
00041 x1 += (*it) * (i+0.5);
00042 x2 += (*it) * (i*i+i+1./3);
00043 }
00044 return (x2 - x1*x1/q ) / q;
00045 }
00046
00047 std::vector<float> SiStripClusterInfo::
00048 stripNoisesRescaledByGain() const {
00049 std::vector<float> noises = stripNoises();
00050 std::vector<float> gains = stripGains();
00051 transform(noises.begin(), noises.end(), gains.begin(),
00052 noises.begin(),
00053 std::divides<double>());
00054 return noises;
00055 }
00056
00057 std::vector<float> SiStripClusterInfo::
00058 stripNoises() const {
00059 SiStripNoises::Range detNoiseRange = noiseHandle->getRange(cluster()->geographicalId());
00060
00061 std::vector<float> noises;
00062 for(size_t i=0; i < width(); i++){
00063 noises.push_back( noiseHandle->getNoise( firstStrip()+i, detNoiseRange) );
00064 }
00065 return noises;
00066 }
00067
00068 std::vector<float> SiStripClusterInfo::
00069 stripGains() const {
00070 SiStripApvGain::Range detGainRange = gainHandle->getRange(cluster()->geographicalId());
00071
00072 std::vector<float> gains;
00073 for(size_t i=0; i< width(); i++){
00074 gains.push_back( gainHandle->getStripGain( firstStrip()+i, detGainRange) );
00075 }
00076 return gains;
00077 }
00078
00079 std::vector<bool> SiStripClusterInfo::
00080 stripQualitiesBad() const {
00081 std::vector<bool> isBad;
00082 for(int i=0; i< width(); i++) {
00083 isBad.push_back( qualityHandle->IsStripBad( cluster()->geographicalId(),
00084 firstStrip()+i) );
00085 }
00086 return isBad;
00087 }
00088
00089 float SiStripClusterInfo::
00090 calculate_noise(const std::vector<float>& noise) const {
00091 float noiseSumInQuadrature = 0;
00092 int numberStripsOverThreshold = 0;
00093 for(int i=0;i<width();i++) {
00094 if(stripCharges().at(i)!=0) {
00095 noiseSumInQuadrature += noise.at(i) * noise.at(i);
00096 numberStripsOverThreshold++;
00097 }
00098 }
00099 return std::sqrt( noiseSumInQuadrature / numberStripsOverThreshold );
00100 }
00101
00102
00103 bool SiStripClusterInfo::
00104 IsAnythingBad() const {
00105 std::vector<bool> stripBad = stripQualitiesBad();
00106 return
00107 IsApvBad() ||
00108 IsFiberBad() ||
00109 IsModuleBad() ||
00110 accumulate(stripBad.begin(), stripBad.end(),
00111 false,
00112 std::logical_or<bool>());
00113 }
00114
00115 bool SiStripClusterInfo::
00116 IsApvBad() const {
00117 return
00118 qualityHandle->IsApvBad( cluster()->geographicalId(), firstStrip()/128 ) ||
00119 qualityHandle->IsApvBad( cluster()->geographicalId(), (firstStrip()+width())/128 ) ;
00120 }
00121
00122 bool SiStripClusterInfo::
00123 IsFiberBad() const {
00124 return
00125 qualityHandle->IsFiberBad( cluster()->geographicalId(), firstStrip()/256 ) ||
00126 qualityHandle->IsFiberBad( cluster()->geographicalId(), (firstStrip()+width())/256 ) ;
00127 }
00128
00129 bool SiStripClusterInfo::
00130 IsModuleBad() const {
00131 return qualityHandle->IsModuleBad( cluster()->geographicalId() );
00132 }
00133
00134 bool SiStripClusterInfo::
00135 IsModuleUsable() const {
00136 return qualityHandle->IsModuleUsable( cluster()->geographicalId() );
00137 }
00138
00139 std::vector<SiStripCluster> SiStripClusterInfo::
00140 reclusterize(const edm::ParameterSet& conf) const {
00141
00142 std::vector<SiStripCluster> clusters;
00143
00144 std::vector<uint8_t> charges = stripCharges();
00145 std::vector<float> gains = stripGains();
00146 for(unsigned i=0; i < charges.size(); i++)
00147 charges[i] = (charges[i] < 254)
00148 ? static_cast<uint8_t>(charges[i] * gains[i])
00149 : charges[i];
00150
00151 std::auto_ptr<StripClusterizerAlgorithm>
00152 algorithm = StripClusterizerAlgorithmFactory::create(conf);
00153 algorithm->initialize(es);
00154
00155 if( algorithm->stripByStripBegin( detId() )) {
00156 for(unsigned i = 0; i<width(); i++)
00157 algorithm->stripByStripAdd( firstStrip()+i, charges[i], clusters );
00158 algorithm->stripByStripEnd( clusters );
00159 }
00160
00161 return clusters;
00162 }
00163