CH-Earth / summa

Structure for Unifying Multiple Modeling Alternatives:
http://www.ral.ucar.edu/projects/summa
GNU General Public License v3.0
79 stars 103 forks source link

revise Gamma distribution used for sub-grid routing #470

Closed h294liu closed 2 years ago

h294liu commented 3 years ago

Make sure all the relevant boxes are checked (and only check the box if you actually completed the step):

The previous code has an error in calculating the cumulative probability "cumProb" of the Gamma distribution. Here the error is fixed by updating "sumProb" based on Equation 13 of Clark at el. (2007, WRR); removing variable "aLambda"; introducing variable "routingGammaMean".

The attached figures show the sensitivity of the updated Gamma distribution to the summa parameters routingGammaShape and routingGammaScale, respectively. shape_sensitivity scale_sensitivity

h294liu commented 3 years ago

Hi Andy, using "tFuture/routingGammaScale" is indeed simplier and correct as well. routingGammaMean is introduced because it appears in Equation 13 of Clark et al. (2007, WRR). In this way, the code stays consistent with the paper.

Regarding routingGammaScale, I like the idea of using "routingGammaMean <= 0._rkind ". In fact, I think both scale and shape parameters need to be positive, not zero.

andywood commented 3 years ago

Hi Hongli, sounds good.

I agree that routingGammaScale needs to be >0. The 'bad values' check would best trap a user trying to enter 0 (perhaps by mistake, perhaps because their optimization wasn't good, etc) to avoid the code crashing on a divide-by-zero error a few lines later.

Btw, gamma is a common statistical distribution, even a fairly common representation of a unit hydrograph - the implementation could probably be referenced to any textbook (or not at all) so we're not tied to a specific paper. But it's fine to leave it as you did.

Cheers, -Andy

On Wed, Jun 16, 2021 at 10:44 PM Hongli Liu @.***> wrote:

Hi Andy, using "tFuture/routingGammaScale" is indeed simplier and correct as well. routingGammaMean is introduced because it appears in Equation 13 of Clark et al. (2007, WRR). In this way, the code stays consistent with the paper. Regarding routingGammaScale, I think mathematically it needs to be positive, cannot be zero. Does this make sense? In summa localParam.txt, its lower and upper limits are 1 and 5000000., respectively.

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/NCAR/summa/pull/470#issuecomment-862919858, or unsubscribe https://github.com/notifications/unsubscribe-auth/ABIKARNEUS35GPAT44Y4LRTTTF4TLANCNFSM462X5VAA .

h294liu commented 3 years ago

Hi repo managers, I'd like to add a bit more explanations and demos for the changes I made. Hope this eases understanding.

The previous code ACTUALLY has an error in calculating the x value of the Gamma distribution. The reason behind is that, when x value is employed in cumprob = gammp(shape, x), the gammp function requires the x being standardized. Given the two-parameter Gamma distribution, the time delay "tFuture" should be standardized by "tFuture/scale". The old summa code has an error in doing the standardization. Here this issue is fixed by removing variable "aLambda" and updating x value for gammp.

In addition, I add two the figures to show the impacts of code changes on routed runoff. This 2-year simulation is based on the Bow at Banff case study. The first figure shows the summa averageRoutedRunoff of the 1st GRU. The second figure shows the mizuRoute IRFroutedRunoff at the basin outlet. You can see that the new code brings down the peak runoff more or less. SUMMA_output mizuRoute_output