diff --git a/src/core_atmosphere/Registry.xml b/src/core_atmosphere/Registry.xml
index 087e8fb620..01be80e53b 100644
--- a/src/core_atmosphere/Registry.xml
+++ b/src/core_atmosphere/Registry.xml
@@ -490,6 +490,7 @@
+
#ifdef DO_PHYSICS
@@ -857,6 +858,7 @@
+
#endif
@@ -1011,38 +1013,81 @@
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
++
++
++
++
++
++
++
++
++
++
++
++
++
++
++
@@ -1079,6 +1124,7 @@
+
#endif
@@ -1099,6 +1145,15 @@
+
+
+
+
+
+
+
+
+
#ifdef DO_PHYSICS
diff --git a/src/core_atmosphere/diagnostics/Registry_convective.xml b/src/core_atmosphere/diagnostics/Registry_convective.xml
index 59ad1c2f58..f8d98f7e12 100644
--- a/src/core_atmosphere/diagnostics/Registry_convective.xml
+++ b/src/core_atmosphere/diagnostics/Registry_convective.xml
@@ -55,5 +55,28 @@
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
diff --git a/src/core_atmosphere/diagnostics/Registry_isobaric.xml b/src/core_atmosphere/diagnostics/Registry_isobaric.xml
index 853be6cde3..97b2be01bc 100644
--- a/src/core_atmosphere/diagnostics/Registry_isobaric.xml
+++ b/src/core_atmosphere/diagnostics/Registry_isobaric.xml
@@ -16,78 +16,102 @@
-
+
+
-
-
-
-
+
+
+
+
-
+
+
-
+
+
-
-
+
+
+
+
-
-
-
-
-
-
-
-
-
-
-
-
-
+
-
-
+
+
@@ -100,114 +124,219 @@
+
+
+
+
+
+
+
+
-
+
+
+
+
+
+
-
-
-
-
+
+
+
+
+
-
+
+
-
+
+
-
-
+
+
+
+
+
+
-
-
-
-
+
+
+
+
-
+
+
-
+
+
-
-
+
+
+
+
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
diff --git a/src/core_atmosphere/diagnostics/mpas_convective_diagnostics.F b/src/core_atmosphere/diagnostics/mpas_convective_diagnostics.F
index 163ee3774f..6d7b0bd4b1 100644
--- a/src/core_atmosphere/diagnostics/mpas_convective_diagnostics.F
+++ b/src/core_atmosphere/diagnostics/mpas_convective_diagnostics.F
@@ -14,6 +14,8 @@ module mpas_convective_diagnostics
type (MPAS_pool_type), pointer :: mesh
type (MPAS_pool_type), pointer :: state
type (MPAS_pool_type), pointer :: diag
+ type (MPAS_pool_type), pointer :: diag_physics
+ type (MPAS_pool_type), pointer :: sfc_input
type (MPAS_clock_type), pointer :: clock
@@ -65,6 +67,8 @@ subroutine convective_diagnostics_setup(all_pools, simulation_clock)
call mpas_pool_get_subpool(all_pools, 'mesh', mesh)
call mpas_pool_get_subpool(all_pools, 'state', state)
call mpas_pool_get_subpool(all_pools, 'diag', diag)
+ call mpas_pool_get_subpool(all_pools, 'diag_physics', diag_physics)
+ call mpas_pool_get_subpool(all_pools, 'sfc_input', sfc_input)
clock => simulation_clock
@@ -240,6 +244,17 @@ subroutine convective_diagnostics_compute()
real (kind=RKIND), dimension(:), pointer :: cin
real (kind=RKIND), dimension(:), pointer :: lcl
real (kind=RKIND), dimension(:), pointer :: lfc
+ real (kind=RKIND), dimension(:), pointer :: fzlev
+ real (kind=RKIND), dimension(:), pointer :: speed_surface
+ real (kind=RKIND), dimension(:), pointer :: heatidx
+ real (kind=RKIND), dimension(:), pointer :: ivt
+ real (kind=RKIND), dimension(:), pointer :: ztd
+ real (kind=RKIND), dimension(:), pointer :: latCell
+ real (kind=RKIND), dimension(:), pointer :: terrain
+ real (kind=RKIND), dimension(:), pointer :: gust_10m
+ real (kind=RKIND), dimension(:), pointer :: ust
+ real (kind=RKIND), dimension(:), pointer :: v10
+ real (kind=RKIND), dimension(:), pointer :: u10
real (kind=RKIND), dimension(:), pointer :: srh_0_1km
real (kind=RKIND), dimension(:), pointer :: srh_0_3km
real (kind=RKIND), dimension(:), pointer :: uzonal_surface
@@ -271,13 +286,16 @@ subroutine convective_diagnostics_compute()
real (kind=RKIND), dimension(:), allocatable :: dudz, dvdz, zp
real (kind=RKIND), dimension(:), allocatable :: zrel, srh
real (kind=RKIND), dimension(:), allocatable :: p_in, t_in, td_in
+ real (kind=RKIND) :: part1, part2
+ real (kind=RKIND) :: zf, zhd, zdiff, zwd_sfc, term1, term2
real (kind=RKIND), dimension(:,:), allocatable :: temperature, dewpoint
- real (kind=RKIND) :: evp
+ real (kind=RKIND) :: evp, tf2, tf, rh2, pressure_in
logical :: need_cape, need_cin, need_lcl, need_lfc, need_srh_01km, need_srh_03km, need_uzonal_sfc, need_uzonal_1km, &
- need_uzonal_6km, need_umerid_sfc, need_umerid_1km, need_umerid_6km, need_tsfc, need_tdsfc
+ need_uzonal_6km, need_umerid_sfc, need_umerid_1km, need_umerid_6km, need_tsfc, need_tdsfc, need_fzlev , &
+ need_speed_surface, need_heatidx, need_temp_wet, need_ivt, need_ztd
logical :: need_any_diags
need_any_diags = .false.
@@ -311,17 +329,46 @@ subroutine convective_diagnostics_compute()
need_tdsfc = MPAS_field_will_be_written('dewpoint_surface')
need_any_diags = need_any_diags .or. need_tdsfc
+ need_fzlev = MPAS_field_will_be_written('fzlev')
+ need_any_diags = need_any_diags .or. need_fzlev
+ need_speed_surface = MPAS_field_will_be_written('speed_surface')
+ need_any_diags = need_any_diags .or. need_speed_surface
+ need_heatidx = MPAS_field_will_be_written('heatidx')
+ need_any_diags = need_any_diags .or. need_heatidx
+ need_temp_wet = MPAS_field_will_be_written('temp_wet')
+ need_any_diags = need_any_diags .or. need_temp_wet
+ need_ivt = MPAS_field_will_be_written('ivt')
+ need_any_diags = need_any_diags .or. need_ivt
+ need_ztd = MPAS_field_will_be_written('ztd')
+ need_any_diags = need_any_diags .or. need_ztd
+
if (need_any_diags) then
call mpas_pool_get_dimension(mesh, 'nCells', nCells)
call mpas_pool_get_dimension(mesh, 'nCellsSolve', nCellsSolve)
call mpas_pool_get_dimension(mesh, 'nVertLevels', nVertLevels)
call mpas_pool_get_dimension(mesh, 'nVertLevelsP1', nVertLevelsP1)
+ call mpas_pool_get_array(mesh, 'latCell', latCell)
call mpas_pool_get_array(diag, 'cape', cape)
call mpas_pool_get_array(diag, 'cin', cin)
call mpas_pool_get_array(diag, 'lcl', lcl)
call mpas_pool_get_array(diag, 'lfc', lfc)
+ call mpas_pool_get_array(diag, 'fzlev', fzlev)
+ call mpas_pool_get_array(diag, 'speed_surface', speed_surface)
+ call mpas_pool_get_array(diag, 'heatidx', heatidx)
+ call mpas_pool_get_array(diag_physics, 't2m', t2m)
+ call mpas_pool_get_array(diag_physics, 'q2', q2)
+ call mpas_pool_get_array(diag_physics, 'ust', ust)
+ call mpas_pool_get_array(diag_physics, 'v10', v10)
+ call mpas_pool_get_array(diag_physics, 'u10', u10)
+ call mpas_pool_get_array(diag, 'surface_pressure', surface_pressure)
+ call mpas_pool_get_array(diag, 'relhum_2m', relhum_2m)
+ call mpas_pool_get_array(diag, 'temp_wet', temp_wet)
+ call mpas_pool_get_array(diag, 'ivt', ivt)
+ call mpas_pool_get_array(diag, 'ztd', ztd)
+ call mpas_pool_get_array(sfc_input, 'ter', terrain)
+ call mpas_pool_get_array(diag, 'gust_10m', gust_10m)
call mpas_pool_get_array(diag, 'srh_0_1km', srh_0_1km)
call mpas_pool_get_array(diag, 'srh_0_3km', srh_0_3km)
call mpas_pool_get_array(diag, 'uzonal_surface', uzonal_surface)
@@ -358,6 +405,18 @@ subroutine convective_diagnostics_compute()
allocate(t_in(nVertLevels))
allocate(td_in(nVertLevels))
+ ctt = 0._RKIND
+ ctpres = 0._RKIND
+ ivt = 0._RKIND
+ ztd = 0._RKIND
+ cape = 0._RKIND
+ cin = 0._RKIND
+ gust_10m = 0._RKIND
+ fzlev = 0._RKIND
+ relhum_2m = 0._RKIND
+ heatidx = 0._RKIND
+ temp_wet = 0._RKIND
+
do iCell = 1,nCells
do k = 1,nVertLevels
temperature(k,iCell) = (theta_m(k,iCell)/(1._RKIND+rvord*scalars(index_qv,k,iCell)))*exner(k,iCell)
@@ -384,6 +443,7 @@ subroutine convective_diagnostics_compute()
umeridional_surface(iCell) = umeridional(1,iCell)
temperature_surface(iCell) = temperature(1,iCell)
dewpoint_surface(iCell) = dewpoint(1,iCell)
+ speed_surface(iCell) = sqrt(uzonal_surface(iCell)**2. + umeridional_surface(iCell)**2.)
if (need_uzonal_1km) then
uzonal_1km(iCell) = column_height_value(uzonal(1:nVertLevels,iCell), zp, 1000.0_RKIND, nVertLevels)
end if
@@ -471,6 +531,111 @@ subroutine convective_diagnostics_compute()
end do
end if
+ if (need_ivt .or. need_ztd) then
+ do iCell=1, nCellsSolve
+ p_in(1:nVertLevels) = (pressure_p(1:nVertLevels,iCell) + pressure_base(1:nVertLevels,iCell))
+ t_in(1:nVertLevels) = temperature(1:nVertLevels,iCell)
+ qvapor_in(1:nVertLevels) = scalars(index_qv,1:nVertLevels,iCell)
+ !for IVT
+ part1=0._RKIND
+ part2=0._RKIND
+ do k = 1,nVertLevels-1
+ if(p_in(k+1) >= 25000.0_RKIND)THEN
+ part1 = part1 + uzonal(k,iCell)*scalars(index_qv,k,iCell) * (p_in(k+1) - p_in(k))
+ part2 = part2 + umeridional(k,iCell)*scalars(index_qv,k,iCell) * (p_in(k+1) - p_in(k))
+ ivt(iCell) = sqrt((part1/9.806)**2+(part2/9.806)**2)
+ endif
+ enddo
+
+ !for ZTD
+ part1=0._RKIND
+ part2=0._RKIND
+ zwd_sfc = 0._RKIND
+! Saastamoinen
+ zf = 1. - 2.66e-3*cos(2.*latCell(iCell))-2.8e-7*terrain(iCell)
+ zhd = 2.2768e-5*surface_pressure(iCell)/zf
+
+ do k = 1,nVertLevels-1
+ zdiff = (height(k+1,iCell) - height(k,iCell))/0.622_RKIND
+ part1 = p_in(k)*qvapor_in(k)/t_in(k)
+ term1 = part1 * zdiff * 2.21e-7
+ term2 = part1 * zdiff * 3.739e-3 / t_in(k)
+ zwd_sfc = zwd_sfc + term1+term2
+ enddo
+ ztd(iCell) = (zwd_sfc + zhd)*100.0_RKIND
+ gust_10m(iCell) = sqrt(u10(iCell)*u10(iCell)+v10(iCell)*v10(iCell))+ (3.*2.4*ust(iCell))
+ end do ! cells
+ end if
+
+! from WRF
+ IF (need_fzlev) then
+ fzlev = -999._RKIND
+ do iCell = 1,nCells
+ do k = 1,nVertLevels
+ IF (temperature(k,icell) .lt. 273.15_RKIND) THEN
+ ! Freezing level calculation
+ ! --------------------------
+ IF (fzlev (icell) .eq. -999._RKIND) THEN ! At freezing level
+ IF (k .ne. 1) THEN ! If not at surface, interpolate.
+ fzlev (icell) = height (k-1,icell) + &
+ ((273.15_RKIND - temperature (k-1,icell)) &
+ /(temperature (k,icell) - temperature (k-1,icell))) &
+ *(height (k,icell) - height (k-1,icell))
+ ELSE ! If at surface, use first level.
+ fzlev(icell) = height (k,icell)
+ ENDIF
+ ENDIF
+ ENDIF
+ enddo
+ enddo
+ ENDIF
+
+ DO iCell = 1, nCells
+ IF(t2m(iCell) > 273.15_RKIND )THEN
+ relhum_2m(iCell) = 0.263*surface_pressure(iCell)*q2(iCell)*(exp(17.625*(t2m(iCell)-273.15_RKIND)/(t2m(iCell)-30.1)))**(-1.)
+ ELSE IF( t2m(iCell) .le. 273.15_RKIND) THEN
+! https://journals.ametsoc.org/view/journals/apme/57/6/jamc-d-17-0334.1.xml
+ relhum_2m(iCell) = 0.263*surface_pressure(iCell)*q2(iCell)*(exp(22.587*(t2m(iCell)-273.15_RKIND)/(t2m(iCell)+0.71)))**(-1.)
+ ENDIF
+ ENDDO
+ where(relhum_2m > 100.) relhum_2m = 100._RKIND
+
+ IF (need_heatidx) then
+ do iCell = 1,nCells
+ IF ( t2m(icell) > 294.26111) THEN
+ tF = ( (t2m(icell) - 273.15_RKIND) * (REAL (9)/REAL (5)) ) + REAL ( 32 )
+ tF2 = tF ** 2
+ rh2 = relhum_2m(iCell)/100. ** 2
+ heatidx(icell) = -42.379 &
+ + ( 2.04901523 * tF ) &
+ + ( 10.14333127 * relhum_2m(iCell)/100. ) &
+ - ( 0.22475541 * tF * relhum_2m(iCell)/100. ) &
+ - ( 6.83783E-03 * tF2 ) &
+ - ( 5.481717E-02 * rh2 ) &
+ + ( 1.22874E-03 * tF2 * relhum_2m(iCell)/100. ) &
+ + ( 8.5282E-04 * tF * rh2 ) &
+ - ( 1.99E-06 * tF2 * rh2 )
+ heatidx(icell) = ((heatidx(icell) - REAL (32)) * (REAL (5)/REAL (9))) + 273.15_RKIND
+ ELSE
+ heatidx(icell) = t2m(icell)
+ END IF
+ enddo
+ ENDIF
+
+ IF(need_temp_wet)THEN
+! Stull 2011: doi:10.1175/JAMC-D-11-0143.1
+ DO iCell =1, nCells
+ IF(t2m(iCell) > 253.15 .AND. relhum_2m(iCell) .ge. 5. .AND. relhum_2m(iCell) .le. 99.)THEN
+ temp_wet(iCell) = (t2m(iCell)-273.15_RKIND)*atan(0.151977*sqrt(relhum_2m(iCell)+8.313659)) &
+ + atan(t2m(iCell)-273.15_RKIND+relhum_2m(iCell))&
+ - atan(relhum_2m(iCell)-1.676331)+0.00391838*((relhum_2m(iCell))**(1.5))*atan(0.023101*relhum_2m(iCell)) &
+ - 4.686035 +273.15_RKIND
+ ELSE
+ temp_wet(iCell) = -9999._RKIND
+ ENDIF
+ ENDDO
+ ENDIF
+
deallocate(temperature)
deallocate(dewpoint)
diff --git a/src/core_atmosphere/diagnostics/mpas_isobaric_diagnostics.F b/src/core_atmosphere/diagnostics/mpas_isobaric_diagnostics.F
index e52c71b125..ef887a2be5 100644
--- a/src/core_atmosphere/diagnostics/mpas_isobaric_diagnostics.F
+++ b/src/core_atmosphere/diagnostics/mpas_isobaric_diagnostics.F
@@ -35,7 +35,14 @@ module mpas_isobaric_diagnostics
need_w_50, need_w_100, need_w_200, need_w_250, need_w_500, need_w_700, need_w_850, need_w_925, &
need_vorticity_50, need_vorticity_100, need_vorticity_200, need_vorticity_250, need_vorticity_500, need_vorticity_700, need_vorticity_850, need_vorticity_925, &
need_t_isobaric, need_z_isobaric, need_meanT_500_300
- logical :: need_temp, need_relhum, need_dewpoint, need_w, need_uzonal, need_umeridional, need_vorticity, need_height
+ LOGICAL :: need_qv_200, need_qv_250, need_qv_500, need_qv_700, need_qv_850, need_qv_925, need_qv, need_temp_300, need_umeridional_300, &
+ need_uzonal_300 , need_qv_300, need_temp_400, need_umeridional_400,need_uzonal_400 , need_qv_400, &
+ need_temp_600, need_umeridional_600,need_uzonal_600 , need_qv_600, need_height_300, need_height_400, need_height_600, &
+ need_temp_1000, need_umeridional_1000,need_uzonal_1000 , need_qv_1000, need_height_1000, &
+ need_temp_20, need_umeridional_20,need_uzonal_20 , need_qv_20, need_height_20, &
+ need_temp_800, need_umeridional_800,need_uzonal_800 , need_qv_800, need_height_800, &
+ need_temp_950, need_umeridional_950,need_uzonal_950 , need_qv_950, need_height_950
+ logical :: need_temp, need_relhum, need_dewpoint, need_w, need_uzonal, need_umeridional, need_vorticity, need_height, need_qv
contains
@@ -100,6 +107,7 @@ subroutine isobaric_diagnostics_compute()
need_umeridional = .false.
need_vorticity = .false.
need_height = .false.
+ need_qv = .false.
need_mslp = MPAS_field_will_be_written('mslp')
need_any_diags = need_any_diags .or. need_mslp
@@ -302,6 +310,135 @@ subroutine isobaric_diagnostics_compute()
need_meanT_500_300 = MPAS_field_will_be_written('meanT_500_300')
need_any_diags = need_any_diags .or. need_meanT_500_300
+! Additional pressure levels
+
+ need_qv_200 = MPAS_field_will_be_written('qv_200hPa')
+ need_qv = need_qv .or. need_qv_200
+ need_any_diags = need_any_diags .or. need_qv_200
+ need_qv_250 = MPAS_field_will_be_written('qv_250hPa')
+ need_qv = need_qv .or. need_qv_250
+ need_any_diags = need_any_diags .or. need_qv_250
+ need_qv_300 = MPAS_field_will_be_written('qv_300hPa')
+ need_qv = need_qv .or. need_qv_300
+ need_any_diags = need_any_diags .or. need_qv_300
+ need_qv_500 = MPAS_field_will_be_written('qv_500hPa')
+ need_qv = need_qv .or. need_qv_500
+ need_any_diags = need_any_diags .or. need_qv_500
+ need_qv_700 = MPAS_field_will_be_written('qv_700hPa')
+ need_qv = need_qv .or. need_qv_700
+ need_any_diags = need_any_diags .or. need_qv_700
+ need_qv_850 = MPAS_field_will_be_written('qv_850hPa')
+ need_qv = need_qv .or. need_qv_850
+ need_any_diags = need_any_diags .or. need_qv_850
+ need_qv_925 = MPAS_field_will_be_written('qv_925hPa')
+ need_qv = need_qv .or. need_qv_925
+ need_any_diags = need_any_diags .or. need_qv_925
+
+ need_temp_20 = MPAS_field_will_be_written('temperature_20hPa')
+ need_temp = need_temp .or. need_temp_20
+ need_any_diags = need_any_diags .or. need_temp_20
+ need_umeridional_20 = MPAS_field_will_be_written('umeridional_20hPa')
+ need_umeridional = need_umeridional .or. need_umeridional_20
+ need_any_diags = need_any_diags .or. need_umeridional_20
+ need_uzonal_20 = MPAS_field_will_be_written('uzonal_20hPa')
+ need_uzonal = need_uzonal .or. need_uzonal_20
+ need_any_diags = need_any_diags .or. need_uzonal_20
+ need_qv_20 = MPAS_field_will_be_written('qv_20hPa')
+ need_qv = need_qv .or. need_qv_20
+ need_any_diags = need_any_diags .or. need_qv_20
+
+ need_temp_300 = MPAS_field_will_be_written('temperature_300hPa')
+ need_temp = need_temp .or. need_temp_300
+ need_any_diags = need_any_diags .or. need_temp_300
+ need_umeridional_300 = MPAS_field_will_be_written('umeridional_300hPa')
+ need_umeridional = need_umeridional .or. need_umeridional_300
+ need_any_diags = need_any_diags .or. need_umeridional_300
+ need_uzonal_300 = MPAS_field_will_be_written('uzonal_300hPa')
+ need_uzonal = need_uzonal .or. need_uzonal_300
+ need_any_diags = need_any_diags .or. need_uzonal_300
+
+ need_temp_400 = MPAS_field_will_be_written('temperature_400hPa')
+ need_temp = need_temp .or. need_temp_400
+ need_any_diags = need_any_diags .or. need_temp_400
+ need_umeridional_400 = MPAS_field_will_be_written('umeridional_400hPa')
+ need_umeridional = need_umeridional .or. need_umeridional_400
+ need_any_diags = need_any_diags .or. need_umeridional_400
+ need_uzonal_400 = MPAS_field_will_be_written('uzonal_400hPa')
+ need_uzonal = need_uzonal .or. need_uzonal_400
+ need_any_diags = need_any_diags .or. need_uzonal_400
+ need_qv_400 = MPAS_field_will_be_written('qv_400hPa')
+ need_qv = need_qv .or. need_qv_400
+ need_any_diags = need_any_diags .or. need_qv_400
+
+ need_temp_600 = MPAS_field_will_be_written('temperature_600hPa')
+ need_temp = need_temp .or. need_temp_600
+ need_any_diags = need_any_diags .or. need_temp_600
+ need_umeridional_600 = MPAS_field_will_be_written('umeridional_600hPa')
+ need_umeridional = need_umeridional .or. need_umeridional_600
+ need_any_diags = need_any_diags .or. need_umeridional_600
+ need_uzonal_600 = MPAS_field_will_be_written('uzonal_600hPa')
+ need_uzonal = need_uzonal .or. need_uzonal_600
+ need_any_diags = need_any_diags .or. need_uzonal_600
+ need_qv_600 = MPAS_field_will_be_written('qv_600hPa')
+ need_qv = need_qv .or. need_qv_600
+ need_any_diags = need_any_diags .or. need_qv_600
+
+ need_temp_800 = MPAS_field_will_be_written('temperature_800hPa')
+ need_temp = need_temp .or. need_temp_800
+ need_any_diags = need_any_diags .or. need_temp_800
+ need_umeridional_800 = MPAS_field_will_be_written('umeridional_800hPa')
+ need_umeridional = need_umeridional .or. need_umeridional_800
+ need_any_diags = need_any_diags .or. need_umeridional_800
+ need_uzonal_800 = MPAS_field_will_be_written('uzonal_800hPa')
+ need_uzonal = need_uzonal .or. need_uzonal_800
+ need_any_diags = need_any_diags .or. need_uzonal_800
+ need_qv_800 = MPAS_field_will_be_written('qv_800hPa')
+ need_qv = need_qv .or. need_qv_800
+ need_any_diags = need_any_diags .or. need_qv_800
+
+ need_temp_950 = MPAS_field_will_be_written('temperature_950hPa')
+ need_temp = need_temp .or. need_temp_950
+ need_any_diags = need_any_diags .or. need_temp_950
+ need_umeridional_950 = MPAS_field_will_be_written('umeridional_950hPa')
+ need_umeridional = need_umeridional .or. need_umeridional_950
+ need_any_diags = need_any_diags .or. need_umeridional_950
+ need_uzonal_950 = MPAS_field_will_be_written('uzonal_950hPa')
+ need_uzonal = need_uzonal .or. need_uzonal_950
+ need_any_diags = need_any_diags .or. need_uzonal_950
+ need_qv_950 = MPAS_field_will_be_written('qv_950hPa')
+ need_qv = need_qv .or. need_qv_950
+ need_any_diags = need_any_diags .or. need_qv_950
+
+ need_temp_1000 = MPAS_field_will_be_written('temperature_1000hPa')
+ need_temp = need_temp .or. need_temp_1000
+ need_any_diags = need_any_diags .or. need_temp_1000
+ need_umeridional_1000 = MPAS_field_will_be_written('umeridional_1000hPa')
+ need_umeridional = need_umeridional .or. need_umeridional_1000
+ need_any_diags = need_any_diags .or. need_umeridional_1000
+ need_uzonal_1000 = MPAS_field_will_be_written('uzonal_1000hPa')
+ need_uzonal = need_uzonal .or. need_uzonal_1000
+ need_any_diags = need_any_diags .or. need_uzonal_1000
+ need_qv_1000 = MPAS_field_will_be_written('qv_1000hPa')
+ need_qv = need_qv .or. need_qv_1000
+ need_any_diags = need_any_diags .or. need_qv_1000
+
+
+ need_height_100 = MPAS_field_will_be_written('height_100hPa')
+ need_height = need_height .or. need_height_100
+ need_any_diags = need_any_diags .or. need_height_100
+ need_height_300 = MPAS_field_will_be_written('height_300hPa')
+ need_height = need_height .or. need_height_300
+ need_any_diags = need_any_diags .or. need_height_300
+ need_height_400 = MPAS_field_will_be_written('height_400hPa')
+ need_height = need_height .or. need_height_400
+ need_any_diags = need_any_diags .or. need_height_400
+ need_height_600 = MPAS_field_will_be_written('height_600hPa')
+ need_height = need_height .or. need_height_600
+ need_any_diags = need_any_diags .or. need_height_600
+ need_height_1000 = MPAS_field_will_be_written('height_1000hPa')
+ need_height = need_height .or. need_height_1000
+ need_any_diags = need_any_diags .or. need_height_1000
+
if (need_any_diags) then
call interp_diagnostics(mesh, state, 1, diag)
end if
@@ -344,15 +481,22 @@ subroutine interp_diagnostics(mesh, state, time_lev, diag)
real (kind=RKIND), dimension(:,:), pointer :: t_isobaric
real (kind=RKIND), dimension(:,:), pointer :: z_isobaric
real (kind=RKIND), dimension(:), pointer :: meanT_500_300
-
+
+ real (kind=RKIND), dimension(:), pointer :: temperature_20hPa
real (kind=RKIND), dimension(:), pointer :: temperature_50hPa
real (kind=RKIND), dimension(:), pointer :: temperature_100hPa
real (kind=RKIND), dimension(:), pointer :: temperature_200hPa
real (kind=RKIND), dimension(:), pointer :: temperature_250hPa
+ real (kind=RKIND), dimension(:), pointer :: temperature_300hPa
+ real (kind=RKIND), dimension(:), pointer :: temperature_400hPa
real (kind=RKIND), dimension(:), pointer :: temperature_500hPa
+ real (kind=RKIND), dimension(:), pointer :: temperature_600hPa
real (kind=RKIND), dimension(:), pointer :: temperature_700hPa
+ real (kind=RKIND), dimension(:), pointer :: temperature_800hPa
real (kind=RKIND), dimension(:), pointer :: temperature_850hPa
real (kind=RKIND), dimension(:), pointer :: temperature_925hPa
+ real (kind=RKIND), dimension(:), pointer :: temperature_950hPa
+ real (kind=RKIND), dimension(:), pointer :: temperature_1000hPa
real (kind=RKIND), dimension(:), pointer :: relhum_50hPa
real (kind=RKIND), dimension(:), pointer :: relhum_100hPa
@@ -371,34 +515,55 @@ subroutine interp_diagnostics(mesh, state, time_lev, diag)
real (kind=RKIND), dimension(:), pointer :: dewpoint_700hPa
real (kind=RKIND), dimension(:), pointer :: dewpoint_850hPa
real (kind=RKIND), dimension(:), pointer :: dewpoint_925hPa
-
+
+ real (kind=RKIND), dimension(:), pointer :: uzonal_20hPa
real (kind=RKIND), dimension(:), pointer :: uzonal_50hPa
real (kind=RKIND), dimension(:), pointer :: uzonal_100hPa
real (kind=RKIND), dimension(:), pointer :: uzonal_200hPa
real (kind=RKIND), dimension(:), pointer :: uzonal_250hPa
+ real (kind=RKIND), dimension(:), pointer :: uzonal_300hPa
+ real (kind=RKIND), dimension(:), pointer :: uzonal_400hPa
real (kind=RKIND), dimension(:), pointer :: uzonal_500hPa
+ real (kind=RKIND), dimension(:), pointer :: uzonal_600hPa
real (kind=RKIND), dimension(:), pointer :: uzonal_700hPa
+ real (kind=RKIND), dimension(:), pointer :: uzonal_800hPa
real (kind=RKIND), dimension(:), pointer :: uzonal_850hPa
real (kind=RKIND), dimension(:), pointer :: uzonal_925hPa
+ real (kind=RKIND), dimension(:), pointer :: uzonal_950hPa
+ real (kind=RKIND), dimension(:), pointer :: uzonal_1000hPa
+ real (kind=RKIND), dimension(:), pointer :: umeridional_20hPa
real (kind=RKIND), dimension(:), pointer :: umeridional_50hPa
real (kind=RKIND), dimension(:), pointer :: umeridional_100hPa
real (kind=RKIND), dimension(:), pointer :: umeridional_200hPa
real (kind=RKIND), dimension(:), pointer :: umeridional_250hPa
+ real (kind=RKIND), dimension(:), pointer :: umeridional_300hPa
+ real (kind=RKIND), dimension(:), pointer :: umeridional_400hPa
real (kind=RKIND), dimension(:), pointer :: umeridional_500hPa
+ real (kind=RKIND), dimension(:), pointer :: umeridional_600hPa
real (kind=RKIND), dimension(:), pointer :: umeridional_700hPa
+ real (kind=RKIND), dimension(:), pointer :: umeridional_800hPa
real (kind=RKIND), dimension(:), pointer :: umeridional_850hPa
real (kind=RKIND), dimension(:), pointer :: umeridional_925hPa
-
+ real (kind=RKIND), dimension(:), pointer :: umeridional_950hPa
+ real (kind=RKIND), dimension(:), pointer :: umeridional_1000hPa
+
+ real (kind=RKIND), dimension(:), pointer :: height_20hPa
real (kind=RKIND), dimension(:), pointer :: height_50hPa
real (kind=RKIND), dimension(:), pointer :: height_100hPa
real (kind=RKIND), dimension(:), pointer :: height_200hPa
real (kind=RKIND), dimension(:), pointer :: height_250hPa
+ real (kind=RKIND), dimension(:), pointer :: height_300hPa
+ real (kind=RKIND), dimension(:), pointer :: height_400hPa
real (kind=RKIND), dimension(:), pointer :: height_500hPa
+ real (kind=RKIND), dimension(:), pointer :: height_600hPa
real (kind=RKIND), dimension(:), pointer :: height_700hPa
+ real (kind=RKIND), dimension(:), pointer :: height_800hPa
real (kind=RKIND), dimension(:), pointer :: height_850hPa
real (kind=RKIND), dimension(:), pointer :: height_925hPa
-
+ real (kind=RKIND), dimension(:), pointer :: height_950hPa
+ real (kind=RKIND), dimension(:), pointer :: height_1000hPa
+
real (kind=RKIND), dimension(:), pointer :: w_50hPa
real (kind=RKIND), dimension(:), pointer :: w_100hPa
real (kind=RKIND), dimension(:), pointer :: w_200hPa
@@ -416,6 +581,22 @@ subroutine interp_diagnostics(mesh, state, time_lev, diag)
real (kind=RKIND), dimension(:), pointer :: vorticity_700hPa
real (kind=RKIND), dimension(:), pointer :: vorticity_850hPa
real (kind=RKIND), dimension(:), pointer :: vorticity_925hPa
+
+ real (kind=RKIND), dimension(:), pointer :: qv_20hPa
+ real (kind=RKIND), dimension(:), pointer :: qv_50hPa
+ real (kind=RKIND), dimension(:), pointer :: qv_100hPa
+ real (kind=RKIND), dimension(:), pointer :: qv_200hPa
+ real (kind=RKIND), dimension(:), pointer :: qv_250hPa
+ real (kind=RKIND), dimension(:), pointer :: qv_300hPa
+ real (kind=RKIND), dimension(:), pointer :: qv_400hPa
+ real (kind=RKIND), dimension(:), pointer :: qv_500hPa
+ real (kind=RKIND), dimension(:), pointer :: qv_600hPa
+ real (kind=RKIND), dimension(:), pointer :: qv_700hPa
+ real (kind=RKIND), dimension(:), pointer :: qv_800hPa
+ real (kind=RKIND), dimension(:), pointer :: qv_850hPa
+ real (kind=RKIND), dimension(:), pointer :: qv_925hPa
+ real (kind=RKIND), dimension(:), pointer :: qv_950hPa
+ real (kind=RKIND), dimension(:), pointer :: qv_1000hPa
real (kind=RKIND) :: evp
@@ -474,16 +655,23 @@ subroutine interp_diagnostics(mesh, state, time_lev, diag)
call mpas_pool_get_array(diag, 't_isobaric', t_isobaric)
call mpas_pool_get_array(diag, 'z_isobaric', z_isobaric)
call mpas_pool_get_array(diag, 'meanT_500_300', meanT_500_300)
-
+
+ call mpas_pool_get_array(diag, 'temperature_20hPa', temperature_20hPa)
call mpas_pool_get_array(diag, 'temperature_50hPa', temperature_50hPa)
call mpas_pool_get_array(diag, 'temperature_100hPa', temperature_100hPa)
call mpas_pool_get_array(diag, 'temperature_200hPa', temperature_200hPa)
call mpas_pool_get_array(diag, 'temperature_250hPa', temperature_250hPa)
+ call mpas_pool_get_array(diag, 'temperature_300hPa', temperature_300hPa)
+ call mpas_pool_get_array(diag, 'temperature_400hPa', temperature_400hPa)
call mpas_pool_get_array(diag, 'temperature_500hPa', temperature_500hPa)
+ call mpas_pool_get_array(diag, 'temperature_600hPa', temperature_600hPa)
call mpas_pool_get_array(diag, 'temperature_700hPa', temperature_700hPa)
+ call mpas_pool_get_array(diag, 'temperature_800hPa', temperature_3800hPa)
call mpas_pool_get_array(diag, 'temperature_850hPa', temperature_850hPa)
call mpas_pool_get_array(diag, 'temperature_925hPa', temperature_925hPa)
-
+ call mpas_pool_get_array(diag, 'temperature_950hPa', temperature_950hPa)
+ call mpas_pool_get_array(diag, 'temperature_1000hPa', temperature_1000hPa)
+
call mpas_pool_get_array(diag, 'relhum_50hPa', relhum_50hPa)
call mpas_pool_get_array(diag, 'relhum_100hPa', relhum_100hPa)
call mpas_pool_get_array(diag, 'relhum_200hPa', relhum_200hPa)
@@ -501,33 +689,54 @@ subroutine interp_diagnostics(mesh, state, time_lev, diag)
call mpas_pool_get_array(diag, 'dewpoint_700hPa', dewpoint_700hPa)
call mpas_pool_get_array(diag, 'dewpoint_850hPa', dewpoint_850hPa)
call mpas_pool_get_array(diag, 'dewpoint_925hPa', dewpoint_925hPa)
-
+
+ call mpas_pool_get_array(diag, 'uzonal_20hPa', uzonal_20hPa)
call mpas_pool_get_array(diag, 'uzonal_50hPa', uzonal_50hPa)
call mpas_pool_get_array(diag, 'uzonal_100hPa', uzonal_100hPa)
call mpas_pool_get_array(diag, 'uzonal_200hPa', uzonal_200hPa)
call mpas_pool_get_array(diag, 'uzonal_250hPa', uzonal_250hPa)
+ call mpas_pool_get_array(diag, 'uzonal_300hPa', uzonal_300hPa)
+ call mpas_pool_get_array(diag, 'uzonal_400hPa', uzonal_400hPa)
call mpas_pool_get_array(diag, 'uzonal_500hPa', uzonal_500hPa)
+ call mpas_pool_get_array(diag, 'uzonal_600hPa', uzonal_600hPa)
call mpas_pool_get_array(diag, 'uzonal_700hPa', uzonal_700hPa)
+ call mpas_pool_get_array(diag, 'uzonal_800hPa', uzonal_800hPa)
call mpas_pool_get_array(diag, 'uzonal_850hPa', uzonal_850hPa)
call mpas_pool_get_array(diag, 'uzonal_925hPa', uzonal_925hPa)
-
+ call mpas_pool_get_array(diag, 'uzonal_950hPa', uzonal_950hPa)
+ call mpas_pool_get_array(diag, 'uzonal_1000hPa', uzonal_1000hPa)
+
+ call mpas_pool_get_array(diag, 'umeridional_20hPa', umeridional_20hPa)
call mpas_pool_get_array(diag, 'umeridional_50hPa', umeridional_50hPa)
call mpas_pool_get_array(diag, 'umeridional_100hPa', umeridional_100hPa)
call mpas_pool_get_array(diag, 'umeridional_200hPa', umeridional_200hPa)
call mpas_pool_get_array(diag, 'umeridional_250hPa', umeridional_250hPa)
+ call mpas_pool_get_array(diag, 'umeridional_300hPa', umeridional_300hPa)
+ call mpas_pool_get_array(diag, 'umeridional_400hPa', umeridional_400hPa)
call mpas_pool_get_array(diag, 'umeridional_500hPa', umeridional_500hPa)
+ call mpas_pool_get_array(diag, 'umeridional_600hPa', umeridional_600hPa)
call mpas_pool_get_array(diag, 'umeridional_700hPa', umeridional_700hPa)
+ call mpas_pool_get_array(diag, 'umeridional_800hPa', umeridional_800hPa)
call mpas_pool_get_array(diag, 'umeridional_850hPa', umeridional_850hPa)
call mpas_pool_get_array(diag, 'umeridional_925hPa', umeridional_925hPa)
-
+ call mpas_pool_get_array(diag, 'umeridional_950hPa', umeridional_950hPa)
+ call mpas_pool_get_array(diag, 'umeridional_1000hPa', umeridional_1000hPa)
+
+ call mpas_pool_get_array(diag, 'height_20hPa', height_20hPa)
call mpas_pool_get_array(diag, 'height_50hPa', height_50hPa)
call mpas_pool_get_array(diag, 'height_100hPa', height_100hPa)
call mpas_pool_get_array(diag, 'height_200hPa', height_200hPa)
call mpas_pool_get_array(diag, 'height_250hPa', height_250hPa)
+ call mpas_pool_get_array(diag, 'height_300hPa', height_300hPa)
+ call mpas_pool_get_array(diag, 'height_400hPa', height_400hPa)
call mpas_pool_get_array(diag, 'height_500hPa', height_500hPa)
+ call mpas_pool_get_array(diag, 'height_600hPa', height_600hPa)
call mpas_pool_get_array(diag, 'height_700hPa', height_700hPa)
+ call mpas_pool_get_array(diag, 'height_800hPa', height_800hPa)
call mpas_pool_get_array(diag, 'height_850hPa', height_850hPa)
call mpas_pool_get_array(diag, 'height_925hPa', height_925hPa)
+ call mpas_pool_get_array(diag, 'height_950hPa', height_950hPa)
+ call mpas_pool_get_array(diag, 'height_1000hPa', height_1000hPa)
call mpas_pool_get_array(diag, 'w_50hPa', w_50hPa)
call mpas_pool_get_array(diag, 'w_100hPa', w_100hPa)
@@ -658,18 +867,25 @@ subroutine interp_diagnostics(mesh, state, time_lev, diag)
end if
!interpolation to fixed pressure levels for fields located at cells centers and at mass points:
- nIntP = 8
+ nIntP = 15
if(.not.allocated(field_interp)) allocate(field_interp(nCells,nIntP) )
if(.not.allocated(press_interp)) allocate(press_interp(nCells,nIntP) )
do iCell = 1, nCells
- press_interp(iCell,1) = 50.0_RKIND
- press_interp(iCell,2) = 100.0_RKIND
- press_interp(iCell,3) = 200.0_RKIND
- press_interp(iCell,4) = 250.0_RKIND
- press_interp(iCell,5) = 500.0_RKIND
- press_interp(iCell,6) = 700.0_RKIND
- press_interp(iCell,7) = 850.0_RKIND
- press_interp(iCell,8) = 925.0_RKIND
+ press_interp(iCell,1) = 20.0_RKIND
+ press_interp(iCell,2) = 50.0_RKIND
+ press_interp(iCell,3) = 100.0_RKIND
+ press_interp(iCell,4) = 200.0_RKIND
+ press_interp(iCell,5) = 250.0_RKIND
+ press_interp(iCell,6) = 300.0_RKIND
+ press_interp(iCell,7) = 400.0_RKIND
+ press_interp(iCell,8) = 500.0_RKIND
+ press_interp(iCell,9) = 600.0_RKIND
+ press_interp(iCell,10) = 700.0_RKIND
+ press_interp(iCell,11) = 800.0_RKIND
+ press_interp(iCell,12) = 850.0_RKIND
+ press_interp(iCell,13) = 925.0_RKIND
+ press_interp(iCell,14) = 950.0_RKIND
+ press_interp(iCell,15) = 1000.0_RKIND
enddo
if(.not.allocated(press_in)) allocate(press_in(nCells,nVertLevels))
@@ -691,14 +907,21 @@ subroutine interp_diagnostics(mesh, state, time_lev, diag)
enddo
enddo
call interp_tofixed_pressure(nCells,nVertLevels,nIntP,press_in,field_in,press_interp,field_interp)
- temperature_50hPa(1:nCells) = field_interp(1:nCells,1)
- temperature_100hPa(1:nCells) = field_interp(1:nCells,2)
- temperature_200hPa(1:nCells) = field_interp(1:nCells,3)
- temperature_250hPa(1:nCells) = field_interp(1:nCells,4)
- temperature_500hPa(1:nCells) = field_interp(1:nCells,5)
- temperature_700hPa(1:nCells) = field_interp(1:nCells,6)
- temperature_850hPa(1:nCells) = field_interp(1:nCells,7)
- temperature_925hPa(1:nCells) = field_interp(1:nCells,8)
+ temperature_20hPa(1:nCells) = field_interp(1:nCells,1)
+ temperature_50hPa(1:nCells) = field_interp(1:nCells,2)
+ temperature_100hPa(1:nCells) = field_interp(1:nCells,3)
+ temperature_200hPa(1:nCells) = field_interp(1:nCells,4)
+ temperature_250hPa(1:nCells) = field_interp(1:nCells,5)
+ temperature_300hPa(1:nCells) = field_interp(1:nCells,6)
+ temperature_400hPa(1:nCells) = field_interp(1:nCells,7)
+ temperature_500hPa(1:nCells) = field_interp(1:nCells,8)
+ temperature_600hPa(1:nCells) = field_interp(1:nCells,9)
+ temperature_700hPa(1:nCells) = field_interp(1:nCells,10)
+ temperature_800hPa(1:nCells) = field_interp(1:nCells,11)
+ temperature_850hPa(1:nCells) = field_interp(1:nCells,12)
+ temperature_925hPa(1:nCells) = field_interp(1:nCells,13)
+ temperature_950hPa(1:nCells) = field_interp(1:nCells,14)
+ temperature_1000hPa(1:nCells) = field_interp(1:nCells,15)
! call mpas_log_write('--- end interpolate temperature:')
end if
@@ -712,14 +935,14 @@ subroutine interp_diagnostics(mesh, state, time_lev, diag)
enddo
enddo
call interp_tofixed_pressure(nCells,nVertLevels,nIntP,press_in,field_in,press_interp,field_interp)
- relhum_50hPa(1:nCells) = field_interp(1:nCells,1)
- relhum_100hPa(1:nCells) = field_interp(1:nCells,2)
- relhum_200hPa(1:nCells) = field_interp(1:nCells,3)
- relhum_250hPa(1:nCells) = field_interp(1:nCells,4)
- relhum_500hPa(1:nCells) = field_interp(1:nCells,5)
- relhum_700hPa(1:nCells) = field_interp(1:nCells,6)
- relhum_850hPa(1:nCells) = field_interp(1:nCells,7)
- relhum_925hPa(1:nCells) = field_interp(1:nCells,8)
+ relhum_50hPa(1:nCells) = field_interp(1:nCells,2)
+ relhum_100hPa(1:nCells) = field_interp(1:nCells,3)
+ relhum_200hPa(1:nCells) = field_interp(1:nCells,4)
+ relhum_250hPa(1:nCells) = field_interp(1:nCells,5)
+ relhum_500hPa(1:nCells) = field_interp(1:nCells,8)
+ relhum_700hPa(1:nCells) = field_interp(1:nCells,10)
+ relhum_850hPa(1:nCells) = field_interp(1:nCells,12)
+ relhum_925hPa(1:nCells) = field_interp(1:nCells,13)
! call mpas_log_write('--- end interpolate relative humidity:')
end if
@@ -732,16 +955,43 @@ subroutine interp_diagnostics(mesh, state, time_lev, diag)
enddo
enddo
call interp_tofixed_pressure(nCells,nVertLevels,nIntP,press_in,field_in,press_interp,field_interp)
- dewpoint_50hPa(1:nCells) = field_interp(1:nCells,1)
- dewpoint_100hPa(1:nCells) = field_interp(1:nCells,2)
- dewpoint_200hPa(1:nCells) = field_interp(1:nCells,3)
- dewpoint_250hPa(1:nCells) = field_interp(1:nCells,4)
- dewpoint_500hPa(1:nCells) = field_interp(1:nCells,5)
- dewpoint_700hPa(1:nCells) = field_interp(1:nCells,6)
- dewpoint_850hPa(1:nCells) = field_interp(1:nCells,7)
- dewpoint_925hPa(1:nCells) = field_interp(1:nCells,8)
+ dewpoint_50hPa(1:nCells) = field_interp(1:nCells,2)
+ dewpoint_100hPa(1:nCells) = field_interp(1:nCells,3)
+ dewpoint_200hPa(1:nCells) = field_interp(1:nCells,4)
+ dewpoint_250hPa(1:nCells) = field_interp(1:nCells,5)
+ dewpoint_500hPa(1:nCells) = field_interp(1:nCells,8)
+ dewpoint_700hPa(1:nCells) = field_interp(1:nCells,10)
+ dewpoint_850hPa(1:nCells) = field_interp(1:nCells,12)
+ dewpoint_925hPa(1:nCells) = field_interp(1:nCells,13)
! call mpas_log_write('--- end interpolate relative humidity:')
end if
+
+ if (NEED_QV) then
+ !... QVAPOR
+ do iCell = 1, nCells
+ do k = 1, nVertLevels
+ kk = nVertLevels+1-k
+ field_in(iCell,kk) = scalars(index_qv,k,iCell)
+ enddo
+ enddo
+ call interp_tofixed_pressure(nCells,nVertLevels,nIntP,press_in,field_in,press_interp,field_interp)
+ qv_20hPa(1:nCells) = field_interp(1:nCells,1)
+ qv_50hPa(1:nCells) = field_interp(1:nCells,2)
+ qv_100hPa(1:nCells) = field_interp(1:nCells,3)
+ qv_200hPa(1:nCells) = field_interp(1:nCells,4)
+ qv_250hPa(1:nCells) = field_interp(1:nCells,5)
+ qv_300hPa(1:nCells) = field_interp(1:nCells,6)
+ qv_400hPa(1:nCells) = field_interp(1:nCells,7)
+ qv_500hPa(1:nCells) = field_interp(1:nCells,8)
+ qv_600hPa(1:nCells) = field_interp(1:nCells,9)
+ qv_700hPa(1:nCells) = field_interp(1:nCells,10)
+ qv_800hPa(1:nCells) = field_interp(1:nCells,11)
+ qv_850hPa(1:nCells) = field_interp(1:nCells,12)
+ qv_925hPa(1:nCells) = field_interp(1:nCells,13)
+ qv_950hPa(1:nCells) = field_interp(1:nCells,14)
+ qv_1000hPa(1:nCells) = field_interp(1:nCells,15)
+! call mpas_log_write('--- end interpolate qv:')
+ end if
if (NEED_UZONAL) then
!... u zonal wind:
@@ -752,14 +1002,21 @@ subroutine interp_diagnostics(mesh, state, time_lev, diag)
enddo
enddo
call interp_tofixed_pressure(nCells,nVertLevels,nIntP,press_in,field_in,press_interp,field_interp)
- uzonal_50hPa(1:nCells) = field_interp(1:nCells,1)
- uzonal_100hPa(1:nCells) = field_interp(1:nCells,2)
- uzonal_200hPa(1:nCells) = field_interp(1:nCells,3)
- uzonal_250hPa(1:nCells) = field_interp(1:nCells,4)
- uzonal_500hPa(1:nCells) = field_interp(1:nCells,5)
- uzonal_700hPa(1:nCells) = field_interp(1:nCells,6)
- uzonal_850hPa(1:nCells) = field_interp(1:nCells,7)
- uzonal_925hPa(1:nCells) = field_interp(1:nCells,8)
+ uzonal_20hPa(1:nCells) = field_interp(1:nCells,1)
+ uzonal_50hPa(1:nCells) = field_interp(1:nCells,2)
+ uzonal_100hPa(1:nCells) = field_interp(1:nCells,3)
+ uzonal_200hPa(1:nCells) = field_interp(1:nCells,4)
+ uzonal_250hPa(1:nCells) = field_interp(1:nCells,5)
+ uzonal_300hPa(1:nCells) = field_interp(1:nCells,6)
+ uzonal_400hPa(1:nCells) = field_interp(1:nCells,7)
+ uzonal_500hPa(1:nCells) = field_interp(1:nCells,8)
+ uzonal_600hPa(1:nCells) = field_interp(1:nCells,9)
+ uzonal_700hPa(1:nCells) = field_interp(1:nCells,10)
+ uzonal_800hPa(1:nCells) = field_interp(1:nCells,11)
+ uzonal_850hPa(1:nCells) = field_interp(1:nCells,12)
+ uzonal_925hPa(1:nCells) = field_interp(1:nCells,13)
+ uzonal_950hPa(1:nCells) = field_interp(1:nCells,14)
+ uzonal_1000hPa(1:nCells) = field_interp(1:nCells,15)
! call mpas_log_write('--- end interpolate zonal wind:')
end if
@@ -772,14 +1029,21 @@ subroutine interp_diagnostics(mesh, state, time_lev, diag)
enddo
enddo
call interp_tofixed_pressure(nCells,nVertLevels,nIntP,press_in,field_in,press_interp,field_interp)
- umeridional_50hPa(1:nCells) = field_interp(1:nCells,1)
- umeridional_100hPa(1:nCells) = field_interp(1:nCells,2)
- umeridional_200hPa(1:nCells) = field_interp(1:nCells,3)
- umeridional_250hPa(1:nCells) = field_interp(1:nCells,4)
- umeridional_500hPa(1:nCells) = field_interp(1:nCells,5)
- umeridional_700hPa(1:nCells) = field_interp(1:nCells,6)
- umeridional_850hPa(1:nCells) = field_interp(1:nCells,7)
- umeridional_925hPa(1:nCells) = field_interp(1:nCells,8)
+ umeridional_20hPa(1:nCells) = field_interp(1:nCells,1)
+ umeridional_50hPa(1:nCells) = field_interp(1:nCells,2)
+ umeridional_100hPa(1:nCells) = field_interp(1:nCells,3)
+ umeridional_200hPa(1:nCells) = field_interp(1:nCells,4)
+ umeridional_250hPa(1:nCells) = field_interp(1:nCells,5)
+ umeridional_300hPa(1:nCells) = field_interp(1:nCells,6)
+ umeridional_400hPa(1:nCells) = field_interp(1:nCells,7)
+ umeridional_500hPa(1:nCells) = field_interp(1:nCells,8)
+ umeridional_600hPa(1:nCells) = field_interp(1:nCells,9)
+ umeridional_700hPa(1:nCells) = field_interp(1:nCells,10)
+ umeridional_800hPa(1:nCells) = field_interp(1:nCells,11)
+ umeridional_850hPa(1:nCells) = field_interp(1:nCells,12)
+ umeridional_925hPa(1:nCells) = field_interp(1:nCells,13)
+ umeridional_950hPa(1:nCells) = field_interp(1:nCells,14)
+ umeridional_1000hPa(1:nCells) = field_interp(1:nCells,15)
! call mpas_log_write('--- end interpolate meridional wind:')
end if
@@ -806,14 +1070,21 @@ subroutine interp_diagnostics(mesh, state, time_lev, diag)
enddo
enddo
call interp_tofixed_pressure(nCells,nVertLevelsP1,nIntP,press_in,field_in,press_interp,field_interp)
- height_50hPa(1:nCells) = field_interp(1:nCells,1)
- height_100hPa(1:nCells) = field_interp(1:nCells,2)
- height_200hPa(1:nCells) = field_interp(1:nCells,3)
- height_250hPa(1:nCells) = field_interp(1:nCells,4)
- height_500hPa(1:nCells) = field_interp(1:nCells,5)
- height_700hPa(1:nCells) = field_interp(1:nCells,6)
- height_850hPa(1:nCells) = field_interp(1:nCells,7)
- height_925hPa(1:nCells) = field_interp(1:nCells,8)
+ height_20hPa(1:nCells) = field_interp(1:nCells,1)
+ height_50hPa(1:nCells) = field_interp(1:nCells,2)
+ height_100hPa(1:nCells) = field_interp(1:nCells,3)
+ height_200hPa(1:nCells) = field_interp(1:nCells,4)
+ height_250hPa(1:nCells) = field_interp(1:nCells,5)
+ height_300hPa(1:nCells) = field_interp(1:nCells,6)
+ height_400hPa(1:nCells) = field_interp(1:nCells,7)
+ height_500hPa(1:nCells) = field_interp(1:nCells,8)
+ height_600hPa(1:nCells) = field_interp(1:nCells,9)
+ height_700hPa(1:nCells) = field_interp(1:nCells,10)
+ height_800hPa(1:nCells) = field_interp(1:nCells,11)
+ height_850hPa(1:nCells) = field_interp(1:nCells,12)
+ height_925hPa(1:nCells) = field_interp(1:nCells,13)
+ height_950hPa(1:nCells) = field_interp(1:nCells,14)
+ height_1000hPa(1:nCells) = field_interp(1:nCells,15)
! call mpas_log_write('--- end interpolate height:')
!... vertical velocity
@@ -824,14 +1095,14 @@ subroutine interp_diagnostics(mesh, state, time_lev, diag)
enddo
enddo
call interp_tofixed_pressure(nCells,nVertLevelsP1,nIntP,press_in,field_in,press_interp,field_interp)
- w_50hPa(1:nCells) = field_interp(1:nCells,1)
- w_100hPa(1:nCells) = field_interp(1:nCells,2)
- w_200hPa(1:nCells) = field_interp(1:nCells,3)
- w_250hPa(1:nCells) = field_interp(1:nCells,4)
- w_500hPa(1:nCells) = field_interp(1:nCells,5)
- w_700hPa(1:nCells) = field_interp(1:nCells,6)
- w_850hPa(1:nCells) = field_interp(1:nCells,7)
- w_925hPa(1:nCells) = field_interp(1:nCells,8)
+ w_50hPa(1:nCells) = field_interp(1:nCells,2)
+ w_100hPa(1:nCells) = field_interp(1:nCells,3)
+ w_200hPa(1:nCells) = field_interp(1:nCells,4)
+ w_250hPa(1:nCells) = field_interp(1:nCells,5)
+ w_500hPa(1:nCells) = field_interp(1:nCells,8)
+ w_700hPa(1:nCells) = field_interp(1:nCells,10)
+ w_850hPa(1:nCells) = field_interp(1:nCells,12)
+ w_925hPa(1:nCells) = field_interp(1:nCells,13)
if(allocated(field_in)) deallocate(field_in)
if(allocated(press_in)) deallocate(press_in)
@@ -893,6 +1164,64 @@ subroutine interp_diagnostics(mesh, state, time_lev, diag)
if(allocated(pressureCp1) ) deallocate(pressureCp1 )
if(allocated(pressure_v) ) deallocate(pressure_v )
+
+! if levels are below surface, set them to missing value (defined in the registry)
+ DO iCell = 1, nCells
+ IF(height_1000hPa(icell) > 360.) qv_1000hPa(icell) = -9999.
+ IF(height_950hPa(icell) > 700. ) qv_950hPa(icell) = -9999.
+ IF(height_925hPa(icell) > 1000.) qv_925hPa(icell) = -9999.
+ IF(height_850hPa(icell) > 1700.) qv_850hPa(icell) = -9999.
+ IF(height_800hPa(icell) > 2150.) qv_800hPa(icell) = -9999.
+ IF(height_700hPa(icell) > 3300.) qv_700hPa(icell) = -9999.
+ IF(height_600hPa(icell) > 4550.) qv_600hPa(icell) = -9999.
+ IF(height_500hPa(icell) > 6000.) qv_500hPa(icell) = -9999.
+ IF(height_400hPa(icell) > 7900.) qv_400hPa(icell) = -9999.
+
+ IF(height_1000hPa(icell) > 360.) temperature_1000hPa(icell) = -9999.
+ IF(height_950hPa(icell) > 700. ) temperature_950hPa(icell) = -9999.
+ IF(height_925hPa(icell) > 1000.) temperature_925hPa(icell) = -9999.
+ IF(height_850hPa(icell) > 1700.) temperature_850hPa(icell) = -9999.
+ IF(height_800hPa(icell) > 2150.) temperature_800hPa(icell) = -9999.
+ IF(height_700hPa(icell) > 3300.) temperature_700hPa(icell) = -9999.
+ IF(height_600hPa(icell) > 4550.) temperature_600hPa(icell) = -9999.
+ IF(height_500hPa(icell) > 6000.) temperature_500hPa(icell) = -9999.
+ IF(height_400hPa(icell) > 7900.) temperature_400hPa(icell) = -9999.
+
+ IF(height_1000hPa(icell) > 360.) umeridional_1000hPa(icell) = -9999.
+ IF(height_950hPa(icell) > 700. ) umeridional_950hPa(icell) = -9999.
+ IF(height_925hPa(icell) > 1000.) umeridional_925hPa(icell) = -9999.
+ IF(height_850hPa(icell) > 1700.) umeridional_850hPa(icell) = -9999.
+ IF(height_800hPa(icell) > 2150.) umeridional_800hPa(icell) = -9999.
+ IF(height_700hPa(icell) > 3300.) umeridional_700hPa(icell) = -9999.
+ IF(height_600hPa(icell) > 4550.) umeridional_600hPa(icell) =-9999.
+ IF(height_500hPa(icell) > 6000.) umeridional_500hPa(icell) = -9999.
+ IF(height_400hPa(icell) > 7900.) umeridional_400hPa(icell) = -9999.
+
+ IF(height_1000hPa(icell) > 360.) uzonal_1000hPa(icell) = -9999.
+ IF(height_950hPa(icell) > 700. ) uzonal_950hPa(icell) = -9999.
+ IF(height_925hPa(icell) > 1000.) uzonal_925hPa(icell) = -9999.
+ IF(height_850hPa(icell) > 1700.) uzonal_850hPa(icell) = -9999.
+ IF(height_800hPa(icell) > 2150.) uzonal_800hPa(icell) = -9999.
+ IF(height_700hPa(icell) > 3300.) uzonal_700hPa(icell) = -9999.
+ IF(height_600hPa(icell) > 4550.) uzonal_600hPa(icell) = -9999.
+ IF(height_500hPa(icell) > 6000.) uzonal_500hPa(icell) = -9999.
+ IF(height_400hPa(icell) > 7900.) uzonal_400hPa(icell) = -9999.
+
+ IF(height_925hPa(icell) > 1000.) w_925hPa(icell) = -9999.
+ IF(height_850hPa(icell) > 1700.) w_850hPa(icell) = -9999.
+ IF(height_700hPa(icell) > 3300.) w_700hPa(icell) = -9999.
+ IF(height_500hPa(icell) > 6000.) w_500hPa(icell) = -9999.
+
+ IF(height_925hPa(icell) > 1000.) vorticity_925hPa(icell) = -9999.
+ IF(height_850hPa(icell) > 1700.) vorticity_850hPa(icell) = -9999.
+ IF(height_700hPa(icell) > 3300.) vorticity_700hPa(icell) = -9999.
+ IF(height_500hPa(icell) > 6000.) vorticity_500hPa(icell) = -9999.
+
+ IF(height_925hPa(icell) > 1000.) dewpoint_925hPa(icell) = -9999.
+ IF(height_850hPa(icell) > 1700.) dewpoint_850hPa(icell) = -9999.
+ IF(height_700hPa(icell) > 3300.) dewpoint_700hPa(icell) = -9999.
+ IF(height_500hPa(icell) > 6000.) dewpoint_500hPa(icell) = -9999.
+ ENDDO
if (need_mslp) then
!... compute SLP (requires temp, height, pressure, qvapor)