Iceland.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "82beda27-ce43-44ca-a7b4-d705469916fa",
   "metadata": {},
   "outputs": [],
   "source": [
    "import glob\n",
    "import os\n",
    "import io\n",
    "import pymap3d\n",
    "\n",
    "import cartopy.feature as cf\n",
    "import geopandas as gpd\n",
    "import holoviews as hv\n",
    "import hvplot.pandas\n",
    "import hvplot.xarray\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import shapely\n",
    "import thalassa\n",
    "import xarray as xr\n",
    "\n",
    "import searvey\n",
    "\n",
    "import pyposeidon\n",
    "import pyposeidon.meteo as pmeteo\n",
    "import pyposeidon.dem as pdem\n",
    "import pyposeidon.boundary as pbound\n",
    "import pyposeidon.mesh as pmesh\n",
    "import pyposeidon.model as pmodel\n",
    "from pyposeidon.utils import cast\n",
    "\n",
    "import pyposeidon.utils.hplot as hplot\n",
    "import pyposeidon.utils.pplot as pplot\n",
    "\n",
    "hv.extension(\"bokeh\")\n",
    "\n",
    "!mkdir -p data/iceland"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "472bb041-4b91-471f-a09f-45eb2b4acd78",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2024-07-10T14:25:52.300893Z",
     "iopub.status.busy": "2024-07-10T14:25:52.300644Z",
     "iopub.status.idle": "2024-07-10T14:25:52.332912Z",
     "shell.execute_reply": "2024-07-10T14:25:52.332032Z",
     "shell.execute_reply.started": "2024-07-10T14:25:52.300869Z"
    }
   },
   "source": [
    "# Geometry"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7b7d8fa2-d30c-449d-b000-f4bf316ca909",
   "metadata": {},
   "outputs": [],
   "source": [
    "lon_min = -28.0\n",
    "lon_max = -11.1\n",
    "lat_min =  62.0\n",
    "lat_max =  68.0\n",
    "\n",
    "bbox = shapely.box(lon_min, lat_min, lon_max, lat_max)\n",
    "geometry = dict(lon_min=lon_min, lon_max=lon_max, lat_min=lat_min, lat_max=lat_max)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6a61f1c2-171a-48c3-9af0-b346cf561237",
   "metadata": {},
   "source": [
    "# Coastlines"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4709fe5f",
   "metadata": {},
   "outputs": [],
   "source": [
    "OSM_FOLDER = \"/home/tomsail/work/python/seareport_org/coastlines/raw/osm/land-polygons-complete-4326\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "46e368f3-cd56-4323-8197-30c3c3d26553",
   "metadata": {},
   "outputs": [],
   "source": [
    "coastlines = gpd.read_file(OSM_FOLDER + \"/land_polygons.shp\", bbox=bbox)\n",
    "coastlines"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c0ce789c-a522-4a49-9d45-fd487d27c830",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2024-07-10T14:26:41.790292Z",
     "iopub.status.busy": "2024-07-10T14:26:41.789720Z",
     "iopub.status.idle": "2024-07-10T14:26:41.818903Z",
     "shell.execute_reply": "2024-07-10T14:26:41.818384Z",
     "shell.execute_reply.started": "2024-07-10T14:26:41.790214Z"
    }
   },
   "source": [
    "# Boundaries"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b78eb9c5-22b6-4f05-9f45-6ffd684af9b6",
   "metadata": {},
   "outputs": [],
   "source": [
    "boundary = pbound.Boundary(geometry=geometry, coastlines=coastlines)\n",
    "boundary.contours.head()\n",
    "len(boundary.contours)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "79adcfc1-f974-4d37-9b0b-88e06e544943",
   "metadata": {},
   "outputs": [],
   "source": [
    "boundary.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f735e7ee-c001-40f0-aec8-632a1f0bbd6a",
   "metadata": {},
   "outputs": [],
   "source": [
    "# boundary.contours.hvplot(geo=True, tiles=True, frame_height=500)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0a581850-d139-4f52-9dcb-b8b866c40905",
   "metadata": {},
   "source": [
    "# DEM"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d4c6e7a4-e435-443c-8412-65329d7a928f",
   "metadata": {},
   "outputs": [],
   "source": [
    "# url = \"https://coastwatch.pfeg.noaa.gov/erddap/griddap/srtm30plus\"\n",
    "url = \"data/iceland.nc\"\n",
    "dem = pdem.Dem(dem_source=url, **geometry)\n",
    "dem.Dataset.load()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5fcdf1f5",
   "metadata": {},
   "outputs": [],
   "source": [
    "dem.Dataset.to_netcdf(\"./data/iceland/dem.nc\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c8e13b4f-5e85-459a-b86c-45a30a1f8aa2",
   "metadata": {},
   "source": [
    "# Mesh"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d9665b26-4b56-497f-9461-7adfa8898b12",
   "metadata": {},
   "outputs": [],
   "source": [
    "mesh = pmesh.set(\n",
    "    mesh_generator='oceanmesh',\n",
    "    bgmesh = \"om\",\n",
    "    dem_source=\"./data/iceland/dem.nc\",\n",
    "    type='tri2d',\n",
    "    geometry=geometry,\n",
    "    coastlines=coastlines,\n",
    "    grad = 0.15,\n",
    "    bathy_gradient= True,\n",
    "    resolution_min=0.02,\n",
    "    resolution_max=1.50,\n",
    "    alpha_wavelength= 100,  # number of element to resolve WL\n",
    "    alpha_slope= 10,  # number of element to resolve bathy gradient\n",
    ")\n",
    "mesh.Dataset.dims"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5386162e",
   "metadata": {},
   "outputs": [],
   "source": [
    "mesh.Dataset.pplot.mesh()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "cc9478be-2bfe-4060-9e8c-c91d1230a5ac",
   "metadata": {},
   "outputs": [],
   "source": [
    "mesh.to_file('./data/iceland/hgrid.gr3')\n",
    "!head ./data/iceland/hgrid.gr3"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "88223b90-0fa6-43d3-b739-61f30d9f4b0e",
   "metadata": {},
   "source": [
    "### Fix bathy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7e91ddf8-72be-484f-ad20-ccc3bd007d7d",
   "metadata": {},
   "outputs": [],
   "source": [
    "#define in a dictionary the properties of the model..\n",
    "model_parameters = {\n",
    "    \"solver_name\": \"telemac\",\n",
    "    \"rpath\": \"./data/iceland/\",\n",
    "    \"dem_source\": \"./data/iceland/dem.nc\",\n",
    "    \"mesh_file\": \"./data/iceland/hgrid.gr3\",\n",
    "    \"update\": [\"dem\"], #set which component should be updated  (meteo,dem,model)\n",
    "    \"start_date\": \"2018-10-01T00:00:00\",\n",
    "    \"time_frame\": \"1d\",\n",
    "    \"global\": False,\n",
    "}\n",
    "model_parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "67ac2c07-bb87-4338-a155-de5a362fde9a",
   "metadata": {},
   "outputs": [],
   "source": [
    "model = pmodel.set(**model_parameters)\n",
    "model.create()\n",
    "model.mesh.to_file('./data/iceland/hgrid.gr3')\n",
    "!head ./data/iceland/hgrid.gr3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "113fe1f8-f5c0-44e1-a916-6b431c786145",
   "metadata": {},
   "outputs": [],
   "source": [
    "def parse_hgrid_nodes(path: os.PathLike[str] | str) -> pd.DataFrame:\n",
    "    with open(path, \"rb\") as fd:\n",
    "        _ = fd.readline()\n",
    "        _, no_points = map(int, fd.readline().strip().split(b\" \"))\n",
    "        content = io.BytesIO(b''.join(next(fd) for _ in range(no_points)))\n",
    "        nodes = pd.read_csv(\n",
    "            content,\n",
    "            engine=\"pyarrow\",\n",
    "            sep=\"\\t\",\n",
    "            header=None,\n",
    "            names=[\"lon\", \"lat\", \"depth\"],\n",
    "            index_col=0\n",
    "        )\n",
    "    nodes = nodes.reset_index(drop=True)\n",
    "    return nodes\n",
    "    \n",
    "def parse_hgrid_elements3(path: os.PathLike[str] | str) -> pd.DataFrame:\n",
    "    with open(path, \"rb\") as fd:\n",
    "        _ = fd.readline()\n",
    "        no_elements, no_points = map(int, fd.readline().strip().split(b\" \"))\n",
    "        for _ in range(no_points):\n",
    "            next(fd) \n",
    "        content = io.BytesIO(b''.join(next(fd) for _ in range(no_elements)))\n",
    "        elements = pd.read_csv(\n",
    "            content,\n",
    "            engine=\"pyarrow\",\n",
    "            sep=\"\\t\",\n",
    "            header=None,\n",
    "            names=[\"no_nodes\", \"n1\", \"n2\", \"n3\"],\n",
    "            index_col=0\n",
    "        )\n",
    "    elements = elements.assign(\n",
    "        n1=elements.n1 - 1,\n",
    "        n2=elements.n2 - 1,\n",
    "        n3=elements.n3 - 1,\n",
    "    ).reset_index(drop=True)\n",
    "    return elements\n",
    "\n",
    "def get_skews_and_base_cfls(lons, lats, depths) -> np.ndarray:\n",
    "    # The shape of each one of the input arrays needs to be (3, <no_triangles>)\n",
    "    #ell = pymap3d.Ellipsoid.from_name(\"wgs84\")\n",
    "    ell = pymap3d.Ellipsoid(6378206.4, 6378206.4, \"schism\", \"schism\")\n",
    "    local_x, local_y, _ = pymap3d.geodetic2enu(lats, lons, depths, lats[0], lons[0], depths[0], ell=ell)\n",
    "    areas = (local_x[1] * local_y[2] - local_x[2] * local_y[1]) * 0.5\n",
    "    rhos = np.sqrt(areas / np.pi)\n",
    "    max_sides = np.maximum(\n",
    "        np.sqrt(local_x[1] ** 2 + local_y[1] ** 2),\n",
    "        np.sqrt(local_x[2] ** 2 + local_y[2] ** 2),\n",
    "        np.sqrt((local_x[2] - local_x[1]) ** 2 + (local_y[2] - local_y[1]) ** 2),\n",
    "    )\n",
    "    skews = max_sides / rhos\n",
    "    base_cfls = np.sqrt(9.81 * np.maximum(0.1, depths.mean(axis=0))) / rhos / 2\n",
    "    return skews, base_cfls\n",
    "\n",
    "def get_skews_and_base_cfls_from_path(path: os.PathLike[str] | str) -> np.ndarray:\n",
    "    nodes = parse_hgrid_nodes(path)\n",
    "    elements = parse_hgrid_elements3(path)\n",
    "    tri = elements[[\"n1\", \"n2\", \"n3\"]].values\n",
    "    lons = nodes.lon.values[tri].T\n",
    "    lats = nodes.lat.values[tri].T\n",
    "    depths = nodes.depth.values[tri].T\n",
    "    skews, base_cfls = get_skews_and_base_cfls(lons=lons, lats=lats, depths=depths)\n",
    "    return skews, base_cfls\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7e7a4f00-01da-4a80-9825-23abf9df1c99",
   "metadata": {},
   "outputs": [],
   "source": [
    "skews, base_cfls = get_skews_and_base_cfls_from_path(\"./data/iceland/hgrid.gr3\")\n",
    "CFL_THRESHOLD = 0.4\n",
    "print(f\"elements violating CFL threshold < {CFL_THRESHOLD}\")\n",
    "print(\"time            N         %\")\n",
    "for dt in (1, 50, 75, 100, 120, 150, 200, 300, 400, 600, 900, 1200, 1800, 3600):\n",
    "    violations = (base_cfls * dt < CFL_THRESHOLD).sum()\n",
    "    print(f\"{dt:>4d} {violations:>12d} {violations / len(base_cfls) * 100:>8.2f}%\")\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d50c240d-24f4-4b30-82b7-6f6d7b100fd0",
   "metadata": {},
   "outputs": [],
   "source": [
    "pd.DataFrame({\"skew\": skews}).describe([0.25, 0.5, 0.75, 0.9, 0.95, 0.99, 0.995, 0.999])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "eaae98ef",
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_meta() -> gpd.GeoDataFrame:\n",
    "    meta_web = searvey.get_ioc_stations().drop(columns=[\"lon\", \"lat\"])\n",
    "    meta_api = (\n",
    "        pd.read_json(\n",
    "            \"http://www.ioc-sealevelmonitoring.org/service.php?query=stationlist&showall=all\"\n",
    "        )\n",
    "        .drop_duplicates()\n",
    "        .drop(columns=[\"lon\", \"lat\"])\n",
    "        .rename(columns={\"Code\": \"ioc_code\", \"Lon\": \"longitude\", \"Lat\": \"latitude\"})\n",
    "    )\n",
    "    merged = pd.merge(\n",
    "        meta_web,\n",
    "        meta_api[[\"ioc_code\", \"longitude\", \"latitude\"]].drop_duplicates(),\n",
    "        on=[\"ioc_code\"],\n",
    "    )\n",
    "    updated = merged.assign(\n",
    "        geometry=gpd.points_from_xy(merged.longitude, merged.latitude, crs=\"EPSG:4326\")\n",
    "    )\n",
    "    return updated\n",
    "\n",
    "ioc_ = get_meta()\n",
    "ioc_[bbox.contains(ioc_.geometry)].to_csv('data/iceland/stations.csv')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "fb76718f-7816-46f1-8f54-72ba32433239",
   "metadata": {},
   "outputs": [],
   "source": [
    "#define in a dictionary the properties of the model..\n",
    "model_parameters = {\n",
    "    \"solver_name\": \"telemac\",\n",
    "    \"tag\": \"telemac2d\",\n",
    "    \"rpath\": \"./data/iceland/20181001\",\n",
    "    \"mesh_file\": \"./data/iceland/hgrid.gr3\",\n",
    "    \"update\": [\"all\"], #set which component should be updated  (meteo,dem,model)\n",
    "    \"meteo_source\": glob.glob(\"data/uvp_*.grib\"),\n",
    "    \"meteo_merge\": \"last\",  # combine meteo\n",
    "    \"meteo_combine_by\": \"nested\",\n",
    "    \"meteo_xr_kwargs\": {\"concat_dim\": \"step\"},\n",
    "    \"start_date\": \"2018-10-01T00:00:00\",\n",
    "    \"time_frame\": \"2d\",\n",
    "    \"obs\": \"data/iceland/stations.csv\",\n",
    "    \"monitor\": True,\n",
    "    \"parameters\": {\n",
    "        \"dt\": 100\n",
    "    }\n",
    "}\n",
    "model_parameters"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "05b44459",
   "metadata": {},
   "source": [
    "## run days 1&2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0ee58df7",
   "metadata": {},
   "outputs": [],
   "source": [
    "a = pmodel.set(**model_parameters)\n",
    "a.create()\n",
    "# IMPORTANT! Here, for a simple surge application, \n",
    "# we will need close all boundaries, otherwise the \n",
    "# model will run out of water\n",
    "a.mesh.Dataset.type[:] = 'closed' # it will create a cli file with all boundaries closed (this can be done only once)\n",
    "# we also need to drop some meteo variables, it is necesarry for zarr export\n",
    "a.output(**{\"global\": False})\n",
    "a.set_obs()\n",
    "a.save() # saves the json model file\n",
    "a.run() # runs the model"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c13a8dfb",
   "metadata": {},
   "source": [
    "## run days 3&4 from hotstart"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a64f3563",
   "metadata": {},
   "outputs": [],
   "source": [
    "# restart model\n",
    "prev_ = pd.Timestamp('2018-10-01')\n",
    "next_ = pd.Timestamp('2018-10-03')\n",
    "end_ = pd.Timestamp('2018-10-05')\n",
    "ppath = os.path.join('data/iceland', prev_.strftime(\"%Y%m%d\"))\n",
    "npath = os.path.join('data/iceland', next_.strftime(\"%Y%m%d\"))\n",
    "m = pyposeidon.model.read(os.path.join(ppath, \"telemac2d_model.json\"))\n",
    "meteo = pmeteo.Meteo(glob.glob(\"data/uvp_*.grib\"),meteo_merge= \"last\", meteo_combine_by= \"nested\", meteo_xr_kwargs= {\"concat_dim\": \"step\"},)\n",
    "rs = cast.set(\n",
    "    solver_name=\"telemac\",\n",
    "    model=m,\n",
    "    ppath=ppath,  # old path\n",
    "    cpath=npath,  # new path\n",
    "    meteo=meteo.Dataset.sel(time=slice(next_, end_)).compute(),\n",
    "    sdate=next_,  # new start date\n",
    "    end_date=end_,  # new end date\n",
    "    start=next_,  # start\n",
    "    copy=True,\n",
    ")\n",
    "b = rs.run(execute=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a6405f1b",
   "metadata": {},
   "outputs": [],
   "source": [
    "b.run()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "af15ee7b",
   "metadata": {},
   "source": [
    "## run 4 days model - check"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9624eb1a",
   "metadata": {},
   "outputs": [],
   "source": [
    "#define in a dictionary the properties of the model..\n",
    "model_parameters = {\n",
    "    \"solver_name\": \"telemac\",\n",
    "    \"tag\": \"telemac2d\",\n",
    "    \"rpath\": \"./data/iceland/20181001-04\",\n",
    "    \"mesh_file\": \"./data/iceland/hgrid.gr3\",\n",
    "    \"update\": [\"all\"], #set which component should be updated  (meteo,dem,model)\n",
    "    \"meteo_source\": glob.glob(\"data/uvp_*.grib\"),\n",
    "    \"meteo_merge\": \"last\",  # combine meteo\n",
    "    \"meteo_combine_by\": \"nested\",\n",
    "    \"meteo_xr_kwargs\": {\"concat_dim\": \"step\"},\n",
    "    \"start_date\": \"2018-10-01T00:00:00\",\n",
    "    \"time_frame\": \"4d\",\n",
    "    \"global\": False,\n",
    "    \"obs\": \"data/iceland/stations.csv\",\n",
    "    \"monitor\": True,\n",
    "    \"parameters\": {\n",
    "        \"dt\": 100\n",
    "    }\n",
    "}\n",
    "model_parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9787559b",
   "metadata": {},
   "outputs": [],
   "source": [
    "c = pmodel.set(**model_parameters)\n",
    "c.create()\n",
    "# IMPORTANT! Here, for a simple surge application, \n",
    "# we will need close all boundaries, otherwise the \n",
    "# model will run out of water\n",
    "c.mesh.Dataset.type[:] = 'closed' # it will create a cli file with all boundaries closed (this can be done only once)\n",
    "# we also need to drop some meteo variables, it is necesarry for zarr export\n",
    "c.meteo.Dataset = c.meteo.Dataset.compute()\n",
    "c.output()\n",
    "c.set_obs()\n",
    "c.save() # saves the json model file\n",
    "c.run() # runs the model"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5da57efc",
   "metadata": {},
   "source": [
    "## compare results"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "bce6f374",
   "metadata": {},
   "outputs": [],
   "source": [
    "res_2days = xr.open_dataset(\"data/iceland/20181001-04/results_2D.slf\")\n",
    "res_day1 = xr.open_dataset(\"data/iceland/20181001/results_2D.slf\")\n",
    "res_day2 = xr.open_dataset(\"data/iceland/20181003/results_2D.slf\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "85bff31a",
   "metadata": {},
   "outputs": [],
   "source": [
    "node = 15000\n",
    "p1 = res_day1.S.isel(node = node).hvplot(label = \"Day 1&2\")\n",
    "p2 = res_day2.S.isel(node = node).hvplot(label = \"Day 3&4\")\n",
    "p3 = res_2days.S.isel(node = node).hvplot(label = \"4 Days\", line_dash='dashed')\n",
    "p1 * p2 * p3"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4692bce5",
   "metadata": {},
   "source": [
    "## compare with observations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "32e3a069",
   "metadata": {},
   "outputs": [],
   "source": [
    "res_2days = xr.open_dataset(\"data/iceland/20181001-04/results_1D.slf\")\n",
    "res_day1 = xr.open_dataset(\"data/iceland/20181001/results_1D.slf\")\n",
    "res_day2 = xr.open_dataset(\"data/iceland/20181003/results_1D.slf\")\n",
    "res_day2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8eb42188",
   "metadata": {},
   "outputs": [],
   "source": [
    "stations = pd.read_csv(\"data/iceland/stations.csv\")\n",
    "stations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c6e13834",
   "metadata": {},
   "outputs": [],
   "source": [
    "from searvey import ioc\n",
    "data = ioc.get_ioc_station_data('reyk', endtime=\"2018-10-05\", period=30)\n",
    "data.index = data['time']\n",
    "data = data.drop(columns=['time'])\n",
    "data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "43b9e48a",
   "metadata": {},
   "outputs": [],
   "source": [
    "# detide \n",
    "from analysea.tide import detide\n",
    "surge = detide(data[\"prs\"],lat = 64.15)\n",
    "surge"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e487c500",
   "metadata": {},
   "outputs": [],
   "source": [
    "p1 = res_day1.S.isel(node = 0).hvplot(label = \"Day 1\")\n",
    "p2 = res_day2.S.isel(node = 0).hvplot(label = \"Day 2\")\n",
    "p3 = res_2days.S.isel(node = 0).hvplot(label = \"2 Days\", line_dash='dashed')\n",
    "obs_ = surge.loc[\"2018-10-01\":\"2018-10-05\"].hvplot(label = \"Observations\", color = 'k', line_dash='dotted')\n",
    "p1 * p2  * p3 * obs_"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "pos_test",
   "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.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}