diff --git a/metric_check.ipynb b/metric_check.ipynb
new file mode 100644
index 0000000..3554ced
--- /dev/null
+++ b/metric_check.ipynb
@@ -0,0 +1,2050 @@
+{
+ "cells": [
+ {
+ "cell_type": "code",
+ "execution_count": 1,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Python 3.11.9\n"
+ ]
+ }
+ ],
+ "source": [
+ "!python --version\n"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 2,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Libraries\n",
+ "\n",
+ "import torch\n",
+ "import numpy as np\n",
+ "import matplotlib.pyplot as plt\n",
+ "import einops\n",
+ "from torch.nn.functional import interpolate\n",
+ "from glob import glob\n",
+ "import xarray as xr\n",
+ "import numpy as np\n",
+ "import netCDF4"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 3,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "def get_deptht(restart, mask):\n",
+ " \"\"\"\n",
+ " Calculate the depth of each vertical level on grid T in the 3D grid.\n",
+ "\n",
+ " Parameters:\n",
+ " restart (xarray.Dataset) : The dataset containing ocean model variables.\n",
+ " mask (xarray.Dataset) : The dataset containing mask variables.\n",
+ "\n",
+ " Returns:\n",
+ " deptht (numpy.array) : The depth of each vertical level.\n",
+ " \"\"\"\n",
+ " ssh = restart.sshn.squeeze()\n",
+ " e3w_0 = mask.e3w_0.squeeze() # initial z axis cell's thickness on grid W - (t,z,y,x)\n",
+ " e3t_0 = mask.e3t_0.squeeze() # initial z axis cell's thickness on grid T - (t,z,y,x)\n",
+ " tmask = mask.tmask.squeeze() # grid T continent mask - (t,z,y,x)\n",
+ " ssmask = tmask[:, 0] # bathymetry - (t,y,x)\n",
+ " bathy = e3t_0.sum(\n",
+ " dim=\"depth\"\n",
+ " ) # initial condition depth 0 - (t,z,y,x)\n",
+ " depth_0 = e3w_0.copy().squeeze()\n",
+ " depth_0[:, 0] = 0.5 * e3w_0[:, 0]\n",
+ " depth_0[:, 1:] = depth_0[:, 0:1].data + e3w_0[:, 1:].cumsum(dim=\"depth\")\n",
+ " deptht = depth_0 * (1 + ssh / (bathy + 1 - ssmask)) * tmask\n",
+ " return deptht\n",
+ "\n",
+ "\n",
+ "def get_density(thetao, so, depth, tmask):\n",
+ " \"\"\"\n",
+ " Compute potential density referenced at the surface and density anomaly.\n",
+ "\n",
+ " Parameters:\n",
+ " thetao (numpy.array) : Temperature array - (t,z,y,x).\n",
+ " so (numpy.array) : Salinity array - (t,z,y,x).\n",
+ " depth (numpy.array) : Depth array - (t,z,y,x).\n",
+ " tmask (numpy.array) : Mask array - (t,z,y,x).\n",
+ "\n",
+ " Returns:\n",
+ " tuple: A tuple containing:\n",
+ " array: Potential density referenced at the surface.\n",
+ " array: Density anomaly.\n",
+ " \"\"\"\n",
+ " rdeltaS = 32.0\n",
+ " r1_S0 = 0.875 / 35.16504\n",
+ " r1_T0 = 1.0 / 40.0\n",
+ " r1_Z0 = 1.0e-4\n",
+ "\n",
+ " EOS000 = 8.0189615746e02\n",
+ " EOS100 = 8.6672408165e02\n",
+ " EOS200 = -1.7864682637e03\n",
+ " EOS300 = 2.0375295546e03\n",
+ " EOS400 = -1.2849161071e03\n",
+ " EOS500 = 4.3227585684e02\n",
+ " EOS600 = -6.0579916612e01\n",
+ " EOS010 = 2.6010145068e01\n",
+ " EOS110 = -6.5281885265e01\n",
+ " EOS210 = 8.1770425108e01\n",
+ " EOS310 = -5.6888046321e01\n",
+ " EOS410 = 1.7681814114e01\n",
+ " EOS510 = -1.9193502195\n",
+ " EOS020 = -3.7074170417e01\n",
+ " EOS120 = 6.1548258127e01\n",
+ " EOS220 = -6.0362551501e01\n",
+ " EOS320 = 2.9130021253e01\n",
+ " EOS420 = -5.4723692739\n",
+ " EOS030 = 2.1661789529e01\n",
+ " EOS130 = -3.3449108469e01\n",
+ " EOS230 = 1.9717078466e01\n",
+ " EOS330 = -3.1742946532\n",
+ " EOS040 = -8.3627885467\n",
+ " EOS140 = 1.1311538584e01\n",
+ " EOS240 = -5.3563304045\n",
+ " EOS050 = 5.4048723791e-01\n",
+ " EOS150 = 4.8169980163e-01\n",
+ " EOS060 = -1.9083568888e-01\n",
+ " EOS001 = 1.9681925209e01\n",
+ " EOS101 = -4.2549998214e01\n",
+ " EOS201 = 5.0774768218e01\n",
+ " EOS301 = -3.0938076334e01\n",
+ " EOS401 = 6.6051753097\n",
+ " EOS011 = -1.3336301113e01\n",
+ " EOS111 = -4.4870114575\n",
+ " EOS211 = 5.0042598061\n",
+ " EOS311 = -6.5399043664e-01\n",
+ " EOS021 = 6.7080479603\n",
+ " EOS121 = 3.5063081279\n",
+ " EOS221 = -1.8795372996\n",
+ " EOS031 = -2.4649669534\n",
+ " EOS131 = -5.5077101279e-01\n",
+ " EOS041 = 5.5927935970e-01\n",
+ " EOS002 = 2.0660924175\n",
+ " EOS102 = -4.9527603989\n",
+ " EOS202 = 2.5019633244\n",
+ " EOS012 = 2.0564311499\n",
+ " EOS112 = -2.1311365518e-01\n",
+ " EOS022 = -1.2419983026\n",
+ " EOS003 = -2.3342758797e-02\n",
+ " EOS103 = -1.8507636718e-02\n",
+ " EOS013 = 3.7969820455e-01\n",
+ "\n",
+ " zh = depth * r1_Z0 # depth\n",
+ " zt = thetao * r1_T0 # temperature\n",
+ " zs = np.sqrt(np.abs(so + rdeltaS) * r1_S0) # square root salinity\n",
+ " ztm = tmask.squeeze()\n",
+ "\n",
+ " zn3 = EOS013 * zt + EOS103 * zs + EOS003\n",
+ " zn2 = (\n",
+ " (EOS022 * zt + EOS112 * zs + EOS012) * zt + (EOS202 * zs + EOS102) * zs + EOS002\n",
+ " )\n",
+ " zn1 = (\n",
+ " (\n",
+ " (\n",
+ " (EOS041 * zt + EOS131 * zs + EOS031) * zt\n",
+ " + (EOS221 * zs + EOS121) * zs\n",
+ " + EOS021\n",
+ " )\n",
+ " * zt\n",
+ " + ((EOS311 * zs + EOS211) * zs + EOS111) * zs\n",
+ " + EOS011\n",
+ " )\n",
+ " * zt\n",
+ " + (((EOS401 * zs + EOS301) * zs + EOS201) * zs + EOS101) * zs\n",
+ " + EOS001\n",
+ " )\n",
+ " zn0 = (\n",
+ " (\n",
+ " (\n",
+ " (\n",
+ " (\n",
+ " (EOS060 * zt + EOS150 * zs + EOS050) * zt\n",
+ " + (EOS240 * zs + EOS140) * zs\n",
+ " + EOS040\n",
+ " )\n",
+ " * zt\n",
+ " + ((EOS330 * zs + EOS230) * zs + EOS130) * zs\n",
+ " + EOS030\n",
+ " )\n",
+ " * zt\n",
+ " + (((EOS420 * zs + EOS320) * zs + EOS220) * zs + EOS120) * zs\n",
+ " + EOS020\n",
+ " )\n",
+ " * zt\n",
+ " + ((((EOS510 * zs + EOS410) * zs + EOS310) * zs + EOS210) * zs + EOS110)\n",
+ " * zs\n",
+ " + EOS010\n",
+ " )\n",
+ " * zt\n",
+ " + (\n",
+ " ((((EOS600 * zs + EOS500) * zs + EOS400) * zs + EOS300) * zs + EOS200) * zs\n",
+ " + EOS100\n",
+ " )\n",
+ " * zs\n",
+ " + EOS000\n",
+ " )\n",
+ "\n",
+ " zn = ((zn3 * zh + zn2) * zh + zn1) * zh + zn0\n",
+ "\n",
+ " rhop = zn0 * ztm # potential density referenced at the surface\n",
+ " rho_insitu = zn * ztm # density anomaly (masked)\n",
+ " return rhop, rho_insitu\n"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 4,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "def check_density(density, epsilon=1e-5) :\n",
+ " \"\"\"\n",
+ " args :\n",
+ " density (xarray) : DataArray (t, depth, lat, lon) with density value for each point of the grid.\n",
+ " return :\n",
+ " (float) proportion of points not respecting density decreasing constraint\n",
+ " \"\"\"\n",
+ " density=density.where(density!=0)\n",
+ " diff = density - density.shift(depth=-1)\n",
+ " return (diff > epsilon).mean().data # Proportion of points not respecting decreasing density\n",
+ "\n",
+ "\n",
+ "\n",
+ "\n",
+ "\n",
+ "def temperature_500m_30NS_metric(temperature, file_mask):\n",
+ " \"\"\"\n",
+ " Metric Extraction Function :\n",
+ " Average Temperature at 500m depth between 30N and 30S.\n",
+ " Unit : °C \n",
+ " \n",
+ " \n",
+ " Input : \n",
+ " - thetao : xarray.DataArray\n",
+ " - file_mask : xarray.Dataset\n",
+ " Output : \n",
+ " - np.float32 or np.float64 depending on recording precision of simulation files \n",
+ " \n",
+ " \"\"\"\n",
+ " \n",
+ " # Taking Temperature At 500m depth and between 30N and 30S.\n",
+ "\n",
+ " t500_30NS=temperature.sel(depth=500,method='nearest').where(abs(temperature.nav_lat)<30,drop=False)\n",
+ "\n",
+ " # Computing Area Weights from Mask over 30N-30S latitude zone and @500m depth\n",
+ " e1t=file_mask.e1t.squeeze()\n",
+ " e2t=file_mask.e2t.squeeze()\n",
+ " tmask=file_mask.tmask.squeeze()\n",
+ " area_500m_30NS=e1t*e2t*tmask.sel(depth=500,method='nearest').where(abs(temperature.nav_lat)<30,drop=False)\n",
+ "\n",
+ " #Returning Average Temperature at 500m depth as a numpy scalar\n",
+ " return ((t500_30NS*area_500m_30NS).sum(dim=[\"nav_lat\",\"nav_lon\"])/area_500m_30NS.sum(dim=[\"nav_lat\",\"nav_lon\"]))\n",
+ "\n",
+ "\n",
+ "\n",
+ "def temperature_BWbox_metric(thetao, file_mask):\n",
+ " \"\"\"\n",
+ " Metric Extraction Function :\n",
+ " Average Temperature in a U-shaped \"Bottom Water\" box corresponding to waters below 3000m or beyond 30 degrees of latitude North and South.\n",
+ " \n",
+ " ________________________________________________ _Surface\n",
+ " | . . . . |__________________________| . . . . |_500m\n",
+ " | . . . . | | . . . . |\n",
+ " | . . . . | Deep Water | . . . . |\n",
+ " | . . . . |__________________________| . . . . |_3000m\n",
+ " | . . . . . . . . Bottom Water . . . . . . . . |\n",
+ " |______________________________________________|_Bottom\n",
+ " S 30S Eq. 30N N\n",
+ " \n",
+ " Figure : Schematic Representation of the Bottom Water box used in this metric.\n",
+ "\n",
+ " Unit : °C \n",
+ " \n",
+ " Input : \n",
+ " - thetao : xarray.DataArray\n",
+ " - file_mask : xarray.Dataset \n",
+ " Output : \n",
+ " - np.float32 or np.float64 depending on recording precision of simulation files \n",
+ " \n",
+ " \"\"\"\n",
+ "\n",
+ " t_BW=thetao.where(1-(thetao.depth<3000)*(abs(thetao.nav_lat)<30))\n",
+ "\n",
+ " # Computing Area Weights from Mask over Box\n",
+ " e1t=file_mask.e1t.squeeze()\n",
+ " e2t=file_mask.e2t.squeeze()\n",
+ " tmask=file_mask.tmask.squeeze()\n",
+ " area_BW=e1t*e2t*tmask.where(1-(thetao.depth<3000)*(abs(thetao.nav_lat)<30))\n",
+ "\n",
+ " #Returning Average Temperature on Box\n",
+ " return ((t_BW*area_BW).sum(dim=[\"nav_lat\",\"nav_lon\",\"depth\"])/area_BW.sum(dim=[\"nav_lat\",\"nav_lon\",\"depth\"]))\n",
+ "\n",
+ "\n",
+ "\n",
+ "def temperature_DWbox_metric(thetao, file_mask):\n",
+ " \"\"\"\n",
+ " Metric Extraction Function :\n",
+ " Average Temperature in a \"Deep Water\" box corresponding to waters between 500m and 3000m depth and 30°N and 30°S.\n",
+ " \n",
+ " ________________________________________________ _Surface\n",
+ " | |__________________________| |_500m\n",
+ " | | . . . . . . . . . . . . .| |\n",
+ " | | . . . .Deep Water . . . .| |\n",
+ " | |__________________________| |_3000m\n",
+ " | Bottom Water |\n",
+ " |______________________________________________|_Bottom\n",
+ " S 30S Eq. 30N N\n",
+ " \n",
+ " Figure : Schematic Representation of the Deep Water box used in this metric.\n",
+ "\n",
+ " Unit : °C \n",
+ " \n",
+ " Input : \n",
+ " - thetao : xarray.DataArray\n",
+ " - file_mask : xarray.Dataset\n",
+ " Output : \n",
+ " - np.float32 or np.float64 depending on recording precision of simulation files \n",
+ " \n",
+ " \"\"\"\n",
+ " e1t=file_mask.e1t.squeeze()\n",
+ " e2t=file_mask.e2t.squeeze()\n",
+ " tmask=file_mask.tmask.squeeze()\n",
+ " t_DW=thetao.where(abs((thetao.depth-1750)<1250)*(abs(thetao.nav_lat)<30))\n",
+ "\n",
+ " # Computing Area Weights from Mask over Box\n",
+ " area_DW=e1t*e2t*tmask.where(abs((thetao.depth-1750)<1250)*(abs(thetao.nav_lat)<30))\n",
+ "\n",
+ "\n",
+ " #Returning Average Temperature on Box\n",
+ " return ((t_DW*area_DW).sum(dim=[\"nav_lat\",\"nav_lon\",\"depth\"])/area_DW.sum(dim=[\"nav_lat\",\"nav_lon\",\"depth\"]))\n",
+ "\n",
+ "\n",
+ "\n",
+ "\n",
+ "## Version 1 of ACC metric : Computes the flux assuming rigid lid (as if ssh didn’t change)\n",
+ "\n",
+ "def ACC_Drake_metric(uo, file_mask):\n",
+ " \"\"\"\n",
+ " Metric Extraction Function :\n",
+ " Antarctic Circumpolar Current Transport at the DINO equivalent of the Drake Passage (x=0).\n",
+ " Unit : Sv\n",
+ "\n",
+ "\n",
+ " Input :\n",
+ " - uo : xarray.DataArray\n",
+ " - file_mask : xarray.Dataset\n",
+ " Output : \n",
+ " - np.float32 or np.float64 depending on recording precision of simulation files\n",
+ " \n",
+ " \n",
+ " \"\"\"\n",
+ " \n",
+ " umask_Drake=file_mask.umask.isel(nav_lon=0).squeeze()\n",
+ " e3u=file_mask.e3u_0.squeeze()\n",
+ " e2u=file_mask.e2u.squeeze()\n",
+ "\n",
+ " # Masking the variables onto the Drake Passage\n",
+ " \n",
+ " u_masked=uo.isel(nav_lon=0)*umask_Drake\n",
+ " e3u_masked=e3u.isel(nav_lon=0)*umask_Drake\n",
+ " e2u_masked=e2u.isel(nav_lon=0)*umask_Drake\n",
+ "\n",
+ " # Multiplying zonal velocity by the sectional areas (e2u*e3u)\n",
+ " \n",
+ " ubar=(u_masked*e3u_masked)\n",
+ " flux=(e2u_masked*ubar).sum(dim=[\"nav_lat\",'depth'])\n",
+ "\n",
+ " #Returning Total Transport across Drake passage as a numpy scalar (unit : Sv)\n",
+ " return flux/1e6\n",
+ "\n",
+ "\n",
+ "\n",
+ "### Version 2 of ACC metric : Computes the flux assuming varying ssh, thus needing to recompute e3u variable from e3u_0\n",
+ "\n",
+ "\n",
+ "def ACC_Drake_metric_2(uo, ssh, file_mask):\n",
+ " \"\"\"\n",
+ " Metric Extraction Function :\n",
+ " Antarctic Circumpolar Current Transport at the DINO equivalent of the Drake Passage (x=0).\n",
+ " Unit : Sv\n",
+ "\n",
+ "\n",
+ " Input :\n",
+ " - uo : xarray.DataArray\n",
+ " - ssh : xarray.DataArray\n",
+ " - file_mask : xarray.Dataset\n",
+ " Output : \n",
+ " - np.float32 or np.float64 depending on recording precision of simulation files\n",
+ " \n",
+ " \n",
+ " \"\"\"\n",
+ "\n",
+ "\n",
+ " e3u_0=file_mask.e3u_0\n",
+ " e2u=file_mask.e2u\n",
+ " umask_Drake=file_mask.umask.isel(nav_lon=0)\n",
+ "\n",
+ " # Recomputing e3u, using ssh to refactor the original e3u_0 cell heights)\n",
+ " \n",
+ " ssh_u = (ssh + ssh.roll(_nav_lon=-1))/2\n",
+ " bathy_u = e3u_0.sum(dim=\"depth\")\n",
+ " ssumask = umask_Drake[:,0]\n",
+ " e3u = e3u_0*(1+ssh_u*ssumask/(bathy_u+1-ssumask))\n",
+ " \n",
+ " # Masking the variables onto the Drake Passage\n",
+ " \n",
+ " u_masked=uo.isel(nav_lon=0)*umask_Drake\n",
+ " e3u_masked=e3u.isel(nav_lon=0)*umask_Drake\n",
+ " e2u_masked=e2u.isel(nav_lon=0)*umask_Drake\n",
+ "\n",
+ " # Multiplying zonal velocity by the sectional areas (e2u*e3u)\n",
+ " \n",
+ " ubar=(u_masked*e3u_masked).sum(dim='depth')\n",
+ " flux=(e2u_masked*ubar).sum()\n",
+ "\n",
+ " #Returning Total Transport across Drake passage as a numpy scalar (unit : Sv)\n",
+ " return flux.data/1e6\n",
+ "\n",
+ "\n",
+ "\n",
+ "### Intensity of the North-Atlantic SubTropical Gyre (NASTG) computed from the local maximum of the Barotropic Stream Function (BSF)\n",
+ "\n",
+ "def NASTG_BSF_max(vo,ssh,file_mask):\n",
+ " \"\"\"\n",
+ " Metric Extraction Function :\n",
+ " Intensity of the North-Atlantic SubTropical Gyre (NASTG) which contains the Gulf-Stream Current.\n",
+ " Computed using the Barotropic Stream Function (BSF).\n",
+ " Unit : Sv\n",
+ " \n",
+ " \n",
+ " Input :\n",
+ " - vo : xarray.DataArray\n",
+ " - file_mask : xarray.Dataset\n",
+ " Output : \n",
+ " - np.float32 or np.float64 depending on recording precision of simulation files \n",
+ " \"\"\"\n",
+ "\n",
+ " e3v_0=file_mask.e3v_0.squeeze()\n",
+ " e1v=file_mask.e1v.squeeze()\n",
+ " vmask=file_mask.vmask.squeeze()\n",
+ " # Updating e3v from e3v_0 and SSH\n",
+ " ssh_v = (ssh + ssh.roll(nav_lat=-1))/2\n",
+ " bathy_v = e3v_0.sum(dim=\"depth\")\n",
+ " ssvmask = vmask.isel(depth=0)\n",
+ " e3v = (e3v_0*(1+ssh_v*ssvmask/(bathy_v+1-ssvmask)))\n",
+ "\n",
+ " # Integrating Meridional Transport (e3v*e1v*vo) along depth and X from Western boundary eastward \n",
+ " # to get Barotropic Stream Function with the \"American continent\" as reference point (BSF=0)\n",
+ " V = (vo * e3v).sum(dim='depth') # == \"Barotropic Velocity\" * Bathymetry\n",
+ " BSF=(V * e1v * ssvmask).cumsum(dim='nav_lon') / 1e6 # Integrating from the West, and converting from m³/s to Sv\n",
+ " # Selecting 0N-40N window where to search for the maximum, which will correspond to the center of rotation for the gyre\n",
+ " BSF_NASPG=BSF.where(abs(BSF.nav_lat-20)<20)\n",
+ "\n",
+ " return BSF_NASPG.max(dim=[\"nav_lat\",\"nav_lon\"])\n",
+ "\n",
+ "\n",
+ " "
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 5,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Successfully loaded dataset from ../nc_files/nc_files/DINO_00576000_restart.nc\n",
+ "Successfully loaded dataset from ../nc_files/nc_files/DINO_1y_grid_T.nc\n",
+ "Successfully loaded dataset from ../nc_files/nc_files/DINO_1y_grid_U.nc\n",
+ "Successfully loaded dataset from ../nc_files/nc_files/DINO_1y_grid_V.nc\n",
+ "Successfully loaded dataset from ../nc_files/nc_files/DINO_1y_grid_W.nc\n",
+ "Successfully loaded dataset from ../nc_files/nc_files/DINO_1m_grid_T.nc\n",
+ "Successfully loaded dataset from ../nc_files/nc_files/DINO_1m_To_1y_grid_T.nc\n",
+ "Successfully loaded dataset from ../nc_files/nc_files/mesh_mask.nc\n",
+ "\n"
+ ]
+ }
+ ],
+ "source": [
+ "# Read restart file\n",
+ "DINO_restart_path = \"../nc_files/nc_files/DINO_00576000_restart.nc\"\n",
+ "DINO_restart_file=xr.open_dataset(DINO_restart_path).rename({\"nav_lev\":\"depth\",\"y\":\"nav_lat\",\"x\":\"nav_lon\"})\n",
+ "print(f\"Successfully loaded dataset from {DINO_restart_path}\")\n",
+ "\n",
+ "#Read grid_T file\n",
+ "DINO_1y_grid_T_path = \"../nc_files/nc_files/DINO_1y_grid_T.nc\"\n",
+ "DINO_1y_grid_T=xr.open_dataset(DINO_1y_grid_T_path,decode_cf=False).rename({\"deptht\":\"depth\",\"y\":\"nav_lat\",\"x\":\"nav_lon\"})\n",
+ "print(f\"Successfully loaded dataset from {DINO_1y_grid_T_path}\")\n",
+ "\n",
+ "#Read grid_U file\n",
+ "DINO_1y_grid_U_path = \"../nc_files/nc_files/DINO_1y_grid_U.nc\"\n",
+ "DINO_1y_grid_U=xr.open_dataset(DINO_1y_grid_U_path,decode_cf=False).rename({\"depthu\":\"depth\",\"y\":\"nav_lat\",\"x\":\"nav_lon\"})\n",
+ "print(f\"Successfully loaded dataset from {DINO_1y_grid_U_path}\")\n",
+ "\n",
+ "#Read grid_V file\n",
+ "DINO_1y_grid_V_path = \"../nc_files/nc_files/DINO_1y_grid_V.nc\"\n",
+ "DINO_1y_grid_V=xr.open_dataset(DINO_1y_grid_V_path,decode_cf=False).rename({\"depthv\":\"depth\",\"y\":\"nav_lat\",\"x\":\"nav_lon\"})\n",
+ "print(f\"Successfully loaded dataset from {DINO_1y_grid_V_path}\")\n",
+ "\n",
+ "#Read grid_W file\n",
+ "DINO_1y_grid_W_path = \"../nc_files/nc_files/DINO_1y_grid_W.nc\"\n",
+ "DINO_1y_grid_W=xr.open_dataset(DINO_1y_grid_W_path,decode_cf=False).rename({\"depthw\":\"depth\",\"y\":\"nav_lat\",\"x\":\"nav_lon\"})\n",
+ "print(f\"Successfully loaded dataset from {DINO_1y_grid_W_path}\")\n",
+ "\n",
+ "#Read grid_T file\n",
+ "DINO_1m_grid_T_path = \"../nc_files/nc_files/DINO_1m_grid_T.nc\"\n",
+ "DINO_1m_grid_T=xr.open_dataset(DINO_1m_grid_T_path,decode_cf=False).rename({\"y\":\"nav_lat\",\"x\":\"nav_lon\"})\n",
+ "print(f\"Successfully loaded dataset from {DINO_1m_grid_T_path}\")\n",
+ "\n",
+ "#Read grid_T file\n",
+ "DINO_1m_To_1y_grid_T_path = \"../nc_files/nc_files/DINO_1m_To_1y_grid_T.nc\"\n",
+ "DINO_1m_To_1y_grid_T=xr.open_dataset(DINO_1m_To_1y_grid_T_path).rename({\"y\":\"nav_lat\",\"x\":\"nav_lon\"})\n",
+ "print(f\"Successfully loaded dataset from {DINO_1m_To_1y_grid_T_path}\")\n",
+ "\n",
+ "DINO_1m_To_1y_grid_T[\"time_counter\"] = DINO_1y_grid_T[\"time_counter\"]\n",
+ "\n",
+ "# Read the corresponding mesh mask file\n",
+ "DINO_mask_path = \"../nc_files/nc_files/mesh_mask.nc\"\n",
+ "DINO_mask_file=xr.open_dataset(DINO_mask_path).rename({\"nav_lev\":\"depth\",\"y\":\"nav_lat\",\"x\":\"nav_lon\"})\n",
+ "print(f\"Successfully loaded dataset from {DINO_mask_path}\\n\")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": []
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 6,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/html": [
+ "
\n",
+ "
<xarray.DataArray (time_counter: 50)> Size: 400B\n",
+ "array([37.37823257, 50.26876607, 53.99026863, 53.63378683, 52.25463586,\n",
+ " 50.67156167, 49.14597273, 47.70837255, 46.46927904, 45.37983725,\n",
+ " 44.3223467 , 43.52822327, 42.72478865, 41.94973811, 41.33527703,\n",
+ " 40.7860674 , 40.32864996, 39.95237118, 39.63370376, 39.34591195,\n",
+ " 39.07324367, 38.8113178 , 38.5579836 , 38.31391372, 38.08504955,\n",
+ " 37.87343851, 37.67977213, 37.50365954, 37.34350209, 37.19826026,\n",
+ " 37.0667913 , 36.94998701, 36.84391707, 36.74824769, 36.65981263,\n",
+ " 36.57487045, 36.49228322, 36.4110251 , 36.33137004, 36.25212382,\n",
+ " 36.17376929, 36.09568445, 36.01735414, 35.93781489, 35.8564326 ,\n",
+ " 35.77270767, 35.68743818, 35.60379858, 35.52292313, 35.44562478])\n",
+ "Coordinates:\n",
+ " * time_counter (time_counter) float64 400B -5.905e+10 ... -5.753e+10\n",
+ " depth float32 4B 5.034
37.38 50.27 53.99 53.63 52.25 50.67 ... 35.77 35.69 35.6 35.52 35.45
array([37.37823257, 50.26876607, 53.99026863, 53.63378683, 52.25463586,\n",
+ " 50.67156167, 49.14597273, 47.70837255, 46.46927904, 45.37983725,\n",
+ " 44.3223467 , 43.52822327, 42.72478865, 41.94973811, 41.33527703,\n",
+ " 40.7860674 , 40.32864996, 39.95237118, 39.63370376, 39.34591195,\n",
+ " 39.07324367, 38.8113178 , 38.5579836 , 38.31391372, 38.08504955,\n",
+ " 37.87343851, 37.67977213, 37.50365954, 37.34350209, 37.19826026,\n",
+ " 37.0667913 , 36.94998701, 36.84391707, 36.74824769, 36.65981263,\n",
+ " 36.57487045, 36.49228322, 36.4110251 , 36.33137004, 36.25212382,\n",
+ " 36.17376929, 36.09568445, 36.01735414, 35.93781489, 35.8564326 ,\n",
+ " 35.77270767, 35.68743818, 35.60379858, 35.52292313, 35.44562478])
PandasIndex
PandasIndex(Index([-59050944000.0, -59019840000.0, -58988736000.0, -58957632000.0,\n",
+ " -58926528000.0, -58895424000.0, -58864320000.0, -58833216000.0,\n",
+ " -58802112000.0, -58771008000.0, -58739904000.0, -58708800000.0,\n",
+ " -58677696000.0, -58646592000.0, -58615488000.0, -58584384000.0,\n",
+ " -58553280000.0, -58522176000.0, -58491072000.0, -58459968000.0,\n",
+ " -58428864000.0, -58397760000.0, -58366656000.0, -58335552000.0,\n",
+ " -58304448000.0, -58273344000.0, -58242240000.0, -58211136000.0,\n",
+ " -58180032000.0, -58148928000.0, -58117824000.0, -58086720000.0,\n",
+ " -58055616000.0, -58024512000.0, -57993408000.0, -57962304000.0,\n",
+ " -57931200000.0, -57900096000.0, -57868992000.0, -57837888000.0,\n",
+ " -57806784000.0, -57775680000.0, -57744576000.0, -57713472000.0,\n",
+ " -57682368000.0, -57651264000.0, -57620160000.0, -57589056000.0,\n",
+ " -57557952000.0, -57526848000.0],\n",
+ " dtype='float64', name='time_counter'))
"
+ ],
+ "text/plain": [
+ " Size: 400B\n",
+ "array([37.37823257, 50.26876607, 53.99026863, 53.63378683, 52.25463586,\n",
+ " 50.67156167, 49.14597273, 47.70837255, 46.46927904, 45.37983725,\n",
+ " 44.3223467 , 43.52822327, 42.72478865, 41.94973811, 41.33527703,\n",
+ " 40.7860674 , 40.32864996, 39.95237118, 39.63370376, 39.34591195,\n",
+ " 39.07324367, 38.8113178 , 38.5579836 , 38.31391372, 38.08504955,\n",
+ " 37.87343851, 37.67977213, 37.50365954, 37.34350209, 37.19826026,\n",
+ " 37.0667913 , 36.94998701, 36.84391707, 36.74824769, 36.65981263,\n",
+ " 36.57487045, 36.49228322, 36.4110251 , 36.33137004, 36.25212382,\n",
+ " 36.17376929, 36.09568445, 36.01735414, 35.93781489, 35.8564326 ,\n",
+ " 35.77270767, 35.68743818, 35.60379858, 35.52292313, 35.44562478])\n",
+ "Coordinates:\n",
+ " * time_counter (time_counter) float64 400B -5.905e+10 ... -5.753e+10\n",
+ " depth float32 4B 5.034"
+ ]
+ },
+ "execution_count": 6,
+ "metadata": {},
+ "output_type": "execute_result"
+ },
+ {
+ "data": {
+ "image/png": "iVBORw0KGgoAAAANSUhEUgAAAiwAAAHXCAYAAACfygGuAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjkuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8hTgPZAAAACXBIWXMAAA9hAAAPYQGoP6dpAABc9klEQVR4nO3dB3gU1doH8De9954QQkhCQktC74L0otIU4VJEUVBQ4Vrgw0JHqldFr4AgTRG4VEWa9CK9h04KJCENQgpJSJ/veQ/supu6STbZkv/veQZ2Z2fPzsxudt495z3nGEiSJBEAAACAFjPU9A4AAAAAlAcBCwAAAGg9BCwAAACg9RCwAAAAgNZDwAIAAABaDwELAAAAaD0ELAAAAKD1ELAAAACA1kPAAgAAAFoPAQuAFpkxYwYZGBho7PW7dOlCTZo00djr67MjR46I91a2nD9/Xm1l79ixo9rKBtAWCFgAapm4uDgRGF2+fJm02ejRo5UuwrIlKCio2LZz586lV155hdzc3MQ2fHwl2b59O/Xq1Ys8PT3JzMyM6tSpQ6+++ipdu3ZN5f364YcfqGHDhuL5Xl5e9NFHH1FmZqbKz//ss8/ol19+ofr165O6tGzZUpQ5duxYtZUJoG2MNb0DAFDzAcvMmTOpXr16FBoaqtWnn4OClStXKq2zs7Mrtt0XX3xB7u7u1KxZM9q3b1+p5YWFhZGDgwNNnDiRnJ2dKSEhgVatWkWtW7emU6dOUUhISJn7M2XKFFq4cKEIcriMGzdu0Pfff0/Xr18v83UV9ejRQ9RkqRMHXiNGjKD8/Hz66aef1Fo2gLZAwAIAWsvY2FhciMsTFRUlArBHjx6Ri4tLqdtNmzat2Lq3335bXPCXLl1Ky5YtK/W58fHx9J///IdGjhxJ69atk69v0KABffDBB7Rz5056+eWXVTouAKg4NAkBaMiJEyeoVatWZG5uTn5+frR8+fJSt/3111+pRYsWZGFhQY6OjjR06FCKiYkpMf/kwoUL1L59e7Gtr6+v0kWY8yj4Ndmbb74pb2ZZs2aNUllcc/Diiy+SpaWlaPbgWgVNKSgooPT09DK34WClslxdXcVxpqamlrkd18BwDQafe0Wy+xs3bqxS85e1tTVFR0fTSy+9JG7zef/vf/8rrxnq2rUrWVlZkY+PD/3222+Vfi0AXYWABUAD+ALUs2dPSkpKEvkWHDxMnz5d5FiUlJ8xatQoCggIEL/wJ02aRAcPHqQXXnih2EU2JSWF+vbtK4IbDjK45uC9994TzR6Mcy9mzZolbnO+A+c98MJlKZbRu3dv0Tzy9ddfi5wRbgrZs2dPuceVlpYmajnKWzIyMlQ6T1lZWWRrayuagThQmzBhgsrPLQuft4cPH4r3gWtYOCDq1q1bmc/JyckR/3MgqIiDHcaBYlUDsz59+pC3t7d47zgIe//990Uwye8H56ksWLCAbGxsxOeBa5UAahUJAGrcgAEDJHNzc+n+/fvydTdu3JCMjIwkxT/Le/fuiXVz585Ven5YWJhkbGystL5z587iuV9//bV8XU5OjhQaGiq5urpKubm5Yt25c+fEdqtXry62X7Iy1q1bp1SGu7u7NHjw4HKPS/b88pY33nij3LL+7//+T5oyZYq0adMmacOGDeI5/NwOHTpIeXl5JT7n4cOHYpvp06eXWXZgYKB8X6ytraUvvvhCKigoKPM5Fy5cENvPnj1baf3evXvl5ZTl8OHDYjv+vyjZsX311VfydSkpKZKFhYVkYGAgbdy4Ub7+1q1bpR4jv6f8GL/HAPoGOSwANYx/SXOC5oABA6hu3bry9Vz7wT1Ydu/eLV+3bds2KiwspCFDhoiaCRlOMOUal8OHD4teJ4o5H+PGjZPfNzU1Ffe5loVrANq2bVvu/nFzhGLeCJfBSamRkZHlPpdrZLiGpjzcS6c88+bNK9b0wvkin3/+OW3ZsqVY00xFrF69WtSq8DHx7adPn4r3xdCw9Ern5s2bU5s2bUQtBzfXcJPZzZs3xbk1MTERZVQV1/bI2NvbU2BgIIWHh4v3X4bX8WOqvB8A+gQBC0AN46YIvrhxwFEUX4wUA5a7d+9ydUuJ2zK+UBYNBDjPQRFf5Nm9e/dUCli4GanoWDDcs+bq1avlPpeboqrTv//9b/ryyy/pwIEDVQpY2rVrJ7/N5XCwyBYvXlzm87Zu3Uqvv/46vfXWW+K+kZGR6NZ89OhRun37NlUF5zIVTRjmprCS3g9er0pgCKBPELAAaDGuXeGLFeeP8MWxpNoQdSvpdRgHTuV5/Pgx5ebmlrsd54GU1D1Zlec5OTmJ11EXDsY4oXX9+vXlBixcs8LJ0hxIcpdoDiS5tosDRVlgqO7zXpX3A0CfIGABqGH8K5ovvHzRK6ror3TuPcQXJu7to8oFkcdY4UHMFGtZ7ty5o9STpjpH0h00aJCobSjPG2+8UaxnkiqePHlSbtflyuAaL04YVhUHKrJaL+5RxV2euacPAFQfBCwANYx/MXOuCg+nzt1YZXksnA9RdPAxDgCmTp0qBnrjrs2KwQYHMlzTwDUOMtztlrtHczMF49oOvs8XeFlzjSyYKa8bb2WoK4clOzub8vLyRI8YRbNnzxbHzb1mKoN7ZXE3ZkXcVMa9rrgXjqKIiAh50FhWDdjkyZNFT6F33323UvsEAKpBwAKgARyA7N27lzp16kTjx48XgQaPmNq4cWOlXBG+WM6ZM0cELXxh5URdvohzl1buAs1dkz/55BOlQICTQnlbrpHZtGmTGIKfRz+V5btwmZy0yeOzcFkcwHAyKdfiVJW6cli4uYVHrR02bJh8KH4O5ji/h4OV/v37K23PXbPv378vukGzY8eOifPGeKA3HruENW3aVHRf5hF+uSmIa7l+/vlnERzNnz9fqUxZN2c+lzI8ui0HU/x8fg6Ph3L27Flau3atUgI1AFQDTXdTAqitjh49KrVo0UIyNTWV6tevLy1btkx0VS3pz3Lr1q1Sx44dJSsrK7EEBQVJEyZMkG7fvq3Upbhx48bS+fPnpXbt2olu0z4+PtIPP/xQrLzff/9datSokegardjFWVZGSd1uuayawl16R4wYIfn7+0uWlpaSmZmZ2C/u9ivrnq1qd2rFbsR8flu2bCk5ODiIY/f09JSGDh0qXb16tViZfLxFj5nPU0hIiHgPbGxspG7dukmHDh1S6ZjK69bMZZZ0XCW9H7xf/fr1K7Ye3ZpBnxnwP9URCAFAzeKRbjm/oyIT+UHN4VGGuSs0NwV26NBB1HJxN3R14KY/7qbNo+3yNAHnzp0r1sQFoOvQJAQAUIO4WY+pM6jgprKBAweqpSwAbYWABQCgBvBUB/v371cac0dduMamusoG0BZoEgLQE2gSAgB9hoAFAAAAtB5mawYAAACth4AFAAAAtJ5eJN3yaJM8JDkPglWdw44DAACA+vDIKjzlBg96WdZs6XoTsHCw4u3trendAAAAgEqIiYkRM5PrfcAim2+ED9jW1lbTuwMAAAAq4AEPucKh6LxhehuwyJqBOFhBwAIAAKBbVEnnQNItAAAAaD0ELAAAAKD1ELAAAACA1kPAAgAAAFoPAQsAAABoPQQsAAAAoPUQsAAAAIDWQ8ACAAAAWg8BCwAAAGg9BCwAAACg9RCwAAAAgNZDwAIAAABaDwFLGTJy8unYnYe0Oyy+5t4RAAAA0M/ZmqtL1MNMGrXqLDlZmVLfph6a3h0AAIBaCzUsZfB1sRL/J2fmUlpWXk29JwAAAFAEApYyWJsZk6uNmbgdlZxZ1qYAAABQjRCwlKP+81qWyIcZ1fk+AAAAQBkQsJTD19la/B/1CDUsAAAAmoKApRx+8hoWBCwAAACagoClHL7OzwMW1LAAAABoDAIWFQOWqEcZVFgo1cR7AgAAAEUgYCmHt6MlGRsaUHZeISWkZ5e3OQAAAFQDBCzlMDEypLqOluI2Em8BAAA0AwGLCtC1GQAAQLMQsKgAibcAAACahYBFBfVdno3Fgq7NAAAAmoGApUI9hTAWCwAAgCYgYFFB/ecBS2xKFuXkF1T3ewIAAABFIGBRgYuNmZgIkYdhiU7OUuUpAAAAoEYIWFRgYGCAxFsAAAANQsBS4a7NyGMBAACoaQhYKjFEPwAAANQsBCwqQtdmAAAAzUHAUsGeQujaDAAAUPMQsKio3vOAJTkzl9Ky8qrzPQEAAIAiELCoiLs1u9maiduRyGMBAACoUQhYKgAj3gIAAGgGApYKQOItAACAZiBgqQAk3gIAAGgGApZKDB4X8RBjsQAAANQkBCwV4OtsLf6/l5xJhTyxEAAAAGh/wDJ//nwxz86kSZNK3aZLly5im6JLv3795NuMHj262OO9e/cmbVPHwYKMDQ0oO6+QEtKzNb07AAAAtYZxZZ947tw5Wr58OQUHB5e53bZt2yg3N1d+Pzk5mUJCQui1115T2o4DlNWrV8vvm5k960KsTUyMDKmuk6WYT4gXT3sLTe8SAABArVCpGpaMjAwaPnw4rVixghwcHMrc1tHRkdzd3eXL/v37ydLSsljAwgGK4nbllav5xFvksQAAAGh1wDJhwgTRpNO9e/cKP/fnn3+moUOHkpXVswu/zJEjR8jV1ZUCAwPpvffeEzUxpcnJyaH09HSlpaa7Nkdg1mYAAADtbRLauHEjXbx4UTQJVdTZs2fp2rVrImgp2hw0aNAg8vX1pYiICPrss8+oT58+dOrUKTIyMipWzrx582jmzJmkCRg8DgAAQMsDlpiYGJo4caJo1jE3N6/wi3Gg0rRpU2rdurXSeq5xkeHHOS/Gz89P1Lp069atWDlTp06ljz76SH6fa1i8vb2pJpuEMDw/AACAljYJXbhwgZKSkqh58+ZkbGwslqNHj9KSJUvE7YKCglKfm5mZKWpnxowZU+7r1K9fn5ydnSk8PLzExznfxdbWVmmpKb7Px2KJTXlKOfmlHy8AAABoqIaFazvCwsKU1r355psUFBREU6ZMKbH5Rmbz5s0i92TEiBHlvk5sbKzIYfHw8CBt42JtJiZCzMjJp+jkLApws9H0LgEAAOi9CgUsNjY21KRJE6V1nDzr5OQkXz9q1Cjy8vISeSZFm4MGDBggti3a44jzUQYPHix6B3EOy+TJk8nf35969epF2obHiOERb6/GponEWwQsAAAAWjwOS2mio6PJ0FC5pen27dt04sQJ+uuvv4ptz7UyV69epbVr11Jqaip5enpSz549afbs2Vo5Foss8ZYDlqhHmZreFQAAgFqhygELJ8aWdZ9xV2VJKnkoewsLC9q3bx/pkvrPh+iPxJxCAAAANQJzCVUh8RY1LAAAADUDAUuVRrtFkxAAAEBNQMBShcHjkjNzKS0rT93vCQAAABSBgKUSrMyMyc32WUIwBpADAACofghYqpx4i2YhAACA6oaApZKQeAsAAFBzELBUEuYUAgAAqDkIWCqJR7tlaBICAACofghYqpjDci85kwoLSx4UDwAAANQDAUsl1XGwIGNDA8rOK6T49Gw1vR0AAABQEgQslWRsZEh1nSzF7Sj0FAIAAKhWCFjU0bX5UYa63g8AAAAoAQKWKkDiLQAAQM1AwKKWrs0YPA4AAKA6IWBRw5xCUWgSAgAAqFYIWKqgvsuzHJbYlKeUk1+grvcEAAAAikDAUgXO1qZkY2ZMkkR0PzkLHy4AAIBqgoClCgwMDORzCmHEWwAAgOqDgKWKMKcQAABA9UPAUkW+z8diweBxAAAA1QcBi5rGYolC12YAAIBqg4BFTV2bMRYLAABA9UHAoqaA5XFmLqVk5qrjPQEAAIAiELBUkZWZsbxZ6Njdh1UtDgAAAEqAgEUN+jbxEP/vuhqvjuIAAACgCAQsavBSyLOA5cidh/QkO08dRQIAAIACBCxqEOhmQ34uVpSbX0gHbiaqo0gAAABQgIBFTSPe9gv2FLfRLAQAAKB+CFjU5KXgZ81Cx+48orSnaBYCAABQJwQsatLAzYYCXK0pt6CQ9t9AsxAAAIA6IWBRo5fkzUJx6iwWAACg1kPAokb9gt3F/8fvPqK0LDQLAQAAqAsCFjXyd7WhIHcbyi+UaN+NBHUWDQAAUKshYFGzfk0xiBwAAIC6IWBRs37Pewv9Hf4IcwsBAACoCQIWNavvYk2NPGyfNQtdR7MQAACAOiBgqcZall1hmFsIAABAHRCwVGMey8mIZErOyKmOlwAAAKhVELBUg3rOVtTEy5YKRLMQBpEDAACoKgQs1TyI3J8YRA4AAECzAcv8+fPFxH+TJk0qdZs1a9aIbRQXc3NzpW0kSaJp06aRh4cHWVhYUPfu3enu3bukD81CpyOT6eETNAsBAABoJGA5d+4cLV++nIKDg8vd1tbWluLj4+XL/fv3lR5fuHAhLVmyhJYtW0ZnzpwhKysr6tWrF2VnZ5Ou8na0pJA6dlQoEe1FbyEAAICaD1gyMjJo+PDhtGLFCnJwcCh3e65VcXd3ly9ubm5KtSvffvstffHFF9S/f38RAK1bt47i4uJox44dpBe9hdAsBAAAUPMBy4QJE6hfv36i6UbVAMfHx4e8vb1FUHL9+nX5Y1FRUZSQkKBUlp2dHbVp04ZOnTpVYnk5OTmUnp6utGijvs+bhc5EPaakdN2tLQIAANC5gGXjxo108eJFmjdvnkrbBwYG0qpVq+j333+nX3/9lQoLC6l9+/YUGxsrHudghSnWusjuyx4ril+bgxrZwoGQNqrjYEnN6tqTJBHtuYZB5AAAAGokYImJiaGJEyfS+vXriyXOlqZdu3Y0atQoCg0Npc6dO9O2bdvIxcVF5L9U1tSpUyktLU2+8H5pK8wtBAAAUMMBy4ULFygpKYmaN29OxsbGYjl69KhImOXbBQUF5ZZhYmJCzZo1o/DwcHGfc1pYYqLyeCV8X/ZYUWZmZiKRV3HRVrJmoXP3H1NCGpqFAAAAqj1g6datG4WFhdHly5flS8uWLUUCLt82MjIqtwwOargM7sLMfH19RWBy8OBB+Tack8K9hbh2Rtd52ltQCx+H581CGKofAACgMowrsrGNjQ01adJEaR13QXZycpKv5+YfLy8veY7LrFmzqG3btuTv70+pqam0aNEi0a357bffFo/LxnGZM2cOBQQEiADmyy+/JE9PTxowYADpA24WunA/hf68Gk9vdvDV9O4AAADod8CiiujoaDI0/KfiJiUlhd555x2RQMtdoFu0aEEnT56kRo0aybeZPHkyZWZm0tixY0VQ07FjR9q7d6/KeTLajpuFZu+6IYKWuNSnotYFAAAAVGcg8UAoOo6bkLi3ECfgams+y5Dlp+hs1GMa3b4ezXilsaZ3BwAAQKeu35hLqIZ82DVA/P/r6fsU+TCjpl4WAABALyBgqSEdA5zpxUAXyi+UaP6eWzX1sgAAAHoBAUsN+qxvQzIyNKC/biTSmcjkmnxpAAAAnYaApQYFuNnQ0FbPRuWdu/smFfLMiAAAAFAuBCw1bFL3BmRlakRXY9PojytxNf3yAAAAOgkBSw1zsTGj8S/6i9uL9t2m7LzyRwcGAACo7RCwaMCYjr7kYWdOD1Kf0qq/ozSxCwAAADoFAYsGmJsY0ae9AsXtHw9H0KOMHE3sBgAAgM5AwKIhA0K9qKmXHWXk5NO3B+5oajcAAAB0AgIWTZ14QwP6vF9DcXvD2RgKT3qiqV0BAADQeghYNKhtfSfq0ciNCgolmrcbg8kBAACUBgGLhv1fnyAyNjSgg7eS6GT4I03vDgAAgFZCwKJhfi7WNLxNXXF7zi4MJgcAAFASBCxaYGL3BmRjbkw34tNp26UHmt4dAAAArYOARQs4WpnS+88Hk1u87zY9zcVgcgAAAIoQsGiJN9rXIy97C0pIz6YfDt/V9O4AAABoFQQsWjSY3JcvPevmvOxoJF2OSdX0LgEAAGgNBCxapHcTD3olxFN0c/7of5fRNAQAAPAcAhYtM6t/Y3KzNaPIh5m0YC/GZgEAAGAIWLSMvaUpLXw1RNxec/Ie/Y2xWQAAABCwaKPODVzkY7N8uvkKpWfnaXqXAAAANAo1LFrqs74NycfJkuLSsmnmHzc0vTsAAAAahYBFS1mZGdPXr4WQgQHR1ouxtO96gqZ3CQAAQGMQsGixlvUcadwLfuL2Z9vC6FFGjqZ3CQAAQCMQsGi5f/cIoCB3G0rOzKXPt4eRJEma3iUAAIAah4BFy5kZG9HXQ0LIxMiA9l1PpG0XMdcQAADUPghYdEBjTzua1L2BuD3jj+sUl/pU07sEAABQoxCw6IhxL9SnZnXt6UlOPn265QoVFqJpCAAAag8ELDrC2MiQ/jMklCxMjOjv8GRae+qepncJAACgxiBg0SG+zlb0Wd8gcXve7luYIBEAAGoNBCw6ZkRbH+rd2J1yCwppwvqL9DgzV9O7BAAAUO0QsOgYAwMDWvhasKhteZD6lCZuvCRmdwYAANBnCFh0kK25CS0d0ZzMTQzp+N1HtOTgXU3vEgAAQLVCwKKjgtxtad6gpuL2kkN36cjtJE3vEgAAQLVBwKLDBjarI2Z15sFvJ226TLEpWZreJQAAgGqBgEXHTXu5EQXXsaPUrDwav/4i5eQXaHqXAAAA1A4Bix4M3f/j8OZkb2lCV2PTaNbOG5reJQAAALVDwKIH6jhY0revh5KBAdH6M9G07WKspncJAABArRCw6Ikuga70YdcAcfuz7WF0KyFd07sEAACgHQHL/PnzxbggkyZNKnWbFStWUKdOncjBwUEs3bt3p7NnzyptM3r0aFGO4tK7d++q7Fqt9GG3AHqhgQtl5xXSe79epPTsPE3vEgAAgGYDlnPnztHy5cspODi4zO2OHDlCw4YNo8OHD9OpU6fI29ubevbsSQ8ePFDajgOU+Ph4+bJhw4bK7lqtZWRoIJqGvOwtKOpRJn38P0ySCAAAtThgycjIoOHDh4vaE641Kcv69etp/PjxFBoaSkFBQbRy5UoqLCykgwcPKm1nZmZG7u7u8qW8cqFkjlamIgnX1MiQ9t9IpFl/3iCJ+z0DAADUtoBlwoQJ1K9fP9G8U1FZWVmUl5dHjo6OxWpiXF1dKTAwkN577z1KTk4utYycnBxKT09XWuAfId72tHhIiLi95uQ9WnE8EqcHAABqV8CyceNGunjxIs2bN69SLzhlyhTy9PRUCna4OWjdunWi1mXBggV09OhR6tOnDxUUlDymCL+2nZ2dfOFmJlD2SognfdGvobj91e5b9Ptl5SY4AAAAXWIgVaC9ICYmhlq2bEn79++X56506dJFNPd8++23KiXpLly4UNSmlJX7EhkZSX5+fnTgwAHq1q1biTUsvMhwDQsHLWlpaWRra6vq4dQKPC7Lqr+jyMTIgNa+2Zra+ztrepcAAADk12+ueFDl+l2hGpYLFy5QUlISNW/enIyNjcXCtSFLliwRt0urEWGLFy8WActff/1VbqJu/fr1ydnZmcLDw0t8nPNd+MAUFygZ17L0C/agvAKJxv1ygW7Go/kMAAB0T4UCFq7tCAsLo8uXL8sXrnHhBFy+bWRkVOLzuFZl9uzZtHfvXrF9eWJjY0UOi4eHR0V2D0pgaGhAX78WQq19HelJTj6NXn2WHqQ+xbkCAAD9DVhsbGyoSZMmSouVlRU5OTmJ22zUqFE0depU+XM4J+XLL7+kVatWUb169SghIUEs3NOI8f+ffvopnT59mu7duyfyWPr370/+/v7Uq1cvdR9vrWRuYkQrRrakBm7WlJieQ2+sOktpWRijBQAAavFIt9HR0WIcFZmlS5dSbm4uvfrqq6LGRLZwExHjWpmrV6/SK6+8Qg0aNKAxY8ZQixYt6Pjx46LpB9TDztKE1rzZmtxtzSk8KYPeWXeesvMwUSIAAOhh0q0+JO3Udjxk/2tLT4nmob5N3emHYc1FsxEAAIDeJN2C7gtyt6Xlo1qIgeV2hyVgYDkAANAJCFhqofZ+zkoDyy3Yexuj4QIAgFZDwFKLB5ab1b+xuL3saATN33sLQQsAAGgtBCy12Kh29eRBy/KjkTRvD4IWAADQTghYajkOWmY/D1p+OhZJX+2+iZoWAADQOghYgEa2q0dzBjwbR2fF8SiaswtBCwAAaBcELCCMaOtDcwc+C1p+PhFFs/9E0AIAANoDAQvIDW/jQ18NbCpu84SJs/68geYhAADQCghYQMm/2tSl+YOeBS2r/75HM3ciaAEAAM1DwALFDG1dlxYMbkoGBs/GaZn+x3XUtAAAgEYhYIESvd6qLi0YFCyClnWn7tN/D4fjTAEAgMYgYIFSDWnlTbP7P0vE/Xr/HTp65yHOFgAAaAQCFii399DQVt7EU2RO3HiJYh5n4YwBAECNQ8AC5ZrxSmMKrmNHqVl59O6vFyg7rwBnDQAAahQCFiiXuYkRLR3RghytTOl6XDp9seMaknABAKBGIWABlXjZW9D3w5qRoQHRlguxtP5MNM4cAADUGAQsoLIO/s40uXeQuD1z53W6GJ2CswcAADUCAQtUyLgX6lOfJu6UVyDR+F8v0sMnOTiDAABQ7RCwQIUYGBjQotdCyM/FihLSs+n93y5SfkEhziIAAFQrBCxQYdZmxrR8ZEuyMjWiM1GPaf6eWziLAABQrRCwQKX4u1rT4tdCxO2VJ6Jo55U4nEkAAKg2CFig0vo09aBxneuL21O2XqU7iU9wNgEAoFogYIEq+bRnIHXwd6Ks3AJ695cL9CQ7D2cUAADUDgELVImxkSEtGdqMPO3MKfJRJn2y+QoGlQMAALVDwAJV5mRtRj+OaEGmRoa073oiLT8WibMKAABqhYAF1CLU217MOcQW7r1FJ8Mf4cwCAIDaIGABtRnW2ptebVGHCiWiDzZcorjUpzi7AACgFghYQK2Dys0Z0IQaedhScmYujV9/kXLyMbMzAABUHQIWUPvMzstGtCA7CxO6HJNKs/+8gTMMAABVhoAF1K6ukyV9OzSUDAyIfj0dLWZ3BgAAqAoELFAtXgx0pYndAsTtz7eH0fW4NJxpAACoNAQsUG0+7BpALwa6UE5+Ib376wVKy8KgcgAAUDkIWKDaGBoa0Devh5K3owXFPH5KkzZdokLuQgQAAFBBCFigWtlbmtLS4S3IzNiQDt9+SN8fCscZBwCACkPAAtWuiZcdzR3YVNz+9uAdOnQrEWcdAAAqBAEL1AgeUG5kWx+SJKKJGy9T1KNMnHkAAFAZAhaoMV++1Iha+jjQk+x8GvfLecrMycfZBwAAlSBggRpjamxIPw5vTi42ZnQnMYMmb72KmZ0BAEAlCFigRrnamtPS4c3J2NCAdl2Np58wszMAAKgAAQvUuJb1HGn6y43E7QV7b9HfmNkZAACqM2CZP3++mPBu0qRJZW63efNmCgoKInNzc2ratCnt3r1b6XFJkmjatGnk4eFBFhYW1L17d7p7925Vdg203Ii2PvKZnd//7SLFpmRpepcAAEAfA5Zz587R8uXLKTg4uMztTp48ScOGDaMxY8bQpUuXaMCAAWK5du2afJuFCxfSkiVLaNmyZXTmzBmysrKiXr16UXZ2dmV3D3RkZuemXnaUkpUnRsLNzsPMzgAAoMaAJSMjg4YPH04rVqwgBweHMrf97rvvqHfv3vTpp59Sw4YNafbs2dS8eXP64Ycf5LUr3377LX3xxRfUv39/EQCtW7eO4uLiaMeOHZXZPdClmZ1HtiBHK1O69iCdPt9+DUm4AACgvoBlwoQJ1K9fP9F0U55Tp04V245rT3g9i4qKooSEBKVt7OzsqE2bNvJtisrJyaH09HSlBXSTl70F/TCsGRkaEG29GEu/nL6v6V0CAAB9CFg2btxIFy9epHnz5qm0PQcjbm5uSuv4Pq+XPS5bV9o2RfFrc1AjW7y9vSt6GKBF2vs709Q+DcXtWTtv0Ll7jzW9SwAAoMsBS0xMDE2cOJHWr18vEmg1ZerUqZSWliZfeL9At73dyZdeDvGk/EKJ3vv1At1Pxki4AABQyYDlwoULlJSUJHJQjI2NxXL06FGRMMu3CwqKJ026u7tTYqLy3DF8n9fLHpetK22boszMzMjW1lZpAd1Pwl0wuCk18rClRxm59Maqs/QoI0fTuwUAALoYsHTr1o3CwsLo8uXL8qVly5YiAZdvGxkZFXtOu3bt6ODBg0rr9u/fL9YzX19fEZgobsM5KdxbSLYN1A6Wpsa05s1WVMfBgu4lZ9Fba85h+H4AAKh4wGJjY0NNmjRRWrgLspOTk7jNRo0aJZpsZLgJae/evfT111/TrVu3aMaMGXT+/Hl6//33xeOycVzmzJlDf/zxhwiIuAxPT0/R/Rlq30i4695qTQ6WJnQ1No3Gr79IeQWFmt4tAADQt5Fuo6OjKT4+Xn6/ffv29Ntvv9FPP/1EISEhtGXLFtFdWRbgsMmTJ9MHH3xAY8eOpVatWolu0xzkaDJPBjSnvos1rRrdiixMjOjonYc0BXMOAQDUegYSD4Si47gJiXsLcQIu8ln0x6FbifTOugtUUCjR+C5+NLl3kKZ3CQAANHT9xlxCoLW6BrnRvIFNxe0fj0TQ2pP3NL1LAACgIQhYQKsNaeVNH/doIG7P2Hmddof909wIAAC1BwIW0Hrvd/Wn4W3qEjdeTtp4mU5HJmt6lwAAoIYhYAGtxz3JZvVvQj0buVFuQSG9s+483UrAdAwAALUJAhbQCUaGBrRkWDNq6eNAT7LzxcBydxKfaHq3AACghiBgAZ2a3XnlGy2pgZs1Jabn0KtLT6J5CACglkDAAjrF3tKUNo1tRy18HCg9O59G/XyW/rwap+ndAgCAaoaABXSOg5UprX+7DfVq/Cyn5f3fLtHK45Ga3i0AAKhGCFhAZ5uHfhzegka3ryfuz9l1k2btvEGFhTo/DiIAAJQAAQvodCLu9Jcb0Wd9n42Au+rvKHp/w0XKzis+azgAAOg2BCyg812ex77gR98NDSUTIwPaHZYg8lpSs3I1vWsAAKBGCFhAL/QP9aK1b7UmG3NjOnvvMb267BTFpmRpercAAEBNELCA3mjv50yb321HHnbmFJ6UQYN+PEmXolM0vVsAAKAGCFhArwS529K28e0pyN2Gkp7k0OvLT9OGs9Ga3i0AAKgiBCygdzzsLERNi6zb89RtYfR/W69STj6ScQEAdBUCFtBLNuYmtGxEC/q0VyAZGBBtPBdDQ5adorjUp5reNQAAqAQELKDXPYgmvOhPa99sTfaWJnQlNo1e/v4EnYx4pOldAwCACkLAAnrvhQYutPP9jtTY05aSM3Np5M9nacWxSJIkDDIHAKArELBAreDtaElb32tPg5p5UUGhRHN336T3N1yizJx8Te8aAACoAAEL1Krh/L8eEkKz+jcmY0MD2nU1ngb++LfoAg0AANoNAQvUuryWUe3q0caxbcnFxozuJGbQS98fp3Wn7qGJCABAiyFggVqpZT1H2vVBR+oU4EzZeYU07ffr9Mbqc5SYnq3pXQMAgBIgYIFay9XWXPQgmvFyIzIzNqRjdx5Sr2+P0e6weE3vGgAAFIGABWo1Q0MDGt3Bl3Z92JGaeNlSalYejV9/kT7adJnSs/M0vXsAAPAcAhYAIvJ3taFt73Wg91/0J0MDom2XHlCfb4/T6chknB8AAC2AgAXgOVNjQ/qkV6AY1r+uoyU9SH1Kw1acpq9238Sw/gAAGoaABaCIFj6OtGdiJxraypt4bLmfjkVSr2+O0b7rCehJBACgIQhYAEpgZWZM8wcH04pRLUX353vJWTTulws09KfTdO1BGs4ZAEANM5D0YHzy9PR0srOzo7S0NLK1tdX07oCeycjJp2VHImjF8UjKyS8UkykOalZHTKzobmeu6d0DAKgV128ELAAq4pyWRXtv0Y7LceK+hYkRjX2hPo3rXJ8sTY1xHgEAKggBC0A1uhyTSnP+vEHn76eI+262ZvRpryAxTxF3kwYAANUgYAGoZtySujssgebtuUmxKU/FukA3GxrT0ZdeCfUU8xYBAEDZELAA1JDsvAJae/Ie/XAonJ48n/nZ2dqUhrfxoRFtfUTCLgAAlAwBC0ANS3uaRxvPRovgJS7t2XxEpkaG1D/Uk8Z08qUgdySDAwAUhYAFQEPyCgpp77UE+vlElMh1keng7ySai7o0cEWeCwDAcwhYALTAhfsptOpEFO25Fk+FzwcP8HW2Esm5A5p5kbejpaZ3EQBAoxCwAGiR2JQs0VS08WyMPM+FtfZ1FMFL32APsjU30eg+AgBoAgIWAC2UmZNPe64l0PZLsXQyIlkM+8/MjA2peyM3Eby80MCFTIwwADUA1A7pGDgOQLvFpz2lHZfiaNvFWLqblCFf72RlSi+HeIqu0c287cmAh9UFANBTFQlYKvRTbunSpRQcHCwK5aVdu3a0Z8+eUrfv0qWL+MItuvTr10++zejRo4s93rt374rsFoDO8bCzoPe6+NFf/36B/vygI73VwVd0h07OzKU1J+/RoB9PUscFh2ne7pti7iI9mEEDAKBKKjQ0/86dO8nIyIgCAgLEF+jatWtp0aJFdOnSJWrcuHGx7R8/fky5ubny+8nJyRQSEkIrV64UgQrj/xMTE2n16tXy7czMzMjBwUHlg8BcQqAP8gsK6fjdR/T75Qe0/0YiZeYWyB+r52RJLwV7itqXQHcbje4nAIBONgk5OjqKoGXMmDHlbvvtt9/StGnTKD4+nqysrOQBS2pqKu3YsaPS+4CABfRxQLrDt5Loz6vxdPBWImXnFcofC3C1ftZsFOJJ9Zyf/R0BAOiiily/Kz1jW0FBAW3evJkyMzNF05Aqfv75Zxo6dKg8WJE5cuQIubq6ilqVrl270pw5c8jJyanUcnJycsSieMAA+oSH9u/T1EMsnKx74GaiCF6O3n4ocl7+s/+OWJrVtaeBzbyoX1MPcrLGqLoAoL8qXMMSFhYmApTs7Gyytram3377jfr27Vvu886ePUtt2rShM2fOUOvWreXrN27cSJaWluTr60sRERH02WefiXJPnTolmp9KMmPGDJo5c2ax9apEaAC6PqIuNxf9cSWOTtx9KB/fxdjQgDo3cBHju3Rv6EYWppjLCABqeZMQ56RER0eLwrds2SLyUY4ePUqNGjUq83njxo0TQcjVq1fL3C4yMpL8/PzowIED1K1bN5VrWLy9vRGwQK2S9CSb/rwSTzsuP6CrsWny9VamRtS7iYeoeWnn50RGmEEaALRUjeawdO/eXQQYy5cvL3Ubbjby9PSkWbNm0cSJE8st08XFRTQLcZCjCuSwQG0XnpQhknW3X3ognz2aedqZ0+ut6tKQVnWIeyYBANS6HBaZwsJCpdqOknCuC28zYsSIcsuLjY0VvYk8PDyqumsAtYa/qzV93DOQPurRQEwJwIEL57zwRIzfHLhD3x28Qy8GutKw1nWpS6ALGWNwOgDQMRWqYZk6dSr16dOH6tatS0+ePBH5KwsWLKB9+/ZRjx49aNSoUeTl5UXz5s1Tel6nTp3Ees5XUZSRkSFyUQYPHkzu7u4ih2Xy5MmibM6V4e7NqkANC0DJPY32XU+g385E05mox/L1brZm9HpLbxrSypvqOGA+IwDQwxqWpKQkEZRwt2R+AR5EThasMM5tMTRUHovu9u3bdOLECfrrr7+KlcdJtZzTwuO5cNdmbjbq2bMnzZ49W+VgBQBK72nUP9RLLBEPM2jTuRjaciGWEtNzaMmhcPr+cDi9EOBCw9vUpW4N3ZDrAgBarco5LNoANSwAqsnJLxC9jDacjaa/w5Pl632cLGl0+3r0WktvsjarcksxAIBKMJcQAJTr3qNM2nAuWswizd2lmY25MQ1t5U1vtK+H5iIAqHYIWABAZVm5+bT14gNafSKKIh9linXcE7p3E3ca09GXmtd1wCSMAFAtELAAQIUVFkp05E4S/XwiSqm5KMTbnt7u6Et9mrijdxEAqBUCFgCoklsJ6bTqRBTtuBxHufmF8jyXcS/40eAWXmRmjJF0AaDqELAAgFo8ysihX0/fp7Un71FK1rM8F1cbM3q7ky/9q40PEnQBoEoQsACA2vNcNpyNoZXHIyk+LVuss7MwoTfa+dDoDr7kaGWKMw4AFYaABQCqBTcP7bj0gJYdjZAn6FqYGIkRdN95wZcw/D8AVAQCFgCoVgWFkhhF98cj4XTtQbpYZ2JkQIOb16H3uviRj5MV3gEAKBcCFgCoETzu5PG7j+i/h8Plw/9zl2geXXd8Fz8KcLPBOwEApULAAgA17vy9x/TD4XA6cvuhuG9gQKIr9IQX/amxpx3eEQAoBgELAGhMWGwa/XD4Lu27nihf1zXIVQQuLXwc8M4AgBwCFgDQuNsJT0SOy84rcVT4fMay9n5OInDh/w24CgYAarX0CszWjMkPAaDa5yxaeiSCtl6MpfznkUtIHTuRnNuzkTsZctILANRK6QhYAEDbPEh9Sj8djaBN52MoO+/Z6Ln1Xazo3Rf8aEAzLzI1NtT0LgJADUPAAgBaKzkjh9acvCdGz03Pzhfr3G3Nxei5PJ6LlZmxpncRAGoIAhYA0HoZOfm04Uw0rTwRSYnpOf+Mntu+nhhB18naTNO7CADVDAELAOiMnPwC2n7xAS0/FklRz0fPNTM2pEHN69CYjvXI3xVjuQDoKwQsAKCzo+fysP9XY9Pk67sEutCYjr7U0d8ZPYsA9AwCFgDQ6dFzz91LERMt7r+ZSNLzLtFB7jb0VkdfeiXEk8xNjDS9mwCgBghYAEAv3E/OpNV/36P/nY+hrNwCsc7Z2pRGtq1Hw9vWJWfkuQDoNAQsAKBX0p7m0caz0aJ3UXxatlhnamRI/YI9aGQ7H2rmbY/mIgAdhIAFAPRSXkEh7bmWQD8fj6QrCnkuTbxsaVTbevRKKJqLAHQJAhYA0HtXYlJp3an7tPNqHOXmF8q7RQ9pWYdGtPUhHycrTe8iAJQDAQsA1BqPM3Np8/kY+vXMfYp5/FSs42mKOjdwoZFtfahLoCsZYfh/AK2EgAUAamW36KN3kkSty9E7D+W9izztzOn1VnVpSKs65GFnoeG9BABFCFgAgGp776JfT9+nLRdiKSUrT6zjSpauQa70rzZ1qXMD1LoAaAMELAAARJSdVyAGo/vtTDSdiXosPydc6zKklTe93sqbUOsCoDkIWAAAioh4mCHmLtp6UbnW5cVAVxraui69GOhCxkaYMRqgJiFgAQCoYK2Li40ZvdqiDg1p6U2+zuhhBFATELAAAKggPClDjKK79UIsJWfmyte38XWkoa29qU8TD0wDAFCNELAAAFQAj+Ny6FYibToXI3oYFT7vYWRjbkwDQr1ErksTLzucUwA1Q8ACAFBJcalPRe8irnmJTXk2rgtr7Gkrmov6h3qSvaUpzi+AGiBgAQCoosJCiU5GJNOm8zG071oC5RY8G03X1NiQejV2FyPqdvBzJkMMSgdQaQhYAADUKCUzl36//IA2nY+lm/Hp8vVe9hb0Wss6Ilm3joMlzjlABSFgAQCoBpIk0fW4dJHrsuPyA3qSnS+fCqCjvzO91tKbejZyQ6IugIoQsAAA1FD3aM51+Ts8Wb7elhN1m3mJfBck6gKUDQELAEANinmcJSZg5GTduLRs+fpGHpyoW4f6h3qRgxUSdQGKQsACAKChCRhPRjwSTUZ/XU/8J1HXyJB6NHYTtS7cdITZowGeQcACAKBhqVmcqBsnmow470XGw85cJOny4uOEEXWhdktPTyc7OztKS0sjW1vbMrc1kDiLrBYdMABATbv2IE00F22/9IDSnj6bx4i1re8oal14RF0LUyO8MVDrpFfg+l2hmb6WLl1KwcHBolBe2rVrR3v27Cl1+zVr1pCBgYHSYm5urrQNx0vTpk0jDw8PsrCwoO7du9Pdu3crslsAAFqNk29nvNKYznzWjX74VzN6oYGL6Fl0OvIxffS/K9Rq7gGaui2MLkWniO9EACjOmCqgTp06NH/+fAoICBB/VGvXrqX+/fvTpUuXqHHjxiU+hwOb27dvy+9z0KJo4cKFtGTJElGWr68vffnll9SrVy+6ceNGseAGAECXmZsY0UvBnmLhEXV5DqPNF2Ip+nEWbTgbLZYAV2sxtsvAZnXEhIwAoKYmIUdHR1q0aBGNGTOmxBqWSZMmUWpqaonP5Zf29PSkjz/+mD755BOxjquF3NzcxHOHDh2q0j6gSQgAdHlEXZ41mnsZ7b4WT9l5zxJ1OTH3xUBX0cvoxSBXMjGqUIU4QO1uElJUUFBAGzdupMzMTNE0VJqMjAzy8fEhb29vURtz/fp1+WNRUVGUkJAgmoFkeMfbtGlDp06dKrXMnJwccZCKCwCALuKh/dv5OdF/Xg+ls593p68GNqVmde1Fj6MDNxNp7C8XqN28gzR31w26k/hE07sLoBtNQiwsLEwEKNnZ2WRtbU3bt2+nRo0albhtYGAgrVq1SuS9cPS0ePFiat++vQhauHmJgxXGNSqK+L7ssZLMmzePZs6cWdFdBwDQarbmJvSvNnXFcjfxiUjU3XrxAT3KyKEVx6PEEuJtT6+1qEMvh3iSnYWJpncZQHubhHJzcyk6OloEIFu2bKGVK1fS0aNHSw1aFOXl5VHDhg1p2LBhNHv2bDp58iR16NCB4uLiRNKtzJAhQ0Suy6ZNm0qtYeFFhmtYuAYHvYQAQN/kFRTS0dsPRffoQ7eSKL/w2Ve22fNJGDnfpb0fxnYB/W8SqnANi6mpKfn7+4vbLVq0oHPnztF3331Hy5cvL/e5JiYm1KxZMwoPDxf33d3dxf+JiYlKAQvfDw0NLbUcMzMzsQAA6DvOXeneyE0sXNOy49IDEbzcScygP67EicXTzpwGNX82tks9Z4ztAvqpyllchYWFSrUd5eW9cJOSLDjhXkEctBw8eFAp2jpz5kyZeTEAALWRs7UZvd2pPu2b9AL98X4HGtnWR8xdxNMB/HA4nLosPkKvLTtJ/zsXQxk5zyZmBNAXFaphmTp1KvXp04fq1q1LT548od9++42OHDlC+/btE4+PGjWKvLy8RI4JmzVrFrVt21bUyHBPIe5NdP/+fXr77bfF49zsw72I5syZI7pKy7o1c8+hAQMGVMfxAgDoPP7uDK5jL5bP+zUUybmbz8fS8bsP6dy9FLFM/+M69WnqLmpd2vo6ieRegFoTsCQlJYmgJD4+XrQ5cTItBys9evQQj3Nui6HhP5U2KSkp9M4774gEWgcHB9GExHkrivkukydPFj2Nxo4dK4Kajh070t69ezEGCwBABcd2SUjLFqPpbr4QQ5EPM2nbxQdi8bK3oMHNvUSzEZqMQFdhaH4AAD3DfSkuxaSKWpc/r8bRk+x/moda1XMQtS59m3qQjTl6GYFmYS4hAAAQsvMK6K8biWJUXW4yet7JiMxNDKl3Y3ca3AK9jEBzELAAAEAxsiajrRdjKTwpQ77e3dacBjTzEs1GAW42OHNQYxCwAABAmU1GV2LTRK0Ld4tWnEG6qZcdDWzmRa+EeopeSQDVCQELAACoJCe/gA7fShIj6h5WGJiO5zLq0sCFBjb3ou4N3URyL4C6IWABAIAKe5yZK5J0OXi5EvPPpLU25sbUr6kH9Q/1oja+jugiDWqDgAUAAKqEc1x4VF3OeXmQ+lS+3sPOnF4J8RTBS0MPGzEmDEBlIWABAAC1KCyU6EzUY/r98gPaFRav1EW6gZu1CFz6h3pSHQdLnHGoMAQsAABQTfkuD0XwcvBWEuXmF8ofa+njQP2beVHfJu7khGRdUBECFgAAqFbcs2jftQTRZHQ6Kpmk5+O7cLJuR39n0WzUs7EbBqeDMiFgAQCAGh3fZeeVOPr9ygO69iBdvt7U2JC6BrqKLtJdg1zR0wiKQcACAAAaEfkwg3Zeiac/rjygiIeZ8vXWZsbUs5EbvRziSR38nUUwA5Ceni7mJkxLSyNbW9syTwjmEgIAgGoZnO5m/BMxMB3Xvij2NLKzMBHTAvQL9qD2fk5kbITgpbZKR8ACAADaFLxcjE4VgcufV+PpUUaO/DEHSxPq3cSDXgr2EGO8IHipXdIRsAAAgDYqKJTobNRjMUDd3msJlJyZK3/M2dqUejdxp35NPam1r6NI4AX9lo6ABQAAtF1+QSGdjnxMu8LiaM+1BErN+mdOI57HqHcTN+rblGtenBC86CkELAAAoFPyCgrpZEQy7Xpe85KuMEAd17z05JwXEbyg2UifIGABAACdxQPSnYx4RLvD4umvG4lKNS9OVs+Cl75N3aldfSTs6joELAAAoDc1L6cikkXwsu96AqUoBC/2libUo6GbyHvpGOBMZsaYUVrXIGABAAC9DF7OiJyXZ8ELzy6tOM4LD07HwUuXQBeyNDXW6L6CahCwAACA3ifsnruXIgIXznlJSM+WP2ZmbEidG7hQn6bu1DXITYz7AtoJAQsAANSqGaUvx6aKuY24t1H04yz5Y8aGBtTOz4l6NXYXI+262pprdF9BGQIWAACo1SPs7r0WL4KXu0kZ8scMDIiaeduL4IWXes5WGt1XIAQsAAAAsrmN9l1PFE1Hl2NSlU5KoJuNmFGag5fGnrZkwBEN1CjUsAAAAJQwq/T+GwkigDkdmUz5hZL8MS97C+re0JV6NHKnNvUdyQTzG9UIBCwAAABlSMvKo4O3ntW8HLvziJ7mFcgfszE3phcDOXhxEz2ObMyRtFtdELAAAACoKDuvgE7cfUT7bySKIOZRxj/dpU2MDKhtfSeRsNutoRt52lvgvKoRAhYAAIBKTs54OSZFjLDLAUzkw0ylxxt52IqmIw5emnrZkSEmaKwSBCwAAABqEPEwQwQuvFyMTiHpn7QXcrUxo24cvAS5UQd/Z7IwxUi7FYWABQAAQM2SM3Lo8O2HdPBmIh2785Ayc//JezE3MaQOfs6i5uXFIBfysEPTkSoQsAAAAFSjnPwCOh35WAQvB24kUlzaPyPtsoYettQtyJVeDHKlUG97MkLTUYkQsAAAANTwYHUcvBy6nSTGe1FsOnK0MqUuDVxE8PJCAxdMFaAAAQsAAIAGm46O3H4oghduOnqSnS9/jGtaWvg4iO7S3HU6yN2mVg9Yl56eTnZ2dpSWlka2trZlbmsgcWhYiw4YAACgJmeYvnA/hQ7dShJLuMJUAczd1lwEL7xw4m5tG/MlHQELAACA9ol5nEVHbieJ5N2TEY8oO69QaaLGlvW49sVVBDA8dYC+176kI2ABAADQ/gHrzkQ9psO3kujonYcU9Uh5zBc3WzPq3MCFOjdwpY7+zmRnqX+1LwhYAAAAdMy9R5ny2hee6ygn/5/aF0MDEr2NOHjpHOgiBq3Th55HCFgAAAB0vPblbNRjUfPCS3iR3BcHSxPqGOBCnQKcxaKr474gYAEAANAjD1Kfih5HR28/pL/DH9GTnH96HrEAV2vqxAFMA2dq4+tIlqbGpAsQsAAAAOhxz6PLMal0/M5DOnb3EV2NTaVChf6+pkaGInlXBDABzmL+I22d8wgBCwAAQC2RmpVLJyOS6fjdh3TsziNRG6OIB65r7+ckEne567S3oyXpYsBiWJGCly5dSsHBwaJQXtq1a0d79uwpdfsVK1ZQp06dyMHBQSzdu3ens2fPKm0zevRo0W1Lcendu3dFdgsAAKDWsrc0pb5NPWjeoGA6MeVFOvRxZ5r5SmMxq7SVqRE9zsylP6/G0/9tC6NOCw9T50WH6fPtYbQnLF4EO7qiQgPH7dy5k4yMjCggIEAMRbx27VpatGgRXbp0iRo3blxs++HDh1OHDh2offv2ZG5uTgsWLKDt27fT9evXycvLSx6wJCYm0urVq+XPMzMzEwGOqjBwHAAAQOnNRyfuPhK5L5diUqlAof2Ih3nhHkft/bj2xYla1XMkcxMj/WwScnR0FEHLmDFjyt22oKBABCI//PADjRo1Sh6wpKam0o4dOyq9DwhYAAAAyvckO4/ORD6mE+GPxFK09xHnvzT3sRczT7f3d6aQOnZkbFShxphqu35XOo2Yg4/NmzdTZmamaBpSRVZWFuXl5YkgR9GRI0fI1dVVBDNdu3alOXPmkJOTU6nl5OTkiEXxgAEAAKBsPPR/90ZuYmEJadmi5uXviEd0MjyZEtKzxSzUvHy9/w5ZmxmLXkccvHANjCZH361wDUtYWJgIULKzs8na2pp+++036tu3r0rPHT9+PO3bt080CXETEdu4cSNZWlqSr68vRURE0GeffSbKPXXqlGh+KsmMGTNo5syZxdZjLiEAAIDK4XAg8lEmneQAJjyZTkUmU9rTPPnjdhYmdPHLHmodsK5am4Ryc3MpOjpaFL5lyxZauXIlHT16lBo1alTm8+bPn08LFy4UtSmcuFuayMhI8vPzowMHDlC3bt1UrmHx9vZGwAIAAKAmnOtyMz79eQ1MMjlbmdJ/Xg8lnc1h4Z4/HGAsX7681G0WL14smnk4CGnZsmW5Zbq4uIjtx40bp9I+IIcFAABA99RIDotMYWGhUm1HUVyrMnfuXNEUpEqwEhsbS8nJyeTh4VHVXQMAAAA9UaGAZerUqdSnTx+qW7cuPXnyROSvcBMPByOMe/5wd+V58+aJ+9yNedq0aWK7evXqUUJCgljPOSq8ZGRkiFyUwYMHk7u7u8hhmTx5Mvn7+1OvXr2q43gBAABA3wOWpKQkEZTEx8eLKhzOReFgpUePHuJxzm0xNDRUGmiOc15effVVpXKmT58uEmc5qfbq1atiPBfu2uzp6Uk9e/ak2bNni7FYAAAAANSSw6INkMMCAACge6ptaH4AAAAATUDAAgAAAFoPAQsAAABoPQQsAAAAoPUQsAAAAIDWQ8ACAAAAWg8BCwAAAGg9BCwAAACg9RCwAAAAgNZDwAIAAABar8qzNWsD2ewCPMQvAAAA6AbZdVuVWYL0ImDhmaOZt7e3pncFAAAAKnEd5zmF9H7yw8LCQoqLiyMbGxsyMDAgbYocOYiKiYkpd1InXYdj1U+16X2tbceLY9VP6Tr2GeYQhIMVT09PMjQ01P8aFj7IOnXqkLbiD40ufHDUAceqn2rT+1rbjhfHqp9sdegzXF7NigySbgEAAEDrIWABAAAArYeApRqZmZnR9OnTxf/6Dseqn2rT+1rbjhfHqp/M9PgzrBdJtwAAAKDfUMMCAAAAWg8BCwAAAGg9BCwAAACg9RCwAAAAgNZDwFIF9erVEyPrKi7z588v8zkRERE0cOBAcnFxEYP6DBkyhBITE5W2efz4MQ0fPlw8bm9vT2PGjKGMjAzStWNNSEigkSNHkru7O1lZWVHz5s1p69atStvcuXOH+vfvT87OzuJ4O3bsSIcPHyZNq67jZbt27aI2bdqQhYUFOTg40IABA0hfj5Xl5ORQaGioKPfy5cukb8d679498Tfq6+sr3lM/Pz/RSyM3N5f08X3Vh+8nfs+Kbi9bNm/eLN/u3Llz1K1bN3Gc/Lfaq1cvunLlCunjsbI1a9ZQcHAwmZubk6urK02YMIG0CvcSgsrx8fGRZs2aJcXHx8uXjIyMUrfnx+rXry8NHDhQunr1qlj69+8vtWrVSiooKJBv17t3bykkJEQ6ffq0dPz4ccnf318aNmyYTh0r69Gjhzi2M2fOSBEREdLs2bMlQ0ND6eLFi/JtAgICpL59+0pXrlyR7ty5I40fP16ytLQU5evj8W7ZskVycHCQli5dKt2+fVu6fv26tGnTJkkfj1Xmww8/lPr06cO9EaVLly5J+nase/bskUaPHi3t27dPPP77779Lrq6u0scffyzp4/uqD99P+fn5StvyMnPmTMna2lp68uSJ2Ib/d3R0FO/trVu3pGvXrkmDBw+W3NzcpNzcXEmfjpV9/fXXkqenp7R+/XopPDxcfCfzZ1mbIGCp4gfnm2++UXl7/kLjP/60tDT5utTUVMnAwEDav3+/uH/jxg3xxX7u3Dn5NvyFyNs8ePBA0pVjZVZWVtK6deuU1vEXwIoVK8Tthw8fimM9duyY/PH09HSxTnY+9Ol48/LyJC8vL2nlypWSNqmOY5XZvXu3FBQUJAIzbQlYqutYFS1cuFDy9fWV9O1Y9en7qajQ0FDprbfekt/nY+RjjY6Olq/jH5m87u7du5I+Hevjx48lCwsL6cCBA5I2Q5NQFXFVnJOTEzVr1owWLVpE+fn5pW7LVeNcDac4oA9XvfFcSCdOnBD3T506JaofW7ZsKd+me/fuYpszZ86Qrhwra9++PW3atElUIfMElRs3bqTs7Gzq0qWLeJzLCgwMpHXr1lFmZqYob/ny5aIqskWLFqRp6j7eixcv0oMHD8R7yWV6eHhQnz596Nq1a6Rvx8q4qfOdd96hX375hSwtLUlbVMexFpWWlkaOjo6kb8eqT99Pii5cuCCaK7l5S4a/m7i8n3/+WTTvPX36VNxu2LChaJbRp2Pdv3+/eL/5+4mPj+fm43QFnkBRq2g6YtJlXIV2+PBhUXXGVfz29vbSv//971K3T0pKkmxtbaWJEydKmZmZohrv/fffFxH72LFjxTZz586VGjRoUOy5Li4u0o8//ijpyrGylJQUqWfPnuL4jI2NxbFzLZOimJgYqUWLFuIXmpGRkeTh4VFis4I+HO+GDRvEY3Xr1hVNQ+fPnxdV6U5OTlJycrKkT8daWFgomg64SYFFRUVpRQ1LdX2OFfGvb97mp59+kvTtWPXp+0nRe++9JzVs2LDY+rCwMMnPz0/UjPMSGBgo3bt3T9Kkr6vhWOfNmyeZmJiI49u7d6906tQpqVu3buJ+Tk6OpC0QsBQxZcoU8cda1nLz5s0ST+bPP/8s/sizs7NLPeH8x895LLIL9IgRI6TmzZtL7777bo1/IVT3sXIw1rp1a1HNePnyZWnGjBmSnZ2dqFaVXdReeeUVkd9w4sQJ6cKFC+KPiZtN4uLi1Hqs2nC83DbMr7F8+XL5c7g8Z2dnadmyZXp1rN99953UoUMH0X5e3QGLpo9VUWxsrLjAjRkzRqoOmj5Wffp+ksnKyhLHuHjx4mLr+VyMGjVKOnv2rLiIcw5L48aNxWP6dKxz584Vr6EYnPIPbA7SOIDRFghYiuA3iT8YZS2lRZyclMVvOidolYfzN/jXDOMkLm7zln34OGJWxLkPHNxs27ZN0pVj5aQtfpy3U8RR+7hx48Rt/lIsmtPDOImPI3510/TxHjp0SGzDiYqK+Evxs88+k/TpWDmZnN9b/tzKFn4O/88XAH06VhnO4eAk8pEjRyol0evTserL95Miztnh2gV+PUWca8bJ04rvJb8edwrg2lJ9OtZVq1aJMrjGWxEfv6ZrChUZa7pJSttwd2NeKoPbBbktl3MwysPdeNmhQ4coKSmJXnnlFXG/Xbt2lJqaKtoZZXkcvA23L3JXWF051qysLPE/b6PIyMhIHEtZ2/B92Tb6dLz8fnL+0u3bt0X3bZaXlye6Hfr4+JA+HeuSJUtozpw58sfi4uJEl1DOj9C3zzHjtv8XX3xRvMerV68utr2+HKu+fD8p4rwU/v4t+lp8PrgMzjuUkd1X9/eTpo+1Q4cO4n/+buL8FcZ5TI8ePVL7d1OVaDpi0lUnT54Umdpcbcrd/3799VdRLar465Grh7kNkLsIKkayXLXIv2Z++eUXkYH/0UcfKZXNbf/NmjUTz+OmEv7Vpslug5U5Vu72xzUlnTp1Euv4eLkakpvCdu3aJa9l4vyNQYMGibK5m+8nn3wifgHwfX07Xsb5S9zkxVWv/IuImw74Vwxn6evbsSrShhyW6jpWfg5vwzURfFux+6i+Hau+fD8p5hvx8XFPp6K4VsPMzEw0U3PvKK7J4CZ8blKpjiZrTR6rrFaUm7v+/vtvkbvz0ksvSY0aNdJoF+6iELBUEudbtGnTRnx4zc3NRRLTV199pdSOKPuS5gQpxbZKbgLiizL/oXMCFedyKOIETP4C4H7ynPT25ptvKvWX15Vj5XFVOBjhCzJXowYHBxfrMsldBznJjwM3GxsbqW3btqIrrCZV5/HyHz+Pz8Hb8PF27969WBW8vhyrtgUs1XWsq1evLjXvQB/fV335fmJTp06VvL29S23C++uvv0QuFpfN4yd17dpV/ODUx2NNS0sTXZ25yY+/j3m8MMUu3drAgP/RdC0PAAAAQFkwDgsAAABoPQQsAAAAoPUQsAAAAIDWQ8ACAAAAWg8BCwAAAGg9BCwAAACg9RCwAAAAgNZDwAJQw0aPHk0DBgzQyfM+Y8YMCg0NrdRzeQoCHtacl8qWoS1kx8LDoqubLn8+QPccO3aMXn75ZfL09BSf6R07dlTo+dnZ2eIz27RpUzI2Ni71s3vkyBFq3ry5mJ7E39+f1qxZU+F9xVxCAGqkOO9ISaZPn07fffcdD4Gqk+f9k08+oQ8++KBKZRw4cEDnA5bqpMufD9A9mZmZFBISQm+99RYNGjSows8vKCggCwsL+vDDD2nr1q0lbhMVFUX9+vWjd999l9avX08HDx6kt99+mzw8PMQcY6pCwAKgRvHx8fLbPNHftGnTxIRiMtbW1mLRVerYfycnJ7FAyezs7HBqoMb06dNHLKXJycmhzz//nDZs2CAmvmzSpAktWLCAunTpIh63srKipUuXitt///232KaoZcuWka+vL3399dfifsOGDenEiRP0zTffVChgQZMQgBq5u7vLF77wcI2L4jq+2Bet8uc/fK61mDRpEjk4OJCbmxutWLFC/PJ58803ycbGRlSh7tmzR+m1rl27Jr5ouEx+zsiRI8XsqqVJTk6mYcOGkZeXF1laWooqXP4Sknn48KHYx6+++kq+7uTJk2Rqaip+EZXUJMTVvK1btxZfWvb29mLW1/v371fonJVXxu+//y6qks3Nzal+/fo0c+ZMys/Plz/OX5Djxo0T54C34S/UP//8U/44/+pr3LixqIquV6+e/EtThtfxMfMvTD7XdevWpZ9++klpm7Nnz1KzZs1E+S1btqRLly4pPZ6SkkLDhw8Xs+Dyr82AgAAxa3NptmzZIs4/b8vBW/fu3cX7zUr6fPCv18mTJ5Ojo6N4j/h9UFTeOeCLQ6dOncTreXt7i/JkrwdQlvfff59OnTpFGzdupKtXr9Jrr71GvXv3prt375Kq+Pn8GVfEgQqvrwgELABaYO3ateTs7CwujBy8vPfee+KLoX379nTx4kXq2bOnCEh4ynvZBapr167iInr+/Hnau3cvJSYm0pAhQ8psa27RogXt2rVLBDtjx44VZfJrMr7Yrlq1SlwMucwnT56Ix/kLq1u3bsXK46CBL6ydO3cWX2T85cNlltcsVpEyjh8/TqNGjaKJEyfSjRs3aPny5aLte+7cueLxwsJCEbTxL7tff/1VbDN//nwyMjISj1+4cEGck6FDh1JYWJg4ti+//LJY+zkHMbJAZPz48eL8y2rGMjIy6KWXXqJGjRqJ8rgMbhpTxGXya3NQefPmTfGLk9/P0mrhOHDkAIm35YCNq+LLagbizwcHdGfOnKGFCxfSrFmzaP/+/Sqdg4iICHGBGTx4sDjHXPPHAQy/rwBliY6OFoH35s2bRcDr5+cnPvsdO3YsMyAvKiEhQQTTivh+eno6PX36VOVyMFszQDXhWXx5VtWi3njjDTGVu0znzp2ljh07yu/n5+dLVlZW0siRI+Xr4uPjxQysspliZ8+eLWa5VhQTEyO2uX37tsr72K9fPzF7tKLx48dLDRo0kP71r39JTZs2VZoJdvr06VJISIh81l5+vSNHjqj0WiXN2FxeGd26dROz0Sr65ZdfJA8PD3F73759kqGhYanHzMfQo0cPpXWffvqp1KhRI/l9Hx8facSIEfL7PHs6z1a8dOlScX/58uWSk5OT9PTpU/k2/Jjisbz88sti1mJVZ9zl5967d6/Ex8v7fLBWrVqJmd9VOQdjxoyRxo4dq7Tu+PHj4jmKxwRARNL27dvlJ+LPP/8U6/j7SHExNjaWhgwZUu5nVyYgIKDY3/GuXbtE2VlZWSqfeOSwAGiB4OBg+W3+ZczNBNxkICP7dZKUlCT+v3LlCh0+fLjEfBL+Rd2gQYMSk+O46eN///sfPXjwgHJzc0X7NDcPKVq8eLFoUuBfVVyjwE0pJeHmCW6+4KrdHj16iCpfrs3gRDpVlVcGHyfXHMhqVGTHwbVFXNvEvXTq1KlT4vEyrsHo37+/0jpucvr2229FObJaCMXzL2vGk51rLoMf56YWmXbt2imVyTUyXIMhqw3jWiOuHSsJJzhyjRW/v3zcvP2rr74qmgNLo7h/jM+PbP/KOwd8DrlmhZMdZbg2h2tmOBmS8wkASsK1i/w3wt8Dsr8VmYrksvHfE9cAK+L7tra2oplSVWgSAtACJiYmSvf5oqm4TtZEwhcZ2RcJd0Xki5Xiwu3KL7zwQomvsWjRItEDZcqUKSLY4e35gsmBS9GAJy4uTrwWd98tC1cLczMOX5y5qYEvmqdPn67QsZdVBh8n56woHiM37fBxcgBRkS+7ip5/2blWBTfJcN7Nv//9b3HuOCAp2mwkw1/83JzDzUfczPT9999TYGCgCB4qs3/lnQM+h5zfongOOYjhc8hV/ACl4SZnDuw5OOY8OsWFgxBVcYAvy4OT4b+BooF/eVDDAqCDOAmVk0k5YZTHPlAF11RwbcOIESPEfb7g3blzR1w0ZTh44cdff/11cRHlroccILi6upb5pcbL1KlTxRfQb7/9Rm3btq3Q8ZRWBh8n55LwF2RpNQ+xsbHiOEqqYeDaAz7uoueBty36i7E0XMYvv/wianVktSwlBWWcA/TGG2+Ihdv7P/30U1FbVRIOOLimhxfuSebj40Pbt2+njz76iCqqvHPA55DzWko7h1C7ZWRkUHh4uPw+B84c1HLtJ3+eOJmc88g4z4v/Rjk5n4MP/txxV2XGny/+7nj8+LHIfZONTyRL0OfuzD/88INIHOfcrUOHDomaXs6nqwjUsADooAkTJogvB07ePHfunKgV2bdvn+hVxL+ISsI9V/hXDff84WYO/tVdtJqWuy+mpaXRkiVLRE0Mf2HxF0xJ+IuNAwyuHeHahb/++kv8aq9IE0N5ZfDFfN26daKW5fr162K/ubfCF198IR7nZF2uUeLmGD42Lo9rLjgJmX388cfiy3X27Nnigs7Jq/zFWVrtR0n+9a9/iQDjnXfeEV/Mu3fvLhaI8H5ybyb+4uf95B46pZ0HTpzlpjlObOakxm3btomLQGWbZso7B/w+8nvOSbayWjjeVyTdAuPPoewHA+OgmW/zZ1pWA8oBC/8t8Y8Ybu7k7xzuTSfTt29f8ZydO3eKJHLF8hh3aebghD+f3CTKwc/KlSsr1KVZQJoRgOaTbidOnKi0DSeCfvPNN2UmxN25c0caOHCgZG9vL1lYWEhBQUHSpEmTRNJoSTjBlV/X2tpaJJV+8cUX0qhRo+T7cvjwYZFMxwmZiomytra20o8//lgs6TYhIUEaMGCASIA1NTUV+zxt2jSpoKBA5aRbVcrYu3ev1L59e3GMvC+tW7eWfvrpJ6Xj4oRXTow1NzeXmjRpIpIFZbZs2SKSbE1MTKS6detKixYtKvdc8zHyscpwsjOv430MDQ2Vtm7dqnQsnATdsGFDsY+Ojo7inEZGRpZ4Hm7cuCH16tVLcnFxkczMzESC8/fff1+hzwc/ztupeg7Onj0rko/5veekyeDgYGnu3Lkl7h+AtjLgfyoW4gAAVBznw/AvLe46jJFuAaCi0CQEADWKk2tL60EDAFAaJN0CQI3grrey0TFL6yoNAFAaNAkBAACA1kOTEAAAAGg9BCwAAACg9RCwAADoMB73gseJ4UVxlmcAfYOABUBDjh07JobX9/T0FBebHTt2FNuGB3bjuXZ4G57zp6Rp3XkEVh5Ijucf4vk9eACxogPC8QBlPColl8Gj1vIorDxTcnl4PqGgoCAxwivPfcODpiniQc94Lhx+bT4G2QiX5eFB73gETZ5LxN7ensaMGSNG3FQ8Jj5ufk0eyVfVC3F55TKeV4dHouVj8vb2FrMfl6e6zjGPKsEDdPHcQDzEPs+lVPT9Le+YuMcVzwBd1kzdAPoAAQuAhmRmZopRH//73/+WejHjC3VkZKQYmZTHL+Eh3Pmixs+V4flreIRJDi6OHj0q5rIZNGiQ/HEe+ZYvpDx0No94yqO9rlmzRj6SZWl4Wx5Jly+Q/Nq8L7xcu3ZN6Rh4qvkFCxZU6Nj5AswjwvLIlzwqLAdvY8eOVdpnvoB/+OGH4njVVS5PZ88BFp9HntCN51eaMWMG/fTTT2WWW13nmIMlHlV42bJlYgRcKysrMfonB0iqHpOpqamY10Vd8yoBaC1Nj1wHAMVHsWW3b98W669duyZfxyPA8gipK1asEPdTU1PFCK6bN2+Wb3Pz5k3xPB6dle3evVsyNDQUo8rKLF26VIwam5OTU+rp5+nj+/Xrp7SuTZs20rhx41QaxbY0PNIrb3vu3Dn5uj179kgGBgbSgwcPVJ6yvjLl8oi9Dg4OSsc9ZcoUKTAwsNRyq+sc84jE7u7uSiPv8mvx6LcbNmyo8LlS9TwB6CrUsABoqZycHPG/bMI9ZmhoKMYwOXHihLjPtQR5eXlKtRDchMPzfPD8PIz/56YVNzc3+Tb8K55rG/iXe2n4eUVrN/h5snIri5/PTRstW7aUr+PX4WPjWobqLJe34Xl3uFZC8Zh4gsWUlBSlnBDZTNXqOsdcHpfL5TOe8ychIUGpXDs7O2rTpo1SudVxrgB0EQIWAC0luyjy5IB8MeXmBm564Zl5OWeB8QWPL758UVPEF05+TLaN4oVU9rjssdKU9ryynqMKfn7R2Z85T4Vnh61K2aqUq8q54BwUnuTNxMREreeYy+NyuXzF9WWd4+o6VwC6CAELgJbiCxwntfIsw3yB4gvd4cOHqU+fPuIXtrpwsignksoWnklYXXhaecWydUHr1q3p1q1b5OXlpdZyuTwul8sHgIrD0PwAWqxFixai501aWpqoYXFxcRFNBrImAk625PWpqalKNQDcg4Ufk21z9uxZpXJlPVz4Me6BpNi7h4Mj2WNFe8IolquKWbNm0SeffKK0jp+flJSktI5703BvmIqUXZQq5ZZ2TLLHSiu3que4tHJl23EvIcXnySaHrK5zBaCLUMMCoAM4t4GDFe7yev78eerfv788oOGamIMHD8q35XwMrjVp166duM//h4WFKV34uMcJd5Nt1KiRaGLw9/eXL7KAhZ+nWK7sebJyVcHNGYply8rliz/nhsgcOnSICgsLRTBWWaqUy9twLxvOSVE8Jm6qcXBwKLFcdZzjkvDM1Rx0KJbLOS+cm6JYbnWcKwCdpOmsX4Da6smTJ6JXDS/8p/if//xH3L5//758m//973/S4cOHpYiICGnHjh2Sj4+PNGjQIKVy3n33Xalu3brSoUOHpPPnz0vt2rUTi0x+fr7UpEkTqWfPntLly5elvXv3ip5GU6dOLXP//v77b8nY2FhavHix6BUzffp00VsmLCxMvk1ycrLY5127dolj2Lhxo7gfHx9fZtm9e/eWmjVrJp05c0Y6ceKEFBAQIA0bNkxpm+vXr4uyXn75ZalLly7yc1WVcrkXjpubmzRy5EjR+4r319LSUlq+fLl8G34u9xqKjY1V6znm8rhcLl9m/vz5kr29vfT7779LV69eFb18fH19padPn1boXDH0EgJ9h4AFQEM4EOGLfNGFLzwy3333nVSnTh0RKPAF84svvijWTZYvbuPHjxfddfniO3DgwGIBw71796Q+ffpIFhYWkrOzs/Txxx9LeXl55e4jB0wNGjSQTE1NpcaNG4vARNHq1atLPAYObsrCgQ5fdK2trUXX3zfffFMEcIo4OCup7KqWe+XKFaljx46i+7CXl5cIGkp6X7irtjrPsazrN5ev2LX5yy+/FEEU70+3bt1Ed/aKHhNDwAL6DrM1AwDoAR4ZmJuPShoxGUAfIIcFAECHHT9+XPTAWr9+vaZ3BaBaoYYFAECHPX36lB48eCBuc+CC3kOgrxCwAAAAgNZDkxAAAABoPQQsAAAAoPUQsAAAAIDWQ8ACAAAAWg8BCwAAAGg9BCwAAACg9RCwAAAAgNZDwAIAAACk7f4fHvXgHTjeVTEAAAAASUVORK5CYII=",
+ "text/plain": [
+ ""
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ }
+ ],
+ "source": [
+ "check_density(DINO_1y_grid_T.rhop)\n",
+ "check_density(get_density(DINO_1y_grid_T.toce,DINO_1y_grid_T.soce,get_deptht(DINO_restart_file,DINO_mask_file),DINO_mask_file.tmask)[0])\n",
+ "\n",
+ "temperature_500m_30NS_metric(DINO_1y_grid_T.toce,DINO_mask_file).plot()\n",
+ "temperature_BWbox_metric(DINO_1y_grid_T.toce,DINO_mask_file)\n",
+ "temperature_DWbox_metric(DINO_1y_grid_U.uoce,DINO_mask_file)\n",
+ "ACC_Drake_metric(DINO_1y_grid_U.uoce,DINO_mask_file)\n",
+ "NASTG_BSF_max(DINO_1y_grid_V.voce,DINO_1m_To_1y_grid_T.ssh,DINO_mask_file)"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 7,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "array(0.00724717)"
+ ]
+ },
+ "execution_count": 7,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "check_density(get_density(DINO_1y_grid_T.toce,DINO_1y_grid_T.soce,get_deptht(DINO_restart_file,DINO_mask_file),DINO_mask_file.tmask)[0])\n"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 8,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# DINO_restart_file"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 9,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Successfully loaded dataset from so_pred\n",
+ "Successfully loaded dataset from thatao_pred\n",
+ "Successfully loaded dataset from zos_pred\n"
+ ]
+ }
+ ],
+ "source": [
+ "# Read prediction file\n",
+ "so_pred_path = \"../nc_files/nc_files/pred_so.npy\"\n",
+ "so_pred = np.load(so_pred_path)\n",
+ "print(f\"Successfully loaded dataset from so_pred\")\n",
+ "\n",
+ "thetao_pred_path = \"../nc_files/nc_files/pred_thetao.npy\"\n",
+ "thetao_pred = np.load(thetao_pred_path)\n",
+ "print(f\"Successfully loaded dataset from thatao_pred\")\n",
+ "\n",
+ "zos_pred_path = \"../nc_files/nc_files/pred_zos.npy\"\n",
+ "zos_pred = np.load(zos_pred_path)\n",
+ "print(f\"Successfully loaded dataset from zos_pred\")\n"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 29,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "soce_array = so_pred\n",
+ "soce_array = soce_array[-1,:,:,:][np.newaxis]\n",
+ "\n",
+ "toce_array = thetao_pred\n",
+ "toce_array = toce_array[-1,:,:,:][np.newaxis] \n",
+ "\n",
+ "\n",
+ "DINO_restart_file.sn[:]=soce_array #salinity(new)\n",
+ "DINO_restart_file.tn[:]=toce_array\n",
+ "\n",
+ "DINO_restart_file.sb[:]=soce_array #salinity(begining)\n",
+ "DINO_restart_file.tb[:]=toce_array"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 30,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "DINO_restart_file.to_netcdf(\"../nc_files/nc_files/restart_updated.nc\")"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 25,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "soce_array shape: (1, 36, 199, 62)\n",
+ "DINO_restart_file.sn shape: (1, 36, 199, 62)\n"
+ ]
+ }
+ ],
+ "source": [
+ "print(\"soce_array shape:\", soce_array.shape)\n",
+ "print(\"DINO_restart_file.sn shape:\", DINO_restart_file.sn.shape)"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 31,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "restart_updated = xr.open_dataset(\"../nc_files/nc_files/restart_updated.nc\")"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 34,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/html": [
+ "\n",
+ "
<xarray.Dataset> Size: 51MB\n",
+ "Dimensions: (nav_lat: 199, nav_lon: 62, depth: 36, time_counter: 1)\n",
+ "Coordinates:\n",
+ " nav_lon (nav_lat, nav_lon) float32 49kB ...\n",
+ " nav_lat (nav_lat, nav_lon) float32 49kB ...\n",
+ " * depth (depth) float32 144B 5.034 15.32 25.96 ... 3.757e+03 4.253e+03\n",
+ " * time_counter (time_counter) float32 4B 5.76e+05\n",
+ "Data variables: (12/29)\n",
+ " kt float64 8B ...\n",
+ " ndastp float64 8B ...\n",
+ " adatrj float64 8B ...\n",
+ " ntime float64 8B ...\n",
+ " utau_b (time_counter, nav_lat, nav_lon) float64 99kB ...\n",
+ " vtau_b (time_counter, nav_lat, nav_lon) float64 99kB ...\n",
+ " ... ...\n",
+ " sshn (time_counter, nav_lat, nav_lon) float64 99kB ...\n",
+ " un (time_counter, depth, nav_lat, nav_lon) float64 4MB ...\n",
+ " vn (time_counter, depth, nav_lat, nav_lon) float64 4MB ...\n",
+ " tn (time_counter, depth, nav_lat, nav_lon) float64 4MB ...\n",
+ " sn (time_counter, depth, nav_lat, nav_lon) float64 4MB ...\n",
+ " rhop (time_counter, depth, nav_lat, nav_lon) float64 4MB 0.0 ......\n",
+ "Attributes:\n",
+ " file_name: DINO_00576000_restart.nc\n",
+ " TimeStamp: 13/11/2024 17:53:42 -0000
- nav_lat: 199
- nav_lon: 62
- depth: 36
- time_counter: 1
kt
()
float64
...
[1 values with dtype=float64]
ndastp
()
float64
...
[1 values with dtype=float64]
adatrj
()
float64
...
[1 values with dtype=float64]
ntime
()
float64
...
[1 values with dtype=float64]
utau_b
(time_counter, nav_lat, nav_lon)
float64
...
[12338 values with dtype=float64]
vtau_b
(time_counter, nav_lat, nav_lon)
float64
...
[12338 values with dtype=float64]
qns_b
(time_counter, nav_lat, nav_lon)
float64
...
[12338 values with dtype=float64]
emp_b
(time_counter, nav_lat, nav_lon)
float64
...
[12338 values with dtype=float64]
sfx_b
(time_counter, nav_lat, nav_lon)
float64
...
[12338 values with dtype=float64]
en
(time_counter, depth, nav_lat, nav_lon)
float64
...
[444168 values with dtype=float64]
avt_k
(time_counter, depth, nav_lat, nav_lon)
float64
...
[444168 values with dtype=float64]
avm_k
(time_counter, depth, nav_lat, nav_lon)
float64
...
[444168 values with dtype=float64]
dissl
(time_counter, depth, nav_lat, nav_lon)
float64
...
[444168 values with dtype=float64]
sbc_hc_b
(time_counter, nav_lat, nav_lon)
float64
...
[12338 values with dtype=float64]
sbc_sc_b
(time_counter, nav_lat, nav_lon)
float64
...
[12338 values with dtype=float64]
qsr_hc_b
(time_counter, depth, nav_lat, nav_lon)
float64
...
[444168 values with dtype=float64]
fraqsr_1lev
(time_counter, nav_lat, nav_lon)
float64
...
[12338 values with dtype=float64]
rdt
()
float64
...
[1 values with dtype=float64]
sshb
(time_counter, nav_lat, nav_lon)
float64
...
[12338 values with dtype=float64]
ub
(time_counter, depth, nav_lat, nav_lon)
float64
...
[444168 values with dtype=float64]
vb
(time_counter, depth, nav_lat, nav_lon)
float64
...
[444168 values with dtype=float64]
tb
(time_counter, depth, nav_lat, nav_lon)
float64
...
[444168 values with dtype=float64]
sb
(time_counter, depth, nav_lat, nav_lon)
float64
...
[444168 values with dtype=float64]
sshn
(time_counter, nav_lat, nav_lon)
float64
...
[12338 values with dtype=float64]
un
(time_counter, depth, nav_lat, nav_lon)
float64
...
[444168 values with dtype=float64]
vn
(time_counter, depth, nav_lat, nav_lon)
float64
...
[444168 values with dtype=float64]
tn
(time_counter, depth, nav_lat, nav_lon)
float64
...
[444168 values with dtype=float64]
sn
(time_counter, depth, nav_lat, nav_lon)
float64
...
[444168 values with dtype=float64]
rhop
(time_counter, depth, nav_lat, nav_lon)
float64
0.0 0.0 0.0 0.0 ... 0.0 0.0 0.0 0.0
array([[[[0., ..., 0.],\n",
+ " ...,\n",
+ " [0., ..., 0.]],\n",
+ "\n",
+ " ...,\n",
+ "\n",
+ " [[0., ..., 0.],\n",
+ " ...,\n",
+ " [0., ..., 0.]]]])
PandasIndex
PandasIndex(Index([ 5.033581733703613, 15.322633743286133, 25.958648681640625,\n",
+ " 37.014060974121094, 48.57634735107422, 60.75108337402344,\n",
+ " 73.66564178466797, 87.47360229492188, 102.36001586914062,\n",
+ " 118.54764556884766, 136.30438232421875, 155.9520263671875,\n",
+ " 177.87646484375, 202.53965759277344, 230.49343872070312,\n",
+ " 262.3951721191406, 299.0254821777344, 341.30792236328125,\n",
+ " 390.3298034667969, 447.3644104003906, 513.892822265625,\n",
+ " 591.623779296875, 682.51025390625, 788.75830078125,\n",
+ " 912.8260498046875, 1057.4073486328125, 1225.396240234375,\n",
+ " 1419.828369140625, 1643.7969970703125, 1900.3427734375,\n",
+ " 2192.322265625, 2522.26416015625, 2892.223876953125,\n",
+ " 3303.654541015625, 3757.30908203125, 4253.1875],\n",
+ " dtype='float32', name='depth'))
PandasIndex
PandasIndex(Index([576000.0], dtype='float32', name='time_counter'))
- file_name :
- DINO_00576000_restart.nc
- TimeStamp :
- 13/11/2024 17:53:42 -0000
"
+ ],
+ "text/plain": [
+ " Size: 51MB\n",
+ "Dimensions: (nav_lat: 199, nav_lon: 62, depth: 36, time_counter: 1)\n",
+ "Coordinates:\n",
+ " nav_lon (nav_lat, nav_lon) float32 49kB ...\n",
+ " nav_lat (nav_lat, nav_lon) float32 49kB ...\n",
+ " * depth (depth) float32 144B 5.034 15.32 25.96 ... 3.757e+03 4.253e+03\n",
+ " * time_counter (time_counter) float32 4B 5.76e+05\n",
+ "Data variables: (12/29)\n",
+ " kt float64 8B ...\n",
+ " ndastp float64 8B ...\n",
+ " adatrj float64 8B ...\n",
+ " ntime float64 8B ...\n",
+ " utau_b (time_counter, nav_lat, nav_lon) float64 99kB ...\n",
+ " vtau_b (time_counter, nav_lat, nav_lon) float64 99kB ...\n",
+ " ... ...\n",
+ " sshn (time_counter, nav_lat, nav_lon) float64 99kB ...\n",
+ " un (time_counter, depth, nav_lat, nav_lon) float64 4MB ...\n",
+ " vn (time_counter, depth, nav_lat, nav_lon) float64 4MB ...\n",
+ " tn (time_counter, depth, nav_lat, nav_lon) float64 4MB ...\n",
+ " sn (time_counter, depth, nav_lat, nav_lon) float64 4MB ...\n",
+ " rhop (time_counter, depth, nav_lat, nav_lon) float64 4MB 0.0 ......\n",
+ "Attributes:\n",
+ " file_name: DINO_00576000_restart.nc\n",
+ " TimeStamp: 13/11/2024 17:53:42 -0000"
+ ]
+ },
+ "execution_count": 34,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "restart_updated"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/html": [
+ "\n",
+ "
<xarray.DataArray (time_counter: 1)> Size: 8B\n",
+ "array([35.51578141])\n",
+ "Coordinates:\n",
+ " * time_counter (time_counter) float32 4B 5.76e+05\n",
+ " depth float32 4B 5.034
"
+ ],
+ "text/plain": [
+ " Size: 8B\n",
+ "array([35.51578141])\n",
+ "Coordinates:\n",
+ " * time_counter (time_counter) float32 4B 5.76e+05\n",
+ " depth float32 4B 5.034"
+ ]
+ },
+ "execution_count": 53,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "check_density(restart_updated.rhop)\n",
+ "check_density(get_density(restart_updated.tn,restart_updated.sn,get_deptht(restart_updated,DINO_mask_file),DINO_mask_file.tmask)[0])\n",
+ "temperature_500m_30NS_metric(restart_updated.tn,DINO_mask_file)\n",
+ "temperature_BWbox_metric(restart_updated.tn,DINO_mask_file)\n",
+ "temperature_DWbox_metric(restart_updated.un,DINO_mask_file)\n",
+ "ACC_Drake_metric(restart_updated.un,DINO_mask_file)\n",
+ "NASTG_BSF_max(restart_updated.vn,restart_updated.sshn,DINO_mask_file)\n",
+ "\n"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 52,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ " Size: 4B\n",
+ "array([576000.], dtype=float32)\n",
+ "Coordinates:\n",
+ " * time_counter (time_counter) float32 4B 5.76e+05\n",
+ "Only one time step available in the restart file.\n"
+ ]
+ }
+ ],
+ "source": [
+ "# Check time variable (usually named \"time_counter\" or similar)\n",
+ "print(restart_updated[\"time_counter\"]) # Replace with actual time variable if different\n",
+ "\n",
+ "# Get time resolution\n",
+ "time_values = restart_updated[\"time_counter\"].values\n",
+ "if len(time_values) > 1:\n",
+ " time_resolution = time_values[1] - time_values[0]\n",
+ " print(f\"Time resolution: {time_resolution} (time units)\")\n",
+ "else:\n",
+ " print(\"Only one time step available in the restart file.\")\n",
+ "\n"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 41,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ " Size: 4MB\n",
+ "[444168 values with dtype=float64]\n",
+ "Coordinates:\n",
+ " nav_lon (nav_lat, nav_lon) float32 49kB ...\n",
+ " nav_lat (nav_lat, nav_lon) float32 49kB -69.85 -69.85 ... 69.85 69.85\n",
+ " * depth (depth) float32 144B 5.034 15.32 25.96 ... 3.757e+03 4.253e+03\n",
+ " * time_counter (time_counter) float32 4B 5.76e+05\n"
+ ]
+ }
+ ],
+ "source": [
+ "print(restart_updated[\"sn\"].isel(time_counter=slice(5))) # Replace \"your_variable\" with the actual variable name\n",
+ "\n"
+ ]
+ }
+ ],
+ "metadata": {
+ "kernelspec": {
+ "display_name": "myenv",
+ "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.9"
+ }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 2
+}
diff --git a/read_metric.py b/read_metric.py
new file mode 100644
index 0000000..8a738e7
--- /dev/null
+++ b/read_metric.py
@@ -0,0 +1,66 @@
+### Clone metrics.py file on CSD3 and read it here
+### Read density.py file
+
+import importlib
+import xarray as xr
+import sys
+from pathlib import Path
+
+# Define paths to the .nc file and metric file
+nc_file_path = "../nemo_4.2.1/tests/DINO/EXP00/DINO_1y_grid_T.nc" # Replace with your .nc file
+metric_file_path = "../Metrics-Ocean/metrics.py" # Path to your metrics file
+density_file_path = 'lib/density.py'
+
+# Import the metrics module dynamically
+metric_module = Path(metric_file_path).stem # Extract module name from file
+spec = importlib.util.spec_from_file_location(metric_module, metric_file_path)
+metrics = importlib.util.module_from_spec(spec)
+sys.modules[metric_module] = metrics
+spec.loader.exec_module(metrics)
+
+# Open the .nc file using xarray
+try:
+ dataset = xr.open_dataset(nc_file_path)
+ print(f"Successfully loaded dataset from {nc_file_path}")
+except FileNotFoundError:
+ print(f"Error: File {nc_file_path} not found.")
+ sys.exit(1)
+
+
+# List all variables in the dataset
+print("\nAvailable variables in the dataset:")
+print(dataset.data_vars)
+
+# Define variable mapping (dataset name -> metric file name)
+variable_mapping = {
+ "rhop": "density",
+ "deptht": "depth",
+}
+
+# Define variables and metrics to apply
+variables_to_check = ["rhop"] # Replace with your variable names
+metrics_to_apply = ["check_density"] # Replace with your metric functions
+
+results = {}
+for var in variables_to_check:
+ if var in dataset:
+ print(f"\nApplying metrics to variable: {var}")
+ variable_data = dataset[var]
+ metric_var_name = variable_mapping.get(var, var)
+ results[var] = {}
+
+ # Apply metrics from metrics.py
+ for metric in metrics_to_apply:
+ if hasattr(metrics, metric):
+ metric_function = getattr(metrics, metric)
+ try:
+ result = metric_function(variable_data, metric_var_name)
+ print(result)
+ results[var][metric] = result
+ print(f"Metric '{metric}' applied to '{metric_var_name}'. Result: {result}")
+ except Exception as e:
+ print(f"Error applying metric '{metric}' to '{metric_var_name}': {e}")
+ else:
+ print(f"Metric '{metric}' not found in {metric_file_path}")
+
+
diff --git a/requirements.txt b/requirements.txt
index 8aea04a..009271b 100644
--- a/requirements.txt
+++ b/requirements.txt
@@ -1,9 +1,131 @@
+asttokens==3.0.0
+attrs==25.1.0
+bokeh==3.6.3
+Bottleneck==1.4.2
+certifi==2025.1.31
+cfgv==3.4.0
+cftime==1.6.4.post1
+charset-normalizer==3.4.1
+click==8.1.8
+cloudpickle==3.1.1
+comm==0.2.2
+contourpy==1.3.1
+coverage==7.6.10
+crc32c==2.7.1
+cycler==0.12.1
+dask==2025.1.0
+debugpy==1.8.12
+decorator==5.1.1
+Deprecated==1.2.18
+distlib==0.3.9
+distributed==2025.1.0
+donfig==0.8.1.post1
+einops==0.8.0
+execnet==2.1.1
+executing==2.2.0
+filelock==3.17.0
+flox==0.10.0
+fonttools==4.55.8
+fsspec==2025.2.0
+h5netcdf==1.5.0
+h5py==3.12.1
+hypothesis==6.125.2
+identify==2.6.6
+idna==3.10
+importlib_metadata==8.6.1
+iniconfig==2.0.0
+ipykernel==6.29.5
+ipython==8.32.0
+jedi==0.19.2
+Jinja2==3.1.5
+joblib==1.4.2
+jupyter_client==8.6.3
+jupyter_core==5.7.2
+kiwisolver==1.4.8
+llvmlite==0.44.0
+locket==1.0.0
+lz4==4.4.3
+MarkupSafe==3.0.2
+matplotlib==3.9.2
+matplotlib-inline==0.1.7
+mpmath==1.3.0
+msgpack==1.1.0
+mypy==1.15.0
+mypy-extensions==1.0.0
+nc-time-axis==1.4.1
+nest-asyncio==1.6.0
+netCDF4==1.7.2
+networkx==3.4.2
+nodeenv==1.9.1
+numba==0.61.0
+numbagg==0.8.2
+numcodecs==0.15.0
numpy==2.0.2
+numpy-groupies==0.11.2
+nvidia-cublas-cu12==12.4.5.8
+nvidia-cuda-cupti-cu12==12.4.127
+nvidia-cuda-nvrtc-cu12==12.4.127
+nvidia-cuda-runtime-cu12==12.4.127
+nvidia-cudnn-cu12==9.1.0.70
+nvidia-cufft-cu12==11.2.1.3
+nvidia-curand-cu12==10.3.5.147
+nvidia-cusolver-cu12==11.6.1.9
+nvidia-cusparse-cu12==12.3.1.170
+nvidia-cusparselt-cu12==0.6.2
+nvidia-nccl-cu12==2.21.5
+nvidia-nvjitlink-cu12==12.4.127
+nvidia-nvtx-cu12==12.4.127
+opt_einsum==3.4.0
+packaging==24.2
pandas==2.2.2
+parso==0.8.4
+partd==1.4.2
+pexpect==4.9.0
+pillow==11.1.0
+platformdirs==4.3.6
+pluggy==1.5.0
+pooch==1.8.2
+pre_commit==4.1.0
+prompt_toolkit==3.0.50
+psutil==6.1.1
+ptyprocess==0.7.0
+pure_eval==0.2.3
+pyarrow==19.0.0
+Pygments==2.19.1
+pyparsing==3.2.1
+pytest==8.3.4
+pytest-cov==6.0.0
+pytest-env==1.1.5
+pytest-timeout==2.3.1
+pytest-xdist==3.6.1
python-dateutil==2.9.0.post0
pytz==2024.1
-six==1.16.0
-matplotlib==3.9.2
-scipy==1.13.1
-xarray[complete]==2024.7.0
+PyYAML==6.0.2
+pyzmq==26.2.1
+requests==2.32.3
+ruff==0.9.4
scikit-learn==1.5.1
+scipy==1.13.1
+seaborn==0.13.2
+six==1.16.0
+sortedcontainers==2.4.0
+stack-data==0.6.3
+sympy==1.13.1
+tblib==3.0.0
+threadpoolctl==3.5.0
+toolz==1.0.0
+torch==2.6.0
+tornado==6.4.2
+traitlets==5.14.3
+triton==3.2.0
+typing_extensions==4.12.2
+tzdata==2025.1
+urllib3==2.3.0
+virtualenv==20.29.1
+wcwidth==0.2.13
+wrapt==1.17.2
+xarray==2024.7.0
+xyzservices==2025.1.0
+zarr==3.0.2
+zict==3.0.0
+zipp==3.21.0