Diagnostics downsampling, fix the issue for symmetric case

- There was a unnecessary check for symmetric case
NOAA-GFDL · Oct 26, 2018 · 16fef4f · 16fef4f
1 parent fe0c041
commit 16fef4f
Show file tree

Hide file tree

Showing 2 changed files with 24 additions and 18 deletions.
diff --git a/src/core/MOM_grid.F90 b/src/core/MOM_grid.F90
@@ -21,7 +21,7 @@ module MOM_grid
   type(MOM_domain_type), pointer :: Domain => NULL() !< Ocean model domain
   type(MOM_domain_type), pointer :: Domain_aux => NULL() !< A non-symmetric auxiliary domain type.
   type(hor_index_type) :: HI !< Horizontal index ranges
-  type(hor_index_type) :: HId2 !< Horizontal index ranges for level-2-downsampling 
+  type(hor_index_type) :: HId2 !< Horizontal index ranges for level-2-downsampling
 
   integer :: isc !< The start i-index of cell centers within the computational domain
   integer :: iec !< The end i-index of cell centers within the computational domain
@@ -348,8 +348,9 @@ subroutine MOM_grid_init(G, param_file, HI, global_indexing, bathymetry_at_vel)
   call get_domain_extent_dsamp2(G%Domain, G%HId2%isc, G%HId2%iec, G%HId2%jsc, G%HId2%jec,&
                                           G%HId2%isd, G%HId2%ied, G%HId2%jsd, G%HId2%jed,&
                                           G%HId2%isg, G%HId2%ieg, G%HId2%jsg, G%HId2%jeg)
-  
+
   ! Set array sizes for fields that are discretized at tracer cell boundaries.
+  G%HId2%IegB = G%HId2%ieg ; G%HId2%JegB = G%HId2%jeg
   G%HId2%IsdB = G%HId2%isd ; G%HId2%JsdB = G%HId2%jsd
   if (G%symmetric) then
     G%HId2%IsdB = G%HId2%isd-1 ; G%HId2%JsdB = G%HId2%jsd-1

diff --git a/src/framework/MOM_diag_mediator.F90 b/src/framework/MOM_diag_mediator.F90
@@ -517,36 +517,41 @@ subroutine set_axes_info_dsamp(G, GV, param_file, diag_cs, id_zl_native, id_zi_n
   integer :: i, j, k, nz, dl
   real, dimension(:), pointer :: gridLonT_dsamp =>NULL()
   real, dimension(:), pointer :: gridLatT_dsamp =>NULL()
+  real, dimension(:), pointer :: gridLonB_dsamp =>NULL()
+  real, dimension(:), pointer :: gridLatB_dsamp =>NULL()
 
   id_zl = id_zl_native ; id_zi = id_zi_native
   !Axes group for native downsampled diagnostics
   do dl=2,MAX_DSAMP_LEV
      if(dl .ne. 2) call MOM_error(FATAL, "set_axes_info_dsamp: Downsample level other than 2 is not supported yet!")
      allocate(gridLonT_dsamp(diag_cs%dsamp(dl)%isg:diag_cs%dsamp(dl)%ieg))
      allocate(gridLatT_dsamp(diag_cs%dsamp(dl)%jsg:diag_cs%dsamp(dl)%jeg))
-
      do i=diag_cs%dsamp(dl)%isg,diag_cs%dsamp(dl)%ieg;  gridLonT_dsamp(i) = G%gridLonT(G%isg+dl*i-2); enddo
      do j=diag_cs%dsamp(dl)%jsg,diag_cs%dsamp(dl)%jeg;  gridLatT_dsamp(j) = G%gridLatT(G%jsg+dl*j-2); enddo
-
-     if (G%symmetric) then
-        call MOM_error(FATAL, "set_axes_info_dsamp: Downsample of symmetric case is not supported yet!")
-     !   id_xq = diag_axis_init('xq', gridLonB_dsamp(G%isgB:G%iegB), G%x_axis_units, 'x', &
-     !             'q point nominal longitude', Domain2=G%Domain%mpp_domain_d2)
-     !   id_yq = diag_axis_init('yq', gridLatB_dsamp(G%jsgB:G%jegB), G%y_axis_units, 'y', &
-     !             'q point nominal latitude', Domain2=G%Domain%mpp_domain_d2)
-     else
-        id_xq = diag_axis_init('xq', gridLonT_dsamp, G%x_axis_units, 'x', &
-                  'q point nominal longitude', Domain2=G%Domain%mpp_domain_d2)
-        id_yq = diag_axis_init('yq', gridLatT_dsamp, G%y_axis_units, 'y', &
-                  'q point nominal latitude', Domain2=G%Domain%mpp_domain_d2)
-     endif
+     allocate(gridLonB_dsamp(diag_cs%dsamp(dl)%isg:diag_cs%dsamp(dl)%ieg))
+     allocate(gridLatB_dsamp(diag_cs%dsamp(dl)%jsg:diag_cs%dsamp(dl)%jeg))
+     do i=diag_cs%dsamp(dl)%isg,diag_cs%dsamp(dl)%ieg;  gridLonB_dsamp(i) = G%gridLonB(G%isg+dl*i-2); enddo
+     do j=diag_cs%dsamp(dl)%jsg,diag_cs%dsamp(dl)%jeg;  gridLatB_dsamp(j) = G%gridLatB(G%jsg+dl*j-2); enddo
+
+!I don't see a need for this since isgB=isg and iegB=ieg
+!     if (G%symmetric) then
+!        id_xq = diag_axis_init('xq', gridLonB_dsamp(G%isgB:G%iegB), G%x_axis_units, 'x', &
+!                  'q point nominal longitude', Domain2=G%Domain%mpp_domain_d2)
+!        id_yq = diag_axis_init('yq', gridLatB_dsamp(G%jsgB:G%jegB), G%y_axis_units, 'y', &
+!                  'q point nominal latitude', Domain2=G%Domain%mpp_domain_d2)
+!     else
+     id_xq = diag_axis_init('xq', gridLonB_dsamp, G%x_axis_units, 'x', &
+               'q point nominal longitude', Domain2=G%Domain%mpp_domain_d2)
+     id_yq = diag_axis_init('yq', gridLatB_dsamp, G%y_axis_units, 'y', &
+               'q point nominal latitude', Domain2=G%Domain%mpp_domain_d2)
+!     endif
      id_xh = diag_axis_init('xh', gridLonT_dsamp, G%x_axis_units, 'x', &
           'h point nominal longitude', Domain2=G%Domain%mpp_domain_d2)
      id_yh = diag_axis_init('yh', gridLatT_dsamp, G%y_axis_units, 'y', &
           'h point nominal latitude', Domain2=G%Domain%mpp_domain_d2)
 
-     deallocate(gridLonT_dsamp)
-     deallocate(gridLatT_dsamp)
+     deallocate(gridLonT_dsamp,gridLatT_dsamp)
+     deallocate(gridLonB_dsamp,gridLatB_dsamp)
 
      ! Axis groupings for the model layers
      call define_axes_group_dsamp(diag_cs, (/ id_xh, id_yh, id_zL /), diag_cs%dsamp(dl)%axesTL, dl, &