00001 #include "CalibTracker/SiStripCommon/interface/ShallowClustersProducer.h"
00002
00003 #include "FWCore/Framework/interface/Event.h"
00004 #include "FWCore/Framework/interface/ESHandle.h"
00005 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00006 #include "RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterInfo.h"
00007 #include "DataFormats/SiStripDigi/interface/SiStripProcessedRawDigi.h"
00008 #include "DataFormats/Common/interface/DetSetVectorNew.h"
00009 #include "DataFormats/SiStripDetId/interface/TECDetId.h"
00010 #include "DataFormats/SiStripDetId/interface/TIBDetId.h"
00011 #include "DataFormats/SiStripDetId/interface/TIDDetId.h"
00012 #include "DataFormats/SiStripDetId/interface/TOBDetId.h"
00013 #include "boost/foreach.hpp"
00014
00015 ShallowClustersProducer::ShallowClustersProducer(const edm::ParameterSet& iConfig)
00016 : theClustersLabel(iConfig.getParameter<edm::InputTag>("Clusters")),
00017 Prefix(iConfig.getParameter<std::string>("Prefix") )
00018 {
00019 produces <std::vector<unsigned> > ( Prefix + "number" );
00020 produces <std::vector<unsigned> > ( Prefix + "width" );
00021 produces <std::vector<float> > ( Prefix + "variance" );
00022 produces <std::vector<float> > ( Prefix + "barystrip" );
00023 produces <std::vector<float> > ( Prefix + "middlestrip" );
00024 produces <std::vector<unsigned> > ( Prefix + "charge" );
00025 produces <std::vector<float> > ( Prefix + "noise" );
00026 produces <std::vector<float> > ( Prefix + "ston" );
00027 produces <std::vector<unsigned> > ( Prefix + "seedstrip" );
00028 produces <std::vector<unsigned> > ( Prefix + "seedindex" );
00029 produces <std::vector<unsigned> > ( Prefix + "seedcharge" );
00030 produces <std::vector<float> > ( Prefix + "seednoise" );
00031 produces <std::vector<float> > ( Prefix + "seedgain" );
00032 produces <std::vector<unsigned> > ( Prefix + "qualityisbad" );
00033
00034 produces <std::vector<float> > ( Prefix + "rawchargeC" );
00035 produces <std::vector<float> > ( Prefix + "rawchargeL" );
00036 produces <std::vector<float> > ( Prefix + "rawchargeR" );
00037 produces <std::vector<float> > ( Prefix + "rawchargeLL" );
00038 produces <std::vector<float> > ( Prefix + "rawchargeRR" );
00039 produces <std::vector<float> > ( Prefix + "eta" );
00040 produces <std::vector<float> > ( Prefix + "foldedeta" );
00041 produces <std::vector<float> > ( Prefix + "etaX" );
00042 produces <std::vector<float> > ( Prefix + "etaasymm" );
00043 produces <std::vector<float> > ( Prefix + "outsideasymm");
00044 produces <std::vector<float> > ( Prefix + "neweta");
00045 produces <std::vector<float> > ( Prefix + "newetaerr");
00046
00047 produces <std::vector<unsigned> > ( Prefix + "detid" );
00048 produces <std::vector<int> > ( Prefix + "subdetid" );
00049 produces <std::vector<int> > ( Prefix + "module" );
00050 produces <std::vector<int> > ( Prefix + "side" );
00051 produces <std::vector<int> > ( Prefix + "layerwheel" );
00052 produces <std::vector<int> > ( Prefix + "stringringrod" );
00053 produces <std::vector<int> > ( Prefix + "petal" );
00054 produces <std::vector<int> > ( Prefix + "stereo" );
00055
00056 }
00057
00058 void ShallowClustersProducer::
00059 produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
00060 std::auto_ptr<std::vector<unsigned> > number ( new std::vector<unsigned>(7,0) );
00061 std::auto_ptr<std::vector<unsigned> > width ( new std::vector<unsigned>() );
00062 std::auto_ptr<std::vector<float> > variance ( new std::vector<float>() );
00063 std::auto_ptr<std::vector<float> > barystrip ( new std::vector<float>() );
00064 std::auto_ptr<std::vector<float> > middlestrip ( new std::vector<float>() );
00065 std::auto_ptr<std::vector<unsigned> > charge ( new std::vector<unsigned>() );
00066 std::auto_ptr<std::vector<float> > noise ( new std::vector<float>() );
00067 std::auto_ptr<std::vector<float> > ston ( new std::vector<float>() );
00068 std::auto_ptr<std::vector<unsigned> > seedstrip ( new std::vector<unsigned>() );
00069 std::auto_ptr<std::vector<unsigned> > seedindex ( new std::vector<unsigned>() );
00070 std::auto_ptr<std::vector<unsigned> > seedcharge ( new std::vector<unsigned>() );
00071 std::auto_ptr<std::vector<float> > seednoise ( new std::vector<float>() );
00072 std::auto_ptr<std::vector<float> > seedgain ( new std::vector<float>() );
00073 std::auto_ptr<std::vector<unsigned> > qualityisbad ( new std::vector<unsigned>() );
00074
00075 std::auto_ptr<std::vector<float> > rawchargeC ( new std::vector<float>() );
00076 std::auto_ptr<std::vector<float> > rawchargeL ( new std::vector<float>() );
00077 std::auto_ptr<std::vector<float> > rawchargeR ( new std::vector<float>() );
00078 std::auto_ptr<std::vector<float> > rawchargeLL ( new std::vector<float>() );
00079 std::auto_ptr<std::vector<float> > rawchargeRR ( new std::vector<float>() );
00080 std::auto_ptr<std::vector<float> > etaX ( new std::vector<float>() );
00081 std::auto_ptr<std::vector<float> > eta ( new std::vector<float>() );
00082 std::auto_ptr<std::vector<float> > foldedeta ( new std::vector<float>() );
00083 std::auto_ptr<std::vector<float> > etaasymm ( new std::vector<float>() );
00084 std::auto_ptr<std::vector<float> > outsideasymm ( new std::vector<float>() );
00085 std::auto_ptr<std::vector<float> > neweta ( new std::vector<float>() );
00086 std::auto_ptr<std::vector<float> > newetaerr ( new std::vector<float>() );
00087
00088 std::auto_ptr<std::vector<unsigned> > detid ( new std::vector<unsigned>() );
00089 std::auto_ptr<std::vector<int> > subdetid ( new std::vector<int>() );
00090 std::auto_ptr<std::vector<int> > side ( new std::vector<int>() );
00091 std::auto_ptr<std::vector<int> > module ( new std::vector<int>() );
00092 std::auto_ptr<std::vector<int> > layerwheel ( new std::vector<int>() );
00093 std::auto_ptr<std::vector<int> > stringringrod ( new std::vector<int>() );
00094 std::auto_ptr<std::vector<int> > petal ( new std::vector<int>() );
00095 std::auto_ptr<std::vector<int> > stereo ( new std::vector<int>());
00096
00097 edm::Handle<edmNew::DetSetVector<SiStripCluster> > clusters;
00098 iEvent.getByLabel(theClustersLabel, clusters);
00099
00100 edm::Handle<edm::DetSetVector<SiStripProcessedRawDigi> > rawProcessedDigis;
00101 iEvent.getByLabel("siStripProcessedRawDigis", "", rawProcessedDigis);
00102
00103 BOOST_FOREACH( const edmNew::DetSet<SiStripCluster> ds, *clusters) {
00104 const moduleVars moduleV(ds.detId());
00105 BOOST_FOREACH( const SiStripCluster cluster, ds) {
00106 const SiStripClusterInfo info(cluster, iSetup);
00107 const NearDigis digis = rawProcessedDigis.isValid() ? NearDigis(info, *rawProcessedDigis) : NearDigis(info);
00108
00109 (number->at(0))++;
00110 (number->at(moduleV.subdetid))++;
00111 width->push_back( cluster.amplitudes().size() );
00112 barystrip->push_back( cluster.barycenter() );
00113 variance->push_back( info.variance() );
00114 middlestrip->push_back( info.firstStrip() + info.width()/2.0 );
00115 charge->push_back( info.charge() );
00116 noise->push_back( info.noiseRescaledByGain() );
00117 ston->push_back( info.signalOverNoise() );
00118 seedstrip->push_back( info.maxStrip() );
00119 seedindex->push_back( info.maxIndex() );
00120 seedcharge->push_back( info.maxCharge() );
00121 seednoise->push_back( info.stripNoisesRescaledByGain().at(info.maxIndex()) );
00122 seedgain->push_back( info.stripGains().at(info.maxIndex()) );
00123 qualityisbad->push_back( info.IsAnythingBad() );
00124
00125 rawchargeC->push_back( digis.max );
00126 rawchargeL->push_back( digis.left );
00127 rawchargeR->push_back( digis.right );
00128 rawchargeLL->push_back( digis.Lleft );
00129 rawchargeRR->push_back( digis.Rright );
00130 etaX->push_back( digis.etaX() );
00131 eta->push_back( digis.eta() );
00132 etaasymm->push_back( digis.etaasymm() );
00133 outsideasymm->push_back( digis.outsideasymm() );
00134 neweta->push_back( (digis.last-digis.first)/info.charge() );
00135 newetaerr->push_back( (sqrt(digis.last+digis.first))/pow(info.charge(),1.5) );
00136
00137 detid->push_back( ds.detId() );
00138 subdetid->push_back( moduleV.subdetid );
00139 side->push_back( moduleV.side );
00140 module->push_back( moduleV.module );
00141 layerwheel->push_back( moduleV.layerwheel );
00142 stringringrod->push_back( moduleV.stringringrod );
00143 petal->push_back( moduleV.petal );
00144 stereo->push_back( moduleV.stereo );
00145 }
00146 }
00147
00148 iEvent.put( number, Prefix + "number" );
00149 iEvent.put( width, Prefix + "width" );
00150 iEvent.put( variance, Prefix + "variance" );
00151 iEvent.put( barystrip, Prefix + "barystrip" );
00152 iEvent.put( middlestrip, Prefix + "middlestrip" );
00153 iEvent.put( charge, Prefix + "charge" );
00154 iEvent.put( noise, Prefix + "noise" );
00155 iEvent.put( ston, Prefix + "ston" );
00156 iEvent.put( seedstrip, Prefix + "seedstrip" );
00157 iEvent.put( seedindex, Prefix + "seedindex" );
00158 iEvent.put( seedcharge, Prefix + "seedcharge" );
00159 iEvent.put( seednoise, Prefix + "seednoise" );
00160 iEvent.put( seedgain, Prefix + "seedgain" );
00161 iEvent.put( qualityisbad, Prefix + "qualityisbad" );
00162
00163 iEvent.put( rawchargeC, Prefix + "rawchargeC" );
00164 iEvent.put( rawchargeL, Prefix + "rawchargeL" );
00165 iEvent.put( rawchargeR, Prefix + "rawchargeR" );
00166 iEvent.put( rawchargeLL, Prefix + "rawchargeLL" );
00167 iEvent.put( rawchargeRR, Prefix + "rawchargeRR" );
00168 iEvent.put( etaX, Prefix + "etaX" );
00169 iEvent.put( eta, Prefix + "eta" );
00170 iEvent.put( foldedeta, Prefix + "foldedeta" );
00171 iEvent.put( etaasymm, Prefix + "etaasymm" );
00172 iEvent.put( outsideasymm, Prefix + "outsideasymm" );
00173 iEvent.put( neweta, Prefix + "neweta" );
00174 iEvent.put( newetaerr, Prefix + "newetaerr" );
00175
00176 iEvent.put( detid, Prefix + "detid" );
00177 iEvent.put( subdetid, Prefix + "subdetid" );
00178 iEvent.put( module, Prefix + "module" );
00179 iEvent.put( side, Prefix + "side" );
00180 iEvent.put( layerwheel, Prefix + "layerwheel" );
00181 iEvent.put( stringringrod, Prefix + "stringringrod" );
00182 iEvent.put( petal, Prefix + "petal" );
00183 iEvent.put( stereo, Prefix + "stereo" );
00184
00185 }
00186
00187 ShallowClustersProducer::NearDigis::
00188 NearDigis(const SiStripClusterInfo& info) {
00189 max = info.maxCharge();
00190 left = info.maxIndex() > uint16_t(0) ? info.stripCharges().at(info.maxIndex()-1) : 0 ;
00191 Lleft = info.maxIndex() > uint16_t(1) ? info.stripCharges().at(info.maxIndex()-2) : 0 ;
00192 right= unsigned(info.maxIndex()+1) < info.stripCharges().size() ? info.stripCharges().at(info.maxIndex()+1) : 0 ;
00193 Rright= unsigned(info.maxIndex()+2) < info.stripCharges().size() ? info.stripCharges().at(info.maxIndex()+2) : 0 ;
00194 first = info.stripCharges().at(0);
00195 last = info.stripCharges().at(info.width()-1);
00196 }
00197
00198 ShallowClustersProducer::NearDigis::
00199 NearDigis(const SiStripClusterInfo& info, const edm::DetSetVector<SiStripProcessedRawDigi>& rawProcessedDigis) {
00200 edm::DetSetVector<SiStripProcessedRawDigi>::const_iterator digiframe = rawProcessedDigis.find(info.detId());
00201 if( digiframe != rawProcessedDigis.end()) {
00202 max = digiframe->data.at(info.maxStrip()).adc() ;
00203 left = info.maxStrip() > uint16_t(0) ? digiframe->data.at(info.maxStrip()-1).adc() : 0 ;
00204 Lleft = info.maxStrip() > uint16_t(1) ? digiframe->data.at(info.maxStrip()-2).adc() : 0 ;
00205 right = unsigned(info.maxStrip()+1) < digiframe->data.size() ? digiframe->data.at(info.maxStrip()+1).adc() : 0 ;
00206 Rright = unsigned(info.maxStrip()+2) < digiframe->data.size() ? digiframe->data.at(info.maxStrip()+2).adc() : 0 ;
00207 first = digiframe->data.at(info.firstStrip()).adc();
00208 last = digiframe->data.at(info.firstStrip()+info.width() - 1).adc();
00209 } else {
00210 *this = NearDigis(info);
00211 }
00212 }
00213
00214 ShallowClustersProducer::moduleVars::
00215 moduleVars(uint32_t detid) {
00216 SiStripDetId subdet(detid);
00217 subdetid = subdet.subDetector();
00218 if( SiStripDetId::TIB == subdetid ) {
00219 TIBDetId tib(detid);
00220 module = tib.module();
00221 side = (tib.isZMinusSide())?-1:1;
00222 layerwheel = tib.layer();
00223 stringringrod = tib.stringNumber();
00224 stereo = tib.isStereo() ? 1 : 0;
00225 } else
00226 if( SiStripDetId::TID == subdetid ) {
00227 TIDDetId tid(detid);
00228 module = tid.moduleNumber();
00229 side = (tid.isZMinusSide())?-1:1;
00230 layerwheel = tid.wheel();
00231 stringringrod = tid.ringNumber();
00232 stereo = tid.isStereo() ? 1 : 0;
00233 } else
00234 if( SiStripDetId::TOB == subdetid ) {
00235 TOBDetId tob(detid);
00236 module = tob.module();
00237 side = (tob.isZMinusSide())?-1:1;
00238 layerwheel = tob.layer();
00239 stringringrod = tob.rodNumber();
00240 stereo = tob.isStereo() ? 1 : 0;
00241 } else
00242 if( SiStripDetId::TEC == subdetid ) {
00243 TECDetId tec(detid);
00244 module = tec.module();
00245 side = (tec.isZMinusSide())?-1:1;
00246 layerwheel = tec.wheel();
00247 stringringrod = tec.ringNumber();
00248 petal = tec.petalNumber();
00249 stereo = tec.isStereo() ? 1 : 0;
00250 } else {
00251 module = 0;
00252 side = 0;
00253 layerwheel=-1;
00254 stringringrod = -1;
00255 petal=-1;
00256 }
00257 }
1.7.3