From faccfaaaa50021106469590434c4863b9e7df193 Mon Sep 17 00:00:00 2001 From: Carlotta Date: Wed, 5 Feb 2025 15:30:38 +0100 Subject: [PATCH 1/5] updated history in netCDF files --- Untitled1.ipynb | 1143 ++++++++++++++++++++++++++++++++++++++++++++++ mapies/viirs.py | 52 ++- run/run_viirs.py | 10 +- 3 files changed, 1190 insertions(+), 15 deletions(-) create mode 100644 Untitled1.ipynb diff --git a/Untitled1.ipynb b/Untitled1.ipynb new file mode 100644 index 00000000..e9f0d743 --- /dev/null +++ b/Untitled1.ipynb @@ -0,0 +1,1143 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": 9, + "id": "2eea51dd-3297-4f9a-b9f3-f171190563e0", + "metadata": {}, + "outputs": [], + "source": [ + "import warnings\n", + "warnings.filterwarnings(\"ignore\")\n", + "\n", + "from functools import lru_cache\n", + "from pathlib import Path\n", + "from numpy import array\n", + "from typing import List\n", + "from pandas import DataFrame, Timestamp, concat\n", + "from tqdm.auto import tqdm\n", + "import xarray as xr\n", + "\n", + "\n", + "# Get the file list using a cached function. Maybe overkill, but listing *all* the\n", + "# files on disk at each timestep might be heavy on the file system, if there are lots of files?\n", + "@lru_cache\n", + "def get_file_list(path: Path, pattern: str):\n", + " return array(list((path.glob('*/*/'+pattern))))\n", + "\n", + "\n", + "def read_obs(\n", + " files: List[Path],\n", + " start: Timestamp | None = None,\n", + " end: Timestamp | None = None,\n", + " varname: str = 'Aerosol_Optical_Thickness_550_Land_Ocean_Best_Estimate',\n", + " tref: Timestamp = Timestamp(1993, 1, 1), # TAI93\n", + " progress: bool = False\n", + ") -> DataFrame:\n", + " \"\"\"\n", + " Read a bunch of files, extract the variable given by the `varname` argument, and return the data as a single dataframe.\n", + " The optional `start` and `end` arguments can be used to limit the temporal extent of the data\n", + " \"\"\"\n", + " dfs = []\n", + " for filename in tqdm(files, disable=not progress, leave=False):\n", + " #logger.debug(f'Read data from {filename}')\n", + " ds = xr.open_dataset(filename)\n", + " df = ds[[varname, 'Scan_Start_Time', 'Aerosol_Type_Land_Ocean', 'Aerosol_Optical_Thickness_QA_Flag_Ocean', 'Aerosol_Optical_Thickness_QA_Flag_Land']].to_dataframe().dropna(subset=varname).reset_index()\n", + " df.loc[:, 'time'] = (tref + df.Scan_Start_Time)\n", + " dfs.append(df.rename(columns={varname: 'obs', 'Longitude': 'lon', 'Latitude': 'lat'}))\n", + " df = concat(dfs)\n", + " if start is not None:\n", + " df = df.loc[df.time >= start]\n", + " if end is not None:\n", + " df = df.loc[df.time < end]\n", + " return df" + ] + }, + { + "cell_type": "code", + "execution_count": 10, + "id": "071b8d5e-4e24-4948-bef5-8f58b017704c", + "metadata": {}, + "outputs": [], + "source": [ + "from numpy import ndarray\n", + "from typing import Tuple\n", + "from pint import Quantity\n", + "from functools import partial\n", + "from numpy import cos, sin, arctan, pi\n", + "\n", + "\n", + "def geo_to_rot(lats: ndarray, lons: ndarray, center_lat: float, center_lon: float) -> Tuple[array, array]:\n", + "\n", + " # Convert to pint quantities, so we don't have to handle the unit conversions\n", + " lats = Quantity(lats, 'deg')\n", + " lons = Quantity(lons, 'deg')\n", + " center_lat = Quantity(center_lat, 'deg')\n", + " center_lon = Quantity(center_lon, 'deg')\n", + "\n", + " # Compute the rotation\n", + " x = cos(center_lat) * cos(lats) * cos(lons - center_lon) + sin(center_lat) * sin(lats)\n", + " y = cos(lats) * sin(lons - center_lon)\n", + " z = - sin(center_lat) * cos(lats) * cos(lons - center_lon) + cos(center_lat) * sin(lats)\n", + " rlat = arctan(z / ((x ** 2 + y ** 2) ** .5))\n", + " rlon = arctan(y / x)\n", + " rlon[x < 0] += pi\n", + "\n", + " # Return as array, in degrees\n", + " return rlat.to('deg').m, rlon.to('deg').m" + ] + }, + { + "cell_type": "code", + "execution_count": 11, + "id": "2e9e801c-372f-4af5-bb93-1f60f8f907e4", + "metadata": {}, + "outputs": [], + "source": [ + "from pandas import Interval, DatetimeIndex\n", + "\n", + "def import_one_chunk_of_data(interval: Interval, file_list: List[str], times: DatetimeIndex, grid):\n", + " \n", + " # Filter the files, adding a bit of safety margin:\n", + " files = file_list[(times >= interval.left - Timedelta(minutes=15)) & (times <= interval.right + Timedelta(minutes=30))]\n", + "\n", + " # Read all the obs of this chunk\n", + " df = read_obs(files, varname='Aerosol_Optical_Thickness_550_Land_Ocean_Best_Estimate', start=interval.left, end=interval.right, progress=False)\n", + "\n", + " # Calculate the rotated coordinates\n", + " df['rotated_lat'], df['rotated_lon'] = geo_to_rot(df.lat.values, df.lon.values, center_lat=grid.center_lat, center_lon=grid.center_lon)\n", + " \n", + " # Filter out the data outside the domain:\n", + " df = df.loc[(df.rotated_lat >= grid.south) & (df.rotated_lat <= -grid.south) & (df.rotated_lon >= grid.west) & (df.rotated_lon <= -grid.west)]\n", + " \n", + " # Calculate the grid indices\n", + " df['ilat'] = ((df['rotated_lat'] - grid.south) / grid.dlat).astype(int)\n", + " df['ilon'] = ((df['rotated_lon'] - grid.west) / grid.dlon).astype(int)\n", + " \n", + " # Aggregate:\n", + " df['nobs'] = 1.\n", + " nobs = df[['ilon', 'ilat', 'nobs']].groupby(['ilon', 'ilat']).sum().reset_index()\n", + " data = df.groupby(['ilon', 'ilat']).mean().reset_index()\n", + " data['nobs'] = nobs.nobs\n", + "\n", + " return data" + ] + }, + { + "cell_type": "code", + "execution_count": 12, + "id": "a47937e0-4c07-4773-a820-a054647a5d9f", + "metadata": {}, + "outputs": [ + { + "data": { + "application/vnd.jupyter.widget-view+json": { + "model_id": "14ef111af95f412581f35b3161f9a8f8", + "version_major": 2, + "version_minor": 0 + }, + "text/plain": [ + " 0%| | 0/31 [00:00 DataFrame:\n", + " chunk_processor = partial(import_one_chunk_of_data, file_list=file_list, times=times, grid=grid)\n", + "\n", + " with Pool() as pp:\n", + " all_data = concat(pp.imap(chunk_processor, tqdm(interval_range(start, end, freq=chunksize))))\n", + "\n", + " return all_data\n", + "\n", + "\n", + "def compute_monthly_data_sp(files, times, grid, start, end, chunksize) -> DataFrame:\n", + " all_data = []\n", + " for trange in tqdm(interval_range(start, end, freq=chunksize)):\n", + " all_data.append(import_one_chunk_of_data(interval=trange, file_list=files, times=times, grid=grid))\n", + " return concat(all_data)\n", + "\n", + "\n", + "# Settings\n", + "start = Timestamp(2024, 1, 1)\n", + "end = Timestamp(2024, 2, 1)\n", + "path_originals : Path = Path('/home/guillaume/BSC/mapies/run/data/VIIRS/original_files/AERDB_L2_VIIRS_NOAA20')\n", + "chunksize = Timedelta('1d')\n", + "grid = SimpleNamespace(center_lat=35, center_lon=15, dlat=.2, dlon=.2, south=-27, west=-35)\n", + "parallelize = True\n", + "\n", + "\n", + "# get the list of files on disk:\n", + "file_list = get_file_list(path_originals, pattern='AERDB_L2_VIIRS_NOAA20.*.nc')\n", + "\n", + "# deduce the corresponding obs times\n", + "times = [f.name.split('.')[1:3] for f in file_list]\n", + "times = array([datetime.strptime(d+h, 'A%Y%j%H%M') for (d, h) in times])\n", + "\n", + "# Choose which function we are going to use, depending on the \"parallelize\" flag\n", + "# compute_monthly_data = {True:compute_monthly_data_mp, False:compute_monthly_data_sp}[parallelize]\n", + "\n", + "all_data = compute_monthly_data_sp(file_list, times, grid, start, end, chunksize)\n", + "all_data = compute_monthly_data_mp(file_list, times, grid, start, end, chunksize)\n" + ] + }, + { + "cell_type": "code", + "execution_count": 17, + "id": "20d7fa95-471e-481c-9d8f-ae9c657b3e09", + "metadata": {}, + "outputs": [], + "source": [ + "# Plot the raw data (just one chunk because the rest doesn't fit in memory.\n", + "\n", + "# import hvplot.pandas\n", + "# import geoviews.feature as gf\n", + "# import geoviews as gv\n", + "# gv.extension('bokeh')\n", + "# df.hvplot.scatter(x='lon', y='lat', c='obs', s=1, global_extent=True, coastline=True, title='raw obs')" + ] + }, + { + "cell_type": "code", + "execution_count": 16, + "id": "50b7c205-052d-4034-b138-0b069d2dac02", + "metadata": {}, + "outputs": [ + { + "data": { + "text/html": [ + "" + ] + }, + "metadata": {}, + "output_type": "display_data" + }, + { + "data": { + "application/javascript": [ + "(function(root) {\n", + " function now() {\n", + " return new Date();\n", + " }\n", + "\n", + " const force = true;\n", + " const py_version = '3.6.2'.replace('rc', '-rc.').replace('.dev', '-dev.');\n", + " const reloading = false;\n", + " const Bokeh = root.Bokeh;\n", + "\n", + " // Set a timeout for this load but only if we are not already initializing\n", + " if (typeof (root._bokeh_timeout) === \"undefined\" || (force || !root._bokeh_is_initializing)) {\n", + " root._bokeh_timeout = Date.now() + 5000;\n", + " root._bokeh_failed_load = false;\n", + " }\n", + "\n", + " function run_callbacks() {\n", + " try {\n", + " root._bokeh_onload_callbacks.forEach(function(callback) {\n", + " if (callback != null)\n", + " callback();\n", + " });\n", + " } finally {\n", + " delete root._bokeh_onload_callbacks;\n", + " }\n", + " console.debug(\"Bokeh: all callbacks have finished\");\n", + " }\n", + "\n", + " function load_libs(css_urls, js_urls, js_modules, js_exports, callback) {\n", + " if (css_urls == null) css_urls = [];\n", + " if (js_urls == null) js_urls = [];\n", + " if (js_modules == null) js_modules = [];\n", + " if (js_exports == null) js_exports = {};\n", + "\n", + " root._bokeh_onload_callbacks.push(callback);\n", + "\n", + " if (root._bokeh_is_loading > 0) {\n", + " // Don't load bokeh if it is still initializing\n", + " console.debug(\"Bokeh: BokehJS is being loaded, scheduling callback at\", now());\n", + " return null;\n", + " } else if (js_urls.length === 0 && js_modules.length === 0 && Object.keys(js_exports).length === 0) {\n", + " // There is nothing to load\n", + " run_callbacks();\n", + " return null;\n", + " }\n", + "\n", + " function on_load() {\n", + " root._bokeh_is_loading--;\n", + " if (root._bokeh_is_loading === 0) {\n", + " console.debug(\"Bokeh: all BokehJS libraries/stylesheets loaded\");\n", + " run_callbacks()\n", + " }\n", + " }\n", + " window._bokeh_on_load = on_load\n", + "\n", + " function on_error(e) {\n", + " const src_el = e.srcElement\n", + " console.error(\"failed to load \" + (src_el.href || src_el.src));\n", + " }\n", + "\n", + " const skip = [];\n", + " if (window.requirejs) {\n", + " window.requirejs.config({'packages': {}, 'paths': {}, 'shim': {}});\n", + " root._bokeh_is_loading = css_urls.length + 0;\n", + " } else {\n", + " root._bokeh_is_loading = css_urls.length + js_urls.length + js_modules.length + Object.keys(js_exports).length;\n", + " }\n", + "\n", + " const existing_stylesheets = []\n", + " const links = document.getElementsByTagName('link')\n", + " for (let i = 0; i < links.length; i++) {\n", + " const link = links[i]\n", + " if (link.href != null) {\n", + " existing_stylesheets.push(link.href)\n", + " }\n", + " }\n", + " for (let i = 0; i < css_urls.length; i++) {\n", + " const url = css_urls[i];\n", + " const escaped = encodeURI(url)\n", + " if (existing_stylesheets.indexOf(escaped) !== -1) {\n", + " on_load()\n", + " continue;\n", + " }\n", + " const element = document.createElement(\"link\");\n", + " element.onload = on_load;\n", + " element.onerror = on_error;\n", + " element.rel = \"stylesheet\";\n", + " element.type = \"text/css\";\n", + " element.href = url;\n", + " console.debug(\"Bokeh: injecting link tag for BokehJS stylesheet: \", url);\n", + " document.body.appendChild(element);\n", + " } var existing_scripts = []\n", + " const scripts = document.getElementsByTagName('script')\n", + " for (let i = 0; i < scripts.length; i++) {\n", + " var script = scripts[i]\n", + " if (script.src != null) {\n", + " existing_scripts.push(script.src)\n", + " }\n", + " }\n", + " for (let i = 0; i < js_urls.length; i++) {\n", + " const url = js_urls[i];\n", + " const escaped = encodeURI(url)\n", + " if (skip.indexOf(escaped) !== -1 || existing_scripts.indexOf(escaped) !== -1) {\n", + " if (!window.requirejs) {\n", + " on_load();\n", + " }\n", + " continue;\n", + " }\n", + " const element = document.createElement('script');\n", + " element.onload = on_load;\n", + " element.onerror = on_error;\n", + " element.async = false;\n", + " element.src = url;\n", + " console.debug(\"Bokeh: injecting script tag for BokehJS library: \", url);\n", + " document.head.appendChild(element);\n", + " }\n", + " for (let i = 0; i < js_modules.length; i++) {\n", + " const url = js_modules[i];\n", + " const escaped = encodeURI(url)\n", + " if (skip.indexOf(escaped) !== -1 || existing_scripts.indexOf(escaped) !== -1) {\n", + " if (!window.requirejs) {\n", + " on_load();\n", + " }\n", + " continue;\n", + " }\n", + " var element = document.createElement('script');\n", + " element.onload = on_load;\n", + " element.onerror = on_error;\n", + " element.async = false;\n", + " element.src = url;\n", + " element.type = \"module\";\n", + " console.debug(\"Bokeh: injecting script tag for BokehJS library: \", url);\n", + " document.head.appendChild(element);\n", + " }\n", + " for (const name in js_exports) {\n", + " const url = js_exports[name];\n", + " const escaped = encodeURI(url)\n", + " if (skip.indexOf(escaped) >= 0 || root[name] != null) {\n", + " if (!window.requirejs) {\n", + " on_load();\n", + " }\n", + " continue;\n", + " }\n", + " var element = document.createElement('script');\n", + " element.onerror = on_error;\n", + " element.async = false;\n", + " element.type = \"module\";\n", + " console.debug(\"Bokeh: injecting script tag for BokehJS library: \", url);\n", + " element.textContent = `\n", + " import ${name} from \"${url}\"\n", + " window.${name} = ${name}\n", + " window._bokeh_on_load()\n", + " `\n", + " document.head.appendChild(element);\n", + " }\n", + " if (!js_urls.length && !js_modules.length) {\n", + " on_load()\n", + " }\n", + " };\n", + "\n", + " function inject_raw_css(css) {\n", + " const element = document.createElement(\"style\");\n", + " element.appendChild(document.createTextNode(css));\n", + " document.body.appendChild(element);\n", + " }\n", + "\n", + " const js_urls = [\"https://cdn.holoviz.org/panel/1.6.0/dist/bundled/reactiveesm/es-module-shims@^1.10.0/dist/es-module-shims.min.js\", \"https://cdn.bokeh.org/bokeh/release/bokeh-3.6.2.min.js\", \"https://cdn.bokeh.org/bokeh/release/bokeh-gl-3.6.2.min.js\", \"https://cdn.bokeh.org/bokeh/release/bokeh-widgets-3.6.2.min.js\", \"https://cdn.bokeh.org/bokeh/release/bokeh-tables-3.6.2.min.js\", \"https://cdn.holoviz.org/panel/1.6.0/dist/panel.min.js\"];\n", + " const js_modules = [];\n", + " const js_exports = {};\n", + " const css_urls = [];\n", + " const inline_js = [ function(Bokeh) {\n", + " Bokeh.set_log_level(\"info\");\n", + " },\n", + "function(Bokeh) {} // ensure no trailing comma for IE\n", + " ];\n", + "\n", + " function run_inline_js() {\n", + " if ((root.Bokeh !== undefined) || (force === true)) {\n", + " for (let i = 0; i < inline_js.length; i++) {\n", + " try {\n", + " inline_js[i].call(root, root.Bokeh);\n", + " } catch(e) {\n", + " if (!reloading) {\n", + " throw e;\n", + " }\n", + " }\n", + " }\n", + " // Cache old bokeh versions\n", + " if (Bokeh != undefined && !reloading) {\n", + " var NewBokeh = root.Bokeh;\n", + " if (Bokeh.versions === undefined) {\n", + " Bokeh.versions = new Map();\n", + " }\n", + " if (NewBokeh.version !== Bokeh.version) {\n", + " Bokeh.versions.set(NewBokeh.version, NewBokeh)\n", + " }\n", + " root.Bokeh = Bokeh;\n", + " }\n", + " } else if (Date.now() < root._bokeh_timeout) {\n", + " setTimeout(run_inline_js, 100);\n", + " } else if (!root._bokeh_failed_load) {\n", + " console.log(\"Bokeh: BokehJS failed to load within specified timeout.\");\n", + " root._bokeh_failed_load = true;\n", + " }\n", + " root._bokeh_is_initializing = false\n", + " }\n", + "\n", + " function load_or_wait() {\n", + " // Implement a backoff loop that tries to ensure we do not load multiple\n", + " // versions of Bokeh and its dependencies at the same time.\n", + " // In recent versions we use the root._bokeh_is_initializing flag\n", + " // to determine whether there is an ongoing attempt to initialize\n", + " // bokeh, however for backward compatibility we also try to ensure\n", + " // that we do not start loading a newer (Panel>=1.0 and Bokeh>3) version\n", + " // before older versions are fully initialized.\n", + " if (root._bokeh_is_initializing && Date.now() > root._bokeh_timeout) {\n", + " // If the timeout and bokeh was not successfully loaded we reset\n", + " // everything and try loading again\n", + " root._bokeh_timeout = Date.now() + 5000;\n", + " root._bokeh_is_initializing = false;\n", + " root._bokeh_onload_callbacks = undefined;\n", + " root._bokeh_is_loading = 0\n", + " console.log(\"Bokeh: BokehJS was loaded multiple times but one version failed to initialize.\");\n", + " load_or_wait();\n", + " } else if (root._bokeh_is_initializing || (typeof root._bokeh_is_initializing === \"undefined\" && root._bokeh_onload_callbacks !== undefined)) {\n", + " setTimeout(load_or_wait, 100);\n", + " } else {\n", + " root._bokeh_is_initializing = true\n", + " root._bokeh_onload_callbacks = []\n", + " const bokeh_loaded = root.Bokeh != null && (root.Bokeh.version === py_version || (root.Bokeh.versions !== undefined && root.Bokeh.versions.has(py_version)));\n", + " if (!reloading && !bokeh_loaded) {\n", + " if (root.Bokeh) {\n", + " root.Bokeh = undefined;\n", + " }\n", + " console.debug(\"Bokeh: BokehJS not loaded, scheduling load and callback at\", now());\n", + " }\n", + " load_libs(css_urls, js_urls, js_modules, js_exports, function() {\n", + " console.debug(\"Bokeh: BokehJS plotting callback run at\", now());\n", + " run_inline_js();\n", + " });\n", + " }\n", + " }\n", + " // Give older versions of the autoload script a head-start to ensure\n", + " // they initialize before we start loading newer version.\n", + " setTimeout(load_or_wait, 100)\n", + "}(window));" + ], + "application/vnd.holoviews_load.v0+json": "(function(root) {\n function now() {\n return new Date();\n }\n\n const force = true;\n const py_version = '3.6.2'.replace('rc', '-rc.').replace('.dev', '-dev.');\n const reloading = false;\n const Bokeh = root.Bokeh;\n\n // Set a timeout for this load but only if we are not already initializing\n if (typeof (root._bokeh_timeout) === \"undefined\" || (force || !root._bokeh_is_initializing)) {\n root._bokeh_timeout = Date.now() + 5000;\n root._bokeh_failed_load = false;\n }\n\n function run_callbacks() {\n try {\n root._bokeh_onload_callbacks.forEach(function(callback) {\n if (callback != null)\n callback();\n });\n } finally {\n delete root._bokeh_onload_callbacks;\n }\n console.debug(\"Bokeh: all callbacks have finished\");\n }\n\n function load_libs(css_urls, js_urls, js_modules, js_exports, callback) {\n if (css_urls == null) css_urls = [];\n if (js_urls == null) js_urls = [];\n if (js_modules == null) js_modules = [];\n if (js_exports == null) js_exports = {};\n\n root._bokeh_onload_callbacks.push(callback);\n\n if (root._bokeh_is_loading > 0) {\n // Don't load bokeh if it is still initializing\n console.debug(\"Bokeh: BokehJS is being loaded, scheduling callback at\", now());\n return null;\n } else if (js_urls.length === 0 && js_modules.length === 0 && Object.keys(js_exports).length === 0) {\n // There is nothing to load\n run_callbacks();\n return null;\n }\n\n function on_load() {\n root._bokeh_is_loading--;\n if (root._bokeh_is_loading === 0) {\n console.debug(\"Bokeh: all BokehJS libraries/stylesheets loaded\");\n run_callbacks()\n }\n }\n window._bokeh_on_load = on_load\n\n function on_error(e) {\n const src_el = e.srcElement\n console.error(\"failed to load \" + (src_el.href || src_el.src));\n }\n\n const skip = [];\n if (window.requirejs) {\n window.requirejs.config({'packages': {}, 'paths': {}, 'shim': {}});\n root._bokeh_is_loading = css_urls.length + 0;\n } else {\n root._bokeh_is_loading = css_urls.length + js_urls.length + js_modules.length + Object.keys(js_exports).length;\n }\n\n const existing_stylesheets = []\n const links = document.getElementsByTagName('link')\n for (let i = 0; i < links.length; i++) {\n const link = links[i]\n if (link.href != null) {\n existing_stylesheets.push(link.href)\n }\n }\n for (let i = 0; i < css_urls.length; i++) {\n const url = css_urls[i];\n const escaped = encodeURI(url)\n if (existing_stylesheets.indexOf(escaped) !== -1) {\n on_load()\n continue;\n }\n const element = document.createElement(\"link\");\n element.onload = on_load;\n element.onerror = on_error;\n element.rel = \"stylesheet\";\n element.type = \"text/css\";\n element.href = url;\n console.debug(\"Bokeh: injecting link tag for BokehJS stylesheet: \", url);\n document.body.appendChild(element);\n } var existing_scripts = []\n const scripts = document.getElementsByTagName('script')\n for (let i = 0; i < scripts.length; i++) {\n var script = scripts[i]\n if (script.src != null) {\n existing_scripts.push(script.src)\n }\n }\n for (let i = 0; i < js_urls.length; i++) {\n const url = js_urls[i];\n const escaped = encodeURI(url)\n if (skip.indexOf(escaped) !== -1 || existing_scripts.indexOf(escaped) !== -1) {\n if (!window.requirejs) {\n on_load();\n }\n continue;\n }\n const element = document.createElement('script');\n element.onload = on_load;\n element.onerror = on_error;\n element.async = false;\n element.src = url;\n console.debug(\"Bokeh: injecting script tag for BokehJS library: \", url);\n document.head.appendChild(element);\n }\n for (let i = 0; i < js_modules.length; i++) {\n const url = js_modules[i];\n const escaped = encodeURI(url)\n if (skip.indexOf(escaped) !== -1 || existing_scripts.indexOf(escaped) !== -1) {\n if (!window.requirejs) {\n on_load();\n }\n continue;\n }\n var element = document.createElement('script');\n element.onload = on_load;\n element.onerror = on_error;\n element.async = false;\n element.src = url;\n element.type = \"module\";\n console.debug(\"Bokeh: injecting script tag for BokehJS library: \", url);\n document.head.appendChild(element);\n }\n for (const name in js_exports) {\n const url = js_exports[name];\n const escaped = encodeURI(url)\n if (skip.indexOf(escaped) >= 0 || root[name] != null) {\n if (!window.requirejs) {\n on_load();\n }\n continue;\n }\n var element = document.createElement('script');\n element.onerror = on_error;\n element.async = false;\n element.type = \"module\";\n console.debug(\"Bokeh: injecting script tag for BokehJS library: \", url);\n element.textContent = `\n import ${name} from \"${url}\"\n window.${name} = ${name}\n window._bokeh_on_load()\n `\n document.head.appendChild(element);\n }\n if (!js_urls.length && !js_modules.length) {\n on_load()\n }\n };\n\n function inject_raw_css(css) {\n const element = document.createElement(\"style\");\n element.appendChild(document.createTextNode(css));\n document.body.appendChild(element);\n }\n\n const js_urls = [\"https://cdn.holoviz.org/panel/1.6.0/dist/bundled/reactiveesm/es-module-shims@^1.10.0/dist/es-module-shims.min.js\", \"https://cdn.bokeh.org/bokeh/release/bokeh-3.6.2.min.js\", \"https://cdn.bokeh.org/bokeh/release/bokeh-gl-3.6.2.min.js\", \"https://cdn.bokeh.org/bokeh/release/bokeh-widgets-3.6.2.min.js\", \"https://cdn.bokeh.org/bokeh/release/bokeh-tables-3.6.2.min.js\", \"https://cdn.holoviz.org/panel/1.6.0/dist/panel.min.js\"];\n const js_modules = [];\n const js_exports = {};\n const css_urls = [];\n const inline_js = [ function(Bokeh) {\n Bokeh.set_log_level(\"info\");\n },\nfunction(Bokeh) {} // ensure no trailing comma for IE\n ];\n\n function run_inline_js() {\n if ((root.Bokeh !== undefined) || (force === true)) {\n for (let i = 0; i < inline_js.length; i++) {\n try {\n inline_js[i].call(root, root.Bokeh);\n } catch(e) {\n if (!reloading) {\n throw e;\n }\n }\n }\n // Cache old bokeh versions\n if (Bokeh != undefined && !reloading) {\n var NewBokeh = root.Bokeh;\n if (Bokeh.versions === undefined) {\n Bokeh.versions = new Map();\n }\n if (NewBokeh.version !== Bokeh.version) {\n Bokeh.versions.set(NewBokeh.version, NewBokeh)\n }\n root.Bokeh = Bokeh;\n }\n } else if (Date.now() < root._bokeh_timeout) {\n setTimeout(run_inline_js, 100);\n } else if (!root._bokeh_failed_load) {\n console.log(\"Bokeh: BokehJS failed to load within specified timeout.\");\n root._bokeh_failed_load = true;\n }\n root._bokeh_is_initializing = false\n }\n\n function load_or_wait() {\n // Implement a backoff loop that tries to ensure we do not load multiple\n // versions of Bokeh and its dependencies at the same time.\n // In recent versions we use the root._bokeh_is_initializing flag\n // to determine whether there is an ongoing attempt to initialize\n // bokeh, however for backward compatibility we also try to ensure\n // that we do not start loading a newer (Panel>=1.0 and Bokeh>3) version\n // before older versions are fully initialized.\n if (root._bokeh_is_initializing && Date.now() > root._bokeh_timeout) {\n // If the timeout and bokeh was not successfully loaded we reset\n // everything and try loading again\n root._bokeh_timeout = Date.now() + 5000;\n root._bokeh_is_initializing = false;\n root._bokeh_onload_callbacks = undefined;\n root._bokeh_is_loading = 0\n console.log(\"Bokeh: BokehJS was loaded multiple times but one version failed to initialize.\");\n load_or_wait();\n } else if (root._bokeh_is_initializing || (typeof root._bokeh_is_initializing === \"undefined\" && root._bokeh_onload_callbacks !== undefined)) {\n setTimeout(load_or_wait, 100);\n } else {\n root._bokeh_is_initializing = true\n root._bokeh_onload_callbacks = []\n const bokeh_loaded = root.Bokeh != null && (root.Bokeh.version === py_version || (root.Bokeh.versions !== undefined && root.Bokeh.versions.has(py_version)));\n if (!reloading && !bokeh_loaded) {\n if (root.Bokeh) {\n root.Bokeh = undefined;\n }\n console.debug(\"Bokeh: BokehJS not loaded, scheduling load and callback at\", now());\n }\n load_libs(css_urls, js_urls, js_modules, js_exports, function() {\n console.debug(\"Bokeh: BokehJS plotting callback run at\", now());\n run_inline_js();\n });\n }\n }\n // Give older versions of the autoload script a head-start to ensure\n // they initialize before we start loading newer version.\n setTimeout(load_or_wait, 100)\n}(window));" + }, + "metadata": {}, + "output_type": "display_data" + }, + { + "data": { + "application/javascript": [ + "\n", + "if ((window.PyViz === undefined) || (window.PyViz instanceof HTMLElement)) {\n", + " window.PyViz = {comms: {}, comm_status:{}, kernels:{}, receivers: {}, plot_index: []}\n", + "}\n", + "\n", + "\n", + " function JupyterCommManager() {\n", + " }\n", + "\n", + " JupyterCommManager.prototype.register_target = function(plot_id, comm_id, msg_handler) {\n", + " if (window.comm_manager || ((window.Jupyter !== undefined) && (Jupyter.notebook.kernel != null))) {\n", + " var comm_manager = window.comm_manager || Jupyter.notebook.kernel.comm_manager;\n", + " comm_manager.register_target(comm_id, function(comm) {\n", + " comm.on_msg(msg_handler);\n", + " });\n", + " } else if ((plot_id in window.PyViz.kernels) && (window.PyViz.kernels[plot_id])) {\n", + " window.PyViz.kernels[plot_id].registerCommTarget(comm_id, function(comm) {\n", + " comm.onMsg = msg_handler;\n", + " });\n", + " } else if (typeof google != 'undefined' && google.colab.kernel != null) {\n", + " google.colab.kernel.comms.registerTarget(comm_id, (comm) => {\n", + " var messages = comm.messages[Symbol.asyncIterator]();\n", + " function processIteratorResult(result) {\n", + " var message = result.value;\n", + " console.log(message)\n", + " var content = {data: message.data, comm_id};\n", + " var buffers = []\n", + " for (var buffer of message.buffers || []) {\n", + " buffers.push(new DataView(buffer))\n", + " }\n", + " var metadata = message.metadata || {};\n", + " var msg = {content, buffers, metadata}\n", + " msg_handler(msg);\n", + " return messages.next().then(processIteratorResult);\n", + " }\n", + " return messages.next().then(processIteratorResult);\n", + " })\n", + " }\n", + " }\n", + "\n", + " JupyterCommManager.prototype.get_client_comm = function(plot_id, comm_id, msg_handler) {\n", + " if (comm_id in window.PyViz.comms) {\n", + " return window.PyViz.comms[comm_id];\n", + " } else if (window.comm_manager || ((window.Jupyter !== undefined) && (Jupyter.notebook.kernel != null))) {\n", + " var comm_manager = window.comm_manager || Jupyter.notebook.kernel.comm_manager;\n", + " var comm = comm_manager.new_comm(comm_id, {}, {}, {}, comm_id);\n", + " if (msg_handler) {\n", + " comm.on_msg(msg_handler);\n", + " }\n", + " } else if ((plot_id in window.PyViz.kernels) && (window.PyViz.kernels[plot_id])) {\n", + " var comm = window.PyViz.kernels[plot_id].connectToComm(comm_id);\n", + " comm.open();\n", + " if (msg_handler) {\n", + " comm.onMsg = msg_handler;\n", + " }\n", + " } else if (typeof google != 'undefined' && google.colab.kernel != null) {\n", + " var comm_promise = google.colab.kernel.comms.open(comm_id)\n", + " comm_promise.then((comm) => {\n", + " window.PyViz.comms[comm_id] = comm;\n", + " if (msg_handler) {\n", + " var messages = comm.messages[Symbol.asyncIterator]();\n", + " function processIteratorResult(result) {\n", + " var message = result.value;\n", + " var content = {data: message.data};\n", + " var metadata = message.metadata || {comm_id};\n", + " var msg = {content, metadata}\n", + " msg_handler(msg);\n", + " return messages.next().then(processIteratorResult);\n", + " }\n", + " return messages.next().then(processIteratorResult);\n", + " }\n", + " })\n", + " var sendClosure = (data, metadata, buffers, disposeOnDone) => {\n", + " return comm_promise.then((comm) => {\n", + " comm.send(data, metadata, buffers, disposeOnDone);\n", + " });\n", + " };\n", + " var comm = {\n", + " send: sendClosure\n", + " };\n", + " }\n", + " window.PyViz.comms[comm_id] = comm;\n", + " return comm;\n", + " }\n", + " window.PyViz.comm_manager = new JupyterCommManager();\n", + " \n", + "\n", + "\n", + "var JS_MIME_TYPE = 'application/javascript';\n", + "var HTML_MIME_TYPE = 'text/html';\n", + "var EXEC_MIME_TYPE = 'application/vnd.holoviews_exec.v0+json';\n", + "var CLASS_NAME = 'output';\n", + "\n", + "/**\n", + " * Render data to the DOM node\n", + " */\n", + "function render(props, node) {\n", + " var div = document.createElement(\"div\");\n", + " var script = document.createElement(\"script\");\n", + " node.appendChild(div);\n", + " node.appendChild(script);\n", + "}\n", + "\n", + "/**\n", + " * Handle when a new output is added\n", + " */\n", + "function handle_add_output(event, handle) {\n", + " var output_area = handle.output_area;\n", + " var output = handle.output;\n", + " if ((output.data == undefined) || (!output.data.hasOwnProperty(EXEC_MIME_TYPE))) {\n", + " return\n", + " }\n", + " var id = output.metadata[EXEC_MIME_TYPE][\"id\"];\n", + " var toinsert = output_area.element.find(\".\" + CLASS_NAME.split(' ')[0]);\n", + " if (id !== undefined) {\n", + " var nchildren = toinsert.length;\n", + " var html_node = toinsert[nchildren-1].children[0];\n", + " html_node.innerHTML = output.data[HTML_MIME_TYPE];\n", + " var scripts = [];\n", + " var nodelist = html_node.querySelectorAll(\"script\");\n", + " for (var i in nodelist) {\n", + " if (nodelist.hasOwnProperty(i)) {\n", + " scripts.push(nodelist[i])\n", + " }\n", + " }\n", + "\n", + " scripts.forEach( function (oldScript) {\n", + " var newScript = document.createElement(\"script\");\n", + " var attrs = [];\n", + " var nodemap = oldScript.attributes;\n", + " for (var j in nodemap) {\n", + " if (nodemap.hasOwnProperty(j)) {\n", + " attrs.push(nodemap[j])\n", + " }\n", + " }\n", + " attrs.forEach(function(attr) { newScript.setAttribute(attr.name, attr.value) });\n", + " newScript.appendChild(document.createTextNode(oldScript.innerHTML));\n", + " oldScript.parentNode.replaceChild(newScript, oldScript);\n", + " });\n", + " if (JS_MIME_TYPE in output.data) {\n", + " toinsert[nchildren-1].children[1].textContent = output.data[JS_MIME_TYPE];\n", + " }\n", + " output_area._hv_plot_id = id;\n", + " if ((window.Bokeh !== undefined) && (id in Bokeh.index)) {\n", + " window.PyViz.plot_index[id] = Bokeh.index[id];\n", + " } else {\n", + " window.PyViz.plot_index[id] = null;\n", + " }\n", + " } else if (output.metadata[EXEC_MIME_TYPE][\"server_id\"] !== undefined) {\n", + " var bk_div = document.createElement(\"div\");\n", + " bk_div.innerHTML = output.data[HTML_MIME_TYPE];\n", + " var script_attrs = bk_div.children[0].attributes;\n", + " for (var i = 0; i < script_attrs.length; i++) {\n", + " toinsert[toinsert.length - 1].childNodes[1].setAttribute(script_attrs[i].name, script_attrs[i].value);\n", + " }\n", + " // store reference to server id on output_area\n", + " output_area._bokeh_server_id = output.metadata[EXEC_MIME_TYPE][\"server_id\"];\n", + " }\n", + "}\n", + "\n", + "/**\n", + " * Handle when an output is cleared or removed\n", + " */\n", + "function handle_clear_output(event, handle) {\n", + " var id = handle.cell.output_area._hv_plot_id;\n", + " var server_id = handle.cell.output_area._bokeh_server_id;\n", + " if (((id === undefined) || !(id in PyViz.plot_index)) && (server_id !== undefined)) { return; }\n", + " var comm = window.PyViz.comm_manager.get_client_comm(\"hv-extension-comm\", \"hv-extension-comm\", function () {});\n", + " if (server_id !== null) {\n", + " comm.send({event_type: 'server_delete', 'id': server_id});\n", + " return;\n", + " } else if (comm !== null) {\n", + " comm.send({event_type: 'delete', 'id': id});\n", + " }\n", + " delete PyViz.plot_index[id];\n", + " if ((window.Bokeh !== undefined) & (id in window.Bokeh.index)) {\n", + " var doc = window.Bokeh.index[id].model.document\n", + " doc.clear();\n", + " const i = window.Bokeh.documents.indexOf(doc);\n", + " if (i > -1) {\n", + " window.Bokeh.documents.splice(i, 1);\n", + " }\n", + " }\n", + "}\n", + "\n", + "/**\n", + " * Handle kernel restart event\n", + " */\n", + "function handle_kernel_cleanup(event, handle) {\n", + " delete PyViz.comms[\"hv-extension-comm\"];\n", + " window.PyViz.plot_index = {}\n", + "}\n", + "\n", + "/**\n", + " * Handle update_display_data messages\n", + " */\n", + "function handle_update_output(event, handle) {\n", + " handle_clear_output(event, {cell: {output_area: handle.output_area}})\n", + " handle_add_output(event, handle)\n", + "}\n", + "\n", + "function register_renderer(events, OutputArea) {\n", + " function append_mime(data, metadata, element) {\n", + " // create a DOM node to render to\n", + " var toinsert = this.create_output_subarea(\n", + " metadata,\n", + " CLASS_NAME,\n", + " EXEC_MIME_TYPE\n", + " );\n", + " this.keyboard_manager.register_events(toinsert);\n", + " // Render to node\n", + " var props = {data: data, metadata: metadata[EXEC_MIME_TYPE]};\n", + " render(props, toinsert[0]);\n", + " element.append(toinsert);\n", + " return toinsert\n", + " }\n", + "\n", + " events.on('output_added.OutputArea', handle_add_output);\n", + " events.on('output_updated.OutputArea', handle_update_output);\n", + " events.on('clear_output.CodeCell', handle_clear_output);\n", + " events.on('delete.Cell', handle_clear_output);\n", + " events.on('kernel_ready.Kernel', handle_kernel_cleanup);\n", + "\n", + " OutputArea.prototype.register_mime_type(EXEC_MIME_TYPE, append_mime, {\n", + " safe: true,\n", + " index: 0\n", + " });\n", + "}\n", + "\n", + "if (window.Jupyter !== undefined) {\n", + " try {\n", + " var events = require('base/js/events');\n", + " var OutputArea = require('notebook/js/outputarea').OutputArea;\n", + " if (OutputArea.prototype.mime_types().indexOf(EXEC_MIME_TYPE) == -1) {\n", + " register_renderer(events, OutputArea);\n", + " }\n", + " } catch(err) {\n", + " }\n", + "}\n" + ], + "application/vnd.holoviews_load.v0+json": "\nif ((window.PyViz === undefined) || (window.PyViz instanceof HTMLElement)) {\n window.PyViz = {comms: {}, comm_status:{}, kernels:{}, receivers: {}, plot_index: []}\n}\n\n\n function JupyterCommManager() {\n }\n\n JupyterCommManager.prototype.register_target = function(plot_id, comm_id, msg_handler) {\n if (window.comm_manager || ((window.Jupyter !== undefined) && (Jupyter.notebook.kernel != null))) {\n var comm_manager = window.comm_manager || Jupyter.notebook.kernel.comm_manager;\n comm_manager.register_target(comm_id, function(comm) {\n comm.on_msg(msg_handler);\n });\n } else if ((plot_id in window.PyViz.kernels) && (window.PyViz.kernels[plot_id])) {\n window.PyViz.kernels[plot_id].registerCommTarget(comm_id, function(comm) {\n comm.onMsg = msg_handler;\n });\n } else if (typeof google != 'undefined' && google.colab.kernel != null) {\n google.colab.kernel.comms.registerTarget(comm_id, (comm) => {\n var messages = comm.messages[Symbol.asyncIterator]();\n function processIteratorResult(result) {\n var message = result.value;\n console.log(message)\n var content = {data: message.data, comm_id};\n var buffers = []\n for (var buffer of message.buffers || []) {\n buffers.push(new DataView(buffer))\n }\n var metadata = message.metadata || {};\n var msg = {content, buffers, metadata}\n msg_handler(msg);\n return messages.next().then(processIteratorResult);\n }\n return messages.next().then(processIteratorResult);\n })\n }\n }\n\n JupyterCommManager.prototype.get_client_comm = function(plot_id, comm_id, msg_handler) {\n if (comm_id in window.PyViz.comms) {\n return window.PyViz.comms[comm_id];\n } else if (window.comm_manager || ((window.Jupyter !== undefined) && (Jupyter.notebook.kernel != null))) {\n var comm_manager = window.comm_manager || Jupyter.notebook.kernel.comm_manager;\n var comm = comm_manager.new_comm(comm_id, {}, {}, {}, comm_id);\n if (msg_handler) {\n comm.on_msg(msg_handler);\n }\n } else if ((plot_id in window.PyViz.kernels) && (window.PyViz.kernels[plot_id])) {\n var comm = window.PyViz.kernels[plot_id].connectToComm(comm_id);\n comm.open();\n if (msg_handler) {\n comm.onMsg = msg_handler;\n }\n } else if (typeof google != 'undefined' && google.colab.kernel != null) {\n var comm_promise = google.colab.kernel.comms.open(comm_id)\n comm_promise.then((comm) => {\n window.PyViz.comms[comm_id] = comm;\n if (msg_handler) {\n var messages = comm.messages[Symbol.asyncIterator]();\n function processIteratorResult(result) {\n var message = result.value;\n var content = {data: message.data};\n var metadata = message.metadata || {comm_id};\n var msg = {content, metadata}\n msg_handler(msg);\n return messages.next().then(processIteratorResult);\n }\n return messages.next().then(processIteratorResult);\n }\n })\n var sendClosure = (data, metadata, buffers, disposeOnDone) => {\n return comm_promise.then((comm) => {\n comm.send(data, metadata, buffers, disposeOnDone);\n });\n };\n var comm = {\n send: sendClosure\n };\n }\n window.PyViz.comms[comm_id] = comm;\n return comm;\n }\n window.PyViz.comm_manager = new JupyterCommManager();\n \n\n\nvar JS_MIME_TYPE = 'application/javascript';\nvar HTML_MIME_TYPE = 'text/html';\nvar EXEC_MIME_TYPE = 'application/vnd.holoviews_exec.v0+json';\nvar CLASS_NAME = 'output';\n\n/**\n * Render data to the DOM node\n */\nfunction render(props, node) {\n var div = document.createElement(\"div\");\n var script = document.createElement(\"script\");\n node.appendChild(div);\n node.appendChild(script);\n}\n\n/**\n * Handle when a new output is added\n */\nfunction handle_add_output(event, handle) {\n var output_area = handle.output_area;\n var output = handle.output;\n if ((output.data == undefined) || (!output.data.hasOwnProperty(EXEC_MIME_TYPE))) {\n return\n }\n var id = output.metadata[EXEC_MIME_TYPE][\"id\"];\n var toinsert = output_area.element.find(\".\" + CLASS_NAME.split(' ')[0]);\n if (id !== undefined) {\n var nchildren = toinsert.length;\n var html_node = toinsert[nchildren-1].children[0];\n html_node.innerHTML = output.data[HTML_MIME_TYPE];\n var scripts = [];\n var nodelist = html_node.querySelectorAll(\"script\");\n for (var i in nodelist) {\n if (nodelist.hasOwnProperty(i)) {\n scripts.push(nodelist[i])\n }\n }\n\n scripts.forEach( function (oldScript) {\n var newScript = document.createElement(\"script\");\n var attrs = [];\n var nodemap = oldScript.attributes;\n for (var j in nodemap) {\n if (nodemap.hasOwnProperty(j)) {\n attrs.push(nodemap[j])\n }\n }\n attrs.forEach(function(attr) { newScript.setAttribute(attr.name, attr.value) });\n newScript.appendChild(document.createTextNode(oldScript.innerHTML));\n oldScript.parentNode.replaceChild(newScript, oldScript);\n });\n if (JS_MIME_TYPE in output.data) {\n toinsert[nchildren-1].children[1].textContent = output.data[JS_MIME_TYPE];\n }\n output_area._hv_plot_id = id;\n if ((window.Bokeh !== undefined) && (id in Bokeh.index)) {\n window.PyViz.plot_index[id] = Bokeh.index[id];\n } else {\n window.PyViz.plot_index[id] = null;\n }\n } else if (output.metadata[EXEC_MIME_TYPE][\"server_id\"] !== undefined) {\n var bk_div = document.createElement(\"div\");\n bk_div.innerHTML = output.data[HTML_MIME_TYPE];\n var script_attrs = bk_div.children[0].attributes;\n for (var i = 0; i < script_attrs.length; i++) {\n toinsert[toinsert.length - 1].childNodes[1].setAttribute(script_attrs[i].name, script_attrs[i].value);\n }\n // store reference to server id on output_area\n output_area._bokeh_server_id = output.metadata[EXEC_MIME_TYPE][\"server_id\"];\n }\n}\n\n/**\n * Handle when an output is cleared or removed\n */\nfunction handle_clear_output(event, handle) {\n var id = handle.cell.output_area._hv_plot_id;\n var server_id = handle.cell.output_area._bokeh_server_id;\n if (((id === undefined) || !(id in PyViz.plot_index)) && (server_id !== undefined)) { return; }\n var comm = window.PyViz.comm_manager.get_client_comm(\"hv-extension-comm\", \"hv-extension-comm\", function () {});\n if (server_id !== null) {\n comm.send({event_type: 'server_delete', 'id': server_id});\n return;\n } else if (comm !== null) {\n comm.send({event_type: 'delete', 'id': id});\n }\n delete PyViz.plot_index[id];\n if ((window.Bokeh !== undefined) & (id in window.Bokeh.index)) {\n var doc = window.Bokeh.index[id].model.document\n doc.clear();\n const i = window.Bokeh.documents.indexOf(doc);\n if (i > -1) {\n window.Bokeh.documents.splice(i, 1);\n }\n }\n}\n\n/**\n * Handle kernel restart event\n */\nfunction handle_kernel_cleanup(event, handle) {\n delete PyViz.comms[\"hv-extension-comm\"];\n window.PyViz.plot_index = {}\n}\n\n/**\n * Handle update_display_data messages\n */\nfunction handle_update_output(event, handle) {\n handle_clear_output(event, {cell: {output_area: handle.output_area}})\n handle_add_output(event, handle)\n}\n\nfunction register_renderer(events, OutputArea) {\n function append_mime(data, metadata, element) {\n // create a DOM node to render to\n var toinsert = this.create_output_subarea(\n metadata,\n CLASS_NAME,\n EXEC_MIME_TYPE\n );\n this.keyboard_manager.register_events(toinsert);\n // Render to node\n var props = {data: data, metadata: metadata[EXEC_MIME_TYPE]};\n render(props, toinsert[0]);\n element.append(toinsert);\n return toinsert\n }\n\n events.on('output_added.OutputArea', handle_add_output);\n events.on('output_updated.OutputArea', handle_update_output);\n events.on('clear_output.CodeCell', handle_clear_output);\n events.on('delete.Cell', handle_clear_output);\n events.on('kernel_ready.Kernel', handle_kernel_cleanup);\n\n OutputArea.prototype.register_mime_type(EXEC_MIME_TYPE, append_mime, {\n safe: true,\n index: 0\n });\n}\n\nif (window.Jupyter !== undefined) {\n try {\n var events = require('base/js/events');\n var OutputArea = require('notebook/js/outputarea').OutputArea;\n if (OutputArea.prototype.mime_types().indexOf(EXEC_MIME_TYPE) == -1) {\n register_renderer(events, OutputArea);\n }\n } catch(err) {\n }\n}\n" + }, + "metadata": {}, + "output_type": "display_data" + }, + { + "data": { + "application/vnd.holoviews_exec.v0+json": "", + "text/html": [ + "
\n", + "
\n", + "
\n", + "" + ] + }, + "metadata": { + "application/vnd.holoviews_exec.v0+json": { + "id": "44faedca-a5ab-4a8e-a0e4-b175acdffb67" + } + }, + "output_type": "display_data" + }, + { + "name": "stderr", + "output_type": "stream", + "text": [ + "WARNING:param.main: geo option cannot be used with kind='scatter' plot type. Geographic plots are only supported for following plot types: ['bivariate', 'contour', 'contourf', 'dataset', 'hexbin', 'image', 'labels', 'paths', 'points', 'points', 'polygons', 'quadmesh', 'rgb', 'vectorfield']\n", + "WARNING:param.main: geo option cannot be used with kind='scatter' plot type. Geographic plots are only supported for following plot types: ['bivariate', 'contour', 'contourf', 'dataset', 'hexbin', 'image', 'labels', 'paths', 'points', 'points', 'polygons', 'quadmesh', 'rgb', 'vectorfield']\n" + ] + }, + { + "data": {}, + "metadata": {}, + "output_type": "display_data" + }, + { + "data": { + "application/vnd.holoviews_exec.v0+json": "", + "text/html": [ + "
\n", + "
\n", + "
\n", + "" + ], + "text/plain": [ + ":Layout\n", + " .Overlay.I :Overlay\n", + " .Scatter.I :Scatter [lon] (lat,obs)\n", + " .Coastline.I :Feature [Longitude,Latitude]\n", + " .Overlay.II :Overlay\n", + " .Scatter.I :Scatter [lon] (lat,nobs)\n", + " .Coastline.I :Feature [Longitude,Latitude]" + ] + }, + "execution_count": 16, + "metadata": { + "application/vnd.holoviews_exec.v0+json": { + "id": "6e745ea5-6980-4d0c-837d-868b40470ebe" + } + }, + "output_type": "execute_result" + } + ], + "source": [ + "import hvplot.pandas\n", + "import geoviews.feature as gf\n", + "import geoviews as gv\n", + "(\n", + " data.hvplot.scatter(x='lon', y='lat', c='obs', s=1, global_extent=True, coastline=True, title='aggregated obs') + \n", + " data.hvplot.scatter(x='lon', y='lat', c='nobs', s=1, global_extent=True, coastline=True, title='nobs')\n", + ").cols(1)" + ] + }, + { + "cell_type": "code", + "execution_count": 18, + "id": "15159bb2-a2db-4f17-819c-21443031ff06", + "metadata": {}, + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "WARNING:param.main: geo option cannot be used with kind='scatter' plot type. Geographic plots are only supported for following plot types: ['bivariate', 'contour', 'contourf', 'dataset', 'hexbin', 'image', 'labels', 'paths', 'points', 'points', 'polygons', 'quadmesh', 'rgb', 'vectorfield']\n" + ] + }, + { + "data": {}, + "metadata": {}, + "output_type": "display_data" + }, + { + "data": { + "application/vnd.holoviews_exec.v0+json": "", + "text/html": [ + "
\n", + "
\n", + "
\n", + "" + ], + "text/plain": [ + ":Overlay\n", + " .Scatter.I :Scatter [rotated_lon] (rotated_lat,obs)\n", + " .Coastline.I :Feature [Longitude,Latitude]" + ] + }, + "execution_count": 18, + "metadata": { + "application/vnd.holoviews_exec.v0+json": { + "id": "d291f6de-424b-4be4-a8b0-b140088953a5" + } + }, + "output_type": "execute_result" + } + ], + "source": [ + "# Check that the rotation worked\n", + "\n", + "data.hvplot.scatter(x='rotated_lon', y='rotated_lat', c='obs', s=1, global_extent=True, coastline=True, title='aggregated obs, rotated (coastline not rotated ...)')" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "d356b9e7-50c1-4b67-908b-0349c880eaa2", + "metadata": {}, + "outputs": [], + "source": [ + "from dataclasses import dataclass\n", + "\n", + "\n", + "def my_actually_costly_but_externalized_operation(arg):\n", + " ...\n", + "\n", + "\n", + "@dataclass(kw_only=True\n", + "class MyClass:\n", + " parallelization : str | None = None\n", + "\n", + " \n", + " def my_costly_operation(self, *args) -> Any:\n", + " if self.parallelization == 'mp':\n", + " return self._my_costly_but_parallelized_operation(*args)\n", + " if self.parallelization == 'mpi':\n", + " return self._my_costly_but_massively_parallelized_operation(*args)\n", + " return self._my_costly_operation()\n", + "\n", + " def _my_costly_operation(self, *args):\n", + " \"\"\"\n", + " Serial, good for developing and when multiprocessing isn't needed\n", + " \"\"\"\n", + " out = []\n", + " for arg in some_iterator:\n", + " out.append(my_actually_costly_but_externalized_operation(arg))\n", + " return out\n", + "\n", + " def _my_costly_but_parallelized_operation(self, *args):\n", + " \"\"\"\n", + " multiprocessing. Good for CPU intensive tasks\n", + " \"\"\"\n", + " import multiprocessing as mp\n", + " with mp.Pool() as pool:\n", + " return pool.map(my_actually_costly_but_externalized_operation, some_iterator)\n", + "\n", + " def _my_costly_but_massively_parallelized_operation(self, *args):\n", + " \"\"\"\n", + " MPI. Good when it really doesn't fit on one node\n", + " \"\"\"\n", + " raise NotImplementedError('placeholder for MPI parallelization') # I do that trick a lot for things I might implement later if needed" + ] + } + ], + "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 +} diff --git a/mapies/viirs.py b/mapies/viirs.py index 6feb60e2..cf1f54c5 100644 --- a/mapies/viirs.py +++ b/mapies/viirs.py @@ -11,6 +11,7 @@ import multiprocessing import gc import dask.array as da from glob import glob as glob_module +import sys import numpy as np import pandas as pd @@ -43,7 +44,7 @@ def process_single_file(r, file, obs_var, lat_var, lon_var, time_var, start_date """ try: ds = xr.open_dataset(file, engine="h5netcdf") - + obs, lat, lon, time_values, time_values_index = new_preprocess_vars(ds, obs_var, lat_var, lon_var, time_var, start_date, end_date) # obs = ds[obs_var].values.flatten() # lat = ds[lat_var].values.flatten() @@ -62,6 +63,9 @@ def process_single_file(r, file, obs_var, lat_var, lon_var, time_var, start_date else: ds.close() return obs, lat, lon, time_values + + + except Exception as e: print(f"Error processing file {file}: {e}") @@ -340,10 +344,12 @@ class VIIRS(MAPIES): self.lat_values = final_lat self.obs = final_obs self.count_obs = cumulative_count + except Exception as e: print(f"Error processing data for the time period from {self.start_date} to {self.end_date}: {e}") return + @timeit def process_monthly_data(self, batch_size, month=None): """ @@ -436,6 +442,7 @@ class VIIRS(MAPIES): return final_lon, final_lat, final_obs, cumulative_count + @timeit def yearly_average(self): """ @@ -554,25 +561,50 @@ class VIIRS(MAPIES): continue try: - # Save directly to NetCDF + # Convert time to seconds since TIME_ORIGIN + final_time_seconds = (final_time - TIME_ORIGIN) / np.timedelta64(1, "s") + final_time_seconds = np.round(final_time_seconds, 0).astype(int) + + print(f'Final time in seconds unique values: {np.unique(final_time_seconds)}') + + # Create xarray Dataset ds = xr.Dataset( - coords={"time": final_time}, + coords={"time": ("time", final_time_seconds, { + "units": "seconds since 1900-01-01 00:00:00", + "calendar": "gregorian" + })}, data_vars={ - "lon": (["time"], final_lon), - "lat": (["time"], final_lat), - "obs": (["time"], final_obs), + "lon": (["time"], final_lon, {"description": "Longitude of observations", "units": "degrees_east"}), + "lat": (["time"], final_lat, {"description": "Latitude of observations", "units": "degrees_north"}), + "obs": (["time"], final_obs, {"description": "Observed values", "units": "unitless"}), + }, + attrs={ + "title": f"VIIRS daily data for day {day}", + "institution": "Barcelona Supercomputing Center", + "grid": f"{self.grid_repr}: {self.grid_config}", + "Developers": "MAPIES team", + "QA": f"Quality assurance applied: {self.qualityFlags}", + "history": f"Created {datetime.now()}", }, ) + + # Define encoding for compact storage + encoding = { + "time": {"dtype": "float32"}, + "lon": {"dtype": "float32"}, + "lat": {"dtype": "float32"}, + "obs": {"dtype": "float32"}, + } + + # Save to NetCDF file filename = f"{self.dest}/daily_{day}.nc" - ds.to_netcdf(filename, encoding={}) + ds.to_netcdf(filename, mode='w', encoding=encoding) print(f"Saved daily data for day {day} to {filename}") - + except Exception as e: print(f"Error saving daily data for day {day}: {e}") continue - - # ============================================================================= # Supporting functions diff --git a/run/run_viirs.py b/run/run_viirs.py index bd8b1504..3f4e634e 100644 --- a/run/run_viirs.py +++ b/run/run_viirs.py @@ -28,8 +28,8 @@ if __name__ == "__main__": # indir = "/home/cgile/bscearth000/esarchive/obs/nasa/viirs_noaa20_aerdb_l2/original_files/VIIRS" indir = "/home/cgile/Documents/mapies/VIIRS" - filename_obs = 'carlotta' - filename_num_obs = 'carlotta' + filename_obs = 'January_1day_obs_centre_lat:35centre_lon:20west:-51south:-35dlon:0.2dlat:0.1_MP' + filename_num_obs = 'January_1day_num_obs_obs_centre_lat:35centre_lon:20west:-51south:-35dlon:0.2dlat:0.1_MP' start_time = time.time() @@ -37,9 +37,9 @@ if __name__ == "__main__": c.read_nc() c.process_data(monthly_avg=False, batch_size = 100) # c.yearly_average() - c.plot_2D_observations(months=[0], filename=filename_obs, outdir=outdir) - c.plot_2D_num_obs(months=[0], filename=filename_num_obs, outdir=outdir) - # c.process_lazy_data() + # c.plot_2D_observations(months=[0], filename=filename_obs, outdir=outdir) + # c.plot_2D_num_obs(months=[0], filename=filename_num_obs, outdir=outdir) + c.process_lazy_data() end_time = time.time() elapsed_time = end_time - start_time print(f"Script executed in {elapsed_time:.2f} seconds.") -- GitLab From ea8fddc7e8354ecb3efc03c08a8bc19b2953efd3 Mon Sep 17 00:00:00 2001 From: Carlotta Date: Wed, 5 Feb 2025 16:29:29 +0100 Subject: [PATCH 2/5] changes --- mapies/viirs.py | 10 ++++++---- 1 file changed, 6 insertions(+), 4 deletions(-) diff --git a/mapies/viirs.py b/mapies/viirs.py index cf1f54c5..5a2fa947 100644 --- a/mapies/viirs.py +++ b/mapies/viirs.py @@ -27,6 +27,7 @@ from mapies.util.func_tools import * from mapies.grids.monarch import RotatedGrid, IrregularRotatedGrid, RegularGrid, regridding_function from pathlib import Path from glob import glob as glob_module +import json os.environ["OMP_NUM_THREADS"] = "8" @@ -565,7 +566,8 @@ class VIIRS(MAPIES): final_time_seconds = (final_time - TIME_ORIGIN) / np.timedelta64(1, "s") final_time_seconds = np.round(final_time_seconds, 0).astype(int) - print(f'Final time in seconds unique values: {np.unique(final_time_seconds)}') + grid_representation_str = json.dumps(self.grid_config, ensure_ascii=False) + qa_flags_str = json.dumps(self.qualityFlags, ensure_ascii=False) # Create xarray Dataset ds = xr.Dataset( @@ -581,16 +583,16 @@ class VIIRS(MAPIES): attrs={ "title": f"VIIRS daily data for day {day}", "institution": "Barcelona Supercomputing Center", - "grid": f"{self.grid_repr}: {self.grid_config}", + "grid": f"{self.grid_repr}: {grid_representation_str}", "Developers": "MAPIES team", - "QA": f"Quality assurance applied: {self.qualityFlags}", + "QA": f"Quality assurance applied: {qa_flags_str}", "history": f"Created {datetime.now()}", }, ) # Define encoding for compact storage encoding = { - "time": {"dtype": "float32"}, + "time": {"dtype": "int32"}, "lon": {"dtype": "float32"}, "lat": {"dtype": "float32"}, "obs": {"dtype": "float32"}, -- GitLab From 5a30d46421de069c456b287881ff64e758c1e26d Mon Sep 17 00:00:00 2001 From: Carlotta Date: Wed, 5 Feb 2025 17:07:49 +0100 Subject: [PATCH 3/5] updates --- mapies/viirs.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/mapies/viirs.py b/mapies/viirs.py index 5a2fa947..f12d2742 100644 --- a/mapies/viirs.py +++ b/mapies/viirs.py @@ -592,7 +592,7 @@ class VIIRS(MAPIES): # Define encoding for compact storage encoding = { - "time": {"dtype": "int32"}, + "time": {"dtype": "float32"}, "lon": {"dtype": "float32"}, "lat": {"dtype": "float32"}, "obs": {"dtype": "float32"}, -- GitLab From 66cfe469908f31290e8d31b089632d160dc63e66 Mon Sep 17 00:00:00 2001 From: Carlotta Date: Thu, 6 Feb 2025 17:47:19 +0100 Subject: [PATCH 4/5] new to_netCDF method --- mapies/util/variables_names.py | 6 +++ mapies/viirs.py | 95 +++++++++++++++++++--------------- 2 files changed, 58 insertions(+), 43 deletions(-) create mode 100644 mapies/util/variables_names.py diff --git a/mapies/util/variables_names.py b/mapies/util/variables_names.py new file mode 100644 index 00000000..12be66d7 --- /dev/null +++ b/mapies/util/variables_names.py @@ -0,0 +1,6 @@ +VARS_ATTRIBUTES={ + "lon": {"description": "Longitude of observations", "units": "degrees_east"}, + "lat": {"description": "Latitude of observations", "units": "degrees_north"}, + "Aerosol_Optical_Thickness_550_Land_Ocean_Best_Estimate": {"description": "Aerosol_Optical_Thickness_550_Land_Ocean_Best_Estimate values", "units": "unitless"}, + "title": {"description": "VIIRS data"}, + } \ No newline at end of file diff --git a/mapies/viirs.py b/mapies/viirs.py index f12d2742..26e7c928 100644 --- a/mapies/viirs.py +++ b/mapies/viirs.py @@ -25,6 +25,7 @@ from .mapies import MAPIES from mapies.util.func_tools import * #timeit, get_file_list, time_converter, frequency_int, error_estimation, time_domain_selection, geo_to_rot from mapies.grids.monarch import RotatedGrid, IrregularRotatedGrid, RegularGrid, regridding_function +from mapies.util.variables_names import VARS_ATTRIBUTES from pathlib import Path from glob import glob as glob_module import json @@ -290,7 +291,7 @@ class VIIRS(MAPIES): @timeit - def process_data(self, monthly_avg=False, batch_size=100): + def process_data(self, monthly_avg=False, batch_size=100, save=False): """ Process the data for the specified year and months. @@ -346,6 +347,9 @@ class VIIRS(MAPIES): self.obs = final_obs self.count_obs = cumulative_count + if save: + self.to_netCDF(self.obs, self.lat_values, self.lon_values, self.start_date, self.end_date) + except Exception as e: print(f"Error processing data for the time period from {self.start_date} to {self.end_date}: {e}") return @@ -561,52 +565,57 @@ class VIIRS(MAPIES): print(f"Error processing data for day {day}: {e}") continue - try: - # Convert time to seconds since TIME_ORIGIN - final_time_seconds = (final_time - TIME_ORIGIN) / np.timedelta64(1, "s") - final_time_seconds = np.round(final_time_seconds, 0).astype(int) - - grid_representation_str = json.dumps(self.grid_config, ensure_ascii=False) - qa_flags_str = json.dumps(self.qualityFlags, ensure_ascii=False) - - # Create xarray Dataset - ds = xr.Dataset( - coords={"time": ("time", final_time_seconds, { - "units": "seconds since 1900-01-01 00:00:00", - "calendar": "gregorian" - })}, - data_vars={ - "lon": (["time"], final_lon, {"description": "Longitude of observations", "units": "degrees_east"}), - "lat": (["time"], final_lat, {"description": "Latitude of observations", "units": "degrees_north"}), - "obs": (["time"], final_obs, {"description": "Observed values", "units": "unitless"}), - }, - attrs={ - "title": f"VIIRS daily data for day {day}", - "institution": "Barcelona Supercomputing Center", - "grid": f"{self.grid_repr}: {grid_representation_str}", - "Developers": "MAPIES team", - "QA": f"Quality assurance applied: {qa_flags_str}", - "history": f"Created {datetime.now()}", - }, - ) + self.to_netCDF(final_obs, final_lat, final_lon, final_time, self.start_date, self.end_date) - # Define encoding for compact storage - encoding = { - "time": {"dtype": "float32"}, - "lon": {"dtype": "float32"}, - "lat": {"dtype": "float32"}, - "obs": {"dtype": "float32"}, - } + def to_netCDF(self, obs, lat, lon, time, start_time, end_time): + try: + # Convert time to seconds since TIME_ORIGIN + final_time_seconds = (time - TIME_ORIGIN) / np.timedelta64(1, "s") + final_time_seconds = np.round(final_time_seconds, 0) - # Save to NetCDF file - filename = f"{self.dest}/daily_{day}.nc" - ds.to_netcdf(filename, mode='w', encoding=encoding) + grid_representation_str = json.dumps(self.grid_config, ensure_ascii=False) + qa_flags_str = json.dumps(self.qualityFlags, ensure_ascii=False) - print(f"Saved daily data for day {day} to {filename}") + # Create xarray Dataset + ds = xr.Dataset( + coords={"time": ("time", final_time_seconds, { + "units": "seconds since 1900-01-01 00:00:00", + "calendar": "gregorian" + })}, + + data_vars={ + "lon": (["time"], lon, VARS_ATTRIBUTES['lon'].copy()), + "lat": (["time"], lat, VARS_ATTRIBUTES['lat'].copy()), + self.obs_var: (["time"], obs, VARS_ATTRIBUTES[self.obs_var].copy()), + }, + attrs={ + "title": f"{VARS_ATTRIBUTES['title']['description']} from {start_time} to {end_time}", + "institution": "Barcelona Supercomputing Center", + "grid": f"{self.grid_repr}: {grid_representation_str}", + "Developers": "MAPIES team", + "QA": f"Quality assurance applied: {qa_flags_str}", + "history": f"Created {datetime.now()}", + }, + ) + + # Define encoding for compact storage + encoding = { + "time": {"dtype": "float32"}, + "lon": {"dtype": "float32"}, + "lat": {"dtype": "float32"}, + self.obs_var: {"dtype": "float32"}, + } + + # Save to NetCDF file + filename = f"{self.dest}/Data_{start_time}_{end_time}.nc" + ds.to_netcdf(filename, mode='w', encoding=encoding) + + print(f"Saved data to {filename}") + + except Exception as e: + print(f"Error saving data into netCDF: {e}") + return - except Exception as e: - print(f"Error saving daily data for day {day}: {e}") - continue # ============================================================================= # Supporting functions -- GitLab From 4770746d1aba85d1e645215daea34cd3d62dce0f Mon Sep 17 00:00:00 2001 From: Carlotta Date: Mon, 10 Feb 2025 12:09:43 +0100 Subject: [PATCH 5/5] updated lon, lat names --- mapies/viirs.py | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/mapies/viirs.py b/mapies/viirs.py index 26e7c928..c6aef82e 100644 --- a/mapies/viirs.py +++ b/mapies/viirs.py @@ -342,13 +342,13 @@ class VIIRS(MAPIES): final_lon, final_lat, final_obs, cumulative_count = self.process_monthly_data( batch_size=batch_size, month=None ) - self.lon_values = final_lon - self.lat_values = final_lat + self.lon = final_lon + self.lat = final_lat self.obs = final_obs self.count_obs = cumulative_count if save: - self.to_netCDF(self.obs, self.lat_values, self.lon_values, self.start_date, self.end_date) + self.to_netCDF(self.obs, self.lat, self.lon, self.start_date, self.end_date) except Exception as e: print(f"Error processing data for the time period from {self.start_date} to {self.end_date}: {e}") @@ -725,8 +725,8 @@ class VIIRS(MAPIES): if months == [0]: try: - lon_values = self.lon_values - lat_values = self.lat_values + lon_values = self.lon + lat_values = self.lat obs_values = self.obs except Exception as e: print(f"Error plotting all data: {e}") @@ -785,8 +785,8 @@ class VIIRS(MAPIES): if months == [0]: try: - lon_values = self.lon_values - lat_values = self.lat_values + lon_values = self.lon + lat_values = self.lat cum_count = self.count_obs except Exception as e: print(f"Error plotting all data: {e}") -- GitLab