Xarray-exported nc to GMT

Q&A
1

hi,

i want to use some magic from GMT on data exported from xarray.
it is a cube having time, lat, lon and several layers (dataarrays).
grdinfo shows only one of the layers (dataarray ‘vel’) and actions like grdblend (to merge two such datacubes) generates the output without errors, but again, only ‘vel’ is processed.
note that the netcdf is read normally in QGIS etc.
exporting to NetCDF3 or ‘classic’ didn’t help.

does anyone have experience with such? i will appreciate some comments on proper export from xarray to ‘standard GMT-readible format’.
the structure is:

(i use only excerpt here of 5x5x5 points)

<xarray.Dataset>
Dimensions:      (lat: 5, lon: 5, time: 5)
Coordinates:
  * lat          (lat) float64 49.76 49.76 49.76 49.76 49.76
  * lon          (lon) float64 18.4 18.4 18.4 18.4 18.4
  * time         (time) datetime64[ns] 2017-12-21 2017-12-27 ... 2018-01-14
    spatial_ref  int64 0
Data variables:
    cum          (time, lat, lon) float32 17.702633 -8.929522 ... -21.63409
    vel          (lat, lon) float64 -5.039 -7.906 -19.88 ... -17.61 -14.55
    vintercept   (lat, lon) float64 4.582 -6.325 7.294 ... 16.18 17.77 6.928
    coh          (lat, lon) float64 0.3168 0.3031 0.3074 ... 0.3187 0.3245
    rms          (lat, lon) float64 10.46 9.867 9.351 9.69 ... 8.953 8.673 8.515
    vstd         (lat, lon) float64 0.4425 0.5942 0.3812 ... 0.422 0.5334 0.4843
    vel_filt     (lat, lon) float64 -5.6 -7.543 -14.49 ... -16.44 -16.4 -14.46
Attributes:
    inc_angle:     33.518
    grid_mapping:  spatial_ref

output of the full file from grdinfo:

gmt grdinfo test.nc 
test.nc: Title: 
test.nc: Command: 
test.nc: Remark: 
test.nc: Gridline node registration used [Cartesian grid]
test.nc: Grid file format: nf = GMT netCDF format (32-bit float), CF-1.7
test.nc: x_min: 18.3719179798 x_max: 18.772457238 x_inc: 0.000515494540738 name: lon n_columns: 778
test.nc: y_min: 49.8242542261 y_max: 49.9711701702 y_inc: 0.000515494540738 name: lat n_rows: 286
test.nc: z_min: 0 z_max: 0 name: vel
test.nc: scale_factor: 1 add_offset: 0
test.nc: format: netCDF-4 chunk_size: 0,0 shuffle: off deflation_level: 0

thank you in advance!

M.

2

See the CookBook on how to access multi-layered nc grids. Basically append ?[layer_num] to the name and where layer_num is the layer number start counting on 0.

3

thank you for fast reaction Joaquim.
1, i would expect that without subsetting by ?[num], all layers should be processed by GMT
2, i just tried applying the subsetting (used 0, 1, 2) - the output grds are all same = only ‘vel’ is exported every time

4

Did you read the docs? I’m pretty sure GMT can read any layer of any variable but I don’t have the answer on top of my head.

5

you’re right - the cookbook is written very well.
i have checked it now, and indeed, providing .nc?coh (or other layer) works fine (except of layers with time coords, but i guess i can find answer in docs)
so the problem here is that i expected grdblend to merge all layers at once. instead, it seems it works only for single layers… i sincerely hope i am wrong and there is a way to have multilayer functionality in GMT codes.

6

GMT works on single layers. I don’t even see the sense of blending multiple layers that, when 4rth dim is time, all cover the same region.
And yes, it works also to access individual time layers.

7

many thanks Joaquim.
perhaps my last note on this - say i have A.nc and B.nc - they have both 5 layers (same names, same order). they cover a bit different areas spatially, they have a small overlap.
by blending to C.nc, I want to have all 5 layers (resampled) in C.nc.

at this moment, i understand that i should do this separately for each layer (the one with time dimensions was not working, but it may be some minor problem now) and then combine the layers to one C.nc. actually i do not see the sense of doing this one by one

8

easy :smile: write a code to do it (quite easy in Juli for example). But there will be another problem. GMT does not have yet a tool to combine the several layer in one single nc file. You’ll have to use another program to do that.