Le jeudi 30 août 2012 22:52:00 topodom a écrit :
> Mais ceci dit, je trouves déjà bien compliqué d'extraire une simple couche
> point de la base osm ! Il faut quand même sortir l'artillerie lourde juste
> pour sélectionner les épicerie et les campings dans un tampon de 10 kms !
Oulà, t'es en train de me donner des idées toi ;-)
Allez, voilà de quoi démarrer. En PJ un script python qui a besoin d'une URL
xapi et d'une projection et te génère le SHP avec le tag name de chaque point
dans la table d'attributs.
C'est du vite fait pour rigoler (si les objets ne sont pas des points ça va
planter), mais si ça peut intéresser du monde, je pourrai voir pour ajouter un
export SHP au xapiviewer.
Cela nécessite gdal et OsmApi.py :
http://svn.openstreetmap.org/applications/utils/python_lib/OsmApi/OsmApi.py
Pour la suite, on peut passer sur dev-fr
Bon app'
--
Nicolas Dumoulin
http://wiki.openstreetmap.org/wiki/User:NicolasDumoulin
# -*- coding: utf-8 -*-
import sys
import ogr, osr
import string
import OsmApi
import urllib2
# Settings
destEPSG = 2154
xapiURL = "http://www.overpass-api.de/api/xapi?*[natural=peak][bbox=3.0486739730835017,45.73113511552445,3.102918968200614,45.759466337869895]";
# init shapefile
driverName = "ESRI Shapefile"
drv = ogr.GetDriverByName( driverName )
if drv is None:
print "%s driver not available.\n" % driverName
sys.exit( 1 )
ds = drv.CreateDataSource( "point_out.shp" )
if ds is None:
print "Creation of output file failed.\n"
sys.exit( 1 )
# init projections
osmSR = osr.SpatialReference()
osmSR.ImportFromEPSG(4326)
destSR = osr.SpatialReference()
destSR.ImportFromEPSG(destEPSG)
srTrans = osr.CoordinateTransformation(osmSR,destSR)
# create data layer in shapefile
lyr = ds.CreateLayer( "point_out", destSR, ogr.wkbPoint )
if lyr is None:
print "Layer creation failed.\n"
sys.exit( 1 )
# load OSM data
osmApi = OsmApi.OsmApi();
osmXml = open("sample.osm","r").read();
osmData = urllib2.urlopen(xapiURL).read
osmData = osmApi.ParseOsm(osmXml);
# TODO find all tags present in the data and create an attribute for each tag found
# add an attribute for the name in the layer
field_defn = ogr.FieldDefn( "Name", ogr.OFTString )
field_defn.SetWidth( 32 )
if lyr.CreateField ( field_defn ) != 0:
print "Creating Name field failed.\n"
sys.exit( 1 )
# write the data
for entity in osmData:
point = srTrans.TransformPoint(entity['data']['lon'], entity['data']['lat'])
#print str(entity['data']['lat']) + " - " + str(entity['data']['lon'])
#print str(point)
x = point[0]
y = point[1]
name = entity['data']['tag']['name']
feat = ogr.Feature( lyr.GetLayerDefn())
feat.SetField( "Name", name )
pt = ogr.Geometry(ogr.wkbPoint)
pt.SetPoint_2D(0, x, y)
feat.SetGeometry(pt)
if lyr.CreateFeature(feat) != 0:
print "Failed to create feature in shapefile.\n"
sys.exit( 1 )
feat.Destroy()
ds = None_______________________________________________
Talk-fr mailing list
Talk-fr@openstreetmap.org
http://lists.openstreetmap.org/listinfo/talk-fr
