ecaldqm::MLClient Class Reference

#include <MLClient.h>

Inheritance diagram for ecaldqm::MLClient:

Public Member Functions
	MLClient ()
void	producePlots (ProcessType) override
	~MLClient () override
Public Member Functions inherited from ecaldqm::DQWorkerClient
void	bookMEs (DQMStore::IBooker &) override
	DQWorkerClient ()
void	endLuminosityBlock (edm::LuminosityBlock const &, edm::EventSetup const &) override
void	releaseMEs () override
void	releaseSource ()
virtual void	resetMEs ()
void	resetPerLumi ()
bool	retrieveSource (DQMStore::IGetter &, ProcessType)
bool	runsOn (ProcessType _type) const
void	setStatusManager (StatusManager const &_manager)
virtual void	setTokens (edm::ConsumesCollector &)
	~DQWorkerClient () override
Public Member Functions inherited from ecaldqm::DQWorker
virtual void	beginLuminosityBlock (edm::LuminosityBlock const &, edm::EventSetup const &)
virtual void	beginRun (edm::Run const &, edm::EventSetup const &)
bool	checkElectronicsMap (bool=true)
bool	checkGeometry (bool=true)
bool	checkTopology (bool=true)
bool	checkTrigTowerMap (bool=true)
	DQWorker ()
virtual void	endRun (edm::Run const &, edm::EventSetup const &)
EcalDQMSetupObjects const	getEcalDQMSetupObjects ()
EcalElectronicsMapping const *	GetElectronicsMap ()
CaloGeometry const *	GetGeometry ()
std::string const &	getName () const
CaloTopology const *	GetTopology ()
EcalTrigTowerConstituentsMap const *	GetTrigTowerMap ()
bool	onlineMode () const
void	setEventNumber (edm::EventNumber_t _e)
void	setLumiNumber (edm::LuminosityBlockNumber_t _l)
void	setRunNumber (edm::RunNumber_t _r)
void	setSetupObjects (edm::EventSetup const &)
void	setSetupObjectsEndLumi (edm::EventSetup const &)
void	setTime (time_t _t)
void	setTokens (edm::ConsumesCollector &)
virtual	~DQWorker () noexcept(false)

Private Member Functions
void	setParams (edm::ParameterSet const &) override

Private Attributes
float	EB_PUcorr_intercept_
float	EB_PUcorr_slope_
std::vector< double >	EBavgOcc
std::deque< std::valarray< std::valarray< float > > >	EBlossMap2dQ
std::deque< std::valarray< float > >	ebOccMap1dQ
float	EBThreshold_
float	EEm_PUcorr_intercept_
float	EEm_PUcorr_slope_
std::vector< double >	EEmavgOcc
std::deque< std::valarray< std::valarray< float > > >	EEmlossMap2dQ
std::deque< std::valarray< float > >	eemOccMap1dQ
float	EEmThreshold_
float	EEp_PUcorr_intercept_
float	EEp_PUcorr_slope_
std::vector< double >	EEpavgOcc
std::deque< std::valarray< std::valarray< float > > >	EEplossMap2dQ
std::deque< std::valarray< float > >	eepOccMap1dQ
float	EEpThreshold_
int	LScount = 0
int	nbadtowerEB
int	nbadtowerEE
std::deque< int >	NEventQ
size_t	nLS = 4
size_t	nLSloss = 6

Static Private Attributes
static const int	nEBEtaTowers = 34
static const int	nEBEtaTowersPad = 36
static const int	nEBPhiTowers = 72
static const int	nEEEtaTowers = 20
static const int	nEEPhiTowers = 20
static const int	nEETowersPad = 22

Additional Inherited Members
Public Types inherited from ecaldqm::DQWorkerClient
enum	ProcessType { kLumi, kJob, nProcessType }
enum	Quality {   kBad = 0, kGood = 1, kUnknown = 2, kMBad = 3,   kMGood = 4, kMUnknown = 5 }
Static Public Member Functions inherited from ecaldqm::DQWorkerClient
static void	fillDescriptions (edm::ParameterSetDescription &)
Static Public Member Functions inherited from ecaldqm::DQWorker
static void	fillDescriptions (edm::ParameterSetDescription &_desc)
Public Attributes inherited from ecaldqm::DQWorker
edm::ESGetToken< EcalElectronicsMapping, EcalMappingRcd >	elecMapHandle
edm::ESGetToken< EcalElectronicsMapping, EcalMappingRcd >	elecMapHandleEndLumi
edm::ESGetToken< CaloGeometry, CaloGeometryRecord >	geomHandle
edm::ESGetToken< CaloGeometry, CaloGeometryRecord >	geomHandleEndLumi
edm::ESGetToken< CaloTopology, CaloTopologyRecord >	topoHandle
edm::ESGetToken< CaloTopology, CaloTopologyRecord >	topoHandleEndLumi
edm::ESGetToken< EcalTrigTowerConstituentsMap, IdealGeometryRecord >	ttMapHandle
edm::ESGetToken< EcalTrigTowerConstituentsMap, IdealGeometryRecord >	ttMapHandleEndLumi
Protected Types inherited from ecaldqm::DQWorker
typedef dqm::legacy::DQMStore	DQMStore
typedef dqm::legacy::MonitorElement	MonitorElement
Protected Member Functions inherited from ecaldqm::DQWorkerClient
void	setME (edm::ParameterSet const &_ps) final
void	setSource (edm::ParameterSet const &) override
void	towerAverage_ (MESet &, MESet const &, float)
bool	using_ (std::string const &_name, ProcessType _type=kJob) const
Protected Member Functions inherited from ecaldqm::DQWorker
void	initialize (std::string const &_name, edm::ParameterSet const &)
void	print_ (std::string const &, int=0) const
void	setVerbosity (int _verbosity)
Protected Attributes inherited from ecaldqm::DQWorkerClient
bool	hasLumiPlots_
std::set< std::string >	qualitySummaries_
MESetCollection	sources_
StatusManager const *	statusManager_
Protected Attributes inherited from ecaldqm::DQWorker
bool	booked_
MESetCollection	MEs_
std::string	name_
bool	onlineMode_
Timestamp	timestamp_
int	verbosity_
bool	willConvertToEDM_

Detailed Description

Definition at line 11 of file MLClient.h.

Constructor & Destructor Documentation

MLClient()

ecaldqm::MLClient::MLClient

(

)

Definition at line 19 of file MLClient.cc.

References ecaldqm::DQWorkerClient::qualitySummaries_.

: DQWorkerClient() { qualitySummaries_.insert("MLQualitySummary"); }

~MLClient()

ecaldqm::MLClient::~MLClient ( )

inlineoverride

Definition at line 14 of file MLClient.h.

{}

Member Function Documentation

producePlots()

void ecaldqm::MLClient::producePlots ( ProcessType  )

overridevirtual

so that real anomalies which persist with time are enhanced and fluctuations are suppressed.

ONNX model running///

Inference on the output from the model///

Implements ecaldqm::DQWorkerClient.

Definition at line 44 of file MLClient.cc.

References ecaldqm::MESetCollection::at(), newFWLiteAna::bin, Skims_PA_cff::content, EB_PUcorr_intercept_, EB_PUcorr_slope_, EBavgOcc, EBlossMap2dQ, ebOccMap1dQ, EBThreshold_, EcalBarrel, EcalEndcap, EEm_PUcorr_intercept_, EEm_PUcorr_slope_, EEmavgOcc, EEmlossMap2dQ, eemOccMap1dQ, EEmThreshold_, EEp_PUcorr_intercept_, EEp_PUcorr_slope_, EEpavgOcc, EEplossMap2dQ, eepOccMap1dQ, EEpThreshold_, web.browse_db::env, edm::FileInPath::fullPath(), ecaldqm::DQWorker::getEcalDQMSetupObjects(), ecaldqm::DQWorker::GetElectronicsMap(), ecaldqm::MESetNonObject::getFloatValue(), ecaldqm::DQWorker::GetTrigTowerMap(), mps_fire::i, l1ctLayer2EG_cff::id, ecaldqm::DQWorker::Timestamp::iLumi, testHGCalDigi_cfg::inFile, dqmiolumiharvest::j, dqmdumpme::k, ecaldqm::DQWorkerClient::kBad, ecaldqm::DQWorkerClient::kGood, ecaldqm::DQWorkerClient::kMBad, ecaldqm::DQWorkerClient::kMGood, LScount, ALPAKA_ACCELERATOR_NAMESPACE::pixelClustering::pixelStatus::mask, ecaldqm::DQWorker::MEs_, nbadtowerEB, nbadtowerEE, nEBEtaTowers, nEBEtaTowersPad, nEBPhiTowers, nEEEtaTowers, nEEPhiTowers, nEETowersPad, NEventQ, nLS, nLSloss, ecaldqm::DQWorker::onlineMode_, makeMuonMisalignmentScenario::outputName, jetsAK4_CHS_cff::outputNames, EcalDQMStatusHelper::PEDESTAL_ONLINE_HIGH_GAIN_RMS_ERROR, EcalDQMStatusHelper::PHYSICS_BAD_CHANNEL_ERROR, EcalDQMStatusHelper::PHYSICS_BAD_CHANNEL_WARNING, funct::pow(), displacedMuons_cfi::pu, quality, ecaldqm::DQWorkerClient::sources_, ecaldqm::DQWorkerClient::statusManager_, AlCaHLTBitMon_QueryRunRegistry::string, ecaldqm::DQWorker::timestamp_, ecaldqm::MESet::const_iterator::toNextChannel(), heppy_batch::val, and EEDetId::zside().

                                         {
    if (!onlineMode_)
      return;
    nbadtowerEB = 0;
    nbadtowerEE = 0;

    using namespace std;
    MESet& meMLQualitySummary(MEs_.at("MLQualitySummary"));
    MESet& meEventsperMLImage(MEs_.at("EventsperMLImage"));

    MESetNonObject const& sPU(static_cast<MESetNonObject&>(sources_.at("PU")));
    MESetNonObject const& sNumEvents(static_cast<MESetNonObject&>(sources_.at("NumEvents")));

    //Get the no.of events and the PU per LS calculated in OccupancyTask
    int nEv = sNumEvents.getFloatValue();
    double pu = sPU.getFloatValue();

    //Do not compute ML quality if PU is non existent.
    if (pu <= 0.) {
      return;
    }
    uint32_t mask(1 << EcalDQMStatusHelper::PEDESTAL_ONLINE_HIGH_GAIN_RMS_ERROR |
                  1 << EcalDQMStatusHelper::PHYSICS_BAD_CHANNEL_WARNING |
                  1 << EcalDQMStatusHelper::PHYSICS_BAD_CHANNEL_ERROR);

    //Inorder to feed the data into the ML model we apply some preprocessing.
    //We use the Digi Occupancy per Lumisection as the input source.
    //The model was trained on each occupancy plot having 500 events.
    //In apprehension of the low luminosity in the beginning of Run3, where in online DQM
    //the no.of events per LS could be lower than 500, we sum the occupancies over a fixed no.of lumisections as a running sum,
    //and require that the total no.of events on this summed occupancy to be atleast 200.
    //(This no.of LS and the no.of events are parameters which would require tuning later)
    //This summed occupancy is now the input image, which is then corrected for PileUp(PU) dependence and
    //change in no.of events, which are derived from training.
    //The input image is also padded by replicating the top and bottom rows so as to prevent the "edge effect"
    //wherein the ML model's learning degrades near the edge of the data set it sees.
    //This padding is then removed during inference on the model output.

    //Get the histogram of the input digi occupancy per lumisection.
    TH2F* hEEmDigiMap((sources_.at("DigiAllByLumi")).getME(0)->getTH2F());
    TH2F* hEbDigiMap((sources_.at("DigiAllByLumi")).getME(1)->getTH2F());
    TH2F* hEEpDigiMap((sources_.at("DigiAllByLumi")).getME(2)->getTH2F());

    size_t nEBTowers = nEBEtaTowers * nEBPhiTowers;  //Each EB occupancy map is of size 34x72 towers
    size_t nEETowers = nEEEtaTowers * nEEPhiTowers;  //Each EE occupancy map is of size 20x20 towers

    //Vectors to feed into the ML network
    std::vector<float> ebOccMap1dCumulPad;
    std::vector<float> eemOccMap1dCumulPad;
    std::vector<float> eepOccMap1dCumulPad;

    //Array to store occupancy maps
    std::valarray<float> ebOccMap1d(nEBTowers);
    std::valarray<float> eemOccMap1d(nEETowers);
    std::valarray<float> eepOccMap1d(nEETowers);

    //Store the values from the input histogram into the array
    //to do preprocessing
    for (int i = 0; i < hEbDigiMap->GetNbinsY(); i++) {  //NbinsY = 34, NbinsX = 72
      for (int j = 0; j < hEbDigiMap->GetNbinsX(); j++) {
        int bin = hEbDigiMap->GetBin(j + 1, i + 1);
        int k = (i * nEBPhiTowers) + j;
        ebOccMap1d[k] = hEbDigiMap->GetBinContent(bin);
      }
    }
    ebOccMap1dQ.push_back(ebOccMap1d);  //Queue which stores input occupancy maps for nLS lumis

    for (int i = 0; i < hEEpDigiMap->GetNbinsY(); i++) {  //NbinsY = 20, NbinsX = 20
      for (int j = 0; j < hEEpDigiMap->GetNbinsX(); j++) {
        int bin = hEEpDigiMap->GetBin(j + 1, i + 1);
        int k = (i * nEEPhiTowers) + j;
        eemOccMap1d[k] = hEEmDigiMap->GetBinContent(bin);
        eepOccMap1d[k] = hEEpDigiMap->GetBinContent(bin);
      }
    }

    //Queue which stores input occupancy maps for nLS lumis
    eemOccMap1dQ.push_back(eemOccMap1d);
    eepOccMap1dQ.push_back(eepOccMap1d);

    NEventQ.push_back(nEv);  //Queue which stores the no.of events per LS for nLS lumis

    if (NEventQ.size() < nLS) {
      return;  //Should have nLS lumis to add the occupancy over.
    }
    if (NEventQ.size() > nLS) {
      NEventQ.pop_front();  //Keep only nLS consecutive LS. Pop the first one if size greater than nLS
    }
    if (ebOccMap1dQ.size() > nLS) {
      ebOccMap1dQ.pop_front();  //Same conditon for the input occupancy maps.
      eemOccMap1dQ.pop_front();
      eepOccMap1dQ.pop_front();
    }

    int TNum = 0;
    for (size_t i = 0; i < nLS; i++) {
      TNum += NEventQ[i];  //Total no.of events over nLS lumis
    }

    if (TNum < 400) {
      return;  //The total no.of events should be atleast 400 over nLS for meaningful statistics
    }
    //Fill the ME to monitor the trend of the total no.of events in each input image to the ML model
    meEventsperMLImage.fill(getEcalDQMSetupObjects(), EcalBarrel, double(timestamp_.iLumi), double(TNum));

    //Array to hold the sum of inputs, which make atleast 400 events.
    std::valarray<float> ebOccMap1dCumul(0., nEBTowers);
    std::valarray<float> eemOccMap1dCumul(0., nEETowers);
    std::valarray<float> eepOccMap1dCumul(0., nEETowers);

    //Sum the input arrays of nLS.
    for (size_t i = 0; i < ebOccMap1dQ.size(); i++) {
      ebOccMap1dCumul += ebOccMap1dQ[i];
      eemOccMap1dCumul += eemOccMap1dQ[i];
      eepOccMap1dCumul += eepOccMap1dQ[i];
    }

    //Applying PU correction derived from training
    ebOccMap1dCumul = ebOccMap1dCumul / (EB_PUcorr_slope_ * pu + EB_PUcorr_intercept_);
    eemOccMap1dCumul = eemOccMap1dCumul / (EEm_PUcorr_slope_ * pu + EEm_PUcorr_intercept_);
    eepOccMap1dCumul = eepOccMap1dCumul / (EEp_PUcorr_slope_ * pu + EEp_PUcorr_intercept_);

    //Scaling up to match input dimensions.
    ebOccMap1dCumul = ebOccMap1dCumul * (nEBEtaTowers * nEBPhiTowers);
    eemOccMap1dCumul = eemOccMap1dCumul * nEEEtaTowers * nEEPhiTowers;  //(nEETowersPad * nEETowersPad);
    eepOccMap1dCumul = eepOccMap1dCumul * nEEEtaTowers * nEEPhiTowers;  //(nEETowersPad * nEETowersPad);

    //Correction for no.of events in each input image as originally model trained with 500 events per image
    ebOccMap1dCumul = ebOccMap1dCumul * (500. / TNum);
    eemOccMap1dCumul = eemOccMap1dCumul * (500. / TNum);
    eepOccMap1dCumul = eepOccMap1dCumul * (500. / TNum);

    std::vector<std::vector<float>> ebOccMap2dCumul(nEBEtaTowers, std::vector<float>(nEBPhiTowers, 0.));
    //Convert 1dCumul to 2d
    for (size_t i = 0; i < nEBEtaTowers; i++) {
      for (size_t j = 0; j < nEBPhiTowers; j++) {
        int k = (i * nEBPhiTowers) + j;
        ebOccMap2dCumul[i][j] = ebOccMap1dCumul[k];
      }
    }

    std::vector<float> pad_top;
    std::vector<float> pad_bottom;
    std::vector<float> pad_left;
    std::vector<float> pad_right;

    pad_top = ebOccMap2dCumul[0];
    pad_bottom = ebOccMap2dCumul[ebOccMap2dCumul.size() - 1];

    ebOccMap2dCumul.insert(ebOccMap2dCumul.begin(), pad_top);
    ebOccMap2dCumul.push_back(pad_bottom);

    std::vector<std::vector<float>> eemOccMap2dCumul(nEEEtaTowers, std::vector<float>(nEEPhiTowers, 0.));
    std::vector<std::vector<float>> eepOccMap2dCumul(nEEEtaTowers, std::vector<float>(nEEPhiTowers, 0.));

    for (size_t i = 0; i < nEEEtaTowers; i++) {
      for (size_t j = 0; j < nEEPhiTowers; j++) {
        int k = (i * nEEPhiTowers) + j;
        eemOccMap2dCumul[i][j] = eemOccMap1dCumul[k];
        eepOccMap2dCumul[i][j] = eepOccMap1dCumul[k];
      }
    }

    // EE - //
    pad_top.clear();
    pad_bottom.clear();
    pad_left.clear();
    pad_right.clear();

    pad_top = eemOccMap2dCumul[0];
    pad_bottom = eemOccMap2dCumul[eemOccMap2dCumul.size() - 1];

    eemOccMap2dCumul.insert(eemOccMap2dCumul.begin(), pad_top);
    eemOccMap2dCumul.push_back(pad_bottom);

    for (auto& row : eemOccMap2dCumul) {
      pad_left.push_back(row[0]);
      pad_right.push_back(row[row.size() - 1]);
    }

    std::size_t Lindex = 0;
    std::size_t Rindex = 0;

    for (auto& row : eemOccMap2dCumul) {
      row.insert(row.begin(), pad_left[Lindex++]);
      row.insert(row.end(), pad_right[Rindex++]);
    }

    // EE + //
    pad_top.clear();
    pad_bottom.clear();

    pad_top = eepOccMap2dCumul[0];
    pad_bottom = eepOccMap2dCumul[eepOccMap2dCumul.size() - 1];

    eepOccMap2dCumul.insert(eepOccMap2dCumul.begin(), pad_top);
    eepOccMap2dCumul.push_back(pad_bottom);

    for (auto& row : eepOccMap2dCumul) {
      pad_left.push_back(row[0]);
      pad_right.push_back(row[row.size() - 1]);
    }

    Lindex = 0;
    Rindex = 0;

    for (auto& row : eepOccMap2dCumul) {
      row.insert(row.begin(), pad_left[Lindex++]);
      row.insert(row.end(), pad_right[Rindex++]);
    }

    //The pre-processed input is now fed into the 1D input tensor vector which will go into the ML model
    for (auto& row : ebOccMap2dCumul) {
      ebOccMap1dCumulPad.insert(ebOccMap1dCumulPad.end(), row.begin(), row.end());
    }

    for (auto& row : eemOccMap2dCumul) {
      eemOccMap1dCumulPad.insert(eemOccMap1dCumulPad.end(), row.begin(), row.end());
    }

    for (auto& row : eepOccMap2dCumul) {
      eepOccMap1dCumulPad.insert(eepOccMap1dCumulPad.end(), row.begin(), row.end());
    }

    //An Autoencoder (AE) network with resnet architecture is used here which is trained on
    //certified good data (EB digi occupancy) from Run 2018 data.
    //On giving an input occupancy map, the encoder part of the AE compresses and reduces the input data, learning its features,
    //and the decoder reconstructs the data from the encoded form into a representation as close to the original input as possible.
    //We then compute the Mean squared error (MSE) between the input and output image, also called the Reconstruction Loss,
    //calculated at a tower by tower basis.
    //Thus, given an anomalous tower the loss should be significantly higher than the loss with respect to good towers, which the model
    //has already seen --> anomaly detection.
    //When calculating the loss we also apply a response correction by dividing each input and output image with the average occupancy from
    //all 2018 data (also to be tuned),to accommodate the difference in response of crystals in different regions of the Ecal Barrel
    //Further each loss map from each input image is then multiplied by the last N loss maps,
    //A quality threshold is then applied on this time multiplied loss map, to mark them as GOOD or BAD,
    //after which it is stored as a quality summary ME.

    std::string instanceName{"AE-DQM-inference"};
    std::string modelFilepath = edm::FileInPath("DQM/EcalMonitorClient/data/onnxModels/resnet.onnx").fullPath();

    Ort::SessionOptions sessionOptions;
    sessionOptions.SetIntraOpNumThreads(1);
    Ort::Env env(OrtLoggingLevel::ORT_LOGGING_LEVEL_WARNING, instanceName.c_str());
    Ort::Session session(env, modelFilepath.c_str(), sessionOptions);

    Ort::AllocatorWithDefaultOptions allocator;

    // Strings returned by session.GetInputNameAllocated are temporary, need to copy them before they are deallocated
    std::string inputName{session.GetInputNameAllocated(0, allocator).get()};

    Ort::TypeInfo inputTypeInfo = session.GetInputTypeInfo(0);
    auto inputTensorInfo = inputTypeInfo.GetTensorTypeAndShapeInfo();

    std::vector<int64_t> inputDims = inputTensorInfo.GetShape();

    std::string outputName{session.GetOutputNameAllocated(0, allocator).get()};

    Ort::TypeInfo outputTypeInfo = session.GetOutputTypeInfo(0);
    auto outputTensorInfo = outputTypeInfo.GetTensorTypeAndShapeInfo();

    std::vector<int64_t> outputDims = outputTensorInfo.GetShape();

    size_t TensorSize = nEBEtaTowersPad * nEBPhiTowers;
    std::vector<float> ebRecoOccMap1dPad(TensorSize);  //To store the output reconstructed occupancy

    std::vector<const char*> inputNames{inputName.c_str()};
    std::vector<const char*> outputNames{outputName.c_str()};
    std::vector<Ort::Value> inputTensors;
    std::vector<Ort::Value> outputTensors;

    Ort::MemoryInfo memoryInfo =
        Ort::MemoryInfo::CreateCpu(OrtAllocatorType::OrtArenaAllocator, OrtMemType::OrtMemTypeDefault);
    inputTensors.push_back(Ort::Value::CreateTensor<float>(
        memoryInfo, ebOccMap1dCumulPad.data(), TensorSize, inputDims.data(), inputDims.size()));

    outputTensors.push_back(Ort::Value::CreateTensor<float>(
        memoryInfo, ebRecoOccMap1dPad.data(), TensorSize, outputDims.data(), outputDims.size()));

    session.Run(Ort::RunOptions{nullptr},
                inputNames.data(),
                inputTensors.data(),
                1,
                outputNames.data(),
                outputTensors.data(),
                1);

    //Endcaps
    // EE- //

    inputDims.clear();
    outputDims.clear();
    inputNames.clear();
    outputNames.clear();
    inputTensors.clear();
    outputTensors.clear();

    modelFilepath = edm::FileInPath("DQM/EcalMonitorClient/data/onnxModels/EEm_resnet2018.onnx").fullPath();

    Ort::Session EEm_session(env, modelFilepath.c_str(), sessionOptions);

    inputName = EEm_session.GetInputNameAllocated(0, allocator).get();

    inputTypeInfo = EEm_session.GetInputTypeInfo(0);
    auto EEm_inputTensorInfo = inputTypeInfo.GetTensorTypeAndShapeInfo();

    inputDims = EEm_inputTensorInfo.GetShape();

    outputName = EEm_session.GetOutputNameAllocated(0, allocator).get();

    //Ort::TypeInfo
    outputTypeInfo = EEm_session.GetOutputTypeInfo(0);
    auto EEm_outputTensorInfo = outputTypeInfo.GetTensorTypeAndShapeInfo();

    outputDims = EEm_outputTensorInfo.GetShape();

    size_t EE_TensorSize = nEETowersPad * nEETowersPad;
    std::vector<float> eemRecoOccMap1dPad(EE_TensorSize);  //To store the output reconstructed occupancy

    inputNames.push_back(inputName.c_str());
    outputNames.push_back(outputName.c_str());

    //Ort::MemoryInfo
    memoryInfo = Ort::MemoryInfo::CreateCpu(OrtAllocatorType::OrtArenaAllocator, OrtMemType::OrtMemTypeDefault);
    inputTensors.push_back(Ort::Value::CreateTensor<float>(
        memoryInfo, eemOccMap1dCumulPad.data(), EE_TensorSize, inputDims.data(), inputDims.size()));

    outputTensors.push_back(Ort::Value::CreateTensor<float>(
        memoryInfo, eemRecoOccMap1dPad.data(), EE_TensorSize, outputDims.data(), outputDims.size()));

    EEm_session.Run(Ort::RunOptions{nullptr},
                    inputNames.data(),
                    inputTensors.data(),
                    1,
                    outputNames.data(),
                    outputTensors.data(),
                    1);

    // EE+ //
    inputDims.clear();
    outputDims.clear();
    inputNames.clear();
    outputNames.clear();
    inputTensors.clear();
    outputTensors.clear();

    modelFilepath = edm::FileInPath("DQM/EcalMonitorClient/data/onnxModels/EEp_resnet2018.onnx").fullPath();

    Ort::Session EEp_session(env, modelFilepath.c_str(), sessionOptions);

    inputName = EEp_session.GetInputNameAllocated(0, allocator).get();

    inputTypeInfo = EEp_session.GetInputTypeInfo(0);
    auto EEp_inputTensorInfo = inputTypeInfo.GetTensorTypeAndShapeInfo();

    inputDims = EEp_inputTensorInfo.GetShape();

    outputName = EEp_session.GetOutputNameAllocated(0, allocator).get();

    outputTypeInfo = EEp_session.GetOutputTypeInfo(0);

    auto EEp_outputTensorInfo = outputTypeInfo.GetTensorTypeAndShapeInfo();

    outputDims = EEp_outputTensorInfo.GetShape();

    std::vector<float> eepRecoOccMap1dPad(EE_TensorSize);  //To store the output reconstructed occupancy

    inputNames.push_back(inputName.c_str());
    outputNames.push_back(outputName.c_str());

    //Ort::MemoryInfo
    memoryInfo = Ort::MemoryInfo::CreateCpu(OrtAllocatorType::OrtArenaAllocator, OrtMemType::OrtMemTypeDefault);
    inputTensors.push_back(Ort::Value::CreateTensor<float>(
        memoryInfo, eepOccMap1dCumulPad.data(), EE_TensorSize, inputDims.data(), inputDims.size()));

    outputTensors.push_back(Ort::Value::CreateTensor<float>(
        memoryInfo, eepRecoOccMap1dPad.data(), EE_TensorSize, outputDims.data(), outputDims.size()));

    EEp_session.Run(Ort::RunOptions{nullptr},
                    inputNames.data(),
                    inputTensors.data(),
                    1,
                    outputNames.data(),
                    outputTensors.data(),
                    1);

    //2D Loss map to store tower by tower loss between the output (reconstructed) and input occupancies,
    //Have same dimensions as the occupancy plot
    std::valarray<std::valarray<float>> EBlossMap2d(std::valarray<float>(nEBPhiTowers), nEBEtaTowers);
    std::valarray<std::valarray<float>> EEmlossMap2d(std::valarray<float>(nEEPhiTowers), nEEEtaTowers);
    std::valarray<std::valarray<float>> EEplossMap2d(std::valarray<float>(nEEPhiTowers), nEEEtaTowers);

    //1D val arrays to store row wise information corresponding to the reconstructed, input and average occupancies, and loss.
    //and to do element wise (tower wise) operations on them to calculate the MSE loss between the reco and input occupancy.
    std::valarray<float> EBrecoOcc1d(0., nEBPhiTowers);
    std::valarray<float> EBinputOcc1d(0., nEBPhiTowers);
    std::valarray<float> EBavgOcc1d(0., nEBPhiTowers);
    std::valarray<float> EBloss_;

    std::valarray<float> EEmrecoOcc1d(0., nEEPhiTowers);
    std::valarray<float> EEminputOcc1d(0., nEEPhiTowers);
    std::valarray<float> EEmavgOcc1d(0., nEEPhiTowers);
    std::valarray<float> EEmloss_;

    std::valarray<float> EEprecoOcc1d(0., nEEPhiTowers);
    std::valarray<float> EEpinputOcc1d(0., nEEPhiTowers);
    std::valarray<float> EEpavgOcc1d(0., nEEPhiTowers);
    std::valarray<float> EEploss_;

    std::string EBOccpath =
        edm::FileInPath("DQM/EcalMonitorClient/data/MLAvgOccupancy/EB_avgocc_Run2022_500ev.dat").fullPath();
    std::ifstream inFile;
    double val;
    inFile.open((EBOccpath).c_str());
    while (inFile) {
      inFile >> val;
      if (inFile.eof())
        break;
      EBavgOcc.push_back(val);
    }
    inFile.close();

    std::string EEmOccpath =
        edm::FileInPath("DQM/EcalMonitorClient/data/MLAvgOccupancy/EEm_avgocc_Run2022_500ev.dat").fullPath();
    inFile.open((EEmOccpath).c_str());
    while (inFile) {
      inFile >> val;
      if (inFile.eof())
        break;
      EEmavgOcc.push_back(val);
    }
    inFile.close();

    std::string EEpOccpath =
        edm::FileInPath("DQM/EcalMonitorClient/data/MLAvgOccupancy/EEp_avgocc_Run2022_500ev.dat").fullPath();
    inFile.open((EEpOccpath).c_str());
    while (inFile) {
      inFile >> val;
      if (inFile.eof())
        break;
      EEpavgOcc.push_back(val);
    }
    inFile.close();

    //Loss calculation
    //Ignore the top and bottom replicated padded rows when doing inference
    //by making index i run over (1,35) instead of (0,36) for EB, and over (1,21) for EE

    MESet const& sAEReco(sources_.at("AEReco"));
    TH2F* hEBRecoMap2d(sAEReco.getME(1)->getTH2F());

    for (int i = 1; i < nEBEtaTowersPad - 1; i++) {
      for (int j = 0; j < nEBPhiTowers; j++) {
        int k = (i * nEBPhiTowers) + j;
        int bin_ = hEBRecoMap2d->GetBin(j + 1, i);
        EBrecoOcc1d[j] = ebRecoOccMap1dPad[k];
        EBinputOcc1d[j] = ebOccMap1dCumulPad[k];
        EBavgOcc1d[j] = EBavgOcc[k];
        double content = ebRecoOccMap1dPad[k];
        hEBRecoMap2d->SetBinContent(bin_, content);
      }
      //Calculate the MSE loss = (output-input)^2, with avg response correction
      EBloss_ = std::pow((EBrecoOcc1d / EBavgOcc1d - EBinputOcc1d / EBavgOcc1d), 2);
      EBlossMap2d[i - 1] = (EBloss_);
    }

    TH2F* hEEmRecoMap2d(sAEReco.getME(0)->getTH2F());
    TH2F* hEEpRecoMap2d(sAEReco.getME(2)->getTH2F());

    for (int i = 1; i < nEETowersPad - 1; i++) {
      for (int j = 0; j < nEEPhiTowers; j++) {
        int k = (i * nEETowersPad) + j + 1;
        int bin_ = hEEmRecoMap2d->GetBin(j + 1, i);

        EEmrecoOcc1d[j] = eemRecoOccMap1dPad[k];
        EEminputOcc1d[j] = eemOccMap1dCumulPad[k];
        EEmavgOcc1d[j] = EEmavgOcc[k];
        double EEmcontent = eemRecoOccMap1dPad[k];
        hEEmRecoMap2d->SetBinContent(bin_, EEmcontent);

        EEprecoOcc1d[j] = eepRecoOccMap1dPad[k];
        EEpinputOcc1d[j] = eepOccMap1dCumulPad[k];
        EEpavgOcc1d[j] = EEpavgOcc[k];
        double EEpcontent = eepRecoOccMap1dPad[k];
        hEEpRecoMap2d->SetBinContent(bin_, EEpcontent);
      }
      //Calculate the MSE loss = (output-input)^2, with avg response correction
      EEmloss_ = std::pow((EEmrecoOcc1d / EEmavgOcc1d - EEminputOcc1d / EEmavgOcc1d), 2);
      EEmlossMap2d[i - 1] = (EEmloss_);

      EEploss_ = std::pow((EEprecoOcc1d / EEpavgOcc1d - EEpinputOcc1d / EEpavgOcc1d), 2);
      EEplossMap2d[i - 1] = (EEploss_);
    }

    //Store each loss map from the output in the queue
    EBlossMap2dQ.push_back(EBlossMap2d);
    EEmlossMap2dQ.push_back(EEmlossMap2d);
    EEplossMap2dQ.push_back(EEplossMap2d);

    //Keep exactly nLSloss loss maps to multiply
    if (EBlossMap2dQ.size() > nLSloss) {
      EBlossMap2dQ.pop_front();
      EEmlossMap2dQ.pop_front();
      EEplossMap2dQ.pop_front();
    }
    if (EBlossMap2dQ.size() < nLSloss) {  //Exit if there are not nLSloss loss maps
      return;
    }

    //To hold the final multiplied loss
    std::valarray<std::valarray<float>> EBlossMap2dMult(std::valarray<float>(1., nEBPhiTowers), nEBEtaTowers);
    std::valarray<std::valarray<float>> EEmlossMap2dMult(std::valarray<float>(1., nEEPhiTowers), nEEEtaTowers);
    std::valarray<std::valarray<float>> EEplossMap2dMult(std::valarray<float>(1., nEEPhiTowers), nEEEtaTowers);

    //Multiply together the last nLSloss loss maps
    //So that real anomalies which persist with time are enhanced and fluctuations are suppressed.
    for (size_t i = 0; i < EBlossMap2dQ.size(); i++) {
      EBlossMap2dMult *= EBlossMap2dQ[i];
      EEmlossMap2dMult *= EEmlossMap2dQ[i];
      EEplossMap2dMult *= EEplossMap2dQ[i];
    }

    //Fill the AELoss ME with the values of this time multiplied loss map
    //MESet const& sAELoss(sources_.at("AELoss"));
    MESet& sAELoss(sources_.at("AELoss"));

    TH2F* hEBLossMap2dMult(sAELoss.getME(1)->getTH2F());

    for (int i = 0; i < hEBLossMap2dMult->GetNbinsY(); i++) {
      for (int j = 0; j < hEBLossMap2dMult->GetNbinsX(); j++) {
        int bin_ = hEBLossMap2dMult->GetBin(j + 1, i + 1);
        double content = EBlossMap2dMult[i][j];
        hEBLossMap2dMult->SetBinContent(bin_, content);
      }
    }

    TH2F* hEEmLossMap2dMult(sAELoss.getME(0)->getTH2F());
    TH2F* hEEpLossMap2dMult(sAELoss.getME(2)->getTH2F());

    for (int i = 0; i < hEEmLossMap2dMult->GetNbinsY(); i++) {
      for (int j = 0; j < hEEmLossMap2dMult->GetNbinsX(); j++) {
        int bin_ = hEEmLossMap2dMult->GetBin(j + 1, i + 1);

        double EEmcontent = EEmlossMap2dMult[i][j];
        hEEmLossMap2dMult->SetBinContent(bin_, EEmcontent);

        double EEpcontent = EEplossMap2dMult[i][j];
        hEEpLossMap2dMult->SetBinContent(bin_, EEpcontent);
      }
    }

    //Apply the quality threshold on the time multiplied loss map stored in the ME AELoss
    //If anomalous, the tower entry will have a large loss value. If good, the value will be close to zero.

    MESet& meBadTowerCount(sources_.at("BadTowerCount"));
    MESet& meBadTowerCountNorm(sources_.at("BadTowerCountNorm"));
    MESet& meTrendMLBadTower(MEs_.at("TrendMLBadTower"));

    LScount++;

    MESet::const_iterator dAEnd(sAELoss.end(GetElectronicsMap()));
    for (MESet::const_iterator dItr(sAELoss.beginChannel(GetElectronicsMap())); dItr != dAEnd;
         dItr.toNextChannel(GetElectronicsMap())) {
      DetId id(dItr->getId());

      bool doMaskML(meMLQualitySummary.maskMatches(id, mask, statusManager_, GetTrigTowerMap()));

      float entries(dItr->getBinContent());

      int quality(doMaskML ? kMGood : kGood);
      float MLThreshold;

      if (id.subdetId() == EcalEndcap) {
        EEDetId eeid(id);
        if (eeid.zside() > 0)
          MLThreshold = EEpThreshold_;
        else
          MLThreshold = EEmThreshold_;
      } else {
        MLThreshold = EBThreshold_;
      }

      //If a trigger tower entry is greater than the ML threshold, set it to Bad quality, otherwise Good.
      if (entries > MLThreshold) {
        quality = doMaskML ? kMBad : kBad;
        meBadTowerCount.fill(getEcalDQMSetupObjects(), id);
        if (id.subdetId() == EcalEndcap)
          nbadtowerEE++;
        else
          nbadtowerEB++;
      }
      //Fill the quality summary with the quality of the given tower id.
      meMLQualitySummary.setBinContent(getEcalDQMSetupObjects(), id, double(quality));

      double badtowcount(meBadTowerCount.getBinContent(getEcalDQMSetupObjects(), id));
      meBadTowerCountNorm.setBinContent(getEcalDQMSetupObjects(), id, double(badtowcount / LScount));
    }  // ML Quality Summary

    meTrendMLBadTower.fill(getEcalDQMSetupObjects(), EcalBarrel, double(timestamp_.iLumi), double(nbadtowerEB));
    meTrendMLBadTower.fill(getEcalDQMSetupObjects(), EcalEndcap, double(timestamp_.iLumi), double(nbadtowerEE));

  }  // producePlots()

setParams()

void ecaldqm::MLClient::setParams ( edm::ParameterSet const &  _params )

overrideprivatevirtual

Reimplemented from ecaldqm::DQWorker.

Definition at line 21 of file MLClient.cc.

References l1extraParticles_cfi::_params, EB_PUcorr_intercept_, EB_PUcorr_slope_, EBThreshold_, EEm_PUcorr_intercept_, EEm_PUcorr_slope_, EEmThreshold_, EEp_PUcorr_intercept_, EEp_PUcorr_slope_, EEpThreshold_, ecaldqm::MESetCollection::erase(), ecaldqm::DQWorker::MEs_, ecaldqm::DQWorker::onlineMode_, ecaldqm::DQWorkerClient::sources_, and AlCaHLTBitMon_QueryRunRegistry::string.

                                                         {
    EBThreshold_ = _params.getUntrackedParameter<double>("EBThreshold");
    EEpThreshold_ = _params.getUntrackedParameter<double>("EEpThreshold");
    EEmThreshold_ = _params.getUntrackedParameter<double>("EEmThreshold");
    EB_PUcorr_slope_ = _params.getUntrackedParameter<double>("EB_PUcorr_slope");
    EB_PUcorr_intercept_ = _params.getUntrackedParameter<double>("EB_PUcorr_intercept");
    EEp_PUcorr_slope_ = _params.getUntrackedParameter<double>("EEp_PUcorr_slope");
    EEp_PUcorr_intercept_ = _params.getUntrackedParameter<double>("EEp_PUcorr_intercept");
    EEm_PUcorr_slope_ = _params.getUntrackedParameter<double>("EEm_PUcorr_slope");
    EEm_PUcorr_intercept_ = _params.getUntrackedParameter<double>("EEm_PUcorr_intercept");

    if (!onlineMode_) {
      MEs_.erase(std::string("MLQualitySummary"));
      MEs_.erase(std::string("EventsperMLImage"));
      sources_.erase(std::string("PU"));
      sources_.erase(std::string("NumEvents"));
      sources_.erase(std::string("DigiAllByLumi"));
      sources_.erase(std::string("AELoss"));
      sources_.erase(std::string("BadTowerCount"));
      sources_.erase(std::string("BadTowerCountNorm"));
    }
  }

Member Data Documentation

EB_PUcorr_intercept_

float ecaldqm::MLClient::EB_PUcorr_intercept_

private

Definition at line 39 of file MLClient.h.

Referenced by producePlots(), and setParams().

EB_PUcorr_slope_

float ecaldqm::MLClient::EB_PUcorr_slope_

private

Definition at line 36 of file MLClient.h.

Referenced by producePlots(), and setParams().

EBavgOcc

std::vector<double> ecaldqm::MLClient::EBavgOcc

private

Definition at line 56 of file MLClient.h.

Referenced by producePlots().

EBlossMap2dQ

std::deque<std::valarray<std::valarray<float> > > ecaldqm::MLClient::EBlossMap2dQ

private

Definition at line 60 of file MLClient.h.

Referenced by producePlots().

ebOccMap1dQ

std::deque<std::valarray<float> > ecaldqm::MLClient::ebOccMap1dQ

private

Definition at line 52 of file MLClient.h.

Referenced by producePlots().

EBThreshold_

float ecaldqm::MLClient::EBThreshold_

private

Definition at line 33 of file MLClient.h.

Referenced by producePlots(), and setParams().

EEm_PUcorr_intercept_

float ecaldqm::MLClient::EEm_PUcorr_intercept_

private

Definition at line 41 of file MLClient.h.

Referenced by producePlots(), and setParams().

EEm_PUcorr_slope_

float ecaldqm::MLClient::EEm_PUcorr_slope_

private

Definition at line 38 of file MLClient.h.

Referenced by producePlots(), and setParams().

EEmavgOcc

std::vector<double> ecaldqm::MLClient::EEmavgOcc

private

Definition at line 58 of file MLClient.h.

Referenced by producePlots().

EEmlossMap2dQ

std::deque<std::valarray<std::valarray<float> > > ecaldqm::MLClient::EEmlossMap2dQ

private

Definition at line 62 of file MLClient.h.

Referenced by producePlots().

eemOccMap1dQ

std::deque<std::valarray<float> > ecaldqm::MLClient::eemOccMap1dQ

private

Definition at line 54 of file MLClient.h.

Referenced by producePlots().

EEmThreshold_

float ecaldqm::MLClient::EEmThreshold_

private

Definition at line 35 of file MLClient.h.

Referenced by producePlots(), and setParams().

EEp_PUcorr_intercept_

float ecaldqm::MLClient::EEp_PUcorr_intercept_

private

Definition at line 40 of file MLClient.h.

Referenced by producePlots(), and setParams().

EEp_PUcorr_slope_

float ecaldqm::MLClient::EEp_PUcorr_slope_

private

Definition at line 37 of file MLClient.h.

Referenced by producePlots(), and setParams().

EEpavgOcc

std::vector<double> ecaldqm::MLClient::EEpavgOcc

private

Definition at line 57 of file MLClient.h.

Referenced by producePlots().

EEplossMap2dQ

std::deque<std::valarray<std::valarray<float> > > ecaldqm::MLClient::EEplossMap2dQ

private

Definition at line 61 of file MLClient.h.

Referenced by producePlots().

eepOccMap1dQ

std::deque<std::valarray<float> > ecaldqm::MLClient::eepOccMap1dQ

private

Definition at line 53 of file MLClient.h.

Referenced by producePlots().

EEpThreshold_

float ecaldqm::MLClient::EEpThreshold_

private

Definition at line 34 of file MLClient.h.

Referenced by producePlots(), and setParams().

LScount

int ecaldqm::MLClient::LScount = 0

private

Definition at line 47 of file MLClient.h.

Referenced by producePlots().

nbadtowerEB

int ecaldqm::MLClient::nbadtowerEB

private

Definition at line 45 of file MLClient.h.

Referenced by producePlots().

nbadtowerEE

int ecaldqm::MLClient::nbadtowerEE

private

Definition at line 46 of file MLClient.h.

Referenced by producePlots().

nEBEtaTowers

const int ecaldqm::MLClient::nEBEtaTowers = 34

staticprivate

Definition at line 23 of file MLClient.h.

Referenced by producePlots().

nEBEtaTowersPad

const int ecaldqm::MLClient::nEBEtaTowersPad = 36

staticprivate

Definition at line 31 of file MLClient.h.

Referenced by producePlots().

nEBPhiTowers

const int ecaldqm::MLClient::nEBPhiTowers = 72

staticprivate

Definition at line 24 of file MLClient.h.

Referenced by producePlots().

nEEEtaTowers

const int ecaldqm::MLClient::nEEEtaTowers = 20

staticprivate

Definition at line 26 of file MLClient.h.

Referenced by producePlots().

nEEPhiTowers

const int ecaldqm::MLClient::nEEPhiTowers = 20

staticprivate

Definition at line 27 of file MLClient.h.

Referenced by producePlots().

nEETowersPad

const int ecaldqm::MLClient::nEETowersPad = 22

staticprivate

Definition at line 32 of file MLClient.h.

Referenced by producePlots().

NEventQ

std::deque<int> ecaldqm::MLClient::NEventQ

private

Definition at line 49 of file MLClient.h.

Referenced by producePlots().

nLS

size_t ecaldqm::MLClient::nLS = 4

private

Definition at line 43 of file MLClient.h.

Referenced by producePlots().

nLSloss

size_t ecaldqm::MLClient::nLSloss = 6

private

Definition at line 44 of file MLClient.h.

Referenced by producePlots().