Skip to content

Commit

Permalink
Merge pull request #5916 from fratnikov/he_darkrning_v1
Browse files Browse the repository at this point in the history 
introduce different scenarios for HE darkening
  • Loading branch information
@cmsbuild
cmsbuild committed Oct 22, 2014
2 parents a01ad1f + c00ea7b commit e243c3d
Show file tree
Hide file tree
Showing 7 changed files with 69 additions and 41 deletions.
24 changes: 17 additions & 7 deletions DataFormats/HcalCalibObjects/interface/HEDarkening.h
View file
Expand Up @@ -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];

};

Expand Down
69 changes: 43 additions & 26 deletions DataFormats/HcalCalibObjects/src/HEDarkening.cc
View file
Expand Up @@ -10,61 +10,78 @@

//#define DebugLog

HEDarkening::HEDarkening() {
HEDarkening::HEDarkening(unsigned int scenario) {
//HE starts at tower 16
ieta_shift = 16;

//scale parameter in fb-1 for exponential darkening from radiation damage for each tile
//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.;

//shift ieta tower index
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]));
}
}

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";
}
4 changes: 2 additions & 2 deletions SLHCUpgradeSimulations/Configuration/python/aging.py
View file
Expand Up @@ -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
Expand Down Expand Up @@ -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'):
Expand Down
2 changes: 1 addition & 1 deletion SimCalorimetry/HcalSimProducers/python/hcalUnsuppressedDigis_cfi.py
View file
Expand Up @@ -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)
)

Expand Down
4 changes: 2 additions & 2 deletions SimCalorimetry/HcalSimProducers/src/HcalDigitizer.cc
View file
Expand Up @@ -148,7 +148,7 @@ HcalDigitizer::HcalDigitizer(const edm::ParameterSet& ps) :
bool doHBHEUpgrade = ps.getParameter<bool>("HBHEUpgradeQIE");
bool doHFUpgrade = ps.getParameter<bool>("HFUpgradeQIE");
deliveredLumi = ps.getParameter<double>("DelivLuminosity");
bool agingFlagHE = ps.getParameter<bool>("HEDarkening");
unsigned int agingFlagHE = ps.getParameter<unsigned>("HEDarkening");
bool agingFlagHF = ps.getParameter<bool>("HFDarkening");

// need to make copies, because they might get different noise generators
Expand Down Expand Up @@ -327,7 +327,7 @@ HcalDigitizer::HcalDigitizer(const edm::ParameterSet& ps) :

hitsProducer_ = ps.getParameter<std::string>("hitsProducer");

if(agingFlagHE) m_HEDarkening = new HEDarkening();
if(agingFlagHE) m_HEDarkening = new HEDarkening(agingFlagHE);
if(agingFlagHF) m_HFRecalibration = new HFRecalibration();
}

Expand Down
5 changes: 3 additions & 2 deletions SimG4CMS/Calo/src/HCalSD.cc
View file
Expand Up @@ -64,7 +64,7 @@ HCalSD::HCalSD(G4String name, const DDCompactView & cpv,
eminHitHF = m_HC.getParameter<double>("EminHitHF")*MeV;
useFibreBundle = m_HC.getParameter<bool>("UseFibreBundleHits");
deliveredLumi = m_HC.getParameter<double>("DelivLuminosity");
bool ageingFlagHE= m_HC.getParameter<bool>("HEDarkening");
unsigned int ageingFlagHE= m_HC.getParameter<unsigned>("HEDarkening");
bool ageingFlagHF= m_HC.getParameter<bool>("HFDarkening");
useHF = m_HC.getUntrackedParameter<bool>("UseHF",true);
bool forTBH2 = m_HC.getUntrackedParameter<bool>("ForTBH2",false);
Expand Down Expand Up @@ -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;

Expand Down Expand Up @@ -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<TFileService> tfile;
Expand Down
2 changes: 1 addition & 1 deletion SimG4Core/Application/python/g4SimHits_cfi.py
View file
Expand Up @@ -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),
Expand Down

0 comments on commit e243c3d

Please sign in to comment.