Berdasarkan isu #89: pemodelan NRECA
Referensi isu:
Deskripsi Permasalahan:
MSTOR
($mm$), Groundwater Storage GSTOR
($mm$), PSUB
, GWF
, Crop Factor CF
, C
, Area AREA
($m^2$),Strategi Penyelesaian:
_
dalam penamaan fungsi).model_NRECA()
dengan kejelasan parameter sehingga bisa dikembangkan lebih lanjut untuk kalibrasi model atau mencari parameter terbaik.Catatan:
PREP
($mm$) dan evapotranspirasi PET
($mm$) setiap bulannya.import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from google.colab import files
try:
pd.read_excel('NRECA_sample.xlsx')
except:
dataset_path = list(files.upload().keys())[0]
dataset = pd.read_excel(dataset_path, header=0, index_col=0, parse_dates=True)
dataset.info()
dataset.head()
<class 'pandas.core.frame.DataFrame'> DatetimeIndex: 120 entries, 1999-01-01 to 2008-12-01 Data columns (total 2 columns): PRECIP 120 non-null float64 PET 120 non-null float64 dtypes: float64(2) memory usage: 2.8 KB
PRECIP | PET | |
---|---|---|
1999-01-01 | 507.000000 | 142.63 |
1999-02-01 | 374.228918 | 128.84 |
1999-03-01 | 211.762683 | 138.19 |
1999-04-01 | 219.793874 | 138.32 |
1999-05-01 | 132.121255 | 125.36 |
# FUNGSI NRECA
def _STORAT(STORAGE, NOMINAL):
return STORAGE/NOMINAL
def _STORAGE(STORAGE, DELSTOR):
return STORAGE + DELSTOR
def _PRERAT(PRECIP, PET):
return PRECIP/PET
def _ETRAT(STORAT, PRERAT):
VAL = (STORAT / 2) + ((1 - (STORAT / 2)) * PRERAT)
return VAL if VAL < 1 else 1
def _AET(ETRAT, PET, CF=1):
return ETRAT * PET * CF
def _WATBAL(PRECIP, AET):
return PRECIP - AET
def _EXMRAT(WATBAL, STORAT):
if WATBAL < 0:
return 0
if STORAT > 1:
return 1 - (0.5 * (2 - STORAT)**2)
else:
return 0.5 * (STORAT**2)
def _EXMST(EXMRAT, WATBAL):
return EXMRAT * WATBAL
def _DELSTOR(WATBAL, EXMST):
return WATBAL - EXMST
def _GWRECH(PSUB, EXMST):
return PSUB * EXMST
def _GWSTOR2(GWSTOR_1, GWRECH):
return GWSTOR_1 + GWRECH
def _GWSTOR1(GWSTOR_2, GWFLOW):
return GWSTOR_2 - GWFLOW
def _GWFLOW(GWRAT, GWSTOR_2):
return GWRAT * GWSTOR_2
def _DFLOW(EXMST, GWRECH):
return EXMST - GWRECH
def _FLOW(GWFLOW, DFLOW):
return GWFLOW + DFLOW
def _DISCHARGE(FLOW, AREA, DAYS):
return (FLOW/1000) * AREA / (DAYS * 24 * 60 * 60)
def _NOMINAL(C, PRECIP_MEAN_ANNUAL):
return 100 + C*PRECIP_MEAN_ANNUAL
def model_NRECA(df, precip_col, pet_col,
MSTOR, GSTOR, PSUB, GWF, CF, C, AREA,
as_df=True, report='discharge'):
# sub_df
data = df.loc[:, [precip_col, pet_col]]
# info df
nrows = data.shape[0]
# initialization
storage, storat, prerat = (np.zeros(nrows) for _ in range(3))
etrat, aet, watbal = (np.zeros(nrows) for _ in range(3))
exmst, exmrat, delstor, gwrech = (np.zeros(nrows) for _ in range(4))
gwstor1, gwstor2, gwflow = (np.zeros(nrows) for _ in range(3))
dflow, flow, discharge = (np.zeros(nrows) for _ in range(3))
# calculation
precip_mean_annual = (data[precip_col].groupby(by=data.index.year)
.sum()
.mean())
nominal = _NOMINAL(C, precip_mean_annual)
days = data.index.days_in_month
precip = data.iloc[:, 0].values
pet = data.iloc[:, 1].values
for i in range(nrows):
if i != 0:
storage[i] = _STORAGE(storage[i - 1], delstor[i - 1])
else:
storage[i] = _STORAGE(MSTOR, 0)
storat[i] = _STORAT(storage[i], nominal)
prerat[i] = _PRERAT(precip[i], pet[i])
etrat[i] = _ETRAT(storat[i], prerat[i])
aet[i] = _AET(etrat[i], pet[i], CF=CF)
watbal[i] = _WATBAL(precip[i], aet[i])
exmrat[i] = _EXMRAT(watbal[i], storat[i])
exmst[i] = _EXMST(exmrat[i], watbal[i])
delstor[i] = _DELSTOR(watbal[i], exmst[i])
gwrech[i] = _GWRECH(PSUB, exmst[i])
if i != 0:
gwstor1[i] = _GWSTOR1(gwstor2[i - 1], gwflow[i - 1])
else:
gwstor1[i] = _GWSTOR1(GSTOR, 0)
gwstor2[i] = _GWSTOR2(gwstor1[i], gwrech[i])
gwflow[i] = _GWFLOW(GWF, gwstor2[i])
dflow[i] = _DFLOW(exmst[i], gwrech[i])
flow[i] = _FLOW(gwflow[i], dflow[i])
discharge[i] = _DISCHARGE(flow[i], AREA, days[i])
# results
if report.lower() == 'full':
results = np.stack((
days, precip, pet, storage, storat, prerat, etrat, aet, watbal,
exmrat, delstor, gwrech, gwstor1, gwstor2, gwflow, dflow, flow,
discharge
), axis=1)
columns_name = [
'DAYS', 'PRECIP', 'PET', 'STORAGE', 'STORAT', 'PRERAT', 'ETRAT',
'AET', 'WATBAL', 'EXMRAT', 'DELSTOR', 'GWRECH', 'GWSTOR1',
'GWSTOR2', 'GWFLOW', 'DFLOW', 'FLOW', 'DISCHARGE'
]
elif report.lower() == 'partial':
results = np.stack((
precip, pet, storage, gwstor2, flow, discharge), axis=1)
columns_name = ['PRECIP', 'PET', 'STORAGE',
'GWSTOR2', 'FLOW', 'DISCHARGE']
elif report.lower() == 'flow':
results = flow
columns_name = ['FLOW']
elif report.lower() == 'discharge':
results = discharge
columns_name = ['DISCHARGE']
else:
raise ValueError(
str(report) + ' not identified. ' +
'Use full / partial / flow / discharge.'
)
if as_df:
return pd.DataFrame(
data=results, index=data.index, columns=columns_name
)
else:
return results
model_NRECA()
¶Fungsi model_NRECA()
membangkitkan nilai debit menggunakan model yang dikembangkan oleh Crawford & Thurin (1981). Untuk mengetahui masing-masing batasan parameter yang digunakan dalam model ini, harap melihat pada makalah terkait. Dalam notebook ini akan lebih fokus dalam penggunaan fungsi.
Fungsi memiliki 10 argumen yang harus diisi, dan 2 argumen opsional. Dibagi menjadi 3 bagian argumen yaitu:
df
: dataset dalam bentuk pandas.DataFrame
.precip_col
: nama kolom untuk curah hujan (presipitasi).pet_col
: nama kolom untuk evapotranspirasi potensial.MSTOR
: Moist Storage ($mm$).GSTOR
: Groundwater Storage ($mm$).PSUB
: Percent Discharge to Subbase. ($0.3-0.8$)GWF
: Groundwater Fraction. ($0.2-0.9$)CF
: Crop Factor.C
: Koefisien karakteristik DAS. ($\{0.2,0.25\}$)AREA
: Luas DAS. ($m^2$)as_df
: True
(default), keluaran berupa pandas.DataFrame
. Keluaran berupa numpy.ndarray
jika False
.report
: discharge
(default). Terdapat beberapa nilai untuk report
:full
: keluaran akan menyertakan seluruh peubah yang dihitung dalam model seperti WATBAL
, ETRAT
, dll.partial
: keluaran hanya menyertakan kolom PRECIP
, PET
, STORAGE
, GWSTOR2
, FLOW
, dan DISCHARGE
.flow
: keluaran hanya menyertakan kolom FLOW
.discharge
: keluaran hanya menyertakan kolom DISCHARGE
.Berikut informasi satuan keluaran model:
PRECIP
, PET
, STORAGE
, AET
, WATBAL
, EXMST
, DELSTOR
, GWRECH
, GWSTOR1
, GWSTOR2
, GWFLOW
, DFLOW
, FLOW
bersatuan $mm/bulan$.STORAGE
, STORAT
, PRERAT
, ETRAT
, EXMRAT
merupakan rasio.as_df=True
, report='discharge'
(default)¶Keluaran berupa pandas.DataFrame
.
model_NRECA(df=dataset, precip_col='PRECIP', pet_col='PET',
MSTOR=1000, GSTOR=100, PSUB=0.4, GWF=0.2, CF=0.6, C=0.25,
AREA=1450.6e6).head()
DISCHARGE | |
---|---|
1999-01-01 | 130.474329 |
1999-02-01 | 124.872153 |
1999-03-01 | 66.438083 |
1999-04-01 | 70.787505 |
1999-05-01 | 41.998656 |
Sebagian argumen bisa disimpan dalam bentuk dictionary.
parameter = {
'MSTOR': 1000,
'GSTOR': 100,
'PSUB': 0.4,
'GWF': 0.2,
'CF': 0.6,
'C': 0.25,
'AREA': 1450.6e6
}
model_NRECA(df=dataset, precip_col='PRECIP', pet_col='PET', **parameter).head()
DISCHARGE | |
---|---|
1999-01-01 | 130.474329 |
1999-02-01 | 124.872153 |
1999-03-01 | 66.438083 |
1999-04-01 | 70.787505 |
1999-05-01 | 41.998656 |
as_df=False
¶Keluaran berupa numpy.array
model_NRECA(df=dataset, precip_col='PRECIP', pet_col='PET', **parameter,
as_df=False)
array([130.47432853, 124.87215318, 66.43808331, 70.78750455, 41.99865614, 28.45640939, 24.28274088, 19.7163445 , 12.90107429, 44.81356496, 53.93315909, 90.06341788, 95.31267461, 145.23779804, 117.57476262, 78.75685187, 91.24644372, 57.71641073, 37.27870151, 22.07183021, 24.97489079, 22.58939002, 35.79315297, 26.23750108, 138.0338475 , 235.17042934, 116.11823068, 87.70033456, 104.97934937, 109.51542424, 116.47710193, 34.51540744, 71.42737681, 33.81058852, 34.87405532, 36.58594762, 110.40885462, 112.85425887, 66.53789581, 97.93354457, 35.04151568, 22.39154574, 47.55438496, 16.65090467, 13.76474786, 10.65657899, 119.65349912, 57.66475933, 137.39388673, 199.12246399, 159.01388739, 172.29628498, 133.74151565, 44.77420904, 34.66390377, 27.73112302, 36.38175359, 120.31395471, 77.61700498, 126.38336198, 114.21457768, 123.66605913, 174.79252429, 195.67904295, 99.25002079, 61.22475824, 71.16772612, 31.42972227, 57.45303916, 44.65959205, 50.75185708, 91.32761479, 141.91241982, 103.15499701, 96.64076894, 63.34595895, 91.80413108, 133.77947175, 81.64626021, 31.2233851 , 56.84185425, 53.0918533 , 111.44781119, 124.64417381, 133.68443001, 82.44572213, 103.95773659, 78.36551337, 70.09140577, 28.12253459, 21.77228484, 17.41782788, 14.39873771, 11.14740984, 26.2002207 , 71.83216778, 77.60770993, 94.77473292, 112.44031817, 75.38748734, 52.32065854, 67.14854195, 35.46514257, 18.59214774, 15.3695088 , 21.20863532, 10.74189578, 49.26924462, 91.25050052, 170.94322505, 77.98879732, 81.89122336, 59.30935231, 27.5965865 , 19.29362528, 28.03710337, 13.98514355, 65.98636155, 69.40576331, 54.21634462])
report
¶report='full'
¶model_NRECA(df=dataset, precip_col='PRECIP', pet_col='PET', **parameter,
report='full').head()
DAYS | PRECIP | PET | STORAGE | STORAT | PRERAT | ETRAT | AET | WATBAL | EXMRAT | DELSTOR | GWRECH | GWSTOR1 | GWSTOR2 | GWFLOW | DFLOW | FLOW | DISCHARGE | |
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
1999-01-01 | 31.0 | 507.000000 | 142.63 | 1000.000000 | 1.323067 | 3.554652 | 1.0 | 85.578 | 421.422000 | 0.770881 | 96.555980 | 129.946408 | 100.000000 | 229.946408 | 45.989282 | 194.919612 | 240.908894 | 130.474329 |
1999-02-01 | 28.0 | 374.228918 | 128.84 | 1096.555980 | 1.450817 | 2.904602 | 1.0 | 77.304 | 296.924918 | 0.849199 | 44.776648 | 100.859308 | 183.957127 | 284.816435 | 56.963287 | 151.288962 | 208.252249 | 124.872153 |
1999-03-01 | 31.0 | 211.762683 | 138.19 | 1141.332627 | 1.510059 | 1.532402 | 1.0 | 82.914 | 128.848683 | 0.879979 | 15.464559 | 45.353650 | 227.853148 | 273.206798 | 54.641360 | 68.030474 | 122.671834 | 66.438083 |
1999-04-01 | 30.0 | 219.793874 | 138.32 | 1156.797186 | 1.530520 | 1.589025 | 1.0 | 82.992 | 136.801874 | 0.889794 | 15.076374 | 48.690200 | 218.565438 | 267.255638 | 53.451128 | 73.035300 | 126.486428 | 70.787505 |
1999-05-01 | 31.0 | 132.121255 | 125.36 | 1171.873560 | 1.550467 | 1.053935 | 1.0 | 75.216 | 56.905255 | 0.898960 | 5.749712 | 20.462217 | 213.804510 | 234.266727 | 46.853345 | 30.693325 | 77.546671 | 41.998656 |
report='partial'
¶model_NRECA(df=dataset, precip_col='PRECIP', pet_col='PET', **parameter,
report='partial').head()
PRECIP | PET | STORAGE | GWSTOR2 | FLOW | DISCHARGE | |
---|---|---|---|---|---|---|
1999-01-01 | 507.000000 | 142.63 | 1000.000000 | 229.946408 | 240.908894 | 130.474329 |
1999-02-01 | 374.228918 | 128.84 | 1096.555980 | 284.816435 | 208.252249 | 124.872153 |
1999-03-01 | 211.762683 | 138.19 | 1141.332627 | 273.206798 | 122.671834 | 66.438083 |
1999-04-01 | 219.793874 | 138.32 | 1156.797186 | 267.255638 | 126.486428 | 70.787505 |
1999-05-01 | 132.121255 | 125.36 | 1171.873560 | 234.266727 | 77.546671 | 41.998656 |
report='flow'
¶model_NRECA(df=dataset, precip_col='PRECIP', pet_col='PET', **parameter,
report='flow').head()
FLOW | |
---|---|
1999-01-01 | 240.908894 |
1999-02-01 | 208.252249 |
1999-03-01 | 122.671834 |
1999-04-01 | 126.486428 |
1999-05-01 | 77.546671 |
report='discharge'
(default)¶model_NRECA(df=dataset, precip_col='PRECIP', pet_col='PET', **parameter,
report='discharge').head()
DISCHARGE | |
---|---|
1999-01-01 | 130.474329 |
1999-02-01 | 124.872153 |
1999-03-01 | 66.438083 |
1999-04-01 | 70.787505 |
1999-05-01 | 41.998656 |
- 20191214 - 1.0.0 - Initial
Source code in this notebook is licensed under a MIT License. Data in this notebook is licensed under a Creative Common Attribution 4.0 International.