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