Skip to content

Commit

Permalink
Give the two convert wind routines unique names.
Browse files Browse the repository at this point in the history
Update unit tests for new specification of wind fields.

Fixes ufs-community#633.
  • Loading branch information
GeorgeGayno-NOAA committed Jan 30, 2023
1 parent 5e6bcf8 commit d45f3c8
Show file tree
Hide file tree
Showing 5 changed files with 77 additions and 51 deletions.
23 changes: 10 additions & 13 deletions sorc/chgres_cube.fd/atm_input_data.F90
Original file line number Diff line number Diff line change
Expand Up @@ -66,7 +66,7 @@ module atm_input_data

public :: read_input_atm_data
public :: cleanup_input_atm_data
public :: convert_winds
public :: convert_winds_to_xyz

contains

Expand Down Expand Up @@ -428,7 +428,7 @@ subroutine read_input_atm_gfs_sigio_file(localpet)
! Convert from 2-d to 3-d component winds.
!---------------------------------------------------------------------------

call convert_winds
call convert_winds_to_xyz

!---------------------------------------------------------------------------
! Compute 3-d pressure from 'ak' and 'bk'.
Expand Down Expand Up @@ -688,7 +688,7 @@ subroutine read_input_atm_gfs_gaussian_nemsio_file(localpet)
! Convert from 2-d to 3-d component winds.
!---------------------------------------------------------------------------

call convert_winds
call convert_winds_to_xyz

!---------------------------------------------------------------------------
! Compute 3-d pressure from 'ak' and 'bk'.
Expand Down Expand Up @@ -953,7 +953,7 @@ subroutine read_input_atm_gaussian_nemsio_file(localpet)
! Convert from 2-d to 3-d component winds.
!---------------------------------------------------------------------------

call convert_winds
call convert_winds_to_xyz

!---------------------------------------------------------------------------
! Compute 3-d pressure. Mid-layer and surface pressure are computed
Expand Down Expand Up @@ -1237,7 +1237,7 @@ subroutine read_input_atm_restart_file(localpet)
! Convert from 2-d to 3-d cartesian winds.
!---------------------------------------------------------------------------

call convert_winds
call convert_winds_to_xyz

!---------------------------------------------------------------------------
! Compute pressures
Expand Down Expand Up @@ -1607,7 +1607,7 @@ subroutine read_input_atm_gaussian_netcdf_file(localpet)
! Convert from 2-d to 3-d cartesian winds.
!---------------------------------------------------------------------------

call convert_winds
call convert_winds_to_xyz

!---------------------------------------------------------------------------
! Compute pressure.
Expand Down Expand Up @@ -1913,7 +1913,7 @@ subroutine read_input_atm_tiled_history_file(localpet)
! Convert from 2-d to 3-d cartesian winds.
!---------------------------------------------------------------------------

call convert_winds
call convert_winds_to_xyz

!---------------------------------------------------------------------------
! Compute pressure.
Expand Down Expand Up @@ -2868,7 +2868,7 @@ subroutine read_input_atm_grib2_file(localpet)
! Convert from 2-d to 3-d component winds.
!---------------------------------------------------------------------------

call convert_winds
call convert_winds_to_xyz

!---------------------------------------------------------------------------
! Convert dpdt to dzdt if needed
Expand Down Expand Up @@ -3103,7 +3103,7 @@ end subroutine read_winds
!> Convert winds from 2-d to 3-d components.
!!
!! @author George Gayno NCEP/EMC
subroutine convert_winds
subroutine convert_winds_to_xyz

implicit none

Expand Down Expand Up @@ -3168,11 +3168,8 @@ subroutine convert_winds
latrad = latptr(i,j) * acos(-1.) / 180.0
lonrad = lonptr(i,j) * acos(-1.) / 180.0
do k = clb(3), cub(3)
! windptr(i,j,k,1) = uptr(i,j,k) * cos(lonrad) - vptr(i,j,k) * sin(latrad) * sin(lonrad)
xptr(i,j,k) = uptr(i,j,k) * cos(lonrad) - vptr(i,j,k) * sin(latrad) * sin(lonrad)
! windptr(i,j,k,2) = uptr(i,j,k) * sin(lonrad) + vptr(i,j,k) * sin(latrad) * cos(lonrad)
yptr(i,j,k) = uptr(i,j,k) * sin(lonrad) + vptr(i,j,k) * sin(latrad) * cos(lonrad)
! windptr(i,j,k,3) = vptr(i,j,k) * cos(latrad)
zptr(i,j,k) = vptr(i,j,k) * cos(latrad)
enddo
enddo
Expand All @@ -3181,7 +3178,7 @@ subroutine convert_winds
call ESMF_FieldDestroy(u_input_grid, rc=rc)
call ESMF_FieldDestroy(v_input_grid, rc=rc)

end subroutine convert_winds
end subroutine convert_winds_to_xyz

!> Compute grid rotation angle for non-latlon grids.
!!
Expand Down
6 changes: 3 additions & 3 deletions sorc/chgres_cube.fd/atmosphere.F90
Original file line number Diff line number Diff line change
Expand Up @@ -458,7 +458,7 @@ subroutine atmosphere_driver(localpet)
! Convert from 3-d to 2-d cartesian winds.
!-----------------------------------------------------------------------------------

call convert_winds
call convert_winds_to_uv

!-----------------------------------------------------------------------------------
! If selected, process thompson microphysics climatological fields.
Expand Down Expand Up @@ -776,7 +776,7 @@ end subroutine create_atm_esmf_fields
!> Convert 3-d component winds to u and v.
!!
!! @author George Gayno
subroutine convert_winds
subroutine convert_winds_to_uv

implicit none

Expand Down Expand Up @@ -912,7 +912,7 @@ subroutine convert_winds
enddo
enddo

end subroutine convert_winds
end subroutine convert_winds_to_uv

!> Computes 3-D pressure given an adjusted surface pressure.
!!
Expand Down
52 changes: 39 additions & 13 deletions tests/chgres_cube/ftst_convert_winds.F90
Original file line number Diff line number Diff line change
Expand Up @@ -14,8 +14,10 @@ program winds
latitude_input_grid, &
longitude_input_grid

use atm_input_data, only : lev_input, convert_winds, &
wind_input_grid, &
use atm_input_data, only : lev_input, convert_winds_to_xyz, &
xwind_input_grid, &
ywind_input_grid, &
zwind_input_grid, &
u_input_grid, &
v_input_grid

Expand All @@ -26,15 +28,17 @@ program winds

real, parameter :: EPSILON=0.0001

integer :: clb(4), cub(4)
integer :: clb(3), cub(3)
integer :: ierr, localpet, npets, rc
integer :: i, j, k

real(esmf_kind_r8), allocatable :: latitude(:,:)
real(esmf_kind_r8), allocatable :: longitude(:,:)
real(esmf_kind_r8), allocatable :: u_wind(:,:,:)
real(esmf_kind_r8), allocatable :: v_wind(:,:,:)
real(esmf_kind_r8), pointer :: windptr(:,:,:,:)
real(esmf_kind_r8), pointer :: xwindptr(:,:,:)
real(esmf_kind_r8), pointer :: ywindptr(:,:,:)
real(esmf_kind_r8), pointer :: zwindptr(:,:,:)

real :: expected_x_component(IPTS,JPTS)
real :: expected_y_component(IPTS,JPTS)
Expand Down Expand Up @@ -95,11 +99,23 @@ program winds
name="input_grid_longitude", &
rc=rc)

wind_input_grid = ESMF_FieldCreate(input_grid, &
xwind_input_grid = ESMF_FieldCreate(input_grid, &
typekind=ESMF_TYPEKIND_R8, &
staggerloc=ESMF_STAGGERLOC_CENTER, &
ungriddedLBound=(/1,1/), &
ungriddedUBound=(/lev_input,3/), rc=rc)
ungriddedLBound=(/1/), &
ungriddedUBound=(/lev_input/), rc=rc)

ywind_input_grid = ESMF_FieldCreate(input_grid, &
typekind=ESMF_TYPEKIND_R8, &
staggerloc=ESMF_STAGGERLOC_CENTER, &
ungriddedLBound=(/1/), &
ungriddedUBound=(/lev_input/), rc=rc)

zwind_input_grid = ESMF_FieldCreate(input_grid, &
typekind=ESMF_TYPEKIND_R8, &
staggerloc=ESMF_STAGGERLOC_CENTER, &
ungriddedLBound=(/1/), &
ungriddedUBound=(/lev_input/), rc=rc)

u_input_grid = ESMF_FieldCreate(input_grid, &
typekind=ESMF_TYPEKIND_R8, &
Expand Down Expand Up @@ -210,21 +226,31 @@ program winds

! Call the routine to unit test.

call convert_winds
call convert_winds_to_xyz

call ESMF_FieldGet(xwind_input_grid, &
computationalLBound=clb, &
computationalUBound=cub, &
farrayPtr=xwindptr, rc=rc)

call ESMF_FieldGet(ywind_input_grid, &
computationalLBound=clb, &
computationalUBound=cub, &
farrayPtr=ywindptr, rc=rc)

call ESMF_FieldGet(wind_input_grid, &
call ESMF_FieldGet(zwind_input_grid, &
computationalLBound=clb, &
computationalUBound=cub, &
farrayPtr=windptr, rc=rc)
farrayPtr=zwindptr, rc=rc)

print*,"Check results."

do j = clb(2), cub(2)
do i = clb(1), cub(1)
do k = clb(3), cub(3)
if (abs(windptr(i,j,k,1) - expected_x_component(i,j)) > EPSILON) stop 2
if (abs(windptr(i,j,k,2) - expected_y_component(i,j)) > EPSILON) stop 3
if (abs(windptr(i,j,k,3) - expected_z_component(i,j)) > EPSILON) stop 4
if (abs(xwindptr(i,j,k) - expected_x_component(i,j)) > EPSILON) stop 2
if (abs(ywindptr(i,j,k) - expected_y_component(i,j)) > EPSILON) stop 3
if (abs(zwindptr(i,j,k) - expected_z_component(i,j)) > EPSILON) stop 4
enddo
enddo
enddo
Expand Down
24 changes: 13 additions & 11 deletions tests/chgres_cube/ftst_read_atm_gaussian_netcdf.F90
Original file line number Diff line number Diff line change
Expand Up @@ -26,7 +26,9 @@ program read_atm_gaussian_netcdf
dzdt_input_grid, &
ps_input_grid, &
pres_input_grid, &
wind_input_grid, &
xwind_input_grid, &
ywind_input_grid, &
zwind_input_grid, &
terrain_input_grid, &
tracers_input_grid

Expand Down Expand Up @@ -67,7 +69,6 @@ program read_atm_gaussian_netcdf

real(esmf_kind_r8), allocatable :: data_one_tile(:,:)
real(esmf_kind_r8), allocatable :: data3d_one_tile(:,:,:)
real(esmf_kind_r8), allocatable :: data4d_one_tile(:,:,:,:)
real(esmf_kind_r8), allocatable :: latitude(:,:)
real(esmf_kind_r8), allocatable :: longitude(:,:)

Expand Down Expand Up @@ -151,7 +152,6 @@ program read_atm_gaussian_netcdf
if (lev_input /= EXPECTED_LEV_INPUT) stop 2
if (levp1_input /= EXPECTED_LEVP1_INPUT) stop 3

allocate(data4d_one_tile(i_input,j_input,lev_input,3))
allocate(data3d_one_tile(i_input,j_input,lev_input))
allocate(data_one_tile(i_input,j_input))

Expand Down Expand Up @@ -187,21 +187,23 @@ program read_atm_gaussian_netcdf
if (abs(data3d_one_tile(1,1,1) - expected_values_pres(1)) > EPSILON) stop 18
if (abs(data3d_one_tile(i_input,j_input,lev_input) - expected_values_pres(2)) > EPSILON) stop 19

call ESMF_FieldGather(wind_input_grid, data4d_one_tile, rootPet=0, rc=rc)
if (abs(data4d_one_tile(1,1,1,1) - expected_values_xwind(1)) > EPSILON) stop 20
if (abs(data4d_one_tile(1,1,lev_input,1) - expected_values_xwind(2)) > EPSILON) stop 21
if (abs(data4d_one_tile(1,1,1,2) - expected_values_ywind(1)) > EPSILON) stop 22
if (abs(data4d_one_tile(1,1,lev_input,2) - expected_values_ywind(2)) > EPSILON) stop 23
if (abs(data4d_one_tile(1,1,1,3) - expected_values_zwind(1)) > EPSILON) stop 24
if (abs(data4d_one_tile(1,1,lev_input,3) - expected_values_zwind(2)) > EPSILON) stop 25
call ESMF_FieldGather(xwind_input_grid, data3d_one_tile, rootPet=0, rc=rc)
if (abs(data3d_one_tile(1,1,1) - expected_values_xwind(1)) > EPSILON) stop 20
if (abs(data3d_one_tile(1,1,lev_input) - expected_values_xwind(2)) > EPSILON) stop 21
call ESMF_FieldGather(ywind_input_grid, data3d_one_tile, rootPet=0, rc=rc)
if (abs(data3d_one_tile(1,1,1) - expected_values_ywind(1)) > EPSILON) stop 22
if (abs(data3d_one_tile(1,1,lev_input) - expected_values_ywind(2)) > EPSILON) stop 23
call ESMF_FieldGather(zwind_input_grid, data3d_one_tile, rootPet=0, rc=rc)
if (abs(data3d_one_tile(1,1,1) - expected_values_zwind(1)) > EPSILON) stop 24
if (abs(data3d_one_tile(1,1,lev_input) - expected_values_zwind(2)) > EPSILON) stop 25

call ESMF_FieldGather(terrain_input_grid, data_one_tile, rootPet=0, rc=rc)
if (abs(data_one_tile(1,1) - expected_values_terrain(1)) > EPSILON) stop 26
if (abs(data_one_tile(i_input,j_input) - expected_values_terrain(2)) > EPSILON) stop 27

print*,"OK"

deallocate(latitude, longitude, data4d_one_tile, data3d_one_tile, data_one_tile)
deallocate(latitude, longitude, data3d_one_tile, data_one_tile)

call ESMF_finalize(endflag=ESMF_END_KEEPMPI)

Expand Down
23 changes: 12 additions & 11 deletions tests/chgres_cube/ftst_read_atm_grib2.F90
Original file line number Diff line number Diff line change
Expand Up @@ -16,7 +16,8 @@ program read_atm_grib2
levp1_input, &
temp_input_grid, tracers_input_grid, &
dzdt_input_grid, pres_input_grid, &
ps_input_grid, wind_input_grid, &
ps_input_grid, xwind_input_grid, &
ywind_input_grid, zwind_input_grid, &
terrain_input_grid

use program_setup, only : input_type, data_dir_input_grid, &
Expand Down Expand Up @@ -52,7 +53,6 @@ program read_atm_grib2
real(esmf_kind_r8), allocatable :: longitude(:,:)
real(esmf_kind_r8), allocatable :: data_one_tile(:,:)
real(esmf_kind_r8), allocatable :: data3d_one_tile(:,:,:)
real(esmf_kind_r8), allocatable :: data4d_one_tile(:,:,:,:)

real :: expected_values_tmp(NUM_VALUES)
real :: expected_values_sphum(NUM_VALUES)
Expand Down Expand Up @@ -174,7 +174,6 @@ program read_atm_grib2
if (levp1_input /= EXPECTED_LEVP1_INPUT) stop 3

allocate(data3d_one_tile(i_input,j_input,lev_input))
allocate(data4d_one_tile(i_input,j_input,lev_input,3))
allocate(data_one_tile(i_input,j_input))

! The i/j/k of the points to be checked.
Expand Down Expand Up @@ -231,21 +230,23 @@ program read_atm_grib2
if (abs(data3d_one_tile(i_check(1),j_check(1),k_check(1)) - expected_values_pres(1)) > EPSILON) stop 22
if (abs(data3d_one_tile(i_check(2),j_check(2),k_check(2)) - expected_values_pres(2)) > EPSILON) stop 23

call ESMF_FieldGather(wind_input_grid, data4d_one_tile, rootPet=0, rc=rc)
if (abs(data4d_one_tile(i_check(1),j_check(1),k_check(1),1) - expected_values_xwind(1)) > EPSILON) stop 24
if (abs(data4d_one_tile(i_check(2),j_check(2),k_check(2),1) - expected_values_xwind(2)) > EPSILON) stop 25
if (abs(data4d_one_tile(i_check(1),j_check(1),k_check(1),2) - expected_values_ywind(1)) > EPSILON) stop 26
if (abs(data4d_one_tile(i_check(2),j_check(2),k_check(2),2) - expected_values_ywind(2)) > EPSILON) stop 27
if (abs(data4d_one_tile(i_check(1),j_check(1),k_check(1),3) - expected_values_zwind(1)) > EPSILON) stop 28
if (abs(data4d_one_tile(i_check(2),j_check(2),k_check(2),3) - expected_values_zwind(2)) > EPSILON) stop 29
call ESMF_FieldGather(xwind_input_grid, data3d_one_tile, rootPet=0, rc=rc)
if (abs(data3d_one_tile(i_check(1),j_check(1),k_check(1)) - expected_values_xwind(1)) > EPSILON) stop 24
if (abs(data3d_one_tile(i_check(2),j_check(2),k_check(2)) - expected_values_xwind(2)) > EPSILON) stop 25
call ESMF_FieldGather(ywind_input_grid, data3d_one_tile, rootPet=0, rc=rc)
if (abs(data3d_one_tile(i_check(1),j_check(1),k_check(1)) - expected_values_ywind(1)) > EPSILON) stop 26
if (abs(data3d_one_tile(i_check(2),j_check(2),k_check(2)) - expected_values_ywind(2)) > EPSILON) stop 27
call ESMF_FieldGather(zwind_input_grid, data3d_one_tile, rootPet=0, rc=rc)
if (abs(data3d_one_tile(i_check(1),j_check(1),k_check(1)) - expected_values_zwind(1)) > EPSILON) stop 28
if (abs(data3d_one_tile(i_check(2),j_check(2),k_check(2)) - expected_values_zwind(2)) > EPSILON) stop 29

call ESMF_FieldGather(ps_input_grid, data_one_tile, rootPet=0, rc=rc)
if (abs(data_one_tile(i_check(1),j_check(1)) - expected_values_ps(1)) > EPSILON) stop 32

call ESMF_FieldGather(terrain_input_grid, data_one_tile, rootPet=0, rc=rc)
if (abs(data_one_tile(i_check(1),j_check(1)) - expected_values_terrain(1)) > EPSILON) stop 34

deallocate(latitude, longitude, data3d_one_tile, data4d_one_tile, data_one_tile)
deallocate(latitude, longitude, data3d_one_tile, data_one_tile)

call ESMF_finalize(endflag=ESMF_END_KEEPMPI)

Expand Down

0 comments on commit d45f3c8

Please sign in to comment.