! Copyright (c) 2004-2024 Lars Nerger ! ! This file is part of PDAF. ! ! PDAF is free software: you can redistribute it and/or modify ! it under the terms of the GNU Lesser General Public License ! as published by the Free Software Foundation, either version ! 3 of the License, or (at your option) any later version. ! ! PDAF is distributed in the hope that it will be useful, ! but WITHOUT ANY WARRANTY; without even the implied warranty of ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ! GNU Lesser General Public License for more details. ! ! You should have received a copy of the GNU Lesser General Public ! License along with PDAF. If not, see <http://www.gnu.org/licenses/>. ! !$Id: PDAFomi_callback.F90 1147 2023-03-12 16:14:34Z lnerger $ !> Generic PDAFomi callback routines !! !! This file provides generic call-back routines for OMI. The OMI structure !! allow us to provide several of the observation-related call-back !! routines for PDAF in a generic form. These routines are collected in !! this file. !! !! The routines here are mainly pure pass-through routines. Thus they !! simply call one of the routines from PDAF-OMI. Partly some addtional !! variable is required, e.g. to specify the offset of an observation !! in the observation vector containing all observation types. These !! cases are described in the routines. !! !! __Revision history:__ !! * 2020-06 - Lars Nerger - Initial code splitting callback_obs_pdafomi.F90 !! * Later revisions - see repository log !! !------------------------------------------------------------------------------- !> Call-back routine for init_obs_f !! !! This routine calls the routine PDAFomi_init_obs_f !! for each observation type !! SUBROUTINE PDAFomi_init_obs_f_cb(step, dim_obs_f, observation_f) ! Include overall pointer to observation variables use PDAFomi, only: n_obstypes, obs_f_all ! Include PDAFomi function USE PDAFomi, ONLY: PDAFomi_init_obs_f IMPLICIT NONE ! *** Arguments *** INTEGER, INTENT(in) :: step !< Current time step INTEGER, INTENT(in) :: dim_obs_f !< Dimension of full observation vector REAL, INTENT(out) :: observation_f(dim_obs_f) !< Full observation vector ! *** local variables *** INTEGER :: i ! Loop counter INTEGER :: offset_obs_f ! Count offset of an observation type in full obs. vector INTEGER :: idummy ! Dummy to prevent compiler warning ! ****************************************** ! *** Initialize full observation vector *** ! ****************************************** ! Initialize dummy to prevent compiler warning idummy = step ! Initialize offset (it will be incremented in PDAFomi_init_obs_f) offset_obs_f = 0 ! The order of the calls has to be consistent with those in obs_op_f_pdafomi DO i=1, n_obstypes CALL PDAFomi_init_obs_f(obs_f_all(i)%ptr, dim_obs_f, observation_f, offset_obs_f) END DO END SUBROUTINE PDAFomi_init_obs_f_cb !------------------------------------------------------------------------------- !> Call-back routine for init_obsvar !! !! This routine calls the routine PDAFomi_init_obsvar_f !! for each observation type !! SUBROUTINE PDAFomi_init_obsvar_cb(step, dim_obs_p, obs_p, meanvar) ! Include overall pointer to observation variables use PDAFomi, only: n_obstypes, obs_f_all ! Include PDAFomi function USE PDAFomi, ONLY: PDAFomi_init_obsvar_f IMPLICIT NONE ! *** Arguments *** INTEGER, INTENT(in) :: step !< Current time step INTEGER, INTENT(in) :: dim_obs_p !< PE-local dimension of observation vector REAL, INTENT(in) :: obs_p(dim_obs_p) !< PE-local observation vector REAL, INTENT(out) :: meanvar !< Mean observation error variance ! *** Local variables *** INTEGER :: i ! Loop counter INTEGER :: cnt_obs_f ! Count observations for offset INTEGER :: idummy ! Dummy to prevent compiler warning REAL :: rdummy ! Dummy to access obs ! ***************************** ! *** Compute mean variance *** ! ***************************** ! Initialize dummy to prevent compiler warning idummy = step rdummy = obs_p(1) ! Initialize observation counter (it will be incremented in PDAFomi_init_obsvar_f) cnt_obs_f = 0 DO i=1, n_obstypes CALL PDAFomi_init_obsvar_f(obs_f_all(i)%ptr, meanvar, cnt_obs_f) END DO END SUBROUTINE PDAFomi_init_obsvar_cb !------------------------------------------------------------------------------- !> Call-back routine for g2l_obs !! !! This routine calls the routine PDAFomi_g2l_obs !! for each observation type !! SUBROUTINE PDAFomi_g2l_obs_cb(domain_p, step, dim_obs_f, dim_obs_l, ostate_f, & ostate_l) ! Include overall pointer to observation variables use PDAFomi, only: n_obstypes, obs_f_all, obs_l_all ! Include PDAFomi function USE PDAFomi, ONLY: PDAFomi_g2l_obs IMPLICIT NONE ! *** Arguments *** INTEGER, INTENT(in) :: domain_p !< Index of current local analysis domain INTEGER, INTENT(in) :: step !< Current time step INTEGER, INTENT(in) :: dim_obs_f !< Dimension of full PE-local observation vector INTEGER, INTENT(in) :: dim_obs_l !< Dimension of local observation vector REAL, INTENT(in) :: ostate_f(dim_obs_f) !< Full PE-local obs.ervation vector REAL, INTENT(out) :: ostate_l(dim_obs_l) !< Observation vector on local domain ! *** local variables *** INTEGER :: i ! Loop counter INTEGER :: idummy ! Dummy to prevent compiler warning ! Initialize dummy to prevent compiler warning idummy = step idummy = domain_p ! ******************************************************* ! *** Perform localization of some observation vector *** ! *** to the current local analysis domain. *** ! ******************************************************* DO i=1, n_obstypes CALL PDAFomi_g2l_obs(obs_l_all(i)%ptr, obs_f_all(i)%ptr, ostate_f, ostate_l) END DO END SUBROUTINE PDAFomi_g2l_obs_cb !------------------------------------------------------------------------------- !> Call-back routine for init_obs_l !! !! This routine calls the routine PDAFomi_init_obs_l !! for each observation type !! SUBROUTINE PDAFomi_init_obs_l_cb(domain_p, step, dim_obs_l, observation_l) ! Include overall pointer to observation variables use PDAFomi, only: n_obstypes, obs_f_all, obs_l_all ! Include PDAFomi function USE PDAFomi, ONLY: PDAFomi_init_obs_l IMPLICIT NONE ! *** Arguments *** INTEGER, INTENT(in) :: domain_p !< Index of current local analysis domain index INTEGER, INTENT(in) :: step !< Current time step INTEGER, INTENT(in) :: dim_obs_l !< Local dimension of observation vector REAL, INTENT(out) :: observation_l(dim_obs_l) !< Local observation vector ! *** local variables *** INTEGER :: i ! Loop counter INTEGER :: idummy ! Dummy to prevent compiler warning ! Initialize dummy to prevent compiler warning idummy = step idummy = domain_p ! ******************************************* ! *** Initialize local observation vector *** ! ******************************************* DO i=1, n_obstypes CALL PDAFomi_init_obs_l(obs_l_all(i)%ptr, obs_f_all(i)%ptr, observation_l) END DO END SUBROUTINE PDAFomi_init_obs_l_cb !------------------------------------------------------------------------------- !> Call-back routine for init_obsvar_l !! !! This routine calls the routine PDAFomi_init_obsvar_l !! for each observation type !! SUBROUTINE PDAFomi_init_obsvar_l_cb(domain_p, step, dim_obs_l, obs_l, meanvar_l) ! Include overall pointer to observation variables use PDAFomi, only: n_obstypes, obs_f_all, obs_l_all ! Include PDAFomi function USE PDAFomi, ONLY: PDAFomi_init_obsvar_l IMPLICIT NONE ! *** Arguments *** INTEGER, INTENT(in) :: domain_p !< Index of current local analysis domain INTEGER, INTENT(in) :: step !< Current time step INTEGER, INTENT(in) :: dim_obs_l !< Local dimension of observation vector REAL, INTENT(in) :: obs_l(dim_obs_l) !< Local observation vector REAL, INTENT(out) :: meanvar_l !< Mean local observation error variance ! *** Local variables *** INTEGER :: i ! Loop counter INTEGER :: cnt_obs_l ! Count local observations for offset INTEGER :: idummy ! Dummy to prevent compiler warning REAL :: rdummy ! Dummy to prevent compiler warning ! Initialize dummy to prevent compiler warning idummy = step idummy = domain_p rdummy = obs_l(1) ! *********************************** ! *** Compute local mean variance *** ! *********************************** ! Initialize observation counter (it will be incremented in PDAFomi_init_obsvar_f) cnt_obs_l = 0 ! The order of the calls is not relevant DO i=1, n_obstypes CALL PDAFomi_init_obsvar_l(obs_l_all(i)%ptr, obs_f_all(i)%ptr, meanvar_l, cnt_obs_l) END DO END SUBROUTINE PDAFomi_init_obsvar_l_cb !------------------------------------------------------------------------------- !> Call-back routine for prodRinvA_l !! !! This routine calls the routine PDAFomi_prodRinvA_l !! for each observation type !! SUBROUTINE PDAFomi_prodRinvA_l_cb(domain_p, step, dim_obs_l, rank, obs_l, A_l, C_l) ! Include overall pointer to observation variables use PDAFomi, only: n_obstypes, obs_f_all, obs_l_all ! Include PDAFomi function USE PDAFomi, ONLY: PDAFomi_prodRinvA_l ! Include filter process rank USE PDAF_mod_filterMPI, ONLY: mype_filter ! Include verbosity information USE PDAF_mod_filter, ONLY: screen #if defined (_OPENMP) ! Include OpenMP function to determine verbosity for OpenMP USE omp_lib, ONLY: omp_get_thread_num #endif IMPLICIT NONE ! *** Arguments *** INTEGER, INTENT(in) :: domain_p !< Index of current local analysis domain INTEGER, INTENT(in) :: step !< Current time step INTEGER, INTENT(in) :: dim_obs_l !< Dimension of local observation vector INTEGER, INTENT(in) :: rank !< Rank of initial covariance matrix REAL, INTENT(in) :: obs_l(dim_obs_l) !< Local vector of observations REAL, INTENT(inout) :: A_l(dim_obs_l, rank) !< Input matrix REAL, INTENT(out) :: C_l(dim_obs_l, rank) !< Output matrix ! *** local variables *** INTEGER :: i ! Loop counter INTEGER :: verbose ! Verbosity flag INTEGER, SAVE :: domain_save = -1 ! Save previous domain index INTEGER, SAVE :: mythread ! Thread variable for OpenMP INTEGER :: idummy ! Dummy to prevent compiler warning REAL :: rdummy ! Dummy to prevent compiler warning !$OMP THREADPRIVATE(mythread, domain_save) ! Initialize dummy to prevent compiler warning idummy = step rdummy = obs_l(1) ! ********************** ! *** INITIALIZATION *** ! ********************** ! For OpenMP parallelization, determine the thread index #if defined (_OPENMP) mythread = omp_get_thread_num() #else mythread = 0 #endif ! Set verbosity flag (Screen output for first analysis domain) IF (screen > 0) THEN IF ((domain_p <= domain_save .OR. domain_save < 0) .AND. mype_filter==0) THEN verbose = 1 ! In case of OpenMP, let only thread 0 write output to the screen IF (mythread>0) verbose = 0 ELSE verbose = 0 END IF ELSE verbose = 0 END IF domain_save = domain_p ! ************************************* ! *** Compute *** ! *** -1 *** ! *** C = W R A *** ! *** *** ! *** where W are the localization *** ! *** weights. *** ! ************************************* DO i=1, n_obstypes CALL PDAFomi_prodRinvA_l(obs_l_all(i)%ptr, obs_f_all(i)%ptr, dim_obs_l, rank, & A_l, C_l, verbose) END DO END SUBROUTINE PDAFomi_prodRinvA_l_cb !------------------------------------------------------------------------------- !> Call-back routine for likelihood_l !! !! This routine calls the routine PDAFomi_likelihood_l !! for each observation type !! SUBROUTINE PDAFomi_likelihood_l_cb(domain_p, step, dim_obs_l, obs_l, resid_l, lhood_l) ! Include overall pointer to observation variables use PDAFomi, only: n_obstypes, obs_f_all, obs_l_all ! Include PDAFomi function USE PDAFomi, ONLY: PDAFomi_likelihood_l ! Include filter process rank USE PDAF_mod_filterMPI, ONLY: mype_filter ! Include verbosity information USE PDAF_mod_filter, ONLY: screen #if defined (_OPENMP) ! Include OpenMP function to determine verbosity for OpenMP USE omp_lib, ONLY: omp_get_thread_num #endif IMPLICIT NONE ! *** Arguments *** INTEGER, INTENT(in) :: domain_p ! Current local analysis domain INTEGER, INTENT(in) :: step !< Current time step INTEGER, INTENT(in) :: dim_obs_l !< PE-local dimension of obs. vector REAL, INTENT(in) :: obs_l(dim_obs_l) !< PE-local vector of observations REAL, INTENT(inout) :: resid_l(dim_obs_l) !< Input vector of residuum REAL, INTENT(out) :: lhood_l !< Output vector - log likelihood ! *** local variables *** INTEGER :: i ! Loop counter INTEGER :: verbose ! Verbosity flag INTEGER, SAVE :: domain_save = -1 ! Save previous domain index INTEGER, SAVE :: mythread ! Thread variable for OpenMP INTEGER :: idummy ! Dummy to prevent compiler warning REAL :: rdummy ! Dummy to prevent compiler warning !$OMP THREADPRIVATE(mythread, domain_save) ! Initialize dummy to prevent compiler warning idummy = step rdummy = obs_l(1) ! ********************** ! *** INITIALIZATION *** ! ********************** ! For OpenMP parallelization, determine the thread index #if defined (_OPENMP) mythread = omp_get_thread_num() #else mythread = 0 #endif ! Set verbosity flag (Screen output for first analysis domain) IF (screen > 0) THEN IF ((domain_p < domain_save .OR. domain_save < 0) .AND. mype_filter==0) THEN verbose = 1 ! In case of OpenMP, let only thread 0 write output to the screen IF (mythread>0) verbose = 0 ELSE verbose = 0 END IF ELSE verbose = 0 END IF domain_save = domain_p ! ******************************** ! *** Compute local likelihood *** ! ******************************** ! Initialize likelihood value before starting computation lhood_l = 0.0 ! Increment likelihood DO i=1, n_obstypes CALL PDAFomi_likelihood_l(obs_l_all(i)%ptr, obs_f_all(i)%ptr, resid_l, & lhood_l, verbose) END DO END SUBROUTINE PDAFomi_likelihood_l_cb !------------------------------------------------------------------------------- !> Call-back routine for prodRinvA_hyb_l !! !! This routine calls the routine PDAFomi_prodRinvA_hyb_l !! for each observation type !! SUBROUTINE PDAFomi_prodRinvA_hyb_l_cb(domain_p, step, dim_obs_l, rank, obs_l, alpha, A_l, C_l) ! Include overall pointer to observation variables use PDAFomi, only: n_obstypes, obs_f_all, obs_l_all ! Include PDAFomi function USE PDAFomi, ONLY: PDAFomi_prodRinvA_hyb_l ! Include filter process rank USE PDAF_mod_filterMPI, ONLY: mype_filter ! Include verbosity information USE PDAF_mod_filter, ONLY: screen #if defined (_OPENMP) ! Include OpenMP function to determine verbosity for OpenMP USE omp_lib, ONLY: omp_get_thread_num #endif IMPLICIT NONE ! *** Arguments *** INTEGER, INTENT(in) :: domain_p !< Index of current local analysis domain INTEGER, INTENT(in) :: step !< Current time step INTEGER, INTENT(in) :: dim_obs_l !< Dimension of local observation vector INTEGER, INTENT(in) :: rank !< Rank of initial covariance matrix REAL, INTENT(in) :: obs_l(dim_obs_l) !< Local vector of observations REAL, INTENT(in) :: alpha !< Hybrid weight REAL, INTENT(inout) :: A_l(dim_obs_l, rank) !< Input matrix REAL, INTENT(out) :: C_l(dim_obs_l, rank) !< Output matrix ! *** local variables *** INTEGER :: i ! Loop counter INTEGER :: verbose ! Verbosity flag INTEGER, SAVE :: domain_save = -1 ! Save previous domain index INTEGER, SAVE :: mythread ! Thread variable for OpenMP INTEGER :: idummy ! Dummy to prevent compiler warning REAL :: rdummy ! Dummy to prevent compiler warning !$OMP THREADPRIVATE(mythread, domain_save) ! Initialize dummy to prevent compiler warning idummy = step rdummy = obs_l(1) ! ********************** ! *** INITIALIZATION *** ! ********************** ! For OpenMP parallelization, determine the thread index #if defined (_OPENMP) mythread = omp_get_thread_num() #else mythread = 0 #endif ! Set verbosity flag (Screen output for first analysis domain) IF (screen > 0) THEN IF ((domain_p <= domain_save .OR. domain_save < 0) .AND. mype_filter==0) THEN verbose = 1 ! In case of OpenMP, let only thread 0 write output to the screen IF (mythread>0) verbose = 0 ELSE verbose = 0 END IF ELSE verbose = 0 END IF domain_save = domain_p ! ************************************* ! *** Compute *** ! *** -1 *** ! *** C = alpha W R A *** ! *** *** ! *** where W are the localization *** ! *** weights. *** ! ************************************* DO i=1, n_obstypes CALL PDAFomi_prodRinvA_hyb_l(obs_l_all(i)%ptr, obs_f_all(i)%ptr, dim_obs_l, rank, & alpha, A_l, C_l, verbose) END DO END SUBROUTINE PDAFomi_prodRinvA_hyb_l_cb !------------------------------------------------------------------------------- !> Call-back routine for likelihood_hyb_l !! !! This routine calls the routine PDAFomi_likelihood_hyb_l !! for each observation type !! SUBROUTINE PDAFomi_likelihood_hyb_l_cb(domain_p, step, dim_obs_l, obs_l, resid_l, alpha, lhood_l) ! Include overall pointer to observation variables use PDAFomi, only: n_obstypes, obs_f_all, obs_l_all ! Include PDAFomi function USE PDAFomi, ONLY: PDAFomi_likelihood_hyb_l ! Include filter process rank USE PDAF_mod_filterMPI, ONLY: mype_filter ! Include verbosity information USE PDAF_mod_filter, ONLY: screen #if defined (_OPENMP) ! Include OpenMP function to determine verbosity for OpenMP USE omp_lib, ONLY: omp_get_thread_num #endif IMPLICIT NONE ! *** Arguments *** INTEGER, INTENT(in) :: domain_p !< Current local analysis domain INTEGER, INTENT(in) :: step !< Current time step INTEGER, INTENT(in) :: dim_obs_l !< PE-local dimension of obs. vector REAL, INTENT(in) :: obs_l(dim_obs_l) !< PE-local vector of observations REAL, INTENT(inout) :: resid_l(dim_obs_l) !< Input vector of residuum REAL, INTENT(in) :: alpha !< Hybrid weight REAL, INTENT(out) :: lhood_l !< Output vector - log likelihood ! *** local variables *** INTEGER :: i ! Loop counter INTEGER :: verbose ! Verbosity flag INTEGER, SAVE :: domain_save = -1 ! Save previous domain index INTEGER, SAVE :: mythread ! Thread variable for OpenMP INTEGER :: idummy ! Dummy to prevent compiler warning REAL :: rdummy ! Dummy to prevent compiler warning !$OMP THREADPRIVATE(mythread, domain_save) ! Initialize dummy to prevent compiler warning idummy = step rdummy = obs_l(1) ! ********************** ! *** INITIALIZATION *** ! ********************** ! For OpenMP parallelization, determine the thread index #if defined (_OPENMP) mythread = omp_get_thread_num() #else mythread = 0 #endif ! Set verbosity flag (Screen output for first analysis domain) IF (screen > 0) THEN IF ((domain_p < domain_save .OR. domain_save < 0) .AND. mype_filter==0) THEN verbose = 1 ! In case of OpenMP, let only thread 0 write output to the screen IF (mythread>0) verbose = 0 ELSE verbose = 0 END IF ELSE verbose = 0 END IF domain_save = domain_p ! ******************************** ! *** Compute local likelihood *** ! ******************************** ! Initialize likelihood value before starting computation lhood_l = 0.0 ! Increment likelihood DO i=1, n_obstypes CALL PDAFomi_likelihood_hyb_l(obs_l_all(i)%ptr, obs_f_all(i)%ptr, resid_l, & alpha, lhood_l, verbose) END DO END SUBROUTINE PDAFomi_likelihood_hyb_l_cb !------------------------------------------------------------------------------- !> Call-back routine for prodRinvA !! !! This routine calls the routine PDAFomi_prodRinvA !! for each observation type !! SUBROUTINE PDAFomi_prodRinvA_cb(step, dim_obs_p, ncol, obs_p, A_p, C_p) ! Include overall pointer to observation variables use PDAFomi, only: n_obstypes, obs_f_all ! Include PDAFomi function USE PDAFomi, ONLY: PDAFomi_prodRinvA IMPLICIT NONE ! *** Arguments *** INTEGER, INTENT(in) :: step !< Current time step INTEGER, INTENT(in) :: dim_obs_p !< Dimension of PE-local observation vector INTEGER, INTENT(in) :: ncol !< Number of columns in A_p and C_p REAL, INTENT(in) :: obs_p(dim_obs_p) !< PE-local vector of observations REAL, INTENT(in) :: A_p(dim_obs_p, ncol) !< Input matrix REAL, INTENT(out) :: C_p(dim_obs_p, ncol) !< Output matrix ! *** local variables *** INTEGER :: i ! Loop counter INTEGER :: idummy ! Dummy to prevent compiler warning REAL :: rdummy ! Dummy to prevent compiler warning ! Initialize dummy to prevent compiler warning idummy = step rdummy = obs_p(1) ! ************************************* ! *** Compute *** ! *** -1 *** ! *** C = R A *** ! ************************************* DO i=1, n_obstypes CALL PDAFomi_prodRinvA(obs_f_all(i)%ptr, ncol, A_p, C_p) END DO END SUBROUTINE PDAFomi_prodRinvA_cb !------------------------------------------------------------------------------- !> Call-back routine for likelihood !! !! This routine calls the routine PDAFomi_likelihood !! for each observation type !! SUBROUTINE PDAFomi_likelihood_cb(step, dim_obs, obs, resid, lhood) ! Include overall pointer to observation variables use PDAFomi, only: n_obstypes, obs_f_all ! Include PDAFomi function USE PDAFomi, ONLY: PDAFomi_likelihood IMPLICIT NONE ! *** Arguments *** INTEGER, INTENT(in) :: step !< Current time step INTEGER, INTENT(in) :: dim_obs !< PE-local dimension of obs. vector REAL, INTENT(in) :: obs(dim_obs) !< PE-local vector of observations REAL, INTENT(in) :: resid(dim_obs) !< Input vector of residuum REAL, INTENT(out) :: lhood !< Output vector - log likelihood ! *** local variables *** INTEGER :: i ! Loop counter REAL :: rdummy ! Dummy to access obs INTEGER :: idummy ! Dummy to prevent compiler warning ! Initialize dummy to prevent compiler warning idummy = step rdummy = obs(1) ! ************************** ! *** Compute likelihood *** ! ************************** ! Initialize likelihood value before starting computation lhood = 0.0 ! Increment likelihood DO i=1, n_obstypes CALL PDAFomi_likelihood(obs_f_all(i)%ptr, resid, lhood) END DO END SUBROUTINE PDAFomi_likelihood_cb !------------------------------------------------------------------------------- !> Call-back routine for add_obs_error !! !! This routine calls the routine PDAFomi_add_obs_error !! for each observation type !! SUBROUTINE PDAFomi_add_obs_error_cb(step, dim_obs_p, C_p) ! Include overall pointer to observation variables use PDAFomi, only: n_obstypes, obs_f_all ! Include PDAFomi function USE PDAFomi, ONLY: PDAFomi_add_obs_error IMPLICIT NONE ! *** Arguments *** INTEGER, INTENT(in) :: step !< Current time step INTEGER, INTENT(in) :: dim_obs_p !< Dimension of PE-local observation vector REAL, INTENT(inout) :: C_p(dim_obs_p,dim_obs_p) ! Matrix to which R is added ! *** local variables *** INTEGER :: i ! Loop counter INTEGER :: idummy ! Dummy to prevent compiler warning ! Initialize dummy to prevent compiler warning idummy = step ! ************************************* ! *** Add observation error *** ! *** *** ! *** Measurements are uncorrelated *** ! *** here, thus R is diagonal *** ! ************************************* DO i=1, n_obstypes CALL PDAFomi_add_obs_error(obs_f_all(i)%ptr, dim_obs_p, C_p) END DO END SUBROUTINE PDAFomi_add_obs_error_cb !------------------------------------------------------------------------------- !> Call-back routine for init_obscovar !! !! This routine calls the routine PDAFomi_init_obscovar !! for each observation type !! SUBROUTINE PDAFomi_init_obscovar_cb(step, dim_obs, dim_obs_p, covar, m_state_p, & isdiag) ! Include overall pointer to observation variables use PDAFomi, only: n_obstypes, obs_f_all ! Include PDAFomi function USE PDAFomi, ONLY: PDAFomi_init_obscovar IMPLICIT NONE ! *** Arguments *** INTEGER, INTENT(in) :: step !< Current time step INTEGER, INTENT(in) :: dim_obs !< Dimension of observation vector INTEGER, INTENT(in) :: dim_obs_p !< PE-local dimension of obs. vector REAL, INTENT(out) :: covar(dim_obs,dim_obs) !< Observation error covar. matrix REAL, INTENT(in) :: m_state_p(dim_obs_p) !< Observation vector LOGICAL, INTENT(out) :: isdiag !< Whether matrix R is diagonal ! *** local variables *** INTEGER :: i ! Loop counter INTEGER :: idummy ! Dummy to prevent compiler warning REAL :: rdummy ! Dummy to prevent compiler warning ! Initialize dummy to prevent compiler warning idummy = step rdummy = m_state_p(1) ! ************************************* ! *** Initialize covariances *** ! ************************************* DO i=1, n_obstypes CALL PDAFomi_init_obscovar(obs_f_all(i)%ptr, dim_obs_p, covar, isdiag) END DO END SUBROUTINE PDAFomi_init_obscovar_cb !------------------------------------------------------------------------------- !> Call-back routine for init_obserr_f_pdaf !! !! This routine calls the routine PDAFomi_init_obserr_f !! for each observation type !! SUBROUTINE PDAFomi_init_obserr_f_cb(step, dim_obs_f, obs_f, obserr_f) ! Include overall pointer to observation variables use PDAFomi, only: n_obstypes, obs_f_all ! Include PDAFomi function USE PDAFomi, ONLY: PDAFomi_init_obserr_f IMPLICIT NONE ! !ARGUMENTS: INTEGER, INTENT(in) :: step !< Current time step INTEGER, INTENT(in) :: dim_obs_f !< Full dimension of observation vector REAL, INTENT(in) :: obs_f(dim_obs_f) !< Full observation vector REAL, INTENT(out) :: obserr_f(dim_obs_f) !< Full observation error stddev ! *** local variables *** INTEGER :: i ! Loop counter INTEGER :: idummy ! Dummy to prevent compiler warning REAL :: rdummy ! Dummy to prevent compiler warning ! Initialize dummy to prevent compiler warning idummy = step rdummy = obs_f(1) ! ***************************************************************************** ! *** Initialize vector of observation errors for generating synthetic obs. *** ! ***************************************************************************** DO i=1, n_obstypes CALL PDAFomi_init_obserr_f(obs_f_all(i)%ptr, obserr_f) END DO END SUBROUTINE PDAFomi_init_obserr_f_cb !------------------------------------------------------------------------------- !> Call-back routine for omit_by_inno_l !! !! This routine calls the routine PDAFomi_omit_by_inno_l !! for each observation type !! SUBROUTINE PDAFomi_omit_by_inno_l_cb(domain_p, dim_obs_l, resid_l, obs_l) ! Include overall pointer to observation variables USE PDAFomi, ONLY: n_obstypes, obs_f_all, obs_l_all ! Include PDAFomi function USE PDAFomi, ONLY: PDAFomi_omit_by_inno_l ! Include array for statistics on omitted observations USE PDAFomi, ONLY: ostats_omit ! Include filter process rank USE PDAF_mod_filterMPI, ONLY: mype_filter ! Include verbosity information USE PDAF_mod_filter, ONLY: screen #if defined (_OPENMP) ! Include OpenMP function to determine verbosity for OpenMP USE omp_lib, ONLY: omp_get_thread_num #endif IMPLICIT NONE ! !ARGUMENTS: INTEGER, INTENT(in) :: domain_p !< Current local analysis domain INTEGER, INTENT(in) :: dim_obs_l !< PE-local dimension of obs. vector REAL, INTENT(inout) :: resid_l(dim_obs_l) !< Input vector of residuum REAL, INTENT(inout) :: obs_l(dim_obs_l) !< Input vector of local observations ! *** local variables *** INTEGER :: i ! Loop counter INTEGER, SAVE :: domain_save = -1 ! Save previous domain index INTEGER, SAVE :: mythread ! Thread variable for OpenMP INTEGER :: verbose ! Verbosity flag INTEGER :: cnt_omit ! Count of omitted observations ! ********************** ! *** INITIALIZATION *** ! ********************** ! For OpenMP parallelization, determine the thread index #if defined (_OPENMP) mythread = omp_get_thread_num() #else mythread = 0 #endif ! Set verbosity flag (Screen output for first analysis domain) IF (screen > 0) THEN IF ((domain_p < domain_save .OR. domain_save < 0) .AND. mype_filter==0) THEN verbose = 1 ! In case of OpenMP, let only thread 0 write output to the screen IF (mythread>0) verbose = 0 ELSE verbose = 0 END IF ELSE verbose = 0 END IF domain_save = domain_p ! ***************************************************************************** ! *** Initialize vector of observation errors for generating synthetic obs. *** ! ***************************************************************************** cnt_omit = 0 DO i=1, n_obstypes CALL PDAFomi_omit_by_inno_l(obs_l_all(i)%ptr, obs_f_all(i)%ptr, resid_l, obs_l, i, cnt_omit, verbose) END DO !$OMP CRITICAL ! Maximum number of excluded obs. IF (cnt_omit > ostats_omit(6)) ostats_omit(6) = cnt_omit ! Maximum number of use obs. after exclusion IF (dim_obs_l - cnt_omit > ostats_omit(7)) ostats_omit(7) = dim_obs_l - cnt_omit IF (cnt_omit > 0) THEN ostats_omit(1) = ostats_omit(1) + 1 ! Number of domains with excluded obs. ostats_omit(3) = ostats_omit(3) + cnt_omit ! Sum of all excluded observations for all domains ostats_omit(4) = ostats_omit(4) + dim_obs_l - cnt_omit ! Sum of all used observations ostats_omit(5) = ostats_omit(5) + dim_obs_l - cnt_omit ! Sum of all used observations in domain with omissions ELSE ostats_omit(2) = ostats_omit(2) + 1 ! Number of domains without exclusions ostats_omit(4) = ostats_omit(4) + dim_obs_l ! Sum of all used observations END IF !$OMP END CRITICAL END SUBROUTINE PDAFomi_omit_by_inno_l_cb !------------------------------------------------------------------------------- !> Call-back routine for omit_by_inno !! !! This routine calls the routine PDAFomi_omit_by_inno !! for each observation type !! SUBROUTINE PDAFomi_omit_by_inno_cb(dim_obs_f, resid_f, obs_f) ! Include overall pointer to observation variables USE PDAFomi, ONLY: n_obstypes, obs_f_all ! Include PDAFomi functions USE PDAFomi, ONLY: PDAFomi_omit_by_inno, PDAFomi_obsstats ! Include array for statistics on omitted observations USE PDAFomi, ONLY: ostats_omit ! Include verbosity information USE PDAF_mod_filter, ONLY: screen IMPLICIT NONE ! !ARGUMENTS: INTEGER, INTENT(in) :: dim_obs_f !< Full dimension of obs. vector REAL, INTENT(inout) :: resid_f(dim_obs_f) !< Input vector of residuum REAL, INTENT(inout) :: obs_f(dim_obs_f) !< Input vector of full observations ! *** local variables *** INTEGER :: i ! Loop counter INTEGER :: cnt_omit ! Count of omitted observations ! ***************************************************************************** ! *** Initialize vector of observation errors for generating synthetic obs. *** ! ***************************************************************************** cnt_omit = 0 DO i=1, n_obstypes CALL PDAFomi_omit_by_inno(obs_f_all(i)%ptr, resid_f, obs_f, i, cnt_omit) END DO !$OMP CRITICAL ! Maximum number of excluded obs. IF (cnt_omit > ostats_omit(6)) ostats_omit(6) = cnt_omit ! Maximum number of use obs. after exclusion IF (dim_obs_f - cnt_omit > ostats_omit(7)) ostats_omit(7) = dim_obs_f - cnt_omit IF (cnt_omit > 0) THEN ostats_omit(1) = ostats_omit(1) + 1 ! Number of domains with excluded obs. ostats_omit(3) = cnt_omit ! Sum of all excluded observations for all domains ostats_omit(4) = dim_obs_f - cnt_omit ! Sum of all used observations ostats_omit(5) = ostats_omit(5) + dim_obs_f - cnt_omit ! Sum of all used observations in domain with omissions ELSE ostats_omit(2) = ostats_omit(2) + 1 ! Number of domains without exclusions ostats_omit(4) = ostats_omit(4) + dim_obs_f ! Sum of all used observations END IF !$OMP END CRITICAL ! Print observation statistics CALL PDAFomi_obsstats(screen) END SUBROUTINE PDAFomi_omit_by_inno_cb