Spaces:
Runtime error
Runtime error
File size: 8,191 Bytes
c569dd7 186b0eb c569dd7 186b0eb c569dd7 186b0eb c569dd7 186b0eb c569dd7 29b5892 c569dd7 186b0eb 29b5892 186b0eb 29b5892 186b0eb c569dd7 29b5892 186b0eb 29b5892 c569dd7 29b5892 47d59cb 186b0eb 47d59cb 186b0eb 0b6e435 47d59cb 29b5892 186b0eb 29b5892 c569dd7 0b6e435 29b5892 0b6e435 29b5892 c569dd7 186b0eb c569dd7 0b6e435 c569dd7 29b5892 c569dd7 |
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 |
{
"cells": [
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"import leafmap\n",
"import solara\n",
"import pystac_client\n",
"import planetary_computer\n",
"import odc.stac\n",
"import geopandas as gpd\n",
"import dask.distributed\n",
"import matplotlib.pyplot as plt\n",
"import rioxarray\n",
"from osgeo import gdal\n"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"\n",
"# Stashed public copies of NPS polygons and CalFire polygons\n",
"\n",
"\n",
"zoom = solara.reactive(14)\n",
"center = solara.reactive((34, -116))\n",
"nps = gpd.read_file(\"/vsicurl/https://huggingface.co/datasets/cboettig/biodiversity/resolve/main/data/NPS.gdb\")\n",
"calfire = gpd.read_file(\"/vsicurl/https://huggingface.co/datasets/cboettig/biodiversity/resolve/main/data/fire22_1.gdb\", layer = \"firep22_1\")\n"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"\n",
"\n",
"\n",
"# fire = gpd.read_file(\"/vsizip/vsicurl/https://edcintl.cr.usgs.gov/downloads/sciweb1/shared/MTBS_Fire/data/composite_data/burned_area_extent_shapefile/mtbs_perimeter_data.zip\"\n",
"\n",
"# extract and reproject the Joshua Tree NP Polygon\n",
"jtree = nps[nps.PARKNAME == \"Joshua Tree\"].to_crs(calfire.crs)\n",
"# All Fires in the DB that intersect the Park\n",
"jtree_fires = jtree.overlay(calfire, how=\"intersection\")\n",
"\n",
"# Extract a polygon if interest. > 2015 for Sentinel, otherwise we can use LandSat\n",
"recent = jtree_fires[jtree_fires.YEAR_ > \"2015\"]\n",
"big = recent[recent.Shape_Area == recent.Shape_Area.max()].to_crs(\"EPSG:4326\")\n",
"box = big.buffer(0.01).bounds.to_numpy()[0] # Fire bbox + buffer\n",
"#box = jtree.to_crs(\"EPSG:4326\").bounds.to_numpy()[0] # Park bbox\n"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [],
"source": [
"from datetime import datetime, timedelta\n",
"alarm_date = datetime.strptime(big.ALARM_DATE.item(), \"%Y-%m-%dT%H:%M:%S+00:00\") \n",
"before_date = alarm_date - timedelta(days=14)\n",
"after_date = alarm_date + timedelta(days=14)"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [],
"source": [
"search_dates = big.ALARM_DATE.item() + \"/\" + big.CONT_DATE.item()\n",
"search_dates = before_date.strftime(\"%Y-%m-%d\") + \"/\" + after_date.strftime(\"%Y-%m-%d\")"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [],
"source": [
"\n",
"def stac_search(box, datetime): \n",
" # STAC Search for this imagery in space/time window\n",
" items = (\n",
" pystac_client.Client.\n",
" open(\"https://planetarycomputer.microsoft.com/api/stac/v1\",\n",
" modifier=planetary_computer.sign_inplace).\n",
" search(collections=[\"sentinel-2-l2a\"],\n",
" bbox=box,\n",
" datetime=datetime,\n",
" query={\"eo:cloud_cover\": {\"lt\": 10}}).\n",
" item_collection())\n",
" return items\n",
"\n",
"def compute_nbs(items, box):\n",
" # Time to compute:\n",
" client = dask.distributed.Client()\n",
" # landsat_bands = [\"nir08\", \"swir16\"]\n",
" sentinel_bands = [\"B08\", \"B12\", \"SCL\"] # NIR, SWIR, and Cloud Mask\n",
"\n",
" # The magic of gdalwarper. Can also resample, reproject, and aggregate on the fly\n",
" data = odc.stac.load(items,\n",
" bands=sentinel_bands,\n",
" bbox=box\n",
" )\n",
" # Compute the Normalized Burn Ratio, must be float\n",
" swir = data[\"B12\"].astype(\"float\")\n",
" nir = data[\"B08\"].astype(\"float\")\n",
" # can resample and aggregate in xarray. compute with dask\n",
" nbs = (((nir - swir) / (nir + swir)).\n",
" # resample(time=\"MS\").\n",
" # median(\"time\", keep_attrs=True).\n",
" compute()\n",
" )\n",
" return nbs\n",
"\n",
"items = stac_search(box, search_dates)\n",
"nbs = compute_nbs(items, box)\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [],
"source": [
"\n",
"nbs.isel(time=0).rio.to_raster(raster_path=\"before.tif\", driver=\"COG\")\n",
"nbs.isel(time=(nbs.time.size-1)).rio.to_raster(raster_path=\"after.tif\", driver=\"COG\")\n"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [],
"source": [
"\n",
"before_url = \"https://huggingface.co/datasets/cboettig/solara-data/resolve/main/before.tif\"\n",
"after_url = \"https://huggingface.co/datasets/cboettig/solara-data/resolve/main/after.tif\"\n",
"\n",
"style = {\n",
" \"stroke\": False,\n",
" \"fill\": True,\n",
" \"fillColor\": \"#ff6666\",\n",
" \"fillOpacity\": 0.5,\n",
"}"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [],
"source": [
"class Map(leafmap.Map):\n",
" def __init__(self, **kwargs):\n",
" super().__init__(**kwargs)\n",
" # Add what you want below\n",
" # self.add_gdf(jtree, layer_name = \"Joshua Tree NP\")\n",
" self.add_gdf(jtree_fires, layer_name = \"All Fires\", style=style)\n",
" self.add_gdf(big, layer_name = big.FIRE_NAME.item())\n",
" #self.add_raster(\"before.tif\", layer_name = \"before\", colormap=\"viridis\")\n",
" #self.add_raster(\"after.tif\", layer_name = \"after\", colormap=\"viridis\")\n",
" self.split_map(before_url, after_url, \n",
" left_label= \"before fire\", \n",
" right_label = \"after fire\")\n",
" #self.add_stac_gui()\n",
"\n",
"\n",
"@solara.component\n",
"def Page():\n",
" with solara.Column(style={\"min-width\": \"500px\"}):\n",
" # solara components support reactive variables\n",
" # solara.SliderInt(label=\"Zoom level\", value=zoom, min=1, max=20)\n",
" # using 3rd party widget library require wiring up the events manually\n",
" # using zoom.value and zoom.set\n",
" Map.element( # type: ignore\n",
" zoom=zoom.value,\n",
" on_zoom=zoom.set,\n",
" center=center.value,\n",
" on_center=center.set,\n",
" scroll_wheel_zoom=True,\n",
" toolbar_ctrl=False,\n",
" data_ctrl=False,\n",
" height=\"780px\",\n",
" )\n",
" solara.Text(f\"Zoom: {zoom.value}\")\n",
" solara.Text(f\"Center: {center.value}\")\n"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {},
"outputs": [
{
"data": {
"application/vnd.jupyter.widget-view+json": {
"model_id": "822080e3202b494888a05630fb5bf798",
"version_major": 2,
"version_minor": 0
},
"text/html": [
"Cannot show widget. You probably want to rerun the code cell above (<i>Click in the code cell, and press Shift+Enter <kbd>⇧</kbd>+<kbd>↩</kbd></i>)."
],
"text/plain": [
"Cannot show ipywidgets in text"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"Page()"
]
}
],
"metadata": {
"kernelspec": {
"display_name": ".venv",
"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.6"
}
},
"nbformat": 4,
"nbformat_minor": 2
}
|