Take the 2-minute tour ×
Geographic Information Systems Stack Exchange is a question and answer site for cartographers, geographers and GIS professionals. It's 100% free, no registration required.

I asked a question earlier about how to calculate a midpoint. I was wondering if there is a Python library that does this alreay?

share|improve this question
    
arcpy from ArcGIS software installation, but I'm assuming your are looking for something open-source. There's probably a module from QGIS if that helps. –  GIS Danny Mar 22 '13 at 17:33
add comment

3 Answers 3

This python package

http://pypi.python.org/pypi/geographiclib

solves the direct and inverse geodesic problems for an ellipsoid. In your case, solve the inverse problem to get a distance and an azimuth at the first point and then solve the direct problem from the first point with the computed azimuth and half the distance. For a quick illustration of using this package, see

http://geographiclib.sourceforge.net/html/other.html#python

Here, for example, is the calculation of the point midway between JFK and SIN airports

# Compute point midway between JFK and SIN airports
import sys
sys.path.append("/usr/local/lib/python/site-packages")
from geographiclib.geodesic import Geodesic

# Coordinates of airports
lat1,lon1 = 40.640,-73.779 # JFK
lat2,lon2 =  1.359,103.989 # SIN

# Compute path from 1 to 2
g = Geodesic.WGS84.Inverse(lat1, lon1, lat2, lon2);

# Compute midpoint starting at 1
h1 = Geodesic.WGS84.Direct(lat1, lon1, g['azi1'], g['s12']/2);
print(h1['lat2'],h1['lon2']);

# Alternatively, compute midpoint starting at 2
h2 = Geodesic.WGS84.Direct(lat2, lon2, g['azi2'], -g['s12']/2);
print(h2['lat2'],h2['lon2']);

On my system, this gives

(70.34117458722292, 97.02347775257729)
(70.34117458722295, 97.02347775257725)
share|improve this answer
add comment

You could try the Angles library, which does bearing and distance between points on the unit sphere, but I think you'd still have to calculate the midpoint yourself.

share|improve this answer
add comment

GDAL/OGR (Geospatial Data Abstraction Library) has a centroid function.

The following script reads in a line shapefile and writes the midpoint out to a point shapefile (source):

#script that finds the midpoint of transect.shp
#script then creates a point shapefile from the line midpoint

from osgeo import ogr
myDriver = ogr.GetDriverByName('ESRI Shapefile')

#create layer from transect shapefile
dataSource = myDriver.Open('transect.shp',1)
lineLayer = dataSource.GetLayer()

#get geometry from first line feature in layer:
lineFeature = lineLayer.GetFeature(0)
lineGeom = lineFeature.geometry()

#get midpoint of line and get x/y coords
midPoint = lineGeom.Centroid()
xMid = midPoint.GetX()
yMid = midPoint.GetY()

#Create a point shapefile for the midpoint
mySource = myDriver.CreateDataSource('MidPoint.shp')
myLayer = mySource.CreateLayer('layer1',geom_type=ogr.wkbPoint)
featObj = ogr.Feature(myLayer.GetLayerDefn())
point_object = ogr.Geometry(ogr.wkbPoint)
point_object.AddPoint(xMid,yMid)
featObj.SetGeometry(point_object)
myLayer.CreateFeature(featObj)

#close both shapefiles
mySource.Destroy()
dataSource.Destroy()
share|improve this answer
add comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.