diff --git a/assimilation_code/modules/assimilation/bnrh_distribution_mod.f90 b/assimilation_code/modules/assimilation/bnrh_distribution_mod.f90 index 57b1590275..97f7a56bdc 100644 --- a/assimilation_code/modules/assimilation/bnrh_distribution_mod.f90 +++ b/assimilation_code/modules/assimilation/bnrh_distribution_mod.f90 @@ -409,24 +409,19 @@ subroutine inv_bnrh_cdf(quantiles, ens_size, sort_ens, & x(i) = lower_bound + (curr_q / q(1)) * (upper_state - lower_bound) else ! Find the mass at the lower bound (which could be unbounded) - ! NOTE: The amplitude here should be one since there is no likelihood. However, there is - ! round-off error that occurs in the statement below. In the long term, amp_adj should be - ! removed from this code block and the onn for the upper region. However, removing it now - ! would require resetting the baseline for the large number of baseline archived experiments. - amp_adj = q(1) / del_q if(bounded_below) then - lower_mass = amp_adj * tail_amp_left * & + lower_mass = tail_amp_left * & normal_cdf(lower_bound, tail_mean_left, tail_sd_left) else lower_mass = 0.0_r8 endif ! Find the mass at the upper bound (ensemble member 1) - upper_mass = amp_adj * tail_amp_left * & + upper_mass = tail_amp_left * & normal_cdf(sort_ens(1), tail_mean_left, tail_sd_left) ! What fraction of this mass difference should we go? fract = curr_q / q(1) target_mass = lower_mass + fract * (upper_mass - lower_mass) - x(i) = inv_weighted_normal_cdf(amp_adj*tail_amp_left, tail_mean_left, & + x(i) = inv_weighted_normal_cdf(tail_amp_left, tail_mean_left, & tail_sd_left, target_mass) endif @@ -441,20 +436,19 @@ subroutine inv_bnrh_cdf(quantiles, ens_size, sort_ens, & else ! Upper tail is (bounded) normal ! Find the mass at the upper bound (which could be unbounded) - amp_adj = (1.0_r8 - q(ens_size)) / del_q if(bounded_above) then - upper_mass = amp_adj * tail_amp_right * & + upper_mass = tail_amp_right * & normal_cdf(upper_bound, tail_mean_right, tail_sd_right) else - upper_mass = amp_adj * 1.0_r8 + upper_mass = 1.0_r8 endif ! Find the mass at the lower edge of the region (ensemble member n) - lower_mass = amp_adj * tail_amp_right * & + lower_mass = tail_amp_right * & normal_cdf(sort_ens(ens_size), tail_mean_right, tail_sd_right) ! What fraction of the last interval do we need to move fract = (curr_q - q(ens_size)) / (1.0_r8 - q(ens_size)) target_mass = lower_mass + fract * (upper_mass - lower_mass) - x(i) = inv_weighted_normal_cdf(amp_adj * tail_amp_right, tail_mean_right, & + x(i) = inv_weighted_normal_cdf(tail_amp_right, tail_mean_right, & tail_sd_right, target_mass) endif