gdapi.c

GIS系统支持库Geospatial Data Abstraction Library代码.GDAL is a translator library for raster geospatial dat

/*
Copyright (C) 1996 Hughes and Applied Research Corporation

Permission to use, modify, and distribute this software and its documentation 
for any purpose without fee is hereby granted, provided that the above 
copyright notice appear in all copies and that both that copyright notice and 
this permission notice appear in supporting documentation.
*/
/*****************************************************************************
REVISIONS:

Aug 31, 1999  Abe Taaheri    Changed memory allocation for utility strings to
                             the size of UTLSTR_MAX_SIZE.
			     Added error check for memory unavailibilty in 
			     several functions.
			     Added check for NULL metabuf returned from 
			     EHmeta... functions. NULL pointer retruned from 
			     EHmeta... functions indicate that memory could not
			     be allocated for metabuf.
Jun  27, 2000  Abe Taaheri   Added support for EASE grid that uses 
                             Behrmann Cylinderical Equal Area (BCEA) projection
Oct  23, 2000  Abe Taaheri   Updated for ISINUS projection, so that both codes
                             31 and 99 can be used for this projection.
Jan  15, 2003   Abe Taaheri  Modified for generalization of EASE Grid. 

Jun  05, 2003 Bruce Beaumont / Abe Taaheri

                             Fixed SQUARE definition.
                             Added static projection number/name translation 
                             Added projection table lookup in GDdefproj.
                             Removed projection table from GDdefproj
                             Added projection table lookup in GDprojinfo
                             Removed projection table from GDprojinfo
                             Added cast for compcode in call to SDsetcompress
                                in GDdeffield to avoid compiler errors
                             Removed declaration for unused variable endptr 
                                in GDSDfldsrch
                             Removed initialization code for unused variables
                                in GDSDfldsrch
                             Removed declarations for unused variables 
                                BCEA_scale, r0, s0, xMtr0, xMtr1, yMtr0, 
                                and yMtr1 in GDll2ij
                             Removed initialization code for unused variables
                                in GDll2ij
                             Added code in GEO projection handling to allow 
                                map to span dateline in GDll2ij
                             Changed "for each point" loop in GDll2ij to 
                                return -2147483648.0 for xVal and yVal if 
                                for_trans returned an error instead of 
                                returning an error to the caller
                                (Note: MAXLONG is defined as 2147483647.0 in
                                 function cproj.c of GCTP)
                             Added code in GDij2ll to use for_trans to 
                                translate the BCEA corner points from packed 
                                degrees to meters
                             Removed declarations for unused variables 
                                BCEA_scale, r0, s0, xMtr, yMtr, epsilon, 
                                beta, qp_cea, kz_cea, eccen, eccen_sq, 
                                phi1, sinphi1, cosphi1, lon, lat, xcor, 
                                ycor, and nlatlon from GDij2ll
                             Removed initialization code for unused variables
                                in GDij2ll
                             Added declarations for xMtr0, yMtr0, xMtr1, and 
                                yMtr1 in GDij2ll
                             Added special-case code for BCEA
                             Changed "for each point" loop in GDij2ll to 
                                return PGSd_GCT_IN_ERROR (1.0e51) for 
                                longitude and latitude values if inv_trans 
                                returned an error instead of return an error 
                                to the caller
                             Removed declaration for unused variable ii in 
                                GDgetpixvalues
                             Removed declaration for unused variable 
                                numTileDims in GDtileinfo
                             Added error message and error return at the 
                                end of GDll2mm_cea
                             Added return statement to GDll2mm_cea
******************************************************************************/
#include "stdio.h"
#include "mfhdf.h"
#include "hcomp.h"
#include <math.h>
#include "cfortHdf.h"
#include "HdfEosDef.h"

#define	GDIDOFFSET 4194304
#define SQUARE(x)       ((x) * (x))   /* x**2 */
#define M_PI		3.14159265358979323846

int32 GDXSDcomb[512*5];
char  GDXSDname[HDFE_NAMBUFSIZE];
char  GDXSDdims[HDFE_DIMBUFSIZE];


#define NGRID 200
/* Grid Structure External Arrays */
struct gridStructure 
{
    int32 active;
    int32 IDTable;
    int32 VIDTable[2];
    int32 fid;
    int32 nSDS;
    int32 *sdsID;
    int32 compcode;
    intn  compparm[5];
    int32 tilecode;
    int32 tilerank;
    int32 tiledims[8];
};
struct gridStructure GDXGrid[NGRID];



#define NGRIDREGN 256
struct gridRegion
{
    int32 fid;
    int32 gridID;
    int32 xStart;
    int32 xCount;
    int32 yStart;
    int32 yCount;
    int32 somStart;
    int32 somCount;
    float64 upleftpt[2];
    float64 lowrightpt[2];
    int32 StartVertical[8];
    int32 StopVertical[8];
    char *DimNamePtr[8];
};
struct gridRegion *GDXRegion[NGRIDREGN];

/* define a macro for the string size of the utility strings and some dimension
   list strings. The value of 80 in the previous version of this code 
   may not be enough in some cases. The length now is 512 which seems to 
   be more than enough to hold larger strings. */
   
#define UTLSTR_MAX_SIZE 512

/* Static projection table */
static struct {
    int32 projcode;
    char *projname;
} Projections[] = {
    {GCTP_GEO,	   "GCTP_GEO"},
    {GCTP_UTM,	   "GCTP_UTM"},
    {GCTP_SPCS,	   "GCTP_SPCS"},
    {GCTP_ALBERS,  "GCTP_ALBERS"},
    {GCTP_LAMCC,   "GCTP_LAMCC"},
    {GCTP_MERCAT,  "GCTP_MERCAT"},
    {GCTP_PS,	   "GCTP_PS"},
    {GCTP_POLYC,   "GCTP_POLYC"},
    {GCTP_EQUIDC,  "GCTP_EQUIDC"},
    {GCTP_TM,	   "GCTP_TM"},
    {GCTP_STEREO,  "GCTP_STEREO"},
    {GCTP_LAMAZ,   "GCTP_LAMAZ"},
    {GCTP_AZMEQD,  "GCTP_AZMEQD"},
    {GCTP_GNOMON,  "GCTP_GNOMON"},
    {GCTP_ORTHO,   "GCTP_ORTHO"},
    {GCTP_GVNSP,   "GCTP_GVNSP"},
    {GCTP_SNSOID,  "GCTP_SNSOID"},
    {GCTP_EQRECT,  "GCTP_EQRECT"},
    {GCTP_MILLER,  "GCTP_MILLER"},
    {GCTP_VGRINT,  "GCTP_VGRINT"},
    {GCTP_HOM,	   "GCTP_HOM"},
    {GCTP_ROBIN,   "GCTP_ROBIN"},
    {GCTP_SOM,	   "GCTP_SOM"},
    {GCTP_ALASKA,  "GCTP_ALASKA"},
    {GCTP_GOOD,	   "GCTP_GOOD"},
    {GCTP_MOLL,	   "GCTP_MOLL"},
    {GCTP_IMOLL,   "GCTP_IMOLL"},
    {GCTP_HAMMER,  "GCTP_HAMMER"},
    {GCTP_WAGIV,   "GCTP_WAGIV"},
    {GCTP_WAGVII,  "GCTP_WAGVII"},
    {GCTP_OBLEQA,  "GCTP_OBLEQA"},
    {GCTP_ISINUS1, "GCTP_ISINUS1"},
    {GCTP_CEA,	   "GCTP_CEA"},
    {GCTP_BCEA,	   "GCTP_BCEA"},
    {GCTP_ISINUS,  "GCTP_ISINUS"},
    {-1,	   NULL}
};


/* Grid Function Prototypes (internal routines) */
intn GDchkgdid(int32, char *, int32 *, int32 *, int32 *);
intn GDfldinfo(int32, char *, int32 *, int32 [], int32 *, char *);
intn GDdeffld(int32, char *, char *, int32, int32);
intn GDwrmeta(int32, char *, char *, int32);
intn GDSDfldsrch(int32, int32, char *, int32 *, int32 *, 
                 int32 *, int32 *, int32 [], int32 *);
intn GDwrrdfield(int32, char *, char *,
	    int32 [], int32 [], int32 [], VOIDP datbuf);
intn GDwrfld(int32, char *, int32 [], int32 [], int32 [], VOIDP);
intn GDrdfld(int32, char *, int32 [], int32 [], int32 [], VOIDP);
intn GDwrrdattr(int32, char *, int32, int32, char *, VOIDP);
intn GDll2ij(int32, int32, float64 [], int32, int32, int32, float64 [], 
             float64 [], int32, float64 [], float64 [],	int32 [], int32 [],
	     float64 [], float64 []);
intn GDij2ll(int32, int32, float64 [], int32, int32, int32, 
	     float64 [], float64 [], int32, int32 [], int32 [],
	     float64 [], float64 [], int32, int32);
intn GDreginfo(int32, int32, char *, int32 *, int32 *, int32 [], int32 *,
	       float64 [], float64 []);
intn  GDgetdefaults(int32, int32, float64 [], int32, float64 [], float64 []);
int32 GDdefvrtreg(int32, int32, char *, float64 []);
intn GDgetpix(int32, int32, float64 [], float64 [], int32 [], int32 []);
int32 GDgetpixval(int32, int32, int32 [], int32 [], char *, VOIDP);
intn GDtangentpnts(int32, float64 [], float64 [], float64 [], float64 [], 
		   float64 [], int32 *);
intn GDwrrdtile(int32, char *, char *, int32 [], VOIDP);
intn GDdeftle(int32, int32, int32, int32 []);
intn GDtleinfo(int32, char *, int32 *, int32 *, int32 []);
intn GDwrtle(int32, char *, int32 [], VOIDP);
intn GDrdtle(int32, char *, int32 [],  VOIDP);
intn GDll2mm_cea(int32 projcode,int32 zonecode, int32 spherecode,
		 float64 projparm[],
		 int32 xdimsize, int32 ydimsize,
		 float64 upleftpt[], float64 lowrightpt[], int32 npnts,
		 float64 lon[],float64 lat[],
		 float64 x[],float64 y[], float64 *scaleX,float64 *scaleY);

intn GDmm2ll_cea(int32 projcode,int32 zonecode, int32 spherecode,
		 float64 projparm[],
		 int32 xdimsize, int32 ydimsize,
		 float64 upleftpt[], float64 lowrightpt[], int32 npnts,
		 float64 x[], float64 y[], 
		 float64 lon[], float64 lat[]);
/*----------------------------------------------------------------------------|
|  BEGIN_PROLOG                                                               |
|                                                                             |
|  FUNCTION: GDopen                                                           |
|                                                                             |
|  DESCRIPTION: Opens or creates HDF file in order to create, read, or write  |
|                a grid.                                                      |
|                                                                             |
|                                                                             |
|  Return Value    Type     Units     Description                             |
|  ============   ======  =========   =====================================   |
|  fid            int32               HDF-EOS file ID                         |
|                                                                             |
|  INPUTS:                                                                    |
|  filename       char                Filename                                |
|  access         intn                HDF access code                         |
|                                                                             |
|                                                                             |
|  OUTPUTS:                                                                   |
|             None                                                            |
|                                                                             |
|  NOTES:                                                                     |
|                                                                             |
|                                                                             |
|   Date     Programmer   Description                                         |
|  ======   ============  =================================================   |
|  Jun 96   Joel Gales    Original Programmer                                 |
|                                                                             |
|  END_PROLOG                                                                 |
-----------------------------------------------------------------------------*/
int32
GDopen(char *filename, intn access)

{
    int32           fid /* HDF-EOS file ID */ ;

    /* Call EHopen to perform file access */
    /* ---------------------------------- */
    fid = EHopen(filename, access);

    return (fid);

}



/*----------------------------------------------------------------------------|
|  BEGIN_PROLOG                                                               |
|                                                                             |
|  FUNCTION: GDcreate                                                         |
|                                                                             |
|  DESCRIPTION: Creates a grid within the file.                               |
|                                                                             |
|                                                                             |
|  Return Value    Type     Units     Description                             |
|  ============   ======  =========   =====================================   |
|  gridID         int32               Grid structure ID                       |
|                                                                             |
|  INPUTS:                                                                    |
|  fid            int32               File ID                                 |
|  gridname       char                Grid structure name                     |
|  xdimsize       int32               Number of columns in grid               |
|  ydimsize       int32               Number of rows in grid                  |
|  upleftpt       float64             Location (m/deg) of upper left corner   |
|  lowrightpt     float64             Location (m/deg) of lower right corner  |
|                                                                             |
|                                                                             |
|  OUTPUTS:                                                                   |
|             None                                                            |
|                                                                             |
|  NOTES:                                                                     |
|                                                                             |
|                                                                             |
|   Date     Programmer   Description                                         |
|  ======   ============  =================================================   |
|  Jun 96   Joel Gales    Original Programmer                                 |
|  Aug 96   Joel Gales    Make metadata ODL compliant                         |
|  Aug 96   Joel Gales    Check grid name for ODL compliance                  |
|                                                                             |
|  END_PROLOG                                                                 |
-----------------------------------------------------------------------------*/
int32
GDcreate(int32 fid, char *gridname, int32 xdimsize, int32 ydimsize,
	 float64 upleftpt[], float64 lowrightpt[])
{
    intn            i;		/* Loop index */
    intn            ngridopen = 0;	/* # of grid structures open */
    intn            status = 0;	/* routine return status variable */

    uint8           access;	/* Read/Write file access code */

    int32           HDFfid;	/* HDF file id */
    int32           vgRef;	/* Vgroup reference number */
    int32           vgid[3];	/* Vgroup ID array */
    int32           gridID = -1;/* HDF-EOS grid ID */