Computing intersections of a Polyline and a line
In getIntersectLon, I didn't find the docstring very helpful. For y = m*x + b, I think of m as rise over run, but in (a[0] - b[0]) / (a[1] - b[1]) I think I'm seeing run over rise. Ok. So we're multiplying 1/m by some delta_latitude, then we subtract, which seems like it could be related to the docstring. Maybe the docstring should talk about how the line described by two points intersects a line of latitude. Given that (lat,lng) are spherical coordinates, the approximation seems to have accuracy that falls off with increasing distance. Maybe this function should be named getLongitudeAt, with 1st parameter being the line of latitude.

It's not obvious to me that a[1] - b[1] shall be non-zero.

The latitudes == 0 test is maybe more stringent than you would like. Consider using abs(latitudes) &lt; epsilon.

If you're mostly interested in sign changes, you might want to define sgn() and then work with small signum values.

def sgn(n):
    return math.copysign(1, n)


The comments "THIS IS THE {|END OF} THE CORE ROUTINE" suggest that you may want to break out a helper routine.