From c3c56c95f88e43237cb9cfe3b381d9b896fd250c Mon Sep 17 00:00:00 2001 From: laurenchilutti <60401591+laurenchilutti@users.noreply.github.com> Date: Fri, 23 Dec 2022 09:13:36 -0500 Subject: [PATCH] fix: Bug Flag cleanup (#1089) BREAKING CHANGE: The following namelist flags are deprecated and will now result in a fatal error if used: grid_center_bug, read_data_bug, retain_cm3_bug, and reproduce_siena --- data_override/data_override.F90 | 15 +++++--- data_override/get_grid_version.F90 | 37 +++++-------------- fms/fms_io.F90 | 6 ++++ fms/read_data_2d.inc | 2 +- fms/read_data_3d.inc | 2 +- fms/read_data_4d.inc | 2 +- horiz_interp/horiz_interp.F90 | 9 ++--- interpolator/interpolator.F90 | 57 +++++++++++------------------- mosaic/mosaic_util.c | 34 ++++-------------- mosaic/mosaic_util.h | 5 --- 10 files changed, 59 insertions(+), 110 deletions(-) diff --git a/data_override/data_override.F90 b/data_override/data_override.F90 index 0dcee39706..e95d5943d9 100644 --- a/data_override/data_override.F90 +++ b/data_override/data_override.F90 @@ -201,6 +201,13 @@ subroutine data_override_init(Atm_domain_in, Ocean_domain_in, Ice_domain_in, Lan unit = stdlog() write(unit, data_override_nml) +! grid_center_bug is no longer supported. +if (grid_center_bug) then + call mpp_error(FATAL, "data_override_init: You have overridden the default value of " // & + "grid_center_bug and set it to .true. in data_override_nml. This was a temporary workaround " // & + "that is no longer supported. Please remove this namelist variable.") +endif + ! if(module_is_initialized) return atm_on = PRESENT(Atm_domain_in) @@ -290,27 +297,27 @@ subroutine data_override_init(Atm_domain_in, Ocean_domain_in, Ice_domain_in, Lan call mpp_get_compute_domain( atm_domain,is,ie,js,je) allocate(lon_local_atm(is:ie,js:je), lat_local_atm(is:ie,js:je)) call get_grid_version_1(grid_file, 'atm', atm_domain, is, ie, js, je, lon_local_atm, lat_local_atm, & - min_glo_lon_atm, max_glo_lon_atm, grid_center_bug ) + min_glo_lon_atm, max_glo_lon_atm ) endif if (ocn_on .and. .not. allocated(lon_local_ocn) ) then call mpp_get_compute_domain( ocn_domain,is,ie,js,je) allocate(lon_local_ocn(is:ie,js:je), lat_local_ocn(is:ie,js:je)) call get_grid_version_1(grid_file, 'ocn', ocn_domain, is, ie, js, je, lon_local_ocn, lat_local_ocn, & - min_glo_lon_ocn, max_glo_lon_ocn, grid_center_bug ) + min_glo_lon_ocn, max_glo_lon_ocn ) endif if (lnd_on .and. .not. allocated(lon_local_lnd) ) then call mpp_get_compute_domain( lnd_domain,is,ie,js,je) allocate(lon_local_lnd(is:ie,js:je), lat_local_lnd(is:ie,js:je)) call get_grid_version_1(grid_file, 'lnd', lnd_domain, is, ie, js, je, lon_local_lnd, lat_local_lnd, & - min_glo_lon_lnd, max_glo_lon_lnd, grid_center_bug ) + min_glo_lon_lnd, max_glo_lon_lnd ) endif if (ice_on .and. .not. allocated(lon_local_ice) ) then call mpp_get_compute_domain( ice_domain,is,ie,js,je) allocate(lon_local_ice(is:ie,js:je), lat_local_ice(is:ie,js:je)) call get_grid_version_1(grid_file, 'ice', ice_domain, is, ie, js, je, lon_local_ice, lat_local_ice, & - min_glo_lon_ice, max_glo_lon_ice, grid_center_bug ) + min_glo_lon_ice, max_glo_lon_ice ) endif else if (atm_on .and. .not. allocated(lon_local_atm) ) then diff --git a/data_override/get_grid_version.F90 b/data_override/get_grid_version.F90 index a5b5f2370a..451c570d43 100644 --- a/data_override/get_grid_version.F90 +++ b/data_override/get_grid_version.F90 @@ -61,15 +61,13 @@ subroutine check_grid_sizes(domain_name, Domain, nlon, nlat) end subroutine check_grid_sizes !> Get global lon and lat of three model (target) grids, with a given file name -subroutine get_grid_version_1(grid_file, mod_name, domain, isc, iec, jsc, jec, lon, lat, min_lon, max_lon, & - & grid_center_bug) +subroutine get_grid_version_1(grid_file, mod_name, domain, isc, iec, jsc, jec, lon, lat, min_lon, max_lon) character(len=*), intent(in) :: grid_file !< name of grid file character(len=*), intent(in) :: mod_name !< module name type(domain2d), intent(in) :: domain !< 2D domain integer, intent(in) :: isc, iec, jsc, jec real, dimension(isc:,jsc:), intent(out) :: lon, lat real, intent(out) :: min_lon, max_lon - logical, intent(in), optional :: grid_center_bug !< Enables legacy behaviour integer :: i, j, siz(4) integer :: nlon, nlat !< size of global lon and lat @@ -83,7 +81,6 @@ subroutine get_grid_version_1(grid_file, mod_name, domain, isc, iec, jsc, jec, l integer :: start(2), nread(2) type(FmsNetcdfDomainFile_t) :: fileobj integer :: ndims !< Number of dimensions - logical :: gc_bug !< local grid_center_bug variable, default is .false. if(.not. open_file(fileobj, grid_file, 'read', domain )) then call mpp_error(FATAL, 'data_override_mod(get_grid_version_1): Error in opening file '//trim(grid_file)) @@ -118,15 +115,6 @@ subroutine get_grid_version_1(grid_file, mod_name, domain, isc, iec, jsc, jec, l lat(:,:) = (lat_vert(:,:,1) + lat_vert(:,:,2) + lat_vert(:,:,3) + lat_vert(:,:,4))*0.25 else - if (present(grid_center_bug)) then - gc_bug = grid_center_bug - else - gc_bug = .false. - endif - - if(gc_bug) call mpp_error(NOTE, & - 'data_override: grid_center_bug is set to true, the grid center location may be incorrect') - ndims = get_variable_num_dimensions(fileobj, 'geolon_vert_t') call get_variable_size(fileobj, 'geolon_vert_t', siz(1:ndims)) nlon = siz(1) - 1; nlat = siz(2) - 1; @@ -141,23 +129,14 @@ subroutine get_grid_version_1(grid_file, mod_name, domain, isc, iec, jsc, jec, l call read_data(fileobj, 'geolon_vert_t', lon_vert(:,:,1), corner=start, edge_lengths=nread) call read_data(fileobj, 'geolat_vert_t', lat_vert(:,:,1), corner=start, edge_lengths=nread) - if(gc_bug) then - do j = jsc, jec - do i = isc, iec - lon(i,j) = (lon_vert(i,j,1) + lon_vert(i+1,j,1))/2. - lat(i,j) = (lat_vert(i,j,1) + lat_vert(i,j+1,1))/2. - enddo - enddo - else - do j = jsc, jec - do i = isc, iec - lon(i,j) = (lon_vert(i,j,1) + lon_vert(i+1,j,1) + & - lon_vert(i+1,j+1,1) + lon_vert(i,j+1,1))*0.25 - lat(i,j) = (lat_vert(i,j,1) + lat_vert(i+1,j,1) + & - lat_vert(i+1,j+1,1) + lat_vert(i,j+1,1))*0.25 - enddo + do j = jsc, jec + do i = isc, iec + lon(i,j) = (lon_vert(i,j,1) + lon_vert(i+1,j,1) + & + lon_vert(i+1,j+1,1) + lon_vert(i,j+1,1))*0.25 + lat(i,j) = (lat_vert(i,j,1) + lat_vert(i+1,j,1) + & + lat_vert(i+1,j+1,1) + lat_vert(i,j+1,1))*0.25 enddo - end if + enddo endif deallocate(lon_vert) deallocate(lat_vert) diff --git a/fms/fms_io.F90 b/fms/fms_io.F90 index 0c900a1859..ce1069948e 100644 --- a/fms/fms_io.F90 +++ b/fms/fms_io.F90 @@ -693,6 +693,12 @@ subroutine fms_io_init() endif ! take namelist options if present +! read_data_bug is no longer supported. + if (read_data_bug) then + call mpp_error(FATAL, "fms_io_init: You have overridden the default value of " // & + "read_data_bug and set it to .true. in fms_io_nml. This was a temporary workaround " // & + "that is no longer supported. Please remove this namelist variable.") + endif ! determine packsize pack_size = size(transfer(doubledata, realarray)) diff --git a/fms/read_data_2d.inc b/fms/read_data_2d.inc index 42d931ddec..c49b76e85c 100644 --- a/fms/read_data_2d.inc +++ b/fms/read_data_2d.inc @@ -43,7 +43,7 @@ ! gdata, len, 0 ) call mpp_broadcast ( gdata, len, mpp_root_pe() ) endif - if(.not. no_halo .or. read_data_bug) then + if(.not. no_halo) then ! data dimensioned on data domain data(is:ie,js:je) = gdata(is:ie,js:je) else diff --git a/fms/read_data_3d.inc b/fms/read_data_3d.inc index f89d584fc2..6b0ebd3b73 100644 --- a/fms/read_data_3d.inc +++ b/fms/read_data_3d.inc @@ -42,7 +42,7 @@ call mpp_broadcast ( gdata, len, mpp_root_pe() ) endif ! return data for compute domain - if(.not. no_halo .or. read_data_bug) then + if(.not. no_halo) then ! data defined on data domain data(is:ie,js:je,:) = gdata(is:ie,js:je,:) diff --git a/fms/read_data_4d.inc b/fms/read_data_4d.inc index 0ed8e79034..d2abe79132 100644 --- a/fms/read_data_4d.inc +++ b/fms/read_data_4d.inc @@ -42,7 +42,7 @@ call mpp_broadcast ( gdata, len, mpp_root_pe() ) endif ! return data for compute domain - if(.not. no_halo .or. read_data_bug) then + if(.not. no_halo) then ! data defined on data domain data(is:ie,js:je,:,:) = gdata(is:ie,js:je,:,:) else diff --git a/horiz_interp/horiz_interp.F90 b/horiz_interp/horiz_interp.F90 index 07e118f6af..9d694f4d21 100644 --- a/horiz_interp/horiz_interp.F90 +++ b/horiz_interp/horiz_interp.F90 @@ -231,10 +231,11 @@ subroutine horiz_interp_init write (unit, nml=horiz_interp_nml) endif - if( reproduce_siena ) then - call mpp_error(FATAL, "horiz_interp_mod: You have overridden the default value of reproduce_siena " // & - "and set it to .true. in horiz_interp_nml. This is a temporary workaround to " // & - "allow for consistency in continuing experiments. Please remove this namelist " ) + if (reproduce_siena) then + call mpp_error(FATAL, "horiz_interp_mod: You have overridden the default value of " // & + "reproduce_siena and set it to .true. in horiz_interp_nml. This was a temporary workaround to " // & + "allow for consistency in continuing experiments and is no longer supported. " // & + "Please remove this namelist.") endif call horiz_interp_conserve_init diff --git a/interpolator/interpolator.F90 b/interpolator/interpolator.F90 index 785ee68835..4a5355040a 100644 --- a/interpolator/interpolator.F90 +++ b/interpolator/interpolator.F90 @@ -321,7 +321,7 @@ module interpolator_mod logical :: read_all_on_init = .false. !< No description integer :: verbose = 0 !< No description logical :: conservative_interp = .true. !< No description -logical :: retain_cm3_bug = .true. !< No description +logical :: retain_cm3_bug = .false. !< No description logical :: use_mpp_io = .false. !< Set to true to use mpp_io, otherwise fms2io is used namelist /interpolator_nml/ & @@ -445,6 +445,13 @@ subroutine interpolator_init( clim_type, file_name, lonb_mod, latb_mod, & read (input_nml_file, nml=interpolator_nml, iostat=io) ierr = check_nml_error(io,'interpolator_nml') + ! retain_cm3_bug is no longer supported. + if (retain_cm3_bug) then + call mpp_error(FATAL, "interpolator_init: You have overridden the default value of " // & + "retain_cm3_bug and set it to .true. in interpolator_nml. This was a temporary workaround " // & + "that is no longer supported. Please remove this namelist variable.") + endif + !--------------------------------------------------------------------- ! write version number and namelist to logfile. !--------------------------------------------------------------------- @@ -1581,9 +1588,7 @@ subroutine obtain_interpolator_time_slices (clim_type, Time) call time_interp(Time, month, clim_type%tweight3, taum, taup, err_msg=err_msg ) ! tweight3 is !! the time weight between the months. - if ( .not. retain_cm3_bug ) then - if (taum==2 .and. taup==2) clim_type%tweight3 = 1. ! protect against post-perth time_interp behavior - end if + if (taum==2 .and. taup==2) clim_type%tweight3 = 1. ! protect against post-perth time_interp behavior if(trim(err_msg) /= '') then call mpp_error(FATAL,'interpolator_timeslice 3: '//trim(err_msg)//' file='//trim(clim_type%file_name), FATAL) endif @@ -1594,9 +1599,7 @@ subroutine obtain_interpolator_time_slices (clim_type, Time) t_prev = set_date(yearm, climmonth, climday, climhour, climminute, climsecond) call time_interp(t_prev, month, clim_type%tweight1, taum, taup, err_msg=err_msg ) !tweight1 is !! the time weight between the climatology years. - if ( .not. retain_cm3_bug ) then - if (taum==2 .and. taup==2) clim_type%tweight1 = 1. ! protect against post-perth time_interp behavior - end if + if (taum==2 .and. taup==2) clim_type%tweight1 = 1. ! protect against post-perth time_interp behavior if(trim(err_msg) /= '') then call mpp_error(FATAL,'interpolator_timeslice 4: '// & & trim(err_msg)//' file='//trim(clim_type%file_name), FATAL) @@ -1608,9 +1611,7 @@ subroutine obtain_interpolator_time_slices (clim_type, Time) t_next = set_date(yearp, climmonth, climday, climhour, climminute, climsecond) call time_interp(t_next, month, clim_type%tweight2, taum, taup, err_msg=err_msg ) !tweight1 is !! the time weight between the climatology years. - if ( .not. retain_cm3_bug ) then - if (taum==2 .and. taup==2) clim_type%tweight2 = 1. ! protect against post-perth time_interp behavior - end if + if (taum==2 .and. taup==2) clim_type%tweight2 = 1. ! protect against post-perth time_interp behavior if(trim(err_msg) /= '') then call mpp_error(FATAL,'interpolator_timeslice 5: '// & & trim(err_msg)//' file='//trim(clim_type%file_name), FATAL) @@ -1949,9 +1950,7 @@ subroutine interpolator_4D(clim_type, Time, phalf, interp_data, & call time_interp(Time, month, clim_type%tweight3, taum, taup, err_msg=err_msg ) ! tweight3 is !! the time weight between the months. - if ( .not. retain_cm3_bug ) then - if (taum==2 .and. taup==2) clim_type%tweight3 = 1. ! protect against post-perth time_interp behavior - end if + if (taum==2 .and. taup==2) clim_type%tweight3 = 1. ! protect against post-perth time_interp behavior if(trim(err_msg) /= '') then call mpp_error(FATAL,'interpolator_4D 3: '//trim(err_msg)//' file='//trim(clim_type%file_name), FATAL) endif @@ -1962,9 +1961,7 @@ subroutine interpolator_4D(clim_type, Time, phalf, interp_data, & t_prev = set_date(yearm, climmonth, climday, climhour, climminute, climsecond) call time_interp(t_prev, month, clim_type%tweight1, taum, taup, err_msg=err_msg ) !tweight1 is !! the time weight between the climatology years. - if ( .not. retain_cm3_bug ) then - if (taum==2 .and. taup==2) clim_type%tweight1 = 1. ! protect against post-perth time_interp behavior - end if + if (taum==2 .and. taup==2) clim_type%tweight1 = 1. ! protect against post-perth time_interp behavior if(trim(err_msg) /= '') then call mpp_error(FATAL,'interpolator_4D 4: '//trim(err_msg)//' file='//trim(clim_type%file_name), FATAL) endif @@ -1974,9 +1971,7 @@ subroutine interpolator_4D(clim_type, Time, phalf, interp_data, & t_next = set_date(yearp, climmonth, climday, climhour, climminute, climsecond) call time_interp(t_next, month, clim_type%tweight2, taum, taup, err_msg=err_msg ) !tweight1 is !! the time weight between the climatology years. - if ( .not. retain_cm3_bug ) then - if (taum==2 .and. taup==2) clim_type%tweight2 = 1. ! protect against post-perth time_interp behavior - end if + if (taum==2 .and. taup==2) clim_type%tweight2 = 1. ! protect against post-perth time_interp behavior if(trim(err_msg) /= '') then call mpp_error(FATAL,'interpolator_4D 5: '//trim(err_msg)//' file='//trim(clim_type%file_name), FATAL) endif @@ -2405,9 +2400,7 @@ subroutine interpolator_3D(clim_type, Time, phalf, interp_data,field_name, is,js call time_interp(Time, month, clim_type%tweight3, taum, taup, err_msg=err_msg ) ! tweight3 is !! the time weight between the months. - if ( .not. retain_cm3_bug ) then - if (taum==2 .and. taup==2) clim_type%tweight3 = 1. ! protect against post-perth time_interp behavior - end if + if (taum==2 .and. taup==2) clim_type%tweight3 = 1. ! protect against post-perth time_interp behavior if(trim(err_msg) /= '') then call mpp_error(FATAL,'interpolator_3D 3: '//trim(err_msg)//' file='//trim(clim_type%file_name), FATAL) endif @@ -2418,9 +2411,7 @@ subroutine interpolator_3D(clim_type, Time, phalf, interp_data,field_name, is,js t_prev = set_date(yearm, climmonth, climday, climhour, climminute, climsecond) call time_interp(t_prev, month, clim_type%tweight1, taum, taup, err_msg=err_msg ) !tweight1 is !! the time weight between the climatology years. - if ( .not. retain_cm3_bug ) then - if (taum==2 .and. taup==2) clim_type%tweight1 = 1. ! protect against post-perth time_interp behavior - end if + if (taum==2 .and. taup==2) clim_type%tweight1 = 1. ! protect against post-perth time_interp behavior if(trim(err_msg) /= '') then call mpp_error(FATAL,'interpolator_3D 4: '//trim(err_msg)//' file='//trim(clim_type%file_name), FATAL) endif @@ -2431,9 +2422,7 @@ subroutine interpolator_3D(clim_type, Time, phalf, interp_data,field_name, is,js t_next = set_date(yearp, climmonth, climday, climhour, climminute, climsecond) call time_interp(t_next, month, clim_type%tweight2, taum, taup, err_msg=err_msg ) !tweight1 is !! the time weight between the climatology years. - if ( .not. retain_cm3_bug ) then - if (taum==2 .and. taup==2) clim_type%tweight2 = 1. ! protect against post-perth time_interp behavior - end if + if (taum==2 .and. taup==2) clim_type%tweight2 = 1. ! protect against post-perth time_interp behavior if(trim(err_msg) /= '') then call mpp_error(FATAL,'interpolator_3D 5: '//trim(err_msg)//' file='//trim(clim_type%file_name), FATAL) endif @@ -2853,9 +2842,7 @@ subroutine interpolator_2D(clim_type, Time, interp_data, field_name, is, js, cli call time_interp(Time, month, clim_type%tweight3, taum, taup, err_msg=err_msg ) ! tweight3 is !! the time weight between the months. - if ( .not. retain_cm3_bug ) then - if (taum==2 .and. taup==2) clim_type%tweight3 = 1. ! protect against post-perth time_interp behavior - end if + if (taum==2 .and. taup==2) clim_type%tweight3 = 1. ! protect against post-perth time_interp behavior if(trim(err_msg) /= '') then call mpp_error(FATAL,'interpolator_2D 3: '//trim(err_msg)//' file='//trim(clim_type%file_name), FATAL) endif @@ -2866,9 +2853,7 @@ subroutine interpolator_2D(clim_type, Time, interp_data, field_name, is, js, cli t_prev = set_date(yearm, climmonth, climday, climhour, climminute, climsecond) call time_interp(t_prev, month, clim_type%tweight1, taum, taup, err_msg=err_msg ) !tweight1 is !! the time weight between the climatology years. - if ( .not. retain_cm3_bug ) then - if (taum==2 .and. taup==2) clim_type%tweight1 = 1. ! protect against post-perth time_interp behavior - end if + if (taum==2 .and. taup==2) clim_type%tweight1 = 1. ! protect against post-perth time_interp behavior if(trim(err_msg) /= '') then call mpp_error(FATAL,'interpolator_2D 4: '//trim(err_msg)//' file='//trim(clim_type%file_name), FATAL) endif @@ -2879,9 +2864,7 @@ subroutine interpolator_2D(clim_type, Time, interp_data, field_name, is, js, cli t_next = set_date(yearp, climmonth, climday, climhour, climminute, climsecond) call time_interp(t_next, month, clim_type%tweight2, taum, taup, err_msg=err_msg ) !tweight1 is !! the time weight between the climatology years. - if ( .not. retain_cm3_bug ) then - if (taum==2 .and. taup==2) clim_type%tweight2 = 1. ! protect against post-perth time_interp behavior - end if + if (taum==2 .and. taup==2) clim_type%tweight2 = 1. ! protect against post-perth time_interp behavior if(trim(err_msg) /= '') then call mpp_error(FATAL,'interpolator_2D 5: '//trim(err_msg)//' file='//trim(clim_type%file_name), FATAL) endif diff --git a/mosaic/mosaic_util.c b/mosaic/mosaic_util.c index 103e882dd2..c37f799f18 100644 --- a/mosaic/mosaic_util.c +++ b/mosaic/mosaic_util.c @@ -43,18 +43,6 @@ void error_handler(char *str) error handler: will print out error message and then abort ***********************************************************/ -int reproduce_siena = 0; - -void set_reproduce_siena_true(void) -{ - reproduce_siena = 1; -} - -void set_reproduce_siena_true_(void) -{ - reproduce_siena = 1; -} - void error_handler(const char *msg) { @@ -298,14 +286,9 @@ double poly_area_dimensionless(const double x[], const double y[], int n) if ( fabs(lat1-lat2) < SMALL_VALUE) /* cheap area calculation along latitude */ area -= dx*sin(0.5*(lat1+lat2)); else { - if(reproduce_siena) { - area += dx*(cos(lat1)-cos(lat2))/(lat1-lat2); - } - else { - dy = 0.5*(lat1-lat2); - dat = sin(dy)/dy; - area -= dx*sin(0.5*(lat1+lat2))*dat; - } + dy = 0.5*(lat1-lat2); + dat = sin(dy)/dy; + area -= dx*sin(0.5*(lat1+lat2))*dat; } } if(area < 0) @@ -335,14 +318,9 @@ double poly_area(const double x[], const double y[], int n) if ( fabs(lat1-lat2) < SMALL_VALUE) /* cheap area calculation along latitude */ area -= dx*sin(0.5*(lat1+lat2)); else { - if(reproduce_siena) { - area += dx*(cos(lat1)-cos(lat2))/(lat1-lat2); - } - else { - dy = 0.5*(lat1-lat2); - dat = sin(dy)/dy; - area -= dx*sin(0.5*(lat1+lat2))*dat; - } + dy = 0.5*(lat1-lat2); + dat = sin(dy)/dy; + area -= dx*sin(0.5*(lat1+lat2))*dat; } } if(area < 0) diff --git a/mosaic/mosaic_util.h b/mosaic/mosaic_util.h index 13700bbb69..c12eb08d03 100644 --- a/mosaic/mosaic_util.h +++ b/mosaic/mosaic_util.h @@ -165,11 +165,6 @@ void setInbound(struct Node *interList, struct Node *list); int isInside(struct Node *node); -void set_reproduce_siena_true(void); - - -void set_reproduce_siena_true_(void); - int inside_a_polygon_(double *lon1, double *lat1, int *npts, double *lon2, double *lat2); #endif