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'm trying to reduce a set of lat/lons (each with it's own weight) into one lat/lon that is the center of mass of all of the points.

however when I map the results the average point falls no where near where I think it should be (I've tested on multiple sets of points).

here's my code, in python, based on these two answers:
How to find a point half way between two other points?
Computing an averaged latitude and longitude coordinates

import math
import nummpy
import math

def toCartesian(t):
        latD,longD = t
        latR = math.radians(latD)
        longR = math.radians(longD)
        return (
            math.cos(latR)*math.cos(longR), 
            math.cos(latR)*math.sin(longR), 
            math.sin(latR)
            )

def toSpherical(t):
    x,y,z = t

    r = math.hypot(x,y)
    if r == 0:
        if z > 0:
            return (90,0)
        elif z< 0:
            return (-90,0)
        else:
            return None
    else:
        return (math.degrees(math.atan2(z, r)), math.degrees(math.atan2(y,x)))


xyz = numpy.asarray([0.0,0.0,0.0])
total = 0
for p in points:
        weight = p["weight"]
        total += weight
        xyz += numpy.asarray(toCartesian((p["lat"],p["long"])))*weight

avgXYZ = xyz/total
avgLat, avgLong = toSpherical(avgXYZ)
print avgLat,avgLong
share|improve this question

closed as too localized by whuber Jun 14 '13 at 17:28

This question is unlikely to help any future visitors; it is only relevant to a small geographic area, a specific moment in time, or an extraordinarily narrow situation that is not generally applicable to the worldwide audience of the internet. For help making this question more broadly applicable, visit the help center.If this question can be reworded to fit the rules in the help center, please edit the question.

2  
Do you mean weighted mean center? –  gm70560 Jun 14 '13 at 14:35
2  
Just glancing over the code for obvious problems, it strikes me that atan2(r,z) ought to be atan2(z,r). With numerical code like this, it's important to write little testing modules for your functions (like toSpherical and toCartesian) and exercise them well so that you know you can trust them. Why don't you do that and see whether your problems go away? –  whuber Jun 14 '13 at 14:35
    
@gm70560 The ArcGIS approach you refer to requires (a) data to be projected and (b) that they not extend very widely over the earth's surface. The approach shown in the code here is intended for global applications where the ArcGIS "Mean Center" will fail or give incorrect results. –  whuber Jun 14 '13 at 14:37
    
@whuber I am confused as to the use of "center of gravity" in the question, that is all. –  gm70560 Jun 14 '13 at 14:51
1  
I ran a whole bunch of tests on this and found that my data parsing up until this point was the problem, not the calculation of the weighted average...the above is correct –  Yanofsky Jun 14 '13 at 17:09
show 2 more comments