Créer des lignes à partir de coordonnées de paires de points avec ArcPy?

11

J'ai des coordonnées de paires de points (points de début et de fin) que je dois transformer en lignes. Jusqu'à présent, j'utilisais un appendice des deux coordonnées dans a pippo.Point(), a pippo.CalculateGeometry()pour définir la géométrie de chaque piont, et pippo.append(defined geometry)pour identifier la paire de points, puis PointsToLine pour obtenir ma ligne. C'est assez cher à faire pour des centaines de lignes.

Existe-t-il un moyen plus court de procéder?

Par exemple, placez le point de départ et le point d'arrivée de chaque ligne dans différents champs d'une même table et importez directement les lignes sans passer pour la géométrie des points.

Annalisa Minelli
la source

Réponses:

8

Cela lit un tableau (feuille Excel dans ce cas, mais pourrait être n'importe quel type de tableau) qui ressemble à ceci:

entrez la description de l'image ici

S_X est le point X de départ, le point X de fin E_X, de même pour les Y. Nous parcourons le tableau d'entrée, puis pour chaque ligne, définissons les X / Y de début / fin en un point, ajoutons ce point à un tableau, puis créons une polyligne à partir du tableau de deux points. Insérez ensuite dans la classe de fonctions. Rincez et répétez.

import arcpy

in_rows = arcpy.SearchCursor(r"D:\Temp\Lines.xls\Sheet1$")

point = arcpy.Point()
array = arcpy.Array()

featureList = []
cursor = arcpy.InsertCursor(r"D:\Temp\Lines.shp")
feat = cursor.newRow()

for in_row in in_rows:
    # Set X and Y for start and end points
    point.X = in_row.S_X
    point.Y = in_row.S_Y
    array.add(point)
    point.X = in_row.E_X
    point.Y = in_row.E_Y
    array.add(point)   
    # Create a Polyline object based on the array of points
    polyline = arcpy.Polyline(array)
    # Clear the array for future use
    array.removeAll()
    # Append to the list of Polyline objects
    featureList.append(polyline)
    # Insert the feature
    feat.shape = polyline
    cursor.insertRow(feat)
del feat
del cursor

Et vous obtenez vos lignes:

entrez la description de l'image ici

Chad Cooper
la source
Merci, je vais essayer d'estimer la durée de mon analyse .. c'est exactement ce que j'essayais de faire :-)
Annalisa Minelli
Pour la ligne point.X = in_row.S_X, elle renvoie une erreur indiquant que la valeur d'entrée n'est pas numérique. J'ai essayé de le faire int ou float ou même numérique mais ne fonctionne pas car le champ n'est pas un nombre est Nonetype. De l'aide?
Federico Gómez
5

J'ai créé un script python la semaine dernière (sans utiliser ArcPy cependant), qui prend des points qui créent la géométrie des lignes de bus (un point shp) en fonction d'un champ numérique séquentiel ("SEQ"). Vous pouvez facilement le modifier pour prendre les coordonnées d'un champ de la même entité (en utilisant la valeur du champ au lieu de la géométrie).

# -*- coding: utf-8 -*-
###############################################################################
from sys import argv
import osgeo.ogr
import os, os.path
###############################################################################

script, srcSHP = argv

#-- Open source shapefile
shapefile = osgeo.ogr.Open(srcSHP)
layer = shapefile.GetLayer(0)
spatialRef = layer.GetSpatialRef()

#-- Output directory
outDir = os.path.dirname(srcSHP)
outDirName = os.path.basename(outDir)

driver = osgeo.ogr.GetDriverByName("ESRI Shapefile")
outFile = driver.CreateDataSource(os.path.join(outDir,outDirName + "_lines.shp"))
outLayer = outFile.CreateLayer("layer", spatialRef)

#-- Adding fields to the output shapefile
fieldDef = osgeo.ogr.FieldDefn("line_no", osgeo.ogr.OFTString)
fieldDef.SetWidth(12)
outLayer.CreateField(fieldDef)

fieldDef = osgeo.ogr.FieldDefn("From_SEQ", osgeo.ogr.OFTReal)
outLayer.CreateField(fieldDef)

fieldDef = osgeo.ogr.FieldDefn("To_SEQ", osgeo.ogr.OFTReal)
outLayer.CreateField(fieldDef)

#-- Going through each feature, one by one
#-- The last point is the end of the line so I don't want to iterate through that one
for i in range(layer.GetFeatureCount()-1):
    lString = osgeo.ogr.Geometry(osgeo.ogr.wkbLineString)  

    feature1 = layer.GetFeature(i)
    feature2 = layer.GetFeature(i+1)

    # When it's a new line, the sequential number restart to 1, so we don't want that line
    if feature1.GetField("SEQ") < feature2.GetField("SEQ"):
        geom1 = feature1.GetGeometryRef()
        geom2 = feature2.GetGeometryRef()

        geom1x = geom1.GetX()
        geom1y = geom1.GetY()
        geom2x = geom2.GetX()
        geom2y = geom2.GetY()

        lString.AddPoint(geom1x, geom1y)
        lString.AddPoint(geom2x, geom2y)     # Adding the destination point

        #-- Adding the information from the source file to the output
        feat = osgeo.ogr.Feature(outLayer.GetLayerDefn())
        feat.SetGeometry(lString)
        feat.SetField("line_no", feature1.GetField("line_no"))
        feat.SetField("From_SEQ", feature1.GetField("SEQ"))
        feat.SetField("To_SEQ", feature2.GetField("SEQ"))
        outLayer.CreateFeature(feat)

print "The End"

Chaque paire de points créera une seule ligne. Il peut y avoir une façon plus élégante de le faire, mais il a créé 3900 lignes en environ 15 secondes, donc cela fonctionne pour moi ...

fgcartographix
la source
Merci, cela ressemble exactement à une élaboration massive .. qui devrait être vraiment utile pour moi. J'essaierai, puis les commentaires. merci pour l'instant.
Annalisa Minelli
3

vous pouvez utiliser ces deux outils pour créer une couche d'événements XY et des points à aligner , en voyant les paramètres nécessaires en points à aligner (champ de ligne, trier les points) et mettre à jour les données de la table d'entrée, la tâche pourrait être plus simple

geogeek
la source
1

ce n'est qu'une mise à jour de la réponse de @ ChadCooper, car les curseurs "da" remplacent désormais avantageusement les curseurs précédents:

with arcpy.da.SearchCursor(input_table,[orig_namefield,x1,y1,x2,y2] ) as in_rows:
    with arcpy.da.InsertCursor(output_lines,["SHAPE@",name_field]) as cursor:
        for row in in_rows:
            # build array for line segment
            array = arcpy.Array([arcpy.Point(row[1],row[2]),arcpy.Point(row[3],row[4])])
            # Create a Polyline object based on the array of points
            polyline = arcpy.Polyline(array)
            # Insert the feature
            cursor.insertRow([polyline,row[0]])
radouxju
la source