program ql21sun
c
c      read CMR5 data on a SUN computer system
c      because UNIX doesn't use record stucture within files
c      this program simulates it using direct access operations
c
      IMPLICIT INTEGER(A-E,G-Z)
C
      PARAMETER(MAX=10,FMISS=-999.9,INDEXCK=35,BPR=192,
     * ID=0,BPW=32,dimbuf=(bpr/bpw+6),dimpk=(BPR-1)/BPW+1,
     * dimun=35,unit=10)
C
      CHARACTER FNAME*64

      COMMON /CMR5/FIELD(35),FTRUEL(35),FTRUEU(35),FUNITS(35)
     +,FBASE(35),BITS(35),OFFSET(35)
C
      DIMENSION BUF(dimbuf),PK(dimpk),UN(dimun),FTRUE(dimun)
C
      DATA FIELD/8HBOX10   ,8HMONTH   ,8HBOX2    ,8HYEAR    ,8HDAY     ,
     +8HHOUR    ,8HX       ,8HY       ,8HS       ,8HBI      ,8HA       ,
     +8HDP      ,8HTI      ,8HU       ,8HV       ,8HDI      ,8HWI      ,
     +8HP       ,8HC       ,8HNH      ,8HCL      ,8HH       ,8HHI      ,
     +8HCM      ,8HCH      ,8HST      ,8HPW      ,8HCD      ,8HLF      ,
     +8HSF      ,8HAF      ,8HRF      ,8HWF      ,8HPF      ,8HCK      /
C
      DATA FTRUEL/3*1.,1800.,1.,3*0.,-5.,0.,-88.,2*0.,2*-102.2,2*0.,870.
     +,17*0./
C
      DATA FTRUEU/648.,12.,16202.,2054.,31.,23.,2*2.,40.,2.,58.,70.,5.
     +,2*102.2,5.,1.,1074.6,2*9.,2*10.,1.,2*10.,7.,99.,999.,0.,5*2.,30./
C
      DATA FUNITS/6*1.,3*.1,1.,2*.1,1.,2*.1,2*1.,.1,17*1./
C
      DATA FBASE/3*0,1799,0,3*-1,-51,-1,-881,2*-1,2*-1023,2*-1,8699
     +,16*-1,0/
C
      DATA BITS/10,4,14,8,4*5,9,2,11,10,3,2*11,3,2,11,4*4,2,3*4,7,10
     +,1,5*2,5/
C
      DATA OFFSET/
     +   0, 10, 14, 28, 36, 41, 46, 51, 56, 65, 67, 78, 88, 91,102,113
     +,116,118,129,133,137,141,145,147,151,155,159,166,176,177,179,181
     +,183,185,187/

      DATA BUF/dimbuf*0/
C
c     open UNIX data file - direct access mode
c     each read will then access one CMR record
      WRITE(*,6) 
  6   FORMAT('  ENTER INPUT FILENAME : ',$)     
      READ(*,'(A)') FNAME

      open(unit=unit,file=fname,
     *       form='unformatted',access='direct',recl=24)
c
      write(*,2)
2     format('1  ****** beginning output for ql21sun *****')
c
C
100   CALL GETRPT(unit,FMISS,FUNITS,FBASE,BITS,OFFSET,INDEXCK,ID
     +,BPR,BPW,BUF,dimbuf,PK,dimpk,UN,dimun,FTRUE,JEOF)
      IF(JEOF.NE.0)GOTO 900
C
c     BUF(2) is an incremental counter of the number of records processed
c     when BUF(2) is equal MAX (10) stop execution.  For complete
c     file processing the code segment below must be modified.

      PRINT 300,(FIELD(I),FTRUE(I),I=1,dimun)
300   FORMAT(6(1X,A5,F7.1))
      IF(BUF(2).LT.MAX)GOTO 100
C
900   PRINT *,' REPORTS ',BUF(2),', EOF ',JEOF
      END
      SUBROUTINE GETRPT(unit,FMISS,FUNITS,FBASE,BITS,OFFSET,INDEXCK,ID
     +,BPR,BPW,BUF,dimbuf,PK,dimpk,UN,dimun,FTRUE,JEOF)
C
C-----RETURN FLOATING POINT VALUES IN FTRUE
C
C     INPUT
C         unit - fortran input unit number
C         FMISS - MISSING VALUE
C         FUNITS(dimun) - UNITS FOR UNCODING
C         FBASE(dimun) - BASE FOR UNCODING
C         BITS(dimun) - BITS FOR UNPACKING
C         OFFSET(dimun) - OFFSET FOR UNPACKING
C         INDEXCK - UN(INDEXCK) = CHECKSUM
C         ID - GROUP NUMBER FOR IDENTIFICATION CHECKSUM
C         BPR - BITS PER REPORT
C         BPW - BITS PER WORD
C     OUTPUT
C         BUF(dimbuf) - RCDIN BUFFER
C         PK(dimpk) - PACKED REPORT
C         UN(dimun) - UNPACKED REPORT
C         FTRUE(dimun) - TRUE VALUES
C         JEOF = end of file indicator, 0=no, -1=yes
C
      IMPLICIT INTEGER(A-E,G-Z)
      DIMENSION FUNITS(dimun),FBASE(dimun),BITS(dimun),OFFSET(dimun)
     +,BUF(dimbuf),PK(dimpk),UN(dimun),FTRUE(dimun)
C

      CALL RCDIN(unit,BUF,dimbuf,PK,dimpk,BPR,BPW,JEOF)
      if(jeof.eq.-1)return
C
C     unpack with gbyte and rescale values

      CK=ID
      DO 230 I=1,dimun
      CALL GBYTE(PK(OFFSET(I)/BPW+1),UN(I),MOD(OFFSET(I),BPW),BITS(I))
      IF(I.EQ.INDEXCK)GOTO 210
      IF(UN(I).EQ.0)GOTO 220
      FTRUE(I)=(UN(I)+FBASE(I))*FUNITS(I)
      CK=CK+UN(I)
      GOTO 230
210   FTRUE(INDEXCK)=UN(INDEXCK)
      GOTO 230
220   FTRUE(I)=FMISS
230   CONTINUE
      IF(MOD(CK,2**BITS(INDEXCK)-1).EQ.UN(INDEXCK))RETURN
C
C-----ERROR
      PRINT *,' SUBROUTINE GETRPT  --  CHECKSUM ERROR, unit = ',unit
     +,', REPORT = ',BUF(2)
      PRINT *,' FTRUE = ',FTRUE
      STOP
C
      END
C=======================================================================
      SUBROUTINE RCDIN(unit,BUF,dimbuf,RCD,dimrcd,BPR,BPW,JEOF)
C
C-----RETURN ONE LOGICAL RECORD IN RCD
C
C     INPUT
C         unit - BUFFER IN UNIT
C         BPR - BITS PER RECORD
C         BPW - BITS PER WORD
C     OUTPUT
C         BUF(dimbuf) - counters and logical record
C         RCD(dimrcd) - LOGICAL RECORD
C         JEOF = end of file indicator, 0=no, 1=yes
C
C     BUF(1) = gbytes offset, always 0 in this case
C     BUF(2) = record count
C     BUF(3) =
C     BUF(4) =
C     BUF(5) = record LENGTH IN BITS
C     BUF(6) =
C
      IMPLICIT INTEGER(A-E,G-Z)
      DIMENSION BUF(dimbuf),RCD(dimrcd)
      data irdata/1/


c     read from UNIX data file

10    read(unit,rec=irdata,end=900)(buf(i),i=7,12)
      irdata=irdata + 1
      CALL SWAP4(BUF(7),BUF(7),24)
c
c     write some data to check, in hexidecimal
c
c     write(6,202)(buf(i),i=7,12)
c202   format(5z12)
c
      jeof=0
100   buf(1)=0
      buf(5)=192
C
C-----GBYTE
      CALL GBYTES(BUF(6+1),RCD,0,BPW,0,dimrcd)
      BUF(2)=BUF(2)+1
c      call cswap_4(rcd,dimrcd)
C      call swap4(rcd,rcd,dimrcd)

      RETURN


c     end of file has been detected set jeof to -1

900   continue
      jeof = -1

      return
      END
      SUBROUTINE GBYTE (IN,IOUT,ISKIP,NBYTE)
      CALL GBYTES (IN,IOUT,ISKIP,NBYTE,0,1)
      RETURN
      END

      SUBROUTINE SBYTE (IOUT,IN,ISKIP,NBYTE)
      CALL SBYTES (IOUT,IN,ISKIP,NBYTE,0,1)
      RETURN
      END
      SUBROUTINE GBYTES(IN,IOUT,ISKIP,NBYTE,NSKIP,N)
C THIS PROGRAM WRITTEN BY.....
C             DR. ROBERT C. GAMMILL, CONSULTANT
C             NATIONAL CENTER FOR ATMOSPHERIC RESEARCH
C             MAY 1972
C THIS IS THE FORTRAN VERSION OF GBYTES.
      COMMON/MACHIN/NBITSW,NBITSC,MASK0,MASKS(32)
      DIMENSION IN(1),IOUT(1)
C  THE STATEMENTS BETWEEN ASTERISK LINES GIVE ALL NECESSARY MACHINE
C  DEPENDENT INFORMATION
C***********************************************************************
C THIS SPECIFICATION IS FOR IBM PC USING MICROSOFT FORTRAN
      INTEGER RGHTSH,OR,AND
      LEFTSH(M,N)=ISHFT(M,N)
      RGHTSH(M,N)=ISHFT(M,-N)
      OR(M,N)=IOR(M,N)
      AND(M,N)=IAND(M,N)
      DATA NBITSW/32/
      DATA MASK0/16#0/
      DATA MASKS/16#1,16#3,16#7,16#F,16#1F,16#3F,16#7F,16#FF,
     2 16#1FF,16#3FF,16#7FF,16#FFF,16#1FFF,16#3FFF,16#7FFF,16#FFFF,
     3 16#1FFFF,16#3FFFF,16#7FFFF,16#FFFFF,16#1FFFFF,16#3FFFFF,
     4 16#7FFFFF,16#FFFFFF,16#1FFFFFF,16#3FFFFFF,16#7FFFFFF,16#FFFFFFF,
     5 16#1FFFFFFF,16#3FFFFFFF,16#7FFFFFFF,16#FFFFFFFF/
C***********************************************************************
C NBYTE MUST BE LESS THAN OR EQUAL TO NBITSW
      ICON=NBITSW-NBYTE
      IF(ICON.LT.0) RETURN
      MASK=MASKS(NBYTE)
C INDEX TELLS HOW MANY WORDS INTO THE ARRAY 'IN' THE NEXT BYTE APPEARS.
      INDEX=ISKIP/NBITSW
C II TELLS HOW MANY BITS THE BYTE IS FROM THE LEFT SIDE OF THE WORD.
      II=MOD(ISKIP,NBITSW)
C ISTEP IS THE DISTANCE IN BITS FROM THE START OF ONE BYTE TO THE NEXT.
      ISTEP=NBYTE+NSKIP
C IWORDS TELLS HOW MANY WORDS TO SKIP FROM ONE BYTE TO THE NEXT.
      IWORDS=ISTEP/NBITSW
C IBITS TELLS HOW MANY BITS TO SKIP AFTER SKIPPING IWORDS.
      IBITS=MOD(ISTEP,NBITSW)
      DO 6 I=1,N
C MOVER SPECIFIES HOW FAR TO THE RIGHT A BYTE MUST BE MOVED IN ORDER
C    TO BE RIGHT ADJUSTED.
      MOVER=ICON-II
      IF(MOVER) 2,3,4
C
C THE BYTE IS SPLIT ACROSS A WORD BREAK.
    2 MOVEL=-MOVER
      MOVER=NBITSW-MOVEL
C EXPAND STATEMENT SINCE MSFORTRAN SEEMS TO FAIL OTHERWISE
      NP1=LEFTSH(IN(INDEX+1),MOVEL)
      NP2=RGHTSH(IN(INDEX+2),MOVER)
C      NP3=OR(LEFTSH(IN(INDEX+1),MOVEL),RGHTSH(IN(INDEX+2),MOVER))
      IOUT(I)=AND(OR(NP1,NP2),MASK)
C      IOUT(I)=AND(OR(LEFTSH(IN(INDEX+1),MOVEL),RGHTSH(IN(INDEX+2),MOVER)
C     1 ),MASK)
      GO TO 5
C
C THE BYTE IS ALREADY RIGHT ADJUSTED.
    3 IOUT(I)=AND(IN(INDEX+1),MASK)
      GO TO 5
C
C RIGHT ADJUST THE BYTE.
    4 IOUT(I)=AND(RGHTSH(IN(INDEX+1),MOVER),MASK)
C
C INCREMENT II AND INDEX.
    5 II=II+IBITS
      INDEX=INDEX+IWORDS
      IF(II.LT.NBITSW) GO TO 6
      II=II-NBITSW
      INDEX=INDEX+1
    6 CONTINUE
      RETURN
      END
      SUBROUTINE SBYTES(IOUT,IN,ISKIP,NBYTE,NSKIP,N)
C THIS PROGRAM WRITTEN BY.....
C             DR. ROBERT C. GAMMILL, CONSULTANT
C             NATIONAL CENTER FOR ATMOSPHERIC RESEARCH
C             JULY 1972
C THIS IS THE FORTRAN VERSIONS OF SBYTES.
      COMMON/MACHIN/NBITSW,NBITSC,MASK0,MASKS(32)
      DIMENSION IN(1),IOUT(1)
C
C  THE STATEMENTS BETWEEN ASTERISK LINES GIVE ALL NECESSARY MACHINE
C  DEPENDENT INFORMATION
C***********************************************************************
C THIS SPECIFICATION IS SUITABLE FOR CDC 6000 SERIES COMPUTERS.
      INTEGER RGHTSH,OR,AND
      OR(M,N)=IOR(M,N)
      AND(M,N)=IAND(M,N)
C
      LEFTSH(M,N)=ISHFT(M,N)
      RGHTSH(M,N)=ISHFT(M,-N)
C***********************************************************************
C
C NBYTE MUST BE LESS THAN OR EQUAL TO NBITSW
      ICON=NBITSW-NBYTE
      IF(ICON.LT.0) RETURN
      MASK=MASKS(NBYTE)
C INDEX TELLS HOW MANY WORDS INTO IOUT THE NEXT BYTE IS TO BE STORED.
      INDEX=ISKIP/NBITSW
C II TELLS HOW MANY BITS IN FROM THE LEFT SIDE OF THE WORD TO STORE IT.
      II=MOD(ISKIP,NBITSW)
C ISTEP IS THE DISTANCE IN BITS FROM ONE BYTE POSITION TO THE NEXT.
      ISTEP=NBYTE+NSKIP
C IWORDS TELLS HOW MANY WORDS TO SKIP FROM ONE BYTE TO THE NEXT.
      IWORDS=ISTEP/NBITSW
C IBITS TELLS HOW MANY BITS TO SKIP AFTER SKIPPING IWORDS.
      IBITS=MOD(ISTEP,NBITSW)
      DO 6 I=1,N
      J=AND(MASK,IN(I))
      MOVEL=ICON-II
      IF(MOVEL) 2,3,4
C
C THE BYTE IS TO BE SPLIT ACROSS A WORD BREAK.
2     MSK=MASKS(NBYTE+MOVEL)
      IOUT(INDEX+1)=OR(AND(NOT(MSK),IOUT(INDEX+1)),RGHTSH(J,-MOVEL))
      ITEMP=AND(MASKS(NBITSW+MOVEL),IOUT(INDEX+2))
      IOUT(INDEX+2)=OR(ITEMP,LEFTSH(J,NBITSW+MOVEL))
      GO TO 5
C
C BYTE IS TO BE STORED RIGHT-ADJUSTED.
3     IOUT(INDEX+1)=OR(AND(NOT(MASK),IOUT(INDEX+1)),J)
      GO TO 5
C
C BYTE IS TO BE STORED IN MIDDLE OF WORD.  SHIFT LEFT.
4     MSK=LEFTSH(MASK,MOVEL)
      IOUT(INDEX+1)=OR(AND(NOT(MSK),IOUT(INDEX+1)),LEFTSH(J,MOVEL))
5     II=II+IBITS
      INDEX=INDEX+IWORDS
      IF(II.LT.NBITSW) GO TO 6
      II=II-NBITSW
      INDEX=INDEX+1
6     CONTINUE
      RETURN
      END
      SUBROUTINE SWAP4(IN,IO,NN)
      LOGICAL*1 IN(1),IO(1),IH
      DO 10 I=1,NN,4
      IH=IN(I)
      IO(I)=IN(I+3)
      IO(I+3)=IH
      IH=IN(I+1)
      IO(I+1)=IN(I+2)
      IO(I+2)=IH
  10  CONTINUE
      RETURN
      END