Scripts working freeboard and thickness

EnKF-MPI-TOPAZ/m_obs.F90

@@ -81,6 +81,10 @@ module m_obs
   real, parameter, private :: SSH_MAX = 3.0d0
   real, parameter, private :: ICEC_MIN = 0.0d0
   real, parameter, private :: ICEC_MAX = 0.999d0        ! [FM] Changed from 0.996 to 0.999
+  real, parameter, private :: RFB_MIN  = 0.0d0          ! FM 2020
+  real, parameter, private :: RFB_MAX  = 10.0d0
+  real, parameter, private :: VT_I_MIN  = 0.0d0          ! FM 2020
+  real, parameter, private :: VT_I_MAX  = 10.0d0
   real, parameter, private :: UVICE_MIN = -100.0
   real, parameter, private :: UVICE_MAX = 100.0
 
@@ -100,28 +104,56 @@ contains
     integer :: rsize
     integer :: ios
     integer :: o
+    CHARACTER(LEN=*), PARAMETER  :: &
+ FMT2 = "(f8.4,X,f8.4,X,a8,X,2(f10.5,X),f4.2,X,2(I3,X),I1,X,4(f5.2,X),L,X,2(I3,X),f5.2,X,I8,X,I1)"
+    real :: myX
+    real :: myY
+
+
+
+!==========  TEST
+!    inquire(iolength = rsize) record
+!    !open(10, file = 'test.txt', form = 'unformatted',&
+!    !     access = 'direct', recl = rsize, status = 'old')
+!    allocate(obs(2))
+!    open(10, file = 'observations.txt')!, form = 'unformatted',&
+!          !access = 'direct', recl = rsize, status = 'old')
+!    !read(10, *) obs(1)
+!
+!    do o = 1, 2
+!      read(10, *) obs(o)
+!      PRINT *, obs(o)
+!    end do
+!    close(10)
+!    stop
+!==========
 
     if (nobs >= 0) then
        return
     end if
 
-    inquire(file = 'observations.uf', exist = exists)
+    ! Testing existence of file
+    inquire(file = 'observations.txt', exist = exists)
+    !inquire(file = 'observations.uf', exist = exists)
     if (.not. exists) then
        if (master) then
-          print *, 'ERROR: obs_getnobs(): file "observations.uf" does not exist'
+          print *, 'ERROR: obs_getnobs(): file "observations.txt" does not exist'
        end if
        stop
     end if
     inquire(iolength = rsize) record
-    open(10, file = 'observations.uf', form = 'unformatted',&
-         access = 'direct', recl = rsize, status = 'old')
-
+    open(10, file = 'observations.txt')!, form = 'unformatted',&
+    ! EXPERIMENTAL
+    !open(10, file = 'observations.uf', form = 'unformatted',&
+    !         access = 'direct', recl = rsize, status = 'old')!, form = 'unformatted',&
+         !access = 'direct', recl = rsize, status = 'old')
+    ! END EXPERIMENTAL
     ! I guess there is no other way to work out the length other than read the
     ! file in fortran - PS
     !
     o = 1
     do while (.true.)
-       read(10, rec = o, iostat = ios) record
+       read(10, *, iostat = ios) record
        if (ios /= 0) then
           nobs = o - 1
           exit
@@ -135,11 +167,17 @@ contains
     ! "Cannot REWIND a file opened for DIRECT access". Therefore reopen.
     !
     close(10)
-    open(10, file = 'observations.uf', form = 'unformatted',&
-         access = 'direct', recl = rsize, status = 'old')
+    open(10, file = 'observations.txt')!, form = 'unformatted',&
+
+    ! BEGIN EXPERIMENTAL
+    !open(10, file = 'observations.uf', form = 'unformatted',&
+    !         access = 'direct', recl = rsize, status = 'old')
+    ! -- END EXPERIMENTAL
     do o = 1, nobs
-       read(10, rec = o) obs(o)
+       read(10, *) obs(o)
+       
        call ucase(obs(o) % id)
+       !PRINT *, obs(o)
     enddo
     close(10)
 
@@ -152,6 +190,8 @@ contains
    call  uobs_get(obs % id, nobs, master)
    
    call obs_testrange
+
+
   end subroutine obs_readobs
 
 
@@ -182,6 +222,12 @@ contains
        elseif (trim(unique_obs(uo)) == 'AT_I') then     ! [FM] Added as we assimilate total ice conc. (opposed to indiv. category
           dmin = ICEC_MIN
           dmax = ICEC_MAX
+       elseif (trim(unique_obs(uo)) == 'RFB') then      ! FM added 2020
+          dmin = RFB_MIN
+          dmax = RFB_MAX
+       elseif (trim(unique_obs(uo)) == 'VT_I') then      ! FM added 2021
+          dmin = VT_I_MIN
+          dmax = VT_I_MAX 
        elseif (trim(unique_obs(uo)) == 'V_ICE'&
             .or. trim(unique_obs(uo)) == 'U_ICE') then
           dmin = UVICE_MIN

EnKF-MPI-TOPAZ/m_prep_4_EnKF.F90

@@ -72,6 +72,8 @@ contains
     integer :: i, j, m, iens
     real*4, dimension(nx,ny) :: fldr4
     real :: readfld(nx, ny), ai1(nx,ny), ai2(nx,ny), ai3(nx,ny), ai4(nx,ny), ai5(nx,ny), uice(nx,ny), vice(nx,ny)
+    real :: vi1(nx,ny), vi2(nx,ny), vi3(nx,ny), vi4(nx,ny), vi5(nx,ny)
+    real :: vs1(nx,ny), vs2(nx,ny), vs3(nx,ny), vs4(nx,ny), vs5(nx,ny)
 
     ! hard-coded for now
     integer, parameter :: drnx = 152, drny = 132
@@ -168,12 +170,119 @@ contains
                 end if
                 stop
              end if
-             readfld=ai1+ai2+ai3+ai4+ai5
+             ! Multipply by 100 to match obs conventions
+             readfld=(ai1+ai2+ai3+ai4+ai5) * 100.0
 
              call Generate_element_Si(S(:, iens), unique_obs(iuobs),&
                   readfld, depths, nx, ny, nz, 0)
           end do
 
+       ! freeboard
+       elseif(trim(unique_obs(iuobs)) == 'VT_I') then
+           do iens = 1, nrens
+             write(cmem, '(i3.3)') iens
+             tlevel = 1
+             call io_mod_fld(ai1,iens,enslist, &
+                 'a_i_htc1', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+             call io_mod_fld(ai2,iens,enslist, &
+                 'a_i_htc2', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+             call io_mod_fld(ai3,iens,enslist, &
+                 'a_i_htc3', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+             call io_mod_fld(ai4,iens,enslist, &
+                 'a_i_htc4', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+             call io_mod_fld(ai5,iens,enslist, &
+                 'a_i_htc5', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+
+             call io_mod_fld(vi1,iens,enslist, &
+                 'v_i_htc1', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+             call io_mod_fld(vi2,iens,enslist, &
+                 'v_i_htc2', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+             call io_mod_fld(vi3,iens,enslist, &
+                 'v_i_htc3', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+             call io_mod_fld(vi4,iens,enslist, &
+                 'v_i_htc4', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+             call io_mod_fld(vi5,iens,enslist, &
+                 'v_i_htc5', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+
+
+             if (tlevel == -1) then
+                 if (master) then
+                   print *, 'ERROR: io_mod_fld_nex(): failed for "SIFB"'
+                 end if
+                 stop
+             end if
+
+
+
+             readfld=(vi1+vi2+vi3+vi4+vi5) 
+             call Generate_element_Si(S(:, iens), unique_obs(iuobs),&
+                  readfld, depths, nx, ny, nz, 0)
+           end do
+
+
+       ! freeboard
+       elseif(trim(unique_obs(iuobs)) == 'RFB') then
+           do iens = 1, nrens
+             write(cmem, '(i3.3)') iens
+             tlevel = 1
+             call io_mod_fld(ai1,iens,enslist, &
+                 'a_i_htc1', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+             call io_mod_fld(ai2,iens,enslist, &
+                 'a_i_htc2', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+             call io_mod_fld(ai3,iens,enslist, &
+                 'a_i_htc3', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+             call io_mod_fld(ai4,iens,enslist, &
+                 'a_i_htc4', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+             call io_mod_fld(ai5,iens,enslist, &
+                 'a_i_htc5', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+
+             call io_mod_fld(vi1,iens,enslist, &
+                 'v_i_htc1', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+             call io_mod_fld(vi2,iens,enslist, &
+                 'v_i_htc2', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+             call io_mod_fld(vi3,iens,enslist, &
+                 'v_i_htc3', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+             call io_mod_fld(vi4,iens,enslist, &
+                 'v_i_htc4', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+             call io_mod_fld(vi5,iens,enslist, &
+                 'v_i_htc5', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+
+             call io_mod_fld(vs1,iens,enslist, &
+                 'v_s_htc1', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+             call io_mod_fld(vs2,iens,enslist, &
+                 'v_s_htc2', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+             call io_mod_fld(vs3,iens,enslist, &
+                 'v_s_htc3', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+             call io_mod_fld(vs4,iens,enslist, &
+                 'v_s_htc4', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+             call io_mod_fld(vs5,iens,enslist, &
+                 'v_s_htc5', 1, 0, 1, nx,ny, 'get',FLOAT(obs(1)%date))
+
+             if (tlevel == -1) then
+                 if (master) then
+                   print *, 'ERROR: io_mod_fld_nex(): failed for "SIFB"'
+                 end if
+                 stop
+             end if
+
+             readfld=(((vi1+vi2+vi3+vi4+vi5) * (1024.0 - 899.5) - 330 * (vs1+vs2+vs3+vs4+vs5)) / &
+                    1024.0-0.25*(vs1 +vs2+vs3+vs4+vs5)) 
+             !readfld=(((vi1+vi2+vi3+vi4+vi5) * (1024.0 - 899.5) - 330 * (vs1+vs2+vs3+vs4+vs5)) / 1024.0 - 0.25 * (vs1+vs2+vs3+vs4+vs5)) / (ai1+ai2+ai3+ai4+ai5)
+             
+             ! Conversion of models' sea ice thickness and snow thickness to
+             ! model's freeboard using fixed densities for snow (330 kg/m3), ice
+             ! (899.5 kg/m3 = average of MYI and FYI from Guerreiro et al. 2017
+             ! and seawater (1024 kg/m3). The model freeboard is then lowered by
+             ! 25% of the snow depth to account for the fact that the radar
+             ! measurement underestimates the actual freeboard due to the lower
+             ! propagation speed of the wave into the snow than in the air.
+             ! Everything is converted from grid cell mean to in situ by
+             ! dividing by concentration (if it is not zero). See exchanges
+             ! e-mail with Sara Fleury 7 December 2020.
+             call Generate_element_Si(S(:, iens), unique_obs(iuobs),&
+                  readfld, depths, nx, ny, nz, 0)
+           end do
+
        elseif (trim(unique_obs(iuobs)) == 'SST') then
           do iens = 1, nrens
              write(cmem,'(i3.3)') iens

EnKF-MPI-TOPAZ/m_uobs.F90

@@ -25,7 +25,7 @@ module m_uobs
 
   public uobs_get
   
-  integer, parameter, private :: MAXNUOBS = 19
+  integer, parameter, private :: MAXNUOBS = 1900
 
   integer, public :: nuobs
   character(OBSTYPESTRLEN), public :: unique_obs(MAXNUOBS)
@@ -52,8 +52,11 @@ contains
     nuobs = 0
     unique_obs = ''
     do o = 1, nrobs
+       !PRINT *, o
+       !PRINT *, '...', nuobs
        obsmatch = .false.
        do uo = 1, nuobs
+          !PRINT *, '-->', uo
           if (trim(tags(o)) == trim(unique_obs(uo))) then
              obsmatch = .true.
              nobseach(uo) = nobseach(uo) + 1
@@ -90,6 +93,8 @@ contains
        do uo = 1, nuobs
           do o = uobs_begin(uo), uobs_end(uo)
              if (trim(tags(o)) /= trim(unique_obs(uo))) then
+                print *, trim(tags(o))
+                print *, trim(unique_obs(uo))
                 print *, 'ERROR: uobs_get(): uinique observations not ',&
                      'continuous in observation array'
                 stop

conversion_uf/make.inc

@@ -1 +1 @@
-../EnKF-MPI-TOPAZ/Config/make.mn4
+../EnKF-MPI-TOPAZ/Config/make.zenobe

conversion_uf/p_prep_obs_ORCA1.F90

@@ -117,9 +117,10 @@ program prep_obs
   read(dummy,'(A)') varname
   call getarg(3, dummy)
   read(dummy,*) rdate
-  if (      trim(varname).ne.'sic' ) then
-    write(*,*) ' (prep_obs) Currently only the total sea ice concentration'
-    write(*,*) '            sic    can be processed                      '
+  if (      trim(varname) .ne.'rfb' .and. (trim(varname) .ne. 'vt_i') ) then
+    write(*,*) ' (prep_obs) Currently only the radar freeboard'
+    write(*,*) '            rfb    can be processed                      '
+    PRINT *, varname
     stop " (prep_obs) Aborted"
   endif
 
@@ -139,11 +140,11 @@ program prep_obs
 
   ! Find VarID of varname
 
-  cvarerror='error_'//varname
+  cvarerror=TRIM(varname)//'_sd'
 
-   ierr = nf90_inq_varid(ncid, cvarerror, varID) ; if (ierr.ne.nf90_noerr) call handle_err(ierr, "inquiring varID")
-   ! get values
-   ierr = nf90_get_var(ncid, varID, errorfld) ; if (ierr.ne.nf90_noerr) call handle_err(ierr, "getting variable")
+  ierr = nf90_inq_varid(ncid, cvarerror, varID) ; if (ierr.ne.nf90_noerr) call handle_err(ierr, "inquiring varID")
+  ! get values
+  ierr = nf90_get_var(ncid, varID, errorfld) ; if (ierr.ne.nf90_noerr) call handle_err(ierr, "getting variable")
 
   ! Close file
   ierr = nf90_close(ncid) ; if (ierr.ne.nf90_noerr) call handle_err(ierr, "closing")
@@ -158,19 +159,22 @@ program prep_obs
     do j=1,size(mask,2)
 
       ! Pick out ocean points where data is existent
-      if ( errorfld(i,j) .GE. 1.0 ) then
-
+      if ( (errorfld(i,j)) ** 2 .GT. 0.0001 .and. (errorfld(i,j)) ** 2 .LT. 1.0 &
+           .and. obsfld(i,j) .LT. 10.0 .and. obsfld(i,j) .GT. 0.0001 ) then
+        
         ! Limit 'obs' to >40°N or <45°S
-        if ( (lats(i,j).ge.40) .or. (lats(i,j).le.-45) ) then
+        !if ( (lats(i,j).ge.40) .or. (lats(i,j).le.-45) ) then
+         if ( (lats(i,j) .ge. 40.0) ) then !.and. (lats(i,j) .le. 80.05) .and. (lons(i,j) .ge. 100.0)) then
+        !if ( (lats(i,j) .ge. 63.0)  .and. (lats(i,j) .le. 65.0))   then
           obscnt = obscnt + 1
-          obs(obscnt)%d    = obsfld(i,j)/100.0
+          obs(obscnt)%d    = obsfld(i,j)
           obs(obscnt)%lon  = lons(i,j)
           obs(obscnt)%lat  = lats(i,j)
           obs(obscnt)%ipiv = i ! the i-point of the grid of the model
           obs(obscnt)%jpiv = j ! the j-point of the grid of the model
 
           ! Put other data into obs type array
-          obs(obscnt)%var    = (errorfld(i,j)/100.0)**2 ! The variance
+          obs(obscnt)%var    = (errorfld(i,j))**2 ! The variance
           obs(obscnt)%id     = TRIM(varname)
           obs(obscnt)%depths = 0
           obs(obscnt)%ns     = 0
@@ -190,6 +194,7 @@ program prep_obs
     end do ! j
   end do ! i
 
+  !obscnt = 3
   !Write data:
   inquire(iolength=i)obs(1)
   open (11, file='observations.uf', status='replace',form='unformatted', access='direct', recl=i)
@@ -198,7 +203,6 @@ program prep_obs
   enddo
   close(11)
 
-PRINT *, rdate
 
 !  !XXX: 
 !  ! Write data to textfile, for visual checking
@@ -206,6 +210,7 @@ PRINT *, rdate
   do j=1,obscnt
      write(12,FMT=53) obs(j)
 53   FORMAT(f8.4,X,f8.4,X,a8,X,2(f10.5,X),f4.2,X,2(I3,X),I1,X,4(f5.2,X),L,X,2(I3,X),f5.2,X,I8,X,I1)
+!      write(12, *) obs(j)
   enddo
   close(12)
 !  !:XXX

enkf_modules_to_load.txt

@@ -3,17 +3,23 @@ set_flags=$-
 #Disabling u and e flags
 set +ue
 
-module purge
-module load intel
-module load impi
-module load MKL
-module load NETCDF/4.4.1.1
+#module purge
+#module load intel
+#module load impi
+#module load MKL
+#module load NETCDF/4.4.1.1
 
 #If u and e flags were set, set it again
 
-if [ $(echo $set_flags | grep -c "u" ) -eq 1 ];then
-   set -u
-fi
-if [ $(echo $set_flags | grep -c "e" ) -eq 1 ];then
-   set -e
-fi
+#if [ $(echo $set_flags | grep -c "u" ) -eq 1 ];then
+#   set -u
+#fi
+#if [ $(echo $set_flags | grep -c "e" ) -eq 1 ];then
+#   set -e
+#fi
+
+
+#====
+
+my_modules
+module load netCDF-Fortran/4.4.4-foss-2016b

sanity_check/p_sanity_check.F90

@@ -1115,23 +1115,27 @@ IF (ierr .NE. nf90_noerr ) THEN
   STOP
 END IF
 
-! 3.C.1.A Sea surface salinity
-! Request variable ID
-ierr = nf90_inq_varid(incid_oce_an_in, csss_m, ivarid)
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file var. inquiry. Abort"
-  WRITE(*,*), TRIM(cfile_analysis_oce)
-  STOP
-END IF
-
-! Get variable
-ierr = nf90_get_var(incid_oce_an_in, ivarid, sss_m_an)
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file get. Abort"
-  WRITE(*,*), TRIM(cfile_analysis_oce)
-  STOP
-END IF
 
+! The following lines commented out because SSS does not appear anymore
+! as restart variable (2020)
+
+!! 3.C.1.A Sea surface salinity
+!! Request variable ID
+!ierr = nf90_inq_varid(incid_oce_an_in, csss_m, ivarid)
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file var. inquiry. Abort"
+!  WRITE(*,*), TRIM(cfile_analysis_oce)
+!  STOP
+!END IF
+!
+!! Get variable
+!ierr = nf90_get_var(incid_oce_an_in, ivarid, sss_m_an)
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file get. Abort"
+!  WRITE(*,*), TRIM(cfile_analysis_oce)
+!  STOP
+!END IF
+!
 ! 3.C.1.B Salinity
 ! Request variable ID
 ierr = nf90_inq_varid(incid_oce_an_in, csn, ivarid)
@@ -1149,22 +1153,24 @@ IF (ierr .NE. nf90_noerr ) THEN
   STOP
 END IF
 
-! 3.C.1.C Sea surface temperature
-! Request variable ID
-ierr = nf90_inq_varid(incid_oce_an_in, csst_m, ivarid)
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file var. inquiry. Abort"
-  WRITE(*,*), TRIM(cfile_analysis_oce)
-  STOP
-END IF
-
-! Get variable
-ierr = nf90_get_var(incid_oce_an_in, ivarid, sst_m_an)
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file get. Abort"
-  WRITE(*,*), TRIM(cfile_analysis_oce)
-  STOP
-END IF
+! The following lines commented out because SST does not longer
+! appear as restart file (2020)
+!! 3.C.1.C Sea surface temperature
+!! Request variable ID
+!ierr = nf90_inq_varid(incid_oce_an_in, csst_m, ivarid)
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file var. inquiry. Abort"
+!  WRITE(*,*), TRIM(cfile_analysis_oce)
+!  STOP
+!END IF
+!
+!! Get variable
+!ierr = nf90_get_var(incid_oce_an_in, ivarid, sst_m_an)
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file get. Abort"
+!  WRITE(*,*), TRIM(cfile_analysis_oce)
+!  STOP
+!END IF
 
 ! 3.C.1.D. Temperature
 ierr = nf90_inq_varid(incid_oce_an_in, ctn, ivarid)
@@ -1246,37 +1252,40 @@ IF (ierr .NE. nf90_noerr ) THEN
   STOP
 END IF
 
-! 3.C.1.I SSU-velocity ("ssu_m")
-ierr = nf90_inq_varid(incid_oce_an_in, cssu_m, ivarid)
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file var. inquiry. Abort"
-  WRITE(*,*), TRIM(cfile_analysis_oce)
-  STOP
-END IF
-  
-! Get variable
-ierr = nf90_get_var(incid_oce_an_in, ivarid, ssu_m_an)
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file get. Abort"
-  WRITE(*,*), TRIM(cfile_analysis_oce)
-  STOP
-END IF
-
-! 3.C.1.J SSV-velocity ("ssv_m")
-ierr = nf90_inq_varid(incid_oce_an_in, cssv_m, ivarid)
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file var. inquiry. Abort"
-  WRITE(*,*), TRIM(cfile_analysis_oce)
-  STOP
-END IF
-  
-! Get variable
-ierr = nf90_get_var(incid_oce_an_in, ivarid, ssv_m_an)
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file get. Abort"
-  WRITE(*,*), TRIM(cfile_analysis_oce)
-  STOP
-END IF
+! The following lines commented out because surface velocities no longer
+! appear as restart variables (2020)
+
+!! 3.C.1.I SSU-velocity ("ssu_m")
+!ierr = nf90_inq_varid(incid_oce_an_in, cssu_m, ivarid)
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file var. inquiry. Abort"
+!  WRITE(*,*), TRIM(cfile_analysis_oce)
+!  STOP
+!END IF
+!  
+!! Get variable
+!ierr = nf90_get_var(incid_oce_an_in, ivarid, ssu_m_an)
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file get. Abort"
+!  WRITE(*,*), TRIM(cfile_analysis_oce)
+!  STOP
+!END IF
+!
+!! 3.C.1.J SSV-velocity ("ssv_m")
+!ierr = nf90_inq_varid(incid_oce_an_in, cssv_m, ivarid)
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file var. inquiry. Abort"
+!  WRITE(*,*), TRIM(cfile_analysis_oce)
+!  STOP
+!END IF
+!  
+!! Get variable
+!ierr = nf90_get_var(incid_oce_an_in, ivarid, ssv_m_an)
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file get. Abort"
+!  WRITE(*,*), TRIM(cfile_analysis_oce)
+!  STOP
+!END IF
 
 ! Close analysis
 ierr = nf90_close(incid_oce_an_in);
@@ -1294,22 +1303,26 @@ IF (ierr .NE. nf90_noerr ) THEN
   STOP
 END IF
 
-! 3.C.2.A Sea surface salinity
-! Request variable ID
-ierr = nf90_inq_varid(incid_oce_fc_in, csss_m, ivarid)
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file var. inquiry. Abort"
-  WRITE(*,*), TRIM(cfile_forecast_oce)
-  STOP
-END IF
 
-! Get variable
-ierr = nf90_get_var(incid_oce_fc_in, ivarid, sss_m_fc)
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file get. Abort"
-  WRITE(*,*), TRIM(cfile_forecast_oce)
-  STOP
-END IF
+! Following lines commented out as the variable no longer appears 
+! in the restart files (2020
+
+!! 3.C.2.A Sea surface salinity
+!! Request variable ID
+!ierr = nf90_inq_varid(incid_oce_fc_in, csss_m, ivarid)
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file var. inquiry. Abort"
+!  WRITE(*,*), TRIM(cfile_forecast_oce)
+!  STOP
+!END IF
+!
+!! Get variable
+!ierr = nf90_get_var(incid_oce_fc_in, ivarid, sss_m_fc)
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file get. Abort"
+!  WRITE(*,*), TRIM(cfile_forecast_oce)
+!  STOP
+!END IF
 
 ! 3.C.2.B. Salinity
 ! Request variable ID
@@ -1328,22 +1341,25 @@ IF (ierr .NE. nf90_noerr ) THEN
   STOP
 END IF
 
-! 3.C.2.C Sea surface temperature
-! Request variable ID
-ierr = nf90_inq_varid(incid_oce_fc_in, csst_m, ivarid)
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file var. inquiry. Abort"
-  WRITE(*,*), TRIM(cfile_forecast_oce)
-  STOP
-END IF
+! Following lines commented out as the variable no longer appears
+! in the restart files (2020)
 
-! Get variable
-ierr = nf90_get_var(incid_oce_fc_in, ivarid, sst_m_fc)
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file get. Abort"
-  WRITE(*,*), TRIM(cfile_forecast_oce)
-  STOP
-END IF
+!! 3.C.2.C Sea surface temperature
+!! Request variable ID
+!ierr = nf90_inq_varid(incid_oce_fc_in, csst_m, ivarid)
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file var. inquiry. Abort"
+!  WRITE(*,*), TRIM(cfile_forecast_oce)
+!  STOP
+!END IF
+!
+!! Get variable
+!ierr = nf90_get_var(incid_oce_fc_in, ivarid, sst_m_fc)
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file get. Abort"
+!  WRITE(*,*), TRIM(cfile_forecast_oce)
+!  STOP
+!END IF
 
 ! 3.C.2.D. Temperature
 ! Request variable ID
@@ -1426,37 +1442,40 @@ IF (ierr .NE. nf90_noerr ) THEN
   STOP
 END IF
 
-! 3.C.2.I SSU-velocity ("ssu_m")
-ierr = nf90_inq_varid(incid_oce_fc_in, cssu_m, ivarid)
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file var. inquiry. Abort"
-  WRITE(*,*), TRIM(cfile_forecast_oce)
-  STOP
-END IF
-  
-! Get variable
-ierr = nf90_get_var(incid_oce_fc_in, ivarid, ssu_m_fc)
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file get. Abort"
-  WRITE(*,*), TRIM(cfile_forecast_oce)
-  STOP
-END IF
-
-! 3.C.2.J SSV-velocity ("ssv_m")
-ierr = nf90_inq_varid(incid_oce_fc_in, cssv_m, ivarid)
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file var. inquiry. Abort"
-  WRITE(*,*), TRIM(cfile_forecast_oce)
-  STOP
-END IF
-  
+! Following lines commented out as the variable no longer appears
+! in the restart files (2020)
+
+!! 3.C.2.I SSU-velocity ("ssu_m")
+!ierr = nf90_inq_varid(incid_oce_fc_in, cssu_m, ivarid)
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file var. inquiry. Abort"
+!  WRITE(*,*), TRIM(cfile_forecast_oce)
+!  STOP
+!END IF
+!  
+!! Get variable
+!ierr = nf90_get_var(incid_oce_fc_in, ivarid, ssu_m_fc)
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file get. Abort"
+!  WRITE(*,*), TRIM(cfile_forecast_oce)
+!  STOP
+!END IF
+
+!! 3.C.2.J SSV-velocity ("ssv_m")
+!ierr = nf90_inq_varid(incid_oce_fc_in, cssv_m, ivarid)
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file var. inquiry. Abort"
+!  WRITE(*,*), TRIM(cfile_forecast_oce)
+!  STOP
+!END IF
+!  
 ! Get variable
-ierr = nf90_get_var(incid_oce_fc_in, ivarid, ssv_m_fc)
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file get. Abort"
-  WRITE(*,*), TRIM(cfile_forecast_oce)
-  STOP
-END IF
+!ierr = nf90_get_var(incid_oce_fc_in, ivarid, ssv_m_fc)
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF Ocean file get. Abort"
+!  WRITE(*,*), TRIM(cfile_forecast_oce)
+!  STOP
+!END IF
 
 
 
@@ -2256,22 +2275,25 @@ IF (ierr .NE. nf90_noerr ) THEN
   STOP
 END IF
 
-! B. sss_m (sea surface salinity)
-cvarname= csss_m
-ierr = nf90_inq_varid(incid_oce_an_out, cvarname, ivarid)
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF oce file var. inquiry. Abort"
-  WRITE(*,*), TRIM(cfileout_oce)
-  STOP
-END IF
+! Following lines commented out as the variable no longer appears
+! in the restart files (2020)
 
-! Put variable 
-ierr = nf90_put_var(incid_oce_an_out, ivarid, sss_m_an(:,:))
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF oce file var. put. Abort"
-  WRITE(*,*), TRIM(cfileout_oce)
-  STOP
-END IF
+!! B. sss_m (sea surface salinity)
+!cvarname= csss_m
+!ierr = nf90_inq_varid(incid_oce_an_out, cvarname, ivarid)
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF oce file var. inquiry. Abort"
+!  WRITE(*,*), TRIM(cfileout_oce)
+!  STOP
+!END IF
+
+!! Put variable 
+!ierr = nf90_put_var(incid_oce_an_out, ivarid, sss_m_an(:,:))
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF oce file var. put. Abort"
+!  WRITE(*,*), TRIM(cfileout_oce)
+!  STOP
+!END IF
 
 ! C. tn (temperature)
 cvarname= ctn
@@ -2291,22 +2313,25 @@ IF (ierr .NE. nf90_noerr ) THEN
   STOP
 END IF
 
-! D. sst_m (sea surface temperature)
-cvarname= csst_m
-ierr = nf90_inq_varid(incid_oce_an_out, cvarname, ivarid)
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF oce file var. inquiry. Abort"
-  WRITE(*,*), TRIM(cfileout_oce)
-  STOP
-END IF
+! Following lines commented out as the variable no longer appears
+! in the restart files (2020)
 
-! Put variable 
-ierr = nf90_put_var(incid_oce_an_out, ivarid, sst_m_an(:,:))
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF oce file var. put. Abort"
-  WRITE(*,*), TRIM(cfileout_oce)
-  STOP
-END IF
+!! D. sst_m (sea surface temperature)
+!cvarname= csst_m
+!ierr = nf90_inq_varid(incid_oce_an_out, cvarname, ivarid)
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF oce file var. inquiry. Abort"
+!  WRITE(*,*), TRIM(cfileout_oce)
+!  STOP
+!END IF
+
+!! Put variable 
+!ierr = nf90_put_var(incid_oce_an_out, ivarid, sst_m_an(:,:))
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF oce file var. put. Abort"
+!  WRITE(*,*), TRIM(cfileout_oce)
+!  STOP
+!END IF
 
 ! E. un (sea velocity, "un")
 cvarname= cun
@@ -2377,39 +2402,42 @@ IF (ierr .NE. nf90_noerr ) THEN
   STOP
 END IF
 
-! I. ssu_m (sea surface velocity, "ssu_m")
-cvarname= cssu_m
-ierr = nf90_inq_varid(incid_oce_an_out, cvarname, ivarid)
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF oce file var. inquiry. Abort"
-  WRITE(*,*), TRIM(cfileout_oce)
-  STOP
-END IF
-
-! Put variable 
-ierr = nf90_put_var(incid_oce_an_out, ivarid, ssu_m_an(:,:))
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF oce file var. put. Abort"
-  WRITE(*,*), TRIM(cfileout_oce)
-  STOP
-END IF
-
-! J. ssv_m (sea surface velocity, "ssv_m")
-cvarname= cssv_m
-ierr = nf90_inq_varid(incid_oce_an_out, cvarname, ivarid)
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF oce file var. inquiry. Abort"
-  WRITE(*,*), TRIM(cfileout_oce)
-  STOP
-END IF
+! Following lines commented out as the variable no longer appears
+! in the restart files (2020)
 
-! Put variable 
-ierr = nf90_put_var(incid_oce_an_out, ivarid, ssv_m_an(:,:))
-IF (ierr .NE. nf90_noerr ) THEN
-  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF oce file var. put. Abort"
-  WRITE(*,*), TRIM(cfileout_oce)
-  STOP
-END IF
+!! I. ssu_m (sea surface velocity, "ssu_m")
+!cvarname= cssu_m
+!ierr = nf90_inq_varid(incid_oce_an_out, cvarname, ivarid)
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF oce file var. inquiry. Abort"
+!  WRITE(*,*), TRIM(cfileout_oce)
+!  STOP
+!END IF
+!
+!! Put variable 
+!ierr = nf90_put_var(incid_oce_an_out, ivarid, ssu_m_an(:,:))
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF oce file var. put. Abort"
+!  WRITE(*,*), TRIM(cfileout_oce)
+!  STOP
+!END IF
+!
+!! J. ssv_m (sea surface velocity, "ssv_m")
+!cvarname= cssv_m
+!ierr = nf90_inq_varid(incid_oce_an_out, cvarname, ivarid)
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF oce file var. inquiry. Abort"
+!  WRITE(*,*), TRIM(cfileout_oce)
+!  STOP
+!END IF
+!
+!! Put variable 
+!ierr = nf90_put_var(incid_oce_an_out, ivarid, ssv_m_an(:,:))
+!IF (ierr .NE. nf90_noerr ) THEN
+!  WRITE(*,*), "(sanity_check_LIM3) Error NetCDF oce file var. put. Abort"
+!  WRITE(*,*), TRIM(cfileout_oce)
+!  STOP
+!END IF
 
 ! Close file
 ierr = nf90_close(incid_oce_an_out)