anchor.diff:    DEALLOCATE(Total_predicted_freq )
     DEALLOCATE( HR_params )
   
  +! Feel free to cut the next five lines - they are just for calculating run time.
  +  CALL  date_and_time(date, time, zone, end_value)
  +  WRITE(Ulog,*) ''
  +  WRITE(Ulog,70) (end_value(5) - start_value(5)),  (end_value(6) - start_value(6)), &
  +       (end_value(7) - start_value(7)),  (end_value(8) - start_value(8))
  +70 FORMAT('Total computation time is: ',I2.2,' hrs, ',I2.2, ' min, ',I2.2, ' s, ',I4.4 ,' ms')  
  +! Close the debug log
  +  CLOSE(Ulog)
     RETURN
  +  
   9998 CONTINUE
     WRITE(0,*) 'End of input stream.'
     ! Close the input/output file units
  
positive.code:   INTEGER :: Uin  ! input file unit number
    INTEGER :: Uout ! output file unit number
  
  
    ! these variables are set to the common block.
    ! They are required by the fitting routine within DQED
    COMMON N_vib, N_rot, nu_low, nu_high, nu_mid, cp_difference, CV_temps
  
  
negative.code:   READ (Uin,*,END=9998) cp_data, data
   ! WRITE (0,*) "cp_data is ",cp_data
   ! WRITE (0,'(A,28I3)') "data is ",data
  
  ! The program may pause at the previous line until it is given more data.
  ! Now that we have some data to start processing,
  ! we can start the clock and open the log file.
  
    ! Feel free to cut this.  I use it to determine how long the code runs.  
    CALL  date_and_time(date, time, zone, start_value)
  
    ! Open the debug log file
    ! use /dev/null to silence
    open (UNIT=Ulog, FILE='frankie.log')
      
    ! call the subroutine which separates the data info and calculates
    ! how many RRHO frequencies it can determine from the structure
    call read_bonds(data, bond_info,  degeneracy)
  