FTN77
$FILES 0,1
      Program SUN7 (1,10)
C     ===========
C
C     Make a map of the sun using the PC-Steer control system.
C
C     MODIFICATIONS :  1992-04-29 JQ
C                    * Now assumes a certain PC-Steer Feed System
C                      arrangement in SysName, NDcode, DScode, TND
C                      Tilted etc. and performs automatic setup
C                      accordingly.
C                    * Checks consistency of Feed System definition
C                      against PC-Steer and also compares clocks.
C                    * Will run to System console by default but can be
C                      set to the Esprit by typing XQ,SUN6,7
C
C                      1992-06-25 JQ
C                    * Now suggests more intelligent defaults for scan
C                      spacing and rate for each Feed system.
C                    * The map need no longer be square and the centre
C                      may also be offset in both Ecliptic Long and Lat.
C                    * PC-Steer Errors now cause the program to suspend
C                      and ask for user intervention.
C
      Implicit  None
C
C     remember to update DATA statement below if MaxSys is changed
C
      Integer   MaxSys, MaxBin
C
      Parameter ( MaxSys=10 , MaxBin=101 )
C
      Logical   Tilted(MaxSys), FilExist, FAULT
C
      Character SysName(MaxSys)*12, Reply*166, JUNK*1, FILNAM*10,
     *          SysonPC*12
C
      Integer*2 IPRAM(5), ITIME(5), IYEAR, UThrs, UTmin, UTsec
C
      Integer   Ad, Am, As, UNcmd(2), Lbuf(500), AntennaLU, FLukeLU,
     *          PlotLU, FeedSys, Nscan, iscan, Count(MaxBin), Nsys,
     *          NDcode(MaxSys), DScode(MaxSys), I, Nfile, iSrate,
     *          ISWIT, SWIN, Istat, Nbins
C
      Real      Atoday, Atomorrow, ScanRate, Width, A0sun, A1sun,
     *          DeltaT, Tscale, Zscale, Out(MaxBin), DVMmax, Base,
     *          FOCUS, TILT, Rdvm, TND(MaxSys), FracDay, SizeLong,
     *          SizeLat, A0off, B0off, SpacDef(MaxSys)
C
      Double Precision MJDPC, MJDHP
C
      Data      AntennaLU /39/, FlukeLU /35/, PlotLU /31/
C
C            Feed System names : must agree with those in PC-STEER
C
      Data      SysName / '2.5cm single',      '3.6cm dicke ',
     *                    '3.6cm single',      '4.5cm single',
     *                    '6cm dicke   ',      '6cm single  ',
     *                    '13cm single ',      '18cm single ',
     *                    '18cm untilt ',      '            ' /
C
C            Noise Diode select codes to be output to DIO channel 3
C
      Data      NDcode  / 2, 1, 1, 2, 3, 2, 4, 5, 5, 0 /
C
C            DVM select codes to be output to DIO channal 5
C
      Data      DScode  / 7, 1, 1, 7, 3, 3, 5, 7, 7, 0 /
C
C            Equivalent Temperature of Noise Diodes in Kelvin
C
      Data      TND     / 158.0,  10.0,  10.0,  82.0,  1.54,
     *                     95.0,  3.96,  33.0,  33.0,   0.0 /
C
C            Required Sub-Reflector position i.e. Tilted ? Y/N
C
      Data      Tilted  / .FALSE., .FALSE., .FALSE., .FALSE.,
     *                    .FALSE., .FALSE., .FALSE., .FALSE.,
     *                     .TRUE., .FALSE. /
C
C            Suitable scan spacing ( approximately FWHM/2 of beam )
C
      Data      SpacDef / 0.030, 0.050, 0.050, 0.060, 0.080,
     *                    0.080, 0.160, 0.250, 0.250, 0.000 /
C
C            Total Number of usable Feed Systems
C
      Data      Nsys    / 9 /
C
C            Default output filename ( 000 replaced by Day Number )
C
      Data      FILNAM / 'SU000A::15' /
C
      Call LGBUF (Lbuf,250)
C
C     recall system console LU from parameter list. Defaults to 1
C
      Call RMPAR (IPRAM)
      Call FSYSU (IPRAM(1),IPRAM(1))
      Print '(a)',CHAR(12)
C
C     check PC vs HP clock and obtain Day Number and Year
C
      Call EXEC (11,ITIME,IYEAR)
      Read (AntennaLU,'(27x,f13.6)') MJDPC
      Call UNLT
      Call MJDAY (IYEAR,ITIME,MJDHP)
      IF ( DABS(MJDPC-MJDHP) .gt. 0.58D-5 ) then
C
C        clocks differ by more than half a second - warn user
C
         Print *,' '
         Print *,'!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!'
         Print *,'!!! PC-STEER & HP1000 Clocks do not agree !!!'
         Print *,'!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!'
         Print *,' '
      endif
C
C     set up default file name using Day Number
C
      Write (FILNAM(3:5),'(I3.3)') ITIME(5)
C
      Print '(16X,"Sun Mapping Program Day ",I3.3,", ",I4)',
     *    ITIME(5), IYEAR
      Print *, '               ================================='
      Print *, ' '
      Print *, 'The HP9872S Plotter MUST be connected to the HPIB'
      Print *, '    of the CONTROL computer during this program.'
      Print *, ' '
      Print *, 'Ephemeris Ref : "The Astronomical Almanac" Section C '
      Print *, ' '
C
      Print *, 'Enter Eclip. Long. of sun for 0h DT today (dd mm ss): _'
      Read *, Ad, Am, As
      Atoday = Ad + Am / 60.0 + As / 3600.0
C
      Print *, ' '
      Print *, 'Enter Eclip. Long. of sun for 0h DT tomorrow: _'
      Read *, Ad, Am, As
      Atomorrow = Ad + Am / 60.0 + As / 3600.0
C
      Print *, ' '
      Print *, 'Enter the difference TAI-UTC from page K9: _'
      Read *, DeltaT
      DeltaT = DeltaT + 32.184
C
10101 FeedSys = 1
      Print *, ' '
      Print *, 'PC-Steer Feed systems as implemented here :'
      Print *, ' '
      DO I=1,Nsys/2*2,2
          Print *,I,' : ',SysName(I),'     ',I+1,' : ',SysName(I+1)
      END DO
      IF (Nsys/2*2 .lt. Nsys)
     *    Print *,Nsys,' : ',SysName(Nsys)
      Print *, ' '
      Print *, 'Enter Feed System (1-8 / =',FeedSys,'): _'
      Read *, FeedSys
      IF (FeedSys .le. 0 .or. FeedSys .gt. Nsys) go to 10101
C
C     check subreflector position
C
      FAULT = .FALSE.
      CALL RHYP (FOCUS, TILT)
      IF (.not. Tilted(FeedSys)) then
         IF ( FOCUS .lt. 7.4 ) then
            PRINT *,' '
            PRINT *,'Hyperbola focus = ',FOCUS,' not 7.5'
            FAULT = .TRUE.
         endif
         IF ( TILT .gt. -0.5 ) then
            PRINT *,' '
            PRINT *,'Hyperbola tilt = ',TILT,' not -0.54'
           FAULT = .TRUE.
         endif
      else
         IF ( FOCUS .gt. 5.1 ) then
            PRINT *,' '
            PRINT *,'Hyperbola focus = ',FOCUS,' not 5.0'
            FAULT = .TRUE.
         endif
         IF ( TILT .lt.  9.7 ) then
            PRINT *,' '
            PRINT *,'Hyperbola tilt = ',TILT,' not 9.80'
            FAULT = .TRUE.
         endif
      endif
      IF (FAULT) go to 10101
C
      SizeLong = 1.0
      Print *, ' '
      Print '(1X,"Enter Map Size in Eclip. Long (deg / = ",F4.2,
     *           "): _")', SizeLong
      Read *, SizeLong
C
      SizeLat = SizeLong
      Print *, ' '
      Print '(1X,"Enter Map Size in Eclip. Lat. (deg / = ",F4.2,
     *           "): _")', SizeLat
      Read *, SizeLat
C
      A0off = 0.0
      Print *, ' '
      Print '(1X,"Enter Map Centre Long offset  (deg / = ",F4.2,
     *           "): _")', A0off
      Read *, A0off
C
      B0off = 0.0
      Print *, ' '
      Print '(1X,"Enter Map Centre Lat. offset  (deg / = ",F4.2,
     *           "): _")', B0off
      Read *, B0off
C
      Width = SpacDef(FeedSys)
      Print *, ' '
      Print '(1X,"Enter Scan Spacing (deg / = ",F4.3,"): _")',Width
      Read *, Width
C
      ScanRate = Width * 1.2
      Print *, ' '
      Print '(1X,"Enter Scan Rate (deg/sec / = ",F4.3,
     *           "): _")',ScanRate
      Read *, ScanRate
C
      Zscale = 7.5
      Print *, ' '
      Print *, 'Enter central height (cm / =',Zscale,'): _'
      Read *, Zscale
C
C     convert to Plot-Units ( normalised to sun after calibration )
C
      Zscale = Zscale * 400.0
C
C     set Fluke to required sample rate
C
      iSrate = 1
      Print *, ' '
      Print *,'Enter DVM sampling rate ( 0=S, 1=M, 2=F, / =',iSrate,
     *        '): _'
      Read *, iSrate
      Write (FlukeLU,'(a,i1)') '*S',iSrate
C
C     select DVM inputs according to FeedSys
C
      Call EXEC (2, 113B, DScode(FeedSys)*11B, 1, 5, 5)
      Call EXEC (12,0,2,0,-1)
C
C     kill the TEMPS program if it is running
C
      Call MESSS (10HOF,TEMPS,1,10)
C
      Nscan = (Nint(SizeLong / Width) / 2) * 2 + 1
      SizeLong = Width * (Nscan - 1)
      Nbins = (Nint(SizeLat / Width) / 2) * 2 + 1
      if (Nbins .gt. MaxBin) then
         Print *,' '
         Print *,'!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!'
         Print *,'!!! Nbins >',MaxBin,' : Latitude size corrected !!'
         Print *,'!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!'
         Print *,' '
         Nbins=MaxBin
      endif
      SizeLat = Width * (Nbins - 1)
C
C     calculate co-ordinates of Sun
C
      Read (AntennaLU,'(27x,f13.6)') MJDPC
      Call UNLT
      FracDay = Dmod(MJDPC,1.0D0) + DeltaT/86400.0
      A1sun = (Atomorrow - Atoday) / 86400.0
      A0sun = Atoday + FracDay * A1sun * 86400.0
C
C     drive to calculated position of Sun ( Note : no Centre offset )
C
      Write (AntennaLU,'(a)') 'Begin FindSun'
      Write (AntennaLU,'(a)') 'Clear_Fifo'
      Write (AntennaLU,'(a)') 'Duration 900'
      Write (AntennaLU,'(a)') 'Coordinate_System ECLIPTIC'
      Write (AntennaLU,'(a)') 'Pointing On'
      Write (AntennaLU,'(a,i2)') 'Feed_System ', FeedSys
      Write (AntennaLU,'(a)') 'A0 0.0'
      Write (AntennaLU,'(a)') 'A1 0.0'
      Write (AntennaLU,'(a)') 'B0 0.0'
      Write (AntennaLU,'(a)') 'B1 0.0'
      Write (AntennaLU,'(a,f7.3)') 'Psi ', A0sun
      Write (AntennaLU,'(a)') 'Phi 0.0'
      Write (AntennaLU,'(a)') 'Theta 0.0'
      Write (AntennaLU,'(a)') 'Parallax 0.0'
      Write (AntennaLU,'(a)') 'End'
      Call UNLT
C
      Print *, ' '
      Print *, 'Driving to Sun .....'
C
C     check PC Selected Feed System Name against internal name
C
      Read (AntennaLU,'(40x,a12)') SysonPC
      Call UNLT
      IF ( SysonPC .ne. SysName(FeedSys) ) then
C
C        PC-STEER Feed name differs from internal - warn user
C
         Print *,' '
         Print *,'!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!'
         Print *,'!!! PC-STEER Feed System list inconsistent !!!'
         Print *,'!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!'
         Print *,' '
      endif
C
C        dither till PC-STEER is 'On Source' or an error occurs
C
1        Call EXEC (12,0,2,0,-3)
         Read (AntennaLU,'(a)') Reply
         Call UNLT
         if (Reply(1:1) .eq. '1') goto 1
C
      if (Reply(1:1) .ne. '0') then
         read (Reply,'(I1)') Istat
         call ERROR (Istat,FAULT)
         if (FAULT) goto 999
         goto 1
      endif
C
      Print *, ' '
      Print *, 'Tracking Sun, Adjust radiometer attenuation'
      Print *, 'Press RETURN when ready: _'
      Read (*,'(a)') junk
C
C     ensure that all noise diodes are turned off
C
      Call EXEC(2,113B,0,1,3,3)
      Call EXEC(12,0,2,0,-2)
C
C     sample signal peak on the sun for scaling purposes
C
      Read (FlukeLU,*) DVMmax
C
      Do iscan = -(Nscan/2) , +(Nscan/2)
C
C       recalculate co-ordinates of Sun for present time
C
        Read (AntennaLU,'(27x,f13.6)') MJDPC
        Call UNLT
        FracDay = Dmod(MJDPC,1.0D0) + DeltaT/86400.0
        UThrs = INT(  Dmod(MJDPC,1.0D0)*24.0  )
        UTmin = INT( (Dmod(MJDPC,1.0D0)*24.0-UThrs)*60.0 )
        UTsec = INT(((Dmod(MJDPC,1.0D0)*24.0-UThrs)*60.0-UTmin)*60.0)
        A0sun = Atoday + FracDay * A1sun * 86400.0
C
C       drive to start of next scan ( Map centre offset from Sun )
C
        Write (AntennaLU,'(a,i2.2)') 'Begin Start', iscan+Nscan/2+1
        Write (AntennaLU,'(a)') 'Duration 400'
        Write (AntennaLU,'(a,f7.3)') 'Psi ', A0sun + A0off
        Write (AntennaLU,'(a,f7.3)') 'Phi ',       - B0off
        Write (AntennaLU,'(a,f7.3)') 'A0 ', Width*iscan
        Write (AntennaLU,'(a,f7.3)') 'B0 ', -SizeLat / 2.0
        Write (antennaLU,'(a)') 'B1 0.0'
        Write (AntennaLU,'(a)') 'End'
        Call UNLT
        Print *,' '
        Print *,'Driving to scan',iscan+Nscan/2+1,'_'
        Print '(3X,I2,2(":",I2.2)," ... _" )', UThrs, UTmin, UTsec
C
        If (iscan .ne. (-Nscan/2)) then
C
C          plot and write to disc the previous scan
C
           CALL OutScan (iscan-1, Nscan, Nbins, Tscale, Zscale,
     *                   Out, Count)
        else
C
C          set up the plotter
C
           Write (PlotLU,'("in;")')
           Write (PlotLU,'("sp 1;pu;pa 4000,10000")')
           Write (PlotLU,'("si 1.5,1.5;lb The SUN",a)') Char(3)
           Write (PlotLU,'(";pu;pa 400,2500")')
           Write (PlotLU,'("si 0.2,0.3;lb Using ",4a)')
     *           SysName(FeedSys), Char(10), Char(13), Char(3)
           Write (PlotLU,'("lb Map Size ",F3.1," by ",F3.1," deg",3a)')
     *           SizeLong, SizeLat, Char(10), Char(13), Char(3)
           Write (PlotLU,'("lb Scan Spacing ",F5.3," deg",3a)')
     *           Width, Char(10), Char(13), Char(3)
           Write (PlotLU,'("lb Scan Rate ",F5.3," deg/s",3a)')
     *           ScanRate, Char(10), Char(13), Char(3)
           IF (UThrs .lt. 10) then
            Write (PlotLU,'("lb ",I1,2(":",I2.2)," Day ",I3.3,I5,3a)')
     *            UThrs, UTmin, UTsec, ITIME(5), IYEAR, Char(10),
     *            Char(13), Char(3)
           else
            Write (PlotLU,'("lb ",I2,2(":",I2.2)," Day ",I3.3,I5,3a)')
     *            UThrs, UTmin, UTsec, ITIME(5), IYEAR, Char(10),
     *            Char(13), Char(3)
           endif
C
C          set up the standard output disc file
C
           FilExist=.TRUE.
           Nfile=0
           DO WHILE (FilExist)
               FILNAM(6:6)=CHAR(ICHAR('A')+Nfile)
               Inquire (FILE=FILNAM, EXIST=FilExist)
               Nfile = Nfile+1
           END DO
           Print *,' '
           Print *,' '
           Print *,'Creating results file ',FILNAM,' ...'
           Open (FILE=FILNAM, UNIT=20, STATUS='NEW')
           Write (PlotLU,'("si 0.2,0.3;lb Output File ",4a)')
     *           FILNAM, Char(10), Char(13), Char(3)
           Write (PlotLU,'("sp 0;")')
        endif
C
C       wait here till drive to start of scan is complete
C
2          Call EXEC (12,0,2,0,-2)
           Read (AntennaLU,'(a)') Reply
           Call UNLT
           If (Reply(1:1) .eq. '1') goto 2
C
        if (Reply(1:1) .ne. '0') then
           read (Reply,'(I1)') Istat
           call ERROR (Istat,FAULT)
           if (FAULT) goto 999
           goto 2
        endif
C
C
        If (iscan .eq. (-Nscan/2)) then
C
C          calibrate on the first scan
C
           Print *, ' '
           Print *, 'Ready to Calibrate : Remove 20dB of attenuation'
           Print *, 'Press RETURN when ready: _'
           Read (*,'(a)') junk
C
           Print *,' '
           CALL CALIBRATE (TND(FeedSys), NDcode(FeedSys), Tscale, FAULT)
           Print *, ' '
           IF (FAULT) goto 999
C
           Print *, 'Done ! Put back 20dB attenuation'
           Print *, 'Press RETURN when ready: _'
           Read (*,'(a)') junk
           Print *, ' '
C
C          sample signal at map corner for base level
C
           Base = 0.0
           DO I=1,20
               Read (FlukeLU,*) Rdvm
           Base = Base + Rdvm
           END DO
           Base = Base / 20.0
C
C          Normalise central height by signal on sun
C
           Zscale = Zscale / (DVMmax - Base)
C
        endif
C
C       reselect Fluke DVM inputs
C
        Call EXEC (2, 113B, DScode(FeedSys)*11B, 1, 5, 5)
        Call EXEC (12,0,2,0,-1)
C
C       duration extended by 5 Seconds to ensure Scan is completed
C
        Write (AntennaLU,'(a,i2.2)') 'Begin SunScn', iscan+Nscan/2+1
        Write (AntennaLU,'(a,f8.3)') 'Duration ', SizeLat/ScanRate+5
        Write (AntennaLU,'(a,f7.3)') 'B1 ', ScanRate
        Write (AntennaLU,'(a)') 'End'
        Call UNLT
        Print *,'Doing scan', iscan+Nscan/2+1,' of',Nscan,' _'
        Call DoScan (Nbins, SizeLat, Base, Out, Count, FAULT)
        if (FAULT) goto 999
      End Do
C
      call OutScan (iscan-1, Nscan, Nbins, Tscale, Zscale, Out, Count)
C
      Close (20,err=900)
C
  900 Print *, ' '
      Print *, '|======================|'
      Print *, '|  SUN scans finished  |'
      Print *, '|======================|'
C
C     Chart Advance on the Plotter
C
      Write (PlotLU,'("af;")')
C
C     Send telescope to Zenith if bit 11 is set
C
      if (ISSW(11) .lt. 0) then
         ISWIT = SWIN(0)
         ISWIT = ISWIT .AND. 173777B
         CALL SWOUT(ISWIT)
         Write (AntennaLU,'(a)') 'Begin ZENITH'
         Write (AntennaLU,'(a)') 'Clear_Fifo'
         Write (AntennaLU,'(a)') 'Duration 9999'
         Write (AntennaLU,'(a)') 'Coordinate_System HORIZON'
         Write (AntennaLU,'(a)') 'Pointing Off'
         Write (AntennaLU,'(a)') 'Feed_System 0'
         Write (AntennaLU,'(a)') 'A0 0.0'
         Write (AntennaLU,'(a)') 'A1 0.0'
         Write (AntennaLU,'(a)') 'B0 90.0'
         Write (AntennaLU,'(a)') 'B1 0.0'
         Write (AntennaLU,'(a)') 'Psi 0.0'
         Write (AntennaLU,'(a)') 'Phi 0.0'
         Write (AntennaLU,'(a)') 'Theta 0.0'
         Write (AntennaLU,'(a)') 'Parallax 0.0'
         Write (AntennaLU,'(a)') 'End'
         Call UNLT
      endif
C
  999 Print *,'Bye !'
      Stop
      End
C
C
      SUBROUTINE CALIBRATE (TND, NDcode, Tscale, FAULT)
C     ====================
C     To obtain the Kelvins / Volt calibration factor
C     Modified to take out linear slope in sampling
C
C     inputs :    TND     real        Cal Noise Diode in Kelvins
C                 NDcode  integer     Cal Noise Diode select code in DIO
C
C     outputs:    Tscale  real        Kelvins / Volt calibration factor
C                 FAULT   logical     .TRUE. if unable to calibrate
C
      Implicit None
C
      Logical          FAULT
C
      Integer          I, J, NVALID(4), NDcode, FlukeLU, Retries, Number
C
      Real             TND, Tscale, DTScale
C
      Double Precision DVM(4), DVMSQD(4), ADVM(4), DMS(4), RMS(4),
     *                 RDVM, RMS_SUM, DMS_SUM, cal_COUNT
C
      Data             FlukeLU / 35 /
C
      Retries = 0
      Number = 150
      FAULT = .FALSE.
      Print *,' System calibration :'
      Print *,' '
C
C         restart point after error  - allow for up to 10 failures
C
   10     Retries = Retries + 1
          DMS_SUM = 0D0
          RMS_SUM = 0D0
C
C         sample with noise cal off (J = 1,4) or on (J = 2,3)
C
          DO J = 1, 4
              DVM(J)  = 0D0
              DVMSQD(J) = 0D0
              NVALID(J)  = 0
              IF (J .EQ. 2) then
C
C                 turn on the noise diode
C
                  Call EXEC (2,113B,NDcode,1,3,3)
              end if
              IF (J .EQ. 1 .OR. J .EQ. 4) then
C
C                 turn off the noise diode
C
                  Call EXEC (2,113B,0,1,3,3)
              end if
              IF (J .NE. 3) then
C
C                 let the system settle
C
                  Call EXEC (12,0,2,0,-2)
              end if
C
C             sample 0.1 sec filtered 'direct' radiometer output
C
              DO I = 1,Number
                  read (FlukeLU,*) RDVM
C
C                 Store DVM and DVM SQUARED readings
C
                  NVALID(J) = NVALID(J) + 1
                  DVM(J) = DVM(J) + RDVM
                  DVMSQD(J) = DVMSQD(J) + RDVM * RDVM
              end do
C
C             sampling done : calculate average DVM reading and rms error
C
              IF (NVALID(J) .gt. 2) then
C
C                 average DVM reading
C
                  ADVM(J) = DVM(J) / NVALID(J)
                  DMS(J) = (DVMSQD(J) - DVM(J)*DVM(J)/NVALID(J))
     *                      / (NVALID(J) - 1.0)
                  DMS_SUM = DMS_SUM + DMS(J)
C
C                 rms error in reading
C
                  RMS(J) = DSQRT(DMS(J))
                  RMS_SUM = RMS_SUM + RMS(J)
C
                  IF (J .EQ. 1 .OR. J .EQ. 4) then
                      Print *,'   Noise cal. off : _'
                  endif
                  IF (J .EQ. 2 .OR. J .EQ. 3) then
                      Print *,'   Noise cal.  on : _'
                  endif
                  Print '(F12.6,'' +-'',F12.6,'' V'')',ADVM(J),RMS(J)
              endif
          end do
C
C         error trap
C
          IF (Retries .ge. 10) then
              FAULT = .TRUE.
              Print *,' '
              Print *,' Abandoning calibration'
              RETURN
          endif
C
C         check noise tube > 4 * average rms noise
C
          cal_COUNT = ((ADVM(2) + ADVM(3)) - (ADVM(1) + ADVM(4)))/2
          if (DABS(cal_COUNT) .lt. RMS_SUM) then
              Print *,' '
              Print *,' Noise cal not seen - trying again'
              Print *,' '
              go to 10
          endif
C
C     calculate Kelvins per DVM unit conversion factor
C
      Tscale = TND / cal_COUNT
C
C     fractional error in Tscale is half of the square root of the sum of
C     the squares of the fractional errors in the four DVM averages
C
      DTscale = Tscale * DSQRT (DMS_SUM / (2*cal_COUNT)**2)
      DTscale = ABS (DTscale)
      Print *,' '
      Print '(1x,"Calibration gives",F12.6," +-",F8.6,
     *           " x 100 Kelvins / Volt")', Tscale, DTscale
C
C     allow for 20dB attenuation to be introduced i.e. multiply by 100
C
      Tscale = Tscale*100.0
C
      Return
      End
C
C
      SUBROUTINE DOSCAN (Nbins, SizeLat, Base, Out, Count, FAULT)
C     =================
C     To sample the DVM and PC-Steer co-ordinates during a scan
C     and then bin and average them as required.
C
C     inputs :    Nbins    integer     No of output bins required
C                 Sizelat  real        Eclip. Lat. length of scan
C                 Base     real        DVM output reading base level
C
C     outputs:    Out(*)   real        Binned averaged DVM output
C                 Count(*) integer     No of values averaged in bin
C                 FAULT    logical     Unrecoverable PC-Steer Error
C
      Implicit None
C
      Integer  Nbins, Count(*), AntennaLU, FlukeLU, i, ii,
     *         Istat, Isam, Nempty
C
      Real     X, DVM(1000), Y(1000), Out(*), SizeLat, Base
C
      Logical  FAULT
C
      DATA     AntennaLU /39/, FlukeLU /35/
C
      Isam = 0
1       Isam = Isam + 1
        Read (AntennaLU,'(i1,157x,f8.3)') Istat, Y(Isam)
        Call UNLT
        Read (FlukeLU,*) DVM(Isam)
        If (Istat .lt. 2) go to 1
C
      If (Istat .ne. 2) then
         call ERROR (Istat,FAULT)
         if (FAULT) goto 999
      endif
C
      Do i = 1 , Nbins
        Count(i) = 0
        Out(i) = 0.0
      End Do
C
      Do i = 1 , Isam
        ii = Nint ((SizeLat/2.0 + Y(i)) / SizeLat * Nbins)
        If (ii .ge. 1 .and. ii .le. Nbins) then
           Count(ii) = Count(ii) + 1
           Out(ii) = Out(ii) + DVM(i)-Base
        endif
      End Do
C
      Nempty=0
      Do i = 1 , Nbins
        If (Count(i) .ne. 0) then
           Out(i) = Out(i) / Count(i)
        else
           Nempty=Nempty+1
        endif
      End do
      Print *,'( ',Nempty,' empty bins )'
C
  999 Return
      End
C
C
      SUBROUTINE OutScan (iscan, Nscan, Nbins, Tscale, Zscale,
     *                    Out, Count)
C     ===================
C     To output the binned data for each scan both into the output
C     file and to the plotter in the form of a continuous line.
C
C     inputs :    iscan    integer     Present scan number
C                 Nscan    integer     Total no of scans in map
C                 Nbins    integer     No of output bins in scan
C                 Tscale   real        Kelvin/Volt calibration factor
C                 Zscale   real        Plot-Unit/Volt scale factor
C                 Out(*)   real        Binned averaged DVM output
C                 Count(*) integer     No of values in bin
C
C     outputs :   None except to file and plotter
C
      Implicit None
C
      Integer  iscan, Nscan, Nbins, Count(*), PlotLU, i
C
      Real     Tscale, Zscale, Out(*)
C
      DATA     PlotLU /31/
C
      Write (PlotLU,'("pu;")')
      If (iscan .lt. 0) Write (PlotLU,'("sp 1;")')
      If (iscan .eq. 0) Write (PlotLU,'("sp 2;")')
      If (iscan .gt. 0) Write (PlotLU,'("sp 3;")')
C
      Do i = 1, Nbins
        If (Count(i) .gt. 0) then
           Write (PlotLU,'("pa ",i6,",",i6,";pd;")')
     *        Nint(FLOAT(i)/Nbins*9000+3000),
     *        Nint(Out(i)*Zscale) +
     *        FLOAT(Iscan+Nscan/2)/Nscan*6000 + 500
           Write (20,'(I3,'','',I3,'','',F12.5,'','',I5)')
     *          -iscan, i-Nbins/2-1, Out(i)*Tscale, Count(i)
        else
           Write (20,'(I3,'','',I3,'','',12X,'','',I5)')
     *          -iscan, i-Nbins/2-1, Count(i)
        end if
      end do
      Write (PlotLU,'("pu;")')
      Write (PlotLU,'("sp 0;")')
C
      Return
      End
C
C
      Subroutine UNLT
C     ===============
C     To reset the PC-Steer HPIB card to UNLISTEN mode
C
C     inputs :    None
C
C     outputs:    None except on HPIB
C
      Integer Cmd(2)
C
      DATA Cmd /2, 2H_?/
C
      Call CMDR (30,Cmd,0)
C
      Return
      End
C
C
      Subroutine MJDAY (IYR, ITIME, MJD)
C     ================
C     To calculate the Modified Julian Day, as defined in the
C     Astronomical Almanac, for a given year, day & time
C
C     inputs :    IYR      Integer*2         Year ( from 1983 - )
C                 ITIME(5) Integer*2         Day No, Hr, Min, Sec, 100th
C
C     outputs:    MJD      Double Precision  Modified Julian Day Number
C
C
      Implicit None
C
      Integer*2         IYR, ITIME(5)
C
      Integer           JYR
C
      Double Precision  MJD
C
      Logical           LEAPYR
C
      MJD = ( ( (DBLE(ITIME(1))/100D0 + DBLE(ITIME(2)) ) / 60.0D0
     *                                  + DBLE(ITIME(3)) ) / 60.0D0
     *                                    + DBLE(ITIME(4)) ) / 24.0D0
     *                                      + DBLE(ITIME(5))
C
      IF (IYR .lt. 1983)  then
C
        Print '(I4," - year not catered for in MJD routine MJDAY")',
     *    IYR
        MJD = 0.0D0
C
      else
C
        MJD = MJD + 45334.0D0      ! 1983 offset
        DO 100 JYR = 1984, IYR
          MJD = MJD + 365.0D0
          IF (LEAPYR(JYR-1))  MJD = MJD + 1.0D0
  100   CONTINUE
C
      endif
C
      Return
      End
C
C
      Logical Function LEAPYR (IYEAR)
C     =======================
C     To determine if any given year is a leap year
C
C     inputs :    IYEAR   Integer  Year to be assessed
C
C     outputs:    LEAPYR  Logical  .TRUE. if IYEAR is a leap year
C
      Implicit None
C
      Integer IYEAR
C
      LEAPYR = (((MOD(IYEAR,4) .EQ. 0) .AND. (MOD(IYEAR,100) .NE. 0))
     *  .OR. (MOD(IYEAR,400) .EQ. 0))
C
      Return
      End
C
C
      Subroutine ERROR (Istat,FAULT)
C     ================
C     To report PC-Steer Errors and request User intervention
C
C     inputs :    Istat   Integer  PC-Steer Error Number
C
C     outputs:    FAULT   Logical  Condition Fatal so terminate
C
      Implicit None
C
      Integer  Istat
C
      Character CHAR
  
      Logical  FAULT
C
      Print *,' '
      Print *,'!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!'
      Print *,'!!! PC-Steer Error',Istat,' has been detected !!'
      Print *,'!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!'
      Print *,' '
C
      FAULT = .FALSE.
1     Print *,'Should Map be terminated  ( <cr> = ',FAULT,' ) ? _'
      Read (*,'(a)') CHAR
      Call Upcase(CHAR)
      If (CHAR .eq. 'Y' .or. CHAR .eq. 'T' ) then
         FAULT =.TRUE.
         goto 999
      endif
      If (CHAR .ne. 'N' .and. CHAR .ne. 'F' .and. CHAR .ne. '/' .and.
     *    CHAR .ne. ' ' ) goto 1
C
  999 Return
      End
C
C=======================
      INCLUDE &RHYP::16
      INCLUDE UPCASE::16
C=======================
