1 /* wndifd.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 WNDIFD ( Difference two DP windows ) */
wndifd_(doublereal * a,doublereal * b,doublereal * c__)9 /* Subroutine */ int wndifd_(doublereal *a, doublereal *b, doublereal *c__)
10 {
11 /* System generated locals */
12 integer i__1;
13
14 /* Local variables */
15 logical keep;
16 integer over;
17 doublereal f;
18 integer acard, bcard;
19 doublereal l;
20 extern integer cardd_(doublereal *);
21 extern /* Subroutine */ int chkin_(char *, ftnlen);
22 integer csize;
23 extern integer sized_(doublereal *);
24 extern /* Subroutine */ int copyd_(doublereal *, doublereal *);
25 integer needed;
26 extern /* Subroutine */ int scardd_(integer *, doublereal *), sigerr_(
27 char *, ftnlen), chkout_(char *, ftnlen), ssized_(integer *,
28 doublereal *), setmsg_(char *, ftnlen), errint_(char *, integer *,
29 ftnlen);
30 extern logical return_(void);
31 logical unrslv;
32 integer apb, bpb, ape, bpe, put;
33
34 /* $ Abstract */
35
36 /* Place the difference of two double precision windows into */
37 /* a third window. */
38
39 /* $ Disclaimer */
40
41 /* THIS SOFTWARE AND ANY RELATED MATERIALS WERE CREATED BY THE */
42 /* CALIFORNIA INSTITUTE OF TECHNOLOGY (CALTECH) UNDER A U.S. */
43 /* GOVERNMENT CONTRACT WITH THE NATIONAL AERONAUTICS AND SPACE */
44 /* ADMINISTRATION (NASA). THE SOFTWARE IS TECHNOLOGY AND SOFTWARE */
45 /* PUBLICLY AVAILABLE UNDER U.S. EXPORT LAWS AND IS PROVIDED "AS-IS" */
46 /* TO THE RECIPIENT WITHOUT WARRANTY OF ANY KIND, INCLUDING ANY */
47 /* WARRANTIES OF PERFORMANCE OR MERCHANTABILITY OR FITNESS FOR A */
48 /* PARTICULAR USE OR PURPOSE (AS SET FORTH IN UNITED STATES UCC */
49 /* SECTIONS 2312-2313) OR FOR ANY PURPOSE WHATSOEVER, FOR THE */
50 /* SOFTWARE AND RELATED MATERIALS, HOWEVER USED. */
51
52 /* IN NO EVENT SHALL CALTECH, ITS JET PROPULSION LABORATORY, OR NASA */
53 /* BE LIABLE FOR ANY DAMAGES AND/OR COSTS, INCLUDING, BUT NOT */
54 /* LIMITED TO, INCIDENTAL OR CONSEQUENTIAL DAMAGES OF ANY KIND, */
55 /* INCLUDING ECONOMIC DAMAGE OR INJURY TO PROPERTY AND LOST PROFITS, */
56 /* REGARDLESS OF WHETHER CALTECH, JPL, OR NASA BE ADVISED, HAVE */
57 /* REASON TO KNOW, OR, IN FACT, SHALL KNOW OF THE POSSIBILITY. */
58
59 /* RECIPIENT BEARS ALL RISK RELATING TO QUALITY AND PERFORMANCE OF */
60 /* THE SOFTWARE AND ANY RELATED MATERIALS, AND AGREES TO INDEMNIFY */
61 /* CALTECH AND NASA FOR ALL THIRD-PARTY CLAIMS RESULTING FROM THE */
62 /* ACTIONS OF RECIPIENT IN THE USE OF THE SOFTWARE. */
63
64 /* $ Required_Reading */
65
66 /* WINDOWS */
67
68 /* $ Keywords */
69
70 /* WINDOWS */
71
72 /* $ Declarations */
73 /* $ Brief_I/O */
74
75 /* VARIABLE I/O DESCRIPTION */
76 /* -------- --- -------------------------------------------------- */
77 /* A, */
78 /* B I Input windows. */
79 /* C I Difference of A and B. */
80
81 /* $ Detailed_Input */
82
83 /* A, */
84 /* B are windows, each of which contains zero or more */
85 /* intervals. */
86
87 /* $ Detailed_Output */
88
89 /* C is the output window, containing the difference */
90 /* of A and B---every point contained in A, but not */
91 /* contained in B. */
92
93 /* C must be distinct from both A and B. */
94
95 /* $ Parameters */
96
97 /* None. */
98
99 /* $ Exceptions */
100
101 /* 1. If the difference of the two windows results in an excess of */
102 /* elements, the error SPICE(WINDOWEXCESS) is signalled. */
103
104 /* $ Files */
105
106 /* None. */
107
108 /* $ Particulars */
109
110 /* Mathematically, the difference of two windows contains every */
111 /* point contained in the first window but not contained in the */
112 /* second window. */
113
114 /* Fortran offers no satisfactory floating point representation */
115 /* of open intervals. Thus, for floating point windows we must */
116 /* return the closure of the set theoretical difference: that is, */
117 /* the difference plus the endpoints of the first window that are */
118 /* contained in the second window. */
119
120 /* $ Examples */
121
122 /* Let A contain the intervals */
123
124 /* [ 1, 3 ] [ 7, 11 ] [ 23, 27 ] */
125
126 /* and B contain the intervals */
127
128 /* [ 2, 4 ] [ 8, 10 ] [ 16, 18 ] */
129
130 /* Then the difference of A and B contains the intervals */
131
132 /* [ 1, 2 ] [ 7, 8 ] [ 10, 11 ] [ 23, 27 ] */
133
134 /* $ Restrictions */
135
136 /* None. */
137
138 /* $ Literature_References */
139
140 /* None. */
141
142 /* $ Author_and_Institution */
143
144 /* H.A. Neilan (JPL) */
145 /* W.L. Taber (JPL) */
146 /* I.M. Underwood (JPL) */
147
148 /* $ Version */
149
150 /* - SPICELIB Version 2.0.0, 16-SEP-1998 (WLT) */
151
152 /* The previous version did not work when removing */
153 /* singletons. This has been corrected. */
154
155 /* - SPICELIB Version 1.0.1, 10-MAR-1992 (WLT) */
156
157 /* Comment section for permuted index source lines was added */
158 /* following the header. */
159
160 /* - SPICELIB Version 1.0.0, 31-JAN-1990 (WLT) (IMU) */
161
162 /* -& */
163 /* $ Index_Entries */
164
165 /* difference two d.p. windows */
166
167 /* -& */
168 /* $ Revisions */
169
170 /* - Beta Version 1.1.0, 27-FEB-1989 (HAN) */
171
172 /* Due to the calling sequence and functionality changes */
173 /* in the routine EXCESS, the method of signalling an */
174 /* excess of elements needed to be changed. */
175
176 /* -& */
177
178 /* SPICELIB functions */
179
180
181 /* Local variables */
182
183
184 /* Standard SPICE error handling. */
185
186 if (return_()) {
187 return 0;
188 }
189 chkin_("WNDIFD", (ftnlen)6);
190
191 /* Find the cardinality of the input windows, and the allowed size */
192 /* of the output window. Also, save the size of the second window. */
193
194 acard = cardd_(a);
195 bcard = cardd_(b);
196 csize = sized_(c__);
197 over = 0;
198
199 /* Empty out the output window. */
200
201 ssized_(&csize, c__);
202
203 /* Let's handle the pathological cases first. */
204
205 if (bcard == 0) {
206 copyd_(a, c__);
207 chkout_("WNDIFD", (ftnlen)6);
208 return 0;
209 } else if (acard == 0) {
210 chkout_("WNDIFD", (ftnlen)6);
211 return 0;
212 }
213
214 /* Now get pointers to the first intervals of A and B. */
215
216 apb = 1;
217 ape = 2;
218 bpb = 1;
219 bpe = 2;
220 put = 1;
221
222 /* As long as the endpointer for A is less than the cardinality */
223 /* of A we need to examine intervals and decide how much of */
224 /* them to keep in C. */
225
226 while(ape <= acard) {
227
228 /* We will work with the interval [F,L] which starts out */
229 /* as the next interval of A. We modify it below as required */
230 /* when subtracting out intervals of B. */
231
232 f = a[apb + 5];
233 l = a[ape + 5];
234
235 /* Right now we have not resolved whether to keep the interval */
236 /* [F,L], but until we know better we assume it is a keeper. */
237
238 unrslv = bpe <= bcard;
239 keep = TRUE_;
240 while(unrslv) {
241 if (l < b[bpb + 5]) {
242
243 /* The interval [F,L] is before the next interval of B, we */
244 /* have resolved what to do with this one. It is a */
245 /* keeper. */
246
247 unrslv = FALSE_;
248 } else if (f > b[bpe + 5]) {
249
250 /* [F,L] is after the end of the current interval in B, */
251 /* we need to look at the next interval of B */
252
253 bpb += 2;
254 bpe += 2;
255 unrslv = bpe <= bcard;
256 } else {
257
258 /* There is some overlap between the current interval */
259 /* of B and the current interval of A. There are */
260 /* several possibilities */
261
262 /* 1) The current interval of A is contained in the */
263 /* current interval of B (This includes singleton */
264 /* intervals in A). We just mark [F,L] so that it */
265 /* won't be kept. We have fully resolved what to */
266 /* do with [F,L]. */
267
268 /* 2) The interval from B overlaps at the beginning */
269 /* of the interval of A */
270
271 /* B interval [......] */
272 /* A interval [............] */
273 /* result of A-B [.........] */
274
275 /* In this case we need to shrink the interval [F,L] */
276 /* but we have not resolved how much of the result */
277 /* to keep. */
278
279 /* 3) The interval from B falls inside the current */
280 /* interval [F,L] */
281
282 /* B interval [......] */
283 /* A interval [............] */
284 /* result of A-B [..] [..] */
285
286 /* If the interval from B is not a singleton, we store */
287 /* the first part of [F,L] in C and then set [F,L] to */
288 /* be the right interval which is still not resolved. */
289
290 /* If the B interval is a singleton we can ignore ignore */
291 /* it. But we have not resolved what to do about */
292 /* [F,L], we need to look at the next interval of B. */
293
294
295 /* 4) The interval from B overlaps at the ending */
296 /* of the interval of A */
297
298 /* B interval [......] */
299 /* A interval [......] */
300 /* result of A-B [....] */
301
302 /* We need to shrink [F,L]. In this case we know we can */
303 /* keep all of what's left because all other intervals */
304 /* of B are to the right of [F,L] */
305
306 if (b[bpb + 5] <= f && l <= b[bpe + 5]) {
307
308 /* Case 1 above */
309
310 keep = FALSE_;
311 unrslv = FALSE_;
312 } else if (b[bpb + 5] <= f) {
313
314 /* Case 2 above */
315
316 f = b[bpe + 5];
317 bpb += 2;
318 bpe += 2;
319 unrslv = bpe <= bcard;
320 } else if (f <= b[bpb + 5] && l >= b[bpe + 5] && b[bpb + 5] <
321 b[bpe + 5]) {
322
323 /* Case 3 above (non-singleton interval of B). */
324
325 if (put < csize) {
326 c__[put + 5] = f;
327 c__[put + 6] = b[bpb + 5];
328 i__1 = put + 1;
329 scardd_(&i__1, c__);
330 put += 2;
331 } else {
332 over += 2;
333 }
334 f = b[bpe + 5];
335
336 /* If the interval from B contained L, we will not */
337 /* want to be keeping the singleton [F,L]. */
338
339 if (f == l) {
340 keep = FALSE_;
341 unrslv = FALSE_;
342 }
343 bpb += 2;
344 bpe += 2;
345 unrslv = unrslv && bpe <= bcard;
346 } else if (f <= b[bpb + 5] && l >= b[bpe + 5] && b[bpb + 5] ==
347 b[bpe + 5]) {
348
349 /* Case 3 above (singleton interval of B). */
350
351 bpb += 2;
352 bpe += 2;
353 unrslv = bpe <= bcard;
354 } else {
355
356 /* Case 4 above */
357
358 l = b[bpb + 5];
359 unrslv = FALSE_;
360 }
361 }
362 }
363
364 /* If there is anything to keep in C, put it there. */
365
366 if (keep) {
367
368 /* Make sure there is sufficient room to do the putting. */
369
370 if (put < csize) {
371 c__[put + 5] = f;
372 c__[put + 6] = l;
373 i__1 = put + 1;
374 scardd_(&i__1, c__);
375 put += 2;
376 } else {
377 over += 2;
378 }
379 }
380
381 /* Move the pointers in A to the next interval. */
382
383 apb += 2;
384 ape += 2;
385 }
386
387 /* We've examined all of the intervals of A and B, but if we */
388 /* didn't actually store all of the difference, signal an error. */
389
390 if (over > 0) {
391 needed = over + csize;
392 setmsg_("The output window did not have sufficient room to contain t"
393 "he result of the window difference. It has room for # endpo"
394 "ints, but # were needed to describe the difference. ", (
395 ftnlen)171);
396 errint_("#", &csize, (ftnlen)1);
397 errint_("#", &needed, (ftnlen)1);
398 sigerr_("SPICE(WINDOWEXCESS)", (ftnlen)19);
399 }
400 chkout_("WNDIFD", (ftnlen)6);
401 return 0;
402 } /* wndifd_ */
403
404