implement support for ctable2 format (read/write)

git-svn-id: http://svn.osgeo.org/metacrs/proj/trunk@2121 4e78687f-474d-0410-85f9-8d5e500ac6b2

4 files changed, 278 insertions, 20 deletions

diff --git a/src/nad2bin.c b/src/nad2bin.c
index 6db9fe67..807859ba 100644
--- a/src/nad2bin.c
+++ b/src/nad2bin.c
@@ -6,8 +6,6 @@
 #define PJ_LIB__
 #include <projects.h>
 #define U_SEC_TO_RAD 4.848136811095359935899141023e-12
-	static char
-*usage = "<ASCII_dist_table local_bin_table";
 
 /************************************************************************/
 /*                             swap_words()                             */
@@ -18,10 +16,11 @@
 static int  byte_order_test = 1;
 #define IS_LSB	(((unsigned char *) (&byte_order_test))[0] == 1)
 
-static void swap_words( unsigned char *data, int word_size, int word_count )
+static void swap_words( void *data_in, int word_size, int word_count )
 
 {
     int	word;
+    unsigned char *data = (unsigned char *) data_in;
 
     for( word = 0; word < word_count; word++ )
     {
@@ -41,6 +40,17 @@ static void swap_words( unsigned char *data, int word_size, int word_count )
 }
 
 /************************************************************************/
+/*                               Usage()                                */
+/************************************************************************/
+
+static void Usage()
+{
+    fprintf(stderr,
+            "usage: nad2bin [-f ctable/ctable2/ntv2] binary_output < ascii_source\n" );
+    exit(1);
+}
+
+/************************************************************************/
 /*                                main()                                */
 /************************************************************************/
 int main(int argc, char **argv) {
@@ -51,7 +61,9 @@ int main(int argc, char **argv) {
     long lam, laml, phi, phil;
     FILE *fp;
 
-    const char *format   = "ntv2";
+    const char *output_file = NULL;
+
+    const char *format   = "ctable2";
     const char *GS_TYPE  = "SECONDS";
     const char *VERSION  = "";
     const char *SYSTEM_F = "NAD27";
@@ -60,11 +72,28 @@ int main(int argc, char **argv) {
     const char *CREATED  = "";
     const char *UPDATED  = "";
 
-    if (argc != 2) {
-        fprintf(stderr,"usage: %s %s\n", argv[0], usage);
-        exit(1);
+/* ==================================================================== */
+/*      Process arguments.                                              */
+/* ==================================================================== */
+    for( i = 1; i < argc; i++ )
+    {
+        if( strcmp(argv[i],"-f") && i < argc-1 ) 
+        {
+            format = argv[++i];
+        }
+        else if( output_file == NULL )
+        {
+            output_file = argv[i];
+        }
+        else
+            Usage();
     }
 
+    if( output_file == NULL )
+        Usage();
+
+    fprintf( stdout, "Output Binary File Format: %s\n", format );
+
 /* ==================================================================== */
 /*      Read the ASCII Table                                            */
 /* ==================================================================== */
@@ -109,8 +138,8 @@ int main(int argc, char **argv) {
 /* ==================================================================== */
     if( strcmp(format,"ctable") == 0 ) 
     {
-	if (!(fp = fopen(argv[1], "wb"))) {
-            perror(argv[1]);
+	if (!(fp = fopen(output_file, "wb"))) {
+            perror(output_file);
             exit(2);
 	}
 	if (fwrite(&ct, sizeof(ct), 1, fp) != 1 ||
@@ -123,13 +152,62 @@ int main(int argc, char **argv) {
     }
 
 /* ==================================================================== */
+/*      Write out the old ctable format - this is machine and byte      */
+/*      order specific.                                                 */
+/* ==================================================================== */
+    if( strcmp(format,"ctable2") == 0 ) 
+    {
+        char header[160];
+
+	if (!(fp = fopen(output_file, "wb"))) {
+            perror(output_file);
+            exit(2);
+	}
+
+        assert( MAX_TAB_ID == 80 );
+        assert( sizeof(int) == 4 ); /* for ct.lim.lam/phi */
+
+        memset( header, 0, sizeof(header) );
+
+        memcpy( header +   0, "CTABLE V2.0     ", 16 );
+        memcpy( header +  16, ct.id, 80 );
+        memcpy( header +  96, &ct.ll.lam, 8 );
+        memcpy( header + 104, &ct.ll.phi, 8 );
+        memcpy( header + 112, &ct.del.lam, 8 );
+        memcpy( header + 120, &ct.del.phi, 8 );
+        memcpy( header + 128, &ct.lim.lam, 4 );
+        memcpy( header + 132, &ct.lim.phi, 4 );
+
+        /* force into LSB format */
+        if( !IS_LSB ) 
+        {
+            swap_words( header +  96, 8, 4 );
+            swap_words( header + 128, 4, 2 );
+            swap_words( ct.cvs, 4, ct.lim.lam * ct.lim.phi );
+        }
+
+        if( fwrite( header, sizeof(header), 1, fp ) != 1 ) {
+            perror( "fwrite" );
+            exit( 2 );
+        }
+
+	if (fwrite(ct.cvs, tsize, 1, fp) != 1) {
+            perror( "fwrite" );
+            exit(2);
+	}
+
+        fclose( fp );
+	exit(0); /* normal completion */
+    }
+
+/* ==================================================================== */
 /*      Write out the NTv2 format grid shift file.                      */
 /* ==================================================================== */
     if( strcmp(format,"ntv2") == 0 ) 
     {
-        if (!(fp = fopen(argv[1], "wb"))) 
+        if (!(fp = fopen(output_file, "wb"))) 
         {
-            perror(argv[1]);
+            perror(output_file);
             exit(2);
         }
         
@@ -266,7 +344,7 @@ int main(int argc, char **argv) {
                 }
 
                 if( !IS_LSB )
-                    swap_words( (unsigned char *) row_buf, 4, ct.lim.lam * 4 );
+                    swap_words( row_buf, 4, ct.lim.lam * 4 );
 
                 if( fwrite( row_buf, sizeof(float), ct.lim.lam*4, fp ) 
                     != 4 * ct.lim.lam )
@@ -280,4 +358,7 @@ int main(int argc, char **argv) {
         fclose( fp );
         exit(0); /* normal completion */
     }
+
+    fprintf( stderr, "Unsupported format, nothing written.\n" );
+    exit( 3 );
 }
diff --git a/src/nad_init.c b/src/nad_init.c
index da10dd2b..20f1adb6 100644
--- a/src/nad_init.c
+++ b/src/nad_init.c
@@ -45,6 +45,38 @@
 #endif /* _WIN32_WCE */
 
 /************************************************************************/
+/*                             swap_words()                             */
+/*                                                                      */
+/*      Convert the byte order of the given word(s) in place.           */
+/************************************************************************/
+
+static int  byte_order_test = 1;
+#define IS_LSB	(((unsigned char *) (&byte_order_test))[0] == 1)
+
+static void swap_words( void *data_in, int word_size, int word_count )
+
+{
+    int	word;
+    unsigned char *data = (unsigned char *) data_in;
+
+    for( word = 0; word < word_count; word++ )
+    {
+        int	i;
+        
+        for( i = 0; i < word_size/2; i++ )
+        {
+            int	t;
+            
+            t = data[i];
+            data[i] = data[word_size-i-1];
+            data[word_size-i-1] = t;
+        }
+        
+        data += word_size;
+    }
+}
+
+/************************************************************************/
 /*                          nad_ctable_load()                           */
 /*                                                                      */
 /*      Load the data portion of a ctable formatted grid.               */
@@ -66,6 +98,79 @@ int nad_ctable_load( projCtx ctx, struct CTABLE *ct, FILE *fid )
         pj_dalloc( ct->cvs );
         ct->cvs = NULL;
 
+        pj_log( ctx, PJ_LOG_ERROR, 
+                "ctable loading failed on fread() - binary incompatible?\n" );
+        pj_ctx_set_errno( ctx, -38 );
+        return 0;
+    }
+
+    return 1;
+} 
+
+/************************************************************************/
+/*                          nad_ctable_init()                           */
+/*                                                                      */
+/*      Read the header portion of a "ctable" format grid.              */
+/************************************************************************/
+
+struct CTABLE *nad_ctable_init( projCtx ctx, FILE * fid )
+{
+    struct CTABLE *ct;
+    int		id_end;
+
+    /* read the table header */
+    ct = (struct CTABLE *) pj_malloc(sizeof(struct CTABLE));
+    if( ct == NULL 
+        || fread( ct, sizeof(struct CTABLE), 1, fid ) != 1 )
+    {
+        pj_ctx_set_errno( ctx, -38 );
+        return NULL;
+    }
+
+    /* do some minimal validation to ensure the structure isn't corrupt */
+    if( ct->lim.lam < 1 || ct->lim.lam > 100000 
+        || ct->lim.phi < 1 || ct->lim.phi > 100000 )
+    {
+        pj_ctx_set_errno( ctx, -38 );
+        return NULL;
+    }
+    
+    /* trim white space and newlines off id */
+    for( id_end = strlen(ct->id)-1; id_end > 0; id_end-- )
+    {
+        if( ct->id[id_end] == '\n' || ct->id[id_end] == ' ' )
+            ct->id[id_end] = '\0';
+        else
+            break;
+    }
+
+    ct->cvs = NULL;
+
+    return ct;
+}
+
+/************************************************************************/
+/*                          nad_ctable2_load()                          */
+/*                                                                      */
+/*      Load the data portion of a ctable2 formatted grid.              */
+/************************************************************************/
+
+int nad_ctable2_load( projCtx ctx, struct CTABLE *ct, FILE *fid )
+
+{
+    int  a_size;
+
+    fseek( fid, 160, SEEK_SET );
+
+    /* read all the actual shift values */
+    a_size = ct->lim.lam * ct->lim.phi;
+    ct->cvs = (FLP *) pj_malloc(sizeof(FLP) * a_size);
+    if( ct->cvs == NULL 
+        || fread(ct->cvs, sizeof(FLP), a_size, fid) != a_size )
+    {
+        pj_dalloc( ct->cvs );
+        ct->cvs = NULL;
+
         if( getenv("PROJ_DEBUG") != NULL )
         {
             fprintf( stderr, 
@@ -76,29 +181,61 @@ int nad_ctable_load( projCtx ctx, struct CTABLE *ct, FILE *fid )
         return 0;
     }
 
+    if( !IS_LSB )
+    {
+        swap_words( ct->cvs, 4, a_size * 2 );
+    }
+
     return 1;
 } 
 
 /************************************************************************/
-/*                          nad_ctable_init()                           */
+/*                          nad_ctable2_init()                          */
 /*                                                                      */
-/*      Read the header portion of a "ctable" format grid.              */
+/*      Read the header portion of a "ctable2" format grid.             */
 /************************************************************************/
 
-struct CTABLE *nad_ctable_init( projCtx ctx, FILE * fid )
+struct CTABLE *nad_ctable2_init( projCtx ctx, FILE * fid )
 {
     struct CTABLE *ct;
     int		id_end;
+    char        header[160];
+
+    if( fread( header, sizeof(header), 1, fid ) != 1 )
+    {
+        pj_ctx_set_errno( ctx, -38 );
+        return NULL;
+    }
+
+    if( !IS_LSB )
+    {
+        swap_words( header +  96, 8, 4 );
+        swap_words( header + 128, 4, 2 );
+    }
+
+    if( strncmp(header,"CTABLE V2",9) != 0 )
+    {
+        pj_log( ctx, PJ_LOG_ERROR, "ctable2 - wrong header!" );
+        pj_ctx_set_errno( ctx, -38 );
+        return NULL;
+    }
 
     /* read the table header */
     ct = (struct CTABLE *) pj_malloc(sizeof(struct CTABLE));
-    if( ct == NULL 
-        || fread( ct, sizeof(struct CTABLE), 1, fid ) != 1 )
+    if( ct == NULL )
     {
         pj_ctx_set_errno( ctx, -38 );
         return NULL;
     }
 
+    memcpy( ct->id,       header +  16, 80 );
+    memcpy( &ct->ll.lam,  header +  96, 8 );
+    memcpy( &ct->ll.phi,  header + 104, 8 );
+    memcpy( &ct->del.lam, header + 112, 8 );
+    memcpy( &ct->del.phi, header + 120, 8 );
+    memcpy( &ct->lim.lam, header + 128, 4 );
+    memcpy( &ct->lim.phi, header + 132, 4 );
+
     /* do some minimal validation to ensure the structure isn't corrupt */
     if( ct->lim.lam < 1 || ct->lim.lam > 100000 
         || ct->lim.phi < 1 || ct->lim.phi > 100000 )
@@ -142,7 +279,7 @@ struct CTABLE *nad_init(projCtx ctx, char *name)
     if (!(fid = pj_open_lib(ctx, fname, "rb"))) {
         return 0;
     }
-    
+
     ct = nad_ctable_init( ctx, fid );
     if( ct != NULL )
     {
diff --git a/src/pj_gridinfo.c b/src/pj_gridinfo.c
index 34e65fb2..cba1e672 100644
--- a/src/pj_gridinfo.c
+++ b/src/pj_gridinfo.c
@@ -122,8 +122,7 @@ int pj_gridinfo_load( projCtx ctx, PJ_GRIDINFO *gi )
         return 0;
 
 /* -------------------------------------------------------------------- */
-/*      ctable is currently loaded on initialization though there is    */
-/*      no real reason not to support delayed loading for it as well.   */
+/*      Original platform specific CTable format.                       */
 /* -------------------------------------------------------------------- */
     if( strcmp(gi->format,"ctable") == 0 )
     {
@@ -146,6 +145,29 @@ int pj_gridinfo_load( projCtx ctx, PJ_GRIDINFO *gi )
     }
 
 /* -------------------------------------------------------------------- */
+/*      CTable2 format.                                                 */
+/* -------------------------------------------------------------------- */
+    else if( strcmp(gi->format,"ctable2") == 0 )
+    {
+        FILE *fid;
+        int result;
+
+        fid = pj_open_lib( ctx, gi->filename, "rb" );
+        
+        if( fid == NULL )
+        {
+            pj_ctx_set_errno( ctx, -38 );
+            return 0;
+        }
+
+        result = nad_ctable2_load( ctx, gi->ct, fid );
+
+        fclose( fid );
+
+        return result;
+    }
+
+/* -------------------------------------------------------------------- */
 /*      NTv1 format.                                                    */
 /*      We process one line at a time.  Note that the array storage     */
 /*      direction (e-w) is different in the NTv1 file and what          */
@@ -810,6 +832,22 @@ PJ_GRIDINFO *pj_gridinfo_init( projCtx ctx, const char *gridname )
         pj_gridinfo_init_gtx( ctx, fp, gilist );
     }
 
+    else if( strncmp(header+0,"CTABLE V2",9) == 0 )
+    {
+        struct CTABLE *ct = nad_ctable2_init( ctx, fp );
+
+        gilist->format = "ctable2";
+        gilist->ct = ct;
+
+        pj_log( ctx, PJ_LOG_DEBUG_MAJOR, 
+                "Ctable2 %s %dx%d: LL=(%.9g,%.9g) UR=(%.9g,%.9g)\n",
+                ct->id, 
+                ct->lim.lam, ct->lim.phi,
+                ct->ll.lam * RAD_TO_DEG, ct->ll.phi * RAD_TO_DEG,
+                (ct->ll.lam + (ct->lim.lam-1)*ct->del.lam) * RAD_TO_DEG, 
+                (ct->ll.phi + (ct->lim.phi-1)*ct->del.phi) * RAD_TO_DEG );
+    }
+
     else
     {
         struct CTABLE *ct = nad_ctable_init( ctx, fp );
diff --git a/src/projects.h b/src/projects.h
index 7c41242b..fb49652d 100644
--- a/src/projects.h
+++ b/src/projects.h
@@ -416,6 +416,8 @@ LP nad_cvt(LP, int, struct CTABLE *);
 struct CTABLE *nad_init(projCtx ctx, char *);
 struct CTABLE *nad_ctable_init( projCtx ctx, FILE * fid );
 int nad_ctable_load( projCtx ctx, struct CTABLE *, FILE * fid );
+struct CTABLE *nad_ctable2_init( projCtx ctx, FILE * fid );
+int nad_ctable2_load( projCtx ctx, struct CTABLE *, FILE * fid );
 void nad_free(struct CTABLE *);
 
 /* higher level handling of datum grid shift files */