How can I create a simple LoD1 citygml-model using Python?
[GIS] Creating simple citygml 3D models based on 2D shapefiles (alkis) shp2gml using Python
gmlpythonshapefile
Related Solutions
I believe this was 'cross-posted' on another site... and the problem was answered with code below. (guidance by Mathew Coyle and Dan Patterson)
import arcpy
# input fc
fc = "D:\\2010.shp"
# the yr 'breaks'...hard-coded as a Python list.
yrs = [1991,1992,1996,1998,1999,2000,2003,2004,2005,2006,2008,2010]
# initialize variable (part of SQL qry) to remove previous selections
subtract = ''
for yr in yrs:
where = '"Year" <= ' + str(yr) + subtract
subtract = ' AND "Year" > ' + str(yr)
filename = str(yr) + '.shp'
arcpy.FeatureClassToFeatureClass_conversion(fc, "D:\\test output", filename, where)
You must first understand how PyQGIS handles geometry (Geometry Handling)
The most important element is the point:
QgsPoint(x,y)
and a line or a segment of line are composed of two points:
QgsGeometry.fromPolyline([QgsPoint(x1,y1),QgsPoint(x2,y2)]));
So to construct a line:
line_start = QgsPoint(50,50)
line_end = QgsPoint(100,150)
line = QgsGeometry.fromPolyline([line_start,line_end])
and with a memory layer (geometry only, without the attributes):
# create a new memory layer
v_layer = QgsVectorLayer("LineString", "line", "memory")
pr = v_layer.dataProvider()
# create a new feature
seg = QgsFeature()
# add the geometry to the feature,
seg.setGeometry(QgsGeometry.fromPolyline([line_start, line_end]))
# ...it was here that you can add attributes, after having defined....
# add the geometry to the layer
pr.addFeatures( [ seg ] )
# update extent of the layer (not necessary)
v_layer.updateExtents()
# show the line
QgsMapLayerRegistry.instance().addMapLayers([v_layer])
the result is:
with 3 points, just add it as a new feature:
newpoint = QgsPoint(143,125)
v_layer = QgsVectorLayer("LineString", "line_3pt", "memory")
pr = v_layer.dataProvider()
seg = QgsFeature()
seg.setGeometry(QgsGeometry.fromPolyline([line_start, line_end]))
# new feature: line from line_end to newpoint
seg = QgsFeature()
seg.setGeometry(QgsGeometry.fromPolyline([line_end, newpoint]))
pr.addFeatures( [ seg ] )
v_layer.updateExtents()
# add the line to
QgsMapLayerRegistry.instance().addMapLayers([v_layer])
and the result is:
And with a for loop you can create a line with many segments:
Best Answer
I build my citygmls based on german alkis-files (shape-files), but it can be used also with other building information sources after changing the code.
Main steps
this file is published to help others to find a accelerated start in creating citygml by using alkis-data (shapefile) or creating citygml at all. this file is not well commented or tested.
Main script:
sub-script: iteration for adding buildings
sub-script: reading building inits
sub-script: read shapefile
sub-script: calculating the building 3D polygons
sub-script: find lower and upper corner of all buildings
call the main script
example result gml for 2 buildings: