Modifying a coordinate

[JohnBrouillet]

Hi!

Sorry for not being able to totally respect #5404 about creating the data within the example since I don't know how to create two datasets that give the results exposed below except one.

We find a strange behavior while working with this file from ucar and another file that has been produced by a coworker. This has been tested with these versions : 0.16.1, 0.20.1 and 0.21.1.

ds = xr.open_dataset('SOI.nc', decode_times=False)  
print(ds['time'].values[:3])  
ds['time'].values[:3] = 1  
print(ds['time'].values[:3])

[186601 186602 186603]
[186601 186602 186603]

The other file gives:

[525972. 526044. 526116.]
[1. 1. 1.]

With another file downloaded from CMEMS, we also observe no modification of time.

[1293969600 1294056000 1294142400]
[1293969600 1294056000 1294142400]

But after modifying it with cdo remapbil, we get:

[1.2939696e+09 1.2940560e+09 1.2941424e+09]
[1. 1. 1.]

Otherwise, casting the time coordinate seems to make the job:

ds = xr.open_dataset('SOI.nc', decode_times=False)  
print(ds['time'].values[:3])  
ds['time'] = ds['time'].astype(np.int64)
ds['time'].values[:3] = 1  
print(ds['time'].values[:3])

[186601 186602 186603]
[1 1 1]

With the dataset example from xarray doc and time adapted:

np.random.seed(0)
temperature = 15 + 8 * np.random.randn(2, 2, 3)
precipitation = 10 * np.random.rand(2, 2, 3)
lon = [[-99.83, -99.32], [-99.79, -99.23]]
lat = [[42.25, 42.21], [42.63, 42.59]]
time = [186601, 186602, 186603]

ds = xr.Dataset(
    data_vars=dict(
        temperature=(["x", "y", "time"], temperature),
        precipitation=(["x", "y", "time"], precipitation),
    ),
    coords=dict(
        lon=(["x", "y"], lon),
        lat=(["x", "y"], lat),
        time=time,
    ),
   attrs=dict(description="Weather related data."),
)
print(ds['time'].values[:3])
ds['time'].values[:3] = 1.0
print(ds['time'].values[:3])

[186601 186602 186603]
[1 1 1]

To sum up, for some files we can change some values of its coordinates with .values[a:b] but for others, we can't. Altering (like a cast) the all coordinate changes its behavior regarding .values[a:b]. For any other variables, I got [1 1 1].

Is this a known issue? It seems to exist many ways to change a coordinate. If I understand correctly a coordinate is handled by a PandasIndex that is immutable. Thus using .values should not be allowed and besides it raises the question of assign_coords usage since we can just do ds['mycoord'].values = np.array([....]) like for every variables. What do you think about it? Is it possible that a netcdf file has some kind of a memory lock for the coordinates?

I take the freedom to use this discussion to ask you if it exists an ULM diagram or other kinds of graphs to see xarray internal architecture. I am really interested to better understand how it works.

Thank you very much !

[JohnBrouillet]

After trying to dig up in the code with the debugger, I looked closer at the solution with astype and I found something. There is something wrong with int32.

Look at this test made with SOI.nc from ucar. time is a int32, I cast it to a int64 and it works. Casting back to int32, it does not work anymore.

import xarray as xr
import numpy as np

ds = xr.open_dataset('SOI.nc')
print('SOI.nc')
print(ds['time'].values[:3])
ds['time'].values[:3] = 1
print(ds['time'].values[:3])

ds = xr.open_dataset('SOI.nc')
ds['time'] = ds['time'].astype(np.int64)
ds.to_netcdf('SOI_int64.nc')

vs = xr.open_dataset('SOI_int64.nc')
print('SOI_int64.nc')
print(vs['time'].values[:3])
vs['time'].values[:3] = 1
print(vs['time'].values[:3])

vs = xr.open_dataset('SOI_int64.nc')
vs['time'] = vs['time'].astype(np.int32)
vs.to_netcdf('SOI_backtoInt32.nc')

ts = xr.open_dataset('SOI_backtoInt32.nc')
print('SOI_backtoInt32.nc')
print(ts['time'].values[:3])
ts['time'].values[:3] = 1
print(ts['time'].values[:3])

The result is:

SOI.nc
[186601 186602 186603]
[186601 186602 186603]

SOI_int64.nc
[186601 186602 186603]
[1 1 1]

SOI_backtoInt32.nc
[186601 186602 186603]
[186601 186602 186603]

I tried again with a CMEMS file. Time coordinates is a datetime64[ns], with decode_times=False, we get an int32.

vs = xr.open_dataset('METOFFICE-GLO-SST-L4-NRT-OBS-SST-MON-V2_1593889554066.nc', decode_times=False)
print(vs['time'].values[:3])
vs['time'].values[:3] = 1
print(vs['time'].values[:3])

[1132056000 1134648000 1137326400]
[1132056000 1134648000 1137326400]

vs = xr.open_dataset('METOFFICE-GLO-SST-L4-NRT-OBS-SST-MON-V2_1593889554066.nc', decode_times=True)
print(vs['time'].values[:3])
vs['time'].values[:3] = 1
print(vs['time'].values[:3])

['2016-11-15T12:00:00.000000000' '2016-12-15T12:00:00.000000000'
 '2017-01-15T12:00:00.000000000']
['1970-01-01T00:00:00.000000001' '1970-01-01T00:00:00.000000001'
 '1970-01-01T00:00:00.000000001']

We can see that time is modified when it is a datetime64[ns] but not when it is an int32.

[JohnBrouillet]

I tried every numpy dtypes with the dataset example from the doc. A few values of a coords can't be changed if its type size is less than 64 bits. Unicode strings coordinates can't be changed.

import xarray as xr
import numpy as np

np.random.seed(0)
temperature = 15 + 8 * np.random.randn(2, 2, 3)
lon = [[-99.83, -99.32], [-99.79, -99.23]]
lat = [[42.25, 42.21], [42.63, 42.59]]

tt = [np.int8, np.int16, np.int32, np.int64, np.float16, np.float32, np.float64, np.uint8, np.uint16, np.uint32, np.uint64]

for t in tt:
    time = np.array([100, 55, 60], dtype=t)

    ds = xr.Dataset(
        data_vars=dict(
            temperature=(["x", "y", "time"], temperature),
        ),
        coords=dict(
            lon=(["x", "y"], lon),
            lat=(["x", "y"], lat),
            time=time,
        )
    )
    print(str(t))
    print(ds['time'].values[:3])
    ds['time'].values[:3] = 1
    print(ds['time'].values[:3])

<class 'numpy.int8'>
[100  55  60]
[100  55  60]

<class 'numpy.int16'>
[100  55  60]
[100  55  60]

<class 'numpy.int32'>
[100  55  60]
[100  55  60]

<class 'numpy.int64'>
[100  55  60]
[1 1 1]

<class 'numpy.float16'>
[100.  55.  60.]
[100.  55.  60.]

<class 'numpy.float32'>
[100.  55.  60.]
[100.  55.  60.]

<class 'numpy.float64'>
[100.  55.  60.]
[1. 1. 1.]

<class 'numpy.uint8'>
[100  55  60]
[100  55  60]

<class 'numpy.uint16'>
[100  55  60]
[100  55  60]

<class 'numpy.uint32'>
[100  55  60]
[100  55  60]

<class 'numpy.uint64'>
[100  55  60]
[1 1 1]