From e204fb826412288f8f789b2df3a71c9c78af2099 Mon Sep 17 00:00:00 2001 From: Tony Craig Date: Thu, 25 Mar 2021 15:35:54 -0700 Subject: [PATCH] Add JRA55 dataset tool and documentation (#582) * Add JRA55 dataset tool Add configuration/tools/jra55_datasets directory and scripts Move convert_restarts to configuration/tools/cice4_restart_conversio Add tools section to the developer guide Update tested compiler versions Change "compliance" to "validation" in documentation when discussing QC testing. * update documentation --- .../convert_restarts.f90 | 0 .../interp_jra55_ncdf_bilinear.py | 441 ++++++++++++++++++ .../tools/jra55_datasets/make_forcing.csh | 49 ++ doc/source/developer_guide/dg_scripts.rst | 24 +- doc/source/developer_guide/dg_tools.rst | 150 ++++++ doc/source/developer_guide/index.rst | 1 + doc/source/user_guide/ug_running.rst | 24 +- doc/source/user_guide/ug_testing.rst | 18 +- 8 files changed, 683 insertions(+), 24 deletions(-) rename configuration/tools/{ => cice4_restart_conversion}/convert_restarts.f90 (100%) create mode 100755 configuration/tools/jra55_datasets/interp_jra55_ncdf_bilinear.py create mode 100755 configuration/tools/jra55_datasets/make_forcing.csh create mode 100644 doc/source/developer_guide/dg_tools.rst diff --git a/configuration/tools/convert_restarts.f90 b/configuration/tools/cice4_restart_conversion/convert_restarts.f90 similarity index 100% rename from configuration/tools/convert_restarts.f90 rename to configuration/tools/cice4_restart_conversion/convert_restarts.f90 diff --git a/configuration/tools/jra55_datasets/interp_jra55_ncdf_bilinear.py b/configuration/tools/jra55_datasets/interp_jra55_ncdf_bilinear.py new file mode 100755 index 000000000..6cc796481 --- /dev/null +++ b/configuration/tools/jra55_datasets/interp_jra55_ncdf_bilinear.py @@ -0,0 +1,441 @@ +#! /usr/bin/env python3 + +import xesmf as xe +from netCDF4 import Dataset +import argparse +import os +import numpy as np +from datetime import datetime + +###################################################### +###################################################### +def make_regridder(lon1, lat1, lon2, lat2, method, periodic, grdname, + lon1_b=None, lat1_b=None, lon2_b=None, lat2_b=None): + ''' + make nearest neighbor xESMF regridder object. + input: + lon1: source longitudes (degrees) + lat1: source latitudes (degrees) + lon2: target longitudes (degrees) + lat2: target latitudes (degrees) + method: regridding method (bilinear, patch, conservative, nearest_s2d) + periodic: True if periodic longitudes, false if not + grdname: filename for regridder (Ugrid or Tgrid) + ''' + if method != "conservative": + # define grids for regridder + grid1 = {'lon' : lon1, 'lat' : lat1} + grid2 = {'lon' : lon2, 'lat' : lat2} + + + else: + # conservative needs boundary lon/lat + grid1 = {'lon' : lon1, 'lat' : lat1, + 'lon_b' : lon1_b, 'lat_b' : lat1_b} + + grid2 = {'lon' : lon2, 'lat' : lat2, + 'lon_b' : lon2_b, 'lat_b' : lat2_b} + + # make regridder + # here specify reuse_weights=False to re-generate weight file. + # if wanted to reuse file inteas of making int, + # check if file exists and change use_file_weights=True. + # see commented out example below + use_file_weights=False + + # check if regrid file exists. + # If so, reuse file instead of regenerating. + # if (os.path.isfile(blin_grid_name)): + # use_file_weights = True + + regridder = xe.Regridder(ds_in=grid1,ds_out=grid2, + method=method, + periodic=periodic, + filename=grdname, + reuse_weights=use_file_weights) + + + return regridder + +######################################### +######################################### +def halo_extrapolate(a,ew_bndy_type,ns_bndy_type): + ''' + Extrapolate to 'halo' cell as in CICE code + ice_boundary.F90:ice_HaloExtrapolate. + inputs: + a: array nx+1, ny+1 (nghost/nhalo hard-coded as 1 for now) + ew_bndy_type: east/west boundary type (cyclic, regional, etc) + ns_bndy_type: norh/south boundary type (cyclic, regional, etc) + + return: a with halo applied + ''' + + # get dimension of a + # expected to be 0:nx+nghost, 0:ny+nghost + nj, ni = a.shape # note with Python NetCDF is nj, ni order + # W/E edges + if ew_bndy_type == 'cyclic': + a[: ,0] = a[:,-2] # -2, since -1 is ghost cell + a[:,-1] = a[:, 1] # 1, since 0 is ghost cell + else: # if (trim(ew_bndy_type) /= 'cyclic') then + a[:, 0] = 2.0*a[:, 1] - a[:, 2] + a[:,-1] = 2.0*a[:,-2] - a[:,-3] + + # south edge + if ns_bndy_type == 'cyclic': + a[0,:] = a[-2,:] # -2, since -1 is ghost cell + else: + a[0,:] = 2.0*a[1,:] - a[2,:] + + # north edge treated a little different, depending + # on if bndy type is tripole + if ns_bndy_type == 'cyclic': + a[-1,:] = a[1,:] # 1, since 0 is ghost cell + + elif (ns_bndy_type != 'cyclic' and + ns_bndy_type != 'tripole' and + ns_bndy_type != 'tripoleT'): + + a[-1,:] = 2.0*a[-2,:] - a[-3,:] + + else: + pass # do nothing + + # return array with halo upated + return a + +######################################### +######################################### + +def Tlatlon(ulat,ulon,ew_bndy_type,ns_bndy_type): + ''' + Make TLAT/TLON from ULAT/ULON. + see ice_grid.F90:Tlatlon for method + Inputs: + ulat: U grid latitude in degrees + ulon: U grid longitude in degrees + + output: + tlat in degrees + tlon in degrees + ''' + + # method obtained from ice_grid.F90: subroutine Tlatlon + ulatcos = np.cos(np.deg2rad(ulat)) + ulatsin = np.sin(np.deg2rad(ulat)) + + uloncos = np.cos(np.deg2rad(ulon)) + ulonsin = np.sin(np.deg2rad(ulon)) + + # initialize array with nghost=1 on each side + nj, ni = ulatcos.shape # note: Python NetCDF is nj, ni order + print("Tlatlon nj, ni", nj, ni) + + nghost = 1 + workdims = (nj+2*nghost,ni+2*nghost) + #print("Tlatlon workdims", workdims) + + ulatcos1 = np.zeros(workdims,dtype='f8') + ulatsin1 = np.zeros(workdims,dtype='f8') + uloncos1 = np.zeros(workdims,dtype='f8') + ulonsin1 = np.zeros(workdims,dtype='f8') + + # fill middle of work arrays + ulatcos1[1:nj+1,1:ni+1] = ulatcos + ulatsin1[1:nj+1,1:ni+1] = ulatsin + + # fill middle of work arrays + ulatcos1[1:nj+1,1:ni+1] = ulatcos + ulatsin1[1:nj+1,1:ni+1] = ulatsin + + uloncos1[1:nj+1,1:ni+1] = uloncos + ulonsin1[1:nj+1,1:ni+1] = ulonsin + + # fill halos + ulatcos1 = halo_extrapolate(ulatcos1,ew_bndy_type,ns_bndy_type) + ulatsin1 = halo_extrapolate(ulatsin1,ew_bndy_type,ns_bndy_type) + uloncos1 = halo_extrapolate(uloncos1,ew_bndy_type,ns_bndy_type) + ulonsin1 = halo_extrapolate(ulonsin1,ew_bndy_type,ns_bndy_type) + + # now do computations as in ice_grid.F90:Tlatlon + + # x, y, z are full 2d + x = uloncos1 * ulatcos1 + y = ulonsin1 * ulatcos1 + z = ulatsin1 + + tx = 0.25*(x[0:nj, 0:ni ] + # x1 + x[0:nj, 1:ni+1] + # x2 + x[1:nj+1,0:ni ] + # x3 + x[1:nj+1,1:ni+1]) # x4 + + #print("Tlonlat: x.shape", x.shape) + #print("Tlonlat: tx.shape", tx.shape) + + + ty = 0.25*(y[0:nj, 0:ni ] + # y1 + y[0:nj, 1:ni+1] + # y2 + y[1:nj+1,0:ni ] + # y3 + y[1:nj+1,1:ni+1]) # y4 + + + tz = 0.25*(z[0:nj, 0:ni ] + # z1 + z[0:nj, 1:ni+1] + # z2 + z[1:nj+1,0:ni ] + # z3 + z[1:nj+1,1:ni+1]) # z4 + + da = np.sqrt(tx*tx + ty*ty + tz*tz) + + tz = tz/da + + tlon = np.arctan2(ty,tx) + tlat = np.arcsin(tz) + + # returnd tlat, tlon in degrees + return np.rad2deg(tlat), np.rad2deg(tlon) + +########################## +########################## + +def get_command_line_args(): + ''' + argument parser for command line arguments + ''' + + dstr = "Interplate JRA55 data" + parser = argparse.ArgumentParser(description=dstr) + + # add arguments + parser.add_argument("JRADTG", type=str, help="JRA55 file date time group") + parser.add_argument("dstgrid", type=str, help="Destination grid file (NetCDF)") + parser.add_argument("ncout", type=str, help="Output file name (NetCDF)") + + + # get the arguments + args = parser.parse_args() + + # return values + return args.JRADTG, args.dstgrid, args.ncout + + +################################ +################################ + +def get_jra55_nc_dict(): + ''' + Create dictionary that links the NetCDF variable name + with the file prefix. The file prefix is appended by + JRADTG from command line + ''' + # specify dictionary with dataset prefix names + jra55dict = {"TPRAT_GDS4_SFC_ave3h" : "fcst_phy2m.061_tprat.reg_tl319", # precip + "DSWRF_GDS4_SFC_ave3h" : "fcst_phy2m.204_dswrf.reg_tl319", # downward shortwave + "DLWRF_GDS4_SFC_ave3h" : "fcst_phy2m.205_dlwrf.reg_tl319", # downward longwave + "TMP_GDS4_HTGL" : "fcst_surf.011_tmp.reg_tl319" , # air temp + "UGRD_GDS4_HTGL" : "fcst_surf.033_ugrd.reg_tl319" , # u velocity + "VGRD_GDS4_HTGL" : "fcst_surf.034_vgrd.reg_tl319" , # v velocity + "SPFH_GDS4_HTGL" : "fcst_surf.051_spfh.reg_tl319"} # specify humidity + + + return jra55dict + +################################ +################################ + +def get_jra55_cice_var(): + ''' + Make dictionary relating JRA55 NetCDF variables + to CICE variables. + ''' + + # specify output variable names + # This is for current CICE expected names + # it might be better to change CICE in long run + cice_var = {"TPRAT_GDS4_SFC_ave3h" : "ttlpcp", + "DSWRF_GDS4_SFC_ave3h" : "glbrad", + "DLWRF_GDS4_SFC_ave3h" : "dlwsfc", + "TMP_GDS4_HTGL" : "airtmp", + "UGRD_GDS4_HTGL" : "wndewd", + "VGRD_GDS4_HTGL" : "wndnwd", + "SPFH_GDS4_HTGL" : "spchmd"} + + return cice_var + +################################ +################################ + +def init_ncout(ncout,nc1,llat,llon): + + ''' + Initialize output NetCDF file + with proper units and dimensions. + ''' + + dsout = Dataset(ncout,'w',format='NETCDF3_64BIT_OFFSET') + + # get dimensions from size of lat + (nlat,nlon) = llat.shape + + # create dimensions + time = dsout.createDimension('time',None) # unlimited + dim_j = dsout.createDimension('dim_j',nlat) + dim_i = dsout.createDimension('dim_i',nlon) + + # create time variable. + # note is defined as 'times' (with and s) to not conflict + # with dimension 'time' + times = dsout.createVariable('time','f8',('time',)) + times.units = nc1['initial_time0_hours'].units + times.calendar = 'gregorian' + + # loop over nc1 times + dates = [] + dates.append(nc1['initial_time0_hours'][0] + nc1['forecast_time1'][1]) + # loop over remaining + for h in nc1['initial_time0_hours'][1:-1]: + for ft in nc1['forecast_time1'][:]: + dates.append(h + ft) + + # include only first forecast_time of last initial time + dates.append(nc1['initial_time0_hours'][-1] + nc1['forecast_time1'][0]) + + # write dates to file + times[:] = dates + + # create LON/LAT variables + LON = dsout.createVariable('LON','f8',('dim_j','dim_i',)) + LON.units = 'degrees_east' + + LAT = dsout.createVariable('LAT','f8',('dim_j','dim_i',)) + LAT.units = 'degrees_north' + + # write LON, LAT to file + LON[:] = llon[:,:] + LAT[:] = llat[:,:] + + + return dsout + +################################ +################################ + + +# main subroutine +if __name__ == "__main__": + + # get jra dtg and ncout from command line + JRADTG, dstgrid, ncout = get_command_line_args() + + # get jra55 variable/filename prefix dictionary + jra55dict = get_jra55_nc_dict() + + # get dictionary linking jra55 variables names + # and CICE forcing varible names + cice_var = get_jra55_cice_var() + + # read input grid. + # use one of the jra55 files. + # it is assumed all JRA data are the same grid for later + fname = f"{jra55dict['TMP_GDS4_HTGL']:s}.{JRADTG:s}.nc" + print("opening dataset ", fname) + grid1_ds = Dataset(fname,'r',format='NETCDF3_64BIT_OFFSET') + lon1 = grid1_ds['g4_lon_3'][:] # 1D + lat1 = grid1_ds['g4_lat_2'][:] # 1D + + # open destination grid + # here it is assumed a CICE NetCDF file. + # the user can update as appropriate + print("Opening ", dstgrid) + grid2_ds = Dataset(dstgrid,'r',format='NETCDF3_64BIT_OFFSET') + ulon2 = grid2_ds["lon"][:,:] # 2D. Assumed ULON in degrees + ulat2 = grid2_ds["lat"][:,:] # 2D. Assumed ULAT in degrees + if np.max(np.abs(ulat2)) < 10. : + ulon2 = np.rad2deg(ulon2) + ulat2 = np.rad2deg(ulat2) + + # make tgrid from ugrid + ew_bndy_type = 'cyclic' + ns_bndy_type = 'open' + tlat2, tlon2 = Tlatlon(ulat2,ulon2,ew_bndy_type,ns_bndy_type) + + # make regridders + print("making bilinear regridder") + method = 'bilinear' + periodic = True + blin_grid_name = 'bilinear_jra55_gx3.nc' + rgrd_bilinear = make_regridder(lon1,lat1,tlon2,tlat2, + method,periodic,blin_grid_name) + + # setup output dataset by adding lat/lon + dsout = init_ncout(ncout,grid1_ds,tlat2,tlon2) + + # no longer need grid1, grid2 + grid1_ds.close() + grid2_ds.close() + + # do the regridding + # Loop over all the files using regridder from above + # and add to dataout + for var, fprefix in jra55dict.items(): + fname = f"{fprefix:s}.{JRADTG:s}.nc" + print("reading ", fname) + jra_ds = Dataset(fname,'r',format='NETCDF3_CLASSIC') + + # make output variable + data = dsout.createVariable(cice_var[var],'f4',('time','dim_j','dim_i')) + + # do interpolation + print("Interpolating ", var) + + if var.find('ave3h') > 0: # ave3r in var + # use bilinear here + d = rgrd_bilinear(jra_ds[var][:,:,:,:]) + + # write to file in correct time order + for t in range(d.shape[0]): + for n in range(d.shape[1]): + #print('indx (2*t)+n = ', (2*t)+n) + data[(2*t)+n,:,:] = d[t,n,:,:] + + else: + # instantaneous use bilinear + d = rgrd_bilinear(jra_ds[var][:,:,:,:]) + + # write to file in correct time order. + # note need to write 2nd forecast_time first. + # in this case first forecast_time is NaN + data[0,:,:] = d[0,1,:,:] + for t in range(1,d.shape[0]-1): + for n in range(d.shape[1]): + #print('indx (2*t)+n-1 = ', (2*t)+n-1) + data[(2*t)+n-1,:,:] = d[t,n,:,:] + + # write first forecast time of last initial time + # second forecast time is NAN + data[-1,:,:] = d[-1,0,:,:] + + # add coordinates attribute + data.coordinates = "LON LAT time" + data.long_name = jra_ds[var].long_name + data.units = jra_ds[var].units + + precip_factor = 1. / 86400. + + # Convert mm / day to kg/m^2/s. + if var.find('PRAT') > 0: + data[:] = data[:] * precip_factor + data.units = 'kg/m2/s' + else: + data.units = jra_ds[var].units + + # close jra55 file + jra_ds.close() + + + # write tou output file + # close output file + dsout.close() + + print("Done") + diff --git a/configuration/tools/jra55_datasets/make_forcing.csh b/configuration/tools/jra55_datasets/make_forcing.csh new file mode 100755 index 000000000..c57871a25 --- /dev/null +++ b/configuration/tools/jra55_datasets/make_forcing.csh @@ -0,0 +1,49 @@ +#!/bin/csh +# ----- +# This is a script that worked on NCAR's cheyenne in March, 2021. +# It converts raw JRA55 datasets to a format that CICE can use. +# This tools is documented in the CICE user guide. The +# tool interpolates to a CICE grid and does things like convert units. +# ----- +# The interp_jra55_ncdf_bilinar.py script was placed in "scripts_dir" +# The raw JRA55 datasets were placed in "jra55_data_dir" +# The CICE grid files were places in "jra55_data_dir" +# The model output was created in "output_data_dir" +# ----- +#PBS -N make_forcing +#PBS -q regular +#PBS -l select=1:ncpus=4:mpiprocs=4 +#PBS -l walltime=06:00:00 +#PBS -A P93300665 + +set scripts_dir = "/glade/work/tcraig/cice-consortium/cice.jra55_tool/configuration/tools/jra55_datasets" +set jra55_data_dir = "/glade/scratch/dbailey/JRA_DATA/" +set output_data_dir = "/glade/scratch/tcraig/JRA_DATA_output" +set grid = "gx3" +set cice_grid_file = "grid_gx3.nc" + +module load python/3.7.9 +source /glade/u/apps/opt/ncar_pylib/ncar_pylib.csh default +module load nco + +mkdir -p ${output_data_dir} +cd ${output_data_dir} + +ln -s ${jra55_data_dir}/fcst_*.nc . +ln -s ${jra55_data_dir}/grid_*.nc . + +ln -s ${scripts_dir}/interp_jra55_ncdf_bilinear.py . + +#foreach year ( 1995 1996 1997 1998 1999 2000 2001 2002 2003 2004 ) +foreach year ( 1997 ) + +./interp_jra55_ncdf_bilinear.py ${year}010100_${year}033121 ${cice_grid_file} JRA55_${grid}_03hr_forcing_${year}-q1.nc +./interp_jra55_ncdf_bilinear.py ${year}040100_${year}063021 ${cice_grid_file} JRA55_${grid}_03hr_forcing_${year}-q2.nc +./interp_jra55_ncdf_bilinear.py ${year}070100_${year}093021 ${cice_grid_file} JRA55_${grid}_03hr_forcing_${year}-q3.nc +./interp_jra55_ncdf_bilinear.py ${year}100100_${year}123121 ${cice_grid_file} JRA55_${grid}_03hr_forcing_${year}-q4.nc + +ncrcat JRA55_${grid}_03hr_forcing_${year}-??.nc JRA55_${grid}_03hr_forcing_${year}.nc + +/bin/rm -f ${jra55_data_dir}/JRA55_${grid}_03hr_forcing_${year}-??.nc + +end diff --git a/doc/source/developer_guide/dg_scripts.rst b/doc/source/developer_guide/dg_scripts.rst index da5ef7d24..50853b3ea 100644 --- a/doc/source/developer_guide/dg_scripts.rst +++ b/doc/source/developer_guide/dg_scripts.rst @@ -161,25 +161,25 @@ To add a new test (for example newtest), several files may be needed, Generating a new test, particularly the **test_newtest.script** usually takes some iteration before it's working properly. -.. _dev_compliance: +.. _dev_validation: -Code Compliance Script +Code Validation Script ---------------------- -The code compliance test validates non bit-for-bit model changes. The directory -**configuration/scripts/tests/QC** contains scripts related to the compliance testing, -and this process is described in :ref:`compliance`. This section will describe a set -of scripts that test and validate the code compliance process. This should be done -when the compliance test or compliance test scripts (i.e., ``cice.t-test.py``) are modified. -Again, this section **documents a validation process for the compliance scripts**; it does not -describe to how run the compliance test itself. +The code validation (aka QC or quality control) test validates non bit-for-bit model changes. The directory +**configuration/scripts/tests/QC** contains scripts related to the validation testing, +and this process is described in :ref:`validation`. This section will describe a set +of scripts that test and validate the QC process. This should be done +when the QC test or QC test scripts (i.e., ``cice.t-test.py``) are modified. +Again, this section **documents a validation process for the QC scripts**; it does not +describe to how run the validation test itself. -Two scripts have been created to automatically validate the code compliance script. +Two scripts have been created to automatically validate the QC script. These scripts are: * ``gen_qc_cases.csh``, which creates the 4 test cases required for validation, builds the executable, and submits to the queue. -* ``compare_qc_cases.csh``, which runs the code compliance script on three combinations +* ``compare_qc_cases.csh``, which runs the QC script on three combinations of the 4 test cases and outputs whether or not the correct response was received. The ``gen_qc_cases.csh`` script allows users to pass some arguments similar @@ -216,7 +216,7 @@ To install the necessary Python packages, the ``pip`` Python utility can be used check to see if there is any Python module (``module avail python``) that you might need to load prior to using ``pip``. -To perform the validation, execute the following commands. +To perform the QC validation, execute the following commands. .. code-block:: bash diff --git a/doc/source/developer_guide/dg_tools.rst b/doc/source/developer_guide/dg_tools.rst new file mode 100644 index 000000000..ba29e0184 --- /dev/null +++ b/doc/source/developer_guide/dg_tools.rst @@ -0,0 +1,150 @@ +:tocdepth: 3 + +.. _tools: + +Tools +============= + + +.. _cice4restart: + +CICE4 restart conversion +------------------------- + +There is a Fortran program in **configuration/tools/cice4_restart_conversion** +that will help convert a CICE4 restart file into a CICE5 restart file. +There is a bit of documentation contained in that source code about how +to build, use, and run the tool. A few prognostic variables were changed +from CICE4 to CICE5 which fundamentally altered the fields saved to +the restart file. See +**configuration/tools/cice4_restart_conversion/convert_restarts.f90** +for additional information. + + +.. _jra55datasettool: + +JRA55 forcing datasets +------------------------ + +This section describes how to generate JRA55 forcing data for the CICE model. +Raw JRA55 files have to be interpolated and processed into input files specifically +for the CICE model. A tool exists in **configuration/tools/jra55_datasets** +to support that process. +The raw JRA55 data is obtained from the NCAR/UCAR Research Data Archive and +the conversion tools are written in python. + +Requirements +********************* + +Python3 is required, and the following +python packages are required with the tested version number in parenthesis. These +versions are not necessarily the only versions that work, they just indicate what +versions were used when the script was recently run. + +- python3 (python3.7.9) +- numpy (1.18.5) +- netCDF4 (1.5.5) +- ESMPy (8.0.0) +- xesmf (0.3.0) + +NCO is required for aggregating the output files into yearly files. + +- netcdf (4.7.4) +- nco (4.9.5) + +Raw JRA55 forcing data +************************* + +The raw JRA55 forcing data is obtained from the UCAR/NCAR Research Data Archive, +https://rda.ucar.edu/. You must first register (free) and then sign in. The +"JRA-55 Reanalysis Daily 3-Hourly and 6-Hourly Data" is ds628.0 and can be found here, +https://rda.ucar.edu/datasets/ds628.0. + +The "Data access" tabs will provide a list of product categories. +The JRA55 data of interest are located in 2 separate products. Winds, air +temperature, and specific humidity fields are included in "JRA-55 +3-Hourly Model Resolution 2-Dimensional Instantaneous Diagnostic Fields". +Precipitation and downward radiation fluxes are found in "JRA-55 3-Hourly +Model Resolution 2-Dimensional Average Diagnostic Fields". (Note the +difference between instantaneous and averaged data products. There are several +JRA55 datasets available, you will likely have to scroll down the page to find +these datasets.) Data are also available on a coarser 1.25° grid, but the tools +are best used with the native TL319 JRA55 grid. + +The fields needed for CICE are + +- specific humidity (3-hourly instantaneous), Qa +- temperature (3-hourly instantaneous), Tair +- u-component of wind (3-hourly instantaneous), uatm +- v-component of wind(3-hourly instantaneous), vatm +- downward longwave radiation flux (3 hourly average), flw +- downward solar radiation flux (3 hourly average), fsw +- total precipitation (3 hourly average), fsnow + +To customize the dataset for download, choose the “Get a Subset” option. Select +the desired times in the “Temporal Selection” section, then click on desired parameters +(see list above). After clicking continue, select Output Format "Converted to NetCDF". + +Once the data request is made, an email notification will be sent with a dedicated +URL that will provide a variety of options for downloading the data remotely. +The data will be available to download for 5 days. +The raw data consists of multiple files, each containing three months of data for +one field. + + +Data conversion +************************* + +The script, **configuration/tools/jra55_datasets/interp_jra55_ncdf_bilinear.py**, +converts the raw data to CICE input files. + +The script uses a bilinear regridding algorithm to regrid from the JRA55 grid to +the CICE grid. The scripts use the Python package ‘xesmf’ to generate bilinear +regridding weights, and these regridding weights are written to the file defined by +the variable "blin_grid_name" in **interp_jra55_ncdf_bilinear.py**. This filename +can be modified by editing **interp_jra55_ncdf_bilinear.py**. +The weights file can be re-used if interpolating different data on the same grid. +Although not tested in this version of the scripts, additional regridding options +are available by xesmf, including ‘conservative’ and ‘nearest neighbor’. These +methods have not been tested in the current version of the scripts. The reader +is referred to the xESMF web page for further documentation +(https://xesmf.readthedocs.io/en/latest/ last accessed 5 NOV 2020). + +To use the **interp_jra55_ncdf_bilinear** script, do :: + + python3 interp_jra55_ncdf_bilinear.py –h + +to see the latest interface information :: + + usage: interp_jra55_ncdf_bilinear.py [-h] JRADTG gridout ncout + + Interpolate JRA55 data to CICE grid + + positional arguments: + JRADTG JRA55 input file date time group + gridout CICE grid file (NetCDF) + ncout Output NetCDF filename + + optional arguments: + -h, --help show this help message and exit + +Sample usage is :: + + ./interp_jra55_ncdf_bilinear.py 1996010100_1996033121 grid_gx3.nc JRA55_gx3_03hr_forcing_1996-q1.nc + ./interp_jra55_ncdf_bilinear.py 1996040100_1996063021 grid_gx3.nc JRA55_gx3_03hr_forcing_1996-q2.nc + ./interp_jra55_ncdf_bilinear.py 1996070100_1996093021 grid_gx3.nc JRA55_gx3_03hr_forcing_1996-q3.nc + ./interp_jra55_ncdf_bilinear.py 1996100100_1996123121 grid_gx3.nc JRA55_gx3_03hr_forcing_1996-q4.nc + +In this case, the 4 quarters of 1996 JRA55 data is going to be interpolated to the gx3 grid. +NCO can be used to aggregate these files into a single file :: + + ncrcat JRA55_gx3_03hr_forcing_1996-??.nc JRA55_${grid}_03hr_forcing_1996.nc + +NOTES + +- The scripts are designed to read a CICE grid file in netCDF format. This is the "grid_gx3.nc" file above. The NetCDF grid names are hardcoded in **interp_jra55_ncdf_bilinear.py**. If you are using a different grid file with different variable names, this subroutine needs to be updated. +- All files should be placed in a common directory. This includes the raw JRA55 input files, the CICE grid file, and **interp_jra55_ncdf_bilinear.py**. The output files will be written to the same directory. +- The script **configuration/tools/jra55_datasets/make_forcing.csh** was used on the NCAR cheyenne machine in March, 2021 to generate CICE forcing data. It assumes the raw JRA55 is downloaded, but then sets up the python environment, links all the data in a common directory, runs **interp_jra55_ncdf_bilinear.py** and then aggregates the quarterly data using NCO. +- The new forcing files can then be defined in the **ice_in** namelist file using the input variables, ``atm_data_type``, ``atm_data_format``, ``atm_data_dir``, ``fyear_init``, and ``ycycle``. See :ref:`forcing` for more information. +- The total precipitation field is mm/day in JRA55. This field is initially read in as snow, but prepare_forcing in **ice_forcing.F90** splits that into rain or snow forcing depending on the air temperature. + diff --git a/doc/source/developer_guide/index.rst b/doc/source/developer_guide/index.rst index ab5b2d1e6..6fc3356f4 100644 --- a/doc/source/developer_guide/index.rst +++ b/doc/source/developer_guide/index.rst @@ -17,5 +17,6 @@ Developer Guide dg_forcing.rst dg_icepack.rst dg_scripts.rst + dg_tools.rst dg_other.rst diff --git a/doc/source/user_guide/ug_running.rst b/doc/source/user_guide/ug_running.rst index 541fa81a4..a79e2ccb2 100644 --- a/doc/source/user_guide/ug_running.rst +++ b/doc/source/user_guide/ug_running.rst @@ -36,12 +36,19 @@ The Consortium has tested the following compilers at some point, - Intel 17.0.2.174 - Intel 17.0.5.239 - Intel 18.0.1.163 +- Intel 18.0.5 - Intel 19.0.2 - Intel 19.0.3.199 +- Intel 19.1.0.166 +- Intel 19.1.1.217 - PGI 16.10.0 +- PGI 19.9-0 +- PGI 20.1-0 - GNU 6.3.0 - GNU 7.2.0 - GNU 7.3.0 +- GNU 8.3.0 +- GNU 9.3.0 - Cray 8.5.8 - Cray 8.6.4 - NAG 6.2 @@ -54,22 +61,33 @@ The Consortium has tested the following mpi versions, - MPICH 7.6.3 - MPICH 7.7.6 - Intel MPI 18.0.1 +- Intel MPI 18.0.4 +- Intel MPI 2019 Update 6 - MPT 2.14 - MPT 2.17 - MPT 2.18 - MPT 2.19 +- MPT 2.20 +- MPT 2.21 +- mvapich2-2.3.3 - OpenMPI 1.6.5 +- OpenMPI 4.0.2 The NetCDF implementation is relatively general and should work with any version of NetCDF 3 or 4. The Consortium has tested - NetCDF 4.3.0 - NetCDF 4.3.2 - NetCDF 4.4.0 -- NetCDF 4.4.1.1.32 +- NetCDF 4.4.1.1.3 - NetCDF 4.4.1.1 - NetCDF 4.4.2 - NetCDF 4.5.0 +- NetCDF 4.5.2 - NetCDF 4.6.1.3 +- NetCDF 4.6.3 +- NetCDF 4.6.3.2 +- NetCDF 4.7.2 +- NetCDF 4.7.4 Please email the Consortium if this list can be extended. @@ -772,7 +790,7 @@ A few notes about the conda configuration: - It is not recommeded to run other test suites than ``quick_suite`` or ``travis_suite`` on a personal computer. - The conda environment is automatically activated when compiling or running the model using the ``./cice.build`` and ``./cice.run`` scripts in the case directory. These scripts source the file ``env.conda_{linux.macos}``, which calls ``conda activate cice``. -- To use the "cice" conda environment with the Python plotting (see :ref:`timeseries`) and quality control scripts (see :ref:`CodeCompliance`), you must manually activate the environment: +- To use the "cice" conda environment with the Python plotting (see :ref:`timeseries`) and quality control (QC) scripts (see :ref:`CodeValidation`), you must manually activate the environment: .. code-block:: bash @@ -897,7 +915,7 @@ To use the ``timeseries.py`` script, the following requirements must be met: * matplotlib Python package * datetime Python package -See :ref:`CodeCompliance` for additional information about how to setup the Python +See :ref:`CodeValidation` for additional information about how to setup the Python environment, but we recommend using ``pip`` as follows: :: pip install --user numpy diff --git a/doc/source/user_guide/ug_testing.rst b/doc/source/user_guide/ug_testing.rst index 4eb7b6795..4a4ced555 100644 --- a/doc/source/user_guide/ug_testing.rst +++ b/doc/source/user_guide/ug_testing.rst @@ -525,7 +525,7 @@ Test Suite Examples This will compare to results saved in the baseline [bdir] directory under the subdirectory cice.v01a. With the ``--bcmp`` option, the results will be tested against prior baselines to verify bit-for-bit, which is an important step prior - to approval of many (not all, see :ref:`compliance`) Pull Requests to incorporate code into + to approval of many (not all, see :ref:`validation`) Pull Requests to incorporate code into the CICE Consortium master code. You can use other regression options as well. (``--bdir`` and ``--bgen``) @@ -771,9 +771,9 @@ assess test coverage. ..in the future. -.. _compliance: +.. _validation: -Code Compliance Test (non bit-for-bit validation) +Code Validation Test (non bit-for-bit validation) ---------------------------------------------------- A core tenet of CICE dycore and CICE innovations is that they must not change @@ -898,7 +898,7 @@ autocorrelation :math:`r_1`. .. _quadratic: -Quadratic Skill Compliance Test +Quadratic Skill Validation Test ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ In addition to the two-stage test of mean sea ice thickness, we also @@ -982,12 +982,12 @@ hemispheres, and must exceed a critical value nominally set to test and the Two-Stage test described in the previous section are provided in :cite:`Hunke18`. -.. _CodeCompliance: +.. _CodeValidation: -Code Compliance Testing Procedure +Code Validation Testing Procedure ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ -The CICE code compliance test is performed by running a python script +The CICE code validation (QC) test is performed by running a python script (**configurations/scripts/tests/QC/cice.t-test.py**). In order to run the script, the following requirements must be met: @@ -1001,7 +1001,7 @@ QC testing should be carried out using configurations (ie. namelist settings) th exercise the active code modifications. Multiple configurations may need to be tested in some cases. Developers can contact the Consortium for guidance or if there are questions. -In order to generate the files necessary for the compliance test, test cases should be +In order to generate the files necessary for the validation test, test cases should be created with the ``qc`` option (i.e., ``--set qc``) when running cice.setup. This option results in daily, non-averaged history files being written for a 5 year simulation. @@ -1013,7 +1013,7 @@ To install the necessary Python packages, the ``pip`` Python utility can be used pip install --user numpy pip install --user matplotlib -To run the compliance test, setup a baseline run with the original baseline model and then +To run the validation test, setup a baseline run with the original baseline model and then a perturbation run based on recent model changes. Use ``--set qc`` in both runs in addition to other settings needed. Then use the QC script to compare history output,