File indexing completed on 2024-04-06 11:57:57
0001
0002
0003
0004
0005
0006 #include <memory>
0007
0008
0009 #include "FWCore/Framework/interface/global/EDAnalyzer.h"
0010 #include "FWCore/Framework/interface/MakerMacros.h"
0011 #include "FWCore/Framework/interface/EventSetup.h"
0012
0013 #include "CalibCalorimetry/EcalLaserCorrection/interface/EcalLaserDbService.h"
0014 #include "CalibCalorimetry/EcalLaserCorrection/interface/EcalLaserDbRecord.h"
0015
0016 #include "CondFormats/EcalObjects/interface/EcalLaserAlphas.h"
0017 #include "CondFormats/EcalObjects/interface/EcalLaserAPDPNRatiosRef.h"
0018 #include "CondFormats/EcalObjects/interface/EcalLaserAPDPNRatios.h"
0019 #include "CondFormats/EcalObjects/interface/EcalLinearCorrections.h"
0020
0021 #include "DataFormats/EcalDetId/interface/EBDetId.h"
0022
0023 class EcalLaserDbAnalyzer : public edm::global::EDAnalyzer<> {
0024 public:
0025 explicit EcalLaserDbAnalyzer(const edm::ParameterSet&);
0026 ~EcalLaserDbAnalyzer() override = default;
0027
0028 void analyze(edm::StreamID, edm::Event const&, edm::EventSetup const&) const override;
0029
0030 private:
0031
0032 const edm::ESGetToken<EcalLaserDbService, EcalLaserDbRecord> laserDbToken_;
0033 };
0034
0035
0036
0037
0038
0039
0040
0041
0042
0043
0044
0045
0046 EcalLaserDbAnalyzer::EcalLaserDbAnalyzer(const edm::ParameterSet& iConfig) : laserDbToken_(esConsumes()) {}
0047
0048
0049
0050
0051
0052
0053 void EcalLaserDbAnalyzer::analyze(edm::StreamID, edm::Event const& iEvent, edm::EventSetup const& iSetup) const {
0054
0055 const auto& setup = iSetup.getData(laserDbToken_);
0056 edm::LogInfo("EcalLaserDbService") << "EcalLaserDbAnalyzer::analyze-> got EcalLaserDbRecord:";
0057
0058 EcalLaserAPDPNRatios::EcalLaserAPDPNpair apdpnpair;
0059 const EcalLaserAPDPNRatios* myapdpn = setup.getAPDPNRatios();
0060 const EcalLaserAPDPNRatios::EcalLaserAPDPNRatiosMap& laserRatiosMap = myapdpn->getLaserMap();
0061
0062
0063
0064
0065 EcalLaserAPDPNref apdpnref;
0066 const EcalLaserAPDPNRatiosRef* myapdpnref = setup.getAPDPNRatiosRef();
0067 const EcalLaserAPDPNRatiosRefMap& laserRefMap = myapdpnref->getMap();
0068
0069 EcalLaserAlpha alpha;
0070 const EcalLaserAlphas* myalpha = setup.getAlphas();
0071 const EcalLaserAlphaMap& laserAlphaMap = myalpha->getMap();
0072
0073 EcalLinearCorrections::Values linValues;
0074 const EcalLinearCorrections* mylinear = setup.getLinearCorrections();
0075 const EcalLinearCorrections::EcalValueMap& linearValueMap = mylinear->getValueMap();
0076
0077 for (int ieta = -EBDetId::MAX_IETA; ieta <= EBDetId::MAX_IETA; ++ieta) {
0078 if (ieta == 0)
0079 continue;
0080 for (int iphi = EBDetId::MIN_IPHI; iphi <= EBDetId::MAX_IPHI; ++iphi) {
0081 EBDetId ebdetid(ieta, iphi);
0082
0083 edm::LogVerbatim("EcalLaserDbService") << ebdetid << " " << ebdetid.ietaSM() << " " << ebdetid.iphiSM() << "\n";
0084
0085 EcalLaserAPDPNRatios::EcalLaserAPDPNRatiosMap::const_iterator itratio = laserRatiosMap.find(ebdetid);
0086 if (itratio != laserRatiosMap.end()) {
0087 apdpnpair = (*itratio);
0088 edm::LogVerbatim("EcalLaserDbService")
0089 << " APDPN pair = " << apdpnpair.p1 << " , " << apdpnpair.p2 << " , " << apdpnpair.p3 << "\n";
0090 } else {
0091 edm::LogError("EcalLaserDbService") << "error with laserRatiosMap!";
0092 }
0093
0094 EcalLinearCorrections::EcalValueMap::const_iterator itlin = linearValueMap.find(ebdetid);
0095 if (itlin != linearValueMap.end()) {
0096 linValues = (*itlin);
0097 edm::LogVerbatim("EcalLaserDbService")
0098 << " APDPN pair = " << linValues.p1 << " , " << linValues.p2 << " , " << linValues.p3 << "\n";
0099 } else {
0100 edm::LogError("EcalLaserDbService") << "error with linearValuesMap!";
0101 }
0102
0103 EcalLaserAPDPNRatiosRefMap::const_iterator itref = laserRefMap.find(ebdetid);
0104 if (itref != laserRefMap.end()) {
0105 apdpnref = (*itref);
0106 edm::LogVerbatim("EcalLaserDbService") << " APDPN ref = " << apdpnref << "\n";
0107 } else {
0108 edm::LogError("EcalLaserDbService") << "error with laserRefMap!";
0109 }
0110
0111 EcalLaserAlphaMap::const_iterator italpha = laserAlphaMap.find(ebdetid);
0112 if (italpha != laserAlphaMap.end()) {
0113 alpha = (*italpha);
0114 edm::LogVerbatim("EcalLaserDbService") << " ALPHA = " << alpha << "\n";
0115 } else {
0116 edm::LogError("EcalLaserDbService") << "error with laserAlphaMap!";
0117 }
0118 }
0119 }
0120 }
0121
0122
0123 DEFINE_FWK_MODULE(EcalLaserDbAnalyzer);