Skip to content

Commit

Permalink
add option to include CRs and VRs in the full likelihood
Browse files Browse the repository at this point in the history
Option introduced in the `.info` file of the analysis e.g.
```
  <pyhf id="RegionC">
    <name>atlas_susy_2018_31_SRC.json</name>
    <regions>
      <channel name="CRtt_cuts"> </channel>
      <channel name="SR_metsigST">
      SRC_22
      SRC_24
      SRC_26
      SRC_28
      </channel>
      <channel name="VRzAlt_cuts" is_included="True"> </channel>
      <channel name="CRz_cuts" is_included="True"> </channel>
      <channel name="VRtt_cuts"> </channel>
      <channel name="VRz_cuts"> </channel>
    </regions>
  </pyhf>
```
If there is no data present in the channel code will look for `is_included`
tag which can be "True", "yes", "1". If true the channel will be included
in the likelihood calculation if not it will be removed. If `is_included`
is not present and channel does not have any data, it will be removed by
default.
  • Loading branch information
jackaraz committed Jan 11, 2022
1 parent e5fa95d commit 37093dd
Show file tree
Hide file tree
Showing 2 changed files with 27 additions and 20 deletions.
28 changes: 15 additions & 13 deletions madanalysis/misc/histfactory_reader.py
Original file line number Diff line number Diff line change
Expand Up @@ -34,19 +34,15 @@ def __init__(self,pyhf_config):
self.path = pyhf_config.get('path', 'missing_path')
self.name = pyhf_config.get('name', 'missing_name')
self.logger = logging.getLogger('MA5')
if type(self) == HF_Background:
if isinstance(self, HF_Background):
self.hf = {}
self.global_config = self.pyhf_config
elif type(self) == HF_Signal:
elif isinstance(self, HF_Signal):
self.hf = []

def __call__(self,lumi):
return self.extrapolate(lumi)

@classmethod
def __type__(self):
return self.__name__

def extrapolate(self,lumi):
""" To calculate the HL variables HF needs to be extrapolated. Expected
observables will be extrapolated and summed, summation is superseeded
Expand Down Expand Up @@ -107,7 +103,7 @@ def extrapolate(self,lumi):
if item != []:
HF['observations'][iobs]['data'] = item

elif isinstance(self,HF_Signal):#type(self) == HF_Signal:
elif isinstance(self, HF_Signal):#type(self) == HF_Signal:
# Signal extrapolation
for i in range(len(HF)):
if HF[i]['op'] == 'remove':
Expand Down Expand Up @@ -214,23 +210,27 @@ class HF_Signal(HistFactory):
is failed and correct pyhf_config is needed.
"""
def __init__(self,pyhf_config, regiondata, xsection=-1, **kwargs):
super(HF_Signal,self).__init__(pyhf_config)
super(HF_Signal, self).__init__(pyhf_config)
self.signal_config = {}

with open(os.path.join(self.path, self.name), 'r') as json_file:
tmp_bkg = json.load(json_file)

bin_sizes = [len(x.get('data', [])) for x in tmp_bkg.get('observations',[])]

for key, item in self.pyhf_config.items():
if key != 'lumi':
self.signal_config[key] = {}
if item['data'] == []:
if not item['is_included']:
self.signal_config[key]['op'] = 'remove'
self.signal_config[key]["path"] = '/channels/' + str(item['channels'])
else:
self.signal_config[key]['op'] = 'add'
self.signal_config[key]["path"] = \
'/channels/' + str(item['channels']) + '/samples/' + \
str(len(tmp_bkg["channels"][int(item['channels'])]["samples"])-1)
self.signal_config[key]["bin_size"] = \
bin_sizes[int(self.signal_config[key]["path"].split('/')[2])]

self.signal_config[key]['data'] = []
for SRname in item['data']:
Expand All @@ -247,7 +247,7 @@ def __init__(self,pyhf_config, regiondata, xsection=-1, **kwargs):
add_normsys = kwargs.get('add_normsys', []),
add_histosys = kwargs.get('add_histosys',[]),)

def set_HF(self,xsection,**kwargs):
def set_HF(self, xsection, **kwargs):
HF = []
if xsection<=0.:
return HF
Expand All @@ -269,6 +269,8 @@ def set_HF(self,xsection,**kwargs):
]
}
}
if len(SR_tmp["value"]["data"]) == 0:
SR_tmp["value"]["data"] = [0.0] * self.signal_config[SR]['bin_size']
HF.append(SR_tmp)
else:
toRemove.append(
Expand All @@ -287,13 +289,13 @@ def set_HF(self,xsection,**kwargs):
HF = self.add_normsys(HF,sys['hi_data'],sys['lo_data'],sys['name'])

background = kwargs.get('background',{})
if type(background) == HF_Background:
if not self.validate_bins(background,HF):
if isinstance(background, HF_Background):
if not self.validate_bins(background, HF):
self.logger.warning('Signal HistFactory validation failed.')
return []
return HF

def validate_bins(self,background,HF=[]):
def validate_bins(self, background, HF = []):
if HF == []:
HF = self.hf
bkg_bins = background.size()
Expand Down
19 changes: 12 additions & 7 deletions madanalysis/misc/run_recast.py
Original file line number Diff line number Diff line change
Expand Up @@ -923,11 +923,10 @@ def parse_info_file(self, etree, analysis, extrapolated_lumi):
return -1,-1, -1
## Getting the XML information
try:
info_input = open(filename)
info_tree = etree.parse(info_input)
info_input.close()
with open(filename, "r") as info_input:
info_tree = etree.parse(info_input)
except Exception as err:
self.logger.warning("Error during HMTL parsing: "+str(err))
self.logger.warning("Error during XML parsing: "+str(err))
self.logger.warning('Cannot parse the info file')
return -1,-1, -1

Expand Down Expand Up @@ -1251,9 +1250,15 @@ def pyhf_info_file(self,info_root):
if channel.text != None:
data = channel.text.split()
pyhf_config[likelihood_profile]['SR'][channel.attrib['name']] = {
'channels' : channel.attrib.get('id',-1),
'data' : data
"channels" : channel.get('id', default = -1),
"data" : data,
}
is_included = (
channel.get("is_included", default = 0) in ["True", "1", "yes"]
) if len(data) == 0 else True
pyhf_config[likelihood_profile]['SR'][channel.attrib['name']].update(
{"is_included" : is_included}
)
if pyhf_config[likelihood_profile]['SR'][channel.attrib['name']]['channels'] == -1:
file = os.path.join(
pyhf_config[likelihood_profile]['path'],
Expand All @@ -1276,7 +1281,7 @@ def pyhf_info_file(self,info_root):
continue
# validat pyhf config
background = HF_Background(config)
signal = HF_Signal(config,{},xsection=1.,
signal = HF_Signal(config, {}, xsection=1.,
background = background,
validate = True)

Expand Down

0 comments on commit 37093dd

Please sign in to comment.