1 : /* pack_gp.f -- translated by f2c (version 20031025).
2 : You must link the resulting object file with libf2c:
3 : on Microsoft Windows system, link with libf2c.lib;
4 : on Linux or Unix systems, link with .../path/to/libf2c.a -lm
5 : or, if you install libf2c.a in a standard place, with -lf2c -lm
6 : -- in that order, at the end of the command line, as in
7 : cc *.o -lf2c -lm
8 : Source for libf2c is in /netlib/f2c/libf2c.zip, e.g.,
9 :
10 : http://www.netlib.org/f2c/libf2c.zip
11 : */
12 :
13 : /*#include "f2c.h"*/
14 : #include <stdlib.h>
15 : #include "grib2.h"
16 : typedef g2int integer;
17 : typedef g2int logical;
18 : #define TRUE_ (1)
19 : #define FALSE_ (0)
20 :
21 0 : /* Subroutine */ int pack_gp(integer *kfildo, integer *ic, integer *nxy,
22 : integer *is523, integer *minpk, integer *inc, integer *missp, integer
23 : *misss, integer *jmin, integer *jmax, integer *lbit, integer *nov,
24 : integer *ndg, integer *lx, integer *ibit, integer *jbit, integer *
25 : kbit, integer *novref, integer *lbitref, integer *ier)
26 : {
27 : /* Initialized data */
28 :
29 0 : const integer mallow = 1073741825; /* MALLOW=2**30+1 */
30 : static integer ifeed = 12;
31 : static integer ifirst = 0;
32 :
33 : /* System generated locals */
34 : integer i__1, i__2, i__3;
35 :
36 : /* Local variables */
37 : static integer j, k, l;
38 : static logical adda;
39 : static integer ired, kinc, mina, maxa, minb, maxb, minc, maxc, ibxx2[31];
40 : static char cfeed[1];
41 : static integer nenda, nendb, ibita, ibitb, minak, minbk, maxak, maxbk,
42 : minck, maxck, nouta, lmiss, itest, nount;
43 : extern /* Subroutine */ int reduce(integer *, integer *, integer *,
44 : integer *, integer *, integer *, integer *, integer *, integer *,
45 : integer *, integer *, integer *, integer *);
46 : static integer ibitbs, mislla, misllb, misllc, iersav, lminpk, ktotal,
47 : kounta, kountb, kstart, mstart, mintst, maxtst,
48 : kounts, mintstk, maxtstk;
49 : integer *misslx;
50 :
51 :
52 : /* FEBRUARY 1994 GLAHN TDL MOS-2000 */
53 : /* JUNE 1995 GLAHN MODIFIED FOR LMISS ERROR. */
54 : /* JULY 1996 GLAHN ADDED MISSS */
55 : /* FEBRUARY 1997 GLAHN REMOVED 4 REDUNDANT TESTS FOR */
56 : /* MISSP.EQ.0; INSERTED A TEST TO BETTER */
57 : /* HANDLE A STRING OF 9999'S */
58 : /* FEBRUARY 1997 GLAHN ADDED LOOPS TO ELIMINATE TEST FOR */
59 : /* MISSS WHEN MISSS = 0 */
60 : /* MARCH 1997 GLAHN CORRECTED FOR SECONDARY MISSING VALUE */
61 : /* MARCH 1997 GLAHN CORRECTED FOR USE OF LOCAL VALUE */
62 : /* OF MINPK */
63 : /* MARCH 1997 GLAHN CORRECTED FOR SECONDARY MISSING VALUE */
64 : /* MARCH 1997 GLAHN CHANGED CALCULATING NUMBER OF BITS */
65 : /* THROUGH EXPONENTS TO AN ARRAY (IMPROVED */
66 : /* OVERALL PACKING PERFORMANCE BY ABOUT */
67 : /* 35 PERCENT!). ALLOWED 0 BITS FOR */
68 : /* PACKING JMIN( ), LBIT( ), AND NOV( ). */
69 : /* MAY 1997 GLAHN A NUMBER OF CHANGES FOR EFFICIENCY. */
70 : /* MOD FUNCTIONS ELIMINATED AND ONE */
71 : /* IFTHEN ADDED. JOUNT REMOVED. */
72 : /* RECOMPUTATION OF BITS NOT MADE UNLESS */
73 : /* NECESSARY AFTER MOVING POINTS FROM */
74 : /* ONE GROUP TO ANOTHER. NENDB ADJUSTED */
75 : /* TO ELIMINATE POSSIBILITY OF VERY */
76 : /* SMALL GROUP AT THE END. */
77 : /* ABOUT 8 PERCENT IMPROVEMENT IN */
78 : /* OVERALL PACKING. ISKIPA REMOVED; */
79 : /* THERE IS ALWAYS A GROUP B THAT CAN */
80 : /* BECOME GROUP A. CONTROL ON SIZE */
81 : /* OF GROUP B (STATEMENT BELOW 150) */
82 : /* ADDED. ADDED ADDA, AND USE */
83 : /* OF GE AND LE INSTEAD OF GT AND LT */
84 : /* IN LOOPS BETWEEN 150 AND 160. */
85 : /* IBITBS ADDED TO SHORTEN TRIPS */
86 : /* THROUGH LOOP. */
87 : /* MARCH 2000 GLAHN MODIFIED FOR GRIB2; CHANGED NAME FROM */
88 : /* PACKGP */
89 : /* JANUARY 2001 GLAHN COMMENTS; IER = 706 SUBSTITUTED FOR */
90 : /* STOPS; ADDED RETURN1; REMOVED STATEMENT */
91 : /* NUMBER 110; ADDED IER AND * RETURN */
92 : /* NOVEMBER 2001 GLAHN CHANGED SOME DIAGNOSTIC FORMATS TO */
93 : /* ALLOW PRINTING LARGER NUMBERS */
94 : /* NOVEMBER 2001 GLAHN ADDED MISSLX( ) TO PUT MAXIMUM VALUE */
95 : /* INTO JMIN( ) WHEN ALL VALUES MISSING */
96 : /* TO AGREE WITH GRIB STANDARD. */
97 : /* NOVEMBER 2001 GLAHN CHANGED TWO TESTS ON MISSP AND MISSS */
98 : /* EQ 0 TO TESTS ON IS523. HOWEVER, */
99 : /* MISSP AND MISSS CANNOT IN GENERAL BE */
100 : /* = 0. */
101 : /* NOVEMBER 2001 GLAHN ADDED CALL TO REDUCE; DEFINED ITEST */
102 : /* BEFORE LOOPS TO REDUCE COMPUTATION; */
103 : /* STARTED LARGE GROUP WHEN ALL SAME */
104 : /* VALUE */
105 : /* DECEMBER 2001 GLAHN MODIFIED AND ADDED A FEW COMMENTS */
106 : /* JANUARY 2002 GLAHN REMOVED LOOP BEFORE 150 TO DETERMINE */
107 : /* A GROUP OF ALL SAME VALUE */
108 : /* JANUARY 2002 GLAHN CHANGED MALLOW FROM 9999999 TO 2**30+1, */
109 : /* AND MADE IT A PARAMETER */
110 : /* MARCH 2002 GLAHN ADDED NON FATAL IER = 716, 717; */
111 : /* REMOVED NENDB=NXY ABOVE 150; */
112 : /* ADDED IERSAV=0; COMMENTS */
113 :
114 : /* PURPOSE */
115 : /* DETERMINES GROUPS OF VARIABLE SIZE, BUT AT LEAST OF */
116 : /* SIZE MINPK, THE ASSOCIATED MAX (JMAX( )) AND MIN (JMIN( )), */
117 : /* THE NUMBER OF BITS NECESSARY TO HOLD THE VALUES IN EACH */
118 : /* GROUP (LBIT( )), THE NUMBER OF VALUES IN EACH GROUP */
119 : /* (NOV( )), THE NUMBER OF BITS NECESSARY TO PACK THE JMIN( ) */
120 : /* VALUES (IBIT), THE NUMBER OF BITS NECESSARY TO PACK THE */
121 : /* LBIT( ) VALUES (JBIT), AND THE NUMBER OF BITS NECESSARY */
122 : /* TO PACK THE NOV( ) VALUES (KBIT). THE ROUTINE IS DESIGNED */
123 : /* TO DETERMINE THE GROUPS SUCH THAT A SMALL NUMBER OF BITS */
124 : /* IS NECESSARY TO PACK THE DATA WITHOUT EXCESSIVE */
125 : /* COMPUTATIONS. IF ALL VALUES IN THE GROUP ARE ZERO, THE */
126 : /* NUMBER OF BITS TO USE IN PACKING IS DEFINED AS ZERO WHEN */
127 : /* THERE CAN BE NO MISSING VALUES; WHEN THERE CAN BE MISSING */
128 : /* VALUES, THE NUMBER OF BITS MUST BE AT LEAST 1 TO HAVE */
129 : /* THE CAPABILITY TO RECOGNIZE THE MISSING VALUE. HOWEVER, */
130 : /* IF ALL VALUES IN A GROUP ARE MISSING, THE NUMBER OF BITS */
131 : /* NEEDED IS 0, AND THE UNPACKER RECOGNIZES THIS. */
132 : /* ALL VARIABLES ARE INTEGER. EVEN THOUGH THE GROUPS ARE */
133 : /* INITIALLY OF SIZE MINPK OR LARGER, AN ADJUSTMENT BETWEEN */
134 : /* TWO GROUPS (THE LOOKBACK PROCEDURE) MAY MAKE A GROUP */
135 : /* SMALLER THAN MINPK. THE CONTROL ON GROUP SIZE IS THAT */
136 : /* THE SUM OF THE SIZES OF THE TWO CONSECUTIVE GROUPS, EACH OF */
137 : /* SIZE MINPK OR LARGER, IS NOT DECREASED. WHEN DETERMINING */
138 : /* THE NUMBER OF BITS NECESSARY FOR PACKING, THE LARGEST */
139 : /* VALUE THAT CAN BE ACCOMMODATED IN, SAY, MBITS, IS */
140 : /* 2**MBITS-1; THIS LARGEST VALUE (AND THE NEXT SMALLEST */
141 : /* VALUE) IS RESERVED FOR THE MISSING VALUE INDICATOR (ONLY) */
142 : /* WHEN IS523 NE 0. IF THE DIMENSION NDG */
143 : /* IS NOT LARGE ENOUGH TO HOLD ALL THE GROUPS, THE LOCAL VALUE */
144 : /* OF MINPK IS INCREASED BY 50 PERCENT. THIS IS REPEATED */
145 : /* UNTIL NDG WILL SUFFICE. A DIAGNOSTIC IS PRINTED WHENEVER */
146 : /* THIS HAPPENS, WHICH SHOULD BE VERY RARELY. IF IT HAPPENS */
147 : /* OFTEN, NDG IN SUBROUTINE PACK SHOULD BE INCREASED AND */
148 : /* A CORRESPONDING INCREASE IN SUBROUTINE UNPACK MADE. */
149 : /* CONSIDERABLE CODE IS PROVIDED SO THAT NO MORE CHECKING */
150 : /* FOR MISSING VALUES WITHIN LOOPS IS DONE THAN NECESSARY; */
151 : /* THE ADDED EFFICIENCY OF THIS IS RELATIVELY MINOR, */
152 : /* BUT DOES NO HARM. FOR GRIB2, THE REFERENCE VALUE FOR */
153 : /* THE LENGTH OF GROUPS IN NOV( ) AND FOR THE NUMBER OF */
154 : /* BITS NECESSARY TO PACK GROUP VALUES ARE DETERMINED, */
155 : /* AND SUBTRACTED BEFORE JBIT AND KBIT ARE DETERMINED. */
156 :
157 : /* WHEN 1 OR MORE GROUPS ARE LARGE COMPARED TO THE OTHERS, */
158 : /* THE WIDTH OF ALL GROUPS MUST BE AS LARGE AS THE LARGEST. */
159 : /* A SUBROUTINE REDUCE BREAKS UP LARGE GROUPS INTO 2 OR */
160 : /* MORE TO REDUCE TOTAL BITS REQUIRED. IF REDUCE SHOULD */
161 : /* ABORT, PACK_GP WILL BE EXECUTED AGAIN WITHOUT THE CALL */
162 : /* TO REDUCE. */
163 :
164 : /* DATA SET USE */
165 : /* KFILDO - UNIT NUMBER FOR OUTPUT (PRINT) FILE. (OUTPUT) */
166 :
167 : /* VARIABLES IN CALL SEQUENCE */
168 : /* KFILDO = UNIT NUMBER FOR OUTPUT (PRINT) FILE. (INPUT) */
169 : /* IC( ) = ARRAY TO HOLD DATA FOR PACKING. THE VALUES */
170 : /* DO NOT HAVE TO BE POSITIVE AT THIS POINT, BUT */
171 : /* MUST BE IN THE RANGE -2**30 TO +2**30 (THE */
172 : /* THE VALUE OF MALLOW). THESE INTEGER VALUES */
173 : /* WILL BE RETAINED EXACTLY THROUGH PACKING AND */
174 : /* UNPACKING. (INPUT) */
175 : /* NXY = NUMBER OF VALUES IN IC( ). ALSO TREATED */
176 : /* AS ITS DIMENSION. (INPUT) */
177 : /* IS523 = missing value management */
178 : /* 0=data contains no missing values */
179 : /* 1=data contains Primary missing values */
180 : /* 2=data contains Primary and secondary missing values */
181 : /* (INPUT) */
182 : /* MINPK = THE MINIMUM SIZE OF EACH GROUP, EXCEPT POSSIBLY */
183 : /* THE LAST ONE. (INPUT) */
184 : /* INC = THE NUMBER OF VALUES TO ADD TO AN ALREADY */
185 : /* EXISTING GROUP IN DETERMINING WHETHER OR NOT */
186 : /* TO START A NEW GROUP. IDEALLY, THIS WOULD BE */
187 : /* 1, BUT EACH TIME INC VALUES ARE ATTEMPTED, THE */
188 : /* MAX AND MIN OF THE NEXT MINPK VALUES MUST BE */
189 : /* FOUND. THIS IS "A LOOP WITHIN A LOOP," AND */
190 : /* A SLIGHTLY LARGER VALUE MAY GIVE ABOUT AS GOOD */
191 : /* RESULTS WITH SLIGHTLY LESS COMPUTATIONAL TIME. */
192 : /* IF INC IS LE 0, 1 IS USED, AND A DIAGNOSTIC IS */
193 : /* OUTPUT. NOTE: IT IS EXPECTED THAT INC WILL */
194 : /* EQUAL 1. THE CODE USES INC PRIMARILY IN THE */
195 : /* LOOPS STARTING AT STATEMENT 180. IF INC */
196 : /* WERE 1, THERE WOULD NOT NEED TO BE LOOPS */
197 : /* AS SUCH. HOWEVER, KINC (THE LOCAL VALUE OF */
198 : /* INC) IS SET GE 1 WHEN NEAR THE END OF THE DATA */
199 : /* TO FORESTALL A VERY SMALL GROUP AT THE END. */
200 : /* (INPUT) */
201 : /* MISSP = WHEN MISSING POINTS CAN BE PRESENT IN THE DATA, */
202 : /* THEY WILL HAVE THE VALUE MISSP OR MISSS. */
203 : /* MISSP IS THE PRIMARY MISSING VALUE AND MISSS */
204 : /* IS THE SECONDARY MISSING VALUE . THESE MUST */
205 : /* NOT BE VALUES THAT WOULD OCCUR WITH SUBTRACTING */
206 : /* THE MINIMUM (REFERENCE) VALUE OR SCALING. */
207 : /* FOR EXAMPLE, MISSP = 0 WOULD NOT BE ADVISABLE. */
208 : /* (INPUT) */
209 : /* MISSS = SECONDARY MISSING VALUE INDICATOR (SEE MISSP). */
210 : /* (INPUT) */
211 : /* JMIN(J) = THE MINIMUM OF EACH GROUP (J=1,LX). (OUTPUT) */
212 : /* JMAX(J) = THE MAXIMUM OF EACH GROUP (J=1,LX). THIS IS */
213 : /* NOT REALLY NEEDED, BUT SINCE THE MAX OF EACH */
214 : /* GROUP MUST BE FOUND, SAVING IT HERE IS CHEAP */
215 : /* IN CASE THE USER WANTS IT. (OUTPUT) */
216 : /* LBIT(J) = THE NUMBER OF BITS NECESSARY TO PACK EACH GROUP */
217 : /* (J=1,LX). IT IS ASSUMED THE MINIMUM OF EACH */
218 : /* GROUP WILL BE REMOVED BEFORE PACKING, AND THE */
219 : /* VALUES TO PACK WILL, THEREFORE, ALL BE POSITIVE. */
220 : /* HOWEVER, IC( ) DOES NOT NECESSARILY CONTAIN */
221 : /* ALL POSITIVE VALUES. IF THE OVERALL MINIMUM */
222 : /* HAS BEEN REMOVED (THE USUAL CASE), THEN IC( ) */
223 : /* WILL CONTAIN ONLY POSITIVE VALUES. (OUTPUT) */
224 : /* NOV(J) = THE NUMBER OF VALUES IN EACH GROUP (J=1,LX). */
225 : /* (OUTPUT) */
226 : /* NDG = THE DIMENSION OF JMIN( ), JMAX( ), LBIT( ), AND */
227 : /* NOV( ). (INPUT) */
228 : /* LX = THE NUMBER OF GROUPS DETERMINED. (OUTPUT) */
229 : /* IBIT = THE NUMBER OF BITS NECESSARY TO PACK THE JMIN(J) */
230 : /* VALUES, J=1,LX. (OUTPUT) */
231 : /* JBIT = THE NUMBER OF BITS NECESSARY TO PACK THE LBIT(J) */
232 : /* VALUES, J=1,LX. (OUTPUT) */
233 : /* KBIT = THE NUMBER OF BITS NECESSARY TO PACK THE NOV(J) */
234 : /* VALUES, J=1,LX. (OUTPUT) */
235 : /* NOVREF = REFERENCE VALUE FOR NOV( ). (OUTPUT) */
236 : /* LBITREF = REFERENCE VALUE FOR LBIT( ). (OUTPUT) */
237 : /* IER = ERROR RETURN. */
238 : /* 706 = VALUE WILL NOT PACK IN 30 BITS--FATAL */
239 : /* 714 = ERROR IN REDUCE--NON-FATAL */
240 : /* 715 = NGP NOT LARGE ENOUGH IN REDUCE--NON-FATAL */
241 : /* 716 = MINPK INCEASED--NON-FATAL */
242 : /* 717 = INC SET = 1--NON-FATAL */
243 : /* (OUTPUT) */
244 : /* * = ALTERNATE RETURN WHEN IER NE 0 AND FATAL ERROR. */
245 :
246 : /* INTERNAL VARIABLES */
247 : /* CFEED = CONTAINS THE CHARACTER REPRESENTATION */
248 : /* OF A PRINTER FORM FEED. */
249 : /* IFEED = CONTAINS THE INTEGER VALUE OF A PRINTER */
250 : /* FORM FEED. */
251 : /* KINC = WORKING COPY OF INC. MAY BE MODIFIED. */
252 : /* MINA = MINIMUM VALUE IN GROUP A. */
253 : /* MAXA = MAXIMUM VALUE IN GROUP A. */
254 : /* NENDA = THE PLACE IN IC( ) WHERE GROUP A ENDS. */
255 : /* KSTART = THE PLACE IN IC( ) WHERE GROUP A STARTS. */
256 : /* IBITA = NUMBER OF BITS NEEDED TO HOLD VALUES IN GROUP A. */
257 : /* MINB = MINIMUM VALUE IN GROUP B. */
258 : /* MAXB = MAXIMUM VALUE IN GROUP B. */
259 : /* NENDB = THE PLACE IN IC( ) WHERE GROUP B ENDS. */
260 : /* IBITB = NUMBER OF BITS NEEDED TO HOLD VALUES IN GROUP B. */
261 : /* MINC = MINIMUM VALUE IN GROUP C. */
262 : /* MAXC = MAXIMUM VALUE IN GROUP C. */
263 : /* KTOTAL = COUNT OF NUMBER OF VALUES IN IC( ) PROCESSED. */
264 : /* NOUNT = NUMBER OF VALUES ADDED TO GROUP A. */
265 : /* LMISS = 0 WHEN IS523 = 0. WHEN PACKING INTO A */
266 : /* SPECIFIC NUMBER OF BITS, SAY MBITS, */
267 : /* THE MAXIMUM VALUE THAT CAN BE HANDLED IS */
268 : /* 2**MBITS-1. WHEN IS523 = 1, INDICATING */
269 : /* PRIMARY MISSING VALUES, THIS MAXIMUM VALUE */
270 : /* IS RESERVED TO HOLD THE PRIMARY MISSING VALUE */
271 : /* INDICATOR AND LMISS = 1. WHEN IS523 = 2, */
272 : /* THE VALUE JUST BELOW THE MAXIMUM (I.E., */
273 : /* 2**MBITS-2) IS RESERVED TO HOLD THE SECONDARY */
274 : /* MISSING VALUE INDICATOR AND LMISS = 2. */
275 : /* LMINPK = LOCAL VALUE OF MINPK. THIS WILL BE ADJUSTED */
276 : /* UPWARD WHENEVER NDG IS NOT LARGE ENOUGH TO HOLD */
277 : /* ALL THE GROUPS. */
278 : /* MALLOW = THE LARGEST ALLOWABLE VALUE FOR PACKING. */
279 : /* MISLLA = SET TO 1 WHEN ALL VALUES IN GROUP A ARE MISSING. */
280 : /* THIS IS USED TO DISTINGUISH BETWEEN A REAL */
281 : /* MINIMUM WHEN ALL VALUES ARE NOT MISSING */
282 : /* AND A MINIMUM THAT HAS BEEN SET TO ZERO WHEN */
283 : /* ALL VALUES ARE MISSING. 0 OTHERWISE. */
284 : /* NOTE THAT THIS DOES NOT DISTINGUISH BETWEEN */
285 : /* PRIMARY AND SECONDARY MISSINGS WHEN SECONDARY */
286 : /* MISSINGS ARE PRESENT. THIS MEANS THAT */
287 : /* LBIT( ) WILL NOT BE ZERO WITH THE RESULTING */
288 : /* COMPRESSION EFFICIENCY WHEN SECONDARY MISSINGS */
289 : /* ARE PRESENT. ALSO NOTE THAT A CHECK HAS BEEN */
290 : /* MADE EARLIER TO DETERMINE THAT SECONDARY */
291 : /* MISSINGS ARE REALLY THERE. */
292 : /* MISLLB = SET TO 1 WHEN ALL VALUES IN GROUP B ARE MISSING. */
293 : /* THIS IS USED TO DISTINGUISH BETWEEN A REAL */
294 : /* MINIMUM WHEN ALL VALUES ARE NOT MISSING */
295 : /* AND A MINIMUM THAT HAS BEEN SET TO ZERO WHEN */
296 : /* ALL VALUES ARE MISSING. 0 OTHERWISE. */
297 : /* MISLLC = PERFORMS THE SAME FUNCTION FOR GROUP C THAT */
298 : /* MISLLA AND MISLLB DO FOR GROUPS B AND C, */
299 : /* RESPECTIVELY. */
300 : /* IBXX2(J) = AN ARRAY THAT WHEN THIS ROUTINE IS FIRST ENTERED */
301 : /* IS SET TO 2**J, J=0,30. IBXX2(30) = 2**30, WHICH */
302 : /* IS THE LARGEST VALUE PACKABLE, BECAUSE 2**31 */
303 : /* IS LARGER THAN THE INTEGER WORD SIZE. */
304 : /* IFIRST = SET BY DATA STATEMENT TO 0. CHANGED TO 1 ON */
305 : /* FIRST */
306 : /* ENTRY WHEN IBXX2( ) IS FILLED. */
307 : /* MINAK = KEEPS TRACK OF THE LOCATION IN IC( ) WHERE THE */
308 : /* MINIMUM VALUE IN GROUP A IS LOCATED. */
309 : /* MAXAK = DOES THE SAME AS MINAK, EXCEPT FOR THE MAXIMUM. */
310 : /* MINBK = THE SAME AS MINAK FOR GROUP B. */
311 : /* MAXBK = THE SAME AS MAXAK FOR GROUP B. */
312 : /* MINCK = THE SAME AS MINAK FOR GROUP C. */
313 : /* MAXCK = THE SAME AS MAXAK FOR GROUP C. */
314 : /* ADDA = KEEPS TRACK WHETHER OR NOT AN ATTEMPT TO ADD */
315 : /* POINTS TO GROUP A WAS MADE. IF SO, THEN ADDA */
316 : /* KEEPS FROM TRYING TO PUT ONE BACK INTO B. */
317 : /* (LOGICAL) */
318 : /* IBITBS = KEEPS CURRENT VALUE IF IBITB SO THAT LOOP */
319 : /* ENDING AT 166 DOESN'T HAVE TO START AT */
320 : /* IBITB = 0 EVERY TIME. */
321 : /* MISSLX(J) = MALLOW EXCEPT WHEN A GROUP IS ALL ONE VALUE (AND */
322 : /* LBIT(J) = 0) AND THAT VALUE IS MISSING. IN */
323 : /* THAT CASE, MISSLX(J) IS MISSP OR MISSS. THIS */
324 : /* GETS INSERTED INTO JMIN(J) LATER AS THE */
325 : /* MISSING INDICATOR; IT CAN'T BE PUT IN UNTIL */
326 : /* THE END, BECAUSE JMIN( ) IS USED TO CALCULATE */
327 : /* THE MAXIMUM NUMBER OF BITS (IBITS) NEEDED TO */
328 : /* PACK JMIN( ). */
329 : /* 1 2 3 4 5 6 7 X */
330 :
331 : /* NON SYSTEM SUBROUTINES CALLED */
332 : /* NONE */
333 :
334 :
335 :
336 : /* MISSLX( ) was AN AUTOMATIC ARRAY. */
337 0 : misslx = (integer *)calloc(*ndg,sizeof(integer));
338 :
339 :
340 : /* Parameter adjustments */
341 0 : --ic;
342 0 : --nov;
343 0 : --lbit;
344 0 : --jmax;
345 0 : --jmin;
346 :
347 : /* Function Body */
348 :
349 0 : *ier = 0;
350 0 : iersav = 0;
351 : /* CALL TIMPR(KFILDO,KFILDO,'START PACK_GP ') */
352 0 : *(unsigned char *)cfeed = (char) ifeed;
353 :
354 0 : ired = 0;
355 : /* IRED IS A FLAG. WHEN ZERO, REDUCE WILL BE CALLED. */
356 : /* IF REDUCE ABORTS, IRED = 1 AND IS NOT CALLED. IN */
357 : /* THIS CASE PACK_GP EXECUTES AGAIN EXCEPT FOR REDUCE. */
358 :
359 0 : if (*inc <= 0) {
360 0 : iersav = 717;
361 : /* WRITE(KFILDO,101)INC */
362 : /* 101 FORMAT(/' ****INC ='I8,' NOT CORRECT IN PACK_GP. 1 IS USED.') */
363 : }
364 :
365 : /* THERE WILL BE A RESTART OF PACK_GP IF SUBROUTINE REDUCE */
366 : /* ABORTS. THIS SHOULD NOT HAPPEN, BUT IF IT DOES, PACK_GP */
367 : /* WILL COMPLETE WITHOUT SUBROUTINE REDUCE. A NON FATAL */
368 : /* DIAGNOSTIC RETURN IS PROVIDED. */
369 :
370 : L102:
371 : /*kinc = max(*inc,1);*/
372 0 : kinc = (*inc > 1) ? *inc : 1;
373 0 : lminpk = *minpk;
374 :
375 : /* CALCULATE THE POWERS OF 2 THE FIRST TIME ENTERED. */
376 :
377 0 : if (ifirst == 0) {
378 0 : ifirst = 1;
379 0 : ibxx2[0] = 1;
380 :
381 0 : for (j = 1; j <= 30; ++j) {
382 0 : ibxx2[j] = ibxx2[j - 1] << 1;
383 : /* L104: */
384 : }
385 :
386 : }
387 :
388 : /* THERE WILL BE A RESTART AT 105 IS NDG IS NOT LARGE ENOUGH. */
389 : /* A NON FATAL DIAGNOSTIC RETURN IS PROVIDED. */
390 :
391 : L105:
392 0 : kstart = 1;
393 0 : ktotal = 0;
394 0 : *lx = 0;
395 0 : adda = FALSE_;
396 0 : lmiss = 0;
397 0 : if (*is523 == 1) {
398 0 : lmiss = 1;
399 : }
400 0 : if (*is523 == 2) {
401 0 : lmiss = 2;
402 : }
403 :
404 : /* ************************************* */
405 :
406 : /* THIS SECTION COMPUTES STATISTICS FOR GROUP A. GROUP A IS */
407 : /* A GROUP OF SIZE LMINPK. */
408 :
409 : /* ************************************* */
410 :
411 0 : ibita = 0;
412 0 : mina = mallow;
413 0 : maxa = -mallow;
414 0 : minak = mallow;
415 0 : maxak = -mallow;
416 :
417 : /* FIND THE MIN AND MAX OF GROUP A. THIS WILL INITIALLY BE OF */
418 : /* SIZE LMINPK (IF THERE ARE STILL LMINPK VALUES IN IC( )), BUT */
419 : /* WILL INCREASE IN SIZE IN INCREMENTS OF INC UNTIL A NEW */
420 : /* GROUP IS STARTED. THE DEFINITION OF GROUP A IS DONE HERE */
421 : /* ONLY ONCE (UPON INITIAL ENTRY), BECAUSE A GROUP B CAN ALWAYS */
422 : /* BECOME A NEW GROUP A AFTER A IS PACKED, EXCEPT IF LMINPK */
423 : /* HAS TO BE INCREASED BECAUSE NDG IS TOO SMALL. THEREFORE, */
424 : /* THE SEPARATE LOOPS FOR MISSING AND NON-MISSING HERE BUYS */
425 : /* ALMOST NOTHING. */
426 :
427 : /* Computing MIN */
428 0 : i__1 = kstart + lminpk - 1;
429 : /*nenda = min(i__1,*nxy);*/
430 0 : nenda = (i__1 < *nxy) ? i__1 : *nxy;
431 0 : if (*nxy - nenda <= lminpk / 2) {
432 0 : nenda = *nxy;
433 : }
434 : /* ABOVE STATEMENT GUARANTEES THE LAST GROUP IS GT LMINPK/2 BY */
435 : /* MAKING THE ACTUAL GROUP LARGER. IF A PROVISION LIKE THIS IS */
436 : /* NOT INCLUDED, THERE WILL MANY TIMES BE A VERY SMALL GROUP */
437 : /* AT THE END. USE SEPARATE LOOPS FOR MISSING AND NO MISSING */
438 : /* VALUES FOR EFFICIENCY. */
439 :
440 : /* DETERMINE WHETHER THERE IS A LONG STRING OF THE SAME VALUE */
441 : /* UNLESS NENDA = NXY. THIS MAY ALLOW A LARGE GROUP A TO */
442 : /* START WITH, AS WITH MISSING VALUES. SEPARATE LOOPS FOR */
443 : /* MISSING OPTIONS. THIS SECTION IS ONLY EXECUTED ONCE, */
444 : /* IN DETERMINING THE FIRST GROUP. IT HELPS FOR AN ARRAY */
445 : /* OF MOSTLY MISSING VALUES OR OF ONE VALUE, SUCH AS */
446 : /* RADAR OR PRECIP DATA. */
447 :
448 0 : if (nenda != *nxy && ic[kstart] == ic[kstart + 1]) {
449 : /* NO NEED TO EXECUTE IF FIRST TWO VALUES ARE NOT EQUAL. */
450 :
451 0 : if (*is523 == 0) {
452 : /* THIS LOOP IS FOR NO MISSING VALUES. */
453 :
454 0 : i__1 = *nxy;
455 0 : for (k = kstart + 1; k <= i__1; ++k) {
456 :
457 0 : if (ic[k] != ic[kstart]) {
458 : /* Computing MAX */
459 0 : i__2 = nenda, i__3 = k - 1;
460 : /*nenda = max(i__2,i__3);*/
461 0 : nenda = (i__2 > i__3) ? i__2 : i__3;
462 0 : goto L114;
463 : }
464 :
465 : /* L111: */
466 : }
467 :
468 0 : nenda = *nxy;
469 : /* FALL THROUGH THE LOOP MEANS ALL VALUES ARE THE SAME. */
470 :
471 0 : } else if (*is523 == 1) {
472 : /* THIS LOOP IS FOR PRIMARY MISSING VALUES ONLY. */
473 :
474 0 : i__1 = *nxy;
475 0 : for (k = kstart + 1; k <= i__1; ++k) {
476 :
477 0 : if (ic[k] != *missp) {
478 :
479 0 : if (ic[k] != ic[kstart]) {
480 : /* Computing MAX */
481 0 : i__2 = nenda, i__3 = k - 1;
482 : /*nenda = max(i__2,i__3);*/
483 0 : nenda = (i__2 > i__3) ? i__2 : i__3;
484 0 : goto L114;
485 : }
486 :
487 : }
488 :
489 : /* L112: */
490 : }
491 :
492 0 : nenda = *nxy;
493 : /* FALL THROUGH THE LOOP MEANS ALL VALUES ARE THE SAME. */
494 :
495 : } else {
496 : /* THIS LOOP IS FOR PRIMARY AND SECONDARY MISSING VALUES. */
497 :
498 0 : i__1 = *nxy;
499 0 : for (k = kstart + 1; k <= i__1; ++k) {
500 :
501 0 : if (ic[k] != *missp && ic[k] != *misss) {
502 :
503 0 : if (ic[k] != ic[kstart]) {
504 : /* Computing MAX */
505 0 : i__2 = nenda, i__3 = k - 1;
506 : /*nenda = max(i__2,i__3);*/
507 0 : nenda = (i__2 > i__3) ? i__2 : i__3;
508 0 : goto L114;
509 : }
510 :
511 : }
512 :
513 : /* L113: */
514 : }
515 :
516 0 : nenda = *nxy;
517 : /* FALL THROUGH THE LOOP MEANS ALL VALUES ARE THE SAME. */
518 : }
519 :
520 : }
521 :
522 : L114:
523 0 : if (*is523 == 0) {
524 :
525 0 : i__1 = nenda;
526 0 : for (k = kstart; k <= i__1; ++k) {
527 0 : if (ic[k] < mina) {
528 0 : mina = ic[k];
529 0 : minak = k;
530 : }
531 0 : if (ic[k] > maxa) {
532 0 : maxa = ic[k];
533 0 : maxak = k;
534 : }
535 : /* L115: */
536 : }
537 :
538 0 : } else if (*is523 == 1) {
539 :
540 0 : i__1 = nenda;
541 0 : for (k = kstart; k <= i__1; ++k) {
542 0 : if (ic[k] == *missp) {
543 0 : goto L117;
544 : }
545 0 : if (ic[k] < mina) {
546 0 : mina = ic[k];
547 0 : minak = k;
548 : }
549 0 : if (ic[k] > maxa) {
550 0 : maxa = ic[k];
551 0 : maxak = k;
552 : }
553 : L117:
554 : ;
555 : }
556 :
557 : } else {
558 :
559 0 : i__1 = nenda;
560 0 : for (k = kstart; k <= i__1; ++k) {
561 0 : if (ic[k] == *missp || ic[k] == *misss) {
562 : goto L120;
563 : }
564 0 : if (ic[k] < mina) {
565 0 : mina = ic[k];
566 0 : minak = k;
567 : }
568 0 : if (ic[k] > maxa) {
569 0 : maxa = ic[k];
570 0 : maxak = k;
571 : }
572 : L120:
573 : ;
574 : }
575 :
576 : }
577 :
578 0 : kounta = nenda - kstart + 1;
579 :
580 : /* INCREMENT KTOTAL AND FIND THE BITS NEEDED TO PACK THE A GROUP. */
581 :
582 0 : ktotal += kounta;
583 0 : mislla = 0;
584 0 : if (mina != mallow) {
585 0 : goto L125;
586 : }
587 : /* ALL MISSING VALUES MUST BE ACCOMMODATED. */
588 0 : mina = 0;
589 0 : maxa = 0;
590 0 : mislla = 1;
591 0 : ibitb = 0;
592 0 : if (*is523 != 2) {
593 0 : goto L130;
594 : }
595 : /* WHEN ALL VALUES ARE MISSING AND THERE ARE NO */
596 : /* SECONDARY MISSING VALUES, IBITA = 0. */
597 : /* OTHERWISE, IBITA MUST BE CALCULATED. */
598 :
599 : L125:
600 0 : itest = maxa - mina + lmiss;
601 :
602 0 : for (ibita = 0; ibita <= 30; ++ibita) {
603 0 : if (itest < ibxx2[ibita]) {
604 0 : goto L130;
605 : }
606 : /* *** THIS TEST IS THE SAME AS: */
607 : /* *** IF(MAXA-MINA.LT.IBXX2(IBITA)-LMISS)GO TO 130 */
608 : /* L126: */
609 : }
610 :
611 : /* WRITE(KFILDO,127)MAXA,MINA */
612 : /* 127 FORMAT(' ****ERROR IN PACK_GP. VALUE WILL NOT PACK IN 30 BITS.', */
613 : /* 1 ' MAXA ='I13,' MINA ='I13,'. ERROR AT 127.') */
614 0 : *ier = 706;
615 0 : goto L900;
616 :
617 : L130:
618 :
619 : /* ***D WRITE(KFILDO,131)KOUNTA,KTOTAL,MINA,MAXA,IBITA,MISLLA */
620 : /* ***D131 FORMAT(' AT 130, KOUNTA ='I8,' KTOTAL ='I8,' MINA ='I8, */
621 : /* ***D 1 ' MAXA ='I8,' IBITA ='I3,' MISLLA ='I3) */
622 :
623 : L133:
624 0 : if (ktotal >= *nxy) {
625 0 : goto L200;
626 : }
627 :
628 : /* ************************************* */
629 :
630 : /* THIS SECTION COMPUTES STATISTICS FOR GROUP B. GROUP B IS A */
631 : /* GROUP OF SIZE LMINPK IMMEDIATELY FOLLOWING GROUP A. */
632 :
633 : /* ************************************* */
634 :
635 : L140:
636 0 : minb = mallow;
637 0 : maxb = -mallow;
638 0 : minbk = mallow;
639 0 : maxbk = -mallow;
640 0 : ibitbs = 0;
641 0 : mstart = ktotal + 1;
642 :
643 : /* DETERMINE WHETHER THERE IS A LONG STRING OF THE SAME VALUE. */
644 : /* THIS WORKS WHEN THERE ARE NO MISSING VALUES. */
645 :
646 0 : nendb = 1;
647 :
648 0 : if (mstart < *nxy) {
649 :
650 0 : if (*is523 == 0) {
651 : /* THIS LOOP IS FOR NO MISSING VALUES. */
652 :
653 0 : i__1 = *nxy;
654 0 : for (k = mstart + 1; k <= i__1; ++k) {
655 :
656 0 : if (ic[k] != ic[mstart]) {
657 0 : nendb = k - 1;
658 0 : goto L150;
659 : }
660 :
661 : /* L145: */
662 : }
663 :
664 0 : nendb = *nxy;
665 : /* FALL THROUGH THE LOOP MEANS ALL REMAINING VALUES */
666 : /* ARE THE SAME. */
667 : }
668 :
669 : }
670 :
671 : L150:
672 : /* Computing MAX */
673 : /* Computing MIN */
674 0 : i__3 = ktotal + lminpk;
675 : /*i__1 = nendb, i__2 = min(i__3,*nxy);*/
676 0 : i__1 = nendb, i__2 = (i__3 < *nxy) ? i__3 : *nxy;
677 : /*nendb = max(i__1,i__2);*/
678 0 : nendb = (i__1 > i__2) ? i__1 : i__2;
679 : /* **** 150 NENDB=MIN(KTOTAL+LMINPK,NXY) */
680 :
681 0 : if (*nxy - nendb <= lminpk / 2) {
682 0 : nendb = *nxy;
683 : }
684 : /* ABOVE STATEMENT GUARANTEES THE LAST GROUP IS GT LMINPK/2 BY */
685 : /* MAKING THE ACTUAL GROUP LARGER. IF A PROVISION LIKE THIS IS */
686 : /* NOT INCLUDED, THERE WILL MANY TIMES BE A VERY SMALL GROUP */
687 : /* AT THE END. USE SEPARATE LOOPS FOR MISSING AND NO MISSING */
688 :
689 : /* USE SEPARATE LOOPS FOR MISSING AND NO MISSING VALUES */
690 : /* FOR EFFICIENCY. */
691 :
692 0 : if (*is523 == 0) {
693 :
694 0 : i__1 = nendb;
695 0 : for (k = mstart; k <= i__1; ++k) {
696 0 : if (ic[k] <= minb) {
697 0 : minb = ic[k];
698 : /* NOTE LE, NOT LT. LT COULD BE USED BUT THEN A */
699 : /* RECOMPUTE OVER THE WHOLE GROUP WOULD BE NEEDED */
700 : /* MORE OFTEN. SAME REASONING FOR GE AND OTHER */
701 : /* LOOPS BELOW. */
702 0 : minbk = k;
703 : }
704 0 : if (ic[k] >= maxb) {
705 0 : maxb = ic[k];
706 0 : maxbk = k;
707 : }
708 : /* L155: */
709 : }
710 :
711 0 : } else if (*is523 == 1) {
712 :
713 0 : i__1 = nendb;
714 0 : for (k = mstart; k <= i__1; ++k) {
715 0 : if (ic[k] == *missp) {
716 0 : goto L157;
717 : }
718 0 : if (ic[k] <= minb) {
719 0 : minb = ic[k];
720 0 : minbk = k;
721 : }
722 0 : if (ic[k] >= maxb) {
723 0 : maxb = ic[k];
724 0 : maxbk = k;
725 : }
726 : L157:
727 : ;
728 : }
729 :
730 : } else {
731 :
732 0 : i__1 = nendb;
733 0 : for (k = mstart; k <= i__1; ++k) {
734 0 : if (ic[k] == *missp || ic[k] == *misss) {
735 : goto L160;
736 : }
737 0 : if (ic[k] <= minb) {
738 0 : minb = ic[k];
739 0 : minbk = k;
740 : }
741 0 : if (ic[k] >= maxb) {
742 0 : maxb = ic[k];
743 0 : maxbk = k;
744 : }
745 : L160:
746 : ;
747 : }
748 :
749 : }
750 :
751 0 : kountb = nendb - ktotal;
752 0 : misllb = 0;
753 0 : if (minb != mallow) {
754 0 : goto L165;
755 : }
756 : /* ALL MISSING VALUES MUST BE ACCOMMODATED. */
757 0 : minb = 0;
758 0 : maxb = 0;
759 0 : misllb = 1;
760 0 : ibitb = 0;
761 :
762 0 : if (*is523 != 2) {
763 0 : goto L170;
764 : }
765 : /* WHEN ALL VALUES ARE MISSING AND THERE ARE NO SECONDARY */
766 : /* MISSING VALUES, IBITB = 0. OTHERWISE, IBITB MUST BE */
767 : /* CALCULATED. */
768 :
769 : L165:
770 0 : for (ibitb = ibitbs; ibitb <= 30; ++ibitb) {
771 0 : if (maxb - minb < ibxx2[ibitb] - lmiss) {
772 0 : goto L170;
773 : }
774 : /* L166: */
775 : }
776 :
777 : /* WRITE(KFILDO,167)MAXB,MINB */
778 : /* 167 FORMAT(' ****ERROR IN PACK_GP. VALUE WILL NOT PACK IN 30 BITS.', */
779 : /* 1 ' MAXB ='I13,' MINB ='I13,'. ERROR AT 167.') */
780 0 : *ier = 706;
781 0 : goto L900;
782 :
783 : /* COMPARE THE BITS NEEDED TO PACK GROUP B WITH THOSE NEEDED */
784 : /* TO PACK GROUP A. IF IBITB GE IBITA, TRY TO ADD TO GROUP A. */
785 : /* IF NOT, TRY TO ADD A'S POINTS TO B, UNLESS ADDITION TO A */
786 : /* HAS BEEN DONE. THIS LATTER IS CONTROLLED WITH ADDA. */
787 :
788 : L170:
789 :
790 : /* ***D WRITE(KFILDO,171)KOUNTA,KTOTAL,MINA,MAXA,IBITA,MISLLA, */
791 : /* ***D 1 MINB,MAXB,IBITB,MISLLB */
792 : /* ***D171 FORMAT(' AT 171, KOUNTA ='I8,' KTOTAL ='I8,' MINA ='I8, */
793 : /* ***D 1 ' MAXA ='I8,' IBITA ='I3,' MISLLA ='I3, */
794 : /* ***D 2 ' MINB ='I8,' MAXB ='I8,' IBITB ='I3,' MISLLB ='I3) */
795 :
796 0 : if (ibitb >= ibita) {
797 0 : goto L180;
798 : }
799 0 : if (adda) {
800 0 : goto L200;
801 : }
802 :
803 : /* ************************************* */
804 :
805 : /* GROUP B REQUIRES LESS BITS THAN GROUP A. PUT AS MANY OF A'S */
806 : /* POINTS INTO B AS POSSIBLE WITHOUT EXCEEDING THE NUMBER OF */
807 : /* BITS NECESSARY TO PACK GROUP B. */
808 :
809 : /* ************************************* */
810 :
811 0 : kounts = kounta;
812 : /* KOUNTA REFERS TO THE PRESENT GROUP A. */
813 0 : mintst = minb;
814 0 : maxtst = maxb;
815 0 : mintstk = minbk;
816 0 : maxtstk = maxbk;
817 :
818 : /* USE SEPARATE LOOPS FOR MISSING AND NO MISSING VALUES */
819 : /* FOR EFFICIENCY. */
820 :
821 0 : if (*is523 == 0) {
822 :
823 0 : i__1 = kstart;
824 0 : for (k = ktotal; k >= i__1; --k) {
825 : /* START WITH THE END OF THE GROUP AND WORK BACKWARDS. */
826 0 : if (ic[k] < minb) {
827 0 : mintst = ic[k];
828 0 : mintstk = k;
829 0 : } else if (ic[k] > maxb) {
830 0 : maxtst = ic[k];
831 0 : maxtstk = k;
832 : }
833 0 : if (maxtst - mintst >= ibxx2[ibitb]) {
834 0 : goto L174;
835 : }
836 : /* NOTE THAT FOR THIS LOOP, LMISS = 0. */
837 0 : minb = mintst;
838 0 : maxb = maxtst;
839 0 : minbk = mintstk;
840 0 : maxbk = maxtstk;
841 0 : --kounta;
842 : /* THERE IS ONE LESS POINT NOW IN A. */
843 : /* L1715: */
844 : }
845 :
846 0 : } else if (*is523 == 1) {
847 :
848 0 : i__1 = kstart;
849 0 : for (k = ktotal; k >= i__1; --k) {
850 : /* START WITH THE END OF THE GROUP AND WORK BACKWARDS. */
851 0 : if (ic[k] == *missp) {
852 0 : goto L1718;
853 : }
854 0 : if (ic[k] < minb) {
855 0 : mintst = ic[k];
856 0 : mintstk = k;
857 0 : } else if (ic[k] > maxb) {
858 0 : maxtst = ic[k];
859 0 : maxtstk = k;
860 : }
861 0 : if (maxtst - mintst >= ibxx2[ibitb] - lmiss) {
862 0 : goto L174;
863 : }
864 : /* FOR THIS LOOP, LMISS = 1. */
865 0 : minb = mintst;
866 0 : maxb = maxtst;
867 0 : minbk = mintstk;
868 0 : maxbk = maxtstk;
869 0 : misllb = 0;
870 : /* WHEN THE POINT IS NON MISSING, MISLLB SET = 0. */
871 : L1718:
872 0 : --kounta;
873 : /* THERE IS ONE LESS POINT NOW IN A. */
874 : /* L1719: */
875 : }
876 :
877 : } else {
878 :
879 0 : i__1 = kstart;
880 0 : for (k = ktotal; k >= i__1; --k) {
881 : /* START WITH THE END OF THE GROUP AND WORK BACKWARDS. */
882 0 : if (ic[k] == *missp || ic[k] == *misss) {
883 : goto L1729;
884 : }
885 0 : if (ic[k] < minb) {
886 0 : mintst = ic[k];
887 0 : mintstk = k;
888 0 : } else if (ic[k] > maxb) {
889 0 : maxtst = ic[k];
890 0 : maxtstk = k;
891 : }
892 0 : if (maxtst - mintst >= ibxx2[ibitb] - lmiss) {
893 0 : goto L174;
894 : }
895 : /* FOR THIS LOOP, LMISS = 2. */
896 0 : minb = mintst;
897 0 : maxb = maxtst;
898 0 : minbk = mintstk;
899 0 : maxbk = maxtstk;
900 0 : misllb = 0;
901 : /* WHEN THE POINT IS NON MISSING, MISLLB SET = 0. */
902 : L1729:
903 0 : --kounta;
904 : /* THERE IS ONE LESS POINT NOW IN A. */
905 : /* L173: */
906 : }
907 :
908 : }
909 :
910 : /* AT THIS POINT, KOUNTA CONTAINS THE NUMBER OF POINTS TO CLOSE */
911 : /* OUT GROUP A WITH. GROUP B NOW STARTS WITH KSTART+KOUNTA AND */
912 : /* ENDS WITH NENDB. MINB AND MAXB HAVE BEEN ADJUSTED AS */
913 : /* NECESSARY TO REFLECT GROUP B (EVEN THOUGH THE NUMBER OF BITS */
914 : /* NEEDED TO PACK GROUP B HAVE NOT INCREASED, THE END POINTS */
915 : /* OF THE RANGE MAY HAVE). */
916 :
917 : L174:
918 0 : if (kounta == kounts) {
919 0 : goto L200;
920 : }
921 : /* ON TRANSFER, GROUP A WAS NOT CHANGED. CLOSE IT OUT. */
922 :
923 : /* ONE OR MORE POINTS WERE TAKEN OUT OF A. RANGE AND IBITA */
924 : /* MAY HAVE TO BE RECOMPUTED; IBITA COULD BE LESS THAN */
925 : /* ORIGINALLY COMPUTED. IN FACT, GROUP A CAN NOW CONTAIN */
926 : /* ONLY ONE POINT AND BE PACKED WITH ZERO BITS */
927 : /* (UNLESS MISSS NE 0). */
928 :
929 0 : nouta = kounts - kounta;
930 0 : ktotal -= nouta;
931 0 : kountb += nouta;
932 0 : if (nenda - nouta > minak && nenda - nouta > maxak) {
933 0 : goto L200;
934 : }
935 : /* WHEN THE ABOVE TEST IS MET, THE MIN AND MAX OF THE */
936 : /* CURRENT GROUP A WERE WITHIN THE OLD GROUP A, SO THE */
937 : /* RANGE AND IBITA DO NOT NEED TO BE RECOMPUTED. */
938 : /* NOTE THAT MINAK AND MAXAK ARE NO LONGER NEEDED. */
939 0 : ibita = 0;
940 0 : mina = mallow;
941 0 : maxa = -mallow;
942 :
943 : /* USE SEPARATE LOOPS FOR MISSING AND NO MISSING VALUES */
944 : /* FOR EFFICIENCY. */
945 :
946 0 : if (*is523 == 0) {
947 :
948 0 : i__1 = nenda - nouta;
949 0 : for (k = kstart; k <= i__1; ++k) {
950 0 : if (ic[k] < mina) {
951 0 : mina = ic[k];
952 : }
953 0 : if (ic[k] > maxa) {
954 0 : maxa = ic[k];
955 : }
956 : /* L1742: */
957 : }
958 :
959 0 : } else if (*is523 == 1) {
960 :
961 0 : i__1 = nenda - nouta;
962 0 : for (k = kstart; k <= i__1; ++k) {
963 0 : if (ic[k] == *missp) {
964 0 : goto L1744;
965 : }
966 0 : if (ic[k] < mina) {
967 0 : mina = ic[k];
968 : }
969 0 : if (ic[k] > maxa) {
970 0 : maxa = ic[k];
971 : }
972 : L1744:
973 : ;
974 : }
975 :
976 : } else {
977 :
978 0 : i__1 = nenda - nouta;
979 0 : for (k = kstart; k <= i__1; ++k) {
980 0 : if (ic[k] == *missp || ic[k] == *misss) {
981 : goto L175;
982 : }
983 0 : if (ic[k] < mina) {
984 0 : mina = ic[k];
985 : }
986 0 : if (ic[k] > maxa) {
987 0 : maxa = ic[k];
988 : }
989 : L175:
990 : ;
991 : }
992 :
993 : }
994 :
995 0 : mislla = 0;
996 0 : if (mina != mallow) {
997 0 : goto L1750;
998 : }
999 : /* ALL MISSING VALUES MUST BE ACCOMMODATED. */
1000 0 : mina = 0;
1001 0 : maxa = 0;
1002 0 : mislla = 1;
1003 0 : if (*is523 != 2) {
1004 0 : goto L177;
1005 : }
1006 : /* WHEN ALL VALUES ARE MISSING AND THERE ARE NO SECONDARY */
1007 : /* MISSING VALUES IBITA = 0 AS ORIGINALLY SET. OTHERWISE, */
1008 : /* IBITA MUST BE CALCULATED. */
1009 :
1010 : L1750:
1011 0 : itest = maxa - mina + lmiss;
1012 :
1013 0 : for (ibita = 0; ibita <= 30; ++ibita) {
1014 0 : if (itest < ibxx2[ibita]) {
1015 0 : goto L177;
1016 : }
1017 : /* *** THIS TEST IS THE SAME AS: */
1018 : /* *** IF(MAXA-MINA.LT.IBXX2(IBITA)-LMISS)GO TO 177 */
1019 : /* L176: */
1020 : }
1021 :
1022 : /* WRITE(KFILDO,1760)MAXA,MINA */
1023 : /* 1760 FORMAT(' ****ERROR IN PACK_GP. VALUE WILL NOT PACK IN 30 BITS.', */
1024 : /* 1 ' MAXA ='I13,' MINA ='I13,'. ERROR AT 1760.') */
1025 0 : *ier = 706;
1026 0 : goto L900;
1027 :
1028 : L177:
1029 0 : goto L200;
1030 :
1031 : /* ************************************* */
1032 :
1033 : /* AT THIS POINT, GROUP B REQUIRES AS MANY BITS TO PACK AS GROUPA. */
1034 : /* THEREFORE, TRY TO ADD INC POINTS TO GROUP A WITHOUT INCREASING */
1035 : /* IBITA. THIS AUGMENTED GROUP IS CALLED GROUP C. */
1036 :
1037 : /* ************************************* */
1038 :
1039 : L180:
1040 0 : if (mislla == 1) {
1041 0 : minc = mallow;
1042 0 : minck = mallow;
1043 0 : maxc = -mallow;
1044 0 : maxck = -mallow;
1045 : } else {
1046 0 : minc = mina;
1047 0 : maxc = maxa;
1048 0 : minck = minak;
1049 0 : maxck = minak;
1050 : }
1051 :
1052 0 : nount = 0;
1053 0 : if (*nxy - (ktotal + kinc) <= lminpk / 2) {
1054 0 : kinc = *nxy - ktotal;
1055 : }
1056 : /* ABOVE STATEMENT CONSTRAINS THE LAST GROUP TO BE NOT LESS THAN */
1057 : /* LMINPK/2 IN SIZE. IF A PROVISION LIKE THIS IS NOT INCLUDED, */
1058 : /* THERE WILL MANY TIMES BE A VERY SMALL GROUP AT THE END. */
1059 :
1060 : /* USE SEPARATE LOOPS FOR MISSING AND NO MISSING VALUES */
1061 : /* FOR EFFICIENCY. SINCE KINC IS USUALLY 1, USING SEPARATE */
1062 : /* LOOPS HERE DOESN'T BUY MUCH. A MISSING VALUE WILL ALWAYS */
1063 : /* TRANSFER BACK TO GROUP A. */
1064 :
1065 0 : if (*is523 == 0) {
1066 :
1067 : /* Computing MIN */
1068 0 : i__2 = ktotal + kinc;
1069 : /*i__1 = min(i__2,*nxy);*/
1070 0 : i__1 = (i__2 < *nxy) ? i__2 : *nxy;
1071 0 : for (k = ktotal + 1; k <= i__1; ++k) {
1072 0 : if (ic[k] < minc) {
1073 0 : minc = ic[k];
1074 0 : minck = k;
1075 : }
1076 0 : if (ic[k] > maxc) {
1077 0 : maxc = ic[k];
1078 0 : maxck = k;
1079 : }
1080 0 : ++nount;
1081 : /* L185: */
1082 : }
1083 :
1084 0 : } else if (*is523 == 1) {
1085 :
1086 : /* Computing MIN */
1087 0 : i__2 = ktotal + kinc;
1088 : /*i__1 = min(i__2,*nxy);*/
1089 0 : i__1 = (i__2 < *nxy) ? i__2 : *nxy;
1090 0 : for (k = ktotal + 1; k <= i__1; ++k) {
1091 0 : if (ic[k] == *missp) {
1092 0 : goto L186;
1093 : }
1094 0 : if (ic[k] < minc) {
1095 0 : minc = ic[k];
1096 0 : minck = k;
1097 : }
1098 0 : if (ic[k] > maxc) {
1099 0 : maxc = ic[k];
1100 0 : maxck = k;
1101 : }
1102 : L186:
1103 0 : ++nount;
1104 : /* L187: */
1105 : }
1106 :
1107 : } else {
1108 :
1109 : /* Computing MIN */
1110 0 : i__2 = ktotal + kinc;
1111 : /*i__1 = min(i__2,*nxy);*/
1112 0 : i__1 = (i__2 < *nxy) ? i__2 : *nxy;
1113 0 : for (k = ktotal + 1; k <= i__1; ++k) {
1114 0 : if (ic[k] == *missp || ic[k] == *misss) {
1115 : goto L189;
1116 : }
1117 0 : if (ic[k] < minc) {
1118 0 : minc = ic[k];
1119 0 : minck = k;
1120 : }
1121 0 : if (ic[k] > maxc) {
1122 0 : maxc = ic[k];
1123 0 : maxck = k;
1124 : }
1125 : L189:
1126 0 : ++nount;
1127 : /* L190: */
1128 : }
1129 :
1130 : }
1131 :
1132 : /* ***D WRITE(KFILDO,191)KOUNTA,KTOTAL,MINA,MAXA,IBITA,MISLLA, */
1133 : /* ***D 1 MINC,MAXC,NOUNT,IC(KTOTAL),IC(KTOTAL+1) */
1134 : /* ***D191 FORMAT(' AT 191, KOUNTA ='I8,' KTOTAL ='I8,' MINA ='I8, */
1135 : /* ***D 1 ' MAXA ='I8,' IBITA ='I3,' MISLLA ='I3, */
1136 : /* ***D 2 ' MINC ='I8,' MAXC ='I8, */
1137 : /* ***D 3 ' NOUNT ='I5,' IC(KTOTAL) ='I9,' IC(KTOTAL+1) =',I9) */
1138 :
1139 : /* IF THE NUMBER OF BITS NEEDED FOR GROUP C IS GT IBITA, */
1140 : /* THEN THIS GROUP A IS A GROUP TO PACK. */
1141 :
1142 0 : if (minc == mallow) {
1143 0 : minc = mina;
1144 0 : maxc = maxa;
1145 0 : minck = minak;
1146 0 : maxck = maxak;
1147 0 : misllc = 1;
1148 0 : goto L195;
1149 : /* WHEN THE NEW VALUE(S) ARE MISSING, THEY CAN ALWAYS */
1150 : /* BE ADDED. */
1151 :
1152 : } else {
1153 0 : misllc = 0;
1154 : }
1155 :
1156 0 : if (maxc - minc >= ibxx2[ibita] - lmiss) {
1157 0 : goto L200;
1158 : }
1159 :
1160 : /* THE BITS NECESSARY FOR GROUP C HAS NOT INCREASED FROM THE */
1161 : /* BITS NECESSARY FOR GROUP A. ADD THIS POINT(S) TO GROUP A. */
1162 : /* COMPUTE THE NEXT GROUP B, ETC., UNLESS ALL POINTS HAVE BEEN */
1163 : /* USED. */
1164 :
1165 : L195:
1166 0 : ktotal += nount;
1167 0 : kounta += nount;
1168 0 : mina = minc;
1169 0 : maxa = maxc;
1170 0 : minak = minck;
1171 0 : maxak = maxck;
1172 0 : mislla = misllc;
1173 0 : adda = TRUE_;
1174 0 : if (ktotal >= *nxy) {
1175 0 : goto L200;
1176 : }
1177 :
1178 0 : if (minbk > ktotal && maxbk > ktotal) {
1179 0 : mstart = nendb + 1;
1180 : /* THE MAX AND MIN OF GROUP B WERE NOT FROM THE POINTS */
1181 : /* REMOVED, SO THE WHOLE GROUP DOES NOT HAVE TO BE LOOKED */
1182 : /* AT TO DETERMINE THE NEW MAX AND MIN. RATHER START */
1183 : /* JUST BEYOND THE OLD NENDB. */
1184 0 : ibitbs = ibitb;
1185 0 : nendb = 1;
1186 0 : goto L150;
1187 : } else {
1188 : goto L140;
1189 : }
1190 :
1191 : /* ************************************* */
1192 :
1193 : /* GROUP A IS TO BE PACKED. STORE VALUES IN JMIN( ), JMAX( ), */
1194 : /* LBIT( ), AND NOV( ). */
1195 :
1196 : /* ************************************* */
1197 :
1198 : L200:
1199 0 : ++(*lx);
1200 0 : if (*lx <= *ndg) {
1201 0 : goto L205;
1202 : }
1203 0 : lminpk += lminpk / 2;
1204 : /* WRITE(KFILDO,201)NDG,LMINPK,LX */
1205 : /* 201 FORMAT(' ****NDG ='I5,' NOT LARGE ENOUGH.', */
1206 : /* 1 ' LMINPK IS INCREASED TO 'I3,' FOR THIS FIELD.'/ */
1207 : /* 2 ' LX = 'I10) */
1208 0 : iersav = 716;
1209 0 : goto L105;
1210 :
1211 : L205:
1212 0 : jmin[*lx] = mina;
1213 0 : jmax[*lx] = maxa;
1214 0 : lbit[*lx] = ibita;
1215 0 : nov[*lx] = kounta;
1216 0 : kstart = ktotal + 1;
1217 :
1218 0 : if (mislla == 0) {
1219 0 : misslx[*lx - 1] = mallow;
1220 : } else {
1221 0 : misslx[*lx - 1] = ic[ktotal];
1222 : /* IC(KTOTAL) WAS THE LAST VALUE PROCESSED. IF MISLLA NE 0, */
1223 : /* THIS MUST BE THE MISSING VALUE FOR THIS GROUP. */
1224 : }
1225 :
1226 : /* ***D WRITE(KFILDO,206)MISLLA,IC(KTOTAL),KTOTAL,LX,JMIN(LX),JMAX(LX), */
1227 : /* ***D 1 LBIT(LX),NOV(LX),MISSLX(LX) */
1228 : /* ***D206 FORMAT(' AT 206, MISLLA ='I2,' IC(KTOTAL) ='I5,' KTOTAL ='I8, */
1229 : /* ***D 1 ' LX ='I6,' JMIN(LX) ='I8,' JMAX(LX) ='I8, */
1230 : /* ***D 2 ' LBIT(LX) ='I5,' NOV(LX) ='I8,' MISSLX(LX) =',I7) */
1231 :
1232 0 : if (ktotal >= *nxy) {
1233 0 : goto L209;
1234 : }
1235 :
1236 : /* THE NEW GROUP A WILL BE THE PREVIOUS GROUP B. SET LIMITS, ETC. */
1237 :
1238 0 : ibita = ibitb;
1239 0 : mina = minb;
1240 0 : maxa = maxb;
1241 0 : minak = minbk;
1242 0 : maxak = maxbk;
1243 0 : mislla = misllb;
1244 0 : nenda = nendb;
1245 0 : kounta = kountb;
1246 0 : ktotal += kounta;
1247 0 : adda = FALSE_;
1248 0 : goto L133;
1249 :
1250 : /* ************************************* */
1251 :
1252 : /* CALCULATE IBIT, THE NUMBER OF BITS NEEDED TO HOLD THE GROUP */
1253 : /* MINIMUM VALUES. */
1254 :
1255 : /* ************************************* */
1256 :
1257 : L209:
1258 0 : *ibit = 0;
1259 :
1260 0 : i__1 = *lx;
1261 0 : for (l = 1; l <= i__1; ++l) {
1262 : L210:
1263 0 : if (jmin[l] < ibxx2[*ibit]) {
1264 0 : goto L220;
1265 : }
1266 0 : ++(*ibit);
1267 0 : goto L210;
1268 : L220:
1269 : ;
1270 : }
1271 :
1272 : /* INSERT THE VALUE IN JMIN( ) TO BE USED FOR ALL MISSING */
1273 : /* VALUES WHEN LBIT( ) = 0. WHEN SECONDARY MISSING */
1274 : /* VALUES CAN BE PRESENT, LBIT(L) WILL NOT = 0. */
1275 :
1276 0 : if (*is523 == 1) {
1277 :
1278 0 : i__1 = *lx;
1279 0 : for (l = 1; l <= i__1; ++l) {
1280 :
1281 0 : if (lbit[l] == 0) {
1282 :
1283 0 : if (misslx[l - 1] == *missp) {
1284 0 : jmin[l] = ibxx2[*ibit] - 1;
1285 : }
1286 :
1287 : }
1288 :
1289 : /* L226: */
1290 : }
1291 :
1292 : }
1293 :
1294 : /* ************************************* */
1295 :
1296 : /* CALCULATE JBIT, THE NUMBER OF BITS NEEDED TO HOLD THE BITS */
1297 : /* NEEDED TO PACK THE VALUES IN THE GROUPS. BUT FIND AND */
1298 : /* REMOVE THE REFERENCE VALUE FIRST. */
1299 :
1300 : /* ************************************* */
1301 :
1302 : /* WRITE(KFILDO,228)CFEED,LX */
1303 : /* 228 FORMAT(A1,/' *****************************************' */
1304 : /* 1 /' THE GROUP WIDTHS LBIT( ) FOR ',I8,' GROUPS' */
1305 : /* 2 /' *****************************************') */
1306 : /* WRITE(KFILDO,229) (LBIT(J),J=1,MIN(LX,100)) */
1307 : /* 229 FORMAT(/' '20I6) */
1308 :
1309 0 : *lbitref = lbit[1];
1310 :
1311 0 : i__1 = *lx;
1312 0 : for (k = 1; k <= i__1; ++k) {
1313 0 : if (lbit[k] < *lbitref) {
1314 0 : *lbitref = lbit[k];
1315 : }
1316 : /* L230: */
1317 : }
1318 :
1319 0 : if (*lbitref != 0) {
1320 :
1321 0 : i__1 = *lx;
1322 0 : for (k = 1; k <= i__1; ++k) {
1323 0 : lbit[k] -= *lbitref;
1324 : /* L240: */
1325 : }
1326 :
1327 : }
1328 :
1329 : /* WRITE(KFILDO,241)CFEED,LBITREF */
1330 : /* 241 FORMAT(A1,/' *****************************************' */
1331 : /* 1 /' THE GROUP WIDTHS LBIT( ) AFTER REMOVING REFERENCE ', */
1332 : /* 2 I8, */
1333 : /* 3 /' *****************************************') */
1334 : /* WRITE(KFILDO,242) (LBIT(J),J=1,MIN(LX,100)) */
1335 : /* 242 FORMAT(/' '20I6) */
1336 :
1337 0 : *jbit = 0;
1338 :
1339 0 : i__1 = *lx;
1340 0 : for (k = 1; k <= i__1; ++k) {
1341 : L310:
1342 0 : if (lbit[k] < ibxx2[*jbit]) {
1343 0 : goto L320;
1344 : }
1345 0 : ++(*jbit);
1346 0 : goto L310;
1347 : L320:
1348 : ;
1349 : }
1350 :
1351 : /* ************************************* */
1352 :
1353 : /* CALCULATE KBIT, THE NUMBER OF BITS NEEDED TO HOLD THE NUMBER */
1354 : /* OF VALUES IN THE GROUPS. BUT FIND AND REMOVE THE */
1355 : /* REFERENCE FIRST. */
1356 :
1357 : /* ************************************* */
1358 :
1359 : /* WRITE(KFILDO,321)CFEED,LX */
1360 : /* 321 FORMAT(A1,/' *****************************************' */
1361 : /* 1 /' THE GROUP SIZES NOV( ) FOR ',I8,' GROUPS' */
1362 : /* 2 /' *****************************************') */
1363 : /* WRITE(KFILDO,322) (NOV(J),J=1,MIN(LX,100)) */
1364 : /* 322 FORMAT(/' '20I6) */
1365 :
1366 0 : *novref = nov[1];
1367 :
1368 0 : i__1 = *lx;
1369 0 : for (k = 1; k <= i__1; ++k) {
1370 0 : if (nov[k] < *novref) {
1371 0 : *novref = nov[k];
1372 : }
1373 : /* L400: */
1374 : }
1375 :
1376 0 : if (*novref > 0) {
1377 :
1378 0 : i__1 = *lx;
1379 0 : for (k = 1; k <= i__1; ++k) {
1380 0 : nov[k] -= *novref;
1381 : /* L405: */
1382 : }
1383 :
1384 : }
1385 :
1386 : /* WRITE(KFILDO,406)CFEED,NOVREF */
1387 : /* 406 FORMAT(A1,/' *****************************************' */
1388 : /* 1 /' THE GROUP SIZES NOV( ) AFTER REMOVING REFERENCE ',I8, */
1389 : /* 2 /' *****************************************') */
1390 : /* WRITE(KFILDO,407) (NOV(J),J=1,MIN(LX,100)) */
1391 : /* 407 FORMAT(/' '20I6) */
1392 : /* WRITE(KFILDO,408)CFEED */
1393 : /* 408 FORMAT(A1,/' *****************************************' */
1394 : /* 1 /' THE GROUP REFERENCES JMIN( )' */
1395 : /* 2 /' *****************************************') */
1396 : /* WRITE(KFILDO,409) (JMIN(J),J=1,MIN(LX,100)) */
1397 : /* 409 FORMAT(/' '20I6) */
1398 :
1399 0 : *kbit = 0;
1400 :
1401 0 : i__1 = *lx;
1402 0 : for (k = 1; k <= i__1; ++k) {
1403 : L410:
1404 0 : if (nov[k] < ibxx2[*kbit]) {
1405 0 : goto L420;
1406 : }
1407 0 : ++(*kbit);
1408 0 : goto L410;
1409 : L420:
1410 : ;
1411 : }
1412 :
1413 : /* DETERMINE WHETHER THE GROUP SIZES SHOULD BE REDUCED */
1414 : /* FOR SPACE EFFICIENCY. */
1415 :
1416 0 : if (ired == 0) {
1417 0 : reduce(kfildo, &jmin[1], &jmax[1], &lbit[1], &nov[1], lx, ndg, ibit,
1418 : jbit, kbit, novref, ibxx2, ier);
1419 :
1420 0 : if (*ier == 714 || *ier == 715) {
1421 : /* REDUCE HAS ABORTED. REEXECUTE PACK_GP WITHOUT REDUCE. */
1422 : /* PROVIDE FOR A NON FATAL RETURN FROM REDUCE. */
1423 0 : iersav = *ier;
1424 0 : ired = 1;
1425 0 : *ier = 0;
1426 0 : goto L102;
1427 : }
1428 :
1429 : }
1430 :
1431 0 : if ( misslx != 0 ) {
1432 0 : free(misslx);
1433 0 : misslx=0;
1434 : }
1435 : /* CALL TIMPR(KFILDO,KFILDO,'END PACK_GP ') */
1436 0 : if (iersav != 0) {
1437 0 : *ier = iersav;
1438 0 : return 0;
1439 : }
1440 :
1441 : /* 900 IF(IER.NE.0)RETURN1 */
1442 :
1443 : L900:
1444 0 : if ( misslx != 0 ) free(misslx);
1445 0 : return 0;
1446 : } /* pack_gp__ */
1447 :
|