Comment 2 for bug 891990

# geo kludge
# I take no responsibility for writing this

import jpype as jp
import shapely

from django.contrib.gis.geos.collections import Polygon

JAR_PATH = '/usr/share/java'
jp.startJVM(jp.getDefaultJVMPath(), '-Djava.ext.dirs=%s' % JAR_PATH)

FT = jp.JClass('org.opensha.sha.faultSurface.FaultTrace')
LOC = jp.JClass('org.opensha.commons.geo.Location')
LOC_LIST = jp.JClass('org.opensha.commons.geo.LocationList')
SGS = jp.JClass('org.opensha.sha.faultSurface.StirlingGriddedSurface')

def fault_poly_from_mls(fault_source_geom, dip, upp_seis_depth, low_seis_depth,
           grid_spacing=1.0):
    """Given a fault source geometry (as a MultiLineString), dip, upper
    seismogenic depth, lower seismogenic depth, and grid spacing (in km),
    create a 3D polygon to represent the fault.

    :param fault_source_geom:
        :class:`django.contrib.gis.geos.collections.MultiLineString`
    :param float dip:
        Angle of dip, from 0.0 to 90.0 degrees (inclusive)
    :param float upp_seis_depth:
        Upper seismogenic depth
    :param float low_seis_depth:
        Lower seismogenic depth
    :param float grid_spacing:
        Grid spacing in kilometers. Default to 1.0.

    :returns:
        3D polygon representing the complete fault geometry
    :rtype:
        :class:`django.contrib.gis.geos.collections.Polygon`
    """

    coords = fault_source_geom.coords

    fault_trace = FT('')
    for line_str in coords:
        for lon, lat in line_str:
            # warning: the ordering of lat/lon is switched here
            # be careful
            loc = LOC(lat, lon)
            fault_trace.add(loc)

    surface = SGS(fault_trace, dip, upp_seis_depth, low_seis_depth,
                  grid_spacing)

    # now we make a polygon with the perimeter coords:
    poly_coords = []
    for per_loc in surface.getSurfacePerimeterLocsList():
        lon = per_loc.getLongitude()
        lat = per_loc.getLatitude()
        depth = per_loc.getDepth()

        poly_coords.append((lon, lat, depth))

    return Polygon(poly_coords)

if __name__ == '__main__':
    ls1 = LineString((0, 0), (1, 1))
    ls2 = LineString((2, 2), (3, 3))
    mls1 = MultiLineString(ls1, ls2)

    p = kludge(mls1, 45.0, 10.0, 40.0)

    print p.coords