File indexing completed on 2023-10-25 09:35:37
0001 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0002 #include "FWCore/ParameterSet/interface/FileInPath.h"
0003 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0004 #include "FWCore/Utilities/interface/Exception.h"
0005
0006 #include "CommonTools/ConditionDBWriter/interface/ConditionDBWriter.h"
0007 #include "CondFormats/SiStripObjects/interface/SiStripBadStrip.h"
0008 #include "CalibTracker/SiStripCommon/interface/SiStripDetInfoFileReader.h"
0009
0010 #include <vector>
0011 #include <iostream>
0012 #include <fstream>
0013 #include <ext/hash_map>
0014
0015 class SiStripBadModuleByHandBuilder : public ConditionDBWriter<SiStripBadStrip> {
0016 public:
0017 explicit SiStripBadModuleByHandBuilder(const edm::ParameterSet&);
0018 ~SiStripBadModuleByHandBuilder() override = default;
0019
0020 private:
0021 std::unique_ptr<SiStripBadStrip> getNewObject() override;
0022
0023 private:
0024 edm::FileInPath fp_;
0025 bool printdebug_;
0026 std::vector<uint32_t> BadModuleList_;
0027 };
0028
0029 SiStripBadModuleByHandBuilder::SiStripBadModuleByHandBuilder(const edm::ParameterSet& iConfig)
0030 : ConditionDBWriter<SiStripBadStrip>(iConfig) {
0031 fp_ = iConfig.getUntrackedParameter<edm::FileInPath>("file", edm::FileInPath(SiStripDetInfoFileReader::kDefaultFile));
0032 BadModuleList_ = iConfig.getUntrackedParameter<std::vector<uint32_t> >("BadModuleList");
0033 printdebug_ = iConfig.getUntrackedParameter<bool>("printDebug", false);
0034 }
0035
0036 std::unique_ptr<SiStripBadStrip> SiStripBadModuleByHandBuilder::getNewObject() {
0037 const auto detInfo = SiStripDetInfoFileReader::read(fp_.fullPath());
0038
0039 auto obj = std::make_unique<SiStripBadStrip>();
0040
0041 unsigned int firstBadStrip = 0;
0042 unsigned short NconsecutiveBadStrips;
0043 unsigned int theBadStripRange;
0044
0045 for (std::vector<uint32_t>::const_iterator it = BadModuleList_.begin(); it != BadModuleList_.end(); ++it) {
0046 std::vector<unsigned int> theSiStripVector;
0047
0048 NconsecutiveBadStrips = detInfo.getNumberOfApvsAndStripLength(*it).first * 128;
0049 theBadStripRange = obj->encode(firstBadStrip, NconsecutiveBadStrips);
0050 if (printdebug_)
0051 edm::LogInfo("SiStripBadModuleByHandBuilder")
0052 << " BadModule " << *it << " \t"
0053 << " firstBadStrip " << firstBadStrip << "\t "
0054 << " NconsecutiveBadStrips " << NconsecutiveBadStrips << "\t "
0055 << " packed integer " << std::hex << theBadStripRange << std::dec << std::endl;
0056
0057 theSiStripVector.push_back(theBadStripRange);
0058 SiStripBadStrip::Range range(theSiStripVector.begin(), theSiStripVector.end());
0059 if (!obj->put(*it, range))
0060 edm::LogError("SiStripBadModuleByHandBuilder")
0061 << "[SiStripBadModuleByHandBuilder::analyze] detid already exists" << std::endl;
0062 }
0063 return obj;
0064 }
0065
0066 #include "FWCore/Framework/interface/MakerMacros.h"
0067 DEFINE_FWK_MODULE(SiStripBadModuleByHandBuilder);