How to use the salem.Grid function in salem

try:
                        data = imresize(data.filled(np.NaN),
                                        (self.grid.ny, self.grid.nx),
                                        order=interp, mode='edge',
                                        anti_aliasing=True)
                    except RuntimeError:
                        # For some order anti_aliasing doesnt work with 'edge'
                        data = imresize(data.filled(np.NaN),
                                        (self.grid.ny, self.grid.nx),
                                        order=interp, mode='edge',
                                        anti_aliasing=False)

            return data

        crs = gis.check_crs(crs, raise_on_error=True)
        if isinstance(crs, Grid):
            # Remap
            if overplot:
                data = self.grid.map_gridded_data(data, crs, interp=interp,
                                                  out=self.data)
            else:
                data = self.grid.map_gridded_data(data, crs, interp=interp)
        else:
            raise ValueError('crs should be a grid, not a proj')

        return data

entity.RGIId = 'RGI50-00.00000'
    entity.O1Region = '00'
    entity.O2Region = '0'
    gdir = GlacierDirectory(entity, base_dir=base_dir, reset=reset)
    gdir.write_shapefile(gpd.GeoDataFrame([entity]), 'outlines')

    # Idealized flowline
    coords = np.arange(0, len(surface_h) - 0.5, 1)
    line = shpg.LineString(np.vstack([coords, coords * 0.]).T)
    fl = Centerline(line, dx=flowline_dx, surface_h=surface_h)
    fl.widths = widths_m / map_dx
    fl.is_rectangular = np.ones(fl.nx).astype(np.bool)
    gdir.write_pickle([fl], 'inversion_flowlines')

    # Idealized map
    grid = salem.Grid(nxny=(1, 1), dxdy=(map_dx, map_dx), x0y0=(0, 0))
    grid.to_json(gdir.get_filepath('glacier_grid'))

    return gdir

ny = nx * yy / xx
        if ny is not None:
            nx = ny * xx / yy
    nx = np.rint(nx)
    ny = np.rint(ny)

    e, n = pyproj.transform(wgs84, projloc, lon, lat)

    if order=='ul':
        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, corner=corner, nxny=(nx, ny), dxdy=dxdy,
                      pixel_ref='corner')

Notes
        -----
        To obtain the exact domain specified in `x` and `y` you may have to
        play with the `size_x` and `size_y` kwargs.
        """

        global API_KEY

        if 'zoom' in kwargs or 'center_ll' in kwargs:
            raise ValueError('incompatible kwargs.')

        # Transform to lonlat
        crs = gis.check_crs(crs)
        if isinstance(crs, pyproj.Proj):
            lon, lat = gis.transform_proj(crs, wgs84, x, y)
        elif isinstance(crs, Grid):
            lon, lat = crs.ij_to_crs(x, y, crs=wgs84)
        else:
            raise NotImplementedError()

        # surely not the smartest way to do but should be enough for now
        mc = (np.mean(lon), np.mean(lat))
        zoom = 20
        while zoom >= 0:
            grid = gis.googlestatic_mercator_grid(center_ll=mc, nx=size_x,
                                                  ny=size_y, zoom=zoom,
                                                  scale=scale)
            dx, dy = grid.transform(lon, lat, maskout=True)
            if np.any(dx.mask):
                zoom -= 1
            else:
                break

# For tidewater glaciers we force border to 10
    if gdir.is_tidewater and cfg.PARAMS['clip_tidewater_border']:
        border = 10

    # Corners, incl. a buffer of N pix
    ulx = np.min(xx) - border * dx
    lrx = np.max(xx) + border * dx
    uly = np.max(yy) + border * dx
    lry = np.min(yy) - border * dx
    # n pixels
    nx = np.int((lrx - ulx) / dx)
    ny = np.int((uly - lry) / dx)

    # Back to lon, lat for DEM download/preparation
    tmp_grid = salem.Grid(proj=proj_out, nxny=(nx, ny), x0y0=(ulx, uly),
                          dxdy=(dx, -dx), pixel_ref='corner')
    minlon, maxlon, minlat, maxlat = tmp_grid.extent_in_crs(crs=salem.wgs84)

    # Open DEM
    source = entity.DEM_SOURCE if hasattr(entity, 'DEM_SOURCE') else None
    if not is_dem_source_available(source,
                                   (minlon, maxlon),
                                   (minlat, maxlat)):
        raise InvalidDEMError('Source: {} not available for glacier {}'
                              .format(source, gdir.rgi_id))
    dem_list, dem_source = get_topo_file((minlon, maxlon), (minlat, maxlat),
                                         rgi_region=gdir.rgi_region,
                                         rgi_subregion=gdir.rgi_subregion,
                                         dx_meter=dx,
                                         source=source)
    log.debug('(%s) DEM source: %s', gdir.rgi_id, dem_source)

raise ValueError(file + ' does not seem to be an ITMIX file.')
        zone = int(bname[pok+3:])
        south = False
        if zone < 0:
            south = True
            zone = -zone
        proj = pyproj.Proj(proj='utm', zone=zone, ellps='WGS84',
                           south=south)

        # brutally efficient
        with rasterio.Env():
            with rasterio.open(file) as src:
                nxny = (src.width, src.height)
                ul_corner = (src.bounds.left, src.bounds.top)
                dxdy = (src.res[0], -src.res[1])
                grid = Grid(x0y0=ul_corner, nxny=nxny, dxdy=dxdy,
                            pixel_ref='corner', proj=proj)
        # done
        self.file = file
        GeoDataset.__init__(self, grid)

# Make a local proj
    lon, lat = center_ll
    proj_params = dict(proj='merc', datum='WGS84')
    projloc = pyproj.Proj(proj_params)

    # Meter per pixel
    mpix = (2 * np.pi * google_earth_radius) / google_pix / (2**zoom)
    xx = nx * mpix
    yy = ny * mpix

    e, n = pyproj.transform(wgs84, projloc, lon, lat)
    corner = (-xx / 2. + e, yy / 2. + n)
    dxdy = (xx / nx, - yy / ny)

    return Grid(proj=projloc, corner=corner, nxny=(nx, ny), dxdy=dxdy,
                pixel_ref='corner')

interp : str, default 'nearest'
            'nearest', 'linear', or 'spline'
        natural_earth : str
           'lr', 'mr' or 'hr' (low res, medium or high res)
           natural earth background img
        """

        if natural_earth is not None:
            from matplotlib.image import imread
            with warnings.catch_warnings():
                # DecompressionBombWarning
                warnings.simplefilter("ignore")
                img = imread(utils.get_natural_earth_file(natural_earth))
            ny, nx = img.shape[0], img.shape[1]
            dx, dy = 360. / nx, 180. / ny
            grid = Grid(nxny=(nx, ny), dxdy=(dx, -dy), x0y0=(-180., 90.),
                        pixel_ref='corner').center_grid
            return self.set_rgb(img, grid, interp='linear')

        if (len(img.shape) != 3) or (img.shape[-1] not in [3, 4]):
            raise ValueError('img should be of shape (y, x, 3) or (y, x, 4)')

        # Unefficient but by far easiest right now
        out = []
        for i in range(img.shape[-1]):
            out.append(self._check_data(img[..., i], crs=crs, interp=interp))
        self._rgb = np.dstack(out)

file = os.path.join(os.path.abspath(os.path.dirname(__file__)),
                        'data', 'demo_glaciers.csv')
    DATA['demo_glaciers'] = pd.read_csv(file, index_col=0)

    # Add other things
    if 'dem_grids' not in DATA:
        grids = {}
        for grid_json in ['gimpdem_90m_v01.1.json',
                          'arcticdem_mosaic_100m_v3.0.json',
                          'AntarcticDEM_wgs84.json',
                          'REMA_100m_dem.json']:
            if grid_json not in grids:
                fp = os.path.join(os.path.abspath(os.path.dirname(__file__)),
                                  'data', grid_json)
                try:
                    grids[grid_json] = salem.Grid.from_json(fp)
                except NameError:
                    pass
        DATA['dem_grids'] = grids