Methods
Generic methods
- nes.methods.cell_measures.calculate_cell_area(grid_corner_lon, grid_corner_lat, earth_radius_minor_axis=6356752.3142, earth_radius_major_axis=6378137.0)[source]
Calculate the area of each cell of a grid.
Parameters
- grid_corner_lonnp.array
Array with longitude bounds of grid.
- grid_corner_latnp.array
Array with longitude bounds of grid.
- earth_radius_minor_axisfloat
Radius of the minor axis of the Earth.
- earth_radius_major_axisfloat
Radius of the major axis of the Earth.
- nes.methods.cell_measures.calculate_geometry_area(geometry_list, earth_radius_minor_axis=6356752.3142, earth_radius_major_axis=6378137.0)[source]
Get coordinate bounds and call function to calculate the area of each cell of a set of geometries.
Parameters
- geometry_listList
List with polygon geometries.
- earth_radius_minor_axisfloat
Radius of the minor axis of the Earth.
- earth_radius_major_axisfloat
Radius of the major axis of the Earth.
- nes.methods.cell_measures.calculate_grid_area(self)[source]
Get coordinate bounds and call function to calculate the area of each cell of a grid.
Parameters
- selfnes.Nes
Source projection Nes Object.
- nes.methods.cell_measures.cross_product(a, b)[source]
Calculate cross product between two points.
Parameters
- anp.array
Position of point a in cartesian coordinates.
- bnp.array
Position of point b in cartesian coordinates.
- nes.methods.cell_measures.lon_lat_to_cartesian(lon, lat, earth_radius_major_axis=6378137.0)[source]
Calculate lon, lat coordinates of a point on a sphere.
Parameters
- lonnp.array
Longitude values.
- latnp.array
Latitude values.
- earth_radius_major_axisfloat
Radius of the major axis of the Earth.
- nes.methods.cell_measures.mod_huiliers_area(cell_corner_lon, cell_corner_lat)[source]
Calculate the area of each cell according to Huilier’s theorem. Reference: CDO (https://earth.bsc.es/gitlab/ces/cdo/).
Parameters
- cell_corner_lonnp.array
Longitude boundaries of each cell.
- cell_corner_latnp.array
Latitude boundaries of each cell.
- nes.methods.cell_measures.norm(cp)[source]
Normalize the result of the cross product operation.
Parameters
- cpnp.array
Cross product between two points.
- nes.methods.cell_measures.tri_area(point_0, point_1, point_2)[source]
Calculate area between three points that form a triangle. Reference: CDO (https://earth.bsc.es/gitlab/ces/cdo/).
Parameters
- point_0np.array
Position of first point in cartesian coordinates.
- point_1np.array
Position of second point in cartesian coordinates.
- point_2np.array
Position of third point in cartesian coordinates.
Horizontal interpolation
- nes.methods.horizontal_interpolation.create_area_conservative_weight_matrix(self, dst_nes, wm_path=None, flux=False, info=False)[source]
To create the weight matrix with the area conservative method.
Parameters
- selfnes.Nes
Source projection Nes Object.
- dst_nesnes.Nes
Final projection Nes object.
- wm_pathstr
Path where write the weight matrix.
- fluxbool
Indicates if you want to calculate the weight matrix for flux variables.
- info: bool
Indicates if you want to print extra info during the methods process.
Returns
- nes.Nes
Weight matrix.
- nes.methods.horizontal_interpolation.create_nn_weight_matrix(self, dst_grid, n_neighbours=4, wm_path=None, info=False)[source]
To create the weight matrix with the nearest neighbours method.
Parameters
- selfnes.Nes
Source projection Nes Object.
- dst_gridnes.Nes
Final projection Nes object.
- n_neighboursint
Used if kind == NearestNeighbour. Number of nearest neighbours to interpolate. Default: 4.
- wm_pathstr
Path where write the weight matrix.
- info: bool
Indicates if you want to print extra info during the methods process.
Returns
- nes.Nes
Weight matrix.
- nes.methods.horizontal_interpolation.get_src_data(comm, var_data, idx, parallel_method)[source]
To obtain the needed src data to interpolate.
Parameters
- commMPI.Communicator.
MPI communicator.
- var_datanp.array
Rank source data.
- idxnp.array
Index of the needed data in a 2D flatten way.
- parallel_method: str
Source parallel method.
Returns
- np.array
Flatten source needed data.
- nes.methods.horizontal_interpolation.get_weights_idx_t_axis(self, dst_grid, weight_matrix_path, kind, n_neighbours, only_create, wm, flux)[source]
To obtain the weights and source data index through the T axis.
Parameters
- selfnes.Nes
Source projection Nes Object.
- dst_gridnes.Nes
Final projection Nes object.
- weight_matrix_pathstr, None
Path to the weight matrix to read/create.
- kindstr
Kind of horizontal interpolation. Accepted values: [‘NearestNeighbour’, ‘Conservative’].
- n_neighboursint
Used if kind == NearestNeighbour. Number of nearest neighbours to interpolate. Default: 4.
- only_createbool
Indicates if you want to only create the Weight Matrix.
- wmNes
Weight matrix Nes File.
- fluxbool
Indicates if you want to calculate the weight matrix for flux variables.
Returns
- tuple
Weights and source data index.
- nes.methods.horizontal_interpolation.get_weights_idx_xy_axis(self, dst_grid, weight_matrix_path, kind, n_neighbours, only_create, wm, flux)[source]
To obtain the weights and source data index through the X or Y axis.
Parameters
- selfnes.Nes
Source projection Nes Object.
- dst_gridnes.Nes
Final projection Nes object.
- weight_matrix_pathstr, None
Path to the weight matrix to read/create.
- kindstr
Kind of horizontal interpolation. Accepted values: [‘NearestNeighbour’, ‘Conservative’].
- n_neighboursint
Used if kind == NearestNeighbour. Number of nearest neighbours to interpolate. Default: 4.
- only_createbool
Indicates if you want to only create the Weight Matrix.
- wmNes
Weight matrix Nes File.
- fluxbool
Indicates if you want to calculate the weight matrix for flux variables.
Returns
- tuple
Weights and source data index.
- nes.methods.horizontal_interpolation.interpolate_horizontal(self, dst_grid, weight_matrix_path=None, kind='NearestNeighbour', n_neighbours=4, info=False, to_providentia=False, only_create_wm=False, wm=None, flux=False)[source]
Horizontal methods from one grid to another one.
Parameters
- selfnes.Nes
Source projection Nes Object.
- dst_gridnes.Nes
Final projection Nes object.
- weight_matrix_pathstr, None
Path to the weight matrix to read/create.
- kindstr
Kind of horizontal interpolation. Accepted values: [‘NearestNeighbour’, ‘Conservative’].
- n_neighboursint
Used if kind == NearestNeighbour. Number of nearest neighbours to interpolate. Default: 4.
- infobool
Indicates if you want to print extra info during the methods process.
- to_providentiabool
Indicates if we want the interpolated grid in Providentia format.
- only_create_wmbool
Indicates if you want to only create the Weight Matrix.
- wmNes
Weight matrix Nes File.
- fluxbool
Indicates if you want to calculate the weight matrix for flux variables.
- nes.methods.horizontal_interpolation.lon_lat_to_cartesian(lon, lat, radius=6378137.0)[source]
Calculate lon, lat coordinates of a point on a sphere.
DEPRECATED!!!!
Parameters
- lonnp.array
Longitude values.
- latnp.array
Latitude values.
- radiusfloat
Radius of the sphere to get the distances.
- nes.methods.horizontal_interpolation.lon_lat_to_cartesian_ecef(lon, lat)[source]
Convert observational/model geographic longitude/latitude coordinates to cartesian ECEF (Earth Centred, Earth Fixed) coordinates, assuming WGS84 datum and ellipsoid, and that all heights = 0. ECEF coordiantes represent positions (in meters) as X, Y, Z coordinates, approximating the earth surface as an ellipsoid of revolution. This conversion is for the subsequent calculation of euclidean distances of the model gridcell centres from each observational station. Defining the distance between two points on the earth’s surface as simply the euclidean distance between the two lat/lon pairs could lead to inaccurate results depending on the distance between two points (i.e. 1 deg. of longitude varies with latitude).
Parameters
- lonnp.array
Longitude values.
- latnp.array
Latitude values.
- nes.methods.horizontal_interpolation.read_weight_matrix(weight_matrix_path, comm=None, parallel_method='T')[source]
Read weight matrix.
Parameters
- weight_matrix_pathstr
Path of the weight matrix.
- commMPI.Communicator
Communicator to read the weight matrix.
- parallel_methodstr
Nes parallel method to read the weight matrix.
Returns
- nes.Nes
Weight matrix.
Spatial join
- nes.methods.spatial_join.get_bbox(self)[source]
Obtain the bounding box of the rank data (lon_min, lat_min, lon_max, lat_max).
Parameters
- selfnes.Nes
Nes Object.
Returns
- tuple
Bounding box
- nes.methods.spatial_join.prepare_external_shapefile(self, ext_shp, var_list, info=False, apply_bbox=True)[source]
Prepare the external shapefile.
It is high recommended to pass ext_shp parameter as string because it will clip the external shapefile to the rank.
Read if it is not already read
Filter variables list
Standardize projections
Parameters
- selfnes.Nes
Nes Object.
- ext_shpGeoDataFrame or str
External shapefile or path to it.
- var_listList[str] or None
External shapefile variables to be computed.
- infobool
Indicates if you want to print the information.
- apply_bboxbool
Indicates if you want to reduce the shapefile to a bbox.
Returns
- GeoDataFrame
External shapefile.
- nes.methods.spatial_join.spatial_join(self, ext_shp, method=None, var_list=None, info=False, apply_bbox=True)[source]
Compute overlay intersection of two GeoPandasDataFrames.
Parameters
- selfnes.Nes
Nes Object.
- ext_shpGeoPandasDataFrame or str
File or path from where the data will be obtained on the intersection.
- methodstr
Overlay method. Accepted values: [‘nearest’, ‘intersection’, ‘centroid’].
- var_listList or None or str
Variables that will be included in the resulting shapefile.
- infobool
Indicates if you want to print the process info.
- apply_bboxbool
Indicates if you want to reduce the shapefile to a bbox.
- nes.methods.spatial_join.spatial_join_centroid(self, ext_shp, info=False)[source]
Perform the spatial join using the centroid method.
Parameters
- selfnes.Nes
Nes Object.
- ext_shpGeoDataFrame
External shapefile.
- infobool
Indicates if you want to print the information.
Vertical interpolation
- nes.methods.vertical_interpolation.add_4d_vertical_info(self, info_to_add)[source]
To add the vertical information from other source.
Parameters
- selfnes.Nes
Source Nes object.
- info_to_addnes.Nes, str
Nes object with the vertical information as variable or str with the path to the NetCDF file that contains the vertical data.
- nes.methods.vertical_interpolation.interpolate_vertical(self, new_levels, new_src_vertical=None, kind='linear', extrapolate=None, info=None, overwrite=False)[source]
Vertical interpolation.
Parameters
- selfNes
Source Nes object.
- new_levelsList
List of new vertical levels.
- new_src_verticalnes.Nes, str
Nes object with the vertical information as variable or str with the path to the NetCDF file that contains the vertical data.
- kindstr
Vertical methods type.
- extrapolateNone, tuple, str
Extrapolate method (for non linear operations).
- info: None, bool
Indicates if you want to print extra information.
- overwrite: bool
Indicates if you want to compute the vertical interpolation in the same object or not.