*  subroutines used for UCAC2 release utility programs
*
*  open_zfile    : open direct access, unformatted zone file
*  read_u2line   : read single record = all items for a star
*  cat_id        : identify proper motion catalogs from flag
*  flip2         : flip 2 byte integer
*  flip4         : flip 4 byte integer
*  valid_range   : restrict R*8 data item to given min,max 
*  get_zone_range: declination range --> required zone numbers
*  get_ra_ragne  : RA range --> required index for 0.1h bins
*  chk_byte_flip : read first record of z001, is a byte flip required?
*  nx_byte_flip  : check / apply byte flip on index array
*  count_id      : count proper motion catalogs from flag
*  as2hms        : convert arcsec (RA,Dec) into hms format
*
*  030528 change "rflg" to "epos", remove "READONLY" option in OPEN
*  030529 fix RA range 24/0h
*
************************************************************************

      SUBROUTINE open_zfile (pathz,un,zn,only_rd)
C
C  input : pathz = path name for zone files
C          un    = Fortran file unit number
C          zn    = zone number = 1, 288
C          only_rd = .TRUE. if only read access

      IMPLICIT NONE
      CHARACTER*(*) pathz
      INTEGER  un,zn, jp
      LOGICAL  only_rd, ifex
      CHARACTER*80 fnzone, answer

      IF (zn.LT.1.OR.zn.GT.288) THEN
        WRITE (*,'(a,i5)') '<open_zf> invalid zone number = ',zn
        STOP
      ENDIF

      IF (only_rd) THEN               ! read file, check existence
 51     jp = INDEX (pathz,' ') - 1
        WRITE (fnzone,'(a,a,i3.3)') pathz(1:jp),'z',zn

        INQUIRE (FILE=fnzone,EXIST=ifex)

        IF (.NOT.ifex) THEN
          WRITE (*,'(/a)') 'can not find the file:'
          WRITE (*,'(a)') fnzone
          WRITE (*,'(a)') 'please use correct CD or enter new path:'
          WRITE (*,'(a)') '(exit with "x")'
          READ (*,'(a)') answer
          IF (answer.NE.' ') pathz = answer
          IF (pathz(1:1).EQ.'x'.AND.pathz(2:2).EQ.' ') THEN
            STOP
          ELSE
            GOTO 51
          ENDIF
        ENDIF        ! file does not exist
 
CC      OPEN (un,FILE=fnzone,ACCESS='direct',RECL=44,READONLY)
        OPEN (un,FILE=fnzone,ACCESS='direct',RECL=44) 

      ELSE
        jp = INDEX (pathz,' ') - 1
        WRITE (fnzone,'(a,a,i3.3)') pathz(1:jp),'z',zn
        OPEN (un,FILE=fnzone,ACCESS='direct',RECL=44)
      ENDIF ! read or write access

      END   ! subr. <open_zfile>

************************************************************************

      SUBROUTINE read_u2line (un,recn,bf,idat,errflg)
C
C  read a single record of UCAC2 data = 1 star
C  input:
C    un   = unit number of file  (assumed to be open)
C    recn = record number on that file
C    bf   = .TRUE. if byte flip required
C  output:
C    idat = integer*4 vector of 23 items (see readme2.txt)
C    errflg = true  if error occured (like end of file)
C             else false

      IMPLICIT NONE
      INTEGER  un,recn, idat(25)  ! item #24,25 = star ID options
      LOGICAL  bf, errflg 

      INTEGER   ra2000, dc2000, pmx,pmy, id2m, u2id,r11
      INTEGER*2 mag, cepx,cepy, j2m,h2m,k2m
      BYTE      sigx,sigy,nobs,epos,ncat,cflg    ! INTEGER*1
     .         ,spmx,spmy, rx,ry, ph,cc          ! signed integer

      errflg = .FALSE.     ! default

CC    READ (un,REC=recn,ERR=99) ra2000,dc2000
CC   .  ,mag,sigx,sigy, nobs,rflg,ncat,cflg
CC   .  ,cepx,cepy, pmx,pmy, spmx,spmy, rx,ry
CC   .  ,id2m, j2m,h2m,k2m, ph,cc, u2id,r11   ! incl. ID numbers

      READ (un,REC=recn,ERR=99) ra2000,dc2000
     .  ,mag,sigx,sigy, nobs,epos,ncat,cflg
     .  ,cepx,cepy, pmx,pmy, spmx,spmy, rx,ry
     .  ,id2m, j2m,h2m,k2m, ph,cc

      IF (bf) THEN
        CALL flip4 (ra2000)
        CALL flip4 (dc2000)
        CALL flip2 (mag)
        CALL flip2 (cepx)
        CALL flip2 (cepy)
        CALL flip4 (pmx)
        CALL flip4 (pmy)
        CALL flip4 (id2m)
        CALL flip2 (j2m)
        CALL flip2 (h2m)
        CALL flip2 (k2m)
CC      CALL flip4 (u2id)
CC      CALL flip4 (r11)
      ENDIF

* note: first assign I*1 to idat(I*4), 
*       then add 127 to avoid overflow

      idat ( 1) = ra2000
      idat ( 2) = dc2000
      idat ( 3) = mag
      idat ( 4) = sigx
      idat ( 4) = idat ( 4) + 127
      idat ( 5) = sigy
      idat ( 5) = idat ( 5) + 127
      idat ( 6) = nobs
      idat ( 7) = epos 
      idat ( 7) = idat ( 7) + 127
      idat ( 8) = ncat
      idat ( 9) = cflg
      idat (10) = cepx
      idat (11) = cepy
      idat (12) = pmx
      idat (13) = pmy
      idat (14) = spmx
      idat (14) = idat (14) + 127
      idat (15) = spmy
      idat (15) = idat (15) + 127
      idat (16) = rx
      idat (16) = idat (16) + 127
      idat (17) = ry
      idat (17) = idat (17) + 127
      idat (18) = id2m
      idat (19) = j2m
      idat (20) = h2m
      idat (21) = k2m
      idat (22) = ph 
      idat (22) = idat (22) + 127
      idat (23) = cc
      idat (23) = idat (23) + 127
CC    idat (24) = u2id  ! option for test runs
CC    idat (25) = r11   ! including cross references
      idat (24) = 0     ! here don't use item 24,25
      idat (25) = 0     ! but keep data structure
      RETURN

 99   errflg = .TRUE.
      END   ! subr. <read_u2line>

************************************************************************

      SUBROUTINE cat_id (cflg,icat)
C
C  input : cflg   = combined flag for catalog ID's
C  output: icat(7)= 1 if catalog (1 to 7) is included, else 0

      IMPLICIT NONE
      INTEGER  cflg, icat(7), cc, j

      cc = cflg

      DO j=1,7
        icat(j) = 0
      ENDDO

      IF (cc.GE.64) THEN         ! USNO-A2
        icat(7) = 1
        cc = cc - 64
      ENDIF

      IF (cc.GE.32) THEN         ! NLTT
        icat(6) = 1
        cc = cc - 32
      ENDIF

      IF (cc.GE.16) THEN         ! Hipparcos
        icat(5) = 1
        cc = cc - 16
      ENDIF

      IF (cc.GE.8) THEN          ! AGK2 
        icat(4) = 1
        cc = cc - 8
      ENDIF

      IF (cc.GE.4) THEN          ! Tycho-2
        icat(3) = 1
        cc = cc - 4
      ENDIF

      IF (cc.GE.2) THEN          ! AC2000
        icat(2) = 1
        cc = cc - 2
      ENDIF

      IF (cc.EQ.1) icat(1) = 1   ! YS

      END  ! subr. <cat_id>

************************************************************************

      SUBROUTINE flip2 (i2)
C
C input:  Integer*2 value i2
C output: same with byte fliped

      IMPLICIT NONE
      INTEGER*2  i2, in, out
      BYTE       a(2), b(2)
      EQUIVALENCE (in,a)
      EQUIVALENCE (out,b)

      in = i2
      b(1) = a(2)
      b(2) = a(1)
      i2 = out

      END    ! subr. <flip2>
 
************************************************************************

      SUBROUTINE flip4 (i4)
C
C input:  Integer*4 value i4
C output: same with byte fliped

      IMPLICIT NONE
      INTEGER*4  i4, in, out
      BYTE       a(4), b(4)
      EQUIVALENCE (in,a)
      EQUIVALENCE (out,b)

      in = i4
      b(1) = a(4)
      b(2) = a(3)
      b(3) = a(2)
      b(4) = a(1)
      i4 = out

      END    ! subr. <flip4>

************************************************************************

      SUBROUTINE valid_range (data,dmin,dmax)
C
      IMPLICIT NONE
      REAL*8   data, dmin,dmax

      IF (data.LT.dmin) data = dmin
      IF (data.GT.dmax) data = dmax

      END    ! subr. <valid_range>

************************************************************************

      SUBROUTINE get_zone_range (dc1,dc2,zmax, d1m,d2m,z1,z2,nz)
C
C input:  dc1,dc2 = declination range (degree)
C         zmax    = largest zone number available
C output: d1m,d2m = declination range in mas
C         z1, z2  = req. range of zone numbers (0.5 deg steps)
C         nz      = number of zones, or 0 if out of range

      IMPLICIT NONE
      REAL*8  dc1,dc2
      INTEGER zmax, d1m,d2m, z1,z2, nz
      REAL*8  dcx 

      IF (dc1.LT.-90.0d0.AND.dc2.LT.-90.0d0) THEN
        nz = 0
        z1 = 1
        z2 = 0
        RETURN
      ENDIF

      CALL valid_range (dc1,-90.0d0,90.0d0)
      CALL valid_range (dc2,-90.0d0,90.0d0)

      IF (dc1.GT.dc2) THEN     ! flip range
        dcx = dc1
        dc1 = dc2
        dc2 = dcx
      ENDIF

      d1m = IDNINT (dc1 * 3.6d6)    ! declination (mas)
      d2m = IDNINT (dc2 * 3.6d6)

      z1 = (d1m + 324000000) / 1800000 + 1
      z2 = (d2m + 323999999) / 1800000 + 1

      IF (z2.GT.zmax) z2 = zmax

      IF (z1.GT.zmax) THEN      ! out of available zone range
        z1 = zmax + 1
        nz = 0  
      ELSE
        nz = z2 - z1 + 1
      ENDIF

      END   !  subr. <get_zone_range>

************************************************************************

      SUBROUTINE get_ra_range (ra1,ra2, ralo,rahi,i1,i2,nr)
C
C  input:  ra1,ra2   = RA range (hour)
C  output: ralo,rahi = range of RA in mas (1 or 2)
C          i1, i2    = range in index for 0.1 h boxes
C          nr        = number of ranges = 1 or 2
C    2 ranges possible, if ra1 > ra2  (e.g. 23.0, 1.0)
C    assume cross over 24/0 hour in RA --> 2 ranges
C    (like  23.0 ... 24.0  and  0.0 ... 1.0 hour for output)

      IMPLICIT NONE
      REAL*8  ra1,ra2
      INTEGER ralo(2),rahi(2), i1(2),i2(2), nr
      INTEGER r1m,r2m
      REAL*8  rax

      CALL valid_range (ra1, 0.0d0,24.0d0)
      CALL valid_range (ra2, 0.0d0,24.0d0)

      r1m = IDNINT (ra1 * 5.4d7)   ! RA in mas
      r2m = IDNINT (ra2 * 5.4d7)

      IF (r1m.LE.r2m) THEN         ! normal case
        nr = 1
        i1(1) =  r1m    / 5400000 + 1
        i2(1) = (r2m-1) / 5400000 + 1
        i1(2) = 1
        i2(2) = 0
        ralo(1) = r1m
        rahi(1) = r2m
        ralo(2) = 0
        rahi(2) = 0

      ELSE                         ! cross over 24/0
        nr = 2
        i1(1) =  r1m    / 5400000 + 1
        i2(1) = 240
        i1(2) =   1
        i2(2) = (r2m-1) / 5400000 + 1
        ralo(1) = r1m
        rahi(1) = 1296000000       ! 24 hour in mas
        ralo(2) = 0
        rahi(2) = r2m
      ENDIF

      END   ! subr. <get_ra_range>

************************************************************************

      SUBROUTINE chk_byte_flip (pathz,un,bf)
C
C  input : pathz  = path for zone files
C          un     = Fortran unit number for zone file
C  output: bf     = .TRUE. if byte flip is required

      IMPLICIT NONE
      CHARACTER*(*) pathz
      INTEGER       un, zn, idat(25)
      LOGICAL       bf, bft, errflg, only_rd
      CHARACTER*40 fnz, a1*1
      INTEGER*2    mag

      only_rd = .TRUE.
      zn = 1

      WRITE (*,'(/a)') 'open first zone, read first record' 

      CALL open_zfile (pathz,un,zn,only_rd)

      bft = .FALSE.                 ! first test with no byte flip

      CALL read_u2line (un,1,bft,idat,errflg)
      CLOSE (un)

      WRITE (*,'(a,i6)') 'mag of first star = ',idat(3)

      IF (idat(3).EQ.1591) THEN     ! magnitude of first star
        bf = .FALSE.
        WRITE (*,'(/a)') '-- no byte flip required'

      ELSE
        mag = idat(3)
        CALL flip2 (mag)
        IF (mag.EQ.1591) THEN
          WRITE (*,'(/a)') '** byte flip is required, will do'
          bf = .TRUE.
        ELSE
          WRITE (*,'(/a)') '** WARNING: byte flip test inconclusive'
          bf = .FALSE.
        ENDIF
      ENDIF

      WRITE (*,'(a,$)') 'hit "enter" to continue '
      READ (*,'(a)') a1

      END   ! subr. <chk_byte_flip>

************************************************************************

      SUBROUTINE nx_byte_flip (nx,zmax,bf)
C
C input : nx  = array with index
C         zmax= dimension of nx, max. number of zones
C output: nx  = same with byte flip applied (if required)
C         bf  = .TRUE. if byte flip was applied

      INTEGER zmax
      INTEGER nx (zmax,240)
      LOGICAL bf
      INTEGER zn,j, i4
      CHARACTER*1 a1

      IF (nx(1,1).EQ.2) THEN
        WRITE (*,'(a)') 'index array: no byte flip required'
        bf = .FALSE.

      ELSEIF (nx(1,1).EQ.33554432) THEN
        WRITE (*,'(a)') 'index array: byte flip is required'
        bf = .TRUE. 
        DO zn= 1,zmax
        DO j = 1,240
          CALL flip4 (nx(zn,j))
        ENDDO
        ENDDO
        WRITE (*,'(a)') 'index array: byte flip done' 
 
      ELSE
        WRITE (*,'(a)') 
     .   'WARNING: index array: byte flip inconclusive'
        bf = .FALSE.
      ENDIF

      WRITE (*,'(a,$)') 'hit "enter" to continue '
      READ (*,'(a)') a1

      END   ! subr. <nx_byte_flip>

************************************************************************

      SUBROUTINE count_id (cflg,ncat)
C
C  input : cflg   = combined flag for catalog ID's
C  output: ncat(7)= count for 7 catalogs coded in cflg

      IMPLICIT NONE
      INTEGER  cflg, ncat(7), cc

      cc = cflg

      IF (cc.GE.64) THEN         ! USNO-A2
        ncat(7) = ncat(7) + 1
        cc = cc - 64
      ENDIF

      IF (cc.GE.32) THEN         ! NLTT
        ncat(6) = ncat(6) + 1
        cc = cc - 32
      ENDIF

      IF (cc.GE.16) THEN         ! Hipparcos
        ncat(5) = ncat(5) + 1
        cc = cc - 16
      ENDIF

      IF (cc.GE.8) THEN          ! AGK2 
        ncat(4) = ncat(4) + 1
        cc = cc - 8
      ENDIF

      IF (cc.GE.4) THEN          ! Tycho-2
        ncat(3) = ncat(3) + 1
        cc = cc - 4
      ENDIF

      IF (cc.GE.2) THEN          ! AC2000
        ncat(2) = ncat(2) + 1
        cc = cc - 2
      ENDIF

      IF (cc.EQ.1) ncat(1) = ncat(1) + 1   ! YS

      END  ! subr. <count_id>

************************************************************************

      SUBROUTINE as2hms (ra,dk,crekt,cdekl)
C
C convert R*8 RA, DC (arcsec) to hms, dms strings
C
C 940725 NZ update to CHARACTER*13 to 1/1000 arcsec

        IMPLICIT REAL*8 (A-H,L-Z)

        REAL*8     RA, DK                            ! added 2-3-93
        INTEGER*4  IRASTD, IRAMIN, IDKGRD,IDKMIN     ! added 2-3-93
        CHARACTER*1  CVZ
        CHARACTER*13 CREKT,CDEKL


C     03. TRANSFORMATIONEN

        IF (RA.GT.1296000.D0)  RA = RA - 1296000.D0
        IF (RA.LT.      0.D0)  RA = RA + 1296000.D0
        IF (RA.GT.1296000.D0.OR.RA.LT.0.D0)  THEN
          WRITE (90,'(1X//1X,A,F13.3/)')  'RA > 24 hours or < 0 :',RA
          RETURN
        END IF
        RASTD = RA/(3600.D0*15.D0)
        IRASTD= IDINT(RASTD)
        RAREST= RASTD-DFLOAT(IRASTD)
        IRAMIN= IDINT(RAREST*60.D0)
        RASEC = RAREST*3600.D0-DFLOAT(IRAMIN)*60.D0
        IF (DABS(RASEC-60.D0).LT.0.001D0)  THEN
          RASEC = 0.D0
          IRAMIN= IRAMIN+1
          IF (IRAMIN.EQ.60)  THEN
            IRAMIN= 0
            IRASTD= IRASTD+1
            IF (IRASTD.EQ.24)  IRASTD= 0
          END IF
        END IF

        IF (DABS(DK).GT.324000.D0)  THEN
          WRITE (90,'(1X//1X,A,F13.3/)')  'abs (DC)  > 90 deg :',DK
          RETURN
        END IF
        DKGRD = DK/3600.D0
        CVZ= '+'
        IF (DK.LT.0.D0)  THEN
          CVZ= '-'
          DKGRD= -DKGRD
        END IF
        IDKGRD= IDINT(DKGRD)
        DKREST= DKGRD-DFLOAT(IDKGRD)
        IDKMIN= IDINT(DKREST*60.D0)
        DKSEC = DKREST*3600.D0-DFLOAT(IDKMIN)*60.D0
        IF (DABS(DKSEC-60.D0).LT.0.01D0)  THEN
          DKSEC = 0.D0
          IDKMIN= IDKMIN+1
          IF (IDKMIN.EQ.60)  THEN
            IDKMIN= 0
            IDKGRD= IDKGRD+1
          END IF
        END IF

        WRITE (CREKT,'(   I2.2,1X,I2.2,1X,F7.4)')
     A    IRASTD,IRAMIN,RASEC
        IF (CREKT(7:7).EQ.' ')  CREKT(7:7)= '0'
        WRITE (CDEKL,'(A1,I2.2,1X,I2.2,1X,F6.3)')
     A    CVZ,IDKGRD,IDKMIN,DKSEC
        IF (CDEKL(8:8).EQ.' ')  CDEKL(8:8)= '0'

        RETURN
      END
************************************************************************
