SUBROUTINE VREAD C C Subroutine to read data from vis file, including calibration tables C C VERSION 3.1 14 OCTOBER 1988 C C 25 OCT 1988 corrected error checking for scans with no data C IMPLICIT UNDEFINED (A-Z) SAVE CHARACTER *11 VERSN, VERSD PARAMETER ( VERSN = 'VERSION 3.1' ) PARAMETER ( VERSD = '23 MAR 1989' ) INTEGER *2 I2SCAN, I2OK INTEGER *4 NSUM, IERR, N, I, J, K, L, VCNT, NPTS, ILAST, IOK INTEGER *4 NDATA(0:4), NDARK, NAME(4), SPECTRUM, ILEN, IFIRST INTEGER *4 ISIZE, IFMT REAL *4 DEN, PAR, RV, PMRA, PMDEC, COSMA, TEMP(16), NEWDIAM REAL *4 NEWMAG, NEWBV REAL *8 RAM, DECM INCLUDE 'VPLOT.INC' CHARACTER * 1 ANS CHARACTER * 3 MNTH(12) LOGICAL NOPCAL, NOSCAL, INSCAN C C The input data format C INTEGER *2 IA(38), ICODE, ISTAR, IBASE, TAU0, DTAU0 INTEGER *2 TSAMPLE, IMON, IDAY, IYEAR, JITTR, ILOCK INTEGER *4 IDARK(4), IFILT(4), DJTR REAL *8 TIME0, DLAS REAL *4 ANUM(6), STROKE, ADEN(4), AVAR(4), HA0, DEC0, ZD0 REAL *4 TFIRST, TLAST, PHI12, PHIRMS C C DATA REC TIME REC DARK REC EQUIVALENCE (IA( 1), TIME0 ), $ (IA( 5), ANUM(1), ICODE ), $ (IA( 6), ISTAR ), $ (IA( 7), IBASE ), $ (IA( 8), TAU0, DTAU0 ), $ (IA( 9), TSAMPLE, DLAS ), $ (IA(10), IMON ), $ (IA(11), IDAY ), $ (IA(12), IYEAR ), $ (IA(13), STROKE, IDARK(1)), $ (IA(15), IFILT(1) ) EQUIVALENCE (IA(17), ADEN(1) ), $ (IA(21), HA0 ), $ (IA(23), DEC0 ), $ (IA(25), AVAR(1), ZD0 ), $ (IA(27), DJTR ), $ (IA(29), TFIRST ), $ (IA(31), TLAST ), $ (IA(33), JITTR, ILOCK ), $ (IA(35), PHIRMS ), $ (IA(37), PHI12 ) C********************************************************************** DATA MNTH / 'JAN', 'FEB', 'MAR', 'APR', 'MAY', 'JUN', $ 'JUL', 'AUG', 'SEP', 'OCT', 'NOV', 'DEC' / C---------------------------------------------------------------------- C Ask user for and open the visibility file. C Decide if the calibration files exist C---------------------------------------------------------------------- 10 CONTINUE WRITE (6,5004) 5004 FORMAT ( ' Filename (must be of form YYMMDD[.vis]): ', $ ) READ (5,'(A)',IOSTAT=IERR) FNAME IF ( IERR .NE. 0 ) THEN WRITE(6,*) ' INPUT ERROR ' FNAME = ' ' RETURN END IF ILEN = LEN(FNAME) IROOT = ILEN ILAST = ILEN DO 50 I = ILEN, 1, -1 IF ( FNAME(I:I) .EQ. '.' ) IROOT = I - 1 IF ( FNAME(I:I) .EQ. ' ' ) ILAST = I - 1 50 CONTINUE IROOT = MIN ( IROOT, ILAST ) IF ( ILAST .EQ. 0 ) THEN WRITE(6,*) ' FILE NAME MUST NOT START WITH A BLANK ' RETURN ELSE IF ( IROOT + 4 .GT. ILEN ) THEN WRITE(6,*) ' FILE NAME IS TOO LONG ' GO TO 10 END IF IF ( IROOT .EQ. ILAST ) THEN ILAST = IROOT + 4 FNAME(IROOT+1:ILAST) = '.VIS' END IF OPEN (UNIT=4, FILE=FNAME(1:ILAST), STATUS='OLD', $ ACCESS='DIRECT', FORM='UNFORMATTED', RECL=76, $ IOSTAT=IERR ) IF ( IERR .NE. 0 ) THEN WRITE(6,'(A,I4,2A/A)') $ ' Error ', IERR, ' opening file ', FNAME, $ ' file not found. ' RETURN ELSE WRITE(6,*) ' data file ', FNAME(1:ILAST), ' opened ' END IF chummel: remove path from name to write to current directory fname=fname(iroot-5:iroot) iroot=6 SNAME = FNAME(1:IROOT) // '.SCN' PNAME = FNAME(1:IROOT) // '.PNT' C C open scan calibration table. C OPEN (UNIT=3, FILE=SNAME(1:ILAST), STATUS='OLD', $ ACCESS='DIRECT', FORM='UNFORMATTED', RECL=64, $ IOSTAT=IERR ) IF ( IERR .NE. 0 ) THEN WRITE(6,'(A)') ' Cannot find scan calibration. ' WRITE(6,'(A)') ' Scan calibration table will be generated ' OPEN (UNIT=3, FILE=SNAME(1:ILAST), STATUS='NEW', $ ACCESS='DIRECT', FORM='UNFORMATTED', RECL=64, $ IOSTAT=IERR ) NOSCAL = .TRUE. IF ( IERR .NE. 0 ) THEN WRITE(6,*) ' Cannot create new file for Calibration', $ ' table.' WRITE(6,*) ' Disk is probably full or write protected' STOP END IF ELSE NOSCAL = .FALSE. WRITE(6,'(A)') ' Scan calibration table opened ' END IF C C open point calibration table. C OPEN (UNIT=2, FILE=PNAME(1:ILAST), STATUS='OLD', $ ACCESS='DIRECT', FORM='UNFORMATTED', RECL=4, $ IOSTAT=IERR ) IF ( IERR .EQ. 0 ) THEN WRITE(6,'(A)') ' Point calibration table opened ' NOPCAL = .FALSE. ELSE WRITE(6,'(A)') ' Cannot find point calibration. ' WRITE(6,'(A)') ' Point calibration table will be generated ' OPEN (UNIT=2, FILE=PNAME(1:ILAST), STATUS='NEW', $ ACCESS='DIRECT', FORM='UNFORMATTED', RECL=4, $ IOSTAT=IERR ) NOPCAL = .TRUE. IF ( IERR .NE. 0 ) THEN WRITE(6,*) ' Cannot create new file for Calibration', $ ' table.' WRITE(6,*) ' Disk is probably full or write protected' STOP END IF END IF C C Reset a few variables, so that another file can be read. C DAY = 0 MONTH = 0 YEAR = 0 C----------------------------------------------------------------------- C Generate the SCAN calibration table, if necessary. C----------------------------------------------------------------------- IF ( NOSCAL ) THEN C C Read the data header to determine the coherent int time. C READ ( 4, REC=1, IOSTAT=IERR) IA IF ( IERR .NE. 0 ) THEN WRITE(6,*) ' Read error record 1 of data file ' STOP END IF IF ( TIME0 .NE. 0.D0 ) THEN WRITE(6,*) ' First record of data file is not header ' STOP END IF IFMT = IA(5) IF ( IFMT .EQ. 2 ) THEN COHINT = IA(8) ELSE COHINT = 4 END IF WRITE(6,'(A,I5,A,I5)') ' VIS FORMAT = ', IFMT, $ ' COHERENT INT TIME = ', COHINT DO I = 1, 16 TEMP(I) = 0 END DO DO I = 1, 4 NDFRATIO(I) = 1.0 END DO FORMATID = 1 WRITE(3,REC=1) (TEMP(J),J=1,8), FORMATID, NDFRATIO, $ COHINT, (TEMP(J),J=15,16) DO 85 I = 1, NFIT DO 80 J = 1, 4 TCAL(I,J) = 1. HCAL(I,J) = 1. ZCAL(I,J) = 1. MCAL(I,J) = 1. 80 CONTINUE 85 CONTINUE N = 0 ISCAN = 0 ILAST = 0 NDARK = 0 DO 90, J = 1, 4 SCAL(1,J) = 1. FILT(J) = 0 90 CONTINUE INSCAN = .FALSE. WRITE(6,'(A)') ' Generating scan calibration table' 100 CONTINUE N = N + 1 READ ( 4, REC=N, IOSTAT=IERR) IA IF ( IERR .NE. 0 ) THEN WRITE(6,*) ' Error ', IERR, ' encountered on read ' WRITE(6,*) N, ' records read ' IF (INSCAN .AND. (ILAST .EQ. 0) ) THEN ISCAN = ISCAN - 1 ELSE IF ( INSCAN ) THEN I2SCAN = ISCAN I2OK = 0 WRITE(3,REC=ISCAN+1) I2SCAN, I2OK, IFIRST, ILAST, $ NDARK, 0., 0., 0., 0., (SCAL(1,J),J=1,4), $ 0., 0., 0., 0. END IF GO TO 120 END IF IF ( NOPCAL ) THEN WRITE(2,REC=N) 0 END IF C Header record IF ( TIME0.EQ.0.D0 ) THEN C WRITE(6,*) ' Record ', N, ' is header record' IF (INSCAN .AND. (ILAST .EQ. 0) ) THEN ISCAN = ISCAN - 1 ELSE IF ( INSCAN ) THEN I2SCAN = ISCAN I2OK = 0 WRITE(3,REC=ISCAN+1) I2SCAN, I2OK, IFIRST, ILAST, $ NDARK, 0., 0., 0., 0., (SCAL(1,J),J=1,4), $ 0., 0., 0., 0. END IF INSCAN = .FALSE. IFMT = IA(5) ISIZE = IA(6) IF ( IFMT .EQ. 2 ) THEN COHINT = IA(8) ELSE COHINT = 4 END IF WRITE(6,1120) N, COHINT 1120 FORMAT ( ' Record ', I5, ' is header record. ', $ ' COHERENT INTEG. TIME =', I5, ' ms ' ) C Data record ELSE IF ( TIME0 .GT. 0.D0 ) THEN IF ( INSCAN ) THEN ILAST = N ELSE WRITE(6,*) ' Data record ', N, ' is not in a scan' END IF C Time/new scan block ELSE IF ( ICODE .EQ. 0 ) THEN IF ( INSCAN .AND. (ILAST.EQ.0) ) THEN ISCAN = ISCAN - 1 ELSE IF ( INSCAN ) THEN I2SCAN = ISCAN I2OK = 0 WRITE(3,REC=ISCAN+1) I2SCAN, I2OK, IFIRST, ILAST, $ NDARK, 0., 0., 0., 0., (SCAL(1,J),J=1,4), $ 0., 0., 0., 0. END IF INSCAN = .TRUE. ISCAN = ISCAN + 1 IFIRST = N ILAST = 0 NDARK = 0 IF ((DAY.EQ.0).OR.(MONTH.EQ.0).OR.(YEAR.EQ.0)) THEN DAY = IDAY MONTH = IMON YEAR = IYEAR ELSE IF ((DAY.NE.IDAY).OR.(MONTH.NE.IMON) $ .OR.( YEAR.NE.IYEAR) ) THEN WRITE(6,*) ' Date changes in mid file ' WRITE(6,*) ' Was ', YEAR, MONTH, DAY, $ ' is ', IYEAR, IMON, IDAY WRITE(6,*) ' Check file integrity' ISCAN = ISCAN - 1 WRITE(6,*) ' Data past scan number ', ISCAN, $ ' ignored ' GO TO 120 END IF DO 105 J = 1, 4 IF ( IFILT(J) .EQ. 0 ) GO TO 105 IFILT(J) = MAX(0, IFILT(J) ) IF ( FILT(J) .EQ. 0 ) FILT(J) = IFILT(J) IF ( FILT(J) .NE. IFILT(J) ) THEN 104 CONTINUE WRITE(6,1200) FILT, IFILT, N READ ( 5,'(A)') ANS CALL CAPS ( ANS ) IF ( ANS .EQ. 'O' ) GO TO 106 IF ( ANS .EQ. 'N' ) THEN FILT(1) = IFILT(1) FILT(2) = IFILT(2) FILT(3) = IFILT(3) FILT(4) = IFILT(4) GO TO 106 END IF GO TO 104 END IF 105 CONTINUE 106 CONTINUE 1200 FORMAT ( ' filters changed in mid file ' / $ ' old filter list = ', 4I4 / $ ' new filter list = ', 4I4 / $ ' problem is with record number ', I5, // $ ' Should I use new or old list (N/O)?' ) C Dark data block ELSE IF ( ICODE .EQ. 1 ) THEN IF ( INSCAN .AND. (ILAST.EQ.0) ) THEN ISCAN = ISCAN - 1 ELSE IF ( INSCAN ) THEN NDARK = N C C DTAU0 is integration time for the dark counts in ms. The dark C count is converted to the 4 ms integration times of the fringe data. C DO 110 J = 1, 4 D(ISCAN,J) = COHINT * IDARK(J) / FLOAT(DTAU0) 110 CONTINUE I2SCAN = ISCAN I2OK = 0 WRITE(3,REC=ISCAN+1) I2SCAN, I2OK, IFIRST, ILAST, $ NDARK, (D(ISCAN,J),J=1,4), (SCAL(1,J),J=1,4), $ 0., 0., 0., 0. ELSE WRITE(6,*) ' Extra dark data record ', N END IF INSCAN = .FALSE. ELSE WRITE(6,*) ' Record ', N, ' is illegal format ' WRITE(6,*) ' *** Bozo city *** ' STOP END IF GO TO 100 120 CONTINUE C C Include header record with date and filter numbers C I2SCAN = ISCAN I2OK = 0 FORMATID = 1 WRITE(3,REC=1) I2SCAN, I2OK, YEAR, MONTH, DAY, $ (FILT(J),J=1,4), FORMATID, $ (NDFRATIO(J),J=1,4), 0., 0., 0. WRITE(6,*) ' FILTER NUMBERS ARE ', FILT C C Add the information records C C WRITE(3,REC=ISCAN+2) NFIT, TCAL0, TCAL1, HCAL0, HCAL1, $ ZCAL0, ZCAL1, MCAL0, MCAL1, $ 0., 0., 0., 0., 0., 0., 0. WRITE(3, REC=ISCAN+3) 0,0,0,0, 0.,100.,0,0, 0,0,0,0, 0,0,0,0 WRITE(3, REC=ISCAN+4) 0,0,0,0, 0., 0.,0,0, 0,0,0,0, 0,0,0,0 WRITE(3, REC=ISCAN+5) 0,0,0,0, 0., 0.,0,0, 0,0,0,0, 0,0,0,0 C C Add the calibration arrays C DO 130 I = 1, NFIT WRITE(3,REC=I+5+ISCAN, IOSTAT=IERR) $ ( TCAL(I,J),J=1,4 ), ( HCAL(I,J),J=1,4 ), $ ( ZCAL(I,J),J=1,4 ), ( MCAL(I,J),J=1,4 ) IF ( IERR .NE. 0 ) THEN WRITE(6,*) ' FATAL ERROR ', IERR, ' WRITING SCAN ', $ ' CALIBRATION TABLE ' WRITE(6,*) ' FILE CREATION NOT SUCESSFUL ' END IF 130 CONTINUE C----------------------------------------------------------------------- C Generate just the point calibration table C----------------------------------------------------------------------- ELSE IF ( NOPCAL ) THEN N = 0 133 CONTINUE N = N + 1 READ ( 4, REC=N, IOSTAT=IERR) IA IF ( IERR .NE. 0 ) GO TO 134 WRITE(2,REC=N) 0 GO TO 133 134 CONTINUE END IF C----------------------------------------------------------------------- C Read in data, using the calibration tables. C The first record is a header record. C----------------------------------------------------------------------- CALL SCANCALC C----------------------------------------------------------------------- C Create a list of observed stars and baselines. Look up each star C in the catalog to determine the magnitude and angular diameter. C Then compute the estimated visibility for each observation. C C Note that the uniform disk angular diameter is a function of C wavelength. The values in DIAMETER.DAT are for 800 nm. C QCAL scales the diameter from 800 nm to the appropriate wavelength. C----------------------------------------------------------------------- NSTAR = 1 NBASE = 1 BLIST(NBASE) = BASE(1) SLIST(NSTAR) = STAR(1) CALL CATALOG( SLIST(NSTAR), NAME, RAM, DECM, PMRA, PMDEC, $ PAR, RV, MAG(NSTAR), BV(NSTAR), SPECTRUM, $ ANGDIA(NSTAR) ) WRITE(STARNAME(NSTAR), '(2A4)') NAME(1), NAME(2) WRITE(6,2050) NSTAR, SLIST(NSTAR), STARNAME(NSTAR) IF (ANGDIA(NSTAR) .GT. 1. ) THEN ANGDIA(NSTAR) = 0. ELSE ANGDIA(NSTAR) = 1000. * ANGDIA(NSTAR) END IF DO 200 J = 1, 4 CALL QCAL ( BASE(1), HA(1), DEC(1), WAVE(J), $ ANGDIA(1), VEST(1,J) ) 200 CONTINUE DO 250 I = 1, ISCAN DO 210 J = 1, NBASE IF ( BASE(I) .EQ. BLIST(J) ) GO TO 220 210 CONTINUE NBASE = NBASE + 1 BLIST(NBASE) = BASE(I) 220 CONTINUE DO 230 J = 1, NSTAR K = J IF ( STAR(I) .EQ. SLIST(J) ) GO TO 235 230 CONTINUE NSTAR = NSTAR + 1 SLIST(NSTAR) = STAR(I) K = NSTAR CALL CATALOG( SLIST(NSTAR), NAME, RAM, DECM, PMRA, PMDEC, $ PAR, RV, MAG(NSTAR), BV(NSTAR), SPECTRUM, $ ANGDIA(NSTAR) ) WRITE(STARNAME(NSTAR), '(2A4)') NAME(1), NAME(2) WRITE(6,2050) NSTAR, SLIST(NSTAR), STARNAME(NSTAR) IF (ANGDIA(NSTAR) .GT. 1. ) THEN ANGDIA(NSTAR) = 0. ELSE ANGDIA(NSTAR) = 1000. * ANGDIA(NSTAR) END IF 235 CONTINUE DO 240 J = 1, 4 CALL QCAL ( BASE(I), HA(I), DEC(I), WAVE(J), $ ANGDIA(K), VEST(I,J) ) 240 CONTINUE 250 CONTINUE C C Update the diameters read from the catalog with the C more current compilation in DIAMETER.DAT C C Note that the uniform disk angular diameter is a function of C wavelength. The values in DIAMETER.DAT are for 800 nm. C QCAL scales the diameter from 800 nm to the appropriate wavelength. C OPEN(UNIT=59,FILE='X:\\DIAMETER.DAT',STATUS='OLD', IOSTAT=IERR ) IF (IERR .NE. 0 ) THEN WRITE(6,*) ' DIAMETER FILE DOES NOT EXIST ' GO TO 400 END IF 325 CONTINUE READ (59,1300,END=360) I, NEWDIAM 1300 FORMAT ( 1X, I5, F7.0 ) DO J = 1, NSTAR IF ( I .EQ. SLIST(J) ) THEN K = J GO TO 335 END IF END DO GO TO 325 335 CONTINUE WRITE(6,*) ' STAR = ', I, ' DIAMETER = ', NEWDIAM ANGDIA(K) = NEWDIAM C WRITE(6,1700) WAVE DO 350 I = 1, ISCAN IF ( STAR(I) .NE. SLIST(K) ) GO TO 350 DO J = 1, 4 CALL QCAL ( BASE(I), HA(I), DEC(I), WAVE(J), $ ANGDIA(K), VEST(I,J) ) END DO 350 CONTINUE GO TO 325 360 CONTINUE CLOSE(UNIT=59) 400 CONTINUE C C Update the magnitudes and colors from the FK5 catalog with the C more current compilation in STARFLUX.DAT C OPEN(UNIT=59,FILE='X:\\STARFLUX.DAT',STATUS='OLD', IOSTAT=IERR ) IF (IERR .NE. 0 ) THEN WRITE(6,*) ' STARFLUX FILE DOES NOT EXIST ' GO TO 500 END IF 425 CONTINUE READ (59,*,END=460) I, NEWMAG, NEWBV C WRITE(6,*) ' STAR = ', I, ' MAGNITUDE = ', NEWMAG, C $ ' COLOR = ', NEWBV DO 430 J = 1, NSTAR IF ( I .EQ. SLIST(J) ) THEN MAG(J) = NEWMAG BV(J) = NEWBV END IF 430 CONTINUE GO TO 425 460 CONTINUE CLOSE(UNIT=59) 500 CONTINUE WRITE(6,2020) NSTAR WRITE(6,2021) NBASE WRITE(6,2025) ( BLIST(J),J=1,NBASE ) 950 CONTINUE RETURN 1700 FORMAT ( 10X, ' wavelengths = ', 4F8.0, $ / ' STAR B HA DEC DIA CH 1 CH 2 CH 3 CH 4') 1705 FORMAT ( I5, I3, 2F6.2, F6.1, 4F8.3 ) 2020 FORMAT( I7, ' OBSERVED STARS' ) 2021 FORMAT( ' THE', I3, ' OBSERVED BASELINES ARE' ) 2025 FORMAT( 5I10 ) 2050 FORMAT ( 2I5, 2X, 2A8 ) END