Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension


Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
10 changes: 5 additions & 5 deletions .github/workflows/main.yml
Original file line number Diff line number Diff line change
Expand Up @@ -207,7 +207,7 @@ jobs:
-t 0 --unblind --doImpacts --globalImpacts
--saveHists --saveHistsPerProcess --computeHistErrors --computeHistErrorsPerProcess
--computeHistCov --computeHistImpacts --computeVariations
-m Basemodel
-m BaseMapping
-m Select ch0 'x:rebin(-5,-2,2,5)'
-m Project ch1 a
-m Project ch1 b
Expand All @@ -217,15 +217,15 @@ jobs:
-m Ratio ch0 ch1 'None:None' 'b:slice(0,2),b:sum'
-m Ratio ch0 ch1 'x:slice(None,None,2),x:sum' 'a:sum,b:sum'
-m Normratio ch1 ch1 sig sig,bkg 'b:sum' 'b:2'
-m tests.param_model.Param 'slope_signal' 1 0
-m tests.param_mapping.Param 'slope_signal' 1 0

- name: nominal fit
run: >-
rabbit_fit.py $RABBIT_OUTDIR/test_tensor.hdf5 -o $RABBIT_OUTDIR/ --postfix composite
-t 0 --unblind '.*' --doImpacts --globalImpacts
--saveHists --saveHistsPerProcess --computeHistErrors --computeHistErrorsPerProcess
--computeHistCov --computeHistImpacts --computeVariations
--compositeModel -m Project ch1 a -m Project ch1 b
--compositeMapping -m Project ch1 a -m Project ch1 b

- name: nominal fit blinded
run: >-
Expand Down Expand Up @@ -320,13 +320,13 @@ jobs:
- name: plot prefit distributions
run: >-
rabbit_plot_hists.py $RABBIT_OUTDIR/fitresults.hdf5 -o $WEB_DIR/$PLOT_DIR
--extraTextLoc '0.05' '0.7' --legCols 1 -m Basemodel -m Project ch1 a -m Project ch1 b --yscale '1.2'
--extraTextLoc '0.05' '0.7' --legCols 1 -m BaseMapping -m Project ch1 a -m Project ch1 b --yscale '1.2'
--subtitle Preliminary --prefit --config tests/style_config.py

- name: plot postfit distributions
run: >-
rabbit_plot_hists.py $RABBIT_OUTDIR/fitresults.hdf5 -o $WEB_DIR/$PLOT_DIR
--extraTextLoc '0.05' '0.7' --legCols 1 -m Basemodel -m Project ch1 a -m Project ch1 b -m Normalize ch0 x --yscale '1.2'
--extraTextLoc '0.05' '0.7' --legCols 1 -m BaseMapping -m Project ch1 a -m Project ch1 b -m Normalize ch0 x --yscale '1.2'
--subtitle Preliminary --varNames 'slope_signal' --varColors orange --config tests/style_config.py --rrange 0.8 1.2

- name: plot histogram uncertainties
Expand Down
18 changes: 10 additions & 8 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -122,19 +122,19 @@ rabbit_fit test_tensor.hdf5 -o results/fitresult.hdf5 -t 0 --doImpacts --globalI
### Bin-by-bin statistical uncertainties
Bin-by-bin statistical uncertainties on the templates are added by default and can be disabled at runtime using the `--noBinByBinStat` option. The Barlow-Beeston lite method is used to add implicit nuisance parameters for each template bin. By default this is implemented using a gamma distribution for the probability density, but Gaussian uncertainties can also be used with `--binByBinStatType normal`.

### Physics models
Physics models are used to perform transformation on the parameters and observables (the histogram bins in the (masked) channels).
Baseline models are defined in `rabbit/physicsmodels/` and can be called in `rabbit_fit` with the `--PhysicsModel` or `-m` option e.g. `-m Select ch0 -m Project ch1 b`.
The first argument is the physics model name followed by arguments passed into the physics model.
### Mappings
Perform mappings on the parameters and observables (the histogram bins in the (masked) channels).
Baseline mappings are defined in `rabbit/mappings/` and can be called in `rabbit_fit` with the `--mapping` or `-m` option e.g. `-m Select ch0 -m Project ch1 b`.
The first argument is the mapping name followed by arguments passed into the mapping.
Available physics models are
* `Basemodel`: Compute histograms in all bins and all channels.
* `BaseMapping`: Compute histograms in all bins and all channels.
* `Select`: To select histograms of a channel, and perform a selection of processes and bins, supporting rebinning.
* `Project`: To project histograms to lower dimensions, respecting the covariance matrix across bins.
* `Normalize`: To normalize histograms to their sum (and project them) e.g. to compute normalized differential cross sections.
* `Ratio`: To compute the ratio between channels, processes, or histogram bins.
* `Normratio`: To compute the ratio of normalized histograms.

Models can be specified in the comand line and can feature different parsing syntax.
Mapping can be specified in the comand line and can feature different parsing syntax.
A convension is set up for parsing process and axes selections (e.g. in the `Select` and `Ratio` models). For selecting processes a comma separated list, e.g. <process_0>,<process_1>...
and for axes selecitons <axis_name_0>:<selection_0>,<axis_name_1>:<selection_1>,... i.e. a comma separated list of axis names and selections separated by ":".
Selections can be
Expand All @@ -145,10 +145,12 @@ Selections can be
- `None:None` for whch `None` is returned, indicating no selection
Multiple selection per axis can be specified, e.g. `x:slice(2,8),x:sum`.

Custom physics models can be used to make the desired transformation.
They can be specified with the full path to the custom model e.g. `-m custom_models.MyCustomModel`.
Custom mappings can be defined.
They can be specified with the full path to the custom mapping e.g. `-m custom_mapping.MyCustomMapping`.
The path must be accessable from your `$PYTHONPATH` variable and an `__ini__.py` file must be in the directory.

### Physics models
TBD

## Fit diagnostics

Expand Down
86 changes: 43 additions & 43 deletions bin/rabbit_fit.py
Original file line number Diff line number Diff line change
Expand Up @@ -12,8 +12,8 @@
import scipy

from rabbit import fitter, inputdata, io_tools, workspace
from rabbit.physicsmodels import helpers as ph
from rabbit.physicsmodels import physicsmodel as pm
from rabbit.mappings import helpers as mh
from rabbit.mappings import mapping as mp
from rabbit.tfhelpers import edmval_cov

from wums import output_tools, logging # isort: skip
Expand Down Expand Up @@ -282,22 +282,22 @@ def make_parser():
)
parser.add_argument(
"-m",
"--physicsModel",
"--mapping",
nargs="+",
action="append",
default=[],
help="""
add physics model to perform transformations on observables for the prefit and postfit histograms,
specifying the model defined in rabbit/physicsmodels/ followed by arguments passed in the model __init__,
perform mappings on observables or parameters for the prefit and postfit histograms,
specifying the mapping defined in rabbit/mappings/ followed by arguments passed in the mapping __init__,
e.g. '-m Project ch0 eta pt' to get a 2D projection to eta-pt or '-m Project ch0' to get the total yield.
This argument can be called multiple times.
Custom models can be specified with the full path to the custom model e.g. '-m custom_models.MyCustomModel'.
Custom mappings can be specified with the full path to the custom mapping e.g. '-m custom_mappings.MyCustomMapping'.
""",
)
parser.add_argument(
"--compositeModel",
"--compositeMapping",
action="store_true",
help="Make a composite model and compute the covariance matrix across all physics models.",
help="Make a composite mapping and compute the covariance matrix across all mappings.",
)
parser.add_argument(
"--doImpacts",
Expand Down Expand Up @@ -368,42 +368,42 @@ def make_parser():
return parser.parse_args()


def save_observed_hists(args, models, fitter, ws):
for model in models:
if not model.has_data:
def save_observed_hists(args, mappings, fitter, ws):
for mapping in mappings:
if not mapping.has_data:
continue

print(f"Save data histogram for {model.key}")
print(f"Save data histogram for {mapping.key}")
ws.add_observed_hists(
model,
mapping,
fitter.indata.data_obs,
fitter.nobs.value(),
fitter.indata.data_cov_inv,
fitter.data_cov_inv,
)


def save_hists(args, models, fitter, ws, prefit=True, profile=False):
def save_hists(args, mappings, fitter, ws, prefit=True, profile=False):

for model in models:
logger.info(f"Save inclusive histogram for {model.key}")
for mapping in mappings:
logger.info(f"Save inclusive histogram for {mapping.key}")

if prefit and model.skip_prefit:
if prefit and mapping.skip_prefit:
continue

if not model.skip_incusive:
if not mapping.skip_incusive:
exp, aux = fitter.expected_events(
model,
mapping,
inclusive=True,
compute_variance=args.computeHistErrors,
compute_cov=args.computeHistCov,
compute_chi2=not args.noChi2 and model.has_data,
compute_chi2=not args.noChi2 and mapping.has_data,
compute_global_impacts=args.computeHistImpacts,
profile=profile,
)

ws.add_expected_hists(
model,
mapping,
exp,
var=aux[0],
cov=aux[1],
Expand All @@ -413,20 +413,20 @@ def save_hists(args, models, fitter, ws, prefit=True, profile=False):
)

if aux[4] is not None:
ws.add_chi2(aux[4], aux[5], prefit, model)
ws.add_chi2(aux[4], aux[5], prefit, mapping)

if args.saveHistsPerProcess and not model.skip_per_process:
logger.info(f"Save processes histogram for {model.key}")
if args.saveHistsPerProcess and not mapping.skip_per_process:
logger.info(f"Save processes histogram for {mapping.key}")

exp, aux = fitter.expected_events(
model,
mapping,
inclusive=False,
compute_variance=args.computeHistErrorsPerProcess,
profile=profile,
)

ws.add_expected_hists(
model,
mapping,
exp,
var=aux[0],
process_axis=True,
Expand All @@ -439,15 +439,15 @@ def save_hists(args, models, fitter, ws, prefit=True, profile=False):
fitter.cov.assign(fitter.prefit_covariance(unconstrained_err=1.0))

exp, aux = fitter.expected_events(
model,
mapping,
inclusive=True,
compute_variance=False,
compute_variations=True,
profile=profile,
)

ws.add_expected_hists(
model,
mapping,
exp,
var=aux[0],
variations=True,
Expand Down Expand Up @@ -655,18 +655,18 @@ def main():
indata = inputdata.FitInputData(args.filename, args.pseudoData)
ifitter = fitter.Fitter(indata, args, do_blinding=any(blinded_fits))

# physics models for observables and transformations
if len(args.physicsModel) == 0 and args.saveHists:
# if no model is explicitly added and --saveHists is specified, fall back to Basemodel
args.physicsModel = [["Basemodel"]]
models = []
for margs in args.physicsModel:
model = ph.instance_from_class(margs[0], indata, *margs[1:])
models.append(model)

if args.compositeModel:
models = [
pm.CompositeModel(models),
# mappings for observables and parameters
if len(args.mapping) == 0 and args.saveHists:
# if no mapping is explicitly added and --saveHists is specified, fall back to BaseMapping
args.mapping = [["BaseMapping"]]
mappings = []
for margs in args.mapping:
mapping = mh.instance_from_class(margs[0], indata, *margs[1:])
mappings.append(mapping)

if args.compositeMapping:
mappings = [
mp.CompositeMapping(mappings),
]

np.random.seed(args.seed)
Expand Down Expand Up @@ -746,8 +746,8 @@ def main():
)

if args.saveHists:
save_observed_hists(args, models, ifitter, ws)
save_hists(args, models, ifitter, ws, prefit=True)
save_observed_hists(args, mappings, ifitter, ws)
save_hists(args, mappings, ifitter, ws, prefit=True)
prefit_time.append(time.time())

if not args.prefitOnly:
Expand All @@ -758,7 +758,7 @@ def main():
if args.saveHists:
save_hists(
args,
models,
mappings,
ifitter,
ws,
prefit=False,
Expand Down
32 changes: 16 additions & 16 deletions bin/rabbit_plot_hists.py
Original file line number Diff line number Diff line change
Expand Up @@ -203,13 +203,13 @@ def parseArgs():
)
parser.add_argument(
"-m",
"--physicsModel",
"--mapping",
nargs="+",
action="append",
default=[],
help="""
Make plot of physics model prefit and postfit histograms. Loop over all by deault.
Can also specify the model name, followed by the arguments, e.g. "-m Project ch0 eta pt".
Make plot of mapping prefit and postfit histograms. Loop over all by deault.
Can also specify the mapping name, followed by the arguments, e.g. "-m Project ch0 eta pt".
This argument can be called multiple times.
""",
)
Expand Down Expand Up @@ -258,7 +258,7 @@ def parseArgs():
type=str,
default="automatic",
choices=["automatic", "saturated", "linear", " ", "none", None],
help="Type of chi2 to print on plot (saturated from fit likelihood. linear from observables, or none) 'automatic' means pick saturated for basemodel and otherwise linear",
help="Type of chi2 to print on plot (saturated from fit likelihood. linear from observables, or none) 'automatic' means pick saturated for basemapping and otherwise linear",
)
parser.add_argument(
"--dataName", type=str, default="Data", help="Data name for plot labeling"
Expand Down Expand Up @@ -1474,16 +1474,16 @@ def main():
varMarkers=args.varMarkers,
)

results = fitresult["physics_models"]
for margs in args.physicsModel:
results = fitresult.get("mappings", fitresult.get("physics_models"))
for margs in args.mapping:
if margs == []:
instance_keys = results.keys()
else:
model_key = " ".join(margs)
instance_keys = [k for k in results.keys() if k.startswith(model_key)]
mapping_key = " ".join(margs)
instance_keys = [k for k in results.keys() if k.startswith(mapping_key)]
if len(instance_keys) == 0:
raise ValueError(
f"No model found under {model_key}. Available models: {results.keys()}."
f"No mapping found under {mapping_key}. Available mappings: {results.keys()}."
)

for instance_key in instance_keys:
Expand All @@ -1499,7 +1499,7 @@ def main():
fitresult
if fittype == "postfit"
and (
(instance_key == "Basemodel" and args.chisq != "linear")
(instance_key == "BaseMapping" and args.chisq != "linear")
or args.chisq == "saturated"
)
else instance
Expand All @@ -1513,7 +1513,7 @@ def main():
continue
logger.info(f"Make plot for {instance_key} in channel {channel}")

if instance_key == "CompositeModel":
if instance_key == "CompositeMapping":
info = channel_info.get(" ".join(channel.split(" ")[-1:]), {})
else:
info = channel_info.get(channel, {})
Expand All @@ -1538,11 +1538,11 @@ def main():
if any(
instance_key.startswith(x)
for x in [
"Basemodel",
"BaseMapping",
"Project",
"Select",
"Norm",
"CompositeModel",
"CompositeMapping",
]
)
and not args.noBinWidthNorm
Expand All @@ -1552,9 +1552,9 @@ def main():
opts["counts"] = counts

varResults = [
r["physics_models"][instance_key.replace("_masked", "")][
"channels"
][channel.replace("_masked", "")]
r.get("mappings", r.get("physics_models"))[
instance_key.replace("_masked", "")
]["channels"][channel.replace("_masked", "")]
for r in varFitresults
]

Expand Down
Loading