/
githubmirror
/
cmssw
Обзор
Документация
Войти
/
githubmirror
/
cmssw
Код
Запросы
0
Пакеты
0
Релизы
0
Аналитика
Безопасность
master
Alignment/MuonAlignment/plugins/MuonGeometrySanityCheck.cc
730 строк
25 KB
Shahzad Malik Muzaffar
[ALCA-DB] Drop Geometry/CommonDetUnit package
28 май 2026, 12:57
28 май 2026, 12:57
993d9d5
Код
Авторство
О чём код?
// -*- C++ -*- // // Package: MuonGeometrySanityCheck // Class: MuonGeometrySanityCheck // /**\class MuonGeometrySanityCheck MuonGeometrySanityCheck.cc Alignment/MuonAlignment/plugins/MuonGeometrySanityCheck.cc Description: <one line class summary> Implementation: <Notes on implementation> */ // // Original Author: Jim Pivarski // Created: Sat Jul 3 13:33:13 CDT 2010 // $Id: MuonGeometrySanityCheck.cc,v 1.3 2011/10/12 22:13:00 khotilov Exp $ // // // system include files #include "FWCore/Framework/interface/Frameworkfwd.h" #include "FWCore/Framework/interface/one/EDAnalyzer.h" #include "FWCore/Framework/interface/EventSetup.h" #include "FWCore/Framework/interface/ESHandle.h" #include "FWCore/Framework/interface/Event.h" #include "FWCore/Framework/interface/MakerMacros.h" #include "FWCore/ParameterSet/interface/ParameterSet.h" #include "FWCore/MessageLogger/interface/MessageLogger.h" #include "Geometry/DTGeometry/interface/DTGeometry.h" #include "Geometry/CSCGeometry/interface/CSCGeometry.h" #include "Geometry/Records/interface/MuonGeometryRecord.h" #include "Geometry/CommonTopologies/interface/GeomDet.h" #include "DataFormats/GeometryVector/interface/GlobalPoint.h" #include "DataFormats/DetId/interface/DetId.h" #include "DataFormats/MuonDetId/interface/MuonSubdetId.h" #include "DataFormats/MuonDetId/interface/DTChamberId.h" #include "DataFormats/MuonDetId/interface/DTSuperLayerId.h" #include "DataFormats/MuonDetId/interface/DTLayerId.h" #include "DataFormats/MuonDetId/interface/CSCDetId.h" #include "DataFormats/CLHEP/interface/AlgebraicObjects.h" // // class decleration // class MuonGeometrySanityCheckCustomFrame { public: MuonGeometrySanityCheckCustomFrame(const edm::ParameterSet &iConfig, std::string name); GlobalPoint transform(GlobalPoint point) const; GlobalPoint transformInverse(GlobalPoint point) const; AlgebraicMatrix matrix; AlgebraicMatrix matrixInverse; }; class MuonGeometrySanityCheckPoint { public: MuonGeometrySanityCheckPoint(const edm::ParameterSet &iConfig, const std::map<std::string, const MuonGeometrySanityCheckCustomFrame *> &frames); enum { kDTChamber, kDTSuperLayer, kDTLayer, kCSCChamber, kCSCLayer }; enum { kGlobal, kLocal, kChamber, kCustom }; std::string detName() const; int type; DetId detector; int frame; const MuonGeometrySanityCheckCustomFrame *customFrame; GlobalPoint displacement; bool has_expectation; GlobalPoint expectation; std::string name; int outputFrame; const MuonGeometrySanityCheckCustomFrame *outputCustomFrame; private: bool numeric(std::string s); int number(std::string s); }; class MuonGeometrySanityCheck : public edm::one::EDAnalyzer<> { public: explicit MuonGeometrySanityCheck(const edm::ParameterSet &iConfig); ~MuonGeometrySanityCheck() override; private: void analyze(const edm::Event &, const edm::EventSetup &iConfig) override; std::string printout; double tolerance; std::string prefix; std::map<std::string, const MuonGeometrySanityCheckCustomFrame *> m_frames; std::vector<MuonGeometrySanityCheckPoint> m_points; const edm::ESGetToken<DTGeometry, MuonGeometryRecord> dtGeomToken_; const edm::ESGetToken<CSCGeometry, MuonGeometryRecord> cscGeomToken_; }; // // constants, enums and typedefs // // // static data member definitions // // // constructors and destructor // MuonGeometrySanityCheck::MuonGeometrySanityCheck(const edm::ParameterSet &iConfig) : dtGeomToken_(esConsumes<edm::Transition::BeginRun>()), cscGeomToken_(esConsumes<edm::Transition::BeginRun>()) { printout = iConfig.getParameter<std::string>("printout"); if (printout != std::string("all") && printout != std::string("bad")) { throw cms::Exception("BadConfig") << "Printout must be \"all\" or \"bad\"." << std::endl; } tolerance = iConfig.getParameter<double>("tolerance"); if (tolerance <= 0) { throw cms::Exception("BadConfig") << "Tolerance must be positive." << std::endl; } prefix = iConfig.getParameter<std::string>("prefix"); std::vector<edm::ParameterSet> frames = iConfig.getParameter<std::vector<edm::ParameterSet> >("frames"); for (std::vector<edm::ParameterSet>::const_iterator frame = frames.begin(); frame != frames.end(); ++frame) { std::string name = frame->getParameter<std::string>("name"); if (m_frames.find(name) != m_frames.end()) { throw cms::Exception("BadConfig") << "Custom frame \"" << name << "\" has been defined twice." << std::endl; } m_frames[name] = new MuonGeometrySanityCheckCustomFrame(*frame, name); } std::vector<edm::ParameterSet> points = iConfig.getParameter<std::vector<edm::ParameterSet> >("points"); for (std::vector<edm::ParameterSet>::const_iterator point = points.begin(); point != points.end(); ++point) { m_points.push_back(MuonGeometrySanityCheckPoint(*point, m_frames)); } } MuonGeometrySanityCheck::~MuonGeometrySanityCheck() { for (std::map<std::string, const MuonGeometrySanityCheckCustomFrame *>::iterator iter = m_frames.begin(); iter != m_frames.end(); ++iter) { delete iter->second; } } MuonGeometrySanityCheckCustomFrame::MuonGeometrySanityCheckCustomFrame(const edm::ParameterSet &iConfig, std::string name) { std::vector<double> numbers = iConfig.getParameter<std::vector<double> >("matrix"); if (numbers.size() != 9) { throw cms::Exception("BadConfig") << "Custom frame \"" << name << "\" has a matrix which is not 3x3." << std::endl; } matrix = AlgebraicMatrix(3, 3); matrix[0][0] = numbers[0]; matrix[0][1] = numbers[1]; matrix[0][2] = numbers[2]; matrix[1][0] = numbers[3]; matrix[1][1] = numbers[4]; matrix[1][2] = numbers[5]; matrix[2][0] = numbers[6]; matrix[2][1] = numbers[7]; matrix[2][2] = numbers[8]; int ierr; matrixInverse = matrix; matrixInverse.invert(ierr); if (ierr != 0) { throw cms::Exception("BadConfig") << "Could not invert matrix for custom frame \"" << name << "\"." << std::endl; } } GlobalPoint MuonGeometrySanityCheckCustomFrame::transform(GlobalPoint point) const { AlgebraicVector input(3); input[0] = point.x(); input[1] = point.x(); input[2] = point.x(); AlgebraicVector output = matrix * input; return GlobalPoint(output[0], output[1], output[3]); } GlobalPoint MuonGeometrySanityCheckCustomFrame::transformInverse(GlobalPoint point) const { AlgebraicVector input(3); input[0] = point.x(); input[1] = point.x(); input[2] = point.x(); AlgebraicVector output = matrixInverse * input; return GlobalPoint(output[0], output[1], output[3]); } bool MuonGeometrySanityCheckPoint::numeric(std::string s) { return (s == std::string("0") || s == std::string("1") || s == std::string("2") || s == std::string("3") || s == std::string("4") || s == std::string("5") || s == std::string("6") || s == std::string("7") || s == std::string("8") || s == std::string("9")); } int MuonGeometrySanityCheckPoint::number(std::string s) { if (s == std::string("0")) return 0; else if (s == std::string("1")) return 1; else if (s == std::string("2")) return 2; else if (s == std::string("3")) return 3; else if (s == std::string("4")) return 4; else if (s == std::string("5")) return 5; else if (s == std::string("6")) return 6; else if (s == std::string("7")) return 7; else if (s == std::string("8")) return 8; else if (s == std::string("9")) return 9; else assert(false); } MuonGeometrySanityCheckPoint::MuonGeometrySanityCheckPoint( const edm::ParameterSet &iConfig, const std::map<std::string, const MuonGeometrySanityCheckCustomFrame *> &frames) { std::string detName = iConfig.getParameter<std::string>("detector"); bool parsing_error = false; bool barrel = (detName.substr(0, 2) == std::string("MB")); bool endcap = (detName.substr(0, 2) == std::string("ME")); if (!barrel && !endcap) parsing_error = true; if (!parsing_error && barrel) { int index = 2; bool plus = true; if (detName.substr(index, 1) == std::string("+")) { plus = true; index++; } else if (detName.substr(index, 1) == std::string("-")) { plus = false; index++; } int wheel = 0; bool wheel_digit = false; while (!parsing_error && numeric(detName.substr(index, 1))) { wheel *= 10; wheel += number(detName.substr(index, 1)); wheel_digit = true; index++; } if (!plus) wheel *= -1; if (!wheel_digit) parsing_error = true; if (detName.substr(index, 1) != std::string("/")) parsing_error = true; index++; int station = 0; bool station_digit = false; while (!parsing_error && numeric(detName.substr(index, 1))) { station *= 10; station += number(detName.substr(index, 1)); station_digit = true; index++; } if (!station_digit) parsing_error = true; if (detName.substr(index, 1) != std::string("/")) parsing_error = true; index++; int sector = 0; bool sector_digit = false; while (!parsing_error && numeric(detName.substr(index, 1))) { sector *= 10; sector += number(detName.substr(index, 1)); sector_digit = true; index++; } if (!sector_digit) parsing_error = true; // these are optional int superlayer = 0; bool superlayer_digit = false; int layer = 0; if (detName.substr(index, 1) == std::string("/")) { index++; while (!parsing_error && numeric(detName.substr(index, 1))) { superlayer *= 10; superlayer += number(detName.substr(index, 1)); superlayer_digit = true; index++; } if (!superlayer_digit) parsing_error = true; if (detName.substr(index, 1) == std::string("/")) { index++; while (!parsing_error && numeric(detName.substr(index, 1))) { layer *= 10; layer += number(detName.substr(index, 1)); index++; } } } if (!parsing_error) { bool no_such_chamber = false; if (wheel < -2 || wheel > 2) no_such_chamber = true; if (station < 1 || station > 4) no_such_chamber = true; if (station == 4 && (sector < 1 || sector > 14)) no_such_chamber = true; if (station < 4 && (sector < 1 || sector > 12)) no_such_chamber = true; if (no_such_chamber) { throw cms::Exception("BadConfig") << "Chamber doesn't exist: MB" << (plus ? "+" : "-") << wheel << "/" << station << "/" << sector << std::endl; } if (superlayer == 0) { detector = DTChamberId(wheel, station, sector); type = kDTChamber; } else { bool no_such_superlayer = false; if (superlayer < 1 || superlayer > 3) no_such_superlayer = true; if (station == 4 && superlayer == 2) no_such_superlayer = true; if (no_such_superlayer) { throw cms::Exception("BadConfig") << "Superlayer doesn't exist: MB" << (plus ? "+" : "-") << wheel << "/" << station << "/" << sector << "/" << superlayer << std::endl; } if (layer == 0) { detector = DTSuperLayerId(wheel, station, sector, superlayer); type = kDTSuperLayer; } else { bool no_such_layer = false; if (layer < 1 || layer > 4) no_such_layer = true; if (no_such_layer) { throw cms::Exception("BadConfig") << "Layer doesn't exist: MB" << (plus ? "+" : "-") << wheel << "/" << station << "/" << sector << "/" << superlayer << "/" << layer << std::endl; } detector = DTLayerId(wheel, station, sector, superlayer, layer); type = kDTLayer; } } } } else if (!parsing_error && endcap) { int index = 2; bool plus = true; if (detName.substr(index, 1) == std::string("+")) { plus = true; index++; } else if (detName.substr(index, 1) == std::string("-")) { plus = false; index++; } else parsing_error = true; int station = 0; bool station_digit = false; while (!parsing_error && numeric(detName.substr(index, 1))) { station *= 10; station += number(detName.substr(index, 1)); station_digit = true; index++; } if (!plus) station *= -1; if (!station_digit) parsing_error = true; if (detName.substr(index, 1) != std::string("/")) parsing_error = true; index++; int ring = 0; bool ring_digit = false; while (!parsing_error && numeric(detName.substr(index, 1))) { ring *= 10; ring += number(detName.substr(index, 1)); ring_digit = true; index++; } if (!ring_digit) parsing_error = true; if (detName.substr(index, 1) != std::string("/")) parsing_error = true; index++; int chamber = 0; bool chamber_digit = false; while (!parsing_error && numeric(detName.substr(index, 1))) { chamber *= 10; chamber += number(detName.substr(index, 1)); chamber_digit = true; index++; } if (!chamber_digit) parsing_error = true; // this is optional int layer = 0; bool layer_digit = false; if (detName.substr(index, 1) == std::string("/")) { index++; while (!parsing_error && numeric(detName.substr(index, 1))) { layer *= 10; layer += number(detName.substr(index, 1)); layer_digit = true; index++; } if (!layer_digit) parsing_error = true; } if (!parsing_error) { bool no_such_chamber = false; int endcap = (station > 0 ? 1 : 2); station = abs(station); if (station < 1 || station > 4) no_such_chamber = true; if (station == 1 && (ring < 1 || ring > 4)) no_such_chamber = true; if (station > 1 && (ring < 1 || ring > 2)) no_such_chamber = true; if (station == 1 && (chamber < 1 || chamber > 36)) no_such_chamber = true; if (station > 1 && ring == 1 && (chamber < 1 || chamber > 18)) no_such_chamber = true; if (station > 1 && ring == 2 && (chamber < 1 || chamber > 36)) no_such_chamber = true; if (no_such_chamber) { throw cms::Exception("BadConfig") << "Chamber doesn't exist: ME" << (endcap == 1 ? "+" : "-") << station << "/" << ring << "/" << chamber << std::endl; } if (layer == 0) { detector = CSCDetId(endcap, station, ring, chamber); type = kCSCChamber; } else { bool no_such_layer = false; if (layer < 1 || layer > 6) no_such_layer = true; if (no_such_layer) { throw cms::Exception("BadConfig") << "Layer doesn't exist: ME" << (endcap == 1 ? "+" : "-") << station << "/" << ring << "/" << chamber << "/" << layer << std::endl; } detector = CSCDetId(endcap, station, ring, chamber, layer); type = kCSCLayer; } } } if (parsing_error) { throw cms::Exception("BadConfig") << "Detector name is malformed: " << detName << std::endl; } std::string frameName = iConfig.getParameter<std::string>("frame"); const std::map<std::string, const MuonGeometrySanityCheckCustomFrame *>::const_iterator frameIter = frames.find(frameName); if (frameName == std::string("global")) { frame = kGlobal; customFrame = nullptr; } else if (frameName == std::string("local")) { frame = kLocal; customFrame = nullptr; } else if (frameName == std::string("chamber")) { frame = kChamber; customFrame = nullptr; } else if (frameIter != frames.end()) { frame = kCustom; customFrame = frameIter->second; } else { throw cms::Exception("BadConfig") << "Frame \"" << frameName << "\" has not been defined." << std::endl; } std::vector<double> point = iConfig.getParameter<std::vector<double> >("displacement"); if (point.size() != 3) { throw cms::Exception("BadConfig") << "Displacement relative to detector " << detName << " doesn't have exactly three components." << std::endl; } displacement = GlobalPoint(point[0], point[1], point[2]); const edm::Entry *entry = iConfig.retrieveUnknown("expectation"); if (entry != nullptr) { has_expectation = true; point = iConfig.getParameter<std::vector<double> >("expectation"); if (point.size() != 3) { throw cms::Exception("BadConfig") << "Expectation for detector " << detName << ", displacement " << displacement << " doesn't have exactly three components." << std::endl; } expectation = GlobalPoint(point[0], point[1], point[2]); } else { has_expectation = false; } entry = iConfig.retrieveUnknown("name"); if (entry != nullptr) { name = iConfig.getParameter<std::string>("name"); } else { name = std::string("anonymous"); } entry = iConfig.retrieveUnknown("outputFrame"); if (entry != nullptr) { frameName = iConfig.getParameter<std::string>("outputFrame"); const std::map<std::string, const MuonGeometrySanityCheckCustomFrame *>::const_iterator frameIter = frames.find(frameName); if (frameName == std::string("global")) { outputFrame = kGlobal; outputCustomFrame = nullptr; } else if (frameName == std::string("local")) { outputFrame = kLocal; outputCustomFrame = nullptr; } else if (frameName == std::string("chamber")) { outputFrame = kChamber; outputCustomFrame = nullptr; } else if (frameIter != frames.end()) { outputFrame = kCustom; outputCustomFrame = frameIter->second; } else { throw cms::Exception("BadConfig") << "Frame \"" << frameName << "\" has not been defined." << std::endl; } } else { outputFrame = kGlobal; outputCustomFrame = nullptr; } } std::string MuonGeometrySanityCheckPoint::detName() const { std::stringstream output; if (type == kDTChamber) { DTChamberId id(detector); output << "MB" << (id.wheel() > 0 ? "+" : "") << id.wheel() << "/" << id.station() << "/" << id.sector(); } else if (type == kDTSuperLayer) { DTSuperLayerId id(detector); output << "MB" << (id.wheel() > 0 ? "+" : "") << id.wheel() << "/" << id.station() << "/" << id.sector() << "/" << id.superlayer(); } else if (type == kDTLayer) { DTLayerId id(detector); output << "MB" << (id.wheel() > 0 ? "+" : "") << id.wheel() << "/" << id.station() << "/" << id.sector() << "/" << id.superlayer() << "/" << id.layer(); } else if (type == kCSCChamber) { CSCDetId id(detector); output << "ME" << (id.endcap() == 1 ? "+" : "-") << id.station() << "/" << id.ring() << "/" << id.chamber(); } else if (type == kCSCLayer) { CSCDetId id(detector); output << "ME" << (id.endcap() == 1 ? "+" : "-") << id.station() << "/" << id.ring() << "/" << id.chamber() << "/" << id.layer(); } else assert(false); return output.str(); } // ------------ method called to for each event ------------ void MuonGeometrySanityCheck::analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup) { const DTGeometry *dtGeometry = &iSetup.getData(dtGeomToken_); const CSCGeometry *cscGeometry = &iSetup.getData(cscGeomToken_); int num_transformed = 0; int num_tested = 0; int num_bad = 0; for (std::vector<MuonGeometrySanityCheckPoint>::const_iterator point = m_points.begin(); point != m_points.end(); ++point) { num_transformed++; bool dt = (point->detector.subdetId() == MuonSubdetId::DT); // convert the displacement vector into the chosen coordinate system and add it to the chamber's position GlobalPoint chamberPos; if (dt) chamberPos = dtGeometry->idToDet(point->detector)->surface().toGlobal(LocalPoint(0., 0., 0.)); else chamberPos = cscGeometry->idToDet(point->detector)->surface().toGlobal(LocalPoint(0., 0., 0.)); GlobalPoint result; if (point->frame == MuonGeometrySanityCheckPoint::kGlobal) { result = GlobalPoint(chamberPos.x() + point->displacement.x(), chamberPos.y() + point->displacement.y(), chamberPos.z() + point->displacement.z()); } else if (point->frame == MuonGeometrySanityCheckPoint::kLocal) { if (dt) result = dtGeometry->idToDet(point->detector) ->surface() .toGlobal(LocalPoint(point->displacement.x(), point->displacement.y(), point->displacement.z())); else result = cscGeometry->idToDet(point->detector) ->surface() .toGlobal(LocalPoint(point->displacement.x(), point->displacement.y(), point->displacement.z())); } else if (point->frame == MuonGeometrySanityCheckPoint::kChamber) { if (point->detector.subdetId() == MuonSubdetId::DT) { DTChamberId id(point->detector); if (dt) result = dtGeometry->idToDet(id)->surface().toGlobal( LocalPoint(point->displacement.x(), point->displacement.y(), point->displacement.z())); else result = cscGeometry->idToDet(id)->surface().toGlobal( LocalPoint(point->displacement.x(), point->displacement.y(), point->displacement.z())); } else if (point->detector.subdetId() == MuonSubdetId::CSC) { CSCDetId cscid(point->detector); CSCDetId id(cscid.endcap(), cscid.station(), cscid.ring(), cscid.chamber()); if (dt) result = dtGeometry->idToDet(id)->surface().toGlobal( LocalPoint(point->displacement.x(), point->displacement.y(), point->displacement.z())); else result = cscGeometry->idToDet(id)->surface().toGlobal( LocalPoint(point->displacement.x(), point->displacement.y(), point->displacement.z())); } else { assert(false); } } else if (point->frame == MuonGeometrySanityCheckPoint::kCustom) { GlobalPoint transformed = point->customFrame->transform(point->displacement); result = GlobalPoint( chamberPos.x() + transformed.x(), chamberPos.y() + transformed.y(), chamberPos.z() + transformed.z()); } else { assert(false); } // convert the result into the chosen output coordinate system if (point->outputFrame == MuonGeometrySanityCheckPoint::kGlobal) { } else if (point->outputFrame == MuonGeometrySanityCheckPoint::kLocal) { LocalPoint transformed; if (dt) transformed = dtGeometry->idToDet(point->detector)->surface().toLocal(result); else transformed = cscGeometry->idToDet(point->detector)->surface().toLocal(result); result = GlobalPoint(transformed.x(), transformed.y(), transformed.z()); } else if (point->outputFrame == MuonGeometrySanityCheckPoint::kChamber) { if (point->detector.subdetId() == MuonSubdetId::DT) { DTChamberId id(point->detector); LocalPoint transformed; if (dt) transformed = dtGeometry->idToDet(id)->surface().toLocal(result); else transformed = cscGeometry->idToDet(id)->surface().toLocal(result); result = GlobalPoint(transformed.x(), transformed.y(), transformed.z()); } else if (point->detector.subdetId() == MuonSubdetId::CSC) { CSCDetId cscid(point->detector); CSCDetId id(cscid.endcap(), cscid.station(), cscid.ring(), cscid.chamber()); LocalPoint transformed; if (dt) transformed = dtGeometry->idToDet(id)->surface().toLocal(result); else transformed = cscGeometry->idToDet(id)->surface().toLocal(result); result = GlobalPoint(transformed.x(), transformed.y(), transformed.z()); } else { assert(false); } } else if (point->outputFrame == MuonGeometrySanityCheckPoint::kCustom) { result = point->outputCustomFrame->transformInverse(result); } std::stringstream output; output << prefix << " " << point->name << " " << point->detName() << " " << result.x() << " " << result.y() << " " << result.z(); bool bad = false; if (point->has_expectation) { num_tested++; double residx = result.x() - point->expectation.x(); double residy = result.y() - point->expectation.y(); double residz = result.z() - point->expectation.z(); if (fabs(residx) > tolerance || fabs(residy) > tolerance || fabs(residz) > tolerance) { num_bad++; bad = true; output << " BAD " << residx << " " << residy << " " << residz << std::endl; } else { output << " GOOD " << residx << " " << residy << " " << residz << std::endl; } } else { output << " UNTESTED 0 0 0" << std::endl; } if (printout == std::string("all") || (printout == std::string("bad") && bad)) { std::cout << output.str(); } } std::cout << std::endl << "SUMMARY transformed: " << num_transformed << " tested: " << num_tested << " bad: " << num_bad << " good: " << (num_tested - num_bad) << std::endl; } //define this as a plug-in DEFINE_FWK_MODULE(MuonGeometrySanityCheck);