1 : /******************************************************************************
2 : * $Id: gctp_wrap.c 22751 2011-07-18 15:25:05Z winkey $
3 : *
4 : * Project: Hierarchical Data Format Release 4 (HDF4)
5 : * Purpose: This is the wrapper code to use OGR Coordinate Transformation
6 : * services instead of GCTP library
7 : * Author: Andrey Kiselev, dron@ak4719.spb.edu
8 : *
9 : ******************************************************************************
10 : * Copyright (c) 2004, Andrey Kiselev <dron@ak4719.spb.edu>
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 : #include "ogr_srs_api.h"
32 : #include <stdlib.h>
33 :
34 : #include "mfhdf.h"
35 :
36 : #include <math.h>
37 :
38 : #ifndef PI
39 : #ifndef M_PI
40 : #define PI (3.141592653589793238)
41 : #else
42 : #define PI (M_PI)
43 : #endif
44 : #endif
45 :
46 : #define DEG (180.0 / PI)
47 : #define RAD (PI / 180.0)
48 :
49 : /***** static vars to store the transformers in *****/
50 : /***** this is not thread safe *****/
51 :
52 : static OGRCoordinateTransformationH hForCT, hInvCT;
53 :
54 : /******************************************************************************
55 : function for forward gctp transformation
56 :
57 : gctp expects Longitude and Latitude values to be in radians
58 : ******************************************************************************/
59 :
60 0 : int32 osr_for(
61 : double lon, /* (I) Longitude */
62 : double lat, /* (I) Latitude */
63 : double *x, /* (O) X projection coordinate */
64 : double *y) /* (O) Y projection coordinate */
65 : {
66 :
67 0 : double dfX, dfY, dfZ = 0.0;
68 :
69 0 : dfX = lon * DEG;
70 0 : dfY = lat * DEG;
71 :
72 0 : OCTTransform( hForCT, 1, &dfX, &dfY, &dfZ);
73 :
74 0 : *x = dfX;
75 0 : *y = dfY;
76 :
77 0 : return 0;
78 : }
79 :
80 : /******************************************************************************
81 : function to init a forward gctp transformer
82 : ******************************************************************************/
83 :
84 0 : void for_init(
85 : int32 outsys, /* output system code */
86 : int32 outzone, /* output zone number */
87 : float64 *outparm, /* output array of projection parameters */
88 : int32 outdatum, /* output datum */
89 : char *fn27, /* NAD 1927 parameter file */
90 : char *fn83, /* NAD 1983 parameter file */
91 : int32 *iflg, /* status flag */
92 : int32 (*for_trans[])(double, double, double *, double *))
93 : /* forward function pointer */
94 : {
95 0 : OGRSpatialReferenceH hOutSourceSRS, hLatLong = NULL;
96 :
97 0 : *iflg = 0;
98 :
99 0 : hOutSourceSRS = OSRNewSpatialReference( NULL );
100 0 : OSRImportFromUSGS( hOutSourceSRS, outsys, outzone, outparm, outdatum );
101 0 : hLatLong = OSRNewSpatialReference ( SRS_WKT_WGS84 );
102 :
103 0 : hForCT = OCTNewCoordinateTransformation( hLatLong, hOutSourceSRS );
104 :
105 0 : OSRDestroySpatialReference( hOutSourceSRS );
106 0 : OSRDestroySpatialReference( hLatLong );
107 :
108 0 : for_trans[outsys] = osr_for;
109 0 : }
110 :
111 : /******************************************************************************
112 : function for inverse gctp transformation
113 :
114 : gctp returns Longitude and Latitude values in radians
115 : ******************************************************************************/
116 :
117 0 : int32 osr_inv(
118 : double x, /* (I) X projection coordinate */
119 : double y, /* (I) Y projection coordinate */
120 : double *lon, /* (O) Longitude */
121 : double *lat) /* (O) Latitude */
122 : {
123 :
124 0 : double dfX, dfY, dfZ = 0.0;
125 :
126 0 : dfX = x;
127 0 : dfY = y;
128 :
129 0 : OCTTransform( hInvCT, 1, &dfX, &dfY, &dfZ );
130 :
131 0 : *lon = dfX * RAD;
132 0 : *lat = dfY * RAD;
133 :
134 0 : return 0;
135 : }
136 :
137 : /******************************************************************************
138 : function to init a inverse gctp transformer
139 : ******************************************************************************/
140 :
141 0 : void inv_init(
142 : int32 insys, /* input system code */
143 : int32 inzone, /* input zone number */
144 : float64 *inparm, /* input array of projection parameters */
145 : int32 indatum, /* input datum code */
146 : char *fn27, /* NAD 1927 parameter file */
147 : char *fn83, /* NAD 1983 parameter file */
148 : int32 *iflg, /* status flag */
149 : int32 (*inv_trans[])(double, double, double*, double*))
150 : /* inverse function pointer */
151 : {
152 :
153 0 : OGRSpatialReferenceH hInSourceSRS, hLatLong = NULL;
154 0 : *iflg = 0;
155 :
156 0 : hInSourceSRS = OSRNewSpatialReference( NULL );
157 0 : OSRImportFromUSGS( hInSourceSRS, insys, inzone, inparm, indatum );
158 :
159 0 : hLatLong = OSRNewSpatialReference ( SRS_WKT_WGS84 );
160 :
161 0 : hInvCT = OCTNewCoordinateTransformation( hInSourceSRS, hLatLong );
162 :
163 0 : OSRDestroySpatialReference( hInSourceSRS );
164 0 : OSRDestroySpatialReference( hLatLong );
165 :
166 0 : inv_trans[insys] = osr_inv;
167 0 : }
168 :
169 : /******************************************************************************
170 : function to cleanup the transformers
171 :
172 : note: gctp does not have a function that does this
173 : ******************************************************************************/
174 :
175 0 : void gctp_destroy(void) {
176 0 : OCTDestroyCoordinateTransformation ( hForCT );
177 0 : OCTDestroyCoordinateTransformation ( hInvCT );
178 0 : }
|
