ll-geo / antbgm Goto Github PK

def fetch_gravity(
    type: str, plot: bool = False, info: bool = False, region=None, spacing=10e3
) -> xr.DataArray:
    """
    Loads an Antarctic gravity grid
    Preliminary compilation of Antarctica gravity and gravity gradient data.
    Updates on 2016 AntGG compilation.
    Accessed from https://ftp.space.dtu.dk/pub/RF/4D-ANTARCTICA/.
    Parameters
    ----------
    type : str
        either 'FA' or 'BA', for free-air and bouguer anomalies, respectively.
    plot : bool, optional
        choose to plot grid, by default False
    info : bool, optional
        choose to print info on grid, by default False
    region : str or np.ndarray, optional
        GMT-format region to clip the loaded grid to, by default doesn't clip
    spacing : str or int, optional
        grid spacing to resample the loaded grid to, by default 10e3
    Returns
    -------
    xr.DataArray
        Returns a loaded, and optional clip/resampled grid of either free-air or
        Bouguer gravity anomalies.
    """

    if region is None:
        region = (-3330000, 3330000, -3330000, 3330000)
    path = pooch.retrieve(
        url="https://ftp.space.dtu.dk/pub/RF/4D-ANTARCTICA/ant4d_gravity.zip",
        known_hash=None,
        processor=pooch.Unzip(),
        progressbar=True,
    )
    file = [p for p in path if p.endswith(".dat")][0]
    df = pd.read_csv(
        file,
        delim_whitespace=True,
        skiprows=3,
        names=["id", "lat", "lon", "FA", "Err", "DG", "BA"],
    )
    transformer = Transformer.from_crs("epsg:4326", "epsg:3031")
    df["x"], df["y"] = transformer.transform(df.lat.tolist(), df.lon.tolist())
    df = pygmt.blockmedian(
        df[["x", "y", type]], spacing=spacing, region=region, verbose="q"
    )
    grd = pygmt.surface(
        data=df[["x", "y", type]], spacing=spacing, region=region, M="2c", verbose="q"
    )
    if plot is True:
        grd.plot(robust=True)
    if info is True:
        print(pygmt.grdinfo(grd))
    return grd