From b8411d7a10e31fe6125aa9981de2dee10f9663cc Mon Sep 17 00:00:00 2001
From: sanAkel <santha.akella@nasa.gov>
Date: Sat, 21 Oct 2023 15:43:19 -0400
Subject: [PATCH 01/20] 1. Add main program: sst_ice_helpers.F90 - to calculate
 daily clim, 2. Adjust CMakeLists, 3. Add scanner/reader of DAO OPS binary
 file.

---
 pre/prepare_ocnExtData/CMakeLists.txt         |   6 ++
 .../daily_clim_SST_FRACI_eight.F90            |  43 ++++++++
 pre/prepare_ocnExtData/sst_ice_helpers.F90    | 100 ++++++++++++++++++
 3 files changed, 149 insertions(+)
 create mode 100644 pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90

diff --git a/pre/prepare_ocnExtData/CMakeLists.txt b/pre/prepare_ocnExtData/CMakeLists.txt
index 19b827c0..63db5721 100644
--- a/pre/prepare_ocnExtData/CMakeLists.txt
+++ b/pre/prepare_ocnExtData/CMakeLists.txt
@@ -19,6 +19,12 @@ ecbuild_add_executable (
   LIBS ${this})
 set_target_properties(sst_sic_eight.x PROPERTIES Fortran_MODULE_DIRECTORY ${include_${this}})
 
+ecbuild_add_executable (
+  TARGET cal_daily_clim.x
+  SOURCES daily_clim_SST_FRACI_eight.F90
+  LIBS ${this})
+set_target_properties(cal_daily_clim.x PROPERTIES Fortran_MODULE_DIRECTORY ${include_${this}})
+
 if (USE_F2PY)
    find_package(F2PY2)
    if (F2PY2_FOUND)
diff --git a/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90 b/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90
new file mode 100644
index 00000000..71ec04a3
--- /dev/null
+++ b/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90
@@ -0,0 +1,43 @@
+!
+PROGRAM daily_clim_SST_FRACI_eight
+!---------------------------------------------------------------------------
+  USE sst_ice_helpers, only: read_bin_SST_ICE, &
+                             write_bin,        &
+                             write_netcdf
+
+  IMPLICIT NONE
+
+  CHARACTER (LEN = 200)   :: sst_file, ice_file
+ 
+! CHARACTER(LEN = 4)      :: Year
+! CHARACTER(LEN = 2)      :: Mon, Day
+
+  INTEGER                 :: nymd_in
+  INTEGER                 :: nymd_out, NLON, NLAT
+  REAL                    :: LON(2880), LAT(1440)
+
+  REAL                    :: sst(1440, 2880), ice(1440, 2880)
+  
+  sst_file = 'dataoceanfile_OSTIA_REYNOLDS_SST.2880x1440.2023.data'
+  ice_file = 'dataoceanfile_OSTIA_REYNOLDS_ICE.2880x1440.2023.data'
+
+  nymd_in  = 20230101
+
+  CALL read_bin_SST_ICE( sst_file, nymd_in, nymd_out, NLON, NLAT, LON, LAT, sst)
+  CALL read_bin_SST_ICE( ice_file, nymd_in, nymd_out, NLON, NLAT, LON, LAT, ice)
+
+  print *, " "
+  print *, "Read: ", sst_file
+  print *, nymd_out
+  print *, " "
+  print *, "Write out the data..."
+
+  CALL write_bin( 2023, 1, 1, 2023, 1, 2, '20230101',  1440, 2880, transpose(sst), transpose(ice))
+  CALL write_netcdf( 2023, 1, 1, '20230101',  1440, 2880, transpose(sst), transpose(ice))
+
+  print *, "Done."
+  print *, " "
+
+!---------------------------------------------------------------------------
+END PROGRAM daily_clim_SST_FRACI_eight
+!
diff --git a/pre/prepare_ocnExtData/sst_ice_helpers.F90 b/pre/prepare_ocnExtData/sst_ice_helpers.F90
index 4304d275..1a408351 100644
--- a/pre/prepare_ocnExtData/sst_ice_helpers.F90
+++ b/pre/prepare_ocnExtData/sst_ice_helpers.F90
@@ -17,6 +17,7 @@ module sst_ice_helpers
 public read_input_quart
 public write_bin
 public write_netcdf
+public read_bin_SST_ICE
 public check
 
 contains
@@ -1150,5 +1151,104 @@ SUBROUTINE write_netcdf( today_year, today_mon, today_day, &
 
     END SUBROUTINE write_netcdf
 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+!
+! Input:  
+!   - File name containing (binary formatted) SST or ICE at (1440, 2880) resolution.
+!   - Date (format): YYYYMMDD (4 digit year, 2 digit month, 2 digit day).
+!
+! Output: 
+!   - Date (= above input date).
+!   - Array of lat (1440), lon (2880), field (1440, 2880) for user input date; field = SST or ICE.
+!
+! .....................................................................
+   SUBROUTINE read_bin_SST_ICE( fileName, nymd_in, nymd_out, NLON, NLAT, LON, LAT, bcs_field)
+
+   IMPLICIT NONE
+
+   CHARACTER (LEN = *), INTENT(IN)  :: fileName
+   INTEGER,             INTENT(IN)  :: nymd_in
+
+   INTEGER,             INTENT(OUT) :: nymd_out
+   INTEGER,             INTENT(OUT) :: NLON
+   INTEGER,             INTENT(OUT) :: NLAT
+
+   REAL,                INTENT(OUT) :: LON(2880)
+   REAL,                INTENT(OUT) :: LAT(1440)
+   REAL,                INTENT(OUT) :: bcs_field(1440, 2880)
+
+! LOCAL VARS
+   real    year1,month1,day1,hour1,min1,sec1
+   real    year2,month2,day2,hour2,min2,sec2
+   real    dum1,dum2
+   real    dLat,dLon
+   real    field(2880,1440)
+
+   integer nymd1,nhms1
+   integer nymd2,nhms2
+   
+   integer rc, ix
+!  ....................................................................
+
+      print *, 'Input date: ', nymd_in
+
+      open (10,file=fileName,form='unformatted',access='sequential', STATUS = 'old')
+!     ....................................................................
+      rc = 0
+      do while (rc.eq.0)
+         ! read header
+         read (10,iostat=rc)  year1,month1,day1,hour1,min1,sec1,                            &
+                              year2,month2,day2,hour2,min2,sec2,dum1,dum2
+         if( rc.eq.0 ) then
+           read (10,iostat=rc)  field
+
+           NLON  = nint(dum1)
+           NLAT  = nint(dum2)
+           IF( (NLON /= 2880) .or. (NLAT /= 1440)) THEN
+              PRINT *, 'ERROR in LAT/LON dimension in file: ', fileName, 'NLON=', NLON, 'NLAT=', NLAT
+              EXIT
+           END IF
+           dLat = 180./NLAT
+           dLon = 360./NLON
+
+           ! start date
+           nymd1 = nint( year1*10000 )  + nint (month1*100)  +  nint( day1 )
+           nhms1 = nint( hour1*10000 )  + nint (  min1*100)  +  nint( sec1 )
+
+           ! end date
+           nymd2 = nint( year2*10000 )  + nint (month2*100)  +  nint( day2 )
+           nhms2 = nint( hour2*10000 )  + nint (  min2*100)  +  nint( sec2 )
+         end if
+
+         !print *, nymd1, nymd2
+
+         if (nymd1 .eq. nymd_in) then ! Found the date for which data is requested.
+            nymd_out = nymd1
+
+            print *, 'Yay! Found data for: ', nymd_out
+!!          print *, bcs_field
+
+            LAT(1) = -90. + dLat/2.
+            DO ix = 2, NLAT
+              LAT(ix) = LAT(ix-1) + dLat
+            END DO
+
+            LON(1) = -180. + dLon/2.
+            DO ix = 2, NLON
+              LON(ix) = LON(ix-1) + dLon
+            END DO
+
+            bcs_field = TRANSPOSE(field)
+            close(10)
+            RETURN 
+         else
+            print *, "Date in sequential file: ", nymd1, "Continue looking: ", nymd_in
+         end if
+      end do 
+!     ....................................................................
+      close(10)
+!---------------------------------------------------------------------------
+
+    END SUBROUTINE read_bin_SST_ICE
+!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
 
 end module sst_ice_helpers

From e02b65dd4611a107dbb68c2b0230682a55997309 Mon Sep 17 00:00:00 2001
From: sanAkel <santha.akella@nasa.gov>
Date: Sun, 22 Oct 2023 22:23:48 -0400
Subject: [PATCH 02/20] to generate conservative regridding weights to map
 original OSTIA data to GEOS 1/8-deg grid: use xesmf

---
 .../gen_wts_regrid_ostia_geos.py              | 22 +++++++++++++++++++
 1 file changed, 22 insertions(+)
 create mode 100755 pre/prepare_ocnExtData/gen_wts_regrid_ostia_geos.py

diff --git a/pre/prepare_ocnExtData/gen_wts_regrid_ostia_geos.py b/pre/prepare_ocnExtData/gen_wts_regrid_ostia_geos.py
new file mode 100755
index 00000000..c3456c3e
--- /dev/null
+++ b/pre/prepare_ocnExtData/gen_wts_regrid_ostia_geos.py
@@ -0,0 +1,22 @@
+#!/usr/bin/env python3
+
+import xarray as xr
+import xesmf as xe
+
+ostia_path = "/discover/nobackup/sakella/projects/bcs_pert/data/data_OSTIA"
+geos_path  = "/discover/nobackup/sakella/projects/bcs_pert/tests"
+
+ds_ostia = xr.open_dataset(ostia_path + "/" + "20230101-UKMO-L4HRfnd-GLOB-v01-fv02-OSTIA.nc").squeeze()
+ds_geos  = xr.open_dataset(geos_path  + "/" + "sst_fraci_20230101.nc4").squeeze()
+
+print("Generating weights...")
+regridder = xe.Regridder(ds_ostia, ds_geos, "conservative", periodic=True)
+print("Done.")
+
+regridder.filename = 'wts_cons_ostia_to_geos_eight.nc'
+regridder.to_netcdf()
+print("Saved weights for reuse to: ", regridder.filename)
+
+# use generated weights
+# regridder = xe.Regridder(ds_ostia, ds_geos, 'conservative', weights='wts_cons_ostia_to_geos_eight.nc')
+#ds_ostia_on_geos_grid = regridder(ds_ostia, keep_attrs=True)

From 407873054ece0c83396d66864b1062894e26cb2f Mon Sep 17 00:00:00 2001
From: sanAkel <santha.akella@nasa.gov>
Date: Sun, 22 Oct 2023 22:24:58 -0400
Subject: [PATCH 03/20] Plot difference between conservatively remapped and
 binned (GEOS FP binary file)

---
 ...plot_diff_ostia_native_regridded_binned.py | 48 +++++++++++++++++++
 1 file changed, 48 insertions(+)
 create mode 100755 pre/prepare_ocnExtData/plot_diff_ostia_native_regridded_binned.py

diff --git a/pre/prepare_ocnExtData/plot_diff_ostia_native_regridded_binned.py b/pre/prepare_ocnExtData/plot_diff_ostia_native_regridded_binned.py
new file mode 100755
index 00000000..4d569cf8
--- /dev/null
+++ b/pre/prepare_ocnExtData/plot_diff_ostia_native_regridded_binned.py
@@ -0,0 +1,48 @@
+#!/usr/bin/env python3
+
+import xarray as xr
+import xesmf as xe
+
+import matplotlib.pyplot as plt
+
+import cartopy.crs as ccrs
+import cartopy.feature as cfeature
+import sys
+
+ostia_path = "/discover/nobackup/sakella/projects/bcs_pert/data/data_OSTIA"
+geos_path  = "/discover/nobackup/sakella/projects/bcs_pert/tests"
+
+date_str = sys.argv[1]#'20230101'
+ds_ostia = xr.open_dataset(ostia_path + "/" + date_str + "-UKMO-L4HRfnd-GLOB-v01-fv02-OSTIA.nc").squeeze()
+ds_geos  = xr.open_dataset(geos_path  + "/" + "sst_fraci_" + date_str + ".nc4").squeeze()
+
+regridder = xe.Regridder(ds_ostia, ds_geos, 'conservative', weights='wts_cons_ostia_to_geos_eight.nc')
+ds_ostia_on_geos_grid = regridder(ds_ostia, keep_attrs=True)
+#
+
+fig = plt.figure(figsize=(16,4))
+
+ax1 = plt.subplot(121, projection=ccrs.PlateCarree( central_longitude=0))
+ax1.coastlines()
+ax1.add_feature(cfeature.LAND,  facecolor='white', zorder=100)
+ax1.add_feature(cfeature.LAKES, edgecolor='gray',  zorder=50)
+im1 = ax1.pcolormesh(ds_geos['lon'], ds_geos['lat'], ds_ostia_on_geos_grid.analysed_sst-ds_geos.SST,  transform=ccrs.PlateCarree(), cmap=plt.cm.bwr, vmin=-0.25, vmax=0.25)
+gl = ax1.gridlines(crs=ccrs.PlateCarree(), linewidth=1, color='0.85', alpha=0.5, linestyle='-', draw_labels=True)
+gl.top_labels = False
+cbar=fig.colorbar(im1, extend='both', shrink=0.5, ax=ax1, pad=0.1)
+cbar.set_label(r'diff SST')
+ax1.set_title(r'%s'%(date_str))
+
+ax2 = plt.subplot(122, projection=ccrs.PlateCarree())
+ax2.coastlines()
+ax2.add_feature(cfeature.LAND,  facecolor='white', zorder=100)
+ax2.add_feature(cfeature.LAKES, edgecolor='gray',  zorder=50)
+im2 = ax2.pcolormesh(ds_geos['lon'], ds_geos['lat'], ds_ostia_on_geos_grid.sea_ice_fraction-ds_geos.FRACI,  transform=ccrs.PlateCarree(), cmap=plt.cm.bwr, vmin=-0.1, vmax=0.1)
+gl = ax2.gridlines(crs=ccrs.PlateCarree(), linewidth=1, color='0.85', alpha=0.5, linestyle='-', draw_labels=True)
+gl.top_labels = False
+cbar=fig.colorbar(im2, extend='both', shrink=0.5, ax=ax2, pad=0.1)
+cbar.set_label(r'diff FRACI')
+ax2.set_title(r'%s'%(date_str))
+
+fName_fig = 'consr_binned_sst_ice_diff_' + date_str + '.png'
+plt.savefig(fName_fig, dpi=100)

From b8c592bfb03114b3e955a1867ce7eae6ce43272f Mon Sep 17 00:00:00 2001
From: sanAkel <santha.akella@nasa.gov>
Date: Sun, 22 Oct 2023 22:27:42 -0400
Subject: [PATCH 04/20] Add utility to extract SST and Ice concentration for
 any given day: extract_daily_sst_ice.x 2023 1 5

---
 pre/prepare_ocnExtData/CMakeLists.txt         |  6 ++
 .../daily_clim_SST_FRACI_eight.F90            |  2 -
 .../extract_day_SST_FRACI_eight.F90           | 71 +++++++++++++++++++
 3 files changed, 77 insertions(+), 2 deletions(-)
 create mode 100644 pre/prepare_ocnExtData/extract_day_SST_FRACI_eight.F90

diff --git a/pre/prepare_ocnExtData/CMakeLists.txt b/pre/prepare_ocnExtData/CMakeLists.txt
index 63db5721..c6b17def 100644
--- a/pre/prepare_ocnExtData/CMakeLists.txt
+++ b/pre/prepare_ocnExtData/CMakeLists.txt
@@ -25,6 +25,12 @@ ecbuild_add_executable (
   LIBS ${this})
 set_target_properties(cal_daily_clim.x PROPERTIES Fortran_MODULE_DIRECTORY ${include_${this}})
 
+ecbuild_add_executable (
+  TARGET extract_daily_sst_ice.x
+  SOURCES extract_day_SST_FRACI_eight.F90
+  LIBS ${this})
+set_target_properties(extract_daily_sst_ice.x PROPERTIES Fortran_MODULE_DIRECTORY ${include_${this}})
+
 if (USE_F2PY)
    find_package(F2PY2)
    if (F2PY2_FOUND)
diff --git a/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90 b/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90
index 71ec04a3..fc2e9092 100644
--- a/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90
+++ b/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90
@@ -1,4 +1,3 @@
-!
 PROGRAM daily_clim_SST_FRACI_eight
 !---------------------------------------------------------------------------
   USE sst_ice_helpers, only: read_bin_SST_ICE, &
@@ -40,4 +39,3 @@ PROGRAM daily_clim_SST_FRACI_eight
 
 !---------------------------------------------------------------------------
 END PROGRAM daily_clim_SST_FRACI_eight
-!
diff --git a/pre/prepare_ocnExtData/extract_day_SST_FRACI_eight.F90 b/pre/prepare_ocnExtData/extract_day_SST_FRACI_eight.F90
new file mode 100644
index 00000000..2d0ff482
--- /dev/null
+++ b/pre/prepare_ocnExtData/extract_day_SST_FRACI_eight.F90
@@ -0,0 +1,71 @@
+PROGRAM extract_day_SST_FRACI_eight
+!---------------------------------------------------------------------------
+  USE sst_ice_helpers, only: read_bin_SST_ICE, &
+!                            write_bin,        &
+                             write_netcdf
+
+  IMPLICIT NONE
+
+  integer, parameter :: nlon = 2880
+  integer, parameter :: nlat = 1440
+  real               :: lon(nlon), lat(nlat)
+  real               :: sst(nlat, nlon), ice(nlat, nlon)
+
+  integer               :: extract_date
+  integer               :: date_out, nlon_out, nlat_out
+
+  character (len=8)     :: date_str, format_str = "(I8)"
+  character (len = 200) :: sst_file, ice_file
+
+  character (len=4) :: arg
+  integer :: numArgs
+  integer :: year_s, month_s, day_s
+!----
+
+  numArgs = iargc()!nargs() for ifort
+  if (numArgs < 3)  then
+    print *, "Need 3 inputs, in following format. Try again!"
+    print *, "year month day"
+    STOP
+  end if
+
+  !-- read user inputs
+  call getarg(1, arg)
+  read(arg, '(I14)') year_s ! default
+  call getarg(2, arg)
+  read(arg, '(I14)') month_s
+  call getarg(3, arg)
+  read(arg, '(I14)') day_s
+  
+  !-- form a string of input date
+  extract_date = year_s*10000+month_s*100+day_s
+  write (date_str, format_str) extract_date
+
+  sst_file = 'dataoceanfile_OSTIA_REYNOLDS_SST.2880x1440.'//trim(date_str(1:4))//'.data'
+  ice_file = 'dataoceanfile_OSTIA_REYNOLDS_ICE.2880x1440.'//trim(date_str(1:4))//'.data'
+
+  print *, " "
+  print *, "Read: ", sst_file, ice_file
+  print *, " "
+
+  ! read sst
+  CALL read_bin_SST_ICE( sst_file, extract_date, date_out, nlon_out, nlat_out, lon, lat, sst)
+  ! read sea ice
+  CALL read_bin_SST_ICE( ice_file, extract_date, date_out, nlon_out, nlat_out, lon, lat, ice)
+
+  print *, date_out
+  print *, " "
+  print *, "Write out the data..."
+
+  ! binary write
+  ! CALL write_bin( 2023, 1, 1, 2023, 1, 2, '20230101',  1440, 2880, transpose(sst), transpose(ice))
+
+  ! netcdf
+  !CALL write_netcdf( 2023, 1, 1, '20230101',  1440, 2880, transpose(sst), transpose(ice))
+  CALL write_netcdf( year_s, month_s, day_s, date_str, nlat, nlon, transpose(sst), transpose(ice))
+
+  print *, "Done."
+  print *, " "
+
+!---------------------------------------------------------------------------
+END PROGRAM extract_day_SST_FRACI_eight

From c668923d4968372708003fbbdce4db5befecd4cc Mon Sep 17 00:00:00 2001
From: sanAkel <santha.akella@nasa.gov>
Date: Sun, 22 Oct 2023 22:33:03 -0400
Subject: [PATCH 05/20] get rid of commented line

---
 pre/prepare_ocnExtData/extract_day_SST_FRACI_eight.F90 | 1 -
 1 file changed, 1 deletion(-)

diff --git a/pre/prepare_ocnExtData/extract_day_SST_FRACI_eight.F90 b/pre/prepare_ocnExtData/extract_day_SST_FRACI_eight.F90
index 2d0ff482..0560fa6c 100644
--- a/pre/prepare_ocnExtData/extract_day_SST_FRACI_eight.F90
+++ b/pre/prepare_ocnExtData/extract_day_SST_FRACI_eight.F90
@@ -61,7 +61,6 @@ PROGRAM extract_day_SST_FRACI_eight
   ! CALL write_bin( 2023, 1, 1, 2023, 1, 2, '20230101',  1440, 2880, transpose(sst), transpose(ice))
 
   ! netcdf
-  !CALL write_netcdf( 2023, 1, 1, '20230101',  1440, 2880, transpose(sst), transpose(ice))
   CALL write_netcdf( year_s, month_s, day_s, date_str, nlat, nlon, transpose(sst), transpose(ice))
 
   print *, "Done."

From 026c86edeb17e2b080ef07262e2cac07af1ac366 Mon Sep 17 00:00:00 2001
From: sanAkel <santha.akella@nasa.gov>
Date: Wed, 8 Nov 2023 20:49:32 -0500
Subject: [PATCH 06/20] Added time ticker from @rtodling

---
 pre/prepare_ocnExtData/CMakeLists.txt         |   1 +
 .../daily_clim_SST_FRACI_eight.F90            |   9 ++
 pre/prepare_ocnExtData/hermes_m_tick.f90      | 146 ++++++++++++++++++
 3 files changed, 156 insertions(+)
 create mode 100644 pre/prepare_ocnExtData/hermes_m_tick.f90

diff --git a/pre/prepare_ocnExtData/CMakeLists.txt b/pre/prepare_ocnExtData/CMakeLists.txt
index c6b17def..9ef9caa8 100644
--- a/pre/prepare_ocnExtData/CMakeLists.txt
+++ b/pre/prepare_ocnExtData/CMakeLists.txt
@@ -2,6 +2,7 @@ esma_set_this()
 
 set ( srcs
   sst_ice_helpers.F90
+  hermes_m_tick.f90
   )
 
 esma_add_library (${this} SRCS ${srcs} DEPENDENCIES MAPL NetCDF::NetCDF_Fortran)
diff --git a/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90 b/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90
index fc2e9092..a0d05b82 100644
--- a/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90
+++ b/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90
@@ -1,5 +1,6 @@
 PROGRAM daily_clim_SST_FRACI_eight
 !---------------------------------------------------------------------------
+  USE m_tick, only: INCYMD
   USE sst_ice_helpers, only: read_bin_SST_ICE, &
                              write_bin,        &
                              write_netcdf
@@ -22,6 +23,14 @@ PROGRAM daily_clim_SST_FRACI_eight
 
   nymd_in  = 20230101
 
+  print *, "----------------"
+  print *, INCYMD (nymd_in,1)
+  print *, INCYMD (nymd_in,2)
+  print *, INCYMD (nymd_in,3)
+  print *, INCYMD (nymd_in,4)
+  print *, INCYMD (nymd_in,5)
+  print *, "----------------"
+
   CALL read_bin_SST_ICE( sst_file, nymd_in, nymd_out, NLON, NLAT, LON, LAT, sst)
   CALL read_bin_SST_ICE( ice_file, nymd_in, nymd_out, NLON, NLAT, LON, LAT, ice)
 
diff --git a/pre/prepare_ocnExtData/hermes_m_tick.f90 b/pre/prepare_ocnExtData/hermes_m_tick.f90
new file mode 100644
index 00000000..dd8b741b
--- /dev/null
+++ b/pre/prepare_ocnExtData/hermes_m_tick.f90
@@ -0,0 +1,146 @@
+! --
+! Copied from and renamed
+! GMAO_Shared/main/GMAO_hermes/m_tick.f90
+! --
+
+module m_tick
+private
+public tick
+public INCYMD
+contains
+      subroutine tick (nymd,nhms,ndt)
+!***********************************************************************
+!  Purpose
+!     Tick the Date (nymd) and Time (nhms) by NDT (seconds)
+!
+!***********************************************************************
+!*                  GODDARD LABORATORY FOR ATMOSPHERES                 *
+!***********************************************************************
+
+      IF(NDT.NE.0) THEN
+      NSEC = NSECF(NHMS) + NDT
+
+      IF (NSEC.GT.86400)  THEN
+      DO WHILE (NSEC.GT.86400)
+      NSEC = NSEC - 86400
+      NYMD = INCYMD (NYMD,1)
+      ENDDO
+      ENDIF  
+              
+      IF (NSEC.EQ.86400)  THEN
+      NSEC = 0
+      NYMD = INCYMD (NYMD,1)
+      ENDIF  
+              
+      IF (NSEC.LT.00000)  THEN
+      DO WHILE (NSEC.LT.0)
+      NSEC = 86400 + NSEC
+      NYMD = INCYMD (NYMD,-1)
+      ENDDO
+      ENDIF  
+              
+      NHMS = NHMSF (NSEC)
+      ENDIF
+
+      RETURN
+      END
+
+      FUNCTION INCYMD (NYMD,M)
+!***********************************************************************
+!  PURPOSE
+!     INCYMD:  NYMD CHANGED BY ONE DAY
+!     MODYMD:  NYMD CONVERTED TO JULIAN DATE
+!  DESCRIPTION OF PARAMETERS
+!     NYMD     CURRENT DATE IN YYMMDD FORMAT
+!     M        +/- 1 (DAY ADJUSTMENT)
+!
+!***********************************************************************
+!*                  GODDARD LABORATORY FOR ATMOSPHERES                 *
+!***********************************************************************
+
+      INTEGER NDPM(12)
+      DATA    NDPM /31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31/
+      LOGICAL LEAP
+      DATA    NY00     / 1900 /
+      LEAP(NY) = MOD(NY,4).EQ.0 .AND. (NY.NE.0 .OR. MOD(NY00,400).EQ.0)
+
+!***********************************************************************
+!
+      NY = NYMD / 10000
+      NM = MOD(NYMD,10000) / 100
+      ND = MOD(NYMD,100) + M
+
+      IF (ND.EQ.0) THEN
+      NM = NM - 1
+      IF (NM.EQ.0) THEN
+          NM = 12
+          NY = NY - 1
+      ENDIF
+      ND = NDPM(NM)
+      IF (NM.EQ.2 .AND. LEAP(NY))  ND = 29
+      ENDIF
+
+      IF (ND.EQ.29 .AND. NM.EQ.2 .AND. LEAP(NY))  GO TO 20
+
+      IF (ND.GT.NDPM(NM)) THEN
+      ND = 1
+      NM = NM + 1
+      IF (NM.GT.12) THEN
+          NM = 1
+          NY = NY + 1
+      ENDIF
+      ENDIF
+
+   20 CONTINUE
+      INCYMD = NY*10000 + NM*100 + ND
+      RETURN
+
+!***********************************************************************
+!                      E N T R Y    M O D Y M D
+!***********************************************************************
+
+      ENTRY MODYMD (NYMD)
+      NY = NYMD / 10000
+      NM = MOD(NYMD,10000) / 100
+      ND = MOD(NYMD,100)
+
+   40 CONTINUE
+      IF (NM.LE.1)  GO TO 60
+      NM = NM - 1
+      ND = ND + NDPM(NM)
+      IF (NM.EQ.2 .AND. LEAP(NY))  ND = ND + 1
+      GO TO 40
+
+   60 CONTINUE
+      MODYMD = ND
+      RETURN
+      END
+
+      function nhmsf (nsec)
+!***********************************************************************
+!  Purpose
+!     Converts Total Seconds to NHMS format
+!
+!***********************************************************************
+!*                  GODDARD LABORATORY FOR ATMOSPHERES                 *
+!***********************************************************************
+      implicit none
+      integer  nhmsf, nsec
+      nhmsf =  nsec/3600*10000 + mod(nsec,3600)/60*100 + mod(nsec,60)
+      return
+      end
+
+      function nsecf (nhms)
+!***********************************************************************
+!  Purpose
+!     Converts NHMS format to Total Seconds
+!
+!***********************************************************************
+!*                  GODDARD LABORATORY FOR ATMOSPHERES                 *
+!***********************************************************************
+      implicit none
+      integer  nhms, nsecf
+      nsecf =  nhms/10000*3600 + mod(nhms,10000)/100*60 + mod(nhms,100)
+      return
+      end
+end module m_tick

From 78bc8c12c2ba78f2bf54c83c12365f88cc7059e9 Mon Sep 17 00:00:00 2001
From: sanAkel <santha.akella@nasa.gov>
Date: Tue, 19 Dec 2023 22:25:42 -0500
Subject: [PATCH 07/20] Getting there! Add standard deviation calculation and
 user input args

---
 .../daily_clim_SST_FRACI_eight.F90            | 81 +++++++++++++------
 pre/prepare_ocnExtData/sst_ice_helpers.F90    | 10 +--
 2 files changed, 63 insertions(+), 28 deletions(-)

diff --git a/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90 b/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90
index a0d05b82..b02498e6 100644
--- a/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90
+++ b/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90
@@ -7,42 +7,77 @@ PROGRAM daily_clim_SST_FRACI_eight
 
   IMPLICIT NONE
 
-  CHARACTER (LEN = 200)   :: sst_file, ice_file
+  INTEGER                 :: month, day
+  INTEGER                 :: clim_year, year_start, year_end
+  CHARACTER (LEN = 200)   :: data_path
+  CHARACTER (LEN = 200)   :: sst_file_pref, ice_file_pref
+  CHARACTER (LEN = 200)   :: sst_file_suff, ice_file_suff
  
-! CHARACTER(LEN = 4)      :: Year
-! CHARACTER(LEN = 2)      :: Mon, Day
-
   INTEGER                 :: nymd_in
   INTEGER                 :: nymd_out, NLON, NLAT
   REAL                    :: LON(2880), LAT(1440)
 
-  REAL                    :: sst(1440, 2880), ice(1440, 2880)
+  INTEGER                 :: year
+  INTEGER                 :: clim_tom, clim_tom_year, clim_tom_mon, clim_tom_day
+  CHARACTER (LEN = 200)   :: year_str, sst_file, ice_file
+  CHARACTER (LEN = 8)     :: clim_date_str
+  REAL                    :: mean_sst(1440, 2880), mean_ice(1440, 2880)
+  REAL                    ::      sst(1440, 2880),      ice(1440, 2880)
   
-  sst_file = 'dataoceanfile_OSTIA_REYNOLDS_SST.2880x1440.2023.data'
-  ice_file = 'dataoceanfile_OSTIA_REYNOLDS_ICE.2880x1440.2023.data'
+! should be user inputs
+  clim_year  = 1 ! Can't be "0" since GEOS/ESMF does not like that!
+  year_start = 2007
+  year_end   = 2023
+  month      = 1
+  day        = 2
+  data_path= '/discover/nobackup/projects/gmao/share/dao_ops/fvInput/g5gcm/bcs/realtime/OSTIA_REYNOLDS/2880x1440/'
+  sst_file_pref = 'dataoceanfile_OSTIA_REYNOLDS_SST.2880x1440.'
+  ice_file_pref = 'dataoceanfile_OSTIA_REYNOLDS_ICE.2880x1440.'
+  sst_file_suff = '.data'
+  ice_file_suff = '.data'
 
-  nymd_in  = 20230101
+  print *, " "
+! mean over year_start: year_end
+  do year = year_start, year_end
+    write (year_str,'(I4.0)') year
+    sst_file = trim(data_path)//"/"//trim(sst_file_pref)//trim(year_str)//trim(sst_file_suff)
+    ice_file = trim(data_path)//"/"//trim(ice_file_pref)//trim(year_str)//trim(ice_file_suff)
 
-  print *, "----------------"
-  print *, INCYMD (nymd_in,1)
-  print *, INCYMD (nymd_in,2)
-  print *, INCYMD (nymd_in,3)
-  print *, INCYMD (nymd_in,4)
-  print *, INCYMD (nymd_in,5)
-  print *, "----------------"
+    nymd_in = year*10000 + month*100 + day
+    !print *, year, year_str,nymd_in
+    !print *, sst_file, ice_file
 
-  CALL read_bin_SST_ICE( sst_file, nymd_in, nymd_out, NLON, NLAT, LON, LAT, sst)
-  CALL read_bin_SST_ICE( ice_file, nymd_in, nymd_out, NLON, NLAT, LON, LAT, ice)
+    CALL read_bin_SST_ICE( sst_file, nymd_in, nymd_out, NLON, NLAT, LON, LAT, sst)
+    CALL read_bin_SST_ICE( ice_file, nymd_in, nymd_out, NLON, NLAT, LON, LAT, ice)
+    print *, "Read SST and FRACI for: ", nymd_out
 
+    if (year == year_start) then
+      mean_sst = sst
+      mean_ice = ice
+    else
+      mean_sst = mean_sst + sst 
+      mean_ice = mean_ice + ice
+    endif
+    !print *, " "
+  end do
   print *, " "
-  print *, "Read: ", sst_file
-  print *, nymd_out
+  
+  print *, "Number of years= ", year_end-year_start+1
+  mean_sst = mean_sst/(year_end-year_start+1)
+  mean_ice = mean_ice/(year_end-year_start+1)
+! --
+
   print *, " "
+  clim_tom = INCYMD( clim_year*10000+month*100+day, 1)
+  clim_tom_year = int(clim_tom/10000)
+  clim_tom_mon  = int((clim_tom-clim_tom_year*10000)/100)
+  clim_tom_day = clim_tom - (clim_tom_year*10000+clim_tom_mon*100)
+  !print *, "Tomorrow day:", clim_tom_year, clim_tom_mon, clim_tom_day
+  write (clim_date_str,'(I0.8)') clim_year*10000+month*100+day
+  !print *, clim_date_str
   print *, "Write out the data..."
-
-  CALL write_bin( 2023, 1, 1, 2023, 1, 2, '20230101',  1440, 2880, transpose(sst), transpose(ice))
-  CALL write_netcdf( 2023, 1, 1, '20230101',  1440, 2880, transpose(sst), transpose(ice))
-
+  CALL write_bin(    clim_year, month, day, clim_tom_year, clim_tom_mon, clim_tom_day, clim_date_str, 1440, 2880, transpose(mean_sst), transpose(mean_ice))
+  CALL write_netcdf( clim_year, month, day, clim_date_str, 1440, 2880, transpose(mean_sst), transpose(mean_ice))
   print *, "Done."
   print *, " "
 
diff --git a/pre/prepare_ocnExtData/sst_ice_helpers.F90 b/pre/prepare_ocnExtData/sst_ice_helpers.F90
index 1a408351..ff7b8018 100644
--- a/pre/prepare_ocnExtData/sst_ice_helpers.F90
+++ b/pre/prepare_ocnExtData/sst_ice_helpers.F90
@@ -1080,7 +1080,7 @@ SUBROUTINE write_netcdf( today_year, today_mon, today_day, &
     file_hist = "File created on " // trim(date_creation) // trim(time_creation) // ". In yyyymmddHHMMSS.millisec"
 
     write(DTSTR,10) today_year, today_mon, today_day
- 10 FORMAT(I4,'-', I2, '-', I0.2)
+ 10 FORMAT(I0.4,'-', I0.2, '-', I0.2)
     TIME_UNITS = "days since " // trim(DTSTR) //" 12:00:00" 
 
     dlat = 180./real(nlat); dlon = 360./real(nlon)
@@ -1097,7 +1097,7 @@ SUBROUTINE write_netcdf( today_year, today_mon, today_day, &
 
     sst_out(:,:, 1) = sst; fraci_out(:,:, 1) = ice
 
-    fileName = 'sst_fraci_' // today //'.nc4'
+    fileName = 'sst_fraci_' // today //'.nc'
 
     call check( nf90_create(fileName, nf90_clobber, ncid)) ! Create the file. 
 
@@ -1189,7 +1189,7 @@ SUBROUTINE read_bin_SST_ICE( fileName, nymd_in, nymd_out, NLON, NLAT, LON, LAT,
    integer rc, ix
 !  ....................................................................
 
-      print *, 'Input date: ', nymd_in
+      !print *, 'Input date: ', nymd_in
 
       open (10,file=fileName,form='unformatted',access='sequential', STATUS = 'old')
 !     ....................................................................
@@ -1224,7 +1224,7 @@ SUBROUTINE read_bin_SST_ICE( fileName, nymd_in, nymd_out, NLON, NLAT, LON, LAT,
          if (nymd1 .eq. nymd_in) then ! Found the date for which data is requested.
             nymd_out = nymd1
 
-            print *, 'Yay! Found data for: ', nymd_out
+!           print *, 'Yay! Found data for: ', nymd_out
 !!          print *, bcs_field
 
             LAT(1) = -90. + dLat/2.
@@ -1241,7 +1241,7 @@ SUBROUTINE read_bin_SST_ICE( fileName, nymd_in, nymd_out, NLON, NLAT, LON, LAT,
             close(10)
             RETURN 
          else
-            print *, "Date in sequential file: ", nymd1, "Continue looking: ", nymd_in
+!           print *, "Date in sequential file: ", nymd1, "Continue looking: ", nymd_in
          end if
       end do 
 !     ....................................................................

From 0e001901dc8488751ea8ae433ccf4c299e950441 Mon Sep 17 00:00:00 2001
From: sanAkel <santha.akella@nasa.gov>
Date: Fri, 22 Dec 2023 15:12:28 -0500
Subject: [PATCH 08/20] Working daily climatology production- including leap
 year/day (02/29)

---
 .../daily_clim_SST_FRACI_eight.F90            | 143 ++++++++++++++----
 pre/prepare_ocnExtData/sst_ice_helpers.F90    |  15 +-
 2 files changed, 123 insertions(+), 35 deletions(-)

diff --git a/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90 b/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90
index b02498e6..9d467a6e 100644
--- a/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90
+++ b/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90
@@ -13,7 +13,7 @@ PROGRAM daily_clim_SST_FRACI_eight
   CHARACTER (LEN = 200)   :: sst_file_pref, ice_file_pref
   CHARACTER (LEN = 200)   :: sst_file_suff, ice_file_suff
  
-  INTEGER                 :: nymd_in
+  INTEGER                 :: nymd_in, nYears
   INTEGER                 :: nymd_out, NLON, NLAT
   REAL                    :: LON(2880), LAT(1440)
 
@@ -22,49 +22,128 @@ PROGRAM daily_clim_SST_FRACI_eight
   CHARACTER (LEN = 200)   :: year_str, sst_file, ice_file
   CHARACTER (LEN = 8)     :: clim_date_str
   REAL                    :: mean_sst(1440, 2880), mean_ice(1440, 2880)
+  REAL                    :: sdev_sst(1440, 2880), sdev_ice(1440, 2880)
   REAL                    ::      sst(1440, 2880),      ice(1440, 2880)
+  REAL                    ::     dsst(1440, 2880),     dice(1440, 2880)
+
+  LOGICAL                 :: verbose
+
+  character (len=4) :: arg
+  integer :: numArgs
+  numArgs = iargc()
+
+  verbose = .false.
+
+! Inputs
+  if (numArgs < 4)  then
+    print *, "Need 4 inputs, in following format. Try again!"
+    print *, "start_year end_year month day. Last 2 inputs are for the climatology of MM/DD."
+    STOP
+  else ! read user inputs
+    call getarg(1, arg)
+    read(arg, '(I4)') year_start !2007
+    call getarg(2, arg)
+    read(arg, '(I4)') year_end   !2023
+    call getarg(3, arg)
+    read(arg, '(I2)') month
+    call getarg(4, arg)
+    read(arg, '(I2)') day
+  end if
   
-! should be user inputs
-  clim_year  = 1 ! Can't be "0" since GEOS/ESMF does not like that!
-  year_start = 2007
-  year_end   = 2023
-  month      = 1
-  day        = 2
+! - Set inputs
+  clim_year  = 1 ! Set, can't be "0" since GEOS/ESMF does not like that!
+! GMAO OPS has following fixed
   data_path= '/discover/nobackup/projects/gmao/share/dao_ops/fvInput/g5gcm/bcs/realtime/OSTIA_REYNOLDS/2880x1440/'
   sst_file_pref = 'dataoceanfile_OSTIA_REYNOLDS_SST.2880x1440.'
   ice_file_pref = 'dataoceanfile_OSTIA_REYNOLDS_ICE.2880x1440.'
   sst_file_suff = '.data'
   ice_file_suff = '.data'
 
+! Check input range
+  if (year_start < 2007) then
+    print *, "Invalid starting year: ", year_start, "Must be 2007 or later, try again."
+    STOP
+  endif
+  if ( (month < 1) .or. (month > 12)) then
+    print *, "Invalid month: ", month, "Try again."
+    STOP
+  endif
+  if ( (day < 1) .or. (day > 31)) then
+    print *, "Invalid day: ", day, "Try again."
+    STOP
+  endif
+  if(verbose) print *, "User inputs: ", year_start, year_end, month, day
+!--
+
   print *, " "
-! mean over year_start: year_end
+! mean and std dev over year_start: year_end
+! 1. mean
+  nYears = 0
   do year = year_start, year_end
     write (year_str,'(I4.0)') year
     sst_file = trim(data_path)//"/"//trim(sst_file_pref)//trim(year_str)//trim(sst_file_suff)
     ice_file = trim(data_path)//"/"//trim(ice_file_pref)//trim(year_str)//trim(ice_file_suff)
 
     nymd_in = year*10000 + month*100 + day
-    !print *, year, year_str,nymd_in
-    !print *, sst_file, ice_file
+    if(verbose) then
+      print *, year, year_str, nymd_in
+      print *, sst_file, ice_file
+      print *, nymd_in
+    endif
 
-    CALL read_bin_SST_ICE( sst_file, nymd_in, nymd_out, NLON, NLAT, LON, LAT, sst)
+    CALL read_bin_SST_ICE( sst_file, nymd_in, nymd_out, NLON, NLAT, LON, LAT, sst) ! Read GMAO OPS dataset
     CALL read_bin_SST_ICE( ice_file, nymd_in, nymd_out, NLON, NLAT, LON, LAT, ice)
-    print *, "Read SST and FRACI for: ", nymd_out
-
-    if (year == year_start) then
-      mean_sst = sst
-      mean_ice = ice
-    else
-      mean_sst = mean_sst + sst 
-      mean_ice = mean_ice + ice
+
+    if ( (nymd_out-nymd_in) .eq. 0) then ! Found data for matching dates
+      print *, "Read SST and FRACI for: ", nymd_out
+      if (verbose) print *, "Found data in file, number of years =", nYears
+      nYears = nYears + 1
+
+      if (nYears == 1) then
+        mean_sst = sst
+        mean_ice = ice
+      else
+        mean_sst = mean_sst + sst 
+        mean_ice = mean_ice + ice
+      endif
     endif
-    !print *, " "
   end do
+
   print *, " "
-  
-  print *, "Number of years= ", year_end-year_start+1
-  mean_sst = mean_sst/(year_end-year_start+1)
-  mean_ice = mean_ice/(year_end-year_start+1)
+  print *, "Number of years= ", nYears
+
+  mean_sst = mean_sst/REAL(nYears) ! nYears should be _promoted_ to real, being safe!
+  mean_ice = mean_ice/REAL(nYears)
+
+! 2. sdev
+  nYears = 0
+  do year = year_start, year_end
+
+    write (year_str,'(I4.0)') year
+    sst_file = trim(data_path)//"/"//trim(sst_file_pref)//trim(year_str)//trim(sst_file_suff)
+    ice_file = trim(data_path)//"/"//trim(ice_file_pref)//trim(year_str)//trim(ice_file_suff)
+
+    nymd_in = year*10000 + month*100 + day
+    CALL read_bin_SST_ICE( sst_file, nymd_in, nymd_out, NLON, NLAT, LON, LAT, sst)
+    CALL read_bin_SST_ICE( ice_file, nymd_in, nymd_out, NLON, NLAT, LON, LAT, ice)
+
+    if ( (nymd_out-nymd_in) .eq. 0) then
+      nYears = nYears + 1
+
+      if (nYears == 1) then
+        dsst = (sst - mean_sst)**2.
+        dice = (ice - mean_ice)**2.
+      else
+        dsst = dsst + (sst - mean_sst)**2.
+        dice = dice + (ice - mean_ice)**2.
+      endif
+    endif
+  end do
+
+  print *, "Number of years= ", nYears
+! based on a sample not population
+  sdev_sst = sqrt( dsst/REAL(nYears-1))
+  sdev_ice = sqrt( dice/REAL(nYears-1))
 ! --
 
   print *, " "
@@ -72,12 +151,18 @@ PROGRAM daily_clim_SST_FRACI_eight
   clim_tom_year = int(clim_tom/10000)
   clim_tom_mon  = int((clim_tom-clim_tom_year*10000)/100)
   clim_tom_day = clim_tom - (clim_tom_year*10000+clim_tom_mon*100)
-  !print *, "Tomorrow day:", clim_tom_year, clim_tom_mon, clim_tom_day
+  if(verbose) print *, "Tomorrow day:", clim_tom_year, clim_tom_mon, clim_tom_day
   write (clim_date_str,'(I0.8)') clim_year*10000+month*100+day
-  !print *, clim_date_str
-  print *, "Write out the data..."
-  CALL write_bin(    clim_year, month, day, clim_tom_year, clim_tom_mon, clim_tom_day, clim_date_str, 1440, 2880, transpose(mean_sst), transpose(mean_ice))
-  CALL write_netcdf( clim_year, month, day, clim_date_str, 1440, 2880, transpose(mean_sst), transpose(mean_ice))
+  if(verbose) print *, clim_date_str
+
+! Only write mean in binary format- good enough.
+  print *, "Write output ..."
+  CALL write_bin(    'daily_clim_mean', clim_year, month, day, clim_tom_year, clim_tom_mon, clim_tom_day, clim_date_str, 1440, 2880, transpose(mean_sst), transpose(mean_ice))
+
+! mean
+  CALL write_netcdf( 'daily_clim_mean_sst_fraci', clim_year, month, day, clim_date_str, 1440, 2880, transpose(mean_sst), transpose(mean_ice))
+! std dev
+  CALL write_netcdf( 'daily_clim_sdev_sst_fraci', clim_year, month, day, clim_date_str, 1440, 2880, transpose(sdev_sst), transpose(sdev_ice))
   print *, "Done."
   print *, " "
 
diff --git a/pre/prepare_ocnExtData/sst_ice_helpers.F90 b/pre/prepare_ocnExtData/sst_ice_helpers.F90
index ff7b8018..e0ebc3b7 100644
--- a/pre/prepare_ocnExtData/sst_ice_helpers.F90
+++ b/pre/prepare_ocnExtData/sst_ice_helpers.F90
@@ -975,7 +975,7 @@ END SUBROUTINE check
 !
 ! Write out a binary file with `HEADER` and data, see below for its format.
 !
-    SUBROUTINE write_bin( today_year, today_mon, today_day, &
+    SUBROUTINE write_bin( fileName_pref, today_year, today_mon, today_day, &
                           tomrw_year, tomrw_mon, tomrw_day, &
                           today, nlat, nlon, sst, ice)
 !---------------------------------------------------------------------------
@@ -984,6 +984,8 @@ SUBROUTINE write_bin( today_year, today_mon, today_day, &
     INTEGER, INTENT(IN)     :: today_year, today_mon, today_day, &
                                tomrw_year, tomrw_mon, tomrw_day, &
                                nlat, nlon
+
+    CHARACTER (LEN=*),INTENT(IN) :: fileName_pref
     CHARACTER (LEN=*),INTENT(IN) :: today
 
     REAL, INTENT(IN)        :: sst(nlon, nlat), &
@@ -1002,8 +1004,8 @@ SUBROUTINE write_bin( today_year, today_mon, today_day, &
 !---------------------------------------------------------------------------
 !   Write out SST and ICE concentration data
 !
-    fileName_sst = 'Ostia_sst_' // today //'.bin'
-    fileName_ice = 'Ostia_ice_' // today //'.bin'
+    fileName_sst = fileName_pref // '_sst_' // today //'.bin'
+    fileName_ice = fileName_pref // '_ice_' // today //'.bin'
 
     OPEN (UNIT = 991, FILE = fileName_sst, FORM = 'unformatted', STATUS = 'new')
     OPEN (UNIT = 992, FILE = fileName_ice, FORM = 'unformatted', STATUS = 'new')
@@ -1024,7 +1026,7 @@ END SUBROUTINE write_bin
 !
 ! Write out a netcdf file, see below for its format.
 !
-    SUBROUTINE write_netcdf( today_year, today_mon, today_day, &
+    SUBROUTINE write_netcdf( fileName_pref, today_year, today_mon, today_day, &
                              today, nlat, nlon, sst, ice)
 !---------------------------------------------------------------------------
     IMPLICIT NONE
@@ -1036,7 +1038,8 @@ SUBROUTINE write_netcdf( today_year, today_mon, today_day, &
     REAL, INTENT(IN)        :: sst(nlon, nlat), &
                                ice(nlon, nlat)
 
-    CHARACTER (LEN = 40)    :: fileName
+    CHARACTER (LEN=*),INTENT(IN)   :: fileName_pref
+    CHARACTER (LEN=60)             :: fileName
 
     integer                        :: ncid
     integer, parameter             :: ndims     = 3
@@ -1097,7 +1100,7 @@ SUBROUTINE write_netcdf( today_year, today_mon, today_day, &
 
     sst_out(:,:, 1) = sst; fraci_out(:,:, 1) = ice
 
-    fileName = 'sst_fraci_' // today //'.nc'
+    fileName = fileName_pref // '_' // today //'.nc'
 
     call check( nf90_create(fileName, nf90_clobber, ncid)) ! Create the file. 
 

From c9c463035a38234877a3f698e42d002d77d1d51e Mon Sep 17 00:00:00 2001
From: sanAkel <santha.akella@nasa.gov>
Date: Fri, 29 Dec 2023 04:56:05 -0500
Subject: [PATCH 09/20] remove hard coded file prefix/suffix

---
 .../extract_day_SST_FRACI_eight.F90           | 19 +++++++++++++++----
 1 file changed, 15 insertions(+), 4 deletions(-)

diff --git a/pre/prepare_ocnExtData/extract_day_SST_FRACI_eight.F90 b/pre/prepare_ocnExtData/extract_day_SST_FRACI_eight.F90
index 0560fa6c..105cc820 100644
--- a/pre/prepare_ocnExtData/extract_day_SST_FRACI_eight.F90
+++ b/pre/prepare_ocnExtData/extract_day_SST_FRACI_eight.F90
@@ -14,6 +14,10 @@ PROGRAM extract_day_SST_FRACI_eight
   integer               :: extract_date
   integer               :: date_out, nlon_out, nlat_out
 
+  character (len = 200) :: data_path
+  character (len = 200) :: sst_file_pref, ice_file_pref
+  character (len = 200) :: sst_file_suff, ice_file_suff
+
   character (len=8)     :: date_str, format_str = "(I8)"
   character (len = 200) :: sst_file, ice_file
 
@@ -41,8 +45,15 @@ PROGRAM extract_day_SST_FRACI_eight
   extract_date = year_s*10000+month_s*100+day_s
   write (date_str, format_str) extract_date
 
-  sst_file = 'dataoceanfile_OSTIA_REYNOLDS_SST.2880x1440.'//trim(date_str(1:4))//'.data'
-  ice_file = 'dataoceanfile_OSTIA_REYNOLDS_ICE.2880x1440.'//trim(date_str(1:4))//'.data'
+  ! GMAO OPS has following fixed
+  data_path= '/discover/nobackup/projects/gmao/share/dao_ops/fvInput/g5gcm/bcs/realtime/OSTIA_REYNOLDS/2880x1440/'
+  sst_file_pref = 'dataoceanfile_OSTIA_REYNOLDS_SST.2880x1440.'
+  ice_file_pref = 'dataoceanfile_OSTIA_REYNOLDS_ICE.2880x1440.'
+  sst_file_suff = '.data'
+  ice_file_suff = '.data'
+
+  sst_file = trim(data_path)//"/"//trim(sst_file_pref)//trim(date_str(1:4))//trim(sst_file_suff)
+  ice_file = trim(data_path)//"/"//trim(ice_file_pref)//trim(date_str(1:4))//trim(ice_file_suff)
 
   print *, " "
   print *, "Read: ", sst_file, ice_file
@@ -58,10 +69,10 @@ PROGRAM extract_day_SST_FRACI_eight
   print *, "Write out the data..."
 
   ! binary write
-  ! CALL write_bin( 2023, 1, 1, 2023, 1, 2, '20230101',  1440, 2880, transpose(sst), transpose(ice))
+  ! CALL write_bin( 'daily_', 2023, 1, 1, 2023, 1, 2, '20230101',  1440, 2880, transpose(sst), transpose(ice))
 
   ! netcdf
-  CALL write_netcdf( year_s, month_s, day_s, date_str, nlat, nlon, transpose(sst), transpose(ice))
+  CALL write_netcdf( 'sst_ice', year_s, month_s, day_s, date_str, nlat, nlon, transpose(sst), transpose(ice))
 
   print *, "Done."
   print *, " "

From 860920ceaa4d4f254de63a51de989eb607881ee8 Mon Sep 17 00:00:00 2001
From: sanAkel <santha.akella@nasa.gov>
Date: Wed, 3 Jan 2024 23:17:22 -0500
Subject: [PATCH 10/20] Improve readability

---
 pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90 | 1 +
 1 file changed, 1 insertion(+)

diff --git a/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90 b/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90
index 9d467a6e..22ca2762 100644
--- a/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90
+++ b/pre/prepare_ocnExtData/daily_clim_SST_FRACI_eight.F90
@@ -152,6 +152,7 @@ PROGRAM daily_clim_SST_FRACI_eight
   clim_tom_mon  = int((clim_tom-clim_tom_year*10000)/100)
   clim_tom_day = clim_tom - (clim_tom_year*10000+clim_tom_mon*100)
   if(verbose) print *, "Tomorrow day:", clim_tom_year, clim_tom_mon, clim_tom_day
+
   write (clim_date_str,'(I0.8)') clim_year*10000+month*100+day
   if(verbose) print *, clim_date_str
 

From 90a3430f4d5d93119db0e1075b5f7173259d82ee Mon Sep 17 00:00:00 2001
From: sanAkel <santha.akella@nasa.gov>
Date: Wed, 3 Jan 2024 23:17:44 -0500
Subject: [PATCH 11/20] Add option to save binary file

---
 .../extract_day_SST_FRACI_eight.F90           | 26 ++++++++++++++-----
 1 file changed, 20 insertions(+), 6 deletions(-)

diff --git a/pre/prepare_ocnExtData/extract_day_SST_FRACI_eight.F90 b/pre/prepare_ocnExtData/extract_day_SST_FRACI_eight.F90
index 105cc820..67359777 100644
--- a/pre/prepare_ocnExtData/extract_day_SST_FRACI_eight.F90
+++ b/pre/prepare_ocnExtData/extract_day_SST_FRACI_eight.F90
@@ -1,7 +1,8 @@
 PROGRAM extract_day_SST_FRACI_eight
 !---------------------------------------------------------------------------
+  USE m_tick, only: INCYMD
   USE sst_ice_helpers, only: read_bin_SST_ICE, &
-!                            write_bin,        &
+                             write_bin,        &
                              write_netcdf
 
   IMPLICIT NONE
@@ -24,12 +25,15 @@ PROGRAM extract_day_SST_FRACI_eight
   character (len=4) :: arg
   integer :: numArgs
   integer :: year_s, month_s, day_s
+  integer :: tom, tom_year, tom_mon, tom_day
+
+  logical :: save_bin
 !----
 
   numArgs = iargc()!nargs() for ifort
-  if (numArgs < 3)  then
-    print *, "Need 3 inputs, in following format. Try again!"
-    print *, "year month day"
+  if (numArgs < 4)  then
+    print *, "Need minmum of 4 inputs, in following format (example). Try again!"
+    print *, "year month day .false."
     STOP
   end if
 
@@ -40,6 +44,8 @@ PROGRAM extract_day_SST_FRACI_eight
   read(arg, '(I14)') month_s
   call getarg(3, arg)
   read(arg, '(I14)') day_s
+  call getarg(4, arg)
+  read(arg, *) save_bin
   
   !-- form a string of input date
   extract_date = year_s*10000+month_s*100+day_s
@@ -69,8 +75,16 @@ PROGRAM extract_day_SST_FRACI_eight
   print *, "Write out the data..."
 
   ! binary write
-  ! CALL write_bin( 'daily_', 2023, 1, 1, 2023, 1, 2, '20230101',  1440, 2880, transpose(sst), transpose(ice))
-
+  if (save_bin) then
+    tom = INCYMD( year_s*10000+month_s*100+day_s, 1) ! next day
+
+    tom_year = int(tom/10000)
+    tom_mon  = int((tom-tom_year*10000)/100)
+    tom_day  = tom - (tom_year*10000+tom_mon*100)
+    !print *, "Tomorrow day:", tom_year, tom_mon, tom_day
+    CALL write_bin( 'sst_ice', year_s, month_s, day_s, tom_year, tom_mon, tom_day, date_str,  nlat, nlon, transpose(sst), transpose(ice))
+  end if
+  
   ! netcdf
   CALL write_netcdf( 'sst_ice', year_s, month_s, day_s, date_str, nlat, nlon, transpose(sst), transpose(ice))
 

From e4ac9c5d5d8457e79fb3b5627b9c1223060ca592 Mon Sep 17 00:00:00 2001
From: sanAkel <santha.akella@nasa.gov>
Date: Sun, 7 Jan 2024 16:24:04 -0500
Subject: [PATCH 12/20] Make defaults argument(s), as asked by @rtodling for
 input BCS data path. Followed his example code:
 https://github.com/GEOS-ESM/GMAO_Shared/blob/main/GMAO_hermes/dyn_recenter.f90

---
 .../extract_day_SST_FRACI_eight.F90           | 67 ++++++++++++++-----
 1 file changed, 51 insertions(+), 16 deletions(-)

diff --git a/pre/prepare_ocnExtData/extract_day_SST_FRACI_eight.F90 b/pre/prepare_ocnExtData/extract_day_SST_FRACI_eight.F90
index 67359777..93e90ca5 100644
--- a/pre/prepare_ocnExtData/extract_day_SST_FRACI_eight.F90
+++ b/pre/prepare_ocnExtData/extract_day_SST_FRACI_eight.F90
@@ -22,7 +22,8 @@ PROGRAM extract_day_SST_FRACI_eight
   character (len=8)     :: date_str, format_str = "(I8)"
   character (len = 200) :: sst_file, ice_file
 
-  character (len=4) :: arg
+  integer :: i, iarg, argc, iargc
+  character (len=255) :: argv
   integer :: numArgs
   integer :: year_s, month_s, day_s
   integer :: tom, tom_year, tom_mon, tom_day
@@ -30,29 +31,63 @@ PROGRAM extract_day_SST_FRACI_eight
   logical :: save_bin
 !----
 
-  numArgs = iargc()!nargs() for ifort
-  if (numArgs < 4)  then
-    print *, "Need minmum of 4 inputs, in following format (example). Try again!"
-    print *, "year month day .false."
+  argc = iargc()
+  if (argc < 3)  then
+    print *, "Need minmum 3 inputs, see below examples. Try again!"
+    print *, "extract_daily_sst_ice.x -year 2010 -month 1 -day 1"
+    print *, " "
+    print *, "** 2 additional optional inputs are allowed ** "
+    print *, "extract_daily_sst_ice.x -year 2010 -month 1 -day 1 -save_bin"
+    print *, " "
+    print *, "extract_daily_sst_ice.x -year 2010 -month 1 -day 1 -input_data_path xx"
+    print *, " "
+    print *, "extract_daily_sst_ice.x -year 2010 -month 1 -day 1 -save_bin -input_data_path xx"
+    print *, " "
     STOP
   end if
 
+  ! Following 2 variables have defaults that can be overridden by user inputs.
+  save_bin = .false.
+  ! GMAO OPS data path
+  data_path= '/discover/nobackup/projects/gmao/share/dao_ops/fvInput/g5gcm/bcs/realtime/OSTIA_REYNOLDS/2880x1440/'
+
   !-- read user inputs
-  call getarg(1, arg)
-  read(arg, '(I14)') year_s ! default
-  call getarg(2, arg)
-  read(arg, '(I14)') month_s
-  call getarg(3, arg)
-  read(arg, '(I14)') day_s
-  call getarg(4, arg)
-  read(arg, *) save_bin
-  
+  iarg = 0
+  do i = 1, 99
+    iarg = iarg + 1
+    if ( iarg > argc ) exit
+    call getarg(iarg, argv)
+
+    select case (argv)
+      case ("-year")
+        iarg = iarg + 1
+        call getarg ( iarg, argv )
+        read(argv, '(I14)') year_s
+      case ("-month")
+        iarg = iarg + 1
+        call getarg ( iarg, argv )
+        read(argv, '(I14)') month_s
+      case ("-day")
+        iarg = iarg + 1
+        call getarg ( iarg, argv )
+        read(argv, '(I14)') day_s
+      case ("-save_bin")
+        save_bin = .true.
+      case ("-input_data_path")
+        iarg = iarg + 1
+        call getarg ( iarg, argv )
+        read(argv, *) data_path
+        data_path = trim(data_path)
+      case default
+    end select
+  end do 
+  print *, year_s, month_s, day_s, save_bin, data_path
+
   !-- form a string of input date
   extract_date = year_s*10000+month_s*100+day_s
   write (date_str, format_str) extract_date
 
   ! GMAO OPS has following fixed
-  data_path= '/discover/nobackup/projects/gmao/share/dao_ops/fvInput/g5gcm/bcs/realtime/OSTIA_REYNOLDS/2880x1440/'
   sst_file_pref = 'dataoceanfile_OSTIA_REYNOLDS_SST.2880x1440.'
   ice_file_pref = 'dataoceanfile_OSTIA_REYNOLDS_ICE.2880x1440.'
   sst_file_suff = '.data'
@@ -82,7 +117,7 @@ PROGRAM extract_day_SST_FRACI_eight
     tom_mon  = int((tom-tom_year*10000)/100)
     tom_day  = tom - (tom_year*10000+tom_mon*100)
     !print *, "Tomorrow day:", tom_year, tom_mon, tom_day
-    CALL write_bin( 'sst_ice', year_s, month_s, day_s, tom_year, tom_mon, tom_day, date_str,  nlat, nlon, transpose(sst), transpose(ice))
+    CALL write_bin( 'bcs_2880x1440', year_s, month_s, day_s, tom_year, tom_mon, tom_day, date_str,  nlat, nlon, transpose(sst), transpose(ice))
   end if
   
   ! netcdf

From a1c57524c4943080b927ce0f310ac4effa110083 Mon Sep 17 00:00:00 2001
From: sanAkel <santha.akella@nasa.gov>
Date: Sun, 14 Jan 2024 18:20:52 -0500
Subject: [PATCH 13/20] add utility to generate forecast (persistence), clean
 up

---
 pre/prepare_ocnExtData/CMakeLists.txt         |   6 +
 .../forecast_SST_FRACI_eight.F90              | 166 ++++++++++++++++++
 .../proc_SST_FRACI_eight.F90                  |   4 +-
 pre/prepare_ocnExtData/sst_ice_helpers.F90    |   6 +-
 4 files changed, 177 insertions(+), 5 deletions(-)
 create mode 100644 pre/prepare_ocnExtData/forecast_SST_FRACI_eight.F90

diff --git a/pre/prepare_ocnExtData/CMakeLists.txt b/pre/prepare_ocnExtData/CMakeLists.txt
index 9ef9caa8..e67d9938 100644
--- a/pre/prepare_ocnExtData/CMakeLists.txt
+++ b/pre/prepare_ocnExtData/CMakeLists.txt
@@ -32,6 +32,12 @@ ecbuild_add_executable (
   LIBS ${this})
 set_target_properties(extract_daily_sst_ice.x PROPERTIES Fortran_MODULE_DIRECTORY ${include_${this}})
 
+ecbuild_add_executable (
+  TARGET gen_forecast_bcs.x
+  SOURCES forecast_SST_FRACI_eight.F90
+  LIBS ${this})
+set_target_properties(gen_forecast_bcs.x PROPERTIES Fortran_MODULE_DIRECTORY ${include_${this}})
+
 if (USE_F2PY)
    find_package(F2PY2)
    if (F2PY2_FOUND)
diff --git a/pre/prepare_ocnExtData/forecast_SST_FRACI_eight.F90 b/pre/prepare_ocnExtData/forecast_SST_FRACI_eight.F90
new file mode 100644
index 00000000..e10e3ee8
--- /dev/null
+++ b/pre/prepare_ocnExtData/forecast_SST_FRACI_eight.F90
@@ -0,0 +1,166 @@
+PROGRAM forecast_SST_FRACI_eight
+!---------------------------------------------------------------------------
+  USE m_tick, only: INCYMD
+  USE sst_ice_helpers, only: read_bin_SST_ICE, &
+                             write_bin,        &
+                             write_netcdf
+
+  IMPLICIT NONE
+
+  integer, parameter :: nlon = 2880
+  integer, parameter :: nlat = 1440
+  real               :: lon(nlon), lat(nlat)
+  real               :: sst(nlat, nlon), ice(nlat, nlon)
+  real               :: sst0(nlat, nlon), ice0(nlat, nlon)
+
+  integer               :: extract_date
+  integer               :: date_out, nlon_out, nlat_out
+
+  character (len = 200) :: data_path
+  character (len = 200) :: sst_file_pref, ice_file_pref
+  character (len = 200) :: sst_file_suff, ice_file_suff
+
+  character (len=8)     :: format_str = "(I8)"
+  character (len=8)     :: date_str, today_str, tom_str
+  character (len = 200) :: sst_file, ice_file
+
+  integer :: i, iarg, argc, iargc
+  character (len=255) :: argv
+  integer :: numArgs
+  integer :: year_s, month_s, day_s, tom
+  integer :: tom_year, tom_mon, tom_day
+  integer :: today_year, today_mon, today_day
+
+  integer :: fcst_nDays
+
+  character (len=255), parameter :: method = "persistence"
+!----
+
+  print *, " "
+  print *, "-----------------------------------------------------------"
+  print *, "        Generating SST and FRACI for forecasts             "
+  print *, "-----------------------------------------------------------"
+  print *, " "
+
+  argc = iargc()
+  if (argc < 3)  then
+    print *, "Need minmum 3 inputs, see below examples. Try again!"
+    print *, "gen_forecast_bcs.x -year 2010 -month 1 -day 1"
+    print *, " "
+    print *, "** 2 additional optional inputs are allowed ** "
+    print *, "gen_forecast_bcs.x -year 2010 -month 1 -day 1 -fcst_nDays NN"
+    print *, " "
+    print *, "gen_forecast_bcs.x -year 2010 -month 1 -day 1 -input_data_path xx"
+    print *, " "
+    print *, "gen_forecast_bcs.x -year 2010 -month 1 -day 1 -fcst_nDays NN -input_data_path xx"
+    print *, " "
+    STOP
+  end if
+
+  ! Following 2 variables have defaults that can be overridden by user inputs.
+  fcst_nDays = 15 ! default length of each forecast: 15 days.
+  ! GMAO OPS data path
+  data_path= '/discover/nobackup/projects/gmao/share/dao_ops/fvInput/g5gcm/bcs/realtime/OSTIA_REYNOLDS/2880x1440/'
+
+  !-- read user inputs
+  iarg = 0
+  do i = 1, 99
+    iarg = iarg + 1
+    if ( iarg > argc ) exit
+    call getarg(iarg, argv)
+
+    select case (argv)
+      case ("-year")
+        iarg = iarg + 1
+        call getarg ( iarg, argv )
+        read(argv, '(I14)') year_s
+      case ("-month")
+        iarg = iarg + 1
+        call getarg ( iarg, argv )
+        read(argv, '(I14)') month_s
+      case ("-day")
+        iarg = iarg + 1
+        call getarg ( iarg, argv )
+        read(argv, '(I14)') day_s
+      case("-fcst_nDays")
+        iarg = iarg + 1
+        call getarg ( iarg, argv )
+        read(argv, '(I14)') fcst_nDays
+      case ("-input_data_path")
+        iarg = iarg + 1
+        call getarg ( iarg, argv )
+        read(argv, *) data_path
+        data_path = trim(data_path)
+      case default
+    end select
+  end do 
+  print *, "Inputs: "
+  print *, year_s, month_s, day_s, fcst_nDays, data_path
+  print *, " "
+
+  ! --- Read forecast start date data ---
+  ! --
+  !-- form a string of input date
+  extract_date = year_s*10000+month_s*100+day_s
+  write (date_str, format_str) extract_date
+
+  ! GMAO OPS has following fixed
+  sst_file_pref = 'dataoceanfile_OSTIA_REYNOLDS_SST.2880x1440.'
+  ice_file_pref = 'dataoceanfile_OSTIA_REYNOLDS_ICE.2880x1440.'
+  sst_file_suff = '.data'
+  ice_file_suff = '.data'
+
+  sst_file = trim(data_path)//"/"//trim(sst_file_pref)//trim(date_str(1:4))//trim(sst_file_suff)
+  ice_file = trim(data_path)//"/"//trim(ice_file_pref)//trim(date_str(1:4))//trim(ice_file_suff)
+
+  ! read start date sst
+  CALL read_bin_SST_ICE( sst_file, extract_date, date_out, nlon_out, nlat_out, lon, lat, sst0)
+  ! read start date fraci
+  CALL read_bin_SST_ICE( ice_file, extract_date, date_out, nlon_out, nlat_out, lon, lat, ice0)
+  ! --
+
+  print *, " "
+  print *, "For:  ", date_out
+  print *, "Read: ", sst_file, ice_file
+  print *, " "
+
+  ! day 1: [year_s,   month_s,   day_s   --> year_s+1, month_s+1, day_s+1]
+  ! day 2: [year_s+1, month_s+1, day_s+1 --> year_s+2, month_s+2, day_s+2]
+  ! ..
+  print *, "Write out: ", trim(method), " BCs data..."
+  print *, " "
+
+  today_year=year_s; today_mon=month_s; today_day=day_s
+  do i = 1, fcst_nDays
+    
+    tom = INCYMD( today_year*10000+today_mon*100+today_day, 1) ! tomorrow = today + (1 day)
+    tom_year = int(tom/10000)
+    tom_mon  = int((tom-tom_year*10000)/100)
+    tom_day  = tom - (tom_year*10000+tom_mon*100)
+
+    if (method == "persistence") then
+     sst = sst0; ice = ice0
+    end if
+
+    write (today_str, format_str) today_year*10000+today_mon*100+today_day
+    write (tom_str,   format_str) tom_year*10000  +tom_mon*100  +tom_day
+
+    !print *, today_year,"/",today_mon,"/",today_day, ">>", tom_year,"/",tom_mon,"/",tom_day
+    print *, trim(today_str), "---->", trim(tom_str)
+
+    CALL write_bin( 'bcs_2880x1440', today_year, today_mon, today_day, tom_year, tom_mon, tom_day, today_str, &
+    nlat, nlon, transpose(sst), transpose(ice))
+
+    CALL write_netcdf( 'sst_ice', today_year, today_mon, today_day, today_str, &
+    nlat, nlon, transpose(sst), transpose(ice))
+
+    ! done with today; increment it to tomorrow
+    today_year=tom_year; today_mon=tom_mon; today_day=tom_day
+  end do
+
+  print *, " "
+  print *, "Done."
+  print *, " "
+
+!---------------------------------------------------------------------------
+END PROGRAM forecast_SST_FRACI_eight
diff --git a/pre/prepare_ocnExtData/proc_SST_FRACI_eight.F90 b/pre/prepare_ocnExtData/proc_SST_FRACI_eight.F90
index 2de2d495..06f85944 100644
--- a/pre/prepare_ocnExtData/proc_SST_FRACI_eight.F90
+++ b/pre/prepare_ocnExtData/proc_SST_FRACI_eight.F90
@@ -255,14 +255,14 @@ PROGRAM proc_SST_FRACI_eight
         READ( tomrw_Day,  99) tomrw_iDay
 
         if (bin_output == 'yes') then
-          call write_bin(today_iYear, today_iMon, today_iDay, &
+          call write_bin('bcs_2880x1440', today_iYear, today_iMon, today_iDay, &
                          tomrw_iYear, tomrw_iMon, tomrw_iDay, &
                          today,                               &
                          NLAT_out,    NLON_out,               &
                          ostia_SST_eigth, ostia_ICE_eigth)
         endif
 
-        call write_netcdf(today_iYear, today_iMon, today_iDay, &
+        call write_netcdf('sst_ice', today_iYear, today_iMon, today_iDay, &
                           today,                               &
                           NLAT_out,    NLON_out,               &
                           ostia_SST_eigth, ostia_ICE_eigth)
diff --git a/pre/prepare_ocnExtData/sst_ice_helpers.F90 b/pre/prepare_ocnExtData/sst_ice_helpers.F90
index e0ebc3b7..e11bb4bf 100644
--- a/pre/prepare_ocnExtData/sst_ice_helpers.F90
+++ b/pre/prepare_ocnExtData/sst_ice_helpers.F90
@@ -1017,8 +1017,8 @@ SUBROUTINE write_bin( fileName_pref, today_year, today_mon, today_day, &
     CLOSE(991)
     CLOSE(992)
 
-    PRINT*, "Finished writing:"
-    PRINT*, fileName_sst, fileName_ice
+    !PRINT*, "Finished writing:"
+    !PRINT*, fileName_sst, fileName_ice
 !---------------------------------------------------------------------------
 
     END SUBROUTINE write_bin
@@ -1150,7 +1150,7 @@ SUBROUTINE write_netcdf( fileName_pref, today_year, today_mon, today_day, &
     call check( nf90_close(ncid) ) ! Close the file.
 !---------------------------------------------------------------------------
 
-    PRINT*, "Finished writing ... ", fileName
+    !PRINT*, "Finished writing ... ", fileName
 
     END SUBROUTINE write_netcdf
 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

From 34816943e848aa4c3766c9daf0b027652d337cd2 Mon Sep 17 00:00:00 2001
From: sanAkel <santha.akella@nasa.gov>
Date: Sun, 14 Jan 2024 18:50:00 -0500
Subject: [PATCH 14/20] remove commented lines

---
 pre/prepare_ocnExtData/forecast_SST_FRACI_eight.F90 | 2 --
 1 file changed, 2 deletions(-)

diff --git a/pre/prepare_ocnExtData/forecast_SST_FRACI_eight.F90 b/pre/prepare_ocnExtData/forecast_SST_FRACI_eight.F90
index e10e3ee8..1b0958d4 100644
--- a/pre/prepare_ocnExtData/forecast_SST_FRACI_eight.F90
+++ b/pre/prepare_ocnExtData/forecast_SST_FRACI_eight.F90
@@ -144,8 +144,6 @@ PROGRAM forecast_SST_FRACI_eight
 
     write (today_str, format_str) today_year*10000+today_mon*100+today_day
     write (tom_str,   format_str) tom_year*10000  +tom_mon*100  +tom_day
-
-    !print *, today_year,"/",today_mon,"/",today_day, ">>", tom_year,"/",tom_mon,"/",tom_day
     print *, trim(today_str), "---->", trim(tom_str)
 
     CALL write_bin( 'bcs_2880x1440', today_year, today_mon, today_day, tom_year, tom_mon, tom_day, today_str, &

From 69ccffcd51563d819d9f69e335a3ad61f296c945 Mon Sep 17 00:00:00 2001
From: sanAkel <santha.akella@nasa.gov>
Date: Wed, 17 Jan 2024 18:50:11 -0500
Subject: [PATCH 15/20] Add option to generate BCs using persistence of initial
 anomaly

---
 .../forecast_SST_FRACI_eight.F90              | 85 +++++++++++++++++--
 1 file changed, 76 insertions(+), 9 deletions(-)

diff --git a/pre/prepare_ocnExtData/forecast_SST_FRACI_eight.F90 b/pre/prepare_ocnExtData/forecast_SST_FRACI_eight.F90
index 1b0958d4..b754f4fa 100644
--- a/pre/prepare_ocnExtData/forecast_SST_FRACI_eight.F90
+++ b/pre/prepare_ocnExtData/forecast_SST_FRACI_eight.F90
@@ -12,6 +12,8 @@ PROGRAM forecast_SST_FRACI_eight
   real               :: lon(nlon), lat(nlat)
   real               :: sst(nlat, nlon), ice(nlat, nlon)
   real               :: sst0(nlat, nlon), ice0(nlat, nlon)
+  real               :: sst_daily_clim(nlat, nlon), ice_daily_clim(nlat, nlon)
+  real               :: anom_sst0(nlat, nlon), anom_ice0(nlat, nlon)
 
   integer               :: extract_date
   integer               :: date_out, nlon_out, nlat_out
@@ -26,14 +28,13 @@ PROGRAM forecast_SST_FRACI_eight
 
   integer :: i, iarg, argc, iargc
   character (len=255) :: argv
+  character (len=255) :: method
   integer :: numArgs
   integer :: year_s, month_s, day_s, tom
   integer :: tom_year, tom_mon, tom_day
   integer :: today_year, today_mon, today_day
 
   integer :: fcst_nDays
-
-  character (len=255), parameter :: method = "persistence"
 !----
 
   print *, " "
@@ -47,20 +48,23 @@ PROGRAM forecast_SST_FRACI_eight
     print *, "Need minmum 3 inputs, see below examples. Try again!"
     print *, "gen_forecast_bcs.x -year 2010 -month 1 -day 1"
     print *, " "
-    print *, "** 2 additional optional inputs are allowed ** "
+    print *, "** 3 additional optional inputs are allowed ** "
     print *, "gen_forecast_bcs.x -year 2010 -month 1 -day 1 -fcst_nDays NN"
     print *, " "
     print *, "gen_forecast_bcs.x -year 2010 -month 1 -day 1 -input_data_path xx"
     print *, " "
     print *, "gen_forecast_bcs.x -year 2010 -month 1 -day 1 -fcst_nDays NN -input_data_path xx"
     print *, " "
+    print *, "gen_forecast_bcs.x -year 2010 -month 1 -day 1 -fcst_nDays NN -input_data_path xx -method yy"
+    print *, " "
     STOP
   end if
 
-  ! Following 2 variables have defaults that can be overridden by user inputs.
+  ! Following variables have defaults that can be overridden by user inputs.
   fcst_nDays = 15 ! default length of each forecast: 15 days.
   ! GMAO OPS data path
   data_path= '/discover/nobackup/projects/gmao/share/dao_ops/fvInput/g5gcm/bcs/realtime/OSTIA_REYNOLDS/2880x1440/'
+  method   = 'persistence'
 
   !-- read user inputs
   iarg = 0
@@ -91,11 +95,16 @@ PROGRAM forecast_SST_FRACI_eight
         call getarg ( iarg, argv )
         read(argv, *) data_path
         data_path = trim(data_path)
+      case ("-method")
+        iarg = iarg + 1
+        call getarg ( iarg, argv )
+        read(argv, *) method
+        method = trim(method)
       case default
     end select
   end do 
   print *, "Inputs: "
-  print *, year_s, month_s, day_s, fcst_nDays, data_path
+  print *, year_s, month_s, day_s, fcst_nDays, data_path, method
   print *, " "
 
   ! --- Read forecast start date data ---
@@ -138,14 +147,27 @@ PROGRAM forecast_SST_FRACI_eight
     tom_mon  = int((tom-tom_year*10000)/100)
     tom_day  = tom - (tom_year*10000+tom_mon*100)
 
-    if (method == "persistence") then
-     sst = sst0; ice = ice0
-    end if
-
     write (today_str, format_str) today_year*10000+today_mon*100+today_day
     write (tom_str,   format_str) tom_year*10000  +tom_mon*100  +tom_day
     print *, trim(today_str), "---->", trim(tom_str)
 
+    if (method == "persistence") then
+      sst = sst0; ice = ice0
+    else if (method == "persist_anomaly") then
+      call get_daily_clim( today_str, nlat_out, nlon_out, sst_daily_clim, ice_daily_clim)
+      if (i > 1) then
+        sst = sst_daily_clim + anom_sst0
+        ice = ice_daily_clim + anom_ice0
+      else
+        sst = sst0; ice = ice0 ! forecast day-1: just write out the same field(s).
+       anom_sst0 = sst0 - sst_daily_clim
+       anom_ice0 = ice0 - ice_daily_clim
+      end if
+    else
+      print *, "Unknown method: ", method, "Exiting. Fix and try again."
+      STOP
+    end if
+
     CALL write_bin( 'bcs_2880x1440', today_year, today_mon, today_day, tom_year, tom_mon, tom_day, today_str, &
     nlat, nlon, transpose(sst), transpose(ice))
 
@@ -160,5 +182,50 @@ PROGRAM forecast_SST_FRACI_eight
   print *, "Done."
   print *, " "
 
+!---------------------------------------------------------------------------
+  CONTAINS
+  SUBROUTINE get_daily_clim(today_str, nlat_out, nlon_out, sst_daily_clim, ice_daily_clim)
+
+    integer :: date_out, nlon_out, nlat_out
+    real    :: sst_daily_clim(nlat, nlon), ice_daily_clim(nlat, nlon)
+    character (len = 200) :: clim_data_path, sst_clim_pref, ice_clim_pref
+    character (len=8)     :: today_str
+
+    clim_data_path = "/discover/nobackup/projects/gmao/advda/sakella/future_sst_fraci/data/binFiles/"
+    sst_clim_pref  = "daily_clim_mean_sst_0001" !clim_year  = 1 ! from daily_clim_SST_FRACI_eight.F90
+    ice_clim_pref  = "daily_clim_mean_ice_0001"
+
+    sst_file = trim(clim_data_path)//"/"//trim(sst_clim_pref)//trim(today_str(5:8))//".bin"
+    ice_file = trim(clim_data_path)//"/"//trim(ice_clim_pref)//trim(today_str(5:8))//".bin"
+    print *, "Reading daily clim from: "
+    print *, sst_file
+    print *, ice_file
+
+    CALL read_daily_clim(sst_file, sst_daily_clim, nlat, nlon) ! read sst anomaly
+    CALL read_daily_clim(ice_file, ice_daily_clim, nlat, nlon) ! read fraci anomaly
+    print *, nlat, nlon
+  END SUBROUTINE get_daily_clim
+
+  SUBROUTINE read_daily_clim(fName, bcs_field, nlat, nlon)
+    integer, intent(in) :: nlat, nlon
+    character (len=*), intent(in) :: fName
+    real, intent(out) :: bcs_field(nlat, nlon)
+
+    real    year1,month1,day1,hour1,min1,sec1, year2,month2,day2,hour2,min2,sec2
+    real    dum1,dum2, field(nlon,nlat)
+    integer rc
+
+    open (10,file=fName,form='unformatted',access='sequential', STATUS = 'old') ! these files only 1 record (header+data)
+    rc = 0
+    do while (rc.eq.0)
+      ! read header
+      read (10,iostat=rc) year1,month1,day1,hour1,min1,sec1,&
+                          year2,month2,day2,hour2,min2,sec2,dum1,dum2
+      ! read data
+      if( rc.eq.0 ) read (10,iostat=rc) field
+      bcs_field = TRANSPOSE(field)
+      close(10)
+    end do
+  END SUBROUTINE read_daily_clim
 !---------------------------------------------------------------------------
 END PROGRAM forecast_SST_FRACI_eight

From a80a172f41e396ff7bd4cfd920e834db3b09d6e1 Mon Sep 17 00:00:00 2001
From: sanAkel <santha.akella@nasa.gov>
Date: Wed, 17 Jan 2024 18:53:42 -0500
Subject: [PATCH 16/20] update change log

---
 CHANGELOG.md | 4 ++++
 1 file changed, 4 insertions(+)

diff --git a/CHANGELOG.md b/CHANGELOG.md
index 1dc03a07..141ff15a 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -19,6 +19,10 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 ## [2.0.2] - 2023-06-28
 
+### Added
+
+- Addition of SST and FRACI forecast (ocean) boundary conditions generation capability in `pre/prepare_ocnExtData`
+
 ### Fixed
 
 - Updates to fix dates for comparison plots. Updates for plot levels.

From 3210ede3599682cd97dca9bd083044fe6933d605 Mon Sep 17 00:00:00 2001
From: Santha Akella <18061653+sanAkel@users.noreply.github.com>
Date: Thu, 18 Jan 2024 19:28:58 -0500
Subject: [PATCH 17/20] Update CHANGELOG.md

addressing https://github.com/GEOS-ESM/GEOS_Util/pull/47#discussion_r1457443445 from @mathomp4
---
 CHANGELOG.md | 5 +----
 1 file changed, 1 insertion(+), 4 deletions(-)

diff --git a/CHANGELOG.md b/CHANGELOG.md
index e98e2ecd..ef35f96c 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -6,7 +6,7 @@ The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/),
 and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
 
 ## [Unreleased]
-
+- Addition of SST and FRACI forecast (ocean) boundary conditions generation capability in `pre/prepare_ocnExtData`
 ### Added
 
 ### Changed
@@ -52,9 +52,6 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 ## [2.0.2] - 2023-06-28
 
-### Added
-
-- Addition of SST and FRACI forecast (ocean) boundary conditions generation capability in `pre/prepare_ocnExtData`
 
 ### Fixed
 

From 176f65bdb3719b6c01db413dbbf5d37d5829321f Mon Sep 17 00:00:00 2001
From: Santha Akella <18061653+sanAkel@users.noreply.github.com>
Date: Thu, 18 Jan 2024 19:31:24 -0500
Subject: [PATCH 18/20] Update gen_wts_regrid_ostia_geos.py

change python3 directive to python to accommodate GEOSadas which apparently is not caught up with python3 (per @rtodling)
---
 pre/prepare_ocnExtData/gen_wts_regrid_ostia_geos.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/pre/prepare_ocnExtData/gen_wts_regrid_ostia_geos.py b/pre/prepare_ocnExtData/gen_wts_regrid_ostia_geos.py
index c3456c3e..4e4d79e9 100755
--- a/pre/prepare_ocnExtData/gen_wts_regrid_ostia_geos.py
+++ b/pre/prepare_ocnExtData/gen_wts_regrid_ostia_geos.py
@@ -1,4 +1,4 @@
-#!/usr/bin/env python3
+#!/usr/bin/env python
 
 import xarray as xr
 import xesmf as xe

From 2f44d7720547637a5fa28ff36ae2efe0826fde47 Mon Sep 17 00:00:00 2001
From: Santha Akella <18061653+sanAkel@users.noreply.github.com>
Date: Thu, 18 Jan 2024 19:32:18 -0500
Subject: [PATCH 19/20] Update plot_diff_ostia_native_regridded_binned.py

change python3 directive to python to accommodate GEOSadas which apparently is not caught up with python3 (per @rtodling)
---
 .../plot_diff_ostia_native_regridded_binned.py                  | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/pre/prepare_ocnExtData/plot_diff_ostia_native_regridded_binned.py b/pre/prepare_ocnExtData/plot_diff_ostia_native_regridded_binned.py
index 4d569cf8..568f9b4c 100755
--- a/pre/prepare_ocnExtData/plot_diff_ostia_native_regridded_binned.py
+++ b/pre/prepare_ocnExtData/plot_diff_ostia_native_regridded_binned.py
@@ -1,4 +1,4 @@
-#!/usr/bin/env python3
+#!/usr/bin/env python
 
 import xarray as xr
 import xesmf as xe

From 7f64f4fe4cfcd2920cee0d01500c51df168b1a8e Mon Sep 17 00:00:00 2001
From: Matthew Thompson <matthew.thompson@nasa.gov>
Date: Fri, 19 Jan 2024 08:13:58 -0500
Subject: [PATCH 20/20] Update CHANGELOG.md

---
 CHANGELOG.md | 4 +++-
 1 file changed, 3 insertions(+), 1 deletion(-)

diff --git a/CHANGELOG.md b/CHANGELOG.md
index ef35f96c..07da8d5b 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -6,9 +6,11 @@ The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/),
 and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
 
 ## [Unreleased]
-- Addition of SST and FRACI forecast (ocean) boundary conditions generation capability in `pre/prepare_ocnExtData`
+
 ### Added
 
+- Addition of SST and FRACI forecast (ocean) boundary conditions generation capability in `pre/prepare_ocnExtData`
+
 ### Changed
 
 ### Fixed