File indexing completed on 2023-03-17 10:47:52
0001
0002
0003
0004
0005
0006
0007
0008 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0009 #include "FWCore/Framework/interface/ESHandle.h"
0010 #include "FWCore/Framework/interface/Event.h"
0011 #include "FWCore/Framework/interface/MakerMacros.h"
0012 #include <cmath>
0013 #include <fstream>
0014 #include <iostream>
0015 #include <map>
0016 #include <sstream>
0017 #include <stdexcept>
0018 #include <string>
0019
0020 #include "FWCore/Framework/interface/EventSetup.h"
0021 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0022
0023 #include "CondTools/DT/test/validate/DTPerformanceValidateDBRead.h"
0024 #include "CondFormats/DTObjects/interface/DTPerformance.h"
0025 #include "CondFormats/DataRecord/interface/DTPerformanceRcd.h"
0026
0027 DTPerformanceValidateDBRead::DTPerformanceValidateDBRead(edm::ParameterSet const& p)
0028 : dataFileName(p.getParameter<std::string>("chkFile")),
0029 elogFileName(p.getParameter<std::string>("logFile")),
0030 dtperfToken_(esConsumes()) {}
0031
0032 DTPerformanceValidateDBRead::DTPerformanceValidateDBRead(int i) : dtperfToken_(esConsumes()) {}
0033
0034 void DTPerformanceValidateDBRead::analyze(const edm::Event& e, const edm::EventSetup& context) {
0035 using namespace edm::eventsetup;
0036
0037 std::cout << " I AM IN RUN NUMBER " << e.id().run() << std::endl;
0038 std::cout << " ---EVENT NUMBER " << e.id().event() << std::endl;
0039 std::stringstream run_fn;
0040 run_fn << "run" << e.id().run() << dataFileName;
0041 std::ifstream chkFile(run_fn.str().c_str());
0042 std::ofstream logFile(elogFileName.c_str(), std::ios_base::app);
0043 auto mP = context.getHandle(dtperfToken_);
0044 std::cout << mP->version() << std::endl;
0045 std::cout << std::distance(mP->begin(), mP->end()) << " data in the container" << std::endl;
0046 int status;
0047 int whe;
0048 int sta;
0049 int sec;
0050 int qua;
0051 float meanT0;
0052 float meanTtrig;
0053 float meanMtime;
0054 float meanNoise;
0055 float meanAfterPulse;
0056 float meanResolution;
0057 float meanEfficiency;
0058 float ckmeanT0;
0059 float ckmeanTtrig;
0060 float ckmeanMtime;
0061 float ckmeanNoise;
0062 float ckmeanAfterPulse;
0063 float ckmeanResolution;
0064 float ckmeanEfficiency;
0065
0066 DTPerformance::const_iterator iter = mP->begin();
0067 DTPerformance::const_iterator iend = mP->end();
0068 while (iter != iend) {
0069 const DTPerformanceId& mpId = iter->first;
0070 const DTPerformanceData& mpData = iter->second;
0071 status = mP->get(mpId.wheelId,
0072 mpId.stationId,
0073 mpId.sectorId,
0074 mpId.slId,
0075 meanT0,
0076 meanTtrig,
0077 meanMtime,
0078 meanNoise,
0079 meanAfterPulse,
0080 meanResolution,
0081 meanEfficiency,
0082 DTTimeUnits::counts);
0083 if (status)
0084 logFile << "ERROR while getting sl performance " << mpId.wheelId << " " << mpId.stationId << " " << mpId.sectorId
0085 << " " << mpId.slId << " , status = " << status << std::endl;
0086 if ((fabs(mpData.meanT0 - meanT0) > 0.0001) || (fabs(mpData.meanTtrig - meanTtrig) > 0.1) ||
0087 (fabs(mpData.meanMtime - meanMtime) > 0.0001) || (fabs(mpData.meanNoise - meanNoise) > 0.0001) ||
0088 (fabs(mpData.meanAfterPulse - meanAfterPulse) > 0.0001) ||
0089 (fabs(mpData.meanResolution - meanResolution) > 0.0001) ||
0090 (fabs(mpData.meanEfficiency - meanEfficiency) > 0.0001))
0091 logFile << "MISMATCH WHEN READING sl performance " << mpId.wheelId << " " << mpId.stationId << " "
0092 << mpId.sectorId << " " << mpId.slId << " : " << mpData.meanT0 << " " << mpData.meanTtrig << " "
0093 << mpData.meanMtime << " " << mpData.meanNoise << " " << mpData.meanAfterPulse << " "
0094 << mpData.meanResolution << " " << mpData.meanEfficiency << " -> " << meanT0 << " " << meanTtrig << " "
0095 << meanMtime << " " << meanNoise << " " << meanAfterPulse << " " << meanResolution << " "
0096 << meanEfficiency << std::endl;
0097 iter++;
0098 }
0099
0100 while (chkFile >> whe >> sta >> sec >> qua >> ckmeanT0 >> ckmeanTtrig >> ckmeanMtime >> ckmeanNoise >>
0101 ckmeanAfterPulse >> ckmeanResolution >> ckmeanEfficiency) {
0102 status = mP->get(whe,
0103 sta,
0104 sec,
0105 qua,
0106 meanT0,
0107 meanTtrig,
0108 meanMtime,
0109 meanNoise,
0110 meanAfterPulse,
0111 meanResolution,
0112 meanEfficiency,
0113 DTTimeUnits::counts);
0114 if ((fabs(ckmeanT0 - meanT0) > 0.0001) || (fabs(ckmeanTtrig - meanTtrig) > 0.1) ||
0115 (fabs(ckmeanMtime - meanMtime) > 0.01) || (fabs(ckmeanNoise - meanNoise) > 0.0001) ||
0116 (fabs(ckmeanAfterPulse - meanAfterPulse) > 0.0001) || (fabs(ckmeanResolution - meanResolution) > 0.01) ||
0117 (fabs(ckmeanEfficiency - meanEfficiency) > 0.0001))
0118 logFile << "MISMATCH IN WRITING AND READING sl performance " << whe << " " << sta << " " << sec << " " << qua
0119 << " : " << ckmeanT0 << " " << ckmeanTtrig << " " << ckmeanMtime << " " << ckmeanNoise << " "
0120 << ckmeanAfterPulse << " " << ckmeanResolution << " " << ckmeanEfficiency << " -> " << meanT0 << " "
0121 << meanTtrig << " " << meanMtime << " " << meanNoise << " " << meanAfterPulse << " " << meanResolution
0122 << " " << meanEfficiency << std::endl;
0123 }
0124 }
0125
0126 void DTPerformanceValidateDBRead::endJob() {
0127 std::ifstream logFile(elogFileName.c_str());
0128 char* line = new char[1000];
0129 int errors = 0;
0130 std::cout << "Mean performance validation result:" << std::endl;
0131 while (logFile.getline(line, 1000)) {
0132 std::cout << line << std::endl;
0133 errors++;
0134 }
0135 if (!errors) {
0136 std::cout << " ********************************* " << std::endl;
0137 std::cout << " *** *** " << std::endl;
0138 std::cout << " *** NO ERRORS FOUND *** " << std::endl;
0139 std::cout << " *** *** " << std::endl;
0140 std::cout << " ********************************* " << std::endl;
0141 }
0142 return;
0143 }
0144
0145 DEFINE_FWK_MODULE(DTPerformanceValidateDBRead);