#!/usr/bin/env python
#******************************************************************************
#  $Id: epsg_tr.py 12994 2007-11-24 19:27:20Z hobu $
# 
#  Project:  CFS OGC MapServer
#  Purpose:  Script to create WKT and PROJ.4 dictionaries for EPSG GCS/PCS
#            codes.
#  Author:   Frank Warmerdam, warmerdam@pobox.com
# 
#******************************************************************************
#  Copyright (c) 2001, Frank Warmerdam
# 
#  Permission is hereby granted, free of charge, to any person obtaining a
#  copy of this software and associated documentation files (the "Software"),
#  to deal in the Software without restriction, including without limitation
#  the rights to use, copy, modify, merge, publish, distribute, sublicense,
#  and/or sell copies of the Software, and to permit persons to whom the
#  Software is furnished to do so, subject to the following conditions:
# 
#  The above copyright notice and this permission notice shall be included
#  in all copies or substantial portions of the Software.
# 
#  THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
#  OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
#  FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
#  THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
#  LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
#  FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
#  DEALINGS IN THE SOFTWARE.
#******************************************************************************

try:
    from osgeo import osr
    from osgeo import gdal
except ImportError:
    import osr
    import gdal

import sys
import string

# =============================================================================
def Usage():

    print 'Usage: epsg_tr.py [-wkt] [-pretty_wkt] [-proj4] [-xml] [-postgis]'
    print '                  [-skip] [-list filename] [start_code [end_code]]'
    sys.exit(1)

# =============================================================================
def trHandleCode(code, gen_dict_line, report_error, output_format):

    import time

    try:
        err = prj_srs.ImportFromEPSG( code )
    except:
        err = 1

    if err != 0 and report_error:
        print 'Unable to lookup ',code,', either not a valid EPSG'
        print 'code, or it the EPSG csv files are not accessable.'
        sys.exit(2)
    else:
        if output_format == '-pretty_wkt':
            if gen_dict_line:
                print 'EPSG:',code

            print prj_srs.ExportToPrettyWkt()

        if output_format == '-xml':
            print prj_srs.ExportToXML()
            
        if output_format == '-wkt':
            if gen_dict_line:
                print 'EPSG:',code
                    
            print prj_srs.ExportToWkt()
                
        if output_format == '-proj4':
            out_string = prj_srs.ExportToProj4()

            name = prj_srs.GetAttrValue('PROJCS')
            if name is None:
                name = prj_srs.GetAttrValue('GEOGCS')

            if name is None:
                name = 'Unknown'
            
            print '# %s' % name
            if err == 0 and string.find(out_string,'+proj=') > -1:
                print '<%s> %s <>' % (str(code), out_string)
            else:
                print '# Unable to translate coordinate system EPSG:%d into PROJ.4 format.' % code
                print '#'

        if output_format == '-postgis':
            name = prj_srs.GetAttrValue('PROJCS')
            if name is None:
                name = prj_srs.GetAttrValue('GEOGCS')

            try:
                proj4text = prj_srs.ExportToProj4()
            except:
                err = 1
            wkt = prj_srs.ExportToWkt()
            
            print '---'
            print '--- EPSG %d : %s' % (code, name)
            print '---'

            if err:
                print '-- (unable to translate)'
            else:
                wkt = gdal.EscapeString(wkt,scheme=gdal.CPLES_SQL)
                proj4text = gdal.EscapeString(proj4text,scheme=gdal.CPLES_SQL)
                print 'INSERT INTO "spatial_ref_sys" ("srid","auth_name","auth_srid","srtext","proj4text") VALUES (%s,\'EPSG\',%s,\'%s\',\'%s\');' % \
                      (str(code),str(code),wkt,proj4text)
            
# =============================================================================

if __name__ == '__main__':

    start_code = -1
    end_code = -1
    list_file = None
    output_format = '-pretty_wkt'
    report_error = 1
    
    argv = gdal.GeneralCmdLineProcessor( sys.argv )
    if argv is None:
        sys.exit( 0 )
        
    # Parse command line arguments.
    
    i = 1
    while i < len(argv):
        arg = argv[i]

        if arg == '-wkt' or arg == '-pretty_wkt' or arg == '-proj4' \
           or arg == '-postgis' or arg == '-xml':
            output_format = arg

        elif arg[:5] == '-skip':
            report_error = 0
            
        elif arg == '-list' and i < len(argv)-1:
            i = i + 1
            list_file = argv[i]
            
        elif arg[0] == '-':
            Usage()

        elif int(arg) > 0:
            
            if start_code == -1:
                start_code = int(arg)
                end_code = int(arg)
            elif end_code == start_code:
                end_code = int(arg)
            else:
                Usage()
        else:
            Usage()

        i = i + 1

    # Output BEGIN transaction for PostGIS
    if output_format == '-postgis':
        print 'BEGIN;'

    # Do we need to produce a single output definition, or include a
    # dictionary line for each entry?

    gen_dict_line = start_code != end_code

    # loop over all codes to generate output
    
    prj_srs = osr.SpatialReference()

    if start_code != -1:
        for code in range(start_code,end_code+1):
            trHandleCode(code, gen_dict_line, report_error, output_format)

    # loop over codes read from list file.

    elif list_file is not None:

        list_fd = open( list_file )
        line = list_fd.readline()
        while len(line) > 0:
            try:
                c_offset = string.find(line,',')
                if c_offset > 0:
                    line = line[:c_offset]
                    
                code = string.atoi(line)
            except:
                code = -1

            if code <> -1:
                trHandleCode(code, gen_dict_line, report_error, output_format)
                
            line = list_fd.readline()

    else:
        Usage()
        
    # Output COMMIT transaction for PostGIS
    if output_format == '-postgis':
        print 'COMMIT;'
        print 'VACUUM ANALYZE spatial_ref_sys;'