ec-jrc · ecCinziaMazzetti · Jan 25, 2023 · Mar 23, 2023 · Mar 23, 2023 · Mar 23, 2023
diff --git a/VERSION b/VERSION
@@ -1 +1 @@
-4.3.1
+4.3.1 MCT
diff --git a/docs/3_step3_preparing-setting-file/index.md b/docs/3_step3_preparing-setting-file/index.md
@@ -378,6 +378,41 @@ These parameters are all related to the [routing of water in the channels](https
 
 
 
+### Diffusive wave routing parameters
+
+The following parameters are related to the [diffusive wave routing](https://ec-jrc.github.io/lisflood-model/3_14_optLISFLOOD_diffusive-wave/) in river channels. The multiplier *CalChanMan3* can be used to fine-tune the diffusive wave propagation when using the Muskingum-Cunge-Todini (MCT) routing, and it can be defined as either a single value or a map. The map *ChannelsMCT* is a Bolean map with the mask of rivers where MCT wave routing must be used. The parameter *ChanGradMaxMCT* defines the maximum riverbed slope for river grid cells using the MCT wave routing. The parameter is provided as a single number and it is recommended to set it to values < 0.001 and > *ChanGradMin*
+```xml
+	<comment>
+	**************************************************************
+	MUSKINGUM-CUNGE-TODINI ROUTING PARAMETERS
+	**************************************************************
+	</comment>
+	<textvar name="CalChanMan3" value="$(PathParams)/params_CalChanMan3">
+	<comment>
+	default: 3.0 [-]
+	Multiplier [-] applied to Channel Manning's n for MCT routing
+	</comment>
+	</textvar>
+	<textvar name="ChannelsMCT" value="$(PathRoot)/maps/chanmct">
+	<comment>
+	Boolean map with value 1 at channel pixels where MCT is
+	used, and 0 at all other pixels
+	</comment>
+	</textvar>
+	<textvar name="ChanGradMaxMCT" value="0.001">
+	<comment>
+	Maximum channel gradient for channels using MCT routing [-] (for MCT wave: slope cannot be 0)
+	Default: 0.001
+	</comment>
+	</textvar>                                                       
+```
+
+- **CalChanMan3** is a multiplier that is applied to the Manning’s roughness map of the [channel system](https://ec-jrc.github.io/lisflood-model/2_16_stdLISFLOOD_channel-routing/) [-] for the grid cells where MCT routing is used
+- **ChannelsMCT** is Bolean mask including the rivers grid cells using the MCT wave routing [-]
+- **ChanGradMaxMCT** is a upper limit for the channel gradient used in the calculation of the MCT wave routing [m/m]
+
+
+
 ### Parameters related to numerics 
 
 This category only contains one parameter at the moment, which can only be a single value. We strongly recommend keeping this parameter at its default value.
@@ -530,8 +565,8 @@ Note that a couple of variables can be [initialized internally](https://ec-jrc.g
 	</comment>                                                          
 	</textvar>                                                          
 	<comment>                                                           
-**************************************************************               
-	The following variables can also be initialized in the model internally. if you want this to happen set them to bogus value of -9999                                                                 
+	**************************************************************               
+	The following variables can also be initialized in the model internally. If you want this to happen set them to bogus value of -9999                                                                 
 	**************************************************************               
 	</comment>                                                          
 	<textvar name="LZInitValue" value="-9999">                      
@@ -571,7 +606,19 @@ Note that a couple of variables can be [initialized internally](https://ec-jrc.g
 	only needed for lakes reservoirs and transmission loss                
 	-9999: use discharge of half bankfull                                 
 	</comment>                                                          
-	</textvar>                                                          
+	</textvar>
+	<textvar name="PrevCmMCTInitValue" value="-9999">
+	<comment>
+	Courant number at previous step for MCT routing
+	Cold start: -9999: use 1
+	</comment>
+	</textvar>
+	<textvar name="PrevDmMCTInitValue" value="-9999">
+	<comment>
+	Reynolds number at previous step for MCT routing
+	Cold start: -9999: use 0
+	</comment>
+	</textvar>		                                                          
 ```
 
 - **WaterDepthInitValue** is the initial amount of water on the soil surface $[mm]$
@@ -604,6 +651,10 @@ Note that a couple of variables can be [initialized internally](https://ec-jrc.g
 
 - **PrevDischarge** is the initial discharge from previous run $[\frac{m^3} {s}]$ used for lakes, reservoirs and transmission loss (only needed if option is on for lakes or reservoirs or transmission loss). Note that PrevDischarge is discharge as an average over the time step (a flux) . A value of **-9999** sets the initial amount of discharge to equivalent of half bankfull.
 
+- **PrevCmMCTInitValue** is the Courant number at the end of the previous step and it is only used for MCT wave routing [-]. A value of **-9999 ** sets the initial value to 1.
+
+- **PrevDmMCTInitValue** is the Reynols number at the end of the previous step and it is only used for MCT wave routing [-]. A value of **-9999 ** sets the initial value to 0.
+
 ```xml
 	<comment>                                                           
 	**************************************************************               

diff --git a/docs/3_step4_preparing-input-files/index.md b/docs/3_step4_preparing-input-files/index.md
@@ -31,17 +31,17 @@ The section [Static Maps](../4_Static-Maps-introduction) provides detailed guide
 +  [land use maps](../4_Static-Maps_land-use/): fraction of forest; fraction of irrigated crops; fraction of rice crops; fraction of inland water; fraction of sealed surfaces; fraction of other land uses.
 +   [land use depending](../4_Static-Maps_land-use-depending/):crop coefficient; crop group number; Manning/s's surface roughness; soil depth.
 +   [soil hydraulic properties](../4_Static-Maps_soil-hydraulic-properties/): saturated hydraulic conductivity; soil water content at saturation; residual soil water content; parameters alpha and lambda of Van Genuchten's equations.
-+   [channel geometry](../4_Static-Maps_channel-geometry/): channels mask; channels side slope; channels length; channels gradient; Manning's rougheness coefficient of the channels; channels bottom width; floodplain width; bankfull channels depth.
++   [channel geometry](../4_Static-Maps_channel-geometry/): channels mask; channels side slope; channels length; channels gradient; Manning's rougheness coefficient of the channels; channels bottom width; floodplain width; bankfull channels depth; MCT diffusive wave routing channels.
 +   [leaf area index](../4_Static-Maps_leaf-area-index/): evolution of vegetation over time (leaf area index) for land covers forest, irrigated areas, others.
 +   [reservoirs and lakes](../4_Static-Maps_reservoirs-lakes/): lake mask; lakes ID points; reservoirs ID points. These maps are required only upon activation of the [lakes module](https://ec-jrc.github.io/lisflood-model/3_02_optLISFLOOD_lakes/) and/or of the [reservoirs module](https://ec-jrc.github.io/lisflood-model/3_03_optLISFLOOD_reservoirs/).
 +  [rice calendar](../4_Static-Maps_rice-calendar/): rice calendar for planting and harvesting seasons. These maps are required only when activating the [rice irrigation module](https://ec-jrc.github.io/lisflood-model/2_17_stdLISFLOOD_irrigation/)
 +   inflow points: locations and IDs of the points in which LISFLOOD adds an inflow hydrograph, as explained [here](https://ec-jrc.github.io/lisflood-model/3_09_optLISFLOOD_inflow-hydrograph/)
 +  water demand maps: domestic, energetic, livestock, industrial water use. These maps represent the time series of spatially distributed values of water demand for domestic, energetic, livestock, and industrial water use. These maps  are required only when activating the [water use module](https://ec-jrc.github.io/lisflood-model/2_18_stdLISFLOOD_water-use/)
 + outlet points: locations and IDs of the points for which LISFLOOD provides the time series of discharge values.
 
-#### Role of "mask" and "channels" maps 
+#### Role of "mask", "channels" ans "channelsMCT" maps 
 
-The mask map (i.e. "area.map") defines the model domain. In order to avoid unexpected results, **it is vital that all maps that are related to topography, land use and soil are defined** (i.e. don't contain a missing value) for each pixel that is "true" (has a Boolean 1 value) on the mask map. The same applies for all meteorological input and the Leaf Area Index maps. Similarly, all pixels that are "true" on the channels map must have some valid (non-missing) value on each of the channel parameter maps. Undefined pixels can lead to unexpected behaviour of the model, output that is full of missing values, loss of mass balance and possibly even model crashes. Some maps needs to have values in a defined range e.g. the gradient map has to be greater than 0.
+The mask map (i.e. "area.map") defines the model domain. In order to avoid unexpected results, **it is vital that all maps that are related to topography, land use and soil are defined** (i.e. don't contain a missing value) for each pixel that is "true" (has a Boolean 1 value) on the mask map. The same applies for all meteorological input and the Leaf Area Index maps. Similarly, all pixels that are "true" on the channels map must have some valid (non-missing) value on each of the channel parameter maps. At the same time, all pixels that have value "true" in the MCT rivers mask must also belong to the "channels" map. Undefined pixels can lead to unexpected behaviour of the model, output that is full of missing values, loss of mass balance and possibly even model crashes. Some maps needs to have values in a defined range e.g. the gradient map has to be greater than 0.
 
 #### Geometrical properties of the computational grid cell
 

diff --git a/docs/4_Static-Maps_channel-geometry/index.md b/docs/4_Static-Maps_channel-geometry/index.md
@@ -2,6 +2,7 @@
 
 In the LISFLOOD model flow through the channel is simulated using the kinematic wave equations. Channel maps describe the sub grid information of the channel geometry, i.e. the length, slope, width and depth of the main channel inside a grid-cell.  <br/>
 + The **channel mask** map is the Boolean field with '1' for all grid-cells with channels and NoData for all grid-cells with no channels.  <br/>
++ The **MCT channel mask** map is the Boolean field with '1' for all river grid-cells using MCT wave routing and '0' for other channel grid cells.  <br/>
 + The **channel side slope** map (referred as 's' in Figure 41) defines the slope of the channel banks. <br/>
 + The **channel length** map is the length of the river in each grid-cell, and the value can exceed grid-size to account for meandering rivers. <br/>
 + The **channel gradient** (or channel slope) map is the average gradient of the main river inside a cell. <br/>
@@ -23,6 +24,7 @@ Channel characteristics, explained above, are shown in the Figure 41 below.  <br
 | Map name | File name;type | Units; range | Description |
 | :---| :--- | :--- | :--- |
 |Channel mask         | chan.nc; <br> Type: Boolean          | Units: -;<br> Range: NoData or 1         |Boolean map that identifies the channel grid-cells |
+|MCT Channel mask         | chanmct.nc; <br> Type: Boolean          | Units: -;<br> Range: [0-1]         |Boolean map that identifies the channel grid-cells using the MCT diffusive river routing |
 |Side slope        |chans.nc; <br> Type: Float32          | Units: m;<br> Range>0          |Channel side slope|
 |Channel length         |chanlength.nc; <br> Type: Float32          | Units: m;<br> Range>0         |Channel length (value can exceed grid size, to account for meandering rivers)|
 |Channel gradient         |changrad.nc; <br> Type: Float32           |Units: m/m;<br> Range: [0-1]          |Channel longitudinal gradient|
@@ -45,6 +47,9 @@ Channel characteristics, explained above, are shown in the Figure 41 below.  <br
 ### Channel mask (chan)
 The channel mask map is used to identify the cells that have channels. The grid-cells that have a channel length (see chanlength map creation below) above zero are assigned to the Boolean field '1', the grid-cells that have a channel length below or equal to zero are assigned with NoData.
 
+### MCT Channel mask (chanmct)
+The  MCT channel mask map is used to identify the cells using the Muskingum-Cunge-Todini diffusive wave routing. The grid-cells that have a riverbed slope < *ChanGradMaxMCT* (default value 0.001) and a set number of upstream grid cells also meeting the same condition are assigned to the mask. All downstream channel pixels of any of the pixels using MCT wave routing are also added to the mask.
+
 ### Side slope (chans)
 The channel side slope map is calculated by dividing the horizontal distance (referred as 'dx' in Figure 42) by vertical distance (referred as 'dy' in Figure 42); here ‘1’ was assigned to all the grid cells, which correspond to a 45° angle of the side slope.
 

diff --git a/docs/4_annex_input-files/index.md b/docs/4_annex_input-files/index.md
@@ -86,6 +86,7 @@ To this end a dedicated pre-processing application has been developed (LISVAP),
 | Sat2 for forest and other                                 | ksat2.map           | U.: $\frac{mm} {day}$ <br> R.:  map>0           | Saturated conductivity layer 2                               |
 | **CHANNEL GEOMETRY**                                      |                     |                                                        |                                                              |
 | Channels                                                  | chan.map            | U.: [-] <br> R.: 0 or 1                                | Map with Boolean 1 for all channel pixels, and Boolean 0 for all other pixels on MaskMap |
+| MCT Channels                                                  | chanmct.map            | U.: [-] <br> R.: 0 or 1                                | Map with Boolean 1 for channel pixels using MCT diffusive wave routing, and Boolean 0 for all other pixels on MaskMap |
 | ChanGrad                                                  | changrad.map        | U.: $\frac{m} {m}$ <br> R.: map > 0  <br> !!!          | Channel gradient                                             |
 | ChanMan                                                   | chanman.map         | U.: $m^{-1/3} s$ <br> R.: map > 0                      | Manning's roughness coefficient for channels                 |
 | ChanLength                                                | chanleng.map        | U.: $m$ <br> R.: map > 0                               | Channel length (can exceed grid size, to account for meandering rivers) |

diff --git a/docs/4_annex_output-files/index.md b/docs/4_annex_output-files/index.md
@@ -9,7 +9,8 @@ LISFLOOD can generate a wide variety of outputs. Output files can be time series
 | Description                                                  | Units            | File name             |
 | ------------------------------------------------------------ | ---------------- | --------------------- |
 | **RATE VARIABLES AT GAUGES**                                 |                  |                       |
-| $^{1,2}$ channel discharge                                   | $\frac{m^3} {s}$ | dis.tss               |
+| $^{1,2}$ average river discharge                                   | $\frac{m^3} {s}$ | dis.tss               |
+| $^{1,2}$ instantaneous rivers discharge                                   | $\frac{m^3} {s}$ | chanq.tss               |
 | **NUMERICAL CHECKS**                                             |                  |                       |
 | $^2$ cumulative mass balance error                           | $m^3$            | mbError.tss           |
 | $^2$ cumulative mass balance error, expressed as mm water slice (average over catchment) | $mm$             | mbErrorMm.tss         |