Merge branch 'release/MAPL-v3' into feature/tclune/#3163-refactor-sha…

…redio

CHANGELOG.md

@@ -70,6 +70,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 - Allow update offsets of ±timestep in ExtData2G
 - Minor revision (and generalization) of grid-def for GSI purposes
+- Trajectory sampler: fix a bug when group_name does not exist in netCDF file and a bug that omitted the first time point
 
 ### Changed
 

base/Plain_netCDF_Time.F90

@@ -218,20 +218,21 @@ subroutine get_v1d_netcdf_R8(filename, name, array, Xdim, group_name, rc)
     real(REAL64), dimension(Xdim), intent(out) :: array
     integer, optional, intent(out) :: rc
     integer :: status
-    integer :: ncid, varid, ncid2
+    integer :: ncid, varid, ncid2, ncid_sv
 
     call check_nc_status(nf90_open(trim(fileName), NF90_NOWRITE, ncid), _RC)
+    ncid_sv = ncid
+
     if(present(group_name)) then
-       ncid2= ncid
-       call check_nc_status(nf90_inq_ncid(ncid2, group_name, ncid), _RC)
+       if(group_name/='') then
+          ncid2= ncid
+          call check_nc_status(nf90_inq_ncid(ncid2, group_name, ncid), _RC)
+       end if
     end if
     call check_nc_status(nf90_inq_varid(ncid, name, varid), _RC)
     call check_nc_status(nf90_get_var(ncid, varid, array), _RC)
-    if(present(group_name)) then
-       call check_nc_status(nf90_close(ncid2), _RC)
-    else
-       call check_nc_status(nf90_close(ncid), _RC)
-    end if
+
+    call check_nc_status(nf90_close(ncid_sv), _RC)
     _RETURN(_SUCCESS)
 
   end subroutine get_v1d_netcdf_R8
@@ -256,9 +257,11 @@ subroutine get_v1d_netcdf_R8_complete(filename, varname, array, att_name, att_va
     call check_nc_status(nf90_open(trim(fileName), NF90_NOWRITE, ncid), _RC)
     ncid_sv = ncid
     if(present(group_name)) then
-       call check_nc_status(nf90_inq_ncid(ncid, group_name, ncid_grp), _RC)
-       ! mod
-       ncid = ncid_grp
+       if(group_name/='') then
+          call check_nc_status(nf90_inq_ncid(ncid, group_name, ncid_grp), _RC)
+          ! mod
+          ncid = ncid_grp
+       end if
     end if
     call check_nc_status(nf90_inq_varid(ncid, varname, varid), _RC)
     call check_nc_status(nf90_get_var(ncid, varid, array), _RC)
@@ -296,9 +299,11 @@ subroutine get_att_real_netcdf(filename, varname, att_name, att_value, group_nam
     call check_nc_status(nf90_open(trim(fileName), NF90_NOWRITE, ncid), _RC)
     ncid_sv = ncid
     if(present(group_name)) then
-       call check_nc_status(nf90_inq_ncid(ncid, group_name, ncid_grp), _RC)
-       ! overwrite
-       ncid = ncid_grp
+       if(group_name/='') then
+          call check_nc_status(nf90_inq_ncid(ncid, group_name, ncid_grp), _RC)
+          ! overwrite
+          ncid = ncid_grp
+       end if
     end if
     call check_nc_status(nf90_inq_varid(ncid, varname, varid), _RC)
     call check_nc_status(nf90_get_att(ncid, varid, att_name, att_value), _RC)
@@ -324,9 +329,11 @@ subroutine get_att_char_netcdf(filename, varname, att_name, att_value, group_nam
     call check_nc_status(nf90_open(trim(fileName), NF90_NOWRITE, ncid), _RC)
     ncid_sv = ncid
     if(present(group_name)) then
-       call check_nc_status(nf90_inq_ncid(ncid, group_name, ncid_grp), _RC)
-       ! overwrite
-       ncid = ncid_grp
+       if(group_name/='') then
+          call check_nc_status(nf90_inq_ncid(ncid, group_name, ncid_grp), _RC)
+          ! overwrite
+          ncid = ncid_grp
+       end if
     end if
     call check_nc_status(nf90_inq_varid(ncid, varname, varid), _RC)
     call check_nc_status(nf90_get_att(ncid, varid, att_name, att_value), _RC)

generic3g/ComponentSpecParser/parse_geometry_spec.F90

@@ -108,7 +108,8 @@ module function parse_geometry_spec(mapl_cfg, registry, rc) result(geometry_spec
             vertical_grid = FixedLevelsVerticalGrid(standard_name, levels, units)
          case('model')
             num_levels = ESMF_HConfigAsI4(vertical_grid_cfg, keyString='num_levels', _RC)
-            vertical_grid = ModelVerticalGrid(num_levels=num_levels)
+            units = ESMF_HConfigAsString(vertical_grid_cfg, keyString='units', _RC)
+            vertical_grid = ModelVerticalGrid(num_levels=num_levels, units=units)
             short_name = ESMF_HConfigAsString(vertical_grid_cfg, keyString='short_name', _RC)
             select type(vertical_grid)
             type is(ModelVerticalGrid)

generic3g/registry/ExtensionFamily.F90

@@ -122,6 +122,7 @@ function find_closest_extension(family, goal_spec, rc) result(closest_extension)
       type(StateItemExtensionPtr) :: extension_ptr
       type(StateItemExtension), pointer :: primary
       class(StateItemSpec), pointer :: spec
+      logical :: match
 
       closest_extension => null()
       subgroup = family%get_extensions()
@@ -135,7 +136,8 @@ function find_closest_extension(family, goal_spec, rc) result(closest_extension)
             extension_ptr = subgroup%of(j)
             spec => extension_ptr%ptr%get_spec()
             associate (adapter => adapters(i)%adapter)
-              if (adapter%match(spec)) then
+              match = adapter%match(spec, _RC)
+              if (match) then
                  call new_subgroup%push_back(extension_ptr)
               end if
             end associate

generic3g/registry/StateItemExtension.F90

@@ -117,14 +117,16 @@ recursive function make_extension(this, goal, rc) result(extension)
       type(ESMF_GridComp) :: coupler_gridcomp
       type(StateItemAdapterWrapper), allocatable :: adapters(:)
       type(ESMF_Clock) :: fake_clock
+      logical :: match
 
       call this%spec%set_active()
 
       new_spec = this%spec
       adapters = this%spec%make_adapters(goal, _RC)
       do i = 1, size(adapters)
-         if (adapters(i)%adapter%match(new_spec)) cycle
-         call adapters(i)%adapter%adapt(new_spec, action)
+         match = adapters(i)%adapter%match(new_spec, _RC)
+         if (match) cycle
+         call adapters(i)%adapter%adapt(new_spec, action, _RC)
          exit
       end do
 

generic3g/specs/FieldSpec.F90

@@ -360,13 +360,16 @@ subroutine write_formatted(this, unit, iotype, v_list, iostat, iomsg)
          write(unit, "(3x, a, a, a)", iostat=iostat, iomsg=iomsg) "long name:", this%long_name, new_line("a")
       end if
       if (allocated(this%units)) then
-         write(unit, "(3x, a, a, a)", iostat=iostat, iomsg=iomsg) "unit:", this%units, new_line("a")
+         write(unit, "(3x, a, a, a)", iostat=iostat, iomsg=iomsg) "units:", this%units, new_line("a")
       end if
       write(unit, "(3x, dt'g0', a)", iostat=iostat, iomsg=iomsg) this%vertical_dim_spec, new_line("a")
       if (allocated(this%vertical_grid)) then
-         write(unit, "(3x, dt'g0', a)", iostat=iostat, iomsg=iomsg) this%vertical_grid, new_line("a")
+         write(unit, "(3x, dt'g0', a)", iostat=iostat, iomsg=iomsg) this%vertical_grid
       end if
       write(unit, "(a)") ")"
+
+      _UNUSED_DUMMY(iotype)
+      _UNUSED_DUMMY(v_list)
    end subroutine write_formatted
 
    function get_ungridded_bounds(this, rc) result(bounds)
@@ -585,8 +588,8 @@ logical function can_connect_to(this, src_spec, rc)
          can_connect_to = all ([ &
               can_match(this%geom,src_spec%geom), &
               can_match(this%vertical_grid, src_spec%vertical_grid), &
-              match(this%vertical_dim_spec,src_spec%vertical_dim_spec), &
-              match(this%ungridded_dims,src_spec%ungridded_dims), &
+              match(this%vertical_dim_spec, src_spec%vertical_dim_spec), &
+              match(this%ungridded_dims, src_spec%ungridded_dims), &
               includes(this%attributes, src_spec%attributes), &
               can_convert_units &
               ])
@@ -843,13 +846,23 @@ subroutine adapt_vertical_grid(this, spec, action, rc)
       type(GriddedComponentDriver), pointer :: v_in_coupler
       type(GriddedComponentDriver), pointer :: v_out_coupler
       type(ESMF_Field) :: v_in_coord, v_out_coord
+      type(ESMF_TypeKind_Flag) :: typekind_in, typekind_out
       integer :: status
 
       select type (spec)
       type is (FieldSpec)
-         call spec%vertical_grid%get_coordinate_field(v_in_coord, v_in_coupler, &
-              'ignore', spec%geom, spec%typekind, spec%units, spec%vertical_dim_spec, _RC)
-         call this%vertical_grid%get_coordinate_field(v_out_coord, v_out_coupler, &
+         ! TODO: DO WE NEED TO RESTRICT SPEC's VERTICAL GRID TO MODEL?
+         ! NOTE: we cannot import ModelVerticalGrid (circular dependency)
+         _ASSERT(spec%vertical_grid%get_units() == this%vertical_grid%get_units(), 'units must match')
+         _ASSERT(spec%vertical_dim_spec == this%vertical_dim_spec, 'temporary restriction')
+         ! Field (to be regridded) should have the same typekind as the underlying vertical grid
+         ! TODO: Should we add a typekind class variable to VerticalGrid?
+         _ASSERT(spec%typekind == this%typekind, 'typekind must match')
+         call spec%vertical_grid%get_coordinate_field( &
+              v_in_coord, v_in_coupler, & ! output
+              'ignore', spec%geom, spec%typekind, this%vertical_grid%get_units(), spec%vertical_dim_spec, _RC)
+         call this%vertical_grid%get_coordinate_field( &
+              v_out_coord, v_out_coupler, & ! output
               'ignore', this%geom, this%typekind, this%units, this%vertical_dim_spec, _RC)
          action = VerticalRegridAction(v_in_coord, v_out_coupler, v_out_coord, v_out_coupler, this%regrid_method)
          spec%vertical_grid = this%vertical_grid
@@ -975,7 +988,7 @@ subroutine adapt_units(this, spec, action, rc)
       _RETURN(_SUCCESS)
    end subroutine adapt_units
 
-  logical function adapter_match_units(this, spec, rc) result(match)
+   logical function adapter_match_units(this, spec, rc) result(match)
       class(UnitsAdapter), intent(in) :: this
       class(StateItemSpec), intent(in) :: spec
       integer, optional, intent(out) :: rc

generic3g/tests/Test_ModelVerticalGrid.pf

@@ -65,19 +65,19 @@ contains
       rc = 0
       ! Inside user "set_geom" phase.
       geom = make_geom(_RC)
-      vgrid = ModelVerticalGrid(num_levels=LM)
+      vgrid = ModelVerticalGrid(num_levels=LM, units="hPa")
       call vgrid%add_variant(short_name=var_name)
 
       ! inside OuterMeta
-      r = StateRegistry('dyn')
+      r = StateRegistry("dyn")
       call vgrid%set_registry(r) ! MAPL_SetVerticalGrid(...)
 
-      v_pt = VirtualConnectionPt(state_intent='export', short_name=var_name)
+      v_pt = VirtualConnectionPt(state_intent="export", short_name=var_name)
       var_spec = VariableSpec(&
            short_name=var_name, &
            state_intent=ESMF_STATEINTENT_EXPORT, &
-           standard_name='air_pressure', &
-           units='hPa', &
+           standard_name="air_pressure", &
+           units="hPa", &
            vertical_dim_spec=vertical_dim_spec, &
            default_value=3.)
       allocate(fld_spec, source=make_itemSpec(var_spec, r, rc=status))
@@ -96,41 +96,37 @@ contains
    function make_geom(rc) result(geom)
       integer, intent(out) :: rc
       type(ESMF_Geom) :: geom
-      type(ESMF_Grid) :: grid
       integer :: status
       type(ESMF_HConfig) :: hconfig
       type(GeomManager), pointer :: geom_mgr
-      class(GeomSpec), allocatable :: geom_spec
       type(MaplGeom), pointer :: mapl_geom
 
       rc = 0
       geom_mgr => get_geom_manager()
-      hconfig = ESMF_HConfigCreate(content='{class: latlon, im_world: 6, jm_world: 7, pole: PC, dateline: DC}', _RC)
+      hconfig = ESMF_HConfigCreate(content="{class: latlon, im_world: 6, jm_world: 7, pole: PC, dateline: DC}", _RC)
       mapl_geom => geom_mgr%get_mapl_geom(hconfig, _RC)
       geom = mapl_geom%get_geom()
    end function make_geom
 
    @test
    subroutine test_num_levels()
       type(ModelVerticalGrid) :: vgrid
-
       integer :: num_levels
 
       num_levels = 10
-      vgrid = ModelVerticalGrid(num_levels=num_levels)
+      vgrid = ModelVerticalGrid(num_levels=num_levels, units="hPa")
       @assert_that(vgrid%get_num_levels(), is(num_levels))
    end subroutine test_num_levels
 
    @test
    subroutine test_num_variants()
       type(ModelVerticalGrid) :: vgrid
-      integer :: num_variants
 
-      vgrid = ModelVerticalGrid(num_levels=3)
+      vgrid = ModelVerticalGrid(num_levels=3, units="hPa")
       @assert_that(vgrid%get_num_variants(), is(0))
-      call vgrid%add_variant(short_name='PLE')
+      call vgrid%add_variant(short_name="PLE")
       @assert_that(vgrid%get_num_variants(), is(1))
-      call vgrid%add_variant(short_name='ZLE')
+      call vgrid%add_variant(short_name="ZLE")
       @assert_that(vgrid%get_num_variants(), is(2))
    end subroutine test_num_variants
 
@@ -149,17 +145,18 @@ contains
 
       call setup("PLE", vgrid, _RC)
 
-      ple_pt = VirtualConnectionPt(state_intent='export', short_name='PLE')
+      ple_pt = VirtualConnectionPt(state_intent="export", short_name="PLE")
       extension => r%get_primary_extension(ple_pt, _RC)
       spec => extension%get_spec()
 
       multi_state = MultiState()
       call spec%add_to_state(multi_state, ActualConnectionPt(ple_pt), _RC)
-      call ESMF_StateGet(multi_state%exportState, itemName='PLE', field=ple, _RC)
+      call ESMF_StateGet(multi_state%exportState, itemName="PLE", field=ple, _RC)
       call ESMF_FieldGet(ple, rank=rank, _RC)
       allocate(localElementCount(rank))
       call ESMF_FieldGet(ple, localElementCount=localElementCount, _RC)
       @assert_that(localElementCount, is(equal_to([IM,JM,LM+1])))
+      _UNUSED_DUMMY(this)
    end subroutine test_created_fields_have_num_levels
 
    @test(type=ESMF_TestMethod, npes=[1])
@@ -180,16 +177,17 @@ contains
 
       call vgrid%get_coordinate_field( &
            vcoord, coupler, &
-           standard_name='air_pressure', &
+           standard_name="air_pressure", &
            geom=geom, &
            typekind=ESMF_TYPEKIND_R4, &
-           units='hPa', &
+           units="hPa", &
            vertical_dim_spec=VERTICAL_DIM_EDGE, &
            _RC)
       @assert_that(associated(coupler), is(false()))
 
       call ESMF_FieldGet(vcoord, fArrayPtr=a, _RC)
       @assert_that(a, every_item(is(equal_to(3.))))
+      _UNUSED_DUMMY(this)
    end subroutine test_get_coordinate_field_simple
 
    @test(type=ESMF_TestMethod, npes=[1])
@@ -213,10 +211,10 @@ contains
 
       call vgrid%get_coordinate_field( &
            vcoord, coupler, &
-           standard_name='air_pressure', &
+           standard_name="air_pressure", &
            geom=geom, &
            typekind=ESMF_TYPEKIND_R4, &
-           units='Pa', &
+           units="Pa", &
            vertical_dim_spec=VERTICAL_DIM_EDGE, &
            _RC)
       @assert_that(associated(coupler), is(true()))
@@ -234,6 +232,7 @@ contains
       end do
       @assert_that(shape(a), is(equal_to([IM, JM, LM+1])))
       @assert_that(a, every_item(is(equal_to(300.))))
+      _UNUSED_DUMMY(this)
    end subroutine test_get_coordinate_field_change_units_edge
 
    @test(type=ESMF_TestMethod, npes=[1])
@@ -257,9 +256,9 @@ contains
 
       call vgrid%get_coordinate_field( &
            vcoord, coupler, &
-           standard_name='air_pressure', &
+           standard_name="air_pressure", &
            geom=geom, &
-           typekind=ESMF_TYPEKIND_R4, units='Pa', &
+           typekind=ESMF_TYPEKIND_R4, units="Pa", &
            vertical_dim_spec=VERTICAL_DIM_CENTER, &
            _RC)
       @assert_that(associated(coupler), is(true()))
@@ -277,6 +276,7 @@ contains
       end do
       @assert_that(shape(a), is(equal_to([IM, JM, LM])))
       @assert_that(a, every_item(is(equal_to(300.))))
+      _UNUSED_DUMMY(this)
    end subroutine test_get_coordinate_field_change_units_center
 
 end module Test_ModelVerticalGrid

generic3g/tests/Test_Scenarios.pf

@@ -127,9 +127,13 @@ contains
               ScenarioDescription('export_dependency',       'parent.yaml',      check_name, check_stateitem), &
               ScenarioDescription('regrid',                  'cap.yaml',         check_name, check_stateitem), &
               ScenarioDescription('propagate_geom',          'parent.yaml',      check_name, check_stateitem), &
-              ScenarioDescription('vertical_regridding',     'parent.yaml',      check_name, check_stateitem) &
-              ! ScenarioDescription('vertical_regridding_2',   'parent.yaml',      check_name, check_stateitem) &
-      ]
+              ScenarioDescription('vertical_regridding',     'parent.yaml',      check_name, check_stateitem)  &
+#ifndef __GFORTRAN__
+              , &
+              ScenarioDescription('vertical_regridding_2',   'parent.yaml',      check_name, check_stateitem), &
+              ScenarioDescription('vertical_regridding_3',   'parent.yaml',      check_name, check_stateitem)  &
+#endif
+              ]
       end function add_params
 
    end function get_parameters

generic3g/tests/scenarios/vertical_regridding_3/A.yaml

@@ -0,0 +1,28 @@
+mapl:
+
+  geometry:
+    esmf_geom:
+      class: latlon
+      im_world: 12
+      jm_world: 13
+      pole: PC
+      dateline: DC
+    vertical_grid:
+      class: model
+      short_name: PL
+      units: hPa
+      num_levels: 4
+
+  states:
+    import: {}
+    export:
+      PL:
+        standard_name: air_pressure_a
+        units: hPa
+        default_value: 17.
+        vertical_dim_spec: center
+      E_A:
+        standard_name: temperature_a
+        units: K
+        default_value: 17.
+        vertical_dim_spec: center