Index: branches/UKMO/dev_r4650_general_vert_coord_obsoper_surf_bgc/NEMOGCM/NEMO/OPA_SRC/OBS/diaobs.F90
===================================================================
--- branches/UKMO/dev_r4650_general_vert_coord_obsoper_surf_bgc/NEMOGCM/NEMO/OPA_SRC/OBS/diaobs.F90	(revision 6855)
+++ branches/UKMO/dev_r4650_general_vert_coord_obsoper_surf_bgc/NEMOGCM/NEMO/OPA_SRC/OBS/diaobs.F90	(revision 6856)
@@ -32,4 +32,5 @@
    USE obs_read_logchl          ! Reading and allocation of logchl observations
    USE obs_read_spm             ! Reading and allocation of spm observations
+   USE obs_read_fco2            ! Reading and allocation of fco2 observations
    USE obs_prep                 ! Preparation of obs. (grid search etc).
    USE obs_oper                 ! Observation operators
@@ -45,4 +46,5 @@
    USE obs_logchl               ! logchl data storage
    USE obs_spm                  ! spm data storage
+   USE obs_fco2                 ! fco2 data storage
    USE obs_types                ! Definitions for observation types
    USE mpp_map                  ! MPP mapping
@@ -90,4 +92,6 @@
    LOGICAL, PUBLIC :: ln_spm         !: Logical switch for spm
    LOGICAL, PUBLIC :: ln_spmfb       !: Logical switch for spm from feedback files
+   LOGICAL, PUBLIC :: ln_fco2        !: Logical switch for fco2
+   LOGICAL, PUBLIC :: ln_fco2fb      !: Logical switch for fco2 from feedback files
    LOGICAL, PUBLIC :: ln_ssh         !: Logical switch for sea surface height
    LOGICAL, PUBLIC :: ln_sss         !: Logical switch for sea surface salinity
@@ -177,4 +181,6 @@
       CHARACTER(len=128) :: spmfiles(MaxNumFiles)
       CHARACTER(len=128) :: spmfbfiles(MaxNumFiles)
+      CHARACTER(len=128) :: fco2files(MaxNumFiles)
+      CHARACTER(len=128) :: fco2fbfiles(MaxNumFiles)
       CHARACTER(LEN=128) :: reysstname
       CHARACTER(LEN=12)  :: reysstfmt
@@ -205,4 +211,6 @@
          &            ln_spm, ln_spmfb,                               &
          &            spmfiles, spmfbfiles,                           &
+         &            ln_fco2, ln_fco2fb,                             &
+         &            fco2files, fco2fbfiles,                         &
          &            ln_profb_enatim, ln_ignmis, ln_cl4,             &
          &            ln_sstbias, sstbias_files
@@ -230,4 +238,6 @@
       INTEGER :: jnumspm
       INTEGER :: jnumspmfb
+      INTEGER :: jnumfco2
+      INTEGER :: jnumfco2fb
       INTEGER :: ji
       INTEGER :: jset
@@ -243,4 +253,6 @@
       ln_spm      = .FALSE.
       ln_spmfb    = .FALSE.
+      ln_fco2     = .FALSE.
+      ln_fco2fb   = .FALSE.
       
       !Initalise all values in namelist arrays
@@ -267,4 +279,6 @@
       spmfiles(:) = ''
       spmfbfiles(:) = ''
+      fco2files(:) = ''
+      fco2fbfiles(:) = ''
       sstbias_files(:) = ''
       endailyavtypes(:) = -1
@@ -385,4 +399,14 @@
          WHERE (spmfbfiles(:) /= '') lmask(:) = .TRUE.
          jnumspmfb = COUNT(lmask)
+      ENDIF
+      IF (ln_fco2) THEN
+         lmask(:) = .FALSE.
+         WHERE (fco2files(:) /= '') lmask(:) = .TRUE.
+         jnumfco2 = COUNT(lmask)
+      ENDIF
+      IF (ln_fco2fb) THEN
+         lmask(:) = .FALSE.
+         WHERE (fco2fbfiles(:) /= '') lmask(:) = .TRUE.
+         jnumfco2fb = COUNT(lmask)
       ENDIF
       
@@ -420,4 +444,6 @@
          WRITE(numout,*) '             Logical switch for spm observations                ln_spm = ', ln_spm
          WRITE(numout,*) '             Logical switch for feedback spm data             ln_spmfb = ', ln_spmfb
+         WRITE(numout,*) '             Logical switch for fco2 observations              ln_fco2 = ', ln_fco2
+         WRITE(numout,*) '             Logical switch for feedback fco2 data           ln_fco2fb = ', ln_fco2fb
          WRITE(numout,*) '             Global distribtion of observations         ln_grid_global = ',ln_grid_global
          WRITE(numout,*) &
@@ -540,4 +566,16 @@
                WRITE(numout,'(1X,2A)') '             Feedback spm input observation file name   spmfbfiles = ', &
                   TRIM(spmfbfiles(ji))
+            END DO
+         ENDIF
+         IF (ln_fco2) THEN
+            DO ji = 1, jnumfco2
+               WRITE(numout,'(1X,2A)') '             fco2 input observation file name  fco2files = ', &
+                  TRIM(fco2files(ji))
+            END DO
+         ENDIF
+         IF (ln_fco2fb) THEN
+            DO ji = 1, jnumfco2fb
+               WRITE(numout,'(1X,2A)') '             Feedback fco2 input observation file name  fco2fbfiles = ', &
+                  TRIM(fco2fbfiles(ji))
             END DO
          ENDIF
@@ -577,9 +615,9 @@
          & ( .NOT. ln_ssh ).AND.( .NOT. ln_sst ).AND.( .NOT. ln_sss ).AND. &
          & ( .NOT. ln_seaice ).AND.( .NOT. ln_vel3d ).AND.( .NOT. ln_logchl ).AND. &
-         & ( .NOT. ln_spm ) ) THEN
+         & ( .NOT. ln_spm ).AND.( .NOT. ln_fco2 ) ) THEN
          IF(lwp) WRITE(numout,cform_war)
          IF(lwp) WRITE(numout,*) ' key_diaobs is activated but logical flags', &
             &                    ' ln_t3d, ln_s3d, ln_sla, ln_ssh, ln_sst, ln_sss, ln_seaice, ln_vel3d,', &
-            &                    ' ln_logchl, ln_spm are all set to .FALSE.'
+            &                    ' ln_logchl, ln_spm, ln_fco2 are all set to .FALSE.'
          nwarn = nwarn + 1
       ENDIF
@@ -1189,4 +1227,59 @@
  
       ENDIF
+
+      !  - fco2
+      
+      IF ( ln_fco2 ) THEN
+
+         ! Set the number of variables for fco2 to 1
+         nfco2vars = 1
+
+         ! Set the number of extra variables for fco2 to 0
+         nfco2extr = 0
+         
+         IF ( ln_fco2fb ) THEN
+            nfco2sets = jnumfco2fb
+         ELSE
+            nfco2sets = 1
+         ENDIF
+
+         ALLOCATE(fco2data(nfco2sets))
+         ALLOCATE(fco2datqc(nfco2sets))
+         fco2data(:)%nsurf=0
+         fco2datqc(:)%nsurf=0
+
+         nfco2sets = 0
+
+         IF ( ln_fco2fb ) THEN             ! Feedback file format
+
+            DO jset = 1, jnumfco2fb
+            
+               nfco2sets = nfco2sets + 1
+
+               CALL obs_rea_fco2( 0, fco2data(nfco2sets), 1, &
+                  &                 fco2fbfiles(jset:jset), &
+                  &                 nfco2vars, nfco2extr, nitend-nit000+2, &
+                  &                 dobsini, dobsend, ln_ignmis, .FALSE. )
+
+               CALL obs_pre_fco2( fco2data(nfco2sets), fco2datqc(nfco2sets), &
+                  &                 ln_fco2, ln_nea )
+            
+            ENDDO
+
+         ELSE                              ! Original file format
+
+            nfco2sets = nfco2sets + 1
+
+            CALL obs_rea_fco2( 1, fco2data(nfco2sets), jnumfco2, &
+               &                 fco2files(1:jnumfco2), &
+               &                 nfco2vars, nfco2extr, nitend-nit000+2, &
+               &                 dobsini, dobsend, ln_ignmis, .FALSE. )
+
+            CALL obs_pre_fco2( fco2data(nfco2sets), fco2datqc(nfco2sets), &
+               &                 ln_fco2, ln_nea )
+
+         ENDIF
+ 
+      ENDIF
      
    END SUBROUTINE dia_obs_init
@@ -1208,4 +1301,5 @@
       !!               - Sea surface log10(chlorophyll)
       !!               - Sea surface spm
+      !!               - Sea surface fco2
       !!
       !! ** Action  : 
@@ -1246,13 +1340,15 @@
 #endif
 #if defined key_hadocc
-      USE trc, ONLY :  &                ! HadOCC chlorophyll
+      USE trc, ONLY :  &                ! HadOCC chlorophyll and fCO2
          & HADOCC_CHL, &
+         & HADOCC_FCO2, &
          & HADOCC_FILL_FLT
 #elif defined key_medusa && defined key_foam_medusa
-      USE trc, ONLY :  &                ! MEDUSA chlorophyll
+      USE trc, ONLY :  &                ! MEDUSA chlorophyll and fCO2
          & MEDUSA_CHL, &
+         & MEDUSA_FCO2, &
          & MEDUSA_FILL_FLT
 #elif defined key_fabm
-      !USE ???                           ! ERSEM chlorophyll
+      !USE ???                           ! ERSEM chlorophyll and fCO2
 #endif
 #if defined key_spm
@@ -1273,4 +1369,5 @@
       INTEGER :: jlogchlset             ! logchl data set loop variable
       INTEGER :: jspmset                ! spm data set loop variable
+      INTEGER :: jfco2set               ! fco2 data set loop variable
       INTEGER :: jvar                   ! Variable number    
 #if ! defined key_lim2 && ! defined key_lim3
@@ -1284,4 +1381,8 @@
       REAL(wp), DIMENSION(jpi,jpj) :: &
          spm                            ! array for spm
+      REAL(wp), DIMENSION(jpi,jpj) :: &
+         fco2                           ! array for fco2
+      REAL(wp), DIMENSION(jpi,jpj) :: &
+         maskfco2                       ! array for special fco2 mask
       INTEGER :: jn                     ! loop index
       CHARACTER(LEN=20) :: datestr=" ",timestr=" "
@@ -1445,4 +1546,40 @@
       ENDIF
 
+      IF ( ln_fco2 ) THEN
+         maskfco2(:,:) = tmask(:,:,1)         ! create a special mask to exclude certain things
+#if defined key_hadocc
+         fco2(:,:) = HADOCC_FCO2(:,:)    ! fCO2 from HadOCC
+         IF ( ( MINVAL( HADOCC_FCO2 ) == HADOCC_FILL_FLT ).AND.( MAXVAL( HADOCC_FCO2 ) == HADOCC_FILL_FLT ) ) THEN
+            fco2(:,:) = obfillflt
+            maskfco2(:,:) = 0
+            CALL ctl_warn( ' HadOCC fCO2 values masked out for observation operator', &
+               &           ' on timestep ' // TRIM(STR(kstp)),                              &
+               &           ' as HADOCC_FCO2(:,:) == HADOCC_FILL_FLT' )
+         ENDIF
+#elif defined key_medusa && defined key_foam_medusa
+         fco2(:,:) = MEDUSA_FCO2(:,:)    ! fCO2 from MEDUSA
+         IF ( ( MINVAL( MEDUSA_FCO2 ) == MEDUSA_FILL_FLT ).AND.( MAXVAL( MEDUSA_FCO2 ) == MEDUSA_FILL_FLT ) ) THEN
+            fco2(:,:) = obfillflt
+            maskfco2(:,:) = 0
+            CALL ctl_warn( ' MEDUSA fCO2 values masked out for observation operator', &
+               &           ' on timestep ' // TRIM(STR(kstp)),                              &
+               &           ' as MEDUSA_FCO2(:,:) == MEDUSA_FILL_FLT' )
+         ENDIF
+#elif defined key_fabm
+         !fco2(:,:)  =  ???                 ! fCO2 from ERSEM
+         CALL ctl_stop( ' Trying to run fco2 observation operator', &
+            &           ' but not properly implemented for FABM-ERSEM yet' )
+#else
+         CALL ctl_stop( ' Trying to run fco2 observation operator', &
+            &           ' but no biogeochemical model appears to have been defined' )
+#endif
+
+         DO jfco2set = 1, nfco2sets
+             CALL obs_fco2_opt( fco2datqc(jfco2set),                      &
+               &                kstp, jpi, jpj, nit000, fco2(:,:), &
+               &                maskfco2(:,:), n2dint )
+         END DO
+      ENDIF
+
 #if ! defined key_lim2 && ! defined key_lim3
       CALL wrk_dealloc(jpi,jpj,frld) 
@@ -1479,4 +1616,5 @@
       INTEGER :: jlogchlset               ! logchl data set loop variable
       INTEGER :: jspmset                  ! spm data set loop variable
+      INTEGER :: jfco2set                 ! fco2 data set loop variable
       INTEGER :: jset
       INTEGER :: jfbini
@@ -1771,4 +1909,34 @@
             WRITE(cdtmp,'(A,I2.2)')'spmfb_',jspmset
             CALL obs_wri_spm( cdtmp, spmdata(jspmset) )
+
+         END DO
+
+      ENDIF
+
+      !  - fco2
+      IF ( ln_fco2 ) THEN
+
+         ! Copy data from fco2datqc to fco2data structures
+         DO jfco2set = 1, nfco2sets
+
+            CALL obs_surf_decompress( fco2datqc(jfco2set), &
+                 &                    fco2data(jfco2set), .TRUE., numout )
+
+         END DO
+         
+         ! Mark as bad observations with no valid model counterpart due to fCO2 not being in the restart
+         ! Seem to need to set to fill value rather than marking as bad to be effective, so do both
+         DO jfco2set = 1, nfco2sets
+            WHERE ( fco2data(jfco2set)%rmod(:,1) == obfillflt )
+               fco2data(jfco2set)%nqc(:)    = 1
+               fco2data(jfco2set)%robs(:,1) = obfillflt
+            END WHERE
+         END DO
+
+         ! Write the fco2 data
+         DO jfco2set = 1, nfco2sets
+      
+            WRITE(cdtmp,'(A,I2.2)')'fco2fb_',jfco2set
+            CALL obs_wri_fco2( cdtmp, fco2data(jfco2set) )
 
          END DO
