IFT_Step00_unimpaired_flows.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import ulmo\n",
    "#from ulmo import cdec\n",
    "from ulmo.util import convert_datetime\n",
    "\n",
    "import os\n",
    "import shutil\n",
    "import sys\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "import plotly.express as px\n",
    "import time\n",
    "from datetime import datetime"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "outer_dir = os.path.abspath(os.path.join(os.getcwd() ,\"..\"))\n",
    "unimpath = os.path.join(outer_dir, \"IFT_files\", 'Unimpaired Flow')\n",
    "# Check if the directory exists, if not, create it\n",
    "if not os.path.exists(unimpath):\n",
    "    os.makedirs(unimpath)\n",
    "ref_dir = os.path.abspath(os.path.join(os.getcwd() ,\"Reference Files\"))\n",
    "nav_div_csv = os.path.join(ref_dir, 'Indian','Navarro_estimated_diversions.csv')\n",
    "# Current date as text in YYYYMMDD format\n",
    "date = datetime.today().strftime('%Y%m%d')\n",
    "\n",
    "\n",
    "LOI_name = \"Indian\" #Sproul or Indian \n",
    "if LOI_name == \"Sproul\":\n",
    "    site_code = '11476600' # Bull Creek gage https://waterdata.usgs.gov/nwis/inventory/?site_no=11476600&agency_cd=USGS\n",
    "elif LOI_name == \"Indian\":\n",
    "    site_code = '11468000' # Navarro Gage https://waterdata.usgs.gov/nwis/inventory/?site_no=11468000&agency_cd=USGS\n",
    "else:\n",
    "    print(\"LOI not recognized\")\n",
    "\n",
    "file_out = os.path.join(unimpath, LOI_name + \"_\" + date + \".csv\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "def download_daily_gage_data(gage_id):\n",
    "    \"\"\"Function to download stream gage data from USGS API\"\"\"\n",
    "    result = ulmo.usgs.nwis.get_site_data(site_code=gage_id, \n",
    "                                          service='dv',\n",
    "                                          parameter_code='00060',\n",
    "                                          methods='all',\n",
    "                                          start = '1850-1-1',\n",
    "                                          end = '2021-9-30'\n",
    "                                         )\n",
    "    if result == {}:\n",
    "        result = ulmo.usgs.nwis.get_site_data(site_code=str(site_code), \n",
    "                                              service='dv',\n",
    "                                              parameter_code='72137',\n",
    "                                              methods='all',\n",
    "                                              start = '1850-1-1',\n",
    "                                              end = '2021-9-30'\n",
    "                                             )\n",
    "        if result == {}:\n",
    "            print(\"discharge parameters not available for this gage\")\n",
    "            return\n",
    "\n",
    "    first_key = list(result.keys())[0]\n",
    "    # Convert to dataframe\n",
    "    first_key = list(result.keys())[0]\n",
    "    df = pd.DataFrame.from_dict(result[first_key]['values'])\n",
    "    df[\"value\"] = df[\"value\"].replace(['---','BRT','ART','-1', 'nan'],np.NaN)\n",
    "    df[['value']] = df[['value']].apply(pd.to_numeric)\n",
    "    df['datetime'] =  pd.to_datetime(df['datetime'])\n",
    "    df['gage_id'] =  gage_id\n",
    "    return df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "def evaluate_record(df):\n",
    "    \"\"\"Evaluate the gage data completeness\"\"\"\n",
    "    first = df['datetime'].min()\n",
    "    last = df['datetime'].max()\n",
    "    Qmin = df['value'].min()\n",
    "    Qmax = df['value'].max()\n",
    "    delta = last - first\n",
    "    all_records = len(df)\n",
    "    df_valid = df.dropna(subset=['value'])\n",
    "    valid_records = len(df_valid)\n",
    "    pct_missing = all_records / valid_records\n",
    "    pct_valid = valid_records / delta.days\n",
    "    print('First date: {}'.format(first))\n",
    "    print('Last date: {}'.format(last))\n",
    "    print('Min Flow: {} cfs'.format(Qmin))\n",
    "    print('Max Flow: {} cfs'.format(Qmax))\n",
    "    print('Percent of records with missing data: {0:.1%}'.format(1 - pct_missing))\n",
    "    print('Percent of the timeframe with valid records: {0:.1%}'.format(pct_valid))\n",
    "    df_month = df_valid.resample('M', on='datetime').count()\n",
    "    df_month = df_month['value'].reset_index()\n",
    "    df_month.rename(columns={'value': 'valid_count'}, inplace=True)\n",
    "    #Remove months with zero days\n",
    "    df_month = df_month.loc[df_month.valid_count > 0]\n",
    "    df_month['days_in_month'] = df_month['datetime'].dt.daysinmonth\n",
    "    df_month['pct_valid_in_month'] = df_month['valid_count'] / df_month['days_in_month']\n",
    "    pct_partial_month = len(df_month.loc[df_month['pct_valid_in_month'] < 1]) / len(df_month)\n",
    "    print('Percent of partial months: {0:.1%}'.format(pct_partial_month))\n",
    "    \n",
    "    \n",
    "    #Evaluation the number in days in each month\n",
    "    df_day_freq = (df_month.groupby('valid_count')\n",
    "               .agg({'days_in_month':'count'})\n",
    "               .reset_index()\n",
    "               .rename(columns={'days_in_month':'Count_months'})\n",
    "              )\n",
    "    print('Number of unique day counts per month: {}'.format(len(df_day_freq)))\n",
    "    return df_day_freq"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "processing data from request: https://waterservices.usgs.gov/nwis/dv/?format=waterml&site=11468000&parameterCd=00060&startDT=1850-01-01&endDT=2021-09-30\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "First date: 1950-10-01 00:00:00\n",
      "Last date: 2021-09-30 00:00:00\n",
      "Min Flow: 0.0 cfs\n",
      "Max Flow: 45100.0 cfs\n",
      "Percent of records with missing data: 0.0%\n",
      "Percent of the timeframe with valid records: 100.0%\n",
      "Percent of partial months: 0.0%\n",
      "Number of unique day counts per month: 4\n"
     ]
    }
   ],
   "source": [
    "df = download_daily_gage_data(str(site_code))\n",
    "df_eval = evaluate_record(df)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "E:\\Users\\kklausmeyer\\AppData\\Local\\Temp\\192\\ipykernel_62648\\344424719.py:9: FutureWarning: The pandas.datetime class is deprecated and will be removed from pandas in a future version. Import from datetime module instead.\n",
      "  return(pd.datetime(row.datetime.year+1,1,1).year)\n",
      "E:\\Users\\kklausmeyer\\AppData\\Local\\Temp\\192\\ipykernel_62648\\344424719.py:11: FutureWarning: The pandas.datetime class is deprecated and will be removed from pandas in a future version. Import from datetime module instead.\n",
      "  return(pd.datetime(row.datetime.year,1,1).year)\n"
     ]
    }
   ],
   "source": [
    "# Add in month and day fields\n",
    "df['Year'] = df.datetime.dt.year\n",
    "df['Month'] = df.datetime.dt.month\n",
    "df['Day'] = df.datetime.dt.day\n",
    "\n",
    "#Add in water year\n",
    "def assign_wy(row):\n",
    "    if row.datetime.month>=10:\n",
    "        return(pd.datetime(row.datetime.year+1,1,1).year)\n",
    "    else:\n",
    "        return(pd.datetime(row.datetime.year,1,1).year)\n",
    "\n",
    "df['Water Year'] = df.apply(lambda x: assign_wy(x), axis=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>gage_id</th>\n",
       "      <th>Water Year</th>\n",
       "      <th>avg_ann_flow</th>\n",
       "      <th>days_with_data</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>11468000</td>\n",
       "      <td>1951</td>\n",
       "      <td>619.506575</td>\n",
       "      <td>365</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>11468000</td>\n",
       "      <td>1953</td>\n",
       "      <td>603.172603</td>\n",
       "      <td>365</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>11468000</td>\n",
       "      <td>1954</td>\n",
       "      <td>520.427397</td>\n",
       "      <td>365</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>11468000</td>\n",
       "      <td>1955</td>\n",
       "      <td>251.441096</td>\n",
       "      <td>365</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>11468000</td>\n",
       "      <td>1957</td>\n",
       "      <td>320.472877</td>\n",
       "      <td>365</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>53</th>\n",
       "      <td>11468000</td>\n",
       "      <td>2004</td>\n",
       "      <td>474.798661</td>\n",
       "      <td>366</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>57</th>\n",
       "      <td>11468000</td>\n",
       "      <td>2008</td>\n",
       "      <td>325.683169</td>\n",
       "      <td>366</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>61</th>\n",
       "      <td>11468000</td>\n",
       "      <td>2012</td>\n",
       "      <td>220.550437</td>\n",
       "      <td>366</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>65</th>\n",
       "      <td>11468000</td>\n",
       "      <td>2016</td>\n",
       "      <td>507.661995</td>\n",
       "      <td>366</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>69</th>\n",
       "      <td>11468000</td>\n",
       "      <td>2020</td>\n",
       "      <td>76.058142</td>\n",
       "      <td>366</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>71 rows × 4 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "     gage_id  Water Year  avg_ann_flow  days_with_data\n",
       "0   11468000        1951    619.506575             365\n",
       "2   11468000        1953    603.172603             365\n",
       "3   11468000        1954    520.427397             365\n",
       "4   11468000        1955    251.441096             365\n",
       "6   11468000        1957    320.472877             365\n",
       "..       ...         ...           ...             ...\n",
       "53  11468000        2004    474.798661             366\n",
       "57  11468000        2008    325.683169             366\n",
       "61  11468000        2012    220.550437             366\n",
       "65  11468000        2016    507.661995             366\n",
       "69  11468000        2020     76.058142             366\n",
       "\n",
       "[71 rows x 4 columns]"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Check for complete years\n",
    "df_year = (df.groupby(['gage_id', 'Water Year'])\n",
    "           .agg({'value':'mean', 'datetime':'count'})\n",
    "           .reset_index()\n",
    "           .rename(columns={'value':'avg_ann_flow', 'datetime':'days_with_data'})\n",
    "          )\n",
    "df_year.sort_values(by=['days_with_data', 'Water Year'], ascending=[True, True])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Calculate the unimpaired flow  \n",
    "if LOI_name == \"Sproul\":\n",
    "    # formula from \"...instream-flow-management\\Reference Files\\SF_Sproul\\Sproul Hydrographs fr Bull Creek v1.xlsm\"\n",
    "    # Documentation from California Trout and Humboldt State University 2018 Sproul Creek Instream Flow Study -- Final Report, prepared for the North Coast Regional Water Quality Control Board and State Water Resources Control Board, Arcata, California. \n",
    "    df['flow'] = (np.where(df['value']<2.706, 0.218*df['value']**2-0.0144*df['value'], df['value']*0.739))/24 * 5\n",
    "    df['flow'] = np.where(df['flow']<0, 0, df['flow'])\n",
    "elif LOI_name == \"Indian\":\n",
    "    # values from \"...instream-flow-management\\Reference Files\\Indian\\Navarro_to_Indian_Creek_Flow.xlsx\"\"\n",
    "    DApod = 39 #sqmiles\n",
    "    DAgage = 303\n",
    "    Ppod = 40.4 # inches from Philo Rain Gage\n",
    "    Pgage = 40\n",
    "    factor = (DApod/ DAgage) * (Ppod/ Pgage)\n",
    "    df_nav = pd.read_csv(nav_div_csv)\n",
    "    df_nav = df_nav.drop(['Month'], axis=1)\n",
    "    df1 = pd.merge(df, df_nav[['Month.1', 'Daily cfs']], left_on='Month', right_on=\"Month.1\", how='left')\n",
    "    df1[\"nav_unimpaired\"] = df1['value'] + df1['Daily cfs']\n",
    "    #Calc Indian Creek unimpaired\n",
    "    df1['flow'] = df1[\"nav_unimpaired\"] * factor\n",
    "    df = df1.drop(['Month.1', \"nav_unimpaired\", 'Daily cfs'], axis=1)\n",
    "    \n",
    "\n",
    "else:\n",
    "    print(\"LOI not recognized\")\n",
    "\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>value</th>\n",
       "      <th>qualifiers</th>\n",
       "      <th>datetime</th>\n",
       "      <th>gage_id</th>\n",
       "      <th>Year</th>\n",
       "      <th>Month</th>\n",
       "      <th>Day</th>\n",
       "      <th>Water Year</th>\n",
       "      <th>flow</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>260.00</td>\n",
       "      <td>A</td>\n",
       "      <td>1950-10-01</td>\n",
       "      <td>11468000</td>\n",
       "      <td>1950</td>\n",
       "      <td>10</td>\n",
       "      <td>1</td>\n",
       "      <td>1951</td>\n",
       "      <td>33.839921</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>260.00</td>\n",
       "      <td>A</td>\n",
       "      <td>1950-10-02</td>\n",
       "      <td>11468000</td>\n",
       "      <td>1950</td>\n",
       "      <td>10</td>\n",
       "      <td>2</td>\n",
       "      <td>1951</td>\n",
       "      <td>33.839921</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>260.00</td>\n",
       "      <td>A</td>\n",
       "      <td>1950-10-03</td>\n",
       "      <td>11468000</td>\n",
       "      <td>1950</td>\n",
       "      <td>10</td>\n",
       "      <td>3</td>\n",
       "      <td>1951</td>\n",
       "      <td>33.839921</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>260.00</td>\n",
       "      <td>A</td>\n",
       "      <td>1950-10-04</td>\n",
       "      <td>11468000</td>\n",
       "      <td>1950</td>\n",
       "      <td>10</td>\n",
       "      <td>4</td>\n",
       "      <td>1951</td>\n",
       "      <td>33.839921</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>260.00</td>\n",
       "      <td>A</td>\n",
       "      <td>1950-10-05</td>\n",
       "      <td>11468000</td>\n",
       "      <td>1950</td>\n",
       "      <td>10</td>\n",
       "      <td>5</td>\n",
       "      <td>1951</td>\n",
       "      <td>33.839921</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25928</th>\n",
       "      <td>0.01</td>\n",
       "      <td>A</td>\n",
       "      <td>2021-09-26</td>\n",
       "      <td>11468000</td>\n",
       "      <td>2021</td>\n",
       "      <td>9</td>\n",
       "      <td>26</td>\n",
       "      <td>2021</td>\n",
       "      <td>0.588453</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25929</th>\n",
       "      <td>0.09</td>\n",
       "      <td>A</td>\n",
       "      <td>2021-09-27</td>\n",
       "      <td>11468000</td>\n",
       "      <td>2021</td>\n",
       "      <td>9</td>\n",
       "      <td>27</td>\n",
       "      <td>2021</td>\n",
       "      <td>0.598853</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25930</th>\n",
       "      <td>0.36</td>\n",
       "      <td>A</td>\n",
       "      <td>2021-09-28</td>\n",
       "      <td>11468000</td>\n",
       "      <td>2021</td>\n",
       "      <td>9</td>\n",
       "      <td>28</td>\n",
       "      <td>2021</td>\n",
       "      <td>0.633953</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25931</th>\n",
       "      <td>0.27</td>\n",
       "      <td>A</td>\n",
       "      <td>2021-09-29</td>\n",
       "      <td>11468000</td>\n",
       "      <td>2021</td>\n",
       "      <td>9</td>\n",
       "      <td>29</td>\n",
       "      <td>2021</td>\n",
       "      <td>0.622253</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25932</th>\n",
       "      <td>0.14</td>\n",
       "      <td>A</td>\n",
       "      <td>2021-09-30</td>\n",
       "      <td>11468000</td>\n",
       "      <td>2021</td>\n",
       "      <td>9</td>\n",
       "      <td>30</td>\n",
       "      <td>2021</td>\n",
       "      <td>0.605353</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>25933 rows × 9 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "        value qualifiers   datetime   gage_id  Year  Month  Day  Water Year  \\\n",
       "0      260.00          A 1950-10-01  11468000  1950     10    1        1951   \n",
       "1      260.00          A 1950-10-02  11468000  1950     10    2        1951   \n",
       "2      260.00          A 1950-10-03  11468000  1950     10    3        1951   \n",
       "3      260.00          A 1950-10-04  11468000  1950     10    4        1951   \n",
       "4      260.00          A 1950-10-05  11468000  1950     10    5        1951   \n",
       "...       ...        ...        ...       ...   ...    ...  ...         ...   \n",
       "25928    0.01          A 2021-09-26  11468000  2021      9   26        2021   \n",
       "25929    0.09          A 2021-09-27  11468000  2021      9   27        2021   \n",
       "25930    0.36          A 2021-09-28  11468000  2021      9   28        2021   \n",
       "25931    0.27          A 2021-09-29  11468000  2021      9   29        2021   \n",
       "25932    0.14          A 2021-09-30  11468000  2021      9   30        2021   \n",
       "\n",
       "            flow  \n",
       "0      33.839921  \n",
       "1      33.839921  \n",
       "2      33.839921  \n",
       "3      33.839921  \n",
       "4      33.839921  \n",
       "...          ...  \n",
       "25928   0.588453  \n",
       "25929   0.598853  \n",
       "25930   0.633953  \n",
       "25931   0.622253  \n",
       "25932   0.605353  \n",
       "\n",
       "[25933 rows x 9 columns]"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Export data\n",
    "df.rename(columns={'datetime': 'DTTM'}, inplace=True)\n",
    "df.drop(['qualifiers', 'value', 'gage_id'], axis=1, inplace=True)\n",
    "df.to_csv(file_out, encoding='utf-8', index=False)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "instream-flow-management",
   "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.15"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}