diff --git a/.gitignore b/.gitignore
new file mode 100644
index 0000000..3e74b12
--- /dev/null
+++ b/.gitignore
@@ -0,0 +1,141 @@
+# Visual Studio Code project settings
+.vscode
+
+# Byte-compiled / optimized / DLL files
+__pycache__/
+*.py[cod]
+*$py.class
+
+# C extensions
+*.so
+
+# Distribution / packaging
+.Python
+build/
+develop-eggs/
+dist/
+downloads/
+eggs/
+.eggs/
+lib/
+lib64/
+parts/
+sdist/
+var/
+wheels/
+share/python-wheels/
+*.egg-info/
+.installed.cfg
+*.egg
+MANIFEST
+
+# PyInstaller
+# Usually these files are written by a python script from a template
+# before PyInstaller builds the exe, so as to inject date/other infos into it.
+*.manifest
+*.spec
+
+# Installer logs
+pip-log.txt
+pip-delete-this-directory.txt
+
+# Unit test / coverage reports
+htmlcov/
+.tox/
+.nox/
+.coverage
+.coverage.*
+.cache
+nosetests.xml
+coverage.xml
+*.cover
+*.py,cover
+.hypothesis/
+.pytest_cache/
+cover/
+
+# Translations
+*.mo
+*.pot
+
+# Django stuff:
+*.log
+local_settings.py
+db.sqlite3
+db.sqlite3-journal
+
+# Flask stuff:
+instance/
+.webassets-cache
+
+# Scrapy stuff:
+.scrapy
+
+# Sphinx documentation
+docs/_build/
+
+# PyBuilder
+.pybuilder/
+target/
+
+# Jupyter Notebook
+.ipynb_checkpoints
+
+# IPython
+profile_default/
+ipython_config.py
+
+# pyenv
+# For a library or package, you might want to ignore these files since the code is
+# intended to run in multiple environments; otherwise, check them in:
+# .python-version
+
+# pipenv
+# According to pypa/pipenv#598, it is recommended to include Pipfile.lock in version control.
+# However, in case of collaboration, if having platform-specific dependencies or dependencies
+# having no cross-platform support, pipenv may install dependencies that don't work, or not
+# install all needed dependencies.
+#Pipfile.lock
+
+# PEP 582; used by e.g. github.com/David-OConnor/pyflow
+__pypackages__/
+
+# Celery stuff
+celerybeat-schedule
+celerybeat.pid
+
+# SageMath parsed files
+*.sage.py
+
+# Environments
+.env
+.venv
+env/
+venv/
+ENV/
+env.bak/
+venv.bak/
+
+# Spyder project settings
+.spyderproject
+.spyproject
+
+# Rope project settings
+.ropeproject
+
+# mkdocs documentation
+/site
+
+# mypy
+.mypy_cache/
+.dmypy.json
+dmypy.json
+
+# Pyre type checker
+.pyre/
+
+# pytype static type analyzer
+.pytype/
+
+# Cython debug symbols
+cython_debug/
diff --git a/README.md b/README.md
index 7bd2b66..0d741af 100644
--- a/README.md
+++ b/README.md
@@ -1,2 +1,147 @@
# RUBEM
-Projeto teste para upload no repositorio do Qgis
+
+
+
+
+
+ 
+
+
+
RUBEM
+
+
+ Rainfall rUnoff Balance Enhanced Model is a distributed hydrological model to calculate monthly flows with changes in land use over time
+
+ Explore the docs »
+
+
+ View Demo
+ ·
+ Report Bug
+ ·
+ Request Feature
+
+
+
+
+
+
+
+ Table of Contents
+
+ -
+ About The Project
+
+
+ -
+ Getting Started
+
+
+ - Usage
+ - Roadmap
+ - Contributing
+ - License
+ - Contact
+ - Acknowledgements
+
+
+
+
+
+
+## About The Project
+
+[![Product Name Screen Shot][product-screenshot]](https://example.com)
+
+Rainfall rUnoff Balance Enhanced Model (RUBEM) is an improved model of balance between rain and runoff. The distributed hydrological model for transforming precipitation into total flow is based on equations that represent the physical processes of the hydrological cycle, with spatial distribution defined in a grid and monthly time scale. The model was developed based on classic concepts of hydrological processes and equations based mainly on the formulations of the SPHY (TERINK et al., 2015), WEAP (YATES et al., 2005) and WetSpass-M (ABDOLLAHI et al., 2017).
+
+The name is a posthumous tribute to Professor Rubem La Laina Porto, dean of the Department of Hydraulic and Environmental Engineering, of the Polytechnic School of USP, who dedicated his professional life to the study, development and practices in hydrological sciences, contributing to the improvement of the state of art and training of students and professionals working in the area.
+
+### Built With
+
+* []()
+* []()
+* []()
+
+
+
+
+## Getting Started
+
+To get a local copy up and running follow these simple steps.
+
+### Prerequisites
+
+This is an example of how to list things you need to use the software and how to install them.
+* npm
+ ```sh
+ ...
+ ```
+
+### Installation
+
+1. Clone the repo
+ ```sh
+ git clone https://github.com/LabSid-USP/RUBEM.git
+ ```
+2. Install
+ ```sh
+ ...
+ ```
+
+
+
+
+## Usage
+
+Use this space to show useful examples of how a project can be used. Additional screenshots, code examples and demos work well in this space. You may also link to more resources.
+
+_For more examples, please refer to the [Documentation](https://example.com)_
+
+
+
+
+## Roadmap
+
+See the [open issues](https://github.com/LabSid-USP/RUBEM/issues) for a list of proposed features (and known issues).
+
+
+
+
+## Contributing
+
+Contributions are what make the open source community such an amazing place to be learn, inspire, and create. Any contributions you make are **greatly appreciated**.
+
+1. Fork the Project
+2. Create your Feature Branch (`git checkout -b feature/AmazingFeature`)
+3. Commit your Changes (`git commit -m 'Add some AmazingFeature'`)
+4. Push to the Branch (`git push origin feature/AmazingFeature`)
+5. Open a Pull Request
+
+
+
+
+## License
+
+Distributed under the GPLv3 License. See `LICENSE` for more information.
+
+
+
+
+## Contact
+
+Project Link: [https://github.com/LabSid-USP/RUBEM](https://github.com/LabSid-USP/RUBEM)
+
+
+
+
+## Acknowledgements
+
+* [Laboratório de Sistemas de Suporte a Decisões Aplicados à Engenharia Ambiental e de Recursos Hídricos](http://labsid.eng.br/Contato.aspx)
+* [Departamento de Engenharia Hidráulica e Ambiental da Escola Politécnica da Universidade de São Paulo](http://www.pha.poli.usp.br/)
+* [Fundo Patrimonial Amigos da Poli](https://www.amigosdapoli.com.br/)
diff --git a/RainfallRunoff.py b/RainfallRunoff.py
index 5a088f6..11291ad 100644
--- a/RainfallRunoff.py
+++ b/RainfallRunoff.py
@@ -1,30 +1,35 @@
-import time
-t1 = time.time()
-print("Inicio", flush=True)
+# coding=utf-8
-import configparser
import os
-import gdal
+import time
+import datetime
+import calendar
+import configparser
+# importacao GDAL mais recente mantendo compatibilidade
+try:
+ from osgeo import gdal
+except ImportError:
+ import gdal
import numpy as np
-from pcraster.framework import dynamicBase
from pcraster.framework import *
import pcraster as pcr
+
# inportacao de funcoes do modelo chuva vazao
from interception import *
from evapotranspiration import *
from surface_runoff import *
from soil import *
+t1 = time.time()
+print("Inicio", flush=True)
+
# Nome do arquivo de configuraçoes
configFile = 'config.ini'
# Leitura de arquivo config.ini
config = configparser.ConfigParser()
config.read(configFile)
-# Leitura de local de arquivos de entrada
-files = config.get('FILES', 'dir')
-# Nome do arquivo com série de vazao observada nos postos dados - deve ser criado no formato compatível
-obsData = config.get('FILES', 'obsdata')
+
# Data inicial e final da simulacao
startDate = config.get('SIM_TIME', 'start')
endDate = config.get('SIM_TIME', 'end')
@@ -43,7 +48,13 @@
########## Funcoes auxiliares ##########
gdal.UseExceptions()
def getRefInfo(self, sourceTif):
- import gdal
+ """
+ :param sourceTif:
+ :sourceTif type:
+
+ :returns:
+ :rtype:
+ """
gdal.AllRegister()
ds = gdal.OpenEx(sourceTif)
cols = ds.RasterXSize
@@ -55,7 +66,25 @@ def getRefInfo(self, sourceTif):
return(Ref)
def reportTif(self, tifRef, pcrObj, fileName, outpath, din = 0):
- import gdal
+ """
+ :param tifRef:
+ :tifRef type:
+
+ :param pcrObj:
+ :pcrObj type:
+
+ :param fileName:
+ :fileName type:
+
+ :param outpath:
+ :outpath type:
+
+ :param din: Defaults to 0.
+ :din type:
+
+ :returns:
+ :rtype:
+ """
# sourceTif = file to get attibutes from - DEM
# pcrObj = pcraster to export
# fileName = string format
@@ -91,7 +120,16 @@ def reportTif(self, tifRef, pcrObj, fileName, outpath, din = 0):
# Calculo de numero de meses (steps) com base nas datas inicial e final de simulacao
def totalSteps(startDate, endDate):
- import time, datetime
+ """
+ :param startDate:
+ :startDate type:
+
+ :param endDate:
+ :endDate type:
+
+ :returns:
+ :rtype:
+ """
# begin simulation
yBegin = int(datetime.datetime.strptime(startDate ,'%d/%m/%Y').strftime("%Y"))
monthBegin = int(datetime.datetime.strptime(startDate ,'%d/%m/%Y').strftime("%m"))
@@ -107,7 +145,16 @@ def totalSteps(startDate, endDate):
# Calculo de numero de dias no mes a partir do timestep (para conversao de vazao de mm para m3/s)
def daysOfMonth(startDate, timestep):
- import calendar, datetime
+ """
+ :param startDate:
+ :startDate type:
+
+ :param timestep:
+ :timestep type:
+
+ :returns:
+ :rtype:
+ """
sourcedate = datetime.datetime.strptime(startDate,'%d/%m/%Y')
month = sourcedate.month -2 + timestep
year = sourcedate.year + month // 12
@@ -118,298 +165,299 @@ def daysOfMonth(startDate, timestep):
########## Dynamic Mode ##########
#raise SystemExit
class Modelo(DynamicModel):
- def __init__(self):
- DynamicModel.__init__(self)
- print("Lendo arquivos de entrada", flush=True)
- # Read file locations
- config = configparser.ConfigParser()
- config.read(configFile)
- self.inpath = config.get('FILES', 'input')
- self.dem_file = config.get('FILES', 'dem')
- self.demTif = config.get('FILES', 'demtif')
- self.clone_file = config.get('FILES', 'clone')
- self.ldd_file = config.get('FILES', 'lddTif')
- self.etp_path = config.get('FILES', 'etp')
- self.prec_path = config.get('FILES', 'prec')
- self.ndvi_path = config.get('FILES', 'ndvi')
- self.kp_path = config.get('FILES', 'kp')
- self.land_path = config.get('FILES', 'landuse')
- self.soil_path = config.get('FILES', 'solo')
- self.outpath = config.get('FILES', 'output')
- self.sampleLocs = config.get('FILES', 'samples')
-
- # Set clone
- pcr.setclone(self.clone_file)
-
- # Read temporal filenames prefix
- self.etpPrefix = config.get('FILES', 'etpFilePrefix')
- self.precPrefix = config.get('FILES', 'precFilePrefix')
- self.ndviPrefix = config.get('FILES', 'ndviFilePrefix')
- self.ndviMaxFile = config.get('FILES', 'ndvimaxPrefix')
- self.ndviMinFile = config.get('FILES', 'ndviminPrefix')
- self.kpPrefix = config.get('FILES', 'kpFilePrefix')
- self.coverPrefix = config.get('FILES', 'landuseFilePrefix')
-
- # Read text lookuptables from config file
- self.rainyDaysTable = config.get('PARAMETERS', 'rainydays')
- self.aiTable = config.get('PARAMETERS', 'a_i')
- self.aoTable = config.get('PARAMETERS', 'a_o')
- self.asTable = config.get('PARAMETERS', 'a_s')
- self.avTable = config.get('PARAMETERS', 'a_v')
- self.manningTable = config.get('PARAMETERS', 'manning')
- self.dgTable = config.get('PARAMETERS', 'dg')
- self.KrTable = config.get('PARAMETERS', 'kr')
- self.TccTable = config.get('PARAMETERS', 'capCampo')
- self.Tporosidade = config.get('PARAMETERS', 'porosidade')
- self.TsatTable = config.get('PARAMETERS', 'saturacao')
- self.TwTable = config.get('PARAMETERS', 'pontomurcha')
- self.ZrTable = config.get('PARAMETERS', 'zr')
- self.KcminTable = config.get('PARAMETERS', 'kcmin')
- self.KcmaxTable = config.get('PARAMETERS', 'kcmax')
-
- # Area da celula #pendente = calculo automatico da area da celula
- self.A = config.getfloat('GRID', 'grid')
-
- # Read calibration parameters from config file
- self.alfa = config.getfloat('CALIBRATION', 'alfa')
- self.b = config.getfloat('CALIBRATION', 'b')
- self.w1 = config.getfloat('CALIBRATION', 'w1')
- self.w2 = config.getfloat('CALIBRATION', 'w2')
- self.w3 = config.getfloat('CALIBRATION', 'w3')
- self.RCD = config.getfloat('CALIBRATION', 'rcd')
- self.f = config.getfloat('CALIBRATION', 'f')
- self.alfa_gw = config.getfloat('CALIBRATION', 'alfa_gw')
- self.x = config.getfloat('CALIBRATION', 'x')
-
- # Read soil conditions from config file
- # teor de umidade inicial da zona radicular (fracao do teor de saturacao)
- self.ftur_ini = config.getfloat('INITIAL SOIL CONDITIONS', 'ftur_ini')
- # escoamento basico inicial
- self.EBini = scalar(config.getfloat('INITIAL SOIL CONDITIONS', 'eb_ini'))
- # limite para escoamento basico
- self.EBlim = scalar(config.getfloat('INITIAL SOIL CONDITIONS', 'eb_lim'))
- # teor de umidade inicial da zona saturada
- self.Tusini = scalar(config.getfloat('INITIAL SOIL CONDITIONS', 'tus_ini'))
-
- # constantes
- self.fpar_max = config.getfloat('CONSTANT', 'fpar_max')
- self.fpar_min = config.getfloat('CONSTANT', 'fpar_min')
- self.lai_max = config.getfloat('CONSTANT', 'lai_max')
- self.I_i = config.getfloat('CONSTANT', 'i_imp')
-
- # # Initialize time series output
- self.OutTssRun= ('outRun')
-
- # Report file
- import time
- # name
- self.timeStamp = str((time.strftime("%Y%m%d_%H%M%S",time.localtime(t1))))
- # header
- self.t_round = str(time.strftime("%Y %m %d %H:%M:%S",time.localtime(t1)))
-
- # Get Tif file Reference
- self.ref = getRefInfo(self, self.demTif)
-
- def initial(self):
-
- # Read dem file
- self.dem = readmap(self.dem_file)
-
- # Generate the local drain direction map on basis of the elevation map
- # self.ldd = lddcreate(self.dem, 1e31, 1e31, 1e31, 1e31)
- self.ldd = pcr.ldd(readmap(self.ldd_file))
+ """Constructor"""
+ def __init__(self):
+ DynamicModel.__init__(self)
+ print("Lendo arquivos de entrada", flush=True)
+ # Read file locations
+ config = configparser.ConfigParser()
+ config.read(configFile)
+ self.inpath = config.get('FILES', 'input')
+ self.dem_file = config.get('FILES', 'dem')
+ self.demTif = config.get('FILES', 'demtif')
+ self.clone_file = config.get('FILES', 'clone')
+ self.ldd_file = config.get('FILES', 'lddTif')
+ self.etp_path = config.get('FILES', 'etp')
+ self.prec_path = config.get('FILES', 'prec')
+ self.ndvi_path = config.get('FILES', 'ndvi')
+ self.kp_path = config.get('FILES', 'kp')
+ self.land_path = config.get('FILES', 'landuse')
+ self.soil_path = config.get('FILES', 'solo')
+ self.outpath = config.get('FILES', 'output')
+ self.sampleLocs = config.get('FILES', 'samples')
- # Create slope map based on DEM
- self.S = slope(self.dem)
-
- # Create outflow points map based on ldd
- self.pits = pit(self.ldd)
-
- # Creat cacthment area basins
- subbasins = catchment(self.ldd, self.pits)
-
- # Initializa Tss report at sample locations or pits
- self.TssFileRun = TimeoutputTimeseries(self.OutTssRun, self, self.sampleLocs, noHeader=True)
+ # Set clone
+ pcr.setclone(self.clone_file)
- # Read min and max ndvi
- self.ndvi_min = scalar(readmap(self.ndvi_path+self.ndviMinFile))
- self.ndvi_max = scalar(readmap(self.ndvi_path+self.ndviMaxFile))
-
- # Compute min and max sr
- self.sr_min = sr_calc(pcr,self,self.ndvi_min)
- self.sr_max = sr_calc(pcr,self,self.ndvi_max)
-
- # Read soil atributes
- solo = self.readmap((self.soil_path))
- self.Kr = lookupscalar(self.KrTable,solo) #coeficiente de condutividade hidráulica
- self.dg = lookupscalar(self.dgTable,solo) #densidade do solo
- self.Zr = lookupscalar(self.ZrTable,solo) # profundidade da zona radicular [cm]
- self.TUsat = lookupscalar(self.TsatTable,solo)*self.dg*self.Zr*10 # umidade para saturação da primeira camada [mm]
- self.TUr_ini = (self.TUsat)*(self.ftur_ini) # teor de umidade inicial da zona radicular [mm]
- self.TUw = lookupscalar(self.TwTable,solo)*self.dg*self.Zr*10 # ponto de mucrha do solo [mm]
- self.TUcc = lookupscalar(self.TccTable,solo)*self.dg*self.Zr*10 # capacidade de campo [mm]
- self.Tpor = lookupscalar(self.Tporosidade,solo) # porosidade [%]
- self.EB_ini = self.EBini # escoamento basico inicial [mm]
- self.EB_lim = self.EBlim # limite para condicao de escoamento basico [mm]
- self.TUs_ini = self.Tusini # teor de umidade inicial da camada saturada [mm]
-
- # steps
- self.steps = totalSteps(startDate,endDate)
- self.lastStep = steps[1]
-
- # Conditions for t = first loop
- self.TUrprev = self.TUr_ini
- self.TUsprev = self.TUs_ini
- self.EBprev = self.EB_ini
- self.TUr = self.TUr_ini
- self.TUs = self.TUs_ini
- self.EB = self.EB_ini
- self.Qini = scalar(0)
- self.Qprev = self.Qini
-
- # initialize first landuse map
- self.landuse = self.readmap(self.land_path + self.coverPrefix)
- self.landuse_ant = self.landuse
-
- def dynamic(self):
- t = self.currentStep
- #print(t)
- print("Tempo: "+str(t), flush=True)
-
- # Read NDVI
- try:
- NDVI = self.readmap(self.ndvi_path + self.ndviPrefix)
- self.ndvi_ant = NDVI
-
- except RuntimeError:
- NDVI = self.ndvi_ant
-
- # Read Landuse Maps
- try:
- self.landuse = self.readmap(self.land_path + self.coverPrefix)
- self.landuse_ant = self.landuse
-
- except RuntimeError:
- self.landuse = self.landuse_ant
-
- # Read precipitation maps
- precipitation = pcr.scalar(self.readmap(self.prec_path + self.precPrefix))
-
- # Read potential evapotranspiration
- ETp = pcr.scalar(self.readmap(self.etp_path + self.etpPrefix))
-
- # Read Kp
- Kp = pcr.scalar(self.readmap(self.kp_path + self.kpPrefix))
-
- # Number of rainy days
- month = ((t-1)%12)+1
- rainyDays = lookupscalar(self.rainyDaysTable, month)
-
- # Read Landuse attributes
- n_manning = lookupscalar(self.manningTable,self.landuse)
- Av = lookupscalar(self.avTable,self.landuse)
- Ao = lookupscalar(self.aoTable,self.landuse)
- As = lookupscalar(self.asTable,self.landuse)
- Ai = lookupscalar(self.aiTable,self.landuse)
- self.kc_min = lookupscalar(self.KcminTable,self.landuse)
- self.kc_max = lookupscalar(self.KcmaxTable,self.landuse)
-
- ######### compute interception #########
- SR = sr_calc(pcr,self,NDVI)
- FPAR = fpar_calc(pcr, self, self.fpar_min, self.fpar_max, SR, self.sr_min, self.sr_max)
- LAI = lai_function(pcr, self, FPAR, self.fpar_max, self.lai_max)
- Id, Ir, Iv, I = Interception_function(pcr, self, self.alfa, LAI, precipitation, rainyDays, Av)
+ # Read temporal filenames prefix
+ self.etpPrefix = config.get('FILES', 'etpFilePrefix')
+ self.precPrefix = config.get('FILES', 'precFilePrefix')
+ self.ndviPrefix = config.get('FILES', 'ndviFilePrefix')
+ self.ndviMaxFile = config.get('FILES', 'ndvimaxPrefix')
+ self.ndviMinFile = config.get('FILES', 'ndviminPrefix')
+ self.kpPrefix = config.get('FILES', 'kpFilePrefix')
+ self.coverPrefix = config.get('FILES', 'landuseFilePrefix')
+
+ # Read text lookuptables from config file
+ self.rainyDaysTable = config.get('PARAMETERS', 'rainydays')
+ self.aiTable = config.get('PARAMETERS', 'a_i')
+ self.aoTable = config.get('PARAMETERS', 'a_o')
+ self.asTable = config.get('PARAMETERS', 'a_s')
+ self.avTable = config.get('PARAMETERS', 'a_v')
+ self.manningTable = config.get('PARAMETERS', 'manning')
+ self.dgTable = config.get('PARAMETERS', 'dg')
+ self.KrTable = config.get('PARAMETERS', 'kr')
+ self.TccTable = config.get('PARAMETERS', 'capCampo')
+ self.Tporosidade = config.get('PARAMETERS', 'porosidade')
+ self.TsatTable = config.get('PARAMETERS', 'saturacao')
+ self.TwTable = config.get('PARAMETERS', 'pontomurcha')
+ self.ZrTable = config.get('PARAMETERS', 'zr')
+ self.KcminTable = config.get('PARAMETERS', 'kcmin')
+ self.KcmaxTable = config.get('PARAMETERS', 'kcmax')
+
+ # Area da celula #pendente = calculo automatico da area da celula
+ self.A = config.getfloat('GRID', 'grid')
+
+ # Read calibration parameters from config file
+ self.alfa = config.getfloat('CALIBRATION', 'alfa')
+ self.b = config.getfloat('CALIBRATION', 'b')
+ self.w1 = config.getfloat('CALIBRATION', 'w1')
+ self.w2 = config.getfloat('CALIBRATION', 'w2')
+ self.w3 = config.getfloat('CALIBRATION', 'w3')
+ self.RCD = config.getfloat('CALIBRATION', 'rcd')
+ self.f = config.getfloat('CALIBRATION', 'f')
+ self.alfa_gw = config.getfloat('CALIBRATION', 'alfa_gw')
+ self.x = config.getfloat('CALIBRATION', 'x')
+
+ # Read soil conditions from config file
+ # teor de umidade inicial da zona radicular (fracao do teor de saturacao)
+ self.ftur_ini = config.getfloat('INITIAL SOIL CONDITIONS', 'ftur_ini')
+ # escoamento basico inicial
+ self.EBini = scalar(config.getfloat('INITIAL SOIL CONDITIONS', 'eb_ini'))
+ # limite para escoamento basico
+ self.EBlim = scalar(config.getfloat('INITIAL SOIL CONDITIONS', 'eb_lim'))
+ # teor de umidade inicial da zona saturada
+ self.Tusini = scalar(config.getfloat('INITIAL SOIL CONDITIONS', 'tus_ini'))
+
+ # constantes
+ self.fpar_max = config.getfloat('CONSTANT', 'fpar_max')
+ self.fpar_min = config.getfloat('CONSTANT', 'fpar_min')
+ self.lai_max = config.getfloat('CONSTANT', 'lai_max')
+ self.I_i = config.getfloat('CONSTANT', 'i_imp')
+
+ # # Initialize time series output
+ self.OutTssRun= ('outRun')
+
+ # Report file
+ # name
+ self.timeStamp = str((time.strftime("%Y%m%d_%H%M%S",time.localtime(t1))))
+ # header
+ self.t_round = str(time.strftime("%Y %m %d %H:%M:%S",time.localtime(t1)))
+
+ # Get Tif file Reference
+ self.ref = getRefInfo(self, self.demTif)
- print("Interceptação...ok", flush=True)
- ######### Compute Evapotranspiration #########
-
- Kc_1 = kc_calc(pcr, self, NDVI, self.ndvi_min, self.ndvi_max, self.kc_min, self.kc_max)
- # condicao do kc, se NDVI < 1.1NDVI_min, kc = kc_min
- kc_cond1 = scalar(NDVI < 1.1*self.ndvi_min)
- kc_cond2 = scalar(NDVI > 1.1*self.ndvi_min)
- Kc = pcr.scalar((kc_cond2*Kc_1) + (kc_cond1*self.kc_min))
- Ks = pcr.scalar(Ks_calc(pcr, self, self.TUr, self.TUw, self.TUcc))
-
- # Vegetated area
- self.ET_av = ETav_calc(pcr, self, ETp, Kc, Ks)
-
- # Impervious area
- # ET impervious area = Interception of impervious area
- # condicao leva em conta a chuva igual a zero
- # mascara: chuva = 0 -> 0, chuva <> 0 -> 1
- cond = pcr.scalar((precipitation != 0))
- self.ET_ai = (self.I_i*cond)
-
- # Open water
- self.ET_ao = ETao_calc(pcr, self, ETp, Kp, precipitation, Ao)
- #self.report(ET_ao,self.outpath+'ETao')
-
- # Bare soil
- self.ET_as = ETas_calc(pcr, self, ETp, self.kc_min, Ks)
- self.ETr = (Av*self.ET_av) + (Ai*self.ET_ai) + (Ao*self.ET_ao) + (As*self.ET_as)
+ def initial(self):
+ """ """
+ # Read dem file
+ self.dem = readmap(self.dem_file)
- print("Evapotranspiração...ok", flush=True)
-
- ######### Surface Runoff #########
- Pdm = (precipitation/rainyDays)
- Ch = Ch_calc(pcr, self, self.TUr, self.dg, self.Zr, self.Tpor, self.b)
- Cper = Cper_calc(pcr, self, self.TUw, self.dg, self.Zr, self.S, n_manning, self.w1, self.w2, self.w3)
- Aimp, Cimp = Cimp_calc(pcr, self, Ao, Ai)
- Cwp = Cwp_calc(pcr, self, Aimp, Cper, Cimp)
- Csr = Csr_calc(pcr, self, Cwp, Pdm, self.RCD)
-
- self.ES = ES_calc(pcr, self, Csr, Ch, precipitation, I, Ao, self.ET_ao)
+ # Generate the local drain direction map on basis of the elevation map
+ # self.ldd = lddcreate(self.dem, 1e31, 1e31, 1e31, 1e31)
+ self.ldd = pcr.ldd(readmap(self.ldd_file))
- print("Escoamento Superficial...ok", flush=True)
- ######### Lateral Flow #########
- self.LF = LF_calc(pcr, self, self.f, self.Kr, self.TUr, self.TUsat)
+ # Create slope map based on DEM
+ self.S = slope(self.dem)
- print("Fluxo Lateral...ok", flush=True)
- ######### Recharge Flow #########
- self.REC = REC_calc(pcr, self, self.f, self.Kr, self.TUr, self.TUsat)
+ # Create outflow points map based on ldd
+ self.pits = pit(self.ldd)
- print("Recarga...ok", flush=True)
- ######### Base Flow #########
- # reportTif(self, self.ref, self.EBprev, 'EBprev', self.outpath, din = 1)
- # reportTif(self, self.ref, self.TUs, 'TUs2', self.outpath, din = 1)
+ # Creat cacthment area basins
+ subbasins = catchment(self.ldd, self.pits)
- self.EB = EB_calc(pcr, self, self.EBprev, self.alfa_gw, self.REC, self.TUs, self.EB_lim)
- self.EBprev = self.EB
- # reportTif(self, self.ref, self.EB, 'EB', self.outpath, din = 1)
-
- ######### Soil Balance #########
- self.TUr = TUr_calc(pcr, self, self.TUrprev, precipitation, I, self.ES, self.LF, self.REC, self.ETr, Ao, self.TUsat)
- self.TUs = TUs_calc(pcr, self, self.TUsprev, self.REC, self.EB)
- self.TUrprev = self.TUr
-
- self.TUsprev = self.TUs
- print("Balanço hídrico do solo...ok", flush=True)
- ######### Compute Runoff ########
- days = daysOfMonth(startDate,t)
- c = days*24*3600
-
- self.Qtot = ((self.ES + self.LF + self.EB)) # [mm]
- self.Qtotvol = self.Qtot*self.A*0.001/c # [m3/s]
-
- self.Qt = accuflux(self.ldd, self.Qtotvol)
-
- self.runoff = self.x*self.Qprev + (1-self.x)*self.Qt
- self.Qprev = self.runoff
+ # Initializa Tss report at sample locations or pits
+ self.TssFileRun = TimeoutputTimeseries(self.OutTssRun, self, self.sampleLocs, noHeader=True)
- print("Vazão...ok", flush=True)
-
- # # Create tss files
- os.chdir(self.outpath)
-
- # Lista de arquivos para exportar - ordenados igual a lista de entrada
- filesList = [self.EB, self. ES, self.ETr, self.LF, I, self.REC, self.TUr, self.runoff, self.Qtot, self.REC]
- for i in range(len(filesList)):
- if genFilesDic[genFilesList[i]] == 'True':
- reportTif(self, self.ref, filesList[i], str(genFilesList[i]), self.outpath, din = 1)
- else:
- pass
-
- print("Finalizando ciclo "+str(t) + " de "+ str(self.lastStep), flush=True)
+ # Read min and max ndvi
+ self.ndvi_min = scalar(readmap(self.ndvi_path+self.ndviMinFile))
+ self.ndvi_max = scalar(readmap(self.ndvi_path+self.ndviMaxFile))
+
+ # Compute min and max sr
+ self.sr_min = sr_calc(pcr,self,self.ndvi_min)
+ self.sr_max = sr_calc(pcr,self,self.ndvi_max)
+
+ # Read soil atributes
+ solo = self.readmap((self.soil_path))
+ self.Kr = lookupscalar(self.KrTable,solo) #coeficiente de condutividade hidráulica
+ self.dg = lookupscalar(self.dgTable,solo) #densidade do solo
+ self.Zr = lookupscalar(self.ZrTable,solo) # profundidade da zona radicular [cm]
+ self.TUsat = lookupscalar(self.TsatTable,solo)*self.dg*self.Zr*10 # umidade para saturação da primeira camada [mm]
+ self.TUr_ini = (self.TUsat)*(self.ftur_ini) # teor de umidade inicial da zona radicular [mm]
+ self.TUw = lookupscalar(self.TwTable,solo)*self.dg*self.Zr*10 # ponto de mucrha do solo [mm]
+ self.TUcc = lookupscalar(self.TccTable,solo)*self.dg*self.Zr*10 # capacidade de campo [mm]
+ self.Tpor = lookupscalar(self.Tporosidade,solo) # porosidade [%]
+ self.EB_ini = self.EBini # escoamento basico inicial [mm]
+ self.EB_lim = self.EBlim # limite para condicao de escoamento basico [mm]
+ self.TUs_ini = self.Tusini # teor de umidade inicial da camada saturada [mm]
+
+ # steps
+ self.steps = totalSteps(startDate,endDate)
+ self.lastStep = steps[1]
+
+ # Conditions for t = first loop
+ self.TUrprev = self.TUr_ini
+ self.TUsprev = self.TUs_ini
+ self.EBprev = self.EB_ini
+ self.TUr = self.TUr_ini
+ self.TUs = self.TUs_ini
+ self.EB = self.EB_ini
+ self.Qini = scalar(0)
+ self.Qprev = self.Qini
+
+ # initialize first landuse map
+ self.landuse = self.readmap(self.land_path + self.coverPrefix)
+ self.landuse_ant = self.landuse
+
+ def dynamic(self):
+ """ """
+ t = self.currentStep
+ #print(t)
+ print("Tempo: "+str(t), flush=True)
+
+ # Read NDVI
+ try:
+ NDVI = self.readmap(self.ndvi_path + self.ndviPrefix)
+ self.ndvi_ant = NDVI
+
+ except RuntimeError:
+ NDVI = self.ndvi_ant
+
+ # Read Landuse Maps
+ try:
+ self.landuse = self.readmap(self.land_path + self.coverPrefix)
+ self.landuse_ant = self.landuse
+
+ except RuntimeError:
+ self.landuse = self.landuse_ant
+
+ # Read precipitation maps
+ precipitation = pcr.scalar(self.readmap(self.prec_path + self.precPrefix))
+
+ # Read potential evapotranspiration
+ ETp = pcr.scalar(self.readmap(self.etp_path + self.etpPrefix))
+
+ # Read Kp
+ Kp = pcr.scalar(self.readmap(self.kp_path + self.kpPrefix))
+
+ # Number of rainy days
+ month = ((t-1)%12)+1
+ rainyDays = lookupscalar(self.rainyDaysTable, month)
+
+ # Read Landuse attributes
+ n_manning = lookupscalar(self.manningTable,self.landuse)
+ Av = lookupscalar(self.avTable,self.landuse)
+ Ao = lookupscalar(self.aoTable,self.landuse)
+ As = lookupscalar(self.asTable,self.landuse)
+ Ai = lookupscalar(self.aiTable,self.landuse)
+ self.kc_min = lookupscalar(self.KcminTable,self.landuse)
+ self.kc_max = lookupscalar(self.KcmaxTable,self.landuse)
+
+ ######### compute interception #########
+ SR = sr_calc(pcr,self,NDVI)
+ FPAR = fpar_calc(pcr, self, self.fpar_min, self.fpar_max, SR, self.sr_min, self.sr_max)
+ LAI = lai_function(pcr, self, FPAR, self.fpar_max, self.lai_max)
+ Id, Ir, Iv, I = Interception_function(pcr, self, self.alfa, LAI, precipitation, rainyDays, Av)
+
+ print("Interceptação...ok", flush=True)
+ ######### Compute Evapotranspiration #########
+
+ Kc_1 = kc_calc(pcr, self, NDVI, self.ndvi_min, self.ndvi_max, self.kc_min, self.kc_max)
+ # condicao do kc, se NDVI < 1.1NDVI_min, kc = kc_min
+ kc_cond1 = scalar(NDVI < 1.1*self.ndvi_min)
+ kc_cond2 = scalar(NDVI > 1.1*self.ndvi_min)
+ Kc = pcr.scalar((kc_cond2*Kc_1) + (kc_cond1*self.kc_min))
+ Ks = pcr.scalar(Ks_calc(pcr, self, self.TUr, self.TUw, self.TUcc))
+
+ # Vegetated area
+ self.ET_av = ETav_calc(pcr, self, ETp, Kc, Ks)
+
+ # Impervious area
+ # ET impervious area = Interception of impervious area
+ # condicao leva em conta a chuva igual a zero
+ # mascara: chuva = 0 -> 0, chuva <> 0 -> 1
+ cond = pcr.scalar((precipitation != 0))
+ self.ET_ai = (self.I_i*cond)
+
+ # Open water
+ self.ET_ao = ETao_calc(pcr, self, ETp, Kp, precipitation, Ao)
+ #self.report(ET_ao,self.outpath+'ETao')
+
+ # Bare soil
+ self.ET_as = ETas_calc(pcr, self, ETp, self.kc_min, Ks)
+ self.ETr = (Av*self.ET_av) + (Ai*self.ET_ai) + (Ao*self.ET_ao) + (As*self.ET_as)
+
+ print("Evapotranspiração...ok", flush=True)
+
+ ######### Surface Runoff #########
+ Pdm = (precipitation/rainyDays)
+ Ch = Ch_calc(pcr, self, self.TUr, self.dg, self.Zr, self.Tpor, self.b)
+ Cper = Cper_calc(pcr, self, self.TUw, self.dg, self.Zr, self.S, n_manning, self.w1, self.w2, self.w3)
+ Aimp, Cimp = Cimp_calc(pcr, self, Ao, Ai)
+ Cwp = Cwp_calc(pcr, self, Aimp, Cper, Cimp)
+ Csr = Csr_calc(pcr, self, Cwp, Pdm, self.RCD)
+
+ self.ES = ES_calc(pcr, self, Csr, Ch, precipitation, I, Ao, self.ET_ao)
+
+ print("Escoamento Superficial...ok", flush=True)
+ ######### Lateral Flow #########
+ self.LF = LF_calc(pcr, self, self.f, self.Kr, self.TUr, self.TUsat)
+
+ print("Fluxo Lateral...ok", flush=True)
+ ######### Recharge Flow #########
+ self.REC = REC_calc(pcr, self, self.f, self.Kr, self.TUr, self.TUsat)
+
+ print("Recarga...ok", flush=True)
+ ######### Base Flow #########
+ # reportTif(self, self.ref, self.EBprev, 'EBprev', self.outpath, din = 1)
+ # reportTif(self, self.ref, self.TUs, 'TUs2', self.outpath, din = 1)
+
+ self.EB = EB_calc(pcr, self, self.EBprev, self.alfa_gw, self.REC, self.TUs, self.EB_lim)
+ self.EBprev = self.EB
+ # reportTif(self, self.ref, self.EB, 'EB', self.outpath, din = 1)
+
+ ######### Soil Balance #########
+ self.TUr = TUr_calc(pcr, self, self.TUrprev, precipitation, I, self.ES, self.LF, self.REC, self.ETr, Ao, self.TUsat)
+ self.TUs = TUs_calc(pcr, self, self.TUsprev, self.REC, self.EB)
+ self.TUrprev = self.TUr
+
+ self.TUsprev = self.TUs
+ print("Balanço hídrico do solo...ok", flush=True)
+ ######### Compute Runoff ########
+ days = daysOfMonth(startDate,t)
+ c = days*24*3600
+
+ self.Qtot = ((self.ES + self.LF + self.EB)) # [mm]
+ self.Qtotvol = self.Qtot*self.A*0.001/c # [m3/s]
+
+ self.Qt = accuflux(self.ldd, self.Qtotvol)
+
+ self.runoff = self.x*self.Qprev + (1-self.x)*self.Qt
+ self.Qprev = self.runoff
+
+ print("Vazão...ok", flush=True)
+
+ # # Create tss files
+ os.chdir(self.outpath)
+
+ # Lista de arquivos para exportar - ordenados igual a lista de entrada
+ filesList = [self.EB, self. ES, self.ETr, self.LF, I, self.REC, self.TUr, self.runoff, self.Qtot, self.REC]
+ for i in range(len(filesList)):
+ if genFilesDic[genFilesList[i]] == 'True':
+ reportTif(self, self.ref, filesList[i], str(genFilesList[i]), self.outpath, din = 1)
+ else:
+ pass
+
+ print("Finalizando ciclo "+str(t) + " de "+ str(self.lastStep), flush=True)
steps = totalSteps(startDate,endDate)
@@ -418,7 +466,6 @@ def dynamic(self):
myModel = Modelo()
dynamicModel = DynamicFramework(myModel,lastTimeStep=end, firstTimestep=start)
dynamicModel.run()
-import time
tempoExec = time.time() - t1
print("Tempo de execução: {:.2f} segundos".format(tempoExec))
print("Fim", flush=True)
\ No newline at end of file
diff --git a/evapotranspiration.py b/evapotranspiration.py
index 63edb2d..9e00544 100644
--- a/evapotranspiration.py
+++ b/evapotranspiration.py
@@ -1,24 +1,93 @@
+# coding=utf-8
+
########## Evapotranspiration Module ##########
# - Function that returns Ks for evapotranspiration of vegetated area
-def Ks_calc(pcr, self, TUr, TUw, TUcc):
+def Ks_calc(self, pcr, TUr, TUw, TUcc):
+ """
+ :param pcr:
+ :pcr type:
+
+ :param TUr:
+ :TUr type:
+
+ :param TUw:
+ :TUw type:
+
+ :param TUcc:
+ :TUcc type:
+
+ :returns:
+ :rtype:
+ """
ks_cond = pcr.scalar(TUr > TUw) #caso TUr < TUw, (false, ks = 0)
# Multiply (TUr - TUw) to avoid negative ln
Ks = (pcr.ln((TUr - TUw)*ks_cond + 1))/(pcr.ln(TUcc - TUw + 1))
return Ks
# - Function that returns evapotranspiration of vegetated area
-def ETav_calc(pcr, self, ETp, Kc, Ks):
+def ETav_calc(self, pcr, ETp, Kc, Ks):
+ """
+ :param pcr:
+ :pcr type:
+
+ :param ETp:
+ :ETp type:
+
+ :param Kc:
+ :Kc type:
+
+ :param Ks:
+ :Ks type:
+
+ :returns:
+ :rtype:
+ """
ETav = ETp * Kc * Ks
return ETav
# - Function that returns Kp for evapotranspiration of open water area
-def Kp_calc(pcr, self, B, U_2, UR):
+def Kp_calc(self, pcr, B, U_2, UR):
+ """
+ :param pcr:
+ :pcr type:
+
+ :param B:
+ :B type:
+
+ :param U_2:
+ :U_2 type:
+
+ :param UR:
+ :UR type:
+
+ :returns:
+ :rtype:
+ """
Kp = 0.482+0.024*pcr.ln(B)-0.000376*U_2 + 0.0045*UR
return Kp
# - Function that returns Ks for evapotranspiration of water area
-def ETao_calc(pcr, self, ETp, Kp, prec, Ao):
+def ETao_calc(self, pcr, ETp, Kp, prec, Ao):
+ """
+ :param pcr:
+ :pcr type:
+
+ :param ETp:
+ :ETp type:
+
+ :param Kp:
+ :Kp type:
+
+ :param prec:
+ :prec type:
+
+ :param Ao:
+ :Ao type:
+
+ :returns:
+ :rtype:
+ """
# condition for pixel of water
cond1 = pcr.scalar(Ao == 1)
@@ -32,7 +101,23 @@ def ETao_calc(pcr, self, ETp, Kp, prec, Ao):
return ETao
# - Function that returns Ks for evapotranspiration of bare soil area
-def ETas_calc(pcr, self, ETp, kc_min, Ks):
+def ETas_calc(self, pcr, ETp, kc_min, Ks):
+ """
+ :param pcr:
+ :pcr type:
+
+ :param ETp:
+ :ETp type:
+
+ :param kc_min:
+ :kc_min type:
+
+ :param Ks:
+ :Ks type:
+
+ :returns:
+ :rtype:
+ """
cond = 1*pcr.scalar(Ks != 0) #caso ks seja diferente de 0
ETas = (ETp*kc_min*Ks*cond)
return ETas
\ No newline at end of file
diff --git a/interception.py b/interception.py
index 2654788..afe4010 100644
--- a/interception.py
+++ b/interception.py
@@ -1,28 +1,122 @@
+# coding=utf-8
+
########## Interception Module ##########
# - Function that returns srmin and srmax
-def sr_calc(pcr, self, NDVI):
+def sr_calc(self, pcr, NDVI):
+ """
+ :param pcr:
+ :pcr type:
+
+ :param NDVI:
+ :NDVI type:
+
+ :returns:
+ :rtype:
+ """
SR = (1+NDVI)/(1-NDVI)
return SR
# - Function that returns Kc
-def kc_calc(pcr, self, NDVI, ndvi_min, ndvi_max, kc_min, kc_max):
+def kc_calc(self, pcr, NDVI, ndvi_min, ndvi_max, kc_min, kc_max):
+ """
+ :param pcr:
+ :pcr type:
+
+ :param NDVI:
+ :NDVI type:
+
+ :param ndvi_min:
+ :ndvi_min type:
+
+ :param ndvi_max:
+ :ndvi_max type:
+
+ :param kc_min:
+ :kc_min type:
+
+ :param kc_max:
+ :kc_max type:
+
+ :returns:
+ :rtype:
+ """
Kc = kc_min+((kc_max-kc_min)*((NDVI - pcr.scalar(ndvi_min))/(pcr.scalar(ndvi_max) - pcr.scalar(ndvi_min))))
return Kc
# - Function that returns fpar
-def fpar_calc(pcr, self,fpar_min, fpar_max, SR, sr_min, sr_max):
+def fpar_calc(self, pcr, fpar_min, fpar_max, SR, sr_min, sr_max):
+ """
+ :param pcr:
+ :pcr type:
+
+ :param fpar_min:
+ :fpar_min type:
+
+ :param fpar_max:
+ :fpar_max type:
+
+ :param SR:
+ :SR type:
+
+ :param sr_min:
+ :sr_min type:
+
+ :param sr_max:
+ :sr_max type:
+
+ :returns:
+ :rtype:
+ """
fpar_comp = (((SR-sr_min)*(fpar_max-fpar_min)/(sr_max-sr_min))+fpar_min)
FPAR = pcr.min(fpar_comp, fpar_max)
return FPAR
# - Function that returns LAI
-def lai_function(pcr, self, FPAR, fpar_max, lai_max):
+def lai_function(self, pcr, FPAR, fpar_max, lai_max):
+ """
+ :param pcr:
+ :pcr type:
+
+ :param FPAR:
+ :FPAR type:
+
+ :param fpar_max:
+ :fpar_max type:
+
+ :param lai_max:
+ :lai_max type:
+
+ :returns:
+ :rtype:
+ """
LAI = lai_max*((pcr.log10(1-FPAR))/(pcr.log10(1-fpar_max)))
return LAI
# - Function that returns Interception
-def Interception_function(pcr, self, alfa, LAI, precipitation, rainy_days, a_v):
+def Interception_function(self, pcr, alfa, LAI, precipitation, rainy_days, a_v):
+ """
+ :param pcr:
+ :pcr type:
+
+ :param alfa:
+ :alfa type:
+
+ :param LAI:
+ :LAI type:
+
+ :param precipitation:
+ :precipitation type:
+
+ :param rainy_days:
+ :rainy_days type:
+
+ :param a_v:
+ :a_v type:
+
+ :returns:
+ :rtype:
+ """
# condition of precipitation, to divide by non zero number (missing value)
cond1 = pcr.scalar((precipitation != 0))
cond2 = pcr.scalar((precipitation == 0))
diff --git a/soil.py b/soil.py
index 5fef452..b7b765b 100644
--- a/soil.py
+++ b/soil.py
@@ -1,15 +1,74 @@
+# coding=utf-8
+
########## Lateral Flow ##########
-def LF_calc(pcr, self, f, Kr, TUr, TUsat):
+def LF_calc(self, pcr, f, Kr, TUr, TUsat):
+ """
+ :param pcr:
+ :pcr type:
+
+ :param f:
+ :f type:
+
+ :param TUr:
+ :TUr type:
+
+ :param TUsat:
+ :TUsat type:
+
+ :returns:
+ :rtype:
+ """
LF = f*Kr*((TUr/TUsat)**2)
return LF
########## Recharge ##########
-def REC_calc(pcr, self, f, Kr, TUr, TUsat):
+def REC_calc(self, pcr, f, Kr, TUr, TUsat):
+ """
+ :param pcr:
+ :pcr type:
+
+ :param f:
+ :f type:
+
+ :param Kr:
+ :Kr type:
+
+ :param TUr:
+ :TUr type:
+
+ :param TUsat:
+ :TUsat type:
+
+ :returns:
+ :rtype:
+ """
REC = (1-f)*Kr*((TUr/TUsat)**2)
return REC
########## Base Flow ##########
-def EB_calc(pcr, self, EB_prev, alfaS, REC, TUs, EB_lim):
+def EB_calc(self, pcr, EB_prev, alfaS, REC, TUs, EB_lim):
+ """
+ :param pcr:
+ :pcr type:
+
+ :param EB_prev:
+ :EB_prev type:
+
+ :param alfaS:
+ :alfaS type:
+
+ :param REC:
+ :REC type:
+
+ :param TUs:
+ :TUs type:
+
+ :param EB_lim:
+ :EB_lim type:
+
+ :returns:
+ :rtype:
+ """
# limit condition for base flow
cond = pcr.scalar(TUs > EB_lim)
EB = ((EB_prev*((pcr.exp(1))**-alfaS))+(1-((pcr.exp(1))**-alfaS))*REC)*cond
@@ -17,7 +76,41 @@ def EB_calc(pcr, self, EB_prev, alfaS, REC, TUs, EB_lim):
########## Soil Balance ##########
# First soil layer
-def TUr_calc(pcr, self, TUrprev, P, I, ES, LF, REC, ETr, Ao, Tsat):
+def TUr_calc(self, pcr, TUrprev, P, I, ES, LF, REC, ETr, Ao, Tsat):
+ """
+ :param pcr:
+ :pcr type:
+
+ :param TUrprev:
+ :TUrprev type:
+
+ :param P:
+ :P type:
+
+ :param I:
+ :I type:
+
+ :param ES:
+ :ES type:
+
+ :param LF:
+ :LF type:
+
+ :param REC:
+ :REC type:
+
+ :param ETr:
+ :ETr type:
+
+ :param Ao:
+ :Ao type:
+
+ :param Tsat:
+ :Tsat type:
+
+ :returns:
+ :rtype:
+ """
# condition for pixel of water, if Ao different of 1 (not water)
condw1 = pcr.scalar(Ao != 1)
# soil balance
@@ -29,7 +122,23 @@ def TUr_calc(pcr, self, TUrprev, P, I, ES, LF, REC, ETr, Ao, Tsat):
return TUr
# Second soil layer
-def TUs_calc(pcr, self, TUsprev, REC, EB):
+def TUs_calc(self, pcr, TUsprev, REC, EB):
+ """
+ :param pcr:
+ :pcr type:
+
+ :param TUsprev:
+ :TUsprev type:
+
+ :param REC:
+ :REC type:
+
+ :param EB:
+ :EB type:
+
+ :returns:
+ :rtype:
+ """
# soil balance
balance = TUsprev + REC - EB
TUs = balance
diff --git a/surface_runoff.py b/surface_runoff.py
index 695ff9b..e961ba5 100644
--- a/surface_runoff.py
+++ b/surface_runoff.py
@@ -1,35 +1,160 @@
+# coding=utf-8
+
########## Surface runoff ##########
# - Function that returns Ch
-def Ch_calc(pcr, self, TUr, dg, Zr, tpor, b):
+def Ch_calc(self, pcr, TUr, dg, Zr, tpor, b):
+ """
+ :param pcr:
+ :pcr type:
+
+ :param TUr:
+ :TUr type:
+
+ :param dg:
+ :dg type:
+
+ :param Zr:
+ :Zr type:
+
+ :param tpor:
+ :tpor type:
+
+ :param b:
+ :b type:
+
+ :returns:
+ :rtype:
+ """
tur = TUr/(dg*Zr*10) # [%] soil moisture
Ch = (tur/tpor)**b
return Ch
# - Function that returns Cper
-def Cper_calc(pcr, self, TUw, dg, Zr, S, manning, w1, w2, w3):
+def Cper_calc(self, pcr, TUw, dg, Zr, S, manning, w1, w2, w3):
+ """
+ :param pcr:
+ :pcr type:
+
+ :param TUw:
+ :TUw type:
+
+ :param dg:
+ :dg type:
+
+ :param Zr:
+ :Zr type:
+
+ :param S:
+ :S type:
+
+ :param manning:
+ :manning type:
+
+ :param w1:
+ :w1 type:
+
+ :param w2:
+ :w2 type:
+
+ :param w3:
+ :w3 type:
+
+ :returns:
+ :rtype:
+ """
tuw = TUw/(dg*Zr*10) # [%] soil wilting point
Cper = w1*(0.02/manning) + w2*(tuw/(1-tuw)) + w3*((S/(10+S)))
return Cper
# - Function that returns Cimp
-def Cimp_calc(pcr, self, ao, ai):
+def Cimp_calc(self, pcr, ao, ai):
+ """
+ :param pcr:
+ :pcr type:
+
+ :param ao:
+ :ao type:
+
+ :param ai:
+ :ai type:
+
+ :returns:
+ :rtype:
+ """
Aimp = ao + ai
Cimp = (0.09*pcr.exp((2.4*Aimp)))
return Aimp, Cimp
# - Function that returns Cwp
-def Cwp_calc(pcr, self, Aimp, Cper, Cimp):
+def Cwp_calc(self, pcr, Aimp, Cper, Cimp):
+ """
+ :param pcr:
+ :pcr type:
+
+ :param Aimp:
+ :Aimp type:
+
+ :param Cper:
+ :Cper type:
+
+ :param Cimp:
+ :Cimp type:
+
+ :returns:
+ :rtype:
+ """
Cwp = (1-Aimp)*Cper+Aimp*Cimp
return Cwp
# - Function that returns Csr
-def Csr_calc(pcr, self, Cwp, P_24, RCD):
+def Csr_calc(self, pcr, Cwp, P_24, RCD):
+ """
+ :param pcr:
+ :pcr type:
+
+ :param Cwp:
+ :Cwp type:
+
+ :param P_24:
+ :P_24 type:
+
+ :param RCD:
+ :RCD type:
+
+ :returns:
+ :rtype:
+ """
Csr = (Cwp*P_24)/(Cwp*P_24 - RCD*Cwp + RCD)
return Csr
# - Function that returns Surface runoff
-def ES_calc(pcr, self, Csr, Ch, prec, I, Ao, ETao):
+def ES_calc(self, pcr, Csr, Ch, prec, I, Ao, ETao):
+ """
+ :param pcr:
+ :pcr type:
+
+ :param Csr:
+ :Csr type:
+
+ :param Ch:
+ :Ch type:
+
+ :param prec:
+ :prec type:
+
+ :param I:
+ :I type:
+
+ :param Ao:
+ :Ao type:
+
+ :param ETao:
+ :ETao type:
+
+ :returns:
+ :rtype:
+ """
# condition for pixel of water
cond1 = pcr.scalar(Ao == 1)
# condition for positive value for (prec - ETao)