@@ -0,0 +1,2863 @@
+from flyingpigeon import utils%0Afrom flyingpigeon.ocgis_module import call%0A%0A%0Adef get_anomalies(nc_file, frac=0.2, reference=None):%0A    %22%22%22%0A    Anomalisation of data subsets for weather classification by subtracting a smoothed annual cycle%0A%0A    :param nc_file: input netCDF file%0A    :param frac: Number between 0-1 for strength of smoothing%0A               (0 = close to the original data, 1 = flat line)%0A               default = 0.2%0A    :param reference: Period to calculate annual cycle%0A%0A    :returns str: path to output netCDF file%0A    %22%22%22%0A    try:%0A        variable = utils.get_variable(nc_file)%0A        calc = %5B%7B'func': 'mean'%7D%5D  # ,'name': variable%0A        calc_grouping = %5B'day', 'month'%5D%0A        nc_anual_cycle = call(nc_file,%0A                              calc=calc,%0A                              calc_grouping=calc_grouping,%0A                              time_range=reference)%0A        print('annual cycle calculated')%0A    except Exception as e:%0A        msg = 'failed to calcualte annual cycle %25s' %25 e%0A        print(msg)%0A        # raise Exception(msg)%0A    #%0A    # try:%0A    #     # spline for smoothing%0A    #     import statsmodels.api as sm%0A    #     from numpy import tile, empty, linspace%0A    #     from netCDF4 import Dataset%0A    #     from cdo import Cdo%0A    #     cdo = Cdo()%0A    #     # variable = utils.get_variable(nc_file)%0A    #     ds = Dataset(nc_anual_cycle, mode='a')%0A    #     vals = ds.variables%5Bvariable%5D%0A    #     vals_sm = empty(vals.shape)%0A    #     ts = vals.shape%5B0%5D%0A    #     x = linspace(1, ts*3, num=ts*3, endpoint=True)%0A    #     for lat in range(vals.shape%5B1%5D):%0A    #         for lon in range(vals.shape%5B2%5D):%0A    #             try:%0A    #                 y = tile(vals%5B:, lat, lon%5D, 3)%0A    #                 # ys = smooth(y, window_size=91, order=2, deriv=0, rate=1)%5Bts:ts*2%5D%0A    #                 ys = sm.nonparametric.lowess(y, x, frac=frac)%5Bts:ts*2, 1%5D%0A    #                 vals_sm%5B:, lat, lon%5D = ys%0A    #             except:%0A    #                 msg = 'failed for lat %25s lon %25s' %25 (lat, lon)%0A    #                 print(msg)%0A    #                 raise Exception(msg)%0A    #         print('done for %25s - %25s ' %25 (lat, lon))%0A    #     vals%5B:, :, :%5D = vals_sm%5B:, :, :%5D%0A    #     ds.close()%0A    #     print('smothing of annual cycle done')%0A    # except:%0A    #     msg = 'failed smothing of annual cycle'%0A    #     print(msg)%0A    #     raise Exception(msg)%0A    # try:%0A    #     ip, nc_anomal = mkstemp(dir='.', suffix='.nc')%0A    #     nc_anomal = cdo.sub(input=%5Bnc_file, nc_anual_cycle%5D, output=nc_anomal)%0A    #     print('cdo.sub; anomalisation done: %25s ' %25 nc_anomal)%0A    # except:%0A    #     msg = 'failed substraction of annual cycle'%0A    #     print(msg)%0A    #     raise Exception(msg)%0A    return nc_anual_cycle%0A%0A%0Anc_file = '/home/nils/data/409f9cdc-67b0-11e7-9260-9cb6d0d3acd7.nc'%0A%0A%0Anc_an = get_anomalies(nc_file)%0A%0A%0Aprint nc_an%0A

add get_anomalies script