C test program for 'FORTRAN' version of height adjustment component 
C of air/sea boundary layer model 


      IMPLICIT NONE

C  feeds input into the sea state / flux subroutine   
      integer astab, air_moist_prm, count, dyn_in_prm, eqv_neut, 
     +        flux_model, ht_adj, ht_adj_, i, num, sfc_moist_prm, 
     &        ss_prm, sst_prm, warn
      real    dom_phs_spd, h_sig, lhf, q_at_z, q_star, shf, t_at_z, 
     +        tau(2), t_star, u_at_z, 
     +        u_star(2), wave_age, ww_stab, z_over_L, zo_m(2)
      real    air_moist_val, CONVECT, CONV_CRIT, dtr, dyn_in_val, 
     &        dyn_in_val2, 
     +        pressure, Qnet, ref_ht_tq, ref_ht_wind, rel_wind_ang, 
     +        salinity, sfc_moist_val, ss_val, t_air, t_skin, 
     +        tau_x, tau_y, wave_ang, wind_ang, z_wanted

      OPEN( unit=9, file='testdata04.dat', status = 'old' )

C     convergence critereon (fractional change)  []  
      CONV_CRIT = 0.00005

C     convective parameter 
      CONVECT = 1.25

C     warning level (1 warnings; 0 no warnings)
      warn = 0 

C     flag for wind output: 0 = wind speed, 1 = equivalent neutral wind speed
      eqv_neut = 0

C     Height for which winds, potential temperature, and humidity are adjusted.
      z_wanted = 15.0

C     conversion from degrees to radians
      dtr = 3.14159 / 180.0

C     begin loop for series of test_data 
C     skip the header line 
      READ( 9, 10)

   10 FORMAT( ) 
   11 FORMAT( A,A,A ) 
   12 FORMAT( i2, 1x, F5.2, 1x, F6.3, 1x, F6.3, 1x, F6.3, 1x, F7.4, 
     + 1x, F9.6, 1x, F8.5, 1x, F6.2, 1x, F6.2, 1x, F5.2, 1x, F6.3, 1x, 
     + F6.3, 1x, F6.2, 1x, F6.2, 1x, F6.2, 1x, F6.2, 1x, F6.3 )  
 
      WRITE(*,11) 'run  U   |ustar| ustar1 ustar2 tstar   qstar     ',
     + 'zref/L     cp     wa   Hsig   tau1   tau2    shf   lhf ',
     + '   u(z)   t(z)   q(z)'

      DO 100 i=1, 62

        READ(9,*) num, dyn_in_prm, dyn_in_val, dyn_in_val2, wave_ang, 
     +   ss_prm, ss_val, air_moist_prm, air_moist_val, sfc_moist_prm, 
     +   sfc_moist_val, t_skin, t_air, ref_ht_wind, ref_ht_tq, pressure, 
     +   salinity, CONVECT, astab
C23456789012345678901234567890123456789012345678901234567890123456789012

        count = ht_adj( %val(dyn_in_prm), %val(dyn_in_val), 
     &   %val(dyn_in_val2), %val(CONVECT), 
     &   %val(CONV_CRIT), %val(pressure), %val(air_moist_prm), 
     &   %val(air_moist_val), %val(sfc_moist_prm), %val(sfc_moist_val), 
     &   %val(salinity), %val(ss_prm), %val(ss_val), %val(t_air), 
     &   %val(sst_prm), %val(t_skin), %val(ref_ht_wind), 
     &   %val(ref_ht_tq), %val(astab), %val(Qnet), %val(warn), 
     &   shf, lhf, tau, u_star, t_star, q_star, z_over_L, wave_age, 
     &   dom_phs_spd, h_sig, ww_stab, zo_m, %val(eqv_neut), 
     &   %val(z_wanted), u_at_z, 
     &   t_at_z, q_at_z, %val(flux_model) )

C        check for convergence
         if ( count .EQ. 1 ) print*, 'non-convergence '  
         WRITE(*,12) i, dyn_in_val, sqrt( u_star(1) * u_star(1) + 
     +    u_star(2) * u_star(2) ), u_star(1), u_star(2), t_star, 
     +    q_star, z_over_L, dom_phs_spd, wave_age, h_sig, tau(1), 
     +    tau(2), shf, lhf, u_at_z, t_at_z, q_at_z 

  100 CONTINUE 

C     end of test program: ftest_ht_adj99.f
      END