MOM6
mom_oda_driver_mod Module Reference

Detailed Description

Interfaces for MOM6 ensembles and data assimilation.

Data Types

type  oda_cs
 Control structure that contains a transpose of the ocean state across ensemble members. More...
 
type  ptr_mpp_domain
 A structure with a pointer to a domain2d, to allow for the creation of arrays of pointers. More...
 
integer, parameter no_assim = 0
 DA parameters. More...
 
integer, parameter oi_assim =1
 DA parameters. More...
 
integer, parameter eakf_assim =2
 DA parameters. More...
 
subroutine, public init_oda (Time, G, GV, CS)
 initialize First_guess (prior) and Analysis grid information for all ensemble members More...
 
subroutine, public set_prior_tracer (Time, G, GV, h, tv, CS)
 Copy ensemble member tracers to ensemble vector. More...
 
subroutine, public get_posterior_tracer (Time, CS, h, tv, increment)
 Returns posterior adjustments or full state Note that only those PEs associated with an ensemble member receive data. More...
 
subroutine, public oda (Time, CS)
 Gather observations and sall ODA routines. More...
 
subroutine, public oda_end (CS)
 Finalize DA module. More...
 
subroutine init_ocean_ensemble (CS, Grid, GV, ens_size)
 Initialize DA module. More...
 
subroutine, public set_analysis_time (Time, CS)
 Set the next analysis time. More...
 
subroutine, public save_obs_diff (filename, CS)
 Write observation differences to a file. More...
 
subroutine, public apply_oda_tracer_increments (dt, G, tv, h, CS)
 Apply increments to tracers. More...
 

Function/Subroutine Documentation

◆ apply_oda_tracer_increments()

subroutine, public mom_oda_driver_mod::apply_oda_tracer_increments ( real, intent(in)  dt,
type(ocean_grid_type), intent(in)  G,
type(thermo_var_ptrs), intent(inout)  tv,
real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in)  h,
type(oda_cs), intent(inout)  CS 
)

Apply increments to tracers.

Parameters
[in]dtThe tracer timestep [s]
[in]gocean grid structure
[in,out]tvA structure pointing to various thermodynamic variables
[in]hlayer thickness [H ~> m or kg m-2]
[in,out]csthe data assimilation structure

Definition at line 548 of file MOM_oda_driver.F90.

548  real, intent(in) :: dt !< The tracer timestep [s]
549  type(ocean_grid_type), intent(in) :: G !< ocean grid structure
550  type(thermo_var_ptrs), intent(inout) :: tv !< A structure pointing to various thermodynamic variables
551  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
552  intent(in) :: h !< layer thickness [H ~> m or kg m-2]
553  type(ODA_CS), intent(inout) :: CS !< the data assimilation structure
554 

Referenced by mom::step_mom_thermo().

Here is the caller graph for this function:

◆ get_posterior_tracer()

subroutine, public mom_oda_driver_mod::get_posterior_tracer ( type(time_type), intent(in)  Time,
type(oda_cs), pointer  CS,
real, dimension(:,:,:), pointer  h,
type(thermo_var_ptrs), pointer  tv,
logical, intent(in), optional  increment 
)

Returns posterior adjustments or full state Note that only those PEs associated with an ensemble member receive data.

Parameters
[in]timethe current model time
csocean DA control structure
hLayer thicknesses [H ~> m or kg m-2]
tvA structure pointing to various thermodynamic variables
[in]incrementTrue if returning increment only

Definition at line 384 of file MOM_oda_driver.F90.

384  type(time_type), intent(in) :: Time !< the current model time
385  type(ODA_CS), pointer :: CS !< ocean DA control structure
386  real, dimension(:,:,:), pointer :: h !< Layer thicknesses [H ~> m or kg m-2]
387  type(thermo_var_ptrs), pointer :: tv !< A structure pointing to various thermodynamic variables
388  logical, optional, intent(in) :: increment !< True if returning increment only
389 
390  type(ocean_control_struct), pointer :: Ocean_increment=>null()
391  integer :: i, j, m
392  logical :: used, get_inc
393 
394  ! return if not analysis time (retain pointers for h and tv)
395  if (time < cs%Time) return
396 
397 
398  !! switch to global pelist
399  call set_current_pelist(cs%filter_pelist)
400  call mom_mesg('Getting posterior')
401 
402  get_inc = .true.
403  if (present(increment)) get_inc = increment
404 
405  if (get_inc) then
406  allocate(ocean_increment)
407  call init_ocean_ensemble(ocean_increment,cs%Grid,cs%GV,cs%ensemble_size)
408  ocean_increment%T = cs%Ocean_posterior%T - cs%Ocean_prior%T
409  ocean_increment%S = cs%Ocean_posterior%S - cs%Ocean_prior%S
410  endif
411  do m=1,cs%ensemble_size
412  if (get_inc) then
413  call mpp_redistribute(cs%mpp_domain, ocean_increment%T(:,:,:,m), &
414  cs%domains(m)%mpp_domain, cs%tv%T, complete=.true.)
415  call mpp_redistribute(cs%mpp_domain, ocean_increment%S(:,:,:,m), &
416  cs%domains(m)%mpp_domain, cs%tv%S, complete=.true.)
417  else
418  call mpp_redistribute(cs%mpp_domain, cs%Ocean_posterior%T(:,:,:,m), &
419  cs%domains(m)%mpp_domain, cs%tv%T, complete=.true.)
420  call mpp_redistribute(cs%mpp_domain, cs%Ocean_posterior%S(:,:,:,m), &
421  cs%domains(m)%mpp_domain, cs%tv%S, complete=.true.)
422  endif
423  enddo
424 
425  tv => cs%tv
426  h => cs%h
427  !! switch back to ensemble member pelist
428  call set_current_pelist(cs%ensemble_pelist(cs%ensemble_id,:))
429 

References init_ocean_ensemble(), and mom_error_handler::mom_mesg().

Here is the call graph for this function:

◆ init_ocean_ensemble()

subroutine mom_oda_driver_mod::init_ocean_ensemble ( type(ocean_control_struct), pointer  CS,
type(ocean_grid_type), pointer  Grid,
type(verticalgrid_type), pointer  GV,
integer, intent(in)  ens_size 
)
private

Initialize DA module.

Parameters
csPointer to ODA control structure
gridPointer to ocean analysis grid
gvPointer to DA vertical grid
[in]ens_sizeensemble size

Definition at line 467 of file MOM_oda_driver.F90.

467  type(ocean_control_struct), pointer :: CS !< Pointer to ODA control structure
468  type(ocean_grid_type), pointer :: Grid !< Pointer to ocean analysis grid
469  type(verticalGrid_type), pointer :: GV !< Pointer to DA vertical grid
470  integer, intent(in) :: ens_size !< ensemble size
471 
472  integer :: n,is,ie,js,je,nk
473 
474  nk=gv%ke
475  is=grid%isd;ie=grid%ied
476  js=grid%jsd;je=grid%jed
477  cs%ensemble_size=ens_size
478  allocate(cs%T(is:ie,js:je,nk,ens_size))
479  allocate(cs%S(is:ie,js:je,nk,ens_size))
480  allocate(cs%SSH(is:ie,js:je,ens_size))
481  allocate(cs%id_t(ens_size));cs%id_t(:)=-1
482  allocate(cs%id_s(ens_size));cs%id_s(:)=-1
483 ! allocate(CS%U(is:ie,js:je,nk,ens_size))
484 ! allocate(CS%V(is:ie,js:je,nk,ens_size))
485 ! allocate(CS%id_u(ens_size));CS%id_u(:)=-1
486 ! allocate(CS%id_v(ens_size));CS%id_v(:)=-1
487  allocate(cs%id_ssh(ens_size));cs%id_ssh(:)=-1
488 
489  return

Referenced by get_posterior_tracer(), and init_oda().

Here is the caller graph for this function:

◆ init_oda()

subroutine, public mom_oda_driver_mod::init_oda ( type(time_type), intent(in)  Time,
type(ocean_grid_type), pointer  G,
type(verticalgrid_type), intent(in)  GV,
type(oda_cs), intent(inout), pointer  CS 
)

initialize First_guess (prior) and Analysis grid information for all ensemble members

Parameters
[in]timeThe current model time.
gdomain and grid information for ocean model
[in]gvThe ocean's vertical grid structure
[in,out]csThe DA control structure

Definition at line 116 of file MOM_oda_driver.F90.

116 
117  type(time_type), intent(in) :: Time !< The current model time.
118  type(ocean_grid_type), pointer :: G !< domain and grid information for ocean model
119  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
120  type(ODA_CS), pointer, intent(inout) :: CS !< The DA control structure
121 
122 ! Local variables
123  type(thermo_var_ptrs) :: tv_dummy
124  type(dyn_horgrid_type), pointer :: dG=> null()
125  type(hor_index_type), pointer :: HI=> null()
126  type(directories) :: dirs
127 
128  type(grid_type), pointer :: T_grid !< global tracer grid
129  real, dimension(:,:), allocatable :: global2D, global2D_old
130  real, dimension(:), allocatable :: lon1D, lat1D, glon1D, glat1D
131  type(param_file_type) :: PF
132  integer :: n, m, k, i, j, nk
133  integer :: is,ie,js,je,isd,ied,jsd,jed
134  integer :: stdout_unit
135  character(len=32) :: assim_method
136  integer :: npes_pm, ens_info(6), ni, nj
137  character(len=128) :: mesg
138  character(len=32) :: fldnam
139  character(len=30) :: coord_mode
140  character(len=200) :: inputdir, basin_file
141  logical :: reentrant_x, reentrant_y, tripolar_N, symmetric
142 
143  if (associated(cs)) call mpp_error(fatal,'Calling oda_init with associated control structure')
144  allocate(cs)
145 ! Use ens1 parameters , this could be changed at a later time
146 ! if it were desirable to have alternate parameters, e.g. for the grid
147 ! for the analysis
148  call get_mom_input(pf,dirs,ensemble_num=0)
149 
150  call unit_scaling_init(pf, cs%US)
151 
152  call get_param(pf, "MOM", "ASSIM_METHOD", assim_method, &
153  "String which determines the data assimilation method "//&
154  "Valid methods are: \'EAKF\',\'OI\', and \'NO_ASSIM\'", default='NO_ASSIM')
155  call get_param(pf, "MOM", "ASSIM_FREQUENCY", cs%assim_frequency, &
156  "data assimilation frequency in hours")
157  call get_param(pf, "MOM", "USE_REGRIDDING", cs%use_ALE_algorithm , &
158  "If True, use the ALE algorithm (regridding/remapping).\n"//&
159  "If False, use the layered isopycnal algorithm.", default=.false. )
160  call get_param(pf, "MOM", "REENTRANT_X", cs%reentrant_x, &
161  "If true, the domain is zonally reentrant.", default=.true.)
162  call get_param(pf, "MOM", "REENTRANT_Y", cs%reentrant_y, &
163  "If true, the domain is meridionally reentrant.", &
164  default=.false.)
165  call get_param(pf,"MOM", "TRIPOLAR_N", cs%tripolar_N, &
166  "Use tripolar connectivity at the northern edge of the "//&
167  "domain. With TRIPOLAR_N, NIGLOBAL must be even.", &
168  default=.false.)
169  call get_param(pf,"MOM", "NIGLOBAL", cs%ni, &
170  "The total number of thickness grid points in the "//&
171  "x-direction in the physical domain.")
172  call get_param(pf,"MOM", "NJGLOBAL", cs%nj, &
173  "The total number of thickness grid points in the "//&
174  "y-direction in the physical domain.")
175  call get_param(pf, 'MOM', "INPUTDIR", inputdir)
176  inputdir = slasher(inputdir)
177 
178  select case(lowercase(trim(assim_method)))
179  case('eakf')
180  cs%assim_method = eakf_assim
181  case('oi')
182  cs%assim_method = oi_assim
183  case('no_assim')
184  cs%assim_method = no_assim
185  case default
186  call mpp_error(fatal,'Invalid assimilation method provided')
187  end select
188 
189  ens_info = get_ensemble_size()
190  cs%ensemble_size = ens_info(1)
191  npes_pm=ens_info(3)
192  cs%ensemble_id = get_ensemble_id()
193  !! Switch to global pelist
194  allocate(cs%ensemble_pelist(cs%ensemble_size,npes_pm))
195  allocate(cs%filter_pelist(cs%ensemble_size*npes_pm))
196  call get_ensemble_pelist(cs%ensemble_pelist,'ocean')
197  call get_ensemble_filter_pelist(cs%filter_pelist,'ocean')
198 
199  call set_current_pelist(cs%filter_pelist)
200 
201  allocate(cs%domains(cs%ensemble_size))
202  cs%domains(cs%ensemble_id)%mpp_domain => g%Domain%mpp_domain
203  do n=1,cs%ensemble_size
204  if (.not. associated(cs%domains(n)%mpp_domain)) allocate(cs%domains(n)%mpp_domain)
205  call set_root_pe(cs%ensemble_pelist(n,1))
206  call mpp_broadcast_domain(cs%domains(n)%mpp_domain)
207  enddo
208  call set_root_pe(cs%filter_pelist(1))
209  allocate(cs%Grid)
210  ! params NIHALO_ODA, NJHALO_ODA set the DA halo size
211  call mom_domains_init(cs%Grid%Domain,pf,param_suffix='_ODA')
212  allocate(hi)
213  call hor_index_init(cs%Grid%Domain, hi, pf, &
214  local_indexing=.false.) ! Use global indexing for DA
215  call verticalgridinit( pf, cs%GV, cs%US )
216  allocate(dg)
217  call create_dyn_horgrid(dg, hi)
218  call clone_mom_domain(cs%Grid%Domain, dg%Domain,symmetric=.false.)
219  call set_grid_metrics(dg,pf)
220  call mom_initialize_topography(dg%bathyT,dg%max_depth,dg,pf)
221  call mom_initialize_coord(cs%GV, cs%US, pf, .false., &
222  dirs%output_directory, tv_dummy, dg%max_depth)
223  call ale_init(pf, cs%GV, cs%US, dg%max_depth, cs%ALE_CS)
224  call mom_grid_init(cs%Grid, pf, global_indexing=.true.)
225  call ale_updateverticalgridtype(cs%ALE_CS, cs%GV)
226  call copy_dyngrid_to_mom_grid(dg, cs%Grid, cs%US)
227  cs%mpp_domain => cs%Grid%Domain%mpp_domain
228  cs%Grid%ke = cs%GV%ke
229  cs%nk = cs%GV%ke
230  ! initialize storage for prior and posterior
231  allocate(cs%Ocean_prior)
232  call init_ocean_ensemble(cs%Ocean_prior,cs%Grid,cs%GV,cs%ensemble_size)
233  allocate(cs%Ocean_posterior)
234  call init_ocean_ensemble(cs%Ocean_posterior,cs%Grid,cs%GV,cs%ensemble_size)
235  allocate(cs%tv)
236 
237  call get_param(pf, 'oda_driver', "REGRIDDING_COORDINATE_MODE", coord_mode, &
238  "Coordinate mode for vertical regridding.", &
239  default="ZSTAR", fail_if_missing=.false.)
240  call initialize_regridding(cs%regridCS, cs%GV, cs%US, dg%max_depth,pf,'oda_driver',coord_mode,'','')
241  call initialize_remapping(cs%remapCS,'PLM')
242  call set_regrid_params(cs%regridCS, min_thickness=0.)
243  call mpp_get_data_domain(g%Domain%mpp_domain,isd,ied,jsd,jed)
244  if (.not. associated(cs%h)) then
245  allocate(cs%h(isd:ied,jsd:jed,cs%GV%ke)); cs%h(:,:,:)=0.0
246  ! assign thicknesses
247  call ale_initthicknesstocoord(cs%ALE_CS,g,cs%GV,cs%h)
248  endif
249  allocate(cs%tv%T(isd:ied,jsd:jed,cs%GV%ke)); cs%tv%T(:,:,:)=0.0
250  allocate(cs%tv%S(isd:ied,jsd:jed,cs%GV%ke)); cs%tv%S(:,:,:)=0.0
251 
252  call set_axes_info(cs%Grid, cs%GV, cs%US, pf, cs%diag_cs, set_vertical=.true.)
253  do n=1,cs%ensemble_size
254  write(fldnam,'(a,i2.2)') 'temp_prior_',n
255  cs%Ocean_prior%id_t(n)=register_diag_field('ODA',trim(fldnam),cs%diag_cs%axesTL%handles,time, &
256  'ocean potential temperature','degC')
257  write(fldnam,'(a,i2.2)') 'salt_prior_',n
258  cs%Ocean_prior%id_s(n)=register_diag_field('ODA',trim(fldnam),cs%diag_cs%axesTL%handles,time, &
259  'ocean salinity','psu')
260  write(fldnam,'(a,i2.2)') 'temp_posterior_',n
261  cs%Ocean_posterior%id_t(n)=register_diag_field('ODA',trim(fldnam),cs%diag_cs%axesTL%handles,time, &
262  'ocean potential temperature','degC')
263  write(fldnam,'(a,i2.2)') 'salt_posterior_',n
264  cs%Ocean_posterior%id_s(n)=register_diag_field('ODA',trim(fldnam),cs%diag_cs%axesTL%handles,time, &
265  'ocean salinity','psu')
266  enddo
267 
268  call mpp_get_data_domain(cs%mpp_domain,isd,ied,jsd,jed)
269  allocate(cs%oda_grid)
270  cs%oda_grid%x => cs%Grid%geolonT
271  cs%oda_grid%y => cs%Grid%geolatT
272 
273  call get_param(pf, 'oda_driver', "BASIN_FILE", basin_file, &
274  "A file in which to find the basin masks, in variable 'basin'.", &
275  default="basin.nc")
276  basin_file = trim(inputdir) // trim(basin_file)
277  allocate(cs%oda_grid%basin_mask(isd:ied,jsd:jed))
278  cs%oda_grid%basin_mask(:,:) = 0.0
279  call mom_read_data(basin_file,'basin',cs%oda_grid%basin_mask,cs%Grid%domain, timelevel=1)
280 
281 ! get global grid information from ocean_model
282  allocate(t_grid)
283  allocate(t_grid%x(cs%ni,cs%nj))
284  allocate(t_grid%y(cs%ni,cs%nj))
285  allocate(t_grid%basin_mask(cs%ni,cs%nj))
286  call mpp_global_field(cs%mpp_domain, cs%Grid%geolonT, t_grid%x)
287  call mpp_global_field(cs%mpp_domain, cs%Grid%geolatT, t_grid%y)
288  call mpp_global_field(cs%mpp_domain, cs%oda_grid%basin_mask, t_grid%basin_mask)
289  t_grid%ni = cs%ni
290  t_grid%nj = cs%nj
291  t_grid%nk = cs%nk
292  allocate(t_grid%mask(cs%ni,cs%nj,cs%nk))
293  allocate(t_grid%z(cs%ni,cs%nj,cs%nk))
294  allocate(global2d(cs%ni,cs%nj))
295  allocate(global2d_old(cs%ni,cs%nj))
296  t_grid%mask(:,:,:) = 0.0
297  t_grid%z(:,:,:) = 0.0
298 
299  do k = 1, cs%nk
300  call mpp_global_field(g%Domain%mpp_domain, cs%h(:,:,k), global2d)
301  do i=1, cs%ni; do j=1, cs%nj
302  if ( global2d(i,j) > 1 ) then
303  t_grid%mask(i,j,k) = 1.0
304  endif
305  enddo ; enddo
306  if (k == 1) then
307  t_grid%z(:,:,k) = global2d/2
308  else
309  t_grid%z(:,:,k) = t_grid%z(:,:,k-1) + (global2d + global2d_old)/2
310  endif
311  global2d_old = global2d
312  enddo
313 
314  call ocean_da_core_init(cs%mpp_domain, t_grid, cs%Profiles, time)
315 
316  cs%Time=time
317  !! switch back to ensemble member pelist
318  call set_current_pelist(cs%ensemble_pelist(cs%ensemble_id,:))

References mom_ale::ale_init(), mom_ale::ale_initthicknesstocoord(), mom_ale::ale_updateverticalgridtype(), mom_transcribe_grid::copy_dyngrid_to_mom_grid(), mom_dyn_horgrid::create_dyn_horgrid(), eakf_assim, mom_get_input::get_mom_input(), mom_hor_index::hor_index_init(), init_ocean_ensemble(), mom_remapping::initialize_remapping(), mom_string_functions::lowercase(), mom_domains::mom_domains_init(), mom_grid::mom_grid_init(), mom_coord_initialization::mom_initialize_coord(), mom_fixed_initialization::mom_initialize_topography(), no_assim, oi_assim, mom_diag_mediator::set_axes_info(), mom_grid_initialize::set_grid_metrics(), mom_regridding::set_regrid_params(), mom_unit_scaling::unit_scaling_init(), and mom_verticalgrid::verticalgridinit().

Here is the call graph for this function:

◆ oda()

subroutine, public mom_oda_driver_mod::oda ( type(time_type), intent(in)  Time,
type(oda_cs), intent(inout)  CS 
)

Gather observations and sall ODA routines.

Parameters
[in]timethe current model time
[in,out]csthe ocean DA control structure

Definition at line 434 of file MOM_oda_driver.F90.

434  type(time_type), intent(in) :: Time !< the current model time
435  type(oda_CS), intent(inout) :: CS !< the ocean DA control structure
436 
437  integer :: i, j
438  integer :: m
439  integer :: yr, mon, day, hr, min, sec
440 
441  if ( time >= cs%Time ) then
442 
443  !! switch to global pelist
444  call set_current_pelist(cs%filter_pelist)
445 
446  call get_profiles(time, cs%Profiles, cs%CProfiles)
447 #ifdef ENABLE_ECDA
448  call ensemble_filter(cs%Ocean_prior, cs%Ocean_posterior, cs%CProfiles, cs%kdroot, cs%mpp_domain, cs%oda_grid)
449 #endif
450 
451  !! switch back to ensemble member pelist
452  call set_current_pelist(cs%ensemble_pelist(cs%ensemble_id,:))
453 
454  endif
455 
456  return

◆ oda_end()

subroutine, public mom_oda_driver_mod::oda_end ( type(oda_cs), intent(inout)  CS)

Finalize DA module.

Parameters
[in,out]csthe ocean DA control structure

Definition at line 461 of file MOM_oda_driver.F90.

461  type(ODA_CS), intent(inout) :: CS !< the ocean DA control structure
462 

◆ save_obs_diff()

subroutine, public mom_oda_driver_mod::save_obs_diff ( character(len=*), intent(in)  filename,
type(oda_cs), intent(in), pointer  CS 
)

Write observation differences to a file.

Parameters
[in]filenamename of output file
[in]cspointer to DA control structure

Definition at line 521 of file MOM_oda_driver.F90.

521  character(len=*), intent(in) :: filename !< name of output file
522  type(ODA_CS), pointer, intent(in) :: CS !< pointer to DA control structure
523 
524  integer :: fid ! profile file handle
525  type(ocean_profile_type), pointer :: Prof=>null()
526 
527  fid = open_profile_file(trim(filename), nvar=2, thread=mpp_single, fset=mpp_single)
528  prof=>cs%CProfiles
529 
530  !! switch to global pelist
531  !call set_current_pelist(CS%filter_pelist)
532 
533  do while (associated(prof))
534  call write_profile(fid,prof)
535  prof=>prof%cnext
536  enddo
537  call close_profile_file(fid)
538 
539  !! switch back to ensemble member pelist
540  !call set_current_pelist(CS%ensemble_pelist(CS%ensemble_id,:))
541 
542  return

◆ set_analysis_time()

subroutine, public mom_oda_driver_mod::set_analysis_time ( type(time_type), intent(in)  Time,
type(oda_cs), intent(inout), pointer  CS 
)

Set the next analysis time.

Parameters
[in]timethe current model time
[in,out]csthe DA control structure

Definition at line 494 of file MOM_oda_driver.F90.

494  type(time_type), intent(in) :: Time !< the current model time
495  type(ODA_CS), pointer, intent(inout) :: CS !< the DA control structure
496 
497  character(len=160) :: mesg ! The text of an error message
498  integer :: yr, mon, day, hr, min, sec
499 
500  if (time >= cs%Time) then
501  cs%Time=increment_time(cs%Time,cs%assim_frequency*3600)
502 
503  call get_date(time, yr, mon, day, hr, min, sec)
504  write(mesg,*) 'Model Time: ', yr, mon, day, hr, min, sec
505  call mom_mesg("set_analysis_time: "//trim(mesg))
506  call get_date(cs%time, yr, mon, day, hr, min, sec)
507  write(mesg,*) 'Assimilation Time: ', yr, mon, day, hr, min, sec
508  call mom_mesg("set_analysis_time: "//trim(mesg))
509  endif
510  if (cs%Time < time) then
511  call mom_error(fatal, " set_analysis_time: " // &
512  "assimilation interval appears to be shorter than " // &
513  "the model timestep")
514  endif
515  return
516 

References mom_error_handler::mom_error(), and mom_error_handler::mom_mesg().

Referenced by mom::step_mom().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ set_prior_tracer()

subroutine, public mom_oda_driver_mod::set_prior_tracer ( type(time_type), intent(in)  Time,
type(ocean_grid_type), pointer  G,
type(verticalgrid_type), intent(in)  GV,
real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in)  h,
type(thermo_var_ptrs), intent(in)  tv,
type(oda_cs), pointer  CS 
)

Copy ensemble member tracers to ensemble vector.

Parameters
[in]timeThe current model time
gdomain and grid information for ocean model
[in]gvThe ocean's vertical grid structure
[in]hLayer thicknesses [H ~> m or kg m-2]
[in]tvA structure pointing to various thermodynamic variables
csocean DA control structure

Definition at line 323 of file MOM_oda_driver.F90.

323  type(time_type), intent(in) :: Time !< The current model time
324  type(ocean_grid_type), pointer :: G !< domain and grid information for ocean model
325  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
326  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
327  type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various thermodynamic variables
328 
329  type(ODA_CS), pointer :: CS !< ocean DA control structure
330  real, dimension(:,:,:), allocatable :: T, S
331  type(ocean_grid_type), pointer :: Grid=>null()
332  integer :: i,j, m, n, ss
333  integer :: is, ie, js, je
334  integer :: isc, iec, jsc, jec
335  integer :: isd, ied, jsd, jed
336  integer :: id
337  logical :: used
338 
339  ! return if not time for analysis
340  if (time < cs%Time) return
341 
342  if (.not. associated(cs%Grid)) call mom_error(fatal,'ODA_CS ensemble horizontal grid not associated')
343  if (.not. associated(cs%GV)) call mom_error(fatal,'ODA_CS ensemble vertical grid not associated')
344 
345  !! switch to global pelist
346  call set_current_pelist(cs%filter_pelist)
347  call mom_mesg('Setting prior')
348 
349  isc=cs%Grid%isc;iec=cs%Grid%iec;jsc=cs%Grid%jsc;jec=cs%Grid%jec
350  call mpp_get_compute_domain(cs%domains(cs%ensemble_id)%mpp_domain,is,ie,js,je)
351  call mpp_get_data_domain(cs%domains(cs%ensemble_id)%mpp_domain,isd,ied,jsd,jed)
352  allocate(t(isd:ied,jsd:jed,cs%nk))
353  allocate(s(isd:ied,jsd:jed,cs%nk))
354 
355  do j=js,je; do i=is,ie
356  call remapping_core_h(cs%remapCS, gv%ke, h(i,j,:), tv%T(i,j,:), &
357  cs%nk, cs%h(i,j,:), t(i,j,:))
358  call remapping_core_h(cs%remapCS, gv%ke, h(i,j,:), tv%S(i,j,:), &
359  cs%nk, cs%h(i,j,:), s(i,j,:))
360  enddo ; enddo
361 
362  do m=1,cs%ensemble_size
363  call mpp_redistribute(cs%domains(m)%mpp_domain, t,&
364  cs%mpp_domain, cs%Ocean_prior%T(:,:,:,m), complete=.true.)
365  call mpp_redistribute(cs%domains(m)%mpp_domain, s,&
366  cs%mpp_domain, cs%Ocean_prior%S(:,:,:,m), complete=.true.)
367  if (cs%Ocean_prior%id_t(m)>0) &
368  used=send_data(cs%Ocean_prior%id_t(m), cs%Ocean_prior%T(isc:iec,jsc:jec,:,m), cs%Time)
369  if (cs%Ocean_prior%id_s(m)>0) &
370  used=send_data(cs%Ocean_prior%id_s(m), cs%Ocean_prior%S(isc:iec,jsc:jec,:,m), cs%Time)
371  enddo
372  deallocate(t,s)
373 
374  !! switch back to ensemble member pelist
375  call set_current_pelist(cs%ensemble_pelist(cs%ensemble_id,:))
376 
377  return
378 

References mom_error_handler::mom_error(), mom_error_handler::mom_mesg(), and mom_remapping::remapping_core_h().

Referenced by mom::step_mom().

Here is the call graph for this function:
Here is the caller graph for this function:

Variable Documentation

◆ eakf_assim

integer, parameter mom_oda_driver_mod::eakf_assim =2
private

DA parameters.

Definition at line 108 of file MOM_oda_driver.F90.

Referenced by init_oda().

◆ no_assim

integer, parameter mom_oda_driver_mod::no_assim = 0
private

DA parameters.

Definition at line 108 of file MOM_oda_driver.F90.

108 integer, parameter :: NO_ASSIM = 0, oi_assim=1, eakf_assim=2

Referenced by init_oda().

◆ oi_assim

integer, parameter mom_oda_driver_mod::oi_assim =1
private

DA parameters.

Definition at line 108 of file MOM_oda_driver.F90.

Referenced by init_oda().