1 /* zzpdpltc.f -- translated by f2c (version 19980913).
2 You must link the resulting object file with the libraries:
3 -lf2c -lm (in that order)
4 */
5
6 #include "f2c.h"
7
8 /* $Procedure ZZPDPLTC (Planetodetic coordinates, point latitude check) */
zzpdpltc_(doublereal * re,doublereal * f,doublereal * p,doublereal * lat)9 logical zzpdpltc_(doublereal *re, doublereal *f, doublereal *p, doublereal *
10 lat)
11 {
12 /* System generated locals */
13 logical ret_val;
14
15 /* Builtin functions */
16 double sqrt(doublereal);
17
18 /* Local variables */
19 doublereal xxpt, yxpt, a, b;
20 extern /* Subroutine */ int zzelnaxx_(doublereal *, doublereal *,
21 doublereal *, doublereal *, doublereal *);
22 doublereal r__;
23 extern /* Subroutine */ int chkin_(char *, ftnlen), errdp_(char *,
24 doublereal *, ftnlen);
25 doublereal r2;
26 extern logical failed_(void);
27 extern /* Subroutine */ int sigerr_(char *, ftnlen), chkout_(char *,
28 ftnlen), setmsg_(char *, ftnlen);
29 extern logical return_(void);
30
31 /* $ Abstract */
32
33 /* SPICE Private routine intended solely for the support of SPICE */
34 /* routines. Users should not call this routine directly due to the */
35 /* volatile nature of this routine. */
36
37 /* Indicate whether a given point on a planetodetic latitude cone */
38 /* has the correct latitude sign. */
39
40 /* $ Disclaimer */
41
42 /* THIS SOFTWARE AND ANY RELATED MATERIALS WERE CREATED BY THE */
43 /* CALIFORNIA INSTITUTE OF TECHNOLOGY (CALTECH) UNDER A U.S. */
44 /* GOVERNMENT CONTRACT WITH THE NATIONAL AERONAUTICS AND SPACE */
45 /* ADMINISTRATION (NASA). THE SOFTWARE IS TECHNOLOGY AND SOFTWARE */
46 /* PUBLICLY AVAILABLE UNDER U.S. EXPORT LAWS AND IS PROVIDED "AS-IS" */
47 /* TO THE RECIPIENT WITHOUT WARRANTY OF ANY KIND, INCLUDING ANY */
48 /* WARRANTIES OF PERFORMANCE OR MERCHANTABILITY OR FITNESS FOR A */
49 /* PARTICULAR USE OR PURPOSE (AS SET FORTH IN UNITED STATES UCC */
50 /* SECTIONS 2312-2313) OR FOR ANY PURPOSE WHATSOEVER, FOR THE */
51 /* SOFTWARE AND RELATED MATERIALS, HOWEVER USED. */
52
53 /* IN NO EVENT SHALL CALTECH, ITS JET PROPULSION LABORATORY, OR NASA */
54 /* BE LIABLE FOR ANY DAMAGES AND/OR COSTS, INCLUDING, BUT NOT */
55 /* LIMITED TO, INCIDENTAL OR CONSEQUENTIAL DAMAGES OF ANY KIND, */
56 /* INCLUDING ECONOMIC DAMAGE OR INJURY TO PROPERTY AND LOST PROFITS, */
57 /* REGARDLESS OF WHETHER CALTECH, JPL, OR NASA BE ADVISED, HAVE */
58 /* REASON TO KNOW, OR, IN FACT, SHALL KNOW OF THE POSSIBILITY. */
59
60 /* RECIPIENT BEARS ALL RISK RELATING TO QUALITY AND PERFORMANCE OF */
61 /* THE SOFTWARE AND ANY RELATED MATERIALS, AND AGREES TO INDEMNIFY */
62 /* CALTECH AND NASA FOR ALL THIRD-PARTY CLAIMS RESULTING FROM THE */
63 /* ACTIONS OF RECIPIENT IN THE USE OF THE SOFTWARE. */
64
65 /* $ Required_Reading */
66
67 /* None. */
68
69 /* $ Keywords */
70
71 /* COORDINATES */
72 /* GEOMETRY */
73 /* PLANETODETIC */
74 /* LATITUDE */
75 /* MATH */
76
77 /* $ Declarations */
78 /* $ Brief_I/O */
79
80 /* Variable I/O Description */
81 /* -------- --- -------------------------------------------------- */
82 /* RE I Equatorial radius. */
83 /* F I Flattening coefficient. */
84 /* P I Three-dimensional point. */
85 /* LAT I Planetodetic latitude. */
86
87 /* The function returns .TRUE. if the sign of the planetodetic */
88 /* latitude of the input point matches that of the input */
89 /* planetodetic latitude. */
90
91 /* $ Detailed_Input */
92
93 /* RE, */
94 /* F are, respectively, the equatorial radius */
95 /* and flattening coefficient of a biaxial */
96 /* spheroid. */
97
98 /* The polar radius RP of the spheroid is */
99
100 /* RP = RE * ( 1 - F ) */
101
102 /* RP may be less than, equal to, or greater than RE. */
103
104
105 /* P is a point (equivalently, a vector) in */
106 /* three-dimensional space. P is expressed in Cartesian */
107 /* coordinates. P must lie on the planetodetic */
108 /* latitude cone corresponding to LAT. */
109
110 /* The units of P must be consistent with those of RE. */
111
112
113 /* LAT is a planetodetic latitude value defining a cone. */
114 /* Units are radius. */
115
116
117 /* $ Detailed_Output */
118
119 /* The function returns .TRUE. if any of the following */
120 /* are true: */
121
122 /* - The input spheroid is prolate or spherical. */
123
124 /* - The input latitude is zero. */
125
126 /* - The input point is determined to have planetodetic */
127 /* latitude of the same sign as the input latitude. */
128
129 /* Otherwise the function returns .FALSE. */
130
131 /* See Particulars below for details. */
132
133 /* $ Parameters */
134
135 /* None. */
136
137 /* $ Exceptions */
138
139 /* 1) If either the equatorial radius or flattening coefficient is */
140 /* invalid, the error SPICE(VALUEOUTOFRANGE) will be signaled. */
141
142 /* $ Files */
143
144 /* None. */
145
146 /* $ Particulars */
147
148 /* This routine supports computation of the intersection between */
149 /* a ray and a planetodetic volume element. The routine is used */
150 /* to determine whether the intersection of the ray and latitude */
151 /* cone is actually on the element's boundary. */
152
153 /* This function serves as a "macro" that executes a logical test */
154 /* that would be awkward to perform in-line. */
155
156 /* For a given reference spheroid, all points having a given */
157 /* planetodetic latitude lie on a cone. The axis of the cone */
158 /* coincides with the Z axis of the coordinate system. */
159
160 /* The possibility that a point on a given latitude can have a Z */
161 /* coordinate with the wrong sign exists for oblate reference */
162 /* spheroids. For these shapes, the vertex of a positive latitude */
163 /* cone has a negative Z component, and the vertex of a negative */
164 /* latitude cone has a positive Z component. */
165
166 /* The purpose of the function is to determine, for a point that */
167 /* lies on a given planetodetic latitude cone, whether that point is */
168 /* on the correct side of the X-Y plane: that is, the side */
169 /* corresponding to the sign of the input latitude value. */
170
171 /* This check is not as simple as checking the sign of the Z */
172 /* component of the input point. For values of LAT that are non-zero */
173 /* but have very small magnitude, points that are nominally on the */
174 /* corresponding latitude cone may have Z components of the wrong */
175 /* sign due to round-off errors that occurred in the process of */
176 /* computing those points. For such cases, the input point is */
177 /* checked by comparing its distance from the Z axis to the X */
178 /* intercept of a line normal to the reference spheroid at a point */
179 /* having planetodetic latitude LAT and longitude 0. This check can */
180 /* be performed accurately even when the Z component of the input */
181 /* point consists of noise. */
182
183 /* This routine does not test the input point to determine whether */
184 /* it lies on the latitude cone defined by the input argument LAT. */
185 /* The caller must ensure that this is the case. */
186
187 /* $ Examples */
188
189 /* None. */
190
191 /* $ Restrictions */
192
193 /* This routine does not test the input point to determine whether */
194 /* it lies on the latitude cone defined by the input argument LAT. */
195 /* The caller must ensure that this is the case. */
196
197 /* $ Literature_References */
198
199 /* None. */
200
201 /* $ Author_and_Institution */
202
203 /* N.J. Bachman (JPL) */
204
205 /* $ Version */
206
207 /* - SPICELIB Version 1.0.0, 13-JAN-2017 (NJB) */
208
209 /* -& */
210 /* $ Index_Entries */
211
212 /* is point on planetodetic latitude boundary */
213
214 /* -& */
215
216 /* SPICELIB functions */
217
218
219 /* Local parameters */
220
221
222 /* Local variables */
223
224
225 /* Give the function a default value. */
226
227 ret_val = FALSE_;
228 if (return_()) {
229 return ret_val;
230 }
231 chkin_("ZZPDPLTC", (ftnlen)8);
232
233 /* The equatorial radius must be greater than zero. */
234
235 if (*re <= 0.) {
236 setmsg_("Equatorial radius was *.", (ftnlen)24);
237 errdp_("*", re, (ftnlen)1);
238 sigerr_("SPICE(VALUEOUTOFRANGE)", (ftnlen)22);
239 chkout_("ZZPDPLTC", (ftnlen)8);
240 return ret_val;
241 }
242
243 /* If the flattening coefficient is greater than one, the polar */
244 /* radius computed below is negative. If it's equal to one, the */
245 /* polar radius is zero. Either case is a problem, so signal an */
246 /* error and check out. */
247
248 if (*f >= 1.) {
249 setmsg_("Flattening coefficient was *.", (ftnlen)29);
250 errdp_("*", f, (ftnlen)1);
251 sigerr_("SPICE(VALUEOUTOFRANGE)", (ftnlen)22);
252 chkout_("ZZPDPLTC", (ftnlen)8);
253 return ret_val;
254 }
255
256 /* The input point is assumed to be on the cone */
257 /* corresponding to the input latitude. */
258
259 /* If the reference spheroid is prolate or spherical, there's */
260 /* nothing to do: the point is automatically on the correct side of */
261 /* the X-Y plane. */
262
263 if (*f <= 0.) {
264 ret_val = TRUE_;
265 } else {
266
267 /* This is the oblate case. */
268
269
270 /* If the point is on the "correct" side of the X-Y plane--- */
271 /* that is, its Z component as the same sign as LAT, the */
272 /* point is considered to have the correct latitude. */
273
274 /* If the point is on the X-Y plane, or if LAT is zero, the point */
275 /* is considered to have the indicated latitude. We condense */
276 /* these cases by requiring that */
277
278 /* LAT * P(3) >= 0 */
279
280 /* rather than */
281
282 /* LAT * P(3) > 0 */
283
284
285 if (p[2] * *lat >= 0.) {
286 ret_val = TRUE_;
287 } else if (abs(*lat) >= .01) {
288
289 /* Ideally, the input point is considered to have the given */
290 /* latitude if the point is on the side of the X-Y plane */
291 /* corresponding to the sign of the input latitude. The */
292 /* problem with this criterion is that it can't be applied */
293 /* correctly when LAT has very small magnitude. */
294
295 /* If the magnitude of LAT is above the limit, it's ok to */
296 /* use the sign of P(3) to determine whether the point */
297 /* has the given latitude. */
298
299 /* The point has the indicated latitude if both LAT and P(3) */
300 /* have the same sign. In this case, we know they have the */
301 /* opposite sign. */
302
303 ret_val = FALSE_;
304 } else {
305
306 /* At this point we know LAT is non-zero, so the cone */
307 /* corresponding to LAT has its vertex on the opposite side of */
308 /* the X-Y plane from any point having latitude LAT. So it's */
309 /* possible for a point to be on the cone but not have the */
310 /* correct latitude. */
311
312 /* We're in the special case where the point's Z component */
313 /* has the opposite sign as LAT, and the magnitude of LAT */
314 /* is below the limit. We don't automatically reject the */
315 /* point in this case: we'll accept it if it is far enough */
316 /* from the Z axis to be outside the portion of the latitude */
317 /* cone on the wrong side of the X-Y plane. */
318
319 a = *re;
320 b = a * (1. - *f);
321
322 /* Compute the intercepts of a normal vector of a point */
323 /* at latitude LAT, longitude 0, with the X and Y axes. */
324
325 zzelnaxx_(&a, &b, lat, &xxpt, &yxpt);
326 if (failed_()) {
327 chkout_("ZZPDPLTC", (ftnlen)8);
328 return ret_val;
329 }
330
331 /* We check the point's distance from the Z axis. This can be */
332 /* done accurately even when the Z component of P consists of */
333 /* noise. */
334
335 /* The point is considered to have the correct latitude when */
336 /* it is farther from the Z axis than the intercept on the X */
337 /* axis of a normal line passing through a point having */
338 /* latitude LAT and longitude 0. Ideally, a point that is on */
339 /* the latitude cone and that satisfies this criterion must be */
340 /* on the correct side of the X-Y plane. */
341
342 r2 = p[0] * p[0] + p[1] * p[1];
343 r__ = sqrt((max(r2,0.)));
344 ret_val = r__ >= xxpt;
345 }
346 }
347 chkout_("ZZPDPLTC", (ftnlen)8);
348 return ret_val;
349 } /* zzpdpltc_ */
350
351