Houdini GeoJSON importer script (python)

hou_geojson_import.py

"""
This script is will 'import' a custom GeoJSON file into Houdini. It is intended to be used in the Houdini SOP context, either in the default python SOP or a custom python operator. These operators need to have the four custom paramters as listed in the 'node parameters' section. 

There are several steps/limitations with the import that must be taken into account:
- Will handle point, line and polygone features
- Will not process polygons with holes
- Will not handle geographic coordinates properly (latitude and longitude). The features must use a projected (X,Y,Z) datum. This technically means that the input GeoJSON will not adhere to the spec which demands coodrinates on WGS84).
- Will not handle NULL or None type attributes for any feature. All features must have a value for all attributes. If there is a NULL for a string use and empty: "". If there is a null for a number use -9999.
"""

# ---------------------------------------------------------------------------------------
# Geo HOU module links
# ---------------------------------------------------------------------------------------

node = hou.pwd()
geo = node.geometry()

# ---------------------------------------------------------------------------------------
# Node parameters
# ---------------------------------------------------------------------------------------

# The GeoJSON File
geojson_file = node.parm('geojson_file').eval()
# The projected geographic centroid of the GeoJSON file - this is written automatically when the file is loaded
geo_cent = [node.parm('geo_centx').eval(), node.parm('geo_centy').eval()]
# The projected geographic centroid of the scene
ref_geo_cent = [node.parm('ref_geo_centx').eval(), node.parm('ref_geo_centy').eval()]
# Toggle to offset the imported features using the reference centroid
use_ref = node.parm('use_ref').eval()

# ---------------------------------------------------------------------------------------
# Import Modules
# ---------------------------------------------------------------------------------------

import json

# ---------------------------------------------------------------------------------------
# Geom and Attribute Creation Defenitions
# ---------------------------------------------------------------------------------------

def vec2tovec3(xyz):

    if len(xyz) < 3:
        return [xyz[0], xyz[1], 0]
    else:
        return xyz

# ---------------------------------------------------------------------------------------

def geoCentroid(features):

    positions = []

    for feature in features:

        geoType = feature['geometry']['type'] # string
        coords = feature['geometry']['coordinates'] # list

        if geoType == 'Point':

            pos = vec2tovec3(coords)
            positions.append(pos)

        elif geoType == 'MultiPoint':

            for pos in coords:
                
                positions.append(vec2tovec3(pos))

        elif geoType == 'LineString':

            for pos in coords:
                
                positions.append(vec2tovec3(pos))

        elif geoType == 'MultiLineString':

            for part in coords:

                for pos in part:

                    positions.append(vec2tovec3(pos))
                
        elif geoType == 'Polygon':

            for pos in coords[0]:

                positions.append(vec2tovec3(pos))

        elif geoType == 'MultiPolygon':

            for part in coords:

                for pos in part[0]:

                    positions.append(vec2tovec3(pos))

    minEast = min(p[0] for p in positions)
    maxEast = max(p[0] for p in positions)

    minNorth = min(p[1] for p in positions)
    maxNorth = max(p[1] for p in positions)

    minAlt = min(p[2] for p in positions)
    maxAlt = max(p[2] for p in positions)

    centroidEast = (minEast + maxEast) / 2
    centroidNorth = (minNorth + maxNorth) / 2
    centroidAltitude = (minAlt + maxAlt) / 2

    return [centroidEast, centroidNorth, centroidAltitude]

# ---------------------------------------------------------------------------------------

def reproject(coord):

    x = coord[0] - centEast         # EASTING
    y = coord[2]                    # ALTITUDE
    z = -(coord[1] - centNorth)     # NORTHING

    return(x, y, z)

# ---------------------------------------------------------------------------------------

def createPt(coord):

    pos = vec2tovec3(coord)

    p = geo.createPoint()
    p.setPosition(reproject(pos))
    return p

# ---------------------------------------------------------------------------------------

def createPts(coords):

    for coord in coords:
        createPt(coord)

# ---------------------------------------------------------------------------------------

def createPrim(coords, openStatus):

    poly = geo.createPolygon()
    poly.setIsClosed(openStatus)

    for coord in coords:
        pt = createPt(coord)
        poly.addVertex(pt)

    return poly

# ---------------------------------------------------------------------------------------

def createAttribs(geoType, name, value):

    if geoType == 'Point' or geoType == 'MultiPoint':
        geo.addAttrib(hou.attribType.Point, name, value)
    else:
        geo.addAttrib(hou.attribType.Prim, name, value)

# ---------------------------------------------------------------------------------------

def setAttribs(geom, attribs):

    try:
        for k, v in attribs.iteritems():
            try:
                if v == '-9999':
                    geom.setAttribValue(k, '')
                elif v == -9999:
                    geom.setAttribValue(k, 0)
                else:
                    geom.setAttribValue(k, v)
            except hou.OperationFailed:
                pass
    except:
        pass

# ---------------------------------------------------------------------------------------

def createFeature(feature):

        geoType = feature['geometry']['type'] # string
        coords = feature['geometry']['coordinates'] # list
        attribs = feature['properties'] # dictionary

        if geoType == 'Point':

            pt = createPt(coords)
            setAttribs(pt, attribs)

        elif geoType == 'MultiPoint':

            for part in coords:

                pt = createPt(part)
                setAttribs(pt, attribs)

        elif geoType == 'LineString':

            line = createPrim(coords, False)
            setAttribs(line, attribs)

        elif geoType == 'MultiLineString':

            for part in coords:

                line = createPrim(coords, False)
                setAttribs(line, attribs)

        elif geoType == 'Polygon':

            prim = createPrim(coords[0], True) # Need to implement prim holes
            setAttribs(prim, attribs)

        elif geoType == 'MultiPolygon':

            for part in coords:

                prim = createPrim(part[0], True) # Need to implement prim holes
                setAttribs(prim, fields)

# ---------------------------------------------------------------------------------------
# Parse GeoJSON File
# ---------------------------------------------------------------------------------------

with open(geojson_file) as f:
    data = json.load(f)

# ---------------------------------------------------------------------------------------
# Set Geo-Centroid of Node
# ---------------------------------------------------------------------------------------

geoCent = geoCentroid(data['features'])
node.setParms({'geo_centx' : geoCent[0], 'geo_centy' : geoCent[1]})

centEast = geoCent[0]
centNorth = geoCent[1]

if use_ref:
    
    centEast = ref_geo_cent[0]
    centNorth = ref_geo_cent[1]

# ---------------------------------------------------------------------------------------
# Create Feilds
# ---------------------------------------------------------------------------------------

fields = [str(k) for k, v in data['features'][0]['properties'].iteritems()]
globalGeoType = data['features'][0]['geometry']['type']

for field in fields:

    fieldValues = []

    for feature in data['features']:

        for k, v in feature['properties'].iteritems():

            if str(k) == field:

                if v is not None:

                    fieldValues.append(v)

    value = min(v for v in fieldValues)

    if isinstance(value, basestring):

        createAttribs(globalGeoType, field, str(''))

    elif isinstance(value, int):

        createAttribs(globalGeoType, field, int(0))

    elif isinstance(value, long):

        createAttribs(globalGeoType, field, long(0))

    elif isinstance(value, float):

        createAttribs(globalGeoType, field, float(0.0))

    del fieldValues

for feature in data['features']:

    createFeature(feature)