<jupyter_start><jupyter_text>Resample data and create climatologies<jupyter_code>%%time
ds_sst.analysed_sst[0,:1799,:3600].plot()
%%time
ds_sst.mask[0,:1799,:3600].plot()
%%time
sst_masked = ds_sst.where(ds_sst.mask==1)
sst_masked.analysed_sst[0,2500:3500,:3600].plot()
masked_25km = sst_masked.interp({'lon':ds_ccmp.longitude,'lat':ds_ccmp.latitude})
cyr1,cyr2='2010-01-01','2019-12-31'
#make day average
ds_ccmp = ds_ccmp.resample(time='D').mean()
ds_ccmp_clim = ds_ccmp.sel(time=slice(cyr1,cyr2))
ds_ccmp_clim = ds_ccmp_clim.groupby('time.dayofyear').mean('time',keep_attrs=True,skipna=False)

#sst 25km
ds_sst_clim = sst_masked.sel(time=slice(cyr1,cyr2))
ds_sst_clim = ds_sst_clim.groupby('time.dayofyear').mean('time',keep_attrs=True,skipna=False)

    
    #put data into a dictionary
    data_dict={'wnd':ds_ccmp,
               'sst':ds_sst}
    clim_dict={'wnd_clim':ds_ccmp_clim,
               'sst_clim':ds_sst_clim}
  
    return data_dict,clim_dict


def multi_apply_along_axis(func1d, axis, arrs, *args, **kwargs):
    """
    function from: https://climate-cms.org/2019/07/29/multi-apply-along-axis.html
    Given a function `func1d(A, B, C, ..., *args, **kwargs)`  that acts on 
    multiple one dimensional arrays, apply that function to the N-dimensional
    arrays listed by `arrs` along axis `axis`
    
    If `arrs` are one dimensional this is equivalent to::
    
        func1d(*arrs, *args, **kwargs)
    
    If there is only one array in `arrs` this is equivalent to::
    
        numpy.apply_along_axis(func1d, axis, arrs[0], *args, **kwargs)
        
    All arrays in `arrs` must have compatible dimensions to be able to run
    `numpy.concatenate(arrs, axis)`
    
    Arguments:
        func1d:   Function that operates on `len(arrs)` 1 dimensional arrays,
                  with signature `f(*arrs, *args, **kwargs)`
        axis:     Axis of all `arrs` to apply the function along
        arrs:     Iterable of numpy arrays
        *args:    Passed to func1d after array arguments
        **kwargs: Passed to func1d as keyword arguments
    """
    import numpy
    # Concatenate the input arrays along the calculation axis to make one big
    # array that can be passed in to `apply_along_axis`
    carrs = numpy.concatenate(arrs, axis)
    
    # We'll need to split the concatenated arrays up before we apply `func1d`,
    # here's the offsets to split them back into the originals
    offsets=[]
    start=0
    for i in range(len(arrs)-1):
        start += arrs[i].shape[axis]
        offsets.append(start)
            
    # The helper closure splits up the concatenated array back into the components of `arrs`
    # and then runs `func1d` on them
    def helperfunc(a, *args, **kwargs):
        arrs = numpy.split(a, offsets)
        return func1d(*[*arrs, *args], **kwargs)
    
    # Run `apply_along_axis` along the concatenated array
    return numpy.apply_along_axis(helperfunc, axis, carrs, *args, **kwargs)

data,clim=get_data()
data['sst']
data['wnd']
ds_sst, ds_wnd = data['sst'],data['wnd']
ds_wnd['wspd']=(ds_wnd.uwnd**2+ds_wnd.vwnd**2)**.5
#interp doesn't work on chunked dims so rechunk
ds_sst = ds_sst.chunk({'time':1,'lat':901,'lon':1800})
ds_new = ds_sst.interp(lat = ds_wnd.lat,lon=ds_wnd.lon)
#SST
#remove mean, seasonal cycle, trend before analysis
clim = ds_new.groupby('time.dayofyear').mean('time',keep_attrs=True,skipna=False)
ds_new_tem = ds_new.groupby('time.dayofyear')-clim
#detrending doesn't work with nan so fill with nan but will have to mask at end and pixels with nan in timeseries
tem = ds_new_tem.fillna(0)
ds_detrended_sst = signal.detrend(tem.analysed_sst,axis=0)
ds_new_tem.analysed_sst.isel(lon=720,lat=310).plot()
plt.plot(ds_new_tem.time,ds_detrended_sst[:,310,720])
#WSPD
#remove mean, seasonal cycle, trend before analysis
clim = ds_wnd.groupby('time.dayofyear').mean('time',keep_attrs=True,skipna=False)
ds_new_tem = ds_wnd.groupby('time.dayofyear')-clim
#detrending doesn't work with nan so fill with nan but will have to mask at end and pixels with nan in timeseries
tem = ds_new_tem.fillna(0)
ds_detrended_wnd = signal.detrend(tem.wspd,axis=0)
ds_new_tem.wspd.isel(lon=720,lat=310).plot()
plt.plot(ds_new_tem.time,ds_detrended_wnd[:,310,720])
corr = multi_apply_along_axis(pearsonr, 0, [ds_detrended_sst,ds_detrended_wnd])
print(corr.shape)
print(type(corr))
#not sure why by 1440 is blank
acorr = corr[0,:,:]
acorr2 = np.concatenate([acorr[:,720:1439],acorr[:,:1439],acorr[:,:720]],axis=1)
#plt.imshow(acorr[:,:1439],vmin=-.1,vmax=.1,cmap='RdBu')
#plt.imshow(acorr2,vmin=-.1,vmax=.1,cmap='viridis')
data = acorr2
data = np.nan_to_num(data,0)
lowpass = ndimage.gaussian_filter(data, 40)
gauss_highpass = data - lowpass





plt.imshow(lowpass,vmin=-.1,vmax=.1,cmap='RdBu')
plt.imshow(gauss_highpass,vmin=-.1,vmax=.1,cmap='RdBu')<jupyter_output><empty_output>
StarcoderdataJupyterStructuredCleanDedup