Changeset 8723 for branches/2017/dev_r8657_UKMO_OBSoper/NEMOGCM/NEMO/OPA_SRC/OBS/diaobs.F90

Timestamp:

2017-11-16T13:18:57+01:00 (6 years ago)

Author:

mattmartin

Message:

Commit changes to obs operator to work with SST bias correction and with the horizontal averaging operator.

File:

• branches/2017/dev_r8657_UKMO_OBSoper/NEMOGCM/NEMO/OPA_SRC/OBS/diaobs.F90 (modified) (10 diffs)

branches/2017/dev_r8657_UKMO_OBSoper/NEMOGCM/NEMO/OPA_SRC/OBS/diaobs.F90

@@ -48 +48 @@
    LOGICAL, PUBLIC :: ln_diaobs   !: Logical switch for the obs operator
    LOGICAL :: ln_sstnight         !: Logical switch for night mean SST obs
-   
+   LOGICAL :: ln_sla_fp_indegs    !: T=> SLA obs footprint size specified in degrees, F=> in metres
+   LOGICAL :: ln_sst_fp_indegs    !: T=> SST obs footprint size specified in degrees, F=> in metres
+   LOGICAL :: ln_sss_fp_indegs    !: T=> SSS obs footprint size specified in degrees, F=> in metres
+   LOGICAL :: ln_sic_fp_indegs    !: T=> sea-ice obs footprint size specified in degrees, F=> in metres
+
+   REAL(wp) :: rn_sla_avglamscl !: E/W diameter of SLA observation footprint (metres)
+   REAL(wp) :: rn_sla_avgphiscl !: N/S diameter of SLA observation footprint (metres)
+   REAL(wp) :: rn_sst_avglamscl !: E/W diameter of SST observation footprint (metres)
+   REAL(wp) :: rn_sst_avgphiscl !: N/S diameter of SST observation footprint (metres)
+   REAL(wp) :: rn_sss_avglamscl !: E/W diameter of SSS observation footprint (metres)
+   REAL(wp) :: rn_sss_avgphiscl !: N/S diameter of SSS observation footprint (metres)
+   REAL(wp) :: rn_sic_avglamscl !: E/W diameter of sea-ice observation footprint (metres)
+   REAL(wp) :: rn_sic_avgphiscl !: N/S diameter of sea-ice observation footprint (metres)
+
    INTEGER :: nn_1dint       !: Vertical interpolation method
-   INTEGER :: nn_2dint       !: Horizontal interpolation method
+   INTEGER :: nn_2dint       !: Default horizontal interpolation method
+   INTEGER :: nn_2dint_sla   !: SLA horizontal interpolation method 
+   INTEGER :: nn_2dint_sst   !: SST horizontal interpolation method 
+   INTEGER :: nn_2dint_sss   !: SSS horizontal interpolation method 
+   INTEGER :: nn_2dint_sic   !: Seaice horizontal interpolation method 
    INTEGER, DIMENSION(imaxavtypes) :: &
       & nn_profdavtypes      !: Profile data types representing a daily average
@@ -61 +78 @@
       & nextrprof, &         !: Number of profile extra variables
       & nextrsurf            !: Number of surface extra variables
-   INTEGER, PUBLIC, ALLOCATABLE, DIMENSION(:) :: sstbias_type !SST bias type    
+   INTEGER, DIMENSION(:), ALLOCATABLE :: &
+      & n2dintsurf           !: Interpolation option for surface variables
+   REAL(wp), DIMENSION(:), ALLOCATABLE :: &
+      & ravglamscl, &        !: E/W diameter of averaging footprint for surface variables
+      & ravgphiscl           !: N/S diameter of averaging footprint for surface variables
+   LOGICAL, DIMENSION(:), ALLOCATABLE :: &
+      & lfpindegs, &         !: T=> surface obs footprint size specified in degrees, F=> in metres
+      & llnightav            !: Logical for calculating night-time averages
+
    TYPE(obs_surf), PUBLIC, POINTER, DIMENSION(:) :: &
       & surfdata, &          !: Initial surface data
@@ -116 +141 @@
          & cn_profbfiles, &      ! T/S profile input filenames
          & cn_sstfbfiles, &      ! Sea surface temperature input filenames
+         & cn_sssfbfiles, &      ! Sea surface salinity input filenames
          & cn_slafbfiles, &      ! Sea level anomaly input filenames
          & cn_sicfbfiles, &      ! Seaice concentration input filenames
          & cn_velfbfiles, &      ! Velocity profile input filenames
-         & cn_sstbias_files      ! SST bias input filenames
+         & cn_sstbiasfiles      ! SST bias input filenames
       CHARACTER(LEN=128) :: &
          & cn_altbiasfile        ! Altimeter bias input filename
@@ -130 +156 @@
       LOGICAL :: ln_sla          ! Logical switch for sea level anomalies 
       LOGICAL :: ln_sst          ! Logical switch for sea surface temperature
+      LOGICAL :: ln_sss          ! Logical switch for sea surface salinity
       LOGICAL :: ln_sic          ! Logical switch for sea ice concentration
       LOGICAL :: ln_vel3d        ! Logical switch for velocity (u,v) obs
       LOGICAL :: ln_nea          ! Logical switch to remove obs near land
       LOGICAL :: ln_altbias      ! Logical switch for altimeter bias
-      LOGICAL :: ln_sstbias     !: Logical switch for bias corection of SST 
+      LOGICAL :: ln_sstbias      ! Logical switch for bias corection of SST 
       LOGICAL :: ln_ignmis       ! Logical switch for ignoring missing files
       LOGICAL :: ln_s_at_t       ! Logical switch to compute model S at T obs
+      LOGICAL :: ln_bound_reject ! Logical to remove obs near boundaries in LAMs.
       LOGICAL :: llvar1          ! Logical for profile variable 1
       LOGICAL :: llvar2          ! Logical for profile variable 1
-      LOGICAL :: llnightav       ! Logical for calculating night-time averages
       LOGICAL, DIMENSION(jpmaxnfiles) :: lmask ! Used for finding number of sstbias files
 
@@ -189 +216 @@
       ! Read namelist parameters
       !-----------------------------------------------------------------------
-      
-      !Initalise all values in namelist arrays
-      ALLOCATE(sstbias_type(jpmaxnfiles))
       ! Some namelist arrays need initialising
       cn_profbfiles(:) = ''
@@ -199 +223 @@
       cn_velfbfiles(:) = ''
       cn_sssfbfiles(:)    = ''
-      cn_sstbias_files(:) = ''
+      cn_sstbiasfiles(:) = ''
       nn_profdavtypes(:) = -1
 
@@ -262 +286 @@
       IF (ln_sstbias) THEN 
          lmask(:) = .FALSE. 
-         WHERE (cn_sstbias_files(:) /= '') lmask(:) = .TRUE. 
+         WHERE (cn_sstbiasfiles(:) /= '') lmask(:) = .TRUE. 
          jnumsstbias = COUNT(lmask) 
          lmask(:) = .FALSE. 
@@ -463 +487 @@
             nvarssurf(jtype) = 1
             nextrsurf(jtype) = 0
-            llnightav = .FALSE.
+            llnightav(jtype) = .FALSE.
             IF ( TRIM(cobstypessurf(jtype)) == 'sla' ) nextrsurf(jtype) = 2
-            IF ( TRIM(cobstypessurf(jtype)) == 'sst' ) llnightav = ln_sstnight
+            IF ( TRIM(cobstypessurf(jtype)) == 'sst' ) llnightav(jtype) = ln_sstnight
 
             !Read in surface obs types
@@ -581 +605 @@
          & zgphi1,    &            ! Model latitudes for prof variable 1
          & zgphi2                  ! Model latitudes for prof variable 2
-      LOGICAL :: llnightav        ! Logical for calculating night-time average
 
       !Allocate local work arrays
@@ -971 +994 @@
    END SUBROUTINE fin_date
    
+    SUBROUTINE obs_settypefiles( ntypes, jpmaxnfiles, jtype, ctypein, &
+       &                         cfilestype, ifiles, cobstypes, cfiles )
+
+    INTEGER, INTENT(IN) :: ntypes      ! Total number of obs types
+    INTEGER, INTENT(IN) :: jpmaxnfiles ! Maximum number of files allowed for each type
+    INTEGER, INTENT(IN) :: jtype       ! Index of the current type of obs
+    INTEGER, DIMENSION(ntypes), INTENT(INOUT) :: &
+       &                   ifiles      ! Out appended number of files for this type
+
+    CHARACTER(len=6), INTENT(IN) :: ctypein 
+    CHARACTER(len=128), DIMENSION(jpmaxnfiles), INTENT(IN) :: &
+       &                   cfilestype  ! In list of files for this obs type
+    CHARACTER(len=6), DIMENSION(ntypes), INTENT(INOUT) :: &
+       &                   cobstypes   ! Out appended list of obs types
+    CHARACTER(len=128), DIMENSION(ntypes, jpmaxnfiles), INTENT(INOUT) :: &
+       &                   cfiles      ! Out appended list of files for all types
+
+    !Local variables
+    INTEGER :: jfile
+
+    cfiles(jtype,:) = cfilestype(:)
+    cobstypes(jtype) = ctypein
+    ifiles(jtype) = 0
+    DO jfile = 1, jpmaxnfiles
+       IF ( trim(cfiles(jtype,jfile)) /= '' ) &
+                 ifiles(jtype) = ifiles(jtype) + 1
+    END DO
+
+    IF ( ifiles(jtype) == 0 ) THEN
+         CALL ctl_stop( 'Logical for observation type '//TRIM(ctypein)//   &
+            &           ' set to true but no files available to read' )
+    ENDIF
+
+    IF(lwp) THEN    
+       WRITE(numout,*) '             '//cobstypes(jtype)//' input observation file names:'
+       DO jfile = 1, ifiles(jtype)
+          WRITE(numout,*) '                '//TRIM(cfiles(jtype,jfile))
+       END DO
+    ENDIF
+
+    END SUBROUTINE obs_settypefiles
+
+    SUBROUTINE obs_setinterpopts( ntypes, jtype, ctypein,             &
+               &                  n2dint_default, n2dint_type,        &
+               &                  ravglamscl_type, ravgphiscl_type,   &
+               &                  lfp_indegs_type, lavnight_type,     &
+               &                  n2dint, ravglamscl, ravgphiscl,     &
+               &                  lfpindegs, lavnight )
+
+    INTEGER, INTENT(IN)  :: ntypes             ! Total number of obs types
+    INTEGER, INTENT(IN)  :: jtype              ! Index of the current type of obs
+    INTEGER, INTENT(IN)  :: n2dint_default     ! Default option for interpolation type
+    INTEGER, INTENT(IN)  :: n2dint_type        ! Option for interpolation type
+    REAL(wp), INTENT(IN) :: &
+       &                    ravglamscl_type, & !E/W diameter of obs footprint for this type
+       &                    ravgphiscl_type    !N/S diameter of obs footprint for this type
+    LOGICAL, INTENT(IN)  :: lfp_indegs_type    !T=> footprint in degrees, F=> in metres
+    LOGICAL, INTENT(IN)  :: lavnight_type      !T=> obs represent night time average
+    CHARACTER(len=6), INTENT(IN) :: ctypein 
+
+    INTEGER, DIMENSION(ntypes), INTENT(INOUT) :: &
+       &                    n2dint 
+    REAL(wp), DIMENSION(ntypes), INTENT(INOUT) :: &
+       &                    ravglamscl, ravgphiscl
+    LOGICAL, DIMENSION(ntypes), INTENT(INOUT) :: &
+       &                    lfpindegs, lavnight
+
+    lavnight(jtype) = lavnight_type
+
+    IF ( (n2dint_type >= 1) .AND. (n2dint_type <= 6) ) THEN
+       n2dint(jtype) = n2dint_type
+    ELSE
+       n2dint(jtype) = n2dint_default
+    ENDIF
+
+    ! For averaging observation footprints set options for size of footprint 
+    IF ( (n2dint(jtype) > 4) .AND. (n2dint(jtype) <= 6) ) THEN
+       IF ( ravglamscl_type > 0._wp ) THEN
+          ravglamscl(jtype) = ravglamscl_type
+       ELSE
+          CALL ctl_stop( 'Incorrect value set for averaging footprint '// &
+                         'scale (ravglamscl) for observation type '//TRIM(ctypein) )      
+       ENDIF
+
+       IF ( ravgphiscl_type > 0._wp ) THEN
+          ravgphiscl(jtype) = ravgphiscl_type
+       ELSE
+          CALL ctl_stop( 'Incorrect value set for averaging footprint '// &
+                         'scale (ravgphiscl) for observation type '//TRIM(ctypein) )      
+       ENDIF
+
+       lfpindegs(jtype) = lfp_indegs_type 
+
+    ENDIF
+
+    ! Write out info 
+    IF(lwp) THEN
+       IF ( n2dint(jtype) <= 4 ) THEN
+          WRITE(numout,*) '             '//TRIM(ctypein)// &
+             &            ' model counterparts will be interpolated horizontally'
+       ELSE IF ( n2dint(jtype) <= 6 ) THEN
+          WRITE(numout,*) '             '//TRIM(ctypein)// &
+             &            ' model counterparts will be averaged horizontally'
+          WRITE(numout,*) '             '//'    with E/W scale: ',ravglamscl(jtype)
+          WRITE(numout,*) '             '//'    with N/S scale: ',ravgphiscl(jtype)
+          IF ( lfpindegs(jtype) ) THEN
+              WRITE(numout,*) '             '//'    (in degrees)'
+          ELSE
+              WRITE(numout,*) '             '//'    (in metres)'
+          ENDIF
+       ENDIF
+    ENDIF
+
+    END SUBROUTINE obs_setinterpopts
+
 END MODULE diaobs