ufo_gnssroonedvarcheck_setom1_mod.f90

Go to the documentation of this file.

!-------------------------------------------------------------------------------
! (C) Crown copyright Met Office. All rights reserved.
!     Refer to COPYRIGHT.txt of this distribution for details.
!-------------------------------------------------------------------------------
! Find a solution to the GPSRO inverse problem.
!-------------------------------------------------------------------------------
 
module ufo_gnssroonedvarcheck_setom1_mod
 
use kinds
use missing_values_mod
use ufo_roobserror_utils_mod, only: rmatrix_type
 
private
public :: ops_gpsro_setom1
 
contains
 
SUBROUTINE ops_gpsro_setom1 (nobs,      &
                             zobs,      &
                             yobs,      &
                             Rmatrix,   &
                             OSigma,    &
                             OM1,       &
                             OM1_error)
 
USE ufo_utils_mod, ONLY: invertmatrix
 
IMPLICIT NONE
 
! Subroutine arguments:
INTEGER, INTENT(IN)             :: nobs
REAL(kind_real), INTENT(IN)     :: zobs(:)
REAL(kind_real), INTENT(IN)     :: yobs(:)
TYPE (rmatrix_type), INTENT(IN) :: rmatrix
REAL(kind_real), INTENT(OUT)    :: osigma(:)
REAL(kind_real), INTENT(OUT)    :: om1(:,:)
LOGICAL, INTENT(OUT)            :: om1_error
 
! Local declarations:
CHARACTER(len=*), PARAMETER     :: routinename = "Ops_GPSRO_setOM1"
INTEGER                         :: n
INTEGER                         :: i
INTEGER                         :: returncode
REAL(kind_real)                 :: frac_err
REAL(kind_real)                 :: omat(nobs,nobs)
REAL(kind_real), ALLOCATABLE    :: gradient(:)
 
IF (rmatrix % satid <= 0) THEN
  ! Rmatrix has not been set up
  om1_error = .true.
  osigma(:) = 1
  om1(:,:) = 1
ELSE
  ALLOCATE (gradient(rmatrix % num_heights - 1))
 
  DO i = 1, rmatrix % num_heights - 1
 
    ! Calculate the gradient of fractional error with height
 
    gradient(i) = (rmatrix % frac_err(i + 1) - rmatrix % frac_err(i)) / &
                  (rmatrix % height(i + 1) - rmatrix % height(i))
 
  END DO
 
  om1_error = .false.
 
  ! Initialise covariance matrix
 
  omat(:,:) = 0.0
 
  ! Calculate the variance values
 
  DO n = 1, nobs
 
    i = 1
 
    DO
 
      IF (zobs(n) < rmatrix % height(i + 1) .OR. &
          i + 1 >= rmatrix % num_heights) THEN
        EXIT
      END IF
 
      i = i + 1
 
    END DO
 
    ! Fractional error
 
    frac_err = rmatrix % frac_err(i) + &
               gradient(i) * (zobs(n) - rmatrix % height(i))
 
    ! Standard deviation
 
    osigma(n) = max(frac_err * yobs(n), rmatrix % min_error)
 
    ! Variance
 
    omat(n,n) = osigma(n) ** 2
 
  END DO
 
  ! Calculate the covariances
 
  DO n = 1,nobs
 
    DO i = n + 1,nobs
 
      omat(n,i) = sqrt(omat(n,n) * omat(i,i)) * &
                  exp(-rmatrix % clen * abs(zobs(i) - zobs(n)))
 
      omat(i,n) = omat(n,i)
 
    END DO
 
  END DO
 
  ! Invert the matrix
 
  CALL invertmatrix (nobs,   &
                     nobs,   &
                     omat,   &
                     returncode)
 
  ! Set the inverse
 
  om1(:,:) = omat(:,:)
 
  IF (returncode /= 0) om1_error = .true.
 
  DEALLOCATE (gradient)
END IF
 
END SUBROUTINE ops_gpsro_setom1
 
end module ufo_gnssroonedvarcheck_setom1_mod