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 use this python code to calculate a greatcircle: https://code.google.com/p/pyroms/source/browse/trunk/pyroms/greatcircle.py?r=39 After initializing the greatcircle I can calculate points on the line.

I would like to do the same with an arbitrary smallcircle. I would like to define the smallcircle using its center location (in latlon) and one location on the line (also latlon), so not using a diameter. Say as an input I would use a certain point spacing in meters, how would I then proceed to calculate the other locations on the smallcircle line?

Thanks!

share|improve this question
add comment

1 Answer

Python:

# From: http://www.movable-type.co.uk/scripts/latlong.html

def distance(lon1, lat1, lon2, lat2):
    # http://code.activestate.com/recipes/576779-calculating-distance-between-two-geographic-points/
    # http://en.wikipedia.org/wiki/Haversine_formula
    dlat = lat2 - lat1
    dlon = lon2 - lon1
    q = sin(dlat/2)**2 + (cos(lat1) * cos(lat2) * (sin(dlon/2)**2))
    return 2 * atan2(sqrt(q), sqrt(1-q))

from math import *

def SmallCircle(radius,clon,clat,llon,llat,meters):
    # clon = circle central longitude (deg)
    # clat = circle central latitude (deg)
    # llon = longitude of point on circle line (deg)
    # llat = latitude of point on circle line (deg)
    # meters = point spacing (meters). This is not the exact spacing, it will be used to calculate nrpoints.
    # radius = radius of sphere

    dlat = radians(llat-clat)
    dlon = radians(llon-clon)
    clat = radians(clat)
    clon = radians(clon)
    llat = radians(llat)
    llon = radians(llon)

    a = sin(dlat/2) * sin(dlat/2) + sin(dlon/2) * sin(dlon/2) * cos(clat) * cos(llat)
    c = 2 * atan2(sqrt(a), sqrt(1-a))
    d = radius * c

    # bearing
    y = sin(dlon) * cos(llat)
    x = cos(clat)*sin(llat) - sin(clat)*cos(llat)*cos(dlon)
    brng = atan2(y, x)

    # circle
    angle = distance(clon,clat,llon,llat)
    circleradius = sin(angle) * radius
    circlelength = 2*pi*circleradius

    # loop circle
    list = []
    list.append([degrees(llon),degrees(llat)])
    nrpoints = int(round(circlelength / meters))
    for i in range(1, nrpoints):
        brng = brng + (float(i) / float(nrpoints))*2*pi
        lat = asin( sin(clat)*cos(d/radius) + cos(clat)*sin(d/radius)*cos(brng) )
        lon = clon + atan2(sin(brng)*sin(d/radius)*cos(clat), cos(d/radius)-sin(clat)*sin(lat))
        list.append([degrees(lon),degrees(lat)])

    return list
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.