diff --git a/src/framework/MOM_horizontal_regridding.F90 b/src/framework/MOM_horizontal_regridding.F90 index bbb5ae0e15..7a0ace9279 100644 --- a/src/framework/MOM_horizontal_regridding.F90 +++ b/src/framework/MOM_horizontal_regridding.F90 @@ -250,14 +250,13 @@ subroutine fill_miss_2d(aout, good, fill, prev, G, acrit, num_pass, relc, debug, end subroutine fill_miss_2d !> Extrapolate and interpolate from a file record -subroutine horiz_interp_and_extrap_tracer_record(filename, varnam, conversion, recnum, G, tr_z, mask_z, & - z_in, z_edges_in, missing_value, reentrant_x, tripolar_n, & +subroutine horiz_interp_and_extrap_tracer_record(filename, varnam, recnum, G, tr_z, mask_z, & + z_in, z_edges_in, missing_value, scale, & homogenize, m_to_Z, answers_2018, ongrid, tr_iter_tol, answer_date) character(len=*), intent(in) :: filename !< Path to file containing tracer to be !! interpolated. character(len=*), intent(in) :: varnam !< Name of tracer in file. - real, intent(in) :: conversion !< Conversion factor for tracer [CU conc-1 ~> 1] integer, intent(in) :: recnum !< Record number of tracer to be read. type(ocean_grid_type), intent(inout) :: G !< Grid object real, allocatable, dimension(:,:,:), intent(out) :: tr_z @@ -271,10 +270,10 @@ subroutine horiz_interp_and_extrap_tracer_record(filename, varnam, conversion, real, allocatable, dimension(:), intent(out) :: z_edges_in !< Cell grid edge values for input data [Z ~> m] real, intent(out) :: missing_value !< The missing value in the returned array, scaled - !! with conversion to avoid accidentally having valid - !! values match missing values [CU ~> conc] - logical, intent(in) :: reentrant_x !< If true, this grid is reentrant in the x-direction - logical, intent(in) :: tripolar_n !< If true, this is a northern tripolar grid + !! to avoid accidentally having valid values match + !! missing values [CU ~> conc] + real, intent(in) :: scale !< Scaling factor for tracer into the internal + !! units of the model [CU conc-1 ~> 1] logical, optional, intent(in) :: homogenize !< If present and true, horizontally homogenize data !! to produce perfectly "flat" initial conditions real, optional, intent(in) :: m_to_Z !< A conversion factor from meters to the units @@ -311,7 +310,7 @@ subroutine horiz_interp_and_extrap_tracer_record(filename, varnam, conversion, real, dimension(:), allocatable :: lon_in ! The longitudes in the input file [degreesE] then [radians] real, dimension(:), allocatable :: lat_in ! The latitudes in the input file [degreesN] then [radians] real, dimension(:), allocatable :: lat_inp ! The input file latitudes expanded to the pole [degreesN] then [radians] - real :: max_lat ! The maximum latitude on the input grid [degreeN] + real :: max_lat ! The maximum latitude on the input grid [degreesN] real :: pole ! The sum of tracer values at the pole [conc] real :: max_depth ! The maximum depth of the ocean [Z ~> m] real :: npole ! The number of points contributing to the pole value [nondim] @@ -330,7 +329,7 @@ subroutine horiz_interp_and_extrap_tracer_record(filename, varnam, conversion, integer :: isd, ied, jsd, jed ! data domain indices integer :: id_clock_read logical :: debug=.false. - real :: I_scale ! The inverse of the conversion factor for diagnostic output [conc CU-1 ~> 1] + real :: I_scale ! The inverse of the scale factor for diagnostic output [conc CU-1 ~> 1] real :: dtr_iter_stop ! The tolerance for changes in tracer concentrations between smoothing ! iterations that determines when to stop iterating [CU ~> conc] real :: npoints ! The number of points in an average [nondim] @@ -355,10 +354,10 @@ subroutine horiz_interp_and_extrap_tracer_record(filename, varnam, conversion, is_ongrid = .false. if (present(ongrid)) is_ongrid = ongrid - dtr_iter_stop = 1.0e-3*conversion + dtr_iter_stop = 1.0e-3*scale if (present(tr_iter_tol)) dtr_iter_stop = tr_iter_tol - I_scale = 1.0 / conversion + I_scale = 1.0 / scale PI_180 = atan(1.0)/45. @@ -371,6 +370,9 @@ subroutine horiz_interp_and_extrap_tracer_record(filename, varnam, conversion, call cpu_clock_begin(id_clock_read) + ! A note by MJH copied from elsewhere suggests that this code may be using the model connectivity + ! (e.g., reentrant or tripolar) but should use the dataset's connectivity instead. + call get_var_axes_info(trim(filename), trim(varnam), axes_info) if (allocated(z_in)) deallocate(z_in) @@ -418,7 +420,7 @@ subroutine horiz_interp_and_extrap_tracer_record(filename, varnam, conversion, if (.not. found_attr) call MOM_error(FATAL, & "error finding missing value for " // trim(varnam) // & " in file " // trim(filename) // " in hinterp_extrap") - missing_value = conversion * missing_val_in + missing_value = scale * missing_val_in call read_attribute(trim(filename), "scale_factor", scale_factor, & varname=trim(varnam), found=found_attr) @@ -471,7 +473,7 @@ subroutine horiz_interp_and_extrap_tracer_record(filename, varnam, conversion, do i=is,ie if (abs(tr_in(i,j)-missing_val_in) > abs(roundoff*missing_val_in)) then mask_in(i,j) = 1.0 - tr_in(i,j) = (tr_in(i,j)*scale_factor+add_offset) * conversion + tr_in(i,j) = (tr_in(i,j)*scale_factor+add_offset) * scale else tr_in(i,j) = missing_value endif @@ -511,7 +513,7 @@ subroutine horiz_interp_and_extrap_tracer_record(filename, varnam, conversion, do j=1,jdp ; do i=1,id if (abs(tr_inp(i,j)-missing_val_in) > abs(roundoff*missing_val_in)) then mask_in(i,j) = 1.0 - tr_inp(i,j) = (tr_inp(i,j)*scale_factor+add_offset) * conversion + tr_inp(i,j) = (tr_inp(i,j)*scale_factor+add_offset) * scale else tr_inp(i,j) = missing_value endif @@ -596,18 +598,17 @@ subroutine horiz_interp_and_extrap_tracer_record(filename, varnam, conversion, end subroutine horiz_interp_and_extrap_tracer_record !> Extrapolate and interpolate using a FMS time interpolation handle -subroutine horiz_interp_and_extrap_tracer_fms_id(fms_id, Time, conversion, G, tr_z, mask_z, & - z_in, z_edges_in, missing_value, reentrant_x, tripolar_n, & +subroutine horiz_interp_and_extrap_tracer_fms_id(fms_id, Time, G, tr_z, mask_z, & + z_in, z_edges_in, missing_value, scale, & homogenize, spongeOngrid, m_to_Z, & answers_2018, tr_iter_tol, answer_date) integer, intent(in) :: fms_id !< A unique id used by the FMS time interpolator type(time_type), intent(in) :: Time !< A FMS time type - real, intent(in) :: conversion !< Conversion factor for tracer. type(ocean_grid_type), intent(inout) :: G !< Grid object real, allocatable, dimension(:,:,:), intent(out) :: tr_z !< Allocatable tracer array on the horizontal - !! model grid and input-file vertical levels. [CU ~> conc] + !! model grid and input-file vertical levels [CU ~> conc] real, allocatable, dimension(:,:,:), intent(out) :: mask_z !< Allocatable tracer mask array on the horizontal !! model grid and input-file vertical levels [nondim] @@ -616,10 +617,10 @@ subroutine horiz_interp_and_extrap_tracer_fms_id(fms_id, Time, conversion, G, t real, allocatable, dimension(:), intent(out) :: z_edges_in !< Cell grid edge values for input data [Z ~> m] real, intent(out) :: missing_value !< The missing value in the returned array, scaled - !! with conversion to avoid accidentally having valid - !! values match missing values [CU ~> conc] - logical, intent(in) :: reentrant_x !< If true, this grid is reentrant in the x-direction - logical, intent(in) :: tripolar_n !< If true, this is a northern tripolar grid + !! to avoid accidentally having valid values match + !! missing values [CU ~> conc] + real, intent(in) :: scale !< Scaling factor for tracer into the internal + !! units of the model [CU conc-1 ~> 1] logical, optional, intent(in) :: homogenize !< If present and true, horizontally homogenize data !! to produce perfectly "flat" initial conditions logical, optional, intent(in) :: spongeOngrid !< If present and true, the sponge data are on the model grid @@ -655,7 +656,7 @@ subroutine horiz_interp_and_extrap_tracer_fms_id(fms_id, Time, conversion, G, t real, dimension(:), allocatable :: lon_in ! The longitudes in the input file [degreesE] then [radians] real, dimension(:), allocatable :: lat_in ! The latitudes in the input file [degreesN] then [radians] real, dimension(:), allocatable :: lat_inp ! The input file latitudes expanded to the pole [degreesN] then [radians] - real :: max_lat ! The maximum latitude on the input grid [degreeN] + real :: max_lat ! The maximum latitude on the input grid [degreesN] real :: pole ! The sum of tracer values at the pole [conc] real :: max_depth ! The maximum depth of the ocean [Z ~> m] real :: npole ! The number of points contributing to the pole value [nondim] @@ -673,7 +674,7 @@ subroutine horiz_interp_and_extrap_tracer_fms_id(fms_id, Time, conversion, G, t logical :: debug=.false. logical :: is_ongrid integer :: ans_date ! The vintage of the expressions and order of arithmetic to use - real :: I_scale ! The inverse of the conversion factor for diagnostic output [conc CU-1 ~> 1] + real :: I_scale ! The inverse of the scale factor for diagnostic output [conc CU-1 ~> 1] real :: dtr_iter_stop ! The tolerance for changes in tracer concentrations between smoothing ! iterations that determines when to stop iterating [CU ~> conc] real :: npoints ! The number of points in an average [nondim] @@ -699,10 +700,10 @@ subroutine horiz_interp_and_extrap_tracer_fms_id(fms_id, Time, conversion, G, t id_clock_read = cpu_clock_id('(Initialize tracer from Z) read', grain=CLOCK_LOOP) - dtr_iter_stop = 1.0e-3*conversion + dtr_iter_stop = 1.0e-3*scale if (present(tr_iter_tol)) dtr_iter_stop = tr_iter_tol - I_scale = 1.0 / conversion + I_scale = 1.0 / scale PI_180 = atan(1.0)/45. @@ -716,7 +717,7 @@ subroutine horiz_interp_and_extrap_tracer_fms_id(fms_id, Time, conversion, G, t call cpu_clock_begin(id_clock_read) call get_external_field_info(fms_id, size=fld_sz, axes=axes_data, missing=missing_val_in) - missing_value = conversion*missing_val_in + missing_value = scale*missing_val_in verbosity = MOM_get_verbosity() @@ -823,7 +824,7 @@ subroutine horiz_interp_and_extrap_tracer_fms_id(fms_id, Time, conversion, G, t do j=1,jdp ; do i=1,id if (abs(tr_inp(i,j)-missing_val_in) > abs(roundoff*missing_val_in)) then mask_in(i,j) = 1.0 - tr_inp(i,j) = tr_inp(i,j) * conversion + tr_inp(i,j) = tr_inp(i,j) * scale else tr_inp(i,j) = missing_value endif @@ -909,7 +910,7 @@ subroutine horiz_interp_and_extrap_tracer_fms_id(fms_id, Time, conversion, G, t do k=1,kd do j=js,je do i=is,ie - tr_z(i,j,k) = data_in(i,j,k) * conversion + tr_z(i,j,k) = data_in(i,j,k) * scale if (ans_date >= 20190101) mask_z(i,j,k) = 1. if (abs(tr_z(i,j,k)-missing_value) < abs(roundoff*missing_value)) mask_z(i,j,k) = 0. enddo diff --git a/src/initialization/MOM_state_initialization.F90 b/src/initialization/MOM_state_initialization.F90 index 62438cbc7f..c5b57f3d57 100644 --- a/src/initialization/MOM_state_initialization.F90 +++ b/src/initialization/MOM_state_initialization.F90 @@ -2483,8 +2483,12 @@ subroutine MOM_temp_salt_initialize_from_Z(h, tv, depth_tot, G, GV, US, PF, just real :: PI_180 ! for conversion from degrees to radians [radian degree-1] real :: Hmix_default ! The default initial mixed layer depth [Z ~> m]. real :: Hmix_depth ! The mixed layer depth in the initial condition [Z ~> m]. - real :: missing_value_temp ! The missing value in the input temperature field - real :: missing_value_salt ! The missing value in the input salinity field + real :: missing_value_temp ! The missing value in the input temperature field [C ~> degC] + real :: missing_value_salt ! The missing value in the input salinity field [S ~> ppt] + real :: tol_temp ! The tolerance for changes in temperature during the horizontal + ! interpolation from an input dataset [C ~> degC] + real :: tol_sal ! The tolerance for changes in salinity during the horizontal + ! interpolation from an input dataset [S ~> ppt] logical :: correct_thickness real :: h_tolerance ! A parameter that controls the tolerance when adjusting the ! thickness to fit the bathymetry [Z ~> m]. @@ -2494,19 +2498,19 @@ subroutine MOM_temp_salt_initialize_from_Z(h, tv, depth_tot, G, GV, US, PF, just logical :: adjust_temperature = .true. ! fit t/s to target densities real :: temp_land_fill ! A temperature value to use for land points [C ~> degC] real :: salt_land_fill ! A salinity value to use for land points [C ~> degC] - logical :: reentrant_x, tripolar_n ! data arrays - real, dimension(:), allocatable :: z_edges_in, z_in ! Interface heights [Z ~> m] - real, dimension(:), allocatable :: Rb ! Interface densities [R ~> kg m-3] + real, dimension(:), allocatable :: z_edges_in ! Input data interface heights or depths [Z ~> m] + real, dimension(:), allocatable :: z_in ! Input data cell heights or depths [Z ~> m] + real, dimension(:), allocatable :: Rb ! Interface densities [R ~> kg m-3] real, dimension(:,:,:), allocatable, target :: temp_z ! Input temperatures [C ~> degC] real, dimension(:,:,:), allocatable, target :: salt_z ! Input salinities [S ~> ppt] real, dimension(:,:,:), allocatable, target :: mask_z ! 1 for valid data points [nondim] - real, dimension(:,:,:), allocatable :: rho_z ! Densities in Z-space [R ~> kg m-3] + real, dimension(:,:,:), allocatable :: rho_z ! Densities in Z-space [R ~> kg m-3] real, dimension(SZI_(G),SZJ_(G),SZK_(GV)+1) :: zi ! Interface heights [Z ~> m]. - real, dimension(SZI_(G),SZJ_(G)) :: Z_bottom ! The (usually negative) height of the seafloor - ! relative to the surface [Z ~> m]. - integer, dimension(SZI_(G),SZJ_(G)) :: nlevs ! The number of levels in each column with valid data + real, dimension(SZI_(G),SZJ_(G)) :: Z_bottom ! The (usually negative) height of the seafloor + ! relative to the surface [Z ~> m]. + integer, dimension(SZI_(G),SZJ_(G)) :: nlevs ! The number of levels in each column with valid data real, dimension(SZI_(G)) :: press ! Pressures [R L2 T-2 ~> Pa]. ! Local variables for ALE remapping @@ -2569,9 +2573,6 @@ subroutine MOM_temp_salt_initialize_from_Z(h, tv, depth_tot, G, GV, US, PF, just eos => tv%eqn_of_state - reentrant_x = .false. ; call get_param(PF, mdl, "REENTRANT_X", reentrant_x, default=.true.) - tripolar_n = .false. ; call get_param(PF, mdl, "TRIPOLAR_N", tripolar_n, default=.false.) - call get_param(PF, mdl, "TEMP_SALT_Z_INIT_FILE", filename, & "The name of the z-space input file used to initialize "//& "temperatures (T) and salinities (S). If T and S are not "//& @@ -2701,6 +2702,8 @@ subroutine MOM_temp_salt_initialize_from_Z(h, tv, depth_tot, G, GV, US, PF, just !### These hard-coded constants should be made into runtime parameters temp_land_fill = 0.0*US%degC_to_C salt_land_fill = 35.0*US%ppt_to_S + tol_temp = 1.0e-3*US%degC_to_C + tol_sal = 1.0e-3*US%ppt_to_S eps_z = GV%Angstrom_Z eps_rho = 1.0e-10*US%kg_m3_to_R @@ -2720,15 +2723,15 @@ subroutine MOM_temp_salt_initialize_from_Z(h, tv, depth_tot, G, GV, US, PF, just ! to the North/South Pole past the limits of the input data, they are extrapolated using the average ! value at the northernmost/southernmost latitude. - call horiz_interp_and_extrap_tracer(tfilename, potemp_var, US%degC_to_C, 1, & - G, temp_z, mask_z, z_in, z_edges_in, missing_value_temp, reentrant_x, & - tripolar_n, homogenize, m_to_Z=US%m_to_Z, answer_date=hor_regrid_answer_date, & - ongrid=pre_gridded, tr_iter_tol=1.0e-3*US%degC_to_C) + call horiz_interp_and_extrap_tracer(tfilename, potemp_var, 1, & + G, temp_z, mask_z, z_in, z_edges_in, missing_value_temp, & + scale=US%degC_to_C, homogenize=homogenize, m_to_Z=US%m_to_Z, & + answer_date=hor_regrid_answer_date, ongrid=pre_gridded, tr_iter_tol=tol_temp) - call horiz_interp_and_extrap_tracer(sfilename, salin_var, US%ppt_to_S, 1, & - G, salt_z, mask_z, z_in, z_edges_in, missing_value_salt, reentrant_x, & - tripolar_n, homogenize, m_to_Z=US%m_to_Z, answer_date=hor_regrid_answer_date, & - ongrid=pre_gridded, tr_iter_tol=1.0e-3*US%ppt_to_S) + call horiz_interp_and_extrap_tracer(sfilename, salin_var, 1, & + G, salt_z, mask_z, z_in, z_edges_in, missing_value_salt, & + scale=US%ppt_to_S, homogenize=homogenize, m_to_Z=US%m_to_Z, & + answer_date=hor_regrid_answer_date, ongrid=pre_gridded, tr_iter_tol=tol_sal) kd = size(z_in,1) diff --git a/src/initialization/MOM_tracer_initialization_from_Z.F90 b/src/initialization/MOM_tracer_initialization_from_Z.F90 index 04c03a5b43..7c62ea496e 100644 --- a/src/initialization/MOM_tracer_initialization_from_Z.F90 +++ b/src/initialization/MOM_tracer_initialization_from_Z.F90 @@ -41,11 +41,12 @@ subroutine MOM_initialize_tracer_from_Z(h, tr, G, GV, US, PF, src_file, src_var_ type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), & intent(in) :: h !< Layer thickness [H ~> m or kg m-2]. - real, dimension(:,:,:), pointer :: tr !< Pointer to array to be initialized + real, dimension(:,:,:), pointer :: tr !< Pointer to array to be initialized [CU ~> conc] type(param_file_type), intent(in) :: PF !< parameter file character(len=*), intent(in) :: src_file !< source filename character(len=*), intent(in) :: src_var_nam !< variable name in file - real, optional, intent(in) :: src_var_unit_conversion !< optional multiplicative unit conversion + real, optional, intent(in) :: src_var_unit_conversion !< optional multiplicative unit conversion, + !! often used for rescaling into model units [CU conc-1 ~> 1] integer, optional, intent(in) :: src_var_record !< record to read for multiple time-level files logical, optional, intent(in) :: homogenize !< optionally homogenize to mean value logical, optional, intent(in) :: useALEremapping !< to remap or not (optional) @@ -53,11 +54,11 @@ subroutine MOM_initialize_tracer_from_Z(h, tr, G, GV, US, PF, src_file, src_var_ character(len=*), optional, intent(in) :: src_var_gridspec !< Source variable name in a gridspec file. !! This is not implemented yet. ! Local variables - real :: land_fill = 0.0 - real :: convert + real :: land_fill = 0.0 ! A value to use to replace missing values [CU ~> conc] + real :: convert ! A conversion factor into the model's internal units [CU conc-1 ~> 1] integer :: recnum character(len=64) :: remapScheme - logical :: homog,useALE + logical :: homog, useALE ! This include declares and sets the variable "version". # include "version_variable.h" @@ -66,8 +67,12 @@ subroutine MOM_initialize_tracer_from_Z(h, tr, G, GV, US, PF, src_file, src_var_ integer :: is, ie, js, je, nz ! compute domain indices integer :: isd, ied, jsd, jed ! data domain indices integer :: i, j, k, kd - real, allocatable, dimension(:,:,:), target :: tr_z, mask_z - real, allocatable, dimension(:), target :: z_edges_in, z_in + real, allocatable, dimension(:,:,:), target :: tr_z ! Tracer array on the horizontal model grid + ! and input-file vertical levels [CU ~> conc] + real, allocatable, dimension(:,:,:), target :: mask_z ! Missing value mask on the horizontal model grid + ! and input-file vertical levels [nondim] + real, allocatable, dimension(:), target :: z_edges_in ! Cell edge depths for input data [Z ~> m] + real, allocatable, dimension(:), target :: z_in ! Cell center depths for input data [Z ~> m] ! Local variables for ALE remapping real, dimension(:,:,:), allocatable :: hSrc ! Source thicknesses [H ~> m or kg m-2]. @@ -75,8 +80,8 @@ subroutine MOM_initialize_tracer_from_Z(h, tr, G, GV, US, PF, src_file, src_var_ real :: zTopOfCell, zBottomOfCell, z_bathy ! Heights [Z ~> m]. type(remapping_CS) :: remapCS ! Remapping parameters and work arrays - real :: missing_value - integer :: nPoints + real :: missing_value ! A value indicating that there is no valid input data at this point [CU ~> conc] + integer :: nPoints ! The number of valid input data points in a column integer :: id_clock_routine, id_clock_ALE integer :: default_answer_date ! The default setting for the various ANSWER_DATE flags. logical :: default_2018_answers ! The default setting for the various 2018_ANSWERS flags. @@ -94,7 +99,6 @@ subroutine MOM_initialize_tracer_from_Z(h, tr, G, GV, US, PF, src_file, src_var_ ! for horizontal regridding. Values below 20190101 recover the ! answers from 2018, while higher values use expressions that have ! been rearranged for rotational invariance. - logical :: reentrant_x, tripolar_n id_clock_routine = cpu_clock_id('(Initialize tracer from Z)', grain=CLOCK_ROUTINE) id_clock_ALE = cpu_clock_id('(Initialize tracer from Z) ALE', grain=CLOCK_LOOP) @@ -153,22 +157,17 @@ subroutine MOM_initialize_tracer_from_Z(h, tr, G, GV, US, PF, src_file, src_var_ "If both HOR_REGRID_2018_ANSWERS and HOR_REGRID_ANSWER_DATE are specified, the "//& "latter takes precedence.", default=default_hor_reg_ans_date) - ! These are model grid properties, but being applied to the data grid for now. - ! need to revisit this (mjh) - reentrant_x = .false. ; call get_param(PF, mdl, "REENTRANT_X", reentrant_x,default=.true.) - tripolar_n = .false. ; call get_param(PF, mdl, "TRIPOLAR_N", tripolar_n, default=.false.) - if (PRESENT(homogenize)) homog=homogenize if (PRESENT(useALEremapping)) useALE=useALEremapping if (PRESENT(remappingScheme)) remapScheme=remappingScheme - recnum=1 + recnum = 1 if (PRESENT(src_var_record)) recnum = src_var_record - convert=1.0 + convert = 1.0 if (PRESENT(src_var_unit_conversion)) convert = src_var_unit_conversion - call horiz_interp_and_extrap_tracer(src_file, src_var_nam, convert, recnum, & - G, tr_z, mask_z, z_in, z_edges_in, missing_value, reentrant_x, tripolar_n, & - homog, m_to_Z=US%m_to_Z, answer_date=hor_regrid_answer_date) + call horiz_interp_and_extrap_tracer(src_file, src_var_nam, recnum, & + G, tr_z, mask_z, z_in, z_edges_in, missing_value, & + scale=convert, homogenize=homog, m_to_Z=US%m_to_Z, answer_date=hor_regrid_answer_date) kd = size(z_edges_in,1)-1 call pass_var(tr_z,G%Domain) @@ -221,7 +220,7 @@ subroutine MOM_initialize_tracer_from_Z(h, tr, G, GV, US, PF, src_file, src_var_ ! Fill land values do k=1,nz ; do j=js,je ; do i=is,ie if (tr(i,j,k) == missing_value) then - tr(i,j,k)=land_fill + tr(i,j,k) = land_fill endif enddo ; enddo ; enddo diff --git a/src/ocean_data_assim/MOM_oda_driver.F90 b/src/ocean_data_assim/MOM_oda_driver.F90 index 2a1a96168a..4ad11592f9 100644 --- a/src/ocean_data_assim/MOM_oda_driver.F90 +++ b/src/ocean_data_assim/MOM_oda_driver.F90 @@ -559,19 +559,22 @@ subroutine get_bias_correction_tracer(Time, US, CS) type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(ODA_CS), pointer :: CS !< ocean DA control structure - integer :: i,j,k + ! Local variables real, allocatable, dimension(:,:,:) :: T_bias ! Temperature biases [C ~> degC] real, allocatable, dimension(:,:,:) :: S_bias ! Salinity biases [C ~> degC] - real, allocatable, dimension(:,:,:) :: mask_z - real, allocatable, dimension(:), target :: z_in, z_edges_in - real :: missing_value - integer,dimension(3) :: fld_sz + real, allocatable, dimension(:,:,:) :: mask_z ! Missing value mask on the horizontal model grid + ! and input-file vertical levels [nondim] + real, allocatable, dimension(:), target :: z_in ! Cell center depths for input data [Z ~> m] + real, allocatable, dimension(:), target :: z_edges_in ! Cell edge depths for input data [Z ~> m] + real :: missing_value ! A value indicating that there is no valid input data at this point [CU ~> conc] + integer, dimension(3) :: fld_sz + integer :: i,j,k call cpu_clock_begin(id_clock_bias_adjustment) - call horiz_interp_and_extrap_tracer(CS%INC_CS%T_id, Time, US%degC_to_C, CS%G, T_bias, & - mask_z, z_in, z_edges_in, missing_value, .true., .false., .false., .true.) - call horiz_interp_and_extrap_tracer(CS%INC_CS%S_id, Time, US%ppt_to_S, CS%G, S_bias, & - mask_z, z_in, z_edges_in, missing_value, .true., .false., .false., .true.) + call horiz_interp_and_extrap_tracer(CS%INC_CS%T_id, Time, CS%G, T_bias, & + mask_z, z_in, z_edges_in, missing_value, scale=US%degC_to_C, spongeOngrid=.true.) + call horiz_interp_and_extrap_tracer(CS%INC_CS%S_id, Time, CS%G, S_bias, & + mask_z, z_in, z_edges_in, missing_value, scale=US%ppt_to_S, spongeOngrid=.true.) ! This should be replaced to use mask_z instead of the following lines ! which are intended to zero land values using an arbitrary limit. diff --git a/src/parameterizations/vertical/MOM_ALE_sponge.F90 b/src/parameterizations/vertical/MOM_ALE_sponge.F90 index 9f5241bb9a..a121b05d30 100644 --- a/src/parameterizations/vertical/MOM_ALE_sponge.F90 +++ b/src/parameterizations/vertical/MOM_ALE_sponge.F90 @@ -133,9 +133,6 @@ module MOM_ALE_sponge logical :: time_varying_sponges !< True if using newer sponge code logical :: spongeDataOngrid !< True if the sponge data are on the model horizontal grid - logical :: reentrant_x !< grid is reentrant in the x direction - logical :: tripolar_N !< grid is folded at its north edge - !>@{ Diagnostic IDs integer, dimension(MAX_FIELDS_) :: id_sp_tendency !< Diagnostic ids for tracers !! tendency due to sponges @@ -257,11 +254,6 @@ subroutine initialize_ALE_sponge_fixed(Iresttime, G, GV, param_file, CS, data_h, "while later versions add parentheses for rotational symmetry. "//& "If both HOR_REGRID_2018_ANSWERS and HOR_REGRID_ANSWER_DATE are specified, the "//& "latter takes precedence.", default=default_hor_reg_ans_date) - call get_param(param_file, mdl, "REENTRANT_X", CS%reentrant_x, & - "If true, the domain is zonally reentrant.", default=.true.) - call get_param(param_file, mdl, "TRIPOLAR_N", CS%tripolar_N, & - "Use tripolar connectivity at the northern edge of the "//& - "domain. With TRIPOLAR_N, NIGLOBAL must be even.", default=.false.) CS%time_varying_sponges = .false. CS%nz = GV%ke @@ -551,11 +543,6 @@ subroutine initialize_ALE_sponge_varying(Iresttime, G, GV, param_file, CS, Irest "When defined, the incoming sponge data are "//& "assumed to be on the model grid " , & default=.false.) - call get_param(param_file, mdl, "REENTRANT_X", CS%reentrant_x, & - "If true, the domain is zonally reentrant.", default=.true.) - call get_param(param_file, mdl, "TRIPOLAR_N", CS%tripolar_N, & - "Use tripolar connectivity at the northern edge of the "//& - "domain. With TRIPOLAR_N, NIGLOBAL must be even.", default=.false.) CS%time_varying_sponges = .true. CS%nz = GV%ke @@ -985,10 +972,10 @@ subroutine apply_ALE_sponge(h, dt, G, GV, US, CS, Time) if (CS%time_varying_sponges) then do m=1,CS%fldno nz_data = CS%Ref_val(m)%nz_data - call horiz_interp_and_extrap_tracer(CS%Ref_val(m)%id, Time, CS%Ref_val(m)%scale, G, sp_val, & - mask_z, z_in, z_edges_in, missing_value, CS%reentrant_x, CS%tripolar_N, .false., & - spongeOnGrid=CS%SpongeDataOngrid, m_to_Z=US%m_to_Z, & - answer_date=CS%hor_regrid_answer_date) + call horiz_interp_and_extrap_tracer(CS%Ref_val(m)%id, Time, G, sp_val, & + mask_z, z_in, z_edges_in, missing_value, & + scale=CS%Ref_val(m)%scale, spongeOnGrid=CS%SpongeDataOngrid, m_to_Z=US%m_to_Z, & + answer_date=CS%hor_regrid_answer_date) allocate( hsrc(nz_data) ) allocate( tmpT1d(nz_data) ) do c=1,CS%num_col @@ -1069,10 +1056,10 @@ subroutine apply_ALE_sponge(h, dt, G, GV, US, CS, Time) if (CS%time_varying_sponges) then nz_data = CS%Ref_val_u%nz_data ! Interpolate from the external horizontal grid and in time - call horiz_interp_and_extrap_tracer(CS%Ref_val_u%id, Time, CS%Ref_val_u%scale, G, sp_val, & - mask_z, z_in, z_edges_in, missing_value, CS%reentrant_x, CS%tripolar_N, .false., & - spongeOnGrid=CS%SpongeDataOngrid, m_to_Z=US%m_to_Z,& - answer_date=CS%hor_regrid_answer_date) + call horiz_interp_and_extrap_tracer(CS%Ref_val_u%id, Time, G, sp_val, & + mask_z, z_in, z_edges_in, missing_value, & + scale=CS%Ref_val_u%scale, spongeOnGrid=CS%SpongeDataOngrid, m_to_Z=US%m_to_Z, & + answer_date=CS%hor_regrid_answer_date) ! Initialize mask_z halos to zero before pass_var, in case of no update mask_z(G%isc-1, G%jsc:G%jec, :) = 0. @@ -1118,10 +1105,10 @@ subroutine apply_ALE_sponge(h, dt, G, GV, US, CS, Time) deallocate(sp_val, mask_u, mask_z, hsrc) nz_data = CS%Ref_val_v%nz_data ! Interpolate from the external horizontal grid and in time - call horiz_interp_and_extrap_tracer(CS%Ref_val_v%id, Time, CS%Ref_val_v%scale, G, sp_val, & - mask_z, z_in, z_edges_in, missing_value, CS%reentrant_x, CS%tripolar_N, .false., & - spongeOnGrid=CS%SpongeDataOngrid, m_to_Z=US%m_to_Z,& - answer_date=CS%hor_regrid_answer_date) + call horiz_interp_and_extrap_tracer(CS%Ref_val_v%id, Time, G, sp_val, & + mask_z, z_in, z_edges_in, missing_value, & + scale=CS%Ref_val_v%scale, spongeOnGrid=CS%SpongeDataOngrid, m_to_Z=US%m_to_Z,& + answer_date=CS%hor_regrid_answer_date) ! Initialize mask_z halos to zero before pass_var, in case of no update mask_z(G%isc:G%iec, G%jsc-1, :) = 0. mask_z(G%isc:G%iec, G%jec+1, :) = 0.