00001 #include "CalibTracker/SiStripESProducers/interface/SiStripNoisesGenerator.h"
00002 #include <boost/cstdint.hpp>
00003 #include "CalibTracker/SiStripCommon/interface/SiStripDetInfoFileReader.h"
00004 #include "FWCore/ParameterSet/interface/FileInPath.h"
00005 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
00006 #include "DataFormats/SiStripDetId/interface/TIBDetId.h"
00007 #include "DataFormats/SiStripDetId/interface/TIDDetId.h"
00008 #include "DataFormats/SiStripDetId/interface/TOBDetId.h"
00009 #include "DataFormats/SiStripDetId/interface/TECDetId.h"
00010
00011 #include "CLHEP/Random/RandFlat.h"
00012 #include "CLHEP/Random/RandGauss.h"
00013
00014 SiStripNoisesGenerator::SiStripNoisesGenerator(const edm::ParameterSet& iConfig,const edm::ActivityRegistry& aReg):
00015 SiStripCondObjBuilderBase<SiStripNoises>::SiStripCondObjBuilderBase(iConfig),
00016 electronsPerADC_(0.),
00017 minimumPosValue_(0.),
00018 stripLengthMode_(true),
00019 printDebug_(0)
00020 {
00021 edm::LogInfo("SiStripNoisesGenerator") << "[SiStripNoisesGenerator::SiStripNoisesGenerator]";
00022 }
00023
00024 SiStripNoisesGenerator::~SiStripNoisesGenerator()
00025 {
00026 edm::LogInfo("SiStripNoisesGenerator") << "[SiStripNoisesGenerator::~SiStripNoisesGenerator]";
00027 }
00028
00029 void SiStripNoisesGenerator::createObject()
00030 {
00031 obj_ = new SiStripNoises();
00032
00033 stripLengthMode_ = _pset.getParameter<bool>("StripLengthMode");
00034
00035
00036 std::map<int, std::vector<double> > meanNoise;
00037 fillParameters(meanNoise, "MeanNoise");
00038 std::map<int, std::vector<double> > sigmaNoise;
00039 fillParameters(sigmaNoise, "SigmaNoise");
00040 minimumPosValue_ = _pset.getParameter<double>("MinPositiveNoise");
00041
00042
00043 std::map<int, std::vector<double> > noiseStripLengthLinearSlope;
00044 fillParameters(noiseStripLengthLinearSlope, "NoiseStripLengthSlope");
00045 std::map<int, std::vector<double> > noiseStripLengthLinearQuote;
00046 fillParameters(noiseStripLengthLinearQuote, "NoiseStripLengthQuote");
00047 electronsPerADC_ = _pset.getParameter<double>("electronPerAdc");
00048
00049 printDebug_ = _pset.getUntrackedParameter<uint32_t>("printDebug", 5);
00050
00051 uint32_t count = 0;
00052 edm::FileInPath fp_ = _pset.getParameter<edm::FileInPath>("file");
00053 SiStripDetInfoFileReader reader(fp_.fullPath());
00054 const std::map<uint32_t, SiStripDetInfoFileReader::DetInfo > DetInfos = reader.getAllData();
00055
00056 for(std::map<uint32_t, SiStripDetInfoFileReader::DetInfo >::const_iterator it = DetInfos.begin(); it != DetInfos.end(); ++it) {
00057
00058
00059 SiStripNoises::InputVector theSiStripVector;
00060 float noise = 0.;
00061 uint32_t detId = it->first;
00062 std::pair<int, int> sl = subDetAndLayer(detId);
00063 unsigned short nApvs = it->second.nApvs;
00064
00065
00066 if(stripLengthMode_) {
00067
00068 double linearSlope = noiseStripLengthLinearSlope[sl.first][sl.second];
00069 double linearQuote = noiseStripLengthLinearQuote[sl.first][sl.second];
00070 double stripLength = it->second.stripLength;
00071 for( unsigned short j=0; j<128*nApvs; ++j ) {
00072 noise = ( linearSlope*stripLength + linearQuote) / electronsPerADC_;
00073 if( count<printDebug_ ) printLog(detId, j, noise);
00074 obj_->setData(noise, theSiStripVector);
00075 }
00076 }
00077 else {
00078
00079 double meanN = meanNoise[sl.first][sl.second];
00080 double sigmaN = sigmaNoise[sl.first][sl.second];
00081 for( unsigned short j=0; j<128*nApvs; ++j ) {
00082 noise = CLHEP::RandGauss::shoot(meanN, sigmaN);
00083 if( noise<=minimumPosValue_ ) noise = minimumPosValue_;
00084 if( count<printDebug_ ) printLog(detId, j, noise);
00085 obj_->setData(noise, theSiStripVector);
00086 }
00087 }
00088 ++count;
00089
00090 if ( ! obj_->put(it->first,theSiStripVector) ) {
00091 edm::LogError("SiStripNoisesFakeESSource::produce ")<<" detid already exists"<<std::endl;
00092 }
00093 }
00094 }
00095
00096 std::pair<int, int> SiStripNoisesGenerator::subDetAndLayer( const uint32_t detId ) const
00097 {
00098 int layerId = 0;
00099
00100 StripSubdetector subid(detId);
00101 int subId = subid.subdetId();
00102
00103 if( subId == int(StripSubdetector::TIB)) {
00104 TIBDetId theTIBDetId(detId);
00105 layerId = theTIBDetId.layer() - 1;
00106 }
00107 else if(subId == int(StripSubdetector::TOB)) {
00108 TOBDetId theTOBDetId(detId);
00109 layerId = theTOBDetId.layer() - 1;
00110 }
00111 else if(subId == int(StripSubdetector::TID)) {
00112 TIDDetId theTIDDetId(detId);
00113 layerId = theTIDDetId.ring() - 1;
00114 }
00115 if(subId == int(StripSubdetector::TEC)) {
00116 TECDetId theTECDetId = TECDetId(detId);
00117 layerId = theTECDetId.ring() - 1;
00118 }
00119 return std::make_pair(subId, layerId);
00120 }
00121
00122 void SiStripNoisesGenerator::fillParameters(std::map<int, std::vector<double> > & mapToFill, const std::string & parameterName) const
00123 {
00124 int layersTIB = 4;
00125 int ringsTID = 3;
00126 int layersTOB = 6;
00127 int ringsTEC = 7;
00128
00129 fillSubDetParameter( mapToFill, _pset.getParameter<std::vector<double> >(parameterName+"TIB"), int(StripSubdetector::TIB), layersTIB );
00130 fillSubDetParameter( mapToFill, _pset.getParameter<std::vector<double> >(parameterName+"TID"), int(StripSubdetector::TID), ringsTID );
00131 fillSubDetParameter( mapToFill, _pset.getParameter<std::vector<double> >(parameterName+"TOB"), int(StripSubdetector::TOB), layersTOB );
00132 fillSubDetParameter( mapToFill, _pset.getParameter<std::vector<double> >(parameterName+"TEC"), int(StripSubdetector::TEC), ringsTEC );
00133 }
00134
00135 void SiStripNoisesGenerator::fillSubDetParameter(std::map<int, std::vector<double> > & mapToFill, const std::vector<double> & v, const int subDet, const unsigned short layers) const
00136 {
00137 if( v.size() == layers ) {
00138 mapToFill.insert(std::make_pair( subDet, v ));
00139 }
00140 else if( v.size() == 1 ) {
00141 std::vector<double> parV(layers, v[0]);
00142 mapToFill.insert(std::make_pair( subDet, parV ));
00143 }
00144 else {
00145 throw cms::Exception("Configuration") << "ERROR: number of parameters for subDet " << subDet << " are " << v.size() << ". They must be either 1 or " << layers << std::endl;
00146 }
00147 }