diff --git a/src/inversion_scripts/invert.py b/src/inversion_scripts/invert.py index 777dd915..2407a077 100644 --- a/src/inversion_scripts/invert.py +++ b/src/inversion_scripts/invert.py @@ -150,7 +150,10 @@ def do_inversion( for p in obs_GC[:, 4] ] ) + # Define observational errors (diagonal entries of S_o matrix) + obs_error = np.power(obs_error, 2) gP = s_superO_p**2 / s_superO_1**2 + # scale error variance by gP obs_error = gP * obs_err # check to make sure obs_err isn't negative, set 1 as default value @@ -179,9 +182,6 @@ def do_inversion( scaling_matrix = np.tile(scale_factors, (reps, 1)) K *= scaling_matrix - # Define observational errors (diagonal entries of S_o matrix) - obs_error = np.power(obs_error, 2) - # Measurement-model mismatch: TROPOMI columns minus GEOS-Chem virtual TROPOMI columns # This is (y - F(xA)), i.e., (y - (K*xA + c)) or (y - K*xA) in shorthand delta_y = obs_GC[:, 0] - obs_GC[:, 1] # [ppb] diff --git a/src/inversion_scripts/lognormal_invert.py b/src/inversion_scripts/lognormal_invert.py index 4527f2f2..0360d0ce 100644 --- a/src/inversion_scripts/lognormal_invert.py +++ b/src/inversion_scripts/lognormal_invert.py @@ -82,7 +82,7 @@ def lognormal_invert(config, state_vector_filepath, jacobian_sf): # get the So matrix ds = np.load("so_super.npz") - so = ds["so"] ** 2 + so = ds["so"] # Calculate the difference between tropomi and the background # simulation, which has no emissions diff --git a/src/inversion_scripts/merge_partial_k.py b/src/inversion_scripts/merge_partial_k.py index 3de95977..ea2a56d3 100644 --- a/src/inversion_scripts/merge_partial_k.py +++ b/src/inversion_scripts/merge_partial_k.py @@ -86,7 +86,9 @@ def merge_partial_k(satdat_dir, lat_bounds, lon_bounds, obs_err, precomp_K): for p in obs_GC[:, 4] ] ) + # scale error variance by gP value following Chen et al. 2023 gP = s_superO_p**2 / s_superO_1**2 + obs_error = obs_error ** 2 obs_error = gP * obs_error # check to make sure obs_err isn't negative, set 1 as default value