From 70f86872e2c37d6fbbe83d401e9be5726e3d2ebb Mon Sep 17 00:00:00 2001 From: William Xu Date: Tue, 10 Dec 2024 09:03:13 -0400 Subject: [PATCH] Streaming Filter The filters and their target frequencies are no longer hard-coded. Instead, multiple filters with tidal or arbitrary frequencies as their target frequencies can be turned on. The filter names are specified in MOM_input and must consist of two letters/numbers. If the name of a filter is the same as the name of a tidal constituent, then the corresponding tidal frequency will be used as its target frequency. Otherwise, the user must specify the target frequency by "FILTER_${FILTER_NAME}_OMEGA" in MOM_input. The restarting capability has also been implemented. Because the filtering is a point-wise operation, all variables are considered as fields, even if they are velocity components. --- src/core/MOM_barotropic.F90 | 94 +++----- .../lateral/MOM_streaming_filter.F90 | 205 ++++++++++++------ 2 files changed, 171 insertions(+), 128 deletions(-) diff --git a/src/core/MOM_barotropic.F90 b/src/core/MOM_barotropic.F90 index 21d484b9b6..120279b99c 100644 --- a/src/core/MOM_barotropic.F90 +++ b/src/core/MOM_barotropic.F90 @@ -26,8 +26,7 @@ module MOM_barotropic use MOM_restart, only : query_initialized, MOM_restart_CS use MOM_self_attr_load, only : scalar_SAL_sensitivity use MOM_self_attr_load, only : SAL_CS -use MOM_streaming_filter, only : Filt_register, Filt_accum, Filter_CS -use MOM_tidal_forcing, only : tidal_frequency +use MOM_streaming_filter, only : Filt_register, Filt_init, Filt_accum, Filter_CS use MOM_time_manager, only : time_type, real_to_time, operator(+), operator(-) use MOM_unit_scaling, only : unit_scale_type use MOM_variables, only : BT_cont_type, alloc_bt_cont_type @@ -250,10 +249,8 @@ module MOM_barotropic logical :: linearized_BT_PV !< If true, the PV and interface thicknesses used !! in the barotropic Coriolis calculation is time !! invariant and linearized. - logical :: use_filter_m2 !< If true, apply streaming band-pass filter for detecting - !! instantaneous tidal signals. - logical :: use_filter_k1 !< If true, apply streaming band-pass filter for detecting - !! instantaneous tidal signals. + logical :: use_filter !< If true, use streaming band-pass filter to detect the + !! instantaneous tidal signals in the simulation. logical :: use_wide_halos !< If true, use wide halos and march in during the !! barotropic time stepping for efficiency. logical :: clip_velocity !< If true, limit any velocity components that are @@ -297,10 +294,8 @@ module MOM_barotropic type(hor_index_type), pointer :: debug_BT_HI => NULL() !< debugging copy of horizontal index_type type(SAL_CS), pointer :: SAL_CSp => NULL() !< Control structure for SAL type(harmonic_analysis_CS), pointer :: HA_CSp => NULL() !< Control structure for harmonic analysis - type(Filter_CS) :: Filt_CS_um2, & !< Control structures for the M2 streaming filter - Filt_CS_vm2, & !< Control structures for the M2 streaming filter - Filt_CS_uk1, & !< Control structures for the K1 streaming filter - Filt_CS_vk1 !< Control structures for the K1 streaming filter + type(Filter_CS) :: Filt_CS_u, & !< Control structures for the streaming band-pass filter of ubt + Filt_CS_v !< Control structures for the streaming band-pass filter of vbt logical :: module_is_initialized = .false. !< If true, module has been initialized integer :: isdw !< The lower i-memory limit for the wide halo arrays. @@ -608,8 +603,8 @@ subroutine btstep(U_in, V_in, eta_in, dt, bc_accel_u, bc_accel_v, forces, pbce, DCor_v, & ! An averaged total thickness at v points [H ~> m or kg m-2]. Datv ! Basin depth at v-velocity grid points times the x-grid ! spacing [H L ~> m2 or kg m-1]. - real, dimension(:,:), pointer :: um2, uk1, vm2, vk1 - ! M2 and K1 velocities from the output of streaming filters [m s-1] + real, dimension(:,:,:), pointer :: ufilt, vfilt + ! Filtered velocities from the output of streaming filters [m s-1] real, target, dimension(SZIW_(CS),SZJW_(CS)) :: & eta, & ! The barotropic free surface height anomaly or column mass ! anomaly [H ~> m or kg m-2] @@ -1598,15 +1593,9 @@ subroutine btstep(U_in, V_in, eta_in, dt, bc_accel_u, bc_accel_v, forces, pbce, endif ; enddo ; enddo endif - ! Here is an example of how the filter equations are time stepped to determine the M2 and K1 velocities. - ! The filters are initialized and registered in subroutine barotropic_init. - if (CS%use_filter_m2) then - call Filt_accum(ubt, um2, CS%Time, US, CS%Filt_CS_um2) - call Filt_accum(vbt, vm2, CS%Time, US, CS%Filt_CS_vm2) - endif - if (CS%use_filter_k1) then - call Filt_accum(ubt, uk1, CS%Time, US, CS%Filt_CS_uk1) - call Filt_accum(vbt, vk1, CS%Time, US, CS%Filt_CS_vk1) + if (CS%use_filter) then + call Filt_accum(ubt, ufilt, CS%Time, US, CS%Filt_CS_u) + call Filt_accum(vbt, vfilt, CS%Time, US, CS%Filt_CS_v) endif ! Zero out the arrays for various time-averaged quantities. @@ -4984,6 +4973,12 @@ subroutine barotropic_init(u, v, h, eta, Time, G, GV, US, param_file, diag, CS, endif ! len_trim(wave_drag_file) > 0 endif ! CS%linear_wave_drag + ! Initialize streaming band-pass filters + if (CS%use_filter) then + call Filt_init(param_file, US, CS%Filt_CS_u, restart_CS) + call Filt_init(param_file, US, CS%Filt_CS_v, restart_CS) + endif + CS%dtbt_fraction = 0.98 ; if (dtbt_input < 0.0) CS%dtbt_fraction = -dtbt_input dtbt_tmp = -1.0 @@ -5268,9 +5263,8 @@ subroutine register_barotropic_restarts(HI, GV, US, param_file, CS, restart_CS) ! Local variables type(vardesc) :: vd(3) character(len=40) :: mdl = "MOM_barotropic" ! This module's name. + integer :: n_filters !< Number of streaming band-pass filters to be used in the simulation. integer :: isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB - real :: am2, ak1 !< Bandwidth parameters of the M2 and K1 streaming filters [nondim] - real :: om2, ok1 !< Target frequencies of the M2 and K1 streaming filters [rad T-1 ~> rad s-1] isd = HI%isd ; ied = HI%ied ; jsd = HI%jsd ; jed = HI%jed IsdB = HI%IsdB ; IedB = HI%IedB ; JsdB = HI%JsdB ; JedB = HI%JedB @@ -5283,33 +5277,6 @@ subroutine register_barotropic_restarts(HI, GV, US, param_file, CS, restart_CS) "sum(u dh_dt) while also correcting for truncation errors.", & default=.false., do_not_log=.true.) - call get_param(param_file, mdl, "STREAMING_FILTER_M2", CS%use_filter_m2, & - "If true, turn on streaming band-pass filter for detecting "//& - "instantaneous tidal signals.", default=.false.) - call get_param(param_file, mdl, "STREAMING_FILTER_K1", CS%use_filter_k1, & - "If true, turn on streaming band-pass filter for detecting "//& - "instantaneous tidal signals.", default=.false.) - call get_param(param_file, mdl, "FILTER_ALPHA_M2", am2, & - "Bandwidth parameter of the streaming filter targeting the M2 frequency. "//& - "Must be positive. To turn off filtering, set FILTER_ALPHA_M2 <= 0.0.", & - default=0.0, units="nondim", do_not_log=.not.CS%use_filter_m2) - call get_param(param_file, mdl, "FILTER_ALPHA_K1", ak1, & - "Bandwidth parameter of the streaming filter targeting the K1 frequency. "//& - "Must be positive. To turn off filtering, set FILTER_ALPHA_K1 <= 0.0.", & - default=0.0, units="nondim", do_not_log=.not.CS%use_filter_k1) - call get_param(param_file, mdl, "TIDE_M2_FREQ", om2, & - "Frequency of the M2 tidal constituent. "//& - "This is only used if TIDES and TIDE_M2"// & - " are true, or if OBC_TIDE_N_CONSTITUENTS > 0 and M2"// & - " is in OBC_TIDE_CONSTITUENTS.", units="rad s-1", default=tidal_frequency("M2"), & - scale=US%T_to_s, do_not_log=.true.) - call get_param(param_file, mdl, "TIDE_K1_FREQ", ok1, & - "Frequency of the K1 tidal constituent. "//& - "This is only used if TIDES and TIDE_K1"// & - " are true, or if OBC_TIDE_N_CONSTITUENTS > 0 and K1"// & - " is in OBC_TIDE_CONSTITUENTS.", units="rad s-1", default=tidal_frequency("K1"), & - scale=US%T_to_s, do_not_log=.true.) - ALLOC_(CS%ubtav(IsdB:IedB,jsd:jed)) ; CS%ubtav(:,:) = 0.0 ALLOC_(CS%vbtav(isd:ied,JsdB:JedB)) ; CS%vbtav(:,:) = 0.0 if (CS%gradual_BT_ICs) then @@ -5338,22 +5305,17 @@ subroutine register_barotropic_restarts(HI, GV, US, param_file, CS, restart_CS) call register_restart_field(CS%dtbt, "DTBT", .false., restart_CS, & longname="Barotropic timestep", units="seconds", conversion=US%T_to_s) - ! Initialize and register streaming filters - if (CS%use_filter_m2) then - if (am2 > 0.0 .and. om2 > 0.0) then - call Filt_register(am2, om2, 'u', HI, CS%Filt_CS_um2) - call Filt_register(am2, om2, 'v', HI, CS%Filt_CS_vm2) - else - CS%use_filter_m2 = .false. - endif - endif - if (CS%use_filter_k1) then - if (ak1 > 0.0 .and. ok1 > 0.0) then - call Filt_register(ak1, ok1, 'u', HI, CS%Filt_CS_uk1) - call Filt_register(ak1, ok1, 'v', HI, CS%Filt_CS_vk1) - else - CS%use_filter_k1 = .false. - endif + ! Initialize and register streaming band-pass filters + call get_param(param_file, mdl, "USE_FILTER", CS%use_filter, & + "If true, use streaming band-pass filters to detect the "//& + "instantaneous tidal signals in the simulation.", default=.false.) + call get_param(param_file, mdl, "N_FILTERS", n_filters, & + "Number of streaming band-pass filters to be used in the simulation.", & + default=0, do_not_log=.not.CS%use_filter) + if (n_filters<=0) CS%use_filter = .false. + if (CS%use_filter) then + call Filt_register(n_filters, 'ubt', 'u', HI, CS%Filt_CS_u, restart_CS) + call Filt_register(n_filters, 'vbt', 'v', HI, CS%Filt_CS_v, restart_CS) endif end subroutine register_barotropic_restarts diff --git a/src/parameterizations/lateral/MOM_streaming_filter.F90 b/src/parameterizations/lateral/MOM_streaming_filter.F90 index f7a4f736c7..b5e23ed64d 100644 --- a/src/parameterizations/lateral/MOM_streaming_filter.F90 +++ b/src/parameterizations/lateral/MOM_streaming_filter.F90 @@ -1,119 +1,200 @@ !> Streaming band-pass filter for detecting the instantaneous tidal signals in the simulation + module MOM_streaming_filter -use MOM_error_handler, only : MOM_mesg, MOM_error, FATAL +use MOM_error_handler, only : MOM_mesg, MOM_error, NOTE, FATAL +use MOM_file_parser, only : get_param, param_file_type use MOM_hor_index, only : hor_index_type +use MOM_io, only : axis_info, set_axis_info +use MOM_restart, only : register_restart_field, query_initialized, MOM_restart_CS +use MOM_tidal_forcing, only : tidal_frequency use MOM_time_manager, only : time_type, time_type_to_real use MOM_unit_scaling, only : unit_scale_type implicit none ; private -public Filt_register, Filt_accum +public Filt_register, Filt_init, Filt_accum #include -!> The control structure for storing the filter infomation of a particular field +!> Control structure for the MOM_streaming_filter module type, public :: Filter_CS ; private - real :: a, & !< Parameter that determines the bandwidth [nondim] - om, & !< Target frequency of the filter [rad T-1 ~> rad s-1] - old_time = -1.0 !< The time of the previous accumulating step [T ~> s] - real, allocatable, dimension(:,:) :: s1, & !< Dummy variable [A] - u1 !< Filtered data [A] + integer :: nf !< Number of filters to be used in the simulation !>@{ Lower and upper bounds of input data integer :: is, ie, js, je !>@} + character(len=8) :: key !< Identifier of the variable to be filtered + character(len=2), allocatable, dimension(:) :: filter_names !< Names of filters + real, allocatable, dimension(:) :: filter_omega !< Target frequencies of filters [rad T-1 ~> rad s-1] + real, allocatable, dimension(:) :: filter_alpha !< Bandwidth parameters of filters [nondim] + real, allocatable, dimension(:,:,:) :: s1, & !< Dummy variable [A] + u1 !< Filtered data [A] + real :: old_time = -1.0 !< The time of the previous accumulating step [T ~> s] end type Filter_CS contains -!> This subroutine registers each of the fields to be filtered. -subroutine Filt_register(a, om, grid, HI, CS) - real, intent(in) :: a !< Parameter that determines the bandwidth [nondim] - real, intent(in) :: om !< Target frequency of the filter [rad T-1 ~> rad s-1] - character(len=*), intent(in) :: grid !< Horizontal grid location: h, u, or v - type(hor_index_type), intent(in) :: HI !< Horizontal index type structure - type(Filter_CS), intent(out) :: CS !< Control structure for the current field +!> This subroutine registers the filter variables given the number of filters and the grid +subroutine Filt_register(nf, key, grid, HI, CS, restart_CS) + integer, intent(in) :: nf !< Number of filters to be used in the simulation + character(len=*), intent(in) :: key !< Identifier of the variable to be filtered + character(len=*), intent(in) :: grid !< Horizontal grid location: "h", "u", or "v" + type(hor_index_type), intent(in) :: HI !< Horizontal index type structure + type(Filter_CS), intent(out) :: CS !< Control structure of MOM_streaming_filter + type(MOM_restart_CS), intent(inout) :: restart_CS !< MOM restart control structure ! Local variables - integer :: isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB - - if (a <= 0.0) call MOM_error(FATAL, "MOM_streaming_filter: bandwidth <= 0") - if (om <= 0.0) call MOM_error(FATAL, "MOM_streaming_filter: target frequency <= 0") - - CS%a = a - CS%om = om + type(axis_info) :: filter_axis(1) + real, dimension(:), allocatable :: n_filters !< Labels of filters [nondim] + integer :: c - isd = HI%isd ; ied = HI%ied ; jsd = HI%jsd ; jed = HI%jed - IsdB = HI%IsdB ; IedB = HI%IedB ; JsdB = HI%JsdB ; JedB = HI%JedB + CS%nf = nf + CS%key = key select case (trim(grid)) case ('h') - allocate(CS%s1(isd:ied,jsd:jed)) ; CS%s1(:,:) = 0.0 - allocate(CS%u1(isd:ied,jsd:jed)) ; CS%u1(:,:) = 0.0 - CS%is = isd ; CS%ie = ied ; CS%js = jsd ; CS%je = jed + CS%is = HI%isd ; CS%ie = HI%ied ; CS%js = HI%jsd ; CS%je = HI%jed case ('u') - allocate(CS%s1(IsdB:IedB,jsd:jed)) ; CS%s1(:,:) = 0.0 - allocate(CS%u1(IsdB:IedB,jsd:jed)) ; CS%u1(:,:) = 0.0 - CS%is = IsdB ; CS%ie = IedB ; CS%js = jsd ; CS%je = jed + CS%is = HI%IsdB ; CS%ie = HI%IedB ; CS%js = HI%jsd ; CS%je = HI%jed case ('v') - allocate(CS%s1(isd:ied,JsdB:JedB)) ; CS%s1(:,:) = 0.0 - allocate(CS%u1(isd:ied,JsdB:JedB)) ; CS%u1(:,:) = 0.0 - CS%is = isd ; CS%ie = ied ; CS%js = JsdB ; CS%je = JedB + CS%is = HI%isd ; CS%ie = HI%ied ; CS%js = HI%JsdB ; CS%je = HI%JedB case default call MOM_error(FATAL, "MOM_streaming_filter: horizontal grid not supported") end select + allocate(CS%s1(CS%is:CS%ie, CS%js:CS%je, nf)) ; CS%s1(:,:,:) = 0.0 + allocate(CS%u1(CS%is:CS%ie, CS%js:CS%je, nf)) ; CS%u1(:,:,:) = 0.0 + + ! Register restarts for s1 and u1 + allocate(n_filters(nf)) + + do c=1,nf ; n_filters(c) = c ; enddo + + call set_axis_info(filter_axis(1), "n_filters", "", "number of filters", nf, n_filters, "N", 1) + + call register_restart_field(CS%s1(:,:,:), "Filter_"//trim(key)//"_s1", .false., restart_CS, & + longname="Dummy variable for streaming band-pass filter", & + hor_grid=trim(grid), z_grid="1", t_grid="s", extra_axes=filter_axis) + call register_restart_field(CS%u1(:,:,:), "Filter_"//trim(key)//"_u1", .false., restart_CS, & + longname="Output of streaming band-pass filter", & + hor_grid=trim(grid), z_grid="1", t_grid="s", extra_axes=filter_axis) + end subroutine Filt_register -!> This subroutine timesteps the filter equations. It takes model output u at the current time step as the input, -!! and returns tidal signal u1 as the output, which is the solution of a set of two ODEs (the filter equations). +!> This subroutine initializes the filters +subroutine Filt_init(param_file, US, CS, restart_CS) + type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters + type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type + type(Filter_CS), intent(inout) :: CS !< Control structure of MOM_streaming_filter + type(MOM_restart_CS), intent(in) :: restart_CS !< MOM restart control structure + + ! Local variables + character(len=40) :: mdl = "MOM_streaming_filter" !< This module's name + character(len=50) :: filter_name_str !< List of filters to be registered + character(len=200) :: mesg + integer :: c + + call get_param(param_file, mdl, "FILTER_NAMES", filter_name_str, & + "Names of streaming band-pass filters to be used in the simulation.", & + fail_if_missing=.true.) + allocate(CS%filter_names(CS%nf)) + allocate(CS%filter_omega(CS%nf)) + allocate(CS%filter_alpha(CS%nf)) + read(filter_name_str, *) CS%filter_names + + do c=1,CS%nf + ! If filter_name_str consists of tidal constituents, use tidal frequencies. + call get_param(param_file, mdl, "FILTER_"//trim(CS%filter_names(c))//"_OMEGA", & + CS%filter_omega(c), "Target frequency of the "//trim(CS%filter_names(c))//& + " filter. This is used if USE_FILTER is true and "//trim(CS%filter_names(c))//& + " is in FILTER_NAMES.", units="rad s-1", scale=US%T_to_s, default=0.0) + call get_param(param_file, mdl, "FILTER_"//trim(CS%filter_names(c))//"_ALPHA", & + CS%filter_alpha(c), "Bandwidth parameter of the "//trim(CS%filter_names(c))//& + " filter. Must be positive.", units="nondim", fail_if_missing=.true.) + + if (CS%filter_omega(c)<=0.0) CS%filter_omega(c) = tidal_frequency(trim(CS%filter_names(c))) + if (CS%filter_alpha(c)<=0.0) call MOM_error(FATAL, "MOM_streaming_filter: bandwidth <= 0") + + write(mesg,*) "MOM_streaming_filter: ", trim(CS%filter_names(c)), & + " filter registered, target frequency = ", CS%filter_omega(c), & + ", bandwidth = ", CS%filter_alpha(c) + call MOM_error(NOTE, trim(mesg)) + enddo + + if (query_initialized(CS%s1, "Filter_"//trim(CS%key)//"_s1", restart_CS)) then + write(mesg,*) "MOM_streaming_filter: Dummy variable for filter ", trim(CS%key), & + " found in restart files." + else + write(mesg,*) "MOM_streaming_filter: Dummy variable for filter ", trim(CS%key), & + " not found in restart files. The filter will spin up from zeros." + endif + call MOM_error(NOTE, trim(mesg)) + + if (query_initialized(CS%u1, "Filter_"//trim(CS%key)//"_u1", restart_CS)) then + write(mesg,*) "MOM_streaming_filter: Output of filter ", trim(CS%key), & + " found in restart files." + else + write(mesg,*) "MOM_streaming_filter: Output of filter ", trim(CS%key), & + " not found in restart files. The filter will spin up from zeros." + endif + call MOM_error(NOTE, trim(mesg)) + +end subroutine Filt_init + +!> This subroutine timesteps the filter equations. Here, u is the broadband input signal from the model, +!! and u1 is the filtered, narrowband output signal, obtained from the solution of the filter equations. subroutine Filt_accum(u, u1, Time, US, CS) - real, dimension(:,:), pointer, intent(out) :: u1 !< Output of the filter [A] - type(time_type), intent(in) :: Time !< The current model time - type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type - type(Filter_CS), target, intent(inout) :: CS !< Control structure of the MOM_streaming_filter module + real, dimension(:,:,:), pointer, intent(out) :: u1 !< Output of the filter [A] + type(time_type), intent(in) :: Time !< The current model time + type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type + type(Filter_CS), target, intent(inout) :: CS !< Control structure of MOM_streaming_filter real, dimension(CS%is:CS%ie,CS%js:CS%je), intent(in) :: u !< Input into the filter [A] ! Local variables real :: now, & !< The current model time [T ~> s] dt, & !< Time step size for the filter equations [T ~> s] c1, c2 !< Coefficients for the filter equations [nondim] - integer :: i, j, is, ie, js, je + integer :: i, j, k now = US%s_to_T * time_type_to_real(Time) - is = CS%is ; ie = CS%ie ; js = CS%js ; je = CS%je - - ! Initialize u1 - if (CS%old_time < 0.0) then - CS%old_time = now - CS%u1(:,:) = u(:,:) - endif - + if (CS%old_time<0.0) CS%old_time = now dt = now - CS%old_time CS%old_time = now ! Timestepping - c1 = CS%om * dt - c2 = 1.0 - CS%a * c1 - - do j=js,je ; do i=is,ie - CS%s1(i,j) = c1 * CS%u1(i,j) + CS%s1(i,j) - CS%u1(i,j) = -c1 * (CS%s1(i,j) - CS%a * u(i,j)) + c2 * CS%u1(i,j) - enddo; enddo + do k=1,CS%nf + c1 = CS%filter_omega(k) * dt + c2 = 1.0 - CS%filter_alpha(k) * c1 + + do j=CS%js,CS%je ; do i=CS%is,CS%ie + CS%s1(i,j,k) = c1 * CS%u1(i,j,k) + CS%s1(i,j,k) + CS%u1(i,j,k) = -c1 * (CS%s1(i,j,k) - CS%filter_alpha(k) * u(i,j)) + c2 * CS%u1(i,j,k) + enddo; enddo + enddo u1 => CS%u1 end subroutine Filt_accum -!> \namespace streaming_filter +!> \namespace mom_streaming_filter +!! +!! By William Xu (chengzhu.xu@oregonstate.edu) and Ed Zaron +!! +!! The algorithm detects the instantaneous, narrowband tidal signals (u1) from the broadband model +!! output (u) by solving a set of coupled ODEs (the filter equations) at each time step. +!! +!! Major revision on Dec 9, 2024: The filters are no longer hard-coded. Instead, multiple filters +!! with tidal frequencies or arbitrary frequencies as their target frequencies can be turned on. +!! The filter names are specified in MOM_input and must consist of two letters/numbers. If the +!! name of a filter is the same as the name of a tidal constituent, then the corresponding tidal +!! frequency will be used as its target frequency. Otherwise, the user must specify the target +!! frequency. In either case, the target frequency is specified by "FILTER_${FILTER_NAME}_OMEGA". !! -!! This module detects instantaneous tidal signals in the model output using a set of coupled ODEs (the filter -!! equations), given the target frequency (om) and the bandwidth parameter (a) of the filter. At each timestep, -!! the filter takes model output (u) as the input and returns a time series consisting of sinusoidal motions (u1) -!! near its target frequency. The filtered tidal signals can be used to parameterize frequency-dependent drag, or -!! to detide the model output. See Xu & Zaron (2024) for detail. +!! The restarting capability has also been implemented. Because the filtering is a point-wise +!! operation, all variables are considered as fields, even if they are velocity components. !! !! Reference: Xu, C., & Zaron, E. D. (2024). Detecting instantaneous tidal signals in ocean models utilizing -!! streaming band-pass filters. Journal of Advances in Modeling Earth Systems. Under review. +!! streaming band-pass filters. Journal of Advances in Modeling Earth Systems, 16, e2024MS004319. +!! https://doi.org/10.1029/2024MS004319 end module MOM_streaming_filter