rsspice/generated/spicelib/zzinpdt.rs
1//
2// GENERATED FILE
3//
4
5use super::*;
6use f2rust_std::*;
7
8const XFRACT: f64 = 0.0000000001;
9const KEYXFR: i32 = 1;
10const SGREED: f64 = 0.00000001;
11const KEYSGR: i32 = (KEYXFR + 1);
12const SGPADM: f64 = 0.0000000001;
13const KEYSPM: i32 = (KEYSGR + 1);
14const PTMEMM: f64 = 0.0000001;
15const KEYPTM: i32 = (KEYSPM + 1);
16const ANGMRG: f64 = 0.000000000001;
17const KEYAMG: i32 = (KEYPTM + 1);
18const LONALI: f64 = 0.000000000001;
19const KEYLAL: i32 = (KEYAMG + 1);
20pub const NONE: i32 = 0;
21pub const LONIDX: i32 = 1;
22pub const LATIDX: i32 = 2;
23pub const ALTIDX: i32 = 3;
24pub const LATMRG: f64 = 0.00000001;
25const WEST: i32 = 1;
26const EAST: i32 = 2;
27const SOUTH: i32 = 1;
28const NORTH: i32 = 2;
29const LOWER: i32 = 1;
30const UPPER: i32 = 2;
31const LT: i32 = -1;
32const GT: i32 = 1;
33
34struct SaveVars {
35    ALTBDS: StackArray2D<f64, 6>,
36    HPI: f64,
37    PI2: f64,
38    FIRST: bool,
39}
40
41impl SaveInit for SaveVars {
42    fn new() -> Self {
43        let mut ALTBDS = StackArray2D::<f64, 6>::new(1..=2, 1..=3);
44        let mut HPI: f64 = 0.0;
45        let mut PI2: f64 = 0.0;
46        let mut FIRST: bool = false;
47
48        FIRST = true;
49        HPI = -1.0;
50
51        Self {
52            ALTBDS,
53            HPI,
54            PI2,
55            FIRST,
56        }
57    }
58}
59
60//$Procedure ZZINPDT ( DSK, in planetodetic element? )
61pub fn ZZINPDT(
62    P: &[f64],
63    BOUNDS: &[f64],
64    CORPAR: &[f64],
65    MARGIN: f64,
66    EXCLUD: i32,
67    INSIDE: &mut bool,
68    ctx: &mut Context,
69) -> f2rust_std::Result<()> {
70    let save = ctx.get_vars::<SaveVars>();
71    let save = &mut *save.borrow_mut();
72
73    let P = DummyArray::new(P, 1..=3);
74    let BOUNDS = DummyArray2D::new(BOUNDS, 1..=2, 1..=3);
75    let CORPAR = DummyArray::new(CORPAR, 1..);
76    let mut AMNALT: f64 = 0.0;
77    let mut AMNLAT: f64 = 0.0;
78    let mut AMNLON: f64 = 0.0;
79    let mut AMXALT: f64 = 0.0;
80    let mut AMXLAT: f64 = 0.0;
81    let mut AMXLON: f64 = 0.0;
82    let mut DLON: f64 = 0.0;
83    let mut F: f64 = 0.0;
84    let mut LON: f64 = 0.0;
85    let mut LONMRG: f64 = 0.0;
86    let mut MAXALT: f64 = 0.0;
87    let mut MAXLAT: f64 = 0.0;
88    let mut MAXLON: f64 = 0.0;
89    let mut MINALT: f64 = 0.0;
90    let mut MINLAT: f64 = 0.0;
91    let mut MINLON: f64 = 0.0;
92    let mut PCNLAT: f64 = 0.0;
93    let mut R: f64 = 0.0;
94    let mut RE: f64 = 0.0;
95    let mut RELMIN: i32 = 0;
96    let mut RELMAX: i32 = 0;
97    let mut ALPASS: bool = false;
98
99    //
100    // SPICELIB functions
101    //
102
103    //
104    // Local parameters
105    //
106
107    //
108    // Element boundary indices:
109    //
110
111    //
112    // Numeric relation codes returned by ZZPDCMPL:
113    //
114
115    //
116    // The code EQ can be returned by ZZPDCMPL, but we make no
117    // references to it, so it's not declared here. For the
118    // record, EQ is set to 0.
119    //
120
121    //
122    // Local variables
123    //
124
125    //
126    // Saved variables
127    //
128
129    //
130    // Initial values
131    //
132
133    if RETURN(ctx) {
134        return Ok(());
135    }
136
137    CHKIN(b"ZZINPDT", ctx)?;
138
139    if save.FIRST {
140        save.HPI = HALFPI(ctx);
141        save.PI2 = TWOPI(ctx);
142        //
143        // Initialize the local array used for altitude checks.
144        //
145        save.ALTBDS[[SOUTH, LATIDX]] = -HALFPI(ctx);
146        save.ALTBDS[[NORTH, LATIDX]] = HALFPI(ctx);
147        save.ALTBDS[[WEST, LONIDX]] = -PI(ctx);
148        save.ALTBDS[[EAST, LONIDX]] = PI(ctx);
149        save.ALTBDS[[LOWER, ALTIDX]] = 0.0;
150        save.ALTBDS[[UPPER, ALTIDX]] = 0.0;
151
152        save.FIRST = false;
153    }
154
155    if ((EXCLUD < 0) || (EXCLUD > 3)) {
156        SETMSG(b"EXCLUD must be in the range 0:3 but was #.", ctx);
157        ERRINT(b"#", EXCLUD, ctx);
158        SIGERR(b"SPICE(VALUEOUTOFRANGE)", ctx)?;
159        CHKOUT(b"ZZINPDT", ctx)?;
160        return Ok(());
161    }
162
163    //
164    // Get the planetocentric [sic] coordinates of the input point. The
165    // latitude we obtain will be planetocentric. To emphasize this, we
166    // use the name "PCNLAT."
167    //
168    RECLAT(P.as_slice(), &mut R, &mut LON, &mut PCNLAT);
169    //
170    // RECLAT is error free, so we don't call FAILED() here.
171    //
172
173    if (MARGIN == 0.0) {
174        //
175        // ZZINPDT0 contains the logic required to determine whether
176        // the input point is contained in the element.
177        //
178        ZZINPDT0(
179            P.as_slice(),
180            LON,
181            BOUNDS.as_slice(),
182            CORPAR.as_slice(),
183            EXCLUD,
184            INSIDE,
185            ctx,
186        )?;
187
188        CHKOUT(b"ZZINPDT", ctx)?;
189        return Ok(());
190    } else if (MARGIN < 0.0) {
191        SETMSG(b"Margin must be non-negative but was #.", ctx);
192        ERRDP(b"#", MARGIN, ctx);
193        SIGERR(b"SPICE(VALUEOUTOFRANGE)", ctx)?;
194        CHKOUT(b"ZZINPDT", ctx)?;
195        return Ok(());
196    }
197
198    //
199    // At this point a more detailed analysis is needed.
200    //
201    //
202    // Assume the point is outside to start. This allows us
203    // to skip setting INSIDE when we find a boundary test
204    // failure.
205    //
206    *INSIDE = false;
207
208    //
209    // We'll use the shape parameters for latitude and altitude
210    // comparisons.
211    //
212    RE = CORPAR[1];
213    F = CORPAR[2];
214
215    //
216    // Get local copies of the coordinate bounds. Don't normalize the
217    // longitude bounds until we know we need them.
218    //
219    MINLAT = BOUNDS[[SOUTH, LATIDX]];
220    MAXLAT = BOUNDS[[NORTH, LATIDX]];
221
222    MINALT = BOUNDS[[LOWER, ALTIDX]];
223    MAXALT = BOUNDS[[UPPER, ALTIDX]];
224
225    //
226    // Compare coordinates to adjusted latitude boundaries.
227    //
228    if (EXCLUD != LATIDX) {
229        //
230        // Create adjusted latitude bounds.
231        //
232        AMNLAT = intrinsics::DMAX1(&[(-save.HPI - ANGMRG), (MINLAT - MARGIN)]);
233        AMXLAT = intrinsics::DMIN1(&[(save.HPI + ANGMRG), (MAXLAT + MARGIN)]);
234
235        //
236        // Compare the latitude of the input point to the bounds.
237        //
238        ZZPDCMPL(RE, F, P.as_slice(), AMNLAT, &mut RELMIN, ctx)?;
239        ZZPDCMPL(RE, F, P.as_slice(), AMXLAT, &mut RELMAX, ctx)?;
240
241        if FAILED(ctx) {
242            CHKOUT(b"ZZINPDT", ctx)?;
243            return Ok(());
244        }
245
246        if ((RELMIN == LT) || (RELMAX == GT)) {
247            //
248            // The latitude of P is strictly outside of the element's
249            // latitude bounds.
250            //
251            CHKOUT(b"ZZINPDT", ctx)?;
252            return Ok(());
253        }
254    }
255
256    //
257    // Test the point for inclusion in the region between the bounding
258    // ellipsoids that act as proxies for altitude boundaries.
259    //
260    if (EXCLUD != ALTIDX) {
261        //
262        // Extract altitude bounds from the segment descriptor.
263        //
264        MINALT = BOUNDS[[LOWER, ALTIDX]];
265        MAXALT = BOUNDS[[UPPER, ALTIDX]];
266        //
267        // Adjust altitude bounds to account for the margin.
268        //
269        AMNALT = (MINALT - (MARGIN * f64::abs(MINALT)));
270        AMXALT = (MAXALT + (MARGIN * f64::abs(MAXALT)));
271        //
272        // Set up a "boundary" array so that we can use ZZINPDT0
273        // to do the altitude check for us. We'll exclude longitude
274        // tests; the latitude test is set up for an automatic pass
275        // (the latitude range of ALTBDS is [-pi/2, pi/2]).
276        //
277        save.ALTBDS[[LOWER, ALTIDX]] = AMNALT;
278        save.ALTBDS[[UPPER, ALTIDX]] = AMXALT;
279
280        ZZINPDT0(
281            P.as_slice(),
282            LON,
283            save.ALTBDS.as_slice(),
284            CORPAR.as_slice(),
285            LONIDX,
286            &mut ALPASS,
287            ctx,
288        )?;
289
290        if !ALPASS {
291            CHKOUT(b"ZZINPDT", ctx)?;
292            return Ok(());
293        }
294    }
295
296    //
297    // At this point, the input altitude and latitude are within the
298    // adjusted bounds, if their tests haven't been excluded by
299    // the caller.
300    //
301    // Perform longitude tests, unless they're excluded by the
302    // caller.
303    //
304    if (EXCLUD != LONIDX) {
305        //
306        // Start out by normalizing the element's longitude bounds.
307        //
308        ZZNRMLON(
309            BOUNDS[[LOWER, LONIDX]],
310            BOUNDS[[UPPER, LONIDX]],
311            ANGMRG,
312            &mut MINLON,
313            &mut MAXLON,
314            ctx,
315        )?;
316
317        if FAILED(ctx) {
318            CHKOUT(b"ZZINPDT", ctx)?;
319            return Ok(());
320        }
321        //
322        // Set the margin to be used for longitude interval
323        // inclusion tests.
324        //
325        LONMRG = intrinsics::DMAX1(&[f64::abs(ANGMRG), f64::abs(MARGIN)]);
326
327        //
328        // We have a special case for segments that include the poles. If
329        // the latitude and altitude of the input point are within
330        // bounds, and if the latitude of the point is close enough to a
331        // pole we consider the point to be included in the segment,
332        // regardless of the point's longitude. All other points get the
333        // normal longitude test.
334        //
335        // We use planetocentric latitude to determine whether the point
336        // is "close" to a pole.
337        //
338        ZZPDCMPL(RE, F, P.as_slice(), (save.HPI - LATMRG), &mut RELMAX, ctx)?;
339        ZZPDCMPL(RE, F, P.as_slice(), (-save.HPI + LATMRG), &mut RELMIN, ctx)?;
340
341        if FAILED(ctx) {
342            CHKOUT(b"ZZINPDT", ctx)?;
343            return Ok(());
344        }
345
346        if ((RELMAX != GT) && (RELMIN != LT)) {
347            //
348            // This is the usual case: the latitude of the input point is
349            // bounded away from the poles.
350            //
351            // Check the point's longitude against the segment's longitude
352            // bounds.
353            //
354            // We'll scale the longitude margin to compensate for the
355            // latitude of the input point. Note that the division
356            // below is safe; presuming a reasonable value of MARGIN,
357            // we know that
358            //
359            //    DLON << 1
360            //
361            // Note that we use planetocentric latitude for scaling the
362            // longitude margin. This substitution (for planetodetic
363            // latitude) is adequate for this purpose.
364            //
365            DLON = (LONMRG / intrinsics::DMAX1(&[f64::abs(f64::cos(PCNLAT)), LATMRG]));
366
367            AMNLON = (MINLON - DLON);
368            AMXLON = (MAXLON + DLON);
369
370            //
371            // Now move the input point's longitude into range, if
372            // necessary, so we can make a valid comparison against
373            // the longitude bounds.
374            //
375            if (LON < AMNLON) {
376                if (LON < (AMNLON - LONALI)) {
377                    //
378                    // See whether an aliased version of LON is a match.
379                    //
380                    LON = (LON + save.PI2);
381                } else {
382                    //
383                    // Consider LON to be a match with the lower bound.
384                    //
385                    LON = AMNLON;
386                }
387            } else if (LON > AMXLON) {
388                if (LON > (AMXLON + LONALI)) {
389                    //
390                    // See whether an aliased version of LON is a match.
391                    //
392                    LON = (LON - save.PI2);
393                } else {
394                    //
395                    // Consider LON to be a match with the upper bound.
396                    //
397                    LON = AMXLON;
398                }
399            }
400
401            if ((LON < AMNLON) || (LON > AMXLON)) {
402                CHKOUT(b"ZZINPDT", ctx)?;
403                return Ok(());
404            }
405        } else {
406
407            //
408            // The latitude of the input point is close to one of the
409            // poles.
410            //
411            // This is a no-op case.
412            //
413            // The input point has already passed whichever of the
414            // altitude and latitude tests that were not excluded.
415            //
416            // If the element has a non-degenerate latitude boundary
417            // having the same sign as the latitude of the input point,
418            // and if latitude is excluded because the input point is
419            // already nominally on that boundary, then passing the
420            // altitude check implies that the point is close to the
421            // element.
422            //
423            // If the element has a degenerate latitude boundary having
424            // the same sign as the latitude of the input point---namely,
425            // the element contains the pole, and latitude is excluded,
426            // then then passing the altitude check implies that the point
427            // is close to the portion of the Z-axis contained in the
428            // element.
429            //
430            // If the altitude check has been excluded because the point
431            // is already nominally on one of the element's altitude
432            // boundaries, the passing the latitude test implies the point
433            // is close to the element.
434            //
435            // In all cases, as long as EXCLUD has been set appropriately,
436            // the point is close to the element. We consider the point to
437            // be in the expanded element.
438            //
439        }
440    }
441
442    //
443    // All tests that were commanded have been passed. The input
444    // point is considered to be contained in the expanded volume
445    // element.
446    //
447    *INSIDE = true;
448
449    CHKOUT(b"ZZINPDT", ctx)?;
450    Ok(())
451}