.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/13_polygon_flagging.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_auto_examples_13_polygon_flagging.py: Polygon Field Flagging ====================== This example demonstrates end-to-end polygon-based flagging: 1. Build a regular block model. 2. Define a named ``PolygonField`` from Shapely geometry. 3. Visualise the polygon against block centroids in 2D. 4. Evaluate a boolean mask from centroid-in-polygon tests (XY only). 5. Encode the result as a categorical block-model attribute. 6. Plot the final 3D block model by that categorical flag. .. GENERATED FROM PYTHON SOURCE LINES 14-28 .. code-block:: Python import tempfile from pathlib import Path import geopandas as gpd import matplotlib.pyplot as plt import numpy as np import pandas as pd import pyvista as pv from shapely.geometry import Polygon from parq_blockmodel import ParquetBlockModel, PolygonField, RegularGeometry .. GENERATED FROM PYTHON SOURCE LINES 29-32 Create a geometry-backed block model ------------------------------------ We start from geometry only, then materialise a canonical ``.pbm``. .. GENERATED FROM PYTHON SOURCE LINES 32-46 .. code-block:: Python temp_dir = Path(tempfile.gettempdir()) / "polygon_field_flagging_example" temp_dir.mkdir(parents=True, exist_ok=True) pbm_path = temp_dir / "polygon_field_flagging.pbm" geometry = RegularGeometry( corner=(0.0, 0.0, 0.0), block_size=(1.0, 1.0, 1.0), shape=(20, 16, 4), ) pbm = ParquetBlockModel.from_geometry(geometry=geometry, path=pbm_path) blocks = pbm.read(index=None) .. GENERATED FROM PYTHON SOURCE LINES 47-51 Persist and ingest a named polygon field (GeoParquet) ----------------------------------------------------- We write a simple named polygon table and then ingest by name. This demonstrates the persisted polygon interchange pattern. .. GENERATED FROM PYTHON SOURCE LINES 51-72 .. code-block:: Python lease_outline = Polygon( [ (3.0, 2.0), (15.0, 2.5), (17.0, 8.0), (13.0, 13.5), (6.0, 12.0), (2.5, 7.0), ] ) polygon_path = temp_dir / "lease_regions.parquet" regions = gpd.GeoDataFrame( {"name": ["lease_a"], "geometry": [lease_outline]}, crs="EPSG:4326", ) regions.to_parquet(polygon_path) polygon = PolygonField.from_geoparquet(polygon_path, name="lease_a") print("Polygon name:", polygon.name) .. rst-class:: sphx-glr-script-out .. code-block:: none Polygon name: lease_a .. GENERATED FROM PYTHON SOURCE LINES 73-76 Visualise the polygon against centroid locations in 2D ------------------------------------------------------ Evaluation is strictly in world XY; Z is ignored. .. GENERATED FROM PYTHON SOURCE LINES 76-109 .. code-block:: Python mask = polygon.evaluate(pbm.geometry) first_layer = blocks.loc[blocks["k"] == 0, ["block_id", "x", "y"]].copy() first_layer["inside"] = mask[first_layer["block_id"].to_numpy(dtype=np.int64)] fig, ax = plt.subplots(figsize=(7, 5)) ax.scatter( first_layer.loc[~first_layer["inside"], "x"], first_layer.loc[~first_layer["inside"], "y"], s=24, color="#bdbdbd", label="outside", ) ax.scatter( first_layer.loc[first_layer["inside"], "x"], first_layer.loc[first_layer["inside"], "y"], s=28, color="#1b9e77", label="inside", ) poly_x, poly_y = polygon.geometry.exterior.xy ax.plot(poly_x, poly_y, color="#222222", linewidth=2.0, label=f"{polygon.name} outline") ax.fill(poly_x, poly_y, color="#1b9e77", alpha=0.12) ax.set_title("Centroids (k=0) vs PolygonField") ax.set_xlabel("x") ax.set_ylabel("y") ax.set_aspect("equal", adjustable="box") ax.legend(loc="best") fig.tight_layout() plt.show() .. image-sg:: /auto_examples/images/sphx_glr_13_polygon_flagging_001.png :alt: Centroids (k=0) vs PolygonField :srcset: /auto_examples/images/sphx_glr_13_polygon_flagging_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 110-113 Encode the boolean mask as a categorical block-model attribute -------------------------------------------------------------- We map ``True`` to the polygon name and ``False`` to missing (``NA``). .. GENERATED FROM PYTHON SOURCE LINES 113-124 .. code-block:: Python inside_label = polygon.name or "inside" domain = np.full(len(blocks), pd.NA, dtype=object) domain[mask[blocks["block_id"].to_numpy(dtype=np.int64)]] = inside_label blocks["lease_domain"] = pd.Categorical(domain, categories=[inside_label], ordered=False) pbm.write(blocks[["block_id", "lease_domain"]], merge=True) print("lease_domain (NA-aware):") print(pbm.read(columns=["lease_domain"], index=None)["lease_domain"].value_counts(dropna=False)) .. rst-class:: sphx-glr-script-out .. code-block:: none lease_domain (NA-aware): lease_domain NaN 776 lease_a 504 Name: count, dtype: int64 .. GENERATED FROM PYTHON SOURCE LINES 125-127 Plot the final 3D block model using the encoded categorical domain ------------------------------------------------------------------ .. GENERATED FROM PYTHON SOURCE LINES 127-144 .. code-block:: Python plotter: pv.Plotter = pbm.plot(scalar="lease_domain", threshold=False, enable_picking=True, grid_type="image") # Keep perspective, but bias the camera toward an XY-aligned oblique view so # the 3D display is easier to relate to the 2D polygon map above. b = plotter.bounds # xmin, xmax, ymin, ymax, zmin, zmax cx = 0.5 * (b[0] + b[1]) cy = 0.5 * (b[2] + b[3]) cz = 0.5 * (b[4] + b[5]) dx = b[1] - b[0] dy = b[3] - b[2] span_xy = max(dx, dy) plotter.camera_position = [ (cx + 0.15 * dx, cy - 2.0 * dy, cz + 0.9 * span_xy), # camera position (cx, cy, cz), # focal point (0.0, 0.0, 1.0), # up-vector ] plotter.show() .. image-sg:: /auto_examples/images/sphx_glr_13_polygon_flagging_002.png :alt: 13 polygon flagging :srcset: /auto_examples/images/sphx_glr_13_polygon_flagging_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 1.221 seconds) .. _sphx_glr_download_auto_examples_13_polygon_flagging.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: 13_polygon_flagging.ipynb <13_polygon_flagging.ipynb>` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: 13_polygon_flagging.py <13_polygon_flagging.py>` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: 13_polygon_flagging.zip <13_polygon_flagging.zip>` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_