The **GR4J model **has four parameters to optimise during calibration:

: the production store maximal capacity (mm),*X*1

: the catchment water exchange coefficient (mm/day),*X*2

: the one-day maximal capacity of the routing reservoir (mm),*X*3

: the*X*4*HU1*unit hydrograph time base (days).

We denote by *P* (mm/day) the **rainfall amount** and by *E* (mm/day) the **potential evapotranspiration** (PET). *P* is an estimation of the catchment rainfall and *E* can come from a mean interannual PET curve.

The following equations correspond to the equations integrated over a time step.

The first operation is the neutralization of *P* by *E* to determine a **net rainfall** *Pn* and a **net evapotranspiration** *En*, calculated by:

If *P* > *E*, then *Pn* = *P* – *E* and *En* = 0

If *P* < *E*, then *Pn* = 0 and *En* = *E* – *P*

In case *Pn* is different from zero, a fraction *Ps* of *Pn* goes to the production reservoir and is calculated by:

where *X*1 (mm) and *S* are, respectively, the maximum capacity and the **production store **level.

Otherwise, when *En* is different from zero, a part of evaporation *Es* is removed from the production store. It is given by:

The production store level is updated through:

*S* = *S* – *Es* + *Ps*

A **percolation **called *Perc *coming from the production store is then calculated:

The production store level is then again updated:

*S = S – Perc*

The water quantity *Pr* that finally reaches the **routing part of the model **is:

*Pr = Perc + (Pn – Ps)*

*Pr* is divided into two flow components, 90 % being routed by a **unit hydrograph** *HU*1 and a routing store and 10 % by a unique unit hydrograph *HU*2.

*HU*1 and *HU*2 depend on the same parameter *X*4, the time base of *HU*1 expressed in days.

The **hydrograms ordinates** are calculated from the S curves (the accumulation of the proportion of unit rainfall treated by the hydrogram in function of time), respectively named *SH*1 and *SH*2.

*SH*1 is defined in function of time by:

*For t = 0 *

*For 0 < t < X4 *

*For t **> **X4 *

*SH*2 is defined in function of time by:

*For t = 0 *

*For 0 < t < X4 *

*For X4 < t < 2X4*

*For t > 2X4*

**The ordinates of HU1 and HU2 are then obtained from:**

where *j* is an integer.

For each time step *i*, the outputs Q9 and Q1 of the two hydrograms are calculated with:

with l = int(X4)+1 and m = int(2.X4)+1, with int(.) the integer part.

A **groundwater exchange term** (loss or gain) is calculated with:

with *R* the routing store level, *X*3 the one-day maximal capacity of the store and *X*2 the water exchange coefficient, which is positive in case of a gain, and negative in case of a loss, or null.

The level in the** routing store** is updated by adding the *Q*9 output of the hydrogram *HU*1 and *F*:

*R* = max (0 ; *R* + *Q*9 + *F*)

Then, it empties in an output *Qr* given by:

The level in the store becomes: *R* = *R* – *Qr*

The output *Q*1 of the hydrogram *HU*2 goes through the same exchanges to give the flow component *Qd*:

*Qd = max (0 ; Q1+F) *

The **total streamflow ***Q* is finally given by: *Q* = *Qr* + *Qd*

**To know more**: check our publications.