Created a tag for the 2012 HWRF baseline tests.
[WPS-merge.git] / hwrf-baseline-20120103-1354 / ungrib / src / ngl / w3 / w3fi75.f
blobe58cf9c9ed78f36d63f419297ce69eb677f95cca
1 SUBROUTINE W3FI75 (IBITL,ITYPE,ITOSS,FLD,IFLD,IBMAP,IBDSFL,
2 & NPTS,BDS11,IPFLD,PFLD,LEN,LENBDS,IBERR,PDS,IGDS)
3 C $$ SUBPROGRAM DOCUMENTATION BLOCK
4 C . . . .
5 C SUBPROGRAM: W3FI75 GRIB PACK DATA AND FORM BDS OCTETS(1-11)
6 C PRGMMR: FARLEY ORG: NMC421 DATE:94-11-22
8 C ABSTRACT: THIS ROUTINE PACKS A GRIB FIELD AND FORMS OCTETS(1-11)
9 C OF THE BINARY DATA SECTION (BDS).
11 C PROGRAM HISTORY LOG:
12 C 92-07-10 M. FARLEY ORIGINAL AUTHOR
13 C 92-10-01 R.E.JONES CORRECTION FOR FIELD OF CONSTANT DATA
14 C 92-10-16 R.E.JONES GET RID OF ARRAYS FP AND INT
15 C 93-08-06 CAVANAUGH ADDED ROUTINES FI7501, FI7502, FI7503
16 C TO ALLOW SECOND ORDER PACKING IN PDS.
17 C 93-07-21 STACKPOLE ASSORTED REPAIRS TO GET 2ND DIFF PACK IN
18 C 93-10-28 CAVANAUGH COMMENTED OUT NONOPERATIONAL PRINTS AND
19 C WRITE STATEMENTS
20 C 93-12-15 CAVANAUGH CORRECTED LOCATION OF START OF FIRST ORDER
21 C VALUES AND START OF SECOND ORDER VALUES TO
22 C REFLECT A BYTE LOCATION IN THE BDS INSTEAD
23 C OF AN OFFSET IN SUBROUTINE FI7501.
24 C 94-01-27 CAVANAUGH ADDED IGDS AS INPUT ARGUMENT TO THIS ROUTINE
25 C AND ADDED PDS AND IGDS ARRAYS TO THE CALL TO
26 C W3FI82 TO PROVIDE INFORMATION NEEDED FOR
27 C BOUSTROPHEDONIC PROCESSING.
28 C 94-05-25 CAVANAUGH SUBROUTINE FI7503 HAS BEEN ADDED TO PROVIDE
29 C FOR ROW BY ROW OR COLUMN BY COLUMN SECOND
30 C ORDER PACKING. THIS FEATURE CAN BE ACTIVATED
31 C BY SETTING IBDSFL(7) TO ZERO.
32 C 94-07-08 CAVANAUGH COMMENTED OUT PRINT STATEMENTS USED FOR DEBUG
33 C 94-11-22 FARLEY ENLARGED WORK ARRAYS TO HANDLE .5DEGREE GRIDS
34 C 95-06-01 R.E.JONES CORRECTION FOR NUMBER OF UNUSED BITS AT END
35 C OF SECTION 4, IN BDS BYTE 4, BITS 5-8.
36 C 95-10-31 IREDELL REMOVED SAVES AND PRINTS
38 C USAGE: CALL W3FI75 (IBITL,ITYPE,ITOSS,FLD,IFLD,IBMAP,IBDSFL,
39 C & NPTS,BDS11,IPFLD,PFLD,LEN,LENBDS,IBERR,PDS,IGDS)
40 C INPUT ARGUMENT LIST:
41 C IBITL - 0, COMPUTER COMPUTES PACKING LENGTH FROM POWER
42 C OF 2 THAT BEST FITS THE DATA.
43 C 8, 12, ETC. COMPUTER RESCALES DATA TO FIT INTO
44 C SET NUMBER OF BITS.
45 C ITYPE - 0 = IF INPUT DATA IS FLOATING POINT (FLD)
46 C 1 = IF INPUT DATA IS INTEGER (IFLD)
47 C ITOSS - 0 = NO BIT MAP IS INCLUDED (DON'T TOSS DATA)
48 C 1 = TOSS NULL DATA ACCORDING TO IBMAP
49 C FLD - REAL ARRAY OF DATA TO BE PACKED IF ITYPE=0
50 C IFLD - INTEGER ARRAY TO BE PACKED IF ITYPE=1
51 C IBMAP - BIT MAP SUPPLIED FROM USER
52 C IBDSFL - INTEGER ARRAY CONTAINING TABLE 11 FLAG INFO
53 C BDS OCTET 4:
54 C (1) 0 = GRID POINT DATA
55 C 1 = SPHERICAL HARMONIC COEFFICIENTS
56 C (2) 0 = SIMPLE PACKING
57 C 1 = SECOND ORDER PACKING
58 C (3) 0 = ORIGINAL DATA WERE FLOATING POINT VALUES
59 C 1 = ORIGINAL DATA WERE INTEGER VALUES
60 C (4) 0 = NO ADDITIONAL FLAGS AT OCTET 14
61 C 1 = OCTET 14 CONTAINS FLAG BITS 5-12
62 C (5) 0 = RESERVED - ALWAYS SET TO 0
63 C (6) 0 = SINGLE DATUM AT EACH GRID POINT
64 C 1 = MATRIX OF VALUES AT EACH GRID POINT
65 C (7) 0 = NO SECONDARY BIT MAPS
66 C 1 = SECONDARY BIT MAPS PRESENT
67 C (8) 0 = SECOND ORDER VALUES HAVE CONSTANT WIDTH
68 C 1 = SECOND ORDER VALUES HAVE DIFFERENT WIDTHS
69 C NPTS - NUMBER OF GRIDPOINTS IN ARRAY TO BE PACKED
70 C IGDS - ARRAY OF GDS INFORMATION
72 C OUTPUT ARGUMENT LIST:
73 C BDS11 - FIRST 11 OCTETS OF BDS
74 C PFLD - PACKED GRIB FIELD
75 C LEN - LENGTH OF PFLD
76 C LENBDS - LENGTH OF BDS
77 C IBERR - 1, ERROR CONVERTING IEEE F.P. NUMBER TO IBM370 F.P.
79 C REMARKS: SUBPROGRAM CAN BE CALLED FROM A MULTIPROCESSING ENVIRONMENT.
81 C ATTRIBUTES:
82 C LANGUAGE: IBM VS FORTRAN 77, CRAY CFT77 FORTRAN
83 C MACHINE: HDS, CRAY C916/256, Y-MP8/64, Y-MP EL92/256
85 C $$
87 REAL FLD(*)
88 C REAL FWORK(260000)
90 C FWORK CAN USE DYNAMIC ALLOCATION OF MEMORY ON CRAY
92 REAL FWORK(NPTS)
93 REAL RMIN,REFNCE
95 INTEGER IPFLD(*)
96 INTEGER IBDSFL(*)
97 INTEGER IBMAP(*)
98 INTEGER IFLD(*),IGDS(*)
99 C INTEGER IWORK(260000)
101 C IWORK CAN USE DYNAMIC ALLOCATION OF MEMORY ON CRAY
103 INTEGER IWORK(NPTS)
105 LOGICAL CONST
107 CHARACTER * 1 BDS11(11),PDS(*)
108 CHARACTER * 1 PFLD(*)
109 CHARACTER * 1 CIEXP(8)
110 CHARACTER * 1 CIMANT(8)
112 EQUIVALENCE (IEXP,CIEXP(1))
113 EQUIVALENCE (IMANT,CIMANT(1))
115 C 1.0 PACK THE FIELD.
117 C 1.1 TOSS DATA IF BITMAP BEING USED,
118 C MOVING 'DATA' TO WORK AREA...
120 CONST = .FALSE.
121 IBERR = 0
122 IW = 0
124 IF (ITOSS .EQ. 1) THEN
125 IF (ITYPE .EQ. 0) THEN
126 DO 110 IT=1,NPTS
127 IF (IBMAP(IT) .EQ. 1) THEN
128 IW = IW + 1
129 FWORK(IW) = FLD(IT)
130 ENDIF
131 110 CONTINUE
132 NPTS = IW
133 ELSE IF (ITYPE .EQ. 1) THEN
134 DO 111 IT=1,NPTS
135 IF (IBMAP(IT) .EQ. 1) THEN
136 IW = IW + 1
137 IWORK(IW) = IFLD(IT)
138 ENDIF
139 111 CONTINUE
140 NPTS = IW
141 ENDIF
143 C ELSE, JUST MOVE DATA TO WORK ARRAY
145 ELSE IF (ITOSS .EQ. 0) THEN
146 IF (ITYPE .EQ. 0) THEN
147 DO 112 IT=1,NPTS
148 FWORK(IT) = FLD(IT)
149 112 CONTINUE
150 ELSE IF (ITYPE .EQ. 1) THEN
151 DO 113 IT=1,NPTS
152 IWORK(IT) = IFLD(IT)
153 113 CONTINUE
154 ENDIF
155 ENDIF
157 C 1.2 CONVERT DATA IF NEEDED PRIOR TO PACKING.
158 C (INTEGER TO F.P. OR F.P. TO INTEGER)
159 C ITYPE = 0...FLOATING POINT DATA
160 C IBITL = 0...PACK IN LEAST # BITS...CONVERT TO INTEGER
161 C ITYPE = 1...INTEGER DATA
162 C IBITL > 0...PACK IN FIXED # BITS...CONVERT TO FLOATING POINT
164 IF (ITYPE .EQ. 0 .AND. IBITL .EQ. 0) THEN
165 DO 120 IF=1,NPTS
166 IWORK(IF) = NINT(FWORK(IF))
167 120 CONTINUE
168 ELSE IF (ITYPE .EQ. 1 .AND. IBITL .NE. 0) THEN
169 DO 123 IF=1,NPTS
170 FWORK(IF) = FLOAT(IWORK(IF))
171 123 CONTINUE
172 ENDIF
174 C 1.3 PACK THE DATA.
176 IF (IBDSFL(2).NE.0) THEN
177 C SECOND ORDER PACKING
179 C PRINT*,' DOING SECOND ORDER PACKING...'
180 IF (IBITL.EQ.0) THEN
182 C PRINT*,' AND VARIABLE BIT PACKING'
184 C WORKING WITH INTEGER VALUES
185 C SINCE DOING VARIABLE BIT PACKING
187 MAX = IWORK(1)
188 MIN = IWORK(1)
189 DO 300 I = 2, NPTS
190 IF (IWORK(I).LT.MIN) THEN
191 MIN = IWORK(I)
192 ELSE IF (IWORK(I).GT.MAX) THEN
193 MAX = IWORK(I)
194 END IF
195 300 CONTINUE
196 C EXTRACT MINIMA
197 DO 400 I = 1, NPTS
198 C IF (IWORK(I).LT.0) THEN
199 C PRINT *,'MINIMA 400',I,IWORK(I),NPTS
200 C END IF
201 IWORK(I) = IWORK(I) - MIN
202 400 CONTINUE
203 REFNCE = MIN
204 IDIFF = MAX - MIN
205 C PRINT *,'REFERENCE VALUE',REFNCE
207 C WRITE (6,FMT='('' MINIMA REMOVED = '',/,
208 C & 10(3X,10I10,/))') (IWORK(I),I=1,6)
209 C WRITE (6,FMT='('' END OF ARRAY = '',/,
210 C & 10(3X,10I10,/))') (IWORK(I),I=NPTS-5,NPTS)
212 C FIND BIT WIDTH OF IDIFF
214 CALL FI7505 (IDIFF,KWIDE)
215 C PRINT*,' BIT WIDTH FOR ORIGINAL DATA', KWIDE
216 ISCAL2 = 0
218 C MULTIPLICATIVE SCALE FACTOR SET TO 1
219 C IN ANTICIPATION OF POSSIBLE USE IN GLAHN 2DN DIFF
221 SCAL2 = 1.
223 ELSE
225 C PRINT*,' AND FIXED BIT PACKING, IBITL = ', IBITL
226 C FIXED BIT PACKING
227 C - LENGTH OF FIELD IN IBITL
228 C - MUST BE REAL DATA
229 C FLOATING POINT INPUT
231 RMAX = FWORK(1)
232 RMIN = FWORK(1)
233 DO 100 I = 2, NPTS
234 IF (FWORK(I).LT.RMIN) THEN
235 RMIN = FWORK(I)
236 ELSE IF (FWORK(I).GT.RMAX) THEN
237 RMAX = FWORK(I)
238 END IF
239 100 CONTINUE
240 REFNCE = RMIN
241 C PRINT *,'100 REFERENCE',REFNCE
242 C EXTRACT MINIMA
243 DO 200 I = 1, NPTS
244 FWORK(I) = FWORK(I) - RMIN
245 200 CONTINUE
246 C PRINT *,'REFERENCE VALUE',REFNCE
247 C WRITE (6,FMT='('' MINIMA REMOVED = '',/,
248 C & 10(3X,10F8.2,/))') (FWORK(I),I=1,6)
249 C WRITE (6,FMT='('' END OF ARRAY = '',/,
250 C & 10(3X,10F8.2,/))') (FWORK(I),I=NPTS-5,NPTS)
251 C FIND LARGEST DELTA
252 IDELT = NINT(RMAX - RMIN)
253 C DO BINARY SCALING
254 C FIND OUT WHAT BINARY SCALE FACTOR
255 C PERMITS CONTAINMENT OF
256 C LARGEST DELTA
257 CALL FI7505 (IDELT,IWIDE)
259 C BINARY SCALING
261 ISCAL2 = IWIDE - IBITL
262 C PRINT *,'SCALING NEEDED TO FIT =',ISCAL2
263 C PRINT*,' RANGE OF = ',IDELT
265 C EXPAND DATA WITH BINARY SCALING
266 C CONVERT TO INTEGER
267 SCAL2 = 2.0**ISCAL2
268 SCAL2 = 1./ SCAL2
269 DO 600 I = 1, NPTS
270 IWORK(I) = NINT(FWORK(I) * SCAL2)
271 600 CONTINUE
272 KWIDE = IBITL
273 END IF
275 C *****************************************************************
277 C FOLLOWING IS FOR GLAHN SECOND DIFFERENCING
278 C NOT STANDARD GRIB
280 C TEST FOR SECOND DIFFERENCE PACKING
281 C BASED OF SIZE OF PDS - SIZE IN FIRST 3 BYTES
283 CALL GBYTE (PDS,IPDSIZ,0,24)
284 IF (IPDSIZ.EQ.50) THEN
285 C PRINT*,' DO SECOND DIFFERENCE PACKING '
287 C GLAHN PACKING TO 2ND DIFFS
289 C WRITE (6,FMT='('' CALL TO W3FI82 WITH = '',/,
290 C & 10(3X,10I6,/))') (IWORK(I),I=1,NPTS)
292 CALL W3FI82 (IWORK,FVAL1,FDIFF1,NPTS,PDS,IGDS)
294 C PRINT *,'GLAHN',FVAL1,FDIFF1
295 C WRITE (6,FMT='('' OUT FROM W3FI82 WITH = '',/,
296 C & 10(3X,10I6,/))') (IWORK(I),I=1,NPTS)
298 C MUST NOW RE-REMOVE THE MINIMUM VALUE
299 C OF THE SECOND DIFFERENCES TO ASSURE
300 C ALL POSITIVE NUMBERS FOR SECOND ORDER GRIB PACKING
302 C ORIGINAL REFERENCE VALUE ADDED TO FIRST POINT
303 C VALUE FROM THE 2ND DIFF PACKER TO BE ADDED
304 C BACK IN WHEN THE 2ND DIFF VALUES ARE
305 C RECONSTRUCTED BACK TO THE BASIC VALUES
307 C ALSO, THE REFERENCE VALUE IS
308 C POWER-OF-TWO SCALED TO MATCH
309 C FVAL1. ALL OF THIS SCALING
310 C WILL BE REMOVED AFTER THE
311 C GLAHN SECOND DIFFERENCING IS UNDONE.
312 C THE SCALING FACTOR NEEDED TO DO THAT
313 C IS SAVED IN THE PDS AS A SIGNED POSITIVE
314 C TWO BYTE INTEGER
316 C THE SCALING FOR THE 2ND DIF PACKED
317 C VALUES IS PROPERLY SET TO ZERO
319 FVAL1 = FVAL1 + REFNCE*SCAL2
320 C FIRST TEST TO SEE IF
321 C ON 32 OR 64 BIT COMPUTER
322 CALL W3FI01(LW)
323 IF (LW.EQ.4) THEN
324 CALL W3FI76 (FVAL1,IEXP,IMANT,32)
325 ELSE
326 CALL W3FI76 (FVAL1,IEXP,IMANT,64)
327 END IF
328 CALL SBYTE (PDS,IEXP,320,8)
329 CALL SBYTE (PDS,IMANT,328,24)
331 IF (LW.EQ.4) THEN
332 CALL W3FI76 (FDIFF1,IEXP,IMANT,32)
333 ELSE
334 CALL W3FI76 (FDIFF1,IEXP,IMANT,64)
335 END IF
336 CALL SBYTE (PDS,IEXP,352,8)
337 CALL SBYTE (PDS,IMANT,360,24)
339 C TURN ISCAL2 INTO SIGNED POSITIVE INTEGER
340 C AND STORE IN TWO BYTES
342 IF(ISCAL2.GE.0) THEN
343 CALL SBYTE (PDS,ISCAL2,384,16)
344 ELSE
345 CALL SBYTE (PDS,1,384,1)
346 ISCAL2 = - ISCAL2
347 CALL SBYTE( PDS,ISCAL2,385,15)
348 ENDIF
350 MAX = IWORK(1)
351 MIN = IWORK(1)
352 DO 700 I = 2, NPTS
353 IF (IWORK(I).LT.MIN) THEN
354 MIN = IWORK(I)
355 ELSE IF (IWORK(I).GT.MAX) THEN
356 MAX = IWORK(I)
357 END IF
358 700 CONTINUE
359 C EXTRACT MINIMA
360 DO 710 I = 1, NPTS
361 IWORK(I) = IWORK(I) - MIN
362 710 CONTINUE
363 REFNCE = MIN
364 C PRINT *,'710 REFERENCE',REFNCE
365 ISCAL2 = 0
367 C AND RESET VALUE OF KWIDE - THE BIT WIDTH
368 C FOR THE RANGE OF THE VALUES
370 IDIFF = MAX - MIN
371 CALL FI7505 (IDIFF,KWIDE)
373 C PRINT*,'BIT WIDTH (KWIDE) OF 2ND DIFFS', KWIDE
375 C **************************** END OF GLAHN PACKING ************
376 ELSE IF (IBDSFL(2).EQ.1.AND.IBDSFL(7).EQ.0) THEN
377 C HAVE SECOND ORDER PACKING WITH NO SECOND ORDER
378 C BIT MAP. ERGO ROW BY ROW - COL BY COL
379 CALL FI7503 (IWORK,IPFLD,NPTS,IBDSFL,BDS11,
380 * LEN,LENBDS,PDS,REFNCE,ISCAL2,KWIDE,IGDS)
381 RETURN
382 END IF
383 C WRITE (6,FMT='('' CALL TO FI7501 WITH = '',/,
384 C & 10(3X,10I6,/))') (IWORK(I),I=1,NPTS)
385 C WRITE (6,FMT='('' END OF ARRAY = '',/,
386 C & 10(3X,10I6,/))') (IWORK(I),I=NPTS-5,NPTS)
387 C PRINT*,' REFNCE,ISCAL2, KWIDE AT CALL TO FI7501',
388 C & REFNCE, ISCAL2,KWIDE
390 C SECOND ORDER PACKING
392 CALL FI7501 (IWORK,IPFLD,NPTS,IBDSFL,BDS11,
393 * LEN,LENBDS,PDS,REFNCE,ISCAL2,KWIDE)
395 C BDS COMPLETELY ASSEMBLED IN FI7501 FOR SECOND ORDER
396 C PACKING.
398 ELSE
399 C SIMPLE PACKING
401 C PRINT*,' SIMPLE FIRST ORDER PACKING...'
402 IF (IBITL.EQ.0) THEN
403 C PRINT*,' WITH VARIABLE BIT LENGTH'
405 C WITH VARIABLE BIT LENGTH, ADJUSTED
406 C TO ACCOMMODATE LARGEST VALUE
407 C BINARY SCALING ALWAYS = 0
409 CALL W3FI58(IWORK,NPTS,IWORK,PFLD,NBITS,LEN,KMIN)
410 RMIN = KMIN
411 REFNCE = RMIN
412 ISCALE = 0
413 C PRINT*,' BIT LENGTH CAME OUT AT ...',NBITS
415 C SET CONST .TRUE. IF ALL VALUES ARE THE SAME
417 IF (LEN.EQ.0.AND.NBITS.EQ.0) CONST = .TRUE.
419 ELSE
420 C PRINT*,' FIXED BIT LENGTH, IBITL = ', IBITL
422 C FIXED BIT LENGTH PACKING (VARIABLE PRECISION)
423 C VALUES SCALED BY POWER OF 2 (ISCALE) TO
424 C FIT LARGEST VALUE INTO GIVEN BIT LENGTH (IBITL)
426 CALL W3FI59(FWORK,NPTS,IBITL,IWORK,PFLD,ISCALE,LEN,RMIN)
427 REFNCE = RMIN
428 C PRINT *,' SCALING NEEDED TO FIT IS ...', ISCALE
429 NBITS = IBITL
431 C SET CONST .TRUE. IF ALL VALUES ARE THE SAME
433 IF (LEN.EQ.0) THEN
434 CONST = .TRUE.
435 NBITS = 0
436 END IF
437 END IF
439 C COMPUTE LENGTH OF BDS IN OCTETS
441 INUM = NPTS * NBITS + 88
442 C PRINT *,'NUMBER OF BITS BEFORE FILL ADDED',INUM
444 C NUMBER OF FILL BITS
445 NFILL = 0
446 NLEFT = MOD(INUM,16)
447 IF (NLEFT.NE.0) THEN
448 INUM = INUM + 16 - NLEFT
449 NFILL = 16 - NLEFT
450 END IF
451 C PRINT *,'NUMBER OF BITS AFTER FILL ADDED',INUM
452 C LENGTH OF BDS IN BYTES
453 LENBDS = INUM / 8
455 C 2.0 FORM THE BINARY DATA SECTION (BDS).
457 C CONCANTENATE ALL FIELDS FOR BDS
459 C BYTES 1-3
460 CALL SBYTE (BDS11,LENBDS,0,24)
462 C BYTE 4
463 C FLAGS
464 CALL SBYTE (BDS11,IBDSFL(1),24,1)
465 CALL SBYTE (BDS11,IBDSFL(2),25,1)
466 CALL SBYTE (BDS11,IBDSFL(3),26,1)
467 CALL SBYTE (BDS11,IBDSFL(4),27,1)
468 C NR OF FILL BITS
469 CALL SBYTE (BDS11,NFILL,28,4)
471 C FILL OCTETS 5-6 WITH THE SCALE FACTOR.
473 C BYTE 5-6
474 IF (ISCALE.LT.0) THEN
475 CALL SBYTE (BDS11,1,32,1)
476 ISCALE = - ISCALE
477 CALL SBYTE (BDS11,ISCALE,33,15)
478 ELSE
479 CALL SBYTE (BDS11,ISCALE,32,16)
480 END IF
482 C FILL OCTET 7-10 WITH THE REFERENCE VALUE
483 C CONVERT THE FLOATING POINT OF YOUR MACHINE TO IBM370 32 BIT
484 C FLOATING POINT NUMBER
486 C BYTE 7-10
487 C REFERENCE VALUE
488 C FIRST TEST TO SEE IF
489 C ON 32 OR 64 BIT COMPUTER
490 CALL W3FI01(LW)
491 IF (LW.EQ.4) THEN
492 CALL W3FI76 (REFNCE,IEXP,IMANT,32)
493 ELSE
494 CALL W3FI76 (REFNCE,IEXP,IMANT,64)
495 END IF
496 CALL SBYTE (BDS11,IEXP,48,8)
497 CALL SBYTE (BDS11,IMANT,56,24)
500 C FILL OCTET 11 WITH THE NUMBER OF BITS.
502 C BYTE 11
503 CALL SBYTE (BDS11,NBITS,80,8)
504 END IF
506 RETURN
508 SUBROUTINE FI7501 (IWORK,IPFLD,NPTS,IBDSFL,BDS11,
509 * LEN,LENBDS,PDS,REFNCE,ISCAL2,KWIDE)
510 C $$ SUBPROGRAM DOCUMENTATION BLOCK
511 C . . . .
512 C SUBPROGRAM: FI7501 BDS SECOND ORDER PACKING
513 C PRGMMR: CAVANAUGH ORG: W/NMC42 DATE: 93-08-06
515 C ABSTRACT: PERFORM SECONDARY PACKING ON GRID POINT DATA,
516 C GENERATING ALL BDS INFORMATION.
518 C PROGRAM HISTORY LOG:
519 C 93-08-06 CAVANAUGH
520 C 93-12-15 CAVANAUGH CORRECTED LOCATION OF START OF FIRST ORDER
521 C VALUES AND START OF SECOND ORDER VALUES TO
522 C REFLECT A BYTE LOCATION IN THE BDS INSTEAD
523 C OF AN OFFSET.
524 C 95-10-31 IREDELL REMOVED SAVES AND PRINTS
526 C USAGE: CALL FI7501 (IWORK,IPFLD,NPTS,IBDSFL,BDS11,
527 C * LEN,LENBDS,PDS,REFNCE,ISCAL2,KWIDE)
528 C INPUT ARGUMENT LIST:
529 C IWORK - INTEGER SOURCE ARRAY
530 C NPTS - NUMBER OF POINTS IN IWORK
531 C IBDSFL - FLAGS
533 C OUTPUT ARGUMENT LIST: (INCLUDING WORK ARRAYS)
534 C IPFLD - CONTAINS BDS FROM BYTE 12 ON
535 C BDS11 - CONTAINS FIRST 11 BYTES FOR BDS
536 C LEN - NUMBER OF BYTES FROM 12 ON
537 C LENBDS - TOTAL LENGTH OF BDS
539 C REMARKS: SUBPROGRAM CAN BE CALLED FROM A MULTIPROCESSING ENVIRONMENT.
541 C ATTRIBUTES:
542 C LANGUAGE: IBM VS FORTRAN 77, CRAY CFT77 FORTRAN
543 C MACHINE: HDS, CRAY C916/256, Y-MP8/64, Y-MP EL92/256
545 C $$
546 CHARACTER*1 BDS11(*),PDS(*)
548 REAL REFNCE
550 INTEGER ISCAL2,KWIDE
551 INTEGER LENBDS
552 INTEGER IPFLD(*)
553 INTEGER LEN,KBDS(22)
554 INTEGER IWORK(*)
555 C OCTET NUMBER IN SECTION, FIRST ORDER PACKING
556 C INTEGER KBDS(12)
557 C FLAGS
558 INTEGER IBDSFL(*)
559 C EXTENDED FLAGS
560 C INTEGER KBDS(14)
561 C OCTET NUMBER FOR SECOND ORDER PACKING
562 C INTEGER KBDS(15)
563 C NUMBER OF FIRST ORDER VALUES
564 C INTEGER KBDS(17)
565 C NUMBER OF SECOND ORDER PACKED VALUES
566 C INTEGER KBDS(19)
567 C WIDTH OF SECOND ORDER PACKING
568 INTEGER ISOWID(50000)
569 C SECONDARY BIT MAP
570 INTEGER ISOBMP(8200)
571 C FIRST ORDER PACKED VALUES
572 INTEGER IFOVAL(50000)
573 C SECOND ORDER PACKED VALUES
574 INTEGER ISOVAL(100000)
576 C INTEGER KBDS(11)
577 C BIT WIDTH TABLE
578 INTEGER IBITS(31)
580 DATA IBITS/1,3,7,15,31,63,127,255,511,1023,
581 * 2047,4095,8191,16383,32767,65535,131072,
582 * 262143,524287,1048575,2097151,4194303,
583 * 8388607,16777215,33554431,67108863,
584 * 134217727,268435455,536870911,
585 * 1073741823,2147483647/
586 C ----------------------------------
587 C INITIALIZE ARRAYS
588 DO 100 I = 1, 50000
589 ISOWID(I) = 0
590 IFOVAL(I) = 0
591 100 CONTINUE
593 DO 101 I = 1, 8200
594 ISOBMP(I) = 0
595 101 CONTINUE
596 DO 102 I = 1, 100000
597 ISOVAL(I) = 0
598 102 CONTINUE
599 C INITIALIZE POINTERS
600 C SECONDARY BIT WIDTH POINTER
601 IWDPTR = 0
602 C SECONDARY BIT MAP POINTER
603 IBMP2P = 0
604 C FIRST ORDER VALUE POINTER
605 IFOPTR = 0
606 C BYTE POINTER TO START OF 1ST ORDER VALUES
607 KBDS(12) = 0
608 C BYTE POINTER TO START OF 2ND ORDER VALUES
609 KBDS(15) = 0
610 C TO CONTAIN NUMBER OF FIRST ORDER VALUES
611 KBDS(17) = 0
612 C TO CONTAIN NUMBER OF SECOND ORDER VALUES
613 KBDS(19) = 0
614 C SECOND ORDER PACKED VALUE POINTER
615 ISOPTR = 0
616 C =======================================================
618 C DATA IS IN IWORK
620 KBDS(11) = KWIDE
622 C DATA PACKING
624 ITER = 0
625 INEXT = 1
626 ISTART = 1
627 C -----------------------------------------------------------
628 KOUNT = 0
629 C DO 1 I = 1, NPTS, 10
630 C PRINT *,I,(IWORK(K),K=I, I+9)
631 C 1 CONTINUE
632 2000 CONTINUE
633 ITER = ITER + 1
634 C PRINT *,'NEXT ITERATION STARTS AT',ISTART
635 IF (ISTART.GT.NPTS) THEN
636 GO TO 4000
637 ELSE IF (ISTART.EQ.NPTS) THEN
638 KPTS = 1
639 MXDIFF = 0
640 GO TO 2200
641 END IF
643 C LOOK FOR REPITITIONS OF A SINGLE VALUE
644 CALL FI7502 (IWORK,ISTART,NPTS,ISAME)
645 IF (ISAME.GE.15) THEN
646 KOUNT = KOUNT + 1
647 C PRINT *,'FI7501 - FOUND IDENTICAL SET OF ',ISAME
648 MXDIFF = 0
649 KPTS = ISAME
650 ELSE
652 C LOOK FOR SETS OF VALUES IN TREND SELECTED RANGE
653 CALL FI7513 (IWORK,ISTART,NPTS,NMAX,NMIN,INRNGE)
654 C PRINT *,'ISTART ',ISTART,' INRNGE',INRNGE,NMAX,NMIN
655 IEND = ISTART + INRNGE - 1
656 C DO 2199 NM = ISTART, IEND, 10
657 C PRINT *,' ',(IWORK(NM+JK),JK=0,9)
658 C2199 CONTINUE
659 MXDIFF = NMAX - NMIN
660 KPTS = INRNGE
661 END IF
662 2200 CONTINUE
663 C PRINT *,' RANGE ',MXDIFF,' MAX',NMAX,' MIN',NMIN
664 C INCREMENT NUMBER OF FIRST ORDER VALUES
665 KBDS(17) = KBDS(17) + 1
666 C ENTER FIRST ORDER VALUE
667 IF (MXDIFF.GT.0) THEN
668 DO 2220 LK = 0, KPTS-1
669 IWORK(ISTART+LK) = IWORK(ISTART+LK) - NMIN
670 2220 CONTINUE
671 CALL SBYTE (IFOVAL,NMIN,IFOPTR,KBDS(11))
672 ELSE
673 CALL SBYTE (IFOVAL,IWORK(ISTART),IFOPTR,KBDS(11))
674 END IF
675 IFOPTR = IFOPTR + KBDS(11)
676 C PROCESS SECOND ORDER BIT WIDTH
677 IF (MXDIFF.GT.0) THEN
678 DO 2330 KWIDE = 1, 31
679 IF (MXDIFF.LE.IBITS(KWIDE)) THEN
680 GO TO 2331
681 END IF
682 2330 CONTINUE
683 2331 CONTINUE
684 ELSE
685 KWIDE = 0
686 END IF
687 CALL SBYTE (ISOWID,KWIDE,IWDPTR,8)
688 IWDPTR = IWDPTR + 8
689 C PRINT *,KWIDE,' IFOVAL=',NMIN,IWORK(ISTART),KPTS
690 C IF KWIDE NE 0, SAVE SECOND ORDER VALUE
691 IF (KWIDE.GT.0) THEN
692 CALL SBYTES (ISOVAL,IWORK(ISTART),ISOPTR,KWIDE,0,KPTS)
693 ISOPTR = ISOPTR + KPTS * KWIDE
694 KBDS(19) = KBDS(19) + KPTS
695 C PRINT *,' SECOND ORDER VALUES'
696 C PRINT *,(IWORK(ISTART+I),I=0,KPTS-1)
697 END IF
698 C ADD TO SECOND ORDER BITMAP
699 CALL SBYTE (ISOBMP,1,IBMP2P,1)
700 IBMP2P = IBMP2P + KPTS
701 ISTART = ISTART + KPTS
702 GO TO 2000
703 C --------------------------------------------------------------
704 4000 CONTINUE
705 C PRINT *,'THERE WERE ',ITER,' SECOND ORDER GROUPS'
706 C PRINT *,'THERE WERE ',KOUNT,' STRINGS OF CONSTANTS'
707 C CONCANTENATE ALL FIELDS FOR BDS
709 C REMAINDER GOES INTO IPFLD
710 IPTR = 0
711 C BYTES 12-13
712 C VALUE FOR N1
713 C LEAVE SPACE FOR THIS
714 IPTR = IPTR + 16
715 C BYTE 14
716 C EXTENDED FLAGS
717 CALL SBYTE (IPFLD,IBDSFL(5),IPTR,1)
718 IPTR = IPTR + 1
719 CALL SBYTE (IPFLD,IBDSFL(6),IPTR,1)
720 IPTR = IPTR + 1
721 CALL SBYTE (IPFLD,IBDSFL(7),IPTR,1)
722 IPTR = IPTR + 1
723 CALL SBYTE (IPFLD,IBDSFL(8),IPTR,1)
724 IPTR = IPTR + 1
725 CALL SBYTE (IPFLD,IBDSFL(9),IPTR,1)
726 IPTR = IPTR + 1
727 CALL SBYTE (IPFLD,IBDSFL(10),IPTR,1)
728 IPTR = IPTR + 1
729 CALL SBYTE (IPFLD,IBDSFL(11),IPTR,1)
730 IPTR = IPTR + 1
731 CALL SBYTE (IPFLD,IBDSFL(12),IPTR,1)
732 IPTR = IPTR + 1
733 C BYTES 15-16
734 C SKIP OVER VALUE FOR N2
735 IPTR = IPTR + 16
736 C BYTES 17-18
737 C P1
738 CALL SBYTE (IPFLD,KBDS(17),IPTR,16)
739 IPTR = IPTR + 16
740 C BYTES 19-20
741 C P2
742 CALL SBYTE (IPFLD,KBDS(19),IPTR,16)
743 IPTR = IPTR + 16
744 C BYTE 21 - RESERVED LOCATION
745 CALL SBYTE (IPFLD,0,IPTR,8)
746 IPTR = IPTR + 8
747 C BYTES 22 - ?
748 C WIDTHS OF SECOND ORDER PACKING
749 IX = (IWDPTR + 32) / 32
750 CALL SBYTES (IPFLD,ISOWID,IPTR,32,0,IX)
751 IPTR = IPTR + IWDPTR
752 C SECONDARY BIT MAP
753 IJ = (IBMP2P + 32) / 32
754 CALL SBYTES (IPFLD,ISOBMP,IPTR,32,0,IJ)
755 IPTR = IPTR + IBMP2P
756 IF (MOD(IPTR,8).NE.0) THEN
757 IPTR = IPTR + 8 - MOD(IPTR,8)
758 END IF
759 C DETERMINE LOCATION FOR START
760 C OF FIRST ORDER PACKED VALUES
761 KBDS(12) = IPTR / 8 + 12
762 C STORE LOCATION
763 CALL SBYTE (IPFLD,KBDS(12),0,16)
764 C MOVE IN FIRST ORDER PACKED VALUES
765 IPASS = (IFOPTR + 32) / 32
766 CALL SBYTES (IPFLD,IFOVAL,IPTR,32,0,IPASS)
767 IPTR = IPTR + IFOPTR
768 IF (MOD(IPTR,8).NE.0) THEN
769 IPTR = IPTR + 8 - MOD(IPTR,8)
770 END IF
771 C PRINT *,'IFOPTR =',IFOPTR,' ISOPTR =',ISOPTR
772 C DETERMINE LOCATION FOR START
773 C OF SECOND ORDER VALUES
774 KBDS(15) = IPTR / 8 + 12
775 C SAVE LOCATION OF SECOND ORDER VALUES
776 CALL SBYTE (IPFLD,KBDS(15),24,16)
777 C MOVE IN SECOND ORDER PACKED VALUES
778 IX = (ISOPTR + 32) / 32
779 CALL SBYTES (IPFLD,ISOVAL,IPTR,32,0,IX)
780 IPTR = IPTR + ISOPTR
781 NLEFT = MOD(IPTR+88,16)
782 IF (NLEFT.NE.0) THEN
783 NLEFT = 16 - NLEFT
784 IPTR = IPTR + NLEFT
785 END IF
786 C COMPUTE LENGTH OF DATA PORTION
787 LEN = IPTR / 8
788 C COMPUTE LENGTH OF BDS
789 LENBDS = LEN + 11
790 C -----------------------------------
791 C BYTES 1-3
792 C THIS FUNCTION COMPLETED BELOW
793 C WHEN LENGTH OF BDS IS KNOWN
794 CALL SBYTE (BDS11,LENBDS,0,24)
795 C BYTE 4
796 CALL SBYTE (BDS11,IBDSFL(1),24,1)
797 CALL SBYTE (BDS11,IBDSFL(2),25,1)
798 CALL SBYTE (BDS11,IBDSFL(3),26,1)
799 CALL SBYTE (BDS11,IBDSFL(4),27,1)
800 C ENTER NUMBER OF FILL BITS
801 CALL SBYTE (BDS11,NLEFT,28,4)
802 C BYTE 5-6
803 IF (ISCAL2.LT.0) THEN
804 CALL SBYTE (BDS11,1,32,1)
805 ISCAL2 = - ISCAL2
806 ELSE
807 CALL SBYTE (BDS11,0,32,1)
808 END IF
809 CALL SBYTE (BDS11,ISCAL2,33,15)
811 C FILL OCTET 7-10 WITH THE REFERENCE VALUE
812 C CONVERT THE FLOATING POINT OF YOUR MACHINE TO IBM370 32 BIT
813 C FLOATING POINT NUMBER
814 C REFERENCE VALUE
815 C FIRST TEST TO SEE IF
816 C ON 32 OR 64 BIT COMPUTER
817 CALL W3FI01(LW)
818 IF (LW.EQ.4) THEN
819 CALL W3FI76 (REFNCE,IEXP,IMANT,32)
820 ELSE
821 CALL W3FI76 (REFNCE,IEXP,IMANT,64)
822 END IF
823 CALL SBYTE (BDS11,IEXP,48,8)
824 CALL SBYTE (BDS11,IMANT,56,24)
826 C BYTE 11
828 CALL SBYTE (BDS11,KBDS(11),80,8)
830 RETURN
832 SUBROUTINE FI7502 (IWORK,ISTART,NPTS,ISAME)
833 C $$ SUBPROGRAM DOCUMENTATION BLOCK
834 C . . . .
835 C SUBPROGRAM: FI7502 SECOND ORDER SAME VALUE COLLECTION
836 C PRGMMR: CAVANAUGH ORG: W/NMC42 DATE: 93-06-23
838 C ABSTRACT: COLLECT SEQUENTIAL SAME VALUES FOR PROCESSING
839 C AS SECOND ORDER VALUE FOR GRIB MESSAGES.
841 C PROGRAM HISTORY LOG:
842 C 93-06-23 CAVANAUGH
843 C 95-10-31 IREDELL REMOVED SAVES AND PRINTS
845 C USAGE: CALL FI7502 (IWORK,ISTART,NPTS,ISAME)
846 C INPUT ARGUMENT LIST:
847 C IWORK - ARRAY CONTAINING SOURCE DATA
848 C ISTART - STARTING LOCATION FOR THIS TEST
849 C NPTS - NUMBER OF POINTS IN IWORK
851 C OUTPUT ARGUMENT LIST: (INCLUDING WORK ARRAYS)
852 C ISAME - NUMBER OF SEQUENTIAL POINTS HAVING THE SAME VALUE
854 C REMARKS: SUBPROGRAM CAN BE CALLED FROM A MULTIPROCESSING ENVIRONMENT.
856 C ATTRIBUTES:
857 C LANGUAGE: IBM VS FORTRAN 77, CRAY CFT77 FORTRAN
858 C MACHINE: HDS, CRAY C916/256, Y-MP8/64, Y-MP EL92/256
860 C $$
861 INTEGER IWORK(*)
862 INTEGER ISTART
863 INTEGER ISAME
864 INTEGER K
865 INTEGER NPTS
866 C -------------------------------------------------------------
867 ISAME = 0
868 DO 100 K = ISTART, NPTS
869 IF (IWORK(K).NE.IWORK(ISTART)) THEN
870 RETURN
871 END IF
872 ISAME = ISAME + 1
873 100 CONTINUE
874 RETURN
876 SUBROUTINE FI7503 (IWORK,IPFLD,NPTS,IBDSFL,BDS11,
877 * LEN,LENBDS,PDS,REFNCE,ISCAL2,KWIDE,IGDS)
878 C $$ SUBPROGRAM DOCUMENTATION BLOCK
879 C . . . .
880 C SUBPROGRAM: FI7501 ROW BY ROW, COL BY COL PACKING
881 C PRGMMR: CAVANAUGH ORG: W/NMC42 DATE: 94-05-20
883 C ABSTRACT: PERFORM ROW BY ROW OR COLUMN BY COLUMN PACKING
884 C GENERATING ALL BDS INFORMATION.
886 C PROGRAM HISTORY LOG:
887 C 93-08-06 CAVANAUGH
888 C 95-10-31 IREDELL REMOVED SAVES AND PRINTS
890 C USAGE: CALL FI7503 (IWORK,IPFLD,NPTS,IBDSFL,BDS11,
891 C * LEN,LENBDS,PDS,REFNCE,ISCAL2,KWIDE,IGDS)
892 C INPUT ARGUMENT LIST:
893 C IWORK - INTEGER SOURCE ARRAY
894 C NPTS - NUMBER OF POINTS IN IWORK
895 C IBDSFL - FLAGS
897 C OUTPUT ARGUMENT LIST: (INCLUDING WORK ARRAYS)
898 C IPFLD - CONTAINS BDS FROM BYTE 12 ON
899 C BDS11 - CONTAINS FIRST 11 BYTES FOR BDS
900 C LEN - NUMBER OF BYTES FROM 12 ON
901 C LENBDS - TOTAL LENGTH OF BDS
903 C REMARKS: SUBPROGRAM CAN BE CALLED FROM A MULTIPROCESSING ENVIRONMENT.
905 C ATTRIBUTES:
906 C LANGUAGE: IBM VS FORTRAN 77, CRAY CFT77 FORTRAN
907 C MACHINE: HDS, CRAY C916/256, Y-MP8/64, Y-MP EL92/256
909 C $$
910 CHARACTER*1 BDS11(*),PDS(*)
912 REAL REFNCE
914 INTEGER ISCAL2,KWIDE
915 INTEGER LENBDS
916 INTEGER IPFLD(*),IGDS(*)
917 INTEGER LEN,KBDS(22)
918 INTEGER IWORK(*)
919 C OCTET NUMBER IN SECTION, FIRST ORDER PACKING
920 C INTEGER KBDS(12)
921 C FLAGS
922 INTEGER IBDSFL(*)
923 C EXTENDED FLAGS
924 C INTEGER KBDS(14)
925 C OCTET NUMBER FOR SECOND ORDER PACKING
926 C INTEGER KBDS(15)
927 C NUMBER OF FIRST ORDER VALUES
928 C INTEGER KBDS(17)
929 C NUMBER OF SECOND ORDER PACKED VALUES
930 C INTEGER KBDS(19)
931 C WIDTH OF SECOND ORDER PACKING
932 INTEGER ISOWID(50000)
933 C SECONDARY BIT MAP
934 INTEGER ISOBMP(8200)
935 C FIRST ORDER PACKED VALUES
936 INTEGER IFOVAL(50000)
937 C SECOND ORDER PACKED VALUES
938 INTEGER ISOVAL(100000)
940 C INTEGER KBDS(11)
941 C ----------------------------------
942 C INITIALIZE ARRAYS
943 DO 100 I = 1, 50000
944 ISOWID(I) = 0
945 IFOVAL(I) = 0
946 100 CONTINUE
948 DO 101 I = 1, 8200
949 ISOBMP(I) = 0
950 101 CONTINUE
951 DO 102 I = 1, 100000
952 ISOVAL(I) = 0
953 102 CONTINUE
954 C INITIALIZE POINTERS
955 C SECONDARY BIT WIDTH POINTER
956 IWDPTR = 0
957 C SECONDARY BIT MAP POINTER
958 IBMP2P = 0
959 C FIRST ORDER VALUE POINTER
960 IFOPTR = 0
961 C BYTE POINTER TO START OF 1ST ORDER VALUES
962 KBDS(12) = 0
963 C BYTE POINTER TO START OF 2ND ORDER VALUES
964 KBDS(15) = 0
965 C TO CONTAIN NUMBER OF FIRST ORDER VALUES
966 KBDS(17) = 0
967 C TO CONTAIN NUMBER OF SECOND ORDER VALUES
968 KBDS(19) = 0
969 C SECOND ORDER PACKED VALUE POINTER
970 ISOPTR = 0
971 C =======================================================
972 C BUILD SECOND ORDER BIT MAP IN EITHER
973 C ROW BY ROW OR COL BY COL FORMAT
974 IF (IAND(IGDS(13),32).NE.0) THEN
975 C COLUMN BY COLUMN
976 KOUT = IGDS(4)
977 KIN = IGDS(5)
978 C PRINT *,'COLUMN BY COLUMN',KOUT,KIN
979 ELSE
980 C ROW BY ROW
981 KOUT = IGDS(5)
982 KIN = IGDS(4)
983 C PRINT *,'ROW BY ROW',KOUT,KIN
984 END IF
985 KBDS(17) = KOUT
986 KBDS(19) = NPTS
988 C DO 4100 J = 1, NPTS, 53
989 C WRITE (6,4101) (IWORK(K),K=J,J+52)
990 4101 FORMAT (1X,25I4)
991 C PRINT *,' '
992 C4100 CONTINUE
994 C INITIALIZE BIT MAP POINTER
995 IBMP2P = 0
996 C CONSTRUCT WORKING BIT MAP
997 DO 2000 I = 1, KOUT
998 DO 1000 J = 1, KIN
999 IF (J.EQ.1) THEN
1000 CALL SBYTE (ISOBMP,1,IBMP2P,1)
1001 ELSE
1002 CALL SBYTE (ISOBMP,0,IBMP2P,1)
1003 END IF
1004 IBMP2P = IBMP2P + 1
1005 1000 CONTINUE
1006 2000 CONTINUE
1007 LEN = IBMP2P / 32 + 1
1008 C CALL BINARY(ISOBMP,LEN)
1010 C PROCESS OUTER LOOP OF ROW BY ROW OR COL BY COL
1012 KPTR = 1
1013 KBDS(11) = KWIDE
1014 DO 6000 I = 1, KOUT
1015 C IN CURRENT ROW OR COL
1016 C FIND FIRST ORDER VALUE
1017 JPTR = KPTR
1018 LOWEST = IWORK(JPTR)
1019 DO 4000 J = 1, KIN
1020 IF (IWORK(JPTR).LT.LOWEST) THEN
1021 LOWEST = IWORK(JPTR)
1022 END IF
1023 JPTR = JPTR + 1
1024 4000 CONTINUE
1025 C SAVE FIRST ORDER VALUE
1026 CALL SBYTE (IFOVAL,LOWEST,IFOPTR,KWIDE)
1027 IFOPTR = IFOPTR + KWIDE
1028 C PRINT *,'FOVAL',I,LOWEST,KWIDE
1029 C SUBTRACT FIRST ORDER VALUE FROM OTHER VALS
1030 C GETTING SECOND ORDER VALUES
1031 JPTR = KPTR
1032 IBIG = IWORK(JPTR) - LOWEST
1033 DO 4200 J = 1, KIN
1034 IWORK(JPTR) = IWORK(JPTR) - LOWEST
1035 IF (IWORK(JPTR).GT.IBIG) THEN
1036 IBIG = IWORK(JPTR)
1037 END IF
1038 JPTR = JPTR + 1
1039 4200 CONTINUE
1040 C HOW MANY BITS TO CONTAIN LARGEST SECOND
1041 C ORDER VALUE IN SEGMENT
1042 CALL FI7505 (IBIG,NWIDE)
1043 C SAVE BIT WIDTH
1044 CALL SBYTE (ISOWID,NWIDE,IWDPTR,8)
1045 IWDPTR = IWDPTR + 8
1046 C PRINT *,I,'SOVAL',IBIG,' IN',NWIDE,' BITS'
1047 C WRITE (6,4101) (IWORK(K),K=KPTR,KPTR+52)
1048 C SAVE SECOND ORDER VALUES OF THIS SEGMENT
1049 DO 5000 J = 0, KIN-1
1050 CALL SBYTE (ISOVAL,IWORK(KPTR+J),ISOPTR,NWIDE)
1051 ISOPTR = ISOPTR + NWIDE
1052 5000 CONTINUE
1053 KPTR = KPTR + KIN
1054 6000 CONTINUE
1055 C =======================================================
1056 C CONCANTENATE ALL FIELDS FOR BDS
1058 C REMAINDER GOES INTO IPFLD
1059 IPTR = 0
1060 C BYTES 12-13
1061 C VALUE FOR N1
1062 C LEAVE SPACE FOR THIS
1063 IPTR = IPTR + 16
1064 C BYTE 14
1065 C EXTENDED FLAGS
1066 CALL SBYTE (IPFLD,IBDSFL(5),IPTR,1)
1067 IPTR = IPTR + 1
1068 CALL SBYTE (IPFLD,IBDSFL(6),IPTR,1)
1069 IPTR = IPTR + 1
1070 CALL SBYTE (IPFLD,IBDSFL(7),IPTR,1)
1071 IPTR = IPTR + 1
1072 CALL SBYTE (IPFLD,IBDSFL(8),IPTR,1)
1073 IPTR = IPTR + 1
1074 CALL SBYTE (IPFLD,IBDSFL(9),IPTR,1)
1075 IPTR = IPTR + 1
1076 CALL SBYTE (IPFLD,IBDSFL(10),IPTR,1)
1077 IPTR = IPTR + 1
1078 CALL SBYTE (IPFLD,IBDSFL(11),IPTR,1)
1079 IPTR = IPTR + 1
1080 CALL SBYTE (IPFLD,IBDSFL(12),IPTR,1)
1081 IPTR = IPTR + 1
1082 C BYTES 15-16
1083 C SKIP OVER VALUE FOR N2
1084 IPTR = IPTR + 16
1085 C BYTES 17-18
1086 C P1
1087 CALL SBYTE (IPFLD,KBDS(17),IPTR,16)
1088 IPTR = IPTR + 16
1089 C BYTES 19-20
1090 C P2
1091 CALL SBYTE (IPFLD,KBDS(19),IPTR,16)
1092 IPTR = IPTR + 16
1093 C BYTE 21 - RESERVED LOCATION
1094 CALL SBYTE (IPFLD,0,IPTR,8)
1095 IPTR = IPTR + 8
1096 C BYTES 22 - ?
1097 C WIDTHS OF SECOND ORDER PACKING
1098 IX = (IWDPTR + 32) / 32
1099 CALL SBYTES (IPFLD,ISOWID,IPTR,32,0,IX)
1100 IPTR = IPTR + IWDPTR
1101 C PRINT *,'ISOWID',IWDPTR,IX
1102 C CALL BINARY (ISOWID,IX)
1104 C NO SECONDARY BIT MAP
1106 C DETERMINE LOCATION FOR START
1107 C OF FIRST ORDER PACKED VALUES
1108 KBDS(12) = IPTR / 8 + 12
1109 C STORE LOCATION
1110 CALL SBYTE (IPFLD,KBDS(12),0,16)
1111 C MOVE IN FIRST ORDER PACKED VALUES
1112 IPASS = (IFOPTR + 32) / 32
1113 CALL SBYTES (IPFLD,IFOVAL,IPTR,32,0,IPASS)
1114 IPTR = IPTR + IFOPTR
1115 C PRINT *,'IFOVAL',IFOPTR,IPASS,KWIDE
1116 C CALL BINARY (IFOVAL,IPASS)
1117 IF (MOD(IPTR,8).NE.0) THEN
1118 IPTR = IPTR + 8 - MOD(IPTR,8)
1119 END IF
1120 C PRINT *,'IFOPTR =',IFOPTR,' ISOPTR =',ISOPTR
1121 C DETERMINE LOCATION FOR START
1122 C OF SECOND ORDER VALUES
1123 KBDS(15) = IPTR / 8 + 12
1124 C SAVE LOCATION OF SECOND ORDER VALUES
1125 CALL SBYTE (IPFLD,KBDS(15),24,16)
1126 C MOVE IN SECOND ORDER PACKED VALUES
1127 IX = (ISOPTR + 32) / 32
1128 CALL SBYTES (IPFLD,ISOVAL,IPTR,32,0,IX)
1129 IPTR = IPTR + ISOPTR
1130 C PRINT *,'ISOVAL',ISOPTR,IX
1131 C CALL BINARY (ISOVAL,IX)
1132 NLEFT = MOD(IPTR+88,16)
1133 IF (NLEFT.NE.0) THEN
1134 NLEFT = 16 - NLEFT
1135 IPTR = IPTR + NLEFT
1136 END IF
1137 C COMPUTE LENGTH OF DATA PORTION
1138 LEN = IPTR / 8
1139 C COMPUTE LENGTH OF BDS
1140 LENBDS = LEN + 11
1141 C -----------------------------------
1142 C BYTES 1-3
1143 C THIS FUNCTION COMPLETED BELOW
1144 C WHEN LENGTH OF BDS IS KNOWN
1145 CALL SBYTE (BDS11,LENBDS,0,24)
1146 C BYTE 4
1147 CALL SBYTE (BDS11,IBDSFL(1),24,1)
1148 CALL SBYTE (BDS11,IBDSFL(2),25,1)
1149 CALL SBYTE (BDS11,IBDSFL(3),26,1)
1150 CALL SBYTE (BDS11,IBDSFL(4),27,1)
1151 C ENTER NUMBER OF FILL BITS
1152 CALL SBYTE (BDS11,NLEFT,28,4)
1153 C BYTE 5-6
1154 IF (ISCAL2.LT.0) THEN
1155 CALL SBYTE (BDS11,1,32,1)
1156 ISCAL2 = - ISCAL2
1157 ELSE
1158 CALL SBYTE (BDS11,0,32,1)
1159 END IF
1160 CALL SBYTE (BDS11,ISCAL2,33,15)
1162 C FILL OCTET 7-10 WITH THE REFERENCE VALUE
1163 C CONVERT THE FLOATING POINT OF YOUR MACHINE TO IBM370 32 BIT
1164 C FLOATING POINT NUMBER
1165 C REFERENCE VALUE
1166 C FIRST TEST TO SEE IF
1167 C ON 32 OR 64 BIT COMPUTER
1168 CALL W3FI01(LW)
1169 IF (LW.EQ.4) THEN
1170 CALL W3FI76 (REFNCE,IEXP,IMANT,32)
1171 ELSE
1172 CALL W3FI76 (REFNCE,IEXP,IMANT,64)
1173 END IF
1174 CALL SBYTE (BDS11,IEXP,48,8)
1175 CALL SBYTE (BDS11,IMANT,56,24)
1177 C BYTE 11
1179 CALL SBYTE (BDS11,KBDS(11),80,8)
1181 KLEN = LENBDS / 4 + 1
1182 C PRINT *,'BDS11 LISTING',4,LENBDS
1183 C CALL BINARY (BDS11,4)
1184 C PRINT *,'IPFLD LISTING'
1185 C CALL BINARY (IPFLD,KLEN)
1186 RETURN
1188 SUBROUTINE FI7505 (N,NBITS)
1189 C $$ SUBPROGRAM DOCUMENTATION BLOCK
1190 C . . . .
1191 C SUBPROGRAM: FI7505 DETERMINE NUMBER OF BITS TO CONTAIN VALUE
1192 C PRGMMR: CAVANAUGH ORG: W/NMC42 DATE: 93-06-23
1194 C ABSTRACT: CALCULATE NUMBER OF BITS TO CONTAIN VALUE N, WITH A
1195 C MAXIMUM OF 32 BITS.
1197 C PROGRAM HISTORY LOG:
1198 C 93-06-23 CAVANAUGH
1199 C 95-10-31 IREDELL REMOVED SAVES AND PRINTS
1201 C USAGE: CALL FI7505 (N,NBITS)
1202 C INPUT ARGUMENT LIST:
1203 C N - INTEGER VALUE
1205 C OUTPUT ARGUMENT LIST: (INCLUDING WORK ARRAYS)
1206 C NBITS - NUMBER OF BITS TO CONTAIN N
1208 C REMARKS: SUBPROGRAM CAN BE CALLED FROM A MULTIPROCESSING ENVIRONMENT.
1210 C ATTRIBUTES:
1211 C LANGUAGE: IBM VS FORTRAN 77, CRAY CFT77 FORTRAN
1212 C MACHINE: HDS, CRAY C916/256, Y-MP8/64, Y-MP EL92/256
1214 C $$
1215 INTEGER N,NBITS
1216 INTEGER IBITS(31)
1218 DATA IBITS/1,3,7,15,31,63,127,255,511,1023,2047,
1219 * 4095,8191,16383,32767,65535,131071,262143,
1220 * 524287,1048575,2097151,4194303,8388607,
1221 * 16777215,33554431,67108863,134217727,268435455,
1222 * 536870911,1073741823,2147483647/
1223 C ----------------------------------------------------------------
1225 DO 1000 NBITS = 1, 31
1226 IF (N.LE.IBITS(NBITS)) THEN
1227 RETURN
1228 END IF
1229 1000 CONTINUE
1230 RETURN
1232 SUBROUTINE FI7513 (IWORK,ISTART,NPTS,MAX,MIN,INRNGE)
1233 C $$ SUBPROGRAM DOCUMENTATION BLOCK
1234 C . . . .
1235 C SUBPROGRAM: FI7513 SELECT BLOCK OF DATA FOR PACKING
1236 C PRGMMR: CAVANAUGH ORG: W/NMC42 DATE: 94-01-21
1238 C ABSTRACT: SELECT A BLOCK OF DATA FOR PACKING
1240 C PROGRAM HISTORY LOG:
1241 C 94-01-21 CAVANAUGH
1242 C 95-10-31 IREDELL REMOVED SAVES AND PRINTS
1244 C USAGE: CALL FI7513 (IWORK,ISTART,NPTS,MAX,MIN,INRNGE)
1245 C INPUT ARGUMENT LIST:
1246 C * - RETURN ADDRESS IF ENCOUNTER SET OF SAME VALUES
1247 C IWORK -
1248 C ISTART -
1249 C NPTS -
1251 C OUTPUT ARGUMENT LIST: (INCLUDING WORK ARRAYS)
1252 C MAX -
1253 C MIN -
1254 C INRNGE -
1256 C REMARKS: SUBPROGRAM CAN BE CALLED FROM A MULTIPROCESSING ENVIRONMENT.
1258 C ATTRIBUTES:
1259 C LANGUAGE: IBM VS FORTRAN 77, CRAY CFT77 FORTRAN
1260 C MACHINE: HDS, CRAY C916/256, Y-MP8/64, Y-MP EL92/256
1262 C $$
1263 INTEGER IWORK(*),NPTS,ISTART,INRNGE,INRNGA,INRNGB
1264 INTEGER MAX,MIN,MXVAL,MAXB,MINB,MXVALB
1265 INTEGER IBITS(31)
1267 DATA IBITS/1,3,7,15,31,63,127,255,511,1023,2047,
1268 * 4095,8191,16383,32767,65535,131071,262143,
1269 * 524287,1048575,2097151,4194303,8388607,
1270 * 16777215,33554431,67108863,134217727,268435455,
1271 * 536870911,1073741823,2147483647/
1272 C ----------------------------------------------------------------
1273 C IDENTIFY NEXT BLOCK OF DATA FOR PACKING AND
1274 C RETURN TO CALLER
1275 C ********************************************************************
1276 ISTRTA = ISTART
1278 C GET BLOCK A
1279 CALL FI7516 (IWORK,NPTS,INRNGA,ISTRTA,
1280 * MAX,MIN,MXVAL,LWIDE)
1281 C ********************************************************************
1283 ISTRTB = ISTRTA + INRNGA
1284 2000 CONTINUE
1285 C IF HAVE PROCESSED ALL DATA, RETURN
1286 IF (ISTRTB.GT.NPTS) THEN
1287 C NO MORE DATA TO LOOK AT
1288 INRNGE = INRNGA
1289 RETURN
1290 END IF
1291 C GET BLOCK B
1292 CALL FI7502 (IWORK,ISTRTB,NPTS,ISAME)
1293 IF (ISAME.GE.15) THEN
1294 C PRINT *,'BLOCK B HAS ALL IDENTICAL VALUES'
1295 C PRINT *,'BLOCK A HAS INRNGE =',INRNGA
1296 C BLOCK B CONTAINS ALL IDENTICAL VALUES
1297 INRNGE = INRNGA
1298 C EXIT WITH BLOCK A
1299 RETURN
1300 END IF
1301 C GET BLOCK B
1303 ISTRTB = ISTRTA + INRNGA
1304 CALL FI7516 (IWORK,NPTS,INRNGB,ISTRTB,
1305 * MAXB,MINB,MXVALB,LWIDEB)
1306 C PRINT *,'BLOCK A',INRNGA,' BLOCK B',INRNGB
1307 C ********************************************************************
1308 C PERFORM TREND ANALYSIS TO DETERMINE
1309 C IF DATA COLLECTION CAN BE IMPROVED
1311 KTRND = LWIDE - LWIDEB
1312 C PRINT *,'TREND',LWIDE,LWIDEB
1313 IF (KTRND.LE.0) THEN
1314 C PRINT *,'BLOCK A - SMALLER, SHOULD EXTEND INTO BLOCK B'
1315 MXVAL = IBITS(LWIDE)
1317 C IF BLOCK A REQUIRES THE SAME OR FEWER BITS
1318 C LOOK AHEAD
1319 C AND GATHER THOSE DATA POINTS THAT CAN
1320 C BE RETAINED IN BLOCK A
1321 C BECAUSE THIS BLOCK OF DATA
1322 C USES FEWER BITS
1324 CALL FI7518 (IRET,IWORK,NPTS,ISTRTA,INRNGA,INRNGB,
1325 * MAX,MIN,LWIDE,MXVAL)
1326 IF(IRET.EQ.1) GO TO 8000
1327 C PRINT *,'18 INRNGA IS NOW ',INRNGA
1328 IF (INRNGB.LT.20) THEN
1329 RETURN
1330 ELSE
1331 GO TO 2000
1332 END IF
1333 ELSE
1334 C PRINT *,'BLOCK A - LARGER, B SHOULD EXTEND BACK INTO A'
1335 MXVALB = IBITS(LWIDEB)
1337 C IF BLOCK B REQUIRES FEWER BITS
1338 C LOOK BACK
1339 C SHORTEN BLOCK A BECAUSE NEXT BLOCK OF DATA
1340 C USES FEWER BITS
1342 CALL FI7517 (IRET,IWORK,NPTS,ISTRTB,INRNGA,
1343 * MAXB,MINB,LWIDEB,MXVALB)
1344 IF(IRET.EQ.1) GO TO 8000
1345 C PRINT *,'17 INRNGA IS NOW ',INRNGA
1346 END IF
1348 C PACK UP BLOCK A
1349 C UPDATA POINTERS
1350 8000 CONTINUE
1351 INRNGE = INRNGA
1352 C GET NEXT BLOCK A
1353 9000 CONTINUE
1354 RETURN
1356 SUBROUTINE FI7516 (IWORK,NPTS,INRNG,ISTART,MAX,MIN,MXVAL,LWIDTH)
1357 C $$ SUBPROGRAM DOCUMENTATION BLOCK
1358 C . . . .
1359 C SUBPROGRAM: FI7516 SCAN NUMBER OF POINTS
1360 C PRGMMR: CAVANAUGH ORG: W/NMC42 DATE: 94-01-21
1362 C ABSTRACT: SCAN FORWARD FROM CURRENT POSITION. COLLECT POINTS AND
1363 C DETERMINE MAXIMUM AND MINIMUM VALUES AND THE NUMBER
1364 C OF POINTS THAT ARE INCLUDED. FORWARD SEARCH IS TERMINATED
1365 C BY ENCOUNTERING A SET OF IDENTICAL VALUES, BY REACHING
1366 C THE NUMBER OF POINTS SELECTED OR BY REACHING THE END
1367 C OF DATA.
1369 C PROGRAM HISTORY LOG:
1370 C 94-01-21 CAVANAUGH
1371 C 95-10-31 IREDELL REMOVED SAVES AND PRINTS
1373 C USAGE: CALL FI7516 (IWORK,NPTS,INRNG,ISTART,MAX,MIN,MXVAL,LWIDTH)
1374 C INPUT ARGUMENT LIST:
1375 C * - RETURN ADDRESS IF ENCOUNTER SET OF SAME VALUES
1376 C IWORK - DATA ARRAY
1377 C NPTS - NUMBER OF POINTS IN DATA ARRAY
1378 C ISTART - STARTING LOCATION IN DATA
1380 C OUTPUT ARGUMENT LIST: (INCLUDING WORK ARRAYS)
1381 C INRNG - NUMBER OF POINTS SELECTED
1382 C MAX - MAXIMUM VALUE OF POINTS
1383 C MIN - MINIMUM VALUE OF POINTS
1384 C MXVAL - MAXIMUM VALUE THAT CAN BE CONTAINED IN LWIDTH BITS
1385 C LWIDTH - NUMBER OF BITS TO CONTAIN MAX DIFF
1387 C REMARKS: SUBPROGRAM CAN BE CALLED FROM A MULTIPROCESSING ENVIRONMENT.
1389 C ATTRIBUTES:
1390 C LANGUAGE: IBM VS FORTRAN 77, CRAY CFT77 FORTRAN
1391 C MACHINE: HDS, CRAY C916/256, Y-MP8/64, Y-MP EL92/256
1393 C $$
1394 INTEGER IWORK(*),NPTS,ISTART,INRNG,MAX,MIN,LWIDTH,MXVAL
1395 INTEGER IBITS(31)
1397 DATA IBITS/1,3,7,15,31,63,127,255,511,1023,2047,
1398 * 4095,8191,16383,32767,65535,131071,262143,
1399 * 524287,1048575,2097151,4194303,8388607,
1400 * 16777215,33554431,67108863,134217727,268435455,
1401 * 536870911,1073741823,2147483647/
1402 C ----------------------------------------------------------------
1404 INRNG = 1
1405 JQ = ISTART + 19
1406 MAX = IWORK(ISTART)
1407 MIN = IWORK(ISTART)
1408 DO 1000 I = ISTART+1, JQ
1409 CALL FI7502 (IWORK,I,NPTS,ISAME)
1410 IF (ISAME.GE.15) THEN
1411 GO TO 5000
1412 END IF
1413 INRNG = INRNG + 1
1414 IF (IWORK(I).GT.MAX) THEN
1415 MAX = IWORK(I)
1416 ELSE IF (IWORK(I).LT.MIN) THEN
1417 MIN = IWORK(I)
1418 END IF
1419 1000 CONTINUE
1420 5000 CONTINUE
1421 KRNG = MAX - MIN
1423 DO 9000 LWIDTH = 1, 31
1424 IF (KRNG.LE.IBITS(LWIDTH)) THEN
1425 C PRINT *,'RETURNED',INRNG,' VALUES'
1426 RETURN
1427 END IF
1428 9000 CONTINUE
1429 RETURN
1431 SUBROUTINE FI7517 (IRET,IWORK,NPTS,ISTRTB,INRNGA,
1432 * MAXB,MINB,MXVALB,LWIDEB)
1433 C $$ SUBPROGRAM DOCUMENTATION BLOCK
1434 C . . . .
1435 C SUBPROGRAM: FI7517 SCAN BACKWARD
1436 C PRGMMR: CAVANAUGH ORG: W/NMC42 DATE: 94-01-21
1438 C ABSTRACT: SCAN BACKWARDS UNTIL A VALUE EXCEEDS RANGE OF GROUP B
1439 C THIS MAY SHORTEN GROUP A
1441 C PROGRAM HISTORY LOG:
1442 C 94-01-21 CAVANAUGH
1443 C 95-10-31 IREDELL REMOVED SAVES AND PRINTS
1444 C 98-06-17 IREDELL REMOVED ALTERNATE RETURN
1446 C USAGE: CALL FI7517 (IRET,IWORK,NPTS,ISTRTB,INRNGA,
1447 C * MAXB,MINB,MXVALB,LWIDEB)
1448 C INPUT ARGUMENT LIST:
1449 C IWORK -
1450 C ISTRTB -
1451 C NPTS -
1452 C INRNGA -
1454 C OUTPUT ARGUMENT LIST: (INCLUDING WORK ARRAYS)
1455 C IRET -
1456 C JLAST -
1457 C MAXB -
1458 C MINB -
1459 C LWIDTH - NUMBER OF BITS TO CONTAIN MAX DIFF
1461 C REMARKS: SUBPROGRAM CAN BE CALLED FROM A MULTIPROCESSING ENVIRONMENT.
1463 C ATTRIBUTES:
1464 C LANGUAGE: IBM VS FORTRAN 77, CRAY CFT77 FORTRAN
1465 C MACHINE: HDS, CRAY C916/256, Y-MP8/64, Y-MP EL92/256
1467 C $$
1468 INTEGER IWORK(*),NPTS,ISTRTB,INRNGA
1469 INTEGER MAXB,MINB,LWIDEB,MXVALB
1470 INTEGER IBITS(31)
1472 DATA IBITS/1,3,7,15,31,63,127,255,511,1023,2047,
1473 * 4095,8191,16383,32767,65535,131071,262143,
1474 * 524287,1048575,2097151,4194303,8388607,
1475 * 16777215,33554431,67108863,134217727,268435455,
1476 * 536870911,1073741823,2147483647/
1477 C ----------------------------------------------------------------
1478 IRET=0
1479 C PRINT *,' FI7517'
1480 NPOS = ISTRTB - 1
1481 ITST = 0
1482 KSET = INRNGA
1484 1000 CONTINUE
1485 C PRINT *,'TRY NPOS',NPOS,IWORK(NPOS),MAXB,MINB
1486 ITST = ITST + 1
1487 IF (ITST.LE.KSET) THEN
1488 IF (IWORK(NPOS).GT.MAXB) THEN
1489 IF ((IWORK(NPOS)-MINB).GT.MXVALB) THEN
1490 C PRINT *,'WENT OUT OF RANGE AT',NPOS
1491 IRET=1
1492 RETURN
1493 ELSE
1494 MAXB = IWORK(NPOS)
1495 END IF
1496 ELSE IF (IWORK(NPOS).LT.MINB) THEN
1497 IF ((MAXB-IWORK(NPOS)).GT.MXVALB) THEN
1498 C PRINT *,'WENT OUT OF RANGE AT',NPOS
1499 IRET=1
1500 RETURN
1501 ELSE
1502 MINB = IWORK(NPOS)
1503 END IF
1504 END IF
1505 INRNGA = INRNGA - 1
1506 NPOS = NPOS - 1
1507 GO TO 1000
1508 END IF
1509 C ----------------------------------------------------------------
1511 9000 CONTINUE
1512 RETURN
1514 SUBROUTINE FI7518 (IRET,IWORK,NPTS,ISTRTA,INRNGA,INRNGB,
1515 * MAXA,MINA,LWIDEA,MXVALA)
1516 C $$ SUBPROGRAM DOCUMENTATION BLOCK
1517 C . . . .
1518 C SUBPROGRAM: FI7518 SCAN FORWARD
1519 C PRGMMR: CAVANAUGH ORG: W/NMC42 DATE: 94-01-21
1521 C ABSTRACT: SCAN FORWARD FROM START OF BLOCK B TOWARDS END OF BLOCK B
1522 C IF NEXT POINT UNDER TEST FORCES A LARGER MAXVALA THEN
1523 C TERMINATE INDICATING LAST POINT TESTED FOR INCLUSION
1524 C INTO BLOCK A.
1526 C PROGRAM HISTORY LOG:
1527 C 94-01-21 CAVANAUGH
1528 C 95-10-31 IREDELL REMOVED SAVES AND PRINTS
1529 C 98-06-17 IREDELL REMOVED ALTERNATE RETURN
1531 C USAGE: CALL FI7518 (IRET,IWORK,NPTS,ISTRTA,INRNGA,INRNGB,
1532 C * MAXA,MINA,LWIDEA,MXVALA)
1533 C INPUT ARGUMENT LIST:
1534 C IFLD -
1535 C JSTART -
1536 C NPTS -
1538 C OUTPUT ARGUMENT LIST: (INCLUDING WORK ARRAYS)
1539 C IRET -
1540 C JLAST -
1541 C MAX -
1542 C MIN -
1543 C LWIDTH - NUMBER OF BITS TO CONTAIN MAX DIFF
1545 C REMARKS: SUBPROGRAM CAN BE CALLED FROM A MULTIPROCESSING ENVIRONMENT.
1547 C ATTRIBUTES:
1548 C LANGUAGE: IBM VS FORTRAN 77, CRAY CFT77 FORTRAN
1549 C MACHINE: HDS, CRAY C916/256, Y-MP8/64, Y-MP EL92/256
1551 C $$
1552 INTEGER IWORK(*),NPTS,ISTRTA,INRNGA
1553 INTEGER MAXA,MINA,LWIDEA,MXVALA
1554 INTEGER IBITS(31)
1556 DATA IBITS/1,3,7,15,31,63,127,255,511,1023,2047,
1557 * 4095,8191,16383,32767,65535,131071,262143,
1558 * 524287,1048575,2097151,4194303,8388607,
1559 * 16777215,33554431,67108863,134217727,268435455,
1560 * 536870911,1073741823,2147483647/
1561 C ----------------------------------------------------------------
1562 IRET=0
1563 C PRINT *,' FI7518'
1564 NPOS = ISTRTA + INRNGA
1565 ITST = 0
1567 1000 CONTINUE
1568 ITST = ITST + 1
1569 IF (ITST.LE.INRNGB) THEN
1570 C PRINT *,'TRY NPOS',NPOS,IWORK(NPOS),MAXA,MINA
1571 IF (IWORK(NPOS).GT.MAXA) THEN
1572 IF ((IWORK(NPOS)-MINA).GT.MXVALA) THEN
1573 C PRINT *,'FI7518A -',ITST,' RANGE EXCEEDS MAX'
1574 IRET=1
1575 RETURN
1576 ELSE
1577 MAXA = IWORK(NPOS)
1578 END IF
1579 ELSE IF (IWORK(NPOS).LT.MINA) THEN
1580 IF ((MAXA-IWORK(NPOS)).GT.MXVALA) THEN
1581 C PRINT *,'FI7518B -',ITST,' RANGE EXCEEDS MAX'
1582 IRET=1
1583 RETURN
1584 ELSE
1585 MINA = IWORK(NPOS)
1586 END IF
1587 END IF
1588 INRNGA = INRNGA + 1
1589 C PRINT *,' ',ITST,INRNGA
1590 NPOS = NPOS +1
1591 GO TO 1000
1592 END IF
1593 C ----------------------------------------------------------------
1594 9000 CONTINUE
1595 RETURN