Release zones (#127)

* Shorten lines * Tidy module * Make target length argument optional * Remove redundant function * Update cookbook Remove dependency on preprepared shapefile data when demonstrating how to create a set of release zones around a polygon representing an island or another object. * Remove old reference shapefile data * Release zone class Add coordinate system argument to ReleaseZone init method which makes it possible to creeate release zones in cartesian or geographic coordinates. The class handles the conversion to UTM coordinates internally should a set of geographic coordiantes be passed in. * Add PyLagAttributeError exception * Switch to using properties * Shorten line * Expand ReleaseZone's interface * Release zone creation Enable support for creating single release zones or release zones along a cord by passing in either Cartesian or geographic coordiantes. * Release zones Further support for defining and using release zones using either cartesian or geographic coordinates. * Align input module with new changes * Initial positions documentation examples Make examples work given the recent changes to release_zone.py * Return copies of mutable objects * Sync more notebooks * Sync global example notebook * Sync the final set of doc examples * Sync release zone unit tests * Print PyLag git commit in pylag-package-conda-linux.yml * Update pylag-package-conda-linux.yml
pmlmodelling · Nov 22, 2023 · d9a8e13 · d9a8e13
1 parent 3e1791d
commit d9a8e13
Show file tree

Hide file tree

Showing 17 changed files with 879 additions and 618 deletions.
diff --git a/.github/workflows/pylag-package-conda-linux.yml b/.github/workflows/pylag-package-conda-linux.yml
@@ -27,6 +27,11 @@ jobs:
         run: |
           conda info
           conda list
+      - name: Install local changes
+           shell: bash -l {0}
+        run: |
+            pip install .
+            python -c "import pylag; print(f'PyLag Git commit: {pylag.version.git_revision}')"
 #      - name: Lint
 #        shell: bash -l {0}
 #        run: |

diff --git a/doc/source/documentation/cookbook/initial_positions.ipynb b/doc/source/documentation/cookbook/initial_positions.ipynb
diff --git a/doc/source/documentation/resources/ref_shapefile.dbf b/doc/source/documentation/resources/ref_shapefile.dbf
diff --git a/doc/source/documentation/resources/ref_shapefile.prj b/doc/source/documentation/resources/ref_shapefile.prj
diff --git a/doc/source/documentation/resources/ref_shapefile.qpj b/doc/source/documentation/resources/ref_shapefile.qpj
diff --git a/doc/source/documentation/resources/ref_shapefile.shp b/doc/source/documentation/resources/ref_shapefile.shp
diff --git a/doc/source/documentation/resources/ref_shapefile.shx b/doc/source/documentation/resources/ref_shapefile.shx
diff --git a/doc/source/examples/arakawa_a_forward_tracking.ipynb b/doc/source/examples/arakawa_a_forward_tracking.ipynb
@@ -66,10 +66,6 @@
     "lat = 50.32\n",
     "lon = -4.17\n",
     "\n",
-    "# Convert to UTM coordinates\n",
-    "epsg_code = '32630'\n",
-    "easting, northing, _ = utm_from_lonlat([lon], [lat], epsg_code=epsg_code)\n",
-    "\n",
     "# Release zone radius (m)\n",
     "radius = 200.0\n",
     "\n",
@@ -84,31 +80,27 @@
     "# Create the release zone\n",
     "surface_release_zone = create_release_zone(group_id = group_id,\n",
     "                                           radius = radius,\n",
-    "                                           centre = [easting, northing],\n",
+    "                                           centre = [lon, lat],\n",
+    "                                           coordinate_system = 'geographic',\n",
     "                                           n_particles = n_particles_target,\n",
     "                                           depth = depth_below_surface,\n",
     "                                           random = False)\n",
     "\n",
     "# Get the actual number of particles\n",
     "n_particles = surface_release_zone.get_number_of_particles()\n",
     "\n",
-    "# Convert to lat/lon coordinates\n",
-    "lons, lats = lonlat_from_utm(surface_release_zone.get_eastings(),\n",
-    "                             surface_release_zone.get_northings(),\n",
-    "                             epsg_code=epsg_code)\n",
-    "\n",
-    "# Get depths\n",
-    "depths = surface_release_zone.get_depths()\n",
+    "# Get coordinates\n",
+    "lons, lats, depths = surface_release_zone.get_coordinates()\n",
     "\n",
     "# Create input sub-directory\n",
-    "input_dir = '{}/input'.format(simulation_dir)\n",
+    "input_dir = f'{simulation_dir}/input'\n",
     "try:\n",
     "    os.makedirs(input_dir)\n",
     "except FileExistsError:\n",
     "    pass\n",
     "\n",
     "# Output filename\n",
-    "file_name = '{}/initial_positions.dat'.format(input_dir)\n",
+    "file_name = f'{input_dir}/initial_positions.dat'\n",
     "\n",
     "# Write data to file\n",
     "create_initial_positions_file_single_group(file_name,\n",
@@ -246,7 +238,8 @@
     "# Create figure\n",
     "font_size = 15\n",
     "cmap = colourmap('h_r')\n",
-    "fig, ax = create_figure(figure_size=(26., 26.), projection=ccrs.PlateCarree(), font_size=font_size, bg_color='gray')\n",
+    "fig, ax = create_figure(figure_size=(26., 26.), projection=ccrs.PlateCarree(),\n",
+    "                        font_size=font_size, bg_color='gray')\n",
     "\n",
     "# Configure plotter\n",
     "plotter = ArakawaAPlotter(grid_metrics_file_name,\n",
@@ -261,8 +254,9 @@
     "plotter.scatter(ax, [-4.17], [50.25], marker='o', c='r')\n",
     "x, y = ax.projection.transform_point(-4.17, 50.25, src_crs=ccrs.PlateCarree())\n",
     "ax.annotate('L4', xy=(x, y), xytext=(0.2, 0.35), textcoords='axes fraction',\n",
-    "             arrowprops=dict(facecolor='red', edgecolor='red', width=1, headwidth=10, headlength=20, shrink=0.05),\n",
-    "             fontsize=font_size, color='red')\n",
+    "            arrowprops=dict(facecolor='red', edgecolor='red', width=1, headwidth=10,\n",
+    "                            headlength=20, shrink=0.05),\n",
+    "            fontsize=font_size, color='red')\n",
     "plotter.set_title(ax, 'Full grid (without grid lines)')\n",
     "\n",
     "# Grid subset (centered on Rame Head)\n",
@@ -271,7 +265,8 @@
     "# Create figure\n",
     "font_size = 15\n",
     "cmap = colourmap('h_r')\n",
-    "fig, ax = create_figure(figure_size=(26., 26.), projection=ccrs.PlateCarree(), font_size=font_size, bg_color='gray')\n",
+    "fig, ax = create_figure(figure_size=(26., 26.), projection=ccrs.PlateCarree(),\n",
+    "                        font_size=font_size, bg_color='gray')\n",
     "\n",
     "# Configure plotter\n",
     "plotter = ArakawaAPlotter(grid_metrics_file_name,\n",
@@ -496,7 +491,8 @@
     "# Plot extents\n",
     "extents = np.array([-4.3, -4.05, 50.25, 50.4], dtype=float)\n",
     "\n",
-    "fig, ax = create_figure(figure_size=(26., 26.), projection=ccrs.PlateCarree(), font_size=font_size, bg_color='gray')\n",
+    "fig, ax = create_figure(figure_size=(26., 26.), projection=ccrs.PlateCarree(),\n",
+    "                        font_size=font_size, bg_color='gray')\n",
     "\n",
     "# Create plot of release \n",
     "plotter = ArakawaAPlotter(grid_metrics_file_name,\n",
@@ -549,7 +545,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.9.15"
+   "version": "3.11.4"
   }
  },
  "nbformat": 4,

diff --git a/doc/source/examples/fvcom_backward_tracking.ipynb b/doc/source/examples/fvcom_backward_tracking.ipynb
@@ -65,10 +65,6 @@
     "lat = 50.32\n",
     "lon = -4.17\n",
     "\n",
-    "# Convert to UTM coordinates\n",
-    "epsg_code = get_epsg_code(lon, lat)\n",
-    "easting, northing, _ = utm_from_lonlat([lon], [lat], epsg_code=epsg_code)\n",
-    "\n",
     "# Release zone radius (m)\n",
     "radius = 200.0\n",
     "\n",
@@ -83,14 +79,22 @@
     "# Create the release zone\n",
     "surface_release_zone = create_release_zone(group_id = group_id,\n",
     "                                           radius = radius,\n",
-    "                                           centre = [easting, northing],\n",
+    "                                           centre = [lon, lat],\n",
+    "                                           coordinate_system = 'geographic',\n",
     "                                           n_particles = n_particles_target,\n",
     "                                           depth = depth_below_surface,\n",
     "                                           random = False)\n",
     "\n",
     "# Get the actual number of particles\n",
     "n_particles = surface_release_zone.get_number_of_particles()\n",
     "\n",
+    "# Get coordinates of the particles. NB - we're going to run\n",
+    "# the model in cartesian coordinates. Objects of type\n",
+    "# ReleaseZone have a method that allows one to transform\n",
+    "# geographic coordinates into UTM coordinates with units\n",
+    "# of m. We use this here.\n",
+    "eastings, northings, depths, epsg_code = surface_release_zone.get_utm_coordinates()\n",
+    "\n",
     "# Create input sub-directory\n",
     "input_dir = f\"{simulation_dir}/input\"\n",
     "try:\n",
@@ -105,9 +109,9 @@
     "create_initial_positions_file_single_group(file_name,\n",
     "                                           n_particles,\n",
     "                                           group_id,\n",
-    "                                           surface_release_zone.get_eastings(),\n",
-    "                                           surface_release_zone.get_northings(),\n",
-    "                                           surface_release_zone.get_depths())"
+    "                                           eastings,\n",
+    "                                           northings,\n",
+    "                                           depths)"
    ]
   },
   {
@@ -123,18 +127,7 @@
    "cell_type": "code",
    "execution_count": 3,
    "metadata": {},
-   "outputs": [
-    {
-     "data": {
-      "text/plain": [
-       "'/home/jcl/code/git/PyLag/PyLag/doc/source/examples/simulations/fvcom_backward/input/grid_metrics.nc'"
-      ]
-     },
-     "execution_count": 3,
-     "metadata": {},
-     "output_type": "execute_result"
-    }
-   ],
+   "outputs": [],
    "source": [
     "from shutil import copyfile\n",
     "\n",
@@ -210,7 +203,18 @@
    "cell_type": "code",
    "execution_count": 5,
    "metadata": {},
-   "outputs": [],
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "\n",
+      "Starting ensemble member 1 ...\n",
+      "Progress:\n",
+      "100% |###########################################|\n"
+     ]
+    }
+   ],
    "source": [
     "# Change to the run directory\n",
     "os.chdir(f\"{simulation_dir}\")\n",
@@ -282,7 +286,8 @@
     "# Plot extents\n",
     "extents = np.array([-4.3, -4.05, 50.25, 50.4], dtype=float)\n",
     "\n",
-    "fig, ax = create_figure(figure_size=(26., 26.), projection=ccrs.PlateCarree(), font_size=font_size, bg_color='gray')\n",
+    "fig, ax = create_figure(figure_size=(26., 26.), projection=ccrs.PlateCarree(),\n",
+    "                        font_size=font_size, bg_color='gray')\n",
     "\n",
     "# Configure plotter\n",
     "plotter = FVCOMPlotter(grid_metrics_file_name,\n",
@@ -329,7 +334,7 @@
  ],
  "metadata": {
   "kernelspec": {
-   "display_name": "Python 3",
+   "display_name": "Python 3 (ipykernel)",
    "language": "python",
    "name": "python3"
   },
@@ -343,7 +348,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.8.5"
+   "version": "3.11.4"
   }
  },
  "nbformat": 4,

diff --git a/doc/source/examples/fvcom_forward_tracking.ipynb b/doc/source/examples/fvcom_forward_tracking.ipynb
@@ -157,10 +157,6 @@
     "lat = 50.32\n",
     "lon = -4.17\n",
     "\n",
-    "# Convert to UTM coordinates\n",
-    "epsg_code = get_epsg_code(lon, lat)\n",
-    "easting, northing, _ = utm_from_lonlat([lon], [lat], epsg_code=epsg_code)\n",
-    "\n",
     "# Release zone radius (m)\n",
     "radius = 200.0\n",
     "\n",
@@ -175,14 +171,22 @@
     "# Create the release zone\n",
     "surface_release_zone = create_release_zone(group_id = group_id,\n",
     "                                           radius = radius,\n",
-    "                                           centre = [easting, northing],\n",
+    "                                           centre = [lon, lat],\n",
+    "                                           coordinate_system = 'geographic',\n",
     "                                           n_particles = n_particles_target,\n",
     "                                           depth = depth_below_surface,\n",
     "                                           random = False)\n",
     "\n",
     "# Get the actual number of particles\n",
     "n_particles = surface_release_zone.get_number_of_particles()\n",
     "\n",
+    "# Get coordinates of the particles. NB - we're going to run\n",
+    "# the model in cartesian coordinates. Objects of type\n",
+    "# ReleaseZone have a method that allows one to transform\n",
+    "# geographic coordinates into UTM coordinates with units\n",
+    "# of m. We use this here.\n",
+    "eastings, northings, depths, epsg_code = surface_release_zone.get_utm_coordinates() \n",
+    "\n",
     "# Create input sub-directory\n",
     "input_dir = f\"{simulation_dir}/input\"\n",
     "try:\n",
@@ -197,9 +201,9 @@
     "create_initial_positions_file_single_group(file_name,\n",
     "                                           n_particles,\n",
     "                                           group_id,\n",
-    "                                           surface_release_zone.get_eastings(),\n",
-    "                                           surface_release_zone.get_northings(),\n",
-    "                                           surface_release_zone.get_depths())\n"
+    "                                           eastings,\n",
+    "                                           northings,\n",
+    "                                           depths)"
    ]
   },
   {
@@ -226,10 +230,8 @@
     }
    ],
    "source": [
-    "# Convert utm coords to degrees\n",
-    "lons, lats = lonlat_from_utm(surface_release_zone.get_eastings(),\n",
-    "                             surface_release_zone.get_northings(),\n",
-    "                             epsg_code=epsg_code)\n",
+    "# Get lon/lat coordinates\n",
+    "lons, lats, _ = surface_release_zone.get_coordinates()\n",
     "\n",
     "# Create figure\n",
     "font_size = 15\n",
@@ -244,8 +246,8 @@
     "\n",
     "# Plot bathymetry\n",
     "extents = np.array([-4.21, -4.15, 50.30, 50.33], dtype=float)\n",
-    "ax, plot = plotter.plot_field(ax, bathy, extents=extents, add_colour_bar=True, cb_label='Depth (m)',\n",
-    "                              vmin=-60., vmax=0., cmap=cmap)\n",
+    "ax, plot = plotter.plot_field(ax, bathy, extents=extents, add_colour_bar=True,\n",
+    "                              cb_label='Depth (m)', vmin=-60., vmax=0., cmap=cmap)\n",
     "\n",
     "# Overlay grid\n",
     "plotter.draw_grid(ax, linewidth=1.0)\n",
@@ -377,7 +379,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": null,
+   "execution_count": 7,
    "metadata": {},
    "outputs": [],
    "source": [
@@ -473,8 +475,8 @@
     "                       font_size=font_size)\n",
     "\n",
     "# Plot the bathymetry again. We'll overlay pathlines on top of this.\n",
-    "plotter.plot_field(ax, bathy, extents=extents, add_colour_bar=True, cb_label='Depth (m)',\n",
-    "                   vmin=-60., vmax=0., cmap=cmap)\n",
+    "plotter.plot_field(ax, bathy, extents=extents, add_colour_bar=True,\n",
+    "                   cb_label='Depth (m)', vmin=-60., vmax=0., cmap=cmap)\n",
     "\n",
     "# Dataset holding particle positions\n",
     "viewer = Viewer(file_name, time_rounding=900)\n",
@@ -537,7 +539,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.9.15"
+   "version": "3.11.4"
   }
  },
  "nbformat": 4,

diff --git a/doc/source/examples/global_arakawa_a_forward_tracking.ipynb b/doc/source/examples/global_arakawa_a_forward_tracking.ipynb
@@ -74,11 +74,6 @@
     "lat = 28.5\n",
     "lon = -78.7\n",
     "\n",
-    "# Convert to UTM coordinates\n",
-    "easting, northing, epsg_code = utm_from_lonlat([lon], [lat])\n",
-    "easting = easting[0]\n",
-    "northing = northing[0]\n",
-    "\n",
     "# Release zone radius (m)\n",
     "radius = 10000.0\n",
     "\n",
@@ -91,7 +86,8 @@
     "# Create the release zone\n",
     "surface_release_zone = create_release_zone(group_id = group_id,\n",
     "                                           radius = radius,\n",
-    "                                           centre = [easting, northing],\n",
+    "                                           centre = [lon, lat],\n",
+    "                                           coordinate_system = 'geographic',\n",
     "                                           n_particles = n_particles_target,\n",
     "                                           depth = depth_below_surface,\n",
     "                                           random = True)\n",
@@ -100,12 +96,7 @@
     "n_particles = surface_release_zone.get_number_of_particles()\n",
     "\n",
     "# Convert to lat/lon coordinates\n",
-    "lons, lats = lonlat_from_utm(surface_release_zone.get_eastings(),\n",
-    "                             surface_release_zone.get_northings(),\n",
-    "                             epsg_code)\n",
-    "\n",
-    "# Get depths\n",
-    "depths = surface_release_zone.get_depths()\n",
+    "lons, lats, depths = surface_release_zone.get_coordinates()\n",
     "\n",
     "# Output filename\n",
     "file_name = '{}/initial_positions.dat'.format(input_dir)\n",
@@ -412,14 +403,14 @@
     "extents = np.array([-90, -30, 20, 50], dtype=float)\n",
     "\n",
     "# Projection\n",
-    "#projection = ccrs.NearsidePerspective(central_longitude=-60.0, central_latitude=35.0, satellite_height=35785800)\n",
     "projection = ccrs.PlateCarree()\n",
     "\n",
     "# Plot positions after 180 and 360 days\n",
     "for tof in [180, 360]:\n",
     "\n",
     "    # Create figure\n",
-    "    fig, ax1 = create_figure(figure_size=(26., 26.), projection=projection, font_size=font_size, bg_color='white')\n",
+    "    fig, ax1 = create_figure(figure_size=(26., 26.), projection=projection,\n",
+    "                             font_size=font_size, bg_color='white')\n",
     "\n",
     "    # Create plotter \n",
     "    plotter = ArakawaAPlotter(grid_metrics_file_name,\n",
@@ -469,7 +460,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.9.7"
+   "version": "3.11.4"
   }
  },
  "nbformat": 4,