PageRenderTime 84ms CodeModel.GetById 42ms app.highlight 35ms RepoModel.GetById 1ms app.codeStats 0ms

/Proj4/pj_gridinfo.c

http://github.com/route-me/route-me
C | 716 lines | 430 code | 119 blank | 167 comment | 122 complexity | f65adc5c2251229b627b73a48ba1943a MD5 | raw file
  1/******************************************************************************
  2 * $Id: pj_gridinfo.c,v 1.8 2006/11/17 22:16:30 mloskot Exp $
  3 *
  4 * Project:  PROJ.4
  5 * Purpose:  Functions for handling individual PJ_GRIDINFO's.  Includes
  6 *           loaders for all formats but CTABLE (in nad_init.c).
  7 * Author:   Frank Warmerdam, warmerdam@pobox.com
  8 *
  9 ******************************************************************************
 10 * Copyright (c) 2000, Frank Warmerdam <warmerdam@pobox.com>
 11 *
 12 * Permission is hereby granted, free of charge, to any person obtaining a
 13 * copy of this software and associated documentation files (the "Software"),
 14 * to deal in the Software without restriction, including without limitation
 15 * the rights to use, copy, modify, merge, publish, distribute, sublicense,
 16 * and/or sell copies of the Software, and to permit persons to whom the
 17 * Software is furnished to do so, subject to the following conditions:
 18 *
 19 * The above copyright notice and this permission notice shall be included
 20 * in all copies or substantial portions of the Software.
 21 *
 22 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
 23 * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
 24 * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
 25 * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
 26 * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
 27 * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
 28 * DEALINGS IN THE SOFTWARE.
 29 ******************************************************************************
 30 *
 31 * $Log: pj_gridinfo.c,v $
 32 * Revision 1.8  2006/11/17 22:16:30  mloskot
 33 * Uploaded PROJ.4 port for Windows CE.
 34 *
 35 * Revision 1.7  2005/07/07 00:16:03  fwarmerdam
 36 * Fixed debug fprintf syntax per:
 37 * http://bugzilla.remotesensing.org/show_bug.cgi?id=886
 38 *
 39 * Revision 1.6  2004/10/30 04:03:03  fwarmerdam
 40 * fixed reported information in ctable debug message
 41 *
 42 * Revision 1.5  2003/08/20 13:23:58  warmerda
 43 * Avoid unsigned char / char casting issues for VC++.
 44 *
 45 * Revision 1.4  2003/03/19 03:36:41  warmerda
 46 * Fixed so swap_words() works when it should.
 47 *
 48 * Revision 1.3  2003/03/17 19:44:45  warmerda
 49 * improved debugging, reduce header read size
 50 *
 51 * Revision 1.2  2003/03/17 18:56:34  warmerda
 52 * implement heirarchical NTv2 gridinfos
 53 *
 54 * Revision 1.1  2003/03/15 06:01:18  warmerda
 55 * New
 56 *
 57 */
 58
 59#define PJ_LIB__
 60
 61#include "projects.h"
 62#include <string.h>
 63#include <math.h>
 64#include <errno.h>
 65
 66#ifdef _WIN32_WCE
 67/* assert.h includes all Windows API headers and causes 'LP' name clash.
 68 * Here assert we disable assert() for Windows CE.
 69 * TODO - mloskot: re-implement porting friendly assert
 70 */
 71# define assert(exp)	((void)0)
 72#else
 73# include <assert.h>
 74#endif /* _WIN32_WCE */
 75
 76/************************************************************************/
 77/*                             swap_words()                             */
 78/*                                                                      */
 79/*      Convert the byte order of the given word(s) in place.           */
 80/************************************************************************/
 81
 82static int  byte_order_test = 1;
 83#define IS_LSB	(((unsigned char *) (&byte_order_test))[0] == 1)
 84
 85static void swap_words( unsigned char *data, int word_size, int word_count )
 86
 87{
 88    int	word;
 89
 90    for( word = 0; word < word_count; word++ )
 91    {
 92        int	i;
 93        
 94        for( i = 0; i < word_size/2; i++ )
 95        {
 96            int	t;
 97            
 98            t = data[i];
 99            data[i] = data[word_size-i-1];
100            data[word_size-i-1] = t;
101        }
102        
103        data += word_size;
104    }
105}
106
107/************************************************************************/
108/*                          pj_gridinfo_free()                          */
109/************************************************************************/
110
111void pj_gridinfo_free( PJ_GRIDINFO *gi )
112
113{
114    if( gi == NULL )
115        return;
116
117    if( gi->child != NULL )
118    {
119        PJ_GRIDINFO *child, *next;
120
121        for( child = gi->child; child != NULL; child=next)
122        {
123            next=child->next;
124            pj_gridinfo_free( child );
125        }
126    }
127
128    if( gi->ct != NULL )
129        nad_free( gi->ct );
130    
131    free( gi->gridname );
132    if( gi->filename != NULL )
133        free( gi->filename );
134
135    pj_dalloc( gi );
136}
137
138/************************************************************************/
139/*                          pj_gridinfo_load()                          */
140/*                                                                      */
141/*      This function is intended to implement delayed loading of       */
142/*      the data contents of a grid file.  The header and related       */
143/*      stuff are loaded by pj_gridinfo_init().                         */
144/************************************************************************/
145
146int pj_gridinfo_load( PJ_GRIDINFO *gi )
147
148{
149    if( gi == NULL || gi->ct == NULL )
150        return 0;
151
152/* -------------------------------------------------------------------- */
153/*      ctable is currently loaded on initialization though there is    */
154/*      no real reason not to support delayed loading for it as well.   */
155/* -------------------------------------------------------------------- */
156    if( strcmp(gi->format,"ctable") == 0 )
157    {
158        FILE *fid;
159        int result;
160
161        fid = pj_open_lib( gi->filename, "rb" );
162        
163        if( fid == NULL )
164        {
165            pj_errno = -38;
166            return 0;
167        }
168
169        result = nad_ctable_load( gi->ct, fid );
170
171        fclose( fid );
172
173        return result;
174    }
175
176/* -------------------------------------------------------------------- */
177/*      NTv1 format.                                                    */
178/*      We process one line at a time.  Note that the array storage     */
179/*      direction (e-w) is different in the NTv1 file and what          */
180/*      the CTABLE is supposed to have.  The phi/lam are also           */
181/*      reversed, and we have to be aware of byte swapping.             */
182/* -------------------------------------------------------------------- */
183    else if( strcmp(gi->format,"ntv1") == 0 )
184    {
185        double	*row_buf;
186        int	row;
187        FILE *fid;
188
189        fid = pj_open_lib( gi->filename, "rb" );
190        
191        if( fid == NULL )
192        {
193            pj_errno = -38;
194            return 0;
195        }
196
197        fseek( fid, gi->grid_offset, SEEK_SET );
198
199        row_buf = (double *) pj_malloc(gi->ct->lim.lam * sizeof(double) * 2);
200        gi->ct->cvs = (FLP *) pj_malloc(gi->ct->lim.lam*gi->ct->lim.phi*sizeof(FLP));
201        if( row_buf == NULL || gi->ct->cvs == NULL )
202        {
203            pj_errno = -38;
204            return 0;
205        }
206        
207        for( row = 0; row < gi->ct->lim.phi; row++ )
208        {
209            int	    i;
210            FLP     *cvs;
211            double  *diff_seconds;
212
213            if( fread( row_buf, sizeof(double), gi->ct->lim.lam * 2, fid ) 
214                != 2 * gi->ct->lim.lam )
215            {
216                pj_dalloc( row_buf );
217                pj_dalloc( gi->ct->cvs );
218                pj_errno = -38;
219                return 0;
220            }
221
222            if( IS_LSB )
223                swap_words( (unsigned char *) row_buf, 8, gi->ct->lim.lam*2 );
224
225            /* convert seconds to radians */
226            diff_seconds = row_buf;
227
228            for( i = 0; i < gi->ct->lim.lam; i++ )
229            {
230                cvs = gi->ct->cvs + (row) * gi->ct->lim.lam
231                    + (gi->ct->lim.lam - i - 1);
232
233                cvs->phi = *(diff_seconds++) * ((PI/180.0) / 3600.0);
234                cvs->lam = *(diff_seconds++) * ((PI/180.0) / 3600.0);
235            }
236        }
237
238        pj_dalloc( row_buf );
239
240        fclose( fid );
241
242        return 1;
243    }
244
245/* -------------------------------------------------------------------- */
246/*      NTv2 format.                                                    */
247/*      We process one line at a time.  Note that the array storage     */
248/*      direction (e-w) is different in the NTv2 file and what          */
249/*      the CTABLE is supposed to have.  The phi/lam are also           */
250/*      reversed, and we have to be aware of byte swapping.             */
251/* -------------------------------------------------------------------- */
252    else if( strcmp(gi->format,"ntv2") == 0 )
253    {
254        float	*row_buf;
255        int	row;
256        FILE *fid;
257
258        if( getenv("PROJ_DEBUG") != NULL )
259        {
260            fprintf( stderr, "NTv2 - loading grid %s\n", gi->ct->id );
261        }
262
263        fid = pj_open_lib( gi->filename, "rb" );
264        
265        if( fid == NULL )
266        {
267            pj_errno = -38;
268            return 0;
269        }
270
271        fseek( fid, gi->grid_offset, SEEK_SET );
272
273        row_buf = (float *) pj_malloc(gi->ct->lim.lam * sizeof(float) * 4);
274        gi->ct->cvs = (FLP *) pj_malloc(gi->ct->lim.lam*gi->ct->lim.phi*sizeof(FLP));
275        if( row_buf == NULL || gi->ct->cvs == NULL )
276        {
277            pj_errno = -38;
278            return 0;
279        }
280        
281        for( row = 0; row < gi->ct->lim.phi; row++ )
282        {
283            int	    i;
284            FLP     *cvs;
285            float   *diff_seconds;
286
287            if( fread( row_buf, sizeof(float), gi->ct->lim.lam*4, fid ) 
288                != 4 * gi->ct->lim.lam )
289            {
290                pj_dalloc( row_buf );
291                pj_dalloc( gi->ct->cvs );
292                gi->ct->cvs = NULL;
293                pj_errno = -38;
294                return 0;
295            }
296
297            if( !IS_LSB )
298                swap_words( (unsigned char *) row_buf, 4, 
299                            gi->ct->lim.lam*4 );
300
301            /* convert seconds to radians */
302            diff_seconds = row_buf;
303
304            for( i = 0; i < gi->ct->lim.lam; i++ )
305            {
306                cvs = gi->ct->cvs + (row) * gi->ct->lim.lam
307                    + (gi->ct->lim.lam - i - 1);
308
309                cvs->phi = *(diff_seconds++) * ((PI/180.0) / 3600.0);
310                cvs->lam = *(diff_seconds++) * ((PI/180.0) / 3600.0);
311                diff_seconds += 2; /* skip accuracy values */
312            }
313        }
314
315        pj_dalloc( row_buf );
316
317        fclose( fid );
318
319        return 1;
320    }
321
322    else
323    {
324        return 0;
325    }
326}
327
328/************************************************************************/
329/*                       pj_gridinfo_init_ntv2()                        */
330/*                                                                      */
331/*      Load a ntv2 (.gsb) file.                                        */
332/************************************************************************/
333
334static int pj_gridinfo_init_ntv2( FILE *fid, PJ_GRIDINFO *gilist )
335
336{
337    unsigned char header[11*16];
338    int num_subfiles, subfile;
339
340    assert( sizeof(int) == 4 );
341    assert( sizeof(double) == 8 );
342    if( sizeof(int) != 4 || sizeof(double) != 8 )
343    {
344        fprintf( stderr, 
345                 "basic types of inappropraiate size in pj_gridinfo_init_ntv2()\n" );
346        pj_errno = -38;
347        return 0;
348    }
349
350/* -------------------------------------------------------------------- */
351/*      Read the overview header.                                       */
352/* -------------------------------------------------------------------- */
353    if( fread( header, sizeof(header), 1, fid ) != 1 )
354    {
355        pj_errno = -38;
356        return 0;
357    }
358
359/* -------------------------------------------------------------------- */
360/*      Byte swap interesting fields if needed.                         */
361/* -------------------------------------------------------------------- */
362    if( !IS_LSB )
363    {
364        swap_words( header+8, 4, 1 );
365        swap_words( header+8+16, 4, 1 );
366        swap_words( header+8+32, 4, 1 );
367        swap_words( header+8+7*16, 8, 1 );
368        swap_words( header+8+8*16, 8, 1 );
369        swap_words( header+8+9*16, 8, 1 );
370        swap_words( header+8+10*16, 8, 1 );
371    }
372
373/* -------------------------------------------------------------------- */
374/*      Get the subfile count out ... all we really use for now.        */
375/* -------------------------------------------------------------------- */
376    memcpy( &num_subfiles, header+8+32, 4 );
377
378/* ==================================================================== */
379/*      Step through the subfiles, creating a PJ_GRIDINFO for each.     */
380/* ==================================================================== */
381    for( subfile = 0; subfile < num_subfiles; subfile++ )
382    {
383        struct CTABLE *ct;
384        LP ur;
385        int gs_count;
386        PJ_GRIDINFO *gi;
387
388/* -------------------------------------------------------------------- */
389/*      Read header.                                                    */
390/* -------------------------------------------------------------------- */
391        if( fread( header, sizeof(header), 1, fid ) != 1 )
392        {
393            pj_errno = -38;
394            return 0;
395        }
396
397        if( strncmp((const char *) header,"SUB_NAME",8) != 0 )
398        {
399            pj_errno = -38;
400            return 0;
401        }
402        
403/* -------------------------------------------------------------------- */
404/*      Byte swap interesting fields if needed.                         */
405/* -------------------------------------------------------------------- */
406        if( !IS_LSB )
407        {
408            swap_words( header+8+16*4, 8, 1 );
409            swap_words( header+8+16*5, 8, 1 );
410            swap_words( header+8+16*6, 8, 1 );
411            swap_words( header+8+16*7, 8, 1 );
412            swap_words( header+8+16*8, 8, 1 );
413            swap_words( header+8+16*9, 8, 1 );
414            swap_words( header+8+16*10, 4, 1 );
415        }
416        
417/* -------------------------------------------------------------------- */
418/*      Initialize a corresponding "ct" structure.                      */
419/* -------------------------------------------------------------------- */
420        ct = (struct CTABLE *) pj_malloc(sizeof(struct CTABLE));
421        strncpy( ct->id, (const char *) header + 8, 8 );
422        ct->id[8] = '\0';
423
424        ct->ll.lam = - *((double *) (header+7*16+8)); /* W_LONG */
425        ct->ll.phi = *((double *) (header+4*16+8));   /* S_LAT */
426
427        ur.lam = - *((double *) (header+6*16+8));     /* E_LONG */
428        ur.phi = *((double *) (header+5*16+8));       /* N_LAT */
429
430        ct->del.lam = *((double *) (header+9*16+8));
431        ct->del.phi = *((double *) (header+8*16+8));
432
433        ct->lim.lam = (int) (fabs(ur.lam-ct->ll.lam)/ct->del.lam + 0.5) + 1;
434        ct->lim.phi = (int) (fabs(ur.phi-ct->ll.phi)/ct->del.phi + 0.5) + 1;
435
436        if( getenv("PROJ_DEBUG") != NULL )
437            fprintf( stderr, 
438                     "NTv2 %s %dx%d: LL=(%.9g,%.9g) UR=(%.9g,%.9g)\n",
439                     ct->id, 
440                     ct->lim.lam, ct->lim.phi,
441                     ct->ll.lam/3600.0, ct->ll.phi/3600.0,
442                     ur.lam/3600.0, ur.phi/3600.0 );
443
444        ct->ll.lam *= DEG_TO_RAD/3600.0;
445        ct->ll.phi *= DEG_TO_RAD/3600.0;
446        ct->del.lam *= DEG_TO_RAD/3600.0;
447        ct->del.phi *= DEG_TO_RAD/3600.0;
448
449        memcpy( &gs_count, header + 8 + 16*10, 4 );
450        if( gs_count != ct->lim.lam * ct->lim.phi )
451        {
452            fprintf( stderr, 
453                     "GS_COUNT(%d) does not match expected cells (%dx%d=%d)\n",
454                     gs_count, ct->lim.lam, ct->lim.phi, 
455                     ct->lim.lam * ct->lim.phi );
456            pj_errno = -38;
457            return 0;
458        }
459
460        ct->cvs = NULL;
461
462/* -------------------------------------------------------------------- */
463/*      Create a new gridinfo for this if we aren't processing the      */
464/*      1st subfile, and initialize our grid info.                      */
465/* -------------------------------------------------------------------- */
466        if( subfile == 0 )
467            gi = gilist;
468        else
469        {
470            gi = (PJ_GRIDINFO *) pj_malloc(sizeof(PJ_GRIDINFO));
471            memset( gi, 0, sizeof(PJ_GRIDINFO) );
472    
473            gi->gridname = strdup( gilist->gridname );
474            gi->filename = strdup( gilist->filename );
475            gi->next = NULL;
476        }
477
478        gi->ct = ct;
479        gi->format = "ntv2";
480        gi->grid_offset = ftell( fid );
481
482/* -------------------------------------------------------------------- */
483/*      Attach to the correct list or sublist.                          */
484/* -------------------------------------------------------------------- */
485        if( strncmp((const char *)header+24,"NONE",4) == 0 )
486        {
487            if( gi != gilist )
488            {
489                PJ_GRIDINFO *lnk;
490
491                for( lnk = gilist; lnk->next != NULL; lnk = lnk->next ) {}
492                lnk->next = gi;
493            }
494        }
495
496        else
497        {
498            PJ_GRIDINFO *lnk;
499            PJ_GRIDINFO *gp = gilist;
500            
501            while( gp != NULL 
502                   && strncmp(gp->ct->id,(const char*)header+24,8) != 0 )
503                gp = gp->next;
504
505            if( gp == NULL )
506            {
507                if( getenv("PROJ_DEBUG") != NULL )
508                    fprintf( stderr, "pj_gridinfo_init_ntv2(): "
509                             "failed to find parent %8.8s for %s.\n", 
510                             (const char *) header+24, gi->ct->id );
511
512                if (gp) {
513			     for( lnk = gp; lnk->next != NULL; lnk = lnk->next ) {}
514                    lnk->next = gi;
515                }
516            }
517            else if( gp->child == NULL )
518            {
519                gp->child = gi;
520            }
521            else
522            {
523                for( lnk = gp->child; lnk->next != NULL; lnk = lnk->next ) {}
524                lnk->next = gi;
525            }
526        }
527
528/* -------------------------------------------------------------------- */
529/*      Seek past the data.                                             */
530/* -------------------------------------------------------------------- */
531        fseek( fid, gs_count * 16, SEEK_CUR );
532    }
533
534    return 1;
535}
536
537/************************************************************************/
538/*                       pj_gridinfo_init_ntv1()                        */
539/*                                                                      */
540/*      Load an NTv1 style Canadian grid shift file.                    */
541/************************************************************************/
542
543static int pj_gridinfo_init_ntv1( FILE * fid, PJ_GRIDINFO *gi )
544
545{
546    unsigned char header[176];
547    struct CTABLE *ct;
548    LP		ur;
549    
550    assert( sizeof(int) == 4 );
551    assert( sizeof(double) == 8 );
552    if( sizeof(int) != 4 || sizeof(double) != 8 )
553    {
554        fprintf( stderr, 
555                 "basic types of inappropraiate size in nad_load_ntv1()\n" );
556        pj_errno = -38;
557        return 0;
558    }
559
560/* -------------------------------------------------------------------- */
561/*      Read the header.                                                */
562/* -------------------------------------------------------------------- */
563    if( fread( header, sizeof(header), 1, fid ) != 1 )
564    {
565        pj_errno = -38;
566        return 0;
567    }
568
569/* -------------------------------------------------------------------- */
570/*      Regularize fields of interest.                                  */
571/* -------------------------------------------------------------------- */
572    if( IS_LSB )
573    {
574        swap_words( header+8, 4, 1 );
575        swap_words( header+24, 8, 1 );
576        swap_words( header+40, 8, 1 );
577        swap_words( header+56, 8, 1 );
578        swap_words( header+72, 8, 1 );
579        swap_words( header+88, 8, 1 );
580        swap_words( header+104, 8, 1 );
581    }
582
583    if( *((int *) (header+8)) != 12 )
584    {
585        pj_errno = -38;
586        printf("NTv1 grid shift file has wrong record count, corrupt?\n");
587        return 0;
588    }
589
590/* -------------------------------------------------------------------- */
591/*      Fill in CTABLE structure.                                       */
592/* -------------------------------------------------------------------- */
593    ct = (struct CTABLE *) pj_malloc(sizeof(struct CTABLE));
594    strcpy( ct->id, "NTv1 Grid Shift File" );
595
596    ct->ll.lam = - *((double *) (header+72));
597    ct->ll.phi = *((double *) (header+24));
598    ur.lam = - *((double *) (header+56));
599    ur.phi = *((double *) (header+40));
600    ct->del.lam = *((double *) (header+104));
601    ct->del.phi = *((double *) (header+88));
602    ct->lim.lam = (int) (fabs(ur.lam-ct->ll.lam)/ct->del.lam + 0.5) + 1;
603    ct->lim.phi = (int) (fabs(ur.phi-ct->ll.phi)/ct->del.phi + 0.5) + 1;
604
605    if( getenv("PROJ_DEBUG") != NULL )
606        fprintf( stderr, 
607                 "NTv1 %dx%d: LL=(%.9g,%.9g) UR=(%.9g,%.9g)\n",
608                 ct->lim.lam, ct->lim.phi,
609                 ct->ll.lam, ct->ll.phi, ur.lam, ur.phi );
610
611    ct->ll.lam *= DEG_TO_RAD;
612    ct->ll.phi *= DEG_TO_RAD;
613    ct->del.lam *= DEG_TO_RAD;
614    ct->del.phi *= DEG_TO_RAD;
615    ct->cvs = NULL;
616
617    gi->ct = ct;
618    gi->grid_offset = ftell( fid );
619    gi->format = "ntv1";
620
621    return 1;
622}
623
624/************************************************************************/
625/*                          pj_gridinfo_init()                          */
626/*                                                                      */
627/*      Open and parse header details from a datum gridshift file       */
628/*      returning a list of PJ_GRIDINFOs for the grids in that          */
629/*      file.  This superceeds use of nad_init() for modern             */
630/*      applications.                                                   */
631/************************************************************************/
632
633PJ_GRIDINFO *pj_gridinfo_init( const char *gridname )
634
635{
636    char 	fname[MAX_PATH_FILENAME+1];
637    PJ_GRIDINFO *gilist;
638    FILE 	*fp;
639    char	header[160];
640
641    errno = pj_errno = 0;
642
643/* -------------------------------------------------------------------- */
644/*      Initialize a GRIDINFO with stub info we would use if it         */
645/*      cannot be loaded.                                               */
646/* -------------------------------------------------------------------- */
647    gilist = (PJ_GRIDINFO *) pj_malloc(sizeof(PJ_GRIDINFO));
648    memset( gilist, 0, sizeof(PJ_GRIDINFO) );
649    
650    gilist->gridname = strdup( gridname );
651    gilist->filename = NULL;
652    gilist->format = "missing";
653    gilist->grid_offset = 0;
654    gilist->ct = NULL;
655    gilist->next = NULL;
656
657/* -------------------------------------------------------------------- */
658/*      Open the file using the usual search rules.                     */
659/* -------------------------------------------------------------------- */
660    strcpy(fname, gridname);
661    if (!(fp = pj_open_lib(fname, "rb"))) {
662        pj_errno = errno;
663        return gilist;
664    }
665
666    gilist->filename = strdup(fname);
667    
668/* -------------------------------------------------------------------- */
669/*      Load a header, to determine the file type.                      */
670/* -------------------------------------------------------------------- */
671    if( fread( header, sizeof(header), 1, fp ) != 1 )
672    {
673        fclose( fp );
674        pj_errno = -38;
675        return gilist;
676    }
677
678    fseek( fp, SEEK_SET, 0 );
679
680/* -------------------------------------------------------------------- */
681/*      Determine file type.                                            */
682/* -------------------------------------------------------------------- */
683    if( strncmp(header + 0, "HEADER", 6) == 0 
684        && strncmp(header + 96, "W GRID", 6) == 0 
685        && strncmp(header + 144, "TO      NAD83   ", 16) == 0 )
686    {
687        pj_gridinfo_init_ntv1( fp, gilist );
688    }
689    
690    else if( strncmp(header + 0, "NUM_OREC", 8) == 0 
691             && strncmp(header + 48, "GS_TYPE", 7) == 0 )
692    {
693        pj_gridinfo_init_ntv2( fp, gilist );
694    }
695    
696    else
697    {
698        struct CTABLE *ct = nad_ctable_init( fp );
699
700        gilist->format = "ctable";
701        gilist->ct = ct;
702
703        if( getenv("PROJ_DEBUG") != NULL )
704            fprintf( stderr, 
705                     "Ctable %s %dx%d: LL=(%.9g,%.9g) UR=(%.9g,%.9g)\n",
706                     ct->id, 
707                     ct->lim.lam, ct->lim.phi,
708                     ct->ll.lam * RAD_TO_DEG, ct->ll.phi * RAD_TO_DEG,
709                     (ct->ll.lam + (ct->lim.lam-1)*ct->del.lam) * RAD_TO_DEG, 
710                     (ct->ll.phi + (ct->lim.phi-1)*ct->del.phi) * RAD_TO_DEG );
711    }
712
713    fclose(fp);
714
715    return gilist;
716}