________ _______ ________ _________
|\ _____\|\ ___ \ |\ ____\ |\___ ___\
\ \ \__/ \ \ __/| \ \ \___|_ \|___ \ \_|
\ \ __\ \ \ \_|/__ \ \_____ \ \ \ \
\ \ \_| \ \ \_|\ \ \|____|\ \ \ \ \
\ \__\ \ \_______\ ____\_\ \ \ \__\
\|__| \|_______| |\_________\ \|__|
\|_________|
S P A T I A L L Y - D I S T R I B U T E D
H Y D R O L O G I C A L - M O D E L
FeST is a spatially distributed physically based hydrological model. FeST is the acronym of “flash–Flood Event–based Spatially distributed rainfall–runoff Transformation”. It was initially developed by Mancini (1990), as a model oriented to the simulation of rainfall-runoff transformation of single flood events. Later the FeST model was merged with the soil water balance scheme from TOPLATS model (Famiglietti and Wood, 1994), transforming it into a continuous model (Montaldo et al., 2007). Then the FeST code was redesigned and rewritten from scratch while keeping the basic assumptions of the previous release (Rabuffetti et al., 2008). In 2011 the FeST was upgraded with a routine to solve the system of water mass and energy balance in order to better simulate the actual evapotranspiration and interface the model to remotely sensed data (Corbari et al., 2011; Corbari & Mancini, 2014). At the same year, 2011, a new module for simulating groundwater flux and river-groundwater interaction was developed and implemented in the FeST (Ravazzani et al., 2011). In 2013 a new version of the code was released built on top of the MOSAICO library (Ravazzani, 2013). In 2014 the FeST was upgraded with a module for glaciers modelling (Boscarello et al., 2014). In 2021 a forest growth component was implemented in the FeST (Feki et al., 2021).
Boscarello, L., Ravazzani, G., Rabuffetti, D., & Mancini, M.. (2014). Integrating glaciers raster-based modelling in large catchments hydrological balance: the Rhone case study. Hydrological processes, 28, 496–508.
Corbari, C., & Mancini, M. (2014). Calibration and validation of a distributed energy–water balance model using satellite data of land surface temperature and ground discharge measurements. Journal of hydrometeorology, 15(1), 376-392.
Corbari, C., Ravazzani, G., & Mancini, M. (2011). A distributed thermodynamic model for energy and mass balance computation: FEST–EWB. Hydrological Processes, 25(9), 1443-1452.
Famiglietti, J.S., & Wood, E.F. (1994). Multiscale modeling of spatially variable water and energy balanceprocess. Water Resour. Res., 30, 3061–3078, https://doi.org/10.1029/94WR01498.
Feki, M., Ravazzani, G., Ceppi, A., Pellicone, G., & Caloiero, T.. (2021). Integration of forest growth component in the fest-wb distributed hydrological model: the bonis catchment case study. Forests, 12(12).
Mancini, M.. (1990). La modellazione distribuita della risposta idrologica: effetti della variabilità spaziale e della scala di rappresentazione del fenomeno dell'assorbimento. PhD thesis. Politecnico di Milano, Istituto di idraulica
Montaldo, N., Ravazzani, G., & Mancini, M.. (2007). On the prediction of the toce alpine basin floods with distributed hydrologic models. Hydrological processes, 21, 608–621
Rabuffetti, D., Ravazzani, G., Corbari, C., & Mancini, M.. (2008). Verification of operational quantitative discharge forecast (QDF) for a regional warning system – the AMPHORE case studies in the upper Po river. Natural hazards and earth system sciences, 8, 161–173.
Ravazzani, G.. (2013). Mosaico, a library for raster based hydrological applications. Computers & geosciences, 51, 1–6.
Ravazzani, G., Rametta, D., & Mancini, M.. (2011). Macroscopic cellular automata for groundwater modelling: a first approach. Environmental modelling & software, 26(5), 634–643.
license: GNU GPL http://www.gnu.org/licenses/
output ragguagli falda
| Type | Attributes | Name | Initial | |||
|---|---|---|---|---|---|---|
| type(grid_real) | :: | Ge | ||||
| type(grid_real) | :: | Hse | ||||
| type(grid_real) | :: | Rnetta | ||||
| type(grid_real) | :: | Ta_prec | ||||
| type(grid_real) | :: | Ts | ||||
| type(grid_real) | :: | Ts_prec | ||||
| type(grid_real) | :: | Velocita | ||||
| type(grid_real) | :: | Xle | ||||
| character(len=300) | :: | arg |
command line arguments |
|||
| character(len=300) | :: | cron | ||||
| type(CronTab) | :: | cronHotstart | ||||
| integer(kind=short) | :: | dtArray(4) | ||||
| integer(kind=short) | :: | dtCalc | ||||
| integer | :: | dtCalcSediment | ||||
| integer | :: | dtCalcSedimentRouting | ||||
| integer(kind=short) | :: | dtOutBasinAverage | ||||
| integer(kind=short) | :: | dtOutCanopy | ||||
| integer(kind=short) | :: | dtOutGroundwater | ||||
| integer(kind=short) | :: | dtOutIce | ||||
| integer(kind=short) | :: | dtOutIrrigation | ||||
| integer(kind=short) | :: | dtOutMeteo | ||||
| integer(kind=short) | :: | dtOutPlants | ||||
| integer | :: | dtOutSediment | ||||
| integer | :: | dtOutSedimentRouting | ||||
| integer(kind=short) | :: | dtOutSnow | ||||
| integer | :: | dtOutSoilBalance | ||||
| logical | :: | fatta_previ | = | .FALSE. | ||
| type(IniList) | :: | festIni | ||||
| logical | :: | hotstart | ||||
| integer | :: | i | ||||
| logical | :: | iniFound |
says if configuration file has been passed as command argument |
|||
| integer | :: | j | ||||
| integer | :: | k | ||||
| character(len=1000) | :: | pathHotstart | ||||
| character(len=1000) | :: | path_out | ||||
| type(grid_real) | :: | raineff | ||||
| logical | :: | saveLast | ||||
| type(grid_real) | :: | storageChannelSurf | ||||
| character(len=1000) | :: | string | ||||
| type(DateTime) | :: | time | ||||
| type(DateTime) | :: | timeNewBasinAverage | ||||
| type(DateTime) | :: | timeNewCanopyInterception | ||||
| type(DateTime) | :: | timeNewDischargeRouting | ||||
| type(DateTime) | :: | timeNewGroundwater | ||||
| type(DateTime) | :: | timeNewIce | ||||
| type(DateTime) | :: | timeNewIrrigation | ||||
| type(DateTime) | :: | timeNewMeteo | ||||
| type(DateTime) | :: | timeNewPlants | ||||
| type(DateTime) | :: | timeNewRasterExport | ||||
| type(DateTime) | :: | timeNewSnow | ||||
| type(DateTime) | :: | timeNewSoilBalance | ||||
| type(DateTime) | :: | timeOutCanopy | ||||
| type(DateTime) | :: | timeOutGroundwater | ||||
| type(DateTime) | :: | timeOutIce | ||||
| type(DateTime) | :: | timeOutIrrigation | ||||
| type(DateTime) | :: | timeOutMeteo | ||||
| type(DateTime) | :: | timeOutPlants | ||||
| type(DateTime) | :: | timeOutSediment | ||||
| type(DateTime) | :: | timeOutSedimentRouting | ||||
| type(DateTime) | :: | timeOutSnow | ||||
| type(DateTime) | :: | timeOutSoilBalance | ||||
| type(DateTime) | :: | timeStart | ||||
| type(DateTime) | :: | timeStop | ||||
| type(DateTime) | :: | timeUpdateSediment | ||||
| type(DateTime) | :: | timeUpdateSedimentRouting |
Save state variables for starting a new simulation from the current condition
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| type(DateTime), | optional | :: | time |
PROGRAM Fest
USE IniLib, ONLY: &
!Imported derived types:
IniList, &
!Imported routines:
IniOpen, IniClose, &
IniReadInt, IniReadReal, IniReadDouble, &
IniReadString, SectionIsPresent, &
KeyIsPresent
USE HydroNetwork, ONLY : &
!imported definitions:
!NetworkSurfaceFlow,
NetworkGroundwater
USE Reservoirs, ONLY : &
!imported definitions:
Reservoir, &
!Imported variables:
dtReservoir
USE SpatialAverage, ONLY : &
!imported routines:
InitSpatialAverageMeteo, &
InitSpatialAverageBalance, &
InitSpatialAverageSnow, &
InitSpatialAverageGlaciers, &
InitSpatialAverageSediment, &
InitSpatialAverageCanopy, &
InitSpatialAveragePlants, &
ComputeSpatialAverageMeteo,&
ComputeSpatialAverageBalance, &
ComputeSpatialAverageSnow, &
ComputeSpatialAverageGlaciers, &
ComputeSpatialAverageSediment, &
ComputeSpatialAverageCanopy, &
ComputeSpatialAveragePlants, &
ExportSpatialAverageMeteo, &
ExportSpatialAverageBalance, &
ExportSpatialAverageSnow, &
ExportSpatialAverageGlaciers, &
ExportSpatialAverageSediment, &
ExportSpatialAverageCanopy, &
ExportSpatialAveragePlants, &
timeSpatialAverageMeteo, &
timeSpatialAverageBalance,&
timeSpatialAverageSnow, &
timeSpatialAverageIce, &
timeSpatialAverageSediment, &
timeSpatialAverageCanopy, &
timeSpatialAveragePlants
USE BasinAverage, ONLY : &
!Imported routines:
InitBasinAverage, &
ReadPointFileBasinAverage, &
ExportBasinAverage, &
!Imported variables:
dtBasinAverage
USE RasterExport, ONLY : &
!Imported routines:
InitRasterExport, &
ExportMaps
USE Snow, ONLY : &
!imported routines:
SnowInit, SnowUpdate, &
SnowPointInit, &
!Imported variables:
dtSnow, meltCoefficient, &
swe, waterInSnow, snowMelt, &
rainfallRate
USE Glacier, ONLY : &
!imported routines:
GlacierInit, GlacierPointInit, &
GlacierUpdate, &
!imported variables:
dtIce, icewe, waterInIce, iceMelt
USE SoilBalance, ONLY : &
!imported routines
InitSoilBalance, SolveSoilBalance, &
!imported variables:
soilMoisture, soilMoistureRZ, &
soilMoistureTZ, soilSatRZ, &
soilSatTZ, sEff, rainFlag, &
runoff, infilt, percolation, caprise, &
soilDepth, wiltingPoint, fieldCapacity, &
psdi, ksat, et, pet, dtSoilBalance, &
QinSoilSub, soilSat, balanceError, &
interstormDuration, deltaSoilMoisture
USE Infiltration, ONLY : &
!imported variables:
infiltrationModel, SCS_CN, &
PHILIPEQ, GREEN_AMPT, &
sEff, raincum, cuminf
USE Groundwater, ONLY : &
!Imported routines:
GroundwaterInit, &
GroundwaterUpdate, &
GroundwaterPointInit, &
GroundwaterOutputInit, &
GroundwaterOutput, &
GroundwaterRiverUpdate, &
!imported variables:
dtGroundwater, &
vadoseDepth, &
riverGroundwaterInteract, &
riverToGroundwater, &
groundwaterToRiver
USE DischargeRouting, ONLY : &
!imported routines:
InitDischargeRouting, DischargeRoute, &
DischargePointInit, &
!imported variables:
Qin, Qout, Pin, Pout, Plat, &
dtDischargeRouting, &
topWidth, waterDepth
USE Irrigation, ONLY : &
!Imported routines:
IrrigationConfig, IrrigationUpdate, &
!Imported variables:
dtIrrigation, &
irrigationFlux
USE Sediment, ONLY: &
!Imported routines:
SedimentInit, InterrillDetachmentRate, &
SedimentRouting, &
!imported variables:
interrillErosion, routeSediment, QoutSed, deltaSed
USE DomainProperties, ONLY: &
!imported variables:
mask, albedoGround, albedo, &
albedo_loaded, &
!imported routines:
DomainInit
USE MorphologicalProperties, ONLY: &
!imported variables:
dem, dem_loaded, &
flowDirection, &
streamNetwork, &
!imported routines:
MorphologyInit
USE Meteo, ONLY: &
!Imported routines:
MeteoInit, MeteoRead , &
MeteoPointInit, &
!Imported variables:
dtMeteo
USE Precipitation, ONLY: &
!imported variables:
precipitationRate
USE PrecipitationDaily, ONLY: &
!imported variables:
precipitationRateDaily
USE AirTemperature, ONLY: &
!imported variables:
temperatureAir
USE AirTemperatureDailyMean, ONLY: &
!imported variables:
temperatureAirDailyMean
USE AirTemperatureDailyMax, ONLY: &
!imported variables:
temperatureAirDailyMax
USE AirTemperatureDailyMin, ONLY: &
!imported variables:
temperatureAirDailyMin
USE SolarRadiation, ONLY: &
!imported variables:
radiation, netRadiation
USE AirRelativeHumidity, ONLY: &
!imported variables:
relHumidityAir
USE WindFlux, ONLY: &
!imported variables:
windSpeed
USE Plants, ONLY : &
!imported routines:
PlantsConfig, PlantsGrow, &
PlantsParameterUpdate, &
!imported variables:
dtPlants, lai, fvcover, &
gpp, npp, carbonroot, &
carbonstem, carbonleaf, dbh, &
plantsHeight, density, stemyield, &
rsMinLoaded, updatePlantsParameters
USE PlantsInterception, ONLY: &
!imported variables:
dtCanopyInterception, canopyStorage, &
canopyPT, &
!imported routines:
InterceptionInit, Throughfall, &
AdjustPT
USE DataTypeSizes, ONLY: &
!Imported parameters:
short, long, float
USE LogLib, ONLY: &
!Imported routines:
LogInit, LogStop, Catch, &
!Imported variables:
verbose
USE GeoLib, ONLY: &
!Imported routines:
GeoInit
USE Chronos, ONLY: &
!Imported routines:
GetMinute, GetHour, GetMonth, &
GetDay, GetDayOfWeek, &
!Imported definitions:
DateTime, &
!Imported variables:
timeString, &
!Imported operators:
ASSIGNMENT (=), &
OPERATOR (+), &
OPERATOR (-), &
OPERATOR (==), &
OPERATOR (>), &
OPERATOR (<), &
OPERATOR (<=)
USE GridOperations, ONLY : &
!Imported operators:
ASSIGNMENT (=)
USE GridLib, ONLY : &
!imported definition:
grid_real, &
!imported routines:
NewGrid, ExportGrid, &
!imported parameters:
ESRI_BINARY
USE StringManipulation, ONLY : &
!imported routines
StringCompact, StringTokenize, &
StringToLong
USE CronSchedule, ONLY : &
!Imported types:
CronTab, &
!Imported routines:
CronParseString, &
CronIsTime
IMPLICIT NONE
!configuration
TYPE (IniList) :: festIni
CHARACTER (len = 300) :: arg !! command line arguments
CHARACTER (len = 1000) :: string
LOGICAL :: iniFound !!says if configuration file has been passed as command argument
!simulation time
TYPE (DateTime) :: time !current step time
TYPE (DateTime) :: timeStart !start time
TYPE (DateTime) :: timeStop !stop time
!meteo
TYPE (DateTime) :: timeNewMeteo
TYPE (DateTime) :: timeOutMeteo
INTEGER (KIND = short) :: dtOutMeteo
!irrigation
TYPE (DateTime) :: timeNewIrrigation
TYPE (DateTime) :: timeOutIrrigation
INTEGER (KIND = short) :: dtOutIrrigation
!snow
TYPE (DateTime) :: timeNewSnow
TYPE (DateTime) :: timeOutSnow
INTEGER (KIND = short) :: dtOutSnow
!glacier
TYPE(DateTime) :: timeNewIce
TYPE(DateTime) :: timeOutIce
INTEGER (KIND = short) :: dtOutIce
!canopy interception
TYPE(DateTime) :: timeNewCanopyInterception
TYPE(grid_real) :: raineff ! throughfall[m/s]
INTEGER (KIND = short) :: dtOutCanopy
TYPE (DateTime) :: timeOutCanopy
!forest dynamic
TYPE(DateTime) :: timeNewPlants
INTEGER (KIND = short) :: dtOutPlants
TYPE (DateTime) :: timeOutPlants
!soil balance
TYPE(DateTime) timeNewSoilBalance
!bilancio energetico
TYPE (grid_real)Ts !temperatura supeficiale
TYPE (grid_real)Rnetta
TYPE (grid_real)Xle
TYPE (grid_real)Hse
TYPE (grid_real)Ge
TYPE (grid_real)Ta_prec !temperatura aria istante precedente
TYPE (grid_real)Ts_prec !temperatura suolo istante precedente
!output
TYPE(DateTime) :: timeOutSoilBalance
INTEGER :: dtOutSoilBalance
! groundwater
INTEGER (KIND = short) :: dtOutGroundwater
TYPE(DateTime) :: timeNewGroundwater
TYPE(DateTime) :: timeOutGroundwater
!sediment
INTEGER :: dtCalcSediment
TYPE(DateTime) :: timeUpdateSediment
INTEGER :: dtOutSediment
TYPE(DateTime) :: timeOutSediment
!sediment routing
INTEGER :: dtCalcSedimentRouting
TYPE(DateTime) :: timeUpdateSedimentRouting
INTEGER :: dtOutSedimentRouting
TYPE(DateTime) :: timeOutSedimentRouting
! TYPE(rete_osservativa) :: xsOutSedimentRouting
! discharge routing
TYPE (DateTime) :: timeNewDischargeRouting
TYPE(grid_real)Velocita
TYPE (grid_real) storageChannelSurf
!basin average
TYPE (DateTime) :: timeNewBasinAverage
INTEGER (KIND = short) :: dtOutBasinAverage
!raster export
TYPE (DateTime) :: timeNewRasterExport
!hotstart
CHARACTER (len=1000) :: pathHotstart
LOGICAL :: saveLast
LOGICAL :: hotstart
CHARACTER (LEN = 300) :: cron
TYPE (CronTab) :: cronHotstart
!general variables
INTEGER (KIND = short) :: dtCalc !calculation time step
INTEGER (KIND = short) :: dtArray (4)
CHARACTER (len = 1000) :: path_out
INTEGER :: k,i,j
LOGICAL :: fatta_previ = .FALSE.
!--------------------------end of declaration----------------------------------
! splash screen
IF ( verbose ) THEN
CALL Splash ()
END IF
!------------------------------------------------------------------------------
! initialize log
!------------------------------------------------------------------------------
CALL LogInit
!------------------------------------------------------------------------------
! initialize cartographic engine
!------------------------------------------------------------------------------
CALL GeoInit ('GeoLib.ini')
!------------------------------------------------------------------------------
! check command line arguments to read configuration file name
!------------------------------------------------------------------------------
i = 1
iniFound = .FALSE.
DO WHILE ( .not. (arg == '') )
CALL Getarg ( i, arg )
SELECT CASE (arg)
CASE ( '-inifile' )
i = i + 1
CALL Getarg ( i , string )
iniFound = .TRUE.
CASE DEFAULT
!case unknown
END SELECT
i = i + 1
END DO
!------------------------------------------------------------------------------
! read name of configuration file if not yet defined
!------------------------------------------------------------------------------
IF ( .NOT. iniFound) THEN
WRITE(*,*) 'Insert the name of configuration file and press enter'
READ(*,*) string
END IF
CALL Catch ('info', 'Fest', &
'start configuration from file: ', argument = string )
!------------------------------------------------------------------------------
! read configuration file
!------------------------------------------------------------------------------
IF ( verbose ) THEN
WRITE(*,*)
WRITE(*,*) 'Loading configuration...'
END IF
CALL IniOpen (string, festIni)
!------------------------------------------------------------------------------
! configure time
!------------------------------------------------------------------------------
CALL Catch ('info', 'Fest', 'time configuration' )
IF ( verbose ) THEN
WRITE(*,*)
WRITE(*,*) '...configuring time...'
END IF
!check whether time is present
IF ( .NOT. SectionIsPresent ('time', festIni) ) THEN !section is mandatory
CALL Catch ('error', 'Fest', 'time section missing in configuration file' )
END IF
!start time
IF (KeyIsPresent('start', festIni, section = 'time')) THEN
timeString = IniReadString ('start', festIni, section = 'time')
ELSE
CALL Catch ('error', 'Fest', 'start missing in time' )
END IF
timeStart = timeString
CALL Catch ('info', 'Fest', 'start time: ', argument = timeString )
!stop time
IF (KeyIsPresent('stop', festIni, section = 'time')) THEN
timeString = IniReadString ('stop', festIni, section = 'time')
ELSE
CALL Catch ('error', 'Fest', 'stop missing in time' )
END IF
timeStop = timeString
CALL Catch ('info', 'Fest', 'stop time: ', argument = timeString )
!------------------------------------------------------------------------------
! result folder
!------------------------------------------------------------------------------
CALL Catch ('info', 'Fest', 'result folder configuration' )
IF ( verbose ) THEN
WRITE(*,*)
WRITE(*,*) '...setting folder for results...'
END IF
!check whether section is present
IF ( .NOT. SectionIsPresent ('result', festIni) ) THEN !section is mandatory
CALL Catch ('error', 'Fest', 'result section missing in configuration file' )
END IF
!folder
IF (KeyIsPresent('folder', festIni, section = 'result')) THEN
path_out = IniReadString ('folder', festIni, section = 'result')
ELSE
CALL Catch ('error', 'Fest', 'folder missing in result' )
END IF
CALL Catch ('info', 'Fest', 'result folder: ', argument = TRIM(path_out) )
!------------------------------------------------------------------------------
! configure simulation domain
!------------------------------------------------------------------------------
CALL Catch ('info', 'Fest', 'simulation domain configuration' )
IF ( verbose ) THEN
WRITE(*,*)
WRITE(*,*) '...configuring simulation domain...'
END IF
!check whether section is present
IF ( .NOT. SectionIsPresent ('domain', festIni) ) THEN !section is mandatory
CALL Catch ('error', 'Fest', 'domain section missing in configuration file' )
END IF
!find configuration file
IF (KeyIsPresent('conf-file', festIni, section = 'domain')) THEN
string = IniReadString ('conf-file', festIni, section = 'domain')
ELSE
CALL Catch ('error', 'Fest', 'conf-file missing in domain' )
END IF
!read data
CALL DomainInit (string)
!------------------------------------------------------------------------------
! morphological properties
!------------------------------------------------------------------------------
IF ( verbose ) THEN
WRITE(*,*)
WRITE(*,*) '...configuring morphological properties...'
END IF
!check whether section is present. Optional section
IF ( .NOT. SectionIsPresent ('morphology', festIni) ) THEN !section is mandatory
CALL Catch ('info', 'Fest', 'morphological properties configuration' )
END IF
!find configuration file
IF (KeyIsPresent('conf-file', festIni, section = 'morphology')) THEN
string = IniReadString ('conf-file', festIni, section = 'morphology')
!read data
CALL MorphologyInit (string, mask)
ELSE
CALL Catch ('error', 'Fest', 'conf-file missing in morphology' )
END IF
!------------------------------------------------------------------------------
! configure meteo input
!------------------------------------------------------------------------------
CALL Catch ('info', 'Fest', 'meteo input configuration' )
IF ( verbose ) THEN
WRITE (*,*)
WRITE (*,*) '...configuring meteorological forcings...'
END IF
!check whether section is present
IF ( .NOT. SectionIsPresent ('meteo', festIni) ) THEN !section is mandatory
CALL Catch ('error', 'Fest', 'meteo section missing in configuration file')
END IF
!time step
IF (KeyIsPresent('dt', festIni, section = 'meteo')) THEN
dtMeteo = IniReadInt ('dt', festIni, section = 'meteo')
ELSE
CALL Catch ('error', 'Fest', 'dt missing in meteo' )
END IF
!find configuration file
IF (KeyIsPresent('conf-file', festIni, section = 'meteo')) THEN
string = IniReadString ('conf-file', festIni, section = 'meteo')
ELSE
CALL Catch ('error', 'Fest', 'conf-file missing in meteo' )
END IF
!read data
CALL MeteoInit (string, timeStart, mask, dem, dem_loaded, albedo_loaded)
!spatial data out dt
IF (KeyIsPresent('dt-out-spatial', festIni, section = 'meteo')) THEN
dtOutMeteo = IniReadInt ('dt-out-spatial', festIni, section = 'meteo')
IF ( dtOutMeteo > 0 ) THEN
timeOutMeteo = timeStart
ELSE
timeOutMeteo = timeStart + (-1)
END IF
ELSE
CALL Catch ('warning', 'Fest', 'dt-out-spatial missing in meteo' )
dtOutMeteo = - 1
timeOutMeteo = timeStart + (-1)
END IF
!point data out
IF (KeyIsPresent('out-point-file', festIni, section = 'meteo')) THEN
string = IniReadString ('out-point-file', festIni, section = 'meteo')
CALL MeteoPointInit (string, path_out, timeStart)
END IF
timeNewMeteo = timeStart
!------------------------------------------------------------------------------
! irrigation
!------------------------------------------------------------------------------
CALL Catch ('info', 'Fest', 'irrigation configuration' )
!check whether section is present
IF ( SectionIsPresent ('irrigation', festIni) ) THEN !section is optional
IF ( verbose ) THEN
WRITE(*,*)
WRITE(*,*) '...configuring irrigation...'
END IF
!time step
IF (KeyIsPresent('dt', festIni, section = 'irrigation')) THEN
dtIrrigation = IniReadInt ('dt', festIni, section = 'irrigation')
ELSE
CALL Catch ('error', 'Fest', 'dt missing in irrigation' )
END IF
!dt output
IF (KeyIsPresent('dt-out', festIni, section = 'irrigation')) THEN
dtOutIrrigation = IniReadInt ('dt-out', festIni, section = 'irrigation')
IF ( dtOutIrrigation > 0 ) THEN
timeOutIrrigation = timeStart
ELSE
timeOutIrrigation = timeStart + (-1)
END IF
ELSE
CALL Catch ('error', 'Fest', 'dt-out missing in irrigation' )
END IF
!configuration file
IF ( dtIrrigation > 0 ) THEN
IF (KeyIsPresent('conf-file', festIni, section = 'irrigation')) THEN
string = IniReadString ('conf-file', festIni, section = 'irrigation')
CALL IrrigationConfig (string, path_out, mask, dtOutIrrigation )
ELSE
CALL Catch ('error', 'Fest', 'conf-file missing in irrigation' )
END IF
timeNewIrrigation = timeStart
END IF
ELSE
dtIrrigation = 0
END IF
!------------------------------------------------------------------------------
! canopy interception
!------------------------------------------------------------------------------
!_________________________________________
! da implementare che se c'è neve non calcolo intercettazione BO! serve?
!________________________
CALL Catch ('info', 'Fest', 'canopy interception configuration' )
!raineff is used even when canopy interception is not simulated
CALL NewGrid (raineff, mask, 0.)
!check whether section is present
IF ( SectionIsPresent ('canopy-interception', festIni) ) THEN !section is optional
IF ( verbose ) THEN
WRITE(*,*)
WRITE(*,*) '...configuring canopy interception...'
END IF
!dt
IF (KeyIsPresent('dt', festIni, section = 'canopy-interception')) THEN
dtCanopyInterception = IniReadInt ('dt', festIni, section = 'canopy-interception')
ELSE
CALL Catch ('error', 'Fest', 'dt missing in canopy-interception' )
END IF
IF ( dtCanopyInterception > 0 ) THEN
!configuration file
IF (KeyIsPresent('conf-file', festIni, section = 'canopy-interception')) THEN
string = IniReadString ('conf-file', festIni, section = 'canopy-interception')
CALL InterceptionInit ( string, mask )
ELSE
CALL Catch ('error', 'Fest', 'conf-file missing in canopy-interception' )
END IF
!spatial data out dt
IF (KeyIsPresent('dt-out-spatial', festIni, section = 'canopy-interception')) THEN
dtOutCanopy = IniReadInt ('dt-out-spatial', festIni, section = 'canopy-interception')
ELSE
CALL Catch ('error', 'Fest', 'dt-out-spatial missing in canopy-interception' )
END IF
timeNewCanopyInterception = timeStart
IF (dtOutCanopy > 0) THEN
timeOutCanopy = timeStart
ELSE
timeOutCanopy = timeStart + (-1)
END IF
ELSE
dtOutCanopy = 0
timeNewCanopyInterception = timeStart + (-1)
timeOutCanopy = timeStart + (-1)
END IF
ELSE
dtCanopyInterception = 0
dtOutCanopy = 0
timeNewCanopyInterception = timeStart + (-1)
timeOutCanopy = timeStart + (-1)
END IF
!------------------------------------------------------------------------------
! Plants dynamic simulation
!------------------------------------------------------------------------------
CALL Catch ('info', 'Fest', 'plants configuration' )
!check whether section is present
IF ( SectionIsPresent ('plants', festIni) ) THEN !section is optional
IF ( verbose ) THEN
WRITE(*,*)
WRITE(*,*) '...configuring plants simulation...'
END IF
!dt
IF (KeyIsPresent('dt', festIni, section = 'plants')) THEN
dtPlants = IniReadInt ('dt', festIni, section = 'plants')
ELSE
CALL Catch ('error', 'Fest', 'dt missing in plants' )
END IF
!configuration file
IF (KeyIsPresent('conf-file', festIni, section = 'plants')) THEN
string = IniReadString ('conf-file', festIni, section = 'plants')
CALL PlantsConfig (string, mask, timeStart, timeStop)
ELSE
CALL Catch ('error', 'Fest', 'conf-file missing in plants' )
END IF
!spatial data out dt
IF (KeyIsPresent('dt-out-spatial', festIni, section = 'plants')) THEN
dtOutPlants = IniReadInt ('dt-out-spatial', festIni, section = 'plants')
ELSE
CALL Catch ('warning', 'Fest', 'dt-out-spatial missing in plants' )
dtOutPlants = 0
END IF
IF (dtOutPlants > 0) THEN
timeOutPlants = timeStart
END IF
IF ( dtPlants > 0 ) THEN
timeNewPlants = timeStart
ELSE
dtOutPlants = 0
timeNewPlants = timeStart + (-1)
timeOutPlants = timeStart + (-1)
END IF
ELSE
dtPlants = 0
dtOutPlants = 0
timeNewPlants = timeStart + (-1)
timeOutPlants = timeStart + (-1)
END IF
!------------------------------------------------------------------------------
! snow accumulation and melting
!------------------------------------------------------------------------------
CALL Catch ('info', 'Fest', 'snow configuration' )
!rainfallRate is used even when snow is not simulated
CALL NewGrid (rainfallRate, mask, 0.)
!check whether section is present
IF ( SectionIsPresent ('snow', festIni) ) THEN !section is optional
IF ( verbose ) THEN
WRITE(*,*)
WRITE(*,*) '...configuring snow accumulation and melting...'
END IF
!dt
IF (KeyIsPresent('dt', festIni, section = 'snow')) THEN
dtSnow = IniReadInt ('dt', festIni, section = 'snow')
ELSE
CALL Catch ('error', 'Fest', 'dt missing in snow' )
END IF
IF (dtSnow > 0) THEN
!check dt
IF ( dtSnow /= dtMeteo ) THEN
dtSnow=dtMeteo
WRITE(*,*)'WARNING: dtSnow \E8 posto uguale a dt_acquisizione_meteo'
END IF
timeNewSnow = timeStart
!configuration file
IF (KeyIsPresent('conf-file', festIni, section = 'snow')) THEN
string = IniReadString ('conf-file', festIni, section = 'snow')
CALL SnowInit (string, mask, timeStart )
ELSE
CALL Catch ('error', 'Fest', 'conf-file missing in snow' )
END IF
!spatial data out dt
IF (KeyIsPresent('dt-out-spatial', festIni, section = 'snow')) THEN
dtOutSnow = IniReadInt ('dt-out-spatial', festIni, section = 'snow')
ELSE
CALL Catch ('error', 'Fest', 'dt-out-spatial missing in snow' )
END IF
IF ( dtOutSnow > 0 ) THEN
timeOutSnow = timeStart
END IF
!point data out
IF (KeyIsPresent('out-point-file', festIni, section = 'snow')) THEN
string = IniReadString ('out-point-file', festIni, section = 'snow')
CALL SnowPointInit ( string, path_out, timeStart )
END IF
ELSE
dtOutSnow = 0
timenEWSnow = timeStart + (-1)
timeOutSnow = timeStart + (-1)
END IF
ELSE
dtSnow = 0
dtOutSnow = 0
timenEWSnow = timeStart + (-1)
timeOutSnow = timeStart + (-1)
END IF
!------------------------------------------------------------------------------
! glacier
!------------------------------------------------------------------------------
CALL Catch ('info', 'Fest', 'glacer configuration' )
!check whether section is present
IF ( SectionIsPresent ('glacier', festIni) ) THEN !section is optional
IF ( verbose ) THEN
WRITE(*,*)
WRITE(*,*) '...configuring glacier accumulation and ablation...'
END IF
!dt
IF (KeyIsPresent('dt', festIni, section = 'glacier')) THEN
dtIce = IniReadInt ('dt', festIni, section = 'glacier')
ELSE
CALL Catch ('error', 'Fest', 'dt missing in glacier' )
END IF
IF (dtIce > 0) THEN
!check dt
IF ( dtIce /= dtMeteo ) THEN
dtIce = dtMeteo
CALL Catch ('info', 'Fest', 'dtIce set equal to dtMeteo' )
END IF
timeNewIce = timeStart
!configuration file
IF (KeyIsPresent('conf-file', festIni, section = 'glacier')) THEN
string = IniReadString ('conf-file', festIni, section = 'glacier')
CALL GlacierInit (string, mask, timeStart )
ELSE
CALL Catch ('error', 'Fest', 'conf-file missing in glacier' )
END IF
!spatial data out dt
IF (KeyIsPresent('dt-out-spatial', festIni, section = 'glacier')) THEN
dtOutIce = IniReadInt ('dt-out-spatial', festIni, section = 'glacier')
ELSE
CALL Catch ('error', 'Fest', 'dt-out-spatial missing in glacier' )
END IF
IF ( dtOutIce > 0 ) THEN
timeOutIce = timeStart
END IF
!point data out
IF (KeyIsPresent('out-point-file', festIni, section = 'glacier')) THEN
string = IniReadString ('out-point-file', festIni, section = 'glacier')
CALL GlacierPointInit ( string, path_out, timeStart )
END IF
ELSE
dtOutIce = 0
timeNewIce = timeStart + (-1)
timeOutIce = timeStart + (-1)
END IF
ELSE
dtIce = 0
dtOutIce = 0
timeNewIce = timeStart + (-1)
timeOutIce = timeStart + (-1)
END IF
!------------------------------------------------------------------------------
! Soil balance
!------------------------------------------------------------------------------
CALL Catch ('info', 'Fest', 'soil balance configuration' )
!check whether section is present
IF ( SectionIsPresent ('soil-balance', festIni) ) THEN !section is optional
IF ( verbose ) THEN
WRITE(*,*)
WRITE(*,*) '...configuring soil water balance...'
END IF
!dt
IF (KeyIsPresent('dt', festIni, section = 'soil-balance')) THEN
dtSoilBalance = IniReadInt ('dt', festIni, section = 'soil-balance')
ELSE
CALL Catch ('error', 'Fest', 'dt missing in soil-balance' )
END IF
IF ( dtSoilBalance > 0 ) THEN
timeNewSoilBalance = timeStart
!configuration file
IF (KeyIsPresent('conf-file', festIni, section = 'soil-balance')) THEN
string = IniReadString ('conf-file', festIni, section = 'soil-balance')
CALL InitSoilBalance (string, flowDirection, timeStart)
ELSE
CALL Catch ('error', 'Fest', 'conf-file missing in soil-balance' )
END IF
!spatial data out dt
IF (KeyIsPresent('dt-out-spatial', festIni, section = 'soil-balance')) THEN
dtOutSoilBalance = IniReadInt ('dt-out-spatial', festIni, &
section = 'soil-balance')
ELSE
dtOutSoilBalance = 0
END IF
IF (dtOutSoilBalance > 0) THEN
timeOutSoilBalance = timeStart
END IF
ELSE
dtOutSoilBalance = 0
timeNewSoilBalance = timeStart + (-1)
timeOutSoilBalance = timeStart + (-1)
END IF
ELSE
dtSoilBalance = 0
dtOutSoilBalance = 0
timeNewSoilBalance = timeStart + (-1)
timeOutSoilBalance = timeStart + (-1)
END IF
!------------------------------------------------------------------------------
! groundwater
!------------------------------------------------------------------------------
CALL Catch ('info', 'Fest', 'groundwater configuration' )
!check whether section is present
IF ( SectionIsPresent ('groundwater', festIni) ) THEN !section is optional
IF ( verbose ) THEN
WRITE(*,*)
WRITE(*,*) '...configuring groundwater...'
END IF
!dt
IF (KeyIsPresent('dt', festIni, section = 'groundwater')) THEN
dtGroundwater = IniReadInt ('dt', festIni, section = 'groundwater')
ELSE
CALL Catch ('error', 'Fest', 'dt missing in groundwater' )
END IF
IF ( dtGroundwater > 0 ) THEN
timeNewGroundwater = timeStart
!configuration file
IF (KeyIsPresent('conf-file', festIni, section = 'groundwater')) THEN
string = IniReadString ('conf-file', festIni, section = 'groundwater')
CALL GroundwaterInit (string)
ELSE
CALL Catch ('error', 'Fest', 'conf-file missing in groundwater' )
END IF
!spatial data out dt
IF (KeyIsPresent('dt-out-aquifer', festIni, section = 'groundwater')) THEN
dtOutGroundwater = IniReadInt ('dt-out-aquifer', festIni, section = 'groundwater')
IF ( dtOutGroundwater > 0 ) THEN
timeOutGroundwater = timeStart
CALL GroundwaterOutputInit ( path_out)
ELSE
timeOutGroundwater = timeStart + (-1)
END IF
ELSE
CALL Catch ('error', 'Fest', 'dt-out-aquifer missing in groundwater' )
END IF
!point data out
IF (KeyIsPresent('out-point-file', festIni, section = 'groundwater')) THEN
string = IniReadString ('out-point-file', festIni, section = 'groundwater')
CALL GroundwaterPointInit ( string, path_out, timestart )
END IF
ELSE
dtOutGroundwater = 0
timeNewGroundwater = timeStart + (-1)
timeOutGroundwater = timeStart + (-1)
END IF
ELSE
dtGroundwater = 0
dtOutGroundwater = 0
timeNewGroundwater = timeStart + (-1)
timeOutGroundwater = timeStart + (-1)
END IF
!------------------------------------------------------------------------------
! soil erosion and sediment transport
!------------------------------------------------------------------------------
CALL Catch ('info', 'Fest', 'sediment erosion configuration' )
!check whether section is present
IF ( SectionIsPresent ('sediment', festIni) ) THEN !section is optional
IF ( verbose ) THEN
WRITE(*,*)
WRITE(*,*) '...configuring soil erosion and sediment transport...'
END IF
!dt
IF (KeyIsPresent('dt', festIni, section = 'sediment')) THEN
dtCalcSediment = IniReadInt ('dt', festIni, section = 'sediment')
ELSE
CALL Catch ('error', 'Fest', 'dt missing in sediment' )
END IF
IF ( dtCalcSediment > 0 ) THEN
!configuration file
IF (KeyIsPresent('conf-file', festIni, section = 'sediment')) THEN
string = IniReadString ('conf-file', festIni, section = 'sediment')
CALL SedimentInit (string, dtCalcSedimentRouting, string)
timeUpdateSediment = timeStart
timeUpdateSedimentRouting = timeStart
ELSE
CALL Catch ('error', 'Fest', 'conf-file missing in sediment' )
END IF
!spatial data out dt
IF (KeyIsPresent('dt-out-spatial', festIni, section = 'sediment')) THEN
dtOutSediment = IniReadInt ('dt-out-spatial', festIni, section = 'sediment')
ELSE
CALL Catch ('error', 'Fest', 'dt-out-spatial missing in sediment' )
END IF
IF (dtOutSediment > 0) THEN
timeOutSoilBalance = timeStart
END IF
timeOutSediment = timeStart
ELSE
dtOutSediment = 0
timeUpdateSediment = timeStart + (-1)
timeOutSediment = timeStart + (-1)
END IF
ELSE
dtCalcSediment = 0
dtOutSediment = 0
timeUpdateSediment = timeStart + (-1)
timeOutSediment = timeStart + (-1)
END IF
!TODO
! IF (routeSediment) THEN
! !set file for output of sediment routing
! fileunit_out_sediment_routing = GetUnit ()
! OPEN (fileunit_out_sediment_routing,file=string)
! WRITE(*,*)'nome file anagrafica out punti propagazione sedimenti: ',string(1:LEN_TRIM(string))
! CALL leggi_anagrafica(xsOutSedimentRouting,fileunit_out_sediment_routing)
! !check dt out sediment routing
! IF(.not.(mod(xsOutSedimentRouting%cadenza_misure,dtCalcSedimentRouting)==0)) THEN
! xsOutSedimentRouting%cadenza_misure=dtCalcSedimentRouting
! WRITE(*,*)'WARNING: xsOutSedimentRouting%cadenza_misure non multiplo di dtCalcSedimentRouting'
! WRITE(*,*)' si impone: xsOutSedimentRouting%cadenza_misure = dtCalcSedimentRouting'
! END IF
! WRITE(*,*)xsOutSedimentRouting%num_staz,' sezioni per risultati propagazione sedimenti'
! DO k=1,xsOutSedimentRouting%num_staz
! xsOutSedimentRouting%staz_oss(k)%misura=0.
! END DO
! CALL righe_colonne(xsOutSedimentRouting,grid_r=QoutSed)
! CLOSE (fileunit_out_sediment_routing)
!
! !write header in output file
! fileunit_out_sediment_routing = GetUnit ()
! OPEN (fileunit_out_sediment_routing,file=path_out(1:LEN_TRIM(path_out))//'sediment_routing.fest')
! WRITE(fileunit_out_sediment_routing,'(a)')'FEST-EWB: stima dei sedimentogrammi'
! WRITE(fileunit_out_sediment_routing,'(a)')'anagrafica'
! CALL scrivi_anagrafica(xsOutSedimentRouting,fileunit_out_sediment_routing)
! WRITE(fileunit_out_sediment_routing,'(a)')
! WRITE(fileunit_out_sediment_routing,'(a)')
! WRITE(fileunit_out_sediment_routing,'(a)')'dati'
! WRITE(fileunit_out_sediment_routing,'(a," ",\)')'#'
! DO k=1,xsOutSedimentRouting%num_staz-1
! WRITE(fileunit_out_sediment_routing,'(a," ",\)')(xsOutSedimentRouting%staz_oss(k)%denominazione)
! END DO
! WRITE(fileunit_out_sediment_routing,'(a)')(xsOutSedimentRouting%staz_oss(xsOutSedimentRouting%num_staz)%denominazione)
!
! !time
! timeOutSedimentRouting = timeStart
!END IF
!
!------------------------------------------------------------------------------
! spatial average
!------------------------------------------------------------------------------
CALL Catch ('info', 'Fest', 'spatial-average configuration' )
!check whether section is present
IF ( SectionIsPresent ('spatial-average', festIni) ) THEN !section is optional
IF ( verbose ) THEN
WRITE(*,*)
WRITE(*,*) '...configuring spatial average for output...'
END IF
!configuration file
IF (KeyIsPresent('conf-file', festIni, section = 'spatial-average')) THEN
string = IniReadString ('conf-file', festIni, section = 'spatial-average')
CALL InitSpatialAverageMeteo (string, path_out, temperatureAir, &
temperatureAirDailyMean, temperatureAirDailyMax, &
temperatureAirDailyMin, precipitationRate, &
relHumidityAir, radiation,netRadiation, &
windSpeed, precipitationRateDaily, irrigationFlux )
CALL InitSpatialAverageBalance (string, path_out, soilMoisture, &
soilMoistureRZ, soilMoistureRZ, &
runoff, infilt, percolation, &
et, pet, caprise, balanceError)
CALL InitSpatialAverageSnow (string, path_out, rainfallRate, swe, &
meltCoefficient, waterInSnow, snowMelt )
CALL InitSpatialAverageGlaciers (string, path_out, icewe, waterInIce, iceMelt )
CALL InitSpatialAverageSediment (string, path_out, interrillErosion)
CALL InitSpatialAverageCanopy (string, path_out, canopyStorage, raineff, canopyPT )
CALL InitSpatialAveragePlants (string, path_out, lai, gpp, npp, carbonstem, &
carbonroot, carbonleaf, fvcover, dbh, plantsHeight, &
density, stemyield )
ELSE
CALL Catch ('error', 'Fest', 'conf-file missing in spatial-average' )
END IF
!initialize meteo output file
IF ( dtOutMeteo > 0 ) THEN
CALL ComputeSpatialAverageMeteo (dtMeteo, temperatureAir, &
temperatureAirDailyMean, &
temperatureAirDailyMax, temperatureAirDailyMin, &
precipitationRate, relHumidityAir, radiation, &
netRadiation, windSpeed, precipitationRateDaily, &
irrigationFlux)
CALL ExportSpatialAverageMeteo ( init = .TRUE. )
END IF
!initialize snow output file
IF ( dtOutSnow > 0 ) THEN
CALL ComputeSpatialAverageSnow ( dtOutSnow, rainfallRate, &
swe, meltCoefficient, waterInSnow, &
snowMelt )
CALL ExportSpatialAverageSnow (init = .TRUE.)
END IF
!initialize glaciers output file
IF ( dtOutIce > 0 ) THEN
CALL ComputeSpatialAverageGlaciers (dtOutIce, icewe, &
waterInIce, iceMelt )
CALL ExportSpatialAverageGlaciers (init = .TRUE.)
END IF
!initailize vegetation canopy output file
IF (dtOutCanopy > 0) THEN
CALL ExportSpatialAverageCanopy (init = .TRUE.)
END IF
!initailize plants output file
IF (dtOutPlants > 0) THEN
CALL ExportSpatialAveragePlants (init = .TRUE.)
END IF
!initailize soil balance output file
IF (dtOutSoilBalance>0) THEN
CALL ExportSpatialAverageBalance (init = .TRUE.)
END IF
!initailize sediment output file
IF (dtCalcSediment > 0) THEN
CALL ExportSpatialAverageSediment (init = .TRUE.)
END IF
END IF
!------------------------------------------------------------------------------
! basin average
!------------------------------------------------------------------------------
CALL Catch ('info', 'Fest', 'basin-average configuration' )
!check whether section is present
IF ( SectionIsPresent ('basin-average', festIni) ) THEN !section is optional
IF ( verbose ) THEN
WRITE(*,*)
WRITE(*,*) '...configuring basin average for output...'
END IF
!output file
IF (KeyIsPresent('out-point-file', festIni, section = 'basin-average')) THEN
string = IniReadString ('out-point-file', festIni, section = 'basin-average')
CALL ReadPointFileBasinAverage ( string )
ELSE
CALL Catch ('error', 'Fest', 'out-point-file missing in basin-average' )
END IF
!configuration file
IF (KeyIsPresent('conf-file', festIni, section = 'basin-average')) THEN
string = IniReadString ('conf-file', festIni, section = 'basin-average')
CALL InitBasinAverage (string, path_out, temperatureAir, &
temperatureAirDailyMean, temperatureAirDailyMax, &
temperatureAirDailyMin, precipitationRate, &
relHumidityAir, radiation,netRadiation, &
windSpeed, precipitationRateDaily, irrigationFlux, &
swe, soilMoisture, runoff, infilt, percolation, &
et, pet, deltaSoilMoisture )
!set first export time
timeNewBasinAverage = timeStart
ELSE
CALL Catch ('error', 'Fest', 'conf-file missing in basin-average' )
END IF
ELSE
timeNewBasinAverage = timeStart + (-1)
END IF
!------------------------------------------------------------------------------
! raster export
!------------------------------------------------------------------------------
CALL Catch ('info', 'FeST', 'raster-export configuration' )
!check whether section is present
IF ( SectionIsPresent ('raster-export', festIni) ) THEN !section is optional
IF ( verbose ) THEN
WRITE(*,*)
WRITE(*,*) '...configuring raster export for output...'
END IF
!configuration file
IF (KeyIsPresent('conf-file', festIni, section = 'raster-export')) THEN
string = IniReadString ('conf-file', festIni, section = 'raster-export')
CALL InitRasterExport (string, temperatureAir, precipitationRate, &
relHumidityAir, radiation, netRadiation, &
windSpeed,swe, soilMoisture, runoff, infilt, &
percolation, et, pet )
!set first time
timeNewRasterExport = timeStart
ELSE
CALL Catch ('error', 'Fest', 'conf-file missing in raster-export' )
END IF
ELSE
timeNewRasterExport = timeStart + (-1)
END IF
!------------------------------------------------------------------------------
! discharge routing
!------------------------------------------------------------------------------
CALL Catch ('info', 'Fest', 'surface routing configuration' )
!check whether section is present
IF ( SectionIsPresent ('discharge-routing', festIni) ) THEN !section is optional
IF ( verbose ) THEN
WRITE(*,*)
WRITE(*,*) '...configuring surface routing...'
END IF
!dt
IF (KeyIsPresent('dt', festIni, section = 'discharge-routing')) THEN
dtDischargeRouting = IniReadInt ('dt', festIni, section = 'discharge-routing')
ELSE
CALL Catch ('error', 'Fest', 'dt missing in discharge-routing' )
END IF
IF ( dtDischargeRouting > 0 ) THEN
timeNewDischargeRouting = timeStart
!configuration file
IF (KeyIsPresent('conf-file', festIni, section = 'discharge-routing')) THEN
string = IniReadString ('conf-file', festIni, section = 'discharge-routing')
CALL InitDischargeRouting (fileini = string, time = timeStart, &
path_out = path_out, &
flowdirection = flowDirection, &
domain = mask , &
storage = storageChannelSurf, &
netRainfall = runoff, &
dtrk = dtReservoir )
ELSE
CALL Catch ('error', 'Fest', 'conf-file missing in discharge-routing' )
END IF
!output data file
IF (KeyIsPresent('out-point-file', festIni, section = 'discharge-routing')) THEN
string = IniReadString ('out-point-file', festIni, section = 'discharge-routing')
CALL DischargePointInit (string, path_out, timeStart)
END IF
ELSE
timeNewDischargeRouting = timeStart + (-1)
END IF
ELSE
dtDischargeRouting = 0
timeNewDischargeRouting = timeStart + (-1)
END IF
!------------------------------------------------------------------------------
! hot start
!------------------------------------------------------------------------------
CALL Catch ('info', 'Fest', 'save-hot-start configuration' )
!check whether section is present
IF ( SectionIsPresent ('save-hot-start', festIni) ) THEN !section is optional
hotstart = .TRUE.
IF ( verbose ) THEN
WRITE(*,*)
WRITE(*,*) '...configuring hot start...'
END IF
!time
IF (KeyIsPresent('time', festIni, section = 'save-hot-start')) THEN
cron = IniReadString ('time', festIni, section = 'save-hot-start')
CALL CronParseString (cron, cronHotstart )
!CALL ConfigureHotStart (cron)
ELSE
CALL Catch ('error', 'Fest', 'time missing in save-hot-start' )
END IF
!folder
IF (KeyIsPresent('folder', festIni, section = 'save-hot-start')) THEN
pathHotstart = IniReadString ('folder', festIni, section = 'save-hot-start')
ELSE
CALL Catch ('error', 'Fest', 'folder missing in save-hot-start' )
END IF
!save last time
IF (KeyIsPresent('save-last', festIni, section = 'save-hot-start')) THEN
IF ( IniReadInt ('save-last', festIni, section = 'save-hot-start') == 1 ) THEN
saveLast = .TRUE.
ELSE
saveLast = .FALSE.
END IF
ELSE !set to default
saveLast = .FALSE.
END IF
!IF ( dtHotstart > 0 ) THEN
!
! timeNewHotstart = timeStart
!
! !folder
! IF (KeyIsPresent('folder', festIni, section = 'save-hot-start')) THEN
! pathHotstart = IniReadString ('folder', festIni, section = 'save-hot-start')
! ELSE
! CALL Catch ('error', 'Fest', 'folder missing in save-hot-start' )
! END IF
!
!ELSE
! timeNewHotstart = timeStart + (-1)
!END IF
!ELSE
! dtHotstart = 0
! timeNewHotstart = timeStart + (-1)
ELSE
hotstart = .FALSE.
END IF
IF ( verbose ) THEN
WRITE(*,*)
WRITE(*,*) '....configuration completed!'
END IF
!------------------------------------------------------------------------------
! set calculation time step
!------------------------------------------------------------------------------
dtArray (1) = dtMeteo
dtArray (2) = dtSoilBalance
dtArray (3) = dtGroundwater
dtArray (4) = dtDischargeRouting
dtCalc = HUGE ( dtCalc )
DO i = 1, SIZE (dtArray)
IF ( dtArray (i) > 0 .AND. dtArray (i) < dtCalc ) THEN
dtCalc = dtArray (i)
END IF
END DO
!------------------------------------------------------------------------------
! start simulation
!------------------------------------------------------------------------------
IF ( verbose ) THEN
WRITE(*,*)
WRITE(*,*) 'starting simulation...'
WRITE(*,*)
WRITE(*,'(a,I10,a)') 'computation time step: ', dtCalc, ' second'
END IF
time = timeStart
!DO WHILE( .NOT. ( time < timeStop) )
DO WHILE( time <= timeStop )
timeString = time
IF ( verbose ) THEN
WRITE(*,'(a25)') timeString
END IF
!-----------------------------------------------------------------------------
! read meteorological forcings
!-----------------------------------------------------------------------------
IF (time == timeNewMeteo) THEN
IF ( verbose ) THEN
WRITE(*,*)'read meteorological forcings'
END IF
CALL MeteoRead (time,dem, albedo)
timeNewMeteo = timeNewMeteo + dtMeteo
IF ( dtOutMeteo > 0 ) THEN
timeSpatialAverageMeteo = time
CALL ComputeSpatialAverageMeteo (dtMeteo, temperatureAir, &
temperatureAirDailyMean, &
temperatureAirDailyMax, temperatureAirDailyMin, &
precipitationRate, relHumidityAir, radiation, &
netRadiation, windSpeed, precipitationRateDaily, &
irrigationFlux)
END IF
END IF
!-----------------------------------------------------------------------------
! update irrigation
!-----------------------------------------------------------------------------
IF (time == timeNewIrrigation) THEN
IF ( verbose ) THEN
WRITE(*,*) 'update irrigation'
END IF
CALL IrrigationUpdate (time, soilSat, Qout, dtOutIrrigation, &
timeOutIrrigation)
timeNewIrrigation = timeNewIrrigation + dtIrrigation
END IF
!-----------------------------------------------------------------------------
! snow
!-----------------------------------------------------------------------------
IF (dtSnow>0) THEN
IF (time == timeNewSnow) THEN
IF ( verbose ) THEN
WRITE(*,*)'snow'
END IF
CALL SnowUpdate (time, mask)
timeNewSnow = timeNewSnow + dtSnow
timeSpatialAverageSnow = time
CALL ComputeSpatialAverageSnow (dtOutSnow, rainfallRate, swe, &
meltCoefficient, waterInSnow, &
snowMelt)
END IF
ELSE
rainfallRate = precipitationRate
END IF
!-----------------------------------------------------------------------------
! glacier
!-----------------------------------------------------------------------------
IF ( time == timeNewIce ) THEN
IF ( verbose ) THEN
WRITE(*,*) 'glaciers'
END IF
CALL GlacierUpdate (time, mask, rainfallRate )
timeNewIce = timeNewIce + dtIce
timeSpatialAverageIce = time
CALL ComputeSpatialAverageGlaciers (dtOutIce, icewe, &
waterInIce, iceMelt)
END IF
!-----------------------------------------------------------------------------
!CANOPY INTERCEPTION
!-----------------------------------------------------------------------------
IF ( dtCanopyInterception > 0 ) THEN
IF ( time == timeNewCanopyInterception) THEN
IF ( verbose ) THEN
WRITE(*,*) 'canopy interception'
END IF
CALL Throughfall (rainfallRate , lai, mask, fvcover, raineff)
!CALL AdjustPT(fvcover, mask, tp_oss, dtetp)
timeNewCanopyInterception = timeNewCanopyInterception + dtCanopyInterception
timeSpatialAverageCanopy = time
CALL ComputeSpatialAverageCanopy (dtCanopyInterception, canopyStorage, raineff, canopyPT)
END IF
ELSE !copy rainfallRate into raineff
raineff = rainfallRate
END IF
!-------------------------------------------------------------------------------------------------------------
!FOREST DYNAMIC
!-------------------------------------------------------------------------------------------------------------
!check if plants parameters need update
IF ( updatePlantsParameters ) THEN
CALL PlantsParameterUpdate ( time )
END IF
IF ( dtPlants > 0 ) THEN
IF ( time == timeNewPlants) THEN
IF ( verbose ) THEN
WRITE(*,*) 'plants'
END IF
! CALL PlantsGrow (t, radiation, temperatureAir, soilMoisture, fieldCapacity, wiltingPoint, relHumidityAir, co2 % npd (1) % value)
timeNewPlants = timeNewPlants + dtPlants
timeSpatialAveragePlants = time
CALL ComputeSpatialAveragePlants (dtPlants, lai, gpp, npp, carbonstem, carbonroot, carbonleaf, fvcover, dbh, plantsHeight, density, stemyield)
END IF
END IF
!-----------------------------------------------------------------------------
! soil water balance
!-----------------------------------------------------------------------------
IF ( dtSoilBalance > 0 ) THEN
IF ( time == timeNewSoilBalance ) THEN
IF ( verbose ) THEN
WRITE(*,*)'soil water balance'
END IF
CALL SolveSoilBalance ( time, raineff, irrigationFlux, &
flowDirection, fvcover, vadoseDepth)
timeSpatialAverageBalance = time
CALL ComputeSpatialAverageBalance (dtSoilBalance, soilMoisture, &
soilMoistureRZ, soilMoistureTZ, &
runoff, infilt, percolation, &
et, pet, caprise, balanceError)
timeNewSoilBalance = timeNewSoilBalance + dtSoilBalance
END IF
ELSE
IF ( ALLOCATED (runoff % mat) ) THEN
DO i = 1, raineff % idim
DO j = 1, raineff % jdim
runoff % mat(i,j) = raineff % mat(i,j)
!percolation%mat(i,j)=0.0
END DO
END DO
END IF
END IF
!--------------------------------------------------------------------------
! groundwater
!--------------------------------------------------------------------------
IF (time == timeNewGroundwater) THEN
IF ( verbose ) THEN
WRITE(*,*) 'groundwater'
END IF
!update river groundwater interaction discharge
IF (riverGroundwaterInteract ) THEN
CALL GroundwaterRiverUpdate ( waterDepth, topWidth )
END IF
!update groundwater head
CALL GroundwaterUpdate ( time, QinSoilSub, percolation, caprise )
IF ( time == timeOutGroundwater) THEN
CALL GroundwaterOutput (time)
timeOutGroundwater = timeOutGroundwater + dtOutGroundwater
END IF
timeNewGroundwater = timeNewGroundwater + dtGroundwater
ENDIF
!-----------------------------------------------------------------------------
! discharge routing
!-----------------------------------------------------------------------------
IF ( time == timeNewDischargeRouting ) THEN
IF ( verbose ) THEN
WRITE(*,*)'surface routing'
END IF
CALL DischargeRoute (dt = dtDischargeRouting, time = time, &
flowdir = flowDirection, &
runoff = runoff, dtrk = dtReservoir, &
riverToGroundwater = riverToGroundwater, &
groundwaterToriver = groundwaterToRiver, &
storage = storageChannelSurf)
timeNewDischargeRouting = timeNewDischargeRouting + dtDischargeRouting
END IF
!-----------------------------------------------------------------------------
! sediment detachment
!-----------------------------------------------------------------------------
IF (time == timeUpdateSediment) THEN
CALL Catch ('info', 'FEST-EWB', 'update sediment erosion and deposition')
CALL InterrillDetachmentRate (raineff)
timeUpdateSediment = timeUpdateSediment + dtCalcSediment
timeSpatialAverageSediment = time
CALL ComputeSpatialAverageSediment (dtCalcSediment, interrillErosion)
END IF
!sediment routing
IF (routeSediment) THEN
IF (time == timeUpdateSedimentRouting) THEN
CALL SedimentRouting (dtCalcSedimentRouting, Qin, Velocita, waterDepth)
!update timeUpdateSedimentRouting
timeUpdateSedimentRouting = timeUpdateSedimentRouting + dtCalcSedimentRouting
END IF
END IF
!-----------------------------------------------------------------------------
! write spatial average results
!-----------------------------------------------------------------------------
!output meteo spatial average
IF (time == timeOutMeteo) THEN
IF ( verbose ) THEN
WRITE(*,*) 'export meteo'
END IF
timeOutMeteo = timeOutMeteo + dtOutMeteo
CALL ExportSpatialAverageMeteo ()
END IF
!output canopy interception spatial average
IF ( time == timeOutCanopy ) THEN
IF ( verbose ) THEN
WRITE(*,*) 'export canopy interception'
END IF
timeOutCanopy = timeOutCanopy + dtOutCanopy
CALL ExportSpatialAverageCanopy ()
END IF
!output plants spatial average
IF ( time == timeOutPlants ) THEN
IF ( verbose ) THEN
WRITE(*,*) 'export plants'
END IF
timeOutPlants = timeOutPlants + dtOutPlants
CALL ExportSpatialAveragePlants ()
END IF
!output snow spatial average
IF (time==timeOutSnow) THEN
IF ( verbose ) THEN
WRITE(*,*) 'export snow'
END IF
timeOutSnow = timeOutSnow + dtOutSnow
CALL ExportSpatialAverageSnow ()
END IF
!output glaciers spatial average
IF (time == timeOutIce) THEN
IF ( verbose ) THEN
WRITE(*,*) 'export glaciers'
END IF
timeOutIce = timeOutIce + dtOutIce
CALL ExportSpatialAverageGlaciers ()
END IF
!output soil balance spatial average
IF ( time == timeOutSoilBalance ) THEN
IF ( verbose ) THEN
WRITE(*,*) 'export soil balance'
END IF
CALL ExportSpatialAverageBalance ()
timeOutSoilBalance = timeOutSoilBalance + dtOutSoilBalance
END IF
!!output ragguagli falda
!IF ( time == timeOutGroundwater ) THEN
! IF ( verbose ) THEN
! WRITE(*,*) 'export groundwater'
! END IF
! timeOutGroundwater = timeOutGroundwater+dtOutGroundwater
! !CALL esporta_ragguagli_falda(path_out)
!END IF
!output sediment spatial average
IF (time == timeOutSediment) THEN
IF ( verbose ) THEN
WRITE(*,*) 'export sediment'
END IF
timeOutSediment = timeOutSediment + dtOutSediment
CALL ExportSpatialAverageSediment ()
END IF
!output sediment routing
!IF (time == timeOutSedimentRouting) THEN
! xsOutSedimentRouting % istante = time
! CALL celle2stazioni(QoutSed,xsOutSedimentRouting)
! CALL scrivi_dato(xsOutSedimentRouting,fileunit_out_sediment_routing)
! !t
! timeOutSedimentRouting = timeOutSedimentRouting + xsOutSedimentRouting % cadenza_misure
!
!END IF
!-----------------------------------------------------------------------------
! basin average
!-----------------------------------------------------------------------------
IF ( time == timeNewBasinAverage ) THEN
IF ( verbose ) THEN
WRITE (*,*) 'export basin average'
END IF
CALL ExportBasinAverage (time, temperatureAir, &
temperatureAirDailyMean, temperatureAirDailyMax, &
temperatureAirDailyMin, precipitationRate, &
relHumidityAir, radiation,netRadiation, &
windSpeed, precipitationRateDaily, irrigationFlux, &
swe, soilMoisture, runoff, infilt, percolation, &
et, pet, deltaSoilMoisture )
timeNewBasinAverage = timeNewBasinAverage + dtBasinAverage
END IF
!-----------------------------------------------------------------------------
! raster export
!-----------------------------------------------------------------------------
IF ( time == timeNewRasterExport ) THEN
IF ( verbose ) THEN
WRITE (*,*) 'update raster maps'
END IF
CALL ExportMaps (time, dtCalc, temperatureAir, precipitationRate, &
relHumidityAir, radiation, netRadiation, &
windSpeed, swe, soilMoisture, runoff, &
infilt, percolation, et, pet )
timeNewRasterExport = timeNewRasterExport + dtCalc
END IF
!-----------------------------------------------------------------------------
! save status
!-----------------------------------------------------------------------------
IF ( hotstart ) THEN
!IF ( IsTimeHotStart (time) ) THEN
IF ( CronIsTime (time, cronHotstart) ) THEN
IF ( verbose ) THEN
WRITE (*,*) 'save state variables for hot start'
END IF
CALL SaveHotStart (time)
END IF
END IF
!-----------------------------------------------------------------------------
! time forward
!-----------------------------------------------------------------------------
time = time + dtCalc
IF ( time > timeStop ) THEN
IF ( hotstart ) THEN
IF ( saveLast > 0 ) THEN
IF ( verbose ) THEN
WRITE (*,*) 'save final state variables for hot start'
END IF
CALL SaveHotStart
END IF
END IF
END IF
END DO
!final message
IF ( verbose ) THEN
WRITE(*,*) 'simulation finished'
END IF
!-----------------------------------------
! terminate log
!-----------------------------------------
CALL LogStop ()
CONTAINS
!==============================================================================
!| Description:
! Save state variables for starting a new simulation from the current condition
!
SUBROUTINE SaveHotStart &
!
(time)
IMPLICIT NONE
TYPE(DateTime), OPTIONAL :: time
CHARACTER (LEN = 26) :: prefix
CHARACTER (LEN = 300) :: file
!-------------------------end of declarations----------------------------------
IF (PRESENT(time)) THEN
prefix = time
prefix = prefix (1:19) // '_'
prefix (14:14) = '-'
prefix (17:17) = '-'
ELSE
prefix = ' '
END IF
!snow
IF ( dtSnow > 0 ) THEN
!swe
file = pathHotstart (1:LEN_TRIM (pathHotstart)) // &
prefix (1:LEN_TRIM(prefix)) // 'swe'
CALL ExportGrid (layer = swe, fileName = file, &
fileFormat = ESRI_BINARY)
!water in snow
file = pathHotstart (1:LEN_TRIM (pathHotstart)) // &
prefix (1:LEN_TRIM(prefix)) // 'water-snow'
CALL ExportGrid (layer = waterInSnow, fileName = file, &
fileFormat = ESRI_BINARY)
END IF
!soil balance
IF ( dtSoilBalance > 0 ) THEN
!soil saturation root zone
file = pathHotstart (1:LEN_TRIM (pathHotstart)) // &
prefix (1:LEN_TRIM(prefix)) // 'sat-rz'
CALL ExportGrid (layer = soilSatRZ, fileName = file, &
fileFormat = ESRI_BINARY)
!soil saturation transmission zone
file = pathHotstart (1:LEN_TRIM (pathHotstart)) // &
prefix (1:LEN_TRIM(prefix)) // 'sat-tz'
CALL ExportGrid (layer = soilSatTZ, fileName = file, &
fileFormat = ESRI_BINARY)
!precipitation status
file = pathHotstart (1:LEN_TRIM (pathHotstart)) // &
prefix (1:LEN_TRIM(prefix)) // 'rainflag'
CALL ExportGrid (layer = rainFlag, fileName = file, &
fileFormat = ESRI_BINARY)
!interstorm duration
file = pathHotstart (1:LEN_TRIM (pathHotstart)) // &
prefix (1:LEN_TRIM(prefix)) // 'interstorm'
CALL ExportGrid (layer = interstormDuration, fileName = file, &
fileFormat = ESRI_BINARY)
IF ( infiltrationModel == SCS_CN ) THEN
!soil retention
file = pathHotstart (1:LEN_TRIM (pathHotstart)) // &
prefix (1:LEN_TRIM(prefix)) // 'seff'
CALL ExportGrid (layer = sEff, fileName = file, &
fileFormat = ESRI_BINARY)
!cumulative precipitation
file = pathHotstart (1:LEN_TRIM (pathHotstart)) // &
prefix (1:LEN_TRIM(prefix)) // 'raincum'
CALL ExportGrid (layer = raincum, fileName = file, &
fileFormat = ESRI_BINARY)
END IF
IF ( infiltrationModel == PHILIPEQ .OR. &
infiltrationModel == GREEN_AMPT ) THEN
!cumulative infiltration
file = pathHotstart (1:LEN_TRIM (pathHotstart)) // &
prefix (1:LEN_TRIM(prefix)) // 'cuminf'
CALL ExportGrid (layer = cuminf, fileName = file, &
fileFormat = ESRI_BINARY)
END IF
END IF
! discharge routing
IF ( dtDischargeRouting > 0 ) THEN
!Qin
file = pathHotstart (1:LEN_TRIM (pathHotstart)) // &
prefix (1:LEN_TRIM(prefix)) // 'Qin'
CALL ExportGrid (layer = Qin, fileName = file, &
fileFormat = ESRI_BINARY)
!Qout
file = pathHotstart (1:LEN_TRIM (pathHotstart)) // &
prefix (1:LEN_TRIM(prefix)) // 'Qout'
CALL ExportGrid (layer = Qout, fileName = file, &
fileFormat = ESRI_BINARY)
!Qlat
file = pathHotstart (1:LEN_TRIM (pathHotstart)) // &
prefix (1:LEN_TRIM(prefix)) // 'Qlat'
CALL ExportGrid (layer = Plat, fileName = file, &
fileFormat = ESRI_BINARY)
!TODO lakes
END IF
! !propagazione sup
! IF (dtcalcolo_propagazione > 0) THEN
! CALL write_raster_bin(path_hotstart(1:LEN_TRIM(path_hotstart))// &
! pref_s(1:LEN_TRIM(pref_s))//'Q_in',matrice_real=Q_in)
! CALL write_raster_bin(path_hotstart(1:LEN_TRIM(path_hotstart))// &
! pref_s(1:LEN_TRIM(pref_s))//'Q_out',matrice_real=Q_out)
! END IF
!! !invasi
! IF ( dtcalcolo_invasi > 0 ) THEN
! CALL esporta_livello_invasi(path_hotstart(1:LEN_TRIM(path_hotstart))// &
! pref_s(1:LEN_TRIM(pref_s))//'livello_invasi.txt',laghi,time)
! ENDIF
! !Falda
!IF (dtGroundwater > 0) THEN
! DO I=1,AltezzaPiezometrica%NumeroStrati
!
! CALL write_raster_bin(path_hotstart(1:LEN_TRIM(path_hotstart))// &
! pref_s(1:LEN_TRIM(pref_s))//'strato'//CHAR(64+I), &
! matrice_real=AltezzaPiezometrica%piezometria(I))
!
! IF (interazione_falda_fiume ) THEN
!
! CALL write_raster_bin(path_hotstart(1:LEN_TRIM(path_hotstart))// &
! pref_s(1:LEN_TRIM(pref_s))//'baseflow_falda', &
! matrice_real=BaseFlow_falda)
!
! CALL write_raster_bin(path_hotstart(1:LEN_TRIM(path_hotstart))// &
! pref_s(1:LEN_TRIM(pref_s))//'ricarica_falda', &
! matrice_real=ricarica_falda_fiume)
!
! END IF
! ENDDO
!END IF
!sediment
!IF (dtCalcSediment > 0) THEN
! CALL write_raster_bin(path_hotstart(1:LEN_TRIM(path_hotstart))// &
! pref_s(1:LEN_TRIM(pref_s))//'deltaSed',matrice_real=deltaSed)
!END IF
!
!IF (dtCalcSedimentRouting > 0) THEN
!
!!to do
!
!END IF
END SUBROUTINE SaveHotStart
SUBROUTINE Splash ()
WRITE (*,*)
!WRITE (*,'(a)') ' ______ ______ _____ _______ '
!WRITE (*,'(a)') ' | ____| ____|/ ____|__ __|'
!WRITE (*,'(a)') ' | |__ | |__ | (___ | | '
!WRITE (*,'(a)') ' | __| | __| \___ \ | | '
!WRITE (*,'(a)') ' | | | |____ ____) | | | '
!WRITE (*,'(a)') ' |_| |______|_____/ |_| '
!WRITE (*,'(a)')
!WRITE (*,'(a)') '-o--o o--o o-o o--O-o '
!WRITE (*,'(a)') ' | | | | '
!WRITE (*,'(a)') ' O-o O-o o-o | '
!WRITE (*,'(a)') ' | | | | '
!WRITE (*,'(a)') ' o o--o o--o o '
!WRITE (*,'(a)')
WRITE (*,'(a)') ' ________ _______ ________ _________ '
WRITE (*,'(a)') '|\ _____\|\ ___ \ |\ ____\ |\___ ___\ '
WRITE (*,'(a)') '\ \ \__/ \ \ __/| \ \ \___|_ \|___ \ \_| '
WRITE (*,'(a)') ' \ \ __\ \ \ \_|/__ \ \_____ \ \ \ \ '
WRITE (*,'(a)') ' \ \ \_| \ \ \_|\ \ \|____|\ \ \ \ \ '
WRITE (*,'(a)') ' \ \__\ \ \_______\ ____\_\ \ \ \__\ '
WRITE (*,'(a)') ' \|__| \|_______| |\_________\ \|__| '
WRITE (*,'(a)') ' \|_________| '
WRITE (*,'(a)') ' '
! (*,'(a)') ' E N E R G Y - W A T E R - B A L A N C E '
WRITE (*,'(a)') ' S P A T I A L L Y - D I S T R I B U T E D '
WRITE (*,'(a)') ' H Y D R O L O G I C A L - M O D E L '
!WRITE (*,'(a)') ' Flood Event Spatially distributed rainfall-runoff Transformation model'
WRITE (*,'(a)')
RETURN
END SUBROUTINE Splash
!-------------------------------------------------------------------------------------------------------------
!-------------------------------------------------------------------------------------------------------------
END PROGRAM Fest