[GRASS-SVN] r59141 - in grass-addons/grass7/vector: . v.civil v.civil/v.civil.river v.civil/v.civil.road v.civil/v.civil.tools v.civil/v.civil.topo
svn_grass at osgeo.org
svn_grass at osgeo.org
Tue Feb 25 16:37:16 PST 2014
Author: jfc
Date: 2014-02-25 16:37:16 -0800 (Tue, 25 Feb 2014)
New Revision: 59141
Added:
grass-addons/grass7/vector/v.civil/
grass-addons/grass7/vector/v.civil/v.civil.river/
grass-addons/grass7/vector/v.civil/v.civil.river/Makefile
grass-addons/grass7/vector/v.civil/v.civil.river/v.civil.river.html
grass-addons/grass7/vector/v.civil/v.civil.river/v.civil.river.py
grass-addons/grass7/vector/v.civil/v.civil.road/
grass-addons/grass7/vector/v.civil/v.civil.road/Makefile
grass-addons/grass7/vector/v.civil/v.civil.road/v.civil.road.html
grass-addons/grass7/vector/v.civil/v.civil.road/v.civil.road.py
grass-addons/grass7/vector/v.civil/v.civil.tools/
grass-addons/grass7/vector/v.civil/v.civil.tools/Makefile
grass-addons/grass7/vector/v.civil/v.civil.tools/v.civil.tools.html
grass-addons/grass7/vector/v.civil/v.civil.tools/v.civil.tools.py
grass-addons/grass7/vector/v.civil/v.civil.topo/
grass-addons/grass7/vector/v.civil/v.civil.topo/Makefile
grass-addons/grass7/vector/v.civil/v.civil.topo/v.civil.topo.html
grass-addons/grass7/vector/v.civil/v.civil.topo/v.civil.topo.py
Log:
v.civil.*: modules for design of roads, channels, ports, etc
Added: grass-addons/grass7/vector/v.civil/v.civil.river/Makefile
===================================================================
--- grass-addons/grass7/vector/v.civil/v.civil.river/Makefile (rev 0)
+++ grass-addons/grass7/vector/v.civil/v.civil.river/Makefile 2014-02-26 00:37:16 UTC (rev 59141)
@@ -0,0 +1,9 @@
+# fix this relative to include/
+# or use absolute path to the GRASS source code
+MODULE_TOPDIR = ../..
+
+PGM = v.civil.river
+
+include $(MODULE_TOPDIR)/include/Make/Script.make
+
+default: script
Property changes on: grass-addons/grass7/vector/v.civil/v.civil.river/Makefile
___________________________________________________________________
Added: svn:mime-type
+ text/x-makefile
Added: svn:eol-style
+ native
Added: grass-addons/grass7/vector/v.civil/v.civil.river/v.civil.river.html
===================================================================
--- grass-addons/grass7/vector/v.civil/v.civil.river/v.civil.river.html (rev 0)
+++ grass-addons/grass7/vector/v.civil/v.civil.river/v.civil.river.html 2014-02-26 00:37:16 UTC (rev 59141)
@@ -0,0 +1,47 @@
+<h2>DESCRIPTION</h2>
+
+<em>v.civil.river</em> Export - Import geometry data to/from an 1-D hydrodynamic model.
+
+<h3>Export</h3>
+
+To export to an 1-D hydrodynamic model a cross section and cross-displ maps from the output of v.civil.road is required.
+And a raster dem too.
+<p>
+You must give a prefix name for output files.
+<p>
+There are two format files, RS,X,Y,Z format and RS,Station,Elevation:
+
+<ul>
+ <li><em>PrefixName</em>_format_xyz.csv</li>
+ <li><em>PrefixName</em>_format_stations.csv</li>
+</ul>
+<p>
+The schematic of the edge and the schematic of the river (lowest point of each cross-section)
+
+<ul>
+ <li><em>PrefixName</em>_schematic.txt</li>
+ <li><em>PrefixName</em>_schematic_real.txt</li>
+</ul>
+<p>
+Tree files to populate geometry tables
+<ul>
+ <li><em>PrefixName</em>_banks.txt</li>
+ <li><em>PrefixName</em>_length.txt</li>
+ <li><em>PrefixName</em>_manning.txt</li>
+</ul>
+
+<h3>Import</h3>
+<p>
+To import from 1-D hydrodynamic model a <em>.sdf</em> file (from export gis data) is required with only one profile and with
+water surface extents.
+<p>
+The output is a polygon with the water surface extents in each cross-section.
+<p>
+The <em>-d</em> flag create ASCII file describing the water surface profile along the channel axis for use within r.inund.fluv
+addon-module.
+
+<h2> AUTHOR</h2>
+<p>
+Jesus Fernandez-Capel Rosillo<br >
+Civil Engineer, Spain<br >
+jfc at alcd net<br >
\ No newline at end of file
Property changes on: grass-addons/grass7/vector/v.civil/v.civil.river/v.civil.river.html
___________________________________________________________________
Added: svn:mime-type
+ text/html
Added: svn:keywords
+ Author Date Id
Added: svn:eol-style
+ native
Added: grass-addons/grass7/vector/v.civil/v.civil.river/v.civil.river.py
===================================================================
--- grass-addons/grass7/vector/v.civil/v.civil.river/v.civil.river.py (rev 0)
+++ grass-addons/grass7/vector/v.civil/v.civil.river/v.civil.river.py 2014-02-26 00:37:16 UTC (rev 59141)
@@ -0,0 +1,447 @@
+#!/usr/bin/env python
+# -*- coding: utf-8
+############################################################################
+#
+# MODULE: v.civil.river, v0.5.0
+#
+# AUTHOR: Jesús Fernández-Capel Rosillo
+# Civil Engineer, Spain
+# jfc at alcd net
+#
+# PURPOSE: Export - Import geometry data to/from an 1-D hydrodynamic model.
+#
+# COPYRIGHT: (c) 2014 Jesús Fernández-Capel Rosillo.
+#
+# This program is free software under the GNU General Public
+# License (>=v2). Read the file COPYING that comes with GRASS
+# for details.
+#
+#############################################################################
+
+#%Module
+#% description: Export - Import geometry data to/from an 1-D hydrodynamic model.
+#% keywords: vector, Hec-Ras, hydraulic, river
+#%End
+
+#### Required section ####
+
+#%option
+#% key: ini
+#% type: string
+#% key_desc: Iniciar
+#% description: Title of project
+#% required: yes
+#%end
+
+
+#### Export section ####
+
+#%flag
+#% key: e
+#% description: Export to Hec-Ras
+#% guisection: Export
+#%end
+
+#%option
+#% key: mapcross
+#% type: string
+#% gisprompt: old,vector,vector
+#% description: Name for croos section map
+#% required: no
+#% guisection: Export
+#%end
+
+#%option
+#% key: dem
+#% type: string
+#% gisprompt: old,cell,raster
+#% key_desc: raster dem
+#% description: Name of DEM raster map
+#% guisection: Export
+#%end
+
+#%option
+#% key: mapcrossdispl
+#% type: string
+#% gisprompt: old,vector,vector
+#% description: Name for cross-displ map
+#% required: no
+#% guisection: Export
+#%end
+
+#%option
+#% key: manning
+#% type: string
+#% description: Values for Manning three columns ej: 0.050,0.045,0.050
+#% required: no
+#% answer: 0.050,0.040,0.050
+#% guisection: Export
+#%end
+
+#%option G_OPT_F_OUTPUT
+#% description: Prefix name for output geometry csv files
+#% required: no
+#% guisection: Export
+#%end
+
+#### Import section ####
+
+#%flag
+#% key: i
+#% description: Import from Hec-Ras
+#% guisection: Import
+#%end
+
+#%flag
+#% key: d
+#% description: Create ASCII file describing the water surface profile along the channel axis (r.inund.fluv)
+#% guisection: Import
+#%end
+
+#%option G_OPT_F_INPUT
+#% key: inputfile
+#% description: Name for input hec-ras file
+#% required: no
+#% guisection: Import
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: outwaterfile
+#% description: Name for output water surface polygon
+#% required: no
+#% guisection: Import
+#%end
+
+#%option G_OPT_F_OUTPUT
+#% key: outprofilefile
+#% description: Name for output water surface profile file (r.inund.fluv)
+#% required: no
+#% guisection: Import
+#%end
+
+
+
+
+
+import os, sys
+from math import *
+import grass.script as grass
+from grass.lib.gis import *
+from grass.lib.vector import *
+from grass.lib.raster import *
+import grass.pygrass.raster as raster
+import grass.pygrass.vector as Vector
+from grass.pygrass.gis.region import Region
+
+###############################################################################
+
+def write_Polygon(puntos,name):
+
+ # Write Polygon
+ sal_linea="B "+str(len(puntos)+1)+" 1\n"
+ for pp in puntos:
+ sal_linea+=''+str(pp[0])+" "+str(pp[1])+"\n"
+ sal_linea+=''+str(puntos[0][0])+" "+str(puntos[0][1])+"\n"
+ sal_linea+="1 1\n"
+ #print sal_linea
+ grass.write_command('v.in.ascii', flags='nz', output=name, stdin=sal_linea, input='-', format='standard', overwrite=True, quiet=True)
+ return 0
+
+
+###############################################################################
+
+def azimut(a,b,c,d):
+
+ Az=atan((c-a)/(d-b))
+ #if a < c and b < d: Az=Az*1 # Az>0 -> Az > 0
+ if a < c and b > d: Az=Az+pi # Az<0 -> Az > 100
+ elif a > c and b > d: Az=Az+pi # Az>0 -> Az > 200
+ elif a > c and b < d: Az=Az+2*pi # Az<0 -> Az > 300
+ return Az
+
+
+def discr_Lines(lines,discr):
+
+ disLines=[]
+ for pto in lines:
+ resto,acum=0,0
+ cat_r=1
+ disptos=[]
+ for i in range(len(pto[:-1])):
+
+ disptos.append(pto[i])
+ Az=azimut(pto[i][0],pto[i][1],pto[i+1][0],pto[i+1][1])
+ Lrecta=sqrt((pto[i+1][0]-pto[i][0])**2+(pto[i+1][1]-pto[i][1])**2)
+ Ptos_recta,resto,acum,cat_r=get_pts_Recta(pto[i][0],pto[i][1],0,Az,discr-resto,Lrecta,discr,acum,cat_r,i)
+ disptos.extend(Ptos_recta)
+ disptos.append(pto[i+1])
+ disLines.append(disptos)
+ return disLines
+
+def get_pts_Recta(xo,yo,zo,Az,Ini,Fin,L_int,Lacum,cat,ali):
+ # Return matriz de puntos, resto
+ if L_int==0: L_int=1
+ M=[]
+ while Ini <= Fin:
+ x2=xo+Ini*sin(Az)
+ y2=yo+Ini*cos(Az)
+ Lacum+=L_int
+ M.append([x2,y2,0,cat,Lacum,Az,'Recta',ali])
+ cat=cat+1
+ Ini+=L_int
+ return M,Fin-(Ini-L_int),Lacum,cat
+
+
+def drape_LinesPoints(lines,dem):
+
+ elev = raster.RasterRowIO(dem)
+ elev.open('r')
+ region=Region()
+ salida=[]
+ for i,line in enumerate(lines):
+
+ ptos=[]
+ for j,pto in enumerate(line):
+ #print j
+ pto_col = int((pto[0] - region.west)/region.ewres)
+ pto_row = int((region.north - pto[1])/region.nsres)
+ ptos.append(lines[i][j][:2]+[elev[pto_row][pto_col]]+lines[i][j][3:])
+
+ salida.append(ptos)
+ elev.close()
+ return salida
+
+###############################################################################
+
+def read_MapTrans(MapTrans):
+
+ verti=grass.read_command('v.out.ascii', input=MapTrans, output='-', format='standard', layer=1)
+ verti = re.findall(r'L [0-9]* [0-9]*\n(.*?)1 ([0-9]*)',verti,re.S|re.I)
+ verti.sort(key=lambda x: float(x[1]))
+ lineas,cats=[],[]
+ for line in verti:
+ lineas.append([d.split() for d in line[0].splitlines(0)])
+ cats.append(line[1])
+
+ for i,line in enumerate(lineas):
+ lineas[i]=[[float(p) for p in lin] for lin in line[:-1]]
+
+ return lineas,cats
+
+
+def get_Banks(MapTrans):
+
+ lineas,cats=read_MapTrans(MapTrans)
+
+ banks=[]
+ for i,cat in enumerate(cats):
+
+ dx1 = sqrt((lineas[i][0][0]-lineas[i][2][0])**2+(lineas[i][0][1]-lineas[i][2][1])**2)
+ dx2 = sqrt((lineas[i][0][0]-lineas[i][4][0])**2+(lineas[i][0][1]-lineas[i][4][1])**2)
+
+ banks.append([round(dx1,2),round(dx2,2)])
+
+ return banks
+
+
+
+###############################################################################
+
+def write_HecRasGeom2(MapTrans,Dem,namesal):
+
+ # puntos terreno transversales
+ lineas,cats = read_MapTrans(MapTrans)
+ lineas=discr_Lines(lineas,1)
+ lineas = drape_LinesPoints(lineas,Dem)
+
+ # Primer filtro. Eliminar ptos repetidos
+ lineas2,distcota=[],[]
+ for i,line in enumerate(lineas):
+ line2=[]
+ distcota.append([])
+ for j,pto in enumerate(line[:-1]):
+
+ dx1 = sqrt((pto[0]-line[j+1][0])**2+(pto[1]-line[j+1][1])**2)
+ dx2 = sqrt((pto[0]-line[0][0])**2+(pto[1]-line[0][1])**2)
+ if dx1 != 0:
+ distcota[-1].append([round(dx2,2),round(pto[2],2)])
+ line2.append([round(pto[0],2),round(pto[1],2),round(pto[2],2)])
+
+ line2.append([round(p,2) for p in line[j+1] ])
+ lineas2.append(line2)
+
+ # Pto mas bajo de cada seccion. cauce real
+ creal=[]
+ for j,line in enumerate(lineas2):
+ cotaline=[p[2] for p in line]
+ minline=min(cotaline)
+ ind = cotaline.index(minline)
+ creal.append(lineas2[j][ind])
+ sal4=""
+ for j,pto in enumerate(creal):
+ sal4+=str(pto[0])+'\t'+str(pto[1])+"\r\n"
+
+ # Ptos del eje de cada seccion
+ eje = [p[4] for p in lineas]
+ sal3=""
+ for j,pto in enumerate(eje):
+ sal3+=str(pto[0])+'\t'+str(pto[1])+"\r\n"
+
+
+ sal="RS,X,Y,Z\r\n"
+ for j,line in enumerate(lineas2):
+ for pto in line:
+ sal+=str(int(cats[-1])-int(cats[j]))+","+str(pto[0])+","+str(pto[1])+","+str(pto[2])+"\r\n"
+ sal2="RS,ST,Z\r\n"
+ for j,line in enumerate(distcota):
+ for pto in line:
+ sal2+=str(int(cats[-1])-int(cats[j]))+","+str(pto[0])+","+str(pto[1])+"\r\n"
+
+
+ with open(namesal+'_format_xyz.csv', 'wb') as f:
+ f.write(sal)
+ f.closed
+ with open(namesal+'_format_stations.csv', 'wb') as f:
+ f.write(sal2)
+ f.closed
+ with open(namesal+'_schematic.txt', 'wb') as f:
+ f.write(sal3)
+ f.closed
+ with open(namesal+'_schematic_real.txt', 'wb') as f:
+ f.write(sal4)
+ f.closed
+ #print sal
+ return 0
+
+
+
+
+##############################################################################
+
+def write_HecRasTables2(MapTrans,MapCrossDispl,manning,namesal):
+
+ # puntos desplazados
+ verti=grass.read_command('v.out.ascii', input=MapCrossDispl, output='-', format='point', columns='Type,PK', layer=1)
+ verti=[p.split('|') for p in verti.splitlines(0)]
+
+ linea_ptos=[]
+ for pto in verti:
+ if float(pto[-1]) == 0: linea_ptos.append([])
+ linea_ptos[-1].append([float(p) for p in pto[:3]]+[float(pto[4][1:])]+[float(pto[5])])
+
+ banks = get_Banks(MapTrans)
+ bankstr = ""
+ for i,bank in enumerate(banks):
+
+ bankstr += str(bank[0])+'\t'+str(bank[1])+"\r\n"
+
+ # Pks transversales
+ pks_trans = grass.read_command('v.db.select', map=MapTrans, columns='pk')
+ pks_trans = pks_trans.splitlines(0)
+ pks_trans = [float(p.replace('+','')) for p in pks_trans[1:]]
+
+ manning=manning.replace(',',' \t')
+ manningsal=manning+'\r\n'
+ distran = ""
+ for i,punt in enumerate(linea_ptos[0][:-1]):
+ manningsal+= manning+'\r\n'
+ distran += str(linea_ptos[0][i+1][4]-linea_ptos[0][i][4])+'\t'+str(pks_trans[i+1]-pks_trans[i])+'\t'+str(linea_ptos[-1][i+1][4]-linea_ptos[-1][i][4])+'\r\n'
+
+ with open(namesal+'_lengths.txt', 'wb') as f:
+ f.write(distran)
+ f.closed
+ with open(namesal+'_manning.txt', 'wb') as f:
+ f.write(manningsal)
+ f.closed
+ with open(namesal+'_banks.txt', 'wb') as f:
+ f.write(bankstr)
+ f.closed
+
+ return 0
+
+###############################################################################
+
+
+def import_HecRas(nament,namesal):
+
+ with open(nament, 'r') as f:
+ fileH=f.read()
+ f.closed
+
+ eje = re.findall(r'CENTERLINE:(.*?)END:',fileH,re.S|re.I)
+ eje = eje[0].split(', ,')
+ eje = [d.split(',') for d in eje]
+ eje = eje[0:-1]
+ for j,punt in enumerate(eje):
+ eje[j]= [d.strip() for d in punt]
+ #print eje
+
+ cotas = re.findall(r'WATER ELEVATION:(.*)',fileH)
+ cotas = [d.strip() for d in cotas]
+
+ inund=''
+ for i,p in enumerate(eje):
+ inund+= p[0]+' '+p[1]+' '+cotas[i]+' '+str(len(eje)-i)+'\n'
+ inund+= eje[-1][0]+' '+eje[-1][0]+' '+cotas[-1]+' 1\n'
+
+ with open(options['outprofilefile'], 'w') as f:
+ f.write(inund)
+ f.closed
+
+ coor = re.findall(r'WATER SURFACE EXTENTS:(.*?)END:',fileH,re.S|re.I)
+ coor = [d.split(',') for d in coor]
+ for j,punt in enumerate(coor):
+ coor[j]= [d.strip() for d in punt]
+
+ polygon=[]
+ polygon2=[]
+ for j,punt in enumerate(coor):
+
+ polygon.append(coor[j][0:2])
+ polygon2.append(coor[j][2:])
+
+ polygon2.reverse()
+ for j,punt in enumerate(polygon2):
+ polygon.append(punt)
+ for j,punt in enumerate(polygon):
+ polygon[j].append('0.0')
+ #print polygon
+ #print coor
+ write_Polygon(polygon,namesal)
+
+ return 0
+
+
+
+
+# ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### #
+# ### Main
+# ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### #
+
+def main():
+
+ MapTrans=options['mapcross']
+ MapCrossDispl=options['mapcrossdispl']
+ Dem=options['dem']
+
+ #### Hec-Ras section ####
+
+ if flags['e']:
+
+ write_HecRasGeom2(MapTrans,Dem,options['output'])
+
+ write_HecRasTables2(MapTrans,MapCrossDispl,options['manning'],options['output'])
+
+
+
+ if flags['i']:
+
+ infile = options['inputfile']
+ import_HecRas(infile,options['outwaterfile'])
+
+
+ sys.exit(0)
+
+if __name__ == "__main__":
+ options, flags = grass.parser()
+ main()
Property changes on: grass-addons/grass7/vector/v.civil/v.civil.river/v.civil.river.py
___________________________________________________________________
Added: svn:mime-type
+ text/x-python
Added: svn:eol-style
+ native
Added: grass-addons/grass7/vector/v.civil/v.civil.road/Makefile
===================================================================
--- grass-addons/grass7/vector/v.civil/v.civil.road/Makefile (rev 0)
+++ grass-addons/grass7/vector/v.civil/v.civil.road/Makefile 2014-02-26 00:37:16 UTC (rev 59141)
@@ -0,0 +1,9 @@
+# fix this relative to include/
+# or use absolute path to the GRASS source code
+MODULE_TOPDIR = ../..
+
+PGM = v.civil.road
+
+include $(MODULE_TOPDIR)/include/Make/Script.make
+
+default: script
Property changes on: grass-addons/grass7/vector/v.civil/v.civil.road/Makefile
___________________________________________________________________
Added: svn:mime-type
+ text/x-makefile
Added: svn:eol-style
+ native
Added: grass-addons/grass7/vector/v.civil/v.civil.road/v.civil.road.html
===================================================================
--- grass-addons/grass7/vector/v.civil/v.civil.road/v.civil.road.html (rev 0)
+++ grass-addons/grass7/vector/v.civil/v.civil.road/v.civil.road.html 2014-02-26 00:37:16 UTC (rev 59141)
@@ -0,0 +1,315 @@
+<h2>DESCRIPTION</h2>
+
+<em>v.civil.road</em> generate a topography edge for design roads, channels
+and ports.
+<p>
+Applications:
+
+<ul>
+ <li>Roads, streets, railroads, walls,...</li>
+ <li>Channels, reservoirs, supply and waste water conduction, ...</li>
+ <li>Ports, </li>
+</ul>
+
+<h2>USAGE</h2>
+
+For run the program, a map with one polyline (horizontal polygon or edge), is
+required. If a table in layer 1 exist, this will no be overwrite, this layer is
+free for store information about the edge by the user.
+
+
+<p>
+There are three options for run the program, Create/Update edge, Update
+solution and Insert point in edge layer.
+
+
+<p>
+The <em>-n</em> flag (Create/Update new edge), will convert the polyline
+in a edge, i.e., four tables will be add to the polyline layers map, attached
+to points in its respective layers, where the parameter of the road will be
+stored . This flag must be selected for create a edge from a polyline or if
+the vertices of the edge are moved or added another ones.
+
+<p>
+The <em>-u</em> flag "Update solution" will write the selected maps,
+with its additional parameters for some of them, and with the parameters
+included in the tables of the edge.
+
+<p>
+The name of the outs maps, by default are selected. If the out name map is
+preceded with character "_", this name will be added to the edge name map.
+This is useful for change of edges quickly.
+
+<p>
+The <em>-i</em> flag (Insert points), new points for defining the road,
+can be added to the layers _Vertical, _Section and _Trans, referencing by pk.
+
+<p>
+The user must select the layer and give a list with the pks. The options
+of the fist pk where the new will be inserted, will be copied in the new.
+
+<p>
+The other options are described below.
+
+
+<h3>Plan</h3>
+
+
+<p>
+The <em>-y</em> flag write a 3D polyline edge, in horizontal alignment, with
+the name in <em>plantpoly</em> option (default: <em>NameEdge</em>_plantpoly).
+
+<p>
+The <em>-h</em> flag write the edge in horizontal alignment, with the name
+in <em>plant</em> option. The map created, contains the segments of the
+alignment and its information is stored in three tables.(default: <em>NameEdge</em>_plan).
+
+<p>
+The <em>-k</em> flag write pk marks in horizontal, with the name in
+<em>pks</em> option.(default: <em>NameEdge</em>_Pks).
+
+<p>
+The <em>pkopt</em> option define the distance between normal marks (npk),
+distance between principal marks (mpk), longitude of normal marks and
+longitude of princicipal marks.
+
+<p>
+The <em>-d</em> flag write the displaced platform lines in horizontal,
+with the name in <em>displ</em> option.
+
+<p>
+The <em>-a</em> flag write an area map with displaced platform lines in
+horizontal, with the name in <em>displ_area</em> option.
+
+
+
+<h3>Vertical</h3>
+
+
+<p>
+The <em>-v</em> flag write the vertical alignment in horizontal, with the name
+in <em>raised</em> option. The map created, contains the segments of the
+alignment and its information is stored in three tables.(default: <em>NameEdge</em>_Vert)
+
+<p>
+The <em>-l</em> flag write the raised of the edge (vertical alignment), the
+terrain profile in a longitudinal profile map, and vertical polygon with the
+names given in <em>LPras</em> ,<em>LPterrain</em> and <em>LPejeref</em> options
+(default: <em>NameEdge</em>_LP_Ras, <em>NameEdge</em>_LP_terrain and
+<em>NameEdge</em>_LP_Edge).
+
+<p>
+The <em>-m</em> flag write the coordinates edges and the profile guitar in the
+longitudinal profile with its marks, with the names given in <em>LPedgeX</em>,
+<em>LPedgeX_marks</em>, <em>LPedgeY</em> and <em>LPedgeY_marks</em> options
+(default: <em>NameEdge</em>_LP_EdgeX, <em>NameEdge</em>_LP_EdgeXmarks,
+<em>NameEdge</em>_LP_EdgeY and <em>NameEdge</em>_LP_EdgeYmarks).
+
+<p>
+With the option <em>LPScale</em>, the vertical scale of the longitudinal profile can be
+given.
+
+<p>
+The <em>LPopt</em> option give the longitude of marks, distance between marks in edge x
+and y, and distance between lines of the guitar.
+
+
+
+
+
+<h3>Cross</h3>
+
+
+
+<p>
+The <em>-c</em> flag write the projection of cross section in horizontal, with the name
+in <em>cross</em> option, and with the lines selected in <em>cross_opt</em> option.
+
+<p>
+The <em>-r</em> flag write a map with points of intersection of displaced lines with
+cross-section lines, in horizontal, with the name in <em>crossdispl</em> option. The
+lines that will be crossed are selected in <em>displ_opt</em>
+
+
+
+<p>
+The <em>-f</em> flag write the raised of cross-section (cross alignment) and the
+terrain profile in a traversal profile map with the
+names given in <em>LTras</em> and <em>LTterrain</em> options
+(default: <em>NameEdge</em>_TP_Ras and <em>NameEdge</em>_TP_terrain
+
+<p>
+The <em>-m</em> flag write the coordinates edges and the profile guitar in the
+longitudinal profile with its marks, with the names given in <em>LTedgeX</em> and
+<em>LTedgeY</em> options (default: <em>NameEdge</em>_TP_EdgeX and
+<em>NameEdge</em>_TP_EdgeY ).
+
+<p>
+With the option <em>LTScale</em>, the vertical scale of the longitudinal profile can be
+given.
+
+<p>
+The <em>LTopt</em> option give the longitude of marks, distance between marks in edge x
+and y, and distance between lines of the guitar.
+
+<p>
+The <em>LTopt2</em> option give the number of rows of cross-section to display,
+distance between cross-section in edge x and y.
+
+
+<h3>Terrain</h3>
+
+
+<p>
+The <em>-t</em> flag write the projection of edge on a DEM (3D line of terrain)
+in horizontal, with the name in <em>outtlong</em> option.
+
+<p>
+The <em>-q</em> flag write the projection of transects lines on a DEM (3D lines of terrain)
+in horizontal, with the name in <em>outtcross</em> option.
+
+<p>
+The <em>-s</em> flag write the slope soil lines calculated with the last
+displaced lines (perpendicular to the edge) in horizontal, with the name
+in <em>outslope</em> option
+
+<p>
+The <em>-e</em> flag write the areas map with slope soil lines and the last
+displaced lines in horizontal, with the name in <em>outslopeareas</em> option
+
+<p>
+The <em>-p</em> flag write a 3D points map of lines selected in <em>pts_opt</em> option
+in horizontal, with the name in <em>outpoints</em> option
+
+<p>
+The <em>-b</em> flag write a 3D lines map of lines selected in <em>break_opt</em> option
+in horizontal, with the name in <em>outbreak</em> option
+
+<p>
+The <em>-o</em> flag write a hull area map of lines selected in <em>hull_opt</em> option
+in horizontal, with the name in <em>outhull</em> option
+
+
+<h2>Tables</h2>
+
+<p>
+Geometric roadway design can be broken into four main parts: horizontal
+alignment, profile (vertical alignment), platforms definition and cross-section.
+
+There are maps in horizontal
+alignment for display in plan, maps in vertical alignment that requires a
+new display for visualization in profile view, and maps in cross-section that
+requires a new display for visualization in profile cross-section view.
+
+<p>
+The tables stored in the edge map are:
+
+<h3>Layer 1:</h3>
+<p>
+Created with only cat column (if not exist), is for free use to add columns with information
+of the edge.
+
+
+
+<h3>Layer 2:</h3>
+<p>
+<b>"NameEdge"_Horizontal</b>, for insert the parameters of the horizontal
+alignment. This layer have points in all vertices of the horizontal polygon.
+No more points can be added to this layer.
+
+<p>
+The columns for editing by the user are:
+
+ <ul>
+ <li><b>radio</b>: radio of curve
+ <p>+ for clockwise</p>
+ <p>- for anticlockwise</p>
+ </li>
+ <li><b>a_in</b>: Parameter A of input Clothoid</li>
+ <li><b>a_out</b>: Parameter A of output Clothoid</li>
+ <li><b>widening</b>: Widening of curve (where this widening is growing
+ in the clothoid, only in mode exact)</li>
+ </ul>
+
+
+
+<h3>Layer 3:</h3>
+<p>
+<b>"NameEdge"_Vertical</b>, for insert the parameters of the vertical
+alignment.This layer are created with the first an last vertices of the
+horizontal polygon. New point can be added by the flag i (Insert point),
+or with edit gui. The added points, are inserted in the edge referencing by pk.
+
+<p>
+The columns for editing by the user are:
+
+ <ul>
+ <li><b>pk</b>: kilometric point of the edge</li>
+ <li><b>elev</b>: Elevation of the vertice of the vertical alignment</li>
+ <li><b>kv</b>: Parameter Kv of the vertical alignment</li>
+ <li><b>l</b>: Leng of parabolic curve (no yet implemented)</li>
+ <li><b>B</b>: Height of vertice of vertical polygon to the parabolic
+ curve (no yet implemented)</li>
+ </ul>
+
+
+
+
+<h3>Layer 4: </h3>
+
+<b>"NameEdge"_Section</b>, for insert the parameters of the platform.
+This layer are created with the first an last vertices of the horizontal
+polygon. New point can be added by the flag i (Insert point), or with edit
+gui. The added points, are inserted in the edge referencing by pk.
+<p>
+The columns for editing by the user are:
+
+ <ul>
+ <li><b>pk</b>: kilometric point of the edge</li>
+ <li><b>sec_left, sec_right</b>: For defining left and right platform
+ lines (distance and height to the edge), separated by ";" </li>
+ <li><b>type_left, type_right</b>: Types of left and right lines
+ (types are: l,e,r,0). You can enter the type of multiple lines
+ separated by ";"
+ <b>l</b>: linear approximation between definition points.
+ <b>e</b>: displaced edge, mode exact.
+ <b>c</b>: ellipse between definition points.
+ <b>rR,A</b>: circle between definition points ((x-A)^2+y^2=R^2))</li>
+
+ <li><b>cut_left, cut_right</b>: Cut slope left and right </li>
+ <li><b>fill_left, fill_right</b>: Fill slope left and right </li>
+ <li><b>superelev_left,superelev_right</b>: Superelevation left and right
+ lanes (no yet implemented)</li>
+ </ul>
+<p>
+Between two points, if the second input distance equal zero, the line will be
+stopped. If the second input distance equal -1, this point will not be processed
+and the next one will be considered.
+
+
+
+<h3>Layer 5:</h3>
+<p>
+<b>"NameEdge"_Tran</b>s, for insert the parameters of transects to the edge. This layer
+are created with the first an last vertices of the horizontal polygon. New point can be added by the flag
+i (Insert point), or with edit gui. The added points, are inserted in the edge referencing by pk
+<p>
+The columns for editing by the user are:
+
+ <ul>
+ <li><b>pk</b>: kilometric point of the edge</li>
+ <li><b>dist_left,dist_right</b>: Distance left and right from the edge</li>
+ <li><b>npk</b>: Distance between trans</li>
+ </ul>
+
+
+
+
+<h2> AUTHOR</h2>
+<p>
+Jesus Fernandez-Capel Rosillo<br >
+Civil Engineer, Spain<br >
+jfc at alcd net<br >
+
+
+
Property changes on: grass-addons/grass7/vector/v.civil/v.civil.road/v.civil.road.html
___________________________________________________________________
Added: svn:mime-type
+ text/html
Added: svn:keywords
+ Author Date Id
Added: svn:eol-style
+ native
Added: grass-addons/grass7/vector/v.civil/v.civil.road/v.civil.road.py
===================================================================
--- grass-addons/grass7/vector/v.civil/v.civil.road/v.civil.road.py (rev 0)
+++ grass-addons/grass7/vector/v.civil/v.civil.road/v.civil.road.py 2014-02-26 00:37:16 UTC (rev 59141)
@@ -0,0 +1,3307 @@
+#!/usr/bin/env python
+# -*- coding: utf-8
+############################################################################
+#
+# MODULE: v.civil.road, v0.5.9
+#
+# AUTHOR: Jesús Fernández-Capel Rosillo
+# Civil Engineer, Spain
+# jfc at alcd net
+#
+# PURPOSE: Desing roads, channels, ports...
+#
+# COPYRIGHT: (c) 2014 Jesús Fernández-Capel Rosillo.
+#
+# This program is free software under the GNU General Public
+# License (>=v2). Read the file COPYING that comes with GRASS
+# for details.
+#
+#############################################################################
+
+#%Module
+#% description: Generate a topography edge for desing roads, channels, ports...
+#% keywords: vector, ROADS, CHANNELS, PORTS.
+#%End
+
+#### Required section ####
+
+
+#%option G_OPT_V_INPUT
+#% key: edge
+#% type: string
+#% gisprompt: old,vector,vector
+#% description: Name for edge (horizontal polygon)
+#% required: yes
+#%end
+
+
+#### Plant section ####
+
+#%flag
+#% key: y
+#% description: Write edge polyline
+#% guisection: Plan
+#%end
+
+#%flag
+#% key: h
+#% description: Write horizontal alignment
+#% guisection: Plan
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: plantpoly
+#% description: Name for output edge polyline
+#% required: no
+#% answer: _Poly
+#% guisection: Plan
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: plant
+#% description: Name for output horizontal alignment
+#% required: no
+#% answer: _Plant
+#% guisection: Plan
+#%end
+
+####
+
+#%flag
+#% key: k
+#% description: Write Pks marks
+#% guisection: Plan
+#%end
+
+##%option
+##% key: numpk
+##% type: integer
+##% description: Specific pk to draw
+##% required: no
+##% answer: 0
+##% guisection: Plan
+##%end
+
+#%option
+#% key: pkopt
+#% type: string
+#% description: Pks marks options values. (npk,mpk,dist,m)
+#% required: no
+#% answer: 20,100,2,1
+#% guisection: Plan
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: pks
+#% description: Name of Pks marks
+#% required: no
+#% answer: _Pks
+#% guisection: Plan
+#%end
+
+####
+
+#%flag
+#% key: d
+#% description: Write plataform displaced lines
+#% guisection: Plan
+#%end
+
+#%flag
+#% key: a
+#% description: Write plataform areas
+#% guisection: Plan
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: displ
+#% description: Name of plataform displaced lines
+#% required: no
+#% answer: _Displ
+#% guisection: Plan
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: displ_area
+#% description: Name of plataform areas
+#% required: no
+#% answer: _Displ_Area
+#% guisection: Plan
+#%end
+
+
+
+#### Alz section ####
+
+#%flag
+#% key: v
+#% description: Write vertical alignment in horizontal
+#% guisection: Vertical
+#%end
+
+#%flag
+#% key: l
+#% description: Draw longitudinal profile (dem required)
+#% guisection: Vertical
+#%end
+
+#%flag
+#% key: m
+#% description: Draw longitudinal profile coord edges (dem required)
+#% guisection: Vertical
+#%end
+
+
+#%option
+#% key: LPScale
+#% type: double
+#% description: Long profile vertical scale (V/H, V/1)
+#% required: no
+#% answer: 8
+#% guisection: Vertical
+#%end
+
+#%option
+#% key: LPopt
+#% type: string
+#% description: Long profile Values for Longmark,distMark_x,distMark_y,DistGitarr.
+#% required: no
+#% answer: 2,20,1,20
+#% guisection: Vertical
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: raised
+#% description: Name of vertical alignment
+#% required: no
+#% answer: _Vertical
+#% guisection: Vertical
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: LPedgeX
+#% description: Name of long profile edge x
+#% required: no
+#% answer: _LP_EdgeX
+#% guisection: Vertical
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: LPedgeXmarks
+#% description: Name of long profile edge x marks
+#% required: no
+#% answer: _LP_EdgeXmarks
+#% guisection: Vertical
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: LPedgeY
+#% description: Name of long profile edge y
+#% required: no
+#% answer: _LP_EdgeY
+#% guisection: Vertical
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: LPedgeYmarks
+#% description: Name of long profile edge y marks
+#% required: no
+#% answer: _LP_EdgeYmarks
+#% guisection: Vertical
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: LPterrain
+#% description: Name of long profile terrain
+#% required: no
+#% answer: _LP_Terr
+#% guisection: Vertical
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: LPras
+#% description: Name of long profile vertical alignment
+#% required: no
+#% answer: _LP_Ras
+#% guisection: Vertical
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: LPejeref
+#% description: Name of long profile vertical polygon
+#% required: no
+#% answer: _LP_Edge
+#% guisection: Vertical
+#%end
+
+#### Section section ####
+
+
+
+#### Cross section ####
+
+#%flag
+#% key: c
+#% description: Write cross section
+#% guisection: Cross
+#%end
+
+#%flag
+#% key: r
+#% description: Write intersection cross section-displaced
+#% guisection: Cross
+#%end
+
+
+
+#%flag
+#% key: f
+#% description: Draw cross section (dem required)
+#% guisection: Cross
+#%end
+
+#%flag
+#% key: g
+#% description: Draw cross section coord edges (dem required)
+#% guisection: Cross
+#%end
+
+
+
+#%option
+#% key: LTScale
+#% type: double
+#% description: Cross section vertical scale (V/H, V/1)
+#% required: no
+#% answer: 2
+#% guisection: Cross
+#%end
+
+#%option
+#% key: LTopt
+#% type: string
+#% description: Cross section options values for Longmark,distMark_x,distMark_y.
+#% required: no
+#% answer: 1,5,2
+#% guisection: Cross
+#%end
+
+#%option
+#% key: LTopt2
+#% type: string
+#% description: Cross section options values for nrows,distTP_x,distTP_y.
+#% required: no
+#% answer: 5,20,20
+#% guisection: Cross
+#%end
+
+#%option G_OPT_V_TYPE
+#% key: cross_opt
+#% options: slope_left,displ_left,edge,displ_rigth,slope_rigth
+#% answer: slope_left,displ_left,edge,displ_rigth,slope_rigth
+#% required: no
+#% description: Add lines to the cross section
+#% guisection: Cross
+#%end
+
+#%option G_OPT_V_TYPE
+#% key: displ_opt
+#% options: displ_left0,displ_left,displ_left-1,displ_rigth0,displ_rigth,displ_rigth-1
+#% required: no
+#% answer: displ_left0,displ_rigth-1
+#% description: Cross points section-Displaced
+#% guisection: Cross
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: cross
+#% description: Name of cross section
+#% required: no
+#% answer: _Cross
+#% guisection: Cross
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: crossdispl
+#% description: Name of intersection cross section - displaced
+#% required: no
+#% answer: _CrossDispl
+#% guisection: Cross
+#%end
+
+
+
+#%option G_OPT_V_OUTPUT
+#% key: LTedgeX
+#% description: Name of cross section coord edge x
+#% required: no
+#% answer: _TP_EdgeX
+#% guisection: Cross
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: LTedgeXmarks
+#% description: Name of cross section coord edge x
+#% required: no
+#% answer: _TP_EdgeXmarks
+#% guisection: Cross
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: LTedgeY
+#% description: Name of cross section coord edge y
+#% required: no
+#% answer: _TP_EdgeY
+#% guisection: Cross
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: LTedgeYmarks
+#% description: Name of cross section coord edge y
+#% required: no
+#% answer: _TP_EdgeYmarks
+#% guisection: Cross
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: LTterrain
+#% description: Name of cross sections terrain
+#% required: no
+#% answer: _TP_Terr
+#% guisection: Cross
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: LTras
+#% description: Name of cross section plataform
+#% required: no
+#% answer: _TP_Ras
+#% guisection: Cross
+#%end
+
+
+#### Terrain section ####
+
+#%flag
+#% key: t
+#% description: Write terrain edge
+#% guisection: Terr
+#%end
+
+#%flag
+#% key: q
+#% description: Write terrain cross section
+#% guisection: Terr
+#%end
+
+#%flag
+#% key: s
+#% description: Write slope soil
+#% guisection: Terr
+#%end
+
+#%flag
+#% key: e
+#% description: Write slope soil areas
+#% guisection: Terr
+#%end
+
+#%flag
+#% key: p
+#% description: Write points
+#% guisection: Terr
+#%end
+
+#%flag
+#% key: b
+#% description: Write breaklines
+#% guisection: Terr
+#%end
+
+#%flag
+#% key: o
+#% description: Write hull
+#% guisection: Terr
+#%end
+
+#%option G_OPT_V_TYPE
+#% key: pts_opt
+#% options: slope_left,displ_left,edge,displ_rigth,slope_rigth
+#% answer: slope_left,displ_left,edge,displ_rigth,slope_rigth
+#% required: no
+#% description: Points
+#% guisection: Terr
+#%end
+
+#%option G_OPT_V_TYPE
+#% key: break_opt
+#% options: slope_left,displ_left,edge,displ_rigth,slope_rigth
+#% answer: slope_left,displ_left,edge,displ_rigth,slope_rigth
+#% required: no
+#% description: Breaklines
+#% guisection: Terr
+#%end
+
+#%option G_OPT_V_TYPE
+#% key: hull_opt
+#% options: slope_left,slope_rigth
+#% answer: slope_left,slope_rigth
+#% required: no
+#% description: Hull
+#% guisection: Terr
+#%end
+
+#%option
+#% key: dem
+#% type: string
+#% gisprompt: old,cell,raster
+#% key_desc: raster dem
+#% description: Name of DEM raster
+#% guisection: Terr
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: outtlong
+#% description: Name of terrain edge
+#% required: no
+#% answer: _LongTerr
+#% guisection: Terr
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: outtcross
+#% description: Name of terrain cross section
+#% required: no
+#% answer: _CrossTerr
+#% guisection: Terr
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: outslope
+#% description: Name of slope soil
+#% required: no
+#% answer: _Slope
+#% guisection: Terr
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: outslopeareas
+#% description: Name of slope soil areas
+#% required: no
+#% answer: _SlopeAreas
+#% guisection: Terr
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: outpoints
+#% description: Name of points
+#% required: no
+#% answer: _Points
+#% guisection: Terr
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: outbreak
+#% description: Name of breaklines
+#% required: no
+#% answer: _Breaklines
+#% guisection: Terr
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: outhull
+#% description: Name of hull
+#% required: no
+#% answer: _Hull
+#% guisection: Terr
+#%end
+
+
+#### Options section ####
+
+#%flag
+#% key: n
+#% description: Create/Update new edge
+#%end
+
+#%flag
+#% key: u
+#% description: Update solution
+#%end
+
+#%flag
+#% key: i
+#% description: Insert point in edge layers (vertical, section an transversal)
+#%end
+
+#%option
+#% key: pklayer
+#% type: string
+#% label: Layer to add Pks points
+#% options: Vertical,Section,Trans
+#% required: no
+#%end
+
+#%option
+#% key: pklist
+#% type: string
+#% description: List of Pks points.
+#% required: no
+#%end
+
+#%option
+#% key: intervR
+#% type: integer
+#% description: Interval calculation
+#% required: no
+#% answer: 1
+#%end
+
+#%option
+#% key: intervC
+#% type: integer
+#% description: Interval calculation in curves
+#% required: no
+#% answer: 1
+#%end
+
+import os, sys
+
+from math import *
+import grass.script as g
+from grass.lib.gis import *
+from grass.lib.vector import *
+from grass.lib.raster import *
+import numpy as np
+#from grass.pygrass.raster import RasterNumpy
+#from grass.pygrass.vector import VectorTopo
+#from grass.pygrass.vector.geometry import Point
+#from grass.pygrass.functions import pixel2coor
+import grass.pygrass.raster as raster
+#import grass.pygrass.gis as gis
+from grass.pygrass.gis.region import Region
+from grass.script import array as garray
+from itertools import groupby
+import time
+#import Gnuplot
+
+#### Alineaciones
+
+def aprox_coord(L, Tau):
+
+ n_iter=10;x=0;y=0
+ for n in range(n_iter):
+ x+=(((-1)**n*Tau**(2*n))/((4*n+1)*factorial(2*n)))
+ y+=(((-1)**n*Tau**(2*n+1))/((4*n+3)*factorial(2*n+1)))
+ x=x*L
+ y=y*L
+ return [x,y]
+
+def aprox_coord2(R, Tau):
+
+ n_iter=10;x=0;y=0
+
+ for n in range(n_iter):
+ x+=((-1)**n*(2*Tau**(2*n+1)/((4*n+1)*factorial(2*n))))-(-1)**n*(Tau**(2*n+1)/factorial(2*n+1))
+ y+=((-1)**n*(2*Tau**(2*n+2)/((4*n+3)*factorial(2*n+1))))+((-1)**n*(Tau**(2*n)/factorial(2*n)))
+ x=x*R
+ y=y*R-R
+ return [x,y]
+
+def cambio_coord(x,y,Az,c,d,R):
+
+ if R>0:
+ x1=c-x*sin(Az)+y*cos(Az)
+ y1=d-x*cos(Az)-y*sin(Az)
+ else:
+ x1=c-x*sin(Az)+y*cos(Az)
+ y1=d-x*cos(Az)-y*sin(Az)
+ return [x1,y1]
+
+def azimut(a,b,c,d):
+
+ if a > c and b == d: Az=3*pi/2
+ elif a < c and b == d: Az=pi/2
+ elif a == c and b > d: Az=pi
+ elif a == c and b < d: Az=2*pi
+ elif a == c and b == d: Az=0
+ else: Az=atan((c-a)/(d-b))
+ #if a > c and b == d: Az=Az+pi # Az>0 -> Az > 0
+ if a < c and b > d: Az=Az+pi # Az<0 -> Az > 100
+ elif a > c and b > d: Az=Az+pi # Az>0 -> Az > 200
+ elif a > c and b < d: Az=Az+2*pi # Az<0 -> Az > 300
+ return Az
+
+
+def format_Pk(pk):
+
+ if str(pk).find(".")==-1:
+ integer = str(pk)
+ decimal= "0"
+ else:
+ integer, decimal = str(pk).split(".")
+ integer = re.sub(r"\B(?=(?:\d{3})+$)", "+", "{:0>4}".format(integer))
+ return "'"+integer + "." + decimal + "'"
+
+
+def get_pts_Clot_ent(A,xad,yad,Az,Ini,Fin,L_int,R,Lacum,cat,ali):
+ # Return matriz de puntos, resto
+ M = []
+ while Ini <= Fin:
+ Rclo=A**2/(Ini)
+ Tau_clo=(Ini)/(2*Rclo)
+ xclo,yclo=aprox_coord((Ini),Tau_clo)
+ if R>0:
+ x_clo=xad-xclo*sin(-Az)+yclo*cos(-Az)
+ y_clo=yad+xclo*cos(-Az)+yclo*sin(-Az)
+ Az1=Az+Tau_clo
+ elif R<0:
+ x_clo=xad+xclo*sin(Az)-yclo*cos(Az)
+ y_clo=yad+xclo*cos(Az)+yclo*sin(Az)
+ Az1=Az-Tau_clo
+ Lacum+=L_int
+ M.append([x_clo,y_clo,0,cat,Lacum,Az1,'Clot_in',ali])
+ cat=cat+1
+ Ini=Ini+L_int
+ return M,Fin-(Ini-L_int),Lacum,cat
+
+def get_pts_Clot_sal(A,xda,yda,Az,Ini,Fin,L_int,R,Lacum,cat,ali):
+ # Return matriz de puntos, resto
+ M = []
+ while Ini <= Fin:
+ #if Fin-Ini==0: Fin=0; Ini=Ini*-1
+ Rclo=A**2/(Fin)
+ Tau_clo=(Fin)/(2*Rclo)
+ xclo,yclo=aprox_coord((Fin),Tau_clo)
+ if R>0:
+ x_clo=xda-xclo*sin(Az)+yclo*cos(Az)
+ y_clo=yda-xclo*cos(Az)-yclo*sin(Az)
+ Az1=Az-Tau_clo
+ elif R<0:
+ x_clo=xda+xclo*sin(-Az)-yclo*cos(-Az)
+ y_clo=yda-xclo*cos(-Az)-yclo*sin(-Az)
+ Az1=Az+Tau_clo
+ Lacum+=L_int
+ M.append([x_clo,y_clo,0,cat,Lacum,Az1,'Clot_out',ali])
+ cat=cat+1
+ Fin=Fin-L_int
+ #Ini=Ini+L_int
+ #print M,Ini-(Fin+L_int),Lacum,cat,Ini,Fin,L_int
+ return M,Ini+(Fin+L_int),Lacum,cat
+
+def get_pts_Curva(xc,yc,AzIni,AzFin,incremento,R,Lacum,dc,cat,ali):
+
+ M = []
+ incremento=incremento/abs(R)
+ if R>0:
+ if AzIni>AzFin: AzIni=AzIni-2*pi
+ ii = AzIni+incremento
+ while ii <= AzFin:
+ x1=xc+R*sin(ii)
+ y1=yc+R*cos(ii)
+ Az1=(ii+pi/2)
+ ii += incremento
+ Lacum+=incremento*abs(R)
+ if Az1>2*pi: Az1=Az1-2*pi
+ M.append([x1,y1,0,cat,Lacum,Az1,'Curve',ali])
+ cat=cat+1
+ resto=(AzFin-(ii-incremento))*abs(R)
+ elif R<0:
+ if AzIni<AzFin: AzFin=AzFin-2*pi
+ ii = AzIni-incremento
+ while ii >= AzFin:
+ x1=xc-R*sin(ii)
+ y1=yc-R*cos(ii)
+ Az1=(ii-pi/2)
+ ii -= incremento
+ Lacum+=incremento*abs(R)
+ if Az1<0: Az1=Az1+2*pi
+ M.append([x1,y1,0,cat,Lacum,Az1,'Curve',ali])
+ cat=cat+1
+ resto= ((ii+incremento)-AzFin)*abs(R)
+ #print M,resto,Lacum,cat
+ return M,resto,Lacum,cat
+
+def get_pts_Recta(xo,yo,zo,Az,Ini,Fin,L_int,Lacum,cat,ali):
+ # Return matriz de puntos, resto
+ if L_int==0: L_int=1
+ M=[]
+ while Ini <= Fin:
+ x2=xo+Ini*sin(Az)
+ y2=yo+Ini*cos(Az)
+ Lacum+=L_int
+ M.append([x2,y2,zo,cat,Lacum,Az,'Line',ali])
+ cat=cat+1
+ Ini+=L_int
+ return M,Fin-(Ini-L_int),Lacum,cat
+
+
+def get_PlantaXY(PK,puntos_eje):
+
+ acum=0
+ for i in range(0,len(puntos_eje)-1,1):
+
+ x_ini,y_ini,z_ini,Lrecta,LrAcum,Az_ent=puntos_eje[i][0]
+ Aent,Az_ent,xad,yad,xar,yar,Lent,LeAcum,R,Lini_e,xadp,yadp=puntos_eje[i][1]
+ R,alpha,xc,yc,Az_ini,Az_fin,Dc,LcAcum=puntos_eje[i][2]
+ Asal,Az_sal,xda,yda,xra,yra,Lsal,LsAcum,R,Lini_s,xdap,ydap=puntos_eje[i][3]
+ M=[[0,0]]
+
+ if 0 <= PK and PK <= LrAcum:
+ M,rest,acum,cat=get_pts_Recta(x_ini,y_ini,z_ini,Az_ent,PK-acum,PK-acum,1,PK-1,1,i); break
+
+ elif LrAcum < PK and PK <= LeAcum:
+ M,rest,acum,cat=get_pts_Clot_ent(Aent,xad,yad,Az_ent,PK-LrAcum,PK-LrAcum,1,R,PK-1,1,i); break
+
+ elif LeAcum < PK and PK <= LcAcum:
+ a=Az_ini+(PK-LeAcum)/R
+ M,rest,acum,cat=get_pts_Curva(xc,yc,a-1/R,a,1,R,PK-1,Dc,1,i); break
+
+ elif LcAcum < PK and PK <= LsAcum:
+ if LsAcum-PK == 0: distcond = PK
+ else: distcond = LsAcum-PK
+ M,rest,acum,cat=get_pts_Clot_sal(Asal,xda,yda,Az_sal,distcond,distcond,1,R,PK-1,1,i); break
+
+ acum=LsAcum
+
+ if PK > LsAcum:
+ M,rest,acum,cat=get_pts_Recta(xda,yda,z_ini,Az_sal,PK-LsAcum,PK-LsAcum,1,PK-1,1,i)
+
+ return M[0]
+
+
+#### ############## Planta #### ##################
+
+def read_Map(planta,layer):
+
+ # Get Plant Map
+ sal=g.read_command('v.report', map=planta, layer=layer, option='coor',
+ units='me', quiet=True)
+ sal = [d.split('|') for d in sal.splitlines(0)]
+ del sal[0]
+ return sal
+
+def clotoide_Locales(A,R):
+
+ if R==0:
+ L,Tau=0,0
+ else:
+ L=A**2/abs(R)
+ Tau=L/(2*R)
+ xe,ye=aprox_coord(L,Tau)
+ xo=xe-R*sin(Tau)
+ yo=ye+R*cos(Tau)-R
+ return xo,yo,Tau,L,xe,ye
+
+def angulos_Alings(a,b,c,d,e,f):
+
+ # Angulo entre las dos rectas
+ Az_ent= azimut(a,b,c,d)
+ Az_sal=azimut(c,d,e,f)
+
+ if ((a <= c and b <= d) or (a <= c and b >= d)): #1er cuadrante o 2do cuadrante
+ #print "1er o 2do cuadrante"
+ if (Az_ent < Az_sal and Az_ent+pi < Az_sal):
+ w=2*pi-abs(Az_ent-Az_sal)
+ else:
+ w=abs(Az_ent-Az_sal)
+
+ if ((a >= c and b >= d) or (a >= c and b <= d)): #3er cuadrante o 4to cuadrante
+ #print "3er o 4to cuadrante"
+ if (Az_ent > Az_sal and Az_ent-pi > Az_sal):
+ w=2*pi-abs(Az_ent-Az_sal)
+ else:
+ w=abs(Az_ent-Az_sal)
+ #print Az_ent*200/pi,Az_sal*200/pi,w
+ return Az_ent,Az_sal,w
+
+def pto_corte_2_rectas(x1,y1,x2,y2,x3,y3,x4,y4):
+
+ if x2 == x1: m11=(y2-y1)/0.0000001
+ else: m11=(y2-y1)/(x2-x1)
+
+ if x3 == x4: m22=(y4-y3)/0.0000001
+ else: m22=(y4-y3)/(x4-x3)
+
+ x=(m11*x1-m22*x3-y1+y3)/(m11-m22)
+ y=m11*(x-x1)+y1
+
+ return x,y
+
+
+def get_PtosEjePlanta(table_plant):
+
+ LAcum=0
+ x_ini,y_ini = table_plant[0][0],table_plant[0][1]
+ puntos_eje=[]
+ centros = []
+ for i in range(1,len(table_plant)-1,1):
+
+ a,b = table_plant[i-1][0],table_plant[i-1][1]
+ c,d = table_plant[i][0],table_plant[i][1]
+ e,f = table_plant[i+1][0],table_plant[i+1][1]
+
+ Az_ent,Az_sal,w=angulos_Alings(a,b,c,d,e,f)
+
+ R,Ae,As = table_plant[i][5],table_plant[i][6],table_plant[i][7]
+
+ if R == 0:
+ Lrecta=sqrt((c-x_ini)**2+(d-y_ini)**2)
+ LAcum+=Lrecta
+ Lini_e,Lini_s=0,0
+ xadp,yadp = c,d
+ xdap,ydap = c,d
+ puntos_eje.append([[x_ini,y_ini,0,Lrecta,LAcum,Az_ent],
+ [0,Az_ent,c,d,c,d,0,LAcum,0,Lini_e,xadp,yadp],
+ [0,0,0,0,0,0,0,LAcum],
+ [0,Az_sal,c,d,c,d,0,LAcum,0,Lini_s,xdap,ydap]])
+ x_ini,y_ini=c,d
+ continue
+
+ xo_e,yo_e,Tau_e,Le,xe,ye=clotoide_Locales(Ae,R)
+ xo_s,yo_s,Tau_s,Ls,xs,ys=clotoide_Locales(As,R)
+
+ # Centros de las curvas dados por el corte de las rectas paralelas a R+AR
+ if R > 0:
+ g90 = pi/2
+ alpha = abs(w - Tau_e - Tau_s)
+ else:
+ g90 = -pi/2
+ alpha = abs(w + Tau_e + Tau_s)
+
+ ac=a+abs(R+yo_e)*sin(Az_ent+g90)
+ bc=b+abs(R+yo_e)*cos(Az_ent+g90)
+
+ cc1=c+abs(R+yo_e)*sin(Az_ent+g90)
+ dc1=d+abs(R+yo_e)*cos(Az_ent+g90)
+
+ cc2=c+abs(R+yo_s)*sin(Az_sal+g90)
+ dc2=d+abs(R+yo_s)*cos(Az_sal+g90)
+
+ ec=e+abs(R+yo_s)*sin(Az_sal+g90)
+ fc=f+abs(R+yo_s)*cos(Az_sal+g90)
+
+ xc,yc=pto_corte_2_rectas(ac,bc,cc1,dc1,cc2,dc2,ec,fc)
+ centros.append([xc,yc,0])
+
+ Lini_e,Lini_s=0,0
+ if Ae <= 0:
+ xar = xc+abs(R)*sin(Az_ent-g90-Tau_e)
+ yar = yc+abs(R)*cos(Az_ent-g90-Tau_e)
+ xad,yad = xar,yar
+
+ Lrecta = sqrt((xad-x_ini)**2+(yad-y_ini)**2)
+ if Ae < 0:
+ x_ini,y_ini = xad,yad
+ Lrecta = 0
+ Le = 0
+ xadp,yadp = xad,yad
+
+ else:
+ xt1= xc+abs(R+yo_e)*sin(Az_ent-g90)
+ yt1= yc+abs(R+yo_e)*cos(Az_ent-g90)
+
+ xad = xt1+xo_e*sin(Az_ent+pi)
+ yad = yt1+xo_e*cos(Az_ent+pi)
+
+ xar = xc+abs(R)*sin(Az_ent-g90+Tau_e)
+ yar = yc+abs(R)*cos(Az_ent-g90+Tau_e)
+
+ Lrecta = sqrt((xad-x_ini)**2+(yad-y_ini)**2)
+ R1,Ae1,As1 = table_plant[i-1][5],table_plant[i-1][6],table_plant[i-1][7]
+ xadp,yadp = xad,yad
+ if As1 < 0:
+ xo_s1,yo_s1,Tau_s1,Ls1,xs1,ys1=clotoide_Locales(As1,R1)
+ Lini_e = Ls1
+ Ptos_ent,resto,acum,cat=get_pts_Clot_ent(Ae,xad,yad,Az_ent,Lini_e,Lini_e,1,R,0,0,0)
+ xadp,yadp = Ptos_ent[0][0],Ptos_ent[0][1]
+ x_ini,y_ini = xadp,yadp
+ Lrecta = 0
+
+ Dc = abs(R)*alpha
+
+ if As <= 0:
+ xra = xc+abs(R)*sin(Az_sal-g90+Tau_s)
+ yra = yc+abs(R)*cos(Az_sal-g90+Tau_s)
+ xda,yda = xra,yra
+ Ls = 0
+ xdap,ydap = xda,yda
+
+ else:
+ xt1= xc+abs(R+yo_s)*sin(Az_sal-g90)
+ yt1= yc+abs(R+yo_s)*cos(Az_sal-g90)
+
+ xda = xt1+xo_s*sin(Az_sal)
+ yda = yt1+xo_s*cos(Az_sal)
+
+ xra = xc+abs(R)*sin(Az_sal-g90-Tau_s)
+ yra = yc+abs(R)*cos(Az_sal-g90-Tau_s)
+
+ R2,Ae2,As2 = table_plant[i+1][5],table_plant[i+1][6],table_plant[i+1][7]
+ xdap,ydap = xda,yda
+ if Ae2 < 0:
+ xo_e2,yo_e2,Tau_e2,Le2,xe2,ye2=clotoide_Locales(Ae2,R2)
+ Lini_s = Le2
+ Ptos_sal,resto,acum,cat_s=get_pts_Clot_sal(As,xda,yda,Az_sal,Lini_s,Lini_s,1,R,0,0,0)
+ xdap,ydap = Ptos_sal[0][0],Ptos_sal[0][1]
+
+ #Az_ini = Az_ent-g90+Tau_e
+ #Az_fin = Az_sal-g90-Tau_s
+ Az_ini=azimut(xc,yc,xar,yar)
+ Az_fin=azimut(xc,yc,xra,yra)
+
+ puntos_eje.append([[x_ini,y_ini,0,Lrecta,LAcum+Lrecta,Az_ent],
+ [Ae,Az_ent,xad,yad,xar,yar,(Le-Lini_e),LAcum+Lrecta+(Le-Lini_e),R,Lini_e,xadp,yadp],
+ [R,alpha,xc,yc,Az_ini,Az_fin,Dc,LAcum+Lrecta+(Le-Lini_e)+Dc],
+ [As,Az_sal,xda,yda,xra,yra,(Ls-Lini_s),LAcum+Lrecta+(Le-Lini_e)+Dc+(Ls-Lini_s),R,Lini_s,xdap,ydap]])
+
+ x_ini,y_ini = xdap,ydap
+ LAcum = LAcum+Lrecta+(Le-Lini_e)+Dc+(Ls-Lini_s)
+
+ Lrecta = sqrt((e-x_ini)**2+(f-y_ini)**2)
+ puntos_eje.append([[x_ini,y_ini,0,Lrecta,LAcum+Lrecta,Az_sal],
+ [0,Az_sal,x_ini,y_ini,0,0,0,LAcum+Lrecta,0,Lini_e,xadp,yadp],
+ [0,0,0,0,0,0,0,LAcum+Lrecta],
+ [0,Az_sal,x_ini,y_ini,0,0,0,LAcum+Lrecta,0,Lini_s,xdap,ydap]])
+ #for jj in puntos_eje:
+ #for pp in jj: print pp
+ #print ''
+ return puntos_eje
+
+
+def generate_Pts(puntos_eje,a,b,c,interv,intervC):
+
+ x_ini,y_ini,z_ini,Lrecta,LAcum,Az_ent=puntos_eje[0][0]
+
+ resto=0
+ cat=2
+ puntos,seg,puntos_caract,puntos_centros=[],[],[],[]
+ if a == 1: puntos=[[x_ini,y_ini,z_ini,1,0,Az_ent,'Ini',0]]
+ if b == 1: puntos_caract=[[x_ini,y_ini,z_ini,1,0,Az_ent,"'Ini'",0]]
+
+ Ini=[[x_ini,y_ini,z_ini,cat,0,Az_ent,"'Ini'",0]]
+
+ acum=0
+ h=1
+ for i in range(0,len(puntos_eje)-1,1):
+ x_ini,y_ini,z_ini,Lrecta,LrAcum,Az_ent=puntos_eje[i][0]
+ Aent,Az_ent,xad,yad,xar,yar,Lent,LeAcum,R,Lini_e,xadp,yadp=puntos_eje[i][1]
+ R,alpha,xc,yc,Az_ini,Az_fin,Dc,LcAcum=puntos_eje[i][2]
+ Asal,Az_sal,xda,yda,xra,yra,Lsal,LsAcum,R,Lini_s,xdap,ydap=puntos_eje[i][3]
+
+ Ptos_recta,resto,acum,cat_r=get_pts_Recta(x_ini,y_ini,z_ini,Az_ent,interv-resto,Lrecta,interv,acum,cat,i)
+ if R!=0:
+ if Aent > 0:
+ Ptos_ent,resto,acum,cat_e=get_pts_Clot_ent(Aent,xad,yad,Az_ent,interv-resto+Lini_e,Lent,intervC,R,acum,cat_r,i)
+ else:
+ Ptos_ent=[]
+ cat_e=cat_r
+ #print xc,yc,Az_ini-resto/R,Az_fin,interv,R,acum,Dc,cat_e,i
+ Ptos_curva,resto,acum,cat_c=get_pts_Curva(xc,yc,Az_ini-resto/R,Az_fin,intervC,R,acum,Dc,cat_e,i)
+ if Asal > 0:
+ #print Asal,xda,yda,Az_sal,interv-resto+Lini_s,Lsal,interv,R,acum,cat_c,i
+ Ptos_sal,resto,acum,cat_s=get_pts_Clot_sal(Asal,xda,yda,Az_sal,0,Lsal-(intervC-resto+Lini_s),intervC,R,acum,cat_c,i)
+ else:
+ if Lini_s != 0: xda,yda = xdap,ydap
+ Ptos_sal=[]
+ cat_s=cat_c
+ else:
+ Ptos_ent,Ptos_curva,Ptos_sal=[[],[],[]]
+ cat = cat_r
+ cat_r,cat_e,cat_c,cat_s=1,2,3,4
+
+ if a == 1:
+ puntos.extend(Ptos_recta+Ptos_ent+Ptos_curva+Ptos_sal)
+ if c == 1:
+ puntos_centros.append([[xc,yc,0,i,"'C'",i],[0,0,0,h,LrAcum,"'Line'",Lrecta,0,Az_ent*200/pi],
+ [0,0,0,h+1,LeAcum,"'Clot_in'",Lent-Lini_e,Aent,Az_ent*200/pi],
+ [0,0,0,h+2,LcAcum,"'Curve'",abs(Dc),R,Az_ini*200/pi],
+ [0,0,0,h+3,LsAcum,"'Clot_out'",Lsal-Lini_s,Asal,Az_sal*200/pi]])
+ if b == 1:
+ puntos_caract.extend([[xadp,yadp,0,h+1,LrAcum,Az_ent*200/pi,"'Te'",i],
+ [xar,yar,0,h+2,LeAcum,Az_ini*200/pi,"'Tec'",i],
+ [xra,yra,0,h+3,LcAcum,Az_fin*200/pi,"'Tsc'",i],
+ [xdap,ydap,0,h+4,LsAcum,Az_sal*200/pi,"'Ts'",i]])
+ if c == 1:
+ seg.append(Ini+Ptos_recta+[[xadp,yadp,0,cat_r,LrAcum,"'Te'",i]])
+ seg.append([[xadp,yadp,0,cat_r,LrAcum,"'Te'",i]]+Ptos_ent+[[xar,yar,0,cat_e,LeAcum,"'Tec'",i]])
+ seg.append([[xar,yar,0,cat_e,LeAcum,"'Tec'",i]]+Ptos_curva+[[xra,yra,0,cat_c,LcAcum,"'Tsc'",i]])
+ seg.append([[xra,yra,0,cat_c,LcAcum,"'Tsc'",i]]+Ptos_sal+[[xdap,ydap,0,cat_s,LsAcum,"'Ts'",i]])
+ Ini=[[xdap,ydap,0,cat_s,LsAcum,"'Ts'",i]]
+ cat = cat_s
+ h=h+4
+
+ x_ini,y_ini,z_ini,Lrecta,LAcum,Az_sal=puntos_eje[-1][0]
+ Ptos_recta,resto,acum,cat=get_pts_Recta(x_ini,y_ini,z_ini,Az_sal,interv-resto,Lrecta,interv,acum,cat_s,i+1)
+ x_fin,y_fin,z_fin=get_pts_Recta(x_ini,y_ini,z_ini,Az_sal,Lrecta,Lrecta,interv,acum,cat_s,i+1)[0][0][:3]
+
+ if a == 1: puntos.extend(Ptos_recta)
+ if a == 1: puntos.append([x_fin,y_fin,z_fin,cat,LAcum,Az_sal,'End',i+1])
+ if b == 1: puntos_caract.append([x_fin,y_fin,z_fin,h+1,LAcum,Az_sal,"'End'",i+1])
+ if c == 1: puntos_centros.append([[],[x_fin,y_fin,z_fin,h,LAcum,"'Line'",Lrecta,0,Az_sal*200/pi]])
+ if c == 1: seg.append(Ini+Ptos_recta+[[x_fin,y_fin,z_fin,cat,LAcum,"'End'",i+1]])
+
+ #for jj in puntos: print jj
+ return puntos,seg,puntos_caract,puntos_centros
+
+
+#### ############## Alzado #### ##################
+
+def get_Cota(puntos,puntos_eje_alz):
+
+ ins=len(puntos[0])
+ i=0
+ for pto in puntos:
+
+ pk=float(pto[4])
+ for h,pkalz in enumerate(puntos_eje_alz):
+
+ pk_ini,zini,pke,ze,pks,zs,pkv,zv,Kv,L,pe,ps=pkalz
+
+ if Kv == 0 and pk_ini <= pk and pk <= pkv:
+ cota=zini+pe*(pk-pk_ini)
+ pend=type_Pend(pe)
+ aling=h; continue
+ elif pk_ini <= pk and pk < pke:
+ cota=zini+pe*(pk-pk_ini)
+ pend=type_Pend(pe)
+ aling=h; continue
+ elif pke <= pk and pk <= pkv:
+ cota=ze+pe*(pk-pke)-(1/(2*Kv))*(pk-pke)**2
+ pend=type_Acuerdo(pe,ps)
+ aling=h; continue
+ elif pkv < pk and pk <= pks:
+ #cota=ze+ps*(pk-pkv)-(1/(2*Kv))*(pk-pkv)**2
+ cota=ze+pe*(pk-pke)-(1/(2*Kv))*(pk-pke)**2
+ pend=type_Acuerdo(pe,ps)
+ aling=h; continue
+ if pks < pk:
+ cota=puntos_eje_alz[-1][7]
+ pend=puntos_eje_alz[-1][-1]
+ aling=h
+ puntos[i][2]=cota
+ if ins==8: puntos[i].extend([pend,aling])
+ else: puntos[i][8],puntos[i][9]=pend,aling
+ i=i+1
+ #for jj in puntos: print jj
+ return puntos
+
+
+def get_PtosEjeAlzado(alz):
+
+ pk_ini,zini=alz[0][4],alz[0][5]
+ pkv,zv=[pk_ini,zini]
+ puntos_eje_alz=[]
+
+ for i in range(1,len(alz)-1,1):
+
+ pk1,z1=alz[i-1][4],alz[i-1][5]
+ pkv,zv=alz[i][4],alz[i][5]
+ pk2,z2=alz[i+1][4],alz[i+1][5]
+
+ pe=(zv-z1)/(pkv-pk1)
+ if pk2 == pkv: ps=0
+ else: ps=(z2-zv)/(pk2-pkv)
+
+ if alz[i][7]!='':Kv=alz[i][7]
+ else:
+ puntos_eje_alz.append([pk_ini,zini,pkv,zv,pkv,zv,pkv,zv,0,0,pe,ps])
+ pk_ini=pkv; zini=zv
+ continue
+ phi=(pe-ps)
+
+ T=Kv*phi/2
+ L=2*T
+ B=Kv*phi**2/8
+ pke=pkv-T
+ pks=pkv+T
+ ze=zv-pe*T
+ zs=zv+ps*T
+
+ puntos_eje_alz.append([pk_ini,zini,pke,ze,pks,zs,pkv,zv,Kv,L,pe,ps])
+ pk_ini=pks; zini=zs
+
+ pk_fin,zfin=alz[-1][4],alz[-1][5]
+ pk2,z2=[pk_fin,zfin]
+ if pk2 == pkv: ps=0
+ else: ps=(z2-zv)/(pk2-pkv)
+ puntos_eje_alz.append([pk_ini,zini,pk_fin,zfin,pk_fin,zfin,pk_fin,zfin,0,0,ps,0])
+ #for jj in puntos_eje_alz: print jj
+ return puntos_eje_alz
+
+def type_Pend(pe):
+
+ if pe>0: typep="'Up'"
+ else: typep="'Down'"
+ return typep
+
+def type_Acuerdo(pe,ps):
+
+ if ((pe>0 and ps<0) or abs(pe)<abs(ps)): typeac="'Convex'"
+ else: typeac="'Concave'"
+ return typeac
+
+
+def generate_PtsAlz(puntos,puntos_eje,puntos_alz):
+
+ p_seg,p_vert=[],[]
+ ps=0
+ pk_ini,zini,pks,pe=puntos_alz[0][0],puntos_alz[0][1],0,puntos_alz[0][10]
+ x0,y0=get_PlantaXY(float(pk_ini),puntos_eje)[:2]
+ p_caract=[[x0,y0,zini,1,pk_ini,"'Ini'",1]]
+ Ini=[[x0,y0,zini,1,pk_ini,"'Ini'",1]]
+ i,h,k=1,1,1
+
+ for i in range(0,len(puntos_alz)-1,1):
+
+ pk_ini,zini,pke,ze,pks,zs,pkv,zv,Kv,L,pe,ps=puntos_alz[i]
+ if Kv==0: cota=zv
+ else: cota=ze+pe*(pkv-pke)-(1/(2*Kv))*(pkv-pke)**2
+ x1,y1=get_PlantaXY(float(pke),puntos_eje)[:2]
+ x2,y2=get_PlantaXY(float(pkv),puntos_eje)[:2]
+ x3,y3=get_PlantaXY(float(pks),puntos_eje)[:2]
+
+ p_caract.append([x1,y1,ze,k+1,pke,"'Ae'",i+1])
+ p_caract.append([x2,y2,cota,k+2,pkv,"'v'",i+1])
+ p_caract.append([x3,y3,zs,k+3,pks,"'As'",i+1])
+
+ p_vert.append([[x2,y2,zv,i,"'V'",i],
+ [0,0,0,h,pk_ini,type_Pend(pe),pke-pk_ini,pe],
+ [0,0,0,h+1,pke,type_Acuerdo(pe,ps),L/2,Kv],
+ [0,0,0,h+2,pkv,type_Acuerdo(pe,ps),L/2,Kv]])
+
+ p_seg.append(Ini+[p for p in puntos if (p[4]> pk_ini and p[4]< pke)]+
+ [[x1,y1,ze,1,pke,"'Ae'",i+1]])
+ p_seg.append([[x1,y1,ze,1,pke,"'Ae'",i+1]]+
+ [p for p in puntos if (p[4]> pke and p[4]< pkv)]+
+ [[x2,y2,cota,1,pkv,"'V'",i+1]])
+ p_seg.append([[x2,y2,cota,1,pkv,"'V'",i+1]]+
+ [p for p in puntos if (p[4]> pkv and p[4]< pks)]+
+ [[x3,y3,zs,1,pks,"'As'",i+1]])
+ Ini=[[x3,y3,zs,i,pks,"'As'",i+1]]
+
+ h=h+3
+ k=k+3
+
+ pk_fin,zfin=puntos_alz[-1][2],puntos_alz[-1][3]
+
+ x4,y4=get_PlantaXY(float(pk_fin),puntos_eje)[:2]
+
+ p_caract.append([x4,y4,zfin,k+1,pk_fin,"'End'",i+1])
+ p_seg.append(Ini+[p for p in puntos if (p[4]> pks and p[4]< pk_fin)]+
+ [[x4,y4,zfin,1,pk_fin,"'End'",i+1]])
+ p_vert.append([[x4,y4,zfin,h,0,"'V'",i],
+ [x4,x4,zfin,h,pk_fin,type_Pend(pe),pk_fin-pks,ps]])
+ #for jj in p_seg: print jj
+ return p_seg,p_caract,p_vert
+
+
+#### ############## Section #### ##################
+
+def read_ColumnSecc(table_secc,column):
+
+ Sec,Cota,pks,type1=[],[],[],[]
+ for i in range(0,len(table_secc),1):
+
+ if table_secc[i][column]!='':
+
+ if table_secc[i][column].find(';')==-1: # not found
+ Sec.append([table_secc[i][column]])
+ else:
+ Sec.append(table_secc[i][column].split(';'))
+ pks.append(table_secc[i][4])
+ if table_secc[i][column+2] != '':
+ type1.append(table_secc[i][column+2].split(';'))
+ else:
+ type1.append([ 'l' for p in Sec[i]])
+
+ if Sec!=[]:
+ Sec = map(lambda *row: [elem or '0 0' for elem in row], *Sec) # transpuesta
+ for i,line in enumerate(Sec):
+ Sec[i]=[float(p.split(' ')[0]) for p in line]
+ Cota.append([float(p.split(' ')[1]) for p in line])
+ type1 = map(lambda *row: [elem or 'l' for elem in row], *type1) # transpuesta
+
+ return Sec,Cota,pks,type1
+
+
+def get_Desplaz(table_plant,Sec,Cota,pks,type1,Puntos,izq):
+
+ M_ptos,M_ptos_caract=[],[]
+ type2=[]
+ nume=0
+ #incremento de sobreancho a cada desplazado
+ if izq==1:
+ for tt in type1:
+ if "e" in tt: nume=nume+1
+ sobrei=nume+1
+ else: sobrei=0
+
+ for i,line in enumerate(Sec):
+
+ Puntos_Despl=[[]]
+ Pc_izq=[]
+ pkant=0.0
+ npks=pks[:]
+ line2=line[:]
+ type2=type1[i][:]
+
+ for j,dist in enumerate(line2[:-1]):
+
+ if line2[j]!=-1 and line2[j+1]==-1:
+
+ line2[j+1]=line2[j]
+ npks[j+1]=npks[j]
+ type2[j+1]=type2[j]
+ continue
+
+ #incremento de sobreancho a cada desplazado
+ if type2[j]=="e" and izq==1: sobrei=sobrei-1
+ elif type2[j]=="e" and izq==0: sobrei=sobrei+1
+
+ #if npks[j+1]==npks[-1] or line2[j+1]==0: npks[j+1]=npks[j+1]
+
+ ptosDes,Pc=get_PtosDesplazdos(Puntos,table_plant,line2[j],Cota[i][j],line2[j+1],Cota[i][j+1],sobrei,izq,npks[j],npks[j+1],type2[j])
+
+ #Corregimos los pks de los desplazados
+ ptosDes2=[]
+ for h,row in enumerate(ptosDes):
+ if row!=[]:
+ ptosDes2.append(row)
+
+ if ptosDes2 != []:
+ if Puntos_Despl[-1]==[]:
+ Puntos_Despl=[]
+ x_ant=ptosDes2[0][0]
+ y_ant=ptosDes2[0][1]
+ pkant=0
+ elif ptosDes2[-1]!=[]:
+ x_ant=Puntos_Despl[-1][0]
+ y_ant=Puntos_Despl[-1][1]
+ pkant=Puntos_Despl[-1][4]
+
+ x2=ptosDes2[0][0]
+ y2=ptosDes2[0][1]
+
+ l=sqrt((x2-x_ant)**2+(y2-y_ant)**2)
+
+ for h,row in enumerate(ptosDes):
+ if row!=[]:
+ ptosDes[h][4]=row[4]+pkant+l
+
+ Pc_izq.extend(Pc)
+ Puntos_Despl.extend(ptosDes)
+
+ M_ptos.append(Puntos_Despl)
+ if Pc_izq!=[]: M_ptos_caract.append(Pc_izq)
+
+ return M_ptos,M_ptos_caract
+
+
+
+def generate_Desplaz(table_plant,table_secc,Puntos):
+
+ M_ptos_izq,M_ptos_der,M_ptos_caract=[],[],[]
+ SecIzq,CotaIzq,pks,typel=read_ColumnSecc(table_secc,5)
+ SecDer,CotaDer,pks,typer=read_ColumnSecc(table_secc,6)
+
+ M_ptos_izq,M_ptos_caract=get_Desplaz(table_plant,SecIzq,CotaIzq,pks,typel,Puntos,1)
+ M_ptos_der,M_ptos_caract2=get_Desplaz(table_plant,SecDer,CotaDer,pks,typer,Puntos,0)
+
+ M_ptos_caract=M_ptos_caract+M_ptos_caract2
+
+ M_ptos_caract = [m for p in M_ptos_caract for m in p]
+
+ return M_ptos_izq,M_ptos_der,M_ptos_caract
+
+
+
+def clotoide_GetA(radio,yo_ent,sobreancho,izq):
+
+ h=1*pi/200
+ Tau=0.0001*pi/200
+ xo,yo=aprox_coord2(abs(radio),Tau)
+ if (izq and radio>0) or (izq==0 and radio<0): sobreancho=sobreancho
+ elif (izq and radio<0) or (izq==0 and radio>0): sobreancho=sobreancho*(-1)
+
+ #Comprobar que el sobreancho sea mayor que el retranqueo
+ while abs(yo-(yo_ent-sobreancho))>0.000001:
+
+ xo,yo=aprox_coord2(abs(radio),Tau)
+ if yo>yo_ent-sobreancho:
+ Tau=Tau-h
+ h=h/10
+ else:
+ Tau=Tau+h
+
+ L=Tau*2*abs(radio)
+ A=sqrt(L*abs(radio))
+ return A
+
+
+def get_Paralles(table_plant,dist,sobreancho,izq):
+
+ a,b,z1=float(table_plant[0][0]),float(table_plant[0][1]),float(table_plant[0][2])
+ c,d,z2=float(table_plant[1][0]),float(table_plant[1][1]),float(table_plant[1][2])
+ r=sqrt((c-a)**2+(d-b)**2)
+ v=[(d-b)/r,-(c-a)/r]
+ if izq: v[0]=v[0]*-1; v[1]=v[1]*-1
+ salida=[]
+ salida.append([v[0]*dist+a,v[1]*dist+b,table_plant[0][2],table_plant[0][3],0,0.0,0.0,0.0,0.0])
+ Lrecta=0
+ for ii in range(1,len(table_plant)-1,1):
+
+ sobreanchoTable=float(table_plant[ii][8])
+ a,b,z1=float(table_plant[ii-1][0]),float(table_plant[ii-1][1]),float(table_plant[ii-1][2])
+ c,d,z2=float(table_plant[ii][0]),float(table_plant[ii][1]),float(table_plant[ii][2])
+ e,f,z3=float(table_plant[ii+1][0]),float(table_plant[ii+1][1]),float(table_plant[ii+1][2])
+
+ r=sqrt((c-a)**2+(d-b)**2)
+ v=[(d-b)/r,-(c-a)/r]
+ if izq: v[0]=v[0]*-1; v[1]=v[1]*-1
+ a1=v[0]*dist+a
+ b1=v[1]*dist+b
+ c1=v[0]*dist+c
+ d1=v[1]*dist+d
+
+ s=sqrt((e-c)**2+(f-d)**2)
+ w=[(f-d)/s,-(e-c)/s]
+ if izq: w[0]=w[0]*-1; w[1]=w[1]*-1
+ e1=w[0]*dist+e
+ f1=w[1]*dist+f
+ c2=w[0]*dist+c
+ d2=w[1]*dist+d
+
+ c3,d3=pto_corte_2_rectas(a1,b1,c1,d1,e1,f1,c2,d2)
+
+ Lrecta+=sqrt((c3-a1)**2+(d3-b1)**2)
+ R=float(table_plant[ii][5])
+ if R==0:
+ salida.append([c3,d3,table_plant[ii][2],table_plant[ii][3],Lrecta,0.0,0.0,0.0,0.0])
+ continue
+
+ Aent,Asal=table_plant[ii][6],table_plant[ii][7]
+
+ # Clotoide de entrada en locales
+ xo_ent,yo_ent,Tau_ent,Lent,xe,ye=clotoide_Locales(Aent,abs(R))
+ # Clotoide de salida en locales
+ xo_sal,yo_sal,Tau_sal,Lsal,xs,ys=clotoide_Locales(Asal,abs(R))
+
+ if izq!=0 and R<0: radio=R+dist+(sobreancho*sobreanchoTable)
+ elif izq!=0 and R>0: radio=R+dist-(sobreancho*sobreanchoTable)
+ elif izq==0 and R>0: radio=R-dist+(sobreancho*sobreanchoTable)
+ elif izq==0 and R<0: radio=R-dist-(sobreancho*sobreanchoTable)
+
+ # Parametro A de entrada
+ if Aent==0: Aent1=0
+ else: Aent1=clotoide_GetA(radio,yo_ent,sobreancho*sobreanchoTable,izq)
+ # Parametro A de salida
+ if Asal==0: Asal1=0
+ else: Asal1=clotoide_GetA(radio,yo_sal,sobreancho*sobreanchoTable,izq)
+
+ salida.append([c3,d3,table_plant[ii][2],table_plant[ii][3],Lrecta,radio,Aent1,Asal1,0])
+
+ a=float(table_plant[-2][0]); b=float(table_plant[-2][1]); z1=float(table_plant[-2][2])
+ c=float(table_plant[-1][0]); d=float(table_plant[-1][1]); z2=float(table_plant[-1][2])
+ r=sqrt((c-a)**2+(d-b)**2)
+ v=[(d-b)/r,-(c-a)/r]
+ if izq: v[0]=v[0]*-1; v[1]=v[1]*-1
+ Lrecta+=sqrt((c-a)**2+(d-b)**2)
+ salida.append([v[0]*dist+c,v[1]*dist+d,table_plant[-1][2],table_plant[-1][3],Lrecta,0.0,0.0,0.0,0.0])
+ #print salida
+ return salida
+
+#### ############## Transversales #### ##################
+
+#Pto de corte entre la perpendicular a un pto de una recta y una recta
+def get_PtoCorte_Perpend_Recta(Az,x,y,xref_d,yref_d,Az_d,Lrecta_d):
+
+ h=Lrecta_d/2
+ Li=-0.00001
+ eq=-0.001
+ while abs(eq)>0.00001 and Li<=Lrecta_d:
+
+ eq_ant=eq
+ x1=xref_d+Li*sin(Az_d)
+ y1=yref_d+Li*cos(Az_d)
+ eq=y1-(y-tan((Az))*(x1-x))
+ if (eq_ant>0 and eq<0) or (eq_ant<0 and eq>0):
+ Li=Li-h
+ h=h/2
+ else:
+ Li=Li+h
+ if abs(eq)>0.001: x1,y1=[0,0]
+ return x1,y1,Li
+
+#Pto de corte entre la perpendicular a un pto de un circulo y un circulo
+def get_PtoCorte_Perpend_Circulo(Az,x,y,xref_d,yref_d,Az_ini_d,Az_fin_d,Dc_d,R_d):
+
+ h=Dc_d/2
+ Li=-0.000001
+ eq=-0.001
+ cond=True
+ if R_d>0:
+ if Az_ini_d>Az_fin_d: Az_ini_d=Az_ini_d-2*pi
+ while abs(eq)>0.0001 and cond:
+ eq_ant=eq
+ x1=xref_d+R_d*sin(Az_ini_d+abs(Li/R_d))
+ y1=yref_d+R_d*cos(Az_ini_d+abs(Li/R_d))
+ eq=y1-(y-tan((Az))*(x1-x))
+ if (eq_ant>0 and eq<0) or (eq_ant<0 and eq>0):
+ Li=Li-h
+ h=h/2
+ else:
+ Li=Li+h
+ cond=(Az_ini_d+(Li/R_d)<=Az_fin_d)
+ elif R_d<0:
+ if Az_ini_d<Az_fin_d: Az_fin_d=Az_fin_d-2*pi
+ while abs(eq)>0.0001 and cond:
+ eq_ant=eq
+ #print (Az_ini_d-abs(Li/R_d))*200/pi
+ x1=xref_d-R_d*sin(Az_ini_d-abs(Li/R_d))
+ y1=yref_d-R_d*cos(Az_ini_d-abs(Li/R_d))
+ eq=y1-(y-tan((Az))*(x1-x))
+ if (eq_ant>0 and eq<0) or (eq_ant<0 and eq>0):
+ Li=Li-h
+ h=h/2
+ else:
+ Li=Li+h
+ cond=(Az_ini_d-abs(Li/R_d)>=Az_fin_d)
+ if abs(eq)>0.001: x1,y1=[0,0]
+ return x1,y1,Li
+
+#Pto de corte entre la perpendicular de una clotoide y la clotoide desplazada
+def get_PtoCorte_Perpend_Clot(Az,x,y,xref_d,yref_d,Az_d,A_d,L_d,R_d,tipo):
+
+ h=L_d/2
+ Li=-0.000001
+ eq=-0.001
+ while abs(eq)>0.0001 and Li<=L_d:
+ eq_ant=eq
+ Ri=A_d**2/Li
+ Tau_i=Li/(2*Ri)
+ xo,yo=aprox_coord(Li,Tau_i)
+ if tipo=="ent":
+ if R_d>0:
+ x1=xref_d-xo*sin(-Az_d)+yo*cos(-Az_d)
+ y1=yref_d+xo*cos(-Az_d)+yo*sin(-Az_d)
+ elif R_d<0:
+ x1=xref_d+xo*sin(Az_d)-yo*cos(Az_d)
+ y1=yref_d+xo*cos(Az_d)+yo*sin(Az_d)
+ elif tipo=="sal":
+ if R_d>0:
+ x1=xref_d-xo*sin(Az_d)+yo*cos(Az_d)
+ y1=yref_d-xo*cos(Az_d)-yo*sin(Az_d)
+ elif R_d<0:
+ x1=xref_d+xo*sin(-Az_d)-yo*cos(-Az_d)
+ y1=yref_d-xo*cos(-Az_d)-yo*sin(-Az_d)
+ eq=y1-(y-tan((Az))*(x1-x))
+ if (eq_ant>0 and eq<0) or (eq_ant<0 and eq>0):
+ Li=Li-h
+ h=h/2
+ else:
+ Li=Li+h
+ if abs(eq)>0.001: x1,y1=[0,0]
+ return x1,y1,Li
+
+
+# busca desde la alineacion anterior a la posterior
+def get_ptoDesplaz(x,y,PK,Az,tipo,align,puntos_DesplEje):
+
+ x1,y1,PK2=0,0,0
+ x_ini_d,y_ini_d,z_ini_d,Lrecta_d,LrAcum_d,Az_ent_d=puntos_DesplEje[align][0]
+ Aent_d,Az_ent_d,xad_d,yad_d,xar_d,yar_d,Lent_d,LeAcum_d=puntos_DesplEje[align][1][0:8]
+ R_d,alpha_d,xc_d,yc_d,Az_ini_d,Az_fin_d,Dc_d,LcAcum_d=puntos_DesplEje[align][2][0:8]
+ Asal_d,Az_sal_d,xda_d,yda_d,xra_d,yra_d,Lsal_d,LsAcum_d=puntos_DesplEje[align][3][0:8]
+
+ if tipo=='Line' or tipo=='Ini' or tipo=='End':
+ if align < len(puntos_DesplEje)-1:
+ x1,y1,Li=get_PtoCorte_Perpend_Recta(Az,x,y,x_ini_d,y_ini_d,Az_ent_d,Lrecta_d)
+ else:
+ x1,y1,Li=get_PtoCorte_Perpend_Recta(Az,x,y,xda_d,yda_d,Az_sal_d,Lrecta_d)
+ PK2=LrAcum_d-Lrecta_d+Li
+ if x1==0 and y1==0 and R_d!=0:
+ x1,y1,Li=get_PtoCorte_Perpend_Clot(Az,x,y,xad_d,yad_d,Az_ent_d,Aent_d,Lent_d,R_d,"ent")
+ PK2=LeAcum_d-Lent_d+Li
+ if x1==0 and y1==0 and R_d!=0:
+ R_d=puntos_DesplEje[align-1][2][0]
+ Asal_d,Az_sal_d,xda_d,yda_d,xra_d,yra_d,Lsal_d,LsAcum_d=puntos_DesplEje[align-1][3][0:8]
+ x1,y1,Li=get_PtoCorte_Perpend_Clot(Az,x,y,xda_d,yda_d,Az_sal_d,Asal_d,Lsal_d,R_d,"sal")
+ PK2=LsAcum_d-Lsal_d+Li
+
+ elif tipo=='Clot_in':
+ x1,y1,Li=get_PtoCorte_Perpend_Clot(Az,x,y,xad_d,yad_d,Az_ent_d,Aent_d,Lent_d,R_d,"ent")
+ PK2=LeAcum_d-Lent_d+Li
+ if x1==0 and y1==0:
+ x1,y1,Li=get_PtoCorte_Perpend_Circulo(Az,x,y,xc_d,yc_d,Az_ini_d,Az_fin_d,Dc_d,R_d)
+ PK2=LcAcum_d-Dc_d+abs(Li)
+ if x1==0 and y1==0:
+ x1,y1,Li=get_PtoCorte_Perpend_Recta(Az,x,y,x_ini_d,y_ini_d,Az_ent_d,Lrecta_d)
+ PK2=LrAcum_d-Lrecta_d+Li
+
+ elif tipo=='Curve':
+ x1,y1,Li=get_PtoCorte_Perpend_Circulo(Az,x,y,xc_d,yc_d,Az_ini_d,Az_fin_d,Dc_d,R_d)
+ PK2=LcAcum_d-Dc_d+abs(Li)
+ if x1==0 and y1==0:
+ x1,y1,Li=get_PtoCorte_Perpend_Clot(Az,x,y,xad_d,yad_d,Az_ent_d,Aent_d,Lent_d,R_d,"ent")
+ PK2=LeAcum_d-Lent_d+Li
+ if x1==0 and y1==0:
+ x1,y1,Li=get_PtoCorte_Perpend_Clot(Az,x,y,xda_d,yda_d,Az_sal_d,Asal_d,Lsal_d,R_d,"sal")
+ PK2=LsAcum_d-Lsal_d+Li
+
+ elif tipo=='Clot_out':
+ x1,y1,Li=get_PtoCorte_Perpend_Clot(Az,x,y,xda_d,yda_d,Az_sal_d,Asal_d,Lsal_d,R_d,"sal")
+ PK2=LsAcum_d-Lsal_d+Li
+ if x1==0 and y1==0:
+ x1,y1,Li=get_PtoCorte_Perpend_Circulo(Az,x,y,xc_d,yc_d,Az_ini_d,Az_fin_d,Dc_d,R_d)
+ PK2=LsAcum_d-Dc_d+abs(Li)
+ if x1==0 and y1==0:
+ if align < len(puntos_DesplEje)-2:
+ x_ini_d,y_ini_d,z_ini_d,Lrecta_d,LrAcum_d,Az_ent_d=puntos_DesplEje[align+1][0]
+ Az_ent_d=puntos_DesplEje[align+1][1][1]
+ x1,y1,Li=get_PtoCorte_Perpend_Recta(Az,x,y,x_ini_d,y_ini_d,Az_ent_d,Lrecta_d)
+ else:
+ Lrecta_d=puntos_DesplEje[align+1][0][3]
+ Asal_d,Az_sal_d,xda_d,yda_d,xra_d,yra_d,Lsal_d,LsAcum_d=puntos_DesplEje[align+1][3][0:8]
+ x1,y1,Li=get_PtoCorte_Perpend_Recta(Az,x,y,xda_d,yda_d,Az_sal_d,Lrecta_d)
+ PK2=LsAcum_d-Lrecta_d+Li
+
+ if x1==0 and y1==0: print 'Point not found ',PK,PK2,tipo
+
+ return x1,y1,PK2
+
+
+
+def get_PtosDesplazdos(puntos,table_plant,dist,cota,dist2,cota2,sobrea,izq,pkini,pkfin,type1):
+
+ lin,Pc=[],[]
+ ptos=[]
+
+ for pp in puntos:
+ if pkini <= pp[4] and pp[4] <= pkfin:
+ ptos.append(pp)
+
+ if round(pkfin - puntos[-1][4],6) == 0:
+ ptos.append(puntos[-1])
+
+ if dist!=0 and dist2==0: dist=0
+ if dist==dist2 and dist!=0 and type1=='e':
+
+ desplaz=get_Paralles(table_plant,dist,sobrea,izq)
+ puntos_DesplEje=get_PtosEjePlanta(desplaz)
+
+ PK2_ini=puntos[int(pkini)][4]
+
+ for pt in ptos:
+
+ x,y,z,cat,PK,Az,tipo,align=pt[0:8]
+ #print x,y,z,cat,PK,Az,tipo,align
+ cota0=cota+((pt[4]-pkini)*(cota2-cota))/(pkfin-pkini)
+ x1,y1,PK2=get_ptoDesplaz(x,y,PK,Az,tipo,align,puntos_DesplEje)
+ if x1!=0 and y1!=0:
+ lin.append([x1,y1,z+cota0,cat,PK2-PK2_ini,Az,'d'+str(dist),PK])
+ else:
+ lin.append([])
+
+ a,b,Pc,c=generate_Pts(puntos_DesplEje,0,1,0,1,1)
+ Pc=[row for row in Pc if row[4]>pkini and row[4]<pkfin]
+
+ elif dist!=0:
+
+ l_ant=0
+ for i,p in enumerate(puntos):
+ if p[4]==pkini:
+ ii=p[3]
+ break
+
+ x,y,z,cat,PK,Az,tipo,align=puntos[ii][0:8]
+
+ #if izq==1: Az1=Az-pi/2
+ #else: Az1=Az+pi/2
+ #x_ant=x+(dist)*sin(Az1)
+ #y_ant=y+(dist)*cos(Az1)
+ x_ant=0
+ y_ant=0
+ #lin.append([x_ant,y_ant,z+cota,cat,0,Az,'d'+str(dist),PK])
+
+ for pt in ptos:
+
+ x,y,z,cat,PK,Az,tipo,align=pt[0:8]
+
+ if type1 == 'l' or type1 == '':
+ distx=dist+((pt[4]-pkini)*(dist2-dist))/(pkfin-pkini)
+
+ elif type1 == 'c':
+ if dist2 < dist:
+ distx=(dist)-sqrt((dist2-dist)**2*(1-cos(acos((pkfin-pt[4])/(pkfin-pkini)))**2))
+ else:
+ distx=(dist2)-sqrt((dist2-dist)**2*(1-cos(acos((pt[4]-pkini)/(pkfin-pkini)))**2))
+
+ elif re.search(r'^r', type1):
+ rad,center=type1[1:].split(",")
+ rad=float(rad)
+ center=float(center)
+
+ if rad > 0:
+ if (rad**2-((pt[4]-pkini)-center)**2) > 0:
+ distx=(rad+dist)-sqrt(rad**2-((pt[4]-pkini)-center)**2)
+ elif (pt[4]-pkini)==0:
+ distx=(rad+dist)
+ else:
+ distx=dist2
+ else:
+ if (rad**2-((pt[4]-pkini)-center)**2) > 0:
+ distx=(rad+dist)+sqrt(rad**2-((pt[4]-pkini)-center)**2)
+ #distx=(dist)-sqrt(rad**2-((pt[4]-pkini)+rad)**2)
+ else:
+ distx=(rad+dist)
+
+ cota0=cota+((pt[4]-pkini)*(cota2-cota))/(pkfin-pkini)
+
+ if izq==1: Az1=Az-pi/2
+ else: Az1=Az+pi/2
+
+ x2=x+(distx)*sin(Az1)
+ y2=y+(distx)*cos(Az1)
+
+ if x_ant == 0: x_ant = x2
+ if y_ant == 0: y_ant = y2
+
+ l=sqrt((x2-x_ant)**2+(y2-y_ant)**2)
+ PK2=l+l_ant
+ lin.append([x2,y2,z+cota0,cat,PK2,Az,'d'+str(round(distx,6)),PK])
+
+ x_ant=x2
+ y_ant=y2
+ l_ant=PK2
+
+ elif dist==0:
+
+ for pt in ptos:
+ lin.append([])
+
+ return lin,Pc
+
+
+def get_Trans(puntos,npk,mpk,m,perpizq,perpder,perpizq2,perpder2,discr,pkini,pkfin):
+
+ trans,pklist=[],[]
+
+ ptos=[]
+ for pp in puntos:
+ if pkini <= pp[4] and pp[4] < pkfin and pp[4]%npk==0:
+ ptos.append(pp)
+ if pkfin == puntos[-1][4]: ptos.append(puntos[-1])
+
+ for pt in ptos:
+
+ if mpk!=0 and (pt[4]%mpk==0):
+ h1=m*perpizq
+ h2=m*perpizq
+ ini_izq,ini_der=perpizq+h1,perpder+h2
+ elif mpk!=0 :
+ h1,h2=0,0
+ ini_izq,ini_der=perpizq,perpder
+ elif discr==0:
+ h1=((pt[4]-pkini)*(perpizq2-perpizq))/(pkfin-pkini)
+ h2=((pt[4]-pkini)*(perpder2-perpder))/(pkfin-pkini)
+ ini_izq,ini_der=perpizq+h1,perpder+h2
+ else:
+ h1=((pt[4]-pkini)*(perpizq2-perpizq))/(pkfin-pkini)
+ h2=((pt[4]-pkini)*(perpder2-perpder))/(pkfin-pkini)
+ ini_izq,ini_der=discr,discr
+
+ x,y,z,cat,PK,Az,tipo,align=pt[0:8]
+
+ Ptos_recta_izq,resto,acum,cat=get_pts_Recta(x,y,z,Az-pi/2,ini_izq,perpizq+h1,discr,0,cat,align)
+ Ptos_recta_der,resto,acum,cat=get_pts_Recta(x,y,z,Az+pi/2,ini_der,perpder+h2,discr,0,cat-1,align)
+ trans.append(Ptos_recta_izq[::-1]+[pt]+Ptos_recta_der)
+ pklist.append(pt)
+ #for jj in trans: print jj
+ return trans,pklist
+
+
+def generate_Transv(puntos,table_transv): #Para modificarlos deben tener cota cero
+
+ tra,trans_pklist=[],[]
+ #table_transv.append([0,0,0,0,table_transv[-1][4]+1,0,0,1])
+ for j,line in enumerate(table_transv[:-1]):
+
+ pkini,pkfin=line[4],table_transv[j+1][4]
+ if line[7]==0.0: continue
+
+ trans,pklist=get_Trans(puntos,int(line[7]),0,0,line[5],line[6],table_transv[j+1][5],table_transv[j+1][6],0,pkini,pkfin)
+ tra.extend(trans)
+ trans_pklist.extend(pklist)
+
+ return tra,trans_pklist
+
+
+#def get_Terrain(puntos,tinmap):
+
+ ##
+ #map_info = pointer(Map_info())
+ #Vect_set_open_level(2)
+ #Vect_open_old(map_info, tinmap, MapSet)
+ #z = c_double()
+ #puntos_terr=[]
+ #for i,punt in enumerate(puntos):
+
+ #Vect_tin_get_z(map_info, float(punt[0]), float(punt[1]), byref(z), None, None)
+
+ #puntos_terr.append([punt[0], punt[1],z.value, punt[3], punt[4]])
+ #Vect_close(map_info)
+ #return puntos_terr
+
+
+###############################################################################
+
+def p_Type(data_type):
+ ptype = POINTER(c_int)
+ if data_type == 0:
+ ptype = POINTER(c_int)
+ elif data_type == 1:
+ ptype = POINTER(c_float)
+ elif data_type == 2:
+ ptype = POINTER(c_double)
+ return ptype
+
+
+
+def get_Terrain(demmap):
+
+ # Get Dem map and return array dem
+ global xref, yref, xres, yres
+ fd = Rast_open_old(demmap, MapSet)
+ data_type = Rast_get_map_type(fd)
+
+ #data_type = G_raster_map_type(demmap, MapSet)
+ ptype = p_Type(data_type)
+ #fd = G_open_cell_old(demmap, MapSet)
+ rast = Rast_allocate_buf(data_type)
+ rast = cast(c_void_p(rast), ptype)
+ window = pointer(Cell_head())
+ G_get_window(window)
+ nrows = window.contents.rows
+ ncols = window.contents.cols
+ xref = window.contents.west
+ yref = window.contents.north
+ xres = window.contents.ew_res
+ yres = window.contents.ns_res
+
+ dem=[[0 for c in range(ncols)] for r in range(nrows)]
+ for i in range(nrows):
+ Rast_get_row(fd, rast, i, data_type)
+ for j in range(ncols):
+ dem[i][j]= rast[j]
+
+ Rast_close(fd)
+ return dem
+
+def drape_Points(puntos,dem):
+
+ salida=[]
+ for i,punt in enumerate(puntos):
+ pto_col = int((punt[0] - xref)/xres)
+ pto_row = int((yref - punt[1])/yres)
+ salida.append(punt[:2]+[dem[pto_row][pto_col]]+punt[3:])
+ return salida
+
+
+def drape_LinesPoints(lines,dem):
+
+ salida=[]
+ for i,line in enumerate(lines):
+ ptos=[]
+ for j,pto in enumerate(line):
+ pto_col = int((pto[0] - xref)/xres)
+ pto_row = int((yref - pto[1])/yres)
+ ptos.append(lines[i][j][:2]+[dem[pto_row][pto_col]]+lines[i][j][3:])
+ salida.append(ptos)
+ return salida
+
+def drape_Points2(puntos,dem):
+
+ elev = raster.RasterRow('elev1')
+ elev.open('r')
+ region=Region()
+ salida=[]
+ for i,punt in enumerate(puntos):
+ pto_col = int((punt[0] - region.west)/region.ewres)
+ pto_row = int((region.north - punt[1])/region.nsres)
+ salida.append(punt[:2]+[elev[pto_row][pto_col]]+punt[3:])
+ elev.close()
+ return salida
+
+
+def drape_LinesPoints2(lines,dem):
+
+ elev = raster.RasterRow('elev1')
+ elev.open('r')
+ region=Region()
+ salida=[]
+ for i,line in enumerate(lines):
+ ptos=[]
+ for j,pto in enumerate(line):
+ pto_col = int((pto[0] - region.west)/region.ewres)
+ pto_row = int((region.north - pto[1])/region.nsres)
+ ptos.append(lines[i][j][:2]+[elev[pto_row][pto_col]]+lines[i][j][3:])
+
+ salida.append(ptos)
+ elev.close()
+ return salida
+
+
+def get_Taludes(puntos,puntos_Despl,des,ter,des2,ter2,dem,pkini,pkfin):
+
+ puntos_talud=[]
+ ptos=[]
+ for pp in puntos:
+ if pkini <= pp[4] and pp[4] < pkfin:
+ ptos.append(pp)
+ if pkfin == puntos[-1][4]: ptos.append(puntos[-1])
+ pta=[]
+ for pt in ptos:
+
+ for ptoD in puntos_Despl:
+ if ptoD!=[] and pt[4] == ptoD[-1]:
+ pta = ptoD
+ continue
+
+ if pta!=[] and des!=0 and ter!=0:
+ a=pt[0]; b=pt[1]; z=pt[2]
+ c=pta[0]; d=pta[1]; z1=pta[2]
+
+ if c!=0 and d!=0:
+ Az=azimut(a,b,c,d)
+ pto_col = int((c - xref)/xres)
+ pto_row = int((yref - d)/yres)
+ zt=dem[pto_row][pto_col]
+ h=1.0
+ Li=0.0
+ zini=z1
+ Li_ant=Li
+ p=0
+ if z1>zt: #Terraplen
+ tipo="Terraplen"
+ talud=ter+((pt[4]-pkini)*(ter2-ter))/(pkfin-pkini)
+ while abs(z1-zt)>0.001:
+ x1=c+Li*sin(Az)
+ y1=d+Li*cos(Az)
+ pto_col = int((x1 - xref)/xres)
+ pto_row = int((yref - y1)/yres)
+ zt=dem[pto_row][pto_col]
+ if p and zt!=zt_ant: zt=zt_ant
+ z1=zini-Li*talud
+ if z1<zt:
+ p=1
+ Li=Li_ant
+ h=h/2
+ else:
+ Li_ant=Li
+ zt_ant=zt
+ Li=Li+h
+
+ elif z1<zt: #Desmonte
+ tipo="Desmonte"
+ talud=des+((pt[4]-pkini)*(des2-des))/float(pkfin-pkini)
+ while abs(zt-z1)>0.001:
+ x1=c+Li*sin(Az)
+ y1=d+Li*cos(Az)
+ pto_col = int((x1 - xref)/xres)
+ pto_row = int((yref - y1)/yres)
+ zt=dem[pto_row][pto_col]
+ if p and zt!=zt_ant: zt=zt_ant #Por si cambia zt en el ultimo instante
+ z1=zini+Li*talud
+ if z1>zt:
+ p=1
+ Li=Li_ant
+ h=h/2
+ else:
+ Li_ant=Li
+ zt_ant=zt
+ Li=Li+h
+
+ puntos_talud.append([x1,y1,z1]+pt[3:]+[tipo])
+ else:
+ puntos_talud.append([])
+
+ return puntos_talud
+
+
+def generate_Taludes(puntos,puntos_Despl_izq,puntos_Despl_der,table_secc,dem):
+
+ tal_izq,tal_der=[],[]
+
+ Despl_izq0,Despl_der1=[],[]
+ if puntos_Despl_izq != [] and puntos_Despl_der != []:
+
+ for j,line in enumerate(table_secc[:-1]):
+
+ pkini,pkfin=line[4],table_secc[j+1][4]
+
+ #if table_secc[j+1]==table_secc[-1]: pkfin=pkfin+1
+
+ taludes=get_Taludes(puntos,puntos_Despl_izq[0],float(line[9]),float(line[11]),
+ float(table_secc[j+1][9]),float(table_secc[j+1][11]),dem,pkini,pkfin)
+ tal_izq.extend(taludes)
+ taludes=get_Taludes(puntos,puntos_Despl_der[-1],float(line[10]),float(line[12]),
+ float(table_secc[j+1][10]),float(table_secc[j+1][12]),dem,pkini,pkfin)
+ tal_der.extend(taludes)
+
+ return tal_izq,tal_der
+
+
+def generate_TaludesAreas(Puntos_Talud_izq,Desplazados_izq,Desplazados_der,Puntos_Talud_der):
+
+ # despl_izq[0] [line1,[],line2,...] --> [[ line1,linet1 ],[ line12,linet2 ],...]
+ # taludes [linet1,[],linet2,...]
+ lines=[]
+ if Puntos_Talud_izq != [] or Puntos_Talud_der != []:
+
+ split_Talud_izq=[list(group) for k, group in groupby(Puntos_Talud_izq, lambda x: x == []) if not k] # split list
+ split_Despl_izq=[list(group) for k, group in groupby(Desplazados_izq[0], lambda x: x == []) if not k] # split list
+ for i,line in enumerate(split_Talud_izq):
+
+ lines.append(split_Talud_izq[i]+split_Despl_izq[i][::-1] )
+
+ split_Talud_der=[list(group) for k, group in groupby(Puntos_Talud_der, lambda x: x == []) if not k] # split list
+ split_Despl_der=[list(group) for k, group in groupby(Desplazados_der[-1], lambda x: x == []) if not k] # split list
+ for i,line in enumerate(split_Despl_der):
+
+ lines.append(split_Despl_der[i]+split_Talud_der[i][::-1] )
+
+ return lines
+
+
+def rellenar_linea(puntos,talud,despl):
+
+ # talud [pto,[],pto,[],...] --> [pto,pto2,pto,pto1,...]
+ # despl [pto1,[],[],pto1,...]
+ # puntos [pto2,pto2,pto2,pto2,...]
+ new=despl+[puntos]
+ for i,pto in enumerate(talud):
+ if pto == []:
+ h=0
+ while new[h][i] == []: h=h+1
+ talud[i] = new[h][i]
+
+ return talud
+
+
+def generate_ContornoAreas(puntos,talud_izq,despl_izq,despl_der,talud_der):
+
+ talud_izq_rell=rellenar_linea(puntos,talud_izq,despl_izq)
+ talud_der_rell=rellenar_linea(puntos,talud_der,despl_der)
+
+ uniq = []
+ for i in talud_izq_rell+talud_der_rell[::-1]:
+ if i not in uniq:
+ uniq.append(i)
+ return uniq
+
+
+##############################################################################
+
+def discr_Lines(lines,discr):
+
+ disLines=[]
+ for pto in lines:
+ resto,acum=0,0
+ cat_r=1
+ disptos=[]
+ for i in range(len(pto[:-1])):
+
+ disptos.append(pto[i])
+ Az=azimut(pto[i][0],pto[i][1],pto[i+1][0],pto[i+1][1])
+ Lrecta=sqrt((pto[i+1][0]-pto[i][0])**2+(pto[i+1][1]-pto[i][1])**2)
+ Ptos_recta,resto,acum,cat_r=get_pts_Recta(pto[i][0],pto[i][1],0,Az,discr-resto,Lrecta,discr,acum,cat_r,i)
+ disptos.extend(Ptos_recta)
+ disptos.append(pto[i+1])
+ disLines.append(disptos)
+ return disLines
+
+
+def get_TransDespl(Trans,desplaz_izq,desplaz_der):
+ # Añade los ptos de corte de los desplazados con los trasversales a los mismos.
+ # Cota cero para poder mover los ptos del transv
+ TransDespl=[]
+ for i,line in enumerate(Trans):
+
+ # Cuidado con la asignacion de listas que se pasan por referencia/mutables
+ tmpline=[line[0][:2]+[0]+line[0][3:]+[0]]
+ for line_desp in desplaz_izq:
+ for i_izq,pto in enumerate(line_desp):
+ if pto!=[] and line[1][4] == pto[-1]:
+ dist_izq=sqrt((line[1][0]-pto[0])**2+(line[1][1]-pto[1])**2)
+ tmpline.append(pto[:2]+[0]+pto[3:]+[round(dist_izq,6)]) # Cota cero para poder mover los ptos del transv
+ continue
+
+ tmpline.append(line[1][:2]+[0]+line[1][3:]+[0])
+
+ for line_desp in desplaz_der:
+ for i_der,pto in enumerate(line_desp):
+ if pto!=[] and line[1][4] == pto[-1]:
+ dist_der=sqrt((line[1][0]-pto[0])**2+(line[1][1]-pto[1])**2)
+ tmpline.append(pto[:2]+[0]+pto[3:]+[round(dist_der,6)])
+ continue
+ tmpline.append(line[2][:2]+[0]+line[2][3:]+[0])
+
+ TransDespl.append(tmpline)
+
+ return TransDespl
+
+
+def get_SeccTerr(Trans,Trans_Pklist,Desplaz_izq,Desplaz_der,Puntos_Talud_izq,Puntos_Talud_der):
+
+ talud_izq=[]
+ for ptoss in Trans_Pklist:
+ esta=0
+ for ptot in Puntos_Talud_izq:
+ if ptot!=[] and ptoss[4] == ptot[4]:
+ talud_izq.append(ptot)
+ esta=1
+ continue
+ if esta==0:
+ talud_izq.append([])
+
+ talud_der=[]
+ for ptoss in Trans_Pklist:
+ esta=0
+ for ptot in Puntos_Talud_der:
+ if ptot!=[] and ptoss[4] == ptot[4]:
+ talud_der.append(ptot)
+ esta=1
+ continue
+ if esta==0:
+ talud_der.append([])
+
+ secc_izq=[talud_izq]
+ for desp in Desplaz_izq:
+ desp_izq=[]
+ for ptoss in Trans_Pklist:
+ esta=0
+ for ptot in desp:
+ if ptot!=[] and ptoss[4] == ptot[-1]:
+ desp_izq.append(ptot)
+ esta=1
+ continue
+ if esta==0:
+ desp_izq.append([])
+ secc_izq.append(desp_izq)
+
+ secc_der=[]
+ for desp in Desplaz_der:
+ desp_der=[]
+ for ptoss in Trans_Pklist:
+ esta=0
+ for ptot in desp:
+ if ptot!=[] and ptoss[4] == ptot[-1]:
+ desp_der.append(ptot)
+ esta=1
+ continue
+ if esta==0:
+ desp_der.append([])
+ secc_der.append(desp_der)
+
+ secc_der.append(talud_der)
+ secc = secc_izq+[Trans_Pklist]+secc_der
+
+ return secc
+
+
+# ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### #
+
+
+
+def gen_LongProfileGuitarr(ASegmentos,APuntos_caract,puntos,puntos_terr,table_alz,escala,opt1):
+
+ mark_lon,mark_x_dist,mark_y_dist,DistGitarr=opt1.split(',')
+ mark_lon,mark_x_dist,mark_y_dist,DistGitarr=int(mark_lon),int(mark_x_dist),int(mark_y_dist),float(DistGitarr)
+
+ dist_ejes_x=DistGitarr
+
+ alto_guit=dist_ejes_x*6
+ cerox=100
+ ceroy=100+alto_guit
+
+ cotaMax_eje=max([p[2] for p in puntos])
+ cotaMin_eje=min([p[2] for p in puntos])
+ cotaMax_terr=max([p[2] for p in puntos_terr])
+ cotaMin_terr=min([p[2] for p in puntos_terr])
+ cotaMax=max([cotaMax_eje,cotaMax_terr])
+ cotaMin=min([cotaMin_eje,cotaMin_terr])
+
+ dif_y_ref=(cotaMin)-int(cotaMin)
+ cero_y_ref=(cotaMin*escala-ceroy)-dif_y_ref*escala
+
+ # Perfil del eje
+ ASeg_ref=[]
+ for seg in ASegmentos:
+
+ ASeg_ref.append([[cerox+p[4],p[2]*escala-cero_y_ref,0,p[3:]] for p in seg])
+
+ APtos_caract_ref=[[cerox+p[4],p[2]*escala-cero_y_ref,0,p[3],p[4],p[5],p[6],] for p in APuntos_caract]
+ #ptos_ref=[[cerox+p[4],p[2]*escala-cero_y_ref,0] for p in puntos]
+
+ # Perfil del terreno
+ ptos_terr_ref=[[cerox+p[4],p[2]*escala-cero_y_ref,0] for p in puntos_terr]
+
+ dist_orig,dist_par,cras,cterr,croja=[],[],[],[],[]
+ d_ant=0
+ for i in range(0,int(puntos[-1][4]),mark_x_dist):
+ dist_orig.append(puntos[i][4])
+ dist_par.append(puntos[i][4]-d_ant)
+ d_ant=puntos[i][4]
+ dist_orig.append(puntos[-1][4])
+ dist_par.append(puntos[-1][4]-d_ant)
+
+ cras=[p[2] for p in puntos if p[4] in dist_orig ]
+ cterr=[p[2] for p in puntos_terr if p[4] in dist_orig ]
+ croja=[cras[i]-cterr[i] for i in range(len(cras))]
+
+ # Ejes de coordenadas
+ # Eje y
+ eje_y=[[[cerox,ceroy,0,1,"'Elev'"],[cerox,(cotaMax)*escala-cero_y_ref,0,1]]]
+
+ # Marcas eje y
+ mark_y=[]
+ cat=1
+ for i in range(int(ceroy),int((cotaMax)*escala-cero_y_ref),int(mark_y_dist*escala)):
+ labely=(cotaMin+(i-ceroy)/escala)-dif_y_ref
+ mark_y.append([[cerox-mark_lon,i,0,cat,labely],[cerox+mark_lon,i,0,cat,labely]])
+ cat=cat+1
+ mark_y.append([[cerox-mark_lon,cotaMax*escala-cero_y_ref,0,cat,round(cotaMax,2)],
+ [cerox+mark_lon,cotaMax*escala-cero_y_ref,0,cat,cotaMax]])
+
+ # Eje x
+ eje_x=[[[cerox,ceroy,0,1,"'Origin distance'"],[cerox+puntos[-1][4],ceroy,0,1]]]
+
+ # Ejes x de distancias y cotas
+ eje_x_labels=["'Origin distance'","'Partial distance'","'Aling elev'","'Ground elev'","'Cota roja'","'Nada'"]
+ for i in range(1,6):
+ eje_x.append([[cerox,ceroy-i*dist_ejes_x,0,i+1,eje_x_labels[i]],[cerox+puntos[-1][4],ceroy-i*dist_ejes_x,0,i+1]])
+
+ # Marcas eje x
+ mark_x=[]
+ cat=1
+ for j in range(0,6):
+ t=0
+ for i in range(0,int(puntos[-1][4]),mark_x_dist):
+ if j==0: label=round(dist_orig[t],2)
+ elif j==1: label=round(dist_par[t],2)
+ elif j==2: label=round(cras[t],2)
+ elif j==3: label=round(cterr[t],2)
+ elif j==4: label=round(croja[t],2)
+ else: label=0
+ t=t+1
+ mark_x.append([[cerox+i,ceroy-mark_lon-j*dist_ejes_x,0,cat,label],[cerox+i,ceroy+mark_lon-j*dist_ejes_x,0,cat,label]])
+ cat=cat+1
+ if j==0: label=round(dist_orig[-1],2)
+ elif j==1: label=round(dist_par[-1],2)
+ elif j==2: label=round(cras[-1],2)
+ elif j==3: label=round(cterr[-1],2)
+ elif j==4: label=round(croja[-1],2)
+ else: label=0
+ mark_x.append([[cerox+puntos[-1][4],ceroy-mark_lon-j*dist_ejes_x,0,cat,label],[cerox+puntos[-1][4],ceroy+mark_lon-j*dist_ejes_x,0,cat,label]])
+ cat=cat+1
+
+ ptos_eje=[[cerox+p[4],p[5]*escala-cero_y_ref,0,p[3],p[4],p[6],p[7]] for p in table_alz]
+
+ return eje_x,eje_y,mark_x,mark_y,ptos_terr_ref,ptos_eje,ASeg_ref,APtos_caract_ref
+
+
+def gen_TransProfile(Trans,Trans_Terr,Trans_Pklist,Puntos_Talud_izq,Puntos_Talud_der,Desplaz_izq,Desplaz_der,escala,opt1,opt2):
+
+ secc = get_SeccTerr(Trans,Trans_Pklist,Desplaz_izq,Desplaz_der,Puntos_Talud_izq,Puntos_Talud_der)
+ secc = [[row[i] for row in secc] for i in range(len(secc[0]))] # transpuesta
+
+ cerox=100
+ ceroy=100
+
+ mark_lon,mark_x_dist,mark_y_dist=opt1.split(',')
+ mark_lon,mark_x_dist,mark_y_dist=int(mark_lon),int(mark_x_dist),int(mark_y_dist)
+
+ filas,sep_x,sep_y=opt2.split(',')
+ filas,sep_x,sep_y=int(filas),float(sep_x),float(sep_y)
+
+ columnas=len(Trans_Pklist)/filas
+ # Ancho y alto de cada fila/columna
+ h=0
+ ancho_colum,centro_secc,max_filas,min_filas,dif_filas=[],[],[],[],[]
+ for j in range (0,columnas):
+ anchos_colum,centros_colum,maxfila,minfila,dif_fila=[],[],[],[],[]
+ for i in range(0,filas):
+ ancho_secc=round(sqrt((Trans[h][0][0]-Trans[h][-1][0])**2+(Trans[h][0][1]-Trans[h][-1][1])**2),6)
+ centro_s=round(sqrt((Trans[h][0][0]-Trans[h][1][0])**2+(Trans[h][0][1]-Trans[h][1][1])**2),6)
+ anchos_colum.append(ancho_secc)
+ centros_colum.append(centro_s)
+ max_trans=max([p[2] for p in secc[h]+Trans_Terr[h] if p!=[]])
+ maxfila.append(max_trans)
+ min_trans=min([p[2] for p in secc[h]+Trans_Terr[h] if p!=[]])
+ minfila.append(min_trans)
+ dif_fila.append((max_trans-min_trans)*escala)
+ h=h+1
+ ancho_colum.append(max(anchos_colum))
+ centro_secc.append(max(centros_colum))
+ max_filas.append(maxfila)
+ min_filas.append(minfila)
+ dif_filas.append(dif_fila)
+
+ max_filas2 = [[row[i] for row in max_filas] for i in range(len(max_filas[0]))] # transpuesta
+ max_filas2 = [max(row) for row in max_filas2]
+ min_filas2 = [[row[i] for row in min_filas] for i in range(len(min_filas[0]))] # transpuesta
+ min_filas2 = [min(row) for row in min_filas2]
+ alto_filas = [(max_filas2[i]-min_filas2[i])*escala for i in range(len(max_filas2)) ]
+
+ hy_tot = sum(alto_filas)+sep_y*filas
+ ceroy = ceroy + hy_tot
+ puntos_terr_ref,ptos_eje_ref=[],[]
+ ejes_x,ejes_y,mark_x,mark_y=[],[],[],[]
+ h,q,w,t=0,0,0,0
+ for j in range (0,columnas):
+
+ sep_eje_x=sep_x/10
+ # origen, final y centro del eje x
+ orig_x=cerox+(sep_x*j)+sum(ancho_colum[0:j])+sep_eje_x
+ final_x=cerox+(sep_x*j)+sum(ancho_colum[0:j])+ancho_colum[j]+sep_eje_x
+ centro_x=orig_x+centro_secc[j]
+ orig_y=ceroy
+
+ for i in range(0,filas):
+
+ cotaMax=max_filas[j][i]
+ cotaMin=min_filas[j][i]
+ dif_y_ref=(cotaMin)-int(cotaMin)
+
+ # origen y final del eje y
+ orig_y=orig_y-(alto_filas[i])
+ final_y=orig_y+dif_filas[j][i]+dif_y_ref
+
+ cero_y_ref=(cotaMin-orig_y)-dif_y_ref
+ ptos_eje,ptos_terr=[],[]
+ b=0
+ for pto in Trans_Terr[h]:
+ dist=sqrt((Trans_Pklist[h][0]-pto[0])**2+(Trans_Pklist[h][1]-pto[1])**2)
+ if dist==0:b=1
+ if b==0: dist=-dist
+ ptos_terr.append([centro_x+dist,orig_y+(pto[2]-cotaMin)*escala+dif_y_ref,0,h+1,Trans_Pklist[h][4]])
+ b=0
+ for pto in secc[h]:
+ if pto !=[]:
+ dist=sqrt((Trans_Pklist[h][0]-pto[0])**2+(Trans_Pklist[h][1]-pto[1])**2)
+ if dist==0:b=1
+ if b==0: dist=-dist
+ cota_ras=orig_y+(pto[2]-cotaMin)*escala+dif_y_ref
+ t=t+1
+ ptos_eje.append([centro_x+dist,cota_ras,0,t,h+1,Trans_Pklist[h][4],round(cota_ras,2),round(dist,2)])
+
+ # Eje y
+ ejes_y.append([[orig_x-sep_eje_x,final_y,0,h+1,"'Y'"],[orig_x-sep_eje_x,orig_y,0,h+1,"'Y'"]])
+ # Eje x
+ ejes_x.append([[orig_x-sep_eje_x,orig_y,0,h+1,"'X'",format_Pk(Trans_Pklist[h][4])],[final_x,orig_y,0,h+1,"'X'",Trans_Pklist[h][4]]])
+
+ # Marcas eje x desde el centro
+ mark2_x=[]
+ mark2_x.append([[orig_x,orig_y-2*mark_lon,0,q+1,-(centro_x-orig_x)],[orig_x,orig_y+2*mark_lon,0]])
+ q=q+1
+ for k in range(mark_x_dist,int(centro_x-orig_x),mark_x_dist):
+ mark2_x.append([[centro_x-k,orig_y-mark_lon,0,q+1,-((centro_x-orig_x)-((centro_x-orig_x) % mark_x_dist)-k+mark_x_dist)],
+ [centro_x-k,orig_y+mark_lon,0]])
+ q=q+1
+ mark2_x.append([[centro_x,orig_y-2*mark_lon,0,q+1,0.0],[centro_x,orig_y+2*mark_lon,0]])
+ q=q+1
+ for k in range(mark_x_dist,int(final_x-centro_x),mark_x_dist):
+ mark2_x.append([[centro_x+k,orig_y-mark_lon,0,q+1,float(k)],[centro_x+k,orig_y+mark_lon,0]])
+ q=q+1
+ mark2_x.append([[final_x,orig_y-2*mark_lon,0,q+1,final_x-centro_x],[final_x,orig_y+2*mark_lon,0]])
+ q=q+1
+
+ # Marcas eje y
+ for k in range(0,int(final_y)-int(orig_y),int(mark_y_dist*escala)):
+ w=w+1
+ mark_y.append([[orig_x-mark_lon-sep_eje_x,k+orig_y,0,w,k/escala+cotaMin-dif_y_ref],[orig_x+mark_lon-sep_eje_x,k+orig_y,0]])
+ w=w+1
+ mark_y.append([[orig_x-mark_lon-sep_eje_x,final_y,0,w,round(cotaMax-dif_y_ref,2)],[orig_x+mark_lon-sep_eje_x,final_y,0]])
+ mark_x.extend(mark2_x)
+ h=h+1
+ puntos_terr_ref.append(ptos_terr)
+ ptos_eje_ref.append(ptos_eje)
+
+
+ return ejes_x,ejes_y,mark_x,mark_y,puntos_terr_ref,ptos_eje_ref
+
+
+# ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### #
+
+
+def read_Table(EjeMap,layer,columns):
+
+ table=g.read_command('v.out.ascii', input=EjeMap, output='-',
+ format='point', layer=layer, columns=columns, quiet=True)
+ table = [d.split('|') for d in table.splitlines(0)]
+ if len(table[0])<len(columns.split(','))+4:
+ for i in range(len(table)):
+ table[i].insert(2,0.0)
+
+ return table
+
+def read_TablePlant(EjeMap):
+
+ plant = read_Table(EjeMap,2,'pk_eje,radio,a_in,a_out,widening')
+ plant = float_List(plant)
+ return plant
+
+def read_TableAlz(EjeMap):
+
+ alzado = read_Table(EjeMap,3,'pk,elev,slope,kv,l,b')
+ alzado = float_List(alzado)
+ return alzado
+
+def read_TableSection(EjeMap):
+
+ section = read_Table(EjeMap,4,'pk,sec_left,sec_right,type_left,type_right,cut_left,cut_right,fill_left,fill_right,superelev_left,superelev_right')
+ for i in range(len(section)):
+ section[i][:5]=[float(p) for p in section[i][:5]]
+ return section
+
+def read_TableTransv(EjeMap):
+
+ trans = read_Table(EjeMap,5,'pk,dist_left,dist_right,npk')
+ trans = float_List(trans)
+ return trans
+
+#-----------------------------------------------------------------------------
+
+def remove_Layer(EjeMap,ext,layer):
+
+ g.write_command('db.execute', database = database1, driver = 'sqlite',
+ stdin="DELETE FROM "+EjeMap+ext+" WHERE cat"+str(layer)+">=0", input='-', quiet=True)
+ g.run_command('v.edit', map=EjeMap, layer=layer, tool='delete', cats='0-100000', quiet=True)
+ return 0
+
+def remove_Plant(EjeMap):
+
+ remove_Layer(EjeMap,'_Plan',2)
+ return 0
+
+def remove_Alz(EjeMap):
+
+ remove_Layer(EjeMap,'_Vertical',3)
+ return 0
+
+def remove_Section(EjeMap):
+
+ remove_Layer(EjeMap,'_Section',4)
+ return 0
+
+def remove_Transv(EjeMap):
+
+ remove_Layer(EjeMap,'_Transv',5)
+ return 0
+
+#-----------------------------------------------------------------------------
+
+def update_Table(EjeMap,ext,layer,ptsList,columns):
+
+ input_Points(EjeMap,layer,ptsList)
+ update_Layer(EjeMap,ext,layer,ptsList,columns)
+ return 0
+
+def input_Points(EjeMap,layer,ptsList):
+
+ sal=''
+ for i,coord in enumerate(ptsList):
+ sal+='P 1 1'+'\n'
+ sal+=str(coord[0])+' '+str(coord[1])+' '+str(coord[2])+'\n'
+ sal+=str(layer)+' '+str(coord[3])+'\n'
+ #print sal
+ os.system('echo "'+sal+'" | v.edit -n tool=add map='+EjeMap+' input=- --quiet')
+ g.run_command('v.to.db', map=EjeMap, layer=layer, type='point', option='cat', columns='cat'+str(layer), quiet=True)
+ return 0
+
+def update_Layer(EjeMap,ext,layer,ptsList,columns):
+
+ sql=''
+ columns=columns.split(',')
+ for i in range(len(ptsList)):
+ sql+="UPDATE "+EjeMap+ext+" SET "
+ sql+=', '.join(a + "=" +str(b) for a,b in zip(columns,ptsList[i][4:]))
+ sql+=" WHERE cat"+str(layer)+"="+str(ptsList[i][3])+";\n"
+ #print sql
+ g.write_command('db.execute', database = database1, driver = 'sqlite', stdin = sql, input='-', quiet=True)
+ return 0
+
+
+def update_TablePlan(EjeMap,ptsList):
+
+ update_Table(EjeMap,'_Plan',2,ptsList,'pk_eje,radio,a_in,a_out,widening')
+ return 0
+
+def update_TableAlz(EjeMap,ptsList):
+
+ update_Table(EjeMap,'_Vertical',3,ptsList,'pk,elev,slope,kv,l,b')
+ return 0
+
+def update_TableSection(EjeMap,ptsList):
+
+ for i,pts in enumerate(ptsList):
+ ptsList[i][5:]=["'"+str(p)+"'" for p in ptsList[i][5:] if str(p).find("'")==-1]
+ update_Table(EjeMap,'_Section',4,ptsList,'pk,sec_left,sec_right,type_left,type_right,cut_left,cut_right,fill_left,fill_right,superelev_left,superelev_right')
+ return 0
+
+def update_TableTransv(EjeMap,ptsList):
+
+ update_Table(EjeMap,'_Transv',5,ptsList,'pk,dist_left,dist_right,npk')
+ return 0
+
+#-----------------------------------------------------------------------------
+
+def corrige_Alzado(puntos_eje,alz,EjeMap):
+
+ alz.sort(key=lambda x: x[4]) # alz=[x,y,z,cat,Pk,Cota,Pend,Kv,L,B]
+
+ for i in range(1,len(alz),1):
+ if i < len(alz)-1:
+ alz[i][0],alz[i][1]=get_PlantaXY(alz[i][4],puntos_eje)[:2]
+ alz[i][3]=i+1
+ alz[i][6]=(float(alz[i][5])-alz[i-1][5])/(float(alz[i][4])-alz[i-1][4])
+ alz[-1][4]=puntos_eje[-1][-1][7]
+ alz[-1][0],alz[-1][1],alz[-1][2]=puntos_eje[-1][0][:3]
+ remove_Alz(EjeMap)
+ update_TableAlz(EjeMap,alz)
+
+ return 0
+
+def corrige_Section(puntos_eje,secc,EjeMap):
+
+ secc.sort(key=lambda x: float(x[4]))
+ for i in range(1,len(secc),1):
+ if i < len(secc)-1:
+ secc[i][0],secc[i][1]=get_PlantaXY(float(secc[i][4]),puntos_eje)[:2]
+ secc[i][3]=i+1
+ secc[-1][4]=puntos_eje[-1][-1][7]
+ secc[-1][0],secc[-1][1],secc[-1][2]=puntos_eje[-1][0][:3]
+ remove_Section(EjeMap)
+ update_TableSection(EjeMap,secc)
+ return 0
+
+def corrige_Transv(puntos_eje,trans,EjeMap):
+
+ trans.sort(key=lambda x: float(x[4]))
+ for i in range(1,len(trans),1):
+ if i < len(trans)-1:
+ trans[i][0],trans[i][1]=get_PlantaXY(float(trans[i][4]),puntos_eje)[:2]
+ trans[i][3]=i+1
+ trans[-1][4]=puntos_eje[-1][-1][7]
+ trans[-1][0],trans[-1][1],trans[-1][2]=puntos_eje[-1][0][:3]
+ remove_Transv(EjeMap)
+ update_TableTransv(EjeMap,trans)
+ return 0
+
+#-----------------------------------------------------------------------------
+
+def float_List(list):
+ for j,punt in enumerate(list):
+ for i,elem in enumerate(punt):
+
+ if list[j][i]=='':
+ list[j][i] = 0.0
+ else:
+ list[j][i] = float(list[j][i])
+ return list
+
+
+def update_EdgeMap(EjeMap):
+
+ g.message("Reading polygon vertices")
+ verti=g.read_command('v.out.ascii', input=EjeMap, output='-', format='standard', layer=1)
+ verti = [d.strip().split() for d in verti.splitlines(0)]
+ verti = verti[11:-1]
+
+ if len(verti[0])==2:
+ for i in range(len(verti)):
+ verti[i].append('0.0')
+ verti=float_List(verti)
+ pk_eje=[0.0]
+
+ for i in range(len(verti)-1):
+ pk_eje.append(sqrt((verti[i+1][0]-verti[i][0])**2+(verti[i+1][1]-verti[i][1])**2)+pk_eje[-1])
+
+
+ dbs=g.vector_db(EjeMap)
+ #print (dbs)
+ if len(dbs) == 5 : # if layers exist
+
+ g.message("Reading old configuration")
+ planta=read_TablePlant(EjeMap)
+ alzado=read_TableAlz(EjeMap)
+ seccion=read_TableSection(EjeMap)
+ transv=read_TableTransv(EjeMap)
+
+ g.message("Deleting old tables")
+ remove_Plant(EjeMap)
+ remove_Alz(EjeMap)
+ remove_Section(EjeMap)
+ remove_Transv(EjeMap)
+
+ g.message("Updating new tables")
+ for i in range(len(pk_eje)):
+ planta[i][:5]=verti[i][0],verti[i][1],verti[i][2],i+1,pk_eje[i]
+ alzado[0][:5]=verti[0][0],verti[0][1],verti[0][2],1,pk_eje[0]
+ alzado[-1][:5]=verti[-1][0],verti[-1][1],verti[-1][2],len(alzado),pk_eje[-1]
+ seccion[0][:5]=verti[0][0],verti[0][1],verti[0][2],1,pk_eje[0]
+ seccion[-1][:5]=verti[-1][0],verti[-1][1],verti[-1][2],len(seccion),pk_eje[-1]
+ transv[0][:5]=verti[0][0],verti[0][1],verti[0][2],1,pk_eje[0]
+ transv[-1][:5]=verti[-1][0],verti[-1][1],verti[-1][2],len(transv),pk_eje[-1]
+
+ update_TablePlan(EjeMap,planta)
+ update_TableAlz(EjeMap,alzado)
+ update_TableSection(EjeMap,seccion)
+ update_TableTransv(EjeMap,transv)
+
+ else:
+
+ #g.message("Deleting old tables")
+ #remove_Plant(EjeMap)
+ #remove_Alz(EjeMap)
+ #remove_Section(EjeMap)
+ #remove_Transv(EjeMap)
+
+ g.message("Adding tables")
+ g.run_command('v.db.addtable', map=EjeMap, layer=2, key='cat2', table=EjeMap+'_Plan',
+ columns='pk_eje double, radio double, a_in double, \
+ a_out double, widening double', quiet=True)
+ g.run_command('v.db.addtable', map=EjeMap, layer=3, key='cat3', table=EjeMap+'_Vertical',
+ columns='pk double, elev double, slope double, \
+ kv double, l double, b double', quiet=True)
+ g.run_command('v.db.addtable', map=EjeMap, layer=4, key='cat4', table=EjeMap+'_Section',
+ columns='pk double, sec_left varchar(25), sec_right varchar(25), \
+ type_left varchar(25), type_right varchar(25), \
+ cut_left varchar(25), cut_right varchar(25), fill_left varchar(25), fill_right varchar(25), \
+ superelev_left varchar(25), superelev_right varchar(25)', quiet=True)
+ g.run_command('v.db.addtable', map=EjeMap, layer=5, key='cat5', table=EjeMap+'_Transv',
+ columns='pk double, dist_left double, dist_right double, npk double', quiet=True)
+
+ g.message("Updating tables")
+ planta=[]
+ for i in range(len(pk_eje)):
+ planta.append([verti[i][0],verti[i][1],verti[i][2],i+1,pk_eje[i],0.0,0.0,0.0,0.0,0.0,0.0])
+ alzado=[[verti[0][0],verti[0][1],verti[0][2],1,pk_eje[0],0.0,0.0,0.0,0.0,0.0],
+ [verti[-1][0],verti[-1][1],verti[-1][2],2,pk_eje[-1],0.0,0.0,0.0,0.0,0.0]]
+ seccion=[[verti[0][0],verti[0][1],verti[0][2],1,pk_eje[0],'','','','','','','','',''],
+ [verti[-1][0],verti[-1][1],verti[-1][2],2,pk_eje[-1],'','','','','','','','','']]
+ transv=[[verti[0][0],verti[0][1],verti[0][2],1,pk_eje[0],0.0,0.0,0.0],
+ [verti[-1][0],verti[-1][1],verti[-1][2],2,pk_eje[-1],0.0,0.0,0.0]]
+
+ update_TablePlan(EjeMap,planta)
+ update_TableAlz(EjeMap,alzado)
+ update_TableSection(EjeMap,seccion)
+ update_TableTransv(EjeMap,transv)
+
+#-----------------------------------------------------------------------------
+
+
+def write_Plant(segmentos,puntos_caract,puntos_centros,EjeMap,ext):
+
+ write_Polylines(segmentos,EjeMap+ext,1)
+ g.run_command('v.db.addtable', map=EjeMap+ext, layer=1, key='cat', table=EjeMap+ext,
+ columns='pk double, type varchar(25), long double, \
+ param double, azimut double, GRASSRGB varchar(12)', quiet=True)
+ puntos_s=[m for p in puntos_centros for m in p[1:]]
+ update_Layer(EjeMap,ext,'',puntos_s,'pk,type,long,param,azimut')
+
+ g.run_command('v.db.addtable', map=EjeMap+ext, layer=2, key='cat2', table=EjeMap+ext+'_PC',
+ columns='pk double, azimut double, type varchar(25), scat int', quiet=True)
+ update_Table(EjeMap+ext,'_PC',2,puntos_caract,'pk,azimut,type,scat')
+
+ g.run_command('v.db.addtable', map=EjeMap+ext, layer=3, key='cat3', table=EjeMap+ext+'_Centers',
+ columns='type varchar(25), scat int', quiet=True)
+ puntos_c=[p[0] for p in puntos_centros[:-1]]
+ update_Table(EjeMap+ext,'_Centers',3,puntos_c,'type,scat')
+ return 0
+
+
+def write_Alz(asegmentos,p_caract,p_vert,EjeMap,ext):
+
+ write_Polylines(asegmentos,EjeMap+ext,1)
+
+ g.run_command('v.db.addtable', map=EjeMap+ext, layer=1, key='cat', table=EjeMap+ext,
+ columns='pk double, type varchar(25), long double, \
+ param double, GRASSRGB varchar(12)', quiet=True)
+ puntos_s=[m for p in p_vert for m in p[1:]]
+ update_Layer(EjeMap,ext,'',puntos_s,'pk,type,long,param')
+
+ g.run_command('v.db.addtable', map=EjeMap+ext, layer=2, key='cat2', table=EjeMap+ext+'_PC',
+ columns='pk double, type varchar(25), scat int', quiet=True)
+ update_Table(EjeMap+ext,'_PC',2,p_caract,'pk,type,scat')
+
+ g.run_command('v.db.addtable', map=EjeMap+ext, layer=3, key='cat3', table=EjeMap+ext+'_Vert',
+ columns='type varchar(25), scat int', quiet=True)
+ puntos_c=[p[0] for p in p_vert[:-1]]
+
+ update_Table(EjeMap+ext,'_Vert',3,puntos_c,'type,scat')
+ return 0
+
+def write_Despl(dlines_izq,dlines_der,p_caract,p_vert,EjeMap,ext):
+
+ dlines=dlines_izq[:]
+ dlines.extend(dlines_der)
+
+ write_Polylines(dlines,EjeMap+ext,1)
+ sdlines,cats=splitdlines(dlines)
+
+ g.run_command('v.db.addtable', map=EjeMap+ext, layer=1, key='cat', table=EjeMap+ext,
+ columns='long double, param double, GRASSRGB varchar(12)', quiet=True)
+ puntos_s=[p[-1][:3]+[cats[i]]+[p[-1][4]] for i,p in enumerate(sdlines)]
+ update_Layer(EjeMap,ext,'',puntos_s,'long')
+
+ g.run_command('v.db.addtable', map=EjeMap+ext, layer=2, key='cat2', table=EjeMap+ext+'_PC',
+ columns='pk double, azimut double, type varchar(25), scat int', quiet=True)
+ update_Table(EjeMap+ext,'_PC',2,p_caract,'pk,azimut,type,scat')
+
+ return 0
+
+
+def write_Transv(lines,trans_Pklist,p_vert,EjeMap,ext):
+
+ write_Polylines(lines,EjeMap+ext,1)
+ g.run_command('v.db.addtable', map=EjeMap+ext, layer=1, key='cat', table=EjeMap+ext,
+ columns='pk varchar(25), azimut double, type varchar(25),x double, y double, z double, GRASSRGB varchar(12)', quiet=True)
+ #puntos_s=[pto[:9] for line in lines for pto in line if len(pto)>8]
+ puntos_s=[pto[:3]+[i+1]+[format_Pk(pto[4])]+[pto[5]]+["'"+pto[6]+"'"]+pto[:3] for i,pto in enumerate(trans_Pklist)]
+
+ update_Layer(EjeMap,ext,'',puntos_s,'pk,azimut,type,x,y,z')
+
+ return 0
+
+def splitdlines(lines):
+
+ # line = [[],[],...]
+ #[[line11,[],line12,...],line3,line4,...] --> [line11,line12,line3,line4,...]
+ # --> [ 1 , 1 , 2 , 3 ,...]
+ tolines,cats=[],[]
+ for j,line in enumerate(lines):
+ if [] in line:
+ splitlist=[list(group) for k, group in groupby(line, lambda x: x == []) if not k] # split list
+ tolines.extend(splitlist)
+ for i in range(len(splitlist)):
+ cats.append(j+1)
+ else:
+ tolines.append(line)
+ cats.append(j+1)
+
+ return tolines,cats
+
+
+#### ############## Out maps #### ##################
+
+def write_Polylines(lines,name,layer):
+
+ sal_linea=""
+ if layer > 1: tool="add"
+ else: tool="create"
+ tolines,cats=splitdlines(lines)
+
+ for j,line in enumerate(tolines):
+ sal_linea+="L "+str(len(line))+" 1\n"
+ for i,pto in enumerate(line):
+ sal_linea+=str(pto[0])+" "+str(pto[1])+" "+str(pto[2])+"\n"
+ sal_linea+=str(layer)+" "+str(cats[j])+"\n"
+
+ g.write_command('v.in.ascii', flags='nz', output=name, stdin=sal_linea,
+ input='-', format='standard', overwrite=True, quiet=True)
+ #g.write_command('v.edit', flags='n', tool=tool, map=name, input='-', stdin=sal_linea, overwrite=True, quiet=True)
+ return 0
+
+def write_PtosPoints(lines,name):
+
+ sal_puntos=""
+ i=0
+ for j,line in enumerate(lines):
+ for pp in line:
+ if pp==[]: continue
+ sal_puntos+=str(i)+"|"+str(pp[4])+"|"+str(pp[5])+"|"+str(pp[6])+"|"+str(pp[0])+"|"+str(pp[1])+"|"+str(pp[2])+"\n"
+ i=i+1
+ g.write_command('v.in.ascii', flags='nz', output=name, stdin=sal_puntos, input='-',
+ columns='cat int, PK double, Az double, Type varchar(25), x double, y double, z double',
+ x=5, y=6, z=7, cat=1, overwrite=True, quiet=True)
+ return 0
+
+def write_Points(puntos,name):
+
+ sal_puntos=""
+ for i,pp in enumerate(puntos):
+ #if len(pp)==6: pp[4]=pp[5]
+ sal_puntos+=str(i)+"|"+str(pp[4])+"|"+str(pp[5])+"|"+str(pp[6])+"|"+str(pp[0])+"|"+str(pp[1])+"|"+str(pp[2])+"\n"
+ g.write_command('v.in.ascii', flags='nz', output=name, stdin=sal_puntos, input='-',
+ columns='cat int, PK double, Az double, Type varchar(25), x double, y double, z double',
+ x=5, y=6, z=7, cat=1, overwrite=True, quiet=True)
+ return 0
+
+def write_Polyline(puntos,name):
+
+ # Write Polyline
+ sal_linea="L "+str(len(puntos))+" 1\n"
+ for i,pp in enumerate(puntos):
+ sal_linea+=str(pp[0])+" "+str(pp[1])+" "+str(pp[2])+"\n"
+ sal_linea+="1 1\n"
+ g.write_command('v.in.ascii', flags='nz', output=name, stdin=sal_linea,
+ input='-', format='standard', overwrite=True, quiet=True)
+ return 0
+
+def write_Polygon(puntos,name):
+
+ # Write Polygon
+ sal_linea="B "+str(len(puntos)+1)+" 1\n"
+ h=0
+ for pp in puntos:
+ if pp==[]: h=h+1; continue
+ sal_linea+=str(pp[0])+" "+str(pp[1])+"\n"
+
+ sal_linea+=str(puntos[0][0])+" "+str(puntos[0][1])+"\n"
+ sal_linea+="1 1\n"
+ #print sal_linea
+ g.write_command('v.in.ascii', flags='nz', output=name, stdin=sal_linea,
+ input='-', format='standard', overwrite=True, quiet=True)
+ return 0
+
+def write_Polygonos(lines,name):
+
+ sal_linea=""
+ for j,line in enumerate(lines):
+ if [] in line: continue
+ sal_linea+="B "+str(len(line)+1)+" 1\n"
+ for i,pp in enumerate(line):
+
+ sal_linea+=str(pp[0])+" "+str(pp[1])+"\n"
+
+ sal_linea+=str(line[0][0])+" "+str(line[0][1])+"\n"
+ sal_linea+="1 "+str(j+1)+"\n"
+ #print sal_linea
+ g.write_command('v.in.ascii', flags='nz', output=name, stdin=sal_linea,
+ input='-', format='standard', overwrite=True, quiet=True)
+ return 0
+
+
+def concatenate_lines(lines):
+
+ # line = [[pto],[pto],...]
+ # despl_izq [[line11,[],line12,...],line3,line4,...] + puntos -->
+ # [[ line11,line3 ],[ line12,line3 ],[ line3,line4 ],...
+ line3=[]
+ for j,line in enumerate(lines[:-1]):
+ if [] in line:
+ splitlist=[list(group) for k, group in groupby(line, lambda x: x == []) if not k] # split list
+ for line2 in splitlist:
+ list2=[]
+ for i,pto in enumerate(line2):
+ list2.append(lines[j+1][int(pto[3])])
+ line3.append(line2+list2[::-1])
+ else:
+ tline=lines[j+1]
+ line3.append(line+tline[::-1])
+
+ return line3
+
+def generate_DesplazAreas(despl_izq,puntos,despl_der):
+
+ new_izq=despl_izq[:]
+ new_izq.append(puntos)
+ new_izq=concatenate_lines(new_izq)
+
+ new_der=despl_der[::-1]
+ new_der.append(puntos)
+ new_der=concatenate_lines(new_der)
+
+ lines=new_izq+new_der[::-1]
+ for i,line in enumerate(lines):
+ if line[0] == line[-1]: del lines[i][0]
+
+ return lines
+
+
+def lista_PksEje(Puntos_Eje,Puntos_EjeAlz,puntos,table_alz,table_secc,table_transv):
+
+ pkpuntos=[round(p[4],6) for p in puntos]
+ for alz in table_alz+table_secc+table_transv:
+ if round(alz[4],6) not in pkpuntos:
+ pkpuntos.append(alz[4])
+ x,y,z,cat,PK,Az,tipo,align=get_PlantaXY(round(alz[4],6),Puntos_Eje)[:8]
+ z=get_Cota([[x,y,z,cat,PK,Az,tipo,align]],Puntos_EjeAlz)[0][2]
+ puntos.append([x,y,z,cat,PK,Az,tipo,align])
+
+ puntos.sort(key=lambda x: x[4])
+ for i,pp in enumerate(puntos):
+ puntos[i][3]=i
+
+ return puntos
+
+
+# ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### #
+# ### Main
+# ###
+
+def main():
+
+ global database1,MapSet,driver
+
+ MapSet = g.gisenv()['MAPSET']
+
+ EjeMap=options['edge']
+ if '@' in EjeMap:
+ NameMap,MapSet=EjeMap.split('@')
+ else:
+ NameMap=EjeMap
+
+
+ if g.find_file(EjeMap, element = 'vector', mapset = MapSet)['file'] is '':
+ g.fatal(_("Vector map %s not found in current mapset") % options['edge'])
+ try:
+ f = g.vector_db(options['edge'])[1]
+ except KeyError:
+ #g.fatal(_("There is no table connected to this map. \
+ #Run v.db.connect or v.db.addtable first."))
+ g.message("There is no table connected to this map. Adding table in layer 1.")
+ g.run_command('v.db.addtable', map=EjeMap, quiet=True)
+ f = g.vector_db(options['edge'])[1]
+
+ table = f['table']
+ database1 = f['database']
+ driver = f['driver']
+
+ #g.message("Creating backup")
+ #g.run_command('g.rename', vect=NameMap+"_old"+","+NameMap+"_old2")
+ #g.run_command('g.copy', vect=NameMap+","+NameMap+"_old")
+ #g.run_command('g.remove', flags='f', vect=NameMap+"_old2")
+ #g.message("Finish backup")
+
+
+#### Edges section ####
+
+ if flags['i']:
+
+ pklist=options["pklist"].split(",")
+ pklist=[float(p) for p in pklist]
+ pklist.sort()
+
+ if options["pklayer"]=="Vertical":
+ table=read_TableAlz(NameMap)
+ elif options["pklayer"]=="Section":
+ table=read_TableSection(NameMap)
+ elif options["pklayer"]=="Trans":
+ table=read_TableTransv(NameMap)
+ else:
+ g.message("Layer no editable")
+
+ tablepks=[]
+ for i,line in enumerate(table):
+ tablepks.append(line[4])
+ pklist2=[]
+ for p in pklist:
+ if p not in tablepks: pklist2.append(p)
+
+ if pklist2!=[]:
+
+ tablepks=tablepks+pklist2
+ tablepks.sort()
+ posi=tablepks.index(pklist2[0])
+ table2=[]
+ categ=len(table)+1
+ for i,pk in enumerate(pklist2):
+ table2.append([table[posi-1][0]+i,table[posi-1][1]+i,table[posi-1][2]+i,categ+i,pk]+table[posi-1][5:])
+
+ table_plant= read_TablePlant(NameMap)
+ Puntos_Eje=get_PtosEjePlanta(table_plant)
+
+ if options["pklayer"]=="Vertical":
+ update_TableAlz(NameMap,table2)
+ table_alz=read_TableAlz(NameMap)
+ corrige_Alzado(Puntos_Eje,table_alz,NameMap)
+
+ elif options["pklayer"]=="Section":
+ update_TableSection(NameMap,table2)
+ table_secc=read_TableSection(NameMap)
+ corrige_Section(Puntos_Eje,table_secc,NameMap)
+
+ elif options["pklayer"]=="Trans":
+ update_TableTransv(NameMap,table2)
+ table_transv=read_TableTransv(NameMap)
+ corrige_Transv(Puntos_Eje,table_transv,NameMap)
+ else: g.message("Pk exist or list empty")
+
+
+ if flags['n']:
+
+ g.message("Creating/Updating edge map")
+ update_EdgeMap(NameMap)
+
+
+ if flags['u']:
+
+ g.message("Reading edge tables")
+ ### Lectura tablas eje
+
+ table_plant= read_TablePlant(NameMap)
+ Puntos_Eje=get_PtosEjePlanta(table_plant)
+
+
+ table_alz=read_TableAlz(NameMap)
+ corrige_Alzado(Puntos_Eje,table_alz,NameMap)
+ table_secc=read_TableSection(NameMap)
+ corrige_Section(Puntos_Eje,table_secc,NameMap)
+ table_secc=read_TableSection(NameMap)
+
+ table_transv=read_TableTransv(NameMap)
+ corrige_Transv(Puntos_Eje,table_transv,NameMap)
+
+ ######################################################################
+ g.message("Generating alings")
+
+ Puntos,Segmentos,Puntos_caract,Puntos_centros=generate_Pts(Puntos_Eje,1,1,1,int(options['intervR']),int(options['intervC']))
+ Puntos_EjeAlz=get_PtosEjeAlzado(table_alz)
+
+ Puntos=get_Cota(Puntos,Puntos_EjeAlz)
+
+ Puntos_caract=get_Cota(Puntos_caract,Puntos_EjeAlz)
+ #Segmentos=get_Cota(Segmentos,Puntos_EjeAlz,0)
+ Puntos2=lista_PksEje(Puntos_Eje,Puntos_EjeAlz,Puntos[:],table_alz,table_secc,table_transv)
+
+ ASegmentos,APuntos_caract,APuntos_vert=generate_PtsAlz(Puntos2,Puntos_Eje,Puntos_EjeAlz)
+
+ Desplazados_izq,Desplazados_der,DPtos_caract=generate_Desplaz(table_plant,table_secc,Puntos2)
+
+ Desplaz_Areas=generate_DesplazAreas(Desplazados_izq[:],Puntos2[:],Desplazados_der[:])
+
+ Transversales,Trans_Pklist=generate_Transv(Puntos,table_transv)
+
+ Transv_Despl = get_TransDespl(Transversales,Desplazados_izq,Desplazados_der)
+
+ Transv_Discr = discr_Lines(Transv_Despl,1)
+
+ T_Pklist=[p[4] for p in Trans_Pklist]
+
+ ######################################################################
+ g.message("Writing Plant/Raised/Section maps")
+
+ if flags['y']:
+ #write_Points(Puntos,name+"_Ptos")
+ if re.search(r'^_', options['plantpoly']): name1=NameMap+options['plantpoly']
+ else: name1=options['plantpoly']
+ write_Polyline(Puntos,name1)
+
+ if flags['h']:
+ if re.search(r'^_', options['plant']): name1=NameMap+options['plant']
+ else: name1=options['plant']
+ write_Plant(Segmentos,Puntos_caract,Puntos_centros,name1,'')
+
+ if flags['v']:
+ if re.search(r'^_', options['raised']): name1=NameMap+options['raised']
+ else: name1=options['raised']
+ write_Alz(ASegmentos,APuntos_caract,APuntos_vert,name1,'')
+
+ if flags['d']:
+ if re.search(r'^_', options['displ']): name1=NameMap+options['displ']
+ else: name1=options['displ']
+ write_Despl(Desplazados_izq,Desplazados_der,DPtos_caract,[],name1,'')
+
+ if flags['a']:
+ if re.search(r'^_', options['displ_area']): name1=NameMap+options['displ_area']
+ else: name1=options['displ_area']
+ if Desplaz_Areas != []:
+ write_Polygonos(Desplaz_Areas,NameMap+"_tmp2")
+ g.run_command('v.centroids', input=NameMap+"_tmp2", output=name1, overwrite=True, quiet=True)
+ g.run_command('g.remove', vect=NameMap+"_tmp2", quiet=True)
+
+ if flags['c']:
+
+ if re.search(r'^_', options['cross']): name1=NameMap+options['cross']
+ else: name1=options['cross']
+
+ write_Transv(Transv_Despl,Trans_Pklist,[],name1,'')
+
+ if flags['r']:
+
+ conj=[]
+ displ_opt=options['displ_opt'].split(',')
+
+ if Desplazados_izq != []:
+ if 'displ_left0' in displ_opt:
+ conj.append([pto for pto in Desplazados_izq[0] if pto != [] and pto[-1] in T_Pklist])
+
+ if 'displ_left-1' in displ_opt:
+ conj.append([pto for pto in Desplazados_izq[-1] if pto != [] and pto[-1] in T_Pklist])
+
+ if Desplazados_der != []:
+ if 'displ_rigth0' in displ_opt:
+ conj.append([pto for pto in Desplazados_der[0] if pto != [] and pto[-1] in T_Pklist])
+
+ if 'displ_rigth-1' in displ_opt:
+ conj.append([pto for pto in Desplazados_der[-1] if pto != [] and pto[-1] in T_Pklist])
+
+ if re.search(r'^_', options['crossdispl']): name1=NameMap+options['crossdispl']
+ else: name1=options['crossdispl']
+ write_PtosPoints(conj,name1)
+
+
+ if flags['k']:
+
+ npk,mpk,dist,m=options['pkopt'].split(',')
+ npk,mpk,dist,m=int(npk),int(mpk),float(dist),float(m)
+ Pks,Pkslist=get_Trans(Puntos,npk,mpk,m,dist,dist,dist,dist,0,0,int(Puntos[-1][4]))
+
+ if re.search(r'^_', options['pks']): name1=NameMap+options['pks']
+ else: name1=options['pks']
+ write_Transv(Pks,Trans_Pklist,[],name1,'')
+
+ ##################################################################
+ if options['dem']:
+
+ g.message("Reading terrain map")
+ Terreno_Array=get_Terrain(options['dem'])
+
+ ##################################################################
+ g.message("Generating terrain maps")
+
+ Puntos_Talud_izq,Puntos_Talud_der=generate_Taludes(Puntos,Desplazados_izq,Desplazados_der,table_secc,Terreno_Array)
+
+ Taludes_Areas=generate_TaludesAreas(Puntos_Talud_izq,Desplazados_izq,Desplazados_der,Puntos_Talud_der)
+
+ Transversales_Terreno = drape_LinesPoints(Transv_Discr,Terreno_Array)
+ Puntos_Long_Terreno = drape_Points(Puntos,Terreno_Array)
+
+ ##################################################################
+ g.message("Writing terrain maps")
+
+ if flags['t']:
+ if re.search(r'^_', options['outtlong']): name1=NameMap+options['outtlong']
+ else: name1=options['outtlong']
+ write_Polyline(Puntos_Long_Terreno,name1)
+
+ if flags['q']:
+ if re.search(r'^_', options['outtcross']): name1=NameMap+options['outtcross']
+ else: name1=options['outtcross']
+ write_Transv(Transversales_Terreno,Trans_Pklist,[],name1,'')
+
+ if flags['s']:
+ if re.search(r'^_', options['outslope']): name1=NameMap+options['outslope']
+ else: name1=options['outslope']
+ #write_Polyline(Puntos_Talud_izq,NameMap+"_Talud_izq")
+ #write_Polyline(Puntos_Talud_der,NameMap+"_Talud_der")
+ write_Polylines([Puntos_Talud_izq,Puntos_Talud_der],name1,1)
+
+ if flags['e']:
+ if re.search(r'^_', options['outslopeareas']): name1=NameMap+options['outslopeareas']
+ else: name1=options['outslopeareas']
+ if Taludes_Areas != []:
+ write_Polygonos(Taludes_Areas,NameMap+"_tmp3")
+ g.run_command('v.centroids', input=NameMap+"_tmp3", output=name1, overwrite=True, quiet=True)
+ g.run_command('g.remove', vect=NameMap+"_tmp3", quiet=True)
+
+ if flags['p']:
+
+ conj=[]
+ pts_opt=options['pts_opt'].split(',')
+ if 'slope_left' in pts_opt:
+ conj.append(Puntos_Talud_izq)
+ if 'displ_left' in pts_opt:
+ conj.extend(Desplazados_izq)
+ if 'edge' in pts_opt:
+ conj.append(Puntos)
+ if 'displ_rigth' in pts_opt:
+ conj.extend(Desplazados_der)
+ if 'slope_rigth' in pts_opt:
+ conj.append(Puntos_Talud_der)
+
+ if re.search(r'^_', options['outpoints']): name1=NameMap+options['outpoints']
+ else: name1=options['outpoints']
+ write_PtosPoints(conj,name1)
+
+ if flags['b']:
+
+ conj=[]
+ break_opt=options['break_opt'].split(',')
+ if 'slope_left' in break_opt:
+ conj.append(Puntos_Talud_izq)
+ if 'displ_left' in break_opt:
+ conj.extend(Desplazados_izq)
+ if 'edge' in break_opt:
+ conj.append(Puntos)
+ if 'displ_rigth' in break_opt:
+ conj.extend(Desplazados_der)
+ if 'slope_rigth' in break_opt:
+ conj.append(Puntos_Talud_der)
+
+ if re.search(r'^_', options['outbreak']): name1=NameMap+options['outbreak']
+ else: name1=options['outbreak']
+ write_Polylines(conj,name1,1)
+
+ if flags['o']:
+
+ conj=[]
+ hull_opt=options['hull_opt'].split(',')
+ if 'slope_left' in hull_opt and 'slope_rigth' in hull_opt:
+ conj= generate_ContornoAreas(Puntos,Puntos_Talud_izq,Desplazados_izq,Desplazados_der,Puntos_Talud_der)
+ elif 'slope_left' in hull_opt:
+ conj=rellenar_linea(Puntos,Puntos_Talud_izq,Desplazados_izq)
+ elif 'slope_rigth' in hull_opt:
+ conj=rellenar_linea(Puntos,Puntos_Talud_der,Desplazados_der)
+
+ if re.search(r'^_', options['outhull']): name1=NameMap+options['outhull']
+ else: name1=options['outhull']
+ write_Polygon(conj,NameMap+"_tmp1")
+
+ g.run_command('v.centroids', input=NameMap+"_tmp1", output=name1, overwrite=True, quiet=True)
+ #g.run_command('v.to.rast', input=NameMap+"_Contorno", output=options['outhull'], use='val', overwrite=True, quiet=True)
+ g.run_command('g.remove', vect=NameMap+"_tmp1", quiet=True)
+
+ ##################################################################
+
+ if flags['l']:
+
+ scale=float(options['LPScale'])
+ opt1=options['LPopt']
+
+ if re.search(r'^_', options['LPterrain']): nameTerr=NameMap+options['LPterrain']
+ else: nameTerr=options['LPterrain']
+
+ if re.search(r'^_', options['LPras']): nameRas=NameMap+options['LPras']
+ else: nameRas=options['LPras']
+
+ (eje_x,eje_y,mark_x,mark_y,ptos_terr_ref,ptos_eje_ref,ASeg_ref,
+ APtos_caract_ref)=gen_LongProfileGuitarr(ASegmentos,APuntos_caract,Puntos,Puntos_Long_Terreno,table_alz,scale,opt1)
+
+ # Terreno alzado
+ write_Polyline(ptos_terr_ref,nameTerr)
+ #write_Polyline(ptos_ref,nameRas)
+
+ # Eje alzado
+ if re.search(r'^_', options['LPejeref']): nameEdgeRef=NameMap+options['LPejeref']
+ else: nameEdgeRef=options['LPejeref']
+
+ write_Polyline(ptos_eje_ref,nameEdgeRef)
+ g.run_command('v.db.addtable', map=nameEdgeRef, layer=1, key='cat', table=nameEdgeRef+"_Lin", columns='label varchar(25)', quiet=True)
+ g.run_command('v.db.addtable', map=nameEdgeRef, layer=2, key='cat2', table=nameEdgeRef+"_Ptos", columns='pk double,slope double,kv double', quiet=True)
+
+ update_Table(nameEdgeRef,'_Ptos',2,ptos_eje_ref[1:-1],'pk,slope,kv')
+
+ # Segmentos alzado
+ write_Polylines(ASeg_ref,nameRas,1)
+
+ g.run_command('v.db.addtable', map=nameRas, layer=1, key='cat', table=nameRas,
+ columns='pk double, type varchar(25), long double, \
+ param double, GRASSRGB varchar(12)', quiet=True)
+ puntos_s=[m for p in APuntos_vert for m in p[1:]]
+ update_Layer(nameRas,"",'',puntos_s,'pk,type,long,param')
+
+ g.run_command('v.db.addtable', map=nameRas, layer=2, key='cat2', table=nameRas+'_PC',
+ columns='pk double, type varchar(25), scat int', quiet=True)
+ update_Table(nameRas,'_PC',2,APtos_caract_ref,'pk,type,scat')
+
+ if flags['m']:
+
+ if re.search(r'^_', options['LPedgeX']): nameEdgeX=NameMap+options['LPedgeX']
+ else: nameEdgeX=options['LPedgeX']
+
+ if re.search(r'^_', options['LPedgeXmarks']): nameEdgeXmarks=NameMap+options['LPedgeXmarks']
+ else: nameEdgeXmarks=options['LPedgeXmarks']
+
+ if re.search(r'^_', options['LPedgeY']): nameEdgeY=NameMap+options['LPedgeY']
+ else: nameEdgeY=options['LPedgeY']
+
+ if re.search(r'^_', options['LPedgeYmarks']): nameEdgeYmarks=NameMap+options['LPedgeYmarks']
+ else: nameEdgeYmarks=options['LPedgeYmarks']
+
+ write_Polylines(eje_x,nameEdgeX,1)
+ write_Polylines(mark_x,nameEdgeXmarks,1)
+ write_Polylines(eje_y,nameEdgeY,1)
+ write_Polylines(mark_y,nameEdgeYmarks,1)
+
+ g.run_command('v.db.addtable', map=nameEdgeX, layer=1, key='cat', table=nameEdgeX+"_x", columns='label varchar(25)', quiet=True)
+ eje_x2=[p[0] for p in eje_x]
+ update_Layer(nameEdgeX,'_x','',eje_x2,'label')
+
+ g.run_command('v.db.addtable', map=nameEdgeXmarks, layer=1, key='cat', table=nameEdgeXmarks+"_Marks", columns='label double', quiet=True)
+ mark_x2=[p[0] for p in mark_x]
+ update_Layer(nameEdgeXmarks,'_Marks','',mark_x2,'label')
+
+ g.run_command('v.db.addtable', map=nameEdgeY, layer=1, key='cat', table=nameEdgeY+"_y", columns='label varchar(25)', quiet=True)
+ eje_y2=[p[0] for p in eje_y]
+ update_Layer(nameEdgeY,'_y','',eje_y2,'label')
+
+ g.run_command('v.db.addtable', map=nameEdgeYmarks, layer=1, key='cat', table=nameEdgeYmarks+"_Marks", columns='label double', quiet=True)
+ mark_y2=[p[0] for p in mark_y]
+ update_Layer(nameEdgeYmarks,'_Marks','',mark_y2,'label')
+
+
+ if flags['f']:
+
+ scale2=float(options['LTScale'])
+ opt1=options['LTopt']
+ opt2=options['LTopt2']
+
+ if re.search(r'^_', options['LTterrain']): nameTTerr=NameMap+options['LTterrain']
+ else: nameTTerr=options['LTterrain']
+
+ if re.search(r'^_', options['LTras']): nameTRas=NameMap+options['LTras']
+ else: nameTRas=options['LTras']
+
+ (ejes_x,ejes_y,mark_x,mark_y,ptos_terr_ref,
+ ptos_eje)=gen_TransProfile(Transversales,Transversales_Terreno,Trans_Pklist,Puntos_Talud_izq,Puntos_Talud_der,Desplazados_izq,Desplazados_der,scale2,opt1,opt2)
+
+ # Terreno
+ write_Polylines(ptos_terr_ref,nameTTerr,1)
+
+ write_Polylines(ptos_eje,nameTRas,1)
+
+ g.run_command('v.db.addtable', map=nameTRas, layer=1, key='cat', table=nameTRas, columns='pk varchar(25)', quiet=True)
+ ptos_eje2=[p[0][:3]+p[0][4:] for p in ptos_eje]
+ update_Layer(nameTRas,'','',ptos_eje2,'pk')
+
+ g.run_command('v.db.addtable', map=nameTRas, layer=2, key='cat2', table=nameTRas+"_ptos", columns='section varchar(25), pk double, elev double, dist double', quiet=True)
+ ptos_eje2=[p for line in ptos_eje for p in line]
+ update_Table(nameTRas,'_ptos',2,ptos_eje2,'section,pk,elev,dist')
+
+ if flags['g']:
+
+ if re.search(r'^_', options['LTedgeX']): nameTEdgeX=NameMap+options['LTedgeX']
+ else: nameTEdgeX=options['LTedgeX']
+
+ if re.search(r'^_', options['LTedgeXmarks']): nameTEdgeXmarks=NameMap+options['LTedgeXmarks']
+ else: nameTEdgeXmarks=options['LTedgeXmarks']
+
+ if re.search(r'^_', options['LTedgeY']): nameTEdgeY=NameMap+options['LTedgeY']
+ else: nameTEdgeY=options['LTedgeY']
+
+ if re.search(r'^_', options['LTedgeYmarks']): nameTEdgeYmarks=NameMap+options['LTedgeYmarks']
+ else: nameTEdgeYmarks=options['LTedgeYmarks']
+
+ write_Polylines(ejes_x,nameTEdgeX,1)
+ write_Polylines(mark_x,nameTEdgeXmarks,1)
+ write_Polylines(ejes_y,nameTEdgeY,1)
+ write_Polylines(mark_y,nameTEdgeYmarks,1)
+
+ g.run_command('v.db.addtable', map=nameTEdgeX, layer=1, key='cat', table=nameTEdgeX+"_x", columns='label varchar(25),pk varchar(10)', quiet=True)
+ ejes_x2=[p[0] for p in ejes_x]
+ update_Layer(nameTEdgeX,'_x','',ejes_x2,'label,pk')
+
+ g.run_command('v.db.addtable', map=nameTEdgeXmarks, layer=1, key='cat', table=nameTEdgeXmarks+"_Marks", columns='label double', quiet=True)
+ mark_x2=[p[0] for p in mark_x]
+ update_Layer(nameTEdgeXmarks,'_Marks','',mark_x2,'label')
+
+ g.run_command('v.db.addtable', map=nameTEdgeY, layer=1, key='cat', table=nameTEdgeY+"_y", columns='label varchar(25)', quiet=True)
+ eje_y2=[p[0] for p in ejes_y]
+ update_Layer(nameTEdgeY,'_y','',eje_y2,'label')
+
+ g.run_command('v.db.addtable', map=nameTEdgeYmarks, layer=1, key='cat', table=nameTEdgeYmarks+"_Marks", columns='label double', quiet=True)
+ mark_y2=[p[0] for p in mark_y]
+ update_Layer(nameTEdgeYmarks,'_Marks','',mark_y2,'label')
+
+
+###########################################################################
+
+ sys.exit(0)
+
+if __name__ == "__main__":
+ options, flags = g.parser()
+ main()
Property changes on: grass-addons/grass7/vector/v.civil/v.civil.road/v.civil.road.py
___________________________________________________________________
Added: svn:executable
+ *
Added: svn:mime-type
+ text/x-python
Added: svn:eol-style
+ native
Added: grass-addons/grass7/vector/v.civil/v.civil.tools/Makefile
===================================================================
--- grass-addons/grass7/vector/v.civil/v.civil.tools/Makefile (rev 0)
+++ grass-addons/grass7/vector/v.civil/v.civil.tools/Makefile 2014-02-26 00:37:16 UTC (rev 59141)
@@ -0,0 +1,9 @@
+# fix this relative to include/
+# or use absolute path to the GRASS source code
+MODULE_TOPDIR = ../..
+
+PGM = v.civil.tools
+
+include $(MODULE_TOPDIR)/include/Make/Script.make
+
+default: script
Property changes on: grass-addons/grass7/vector/v.civil/v.civil.tools/Makefile
___________________________________________________________________
Added: svn:mime-type
+ text/x-makefile
Added: svn:eol-style
+ native
Added: grass-addons/grass7/vector/v.civil/v.civil.tools/v.civil.tools.html
===================================================================
--- grass-addons/grass7/vector/v.civil/v.civil.tools/v.civil.tools.html (rev 0)
+++ grass-addons/grass7/vector/v.civil/v.civil.tools/v.civil.tools.html 2014-02-26 00:37:16 UTC (rev 59141)
@@ -0,0 +1,94 @@
+<h2>DESCRIPTION</h2>
+
+<em>v.civil.tools</em> tools for generate edges, roundabout and street platforms.
+
+<h3>Align</h3>
+
+<p>
+Generates an alignment from a start point. This point can be defined as a point on a curve, then radio, center and
+a point of the curve (first point of alignment) must be given. Or the point can be defined as a point of straight,
+then the radio will be 0, the center of the curve now is the fist point of the straight and the point of the curve is
+the second point of the straight (first point of alignment).
+
+<p>
+In both cases, the parameter of the clothoid must be given. This clothoid ends in the curve 2.
+
+<p>
+The curve 2 can be a list of curves and straights. It is defined by radio, the same input clothoid above
+(with the sign changed in some cases), longitude of the curve, and the parameter of out clothoid.
+
+<p>
+To close the alignment or connect to another curve with different clockwise, a straight must be defined as curve 2, with
+radio, clothoid in and out equal to 0, and with a required longitude.
+
+<p>
+This option call v.civil.road to create the edge and plan alignment. Example:
+
+<p>
+<pre>v.civil.tools -c --o outmap=b1 cradio1=1083.5 ccenter=551139.061,4083285.68 casal1=-290 cpoint1=552102.873,4083780.316 \
+ cradio2=820,2495,296,0,-244,0,-2607,0 \
+ caent2 =290,-290,156,0,115,0,300,0 \
+ casal2 =290,-156,156,0,115,0,0,0 \
+ len2 =115.5,49,245,0.01,428,435,1,1 </pre>
+
+<pre>v.civil.tools -c --o outmap=b2 cradio1=584 ccenter=551313.918903,4082758.71597 casal1=-141 cpoint1=551342.108908,4083342.37232 \
+ cradio2=320,2500,577,0,50,0 \
+ caent2=141,-141,219,0,40,0 \
+ casal2=141,-219,219,0,40,0 \
+ len2=23.5,10,192,451.5,0.001,0.001</pre>
+
+
+
+<h3>Round</h3>
+
+<p>
+Generates an alignment of a roundabout with mayor and minus radio, azimut, and the center of it.
+
+<p>
+This option call v.civil.road to create the edge and plan alignment with the name given in the required section.
+
+
+
+<h3>Street1</h3>
+<p>
+This section give the parameters to get a curve between the platforms lines of two edges.
+
+<p>
+The first parameter is the radio for the two platforms lines.
+
+<p>
+The rest of the parameters of the two edges can be obtained with their respective plan maps.
+
+<p>
+The output is an auxiliary map with lines representing the center of curve and initial and final pks of curve.
+
+
+<h3>Street2</h3>
+<p>
+This section give the parameters to get a curve between the platforms lines of a edge and the platforms lines of an roundabout.
+
+<p>
+The first parameter is the direction of the edge, if the edge go into or go out from the roundabout.
+
+<p>
+Second parameter is the radio for the two platforms lines.
+
+<p>
+The rest of the parameters of the two edges can be obtained with their respective plan maps.
+
+<p>
+The output is an auxiliary map with lines representing the center of curve and initial and final pks of curve.
+
+<h2> AUTHOR</h2>
+<p>
+Jesus Fernandez-Capel Rosillo<br>
+Civil Engineer, Spain<br>
+jfc at alcd net<br>
+
+
+
+
+
+
+
+
Property changes on: grass-addons/grass7/vector/v.civil/v.civil.tools/v.civil.tools.html
___________________________________________________________________
Added: svn:mime-type
+ text/html
Added: svn:keywords
+ Author Date Id
Added: svn:eol-style
+ native
Added: grass-addons/grass7/vector/v.civil/v.civil.tools/v.civil.tools.py
===================================================================
--- grass-addons/grass7/vector/v.civil/v.civil.tools/v.civil.tools.py (rev 0)
+++ grass-addons/grass7/vector/v.civil/v.civil.tools/v.civil.tools.py 2014-02-26 00:37:16 UTC (rev 59141)
@@ -0,0 +1,1550 @@
+#!/usr/bin/env python
+# -*- coding: utf-8
+############################################################################
+#
+# MODULE: v.civil.tools, v0.5.0
+#
+# AUTHOR: Jesús Fernández-Capel Rosillo
+# Civil Engineer, Spain
+# jfc at alcd net
+#
+# PURPOSE: Road tools for use with v.civil.road
+#
+# COPYRIGHT: (c) 2014 Jesús Fernández-Capel Rosillo.
+#
+# This program is free software under the GNU General Public
+# License (>=v2). Read the file COPYING that comes with GRASS
+# for details.
+#
+#############################################################################
+
+#%module
+#% description: Road tool for use with v.civil.road
+#% keywords: vector, Road, civil
+#%end
+
+#### Input section ####
+
+#%option G_OPT_V_OUTPUT
+#% key: outmap
+#% description: Name for output map
+#% required: yes
+#%end
+
+
+
+#### Curves section ####
+
+#%flag
+#% key: c
+#% description: Generate align
+#% guisection: Align
+#%end
+
+
+#%option
+#% key: cradio1
+#% type: double
+#% description: Radio of curve1
+#% required: no
+#% guisection: Align
+#%end
+
+#%option G_OPT_M_COORDS
+#% key: ccenter
+#% description: Center of curve1
+#% required: no
+#% guisection: Align
+#%end
+
+#%option
+#% key: casal1
+#% type: double
+#% description: Parameter A of clothoid out
+#% required: no
+#% guisection: Align
+#%end
+
+#%option G_OPT_M_COORDS
+#% key: cpoint1
+#% type: string
+#% description: Point of curve1
+#% required: no
+#% guisection: Align
+#%end
+
+#%option
+#% key: cradio2
+#% type: string
+#% description: Radios of curve2
+#% required: no
+#% guisection: Align
+#%end
+
+#%option
+#% key: caent2
+#% type: string
+#% description: Parameter A of clothoid in
+#% required: no
+#% guisection: Align
+#%end
+
+#%option
+#% key: casal2
+#% type: string
+#% description: Parameter A of clothoid out
+#% required: no
+#% guisection: Align
+#%end
+
+#%option
+#% key: len2
+#% type: string
+#% description: Longitud of alignment
+#% required: no
+#% guisection: Align
+#%end
+
+
+#### Round section ####
+
+#%flag
+#% key: r
+#% description: Write roundabout edge
+#% guisection: Round
+#%end
+
+#%option
+#% key: rround1
+#% type: string
+#% description: Minus radio for roundabout edge
+#% required: no
+#% guisection: Round
+#%end
+
+#%option
+#% key: rround2
+#% type: string
+#% description: Mayor radio for roundabout edge
+#% required: no
+#% guisection: Round
+#%end
+
+#%option
+#% key: azround
+#% type: string
+#% description: Azimut for roundabout start point
+#% required: no
+#% guisection: Round
+#%end
+
+#%option
+#% key: cround
+#% type: string
+#% description: Center for roundabout
+#% required: no
+#% guisection: Round
+#%end
+
+
+##### Tools section ####
+
+##%flag
+##% key: o
+##% description: Calculate cross point of two straights
+##% guisection: Tools
+##%end
+
+##%flag
+##% key: t
+##% description: Calculate tangencial point of a straight and a circle
+##% guisection: Tools
+##%end
+
+##%flag
+##% key: l
+##% description: along a straight
+##% guisection: Tools
+##%end
+
+##%option
+##% key: tcenter
+##% type: string
+##% description: Center of circle
+##% required: no
+##% guisection: Tools
+##%end
+
+##%option
+##% key: tradio
+##% type: double
+##% description: Radio of circle
+##% required: no
+##% guisection: Tools
+##%end
+
+##%option
+##% key: tpoint1
+##% type: string
+##% description: Point 1 of straight1
+##% required: no
+##% guisection: Tools
+##%end
+
+##%option
+##% key: tpoint2
+##% type: string
+##% description: Point 2 of straight1
+##% required: no
+##% guisection: Tools
+##%end
+
+##%option
+##% key: tpoint3
+##% type: string
+##% description: Point 1 of straight2
+##% required: no
+##% guisection: Tools
+##%end
+
+##%option
+##% key: tpoint4
+##% type: string
+##% description: Point 2 of straight2
+##% required: no
+##% guisection: Tools
+##%end
+
+##%option
+##% key: dist
+##% type: double
+##% description: Distance to along
+##% required: no
+##% guisection: Tools
+##%end
+
+#### Street section ####
+
+#%flag
+#% key: s
+#% description: Curve between two straights
+#% guisection: Street1
+#%end
+
+#%option
+#% key: sradio
+#% type: double
+#% description: Radio of circle for intersection
+#% required: no
+#% guisection: Street1
+#%end
+
+#%option
+#% key: spoint1
+#% type: string
+#% description: Point 1 of straight1
+#% required: no
+#% guisection: Street1
+#%end
+
+#%option
+#% key: spoint2
+#% type: string
+#% description: Point 2 of straight1
+#% required: no
+#% guisection: Street1
+#%end
+
+#%option
+#% key: lado1
+#% type: string
+#% label: left or right side
+#% options: Izq,Der
+#% required: no
+#% guisection: Street1
+#%end
+
+#%option
+#% key: dist_despl1
+#% type: double
+#% description: Distance displaced
+#% required: no
+#% guisection: Street1
+#%end
+
+#%option
+#% key: pkref1
+#% type: double
+#% description: Fist point pk's of straight1
+#% required: no
+#% answer: 0
+#% guisection: Street1
+#%end
+
+#%option
+#% key: spoint3
+#% type: string
+#% description: Point 1 of straight2
+#% required: no
+#% guisection: Street1
+#%end
+
+#%option
+#% key: spoint4
+#% type: string
+#% description: Point 2 of straight2
+#% required: no
+#% guisection: Street1
+#%end
+
+#%option
+#% key: lado2
+#% type: string
+#% label: left or right side
+#% options: Izq,Der
+#% required: no
+#% guisection: Street1
+#%end
+
+#%option
+#% key: dist_despl2
+#% type: double
+#% description: Distance displaced
+#% required: no
+#% guisection: Street1
+#%end
+
+#%option
+#% key: pkref2
+#% type: double
+#% description: Fist point pk's of straight2
+#% required: no
+#% answer: 0
+#% guisection: Street1
+#%end
+
+#### Street2 section ####
+
+#%flag
+#% key: h
+#% description: Curve between a straight and a roundabout
+#% guisection: Street2
+#%end
+
+#%option
+#% key: inout
+#% type: string
+#% label: In or out from roundabout
+#% options: In,Out
+#% required: no
+#% guisection: Street2
+#%end
+
+#%option
+#% key: hradio
+#% type: double
+#% description: Radio of circle for intersection
+#% required: no
+#% guisection: Street2
+#%end
+
+#%option
+#% key: hpoint1
+#% type: string
+#% description: Point 1 of straight1
+#% required: no
+#% guisection: Street2
+#%end
+
+#%option
+#% key: hpoint2
+#% type: string
+#% description: Point 2 of straight1
+#% required: no
+#% guisection: Street2
+#%end
+
+#%option
+#% key: hlado
+#% type: string
+#% label: left or right side
+#% options: Izq,Der
+#% required: no
+#% guisection: Street2
+#%end
+
+#%option
+#% key: hdist_despl
+#% type: double
+#% description: Distance displaced
+#% required: no
+#% guisection: Street2
+#%end
+
+#%option
+#% key: hpkref
+#% type: double
+#% description: Fist point pk's of straight1
+#% required: no
+#% answer: 0
+#% guisection: Street2
+#%end
+
+#%option
+#% key: hcenter
+#% type: string
+#% description: Center of circle
+#% required: no
+#% guisection: Street2
+#%end
+
+#%option
+#% key: hradio1
+#% type: double
+#% description: Radio of circle
+#% required: no
+#% guisection: Street2
+#%end
+
+
+
+import os, sys
+from math import *
+import grass.script as grass
+
+
+
+###############################################################################
+
+def splitdlines(lines):
+
+ # line = [[],[],...]
+ #[[line11,[],line12,...],line3,line4,...] --> [line11,line12,line3,line4,...]
+ # --> [ 1 , 1 , 2 , 3 ,...]
+ tolines,cats=[],[]
+ for j,line in enumerate(lines):
+ if [] in line:
+ splitlist=[list(group) for k, group in groupby(line, lambda x: x == []) if not k] # split list
+ tolines.extend(splitlist)
+ for i in range(len(splitlist)):
+ cats.append(j+1)
+ else:
+ tolines.append(line)
+ cats.append(j+1)
+
+ return tolines,cats
+
+def write_Polylines(lines,name):
+
+ sal_linea=""
+
+ tolines,cats=splitdlines(lines)
+
+ for j,line in enumerate(tolines):
+
+ sal_linea+="L "+str(len(line))+" 1\n"
+ for i,pto in enumerate(line):
+ sal_linea+=str(pto[0])+" "+str(pto[1])+" "+str(pto[2])+"\n"
+ sal_linea+="1 "+str(cats[j])+"\n"
+ #print sal_linea
+ grass.write_command('v.in.ascii', flags='nz', output=name, stdin=sal_linea,
+ input='-', format='standard', overwrite=True, quiet=True)
+ return 0
+
+def write_Polyline(puntos,name):
+
+ # Write Polyline
+ sal_linea="L "+str(len(puntos))+" 1\n"
+ for i,pp in enumerate(puntos):
+ sal_linea+=str(pp[0])+" "+str(pp[1])+" "+str(pp[2])+"\n"
+ sal_linea+="1 1\n"
+ grass.write_command('v.in.ascii', flags='nz', output=name, stdin=sal_linea,
+ input='-', format='standard', overwrite=True, quiet=True)
+ return 0
+
+def write_edge(polygon,radios,A1,A2,outmap):
+
+ write_Polyline(polygon,outmap)
+
+ #grass.run_command('v.db.addtable', map=outmap, quiet=True)
+ grass.run_command('v.civil.road', flags='n', edge=outmap)
+
+ f = grass.vector_db(outmap)[1]
+ table = f['table']
+ database1 = f['database']
+ driver = f['driver']
+
+ sql = ""
+ for i,ra in enumerate(radios):
+ sql += "UPDATE "+outmap+"_Plan SET "
+ sql += "radio="+str(ra)+", a_in="+str(A1[i])+", a_out="+str(A2[i])
+ sql += " WHERE cat2="+str(i+1)+" ;\n"
+
+ grass.write_command('db.execute', database = database1, driver = driver, stdin = sql, input='-', quiet=True)
+ grass.run_command('v.civil.road', flags='hyu', edge=outmap)
+
+ return 0
+
+###############################################################################
+
+def aprox_coord(L, Tau):
+
+ n_iter=10;x=0;y=0
+ for n in range(n_iter):
+ x+=(((-1)**n*Tau**(2*n))/((4*n+1)*factorial(2*n)))
+ y+=(((-1)**n*Tau**(2*n+1))/((4*n+3)*factorial(2*n+1)))
+ x=x*L
+ y=y*L
+ return [x,y]
+
+def clotoide_Locales(A,R):
+
+ if R==0:
+ L,Tau=0,0
+ else:
+ L=A**2/abs(R)
+ Tau=L/(2*R)
+ xe,ye=aprox_coord(L,Tau)
+ xo=xe-R*sin(Tau)
+ yo=ye+R*cos(Tau)-R
+ return xo,yo,Tau,L,xe,ye
+
+def azimut(a,b,c,d):
+
+ if a > c and b == d: Az=3*pi/2
+ elif a < c and b == d: Az=pi/2
+ elif a == c and b > d: Az=pi
+ elif a == c and b < d: Az=2*pi
+ elif a == c and b == d: Az=0
+ else: Az=atan((c-a)/(d-b))
+ #if a > c and b == d: Az=Az+pi # Az>0 -> Az > 0
+ if a < c and b > d: Az=Az+pi # Az<0 -> Az > 100
+ elif a > c and b > d: Az=Az+pi # Az>0 -> Az > 200
+ elif a > c and b < d: Az=Az+2*pi # Az<0 -> Az > 300
+ return Az
+
+def cambio_coord(x,y,Az,c,d,R):
+
+ if R>0:
+ x1=c-x*sin(Az)+y*cos(Az)
+ y1=d-x*cos(Az)-y*sin(Az)
+ else:
+ x1=c-x*sin(Az)+y*cos(Az)
+ y1=d-x*cos(Az)-y*sin(Az)
+ return [x1,y1]
+
+
+
+def recta_tg_circulo(xo,yo,R,xc,yc):
+
+ a=(R**2-xc**2+2*xc*xo-xo**2)
+ b=2*(xc*yc-xc*yo-xo*yc+xo*yo)
+ c=(R**2-yc**2+2*yc*yo-yo**2)
+
+ m1=(-b+sqrt(b**2-4*a*c))/(2*a)
+ m2=(-b-sqrt(b**2-4*a*c))/(2*a)
+
+ x1=(xc+xo*m1**2+yc*m1-yo*m1)/(m1**2+1)
+ y1=-1/m1*(x1-xc)+yc
+
+ x2=(xc+xo*m2**2+yc*m2-yo*m2)/(m2**2+1)
+ y2=m2*(x2-xo)+yo
+
+ if options['izq']=='0':
+ xt,yt=x1,y1
+ else: xt,yt=x2,y2
+
+ return xt,yt
+
+
+def pto_corte_2_rectas(x1,y1,x2,y2,x3,y3,x4,y4):
+
+ if x2 == x1: m11=(y2-y1)/0.0000001
+ else: m11=(y2-y1)/(x2-x1)
+
+ if x3 == x4: m22=(y4-y3)/0.0000001
+ else: m22=(y4-y3)/(x4-x3)
+
+ x=(m11*x1-m22*x3-y1+y3)/(m11-m22)
+ y=m11*(x-x1)+y1
+
+ return x,y
+
+
+def alargar_recta(x1,y1,x2,y2,d):
+
+ Az=azimut(x1,y1,x2,y2)
+
+ x=x2+d*sin(Az)
+ y=y2+d*cos(Az)
+
+ return x,y,x-x2,y-y2
+
+def angulos_Alings(a,b,c,d,e,f):
+
+ Az_ent = azimut(a,b,c,d)
+ Az_sal = azimut(c,d,e,f)
+
+ if ((a <= c and b <= d) or (a <= c and b >= d)): #1er cuadrante o 2do cuadrante
+
+ if (Az_ent < Az_sal and Az_ent+pi < Az_sal):
+ w=2*pi-abs(Az_ent-Az_sal)
+ else:
+ w=abs(Az_ent-Az_sal)
+
+ if ((a >= c and b >= d) or (a >= c and b <= d)): #3er cuadrante o 4to cuadrante
+
+ if (Az_ent > Az_sal and Az_ent-pi > Az_sal):
+ w=2*pi-abs(Az_ent-Az_sal)
+ else:
+ w=abs(Az_ent-Az_sal)
+
+ return Az_ent,Az_sal,w
+
+###############################################################################
+
+def curva_recta(center,radio,param,izq,point):
+
+ xc,yc=center.split(',')
+ xc=float(xc)
+ yc=float(yc)
+ R=float(radio)
+ A=float(param)
+ xo,yo=point.split(',')
+ xo,yo=float(xo),float(yo)
+
+ xo_ent,yo_ent,Tau_ent,Lent,xe,ye=clotoide_Locales(A,R)
+ #print xo_ent,yo_ent,Tau_ent,Lent,xe,ye,A
+ R1=R+yo_ent
+
+ xt,yt=recta_tg_circulo(xo,yo,R1,xc,yc)
+
+ xl=(R+yo_ent)*tan(Tau_ent)
+ d=sqrt((xo-xt)**2+(yo-yt)**2)
+ #print xl,d
+ Az=azimut(xo,yo,xt,yt)
+
+ xv=xo+(d+xl)*sin(Az)
+ yv=yo+(d+xl)*cos(Az)
+
+ x3=xv+(d+xl)*sin(Az-Tau_ent*2)
+ y3=yv+(d+xl)*cos(Az-Tau_ent*2)
+
+
+ write_edge([[x3,y3,0],[xv,yv,0],[xo,yo,0]],[0,R,0],[0,A,0],[0,A,0],options['outmap'].split('@')[0])
+
+ write_Polyline([[xc,yc,0],[xt,yt,0]],'acc')
+ write_Polyline([[xc,yc,0],[xv,yv,0]],'avv')
+
+ return 0
+
+###############################################################################
+
+###############################################################################
+def get_pts_Clot_ent(A,xad,yad,Az,Ini,Fin,L_int,R,Lacum,cat,ali):
+ # Return matriz de puntos, resto
+ M = []
+ while Ini <= Fin:
+ Rclo=A**2/(Ini)
+ Tau_clo=(Ini)/(2*Rclo)
+ xclo,yclo=aprox_coord((Ini),Tau_clo)
+ if R>0:
+ x_clo=xad-xclo*sin(-Az)+yclo*cos(-Az)
+ y_clo=yad+xclo*cos(-Az)+yclo*sin(-Az)
+ Az1=Az+Tau_clo
+ elif R<0:
+ x_clo=xad+xclo*sin(Az)-yclo*cos(Az)
+ y_clo=yad+xclo*cos(Az)+yclo*sin(Az)
+ Az1=Az-Tau_clo
+ Lacum+=L_int
+ M.append([x_clo,y_clo,0,cat,Lacum,Az1,'Clot_in',ali])
+ cat=cat+1
+ Ini=Ini+L_int
+ return M,Fin-(Ini-L_int),Lacum,cat
+
+def get_pts_Clot_sal(A,xda,yda,Az,Ini,Fin,L_int,R,Lacum,cat,ali):
+ # Return matriz de puntos, resto
+ M = []
+ while Ini <= Fin:
+ if Fin-Ini==0: Fin=0; Ini=Ini*-1
+ Rclo=A**2/(Fin-Ini)
+ Tau_clo=(Fin-Ini)/(2*Rclo)
+ xclo,yclo=aprox_coord((Fin-Ini),Tau_clo)
+ if R>0:
+ x_clo=xda-xclo*sin(Az)+yclo*cos(Az)
+ y_clo=yda-xclo*cos(Az)-yclo*sin(Az)
+ Az1=Az-Tau_clo
+ elif R<0:
+ x_clo=xda+xclo*sin(-Az)-yclo*cos(-Az)
+ y_clo=yda-xclo*cos(-Az)-yclo*sin(-Az)
+ Az1=Az+Tau_clo
+ Lacum+=L_int
+ M.append([x_clo,y_clo,0,cat,Lacum,Az1,'Clot_out',ali])
+ cat=cat+1
+ Ini=Ini+L_int
+ return M,Fin-(Ini-L_int),Lacum,cat
+
+def get_pts_Curva(xc,yc,AzIni,AzFin,incremento,R,Lacum,dc,cat,ali):
+
+ M = []
+ incremento=incremento/abs(R)
+ if R>0:
+ if AzIni>AzFin: AzIni=AzIni-2*pi
+ ii = AzIni+incremento
+ while ii <= AzFin:
+ x1=xc+R*sin(ii)
+ y1=yc+R*cos(ii)
+ Az1=(ii+pi/2)
+ ii += incremento
+ Lacum+=incremento*abs(R)
+ if Az1>2*pi: Az1=Az1-2*pi
+ M.append([x1,y1,0,cat,Lacum,Az1,'Curve',ali])
+ cat=cat+1
+ resto=(AzFin-(ii-incremento))*abs(R)
+ elif R<0:
+ if AzIni<AzFin: AzFin=AzFin-2*pi
+ ii = AzIni-incremento
+ while ii >= AzFin:
+ x1=xc-R*sin(ii)
+ y1=yc-R*cos(ii)
+ Az1=(ii-pi/2)
+ ii -= incremento
+ Lacum+=incremento*abs(R)
+ if Az1<0: Az1=Az1+2*pi
+ M.append([x1,y1,0,cat,Lacum,Az1,'Curve',ali])
+ cat=cat+1
+ resto= ((ii+incremento)-AzFin)*abs(R)
+ return M,resto,Lacum,cat
+
+def get_pts_Recta(xo,yo,zo,Az,Ini,Fin,L_int,Lacum,cat,ali):
+ # Return matriz de puntos, resto
+ if L_int==0: L_int=1
+ M=[]
+ while Ini <= Fin:
+ x2=xo+Ini*sin(Az)
+ y2=yo+Ini*cos(Az)
+ Lacum+=L_int
+ M.append([x2,y2,zo,cat,Lacum,Az,'Line',ali])
+ cat=cat+1
+ Ini+=L_int
+ return M,Fin-(Ini-L_int),Lacum,cat
+
+
+
+
+###############################################################################
+
+def curva_R_r(xc1,yc1,R1,xp1,yp1,As1,R2,Ae2,As2,alpha,h,Lacum):
+
+ # R1>R2
+ #Calculamos los puntos de tangencia de la clotoide con los dos circulos
+ xo1_s,yo1_s,Tau1_s,Ls1,xe1,ye1=clotoide_Locales(As1,R1)
+ xo2_e,yo2_e,Tau2_e,Le2,xe2,ye2=clotoide_Locales(Ae2,R2)
+ Az1=azimut(xc1,yc1,xp1,yp1)
+
+ if R1 > 0 and R2 > 0:
+ g90 = pi/2
+ g180 = pi
+
+ elif R1 < 0 and R2 < 0:
+ alpha=-alpha
+ g90 = -pi/2
+ g180 = -pi
+
+ else:
+ print "Error: For change the radio sing a straight must be between the radios"
+ return 0
+
+ xt1=xc1+abs(R1+yo1_s)*sin(Az1-Tau1_s)
+ yt1=yc1+abs(R1+yo1_s)*cos(Az1-Tau1_s)
+
+ xt2=xt1+(xo2_e-xo1_s)*sin((Az1-Tau1_s)+g90)
+ yt2=yt1+(xo2_e-xo1_s)*cos((Az1-Tau1_s)+g90)
+
+ xad2=xt1-(xo1_s)*sin((Az1-Tau1_s)+g90)
+ yad2=yt1-(xo1_s)*cos((Az1-Tau1_s)+g90)
+
+ xc2=xt2+abs(R2+yo2_e)*sin((Az1-Tau1_s)+g180)
+ yc2=yt2+abs(R2+yo2_e)*cos((Az1-Tau1_s)+g180)
+
+ xar2=xc2+abs(R2)*sin(Az1-Tau1_s+Tau2_e)
+ yar2=yc2+abs(R2)*cos(Az1-Tau1_s+Tau2_e)
+
+ xra2=xc2+abs(R2)*sin((Az1-Tau1_s)+(Tau2_e+alpha))
+ yra2=yc2+abs(R2)*cos((Az1-Tau1_s)+(Tau2_e+alpha))
+
+
+ Azent=azimut(xad2,yad2,xt2,yt2)
+ Azini=azimut(xc2,yc2,xar2,yar2)
+ Azfin=azimut(xc2,yc2,xra2,yra2)
+
+ aligns = [["C_ent",[Ae2,Azent,xad2,yad2,xar2,yar2,Le2,Lacum,R2,Le2-Ls1]]]
+ aligns.append(["Curva",[R2,alpha,xc2,yc2,Azini,Azfin,alpha*abs(R2),Lacum]])
+
+ return xc2,yc2,R2,xra2,yra2,[xad2,yad2,xt1,yt1],Lacum,aligns
+
+####-------------
+
+def curva_r_R(xc1,yc1,R1,xp1,yp1,As1,R2,Ae2,As2,alpha,h,Lacum):
+
+ # R1<R2
+ xo1_s,yo1_s,Tau1_s,Ls1,xe1,ye1=clotoide_Locales(As1,R1)
+ xo2_e,yo2_e,Tau2_e,Le2,xe2,ye2=clotoide_Locales(Ae2,R2)
+ Az1=azimut(xc1,yc1,xp1,yp1)
+
+ if R1 > 0 and R2 > 0:
+ g90 = pi/2
+ g180 = pi
+
+ elif R1 < 0 and R2 < 0:
+ alpha=-alpha
+ g90 = -pi/2
+ g180 = -pi
+
+ else:
+ print "Error: For change the radio sing a straight must be between the radios"
+ return 0
+
+ xt1=xc1+abs(R1+yo1_s)*sin(Az1+Tau1_s)
+ yt1=yc1+abs(R1+yo1_s)*cos(Az1+Tau1_s)
+
+ xt2=xt1+(xo1_s-xo2_e)*sin((Az1+Tau1_s)+g90)
+ yt2=yt1+(xo1_s-xo2_e)*cos((Az1+Tau1_s)+g90)
+
+ xda1=xt1+(xo1_s)*sin((Az1+Tau1_s)+g90)
+ yda1=yt1+(xo1_s)*cos((Az1+Tau1_s)+g90)
+
+ xc2=xt2+abs(R2+yo2_e)*sin((Az1+Tau1_s)+g180)
+ yc2=yt2+abs(R2+yo2_e)*cos((Az1+Tau1_s)+g180)
+
+ xar2=xc2+abs(R2)*sin(Az1+Tau1_s-Tau2_e)
+ yar2=yc2+abs(R2)*cos(Az1+Tau1_s-Tau2_e)
+
+ xra2=xc2+abs(R2)*sin(Az1+Tau1_s-Tau2_e+alpha)
+ yra2=yc2+abs(R2)*cos(Az1+Tau1_s-Tau2_e+alpha)
+
+ Azsal=azimut(xt1,yt1,xda1,yda1)
+ Azini=azimut(xc2,yc2,xar2,yar2)
+ Azfin=azimut(xc2,yc2,xra2,yra2)
+
+ aligns = [["C_sal",[As1,Azsal,xda1,yda1,xp1,yp1,Ls1,Lacum,R1,Ls1-Le2]]]
+ aligns.append(["Curva",[R2,alpha,xc2,yc2,Azini,Azfin,alpha*abs(R2)]])
+
+ return xc2,yc2,R2,xra2,yra2,[xt1,yt1,xda1,yda1],Lacum,aligns
+
+####-------------
+
+def curva_recta2(xc1,yc1,R1,xp1,yp1,As1,R2,Ae2,As2,len2,h,Lacum):
+
+ # R1 != 0 y R2 = 0
+ xo1_s,yo1_s,Tau1_s,Ls1,xe1,ye1=clotoide_Locales(As1,R1)
+ Az1=azimut(xc1,yc1,xp1,yp1)
+
+ xra1=xp1
+ yra1=yp1
+
+ if R1 > 0:
+ g90 = pi/2
+ else:
+ g90 = -pi/2
+
+ xt1=xc1+abs(R1+yo1_s)*sin(Az1+Tau1_s)
+ yt1=yc1+abs(R1+yo1_s)*cos(Az1+Tau1_s)
+
+ xda1=xt1+(xo1_s)*sin((Az1+Tau1_s)+g90)
+ yda1=yt1+(xo1_s)*cos((Az1+Tau1_s)+g90)
+
+ if As1 == 0:
+
+ xda1=xra1
+ yda1=yra1
+
+ xp1=xda1+len2*sin((Az1+Tau1_s)+g90)
+ yp1=yda1+len2*cos((Az1+Tau1_s)+g90)
+
+
+ Az1=azimut(xt1,yt1,xp1,yp1)
+ dist=sqrt((xda1-xp1)**2+(yda1-yp1)**2)
+
+ aligns = [["C_sal",[As1,Az1,xda1,yda1,xra1,yra1,Ls1,Lacum,R2]]]
+ aligns.append(["Recta",[xp1,yp1,0,dist,Lacum,Az1]])
+
+ return xda1,yda1,R2,xp1,yp1,[xt1,yt1,xp1,yp1],Lacum,aligns
+
+####-------------
+
+def recta_curva(xc1,yc1,R1,xp1,yp1,As1,R2,Ae2,As2,alpha,h,Lacum):
+
+ # R1 = 0 y R2 != 0
+ xo2_e,yo2_e,Tau2_e,Le2,xe2,ye2=clotoide_Locales(Ae2,R2)
+
+ Az0=azimut(xc1,yc1,xp1,yp1)
+ xad2=xp1
+ yad2=yp1
+
+ xt2=xad2+(xo2_e)*sin(Az0)
+ yt2=yad2+(xo2_e)*cos(Az0)
+
+ if R2 > 0:
+ g90 = pi/2
+
+ else:
+ alpha=-alpha
+ g90 = -pi/2
+
+ xc2=xt2+abs(R2+yo2_e)*sin(Az0+g90)
+ yc2=yt2+abs(R2+yo2_e)*cos(Az0+g90)
+
+ Az1=Az0-g90
+
+ xar2=xc2+abs(R2)*sin(Az1+Tau2_e)
+ yar2=yc2+abs(R2)*cos(Az1+Tau2_e)
+
+ if Ae2 == 0:
+
+ xar2=xad2
+ yar2=yad2
+
+ xra2=xc2+abs(R2)*sin(Az1+Tau2_e+alpha)
+ yra2=yc2+abs(R2)*cos(Az1+Tau2_e+alpha)
+
+
+ aligns = [["C_ent",[Ae2,Az0,xad2,yad2,xar2,yar2,Le2,Lacum,R2]]]
+ aligns.append(["Curva",[R2,alpha,xc2,yc2,Az1,Az1+alpha,alpha*abs(R2),Lacum]])
+
+ return xc2,yc2,R2,xra2,yra2,[xc1,yc1,xp1,yp1],Lacum,aligns
+
+####-------------
+
+def generate_edge_polygon(center,radio1,param1,point1,radio2,param21,param22,len2):
+
+ xc1,yc1=center.split(',')
+ xc1,yc1=float(xc1),float(yc1)
+ R1=float(radio1)
+ As1=float(param1)
+
+ xp1,yp1=point1.split(',')
+ xp1,yp1=float(xp1),float(yp1)
+
+ R2 = [float(p) for p in radio2.split(',')]+[0,0]
+ Ae2 = [float(p) for p in param21.split(',')]+[0,0]
+ As2 = [float(p) for p in param22.split(',')]+[0,0]
+ len2 = [float(p) for p in len2.split(',')]+[0,0]
+
+ Lacum=0
+ if R1 == 0:
+ dist=sqrt((xc1-xp1)**2+(yc1-yp1)**2)
+ Lacum=Lacum+dist
+
+ else:
+ # corregimos xp1
+ Az1=azimut(xc1,yc1,xp1,yp1)
+ xp1=xc1+(R1)*sin(Az1)
+ yp1=yc1+(R1)*cos(Az1)
+
+ h=1
+
+ tabla=[[R1,0,As1]]
+ ppolyg,peje,aligns=[],[],[]
+ for i,ra in enumerate(R2[:-1]):
+
+ if h == 3: h=1
+
+ if R1 == 0 and R2[i] == 0:
+ #print "R1=0 and R2=0"
+ continue
+
+ elif R1 == 0:
+ xc1,yc1,R1,xp1,yp1,pts_polyg,Lacum,align=recta_curva(xc1,yc1,R1,xp1,yp1,As1,R2[i],Ae2[i],As2[i],len2[i]/abs(R2[i]),h,Lacum)
+ #write_Polyline(ptos_eje,options['outmap']+'recta1')
+
+ elif R2[i] == 0:
+ xc1,yc1,R1,xp1,yp1,pts_polyg,Lacum,align=curva_recta2(xc1,yc1,R1,xp1,yp1,As1,R2[i],Ae2[i],As2[i],len2[i],h,Lacum)
+ #write_Polyline(ptos_eje,options['outmap']+'recta2')
+
+ elif abs(R1) > abs(R2[i]) :
+ xc1,yc1,R1,xp1,yp1,pts_polyg,Lacum,align=curva_R_r(xc1,yc1,R1,xp1,yp1,As1,R2[i],Ae2[i],As2[i],len2[i]/abs(R2[i]),h,Lacum)
+ #write_Polyline(ptos_eje,options['outmap']+'recta3')
+
+ elif abs(R1) < abs(R2[i]):
+ xc1,yc1,R1,xp1,yp1,pts_polyg,Lacum,align=curva_r_R(xc1,yc1,R1,xp1,yp1,As1,R2[i],Ae2[i],As2[i],len2[i]/abs(R2[i]),h,Lacum)
+ #write_Polyline(ptos_eje,options['outmap']+'recta4')
+
+ h=h+1
+
+ aligns.extend(align)
+ ppolyg.append(pts_polyg)
+
+ if R2[i] != 0:
+ tabla.append([R2[i],Ae2[i],As2[i]])
+
+ #peje.append(ptos_eje)
+
+ As1=As2[i]
+
+ tabla.append([0,0,0])
+
+ polygon = [ppolyg[0][:2]+[0]]
+
+ for i,pto in enumerate(ppolyg[:-1]):
+
+ Az_1 = azimut(ppolyg[i][0],ppolyg[i][1],ppolyg[i][2],ppolyg[i][3])
+ Az_2 = azimut(ppolyg[i+1][0],ppolyg[i+1][1],ppolyg[i+1][2],ppolyg[i+1][3])
+
+ if round(Az_1,8) == round(Az_2,8): continue
+
+ xv,yv=pto_corte_2_rectas(ppolyg[i][0],ppolyg[i][1],ppolyg[i][2],ppolyg[i][3],ppolyg[i+1][0],ppolyg[i+1][1],ppolyg[i+1][2],ppolyg[i+1][3])
+ polygon.append([xv,yv,0])
+
+ polygon.append(ppolyg[-1][2:]+[0])
+
+ #for jj in aligns: print jj
+ return polygon,tabla
+
+###############################################################################
+
+
+
+
+
+###############################################################################
+
+def generate_Aligns(poly,tabla):
+
+ LAcum=0
+ #for jj in tabla: print jj
+ x_ini,y_ini = poly[0][0],poly[0][1]
+ puntos_eje=[]
+ centros = []
+ for i in range(1,len(poly)-1,1):
+
+ a,b = poly[i-1][0],poly[i-1][1]
+ c,d = poly[i][0],poly[i][1]
+ e,f = poly[i+1][0],poly[i+1][1]
+
+ Az_ent,Az_sal,w=angulos_Alings(a,b,c,d,e,f)
+
+ R,Ae,As = tabla[i][0],tabla[i][1],tabla[i][2]
+
+ if R == 0:
+ Lrecta=sqrt((c-x_ini)**2+(d-y_ini)**2)
+ LAcum+=Lrecta
+ Lini_e,Lini_s=0,0
+ puntos_eje.append([[x_ini,y_ini,0,Lrecta,LAcum,Az_ent],
+ [0,Az_ent,c,d,c,d,0,LAcum,0,Lini_e,0,0],
+ [0,0,0,0,0,0,0,LAcum],
+ [0,Az_sal,c,d,c,d,0,LAcum,0,Lini_s,0,0]])
+ x_ini,y_ini=c,d
+ continue
+
+ xo_e,yo_e,Tau_e,Le,xe,ye=clotoide_Locales(Ae,R)
+ xo_s,yo_s,Tau_s,Ls,xs,ys=clotoide_Locales(As,R)
+
+ # Centros de las curvas dados por el corte de las rectas paralelas a R+AR
+ if R > 0:
+ g90 = pi/2
+ alpha = abs(w - Tau_e - Tau_s)
+ else:
+ g90 = -pi/2
+ alpha = abs(w + Tau_e + Tau_s)
+
+ ac=a+abs(R+yo_e)*sin(Az_ent+g90)
+ bc=b+abs(R+yo_e)*cos(Az_ent+g90)
+
+ cc1=c+abs(R+yo_e)*sin(Az_ent+g90)
+ dc1=d+abs(R+yo_e)*cos(Az_ent+g90)
+
+ cc2=c+abs(R+yo_s)*sin(Az_sal+g90)
+ dc2=d+abs(R+yo_s)*cos(Az_sal+g90)
+
+ ec=e+abs(R+yo_s)*sin(Az_sal+g90)
+ fc=f+abs(R+yo_s)*cos(Az_sal+g90)
+
+ xc,yc=pto_corte_2_rectas(ac,bc,cc1,dc1,cc2,dc2,ec,fc)
+ centros.append([xc,yc,0])
+
+ Lini_e,Lini_s=0,0
+ if Ae <= 0:
+
+ xar = xc+abs(R)*sin(Az_ent-g90-Tau_e)
+ yar = yc+abs(R)*cos(Az_ent-g90-Tau_e)
+ xad,yad = xar,yar
+
+ Lrecta = sqrt((xad-x_ini)**2+(yad-y_ini)**2)
+ if Ae < 0:
+ x_ini,y_ini = xad,yad
+ Lrecta = 0
+ Le = 0
+ xadp,yadp = xad,yad
+
+ else:
+
+ xt1= xc+abs(R+yo_e)*sin(Az_ent-g90)
+ yt1= yc+abs(R+yo_e)*cos(Az_ent-g90)
+
+ xad = xt1+xo_e*sin(Az_ent+pi)
+ yad = yt1+xo_e*cos(Az_ent+pi)
+
+ xar = xc+abs(R)*sin(Az_ent-g90+Tau_e)
+ yar = yc+abs(R)*cos(Az_ent-g90+Tau_e)
+
+ Lrecta = sqrt((xad-x_ini)**2+(yad-y_ini)**2)
+ R1,Ae1,As1 = tabla[i-1][0],tabla[i-1][1],tabla[i-1][2]
+ if As1 < 0:
+ xo_s1,yo_s1,Tau_s1,Ls1,xs1,ys1=clotoide_Locales(As1,R1)
+ Lini_e = Ls1
+ Ptos_ent,resto,acum,cat=get_pts_Clot_ent(Ae,xad,yad,Az_ent,Lini_e,Lini_e,1,R,0,0,0)
+ xadp,yadp = Ptos_ent[0][0],Ptos_ent[0][1]
+ x_ini,y_ini = xadp,yadp
+ Lrecta = 0
+
+ Dc = abs(R)*alpha
+
+ if As <= 0:
+
+ xra = xc+abs(R)*sin(Az_sal-g90+Tau_s)
+ yra = yc+abs(R)*cos(Az_sal-g90+Tau_s)
+ xda,yda = xra,yra
+ Ls = 0
+ xdap,ydap = xda,yda
+
+ else:
+ xt1= xc+abs(R+yo_s)*sin(Az_sal-g90)
+ yt1= yc+abs(R+yo_s)*cos(Az_sal-g90)
+
+ xda = xt1+xo_s*sin(Az_sal)
+ yda = yt1+xo_s*cos(Az_sal)
+
+ xra = xc+abs(R)*sin(Az_sal-g90-Tau_s)
+ yra = yc+abs(R)*cos(Az_sal-g90-Tau_s)
+
+ R2,Ae2,As2 = tabla[i+1][0],tabla[i+1][1],tabla[i+1][2]
+ if Ae2 < 0:
+ xo_e2,yo_e2,Tau_e2,Le2,xe2,ye2=clotoide_Locales(Ae2,R2)
+ Lini_s = Le2
+ Ptos_sal,resto,acum,cat_s=get_pts_Clot_sal(As,xda,yda,Az_sal,Lini_s,Lini_s,1,R,0,0,0)
+ xdap,ydap = Ptos_sal[0][0],Ptos_sal[0][1]
+
+ #Az_ini = Az_ent-g90+Tau_e
+ #Az_fin = Az_sal-g90-Tau_s
+ Az_ini=azimut(xc,yc,xar,yar)
+ Az_fin=azimut(xc,yc,xra,yra)
+
+ puntos_eje.append([[x_ini,y_ini,0,Lrecta,LAcum+Lrecta,Az_ent],
+ [Ae,Az_ent,xad,yad,xar,yar,Le,LAcum+Lrecta+(Le-Lini_e),R,Lini_e],
+ [R,alpha,xc,yc,Az_ini,Az_fin,Dc,LAcum+Lrecta+(Le-Lini_e)+Dc],
+ [As,Az_sal,xda,yda,xra,yra,Ls,LAcum+Lrecta+(Le-Lini_e)+Dc+(Ls-Lini_s),R,Lini_s]])
+
+ x_ini,y_ini = xda,yda
+ LAcum = LAcum+Lrecta+(Le-Lini_e)+Dc+(Ls-Lini_s)
+
+ Lrecta = sqrt((e-x_ini)**2+(f-y_ini)**2)
+ puntos_eje.append([[x_ini,y_ini,0,Lrecta,LAcum+Lrecta,Az_sal],
+ [0,Az_sal,x_ini,y_ini,0,0,0,LAcum+Lrecta,0,Lini_e],
+ [0,0,0,0,0,0,0,LAcum+Lrecta],
+ [0,Az_sal,x_ini,y_ini,0,0,0,LAcum+Lrecta,0,Lini_s]])
+ #for jj in puntos_eje:
+ #for pp in jj: print pp
+ return puntos_eje
+
+###############################################################################
+
+
+def generate_Pts(puntos_eje,a):
+
+ #x_ini,y_ini,z_ini,Lrecta,LAcum,Az_ent=puntos_eje[0][0]
+
+ Az_ent=puntos_eje[1][1][1]
+ resto=1
+ interv=1
+ cat=1
+ puntos,seg,puntos_caract,puntos_centros=[],[],[],[]
+
+ acum=0
+ h=1
+ for i in range(0,len(puntos_eje),1):
+
+ Ptos_recta,Ptos_ent,Ptos_curva,Ptos_sal=[],[],[],[]
+
+ x_ini,y_ini,z_ini,Lrecta,LrAcum,Az_ent=puntos_eje[i][0]
+ Aent,Az_ent,xad,yad,xar,yar,Lent,LeAcum,R,Lini_e=puntos_eje[i][1]
+ R,alpha,xc,yc,Az_ini,Az_fin,Dc,LcAcum=puntos_eje[i][2]
+ Asal,Az_sal,xda,yda,xra,yra,Lsal,LsAcum,R,Lini_s=puntos_eje[i][3]
+
+ Ptos_recta,resto,acum,cat=get_pts_Recta(x_ini,y_ini,z_ini,Az_ent,interv-resto,Lrecta,interv,acum,cat,i+1)
+
+ if R != 0:
+ if Aent != 0:
+ Ptos_ent,resto,acum,cat=get_pts_Clot_ent(Aent,xad,yad,Az_ent,interv-resto+Lini_e,Lent,interv,R,acum,cat,i+1)
+
+ Ptos_curva,resto,acum,cat=get_pts_Curva(xc,yc,Az_ini-resto/R,Az_fin,interv,R,acum,Dc,cat,i+1)
+
+ if Asal != 0:
+ Ptos_sal,resto,acum,cat=get_pts_Clot_sal(Asal,xda,yda,Az_sal,interv-resto+Lini_s,Lsal,interv,R,acum,cat,i+1)
+ else:
+ Ptos_ent,Ptos_curva,Ptos_sal=[[],[],[]]
+
+ x_ini,y_ini=xda,yda
+
+ puntos.extend(Ptos_recta+Ptos_ent+Ptos_curva+Ptos_sal)
+
+ h=h+4
+
+ #x_fin,y_fin,z_fin,Lrecta,LAcum=[ float(t) for t in puntos_eje[-1][0]]
+ #puntos.append([x_fin,y_fin,z_fin,cat,LAcum,Az_sal,'End',i+2])
+
+ #for jj in puntos: print jj
+
+ return puntos
+
+###############################################################################
+###############################################################################
+###############################################################################
+###############################################################################
+
+
+
+def generate_PtsRound(radio,radio2,azimut,center):
+
+ xo,yo=center.split(',')
+ radio=float(radio)
+ radio2=float(radio2)
+ Az=float(azimut)
+
+ x1=float(xo)+radio2*sin(Az)
+ y1=float(yo)+radio2*cos(Az)
+
+ x2=x1+radio*sin(Az+pi/2)
+ y2=y1+radio*cos(Az+pi/2)
+
+ x3=x2+2*radio2*sin(Az+pi)
+ y3=y2+2*radio2*cos(Az+pi)
+
+ x4=x3+2*radio*sin(Az+3*pi/2)
+ y4=y3+2*radio*cos(Az+3*pi/2)
+
+ x5=x4+2*radio2*sin(Az)
+ y5=y4+2*radio2*cos(Az)
+
+ x6=x5+(radio)*sin(Az+pi/2)
+ y6=y5+(radio)*cos(Az+pi/2)
+
+ return [[x1,y1,0],[x2,y2,0],[x3,y3,0],[x4,y4,0],[x5,y5,0],[x6,y6,0]]
+
+###############################################################################
+
+
+def edges_intersec(R,point1,point2,izq1,dist1,point3,point4,izq2,dist2):
+
+ R=float(R)
+ dist1=float(dist1)
+ dist2=float(dist2)
+ x1,y1=point1.split(',')
+ x1,y1=float(x1),float(y1)
+ x2,y2=point2.split(',')
+ x2,y2=float(x2),float(y2)
+ x3,y3=point3.split(',')
+ x3,y3=float(x3),float(y3)
+ x4,y4=point4.split(',')
+ x4,y4=float(x4),float(y4)
+
+ Az1=azimut(x1,y1,x2,y2)
+ Az2=azimut(x3,y3,x4,y4)
+ #alpha=abs(Az1-Az2)
+
+ if izq1=="Izq": izq1=-1
+ else: izq1=1
+
+ # Recta paralela a dist1+R
+ x11=x1+(dist1+R)*sin(Az1+izq1*pi/2)
+ y11=y1+(dist1+R)*cos(Az1+izq1*pi/2)
+ x22=x2+(dist1+R)*sin(Az1+izq1*pi/2)
+ y22=y2+(dist1+R)*cos(Az1+izq1*pi/2)
+
+ # Recta paralela a dist1 para poner el vertice en el punto medio
+ xr1=x1+(dist1)*sin(Az1+izq1*pi/2)
+ yr1=y1+(dist1)*cos(Az1+izq1*pi/2)
+ xr2=x2+(dist1)*sin(Az1+izq1*pi/2)
+ yr2=y2+(dist1)*cos(Az1+izq1*pi/2)
+
+ if izq2=="Izq": izq2=-1
+ else: izq2=1
+
+ # Recta paralela a dist2+R
+ x33=x3+(dist2+R)*sin(Az2+izq2*pi/2)
+ y33=y3+(dist2+R)*cos(Az2+izq2*pi/2)
+ x44=x4+(dist2+R)*sin(Az2+izq2*pi/2)
+ y44=y4+(dist2+R)*cos(Az2+izq2*pi/2)
+
+ # Recta paralela a dist1 para poner el vertice en el punto medio
+ xr3=x3+(dist2)*sin(Az2+izq2*pi/2)
+ yr3=y3+(dist2)*cos(Az2+izq2*pi/2)
+ xr4=x4+(dist2)*sin(Az2+izq2*pi/2)
+ yr4=y4+(dist2)*cos(Az2+izq2*pi/2)
+
+ xv,yv=pto_corte_2_rectas(xr1,yr1,xr2,yr2,xr3,yr3,xr4,yr4)
+ xc,yc=pto_corte_2_rectas(x11,y11,x22,y22,x33,y33,x44,y44)
+
+ Az3=azimut(xc,yc,xv,yv)
+
+ # Punto de corte de la recta centro-vertice con el circulo
+ xp=xc+R*sin(Az3)
+ yp=yc+R*cos(Az3)
+
+ # Rectas perpendiculares de xc,yc y xp,yp al eje 1
+ xt1=xc+(dist1+R)*sin(Az1-izq1*pi/2)
+ yt1=yc+(dist1+R)*cos(Az1-izq1*pi/2)
+ xp1=xp+(dist1+R)*sin(Az1-izq1*pi/2)
+ yp1=yp+(dist1+R)*cos(Az1-izq1*pi/2)
+
+ xt2,yt2=pto_corte_2_rectas(x1,y1,x2,y2,xp,yp,xp1,yp1)
+
+ # Rectas perpendiculares de xc,yc y xp,yp al eje 2
+ xt3=xc+(dist2+R)*sin(Az2-izq2*pi/2)
+ yt3=yc+(dist2+R)*cos(Az2-izq2*pi/2)
+ xp2=xp+(dist2+R)*sin(Az2-izq2*pi/2)
+ yp2=yp+(dist2+R)*cos(Az2-izq2*pi/2)
+
+ xt4,yt4=pto_corte_2_rectas(x3,y3,x4,y4,xp,yp,xp2,yp2)
+
+ # Distancias al primer punto de cada eje
+ d1=sqrt((x1-xt1)**2+(y1-yt1)**2)+float(options['pkref1'])
+ d2=sqrt((x1-xt2)**2+(y1-yt2)**2)+float(options['pkref1'])
+ d3=sqrt((x3-xt3)**2+(y3-yt3)**2)+float(options['pkref2'])
+ d4=sqrt((x3-xt4)**2+(y3-yt4)**2)+float(options['pkref2'])
+
+ if d1 > d2:
+ dtmp=d1
+ d1=d2
+ d2=dtmp
+ dif=d2-d1
+ else:
+ dif=0
+
+ grass.message("Edge 1:")
+ grass.message(" pk1="+str(d1)+", dist="+str(dist1)+", radio: r"+str(R)+","+str(dif))
+ grass.message(" pk2="+str(d2)+", dist="+str(dist1)+", radio: ")
+
+ if d3 > d4:
+ dtmp=d3
+ d3=d4
+ d4=dtmp
+ dif=d4-d3
+ else:
+ dif=0
+
+ grass.message("Edge 2:")
+ grass.message(" pk1="+str(d3)+", dist="+str(dist2)+", radio: r"+str(R)+","+str(dif))
+ grass.message(" pk2="+str(d4)+", dist="+str(dist2)+", radio: ")
+
+ write_Polylines([[[x11,y11,0],[x22,y22,0]],
+ [[x33,y33,0],[x44,y44,0]],
+ [[xc,yc,0],[xv,yv,0]],
+ [[xc,yc,0],[xt1,yt1,0]],
+ [[xp,yp,0],[xt2,yt2,0]],
+ [[xc,yc,0],[xt3,yt3,0]],
+ [[xp,yp,0],[xt4,yt4,0]]],options['outmap'])
+
+ return 0
+
+###############################################################################
+
+def edge_circ_intersec(R,point1,point2,izq1,dist1,center,radio1):
+
+ R=float(R)
+ dist1=float(dist1)
+ x1,y1=point1.split(',')
+ x1,y1=float(x1),float(y1)
+ x2,y2=point2.split(',')
+ x2,y2=float(x2),float(y2)
+
+ xc,yc=center.split(',')
+ xc,yc=float(xc),float(yc)
+ radio1=float(radio1)
+
+ if izq1=="Izq": izq1=-1
+ else: izq1=1
+
+ Az1=azimut(x1,y1,x2,y2)
+
+ # Recta paralela a dist1+R
+ x11=x1+(dist1+R)*sin(Az1+izq1*pi/2)
+ y11=y1+(dist1+R)*cos(Az1+izq1*pi/2)
+ x22=x2+(dist1+R)*sin(Az1+izq1*pi/2)
+ y22=y2+(dist1+R)*cos(Az1+izq1*pi/2)
+
+ xc1=xc+radio1*sin(Az1-izq1*pi/2)
+ yc1=yc+radio1*cos(Az1-izq1*pi/2)
+
+ d1=sqrt((xc-xc1)**2+(yc-yc1)**2)
+
+ # Distancia del centro del circulo al eje
+ xc2,yc2=pto_corte_2_rectas(x1,y1,x2,y2,xc,yc,xc1,yc1)
+
+ d2=sqrt((xc-xc2)**2+(yc-yc2)**2)
+
+ if d1 > d2: d2=-d2
+
+ s=(R+dist1)+d2
+
+ alpha=asin(s/(R+radio1))
+
+ c1=radio1*cos(alpha)
+ c2=(R+radio1)*cos(alpha)
+
+ if options['inout']=='Out':
+
+ # Centro del circulo requerido
+ xcc=xc+(R+radio1)*sin(Az1+izq1*alpha)
+ ycc=yc+(R+radio1)*cos(Az1+izq1*alpha)
+
+ # Punto tangencia dos circulos
+ xp=xc+radio1*sin(Az1+izq1*alpha)
+ yp=yc+radio1*cos(Az1+izq1*alpha)
+
+ xt1=xc2+c2*sin(Az1)
+ yt1=yc2+c2*cos(Az1)
+
+ xt2=xc2+c1*sin(Az1)
+ yt2=yc2+c1*cos(Az1)
+
+ else:
+ xcc=xc+(R+radio1)*sin(Az1+izq1*(pi-alpha))
+ ycc=yc+(R+radio1)*cos(Az1+izq1*(pi-alpha))
+
+ # Punto tangencia dos circulos
+ xp=xc+radio1*sin(Az1+izq1*(pi-alpha))
+ yp=yc+radio1*cos(Az1+izq1*(pi-alpha))
+
+ xt1=xc2+c2*sin(Az1+pi)
+ yt1=yc2+c2*cos(Az1+pi)
+
+ xt2=xc2+c1*sin(Az1+pi)
+ yt2=yc2+c1*cos(Az1+pi)
+
+ xpp=xcc+s*sin(Az1-izq1*pi/2)
+ ypp=ycc+s*cos(Az1-izq1*pi/2)
+
+ d3=sqrt((x1-xt1)**2+(y1-yt1)**2)+float(options['hpkref'])
+ d4=sqrt((x1-xt2)**2+(y1-yt2)**2)+float(options['hpkref'])
+
+ if d3 > d4:
+ dtmp=d3
+ d3=d4
+ d4=dtmp
+ dif=d4-d3
+ else:
+ dif=0
+
+ grass.message("Edge 1:")
+ grass.message(" pk1="+str(d3)+", dist="+str(dist1)+", radio: r"+str(R)+","+str(dif))
+ grass.message(" pk2="+str(d4)+", dist="+str(dist1)+", radio: ")
+
+ write_Polylines([[[xc,yc,0],[xc1,yc1,0]],
+ [[xc,yc,0],[xcc,ycc,0]],
+ [[xcc,ycc,0],[xt1,yt1,0]],
+ [[xc,yc,0],[xpp,ypp,0]],
+ [[xp,yp,0],[xt2,yt2,0]]],options['outmap'])
+
+ return 0
+
+
+# ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### #
+# ### Main
+# ###
+
+def main():
+
+
+ #if flags['a']:
+ #curva_recta(options['acenter'],options['aradio'],options['aparam'],options['izq'],options['apoint'])
+
+ if flags['c']:
+
+ polyg,tabla=generate_edge_polygon(options['ccenter'],options['cradio1'],options['casal1'],options['cpoint1'],options['cradio2'],options['caent2'],options['casal2'],options['len2'])
+ pts=generate_Aligns(polyg,tabla)
+ list_ptos=generate_Pts(pts,1)
+
+
+ radios,A1,A2=[],[],[]
+ for p in tabla:
+ radios.append(p[0])
+ A1.append(p[1])
+ A2.append(p[2])
+ write_Polyline(list_ptos,options['outmap'].split('@')[0]+'_poly2')
+ write_edge(polyg,radios,A1,A2,options['outmap'].split('@')[0])
+
+
+
+#### Roudabout ###
+
+ if flags['h']:
+ edge_circ_intersec(options['hradio'],
+ options['hpoint1'],options['hpoint2'],options['hlado'],options['hdist_despl'],
+ options['hcenter'],options['hradio1'])
+
+ if flags['s']:
+ edges_intersec(options['sradio'],
+ options['spoint1'],options['spoint2'],options['lado1'],options['dist_despl1'],
+ options['spoint3'],options['spoint4'],options['lado2'],options['dist_despl2'])
+
+#### Roudabout ###
+
+ if flags['r']:
+ if options['rround2']=='0': options['rround2']=options['rround1']
+
+ az=float(options['azround'])
+ PtsRound=generate_PtsRound(options['rround1'],options['rround2'],az*pi/200,options['cround'])
+
+ radios,A1,A2=[0],[0],[0]
+ for p in PtsRound[1:-1]:
+ radios.append(options['rround1'])
+ A1.append(0)
+ A2.append(0)
+ write_edge(PtsRound,radios,A1,A2,options['outmap'].split('@')[0])
+ #write_Polyline(PtsRound,options['outmap'])
+
+#### Tools ###
+
+ if flags['o']:
+
+ x1,y1=options['tpoint1'].split(',')
+ x1,y1=float(x1),float(y1)
+ x2,y2=options['tpoint2'].split(',')
+ x2,y2=float(x2),float(y2)
+ x3,y3=options['tpoint3'].split(',')
+ x3,y3=float(x3),float(y3)
+ x4,y4=options['tpoint4'].split(',')
+ x4,y4=float(x4),float(y4)
+
+ xx,yy=pto_corte_2_rectas(x1,y1,x2,y2,x3,y3,x4,y4)
+ grass.message("Coord: "+str(xx)+","+str(yy))
+
+ if flags['t']:
+
+ R=float(options['tradio'])
+ x1,y1=options['tpoint1'].split(',')
+ x1,y1=float(x1),float(y1)
+ xc,yc=options['tcenter'].split(',')
+ xc,yc=float(xc),float(yc)
+
+ xx,yy=recta_tg_circulo(x1,y1,R,xc,yc)
+ grass.message("Coord: "+str(xx)+","+str(yy))
+
+ if flags['l']:
+
+ x1,y1=options['tpoint1'].split(',')
+ x1,y1=float(x1),float(y1)
+ x2,y2=options['tpoint2'].split(',')
+ x2,y2=float(x2),float(y2)
+ d=float(options['dist'])
+
+ xx,yy,dx,dy=alargar_recta(x1,y1,x2,y2,d)
+ grass.message("Coord: "+str(xx)+","+str(yy))
+ grass.message("Dist: "+str(dx)+","+str(dy))
+
+ sys.exit(0)
+
+if __name__ == "__main__":
+ options, flags = grass.parser()
+ main()
Property changes on: grass-addons/grass7/vector/v.civil/v.civil.tools/v.civil.tools.py
___________________________________________________________________
Added: svn:mime-type
+ text/x-python
Added: svn:eol-style
+ native
Added: grass-addons/grass7/vector/v.civil/v.civil.topo/Makefile
===================================================================
--- grass-addons/grass7/vector/v.civil/v.civil.topo/Makefile (rev 0)
+++ grass-addons/grass7/vector/v.civil/v.civil.topo/Makefile 2014-02-26 00:37:16 UTC (rev 59141)
@@ -0,0 +1,9 @@
+# fix this relative to include/
+# or use absolute path to the GRASS source code
+MODULE_TOPDIR = ../..
+
+PGM = v.civil.topo
+
+include $(MODULE_TOPDIR)/include/Make/Script.make
+
+default: script
Property changes on: grass-addons/grass7/vector/v.civil/v.civil.topo/Makefile
___________________________________________________________________
Added: svn:mime-type
+ text/x-makefile
Added: svn:eol-style
+ native
Added: grass-addons/grass7/vector/v.civil/v.civil.topo/v.civil.topo.html
===================================================================
--- grass-addons/grass7/vector/v.civil/v.civil.topo/v.civil.topo.html (rev 0)
+++ grass-addons/grass7/vector/v.civil/v.civil.topo/v.civil.topo.html 2014-02-26 00:37:16 UTC (rev 59141)
@@ -0,0 +1,27 @@
+<h2>DESCRIPTION</h2>
+
+<em>v.civil.topo</em> Create and modific Points-Breaklines maps, and triangulate,
+using triangle with v.triangle and nn with r.surf.nnbathy
+
+<h3>Topo</h3>
+<p>
+This section is under development and is still no usable.
+
+
+<h3>Sup</h3>
+This section follows the same philosophy as v.civil.road with the names of maps.
+<p>
+It has two option to triangulate the points and breaklines maps using triangle or nn.
+The output is a dem and countour lines delimited by the input hull.
+<p>
+The others two options patch the output maps from triangle or nn with the input dem map
+generating a new dem map and countours lines map.
+
+<h2> AUTHOR</h2>
+<p>
+Jesus Fernandez-Capel Rosillo<br>
+Civil Engineer, Spain<br>
+jfc at alcd net<br>
+
+
+
Property changes on: grass-addons/grass7/vector/v.civil/v.civil.topo/v.civil.topo.html
___________________________________________________________________
Added: svn:mime-type
+ text/html
Added: svn:keywords
+ Author Date Id
Added: svn:eol-style
+ native
Added: grass-addons/grass7/vector/v.civil/v.civil.topo/v.civil.topo.py
===================================================================
--- grass-addons/grass7/vector/v.civil/v.civil.topo/v.civil.topo.py (rev 0)
+++ grass-addons/grass7/vector/v.civil/v.civil.topo/v.civil.topo.py 2014-02-26 00:37:16 UTC (rev 59141)
@@ -0,0 +1,614 @@
+#!/usr/bin/env python
+# -*- coding: utf-8
+############################################################################
+#
+# MODULE: v.civil.topo, v0.5.0
+#
+# AUTHOR: Jesús Fernández-Capel Rosillo
+# Civil Engineer, Spain
+# jfc at alcd net
+#
+# PURPOSE: Create and modific Points-Breaklines maps, and triangulate,
+# using triangle with v.triangle and nn with r.surf.nnbathy
+#
+# COPYRIGHT: (c) 2014 Jesús Fernández-Capel Rosillo.
+#
+# This program is free software under the GNU General Public
+# License (>=v2). Read the file COPYING that comes with GRASS
+# for details.
+#
+#############################################################################
+
+#%module
+#% description: Create and modific Points-Breaklines maps, and triangulate, with triangle with v.triangle and nn with r.surf.nnbathy
+#% keywords: vector, TOPO, topography
+#%end
+
+#### Input section ####
+
+#%option
+#% key: points
+#% type: string
+#% gisprompt: old,vector,vector
+#% key_desc: name
+#% description: Name of input points vector map
+#% required: yes
+#% guisection: Input
+#%end
+
+#%option
+#% key: breaklines
+#% type: string
+#% gisprompt: old,vector,vector
+#% key_desc: name
+#% description: Name of input break lines vector map
+#% required: no
+#% guisection: Input
+#%end
+
+#%option
+#% key: contorno
+#% type: string
+#% gisprompt: old,vector,vector
+#% key_desc: name
+#% description: Name of area vector map for mask
+#% required: no
+#% guisection: Input
+#%end
+
+
+#### Topo section ####
+
+#%flag
+#% key: e
+#% description: Patch for edit, points and breaklines maps
+#% guisection: Topo
+#%end
+
+
+#%flag
+#% key: u
+#% description: Update topo map
+#% guisection: Topo
+#%end
+
+#%flag
+#% key: s
+#% description: Split maps in points and breaklines maps
+#% guisection: Topo
+#%end
+
+#%flag
+#% key: l
+#% description: Split breaklines in segments
+#% guisection: Topo
+#%end
+
+#%flag
+#% key: p
+#% description: Create points in vertices of polylines
+#% guisection: Topo
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: topomap
+#% gisprompt: new,vector,vector
+#% description: Name for output topo map
+#% required: no
+#% guisection: Topo
+#%end
+
+#%option
+#% key: cats
+#% type: string
+#% description: Cats of lines
+#% required: no
+#% guisection: Topo
+#%end
+
+#%option
+#% key: splitlength
+#% type: double
+#% description: Split segment length
+#% answer: 10
+#% required: no
+#% guisection: Topo
+#%end
+
+
+#### Sup section ####
+
+#%flag
+#% key: t
+#% description: Triangulate, create Dem and Contours using triangle
+#% guisection: Sup
+#%end
+
+#%flag
+#% key: y
+#% description: Create Dem and Contours using nnbathy
+#% guisection: Sup
+#%end
+
+#%flag
+#% key: f
+#% description: Patch tindem with dem and generate contours
+#% guisection: Sup
+#%end
+
+#%flag
+#% key: n
+#% description: Patch nnbathydem with dem and generate contours
+#% guisection: Sup
+#%end
+
+
+#%option
+#% key: step
+#% type: double
+#% description: Increment between contour levels
+#% answer: 1
+#% guisection: Sup
+#% required: no
+#%end
+
+##%option
+##% key: cut
+##% type: integer
+##% description: It acts like a filter omitting spurs single points etc
+##% answer: 0
+##% guisection: Sup
+##% required: no
+##%end
+
+#%option
+#% key: dem
+#% type: string
+#% gisprompt: old,cell,raster
+#% key_desc: raster dem
+#% description: Name of DEM raster map
+#% guisection: Sup
+#% required: no
+#%end
+
+
+#%option G_OPT_V_OUTPUT
+#% key: tin
+#% description: Name of tin map
+#% required: no
+#% answer: _Tin
+#% guisection: Sup
+#%end
+
+#%option G_OPT_R_OUTPUT
+#% key: tindem
+#% description: Name of dem map created with tin
+#% required: no
+#% answer: _TinDem
+#% guisection: Sup
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: tincontour
+#% description: Name of tin dem curves map
+#% required: no
+#% answer: _TinDemCurves
+#% guisection: Sup
+#%end
+
+#%option G_OPT_R_OUTPUT
+#% key: nndem
+#% description: Name of dem map created with nnbathy
+#% required: no
+#% answer: _NNbathyDem
+#% guisection: Sup
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: nncontour
+#% description: Name of nnbathy dem curves map
+#% required: no
+#% answer: _NNbathyDemCurves
+#% guisection: Sup
+#%end
+
+#%option G_OPT_R_OUTPUT
+#% key: tindemDem
+#% description: Name of dem map created with tin
+#% required: no
+#% answer: _TinDemDem
+#% guisection: Sup
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: tinDemcontour
+#% description: Name of tin dem curves map
+#% required: no
+#% answer: _TinDemDemCurves
+#% guisection: Sup
+#%end
+
+#%option G_OPT_R_OUTPUT
+#% key: nndemDem
+#% description: Name of dem map created with nnbathy
+#% required: no
+#% answer: _NNbathyDemDem
+#% guisection: Sup
+#%end
+
+#%option G_OPT_V_OUTPUT
+#% key: nnDemcontour
+#% description: Name of nnbathy dem curves map
+#% required: no
+#% answer: _NNbathyDemDemCurves
+#% guisection: Sup
+#%end
+
+
+import os, sys, re
+from math import *
+import grass.script as grass
+
+###############################################################################
+
+def get_pts_Recta(xo,yo,zo,Az,Ini,Fin,difcotas,L_int,Lacum,cat,ali):
+ # Return matriz de puntos, resto
+ if L_int==0: L_int=1
+ M=[]
+ Lz=0
+ alpha=atan(difcotas/Fin)
+ print Ini,Fin,difcotas,alpha
+ while Ini <= Fin:
+ x2=xo+Ini*sin(Az)
+ y2=yo+Ini*cos(Az)
+ Lacum+=L_int
+ Lz+=L_int
+ print Lz,L_int,Lacum,Ini
+ z2=zo+(Lz)*tan(alpha)
+ print zo,z2
+ M.append([x2,y2,z2,cat,Lacum,Az,'Line',ali])
+ cat=cat+1
+ Ini+=L_int
+ return M,Fin-(Ini-L_int),Lacum,cat
+
+def discr_Lines(lines,discr):
+
+ disLines=[]
+ for pto in lines:
+ resto,acum=0,0
+ cat_r=1
+ disptos=[]
+ for i in range(len(pto[:-1])):
+
+ disptos.append(pto[i])
+ Az=azimut(pto[i][0],pto[i][1],pto[i+1][0],pto[i+1][1])
+ difcotas=pto[i+1][2]-pto[i][2]
+ print difcotas
+ Lrecta=sqrt((pto[i+1][0]-pto[i][0])**2+(pto[i+1][1]-pto[i][1])**2)
+ Ptos_recta,resto,acum,cat_r=get_pts_Recta(pto[i][0],pto[i][1],pto[i][2],Az,discr-resto,Lrecta,difcotas,discr,acum,cat_r,i)
+ disptos.extend(Ptos_recta)
+ disptos.append(pto[i+1])
+ disLines.append(disptos)
+ return disLines
+
+def azimut(a,b,c,d):
+
+ if a > c and b == d: Az=3*pi/2
+ elif a < c and b == d: Az=pi/2
+ elif a == c and b > d: Az=pi
+ elif a == c and b < d: Az=2*pi
+ elif a == c and b == d: Az=0
+ else: Az=atan((c-a)/(d-b))
+ #if a > c and b == d: Az=Az+pi # Az>0 -> Az > 0
+ if a < c and b > d: Az=Az+pi # Az<0 -> Az > 100
+ elif a > c and b > d: Az=Az+pi # Az>0 -> Az > 200
+ elif a > c and b < d: Az=Az+2*pi # Az<0 -> Az > 300
+ return Az
+
+def write_Polylines(lines,name,layer,catt):
+
+ sal_linea=""
+ if layer > 1: tool="add"
+ else: tool="create"
+ tolines,cats=splitdlines(lines)
+
+ for j,line in enumerate(tolines):
+ sal_linea+="L "+str(len(line))+" 1\n"
+ for i,pto in enumerate(line):
+ sal_linea+=str(pto[0])+" "+str(pto[1])+" "+str(pto[2])+"\n"
+ sal_linea+=str(layer)+" "+str(cats[j]+catt-1)+"\n"
+
+ grass.write_command('v.edit', flags='n', tool=tool, map=name, input='-', stdin=sal_linea, overwrite=True, quiet=True)
+ return 0
+
+def splitdlines(lines):
+
+ # line = [[],[],...]
+ #[[line11,[],line12,...],line3,line4,...] --> [line11,line12,line3,line4,...]
+ # --> [ 1 , 1 , 2 , 3 ,...]
+ tolines,cats=[],[]
+ for j,line in enumerate(lines):
+ if [] in line:
+ splitlist=[list(group) for k, group in groupby(line, lambda x: x == []) if not k] # split list
+ tolines.extend(splitlist)
+ for i in range(len(splitlist)):
+ cats.append(j+1)
+ else:
+ tolines.append(line)
+ cats.append(j+1)
+
+ return tolines,cats
+
+# ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### #
+
+def generate_Segments(splitlength,catts,tmpMap):
+
+ lines=grass.read_command('v.out.ascii', input=tmpMap, output='-', type='line', cats=catts, format='standar', layer=2, quiet=True)
+ lines = lines.split('L ')
+
+ lines = [d.splitlines(0) for d in lines[1:]]
+ for i,line in enumerate(lines):
+ lines[i]=[p.split() for p in line[1:-1]]
+ for i,line in enumerate(lines):
+ for j,pto in enumerate(line):
+ lines[i][j]=[float(p) for p in pto]
+
+ lines = discr_Lines(lines,float(splitlength))
+
+ grass.run_command('v.edit', map=tmpMap, layer=2, tool='delete', where='cat2='+str(catts), quiet=True)
+ #grass.run_command('v.edit', map=tmpMap, layer=2, tool='catdel', where='cat2='+str(cats))
+ write_Polylines(lines,tmpMap,2,int(catts))
+ return 0
+
+def generate_Points(catts,tmpMap):
+
+ lines=grass.read_command('v.out.ascii', input=tmpMap, output='-', type='line', cats=catts, format='standar', layer=2, quiet=True)
+ lines = lines.split('L ')
+
+ lines = [d.splitlines(0) for d in lines[1:]]
+ for i,line in enumerate(lines):
+ lines[i]=[p.split() for p in line[1:-1]]
+
+ ptos=grass.read_command('v.out.ascii', input=tmpMap, type='point')
+ ptos = ptos.splitlines(0)
+ i = int(ptos[-1].split("|")[-1])
+
+ sal_puntos=""
+ for j,line in enumerate(lines):
+ for pp in line:
+ if pp==[]: continue
+ sal_puntos+="P 1 1 \n"
+ sal_puntos+=" "+str(pp[0])+" "+str(pp[1])+" "+str(pp[2])+" \n"
+ sal_puntos+=" 1 "+str(i+1)+" \n"
+ i=i+1
+ print sal_puntos
+ grass.write_command('v.edit', tool='add', map=tmpMap, flags='n', type='point', input='-', stdin=sal_puntos, overwrite=True)
+ grass.run_command('v.to.db', map=tmpMap, option='cat', columns='cat')
+
+ return 0
+
+# ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### #
+
+def patch_maps(pointsMap,linesMap,tmpMap):
+
+ grass.run_command('g.copy', vect=pointsMap+','+tmpMap)
+
+ # lee las columnas de mapa de puntos y inserta las columnas si no existen
+ columsPts=grass.read_command('db.columns', table=tmpMap)
+ columsPts=columsPts.splitlines(0)
+ if 'x' not in columsPts:
+ grass.run_command('v.db.addcolumn', map=tmpMap, columns='x double')
+ if 'y' not in columsPts:
+ grass.run_command('v.db.addcolumn', map=tmpMap, columns='y double')
+ if 'z' not in columsPts:
+ grass.run_command('v.db.addcolumn', map=tmpMap, columns='z double')
+ if 'action' not in columsPts:
+ grass.run_command('v.db.addcolumn', map=tmpMap, columns='action varchar(4)')
+
+ # Inserta una nueva tabla en layer 2
+ grass.run_command('v.db.addtable', map=tmpMap, layer=2, key='cat2', table=tmpMap+"_breaklines")
+
+ # lee mapa de lineas e incrementa el layer a 2
+ lines=grass.read_command('v.out.ascii', input=linesMap, output='-', format='standar', layer=1, quiet=True)
+ lines=lines.replace(' 1 ',' 2 ')
+
+ grass.write_command('v.edit', tool='add', map=tmpMap, input='-', stdin=lines, overwrite=True)
+ grass.run_command('v.to.db', map=tmpMap, layer=2, option='cat', columns='cat2')
+
+ # Insertamos las columnas de breaklines si existen
+ tables=grass.read_command('db.tables', flags='p')
+ tables=tables.splitlines(0)
+
+ if linesMap.split('@')[0] in tables:
+ columsLines=grass.read_command('db.columns', table=linesMap.split('@')[0])
+ columsLines=columsLines.splitlines(0)[1:]
+ columsLines=','.join(columsLines)
+
+ grass.run_command('v.db.join', map=tmpMap, layer=2, column='cat2', otable=linesMap.split('@')[0], ocolumn='cat', scolumns=columsLines)
+
+ return 0
+
+###############################################################################
+
+def update_tmpMap(tmpMap):
+
+ # leemos las coordenadas de la tabla y la accion
+ lines = grass.read_command('v.db.select', map=tmpMap, columns='cat,x,y,z,action', quiet=True)
+ lines = [d.split('|') for d in lines.splitlines(0)]
+ del lines[0]
+
+ # leemos las coordenadas originales de los puntos
+ lines_org = grass.read_command('v.to.db', flags='p', map=tmpMap, layer=1, option='coor', columns='x,y,z', quiet=True)
+ lines_org = [d.split('|') for d in lines_org.splitlines(0)]
+
+ # Calculamos las distancias entre las dos coordenadas
+ distances=[]
+ for i,pto in enumerate(lines):
+ if pto[4] != '' and pto[4] != ' ':
+ for j,pto_org in enumerate(lines_org):
+ if pto[0] == pto_org[0]:
+
+ distan=[pto[0]]
+ if pto[1] != '':
+ distan.append(float(pto[1])-float(pto_org[1]))
+ else:
+ distan.append(0)
+ if pto[2] != '':
+ distan.append(float(pto[2])-float(pto_org[2]))
+ else:
+ distan.append(0)
+ if pto[3] != '':
+ distan.append(float(pto[3])-float(pto_org[3]))
+ else:
+ distan.append(0)
+ break
+ distances.append(distan)
+
+ if distances != []:
+ # Movemos los puntos seleccionados en action
+ for dist in distances:
+
+ g.message(dist[0]+','+str(dist[1])+','+str(dist[2])+','+str(dist[3]))
+ grass.run_command('v.edit', map=tmpMap, type='point', tool='move',
+ move=str(dist[1])+','+str(dist[2])+','+str(dist[3]), cats=dist[0], quiet=True)
+ # Actualizamos action en la tabla
+ grass.run_command('v.db.update', map=tmpMap, layer=1, column='action', value=' ', where='cat='+dist[0])
+
+ return 0
+
+###############################################################################
+
+def split_maps(pointsMap,linesMap,tmpMap):
+
+ grass.run_command('v.extract', input=tmpMap, type='point', where='cat>0', output=pointsMap, quiet=True)
+
+ lines=grass.read_command('v.out.ascii', input=tmpMap, output='-', format='standar', layer=2, quiet=True)
+ lines=lines.replace(' 2 ',' 1 ')
+ grass.write_command('v.edit', tool='create', map=linesMap, input='-', stdin=lines, quiet=True)
+ grass.run_command('v.db.addtable', map=linesMap, layer=1, key='cat', quiet=True)
+ grass.run_command('v.to.db', map=linesMap, layer=1, option='cat', columns='cat', quiet=True)
+
+ columsLines=grass.read_command('db.columns', table=tmpMap+"_breaklines", quiet=True)
+
+ if len(columsLines) > 1:
+ columsLines=columsLines.splitlines(0)[1:]
+ columsLines=','.join(columsLines)
+ grass.run_command('v.db.join', map=linesMap, layer=1, column='cat', otable=tmpMap+"_breaklines", ocolumn='cat2', scolumns=columsLines, quiet=True)
+
+ return 0
+
+# ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### #
+# ### Main
+# ###
+
+def main():
+
+ PtosMap=options['points']
+ BreakMap=options['breaklines']
+ HullMap=options['contorno']
+
+ if re.search(r'_Points', options['points']):
+ OutName=PtosMap.replace("_Points","")
+ else:
+ OutName=options['points']
+
+ OutName=OutName.split('@')[0]
+
+ if re.search(r'^_', options['tin']): nameTin=OutName+options['tin']
+ else: nameTin=options['tin']
+
+ if re.search(r'^_', options['tindem']): nameTinDem=OutName+options['tindem']
+ else: nameTinDem=options['tindem']
+
+ if re.search(r'^_', options['tincontour']): nameTinDemC=OutName+options['tincontour']
+ else: nameTinDemC=options['tincontour']
+
+ if re.search(r'^_', options['nndem']): nameNNDem=OutName+options['nndem']
+ else: nameNNDem=options['nndem']
+
+ if re.search(r'^_', options['nncontour']): nameNNDemC=OutName+options['nncontour']
+ else: nameNNDemC=options['nncontour']
+
+
+ if re.search(r'^_', options['tindemDem']): nameTinDemDem=OutName+options['tindemDem']
+ else: nameTinDemDem=options['tindemDem']
+
+ if re.search(r'^_', options['tinDemcontour']): nameTinDemDemC=OutName+options['tinDemcontour']
+ else: nameTinDemDemC=options['tinDemcontour']
+
+ if re.search(r'^_', options['nndemDem']): nameNNDemDem=OutName+options['nndemDem']
+ else: nameNNDemDem=options['nndemDem']
+
+ if re.search(r'^_', options['nnDemcontour']): nameNNDemDemC=OutName+options['nnDemcontour']
+ else: nameNNDemDemC=options['nnDemcontour']
+
+
+
+#### Point section ####
+
+ if flags['e']:
+ grass.message("Patching points and breaklines maps")
+ patch_maps(PtosMap,BreakMap,options['topomap'])
+
+ if flags['u']:
+ grass.message("Updating topo map")
+ update_tmpMap(options['topomap'])
+
+ if flags['s']:
+ grass.message("Spliting in points and breaklines maps")
+ split_maps(PtosMap,BreakMap,options['topomap'])
+
+#### Breaklines section ####
+
+ if flags['p']:
+ grass.message("Generate points in vertices of polylines")
+ generate_Points(int(options['cats']),options['topomap'])
+
+ if flags['l']:
+ grass.message("Spliting breaklines map")
+ generate_Segments(options['splitlength'],int(options['cats']),options['topomap'])
+
+#### Sup section ####
+
+ if flags['t']:
+ grass.message("Triangulating tin")
+ cmd='/mnt/trb/addons/v.triangle/v.triangle.v7 points='+PtosMap
+ if BreakMap != "": cmd += ' lines='+BreakMap
+ cmd += ' tin='+nameTin+' --o --q'
+ os.system(cmd)
+ print cmd
+
+ grass.run_command('v.tin.to.raster', map=nameTin, output=OutName+"_TinDem_borrar", overwrite=True)
+ grass.run_command('v.to.rast', input=HullMap, output=HullMap, use='val', overwrite=True, quiet=True)
+
+ os.system("r.mapcalc '"+nameTinDem+" = if("+HullMap+"==1"+", "+OutName+"_TinDem_borrar"+", null())' --o --q")
+ grass.run_command('r.contour', input=nameTinDem, output=nameTinDemC, step=options['step'], overwrite=True, quiet=True)
+ grass.run_command('g.remove', rast=OutName+"_TinDem_borrar")
+
+ if flags['y']:
+ grass.message("Triangulating nnbathy")
+ grass.run_command('v.to.rast', input=PtosMap, output=PtosMap, use='z', overwrite=True, quiet=True)
+
+ os.system('/mnt/trb/addons/r.surf.nnbathy/r.surf.nnbathy.sh input='+PtosMap+' output='+OutName+"_NNbathyDem_borrar"+' alg=l --o --q')
+
+ grass.run_command('v.to.rast', input=HullMap, output=HullMap, use='val', overwrite=True, quiet=True)
+
+ os.system("r.mapcalc '"+nameNNDem+" = if("+HullMap+"==1"+", "+OutName+"_NNbathyDem_borrar"+", null())' --o --q")
+ grass.run_command('r.contour', input=nameNNDem, output=nameNNDemC, step=options['step'], overwrite=True, quiet=True)
+ grass.run_command('g.remove', rast=OutName+"_NNbathyDem_borrar")
+
+ if flags['f']:
+ grass.message("Patching timdem with dem")
+ grass.run_command('r.patch', input=nameTinDem+","+options['dem'], output=nameTinDemDem, overwrite=True, quiet=True)
+ grass.run_command('r.contour', input=nameTinDemDem, output=nameTinDemDemC, step=1, overwrite=True, quiet=True)
+
+ if flags['n']:
+ grass.message("Patching nnbathydem with dem")
+ grass.run_command('r.patch', input=nameNNDem+","+options['dem'], output=nameNNDemDem, overwrite=True, quiet=True)
+ grass.run_command('r.contour', input=nameNNDemDem, output=nameNNDemDemC, step=1, overwrite=True, quiet=True)
+
+ sys.exit(0)
+
+if __name__ == "__main__":
+ options, flags = grass.parser()
+ main()
Property changes on: grass-addons/grass7/vector/v.civil/v.civil.topo/v.civil.topo.py
___________________________________________________________________
Added: svn:mime-type
+ text/x-python
Added: svn:eol-style
+ native
More information about the grass-commit
mailing list