Подгонка полинома к частям большого xipray CMIP6 NetCDF? - PullRequest
Новая
       132

Подгонка полинома к частям большого xipray CMIP6 NetCDF?

1 голос
/ 30 апреля 2020

Моя цель - рассчитать квадратичное c соответствие по временному измерению данных «zos» (уровня моря) прогонов piControl различных моделей CMIP6 и вычесть их из других прогонов. В этом конкретном примере данные выглядят следующим образом: 

Out[108]: 
<xarray.Dataset>
Dimensions:    (bnds: 2, ncells: 830305, time: 3, vertices: 16)
Coordinates:
  * time       (time) object 2401-01-16 12:00:00 ... 2401-03-16 12:00:00
    lon        (ncells) float64 dask.array<chunksize=(5000,), meta=np.ndarray>
    lat        (ncells) float64 dask.array<chunksize=(5000,), meta=np.ndarray>
Dimensions without coordinates: bnds, ncells, vertices
Data variables:
    time_bnds  (time, bnds) object dask.array<chunksize=(3, 2), meta=np.ndarray>
    lon_bnds   (ncells, vertices) float64 dask.array<chunksize=(5000, 16), meta=np.ndarray>
    lat_bnds   (ncells, vertices) float64 dask.array<chunksize=(5000, 16), meta=np.ndarray>
    zos        (time, ncells) float32 dask.array<chunksize=(3, 5000), meta=np.ndarray>

Я пытался использовать xarray.Dataset.polyfit () , но это приводит к AttributeError: объект Dataset не имеет атрибута 'polyfit' . Предположительно, это связано с тем, что модели CMIP6 имеют нерегулярные календари и года моделей, которые находятся за пределами поддерживаемого диапазона (например, 'позволяет декодировать временную ось в полные объекты numpy .datetime64, продолжая вместо этого использовать объекты cftime.datetime, причина : даты вне диапазона '). По этой же причине я использую encode_times = False во фрагменте ниже.

В качестве альтернативы я пробовал это решение : 

from numpy.polynomial.polynomial import polyval,polyfit
import pandas as pd
import xarray as xr
import dask.array as da
import numpy as np

xr.set_options(display_style="html")  # fancy HTML repr

mycmip6 = (
    xr.open_dataset('/Volumes/Naamloos/PhD_Data/CMIP6/raw_mergetime/zos/AWI-CM-1-1-MR/zos_Omon_AWI-CM-1-1-MR_piControl_r1i1p1f1_gn_240101-290012.nc',decode_times=False,chunks={'ncells':5000})
    .isel(time=slice(3)) #using a small slice for example
    )
mycmip6

def fit_quadratic(data, time):
    pfit = np.polyfit(time, data,2) 
    
    return np.transpose( np.polyval(pfit,time) )


pfit = xr.apply_ufunc(
    fit_quadratic,  # first the function
    mycmip6.zos, #sea level data
    mycmip6.time,  # time
    input_core_dims=[["time"], ["time"]],  # list with one entry per arg
    output_core_dims=[["time"]],  # returned data has one dimension
    vectorize=True,  # loop over non-core dims
    dask="parallelized",
    output_dtypes=[mycmip6.zos.dtype],
)
Это работает нормально, но как только я пытаюсь манипулировать или сохранять pfit в netcdf, 
pfit.compute()

я сталкиваюсь с этой ошибкой: 

TypeError: wrong number of positional arguments: expected 2, got 3

Когда я вообще не использую chunking, я могу использовать 'pfit', но у меня возникают проблемы с памятью для самых больших наборов данных.

Что я здесь не так делаю?

...