Skip to content

Commit 2ef07e3

Browse files
fix: vtec extrapolation calculation
- result of 4 corners was not been stored correctly, overwritting corner11 - domain of ZG was incorrect - units of ZG for need to be m in the vtec calculation fixes #851
1 parent 4eea462 commit 2ef07e3

File tree

1 file changed

+5
-5
lines changed

1 file changed

+5
-5
lines changed

models/tiegcm/model_mod.f90

Lines changed: 5 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -7,7 +7,6 @@
77
! - Nick Dietrich fix_mmr. When do to this?
88
! - model_time
99
! - get_state_meta_data 2D variables
10-
! - test vtec
1110

1211
module model_mod
1312

@@ -379,9 +378,9 @@ subroutine model_interpolate(state_handle, ens_size, location, iqty, obs_val, is
379378
if ( iqty == QTY_VERTICAL_TEC ) then ! extrapolate vtec
380379

381380
call extrapolate_vtec(state_handle, ens_size, lon_below, lat_below, val11)
382-
call extrapolate_vtec(state_handle, ens_size, lon_below, lat_above, val11)
383-
call extrapolate_vtec(state_handle, ens_size, lon_above, lat_below, val11)
384-
call extrapolate_vtec(state_handle, ens_size, lon_above, lat_above, val11)
381+
call extrapolate_vtec(state_handle, ens_size, lon_below, lat_above, val12)
382+
call extrapolate_vtec(state_handle, ens_size, lon_above, lat_below, val21)
383+
call extrapolate_vtec(state_handle, ens_size, lon_above, lat_above, val22)
385384
obs_val(:) = interpolate(ens_size, lon_fract, lat_fract, val11, val12, val21, val22)
386385
istatus(:) = 0
387386

@@ -1116,7 +1115,7 @@ subroutine extrapolate_vtec(state_handle, ens_size, lon_index, lat_index, vTEC)
11161115
enddo
11171116

11181117
! TI (midpoints)
1119-
var_id = get_varid_from_varname(domain_id(RESTART_DOM), 'TI')
1118+
var_id = get_varid_from_varname(domain_id(SECONDARY_DOM), 'TI')
11201119
do i = 1, nlev
11211120
idx = get_dart_vector_index(lon_index,lat_index, i, &
11221121
domain_id(RESTART_DOM), var_id)
@@ -1135,6 +1134,7 @@ subroutine extrapolate_vtec(state_handle, ens_size, lon_index, lat_index, vTEC)
11351134

11361135
earth_radiusm = earth_radius * 1000.0_r8 ! Convert earth_radius in km to m
11371136
NE = NE * 1.0e+6_r8 ! Convert NE in #/cm^3 to #/m^3
1137+
ZG = ZG * 1.0e-2_r8 ! Convert ZG in cm to m
11381138

11391139
! Gravity at the top layer
11401140
GRAVITYtop(:) = gravity * (earth_radiusm / (earth_radiusm + ZG(nilev,:))) ** 2

0 commit comments

Comments
 (0)