From c00ea7ba4d00dead2ccf9d6c83feee85eb957644 Mon Sep 17 00:00:00 2001 From: Fedor Date: Wed, 8 Oct 2014 00:46:52 +0200 Subject: [PATCH] introduce different scenarios for HE darkening --- .../HcalCalibObjects/interface/HEDarkening.h | 24 +++++-- .../HcalCalibObjects/src/HEDarkening.cc | 69 ++++++++++++------- .../Configuration/python/aging.py | 4 +- .../python/hcalUnsuppressedDigis_cfi.py | 2 +- .../HcalSimProducers/src/HcalDigitizer.cc | 4 +- SimG4CMS/Calo/src/HCalSD.cc | 5 +- SimG4Core/Application/python/g4SimHits_cfi.py | 2 +- 7 files changed, 69 insertions(+), 41 deletions(-) diff --git a/DataFormats/HcalCalibObjects/interface/HEDarkening.h b/DataFormats/HcalCalibObjects/interface/HEDarkening.h index 5e095ad4722d4..7d37e03ed9a71 100644 --- a/DataFormats/HcalCalibObjects/interface/HEDarkening.h +++ b/DataFormats/HcalCalibObjects/interface/HEDarkening.h @@ -4,25 +4,35 @@ // Simple class with parameterizaed function to get darkening attenuation // coefficiant for SLHC conditions // = degradation(int_lumi(intlumi) * dose(layer,Radius)), where -// intlumi is integrated luminosity (fb-1), +// intlumi is integrated luminosity (fb-1), // layer is HE layer number (from -1 up// to 17), NB: 1-19 in HcalTestNumbering // Radius is radius from the beam line (cm) // - -#define maxEta 14 -#define maxLay 19 +// scenario = 0 - full replacement of HE scintillators +// scenario = 1 - no replacement, full stage darkening +// scenario = 2 - only replace scintillators with expected light yield < 20% for 500 fb-1 +// scenario = 3 - replace complete megatiles only in the front 4 layers +// class HEDarkening { public: - HEDarkening(); + HEDarkening(unsigned int scenario = 3); ~HEDarkening(); - float degradation(float intlumi, int ieta, int lay); + float degradation(float intlumi, int ieta, int lay) const; + + static const char* scenarioDescription (unsigned int scenario); + + // 2 contsants below are used in CalibCalorimetry/HcalPlugins HERecalibration + // (1) number of HE ieta bins for darkening + static const unsigned int nEtaBins = 14; + // (2) max. number of HE scint. layers + static const unsigned int nScintLayers = 19; private: int ieta_shift; - float lumiscale[maxEta][maxLay]; + float lumiscale[nEtaBins][nScintLayers]; }; diff --git a/DataFormats/HcalCalibObjects/src/HEDarkening.cc b/DataFormats/HcalCalibObjects/src/HEDarkening.cc index 41dbd34e4ce1b..be617fc4e8260 100644 --- a/DataFormats/HcalCalibObjects/src/HEDarkening.cc +++ b/DataFormats/HcalCalibObjects/src/HEDarkening.cc @@ -10,7 +10,7 @@ //#define DebugLog -HEDarkening::HEDarkening() { +HEDarkening::HEDarkening(unsigned int scenario) { //HE starts at tower 16 ieta_shift = 16; @@ -18,38 +18,44 @@ HEDarkening::HEDarkening() { //fits from laser data for L1 and L7 //L0,Lm1 = L1, L2-L6 interpolated, L8-L17 extrapolated //from Vladimir Epshteyn - float _lumiscale[maxEta][maxLay] = { - {1194.9, 1194.9, 1194.9, 1651.5, 2282.7, 3155.2, 4361.0, 6027.8, 8331.5, 11515.7, 15916.8, 22000.0, 30408.2, 42029.8, 58093.1, 80295.6, 110983.5, 153400.1, 212027.7}, - {952.8, 952.8, 952.8, 1293.9, 1757.1, 2386.1, 3240.3, 4400.3, 5975.4, 8114.5, 11019.3, 14963.9, 20320.6, 27594.9, 37473.2, 50887.7, 69104.3, 93841.9, 127435.0}, - {759.8, 759.8, 759.8, 1013.8, 1352.5, 1804.5, 2407.6, 3212.2, 4285.7, 5717.9, 7628.7, 10178.1, 13579.5, 18117.6, 24172.3, 32250.4, 43028.0, 57407.4, 76592.2}, - {605.9, 605.9, 605.9, 794.2, 1041.1, 1364.7, 1788.9, 2344.9, 3073.7, 4029.1, 5281.4, 6922.9, 9074.7, 11895.2, 15592.4, 20438.8, 26791.5, 35118.8, 46034.2}, - {483.2, 483.2, 483.2, 622.3, 801.4, 1032.1, 1329.2, 1711.8, 2204.5, 2839.1, 3656.3, 4708.8, 6064.3, 7809.9, 10058.0, 12953.2, 16681.8, 21483.8, 27667.9}, - {385.3, 385.3, 385.3, 487.5, 616.9, 780.5, 987.6, 1249.6, 1581.1, 2000.6, 2531.3, 3202.8, 4052.5, 5127.6, 6487.9, 8209.2, 10387.0, 13142.6, 16629.3}, - {307.3, 307.3, 307.3, 382.0, 474.8, 590.3, 733.8, 912.2, 1134.0, 1409.7, 1752.4, 2178.5, 2708.1, 3366.6, 4185.1, 5202.6, 6467.5, 8039.9, 9994.7}, - {245.0, 245.0, 245.0, 299.3, 365.5, 446.4, 545.2, 665.9, 813.3, 993.3, 1213.2, 1481.8, 1809.7, 2210.3, 2699.6, 3297.2, 4027.0, 4918.4, 6007.1}, - {195.4, 195.4, 195.4, 234.5, 281.3, 337.6, 405.1, 486.1, 583.3, 700.0, 839.9, 1007.9, 1209.4, 1451.2, 1741.4, 2089.6, 2507.4, 3008.8, 3610.5}, - {155.8, 155.8, 155.8, 183.7, 216.6, 255.3, 301.0, 354.9, 418.4, 493.2, 581.5, 685.5, 808.2, 952.8, 1123.3, 1324.3, 1561.3, 1840.6, 2170.0}, - {124.3, 124.3, 124.3, 143.9, 166.7, 193.1, 223.6, 259.0, 300.1, 347.5, 402.6, 466.3, 540.1, 625.6, 724.6, 839.3, 972.1, 1126.0, 1304.2}, - {99.1, 99.1, 99.1, 112.8, 128.3, 146.0, 166.2, 189.1, 215.2, 244.9, 278.7, 317.2, 360.9, 410.7, 467.4, 531.9, 605.3, 688.8, 783.9}, - {79.0, 79.0, 79.0, 88.3, 98.8, 110.4, 123.5, 138.0, 154.3, 172.6, 192.9, 215.7, 241.2, 269.7, 301.5, 337.1, 376.9, 421.4, 471.1}, - {63.0, 63.0, 63.0, 69.2, 76.0, 83.5, 91.7, 100.8, 110.7, 121.6, 133.6, 146.7, 161.2, 177.0, 194.5, 213.6, 234.7, 257.8, 283.2} - }; + const float _lumiscale[nEtaBins][nScintLayers] = { + {1194.9, 1194.9, 1194.9, 1651.5, 2282.7, 3155.2, 4361.0, 6027.8, 8331.5, 11515.7, 15916.8, 22000.0, 30408.2, 42029.8, 58093.1, 80295.6, 110983.5, 153400.1, 212027.7}, + {952.8, 952.8, 952.8, 1293.9, 1757.1, 2386.1, 3240.3, 4400.3, 5975.4, 8114.5, 11019.3, 14963.9, 20320.6, 27594.9, 37473.2, 50887.7, 69104.3, 93841.9, 127435.0}, + {759.8, 759.8, 759.8, 1013.8, 1352.5, 1804.5, 2407.6, 3212.2, 4285.7, 5717.9, 7628.7, 10178.1, 13579.5, 18117.6, 24172.3, 32250.4, 43028.0, 57407.4, 76592.2}, + {605.9, 605.9, 605.9, 794.2, 1041.1, 1364.7, 1788.9, 2344.9, 3073.7, 4029.1, 5281.4, 6922.9, 9074.7, 11895.2, 15592.4, 20438.8, 26791.5, 35118.8, 46034.2}, + {483.2, 483.2, 483.2, 622.3, 801.4, 1032.1, 1329.2, 1711.8, 2204.5, 2839.1, 3656.3, 4708.8, 6064.3, 7809.9, 10058.0, 12953.2, 16681.8, 21483.8, 27667.9}, + {385.3, 385.3, 385.3, 487.5, 616.9, 780.5, 987.6, 1249.6, 1581.1, 2000.6, 2531.3, 3202.8, 4052.5, 5127.6, 6487.9, 8209.2, 10387.0, 13142.6, 16629.3}, + {307.3, 307.3, 307.3, 382.0, 474.8, 590.3, 733.8, 912.2, 1134.0, 1409.7, 1752.4, 2178.5, 2708.1, 3366.6, 4185.1, 5202.6, 6467.5, 8039.9, 9994.7}, + {245.0, 245.0, 245.0, 299.3, 365.5, 446.4, 545.2, 665.9, 813.3, 993.3, 1213.2, 1481.8, 1809.7, 2210.3, 2699.6, 3297.2, 4027.0, 4918.4, 6007.1}, + {195.4, 195.4, 195.4, 234.5, 281.3, 337.6, 405.1, 486.1, 583.3, 700.0, 839.9, 1007.9, 1209.4, 1451.2, 1741.4, 2089.6, 2507.4, 3008.8, 3610.5}, + {155.8, 155.8, 155.8, 183.7, 216.6, 255.3, 301.0, 354.9, 418.4, 493.2, 581.5, 685.5, 808.2, 952.8, 1123.3, 1324.3, 1561.3, 1840.6, 2170.0}, + {124.3, 124.3, 124.3, 143.9, 166.7, 193.1, 223.6, 259.0, 300.1, 347.5, 402.6, 466.3, 540.1, 625.6, 724.6, 839.3, 972.1, 1126.0, 1304.2}, + {99.1, 99.1, 99.1, 112.8, 128.3, 146.0, 166.2, 189.1, 215.2, 244.9, 278.7, 317.2, 360.9, 410.7, 467.4, 531.9, 605.3, 688.8, 783.9}, + {79.0, 79.0, 79.0, 88.3, 98.8, 110.4, 123.5, 138.0, 154.3, 172.6, 192.9, 215.7, 241.2, 269.7, 301.5, 337.1, 376.9, 421.4, 471.1}, + {63.0, 63.0, 63.0, 69.2, 76.0, 83.5, 91.7, 100.8, 110.7, 121.6, 133.6, 146.7, 161.2, 177.0, 194.5, 213.6, 234.7, 257.8, 283.2} + }; //store array //flux_factor: to account for increased flux at 14 TeV vs 7-8 TeV (approximate) //*divide* lumiscale params by this since increased flux -> faster darkening double flux_factor = 1.2; - for(int j = 0; j < maxEta; j++){ - for(int i = 0; i < maxLay; i++){ - lumiscale[j][i] = _lumiscale[j][i]/flux_factor; - } + for(unsigned int j = 0; j < nEtaBins; j++){ + for(unsigned int i = 0; i < nScintLayers; i++){ + if (scenario == 0 || + (scenario == 2 && exp(-500./_lumiscale[j][i])<0.2) || + (scenario == 3 && i<4)) { + lumiscale[j][i] = 0; + } + else { + lumiscale[j][i] = _lumiscale[j][i]/flux_factor; + } + } } - } HEDarkening::~HEDarkening() { } -float HEDarkening::degradation(float intlumi, int ieta, int lay) { +float HEDarkening::degradation(float intlumi, int ieta, int lay) const { //no lumi, no darkening if(intlumi <= 0) return 1.; @@ -57,14 +63,25 @@ float HEDarkening::degradation(float intlumi, int ieta, int lay) { ieta -= ieta_shift; //if outside eta range, no darkening - if(ieta < 0 || ieta >= maxEta) return 1.; + if(ieta < 0 || ieta >= (int)nEtaBins) return 1.; //shift layer index by 1 to act as array index lay += 1; //if outside layer range, no darkening - if(lay < 0 || lay >= maxLay) return 1.; + if(lay < 0 || lay >= (int)nScintLayers) return 1.; + + // if replaced tile by scenario + if (lumiscale[ieta][lay] == 0) return 1.; //return darkening factor return (exp(-intlumi/lumiscale[ieta][lay])); -} \ No newline at end of file +} + +const char* HEDarkening::scenarioDescription (unsigned int scenario) { + if (scenario == 0) return "full replacement of HE scintillators, no darkening"; + else if (scenario == 1) return "no replacement, full stage darkening"; + else if (scenario == 2) return "only replace scintillators with expected light yield < 20% for 500 fb-1"; + else if (scenario == 3) return "replace complete megatiles only in the front 4 layers"; + return "undefined scenario: assume no replacement, full stage darkening"; +} diff --git a/SLHCUpgradeSimulations/Configuration/python/aging.py b/SLHCUpgradeSimulations/Configuration/python/aging.py index c1429daae5316..4b75a14e11eb1 100644 --- a/SLHCUpgradeSimulations/Configuration/python/aging.py +++ b/SLHCUpgradeSimulations/Configuration/python/aging.py @@ -74,7 +74,7 @@ def ageHcal(process,lumi): if hasattr(process,'mix') and hasattr(process.mix,'digitizers') and hasattr(process.mix.digitizers,'hcal'): process.mix.digitizers.hcal.DelivLuminosity = cms.double(float(lumi)) # integrated lumi in fb-1 - process.mix.digitizers.hcal.HEDarkening = cms.bool(True) + process.mix.digitizers.hcal.HEDarkening = cms.uint32(1) process.mix.digitizers.hcal.HFDarkening = cms.bool(True) #these lines need to be further activated by tuning on 'complete' aging for HF @@ -225,7 +225,7 @@ def ecal_complete_aging(process): def turn_off_HE_aging(process): if hasattr(process,'mix') and hasattr(process.mix,'digitizers') and hasattr(process.mix.digitizers,'hcal'): - process.mix.digitizers.hcal.HEDarkening = cms.bool(False) + process.mix.digitizers.hcal.HEDarkening = cms.uint32(0) if hasattr(process,'es_hardcode'): process.es_hardcode.HERecalibration = cms.bool(False) if hasattr(process,'simHcalDigis'): diff --git a/SimCalorimetry/HcalSimProducers/python/hcalUnsuppressedDigis_cfi.py b/SimCalorimetry/HcalSimProducers/python/hcalUnsuppressedDigis_cfi.py index 1b1824b569112..a2f4603012247 100644 --- a/SimCalorimetry/HcalSimProducers/python/hcalUnsuppressedDigis_cfi.py +++ b/SimCalorimetry/HcalSimProducers/python/hcalUnsuppressedDigis_cfi.py @@ -38,7 +38,7 @@ CorrFactorFile = cms.FileInPath("SimCalorimetry/HcalSimProducers/data/calor_corr01.txt"), HcalReLabel = HcalReLabel, DelivLuminosity = cms.double(0), - HEDarkening = cms.bool(False), + HEDarkening = cms.uint32(0), HFDarkening = cms.bool(False) ) diff --git a/SimCalorimetry/HcalSimProducers/src/HcalDigitizer.cc b/SimCalorimetry/HcalSimProducers/src/HcalDigitizer.cc index 46bf06e671fe7..b8ec01fb824d2 100644 --- a/SimCalorimetry/HcalSimProducers/src/HcalDigitizer.cc +++ b/SimCalorimetry/HcalSimProducers/src/HcalDigitizer.cc @@ -148,7 +148,7 @@ HcalDigitizer::HcalDigitizer(const edm::ParameterSet& ps) : bool doHBHEUpgrade = ps.getParameter("HBHEUpgradeQIE"); bool doHFUpgrade = ps.getParameter("HFUpgradeQIE"); deliveredLumi = ps.getParameter("DelivLuminosity"); - bool agingFlagHE = ps.getParameter("HEDarkening"); + unsigned int agingFlagHE = ps.getParameter("HEDarkening"); bool agingFlagHF = ps.getParameter("HFDarkening"); // need to make copies, because they might get different noise generators @@ -327,7 +327,7 @@ HcalDigitizer::HcalDigitizer(const edm::ParameterSet& ps) : hitsProducer_ = ps.getParameter("hitsProducer"); - if(agingFlagHE) m_HEDarkening = new HEDarkening(); + if(agingFlagHE) m_HEDarkening = new HEDarkening(agingFlagHE); if(agingFlagHF) m_HFRecalibration = new HFRecalibration(); } diff --git a/SimG4CMS/Calo/src/HCalSD.cc b/SimG4CMS/Calo/src/HCalSD.cc index 785108d2fedbd..7d78950d9598c 100644 --- a/SimG4CMS/Calo/src/HCalSD.cc +++ b/SimG4CMS/Calo/src/HCalSD.cc @@ -64,7 +64,7 @@ HCalSD::HCalSD(G4String name, const DDCompactView & cpv, eminHitHF = m_HC.getParameter("EminHitHF")*MeV; useFibreBundle = m_HC.getParameter("UseFibreBundleHits"); deliveredLumi = m_HC.getParameter("DelivLuminosity"); - bool ageingFlagHE= m_HC.getParameter("HEDarkening"); + unsigned int ageingFlagHE= m_HC.getParameter("HEDarkening"); bool ageingFlagHF= m_HC.getParameter("HFDarkening"); useHF = m_HC.getUntrackedParameter("UseHF",true); bool forTBH2 = m_HC.getUntrackedParameter("ForTBH2",false); @@ -101,6 +101,7 @@ HCalSD::HCalSD(G4String name, const DDCompactView & cpv, << eminHitHO << " HF: " << eminHitHF << "\n" << "Delivered luminosity for Darkening " << deliveredLumi << " Flag (HE) " << ageingFlagHE + << " (" << HEDarkening::scenarioDescription(ageingFlagHE) << ")" << " Flag (HF) " << ageingFlagHF << "\n" << "Application of Fiducial Cut " << applyFidCut; @@ -323,7 +324,7 @@ HCalSD::HCalSD(G4String name, const DDCompactView & cpv, for (int i=0; i<9; ++i) hit_[i] = time_[i]= dist_[i] = 0; hzvem = hzvhad = 0; - if (ageingFlagHE) m_HEDarkening = new HEDarkening(); + if (ageingFlagHE) m_HEDarkening = new HEDarkening(ageingFlagHE); if (ageingFlagHF) m_HFDarkening = new HFDarkening(); #ifdef plotDebug edm::Service tfile; diff --git a/SimG4Core/Application/python/g4SimHits_cfi.py b/SimG4Core/Application/python/g4SimHits_cfi.py index fea5e6a8e948e..11569ad8b7edb 100644 --- a/SimG4Core/Application/python/g4SimHits_cfi.py +++ b/SimG4Core/Application/python/g4SimHits_cfi.py @@ -269,7 +269,7 @@ BetaThreshold = cms.double(0.7), TimeSliceUnit = cms.double(1), IgnoreTrackID = cms.bool(False), - HEDarkening = cms.bool(False), + HEDarkening = cms.uint32(0), HFDarkening = cms.bool(False), UseHF = cms.untracked.bool(True), ForTBH2 = cms.untracked.bool(False),