diff --git a/Geometry/HcalCommonData/interface/HcalDDDRecConstants.h b/Geometry/HcalCommonData/interface/HcalDDDRecConstants.h index 1cd7a60e3b103..070152daa23aa 100644 --- a/Geometry/HcalCommonData/interface/HcalDDDRecConstants.h +++ b/Geometry/HcalCommonData/interface/HcalDDDRecConstants.h @@ -61,6 +61,7 @@ class HcalDDDRecConstants { return gcons; } } + int findDepth(const int& det, const int& eta, const int& phi, const int& zside, const int& lay) const; std::vector getDepth(const int& det, const int& phi, const int& zside, const unsigned int& eta) const; std::vector getDepth(const unsigned int& eta, const bool& extra) const; int getDepthEta16(const int& det, const int& iphi, const int& zside) const { diff --git a/Geometry/HcalCommonData/interface/HcalHitRelabeller.h b/Geometry/HcalCommonData/interface/HcalHitRelabeller.h index f4a166dafd891..c13346197a2e2 100644 --- a/Geometry/HcalCommonData/interface/HcalHitRelabeller.h +++ b/Geometry/HcalCommonData/interface/HcalHitRelabeller.h @@ -14,10 +14,9 @@ class HcalHitRelabeller { void setGeometry(const HcalDDDRecConstants*&); DetId relabel(const uint32_t testId) const; static DetId relabel(const uint32_t testId, const HcalDDDRecConstants* theRecNumber); - -private: double energyWt(const uint32_t testId) const; +private: const HcalDDDRecConstants* theRecNumber; bool neutralDensity_; }; diff --git a/Geometry/HcalCommonData/src/HcalDDDRecConstants.cc b/Geometry/HcalCommonData/src/HcalDDDRecConstants.cc index 773b325b54fce..4ee9d2c6b0d0c 100644 --- a/Geometry/HcalCommonData/src/HcalDDDRecConstants.cc +++ b/Geometry/HcalCommonData/src/HcalDDDRecConstants.cc @@ -7,6 +7,7 @@ #include "CLHEP/Units/GlobalSystemOfUnits.h" #include #include +#include //#define EDM_ML_DEBUG using namespace geant_units::operators; @@ -27,6 +28,24 @@ HcalDDDRecConstants::~HcalDDDRecConstants() { #endif } +int HcalDDDRecConstants::findDepth( + const int& det, const int& eta, const int& phi, const int& zside, const int& lay) const { + int depth = hcons.findDepth(det, eta, phi, zside, lay); + if (depth < 0) { + std::vector depths = getDepth(eta, false); + if ((lay > 0) && (lay <= static_cast(depths.size()))) + depth = depths[lay - 1]; +#ifdef EDM_ML_DEBUG + std::ostringstream st1; + st1 << depths.size() << " depths "; + for (const auto& d : depths) + st1 << ": " << d; + edm::LogVerbatim("HCalGeom") << "HcalDDDRecConstants:: " << st1.str() << " for eta = " << eta << " Depth " << depth; +#endif + } + return depth; +} + std::vector HcalDDDRecConstants::getDepth(const unsigned int& eta, const bool& extra) const { if (!extra) { std::vector::const_iterator last = hpar->layerGroupEtaRec.begin(); diff --git a/Geometry/HcalCommonData/src/HcalHitRelabeller.cc b/Geometry/HcalCommonData/src/HcalHitRelabeller.cc index 34929db5fcaff..41fdc015925ed 100644 --- a/Geometry/HcalCommonData/src/HcalHitRelabeller.cc +++ b/Geometry/HcalCommonData/src/HcalHitRelabeller.cc @@ -89,7 +89,9 @@ double HcalHitRelabeller::energyWt(const uint32_t testId) const { int det, z, depth, eta, phi, layer; HcalTestNumbering::unpackHcalIndex(testId, det, z, depth, eta, phi, layer); int zside = (z == 0) ? (-1) : (1); - double wt = (((det == 1) || (det == 2)) && (depth == 1)) ? theRecNumber->getLayer0Wt(det, phi, zside) : 1.0; + double wt = ((((det == 1) && (layer <= 1)) || ((det == 2) && (layer <= 2))) && (depth == 1)) + ? theRecNumber->getLayer0Wt(det, phi, zside) + : 1.0; #ifdef EDM_ML_DEBUG edm::LogVerbatim("HcalSim") << "EnergyWT::det: " << det << " z: " << z << ":" << zside << " depth: " << depth << " ieta: " << eta << " iphi: " << phi << " layer: " << layer << " wt " << wt; diff --git a/Geometry/HcalCommonData/src/HcalLayerDepthMap.cc b/Geometry/HcalCommonData/src/HcalLayerDepthMap.cc index f6d6590b62ce7..79fb1694ff5e3 100644 --- a/Geometry/HcalCommonData/src/HcalLayerDepthMap.cc +++ b/Geometry/HcalCommonData/src/HcalLayerDepthMap.cc @@ -254,7 +254,7 @@ int HcalLayerDepthMap::getMaxDepthLastHE(const int subdet, const int iphi, const double HcalLayerDepthMap::getLayer0Wt(const int subdet, const int iphi, const int zside) const { double wt = isValid(subdet, iphi, zside) ? wtl0C_ : -1.0; #ifdef EDM_ML_DEBUG - edm::LogVerbatim("HCalGeom") << "Debug info -- getLayer0Wt::Input " << subdet << ":" << iphi << ":" << zside + edm::LogVerbatim("HCalGeom") << "HcalLayerDepthMap -- getLayer0Wt::Input " << subdet << ":" << iphi << ":" << zside << " Output " << wt; #endif return wt; diff --git a/Geometry/HcalCommonData/test/HcalHitRelabellerTester.cc b/Geometry/HcalCommonData/test/HcalHitRelabellerTester.cc new file mode 100644 index 0000000000000..3024205714485 --- /dev/null +++ b/Geometry/HcalCommonData/test/HcalHitRelabellerTester.cc @@ -0,0 +1,93 @@ +// -*- C++ -*- +// +// Package: HcalHitRelabellerTester +// Class: HcalHitRelabellerTester +// +/**\class HcalHitRelabellerTester HcalHitRelabellerTester.cc test/HcalHitRelabellerTester.cc + + Description: + + Implementation: + +*/ +// +// Original Author: Sunanda Banerjee +// Created: Mon 2023/06/24 +// + +// system include files +#include +#include +#include + +// user include files +#include "FWCore/Framework/interface/Frameworkfwd.h" +#include "FWCore/Framework/interface/one/EDAnalyzer.h" + +#include "FWCore/Framework/interface/Event.h" +#include "FWCore/Framework/interface/EventSetup.h" +#include "FWCore/Framework/interface/MakerMacros.h" +#include "FWCore/MessageLogger/interface/MessageLogger.h" +#include "FWCore/ParameterSet/interface/ParameterSet.h" +#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" +#include "DataFormats/HcalDetId/interface/HcalTestNumbering.h" +#include "Geometry/Records/interface/IdealGeometryRecord.h" +#include "Geometry/Records/interface/HcalRecNumberingRecord.h" +#include "Geometry/HcalCommonData/interface/HcalDDDRecConstants.h" +#include "Geometry/HcalCommonData/interface/HcalHitRelabeller.h" + +class HcalHitRelabellerTester : public edm::one::EDAnalyzer<> { +public: + explicit HcalHitRelabellerTester(const edm::ParameterSet&); + ~HcalHitRelabellerTester() override = default; + static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); + + void analyze(edm::Event const& iEvent, edm::EventSetup const&) override; + +private: + bool nd_; + edm::ESGetToken token_; +}; + +HcalHitRelabellerTester::HcalHitRelabellerTester(const edm::ParameterSet& ps) + : nd_(ps.getParameter("neutralDensity")), + token_{esConsumes(edm::ESInputTag{})} { + edm::LogVerbatim("HCalGeom") << "Construct HcalHitRelabellerTester with Neutral Density: " << nd_; +} + +void HcalHitRelabellerTester::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + desc.add("neutralDensity", true); + descriptions.add("hcalHitRelabellerTester", desc); +} + +// ------------ method called to produce the data ------------ +void HcalHitRelabellerTester::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) { + const HcalDDDRecConstants* theRecNumber = &iSetup.getData(token_); + std::unique_ptr relabel = std::make_unique(nd_); + relabel->setGeometry(theRecNumber); + std::vector etas = {-29, -26, -22, -19, -16, -13, -10, -7, -4, -1, 2, 5, 8, 11, 14, 17, 20, 23, 26, 29}; + std::vector layers = {1, 2}; + const int iphi = 63; + edm::LogVerbatim("HCalGeom") << "HcalHitRelabellerTester:: Testing " << etas.size() << " eta, " << layers.size() + << " layer values and iphi = " << iphi; + for (const auto& eta : etas) { + int ieta = std::abs(eta); + int det = (ieta <= 16) ? 1 : 2; + int zside = (eta >= 0) ? 1 : -1; + for (const auto& lay : layers) { + if (ieta == 16) + det = (lay <= 3) ? 1 : 2; + int depth = theRecNumber->findDepth(det, ieta, iphi, zside, lay); + if (depth > 0) { + uint32_t id = HcalTestNumbering::packHcalIndex(det, zside, depth, ieta, iphi, lay); + double wt = relabel->energyWt(id); + edm::LogVerbatim("HCalGeom") << "Det " << det << " Z " << zside << " Eta" << ieta << " Phi " << iphi << " Lay " + << lay << " Depth " << depth << " Layer0Wt " << wt; + } + } + } +} + +//define this as a plug-in +DEFINE_FWK_MODULE(HcalHitRelabellerTester); diff --git a/Geometry/HcalCommonData/test/python/testHitRelabeller_cfg.py b/Geometry/HcalCommonData/test/python/testHitRelabeller_cfg.py new file mode 100644 index 0000000000000..b5a09ae09ca82 --- /dev/null +++ b/Geometry/HcalCommonData/test/python/testHitRelabeller_cfg.py @@ -0,0 +1,209 @@ +############################################################################### +# Way to use this: +# cmsRun testHitRelabeller_cfg.py geometry=2021 +# +# Options for geometry 2016, 2017, 2018, 2021, 2023 +# +############################################################################### +import FWCore.ParameterSet.Config as cms +import os, sys, imp, re +import FWCore.ParameterSet.VarParsing as VarParsing + +#################################################################### +### SETUP OPTIONS +options = VarParsing.VarParsing('standard') +options.register('geometry', + "2021", + VarParsing.VarParsing.multiplicity.singleton, + VarParsing.VarParsing.varType.string, + "geometry of operations: 2016, 2017, 2018, 2021, 2023") + +### get and parse the command line arguments +options.parseArguments() + +print(options) + +#################################################################### +# Use the options + +geomName = "Configuration.Geometry.GeometryExtended" + options.geometry + "Reco_cff" + +if (options.geometry == "2016"): + from Configuration.Eras.Era_Run2_2016_cff import Run2_2016 + process = cms.Process('HitRelabeller',Run2_2016) + baseName = "auto:run2_data" +elif (options.geometry == "2017"): + from Configuration.Eras.Era_Run2_2017_cff import Run2_2017 + process = cms.Process('HitRelabeller',Run2_2017) + baseName = "auto:run2_data" +elif (options.geometry == "2018"): + from Configuration.Eras.Era_Run2_2018_cff import Run2_2018 + process = cms.Process('HitRelabeller',Run2_2018) + baseName = "auto:run2_data" +else: + from Configuration.Eras.Era_Run3_DDD_cff import Run3_DDD + process = cms.Process('HitRelabeller',Run3_DDD) + baseName = "auto:run3_data" + +print("Base file Name: ", baseName) +print("Geom file Name: ", geomName) + + +# import of standard configurations +process.load('Configuration.StandardSequences.Services_cff') +process.load('SimGeneral.HepPDTESSource.pythiapdt_cfi') +process.load('FWCore.MessageService.MessageLogger_cfi') +process.load('Configuration.EventContent.EventContent_cff') +process.load('SimGeneral.MixingModule.mixNoPU_cfi') +process.load('Configuration.StandardSequences.MagneticField_cff') +process.load('Configuration.StandardSequences.Generator_cff') +process.load('IOMC.EventVertexGenerators.VtxSmearedRealistic25ns13p6TeVEarly2022Collision_cfi') +process.load('GeneratorInterface.Core.genFilterSummary_cff') +process.load('Configuration.StandardSequences.SimIdeal_cff') +process.load('Configuration.StandardSequences.EndOfProcess_cff') +process.load('Configuration.StandardSequences.FrontierConditions_GlobalTag_cff') +process.load(geomName) + +if hasattr(process,'MessageLogger'): + process.MessageLogger.HCalGeom=dict() + +process.maxEvents = cms.untracked.PSet( + input = cms.untracked.int32(1), + output = cms.optional.untracked.allowed(cms.int32,cms.PSet) +) + +# Input source +process.source = cms.Source("EmptySource") + +process.options = cms.untracked.PSet( + FailPath = cms.untracked.vstring(), + IgnoreCompletely = cms.untracked.vstring(), + Rethrow = cms.untracked.vstring(), + SkipEvent = cms.untracked.vstring(), + allowUnscheduled = cms.obsolete.untracked.bool, + canDeleteEarly = cms.untracked.vstring(), + emptyRunLumiMode = cms.obsolete.untracked.string, + eventSetup = cms.untracked.PSet( + forceNumberOfConcurrentIOVs = cms.untracked.PSet( + + ), + numberOfConcurrentIOVs = cms.untracked.uint32(1) + ), + fileMode = cms.untracked.string('FULLMERGE'), + forceEventSetupCacheClearOnNewRun = cms.untracked.bool(False), + makeTriggerResults = cms.obsolete.untracked.bool, + numberOfConcurrentLuminosityBlocks = cms.untracked.uint32(1), + numberOfConcurrentRuns = cms.untracked.uint32(1), + numberOfStreams = cms.untracked.uint32(0), + numberOfThreads = cms.untracked.uint32(1), + printDependencies = cms.untracked.bool(False), + sizeOfStackForThreadsInKB = cms.optional.untracked.uint32, + throwIfIllegalParameter = cms.untracked.bool(True), + wantSummary = cms.untracked.bool(False) +) + +# Production Info +process.configurationMetadata = cms.untracked.PSet( + annotation = cms.untracked.string('ZMM_14TeV_TuneCP5_cfi nevts:10'), + name = cms.untracked.string('Applications'), + version = cms.untracked.string('$Revision: 1.19 $') +) + +# Other statements +process.genstepfilter.triggerConditions=cms.vstring("generation_step") +from Configuration.AlCa.GlobalTag import GlobalTag +process.GlobalTag = GlobalTag(process.GlobalTag, baseName, '') + +process.generator = cms.EDFilter("Pythia8GeneratorFilter", + pythiaHepMCVerbosity = cms.untracked.bool(False), + maxEventsToPrint = cms.untracked.int32(0), + pythiaPylistVerbosity = cms.untracked.int32(0), + filterEfficiency = cms.untracked.double(1.0), + comEnergy = cms.double(14000.0), + PythiaParameters = cms.PSet( + pythia8CommonSettings = cms.vstring( + 'Tune:preferLHAPDF = 2', + 'Main:timesAllowErrors = 10000', + 'Check:epTolErr = 0.01', + 'Beams:setProductionScalesFromLHEF = off', + 'SLHA:minMassSM = 1000.', + 'ParticleDecays:limitTau0 = on', + 'ParticleDecays:tau0Max = 10', + 'ParticleDecays:allowPhotonRadiation = on', + ), + pythia8CP5Settings = cms.vstring( + 'Tune:pp 14', + 'Tune:ee 7', + 'MultipartonInteractions:ecmPow=0.03344', + 'MultipartonInteractions:bProfile=2', + 'MultipartonInteractions:pT0Ref=1.41', + 'MultipartonInteractions:coreRadius=0.7634', + 'MultipartonInteractions:coreFraction=0.63', + 'ColourReconnection:range=5.176', + 'SigmaTotal:zeroAXB=off', + 'SpaceShower:alphaSorder=2', + 'SpaceShower:alphaSvalue=0.118', + 'SigmaProcess:alphaSvalue=0.118', + 'SigmaProcess:alphaSorder=2', + 'MultipartonInteractions:alphaSvalue=0.118', + 'MultipartonInteractions:alphaSorder=2', + 'TimeShower:alphaSorder=2', + 'TimeShower:alphaSvalue=0.118', + 'SigmaTotal:mode = 0', + 'SigmaTotal:sigmaEl = 21.89', + 'SigmaTotal:sigmaTot = 100.309', + 'PDF:pSet=LHAPDF6:NNPDF31_nnlo_as_0118', + ), + processParameters = cms.vstring( + 'WeakSingleBoson:ffbar2gmZ = on', + '23:onMode = off', + '23:onIfAny = 13', + 'PhaseSpace:mHatMin = 75.', + ), + parameterSets = cms.vstring('pythia8CommonSettings', + 'pythia8CP5Settings', + 'processParameters', + ) + ) +) + +process.mumugenfilter = cms.EDFilter("MCParticlePairFilter", + MaxEta = cms.untracked.vdouble(4.0, 4.0), + MinEta = cms.untracked.vdouble(-4.0, -4.0), + MinPt = cms.untracked.vdouble(2.5, 2.5), + ParticleCharge = cms.untracked.int32(-1), + ParticleID1 = cms.untracked.vint32(13), + ParticleID2 = cms.untracked.vint32(13), + Status = cms.untracked.vint32(1, 1) +) + + +process.ProductionFilterSequence = cms.Sequence(process.generator+process.mumugenfilter) +process.load("Geometry.HcalCommonData.hcalHitRelabellerTester_cfi") + +# Path and EndPath definitions +process.generation_step = cms.Path(process.pgen) +process.simulation_step = cms.Path(process.psim) +process.genfiltersummary_step = cms.EndPath(process.genFilterSummary) +process.endjob_step = cms.EndPath(process.endOfProcess) +process.analysis_step = cms.EndPath(process.hcalHitRelabellerTester) + +# Schedule definition +process.schedule = cms.Schedule(process.generation_step, + process.genfiltersummary_step, + process.simulation_step, + process.endjob_step, + process.analysis_step) +from PhysicsTools.PatAlgos.tools.helpers import associatePatAlgosToolsTask +associatePatAlgosToolsTask(process) +# filter all path with the production filter sequence +for path in process.paths: + getattr(process,path).insert(0, process.ProductionFilterSequence) + + +# Customisation from command line + +# Add early deletion of temporary data products to reduce peak memory need +from Configuration.StandardSequences.earlyDeleteSettings_cff import customiseEarlyDelete +process = customiseEarlyDelete(process) +# End adding early deletion