[Talk-cz] Import of farmland from LPIS

[Pavel Machek]

Hi!

I'd like to start import of LPIS farmland database, as we have very
good coverage of houses, forests and water, but farmland is very good
at places and completely missing at different places.

Import page is at https://wiki.openstreetmap.org/wiki/LPIS ,  import
script is attached to this email.

Best regards,
									Pavel
-- 
(english) http://www.livejournal.com/~pavelmachek
(cesky, pictures) http://atrey.karlin.mff.cuni.cz/~pavel/picture/horses/blog.html
------------- další část ---------------
#!/usr/bin/python
# http://www.lpis.cz/index.html
# Mapa:
# http://eagri.cz/public/app/lpisext/lpis/verejny/

# Vyzadat data na:
# http://eagri.cz/public/app/lpisext/lpis/verejny/exportDat.html
# Vyzadat: Pudni bloky
# Souradny systrem: wgs-84

# Aha, a ty cisla katastralnich uzemi jde dokonce vykoukat z osm, jako ref u boundary=administrative, admin_level=10
# A... nefunguje jim dobre captcha, takze jde stahnout vic uzemi jednim formularem.
# doubek: 631035
# babice: 600601


# for A in *.zip; do unzip $A; done

import sys
import shapefile
from pyproj import *

print "<?xml version='1.0' encoding='UTF-8'?>"
print "<osm version='0.6' generator='pyshp'>"

#p1 = Proj(init='esri:102067')
p1 = Proj(init='EPSG:4326')
p2 = Proj(init='EPSG:4326')

sf = shapefile.Reader(sys.argv[1])

global field
field = {}

def a_field(num, name):
    global field
    i1, i2, i3, i4 = sf.fields[num]
    if name != i1:
        print 'Unexpected field name <', name
    field[name] = num

a_field(1, 'ID_FB')
a_field(3, 'NKODFB')
a_field(6, 'PLATNYOD')
a_field(7, 'VYMERAM')
a_field(8, 'KULTURA')
a_field(9, 'KULTURA_KL')
a_field(10, 'KULTURAOD')
a_field(11, 'EKO')
a_field(21, 'VYSKA')
a_field(22, 'SVAZITOST')
a_field(26, 'KU_KOD')

global nodeid
nodeid = 0

def convert(point):
    lon, lat = point
    lon, lat = transform(p1, p2, lon, lat)
    #lat += 50.013082018919185-50.013834
    #lon += 14.748617781670555-14.749757
    return lon, lat

def write_point(point):
    global nodeid
    lon, lat = convert(point)
    tags = '<tag k="created_by" v="shpupload"/>'
    tags += '<tag k="source" v="lpis"/>'
    nodeid -= 1
    print '<node id="%d" lon="%f" lat="%f\">%s</node>' % ( nodeid, lon, lat, tags )
    return nodeid

def attr(shrec, name):
    return shrec.record[field[name]]

for shrec in sf.shapeRecords():
    shape = shrec.shape

    pts = []
    for point in shape.points:
        pts += [ write_point(point) ]

    nodeid -= 1
    print '<way id="%d">' % nodeid
    print '  <tag k="created_by" v="pyshp"/>'
    print '  <tag k="id_fb" v="%s"/>' % attr(shrec, 'ID_FB')
    print '  <tag k="source" v="lpis"/>'
    print '  <tag k="lpis:kultura" v="%d"/>' % attr(shrec, 'KULTURA')

    kul = attr(shrec, 'KULTURA')
    if kul == 2:    print '  <tag k="landuse" v="farmland"/>'
    elif kul == 3:  print '  <tag k="landuse" v="farmland"/><tag k="crop" v="hop"/>'
    elif kul == 30: print '  <tag k="landuse" v="farmland"/><tag k="crop" v="hop"/>'
    elif kul == 31: print '  <tag k="landuse" v="farmland"/><tag k="crop" v="hop"/>'
    elif kul == 41: print '  <tag k="landuse" v="vineyard"/><tag k="barrier" v="fence"/>'
    elif kul == 61: print '  <tag k="landuse" v="orchard"/><tag k="barrier" v="fence"/>'
    elif kul == 62: print '  <tag k="landuse" v="orchard"/>'
    elif kul == 7:  print '  <tag k="landuse" v="meadow"/><tag k="meadow" v="agricultural"/>'
    elif kul == 71: print '  <tag k="landuse" v="meadow"/><tag k="meadow" v="agricultural"/><tag k="barrier" v="fence"/>'
    elif kul == 72: print '  <tag k="landuse" v="meadow"/><tag k="meadow" v="agricultural"/>'
    elif kul == 91: print '  <tag k="landuse" v="forest"/><tag k="barrier" v="fence"/>'
    elif kul == 92: print '  <tag k="landuse" v="farm"/><tag k="crop" v="vegetables"/>'
    elif kul == 97: print '  <tag k="landuse" v="reservoir"/>'
    elif kul == 98: print '  <tag k="landuse" v="scrub"/><tag k="note" v="rychle rostouci dreviny"/>'
    elif kul == 99: print '  <tag k="landuse" v="forest"/>'
    else:           print '  <tag k="landuse" v="unknown_farmland_%d"/>' % kul

    print '  <tag k="ele" v="%d"/>' % attr(shrec, 'VYSKA')
    print '  <tag k="ref" v="%s"/>' % attr(shrec, 'NKODFB')

    for pt in pts:
        print '  <nd ref="%d"/>' % pt
    print '</way>'

print '</osm>'