      PROGRAM CAD4TAPE
C
C         ROBERT H. BLESSING
C         MEDICAL FOUNDATION OF BUFFALO
C         73 HIGH STREET
C         BUFFALO, NEW YORK 14203, USA
C         TELEPHONE (716)856-9600
C
C         AUGUST 1985
C
C         RUNS INTERACTIVELY.
C
C         BASED ON PROGRAM CAD4DT FROM THE ENRAF-NONIUS STRUCTURE
C         DETERMINATION PROGRAM PACKAGE.
C
C         DECODES DATA TAPE CONTAINING THETA/TWO-THETA OR OMEGA SCAN
C         PROFILE DATA WRITTEN BY THE ENRAF-NONIUS CAD4 DIFFRACTOMETER.
C
C         WRITES A BINARY FILE NAME.RAW CONTAINING THE RAW REFLECTION
C         DATA RECORDS:
C
C             JREF      SERIAL NUMBER
C
C             JH
C             JK        MILLER INDICES
C             JL
C
C             THETA
C             PHIK      DIFFRACTOMETER
C             OMEGAK    SETTING ANGLES
C             KAPPA
C
C             WIDTH     SCAN WIDTH (DEGREES OMEGA)
C             SPEED     SCAN SPEED (DEGREES OMEGA PER MINUTE)
C
C             BGL       LEFT BACKGROUND  (SUM STEPS  1-16)
C             PK        PEAK             (SUM STEPS 17-80)
C             BGR       RIGHT BACKGROUND (SUM STEPS 81-96)
C
C             RNET      NET REFLECTION INTENSITY
C             SIGR      ESD ON COUNTING STATISTICS ALONE
C
C             XTIME     X-RAY EXPOSURE TIME (HOURS)
C
C             JD(96)    ARRAY OF RAW STEP SCAN COUNTS
C
C         JREF IS AN INTEGER*4 WORD; JH, JK, JL, AND JD(96) ARE 
C         INTEGER*2 WORDS; AND ALL THE OTHER WORDS ARE REAL.
C
C         JREF IS NEGATIVE FOR STANDARD REFERENCE INTENSITY REFLECTIONS.
C
C         WIDTH IS NEGATIVE IF THE BEAM ATTENUATOR WAS USED.
C
C         RNET = SPEED*(PK - 2*(BGL + BGR))
C         SIGR = SPEED*(SQRT(PK + 4*(BGL + BGR)))
C
C         NOTE THAT RNET AND SIGR AS READ FROM THE DIFFRACTOMETER TAPE
C         ARE RNET/SPEED AND SIGR/SPEED, WHICH ARE NOT NORMALIZED TO
C         ACCOUNT FOR VARIABLE SCAN SPEEDS OR FOR THE BEAM ATTENUATOR
C         FACTOR.
C
C-----------------------------------------------------------------------
C
C         PROGRAM CAD4TAPE CALLS A NUMBER OF MAGNETIC TAPE INPUT-OUTPUT
C         SUBROUTINES THAT ARE MACHINE SPECIFIC.  THESE CALLS ARE:
C
C             CALL ERRSET (64,.TRUE.,.FALSE.,.FALSE.,.TRUE.,31)
C             CALL READMT (IBUF,NBYTE)
C             CALL WAITRW (ISTAT,ICNT)
C             CALL SKPFMT (1)
C             CALL CLOSE_MT
C             CALL MTERR (ISTAT)
C
C         LISTINGS AND DESCRIPTIONS OF VAX ASSEMBLY LANGUAGE VERSIONS OF
C         THESE ROUTINES FROM THE MFB LIBRARY ARE INCLUDED IN THE
C         ACCOMPANYING FILE MAGIO.MAR.
C                           ---------
C-----------------------------------------------------------------------
C
C         DATA ON CAD4 MAGNETIC TAPE ARE ORGANIZED IN BLOCKS OF SIX
C         RECORDS.  EACH RECORD CONSISTS OF 64 BYTES.
C
C         64 BYTES PER RECORD
C         6 RECORDS PER BLOCK
C         384 BYTES PER BLOCK
C
C         8 BITS PER BYTE
C         1 BYTE PER ASCII CHARACTER
C         2 CHARACTERS PER INTEGER*2 WORD
C
      PARAMETER (NREC=6,NBYT=64,NBYTE=384,NWORD=192)
      INTEGER*2 IBUF(NWORD),MASK
      DATA MASK /"77577/
      CHARACTER*1 CHAR(NBYT,NREC)
      EQUIVALENCE (CHAR(1,1),IBUF(1))
C-----------------------------------------------------------------------
      REAL KAPPA
      CHARACTER OUTFILE*40
      INTEGER*2 JH,JK,JL,JD(96),JT
      DIMENSION ACODE(6)
      TYPE 999
 999  FORMAT ('0TYPE IN NAME FOR OUTPUT FILE.  DEV:[DIR.SUB]FILENAME.EXT
     1')
      ACCEPT 998,OUTFILE
 998  FORMAT (A)
      OPEN (UNIT=1,FILE=OUTFILE,STATUS='NEW',FORM='UNFORMATTED')
C
C         SET ERROR HANDLER FOR INPUT CONVERSION ERROR SO THAT PROGRAM
C         DOES NOT EXIT ON EOF, EOT, OR PARITY ERROR.
C
      CALL ERRSET (64,.TRUE.,.FALSE.,.FALSE.,.TRUE.,31)
C
C         READ IN A BLOCK OF DATA.
C
      NREF=0
 100  IREF=0
 101  CALL READMT (IBUF,NBYTE)
      CALL WAITRW (ISTAT,ICNT)
      IF (ISTAT.EQ.0) GO TO 50
C
C         ISTAT .NE. 0 INDICATES A MAGNETIC TAPE ERROR CONDITION (EOF,
C         EOT, PARITY ERROR, ETC.).
C
      CALL MTERR (ISTAT)
      IF (IREF.GT.0) GO TO 51
C
C         A SHORT FILE OF HEADER RECORDS PRECEEDS EACH REFLECTION DATA
C         FILE.
C
      TYPE 997
 997  FORMAT ('0HEADER RECORDS')
      GO TO 55
C
C         THE LAST REFLECTION IN A FILE.
C
 51   TYPE 996,NREF
 996  FORMAT ('0LAST REFLECTION WRITTEN ON OUTPUT FILE HAD NREF = 'I6)
      TYPE 995,JREF,JH,JK,JL,THETA,PHIK,OMEGAK,KAPPA,RNET,SIGR,XTIME
 995  FORMAT ('      I   H  K  L     TH     PH     OM     KA       INT
     1       SIG      TIME'/ I7,1X,3I3,4F7.2,2F10.0,F10.2)
C
C         CONTINUE TO THE NEXT FILE?
C
 55   TYPE 994
 994  FORMAT ('$CONTINUE?  Y OR N.  ')
      ACCEPT 993,ANS
 993  FORMAT (A1)
      IF (ANS.NE.'Y'.AND.ANS.NE.'N') GO TO 55
      IF (ANS.EQ.'N') GO TO 9
      GO TO 100
C
C         ZERO THE PARITY BIT (BIT 8) IN EACH BYTE.  (BIT REPRESENTATION
C         OF AN ASCII CHARACTER REQUIRES ONLY SEVEN BITS.)
C
 50   DO 59 I=1,NWORD
 59   IBUF(I)=IBUF(I).AND.MASK
C
C         DECODING LOOP FOR ONE DATA BLOCK
C         -------- ---- --- --- ---- -----
C
      DO 300 I=1,NREC
C
C         DETERMINE RECORD TYPE.
C
      JTYPE=0
      DECODE (2,1400,CHAR(3,I),ERR=8) JTYPE
 1400 FORMAT (I2)
      IF (JTYPE.LT. 1) GO TO 300
      IF (JTYPE.GT.10) GO TO 300
      GO TO (401,402,403,403,403,403,403,403,403,403) JTYPE
C
C         TYPE 1 RECORD
C
  401 IERR=0
      DECODE (58,1401,CHAR(5,I),ERR=8) JREF,JH,JK,JL,(ACODE(J),J=1,6),
     1PSI,NPI,JBGL,JPK,JBGR
 1401 FORMAT (I6,3I5,1X,6A1,F7.2,I4,I6,I7,I6)
      GO TO 300
C
C         TYPE 2 RECORD
C
  402 DECODE (58,1402,CHAR(5,I),ERR=8) JREF,THETA,PHIK,OMEGAK,KAPPA,
     1WIDTH,JXTIME,JFRIDL
 1402 FORMAT (I6,F8.3,3F9.3,F7.3,I7,I3)
      GO TO 300
C
C         TYPES 3 TO 10 RECORDS
C
  403 N=12*(JTYPE-2)
      M=N-11
      DECODE (58,1403,CHAR(5,I),ERR=8) JREF,(JD(J),J=M,N)
 1403 FORMAT (I6,4(1X,3I4))
      IF (N.LT.96) GO TO 300
      IF (IERR.NE.0) GO TO 300
C
C         RECORDS FOR ONE ENTIRE REFLECTION ARE DECODED.
C
      IREF=IREF+1
C
C         THE FIRST REFLECTION IN A FILE.
C
      IF (IREF.GT.1) GO TO 69
      TYPE 992
 992  FORMAT ('0FIRST REFLECTION IN THIS FILE HAS:')
      TYPE 995,JREF,JH,JK,JL,THETA,PHIK,OMEGAK,KAPPA,RNET,SIGR,XTIME
C
C         SKIP THIS FILE?
C
 65   TYPE 991
 991  FORMAT ('$SKIP THIS FILE?  Y OR N.  ')
      ACCEPT 993,ANS
      IF (ANS.NE.'Y'.AND.ANS.NE.'N') GO TO 65
      IF (ANS.EQ.'N') GO TO 69
      CALL SKPFMT (1)
      CALL WAITMT (ISTAT,ICNT)
      IF (ISTAT.NE.0) CALL MTERR (ISTAT)
      GO TO 100
 69   CONTINUE
C
C         FLAG STANDARD REFERENCE INTENSITY REFLECTIONS BY SETTING THEIR
C         SERIAL NUMBER NEGATIVE.
C
      IF (ACODE(1).EQ.'I') JREF=-JREF
C
C         ESTIMATE NET INTEGRATED INTENSITY.
C
      BGL=JBGL
      PK=JPK
      BGR=JBGR
      RNET=PK-2*(BGL+BGR)
      VARR=PK+4*(BGL+BGR)
      SIGR=SQRT(VARR)
C
C         EVALUATE SCAN SPEED.
C
      SPEED=20.1166/NPI
C
C         NEGATIVE SCAN SPEED IS A FLAG THAT THE BEAM ATTENUATOR WAS
C         USED.
C
C         NEGATIVE WIDTH IS A FLAG THAT THE SCAN WAS DONE IN THE
C         NEGATIVE, I.E., POSITIVE TO NEGATIVE THETA, DIRECTION.  (NOTE,
C         HOWEVER, THAT THE PROFILE IS ALWAYS WRITTEN IN THE POSITIVE,
C         I.E., NEGATIVE TO POSITIVE THETA, DIRECTION.)
C
C         TRANSFER ATTENUATOR FLAG FROM SCAN SPEED TO SCAN WIDTH.
C
      WIDTH=ABS(WIDTH)
      IF (SPEED.GE.0.0) GO TO 60
      SPEED=-SPEED
      WIDTH=-WIDTH
  60  CONTINUE
C
C         NORMALIZE NET INTENSITY FOR VARIABLE SCAN SPEED.
C
      RNET=SPEED*RNET
      SIGR=SPEED*SIGR
C
C         CONVERT X-RAY EXPOSURE TIME FROM SECONDS TO HOURS.
C
      XTIME=JXTIME/3600.
C
C         REVERSE ORDER OF PROFILE ARRAY FOR REFLECTIONS MEASURED AT
C         NEGATIVE TWO-THETA.
C
      IF (THETA.LT.0) THEN
      DO 5 J=1,48
      K=96-J+1
      JT=JD(J)
      JD(J)=JD(K)
 5    JD(K)=JT
      END IF
C
C         REVERSE SIGNS OF INDICES FOR PSI-SCAN REFLECTIONS MEASURED AS
C         FRIEDEL EQUIVALENTS BUT AT NEGATIVE TWO-THETA.
C
      IF (PSI.NE.0.AND.THETA.LT.0.AND.JFRIDL.LT.0) THEN
      JH=-JH
      JK=-JK
      JL=-JL
      END IF
C
C         WRITE REFLECTION RECORD TO OUTPUT FILE.
C
      WRITE (1) JREF,JH,JK,JL,THETA,PHIK,OMEGAK,KAPPA,WIDTH,SPEED,
     1BGL,PK,BGR,RNET,SIGR,XTIME,JD
      NREF=NREF+1
      GO TO 300
C
C         SET ERROR FLAG.
C
 8    IERR=1
C
C         END OF THE DECODING LOOP FOR ONE DATA BLOCK
C         --- -- --- -------- ---- --- --- ---- -----
C
 300  CONTINUE
      GO TO 101
C
C         END OUTPUT FILE AND EXIT.
C
 9    ENDFILE 1
      CALL CLOSE_MT
      TYPE 990,NREF
 990  FORMAT ('0NREF = 'I5' REFLECTIONS TRANSFERRED TO OUTPUT FILE')
      STOP 'PROGRAM CAD4TAPE FINIS!'
      END
      SUBROUTINE MTERR (E_CODE)
C
C     MAGTAPE ERROR MESSAGE PROCESSING.
C
      INTEGER*2 MESS_L
      INTEGER*4 E_CODE
      CHARACTER E_MESS*256
      IF(E_CODE.EQ.0)RETURN
      CALL SYS$GETMSG(%VAL(E_CODE),MESS_L,E_MESS,%VAL('F'X),)
      TYPE 100,E_MESS(1:MESS_L)
 100  FORMAT('0',A)
      RETURN
      END
