Created
May 16, 2025 13:55
-
-
Save keewis/8e0437c83eac469873f6ab6e501d06ae to your computer and use it in GitHub Desktop.
cell selection with mocs
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
{ | |
"cells": [ | |
{ | |
"cell_type": "code", | |
"execution_count": 1, | |
"id": "0", | |
"metadata": {}, | |
"outputs": [], | |
"source": [ | |
"import numpy as np\n", | |
"import xarray as xr\n", | |
"import xdggs\n", | |
"import cdshealpix.nested\n", | |
"import astropy.coordinates as ac" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 2, | |
"id": "57728701-982a-4cda-9213-713b99c2fa18", | |
"metadata": {}, | |
"outputs": [], | |
"source": [ | |
"import shapely" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 30, | |
"id": "1", | |
"metadata": {}, | |
"outputs": [], | |
"source": [ | |
"ZOOM = 4\n", | |
"N = int(10e6)\n", | |
"da = xr.DataArray(\n", | |
" data=np.arange(N),\n", | |
" dims=(\"cell_ids\",),\n", | |
" coords={\n", | |
" \"cell_ids\": (\n", | |
" \"cell_ids\",\n", | |
" np.arange(N),\n", | |
" {\"grid_name\": \"healpix\", \"level\": ZOOM, \"indexing_scheme\": \"nested\"},\n", | |
" )\n", | |
" },\n", | |
").dggs.decode()\n", | |
"\n", | |
"\n", | |
"def sel_and_decode(da, cell_ids):\n", | |
" da = da.copy()\n", | |
"\n", | |
" return da.sel(\n", | |
" cell_ids=cell_ids,\n", | |
" ).dggs.decode()\n", | |
"\n", | |
"\n", | |
"def get_cell_ids(lon_lat_list: list[tuple[float, float]], zoom: int):\n", | |
" lon = ac.Longitude(lon_lat_list[:, 0], unit=\"deg\")\n", | |
" lat = ac.Latitude(lon_lat_list[:, 1], unit=\"deg\")\n", | |
" cell_ids, _, _ = cdshealpix.nested.polygon_search(lon, lat,\n", | |
" zoom,\n", | |
" flat=True,\n", | |
" )\n", | |
"\n", | |
" return cell_ids\n", | |
"\n", | |
"\n", | |
"def lon_lat_to_xyz(lon_lat: tuple[float, float]) -> tuple[float, float, float]:\n", | |
" # https://en.wikipedia.org/wiki/Spherical_coordinate_system\n", | |
" r = 1\n", | |
" lat, lon = lon_lat\n", | |
" lat = np.radians(lat)\n", | |
" lon = np.radians(lon)\n", | |
" x = r * np.cos(lat) * np.cos(lon)\n", | |
" y = r * np.cos(lat) * np.sin(lon)\n", | |
" z = r * np.sin(lat)\n", | |
" return x, y, z\n", | |
"\n", | |
"\n", | |
"def get_aus_bbox():\n", | |
" return (113.338953078, -43.6345972634, 153.569469029, -10.6681857235)\n", | |
"\n", | |
"polygon = shapely.box(*get_aus_bbox()).segmentize(0.1)\n", | |
"POLYGON = shapely.get_coordinates(polygon)\n", | |
"\n", | |
"cells = get_cell_ids(POLYGON, ZOOM)\n", | |
"cells1 = np.array(\n", | |
" [\n", | |
" 1344,\n", | |
" 1345,\n", | |
" 1347,\n", | |
" 1348,\n", | |
" 1349,\n", | |
" 1350,\n", | |
" 2492,\n", | |
" 2493,\n", | |
" 2495,\n", | |
" 2530,\n", | |
" 2531,\n", | |
" 2536,\n", | |
" 2537,\n", | |
" 2538,\n", | |
" 2539,\n", | |
" 2540,\n", | |
" 2542,\n", | |
" 2543,\n", | |
" ],\n", | |
" dtype=np.uint64,\n", | |
")\n", | |
"\n", | |
"da_truth = sel_and_decode(da, cells)\n", | |
"da1 = sel_and_decode(da, cells1)" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 16, | |
"id": "221ad3ce-648e-4741-87b2-c21c0db2a57b", | |
"metadata": {}, | |
"outputs": [ | |
{ | |
"data": { | |
"text/html": [ | |
"<div><svg style=\"position: absolute; width: 0; height: 0; overflow: hidden\">\n", | |
"<defs>\n", | |
"<symbol id=\"icon-database\" viewBox=\"0 0 32 32\">\n", | |
"<path d=\"M16 0c-8.837 0-16 2.239-16 5v4c0 2.761 7.163 5 16 5s16-2.239 16-5v-4c0-2.761-7.163-5-16-5z\"></path>\n", | |
"<path d=\"M16 17c-8.837 0-16-2.239-16-5v6c0 2.761 7.163 5 16 5s16-2.239 16-5v-6c0 2.761-7.163 5-16 5z\"></path>\n", | |
"<path d=\"M16 26c-8.837 0-16-2.239-16-5v6c0 2.761 7.163 5 16 5s16-2.239 16-5v-6c0 2.761-7.163 5-16 5z\"></path>\n", | |
"</symbol>\n", | |
"<symbol id=\"icon-file-text2\" viewBox=\"0 0 32 32\">\n", | |
"<path d=\"M28.681 7.159c-0.694-0.947-1.662-2.053-2.724-3.116s-2.169-2.030-3.116-2.724c-1.612-1.182-2.393-1.319-2.841-1.319h-15.5c-1.378 0-2.5 1.121-2.5 2.5v27c0 1.378 1.122 2.5 2.5 2.5h23c1.378 0 2.5-1.122 2.5-2.5v-19.5c0-0.448-0.137-1.23-1.319-2.841zM24.543 5.457c0.959 0.959 1.712 1.825 2.268 2.543h-4.811v-4.811c0.718 0.556 1.584 1.309 2.543 2.268zM28 29.5c0 0.271-0.229 0.5-0.5 0.5h-23c-0.271 0-0.5-0.229-0.5-0.5v-27c0-0.271 0.229-0.5 0.5-0.5 0 0 15.499-0 15.5 0v7c0 0.552 0.448 1 1 1h7v19.5z\"></path>\n", | |
"<path d=\"M23 26h-14c-0.552 0-1-0.448-1-1s0.448-1 1-1h14c0.552 0 1 0.448 1 1s-0.448 1-1 1z\"></path>\n", | |
"<path d=\"M23 22h-14c-0.552 0-1-0.448-1-1s0.448-1 1-1h14c0.552 0 1 0.448 1 1s-0.448 1-1 1z\"></path>\n", | |
"<path d=\"M23 18h-14c-0.552 0-1-0.448-1-1s0.448-1 1-1h14c0.552 0 1 0.448 1 1s-0.448 1-1 1z\"></path>\n", | |
"</symbol>\n", | |
"</defs>\n", | |
"</svg>\n", | |
"<style>/* CSS stylesheet for displaying xarray objects in jupyterlab.\n", | |
" *\n", | |
" */\n", | |
"\n", | |
":root {\n", | |
" --xr-font-color0: var(--jp-content-font-color0, rgba(0, 0, 0, 1));\n", | |
" --xr-font-color2: var(--jp-content-font-color2, rgba(0, 0, 0, 0.54));\n", | |
" --xr-font-color3: var(--jp-content-font-color3, rgba(0, 0, 0, 0.38));\n", | |
" --xr-border-color: var(--jp-border-color2, #e0e0e0);\n", | |
" --xr-disabled-color: var(--jp-layout-color3, #bdbdbd);\n", | |
" --xr-background-color: var(--jp-layout-color0, white);\n", | |
" --xr-background-color-row-even: var(--jp-layout-color1, white);\n", | |
" --xr-background-color-row-odd: var(--jp-layout-color2, #eeeeee);\n", | |
"}\n", | |
"\n", | |
"html[theme=\"dark\"],\n", | |
"html[data-theme=\"dark\"],\n", | |
"body[data-theme=\"dark\"],\n", | |
"body.vscode-dark {\n", | |
" --xr-font-color0: rgba(255, 255, 255, 1);\n", | |
" --xr-font-color2: rgba(255, 255, 255, 0.54);\n", | |
" --xr-font-color3: rgba(255, 255, 255, 0.38);\n", | |
" --xr-border-color: #1f1f1f;\n", | |
" --xr-disabled-color: #515151;\n", | |
" --xr-background-color: #111111;\n", | |
" --xr-background-color-row-even: #111111;\n", | |
" --xr-background-color-row-odd: #313131;\n", | |
"}\n", | |
"\n", | |
".xr-wrap {\n", | |
" display: block !important;\n", | |
" min-width: 300px;\n", | |
" max-width: 700px;\n", | |
"}\n", | |
"\n", | |
".xr-text-repr-fallback {\n", | |
" /* fallback to plain text repr when CSS is not injected (untrusted notebook) */\n", | |
" display: none;\n", | |
"}\n", | |
"\n", | |
".xr-header {\n", | |
" padding-top: 6px;\n", | |
" padding-bottom: 6px;\n", | |
" margin-bottom: 4px;\n", | |
" border-bottom: solid 1px var(--xr-border-color);\n", | |
"}\n", | |
"\n", | |
".xr-header > div,\n", | |
".xr-header > ul {\n", | |
" display: inline;\n", | |
" margin-top: 0;\n", | |
" margin-bottom: 0;\n", | |
"}\n", | |
"\n", | |
".xr-obj-type,\n", | |
".xr-array-name {\n", | |
" margin-left: 2px;\n", | |
" margin-right: 10px;\n", | |
"}\n", | |
"\n", | |
".xr-obj-type {\n", | |
" color: var(--xr-font-color2);\n", | |
"}\n", | |
"\n", | |
".xr-sections {\n", | |
" padding-left: 0 !important;\n", | |
" display: grid;\n", | |
" grid-template-columns: 150px auto auto 1fr 0 20px 0 20px;\n", | |
"}\n", | |
"\n", | |
".xr-section-item {\n", | |
" display: contents;\n", | |
"}\n", | |
"\n", | |
".xr-section-item input {\n", | |
" display: inline-block;\n", | |
" opacity: 0;\n", | |
" height: 0;\n", | |
"}\n", | |
"\n", | |
".xr-section-item input + label {\n", | |
" color: var(--xr-disabled-color);\n", | |
"}\n", | |
"\n", | |
".xr-section-item input:enabled + label {\n", | |
" cursor: pointer;\n", | |
" color: var(--xr-font-color2);\n", | |
"}\n", | |
"\n", | |
".xr-section-item input:focus + label {\n", | |
" border: 2px solid var(--xr-font-color0);\n", | |
"}\n", | |
"\n", | |
".xr-section-item input:enabled + label:hover {\n", | |
" color: var(--xr-font-color0);\n", | |
"}\n", | |
"\n", | |
".xr-section-summary {\n", | |
" grid-column: 1;\n", | |
" color: var(--xr-font-color2);\n", | |
" font-weight: 500;\n", | |
"}\n", | |
"\n", | |
".xr-section-summary > span {\n", | |
" display: inline-block;\n", | |
" padding-left: 0.5em;\n", | |
"}\n", | |
"\n", | |
".xr-section-summary-in:disabled + label {\n", | |
" color: var(--xr-font-color2);\n", | |
"}\n", | |
"\n", | |
".xr-section-summary-in + label:before {\n", | |
" display: inline-block;\n", | |
" content: \"►\";\n", | |
" font-size: 11px;\n", | |
" width: 15px;\n", | |
" text-align: center;\n", | |
"}\n", | |
"\n", | |
".xr-section-summary-in:disabled + label:before {\n", | |
" color: var(--xr-disabled-color);\n", | |
"}\n", | |
"\n", | |
".xr-section-summary-in:checked + label:before {\n", | |
" content: \"▼\";\n", | |
"}\n", | |
"\n", | |
".xr-section-summary-in:checked + label > span {\n", | |
" display: none;\n", | |
"}\n", | |
"\n", | |
".xr-section-summary,\n", | |
".xr-section-inline-details {\n", | |
" padding-top: 4px;\n", | |
" padding-bottom: 4px;\n", | |
"}\n", | |
"\n", | |
".xr-section-inline-details {\n", | |
" grid-column: 2 / -1;\n", | |
"}\n", | |
"\n", | |
".xr-section-details {\n", | |
" display: none;\n", | |
" grid-column: 1 / -1;\n", | |
" margin-bottom: 5px;\n", | |
"}\n", | |
"\n", | |
".xr-section-summary-in:checked ~ .xr-section-details {\n", | |
" display: contents;\n", | |
"}\n", | |
"\n", | |
".xr-array-wrap {\n", | |
" grid-column: 1 / -1;\n", | |
" display: grid;\n", | |
" grid-template-columns: 20px auto;\n", | |
"}\n", | |
"\n", | |
".xr-array-wrap > label {\n", | |
" grid-column: 1;\n", | |
" vertical-align: top;\n", | |
"}\n", | |
"\n", | |
".xr-preview {\n", | |
" color: var(--xr-font-color3);\n", | |
"}\n", | |
"\n", | |
".xr-array-preview,\n", | |
".xr-array-data {\n", | |
" padding: 0 5px !important;\n", | |
" grid-column: 2;\n", | |
"}\n", | |
"\n", | |
".xr-array-data,\n", | |
".xr-array-in:checked ~ .xr-array-preview {\n", | |
" display: none;\n", | |
"}\n", | |
"\n", | |
".xr-array-in:checked ~ .xr-array-data,\n", | |
".xr-array-preview {\n", | |
" display: inline-block;\n", | |
"}\n", | |
"\n", | |
".xr-dim-list {\n", | |
" display: inline-block !important;\n", | |
" list-style: none;\n", | |
" padding: 0 !important;\n", | |
" margin: 0;\n", | |
"}\n", | |
"\n", | |
".xr-dim-list li {\n", | |
" display: inline-block;\n", | |
" padding: 0;\n", | |
" margin: 0;\n", | |
"}\n", | |
"\n", | |
".xr-dim-list:before {\n", | |
" content: \"(\";\n", | |
"}\n", | |
"\n", | |
".xr-dim-list:after {\n", | |
" content: \")\";\n", | |
"}\n", | |
"\n", | |
".xr-dim-list li:not(:last-child):after {\n", | |
" content: \",\";\n", | |
" padding-right: 5px;\n", | |
"}\n", | |
"\n", | |
".xr-has-index {\n", | |
" font-weight: bold;\n", | |
"}\n", | |
"\n", | |
".xr-var-list,\n", | |
".xr-var-item {\n", | |
" display: contents;\n", | |
"}\n", | |
"\n", | |
".xr-var-item > div,\n", | |
".xr-var-item label,\n", | |
".xr-var-item > .xr-var-name span {\n", | |
" background-color: var(--xr-background-color-row-even);\n", | |
" margin-bottom: 0;\n", | |
"}\n", | |
"\n", | |
".xr-var-item > .xr-var-name:hover span {\n", | |
" padding-right: 5px;\n", | |
"}\n", | |
"\n", | |
".xr-var-list > li:nth-child(odd) > div,\n", | |
".xr-var-list > li:nth-child(odd) > label,\n", | |
".xr-var-list > li:nth-child(odd) > .xr-var-name span {\n", | |
" background-color: var(--xr-background-color-row-odd);\n", | |
"}\n", | |
"\n", | |
".xr-var-name {\n", | |
" grid-column: 1;\n", | |
"}\n", | |
"\n", | |
".xr-var-dims {\n", | |
" grid-column: 2;\n", | |
"}\n", | |
"\n", | |
".xr-var-dtype {\n", | |
" grid-column: 3;\n", | |
" text-align: right;\n", | |
" color: var(--xr-font-color2);\n", | |
"}\n", | |
"\n", | |
".xr-var-preview {\n", | |
" grid-column: 4;\n", | |
"}\n", | |
"\n", | |
".xr-index-preview {\n", | |
" grid-column: 2 / 5;\n", | |
" color: var(--xr-font-color2);\n", | |
"}\n", | |
"\n", | |
".xr-var-name,\n", | |
".xr-var-dims,\n", | |
".xr-var-dtype,\n", | |
".xr-preview,\n", | |
".xr-attrs dt {\n", | |
" white-space: nowrap;\n", | |
" overflow: hidden;\n", | |
" text-overflow: ellipsis;\n", | |
" padding-right: 10px;\n", | |
"}\n", | |
"\n", | |
".xr-var-name:hover,\n", | |
".xr-var-dims:hover,\n", | |
".xr-var-dtype:hover,\n", | |
".xr-attrs dt:hover {\n", | |
" overflow: visible;\n", | |
" width: auto;\n", | |
" z-index: 1;\n", | |
"}\n", | |
"\n", | |
".xr-var-attrs,\n", | |
".xr-var-data,\n", | |
".xr-index-data {\n", | |
" display: none;\n", | |
" background-color: var(--xr-background-color) !important;\n", | |
" padding-bottom: 5px !important;\n", | |
"}\n", | |
"\n", | |
".xr-var-attrs-in:checked ~ .xr-var-attrs,\n", | |
".xr-var-data-in:checked ~ .xr-var-data,\n", | |
".xr-index-data-in:checked ~ .xr-index-data {\n", | |
" display: block;\n", | |
"}\n", | |
"\n", | |
".xr-var-data > table {\n", | |
" float: right;\n", | |
"}\n", | |
"\n", | |
".xr-var-name span,\n", | |
".xr-var-data,\n", | |
".xr-index-name div,\n", | |
".xr-index-data,\n", | |
".xr-attrs {\n", | |
" padding-left: 25px !important;\n", | |
"}\n", | |
"\n", | |
".xr-attrs,\n", | |
".xr-var-attrs,\n", | |
".xr-var-data,\n", | |
".xr-index-data {\n", | |
" grid-column: 1 / -1;\n", | |
"}\n", | |
"\n", | |
"dl.xr-attrs {\n", | |
" padding: 0;\n", | |
" margin: 0;\n", | |
" display: grid;\n", | |
" grid-template-columns: 125px auto;\n", | |
"}\n", | |
"\n", | |
".xr-attrs dt,\n", | |
".xr-attrs dd {\n", | |
" padding: 0;\n", | |
" margin: 0;\n", | |
" float: left;\n", | |
" padding-right: 10px;\n", | |
" width: auto;\n", | |
"}\n", | |
"\n", | |
".xr-attrs dt {\n", | |
" font-weight: normal;\n", | |
" grid-column: 1;\n", | |
"}\n", | |
"\n", | |
".xr-attrs dt:hover span {\n", | |
" display: inline-block;\n", | |
" background: var(--xr-background-color);\n", | |
" padding-right: 10px;\n", | |
"}\n", | |
"\n", | |
".xr-attrs dd {\n", | |
" grid-column: 2;\n", | |
" white-space: pre-wrap;\n", | |
" word-break: break-all;\n", | |
"}\n", | |
"\n", | |
".xr-icon-database,\n", | |
".xr-icon-file-text2,\n", | |
".xr-no-icon {\n", | |
" display: inline-block;\n", | |
" vertical-align: middle;\n", | |
" width: 1em;\n", | |
" height: 1.5em !important;\n", | |
" stroke-width: 0;\n", | |
" stroke: currentColor;\n", | |
" fill: currentColor;\n", | |
"}\n", | |
"</style><pre class='xr-text-repr-fallback'><xarray.DataArray (cell_ids: 18)> Size: 144B\n", | |
"array([1344, 1345, 1347, 1348, 1349, 1350, 2492, 2493, 2495, 2530, 2531,\n", | |
" 2536, 2537, 2538, 2539, 2540, 2542, 2543])\n", | |
"Coordinates:\n", | |
" * cell_ids (cell_ids) int64 144B 1344 1345 1347 1348 ... 2539 2540 2542 2543\n", | |
"Indexes:\n", | |
" cell_ids HealpixIndex(level=4, indexing_scheme=nested)</pre><div class='xr-wrap' style='display:none'><div class='xr-header'><div class='xr-obj-type'>xarray.DataArray</div><div class='xr-array-name'></div><ul class='xr-dim-list'><li><span class='xr-has-index'>cell_ids</span>: 18</li></ul></div><ul class='xr-sections'><li class='xr-section-item'><div class='xr-array-wrap'><input id='section-10f12661-65ff-4865-a91b-db8336424d1d' class='xr-array-in' type='checkbox' checked><label for='section-10f12661-65ff-4865-a91b-db8336424d1d' title='Show/hide data repr'><svg class='icon xr-icon-database'><use xlink:href='#icon-database'></use></svg></label><div class='xr-array-preview xr-preview'><span>1344 1345 1347 1348 1349 1350 2492 ... 2537 2538 2539 2540 2542 2543</span></div><div class='xr-array-data'><pre>array([1344, 1345, 1347, 1348, 1349, 1350, 2492, 2493, 2495, 2530, 2531,\n", | |
" 2536, 2537, 2538, 2539, 2540, 2542, 2543])</pre></div></div></li><li class='xr-section-item'><input id='section-63c36c07-0ac2-4a32-914f-d5b7a8623d4b' class='xr-section-summary-in' type='checkbox' checked><label for='section-63c36c07-0ac2-4a32-914f-d5b7a8623d4b' class='xr-section-summary' >Coordinates: <span>(1)</span></label><div class='xr-section-inline-details'></div><div class='xr-section-details'><ul class='xr-var-list'><li class='xr-var-item'><div class='xr-var-name'><span class='xr-has-index'>cell_ids</span></div><div class='xr-var-dims'>(cell_ids)</div><div class='xr-var-dtype'>int64</div><div class='xr-var-preview xr-preview'>1344 1345 1347 ... 2540 2542 2543</div><input id='attrs-2795aa39-9539-4380-86c4-1dff8a54fdbf' class='xr-var-attrs-in' type='checkbox' ><label for='attrs-2795aa39-9539-4380-86c4-1dff8a54fdbf' title='Show/Hide attributes'><svg class='icon xr-icon-file-text2'><use xlink:href='#icon-file-text2'></use></svg></label><input id='data-b9b16a96-fd57-4595-8e06-85f9013ca393' class='xr-var-data-in' type='checkbox'><label for='data-b9b16a96-fd57-4595-8e06-85f9013ca393' title='Show/Hide data repr'><svg class='icon xr-icon-database'><use xlink:href='#icon-database'></use></svg></label><div class='xr-var-attrs'><dl class='xr-attrs'><dt><span>grid_name :</span></dt><dd>healpix</dd><dt><span>level :</span></dt><dd>4</dd><dt><span>indexing_scheme :</span></dt><dd>nested</dd></dl></div><div class='xr-var-data'><pre>array([1344, 1345, 1347, 1348, 1349, 1350, 2492, 2493, 2495, 2530, 2531, 2536,\n", | |
" 2537, 2538, 2539, 2540, 2542, 2543])</pre></div></li></ul></div></li><li class='xr-section-item'><input id='section-0f981f53-0358-4b95-96f7-4612dd0f987b' class='xr-section-summary-in' type='checkbox' ><label for='section-0f981f53-0358-4b95-96f7-4612dd0f987b' class='xr-section-summary' >Indexes: <span>(1)</span></label><div class='xr-section-inline-details'></div><div class='xr-section-details'><ul class='xr-var-list'><li class='xr-var-item'><div class='xr-index-name'><div>cell_ids</div></div><div class='xr-index-preview'>HealpixIndex(level=4, indexing_scheme=nested)</div><input type='checkbox' disabled/><label></label><input id='index-a5e5ff5a-ec4d-4b4a-980c-9586b2b3e49b' class='xr-index-data-in' type='checkbox'/><label for='index-a5e5ff5a-ec4d-4b4a-980c-9586b2b3e49b' title='Show/Hide index repr'><svg class='icon xr-icon-database'><use xlink:href='#icon-database'></use></svg></label><div class='xr-index-data'><pre><xdggs.healpix.HealpixIndex object at 0x7cbfabd01b20></pre></div></li></ul></div></li><li class='xr-section-item'><input id='section-fe767012-0e2f-4889-b706-96c4ec09cdc8' class='xr-section-summary-in' type='checkbox' disabled ><label for='section-fe767012-0e2f-4889-b706-96c4ec09cdc8' class='xr-section-summary' title='Expand/collapse section'>Attributes: <span>(0)</span></label><div class='xr-section-inline-details'></div><div class='xr-section-details'><dl class='xr-attrs'></dl></div></li></ul></div></div>" | |
], | |
"text/plain": [ | |
"<xarray.DataArray (cell_ids: 18)> Size: 144B\n", | |
"array([1344, 1345, 1347, 1348, 1349, 1350, 2492, 2493, 2495, 2530, 2531,\n", | |
" 2536, 2537, 2538, 2539, 2540, 2542, 2543])\n", | |
"Coordinates:\n", | |
" * cell_ids (cell_ids) int64 144B 1344 1345 1347 1348 ... 2539 2540 2542 2543\n", | |
"Indexes:\n", | |
" cell_ids HealpixIndex(level=4, indexing_scheme=nested)" | |
] | |
}, | |
"execution_count": 16, | |
"metadata": {}, | |
"output_type": "execute_result" | |
} | |
], | |
"source": [ | |
"da1" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 5, | |
"id": "f23d4d3b-6d84-4647-95f3-993a68c4cfb8", | |
"metadata": {}, | |
"outputs": [], | |
"source": [ | |
"import lonboard\n", | |
"import geoarrow.rust.core as geoarrow\n", | |
"from arro3.core import Array, Schema, Table" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 6, | |
"id": "924b5942-68a2-478a-97b8-b969e8c5fcf3", | |
"metadata": {}, | |
"outputs": [], | |
"source": [ | |
"def add_polygon(map, polygon):\n", | |
" arr = Array.from_arrow(geoarrow.from_shapely(np.array([polygon])))\n", | |
" arrays = {\"geometry\": arr}\n", | |
" fields = [array.field.with_name(name) for name, array in arrays.items()]\n", | |
" schema = Schema(fields)\n", | |
"\n", | |
" table = Table.from_arrays(list(arrays.values()), schema=schema)\n", | |
" layer = lonboard.PolygonLayer(table=table, get_line_width=2, get_line_color=\"red\", get_fill_color=[0, 0, 0, 180])\n", | |
"\n", | |
" map.layers = (layer, *map.layers)\n", | |
" return map\n", | |
" " | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 7, | |
"id": "c7488da0-50b4-4938-b6fe-57bffd45ab4a", | |
"metadata": {}, | |
"outputs": [ | |
{ | |
"name": "stderr", | |
"output_type": "stream", | |
"text": [ | |
"/home/jmagin/.local/opt/mambaforge/envs/cdshealpix/lib/python3.12/site-packages/lonboard/_geoarrow/ops/reproject.py:33: UserWarning: No CRS exists on data. If no data is shown on the map, double check that your CRS is WGS84.\n", | |
" warn(\n" | |
] | |
}, | |
{ | |
"data": { | |
"application/vnd.jupyter.widget-view+json": { | |
"model_id": "77ad6fbbb45e4690b46220f6af75440c", | |
"version_major": 2, | |
"version_minor": 1 | |
}, | |
"text/plain": [ | |
"Map(custom_attribution='', layers=(PolygonLayer(get_fill_color=[0, 0, 0, 180], get_line_color=(255, 0, 0, 255)…" | |
] | |
}, | |
"execution_count": 7, | |
"metadata": {}, | |
"output_type": "execute_result" | |
} | |
], | |
"source": [ | |
"map = da_truth.dggs.explore(alpha=0.7)\n", | |
"add_polygon(map, polygon)" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 8, | |
"id": "3", | |
"metadata": {}, | |
"outputs": [ | |
{ | |
"name": "stderr", | |
"output_type": "stream", | |
"text": [ | |
"/home/jmagin/.local/opt/mambaforge/envs/cdshealpix/lib/python3.12/site-packages/lonboard/_geoarrow/ops/reproject.py:33: UserWarning: No CRS exists on data. If no data is shown on the map, double check that your CRS is WGS84.\n", | |
" warn(\n" | |
] | |
}, | |
{ | |
"data": { | |
"application/vnd.jupyter.widget-view+json": { | |
"model_id": "5b8e2bb7cc5144d2bf2a5fa1806bd224", | |
"version_major": 2, | |
"version_minor": 1 | |
}, | |
"text/plain": [ | |
"Map(custom_attribution='', layers=(PolygonLayer(get_fill_color=[0, 0, 0, 180], get_line_color=(255, 0, 0, 255)…" | |
] | |
}, | |
"execution_count": 8, | |
"metadata": {}, | |
"output_type": "execute_result" | |
} | |
], | |
"source": [ | |
"map = da1.dggs.explore()\n", | |
"add_polygon(map, polygon)" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 9, | |
"id": "4a45abf1-98d6-45c2-b806-cec1165d5763", | |
"metadata": {}, | |
"outputs": [], | |
"source": [ | |
"import healpix_geo" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 18, | |
"id": "39ddafb9-a86c-43bd-80a9-409ee2d58272", | |
"metadata": {}, | |
"outputs": [ | |
{ | |
"data": { | |
"text/plain": [ | |
"array([[153.56946903, -43.63459726],\n", | |
" [153.56946903, -43.53469905],\n", | |
" [153.56946903, -43.43480083],\n", | |
" ...,\n", | |
" [153.36981386, -43.63459726],\n", | |
" [153.46964145, -43.63459726],\n", | |
" [153.56946903, -43.63459726]], shape=(1467, 2))" | |
] | |
}, | |
"execution_count": 18, | |
"metadata": {}, | |
"output_type": "execute_result" | |
} | |
], | |
"source": [ | |
"POLYGON" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 22, | |
"id": "b564072a-0417-4913-b129-a167e2cc96c9", | |
"metadata": {}, | |
"outputs": [ | |
{ | |
"data": { | |
"text/plain": [ | |
"dtype('int64')" | |
] | |
}, | |
"execution_count": 22, | |
"metadata": {}, | |
"output_type": "execute_result" | |
} | |
], | |
"source": [ | |
"da[\"cell_ids\"].data" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 24, | |
"id": "74111a36-77ef-41ff-a566-f025fad3f8ca", | |
"metadata": {}, | |
"outputs": [ | |
{ | |
"data": { | |
"text/plain": [ | |
"array([[153.56946903, -43.63459726],\n", | |
" [153.56946903, -43.53469905],\n", | |
" [153.56946903, -43.43480083],\n", | |
" ...,\n", | |
" [153.36981386, -43.63459726],\n", | |
" [153.46964145, -43.63459726],\n", | |
" [153.56946903, -43.63459726]], shape=(1467, 2))" | |
] | |
}, | |
"execution_count": 24, | |
"metadata": {}, | |
"output_type": "execute_result" | |
} | |
], | |
"source": [ | |
"POLYGON" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 26, | |
"id": "8043174d-3bff-4a3a-b825-9cb8463f9daa", | |
"metadata": {}, | |
"outputs": [], | |
"source": [ | |
"box = shapely.box(*get_aus_bbox())\n", | |
"polygon_ = shapely.get_coordinates(box)" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 28, | |
"id": "b0c25090-725d-444a-98cd-71a6ee58c097", | |
"metadata": {}, | |
"outputs": [ | |
{ | |
"data": { | |
"text/plain": [ | |
"array([[153.56946903, -43.63459726],\n", | |
" [153.56946903, -10.66818572],\n", | |
" [113.33895308, -10.66818572],\n", | |
" [113.33895308, -43.63459726],\n", | |
" [153.56946903, -43.63459726]])" | |
] | |
}, | |
"execution_count": 28, | |
"metadata": {}, | |
"output_type": "execute_result" | |
} | |
], | |
"source": [ | |
"polygon_" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 29, | |
"id": "20b82cd5-3380-48bf-9dfa-a02cb638733a", | |
"metadata": {}, | |
"outputs": [ | |
{ | |
"name": "stderr", | |
"output_type": "stream", | |
"text": [ | |
"\n", | |
"thread '<unnamed>' panicked at /home/jmagin/.cargo/registry/src/index.crates.io-1949cf8c6b5b557f/cdshealpix-0.7.3/src/special_points_finder.rs:459:3:\n", | |
"p1: (153.56946902900003, -43.6345972634); p2: (153.56946902900003, -41.810314895778596)\n" | |
] | |
}, | |
{ | |
"ename": "PanicException", | |
"evalue": "p1: (153.56946902900003, -43.6345972634); p2: (153.56946902900003, -41.810314895778596)", | |
"output_type": "error", | |
"traceback": [ | |
"\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", | |
"\u001b[0;31mPanicException\u001b[0m Traceback (most recent call last)", | |
"Cell \u001b[0;32mIn[29], line 1\u001b[0m\n\u001b[0;32m----> 1\u001b[0m \u001b[43mhealpix_geo\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mselect\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mcells_in_polygon\u001b[49m\u001b[43m(\u001b[49m\u001b[43mda\u001b[49m\u001b[43m[\u001b[49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[38;5;124;43mcell_ids\u001b[39;49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[43m]\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mdata\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mastype\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[38;5;124;43muint64\u001b[39;49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[43m)\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mda\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mdggs\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mgrid_info\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mlevel\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mpolygon_\u001b[49m\u001b[43m[\u001b[49m\u001b[43m:\u001b[49m\u001b[38;5;241;43m-\u001b[39;49m\u001b[38;5;241;43m1\u001b[39;49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43m:\u001b[49m\u001b[43m:\u001b[49m\u001b[38;5;241;43m-\u001b[39;49m\u001b[38;5;241;43m1\u001b[39;49m\u001b[43m]\u001b[49m\u001b[43m)\u001b[49m\n", | |
"\u001b[0;31mPanicException\u001b[0m: p1: (153.56946902900003, -43.6345972634); p2: (153.56946902900003, -41.810314895778596)" | |
] | |
} | |
], | |
"source": [ | |
"healpix_geo.select.cells_in_polygon(da[\"cell_ids\"].data.astype(\"uint64\"), da.dggs.grid_info.level, polygon_[:-1, ::-1])" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": null, | |
"id": "a6fd9a1a-f035-4bc6-ad95-a95437c19f46", | |
"metadata": {}, | |
"outputs": [], | |
"source": [] | |
} | |
], | |
"metadata": { | |
"kernelspec": { | |
"display_name": "Python 3 (ipykernel)", | |
"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.12.8" | |
} | |
}, | |
"nbformat": 4, | |
"nbformat_minor": 5 | |
} |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment