1

I am merging some ERA5-Land and ERA5 (netCDF) time series:

  1. extract individual ERA5 and ERA5-Land maps and apply scale and offset factors -- this by looping the following
for BAND in $(seq 1 8784) ;do
    gdal_translate -unscale -b $BAND era5.nc era5_band_${BAND}.nc
done

and likewise

for BAND in $(seq 1 8784) ;do
    gdal_translate -unscale -b $BAND era5land.nc era5land_band_${BAND}.nc
done

This, however, does not preserve the NETCDF_DIM_time stamp of each map.

  1. merge corresponding maps in one netCDF by copying ERA5-Land over ERA5 data and resampling to the spatial resolution of ERA5-Land data -- this by looping over:
for BAND in $(seq 1 8784) ;do
    gdal_merge.py -ps 0.1 0.1 -o era5land_band_${BAND}_patched.nc era5_band_${BAND}.nc era5land_band_${BAND}.nc 
done

This produces thousands of files (as many as the time steps in the initial time series). A small subset list:

era5land_band_1_patched.nc
era5land_band_2_patched.nc
era5land_band_3_patched.nc
era5land_band_4_patched.nc
era5land_band_5_patched.nc
era5land_band_6_patched.nc
era5land_band_7_patched.nc
era5land_band_8_patched.nc
era5land_band_9_patched.nc

Finally, I would like to stack them back to a single netCDF file. I went on and updated the variable name and long_name using hints from

ncrename -v Band1,t2m era5land_band_1_patched.nc 
ncatted -O -a long_name,t2m,o,c,'2 metre temperature' era5land_band_1_patched.nc

Given both source data ERA5 and ERA5-Land, have the same temporal resolution and time-stamps, how is it recommended building a time series with the merged maps and the original time-stamps?

I would like to stay practical with GDAL and Bash (also due to the availability of GNU Parallel). Are there better alternatives (simplicity first, speed second)?

1 Answer 1

0

Quick UPDATE of the workflow!

I am merging some ERA5-Land and ERA5 (netCDF) time series using GDAL and Bash as follows:

  1. apply scale and offset factors and extract individual ERA5 and ERA5-Land maps and -- this by looping the following
for BAND in $(seq 1 8784) ;do
    gdal_translate -ot Float32 -unscale -b $BAND era5.nc era5_band_${BAND}.nc
done

and likewise

for BAND in $(seq 1 8784) ;do
    gdal_translate -ot Float32 -unscale -b $BAND era5land.nc era5land_band_${BAND}.nc
done

This, however, does not preserve the NETCDF_DIM_time stamp of each map.

  1. merge corresponding maps in one netCDF by copying ERA5-Land over ERA5 data and resample, at the same time, to the spatial resolution of ERA5-Land data -- this by looping over:
for BAND in $(seq 1 8784)
    gdal_merge.py \
      -ot Float32 \
      -ps 0.1 0.1 \
      -o era5land_band_${BAND}_patched.nc \
      era5_band_${BAND}.nc \
      era5land_band_${BAND}.nc 
done

This produces thousands of files (as many as the time steps in the initial time series). A small subset list:

era5land_band_1_patched.nc
era5land_band_2_patched.nc
era5land_band_3_patched.nc
era5land_band_4_patched.nc
era5land_band_5_patched.nc
era5land_band_6_patched.nc
era5land_band_7_patched.nc
era5land_band_8_patched.nc
era5land_band_9_patched.nc

Given both source data ERA5 and ERA5-Land have the same temporal resolution and time-stamps,

  1. Compute effective scale and offset values based on all input maps via the following formula (Python syntax below)
scale_factor = (data_maximum - data_minimum) / (2 ** bits - 1)
add_offset = data_minimum

where minimum and maximum values are the overall minimum and maximum from all maps

  1. Scale
# define $OUT_MAP, something like
OUT_MAP=$(basename $IN_MAP .nc)_scaled.nc
for IN_MAP in $(/usr/bin/ls era5land_band_*_patched.nc) ;do
  EXPRESSION="(A - $ADD_OFFSET) / $SCALE_FACTOR"
  gdal_calc.py \
      --format NETCDF \
      -A $IN_MAP \
      --outfile=$OUT_MAP \
      --calc="$EXPRESSION" \
      --type=Int16 \
      --quiet
  ncatted -O -h -a add_offset,Band1,c,f,$ADD_OFFSET $OUT_MAP
  ncatted -O -h -a scale_factor,Band1,c,f,$SCALE_FACTOR $OUT_MAP

done

also, create a function with the above and use parallel!

  1. Time-stamp and add-metadata for each individual patched map
# define the $OUT, something like:
HOURS_SINCE=$(gdalinfo "$cds_era5_land_2m_temperature_${YEAR}.nc" |grep -m 3 'NETCDF_DIM_time' |tail -n 1 |cut -d '=' -f2)
for BAND in $(/usr/bin/ls era5land_band_*_patched_scaled.nc) ;do
  IN_TIME=$(basename $BAND .nc)_timedimension.nc
  OUT=$(basename $BAND .nc)_timestamped.nc
  ncecat -u time $1 $IN_TIME  # add time dimension
  ncap2 \
     -h \
     -s "time[time]=array($HOURS_SINCE,1,\$time)" \
     -s 'time@units="hours since 1900-01-01 00:00:00.0"' \
     -O $IN_TIME $OUT
  ncrename -v Band1,t2m $OUT
  ncatted -a long_name,t2m,o,c,'2 metre temperature' $OUT
  ncatted -a units,t2m,o,c,K $OUT
  HOURS_SINCE=$(($HOURS_SINCE + 1))
done

Get HOURS_SINCE from the original file, increase by 1 after each time-stamping iteration, for example:



  9. Finally, merge maps to a time series in a single netCDF file

cdo mergetime era5_and_land_t2m_${YEAR}band*patched_scaled_timestamped.nc era5_and_land_t2m${YEAR}.nc


The workflow doesn't seem elegant (mixture of GDAL, NCO and CDO tools).

- Are there better alternatives (simplicity first, speed second)?
- With CDO, an input and an output is always expected (?). How to edit netCDF attributes of a file using CDO tools in-place?
- If best to use CDO/NCO, how to repeat this exercise using exclusively CDO/NCO?


## Useful

Don't create .xml files when computing statistics!

```bash
GDAL_PAM_ENABLED=NO
export GDAL_PAM_ENABLED

and some helper function

function get_number_of_bands() {
    gdalinfo $1 \
    |grep NETCDF_DIM_time_DEF \
    |cut -d'=' -f2 \
    |tr -d '{' \
    |cut -d',' -f1 ;
}
4
  • While the timestamping part works fine for a few files/maps, it takes very long for when the maps to stamp are thousands, like shown in the examples above. What alternative way can be faster to add a timestamp to each map? Mar 17, 2022 at 22:49
  • Is this an update or and answer?
    – Binx
    May 26, 2023 at 18:06
  • 1
    @Binx It works, as a process, albeit not claiming to be an optimal one. I have updated "my" workflow off-line and I can share it if you are interested. It is a mixture of Python, Bash (as a glue and more), Dask for one part, CDO/NCO. However, designed to work with the JRC's high-throughput computing cluster BDAP. May 30, 2023 at 10:53
  • If it solves your original problem, could you mark it as the answer?
    – Binx
    May 30, 2023 at 14:43

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.