pygeogrids package

Submodules

pygeogrids.geodetic_datum module

class pygeogrids.geodetic_datum.GeodeticDatum(ellps, **kwargs)[source]

Bases: object

Class representing a geodetic datum providing transformation and calculation functionality

Parameters

ellString (string) – String of geodetic datum (ellipsoid) as provided in pyproj

EllM(lat)[source]

Method to calculate the radius of the curvature

Parameters

lat (numpy.array, list or float) – Geodatic latitudes of the points in the grid

Returns

r – radius of the curvature

Return type

np.array, float

EllN(lat)[source]

Method to calculate the radius of the prime vertical

Parameters

lat (numpy.array, list or float) – Geodatic latitudes of the points in the grid in degrees

Returns

r – radius of the prime vertical

Return type

np.array, float

GaussianRadi(lat)[source]

Method to calculate the gaussian radius of the curvature

Parameters

lat (numpy.array, list or float) – Geodatic latitudes of the points in the grid

Returns

r – gaussian radius of the curvature

Return type

np.array, float

GeocentricDistance(lon, lat)[source]

Method to calculate the geocentric distance to given points

Parameters
  • lon (numpy.array, list or float) – Geodatic longitude of the points in the grid

  • lat (numpy.array, list or float) – Geodatic latitudes of the points in the grid

Returns

r – Geocentric radius

Return type

np.array, float

GeocentricLat(lat)[source]

Method to calculate the geocentric from the geodatic latitude.

Parameters

lat (numpy.array, list or float) – Geodatic latitudes of the points in the grid

Returns

lat_geocentric – Geocentric latitude

Return type

np.array, float

GeodeticLat(lat)[source]

Method to calculate the geodatic from the geocentric latitude.

Parameters

lat (numpy.array, list or float) – Geocentric latitudes of the points in the grid

Returns

lat_geodatic – Geodatic latitude

Return type

np.array, float

MeridianArcDist(lat1, lat2)[source]

Method to calculate the distance between two parallels (meridian arc distance)

Parameters
  • lat1 (numpy.array, float) – Geodatic latitudes 1

  • lat2 (numpy.array, float) – Geodatic latitudes 2

Returns

dist – Meridian arc distance

Return type

np.array, float

ParallelArcDist(lat, lon1, lon2)[source]

Method to calculate the distance between two points on a given parallel

Parameters
  • lat (float) – Geodatic latitudes of the points in the grid

  • lon1 (float) – Longitude of point 1 at the given parallel

  • lon2 (float) – Longitude of point 2 at the given parallel

Returns

dist – Parallel arc distance

Return type

np.array, float

ParallelRadi(lat)[source]

Method to get the radius the parallel at a given latitude.

Parameters

lat (numpy.array, list or float) – latitudes of the points in the grid

Returns

radius – Radius of parallel

Return type

np.array, float

ReducedLat(lat)[source]

Method to calculate the reduced from the geodatic latitude.

Parameters

lat (numpy.array, list or float) – Geodatic latitudes of the points in the grid

Returns

lat_reduced – Reduced latitude

Return type

np.array, float

getParameter()[source]

Method to transform lon/lat to ECEF (Earth-Centered, Earth-Fixed) coordinates representing a 3d Cartesian coordinate system.

Parameters
  • lon (numpy.array, list or float) – longitudes of the points in the grid

  • lat (numpy.array, list or float) – latitudes of the points in the grid

Returns

x, y, z – 3D cartesian coordinates

Return type

np.array

toECEF(lon, lat)[source]

Method to transform lon/lat to ECEF (Earth-Centered, Earth-Fixed) coordinates representing a 3d Cartesian coordinate system.

Parameters
  • lon (numpy.array, list or float) – longitudes of the points in the grid

  • lat (numpy.array, list or float) – geodatic latitudes of the points in the grid

Returns

x, y, z – 3D cartesian coordinates

Return type

np.array

pygeogrids.grids module

The grids module defines the grid classes.

class pygeogrids.grids.BasicGrid(lon, lat, gpis=None, geodatum='WGS84', subset=None, setup_kdTree=True, shape=None, transform_lon=None, kd_tree_name='pykdtree')[source]

Bases: object

Grid that just has lat,lon coordinates and can find the nearest neighbour. It can also yield the gpi, lat, lon information in order.

Parameters
  • lon (numpy.array) – longitudes of the points in the grid

  • lat (numpy.array) – latitudes of the points in the grid

  • geodatum (basestring) – Name of the geodatic datum associated with the grid

  • gpis (numpy.array, optional) – if the gpi numbers are in a different order than the lon and lat arrays an array containing the gpi numbers can be given if no array is given here the lon lat arrays are given gpi numbers starting at 0

  • subset (numpy.array, optional) – if the active part of the array is only a subset of all the points then the subset array which is a index into lon and lat can be given here.

  • setup_kdTree (boolean, optional) – if set (default) then the kdTree for nearest neighbour search will be built on initialization

  • shape (tuple, optional) – The shape of the grid array in 2-d space. e.g. for a 1x1 degree global regular grid the shape would be (180,360). if given the grid can be reshaped into the given shape this indicates that it is a regular grid and fills the attributes self.lon2d and self.lat2d which define the grid only be the meridian coordinates(self.lon2d) and the coordinates of the circles of latitude(self.lat2d). The shape has to be given as (lat2d, lon2d) It it is not given the shape is set to the length of the input lon and lat arrays.

  • transform_lon (bool or None, optional (default: None)) – Whether to transform longitudes to values between -180 and 180. By default values are transformed, but a warning is issued. To turn off the warning, set this to True, to turn of transformation set this to False.

  • kd_tree_name (str, optional (default: 'pykdtree')) – KDTree engine, ‘pykdtree’ or ‘scipy’

arrlon

1D array of all longitudes of the grid

Type

numpy.array

arrlat

1D array of all latitudes of the grid

Type

numpy.array

n_gpi

number of gpis in the grid

Type

int

gpidirect

if true the gpi number is equal to the index of arrlon and arrlat

Type

boolean

gpis

gpi number for elements in arrlon and arrlat gpi[i] is located at arrlon[i],arrlat[i]

Type

numpy.array

subset

if given then this contains the indices of a subset of the grid. This can be used if only a part of a grid is interesting for a application. e.g. land points, or only a specific country

Type

numpy.array

allpoints

if False only a subset of the grid is active

Type

boolean

activearrlon

array of longitudes that are active, is defined by arrlon[subset] if a subset is given otherwise equal to arrlon

Type

numpy.array

activearrlat

array of latitudes that are active, is defined by arrlat[subset] if a subset is given otherwise equal to arrlat

Type

numpy.array

activegpis

array of gpis that are active, is defined by gpis[subset] if a subset is given otherwise equal to gpis

Type

numpy.array

geodatum

pygeogrids.geodatic_datum object (reference ellipsoid) associated with the grid

Type

object

issplit

if True then the array was split in n parts with the self.split function

Type

boolean

kdTree

grid.nearest_neighbor.findGeoNN object for nearest neighbor search

Type

object

shape

if given during initialization then this is the shape the grid can be reshaped to

Type

tuple, optional

lat2d

if shape is given this attribute contains all latitudes according to the provided 2d-shape that make up the grid

Type

numpy.array, optional

lon2d

if shape is given this attribute contains all longitudes according to the provided 2d-shape that make up the grid

Type

numpy.array, optional

calc_lut(other, max_dist=inf, into_subset=False)[source]

Takes other BasicGrid or CellGrid objects and computes a lookup table between them. The lut will have the size of self.n_gpis and will for every grid point have the nearest index into other.arrlon etc.

Parameters
  • other (grid object) – to which to calculate the lut to

  • max_dist (float, optional) – maximum allowed distance in meters

  • into_subset (boolean, optional) – if set the returned lut will have the index into the subset if the other grid is a subset of a grid. Example: if e.g. ind_l is used for the warp_grid some datasets will be given as arrays with len(ind_l) elements. These datasets can not be indexed with gpi numbers but have to be indexed with indices into the subset

find_k_nearest_gpi(lon, lat, max_dist=inf, k=1)[source]

Find k nearest gpi, builds kdTree if it does not yet exist.

Parameters
  • lon (float or iterable) – Longitude of point.

  • lat (float or iterable) – Latitude of point.

  • max_dist (float, optional) – Maximum distance to consider for search (default: np.Inf).

  • k (int, optional) – The number of nearest neighbors to return (default: 1).

Returns

  • gpi (np.ndarray) – Grid point indices.

  • distance (np.ndarray) – Distance of gpi(s) to given lon, lat. At the moment not on a great circle but in spherical cartesian coordinates.

find_nearest_gpi(lon, lat, max_dist=inf)[source]

Finds nearest gpi, builds kdTree if it does not yet exist.

Parameters
  • lon (float or iterable) – Longitude of point.

  • lat (float or iterable) – Latitude of point.

  • max_dist (float, optional) – Maximum distance [m] to consider for search (default: np.Inf).

Returns

  • gpi (long) – Grid point index. If no point was found within the maximum distance to consider, an empty array is returned.

  • distance (float) – Distance of gpi to given lon, lat. At the moment not on a great circle but in spherical cartesian coordinates. If no point was found within the maximum distance to consider, an empty array is returned.

get_bbox_grid_points(latmin=- 90, latmax=90, lonmin=- 180, lonmax=180, coords=False, both=False)[source]

Returns all grid points located in a submitted geographic box, optional as coordinates.

Parameters
  • latmin (float, optional (default: -90)) – minimum latitude

  • latmax (float, optional (default: 90)) – maximum latitude

  • lonmin (float, optional (default: -180)) – minimum latitude

  • lonmax (float, optional (default: 180)) – maximum latitude

  • coords (boolean, optional (default: False)) – set to True if coordinates should be returned

  • both (boolean, optional (default: False)) – set to True if gpis and coordinates should be returned

Returns

  • gpi (numpy.ndarray) – grid point indices, if coords=False

  • lat (numpy.ndarray) – longitudes of gpis, if coords=True

  • lon (numpy.ndarray) – longitudes of gpis, if coords=True

get_grid_points(*args)[source]

Returns all active grid points.

Parameters

n (int, optional) – if the grid is split in n parts using the split function then this function will only return the nth part of the grid

Returns

  • gpis (numpy.ndarray) – Grid point indices.

  • arrlon (numpy.ndarray) – Longitudes.

  • arrlat (numpy.ndarray) – Latitudes.

get_shp_grid_points(ply)[source]

Returns all grid points located in a submitted shapefile, optinal as coordinates. Currently only works in WGS84.

Parameters

ply (object, OGRGeometryShadow) – the Geometry of the Feature as returned from ogr.GetGeometryRef

Returns

grid – Subgrid.

Return type

BasicGrid

gpi2lonlat(gpi)[source]

Longitude and latitude for given gpi.

Parameters

gpi (int32 or iterable) – Grid point index.

Returns

  • lon (float) – Longitude of gpi.

  • lat (float) – Latitude of gpi

gpi2rowcol(gpi)[source]

If the grid can be reshaped into a sensible 2D shape then this function gives the row(latitude dimension) and column(longitude dimension) indices of the gpi in the 2D grid.

Parameters

gpi (int32) – Grid point index.

Returns

  • row (int) – Row in 2D array.

  • col (int) – Column in 2D array.

grid_points(*args)[source]

Yields all grid points in order

Parameters

n (int, optional) – if the grid is split in n parts using the split function then this iterator will only iterate of the nth part of the grid

Returns

  • gpi (long) – grid point index

  • lon (float) – longitude of gpi

  • lat (float) – longitude of gpi

static issorted(a)[source]
split(n)[source]

Function splits the grid into n parts this changes not function but grid_points() which takes the argument n and will only iterate through this part of the grid.

Parameters

n (int) – Number of parts the grid should be split into

subgrid_from_gpis(gpis)[source]

Generate a subgrid for given gpis.

Parameters

gpis (int, numpy.ndarray) – Grid point indices.

Returns

grid – Subgrid.

Return type

BasicGrid

to_cell_grid(cellsize=5.0, cellsize_lat=None, cellsize_lon=None)[source]

Convert grid to cellgrid with a cell partition of cellsize.

Parameters
  • cellsize (float, optional) – Cell size in degrees

  • cellsize_lon (float, optional) – Cell size in degrees on the longitude axis

  • cellsize_lat (float, optional) – Cell size in degrees on the latitude axis

Returns

cell_grid – Cell grid object.

Return type

CellGrid object

unite()[source]

Unites a split array, so that it can be iterated over as a whole again.

class pygeogrids.grids.CellGrid(lon, lat, cells, gpis=None, geodatum='WGS84', subset=None, setup_kdTree=False, **kwargs)[source]

Bases: pygeogrids.grids.BasicGrid

Grid that has lat,lon coordinates as well as cell informatin. It can find nearest neighbour. It can also yield the gpi, lat, lon, cell information in cell order. This is important if the data on the grid is saved in cell files on disk as we can go through all grid points with optimized IO performance.

Parameters
  • lon (numpy.ndarray) – Longitudes of the points in the grid.

  • lat (numpy.ndarray) – Latitudes of the points in the grid.

  • cells (numpy.ndarray) – Of same shape as lon and lat, containing the cell number of each gpi.

  • gpis (numpy.ndarray, optional) – If the gpi numbers are in a different order than the lon and lat arrays an array containing the gpi numbers can be given.

  • subset (numpy.array, optional) – If the active part of the array is only a subset of all the points then the subset array which is a index into lon, lat and cells can be given here.

arrcell

Array of cell number with same shape as arrlon, arrlat.

Type

numpy.ndarray

activearrcell

Array of longitudes that are active, is defined by arrlon[subset] if a subset is given otherwise equal to arrlon.

Type

numpy.ndarray

get_bbox_grid_points(latmin=- 90, latmax=90, lonmin=- 180, lonmax=180, coords=False, both=False, sort_by_cells=True)[source]

Returns all grid points located in a submitted geographic box, optional as coordinates. The points are grouped by cells!

Parameters
  • latmin (float, optional (default: -90)) – minimum latitude

  • latmax (float, optional (default: 90)) – maximum latitude

  • lonmin (float, optional (default: -180)) – minimum latitude

  • lonmax (float, optional (default: 180)) – maximum latitude

  • coords (boolean, optional (default: False)) – set to True if coordinates should be returned

  • both (boolean, optional (default: False)) – set to True if gpis and coordinates should be returned

  • sort_by_cells (bool, optional (default: True)) – Points in the bbox are grouped by cells (so that e.g. when iterating over them you can read data from a cell that is already in memory).

Returns

  • gpi (numpy.ndarray) – grid point indices, if coords=False

  • lat (numpy.ndarray) – longitudes of gpis, if coords=True

  • lon (numpy.ndarray) – longitudes of gpis, if coords=True

get_cells()[source]

Function to get all cell numbers of the grid.

Returns

cells – Unique cell numbers.

Return type

numpy.ndarray

get_grid_points(*args)[source]

Returns all active grid points.

Parameters

n (int, optional) – If the grid is split in n parts using the split function then this function will only return the nth part of the grid.

Returns

  • gpis (numpy.ndarray) – Grid point indices.

  • arrlon (numpy.ndarray) – Longitudes.

  • arrlat (numpy.ndarray) – Latitudes.

  • cells (numpy.ndarray) – Cell numbers.

gpi2cell(gpi)[source]

Cell for given gpi.

Parameters

gpi (int32 or iterable) – Grid point index.

Returns

cell – Cell number of gpi.

Return type

int or iterable

grid_points_for_cell(cells)[source]

Get all grid points for a given cell number.

Parameters

cell (int, numpy.ndarray) – Cell numbers.

Returns

  • gpis (numpy.ndarray) – Gpis belonging to cell.

  • lons (numpy.array) – Longitudes belonging to the gpis.

  • lats (numpy.array) – Latitudes belonging to the gpis.

split(n)[source]

Function splits the grid into n parts this changes not function but grid_points() which takes the argument n and will only iterate through this part of the grid.

Parameters

n (int) – Number of parts the grid should be split into.

subgrid_from_cells(cells)[source]

Generate a subgrid for given cells.

Parameters

cells (int, numpy.ndarray) – Cell numbers.

Returns

grid – Subgrid.

Return type

CellGrid

subgrid_from_gpis(gpis)[source]

Generate a subgrid for given gpis.

Parameters

gpis (int, numpy.ndarray) – Grid point indices.

Returns

grid – Subgrid.

Return type

BasicGrid

exception pygeogrids.grids.GridDefinitionError[source]

Bases: Exception

exception pygeogrids.grids.GridIterationError[source]

Bases: Exception

pygeogrids.grids.genreg_grid(grd_spc_lat=1, grd_spc_lon=1, minlat=- 90.0, maxlat=90.0, minlon=- 180.0, maxlon=180.0, **kwargs)[source]

Define a global regular lon lat grid which starts in the North Western Corner of minlon, maxlat. The grid points are defined to be in the middle of a grid cell. e.g. the first point on a 1x1 degree grid with minlon -180.0 and maxlat 90.0 will be at -179.5 longitude, 89.5 latitude.

Parameters
  • grd_spc_lat (float, optional) – Grid spacing in latitude direction.

  • grd_spc_lon (float, optional) – Grid spacing in longitude direction.

  • minlat (float, optional) – Minimum latitude of the grid.

  • maxlat (float, optional) – Maximum latitude of the grid.

  • minlon (float, optional) – Minimum longitude of the grid.

  • maxlon (float, optional) – Maximum longitude of the grid.

pygeogrids.grids.gridfromdims(londim, latdim, origin='top', **kwargs)[source]

Defines new grid object from latitude and longitude dimensions. Latitude and longitude dimensions are 1D arrays that give the latitude and longitude values of a 2D latitude-longitude array.

Parameters
  • londim (numpy.ndarray) – longitude dimension

  • latdim (numpy.ndarray) – latitude dimension

  • origin (Literal['bottom', 'top'], optional (default: 'top')) – If bottom is selected, the GPI origin is at (min_lon, min_lat), i.e. in the bottom left corner. If ‘top’ is selected, the origin is at (min_lon, max_lat), i.e. in the top left corner

Returns

grid – New grid object.

Return type

BasicGrid

pygeogrids.grids.lonlat2cell(lon, lat, cellsize=5.0, cellsize_lon=None, cellsize_lat=None)[source]

Partition lon, lat points into cells.

Parameters
  • lat (float64, or numpy.ndarray) – Latitude.

  • lon (float64, or numpy.ndarray) – Longitude.

  • cellsize (float) – Cell size in degrees.

  • cellsize_lon (float, optional) – Cell size in degrees on the longitude axis.

  • cellsize_lat (float, optional) – Cell size in degrees on the latitude axis.

Returns

cell – Cell numbers.

Return type

int32, or numpy.ndarray

pygeogrids.grids.reorder_to_cellsize(grid, cellsize_lat, cellsize_lon)[source]

Reorder grid points in one grid to follow the ordering of differently sized cells. This is useful if e.g. a 10x10 degree CellGrid should be traversed in an order compatible with a 5x5 degree CellGrid.

Parameters
Returns

new_grid – output grid with original cell sizes but different ordering.

Return type

pygeogrids.grids.CellGrid

pygeogrids.nearest_neighbor module

class pygeogrids.nearest_neighbor.findGeoNN(lon, lat, geodatum, grid=False, kd_tree_name='pykdtree')[source]

Bases: object

class that takes lat,lon coordinates, transformes them to cartesian (X,Y,Z) coordinates and provides a interface to scipy.spatial.kdTree as well as pykdtree if installed

Parameters
  • lon (numpy.array or list) – longitudes of the points in the grid

  • lat (numpy.array or list) – latitudes of the points in the grid

  • geodatum (object) – pygeogrids.geodatic_datum.GeodeticDatum object associated with lons/lats coordinates

  • grid (boolean, optional) – if True then lon and lat are assumed to be the coordinates of a grid and will be used in numpy.meshgrid to get coordinates for all grid points

  • kd_tree_name (string, optional) – name of kdTree implementation to use, either ‘pykdtree’ to use pykdtree or ‘scipy’ to use scipy.spatial.kdTree Fallback is always scipy if any other string is given or if pykdtree is not installed. standard is pykdtree since it is faster

geodatum

pygeogrids.geodatic_datum.GeodeticDatum object used for x,y,z coordinates calculations

Type

object

coords

3D array of cartesian x,y,z coordinates

Type

numpy.array

kd_tree_name

name of kdTree implementation to use, either ‘pykdtree’ to use pykdtree or ‘scipy’ to use scipy.spatial.kdTree Fallback is always scipy if any other string is given or if pykdtree is not installed

Type

string

kdtree

kdTree object that is built only once and saved in this attribute

Type

object

find_nearest_index(lon, lat)[source]

finds the nearest neighbor of the given lon,lat coordinates in the lon,lat arrays given during initialization and returns the index of the nearest neighbour in those arrays.

find_nearest_index(lon, lat, max_dist=inf, k=1)[source]

finds nearest index, builds kdTree if it does not yet exist

Parameters
  • lon (float, list or numpy.array) – longitude of point

  • lat (float, list or numpy.array) – latitude of point

  • max_dist (float, optional) – Maximum distance to consider for search (default: np.Inf).

  • k (int, optional) – The number of nearest neighbors to return (default: 1).

Returns

  • d (float, numpy.array) – distances of query coordinates to the nearest grid point, distance is given in cartesian coordinates and is not the great circle distance at the moment. This should be OK for most applications that look for the nearest neighbor which should not be hundreds of kilometers away. If no point was found within the maximum distance to consider, an empty array is returned.

  • ind (int, numpy.array) – If self.grid is False indices of nearest neighbor. If no point was found within the maximum distance to consider, an empty array is returned.

  • index_lon (numpy.array, optional) – If self.grid is True then return index into lon array of grid definition. If no point was found within the maximum distance to consider, an empty array is returned.

  • index_lat (numpy.array, optional) – If self.grid is True then return index into lat array of grid definition. If no point was found within the maximum distance to consider, an empty array is returned.

pygeogrids.netcdf module

Module for saving grid to netCDF.

pygeogrids.netcdf.load_grid(filename, subset_flag='subset_flag', subset_value=1, location_var_name='gpi', **grid_kwargs)[source]

load a grid from netCDF file

Parameters
  • filename (string) – filename

  • subset_flag (string, optional (default: 'subset_flag')) – name of the subset to load.

  • subset_value (int or list, optional (default: 1)) – Value(s) of the subset variable that points are loaded for.

  • location_var_name (string, optional (default: 'gpi')) – variable name under which the grid point locations are stored

  • **grid_kwargs (additional kwargs that are passed to BasicGrid or CellGrid) –

Returns

grid – grid instance initialized with the loaded data

Return type

BasicGrid or CellGrid instance

pygeogrids.netcdf.save_grid(filename, grid, subset_name='subset_flag', subset_value=1.0, subset_meaning='water land', global_attrs=None)[source]

save a BasicGrid or CellGrid to netCDF it is assumed that a subset should be used as land_points

Parameters
  • filename (string) – name of file

  • grid (BasicGrid or CellGrid object) – grid whose definition to save to netCDF

  • subset_name (string, optional (default: 'subset_flag')) – long_name of the netcdf variable if the subset symbolises something other than a land/sea mask

  • subset_value (float, optional (default: 1.)) – Value that that the subgrid is written down as in the file.

  • subset_meaning (string, optional (default: 'water land')) – will be written into flag_meanings metadata of variable ‘subset_name’

  • global_attrs (dict, optional) – if given will be written as global attributes into netCDF file

pygeogrids.netcdf.save_lonlat(filename, arrlon, arrlat, geodatum, arrcell=None, gpis=None, subsets={}, global_attrs=None, format='NETCDF4', zlib=False, complevel=4, shuffle=True)[source]

saves grid information to netCDF file

Parameters
  • filename (string) – name of file

  • arrlon (numpy.array) – array of longitudes

  • arrlat (numpy.array) – array of latitudes

  • geodatum (object) – pygeogrids.geodetic_datum.GeodeticDatum object associated with lon/lat

  • arrcell (numpy.array, optional) – array of cell numbers

  • gpis (numpy.array, optional) – gpi numbers if not index of arrlon, arrlat

  • subsets (dict of dicts, optional) –

    keys : long_name of the netcdf variables values : dict with the following keys: points, meaning e.g. subsets = {‘subset_flag’: {‘points’: numpy.array,

    ’value’: int, ‘meaning’: ‘water, land’}}

  • global_attrs (dict, optional) – if given will be written as global attributs into netCDF file

  • format (string, optional) –

    choose either from one of these NetCDF formats

    ’NETCDF4’ ‘NETCDF4_CLASSIC’ ‘NETCDF3_CLASSIC’ ‘NETCDF3_64BIT_OFFSET’

  • zlib (boolean, optional) – see netCDF documentation

  • shuffle (boolean, optional) – see netCDF documentation

  • complevel (int, opational) – see netCDF documentation

pygeogrids.netcdf.sort_for_netcdf(lons, lats, values)[source]

Sort an 2D array for storage in a netCDF file. This mans that the latitudes are stored from 90 to -90 and the longitudes from -180 to 180. Input arrays have to have shape latdim, londim which would mean for a global 10 degree grid (18, 36).

Parameters
Returns

  • lons (numpy.ndarray) – 2D numpy array of longitudes, sorted

  • lats (numpy.ndarray) – 2D numpy array of latitudes, sorted

  • values (numpy.ndarray) – 2D numpy array of values to sort, sorted

pygeogrids.plotting module

pygeogrids.plotting.plot_cell_grid_partitioning(output, cellsize_lon=5.0, cellsize_lat=5.0, figsize=(12, 6))[source]

Plot an overview of a global cell partitioning.

Parameters

output (string) – output file name

pygeogrids.shapefile module

Module for extracting grid points from global administrative areas

pygeogrids.shapefile.get_gad_grid_points(grid, gadm_shp_path, level, name=None, oid=None)[source]

Returns all grid points located in a administrative area. For this function the files from http://biogeo.ucdavis.edu/data/gadm2.8/gadm28_levels.shp.zip need to be available in the folder gadm_shp_path Optinal as coordinates. Currently only works in WGS84.

Parameters
  • grid (object) –

  • gadm_shp_path (path) – Location to GADM28 shapefiles

  • level (int) – Global Administrative Database Level 0 : country 1 : province/county/state/region/municipality/… 2 : municipality/District/county/…

  • name (str) – name of region at indicated level. For countries the english name

  • oid (int) – OBJECTID of feature. This only works with the correct level shp.

Returns

grid – Subgrid.

Return type

BasicGrid

Raises
  • ValueError – If name or oid are not found in shapefile of given level:

  • ImportError – If gdal or osgeo are not installed:

Module contents