Index: branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/TOP_SRC/FABM/trcini_fabm.F90
===================================================================
--- branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/TOP_SRC/FABM/trcini_fabm.F90	(revision 10236)
+++ branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/TOP_SRC/FABM/trcini_fabm.F90	(revision 10390)
@@ -65,10 +65,49 @@
       jp_fabm_m1=jptra
       jptra = jptra + jp_fabm
-      jpdia2d = jpdia2d + size(model%horizontal_diagnostic_variables)
-      jpdia3d = jpdia3d + size(model%diagnostic_variables)
+      jp_fabm_2d = size(model%horizontal_diagnostic_variables)
+      jp_fabm_3d = size(model%diagnostic_variables)
+      jpdia2d = jpdia2d + jp_fabm_2d
+      jpdia3d = jpdia3d + jp_fabm_3d
       jpdiabio = jpdiabio + jp_fabm
 
       !Initialize input data structures.
       call initialize_inputs
+
+      ! Get indexes for select state variables
+      jp_fabm_chl1 = fabm_state_index( 'P1_Chl' )
+      jp_fabm_chl2 = fabm_state_index( 'P2_Chl' )
+      jp_fabm_chl3 = fabm_state_index( 'P3_Chl' )
+      jp_fabm_chl4 = fabm_state_index( 'P4_Chl' )
+      jp_fabm_p1c  = fabm_state_index( 'P1_c' )
+      jp_fabm_p1n  = fabm_state_index( 'P1_n' )
+      jp_fabm_p1p  = fabm_state_index( 'P1_p' )
+      jp_fabm_p1s  = fabm_state_index( 'P1_s' )
+      jp_fabm_p2c  = fabm_state_index( 'P2_c' )
+      jp_fabm_p2n  = fabm_state_index( 'P2_n' )
+      jp_fabm_p2p  = fabm_state_index( 'P2_p' )
+      jp_fabm_p3c  = fabm_state_index( 'P3_c' )
+      jp_fabm_p3n  = fabm_state_index( 'P3_n' )
+      jp_fabm_p3p  = fabm_state_index( 'P3_p' )
+      jp_fabm_p4c  = fabm_state_index( 'P4_c' )
+      jp_fabm_p4n  = fabm_state_index( 'P4_n' )
+      jp_fabm_p4p  = fabm_state_index( 'P4_p' )
+      jp_fabm_z4c  = fabm_state_index( 'Z4_c' )
+      jp_fabm_z5c  = fabm_state_index( 'Z5_c' )
+      jp_fabm_z5n  = fabm_state_index( 'Z5_n' )
+      jp_fabm_z5p  = fabm_state_index( 'Z5_p' )
+      jp_fabm_z6c  = fabm_state_index( 'Z6_c' )
+      jp_fabm_z6n  = fabm_state_index( 'Z6_n' )
+      jp_fabm_z6p  = fabm_state_index( 'Z6_p' )
+      jp_fabm_n1p  = fabm_state_index( 'N1_p' )
+      jp_fabm_n3n  = fabm_state_index( 'N3_n' )
+      jp_fabm_n4n  = fabm_state_index( 'N4_n' )
+      jp_fabm_n5s  = fabm_state_index( 'N5_s' )
+      jp_fabm_o2o  = fabm_state_index( 'O2_o' )
+      jp_fabm_o3c  = fabm_state_index( 'O3_c' )
+      jp_fabm_o3a  = fabm_state_index( 'O3_bioalk' )
+
+      ! Get indexes for select diagnostic variables
+      jp_fabm_o3ph = fabm_diag_index( 'O3_pH' )
+      jp_fabm_o3pc = fabm_diag_index( 'O3_pCO2' )
 
       IF (lwp) THEN
@@ -84,4 +123,6 @@
             CALL write_trends_xml(xml_unit,model%state_variables(jn))
 #endif
+            CALL write_25hourm_xml(xml_unit,model%state_variables(jn))
+            CALL write_tmb_xml(xml_unit,model%state_variables(jn))
          END DO
          WRITE (xml_unit,1000) ' </field_group>'
@@ -90,7 +131,9 @@
          DO jn=1,jp_fabm_surface
             CALL write_variable_xml(xml_unit,model%surface_state_variables(jn))
+            CALL write_25hourm_xml(xml_unit,model%surface_state_variables(jn))
          END DO
          DO jn=1,jp_fabm_bottom
             CALL write_variable_xml(xml_unit,model%bottom_state_variables(jn))
+            CALL write_25hourm_xml(xml_unit,model%bottom_state_variables(jn))
          END DO
          WRITE (xml_unit,1000) ' </field_group>'
@@ -99,7 +142,10 @@
          DO jn=1,size(model%diagnostic_variables)
             CALL write_variable_xml(xml_unit,model%diagnostic_variables(jn),3)
+            CALL write_25hourm_xml(xml_unit,model%diagnostic_variables(jn),3)
+            CALL write_tmb_xml(xml_unit,model%diagnostic_variables(jn))
          END DO
          DO jn=1,size(model%horizontal_diagnostic_variables)
             CALL write_variable_xml(xml_unit,model%horizontal_diagnostic_variables(jn))
+            CALL write_25hourm_xml(xml_unit,model%horizontal_diagnostic_variables(jn))
          END DO
          WRITE (xml_unit,1000) ' </field_group>'
@@ -168,4 +214,50 @@
 
    END SUBROUTINE write_variable_xml
+
+   SUBROUTINE write_25hourm_xml(xml_unit,variable,flag_grid_ref)
+      INTEGER,INTENT(IN) :: xml_unit
+      INTEGER,INTENT(IN),OPTIONAL :: flag_grid_ref
+      CLASS (type_external_variable),INTENT(IN) :: variable
+
+      CHARACTER(LEN=20) :: missing_value,string_dimensions
+      INTEGER :: number_dimensions
+
+      ! Check variable dimension for grid_ref specificaiton.
+      ! Default is to not specify the grid_ref in the field definition.
+      IF (present(flag_grid_ref)) THEN
+          number_dimensions=flag_grid_ref
+      ELSE
+          number_dimensions=-1 !default, don't specify grid_ref
+      ENDIF
+
+      WRITE (missing_value,'(E9.3)') 1.e+20
+      WRITE (string_dimensions,'(I1)') number_dimensions
+      SELECT CASE (number_dimensions)
+      CASE (3)
+         WRITE (xml_unit,'(A)') '  <field id="'//TRIM(variable%name)//'25h'//'" long_name="'//TRIM(variable%long_name)//' 25-hour mean'//'" unit="'//TRIM(variable%units)//'" default_value="'//TRIM(ADJUSTL(missing_value))//'" grid_ref="grid_T_3D" />'
+      CASE (2)
+         WRITE (xml_unit,'(A)') '  <field id="'//TRIM(variable%name)//'25h'//'" long_name="'//TRIM(variable%long_name)//' 25-hour mean'//'" unit="'//TRIM(variable%units)//'" default_value="'//TRIM(ADJUSTL(missing_value))//'" grid_ref="grid_T_2D"/>'
+      CASE (0)
+         WRITE (xml_unit,'(A)') '  <field id="'//TRIM(variable%name)//'25h'//'" long_name="'//TRIM(variable%long_name)//' 25-hour mean'//'" unit="'//TRIM(variable%units)//'" default_value="'//TRIM(ADJUSTL(missing_value))//'" grid_ref="1point"/>'
+      CASE (-1)
+         WRITE (xml_unit,'(A)') '  <field id="'//TRIM(variable%name)//'25h'//'" long_name="'//TRIM(variable%long_name)//' 25-hour mean'//'" unit="'//TRIM(variable%units)//'" default_value="'//TRIM(ADJUSTL(missing_value))//'" />'
+      CASE default
+         IF(lwp) WRITE(numout,*) ' trc_ini_fabm: Failing to initialise output of variable '//TRIM(variable%name)//'25h'//': Output of '//TRIM(ADJUSTL(string_dimensions))//'-dimensional variables not supported!!!'
+      END SELECT
+
+   END SUBROUTINE write_25hourm_xml
+
+   SUBROUTINE write_tmb_xml(xml_unit,variable)
+      INTEGER,INTENT(IN) :: xml_unit
+      CLASS (type_external_variable),INTENT(IN) :: variable
+
+      CHARACTER(LEN=20) :: missing_value
+
+      WRITE (missing_value,'(E9.3)') 1.e+20
+      WRITE (xml_unit,'(A)') '  <field id="top_'//TRIM(variable%name)//'" long_name="Top-level '//TRIM(variable%long_name)//'" unit="'//TRIM(variable%units)//'" default_value="'//TRIM(ADJUSTL(missing_value))//'" grid_ref="grid_T_2D"/>'
+      WRITE (xml_unit,'(A)') '  <field id="mid_'//TRIM(variable%name)//'" long_name="Middle-level '//TRIM(variable%long_name)//'" unit="'//TRIM(variable%units)//'" default_value="'//TRIM(ADJUSTL(missing_value))//'" grid_ref="grid_T_2D"/>'
+      WRITE (xml_unit,'(A)') '  <field id="bot_'//TRIM(variable%name)//'" long_name="Bottom-level '//TRIM(variable%long_name)//'" unit="'//TRIM(variable%units)//'" default_value="'//TRIM(ADJUSTL(missing_value))//'" grid_ref="grid_T_2D"/>'
+
+   END SUBROUTINE write_tmb_xml
 
    SUBROUTINE write_trends_xml(xml_unit,variable,flag_grid_ref)
@@ -328,4 +420,68 @@
    END SUBROUTINE trc_ini_fabm
 
+   INTEGER FUNCTION fabm_state_index( state_name )
+      !!----------------------------------------------------------------------
+      !!                     ***  fabm_state_index  ***  
+      !!
+      !! ** Purpose :   return index of a given FABM state variable
+      !!
+      !! ** Method  : - loop through state variables until found
+      !!----------------------------------------------------------------------
+      
+      IMPLICIT NONE
+      
+      CHARACTER(LEN=256), INTENT(IN) :: state_name
+      
+      INTEGER                        :: jn
+
+      !!----------------------------------------------------------------------
+      
+      fabm_state_index = -1
+      DO jn=1,jp_fabm
+         IF (TRIM(model%state_variables(jn)%name) == TRIM(state_name)) THEN
+            fabm_state_index = jn
+            EXIT
+         ENDIF
+      END DO
+      IF (fabm_state_index == -1) THEN
+         CALL ctl_stop( 'Could not find '//TRIM(state_name)//' state variable' )
+      ELSE
+         IF (lwp) WRITE(numout,*) 'Index for '//TRIM(state_name)//' is: ', fabm_state_index
+      ENDIF
+   
+   END FUNCTION fabm_state_index
+
+   INTEGER FUNCTION fabm_diag_index( diag_name )
+      !!----------------------------------------------------------------------
+      !!                     ***  fabm_state_index  ***  
+      !!
+      !! ** Purpose :   return index of a given FABM diagnostic variable
+      !!
+      !! ** Method  : - loop through diagnostic variables until found
+      !!----------------------------------------------------------------------
+      
+      IMPLICIT NONE
+      
+      CHARACTER(LEN=256), INTENT(IN) :: diag_name
+      
+      INTEGER                        :: jn
+
+      !!----------------------------------------------------------------------
+      
+      fabm_diag_index = -1
+      DO jn = 1, SIZE(model%diagnostic_variables)
+         IF (TRIM(model%diagnostic_variables(jn)%name) == TRIM(diag_name)) THEN
+            fabm_diag_index = jn
+            EXIT
+         ENDIF
+      END DO
+      IF (fabm_diag_index == -1) THEN
+         CALL ctl_stop( 'Could not find '//TRIM(diag_name)//' diagnostic' )
+      ELSE
+         IF (lwp) WRITE(numout,*) 'Index for '//TRIM(diag_name)//' is: ', fabm_diag_index
+      ENDIF
+   
+   END FUNCTION fabm_diag_index
+
 #else
    !!----------------------------------------------------------------------
