FV3 CDEPS Inline for GL FVCOM Data

README.md

Steps

Get Code and Build

The work is still in a development fork.

• git clone --recursive -b feature/inline https://github.com/uturuncoglu/ufs-weather-model.git
• cd ufs-weather-model/tests
• ./compile.sh "orion" "-DAPP=ATM -DCCPP_SUITES=FV3_GFS_v15p2,FV3_GFS_v16,FV3_GFS_v17_p8,FV3_RRFS_v1beta,FV3_HRRR,FV3_RAP,FV3_GFS_v15_thompson_mynn_lam3km,FV3_WoFS_v0 -D32BIT=ON -DCDEPS_INLINE=ON" rrfs intel NO NO

Note that this is RRFS configuration (25 km CONUS) and built on "orion" (replace this with "derecho" if you need). Once the code is built, the executable fv3_rrfs.exe can be found in the same directory.

Run Configuration

The pre-configured run directories can be found on Orion and Derecho,

• Orion: /work/noaa/nems/tufuk/SRW/expt_dirs/test_community/2023070600_dev_inline
• Derecho: /glade/work/turuncu/GL/2023070600_dev_inline

The FVCOM data is staged on INPUT_DATA/ folder and stream.config includes the namelist options for the CDEPS Inline.

The input.nml (main configuration file for FV3) has extra arguments called as use_inline. This needs to be set to .true. to use CDEPS Inline capability under FV3 along with the -DCDEPS_INLINE=ON build option in the compile step.

To perform out-of-box run (without using CDEPS inline), just set use_inline = .false. in input.nml and compile without -DCDEPS_INLINE=ON option.

The user could copy the run directory and issue following commands to use executable that is created in the previous step. These are example commands on NCAR's Derecho.

cp -r /glade/work/turuncu/GL/2023070600_dev_inline .
cd 2023070600_dev_inline
ln -sf ../ufs-weather-model/tests/fv3_rrfs.exe fv3.exe
qsub job_card (or sbatch job_card on Orion)

Note: Do not forget to change #PBS -A P93300606 line in the job_card (the actual line could be different for Orion)

Prepare Data

This part is optional if you want to run the model for different dates.

• Download data using ./get_data.sh script. Note that https://www.glerl.noaa.gov/emf/OWAQ_fv3 has only recent data. I think the each file found in the link is belong to different 2-days run (needs to be confirmed).
• Subset data over Great Lakes region using subset.sh script.
• Fix input data using fix.sh script. This script fixes the time axis and applies unit conversion. It also creates a mask file which is used by CDEPS inline to update only specific part of the domain.

The following steps are only required if ESMF mesh file for data is not available or region of interest is changed

• Use gen_script.ncl NCL script to create SCRIP grid definition file that represents the subsetted data. Note that ifile argument could be different in your case.
• Use create_mesh.sh file to create ESMF Mesh file that will be used by CDEPS inline to define the domain. Note that this script uses ESMF module from UFS Weather Model dependencies and path of the model needs to be changed.

create_mesh.sh

#!/bin/bash

module use /work/noaa/nems/tufuk/SRW/ufs-weather-model/modulefiles
module load ufs_orion.intel

ESMF_Scrip2Unstruct scrip.nc ESMFmesh.nc 0 ESMF

fix.sh

#!/bin/bash
# This script fixes data files (fixing time axis, unit conversion etc.)
# and create mask file which is used by CDEPS inline to update only
# cerating part of data in CCPP/Physics

module load cdo
module load nco

# loop over files and fix them
lst=`ls -al tsfc_fv3grid_*.nc | grep -v fixed | awk '{print $9}'`
for i in $lst
do
  ofile=${i/.nc/_fixed.nc}
  taxis_start=`cdo -s -w sinfo $i | grep ":[0-9][0-9]:" | grep -v "RefTime" | head -n 1 | awk '{print $1","$2}'`
  echo "$i --> $ofile [${taxis_start}]"

  cdo -L -chname,time,Time \
         -setreftime,1970-01-01,00:00:00,1hour \
         -settaxis,${taxis_start},1hour \
         -setcalendar,standard \
         -aexpr,'tisfc=tisfc+273.15;tsfc=tsfc+273.15;twsfc=twsfc+273.15;' \
         -setattribute,tisfc@units="K",tsfc@units="K",twsfc@units="K" \
         $i $ofile
done

# create mask file
ncap2 -O -s "glmask_new=glmask+array(0.0,0.0,tsfc);" tsfc_fv3grid_202318612_sub_fixed.nc tmp.nc
ncks -O -v glmask_new -d time,0 tmp.nc mask.nc
rm tmp.nc
ncrename -O -v glmask_new,glmask mask.nc

gen_scrip.ncl

begin
  ;--- arguments ---
  ifile = "tsfc_fv3grid_202318612_sub_fixed.nc"
  ofile = "scrip.nc"
  lat_name = ""
  lon_name = ""
  msk_name = "glmask"

  ;--- read file ---
  grd_file = addfile(ifile, "r")

  ;--- list of variables ---
  lst_var = getfilevarnames(grd_file)

  ;--- check variables ---
  do j = 0, dimsizes(lst_var)-1
    var_type = typeof(grd_file->$lst_var(j)$)
    if (isStrSubset(str_lower(lst_var(j)), "lon") .or. \
        isStrSubset(str_lower(lst_var(j)), "xc")) then
      if (var_type .ne. "integer") then
        lon_name = lst_var(j)
      end if
    end if 
    if (isStrSubset(str_lower(lst_var(j)), "lat") .or. \
        isStrSubset(str_lower(lst_var(j)), "yc")) then
      if (var_type .ne. "integer") then
        lat_name = lst_var(j)
      end if
    end if
    if (isStrSubset(str_lower(lst_var(j)), "mask")) then
      msk_name = lst_var(j)
    end if
  end do

  ;--- read variables ---
  if (.not. str_is_blank(lat_name)) then
    lat = grd_file->$lat_name$
    rank = dimsizes(dimsizes(lat))
  end if
  if (.not. str_is_blank(lon_name)) then
    lon = grd_file->$lon_name$
  end if
  if (.not. str_is_blank(msk_name)) then
    msk = grd_file->$msk_name$
  end if
  print(lat_name+" "+lon_name+" "+msk_name+" "+rank)

  ;--- generate SCRIP file ---
  opt = True
  opt@Debug = False
  opt@Testit = False
  opt@ForceOverwrite = True
  opt@PrintTimings = True
  opt@NetCDFType = "netcdf4"
  opt@Title = "input file: "+ifile
  if (.not. str_is_blank(msk_name)) then
    if (min(msk) .ne. max(msk)) then
      opt@GridMask = msk
      print("mask min = "+sprintf("%8.3f", min(msk))+" max = "+sprintf("%8.3f", max(msk)))
    end if
  end if

  ;--- fix for single point data ---
  dims = dimsizes(lat)
  if (dims(0) .eq. 1 .and. dims(1) .eq. 1) then
    rectilinear_to_SCRIP(ofile, lat(0,0), lon(0,0), opt)
  else
    if (rank .eq. 1) then
      rectilinear_to_SCRIP(ofile, lat, lon, opt)
    end if
    if (rank .eq. 2) then
      processed = False
      if (dims(0) .eq. 1) then
        rectilinear_to_SCRIP(ofile, lat(0,:), lon(0,:), opt)
        processed = True
      end if
      if (dims(1) .eq. 1) then
        rectilinear_to_SCRIP(ofile, lat(:,0), lon(:,0), opt)
        processed = True
      end if
      if (.not. processed) then
        curvilinear_to_SCRIP(ofile, lat, lon, opt)
      end if
    end if
  end if

  ;--- add area to SCRIP file ---
  scripFile = addfile(ofile, "w")

  grid_size = dimsizes(scripFile->grid_center_lat)

  do j = 0, grid_size-1
    temp_tlat = (/ scripFile->grid_corner_lat(j,2), \
                   scripFile->grid_corner_lat(j,1), \
                   scripFile->grid_corner_lat(j,0), \
                   scripFile->grid_corner_lat(j,3) /)
    temp_tlon = (/ scripFile->grid_corner_lon(j,2), \
                   scripFile->grid_corner_lon(j,1), \
                   scripFile->grid_corner_lon(j,0), \
                   scripFile->grid_corner_lon(j,3) /)
  end do
end

get_data.sh

#!/bin/bash
# Gets recent data from https://www.glerl.noaa.gov/emf/OWAQ_fv3

# get file list
wget https://www.glerl.noaa.gov/emf/OWAQ_fv3
lst=`cat OWAQ_fv3 | grep tsfc_fv3grid | awk -F\" '{print $8}' | awk '{print $1}'`
rm OWAQ_fv3

# download them
for i in $lst
do
  wget -c https://www.glerl.noaa.gov/emf/OWAQ_fv3/$i
done

plot_2d.ncl

load "$NCARG_ROOT/lib/ncarg/nclscripts/contrib/cd_string.ncl"

begin
  ;--- arguments ---
  ;vname = "tmp2m"
  vname = "tmpsfc"
  tindx = 11

  ;--- open files ---
  tstr = sprinti("%0.3i", tindx)
  print(""+tstr)
  nc1 = addfile("../2023070600/phyf"+tstr+".nc", "r")
  nc2 = addfile("phyf"+tstr+".nc", "r")

  ;--- read data ---
  var1 = nc1->$vname$(0,73:120,123:183)
  var2 = nc2->$vname$(0,73:120,123:183)

  land = nc1->land(0,73:120,123:183)
  var1 = mask(var1, land .gt. 0, False)
  var2 = mask(var2, land .gt. 0, False)

  lat2d = nc1->lat(73:120,123:183)
  lon2d = nc1->lon(73:120,123:183)

  var1@lon2d = lon2d
  var1@lat2d = lat2d
  var2@lon2d = lon2d
  var2@lat2d = lat2d

  ;--- plot data ---
  wks_type = "newpdf"
  wks = gsn_open_wks(wks_type, "plot_"+vname+"_2d_"+tstr)

  res = True
  res@gsnFrame = False
  res@gsnLeftString = ""
  res@gsnRightString = ""
  res@pmTickMarkDisplayMode = "Always"
  res@tiXAxisFontHeightF = 0.01
  res@tiYAxisFontHeightF = 0.01
  res@tmXTOn = False
  res@tmXTLabelsOn = False
  res@tmXBOn = False
  res@tmXBLabelsOn = False
  res@tmYROn = False
  res@tmYRLabelsOn = False
  res@tmYLOn = False
  res@tmYLLabelsOn = False
  res@tmXBLabelFontHeightF = 0.01
  res@tmYLLabelFontHeightF = 0.01
  res@gsnSpreadColors = True
  res@gsnSpreadColorEnd = -3
  res@cnFillOn = True
  res@cnLinesOn = False
  res@cnLineLabelsOn = False
  res@cnInfoLabelOn = False
  res@cnFillMode = "RasterFill"
  res@lbLabelBarOn = False

  yp = 0.3
  res@vpHeightF = yp

  minlat = 40.0
  maxlat = 49.0
  minlon = 269.0
  maxlon = 286.0

  res@mpProjection = "LambertConformal"
  res@mpOutlineDrawOrder = "PostDraw"
  res@mpGridAndLimbOn = False
  onesix = (maxlon-minlon)/6.0
  std1 = minlon+onesix
  std2 = maxlon-onesix
  print("[warning] -- one-six-rule standard parallels ("+\
  sprintf("%6.2f", std1)+","+sprintf("%6.2f", std2)+") are used.")
  res@mpLambertParallel1F = std1
  res@mpLambertParallel2F = std2

  res@mpDataBaseVersion = "HighRes"
  res@mpLimitMode = "Corners"
  res@mpLeftCornerLatF = minlat
  res@mpLeftCornerLonF = minlon
  res@mpRightCornerLatF = maxlat
  res@mpRightCornerLonF = maxlon
  clon = (maxlon+minlon)/2.0
  res@mpLambertMeridianF = clon
  print("[warning] -- default center longtitude ("+sprintf("%6.2f", clon)+") is used. ")
  res@gsnAddCyclic = False

  res@cnLevelSelectionMode = "ExplicitLevels"
  res@cnFillPalette = "NCV_banded"
  if (any(vname .eq. (/ "tmp2m", "tmpsfc" /))) then
    var1 = var1-273.15 ; K -> degC
    var2 = var2-273.15 ; K -> degC
    res@cnLevels = fspan(10, 30, 21)
  end if  

  ;--- plot REF ---
  res@vpXF = 0.1
  res@vpYF = 0.98
  plot1 = gsn_csm_contour_map(wks, var1, res)

  ;--- plot SST ---
  res@lbLabelBarOn = True
  res@pmLabelBarHeightF = 0.04

  res@vpXF = 0.1
  res@vpYF = res@vpYF-yp
  plot2 = gsn_csm_contour_map(wks, var2, res)

  ;--- plot DIFF ---
  delete(res@cnLevels)
  res@cnFillPalette = "NCV_blu_red"
  if (any(vname .eq. (/ "tmp2m", "tmpsfc" /))) then
    res@cnLevels = fspan(-5, 5, 21)
  end if

  res@vpYF = res@vpYF-yp-0.04
  res@vpXF = 0.1
  diff = var1
  diff = var2-var1
  plot3 = gsn_csm_contour_map(wks, diff, res)

  ;--- add sub-plot title ---
  txres = True
  txres@gsnFrame = False
  txres@txFontHeightF = 0.02
  txres@txJust = "CenterRight"
  txres@txFont = "Helvetica-bold"

  txres@txPosXF = 0.6
  txres@txPosYF = 0.96-yp+0.04
  txid = gsn_create_text(wks, "REF", txres)
  draw(txid)

  txres@txPosXF = 0.6
  txres@txPosYF = txres@txPosYF-yp
  txid = gsn_create_text(wks, "SST", txres)
  draw(txid)

  txres@txPosXF = 0.6
  txres@txPosYF = txres@txPosYF-yp-0.04
  txid = gsn_create_text(wks, "SST-REF", txres)
  draw(txid)

  ;--- add plot title ---
  txres@txJust = "CenterLeft"
  txres@txFontHeightF = 0.01
  txres@txPosXF = 0.2
  txres@txPosYF = 0.99 
  tstr = cd_string(nc1->time, vname+" @ %Y-%N-%D %H:%M:%S")
  txid = gsn_create_text(wks, tstr, txres)
  draw(txid)
 
  frame(wks)
end

stream.config

stream_info:             stream01 stream02

taxmode01:               limit
mapalgo01:               bilinear 
tInterpAlgo01:           linear
readMode01:              single
dtlimit01:               1.5
stream_offset01:         0
yearFirst01:             2023
yearLast01:              2023
yearAlign01:             2023
stream_vectors01:        null
stream_mesh_file01:      INPUT_DATA/ESMFmesh.nc
stream_lev_dimname01:    null
stream_data_files01:     "INPUT_DATA/tsfc_fv3grid_202318700_sub_fixed.nc"
stream_data_variables01: "twsfc So_t" "aice Si_ifrac" "vice Si_vice"
stream_dst_mask01:       1

taxmode02:               cycle 
mapalgo02:               consd 
tInterpAlgo02:           lower
readMode02:              single
dtlimit02:               1.5
stream_offset02:         0
yearFirst02:             2023 
yearLast02:              2023
yearAlign02:             2023
stream_vectors02:        null
stream_mesh_file02:      INPUT_DATA/ESMFmesh.nc
stream_lev_dimname02:    null
stream_data_files02:     INPUT_DATA/mask.nc
stream_data_variables02: "glmask So_omask"
stream_dst_mask02:       1

subset.sh

#!/bin/bash
# Subsets data over great lakes region once it is downloaded
module load nco
lst=`ls -al tsfc_fv3grid_*.nc | grep -v sub | awk '{print $9}'`
for i in $lst
do
  ncks -O -d lat,940,1250 -d lon,2490,2960 $i ${i/.nc/_sub.nc}
done