How to use the salem.gis.Grid function in salem

To help you get started, we’ve selected a few salem examples, based on popular ways it is used in public projects.

Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.

github fmaussion / salem / salem / sio.py View on Github external
nx = len(ds.dimensions['west_east'])
        ny = len(ds.dimensions['south_north'])
    except AttributeError:
        # maybe an xarray dataset
        nx = ds.dims['west_east']
        ny = ds.dims['south_north']
    if hasattr(ds, 'PROJ_ENVI_STRING'):
        # HAR
        x0 = ds['west_east'][0]
        y0 = ds['south_north'][0]
    else:
        # Normal WRF file
        e, n = gis.transform_proj(wgs84, proj, cen_lon, cen_lat)
        x0 = -(nx-1) / 2. * dx + e  # DL corner
        y0 = -(ny-1) / 2. * dy + n  # DL corner
    grid = gis.Grid(nxny=(nx, ny), x0y0=(x0, y0), dxdy=(dx, dy), proj=proj)

    if tmp_check_wrf:
        #  Temporary asserts
        if 'XLONG' in ds.variables:
            # Normal WRF
            reflon = ds.variables['XLONG']
            reflat = ds.variables['XLAT']
        elif 'XLONG_M' in ds.variables:
            # geo_em
            reflon = ds.variables['XLONG_M']
            reflat = ds.variables['XLAT_M']
        elif 'lon' in ds.variables:
            # HAR
            reflon = ds.variables['lon']
            reflat = ds.variables['lat']
        else:
github fmaussion / salem / salem / gis.py View on Github external
operations to do with the same grid, set ``lut`` to an existing
            table obtained from a previous call to  :py:meth:`Grid.grid_lookup`
        return_lut : bool, optional
            set to True if you want to return the lookup table for later use.
            in this case, returns a tuple

        Returns
        -------
        An aggregated ndarray of the data, in ``self`` coordinates.
        If ``return_lut==True``, also return the lookup table
        """

        # Input checks
        if grid is None:
            grid = check_crs(data)  # xarray
        if not isinstance(grid, Grid):
            raise ValueError('grid should be a Grid instance')
        if hasattr(data, 'values'):
            data = data.values  # xarray

        # dimensional check
        in_shape = data.shape
        ndims = len(in_shape)
        if (ndims < 2) or (ndims > 4):
            raise ValueError('data dimension not accepted')
        if (in_shape[-1] != grid.nx) or (in_shape[-2] != grid.ny):
            raise ValueError('data dimension not compatible')

        if lut is None:
            lut = self.grid_lookup(grid)

        # Prepare the output
github fmaussion / salem / salem / wrftools.py View on Github external
pwrf = pwrf.format(**pargs)
    else:
        raise NotImplementedError('WRF proj not implemented yet: '
                                  '{}'.format(map_proj))
    pwrf = gis.check_crs(pwrf)

    # get easting and northings from dom center (probably unnecessary here)
    e, n = gis.transform_proj(wgs84, pwrf, pargs['ref_lon'], pargs['lat_0'])

    # LL corner
    nx, ny = e_we[0]-1, e_sn[0]-1
    x0 = -(nx-1) / 2. * dx + e  # -2 because of staggered grid
    y0 = -(ny-1) / 2. * dy + n

    # parent grid
    grid = gis.Grid(nxny=(nx, ny), x0y0=(x0, y0), dxdy=(dx, dy), proj=pwrf)

    # child grids
    out = [grid]
    for ips, jps, pid, ratio, we, sn in zip(i_parent_start, j_parent_start,
                                            parent_id, parent_ratio,
                                            e_we, e_sn):
        if ips == 1:
            continue
        ips -= 1
        jps -= 1
        we -= 1
        sn -= 1
        nx = we / ratio
        ny = sn / ratio
        if nx != (we / ratio):
            raise RuntimeError('e_we and ratios are incompatible: '
github fmaussion / salem / salem / gis.py View on Github external
else:
        ny = nx * yy / xx

    nx = np.rint(nx)
    ny = np.rint(ny)

    e, n = transform_proj(wgs84, projloc, lon, lat)

    if origin== 'upper-left':
        corner = (-xx / 2. + e, yy / 2. + n)
        dxdy = (xx / nx, - yy / ny)
    else:
        corner = (-xx / 2. + e, -yy / 2. + n)
        dxdy = (xx / nx, yy / ny)

    return Grid(proj=projloc, x0y0=corner, nxny=(nx, ny), dxdy=dxdy,
                pixel_ref='corner')
github fmaussion / salem / salem / sio.py View on Github external
lon = ds.variables[x][:]
    lat = ds.variables[y][:]

    # double check for dubious variables
    if not utils.str_in_list([x], utils.valid_names['lon_var']) or \
            not utils.str_in_list([y], utils.valid_names['lat_var']):
        # name not usual. see if at least the range follows some conv
        if (np.max(np.abs(lon)) > 360.1) or (np.max(np.abs(lat)) > 90.1):
            return None

    # Make the grid
    dx = lon[1]-lon[0]
    dy = lat[1]-lat[0]
    args = dict(nxny=(lon.shape[0], lat.shape[0]), proj=wgs84, dxdy=(dx, dy),
                x0y0=(lon[0], lat[0]))
    return gis.Grid(**args)
github fmaussion / salem / salem / wrftools.py View on Github external
we -= 1
        sn -= 1
        nx = we / ratio
        ny = sn / ratio
        if nx != (we / ratio):
            raise RuntimeError('e_we and ratios are incompatible: '
                               '(e_we - 1) / ratio must be integer!')
        if ny != (sn / ratio):
            raise RuntimeError('e_sn and ratios are incompatible: '
                               '(e_sn - 1) / ratio must be integer!')

        prevgrid = out[pid - 1]
        xx, yy = prevgrid.corner_grid.x_coord, prevgrid.corner_grid.y_coord
        dx = prevgrid.dx / ratio
        dy = prevgrid.dy / ratio
        grid = gis.Grid(nxny=(we, sn),
                        x0y0=(xx[ips], yy[jps]),
                        dxdy=(dx, dy),
                        pixel_ref='corner',
                        proj=pwrf)
        out.append(grid.center_grid)

    maps = None
    if do_maps:
        from salem import Map
        import shapely.geometry as shpg

        if map_kwargs is None:
            map_kwargs = {}

        maps = []
        for i, g in enumerate(out):
github fmaussion / salem / salem / grids.py View on Github external
nx = len(nc.dimensions['west_east'])
        ny = len(nc.dimensions['south_north'])
    except AttributeError:
        # maybe an xarray dataset
        nx = nc.dims['west_east']
        ny = nc.dims['south_north']
    if hasattr(nc, 'PROJ_ENVI_STRING'):
        # HAR
        x0 = nc.GRID_X00
        y0 = nc.GRID_Y00
    else:
        # Normal WRF file
        e, n = gis.transform_proj(wgs84, proj, cen_lon, cen_lat)
        x0 = -(nx-1) / 2. * dx + e  # DL corner
        y0 = -(ny-1) / 2. * dy + n  # DL corner
    grid = gis.Grid(nxny=(nx, ny), ll_corner=(x0, y0), dxdy=(dx, dy),
                    proj=proj)

    if tmp_check_wrf:
        #  Temporary asserts
        if 'XLONG' in nc.variables:
            # Normal WRF
            mylon, mylat = grid.ll_coordinates
            reflon = nc.variables['XLONG']
            reflat = nc.variables['XLAT']
            if len(reflon.shape) == 3:
                reflon = reflon[0, :, :]
                reflat = reflat[0, :, :]
            assert np.allclose(reflon, mylon, atol=1e-4)
            assert np.allclose(reflat, mylat, atol=1e-4)
        if 'lon' in nc.variables:
            # HAR