Geopandas auto-fractionné Linestrame Geodataframe d'une manière rapide sans perdre tous les attributs

Nov 24 2020

J'ai geodataframe de nombreux LineStrings. Les LineStrings se croisent mais ne sont pas divisés à ces intersections. Ma solution actuelle pour y parvenir est d' ici :

network = gp.read_file(filenameNetwork)
newNetwork = gp.GeoDataFrame()
for splittedGeom in network.geometry.unary_union:
    part = gp.GeoDataFrame([[splittedGeom]], columns=['geometry'])
    newNetwork = newNetwork.append(part)

Cependant, je perds toutes les colonnes avec cette solution. J'ai aussi essayé ceci, mais cela prend beaucoup trop de temps:

from shapely import ops

streets = streets.reset_index(drop=True)
streets = streets[['geometry', 'costs']]
headers = list(streets.columns)

index = 0
newStreets = gp.GeoDataFrame( columns=['geometry'])
for line in range(len(streets)-1):
    print(line, len(streets))
    linegeom = streets.at[line, 'geometry']
    isNotSplitted = True
    for line2 in range(len(streets)):
        if line2 == line:
            continue
        linegeom2 = streets.at[line2, 'geometry']
        if linegeom2.crosses(linegeom):
            try:
                linegeomsplitted = ops.split(linegeom, linegeom2)
            except:
                continue
            isNotSplitted = False
            for split in range((len(list(linegeomsplitted.geoms)))):
                splittedline = (list(linegeomsplitted.geoms))[split]
                for head in headers:
                    if head == 'geometry':
                        headValue = splittedline
                    else:
                        headValue = streets.at[line, head]
                    newStreets.at[index, head] = headValue
                index += 1
    if isNotSplitted:
        for head in headers:
            headValue = streets.at[line, head]
            newStreets.at[index, head] = headValue
        index += 1
streets = newStreets
streets = streets.drop_duplicates(subset=['geometry'],
                                  keep='first')

Aucune suggestion?

Réponses

4 gene Nov 30 2020 at 11:31

J'ai trouvé une solution.

En utilisant mon exemple :

a) Le fichier de formes d'origine

import geopandas as gpd
df = gpd.read_file("stac-graphe.shp")
df
id   test                geometry
1   test1   LINESTRING (10.244 -273.317, 784.201 -222.924)
2   test2   LINESTRING (210.484 -553.461, 324.991 -4.534)
3   test3   LINESTRING (169.970 -134.276, 126.511 -218.533...
4   test4   LINESTRING (100.000 -433.317, 724.390 -112.341...
5   test5   LINESTRING (232.683 -113.317, 694.146 -445.024...
6   test6   LINESTRING (563.415 -552.341, 559.512 -22.585)

b) Tamponnez la géométrie d'origine pour éviter les problèmes d'arithmétique flottante (en intersectsou within)

df2 = df.copy()
df2.geometry = df2.geometry.buffer(0.01)

c) Utilisez unary_unionpour diviser toutes les LineStrings auto-intersectées

un = df.geometry.unary_union
geom = [i for i in un]
id = [j for j in range(len(geom))]
unary = gpd.GeoDataFrame({"id":id,"geometry":geom})
unary.head()
id                   geometry
0   LINESTRING (10.244 -273.317, 192.920 -261.423)
1   LINESTRING (192.920 -261.423, 272.484 -256.242)
2   LINESTRING (272.484 -256.242, 418.308 -246.748)
3   LINESTRING (418.308 -246.748, 469.403 -243.421)
4   LINESTRING (469.403 -243.421, 561.095 -237.451)

d) Utilisez une jointure spatiale (avec withinou intersect) pour joindre les deux dataframes et récupérer les attributs d'origine

from geopandas.tools import sjoin
result =sjoin(unary, df2, how="inner",op='within')
result.head()
id_left                   geometry               index_right id_right   test
0   LINESTRING (10.244 -273.317, 192.920 -261.423)   0         1       test1
1   LINESTRING (192.920 -261.423, 272.484 -256.242)  0         1       test1
2   LINESTRING (272.484 -256.242, 418.308 -246.748)  0         1       test1
3   LINESTRING (418.308 -246.748, 469.403 -243.421)  0         1       test1
4   LINESTRING (469.403 -243.421, 561.095 -237.451)  0         1       test1

1 GevaertJoep Nov 28 2020 at 20:14

Ce n'est pas une solution, mais cela peut aider: essayer de créer une intersection de l'union donne deux itérables pour chaque attribut où ils se croisent. Cependant, certaines parties de ligne deviennent des points ...:

l'original:

production:

Peut-être qu'avec quelques modifications sur le code ci-dessous, cela pourrait fonctionner?

import geopandas as gp
network = gp.read_file(filepath)
newNetwork = gp.GeoDataFrame()

geom = network.unary_union
newNetwork = gp.GeoDataFrame(columns=network.columns)

for i in range(len(network)):
    for splittedGeom in network.intersection(network.unary_union)[i]:
        part = gp.GeoDataFrame([list(network.loc[i,network.columns[:-1]])+[splittedGeom]],columns=list(network.columns))
        newNetwork = newNetwork.append(part)

newNetwork.plot()