da_get_height.inc

1 subroutine da_get_height( input_file, dim1, dim2, dim3, height)
2 
3    !---------------------------------------------------------------------------
4    ! Purpose: Calculates T, RH from input WRF file.
5    !---------------------------------------------------------------------------
6    
7    implicit none
8    
9    character(len=200), intent(in)  :: input_file       ! NETCDF file nane.
10    integer,            intent(in)  :: dim1, dim2, dim3          ! Dimensions.
11    real,               intent(out) :: height(1:dim1,1:dim2,1:dim3) ! Height.
12 
13    character(len=10) :: var                            ! Variable to search for.
14    integer           :: k                              ! Loop counter.
15    real              :: gravity_inv                    ! 1/gravity.
16    real              :: phb(1:dim1,1:dim2)             ! Base state geopotential.
17    real              :: ph(1:dim1,1:dim2)              ! Perturbation geopotential.
18    real              :: phfull(1:dim1,1:dim2,1:dim3+1) ! Geopotential.
19 
20    if (trace_use) call da_trace_entry("da_get_height")
21 
22    gravity_inv = 1.0 / gravity
23 
24    do k = 1, dim3+1
25       var = "PHB"
26       call da_get_field( input_file, var, 3, dim1, dim2, dim3+1, k, phb)
27       var = "PH"
28       call da_get_field( input_file, var, 3, dim1, dim2, dim3+1, k, ph)
29 
30       phfull(:,:,k) = phb + ph ! Calculate full geopotential on full(w) model levels:
31    end do
32 
33    do k = 1, dim3
34       height(:,:,k) = 0.5 *( phfull(:,:,k+1) + phfull(:,:,k)) * gravity_inv
35    end do
36 
37    if (trace_use) call da_trace_exit("da_get_height")
38 
39 end subroutine da_get_height
40 
41