Thanks for the tip! It works well.

Is it better to raise a warning or something else to remind users there're nan values in the mask?

@aidanheerdegen Thanks for the code. I suppose it's better to mention this method for DataArray in User Guide. @dcherian Should I create a PR for example like this?

``` air = xr.tutorial.open_dataset("air_temperature")['air']

air.isel(lon=10, lat=[19, 21, 22]).plot.line(x="time", marker='o',linewidth=0.,markersize=1) ```

@andersy005 Thanks a lot, I realize that it's mentioned in the comments of Guide.

Hi @TomNicholas, thanks and yes that's the same issue. Shall we close this duplicated one?

Old figure:

New figure:

Ha, thanks. It makes sense now. Shall we close this?

@jklymak Thanks for the explanation.

To get the old behaviour you simply need to do pcolormesh(x, y, Z[:-1, :-1], shading='flat') or make x and y one larger than Z in each dimension and specify the the corners of the quadrilaterals.

Method1: Subset value

This method works for the xarray tutorial data, but not for the TROPOMI polar-orbiting satellite data.

``` %matplotlib inline

import xarray as xr import cartopy.crs as ccrs import matplotlib.pyplot as plt

plt.figure(figsize=(14,6)) ax = plt.axes(projection=ccrs.PlateCarree())

ds = xr.open_dataset('./S5P_OFFL_L2__NO2____20190810T212136_20190810T230306_09456_01_010302_20190816T233944.nc', group='PRODUCT').isel(time=0)

m = ax.pcolormesh(ds['longitude'], ds['latitude'], ds['nitrogendioxide_tropospheric_column'][:-1, :-1], # ds['nitrogendioxide_tropospheric_column'], # shading='auto', transform=ccrs.PlateCarree(), vmin=0, vmax=1e-4, cmap='Spectral_r') ```

(The TROPOMI example data is uploaded to Google Drive)

Method2: bounds

This issue still exists with bounds data: ``` %matplotlib inline import numpy as np import xarray as xr import cartopy.crs as ccrs import matplotlib.pyplot as plt

def prepare_geo(bounds_data): """Prepare lat/lon bounds for pcolormesh. lat/lon bounds are ordered in the following way:: 3----2 | | 0----1 Extend longitudes and latitudes with one element to support "pcolormesh":: (X[i+1, j], Y[i+1, j]) (X[i+1, j+1], Y[i+1, j+1]) +--------+ | C[i,j] | +--------+ (X[i, j], Y[i, j]) (X[i, j+1], Y[i, j+1]) """ # Create the left array left = np.vstack([bounds_data[:, :, 0], bounds_data[-1:, :, 3]]) # Create the right array right = np.vstack([bounds_data[:, -1:, 1], bounds_data[-1:, -1:, 2]]) # Stack horizontally dest = np.hstack([left, right]) # Convert to DataArray dest = xr.DataArray(dest, dims=('y_bounds', 'x_bounds'), attrs=bounds_data.attrs ) return dest

ds = xr.open_dataset('./S5P_OFFL_L2__NO2_20190810T21213620190810T230306_09456_01_010302_20190816T233944.nc', group='PRODUCT').isel(time=0) ds_geo = xr.open_dataset('./S5P_OFFL_L2NO2____20190810T212136_20190810T230306_09456_01_010302_20190816T233944.nc', group='/PRODUCT/SUPPORT_DATA/GEOLOCATIONS').isel(time=0)

lon_bounds = prepare_geo(ds_geo['longitude_bounds']) lat_bounds = prepare_geo(ds_geo['latitude_bounds'])

plt.figure(figsize=(14,6)) ax = plt.axes(projection=ccrs.PlateCarree())

m = ax.pcolormesh(lon_bounds, lat_bounds, ds['nitrogendioxide_tropospheric_column'], # shading='auto', transform=ccrs.PlateCarree(), vmin=0, vmax=1e-4, cmap='Spectral_r') ```

@QuLogic Ha, it looks well with the latest cartopy (0.20.1). Thanks a lot.

@TomNicholas So, is it better to keep this open until the doc is updated?

@TomNicholas I checked the doc and this issue begins from v0.16.1. Note that there're also small spurious lines after v0.10.9. Before v0.10.9, the figure looks fine. It's may be also related to matplotlib ... CC @jklymak and @timhoffm.

BTW, the question on StackOverflow, which was raised by @gerritholl a long time ago, looks similar. I'm not sure whether this is the cartopy issue, CC @QuLogic, and @greglucas.

Sorry for this issue. This actually caused by the missing args like input_core_dims, exclude_dims, etc.

Anyway, this one works well: res = xr.apply_ufunc(reduceat_np, dask_data, bins_reduceat[:5], input_core_dims=[['x'], ['new_x']], exclude_dims=set(('x',)), output_core_dims=[['new_x']], dask="parallelized", output_dtypes=[data.dtype], dask_gufunc_kwargs={'allow_rechunk': True}, ) res.compute()

@dcherian Has this method been improved in dask_groupby? Could you provide a simple example we can follow? I got lost in the dask_groupby documentation ...

@dcherian Thanks! That's simple ;) However, the groupby_bins method is a little different from binned_statistic.

binned_statistic:

All but the last (righthand-most) bin is half-open. In other words, if bins is [1, 2, 3, 4], then the first bin is [1, 2) (including 1, but excluding 2) and the second [2, 3). The last bin, however, is [3, 4], which includes 4.

groupby_bins:

right (bool, default: True) – Indicates whether the bins include the rightmost edge or not. If right == True (the default), then the bins [1,2,3,4] indicate (1,2], (2,3], (3,4].

So, let's check this shorter example: ``` from scipy.stats import binned_statistic import numpy as np import xarray as xr

--- scipy method ---

x = np.arange(10) values = x*5 statistics, _, _ = binned_statistic(x, values, statistic='min', bins=10, range=(0, 10))

--- xarray method ---

x = xr.DataArray(x) values = xr.DataArray(values) bin_res = values.groupby_bins('dim_0', bins=np.linspace(0, 10, 10), right=False, include_lowest=True).min()

print('scipy: \n', statistics) print('xarray: \n', bin_res) ```

Output: ``` scipy: [ 0. 5. 10. 15. 20. 25. 30. 35. 40. 45.]

xarray: <xarray.DataArray (dim_0_bins: 9)> array([ 0, 10, 15, 20, 25, 30, 35, 40, 45]) Coordinates: * dim_0_bins (dim_0_bins) object [0.0, 1.111) ... [8.889, 10.0) ```

The scipy method has one more value ...

Summary

These produce the same results: binned_statistic(x, values, statistic='min', bins=10, range=(0, 10)) values.groupby_bins('dim_0', bins=np.linspace(0, 10, 11), right=False, include_lowest=True).min()

Output: scipy: [ 0. 5. 10. 15. 20. 25. 30. 35. 40. 45.] xarray: <xarray.DataArray (dim_0_bins: 10)> array([ 0, 5, 10, 15, 20, 25, 30, 35, 40, 45]) Coordinates: * dim_0_bins (dim_0_bins) object [0.0, 1.0) [1.0, 2.0) ... [9.0, 10.0)

@markusritschel I tested v0.17.0 and it doesn't work.

data = np.random.rand(4, 3) locs = ["IA", "IL", "IN"] times = pd.date_range("2000-01-01", periods=4) foo = xr.DataArray(data, coords=[times, locs], dims=["time", "space"]) foo.groupby('space').argmax('time')

Error: AttributeError: 'DataArrayGroupBy' object has no attribute 'argmax'

@JavierRuano I find the simpler solution from a similar question in stack overflow.

sort_pair = np.take_along_axis(pair.values, cld.argsort(axis=0), axis=0)

Complete example

``` import xarray as xr import numpy as np

x = 4 y = 2 z = 4

data = np.arange(xyz).reshape(z, y, x)

3d array with coords

cld_1 = xr.DataArray(data, dims=['z', 'y', 'x'], coords={'z': np.arange(z)})

2d array without coords

cld_2 = xr.DataArray(np.arange(xy).reshape(y, x)1.5+1, dims=['y', 'x'])

expand 2d to 3d

cld_2 = cld_2.expand_dims(z=[4])

concat

cld = xr.concat([cld_1, cld_2], dim='z')

paired array

pair = cld.copy(data=np.arange(xy(z+1)).reshape(z+1, y, x))

sort_pair = np.take_along_axis(pair.values, cld.argsort(axis=0), axis=0)

print(cld) print(pair) print(sort_pair) ```

Output: ``` <xarray.DataArray (z: 5, y: 2, x: 4)> array([[[ 0. , 1. , 2. , 3. ], [ 4. , 5. , 6. , 7. ]],

   [[ 8. ,  9. , 10. , 11. ],
    [12. , 13. , 14. , 15. ]],

   [[16. , 17. , 18. , 19. ],
    [20. , 21. , 22. , 23. ]],

   [[24. , 25. , 26. , 27. ],
    [28. , 29. , 30. , 31. ]],

   [[ 1. ,  2.5,  4. ,  5.5],
    [ 7. ,  8.5, 10. , 11.5]]])

Coordinates: * z (z) int64 0 1 2 3 4 Dimensions without coordinates: y, x <xarray.DataArray (z: 5, y: 2, x: 4)> array([[[ 0, 1, 2, 3], [ 4, 5, 6, 7]],

   [[ 8,  9, 10, 11],
    [12, 13, 14, 15]],

   [[16, 17, 18, 19],
    [20, 21, 22, 23]],

   [[24, 25, 26, 27],
    [28, 29, 30, 31]],

   [[32, 33, 34, 35],
    [36, 37, 38, 39]]])

Coordinates: * z (z) int64 0 1 2 3 4 Dimensions without coordinates: y, x [[[ 0 1 2 3] [ 4 5 6 7]]

[[32 33 34 35] [36 37 38 39]]

[[ 8 9 10 11] [12 13 14 15]]

[[16 17 18 19] [20 21 22 23]]

[[24 25 26 27] ```

Note, I have to use pair.values instead of pair in the last sorting step. Otherwise, I will get this error:

IndexError: Unlabeled multi-dimensional array cannot be used for indexing: y

Thanks, it works well.

@keewis Thanks! I will try to apply this method and check the results.

@dcherian Sorry. I made a mistake in the expected result. It should be [[0], [1, 2, 3]]

@JavierRuano The time indexes are same in my real case. Maybe, I have to merge these data if I can't find the solution.

@dcherian Oh, thanks! After scipy is installed, it works: ``` =============================================== test session starts =============================================== platform linux -- Python 3.7.6, pytest-5.4.1, py-1.8.1, pluggy-0.12.0 -- /yin_raid/xin/miniconda3/envs/xarray_dev/bin/python cachedir: .pytest_cache rootdir: /yin_raid/xin/github/xarray, inifile: setup.cfg collected 2 items

test_interp.py::test_nans[True] SKIPPED [ 50%] test_interp.py::test_nans[False] PASSED [100%]

================================================ warnings summary ================================================= xarray/tests/test_interp.py::test_nans[False] xarray/tests/test_interp.py::test_nans[False] /yin_raid/xin/miniconda3/envs/xarray_dev/lib/python3.7/importlib/_bootstrap.py:219: RuntimeWarning: numpy.ufunc size changed, may indicate binary incompatibility. Expected 192 from C header, got 216 from PyObject return f(args, *kwds)

-- Docs: https://docs.pytest.org/en/latest/warnings.html ==================================== 1 passed, 1 skipped, 2 warnings in 0.88s ===================================== ```

I will update the test and pull it soon.

It's my first time to write test for xarray. I tried this method: pytest test_interp.py::test_nans -v, but all tests are skipped: ``` =================================================== test session starts ==================================================== platform linux -- Python 3.7.6, pytest-5.4.1, py-1.8.1, pluggy-0.12.0 -- /yin_raid/xin/miniconda3/envs/xarray_dev/bin/python cachedir: .pytest_cache rootdir: /yin_raid/xin/github/xarray, inifile: setup.cfg collected 2 items

test_interp.py::test_nans[True] SKIPPED [ 50%] test_interp.py::test_nans[False] SKIPPED [100%]

==================================================== 2 skipped in 0.50s ==================================================== ```

How to make the test actually run?

Any update? This issue could result in errors for many functions of xarray, like interp and dot.

I need to use df.index = pd.MultiIndex.from_arrays(.....). See https://github.com/pandas-dev/pandas/issues/33420

@JavierRuano When the dataframe is converted back to dataset, the values aren't changed because of the unchanged Multiindex in dataframe ... I have tried this: df = ds.to_dataframe() new_df = df.sort_values(by=['x', 'y', 'cld']) new_df.index.set_levels(list(np.arange(ds['cld'].sizes['z'])), level='z', inplace=True) But, it doesn't work. Still trying ...

@JavierRuano Nice suggestion! I combine them to dataset, convert it to dataframe and then sort_values. Finally, convert the dataframe back to dataset: ``` ds = cld.to_dataset(name='cld') ds['pair'] = pair

df = ds.to_dataframe() new_ds = df.sort_values(by='cld').to_xarray().transpose() ```

@JavierRuano Thank you very much. This example is a special case. If the order of z is different for each x and y, do we need to create a tmp DataArray to save the result of looping x and y?

@kmuehlbauer Thanks, Nice trick! It works well for this situation.

Thanks, @fujiisoup . @dcherian decided to improve the error message later. So, I will leave this open.

@johnomotani Thanks, it works. Sorry for this simple question ...

@dcherian Sorry for the misunderstanding. I tried again for the 3d array, it works well ;) ``` import xarray as xr import numpy as np

x = 2 y = 4 z = 3 data = np.arange(xyz).reshape(z, x, y)

input array

a = xr.DataArray(data, dims=['z', 'y', 'x'])

start_index array

sindex = xr.DataArray(np.full_like(a[0, ...], 0), dims=['y', 'x'])

end_index array

eindex = xr.DataArray(np.full_like(a[0, ...], 1), dims=['y', 'x'])

zindex = a.z.copy(data=np.arange(a.sizes["z"]))

sub_z = (zindex >= sindex) & (zindex <= eindex) sum_a = a.where(sub_z).sum('z', keepdims=True)

print(a) print(sum_a) ```

``` <xarray.DataArray (z: 3, y: 2, x: 4)> array([[[ 0, 1, 2, 3], [ 4, 5, 6, 7]],

   [[ 8,  9, 10, 11],
    [12, 13, 14, 15]],

   [[16, 17, 18, 19],
    [20, 21, 22, 23]]])

Dimensions without coordinates: z, y, x

<xarray.DataArray (z: 1, y: 2, x: 4)> array([[[ 8., 10., 12., 14.], [16., 18., 20., 22.]]]) Dimensions without coordinates: z, y, x ```

@dcherian Excellent solution!

If we upgrade this to 3d array and sum by z axis, it seems that method isn't suitable:

``` import xarray as xr import numpy as np

x = 2 y = 2 z = 3 data = np.arange(xyz).reshape(z, y, x)

input array

a = xr.DataArray(data, dims=['z', 'y', 'x'])

start_index array

sindex = xr.DataArray(np.full_like(a[0, ...], 0), dims=['y', 'x'])

end_index array

eindex = xr.DataArray(np.full_like(a[0, ...], 1), dims=['y', 'x']) ```

Solution (Boolean)

```

stack indexes

index_list = np.column_stack((sindex, eindex))

all false array

boolean_array = np.zeros(a.shape, dtype=bool)

iterate and assign true

for row in range(len(index_list)): boolean_array[row, np.arange(index_list[row][0], index_list[row][1]+1)] = True

sum_a = a.where(boolean_array).sum(dim='y') ```

@spencerkclark Maybe converting the datetime into number? BTW, Why not forcing numpy >= 1.18.1?

@mathause Thanks! Shall we close this issue?

@mathause For .values, if I delete vectorize=True, I got this error: File "/home/xin/miniconda3/envs/satpy/lib/python3.7/site-packages/scipy/interpolate/interpolate.py", line 455, in __init__ raise ValueError("the x array must have exactly one dimension.") ValueError: the x array must have exactly one dimension. Then, I keep vectorize=True deleted and use the np.interp, I got this error: File "/mnt/d/Github/s5p-wrfchem/s5p_utils.py", line 264, in interp1d_np return np.interp(xi, x, data) File "<__array_function__ internals>", line 6, in interp File "/home/xin/miniconda3/envs/satpy/lib/python3.7/site-packages/numpy/lib/function_base.py", line 1412, in interp return interp_func(x, xp, fp, left, right) ValueError: object too deep for desired array If I let vectorize=True shows again and use the np.interp, I got the error mentioned before: File "/home/xin/miniconda3/envs/satpy/lib/python3.7/site-packages/numpy/lib/function_base.py", line 1830, in _update_dim_sizes % (dim, size, dim_sizes[dim])) ValueError: inconsistent size for core dimension 'dim0': 2 vs 39

For the one without .values, this is the result of repr(s5p['p']): <xarray.DataArray (bottom_top: 25, y: 389, x: 450)> dask.array<where, shape=(25, 389, 450), dtype=float32, chunksize=(25, 389, 450), chunktype=numpy.ndarray> Coordinates: * bottom_top (bottom_top) int32 0 1 2 3 4 5 6 7 8 ... 17 18 19 20 21 22 23 24 vertices int32 0 crs object +proj=latlong +datum=WGS84 +ellps=WGS84 +type=crs Dimensions without coordinates: y, x Attributes: name: p resolution: None calibration: None polarization: None level: None modifiers: () units: hPa After the bottom_up in renamed to new_dim, it works without error for both scipy and numpy interpolation function.

I tested again with a subset of my data: subset_no2 = regrid_vars['no2'].isel(x=277, y=[212, 213]) subset_p = regrid_vars['p'].isel(x=277, y=[212, 213]) subset_interp = s5p['p'].isel(x=277, y=[212, 213]) interped = xr.apply_ufunc( interp1d_np, subset_no2, subset_p, subset_interp, input_core_dims=[["bottom_top"], ["bottom_top"], ["new_dim"]], output_core_dims=[["new_dim"]], exclude_dims=set(("bottom_top",)), vectorize=True, )

Error without .values: File "/home/xin/miniconda3/envs/satpy/lib/python3.7/site-packages/xarray/core/computation.py", line 508, in broadcast_compat_data list(core_dims), missing_core_dims ValueError: operand to apply_ufunc has required core dimensions ['new_dim'], but some of these dimensions are absent on an input variable: ['new_dim']

Error with .values: File "/home/xin/miniconda3/envs/satpy/lib/python3.7/site-packages/numpy/lib/function_base.py", line 1830, in _update_dim_sizes % (dim, size, dim_sizes[dim])) ValueError: inconsistent size for core dimension 'dim0': 2 vs 39

@dcherian If .values is removed, I got this error: File "/home/xin/miniconda3/envs/satpy/lib/python3.7/site-packages/xarray/core/computation.py", line 508, in broadcast_compat_data list(core_dims), missing_core_dims ValueError: operand to apply_ufunc has required core dimensions ['new_dim'], but some of these dimensions are absent on an input variable: ['new_dim']

Here's the information of regrid_vars['no2'], regrid_vars['p'] and s5p['p']: <xarray.DataArray 'no2' (bottom_top: 39, y: 389, x: 450)> <xarray.DataArray (bottom_top: 39, y: 389, x: 450)> <xarray.DataArray (bottom_top: 25, y: 389, x: 450)>

BTW, I have nan values in regrid_vars['no2'] and regrid_vars['p']. I think that wouldn't cause that error.

@mathause Thanks! It works well. Here's the solution:

Code

``` def interp1d_np(data, x, xi): from scipy import interpolate # return np.interp(xi, x, data) f = interpolate.interp1d(x, data, fill_value='extrapolate') return f(xi)

interped = xr.apply_ufunc( interp1d_np, # first the function bottom_up, # now arguments in the order expected by 'interp1_np' pressure.values, # as above interp_p.values, # as above input_core_dims=[["z"], ["z"], ["new_z"]], # list with one entry per arg output_core_dims=[["new_z"]], # returned data has one dimension exclude_dims=set(("z",)), # dimensions allowed to change size. Must be a set! vectorize=True, # loop over non-core dims ) interped = interped.rename({"new_z": "z"})

print(np.testing.assert_allclose(output.values, interped.values)) ```

Result:

None

However, when I apply it to my real data, I got some errors:

Code

``` def interp1d_np(data, x, xi): from scipy import interpolate f = interpolate.interp1d(x, data, fill_value='extrapolate') return f(xi)

interped = xr.apply_ufunc(
    interp1d_np,
    regrid_vars['no2'],
    regrid_vars['p'].values,
    s5p['p'].values,
    input_core_dims=[["bottom_top"], ["bottom_top"], ["new_dim"]],
    output_core_dims=[["new_dim"]],
    exclude_dims=set(("bottom_top",)),
    vectorize=True,
)

```

Error:

File "/home/xin/miniconda3/envs/satpy/lib/python3.7/site-packages/numpy/lib/function_base.py", line 1830, in _update_dim_sizes % (dim, size, dim_sizes[dim])) ValueError: inconsistent size for core dimension 'dim0': 450 vs 39 Here's the output of print(regrid_vars['no2'].shape, regrid_vars['p'].values.shape, s5p['p'].values.shape): (39, 389, 450) (39, 389, 450) (25, 389, 450) The shape looks fine.

Hi @max-sixty, thanks. If this looks well, I'm glad to add the test for it.

@DocOtak Thank you for your explanation! It works well now :)

@shoyer Thanks. It's not datetime64 arrays, this is the result of np.isnat(t): File "/public/software/anaconda/anaconda3/envs/python36/lib/python3.6/site-packages/xarray-0.12.3-py3.6.egg/xarray/core/arithmetic.py", line 69, in __array_ufunc__ dask='allowed') File "/public/software/anaconda/anaconda3/envs/python36/lib/python3.6/site-packages/xarray-0.12.3-py3.6.egg/xarray/core/computation.py", line 969, in apply_ufunc keep_attrs=keep_attrs) File "/public/software/anaconda/anaconda3/envs/python36/lib/python3.6/site-packages/xarray-0.12.3-py3.6.egg/xarray/core/computation.py", line 217, in apply_dataarray_vfunc result_var = func(*data_vars) File "/public/software/anaconda/anaconda3/envs/python36/lib/python3.6/site-packages/xarray-0.12.3-py3.6.egg/xarray/core/computation.py", line 564, in apply_variable_ufunc result_data = func(*input_data) TypeError: ufunc 'isnat' is only defined for datetime and timedelta.

I use pd.isnull(t).all() to check it, it works. Actually it's all nan. There's something wrong with the nc file, I will contact the data center. Thank you for all your help :)

@keewis I tried to using np.isnan(t.values).all() to check whether it's all nan. But, I got this error: print (np.isnan(t.values).all()) TypeError: ufunc 'isnan' not supported for the input types, and the inputs could not be safely coerced to any supported types according to the casting rule ''safe'' This is the type of t.values: <class 'numpy.ndarray'>

@shoyer Thank. It works now. But, I get another question. This is the result of t = ds['time_utc']: <xarray.DataArray 'time_utc' (time: 1, scanline: 357, ground_pixel: 450)> array([[[nan, nan, ..., nan, nan], [nan, nan, ..., nan, nan], ..., [nan, nan, ..., nan, nan], [nan, nan, ..., nan, nan]]], dtype=object) Coordinates: * scanline (scanline) float64 1.0 2.0 3.0 4.0 ... 354.0 355.0 356.0 357.0 * ground_pixel (ground_pixel) float64 1.0 2.0 3.0 4.0 ... 448.0 449.0 450.0 * time (time) datetime64[ns] 2019-08-25 Attributes: long_name: Time of observation as ISO 8601 date-time string If I want to get the minimum value by t.min(skipna=True), I get the strange type: <xarray.DataArray 'time_utc' ()> array(<xarray.core.dtypes.AlwaysGreaterThan object at 0x7f96ac188550>, dtype=object) Can't convert it to string by str(t.min(skipna=True)).

@max-sixty Actually, I'm using numpy = 1.13.1 and I need skipna= True. Don't understand the error it shows.

@dcherian Can't find how to do that by rcParams. Could you give a simple example?

@dcherian Thanks! Figure out now:

``` import xarray as xr import matplotlib.pyplot as plt

airtemps = xr.tutorial.open_dataset('air_temperature') air = airtemps.air - 273.15 air2d = air.isel(time=500)

im = air2d.plot.pcolormesh(add_colorbar=False) cb = plt.colorbar(im, orientation="horizontal", pad=0.15) cb.set_label(label='Temperature ($^{\circ}$C)', size='large', weight='bold') cb.ax.tick_params(labelsize='large') ```

@dcherian It works by netCDF4, but not for xarray: ``` file = Dataset('ds1.nc') print (file.variables['time'],'\n')

with xr.open_dataset('ds1.nc') as f: print (f.time.attrs) Output: <class 'netCDF4._netCDF4.Variable'> int64 time(time) units: hours since 2015-01-01 unlimited dimensions: current shape = (3,) filling on, default _FillValue of -9223372036854775806 used

OrderedDict() ```

What's the difference between units and attrs?

@xylar Thanks! I just found another question similar to this one.

I've tried some operations: with xr.open_dataset('merge.nc') as f: print (f['temperature'],'\n') print ('---------------------------') print (f.mean(dim='time')) print ('---------------------------') print (f['temperature'].loc[:,:,'2015-01-05T04:00:00',]) print ('---------------------------')

It works fine: ``` <xarray.DataArray 'temperature' (x: 2, y: 2, time: 6)> array([[[-0.022611, -1.428088, -0.655508, 0.977389, -0.428088, 0.344492], [ 0.430102, 0.996973, -0.882054, 1.430102, 1.996973, 0.117946]],

   [[ 0.157233, -0.230397, -0.505775,  1.157233,  0.769603,  0.494225],
    [-0.075826, -1.933904, -0.823982,  0.924174, -0.933904,  0.176018]]])

Coordinates: lon (x, y) float64 ... lat (x, y) float64 ... * time (time) datetime64[ns] 2015-01-05T04:00:00 2015-01-05T05:00:00 ... Dimensions without coordinates: x, y

<xarray.Dataset> Dimensions: (x: 2, y: 2) Coordinates: lon (x, y) float64 ... lat (x, y) float64 ... Dimensions without coordinates: x, y Data variables: temperature (x, y) float64 -0.2021 0.6817 0.307 -0.4446

<xarray.DataArray 'temperature' (x: 2, y: 2)> array([[-0.022611, 0.430102], [ 0.157233, -0.075826]]) Coordinates: lon (x, y) float64 ... lat (x, y) float64 ... time datetime64[ns] 2015-01-05T04:00:00 Dimensions without coordinates: x, y

```