In [3]:
from scipy import *
from cmath import *
from numpy import *
import numpy as np
import fileinput
import re
import pandas as pd
from scipy.interpolate import spline
from sklearn.metrics import *
import seaborn as sns
from scipy.interpolate import interp1d
from datetime import timedelta
import datetime
import calendar
import matplotlib as plt
import xml.etree.ElementTree as ET
from os.path import join
#get_ipython().magic(u'pylab inline')
from CasCade import *
%pylab inline
#import pylab as plt
In [4]:
from CasCade import *
import matplotlib.pyplot as plt
In [5]:
filepath= 'W:\\w2\\studenten_praktikanten\\Guarneri\\Model\\Bkup\\Para02_2016_09_21'#'W:\\w2\\studenten_praktikanten\\Guarneri\\Model\\Para02'#'W:\\w2\\studenten_praktikanten\\Guarneri\\Iff-Mainz'#
filename= 'para02.xml'#'Iff-Mainz_P2008-mitRet.xml'##'Iff-Mainz_P2008-mitRet.xml'#
In [6]:
pj = CasCade_Project()
pj.load_xml(filepath,filename)
In [7]:
margins_hight_left = pd.DataFrame(map(lambda i: pj.cross_sections_describe.profile[i].y[pj.cross_sections_describe.profile[i].x<300].max(), range(len(pj.cross_sections_describe))))
margins_hight_left['riverstation'] = pj.cross_sections_describe['riverstation']
margins_hight_right = pd.DataFrame(map(lambda i: pj.cross_sections_describe.profile[i].y[pj.cross_sections_describe.profile[i].x>300].max(), range(len(pj.cross_sections_describe))))
margins_hight_right['riverstation'] = pj.cross_sections_describe['riverstation']
In [8]:
scatter(margins_hight_left.riverstation.astype('float'),margins_hight_left[0],c='r')
scatter(margins_hight_right.riverstation.astype('float'),margins_hight_right[0],c='b')
Out[8]:
In [10]:
for i in range(len(pj.cross_sections_describe.profile)):
print(pj.cross_sections_describe.riverstation[i])
pj.cross_sections_describe.profile[i].plot('x','y',figsize=(10,1))
#scatter(margins_hight_left.x[i],margins_hight_left.y[i])
#scatter(margins_hight_right.x[i],margins_hight_right.y[i])
xlim(0,2500)
show()
In [9]:
pj.load_steady_runfile(filepath,'RunPara02.xml')
Discharge boarder conditions:
In [10]:
pj.runfile.border_condition_q.values
Out[10]:
Water Level boarder conditions:
In [46]:
pj.runfile.border_condition_y.values
Out[46]:
Exporting
In [12]:
#pj.export_runfile(filepath,'RunPara02.xml','Stationary')
Running the model
In [13]:
pj.run(filepath)
Open BoatSim and export .csv data
Read exported csv
In [11]:
expd = pd.read_table(join(filepath,'para02_out.csv'))
Steps
In [12]:
len(expd.Zeit.unique())
Out[12]:
In [13]:
expd.head(5)
Out[13]:
In [17]:
expd.tail(5)
Out[17]:
In [14]:
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
#plt.plot(expd.Ort[expd.Zeit == 0.000012],expd.Wasserspiegel[expd.Zeit == 0.000012])
ax1.set_ylim(74,92)
ax1.set_xlim(0,150)
xticks(range(1,150,11))
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
#y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
#ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
#ax1.scatter(x,y_max,c='orange')
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Q = 50 $m^3/s$',color='g')
ax2.plot(expd.Ort[expd.Zeit == 0.000012],expd.Geschwindigkeit[expd.Zeit == 0.000012],c='r',lw=1)
ax1.plot(expd.Ort[expd.Zeit == 0.005787],expd.Wasserspiegel[expd.Zeit == 0.005787])
plt.savefig(join(filepath,'Revision_Points','result_q50.png'),dpi=600)
The result above shows the result of a very low discharge scenario. With that, it is possible to assess the unconformities of the river slope. There are at least 7 velocity peaks, tha suggest that at this locations cross-sections are working as a weirs. We will have a closer look to understand what is happening in every location.
In [15]:
locations = pd.Series([2.5, 26.,64.,76.,96.5,121.,123.2])
In [16]:
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(15,5))
#missing legend
ax1.set_ylim(70,92)
ax1.set_xlim(0,150)
xticks(range(1,150,11))
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
ax1.scatter(x,y_min,c='g',s=5)
#ax1.scatter(x,y_mean,c='y')
ax1.scatter(x,y_max,c='orange',s=5)
for i in range(len(locations)):
ax1.text(locations[i],90,str(i+1),ha='center',size=12)
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Q = 50 $m^3/s$',color='g')
ax2.plot(expd.Ort[expd.Zeit == 0.000012],expd.Geschwindigkeit[expd.Zeit == 0.000012],c='r',lw=1)
ax1.plot(expd.Ort[expd.Zeit == 0.005787],expd.Wasserspiegel[expd.Zeit == 0.005787])
plt.savefig(join(filepath,'Revision_Points','result_q50_long.png'),dpi=600)
a closer look gives us:
In [17]:
for j in range(len(locations)):
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
ax1.set_ylim(70,92)
ax1.set_xlim(locations[j]-5,locations[j]+5)
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
ax1.scatter(x,y_max,c='orange')
ax1.set_xticks(range(1,150))
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Location = '+ str(locations[j])+u'\n' + r' Q = 50 $m^3/s$')
ax2.set_xlim(locations[j]-5,locations[j]+5)
ax2.plot(expd.Ort[expd.Zeit == 0.000012],expd.Geschwindigkeit[expd.Zeit == 0.000012],c='r',lw=1)
ax1.plot(expd.Ort[expd.Zeit == 0.005787],expd.Wasserspiegel[expd.Zeit == 0.005787])
plt.savefig(join(filepath,'Revision_Points','P'+str(j+1).zfill(2)+'_Initial_Profiles.png'),dpi=200)
plt.show()
In [18]:
j=0
After the alteration it does not show the problem anymore, need to fix it.
In [19]:
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
ax1.set_ylim(70,92)
ax1.set_xlim(locations[j]-5,locations[j]+5)
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
ax1.scatter(x,y_max,c='orange')
ax1.set_xticks(range(1,150))
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Location = '+ str(locations[j])+u'\n' + r' Q = 50 $m^3/s$')
ax2.set_xlim(locations[j]-5,locations[j]+5)
ax2.plot(expd.Ort[expd.Zeit == 0.000012],expd.Geschwindigkeit[expd.Zeit == 0.000012],c='r',lw=1)
ax1.plot(expd.Ort[expd.Zeit == 0.005787],expd.Wasserspiegel[expd.Zeit == 0.005787])
plt.show()
In this location we can see clearly that the increase in velocity and the 'waterstep' are caused by a sudden decrese of available depth.
In [20]:
index_cs = pj.cross_sections_describe.riverstation.astype(float).index[(pj.cross_sections_describe.riverstation.astype(float)/1000>locations[j]-1)&(pj.cross_sections_describe.riverstation.astype(float)/1000<locations[j]+1)]
After the alteration it does not show the problem anymore, need to fix it.
In [22]:
fig, ax = plt.subplots()
sns.set_style("darkgrid", {'axes.grid' : True})
for i in index_cs:
#print(pj.cross_sections_describe.riverstation[i]+' meters')
pj.cross_sections_describe.profile[i].plot('x','y',figsize=(6,3),ax=ax,label=pj.cross_sections_describe.riverstation[i]+' meters')
ax.set_yticks(range(78,88,1))
ylabel('Altitude [m]')
xlabel('Profile [m]')
#scatter(margins_hight_left.x[i],margins_hight_left.y[i])
#scatter(margins_hight_right.x[i],margins_hight_right.y[i])
ax.set_xlim(0,800)
plt.savefig(join(filepath,'Revision_Points','P'+str(j+1).zfill(2)+'_CrossSections_Initial.png'),dpi=200)
fig.show()
With this cross sections from location 1 ploted it is clear that there is a sudden rise of 3 meters between the cross-section 1800 and 2100, this variation is only normalized at the riverstation 3000 (there is also a significant fariation in area). Coinsidently, 3 meters is the firts 'level curve' from the nautical chart from were the depth informations for the cross sections were extracted. We need to have a closer look at the nautical chart.
In [15]:
pj.plot_position() #Not using this plot but the elements self.cross_sections_describe.position it generates.
In [50]:
from matplotlib.pyplot import *
%matplotlib inline
In [51]:
mappath = 'C:\\Users\\guarneri\\Documents\\ipnb\\figures'
In [28]:
index_cs
Out[28]:
In [29]:
easting, northing = (float(pj.cross_sections_describe.position[8][0]),float(pj.cross_sections_describe.position[8][1]))
easting2, northing2 = (float(pj.cross_sections_describe.position[8][2]),float(pj.cross_sections_describe.position[8][3]))
p1 = utmToLatLng(21, easting, northing, northernHemisphere=False)
p2 = utmToLatLng(21, easting2, northing2, northernHemisphere=False)
pm = (mean([p1[0],p2[0]]),mean([p1[1],p2[1]]))
right_margin = map(lambda x: utmToLatLng(21,pj.cross_sections_describe.position_list[0][1][x],pj.cross_sections_describe.position_list[1][1][x], northernHemisphere=False),range(len(pj.cross_sections_describe.position_list[0][0])))
left_margin = map(lambda x: utmToLatLng(21,pj.cross_sections_describe.position_list[0][0][x],pj.cross_sections_describe.position_list[1][0][x], northernHemisphere=False),range(len(pj.cross_sections_describe.position_list[0][0])))
In [30]:
#known projection error.
from urllib2 import urlopen, HTTPError
from pylab import imshow, imread, show
import time
pyplot.clf()
fig = plt.figure()
lon = [pm[0]-0.01,pm[0]+0.01]
lat = [pm[1]-0.02,pm[1]+0.02]
scale = 20000
print "Downloading map... "
tries = 0
url = None
while tries < 5:
tries += 1
print 'Try {}...'.format(tries)
try:
url = urlopen('http://parent.tile.openstreetmap.org/cgi-bin/export?'
'bbox={lat1:.2f},{lon1:.2f},{lat2:.2f},{lon2:.2f}&'
'scale={scale:d}&format=png'.format(lat1=lat[0],
lat2=lat[1],
lon1=lon[0],
lon2=lon[1],
scale=scale))
except HTTPError:
time.sleep(5)
continue
else:
print 'Map successfully downloaded.'
break
if url is None:
print 'Failed to download a map.'
else:
m = imread(url)
imshow(m, extent=lat+lon, aspect='equal')
grid(0)
plot((pd.DataFrame(right_margin)[1][index_cs],pd.DataFrame(left_margin)[1][index_cs]),(pd.DataFrame(right_margin)[0][index_cs],pd.DataFrame(left_margin)[0][index_cs]))
#plot((E1,E2),(N1,N2))
xlim(lat)
ylim(lon)
for i in index_cs:
scatter(pd.DataFrame(right_margin)[1][i],pd.DataFrame(right_margin)[0][i])
text(pd.DataFrame(right_margin)[1][i]+0.001,pd.DataFrame(right_margin)[0][i],pj.cross_sections_describe.riverstation[i])
#show()
savefig(join(filepath,'Revision_Points','P'+str(j+1).zfill(2)+'_MAP.png'),dpi=200,bbox_inches='tight')
In [31]:
import matplotlib.image as mpimg
img = mpimg.imread(join(filepath,'Revision_Points','Nautica_Chart_P01.PNG'))
figure(figsize=(11,5))
imgplot = plt.imshow(img)
grid(0)
plt.savefig(join(filepath,'Revision_Points','P'+str(j+1).zfill(2)+'_CHART.png'),dpi=200,bbox_iches='tight')
It is possible to see that in the river station 2400 the depths are not as low as the cross-section sugests. in fact in this regions they should look preaty much alike. We are going to replace the x,y values of the 2100,2400 and 2700 cross sections by the 3000 cross section. That would give an ok transition formm 1800 to 3300.
In [27]:
pj.cross_sections_describe.profile[6] = pj.cross_sections_describe.profile[9]
pj.cross_sections_describe.profile[7] = pj.cross_sections_describe.profile[9]
pj.cross_sections_describe.profile[8] = pj.cross_sections_describe.profile[9]
Now we export the cross-sections again, and run the model again to check if it corrected the issue.
In [ ]:
pj.cross_update_from_table()
In [31]:
pj.cross_update_from_table()
In [32]:
pj.cross_sections[6].export(join('W:\\w2\\studenten_praktikanten\\Guarneri\\Model\\Para02\\Profile','Rio_Paraguai_Lad_Pesper_7_2.xml'))
pj.cross_sections[7].export(join('W:\\w2\\studenten_praktikanten\\Guarneri\\Model\\Para02\\Profile','Rio_Paraguai_Lad_Pesper_8_2.xml'))
pj.cross_sections[8].export(join('W:\\w2\\studenten_praktikanten\\Guarneri\\Model\\Para02\\Profile','Rio_Paraguai_Lad_Pesper_9_2.xml'))
In [28]:
pj.run(filepath)
In [32]:
r1 = pd.read_table(join(filepath,'para02_out_revision_01.csv'))
In [33]:
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
#plt.plot(expd.Ort[expd.Zeit == 0.000012],expd.Wasserspiegel[expd.Zeit == 0.000012])
ax1.set_ylim(74,92)
ax1.set_xlim(0,150)
xticks(range(1,150,11))
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
#y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
#ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
#ax1.scatter(x,y_max,c='orange')
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Q = 50 $m^3/s$',color='g')
ax2.plot(r1.Ort[r1.Zeit == 0.000012],r1.Geschwindigkeit[r1.Zeit == 0.000012],c='r',lw=1)
ax1.plot(r1.Ort[r1.Zeit == 0.005787],r1.Wasserspiegel[r1.Zeit == 0.005787])
plt.savefig(join(filepath,'Revision_Points','result_q50_r1.png'),dpi=600)
In [34]:
j=0
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
ax1.set_ylim(70,92)
ax1.set_xlim(locations[j]-5,locations[j]+5)
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
ax1.scatter(x,y_max,c='orange')
ax1.set_xticks(range(1,150))
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Location = '+ str(locations[j])+u'\n' + r' Q = 50 $m^3/s$')
ax2.set_xlim(locations[j]-5,locations[j]+5)
ax2.plot(expd.Ort[expd.Zeit == 0.000012],expd.Geschwindigkeit[expd.Zeit == 0.000012],'--r',lw=1,label='V old')
ax1.plot(expd.Ort[expd.Zeit == 0.005787],expd.Wasserspiegel[expd.Zeit == 0.005787],'--b',label='WL old')
ax2.plot(r1.Ort[r1.Zeit == 0.000012],r1.Geschwindigkeit[r1.Zeit == 0.000012],'r',lw=1,label='V new')
ax1.plot(r1.Ort[r1.Zeit == 0.005787],r1.Wasserspiegel[r1.Zeit == 0.005787],'b',label='WL new')
ax1.legend(loc=2)
ax2.legend(loc=1)
plt.savefig(join(filepath,'Revision_Points','P01_Comparison_Profiles.png'),dpi=600)
plt.show()
We can see that the velocity peak and the water step have vanished. We can close this matter and follow to the next.
In [35]:
j=1
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
ax1.set_ylim(70,92)
ax1.set_xlim(locations[j]-5,locations[j]+5)
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
ax1.scatter(x,y_max,c='orange')
ax1.set_xticks(range(1,150))
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Location = '+ str(locations[j])+u'\n' + r' Q = 50 $m^3/s$')
ax2.set_xlim(locations[j]-5,locations[j]+5)
ax2.plot(expd.Ort[expd.Zeit == 0.000012],expd.Geschwindigkeit[expd.Zeit == 0.000012],c='r',lw=1,label='V')
ax1.plot(expd.Ort[expd.Zeit == 0.005787],expd.Wasserspiegel[expd.Zeit == 0.005787],'b',label='WL')
#ax2.plot(r1.Ort[r1.Zeit == 0.000012],r1.Geschwindigkeit[r1.Zeit == 0.000012],'--r',lw=1,label='V new')
#ax1.plot(r1.Ort[r1.Zeit == 0.005787],r1.Wasserspiegel[r1.Zeit == 0.005787],'--b',label='WL new')
ax1.legend(loc=2)
ax2.legend(loc=1)
#plt.savefig(join(filepath,'Revision_Points','P01_Results.png'),dpi=200)
#plt.show()
Out[35]:
In this location, there seems to exist two problems. A very low depth cross-section (in green) causing a waterstep and a peak in velocity, and a very low margin (orange) that is not causing problemns now for a dischange of 50 m3 but it could cause later for larger discharges.
In [52]:
sns.palplot(sns.cubehelix_palette(12, start=.5, rot=-.75))
In [37]:
index_cs
Out[37]:
In [38]:
index_cs = pj.cross_sections_describe.riverstation.astype(float).index[(pj.cross_sections_describe.riverstation.astype(float)/1000>locations[j]-1)&(pj.cross_sections_describe.riverstation.astype(float)/1000<locations[j]+1)]
fig, ax = plt.subplots()
sns.set_style("darkgrid", {'axes.grid' : True})
k=0
for i in index_cs:
#print(pj.cross_sections_describe.riverstation[i]+' meters')
pj.cross_sections_describe.profile[i].plot('x','y',figsize=(6,3),c=sns.cubehelix_palette(11, start=.5, rot=-.75)[k],ax=ax,label=pj.cross_sections_describe.riverstation[i]+' meters')
ax.set_yticks(range(78,88,1))
ylabel('Altitude [m]')
xlabel('Profile [m]')
legend(loc=4)
#scatter(margins_hight_left.x[i],margins_hight_left.y[i])
#scatter(margins_hight_right.x[i],margins_hight_right.y[i])
ax.set_xlim(0,1000)
k+=1
fig.show()
It seems that the problem is related with the presence of an island. In these cases, an island is best represented by a river split. We shall performe the split, run the model and check the results. If the problem persists we may need to have a closer look at the cross-sections depth.
In [44]:
import matplotlib.image as mpimg
img = mpimg.imread(join(filepath,'Revision_Points','P02_MAP.PNG'))
figure(figsize=(20,6))
imgplot = plt.imshow(img)
grid(0)
#plt.savefig(join(filepath,'Revision_Points','P'+str(j+1).zfill(2)+'_MAP.png'),dpi=200,bbox_iches='tight')
In [53]:
import utm
In [109]:
easting, northing = (float(pj.cross_sections_describe.position[index_cs[3]][0]),float(pj.cross_sections_describe.position[index_cs[3]][1]))
easting2, northing2 = (float(pj.cross_sections_describe.position[index_cs[3]][2]),float(pj.cross_sections_describe.position[index_cs[3]][3]))
p1 = utm.to_latlon(easting, northing, 21, 'K')
p2 = utm.to_latlon(easting2, northing2, 21, 'K')
pm = (mean([p1[1],p2[1]]),mean([p1[0],p2[0]]))
In [110]:
mappath = 'figures'
mapname = 'ldpe_r2.html'
central_latitute = pm[1]
central_longitude = pm[0]
In [111]:
pj.plot_folium(central_latitute,central_longitude,mappath,mapname)
In [45]:
%%html
<iframe width='100%' height='500' src='figures\\ldpe_r2.html' ></iframe>
To perform this split, we need to:
Creat one node at the begging of the split and one at the end. This would increase the number of branches from 1 to 4. We need to carefully create the new braches and declare them.
Split the following cross sections from the island in left and right: 25844, 26100, 26400
Relocate cross-sections to each different branches.
In [47]:
pyplot.clf()
Let's look the process for 1 cross section and than do faster for all of them.
In [48]:
cs_25844_1 = Cross_Section()
In [49]:
index_cs[2]
Out[49]:
In [50]:
pj.cross_sections[index_cs[2]].riverstation
Out[50]:
In [51]:
fig = pj.cross_sections[index_cs[2]].plot()
xticks(range(0,800,50))
grid(1)
In [42]:
pj.load_cross_sections()
In [43]:
cs_25844_1 = pj.cross_sections[index_cs[2]]
f = interp1d(cs_25844_1.profile.x.values,cs_25844_1.profile.y.values)
xnew = linspace(cs_25844_1.profile.x.min(),cs_25844_1.profile.x.max(),1000)
plt.plot(cs_25844_1.profile.x.values,cs_25844_1.profile.y.values,'o',xnew,f(xnew),'-')
new_profile = pd.DataFrame([pd.Series(xnew),pd.Series(f(xnew))]).T
new_profile.columns = ['x','y']
cs_25844_1.profile = new_profile
cs_25844_1.split(400)
cs_25844_1.profile = cs_25844_1.right_profile
cs_25844_1.export(join('W:\\w2\\studenten_praktikanten\\Guarneri\\Model\\Para02\\Profile','Rio_Paraguai_Lad_Pesper_'+str(index_cs[2]+1)+'_r'))
cs_25844_1.profile = cs_25844_1.left_profile
cs_25844_1.export(join('W:\\w2\\studenten_praktikanten\\Guarneri\\Model\\Para02\\Profile','Rio_Paraguai_Lad_Pesper_'+str(index_cs[2]+1)+'_l'))
In [44]:
cs_25844_1.right_profile.plot('x','y',figsize=(2,2))
xlim(0,400)
show()
cs_25844_1.left_profile.plot('x','y',figsize=(2,2))
xlim(400,750)
show()
In [45]:
index_cs[3]
Out[45]:
In [46]:
pj.cross_sections[index_cs[3]].riverstation
Out[46]:
In [47]:
fc3 = pj.cross_sections[index_cs[3]].plot()
In [48]:
cs_26100_1 = pj.cross_sections[index_cs[3]]
f = interp1d(cs_26100_1.profile.x.values,cs_26100_1.profile.y.values)
xnew = linspace(cs_26100_1.profile.x.min(),cs_26100_1.profile.x.max(),500)
plt.plot(cs_26100_1.profile.x.values,cs_26100_1.profile.y.values,'o',xnew,f(xnew),'-')
new_profile = pd.DataFrame([pd.Series(xnew),pd.Series(f(xnew))]).T
new_profile.columns = ['x','y']
cs_26100_1.profile = new_profile
In [49]:
cs_26100_1.split(350)
cs_26100_1.profile = cs_26100_1.right_profile
cs_26100_1.export(join('W:\\w2\\studenten_praktikanten\\Guarneri\\Model\\Para02\\Profile','Rio_Paraguai_Lad_Pesper_'+str(index_cs[3]+1)+'_r'))
cs_26100_1.profile = cs_26100_1.left_profile
cs_26100_1.export(join('W:\\w2\\studenten_praktikanten\\Guarneri\\Model\\Para02\\Profile','Rio_Paraguai_Lad_Pesper_'+str(index_cs[3]+1)+'_l'))
In [50]:
cs_26100_1.right_profile.plot('x','y',figsize=(2,2))
xlim(0,300)
show()
cs_26100_1.left_profile.plot('x','y',figsize=(2,2))
xlim(300,750)
show()
26400
In [51]:
print index_cs[4], pj.cross_sections[index_cs[4]].riverstation
In [52]:
fc3 = pj.cross_sections[index_cs[4]].plot()
In [53]:
cs_26400_1 = pj.cross_sections[index_cs[4]]
f = interp1d(cs_26400_1.profile.x.values,cs_26400_1.profile.y.values)
xnew = linspace(cs_26400_1.profile.x.min(),cs_26400_1.profile.x.max(),500)
plt.plot(cs_26400_1.profile.x.values,cs_26400_1.profile.y.values,'o',xnew,f(xnew),'-')
xticks(range(0,600,25),rotation=90)
new_profile = pd.DataFrame([pd.Series(xnew),pd.Series(f(xnew))]).T
new_profile.columns = ['x','y']
cs_26400_1.profile = new_profile
In [54]:
cs_26400_1.split(270)
cs_26400_1.profile = cs_26400_1.right_profile
cs_26400_1.export(join('W:\\w2\\studenten_praktikanten\\Guarneri\\Model\\Para02\\Profile','Rio_Paraguai_Lad_Pesper_'+str(index_cs[4]+1)+'_r'))
cs_26400_1.profile = cs_26400_1.left_profile
cs_26400_1.export(join('W:\\w2\\studenten_praktikanten\\Guarneri\\Model\\Para02\\Profile','Rio_Paraguai_Lad_Pesper_'+str(index_cs[4]+1)+'_l'))
In [55]:
cs_26400_1.right_profile.plot('x','y',figsize=(2,2))
xlim(0,275)
show()
cs_26400_1.left_profile.plot('x','y',figsize=(2,2))
xlim(275,750)
show()
CS spliting STATUS: Done. Let's correct the model file.
Altering Model file
First we need to create de new links
Than declare the new nodes
Than change the information on the branches table.
This is our current model structure.
In [56]:
pj.model.links
Out[56]:
We will now create the new links structure.
In [57]:
# [Node,To,Branch]
link1 = ['Lad_Pesper_Node_1','LE25200','Lad_LE25200']
link2 = ['LE25200','LE26700','LE25200_LE26700_1']
link3 = ['LE25200','LE26700','LE25200_LE26700_2']
link4 = ['LE26700','Lad_Pesper_Node_2','LE26700_Pesper']
In [58]:
links = [link1,link2,link3,link4]
In [59]:
pd.DataFrame(links,columns=pj.model.links.columns)
Out[59]:
In [60]:
new_links = pd.DataFrame(links,columns=pj.model.links.columns)
In [61]:
pj.model.links = new_links
Now we need to declare the new nodes
In [62]:
pj.nodes #old
Out[62]:
In [63]:
pd.Series(['Lad_Pesper_Node_1','LE25200','LE26700','Lad_Pesper_Node_2']) #new
Out[63]:
In [64]:
pj.nodes = pd.Series(['Lad_Pesper_Node_1','LE25200','LE26700','Lad_Pesper_Node_2'])
And now update the branches table, than update the model file
In [65]:
new_links.Branch[0]
Out[65]:
Renaming the branches at the extremes (beggining and end) of the stretch.
In [66]:
pj.branches['branch_names'][pj.branches.index<index_cs[2]] = new_links.Branch[0]
pj.branches['branch_names'][pj.branches.index>index_cs[4]] = new_links.Branch[3]
Renaming the branch name and the reach names for the firt split branch.
In [67]:
pj.branches['branch_names'][(pj.branches.index>=index_cs[2])&(pj.branches.index<=index_cs[4])] = new_links.Branch[1]
In [68]:
pj.branches['reach_names'][(pj.branches.index>=index_cs[2])&(pj.branches.index<=index_cs[4])] = pj.branches['reach_names'][(pj.branches.index>=index_cs[2])&(pj.branches.index<=index_cs[4])]+'_1'
Declaring the names of the new cross-sections files
In [69]:
pj.branches['crosssectionref_sources'][82] = 'Profile/Rio_Paraguai_Lad_Pesper_83_l.xml'
pj.branches['crosssectionref_sources'][83] = 'Profile/Rio_Paraguai_Lad_Pesper_84_l.xml'
pj.branches['crosssectionref_sources'][84] = 'Profile/Rio_Paraguai_Lad_Pesper_85_l.xml'
In [70]:
pj.branches['crosssectionref_names'][82] = 'Rio_Paraguai_Lad_Pesper_83_l'
pj.branches['crosssectionref_names'][83] = 'Rio_Paraguai_Lad_Pesper_84_l'
pj.branches['crosssectionref_names'][84] = 'Rio_Paraguai_Lad_Pesper_85_l'
In [71]:
pj.branches['crosssectionref_names'][84]
Out[71]:
Make a new branch from a copy of the first split
In [72]:
new_branches = pj.branches[(pj.branches.index>=index_cs[2])&(pj.branches.index<=index_cs[4])]
Declare its name.
In [73]:
new_branches['branch_names'] = new_links.Branch[2]
In [74]:
new_branches['branch_names']
Out[74]:
Declare it's reaches names.
In [77]:
new_branches['reach_names'][82] = new_branches['reach_names'][82][:-1]+'2'
new_branches['reach_names'][83] = new_branches['reach_names'][83][:-1]+'2'
new_branches['reach_names'][84] = new_branches['reach_names'][84][:-1]+'2'
In [78]:
new_branches['reach_names'][84]
Out[78]:
In [79]:
new_branches['crosssectionref_sources'][82]
Out[79]:
Declare its new cross_sections file names and location:
In [80]:
new_branches['crosssectionref_sources'][82] = 'Profile/Rio_Paraguai_Lad_Pesper_83_r.xml'
new_branches['crosssectionref_sources'][83] = 'Profile/Rio_Paraguai_Lad_Pesper_84_r.xml'
new_branches['crosssectionref_sources'][84] = 'Profile/Rio_Paraguai_Lad_Pesper_85_r.xml'
In [81]:
new_branches['crosssectionref_names'][82] = 'Rio_Paraguai_Lad_Pesper_83_r'
new_branches['crosssectionref_names'][83] = 'Rio_Paraguai_Lad_Pesper_84_r'
new_branches['crosssectionref_names'][84] = 'Rio_Paraguai_Lad_Pesper_85_r'
In [82]:
new_branches['crosssectionref_names'][84]
Out[82]:
In [83]:
new_branches['crosssectionref_sources'][82]
Out[83]:
In [84]:
pj.branches['branch_names'][(pj.branches.index>=index_cs[2])&(pj.branches.index<=index_cs[4])]
Out[84]:
In [85]:
new_branches['branch_names']
Out[85]:
In [86]:
a = pj.branches
b = new_branches
In [87]:
frames = [a,b]
result = pd.concat(frames)
In [88]:
result.index = range(len(result))
In [89]:
result
Out[89]:
In [90]:
pj.branches = result
In [91]:
pj.export_model(join(filepath,'para02_t'))
In [108]:
pj.run(filepath)
In [105]:
pj.branches.branch_names.unique()
Out[105]:
Partial Resuls
In [110]:
r2.Zeit
Out[110]:
In [62]:
filepath2 = 'W:\\w2\\studenten_praktikanten\\Guarneri\\Model\\Bkup\\Para02_2016_09_22'
In [63]:
filepath2 = 'W:\\w2\\studenten_praktikanten\\Guarneri\\Model\\Bkup\\Para02_2016_09_22'
r2 = pd.read_table(join(filepath2,'para02_out_revision_02.csv'))
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
#plt.plot(expd.Ort[expd.Zeit == 0.000012],expd.Wasserspiegel[expd.Zeit == 0.000012])
ax1.set_ylim(74,92)
ax1.set_xlim(0,150)
xticks(range(1,150,11))
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
#y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
#ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
#ax1.scatter(x,y_max,c='orange')
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Q = 100 $m^3/s$',color='g')
ax2.plot(r2.Ort[r2.Zeit == 0.000012],r2.Geschwindigkeit[r2.Zeit == 0.000012],c='r',lw=1)
ax1.plot(r2.Ort[r2.Zeit == 0.005787],r2.Wasserspiegel[r2.Zeit == 0.005787])
#plt.savefig(join(filepath,'Revision_Points','result_q50_r1.png'),dpi=600)
Out[63]:
The problem is worst. There is something wrong with the spliting concept. Probably related to removin 1 reach. I should ask Rolf:
Test: Remove the cross sections 83, 84 and 85
In [65]:
filepath2 = 'W:\\w2\\studenten_praktikanten\\Guarneri\\Model\\Bkup\\Para02_2016_09_22'
r2 = pd.read_table(join(filepath2,'para02_out_revision_02.csv'))
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
#plt.plot(expd.Ort[expd.Zeit == 0.000012],expd.Wasserspiegel[expd.Zeit == 0.000012])
ax1.set_ylim(74,92)
ax1.set_xlim(0,150)
xticks(range(1,150,11))
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
#y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
#ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
#ax1.scatter(x,y_max,c='orange')
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Q = 100 $m^3/s$',color='g')
ax2.plot(r2.Ort[r2.Zeit == 0.000012],r2.Geschwindigkeit[r2.Zeit == 0.000012],c='r',lw=1)
ax1.plot(r2.Ort[r2.Zeit == 0.005787],r2.Wasserspiegel[r2.Zeit == 0.005787])
#plt.savefig(join(filepath,'Revision_Points','result_q50_r1.png'),dpi=600)
Out[65]:
In [70]:
j=1
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
ax1.set_ylim(70,92)
ax1.set_xlim(locations[j]-5,locations[j]+5)
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
ax1.scatter(x,y_max,c='orange')
ax1.set_xticks(range(1,150))
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Location = '+ str(locations[j])+u'\n' + r' Q = 100 $m^3/s$')
ax2.set_xlim(locations[j]-5,locations[j]+5)
ax2.plot(expd.Ort[expd.Zeit == 0.000012],expd.Geschwindigkeit[expd.Zeit == 0.000012],'--r',lw=1,label='V old')
ax1.plot(expd.Ort[expd.Zeit == 0.005787],expd.Wasserspiegel[expd.Zeit == 0.005787],'--b',label='WL old')
ax2.plot(r2.Ort[r2.Zeit == 0.000012],r2.Geschwindigkeit[r2.Zeit == 0.000012],'r',lw=1,label='V new')
ax1.plot(r2.Ort[r2.Zeit == 0.005787],r2.Wasserspiegel[r2.Zeit == 0.005787],'b',label='WL new')
ax1.legend(loc=2)
ax2.legend(loc=1)
#plt.savefig(join(filepath,'Revision_Points','P01_Comparison_Profiles.png'),dpi=600)
plt.show()
Removing 3 cross-sections solved the problem. It solved the problem. For now lets keep it withoud the 3 cross-sections
In [48]:
j=2
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
ax1.set_ylim(70,92)
ax1.set_xlim(locations[j]-5,locations[j]+5)
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
ax1.scatter(x,y_max,c='orange')
ax1.set_xticks(range(1,150))
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Location = '+ str(locations[j])+u'\n' + r' Q = 50 $m^3/s$')
ax2.set_xlim(locations[j]-5,locations[j]+5)
ax2.plot(expd.Ort[expd.Zeit == 0.000012],expd.Geschwindigkeit[expd.Zeit == 0.000012],c='r',lw=1,label='V')
ax1.plot(expd.Ort[expd.Zeit == 0.005787],expd.Wasserspiegel[expd.Zeit == 0.005787],'b',label='WL')
#ax2.plot(r1.Ort[r1.Zeit == 0.000012],r1.Geschwindigkeit[r1.Zeit == 0.000012],'--r',lw=1,label='V new')
#ax1.plot(r1.Ort[r1.Zeit == 0.005787],r1.Wasserspiegel[r1.Zeit == 0.005787],'--b',label='WL new')
ax1.legend(loc=2)
ax2.legend(loc=1)
#plt.savefig(join(filepath,'Revision_Points','P01_Results.png'),dpi=200)
#plt.show()
Out[48]:
In [51]:
sns.palplot(sns.cubehelix_palette(6, start=.5, rot=-.75))
In [64]:
index_cs = pj.cross_sections_describe.riverstation.astype(float).index[(pj.cross_sections_describe.riverstation.astype(float)/1000>locations[j]-1)&(pj.cross_sections_describe.riverstation.astype(float)/1000<locations[j]+1)]
fig, ax = plt.subplots()
sns.set_style("darkgrid", {'axes.grid' : True})
k=0
for i in index_cs:
#print(pj.cross_sections_describe.riverstation[i]+' meters')
pj.cross_sections_describe.profile[i].plot('x','y',figsize=(6,3),c=sns.cubehelix_palette(6, start=.5, rot=-.75)[k],ax=ax,label=pj.cross_sections_describe.riverstation[i]+' meters'+' - '+str(i))
ax.set_yticks(range(78,88,1))
ylabel('Altitude [m]')
xlabel('Profile [m]')
legend(loc=4)
#scatter(margins_hight_left.x[i],margins_hight_left.y[i])
#scatter(margins_hight_right.x[i],margins_hight_right.y[i])
ax.set_xlim(0,1000)
k+=1
fig.show()
In [63]:
index_cs
Out[63]:
In [57]:
import utm
In [58]:
import utm
easting, northing = (float(pj.cross_sections_describe.position[index_cs[3]][0]),float(pj.cross_sections_describe.position[index_cs[3]][1]))
easting2, northing2 = (float(pj.cross_sections_describe.position[index_cs[3]][2]),float(pj.cross_sections_describe.position[index_cs[3]][3]))
p1 = utm.to_latlon(easting, northing, 21, 'K')
p2 = utm.to_latlon(easting2, northing2, 21, 'K')
pm = (mean([p1[1],p2[1]]),mean([p1[0],p2[0]]))
mappath = 'figures'
mapname = 'ldpe_r3.html'
central_latitute = pm[1]
central_longitude = pm[0]
pj.plot_folium(central_latitute,central_longitude,mappath,mapname)
In [60]:
%%html
<iframe width='100%' height='500' src='figures\\ldpe_r3.html' ></iframe>
In [62]:
import matplotlib.image as mpimg
img = mpimg.imread(join(filepath,'Revision_Points','P03_MAP.PNG'))
figure(figsize=(20,10))
imgplot = plt.imshow(img)
grid(0)
#plt.savefig(join(filepath,'Revision_Points','P'+str(j+1).zfill(2)+'_MAP.png'),dpi=200,bbox_iches='tight')
It seems that the cross-sections are representing poorly the river depths in this sections. CS 64200 and 62500 seem to be the problem. Let's remove them and run the model to check the results.
In [91]:
pj.run(filepath)
In [92]:
r3 = pd.read_table(join(filepath,'para02_out_revision_03.csv'))
In [93]:
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
#plt.plot(expd.Ort[expd.Zeit == 0.000012],expd.Wasserspiegel[expd.Zeit == 0.000012])
ax1.set_ylim(74,92)
ax1.set_xlim(0,150)
xticks(range(1,150,11))
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
#y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
#ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
#ax1.scatter(x,y_max,c='orange')
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Q = 50 $m^3/s$',color='g')
ax2.plot(r1.Ort[r1.Zeit == 0.000012],r1.Geschwindigkeit[r1.Zeit == 0.000012],c='r',lw=1)
ax1.plot(r1.Ort[r1.Zeit == 0.005787],r1.Wasserspiegel[r1.Zeit == 0.005787])
#plt.savefig(join(filepath,'Revision_Points','result_q50_r1.png'),dpi=600)
Out[93]:
In [97]:
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
#plt.plot(expd.Ort[expd.Zeit == 0.000012],expd.Wasserspiegel[expd.Zeit == 0.000012])
ax1.set_ylim(74,92)
ax1.set_xlim(0,150)
xticks(range(1,150,11))
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
#y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
#ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
#ax1.scatter(x,y_max,c='orange')
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Q = 50 $m^3/s$',color='g')
ax2.plot(r3.Ort[r3.Zeit == 0.000012],r3.Geschwindigkeit[r3.Zeit == 0.000012],c='r',lw=1)
ax1.plot(r3.Ort[r3.Zeit == 0.005787],r3.Wasserspiegel[r3.Zeit == 0.005787])
plt.savefig(join(filepath,'Revision_Points','result_q50_r3.png'),dpi=600)
In [96]:
j=2
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
ax1.set_ylim(70,92)
ax1.set_xlim(locations[j]-5,locations[j]+5)
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
ax1.scatter(x,y_max,c='orange')
ax1.set_xticks(range(1,150))
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Location = '+ str(locations[j])+u'\n' + r' Q = 100 $m^3/s$')
ax2.set_xlim(locations[j]-5,locations[j]+5)
ax2.plot(expd.Ort[expd.Zeit == 0.000012],expd.Geschwindigkeit[expd.Zeit == 0.000012],'--r',lw=1,label='V old')
ax1.plot(expd.Ort[expd.Zeit == 0.005787],expd.Wasserspiegel[expd.Zeit == 0.005787],'--b',label='WL old')
ax2.plot(r3.Ort[r3.Zeit == 0.000012],r3.Geschwindigkeit[r3.Zeit == 0.000012],'r',lw=1,label='V new')
ax1.plot(r3.Ort[r3.Zeit == 0.005787],r3.Wasserspiegel[r3.Zeit == 0.005787],'b',label='WL new')
ax1.legend(loc=2)
ax2.legend(loc=1)
plt.savefig(join(filepath,'Revision_Points','P03_Comparison_Profiles.png'),dpi=600)
plt.show()
In [71]:
j=3
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
ax1.set_ylim(70,92)
ax1.set_xlim(locations[j]-5,locations[j]+5)
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
ax1.scatter(x,y_max,c='orange')
ax1.set_xticks(range(1,150))
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Location = '+ str(locations[j])+u'\n' + r' Q = 50 $m^3/s$')
ax2.set_xlim(locations[j]-5,locations[j]+5)
ax2.plot(expd.Ort[expd.Zeit == 0.000012],expd.Geschwindigkeit[expd.Zeit == 0.000012],c='r',lw=1,label='V')
ax1.plot(expd.Ort[expd.Zeit == 0.005787],expd.Wasserspiegel[expd.Zeit == 0.005787],'b',label='WL')
#ax2.plot(r1.Ort[r1.Zeit == 0.000012],r1.Geschwindigkeit[r1.Zeit == 0.000012],'--r',lw=1,label='V new')
#ax1.plot(r1.Ort[r1.Zeit == 0.005787],r1.Wasserspiegel[r1.Zeit == 0.005787],'--b',label='WL new')
ax1.legend(loc=2)
ax2.legend(loc=1)
plt.savefig(join(filepath,'Revision_Points','P04_Results.png'),dpi=200)
#plt.show()
Apperently, one cross-section at the beginning of km 76 is causing the water step and consequente increase in velocity.
In [56]:
index_cs = pj.cross_sections_describe.riverstation.astype(float).index[(pj.cross_sections_describe.riverstation.astype(float)/1000>locations[j]-1)&(pj.cross_sections_describe.riverstation.astype(float)/1000<locations[j]+1)]
fig, ax = plt.subplots()
sns.set_style("darkgrid", {'axes.grid' : True})
k=0
for i in index_cs:
#print(pj.cross_sections_describe.riverstation[i]+' meters')
pj.cross_sections_describe.profile[i].plot('x','y',figsize=(6,3),c=sns.cubehelix_palette(6, start=.5, rot=-.75)[k],ax=ax,label=pj.cross_sections_describe.riverstation[i]+' meters'+' - '+str(i))
ax.set_yticks(range(78,88,1))
ylabel('Altitude [m]')
xlabel('Profile [m]')
legend(loc=4)
#scatter(margins_hight_left.x[i],margins_hight_left.y[i])
#scatter(margins_hight_right.x[i],margins_hight_right.y[i])
ax.set_xlim(0,1700)
k+=1
plt.savefig(join(filepath,'Revision_Points','P04_CS.png'),dpi=600)
fig.show()
In [21]:
%%html
<iframe width='100%' height='500' src='figures\\ldpe_r3.html' ></iframe>
Alterations performed straight in the source code. And now we run the model to check the results.
In [57]:
pj.run(filepath)
In [58]:
r4 = pd.read_table(join(filepath,'para02_out_revision_04.csv'))
In [59]:
r4 = pd.read_table(join(filepath,'para02_out_revision_04.csv'))
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
#plt.plot(expd.Ort[expd.Zeit == 0.000012],expd.Wasserspiegel[expd.Zeit == 0.000012])
ax1.set_ylim(74,92)
ax1.set_xlim(0,150)
xticks(range(1,150,11))
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
#y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
#ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
#ax1.scatter(x,y_max,c='orange')
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Q = 50 $m^3/s$',color='g')
ax2.plot(r4.Ort[r4.Zeit == 0.000012],r4.Geschwindigkeit[r4.Zeit == 0.000012],c='r',lw=1)
ax1.plot(r4.Ort[r4.Zeit == 0.005787],r4.Wasserspiegel[r4.Zeit == 0.005787])
plt.savefig(join(filepath,'Revision_Points','result_q50_r4.png'),dpi=600)
Its noticeable that the change had an impact on the water level along the whole river upstream section. As a future methodologie remark, I conclude that this alterations should be made from downstream to upstream.
In [60]:
j=3
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
ax1.set_ylim(70,92)
ax1.set_xlim(locations[j]-5,locations[j]+5)
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
ax1.scatter(x,y_max,c='orange')
ax1.set_xticks(range(1,150))
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Location = '+ str(locations[j])+u'\n' + r' Q = 100 $m^3/s$')
ax2.set_xlim(locations[j]-5,locations[j]+5)
ax2.plot(expd.Ort[expd.Zeit == 0.000012],expd.Geschwindigkeit[expd.Zeit == 0.000012],'--r',lw=1,label='V old')
ax1.plot(expd.Ort[expd.Zeit == 0.005787],expd.Wasserspiegel[expd.Zeit == 0.005787],'--b',label='WL old')
ax2.plot(r4.Ort[r4.Zeit == 0.000012],r4.Geschwindigkeit[r4.Zeit == 0.000012],'r',lw=1,label='V new')
ax1.plot(r4.Ort[r4.Zeit == 0.005787],r4.Wasserspiegel[r4.Zeit == 0.005787],'b',label='WL new')
ax1.legend(loc=2)
ax2.legend(loc=1)
plt.savefig(join(filepath,'Revision_Points','P04_Comparison_Profiles.png'),dpi=600)
plt.show()
Notice the expresive decrease in waterlevel and velocity.
In [98]:
j=4
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
ax1.set_ylim(70,92)
ax1.set_xlim(locations[j]-5,locations[j]+5)
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
ax1.scatter(x,y_max,c='orange')
ax1.set_xticks(range(1,150))
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Location = '+ str(locations[j])+u'\n' + r' Q = 50 $m^3/s$')
ax2.set_xlim(locations[j]-5,locations[j]+5)
ax2.plot(expd.Ort[expd.Zeit == 0.000012],expd.Geschwindigkeit[expd.Zeit == 0.000012],c='r',lw=1,label='V')
ax1.plot(expd.Ort[expd.Zeit == 0.005787],expd.Wasserspiegel[expd.Zeit == 0.005787],'b',label='WL')
#ax2.plot(r1.Ort[r1.Zeit == 0.000012],r1.Geschwindigkeit[r1.Zeit == 0.000012],'--r',lw=1,label='V new')
#ax1.plot(r1.Ort[r1.Zeit == 0.005787],r1.Wasserspiegel[r1.Zeit == 0.005787],'--b',label='WL new')
ax1.legend(loc=2)
ax2.legend(loc=1)
plt.savefig(join(filepath,'Revision_Points','P05_Results.png'),dpi=200)
#plt.show()
Apperently, one cross-section at the beginning of km 96.5 is causing the water step and consequente increase in velocity.
In [99]:
index_cs = pj.cross_sections_describe.riverstation.astype(float).index[(pj.cross_sections_describe.riverstation.astype(float)/1000>locations[j]-1)&(pj.cross_sections_describe.riverstation.astype(float)/1000<locations[j]+1)]
fig, ax = plt.subplots()
sns.set_style("darkgrid", {'axes.grid' : True})
k=0
for i in index_cs:
#print(pj.cross_sections_describe.riverstation[i]+' meters')
pj.cross_sections_describe.profile[i].plot('x','y',figsize=(6,3),c=sns.cubehelix_palette(6, start=.5, rot=-.75)[k],ax=ax,label=pj.cross_sections_describe.riverstation[i]+' meters'+' - '+str(i))
ax.set_yticks(range(74,88,1))
ylabel('Altitude [m]')
xlabel('Profile [m]')
legend(loc=4)
#scatter(margins_hight_left.x[i],margins_hight_left.y[i])
#scatter(margins_hight_right.x[i],margins_hight_right.y[i])
ax.set_xlim(0,3700)
k+=1
plt.savefig(join(filepath,'Revision_Points','P05_CS.png'),dpi=600)
fig.show()
This is a big island, and should be considered to be splited. For now lets remove the cross-section 96603 and see the results.
In [100]:
pj.run(filepath)
In [103]:
r5 = pd.read_table(join(filepath,'para02_out_revision_05.csv'))
In [109]:
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
#plt.plot(expd.Ort[expd.Zeit == 0.000012],expd.Wasserspiegel[expd.Zeit == 0.000012])
ax1.set_ylim(74,92)
ax1.set_xlim(0,150)
xticks(range(1,150,11))
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
#y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
#ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
#ax1.scatter(x,y_max,c='orange')
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Q = 50 $m^3/s$',color='g')
ax2.plot(r5.Ort[r5.Zeit == 0.000012],r5.Geschwindigkeit[r5.Zeit == 0.000012],c='r',lw=1)
ax1.plot(r5.Ort[r5.Zeit == 0.005787],r5.Wasserspiegel[r5.Zeit == 0.005787])
#plt.savefig(join(filepath,'Revision_Points','result_q50_r5.png'),dpi=600)
Out[109]:
In [105]:
j=4
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
ax1.set_ylim(70,92)
ax1.set_xlim(locations[j]-5,locations[j]+5)
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
ax1.scatter(x,y_max,c='orange')
ax1.set_xticks(range(1,150))
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Location = '+ str(locations[j])+u'\n' + r' Q = 100 $m^3/s$')
ax2.set_xlim(locations[j]-5,locations[j]+5)
ax2.plot(expd.Ort[expd.Zeit == 0.000012],expd.Geschwindigkeit[expd.Zeit == 0.000012],'--r',lw=1,label='V old')
ax1.plot(expd.Ort[expd.Zeit == 0.005787],expd.Wasserspiegel[expd.Zeit == 0.005787],'--b',label='WL old')
ax2.plot(r5.Ort[r5.Zeit == 0.000012],r5.Geschwindigkeit[r5.Zeit == 0.000012],'r',lw=1,label='V new')
ax1.plot(r5.Ort[r5.Zeit == 0.005787],r5.Wasserspiegel[r5.Zeit == 0.005787],'b',label='WL new')
ax1.legend(loc=2)
ax2.legend(loc=1)
plt.savefig(join(filepath,'Revision_Points','P05_Comparison_Profiles.png'),dpi=600)
plt.show()
The velocity peak droped significantly (~33%). There is still a water step. It can be solved by removing another cross_section. Lets remove the cross-section on km 96 and see what happens.
In [106]:
pj.run(filepath)
In [110]:
r5_2 = pd.read_table(join(filepath,'para02_out_revision_05_2.csv'))
In [111]:
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
#plt.plot(expd.Ort[expd.Zeit == 0.000012],expd.Wasserspiegel[expd.Zeit == 0.000012])
ax1.set_ylim(74,92)
ax1.set_xlim(0,150)
xticks(range(1,150,11))
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
#y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
#ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
#ax1.scatter(x,y_max,c='orange')
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Q = 50 $m^3/s$',color='g')
ax2.plot(r5_2.Ort[r5_2.Zeit == 0.000012],r5_2.Geschwindigkeit[r5_2.Zeit == 0.000012],c='r',lw=1)
ax1.plot(r5_2.Ort[r5_2.Zeit == 0.005787],r5_2.Wasserspiegel[r5_2.Zeit == 0.005787])
plt.savefig(join(filepath,'Revision_Points','result_q50_r5_2.png'),dpi=600)
As we seen before, the change introduces new water steps towards upstream. Nontheless we were able by now to eliminate most of the peak velocities that were way over resonable.
In [112]:
j=4
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
ax1.set_ylim(70,92)
ax1.set_xlim(locations[j]-5,locations[j]+5)
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
ax1.scatter(x,y_max,c='orange')
ax1.set_xticks(range(1,150))
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Location = '+ str(locations[j])+u'\n' + r' Q = 100 $m^3/s$')
ax2.set_xlim(locations[j]-5,locations[j]+5)
ax2.plot(expd.Ort[expd.Zeit == 0.000012],expd.Geschwindigkeit[expd.Zeit == 0.000012],'--r',lw=1,label='V old')
ax1.plot(expd.Ort[expd.Zeit == 0.005787],expd.Wasserspiegel[expd.Zeit == 0.005787],'--b',label='WL old')
ax2.plot(r5_2.Ort[r5_2.Zeit == 0.000012],r5_2.Geschwindigkeit[r5_2.Zeit == 0.000012],'r',lw=1,label='V new')
ax1.plot(r5_2.Ort[r5_2.Zeit == 0.005787],r5_2.Wasserspiegel[r5_2.Zeit == 0.005787],'b',label='WL new')
ax1.legend(loc=2)
ax2.legend(loc=1)
plt.savefig(join(filepath,'Revision_Points','P05_2_Comparison_Profiles.png'),dpi=600)
plt.show()
The removal of two cross-sections eliminated complete the water step and peak velocity in this region.
In [80]:
%%html
<iframe width='100%' height='500' src='figures\\ldpe_r3.html' ></iframe>
In [116]:
j=5
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
ax1.set_ylim(74,92)
ax1.set_xlim(locations[j]-5,locations[j]+5)
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
ax1.scatter(x,y_max,c='orange')
ax1.set_xticks(range(1,150))
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Location = '+ str(locations[j])+u'\n' + r' Q = 50 $m^3/s$')
ax2.set_xlim(locations[j]-5,locations[j]+5)
ax2.plot(expd.Ort[expd.Zeit == 0.000012],expd.Geschwindigkeit[expd.Zeit == 0.000012],c='r',lw=1,label='V')
ax1.plot(expd.Ort[expd.Zeit == 0.005787],expd.Wasserspiegel[expd.Zeit == 0.005787],'b',label='WL')
#ax2.plot(r1.Ort[r1.Zeit == 0.000012],r1.Geschwindigkeit[r1.Zeit == 0.000012],'--r',lw=1,label='V new')
#ax1.plot(r1.Ort[r1.Zeit == 0.005787],r1.Wasserspiegel[r1.Zeit == 0.005787],'--b',label='WL new')
ax1.legend(loc=2)
ax2.legend(loc=1)
plt.savefig(join(filepath,'Revision_Points','P06_Results.png'),dpi=200)
#plt.show()
There is a lot going on on this plot. The Location in question is the km 121, not the 123. It seems that with the removal of 1 or 2 cross-sections, we can solve this problem, but before, lets have a look at the cross sections.
In [118]:
index_cs
Out[118]:
In [121]:
index_cs = pj.cross_sections_describe.riverstation.astype(float).index[(pj.cross_sections_describe.riverstation.astype(float)/1000>locations[j]-1)&(pj.cross_sections_describe.riverstation.astype(float)/1000<locations[j]+1)]
fig, ax = plt.subplots()
sns.set_style("darkgrid", {'axes.grid' : True})
k=0
for i in index_cs:
#print(pj.cross_sections_describe.riverstation[i]+' meters')
pj.cross_sections_describe.profile[i].plot('x','y',figsize=(6,3),c=sns.cubehelix_palette(7, start=.5, rot=-.75)[k],ax=ax,label=pj.cross_sections_describe.riverstation[i]+' meters'+' - '+str(i))
ax.set_yticks(range(74,88,1))
ylabel('Altitude [m]')
xlabel('Profile [m]')
legend(loc=4)
#scatter(margins_hight_left.x[i],margins_hight_left.y[i])
#scatter(margins_hight_right.x[i],margins_hight_right.y[i])
ax.set_xlim(0,1900)
k+=1
plt.savefig(join(filepath,'Revision_Points','P06_CS.png'),dpi=600)
fig.show()
In [114]:
%%html
<iframe width='100%' height='500' src='figures\\ldpe_r3.html' ></iframe>
It seems that if we remove cross-section 121500 e 120900, we can solve the problem.
In [122]:
pj.run(filepath)
In [123]:
r6 = pd.read_table(join(filepath,'para02_out_revision_06.csv'))
In [132]:
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
#plt.plot(expd.Ort[expd.Zeit == 0.000012],expd.Wasserspiegel[expd.Zeit == 0.000012])
ax1.set_ylim(74,92)
ax1.set_xlim(0,150)
xticks(range(1,150,11))
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
#y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
#ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
#ax1.scatter(x,y_max,c='orange')
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Q = 50 $m^3/s$',color='g')
ax2.plot(r6.Ort[r6.Zeit == 0.000012],r6.Geschwindigkeit[r6.Zeit == 0.000012],c='r',lw=1,label='V new')
ax1.plot(r6.Ort[r6.Zeit == 0.005787],r6.Wasserspiegel[r6.Zeit == 0.005787],'b',label='WL new')
ax2.plot(expd.Ort[expd.Zeit == 0.000012],expd.Geschwindigkeit[expd.Zeit == 0.000012],'--y',lw=0.5,label='V old')
ax1.plot(expd.Ort[expd.Zeit == 0.005787],expd.Wasserspiegel[expd.Zeit == 0.005787],'--b',lw=0.5,label='WL old')
ax1.legend(loc=2)
ax2.legend(loc=1)
#plt.savefig(join(filepath,'Revision_Points','result_q50_r5.png'),dpi=600)
Out[132]:
In [131]:
j=5
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
ax1.set_ylim(70,92)
ax1.set_xlim(locations[j]-5,locations[j]+5)
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
ax1.scatter(x,y_max,c='orange')
ax1.set_xticks(range(1,150))
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Location = '+ str(locations[j])+u'\n' + r' Q = 100 $m^3/s$')
ax2.set_xlim(locations[j]-5,locations[j]+5)
ax2.plot(expd.Ort[expd.Zeit == 0.000012],expd.Geschwindigkeit[expd.Zeit == 0.000012],'--r',lw=1,label='V old')
ax1.plot(expd.Ort[expd.Zeit == 0.005787],expd.Wasserspiegel[expd.Zeit == 0.005787],'--b',label='WL old')
ax2.plot(r6.Ort[r6.Zeit == 0.000012],r6.Geschwindigkeit[r6.Zeit == 0.000012],'r',lw=1,label='V new')
ax1.plot(r6.Ort[r6.Zeit == 0.005787],r6.Wasserspiegel[r6.Zeit == 0.005787],'b',label='WL new')
ax1.legend(loc=2)
ax2.legend(loc=1)
plt.savefig(join(filepath,'Revision_Points','P06_Comparison_Profiles.png'),dpi=600)
plt.show()
The water step and the velocity peak seem to be solved.
In [138]:
j=6
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
ax1.set_ylim(70,90)
ax1.set_xlim(locations[j]-5,locations[j]+5)
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
ax1.scatter(x,y_max,c='orange')
ax1.set_xticks(range(1,150))
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Location = '+ str(locations[j])+u'\n' + r' Q = 50 $m^3/s$')
ax2.set_xlim(locations[j]-5,locations[j]+5)
ax2.plot(r6.Ort[r6.Zeit == 0.000012],r6.Geschwindigkeit[r6.Zeit == 0.000012],c='r',lw=1,label='V')
ax1.plot(r6.Ort[r6.Zeit == 0.005787],r6.Wasserspiegel[r6.Zeit == 0.005787],'b',label='WL')
#ax2.plot(r1.Ort[r1.Zeit == 0.000012],r1.Geschwindigkeit[r1.Zeit == 0.000012],'--r',lw=1,label='V new')
#ax1.plot(r1.Ort[r1.Zeit == 0.005787],r1.Wasserspiegel[r1.Zeit == 0.005787],'--b',label='WL new')
ax1.legend(loc=2)
ax2.legend(loc=1)
plt.savefig(join(filepath,'Revision_Points','P07_Results.png'),dpi=200)
#plt.show()
This looks terrible, we have a waterfall. A 5 meter fall in less than 1 km. Lets have a look at the cross-sections
In [134]:
index_cs = pj.cross_sections_describe.riverstation.astype(float).index[(pj.cross_sections_describe.riverstation.astype(float)/1000>locations[j]-1)&(pj.cross_sections_describe.riverstation.astype(float)/1000<locations[j]+1)]
fig, ax = plt.subplots()
sns.set_style("darkgrid", {'axes.grid' : True})
k=0
for i in index_cs:
#print(pj.cross_sections_describe.riverstation[i]+' meters')
pj.cross_sections_describe.profile[i].plot('x','y',figsize=(6,3),c=sns.cubehelix_palette(7, start=.5, rot=-.75)[k],ax=ax,label=pj.cross_sections_describe.riverstation[i]+' meters'+' - '+str(i))
ax.set_yticks(range(74,88,1))
ylabel('Altitude [m]')
xlabel('Profile [m]')
legend(loc=4)
#scatter(margins_hight_left.x[i],margins_hight_left.y[i])
#scatter(margins_hight_right.x[i],margins_hight_right.y[i])
ax.set_xlim(0,1900)
k+=1
plt.savefig(join(filepath,'Revision_Points','P07_CS.png'),dpi=600)
fig.show()
Vamos consultar as cartas náuticas
In [141]:
%%html
<iframe width='100%' height='500' src='figures\\ldpe_r3.html' ></iframe>
In [140]:
import matplotlib.image as mpimg
img = mpimg.imread(join(filepath,'Revision_Points','P07_CHART.PNG'))
figure(figsize=(20,10))
imgplot = plt.imshow(img)
grid(0)
#plt.savefig(join(filepath,'Revision_Points','P'+str(j+1).zfill(2)+'_MAP.png'),dpi=200,bbox_iches='tight')
In [142]:
pj.run(filepath)
In [143]:
r6 = pd.read_table(join(filepath,'para02_out_revision_07.csv'))
In [145]:
j=6
sns.set_style("darkgrid", {'axes.grid' : False})
fig, ax1 = plt.subplots(figsize=(5,5))
ax1.set_ylim(70,90)
ax1.set_xlim(locations[j]-5,locations[j]+5)
for i in range(len(pj.cross_sections_describe)):
x = float(pj.cross_sections_describe.riverstation[i])/1000
y_min = pj.cross_sections_describe.profile[i].y.min()
y_mean = pj.cross_sections_describe.profile[i].y.mean()
y_max = pj.cross_sections_describe.profile[i].y.max()
ax1.scatter(x,y_min,c='g')
#ax1.scatter(x,y_mean,c='y')
ax1.scatter(x,y_max,c='orange')
ax1.set_xticks(range(1,150))
ax1.grid('on')
ax2 = ax1.twinx()
ax1.set_ylabel('Water Level [m]', color='b')
ax2.set_ylabel('Velocity [m/s]', color='r')
ax1.set_xlabel('km')
plt.title(r'Location = '+ str(locations[j])+u'\n' + r' Q = 100 $m^3/s$')
ax2.set_xlim(locations[j]-5,locations[j]+5)
ax2.plot(r6.Ort[r6.Zeit == 0.000012],r6.Geschwindigkeit[r6.Zeit == 0.000012],c='r',lw=1,label='V')
ax1.plot(r6.Ort[r6.Zeit == 0.005787],r6.Wasserspiegel[r6.Zeit == 0.005787],'b',label='WL')
#ax2.plot(r1.Ort[r1.Zeit == 0.000012],r1.Geschwindigkeit[r1.Zeit == 0.000012],'--r',lw=1,label='V new')
#ax1.plot(r1.Ort[r1.Zeit == 0.005787],r1.Wasserspiegel[r1.Zeit == 0.005787],'--b',label='WL new')
ax1.legend(loc=2)
ax2.legend(loc=1)
plt.savefig(join(filepath,'Revision_Points','P07_Results.png'),dpi=200)
#plt.show()
In [ ]: