twd97_gpx_to_twd67_gpx.py

#!/usr/bin/python

#rt[100].get('lat')
#rt[100].get('lon')
#rt[100][0].text
#
#
#for child_of_root in rt:
#    if child_of_root.tag == '{http://www.topografix.com/GPX/1/0}wpt': 
#        for cc in child_of_root:
#            if cc.tag =='{http://www.topografix.com/GPX/1/0}name':
#                print "%s %s %s" % (child_of_root.get('lat'),child_of_root.get('lon'), cc.text)
#
#
#for elem in tree.iterfind('wpt/name'):
#    print elem.text
#
#
#
#for child_of_root in rt:
#    if child_of_root.tag == '{http://www.topografix.com/GPX/1/0}wpt': print "%s %s %s" % (child_of_root.get('lat'),child_of_root.get('lon'), child_of_root[0].text)
#
# For XML:
#http://pycoders-weekly-chinese.readthedocs.org/en/latest/issue6/processing-xml-in-python-with-element-tree.html
#
# For Proj
# http://wiki.osgeo.org/wiki/Taiwan_datums
# http://trac.osgeo.org/proj/
# https://pypi.python.org/pypi/pyproj/1.9.0
# Verify:
# http://www.sunriver.com.tw/taiwanmap/grid_tm2_convert.php#a03

from pyproj import Proj
import xml.etree.cElementTree as ET
tree = ET.ElementTree(file='input.gpx')
rt = tree.getroot()

def lonlatTotwd97tm2(lon, lat):
    p = Proj('+proj=tmerc +ellps=GRS80 +lon_0=121 +x_0=250000 +k=0.9999')
    x,y = p(lon, lat)
    return x,y

def twd97tm2Totwd67tm2(x, y):
    # http://gis.thl.ncku.edu.tw/coordtrans/coordtrans.aspx
    #E67 = E97 - 807.8 - A * E97 - B * N97
    #N67 = N97 + 248.6 - A * N97 - B * E97
    A= 0.00001549
    B= 0.000006521
    E67 = x - 807.8 - A * x - B * y
    N67 = y + 248.6 - A * y - B * x
    return E67, N67

def decimanToms(lonlat):
    decimal = lonlat % 1 # yyy.xxx -> 0.xxx
    v = lonlat - decimal
    mid = (1.0/60.0) # m in decimal
    sid = (1.0/3600.0) # s in decimal
    m = (decimal - (decimal % mid))/mid
    s = (decimal - (decimal % sid))/sid - m*60
    print "%i %i %i" % (v, m, s)


for child_of_root in rt:
    if child_of_root.tag == '{http://www.topografix.com/GPX/1/0}wpt': 
        for cc in child_of_root:
            if cc.tag =='{http://www.topografix.com/GPX/1/0}name':
                x, y = lonlatTotwd97tm2(child_of_root.get('lon'), child_of_root.get('lat'))
                x, y = twd97tm2Totwd67tm2(x, y)
                print "%s %s %s" % (x, y, cc.text)