grafs_stacker.py

import glob
import os
from datetime import datetime 
import numpy as np
import xarray as xr
import netCDF4
import argparse

def extract_time(file_name):
    fname = os.path.splitext(os.path.basename(file_name))[0]
    date_str = fname.split("_")[-1]
    
    return datetime.strptime(date_str, '%Y%m%d')

def stack_year(source_path, year):
    cube = None
    t_dim = []

    files = glob.glob("{}/API_analysis_window_{}*.nc".format(source_path, year))
    files.sort()

    for f in files:
        t_dim.append(extract_time(f))
        ds = xr.open_dataset(f)

        if cube is None:
            cube = ds["API"][:]
        else:
            cube = np.dstack((cube, ds["API"][:]))

    return cube, t_dim


def pack_api(stack, dates, dest):
    
    with netCDF4.Dataset(dest, 'w', format='NETCDF4_CLASSIC') as ds:
        with open('grafs_metadata.json') as data_file:
            attrs = json.load(data_file)
            for key in attrs:
                setattr(ds, key, attrs[key])
        setattr(ds, "date_created", datetime.now().strftime("%Y%m%dT%H%M%S"))
        
        t_dim = ds.createDimension("time", len(dates))
        lon_dim = ds.createDimension("lon", stack.shape[1])
        lat_dim = ds.createDimension("lat", stack.shape[0])

        var = ds.createVariable("time", "f8", ("time",))
        var.units = "seconds since 1970-01-01 00:00:00.0"
        var.calendar = "standard"
        var.long_name = "Time, unix time-stamp"
        var.standard_name = "time"
        print(t_dim)
        var[:] = netCDF4.date2num(dates, units="seconds since 1970-01-01 00:00:00.0", calendar="standard")

        var = ds.createVariable("lon", "f8", ("lon",))
        var.units = "degrees east"
        var.long_name = "longitude"
        var.standard_name = "longitude"
        var[:] = np.linspace(-179.95, 179.95, num=3600)
        
        var = ds.createVariable("lat", "f8", ("lat",))
        var.units = "degrees north"
        var.long_name = "latitude"
        var.standard_name = "latitude"
        var[:] = np.linspace(59.95, -59.95, num=1200)
        
        var = ds.createVariable("API", 'f4', ("time", "lat", "lon"), fill_value=-9999.9)
        var.long_name = "API"
        var[:] = stack


if __name__ == "__main__":
    parser = argparse.ArgumentParser(description="""API data packer""")
    parser.add_argument('-y', '--year', type=int, required=True, help="Year to pack [YYYY].")
    parser.add_argument('-src', '--source', type=str, required=True, help="Path to source folder containing API data sets.")
    parser.add_argument('-dst', '--destination', type=str, required=True, help="Path to destination file to write.")
    args = parser.parse_args()

    stack, t_dim = stack_year(args.source, args.year)
    pack_api(stack, t_dim, args.destination)

Example: /g/data/oe9/software/miniconda3/bin/python ~/grafs_stacker.py -y 2017 -src /g/data/fj4/SatelliteSoilMoistureProducts/GRZSMFS/API_ANALYSIS_WINDOW -dst /g/data/ub8/global/GRAFS/GRAFS_2017.nc