1 /* reduce.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
17 #include "cpl_port.h"
18
19 typedef g2int integer;
20 typedef g2float real;
21
reduce(CPL_UNUSED integer * kfildo,integer * jmin,integer * jmax,integer * lbit,integer * nov,integer * lx,integer * ndg,integer * ibit,integer * jbit,integer * kbit,integer * novref,integer * ibxx2,integer * ier)22 /* Subroutine */ int reduce(CPL_UNUSED integer *kfildo, integer *jmin, integer *jmax,
23 integer *lbit, integer *nov, integer *lx, integer *ndg, integer *ibit,
24 integer *jbit, integer *kbit, integer *novref, integer *ibxx2,
25 integer *ier)
26 {
27 /* Initialized data */
28
29 static integer ifeed = 12;
30
31 /* System generated locals */
32 integer i__1, i__2;
33
34 /* Local variables */
35 static integer newboxtp, j, l, m, jj, lxn, left;
36 static real pimp;
37 static integer move, novl;
38 static char cfeed[1];
39 static integer /* nboxj[31], */ lxnkp, iorigb, ibxx2m1, movmin,
40 ntotbt[31], ntotpr, newboxt;
41 integer *newbox, *newboxp;
42
43
44 /* NOVEMBER 2001 GLAHN TDL GRIB2 */
45 /* MARCH 2002 GLAHN COMMENT IER = 715 */
46 /* MARCH 2002 GLAHN MODIFIED TO ACCOMMODATE LX=1 ON ENTRY */
47
48 /* PURPOSE */
49 /* DETERMINES WHETHER THE NUMBER OF GROUPS SHOULD BE */
50 /* INCREASED IN ORDER TO REDUCE THE SIZE OF THE LARGE */
51 /* GROUPS, AND TO MAKE THAT ADJUSTMENT. BY REDUCING THE */
52 /* SIZE OF THE LARGE GROUPS, LESS BITS MAY BE NECESSARY */
53 /* FOR PACKING THE GROUP SIZES AND ALL THE INFORMATION */
54 /* ABOUT THE GROUPS. */
55
56 /* THE REFERENCE FOR NOV( ) WAS REMOVED IN THE CALLING */
57 /* ROUTINE SO THAT KBIT COULD BE DETERMINED. THIS */
58 /* FURNISHES A STARTING POINT FOR THE ITERATIONS IN REDUCE. */
59 /* HOWEVER, THE REFERENCE MUST BE CONSIDERED. */
60
61 /* DATA SET USE */
62 /* KFILDO - UNIT NUMBER FOR OUTPUT (PRINT) FILE. (OUTPUT) */
63
64 /* VARIABLES IN CALL SEQUENCE */
65 /* KFILDO = UNIT NUMBER FOR OUTPUT (PRINT) FILE. (INPUT) */
66 /* JMIN(J) = THE MINIMUM OF EACH GROUP (J=1,LX). IT IS */
67 /* POSSIBLE AFTER SPLITTING THE GROUPS, JMIN( ) */
68 /* WILL NOT BE THE MINIMUM OF THE NEW GROUP. */
69 /* THIS DOESN'T MATTER; JMIN( ) IS REALLY THE */
70 /* GROUP REFERENCE AND DOESN'T HAVE TO BE THE */
71 /* SMALLEST VALUE. (INPUT/OUTPUT) */
72 /* JMAX(J) = THE MAXIMUM OF EACH GROUP (J=1,LX). */
73 /* (INPUT/OUTPUT) */
74 /* LBIT(J) = THE NUMBER OF BITS NECESSARY TO PACK EACH GROUP */
75 /* (J=1,LX). (INPUT/OUTPUT) */
76 /* NOV(J) = THE NUMBER OF VALUES IN EACH GROUP (J=1,LX). */
77 /* (INPUT/OUTPUT) */
78 /* LX = THE NUMBER OF GROUPS. THIS WILL BE INCREASED */
79 /* IF GROUPS ARE SPLIT. (INPUT/OUTPUT) */
80 /* NDG = THE DIMENSION OF JMIN( ), JMAX( ), LBIT( ), AND */
81 /* NOV( ). (INPUT) */
82 /* IBIT = THE NUMBER OF BITS NECESSARY TO PACK THE JMIN(J) */
83 /* VALUES, J=1,LX. (INPUT) */
84 /* JBIT = THE NUMBER OF BITS NECESSARY TO PACK THE LBIT(J) */
85 /* VALUES, J=1,LX. (INPUT) */
86 /* KBIT = THE NUMBER OF BITS NECESSARY TO PACK THE NOV(J) */
87 /* VALUES, J=1,LX. IF THE GROUPS ARE SPLIT, KBIT */
88 /* IS REDUCED. (INPUT/OUTPUT) */
89 /* NOVREF = REFERENCE VALUE FOR NOV( ). (INPUT) */
90 /* IBXX2(J) = 2**J (J=0,30). (INPUT) */
91 /* IER = ERROR RETURN. (OUTPUT) */
92 /* 0 = GOOD RETURN. */
93 /* 714 = PROBLEM IN ALGORITHM. REDUCE ABORTED. */
94 /* 715 = NGP NOT LARGE ENOUGH. REDUCE ABORTED. */
95 /* NTOTBT(J) = THE TOTAL BITS USED FOR THE PACKING BITS J */
96 /* (J=1,30). (INTERNAL) */
97 /* NBOXJ(J) = NEW BOXES NEEDED FOR THE PACKING BITS J */
98 /* (J=1,30). (INTERNAL) */
99 /* NEWBOX(L) = NUMBER OF NEW BOXES (GROUPS) FOR EACH ORIGINAL */
100 /* GROUP (L=1,LX) FOR THE CURRENT J. (AUTOMATIC) */
101 /* (INTERNAL) */
102 /* NEWBOXP(L) = SAME AS NEWBOX( ) BUT FOR THE PREVIOUS J. */
103 /* THIS ELIMINATES RECOMPUTATION. (AUTOMATIC) */
104 /* (INTERNAL) */
105 /* CFEED = CONTAINS THE CHARACTER REPRESENTATION */
106 /* OF A PRINTER FORM FEED. (CHARACTER) (INTERNAL) */
107 /* IFEED = CONTAINS THE INTEGER VALUE OF A PRINTER */
108 /* FORM FEED. (INTERNAL) */
109 /* IORIGB = THE ORIGINAL NUMBER OF BITS NECESSARY */
110 /* FOR THE GROUP VALUES. (INTERNAL) */
111 /* 1 2 3 4 5 6 7 X */
112
113 /* NON SYSTEM SUBROUTINES CALLED */
114 /* NONE */
115
116
117 /* NEWBOX( ) AND NEWBOXP( ) were AUTOMATIC ARRAYS. */
118 newbox = (integer *)calloc(*ndg,sizeof(integer));
119 newboxp = (integer *)calloc(*ndg,sizeof(integer));
120
121 /* Parameter adjustments */
122 --nov;
123 --lbit;
124 --jmax;
125 --jmin;
126
127 /* Function Body */
128
129 *ier = 0;
130 if (*lx == 1) {
131 goto L410;
132 }
133 /* IF THERE IS ONLY ONE GROUP, RETURN. */
134
135 *(unsigned char *)cfeed = (char) ifeed;
136
137 /* INITIALIZE NUMBER OF NEW BOXES PER GROUP TO ZERO. */
138
139 i__1 = *lx;
140 for (l = 1; l <= i__1; ++l) {
141 newbox[l - 1] = 0;
142 /* L110: */
143 }
144
145 /* INITIALIZE NUMBER OF TOTAL NEW BOXES PER J TO ZERO. */
146
147 for (j = 1; j <= 31; ++j) {
148 ntotbt[j - 1] = 999999999;
149 /* nboxj[j - 1] = 0; */
150 /* L112: */
151 }
152
153 iorigb = (*ibit + *jbit + *kbit) * *lx;
154 /* IBIT = BITS TO PACK THE JMIN( ). */
155 /* JBIT = BITS TO PACK THE LBIT( ). */
156 /* KBIT = BITS TO PACK THE NOV( ). */
157 /* LX = NUMBER OF GROUPS. */
158 ntotbt[*kbit - 1] = iorigb;
159 /* THIS IS THE VALUE OF TOTAL BITS FOR THE ORIGINAL LX */
160 /* GROUPS, WHICH REQUIRES KBITS TO PACK THE GROUP */
161 /* LENGHTS. SETTING THIS HERE MAKES ONE LESS LOOPS */
162 /* NECESSARY BELOW. */
163
164 /* COMPUTE BITS NOW USED FOR THE PARAMETERS DEFINED. */
165
166 /* DETERMINE OTHER POSSIBILITES BY INCREASING LX AND DECREASING */
167 /* NOV( ) WITH VALUES GREATER THAN THRESHOLDS. ASSUME A GROUP IS */
168 /* SPLIT INTO 2 OR MORE GROUPS SO THAT KBIT IS REDUCED WITHOUT */
169 /* CHANGING IBIT OR JBIT. */
170
171 jj = 0;
172
173 /* Computing MIN */
174 i__1 = 30, i__2 = *kbit - 1;
175 /*for (j = min(i__1,i__2); j >= 2; --j) {*/
176 for (j = (i__1 < i__2) ? i__1 : i__2; j >= 2; --j) {
177 /* VALUES GE KBIT WILL NOT REQUIRE SPLITS. ONCE THE TOTAL */
178 /* BITS START INCREASING WITH DECREASING J, STOP. ALSO, THE */
179 /* NUMBER OF BITS REQUIRED IS KNOWN FOR KBITS = NTOTBT(KBIT). */
180
181 newboxt = 0;
182
183 i__1 = *lx;
184 for (l = 1; l <= i__1; ++l) {
185
186 if (nov[l] < ibxx2[j]) {
187 newbox[l - 1] = 0;
188 /* NO SPLITS OR NEW BOXES. */
189 goto L190;
190 } else {
191 novl = nov[l];
192
193 m = (nov[l] - 1) / (ibxx2[j] - 1) + 1;
194 /* M IS FOUND BY SOLVING THE EQUATION BELOW FOR M: */
195 /* (NOV(L)+M-1)/M LT IBXX2(J) */
196 /* M GT (NOV(L)-1)/(IBXX2(J)-1) */
197 /* SET M = (NOV(L)-1)/(IBXX2(J)-1)+1 */
198 L130:
199 novl = (nov[l] + m - 1) / m;
200 /* THE +M-1 IS NECESSARY. FOR INSTANCE, 15 WILL FIT */
201 /* INTO A BOX 4 BITS WIDE, BUT WON'T DIVIDE INTO */
202 /* TWO BOXES 3 BITS WIDE EACH. */
203
204 if (novl < ibxx2[j]) {
205 goto L185;
206 } else {
207 ++m;
208 /* *** WRITE(KFILDO,135)L,NOV(L),NOVL,M,J,IBXX2(J) */
209 /* *** 135 FORMAT(/' AT 135--L,NOV(L),NOVL,M,J,IBXX2(J)',6I10) */
210 goto L130;
211 }
212
213 /* THE ABOVE DO LOOP WILL NEVER COMPLETE. */
214 }
215
216 L185:
217 newbox[l - 1] = m - 1;
218 newboxt = newboxt + m - 1;
219 L190:
220 ;
221 }
222
223 /* nboxj[j - 1] = newboxt; */
224 ntotpr = ntotbt[j];
225 ntotbt[j - 1] = (*ibit + *jbit) * (*lx + newboxt) + j * (*lx +
226 newboxt);
227
228 if (ntotbt[j - 1] >= ntotpr) {
229 jj = j + 1;
230 /* THE PLUS IS USED BECAUSE J DECREASES PER ITERATION. */
231 goto L250;
232 } else {
233
234 /* SAVE THE TOTAL NEW BOXES AND NEWBOX( ) IN CASE THIS */
235 /* IS THE J TO USE. */
236
237 newboxtp = newboxt;
238
239 i__1 = *lx;
240 for (l = 1; l <= i__1; ++l) {
241 newboxp[l - 1] = newbox[l - 1];
242 /* L195: */
243 }
244
245 /* WRITE(KFILDO,197)NEWBOXT,IBXX2(J) */
246 /* 197 FORMAT(/' *****************************************' */
247 /* 1 /' THE NUMBER OF NEWBOXES PER GROUP OF THE TOTAL', */
248 /* 2 I10,' FOR GROUP MAXSIZE PLUS 1 ='I10 */
249 /* 3 /' *****************************************') */
250 /* WRITE(KFILDO,198) (NEWBOX(L),L=1,LX) */
251 /* 198 FORMAT(/' '20I6/(' '20I6)) */
252 }
253
254 /* 205 WRITE(KFILDO,209)KBIT,IORIGB */
255 /* 209 FORMAT(/' ORIGINAL BITS WITH KBIT OF',I5,' =',I10) */
256 /* WRITE(KFILDO,210)(N,N=2,10),(IBXX2(N),N=2,10), */
257 /* 1 (NTOTBT(N),N=2,10),(NBOXJ(N),N=2,10), */
258 /* 2 (N,N=11,20),(IBXX2(N),N=11,20), */
259 /* 3 (NTOTBT(N),N=11,20),(NBOXJ(N),N=11,20), */
260 /* 4 (N,N=21,30),(IBXX2(N),N=11,20), */
261 /* 5 (NTOTBT(N),N=21,30),(NBOXJ(N),N=21,30) */
262 /* 210 FORMAT(/' THE TOTAL BYTES FOR MAXIMUM GROUP LENGTHS BY ROW'// */
263 /* 1 ' J = THE NUMBER OF BITS PER GROUP LENGTH'/ */
264 /* 2 ' IBXX2(J) = THE MAXIMUM GROUP LENGTH PLUS 1 FOR THIS J'/ */
265 /* 3 ' NTOTBT(J) = THE TOTAL BITS FOR THIS J'/ */
266 /* 4 ' NBOXJ(J) = THE NEW GROUPS FOR THIS J'/ */
267 /* 5 4(/10X,9I10)/4(/10I10)/4(/10I10)) */
268
269 /* L200: */
270 }
271
272 L250:
273 pimp = (iorigb - ntotbt[jj - 1]) / (real) iorigb * 100.f;
274 /* WRITE(KFILDO,252)PIMP,KBIT,JJ */
275 /* 252 FORMAT(/' PERCENT IMPROVEMENT =',F6.1, */
276 /* 1 ' BY DECREASING GROUP LENGTHS FROM',I4,' TO',I4,' BITS') */
277 if (pimp >= 2.f) {
278
279 /* WRITE(KFILDO,255)CFEED,NEWBOXTP,IBXX2(JJ) */
280 /* 255 FORMAT(A1,/' *****************************************' */
281 /* 1 /' THE NUMBER OF NEWBOXES PER GROUP OF THE TOTAL', */
282 /* 2 I10,' FOR GROUP MAXSIZE PLUS 1 ='I10 */
283 /* 2 /' *****************************************') */
284 /* WRITE(KFILDO,256) (NEWBOXP(L),L=1,LX) */
285 /* 256 FORMAT(/' '20I6) */
286
287 /* ADJUST GROUP LENGTHS FOR MAXIMUM LENGTH OF JJ BITS. */
288 /* THE MIN PER GROUP AND THE NUMBER OF BITS REQUIRED */
289 /* PER GROUP ARE NOT CHANGED. THIS MAY MEAN THAT A */
290 /* GROUP HAS A MIN (OR REFERENCE) THAT IS NOT ZERO. */
291 /* THIS SHOULD NOT MATTER TO THE UNPACKER. */
292
293 lxnkp = *lx + newboxtp;
294 /* LXNKP = THE NEW NUMBER OF BOXES */
295
296 if (lxnkp > *ndg) {
297 /* DIMENSIONS NOT LARGE ENOUGH. PROBABLY AN ERROR */
298 /* OF SOME SORT. ABORT. */
299 /* WRITE(KFILDO,257)NDG,LXNPK */
300 /* 1 2 3 4 5 6 7 X */
301 /* 257 FORMAT(/' DIMENSIONS OF JMIN, ETC. IN REDUCE =',I8, */
302 /* 1 ' NOT LARGE ENOUGH FOR THE EXPANDED NUMBER OF', */
303 /* 2 ' GROUPS =',I8,'. ABORT REDUCE.') */
304 *ier = 715;
305 goto L410;
306 /* AN ABORT CAUSES THE CALLING PROGRAM TO REEXECUTE */
307 /* WITHOUT CALLING REDUCE. */
308 }
309
310 lxn = lxnkp;
311 /* LXN IS THE NUMBER OF THE BOX IN THE NEW SERIES BEING */
312 /* FILLED. IT DECREASES PER ITERATION. */
313 ibxx2m1 = ibxx2[jj] - 1;
314 /* IBXX2M1 IS THE MAXIMUM NUMBER OF VALUES PER GROUP. */
315
316 for (l = *lx; l >= 1; --l) {
317
318 /* THE VALUES IS NOV( ) REPRESENT THOSE VALUES + NOVREF. */
319 /* WHEN VALUES ARE MOVED TO ANOTHER BOX, EACH VALUE */
320 /* MOVED TO A NEW BOX REPRESENTS THAT VALUE + NOVREF. */
321 /* THIS HAS TO BE CONSIDERED IN MOVING VALUES. */
322
323 if (newboxp[l - 1] * (ibxx2m1 + *novref) + *novref > nov[l] + *
324 novref) {
325 /* IF THE ABOVE TEST IS MET, THEN MOVING IBXX2M1 VALUES */
326 /* FOR ALL NEW BOXES WILL LEAVE A NEGATIVE NUMBER FOR */
327 /* THE LAST BOX. NOT A TOLERABLE SITUATION. */
328 movmin = (nov[l] - newboxp[l - 1] * *novref) / newboxp[l - 1];
329 left = nov[l];
330 /* LEFT = THE NUMBER OF VALUES TO MOVE FROM THE ORIGINAL */
331 /* BOX TO EACH NEW BOX EXCEPT THE LAST. LEFT IS THE */
332 /* NUMBER LEFT TO MOVE. */
333 } else {
334 movmin = ibxx2m1;
335 /* MOVMIN VALUES CAN BE MOVED FOR EACH NEW BOX. */
336 left = nov[l];
337 /* LEFT IS THE NUMBER OF VALUES LEFT TO MOVE. */
338 }
339
340 if (newboxp[l - 1] > 0) {
341 if ((movmin + *novref) * newboxp[l - 1] + *novref <= nov[l] +
342 *novref && (movmin + *novref) * (newboxp[l - 1] + 1)
343 >= nov[l] + *novref) {
344 goto L288;
345 } else {
346 /* ***D WRITE(KFILDO,287)L,MOVMIN,NOVREF,NEWBOXP(L),NOV(L) */
347 /* ***D287 FORMAT(/' AT 287 IN REDUCE--L,MOVMIN,NOVREF,', */
348 /* ***D 1 'NEWBOXP(L),NOV(L)',5I12 */
349 /* ***D 2 ' REDUCE ABORTED.') */
350 /* WRITE(KFILDO,2870) */
351 /* 2870 FORMAT(/' AN ERROR IN REDUCE ALGORITHM. ABORT REDUCE.') */
352 *ier = 714;
353 goto L410;
354 /* AN ABORT CAUSES THE CALLING PROGRAM TO REEXECUTE */
355 /* WITHOUT CALLING REDUCE. */
356 }
357
358 }
359
360 L288:
361 i__1 = newboxp[l - 1] + 1;
362 for (j = 1; j <= i__1; ++j) {
363 /*move = min(movmin,left);*/
364 move = (movmin < left) ? movmin : left;
365 jmin[lxn] = jmin[l];
366 jmax[lxn] = jmax[l];
367 lbit[lxn] = lbit[l];
368 nov[lxn] = move;
369 --lxn;
370 left -= move + *novref;
371 /* THE MOVE OF MOVE VALUES REALLY REPRESENTS A MOVE OF */
372 /* MOVE + NOVREF VALUES. */
373 /* L290: */
374 }
375
376 if (left != -(*novref)) {
377 /* *** WRITE(KFILDO,292)L,LXN,MOVE,LXNKP,IBXX2(JJ),LEFT,NOV(L), */
378 /* *** 1 MOVMIN */
379 /* *** 292 FORMAT(' AT 292 IN REDUCE--L,LXN,MOVE,LXNKP,', */
380 /* *** 1 'IBXX2(JJ),LEFT,NOV(L),MOVMIN'/8I12) */
381 }
382
383 /* L300: */
384 }
385
386 *lx = lxnkp;
387 /* LX IS NOW THE NEW NUMBER OF GROUPS. */
388 *kbit = jj;
389 /* KBIT IS NOW THE NEW NUMBER OF BITS REQUIRED FOR PACKING */
390 /* GROUP LENGHTS. */
391 }
392
393 /* WRITE(KFILDO,406)CFEED,LX */
394 /* 406 FORMAT(A1,/' *****************************************' */
395 /* 1 /' THE GROUP SIZES NOV( ) AFTER REDUCTION IN SIZE', */
396 /* 2 ' FOR'I10,' GROUPS', */
397 /* 3 /' *****************************************') */
398 /* WRITE(KFILDO,407) (NOV(J),J=1,LX) */
399 /* 407 FORMAT(/' '20I6) */
400 /* WRITE(KFILDO,408)CFEED,LX */
401 /* 408 FORMAT(A1,/' *****************************************' */
402 /* 1 /' THE GROUP MINIMA JMIN( ) AFTER REDUCTION IN SIZE', */
403 /* 2 ' FOR'I10,' GROUPS', */
404 /* 3 /' *****************************************') */
405 /* WRITE(KFILDO,409) (JMIN(J),J=1,LX) */
406 /* 409 FORMAT(/' '20I6) */
407
408 L410:
409 if ( newbox != 0 ) free(newbox);
410 if ( newboxp != 0 ) free(newboxp);
411 return 0;
412 } /* reduce_ */
413