File indexing completed on 2024-09-11 04:32:41
0001 #include "DQM/SiStripMonitorSummary/interface/SiStripLorentzAngleDQM.h"
0002 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
0003 #include "DataFormats/TrackerCommon/interface/SiStripSubStructure.h"
0004 #include "TCanvas.h"
0005
0006 SiStripLorentzAngleDQM::SiStripLorentzAngleDQM(edm::ESGetToken<SiStripLorentzAngle, SiStripLorentzAngleRcd> token,
0007 edm::RunNumber_t iRun,
0008 edm::ParameterSet const &hPSet,
0009 edm::ParameterSet const &fPSet,
0010 const TrackerTopology *tTopo,
0011 const TkDetMap *tkDetMap)
0012 : SiStripBaseCondObjDQMGet<SiStripLorentzAngle, SiStripLorentzAngleRcd>{token, iRun, hPSet, fPSet, tTopo} {
0013 if (HistoMaps_On_) {
0014 Tk_HM_ = std::make_unique<TkHistoMap>(tkDetMap, "SiStrip/Histo_Map", "LA_TkMap", 0.);
0015 }
0016 }
0017
0018 SiStripLorentzAngleDQM::~SiStripLorentzAngleDQM() {}
0019
0020 void SiStripLorentzAngleDQM::getActiveDetIds(const edm::EventSetup &eSetup) {
0021 getConditionObject(eSetup);
0022
0023 std::map<uint32_t, float>::const_iterator LAMapIter_;
0024 std::map<uint32_t, float> LAMap_ = condObj_->getLorentzAngles();
0025
0026 for (LAMapIter_ = LAMap_.begin(); LAMapIter_ != LAMap_.end(); LAMapIter_++) {
0027 activeDetIds.push_back((*LAMapIter_).first);
0028 }
0029 }
0030
0031 void SiStripLorentzAngleDQM::fillSummaryMEs(const std::vector<uint32_t> &selectedDetIds) {
0032
0033
0034
0035
0036 bool fillNext = true;
0037 for (unsigned int i = 0; i < selectedDetIds.size(); i++) {
0038 int subDetId_ = DetId(selectedDetIds[i]).subdetId();
0039 if (subDetId_ < 3 || subDetId_ > 6) {
0040 edm::LogError("SiStripLorentzAngle")
0041 << "[SiStripLorentzAngle::fillSummaryMEs] WRONG INPUT : no such "
0042 "subdetector type : "
0043 << subDetId_ << " and detId " << selectedDetIds[i] << " therefore no filling!" << std::endl;
0044 } else if (SummaryOnLayerLevel_On_) {
0045 if (fillNext) {
0046 fillMEsForLayer( selectedDetIds[i]);
0047 }
0048 if (getLayerNameAndId(selectedDetIds[i + 1]) == getLayerNameAndId(selectedDetIds[i])) {
0049 fillNext = false;
0050 } else {
0051 fillNext = true;
0052 }
0053 } else if (SummaryOnStringLevel_On_) {
0054 if (fillNext) {
0055 fillMEsForLayer( selectedDetIds[i]);
0056 }
0057 if (getStringNameAndId(selectedDetIds[i + 1]) == getStringNameAndId(selectedDetIds[i])) {
0058 fillNext = false;
0059 } else {
0060 fillNext = true;
0061 }
0062 }
0063 }
0064
0065 for (std::map<uint32_t, ModMEs>::iterator iter = SummaryMEsMap_.begin(); iter != SummaryMEsMap_.end(); iter++) {
0066 ModMEs selME;
0067 selME = iter->second;
0068
0069 if (SummaryOnStringLevel_On_) {
0070 if (fPSet_.getParameter<bool>("OutputSummaryProfileAtLayerLevelAsImage")) {
0071 TCanvas c1("c1");
0072 selME.SummaryOfProfileDistr->getTProfile()->Draw();
0073 std::string name(selME.SummaryOfProfileDistr->getTProfile()->GetTitle());
0074 name += ".png";
0075 c1.Print(name.c_str());
0076 }
0077
0078 if (fPSet_.getParameter<bool>("OutputCumulativeSummaryAtLayerLevelAsImage")) {
0079 TCanvas c2("c2");
0080 selME.SummaryOfCumulDistr->getTH1()->Draw();
0081 std::string name2(selME.SummaryOfCumulDistr->getTitle());
0082 name2 += ".png";
0083 c2.Print(name2.c_str());
0084 }
0085
0086 } else {
0087 if (hPSet_.getParameter<bool>("FillSummaryProfileAtLayerLevel") &&
0088 fPSet_.getParameter<bool>("OutputSummaryProfileAtLayerLevelAsImage")) {
0089 TCanvas c1("c1");
0090 selME.SummaryOfProfileDistr->getTProfile()->Draw();
0091 std::string name(selME.SummaryOfProfileDistr->getTProfile()->GetTitle());
0092 name += ".png";
0093 c1.Print(name.c_str());
0094 }
0095
0096 if (hPSet_.getParameter<bool>("FillCumulativeSummaryAtLayerLevel") &&
0097 fPSet_.getParameter<bool>("OutputCumulativeSummaryAtLayerLevelAsImage")) {
0098 TCanvas c1("c1");
0099 selME.SummaryOfCumulDistr->getTH1()->Draw();
0100 std::string name(selME.SummaryOfCumulDistr->getTitle());
0101 name += ".png";
0102 c1.Print(name.c_str());
0103 }
0104 }
0105 }
0106 }
0107
0108 void SiStripLorentzAngleDQM::fillMEsForLayer(
0109 uint32_t selDetId_) {
0110 SiStripHistoId hidmanager;
0111
0112 std::string hSummaryOfProfile_description;
0113 hSummaryOfProfile_description = hPSet_.getParameter<std::string>("SummaryOfProfile_description");
0114
0115 std::string hSummary_name;
0116
0117 int subDetId_ = DetId(selDetId_).subdetId();
0118
0119 if (subDetId_ < 3 || subDetId_ > 6) {
0120 edm::LogError("SiStripLorentzAngleDQM") << "[SiStripLorentzAngleDQM::fillMEsForLayer] WRONG INPUT : no such "
0121 "subdetector type : "
0122 << subDetId_ << " no folder set!" << std::endl;
0123 return;
0124 }
0125
0126 std::vector<uint32_t> sameLayerDetIds_;
0127 sameLayerDetIds_.clear();
0128
0129 if (SummaryOnStringLevel_On_) {
0130
0131 hSummary_name =
0132 hidmanager.createHistoLayer(hSummaryOfProfile_description, "layer", getStringNameAndId(selDetId_).first, "");
0133 std::map<uint32_t, ModMEs>::iterator selMEsMapIter_ = SummaryMEsMap_.find(getStringNameAndId(selDetId_).second);
0134
0135 ModMEs selME_;
0136 if (selMEsMapIter_ != SummaryMEsMap_.end())
0137 selME_ = selMEsMapIter_->second;
0138
0139 getSummaryMEs(selME_, selDetId_);
0140
0141
0142 sameLayerDetIds_.clear();
0143
0144 switch (DetId(selDetId_).subdetId()) {
0145 case StripSubdetector::TIB:
0146 if (tTopo_->tibIsInternalString(selDetId_)) {
0147 SiStripSubStructure::getTIBDetectors(
0148 activeDetIds, sameLayerDetIds_, tTopo_, tTopo_->tibLayer(selDetId_), 0, 1, tTopo_->tibString(selDetId_));
0149 }
0150 if (tTopo_->tibIsExternalString(selDetId_)) {
0151 SiStripSubStructure::getTIBDetectors(
0152 activeDetIds, sameLayerDetIds_, tTopo_, tTopo_->tibLayer(selDetId_), 0, 2, tTopo_->tibString(selDetId_));
0153 }
0154 break;
0155 case StripSubdetector::TID:
0156 SiStripSubStructure::getTIDDetectors(activeDetIds, sameLayerDetIds_, tTopo_, 0, 0, 0, 0);
0157 break;
0158 case StripSubdetector::TOB:
0159 SiStripSubStructure::getTOBDetectors(
0160 activeDetIds, sameLayerDetIds_, tTopo_, tTopo_->tobLayer(selDetId_), 0, tTopo_->tobRod(selDetId_));
0161 break;
0162 case StripSubdetector::TEC:
0163 SiStripSubStructure::getTECDetectors(activeDetIds, sameLayerDetIds_, tTopo_, 0, 0, 0, 0, 0, 0);
0164 break;
0165 }
0166
0167
0168
0169 for (unsigned int i = 0; i < sameLayerDetIds_.size(); i++) {
0170 selME_.SummaryOfProfileDistr->Fill(i + 1, condObj_->getLorentzAngle(sameLayerDetIds_[i]));
0171
0172
0173 if (HistoMaps_On_)
0174 Tk_HM_->fill(sameLayerDetIds_[i], condObj_->getLorentzAngle(sameLayerDetIds_[i]));
0175
0176 std::cout << sameLayerDetIds_[i] << "\t" << condObj_->getLorentzAngle(sameLayerDetIds_[i]) << std::endl;
0177
0178 if (fPSet_.getParameter<bool>("TkMap_On") || hPSet_.getParameter<bool>("TkMap_On")) {
0179 fillTkMap(sameLayerDetIds_[i], condObj_->getLorentzAngle(sameLayerDetIds_[i]));
0180 }
0181 }
0182
0183 std::string hSummaryOfCumul_description;
0184 hSummaryOfCumul_description = hPSet_.getParameter<std::string>("SummaryOfCumul_description");
0185
0186 std::string hSummaryOfCumul_name;
0187
0188 if (subDetId_ < 3 || subDetId_ > 6) {
0189 edm::LogError("SiStripLorentzAngleDQM") << "[SiStripLorentzAngleDQM::fillMEsForLayer] WRONG INPUT : no such "
0190 "subdetector type : "
0191 << subDetId_ << " no folder set!" << std::endl;
0192 return;
0193 }
0194
0195 hSummaryOfCumul_name =
0196 hidmanager.createHistoLayer(hSummaryOfCumul_description, "layer", getStringNameAndId(selDetId_).first, "");
0197
0198 for (unsigned int i = 0; i < sameLayerDetIds_.size(); i++) {
0199 selME_.SummaryOfCumulDistr->Fill(condObj_->getLorentzAngle(sameLayerDetIds_[i]));
0200 }
0201 }
0202
0203 else {
0204
0205 std::map<uint32_t, ModMEs>::iterator selMEsMapIter_ = SummaryMEsMap_.find(getLayerNameAndId(selDetId_).second);
0206
0207 ModMEs selME_;
0208 if (selMEsMapIter_ != SummaryMEsMap_.end())
0209 selME_ = selMEsMapIter_->second;
0210
0211 getSummaryMEs(selME_, selDetId_);
0212
0213 if (hPSet_.getParameter<bool>("FillSummaryProfileAtLayerLevel")) {
0214 hSummary_name =
0215 hidmanager.createHistoLayer(hSummaryOfProfile_description, "layer", getLayerNameAndId(selDetId_).first, "");
0216
0217
0218 sameLayerDetIds_.clear();
0219
0220 sameLayerDetIds_ = GetSameLayerDetId(activeDetIds, selDetId_);
0221
0222 for (unsigned int i = 0; i < sameLayerDetIds_.size(); i++) {
0223 selME_.SummaryOfProfileDistr->Fill(i + 1, condObj_->getLorentzAngle(sameLayerDetIds_[i]));
0224
0225
0226 if (HistoMaps_On_)
0227 Tk_HM_->fill(sameLayerDetIds_[i], condObj_->getLorentzAngle(sameLayerDetIds_[i]));
0228
0229 if (fPSet_.getParameter<bool>("TkMap_On") || hPSet_.getParameter<bool>("TkMap_On")) {
0230 fillTkMap(sameLayerDetIds_[i], condObj_->getLorentzAngle(sameLayerDetIds_[i]));
0231 }
0232 }
0233 }
0234
0235 if (hPSet_.getParameter<bool>("FillCumulativeSummaryAtLayerLevel")) {
0236 std::string hSummaryOfCumul_description;
0237 hSummaryOfCumul_description = hPSet_.getParameter<std::string>("SummaryOfCumul_description");
0238
0239 std::string hSummaryOfCumul_name;
0240
0241 if (subDetId_ < 3 || subDetId_ > 6) {
0242 edm::LogError("SiStripLorentzAngleDQM") << "[SiStripLorentzAngleDQM::fillMEsForLayer] WRONG INPUT : no "
0243 "such subdetector type : "
0244 << subDetId_ << " no folder set!" << std::endl;
0245 return;
0246 }
0247
0248 hSummaryOfCumul_name =
0249 hidmanager.createHistoLayer(hSummaryOfCumul_description, "layer", getLayerNameAndId(selDetId_).first, "");
0250
0251 for (unsigned int i = 0; i < sameLayerDetIds_.size(); i++) {
0252 selME_.SummaryOfCumulDistr->Fill(condObj_->getLorentzAngle(sameLayerDetIds_[i]));
0253 }
0254 }
0255 }
0256 }