Skip to content
65 changes: 37 additions & 28 deletions validphys2/src/validphys/dataplots.py
Original file line number Diff line number Diff line change
Expand Up @@ -263,10 +263,6 @@ def _plot_fancy_impl(
nkinlabels = len(table.columns)
ndata = len(table)

# Compute shifts due to the correlated part of the exp cov matrix
lcd_wc = loaded_commondata_with_cuts(commondata, cutlist[0])
Comment thread
enocera marked this conversation as resolved.
theory_predictions = results[1].central_value

# This is easier than cheking every time
if labellist is None:
labellist = [None] * len(results)
Expand All @@ -281,36 +277,49 @@ def _plot_fancy_impl(
# We modify the table, so we pass only the label columns
norm_cv, _ = transform_result(cv, err, table.iloc[:, :nkinlabels], info)

cvcols = []

# Compute systematic shifts
# For unknown reasons, `shifts_from_systematics` may
# randomly fails. If a LinAlgError is raised, shifts are not included in
# the final plot.
if with_shift:
try:
shifts, alpha = shifts_from_systematics(lcd_wc, theory_predictions)
except np.linalg.LinAlgError:
log.warning(
"Error occurred in computing systematic shifts for "
f"{info.ds_metadata.name}. These will be disregarded in the plots."
)
with_shift = False

cvcols = []

for i, (result, cuts) in enumerate(zip(results, cutlist)):
# We modify the table, so we pass only the label columns

mask = cut_mask(cuts)
cv = np.full(ndata, np.nan)
err = np.full(ndata, np.nan)
# Shift the theory when with_shift option is True
if i == 1 and with_shift:
cv[mask] = result.central_value - shifts

shifts = None
alpha = None
do_shift = with_shift

# Compute shifts due to the correlated part on the experiemntal ucnertainty
if with_shift:

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think it would be cleaner if this were instead:

Suggested change
if with_shift:
if do_shift:

where

do_shift = with_shift and not isintance(result, DataResult)

(you might need to import DataResults from results.py)

So that you don't enter the if condition with the Data, so you don't have to check if i>=1, which might not always be reliable if things change in validphys. You still need the alpha of course, so make sure that you add:

if with_shift and isintance(result, DataResult):
    err[mask] = compute_alpha(lcd_wc)

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

In my opinion, we have to enter the condition with the Data, also to get the uncorrelated part of the uncertainty, because these are Data with cuts, and cuts may differ.

lcd_wc = loaded_commondata_with_cuts(commondata, cuts)
theory_predictions = result.central_value

# For unknown reasons, `shifts_from_systematics` may randomly fail.
# If a LinAlgError is raised, shifts are not included in the final plot.
try:
shifts, alpha = shifts_from_systematics(lcd_wc, theory_predictions)

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think alpha should not be an output of this function. The complicated/costly part of this function (inverting a matrix, etc) doesn't need to be done for the data part (which is the only reason you might want to enter this loop with the data), right? I think it would be worth it to compute alpha, given solely lcd_wc as a separate function (shifts_from_systematics can also call that same function).

Btw, one thing that might reduce the report time is to make the result hashable (depends only on commondata and on the theory number -if a theory prediction-) such that the amount of shifts that need to be computed for a report are not duplicated... writing but probably an issue to tackle later.

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I have separated the computation of the shifts and the extraction of the uncorrelated and correlated parts part of the uncertainty.

except np.linalg.LinAlgError:
log.warning(
"Error occurred in computing systematic shifts for "
f"{info.ds_metadata.name}. These will be disregarded in the plots."
)
do_shift = False

# Shift theory predictions, but not data
if i >= 1 and do_shift:
cv[mask] = result.central_value - shifts
else:
cv[mask] = result.central_value

# Retain only the uncorrelated part of the data error if shifting the data
if i == 0 and do_shift:
err[mask] = alpha
else:
err[mask] = result.std_error

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
do_shift = False
# Shift theory predictions, but not data
if i >= 1 and do_shift:
cv[mask] = result.central_value - shifts
else:
cv[mask] = result.central_value
# Retain only the uncorrelated part of the data error if shifting the data
if i == 0 and do_shift:
err[mask] = alpha
else:
err[mask] = result.std_error
shifts = 0.0
cv[mask] = result.central_value - shifts
err[mask] = result.std_error

Since when you are dealing with the data, you are not entering in the if condition, this is safe.

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Done.


# No shifts
else:
cv[mask] = result.central_value
# Retain only the uncorrelated part of the error if shifting the data
if i == 0 and with_shift:
err[mask] = alpha
else:
err[mask] = result.std_error

cv, err = transform_result(cv, err, table.iloc[:, :nkinlabels], info)
Expand Down
Loading