Index: branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/nemo.f90
===================================================================
--- branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/nemo.f90	(revision 4102)
+++ branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/nemo.f90	(revision 4106)
@@ -1,4 +1,4 @@
 PROGRAM nemo
-   USE off_oper
+   USE ooomod
    IMPLICIT NONE
    !!----------------------------------------------------------------------
@@ -8,5 +8,5 @@
    !!              location.
    !!----------------------------------------------------------------------
-   CALL off_obs_oper
+   CALL offline_obs_oper
 END PROGRAM nemo
 
Index: branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/nemogcm.F90
===================================================================
--- branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/nemogcm.F90	(revision 4102)
+++ branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/nemogcm.F90	(revision 4106)
@@ -101,112 +101,4 @@
 CONTAINS
 
-   SUBROUTINE nemo_gcm
-      !!----------------------------------------------------------------------
-      !!                     ***  ROUTINE nemo_gcm  ***
-      !!
-      !! ** Purpose :   NEMO solves the primitive equations on an orthogonal
-      !!              curvilinear mesh on the sphere.
-      !!
-      !! ** Method  : - model general initialization
-      !!              - launch the time-stepping (stp routine)
-      !!              - finalize the run by closing files and communications
-      !!
-      !! References : Madec, Delecluse, Imbard, and Levy, 1997:  internal report, IPSL.
-      !!              Madec, 2008, internal report, IPSL.
-      !!----------------------------------------------------------------------
-      INTEGER ::   istp       ! time step index
-      !!----------------------------------------------------------------------
-      !
-#if defined key_agrif
-      CALL Agrif_Init_Grids()      ! AGRIF: set the meshes
-#endif
-
-      !                            !-----------------------!
-      CALL nemo_init               !==  Initialisations  ==!
-      !                            !-----------------------!
-#if defined key_agrif
-      CALL Agrif_Declare_Var_dom   ! AGRIF: set the meshes for DOM
-      CALL Agrif_Declare_Var       !  "      "   "   "      "  DYN/TRA 
-# if defined key_top
-      CALL Agrif_Declare_Var_top   !  "      "   "   "      "  TOP
-# endif
-# if defined key_lim2
-      CALL Agrif_Declare_Var_lim2  !  "      "   "   "      "  LIM
-# endif
-#endif
-      ! check that all process are still there... If some process have an error,
-      ! they will never enter in step and other processes will wait until the end of the cpu time!
-      IF( lk_mpp )   CALL mpp_max( nstop )
-
-      IF(lwp) WRITE(numout,cform_aaa)   ! Flag AAAAAAA
-
-      !                            !-----------------------!
-      !                            !==   time stepping   ==!
-      !                            !-----------------------!
-      istp = nit000
-#if defined key_c1d
-         DO WHILE ( istp <= nitend .AND. nstop == 0 )
-            CALL stp_c1d( istp )
-            istp = istp + 1
-         END DO
-#else
-          IF( lk_asminc ) THEN
-             IF( ln_bkgwri ) CALL asm_bkg_wri( nit000 - 1 )    ! Output background fields
-             IF( ln_asmdin ) THEN                        ! Direct initialization
-                IF( ln_trainc ) CALL tra_asm_inc( nit000 - 1 )    ! Tracers
-                IF( ln_dyninc ) CALL dyn_asm_inc( nit000 - 1 )    ! Dynamics
-                IF( ln_sshinc ) CALL ssh_asm_inc( nit000 - 1 )    ! SSH
-             ENDIF
-          ENDIF
-
-         DO WHILE ( istp <= nitend .AND. nstop == 0 )
-#if defined key_agrif
-            CALL Agrif_Step( stp )           ! AGRIF: time stepping
-#else
-            CALL stp( istp )                 ! standard time stepping
-#endif
-            istp = istp + 1
-            IF( lk_mpp )   CALL mpp_max( nstop )
-         END DO
-#endif
-
-      IF( lk_diaobs   )   CALL dia_obs_wri
-      !
-      IF( ln_icebergs )   CALL icb_end( nitend )
-
-      !                            !------------------------!
-      !                            !==  finalize the run  ==!
-      !                            !------------------------!
-      IF(lwp) WRITE(numout,cform_aaa)   ! Flag AAAAAAA
-      !
-      IF( nstop /= 0 .AND. lwp ) THEN   ! error print
-         WRITE(numout,cform_err)
-         WRITE(numout,*) nstop, ' error have been found'
-      ENDIF
-      !
-#if defined key_agrif
-      CALL Agrif_ParentGrid_To_ChildGrid()
-      IF( lk_diaobs ) CALL dia_obs_wri
-      IF( nn_timing == 1 )   CALL timing_finalize
-      CALL Agrif_ChildGrid_To_ParentGrid()
-#endif
-      IF( nn_timing == 1 )   CALL timing_finalize
-      !
-      CALL nemo_closefile
-#if defined key_iomput
-      CALL xios_finalize                ! end mpp communications with xios
-# if defined key_oasis3 || defined key_oasis4
-      CALL cpl_prism_finalize           ! end coupling and mpp communications with OASIS
-# endif
-#else
-# if defined key_oasis3 || defined key_oasis4
-      CALL cpl_prism_finalize           ! end coupling and mpp communications with OASIS
-# else
-      IF( lk_mpp )   CALL mppstop       ! end mpp communications
-# endif
-#endif
-      !
-   END SUBROUTINE nemo_gcm
-
 
    SUBROUTINE nemo_init
@@ -335,73 +227,10 @@
       IF( ln_ctl        )   CALL prt_ctl_init   ! Print control
 
-      IF( lk_obc        )   CALL     obc_init   ! Open boundaries
-
                             CALL  istate_init   ! ocean initial state (Dynamics and tracers)
 
-      IF( lk_tide       )   CALL tide_init( nit000 )    ! Initialisation of the tidal harmonics
-
-      IF( lk_bdy        )   CALL     bdy_init       ! Open boundaries initialisation
-      IF( lk_bdy        )   CALL     bdy_dta_init   ! Open boundaries initialisation of external data arrays
-      IF( lk_bdy        )   CALL     bdytide_init   ! Open boundaries initialisation of tidal harmonic forcing
-
-                            CALL dyn_nept_init  ! simplified form of Neptune effect
-
-      !                                     ! Ocean physics
-                            CALL     sbc_init   ! Forcings : surface module
-      !                                         ! Vertical physics
-                            CALL     zdf_init      ! namelist read
-                            CALL zdf_bfr_init      ! bottom friction
-      IF( lk_zdfric     )   CALL zdf_ric_init      ! Richardson number dependent Kz
-      IF( lk_zdftke     )   CALL zdf_tke_init      ! TKE closure scheme
-      IF( lk_zdfgls     )   CALL zdf_gls_init      ! GLS closure scheme
-      IF( lk_zdfkpp     )   CALL zdf_kpp_init      ! KPP closure scheme
-      IF( lk_zdftmx     )   CALL zdf_tmx_init      ! tidal vertical mixing
-      IF( lk_zdfddm .AND. .NOT. lk_zdfkpp )   &
-         &                  CALL zdf_ddm_init      ! double diffusive mixing
-      !                                         ! Lateral physics
-                            CALL ldf_tra_init      ! Lateral ocean tracer physics
-                            CALL ldf_dyn_init      ! Lateral ocean momentum physics
-      IF( lk_ldfslp     )   CALL ldf_slp_init      ! slope of lateral mixing
-
-      !                                     ! Active tracers
-                            CALL tra_qsr_init   ! penetrative solar radiation qsr
-                            CALL tra_bbc_init   ! bottom heat flux
-      IF( lk_trabbl     )   CALL tra_bbl_init   ! advective (and/or diffusive) bottom boundary layer scheme
-      IF( ln_tradmp     )   CALL tra_dmp_init   ! internal damping trends
-                            CALL tra_adv_init   ! horizontal & vertical advection
-                            CALL tra_ldf_init   ! lateral mixing
-                            CALL tra_zdf_init   ! vertical mixing and after tracer fields
-
-      !                                     ! Dynamics
-                            CALL dyn_adv_init   ! advection (vector or flux form)
-                            CALL dyn_vor_init   ! vorticity term including Coriolis
-                            CALL dyn_ldf_init   ! lateral mixing
-                            CALL dyn_hpg_init   ! horizontal gradient of Hydrostatic pressure
-                            CALL dyn_zdf_init   ! vertical diffusion
-                            CALL dyn_spg_init   ! surface pressure gradient
-
-      !                                     ! Misc. options
-      IF( nn_cla == 1   )   CALL cla_init       ! Cross Land Advection
-                            CALL icb_init( rdt, nit000)   ! initialise icebergs instance
-      
-#if defined key_top
-      !                                     ! Passive tracers
-                            CALL     trc_init
-#endif
-      !                                     ! Diagnostics
-      IF( lk_floats     )   CALL     flo_init   ! drifting Floats
-      IF( lk_diaar5     )   CALL dia_ar5_init   ! ar5 diag
-                            CALL dia_ptr_init   ! Poleward TRansports initialization
-      IF( lk_diadct     )   CALL dia_dct_init   ! Sections tranports
-                            CALL dia_hsb_init   ! heat content, salt content and volume budgets
-                            CALL trd_mod_init   ! Mixed-layer/Vorticity/Integral constraints trends
       IF( lk_diaobs     ) THEN                  ! Observation & model comparison
                             CALL dia_obs_init            ! Initialize observational data
                             CALL dia_obs( nit000 - 1 )   ! Observation operator for restart
       ENDIF
-      !                                     ! Assimilation increments
-      IF( lk_asminc     )   CALL asm_inc_init   ! Initialize assimilation increments
-      IF(lwp) WRITE(numout,*) 'Euler time step switch is ', neuler
-      !
    END SUBROUTINE nemo_init
 
Index: branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/obs_fbm.F90
===================================================================
--- branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/obs_fbm.F90	(revision 4106)
+++ branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/obs_fbm.F90	(revision 4106)
@@ -0,0 +1,2263 @@
+MODULE obs_fbm
+   !!======================================================================
+   !!                       ***  MODULE obs_fbm  ***
+   !! Observation operators : I/O + tools for feedback files
+   !!======================================================================
+   !! History : 
+   !!             !  08-11  (K. Mogensen) Initial version
+   !!----------------------------------------------------------------------
+
+   !!----------------------------------------------------------------------
+   !!   init_obfbdata     :  Initialize sizes in obfbdata structure
+   !!   alloc_obfbdata    :  Allocate data in an obfbdata structure
+   !!   dealloc_obfbdata  :  Dellocate data in an obfbdata structure
+   !!   copy_obfbdata     :  Copy an obfbdata structure
+   !!   subsamp_obfbdata  :  Sumsample an obfbdata structure
+   !!   merge_obfbdata    :  Merge multiple obfbdata structures into an one.
+   !!   write_obfbdata    :  Write an obfbdata structure into a netCDF file.
+   !!   read_obfbdata     :  Read an obfbdata structure from a netCDF file.
+   !!----------------------------------------------------------------------
+   USE netcdf
+   USE obs_utils      ! Various utilities for observation operators
+
+   IMPLICIT NONE
+   PUBLIC
+
+   ! Type kinds for feedback data.
+
+   INTEGER, PARAMETER :: fbsp = SELECTED_REAL_KIND( 6, 37) !: single precision
+   INTEGER, PARAMETER :: fbdp = SELECTED_REAL_KIND(12,307) !: double precision
+
+   ! Parameters for string lengths.
+
+   INTEGER, PARAMETER    :: ilenwmo  = 8    !: Length of station identifier
+   INTEGER, PARAMETER    :: ilentyp  = 4    !: Length of type
+   INTEGER, PARAMETER    :: ilenname = 8    !: Length of variable names
+   INTEGER, PARAMETER    :: ilengrid = 1    !: Grid (e.g. 'T') length
+   INTEGER, PARAMETER    :: ilenjuld = 14   !: Lenght of reference julian date
+   INTEGER, PARAMETER    :: idefnqcf = 2    !: Default number of words in QC
+                                            !  flags
+   INTEGER, PARAMETER    :: ilenlong = 128  !: Length of long name
+   INTEGER, PARAMETER    :: ilenunit = 32   !: Length of units
+   
+   ! Missinge data indicators 
+   
+   INTEGER, PARAMETER    :: fbimdi = -99999   !: Integers
+   REAL(fbsp), PARAMETER :: fbrmdi =  99999   !: Reals
+
+   ! Output stream choice
+   LOGICAL               :: ln_cl4 = .FALSE.  !: Logical switch for
+                                              !: class 4 file outputs
+ 
+   ! Main data structure for observation feedback data.
+
+   TYPE obfbdata
+      LOGICAL :: lalloc         !: Allocation status for data
+      LOGICAL :: lgrid          !: Include grid search info
+      INTEGER :: nvar           !: Number of variables
+      INTEGER :: nobs           !: Number of observations
+      INTEGER :: nlev           !: Number of levels
+      INTEGER :: nadd           !: Number of additional entries
+      INTEGER :: next           !: Number of extra variables
+      INTEGER :: nqcf           !: Number of words per qc flag
+      CHARACTER(LEN=ilenwmo), DIMENSION(:), POINTER :: &
+         & cdwmo                !: Identifier
+      CHARACTER(LEN=ilentyp), DIMENSION(:), POINTER :: &
+         & cdtyp                !: Instrument type
+      CHARACTER(LEN=ilenjuld) :: &
+         & cdjuldref            !: Julian date reference
+      INTEGER, DIMENSION(:), POINTER :: &
+         & kindex               !: Index of observations in the original file
+      INTEGER, DIMENSION(:), POINTER :: &
+         & ioqc, &              !: Observation QC
+         & ipqc, &              !: Position QC
+         & itqc                 !: Time QC
+      INTEGER, DIMENSION(:,:), POINTER :: &
+         & ioqcf, &             !: Observation QC flags
+         & ipqcf, &             !: Position QC flags
+         & itqcf                !: Time QC flags
+      INTEGER, DIMENSION(:,:), POINTER :: &
+         & idqc                 !: Depth QC
+      INTEGER, DIMENSION(:,:,:), POINTER :: &
+         & idqcf                !: Depth QC flags
+      REAL(KIND=fbdp), DIMENSION(:), POINTER :: &
+         & plam, &              !: Longitude
+         & pphi, &              !: Latitude
+         & ptim                 !: Time
+      REAL(KIND=fbsp), DIMENSION(:,:), POINTER :: &
+         & pdep                 !: Depth
+      CHARACTER(LEN=ilenname), DIMENSION(:), POINTER  :: &
+         & cname                !: Name of variable
+      REAL(fbsp), DIMENSION(:,:,:), POINTER :: &
+         & pob                  !: Observation
+      CHARACTER(LEN=ilenlong), DIMENSION(:), POINTER :: &
+         & coblong              !: Observation long name (for output)
+      CHARACTER(LEN=ilenunit), DIMENSION(:), POINTER :: &
+         & cobunit              !: Observation units (for output)
+      INTEGER, DIMENSION(:,:), POINTER :: &
+         & ivqc                 !: Variable QC
+      INTEGER, DIMENSION(:,:,:), POINTER :: &
+         & ivqcf                !: Variable QC flags
+      INTEGER, DIMENSION(:,:,:), POINTER :: &
+         & ivlqc                !: Variable level QC
+      INTEGER, DIMENSION(:,:,:,:), POINTER :: &
+         & ivlqcf               !: Variable level QC flags
+      INTEGER, DIMENSION(:,:), POINTER :: &
+         & iproc, &             !: Processor of obs (no I/O for this variable).
+         & iobsi, &             !: Global i index
+         & iobsj                !: Global j index
+      INTEGER, DIMENSION(:,:,:), POINTER :: &
+         & iobsk                !: k index
+      CHARACTER(LEN=ilengrid), DIMENSION(:), POINTER  :: &
+         & cgrid                !: Grid for this variable
+      CHARACTER(LEN=ilenname), DIMENSION(:), POINTER :: &
+         & caddname             !: Additional entries names
+      CHARACTER(LEN=ilenlong), DIMENSION(:,:), POINTER :: &
+         & caddlong             !: Additional entries long name (for output)
+      CHARACTER(LEN=ilenunit), DIMENSION(:,:), POINTER :: &
+         & caddunit             !: Additional entries units (for output)
+      REAL(fbsp), DIMENSION(:,:,:,:)   , POINTER :: &
+         & padd                 !: Additional entries
+      CHARACTER(LEN=ilenname), DIMENSION(:), POINTER :: &
+         & cextname             !: Extra variables names
+      CHARACTER(LEN=ilenlong), DIMENSION(:), POINTER :: &
+         & cextlong             !: Extra variables long name (for output)
+      CHARACTER(LEN=ilenunit), DIMENSION(:), POINTER :: &
+         & cextunit             !: Extra variables units (for output)
+      REAL(fbsp), DIMENSION(:,:,:)   , POINTER :: &
+         & pext                 !: Extra variables
+   END TYPE obfbdata
+
+   PRIVATE putvaratt_obfbdata
+
+   !!----------------------------------------------------------------------
+   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
+   !! $Id$
+   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
+   !!----------------------------------------------------------------------
+
+CONTAINS
+
+   SUBROUTINE init_obfbdata( fbdata )
+      !!----------------------------------------------------------------------
+      !!                    ***  ROUTINE init_obfbdata  ***
+      !!
+      !! ** Purpose :   Initialize sizes in obfbdata structure
+      !!
+      !! ** Method  :   
+      !!
+      !! ** Action  : 
+      !!
+      !!----------------------------------------------------------------------
+      !! * Arguments
+      TYPE(obfbdata) :: fbdata      ! obsfbdata structure
+
+      fbdata%nvar   = 0
+      fbdata%nobs   = 0
+      fbdata%nlev   = 0
+      fbdata%nadd   = 0
+      fbdata%next   = 0
+      fbdata%nqcf   = idefnqcf
+      fbdata%lalloc = .FALSE.
+      fbdata%lgrid  = .FALSE.
+
+   END SUBROUTINE init_obfbdata
+   
+   SUBROUTINE alloc_obfbdata( fbdata, kvar, kobs, klev, kadd, kext, lgrid, &
+      &                       kqcf)
+      !!----------------------------------------------------------------------
+      !!                    ***  ROUTINE alloc_obfbdata  ***
+      !!
+      !! ** Purpose :   Allocate data in an obfbdata structure
+      !!
+      !! ** Method  :   
+      !!
+      !! ** Action  : 
+      !!
+      !!----------------------------------------------------------------------
+      !! * Arguments
+      TYPE(obfbdata) ::  fbdata          ! obsfbdata structure to be allocated
+      INTEGER, INTENT(IN) :: kvar        ! Number of variables
+      INTEGER, INTENT(IN) :: kobs        ! Number of observations
+      INTEGER, INTENT(IN) :: klev        ! Number of levels
+      INTEGER, INTENT(IN) :: kadd        ! Number of additional entries
+      INTEGER, INTENT(IN) :: kext        ! Number of extra variables
+      LOGICAL, INTENT(IN) :: lgrid       ! Include grid search information
+      INTEGER, OPTIONAL ::  kqcf         ! Number of words for QC flags
+      !! * Local variables
+      INTEGER :: ji
+      INTEGER :: jv
+
+      ! Check allocation status and deallocate previous allocated structures
+
+      IF ( fbdata%lalloc ) THEN
+         CALL dealloc_obfbdata( fbdata )
+      ENDIF
+
+      ! Set dimensions
+
+      fbdata%lalloc = .TRUE.
+      fbdata%nvar   = kvar
+      fbdata%nobs   = kobs
+      fbdata%nlev   = MAX( klev, 1 )
+      fbdata%nadd   = kadd
+      fbdata%next   = kext
+      IF ( PRESENT(kqcf) ) THEN
+         fbdata%nqcf = kqcf
+      ELSE
+         fbdata%nqcf = idefnqcf
+      ENDIF
+
+      ! Set data not depending on number of observations
+
+      fbdata%cdjuldref  = REPEAT( 'X', ilenjuld )
+
+      ! Allocate and initialize standard data 
+
+      ALLOCATE( &
+         & fbdata%cname(fbdata%nvar),   &
+         & fbdata%coblong(fbdata%nvar), &
+         & fbdata%cobunit(fbdata%nvar)  &
+         & )
+      DO ji = 1, fbdata%nvar
+         WRITE(fbdata%cname(ji),'(A,I2.2)')'V_',ji
+         fbdata%coblong(ji) = REPEAT( ' ', ilenlong )
+         fbdata%cobunit(ji) = REPEAT( ' ', ilenunit )
+      END DO
+
+      ! Optionally also store grid search information
+      
+      IF ( lgrid ) THEN
+         ALLOCATE ( &
+            & fbdata%cgrid(fbdata%nvar) &
+            & )
+         fbdata%cgrid(:)      = REPEAT( 'X', ilengrid )
+         fbdata%lgrid         = .TRUE.
+      ENDIF
+         
+      ! Allocate and initialize additional entries if present
+         
+      IF ( fbdata%nadd > 0 ) THEN
+         ALLOCATE( &
+            & fbdata%caddname(fbdata%nadd),              &
+            & fbdata%caddlong(fbdata%nadd, fbdata%nvar), &
+            & fbdata%caddunit(fbdata%nadd, fbdata%nvar)  &
+            & )
+         DO ji = 1, fbdata%nadd
+            WRITE(fbdata%caddname(ji),'(A,I2.2)')'A',ji
+         END DO
+         DO jv = 1, fbdata%nvar
+            DO ji = 1, fbdata%nadd
+               fbdata%caddlong(ji,jv) = REPEAT( ' ', ilenlong )
+               fbdata%caddunit(ji,jv) = REPEAT( ' ', ilenunit )
+            END DO
+         END DO
+      ENDIF
+         
+      ! Allocate and initialize additional variables if present
+         
+      IF ( fbdata%next > 0 ) THEN
+         ALLOCATE( &
+            & fbdata%cextname(fbdata%next), &
+            & fbdata%cextlong(fbdata%next), &
+            & fbdata%cextunit(fbdata%next)  &
+            & )
+         DO ji = 1, fbdata%next
+            WRITE(fbdata%cextname(ji),'(A,I2.2)')'E_',ji
+            fbdata%cextlong(ji) = REPEAT( ' ', ilenlong )
+            fbdata%cextunit(ji) = REPEAT( ' ', ilenunit )
+         END DO
+      ENDIF
+
+      ! Data depending on number of observations is only allocated if nobs>0
+
+      IF ( fbdata%nobs > 0 ) THEN
+
+         ALLOCATE( &
+            & fbdata%cdwmo(fbdata%nobs),                                      &
+            & fbdata%cdtyp(fbdata%nobs),                                      &
+            & fbdata%ioqc(fbdata%nobs),                                       &
+            & fbdata%ioqcf(fbdata%nqcf,fbdata%nobs),                          &
+            & fbdata%ipqc(fbdata%nobs),                                       &
+            & fbdata%ipqcf(fbdata%nqcf,fbdata%nobs),                          &
+            & fbdata%itqc(fbdata%nobs),                                       &
+            & fbdata%itqcf(fbdata%nqcf,fbdata%nobs),                          &
+            & fbdata%idqc(fbdata%nlev,fbdata%nobs),                           &
+            & fbdata%idqcf(fbdata%nqcf,fbdata%nlev,fbdata%nobs),              &
+            & fbdata%plam(fbdata%nobs),                                       &
+            & fbdata%pphi(fbdata%nobs),                                       &
+            & fbdata%pdep(fbdata%nlev,fbdata%nobs),                           &
+            & fbdata%ptim(fbdata%nobs),                                       &
+            & fbdata%kindex(fbdata%nobs),                                     &
+            & fbdata%ivqc(fbdata%nobs,fbdata%nvar),                           &
+            & fbdata%ivqcf(fbdata%nqcf,fbdata%nobs,fbdata%nvar),              &
+            & fbdata%ivlqc(fbdata%nlev,fbdata%nobs,fbdata%nvar),              &
+            & fbdata%ivlqcf(fbdata%nqcf,fbdata%nlev,fbdata%nobs,fbdata%nvar), &
+            & fbdata%pob(fbdata%nlev,fbdata%nobs,fbdata%nvar)                 &
+            & )
+         fbdata%kindex(:)       = fbimdi
+         fbdata%cdwmo(:)        = REPEAT( 'X', ilenwmo )
+         fbdata%cdtyp(:)        = REPEAT( 'X', ilentyp )
+         fbdata%ioqc(:)         = fbimdi
+         fbdata%ioqcf(:,:)      = fbimdi
+         fbdata%ipqc(:)         = fbimdi
+         fbdata%ipqcf(:,:)      = fbimdi
+         fbdata%itqc(:)         = fbimdi
+         fbdata%itqcf(:,:)      = fbimdi
+         fbdata%idqc(:,:)       = fbimdi
+         fbdata%idqcf(:,:,:)    = fbimdi
+         fbdata%plam(:)         = fbrmdi
+         fbdata%pphi(:)         = fbrmdi
+         fbdata%pdep(:,:)       = fbrmdi
+         fbdata%ptim(:)         = fbrmdi
+         fbdata%ivqc(:,:)       = fbimdi
+         fbdata%ivqcf(:,:,:)    = fbimdi
+         fbdata%ivlqc(:,:,:)    = fbimdi
+         fbdata%ivlqcf(:,:,:,:) = fbimdi
+         fbdata%pob(:,:,:)      = fbrmdi
+         
+         ! Optionally also store grid search information
+         
+         IF ( lgrid ) THEN
+            ALLOCATE ( &
+               & fbdata%iproc(fbdata%nobs,fbdata%nvar),            &
+               & fbdata%iobsi(fbdata%nobs,fbdata%nvar),            &
+               & fbdata%iobsj(fbdata%nobs,fbdata%nvar),            &
+               & fbdata%iobsk(fbdata%nlev,fbdata%nobs,fbdata%nvar) &
+               & )
+            fbdata%iproc(:,:)    = fbimdi
+            fbdata%iobsi(:,:)    = fbimdi
+            fbdata%iobsj(:,:)    = fbimdi
+            fbdata%iobsk(:,:,:)  = fbimdi
+            fbdata%lgrid         = .TRUE.
+         ENDIF
+         
+         ! Allocate and initialize additional entries if present
+         
+         IF ( fbdata%nadd > 0 ) THEN
+            ALLOCATE( &
+               & fbdata%padd(fbdata%nlev,fbdata%nobs,fbdata%nadd,fbdata%nvar) &
+               & )
+            fbdata%padd(:,:,:,:) = fbrmdi
+         ENDIF
+         
+         ! Allocate and initialize additional variables if present
+         
+         IF ( fbdata%next > 0 ) THEN
+            ALLOCATE( &
+               & fbdata%pext(fbdata%nlev,fbdata%nobs,fbdata%next) &
+               & )
+            fbdata%pext(:,:,:) = fbrmdi
+         ENDIF
+
+      ENDIF
+
+   END SUBROUTINE alloc_obfbdata
+
+   SUBROUTINE dealloc_obfbdata( fbdata )
+      !!----------------------------------------------------------------------
+      !!                    ***  ROUTINE dealloc_obfbdata  ***
+      !!
+      !! ** Purpose :   Deallocate data in an obfbdata strucure
+      !!
+      !! ** Method  :   
+      !!
+      !! ** Action  : 
+      !!
+      !!----------------------------------------------------------------------
+      !! * Arguments
+      TYPE(obfbdata) :: fbdata      ! obsfbdata structure
+
+      ! Deallocate data 
+
+      DEALLOCATE( &
+         & fbdata%cname,  &
+         & fbdata%coblong,&
+         & fbdata%cobunit &
+         & )
+
+      ! Deallocate optional grid search information
+      
+      IF ( fbdata%lgrid ) THEN
+         DEALLOCATE ( &
+            & fbdata%cgrid  &
+            & )
+      ENDIF
+
+      ! Deallocate additional entries
+
+      IF ( fbdata%nadd > 0 ) THEN
+         DEALLOCATE( &
+            & fbdata%caddname, &
+            & fbdata%caddlong, &
+            & fbdata%caddunit  &
+            & )
+      ENDIF
+
+      ! Deallocate extra variables
+
+      IF ( fbdata%next > 0 ) THEN
+         DEALLOCATE( &
+            & fbdata%cextname, &
+            & fbdata%cextlong, &
+            & fbdata%cextunit  &
+            & )
+      ENDIF
+
+      ! Deallocate arrays depending on number of obs (if nobs>0 only).
+
+      IF ( fbdata%nobs > 0 ) THEN
+
+         DEALLOCATE( &
+            & fbdata%cdwmo,  &
+            & fbdata%cdtyp,  &
+            & fbdata%ioqc,   &
+            & fbdata%ioqcf,  &
+            & fbdata%ipqc,   &
+            & fbdata%ipqcf,  &
+            & fbdata%itqc,   &
+            & fbdata%itqcf,  &
+            & fbdata%idqc,   &
+            & fbdata%idqcf,  &
+            & fbdata%plam,   &
+            & fbdata%pphi,   &
+            & fbdata%pdep,   &
+            & fbdata%ptim,   &
+            & fbdata%kindex, &
+            & fbdata%ivqc,   &
+            & fbdata%ivqcf,  &
+            & fbdata%ivlqc,  &
+            & fbdata%ivlqcf, &
+            & fbdata%pob     &
+            & )
+
+
+         ! Deallocate optional grid search information
+      
+         IF ( fbdata%lgrid ) THEN
+            DEALLOCATE ( &
+               & fbdata%iproc, &
+               & fbdata%iobsi, &
+               & fbdata%iobsj, &
+               & fbdata%iobsk  & 
+               & )
+         ENDIF
+
+         ! Deallocate additional entries
+
+         IF ( fbdata%nadd > 0 ) THEN
+            DEALLOCATE( &
+               & fbdata%padd       &
+               & )
+         ENDIF
+
+         ! Deallocate extra variables
+
+         IF ( fbdata%next > 0 ) THEN
+            DEALLOCATE( &
+               & fbdata%pext       &
+               & )
+         ENDIF
+
+      ENDIF
+
+      ! Reset arrays sizes
+
+      fbdata%lalloc = .FALSE.
+      fbdata%lgrid  = .FALSE.
+      fbdata%nvar   = 0
+      fbdata%nobs   = 0
+      fbdata%nlev   = 0
+      fbdata%nadd   = 0
+      fbdata%next   = 0
+   
+   END SUBROUTINE dealloc_obfbdata
+
+   SUBROUTINE copy_obfbdata( fbdata1, fbdata2, kadd, kext, lgrid, kqcf )
+      !!----------------------------------------------------------------------
+      !!                    ***  ROUTINE copy_obfbdata  ***
+      !!
+      !! ** Purpose :   Copy an obfbdata structure
+      !!
+      !! ** Method  :   Copy all data from fbdata1 to fbdata2
+      !!                If fbdata2 is allocated it needs to be compliant
+      !!                with fbdata1.
+      !!                Additional entries can be added by setting nadd
+      !!                Additional extra fields can be added by setting next
+      !!                Grid information can be included with lgrid=.true.
+      !!
+      !! ** Action  : 
+      !!
+      !!----------------------------------------------------------------------
+      !! * Arguments
+      TYPE(obfbdata) :: fbdata1               ! Input obsfbdata structure
+      TYPE(obfbdata) :: fbdata2               ! Output obsfbdata structure
+      INTEGER, INTENT(IN), OPTIONAL :: kadd   ! Number of additional entries
+      INTEGER, INTENT(IN), OPTIONAL :: kext   ! Number of extra variables
+      INTEGER, INTENT(IN), OPTIONAL :: kqcf   ! Number of words per qc flags
+      LOGICAL, OPTIONAL :: lgrid              ! Grid info on output file
+
+      !! * Local variables
+      INTEGER :: nadd
+      INTEGER :: next
+      INTEGER :: nqcf
+      LOGICAL :: llgrid
+      INTEGER :: jv
+      INTEGER :: je
+      INTEGER :: ji
+      INTEGER :: jk
+      INTEGER :: jq
+
+      ! Check allocation status of fbdata1
+
+      IF ( .NOT. fbdata1%lalloc ) THEN
+         CALL fatal_error( 'copy_obfbdata: input data not allocated', &
+            &              __LINE__ )
+      ENDIF
+      
+      ! If nadd,next not specified use the ones from fbdata1
+      ! Otherwise check that they have large than the original ones
+      
+      IF ( PRESENT(kadd) ) THEN
+         nadd = kadd
+         IF ( nadd < fbdata1%nadd ) THEN
+            CALL warning    ( 'copy_obfbdata: ' // &
+               &              'nadd smaller than input nadd', __LINE__ )
+         ENDIF
+      ELSE
+         nadd = fbdata1%nadd
+      ENDIF
+      IF ( PRESENT(kext) ) THEN
+         next = kext
+         IF ( next < fbdata1%next ) THEN
+            CALL fatal_error( 'copy_obfbdata: ' // &
+               &              'next smaller than input next', __LINE__ )
+         ENDIF
+      ELSE
+         next = fbdata1%next
+      ENDIF
+      IF ( PRESENT(lgrid) ) THEN
+         llgrid = lgrid
+         IF ( fbdata1%lgrid .AND. (.NOT. llgrid) ) THEN
+            CALL fatal_error( 'copy_obfbdata: ' // &
+               &              'switching off grid info not possible', &
+               &              __LINE__ )
+         ENDIF
+      ELSE
+         llgrid = fbdata1%lgrid
+      ENDIF
+      IF ( PRESENT(kqcf) ) THEN
+         nqcf = kqcf
+         IF ( nqcf < fbdata1%nqcf ) THEN
+            CALL fatal_error( 'copy_obfbdata: ' // &
+               &              'nqcf smaller than input nqcf', __LINE__ )
+         ENDIF
+      ELSE
+         nqcf = fbdata1%nqcf
+      ENDIF
+
+      ! Check allocation status of fbdata2 and 
+      ! a) check that it conforms in size if already allocated
+      ! b) allocate it if not already allocated
+      
+      IF ( fbdata2%lalloc ) THEN
+         IF ( fbdata1%nvar > fbdata2%nvar ) THEN
+            CALL fatal_error( 'copy_obfbdata: ' // &
+               &              'output kvar smaller than input kvar', __LINE__ )
+         ENDIF
+         IF ( fbdata1%nobs > fbdata2%nobs ) THEN
+            CALL fatal_error( 'copy_obfbdata: ' // &
+               &              'output kobs smaller than input kobs', __LINE__ )
+         ENDIF
+         IF ( fbdata1%nlev > fbdata2%nlev ) THEN
+            CALL fatal_error( 'copy_obfbdata: ' // &
+               &              'output klev smaller than input klev', __LINE__ )
+         ENDIF
+         IF ( fbdata1%nadd > fbdata2%nadd ) THEN
+            CALL warning    ( 'copy_obfbdata: ' // &
+               &              'output nadd smaller than input nadd', __LINE__ )
+         ENDIF
+         IF ( fbdata1%next > fbdata2%next ) THEN
+            CALL fatal_error( 'copy_obfbdata: ' // &
+               &              'output next smaller than input next', __LINE__ )
+         ENDIF
+         IF ( fbdata1%lgrid .NEQV. fbdata2%lgrid ) THEN
+            CALL fatal_error( 'copy_obfbdata: ' // &
+               &              'lgrid inconsistent', __LINE__ )
+         ENDIF
+         IF ( fbdata1%next > fbdata2%next ) THEN
+            CALL fatal_error( 'copy_obfbdata: ' // &
+               &              'output next smaller than input next', __LINE__ )
+         ENDIF
+         IF ( fbdata1%nqcf > fbdata2%nqcf ) THEN
+            CALL fatal_error( 'copy_obfbdata: ' // &
+               &              'output  smaller than input kext', __LINE__ )
+         ENDIF
+      ELSE
+         CALL alloc_obfbdata( fbdata2, fbdata1%nvar, fbdata1%nobs, &
+            &                 fbdata1%nlev, nadd, next, llgrid, kqcf = nqcf )
+      ENDIF
+
+      ! Copy the header data
+
+      fbdata2%cdjuldref = fbdata1%cdjuldref
+
+      DO ji = 1, fbdata1%nobs
+         fbdata2%cdwmo(ji)  = fbdata1%cdwmo(ji)
+         fbdata2%cdtyp(ji)  = fbdata1%cdtyp(ji)
+         fbdata2%ioqc(ji)   = fbdata1%ioqc(ji)
+         fbdata2%ipqc(ji)   = fbdata1%ipqc(ji)
+         fbdata2%itqc(ji)   = fbdata1%itqc(ji)
+         fbdata2%plam(ji)   = fbdata1%plam(ji)
+         fbdata2%pphi(ji)   = fbdata1%pphi(ji)
+         fbdata2%ptim(ji)   = fbdata1%ptim(ji)
+         fbdata2%kindex(ji) = fbdata1%kindex(ji)
+         DO jq = 1, fbdata1%nqcf
+            fbdata2%ioqcf(jq,ji)  = fbdata1%ioqcf(jq,ji)
+            fbdata2%ipqcf(jq,ji)  = fbdata1%ipqcf(jq,ji)
+            fbdata2%itqcf(jq,ji)  = fbdata1%itqcf(jq,ji)
+         END DO
+         DO jk = 1, fbdata1%nlev
+            fbdata2%idqc(jk,ji)  = fbdata1%idqc(jk,ji)
+            fbdata2%pdep(jk,ji)  = fbdata1%pdep(jk,ji)
+            DO jq = 1, fbdata1%nqcf
+               fbdata2%idqcf(jq,jk,ji) = fbdata1%idqcf(jq,jk,ji)
+            END DO
+         END DO
+      END DO
+
+      ! Copy the variable data
+
+      DO jv = 1, fbdata1%nvar
+         fbdata2%cname(jv) = fbdata1%cname(jv)
+         fbdata2%coblong(jv) = fbdata1%coblong(jv)
+         fbdata2%cobunit(jv) = fbdata1%cobunit(jv)
+         DO ji = 1, fbdata1%nobs
+            fbdata2%ivqc(ji,jv)  = fbdata1%ivqc(ji,jv)
+            DO jq = 1, fbdata1%nqcf
+               fbdata2%ivqcf(jq,ji,jv) = fbdata1%ivqcf(jq,ji,jv)
+            END DO
+            DO jk = 1, fbdata1%nlev
+               fbdata2%ivlqc(jk,ji,jv)  = fbdata1%ivlqc(jk,ji,jv)
+               fbdata2%pob(jk,ji,jv)    = fbdata1%pob(jk,ji,jv)
+               DO jq = 1, fbdata1%nqcf
+                  fbdata2%ivlqcf(jq,jk,ji,jv) = fbdata1%ivlqcf(jq,jk,ji,jv)
+               END DO
+            END DO
+         END DO
+      END DO
+
+      ! Copy grid information
+      
+      IF ( fbdata1%lgrid ) THEN
+         DO jv = 1, fbdata1%nvar
+            fbdata2%cgrid(jv) = fbdata1%cgrid(jv)
+            DO ji = 1, fbdata1%nobs
+               fbdata2%iproc(ji,jv) = fbdata1%iproc(ji,jv)
+               fbdata2%iobsi(ji,jv) = fbdata1%iobsi(ji,jv)
+               fbdata2%iobsj(ji,jv) = fbdata1%iobsj(ji,jv)
+               DO jk = 1, fbdata1%nlev
+                  fbdata2%iobsk(jk,ji,jv)  = fbdata1%iobsk(jk,ji,jv)
+               END DO
+            END DO
+         END DO
+      ENDIF
+
+      ! Copy additional information
+      
+      DO je = 1, MIN( fbdata1%nadd, fbdata2%nadd )
+         fbdata2%caddname(je) = fbdata1%caddname(je)
+      END DO
+      DO jv = 1, fbdata1%nvar
+         DO je = 1, MIN( fbdata1%nadd, fbdata2%nadd )
+            fbdata2%caddlong(je,jv) = fbdata1%caddlong(je,jv)
+            fbdata2%caddunit(je,jv) = fbdata1%caddunit(je,jv)
+            DO ji = 1, fbdata1%nobs
+               DO jk = 1, fbdata1%nlev
+                  fbdata2%padd(jk,ji,je,jv) = fbdata1%padd(jk,ji,je,jv)
+               END DO
+            END DO
+         END DO
+      END DO
+      
+      ! Copy extra information
+
+      DO je = 1, fbdata1%next
+         fbdata2%cextname(je) = fbdata1%cextname(je)
+         fbdata2%cextlong(je) = fbdata1%cextlong(je)
+         fbdata2%cextunit(je) = fbdata1%cextunit(je)
+      END DO
+      DO je = 1, fbdata1%next
+         DO ji = 1, fbdata1%nobs
+            DO jk = 1, fbdata1%nlev
+               fbdata2%pext(jk,ji,je) = fbdata1%pext(jk,ji,je)
+            END DO
+         END DO
+      END DO
+
+   END SUBROUTINE copy_obfbdata
+
+   SUBROUTINE subsamp_obfbdata( fbdata1, fbdata2, llvalid )
+      !!----------------------------------------------------------------------
+      !!                    ***  ROUTINE susbamp_obfbdata  ***
+      !!
+      !! ** Purpose :   Subsample an obfbdata structure based on the
+      !!                logical mask.
+      !!
+      !! ** Method  :   Copy all data from fbdata1 to fbdata2 if
+      !!                llvalid(obs)==true
+      !!
+      !! ** Action  : 
+      !!
+      !!----------------------------------------------------------------------
+      !! * Arguments
+      TYPE(obfbdata) :: fbdata1           ! Input obsfbdata structure
+      TYPE(obfbdata) :: fbdata2           ! Output obsfbdata structure
+      LOGICAL, DIMENSION(fbdata1%nobs) :: llvalid     ! Grid info on output file
+      !! * Local variables
+      INTEGER :: nobs
+      INTEGER :: jv
+      INTEGER :: je
+      INTEGER :: ji
+      INTEGER :: jk
+      INTEGER :: jq
+      INTEGER :: ij
+
+      ! Check allocation status of fbdata1
+
+      IF ( .NOT. fbdata1%lalloc ) THEN
+         CALL fatal_error( 'copy_obfbdata: input data not allocated', &
+            &              __LINE__ )
+      ENDIF
+      
+      ! Check allocation status of fbdata2 and abort if already allocated
+      
+      IF ( fbdata2%lalloc ) THEN
+         CALL fatal_error( 'subsample_obfbdata: ' // &
+            &              'fbdata2 already allocated', __LINE__ )
+      ENDIF
+      
+      ! Count number of subsampled observations
+
+      nobs = COUNT(llvalid)
+      
+      ! Allocate new data structure
+
+      CALL alloc_obfbdata( fbdata2, fbdata1%nvar, nobs, &
+         &                 fbdata1%nlev, fbdata1%nadd, fbdata1%next, &
+         &                 fbdata1%lgrid, kqcf = fbdata1%nqcf )
+
+      ! Copy the header data
+
+      fbdata2%cdjuldref = fbdata1%cdjuldref
+      
+      ij = 0
+      DO ji = 1, fbdata1%nobs
+         IF ( llvalid(ji) ) THEN
+            ij = ij +1
+            fbdata2%cdwmo(ij)  = fbdata1%cdwmo(ji)
+            fbdata2%cdtyp(ij)  = fbdata1%cdtyp(ji)
+            fbdata2%ioqc(ij)   = fbdata1%ioqc(ji)
+            fbdata2%ipqc(ij)   = fbdata1%ipqc(ji)
+            fbdata2%itqc(ij)   = fbdata1%itqc(ji)
+            fbdata2%plam(ij)   = fbdata1%plam(ji)
+            fbdata2%pphi(ij)   = fbdata1%pphi(ji)
+            fbdata2%ptim(ij)   = fbdata1%ptim(ji)
+            fbdata2%kindex(ij) = fbdata1%kindex(ji)
+            DO jq = 1, fbdata1%nqcf
+               fbdata2%ioqcf(jq,ij)  = fbdata1%ioqcf(jq,ji)
+               fbdata2%ipqcf(jq,ij)  = fbdata1%ipqcf(jq,ji)
+               fbdata2%itqcf(jq,ij)  = fbdata1%itqcf(jq,ji)
+            END DO
+            DO jk = 1, fbdata1%nlev
+               fbdata2%idqc(jk,ij)  = fbdata1%idqc(jk,ji)
+               fbdata2%pdep(jk,ij)  = fbdata1%pdep(jk,ji)
+               DO jq = 1, fbdata1%nqcf
+                  fbdata2%idqcf(jq,jk,ij) = fbdata1%idqcf(jq,jk,ji)
+               END DO
+            END DO
+         ENDIF
+      END DO
+
+      ! Copy the variable data
+
+      DO jv = 1, fbdata1%nvar
+         fbdata2%cname(jv) = fbdata1%cname(jv)
+         fbdata2%coblong(jv) = fbdata1%coblong(jv)
+         fbdata2%cobunit(jv) = fbdata1%cobunit(jv)
+         ij = 0
+         DO ji = 1, fbdata1%nobs
+            IF ( llvalid(ji) ) THEN
+               ij = ij + 1
+               fbdata2%ivqc(ij,jv)  = fbdata1%ivqc(ji,jv)
+               DO jq = 1, fbdata1%nqcf
+                  fbdata2%ivqcf(jq,ij,jv) = fbdata1%ivqcf(jq,ji,jv)
+               END DO
+               DO jk = 1, fbdata1%nlev
+                  fbdata2%ivlqc(jk,ij,jv)  = fbdata1%ivlqc(jk,ji,jv)
+                  fbdata2%pob(jk,ij,jv)    = fbdata1%pob(jk,ji,jv)
+                  DO jq = 1, fbdata1%nqcf
+                     fbdata2%ivlqcf(jq,jk,ij,jv) = fbdata1%ivlqcf(jq,jk,ji,jv)
+                  END DO
+               END DO
+            ENDIF
+         END DO
+      END DO
+
+      ! Copy grid information
+      
+      IF ( fbdata1%lgrid ) THEN
+         DO jv = 1, fbdata1%nvar
+            fbdata2%cgrid(jv) = fbdata1%cgrid(jv)
+            ij = 0
+            DO ji = 1, fbdata1%nobs
+               IF ( llvalid(ji) ) THEN
+                  ij = ij + 1
+                  fbdata2%iproc(ij,jv) = fbdata1%iproc(ji,jv)
+                  fbdata2%iobsi(ij,jv) = fbdata1%iobsi(ji,jv)
+                  fbdata2%iobsj(ij,jv) = fbdata1%iobsj(ji,jv)
+                  DO jk = 1, fbdata1%nlev
+                     fbdata2%iobsk(jk,ij,jv)  = fbdata1%iobsk(jk,ji,jv)
+                  END DO
+               ENDIF
+            END DO
+         END DO
+      ENDIF
+
+      ! Copy additional information
+      
+      DO je = 1, fbdata1%nadd
+         fbdata2%caddname(je) = fbdata1%caddname(je)
+      END DO
+      DO jv = 1, fbdata1%nvar
+         DO je = 1, fbdata1%nadd
+            fbdata2%caddlong(je,jv) = fbdata1%caddlong(je,jv)
+            fbdata2%caddunit(je,jv) = fbdata1%caddunit(je,jv)
+            ij = 0
+            DO ji = 1, fbdata1%nobs
+               IF ( llvalid(ji) ) THEN
+                  ij = ij + 1
+                  DO jk = 1, fbdata1%nlev
+                     fbdata2%padd(jk,ij,je,jv) = fbdata1%padd(jk,ji,je,jv)
+                  END DO
+               ENDIF
+            END DO
+         END DO
+      END DO
+      
+      ! Copy extra information
+
+      DO je = 1, fbdata1%next
+         fbdata2%cextname(je) = fbdata1%cextname(je)
+         fbdata2%cextlong(je) = fbdata1%cextlong(je)
+         fbdata2%cextunit(je) = fbdata1%cextunit(je)
+      END DO
+      DO je = 1, fbdata1%next
+         ij = 0
+         DO ji = 1, fbdata1%nobs
+            IF ( llvalid(ji) ) THEN
+               ij = ij + 1
+               DO jk = 1, fbdata1%nlev
+                  fbdata2%pext(jk,ij,je) = fbdata1%pext(jk,ji,je)
+               END DO
+            ENDIF
+         END DO
+      END DO
+
+   END SUBROUTINE subsamp_obfbdata
+
+   SUBROUTINE merge_obfbdata( nsets, fbdatain, fbdataout, iset, inum, iind )
+      !!----------------------------------------------------------------------
+      !!                    ***  ROUTINE merge_obfbdata  ***
+      !!
+      !! ** Purpose :   Merge multiple obfbdata structures into an one.
+      !!
+      !! ** Method  :   The order of elements is based on the indices in
+      !!                iind.
+      !!                All input data are assumed to be consistent. This
+      !!                is assumed to be checked before calling this routine.
+      !!                Likewise output data is assume to be consistent as 
+      !!                well without error checking.
+      !!
+      !! ** Action  : 
+      !!
+      !!----------------------------------------------------------------------
+      !! * Arguments
+      INTEGER, INTENT(IN):: nsets      ! Number of input data sets 
+      TYPE(obfbdata), DIMENSION(nsets) :: fbdatain  ! Input obsfbdata structure
+      TYPE(obfbdata) :: fbdataout      ! Output obsfbdata structure
+      INTEGER, INTENT(IN), DIMENSION(fbdataout%nobs) :: &
+         & iset                 ! Set number for a given obs.
+      INTEGER, INTENT(IN), DIMENSION(fbdataout%nobs) :: &
+         & inum                 ! Number within set for an obs
+      INTEGER, INTENT(IN), DIMENSION(fbdataout%nobs) :: &
+         & iind                 ! Indices for copying.
+      !! * Local variables
+
+      INTEGER :: js
+      INTEGER :: jo
+      INTEGER :: jv
+      INTEGER :: je
+      INTEGER :: ji
+      INTEGER :: jk
+      INTEGER :: jq
+
+      ! Check allocation status of fbdatain 
+      
+      DO js = 1, nsets
+         IF ( .NOT. fbdatain(js)%lalloc ) THEN
+            CALL fatal_error( 'merge_obfbdata: input data not allocated', &
+               &              __LINE__ )
+         ENDIF
+      END DO
+
+      ! Check allocation status of fbdataout
+      
+      IF ( .NOT.fbdataout%lalloc ) THEN
+         CALL fatal_error( 'merge_obfbdata: output data not allocated', &
+            &              __LINE__ )
+      ENDIF
+
+      ! Merge various names
+
+      DO jv = 1, fbdatain(1)%nvar
+         fbdataout%cname(jv) = fbdatain(1)%cname(jv)
+         fbdataout%coblong(jv) = fbdatain(1)%coblong(jv)
+         fbdataout%cobunit(jv) = fbdatain(1)%cobunit(jv)
+         IF ( fbdatain(1)%lgrid ) THEN
+            fbdataout%cgrid(jv) = fbdatain(1)%cgrid(jv)
+         ENDIF
+      END DO
+      DO jv = 1, fbdatain(1)%nadd
+         fbdataout%caddname(jv) = fbdatain(1)%caddname(jv)
+      END DO
+      DO jv = 1, fbdatain(1)%nvar
+         DO je = 1, fbdatain(1)%nadd
+            fbdataout%caddlong(je,jv) = fbdatain(1)%caddlong(je,jv)
+            fbdataout%caddunit(je,jv) = fbdatain(1)%caddunit(je,jv)
+         END DO
+      END DO
+      DO jv = 1, fbdatain(1)%next
+         fbdataout%cextname(jv) = fbdatain(1)%cextname(jv)
+         fbdataout%cextlong(jv) = fbdatain(1)%cextlong(jv)
+         fbdataout%cextunit(jv) = fbdatain(1)%cextunit(jv)
+      END DO
+      fbdataout%cdjuldref = fbdatain(1)%cdjuldref
+
+      ! Loop over total views
+
+      DO jo = 1, fbdataout%nobs
+
+         js = iset(iind(jo))
+         ji = inum(iind(jo))
+
+         ! Merge the header data
+
+         fbdataout%cdwmo(jo)  = fbdatain(js)%cdwmo(ji)
+         fbdataout%cdtyp(jo)  = fbdatain(js)%cdtyp(ji)
+         fbdataout%ioqc(jo)   = fbdatain(js)%ioqc(ji)
+         fbdataout%ipqc(jo)   = fbdatain(js)%ipqc(ji)
+         fbdataout%itqc(jo)   = fbdatain(js)%itqc(ji)
+         fbdataout%plam(jo)   = fbdatain(js)%plam(ji)
+         fbdataout%pphi(jo)   = fbdatain(js)%pphi(ji)
+         fbdataout%ptim(jo)   = fbdatain(js)%ptim(ji)
+         fbdataout%kindex(jo) = fbdatain(js)%kindex(ji)
+         DO jq = 1, fbdatain(js)%nqcf
+            fbdataout%ioqcf(jq,jo)  = fbdatain(js)%ioqcf(jq,ji)
+            fbdataout%ipqcf(jq,jo)  = fbdatain(js)%ipqcf(jq,ji)
+            fbdataout%itqcf(jq,jo)  = fbdatain(js)%itqcf(jq,ji)
+         END DO
+         DO jk = 1, fbdatain(js)%nlev
+            fbdataout%pdep(jk,jo)  = fbdatain(js)%pdep(jk,ji)
+            fbdataout%idqc(jk,jo)  = fbdatain(js)%idqc(jk,ji)
+            DO jq = 1, fbdatain(js)%nqcf
+               fbdataout%idqcf(jq,jk,jo) = fbdatain(js)%idqcf(jq,jk,ji)
+            END DO
+         END DO
+
+         ! Merge the variable data
+
+         DO jv = 1, fbdatain(js)%nvar
+            fbdataout%ivqc(jo,jv)  = fbdatain(js)%ivqc(ji,jv)
+            DO jq = 1, fbdatain(js)%nqcf
+               fbdataout%ivqcf(jq,jo,jv) = fbdatain(js)%ivqcf(jq,ji,jv)
+            END DO
+            DO jk = 1, fbdatain(js)%nlev
+               fbdataout%ivlqc(jk,jo,jv)  = fbdatain(js)%ivlqc(jk,ji,jv)
+               fbdataout%pob(jk,jo,jv)    = fbdatain(js)%pob(jk,ji,jv)
+               DO jq = 1, fbdatain(js)%nqcf
+                  fbdataout%ivlqcf(jq,jk,jo,jv) = &
+                     &                     fbdatain(js)%ivlqcf(jq,jk,ji,jv)
+               END DO
+            END DO
+         END DO
+
+         ! Merge grid information
+         
+         IF ( fbdatain(js)%lgrid ) THEN
+            DO jv = 1, fbdatain(js)%nvar
+               fbdataout%cgrid(jv) = fbdatain(js)%cgrid(jv)
+               fbdataout%iproc(jo,jv) = fbdatain(js)%iproc(ji,jv)
+               fbdataout%iobsi(jo,jv) = fbdatain(js)%iobsi(ji,jv)
+               fbdataout%iobsj(jo,jv) = fbdatain(js)%iobsj(ji,jv)
+               DO jk = 1, fbdatain(js)%nlev
+                  fbdataout%iobsk(jk,jo,jv)  = fbdatain(js)%iobsk(jk,ji,jv)
+               END DO
+            END DO
+         ENDIF
+
+         ! Merge additional information
+      
+         DO jv = 1, fbdatain(js)%nvar
+            DO je = 1, fbdatain(js)%nadd
+               DO jk = 1, fbdatain(js)%nlev
+                  fbdataout%padd(jk,jo,je,jv) = fbdatain(js)%padd(jk,ji,je,jv)
+               END DO
+            END DO
+         END DO
+         
+         ! Merge extra information
+         
+         DO je = 1, fbdatain(js)%next
+            DO jk = 1, fbdatain(js)%nlev
+               fbdataout%pext(jk,jo,je) = fbdatain(js)%pext(jk,ji,je)
+            END DO
+         END DO
+
+      END DO
+
+   END SUBROUTINE merge_obfbdata
+
+   SUBROUTINE write_obfbdata( cdfilename, fbdata )
+      !!----------------------------------------------------------------------
+      !!                    ***  ROUTINE write_obfbdata  ***
+      !!
+      !! ** Purpose :   Write an obfbdata structure into a netCDF file.
+      !!
+      !! ** Method  :   Decides which output wrapper to use. 
+      !!
+      !! ** Action  : 
+      !!
+      !!----------------------------------------------------------------------
+      !! * Arguments
+      CHARACTER(len=*) :: cdfilename ! Output filename
+      TYPE(obfbdata)   :: fbdata     ! obsfbdata structure
+#if defined key_offobsoper
+      IF (ln_cl4) THEN
+          ! Class 4 file output stream
+          CALL write_obfbdata_cl( cdfilename, fbdata )
+      ELSE
+#endif
+          ! Standard feedback file output stream
+          CALL write_obfbdata_fb( cdfilename, fbdata )
+#if defined key_offobsoper
+      ENDIF
+#endif
+   END SUBROUTINE write_obfbdata
+
+   SUBROUTINE write_obfbdata_fb( cdfilename, fbdata )
+      !!----------------------------------------------------------------------
+      !!                    ***  ROUTINE write_obfbdata  ***
+      !!
+      !! ** Purpose :   Write an obfbdata structure into a netCDF file.
+      !!
+      !! ** Method  :   
+      !!
+      !! ** Action  : 
+      !!
+      !!----------------------------------------------------------------------
+      !! * Arguments
+      CHARACTER(len=*) :: cdfilename ! Output filename
+      TYPE(obfbdata)   :: fbdata     ! obsfbdata structure
+      !! * Local variables
+      CHARACTER(LEN=14), PARAMETER :: cpname = 'write_obfbdata'
+      ! Dimension ids
+      INTEGER :: idfile
+      INTEGER :: idodim
+      INTEGER :: idldim
+      INTEGER :: idvdim
+      INTEGER :: idadim
+      INTEGER :: idedim
+      INTEGER :: idsndim
+      INTEGER :: idsgdim
+      INTEGER :: idswdim
+      INTEGER :: idstdim
+      INTEGER :: idjddim
+      INTEGER :: idqcdim
+      INTEGER :: idvard
+      INTEGER :: idaddd
+      INTEGER :: idextd
+      INTEGER :: idcdwmo
+      INTEGER :: idcdtyp
+      INTEGER :: idplam
+      INTEGER :: idpphi
+      INTEGER :: idpdep
+      INTEGER :: idptim
+      INTEGER :: idptimr
+      INTEGER :: idioqc         
+      INTEGER :: idioqcf         
+      INTEGER :: idipqc
+      INTEGER :: idipqcf
+      INTEGER :: iditqc
+      INTEGER :: iditqcf
+      INTEGER :: ididqc
+      INTEGER :: ididqcf
+      INTEGER :: idkindex
+      INTEGER, DIMENSION(fbdata%nvar) :: &
+         & idpob,    &
+         & idivqc,   &
+         & idivqcf,  &
+         & idivlqc,  &
+         & idivlqcf, &
+         & idiobsi,  &
+         & idiobsj,  &
+         & idiobsk,  &
+         & idcgrid
+      INTEGER, DIMENSION(fbdata%nadd,fbdata%nvar) :: idpadd
+      INTEGER, DIMENSION(fbdata%next) :: idpext
+      INTEGER, DIMENSION(1) :: incdim1
+      INTEGER, DIMENSION(2) :: incdim2
+      INTEGER, DIMENSION(3) :: incdim3
+      INTEGER, DIMENSION(4) :: incdim4
+
+      INTEGER :: jv
+      INTEGER :: je
+      INTEGER :: ioldfill
+      CHARACTER(len=nf90_max_name) :: &
+         & cdtmp
+      CHARACTER(len=16), PARAMETER :: &
+         & cdqcconv = 'q where q =[0,9]'
+      CHARACTER(len=24), PARAMETER :: &
+         & cdqcfconv = 'NEMOVAR flag conventions'
+      CHARACTER(len=ilenlong) :: &
+         & cdltmp
+
+      ! Open output filename
+
+      CALL chkerr( nf90_create( TRIM( cdfilename ), nf90_clobber, idfile ), &
+         &         cpname, __LINE__ )
+      CALL chkerr( nf90_set_fill( idfile, nf90_nofill, ioldfill ), &
+         &         cpname, __LINE__ )
+      CALL chkerr( nf90_put_att( idfile, nf90_global, 'title', &
+         &                       'NEMO observation operator output' ), &
+         &         cpname, __LINE__ )
+      CALL chkerr( nf90_put_att( idfile, nf90_global, 'Convention', &
+         &                       'NEMO unified observation operator output' ),&
+         &         cpname,__LINE__ )
+
+      ! Create the dimensions
+
+      CALL chkerr( nf90_def_dim( idfile, 'N_OBS'  , fbdata%nobs, idodim ),  &
+         &         cpname,__LINE__ )
+      CALL chkerr( nf90_def_dim( idfile, 'N_LEVELS', fbdata%nlev, idldim ), &
+         &         cpname,__LINE__ )
+      CALL chkerr( nf90_def_dim( idfile, 'N_VARS', fbdata%nvar, idvdim ), &
+         &         cpname,__LINE__ )
+      CALL chkerr( nf90_def_dim( idfile, 'N_QCF', fbdata%nqcf, idqcdim ),&
+         &         cpname,__LINE__ )
+      IF ( fbdata%nadd > 0 ) THEN
+         CALL chkerr( nf90_def_dim( idfile, 'N_ENTRIES', fbdata%nadd, idadim ), &
+            &         cpname,__LINE__ )
+      ENDIF
+      IF ( fbdata%next > 0 ) THEN
+         CALL chkerr( nf90_def_dim( idfile, 'N_EXTRA', fbdata%next, idedim ), &
+            &         cpname,__LINE__ )
+      ENDIF
+      CALL chkerr( nf90_def_dim( idfile, 'STRINGNAM', ilenname, idsndim ), &
+         &         cpname,__LINE__ )
+      IF (fbdata%lgrid) THEN
+         CALL chkerr( nf90_def_dim( idfile, 'STRINGGRID', ilengrid, idsgdim ),&
+            &         cpname,__LINE__ )
+      ENDIF
+      CALL chkerr( nf90_def_dim( idfile, 'STRINGWMO', ilenwmo, idswdim ), &
+         &         cpname,__LINE__ )
+      CALL chkerr( nf90_def_dim( idfile, 'STRINGTYP', ilentyp, idstdim ), &
+         &         cpname,__LINE__ )
+      CALL chkerr( nf90_def_dim( idfile, 'STRINGJULD', ilenjuld, idjddim ), &
+         &         cpname,__LINE__ )
+      
+      ! Define netCDF variables for header information
+      
+      incdim2(1) = idsndim
+      incdim2(2) = idvdim
+
+      CALL chkerr( nf90_def_var( idfile, 'VARIABLES', nf90_char, incdim2, &
+         &                       idvard ), cpname, __LINE__ )
+      CALL putvaratt_obfbdata( idfile, idvard, &
+         &                     'List of variables in feedback files' )
+      
+      IF ( fbdata%nadd > 0 ) THEN
+         incdim2(1) = idsndim
+         incdim2(2) = idadim
+         CALL chkerr( nf90_def_var( idfile, 'ENTRIES', nf90_char, incdim2, &
+            &                       idaddd ), cpname, __LINE__ )
+         CALL putvaratt_obfbdata( idfile, idaddd,  &
+            &                     'List of additional entries for each '// &
+            &                     'variable in feedback files' )
+      ENDIF
+   
+      IF ( fbdata%next > 0 ) THEN
+         incdim2(1) = idsndim
+         incdim2(2) = idedim
+         CALL chkerr( nf90_def_var( idfile, 'EXTRA', nf90_char, incdim2, &
+            &                       idextd ), cpname, __LINE__ )
+         CALL putvaratt_obfbdata(  idfile, idextd, &
+            &                      'List of extra variables' )
+      ENDIF
+
+      incdim2(1) = idswdim
+      incdim2(2) = idodim
+      CALL chkerr( nf90_def_var( idfile, 'STATION_IDENTIFIER', &
+         &                       nf90_char, incdim2, &
+         &                       idcdwmo ), cpname, __LINE__ )
+      CALL putvaratt_obfbdata(  idfile, idcdwmo, &
+         &                      'Station identifier' )
+      incdim2(1) = idstdim
+      incdim2(2) = idodim
+      CALL chkerr( nf90_def_var( idfile, 'STATION_TYPE', &
+         &                       nf90_char, incdim2, &
+         &                       idcdtyp ), cpname, __LINE__ )
+      CALL putvaratt_obfbdata(  idfile, idcdtyp, &
+         &                      'Code instrument type' )
+      incdim1(1) = idodim
+      CALL chkerr( nf90_def_var( idfile, 'LONGITUDE', &
+         &                       nf90_double, incdim1, &
+         &                       idplam ), cpname, __LINE__ )
+      CALL putvaratt_obfbdata(  idfile, idplam, &
+         &                      'Longitude', cdunits = 'degrees_east', &
+         &                      rfillvalue = fbrmdi )
+      CALL chkerr( nf90_def_var( idfile, 'LATITUDE', &
+         &                       nf90_double, incdim1, &
+         &                       idpphi ), cpname, __LINE__ )
+      CALL putvaratt_obfbdata(  idfile, idpphi, &
+         &                      'Latitude', cdunits = 'degrees_north', &
+         &                      rfillvalue = fbrmdi )
+      incdim2(1) = idldim
+      incdim2(2) = idodim
+      CALL chkerr( nf90_def_var( idfile, 'DEPTH', &
+         &                       nf90_double, incdim2, &
+         &                       idpdep ), cpname, __LINE__ )
+      CALL putvaratt_obfbdata(  idfile, idpdep, &
+         &                      'Depth', cdunits = 'metre', &
+         &                      rfillvalue = fbrmdi )
+      incdim3(1) = idqcdim
+      incdim3(2) = idldim
+      incdim3(3) = idodim
+      CALL chkerr( nf90_def_var( idfile, 'DEPTH_QC', &
+         &                       nf90_int, incdim2, &
+         &                       ididqc ), cpname, __LINE__ )
+      CALL putvaratt_obfbdata(  idfile, ididqc, &
+         &                      'Quality on depth',  &
+         &                      conventions = cdqcconv, &
+         &                      ifillvalue = 0 )
+      CALL chkerr( nf90_def_var( idfile, 'DEPTH_QC_FLAGS', &
+         &                       nf90_int, incdim3, &
+         &                       ididqcf ), cpname, __LINE__ )
+      CALL putvaratt_obfbdata(  idfile, ididqcf, &
+         &                      'Quality flags on depth',  &
+         &                      conventions = cdqcfconv )
+      CALL chkerr( nf90_def_var( idfile, 'JULD', &
+         &                       nf90_double, incdim1, &
+         &                       idptim ), cpname, __LINE__ )
+      CALL putvaratt_obfbdata(  idfile, idptim, &
+         &                      'Julian day', &
+         &                      cdunits = 'days since JULD_REFERENCE', &
+         &                      conventions = 'relative julian days with '// &
+         &                                 'decimal part (as parts of day)', &
+         &                      rfillvalue = fbrmdi )
+      incdim1(1) = idjddim
+      CALL chkerr( nf90_def_var( idfile, 'JULD_REFERENCE', &
+         &                       nf90_char, incdim1, &
+         &                       idptimr ), cpname, __LINE__ )
+      CALL putvaratt_obfbdata(  idfile, idptimr, &
+         &                      'Date of reference for julian days ', &
+         &                      conventions = 'YYYYMMDDHHMMSS' )
+      incdim1(1) = idodim
+      CALL chkerr( nf90_def_var( idfile, 'OBSERVATION_QC', &
+         &                       nf90_int, incdim1, &
+         &                       idioqc ), cpname, __LINE__ )
+      CALL putvaratt_obfbdata(  idfile, idioqc, &
+         &                      'Quality on observation',  &
+         &                      conventions = cdqcconv, &
+         &                      ifillvalue = 0 )
+      incdim2(1) = idqcdim
+      incdim2(2) = idodim
+      CALL chkerr( nf90_def_var( idfile, 'OBSERVATION_QC_FLAGS', &
+         &                       nf90_int, incdim2, &
+         &                       idioqcf ), cpname, __LINE__ )
+      CALL putvaratt_obfbdata(  idfile, idioqcf, &
+         &                      'Quality flags on observation',  &
+         &                      conventions = cdqcfconv, &
+         &                      ifillvalue = 0 )
+      CALL chkerr( nf90_def_var( idfile, 'POSITION_QC', &
+         &                       nf90_int, incdim1, &
+         &                       idipqc ), cpname, __LINE__ )
+      CALL putvaratt_obfbdata(  idfile, idipqc, &
+         &                      'Quality on position (latitude and longitude)',  &
+         &                      conventions = cdqcconv, &
+         &                      ifillvalue = 0 )
+      CALL chkerr( nf90_def_var( idfile, 'POSITION_QC_FLAGS', &
+         &                       nf90_int, incdim2, &
+         &                       idipqcf ), cpname, __LINE__ )
+      CALL putvaratt_obfbdata(  idfile, idipqcf, &
+         &                      'Quality flags on position',  &
+         &                      conventions = cdqcfconv, &
+         &                      ifillvalue = 0 )
+      CALL chkerr( nf90_def_var( idfile, 'JULD_QC', &
+         &                       nf90_int, incdim1, &
+         &                       iditqc ), cpname, __LINE__ )
+      CALL putvaratt_obfbdata(  idfile, iditqc, &
+         &                      'Quality on date and time',  &
+         &                      conventions = cdqcconv, &
+         &                      ifillvalue = 0 )
+      CALL chkerr( nf90_def_var( idfile, 'JULD_QC_FLAGS', &
+         &                       nf90_int, incdim2, &
+         &                       iditqcf ), cpname, __LINE__ )
+      CALL putvaratt_obfbdata(  idfile, iditqcf, &
+         &                      'Quality flags on date and time',  &
+         &                      conventions = cdqcfconv, &
+         &                      ifillvalue = 0 )
+      CALL chkerr( nf90_def_var( idfile, 'ORIGINAL_FILE_INDEX', &
+         &                       nf90_int, incdim1, &
+         &                       idkindex ), cpname, __LINE__ )
+      CALL putvaratt_obfbdata(  idfile, idkindex, &
+         &                      'Index in original data file',  &
+         &                      ifillvalue = fbimdi )
+
+      ! Define netCDF variables for individual variables
+
+      DO jv = 1, fbdata%nvar
+
+         incdim1(1) = idodim
+         incdim2(1) = idldim
+         incdim2(2) = idodim
+         WRITE(cdtmp,'(2A)') TRIM(fbdata%cname(jv)),'_OBS'
+         CALL chkerr( nf90_def_var( idfile, cdtmp, nf90_float, &
+            &                       incdim2, idpob(jv) ), &
+            &         cpname, __LINE__ )
+         CALL putvaratt_obfbdata(  idfile, idpob(jv), &
+            &                      fbdata%coblong(jv),  &
+            &                      cdunits =  fbdata%cobunit(jv), &
+            &                      rfillvalue = fbrmdi )
+
+         IF ( fbdata%nadd > 0 ) THEN
+            DO je = 1, fbdata%nadd
+               WRITE(cdtmp,'(3A)') TRIM(fbdata%cname(jv)),'_',&
+                  &                TRIM(fbdata%caddname(je))
+               CALL chkerr( nf90_def_var( idfile, cdtmp, nf90_float, &
+                  &                       incdim2, idpadd(je,jv) ), &
+                  &         cpname, __LINE__ )
+               CALL putvaratt_obfbdata(  idfile, idpadd(je,jv), &
+                  &                      fbdata%caddlong(je,jv), &
+                  &                      cdunits =  fbdata%caddunit(je,jv), &
+                  &                      rfillvalue = fbrmdi )
+            END DO
+         ENDIF
+
+         cdltmp = fbdata%coblong(jv)
+         IF (( cdltmp(1:1) >= 'A' ).AND.( cdltmp(1:1) <= 'Z' )) &
+            & cdltmp(1:1) = ACHAR(IACHAR(cdltmp(1:1)) + 32)
+         WRITE(cdtmp,'(2A)') TRIM(fbdata%cname(jv)),'_QC'
+         CALL chkerr( nf90_def_var( idfile, cdtmp, nf90_int, &
+            &                       incdim1, idivqc(jv) ), &
+            &         cpname, __LINE__ )
+         CALL putvaratt_obfbdata(  idfile, idivqc(jv), &
+            &                      'Quality on '//cdltmp,  &
+            &                      conventions = cdqcconv, &
+            &                      ifillvalue = 0 )
+         incdim2(1) = idqcdim
+         incdim2(2) = idodim
+         WRITE(cdtmp,'(2A)') TRIM(fbdata%cname(jv)),'_QC_FLAGS'
+         CALL chkerr( nf90_def_var( idfile, cdtmp, nf90_int, &
+            &                       incdim2, idivqcf(jv) ), &
+            &         cpname, __LINE__ )
+         CALL putvaratt_obfbdata(  idfile, idivqcf(jv), &
+            &                      'Quality flags on '//cdltmp,  &
+            &                      conventions = cdqcfconv, &
+            &                      ifillvalue = 0 )
+         incdim2(1) = idldim
+         incdim2(2) = idodim
+         WRITE(cdtmp,'(2A)') TRIM(fbdata%cname(jv)),'_LEVEL_QC'
+         CALL chkerr( nf90_def_var( idfile, cdtmp, nf90_int, &
+            &                       incdim2, idivlqc(jv) ), &
+            &         cpname, __LINE__ )
+         CALL putvaratt_obfbdata(  idfile, idivlqc(jv), &
+            &                      'Quality for each level on '//cdltmp,  &
+            &                      conventions = cdqcconv, &
+            &                      ifillvalue = 0 )
+         incdim3(1) = idqcdim
+         incdim3(2) = idldim
+         incdim3(3) = idodim
+         WRITE(cdtmp,'(2A)') TRIM(fbdata%cname(jv)),'_LEVEL_QC_FLAGS'
+         CALL chkerr( nf90_def_var( idfile, cdtmp, nf90_int, &
+            &                       incdim3, idivlqcf(jv) ), &
+            &         cpname, __LINE__ )
+         CALL putvaratt_obfbdata(  idfile, idivlqcf(jv), &
+            &                      'Quality flags for each level on '//&
+            &                      cdltmp,  &
+            &                      conventions = cdqcfconv, &
+            &                      ifillvalue = 0 )
+
+         IF (fbdata%lgrid) THEN
+            incdim2(1) = idldim
+            incdim2(2) = idodim
+            WRITE(cdtmp,'(2A)') TRIM(fbdata%cname(jv)),'_IOBSI'
+            CALL chkerr( nf90_def_var( idfile, cdtmp, nf90_int, &
+               &                       incdim1, idiobsi(jv) ), &
+               &         cpname, __LINE__ )
+            CALL putvaratt_obfbdata(  idfile, idiobsi(jv), &
+               &                      'ORCA grid search I coordinate')
+            WRITE(cdtmp,'(2A)') TRIM(fbdata%cname(jv)),'_IOBSJ'
+            CALL chkerr( nf90_def_var( idfile, cdtmp, nf90_int, &
+               &                       incdim1, idiobsj(jv) ), &
+               &         cpname, __LINE__ )
+            CALL putvaratt_obfbdata(  idfile, idiobsj(jv), &
+               &                      'ORCA grid search J coordinate')
+            WRITE(cdtmp,'(2A)') TRIM(fbdata%cname(jv)),'_IOBSK'
+            CALL chkerr( nf90_def_var( idfile, cdtmp, nf90_int, &
+               &                       incdim2, idiobsk(jv) ), &
+               &         cpname, __LINE__ )
+            CALL putvaratt_obfbdata(  idfile, idiobsk(jv), &
+               &                      'ORCA grid search K coordinate')
+            incdim1(1) = idsgdim
+            WRITE(cdtmp,'(2A)') TRIM(fbdata%cname(jv)),'_GRID'
+            CALL chkerr( nf90_def_var( idfile, cdtmp, nf90_char, incdim1, &
+               &                       idcgrid(jv) ), cpname, __LINE__ )
+            CALL putvaratt_obfbdata(  idfile, idcgrid(jv), &
+               &                      'ORCA grid search grid (T,U,V)')
+         ENDIF
+
+      END DO
+
+      IF ( fbdata%next > 0 ) THEN
+         DO je = 1, fbdata%next
+            incdim2(1) = idldim
+            incdim2(2) = idodim
+            WRITE(cdtmp,'(A)') TRIM(fbdata%cextname(je))
+            CALL chkerr( nf90_def_var( idfile, cdtmp, nf90_float, &
+               &                       incdim2, idpext(je) ), &
+               &         cpname, __LINE__ )
+            CALL putvaratt_obfbdata(  idfile, idpext(je), &
+               &                      fbdata%cextlong(je),  &
+               &                      cdunits =  fbdata%cextunit(je), &
+               &                      rfillvalue = fbrmdi )
+         END DO
+      ENDIF
+      
+      ! Stop definitions
+
+      CALL chkerr( nf90_enddef( idfile ), cpname, __LINE__ )
+      
+      ! Write the variables
+      
+      CALL chkerr( nf90_put_var( idfile, idvard, fbdata%cname ), &
+         &         cpname, __LINE__ )
+      
+      IF ( fbdata%nadd > 0 ) THEN
+         CALL chkerr( nf90_put_var( idfile, idaddd, fbdata%caddname ), &
+            &         cpname, __LINE__ )
+      ENDIF
+      
+      IF ( fbdata%next > 0 ) THEN
+         CALL chkerr( nf90_put_var( idfile, idextd, fbdata%cextname ), &
+            &         cpname, __LINE__ )
+      ENDIF
+
+      CALL chkerr( nf90_put_var( idfile, idptimr, fbdata%cdjuldref ), &
+         &         cpname, __LINE__ )
+
+      ! Only write the data if observation is available
+      
+      IF ( fbdata%nobs > 0 ) THEN
+
+         CALL chkerr( nf90_put_var( idfile, idcdwmo, fbdata%cdwmo ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_put_var( idfile, idcdtyp, fbdata%cdtyp ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_put_var( idfile, idplam, fbdata%plam ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_put_var( idfile, idpphi, fbdata%pphi ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_put_var( idfile, idpdep, fbdata%pdep ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_put_var( idfile, idptim, fbdata%ptim ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_put_var( idfile, idioqc, fbdata%ioqc ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_put_var( idfile, idioqcf, fbdata%ioqcf ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_put_var( idfile, idipqc, fbdata%ipqc ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_put_var( idfile, idipqcf, fbdata%ipqcf ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_put_var( idfile, iditqc, fbdata%itqc ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_put_var( idfile, iditqcf, fbdata%itqcf ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_put_var( idfile, ididqc, fbdata%idqc ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_put_var( idfile, ididqcf, fbdata%idqcf ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_put_var( idfile, idkindex, fbdata%kindex ), &
+            &         cpname, __LINE__ )
+
+         DO jv = 1, fbdata%nvar
+            CALL chkerr( nf90_put_var( idfile, idpob(jv), fbdata%pob(:,:,jv) ), &
+               &         cpname, __LINE__ )
+            IF ( fbdata%nadd > 0 ) THEN
+               DO je = 1, fbdata%nadd
+                  CALL chkerr( nf90_put_var( idfile, idpadd(je,jv), &
+                     &                       fbdata%padd(:,:,je,jv) ), &
+                     &         cpname, __LINE__ )
+               END DO
+            ENDIF
+            CALL chkerr( nf90_put_var( idfile, idivqc(jv), &
+               &                       fbdata%ivqc(:,jv) ),&
+               &         cpname, __LINE__ )
+            CALL chkerr( nf90_put_var( idfile, idivqcf(jv), &
+               &                       fbdata%ivqcf(:,:,jv) ),&
+               &         cpname, __LINE__ )
+            CALL chkerr( nf90_put_var( idfile, idivlqc(jv), &
+               &                       fbdata%ivlqc(:,:,jv) ),&
+               &         cpname, __LINE__ )
+            CALL chkerr( nf90_put_var( idfile, idivlqcf(jv), &
+               &                       fbdata%ivlqcf(:,:,:,jv) ),&
+               &         cpname, __LINE__ )
+            IF (fbdata%lgrid) THEN
+               CALL chkerr( nf90_put_var( idfile, idiobsi(jv), &
+                  &                       fbdata%iobsi(:,jv) ),&
+                  &         cpname, __LINE__ )
+               CALL chkerr( nf90_put_var( idfile, idiobsj(jv), &
+                  &                       fbdata%iobsj(:,jv) ),&
+                  &         cpname, __LINE__ )
+               CALL chkerr( nf90_put_var( idfile, idiobsk(jv), &
+                  &                       fbdata%iobsk(:,:,jv) ),&
+                  &         cpname, __LINE__ )
+               CALL chkerr( nf90_put_var( idfile, idcgrid(jv), &
+                  &                       fbdata%cgrid(jv) ), &
+                  &         cpname, __LINE__ )
+            ENDIF
+         END DO
+
+         IF ( fbdata%next > 0 ) THEN
+            DO je = 1, fbdata%next
+               CALL chkerr( nf90_put_var( idfile, idpext(je), &
+                  &                       fbdata%pext(:,:,je) ), &
+                  &         cpname, __LINE__ )
+            END DO
+         ENDIF
+
+      ENDIF
+
+      ! Close the file
+
+      CALL chkerr( nf90_close( idfile ), cpname, __LINE__ )
+
+      
+   END SUBROUTINE write_obfbdata_fb
+
+#if defined key_offobsoper
+   SUBROUTINE write_obfbdata_cl(cdfilename, fbdata)
+      !!----------------------------------------------------------------------
+      !!                    ***  ROUTINE write_obfbdata_cl  ***
+      !!
+      !! ** Purpose : Write an obfbdata structure into a class 4 file.
+      !!
+      !! ** Method  : 1. Allocate memory needed by off_wri
+      !!              2. Map obfbdata into allocated memory
+      !!              3. Pass mapped data to off_wri
+      !!              4. Deallocate memory
+      !!----------------------------------------------------------------------
+      USE dom_oce, ONLY: narea
+      USE ooo_write
+      USE ooo_data
+      IMPLICIT NONE
+      !! * Arguments
+      CHARACTER(len=*) :: cdfilename ! Feedback filename
+      TYPE(obfbdata)   :: fbdata     ! obsfbdata structure
+      !! * Local variables
+      CHARACTER(len=17), PARAMETER :: cpname = 'write_obfbdata_cl'
+      CHARACTER(len=64) :: &
+              & cdate, &   !: class 4 file validity date 
+              & cconf, &   !: model configuration
+              & csys, &    !: model system
+              & ccont, &   !: contact email
+              & cinst, &   !: institution
+              & cversion   !: model version
+      CHARACTER(len=8) :: &
+              & ckind      !: observation kind
+      CHARACTER(len=3) :: cfield
+      INTEGER :: kobs, &   !: number of observations
+              &  kvars, &  !: number of physical variables
+              &  kdeps, &  !: number of observed depths
+              &  kfcst, &  !: number of forecasts
+              &  kifcst, & !: current forecast number
+              &  kproc     !: processor number
+      INTEGER, DIMENSION(:, :, :), ALLOCATABLE :: &
+              &  kqc       !: quality control counterpart
+      INTEGER(KIND=2), DIMENSION(:, :, :), ALLOCATABLE :: &
+              &  k2qc       !: quality control counterpart
+      REAL(kind=fbdp) :: &
+              &  pmodjuld  !: model Julian day
+      REAL(kind=fbdp), DIMENSION(:), ALLOCATABLE :: &
+              &  plead, &  !: forecast lead time
+              &  plam, &   !: longitude of observation
+              &  pphi, &   !: latitude of observation
+              &  ptim      !: time of observation
+      REAL(kind=fbdp), DIMENSION(:, :), ALLOCATABLE :: &
+              &  pdep      !: depths of observations
+      REAL(kind=fbdp), DIMENSION(:, :, :), ALLOCATABLE :: &
+              &  pob, &    !: observation counterpart
+              &  pextra    !: extra field counterpart
+      REAL(kind=fbdp), DIMENSION(:, :, :), ALLOCATABLE :: &
+              &  pmod      !: model counterpart
+      CHARACTER(len=128) :: &
+              &  clfilename  !: class 4 file name
+      CHARACTER(len=128), DIMENSION(:), ALLOCATABLE :: &
+              &  ctype       !: Instrument type
+      CHARACTER(len=nf90_max_name) :: &
+              & cdtmp        !: NetCDF variable name
+      CHARACTER(len=8), DIMENSION(:), ALLOCATABLE :: &
+              &  cwmo, &     !: Instrument WMO ID
+              &  cunit, &    !: Instrument WMO ID
+              &  cvarname    !: Instrument WMO ID
+      INTEGER :: &
+              &  idep, &     !: Loop variable
+              &  ivar, &     !: Loop variable
+              &  iobs, &     !: Loop variable
+              &  ii, &       !: Loop variable
+              &  ij, &       !: Loop variable
+              &  ik, &       !: Loop variable
+              &  il          !: Loop variable
+      cconf = TRIM(cl4_cfg)
+      csys = TRIM(cl4_sys)
+      cversion = TRIM(cl4_vn)
+      ccont = TRIM(cl4_contact)
+      cinst = TRIM(cl4_inst)
+      cdate = TRIM(cl4_date)
+      CALL locate_kind(cdfilename, ckind)
+      kproc = narea
+      kfcst = cl4_fcst_len
+      kobs = fbdata%nobs
+      kdeps = fbdata%nlev
+      kvars = fbdata%nvar
+      IF (kobs .GT. 0) THEN
+         ALLOCATE(plam(kobs), &
+               &  pphi(kobs), &
+               &  ptim(kobs), &
+               &  plead(kfcst), &
+               &  pdep(kdeps, kobs), &
+               &  kqc(kdeps, kvars, kobs), &
+               &  k2qc(kdeps, kvars, kobs), &
+               &  pob(kdeps, kvars, kobs), &
+               &  pmod(kdeps, kvars, kobs), &
+               &  pextra(kdeps, kvars, kobs), &
+               &  ctype(kobs), &
+               &  cwmo(kobs), &
+               &  cunit(kvars), &
+               &  cvarname(kvars))
+         plam(:) = fbdata%plam(:)
+         pphi(:) = fbdata%pphi(:)
+         ptim(:) = fbdata%ptim(:)
+         pdep(:, :) = fbdata%pdep(:, :)
+         kqc(:,:,:) = 1.
+         DO ii = 1, kvars
+            cvarname(ii)  = fbdata%cname(ii)
+            cunit(ii)     = fbdata%cobunit(ii)
+         END DO
+
+         ! Quality control algorithm
+         k2qc(:,:,:) = NF90_FILL_SHORT
+         DO idep = 1,kdeps
+            DO ivar = 1, kvars
+               DO iobs = 1, kobs
+                  ! 1 symbolises good for fbdata
+                  ! fbimdi symbolises that qc has not been set
+                  ! Essentially, if any fbdata flag is not an element of {1, fbimdi}
+                  ! then set the class 4 flag to bad.
+                  ! Note: fbdata%ioqc is marked good if zero.
+                  IF (((fbdata%ioqc(iobs) /= 0) .AND. &
+                            & (fbdata%ioqc(iobs) /= fbimdi)) .OR. &
+                    & ((fbdata%ipqc(iobs) /= 1) .AND. &
+                            & (fbdata%ipqc(iobs) /= fbimdi)) .OR. &
+                    & ((fbdata%idqc(idep,iobs) /= 1) .AND. &
+                            & (fbdata%idqc(idep,iobs) /= fbimdi)) .OR. &
+                    & ((fbdata%ivqc(iobs,ivar) /= 1) .AND. &
+                            & (fbdata%ivqc(iobs,ivar) /= fbimdi)) .OR. &
+                    & ((fbdata%ivlqc(idep,iobs,ivar) /= 1) .AND. &
+                            & (fbdata%ivlqc(idep,iobs,ivar) /= fbimdi)) .OR. &
+                    & ((fbdata%itqc(iobs) /= 1) .AND. &
+                            & (fbdata%itqc(iobs) /= fbimdi))) THEN
+                     ! 1 symbolises bad for class 4 file
+                     k2qc(idep, ivar, iobs) = 1
+                  ELSE
+                     ! 0 symbolises good for class 4 file
+                     k2qc(idep, ivar, iobs) = 0
+                  END IF 
+               END DO
+            END DO
+         END DO
+
+         ! Permute observation dimensions
+         pob(:,:,:) = RESHAPE(fbdata%pob, (/kdeps, kvars, kobs/), &
+                            & ORDER=(/1, 3, 2/))
+
+         ! Explicit model counterpart dimension permutation
+         ! 1,2,3,4 --> 1,4,2,3
+         pmod(:,:,:) = fbrmdi
+         ij = cl4_fcst_idx(jimatch)
+         DO ii = 1,kdeps
+            DO ik = 1, kvars
+               DO il = 1, kobs
+                  pmod(ii,ik,il) = fbdata%padd(ii,il,1,ik)
+               END DO
+            END DO
+         END DO
+
+         ! Extra fields set to missing for now
+         pextra(:,:,:) = fbrmdi
+
+         ! Lead time of class 4 file is a global parameter
+         plead = cl4_leadtime(1:cl4_fcst_len)
+
+         ! Model Julian day
+         pmodjuld = cl4_modjuld
+
+         ! Observation types
+         ctype(:) = 'X'
+         DO ii = 1,kobs
+            ctype(ii) = fbdata%cdtyp(ii)
+         END DO
+
+         ! World Meteorology Organisation codes
+         cwmo(:) = fbdata%cdwmo(:)
+
+         ! Initialise class 4 file
+         CALL ooo_wri_init(cconf, csys, ckind, cversion, ccont, cinst, cdate, &
+                         & kproc, kobs, kvars, kdeps, kfcst, &
+                         & clfilename)
+
+         ! Write standard variables
+         CALL ooo_wri_default(clfilename, kobs, kvars, kfcst, kdeps, &
+                            & ctype, cwmo, cunit, cvarname, &
+                            & plam, pphi, pdep, ptim, pob, plead, &
+                            & k2qc, pmodjuld)
+         !! Write to optional variables
+         cdtmp = cl4_vars(jimatch)
+         IF ( (TRIM(cdtmp) == "forecast") .OR. &
+              (TRIM(cdtmp) == "persistence") ) THEN
+            !! 4D variables
+            CALL ooo_wri_extra(clfilename, TRIM(cdtmp), kdeps, kfcst, &
+                            &  kvars, kobs, (/ 1,ij,1,1 /), (/ kdeps,1,kvars,kobs /), pmod)
+         ELSE
+            !! 3D variables
+            CALL ooo_wri_extra(clfilename, TRIM(cdtmp), kdeps, &
+                            &  kvars, kobs, (/ 1,1,1 /), (/ kdeps,kvars,kobs /), pmod)
+         ENDIF
+
+         DEALLOCATE(plam, pphi, ptim, pdep, plead, kqc, k2qc, &
+                  & pob, pmod, pextra, ctype, cwmo, &
+                  & cunit, cvarname)
+      END IF
+   END SUBROUTINE write_obfbdata_cl
+#endif
+
+#if defined key_offobsoper
+   SUBROUTINE locate_kind(cdfilename, ckind)
+      !!----------------------------------------------------------------------
+      !!                    ***  ROUTINE locate_kind  ***
+      !!
+      !! ** Purpose : Detect which kind of class 4 file is being produced.
+      !!
+      !! ** Method  : 1. Inspect cdfilename for observation kind.
+      !!----------------------------------------------------------------------
+      IMPLICIT NONE
+      CHARACTER(len=*) :: cdfilename ! Feedback filename
+      CHARACTER(len=8) :: ckind
+      IF (cdfilename(1:3) == 'sst') THEN
+         ckind = 'SST'
+      ELSE IF (cdfilename(1:3) == 'sla') THEN
+         ckind = 'SLA'
+      ELSE IF (cdfilename(1:3) == 'pro') THEN
+         ckind = 'profile'
+      ELSE IF (cdfilename(1:3) == 'ena') THEN
+         ckind = 'profile'
+      ELSE IF (cdfilename(1:3) == 'sea') THEN
+         ckind = 'seaice'
+      ELSE
+         ckind = 'unknown'
+      END IF
+   END SUBROUTINE locate_kind
+#endif
+
+   SUBROUTINE putvaratt_obfbdata( idfile, idvar, cdlongname, cdunits, &
+      &                           conventions, cfillvalue, &
+      &                           ifillvalue, rfillvalue )
+      !!----------------------------------------------------------------------
+      !!                    ***  ROUTINE putvaratt_obfbdata  ***
+      !!
+      !! ** Purpose :   Write netcdf attributes for variable
+      !!
+      !! ** Method  :   
+      !!
+      !! ** Action  : 
+      !!
+      !!----------------------------------------------------------------------
+      !! * Arguments
+      INTEGER :: idfile                    ! File netcdf id.
+      INTEGER :: idvar                     ! Variable netcdf id.
+      CHARACTER(len=*) :: cdlongname       ! Long name for variable
+      CHARACTER(len=*), OPTIONAL :: cdunits       ! Units for variable
+      CHARACTER(len=*), OPTIONAL :: cfillvalue    ! Fill value for character variables
+      INTEGER, OPTIONAL :: ifillvalue             ! Fill value for integer variables
+      REAL(kind=fbsp), OPTIONAL :: rfillvalue     ! Fill value for real variables
+      CHARACTER(len=*), OPTIONAL :: conventions   ! Conventions for variable
+      !! * Local variables
+      CHARACTER(LEN=18), PARAMETER :: &
+         & cpname = 'putvaratt_obfbdata'
+
+      CALL chkerr( nf90_put_att( idfile, idvar, 'long_name', &
+         &                       TRIM(cdlongname) ), &
+         &                       cpname, __LINE__ )
+      
+      IF ( PRESENT(cdunits) ) THEN
+
+         CALL chkerr( nf90_put_att( idfile, idvar, 'units', &
+            &                       TRIM(cdunits) ), &
+            &                       cpname, __LINE__ )
+
+      ENDIF
+
+      IF ( PRESENT(conventions) ) THEN
+
+         CALL chkerr( nf90_put_att( idfile, idvar, 'Conventions', &
+            &                       TRIM(conventions) ), &
+            &                       cpname, __LINE__ )
+
+      ENDIF
+
+      IF ( PRESENT(cfillvalue) ) THEN
+
+         CALL chkerr( nf90_put_att( idfile, idvar, '_Fillvalue', &
+            &                       TRIM(cfillvalue) ), &
+            &                       cpname, __LINE__ )
+
+      ENDIF
+
+      IF ( PRESENT(ifillvalue) ) THEN
+
+         CALL chkerr( nf90_put_att( idfile, idvar, '_Fillvalue', &
+            &                       ifillvalue ), &
+            &                       cpname, __LINE__ )
+
+      ENDIF
+
+      IF ( PRESENT(rfillvalue) ) THEN
+
+         CALL chkerr( nf90_put_att( idfile, idvar, '_Fillvalue', &
+            &                       rfillvalue ), &
+            &                       cpname, __LINE__ )
+
+      ENDIF
+
+   END SUBROUTINE putvaratt_obfbdata
+
+   SUBROUTINE read_obfbdata( cdfilename, fbdata, ldgrid )
+      !!----------------------------------------------------------------------
+      !!                    ***  ROUTINE read_obfbdata  ***
+      !!
+      !! ** Purpose :   Read an obfbdata structure from a netCDF file.
+      !!
+      !! ** Method  :   
+      !!
+      !! ** Action  : 
+      !!
+      !!----------------------------------------------------------------------
+      !! * Arguments
+      CHARACTER(len=*) :: cdfilename  ! Input filename
+      TYPE(obfbdata)   :: fbdata      ! obsfbdata structure
+      LOGICAL, OPTIONAL :: ldgrid     ! Allow forcing of grid info
+      !! * Local variables
+      CHARACTER(LEN=14), PARAMETER :: cpname = 'read_obfbdata'
+      INTEGER :: idfile
+      INTEGER :: idodim
+      INTEGER :: idldim
+      INTEGER :: idvdim
+      INTEGER :: idadim
+      INTEGER :: idedim
+      INTEGER :: idgdim
+      INTEGER :: idvard
+      INTEGER :: idaddd
+      INTEGER :: idextd
+      INTEGER :: idcdwmo
+      INTEGER :: idcdtyp
+      INTEGER :: idplam
+      INTEGER :: idpphi
+      INTEGER :: idpdep
+      INTEGER :: idptim
+      INTEGER :: idptimr
+      INTEGER :: idioqc        
+      INTEGER :: idioqcf
+      INTEGER :: idipqc
+      INTEGER :: idipqcf
+      INTEGER :: ididqc
+      INTEGER :: ididqcf
+      INTEGER :: iditqc
+      INTEGER :: iditqcf
+      INTEGER :: idkindex
+      INTEGER, DIMENSION(:), ALLOCATABLE :: &
+         & idpob,    &
+         & idivqc,   &
+         & idivqcf,  &
+         & idivlqc,  &
+         & idivlqcf, &
+         & idiobsi,  &
+         & idiobsj,  &
+         & idiobsk,  &
+         & idcgrid,  &
+         & idpext
+      INTEGER, DIMENSION(:,:), ALLOCATABLE :: &
+         & idpadd
+      INTEGER :: jv
+      INTEGER :: je
+      INTEGER :: nvar
+      INTEGER :: nobs
+      INTEGER :: nlev
+      INTEGER :: nadd
+      INTEGER :: next
+      LOGICAL :: lgrid
+      CHARACTER(len=NF90_MAX_NAME) :: cdtmp
+
+      ! Check allocation status and deallocate previous allocated structures
+
+      IF ( fbdata%lalloc ) THEN
+         CALL dealloc_obfbdata( fbdata )
+      ENDIF
+
+      ! Open input filename
+
+      CALL chkerr( nf90_open( TRIM( cdfilename ), nf90_nowrite, idfile ), &
+         &         cpname, __LINE__ )
+
+      ! Get input dimensions
+
+      CALL chkerr( nf90_inq_dimid( idfile, 'N_OBS'  , idodim ),  &
+         &         cpname,__LINE__ )
+      CALL chkerr( nf90_inquire_dimension( idfile, idodim, len=nobs ), &
+         &         cpname,__LINE__ )
+      CALL chkerr( nf90_inq_dimid( idfile, 'N_LEVELS', idldim ), &
+         &         cpname,__LINE__ )
+      CALL chkerr( nf90_inquire_dimension( idfile, idldim, len=nlev ), &
+         &         cpname,__LINE__ )
+      CALL chkerr( nf90_inq_dimid( idfile, 'N_VARS', idvdim ), &
+         &         cpname,__LINE__ )
+      CALL chkerr( nf90_inquire_dimension( idfile, idvdim, len=nvar ), &
+         &         cpname,__LINE__ )
+      IF ( nf90_inq_dimid( idfile, 'N_ENTRIES',  idadim ) == 0 ) THEN
+         CALL chkerr( nf90_inquire_dimension( idfile, idadim, len=nadd ), &
+            &         cpname,__LINE__ )
+      ELSE
+         nadd = 0
+      ENDIF
+      IF ( nf90_inq_dimid( idfile, 'N_EXTRA',  idedim ) == 0 ) THEN
+         CALL chkerr( nf90_inquire_dimension( idfile, idedim, len=next ), &
+            &         cpname,__LINE__ )
+      ELSE
+         next = 0
+      ENDIF
+      !
+      ! Check if this input file  contains grid search informations
+      !
+      lgrid = ( nf90_inq_dimid( idfile, 'STRINGGRID',  idgdim ) == 0 )
+
+      ! Allocate data structure
+
+      IF ( PRESENT(ldgrid) ) THEN
+         CALL alloc_obfbdata( fbdata, nvar, nobs, nlev, nadd, next, &
+            & lgrid.OR.ldgrid )
+      ELSE
+         CALL alloc_obfbdata( fbdata, nvar, nobs, nlev, nadd, next, &
+            & lgrid )
+      ENDIF
+
+      ! Allocate netcdf identifiers
+
+      ALLOCATE( &
+         & idpob(fbdata%nvar),    &
+         & idivqc(fbdata%nvar),   &
+         & idivqcf(fbdata%nvar),  &
+         & idivlqc(fbdata%nvar),  &
+         & idivlqcf(fbdata%nvar), &
+         & idiobsi(fbdata%nvar),  &
+         & idiobsj(fbdata%nvar),  &
+         & idiobsk(fbdata%nvar),  &
+         & idcgrid(fbdata%nvar)   &
+         & )
+      IF ( fbdata%nadd > 0 ) THEN
+         ALLOCATE( &
+            & idpadd(fbdata%nadd,fbdata%nvar) &
+            & )
+      ENDIF
+      IF ( fbdata%next > 0 ) THEN
+         ALLOCATE( &
+            & idpext(fbdata%next) &
+            & )
+      ENDIF
+
+      ! Read variables for header information
+
+      CALL chkerr( nf90_inq_varid( idfile, 'VARIABLES',idvard ), &
+         &         cpname, __LINE__ )
+      CALL chkerr( nf90_get_var( idfile, idvard, fbdata%cname ), &
+         &         cpname, __LINE__ )
+      IF ( fbdata%nadd > 0 ) THEN
+         CALL chkerr( nf90_inq_varid( idfile, 'ENTRIES', idaddd ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_get_var( idfile, idaddd, fbdata%caddname ), &
+            &         cpname, __LINE__ )
+      ENDIF
+      IF ( fbdata%next > 0 ) THEN
+         CALL chkerr( nf90_inq_varid( idfile, 'EXTRA', idextd ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_get_var( idfile, idextd, fbdata%cextname ), &
+            &         cpname, __LINE__ )
+      ENDIF
+
+      CALL chkerr( nf90_inq_varid( idfile, 'JULD_REFERENCE', idptimr ), &
+         &         cpname, __LINE__ )
+      CALL chkerr( nf90_get_var( idfile, idptimr, fbdata%cdjuldref ), &
+         &         cpname, __LINE__ )
+
+      IF  ( fbdata%nobs > 0 ) THEN
+         
+         CALL chkerr( nf90_inq_varid( idfile, 'STATION_IDENTIFIER', idcdwmo ),&
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_get_var( idfile, idcdwmo, fbdata%cdwmo ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_inq_varid( idfile, 'STATION_TYPE', idcdtyp ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_get_var( idfile, idcdtyp, fbdata%cdtyp), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_inq_varid( idfile, 'LONGITUDE', idplam ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_get_var( idfile, idplam, fbdata%plam ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_inq_varid( idfile, 'LATITUDE', idpphi ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_get_var( idfile, idpphi, fbdata%pphi ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_inq_varid( idfile, 'DEPTH', idpdep ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_get_var( idfile, idpdep, fbdata%pdep ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_inq_varid( idfile, 'JULD', idptim ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_get_var( idfile, idptim, fbdata%ptim ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_inq_varid( idfile, 'OBSERVATION_QC', idioqc ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_get_var( idfile, idioqc, fbdata%ioqc ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_inq_varid( idfile, 'OBSERVATION_QC_FLAGS', idioqcf ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_get_var( idfile, idioqcf, fbdata%ioqcf ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_inq_varid( idfile, 'POSITION_QC', idipqc ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_get_var( idfile, idipqc, fbdata%ipqc ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_inq_varid( idfile, 'POSITION_QC_FLAGS', idipqcf ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_get_var( idfile, idipqcf, fbdata%ipqcf ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_inq_varid( idfile, 'DEPTH_QC', ididqc ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_get_var( idfile, ididqc, fbdata%idqc ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_inq_varid( idfile, 'DEPTH_QC_FLAGS', ididqcf ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_get_var( idfile, ididqcf, fbdata%idqcf ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_inq_varid( idfile, 'JULD_QC', iditqc ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_get_var( idfile, iditqc, fbdata%itqc ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_inq_varid( idfile, 'JULD_QC_FLAGS', iditqcf ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_get_var( idfile, iditqcf, fbdata%itqcf ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_inq_varid( idfile, 'ORIGINAL_FILE_INDEX', idkindex ), &
+            &         cpname, __LINE__ )
+         CALL chkerr( nf90_get_var( idfile, idkindex, fbdata%kindex ), &
+            &         cpname, __LINE__ )
+         
+         ! Read netCDF variables for individual variables
+         
+         DO jv = 1, fbdata%nvar
+            
+            WRITE(cdtmp,'(2A)') TRIM(fbdata%cname(jv)),'_OBS'
+            CALL chkerr( nf90_inq_varid( idfile, cdtmp, idpob(jv) ), &
+               &         cpname, __LINE__ )
+            CALL chkerr( nf90_get_var( idfile, idpob(jv), &
+               &                       fbdata%pob(:,:,jv) ), &
+               &         cpname, __LINE__ )
+            CALL getvaratt_obfbdata( idfile, idpob(jv), &
+               &                     fbdata%coblong(jv), &
+               &                     fbdata%cobunit(jv) )
+            
+            IF ( fbdata%nadd > 0 ) THEN
+               DO je = 1, fbdata%nadd
+                  WRITE(cdtmp,'(3A)') TRIM(fbdata%cname(jv)),'_',&
+                     &                TRIM(fbdata%caddname(je))
+                  CALL chkerr( nf90_inq_varid( idfile, cdtmp, idpadd(je,jv) ), &
+                     &         cpname, __LINE__ )
+                  CALL chkerr( nf90_get_var( idfile, idpadd(je,jv), &
+                     &                       fbdata%padd(:,:,je,jv) ), &
+                     &         cpname, __LINE__ )
+                  CALL getvaratt_obfbdata( idfile, idpadd(je,jv), &
+                     &                     fbdata%caddlong(je,jv), &
+                     &                     fbdata%caddunit(je,jv) )
+               END DO
+            ENDIF
+            
+            WRITE(cdtmp,'(2A)') TRIM(fbdata%cname(jv)),'_QC'
+            CALL chkerr( nf90_inq_varid( idfile, cdtmp, idivqc(jv) ), &
+            &         cpname, __LINE__ )
+            CALL chkerr( nf90_get_var( idfile, idivqc(jv), &
+               &                       fbdata%ivqc(:,jv) ), &
+               &         cpname, __LINE__ )
+            WRITE(cdtmp,'(2A)') TRIM(fbdata%cname(jv)),'_QC_FLAGS'
+            CALL chkerr( nf90_inq_varid( idfile, cdtmp, idivqcf(jv) ), &
+               &         cpname, __LINE__ )
+            CALL chkerr( nf90_get_var( idfile, idivqcf(jv), &
+               &                       fbdata%ivqcf(:,:,jv) ), &
+               &         cpname, __LINE__ )
+            WRITE(cdtmp,'(2A)') TRIM(fbdata%cname(jv)),'_LEVEL_QC'
+            CALL chkerr( nf90_inq_varid( idfile, cdtmp, idivlqc(jv) ), &
+               &         cpname, __LINE__ )
+            CALL chkerr( nf90_get_var( idfile, idivlqc(jv), &
+               &                       fbdata%ivlqc(:,:,jv) ), &
+               &         cpname, __LINE__ )
+            WRITE(cdtmp,'(2A)') TRIM(fbdata%cname(jv)),'_LEVEL_QC_FLAGS'
+            CALL chkerr( nf90_inq_varid( idfile, cdtmp, idivlqcf(jv) ), &
+               &         cpname, __LINE__ )
+            CALL chkerr( nf90_get_var( idfile, idivlqcf(jv), &
+               &                       fbdata%ivlqcf(:,:,:,jv) ), &
+               &         cpname, __LINE__ )
+            IF ( lgrid ) THEN
+               WRITE(cdtmp,'(2A)') TRIM(fbdata%cname(jv)),'_IOBSI'
+               CALL chkerr( nf90_inq_varid( idfile, cdtmp, idiobsi(jv) ), &
+                  &         cpname, __LINE__ )
+               CALL chkerr( nf90_get_var( idfile, idiobsi(jv), &
+                  &                       fbdata%iobsi(:,jv) ), &
+                  &         cpname, __LINE__ )
+               WRITE(cdtmp,'(2A)') TRIM(fbdata%cname(jv)),'_IOBSJ'
+               CALL chkerr( nf90_inq_varid( idfile, cdtmp, idiobsj(jv) ), &
+                  &         cpname, __LINE__ )
+               CALL chkerr( nf90_get_var( idfile, idiobsj(jv), &
+                  &                       fbdata%iobsj(:,jv) ), &
+                  &         cpname, __LINE__ )
+               WRITE(cdtmp,'(2A)') TRIM(fbdata%cname(jv)),'_IOBSK'
+               CALL chkerr( nf90_inq_varid( idfile, cdtmp, idiobsk(jv) ), &
+                  &         cpname, __LINE__ )
+               CALL chkerr( nf90_get_var( idfile, idiobsk(jv), &
+                  &                       fbdata%iobsk(:,:,jv) ), &
+                  &         cpname, __LINE__ )
+               WRITE(cdtmp,'(2A)') TRIM(fbdata%cname(jv)),'_GRID'
+               CALL chkerr( nf90_inq_varid( idfile, cdtmp, idcgrid(jv) ), &
+                  &         cpname, __LINE__ )
+               CALL chkerr( nf90_get_var( idfile, idcgrid(jv), &
+                  &                       fbdata%cgrid(jv) ), &
+                  &         cpname, __LINE__ )
+            ENDIF
+            
+         END DO
+         
+         IF ( fbdata%next > 0 ) THEN
+            DO je = 1, fbdata%next
+               WRITE(cdtmp,'(A)') TRIM(fbdata%cextname(je))
+               CALL chkerr( nf90_inq_varid( idfile, cdtmp, idpext(je) ), &
+                  &         cpname, __LINE__ )
+               CALL chkerr( nf90_get_var( idfile, idpext(je), &
+                  &                       fbdata%pext(:,:,je) ), &
+                  &         cpname, __LINE__ )
+               CALL getvaratt_obfbdata( idfile, idpext(je), &
+                  &                     fbdata%cextlong(je), &
+                  &                     fbdata%cextunit(je) )
+            END DO
+         ENDIF
+
+      ELSE ! if no observations only get attributes
+
+         DO jv = 1, fbdata%nvar            
+
+            WRITE(cdtmp,'(2A)') TRIM(fbdata%cname(jv)),'_OBS'
+            CALL chkerr( nf90_inq_varid( idfile, cdtmp, idpob(jv) ), &
+               &         cpname, __LINE__ )
+            CALL getvaratt_obfbdata( idfile, idpob(jv), &
+               &                     fbdata%coblong(jv), &
+               &                     fbdata%cobunit(jv) )
+            
+            IF ( fbdata%nadd > 0 ) THEN
+               DO je = 1, fbdata%nadd
+                  WRITE(cdtmp,'(3A)') TRIM(fbdata%cname(jv)),'_',&
+                     &                TRIM(fbdata%caddname(je))
+                  CALL chkerr( nf90_inq_varid( idfile, cdtmp, idpadd(je,jv) ), &
+                     &         cpname, __LINE__ )
+                  CALL getvaratt_obfbdata( idfile, idpadd(je,jv), &
+                     &                     fbdata%caddlong(je,jv), &
+                     &                     fbdata%caddunit(je,jv) )
+               END DO
+            ENDIF
+            
+         END DO
+         
+         IF ( fbdata%next > 0 ) THEN
+            DO je = 1, fbdata%next
+               WRITE(cdtmp,'(A)') TRIM(fbdata%cextname(je))
+               CALL chkerr( nf90_inq_varid( idfile, cdtmp, idpext(je) ), &
+                  &         cpname, __LINE__ )
+               CALL getvaratt_obfbdata( idfile, idpext(je), &
+                  &                     fbdata%cextlong(je), &
+                  &                     fbdata%cextunit(je) )
+            END DO
+         ENDIF
+
+      ENDIF
+
+      ! Close the file
+
+      CALL chkerr( nf90_close( idfile ), cpname, __LINE__ )
+
+   END SUBROUTINE read_obfbdata
+
+   SUBROUTINE getvaratt_obfbdata( idfile, idvar, cdlongname, cdunits )
+      !!----------------------------------------------------------------------
+      !!                    ***  ROUTINE putvaratt_obfbdata  ***
+      !!
+      !! ** Purpose :   Read netcdf attributes for variable
+      !!
+      !! ** Method  :   
+      !!
+      !! ** Action  : 
+      !!
+      !!----------------------------------------------------------------------
+      !! * Arguments
+      INTEGER :: idfile      ! File netcdf id.
+      INTEGER :: idvar       ! Variable netcdf id.
+      CHARACTER(len=*) :: cdlongname  ! Long name for variable
+      CHARACTER(len=*) :: cdunits     ! Units for variable
+      !! * Local variables
+      CHARACTER(LEN=18), PARAMETER :: cpname = 'getvaratt_obfbdata'
+
+      CALL chkerr( nf90_get_att( idfile, idvar, 'long_name', &
+         &                       cdlongname ), &
+         &                       cpname, __LINE__ )
+
+      CALL chkerr( nf90_get_att( idfile, idvar, 'units', &
+         &                       cdunits ), &
+         &                       cpname, __LINE__ )
+
+   END SUBROUTINE getvaratt_obfbdata
+
+END MODULE obs_fbm
Index: branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/off_data.F90
===================================================================
--- branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/off_data.F90	(revision 4102)
+++ 	(revision )
@@ -1,136 +1,0 @@
-MODULE off_data
-   USE par_kind, ONLY: lc
-   IMPLICIT NONE
-   INTEGER, PARAMETER :: MaxNumFiles = 1000
-
-   !! Class 4 file settings
-   INTEGER :: &
-           & cl4_fcst_idx(MaxNumFiles), & !: forecast indices
-           & cl4_match_len, &             !: number of match types
-           & cl4_fcst_len                 !: number of forecast days
-   CHARACTER(len=lc) :: &
-           & cl4_vars(MaxNumFiles), &  !: class 4 variables
-           & cl4_sys, &                !: class 4 system
-           & cl4_cfg, &                !: class 4 configuration
-           & cl4_date, &               !: class 4 date
-           & cl4_vn,  &                !: class 4 version
-           & cl4_prefix, &             !: class 4 prefix
-           & cl4_contact, &            !: class 4 contact
-           & cl4_inst                  !: class 4 institute
-   REAL ::   cl4_modjuld               !: model Julian day
-   REAL :: &
-      & cl4_leadtime(MaxNumFiles)      !: Lead time data
-
-   !! Offline obs_oper settings
-   CHARACTER(len=lc) :: &
-      & off_files(MaxNumFiles)         !: model files
-   INTEGER            :: &
-      & jifile, &                      !: current file list index
-      & n_files, &                     !: number of files
-      & jimatch, &                     !: current match
-      & nn_off_idx(MaxNumFiles), &     !: time_counter indices
-      & nn_off_freq                    !: read frequency in time steps
-   CHARACTER(len=128) :: &
-      & alt_file                       !: altimeter file
-CONTAINS
-   SUBROUTINE off_data_init( ld_cl4 )
-      !!----------------------------------------------------------------------
-      !!                    ***  SUBROUTINE off_data_init ***
-      !!
-      !! ** Purpose : To read namelists and initialise offline_oper run.
-      !!
-      !!----------------------------------------------------------------------
-      USE in_out_manager
-      IMPLICIT NONE
-      INTEGER            :: &
-         & jf                           !: file dummy loop index
-      LOGICAL :: lmask(MaxNumFiles)     !: Logical mask used for counting
-      LOGICAL, INTENT(IN) :: ld_cl4     !: Logical class 4 on/off
-
-      ! Standard offline obs_oper information
-      NAMELIST/namoff/off_files, nn_off_idx, nn_off_freq
-
-      ! Class 4 file specifiers
-      NAMELIST/namcl4/cl4_vars, cl4_sys, cl4_cfg, cl4_date, cl4_vn, &
-         &            cl4_prefix, cl4_contact, cl4_inst, cl4_leadtime, &
-         &            cl4_fcst_idx, cl4_fcst_len, cl4_match_len
-
-      ! Standard offline obs_oper initialisation
-      jimatch = 0                   !: match-up iteration variable 
-      jifile = 1                    !: input file iteration variable 
-      n_files = 0                   !: number of files to cycle through
-      off_files(:) = ''             !: list of files to read in
-      nn_off_idx(:) = 0             !: list of indices inside each file
-      nn_off_freq = -1              !: input frequency in time steps
-
-      ! Class 4 initialisation
-      cl4_leadtime(:) = 0           !: Lead time axis value for each file
-      cl4_fcst_len = 0              !: Length of the forecast dimension
-      cl4_match_len = 1             !: Number of match types
-      cl4_fcst_idx(:) = 0           !: output file forecast index
-      cl4_vars(:) = ''              !: output file variable names
-      cl4_sys = ''                  !: output file system
-      cl4_cfg = ''                  !: output file configuration
-      cl4_date = ''                 !: output file date string
-      cl4_vn = ''                   !: output file version
-      cl4_prefix = 'class4'         !: output file prefix
-      cl4_contact = ''              !: output file contact details
-      cl4_inst = ''                 !: output file institution
-
-      ! Standard offline obs_oper settings
-      READ(numnam, namoff)
-
-      ! Read class 4 output settings
-      IF (ld_cl4) THEN
-         READ(numnam, namcl4)
-      ENDIF
-
-      ! count input files
-      lmask(:) = .FALSE.
-      WHERE (off_files(:) /= '') lmask(:) = .TRUE.
-      n_files = COUNT(lmask)
-
-      !! Initialise sub obs window frequency
-      IF (nn_off_freq == -1) THEN
-         !! Run length
-         nn_off_freq = nitend - nit000 + 1
-      ENDIF
-
-      !! Print summary of settings
-      IF(lwp) THEN
-         WRITE(numout,*)
-         WRITE(numout,*) 'offline obs_oper : Initialization'
-         WRITE(numout,*) '~~~~~~~~~~~~~~~~~'
-         WRITE(numout,*) '   Namelist namoff : set offline obs_oper parameters' 
-         DO jf = 1, n_files
-            WRITE(numout,'(1X,2A)') '   Input forecast file name          forecastfile = ', &
-               TRIM(off_files(jf))
-            WRITE(numout,*) '   Input forecast file index        forecastindex = ', &
-               nn_off_idx(jf)
-            WRITE(numout,*) '   Output forecast leadtime index   leadtimeindex = ', &
-               cl4_fcst_idx(jf)
-            WRITE(numout,*) '   Output forecast leadtime value   leadtimevalue = ', &
-               cl4_leadtime(jf)
-            WRITE(numout,'(1X,2A)') '   Input class 4 variable       class 4 parameter = ', &
-               TRIM(cl4_vars(jf))
-         END DO
-         WRITE(numout, '(1X,2A)') '   Input class 4 system            class 4 system = ', &
-            TRIM(cl4_sys)
-         WRITE(numout, '(1X,2A)') '   Input class 4 config            class 4 config = ', &
-            TRIM(cl4_cfg)
-         WRITE(numout, '(1X,2A)') '   Input class 4 date                class 4 date = ', &
-            TRIM(cl4_date)
-         WRITE(numout, '(1X,2A)') '   Input class 4 version          class 4 version = ', &
-            TRIM(cl4_vn)
-         WRITE(numout, '(1X,2A)') '   Input class 4 prefix            class 4 prefix = ', &
-            TRIM(cl4_prefix)
-         WRITE(numout, '(1X,2A)') '   Input class 4 contact          class 4 contact = ', &
-            TRIM(cl4_contact)
-         WRITE(numout, '(1X,2A)') '   Input class 4 institute      class 4 institute = ', &
-            TRIM(cl4_inst)
-      END IF
-
-   END SUBROUTINE off_data_init
-
-END MODULE off_data
-
Index: branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/off_netcdf.F90
===================================================================
--- branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/off_netcdf.F90	(revision 4102)
+++ 	(revision )
@@ -1,104 +1,0 @@
-MODULE off_netcdf
-   USE obs_utils, ONLY: chkerr
-   USE obs_const, ONLY: obfillflt
-   USE par_oce
-   USE obs_types, ONLY: cwmonam1770
-   USE obs_sla, ONLY: nslasets
-
-   ! Define double precision obfillflt
-   REAL(kind=dp), PARAMETER :: obfilldbl=99999.
-
-   CONTAINS
-
-   SUBROUTINE yyyymmdd_to_ref_date(indate, intime, outstring)
-      IMPLICIT NONE
-      !----------------------------------------
-      ! Routine to create reference date string
-      !----------------------------------------
-      ! Routine arguments
-      CHARACTER(len=8),  INTENT(IN)  :: indate    ! yyyymmdd
-      CHARACTER(len=6),  INTENT(IN)  :: intime    ! hhmmss
-      CHARACTER(len=23), INTENT(OUT) :: outstring ! yyyy-mm-dd hh:mm:ss utc
-      ! Local variables
-      CHARACTER(len=4)  :: year      ! yyyy
-      CHARACTER(len=2)  :: month     ! mm
-      CHARACTER(len=2)  :: day       ! dd
-      CHARACTER(len=2)  :: hour      ! hh
-      CHARACTER(len=2)  :: minute    ! mm
-      CHARACTER(len=2)  :: second    ! ss
-
-      year   = indate(1:4)
-      month  = indate(5:6)
-      day    = indate(7:8)
-      hour   = intime(1:2)
-      minute = intime(3:4)
-      second = intime(5:6)
-
-      outstring = year//"-"//month//"-"//day//" "//hour//":"//minute//":"//second//" utc"
-      
-   END SUBROUTINE yyyymmdd_to_ref_date
-
-   SUBROUTINE date_format(date_str)
-      IMPLICIT NONE
-      !---------------------------------------
-      ! Routine to create creation date string
-      !---------------------------------------
-
-      ! Routine arguments
-      CHARACTER(len=*), INTENT(OUT) :: date_str
-      ! Local variables
-      CHARACTER(8)  :: date
-      CHARACTER(10) :: time
-
-      CHARACTER(10)  :: date_part
-      CHARACTER(8)  :: time_part
-
-      CALL date_and_time(DATE=date)
-      CALL date_and_time(TIME=time)
-      date_part  = date(1:4)//'/'//date(5:6)//'/'//date(7:8)
-      time_part  = time(1:2)//':'//time(3:4)//':'//time(5:6)
-
-      date_str   = date_part//' '//time_part
-
-   END SUBROUTINE date_format
-
-   SUBROUTINE inst_converter(wmo_inst_list, nprofs, obs_names)
-      IMPLICIT NONE
-      !---------------------------------------
-      ! Routine to convert WMO_INST_TYPE to
-      ! 'ARGO', 'TESAC', 'XBT', 'BUOYS'
-      !---------------------------------------
-
-      !! Routine arguments
-      INTEGER,                             INTENT(IN) :: nprofs
-      CHARACTER(len=*),DIMENSION(nprofs),  INTENT(IN) :: wmo_inst_list      
-      CHARACTER(len=128),DIMENSION(nprofs),INTENT(OUT):: obs_names
-
-      !! Local variables
-      INTEGER :: inam !: loop over name
-
-      !! Initialise obs_names
-      obs_names(:) = ''
-
-      !! Convert number to text
-      DO inam = 1,nprofs
-         IF (trim(adjustl(wmo_inst_list(inam))) .EQ. '820') THEN
-            obs_names(inam) = 'BUOYS'
-         ELSEIF (trim(adjustl(wmo_inst_list(inam))) .EQ. '401') THEN
-            obs_names(inam) = 'XBT'
-         ELSEIF (trim(adjustl(wmo_inst_list(inam))) .EQ. '741') THEN
-            obs_names(inam) = 'TESAC'
-         ELSEIF (trim(adjustl(wmo_inst_list(inam))) .EQ. '831') THEN
-            obs_names(inam) = 'ARGO'
-         ELSEIF (trim(adjustl(wmo_inst_list(inam))) .EQ. '22') THEN ! Special case for AATSR
-            obs_names(inam) = '22'
-         ELSEIF (trim(adjustl(wmo_inst_list(inam))) .EQ. '') THEN
-            obs_names(inam) = 'UNKNOWN'
-         ELSE 
-            obs_names(inam) = adjustl(wmo_inst_list(inam)) ! For all other cases
-         ENDIF
-      ENDDO
-
-   END SUBROUTINE inst_converter
-
-END MODULE off_netcdf
Index: branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/off_oper.F90
===================================================================
--- branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/off_oper.F90	(revision 4102)
+++ 	(revision )
@@ -1,98 +1,0 @@
-MODULE off_oper
-   !!----------------------------------------------------------------------
-   !!                    ***  MODULE off_oper ***
-   !! ** Purpose : Run NEMO observation operator in offline mode
-   !!----------------------------------------------------------------------
-   !! NEMO modules
-   USE diaobs
-   USE nemogcm
-   USE lib_mpp
-   USE in_out_manager
-   USE obs_fbm, ONLY: ln_cl4
-   !! Offline obs_oper modules
-   USE off_data
-   USE off_read
-   CONTAINS
-      SUBROUTINE off_oper_init
-         !!----------------------------------------------------------------------
-         !!                    ***  SUBROUTINE off_oper_init ***
-         !!
-         !! ** Purpose : To initialise the model as if it were running online.
-         !!
-         !!----------------------------------------------------------------------
-         IMPLICIT NONE
-         !! Initialise NEMO
-         CALL nemo_init
-         !! Initialise Offline obs_oper data
-         CALL off_data_init( ln_cl4 )
-      END SUBROUTINE off_oper_init
-
-      SUBROUTINE off_obs_oper
-         !!----------------------------------------------------------------------
-         !!                    ***  SUBROUTINE off_obs_oper ***
-         !!
-         !! ** Purpose : To use NEMO components to interpolate model fields
-         !!              to observation space.
-         !!
-         !!----------------------------------------------------------------------
-         IMPLICIT NONE
-         !! Initialise offline obs_oper
-         CALL off_oper_init
-         !! Loop over various model counterparts
-         DO jimatch = 1, cl4_match_len
-            IF (jimatch .GT. 1) THEN
-               !! Initialise obs_oper
-               CALL dia_obs_init
-            END IF
-            !! Interpolate to observation space
-            CALL off_oper_interp
-            !! Pipe to output files
-            CALL dia_obs_wri
-            !! Reset the obs_oper
-            CALL dia_obs_dealloc
-         END DO
-         !! Safely stop MPI
-         CALL off_oper_stop
-      END SUBROUTINE off_obs_oper
-
-      SUBROUTINE off_oper_interp
-         !!----------------------------------------------------------------------
-         !!                    ***  SUBROUTINE off_oper_interp ***
-         !!
-         !! ** Purpose : To interpolate the model as if it were running online.
-         !!
-         !! ** Method : 1. Populate model counterparts
-         !!             2. Call dia_obs at appropriate time steps
-         !!----------------------------------------------------------------------
-         IMPLICIT NONE
-         INTEGER :: istp ! time step index
-         !! Loop over entire run
-         istp = nit000
-         nstop = 0
-         DO WHILE ( istp <= nitend .AND. nstop == 0 )
-            IF (jifile <= n_files + 1) THEN
-               IF ( MOD(istp, nn_off_freq) == nit000 ) THEN
-                  !! Read next model counterpart
-                  CALL off_rea_dri(jifile)
-                  jifile = jifile + 1
-               ENDIF
-               !! Interpolate single time step
-               CALL dia_obs(istp)
-            ENDIF
-            !! Increment model step
-            istp = istp + 1
-         END DO
-      END SUBROUTINE off_oper_interp
-
-      SUBROUTINE off_oper_stop
-         !!----------------------------------------------------------------------
-         !!                    ***  SUBROUTINE off_oper_stop ***
-         !!
-         !! ** Purpose : To finalise the model as if it were running online.
-         !!
-         !!----------------------------------------------------------------------
-         IMPLICIT NONE
-         IF(lk_mpp) CALL mppstop  ! end mpp communications
-      END SUBROUTINE off_oper_stop
-
-END MODULE off_oper
Index: branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/off_read.F90
===================================================================
--- branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/off_read.F90	(revision 4102)
+++ 	(revision )
@@ -1,351 +1,0 @@
-
-MODULE off_read
-   !!======================================================================
-   !!                      *** MODULE off_read ***
-   !! Read routines : I/O for offline obs_oper
-   !!======================================================================
-
-USE mppini
-USE lib_mpp
-USE in_out_manager
-USE par_kind, ONLY: lc
-USE netcdf
-USE oce,     ONLY: tsn, sshn
-USE dom_oce, ONLY: nlci, nlcj, nimpp, njmpp, tmask
-USE par_oce, ONLY: jpi, jpj, jpk
-USE obs_fbm, ONLY: fbimdi, fbrmdi, fbsp, fbdp
-
-USE off_data
-!! * Routine accessibility
-PRIVATE
-
-PUBLIC off_rea_dri
-
-CONTAINS
-   SUBROUTINE off_rea_dri(kfile)
-      IMPLICIT NONE
-      !!------------------------------------------------------------------------
-      !!             *** off_rea_dri ***
-      !!
-      !! Purpose : To choose appropriate read method
-      !! Method  : 
-      !!
-      !! Author  : A. Ryan Oct 2013
-      !!
-      !!------------------------------------------------------------------------
-      INTEGER, INTENT(IN) :: &
-              & kfile         !: File number
-      CHARACTER(len=lc) :: &
-              & cdfilename, & !: File name
-              & cmatchname    !: Match name
-      INTEGER :: &
-              & kindex       !: File index to read
- 
-      !! Filename, index and match-up kind
-      cdfilename = TRIM(off_files(kfile))
-      cmatchname = TRIM(cl4_vars(kfile))
-      kindex = nn_off_idx(kfile)
-
-      !! Update model fields
-      !! Class 4 variables: forecast, persistence,
-      !!                    nrt_analysis, best_estimate
-      !! Feedback variables: empty string
-      IF ( (TRIM(cmatchname) == 'forecast') .OR. &
-         & (TRIM(cmatchname) == 'persistence') .OR. &
-         & (TRIM(cmatchname) == 'nrt_analysis') .OR. &
-         & (TRIM(cmatchname) == 'best_estimate').OR. &
-         & (TRIM(cmatchname) == '') ) THEN
-         CALL off_read_file(TRIM(cdfilename), kindex)
-         CALL off_read_juld(TRIM(cdfilename), kindex, cl4_modjuld)
-      ELSE IF (TRIM(cmatchname) == 'climatology') THEN
-         WRITE(numout,*) 'Interpolating climatologies'
-      ELSE IF (TRIM(cmatchname) == 'altimeter') THEN
-         CALL off_read_altbias(TRIM(cdfilename))
-         CALL off_read_juld(TRIM(cdfilename), kindex, cl4_modjuld)
-      END IF
-
-   END SUBROUTINE off_rea_dri
-
-   SUBROUTINE off_read_altbias(filename)
-      IMPLICIT NONE
-      !!------------------------------------------------------------------------
-      !!                      *** off_read_altbias ***
-      !!
-      !! Purpose : To read altimeter bias and set tn,sn to missing values
-      !! Method  : Use subdomain indices to create start and count matrices
-      !!           for netcdf read.
-      !!
-      !! Author  : A. Ryan Sept 2012
-      !!------------------------------------------------------------------------
-      CHARACTER(len=*), INTENT(IN) :: filename
-      INTEGER                      :: ncid, &
-                                    & varid,&
-                                    & istat,&
-                                    & ntimes,&
-                                    & tdim, &
-                                    & xdim, &
-                                    & ydim, &
-                                    & zdim
-      INTEGER                      :: ii, ij, ik
-      INTEGER, DIMENSION(3)        :: start_s, &
-                                    & count_s
-      REAL(fbdp), DIMENSION(:,:),  ALLOCATABLE :: temp_sshn
-      REAL(fbdp)                     :: fill_val
-
-      IF (TRIM(filename) == 'nofile') THEN
-         tsn(:,:,:,:) = fbrmdi
-         sshn(:,:) = fbrmdi 
-      ELSE
-         ! Open Netcdf file to find dimension id
-         istat = nf90_open(trim(filename),nf90_nowrite,ncid)
-         istat = nf90_inq_dimid(ncid,'x',xdim)
-         istat = nf90_inq_dimid(ncid,'y',ydim)
-         istat = nf90_inq_dimid(ncid,'deptht',zdim)
-         istat = nf90_inq_dimid(ncid,'time',tdim)
-         istat = nf90_inquire_dimension(ncid, tdim, len=ntimes)
-         ! Allocate temporary temperature array
-         ALLOCATE(temp_sshn(nlci,nlcj))
-         ! Create start and count arrays
-         start_s = (/ nimpp, njmpp, 1 /)
-         count_s = (/ nlci,  nlcj,  1 /)
-          
-         ! Altimeter bias
-         istat = nf90_inq_varid(ncid,'altbias',varid)
-         istat = nf90_get_att(ncid, varid, '_FillValue', fill_val)
-         istat = nf90_get_var(ncid, varid, temp_sshn, start_s, count_s)
-         WHERE(temp_sshn(:,:) == fill_val) temp_sshn(:,:) = fbrmdi
-   
-         ! Initialise tsn, sshn to fbrmdi
-         tsn(:,:,:,:) = fbrmdi
-         sshn(:,:) = fbrmdi 
-
-         ! Fill sshn with altimeter bias 
-         sshn(1:nlci,1:nlcj) = temp_sshn(:,:) * tmask(1:nlci,1:nlcj,1)
-
-         ! Remove halo from tmask, sshn to prevent double obs counting
-         IF (jpi > nlci) THEN
-             tmask(nlci+1:,:,:) = 0
-             sshn(nlci+1:,:) = 0
-         END IF
-         IF (jpj > nlcj) THEN
-             tmask(:,nlcj+1:,:) = 0
-             sshn(:,nlcj+1:) = 0
-         END IF
-
-         ! Deallocate arrays
-         DEALLOCATE(temp_sshn)
-
-         ! Close netcdf file
-         istat = nf90_close(ncid)
-      END IF
-   
-   END SUBROUTINE off_read_altbias
-
-   SUBROUTINE off_read_file(filename, ifcst)
-      IMPLICIT NONE
-      !!------------------------------------------------------------------------
-      !!             *** off_read_file ***
-      !!
-      !! Purpose : To fill tn and sn with dailymean field from netcdf files
-      !! Method  : Use subdomain indices to create start and count matrices
-      !!           for netcdf read.
-      !!
-      !! Author  : A. Ryan Oct 2010
-      !!------------------------------------------------------------------------
-
-      INTEGER,          INTENT(IN) :: ifcst
-      CHARACTER(len=*), INTENT(IN) :: filename
-      INTEGER                      :: ncid, &
-                                    & varid,&
-                                    & istat,&
-                                    & ntimes,&
-                                    & tdim, &
-                                    & xdim, &
-                                    & ydim, &
-                                    & zdim
-      INTEGER                      :: ii, ij, ik
-      INTEGER, DIMENSION(4)        :: start_n, &
-                                    & count_n
-      INTEGER, DIMENSION(3)        :: start_s, &
-                                    & count_s
-      REAL(fbdp), DIMENSION(:,:,:),ALLOCATABLE :: temp_tn, &
-                                              & temp_sn
-      REAL(fbdp), DIMENSION(:,:),  ALLOCATABLE :: temp_sshn
-      REAL(fbdp)                     :: fill_val
-
-      ! DEBUG
-      INTEGER :: istage
-
-      IF (TRIM(filename) == 'nofile') THEN
-         tsn(:,:,:,:) = fbrmdi
-         sshn(:,:) = fbrmdi 
-      ELSE
-         WRITE(*,*) "Opening :", TRIM(filename)
-         ! Open Netcdf file to find dimension id
-         istat = nf90_open(path=TRIM(filename), mode=nf90_nowrite, ncid=ncid)
-         IF ( istat /= nf90_noerr ) THEN
-             WRITE(*,*) "WARNING: Could not open ", trim(filename)
-             WRITE(*,*) "ERROR: ", nf90_strerror(istat)
-         ENDIF
-         istat = nf90_inq_dimid(ncid,'x',xdim)
-         istat = nf90_inq_dimid(ncid,'y',ydim)
-         istat = nf90_inq_dimid(ncid,'deptht',zdim)
-         istat = nf90_inq_dimid(ncid,'time_counter',tdim)
-         istat = nf90_inquire_dimension(ncid, tdim, len=ntimes)
-         IF (ifcst .LE. ntimes) THEN   
-            ! Allocate temporary temperature array
-            ALLOCATE(temp_tn(nlci,nlcj,jpk))
-            ALLOCATE(temp_sn(nlci,nlcj,jpk))
-            ALLOCATE(temp_sshn(nlci,nlcj))
-      
-            ! Set temp_tn, temp_sn to 0.
-            temp_tn(:,:,:) = fbrmdi
-            temp_sn(:,:,:) = fbrmdi
-            temp_sshn(:,:) = fbrmdi
-      
-            ! Create start and count arrays
-            start_n = (/ nimpp, njmpp, 1,   ifcst /)
-            count_n = (/ nlci,  nlcj,  jpk, 1     /)
-            start_s = (/ nimpp, njmpp, ifcst /)
-            count_s = (/ nlci,  nlcj,  1     /)
-             
-            ! Read information into temporary arrays
-            ! retrieve varid and read in temperature
-            istat = nf90_inq_varid(ncid,'votemper',varid)
-            istat = nf90_get_att(ncid, varid, '_FillValue', fill_val)
-            istat = nf90_get_var(ncid, varid, temp_tn, start_n, count_n)
-            WHERE(temp_tn(:,:,:) == fill_val) temp_tn(:,:,:) = fbrmdi
-      
-            ! retrieve varid and read in salinity
-            istat = nf90_inq_varid(ncid,'vosaline',varid)
-            istat = nf90_get_att(ncid, varid, '_FillValue', fill_val)
-            istat = nf90_get_var(ncid, varid, temp_sn, start_n, count_n)
-            WHERE(temp_sn(:,:,:) == fill_val) temp_sn(:,:,:) = fbrmdi
-      
-            ! retrieve varid and read in SSH
-            istat = nf90_inq_varid(ncid,'sossheig',varid)
-            IF (istat /= nf90_noerr) THEN
-               ! Altimeter bias
-               istat = nf90_inq_varid(ncid,'altbias',varid)
-            END IF
-      
-            istat = nf90_get_att(ncid, varid, '_FillValue', fill_val)
-            istat = nf90_get_var(ncid, varid, temp_sshn, start_s, count_s)
-            WHERE(temp_sshn(:,:) == fill_val) temp_sshn(:,:) = fbrmdi
-   
-            ! Initialise tsn, sshn to fbrmdi
-            tsn(:,:,:,:) = fbrmdi
-            sshn(:,:) = fbrmdi 
-
-            ! Mask out missing data index
-            tsn(1:nlci,1:nlcj,1:jpk,1) = temp_tn(:,:,:) * tmask(1:nlci,1:nlcj,1:jpk)
-            tsn(1:nlci,1:nlcj,1:jpk,2) = temp_sn(:,:,:) * tmask(1:nlci,1:nlcj,1:jpk)
-            sshn(1:nlci,1:nlcj)        = temp_sshn(:,:) * tmask(1:nlci,1:nlcj,1)
-
-            ! Remove halo from tmask, tsn, sshn to prevent double obs counting
-            IF (jpi > nlci) THEN
-                tmask(nlci+1:,:,:) = 0
-                tsn(nlci+1:,:,:,1) = 0
-                tsn(nlci+1:,:,:,2) = 0
-                sshn(nlci+1:,:) = 0
-            END IF
-            IF (jpj > nlcj) THEN
-                tmask(:,nlcj+1:,:) = 0
-                tsn(:,nlcj+1:,:,1) = 0
-                tsn(:,nlcj+1:,:,2) = 0
-                sshn(:,nlcj+1:) = 0
-            END IF
-
-            ! Deallocate arrays
-            DEALLOCATE(temp_tn, temp_sn, temp_sshn)
-         ELSE   
-            ! Mark all as missing data
-            tsn(:,:,:,:) = fbrmdi
-            sshn(:,:) = fbrmdi 
-         ENDIF
-         ! Close netcdf file
-         WRITE(*,*) "Closing :", TRIM(filename)
-         istat = nf90_close(ncid)
-      END IF
-   END SUBROUTINE off_read_file
-
-   SUBROUTINE off_read_juld(filename, ifcst, julian)
-      USE calendar
-      IMPLICIT NONE
-      !!--------------------------------------------------------------------
-      !!                 *** off_read_juld ***
-      !!
-      !!   Purpose : To read model julian day information from file
-      !!   Author  : A. Ryan Nov 2010
-      !!--------------------------------------------------------------------
-
-      !! Routine arguments
-      CHARACTER(len=*), INTENT(IN)  :: filename
-      INTEGER,          INTENT(IN)  :: ifcst
-      REAL,             INTENT(OUT) :: julian    !: Julian day
-
-      !! Local variables
-      INTEGER :: year,  &   !: Date information
-               & month, &
-               & day,   &
-               & hour,  &
-               & minute,&
-               & second
-      INTEGER :: istat, &   !: Netcdf variables
-               & ncid,  &
-               & dimid, &
-               & varid, &
-               & ntime      
-      REAL,DIMENSION(:),ALLOCATABLE :: r_sec     !: Remainder seconds
-      CHARACTER(len=120) :: time_str  !: time string
-
-      time_str=''
-      !! Read in time_counter and remainder seconds
-      istat = nf90_open(trim(filename),nf90_nowrite,ncid)
-      istat = nf90_inq_dimid(ncid,'time_counter',dimid)
-      IF (istat /= nf90_noerr) THEN
-          istat = nf90_inq_dimid(ncid,'time',dimid)
-      ENDIF
-      istat = nf90_inquire_dimension(ncid,dimid,len=ntime)
-      istat = nf90_inq_varid(ncid,'time_counter',varid)
-      IF (istat /= nf90_noerr) THEN
-          istat = nf90_inq_dimid(ncid,'time',dimid)
-      ENDIF
-      istat = nf90_get_att(ncid,varid,'units',time_str)
-      ALLOCATE(r_sec(ntime))
-      istat = nf90_get_var(ncid,varid, r_sec)
-      istat = nf90_close(ncid)
-
-      !! Fill yyyy-mm-dd hh:mm:ss
-      !! format(('seconds since ', I4.4,'-',I2.2,'-',I2.2,' ',I2.2,':',I2.2,':',I2.2))
-      100 format((14x, I4.4,1x,I2.2,1x,I2.2,1x,I2.2,1x,I2.2,1x,I2.2))
-      READ( time_str, 100 ) year, month, day, hour, minute, second
-
-      CALL ymds2ju(year, month, day, r_sec(ifcst), julian)
-
-      !! To take a comment from the ymds2ju subroutine
-
-   !- In the case of the Gregorian calendar we have chosen to use
-   !- the Lilian day numbers. This is the day counter which starts
-   !- on the 15th October 1582.
-   !- This is the day at which Pope Gregory XIII introduced the
-   !- Gregorian calendar.
-   !- Compared to the true Julian calendar, which starts some
-   !- 7980 years ago, the Lilian days are smaler and are dealt with
-   !- easily on 32 bit machines. With the true Julian days you can only
-   !- the fraction of the day in the real part to a precision of
-   !- a 1/4 of a day with 32 bits.
-      
-      !! The obs operator routines prefer to calculate Julian days since 
-      !! 01/01/1950 00:00:00
-      !! In order to convert to the 1950 version we must adjust by the number
-      !! of days between 15th October 1582 and 1st Jan 1950
-
-      julian = julian - 134123.
-      
-      DEALLOCATE(r_sec)
-      
-   END SUBROUTINE off_read_juld
-
-END MODULE off_read 
-
Index: branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/off_write.F90
===================================================================
--- branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/off_write.F90	(revision 4102)
+++ 	(revision )
@@ -1,618 +1,0 @@
-MODULE off_write
-   !!======================================================================
-   !!                       ***  MODULE off_write  ***
-   !!======================================================================
-   USE netcdf
-   USE obs_utils, ONLY: chkerr
-   USE off_netcdf, ONLY: date_format, inst_converter, yyyymmdd_to_ref_date
-   USE off_data
-   IMPLICIT NONE
-   PRIVATE
-
-   PUBLIC off_wri_init
-   PUBLIC off_wri_default
-   PUBLIC off_wri_extra
-
-   ! Type kinds for class 4 data.
-   INTEGER, PARAMETER :: clsp = SELECTED_REAL_KIND( 6, 37) !: single precision
-   INTEGER, PARAMETER :: cldp = SELECTED_REAL_KIND(12,307) !: double precision
-
-   ! Missinge data indicators 
-   INTEGER, PARAMETER    :: climdi = -99999   !: Integers
-   REAL(clsp), PARAMETER :: clrmdi =  99999   !: Reals
-
-   INTERFACE off_wri_extra
-      MODULE PROCEDURE off_wri_extra_3d_index, off_wri_extra_4d, off_wri_extra_4d_index
-   END INTERFACE
-
-   CONTAINS
-
-      SUBROUTINE off_wri_extra_3d_index(cdfilename, cdvarname, ndeps, nvars, &
-                               &  nobs, kstart, kcount, pdata)
-         !!----------------------------------------------------------------------
-         !!                    ***  ROUTINE off_wri_extra_3d  ***
-         !!
-         !! ** Purpose : Write 3d variables to class 4 file.
-         !!
-         !!----------------------------------------------------------------------
-         INTEGER, INTENT(IN) :: &
-                 & nobs, &       !: number of observations/profiles
-                 & nvars, &      !: number of physical parameters
-                 & ndeps         !: number of depths
-         CHARACTER(len=*), INTENT(IN) :: &
-                 & cdfilename, & !: netcdf file name
-                 & cdvarname     !: netcdf variable name
-         INTEGER, DIMENSION(3), INTENT(IN) :: &
-                 & kstart, &     !: start indices
-                 & kcount        !: count indices
-         REAL(KIND=cldp), DIMENSION(ndeps, nvars, nobs), INTENT(IN) :: &
-                 & pdata         !: 3d data
-         INTEGER :: &
-                 & ncid, &       !:
-                 & dimid, &      !:
-                 & varid         !:
-         CHARACTER(len=16), PARAMETER :: cpname = 'off_wri_extra_3d'
-         ! Open netcdf file
-         CALL chkerr(nf90_open(trim(cdfilename), nf90_write, ncid), cpname, __LINE__ )
-         ! Write data
-         CALL chkerr(nf90_inq_varid(ncid,TRIM(cdvarname), varid),cpname, __LINE__ )
-         CALL chkerr(nf90_put_var(ncid, varid, pdata, kstart, kcount),cpname, __LINE__ )
-         ! Close netcdf file
-         CALL chkerr(nf90_close(ncid), cpname, __LINE__ )
-      END SUBROUTINE off_wri_extra_3d_index
-
-      SUBROUTINE off_wri_extra_4d_index(cdfilename, cdvarname, ndeps, nfcst, &
-                               &  nvars, nobs, kstart, kcount, pdata)
-         !!----------------------------------------------------------------------
-         !!                    ***  ROUTINE off_wri_extra_4d  ***
-         !!
-         !! ** Purpose : Write 4d variables to class 4 file.
-         !!
-         !!----------------------------------------------------------------------
-         INTEGER, INTENT(IN) :: &
-                 & nobs, &       !: number of observations/profiles
-                 & nvars, &      !: number of physical parameters
-                 & ndeps, &      !: number of depths
-                 & nfcst         !: number of forecasts
-         CHARACTER(len=*), INTENT(IN) :: &
-                 & cdfilename, & !: netcdf file name
-                 & cdvarname     !: netcdf variable name
-         INTEGER, DIMENSION(4), INTENT(IN) :: &
-                 & kstart, &     !: start indices
-                 & kcount        !: count indices
-         REAL(KIND=cldp), DIMENSION(ndeps, nvars, nobs), INTENT(IN) :: &
-                 & pdata         !: slice of 4d data
-         INTEGER :: &
-                 & ncid, &       !:
-                 & dimid, &      !:
-                 & varid         !:
-         CHARACTER(len=22), PARAMETER :: cpname = 'off_wri_extra_4d_index'
-         ! Open netcdf file
-         CALL chkerr(nf90_open(trim(cdfilename), nf90_write, ncid), cpname, __LINE__ )
-         ! Write data
-         CALL chkerr(nf90_inq_varid(ncid,TRIM(cdvarname), varid),cpname, __LINE__ )
-         CALL chkerr(nf90_put_var(ncid, varid, pdata, kstart, kcount),cpname, __LINE__ )
-         ! Close netcdf file
-         CALL chkerr(nf90_close(ncid), cpname, __LINE__ )
-      END SUBROUTINE off_wri_extra_4d_index
-
-      SUBROUTINE off_wri_extra_4d(cdfilename, cdvarname, ndeps, nfcst, &
-                               &  nvars, nobs, pdata)
-         !!----------------------------------------------------------------------
-         !!                    ***  ROUTINE off_wri_extra_4d  ***
-         !!
-         !! ** Purpose : Write 4d variables to class 4 file.
-         !!
-         !!----------------------------------------------------------------------
-         INTEGER, INTENT(IN) :: &
-                 & nobs, &       !: number of observations/profiles
-                 & nvars, &      !: number of physical parameters
-                 & ndeps, &      !: number of depths
-                 & nfcst         !: number of forecasts
-         CHARACTER(len=*), INTENT(IN) :: &
-                 & cdfilename, & !: netcdf file name
-                 & cdvarname     !: netcdf variable name
-         REAL(KIND=cldp), DIMENSION(ndeps, nfcst, nvars, nobs), INTENT(IN) :: &
-                 & pdata         !: 4d data
-         INTEGER :: &
-                 & ncid, &       !:
-                 & dimid, &      !:
-                 & varid         !:
-         CHARACTER(len=16), PARAMETER :: cpname = 'off_wri_extra_4d'
-         ! Open netcdf file
-         CALL chkerr(nf90_open(trim(cdfilename), nf90_write, ncid), cpname, __LINE__ )
-         ! Write data
-         CALL chkerr(nf90_inq_varid(ncid,TRIM(cdvarname), varid),cpname, __LINE__ )
-         CALL chkerr(nf90_put_var(ncid, varid, pdata),cpname, __LINE__ )
-         ! Close netcdf file
-         CALL chkerr(nf90_close(ncid), cpname, __LINE__ )
-      END SUBROUTINE off_wri_extra_4d
-
-      SUBROUTINE off_wri_default(cdfilename, nobs, nvars, nfcst, ndeps, & 
-                             &   cdtyp, cdwmo, cunit, cvnam, &
-                             &   plam, pphi, pdep, ptim, pob, plead, &
-                             &   kqc, pmjuld)
-         !!----------------------------------------------------------------------
-         !!                    ***  ROUTINE off_wri_default  ***
-         !!
-         !! ** Purpose : Write standard variables to class 4 file.
-         !! ** Method  : Write the following variables;
-         !!                 observation type
-         !!                 observation wmo id
-         !!                 variable units
-         !!                 variable names
-         !!                 observation latitudes
-         !!                 observation longitudes
-         !!                 observation depths
-         !!                 observation times
-         !!                 observation values
-         !!                 observation qc flags
-         !!                 model julian days
-         !!                 
-         !! ** Returns : 
-         !!----------------------------------------------------------------------
-         INTEGER, INTENT(IN) :: &
-                 & nfcst, &      !: number of forecasts
-                 & nobs, &       !: number of observations/profiles
-                 & nvars, &      !: number of physical parameters
-                 & ndeps         !: number of depths
-         CHARACTER(len=*), INTENT(IN) :: &
-                 & cdfilename    !: file name
-         CHARACTER(LEN=128), DIMENSION(nobs), INTENT(IN) :: &
-                 & cdtyp         !: Instrument type
-         CHARACTER(LEN=8), DIMENSION(nobs), INTENT(IN) :: &
-                 & cdwmo         ! WMO number
-         CHARACTER(LEN=128), DIMENSION(nobs) :: &
-                 & cdpad         !: Instrument type (padded)
-         CHARACTER(LEN=8), DIMENSION(nvars), INTENT(IN) :: &
-                 & cunit, &      !: variable units (e.g. psu)
-                 & cvnam         !: variable names (e.g. vosaline)
-         REAL(KIND=cldp), DIMENSION(nfcst), INTENT(IN) :: &
-                 & plead         !: Leadtime
-         REAL(KIND=cldp), DIMENSION(nobs), INTENT(IN) :: &
-                 & pphi, &       !: Latitude
-                 & plam          !: Longitude
-         REAL(KIND=cldp), DIMENSION(ndeps, nobs), INTENT(IN) :: &
-                 & pdep          !: Depth               
-         REAL(KIND=cldp), DIMENSION(ndeps, nvars, nobs), INTENT(IN) :: &
-                 & pob           !: T observation
-         REAL(KIND=cldp), DIMENSION(nobs), INTENT(IN) :: &
-                 & ptim          !: Time
-         INTEGER(KIND=2), DIMENSION(ndeps, nvars, nobs), INTENT(IN) :: &
-                 & kqc           !: Observation QC
-         REAL, INTENT(IN) :: &
-                 & pmjuld        !: Model Julian Day
-         CHARACTER(len=128) :: &
-                 & cvar          !: variable name placeholder
-         INTEGER :: &
-                 & ncid, &  !:
-                 & dimid, & !:
-                 & varid    !:
-         CHARACTER(len=15), PARAMETER :: cpname = 'off_wri_default'
-         ! Open netcdf file
-         CALL chkerr(nf90_open(trim(cdfilename), nf90_write, ncid), cpname, __LINE__ )
-         ! longitude
-         cvar = 'longitude'
-         WRITE(*,*) TRIM(cvar)
-         CALL chkerr(nf90_inq_varid(ncid,TRIM(cvar), varid),cpname, __LINE__ )
-         CALL chkerr(nf90_put_var(ncid, varid, plam),cpname, __LINE__ )
-         ! latitude
-         cvar = 'latitude'
-         WRITE(*,*) TRIM(cvar)
-         CALL chkerr(nf90_inq_varid(ncid,TRIM(cvar), varid),cpname, __LINE__ )
-         CALL chkerr(nf90_put_var(ncid, varid, pphi),cpname, __LINE__ )
-         ! depth
-         cvar = 'depth'
-         WRITE(*,*) TRIM(cvar)
-         CALL chkerr(nf90_inq_varid(ncid,TRIM(cvar), varid),cpname, __LINE__ )
-         CALL chkerr(nf90_put_var(ncid, varid, pdep),cpname, __LINE__ )
-         ! varname
-         cvar = 'varname'
-         WRITE(*,*) TRIM(cvar)
-         CALL chkerr(nf90_inq_varid(ncid,TRIM(cvar), varid),cpname, __LINE__ )
-         CALL chkerr(nf90_put_var(ncid, varid, cvnam),cpname, __LINE__ )
-         ! unitname
-         cvar = 'unitname'
-         WRITE(*,*) TRIM(cvar)
-         CALL chkerr(nf90_inq_varid(ncid,TRIM(cvar), varid),cpname, __LINE__ )
-         CALL chkerr(nf90_put_var(ncid, varid, cunit),cpname, __LINE__ )
-         ! leadtime
-         cvar = 'leadtime'
-         WRITE(*,*) TRIM(cvar)
-         CALL chkerr(nf90_inq_varid(ncid,TRIM(cvar), varid),cpname, __LINE__ )
-         CALL chkerr(nf90_put_var(ncid, varid, plead),cpname, __LINE__ )
-         ! observation
-         cvar = 'observation'
-         WRITE(*,*) TRIM(cvar)
-         CALL chkerr(nf90_inq_varid(ncid,TRIM(cvar), varid),cpname, __LINE__ )
-         CALL chkerr(nf90_put_var(ncid, varid, pob),cpname, __LINE__ )
-         ! qc
-         cvar = 'qc'
-         WRITE(*,*) TRIM(cvar)
-         CALL chkerr(nf90_inq_varid(ncid,TRIM(cvar), varid),cpname, __LINE__ )
-         CALL chkerr(nf90_put_var(ncid, varid, kqc),cpname, __LINE__ )
-         ! juld
-         cvar = 'juld'
-         WRITE(*,*) TRIM(cvar)
-         CALL chkerr(nf90_inq_varid(ncid,TRIM(cvar), varid),cpname, __LINE__ )
-         CALL chkerr(nf90_put_var(ncid, varid, ptim),cpname, __LINE__ )
-         ! type
-         cvar = 'type'
-         WRITE(*,*) TRIM(cvar)
-         CALL inst_converter(cdtyp,nobs,cdpad)
-         CALL chkerr(nf90_inq_varid(ncid,TRIM(cvar), varid),cpname, __LINE__ )
-         CALL chkerr(nf90_put_var(ncid, varid, cdpad),cpname, __LINE__ )
-         ! id
-         cvar = 'id'
-         WRITE(*,*) TRIM(cvar)
-         CALL chkerr(nf90_inq_varid(ncid,TRIM(cvar), varid),cpname, __LINE__ )
-         CALL chkerr(nf90_put_var(ncid, varid, cdwmo),cpname, __LINE__ )
-         ! Close netcdf file
-         CALL chkerr(nf90_close(ncid), cpname, __LINE__ )
-      END SUBROUTINE off_wri_default
-
-      SUBROUTINE off_wri_init(cconf, csys, ckind, cversion, ccont, &
-                            & cinst, cdate, nproc, nobs, nvars, &
-                            & ndeps, nfcst, cdfilename)
-         !!----------------------------------------------------------------------
-         !!                    ***  ROUTINE off_wri_init  ***
-         !!
-         !! ** Purpose : Initialise a class 4 file.
-         !!
-         !! ** Method  :
-         !! ** Returns : cdfilename
-         !!----------------------------------------------------------------------
-         CHARACTER(len=*), INTENT(IN) :: &
-                 & cconf, &      !: model configuration e.g. orca025
-                 & csys, &       !: model system e.g. FOAM
-                 & ckind, &      !: observation kind e.g. profile
-                 & cversion, &   !: model version e.g. 12.0
-                 & ccont, &      !: contact email
-                 & cinst, &      !: institution description
-                 & cdate         !: e.g. yyyymmdd
-         INTEGER, INTENT(IN) :: &
-                 & nproc, &      !: processor number
-                 & nobs, &       !: number of observations/profiles
-                 & nvars, &      !: number of physical parameters
-                 & ndeps, &      !: number of depths
-                 & nfcst         !: number of forecast days
-         CHARACTER(len=128), INTENT(OUT) :: &
-                 & cdfilename    !: e.g. orca025
-         CHARACTER(len=4) :: &
-                 & cproc         !: character representation of nproc
-         CHARACTER(len=128) :: &
-                 & cattnam, &    !: attribute placeholder
-                 & cattval, &    !: attribute placeholder
-                 & cvar, &       !: variable name placeholder
-                 & cqcdes, &     !: quality control description
-                 & cqcval        !: quality control values
-         INTEGER :: &
-                 & istat, &      !: status of netcdf operation
-                 & ncid, &       !: netcdf file id placeholder
-                 & dimid_d, &    !: netcdf dimension id placeholder
-                 & dimid_f, &    !: netcdf dimension id placeholder
-                 & dimid_o, &    !: netcdf dimension id placeholder
-                 & dimid_v, &    !: netcdf dimension id placeholder
-                 & dimid_s, &    !: netcdf dimension id placeholder
-                 & dimid_l, &    !: netcdf dimension id placeholder
-                 & varid         !: netcdf variable id placeholder
-         INTEGER, DIMENSION(2) :: &
-                 & qc_flags      !: quality control flags
-         CHARACTER(len=23) :: &
-                 & ref_date      !: reference date holder
-         INTEGER, DIMENSION(2) :: &
-                 & dim2a, &      !: 2 dimensional settings
-                 & dim2b, &      !:
-                 & dim2c, &      !:
-                 & dim2d         !:
-         INTEGER, DIMENSION(3) :: &
-                 & dim3a         !: 3 dimensional settings
-         INTEGER, DIMENSION(4) :: &
-                 & dim4a         !: 4 dimensional settings
-         CHARACTER(len=12), PARAMETER :: cpname = 'off_wri_init'
-         ! Global att variables
-         CHARACTER(len=40) :: date_str
-         CHARACTER(len=6)  :: noon
-
-         WRITE(cproc, FMT="(I4.4)") nproc
-         cdfilename = 'class4'
-         cdfilename = TRIM(cdfilename)//'_'//TRIM(cdate)
-         cdfilename = TRIM(cdfilename)//'_'//TRIM(csys)
-         cdfilename = TRIM(cdfilename)//'_'//TRIM(cconf)
-         cdfilename = TRIM(cdfilename)//'_'//TRIM(cversion)
-         cdfilename = TRIM(cdfilename)//'_'//TRIM(ckind)
-         cdfilename = TRIM(cdfilename)//'_'//TRIM(cproc)
-         cdfilename = TRIM(cdfilename)//'.nc'
-         ! QC attribute string settings
-         IF (ckind .EQ. 'profile') THEN
-            cqcval = '0 - good data. 9 - bad data.'
-            cqcdes = ''
-            qc_flags = (/0, 9/)
-         ELSE IF (ckind .EQ. 'SST') THEN
-            cqcval = '0 - good data. 1 - bad data.'
-            cqcdes = 'In situ qc flag set to 0 if prob. of gross error < 0.5, &
-                     &1 otherwise. AATSR qc flag set to 0 for best_quality, 1 otherwise.'
-            qc_flags = (/0, 1/)
-         ELSE IF (ckind .EQ. 'SLA') THEN
-            cqcval = '0 - good data. 9 - bad data.'
-            cqcdes = ''
-            qc_flags = (/0, 9/)
-         ELSE
-            cqcval = ''
-            cqcdes = ''
-            qc_flags = (/0, 9/)
-         END IF
-         ! Check cdfilename exists
-         istat = nf90_open(TRIM(cdfilename), nf90_nowrite, ncid)
-         IF (istat /= nf90_noerr) THEN
-            WRITE(*, *) TRIM(cdfilename)
-            CALL chkerr( nf90_create(TRIM(cdfilename), nf90_clobber, ncid), cpname, __LINE__ )
-            !! Global Dimension section
-            !  Partially hardwired for now
-            CALL date_format(date_str)
-            noon = '120000' ! hhmmss
-            ref_date = ''
-            CALL yyyymmdd_to_ref_date(TRIM(cdate), noon, ref_date)
-
-            CALL chkerr( nf90_put_att(ncid, nf90_global, 'title', &
-               & 'Forecast class 4 file'), cpname, __LINE__ )
-            CALL chkerr( nf90_put_att(ncid, nf90_global, 'version', &
-               & TRIM(cversion)), cpname, __LINE__ )
-            CALL chkerr( nf90_put_att(ncid, nf90_global, 'creation_date', &
-               & TRIM(date_str) ), cpname, __LINE__ )
-            CALL chkerr( nf90_put_att(ncid, nf90_global, 'validity_time', &
-               & TRIM(ref_date) ), cpname, __LINE__ )
-            CALL chkerr( nf90_put_att(ncid, nf90_global, 'contact', &
-               & TRIM(ccont) ), cpname, __LINE__ )
-            CALL chkerr( nf90_put_att(ncid, nf90_global, 'obs_type', &
-               & TRIM(ckind) ), cpname, __LINE__ )
-            CALL chkerr( nf90_put_att(ncid, nf90_global, 'system', &
-               & TRIM(csys)), cpname, __LINE__ )
-            CALL chkerr( nf90_put_att(ncid, nf90_global, 'configuration', &
-               & TRIM(cconf) ), cpname, __LINE__ ) 
-            CALL chkerr( nf90_put_att(ncid, nf90_global, 'institution', &
-               & TRIM(cinst) ), cpname, __LINE__ ) 
-
-            !! Define Dimensions
-            CALL chkerr(nf90_def_dim(ncid, 'numdeps', ndeps, dimid_d), cpname, __LINE__ )
-            CALL chkerr(nf90_def_dim(ncid, 'numfcsts', nfcst, dimid_f), cpname, __LINE__ )
-            CALL chkerr(nf90_def_dim(ncid, 'numvars', nvars, dimid_v), cpname, __LINE__ )
-            CALL chkerr(nf90_def_dim(ncid, 'numobs', nobs, dimid_o), cpname, __LINE__ )
-            CALL chkerr(nf90_def_dim(ncid, 'string_length8', 8, dimid_s), cpname, __LINE__ )
-            CALL chkerr(nf90_def_dim(ncid, 'string_length128', 128, dimid_l), cpname, __LINE__ )
-            !! Define possible dimension permutations
-            ! 2d
-            dim2a(:) = (/ dimid_d, dimid_o /) !: (/ ndeps, nobs  /) 
-            dim2b(:) = (/ dimid_s, dimid_o /) !: (/ 8, nobs  /)
-            dim2c(:) = (/ dimid_s, dimid_v /) !: (/ 8, nvars /)
-            dim2d(:) = (/ dimid_l, dimid_o /) !: (/ 128, nobs  /)
-            ! 3d
-            dim3a(:) = (/ dimid_d, dimid_v, dimid_o/) !: (/ ndeps, nvars, nobs /)
-            ! 4d
-            dim4a(:) = (/ dimid_d, dimid_f, dimid_v, dimid_o /) !: (/ ndeps, nfcst, nvars, nobs /)
-
-            !! Create the variables
-            !  Forecast day
-            cvar = 'leadtime'
-            WRITE(*,*) TRIM(cvar)
-            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dimid_f, varid), cpname, __LINE__ )
-            cattnam = 'long_name'
-            cattval = 'Model forecast day offset'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-            cattnam = 'units'
-            cattval = 'Hours'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-            cattnam = 'comment'
-            cattval = 'Hours between forecast production and validity time'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-
-            !  longitude
-            cvar = 'longitude'
-            WRITE(*,*) TRIM(cvar)
-            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dimid_o, varid), cpname, __LINE__ )
-            cattnam = 'long_name'
-            cattval = 'Longitudes'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-            cattnam = 'units'
-            cattval = 'Degrees'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-
-            !  latitude
-            cvar = 'latitude'
-            WRITE(*,*) TRIM(cvar)
-            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dimid_o, varid), cpname, __LINE__ )
-            cattnam = 'long_name'
-            cattval = 'Latitudes'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-            cattnam = 'units'
-            cattval = 'Degrees'
-
-            !  depth
-            cvar = 'depth'
-            WRITE(*,*) TRIM(cvar)
-            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dim2a, varid), cpname, __LINE__ )
-            CALL chkerr(nf90_put_att(ncid, varid, '_FillValue', clrmdi), cpname, __LINE__ )
-            cattnam = 'long_name'
-            cattval = 'Depths'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-            cattnam = 'units'
-            cattval = 'metre'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-
-            !  varname
-            cvar = 'varname'
-            WRITE(*,*) TRIM(cvar)
-            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_char, dim2c, varid), cpname, __LINE__ )
-            cattnam = 'long_name'
-            cattval = 'Variable name'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-
-            !  unitname
-            cvar = 'unitname'
-            WRITE(*,*) TRIM(cvar)
-            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_char, dim2c, varid), cpname, __LINE__ )
-            cattnam = 'long_name'
-            cattval = 'Unit name'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-
-            !  observation
-            cvar = 'observation'
-            WRITE(*,*) TRIM(cvar)
-            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dim3a, varid), cpname, __LINE__ )
-            CALL chkerr(nf90_put_att(ncid, varid, '_FillValue', clrmdi), cpname, __LINE__ )
-            cattnam = 'long_name'
-            cattval = 'Observation value'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-
-            !  forecast
-            cvar = 'forecast'
-            WRITE(*,*) TRIM(cvar)
-            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dim4a, varid), cpname, __LINE__ )
-            CALL chkerr(nf90_put_att(ncid, varid, '_FillValue', clrmdi), cpname, __LINE__ )
-            cattnam = 'long_name'
-            cattval = 'Model forecast counterpart of obs. value'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-            cattnam = 'comment'
-            cattval = 'Model daily mean valid at noon used for calculation'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-
-            !  persistence
-            cvar = 'persistence'
-            WRITE(*,*) TRIM(cvar)
-            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dim4a, varid), cpname, __LINE__ )
-            CALL chkerr(nf90_put_att(ncid, varid, '_FillValue', clrmdi), cpname, __LINE__ )
-            cattnam = 'long_name'
-            cattval = 'Model persistence counterpart of obs. value'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-            cattnam = 'comment'
-            cattval = 'Model daily mean valid at noon used for calculation'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-
-            !  climatology
-            cvar = 'climatology'
-            WRITE(*,*) TRIM(cvar)
-            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dim3a, varid), cpname, __LINE__ )
-            CALL chkerr(nf90_put_att(ncid, varid, '_FillValue', clrmdi), cpname, __LINE__ )
-            cattnam = 'long_name'
-            cattval = 'Climatological value'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-            cattnam = 'comment'
-            cattval = 'Levitus monthly fields interpolated to the correct day'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-
-            !  daym2
-            cvar = 'best_estimate'
-            WRITE(*,*) TRIM(cvar)
-            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dim3a, varid), cpname, __LINE__ )
-            CALL chkerr(nf90_put_att(ncid, varid, '_FillValue', clrmdi), cpname, __LINE__ )
-            cattnam = 'long_name'
-            cattval = 'Best estimate'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-            cattnam = 'comment'
-            cattval = 'FOAM daym2 field'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-
-            !  daym1
-            cvar = 'nrt_analysis'
-            WRITE(*,*) TRIM(cvar)
-            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dim3a, varid), cpname, __LINE__ )
-            CALL chkerr(nf90_put_att(ncid, varid, '_FillValue', clrmdi), cpname, __LINE__ )
-            cattnam = 'long_name'
-            cattval = 'Near real time analysis'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-            cattnam = 'comment'
-            cattval = 'FOAM daym1 field'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-         
-            !  optional fields
-            IF (TRIM(ckind) .EQ. 'SLA') THEN
-               !  mdt
-               cvar = 'mdt_reference'
-               WRITE(*,*) TRIM(cvar)
-               CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dim3a, varid), cpname, __LINE__ )
-               CALL chkerr(nf90_put_att(ncid, varid, '_FillValue', clrmdi), cpname, __LINE__ )
-               cattnam = 'long_name'
-               cattval = 'Mean dynamic topography'
-               CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-               cattnam = 'comment'
-               cattval = 'MDT reference field'
-               CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-
-               !  altbias
-               cvar = 'altimeter_bias'
-               WRITE(*,*) TRIM(cvar)
-               CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dim3a, varid), cpname, __LINE__ )
-               CALL chkerr(nf90_put_att(ncid, varid, '_FillValue', clrmdi), cpname, __LINE__ )
-               cattnam = 'long_name'
-               cattval = 'Altimeter bias'
-               CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-            ENDIF
-         
-            !  qc
-            cvar = 'qc'
-            WRITE(*,*) TRIM(cvar)
-            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_short, dim3a, varid), cpname, __LINE__ )
-            CALL chkerr(nf90_put_att(ncid, varid, '_FillValue', NF90_FILL_SHORT), cpname, __LINE__ )
-            cattnam = 'long_name'
-            cattval = 'Quality flags'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-            cattnam = 'flag_value'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), qc_flags), cpname, __LINE__ )
-            cattnam = 'flag_meaning'
-            cattval = cqcval
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-            cattnam = 'comment'
-            cattval = cqcdes
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-
-!            !  juld
-            cvar = 'juld'
-            WRITE(*,*) TRIM(cvar)
-            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_double, dimid_o, varid), cpname, __LINE__ )
-!            CALL chkerr(nf90_put_att(ncid, varid, '_FillValue', clrmdi), cpname, __LINE__ )
-            cattnam = 'long_name'
-            cattval = 'Observation time in Julian days'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-            cattnam = 'units'
-            CALL yyyymmdd_to_ref_date('19500101', '000000', ref_date)
-            cattval = 'Days since '//ref_date
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-
-            !  modeljuld 
-            cvar = 'modeljuld'
-            WRITE(*,*) TRIM(cvar)
-            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dimid_f, varid), cpname, __LINE__ )
-            cattnam = 'long_name'
-            cattval = 'Model field date in Julian days'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-            cattnam = 'units'
-            cattval = 'Days since '//ref_date
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-
-            !  type
-            cvar = 'type'
-            WRITE(*,*) TRIM(cvar)
-            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_char, dim2d, varid), cpname, __LINE__ )
-            cattnam = 'long_name'
-            cattval = 'Observation type'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-
-            !  id
-            cvar = 'id'
-            WRITE(*,*) TRIM(cvar)
-            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_char, dim2b, varid), cpname, __LINE__ )
-            cattnam = 'long_name'
-            cattval = 'Observation id'
-            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
-
-            ! Close netcdf file
-            CALL chkerr(nf90_close(ncid), cpname, __LINE__ )
-         ELSE
-            WRITE(*, *) cdfilename, 'already exists.'
-            ! Close netcdf file
-            CALL chkerr(nf90_close(ncid), cpname, __LINE__ )
-         END IF
-      END SUBROUTINE off_wri_init
-
-
-END MODULE off_write
Index: branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/ooo_data.F90
===================================================================
--- branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/ooo_data.F90	(revision 4106)
+++ branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/ooo_data.F90	(revision 4106)
@@ -0,0 +1,136 @@
+MODULE ooo_data
+   USE par_kind, ONLY: lc
+   IMPLICIT NONE
+   INTEGER, PARAMETER :: MaxNumFiles = 1000
+
+   !! Class 4 file settings
+   INTEGER :: &
+           & cl4_fcst_idx(MaxNumFiles), & !: forecast indices
+           & cl4_match_len, &             !: number of match types
+           & cl4_fcst_len                 !: number of forecast days
+   CHARACTER(len=lc) :: &
+           & cl4_vars(MaxNumFiles), &  !: class 4 variables
+           & cl4_sys, &                !: class 4 system
+           & cl4_cfg, &                !: class 4 configuration
+           & cl4_date, &               !: class 4 date
+           & cl4_vn,  &                !: class 4 version
+           & cl4_prefix, &             !: class 4 prefix
+           & cl4_contact, &            !: class 4 contact
+           & cl4_inst                  !: class 4 institute
+   REAL ::   cl4_modjuld               !: model Julian day
+   REAL :: &
+      & cl4_leadtime(MaxNumFiles)      !: Lead time data
+
+   !! Offline obs_oper settings
+   CHARACTER(len=lc) :: &
+      & off_files(MaxNumFiles)         !: model files
+   INTEGER            :: &
+      & jifile, &                      !: current file list index
+      & n_files, &                     !: number of files
+      & jimatch, &                     !: current match
+      & nn_off_idx(MaxNumFiles), &     !: time_counter indices
+      & nn_off_freq                    !: read frequency in time steps
+   CHARACTER(len=128) :: &
+      & alt_file                       !: altimeter file
+CONTAINS
+   SUBROUTINE ooo_data_init( ld_cl4 )
+      !!----------------------------------------------------------------------
+      !!                    ***  SUBROUTINE ooo_data_init ***
+      !!
+      !! ** Purpose : To read namelists and initialise offline_oper run.
+      !!
+      !!----------------------------------------------------------------------
+      USE in_out_manager
+      IMPLICIT NONE
+      INTEGER            :: &
+         & jf                           !: file dummy loop index
+      LOGICAL :: lmask(MaxNumFiles)     !: Logical mask used for counting
+      LOGICAL, INTENT(IN) :: ld_cl4     !: Logical class 4 on/off
+
+      ! Standard offline obs_oper information
+      NAMELIST/namoff/off_files, nn_off_idx, nn_off_freq
+
+      ! Class 4 file specifiers
+      NAMELIST/namcl4/cl4_vars, cl4_sys, cl4_cfg, cl4_date, cl4_vn, &
+         &            cl4_prefix, cl4_contact, cl4_inst, cl4_leadtime, &
+         &            cl4_fcst_idx, cl4_fcst_len, cl4_match_len
+
+      ! Standard offline obs_oper initialisation
+      jimatch = 0                   !: match-up iteration variable 
+      jifile = 1                    !: input file iteration variable 
+      n_files = 0                   !: number of files to cycle through
+      off_files(:) = ''             !: list of files to read in
+      nn_off_idx(:) = 0             !: list of indices inside each file
+      nn_off_freq = -1              !: input frequency in time steps
+
+      ! Class 4 initialisation
+      cl4_leadtime(:) = 0           !: Lead time axis value for each file
+      cl4_fcst_len = 0              !: Length of the forecast dimension
+      cl4_match_len = 1             !: Number of match types
+      cl4_fcst_idx(:) = 0           !: output file forecast index
+      cl4_vars(:) = ''              !: output file variable names
+      cl4_sys = ''                  !: output file system
+      cl4_cfg = ''                  !: output file configuration
+      cl4_date = ''                 !: output file date string
+      cl4_vn = ''                   !: output file version
+      cl4_prefix = 'class4'         !: output file prefix
+      cl4_contact = ''              !: output file contact details
+      cl4_inst = ''                 !: output file institution
+
+      ! Standard offline obs_oper settings
+      READ(numnam, namoff)
+
+      ! Read class 4 output settings
+      IF (ld_cl4) THEN
+         READ(numnam, namcl4)
+      ENDIF
+
+      ! count input files
+      lmask(:) = .FALSE.
+      WHERE (off_files(:) /= '') lmask(:) = .TRUE.
+      n_files = COUNT(lmask)
+
+      !! Initialise sub obs window frequency
+      IF (nn_off_freq == -1) THEN
+         !! Run length
+         nn_off_freq = nitend - nit000 + 1
+      ENDIF
+
+      !! Print summary of settings
+      IF(lwp) THEN
+         WRITE(numout,*)
+         WRITE(numout,*) 'offline obs_oper : Initialization'
+         WRITE(numout,*) '~~~~~~~~~~~~~~~~~'
+         WRITE(numout,*) '   Namelist namoff : set offline obs_oper parameters' 
+         DO jf = 1, n_files
+            WRITE(numout,'(1X,2A)') '   Input forecast file name          forecastfile = ', &
+               TRIM(off_files(jf))
+            WRITE(numout,*) '   Input forecast file index        forecastindex = ', &
+               nn_off_idx(jf)
+            WRITE(numout,*) '   Output forecast leadtime index   leadtimeindex = ', &
+               cl4_fcst_idx(jf)
+            WRITE(numout,*) '   Output forecast leadtime value   leadtimevalue = ', &
+               cl4_leadtime(jf)
+            WRITE(numout,'(1X,2A)') '   Input class 4 variable       class 4 parameter = ', &
+               TRIM(cl4_vars(jf))
+         END DO
+         WRITE(numout, '(1X,2A)') '   Input class 4 system            class 4 system = ', &
+            TRIM(cl4_sys)
+         WRITE(numout, '(1X,2A)') '   Input class 4 config            class 4 config = ', &
+            TRIM(cl4_cfg)
+         WRITE(numout, '(1X,2A)') '   Input class 4 date                class 4 date = ', &
+            TRIM(cl4_date)
+         WRITE(numout, '(1X,2A)') '   Input class 4 version          class 4 version = ', &
+            TRIM(cl4_vn)
+         WRITE(numout, '(1X,2A)') '   Input class 4 prefix            class 4 prefix = ', &
+            TRIM(cl4_prefix)
+         WRITE(numout, '(1X,2A)') '   Input class 4 contact          class 4 contact = ', &
+            TRIM(cl4_contact)
+         WRITE(numout, '(1X,2A)') '   Input class 4 institute      class 4 institute = ', &
+            TRIM(cl4_inst)
+      END IF
+
+   END SUBROUTINE ooo_data_init
+
+END MODULE ooo_data
+
Index: branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/ooo_netcdf.F90
===================================================================
--- branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/ooo_netcdf.F90	(revision 4106)
+++ branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/ooo_netcdf.F90	(revision 4106)
@@ -0,0 +1,104 @@
+MODULE ooo_netcdf
+   USE obs_utils, ONLY: chkerr
+   USE obs_const, ONLY: obfillflt
+   USE par_oce
+   USE obs_types, ONLY: cwmonam1770
+   USE obs_sla, ONLY: nslasets
+
+   ! Define double precision obfillflt
+   REAL(kind=dp), PARAMETER :: obfilldbl=99999.
+
+   CONTAINS
+
+   SUBROUTINE yyyymmdd_to_ref_date(indate, intime, outstring)
+      IMPLICIT NONE
+      !----------------------------------------
+      ! Routine to create reference date string
+      !----------------------------------------
+      ! Routine arguments
+      CHARACTER(len=8),  INTENT(IN)  :: indate    ! yyyymmdd
+      CHARACTER(len=6),  INTENT(IN)  :: intime    ! hhmmss
+      CHARACTER(len=23), INTENT(OUT) :: outstring ! yyyy-mm-dd hh:mm:ss utc
+      ! Local variables
+      CHARACTER(len=4)  :: year      ! yyyy
+      CHARACTER(len=2)  :: month     ! mm
+      CHARACTER(len=2)  :: day       ! dd
+      CHARACTER(len=2)  :: hour      ! hh
+      CHARACTER(len=2)  :: minute    ! mm
+      CHARACTER(len=2)  :: second    ! ss
+
+      year   = indate(1:4)
+      month  = indate(5:6)
+      day    = indate(7:8)
+      hour   = intime(1:2)
+      minute = intime(3:4)
+      second = intime(5:6)
+
+      outstring = year//"-"//month//"-"//day//" "//hour//":"//minute//":"//second//" utc"
+      
+   END SUBROUTINE yyyymmdd_to_ref_date
+
+   SUBROUTINE date_format(date_str)
+      IMPLICIT NONE
+      !---------------------------------------
+      ! Routine to create creation date string
+      !---------------------------------------
+
+      ! Routine arguments
+      CHARACTER(len=*), INTENT(OUT) :: date_str
+      ! Local variables
+      CHARACTER(8)  :: date
+      CHARACTER(10) :: time
+
+      CHARACTER(10)  :: date_part
+      CHARACTER(8)  :: time_part
+
+      CALL date_and_time(DATE=date)
+      CALL date_and_time(TIME=time)
+      date_part  = date(1:4)//'/'//date(5:6)//'/'//date(7:8)
+      time_part  = time(1:2)//':'//time(3:4)//':'//time(5:6)
+
+      date_str   = date_part//' '//time_part
+
+   END SUBROUTINE date_format
+
+   SUBROUTINE inst_converter(wmo_inst_list, nprofs, obs_names)
+      IMPLICIT NONE
+      !---------------------------------------
+      ! Routine to convert WMO_INST_TYPE to
+      ! 'ARGO', 'TESAC', 'XBT', 'BUOYS'
+      !---------------------------------------
+
+      !! Routine arguments
+      INTEGER,                             INTENT(IN) :: nprofs
+      CHARACTER(len=*),DIMENSION(nprofs),  INTENT(IN) :: wmo_inst_list      
+      CHARACTER(len=128),DIMENSION(nprofs),INTENT(OUT):: obs_names
+
+      !! Local variables
+      INTEGER :: inam !: loop over name
+
+      !! Initialise obs_names
+      obs_names(:) = ''
+
+      !! Convert number to text
+      DO inam = 1,nprofs
+         IF (trim(adjustl(wmo_inst_list(inam))) .EQ. '820') THEN
+            obs_names(inam) = 'BUOYS'
+         ELSEIF (trim(adjustl(wmo_inst_list(inam))) .EQ. '401') THEN
+            obs_names(inam) = 'XBT'
+         ELSEIF (trim(adjustl(wmo_inst_list(inam))) .EQ. '741') THEN
+            obs_names(inam) = 'TESAC'
+         ELSEIF (trim(adjustl(wmo_inst_list(inam))) .EQ. '831') THEN
+            obs_names(inam) = 'ARGO'
+         ELSEIF (trim(adjustl(wmo_inst_list(inam))) .EQ. '22') THEN ! Special case for AATSR
+            obs_names(inam) = '22'
+         ELSEIF (trim(adjustl(wmo_inst_list(inam))) .EQ. '') THEN
+            obs_names(inam) = 'UNKNOWN'
+         ELSE 
+            obs_names(inam) = adjustl(wmo_inst_list(inam)) ! For all other cases
+         ENDIF
+      ENDDO
+
+   END SUBROUTINE inst_converter
+
+END MODULE ooo_netcdf
Index: branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/ooo_read.F90
===================================================================
--- branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/ooo_read.F90	(revision 4106)
+++ branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/ooo_read.F90	(revision 4106)
@@ -0,0 +1,351 @@
+
+MODULE ooo_read
+   !!======================================================================
+   !!                      *** MODULE ooo_read ***
+   !! Read routines : I/O for offline obs_oper
+   !!======================================================================
+
+USE mppini
+USE lib_mpp
+USE in_out_manager
+USE par_kind, ONLY: lc
+USE netcdf
+USE oce,     ONLY: tsn, sshn
+USE dom_oce, ONLY: nlci, nlcj, nimpp, njmpp, tmask
+USE par_oce, ONLY: jpi, jpj, jpk
+USE obs_fbm, ONLY: fbimdi, fbrmdi, fbsp, fbdp
+
+USE ooo_data
+!! * Routine accessibility
+PRIVATE
+
+PUBLIC ooo_rea_dri
+
+CONTAINS
+   SUBROUTINE ooo_rea_dri(kfile)
+      IMPLICIT NONE
+      !!------------------------------------------------------------------------
+      !!             *** ooo_rea_dri ***
+      !!
+      !! Purpose : To choose appropriate read method
+      !! Method  : 
+      !!
+      !! Author  : A. Ryan Oct 2013
+      !!
+      !!------------------------------------------------------------------------
+      INTEGER, INTENT(IN) :: &
+              & kfile         !: File number
+      CHARACTER(len=lc) :: &
+              & cdfilename, & !: File name
+              & cmatchname    !: Match name
+      INTEGER :: &
+              & kindex       !: File index to read
+ 
+      !! Filename, index and match-up kind
+      cdfilename = TRIM(off_files(kfile))
+      cmatchname = TRIM(cl4_vars(kfile))
+      kindex = nn_off_idx(kfile)
+
+      !! Update model fields
+      !! Class 4 variables: forecast, persistence,
+      !!                    nrt_analysis, best_estimate
+      !! Feedback variables: empty string
+      IF ( (TRIM(cmatchname) == 'forecast') .OR. &
+         & (TRIM(cmatchname) == 'persistence') .OR. &
+         & (TRIM(cmatchname) == 'nrt_analysis') .OR. &
+         & (TRIM(cmatchname) == 'best_estimate').OR. &
+         & (TRIM(cmatchname) == '') ) THEN
+         CALL ooo_read_file(TRIM(cdfilename), kindex)
+         CALL ooo_read_juld(TRIM(cdfilename), kindex, cl4_modjuld)
+      ELSE IF (TRIM(cmatchname) == 'climatology') THEN
+         WRITE(numout,*) 'Interpolating climatologies'
+      ELSE IF (TRIM(cmatchname) == 'altimeter') THEN
+         CALL ooo_read_altbias(TRIM(cdfilename))
+         CALL ooo_read_juld(TRIM(cdfilename), kindex, cl4_modjuld)
+      END IF
+
+   END SUBROUTINE ooo_rea_dri
+
+   SUBROUTINE ooo_read_altbias(filename)
+      IMPLICIT NONE
+      !!------------------------------------------------------------------------
+      !!                      *** ooo_read_altbias ***
+      !!
+      !! Purpose : To read altimeter bias and set tn,sn to missing values
+      !! Method  : Use subdomain indices to create start and count matrices
+      !!           for netcdf read.
+      !!
+      !! Author  : A. Ryan Sept 2012
+      !!------------------------------------------------------------------------
+      CHARACTER(len=*), INTENT(IN) :: filename
+      INTEGER                      :: ncid, &
+                                    & varid,&
+                                    & istat,&
+                                    & ntimes,&
+                                    & tdim, &
+                                    & xdim, &
+                                    & ydim, &
+                                    & zdim
+      INTEGER                      :: ii, ij, ik
+      INTEGER, DIMENSION(3)        :: start_s, &
+                                    & count_s
+      REAL(fbdp), DIMENSION(:,:),  ALLOCATABLE :: temp_sshn
+      REAL(fbdp)                     :: fill_val
+
+      IF (TRIM(filename) == 'nofile') THEN
+         tsn(:,:,:,:) = fbrmdi
+         sshn(:,:) = fbrmdi 
+      ELSE
+         ! Open Netcdf file to find dimension id
+         istat = nf90_open(trim(filename),nf90_nowrite,ncid)
+         istat = nf90_inq_dimid(ncid,'x',xdim)
+         istat = nf90_inq_dimid(ncid,'y',ydim)
+         istat = nf90_inq_dimid(ncid,'deptht',zdim)
+         istat = nf90_inq_dimid(ncid,'time',tdim)
+         istat = nf90_inquire_dimension(ncid, tdim, len=ntimes)
+         ! Allocate temporary temperature array
+         ALLOCATE(temp_sshn(nlci,nlcj))
+         ! Create start and count arrays
+         start_s = (/ nimpp, njmpp, 1 /)
+         count_s = (/ nlci,  nlcj,  1 /)
+          
+         ! Altimeter bias
+         istat = nf90_inq_varid(ncid,'altbias',varid)
+         istat = nf90_get_att(ncid, varid, '_FillValue', fill_val)
+         istat = nf90_get_var(ncid, varid, temp_sshn, start_s, count_s)
+         WHERE(temp_sshn(:,:) == fill_val) temp_sshn(:,:) = fbrmdi
+   
+         ! Initialise tsn, sshn to fbrmdi
+         tsn(:,:,:,:) = fbrmdi
+         sshn(:,:) = fbrmdi 
+
+         ! Fill sshn with altimeter bias 
+         sshn(1:nlci,1:nlcj) = temp_sshn(:,:) * tmask(1:nlci,1:nlcj,1)
+
+         ! Remove halo from tmask, sshn to prevent double obs counting
+         IF (jpi > nlci) THEN
+             tmask(nlci+1:,:,:) = 0
+             sshn(nlci+1:,:) = 0
+         END IF
+         IF (jpj > nlcj) THEN
+             tmask(:,nlcj+1:,:) = 0
+             sshn(:,nlcj+1:) = 0
+         END IF
+
+         ! Deallocate arrays
+         DEALLOCATE(temp_sshn)
+
+         ! Close netcdf file
+         istat = nf90_close(ncid)
+      END IF
+   
+   END SUBROUTINE ooo_read_altbias
+
+   SUBROUTINE ooo_read_file(filename, ifcst)
+      IMPLICIT NONE
+      !!------------------------------------------------------------------------
+      !!             *** ooo_read_file ***
+      !!
+      !! Purpose : To fill tn and sn with dailymean field from netcdf files
+      !! Method  : Use subdomain indices to create start and count matrices
+      !!           for netcdf read.
+      !!
+      !! Author  : A. Ryan Oct 2010
+      !!------------------------------------------------------------------------
+
+      INTEGER,          INTENT(IN) :: ifcst
+      CHARACTER(len=*), INTENT(IN) :: filename
+      INTEGER                      :: ncid, &
+                                    & varid,&
+                                    & istat,&
+                                    & ntimes,&
+                                    & tdim, &
+                                    & xdim, &
+                                    & ydim, &
+                                    & zdim
+      INTEGER                      :: ii, ij, ik
+      INTEGER, DIMENSION(4)        :: start_n, &
+                                    & count_n
+      INTEGER, DIMENSION(3)        :: start_s, &
+                                    & count_s
+      REAL(fbdp), DIMENSION(:,:,:),ALLOCATABLE :: temp_tn, &
+                                              & temp_sn
+      REAL(fbdp), DIMENSION(:,:),  ALLOCATABLE :: temp_sshn
+      REAL(fbdp)                     :: fill_val
+
+      ! DEBUG
+      INTEGER :: istage
+
+      IF (TRIM(filename) == 'nofile') THEN
+         tsn(:,:,:,:) = fbrmdi
+         sshn(:,:) = fbrmdi 
+      ELSE
+         WRITE(*,*) "Opening :", TRIM(filename)
+         ! Open Netcdf file to find dimension id
+         istat = nf90_open(path=TRIM(filename), mode=nf90_nowrite, ncid=ncid)
+         IF ( istat /= nf90_noerr ) THEN
+             WRITE(*,*) "WARNING: Could not open ", trim(filename)
+             WRITE(*,*) "ERROR: ", nf90_strerror(istat)
+         ENDIF
+         istat = nf90_inq_dimid(ncid,'x',xdim)
+         istat = nf90_inq_dimid(ncid,'y',ydim)
+         istat = nf90_inq_dimid(ncid,'deptht',zdim)
+         istat = nf90_inq_dimid(ncid,'time_counter',tdim)
+         istat = nf90_inquire_dimension(ncid, tdim, len=ntimes)
+         IF (ifcst .LE. ntimes) THEN   
+            ! Allocate temporary temperature array
+            ALLOCATE(temp_tn(nlci,nlcj,jpk))
+            ALLOCATE(temp_sn(nlci,nlcj,jpk))
+            ALLOCATE(temp_sshn(nlci,nlcj))
+      
+            ! Set temp_tn, temp_sn to 0.
+            temp_tn(:,:,:) = fbrmdi
+            temp_sn(:,:,:) = fbrmdi
+            temp_sshn(:,:) = fbrmdi
+      
+            ! Create start and count arrays
+            start_n = (/ nimpp, njmpp, 1,   ifcst /)
+            count_n = (/ nlci,  nlcj,  jpk, 1     /)
+            start_s = (/ nimpp, njmpp, ifcst /)
+            count_s = (/ nlci,  nlcj,  1     /)
+             
+            ! Read information into temporary arrays
+            ! retrieve varid and read in temperature
+            istat = nf90_inq_varid(ncid,'votemper',varid)
+            istat = nf90_get_att(ncid, varid, '_FillValue', fill_val)
+            istat = nf90_get_var(ncid, varid, temp_tn, start_n, count_n)
+            WHERE(temp_tn(:,:,:) == fill_val) temp_tn(:,:,:) = fbrmdi
+      
+            ! retrieve varid and read in salinity
+            istat = nf90_inq_varid(ncid,'vosaline',varid)
+            istat = nf90_get_att(ncid, varid, '_FillValue', fill_val)
+            istat = nf90_get_var(ncid, varid, temp_sn, start_n, count_n)
+            WHERE(temp_sn(:,:,:) == fill_val) temp_sn(:,:,:) = fbrmdi
+      
+            ! retrieve varid and read in SSH
+            istat = nf90_inq_varid(ncid,'sossheig',varid)
+            IF (istat /= nf90_noerr) THEN
+               ! Altimeter bias
+               istat = nf90_inq_varid(ncid,'altbias',varid)
+            END IF
+      
+            istat = nf90_get_att(ncid, varid, '_FillValue', fill_val)
+            istat = nf90_get_var(ncid, varid, temp_sshn, start_s, count_s)
+            WHERE(temp_sshn(:,:) == fill_val) temp_sshn(:,:) = fbrmdi
+   
+            ! Initialise tsn, sshn to fbrmdi
+            tsn(:,:,:,:) = fbrmdi
+            sshn(:,:) = fbrmdi 
+
+            ! Mask out missing data index
+            tsn(1:nlci,1:nlcj,1:jpk,1) = temp_tn(:,:,:) * tmask(1:nlci,1:nlcj,1:jpk)
+            tsn(1:nlci,1:nlcj,1:jpk,2) = temp_sn(:,:,:) * tmask(1:nlci,1:nlcj,1:jpk)
+            sshn(1:nlci,1:nlcj)        = temp_sshn(:,:) * tmask(1:nlci,1:nlcj,1)
+
+            ! Remove halo from tmask, tsn, sshn to prevent double obs counting
+            IF (jpi > nlci) THEN
+                tmask(nlci+1:,:,:) = 0
+                tsn(nlci+1:,:,:,1) = 0
+                tsn(nlci+1:,:,:,2) = 0
+                sshn(nlci+1:,:) = 0
+            END IF
+            IF (jpj > nlcj) THEN
+                tmask(:,nlcj+1:,:) = 0
+                tsn(:,nlcj+1:,:,1) = 0
+                tsn(:,nlcj+1:,:,2) = 0
+                sshn(:,nlcj+1:) = 0
+            END IF
+
+            ! Deallocate arrays
+            DEALLOCATE(temp_tn, temp_sn, temp_sshn)
+         ELSE   
+            ! Mark all as missing data
+            tsn(:,:,:,:) = fbrmdi
+            sshn(:,:) = fbrmdi 
+         ENDIF
+         ! Close netcdf file
+         WRITE(*,*) "Closing :", TRIM(filename)
+         istat = nf90_close(ncid)
+      END IF
+   END SUBROUTINE ooo_read_file
+
+   SUBROUTINE ooo_read_juld(filename, ifcst, julian)
+      USE calendar
+      IMPLICIT NONE
+      !!--------------------------------------------------------------------
+      !!                 *** ooo_read_juld ***
+      !!
+      !!   Purpose : To read model julian day information from file
+      !!   Author  : A. Ryan Nov 2010
+      !!--------------------------------------------------------------------
+
+      !! Routine arguments
+      CHARACTER(len=*), INTENT(IN)  :: filename
+      INTEGER,          INTENT(IN)  :: ifcst
+      REAL,             INTENT(OUT) :: julian    !: Julian day
+
+      !! Local variables
+      INTEGER :: year,  &   !: Date information
+               & month, &
+               & day,   &
+               & hour,  &
+               & minute,&
+               & second
+      INTEGER :: istat, &   !: Netcdf variables
+               & ncid,  &
+               & dimid, &
+               & varid, &
+               & ntime      
+      REAL,DIMENSION(:),ALLOCATABLE :: r_sec     !: Remainder seconds
+      CHARACTER(len=120) :: time_str  !: time string
+
+      time_str=''
+      !! Read in time_counter and remainder seconds
+      istat = nf90_open(trim(filename),nf90_nowrite,ncid)
+      istat = nf90_inq_dimid(ncid,'time_counter',dimid)
+      IF (istat /= nf90_noerr) THEN
+          istat = nf90_inq_dimid(ncid,'time',dimid)
+      ENDIF
+      istat = nf90_inquire_dimension(ncid,dimid,len=ntime)
+      istat = nf90_inq_varid(ncid,'time_counter',varid)
+      IF (istat /= nf90_noerr) THEN
+          istat = nf90_inq_dimid(ncid,'time',dimid)
+      ENDIF
+      istat = nf90_get_att(ncid,varid,'units',time_str)
+      ALLOCATE(r_sec(ntime))
+      istat = nf90_get_var(ncid,varid, r_sec)
+      istat = nf90_close(ncid)
+
+      !! Fill yyyy-mm-dd hh:mm:ss
+      !! format(('seconds since ', I4.4,'-',I2.2,'-',I2.2,' ',I2.2,':',I2.2,':',I2.2))
+      100 format((14x, I4.4,1x,I2.2,1x,I2.2,1x,I2.2,1x,I2.2,1x,I2.2))
+      READ( time_str, 100 ) year, month, day, hour, minute, second
+
+      CALL ymds2ju(year, month, day, r_sec(ifcst), julian)
+
+      !! To take a comment from the ymds2ju subroutine
+
+   !- In the case of the Gregorian calendar we have chosen to use
+   !- the Lilian day numbers. This is the day counter which starts
+   !- on the 15th October 1582.
+   !- This is the day at which Pope Gregory XIII introduced the
+   !- Gregorian calendar.
+   !- Compared to the true Julian calendar, which starts some
+   !- 7980 years ago, the Lilian days are smaler and are dealt with
+   !- easily on 32 bit machines. With the true Julian days you can only
+   !- the fraction of the day in the real part to a precision of
+   !- a 1/4 of a day with 32 bits.
+      
+      !! The obs operator routines prefer to calculate Julian days since 
+      !! 01/01/1950 00:00:00
+      !! In order to convert to the 1950 version we must adjust by the number
+      !! of days between 15th October 1582 and 1st Jan 1950
+
+      julian = julian - 134123.
+      
+      DEALLOCATE(r_sec)
+      
+   END SUBROUTINE ooo_read_juld
+
+END MODULE ooo_read 
+
Index: branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/ooo_write.F90
===================================================================
--- branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/ooo_write.F90	(revision 4106)
+++ branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/ooo_write.F90	(revision 4106)
@@ -0,0 +1,618 @@
+MODULE ooo_write
+   !!======================================================================
+   !!                       ***  MODULE ooo_write  ***
+   !!======================================================================
+   USE netcdf
+   USE obs_utils, ONLY: chkerr
+   USE ooo_netcdf, ONLY: date_format, inst_converter, yyyymmdd_to_ref_date
+   USE ooo_data
+   IMPLICIT NONE
+   PRIVATE
+
+   PUBLIC ooo_wri_init
+   PUBLIC ooo_wri_default
+   PUBLIC ooo_wri_extra
+
+   ! Type kinds for class 4 data.
+   INTEGER, PARAMETER :: clsp = SELECTED_REAL_KIND( 6, 37) !: single precision
+   INTEGER, PARAMETER :: cldp = SELECTED_REAL_KIND(12,307) !: double precision
+
+   ! Missinge data indicators 
+   INTEGER, PARAMETER    :: climdi = -99999   !: Integers
+   REAL(clsp), PARAMETER :: clrmdi =  99999   !: Reals
+
+   INTERFACE ooo_wri_extra
+      MODULE PROCEDURE ooo_wri_extra_3d_index, ooo_wri_extra_4d, ooo_wri_extra_4d_index
+   END INTERFACE
+
+   CONTAINS
+
+      SUBROUTINE ooo_wri_extra_3d_index(cdfilename, cdvarname, ndeps, nvars, &
+                               &  nobs, kstart, kcount, pdata)
+         !!----------------------------------------------------------------------
+         !!                    ***  ROUTINE ooo_wri_extra_3d  ***
+         !!
+         !! ** Purpose : Write 3d variables to class 4 file.
+         !!
+         !!----------------------------------------------------------------------
+         INTEGER, INTENT(IN) :: &
+                 & nobs, &       !: number of observations/profiles
+                 & nvars, &      !: number of physical parameters
+                 & ndeps         !: number of depths
+         CHARACTER(len=*), INTENT(IN) :: &
+                 & cdfilename, & !: netcdf file name
+                 & cdvarname     !: netcdf variable name
+         INTEGER, DIMENSION(3), INTENT(IN) :: &
+                 & kstart, &     !: start indices
+                 & kcount        !: count indices
+         REAL(KIND=cldp), DIMENSION(ndeps, nvars, nobs), INTENT(IN) :: &
+                 & pdata         !: 3d data
+         INTEGER :: &
+                 & ncid, &       !:
+                 & dimid, &      !:
+                 & varid         !:
+         CHARACTER(len=16), PARAMETER :: cpname = 'ooo_wri_extra_3d'
+         ! Open netcdf file
+         CALL chkerr(nf90_open(trim(cdfilename), nf90_write, ncid), cpname, __LINE__ )
+         ! Write data
+         CALL chkerr(nf90_inq_varid(ncid,TRIM(cdvarname), varid),cpname, __LINE__ )
+         CALL chkerr(nf90_put_var(ncid, varid, pdata, kstart, kcount),cpname, __LINE__ )
+         ! Close netcdf file
+         CALL chkerr(nf90_close(ncid), cpname, __LINE__ )
+      END SUBROUTINE ooo_wri_extra_3d_index
+
+      SUBROUTINE ooo_wri_extra_4d_index(cdfilename, cdvarname, ndeps, nfcst, &
+                               &  nvars, nobs, kstart, kcount, pdata)
+         !!----------------------------------------------------------------------
+         !!                    ***  ROUTINE ooo_wri_extra_4d  ***
+         !!
+         !! ** Purpose : Write 4d variables to class 4 file.
+         !!
+         !!----------------------------------------------------------------------
+         INTEGER, INTENT(IN) :: &
+                 & nobs, &       !: number of observations/profiles
+                 & nvars, &      !: number of physical parameters
+                 & ndeps, &      !: number of depths
+                 & nfcst         !: number of forecasts
+         CHARACTER(len=*), INTENT(IN) :: &
+                 & cdfilename, & !: netcdf file name
+                 & cdvarname     !: netcdf variable name
+         INTEGER, DIMENSION(4), INTENT(IN) :: &
+                 & kstart, &     !: start indices
+                 & kcount        !: count indices
+         REAL(KIND=cldp), DIMENSION(ndeps, nvars, nobs), INTENT(IN) :: &
+                 & pdata         !: slice of 4d data
+         INTEGER :: &
+                 & ncid, &       !:
+                 & dimid, &      !:
+                 & varid         !:
+         CHARACTER(len=22), PARAMETER :: cpname = 'ooo_wri_extra_4d_index'
+         ! Open netcdf file
+         CALL chkerr(nf90_open(trim(cdfilename), nf90_write, ncid), cpname, __LINE__ )
+         ! Write data
+         CALL chkerr(nf90_inq_varid(ncid,TRIM(cdvarname), varid),cpname, __LINE__ )
+         CALL chkerr(nf90_put_var(ncid, varid, pdata, kstart, kcount),cpname, __LINE__ )
+         ! Close netcdf file
+         CALL chkerr(nf90_close(ncid), cpname, __LINE__ )
+      END SUBROUTINE ooo_wri_extra_4d_index
+
+      SUBROUTINE ooo_wri_extra_4d(cdfilename, cdvarname, ndeps, nfcst, &
+                               &  nvars, nobs, pdata)
+         !!----------------------------------------------------------------------
+         !!                    ***  ROUTINE ooo_wri_extra_4d  ***
+         !!
+         !! ** Purpose : Write 4d variables to class 4 file.
+         !!
+         !!----------------------------------------------------------------------
+         INTEGER, INTENT(IN) :: &
+                 & nobs, &       !: number of observations/profiles
+                 & nvars, &      !: number of physical parameters
+                 & ndeps, &      !: number of depths
+                 & nfcst         !: number of forecasts
+         CHARACTER(len=*), INTENT(IN) :: &
+                 & cdfilename, & !: netcdf file name
+                 & cdvarname     !: netcdf variable name
+         REAL(KIND=cldp), DIMENSION(ndeps, nfcst, nvars, nobs), INTENT(IN) :: &
+                 & pdata         !: 4d data
+         INTEGER :: &
+                 & ncid, &       !:
+                 & dimid, &      !:
+                 & varid         !:
+         CHARACTER(len=16), PARAMETER :: cpname = 'ooo_wri_extra_4d'
+         ! Open netcdf file
+         CALL chkerr(nf90_open(trim(cdfilename), nf90_write, ncid), cpname, __LINE__ )
+         ! Write data
+         CALL chkerr(nf90_inq_varid(ncid,TRIM(cdvarname), varid),cpname, __LINE__ )
+         CALL chkerr(nf90_put_var(ncid, varid, pdata),cpname, __LINE__ )
+         ! Close netcdf file
+         CALL chkerr(nf90_close(ncid), cpname, __LINE__ )
+      END SUBROUTINE ooo_wri_extra_4d
+
+      SUBROUTINE ooo_wri_default(cdfilename, nobs, nvars, nfcst, ndeps, & 
+                             &   cdtyp, cdwmo, cunit, cvnam, &
+                             &   plam, pphi, pdep, ptim, pob, plead, &
+                             &   kqc, pmjuld)
+         !!----------------------------------------------------------------------
+         !!                    ***  ROUTINE ooo_wri_default  ***
+         !!
+         !! ** Purpose : Write standard variables to class 4 file.
+         !! ** Method  : Write the following variables;
+         !!                 observation type
+         !!                 observation wmo id
+         !!                 variable units
+         !!                 variable names
+         !!                 observation latitudes
+         !!                 observation longitudes
+         !!                 observation depths
+         !!                 observation times
+         !!                 observation values
+         !!                 observation qc flags
+         !!                 model julian days
+         !!                 
+         !! ** Returns : 
+         !!----------------------------------------------------------------------
+         INTEGER, INTENT(IN) :: &
+                 & nfcst, &      !: number of forecasts
+                 & nobs, &       !: number of observations/profiles
+                 & nvars, &      !: number of physical parameters
+                 & ndeps         !: number of depths
+         CHARACTER(len=*), INTENT(IN) :: &
+                 & cdfilename    !: file name
+         CHARACTER(LEN=128), DIMENSION(nobs), INTENT(IN) :: &
+                 & cdtyp         !: Instrument type
+         CHARACTER(LEN=8), DIMENSION(nobs), INTENT(IN) :: &
+                 & cdwmo         ! WMO number
+         CHARACTER(LEN=128), DIMENSION(nobs) :: &
+                 & cdpad         !: Instrument type (padded)
+         CHARACTER(LEN=8), DIMENSION(nvars), INTENT(IN) :: &
+                 & cunit, &      !: variable units (e.g. psu)
+                 & cvnam         !: variable names (e.g. vosaline)
+         REAL(KIND=cldp), DIMENSION(nfcst), INTENT(IN) :: &
+                 & plead         !: Leadtime
+         REAL(KIND=cldp), DIMENSION(nobs), INTENT(IN) :: &
+                 & pphi, &       !: Latitude
+                 & plam          !: Longitude
+         REAL(KIND=cldp), DIMENSION(ndeps, nobs), INTENT(IN) :: &
+                 & pdep          !: Depth               
+         REAL(KIND=cldp), DIMENSION(ndeps, nvars, nobs), INTENT(IN) :: &
+                 & pob           !: T observation
+         REAL(KIND=cldp), DIMENSION(nobs), INTENT(IN) :: &
+                 & ptim          !: Time
+         INTEGER(KIND=2), DIMENSION(ndeps, nvars, nobs), INTENT(IN) :: &
+                 & kqc           !: Observation QC
+         REAL, INTENT(IN) :: &
+                 & pmjuld        !: Model Julian Day
+         CHARACTER(len=128) :: &
+                 & cvar          !: variable name placeholder
+         INTEGER :: &
+                 & ncid, &  !:
+                 & dimid, & !:
+                 & varid    !:
+         CHARACTER(len=15), PARAMETER :: cpname = 'ooo_wri_default'
+         ! Open netcdf file
+         CALL chkerr(nf90_open(trim(cdfilename), nf90_write, ncid), cpname, __LINE__ )
+         ! longitude
+         cvar = 'longitude'
+         WRITE(*,*) TRIM(cvar)
+         CALL chkerr(nf90_inq_varid(ncid,TRIM(cvar), varid),cpname, __LINE__ )
+         CALL chkerr(nf90_put_var(ncid, varid, plam),cpname, __LINE__ )
+         ! latitude
+         cvar = 'latitude'
+         WRITE(*,*) TRIM(cvar)
+         CALL chkerr(nf90_inq_varid(ncid,TRIM(cvar), varid),cpname, __LINE__ )
+         CALL chkerr(nf90_put_var(ncid, varid, pphi),cpname, __LINE__ )
+         ! depth
+         cvar = 'depth'
+         WRITE(*,*) TRIM(cvar)
+         CALL chkerr(nf90_inq_varid(ncid,TRIM(cvar), varid),cpname, __LINE__ )
+         CALL chkerr(nf90_put_var(ncid, varid, pdep),cpname, __LINE__ )
+         ! varname
+         cvar = 'varname'
+         WRITE(*,*) TRIM(cvar)
+         CALL chkerr(nf90_inq_varid(ncid,TRIM(cvar), varid),cpname, __LINE__ )
+         CALL chkerr(nf90_put_var(ncid, varid, cvnam),cpname, __LINE__ )
+         ! unitname
+         cvar = 'unitname'
+         WRITE(*,*) TRIM(cvar)
+         CALL chkerr(nf90_inq_varid(ncid,TRIM(cvar), varid),cpname, __LINE__ )
+         CALL chkerr(nf90_put_var(ncid, varid, cunit),cpname, __LINE__ )
+         ! leadtime
+         cvar = 'leadtime'
+         WRITE(*,*) TRIM(cvar)
+         CALL chkerr(nf90_inq_varid(ncid,TRIM(cvar), varid),cpname, __LINE__ )
+         CALL chkerr(nf90_put_var(ncid, varid, plead),cpname, __LINE__ )
+         ! observation
+         cvar = 'observation'
+         WRITE(*,*) TRIM(cvar)
+         CALL chkerr(nf90_inq_varid(ncid,TRIM(cvar), varid),cpname, __LINE__ )
+         CALL chkerr(nf90_put_var(ncid, varid, pob),cpname, __LINE__ )
+         ! qc
+         cvar = 'qc'
+         WRITE(*,*) TRIM(cvar)
+         CALL chkerr(nf90_inq_varid(ncid,TRIM(cvar), varid),cpname, __LINE__ )
+         CALL chkerr(nf90_put_var(ncid, varid, kqc),cpname, __LINE__ )
+         ! juld
+         cvar = 'juld'
+         WRITE(*,*) TRIM(cvar)
+         CALL chkerr(nf90_inq_varid(ncid,TRIM(cvar), varid),cpname, __LINE__ )
+         CALL chkerr(nf90_put_var(ncid, varid, ptim),cpname, __LINE__ )
+         ! type
+         cvar = 'type'
+         WRITE(*,*) TRIM(cvar)
+         CALL inst_converter(cdtyp,nobs,cdpad)
+         CALL chkerr(nf90_inq_varid(ncid,TRIM(cvar), varid),cpname, __LINE__ )
+         CALL chkerr(nf90_put_var(ncid, varid, cdpad),cpname, __LINE__ )
+         ! id
+         cvar = 'id'
+         WRITE(*,*) TRIM(cvar)
+         CALL chkerr(nf90_inq_varid(ncid,TRIM(cvar), varid),cpname, __LINE__ )
+         CALL chkerr(nf90_put_var(ncid, varid, cdwmo),cpname, __LINE__ )
+         ! Close netcdf file
+         CALL chkerr(nf90_close(ncid), cpname, __LINE__ )
+      END SUBROUTINE ooo_wri_default
+
+      SUBROUTINE ooo_wri_init(cconf, csys, ckind, cversion, ccont, &
+                            & cinst, cdate, nproc, nobs, nvars, &
+                            & ndeps, nfcst, cdfilename)
+         !!----------------------------------------------------------------------
+         !!                    ***  ROUTINE ooo_wri_init  ***
+         !!
+         !! ** Purpose : Initialise a class 4 file.
+         !!
+         !! ** Method  :
+         !! ** Returns : cdfilename
+         !!----------------------------------------------------------------------
+         CHARACTER(len=*), INTENT(IN) :: &
+                 & cconf, &      !: model configuration e.g. orca025
+                 & csys, &       !: model system e.g. FOAM
+                 & ckind, &      !: observation kind e.g. profile
+                 & cversion, &   !: model version e.g. 12.0
+                 & ccont, &      !: contact email
+                 & cinst, &      !: institution description
+                 & cdate         !: e.g. yyyymmdd
+         INTEGER, INTENT(IN) :: &
+                 & nproc, &      !: processor number
+                 & nobs, &       !: number of observations/profiles
+                 & nvars, &      !: number of physical parameters
+                 & ndeps, &      !: number of depths
+                 & nfcst         !: number of forecast days
+         CHARACTER(len=128), INTENT(OUT) :: &
+                 & cdfilename    !: e.g. orca025
+         CHARACTER(len=4) :: &
+                 & cproc         !: character representation of nproc
+         CHARACTER(len=128) :: &
+                 & cattnam, &    !: attribute placeholder
+                 & cattval, &    !: attribute placeholder
+                 & cvar, &       !: variable name placeholder
+                 & cqcdes, &     !: quality control description
+                 & cqcval        !: quality control values
+         INTEGER :: &
+                 & istat, &      !: status of netcdf operation
+                 & ncid, &       !: netcdf file id placeholder
+                 & dimid_d, &    !: netcdf dimension id placeholder
+                 & dimid_f, &    !: netcdf dimension id placeholder
+                 & dimid_o, &    !: netcdf dimension id placeholder
+                 & dimid_v, &    !: netcdf dimension id placeholder
+                 & dimid_s, &    !: netcdf dimension id placeholder
+                 & dimid_l, &    !: netcdf dimension id placeholder
+                 & varid         !: netcdf variable id placeholder
+         INTEGER, DIMENSION(2) :: &
+                 & qc_flags      !: quality control flags
+         CHARACTER(len=23) :: &
+                 & ref_date      !: reference date holder
+         INTEGER, DIMENSION(2) :: &
+                 & dim2a, &      !: 2 dimensional settings
+                 & dim2b, &      !:
+                 & dim2c, &      !:
+                 & dim2d         !:
+         INTEGER, DIMENSION(3) :: &
+                 & dim3a         !: 3 dimensional settings
+         INTEGER, DIMENSION(4) :: &
+                 & dim4a         !: 4 dimensional settings
+         CHARACTER(len=12), PARAMETER :: cpname = 'ooo_wri_init'
+         ! Global att variables
+         CHARACTER(len=40) :: date_str
+         CHARACTER(len=6)  :: noon
+
+         WRITE(cproc, FMT="(I4.4)") nproc
+         cdfilename = 'class4'
+         cdfilename = TRIM(cdfilename)//'_'//TRIM(cdate)
+         cdfilename = TRIM(cdfilename)//'_'//TRIM(csys)
+         cdfilename = TRIM(cdfilename)//'_'//TRIM(cconf)
+         cdfilename = TRIM(cdfilename)//'_'//TRIM(cversion)
+         cdfilename = TRIM(cdfilename)//'_'//TRIM(ckind)
+         cdfilename = TRIM(cdfilename)//'_'//TRIM(cproc)
+         cdfilename = TRIM(cdfilename)//'.nc'
+         ! QC attribute string settings
+         IF (ckind .EQ. 'profile') THEN
+            cqcval = '0 - good data. 9 - bad data.'
+            cqcdes = ''
+            qc_flags = (/0, 9/)
+         ELSE IF (ckind .EQ. 'SST') THEN
+            cqcval = '0 - good data. 1 - bad data.'
+            cqcdes = 'In situ qc flag set to 0 if prob. of gross error < 0.5, &
+                     &1 otherwise. AATSR qc flag set to 0 for best_quality, 1 otherwise.'
+            qc_flags = (/0, 1/)
+         ELSE IF (ckind .EQ. 'SLA') THEN
+            cqcval = '0 - good data. 9 - bad data.'
+            cqcdes = ''
+            qc_flags = (/0, 9/)
+         ELSE
+            cqcval = ''
+            cqcdes = ''
+            qc_flags = (/0, 9/)
+         END IF
+         ! Check cdfilename exists
+         istat = nf90_open(TRIM(cdfilename), nf90_nowrite, ncid)
+         IF (istat /= nf90_noerr) THEN
+            WRITE(*, *) TRIM(cdfilename)
+            CALL chkerr( nf90_create(TRIM(cdfilename), nf90_clobber, ncid), cpname, __LINE__ )
+            !! Global Dimension section
+            !  Partially hardwired for now
+            CALL date_format(date_str)
+            noon = '120000' ! hhmmss
+            ref_date = ''
+            CALL yyyymmdd_to_ref_date(TRIM(cdate), noon, ref_date)
+
+            CALL chkerr( nf90_put_att(ncid, nf90_global, 'title', &
+               & 'Forecast class 4 file'), cpname, __LINE__ )
+            CALL chkerr( nf90_put_att(ncid, nf90_global, 'version', &
+               & TRIM(cversion)), cpname, __LINE__ )
+            CALL chkerr( nf90_put_att(ncid, nf90_global, 'creation_date', &
+               & TRIM(date_str) ), cpname, __LINE__ )
+            CALL chkerr( nf90_put_att(ncid, nf90_global, 'validity_time', &
+               & TRIM(ref_date) ), cpname, __LINE__ )
+            CALL chkerr( nf90_put_att(ncid, nf90_global, 'contact', &
+               & TRIM(ccont) ), cpname, __LINE__ )
+            CALL chkerr( nf90_put_att(ncid, nf90_global, 'obs_type', &
+               & TRIM(ckind) ), cpname, __LINE__ )
+            CALL chkerr( nf90_put_att(ncid, nf90_global, 'system', &
+               & TRIM(csys)), cpname, __LINE__ )
+            CALL chkerr( nf90_put_att(ncid, nf90_global, 'configuration', &
+               & TRIM(cconf) ), cpname, __LINE__ ) 
+            CALL chkerr( nf90_put_att(ncid, nf90_global, 'institution', &
+               & TRIM(cinst) ), cpname, __LINE__ ) 
+
+            !! Define Dimensions
+            CALL chkerr(nf90_def_dim(ncid, 'numdeps', ndeps, dimid_d), cpname, __LINE__ )
+            CALL chkerr(nf90_def_dim(ncid, 'numfcsts', nfcst, dimid_f), cpname, __LINE__ )
+            CALL chkerr(nf90_def_dim(ncid, 'numvars', nvars, dimid_v), cpname, __LINE__ )
+            CALL chkerr(nf90_def_dim(ncid, 'numobs', nobs, dimid_o), cpname, __LINE__ )
+            CALL chkerr(nf90_def_dim(ncid, 'string_length8', 8, dimid_s), cpname, __LINE__ )
+            CALL chkerr(nf90_def_dim(ncid, 'string_length128', 128, dimid_l), cpname, __LINE__ )
+            !! Define possible dimension permutations
+            ! 2d
+            dim2a(:) = (/ dimid_d, dimid_o /) !: (/ ndeps, nobs  /) 
+            dim2b(:) = (/ dimid_s, dimid_o /) !: (/ 8, nobs  /)
+            dim2c(:) = (/ dimid_s, dimid_v /) !: (/ 8, nvars /)
+            dim2d(:) = (/ dimid_l, dimid_o /) !: (/ 128, nobs  /)
+            ! 3d
+            dim3a(:) = (/ dimid_d, dimid_v, dimid_o/) !: (/ ndeps, nvars, nobs /)
+            ! 4d
+            dim4a(:) = (/ dimid_d, dimid_f, dimid_v, dimid_o /) !: (/ ndeps, nfcst, nvars, nobs /)
+
+            !! Create the variables
+            !  Forecast day
+            cvar = 'leadtime'
+            WRITE(*,*) TRIM(cvar)
+            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dimid_f, varid), cpname, __LINE__ )
+            cattnam = 'long_name'
+            cattval = 'Model forecast day offset'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+            cattnam = 'units'
+            cattval = 'Hours'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+            cattnam = 'comment'
+            cattval = 'Hours between forecast production and validity time'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+
+            !  longitude
+            cvar = 'longitude'
+            WRITE(*,*) TRIM(cvar)
+            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dimid_o, varid), cpname, __LINE__ )
+            cattnam = 'long_name'
+            cattval = 'Longitudes'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+            cattnam = 'units'
+            cattval = 'Degrees'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+
+            !  latitude
+            cvar = 'latitude'
+            WRITE(*,*) TRIM(cvar)
+            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dimid_o, varid), cpname, __LINE__ )
+            cattnam = 'long_name'
+            cattval = 'Latitudes'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+            cattnam = 'units'
+            cattval = 'Degrees'
+
+            !  depth
+            cvar = 'depth'
+            WRITE(*,*) TRIM(cvar)
+            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dim2a, varid), cpname, __LINE__ )
+            CALL chkerr(nf90_put_att(ncid, varid, '_FillValue', clrmdi), cpname, __LINE__ )
+            cattnam = 'long_name'
+            cattval = 'Depths'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+            cattnam = 'units'
+            cattval = 'metre'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+
+            !  varname
+            cvar = 'varname'
+            WRITE(*,*) TRIM(cvar)
+            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_char, dim2c, varid), cpname, __LINE__ )
+            cattnam = 'long_name'
+            cattval = 'Variable name'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+
+            !  unitname
+            cvar = 'unitname'
+            WRITE(*,*) TRIM(cvar)
+            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_char, dim2c, varid), cpname, __LINE__ )
+            cattnam = 'long_name'
+            cattval = 'Unit name'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+
+            !  observation
+            cvar = 'observation'
+            WRITE(*,*) TRIM(cvar)
+            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dim3a, varid), cpname, __LINE__ )
+            CALL chkerr(nf90_put_att(ncid, varid, '_FillValue', clrmdi), cpname, __LINE__ )
+            cattnam = 'long_name'
+            cattval = 'Observation value'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+
+            !  forecast
+            cvar = 'forecast'
+            WRITE(*,*) TRIM(cvar)
+            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dim4a, varid), cpname, __LINE__ )
+            CALL chkerr(nf90_put_att(ncid, varid, '_FillValue', clrmdi), cpname, __LINE__ )
+            cattnam = 'long_name'
+            cattval = 'Model forecast counterpart of obs. value'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+            cattnam = 'comment'
+            cattval = 'Model daily mean valid at noon used for calculation'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+
+            !  persistence
+            cvar = 'persistence'
+            WRITE(*,*) TRIM(cvar)
+            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dim4a, varid), cpname, __LINE__ )
+            CALL chkerr(nf90_put_att(ncid, varid, '_FillValue', clrmdi), cpname, __LINE__ )
+            cattnam = 'long_name'
+            cattval = 'Model persistence counterpart of obs. value'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+            cattnam = 'comment'
+            cattval = 'Model daily mean valid at noon used for calculation'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+
+            !  climatology
+            cvar = 'climatology'
+            WRITE(*,*) TRIM(cvar)
+            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dim3a, varid), cpname, __LINE__ )
+            CALL chkerr(nf90_put_att(ncid, varid, '_FillValue', clrmdi), cpname, __LINE__ )
+            cattnam = 'long_name'
+            cattval = 'Climatological value'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+            cattnam = 'comment'
+            cattval = 'Levitus monthly fields interpolated to the correct day'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+
+            !  daym2
+            cvar = 'best_estimate'
+            WRITE(*,*) TRIM(cvar)
+            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dim3a, varid), cpname, __LINE__ )
+            CALL chkerr(nf90_put_att(ncid, varid, '_FillValue', clrmdi), cpname, __LINE__ )
+            cattnam = 'long_name'
+            cattval = 'Best estimate'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+            cattnam = 'comment'
+            cattval = 'FOAM daym2 field'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+
+            !  daym1
+            cvar = 'nrt_analysis'
+            WRITE(*,*) TRIM(cvar)
+            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dim3a, varid), cpname, __LINE__ )
+            CALL chkerr(nf90_put_att(ncid, varid, '_FillValue', clrmdi), cpname, __LINE__ )
+            cattnam = 'long_name'
+            cattval = 'Near real time analysis'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+            cattnam = 'comment'
+            cattval = 'FOAM daym1 field'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+         
+            !  optional fields
+            IF (TRIM(ckind) .EQ. 'SLA') THEN
+               !  mdt
+               cvar = 'mdt_reference'
+               WRITE(*,*) TRIM(cvar)
+               CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dim3a, varid), cpname, __LINE__ )
+               CALL chkerr(nf90_put_att(ncid, varid, '_FillValue', clrmdi), cpname, __LINE__ )
+               cattnam = 'long_name'
+               cattval = 'Mean dynamic topography'
+               CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+               cattnam = 'comment'
+               cattval = 'MDT reference field'
+               CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+
+               !  altbias
+               cvar = 'altimeter_bias'
+               WRITE(*,*) TRIM(cvar)
+               CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dim3a, varid), cpname, __LINE__ )
+               CALL chkerr(nf90_put_att(ncid, varid, '_FillValue', clrmdi), cpname, __LINE__ )
+               cattnam = 'long_name'
+               cattval = 'Altimeter bias'
+               CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+            ENDIF
+         
+            !  qc
+            cvar = 'qc'
+            WRITE(*,*) TRIM(cvar)
+            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_short, dim3a, varid), cpname, __LINE__ )
+            CALL chkerr(nf90_put_att(ncid, varid, '_FillValue', NF90_FILL_SHORT), cpname, __LINE__ )
+            cattnam = 'long_name'
+            cattval = 'Quality flags'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+            cattnam = 'flag_value'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), qc_flags), cpname, __LINE__ )
+            cattnam = 'flag_meaning'
+            cattval = cqcval
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+            cattnam = 'comment'
+            cattval = cqcdes
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+
+!            !  juld
+            cvar = 'juld'
+            WRITE(*,*) TRIM(cvar)
+            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_double, dimid_o, varid), cpname, __LINE__ )
+!            CALL chkerr(nf90_put_att(ncid, varid, '_FillValue', clrmdi), cpname, __LINE__ )
+            cattnam = 'long_name'
+            cattval = 'Observation time in Julian days'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+            cattnam = 'units'
+            CALL yyyymmdd_to_ref_date('19500101', '000000', ref_date)
+            cattval = 'Days since '//ref_date
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+
+            !  modeljuld 
+            cvar = 'modeljuld'
+            WRITE(*,*) TRIM(cvar)
+            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_float, dimid_f, varid), cpname, __LINE__ )
+            cattnam = 'long_name'
+            cattval = 'Model field date in Julian days'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+            cattnam = 'units'
+            cattval = 'Days since '//ref_date
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+
+            !  type
+            cvar = 'type'
+            WRITE(*,*) TRIM(cvar)
+            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_char, dim2d, varid), cpname, __LINE__ )
+            cattnam = 'long_name'
+            cattval = 'Observation type'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+
+            !  id
+            cvar = 'id'
+            WRITE(*,*) TRIM(cvar)
+            CALL chkerr(nf90_def_var(ncid, TRIM(cvar), nf90_char, dim2b, varid), cpname, __LINE__ )
+            cattnam = 'long_name'
+            cattval = 'Observation id'
+            CALL chkerr(nf90_put_att(ncid, varid, TRIM(cattnam), TRIM(cattval)), cpname, __LINE__ )
+
+            ! Close netcdf file
+            CALL chkerr(nf90_close(ncid), cpname, __LINE__ )
+         ELSE
+            WRITE(*, *) cdfilename, 'already exists.'
+            ! Close netcdf file
+            CALL chkerr(nf90_close(ncid), cpname, __LINE__ )
+         END IF
+      END SUBROUTINE ooo_wri_init
+
+
+END MODULE ooo_write
Index: branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/ooomod.F90
===================================================================
--- branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/ooomod.F90	(revision 4106)
+++ branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OOO_SRC/ooomod.F90	(revision 4106)
@@ -0,0 +1,98 @@
+MODULE ooomod
+   !!----------------------------------------------------------------------
+   !!                    ***  MODULE ooo ***
+   !! ** Purpose : Run NEMO observation operator in offline mode
+   !!----------------------------------------------------------------------
+   !! NEMO modules
+   USE diaobs
+   USE nemogcm
+   USE lib_mpp
+   USE in_out_manager
+   USE obs_fbm, ONLY: ln_cl4
+   !! Offline obs_oper modules
+   USE ooo_data
+   USE ooo_read
+   CONTAINS
+      SUBROUTINE ooo_init
+         !!----------------------------------------------------------------------
+         !!                    ***  SUBROUTINE ooo_init ***
+         !!
+         !! ** Purpose : To initialise the model as if it were running online.
+         !!
+         !!----------------------------------------------------------------------
+         IMPLICIT NONE
+         !! Initialise NEMO
+         CALL nemo_init
+         !! Initialise Offline obs_oper data
+         CALL ooo_data_init( ln_cl4 )
+      END SUBROUTINE ooo_init
+
+      SUBROUTINE offline_obs_oper
+         !!----------------------------------------------------------------------
+         !!                    ***  SUBROUTINE offline_obs_oper ***
+         !!
+         !! ** Purpose : To use NEMO components to interpolate model fields
+         !!              to observation space.
+         !!
+         !!----------------------------------------------------------------------
+         IMPLICIT NONE
+         !! Initialise offline obs_oper
+         CALL ooo_init
+         !! Loop over various model counterparts
+         DO jimatch = 1, cl4_match_len
+            IF (jimatch .GT. 1) THEN
+               !! Initialise obs_oper
+               CALL dia_obs_init
+            END IF
+            !! Interpolate to observation space
+            CALL ooo_interp
+            !! Pipe to output files
+            CALL dia_obs_wri
+            !! Reset the obs_oper
+            CALL dia_obs_dealloc
+         END DO
+         !! Safely stop MPI
+         CALL ooo_stop
+      END SUBROUTINE offline_obs_oper
+
+      SUBROUTINE ooo_interp
+         !!----------------------------------------------------------------------
+         !!                    ***  SUBROUTINE ooo_interp ***
+         !!
+         !! ** Purpose : To interpolate the model as if it were running online.
+         !!
+         !! ** Method : 1. Populate model counterparts
+         !!             2. Call dia_obs at appropriate time steps
+         !!----------------------------------------------------------------------
+         IMPLICIT NONE
+         INTEGER :: istp ! time step index
+         !! Loop over entire run
+         istp = nit000
+         nstop = 0
+         DO WHILE ( istp <= nitend .AND. nstop == 0 )
+            IF (jifile <= n_files + 1) THEN
+               IF ( MOD(istp, nn_off_freq) == nit000 ) THEN
+                  !! Read next model counterpart
+                  CALL ooo_rea_dri(jifile)
+                  jifile = jifile + 1
+               ENDIF
+               !! Interpolate single time step
+               CALL dia_obs(istp)
+            ENDIF
+            !! Increment model step
+            istp = istp + 1
+         END DO
+      END SUBROUTINE ooo_interp
+
+      SUBROUTINE ooo_stop
+         !!----------------------------------------------------------------------
+         !!                    ***  SUBROUTINE ooo_stop ***
+         !!
+         !! ** Purpose : To finalise the model as if it were running online.
+         !!
+         !!----------------------------------------------------------------------
+         IMPLICIT NONE
+         IF(lk_mpp) CALL mppstop  ! end mpp communications
+      END SUBROUTINE ooo_stop
+
+END MODULE ooomod
Index: branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OPA_SRC/OBS/obs_fbm.F90
===================================================================
--- branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OPA_SRC/OBS/obs_fbm.F90	(revision 4102)
+++ branches/2013/dev_r3987_UKMO4_OBS/NEMOGCM/NEMO/OPA_SRC/OBS/obs_fbm.F90	(revision 4106)
@@ -1570,6 +1570,6 @@
       !!----------------------------------------------------------------------
       USE dom_oce, ONLY: narea
-      USE off_write
-      USE off_data
+      USE ooo_write
+      USE ooo_data
       IMPLICIT NONE
       !! * Arguments
@@ -1734,10 +1734,10 @@
 
          ! Initialise class 4 file
-         CALL off_wri_init(cconf, csys, ckind, cversion, ccont, cinst, cdate, &
+         CALL ooo_wri_init(cconf, csys, ckind, cversion, ccont, cinst, cdate, &
                          & kproc, kobs, kvars, kdeps, kfcst, &
                          & clfilename)
 
          ! Write standard variables
-         CALL off_wri_default(clfilename, kobs, kvars, kfcst, kdeps, &
+         CALL ooo_wri_default(clfilename, kobs, kvars, kfcst, kdeps, &
                             & ctype, cwmo, cunit, cvarname, &
                             & plam, pphi, pdep, ptim, pob, plead, &
@@ -1748,9 +1748,9 @@
               (TRIM(cdtmp) == "persistence") ) THEN
             !! 4D variables
-            CALL off_wri_extra(clfilename, TRIM(cdtmp), kdeps, kfcst, &
+            CALL ooo_wri_extra(clfilename, TRIM(cdtmp), kdeps, kfcst, &
                             &  kvars, kobs, (/ 1,ij,1,1 /), (/ kdeps,1,kvars,kobs /), pmod)
          ELSE
             !! 3D variables
-            CALL off_wri_extra(clfilename, TRIM(cdtmp), kdeps, &
+            CALL ooo_wri_extra(clfilename, TRIM(cdtmp), kdeps, &
                             &  kvars, kobs, (/ 1,1,1 /), (/ kdeps,kvars,kobs /), pmod)
          ENDIF
