PDAF-OMI Debugging Information

Overview

When implementing an observation with PDAF, or when performing the very first implementation of PDAF with a new model, it is useful to check whether the inputs to the PDAF-routines are correctly used. For this porpose, PDAF-OMI provides a debugging functionality. It allows you to activate debugging output e.g. for a single local analysis domain on a single process of a complex application of a local filter like LEKSTF.

Activating Debugging Output

Debugging output is activated by the routine PDAFomi_set_debug_flag. In particular this call can be inserted in any routines contained in callback_obs_pdafomi.F90.

For example to activate debugging in init_dim_obs_l_pdafomi for the local analysis domain domain_p=10 and filter process 0, one uses

SUBROUTINE init_dim_obs_l_pdafomi(domain_p, step, dim_obs, dim_obs_l)

  USE PDAFomi, ONLY: PDAFomi_set_debug_flag
  USE mod_parallel_pdaf, ONLY: mype_filter

  ...

  IF (domain_p==10 .AND. mype_filter==0) THEN
    CALL PDAFomi_set_debug_flag(domain_p)
  ELSE
    CALL PDAFomi_set_debug_flag(0)
  ENDIF

  ...

Understanding the Debugging Output

The debugging output mainly writes information about the different variables contained in the full data type obs_f allocated as thisobs and the local type obs_l allocate as thisobs_l. For reference we list the full declaration of these types. When reading the debugging output one can check for the meaning of the variables.

  TYPE obs_f
     ! ---- Mandatory variables to be set in INIT_DIM_OBS ----
     INTEGER :: doassim=0                 !< Whether to assimilate this observation type
     INTEGER :: disttype                  !< Type of distance computation to use for localization
     INTEGER :: ncoord                    !< Number of coordinates use for distance computation
     INTEGER, ALLOCATABLE :: id_obs_p(:,:) !< Indices of process-local observed field in state vector

     ! ---- Optional variables - they can be set in INIT_DIM_OBS ----
     REAL, ALLOCATABLE :: icoeff_p(:,:)   !< Interpolation coefficients for obs. operator (optional)
     REAL, ALLOCATABLE :: domainsize(:)   !< Size of domain for periodicity (<=0 for no periodicity) (optional)

     ! ---- Variables with predefined values - they can be changed in INIT_DIM_OBS  ----
     INTEGER :: obs_err_type=0            !< Type of observation error: (0) Gauss, (1) Laplace
     INTEGER :: use_global_obs=1          !< Whether to use (1) global full obs. 
                                          !< or (0) obs. restricted to those relevant for a process domain

     ! ----  The following variables are set in the routine PDAFomi_gather_obs ---
     INTEGER :: dim_obs_p                 !< number of PE-local observations
     INTEGER :: dim_obs_f                 !< number of full observations
     INTEGER :: dim_obs_g                 !< global number of observations
     INTEGER :: off_obs_f                 !< Offset of this observation in overall full obs. vector
     INTEGER :: off_obs_g                 !< Offset of this observation in overall global obs. vector
     INTEGER :: obsid                     !< Index of observation over all assimilated observations
     REAL, ALLOCATABLE :: obs_f(:)        !< Full observed field
     REAL, ALLOCATABLE :: ocoord_f(:,:)   !< Coordinates of full observation vector
     REAL, ALLOCATABLE :: ivar_obs_f(:)   !< Inverse variance of full observations
     INTEGER, ALLOCATABLE :: id_obs_f_lim(:) !< Indices of domain-relevant full obs. in global vector of obs.
  END TYPE obs_f

  TYPE obs_l
     INTEGER :: dim_obs_l                 !< number of local observations
     INTEGER :: off_obs_l                 !< Offset of this observation in overall local obs. vector
     INTEGER, ALLOCATABLE :: id_obs_l(:)  !< Indices of local observations in full obs. vector 
     REAL, ALLOCATABLE :: distance_l(:)   !< Distances of local observations
     REAL, ALLOCATABLE :: ivar_obs_l(:)   !< Inverse variance of local observations
     INTEGER :: locweight                 !< Specify localization function
     REAL :: lradius                      !< localization radius
     REAL :: sradius                      !< support radius for localization function
  END TYPE obs_l

Example output

For illustration we insert the above call to PDAFomi_set_debug flag into init_dim_obs_l_pdafomi in tutorial/online_2D_serialmodel_omi/. then we compile and execute the tutorial program with:

  mpirun -np 4 ./model_pdaf -dim_ens 4 -filtertype 7 -local_range 5 -locweight 2 -assim_B .true. -assim_A .false.

With these settings only observation type B is activated, which are just 3 values in the model domain

The first lines of the debugging output looks like this:

 ++ OMI-debug set_debug_flag: mype_filter           0 activate          10
 ++ OMI-debug:           10 PDAFomi_init_dim_obs_l -- START
 ++ OMI-debug:           10    PDAFomi_init_dim_obs_l -- count local observations
 ++ OMI-debug init_dim_obs_l:          10   Re-init dim_obs_l=0
 ++ OMI-debug init_dim_obs_l:          10   coords_l   1.0000000000000000        10.000000000000000     
 ++ OMI-debug cnt_dim_obs_l:           10   thisobs%ncoord           2
 ++ OMI-debug cnt_dim_obs_l:           10   thisobs_l%lradius   5.0000000000000000     
 ++ OMI-debug cnt_dim_obs_l:           10   Check for observations within radius
 ++ OMI-debug comp_dist2:              10   compute Cartesian distance
 ++ OMI-debug cnt_dim_obs_l:           10   valid observation with coordinates   5.0000000000000000        8.0000000000000000     
 ++ OMI-debug:           10    PDAFomi_init_dim_obs_l -- initialize local observation arrays
 ++ OMI-debug comp_dist2:              10   compute Cartesian distance
 ++ OMI-debug init_dim_obs_l:          10   thisobs_l%dim_obs_l           1
 ++ OMI-debug init_dim_obs_l:          10   thisobs_l%id_obs_l           1
 ++ OMI-debug init_dim_obs_l:          10   thisobs_l%distance_l   4.4721359549995796     
 ++ OMI-debug:           10 PDAFomi_init_dim_obs_l -- END

The first line

 ++ OMI-debug set_debug_flag: mype_filter           0 activate          10

is from PDAFomi_set_debug_flag showing that debugging is activates with value 10 (which is values fo domain_p specified in the call)

The next lines are

 ++ OMI-debug:           10 PDAFomi_init_dim_obs_l -- START           
 ++ OMI-debug:           10    PDAFomi_init_dim_obs_l -- count local observations
 ++ OMI-debug init_dim_obs_l:          10   Re-init dim_obs_l=0

show that debugging output for PDAFomi_init_dim_obs_l is shown. Most routines show such a 'START' line. The second line shows that a segment of the orutine started, the counting of local observations. The third line states that dim_obs_l=0 is set. This line, as many others shows the name of the subroutine in short form without PDAFomi at the beginning of the line.

The following lines show variable values

 ++ OMI-debug init_dim_obs_l:          10   coords_l   1.0000000000000000        10.000000000000000     
 ++ OMI-debug cnt_dim_obs_l:           10   thisobs%ncoord           2
 ++ OMI-debug cnt_dim_obs_l:           10   thisobs_l%lradius   5.0000000000000000     

First we see that the coordinates coords_l of the grid point corresponding to domain_p=10 are (1.0, 10.0). Further we have two dimensions (thisobs%ncoord=2) and the localization radius is set to 5.0. In the following lines

 ++ OMI-debug cnt_dim_obs_l:           10   Check for observations within radius
 ++ OMI-debug comp_dist2:              10   compute Cartesian distance
 ++ OMI-debug cnt_dim_obs_l:           10   valid observation with coordinates   5.0000000000000000        8.0000000000000000     

the it is checked which observations lie within the distance thisobs_l%lradius=5.0 from coords_l=(1.0, 5.0). One observation with coordinates (5.0, 8.0) is found. The followign lines

 ++ OMI-debug:           10    PDAFomi_init_dim_obs_l -- initialize local observation arrays
 ++ OMI-debug comp_dist2:              10   compute Cartesian distance
 ++ OMI-debug init_dim_obs_l:          10   thisobs_l%dim_obs_l           1
 ++ OMI-debug init_dim_obs_l:          10   thisobs_l%id_obs_l           1
 ++ OMI-debug init_dim_obs_l:          10   thisobs_l%distance_l   4.4721359549995796     

show that observation arrays are initialized. A Cartesian distance is computerd. thisobs_l%dim_obs_l confirms that one local observation was found and it's the first element of the full observation vector (thisobs_l%id_obs_l=1) and the distance is thisobs_l%distance_l=4.4721359549995796.

One could now compare this with the input information. Are the values of coords_l correct? Are the coordinates of the first observation (5.0, 8.0)?

Later in the debugging output (not shown above) we find, for example,

 ++ OMI-debug:           10    PDAFomi_init_obs_l -- Get local vector of observations
 ++ OMI-debug g2l_obs:                 10   thisobs%id_obs_l           1
 ++ OMI-debug g2l_obs:                 10   obs_l -0.98253999999999997     
 ++ OMI-debug:           10    PDAFomi_init_obs_l -- Get local vector of inverse obs. variances
 ++ OMI-debug g2l_obs:                 10   thisobs%id_obs_l           1
 ++ OMI-debug g2l_obs:                 10   obs_l   4.0000000000000000     

These lines are from the internal routine PDAFomi_init_obs_l, which initializes the local vector of observations. (This is functionality also existing in the 'traditional' form of implementing the observation. The documention page about init_obs_l_pdaf explains what is done here.) Important is that not only the local observation vector (a single value) is initialized, but also the inverse observation variance (4.0 in this example).

Similarly

 ++ OMI-debug:           10 PDAFomi_g2l_obs -- START Get local observed ensemble member           1
 ++ OMI-debug g2l_obs:                 10   thisobs%id_obs_l           1
 ++ OMI-debug g2l_obs:                 10   obs_l  0.12725110911115362     
 ++ OMI-debug:           10 PDAFomi_g2l_obs -- END

shows the localization of ensemble member 1. This is performed in the interval routine PDAFomi_g2l_obs (see the page about g2l_obs_pdaf). There is analogous output for all of the 4 ensemble states.