~ubuntu-branches/ubuntu/raring/geotranz/raring

« back to all changes in this revision

Viewing changes to dt_cc/ups/ups.c

  • Committer: Bazaar Package Importer
  • Author(s): Roberto Lumbreras
  • Date: 2008-10-17 14:43:09 UTC
  • Revision ID: james.westby@ubuntu.com-20081017144309-jb7uzfi1y1lvez8j
Tags: upstream-2.4.2
ImportĀ upstreamĀ versionĀ 2.4.2

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
/********************************************************************/
 
2
/* RSC IDENTIFIER: UPS
 
3
 *
 
4
 *
 
5
 * ABSTRACT
 
6
 *
 
7
 *    This component provides conversions between geodetic (latitude
 
8
 *    and longitude) coordinates and Universal Polar Stereographic (UPS)
 
9
 *    projection (hemisphere, easting, and northing) coordinates.
 
10
 *
 
11
 *
 
12
 * ERROR HANDLING
 
13
 *
 
14
 *    This component checks parameters for valid values.  If an 
 
15
 *    invalid value is found the error code is combined with the 
 
16
 *    current error code using the bitwise or.  This combining allows  
 
17
 *    multiple error codes to be returned. The possible error codes 
 
18
 *    are:
 
19
 *
 
20
 *         UPS_NO_ERROR           : No errors occurred in function
 
21
 *         UPS_LAT_ERROR          : Latitude outside of valid range
 
22
 *                                   (North Pole: 83.5 to 90,
 
23
 *                                    South Pole: -79.5 to -90)
 
24
 *         UPS_LON_ERROR          : Longitude outside of valid range
 
25
 *                                   (-180 to 360 degrees)
 
26
 *         UPS_HEMISPHERE_ERROR   : Invalid hemisphere ('N' or 'S')
 
27
 *         UPS_EASTING_ERROR      : Easting outside of valid range,
 
28
 *                                   (0 to 4,000,000m)
 
29
 *         UPS_NORTHING_ERROR     : Northing outside of valid range,
 
30
 *                                   (0 to 4,000,000m)
 
31
 *         UPS_A_ERROR            : Semi-major axis less than or equal to zero
 
32
 *         UPS_INV_F_ERROR        : Inverse flattening outside of valid range
 
33
 *                                                                                     (250 to 350)
 
34
 *
 
35
 *
 
36
 * REUSE NOTES
 
37
 *
 
38
 *    UPS is intended for reuse by any application that performs a Universal
 
39
 *    Polar Stereographic (UPS) projection.
 
40
 *
 
41
 *
 
42
 * REFERENCES
 
43
 *
 
44
 *    Further information on UPS can be found in the Reuse Manual.
 
45
 *
 
46
 *    UPS originated from :  U.S. Army Topographic Engineering Center
 
47
 *                           Geospatial Information Division
 
48
 *                           7701 Telegraph Road
 
49
 *                           Alexandria, VA  22310-3864
 
50
 *
 
51
 *
 
52
 * LICENSES
 
53
 *
 
54
 *    None apply to this component.
 
55
 *
 
56
 *
 
57
 * RESTRICTIONS
 
58
 *
 
59
 *    UPS has no restrictions.
 
60
 *
 
61
 *
 
62
 * ENVIRONMENT
 
63
 *
 
64
 *    UPS was tested and certified in the following environments:
 
65
 *
 
66
 *    1. Solaris 2.5 with GCC version 2.8.1
 
67
 *    2. Windows 95 with MS Visual C++ version 6
 
68
 *
 
69
 *
 
70
 * MODIFICATIONS
 
71
 *
 
72
 *    Date              Description
 
73
 *    ----              -----------
 
74
 *    06-11-95          Original Code
 
75
 *    03-01-97          Original Code
 
76
 *
 
77
 *
 
78
 */
 
79
 
 
80
 
 
81
/************************************************************************/
 
82
/*
 
83
 *                               INCLUDES
 
84
 */
 
85
 
 
86
#include <math.h>
 
87
#include "polarst.h"
 
88
#include "ups.h"
 
89
/*
 
90
 *    math.h     - Is needed to call the math functions.
 
91
 *    polar.h    - Is used to convert polar stereographic coordinates
 
92
 *    ups.h      - Defines the function prototypes for the ups module.
 
93
 */
 
94
 
 
95
 
 
96
/************************************************************************/
 
97
/*                               GLOBAL DECLARATIONS
 
98
 *
 
99
 */
 
100
 
 
101
#define PI       3.14159265358979323e0  /* PI     */
 
102
#define PI_OVER    (PI/2.0e0)           /* PI over 2 */
 
103
#define MAX_LAT    ((PI * 90)/180.0)    /* 90 degrees in radians */
 
104
#define MAX_ORIGIN_LAT ((81.114528 * PI) / 180.0)
 
105
#define MIN_NORTH_LAT (83.5*PI/180.0)
 
106
#define MIN_SOUTH_LAT (-79.5*PI/180.0)
 
107
#define MIN_EAST_NORTH 0
 
108
#define MAX_EAST_NORTH 4000000
 
109
 
 
110
/* Ellipsoid Parameters, default to WGS 84  */
 
111
static double UPS_a = 6378137.0;          /* Semi-major axis of ellipsoid in meters   */
 
112
static double UPS_f = 1 / 298.257223563;  /* Flattening of ellipsoid  */
 
113
const double UPS_False_Easting = 2000000;
 
114
const double UPS_False_Northing = 2000000;
 
115
static double UPS_Origin_Latitude = MAX_ORIGIN_LAT;  /*set default = North Hemisphere */
 
116
static double UPS_Origin_Longitude = 0.0;
 
117
 
 
118
 
 
119
/************************************************************************/
 
120
/*                              FUNCTIONS
 
121
 *
 
122
 */
 
123
 
 
124
 
 
125
long Set_UPS_Parameters( double a,
 
126
                         double f)
 
127
{
 
128
/*
 
129
 * The function SET_UPS_PARAMETERS receives the ellipsoid parameters and sets
 
130
 * the corresponding state variables. If any errors occur, the error code(s)
 
131
 * are returned by the function, otherwise UPS_NO_ERROR is returned.
 
132
 *
 
133
 *   a     : Semi-major axis of ellipsoid in meters (input)
 
134
 *   f     : Flattening of ellipsoid                                          (input)
 
135
 */
 
136
 
 
137
  double inv_f = 1 / f;
 
138
  long Error_Code = UPS_NO_ERROR;
 
139
 
 
140
  if (a <= 0.0)
 
141
  { /* Semi-major axis must be greater than zero */
 
142
    Error_Code |= UPS_A_ERROR;
 
143
  }
 
144
  if ((inv_f < 250) || (inv_f > 350))
 
145
  { /* Inverse flattening must be between 250 and 350 */
 
146
    Error_Code |= UPS_INV_F_ERROR;
 
147
  }
 
148
 
 
149
  if (!Error_Code)
 
150
  { /* no errors */
 
151
    UPS_a = a;
 
152
    UPS_f = f;
 
153
  }
 
154
  return (Error_Code);
 
155
}  /* END of Set_UPS_Parameters  */
 
156
 
 
157
 
 
158
void Get_UPS_Parameters( double *a,
 
159
                         double *f)
 
160
{
 
161
/*
 
162
 * The function Get_UPS_Parameters returns the current ellipsoid parameters.
 
163
 *
 
164
 *  a      : Semi-major axis of ellipsoid, in meters (output)
 
165
 *  f      : Flattening of ellipsoid                                           (output)
 
166
 */
 
167
 
 
168
  *a = UPS_a;
 
169
  *f = UPS_f;
 
170
  return;
 
171
} /* END OF Get_UPS_Parameters */
 
172
 
 
173
 
 
174
long Convert_Geodetic_To_UPS ( double Latitude,
 
175
                               double Longitude,
 
176
                               char   *Hemisphere,
 
177
                               double *Easting,
 
178
                               double *Northing)
 
179
{
 
180
/*
 
181
 *  The function Convert_Geodetic_To_UPS converts geodetic (latitude and
 
182
 *  longitude) coordinates to UPS (hemisphere, easting, and northing)
 
183
 *  coordinates, according to the current ellipsoid parameters. If any 
 
184
 *  errors occur, the error code(s) are returned by the function, 
 
185
 *  otherwide UPS_NO_ERROR is returned.
 
186
 *
 
187
 *    Latitude      : Latitude in radians                       (input)
 
188
 *    Longitude     : Longitude in radians                      (input)
 
189
 *    Hemisphere    : Hemisphere either 'N' or 'S'              (output)
 
190
 *    Easting       : Easting/X in meters                       (output)
 
191
 *    Northing      : Northing/Y in meters                      (output)
 
192
 */
 
193
 
 
194
  double tempEasting, tempNorthing;
 
195
  long Error_Code = UPS_NO_ERROR;
 
196
 
 
197
  if ((Latitude < -MAX_LAT) || (Latitude > MAX_LAT))
 
198
  {   /* latitude out of range */
 
199
    Error_Code |= UPS_LAT_ERROR;
 
200
  }
 
201
  if ((Latitude < 0) && (Latitude > MIN_SOUTH_LAT))
 
202
    Error_Code |= UPS_LAT_ERROR;
 
203
  if ((Latitude >= 0) && (Latitude < MIN_NORTH_LAT))
 
204
    Error_Code |= UPS_LAT_ERROR;
 
205
  if ((Longitude < -PI) || (Longitude > (2 * PI)))
 
206
  {  /* slam out of range */
 
207
    Error_Code |= UPS_LON_ERROR;
 
208
  }
 
209
 
 
210
  if (!Error_Code)
 
211
  {  /* no errors */
 
212
    if (Latitude < 0)
 
213
    {
 
214
      UPS_Origin_Latitude = -MAX_ORIGIN_LAT; 
 
215
      *Hemisphere = 'S';
 
216
    }
 
217
    else
 
218
    {
 
219
      UPS_Origin_Latitude = MAX_ORIGIN_LAT; 
 
220
      *Hemisphere = 'N';
 
221
    }
 
222
 
 
223
 
 
224
    Set_Polar_Stereographic_Parameters( UPS_a,
 
225
                                        UPS_f,
 
226
                                        UPS_Origin_Latitude,
 
227
                                        UPS_Origin_Longitude,
 
228
                                        UPS_False_Easting,
 
229
                                        UPS_False_Northing);
 
230
 
 
231
    Convert_Geodetic_To_Polar_Stereographic(Latitude,
 
232
                                            Longitude,
 
233
                                            &tempEasting,
 
234
                                            &tempNorthing);
 
235
 
 
236
    *Easting = tempEasting;
 
237
    *Northing = tempNorthing;
 
238
  }  /*  END of if(!Error_Code)   */
 
239
 
 
240
  return Error_Code;
 
241
}  /* END OF Convert_Geodetic_To_UPS  */
 
242
 
 
243
 
 
244
long Convert_UPS_To_Geodetic(char   Hemisphere,
 
245
                             double Easting,
 
246
                             double Northing,
 
247
                             double *Latitude,
 
248
                             double *Longitude)
 
249
{
 
250
/*
 
251
 *  The function Convert_UPS_To_Geodetic converts UPS (hemisphere, easting, 
 
252
 *  and northing) coordinates to geodetic (latitude and longitude) coordinates
 
253
 *  according to the current ellipsoid parameters.  If any errors occur, the 
 
254
 *  error code(s) are returned by the function, otherwise UPS_NO_ERROR is 
 
255
 *  returned.
 
256
 *
 
257
 *    Hemisphere    : Hemisphere either 'N' or 'S'              (input)
 
258
 *    Easting       : Easting/X in meters                       (input)
 
259
 *    Northing      : Northing/Y in meters                      (input)
 
260
 *    Latitude      : Latitude in radians                       (output)
 
261
 *    Longitude     : Longitude in radians                      (output)
 
262
 */
 
263
 
 
264
  long Error_Code = UPS_NO_ERROR;
 
265
 
 
266
  if ((Hemisphere != 'N') && (Hemisphere != 'S'))
 
267
    Error_Code |= UPS_HEMISPHERE_ERROR;
 
268
  if ((Easting < MIN_EAST_NORTH) || (Easting > MAX_EAST_NORTH))
 
269
    Error_Code |= UPS_EASTING_ERROR;
 
270
  if ((Northing < MIN_EAST_NORTH) || (Northing > MAX_EAST_NORTH))
 
271
    Error_Code |= UPS_NORTHING_ERROR;
 
272
 
 
273
  if (Hemisphere =='N')
 
274
  {UPS_Origin_Latitude = MAX_ORIGIN_LAT;}
 
275
  if (Hemisphere =='S')
 
276
  {UPS_Origin_Latitude = -MAX_ORIGIN_LAT;}
 
277
 
 
278
  if (!Error_Code)
 
279
  {   /*  no errors   */
 
280
    Set_Polar_Stereographic_Parameters( UPS_a,
 
281
                                        UPS_f,
 
282
                                        UPS_Origin_Latitude,
 
283
                                        UPS_Origin_Longitude,
 
284
                                        UPS_False_Easting,
 
285
                                        UPS_False_Northing);
 
286
 
 
287
 
 
288
 
 
289
    Convert_Polar_Stereographic_To_Geodetic( Easting,
 
290
                                             Northing,
 
291
                                             Latitude,
 
292
                                             Longitude); 
 
293
 
 
294
 
 
295
    if ((*Latitude < 0) && (*Latitude > MIN_SOUTH_LAT))
 
296
      Error_Code |= UPS_LAT_ERROR;
 
297
    if ((*Latitude >= 0) && (*Latitude < MIN_NORTH_LAT))
 
298
      Error_Code |= UPS_LAT_ERROR;
 
299
  }  /*  END OF if(!Error_Code) */
 
300
  return (Error_Code);
 
301
}  /*  END OF Convert_UPS_To_Geodetic  */ 
 
302