Changeset 15180 for NEMO/branches/UKMO/NEMO_4.0.4_generic_obs

Timestamp:

2021-08-11T13:24:27+02:00 (3 years ago)

Author:

dford

Message:

Further generification, particularly surrounding additional and extra variables.

Location:

NEMO/branches/UKMO/NEMO_4.0.4_generic_obs

Files:

• cfgs/SHARED/namelist_ref (modified) (1 diff)
• src/OCE/OBS/diaobs.F90 (modified) (13 diffs)
• src/OCE/OBS/obs_bias.F90 (moved) (moved from NEMO/branches/UKMO/NEMO_4.0.4_generic_obs/src/OCE/OBS/obs_sstbias.F90) (14 diffs)
• src/OCE/OBS/obs_field.F90 (modified) (16 diffs)
• src/OCE/OBS/obs_oper.F90 (modified) (9 diffs)
• src/OCE/OBS/obs_prep.F90 (modified) (6 diffs)
• src/OCE/OBS/obs_profiles_def.F90 (modified) (26 diffs)
• src/OCE/OBS/obs_read_prof.F90 (modified) (16 diffs)
• src/OCE/OBS/obs_read_surf.F90 (modified) (14 diffs)
• src/OCE/OBS/obs_readmdt.F90 (modified) (3 diffs)
• src/OCE/OBS/obs_surf_def.F90 (modified) (19 diffs)
• src/OCE/OBS/obs_write.F90 (modified) (9 diffs)

NEMO/branches/UKMO/NEMO_4.0.4_generic_obs/cfgs/SHARED/namelist_ref

@@ -1261 +1261 @@
 &namobs_dta    !  observation and model comparison - external data      (see: namobs)
 !-----------------------------------------------------------------------
-   cn_groupname         = ''               ! Name of obs group (for stdout)
-   ln_enabled           = .true.           ! Logical switch for group being processed not ignored
-   ln_prof              = .false.          ! Logical switch for profile data
-   ln_surf              = .false.          ! Logical switch for surface data
-   cn_obsfiles          = ''               ! Observation file names
-   cn_obstypes          = ''               ! Observation types to read from files
-   ln_nea               = .false.          ! Logical switch for rejecting observations near land
-   ln_bound_reject      = .false.          ! Logical switch for rejecting obs near the boundary
-   ln_ignmis            = .true.           ! Logical switch for ignoring missing files
-   nn_2dint             = 0                ! Type of horizontal interpolation method
+   cn_groupname           = ''             ! Name of obs group (output file will be cn_groupname//'fb_????.nc')
+   ln_enabled             = .true.         ! Logical switch for group being processed not ignored
+   ln_prof                = .false.        ! Logical switch for profile data
+   ln_surf                = .false.        ! Logical switch for surface data
+   cn_obsfiles            = ''             ! Observation file names
+   cn_obstypes            = ''             ! Observation types to read from files
+   ln_nea                 = .false.        ! Logical switch for rejecting observations near land
+   ln_bound_reject        = .false.        ! Logical switch for rejecting obs near the boundary
+   ln_ignmis              = .true.         ! Logical switch for ignoring missing files
+   nn_2dint               = 0              ! Type of horizontal interpolation method
                                            ! Relevant if ln_prof = .true.:
-   nn_1dint             = 0                !    Type of vertical interpolation method
-   nn_profdavtypes      = -1               !    Profile data types representing a daily average
+   nn_1dint               = 0              !    Type of vertical interpolation method
+   nn_profdavtypes        = -1             !    Profile data types representing a daily average
+   ln_all_at_all          = .false.        !    Logical switch for computing all model variables at all obs points
                                            ! Relevant if ln_surf = .true.:
-   ln_fp_indegs         = .true.           !    Logical: T=> averaging footprint is in degrees, F=> in metres
-   rn_avglamscl         = 0.               !    E/W diameter of observation footprint (metres/degrees)
-   rn_avgphiscl         = 0.               !    N/S diameter of observation footprint (metres/degrees)
-   ln_night             = .false.          !    Logical switch for calculating night-time average for obs
-                                           ! Relevant if 'SST' in cn_obstypes:
-   ln_sstbias           = .false.          !    Logical switch for SST bias correction
-   cn_sstbiasfiles      = ''               !    SST bias input file names
+   ln_fp_indegs           = .true.         !    Logical: T=> averaging footprint is in degrees, F=> in metres
+   rn_avglamscl           = 0.             !    E/W diameter of observation footprint (metres/degrees)
+   rn_avgphiscl           = 0.             !    N/S diameter of observation footprint (metres/degrees)
+   ln_night               = .false.        !    Logical switch for calculating night-time average for obs
+   ln_obsbias             = .false.        !    Logical switch for bias correction
+   cn_obsbiasfiles        = ''             !    Bias input file names
+   cn_type_to_biascorrect = ''             !    Observation type to bias correct
+   cn_obsbiasfile_varname = ''             !    Bias variable name in input file
                                            ! Relevant if 'SLA' in cn_obstypes:
-   ln_altbias           = .false.          !    Logical switch for altimeter bias correction
-   cn_altbiasfile       = ''               !    Altimeter bias input file name
-   nn_msshc             = 0                !    MSSH correction scheme
-   rn_mdtcorr           = 1.61             !    MDT correction
-   rn_mdtcutoff         = 65.0             !    MDT cutoff for computed correction
-   ln_time_mean_sla_bkg = .false.          !    Logical switch for applying time mean of SLA background to remove tidal signal
-                                           ! Relevant if 'POTM' and/or 'PSAL' in cn_obstypes:
-   ln_s_at_t            = .false.          !    Logical switch for computing model S at T obs if not there
-   ln_output_clim       = .false.          !    Logical switch for writing climatological values to fdbk files
+   ln_altbias             = .false.        !    Logical switch for altimeter bias correction
+   cn_altbiasfile         = ''             !    Altimeter bias input file name
+   nn_msshc               = 0              !    MSSH correction scheme
+   rn_mdtcorr             = 1.61           !    MDT correction
+   rn_mdtcutoff           = 65.0           !    MDT cutoff for computed correction
+!!! NOT YET IMPLEMENTED:
+!!!   OUTPUT CLIMATOLOGY
+!!!   TIME MEAN BACKGROUND
 /
 !-----------------------------------------------------------------------

NEMO/branches/UKMO/NEMO_4.0.4_generic_obs/src/OCE/OBS/diaobs.F90

@@ -28 +28 @@
    USE par_kind       ! Precision variables
    USE in_out_manager ! I/O manager
+   USE timing         ! Timing
    USE par_oce        ! ocean parameter
    USE dom_oce        ! Ocean space and time domain variables
@@ -34 +35 @@
    USE obs_read_prof  ! Reading and allocation of profile obs
    USE obs_read_surf  ! Reading and allocation of surface obs
-   USE obs_sstbias    ! Bias correction routine for SST 
+   USE obs_bias       ! Bias correction routine
    USE obs_readmdt    ! Reading and allocation of MDT for SLA.
    USE obs_prep       ! Preparation of obs. (grid search etc).
@@ -197 +198 @@
                llvar(:) = .TRUE.
                !
-               CALL obs_rea_prof( sobsgroups(jgroup)%sprofdata, &
-                  &               sobsgroups(jgroup)%nobsfiles, &
-                  &               sobsgroups(jgroup)%cobsfiles, &
-                  &               sobsgroups(jgroup)%nobstypes, &
-                  &               sobsgroups(jgroup)%nextvars,  &
-                  &               nitend-nit000+2,              &
-                  &               rn_dobsini,                   &
-                  &               rn_dobsend,                   &
-                  &               llvar,                        &
-                  &               sobsgroups(jgroup)%lignmis,   &
-                  &               sobsgroups(jgroup)%ls_at_t,   &
-                  &               .FALSE.,                      &
-                  &               sobsgroups(jgroup)%cobstypes, &
+               CALL obs_rea_prof( sobsgroups(jgroup)%sprofdata,   &
+                  &               sobsgroups(jgroup)%nobsfiles,   &
+                  &               sobsgroups(jgroup)%cobsfiles,   &
+                  &               sobsgroups(jgroup)%nobstypes,   &
+                  &               sobsgroups(jgroup)%naddvars,    &
+                  &               sobsgroups(jgroup)%nextvars,    &
+                  &               nitend-nit000+2,                &
+                  &               rn_dobsini,                     &
+                  &               rn_dobsend,                     &
+                  &               llvar,                          &
+                  &               sobsgroups(jgroup)%lignmis,     &
+                  &               sobsgroups(jgroup)%lall_at_all, &
+                  &               .FALSE.,                        &
+                  &               sobsgroups(jgroup)%cobstypes,   &
                   &               kdailyavtypes = sobsgroups(jgroup)%nprofdavtypes )
                !
@@ -237 +239 @@
                   &               sobsgroups(jgroup)%cobsfiles, &
                   &               sobsgroups(jgroup)%nobstypes, &
+                  &               sobsgroups(jgroup)%naddvars,  &
                   &               sobsgroups(jgroup)%nextvars,  &
                   &               nitend-nit000+2,              &
@@ -246 +249 @@
                   &               sobsgroups(jgroup)%cobstypes )
                   !
-               CALL obs_pre_surf( sobsgroups(jgroup)%ssurfdata,   &
-                  &               sobsgroups(jgroup)%ssurfdataqc, &
-                  &               sobsgroups(jgroup)%lnea,        &
+
+               CALL obs_pre_surf( sobsgroups(jgroup)%ssurfdata,      &
+                  &               sobsgroups(jgroup)%ssurfdataqc,    &
+                  &               jpi, jpj,                          &
+                  &               sobsgroups(jgroup)%rmask(:,:,1,:), &
+                  &               sobsgroups(jgroup)%rglam,          &
+                  &               sobsgroups(jgroup)%rgphi,          &
+                  &               sobsgroups(jgroup)%lnea,           &
                   &               sobsgroups(jgroup)%lbound_reject )
                !
@@ -261 +269 @@
                ENDIF
                !
-               IF( sobsgroups(jgroup)%lsst .AND. sobsgroups(jgroup)%lsstbias ) THEN
-                  CALL obs_app_sstbias( sobsgroups(jgroup)%ssurfdataqc,   &
-                     &                  sobsgroups(jgroup)%n2dint,        & 
-                     &                  sobsgroups(jgroup)%nsstbiasfiles, &
-                     &                  sobsgroups(jgroup)%csstbiasfiles ) 
+               IF( sobsgroups(jgroup)%lobsbias ) THEN
+                  CALL obs_app_bias( sobsgroups(jgroup)%ssurfdataqc,   &
+                     &               sobsgroups(jgroup)%nbiasvar,      & 
+                     &               sobsgroups(jgroup)%n2dint,        & 
+                     &               sobsgroups(jgroup)%nobsbiasfiles, &
+                     &               sobsgroups(jgroup)%cobsbiasfiles, &
+                     &               sobsgroups(jgroup)%cbiasvarname ) 
                ENDIF
                !
@@ -315 +325 @@
 
       !-----------------------------------------------------------------------
+
+      IF( ln_timing )   CALL timing_start('dia_obs')
 
       IF(lwp) THEN
@@ -405 +417 @@
                      &               kstp, jpi, jpj,                       &
                      &               nit000, idaystp,                      &
-                     &               zsurfvar,                             &
+                     &               jvar, zsurfvar,                       &
                      &               sobsgroups(jgroup)%rmask(:,:,1,jvar), &
                      &               sobsgroups(jgroup)%n2dint,            &
@@ -421 +433 @@
          ENDIF
       END DO
+
+      IF( ln_timing )   CALL timing_stop('dia_obs')
 
    END SUBROUTINE dia_obs
@@ -449 +463 @@
       !! * Local declarations
       INTEGER :: jgroup                   ! Data set loop variable
-      INTEGER :: jo, jvar, jk
+      INTEGER :: jo, jvar, jk, jadd, jext
       REAL(wp), DIMENSION(:), ALLOCATABLE :: &
          & zu, &
          & zv
+      TYPE(obswriinfo) :: sladd, slext
+
+      IF( ln_timing )   CALL timing_start('dia_obs_wri')
 
       !-----------------------------------------------------------------------
@@ -496 +513 @@
                   &                      sobsgroups(jgroup)%sprofdata, .TRUE., numout )
 
-               CALL obs_wri_prof( sobsgroups(jgroup)%sprofdata )
+               sladd%inum = sobsgroups(jgroup)%sprofdata%nadd
+               IF ( sladd%inum > 0 ) THEN
+                  ALLOCATE( sladd%ipoint(sladd%inum),                                   &
+                     &      sladd%cdname(sladd%inum),                                   &
+                     &      sladd%cdlong(sladd%inum,sobsgroups(jgroup)%sprofdata%nvar), &
+                     &      sladd%cdunit(sladd%inum,sobsgroups(jgroup)%sprofdata%nvar) )
+                  DO jadd = 1, sladd%inum
+                     sladd%ipoint(jadd) = jadd
+                     sladd%cdname(jadd) = sobsgroups(jgroup)%sprofdata%caddvars(jadd)
+                     DO jvar = 1, sobsgroups(jgroup)%sprofdata%nvar
+                        sladd%cdlong(jadd,jvar) = sobsgroups(jgroup)%sprofdata%caddlong(jadd,jvar)
+                        sladd%cdunit(jadd,jvar) = sobsgroups(jgroup)%sprofdata%caddunit(jadd,jvar)
+                     END DO
+                  END DO
+               ENDIF
+               slext%inum = sobsgroups(jgroup)%sprofdata%next
+               IF ( slext%inum > 0 ) THEN
+                  ALLOCATE( slext%ipoint(slext%inum),   &
+                     &      slext%cdname(slext%inum),   &
+                     &      slext%cdlong(slext%inum,1), &
+                     &      slext%cdunit(slext%inum,1) )
+                  DO jext = 1, slext%inum
+                     slext%ipoint(jext)   = jext
+                     slext%cdname(jext)   = sobsgroups(jgroup)%sprofdata%cextvars(jext)
+                     slext%cdlong(jext,1) = sobsgroups(jgroup)%sprofdata%cextlong(jext)
+                     slext%cdunit(jext,1) = sobsgroups(jgroup)%sprofdata%cextunit(jext)
+                  END DO
+               ENDIF
+
+               CALL obs_wri_prof( sobsgroups(jgroup)%sprofdata, sobsgroups(jgroup)%cgroupname, sladd, slext )
+
+               IF ( sladd%inum > 0 ) THEN
+                  DEALLOCATE( sladd%ipoint, sladd%cdname, sladd%cdlong, sladd%cdunit )
+               ENDIF
+               IF ( slext%inum > 0 ) THEN
+                  DEALLOCATE( slext%ipoint, slext%cdname, slext%cdlong, slext%cdunit )
+               ENDIF
 
             ELSEIF (sobsgroups(jgroup)%lsurf) THEN
@@ -503 +556 @@
                   &                      sobsgroups(jgroup)%ssurfdata, .TRUE., numout )
 
-               CALL obs_wri_surf( sobsgroups(jgroup)%ssurfdata )
+               sladd%inum = sobsgroups(jgroup)%ssurfdata%nadd
+               IF ( sladd%inum > 0 ) THEN
+                  ALLOCATE( sladd%ipoint(sladd%inum),                                   &
+                     &      sladd%cdname(sladd%inum),                                   &
+                     &      sladd%cdlong(sladd%inum,sobsgroups(jgroup)%ssurfdata%nvar), &
+                     &      sladd%cdunit(sladd%inum,sobsgroups(jgroup)%ssurfdata%nvar) )
+                  DO jadd = 1, sladd%inum
+                     sladd%ipoint(jadd) = jadd
+                     sladd%cdname(jadd) = sobsgroups(jgroup)%ssurfdata%caddvars(jadd)
+                     DO jvar = 1, sobsgroups(jgroup)%ssurfdata%nvar
+                        sladd%cdlong(jadd,jvar) = sobsgroups(jgroup)%ssurfdata%caddlong(jadd,jvar)
+                        sladd%cdunit(jadd,jvar) = sobsgroups(jgroup)%ssurfdata%caddunit(jadd,jvar)
+                     END DO
+                  END DO
+               ENDIF
+               slext%inum = sobsgroups(jgroup)%ssurfdata%nextra
+               IF ( slext%inum > 0 ) THEN
+                  ALLOCATE( slext%ipoint(slext%inum),   &
+                     &      slext%cdname(slext%inum),   &
+                     &      slext%cdlong(slext%inum,1), &
+                     &      slext%cdunit(slext%inum,1) )
+                  DO jext = 1, slext%inum
+                     slext%ipoint(jext)   = jext
+                     slext%cdname(jext)   = sobsgroups(jgroup)%ssurfdata%cextvars(jext)
+                     slext%cdlong(jext,1) = sobsgroups(jgroup)%ssurfdata%cextlong(jext)
+                     slext%cdunit(jext,1) = sobsgroups(jgroup)%ssurfdata%cextunit(jext)
+                  END DO
+               ENDIF
+
+               CALL obs_wri_surf( sobsgroups(jgroup)%ssurfdata, sobsgroups(jgroup)%cgroupname, sladd, slext )
+
+               IF ( sladd%inum > 0 ) THEN
+                  DEALLOCATE( sladd%ipoint, sladd%cdname, sladd%cdlong, sladd%cdunit )
+               ENDIF
+               IF ( slext%inum > 0 ) THEN
+                  DEALLOCATE( slext%ipoint, slext%cdname, slext%cdlong, slext%cdunit )
+               ENDIF
 
             ENDIF
@@ -510 +599 @@
 
       END DO
+
+      IF( ln_timing )   CALL timing_stop('dia_obs_wri')
 
    END SUBROUTINE dia_obs_wri

NEMO/branches/UKMO/NEMO_4.0.4_generic_obs/src/OCE/OBS/obs_bias.F90

@@ -1 @@
-MODULE obs_sstbias
+MODULE obs_bias
    !!======================================================================
-   !!                       ***  MODULE obs_sstbias  ***
-   !! Observation diagnostics: Read the bias for SST data
+   !!                       ***  MODULE obs_bias  ***
+   !! Observation diagnostics: Read the bias for observation data
    !!======================================================================
    !!----------------------------------------------------------------------
-   !!   obs_app_sstbias : Driver for reading and applying the SST bias
+   !!   obs_app_bias : Driver for reading and applying the bias
    !!----------------------------------------------------------------------
    !! * Modules used   
    USE par_kind, ONLY : &       ! Precision variables
-      & wp, &
-      & dp, &
-      & sp
+      & wp
    USE par_oce, ONLY : &        ! Domain parameters
       & jpi, &
-      & jpj, &
-      & jpim1
+      & jpj
    USE in_out_manager, ONLY : & ! I/O manager
       & lwp,    &
@@ -22 +19 @@
    USE dom_oce, ONLY : &        ! Domain variables
       & tmask, &
-      & tmask_i, &
-      & e1t,   &
-      & e2t,   &
       & gphit, &
       & glamt
-   USE oce, ONLY : &           ! Model variables
-      & sshn
    USE obs_inter_h2d
    USE obs_utils               ! Various observation tools
@@ -35 +27 @@
    !! * Routine accessibility
    PRIVATE
-   PUBLIC obs_app_sstbias     ! Read the altimeter bias
+   PUBLIC obs_app_bias     ! Read the observation bias
 CONTAINS
-   SUBROUTINE obs_app_sstbias( sstdata, k2dint, knumtypes, &
-                               cl_bias_files )
+   SUBROUTINE obs_app_bias( obsdata, kvar, k2dint, knumtypes, &
+                            cl_bias_files, cd_biasname )
       !!---------------------------------------------------------------------
       !!
-      !!                   *** ROUTINE obs_app_sstbias ***
-      !!
-      !! ** Purpose : Read SST bias data from files and apply correction to
-      !!             observations
+      !!                   *** ROUTINE obs_app_bias ***
+      !!
+      !! ** Purpose : Read bias data from files and apply correction to
+      !!              observations
       !!
       !! ** Method  :
@@ -54 +46 @@
       !! History : 
       !!      ! :  2014-08 (J. While) Bias correction code for SST obs,
-      !!      !                      based on obs_rea_altbias
+      !!      !                       based on obs_rea_altbias
+      !!      ! :  2021-07 (D. Ford)  Renamed obs_app_bias and made generic
       !!----------------------------------------------------------------------
       !! * Modules used
@@ -61 +54 @@
       !! * Arguments
 
-      TYPE(obs_surf), INTENT(INOUT) :: sstdata       ! SST data
+      TYPE(obs_surf), INTENT(INOUT) :: obsdata       ! Observation data
+      INTEGER, INTENT(IN) :: kvar    ! Index of obs type being bias corrected
       INTEGER, INTENT(IN) :: k2dint
       INTEGER, INTENT(IN) :: knumtypes !number of bias types to read in
       CHARACTER(LEN=128), DIMENSION(knumtypes), INTENT(IN) :: &
                           cl_bias_files !List of files to read
+      CHARACTER(LEN=128), INTENT(IN) :: cd_biasname  !Variable name in file
       !! * Local declarations
       INTEGER :: jobs         ! Obs loop variable
-      INTEGER :: jpisstbias   ! Number of grid point in latitude for the bias
-      INTEGER :: jpjsstbias   ! Number of grid point in longitude for the bias
       INTEGER :: iico         ! Grid point indices
       INTEGER :: ijco
       INTEGER :: jt
-      INTEGER :: i_nx_id      ! Index to read the NetCDF file
-      INTEGER :: i_ny_id      !
-      INTEGER :: i_file_id    !
-      INTEGER :: i_var_id
       INTEGER, DIMENSION(knumtypes) :: &
          & ibiastypes             ! Array of the bias types in each file
       REAL(wp), DIMENSION(jpi,jpj,knumtypes) :: & 
-         & z_sstbias              ! Array to store the SST bias values
+         & z_obsbias              ! Array to store the bias values
       REAL(wp), DIMENSION(jpi,jpj) :: & 
-         & z_sstbias_2d           ! Array to store the SST bias values   
+         & z_obsbias_2d           ! Array to store the bias values   
       REAL(wp), DIMENSION(1) :: &
          & zext, &
@@ -105 +94 @@
          & igrdi_tmp, &
          & igrdj_tmp   
-      INTEGER :: numsstbias
+      INTEGER :: numobsbias
       INTEGER(KIND=NF90_INT) :: ifile_source
      
@@ -113 +102 @@
       INTEGER :: inumtype
       IF(lwp)WRITE(numout,*) 
-      IF(lwp)WRITE(numout,*) 'obs_rea_sstbias : '
+      IF(lwp)WRITE(numout,*) 'obs_app_bias : '
       IF(lwp)WRITE(numout,*) '----------------- '
-      IF(lwp)WRITE(numout,*) 'Read SST bias '
+      IF(lwp)WRITE(numout,*) 'Read observation bias for ', TRIM(obsdata%cvars(kvar))
       ! Open and read the files
-      z_sstbias(:,:,:)=0.0_wp
+      z_obsbias(:,:,:)=0.0_wp
       DO jtype = 1, knumtypes
      
-         numsstbias=0
+         numobsbias=0
          IF(lwp)WRITE(numout,*) 'Opening ',cl_bias_files(jtype)
-         CALL iom_open( cl_bias_files(jtype), numsstbias, ldstop=.FALSE. )       
-         IF (numsstbias > 0) THEN
+         CALL iom_open( cl_bias_files(jtype), numobsbias, ldstop=.FALSE. )       
+         IF (numobsbias > 0) THEN
      
             !Read the bias type from the file
@@ -130 +119 @@
             !routines directly - should be upgraded in the future
             iret=NF90_OPEN(TRIM(cl_bias_files(jtype)), NF90_NOWRITE, incfile)
-            iret=NF90_GET_ATT( incfile, NF90_GLOBAL, "SST_source", &
+            iret=NF90_GET_ATT( incfile, NF90_GLOBAL, TRIM(obsdata%cvars(kvar))//"_source", &
                               ifile_source )
             ibiastypes(jtype) = ifile_source                 
@@ -136 +125 @@
            
             IF ( iret /= 0  ) CALL ctl_stop( &
-               'obs_rea_sstbias : Cannot read bias type from file '// &
+               'obs_app_bias : Cannot read bias type from file '// &
                cl_bias_files(jtype) )
-            ! Get the SST bias data
-            CALL iom_get( numsstbias, jpdom_data, 'tn', z_sstbias_2d(:,:), 1 )
-            z_sstbias(:,:,jtype) = z_sstbias_2d(:,:)       
+            ! Get the bias data
+            CALL iom_get( numobsbias, jpdom_data, TRIM(cd_biasname), z_obsbias_2d(:,:), 1 )
+            z_obsbias(:,:,jtype) = z_obsbias_2d(:,:)       
             ! Close the file
-            CALL iom_close(numsstbias)       
+            CALL iom_close(numobsbias)       
          ELSE
-            CALL ctl_stop('obs_read_sstbias: File '// & 
+            CALL ctl_stop('obs_app_bias: File '// & 
                            TRIM( cl_bias_files(jtype) )//' Not found')
          ENDIF
@@ -151 +140 @@
       ! Interpolate the bias already on the model grid at the observation point
       ALLOCATE( &
-         & igrdi(2,2,sstdata%nsurf), &
-         & igrdj(2,2,sstdata%nsurf), &
-         & zglam(2,2,sstdata%nsurf), &
-         & zgphi(2,2,sstdata%nsurf), &
-         & zmask(2,2,sstdata%nsurf)  )
+         & igrdi(2,2,obsdata%nsurf), &
+         & igrdj(2,2,obsdata%nsurf), &
+         & zglam(2,2,obsdata%nsurf), &
+         & zgphi(2,2,obsdata%nsurf), &
+         & zmask(2,2,obsdata%nsurf)  )
        
-      DO jobs = 1, sstdata%nsurf 
-         igrdi(1,1,jobs) = sstdata%mi(jobs)-1
-         igrdj(1,1,jobs) = sstdata%mj(jobs)-1
-         igrdi(1,2,jobs) = sstdata%mi(jobs)-1
-         igrdj(1,2,jobs) = sstdata%mj(jobs)
-         igrdi(2,1,jobs) = sstdata%mi(jobs)
-         igrdj(2,1,jobs) = sstdata%mj(jobs)-1
-         igrdi(2,2,jobs) = sstdata%mi(jobs)
-         igrdj(2,2,jobs) = sstdata%mj(jobs)
+      DO jobs = 1, obsdata%nsurf 
+         igrdi(1,1,jobs) = obsdata%mi(jobs)-1
+         igrdj(1,1,jobs) = obsdata%mj(jobs)-1
+         igrdi(1,2,jobs) = obsdata%mi(jobs)-1
+         igrdj(1,2,jobs) = obsdata%mj(jobs)
+         igrdi(2,1,jobs) = obsdata%mi(jobs)
+         igrdj(2,1,jobs) = obsdata%mj(jobs)-1
+         igrdi(2,2,jobs) = obsdata%mi(jobs)
+         igrdj(2,2,jobs) = obsdata%mj(jobs)
       END DO
-      CALL obs_int_comm_2d( 2, 2, sstdata%nsurf, jpi, jpj, &
+      CALL obs_int_comm_2d( 2, 2, obsdata%nsurf, jpi, jpj, &
          &                  igrdi, igrdj, glamt, zglam )
-      CALL obs_int_comm_2d( 2, 2, sstdata%nsurf, jpi, jpj, &
+      CALL obs_int_comm_2d( 2, 2, obsdata%nsurf, jpi, jpj, &
          &                  igrdi, igrdj, gphit, zgphi )
-      CALL obs_int_comm_2d( 2, 2, sstdata%nsurf, jpi, jpj, &
+      CALL obs_int_comm_2d( 2, 2, obsdata%nsurf, jpi, jpj, &
          &                  igrdi, igrdj, tmask(:,:,1), zmask )
       DO jtype = 1, knumtypes
          
          !Find the number observations of type and allocate tempory arrays
-         inumtype = COUNT( sstdata%ntyp(:) == ibiastypes(jtype) )
+         inumtype = COUNT( obsdata%ntyp(:) == ibiastypes(jtype) )
          ALLOCATE( &
             & igrdi_tmp(2,2,inumtype), &
@@ -185 +174 @@
             & zbias( 2,2,inumtype ) )
          jt=1
-         DO jobs = 1, sstdata%nsurf 
-            IF ( sstdata%ntyp(jobs) == ibiastypes(jtype) ) THEN
+         DO jobs = 1, obsdata%nsurf 
+            IF ( obsdata%ntyp(jobs) == ibiastypes(jtype) ) THEN
                igrdi_tmp(:,:,jt) = igrdi(:,:,jobs) 
                igrdj_tmp(:,:,jt) = igrdj(:,:,jobs)
@@ -198 +187 @@
          CALL obs_int_comm_2d( 2, 2, inumtype, jpi, jpj, &
                &           igrdi_tmp(:,:,:), igrdj_tmp(:,:,:), &
-               &           z_sstbias(:,:,jtype), zbias(:,:,:) )
+               &           z_obsbias(:,:,jtype), zbias(:,:,:) )
          jt=1
-         DO jobs = 1, sstdata%nsurf
-            IF ( sstdata%ntyp(jobs) == ibiastypes(jtype) ) THEN
-               zlam = sstdata%rlam(jobs)
-               zphi = sstdata%rphi(jobs)
-               iico = sstdata%mi(jobs)
-               ijco = sstdata%mj(jobs)         
+         DO jobs = 1, obsdata%nsurf
+            IF ( obsdata%ntyp(jobs) == ibiastypes(jtype) ) THEN
+               zlam = obsdata%rlam(jobs)
+               zphi = obsdata%rphi(jobs)
+               iico = obsdata%mi(jobs)
+               ijco = obsdata%mj(jobs)         
                CALL obs_int_h2d_init( 1, 1, k2dint, zlam, zphi,         &
                   &                   zglam_tmp(:,:,jt), &
@@ -211 +200 @@
                   &                   zmask_tmp(:,:,jt), zweig, zobsmask )
                CALL obs_int_h2d( 1, 1, zweig, zbias(:,:,jt),  zext )
-               ! adjust sst with bias field
-               sstdata%robs(jobs,1) = sstdata%robs(jobs,1) - zext(1)
+               ! adjust observations with bias field
+               obsdata%robs(jobs,kvar) = obsdata%robs(jobs,kvar) - zext(1)
                jt=jt+1
             ENDIF
@@ -235 +224 @@
       IF(lwp) THEN
          WRITE(numout,*) " "
-         WRITE(numout,*) "SST bias correction applied successfully"
+         WRITE(numout,*) "Bias correction applied successfully"
          WRITE(numout,*) "Obs types: ",ibiastypes(:), &
                               " Have all been bias corrected\n"
       ENDIF
-   END SUBROUTINE obs_app_sstbias
+   END SUBROUTINE obs_app_bias
  
-END MODULE obs_sstbias
+END MODULE obs_bias

NEMO/branches/UKMO/NEMO_4.0.4_generic_obs/src/OCE/OBS/obs_field.F90

@@ -37 +37 @@
    
    ! Expected names for observation types with special behaviours (not needed for all observation types)
-   CHARACTER(LEN=8) :: cobsname_sst    = 'SST'  ! Expected variable name for SST
-   CHARACTER(LEN=8) :: cobsname_temp3d = 'POTM' ! Expected variable name for 3D temperature
-   CHARACTER(LEN=8) :: cobsname_sal3d  = 'PSAL' ! Expected variable name for 3D salinity
    CHARACTER(LEN=8) :: cobsname_uvel3d = 'UVEL' ! Expected variable name for 3D zonal currents
    CHARACTER(LEN=8) :: cobsname_vvel3d = 'VVEL' ! Expected variable name for 3D meridional currents
@@ -47 +44 @@
    TYPE obs_group
       !
-      CHARACTER(LEN=128)                            :: cgroupname    !: Name of obs group (for stdout)
+      CHARACTER(LEN=25)                             :: cgroupname    !: Name of obs group (for stdout)
       CHARACTER(LEN=8),   DIMENSION(:), ALLOCATABLE :: cobstypes     !: Observation types to read from files
       CHARACTER(LEN=128), DIMENSION(:), ALLOCATABLE :: cobsfiles     !: Observation file names
-      CHARACTER(LEN=128), DIMENSION(:), ALLOCATABLE :: csstbiasfiles !: SST bias input file names
+      CHARACTER(LEN=128), DIMENSION(:), ALLOCATABLE :: cobsbiasfiles !: Bias input file names
+      CHARACTER(LEN=128)                            :: cbiasvarname  !: Bias variable name in input file
       CHARACTER(LEN=128)                            :: caltbiasfile  !: Altimeter bias input file name
       !
@@ -57 +55 @@
       INTEGER  :: nobstypes          !: Number of observation types
       INTEGER  :: nobsfiles          !: Number of observation files
-      INTEGER  :: nextvars           !: Number of extra variables to get
-      INTEGER  :: nsstbiasfiles      !: Number of SST bias files
+      INTEGER  :: nobsbiasfiles      !: Number of bias files
+      INTEGER  :: nbiasvar           !: Index of observation type to be bias corrected
       INTEGER  :: navtypes           !: Number of profile data types representing a daily average
+      INTEGER  :: nextvars           !: Number of extra variables in addition to any in input files
+      INTEGER  :: naddvars           !: Number of additional variables in addition to any in input files
       INTEGER  :: n1dint             !: Type of vertical interpolation method
       INTEGER  :: n2dint             !: Type of horizontal interpolation method
@@ -67 +67 @@
       LOGICAL  :: lsurf              !: Logical switch for surface data
       LOGICAL  :: lprof              !: Logical switch for profile data
-      LOGICAL  :: lsst               !: Logical switch for SST data
-      LOGICAL  :: ltemp3d            !: Logical switch for 3D temperature data
-      LOGICAL  :: lsal3d             !: Logical switch for 3D salinity data
       LOGICAL  :: lvel3d             !: Logical switch for 3D velocity data
       LOGICAL  :: lsla               !: Logical switch for SLA data
       LOGICAL  :: laltbias           !: Logical switch for altimeter bias correction
-      LOGICAL  :: lsstbias           !: Logical switch for SST bias correction
+      LOGICAL  :: lobsbias           !: Logical switch for bias correction
       LOGICAL  :: lnea               !: Logical switch for rejecting observations near land
       LOGICAL  :: lbound_reject      !: Logical switch for rejecting obs near the boundary
       LOGICAL  :: lignmis            !: Logical switch for ignoring missing files
-      LOGICAL  :: ls_at_t            !: Logical switch for computing model S at T obs if not there
+      LOGICAL  :: lall_at_all        !: Logical switch for computing all model variables at all obs points
       LOGICAL  :: lnight             !: Logical switch for calculating night-time average
-      LOGICAL  :: loutput_clim       !: Logical switch for writing climatological values to fdbk files
-      LOGICAL  :: ltime_mean_sla_bkg !: Logical switch for applying time mean of SLA background to remove tidal signal
       LOGICAL  :: lfp_indegs         !: Logical: T=> averaging footprint is in degrees, F=> in metres
       !
@@ -118 +113 @@
       ALLOCATE( sdobsgroup%cobstypes    (sdobsgroup%nobstypes            ), &
          &      sdobsgroup%cobsfiles    (sdobsgroup%nobsfiles            ), &
-         &      sdobsgroup%csstbiasfiles(sdobsgroup%nsstbiasfiles        ), &
+         &      sdobsgroup%cobsbiasfiles(sdobsgroup%nobsbiasfiles        ), &
          &      sdobsgroup%nprofdavtypes(sdobsgroup%navtypes             ), &
          &      sdobsgroup%rglam        (jpi,jpj,    sdobsgroup%nobstypes), &
@@ -144 +139 @@
       DEALLOCATE( sdobsgroup%cobstypes,     &
          &        sdobsgroup%cobsfiles,     &
-         &        sdobsgroup%csstbiasfiles, &
+         &        sdobsgroup%cobsbiasfiles, &
          &        sdobsgroup%nprofdavtypes, &
          &        sdobsgroup%rglam,         &
@@ -174 +169 @@
       CHARACTER(LEN=8),   DIMENSION(jpmaxntypes) :: cn_obstypes
       CHARACTER(LEN=128), DIMENSION(jpmaxnfiles) :: cn_obsfiles
-      CHARACTER(LEN=128), DIMENSION(jpmaxnfiles) :: cn_sstbiasfiles
+      CHARACTER(LEN=128), DIMENSION(jpmaxnfiles) :: cn_obsbiasfiles
+      CHARACTER(LEN=128)                         :: cn_type_to_biascorrect
+      CHARACTER(LEN=128)                         :: cn_obsbiasfile_varname
       CHARACTER(LEN=128)                         :: cn_altbiasfile
       INTEGER,            DIMENSION(imaxavtypes) :: nn_profdavtypes
@@ -184 +181 @@
       LOGICAL                                    :: ln_prof
       LOGICAL                                    :: ln_altbias
-      LOGICAL                                    :: ln_sstbias
+      LOGICAL                                    :: ln_obsbias
       LOGICAL                                    :: ln_nea
       LOGICAL                                    :: ln_bound_reject
       LOGICAL                                    :: ln_ignmis
-      LOGICAL                                    :: ln_s_at_t
+      LOGICAL                                    :: ln_all_at_all
       LOGICAL                                    :: ln_night
-      LOGICAL                                    :: ln_output_clim
-      LOGICAL                                    :: ln_time_mean_sla_bkg
       LOGICAL                                    :: ln_fp_indegs
       REAL(wp)                                   :: rn_avglamscl
@@ -201 +196 @@
          &                cn_obsfiles, cn_obstypes, ln_nea, ln_bound_reject,    &
          &                ln_ignmis, nn_2dint, nn_1dint, nn_profdavtypes,       &
-         &                ln_fp_indegs, rn_avglamscl, rn_avgphiscl, ln_sstbias, &
-         &                cn_sstbiasfiles, ln_night, ln_altbias,                &
+         &                ln_fp_indegs, rn_avglamscl, rn_avgphiscl, ln_obsbias, &
+         &                cn_obsbiasfiles, cn_type_to_biascorrect,              &
+         &                cn_obsbiasfile_varname, ln_night, ln_altbias,         &
          &                cn_altbiasfile, nn_msshc, rn_mdtcorr, rn_mdtcutoff,   &
-         &                ln_time_mean_sla_bkg, ln_s_at_t, ln_output_clim
+         &                ln_all_at_all
       !!----------------------------------------------------------------------
 
       cn_obstypes(:)     = ''
       cn_obsfiles(:)     = ''
-      cn_sstbiasfiles(:) = ''
+      cn_obsbiasfiles(:) = ''
       nn_profdavtypes(:) = -1
 
@@ -229 +225 @@
          sdobsgroup%nobstypes     = 0
          sdobsgroup%nobsfiles     = 0
+         sdobsgroup%naddvars      = 0
          sdobsgroup%nextvars      = 0
          sdobsgroup%navtypes      = 0
-         sdobsgroup%nsstbiasfiles = 0
-         sdobsgroup%lsst          = .false.
-         sdobsgroup%ltemp3d       = .false.
-         sdobsgroup%lsal3d        = .false.
+         sdobsgroup%nobsbiasfiles = 0
          sdobsgroup%lvel3d        = .false.
          sdobsgroup%lsla          = .false.
@@ -254 +248 @@
          END DO
          DO jfile = 1, jpmaxnfiles
-            IF ( TRIM(cn_sstbiasfiles(jfile)) /= '' ) THEN
-               sdobsgroup%nsstbiasfiles = sdobsgroup%nsstbiasfiles + 1
+            IF ( TRIM(cn_obsbiasfiles(jfile)) /= '' ) THEN
+               sdobsgroup%nobsbiasfiles = sdobsgroup%nobsbiasfiles + 1
             ENDIF
          END DO
@@ -266 +260 @@
                itype = itype + 1
                sdobsgroup%cobstypes(itype) = TRIM(cn_obstypes(jtype))
-               IF ( TRIM(sdobsgroup%cobstypes(itype)) == cobsname_sst ) THEN
-                  sdobsgroup%lsst = .true.
-               ELSEIF ( TRIM(sdobsgroup%cobstypes(itype)) == cobsname_temp3d ) THEN
-                  sdobsgroup%ltemp3d = .true.
-                  sdobsgroup%nextvars = sdobsgroup%nextvars + 1
-               ELSEIF ( TRIM(sdobsgroup%cobstypes(itype)) == cobsname_sal3d ) THEN
-                  sdobsgroup%lsal3d = .true.
-               ELSEIF ( (TRIM(sdobsgroup%cobstypes(itype)) == cobsname_uvel3d) .OR. &
-                  &     (TRIM(sdobsgroup%cobstypes(itype)) == cobsname_vvel3d) ) THEN
+               IF ( (TRIM(sdobsgroup%cobstypes(itype)) == cobsname_uvel3d) .OR. &
+                  & (TRIM(sdobsgroup%cobstypes(itype)) == cobsname_vvel3d) ) THEN
                   sdobsgroup%lvel3d = .true.
-                  sdobsgroup%nextvars = sdobsgroup%nextvars + 1
                ELSEIF ( TRIM(sdobsgroup%cobstypes(itype)) == cobsname_sla ) THEN
                   sdobsgroup%lsla = .true.
-                  sdobsgroup%nextvars = sdobsgroup%nextvars + 2
+! THESE WILL EACH NEED TO BE 1 (ADD=SSH, EXT=MDT)
+                  sdobsgroup%naddvars = 0
+                  sdobsgroup%nextvars = 0
+! DO THIS FOR FBD TOO
                ENDIF
                !
@@ -313 +302 @@
          ifile = 0
          DO jfile = 1, jpmaxnfiles
-            IF ( TRIM(cn_sstbiasfiles(jfile)) /= '' ) THEN
+            IF ( TRIM(cn_obsbiasfiles(jfile)) /= '' ) THEN
                ifile = ifile + 1
-               sdobsgroup%csstbiasfiles(ifile) = cn_sstbiasfiles(jfile)
-            ENDIF
-         END DO
+               sdobsgroup%cobsbiasfiles(ifile) = cn_obsbiasfiles(jfile)
+            ENDIF
+         END DO
+         IF ( ln_obsbias ) THEN
+            sdobsgroup%nbiasvar = -1
+            DO jtype = 1, sdobsgroup%nobstypes
+               IF ( TRIM(sdobsgroup%cobstypes(itype)) == TRIM(cn_type_to_biascorrect) ) THEN
+                  sdobsgroup%nbiasvar = jtype
+                  EXIT
+               ENDIF
+            ENDDO
+         ENDIF
 
          sdobsgroup%caltbiasfile       = cn_altbiasfile
@@ -326 +324 @@
          sdobsgroup%lprof              = ln_prof
          sdobsgroup%laltbias           = ln_altbias
-         sdobsgroup%lsstbias           = ln_sstbias
+         sdobsgroup%lobsbias           = ln_obsbias
+         sdobsgroup%cbiasvarname       = cn_obsbiasfile_varname
          sdobsgroup%lnea               = ln_nea
          sdobsgroup%lbound_reject      = ln_bound_reject
          sdobsgroup%lignmis            = ln_ignmis
-         sdobsgroup%ls_at_t            = ln_s_at_t
+         sdobsgroup%lall_at_all        = ln_all_at_all
          sdobsgroup%lnight             = ln_night
-         sdobsgroup%loutput_clim       = ln_output_clim
-         sdobsgroup%ltime_mean_sla_bkg = ln_time_mean_sla_bkg
          sdobsgroup%lfp_indegs         = ln_fp_indegs
          sdobsgroup%ravglamscl         = rn_avglamscl
@@ -395 +392 @@
             WRITE(numout,*) '             N/S diameter of obs footprint                      rn_avgphiscl = ', sdobsgroup%ravgphiscl
             WRITE(numout,*) '             Logical switch for night-time average                  ln_night = ', sdobsgroup%lnight
+            WRITE(numout,*) '             Logical switch for bias correction                   ln_obsbias = ', sdobsgroup%lobsbias
+            IF ( sdobsgroup%lobsbias ) THEN
+               WRITE(numout,*) '             Observation type to be bias corrected    cn_type_to_biascorrect = ', TRIM(sdobsgroup%cobstypes(sdobsgroup%nbiasvar))
+               WRITE(numout,*) '             Bias variable name in bias files         cn_obsbiasfile_varname = ', TRIM(sdobsgroup%cbiasvarname)
+               WRITE(numout,*) '             Bias files in group:', sdobsgroup%nobsbiasfiles
+               DO jfile = 1, sdobsgroup%nobsbiasfiles
+                  WRITE(numout,*) '                ', TRIM(sdobsgroup%cobsbiasfiles(jfile))
+               END DO
+            ENDIF
             WRITE(numout,*) '          Settings only for profile data, which is ', sdobsgroup%lprof
             WRITE(numout,*) '             Type of vertical interpolation method                  nn_1dint = ', sdobsgroup%n1dint
             WRITE(numout,*) '             Daily average types                             nn_profdavtypes = ', sdobsgroup%nprofdavtypes
-            WRITE(numout,*) '          Settings only for SST data, which is ', sdobsgroup%lsst
-            WRITE(numout,*) '             Logical switch for sst bias                          ln_sstbias = ', sdobsgroup%lsstbias
-            IF ( sdobsgroup%lsstbias ) THEN
-               WRITE(numout,*) '             SST bias files in group:'
-               DO jfile = 1, sdobsgroup%nsstbiasfiles
-                  WRITE(numout,*) '                ', TRIM(sdobsgroup%csstbiasfiles(jfile))
-               END DO
-            ENDIF
+            WRITE(numout,*) '             Logical switch to compute all vars at all pts     ln_all_at_all = ', sdobsgroup%lall_at_all
             WRITE(numout,*) '          Settings only for SLA data, which is ', sdobsgroup%lsla
             WRITE(numout,*) '             Logical switch for alt bias                          ln_altbias = ', sdobsgroup%laltbias
             WRITE(numout,*) '             Alt bias file name                               cn_altbiasfile = ', TRIM(sdobsgroup%caltbiasfile)
-            WRITE(numout,*) '             Logical switch for time-mean of SLA        ln_time_mean_sla_bkg = ', sdobsgroup%ltime_mean_sla_bkg
             WRITE(numout,*) '             MSSH correction scheme                                 nn_msshc = ', sdobsgroup%nmsshc
             WRITE(numout,*) '             MDT  correction                                      rn_mdtcorr = ', sdobsgroup%rmdtcorr
             WRITE(numout,*) '             MDT cutoff for computed correction                 rn_mdtcutoff = ', sdobsgroup%rmdtcutoff
-            WRITE(numout,*) '          Settings only for 3D temperature/salinity data, temperature is ', sdobsgroup%ltemp3d
-            WRITE(numout,*) '                                                             salinity is ', sdobsgroup%lsal3d
-            WRITE(numout,*) '             Logical switch to compute model S at T obs            ln_s_at_t = ', sdobsgroup%ls_at_t
-            WRITE(numout,*) '             Logical switch for writing climat. at obs points ln_output_clim = ', sdobsgroup%loutput_clim
          ENDIF
 
@@ -432 +426 @@
          ENDIF
 
-         IF ( (sdobsgroup%lsst) .AND. (sdobsgroup%lsstbias) .AND. (sdobsgroup%nsstbiasfiles == 0) ) THEN
-            CALL ctl_stop( ' No SST bias files specified for this observation group' )
+         IF ( (sdobsgroup%lobsbias) .AND. (sdobsgroup%nobsbiasfiles == 0) ) THEN
+            CALL ctl_stop( ' No bias files specified for this observation group' )
          ENDIF
 

NEMO/branches/UKMO/NEMO_4.0.4_generic_obs/src/OCE/OBS/obs_oper.F90

@@ -103 +103 @@
       INTEGER       , INTENT(in   ) ::   k2dint          ! Horizontal interpolation type (see header)
       INTEGER       , INTENT(in   ) ::   kdaystp         ! Number of time steps per day
-      INTEGER       , INTENT(in   ) ::   kvar            ! Number of variables in prodatqc
+      INTEGER       , INTENT(in   ) ::   kvar            ! Index of variable in prodatqc
       REAL(KIND=wp) , INTENT(in   ), DIMENSION(kpi,kpj,kpk) ::   pvar             ! Model field
       REAL(KIND=wp) , INTENT(in   ), DIMENSION(kpi,kpj,kpk) ::   pmask            ! Land-sea mask
@@ -450 +450 @@
    END SUBROUTINE obs_prof_opt
 
-   SUBROUTINE obs_surf_opt( surfdataqc, kt, kpi, kpj,            &
-      &                     kit000, kdaystp, psurf, psurfmask,   &
-      &                     k2dint, ldnightav, plamscl, pphiscl, &
+   SUBROUTINE obs_surf_opt( surfdataqc, kt, kpi, kpj,                &
+      &                     kit000, kvar, kdaystp, psurf, psurfmask, &
+      &                     k2dint, ldnightav, plamscl, pphiscl,     &
       &                     lindegrees )
 
@@ -499 +499 @@
                                        !   (kit000-1 = restart time)
       INTEGER, INTENT(IN) :: kdaystp   ! Number of time steps per day
+      INTEGER, INTENT(IN) :: kvar      ! Index of variable in surfdataqc  
       INTEGER, INTENT(IN) :: k2dint    ! Horizontal interpolation type (see header)
       REAL(wp), INTENT(IN), DIMENSION(kpi,kpj) :: &
@@ -563 +564 @@
       IF ( ldnightav ) THEN
 
-      ! Initialize array for night mean
+         ! Initialize array for night mean
          IF ( kt == 0 ) THEN
             ALLOCATE ( icount_night(kpi,kpj) )
@@ -581 +582 @@
             DO jj = 1, jpj
                DO ji = 1, jpi
-                  surfdataqc%vdmean(ji,jj) = 0.0
+                  surfdataqc%vdmean(ji,jj,:) = 0.0
                   zmeanday(ji,jj) = 0.0
                   icount_night(ji,jj) = 0
@@ -594 +595 @@
          DO jj = 1, jpj
             DO ji = 1, jpi
-               ! Increment the temperature field for computing night mean and counter
-               surfdataqc%vdmean(ji,jj) = surfdataqc%vdmean(ji,jj)  &
-                      &                    + psurf(ji,jj) * REAL( imask_night(ji,jj) )
+               ! Increment the model field for computing night mean and counter
+               surfdataqc%vdmean(ji,jj,kvar) = surfdataqc%vdmean(ji,jj,kvar)  &
+                      &                        + psurf(ji,jj) * REAL( imask_night(ji,jj) )
                zmeanday(ji,jj)          = zmeanday(ji,jj) + psurf(ji,jj)
                icount_night(ji,jj)      = icount_night(ji,jj) + imask_night(ji,jj)
@@ -610 +611 @@
                   ! Test if "no night" point
                   IF ( icount_night(ji,jj) > 0 ) THEN
-                     surfdataqc%vdmean(ji,jj) = surfdataqc%vdmean(ji,jj) &
-                       &                        / REAL( icount_night(ji,jj) )
+                     surfdataqc%vdmean(ji,jj,kvar) = surfdataqc%vdmean(ji,jj,kvar) &
+                       &                             / REAL( icount_night(ji,jj) )
                   ELSE
                      !At locations where there is no night (e.g. poles),
                      ! calculate daily mean instead of night-time mean.
-                     surfdataqc%vdmean(ji,jj) = zmeanday(ji,jj) * zdaystp
+                     surfdataqc%vdmean(ji,jj,kvar) = zmeanday(ji,jj) * zdaystp
                   ENDIF
                END DO
@@ -689 +690 @@
 
          CALL obs_int_comm_2d( imaxifp,imaxjfp, isurf, kpi, kpj, igrdi, igrdj, &
-            &               surfdataqc%vdmean(:,:), zsurfm )
+            &               surfdataqc%vdmean(:,:,kvar), zsurfm )
 
       ENDIF
@@ -750 +751 @@
 
          ENDIF
-
+! WHERE BEST TO DO THIS?
          IF ( TRIM(surfdataqc%cvars(1)) == 'SLA' .AND. surfdataqc%nextra == 2 ) THEN
             ! ... Remove the MDT from the SSH at the observation point to get the SLA

NEMO/branches/UKMO/NEMO_4.0.4_generic_obs/src/OCE/OBS/obs_prep.F90

@@ -42 +42 @@
 CONTAINS
 
-   SUBROUTINE obs_pre_surf( surfdata, surfdataqc, ld_nea, ld_bound_reject, &
-                            kqc_cutoff )
-      !!----------------------------------------------------------------------
-      !!                    ***  ROUTINE obs_pre_sla  ***
+   SUBROUTINE obs_pre_surf( surfdata, surfdataqc, &
+      &                     kpi, kpj, &   
+      &                     zmask, pglam, pgphi, &
+      &                     ld_nea, ld_bound_reject, &
+      &                     kqc_cutoff )
+      !!----------------------------------------------------------------------
+      !!                    ***  ROUTINE obs_pre_surf  ***
       !!
       !! ** Purpose : First level check and screening of surface observations
@@ -65 +68 @@
       !! * Arguments
       TYPE(obs_surf), INTENT(INOUT) :: surfdata    ! Full set of surface data
-      TYPE(obs_surf), INTENT(INOUT) :: surfdataqc   ! Subset of surface data not failing screening
-      LOGICAL, INTENT(IN) :: ld_nea         ! Switch for rejecting observation near land
+      TYPE(obs_surf), INTENT(INOUT) :: surfdataqc  ! Subset of surface data not failing screening
+      INTEGER, INTENT(IN) :: kpi, kpj              ! Local domain sizes      
+      REAL(wp), INTENT(IN), DIMENSION(kpi,kpj,surfdata%nvar) :: &
+         & zmask      
+      REAL(wp), INTENT(IN), DIMENSION(kpi,kpj,surfdata%nvar) :: &
+         & pglam, &
+         & pgphi
+      LOGICAL, INTENT(IN) :: ld_nea                ! Switch for rejecting observation near land
       LOGICAL, INTENT(IN) :: ld_bound_reject       ! Switch for rejecting obs near the boundary
-      INTEGER, INTENT(IN), OPTIONAL :: kqc_cutoff   ! cut off for QC value
+      INTEGER, INTENT(IN), OPTIONAL :: kqc_cutoff  ! cut off for QC value
       !! * Local declarations
       INTEGER :: iqc_cutoff = 255   ! cut off for QC value
@@ -77 +86 @@
       INTEGER :: imin0
       INTEGER :: icycle       ! Current assimilation cycle
-                              ! Counters for observations that
-      INTEGER :: iotdobs      !  - outside time domain
-      INTEGER :: iosdsobs     !  - outside space domain
-      INTEGER :: ilansobs     !  - within a model land cell
-      INTEGER :: inlasobs     !  - close to land
-      INTEGER :: igrdobs      !  - fail the grid search
-      INTEGER :: ibdysobs     !  - close to open boundary
-                              ! Global counters for observations that
-      INTEGER :: iotdobsmpp     !  - outside time domain
-      INTEGER :: iosdsobsmpp    !  - outside space domain
-      INTEGER :: ilansobsmpp    !  - within a model land cell
-      INTEGER :: inlasobsmpp    !  - close to land
-      INTEGER :: igrdobsmpp     !  - fail the grid search
-      INTEGER :: ibdysobsmpp  !  - close to open boundary
+                                                        ! Counters for observations that are
+      INTEGER                           :: iotdobs      !  - outside time domain
+      INTEGER, DIMENSION(surfdata%nvar) :: iosdsobs     !  - outside space domain
+      INTEGER, DIMENSION(surfdata%nvar) :: ilansobs     !  - within a model land cell
+      INTEGER, DIMENSION(surfdata%nvar) :: inlasobs     !  - close to land
+      INTEGER, DIMENSION(surfdata%nvar) :: ibdysobs     !  - close to open boundary
+      INTEGER                           :: igrdobs      !  - fail the grid search      
+                                                        ! Global counters for observations that
+      INTEGER                           :: iotdobsmpp   !  - outside time domain
+      INTEGER, DIMENSION(surfdata%nvar) :: iosdsobsmpp  !  - outside space domain
+      INTEGER, DIMENSION(surfdata%nvar) :: ilansobsmpp  !  - within a model land cell
+      INTEGER, DIMENSION(surfdata%nvar) :: inlasobsmpp  !  - close to land
+      INTEGER, DIMENSION(surfdata%nvar) :: ibdysobsmpp  !  - close to open boundary
+      INTEGER                           :: igrdobsmpp   !  - fail the grid search
+
       LOGICAL, DIMENSION(:), ALLOCATABLE :: & 
          & llvalid            ! SLA data selection
-      INTEGER :: jobs         ! Obs. loop variable
+      INTEGER :: jobs         ! Obs. loop counter
+      INTEGER :: jvar         ! Variable loop counter
       INTEGER :: jstp         ! Time loop variable
       INTEGER :: inrc         ! Time index variable
+      CHARACTER(LEN=256) :: cout1  ! Diagnostic output line
       !!----------------------------------------------------------------------
 
@@ -110 +122 @@
       icycle = nn_no     ! Assimilation cycle
 
-      ! Diagnotics counters for various failures.
+      ! Diagnostic counters for various failures.
 
       iotdobs  = 0
       igrdobs  = 0
-      iosdsobs = 0
-      ilansobs = 0
-      inlasobs = 0
-      ibdysobs = 0 
+      iosdsobs(:) = 0
+      ilansobs(:) = 0
+      inlasobs(:) = 0
+      ibdysobs(:) = 0 
 
       ! Set QC cutoff to optional value if provided
@@ -147 +159 @@
       ! -----------------------------------------------------------------------
 
-      CALL obs_coo_spc_2d( surfdata%nsurf,              &
-         &                 jpi,          jpj,          &
-         &                 surfdata%mi,   surfdata%mj,   & 
-         &                 surfdata%rlam, surfdata%rphi, &
-         &                 glamt,        gphit,        &
-         &                 tmask(:,:,1), surfdata%nqc,  &
-         &                 iosdsobs,     ilansobs,     &
-         &                 inlasobs,     ld_nea,       &
-         &                 ibdysobs,     ld_bound_reject, &
-         &                 iqc_cutoff                     )
-
-      CALL obs_mpp_sum_integer( iosdsobs, iosdsobsmpp )
-      CALL obs_mpp_sum_integer( ilansobs, ilansobsmpp )
-      CALL obs_mpp_sum_integer( inlasobs, inlasobsmpp )
-      CALL obs_mpp_sum_integer( ibdysobs, ibdysobsmpp )
+      DO jvar = 1, surfdata%nvar
+         CALL obs_coo_spc_2d( surfdata%nsurf,                   &
+            &                 jpi,             jpj,             &
+            &                 surfdata%mi,     surfdata%mj,     & 
+            &                 surfdata%rlam,   surfdata%rphi,   &
+            &                 pglam(:,:,jvar), pgphi(:,:,jvar), &
+            &                 zmask(:,:,jvar), surfdata%nqc,    &
+            &                 iosdsobs(jvar),  ilansobs(jvar),  &
+            &                 inlasobs(jvar),  ld_nea,          &
+            &                 ibdysobs(jvar),  ld_bound_reject, &
+            &                 iqc_cutoff                     )
+         CALL obs_mpp_sum_integer( iosdsobs(jvar), iosdsobsmpp(jvar) )
+         CALL obs_mpp_sum_integer( ilansobs(jvar), ilansobsmpp(jvar) )
+         CALL obs_mpp_sum_integer( inlasobs(jvar), inlasobsmpp(jvar) )
+         CALL obs_mpp_sum_integer( ibdysobs(jvar), ibdysobsmpp(jvar) )
+      END DO
 
       ! -----------------------------------------------------------------------
@@ -191 +204 @@
       
       IF(lwp) THEN
+         DO jvar = 1, surfdataqc%nvar       
+            IF ( jvar == 1 ) THEN
+               cout1=TRIM(surfdataqc%cvars(1))                  
+            ELSE
+               WRITE(cout1,'(A,A1,A)') TRIM(cout1), '/', TRIM(surfdataqc%cvars(jvar))            
+            ENDIF
+         END DO
+               
          WRITE(numout,*)
-         WRITE(numout,*) ' '//surfdataqc%cvars(1)//' data outside time domain                  = ', &
+         WRITE(numout,*) ' '//TRIM(cout1)//' data outside time domain                  = ', &
             &            iotdobsmpp
-         WRITE(numout,*) ' Remaining '//surfdataqc%cvars(1)//' data that failed grid search    = ', &
+         WRITE(numout,*) ' Remaining '//TRIM(cout1)//' data that failed grid search    = ', &
             &            igrdobsmpp
-         WRITE(numout,*) ' Remaining '//surfdataqc%cvars(1)//' data outside space domain       = ', &
-            &            iosdsobsmpp
-         WRITE(numout,*) ' Remaining '//surfdataqc%cvars(1)//' data at land points             = ', &
-            &            ilansobsmpp
-         IF (ld_nea) THEN
-            WRITE(numout,*) ' Remaining '//surfdataqc%cvars(1)//' data near land points (removed) = ', &
-               &            inlasobsmpp
-         ELSE
-            WRITE(numout,*) ' Remaining '//surfdataqc%cvars(1)//' data near land points (kept)    = ', &
-               &            inlasobsmpp
-         ENDIF
-         WRITE(numout,*) ' Remaining '//surfdataqc%cvars(1)//' data near open boundary (removed) = ', &
-            &            ibdysobsmpp  
-         WRITE(numout,*) ' '//surfdataqc%cvars(1)//' data accepted                             = ', &
-            &            surfdataqc%nsurfmpp
+
+         DO jvar = 1, surfdataqc%nvar            
+            WRITE(numout,*) ' Remaining '//surfdataqc%cvars(jvar)//' data outside space domain       = ', &
+                &            iosdsobsmpp(jvar)
+             WRITE(numout,*) ' Remaining '//surfdataqc%cvars(jvar)//' data at land points             = ', &
+                &            ilansobsmpp(jvar)
+             IF (ld_nea) THEN
+                WRITE(numout,*) ' Remaining '//surfdataqc%cvars(jvar)//' data near land points (removed) = ', &
+                   &            inlasobsmpp(jvar)
+             ELSE
+                WRITE(numout,*) ' Remaining '//surfdataqc%cvars(jvar)//' data near land points (kept)    = ', &
+                   &            inlasobsmpp(jvar)
+             ENDIF     
+             WRITE(numout,*) ' Remaining '//surfdataqc%cvars(jvar)//' data near open boundary (removed) = ', &
+                &            ibdysobsmpp(jvar)
+          END DO
+          WRITE(numout,*) ' '//TRIM(cout1)//' data accepted                             = ', &
+             &            surfdataqc%nsurfmpp
 
          WRITE(numout,*)
          WRITE(numout,*) ' Number of observations per time step :'
          WRITE(numout,*)
-         WRITE(numout,'(10X,A,10X,A)')'Time step',surfdataqc%cvars(1)
+         WRITE(numout,'(10X,A,10X,A)')'Time step',TRIM(cout1)
          WRITE(numout,'(10X,A,5X,A)')'---------','-----------------'
          CALL FLUSH(numout)

NEMO/branches/UKMO/NEMO_4.0.4_generic_obs/src/OCE/OBS/obs_profiles_def.F90

@@ -75 +75 @@
 
       REAL(KIND=wp), POINTER, DIMENSION(:,:) :: &
-         & vext           !: Extra variables
+         & vadd           !: Additional variables
 
       INTEGER, POINTER, DIMENSION(:) :: &
@@ -94 +94 @@
 
       INTEGER :: nvar     !: Number of variables
-      INTEGER :: next     !: Number of extra fields
+      INTEGER :: next     !: Number of extra variables
+      INTEGER :: nadd     !: Number of additional variables
       INTEGER :: nprof    !: Total number of profiles within window.
       INTEGER :: nstp     !: Number of time steps
@@ -104 +105 @@
       ! Bookkeeping arrays with sizes equal to number of variables
 
-      CHARACTER(len=8), POINTER, DIMENSION(:) :: &
-         & cvars          !: Variable names
+      CHARACTER(len=ilenname), POINTER, DIMENSION(:) :: &
+         & cvars,    &    !: Variable names
+         & cextvars, &    !: Extra variable names
+         & caddvars       !: Additional variable names
+
+      CHARACTER(len=ilenlong), POINTER, DIMENSION(:) :: &
+         & clong,    &    !: Variable long names
+         & cextlong       !: Extra variable long names
+
+      CHARACTER(len=ilenlong), POINTER, DIMENSION(:,:) :: &
+         & caddlong       !: Additional variable long names
+
+      CHARACTER(len=ilenunit), POINTER, DIMENSION(:) :: &
+         & cunit,    &    !: Variable units
+         & cextunit       !: Extra variable units
+
+      CHARACTER(len=ilenunit), POINTER, DIMENSION(:,:) :: &
+         & caddunit       !: Additional variable units
+
+      CHARACTER(len=ilengrid), POINTER, DIMENSION(:) :: &
+         & cgrid          !: Variable grids
 
       INTEGER, POINTER, DIMENSION(:) :: &
@@ -131 +151 @@
          & rphi           !: Latitude coordinate of profile data
 
-      CHARACTER(LEN=8), POINTER, DIMENSION(:) :: &
+      CHARACTER(LEN=ilenwmo), POINTER, DIMENSION(:) :: &
          & cwmo           !: Profile WMO indentifier
       
@@ -160 +180 @@
       TYPE(obs_prof_var), POINTER, DIMENSION(:) :: var
 
+      REAL(KIND=wp), POINTER, DIMENSION(:,:) :: &
+         & vext        !: Extra variables
+
       ! Arrays with size equal to the number of time steps in the window
 
@@ -197 +220 @@
 CONTAINS
    
-   SUBROUTINE obs_prof_alloc( prof,  kvar, kext, kprof,  &
+   SUBROUTINE obs_prof_alloc( prof,  kvar, kadd, kext, kprof,  &
       &                       ko3dt, kstp, kpi, kpj, kpk )
       !!----------------------------------------------------------------------
@@ -214 +237 @@
       INTEGER, INTENT(IN) :: kprof  ! Number of profiles
       INTEGER, INTENT(IN) :: kvar   ! Number of variables
-      INTEGER, INTENT(IN) :: kext   ! Number of extra fields within each variable
+      INTEGER, INTENT(IN) :: kadd   ! Number of additional fields within each variable
+      INTEGER, INTENT(IN) :: kext   ! Number of extra fields
       INTEGER, INTENT(IN), DIMENSION(kvar) :: &
          & ko3dt     ! Number of observations per variables
@@ -223 +247 @@
 
       !!* Local variables
-      INTEGER :: jvar
+      INTEGER :: jvar, jadd, jext
       INTEGER :: ji
 
@@ -229 +253 @@
 
       prof%nvar      = kvar
+      prof%nadd      = kadd
       prof%next      = kext
       prof%nprof     = kprof
@@ -241 +266 @@
       ALLOCATE( &
          & prof%cvars(kvar),    &
+         & prof%clong(kvar),    &
+         & prof%cunit(kvar),    &
+         & prof%cgrid(kvar),    &
          & prof%nvprot(kvar),   &
          & prof%nvprotmpp(kvar) &
@@ -247 +275 @@
       DO jvar = 1, kvar
          prof%cvars    (jvar) = "NotSet"
+         prof%clong    (jvar) = "NotSet"
+         prof%cunit    (jvar) = "NotSet"
+         prof%cgrid    (jvar) = ""
          prof%nvprot   (jvar) = ko3dt(jvar)
          prof%nvprotmpp(jvar) = 0
+      END DO
+
+      ! Allocate additional/extra variable metadata
+
+      ALLOCATE( &
+         & prof%caddvars(kadd),      &
+         & prof%caddlong(kadd,kvar), &
+         & prof%caddunit(kadd,kvar), &
+         & prof%cextvars(kext),      &
+         & prof%cextlong(kext),      &
+         & prof%cextunit(kext)       &
+         )
+         
+      DO jadd = 1, kadd
+         prof%caddvars(jadd) = "NotSet"
+         DO jvar = 1, kvar
+            prof%caddlong(jadd,jvar) = "NotSet"
+            prof%caddunit(jadd,jvar) = "NotSet"
+         END DO
+      END DO
+         
+      DO jext = 1, kext
+         prof%cextvars(jext) = "NotSet"
+         prof%cextlong(jext) = "NotSet"
+         prof%cextunit(jext) = "NotSet"
       END DO
 
@@ -308 +364 @@
 
          IF ( ko3dt(jvar) >= 0 ) THEN
-            CALL obs_prof_alloc_var( prof, jvar, kext, ko3dt(jvar) )
+            CALL obs_prof_alloc_var( prof, jvar, kadd, ko3dt(jvar) )
          ENDIF
          
       END DO
+      
+      ! Allocate extra variables
+      ALLOCATE( &
+         & prof%vext(kprof,kext) &
+         & )
 
       ! Allocate arrays of size number of time step size
@@ -432 +493 @@
          & )
 
+      ! Deallocate extra variables
+      DEALLOCATE( &
+         & prof%vext &
+         & )
+      
       ! Deallocate arrays of size number of time step size
 
@@ -458 +524 @@
       DEALLOCATE( &
          & prof%cvars,    &
+         & prof%clong,    &
+         & prof%cunit,    &
+         & prof%cgrid,    &
          & prof%nvprot,   &
          & prof%nvprotmpp &
          )
 
+      ! Dellocate additional/extra variables metadata
+
+      DEALLOCATE( &
+         & prof%caddvars, &
+         & prof%caddlong, &
+         & prof%caddunit, &
+         & prof%cextvars, &
+         & prof%cextlong, &
+         & prof%cextunit  &
+         )
+
 
    END SUBROUTINE obs_prof_dealloc
 
 
-   SUBROUTINE obs_prof_alloc_var( prof, kvar, kext, kobs )
+   SUBROUTINE obs_prof_alloc_var( prof, kvar, kadd, kobs )
 
       !!----------------------------------------------------------------------
@@ -480 +560 @@
       TYPE(obs_prof), INTENT(INOUT) :: prof   ! Profile data to be allocated
       INTEGER, INTENT(IN) :: kvar   ! Variable number
-      INTEGER, INTENT(IN) :: kext   ! Number of extra fields within each variable
+      INTEGER, INTENT(IN) :: kadd   ! Number of additional fields within each variable
       INTEGER, INTENT(IN) :: kobs   ! Number of observations
       
@@ -498 +578 @@
          & prof%var(kvar)%nvqcf(idefnqcf,kobs)  &
          & )
-      IF (kext>0) THEN
+      IF (kadd>0) THEN
          ALLOCATE( & 
-            & prof%var(kvar)%vext(kobs,kext) &
+            & prof%var(kvar)%vadd(kobs,kadd) &
             & )
       ENDIF
@@ -534 +614 @@
          & prof%var(kvar)%nvqcf   &
          & )
-      IF (prof%next>0) THEN
+      IF (prof%nadd>0) THEN
          DEALLOCATE( & 
-            & prof%var(kvar)%vext  &
+            & prof%var(kvar)%vadd  &
             & )
       ENDIF
@@ -576 +656 @@
          & invpro
       INTEGER :: jvar
+      INTEGER :: jadd
       INTEGER :: jext
       INTEGER :: ji
@@ -627 +708 @@
       IF ( lallocate ) THEN
          CALL obs_prof_alloc( newprof,   prof%nvar, &
-            &                 prof%next,            &
+            &                 prof%nadd, prof%next, &
             &                 inprof,    invpro,    &
             &                 prof%nstp, prof%npi,  &
@@ -670 +751 @@
 
             newprof%mi(inprof,:)  = prof%mi(ji,:)
-            newprof%mj(inprof,:) = prof%mj(ji,:)
+            newprof%mj(inprof,:)  = prof%mj(ji,:)
             newprof%npidx(inprof) = prof%npidx(ji)
             newprof%npfil(inprof) = prof%npfil(ji)
@@ -741 +822 @@
                      newprof%var(jvar)%vmod(invpro(jvar))   = &
                         &                           prof%var(jvar)%vmod(jj)
-                     DO jext = 1, prof%next
-                        newprof%var(jvar)%vext(invpro(jvar),jext) = &
-                           &                      prof%var(jvar)%vext(jj,jext)
+                     DO jadd = 1, prof%nadd
+                        newprof%var(jvar)%vadd(invpro(jvar),jadd) = &
+                           &                      prof%var(jvar)%vadd(jj,jadd)
                      END DO
                   
@@ -756 +837 @@
             END DO
 
+            DO jext = 1, prof%next
+               newprof%vext(inprof,jext) = prof%vext(ji,jext)
+            END DO
+
          ENDIF
 
@@ -771 +856 @@
 
       newprof%nvar     = prof%nvar
+      newprof%nadd     = prof%nadd
       newprof%next     = prof%next
       newprof%nstp     = prof%nstp
@@ -777 +863 @@
       newprof%npk      = prof%npk
       newprof%cvars(:) = prof%cvars(:)
+      newprof%clong(:) = prof%clong(:)
+      newprof%cunit(:) = prof%cunit(:)
+      newprof%cgrid(:) = prof%cgrid(:)
+      newprof%caddvars(:) = prof%caddvars(:)
+      newprof%caddlong(:) = prof%caddlong(:)
+      newprof%caddunit(:) = prof%caddunit(:)
+      newprof%cextvars(:) = prof%cextvars(:)
+      newprof%cextlong(:) = prof%cextlong(:)
+      newprof%cextunit(:) = prof%cextunit(:)
  
       ! Deallocate temporary data
@@ -810 +905 @@
       !!* Local variables
       INTEGER :: jvar
+      INTEGER :: jadd
       INTEGER :: jext
       INTEGER :: ji
@@ -866 +962 @@
                oldprof%var(jvar)%idqcf(:,jl) = prof%var(jvar)%idqcf(:,jj)
                oldprof%var(jvar)%nvqcf(:,jl) = prof%var(jvar)%nvqcf(:,jj)
-               DO jext = 1, prof%next
-                  oldprof%var(jvar)%vext(jl,jext) = &
-                     &                        prof%var(jvar)%vext(jj,jext)
+               DO jadd = 1, prof%nadd
+                  oldprof%var(jvar)%vadd(jl,jadd) = &
+                     &                        prof%var(jvar)%vadd(jj,jadd)
                END DO
                
             END DO
 
+         END DO
+
+         DO jext = 1, prof%next
+            oldprof%vext(jk,jext) = prof%vext(jj,jext)
          END DO
          

NEMO/branches/UKMO/NEMO_4.0.4_generic_obs/src/OCE/OBS/obs_read_prof.F90

@@ -44 +44 @@
 
    SUBROUTINE obs_rea_prof( profdata, knumfiles, cdfilenames, &
-      &                     kvars, kextr, kstp, ddobsini, ddobsend, &
-      &                     ldvar, ldignmis, ldsatt, &
+      &                     kvars, kadd, kextr, kstp, ddobsini, ddobsend, &
+      &                     ldvar, ldignmis, ldallatall, &
       &                     ldmod, cdvars, kdailyavtypes )
       !!---------------------------------------------------------------------
@@ -72 +72 @@
          & cdfilenames(knumfiles)        ! File names to read in
       INTEGER, INTENT(IN) :: kvars      ! Number of variables in profdata
-      INTEGER, INTENT(IN) :: kextr      ! Number of extra fields for each var
+      INTEGER, INTENT(IN) :: kadd       ! Number of additional fields
+                                        !   in addition to those in the input file(s)
+      INTEGER, INTENT(IN) :: kextr      ! Number of extra fields
+                                        !   in addition to those in the input file(s)
       INTEGER, INTENT(IN) :: kstp       ! Ocean time-step index
       LOGICAL, DIMENSION(kvars), INTENT(IN) :: ldvar     ! Observed variables switches
       LOGICAL, INTENT(IN) :: ldignmis   ! Ignore missing files
-      LOGICAL, INTENT(IN) :: ldsatt     ! Compute salinity at all temperature points
+      LOGICAL, INTENT(IN) :: ldallatall     ! Compute salinity at all temperature points
       LOGICAL, INTENT(IN) :: ldmod      ! Initialize model from input data
       REAL(dp), INTENT(IN) :: ddobsini  ! Obs. ini time in YYYYMMDD.HHMMSS
@@ -87 +90 @@
       CHARACTER(LEN=15), PARAMETER :: cpname='obs_rea_prof'
       CHARACTER(len=8) :: clrefdate
-      CHARACTER(len=8), DIMENSION(:), ALLOCATABLE :: clvarsin
-      INTEGER :: jvar
+      CHARACTER(len=ilenname), DIMENSION(:),   ALLOCATABLE :: clvarsin
+      CHARACTER(len=ilenlong), DIMENSION(:),   ALLOCATABLE :: cllongin
+      CHARACTER(len=ilenunit), DIMENSION(:),   ALLOCATABLE :: clunitin
+      CHARACTER(len=ilengrid), DIMENSION(:),   ALLOCATABLE :: clgridin
+      CHARACTER(len=ilenname), DIMENSION(:),   ALLOCATABLE :: claddvarsin
+      CHARACTER(len=ilenlong), DIMENSION(:,:), ALLOCATABLE :: claddlongin
+      CHARACTER(len=ilenunit), DIMENSION(:,:), ALLOCATABLE :: claddunitin
+      CHARACTER(len=ilenname), DIMENSION(:),   ALLOCATABLE :: clextvarsin
+      CHARACTER(len=ilenlong), DIMENSION(:),   ALLOCATABLE :: clextlongin
+      CHARACTER(len=ilenunit), DIMENSION(:),   ALLOCATABLE :: clextunitin
       INTEGER :: ji
       INTEGER :: jj
       INTEGER :: jk
       INTEGER :: ij
+      INTEGER :: jext
+      INTEGER :: jvar
+      INTEGER :: jadd
+      INTEGER :: jadd2
+      INTEGER :: iadd
+      INTEGER :: iaddin
+      INTEGER :: iextr
       INTEGER :: iflag
       INTEGER :: inobf
@@ -166 +184 @@
       ALLOCATE( inpfiles(inobf) )
 
+      iadd  = 0
+      iextr = 0
+
       prof_files : DO jj = 1, inobf
 
@@ -221 +242 @@
             ENDIF
 
+            IF ( (iextr > 0) .AND. (inpfiles(jj)%next /= iextr) ) THEN
+               CALL ctl_stop( 'Number of extra variables not consistent', &
+                  &           ' with previous files for this type' )
+            ELSE
+               iextr = inpfiles(jj)%next
+            ENDIF
+
+            ! Ignore model counterpart
+            iaddin = inpfiles(jj)%nadd
+            DO ji = 1, iaddin
+               IF ( TRIM(inpfiles(jj)%caddname(ji)) == 'Hx' ) THEN
+                  iaddin = iaddin - 1
+                  EXIT
+               ENDIF
+            END DO
+            IF ( ldmod .AND. ( inpfiles(jj)%nadd == iaddin ) ) THEN
+               CALL ctl_stop( 'Model not in input data' )
+            ENDIF
+
+            IF ( (iadd > 0) .AND. (iaddin /= iadd) ) THEN
+               CALL ctl_stop( 'Number of additional variables not consistent', &
+                  &           ' with previous files for this type' )
+            ELSE
+               iadd = iaddin
+            ENDIF
+
             IF ( jj == 1 ) THEN
                ALLOCATE( clvarsin( inpfiles(jj)%nvar ) )
+               ALLOCATE( cllongin( inpfiles(jj)%nvar ) )
+               ALLOCATE( clunitin( inpfiles(jj)%nvar ) )
+               ALLOCATE( clgridin( inpfiles(jj)%nvar ) )
                DO ji = 1, inpfiles(jj)%nvar
                  clvarsin(ji) = inpfiles(jj)%cname(ji)
+                 cllongin(ji) = inpfiles(jj)%coblong(ji)
+                 clunitin(ji) = inpfiles(jj)%cobunit(ji)
+                 clgridin(ji) = inpfiles(jj)%cgrid(ji)
                  IF ( clvarsin(ji) /= cdvars(ji) ) THEN
                     CALL ctl_stop( 'Feedback file variables do not match', &
@@ -230 +283 @@
                  ENDIF
                END DO
+               IF ( iadd > 0 ) THEN
+                  ALLOCATE( claddvarsin( iadd ) )
+                  ALLOCATE( claddlongin( iadd, inpfiles(jj)%nvar ) )
+                  ALLOCATE( claddunitin( iadd, inpfiles(jj)%nvar ) )
+                  jadd = 0
+                  DO ji = 1, inpfiles(jj)%nadd
+                    IF ( TRIM(inpfiles(jj)%caddname(ji)) /= 'Hx' ) THEN
+                       jadd = jadd + 1
+                       claddvarsin(jadd) = inpfiles(jj)%caddname(ji)
+                       DO jk = 1, inpfiles(jj)%nvar
+                          claddlongin(jadd,jk) = inpfiles(jj)%caddlong(ji,jk)
+                          claddunitin(jadd,jk) = inpfiles(jj)%caddunit(ji,jk)
+                       END DO
+                    ENDIF
+                  END DO
+               ENDIF
+               IF ( iextr > 0 ) THEN
+                  ALLOCATE( clextvarsin( iextr ) )
+                  ALLOCATE( clextlongin( iextr ) )
+                  ALLOCATE( clextunitin( iextr ) )
+                  DO ji = 1, iextr
+                    clextvarsin(ji) = inpfiles(jj)%cextname(ji)
+                    clextlongin(ji) = inpfiles(jj)%cextlong(ji)
+                    clextunitin(ji) = inpfiles(jj)%cextunit(ji)
+                  END DO
+               ENDIF
             ELSE
                DO ji = 1, inpfiles(jj)%nvar
@@ -237 +316 @@
                   ENDIF
                END DO
+               IF ( iadd > 0 ) THEN
+                  jadd = 0
+                  DO ji = 1, inpfiles(jj)%nadd
+                     IF ( TRIM(inpfiles(jj)%caddname(ji)) /= 'Hx' ) THEN
+                        jadd = jadd + 1
+                        IF ( inpfiles(jj)%caddname(ji) /= claddvarsin(jadd) ) THEN
+                           CALL ctl_stop( 'Feedback file additional variables not consistent', &
+                              &           ' with previous files for this type' )
+                        ENDIF
+                     ENDIF
+                  END DO
+               ENDIF
+               IF ( iextr > 0 ) THEN
+                  DO ji = 1, iextr
+                     IF ( inpfiles(jj)%cextname(ji) /= clextvarsin(ji) ) THEN
+                        CALL ctl_stop( 'Feedback file extra variables not consistent', &
+                           &           ' with previous files for this type' )
+                     ENDIF
+                  END DO
+               ENDIF
             ENDIF
 
@@ -499 +598 @@
 
       iv3dt(:) = -1
-      IF (ldsatt) THEN
+      IF (ldallatall) THEN
          iv3dt(:) = ip3dt
       ELSE
          iv3dt(:) = ivart0(:)
       ENDIF
-      CALL obs_prof_alloc( profdata, kvars, kextr, iprof, iv3dt, &
+      CALL obs_prof_alloc( profdata, kvars, kadd+iadd, kextr+iextr, iprof, iv3dt, &
          &                 kstp, jpi, jpj, jpk )
 
@@ -512 +611 @@
       profdata%nvprot(:) = 0
       profdata%cvars(:)  = clvarsin(:)
+      profdata%clong(:)  = cllongin(:)
+      profdata%cunit(:)  = clunitin(:)
+      profdata%cgrid(:)  = clgridin(:)
+      IF ( iadd > 0 ) THEN
+         profdata%caddvars(kadd+1:)   = claddvarsin(:)
+         profdata%caddlong(kadd+1:,:) = claddlongin(:,:)
+         profdata%caddunit(kadd+1:,:) = claddunitin(:,:)
+      ENDIF
+      IF ( iextr > 0 ) THEN
+         profdata%cextvars(kextr+1:) = clextvarsin(:)
+         profdata%cextlong(kextr+1:) = clextlongin(:)
+         profdata%cextunit(kextr+1:) = clextunitin(:)
+      ENDIF
       iprof = 0
 
@@ -644 +756 @@
                      & CYCLE
 
-                  IF (ldsatt) THEN
+                  IF (ldallatall) THEN
 
                      DO jvar = 1, kvars
@@ -663 +775 @@
                      IF ( ( .NOT. BTEST(inpfiles(jj)%ivlqc(ij,ji,jvar),2) .AND. &
                        &   .NOT. BTEST(inpfiles(jj)%idqc(ij,ji),2) .AND. &
-                       &    ldvar(jvar) ) .OR. ldsatt ) THEN
-
-                        IF (ldsatt) THEN
+                       &    ldvar(jvar) ) .OR. ldallatall ) THEN
+
+                        IF (ldallatall) THEN
 
                            ivart(jvar) = ip3dt
@@ -698 +810 @@
                            profdata%var(jvar)%vobs(ivart(jvar)) = &
                               &                inpfiles(jj)%pob(ij,ji,jvar)
-                           IF ( ldmod ) THEN
-                              profdata%var(jvar)%vmod(ivart(jvar)) = &
-                                 &                inpfiles(jj)%padd(ij,ji,1,jvar)
-                           ENDIF
                            ! Count number of profile var1 data as function of type
                            itypvar( profdata%ntyp(iprof) + 1, jvar ) = &
@@ -717 +825 @@
                            & inpfiles(jj)%ivlqcf(:,ij,ji,jvar)
 
-                        ! Profile insitu T value
-                        IF ( TRIM( inpfiles(jj)%cname(jvar) ) == 'POTM' ) THEN
-                           profdata%var(jvar)%vext(ivart(jvar),1) = &
-                              &                inpfiles(jj)%pext(ij,ji,1)
+                        ! Additional variables
+                        IF ( iadd > 0 ) THEN
+                           jadd2 = 0
+                           DO jadd = 1, inpfiles(jj)%nadd
+                              IF ( TRIM(inpfiles(jj)%caddname(jadd)) == 'Hx' ) THEN
+                                 IF ( ldmod ) THEN
+                                    profdata%var(jvar)%vmod(ivart(jvar)) = &
+                                       &                inpfiles(jj)%padd(ij,ji,jadd,jvar)
+                                 ENDIF
+                              ELSE
+                                 jadd2 = jadd2 + 1
+                                 profdata%var(jvar)%vadd(ivart(jvar),kadd+jadd2) = &
+                                    &                inpfiles(jj)%padd(ij,ji,jadd,jvar)
+                              ENDIF
+                           END DO
                         ENDIF
 
@@ -726 +845 @@
                   
                   END DO
+                  
+                  ! Extra variables
+                  IF ( iextr > 0 ) THEN
+                     DO jext = 1, iextr
+                        profdata%vext(iprof,kextr+jext) = inpfiles(jj)%pext(ij,ji,jext)
+                     END DO
+                  ENDIF
 
                END DO loop_p
@@ -777 +903 @@
       ENDIF
 
-      IF (ldsatt) THEN
+      IF (ldallatall) THEN
          profdata%nvprot(:)    = ip3dt
          profdata%nvprotmpp(:) = ip3dtmpp
@@ -810 +936 @@
       ! Deallocate temporary data
       !-----------------------------------------------------------------------
-      DEALLOCATE( ifileidx, iprofidx, zdat, clvarsin )
+      DEALLOCATE( ifileidx, iprofidx, zdat, &
+         &        clvarsin, cllongin, clunitin, clgridin )
+      IF ( iadd > 0 ) THEN
+         DEALLOCATE( claddvarsin, claddlongin, claddunitin)
+      ENDIF
+      IF ( iextr > 0 ) THEN
+         DEALLOCATE( clextvarsin, clextlongin, clextunitin )
+      ENDIF
 
       !-----------------------------------------------------------------------

NEMO/branches/UKMO/NEMO_4.0.4_generic_obs/src/OCE/OBS/obs_read_surf.F90

@@ -39 +39 @@
 
    SUBROUTINE obs_rea_surf( surfdata, knumfiles, cdfilenames, &
-      &                     kvars, kextr, kstp, ddobsini, ddobsend, &
+      &                     kvars, kadd, kextr, kstp, ddobsini, ddobsend, &
       &                     ldignmis, ldmod, ldnightav, cdvars )
       !!---------------------------------------------------------------------
@@ -66 +66 @@
          & cdfilenames(knumfiles)       ! File names to read in
       INTEGER, INTENT(IN) :: kvars      ! Number of variables in surfdata
-      INTEGER, INTENT(IN) :: kextr      ! Number of extra fields for each var
+      INTEGER, INTENT(IN) :: kadd       ! Number of additional fields
+                                        !   in addition to those in the input file(s)
+      INTEGER, INTENT(IN) :: kextr      ! Number of extra fields
+                                        !   in addition to those in the input file(s)
       INTEGER, INTENT(IN) :: kstp       ! Ocean time-step index
       LOGICAL, INTENT(IN) :: ldignmis   ! Ignore missing files
@@ -78 +81 @@
       CHARACTER(LEN=11), PARAMETER :: cpname='obs_rea_surf'
       CHARACTER(len=8) :: clrefdate
-      CHARACTER(len=8), DIMENSION(:), ALLOCATABLE :: clvarsin
+      CHARACTER(len=ilenname), DIMENSION(:),   ALLOCATABLE :: clvarsin
+      CHARACTER(len=ilenlong), DIMENSION(:),   ALLOCATABLE :: cllongin
+      CHARACTER(len=ilenunit), DIMENSION(:),   ALLOCATABLE :: clunitin
+      CHARACTER(len=ilengrid), DIMENSION(:),   ALLOCATABLE :: clgridin
+      CHARACTER(len=ilenname), DIMENSION(:),   ALLOCATABLE :: claddvarsin
+      CHARACTER(len=ilenlong), DIMENSION(:,:), ALLOCATABLE :: claddlongin
+      CHARACTER(len=ilenunit), DIMENSION(:,:), ALLOCATABLE :: claddunitin
+      CHARACTER(len=ilenname), DIMENSION(:),   ALLOCATABLE :: clextvarsin
+      CHARACTER(len=ilenlong), DIMENSION(:),   ALLOCATABLE :: clextlongin
+      CHARACTER(len=ilenunit), DIMENSION(:),   ALLOCATABLE :: clextunitin
       INTEGER :: ji
       INTEGER :: jj
       INTEGER :: jk
+      INTEGER :: jvar
+      INTEGER :: jext
+      INTEGER :: jadd
+      INTEGER :: jadd2
+      INTEGER :: iadd
+      INTEGER :: iaddin
+      INTEGER :: iextr
       INTEGER :: iflag
       INTEGER :: inobf
@@ -121 +140 @@
       TYPE(obfbdata), POINTER, DIMENSION(:) :: &
          & inpfiles
+      CHARACTER(LEN=256) :: cout1  ! Diagnostic output line
 
       ! Local initialization
@@ -132 +152 @@
 
       ALLOCATE( inpfiles(inobf) )
+
+      iadd  = 0
+      iextr = 0
 
       surf_files : DO jj = 1, inobf
@@ -189 +212 @@
             ENDIF
 
+            IF ( (iextr > 0) .AND. (inpfiles(jj)%next /= iextr) ) THEN
+               CALL ctl_stop( 'Number of extra variables not consistent', &
+                  &           ' with previous files for this type' )
+            ELSE
+               iextr = inpfiles(jj)%next
+            ENDIF
+
+            ! Ignore model counterpart
+            iaddin = inpfiles(jj)%nadd
+            DO ji = 1, iaddin
+               IF ( TRIM(inpfiles(jj)%caddname(ji)) == 'Hx' ) THEN
+                  iaddin = iaddin - 1
+                  EXIT
+               ENDIF
+            END DO
+            IF ( ldmod .AND. ( inpfiles(jj)%nadd == iaddin ) ) THEN
+               CALL ctl_stop( 'Model not in input data' )
+            ENDIF
+
+            IF ( (iadd > 0) .AND. (iaddin /= iadd) ) THEN
+               CALL ctl_stop( 'Number of additional variables not consistent', &
+                  &           ' with previous files for this type' )
+            ELSE
+               iadd = iaddin
+            ENDIF
+
             IF ( jj == 1 ) THEN
                ALLOCATE( clvarsin( inpfiles(jj)%nvar ) )
+               ALLOCATE( cllongin( inpfiles(jj)%nvar ) )
+               ALLOCATE( clunitin( inpfiles(jj)%nvar ) )
+               ALLOCATE( clgridin( inpfiles(jj)%nvar ) )
                DO ji = 1, inpfiles(jj)%nvar
                  clvarsin(ji) = inpfiles(jj)%cname(ji)
+                 cllongin(ji) = inpfiles(jj)%coblong(ji)
+                 clunitin(ji) = inpfiles(jj)%cobunit(ji)
+                 clgridin(ji) = inpfiles(jj)%cgrid(ji)
                  IF ( clvarsin(ji) /= cdvars(ji) ) THEN
                     CALL ctl_stop( 'Feedback file variables do not match', &
@@ -198 +253 @@
                  ENDIF
                END DO
+               IF ( iadd > 0 ) THEN
+                  ALLOCATE( claddvarsin( iadd ) )
+                  ALLOCATE( claddlongin( iadd, inpfiles(jj)%nvar ) )
+                  ALLOCATE( claddunitin( iadd, inpfiles(jj)%nvar ) )
+                  jadd = 0
+                  DO ji = 1, inpfiles(jj)%nadd
+                    IF ( TRIM(inpfiles(jj)%caddname(ji)) /= 'Hx' ) THEN
+                       jadd = jadd + 1
+                       claddvarsin(jadd) = inpfiles(jj)%caddname(ji)
+                       DO jk = 1, inpfiles(jj)%nvar
+                          claddlongin(jadd,jk) = inpfiles(jj)%caddlong(ji,jk)
+                          claddunitin(jadd,jk) = inpfiles(jj)%caddunit(ji,jk)
+                       END DO
+                    ENDIF
+                  END DO
+               ENDIF
+               IF ( iextr > 0 ) THEN
+                  ALLOCATE( clextvarsin( iextr ) )
+                  ALLOCATE( clextlongin( iextr ) )
+                  ALLOCATE( clextunitin( iextr ) )
+                  DO ji = 1, iextr
+                    clextvarsin(ji) = inpfiles(jj)%cextname(ji)
+                    clextlongin(ji) = inpfiles(jj)%cextlong(ji)
+                    clextunitin(ji) = inpfiles(jj)%cextunit(ji)
+                  END DO
+               ENDIF
             ELSE
                DO ji = 1, inpfiles(jj)%nvar
@@ -205 +286 @@
                   ENDIF
                END DO
+               IF ( iadd > 0 ) THEN
+                  jadd = 0
+                  DO ji = 1, inpfiles(jj)%nadd
+                     IF ( TRIM(inpfiles(jj)%caddname(ji)) /= 'Hx' ) THEN
+                        jadd = jadd + 1
+                        IF ( inpfiles(jj)%caddname(ji) /= claddvarsin(jadd) ) THEN
+                           CALL ctl_stop( 'Feedback file additional variables not consistent', &
+                              &           ' with previous files for this type' )
+                        ENDIF
+                     ENDIF
+                  END DO
+               ENDIF
+               IF ( iextr > 0 ) THEN
+                  DO ji = 1, iextr
+                     IF ( inpfiles(jj)%cextname(ji) /= clextvarsin(ji) ) THEN
+                        CALL ctl_stop( 'Feedback file extra variables not consistent', &
+                           &           ' with previous files for this type' )
+                     ENDIF
+                  END DO
+               ENDIF
+
             ENDIF
 
@@ -351 +453 @@
          &               iindx   )
 
-      CALL obs_surf_alloc( surfdata, iobs, kvars, kextr, kstp, jpi, jpj )
+      CALL obs_surf_alloc( surfdata, iobs, kvars, kadd+iadd, kextr+iextr, kstp, jpi, jpj )
 
       ! Read obs/positions, QC, all variable and assign to surfdata
@@ -358 +460 @@
 
       surfdata%cvars(:)  = clvarsin(:)
+      surfdata%clong(:)  = cllongin(:)
+      surfdata%cunit(:)  = clunitin(:)
+      surfdata%cgrid(:)  = clgridin(:)
+      IF ( iadd > 0 ) THEN
+         surfdata%caddvars(kadd+1:)   = claddvarsin(:)
+         surfdata%caddlong(kadd+1:,:) = claddlongin(:,:)
+         surfdata%caddunit(kadd+1:,:) = claddunitin(:,:)
+      ENDIF
+      IF ( iextr > 0 ) THEN
+         surfdata%cextvars(kextr+1:) = clextvarsin(:)
+         surfdata%cextlong(kextr+1:) = clextlongin(:)
+         surfdata%cextunit(kextr+1:) = clextunitin(:)
+      ENDIF
 
       ityp   (:) = 0
@@ -433 +548 @@
                surfdata%nsfil(iobs) = iindx(jk)
 
-               ! QC flags
-               surfdata%nqc(iobs) = inpfiles(jj)%ivqc(ji,1)
-
-               ! Observed value
-               surfdata%robs(iobs,1) = inpfiles(jj)%pob(1,ji,1)
-
-
-               ! Model and MDT is set to fbrmdi unless read from file
-               IF ( ldmod ) THEN
-                  surfdata%rmod(iobs,1) = inpfiles(jj)%padd(1,ji,1,1)
-                  IF ( TRIM(surfdata%cvars(1)) == 'SLA' ) THEN
-                     surfdata%rext(iobs,1) = inpfiles(jj)%padd(1,ji,2,1)
-                     surfdata%rext(iobs,2) = inpfiles(jj)%pext(1,ji,1)
+               DO jvar = 1, kvars
+
+                  ! QC flags
+! WHY IS THIS NOT A FUNCTION OF NUM VARS?
+                  surfdata%nqc(iobs) = inpfiles(jj)%ivqc(ji,jvar)
+
+                  ! Observed value
+                  surfdata%robs(iobs,jvar) = inpfiles(jj)%pob(1,ji,jvar)
+
+! THIS NEEDS SORTING
+!                  ! Model and MDT is set to fbrmdi unless read from file
+!                  IF ( ldmod ) THEN
+!                     surfdata%rmod(iobs,jvar) = inpfiles(jj)%padd(1,ji,1,1)
+!                     IF ( TRIM(surfdata%cvars(1)) == 'SLA' ) THEN
+!                        surfdata%rext(iobs,1) = inpfiles(jj)%padd(1,ji,2,1)
+!                        surfdata%rext(iobs,2) = inpfiles(jj)%pext(1,ji,1)
+!                     ENDIF
+!                   ELSE
+!                     surfdata%rmod(iobs,jvar) = fbrmdi
+!                     IF ( TRIM(surfdata%cvars(1)) == 'SLA' ) surfdata%rext(iobs,:) = fbrmdi
+!                  ENDIF
+
+                  ! Additional variables
+                  surfdata%rmod(iobs,jvar) = fbrmdi
+                  IF ( iadd > 0 ) THEN
+                     jadd2 = 0
+                     DO jadd = 1, inpfiles(jj)%nadd
+                        IF ( TRIM(inpfiles(jj)%caddname(jadd)) == 'Hx' ) THEN
+                           IF ( ldmod ) THEN
+                              surfdata%rmod(iobs,jvar) = inpfiles(jj)%padd(1,ji,jadd,jvar)
+                           ENDIF
+                        ELSE
+                           jadd2 = jadd2 + 1
+                           surfdata%radd(iobs,kadd+jadd2,jvar) = &
+                              &                inpfiles(jj)%padd(1,ji,jadd,jvar)
+                        ENDIF
+                     END DO
                   ENDIF
-                ELSE
-                  surfdata%rmod(iobs,1) = fbrmdi
-                  IF ( TRIM(surfdata%cvars(1)) == 'SLA' ) surfdata%rext(iobs,:) = fbrmdi
+
+               END DO
+                  
+               ! Extra variables
+               IF ( iextr > 0 ) THEN
+                  DO jext = 1, iextr
+                     surfdata%rext(iobs,kextr+jext) = inpfiles(jj)%pext(1,ji,jext)
+                  END DO
                ENDIF
             ENDIF
@@ -467 +611 @@
       !-----------------------------------------------------------------------
       IF (lwp) THEN
-
+         DO jvar = 1, surfdata%nvar       
+            IF ( jvar == 1 ) THEN
+               cout1=TRIM(surfdata%cvars(1))                  
+            ELSE
+               WRITE(cout1,'(A,A1,A)') TRIM(cout1), '/', TRIM(surfdata%cvars(jvar))            
+            ENDIF
+         END DO
+ 
          WRITE(numout,*)
-         WRITE(numout,'(1X,A)')TRIM( surfdata%cvars(1) )//' data'
+         WRITE(numout,'(1X,A)')TRIM( cout1 )//' data'
          WRITE(numout,'(1X,A)')'--------------'
          DO jj = 1,8
@@ -479 +630 @@
             & '---------------------------------------------------------------'
          WRITE(numout,'(1X,A,I8)') &
-            & 'Total data for variable '//TRIM( surfdata%cvars(1) )// &
+            & 'Total data for variable '//TRIM( cout1 )// &
             & '           = ', iobsmpp
          WRITE(numout,'(1X,A)') &
@@ -490 +641 @@
       ! Deallocate temporary data
       !-----------------------------------------------------------------------
-      DEALLOCATE( ifileidx, isurfidx, zdat, clvarsin )
+      DEALLOCATE( ifileidx, isurfidx, zdat, clvarsin, &
+         &        cllongin, clunitin, clgridin )
+      IF ( iadd > 0 ) THEN
+         DEALLOCATE( claddvarsin, claddlongin, claddunitin)
+      ENDIF
+      IF ( iextr > 0 ) THEN
+         DEALLOCATE( clextvarsin, clextlongin, clextunitin )
+      ENDIF
 
       !-----------------------------------------------------------------------

NEMO/branches/UKMO/NEMO_4.0.4_generic_obs/src/OCE/OBS/obs_readmdt.F90

@@ -108 +108 @@
          & CALL obs_offset_mdt( jpi, jpj, z_mdt, zfill )
 
-      ! Intepolate the MDT already on the model grid at the observation point
-  
+      ! Interpolate the MDT already on the model grid at the observation point
+
       ALLOCATE( &
          & igrdi(2,2,sladata%nsurf), &
@@ -118 +118 @@
          & zmdtl(2,2,sladata%nsurf)  &
          & )
-         
+
       DO jobs = 1, sladata%nsurf
 
@@ -147 +147 @@
             
          CALL obs_int_h2d( 1, 1, zweig, zmdtl(:,:,jobs),  zext )
- 
-         sladata%rext(jobs,2) = zext(1)
+
+! FIGURE OUT THIS ASSIGNMENT 
+!         sladata%rext(jobs,2) = zext(1)
 
 ! mark any masked data with a QC flag
          IF( zobsmask(1) == 0 )   sladata%nqc(jobs) = IBSET(sladata%nqc(jobs),15)
 
-         END DO
-         
+      END DO
+
       DEALLOCATE( &
          & igrdi, &

NEMO/branches/UKMO/NEMO_4.0.4_generic_obs/src/OCE/OBS/obs_surf_def.F90

@@ -23 +23 @@
    USE obs_mpp, ONLY : &  ! MPP tools 
       obs_mpp_sum_integer
+   USE obs_fbm            ! Obs feedback format
 
    IMPLICIT NONE
@@ -45 +46 @@
       INTEGER :: nsurfmpp   !: Global number of surface data within window
       INTEGER :: nvar       !: Number of variables at observation points
+      INTEGER :: nadd       !: Number of additional fields at observation points
       INTEGER :: nextra     !: Number of extra fields at observation points
       INTEGER :: nstp       !: Number of time steps
@@ -69 +71 @@
          & ntyp           !: Type of surface observation product
 
-      CHARACTER(len=8), POINTER, DIMENSION(:) :: &
-         & cvars          !: Variable names
-
-      CHARACTER(LEN=8), POINTER, DIMENSION(:) :: &
+      CHARACTER(len=ilenname), POINTER, DIMENSION(:) :: &
+         & cvars,    &    !: Variable names
+         & cextvars, &    !: Extra variable names
+         & caddvars       !: Additional variable names
+
+      CHARACTER(len=ilenlong), POINTER, DIMENSION(:) :: &
+         & clong,    &    !: Variable long names
+         & cextlong       !: Extra variable long names
+
+      CHARACTER(len=ilenlong), POINTER, DIMENSION(:,:) :: &
+         & caddlong       !: Additional variable long names
+
+      CHARACTER(len=ilenunit), POINTER, DIMENSION(:) :: &
+         & cunit,    &    !: Variable units
+         & cextunit       !: Extra variable units
+
+      CHARACTER(len=ilenunit), POINTER, DIMENSION(:,:) :: &
+         & caddunit       !: Additional variable units
+
+      CHARACTER(len=ilengrid), POINTER, DIMENSION(:) :: &
+         & cgrid          !: Variable grids
+
+      CHARACTER(LEN=ilenwmo), POINTER, DIMENSION(:) :: &
          & cwmo           !: WMO indentifier
          
@@ -86 +107 @@
          & rext           !: Extra fields interpolated to observation points
 
-      REAL(KIND=wp), POINTER, DIMENSION(:,:) :: &
+      REAL(KIND=wp), POINTER, DIMENSION(:,:,:) :: &
+         & radd           !: Additional fields interpolated to observation points
+
+      REAL(KIND=wp), POINTER, DIMENSION(:,:,:) :: &
          & vdmean         !: Time averaged of model field
 
@@ -121 +145 @@
 CONTAINS
    
-   SUBROUTINE obs_surf_alloc( surf, ksurf, kvar, kextra, kstp, kpi, kpj )
+   SUBROUTINE obs_surf_alloc( surf, ksurf, kvar, kadd, kextra, kstp, kpi, kpj )
       !!----------------------------------------------------------------------
       !!                     ***  ROUTINE obs_surf_alloc  ***
@@ -136 +160 @@
       INTEGER, INTENT(IN) :: ksurf   ! Number of surface observations
       INTEGER, INTENT(IN) :: kvar    ! Number of surface variables
+      INTEGER, INTENT(IN) :: kadd    ! Number of additional fields at observation points
       INTEGER, INTENT(IN) :: kextra  ! Number of extra fields at observation points
       INTEGER, INTENT(IN) :: kstp    ! Number of time steps
@@ -143 +168 @@
       !!* Local variables
       INTEGER :: ji
-      INTEGER :: jvar
+      INTEGER :: jvar, jadd, jext
 
       ! Set bookkeeping variables
@@ -149 +174 @@
       surf%nsurf    = ksurf
       surf%nsurfmpp = 0
+      surf%nadd     = kadd
       surf%nextra   = kextra
       surf%nvar     = kvar
@@ -158 +184 @@
 
       ALLOCATE( &
-         & surf%cvars(kvar)    &
+         & surf%cvars(kvar), &
+         & surf%clong(kvar), &
+         & surf%cunit(kvar), &
+         & surf%cgrid(kvar)  &
          & )
 
       DO jvar = 1, kvar
          surf%cvars(jvar) = "NotSet"
+         surf%clong(jvar) = "NotSet"
+         surf%cunit(jvar) = "NotSet"
+         surf%cgrid(jvar) = ""
+      END DO
+
+      ! Allocate additional/extra variable metadata
+
+      ALLOCATE( &
+         & surf%caddvars(kadd),      &
+         & surf%caddlong(kadd,kvar), &
+         & surf%caddunit(kadd,kvar), &
+         & surf%cextvars(kextra),    &
+         & surf%cextlong(kextra),    &
+         & surf%cextunit(kextra)     &
+         )
+         
+      DO jadd = 1, kadd
+         surf%caddvars(jadd) = "NotSet"
+         DO jvar = 1, kvar
+            surf%caddlong(jadd,jvar) = "NotSet"
+            surf%caddunit(jadd,jvar) = "NotSet"
+         END DO
+      END DO
+         
+      DO jext = 1, kextra
+         surf%cextvars(jext) = "NotSet"
+         surf%cextlong(jext) = "NotSet"
+         surf%cextunit(jext) = "NotSet"
       END DO
       
@@ -205 +262 @@
       surf%rext(:,:) = 0.0_wp 
 
+      ! Allocate arrays of number of additional fields at observation points
+
+      ALLOCATE( & 
+         & surf%radd(ksurf,kadd,kvar) &
+         & )
+
+      surf%radd(:,:,:) = 0.0_wp 
+
       ! Allocate arrays of number of time step size
 
@@ -215 +280 @@
 
       ALLOCATE( &
-         & surf%vdmean(kpi,kpj) &
+         & surf%vdmean(kpi,kpj,kvar) &
          & )
 
@@ -291 +356 @@
          & )
 
+      ! Deallocate arrays of number of additional fields at observation points
+
+      DEALLOCATE( & 
+         & surf%radd &
+         & )
+
       ! Deallocate arrays of size number of grid points size times
       ! number of variables
@@ -308 +379 @@
 
       DEALLOCATE( &
-         & surf%cvars     &
-         & )
+         & surf%cvars, &
+         & surf%clong, &
+         & surf%cunit, &
+         & surf%cgrid  &
+         & )
+
+      ! Dellocate additional/extra variables metadata
+
+      DEALLOCATE( &
+         & surf%caddvars, &
+         & surf%caddlong, &
+         & surf%caddunit, &
+         & surf%cextvars, &
+         & surf%cextlong, &
+         & surf%cextunit  &
+         )
 
    END SUBROUTINE obs_surf_dealloc
@@ -343 +428 @@
       INTEGER :: ji
       INTEGER :: jk
+      INTEGER :: jadd
       LOGICAL, DIMENSION(:), ALLOCATABLE :: llvalid
 
@@ -361 +447 @@
 
       IF ( lallocate ) THEN
-         CALL obs_surf_alloc( newsurf,  insurf, surf%nvar, &
+         CALL obs_surf_alloc( newsurf,  insurf, surf%nvar, surf%nadd, &
             & surf%nextra, surf%nstp, surf%npi, surf%npj )
       ENDIF
@@ -410 +496 @@
                newsurf%rmod(insurf,jk)  = surf%rmod(ji,jk)
                
+               DO jadd = 1, surf%nadd
+                  newsurf%radd(insurf,jadd,jk) = surf%radd(ji,jadd,jk)
+               END DO
+               
             END DO
 
@@ -435 +525 @@
       newsurf%nstp     = surf%nstp
       newsurf%cvars(:) = surf%cvars(:)
+      newsurf%clong(:) = surf%clong(:)
+      newsurf%cunit(:) = surf%cunit(:)
+      newsurf%cgrid(:) = surf%cgrid(:)
+      newsurf%caddvars(:) = surf%caddvars(:)
+      newsurf%caddlong(:) = surf%caddlong(:)
+      newsurf%caddunit(:) = surf%caddunit(:)
+      newsurf%cextvars(:) = surf%cextvars(:)
+      newsurf%cextlong(:) = surf%cextlong(:)
+      newsurf%cextunit(:) = surf%cextunit(:)
       
       ! Set gridded stuff
@@ -470 +569 @@
       INTEGER :: jj
       INTEGER :: jk
+      INTEGER :: jadd
 
       ! Copy data from surf to old surf
@@ -504 +604 @@
             oldsurf%robs(jj,jk)  = surf%robs(ji,jk)
             oldsurf%rmod(jj,jk)  = surf%rmod(ji,jk)
+               
+            DO jadd = 1, surf%nadd
+               oldsurf%radd(jj,jadd,jk) = surf%radd(ji,jadd,jk)
+            END DO
 
          END DO

NEMO/branches/UKMO/NEMO_4.0.4_generic_obs/src/OCE/OBS/obs_write.F90

@@ -54 +54 @@
 CONTAINS
 
-   SUBROUTINE obs_wri_prof( profdata, padd, pext )
+   SUBROUTINE obs_wri_prof( profdata, clfiletype, padd, pext )
       !!-----------------------------------------------------------------------
       !!
@@ -77 +77 @@
       !! * Arguments
       TYPE(obs_prof), INTENT(INOUT) :: profdata      ! Full set of profile data
-      TYPE(obswriinfo), OPTIONAL :: padd             ! Additional info for each variable
-      TYPE(obswriinfo), OPTIONAL :: pext             ! Extra info
+      CHARACTER(LEN=25), INTENT(IN) :: clfiletype    ! Base name for file name
+      TYPE(obswriinfo), OPTIONAL    :: padd          ! Additional info for each variable
+      TYPE(obswriinfo), OPTIONAL    :: pext          ! Extra info
 
       !! * Local declarations
       TYPE(obfbdata) :: fbdata
       CHARACTER(LEN=40) :: clfname
-      CHARACTER(LEN=10) :: clfiletype
       CHARACTER(LEN=ilenlong) :: cllongname  ! Long name of variable
       CHARACTER(LEN=ilenunit) :: clunits     ! Units of variable
@@ -120 +120 @@
       END DO
 
-      SELECT CASE ( TRIM(profdata%cvars(1)) )
-      CASE('POTM')
-
-         clfiletype='profb'
-         CALL alloc_obfbdata( fbdata, 2, profdata%nprof, ilevel, &
-            &                 1 + iadd, 1 + iext, .TRUE. )
-         fbdata%cname(1)      = profdata%cvars(1)
-         fbdata%cname(2)      = profdata%cvars(2)
-         fbdata%coblong(1)    = 'Potential temperature'
-         fbdata%coblong(2)    = 'Practical salinity'
-         fbdata%cobunit(1)    = 'Degrees centigrade'
-         fbdata%cobunit(2)    = 'PSU'
-         fbdata%cextname(1)   = 'TEMP'
-         fbdata%cextlong(1)   = 'Insitu temperature'
-         fbdata%cextunit(1)   = 'Degrees centigrade'
-         fbdata%caddlong(1,1) = 'Model interpolated potential temperature'
-         fbdata%caddlong(1,2) = 'Model interpolated practical salinity'
-         fbdata%caddunit(1,1) = 'Degrees centigrade'
-         fbdata%caddunit(1,2) = 'PSU'
-         fbdata%cgrid(:)      = 'T'
-         DO je = 1, iext
-            fbdata%cextname(1+je) = pext%cdname(je)
-            fbdata%cextlong(1+je) = pext%cdlong(je,1)
-            fbdata%cextunit(1+je) = pext%cdunit(je,1)
-         END DO
-         DO ja = 1, iadd
-            fbdata%caddname(1+ja) = padd%cdname(ja)
-            DO jvar = 1, 2
+      CALL alloc_obfbdata( fbdata, profdata%nvar, profdata%nprof, ilevel, &
+            &                 1 + iadd, iext, .TRUE. )
+      fbdata%caddname(1)   = 'Hx'
+      DO jvar = 1, profdata%nvar
+         fbdata%cname(jvar)      = profdata%cvars(jvar)
+         fbdata%coblong(jvar)    = profdata%clong(jvar)
+         fbdata%cobunit(jvar)    = profdata%cunit(jvar)
+         fbdata%cgrid(jvar)      = profdata%cgrid(jvar)
+         fbdata%caddlong(1,jvar) = 'Model interpolated ' // TRIM(profdata%clong(jvar))
+         fbdata%caddunit(1,jvar) = profdata%cunit(jvar)
+         IF (iadd > 0) THEN
+            DO ja = 1, iadd
+               fbdata%caddname(1+ja) = padd%cdname(ja)
                fbdata%caddlong(1+ja,jvar) = padd%cdlong(ja,jvar)
                fbdata%caddunit(1+ja,jvar) = padd%cdunit(ja,jvar)
             END DO
-         END DO
-
-      CASE('UVEL')
-
-         clfiletype='velfb'
-         CALL alloc_obfbdata( fbdata, 2, profdata%nprof, ilevel, 1, 0, .TRUE. )
-         fbdata%cname(1)      = profdata%cvars(1)
-         fbdata%cname(2)      = profdata%cvars(2)
-         fbdata%coblong(1)    = 'Zonal velocity'
-         fbdata%coblong(2)    = 'Meridional velocity'
-         fbdata%cobunit(1)    = 'm/s'
-         fbdata%cobunit(2)    = 'm/s'
-         DO je = 1, iext
-            fbdata%cextname(je) = pext%cdname(je)
-            fbdata%cextlong(je) = pext%cdlong(je,1)
-            fbdata%cextunit(je) = pext%cdunit(je,1)
-         END DO
-         fbdata%caddlong(1,1) = 'Model interpolated zonal velocity'
-         fbdata%caddlong(1,2) = 'Model interpolated meridional velocity'
-         fbdata%caddunit(1,1) = 'm/s'
-         fbdata%caddunit(1,2) = 'm/s'
-         fbdata%cgrid(1)      = 'U' 
-         fbdata%cgrid(2)      = 'V'
-         DO ja = 1, iadd
-            fbdata%caddname(1+ja) = padd%cdname(ja)
-            fbdata%caddlong(1+ja,1) = padd%cdlong(ja,1)
-            fbdata%caddunit(1+ja,1) = padd%cdunit(ja,1)
-         END DO
-
-      END SELECT
-      
-      IF ( ( TRIM(profdata%cvars(1)) /= 'POTM' ) .AND. &
-         & ( TRIM(profdata%cvars(1)) /= 'UVEL' ) ) THEN
-         CALL alloc_obfbdata( fbdata, 1, profdata%nprof, ilevel, &
-            &                 1 + iadd, iext, .TRUE. )
-         fbdata%cname(1)      = profdata%cvars(1)
-         fbdata%coblong(1)    = cllongname
-         fbdata%cobunit(1)    = clunits
-         fbdata%caddlong(1,1) = 'Model interpolated ' // TRIM(cllongname)
-         fbdata%caddunit(1,1) = clunits
-         fbdata%cgrid(:)      = clgrid
-         DO je = 1, iext
-            fbdata%cextname(je) = pext%cdname(je)
-            fbdata%cextlong(je) = pext%cdlong(je,1)
-            fbdata%cextunit(je) = pext%cdunit(je,1)
-         END DO
-         DO ja = 1, iadd
-            fbdata%caddname(1+ja) = padd%cdname(ja)
-            fbdata%caddlong(1+ja,1) = padd%cdlong(ja,1)
-            fbdata%caddunit(1+ja,1) = padd%cdunit(ja,1)
-         END DO
-      ENDIF
-
-      fbdata%caddname(1)   = 'Hx'
-
-      WRITE(clfname, FMT="(A,'_fdbk_',I4.4,'.nc')") TRIM(clfiletype), nproc
+         ENDIF
+         IF (iext > 0) THEN
+            DO je = 1, iext
+               fbdata%cextname(je) = pext%cdname(je)
+               fbdata%cextlong(je) = pext%cdlong(je,1)
+               fbdata%cextunit(je) = pext%cdunit(je,1)
+            END DO
+         ENDIF
+      END DO
+!fbdata%cextname(1)   = 'TEMP'
+!fbdata%cextlong(1)   = 'Insitu temperature'
+!fbdata%cextunit(1)   = 'Degrees centigrade'
+
+      WRITE(clfname, FMT="(A,'fb_fdbk_',I4.4,'.nc')") TRIM(clfiletype), nproc
 
       IF(lwp) THEN
@@ -274 +217 @@
                ENDIF
                fbdata%iobsk(ik,jo,jvar)  = profdata%var(jvar)%mvk(jk)
-               DO ja = 1, iadd
-                  fbdata%padd(ik,jo,1+ja,jvar) = &
-                     & profdata%var(jvar)%vext(jk,padd%ipoint(ja))
-               END DO
-               DO je = 1, iext
-                  fbdata%pext(ik,jo,1+je) = &
-                     & profdata%var(jvar)%vext(jk,pext%ipoint(je))
-               END DO
-               IF ( ( jvar == 1 ) .AND. &
-                  & ( TRIM(profdata%cvars(1)) == 'POTM' ) ) THEN
-                  fbdata%pext(ik,jo,1) = profdata%var(jvar)%vext(jk,1)
-               ENDIF 
+               IF (iadd > 0) THEN
+                  DO ja = 1, iadd
+                     fbdata%padd(ik,jo,1+ja,jvar) = &
+                        & profdata%var(jvar)%vadd(jk,padd%ipoint(ja))
+                  END DO
+               ENDIF
+! MOVE OUTSIDE JVAR LOOP?
+               IF (iext > 0) THEN
+                  DO je = 1, iext
+                     fbdata%pext(ik,jo,je) = &
+                        & profdata%vext(jk,pext%ipoint(je))
+                  END DO
+               ENDIF
+!IF ( ( jvar == 1 ) .AND. &
+!   & ( TRIM(profdata%cvars(1)) == 'POTM' ) ) THEN
+!   fbdata%pext(ik,jo,1) = profdata%var(jvar)%vext(jk,1)
+!ENDIF 
             END DO
          END DO
       END DO
 
-      IF ( TRIM(profdata%cvars(1)) == 'POTM' ) THEN
-         ! Convert insitu temperature to potential temperature using the model
-         ! salinity if no potential temperature
-         DO jo = 1, fbdata%nobs
-            IF ( fbdata%pphi(jo) < 9999.0 ) THEN
-               DO jk = 1, fbdata%nlev
-                  IF ( ( fbdata%pob(jk,jo,1) >= 9999.0 ) .AND. &
-                     & ( fbdata%pdep(jk,jo) < 9999.0 ) .AND. &
-                     & ( fbdata%padd(jk,jo,1,2) < 9999.0 ) .AND. &
-                     & ( fbdata%pext(jk,jo,1) < 9999.0 ) ) THEN
-                     zpres = dep_to_p( REAL(fbdata%pdep(jk,jo),wp), &
-                        &              REAL(fbdata%pphi(jo),wp) )
-                     fbdata%pob(jk,jo,1) = potemp( &
-                        &                     REAL(fbdata%padd(jk,jo,1,2), wp),  &
-                        &                     REAL(fbdata%pext(jk,jo,1), wp), &
-                        &                     zpres, 0.0_wp )
-                  ENDIF
-               END DO
-            ENDIF
-         END DO
-      ENDIF
+!IF ( TRIM(profdata%cvars(1)) == 'POTM' ) THEN
+!   ! Convert insitu temperature to potential temperature using the model
+!   ! salinity if no potential temperature
+!   DO jo = 1, fbdata%nobs
+!      IF ( fbdata%pphi(jo) < 9999.0 ) THEN
+!         DO jk = 1, fbdata%nlev
+!            IF ( ( fbdata%pob(jk,jo,1) >= 9999.0 ) .AND. &
+!               & ( fbdata%pdep(jk,jo) < 9999.0 ) .AND. &
+!               & ( fbdata%padd(jk,jo,1,2) < 9999.0 ) .AND. &
+!               & ( fbdata%pext(jk,jo,1) < 9999.0 ) ) THEN
+!               zpres = dep_to_p( REAL(fbdata%pdep(jk,jo),wp), &
+!                  &              REAL(fbdata%pphi(jo),wp) )
+!               fbdata%pob(jk,jo,1) = potemp( &
+!                  &                     REAL(fbdata%padd(jk,jo,1,2), wp),  &
+!                  &                     REAL(fbdata%pext(jk,jo,1), wp), &
+!                  &                     zpres, 0.0_wp )
+!            ENDIF
+!         END DO
+!      ENDIF
+!   END DO
+!ENDIF
 
       ! Write the obfbdata structure
@@ -322 +270 @@
    END SUBROUTINE obs_wri_prof
 
-   SUBROUTINE obs_wri_surf( surfdata, padd, pext )
+   SUBROUTINE obs_wri_surf( surfdata, clfiletype, padd, pext )
       !!-----------------------------------------------------------------------
       !!
@@ -342 +290 @@
 
       !! * Arguments
-      TYPE(obs_surf), INTENT(INOUT) :: surfdata         ! Full set of surface data
-      TYPE(obswriinfo), OPTIONAL :: padd               ! Additional info for each variable
-      TYPE(obswriinfo), OPTIONAL :: pext               ! Extra info
+      TYPE(obs_surf), INTENT(INOUT) :: surfdata      ! Full set of surface data
+      CHARACTER(LEN=25), INTENT(IN) :: clfiletype    ! Base name for file name
+      TYPE(obswriinfo), OPTIONAL    :: padd          ! Additional info for each variable
+      TYPE(obswriinfo), OPTIONAL    :: pext          ! Extra info
 
       !! * Local declarations
       TYPE(obfbdata) :: fbdata
       CHARACTER(LEN=40) :: clfname         ! netCDF filename
-      CHARACTER(LEN=10) :: clfiletype
       CHARACTER(LEN=ilenlong) :: cllongname  ! Long name of variable
       CHARACTER(LEN=ilenunit) :: clunits     ! Units of variable
@@ -357 +305 @@
       INTEGER :: ja
       INTEGER :: je
+      INTEGER :: jvar
       INTEGER :: iadd
       INTEGER :: iext
@@ -374 +323 @@
       CALL init_obfbdata( fbdata )
 
-      SELECT CASE ( TRIM(surfdata%cvars(1)) )
-      CASE('SLA')
-         
-         ! SLA needs special treatment because of MDT, so is all done here
-         ! Other variables are done more generically
-         ! No climatology for SLA, MDT is our best estimate of that and is already output.
-
-         CALL alloc_obfbdata( fbdata, 1, surfdata%nsurf, 1, &
-            &                 2 + iadd, 1 + iext, .TRUE. )
-
-         clfiletype = 'slafb'
-         fbdata%cname(1)      = surfdata%cvars(1)
-         fbdata%coblong(1)    = 'Sea level anomaly'
-         fbdata%cobunit(1)    = 'Metres'
-         fbdata%cextname(1)   = 'MDT'
-         fbdata%cextlong(1)   = 'Mean dynamic topography'
-         fbdata%cextunit(1)   = 'Metres'
-         DO je = 1, iext
-            fbdata%cextname(je) = pext%cdname(je)
-            fbdata%cextlong(je) = pext%cdlong(je,1)
-            fbdata%cextunit(je) = pext%cdunit(je,1)
-         END DO
-         fbdata%caddlong(1,1) = 'Model interpolated SSH - MDT'
-         fbdata%caddunit(1,1) = 'Metres' 
-         fbdata%caddname(2)   = 'SSH'
-         fbdata%caddlong(2,1) = 'Model Sea surface height'
-         fbdata%caddunit(2,1) = 'Metres'
-         fbdata%cgrid(1)      = 'T'
-         DO ja = 1, iadd
-            fbdata%caddname(2+ja) = padd%cdname(ja)
-            fbdata%caddlong(2+ja,1) = padd%cdlong(ja,1)
-            fbdata%caddunit(2+ja,1) = padd%cdunit(ja,1)
-         END DO
-
-      CASE('SST')
-
-         clfiletype = 'sstfb'
-         cllongname = 'Sea surface temperature'
-         clunits    = 'Degree centigrade'
-         clgrid     = 'T'
-         
-      CASE('ICECONC')
-
-         clfiletype = 'sicfb'
-         cllongname = 'Sea ice concentration'
-         clunits    = 'Fraction'
-         clgrid     = 'T'
-
-      CASE('SSS')
-
-         clfiletype = 'sssfb'
-         cllongname = 'Sea surface salinity'
-         clunits    = 'psu'
-         clgrid     = 'T'
-
-      CASE DEFAULT
-
-         CALL ctl_stop( 'Unknown observation type '//TRIM(surfdata%cvars(1))//' in obs_wri_surf' )
-
-      END SELECT
-
-      ! SLA needs special treatment because of MDT, so is done above
-      ! Remaining variables treated more generically
-
-      IF ( TRIM(surfdata%cvars(1)) /= 'SLA' ) THEN
-      
-         CALL alloc_obfbdata( fbdata, 1, surfdata%nsurf, 1, &
+      CALL alloc_obfbdata( fbdata, surfdata%nvar, surfdata%nsurf, 1, &
             &                 1 + iadd, iext, .TRUE. )
-
-         fbdata%cname(1)      = surfdata%cvars(1)
-         fbdata%coblong(1)    = cllongname
-         fbdata%cobunit(1)    = clunits
-         DO je = 1, iext
-            fbdata%cextname(je) = pext%cdname(je)
-            fbdata%cextlong(je) = pext%cdlong(je,1)
-            fbdata%cextunit(je) = pext%cdunit(je,1)
-         END DO        
-         IF ( TRIM(surfdata%cvars(1)) == 'ICECONC' ) THEN
-            fbdata%caddlong(1,1) = 'Model interpolated ICE'
-         ELSE
-            fbdata%caddlong(1,1) = 'Model interpolated ' // TRIM(surfdata%cvars(1))
-         ENDIF
-         fbdata%caddunit(1,1) = clunits
-         fbdata%cgrid(1)      = clgrid
-         DO ja = 1, iadd
-            fbdata%caddname(1+ja) = padd%cdname(ja)
-            fbdata%caddlong(1+ja,1) = padd%cdlong(ja,1)
-            fbdata%caddunit(1+ja,1) = padd%cdunit(ja,1)
-         END DO
-      ENDIF
-
       fbdata%caddname(1)   = 'Hx'
-
-      WRITE(clfname, FMT="(A,'_fdbk_',I4.4,'.nc')") TRIM(clfiletype), nproc
+      DO jvar = 1, surfdata%nvar
+         fbdata%cname(jvar)      = surfdata%cvars(jvar)
+         fbdata%coblong(jvar)    = surfdata%clong(jvar)
+         fbdata%cobunit(jvar)    = surfdata%cunit(jvar)
+         fbdata%cgrid(jvar)      = surfdata%cgrid(jvar)
+         fbdata%caddlong(1,jvar) = 'Model interpolated ' // TRIM(surfdata%clong(jvar))
+         fbdata%caddunit(1,jvar) = surfdata%cunit(jvar)
+         IF (iadd > 0) THEN
+            DO ja = 1, iadd
+               fbdata%caddname(1+ja) = padd%cdname(ja)
+               fbdata%caddlong(1+ja,jvar) = padd%cdlong(ja,jvar)
+               fbdata%caddunit(1+ja,jvar) = padd%cdunit(ja,jvar)
+            END DO
+         ENDIF
+         IF (iext > 0) THEN
+            DO je = 1, iext
+               fbdata%cextname(je) = pext%cdname(je)
+               fbdata%cextlong(je) = pext%cdlong(je,1)
+               fbdata%cextunit(je) = pext%cdunit(je,1)
+            END DO
+         ENDIF
+      END DO
+!fbdata%cname(1)      = surfdata%cvars(1)
+!fbdata%coblong(1)    = 'Sea level anomaly'
+!fbdata%cobunit(1)    = 'Metres'
+!fbdata%cextname(1)   = 'MDT'
+!fbdata%cextlong(1)   = 'Mean dynamic topography'
+!fbdata%cextunit(1)   = 'Metres'
+!fbdata%caddname(2)   = 'SSH'
+!fbdata%caddlong(2,1) = 'Model Sea surface height'
+!fbdata%caddunit(2,1) = 'Metres'
+!IF ( TRIM(surfdata%cvars(1)) == 'ICECONC' ) THEN
+!   fbdata%caddlong(1,1) = 'Model interpolated ICE'
+!ELSE
+!   fbdata%caddlong(1,1) = 'Model interpolated ' // TRIM(surfdata%cvars(1))
+!ENDIF
+
+      WRITE(clfname, FMT="(A,'fb_fdbk_',I4.4,'.nc')") TRIM(clfiletype), nproc
 
       IF(lwp) THEN
@@ -514 +410 @@
             &           surfdata%nyea(jo), &
             &           fbdata%ptim(jo),   &
-            &           krefdate = 19500101 )
-         fbdata%padd(1,jo,1,1) = surfdata%rmod(jo,1)
-         IF ( TRIM(surfdata%cvars(1)) == 'SLA' ) fbdata%padd(1,jo,2,1) = surfdata%rext(jo,1)
-         fbdata%pob(1,jo,1)    = surfdata%robs(jo,1) 
+            &           krefdate = 19500101 ) 
          fbdata%pdep(1,jo)     = 0.0
          fbdata%idqc(1,jo)     = 0
          fbdata%idqcf(:,1,jo)  = 0
-         IF ( surfdata%nqc(jo) > 255 ) THEN
-            fbdata%ivqc(jo,1)       = 4
-            fbdata%ivlqc(1,jo,1)    = 4
-            fbdata%ivlqcf(1,1,jo,1) = 0
+         DO jvar = 1, surfdata%nvar
+            fbdata%padd(1,jo,1,jvar) = surfdata%rmod(jo,jvar)
+!IF ( TRIM(surfdata%cvars(1)) == 'SLA' ) fbdata%padd(1,jo,2,1) = surfdata%rext(jo,1)
+            fbdata%pob(1,jo,jvar)    = surfdata%robs(jo,jvar)
+            IF ( surfdata%nqc(jo) > 255 ) THEN
+               fbdata%ivqc(jo,jvar)       = 4
+               fbdata%ivlqc(1,jo,jvar)    = 4
+               fbdata%ivlqcf(1,1,jo,jvar) = 0
 !$AGRIF_DO_NOT_TREAT
-            fbdata%ivlqcf(2,1,jo,1) = IAND(surfdata%nqc(jo),b'0000000011111111')
+               fbdata%ivlqcf(2,1,jo,jvar) = IAND(surfdata%nqc(jo),b'0000000011111111')
 !$AGRIF_END_DO_NOT_TREAT
-         ELSE
-            fbdata%ivqc(jo,1)       = surfdata%nqc(jo)
-            fbdata%ivlqc(1,jo,1)    = surfdata%nqc(jo)
-            fbdata%ivlqcf(:,1,jo,1) = 0
-         ENDIF
-         fbdata%iobsk(1,jo,1)  = 0
-         IF ( TRIM(surfdata%cvars(1)) == 'SLA' ) fbdata%pext(1,jo,1) = surfdata%rext(jo,2)
-         DO ja = 1, iadd
-            fbdata%padd(1,jo,2+ja,1) = &
-               & surfdata%rext(jo,padd%ipoint(ja))
+            ELSE
+               fbdata%ivqc(jo,jvar)       = surfdata%nqc(jo)
+               fbdata%ivlqc(1,jo,jvar)    = surfdata%nqc(jo)
+               fbdata%ivlqcf(:,1,jo,jvar) = 0
+            ENDIF
+            fbdata%iobsk(1,jo,jvar)  = 0
+!IF ( TRIM(surfdata%cvars(1)) == 'SLA' ) fbdata%pext(1,jo,1) = surfdata%rext(jo,2)
+            IF (iadd > 0) THEN
+               DO ja = 1, iadd
+                  fbdata%padd(1,jo,1+ja,jvar) = &
+                     & surfdata%radd(jo,padd%ipoint(ja),jvar)
+               END DO
+            ENDIF
          END DO
-         DO je = 1, iext
-            fbdata%pext(1,jo,1+je) = &
-               & surfdata%rext(jo,pext%ipoint(je))
-         END DO
+         IF (iext > 0) THEN
+            DO je = 1, iext
+               fbdata%pext(1,jo,je) = &
+                  & surfdata%rext(jo,pext%ipoint(je))
+            END DO
+         ENDIF
       END DO