SUBROUTINE CGETSTAR C=================================================================== C This subroutine is used by BINFAKE, the program for generating fake C Mt. Wilson data. It gets star number, observation date, baseline, C filters, widths, and uncertainties; tells the user when the star is C up; and gets the start and stop times and interval between C observations. C C J.T. Armstrong 16 Jan 1991 C=================================================================== INCLUDE 'BINFIT.INC' INTEGER*4 NPARS PARAMETER (NPARS=5+4*MXFILT) INTEGER*4 IDAY, IMON, IYEAR, NBODY, EVERUP, IERR, SUN, EARTH INTEGER*4 MODE REAL*8 ZLIMSUN, ZLIMSTAR, ZANGLE REAL*8 UTSTDOWN, UTSTUP, UTSUNDOWN, UTSUNUP, UTSTART, UTSTOP REAL*8 RAM, DECM, PMRA, PMDEC, PAR, RV, RAA, DECA, RA, DEC, HA REAL*8 RASUN, DECSUN, DIST, INTERVAL REAL*8 DJUL, TJD, JD, UT, UT1, UTINCR, GAST0, GAST, DOT, TDTOFF REAL*8 COSZL, COSDS, SINDS, COSL, SINL REAL*8 B(3), U(3), V(3), W(3), BX(MXBAS), BY(MXBAS), BZ(MXBAS) REAL*8 PARS(NPARS), VALS(NPARS), ENDMRK CHARACTER*64 BFILENAME PARAMETER (ZLIMSUN=102.) ! Sunrise/set limit PARAMETER (TDTOFF=55.D0/86400.D0) ! TDT offset for late 1990 EXTERNAL JD, DOT DATA PARS/'DATE','ZENlimit','INTerval','FKno','BLno', 1 MXFILT*'FILTers',MXFILT*'WIDths', 2 MXFILT*'FRACerr',MXFILT*'ABSerr' / DATA VALS/ 0 , 50D0 , 10.0 , 0 , 0 , 1 MXFILT*0D0, MXFILT*0D0, 2 MXFILT*0D0, MXFILT*0D0 / DATA ENDMRK/'/'/ C Get star number, observation date, baseline, filters, widths, and C uncertainties: MODE = 2 100 WRITE(OUTC,'(A,A,A,A)') ' Enter date (YYMMDD), Zen. angle ', 1 'limit, obs. interval (min), FK#, baseline #, ' WRITE(OUTC,'(A,A,A,A)') ' filters (nm), widths (nm), fractional ', 2 'and absolute errors (%): ' CALL KEYIN(PARS,VALS,NPARS,ENDMRK,MODE,INC,OUTC) DATE = INT(VALS(1)) ZLIMSTAR = VALS(2) INTERVAL = VALS(3)/60. STARNO = INT(VALS(4)) BASELINE(1) = INT(VALS(5)) NFILT = 0 DO IG = 1, MXFILT IF (VALS(5+IG).EQ.0.) THEN GO TO 110 ELSE NFILT = NFILT + 1 LAMBDA0(IG) = VALS(5+IG) DLAMBDA(IG) = VALS(5+MXFILT+IG) FRACERR(IG) = VALS(5+2*MXFILT+IG)/100. ABSERR (IG) = VALS(5+3*MXFILT+IG)/100. LAM0 (IG) = LAMBDA0(IG) * 1.0D-9 ! Convert to meters DLAM (IG) = DLAMBDA(IG) * 1.0D-9 NU0 (IG) = C/LAM0(IG) NULO (IG) = C / ( LAM0(IG) + DLAM(IG)/2.0D0 ) NUHI (IG) = C / ( LAM0(IG) - DLAM(IG)/2.0D0 ) DNU (IG) = ( NUHI(IG) - NULO(IG) ) / NUSTEPS END IF END DO 110 IF (NFILT.EQ.0) GO TO 100 IDAY = MOD(DATE, 100) IMON = MOD(DATE/100, 100) IYEAR = MOD(DATE/10000, 100) + 1900 DJUL = JD(IDAY,IMON,IYEAR) ! Julian date for 0 hrs UT D WRITE(OUTC,'(A,F10.1)') ' DJUL = ',DJUL D WRITE(OUTC,'(1X,I3,A,3(1X,F6.1))') NFILT, D 1 ' filters used; wavelengths ',(LAMBDA0(IG),IG=1,3) C Figure out when night starts and ends (to within a couple minutes) SUN = 0 ! Code for sun for APPLAN EARTH = 3 ! Code for earth for APPLAN CALL APPLAN (DJUL, SUN, EARTH, RASUN, DECSUN, DIST) ! Sun pos. at 0 hrs UT D WRITE(OUTC,'(A,F10.4,1X,F10.4)') ' RA, Dec of Sun: ',RASUN,DECSUN CALL SIDTIM (DJUL, 0., 1, GAST0) ! Sid time at 0 hrs UT CALL RISESET (ZLIMSUN, RASUN, DECSUN, LONG, LAT, GAST0, 1 UTSUNDOWN, UTSUNUP, EVERUP) C Get star number; find out when it's up CALL READFK5 (STARNO, STARNAME, RAM, DECM, PMRA, PMDEC, PAR, RV) NBODY=3 CALL APSTAR (DJUL, NBODY, RAM, DECM, PMRA, PMDEC, PAR, RV, 1 RAA, DECA) CALL RISESET (ZLIMSTAR, RAA, DECA, LONG, LAT, GAST0, 1 UTSTDOWN, UTSTUP, EVERUP) UTSTUP = UTSTUP - 24. WRITE(OUTC,'(A,F10.4,1X,F10.4)') ' Star rise, set times: ', 1 UTSTUP,UTSTDOWN C Restrict calculations to time between sunset and starset; if there C is no such time, try another star IF(EVERUP.EQ.1) THEN UTSTART = UTSUNDOWN ! Circumpolar UTSTOP = UTSUNUP ELSE IF(EVERUP.EQ.0) THEN IF(UTSTUP.LT.UTSUNDOWN .AND. UTSTUP.GT.(UTSUNUP-24.)) THEN UTSTART = UTSUNDOWN ! Star rises before sunset IF(UTSTDOWN.LT.UTSUNDOWN) THEN UTSTOP = UTSUNDOWN ! Star sets before sunset, so EVERUP = -1 ! not observable tonight ELSE IF(UTSTDOWN.LT.UTSUNUP) THEN UTSTOP = UTSTDOWN ! Star sets during night ELSE UTSTOP = UTSUNUP ! Star up all night END IF END IF ELSE IF(UTSTUP.GE.UTSUNDOWN .AND. UTSTUP.LT.UTSUNUP) THEN UTSTART = UTSTUP ! Star rises during night IF(UTSTDOWN.LT.UTSUNUP) THEN UTSTOP = UTSTDOWN ! Star sets during night ELSE UTSTOP = UTSUNUP ! Star up rest of night END IF ELSE IF(UTSTUP.GE.UTSUNUP) THEN UTSTART = UTSUNDOWN+24. ! Star rises after sunrise IF(UTSTDOWN.LT.(UTSUNDOWN+24.)) THEN UTSTOP = UTSUNDOWN ! Star sets before sunset EVERUP = -1 ELSE UTSTOP = UTSTDOWN ! Star sets after sunset END IF ELSE WRITE(OUTC,*) ' Not all cases were covered.' END IF END IF IF(EVERUP.EQ.-1) THEN WRITE(OUTC,'(A,A)') ' This star can''t be observed tonight.', 1 ' Try again.' GO TO 100 END IF WRITE(OUTC,'(A,F7.4,1X,F7.4)') ' Start, stop times: ', 1 UTSTART,UTSTOP C Calculate the observation times and use same baseline for all points NDATA = INT((UTSTOP-UTSTART)/INTERVAL) DO ID = 1, NDATA HOURS(ID) = UTSTART + (ID-1)*INTERVAL BASELINE(ID) = BASELINE(1) END DO D WRITE(OUTC,'(A,I4)') ' Number of scans to calculate: ',NDATA C Read in the baseline file. OPEN(UNIT=INMOD,FILE='BASELINE.DAT',STATUS='OLD',IOSTAT=IERR) IF ( IERR .NE. 0 ) THEN WRITE (OUTC,'(A)') 1 ' File BASELINE.DAT not in current directory.' WRITE (OUTC,'(A,$)') ' Provide name of baseline file: ' READ (INC,'(A)',IOSTAT=IERR) BFILENAME OPEN (UNIT=INMOD,FILE=BFILENAME,STATUS='OLD',IOSTAT=IERR) IF (IERR .NE. 0) THEN WRITE(OUTC,'(/,A)') ' Baseline file not found.' STOP END IF END IF READ(INMOD,'(1(/))') ! skip header lines in file DO I = 1, MXBAS READ(INMOD,'(I5,3(F15.6))',IOSTAT=IERR) 1 J, BX(I), BY(I), BZ(I) IF ( (IERR .NE. 0) .OR. (J .NE. I) ) THEN WRITE(OUTC,'(A,$)') ' Baseline file is corrupt.' STOP END IF END DO CLOSE(INMOD) D WRITE(OUTC,'(A)') ' Baseline file successfully read.' C Calculate U, V, W for given baseline. Use limit on W to set C restrictions on observing time. Try another baseline if the star C can't be observed on this one. C================================================================= C Calculate the (u,v) position for each observation. C DO J = 1, NDATA IF ( J.EQ.1 ) THEN IDAY = MOD(DATE, 100) IMON = MOD(DATE/100, 100) IYEAR = MOD(DATE/10000, 100) + 1900 D WRITE(OUTC,1300) J, DATE, IDAY, IMON, IYEAR 1300 FORMAT( ' SCAN = ', I5, ' DATE = ', I10, 3I5,/ ) DJUL = JD ( IDAY, IMON, IYEAR ) B(1) = BX(BASELINE(J)) ! Baseline components in meters B(2) = BY(BASELINE(J)) B(3) = BZ(BASELINE(J)) D WRITE(OUTC,'(A,3(1X,G12.5))') ' B(1,2,3) (m): ', D 1 B(1),B(2),B(3) END IF UT1 = HOURS(J) / 24.D0 TJD(J) = DJUL + UT1 + TDTOFF CALL SIDTIM ( DJUL, UT1, 1, GAST ) CALL DIURN ( TJD, GAST, RAA, DECA, RA, DEC ) HA = 15.D0 *( GAST - RA ) - LONG IF ( HA .LT. -180 ) HA = HA + 360.D0 C Components of u, v, w axes in local coordinate system U(1) = -SIN(HA*DEGRAD) U(2) = -COS(HA*DEGRAD) U(3) = 0.D0 D WRITE(OUTC,'(A,3(1X,G12.5))') ' U(1,2,3) : ',U(1),U(2),U(3) V(1) = SIN(DEC*DEGRAD) * COS(HA*DEGRAD) V(2) = -SIN(DEC*DEGRAD) * SIN(HA*DEGRAD) V(3) = -COS(DEC*DEGRAD) D WRITE(OUTC,'(A,3(1X,G12.5))') ' V(1,2,3) : ',V(1),V(2),V(3) W(1) = -COS(DEC*DEGRAD) * COS(HA*DEGRAD) W(2) = COS(DEC*DEGRAD) * SIN(HA*DEGRAD) W(3) = -SIN(DEC*DEGRAD) D WRITE(OUTC,'(A,3(1X,G12.5))') ' W(1,2,3) : ',W(1),W(2),W(3) DO IG = 1, NFILT UDATA(J,IG) = DOT( B, U ) / (LAMBDA0(IG)*1.0E-9) VDATA(J,IG) = DOT( B, V ) / (LAMBDA0(IG)*1.0E-9) DELAY(J,IG) = DOT( B, W ) / (LAMBDA0(IG)*1.0E-9) D WRITE(OUTC,'(A,I4,1X,I2,3(1X,G12.5))') D 1 ' Scan, Filter, U, V, W (wavel): ', J, IG, D 2 UDATA(J,IG),VDATA(J,IG),DELAY(J,IG) END DO END DO D WRITE(OUTC,*) ' U, V coordinates calculated' RETURN END