#!/usr/bin/env python
# encoding: utf-8
#title :concave_hull12.py
#description :script accompagnant l'article "Sur la création des enveloppes concaves (concave hull) et les divers moyens d'y parvenir".
#auteur :M.Laloux
#date :20140930
#==============================================================================
from shapely.geometry import shape
from shapely.ops import polygonize
import fiona
import numpy as np
import subprocess
# chemin de l'exécutable hull
hull_path = "./hull"
# lecture du fichier shapefile
points = [[shape(point['geometry']).x,shape(point['geometry']).y] for point in fiona.open("points_test.shp")]
# ou directement
points = [(point['geometry']['coordinates'][0],point['geometry']['coordinates'][1]) for point in fiona.open("points_test.shp")]
# création du fichier nécessaire au programme hull
infile='test.txt'
with open(infile, 'w') as f:
np.savetxt(f, points, delimiter=' ', fmt="%0.7f %0.7f\n")
# appel du programme et obtention de l'alpha shape en fonction d'un rayon de recherche
def get_alpha_shape(infile,rayon):
commande = "%s -A -aa %s -r -m10000000 -oN -oFpoints < %s" % (hull_path, str(rayon), infile)
print >> sys.stderr, "Running command: %s" % command
retcode = subprocess.call(commande, shell=True)
# hull sauve automatiquement le résultat dans un fichier nommé hout-alf
results_file = open("hout-alf")
results_file.next()
results_indices = [[int(i) for i in line.rstrip().split()] for line in results_file]
results_file.close()
result = [(points[i], points[j]) for i,j in results_indices]
return result
.....
finding alpha shape; output to hout-alf
alpha=1212.26
rayon = 100
alpha_shape = get_alpha_shape(infile,rayon)
# polygone résultant avec la fonction polygonize de Shapely qu'il est ensuite possible de sauver avrc Fiona
resultat= list(polygonize(alpha_shape))
