Wednesday, January 10, 2018

convert between latitude longitude and georef

World Geographic Reference System (GEOREF)
https://en.wikipedia.org/wiki/World_Geographic_Reference_System



import sys

georef_uppercase = 'ABCDEFGHJKLMNPQRSTUVWXYZ'

def ll_2_georef( lat, lon ):
    f_lat, f_lon = lat+90, lon+180
    lat_0, lat_1 = divmod( int(f_lat), 15 )
    lon_0, lon_1 = divmod( int(f_lon), 15 )
    lat_m, lon_m = 6000*(f_lat-int(f_lat)), 6000*(f_lon-int(f_lon))
    return "{lon_0}{lat_0}{lon_1}{lat_1}{lon_m:04d}{lat_m:04d}". \
        format(
        lon_0= georef_uppercase[lon_0],
        lat_0= georef_uppercase[lat_0],
        lon_1= georef_uppercase[lon_1],
        lat_1= georef_uppercase[lat_1],
        lon_m= int(lon_m),
        lat_m= int(lat_m),
        )

def georef_2_ll( grid ):
    lon_0, lat_0, lon_1, lat_1= grid[:4]

    rest= grid[4:]
    pos= len(rest)//2
    if pos:
        scale= { 4: 100, 3: 10, 2: 1 }[pos]
        lon_frac, lat_frac = float(rest[:pos])/scale,\
            float(rest[pos:])/scale
    else:
        lon_frac, lat_frac = 0, 0
    lat= georef_uppercase.find(lat_0)*15+georef_uppercase.\
        find(lat_1)+lat_frac/60
    lon= georef_uppercase.find(lon_0)*15+georef_uppercase.\
        find(lon_1)+lon_frac/60
    return lat-90, lon-180



#python georeff.py 45.234  89.67
#python georeff.py TKQA40201404

if len(sys.argv) == 2:
    georef1 = sys.argv[1]
    print georef1
    print(georef_2_ll(georef1))
elif len(sys.argv) > 2:
    lat = float(sys.argv[1])
    lon= float(sys.argv[2])
    print(lat, lon)
    print(ll_2_georef(lat, lon))







1- usage
$ python georeff.py 45.234  89.67
$ python georeff.py TKQA40201404

No comments:

Post a Comment