1

I have some GeoJSON files that I converted from some KML someone gave me to analyze.

I have used GPSBabel to convert them from KML to GeoJSON. The reason for doing so is that I am able to easily process them using python 3.5 libraries like json and shapely (my environment does not work with packages like fiona for version incompatibilities).

A sample command is: gpsbabel -w -t -i kml -f input.kml -o geojson -F output.json

  • The KML contains various <Polygon> each one with one <LinearRing> inside.
  • The GeoJSON contains various geometry coordinates of type called LineString (this comes from the GPSBabel automatic conversion).
  • Each point in the KML/GeoJSON is a 3 fields array which should be [longitude,latitude,altitude] (according to this: https://developers.google.com/kml/documentation/kmlreference#linearring).

I've written the following function to establish if an input Point (Shapely geometry) is contained in the LineStrings (Shapely geometry) found in the input:

from shapely.geometry import Point, LineString

def poly_contains(point, poly_js):
    for feature in poly_js['features']:
        lineString = LineString(feature['geometry']['coordinates'])
        if lineString.contains(point):
            return True
    return False

Each input point is created from a CSV containing the coordinate pairs with something like this: Point(row['LatitudeId'], row['LongitudeId'])

When applying the previous function to all these points I can not find any point contained in those shapes, but I know for sure some of those points are indeed contained!

I thought the reason could be the first and the last point of a LineString don't match (i.e. the "polyline" is not a "polygon"), but the following piece of code makes sure I am dealing with "polygons" defined with matching first and last element:

def check_file(file_content):
    for feature in file_content['features']:
        lineString = LineString(feature['geometry']['coordinates'])
        if lineString.coords[-1] != lineString.coords[0]:
            print("Issue!")

When applying the previous function to these GeoJSONs like this:

import json
with open('.../.../myfile.json') as f:
    stuff = json.load(f)

check_file(stuff)

then I don't see any output like Issue!, so all these LineStrings should represent "polygons" right?

I don't have many details about those projections and stuff like that, but I am assuming this is all Web Mercator in WGS84 coordinate system as this data is used to be represented on Google Maps (both the original KMLs and the points from the CSV that I need to check).

What's wrong with what I am doing? Why that function contains applied to LinearString and Point does not work as expected?

Edit:

I may be confusing WGS84 and Web Mercator (EPSG:3857), I am not sure. However I have tried to project all those shapes (the points from the CSV and the LineString from the KMLs/GeoJSONs) using pyproj, then wrapping those shapes in Polygon instead of LineString like this:

from shapely.geometry import Point, Polygon
from pyproj import Proj
p_web = Proj(init='EPSG:3857')

# applying 2D Web Mercator projection
def poly_contains_v2(point, poly_js):
    point_proj = Point(p_web(point.coords[0][0], point.coords[0][1]))
    for feature in poly_js['features']:
        new_coords = []
        for ring in feature['geometry']['coordinates']:
            x2, y2 = p_web(ring[0], ring[1])
            new_coords.append([x2, y2])
        polygon = Polygon(new_coords)
        if polygon.contains(point_proj):
            return True
    return False

but still I can not find any point included in those "polygons". I have also tried to apply just a Polygon wrap instead of a LineString (with no projection) with no luck.

  • Linear ring is just a line. I feel that you somehow believe that the area inside the ring belongs to the feature but it does not. You must construct a polygon first. – user30184 Nov 21 '17 at 15:45
  • Using Polygon (with or without a projection first) does not give any included point, cf. my edit at the question above ^ – TPPZ Nov 21 '17 at 16:20
  • Have you pulled these files into a GIS and confirmed they overlap? As another comment stated: shapely is a geometry library so the shapely.geometry.contains() function is strict and literal. LineString is a 1-d topology, so points have to fall on the line. Floating point differences can apply here as well. If you upload an image with the situation that could help troubleshoot. – Logan Byers Nov 29 '17 at 21:27
  • Thanks @LoganByers but I am not in a position of sharing this image at the moment. I am sure some of these points fall inside some of those polygons. Bare in mind the source of the polygons (KML) is different than the source of the points (database dump), however I've seen both sides as a layer on top of Google Maps (on separate websites). I hoped this was a well known issue, but probably there's something deeper to investigate regarding projections. Unfortunately I can not provide more details regarding the data at the moment. – TPPZ Dec 1 '17 at 14:45

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.