-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathp.C1.river_direction.py
155 lines (148 loc) · 8.22 KB
/
p.C1.river_direction.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
#!/usr/bin/env python
#
############################################################################
#
# MODULE : p.C1.river_direction.py
# AUTHOR(S) : Sanzana P. 01/12/2014
#
# PURPOSE : To update the direction of each segment river (upstream or downstream)
#
# COPYRIGHT : IRSTEA-UC-UCH
# This file is part of GeoPUMMA
#
# This program is free software; you can redistribute it and/or
# modify it under the terms of the GNU General Public License
# as published by the Free Software Foundation; either version 3
# of the License, or (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program; if not, see <http://www.gnu.org/licenses/>.
#
#
#############################################################################
#
#
#
import grass.script as grass
env = grass.gisenv()
print env
vectors = grass.read_command("g.list", type='vect')
print vectors
ditch_t=raw_input("Please enter the name of the ditch : ")
outlet_t=raw_input("Please enter the name of the oulet : ")
d_flow=raw_input("Please enter direction flow (-1 Upstream or 1 Downstream) : ")
columns_ditch=grass.read_command("v.info",map=ditch_t,flags='c')
print columns_ditch
river_id=raw_input("Please enter column name of the river id : ")
columns_output=grass.read_command("v.info",map=outlet_t,flags='c')
print columns_output
outlet_id=raw_input("Please enter column name of the output id : ")
outlet_count = grass.read_command("v.db.select",map=outlet_t,col=outlet_id,flags='c')
for i in outlet_count.splitlines():
river_where=river_id+"="+str(i)
grass.run_command("v.extract",input=ditch_t,output='ditch',where=river_where,overwrite=True)
ditch="ditch"
outlet_where=outlet_id+"="+str(i)
grass.run_command("v.extract",input=outlet_t,output='outlet',where=outlet_where,overwrite=True)
outlet="outlet"
#add start and end coordinates in ditch table
grass.run_command("v.db.addcol",map=ditch,col='start_E double,start_N double,end_E double,end_N double')
grass.run_command("v.db.addcol",map=ditch,col='start_E_0 double,start_N_0 double,end_E_0 double,end_N_0 double')
grass.run_command("v.db.addcol",map=ditch,col='start_E_1 double,start_N_1 double,end_E_1 double,end_N_1 double')
grass.run_command("v.db.addcol",map=ditch,col='Up_Str double')
grass.run_command("v.db.addcol",map=ditch,col='Do_Str double')
grass.run_command("v.to.db",map=ditch,option='start',columns='start_E_0,start_N_0')
grass.run_command("v.to.db",map=ditch,option='end',columns='end_E_0,end_N_0')
#the next stpes are necessary to resolve floating numeric problems
grass.run_command("v.db.update",map=ditch,col='start_E_1', value='start_E_0/10000')
grass.run_command("v.db.update",map=ditch,col='start_E', value='start_E_1*10000')
grass.run_command("v.db.update",map=ditch,col='start_N_1', value='start_N_0/10000')
grass.run_command("v.db.update",map=ditch,col='start_N', value='start_N_1*10000')
grass.run_command("v.db.update",map=ditch,col='end_E_1', value='end_E_0/10000')
grass.run_command("v.db.update",map=ditch,col='end_E', value='end_E_1*10000')
grass.run_command("v.db.update",map=ditch,col='end_N_1', value='end_N_0/10000')
grass.run_command("v.db.update",map=ditch,col='end_N', value='end_N_1*10000')
#extract nodes
grass.run_command("v.to.points",input=ditch,output='temp_point_1',flags='nt',overwrite=True)
#clean duplicated nodes
points = grass.read_command("v.out.ascii",input='temp_point_1',fs=',')
a=open("points_ascii.txt","w")
a.write(points)
a.close()
grass.run_command("v.in.ascii",input='points_ascii.txt',output='temp_point_2',fs=',',overwrite=True, columns='E double, N double,index int')
grass.run_command("v.clean",input='temp_point_2',output='temp_point_3',tool='rmdupl',overwrite=True)
#get only nodes cleaned
points_cleaned = grass.read_command("v.out.ascii",input='temp_point_3',fs=',')
b=open("points_cleaned.txt","w")
b.write(points_cleaned)
b.close()
grass.run_command("v.in.ascii",input='points_cleaned.txt',output='temp_point_4',fs=',',overwrite=True,columns='E_0 double, N_0 double,index int')
grass.run_command("v.db.addcol",map='temp_point_4',col='E_1 double,N_1 double,E double,N double')
#the next stpes are necessary to resolve floating numeric problems
grass.run_command("v.db.update",map='temp_point_4',col='N_1', value='N_0/10000')
grass.run_command("v.db.update",map='temp_point_4',col='N', value='N_1*10000')
grass.run_command("v.db.update",map='temp_point_4',col='E_1', value='E_0/10000')
grass.run_command("v.db.update",map='temp_point_4',col='E', value='E_1*10000')
grass.run_command("v.db.update",map='temp_point_4',col='index')
#get category number of outlet point
copy_outlet = outlet + ",outlet_temp"
grass.run_command("g.copy",vect=copy_outlet,overwrite=True)
grass.run_command("v.db.droptable",flags='f',map='outlet_temp')
grass.run_command("v.db.addtable",map='outlet_temp')
grass.run_command("v.db.addcol",map='outlet_temp',col='cat_out int')
grass.run_command("v.what.vect",vector='outlet_temp',column='cat_out',qvector='temp_point_4',qcolum='cat',dmax='10')
out_cat_1 = grass.read_command("v.db.select",map='outlet_temp',flags='c',col='cat_out')
print "out_cat_1 " + out_cat_1
#delete space with out information
out_cat_2 = out_cat_1.rsplit()
out_cat_3 = out_cat_2[0]
#add nodes to ditch network
grass.run_command("v.net",input=ditch,points='temp_point_4',output='ditch_points',operation='connect',thresh='1.0',overwrite=True)
#list_index
list_cat = grass.read_command("v.db.select",map='temp_point_4',col='cat',layer='1',flags='c')
#calculate distances to outlet from start node and end node each line
for i in list_cat.splitlines():
c=open("goto.txt","w")
text="1 " + str(i)+" "+str(out_cat_3)
print "texto de nodo a nodo "
print text
c.write(text)
c.close()
grass.run_command("v.net.path",input='ditch_points',out='path_temp',file='goto.txt',overwrite=True)
grass.run_command("v.db.addcol",map='path_temp',col='length double')
grass.run_command("v.db.update",map='path_temp',col='length',value='0')
grass.run_command("v.to.db",map='path_temp',col='length',option='length')
length_temp_1 = grass.read_command("v.db.select",map='path_temp',flags='c',col='length')
length_temp_2 = length_temp_1.rsplit()
length_temp_3 = length_temp_2[0]
print "Distancia " + str(length_temp_3)
E1=grass.read_command("v.db.select",map='temp_point_4',flags='c',col='E',where="cat=%s"%i)
E2=E1.rsplit()
E3=E2[0]
N1=grass.read_command("v.db.select",map='temp_point_4',flags='c',col='N',where="cat=%s"%i)
N2=N1.rsplit()
N3=N2[0]
condition1="start_E=%s"%E3+" AND "+"start_N=%s"%N3
print condition1
grass.run_command("v.db.update",map=ditch,col='Up_Str',value=length_temp_3,where=condition1)
condition2="end_E=%s"%E3+" AND "+"end_N=%s"%N3
print condition2
grass.run_command("v.db.update",map=ditch,col='Do_Str',value=length_temp_3,where=condition2)
grass.run_command("v.db.select",map=ditch,col='cat,Up_Str,Do_Str')
#check for correct direction
list_cat_ditch = grass.read_command("v.db.select",map=ditch,col='cat',layer='1',flags='c')
for i in list_cat_ditch.splitlines():
Up_Str1=grass.read_command("v.db.select",map=ditch,flags='c',col='Up_Str',where="cat=%s"%i)
Up_Str2=Up_Str1.rsplit()
Up_Str3=Up_Str2[0]
Do_Str1=grass.read_command("v.db.select",map=ditch,flags='c',col='Do_Str',where="cat=%s"%i)
Do_Str2=Do_Str1.rsplit()
Do_Str3=Do_Str2[0]
if (float(Up_Str3)-float(Do_Str3))*float(d_flow) < 0:
grass.run_command("v.edit",map=ditch_t,tool='flip',where="cat=%s"%i)
grass.run_command("g.remove",vect='ditch,outlet,ditch_points,outlet_temp,path_temp,temp_point_1,temp_point_2,temp_point_3,temp_point_4',overwrite=True)