cgeniepy.grid ============= .. py:module:: cgeniepy.grid Classes ------- .. autoapisummary:: cgeniepy.grid.GridOperation cgeniepy.grid.Interpolator Module Contents --------------- .. py:class:: GridOperation A set of operations on grid/coordinate data .. py:method:: get_genie_lon(N=36, edge=False, offset_start=-260) get GENIE longitude in 10 degree resolution, if edge is False, then return midpoint :param N: number of grid points :param edge: if True, return edge points :param offset_start: the par_grid_lon_offset option in main configuration file .. py:method:: get_genie_lat(N=36, edge=False) return cGENIE latitude in log-sine normally degree resolution, if edge is False, then return midpoint .. py:method:: get_genie_depth(N=16, edge=False, max_depth=5000) calculate cGENIE vertical depth :param N: number of grid points :param edge: if True, return edge points, otherwise return midpoints .. py:method:: get_normal_lon(N=36, edge=False) Normal longitude in 10 degree resolution (default), if edge is False, then return midpoint .. py:method:: lon_n2g(x, grid_lon_offset=-260) Convert normal longitude (-180, 180) to GENIE longitude (-270, 90) :param x: normal longitude :return: GENIE longitude .. py:method:: lon_g2n(x) Convert GENIE longitude to normal longitude. This is independent on the grid_offset_start option :param x: GENIE longitude :return: normal longitude .. py:method:: lon_e2n(x) Convert eastern longitude to normal longitude :param x: longitude in eastern degree :return: normal longitude .. py:method:: lon_n2e(x) Convert normal longitude (-180, 180) to longitude east(0,360) :param x: normal longitude :return: longitude in eastern degree .. py:method:: xr_n2g(data: xarray.Dataset, longitude='lon', *args, **kwargs) -> xarray.Dataset Apply longitude conversion method n2g for the input data (normal to GENIE) :param data: input data :param longitude: longitude coordinate name :returns: xr.Dataset .. py:method:: xr_g2n(data: xarray.Dataset, longitude='lon', *args, **kwargs) -> xarray.Dataset Apply longitude conversion method g2n for the input data (GENIE to normal) :param data: input data :param longitude: longitude coordinate name :returns: xr.Dataset .. py:method:: xr_e2n(data: xarray.Dataset, longitude='lon', *args, **kwargs) -> xarray.Dataset Apply longitude conversion method e2n for the input data (eastern to normal) :param data: input data :param longitude: longitude coordinate name :returns: xr.Dataset .. py:method:: xr_n2e(data: xarray.Dataset, longitude='lon', *args, **kwargs) -> xarray.Dataset Apply longitude conversion method n2e for the input data (normal to eastern) :param data: input data :param longitude: longitude coordinate name :returns: xr.Dataset .. py:method:: mask_Arctic_Med(array, policy='na') mask Arctic and Meditterean Sea in cGENIE modern continent configuration :param array: 36x36 GENIE array .. py:method:: GENIE_grid_mask(base='worjh2', basin='ALL', subbasin='', invert=False) Get a modern GENIE 36x36 mask array from input data. The input array is flipped (left/right flip -> up/down flip) for easy recognition :continent: worjh2, worlg4, worbe2, GIteiiaa, GIteiiaa, p0055c :basin: Atlantic/Pacific/Indian/ALL/Tanzania :subbasin: N/S/ALL, ALL means Southern Ocean section included :returns: GENIE grid array where continent/ice cap is 0 and ocean is 1, default is 'worjh2' .. py:method:: geniebin_lat(x, *args, **kwargs) Categorize into cGENIE grid bins .. py:method:: geniebin_lon(x, *args, **kwargs) Categorize into cGENIE grid bins .. py:method:: geniebin_depth(x, *args, **kwargs) Categorize into cGENIE grid bins .. py:method:: normbin_lon(x, *args, **kwargs) Categorize into normal grid bins .. py:method:: haversine_distance(lat1, lon1, lat2, lon2) Calculate the Haversine distance between corresponding pairs of points. :param lat1: Array of latitudes for points 1 :param lon1: Array of longitudes for points 1 :param lat2: Array of latitudes for points 2 :param lon2: Array of longitudes for points 2 :returns: Array of distances between corresponding pairs of points .. py:method:: geo_dis3d(point1, points2) Calculate the 3D geographical distance between point1 and multiple points2. :param point1: Tuple/list of coordinates (z, lat, lon) or (lat, lon) :param points2: Numpy Array of shape (n, 3) containing coordinates (z, lat, lon) of points :returns: Array of distances between point1 and each point in points2 .. py:method:: geo_dis2d(point1, points2) Calculate the 2D geographical distance between point1 and multiple points2. :param point1: Tuple/list of coordinates (z, lon, lat) or (lon, lat) :param points2: Array of shape (n, 3) containing coordinates (z, lon, lat) of points :returns: Array of distances between point1 and each point in points2 .. py:method:: check_dimension(input) :staticmethod: check the presence of latitude, longitude, depth, and time in the input tuple :param input: tuple of dimension names :return: tuple of boolean values indicating the presence of latitude, longitude, depth, and time Example ----------- >>> input = ('lat', 'lon', 'depth', 'time') >>> check_dimension(input) (True, True, True, True) .. py:method:: dim_order(input) Determine the order of dimensions in the input data array :return: tuple of index in the input data Example ----------- >>> input = ('lat', 'lon', 'depth') >>> dim_order(input) ## always follow time, depth, lat, lon (2, 0, 1) .. py:method:: set_coordinates(obj, index) :staticmethod: Set the coordinates of the input object based on the input index :param obj: object to set coordinates :param index: tuple of dimension names .. py:class:: Interpolator(dims, coordinates, values, grid_number=200, method='r-linear') A univeral ineteprolator for cgeniepy that can be used to interpolate data for both regular and irregular grid .. py:attribute:: dims .. py:attribute:: coords .. py:attribute:: values .. py:attribute:: grid_number :value: 200 .. py:attribute:: gridded_coord :value: [] .. py:attribute:: meshgrid .. py:attribute:: interp_function .. py:attribute:: gridded_data .. py:method:: new_coordinate(n) create new coordinates for regridding :param n: the resolution of the new grid :returns: a list of new coordinates .. py:method:: new_meshgrid(*args, **kwargs) create meshgrid for regridding :returns: numpy meshgrid .. py:method:: interpolate_data(*args, **kwargs) Use the interpolation function to get regridded values :returns: regridded data array .. py:method:: to_xarray() Combine the regridded data with the new coordinates to create a xarray DataArray :returns: xarray DataArray .. py:method:: to_dataframe() Convert the regridded data to a pandas DataFrame :returns: pandas DataFrame