LFRic2UM pseudo-levels

applications/lfricinputs/source/common/lfricinp_add_um_field_to_file_mod.f90

@@ -11,23 +11,23 @@ module lfricinp_add_um_field_to_file_mod
 ! lfricinp modules
 use lfricinp_check_shumlib_status_mod, only: shumlib
 use lfricinp_grid_type_mod,            only: lfricinp_grid_type
-use lfricinp_stashmaster_mod,          only: get_stashmaster_item, grid,   &
-                                             p_points, u_points, v_points, &
-                                             ozone_points,                 &
-                                             ppfc,                         &
-                                             sm_lbvc => lbvc,              &
-                                             cfff,                         &
-                                             levelt,                       &
-                                             rho_levels, theta_levels,     &
-                                             single_level,                 &
-                                             cfll,                         &
+use lfricinp_stashmaster_mod,          only: get_stashmaster_item, grid,     &
+                                             p_points, u_points, v_points,   &
+                                             ozone_points, land_compressed,  &
+                                             ppfc, p_points_values_over_sea, &
+                                             sm_lbvc => lbvc,                &
+                                             cfff,                           &
+                                             levelt,                         &
+                                             rho_levels, theta_levels,       &
+                                             single_level,                   &
+                                             cfll,                           &
                                              datat
 use lfricinp_um_level_codes_mod,       only: lfricinp_get_first_level_num
-use lfricinp_um_parameters_mod,        only: um_imdi, um_rmdi,             &
-                                             rh_polelat, rh_polelong,      &
-                                             ldc_zsea_theta, ldc_zsea_rho, &
-                                             ldc_c_theta, ldc_c_rho,       &
-                                             rh_deltaEW, rh_deltaNS,       &
+use lfricinp_um_parameters_mod,        only: um_imdi, um_rmdi,               &
+                                             rh_polelat, rh_polelong,        &
+                                             ldc_zsea_theta, ldc_zsea_rho,   &
+                                             ldc_c_theta, ldc_c_rho,         &
+                                             rh_deltaEW, rh_deltaNS,         &
                                              ih_model_levels
 
 
@@ -39,13 +39,13 @@ module lfricinp_add_um_field_to_file_mod
 use f_shum_fieldsfile_mod, only: f_shum_fixed_length_header_len
 use f_shum_file_mod, only: shum_file_type
 use f_shum_field_mod, only: shum_field_type
-use f_shum_lookup_indices_mod, only:                                        &
-    lbyr, lbmon, lbdat, lbhr, lbmin, lbday, lbsec, lbyrd, lbmond, lbdatd,   &
-    lbhrd, lbmind, lbdayd, lbsecd, lbtim, lbft, lbcode, lbhem, lbrow,       &
-    lbnpt, lbpack, lbrel, lbfc, lbcfc, lbproc, lbvc, lbrvc, lbtyp, lblev,   &
-    lbrsvd1, lbrsvd2, lbrsvd3, lbrsvd4, lbsrce, lbuser1, lbuser4, lbuser7,  &
-    bulev, bhulev, brsvd3, brsvd4, bdatum, bacc, blev, brlev,               &
-    bhlev, bhrlev, bplat, bplon, bgor, bzy, bdy, bzx, bdx, bmdi, bmks
+use f_shum_lookup_indices_mod, only:                                         &
+    lbyr, lbmon, lbdat, lbhr, lbmin, lbday, lbsec, lbyrd, lbmond, lbdatd,    &
+    lbhrd, lbmind, lbdayd, lbsecd, lbtim, lbft, lbcode, lbhem, lbrow, lbnpt, &
+    lbpack, lbrel, lbfc, lbcfc, lbproc, lbvc, lbrvc, lbtyp, lblev, lbsrce,   &
+    lbrsvd1, lbrsvd2, lbrsvd3, lbrsvd4, lbuser1, lbuser4, lbuser5, lbuser7,  &
+    bulev, bhulev, brsvd3, brsvd4, bdatum, bacc, blev, brlev, bhlev,         &
+    bhrlev, bplat, bplon, bgor, bzy, bdy, bzx, bdx, bmdi, bmks
 
 use f_shum_fixed_length_header_indices_mod, only:                           &
     vert_coord_type, horiz_grid_type, dataset_type, run_identifier,         &
@@ -138,7 +138,7 @@ subroutine lfricinp_add_um_field_to_file(um_file, stashcode, level_index, &
 
 ! Set horiz grid
 select case(grid_type_code)
-case(p_points, ozone_points)
+case(p_points, ozone_points, land_compressed, p_points_values_over_sea)
   lookup_int(lbrow) = um_grid%num_p_points_y
   lookup_int(lbnpt) = um_grid%num_p_points_x
   ! "Zeroth" start lat/lon so subtract one grid spacing
@@ -225,14 +225,52 @@ subroutine lfricinp_add_um_field_to_file(um_file, stashcode, level_index, &
 
 ! Check vertical coordinate type
 if (lookup_int(lbvc) >= 126 .and. lookup_int(lbvc) <= 139 &
-     .or. lookup_int(lbvc) == 5) then
+     .or. lookup_int(lbvc) == 5 .or. lookup_int(lbvc) == 0 .or. &
+     lookup_int(lbvc) == 275 ) then
+  write(log_scratch_space, '(A,I0,A)') &
+     "Vertical coord type ", lookup_int(lbvc), " treated as single layer"
+  call log_event(log_scratch_space, LOG_LEVEL_INFO)
+  ! Pseudo-level number
+  ! lookup_int(lbuser5) = 5_int64
+
   ! Special codes inc single level, set to 0.0
   lookup_real_tmp(blev)=0.0_real64
   lookup_real_tmp(bhlev)=0.0_real64
   lookup_real_tmp(brlev)=0.0_real64
   lookup_real_tmp(bhrlev)=0.0_real64
   lookup_real_tmp(bulev)=0.0_real64
   lookup_real_tmp(bhulev)=0.0_real64
+
+else if (lookup_int(lbvc) == 6) then ! Deep soil levels
+  ! These are hardcoded to the settings in a UM dump file with 4 soil levels as
+  ! that is currently hardcoded in elsewhere in lfric2um. If at some point that
+  ! gets changed to not be hardcoded, this will also need to change
+  ! bulev is the same as brsvd1
+  write(log_scratch_space, '(A,I0,A)') &
+    "Vertical coord type ", lookup_int(lbvc), " treated as soil field"
+  call log_event(log_scratch_space, LOG_LEVEL_INFO)
+  if (level_number == 1) then
+    lookup_real_tmp(bulev) = 0.0_real64
+    lookup_real_tmp(blev)=0.05_real64
+    lookup_real_tmp(brlev)=0.1_real64
+  else if (level_number == 2) then
+    lookup_real_tmp(bulev) = 0.1_real64
+    lookup_real_tmp(blev)=0.225_real64
+    lookup_real_tmp(brlev)=0.35_real64
+  else if (level_number == 3) then
+    lookup_real_tmp(bulev) = 0.35_real64
+    lookup_real_tmp(blev)=0.675_real64
+    lookup_real_tmp(brlev)=1.0_real64
+  else if (level_number == 4) then
+    lookup_real_tmp(bulev) = 1.0_real64
+    lookup_real_tmp(blev)=2.0_real64
+    lookup_real_tmp(brlev)=3.0_real64
+  else
+    write(log_scratch_space, '(A,I0,A)') "Soil level number ", level_number, &
+      " not supported. Only soil fields with 4 levels are supported currently"
+    call log_event(log_scratch_space, LOG_LEVEL_ERROR)
+  end if
+
 else if (lookup_int(lbvc) == 65) then ! Standard hybrid height levels
   ! height of model level k above mean sea level is
   !       z(i,j,k) = Zsea(k) + C(k)*Zorog(i,j)
@@ -241,6 +279,7 @@ subroutine lfricinp_add_um_field_to_file(um_file, stashcode, level_index, &
   ! brlev,bhrlev      zsea,C of lower level boundary
   ! The level here can refer to either a theta or rho level, with
   ! layer boundaries defined by surrounding rho or theta levels.
+
   if (level_code == theta_levels) then ! theta level (& w)
 
     ! When referencing theta arrays need to add 1 to the level number as level

applications/lfricinputs/source/common/lfricinp_create_lfric_fields_mod.f90

@@ -147,7 +147,7 @@ subroutine lfricinp_create_lfric_fields( mesh, twod_mesh,              &
           ndata_64 = lfricinp_get_num_pseudo_levels(um_grid, stashcode)
           write(log_scratch_space, '(A,I0,A)')                                   &
              "This field has ", ndata_64, " pseudo levels"
-          call log_event(log_scratch_space, LOG_LEVEL_DEBUG)
+          call log_event(log_scratch_space, LOG_LEVEL_INFO)
         end if
         ndata_first = .false.
 
@@ -160,7 +160,7 @@ subroutine lfricinp_create_lfric_fields( mesh, twod_mesh,              &
         ndata_64 = lfricinp_get_num_levels(um_file, stashcode)
         write(log_scratch_space, '(A,I0,A)')                                   &
             "This field has ", ndata_64, " levels"
-        call log_event(log_scratch_space, LOG_LEVEL_DEBUG)
+        call log_event(log_scratch_space, LOG_LEVEL_INFO)
         ndata_first = .true.
 
       case DEFAULT

applications/lfricinputs/source/common/lfricinp_um_grid_mod.f90

@@ -9,7 +9,7 @@ module lfricinp_um_grid_mod
 use, intrinsic :: iso_fortran_env, only : int64
 
 ! LFRic modules
-use log_mod,  only: log_event, LOG_LEVEL_INFO, LOG_LEVEL_ERROR
+use log_mod, only: log_event, log_scratch_space, LOG_LEVEL_INFO, LOG_LEVEL_ERROR
 
 ! UM2LFRic modules
 use lfricinp_grid_type_mod, only: lfricinp_grid_type
@@ -21,7 +21,7 @@ module lfricinp_um_grid_mod
 
 private
 
-public :: lfricinp_set_grid_from_file
+public :: lfricinp_set_grid_from_file, lfricinp_set_grid_from_namelist
 
 contains
 
@@ -55,4 +55,22 @@ subroutine lfricinp_set_grid_from_file(um_input_file, num_snow_layers, &
 
 end subroutine lfricinp_set_grid_from_file
 
+!-----------------------------------------------------------
+
+subroutine lfricinp_set_grid_from_namelist(num_snow_layers, num_surface_types)
+! Description:
+! Extracts grid information from lfric2um namelist to populate grid_info object
+
+implicit none
+
+! Arguments
+integer(kind=int64), intent(in) :: num_snow_layers
+integer(kind=int64), intent(in) :: num_surface_types
+
+! Set pseudo level information
+um_grid%num_snow_layers = num_snow_layers
+um_grid%num_surface_types = num_surface_types
+
+end subroutine lfricinp_set_grid_from_namelist
+
 end module lfricinp_um_grid_mod

applications/lfricinputs/source/common/lfricinp_um_level_codes_mod.f90

@@ -17,7 +17,7 @@ module lfricinp_um_level_codes_mod
 
 ! lfric modules
 use log_mod,  only : log_event, log_scratch_space,         &
-                     LOG_LEVEL_ERROR
+                     LOG_LEVEL_DEBUG, LOG_LEVEL_ERROR, LOG_LEVEL_INFO
 
 implicit none
 
@@ -147,6 +147,7 @@ function lfricinp_get_num_pseudo_levels(um_grid, stashcode) &
 ! Result
 integer(kind=int64) :: num_pseudo_levels
 
+
 num_pseudo_levels = lfricinp_get_last_pseudo_level_num(um_grid, stashcode) - &
              lfricinp_get_first_pseudo_level_num(stashcode) + 1
 
@@ -207,14 +208,23 @@ function lfricinp_get_last_pseudo_level_num(um_grid, stashcode) &
      'lfricinp_get_last_pseudo_level_num'
 
 last_pseudo_level_code =  get_stashmaster_item(stashcode, pseudl)
+write(log_scratch_space, '(A,I0)') &
+     "Last pseudo level code is ", last_pseudo_level_code
+call log_event(log_scratch_space, LOG_LEVEL_DEBUG)
 
 select case(last_pseudo_level_code)
 case (7,9) ! ntypes == ntiles (lfricinputs doesn't support aggregate tile)
   last_pseudo_level_num = um_grid % num_surface_types
+  if (last_pseudo_level_num == -32768) then
+    call log_event("Number of surface types is not set", LOG_LEVEL_ERROR)
+  end if
 case (10)
   last_pseudo_level_num = um_grid%num_ice_cats
 case (11)
   last_pseudo_level_num = um_grid%num_snow_layers * um_grid%num_surface_types
+  if ( um_grid%num_snow_layers == -32768 .or. um_grid%num_surface_types == -32768 ) then
+    call log_event("Either the number of surface types or the number of snow layers is not set", LOG_LEVEL_ERROR)
+  end if
 case DEFAULT
   write(log_scratch_space, '(A,I0,A,I0)') &
        "Last pseudo_level code ", last_pseudo_level_code, &
@@ -224,5 +234,4 @@ function lfricinp_get_last_pseudo_level_num(um_grid, stashcode) &
 
 end function lfricinp_get_last_pseudo_level_num
 
-
 end module lfricinp_um_level_codes_mod

applications/lfricinputs/source/lfric2um/lfric2um_initialise_lfric2um_mod.f90

@@ -17,10 +17,14 @@ module lfric2um_initialise_lfric2um_mod
 !> This includes reading namelists and stashmaster file,
 !> and reading in gridding weights
 subroutine lfric2um_initialise_lfric2um()
-use lfric2um_namelists_mod, only: lfric2um_config
-use lfricinp_stashmaster_mod, only: lfricinp_read_stashmaster
+
+use lfricinp_stashmaster_mod,        only: lfricinp_read_stashmaster
 use lfricinp_stash_to_lfric_map_mod, only: lfricinp_init_stash_to_lfric_map
-use lfric2um_regrid_weights_mod, only: lfric2um_regrid_weightsfile_ctl
+use lfricinp_um_grid_mod,            only: lfricinp_set_grid_from_namelist
+
+use lfric2um_namelists_mod,          only: lfric2um_config
+use lfric2um_regrid_weights_mod,     only: lfric2um_regrid_weightsfile_ctl
+
 implicit none
 
 ! Read namelists
@@ -29,6 +33,10 @@ subroutine lfric2um_initialise_lfric2um()
 ! Read in STASHmaster file
 call lfricinp_read_stashmaster(lfric2um_config%stashmaster_file)
 
+! Set um_grid from lfric2um namelist
+call lfricinp_set_grid_from_namelist(lfric2um_config%num_snow_layers, &
+                                     lfric2um_config%num_surface_types)
+
 ! Read in weights files
 call lfric2um_regrid_weightsfile_ctl()
 

applications/lfricinputs/source/lfric2um/lfric2um_initialise_um_mod.f90

@@ -182,7 +182,8 @@ subroutine lfric2um_set_integer_constants(um_output_file, um_grid, lfric_mesh)
 int_constants(ih_rows) = um_grid%num_p_points_y
 int_constants(ih_model_levels) = lfric_mesh%get_nlayers()
 int_constants(ih_wet_levels) = int_constants(ih_model_levels)
-! hardcode for now, no info on soil in lfric dump
+! hardcode to 4 for now, no info on soil in lfric dump.
+! TODO Potentially get from namelists?
 int_constants(ih_soilT_levels) = 4
 int_constants(ih_soilQ_levels) = int_constants(ih_soilT_levels)
 int_constants(ih_cloud_levels) = lfric_mesh%get_nlayers()
@@ -193,10 +194,15 @@ subroutine lfric2um_set_integer_constants(um_output_file, um_grid, lfric_mesh)
 ! hardcode for now - doesn't exist in lfric yet
 int_constants(ih_1_c_rho_level) = 30
 !  hardcode for now - aquaplanet
+! TODO unhardcode - no longer aquaplanet!
 int_constants(ih_land_points) = 0
 int_constants(ih_ozone_levels) = int_constants(ih_model_levels)
 int_constants(ih_convect_levels) = 0
 
+! Set pseudo-level information
+! um_output_file%num_snow_layers = um_grid%num_snow_layers
+! um_output_file%num_surface_types = um_grid%num_surface_types
+
 call shumlib(routinename//'::set_integer_constants', &
      um_output_file%set_integer_constants(int_constants))
 

applications/lfricinputs/source/lfric2um/lfric2um_main_loop_mod.f90

@@ -23,14 +23,18 @@ module lfric2um_main_loop_mod
 use lfricinp_lfric_driver_mod,         only: lfric_fields, local_rank, comm,   &
                                              twod_mesh
 use lfricinp_regrid_weights_type_mod,  only: lfricinp_regrid_weights_type
-use lfricinp_stash_to_lfric_map_mod,   only: get_field_name
+use lfricinp_stash_to_lfric_map_mod,   only: get_field_name,                   &
+                                             get_lfric_field_kind,             &
+                                             w2h_field, w3_field, w3_field_2d, &
+                                             w3_soil_field, wtheta_field
 use lfricinp_stashmaster_mod,          only: get_stashmaster_item, levelt,     &
-                                             rho_levels, single_level,         &
+                                             pseudt, rho_levels, single_level, &
                                              stashcode_exner, stashcode_p,     &
                                              stashcode_q, stashcode_theta,     &
                                              theta_levels
 use lfricinp_um_grid_mod,              only: um_grid
-use lfricinp_um_level_codes_mod,       only: lfricinp_get_num_levels
+use lfricinp_um_level_codes_mod,       only: lfricinp_get_num_levels,          &
+                                             lfricinp_get_num_pseudo_levels
 
 ! lfric modules
 use field_mod,  only: field_type
@@ -53,7 +57,7 @@ subroutine lfric2um_main_loop()
 implicit none
 
 integer(kind=int64) :: i_stash, level, i_field
-integer(kind=int64) :: stashcode, num_levels
+integer(kind=int64) :: stashcode, num_levels, lfric_field_kind
 character(len=*), parameter :: routinename='lfric2um_main_loop'
 type(field_type), pointer :: lfric_field
 type(lfricinp_regrid_weights_type), pointer :: weights
@@ -76,7 +80,43 @@ subroutine lfric2um_main_loop()
   stashcode = lfric2um_config%stash_list(i_stash)
   write(log_scratch_space, '(A,I0)') "Processing stashcode ", stashcode
   call log_event(log_scratch_space, LOG_LEVEL_INFO)
-  num_levels = lfricinp_get_num_levels(um_output_file, stashcode)
+
+  !---------------------------------------------------------------------------
+  ! Identify if field uses pseudo levels
+  !---------------------------------------------------------------------------
+  lfric_field_kind = get_lfric_field_kind(stashcode)
+  select case (lfric_field_kind)
+
+    !---------------------------------------------------------------------------
+    ! Stashcodes that map to W2h (i.e. winds), W3/rho, Wtheta
+    !---------------------------------------------------------------------------
+    case(w2h_field,w3_field,wtheta_field)
+      num_levels = lfricinp_get_num_levels(um_output_file, stashcode)
+
+    !---------------------------------------------------------------------------
+    ! Stashcodes that need 2D mesh, inc non-soil_levels pseudo levels
+    !---------------------------------------------------------------------------
+    case(w3_field_2d)
+      if ( get_stashmaster_item(stashcode, pseudt) == 0 ) then
+        ! Field has no pseudo levels
+        num_levels = lfricinp_get_num_levels(um_output_file, stashcode)
+      else
+        ! Get number of pseudo levels
+        num_levels = lfricinp_get_num_pseudo_levels(um_grid, stashcode)
+      end if
+
+    !---------------------------------------------------------------------------
+    ! Stashcodes that are soil fields
+    !---------------------------------------------------------------------------
+    case(w3_soil_field)
+      num_levels = lfricinp_get_num_levels(um_output_file, stashcode)
+
+    case DEFAULT
+      write(log_scratch_space, '(A,I0,A)')                                     &
+         "LFRic field kind code ", lfric_field_kind, " not recognised"
+      call log_event(log_scratch_space, LOG_LEVEL_ERROR)
+
+  end select
 
   !---------------------------------------------------------------------------
   ! Select appropriate weights
@@ -120,6 +160,8 @@ subroutine lfric2um_main_loop()
       call lfricinp_add_um_field_to_file(um_output_file, stashcode, &
            level, um_grid, lfric2um_config%lbtim_list(i_stash),     &
            lfric2um_config%lbproc_list(i_stash))
+      ! I think this + loop it's in is what is causing the pseudolevel fields to
+      ! be added as 9 individual fields rather than a field with 9 pseudolevels
 
       !-------------------------------------------------------------------------
       ! Adding a 2D UM field to the file increments num_fields by 1 each time