Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Master test #18

Merged
merged 12 commits into from
Feb 21, 2020
147 changes: 135 additions & 12 deletions GFDL_tools/fv_cmip_diag.F90
Original file line number Diff line number Diff line change
Expand Up @@ -28,21 +28,21 @@ module fv_cmip_diag_mod
use fms_io_mod, only: set_domain, nullify_domain, string
use time_manager_mod, only: time_type
use mpp_domains_mod, only: domain2d
use diag_manager_mod, only: register_diag_field, &
use diag_manager_mod, only: register_diag_field, diag_axis_init, &
send_data, get_diag_field_id, &
register_static_field, &
diag_field_add_attribute, &
DIAG_FIELD_NOT_FOUND
use diag_data_mod, only: CMOR_MISSING_VALUE, null_axis_id
use tracer_manager_mod, only: get_tracer_index
use field_manager_mod, only: MODEL_ATMOS
use constants_mod, only: GRAV
use constants_mod, only: GRAV, RDGAS

use fv_mapz_mod, only: E_Flux
use fv_arrays_mod, only: fv_atmos_type
use fv_diagnostics_mod, only: interpolate_vertical, &
get_height_given_pressure, &
rh_calc, get_height_field
rh_calc, get_height_field, get_vorticity

use atmos_cmip_diag_mod, only: register_cmip_diag_field_2d, &
register_cmip_diag_field_3d, &
Expand All @@ -58,7 +58,7 @@ module fv_cmip_diag_mod

public :: fv_cmip_diag_init, fv_cmip_diag, fv_cmip_diag_end

integer :: sphum
integer :: sphum, nql, nqi, nqa

!-----------------------------------------------------------------------
!--- namelist ---
Expand All @@ -71,12 +71,15 @@ module fv_cmip_diag_mod

type(cmip_diag_id_type) :: ID_ta, ID_ua, ID_va, ID_hus, ID_hur, ID_wap, ID_zg, &
ID_u2, ID_v2, ID_t2, ID_wap2, ID_uv, ID_ut, ID_vt, &
ID_uwap, ID_vwap, ID_twap
ID_uwap, ID_vwap, ID_twap, ID_wa, &
ID_cls, ID_clws, ID_clis

integer :: id_ps, id_orog
integer :: id_ua200, id_va200, id_ua850, id_va850, &
id_ta500, id_ta700, id_ta850, id_zg500, &
id_zg100, id_zg10, id_zg1000, &
id_hus850, id_wap500, id_ua10
integer :: id_rv200, id_rv500, id_rv850, id_vortmean

character(len=5) :: mod_name = 'atmos'

Expand Down Expand Up @@ -107,6 +110,7 @@ subroutine fv_cmip_diag_init ( Atm, axes, Time )
integer, dimension(7) :: id_plevels
integer, parameter :: id_p10=1, id_p100=2, id_p200=3, id_p500=4, id_p700=5, id_p850=6, id_p1000=7
character(len=4) :: plabel
integer :: id_pl700, id_pl700_bnds, id_nv
!-----------------------------------------------------------------------

if (module_is_initialized) then
Expand Down Expand Up @@ -151,6 +155,9 @@ subroutine fv_cmip_diag_init ( Atm, axes, Time )
!-----------------------------------------------------------------------

sphum = get_tracer_index (MODEL_ATMOS, 'sphum')
nql = get_tracer_index (MODEL_ATMOS, 'liq_wat')
nqi = get_tracer_index (MODEL_ATMOS, 'ice_wat')
nqa = get_tracer_index (MODEL_ATMOS, 'cld_amt')

!-----------------------------------------------------------------------
! register cmip 3D variables (on model levels and pressure levels)
Expand All @@ -173,6 +180,9 @@ subroutine fv_cmip_diag_init ( Atm, axes, Time )
ID_hur = register_cmip_diag_field_3d (mod_name, 'hur', Time, &
'Relative Humidity', '%', standard_name='relative_humidity')

ID_wa = register_cmip_diag_field_3d (mod_name, 'wa', Time, &
'Upward Air Velocity', 'm s-1', standard_name='upward_air_velocity')

ID_zg = register_cmip_diag_field_3d (mod_name, 'zg', Time, &
'Geopotential Height', 'm', standard_name='geopotential_height', axis='half')

Expand Down Expand Up @@ -215,6 +225,24 @@ subroutine fv_cmip_diag_init ( Atm, axes, Time )
'Air Temperature times Omega', 'K m s-1', &
standard_name='product_of_omega_and_air_temperature')

!-----------------------------------------------------------------------
! stratiform cloud tracers

if (nql > 0) then
ID_clws = register_cmip_diag_field_3d (mod_name, 'clws', Time, &
'Mass Fraction of Stratiform Cloud Liquid Water', '1.0', &
standard_name='mass_fraction_of_stratiform_cloud_liquid_water_in_air')
endif
if (nqi > 0) then
ID_clis = register_cmip_diag_field_3d (mod_name, 'clis', Time, &
'Mass Fraction of Stratiform Cloud Ice', '1.0', &
standard_name='mass_fraction_of_convective_cloud_ice_in_air')
endif
if (nqa > 0) then
ID_cls = register_cmip_diag_field_3d (mod_name, 'cls', Time, &
'Percentage Cover of Stratiform Cloud', '%', &
standard_name='stratiform_cloud_area_fraction_in_atmosphere_layer')
endif
!-----------------------------------------------------------------------
! 2D fields

Expand All @@ -232,14 +260,14 @@ subroutine fv_cmip_diag_init ( Atm, axes, Time )
id_orog = register_static_field (mod_name, 'orog', axes(1:2), &
'Surface Altitude', 'm', &
standard_name='surface_altitude', &
area=area_id)
area=area_id, interp_method='conserve_order1')
if (id_orog > 0) used = send_data (id_orog, Atm(n)%phis(isc:iec,jsc:jec)/GRAV, Time)
#else
!--- for now output this as 'zsurf' from fv_diagnostics ---
! id_orog = register_diag_field (mod_name, 'orog', axes(1:2), Time, &
! 'Surface Altitude', 'm', &
! standard_name='surface_altitude', &
! area=area_id)
! area=area_id, interp_method='conserve_order1')
#endif

!-----------------------------------------------------------------------
Expand All @@ -259,6 +287,24 @@ subroutine fv_cmip_diag_init ( Atm, axes, Time )
endif
enddo

id_pl700 = register_static_field (mod_name, 'pl700', (/null_axis_id/), &
'700 hPa Average', 'Pa', standard_name='air_pressure')
if (id_pl700 > 0) then
call diag_field_add_attribute (id_pl700, 'axis', 'Z')
call diag_field_add_attribute (id_pl700, 'positive', 'down' )
call diag_field_add_attribute (id_pl700, 'comment', 'average at levels 600,700,850 hPa' )
! add bounds
id_nv = diag_axis_init('nv', (/1.,2./), 'none', 'N', 'vertex number', set_name='nv')
id_pl700_bnds = register_static_field (mod_name, 'pl700_bnds', (/id_nv,null_axis_id/), &
'700 hPa boundaries', 'Pa', standard_name='air_pressure')
if (id_pl700_bnds > 0) then
call diag_field_add_attribute (id_pl700, 'bounds', 'pl700_bnds' )
used = send_data (id_pl700_bnds, (/850.e2,600.e2/), Time)
endif
used = send_data (id_pl700, 700.e2, Time)
endif


!---- register field on single pressure levels ----

id_ua10 = register_cmip_diag_field_2d (mod_name, 'ua10', Time, &
Expand Down Expand Up @@ -311,6 +357,30 @@ subroutine fv_cmip_diag_init ( Atm, axes, Time )
if (id_hus850 > 0 .and. id_plevels(id_p850) > 0) &
call diag_field_add_attribute (id_hus850, 'coordinates', 'p850')

!---- relative vorticity at 200, 500, 850 hPa ----
id_rv200 = register_cmip_diag_field_2d (mod_name, 'rv200', Time, &
'Relative Vorticity at 200 hPa', 's-1', standard_name='atmosphere_relative_vorticity')
if (id_rv200 > 0 .and. id_plevels(id_p200) > 0) &
call diag_field_add_attribute (id_rv200, 'coordinates', 'p200')

id_rv500 = register_cmip_diag_field_2d (mod_name, 'rv500', Time, &
'Relative Vorticity at 500 hPa', 's-1', standard_name='atmosphere_relative_vorticity')
if (id_rv500 > 0 .and. id_plevels(id_p500) > 0) &
call diag_field_add_attribute (id_rv500, 'coordinates', 'p500')

id_rv850 = register_cmip_diag_field_2d (mod_name, 'rv850', Time, &
'Relative Vorticity at 850 hPa', 's-1', standard_name='atmosphere_relative_vorticity')
if (id_rv850 > 0 .and. id_plevels(id_p850) > 0) &
call diag_field_add_attribute (id_rv850, 'coordinates', 'p850')

!---- mean relative vorticity 600, 700, 850 hPa ----

id_vortmean = register_cmip_diag_field_2d (mod_name, 'vortmean', Time, &
'Mean Relative Vorticity over 600-850 hPa', 's-1', &
standard_name='atmosphere_relative_vorticity')
if (id_vortmean > 0 .and. id_pl700 > 0) &
call diag_field_add_attribute (id_vortmean, 'coordinates', 'pl700')

!---- omega at 500 hPa ----

id_wap500 = register_cmip_diag_field_2d (mod_name, 'wap500', Time, &
Expand Down Expand Up @@ -357,15 +427,18 @@ subroutine fv_cmip_diag ( Atm, zvir, Time )
integer :: isc, iec, jsc, jec, n, i, j, k, id
integer :: ngc, npz
logical :: used
logical :: compute_wa , compute_rh

real, dimension(Atm(1)%bd%isc:Atm(1)%bd%iec, &
Atm(1)%bd%jsc:Atm(1)%bd%jec) :: pfull, dat2
Atm(1)%bd%jsc:Atm(1)%bd%jec) :: pfull, dat2, &
rv850, rv700, rv600

real, dimension(Atm(1)%bd%isc:Atm(1)%bd%iec, &
Atm(1)%bd%jsc:Atm(1)%bd%jec,1) :: dat3

real, dimension(Atm(1)%bd%isc:Atm(1)%bd%iec, &
Atm(1)%bd%jsc:Atm(1)%bd%jec, &
Atm(1)%npz) :: rhum
Atm(1)%npz) :: rhum, wa, rv

real, dimension(Atm(1)%bd%isc:Atm(1)%bd%iec, &
Atm(1)%bd%jsc:Atm(1)%bd%jec, &
Expand All @@ -384,26 +457,45 @@ subroutine fv_cmip_diag ( Atm, zvir, Time )

call set_domain(Atm(n)%domain)

! set flags for computing quantities
compute_rh = .false.
compute_wa = .false.
if (count(ID_hur%field_id(:)>0) > 0) compute_rh = .true.
if (count(ID_wa%field_id(:)>0) > 0) compute_wa = .true.

! compute relative humidity at model levels (if needed)
if (count(ID_hur%field_id(:)>0) > 0) then
if (compute_rh .or. compute_wa) then
do k=1,npz
do j=jsc,jec
do i=isc,iec
pfull(i,j) = Atm(n)%delp(i,j,k)/(Atm(n)%peln(i,k+1,j)-Atm(n)%peln(i,k,j))
enddo
enddo
call rh_calc (pfull, Atm(n)%pt(isc:iec,jsc:jec,k), &
! compute relative humidity
if (compute_rh) then
call rh_calc (pfull, Atm(n)%pt(isc:iec,jsc:jec,k), &
Atm(n)%q(isc:iec,jsc:jec,k,sphum), rhum(isc:iec,jsc:jec,k), do_cmip=.true.)
endif
! compute vertical velocity
if (compute_wa) then
wa(isc:iec,jsc:jec,k) = -(Atm(n)%omga(isc:iec,jsc:jec,k)*Atm(n)%pt(isc:iec,jsc:jec,k)/ &
pfull(isc:iec,jsc:jec))*(RDGAS/GRAV)
endif
enddo
endif


! height field (wz) if needed
if (count(ID_zg%field_id(:)>0) > 0 .or. any((/id_zg10,id_zg100,id_zg500,id_zg1000/) > 0)) then
call get_height_field(isc, iec, jsc, jec, ngc, npz, Atm(n)%flagstruct%hydrostatic, Atm(n)%delz, &
wz, Atm(n)%pt, Atm(n)%q, Atm(n)%peln, zvir)
endif

! relative vorticity
if (any((/id_rv200,id_rv500,id_rv850,id_vortmean/) > 0)) then
call get_vorticity(isc, iec, jsc, jec, Atm(n)%bd%isd, Atm(n)%bd%ied, Atm(n)%bd%jsd, Atm(n)%bd%jed, npz, &
Atm(n)%u, Atm(n)%v, rv, Atm(n)%gridstruct%dx, Atm(n)%gridstruct%dy, Atm(n)%gridstruct%rarea)
endif

!----------------------------------------------------------------------
! process 2D fields

Expand Down Expand Up @@ -431,6 +523,10 @@ subroutine fv_cmip_diag ( Atm, zvir, Time )
if (query_cmip_diag_id(ID_hur)) &
used = send_cmip_data_3d (ID_hur, rhum(isc:iec,jsc:jec,:), Time, phalf=Atm(n)%peln, opt=1)

! vertical velocity
if (query_cmip_diag_id(ID_wa)) &
used = send_cmip_data_3d (ID_wa, wa(isc:iec,jsc:jec,:), Time, phalf=Atm(n)%peln, opt=1)

! geopotential height
if (query_cmip_diag_id(ID_zg)) &
used = send_cmip_data_3d (ID_zg, wz, Time, phalf=Atm(n)%peln, opt=1, ext=.true.)
Expand Down Expand Up @@ -478,6 +574,13 @@ subroutine fv_cmip_diag ( Atm, zvir, Time )
used = send_cmip_data_3d (ID_twap, Atm(n)%pt (isc:iec,jsc:jec,:)*Atm(n)%omga(isc:iec,jsc:jec,:), &
Time, phalf=Atm(n)%peln, opt=1)

!----------------------------------------------------------------------
! stratiform cloud tracers (only on model levels)

if (query_cmip_diag_id(ID_cls)) used = send_cmip_data_3d (ID_cls, Atm(n)%q(isc:iec,jsc:jec,:,nqa)*100., Time)
if (query_cmip_diag_id(ID_clws)) used = send_cmip_data_3d (ID_clws, Atm(n)%q(isc:iec,jsc:jec,:,nql), Time)
if (query_cmip_diag_id(ID_clis)) used = send_cmip_data_3d (ID_clis, Atm(n)%q(isc:iec,jsc:jec,:,nqi), Time)

!----------------------------------------------------------------------
! process 2D fields on specific pressure levels
!
Expand Down Expand Up @@ -541,6 +644,26 @@ subroutine fv_cmip_diag ( Atm, zvir, Time )
used = send_data (id_wap500, dat2, Time)
endif

if (id_rv200 > 0) then
call interpolate_vertical (isc, iec, jsc, jec, npz, 200.e2, Atm(n)%peln, rv, dat2)
used = send_data (id_rv200, dat2, Time)
endif

if (id_rv500 > 0) then
call interpolate_vertical (isc, iec, jsc, jec, npz, 500.e2, Atm(n)%peln, rv, dat2)
used = send_data (id_rv500, dat2, Time)
endif

if (id_rv850 > 0 .or. id_vortmean > 0 ) then
call interpolate_vertical (isc, iec, jsc, jec, npz, 850.e2, Atm(n)%peln, rv, rv850)
if (id_rv850 > 0) used = send_data (id_rv850, rv850, Time)
if (id_vortmean > 0) then
call interpolate_vertical (isc, iec, jsc, jec, npz, 600.e2, Atm(n)%peln, rv, rv600)
call interpolate_vertical (isc, iec, jsc, jec, npz, 700.e2, Atm(n)%peln, rv, rv700)
used = send_data (id_vortmean, (rv600+rv700+rv850)/3., Time)
endif
endif

if (id_zg10 > 0) then
call get_height_given_pressure (isc, iec, jsc, jec, npz, wz, 1, (/id_zg10/), &
(/log(10.e2)/), Atm(n)%peln, dat3)
Expand Down
24 changes: 12 additions & 12 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -6,26 +6,26 @@ The GFDL Microphysics is also available via this repository.

# Where to find information

See the [FV3 documentation and references](https://www.gfdl.noaa.gov/fv3/fv3-documentation-and-references/) for more information.
See the [FV3 documentation and references](https://www.gfdl.noaa.gov/fv3/fv3-documentation-and-references/)
for more information.

# Proper usage attribution

Cite either Putman and Lin (2007) or Harris and Lin (2013) when describing a model using the FV3 dynamical core.
Cite Chen et al (2013) and Zhou et al (2019) if using the GFDL Microphysics.
Cite Putman and Lin (2007) and Harris and Lin (2013) when describing a model using the FV3 dynamical core.
Cite Chen et al (2013) and Zhou et al (2019) when using the GFDL Microphysics.

# What files are what

The top level directory structure groups source code and input files as follow:

| File/directory | Purpose |
| -------------- | ------- |
| ```LICENSE.md``` | a copy of the Gnu lesser general public license, version 3. |
| ```README.md``` | this file with basic pointers to more information |
| ```model/``` | contains the source code for core of the FV3 dyanmical core |
| ```model_nh/``` | contains the source code for non-hydrostatic extensions |
| ```driver/``` | contains drivers used by different models/modeling systems |
| ```tools/``` | contains source code of tools used within the core |
| ```GFDL_tools/``` | contains source code of tools specific to GFDL models |
| File/directory | Purpose |
| -------------- | ------- |
| ```LICENSE.md``` | a copy of the Gnu lesser general public license, version 3. |
| ```README.md``` | this file with basic pointers to more information |
| ```model/``` | contains the source code for core of the FV3 dyanmical core |
| ```driver/``` | contains drivers used by different models/modeling systems |
| ```tools/``` | contains source code of tools used within the core |
| ```GFDL_tools/``` | contains source code of tools specific to GFDL models |

# Disclaimer

Expand Down