import numpy as np
import pandas as pd
import geopandas as gpd
from numba import jit
from numba.experimental import jitclass
from numba.typed import List
import xarray as xr
import source.riversection as sectorg
import source.s1driverflow as model
def makeSect(classsect, gdfsorted):
def from3Dto2D(pointz, porg=None):
point = pointz[:,:2]
if porg is None : porg = pointz[0]
v = point[-1] - point[0]
e = v/np.linalg.norm(v)
d = point - porg[:2]
L = np.dot(d, e)
Z = pointz[:,2]
return L, Z
channel = List()
for calc, dist in zip( gdfsorted['calc-input'].values, gdfsorted['distancefromDB'].values ):
typed_ps = List()
typed_ns = List()
for i, c in enumerate(calc['point-data']):
p3d = np.array( c['coordinates'] )
n = np.array( c['manning'] )
if len(n) == 1 : n = np.repeat(n, (len(p3d) - 1))
if i == 0 : porg = p3d[0]
L, Z = from3Dto2D(p3d, porg)
p = np.c_[L, Z]
typed_ps.append(p)
typed_ns.append(n)
channel.append( classsect.section(typed_ps, typed_ns, dist) )
return channel
fin = '8909090001CalcData.geojson'
gdfsect = gpd.read_file(fin)
# Since the downstream end is at point 3.000, the data downstream of that point will be deleted.
val = gdfsect[gdfsect['name'] == '3.000']['distancefromDB'].values[0]
gdfsect = gdfsect[gdfsect['distancefromDB'] >= val]
gdfsect = gdfsect.reset_index(drop=True)
gdfsectsorted = gdfsect.sort_values('distancefromDB', ascending=True)
chD2U = makeSect(sectorg, gdfsectsorted)
gdfsectsorted = gdfsect.sort_values('distancefromDB', ascending=False)
chU2D = makeSect(sectorg, gdfsectsorted)
dfhydro = pd.read_csv('hydrologicData.csv',index_col='date',parse_dates=True)
dt = float(10)
dfhydroip = dfhydro.resample(str(dt) + 'S').interpolate()
Qup = dfhydroip['代継橋_時間流量'].values
Hdown = dfhydroip['小島_時間水位'].values
%%time
# nonuniform flow
Qt = np.full(len(chD2U), Qup[0], dtype=np.float64)
Huni = model.NonUniformflow(chD2U, Qt, Hdown[0])
Auni = np.array( [chD2U[i].H2ABS(hh)[0] for i, hh in enumerate(Huni)] )
A0 = Auni[0]
# unsteady flow : 20hr
Q = np.full_like(Auni, Qup[0])
A, H = Auni[::-1], Huni[::-1]
for n in range(1, int(3600*20/dt)):
A, Q, H = model.UnSteadyflow(chU2D, A, Q, H, A0, Qup[0], dt)
Wall time: 9.41 s
%%time
nmax = len(chU2D)
ntmax = len(Qup)
Amat, Qmat, Hmat = np.zeros((ntmax, nmax)), np.zeros((ntmax, nmax)), np.zeros((ntmax, nmax))
Amat[0], Qmat[0], Hmat[0] = A, Q, H
for n in range(1, ntmax):
A0, _, _ = chU2D[-1].H2ABS(Hdown[n])
A, Q, H = model.UnSteadyflow(chU2D, A, Q, H, A0, Qup[n], dt)
Amat[n], Qmat[n], Hmat[n] = A, Q, H
Wall time: 14.4 s
gdfsectsorted = gdfsect.sort_values('distancefromDB', ascending=False)
ds = xr.Dataset( {'A': (['t','x'], Amat)
,'Q': (['t','x'], Qmat)
,'H': (['t','x'], Hmat)
,'zbmin': (['x'], [c.zbmin() for c in chU2D] )
}
, coords={'x':gdfsectsorted['distancefromDB'].values
, 't':dfhydroip.index.values}
,attrs = {'name':gdfsectsorted['name'].values.tolist()}
)
dout = ds.to_netcdf(r'calout.nc')
del dout