{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "5bbd0865",
   "metadata": {},
   "source": [
    "# Accessing USGS data\n",
    "\n",
    "### This notebook does the following:  \n",
    "* Load a single USGS point\n",
    "* Look at the metadata\n",
    "* Get multiple radiation datasets\n",
    "* Calculate albedo from the datasets\n",
    "* Plot the albedo evolution"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f8319b6d",
   "metadata": {},
   "outputs": [],
   "source": [
    "# imports\n",
    "from datetime import datetime\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "\n",
    "from metloom.pointdata import USGSPointData\n",
    "\n",
    "# For rendering in readthedocs\n",
    "import plotly.offline as py\n",
    "py.init_notebook_mode(connected=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a9911803",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define a point known to have solar measurements\n",
    "pt = USGSPointData(\"395709105582701\", \"Blue Ridge Meteorological Station NR Fraser\")\n",
    "\n",
    "# start data and end date\n",
    "start_date = datetime(2024, 1, 1)\n",
    "end_date = datetime(2024, 7, 1)\n",
    "# Define a list of variables we want\n",
    "incoming_sw = pt.ALLOWED_VARIABLES.DOWNSHORTWAVE\n",
    "outgoing_sw = pt.ALLOWED_VARIABLES.UPSHORTWAVE\n",
    "depth = pt.ALLOWED_VARIABLES.SNOWDEPTH\n",
    "variables = [incoming_sw, outgoing_sw, depth]\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8dc795ff-a2d9-4e93-ae1c-e55259b14204",
   "metadata": {},
   "outputs": [],
   "source": [
    "# LETS GET THAT DATA\n",
    "df = pt.get_hourly_data(start_date, end_date, variables)\n",
    "df.head(10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "68ccbe7e",
   "metadata": {
    "editable": true,
    "slideshow": {
     "slide_type": ""
    }
   },
   "outputs": [],
   "source": [
    "# Check out this data!\n",
    "df = df.reset_index().set_index(\"datetime\")\n",
    "var_names = [v.name for v in variables]\n",
    "\n",
    "# Sample to just the vars\n",
    "df_rad = df.loc[:, var_names]\n",
    "\n",
    "df_rad.loc[:, [incoming_sw.name, outgoing_sw.name]].plot()\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9112e047-7f6f-4771-be01-c6f1efdbd661",
   "metadata": {
    "tags": [
     "nbsphinx-thumbnail",
     "nbsphinx-gallery"
    ]
   },
   "outputs": [],
   "source": [
    "# Let's get the data in a more usable state\n",
    "sw_thresh = 10\n",
    "\n",
    "# mask the SW to decent values\n",
    "df_rad[incoming_sw.name] = df_rad[incoming_sw.name].mask(df_rad[incoming_sw.name] < sw_thresh, np.nan)\n",
    "df_rad[outgoing_sw.name] = df_rad[outgoing_sw.name].mask(df_rad[outgoing_sw.name] < sw_thresh, np.nan)\n",
    "\n",
    "# Resample to daily, based on mean values\n",
    "df_rad = df_rad.resample(\"D\").mean()\n",
    "\n",
    "# Plot again\n",
    "df_rad.loc[:, [incoming_sw.name, outgoing_sw.name]].plot()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5515ac70-b7f7-4eaa-b58a-b7c4201e2717",
   "metadata": {},
   "source": [
    "### Now we can think Albedo\n",
    "\n",
    "The daily data plot looks much better. Now we can start to think about albedo.\n",
    "We can use the SW to calcaulate albedo next"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f9c91bbe",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Calculate albedo\n",
    "albedo_var = \"ALBEDO\"\n",
    "df_rad[albedo_var] = df_rad[outgoing_sw.name] / df_rad[incoming_sw.name]\n",
    "df_rad[albedo_var] = df_rad[albedo_var].mask(df_rad[albedo_var] < 0.1, np.nan)\n",
    "df_rad[albedo_var] = df_rad[albedo_var].mask(df_rad[albedo_var] > 1, np.nan)\n",
    "df_rad[albedo_var].plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "61cb8617-dd47-495c-9164-5a69efaded78",
   "metadata": {},
   "source": [
    "### Better plotting\n",
    "Awesome - now we have albedo. We haven't used snowdepth yet,\n",
    "so let's create a better plot that shows how the two relate"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d0fe841b-0b9a-43fa-9df2-e106ecf2483c",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Get plotly for a nicer plot\n",
    "# !pip install plotly\n",
    "import plotly.express as px\n",
    "import plotly.graph_objects as go\n",
    "from plotly.subplots import make_subplots"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "445efe95-75a7-46db-88d0-4bc20027f789",
   "metadata": {
    "editable": true,
    "slideshow": {
     "slide_type": ""
    },
    "tags": []
   },
   "outputs": [],
   "source": [
    "# Create figure with secondary y-axis\n",
    "fig = make_subplots(specs=[[{\"secondary_y\": True}]])\n",
    "\n",
    "# Add traces\n",
    "fig.add_trace(\n",
    "    go.Scatter(x=df_rad.index, y=df_rad[albedo_var], name=albedo_var),\n",
    "    secondary_y=False,\n",
    ")\n",
    "\n",
    "fig.add_trace(\n",
    "    go.Scatter(\n",
    "        x=df_rad.index, y=df_rad[depth.name].diff(), name=f\"{depth.name} signal\", mode=\"lines+markers\"\n",
    "    ), secondary_y=True, \n",
    ")\n",
    "\n",
    "fig.add_trace(\n",
    "    go.Scatter(\n",
    "        x=df_rad.index, y=df_rad[depth.name], name=f\"{depth.name}\", opacity=0.3,\n",
    "    ), secondary_y=True, \n",
    ")\n",
    "\n",
    "fig.update_layout(\n",
    "    # template='plotly_dark',\n",
    "    title=f'{pt.name}',\n",
    "    xaxis=dict(title='Date'),\n",
    "    yaxis=dict(\n",
    "        title=f'Unitless',\n",
    "        titlefont=dict(color='blue'),\n",
    "        tickfont=dict(color='blue'),\n",
    "        tickvals=[.25, .5, .75, 1.0],\n",
    "        range=[.25, 1]\n",
    "    ),\n",
    "    yaxis2=dict(\n",
    "        title=f'[m]',\n",
    "        titlefont=dict(color='red'),\n",
    "        tickfont=dict(color='red'), overlaying='y', side='right'\n",
    "    )\n",
    ")\n",
    "\n",
    "# Show the plot\n",
    "fig.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "eac2af7c-b240-4e95-baa1-e65611347a95",
   "metadata": {
    "editable": true,
    "slideshow": {
     "slide_type": ""
    },
    "tags": []
   },
   "source": [
    "### Summary\n",
    "* We started with a point that we knew had shortwave measurements.\n",
    "* Next, we pulled all the necessary data, cleaned it, and resampled to daily.\n",
    "* We calculated albedo, and plotted it to see how the snow"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "72ac1ea5-3c9c-4b85-8d74-454055599194",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.9.18"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}