Skip to content

Commit 67b0519

Browse files
committed
fixes for issue 3314
1 parent b0fe77b commit 67b0519

4 files changed

Lines changed: 89 additions & 74 deletions

File tree

Lines changed: 14 additions & 14 deletions
Original file line numberDiff line numberDiff line change
@@ -1,16 +1,16 @@
11
esa_mode,incident_E-Step,Observed_E-Step,Cntr_E,Cntr_E_unc,GF_Trpl_H,GF_Trpl_H_unc_minus,GF_Trpl_H_unc_plus
22
# [1],[1],[1],[keV],[keV],[cm^2sr keV/keV],[cm^2sr keV/keV],[cm^2sr keV/keV]
3-
0,1,1,0.01633,0.28,4.45E-05,4.03E-05,4.20E-05
4-
0,2,2,0.03047,0.43,5.02E-05,4.53E-05,4.72E-05
5-
0,3,3,0.05576,0.89,6.16E-05,5.58E-05,5.82E-05
6-
0,4,4,0.10626,1.9,7.12E-05,4.51E-05,4.89E-05
7-
0,5,5,0.20004,3.4,8.89E-05,5.85E-05,6.31E-05
8-
0,6,6,0.40496,7.3,1.12E-04,6.58E-05,7.23E-05
9-
0,7,7,0.78729,22,1.43E-04,8.25E-05,9.09E-05
10-
1,1,1,0.01719,0.22,1.05E-04,8.82E-05,9.26E-05
11-
1,2,2,0.03236,0.36,1.22E-04,1.02E-04,1.07E-04
12-
1,3,3,0.05948,0.77,1.44E-04,1.21E-04,1.27E-04
13-
1,4,4,0.11441,1.3,1.73E-04,1.17E-04,1.26E-04
14-
1,5,5,0.2137,3,2.15E-04,1.37E-04,1.49E-04
15-
1,6,6,0.43736,5.3,2.82E-04,1.65E-04,1.81E-04
16-
1,7,7,0.83888,11.7,3.61E-04,2.08E-04,2.29E-04
3+
0,1,1,0.01633,0.00028,4.45E-05,4.03E-05,4.20E-05
4+
0,2,2,0.03047,0.00043,5.02E-05,4.53E-05,4.72E-05
5+
0,3,3,0.05576,0.00089,6.16E-05,5.58E-05,5.82E-05
6+
0,4,4,0.10626,0.0019,7.12E-05,4.51E-05,4.89E-05
7+
0,5,5,0.20004,0.0034,8.89E-05,5.85E-05,6.31E-05
8+
0,6,6,0.40496,0.0073,1.12E-04,6.58E-05,7.23E-05
9+
0,7,7,0.78729,0.022,1.43E-04,8.25E-05,9.09E-05
10+
1,1,1,0.01719,0.00022,1.05E-04,8.82E-05,9.26E-05
11+
1,2,2,0.03236,0.00036,1.22E-04,1.02E-04,1.07E-04
12+
1,3,3,0.05948,0.00077,1.44E-04,1.21E-04,1.27E-04
13+
1,4,4,0.11441,0.0013,1.73E-04,1.17E-04,1.26E-04
14+
1,5,5,0.2137,0.003,2.15E-04,1.37E-04,1.49E-04
15+
1,6,6,0.43736,0.0053,2.82E-04,1.65E-04,1.81E-04
16+
1,7,7,0.83888,0.0117,3.61E-04,2.08E-04,2.29E-04
Lines changed: 14 additions & 14 deletions
Original file line numberDiff line numberDiff line change
@@ -1,16 +1,16 @@
11
esa_mode,incident_E-Step,Observed_E-Step,Cntr_E,Cntr_E_unc,GF_Trpl_O,GF_Trpl_O_unc_minus,GF_Trpl_O_unc_plus
22
# [1],[1],[1],[keV],[keV],[cm^2sr keV/keV],[cm^2sr keV/keV],[cm^2sr keV/keV]
3-
0,1,1,0.01919,1.372344105,1.98E-05,1.74E-05,1.81E-05
4-
0,2,2,0.03675,2.537963082,3.18E-05,2.39E-05,2.53E-05
5-
0,3,3,0.07121,6.591339559,5.34E-05,4.05E-05,4.29E-05
6-
0,4,4,0.14141,17.43990907,8.64E-05,6.14E-05,6.56E-05
7-
0,5,5,0.274,35.83900763,1.32E-04,9.07E-05,9.72E-05
8-
0,6,6,0.58503,98.06681395,2.46E-04,1.53E-04,1.67E-04
9-
0,7,7,1.13506,238.1076888,3.81E-04,2.23E-04,2.45E-04
10-
1,1,1,0.02043303552,1.575924599,5.02E-05,4.41E-05,4.61E-05
11-
1,2,2,0.03972,2.953020557,8.13E-05,6.09E-05,6.47E-05
12-
1,3,3,0.07648,7.476494389,1.33E-04,1.01E-04,1.07E-04
13-
1,4,4,0.15353,14.87953071,2.38E-04,1.69E-04,1.80E-04
14-
1,5,5,0.29846,41.48950349,3.87E-04,2.67E-04,2.86E-04
15-
1,6,6,0.61524,83.57629594,6.61E-04,4.11E-04,4.47E-04
16-
1,7,7,1.24282,212.0806857,1.02E-03,6.10E-04,6.67E-04
3+
0,1,1,0.01919,0.001372344105,1.98E-05,1.74E-05,1.81E-05
4+
0,2,2,0.03675,0.002537963082,3.18E-05,2.39E-05,2.53E-05
5+
0,3,3,0.07121,0.006591339559,5.34E-05,4.05E-05,4.29E-05
6+
0,4,4,0.14141,0.01743990907,8.64E-05,6.14E-05,6.56E-05
7+
0,5,5,0.274,0.03583900763,1.32E-04,9.07E-05,9.72E-05
8+
0,6,6,0.58503,0.09806681395,2.46E-04,1.53E-04,1.67E-04
9+
0,7,7,1.13506,0.2381076888,3.81E-04,2.23E-04,2.45E-04
10+
1,1,1,0.02043303552,0.001575924599,5.02E-05,4.41E-05,4.61E-05
11+
1,2,2,0.03972,0.002953020557,8.13E-05,6.09E-05,6.47E-05
12+
1,3,3,0.07648,0.007476494389,1.33E-04,1.01E-04,1.07E-04
13+
1,4,4,0.15353,0.01487953071,2.38E-04,1.69E-04,1.80E-04
14+
1,5,5,0.29846,0.04148950349,3.87E-04,2.67E-04,2.86E-04
15+
1,6,6,0.61524,0.08357629594,6.61E-04,4.11E-04,4.47E-04
16+
1,7,7,1.24282,0.2120806857,1.02E-03,6.10E-04,6.67E-04

imap_processing/lo/l2/lo_l2.py

Lines changed: 44 additions & 27 deletions
Original file line numberDiff line numberDiff line change
@@ -1009,16 +1009,24 @@ def calculate_intensities(dataset: xr.Dataset) -> xr.Dataset:
10091009
dataset["counts_over_eff_squared"]
10101010
) / (dataset["geometric_factor"] * dataset["energy"] * dataset["exposure_factor"])
10111011

1012-
for suffix in ("minus", "plus"):
1013-
dataset[f"ena_intensity_sys_err_{suffix}"] = (
1014-
dataset["ena_intensity"]
1015-
* dataset[f"geometric_factor_stat_uncert_{suffix}"]
1016-
/ dataset["geometric_factor"]
1017-
)
1012+
plus_multiplier = dataset["geometric_factor"] / (
1013+
dataset["geometric_factor"] - dataset["geometric_factor_stat_uncert_minus"]
1014+
)
1015+
minus_multiplier = dataset["geometric_factor"] / (
1016+
dataset["geometric_factor"] + dataset["geometric_factor_stat_uncert_plus"]
1017+
)
1018+
1019+
dataset["ena_intensity_sys_err_plus"] = (
1020+
dataset["ena_intensity"] * plus_multiplier
1021+
) - dataset["ena_intensity"]
10181022

1019-
# Symmetric systematic error (mean of the asymmetric minus/plus bounds)
1020-
dataset["ena_intensity_sys_err"] = 0.5 * (
1021-
dataset["ena_intensity_sys_err_minus"] + dataset["ena_intensity_sys_err_plus"]
1023+
dataset["ena_intensity_sys_err_minus"] = dataset["ena_intensity"] - (
1024+
dataset["ena_intensity"] * minus_multiplier
1025+
)
1026+
1027+
# Symmetric systematic error
1028+
dataset["ena_intensity_sys_err"] = np.sqrt(
1029+
dataset["ena_intensity_sys_err_minus"] * dataset["ena_intensity_sys_err_plus"]
10221030
)
10231031

10241032
return dataset
@@ -1048,16 +1056,24 @@ def calculate_backgrounds(dataset: xr.Dataset) -> xr.Dataset:
10481056
/ dataset["exposure_factor"] ** 2
10491057
)
10501058

1051-
for suffix in ("minus", "plus"):
1052-
dataset[f"bg_rate_sys_err_{suffix}"] = (
1053-
dataset["bg_rate"]
1054-
* dataset[f"geometric_factor_stat_uncert_{suffix}"]
1055-
/ dataset["geometric_factor"]
1056-
)
1059+
plus_multiplier = dataset["geometric_factor"] / (
1060+
dataset["geometric_factor"] - dataset["geometric_factor_stat_uncert_minus"]
1061+
)
1062+
minus_multiplier = dataset["geometric_factor"] / (
1063+
dataset["geometric_factor"] + dataset["geometric_factor_stat_uncert_plus"]
1064+
)
1065+
1066+
dataset["bg_rate_sys_err_plus"] = (dataset["bg_rate"] * plus_multiplier) - dataset[
1067+
"bg_rate"
1068+
]
1069+
1070+
dataset["bg_rate_sys_err_minus"] = dataset["bg_rate"] - (
1071+
dataset["bg_rate"] * minus_multiplier
1072+
)
10571073

1058-
# Symmetric systematic error (mean of the asymmetric minus/plus bounds)
1059-
dataset["bg_rate_sys_err"] = 0.5 * (
1060-
dataset["bg_rate_sys_err_minus"] + dataset["bg_rate_sys_err_plus"]
1074+
# Symmetric systematic error
1075+
dataset["bg_rate_sys_err"] = np.sqrt(
1076+
dataset["bg_rate_sys_err_minus"] * dataset["bg_rate_sys_err_plus"]
10611077
)
10621078

10631079
# Background intensity
@@ -1068,16 +1084,17 @@ def calculate_backgrounds(dataset: xr.Dataset) -> xr.Dataset:
10681084
dataset["geometric_factor"] * dataset["energy"]
10691085
)
10701086

1071-
for suffix in ("minus", "plus"):
1072-
dataset[f"bg_intensity_sys_err_{suffix}"] = (
1073-
dataset["bg_intensity"]
1074-
* dataset[f"geometric_factor_stat_uncert_{suffix}"]
1075-
/ dataset["geometric_factor"]
1076-
)
1087+
dataset["bg_intensity_sys_err_plus"] = (
1088+
dataset["bg_intensity"] * plus_multiplier
1089+
) - dataset["bg_intensity"]
1090+
1091+
dataset["bg_intensity_sys_err_minus"] = dataset["bg_intensity"] - (
1092+
dataset["bg_intensity"] * minus_multiplier
1093+
)
10771094

1078-
# Symmetric systematic error (mean of the asymmetric minus/plus bounds)
1079-
dataset["bg_intensity_sys_err"] = 0.5 * (
1080-
dataset["bg_intensity_sys_err_minus"] + dataset["bg_intensity_sys_err_plus"]
1095+
# Symmetric systematic error
1096+
dataset["bg_intensity_sys_err"] = np.sqrt(
1097+
dataset["bg_intensity_sys_err_minus"] * dataset["bg_intensity_sys_err_plus"]
10811098
)
10821099

10831100
return dataset

imap_processing/tests/lo/test_lo_l2.py

Lines changed: 17 additions & 19 deletions
Original file line numberDiff line numberDiff line change
@@ -1134,17 +1134,17 @@ def test_calculate_intensities_basic(self, sample_dataset_with_geometric_factors
11341134
result["ena_intensity_stat_uncert"], expected_stat_uncert
11351135
)
11361136

1137-
# Check systematic uncertainty calculation. The single `_sys_err` is the
1138-
# mean of the asymmetric minus/plus bounds.
1139-
mean_gf_stat_uncert = 0.5 * (
1140-
sample_dataset_with_geometric_factors["geometric_factor_stat_uncert_minus"]
1141-
+ sample_dataset_with_geometric_factors["geometric_factor_stat_uncert_plus"]
1142-
)
1143-
expected_sys_err = (
1144-
result["ena_intensity"]
1145-
* mean_gf_stat_uncert
1146-
/ sample_dataset_with_geometric_factors["geometric_factor"]
1147-
)
1137+
# Check systematic uncertainty calculation
1138+
gf = sample_dataset_with_geometric_factors["geometric_factor"]
1139+
dg_minus = sample_dataset_with_geometric_factors[
1140+
"geometric_factor_stat_uncert_minus"
1141+
]
1142+
dg_plus = sample_dataset_with_geometric_factors[
1143+
"geometric_factor_stat_uncert_plus"
1144+
]
1145+
expected_sys_err_plus = result["ena_intensity"] * dg_minus / (gf - dg_minus)
1146+
expected_sys_err_minus = result["ena_intensity"] * dg_plus / (gf + dg_plus)
1147+
expected_sys_err = np.sqrt(expected_sys_err_minus * expected_sys_err_plus)
11481148
xr.testing.assert_allclose(result["ena_intensity_sys_err"], expected_sys_err)
11491149

11501150
def test_calculate_intensities_missing_variables(self):
@@ -1197,14 +1197,12 @@ def test_calculate_backgrounds_basic(
11971197
xr.testing.assert_allclose(result["bg_rate_stat_uncert"], expected_stat_uncert)
11981198

11991199
# Check systematic uncertainty calculation
1200-
# (mean(geometric_factor_stat_uncert bounds) / geometric_factor) * bg_rate
1201-
mean_gf_stat_uncert = 0.5 * (
1202-
dataset["geometric_factor_stat_uncert_minus"]
1203-
+ dataset["geometric_factor_stat_uncert_plus"]
1204-
)
1205-
expected_sys_err = (
1206-
result["bg_rate"] * mean_gf_stat_uncert / dataset["geometric_factor"]
1207-
)
1200+
gf = dataset["geometric_factor"]
1201+
dg_minus = dataset["geometric_factor_stat_uncert_minus"]
1202+
dg_plus = dataset["geometric_factor_stat_uncert_plus"]
1203+
expected_sys_err_plus = result["bg_rate"] * dg_minus / (gf - dg_minus)
1204+
expected_sys_err_minus = result["bg_rate"] * dg_plus / (gf + dg_plus)
1205+
expected_sys_err = np.sqrt(expected_sys_err_minus * expected_sys_err_plus)
12081206
xr.testing.assert_allclose(result["bg_rate_sys_err"], expected_sys_err)
12091207

12101208
def test_calculate_backgrounds_zero_exposure(self):

0 commit comments

Comments
 (0)