1629 lines
57 KiB
Fortran
1629 lines
57 KiB
Fortran
SUBROUTINE W3FI75 (IBITL,ITYPE,ITOSS,FLD,IFLD,IBMAP,IBDSFL,
|
|
& NPTS,BDS11,IPFLD,PFLD,LEN,LENBDS,IBERR,PDS,IGDS)
|
|
C$$$ SUBPROGRAM DOCUMENTATION BLOCK
|
|
C . . . .
|
|
C SUBPROGRAM: W3FI75 GRIB PACK DATA AND FORM BDS OCTETS(1-11)
|
|
C PRGMMR: FARLEY ORG: NMC421 DATE:94-11-22
|
|
C
|
|
C ABSTRACT: THIS ROUTINE PACKS A GRIB FIELD AND FORMS OCTETS(1-11)
|
|
C OF THE BINARY DATA SECTION (BDS).
|
|
C
|
|
C PROGRAM HISTORY LOG:
|
|
C 92-07-10 M. FARLEY ORIGINAL AUTHOR
|
|
C 92-10-01 R.E.JONES CORRECTION FOR FIELD OF CONSTANT DATA
|
|
C 92-10-16 R.E.JONES GET RID OF ARRAYS FP AND INT
|
|
C 93-08-06 CAVANAUGH ADDED ROUTINES FI7501, FI7502, FI7503
|
|
C TO ALLOW SECOND ORDER PACKING IN PDS.
|
|
C 93-07-21 STACKPOLE ASSORTED REPAIRS TO GET 2ND DIFF PACK IN
|
|
C 93-10-28 CAVANAUGH COMMENTED OUT NONOPERATIONAL PRINTS AND
|
|
C WRITE STATEMENTS
|
|
C 93-12-15 CAVANAUGH CORRECTED LOCATION OF START OF FIRST ORDER
|
|
C VALUES AND START OF SECOND ORDER VALUES TO
|
|
C REFLECT A BYTE LOCATION IN THE BDS INSTEAD
|
|
C OF AN OFFSET IN SUBROUTINE FI7501.
|
|
C 94-01-27 CAVANAUGH ADDED IGDS AS INPUT ARGUMENT TO THIS ROUTINE
|
|
C AND ADDED PDS AND IGDS ARRAYS TO THE CALL TO
|
|
C W3FI82 TO PROVIDE INFORMATION NEEDED FOR
|
|
C BOUSTROPHEDONIC PROCESSING.
|
|
C 94-05-25 CAVANAUGH SUBROUTINE FI7503 HAS BEEN ADDED TO PROVIDE
|
|
C FOR ROW BY ROW OR COLUMN BY COLUMN SECOND
|
|
C ORDER PACKING. THIS FEATURE CAN BE ACTIVATED
|
|
C BY SETTING IBDSFL(7) TO ZERO.
|
|
C 94-07-08 CAVANAUGH COMMENTED OUT PRINT STATEMENTS USED FOR DEBUG
|
|
C 94-11-22 FARLEY ENLARGED WORK ARRAYS TO HANDLE .5DEGREE GRIDS
|
|
C 95-06-01 R.E.JONES CORRECTION FOR NUMBER OF UNUSED BITS AT END
|
|
C OF SECTION 4, IN BDS BYTE 4, BITS 5-8.
|
|
C 95-10-31 IREDELL REMOVED SAVES AND PRINTS
|
|
C 2001-06-06 GILBERT CHanged gbyte/sbyte calls to refer to
|
|
C Wesley Ebisuzaki's endian independent
|
|
C versions gbytec/sbytec.
|
|
C Use f90 standard routine bit_size to get
|
|
C number of bits in an integer instead of w3fi01.
|
|
C
|
|
C USAGE: CALL W3FI75 (IBITL,ITYPE,ITOSS,FLD,IFLD,IBMAP,IBDSFL,
|
|
C & NPTS,BDS11,IPFLD,PFLD,LEN,LENBDS,IBERR,PDS,IGDS)
|
|
C INPUT ARGUMENT LIST:
|
|
C IBITL - 0, COMPUTER COMPUTES PACKING LENGTH FROM POWER
|
|
C OF 2 THAT BEST FITS THE DATA.
|
|
C 8, 12, ETC. COMPUTER RESCALES DATA TO FIT INTO
|
|
C SET NUMBER OF BITS.
|
|
C ITYPE - 0 = IF INPUT DATA IS FLOATING POINT (FLD)
|
|
C 1 = IF INPUT DATA IS INTEGER (IFLD)
|
|
C ITOSS - 0 = NO BIT MAP IS INCLUDED (DON'T TOSS DATA)
|
|
C 1 = TOSS NULL DATA ACCORDING TO IBMAP
|
|
C FLD - REAL ARRAY OF DATA TO BE PACKED IF ITYPE=0
|
|
C IFLD - INTEGER ARRAY TO BE PACKED IF ITYPE=1
|
|
C IBMAP - BIT MAP SUPPLIED FROM USER
|
|
C IBDSFL - INTEGER ARRAY CONTAINING TABLE 11 FLAG INFO
|
|
C BDS OCTET 4:
|
|
C (1) 0 = GRID POINT DATA
|
|
C 1 = SPHERICAL HARMONIC COEFFICIENTS
|
|
C (2) 0 = SIMPLE PACKING
|
|
C 1 = SECOND ORDER PACKING
|
|
C (3) 0 = ORIGINAL DATA WERE FLOATING POINT VALUES
|
|
C 1 = ORIGINAL DATA WERE INTEGER VALUES
|
|
C (4) 0 = NO ADDITIONAL FLAGS AT OCTET 14
|
|
C 1 = OCTET 14 CONTAINS FLAG BITS 5-12
|
|
C (5) 0 = RESERVED - ALWAYS SET TO 0
|
|
C (6) 0 = SINGLE DATUM AT EACH GRID POINT
|
|
C 1 = MATRIX OF VALUES AT EACH GRID POINT
|
|
C (7) 0 = NO SECONDARY BIT MAPS
|
|
C 1 = SECONDARY BIT MAPS PRESENT
|
|
C (8) 0 = SECOND ORDER VALUES HAVE CONSTANT WIDTH
|
|
C 1 = SECOND ORDER VALUES HAVE DIFFERENT WIDTHS
|
|
C NPTS - NUMBER OF GRIDPOINTS IN ARRAY TO BE PACKED
|
|
C IGDS - ARRAY OF GDS INFORMATION
|
|
C
|
|
C OUTPUT ARGUMENT LIST:
|
|
C BDS11 - FIRST 11 OCTETS OF BDS
|
|
C PFLD - PACKED GRIB FIELD
|
|
C LEN - LENGTH OF PFLD
|
|
C LENBDS - LENGTH OF BDS
|
|
C IBERR - 1, ERROR CONVERTING IEEE F.P. NUMBER TO IBM370 F.P.
|
|
C
|
|
C REMARKS: SUBPROGRAM CAN BE CALLED FROM A MULTIPROCESSING ENVIRONMENT.
|
|
C
|
|
C ATTRIBUTES:
|
|
C LANGUAGE: IBM VS FORTRAN 77, CRAY CFT77 FORTRAN
|
|
C MACHINE: HDS, CRAY C916/256, Y-MP8/64, Y-MP EL92/256
|
|
C
|
|
C$$$
|
|
C
|
|
REAL FLD(*)
|
|
C REAL FWORK(260000)
|
|
C
|
|
C FWORK CAN USE DYNAMIC ALLOCATION OF MEMORY ON CRAY
|
|
C
|
|
REAL FWORK(NPTS)
|
|
REAL RMIN,REFNCE
|
|
C
|
|
character(len=1) IPFLD(*)
|
|
INTEGER IBDSFL(*)
|
|
INTEGER IBMAP(*)
|
|
INTEGER IFLD(*),IGDS(*)
|
|
C INTEGER IWORK(260000)
|
|
C
|
|
C IWORK CAN USE DYNAMIC ALLOCATION OF MEMORY ON CRAY
|
|
C
|
|
INTEGER IWORK(NPTS)
|
|
C
|
|
LOGICAL CONST
|
|
C
|
|
CHARACTER * 1 BDS11(11),PDS(*)
|
|
CHARACTER * 1 PFLD(*)
|
|
C
|
|
C ================================= RCS keyword statements ==========
|
|
CHARACTER*68 RCSKW1,RCSKW2
|
|
DATA RCSKW1,RCSKW2 / '
|
|
.$Source: /fs/hseb/ob72/wfo_rfc/precip_proc/source/gribit/src/RCS/w3fi75.f,v $
|
|
. $', '
|
|
.$Id: w3fi75.f,v 1.1 2006/05/03 13:44:00 gsood Exp $
|
|
. $' /
|
|
C ===================================================================
|
|
C
|
|
C
|
|
C 1.0 PACK THE FIELD.
|
|
C
|
|
C 1.1 TOSS DATA IF BITMAP BEING USED,
|
|
C MOVING 'DATA' TO WORK AREA...
|
|
C
|
|
CONST = .FALSE.
|
|
IBERR = 0
|
|
IW = 0
|
|
C
|
|
IF (ITOSS .EQ. 1) THEN
|
|
IF (ITYPE .EQ. 0) THEN
|
|
DO 110 IT=1,NPTS
|
|
IF (IBMAP(IT) .EQ. 1) THEN
|
|
IW = IW + 1
|
|
FWORK(IW) = FLD(IT)
|
|
ENDIF
|
|
110 CONTINUE
|
|
NPTS = IW
|
|
ELSE IF (ITYPE .EQ. 1) THEN
|
|
DO 111 IT=1,NPTS
|
|
IF (IBMAP(IT) .EQ. 1) THEN
|
|
IW = IW + 1
|
|
IWORK(IW) = IFLD(IT)
|
|
ENDIF
|
|
111 CONTINUE
|
|
NPTS = IW
|
|
ENDIF
|
|
C
|
|
C ELSE, JUST MOVE DATA TO WORK ARRAY
|
|
C
|
|
ELSE IF (ITOSS .EQ. 0) THEN
|
|
IF (ITYPE .EQ. 0) THEN
|
|
DO 112 IT=1,NPTS
|
|
FWORK(IT) = FLD(IT)
|
|
112 CONTINUE
|
|
ELSE IF (ITYPE .EQ. 1) THEN
|
|
DO 113 IT=1,NPTS
|
|
IWORK(IT) = IFLD(IT)
|
|
113 CONTINUE
|
|
ENDIF
|
|
ENDIF
|
|
C
|
|
C 1.2 CONVERT DATA IF NEEDED PRIOR TO PACKING.
|
|
C (INTEGER TO F.P. OR F.P. TO INTEGER)
|
|
C ITYPE = 0...FLOATING POINT DATA
|
|
C IBITL = 0...PACK IN LEAST # BITS...CONVERT TO INTEGER
|
|
C ITYPE = 1...INTEGER DATA
|
|
C IBITL > 0...PACK IN FIXED # BITS...CONVERT TO FLOATING POINT
|
|
C
|
|
IF (ITYPE .EQ. 0 .AND. IBITL .EQ. 0) THEN
|
|
DO 120 IF=1,NPTS
|
|
IWORK(IF) = NINT(FWORK(IF))
|
|
120 CONTINUE
|
|
ELSE IF (ITYPE .EQ. 1 .AND. IBITL .NE. 0) THEN
|
|
DO 123 IF=1,NPTS
|
|
FWORK(IF) = FLOAT(IWORK(IF))
|
|
123 CONTINUE
|
|
ENDIF
|
|
C
|
|
C 1.3 PACK THE DATA.
|
|
C
|
|
IF (IBDSFL(2).NE.0) THEN
|
|
C SECOND ORDER PACKING
|
|
C
|
|
C PRINT*,' DOING SECOND ORDER PACKING...'
|
|
IF (IBITL.EQ.0) THEN
|
|
C
|
|
C PRINT*,' AND VARIABLE BIT PACKING'
|
|
C
|
|
C WORKING WITH INTEGER VALUES
|
|
C SINCE DOING VARIABLE BIT PACKING
|
|
C
|
|
MAX = IWORK(1)
|
|
MIN = IWORK(1)
|
|
DO 300 I = 2, NPTS
|
|
IF (IWORK(I).LT.MIN) THEN
|
|
MIN = IWORK(I)
|
|
ELSE IF (IWORK(I).GT.MAX) THEN
|
|
MAX = IWORK(I)
|
|
END IF
|
|
300 CONTINUE
|
|
C EXTRACT MINIMA
|
|
DO 400 I = 1, NPTS
|
|
C IF (IWORK(I).LT.0) THEN
|
|
C PRINT *,'MINIMA 400',I,IWORK(I),NPTS
|
|
C END IF
|
|
IWORK(I) = IWORK(I) - MIN
|
|
400 CONTINUE
|
|
REFNCE = MIN
|
|
IDIFF = MAX - MIN
|
|
C PRINT *,'REFERENCE VALUE',REFNCE
|
|
C
|
|
C WRITE (6,FMT='('' MINIMA REMOVED = '',/,
|
|
C & 10(3X,10I10,/))') (IWORK(I),I=1,6)
|
|
C WRITE (6,FMT='('' END OF ARRAY = '',/,
|
|
C & 10(3X,10I10,/))') (IWORK(I),I=NPTS-5,NPTS)
|
|
C
|
|
C FIND BIT WIDTH OF IDIFF
|
|
C
|
|
CALL FI7505 (IDIFF,KWIDE)
|
|
C PRINT*,' BIT WIDTH FOR ORIGINAL DATA', KWIDE
|
|
ISCAL2 = 0
|
|
C
|
|
C MULTIPLICATIVE SCALE FACTOR SET TO 1
|
|
C IN ANTICIPATION OF POSSIBLE USE IN GLAHN 2DN DIFF
|
|
C
|
|
SCAL2 = 1.
|
|
C
|
|
ELSE
|
|
C
|
|
C PRINT*,' AND FIXED BIT PACKING, IBITL = ', IBITL
|
|
C FIXED BIT PACKING
|
|
C - LENGTH OF FIELD IN IBITL
|
|
C - MUST BE REAL DATA
|
|
C FLOATING POINT INPUT
|
|
C
|
|
RMAX = FWORK(1)
|
|
RMIN = FWORK(1)
|
|
DO 100 I = 2, NPTS
|
|
IF (FWORK(I).LT.RMIN) THEN
|
|
RMIN = FWORK(I)
|
|
ELSE IF (FWORK(I).GT.RMAX) THEN
|
|
RMAX = FWORK(I)
|
|
END IF
|
|
100 CONTINUE
|
|
REFNCE = RMIN
|
|
C PRINT *,'100 REFERENCE',REFNCE
|
|
C EXTRACT MINIMA
|
|
DO 200 I = 1, NPTS
|
|
FWORK(I) = FWORK(I) - RMIN
|
|
200 CONTINUE
|
|
C PRINT *,'REFERENCE VALUE',REFNCE
|
|
C WRITE (6,FMT='('' MINIMA REMOVED = '',/,
|
|
C & 10(3X,10F8.2,/))') (FWORK(I),I=1,6)
|
|
C WRITE (6,FMT='('' END OF ARRAY = '',/,
|
|
C & 10(3X,10F8.2,/))') (FWORK(I),I=NPTS-5,NPTS)
|
|
C FIND LARGEST DELTA
|
|
IDELT = NINT(RMAX - RMIN)
|
|
C DO BINARY SCALING
|
|
C FIND OUT WHAT BINARY SCALE FACTOR
|
|
C PERMITS CONTAINMENT OF
|
|
C LARGEST DELTA
|
|
CALL FI7505 (IDELT,IWIDE)
|
|
C
|
|
C BINARY SCALING
|
|
C
|
|
ISCAL2 = IWIDE - IBITL
|
|
C PRINT *,'SCALING NEEDED TO FIT =',ISCAL2
|
|
C PRINT*,' RANGE OF = ',IDELT
|
|
C
|
|
C EXPAND DATA WITH BINARY SCALING
|
|
C CONVERT TO INTEGER
|
|
SCAL2 = 2.0**ISCAL2
|
|
SCAL2 = 1./ SCAL2
|
|
DO 600 I = 1, NPTS
|
|
IWORK(I) = NINT(FWORK(I) * SCAL2)
|
|
600 CONTINUE
|
|
KWIDE = IBITL
|
|
END IF
|
|
C
|
|
C *****************************************************************
|
|
C
|
|
C FOLLOWING IS FOR GLAHN SECOND DIFFERENCING
|
|
C NOT STANDARD GRIB
|
|
C
|
|
C TEST FOR SECOND DIFFERENCE PACKING
|
|
C BASED OF SIZE OF PDS - SIZE IN FIRST 3 BYTES
|
|
C
|
|
CALL GBYTEC(PDS,IPDSIZ,0,24)
|
|
IF (IPDSIZ.EQ.50) THEN
|
|
C PRINT*,' DO SECOND DIFFERENCE PACKING '
|
|
C
|
|
C GLAHN PACKING TO 2ND DIFFS
|
|
C
|
|
C WRITE (6,FMT='('' CALL TO W3FI82 WITH = '',/,
|
|
C & 10(3X,10I6,/))') (IWORK(I),I=1,NPTS)
|
|
C
|
|
CALL W3FI82 (IWORK,FVAL1,FDIFF1,NPTS,PDS,IGDS)
|
|
C
|
|
C PRINT *,'GLAHN',FVAL1,FDIFF1
|
|
C WRITE (6,FMT='('' OUT FROM W3FI82 WITH = '',/,
|
|
C & 10(3X,10I6,/))') (IWORK(I),I=1,NPTS)
|
|
C
|
|
C MUST NOW RE-REMOVE THE MINIMUM VALUE
|
|
C OF THE SECOND DIFFERENCES TO ASSURE
|
|
C ALL POSITIVE NUMBERS FOR SECOND ORDER GRIB PACKING
|
|
C
|
|
C ORIGINAL REFERENCE VALUE ADDED TO FIRST POINT
|
|
C VALUE FROM THE 2ND DIFF PACKER TO BE ADDED
|
|
C BACK IN WHEN THE 2ND DIFF VALUES ARE
|
|
C RECONSTRUCTED BACK TO THE BASIC VALUES
|
|
C
|
|
C ALSO, THE REFERENCE VALUE IS
|
|
C POWER-OF-TWO SCALED TO MATCH
|
|
C FVAL1. ALL OF THIS SCALING
|
|
C WILL BE REMOVED AFTER THE
|
|
C GLAHN SECOND DIFFERENCING IS UNDONE.
|
|
C THE SCALING FACTOR NEEDED TO DO THAT
|
|
C IS SAVED IN THE PDS AS A SIGNED POSITIVE
|
|
C TWO BYTE INTEGER
|
|
C
|
|
C THE SCALING FOR THE 2ND DIF PACKED
|
|
C VALUES IS PROPERLY SET TO ZERO
|
|
C
|
|
FVAL1 = FVAL1 + REFNCE*SCAL2
|
|
C FIRST TEST TO SEE IF
|
|
C ON 32 OR 64 BIT COMPUTER
|
|
C CALL W3FI01(LW)
|
|
IF (bit_size(LW).EQ.32) THEN
|
|
CALL W3FI76 (FVAL1,IEXP,IMANT,32)
|
|
ELSE
|
|
CALL W3FI76 (FVAL1,IEXP,IMANT,64)
|
|
END IF
|
|
CALL SBYTEC(PDS,IEXP,320,8)
|
|
CALL SBYTEC(PDS,IMANT,328,24)
|
|
C
|
|
IF (bit_size(LW).EQ.32) THEN
|
|
CALL W3FI76 (FDIFF1,IEXP,IMANT,32)
|
|
ELSE
|
|
CALL W3FI76 (FDIFF1,IEXP,IMANT,64)
|
|
END IF
|
|
CALL SBYTEC(PDS,IEXP,352,8)
|
|
CALL SBYTEC(PDS,IMANT,360,24)
|
|
C
|
|
C TURN ISCAL2 INTO SIGNED POSITIVE INTEGER
|
|
C AND STORE IN TWO BYTES
|
|
C
|
|
IF(ISCAL2.GE.0) THEN
|
|
CALL SBYTEC(PDS,ISCAL2,384,16)
|
|
ELSE
|
|
CALL SBYTEC(PDS,1,384,1)
|
|
ISCAL2 = - ISCAL2
|
|
CALL SBYTEC( PDS,ISCAL2,385,15)
|
|
ENDIF
|
|
C
|
|
MAX = IWORK(1)
|
|
MIN = IWORK(1)
|
|
DO 700 I = 2, NPTS
|
|
IF (IWORK(I).LT.MIN) THEN
|
|
MIN = IWORK(I)
|
|
ELSE IF (IWORK(I).GT.MAX) THEN
|
|
MAX = IWORK(I)
|
|
END IF
|
|
700 CONTINUE
|
|
C EXTRACT MINIMA
|
|
DO 710 I = 1, NPTS
|
|
IWORK(I) = IWORK(I) - MIN
|
|
710 CONTINUE
|
|
REFNCE = MIN
|
|
C PRINT *,'710 REFERENCE',REFNCE
|
|
ISCAL2 = 0
|
|
C
|
|
C AND RESET VALUE OF KWIDE - THE BIT WIDTH
|
|
C FOR THE RANGE OF THE VALUES
|
|
C
|
|
IDIFF = MAX - MIN
|
|
CALL FI7505 (IDIFF,KWIDE)
|
|
C
|
|
C PRINT*,'BIT WIDTH (KWIDE) OF 2ND DIFFS', KWIDE
|
|
C
|
|
C **************************** END OF GLAHN PACKING ************
|
|
ELSE IF (IBDSFL(2).EQ.1.AND.IBDSFL(7).EQ.0) THEN
|
|
C HAVE SECOND ORDER PACKING WITH NO SECOND ORDER
|
|
C BIT MAP. ERGO ROW BY ROW - COL BY COL
|
|
CALL FI7503 (IWORK,IPFLD,NPTS,IBDSFL,BDS11,
|
|
* LEN,LENBDS,PDS,REFNCE,ISCAL2,KWIDE,IGDS)
|
|
RETURN
|
|
END IF
|
|
C WRITE (6,FMT='('' CALL TO FI7501 WITH = '',/,
|
|
C & 10(3X,10I6,/))') (IWORK(I),I=1,NPTS)
|
|
C WRITE (6,FMT='('' END OF ARRAY = '',/,
|
|
C & 10(3X,10I6,/))') (IWORK(I),I=NPTS-5,NPTS)
|
|
C PRINT*,' REFNCE,ISCAL2, KWIDE AT CALL TO FI7501',
|
|
C & REFNCE, ISCAL2,KWIDE
|
|
C
|
|
C SECOND ORDER PACKING
|
|
C
|
|
CALL FI7501 (IWORK,IPFLD,NPTS,IBDSFL,BDS11,
|
|
* LEN,LENBDS,PDS,REFNCE,ISCAL2,KWIDE)
|
|
C
|
|
C BDS COMPLETELY ASSEMBLED IN FI7501 FOR SECOND ORDER
|
|
C PACKING.
|
|
C
|
|
ELSE
|
|
C SIMPLE PACKING
|
|
C
|
|
C PRINT*,' SIMPLE FIRST ORDER PACKING...'
|
|
IF (IBITL.EQ.0) THEN
|
|
C PRINT*,' WITH VARIABLE BIT LENGTH'
|
|
C
|
|
C WITH VARIABLE BIT LENGTH, ADJUSTED
|
|
C TO ACCOMMODATE LARGEST VALUE
|
|
C BINARY SCALING ALWAYS = 0
|
|
C
|
|
CALL W3FI58(IWORK,NPTS,IWORK,PFLD,NBITS,LEN,KMIN)
|
|
RMIN = KMIN
|
|
REFNCE = RMIN
|
|
ISCALE = 0
|
|
C PRINT*,' BIT LENGTH CAME OUT AT ...',NBITS
|
|
C
|
|
C SET CONST .TRUE. IF ALL VALUES ARE THE SAME
|
|
C
|
|
IF (LEN.EQ.0.AND.NBITS.EQ.0) CONST = .TRUE.
|
|
C
|
|
ELSE
|
|
C PRINT*,' FIXED BIT LENGTH, IBITL = ', IBITL
|
|
C
|
|
C FIXED BIT LENGTH PACKING (VARIABLE PRECISION)
|
|
C VALUES SCALED BY POWER OF 2 (ISCALE) TO
|
|
C FIT LARGEST VALUE INTO GIVEN BIT LENGTH (IBITL)
|
|
C
|
|
CALL W3FI59(FWORK,NPTS,IBITL,IWORK,PFLD,ISCALE,LEN,RMIN)
|
|
REFNCE = RMIN
|
|
C PRINT *,' SCALING NEEDED TO FIT IS ...', ISCALE
|
|
NBITS = IBITL
|
|
C
|
|
C SET CONST .TRUE. IF ALL VALUES ARE THE SAME
|
|
C
|
|
IF (LEN.EQ.0) THEN
|
|
CONST = .TRUE.
|
|
NBITS = 0
|
|
END IF
|
|
END IF
|
|
C
|
|
C$ COMPUTE LENGTH OF BDS IN OCTETS
|
|
C
|
|
INUM = NPTS * NBITS + 88
|
|
C PRINT *,'NUMBER OF BITS BEFORE FILL ADDED',INUM
|
|
C
|
|
C NUMBER OF FILL BITS
|
|
NFILL = 0
|
|
NLEFT = MOD(INUM,16)
|
|
IF (NLEFT.NE.0) THEN
|
|
INUM = INUM + 16 - NLEFT
|
|
NFILL = 16 - NLEFT
|
|
END IF
|
|
C PRINT *,'NUMBER OF BITS AFTER FILL ADDED',INUM
|
|
C LENGTH OF BDS IN BYTES
|
|
LENBDS = INUM / 8
|
|
C
|
|
C 2.0 FORM THE BINARY DATA SECTION (BDS).
|
|
C
|
|
C CONCANTENATE ALL FIELDS FOR BDS
|
|
C
|
|
C BYTES 1-3
|
|
CALL SBYTEC (BDS11,LENBDS,0,24)
|
|
C
|
|
C BYTE 4
|
|
C FLAGS
|
|
CALL SBYTEC (BDS11,IBDSFL(1),24,1)
|
|
CALL SBYTEC (BDS11,IBDSFL(2),25,1)
|
|
CALL SBYTEC (BDS11,IBDSFL(3),26,1)
|
|
CALL SBYTEC (BDS11,IBDSFL(4),27,1)
|
|
C NR OF FILL BITS
|
|
CALL SBYTEC (BDS11,NFILL,28,4)
|
|
C
|
|
C$ FILL OCTETS 5-6 WITH THE SCALE FACTOR.
|
|
C
|
|
C BYTE 5-6
|
|
IF (ISCALE.LT.0) THEN
|
|
CALL SBYTEC (BDS11,1,32,1)
|
|
ISCALE = - ISCALE
|
|
CALL SBYTEC (BDS11,ISCALE,33,15)
|
|
ELSE
|
|
CALL SBYTEC (BDS11,ISCALE,32,16)
|
|
END IF
|
|
C
|
|
C$ FILL OCTET 7-10 WITH THE REFERENCE VALUE
|
|
C CONVERT THE FLOATING POINT OF YOUR MACHINE TO IBM370 32 BIT
|
|
C FLOATING POINT NUMBER
|
|
C
|
|
C BYTE 7-10
|
|
C REFERENCE VALUE
|
|
C FIRST TEST TO SEE IF
|
|
C ON 32 OR 64 BIT COMPUTER
|
|
C CALL W3FI01(LW)
|
|
IF (bit_size(LW).EQ.32) THEN
|
|
CALL W3FI76 (REFNCE,IEXP,IMANT,32)
|
|
ELSE
|
|
CALL W3FI76 (REFNCE,IEXP,IMANT,64)
|
|
END IF
|
|
CALL SBYTEC (BDS11,IEXP,48,8)
|
|
CALL SBYTEC (BDS11,IMANT,56,24)
|
|
C
|
|
C
|
|
C$ FILL OCTET 11 WITH THE NUMBER OF BITS.
|
|
C
|
|
C BYTE 11
|
|
CALL SBYTEC (BDS11,NBITS,80,8)
|
|
END IF
|
|
C
|
|
RETURN
|
|
END
|
|
SUBROUTINE FI7501 (IWORK,IPFLD,NPTS,IBDSFL,BDS11,
|
|
* LEN,LENBDS,PDS,REFNCE,ISCAL2,KWIDE)
|
|
C$$$ SUBPROGRAM DOCUMENTATION BLOCK
|
|
C . . . .
|
|
C SUBPROGRAM: FI7501 BDS SECOND ORDER PACKING
|
|
C PRGMMR: CAVANAUGH ORG: W/NMC42 DATE: 93-08-06
|
|
C
|
|
C ABSTRACT: PERFORM SECONDARY PACKING ON GRID POINT DATA,
|
|
C GENERATING ALL BDS INFORMATION.
|
|
C
|
|
C PROGRAM HISTORY LOG:
|
|
C 93-08-06 CAVANAUGH
|
|
C 93-12-15 CAVANAUGH CORRECTED LOCATION OF START OF FIRST ORDER
|
|
C VALUES AND START OF SECOND ORDER VALUES TO
|
|
C REFLECT A BYTE LOCATION IN THE BDS INSTEAD
|
|
C OF AN OFFSET.
|
|
C 95-10-31 IREDELL REMOVED SAVES AND PRINTS
|
|
C
|
|
C USAGE: CALL FI7501 (IWORK,IPFLD,NPTS,IBDSFL,BDS11,
|
|
C * LEN,LENBDS,PDS,REFNCE,ISCAL2,KWIDE)
|
|
C INPUT ARGUMENT LIST:
|
|
C IWORK - INTEGER SOURCE ARRAY
|
|
C NPTS - NUMBER OF POINTS IN IWORK
|
|
C IBDSFL - FLAGS
|
|
C
|
|
C OUTPUT ARGUMENT LIST: (INCLUDING WORK ARRAYS)
|
|
C IPFLD - CONTAINS BDS FROM BYTE 12 ON
|
|
C BDS11 - CONTAINS FIRST 11 BYTES FOR BDS
|
|
C LEN - NUMBER OF BYTES FROM 12 ON
|
|
C LENBDS - TOTAL LENGTH OF BDS
|
|
C
|
|
C REMARKS: SUBPROGRAM CAN BE CALLED FROM A MULTIPROCESSING ENVIRONMENT.
|
|
C
|
|
C ATTRIBUTES:
|
|
C LANGUAGE: IBM VS FORTRAN 77, CRAY CFT77 FORTRAN
|
|
C MACHINE: HDS, CRAY C916/256, Y-MP8/64, Y-MP EL92/256
|
|
C
|
|
C$$$
|
|
CHARACTER*1 BDS11(*),PDS(*)
|
|
C
|
|
REAL REFNCE
|
|
C
|
|
INTEGER ISCAL2,KWIDE
|
|
INTEGER LENBDS
|
|
CHARACTER(len=1) IPFLD(*)
|
|
INTEGER LEN,KBDS(22)
|
|
INTEGER IWORK(*)
|
|
C OCTET NUMBER IN SECTION, FIRST ORDER PACKING
|
|
C INTEGER KBDS(12)
|
|
C FLAGS
|
|
INTEGER IBDSFL(*)
|
|
C EXTENDED FLAGS
|
|
C INTEGER KBDS(14)
|
|
C OCTET NUMBER FOR SECOND ORDER PACKING
|
|
C INTEGER KBDS(15)
|
|
C NUMBER OF FIRST ORDER VALUES
|
|
C INTEGER KBDS(17)
|
|
C NUMBER OF SECOND ORDER PACKED VALUES
|
|
C INTEGER KBDS(19)
|
|
C WIDTH OF SECOND ORDER PACKING
|
|
character(len=1) ISOWID(400000)
|
|
C SECONDARY BIT MAP
|
|
character(len=1) ISOBMP(65600)
|
|
C FIRST ORDER PACKED VALUES
|
|
character(len=1) IFOVAL(400000)
|
|
C SECOND ORDER PACKED VALUES
|
|
character(len=1) ISOVAL(800000)
|
|
C
|
|
C INTEGER KBDS(11)
|
|
C BIT WIDTH TABLE
|
|
INTEGER IBITS(31)
|
|
C
|
|
DATA IBITS/1,3,7,15,31,63,127,255,511,1023,
|
|
* 2047,4095,8191,16383,32767,65535,131072,
|
|
* 262143,524287,1048575,2097151,4194303,
|
|
* 8388607,16777215,33554431,67108863,
|
|
* 134217727,268435455,536870911,
|
|
* 1073741823,2147483647/
|
|
C ----------------------------------
|
|
C INITIALIZE ARRAYS
|
|
|
|
DO I = 1, 400000
|
|
IFOVAL(I) = char(0)
|
|
ISOWID(I) = char(0)
|
|
ENDDO
|
|
C
|
|
DO 101 I = 1, 65600
|
|
ISOBMP(I) = char(0)
|
|
101 CONTINUE
|
|
DO 102 I = 1, 800000
|
|
ISOVAL(I) = char(0)
|
|
102 CONTINUE
|
|
C INITIALIZE POINTERS
|
|
C SECONDARY BIT WIDTH POINTER
|
|
IWDPTR = 0
|
|
C SECONDARY BIT MAP POINTER
|
|
IBMP2P = 0
|
|
C FIRST ORDER VALUE POINTER
|
|
IFOPTR = 0
|
|
C BYTE POINTER TO START OF 1ST ORDER VALUES
|
|
KBDS(12) = 0
|
|
C BYTE POINTER TO START OF 2ND ORDER VALUES
|
|
KBDS(15) = 0
|
|
C TO CONTAIN NUMBER OF FIRST ORDER VALUES
|
|
KBDS(17) = 0
|
|
C TO CONTAIN NUMBER OF SECOND ORDER VALUES
|
|
KBDS(19) = 0
|
|
C SECOND ORDER PACKED VALUE POINTER
|
|
ISOPTR = 0
|
|
C =======================================================
|
|
C
|
|
C DATA IS IN IWORK
|
|
C
|
|
KBDS(11) = KWIDE
|
|
C
|
|
C DATA PACKING
|
|
C
|
|
ITER = 0
|
|
INEXT = 1
|
|
ISTART = 1
|
|
C -----------------------------------------------------------
|
|
KOUNT = 0
|
|
C DO 1 I = 1, NPTS, 10
|
|
C PRINT *,I,(IWORK(K),K=I, I+9)
|
|
C 1 CONTINUE
|
|
2000 CONTINUE
|
|
ITER = ITER + 1
|
|
C PRINT *,'NEXT ITERATION STARTS AT',ISTART
|
|
IF (ISTART.GT.NPTS) THEN
|
|
GO TO 4000
|
|
ELSE IF (ISTART.EQ.NPTS) THEN
|
|
KPTS = 1
|
|
MXDIFF = 0
|
|
GO TO 2200
|
|
END IF
|
|
C
|
|
C LOOK FOR REPITITIONS OF A SINGLE VALUE
|
|
CALL FI7502 (IWORK,ISTART,NPTS,ISAME)
|
|
IF (ISAME.GE.15) THEN
|
|
KOUNT = KOUNT + 1
|
|
C PRINT *,'FI7501 - FOUND IDENTICAL SET OF ',ISAME
|
|
MXDIFF = 0
|
|
KPTS = ISAME
|
|
ELSE
|
|
C
|
|
C LOOK FOR SETS OF VALUES IN TREND SELECTED RANGE
|
|
CALL FI7513 (IWORK,ISTART,NPTS,NMAX,NMIN,INRNGE)
|
|
C PRINT *,'ISTART ',ISTART,' INRNGE',INRNGE,NMAX,NMIN
|
|
IEND = ISTART + INRNGE - 1
|
|
C DO 2199 NM = ISTART, IEND, 10
|
|
C PRINT *,' ',(IWORK(NM+JK),JK=0,9)
|
|
C2199 CONTINUE
|
|
MXDIFF = NMAX - NMIN
|
|
KPTS = INRNGE
|
|
END IF
|
|
2200 CONTINUE
|
|
C PRINT *,' RANGE ',MXDIFF,' MAX',NMAX,' MIN',NMIN
|
|
C INCREMENT NUMBER OF FIRST ORDER VALUES
|
|
KBDS(17) = KBDS(17) + 1
|
|
C ENTER FIRST ORDER VALUE
|
|
IF (MXDIFF.GT.0) THEN
|
|
DO 2220 LK = 0, KPTS-1
|
|
IWORK(ISTART+LK) = IWORK(ISTART+LK) - NMIN
|
|
2220 CONTINUE
|
|
CALL SBYTEC (IFOVAL,NMIN,IFOPTR,KBDS(11))
|
|
ELSE
|
|
CALL SBYTEC (IFOVAL,IWORK(ISTART),IFOPTR,KBDS(11))
|
|
END IF
|
|
IFOPTR = IFOPTR + KBDS(11)
|
|
C PROCESS SECOND ORDER BIT WIDTH
|
|
IF (MXDIFF.GT.0) THEN
|
|
DO 2330 KWIDE = 1, 31
|
|
IF (MXDIFF.LE.IBITS(KWIDE)) THEN
|
|
GO TO 2331
|
|
END IF
|
|
2330 CONTINUE
|
|
2331 CONTINUE
|
|
ELSE
|
|
KWIDE = 0
|
|
END IF
|
|
CALL SBYTEC (ISOWID,KWIDE,IWDPTR,8)
|
|
IWDPTR = IWDPTR + 8
|
|
C PRINT *,KWIDE,' IFOVAL=',NMIN,IWORK(ISTART),KPTS
|
|
C IF KWIDE NE 0, SAVE SECOND ORDER VALUE
|
|
IF (KWIDE.GT.0) THEN
|
|
CALL SBYTESC (ISOVAL,IWORK(ISTART),ISOPTR,KWIDE,0,KPTS)
|
|
ISOPTR = ISOPTR + KPTS * KWIDE
|
|
KBDS(19) = KBDS(19) + KPTS
|
|
C PRINT *,' SECOND ORDER VALUES'
|
|
C PRINT *,(IWORK(ISTART+I),I=0,KPTS-1)
|
|
END IF
|
|
C ADD TO SECOND ORDER BITMAP
|
|
CALL SBYTEC (ISOBMP,1,IBMP2P,1)
|
|
IBMP2P = IBMP2P + KPTS
|
|
ISTART = ISTART + KPTS
|
|
GO TO 2000
|
|
C --------------------------------------------------------------
|
|
4000 CONTINUE
|
|
C PRINT *,'THERE WERE ',ITER,' SECOND ORDER GROUPS'
|
|
C PRINT *,'THERE WERE ',KOUNT,' STRINGS OF CONSTANTS'
|
|
C CONCANTENATE ALL FIELDS FOR BDS
|
|
C
|
|
C REMAINDER GOES INTO IPFLD
|
|
IPTR = 0
|
|
C BYTES 12-13
|
|
C VALUE FOR N1
|
|
C LEAVE SPACE FOR THIS
|
|
IPTR = IPTR + 16
|
|
C BYTE 14
|
|
C EXTENDED FLAGS
|
|
CALL SBYTEC (IPFLD,IBDSFL(5),IPTR,1)
|
|
IPTR = IPTR + 1
|
|
CALL SBYTEC (IPFLD,IBDSFL(6),IPTR,1)
|
|
IPTR = IPTR + 1
|
|
CALL SBYTEC (IPFLD,IBDSFL(7),IPTR,1)
|
|
IPTR = IPTR + 1
|
|
CALL SBYTEC (IPFLD,IBDSFL(8),IPTR,1)
|
|
IPTR = IPTR + 1
|
|
CALL SBYTEC (IPFLD,IBDSFL(9),IPTR,1)
|
|
IPTR = IPTR + 1
|
|
CALL SBYTEC (IPFLD,IBDSFL(10),IPTR,1)
|
|
IPTR = IPTR + 1
|
|
CALL SBYTEC (IPFLD,IBDSFL(11),IPTR,1)
|
|
IPTR = IPTR + 1
|
|
CALL SBYTEC (IPFLD,IBDSFL(12),IPTR,1)
|
|
IPTR = IPTR + 1
|
|
C BYTES 15-16
|
|
C SKIP OVER VALUE FOR N2
|
|
IPTR = IPTR + 16
|
|
C BYTES 17-18
|
|
C P1
|
|
CALL SBYTEC (IPFLD,KBDS(17),IPTR,16)
|
|
IPTR = IPTR + 16
|
|
C BYTES 19-20
|
|
C P2
|
|
CALL SBYTEC (IPFLD,KBDS(19),IPTR,16)
|
|
IPTR = IPTR + 16
|
|
C BYTE 21 - RESERVED LOCATION
|
|
CALL SBYTEC (IPFLD,0,IPTR,8)
|
|
IPTR = IPTR + 8
|
|
C BYTES 22 - ?
|
|
C WIDTHS OF SECOND ORDER PACKING
|
|
IX = (IWDPTR + 32) / 32
|
|
C CALL SBYTESC (IPFLD,ISOWID,IPTR,32,0,IX)
|
|
ijk=IWDPTR/8
|
|
jst=(iptr/8)+1
|
|
ipfld(jst:jst+ijk)=ISOWID(1:ijk)
|
|
IPTR = IPTR + IWDPTR
|
|
C SECONDARY BIT MAP
|
|
IJ = (IBMP2P + 32) / 32
|
|
C CALL SBYTESC (IPFLD,ISOBMP,IPTR,32,0,IJ)
|
|
ijk=(IBMP2P/8)+1
|
|
jst=(iptr/8)+1
|
|
ipfld(jst:jst+ijk)=ISOBMP(1:ijk)
|
|
IPTR = IPTR + IBMP2P
|
|
IF (MOD(IPTR,8).NE.0) THEN
|
|
IPTR = IPTR + 8 - MOD(IPTR,8)
|
|
END IF
|
|
C DETERMINE LOCATION FOR START
|
|
C OF FIRST ORDER PACKED VALUES
|
|
KBDS(12) = IPTR / 8 + 12
|
|
C STORE LOCATION
|
|
CALL SBYTEC (IPFLD,KBDS(12),0,16)
|
|
C MOVE IN FIRST ORDER PACKED VALUES
|
|
IPASS = (IFOPTR + 32) / 32
|
|
C CALL SBYTESC (IPFLD,IFOVAL,IPTR,32,0,IPASS)
|
|
ijk=(IFOPTR/8)+1
|
|
jst=(iptr/8)+1
|
|
ipfld(jst:jst+ijk)=ifoval(1:ijk)
|
|
IPTR = IPTR + IFOPTR
|
|
IF (MOD(IPTR,8).NE.0) THEN
|
|
IPTR = IPTR + 8 - MOD(IPTR,8)
|
|
END IF
|
|
C PRINT *,'IFOPTR =',IFOPTR,' ISOPTR =',ISOPTR
|
|
C DETERMINE LOCATION FOR START
|
|
C OF SECOND ORDER VALUES
|
|
KBDS(15) = IPTR / 8 + 12
|
|
C SAVE LOCATION OF SECOND ORDER VALUES
|
|
CALL SBYTEC (IPFLD,KBDS(15),24,16)
|
|
C MOVE IN SECOND ORDER PACKED VALUES
|
|
IX = (ISOPTR + 32) / 32
|
|
c CALL SBYTESC (IPFLD,ISOVAL,IPTR,32,0,IX)
|
|
ijk=(ISOPTR/8)+1
|
|
jst=(iptr/8)+1
|
|
ipfld(jst:jst+ijk)=isoval(1:ijk)
|
|
IPTR = IPTR + ISOPTR
|
|
NLEFT = MOD(IPTR+88,16)
|
|
IF (NLEFT.NE.0) THEN
|
|
NLEFT = 16 - NLEFT
|
|
IPTR = IPTR + NLEFT
|
|
END IF
|
|
C COMPUTE LENGTH OF DATA PORTION
|
|
LEN = IPTR / 8
|
|
C COMPUTE LENGTH OF BDS
|
|
LENBDS = LEN + 11
|
|
C -----------------------------------
|
|
C BYTES 1-3
|
|
C THIS FUNCTION COMPLETED BELOW
|
|
C WHEN LENGTH OF BDS IS KNOWN
|
|
CALL SBYTEC (BDS11,LENBDS,0,24)
|
|
C BYTE 4
|
|
CALL SBYTEC (BDS11,IBDSFL(1),24,1)
|
|
CALL SBYTEC (BDS11,IBDSFL(2),25,1)
|
|
CALL SBYTEC (BDS11,IBDSFL(3),26,1)
|
|
CALL SBYTEC (BDS11,IBDSFL(4),27,1)
|
|
C ENTER NUMBER OF FILL BITS
|
|
CALL SBYTEC (BDS11,NLEFT,28,4)
|
|
C BYTE 5-6
|
|
IF (ISCAL2.LT.0) THEN
|
|
CALL SBYTEC (BDS11,1,32,1)
|
|
ISCAL2 = - ISCAL2
|
|
ELSE
|
|
CALL SBYTEC (BDS11,0,32,1)
|
|
END IF
|
|
CALL SBYTEC (BDS11,ISCAL2,33,15)
|
|
C
|
|
C$ FILL OCTET 7-10 WITH THE REFERENCE VALUE
|
|
C CONVERT THE FLOATING POINT OF YOUR MACHINE TO IBM370 32 BIT
|
|
C FLOATING POINT NUMBER
|
|
C REFERENCE VALUE
|
|
C FIRST TEST TO SEE IF
|
|
C ON 32 OR 64 BIT COMPUTER
|
|
C CALL W3FI01(LW)
|
|
IF (bit_size(LW).EQ.32) THEN
|
|
CALL W3FI76 (REFNCE,IEXP,IMANT,32)
|
|
ELSE
|
|
CALL W3FI76 (REFNCE,IEXP,IMANT,64)
|
|
END IF
|
|
CALL SBYTEC (BDS11,IEXP,48,8)
|
|
CALL SBYTEC (BDS11,IMANT,56,24)
|
|
C
|
|
C BYTE 11
|
|
C
|
|
CALL SBYTEC (BDS11,KBDS(11),80,8)
|
|
C
|
|
RETURN
|
|
END
|
|
SUBROUTINE FI7502 (IWORK,ISTART,NPTS,ISAME)
|
|
C$$$ SUBPROGRAM DOCUMENTATION BLOCK
|
|
C . . . .
|
|
C SUBPROGRAM: FI7502 SECOND ORDER SAME VALUE COLLECTION
|
|
C PRGMMR: CAVANAUGH ORG: W/NMC42 DATE: 93-06-23
|
|
C
|
|
C ABSTRACT: COLLECT SEQUENTIAL SAME VALUES FOR PROCESSING
|
|
C AS SECOND ORDER VALUE FOR GRIB MESSAGES.
|
|
C
|
|
C PROGRAM HISTORY LOG:
|
|
C 93-06-23 CAVANAUGH
|
|
C 95-10-31 IREDELL REMOVED SAVES AND PRINTS
|
|
C
|
|
C USAGE: CALL FI7502 (IWORK,ISTART,NPTS,ISAME)
|
|
C INPUT ARGUMENT LIST:
|
|
C IWORK - ARRAY CONTAINING SOURCE DATA
|
|
C ISTART - STARTING LOCATION FOR THIS TEST
|
|
C NPTS - NUMBER OF POINTS IN IWORK
|
|
C
|
|
C OUTPUT ARGUMENT LIST: (INCLUDING WORK ARRAYS)
|
|
C ISAME - NUMBER OF SEQUENTIAL POINTS HAVING THE SAME VALUE
|
|
C
|
|
C REMARKS: SUBPROGRAM CAN BE CALLED FROM A MULTIPROCESSING ENVIRONMENT.
|
|
C
|
|
C ATTRIBUTES:
|
|
C LANGUAGE: IBM VS FORTRAN 77, CRAY CFT77 FORTRAN
|
|
C MACHINE: HDS, CRAY C916/256, Y-MP8/64, Y-MP EL92/256
|
|
C
|
|
C$$$
|
|
INTEGER IWORK(*)
|
|
INTEGER ISTART
|
|
INTEGER ISAME
|
|
INTEGER K
|
|
INTEGER NPTS
|
|
C -------------------------------------------------------------
|
|
ISAME = 0
|
|
DO 100 K = ISTART, NPTS
|
|
IF (IWORK(K).NE.IWORK(ISTART)) THEN
|
|
RETURN
|
|
END IF
|
|
ISAME = ISAME + 1
|
|
100 CONTINUE
|
|
RETURN
|
|
END
|
|
SUBROUTINE FI7503 (IWORK,IPFLD,NPTS,IBDSFL,BDS11,
|
|
* LEN,LENBDS,PDS,REFNCE,ISCAL2,KWIDE,IGDS)
|
|
C$$$ SUBPROGRAM DOCUMENTATION BLOCK
|
|
C . . . .
|
|
C SUBPROGRAM: FI7501 ROW BY ROW, COL BY COL PACKING
|
|
C PRGMMR: CAVANAUGH ORG: W/NMC42 DATE: 94-05-20
|
|
C
|
|
C ABSTRACT: PERFORM ROW BY ROW OR COLUMN BY COLUMN PACKING
|
|
C GENERATING ALL BDS INFORMATION.
|
|
C
|
|
C PROGRAM HISTORY LOG:
|
|
C 93-08-06 CAVANAUGH
|
|
C 95-10-31 IREDELL REMOVED SAVES AND PRINTS
|
|
C
|
|
C USAGE: CALL FI7503 (IWORK,IPFLD,NPTS,IBDSFL,BDS11,
|
|
C * LEN,LENBDS,PDS,REFNCE,ISCAL2,KWIDE,IGDS)
|
|
C INPUT ARGUMENT LIST:
|
|
C IWORK - INTEGER SOURCE ARRAY
|
|
C NPTS - NUMBER OF POINTS IN IWORK
|
|
C IBDSFL - FLAGS
|
|
C
|
|
C OUTPUT ARGUMENT LIST: (INCLUDING WORK ARRAYS)
|
|
C IPFLD - CONTAINS BDS FROM BYTE 12 ON
|
|
C BDS11 - CONTAINS FIRST 11 BYTES FOR BDS
|
|
C LEN - NUMBER OF BYTES FROM 12 ON
|
|
C LENBDS - TOTAL LENGTH OF BDS
|
|
C
|
|
C REMARKS: SUBPROGRAM CAN BE CALLED FROM A MULTIPROCESSING ENVIRONMENT.
|
|
C
|
|
C ATTRIBUTES:
|
|
C LANGUAGE: IBM VS FORTRAN 77, CRAY CFT77 FORTRAN
|
|
C MACHINE: HDS, CRAY C916/256, Y-MP8/64, Y-MP EL92/256
|
|
C
|
|
C$$$
|
|
CHARACTER*1 BDS11(*),PDS(*),IPFLD(*)
|
|
C
|
|
REAL REFNCE
|
|
C
|
|
INTEGER ISCAL2,KWIDE
|
|
INTEGER LENBDS
|
|
INTEGER IGDS(*)
|
|
INTEGER LEN,KBDS(22)
|
|
INTEGER IWORK(*)
|
|
C OCTET NUMBER IN SECTION, FIRST ORDER PACKING
|
|
C INTEGER KBDS(12)
|
|
C FLAGS
|
|
INTEGER IBDSFL(*)
|
|
C EXTENDED FLAGS
|
|
C INTEGER KBDS(14)
|
|
C OCTET NUMBER FOR SECOND ORDER PACKING
|
|
C INTEGER KBDS(15)
|
|
C NUMBER OF FIRST ORDER VALUES
|
|
C INTEGER KBDS(17)
|
|
C NUMBER OF SECOND ORDER PACKED VALUES
|
|
C INTEGER KBDS(19)
|
|
C WIDTH OF SECOND ORDER PACKING
|
|
character(len=1) ISOWID(400000)
|
|
C SECONDARY BIT MAP
|
|
character(len=1) ISOBMP(65600)
|
|
C FIRST ORDER PACKED VALUES
|
|
character(len=1) IFOVAL(400000)
|
|
C SECOND ORDER PACKED VALUES
|
|
character(len=1) ISOVAL(800000)
|
|
C
|
|
C INTEGER KBDS(11)
|
|
C ----------------------------------
|
|
C INITIALIZE ARRAYS
|
|
C
|
|
DO I = 1, 400000
|
|
IFOVAL(I) = char(0)
|
|
ISOWID(I) = char(0)
|
|
ENDDO
|
|
C
|
|
DO 101 I = 1, 65600
|
|
ISOBMP(I) = char(0)
|
|
101 CONTINUE
|
|
DO 102 I = 1, 800000
|
|
ISOVAL(I) = char(0)
|
|
102 CONTINUE
|
|
C INITIALIZE POINTERS
|
|
C SECONDARY BIT WIDTH POINTER
|
|
IWDPTR = 0
|
|
C SECONDARY BIT MAP POINTER
|
|
IBMP2P = 0
|
|
C FIRST ORDER VALUE POINTER
|
|
IFOPTR = 0
|
|
C BYTE POINTER TO START OF 1ST ORDER VALUES
|
|
KBDS(12) = 0
|
|
C BYTE POINTER TO START OF 2ND ORDER VALUES
|
|
KBDS(15) = 0
|
|
C TO CONTAIN NUMBER OF FIRST ORDER VALUES
|
|
KBDS(17) = 0
|
|
C TO CONTAIN NUMBER OF SECOND ORDER VALUES
|
|
KBDS(19) = 0
|
|
C SECOND ORDER PACKED VALUE POINTER
|
|
ISOPTR = 0
|
|
C =======================================================
|
|
C BUILD SECOND ORDER BIT MAP IN EITHER
|
|
C ROW BY ROW OR COL BY COL FORMAT
|
|
IF (IAND(IGDS(13),32).NE.0) THEN
|
|
C COLUMN BY COLUMN
|
|
KOUT = IGDS(4)
|
|
KIN = IGDS(5)
|
|
C PRINT *,'COLUMN BY COLUMN',KOUT,KIN
|
|
ELSE
|
|
C ROW BY ROW
|
|
KOUT = IGDS(5)
|
|
KIN = IGDS(4)
|
|
C PRINT *,'ROW BY ROW',KOUT,KIN
|
|
END IF
|
|
KBDS(17) = KOUT
|
|
KBDS(19) = NPTS
|
|
C
|
|
C DO 4100 J = 1, NPTS, 53
|
|
C WRITE (6,4101) (IWORK(K),K=J,J+52)
|
|
4101 FORMAT (1X,25I4)
|
|
C PRINT *,' '
|
|
C4100 CONTINUE
|
|
C
|
|
C INITIALIZE BIT MAP POINTER
|
|
IBMP2P = 0
|
|
C CONSTRUCT WORKING BIT MAP
|
|
DO 2000 I = 1, KOUT
|
|
DO 1000 J = 1, KIN
|
|
IF (J.EQ.1) THEN
|
|
CALL SBYTEC (ISOBMP,1,IBMP2P,1)
|
|
ELSE
|
|
CALL SBYTEC (ISOBMP,0,IBMP2P,1)
|
|
END IF
|
|
IBMP2P = IBMP2P + 1
|
|
1000 CONTINUE
|
|
2000 CONTINUE
|
|
LEN = IBMP2P / 32 + 1
|
|
C CALL BINARY(ISOBMP,LEN)
|
|
C
|
|
C PROCESS OUTER LOOP OF ROW BY ROW OR COL BY COL
|
|
C
|
|
KPTR = 1
|
|
KBDS(11) = KWIDE
|
|
DO 6000 I = 1, KOUT
|
|
C IN CURRENT ROW OR COL
|
|
C FIND FIRST ORDER VALUE
|
|
JPTR = KPTR
|
|
LOWEST = IWORK(JPTR)
|
|
DO 4000 J = 1, KIN
|
|
IF (IWORK(JPTR).LT.LOWEST) THEN
|
|
LOWEST = IWORK(JPTR)
|
|
END IF
|
|
JPTR = JPTR + 1
|
|
4000 CONTINUE
|
|
C SAVE FIRST ORDER VALUE
|
|
CALL SBYTEC (IFOVAL,LOWEST,IFOPTR,KWIDE)
|
|
IFOPTR = IFOPTR + KWIDE
|
|
C PRINT *,'FOVAL',I,LOWEST,KWIDE
|
|
C SUBTRACT FIRST ORDER VALUE FROM OTHER VALS
|
|
C GETTING SECOND ORDER VALUES
|
|
JPTR = KPTR
|
|
IBIG = IWORK(JPTR) - LOWEST
|
|
DO 4200 J = 1, KIN
|
|
IWORK(JPTR) = IWORK(JPTR) - LOWEST
|
|
IF (IWORK(JPTR).GT.IBIG) THEN
|
|
IBIG = IWORK(JPTR)
|
|
END IF
|
|
JPTR = JPTR + 1
|
|
4200 CONTINUE
|
|
C HOW MANY BITS TO CONTAIN LARGEST SECOND
|
|
C ORDER VALUE IN SEGMENT
|
|
CALL FI7505 (IBIG,NWIDE)
|
|
C SAVE BIT WIDTH
|
|
CALL SBYTEC (ISOWID,NWIDE,IWDPTR,8)
|
|
IWDPTR = IWDPTR + 8
|
|
C PRINT *,I,'SOVAL',IBIG,' IN',NWIDE,' BITS'
|
|
C WRITE (6,4101) (IWORK(K),K=KPTR,KPTR+52)
|
|
C SAVE SECOND ORDER VALUES OF THIS SEGMENT
|
|
DO 5000 J = 0, KIN-1
|
|
CALL SBYTEC (ISOVAL,IWORK(KPTR+J),ISOPTR,NWIDE)
|
|
ISOPTR = ISOPTR + NWIDE
|
|
5000 CONTINUE
|
|
KPTR = KPTR + KIN
|
|
6000 CONTINUE
|
|
C =======================================================
|
|
C CONCANTENATE ALL FIELDS FOR BDS
|
|
C
|
|
C REMAINDER GOES INTO IPFLD
|
|
IPTR = 0
|
|
C BYTES 12-13
|
|
C VALUE FOR N1
|
|
C LEAVE SPACE FOR THIS
|
|
IPTR = IPTR + 16
|
|
C BYTE 14
|
|
C EXTENDED FLAGS
|
|
CALL SBYTEC (IPFLD,IBDSFL(5),IPTR,1)
|
|
IPTR = IPTR + 1
|
|
CALL SBYTEC (IPFLD,IBDSFL(6),IPTR,1)
|
|
IPTR = IPTR + 1
|
|
CALL SBYTEC (IPFLD,IBDSFL(7),IPTR,1)
|
|
IPTR = IPTR + 1
|
|
CALL SBYTEC (IPFLD,IBDSFL(8),IPTR,1)
|
|
IPTR = IPTR + 1
|
|
CALL SBYTEC (IPFLD,IBDSFL(9),IPTR,1)
|
|
IPTR = IPTR + 1
|
|
CALL SBYTEC (IPFLD,IBDSFL(10),IPTR,1)
|
|
IPTR = IPTR + 1
|
|
CALL SBYTEC (IPFLD,IBDSFL(11),IPTR,1)
|
|
IPTR = IPTR + 1
|
|
CALL SBYTEC (IPFLD,IBDSFL(12),IPTR,1)
|
|
IPTR = IPTR + 1
|
|
C BYTES 15-16
|
|
C SKIP OVER VALUE FOR N2
|
|
IPTR = IPTR + 16
|
|
C BYTES 17-18
|
|
C P1
|
|
CALL SBYTEC (IPFLD,KBDS(17),IPTR,16)
|
|
IPTR = IPTR + 16
|
|
C BYTES 19-20
|
|
C P2
|
|
CALL SBYTEC (IPFLD,KBDS(19),IPTR,16)
|
|
IPTR = IPTR + 16
|
|
C BYTE 21 - RESERVED LOCATION
|
|
CALL SBYTEC (IPFLD,0,IPTR,8)
|
|
IPTR = IPTR + 8
|
|
C BYTES 22 - ?
|
|
C WIDTHS OF SECOND ORDER PACKING
|
|
IX = (IWDPTR + 32) / 32
|
|
C CALL SBYTESC (IPFLD,ISOWID,IPTR,32,0,IX)
|
|
ijk=IWDPTR/8
|
|
jst=(iptr/8)+1
|
|
ipfld(jst:jst+ijk)=ISOWID(1:ijk)
|
|
IPTR = IPTR + IWDPTR
|
|
C PRINT *,'ISOWID',IWDPTR,IX
|
|
C CALL BINARY (ISOWID,IX)
|
|
C
|
|
C NO SECONDARY BIT MAP
|
|
|
|
C DETERMINE LOCATION FOR START
|
|
C OF FIRST ORDER PACKED VALUES
|
|
KBDS(12) = IPTR / 8 + 12
|
|
C STORE LOCATION
|
|
CALL SBYTEC (IPFLD,KBDS(12),0,16)
|
|
C MOVE IN FIRST ORDER PACKED VALUES
|
|
IPASS = (IFOPTR + 32) / 32
|
|
c CALL SBYTESC (IPFLD,IFOVAL,IPTR,32,0,IPASS)
|
|
ijk=(IFOPTR/8)+1
|
|
jst=(iptr/8)+1
|
|
ipfld(jst:jst+ijk)=ifoval(1:ijk)
|
|
IPTR = IPTR + IFOPTR
|
|
C PRINT *,'IFOVAL',IFOPTR,IPASS,KWIDE
|
|
C CALL BINARY (IFOVAL,IPASS)
|
|
IF (MOD(IPTR,8).NE.0) THEN
|
|
IPTR = IPTR + 8 - MOD(IPTR,8)
|
|
END IF
|
|
C PRINT *,'IFOPTR =',IFOPTR,' ISOPTR =',ISOPTR
|
|
C DETERMINE LOCATION FOR START
|
|
C OF SECOND ORDER VALUES
|
|
KBDS(15) = IPTR / 8 + 12
|
|
C SAVE LOCATION OF SECOND ORDER VALUES
|
|
CALL SBYTEC (IPFLD,KBDS(15),24,16)
|
|
C MOVE IN SECOND ORDER PACKED VALUES
|
|
IX = (ISOPTR + 32) / 32
|
|
C CALL SBYTESC (IPFLD,ISOVAL,IPTR,32,0,IX)
|
|
ijk=(ISOPTR/8)+1
|
|
jst=(iptr/8)+1
|
|
ipfld(jst:jst+ijk)=isoval(1:ijk)
|
|
IPTR = IPTR + ISOPTR
|
|
C PRINT *,'ISOVAL',ISOPTR,IX
|
|
C CALL BINARY (ISOVAL,IX)
|
|
NLEFT = MOD(IPTR+88,16)
|
|
IF (NLEFT.NE.0) THEN
|
|
NLEFT = 16 - NLEFT
|
|
IPTR = IPTR + NLEFT
|
|
END IF
|
|
C COMPUTE LENGTH OF DATA PORTION
|
|
LEN = IPTR / 8
|
|
C COMPUTE LENGTH OF BDS
|
|
LENBDS = LEN + 11
|
|
C -----------------------------------
|
|
C BYTES 1-3
|
|
C THIS FUNCTION COMPLETED BELOW
|
|
C WHEN LENGTH OF BDS IS KNOWN
|
|
CALL SBYTEC (BDS11,LENBDS,0,24)
|
|
C BYTE 4
|
|
CALL SBYTEC (BDS11,IBDSFL(1),24,1)
|
|
CALL SBYTEC (BDS11,IBDSFL(2),25,1)
|
|
CALL SBYTEC (BDS11,IBDSFL(3),26,1)
|
|
CALL SBYTEC (BDS11,IBDSFL(4),27,1)
|
|
C ENTER NUMBER OF FILL BITS
|
|
CALL SBYTEC (BDS11,NLEFT,28,4)
|
|
C BYTE 5-6
|
|
IF (ISCAL2.LT.0) THEN
|
|
CALL SBYTEC (BDS11,1,32,1)
|
|
ISCAL2 = - ISCAL2
|
|
ELSE
|
|
CALL SBYTEC (BDS11,0,32,1)
|
|
END IF
|
|
CALL SBYTEC (BDS11,ISCAL2,33,15)
|
|
C
|
|
C$ FILL OCTET 7-10 WITH THE REFERENCE VALUE
|
|
C CONVERT THE FLOATING POINT OF YOUR MACHINE TO IBM370 32 BIT
|
|
C FLOATING POINT NUMBER
|
|
C REFERENCE VALUE
|
|
C FIRST TEST TO SEE IF
|
|
C ON 32 OR 64 BIT COMPUTER
|
|
C CALL W3FI01(LW)
|
|
IF (bit_size(LW).EQ.32) THEN
|
|
CALL W3FI76 (REFNCE,IEXP,IMANT,32)
|
|
ELSE
|
|
CALL W3FI76 (REFNCE,IEXP,IMANT,64)
|
|
END IF
|
|
CALL SBYTEC (BDS11,IEXP,48,8)
|
|
CALL SBYTEC (BDS11,IMANT,56,24)
|
|
C
|
|
C BYTE 11
|
|
C
|
|
CALL SBYTEC (BDS11,KBDS(11),80,8)
|
|
C
|
|
KLEN = LENBDS / 4 + 1
|
|
C PRINT *,'BDS11 LISTING',4,LENBDS
|
|
C CALL BINARY (BDS11,4)
|
|
C PRINT *,'IPFLD LISTING'
|
|
C CALL BINARY (IPFLD,KLEN)
|
|
RETURN
|
|
END
|
|
SUBROUTINE FI7505 (N,NBITS)
|
|
C$$$ SUBPROGRAM DOCUMENTATION BLOCK
|
|
C . . . .
|
|
C SUBPROGRAM: FI7505 DETERMINE NUMBER OF BITS TO CONTAIN VALUE
|
|
C PRGMMR: CAVANAUGH ORG: W/NMC42 DATE: 93-06-23
|
|
C
|
|
C ABSTRACT: CALCULATE NUMBER OF BITS TO CONTAIN VALUE N, WITH A
|
|
C MAXIMUM OF 32 BITS.
|
|
C
|
|
C PROGRAM HISTORY LOG:
|
|
C 93-06-23 CAVANAUGH
|
|
C 95-10-31 IREDELL REMOVED SAVES AND PRINTS
|
|
C
|
|
C USAGE: CALL FI7505 (N,NBITS)
|
|
C INPUT ARGUMENT LIST:
|
|
C N - INTEGER VALUE
|
|
C
|
|
C OUTPUT ARGUMENT LIST: (INCLUDING WORK ARRAYS)
|
|
C NBITS - NUMBER OF BITS TO CONTAIN N
|
|
C
|
|
C REMARKS: SUBPROGRAM CAN BE CALLED FROM A MULTIPROCESSING ENVIRONMENT.
|
|
C
|
|
C ATTRIBUTES:
|
|
C LANGUAGE: IBM VS FORTRAN 77, CRAY CFT77 FORTRAN
|
|
C MACHINE: HDS, CRAY C916/256, Y-MP8/64, Y-MP EL92/256
|
|
C
|
|
C$$$
|
|
INTEGER N,NBITS
|
|
INTEGER IBITS(31)
|
|
C
|
|
DATA IBITS/1,3,7,15,31,63,127,255,511,1023,2047,
|
|
* 4095,8191,16383,32767,65535,131071,262143,
|
|
* 524287,1048575,2097151,4194303,8388607,
|
|
* 16777215,33554431,67108863,134217727,268435455,
|
|
* 536870911,1073741823,2147483647/
|
|
C ----------------------------------------------------------------
|
|
C
|
|
DO 1000 NBITS = 1, 31
|
|
IF (N.LE.IBITS(NBITS)) THEN
|
|
RETURN
|
|
END IF
|
|
1000 CONTINUE
|
|
RETURN
|
|
END
|
|
SUBROUTINE FI7513 (IWORK,ISTART,NPTS,MAX,MIN,INRNGE)
|
|
C$$$ SUBPROGRAM DOCUMENTATION BLOCK
|
|
C . . . .
|
|
C SUBPROGRAM: FI7513 SELECT BLOCK OF DATA FOR PACKING
|
|
C PRGMMR: CAVANAUGH ORG: W/NMC42 DATE: 94-01-21
|
|
C
|
|
C ABSTRACT: SELECT A BLOCK OF DATA FOR PACKING
|
|
C
|
|
C PROGRAM HISTORY LOG:
|
|
C 94-01-21 CAVANAUGH
|
|
C 95-10-31 IREDELL REMOVED SAVES AND PRINTS
|
|
C
|
|
C USAGE: CALL FI7513 (IWORK,ISTART,NPTS,MAX,MIN,INRNGE)
|
|
C INPUT ARGUMENT LIST:
|
|
C * - RETURN ADDRESS IF ENCOUNTER SET OF SAME VALUES
|
|
C IWORK -
|
|
C ISTART -
|
|
C NPTS -
|
|
C
|
|
C OUTPUT ARGUMENT LIST: (INCLUDING WORK ARRAYS)
|
|
C MAX -
|
|
C MIN -
|
|
C INRNGE -
|
|
C
|
|
C REMARKS: SUBPROGRAM CAN BE CALLED FROM A MULTIPROCESSING ENVIRONMENT.
|
|
C
|
|
C ATTRIBUTES:
|
|
C LANGUAGE: IBM VS FORTRAN 77, CRAY CFT77 FORTRAN
|
|
C MACHINE: HDS, CRAY C916/256, Y-MP8/64, Y-MP EL92/256
|
|
C
|
|
C$$$
|
|
INTEGER IWORK(*),NPTS,ISTART,INRNGE,INRNGA,INRNGB
|
|
INTEGER MAX,MIN,MXVAL,MAXB,MINB,MXVALB
|
|
INTEGER IBITS(31)
|
|
C
|
|
DATA IBITS/1,3,7,15,31,63,127,255,511,1023,2047,
|
|
* 4095,8191,16383,32767,65535,131071,262143,
|
|
* 524287,1048575,2097151,4194303,8388607,
|
|
* 16777215,33554431,67108863,134217727,268435455,
|
|
* 536870911,1073741823,2147483647/
|
|
C ----------------------------------------------------------------
|
|
C IDENTIFY NEXT BLOCK OF DATA FOR PACKING AND
|
|
C RETURN TO CALLER
|
|
C ********************************************************************
|
|
ISTRTA = ISTART
|
|
C
|
|
C GET BLOCK A
|
|
CALL FI7516 (IWORK,NPTS,INRNGA,ISTRTA,
|
|
* MAX,MIN,MXVAL,LWIDE)
|
|
C ********************************************************************
|
|
C
|
|
ISTRTB = ISTRTA + INRNGA
|
|
2000 CONTINUE
|
|
C IF HAVE PROCESSED ALL DATA, RETURN
|
|
IF (ISTRTB.GT.NPTS) THEN
|
|
C NO MORE DATA TO LOOK AT
|
|
INRNGE = INRNGA
|
|
RETURN
|
|
END IF
|
|
C GET BLOCK B
|
|
CALL FI7502 (IWORK,ISTRTB,NPTS,ISAME)
|
|
IF (ISAME.GE.15) THEN
|
|
C PRINT *,'BLOCK B HAS ALL IDENTICAL VALUES'
|
|
C PRINT *,'BLOCK A HAS INRNGE =',INRNGA
|
|
C BLOCK B CONTAINS ALL IDENTICAL VALUES
|
|
INRNGE = INRNGA
|
|
C EXIT WITH BLOCK A
|
|
RETURN
|
|
END IF
|
|
C GET BLOCK B
|
|
C
|
|
ISTRTB = ISTRTA + INRNGA
|
|
CALL FI7516 (IWORK,NPTS,INRNGB,ISTRTB,
|
|
* MAXB,MINB,MXVALB,LWIDEB)
|
|
C PRINT *,'BLOCK A',INRNGA,' BLOCK B',INRNGB
|
|
C ********************************************************************
|
|
C PERFORM TREND ANALYSIS TO DETERMINE
|
|
C IF DATA COLLECTION CAN BE IMPROVED
|
|
C
|
|
KTRND = LWIDE - LWIDEB
|
|
C PRINT *,'TREND',LWIDE,LWIDEB
|
|
IF (KTRND.LE.0) THEN
|
|
C PRINT *,'BLOCK A - SMALLER, SHOULD EXTEND INTO BLOCK B'
|
|
MXVAL = IBITS(LWIDE)
|
|
C
|
|
C IF BLOCK A REQUIRES THE SAME OR FEWER BITS
|
|
C LOOK AHEAD
|
|
C AND GATHER THOSE DATA POINTS THAT CAN
|
|
C BE RETAINED IN BLOCK A
|
|
C BECAUSE THIS BLOCK OF DATA
|
|
C USES FEWER BITS
|
|
C
|
|
CALL FI7518 (IRET,IWORK,NPTS,ISTRTA,INRNGA,INRNGB,
|
|
* MAX,MIN,LWIDE,MXVAL)
|
|
IF(IRET.EQ.1) GO TO 8000
|
|
C PRINT *,'18 INRNGA IS NOW ',INRNGA
|
|
IF (INRNGB.LT.20) THEN
|
|
RETURN
|
|
ELSE
|
|
GO TO 2000
|
|
END IF
|
|
ELSE
|
|
C PRINT *,'BLOCK A - LARGER, B SHOULD EXTEND BACK INTO A'
|
|
MXVALB = IBITS(LWIDEB)
|
|
C
|
|
C IF BLOCK B REQUIRES FEWER BITS
|
|
C LOOK BACK
|
|
C SHORTEN BLOCK A BECAUSE NEXT BLOCK OF DATA
|
|
C USES FEWER BITS
|
|
C
|
|
CALL FI7517 (IRET,IWORK,NPTS,ISTRTB,INRNGA,
|
|
* MAXB,MINB,LWIDEB,MXVALB)
|
|
IF(IRET.EQ.1) GO TO 8000
|
|
C PRINT *,'17 INRNGA IS NOW ',INRNGA
|
|
END IF
|
|
C
|
|
C PACK UP BLOCK A
|
|
C UPDATA POINTERS
|
|
8000 CONTINUE
|
|
INRNGE = INRNGA
|
|
C GET NEXT BLOCK A
|
|
9000 CONTINUE
|
|
RETURN
|
|
END
|
|
SUBROUTINE FI7516 (IWORK,NPTS,INRNG,ISTART,MAX,MIN,MXVAL,LWIDTH)
|
|
C$$$ SUBPROGRAM DOCUMENTATION BLOCK
|
|
C . . . .
|
|
C SUBPROGRAM: FI7516 SCAN NUMBER OF POINTS
|
|
C PRGMMR: CAVANAUGH ORG: W/NMC42 DATE: 94-01-21
|
|
C
|
|
C ABSTRACT: SCAN FORWARD FROM CURRENT POSITION. COLLECT POINTS AND
|
|
C DETERMINE MAXIMUM AND MINIMUM VALUES AND THE NUMBER
|
|
C OF POINTS THAT ARE INCLUDED. FORWARD SEARCH IS TERMINATED
|
|
C BY ENCOUNTERING A SET OF IDENTICAL VALUES, BY REACHING
|
|
C THE NUMBER OF POINTS SELECTED OR BY REACHING THE END
|
|
C OF DATA.
|
|
C
|
|
C PROGRAM HISTORY LOG:
|
|
C 94-01-21 CAVANAUGH
|
|
C 95-10-31 IREDELL REMOVED SAVES AND PRINTS
|
|
C
|
|
C USAGE: CALL FI7516 (IWORK,NPTS,INRNG,ISTART,MAX,MIN,MXVAL,LWIDTH)
|
|
C INPUT ARGUMENT LIST:
|
|
C * - RETURN ADDRESS IF ENCOUNTER SET OF SAME VALUES
|
|
C IWORK - DATA ARRAY
|
|
C NPTS - NUMBER OF POINTS IN DATA ARRAY
|
|
C ISTART - STARTING LOCATION IN DATA
|
|
C
|
|
C OUTPUT ARGUMENT LIST: (INCLUDING WORK ARRAYS)
|
|
C INRNG - NUMBER OF POINTS SELECTED
|
|
C MAX - MAXIMUM VALUE OF POINTS
|
|
C MIN - MINIMUM VALUE OF POINTS
|
|
C MXVAL - MAXIMUM VALUE THAT CAN BE CONTAINED IN LWIDTH BITS
|
|
C LWIDTH - NUMBER OF BITS TO CONTAIN MAX DIFF
|
|
C
|
|
C REMARKS: SUBPROGRAM CAN BE CALLED FROM A MULTIPROCESSING ENVIRONMENT.
|
|
C
|
|
C ATTRIBUTES:
|
|
C LANGUAGE: IBM VS FORTRAN 77, CRAY CFT77 FORTRAN
|
|
C MACHINE: HDS, CRAY C916/256, Y-MP8/64, Y-MP EL92/256
|
|
C
|
|
C$$$
|
|
INTEGER IWORK(*),NPTS,ISTART,INRNG,MAX,MIN,LWIDTH,MXVAL
|
|
INTEGER IBITS(31)
|
|
C
|
|
DATA IBITS/1,3,7,15,31,63,127,255,511,1023,2047,
|
|
* 4095,8191,16383,32767,65535,131071,262143,
|
|
* 524287,1048575,2097151,4194303,8388607,
|
|
* 16777215,33554431,67108863,134217727,268435455,
|
|
* 536870911,1073741823,2147483647/
|
|
C ----------------------------------------------------------------
|
|
C
|
|
INRNG = 1
|
|
JQ = ISTART + 19
|
|
MAX = IWORK(ISTART)
|
|
MIN = IWORK(ISTART)
|
|
DO 1000 I = ISTART+1, JQ
|
|
CALL FI7502 (IWORK,I,NPTS,ISAME)
|
|
IF (ISAME.GE.15) THEN
|
|
GO TO 5000
|
|
END IF
|
|
INRNG = INRNG + 1
|
|
IF (IWORK(I).GT.MAX) THEN
|
|
MAX = IWORK(I)
|
|
ELSE IF (IWORK(I).LT.MIN) THEN
|
|
MIN = IWORK(I)
|
|
END IF
|
|
1000 CONTINUE
|
|
5000 CONTINUE
|
|
KRNG = MAX - MIN
|
|
C
|
|
DO 9000 LWIDTH = 1, 31
|
|
IF (KRNG.LE.IBITS(LWIDTH)) THEN
|
|
C PRINT *,'RETURNED',INRNG,' VALUES'
|
|
RETURN
|
|
END IF
|
|
9000 CONTINUE
|
|
RETURN
|
|
END
|
|
SUBROUTINE FI7517 (IRET,IWORK,NPTS,ISTRTB,INRNGA,
|
|
* MAXB,MINB,MXVALB,LWIDEB)
|
|
C$$$ SUBPROGRAM DOCUMENTATION BLOCK
|
|
C . . . .
|
|
C SUBPROGRAM: FI7517 SCAN BACKWARD
|
|
C PRGMMR: CAVANAUGH ORG: W/NMC42 DATE: 94-01-21
|
|
C
|
|
C ABSTRACT: SCAN BACKWARDS UNTIL A VALUE EXCEEDS RANGE OF GROUP B
|
|
C THIS MAY SHORTEN GROUP A
|
|
C
|
|
C PROGRAM HISTORY LOG:
|
|
C 94-01-21 CAVANAUGH
|
|
C 95-10-31 IREDELL REMOVED SAVES AND PRINTS
|
|
C 98-06-17 IREDELL REMOVED ALTERNATE RETURN
|
|
C
|
|
C USAGE: CALL FI7517 (IRET,IWORK,NPTS,ISTRTB,INRNGA,
|
|
C * MAXB,MINB,MXVALB,LWIDEB)
|
|
C INPUT ARGUMENT LIST:
|
|
C IWORK -
|
|
C ISTRTB -
|
|
C NPTS -
|
|
C INRNGA -
|
|
C
|
|
C OUTPUT ARGUMENT LIST: (INCLUDING WORK ARRAYS)
|
|
C IRET -
|
|
C JLAST -
|
|
C MAXB -
|
|
C MINB -
|
|
C LWIDTH - NUMBER OF BITS TO CONTAIN MAX DIFF
|
|
C
|
|
C REMARKS: SUBPROGRAM CAN BE CALLED FROM A MULTIPROCESSING ENVIRONMENT.
|
|
C
|
|
C ATTRIBUTES:
|
|
C LANGUAGE: IBM VS FORTRAN 77, CRAY CFT77 FORTRAN
|
|
C MACHINE: HDS, CRAY C916/256, Y-MP8/64, Y-MP EL92/256
|
|
C
|
|
C$$$
|
|
INTEGER IWORK(*),NPTS,ISTRTB,INRNGA
|
|
INTEGER MAXB,MINB,LWIDEB,MXVALB
|
|
INTEGER IBITS(31)
|
|
C
|
|
DATA IBITS/1,3,7,15,31,63,127,255,511,1023,2047,
|
|
* 4095,8191,16383,32767,65535,131071,262143,
|
|
* 524287,1048575,2097151,4194303,8388607,
|
|
* 16777215,33554431,67108863,134217727,268435455,
|
|
* 536870911,1073741823,2147483647/
|
|
C ----------------------------------------------------------------
|
|
IRET=0
|
|
C PRINT *,' FI7517'
|
|
NPOS = ISTRTB - 1
|
|
ITST = 0
|
|
KSET = INRNGA
|
|
C
|
|
1000 CONTINUE
|
|
C PRINT *,'TRY NPOS',NPOS,IWORK(NPOS),MAXB,MINB
|
|
ITST = ITST + 1
|
|
IF (ITST.LE.KSET) THEN
|
|
IF (IWORK(NPOS).GT.MAXB) THEN
|
|
IF ((IWORK(NPOS)-MINB).GT.MXVALB) THEN
|
|
C PRINT *,'WENT OUT OF RANGE AT',NPOS
|
|
IRET=1
|
|
RETURN
|
|
ELSE
|
|
MAXB = IWORK(NPOS)
|
|
END IF
|
|
ELSE IF (IWORK(NPOS).LT.MINB) THEN
|
|
IF ((MAXB-IWORK(NPOS)).GT.MXVALB) THEN
|
|
C PRINT *,'WENT OUT OF RANGE AT',NPOS
|
|
IRET=1
|
|
RETURN
|
|
ELSE
|
|
MINB = IWORK(NPOS)
|
|
END IF
|
|
END IF
|
|
INRNGA = INRNGA - 1
|
|
NPOS = NPOS - 1
|
|
GO TO 1000
|
|
END IF
|
|
C ----------------------------------------------------------------
|
|
C
|
|
9000 CONTINUE
|
|
RETURN
|
|
END
|
|
SUBROUTINE FI7518 (IRET,IWORK,NPTS,ISTRTA,INRNGA,INRNGB,
|
|
* MAXA,MINA,LWIDEA,MXVALA)
|
|
C$$$ SUBPROGRAM DOCUMENTATION BLOCK
|
|
C . . . .
|
|
C SUBPROGRAM: FI7518 SCAN FORWARD
|
|
C PRGMMR: CAVANAUGH ORG: W/NMC42 DATE: 94-01-21
|
|
C
|
|
C ABSTRACT: SCAN FORWARD FROM START OF BLOCK B TOWARDS END OF BLOCK B
|
|
C IF NEXT POINT UNDER TEST FORCES A LARGER MAXVALA THEN
|
|
C TERMINATE INDICATING LAST POINT TESTED FOR INCLUSION
|
|
C INTO BLOCK A.
|
|
C
|
|
C PROGRAM HISTORY LOG:
|
|
C 94-01-21 CAVANAUGH
|
|
C 95-10-31 IREDELL REMOVED SAVES AND PRINTS
|
|
C 98-06-17 IREDELL REMOVED ALTERNATE RETURN
|
|
C
|
|
C USAGE: CALL FI7518 (IRET,IWORK,NPTS,ISTRTA,INRNGA,INRNGB,
|
|
C * MAXA,MINA,LWIDEA,MXVALA)
|
|
C INPUT ARGUMENT LIST:
|
|
C IFLD -
|
|
C JSTART -
|
|
C NPTS -
|
|
C
|
|
C OUTPUT ARGUMENT LIST: (INCLUDING WORK ARRAYS)
|
|
C IRET -
|
|
C JLAST -
|
|
C MAX -
|
|
C MIN -
|
|
C LWIDTH - NUMBER OF BITS TO CONTAIN MAX DIFF
|
|
C
|
|
C REMARKS: SUBPROGRAM CAN BE CALLED FROM A MULTIPROCESSING ENVIRONMENT.
|
|
C
|
|
C ATTRIBUTES:
|
|
C LANGUAGE: IBM VS FORTRAN 77, CRAY CFT77 FORTRAN
|
|
C MACHINE: HDS, CRAY C916/256, Y-MP8/64, Y-MP EL92/256
|
|
C
|
|
C$$$
|
|
INTEGER IWORK(*),NPTS,ISTRTA,INRNGA
|
|
INTEGER MAXA,MINA,LWIDEA,MXVALA
|
|
INTEGER IBITS(31)
|
|
C
|
|
DATA IBITS/1,3,7,15,31,63,127,255,511,1023,2047,
|
|
* 4095,8191,16383,32767,65535,131071,262143,
|
|
* 524287,1048575,2097151,4194303,8388607,
|
|
* 16777215,33554431,67108863,134217727,268435455,
|
|
* 536870911,1073741823,2147483647/
|
|
C ----------------------------------------------------------------
|
|
IRET=0
|
|
C PRINT *,' FI7518'
|
|
NPOS = ISTRTA + INRNGA
|
|
ITST = 0
|
|
C
|
|
1000 CONTINUE
|
|
ITST = ITST + 1
|
|
IF (ITST.LE.INRNGB) THEN
|
|
C PRINT *,'TRY NPOS',NPOS,IWORK(NPOS),MAXA,MINA
|
|
IF (IWORK(NPOS).GT.MAXA) THEN
|
|
IF ((IWORK(NPOS)-MINA).GT.MXVALA) THEN
|
|
C PRINT *,'FI7518A -',ITST,' RANGE EXCEEDS MAX'
|
|
IRET=1
|
|
RETURN
|
|
ELSE
|
|
MAXA = IWORK(NPOS)
|
|
END IF
|
|
ELSE IF (IWORK(NPOS).LT.MINA) THEN
|
|
IF ((MAXA-IWORK(NPOS)).GT.MXVALA) THEN
|
|
C PRINT *,'FI7518B -',ITST,' RANGE EXCEEDS MAX'
|
|
IRET=1
|
|
RETURN
|
|
ELSE
|
|
MINA = IWORK(NPOS)
|
|
END IF
|
|
END IF
|
|
INRNGA = INRNGA + 1
|
|
C PRINT *,' ',ITST,INRNGA
|
|
NPOS = NPOS +1
|
|
GO TO 1000
|
|
END IF
|
|
C ----------------------------------------------------------------
|
|
9000 CONTINUE
|
|
RETURN
|
|
END
|