{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# GIS interoperability and export (shapefile, raster)\n", "This notebook shows how to export gridmarthe features to common GIS format.\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Script API" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [], "source": [ "# import modules\n", "import numpy as np\n", "import gridmarthe as gm" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Shapefile, step 1: conversion\n", "\n", "``gridmarthe`` use the ``geopandas`` API to convert grids to a geodataframe object." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "/mnt/d/home/gridmarthe/src/gridmarthe/operasem/gis.py:65: FutureWarning: The return type of `Dataset.dims` will be changed to return a set of dimension names in future, in order to be more consistent with `DataArray.dims`. To access a mapping from dimension names to lengths, please use `Dataset.sizes`.\n", " if 'time' in ds.dims.keys():\n" ] }, { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
chargexydxdytimegeometry
zone
19999.0596.752568.750.50.52012-07-01POLYGON ((596.5 2568.5, 597 2568.5, 597 2569, ...
29999.0597.252568.750.50.52012-07-01POLYGON ((597 2568.5, 597.5 2568.5, 597.5 2569...
39999.0597.752568.750.50.52012-07-01POLYGON ((597.5 2568.5, 598 2568.5, 598 2569, ...
49999.0598.252568.750.50.52012-07-01POLYGON ((598 2568.5, 598.5 2568.5, 598.5 2569...
59999.0598.752568.750.50.52012-07-01POLYGON ((598.5 2568.5, 599 2568.5, 599 2569, ...
\n", "
" ], "text/plain": [ " charge x y dx dy time \\\n", "zone \n", "1 9999.0 596.75 2568.75 0.5 0.5 2012-07-01 \n", "2 9999.0 597.25 2568.75 0.5 0.5 2012-07-01 \n", "3 9999.0 597.75 2568.75 0.5 0.5 2012-07-01 \n", "4 9999.0 598.25 2568.75 0.5 0.5 2012-07-01 \n", "5 9999.0 598.75 2568.75 0.5 0.5 2012-07-01 \n", "\n", " geometry \n", "zone \n", "1 POLYGON ((596.5 2568.5, 597 2568.5, 597 2569, ... \n", "2 POLYGON ((597 2568.5, 597.5 2568.5, 597.5 2569... \n", "3 POLYGON ((597.5 2568.5, 598 2568.5, 598 2569, ... \n", "4 POLYGON ((598 2568.5, 598.5 2568.5, 598.5 2569... \n", "5 POLYGON ((598.5 2568.5, 599 2568.5, 599 2569, ... " ] }, "execution_count": 16, "metadata": {}, "output_type": "execute_result" } ], "source": [ "path_to_data = './data'\n", "ds = gm.load_marthe_grid(f'{path_to_data}/chasim_hallue.out', 'CHARGE', fpastp=f'{path_to_data}/hallue.pastp').isel(time=-1)\n", "gdf = gm.to_geodataframe(ds)\n", "gdf.head()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Shapefile, step 2: export" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# gdf.to_file('chasim.shp')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Export to raster\n", "\n", "For one time step, user can export a regular grid as a TIFF file with function ``glm.to_raster``\n", "For all time steps at once, use ``ds.mart.to_raster(filename_tpl='chasim_hallue')``, the filename template\n", "argument must not contain any extension.\n", "\n", "**CAUTION**: this is only for regular grids. For irregular/nested grids,\n", "the QGIS plugin PostMARTHE is preferable to save chasim as .tiff\n", "It will save the results in raster file with the smallest resolution\n", "available in all nested grids." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Help on function to_raster in module gridmarthe.operasem.gis:\n", "\n", "to_raster(da, x_dim='x', y_dim='y', epsg=27572, fout='raster.tiff')\n", " Write a xr.DataArray to a raster file\n", " \n", " need xarray with rioxarray installed.\n", " Only for regular grids.\n", " TODO: add support for irregular grids, using PostMARTHE QGIS plugin code.\n", " \n", " Parameters\n", " ----------\n", " da : xarray.DataArray\n", " The DataArray to write to a raster file (only 2D, i.e. no time dimension, select\n", " a variable and timestep before).\n", " x_dim : str, optional\n", " The name of the x dimension, by default 'x'.\n", " y_dim : str, optional\n", " The name of the y dimension, by default 'y'.\n", " epsg : int, optional\n", " The EPSG code for the coordinate reference system, by default 27572.\n", " fout : str, optional\n", " The output file path for the raster file, by default 'raster.tiff'.\n", " \n", " Returns\n", " -------\n", " None\n", " The function writes the raster file and returns None.\n", "\n" ] } ], "source": [ "help(gm.to_raster) # function doc to write a raster file" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## CLI for gis transformation\n", "\n", "A command line script is also provided to transform a marthe grid file to shapefile or geopackage." ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "usage: martshp [-h] [--output OUTPUT] [--variable VARIABLE]\n", " [--xyfactor XYFACTOR] [--gpkg] [--mask] [--version]\n", " [--wide_fmt]\n", " [grid timesteps ...]\n", "\n", "Convert a Marthe GridFile to shapefile format.\n", "\n", "positional arguments:\n", " grid timesteps Paths to grid [and timesteps if result] files are\n", " expected\n", "\n", "options:\n", " -h, --help show this help message and exit\n", " --output OUTPUT, -o OUTPUT\n", " Output filename. Default is input.nc\n", " --variable VARIABLE, -v VARIABLE\n", " Variable (field) to read, default is None: i.e\n", " variable will be parsed from file and ONLY the first\n", " variable will be read. Pass 'all' to get all\n", " variables.\n", " --xyfactor XYFACTOR, -x XYFACTOR\n", " Transformation factor for coordinates. Optionnal,\n", " default is 1 (no transformation).\n", " --gpkg, -g Use GPKG format instead of shapefile\n", " --mask, -m Only get a mask of active domain\n", " --version, -V Show version and exit\n", " --wide_fmt, -w Use wide format (columns) for time\n", "\n", "gridmarthe Copyright (C) 2024 BRGM.\n", "This program comes with ABSOLUTELY NO WARRANTY.\n", "This is free software, see GNU General Public Lisence v3 for copying.\n" ] } ], "source": [ "!martshp -h" ] } ], "metadata": { "kernelspec": { "display_name": "gm", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.11.2" } }, "nbformat": 4, "nbformat_minor": 2 }