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

Add capability to use RAP/HRRR MASSDEN as smoke IC/LBCs #923

Open
wants to merge 6 commits into
base: develop
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
3 changes: 2 additions & 1 deletion parm/varmap_tables/GSDphys_var_map.txt
Original file line number Diff line number Diff line change
Expand Up @@ -24,4 +24,5 @@ cnwat cnwat set_to_fill 0.0 S
icetk icetk set_to_fill 265.0 S
weasd weasd set_to_fill 0.0 S
snod snod set_to_fill 0.0 S
tprcp tprcp set_to_fill 0.0 S
tprcp tprcp set_to_fill 0.0 S
massden smoke set_to_fill 1E-12 T
84 changes: 77 additions & 7 deletions sorc/chgres_cube.fd/atm_input_data.F90
Original file line number Diff line number Diff line change
Expand Up @@ -1978,7 +1978,7 @@ subroutine read_input_atm_grib2_file(localpet)

integer, intent(in) :: localpet

integer, parameter :: ntrac_max=14
integer, parameter :: ntrac_max=15
integer, parameter :: max_levs=1000

character(len=300) :: the_file
Expand Down Expand Up @@ -2015,7 +2015,8 @@ subroutine read_input_atm_grib2_file(localpet)
real(esmf_kind_r8), allocatable :: rlevs(:)
real(esmf_kind_r4), allocatable :: dummy2d(:,:)
real(esmf_kind_r8), allocatable :: dummy3d(:,:,:), dummy2d_8(:,:),&
u_tmp_3d(:,:,:), v_tmp_3d(:,:,:)
u_tmp_3d(:,:,:), v_tmp_3d(:,:,:),&
dummy3d_pres(:,:,:)
real(esmf_kind_r8), pointer :: presptr(:,:,:), psptr(:,:),tptr(:,:,:), &
qptr(:,:,:), wptr(:,:,:), &
uptr(:,:,:), vptr(:,:,:)
Expand All @@ -2030,18 +2031,19 @@ subroutine read_input_atm_grib2_file(localpet)

tracers(:) = "NULL"

trac_names_oct10 = (/1, 1, 14, 1, 1, 1, 1, 6, 6, 1, 6, 13, 13, 2 /)
trac_names_oct11 = (/0, 22, 192, 23, 24, 25, 32, 1, 29, 100, 28, 193, 192, 2 /)
trac_names_oct10 = (/1, 1, 14, 1, 1, 1, 1, 6, 6, 1, 6, 13, 13, 2, 20 /)
trac_names_oct11 = (/0, 22, 192, 23, 24, 25, 32, 1, 29, 100, 28, 193, 192, 2, 0 /)


trac_names_vmap = (/"sphum ", "liq_wat ", "o3mr ", "ice_wat ", &
"rainwat ", "snowwat ", "graupel ", "cld_amt ", "ice_nc ", &
"rain_nc ", "water_nc", "liq_aero", "ice_aero", &
"sgs_tke "/)
"sgs_tke ", "massden "/)

tracers_default = (/"sphum ", "liq_wat ", "o3mr ", "ice_wat ", &
"rainwat ", "snowwat ", "graupel ", "cld_amt ", "ice_nc ", &
"rain_nc ", "water_nc", "liq_aero", "ice_aero", &
"sgs_tke "/)
"sgs_tke ", "smoke "/)

the_file = trim(data_dir_input_grid) // "/" // trim(grib2_file_input_grid)

Expand Down Expand Up @@ -2393,11 +2395,19 @@ subroutine read_input_atm_grib2_file(localpet)
allocate(dummy2d(i_input,j_input))
allocate(dummy2d_8(i_input,j_input))
allocate(dummy3d(i_input,j_input,lev_input))
if (trim(external_model) .eq. 'RAP' .or. & ! for smoke conversion
trim(external_model) .eq. 'HRRR' ) then
allocate(dummy3d_pres(i_input,j_input,lev_input))
endif
allocate(dum2d_1(i_input,j_input))
else
allocate(dummy2d(0,0))
allocate(dummy2d_8(0,0))
allocate(dummy3d(0,0,0))
if (trim(external_model) .eq. 'RAP' .or. & ! for smoke conversion
trim(external_model) .eq. 'HRRR' ) then
allocate(dummy3d_pres(0,0,0))
endif
allocate(dum2d_1(0,0))
endif

Expand Down Expand Up @@ -2457,12 +2467,61 @@ subroutine read_input_atm_grib2_file(localpet)
call get_var_cond(vname,this_miss_var_method=method, this_miss_var_value=value, &
this_field_var_name=tmpstr,loc=varnum)

if (n==1 .and. .not. hasspfh) then
if (n==1 .and. .not. hasspfh .or. &
( (trim(external_model) .eq. 'RAP' .or. & ! for smoke conversion
trim(external_model) .eq. 'HRRR' ) .and. &
tracers_input_vmap(n) == trac_names_vmap(15) )) then
print*,"- CALL FieldGather TEMPERATURE."
call ESMF_FieldGather(temp_input_grid,dummy3d,rootPet=0, tile=1, rc=rc)
if(ESMF_logFoundError(rcToCheck=rc,msg=ESMF_LOGERR_PASSTHRU,line=__LINE__,file=__FILE__)) &
call error_handler("IN FieldGet", rc)
endif

if ( (trim(external_model) .eq. 'RAP' .or. & ! for smoke conversion
trim(external_model) .eq. 'HRRR' ) .and. &
tracers_input_vmap(n) == trac_names_vmap(15)) then

if (localpet == 0) then

print*,"- READ PRESSURE FOR SMOKE CONVERSION."

jdisc = 0 ! search for discipline - meteorological products
j = 0 ! search at beginning of file.
jpdt = -9999 ! array of values in product definition template, set towildcard
jids = -9999 ! array of values in identification section, set towildcard
jgdt = -9999 ! array of values in grid definition template, set towildcard
jgdtn = -1 ! search for any grid definition number.
jpdtn = pdt_num ! Search for the product definition template number.
jpdt(1) = 3 ! Sect4/oct 10 - parameter category - mass
jpdt(2) = 0 ! Sect4/oct 11 - parameter number - pressure
jpdt(10) = octet_23 ! Sect4/oct 23 - type of level.
unpack=.true.

do vlev = 1, lev_input

jpdt(12) = nint(rlevs(vlev))
call getgb2(lugb, lugi, j, jdisc, jids, jpdtn, jpdt, jgdtn, jgdt, &
unpack, k, gfld, iret)
if (iret /= 0) then
call error_handler("READING IN PRESSURE AT LEVEL"//trim(slevs(vlev)),iret)
endif

dum2d_1 = reshape(gfld%fld, (/i_input,j_input/) )

dummy3d_pres(:,:,vlev) = dum2d_1

enddo

endif ! localpet == 0

endif ! read pressure for smoke conversion

if (tracers_input_vmap(n) == trac_names_vmap(15) .and. &
(trim(external_model) .ne. 'RAP' .and. & ! for smoke conversion
trim(external_model) .ne. 'HRRR' ) ) then
cycle ! Do not process smoke for non RAP/HRRR
endif


if (localpet == 0) then

Expand Down Expand Up @@ -2557,6 +2616,16 @@ subroutine read_input_atm_grib2_file(localpet)
end if
endif

! Convert smoke from mass density (RAP/HRRR = kg/m^3) to mixing ratio (ug/kg)
if ( tracers_input_vmap(n) == trac_names_vmap(15) ) then
do i = 1, i_input
do j = 1, j_input
dummy2d(i,j) = dummy2d(i,j) * 1.0d9 * &
(287.05 * dummy3d(i,j,vlev) / dummy3d_pres(i,j,vlev))
enddo
enddo
endif

dummy3d(:,:,vlev) = real(dummy2d,esmf_kind_r8)

enddo !vlev
Expand Down Expand Up @@ -2869,6 +2938,7 @@ subroutine read_input_atm_grib2_file(localpet)
endif

deallocate(dummy3d, dum2d_1)
if (allocated(dummy3d_pres)) deallocate(dummy3d_pres)

!---------------------------------------------------------------------------
! Convert from 2-d to 3-d component winds.
Expand Down
Loading