ICAOHEIGHT.f

      SUBROUTINE icaoheight(MAXWP,     & !input
                            MAXWICAOZ)   ! output
 
! SUBPROGRAM:    ICAOHEIGHT      
!   PRGRMMR: CHUANG         ORG: W/NP2      DATE: 09-05-08
!     
! ABSTRACT:
!     THIS ROUTINE CONVERT PRESSURE FIELDS TO ICAO HEIGHT
!
 
! Description:
!   Convert pressure (Pa) to height (m) using ICAO standard atmosphere
!   adpt code from uk
! Method:
!
!
! Code Description:
!   Language:           Fortran 90
!   Software Standards: UMDP3 v6
 
use ctlblk_mod, only: jsta, jend, spval, im, ista, iend
use physcons_post, only: con_g, con_rd
IMPLICIT None
 
! Subroutine Arguments:
!REAL, INTENT(IN) :: SPVAL
REAL, INTENT(IN) :: MAXWP(ista:iend,jsta:jend)       !P field for conversion
 
REAL, INTENT(INOUT) :: MAXWICAOZ(ista:iend,jsta:jend)   !ICAO height in m
!INTEGER, INTENT(INOUT) :: ErrorStatus
 
! Local Constants:
REAL, PARAMETER :: G_over_R = con_g / con_rd
REAL, PARAMETER :: Lapse_RateL = 6.5e-03  ! For levels below 11,000 gpm
REAL, PARAMETER :: Lapse_RateU = -1.0e-03 ! For levels above 11,000 gpm
REAL, PARAMETER :: Press_Bot = 101325.    ! ICAO std: surface pressure
REAL, PARAMETER :: Press_Mid = 22632.     !      pressure @ 11,000 gpm
REAL, PARAMETER :: Press_Top = 5474.87    !      pressure @ 20,000 gpm
REAL, PARAMETER :: Temp_Bot = 288.15      ! Surface temperature
REAL, PARAMETER :: Temp_Top = 216.65      ! Temperature of isotherm
REAL, PARAMETER :: Gpm1 = 11000.0  ! Ht limit (gpm) for std lower
                                   ! lapse rate
REAL, PARAMETER :: Gpm2 = 20000.0  ! Ht (gpm) of top of isothermal layer
REAL, PARAMETER :: ZP1 = lapse_ratel/g_over_r ! Exponents used for
REAL, PARAMETER :: ZP2 = lapse_rateu/g_over_r ! calculation
 
! Local Variables:
INTEGER :: i,j      ! loop counters
REAL :: Pressure    ! Local pressure
 
! End of header --------------------------------------------------------
 
!IF ( ErrorStatus /= StatusOK ) THEN
  ! Previous error - do not proceed
!  GO TO 9999
!END IF
 
!IF ( ASSOCIATED( IHField % RData ) ) THEN
!  DEALLOCATE( IHField % RData )
!END IF
 
 
 
DO j=jsta,jend
  DO i=ista,iend
    pressure = maxwp(i,j)
    IF ( (pressure <= 1000.) .AND. (pressure >= 0.) ) THEN
      pressure = 1000.
      print*,'lower ICAO pressure to 10 mb'
    END IF
    IF ( pressure > press_bot .and. pressure<spval) THEN
      pressure = press_bot
    END IF
 
    IF (pressure == spval) THEN
      maxwicaoz(i,j) = spval
    ELSE IF (pressure > press_mid) THEN ! Hts up to 11,000 GPM
      pressure = pressure/press_bot
      pressure = 1.0 - pressure**zp1
      maxwicaoz(i,j) = pressure*temp_bot/lapse_ratel
 
    ELSE IF (pressure > press_top) THEN ! Hts between 11,000
                                        !     and     20,000 GPM
      pressure = pressure/press_mid
      pressure = -alog(pressure)
      maxwicaoz(i,j) = gpm1 + pressure*temp_top/g_over_r
 
    ELSE                                ! Hts above 20,000 GPM
      pressure = pressure/press_top
      pressure = 1.0 - pressure**zp2
      maxwicaoz(i,j) = gpm2 + pressure*temp_top/lapse_rateu
 
    END IF
 
  ENDDO
ENDDO
 
 
!9999 CONTINUE
 
END SUBROUTINE icaoheight