City Warming Strips .ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [],
   "source": [
    "#!/usr/local/bin/python3.7\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import numpy_indexed as npi\n",
    "import os, netCDF4, time\n",
    "from matplotlib.colors import Normalize\n",
    "\n",
    "os.chdir('/Users/hausfath/Desktop/Climate Science/GHCN Monthly/')\n",
    "\n",
    "berkeley_file = 'Complete_TAVG_LatLong1.nc'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 207,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/usr/local/lib/python3.7/site-packages/ipykernel_launcher.py:5: UserWarning: WARNING: valid_min not used since it\n",
      "cannot be safely cast to variable data type\n",
      "  \"\"\"\n",
      "/usr/local/lib/python3.7/site-packages/ipykernel_launcher.py:5: UserWarning: WARNING: valid_max not used since it\n",
      "cannot be safely cast to variable data type\n",
      "  \"\"\"\n"
     ]
    }
   ],
   "source": [
    "def import_berkeley(filename, verbose = False):\n",
    "    nc = netCDF4.Dataset(filename, 'r')\n",
    "    lats = nc.variables['latitude'][:]\n",
    "    lons = nc.variables['longitude'][:]\n",
    "    temps = nc.variables['temperature'][:,:,:]\n",
    "    times = nc.variables['time'][:]\n",
    "    years = times.astype(int)\n",
    "    return {\n",
    "        'lats' : lats,\n",
    "        'lons' : lons,\n",
    "        'temps' : temps,\n",
    "        'times' : times,\n",
    "        'years' : years\n",
    "    }\n",
    "\n",
    "\n",
    "def local_timeseries(data, local_lat, local_lon):\n",
    "    near_lat = find_nearest(data['lats'], local_lat)\n",
    "    near_lon = find_nearest(data['lons'], local_lon)\n",
    "    near_lat_pos = np.where(data['lats'] == near_lat)[0][0]\n",
    "    near_lon_pos = np.where(data['lons'] == near_lon)[0][0]\n",
    "    try:\n",
    "        anoms = np.swapaxes(data['anoms'],0,2)\n",
    "        anoms = np.swapaxes(anoms,0,1)\n",
    "    except:\n",
    "        anoms = np.swapaxes(data['temps'],0,2)\n",
    "        anoms = np.swapaxes(anoms,0,1)\n",
    "    return anoms[near_lat_pos][near_lon_pos]\n",
    "\n",
    "\n",
    "def find_nearest(array, value):\n",
    "    idx = (np.abs(array - value)).argmin()\n",
    "    return array[idx]\n",
    "\n",
    "\n",
    "def berkeley_local_annual(filename, local_lat, local_lon, data):\n",
    "    local_data = local_timeseries(data, local_lat, local_lon)\n",
    "    unique, mean = npi.group_by(data['years']).mean(local_data)\n",
    "    df = pd.DataFrame({'years' : unique,\n",
    "                       'temps' : mean})\n",
    "    return df\n",
    "\n",
    "#results = berkeley_local_annual(berkeley_file, doha_lat, doha_lon)\n",
    "\n",
    "data = import_berkeley(berkeley_file)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 206,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "138 189\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAd4AAAD7CAYAAADEiPXyAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAACUVJREFUeJzt3FuMXWUZx+H/tNOZttAThWIptbQdtQGkYAlQTq1AEFCMaDEqJ4MEjFFjRNEQT6AGQwgSwEACIhhi1YAGY4ByEAqiEMEKaQmFUtsyzLTFlp7PnfHiq4neevFywfNcTWb2Xutb7/rW/GbfTMfg4GAAgBpD3ukFAMC7ifACQCHhBYBCwgsAhYQXAAoJLwAUEl4AKCS8AFBIeAGgkPACQKHOd+rEf1mxbnDued9Mkix7/Oa8um57kqR/884se2tLkuSzRx+SJOndtDO3Lnw9SXLbvA9m5LD298KegfbvLp9auTGbdu5Jkpw6ZWxWbdyZJHnlX1uTJAuXrk1X59Akyezp47Njz94kyaKVG5Ikd57anb6RU5IkL/RtTpKcPW10rn70n0mSaRP2y2UzJyRJunr/kSR5csj7M3fg1STJhc925e5zJydJZl33QpJkQ/9b6Zk1LUlyz8Ufypotu/9nTR+ZPi6jutuatuwayKjOdi0dz96fJFk985PZve/6Dut9Jpt6Tk2SzPnhE0mS578zO529L7U5TJmVX7/SruXM6Qe0dVzxy6y66r1Jkt1HnJ4zb/lrkuSa845MkgwMJqe88VCS5NiHx7djHz85z77UnySZN3dabrixreWO6y5Jkhw3aXQOHNbmfPF9S3PtOTOSJJMevSlJsuCoy/LRnnFJkiFb12XDsPb19j3tOiZ27sj+c69Kknzt2ivz3dPbfB5a9naS5IiD9s9hT96cJLl32kU58b1jkyQ/WrA0SfLiov489v3TkySX3LsoD1x2bJLkx0+0+zTvqENy/ePtnvSu3pI/HdebJLny7aOTJFecOCXL1m1Lkjzyytr2uvXbcv/5U5Mk585/PZ87vs3spt8tTpLcdvnxOWxsd5Lkewtey21Htf1x5iNJktz+maMztb/NtqN7RGbeuTFJ8sI1c9ocNq/Nkr1tvjOW3J/eYz7d1jK/7aPe19ZlyWUjkyRbe07JwpXt/f2b2x4+q2d8ujs7kiSL127LN37xtyTJw99qx1+wbH2OmTg6SfLyvufmrOnj8nx/22eTx3TnE9c9mSQZe9B+SZLRY4fnP/Yf3plLTmh7f+q4EUmSoR0d+f2S1UmSvQODuXTWpCTJhIG2xzYMG5cxHW19Ny9an3EjhiVJTpvW9t4h+3Xmub52/gNHdmX99rb3f/X3dj+279qbq07rSZJ8YMfy/GHLwUmSkya361i9dU862iWnZ1x3vv3Qa23+Z0xPktzw1IpcO6vdk7yxpP1sbU+u/nC7jz99ZlU+deTEdvysSZL0dU/MxN1vJUnuWjEkXzi07YMvPdOu48snT838RW8mSX6w36K8ccTHkyRj9j2jYwa3Ze3eNreNOwdywS1/TpJM72n39oTp43PlxHauZzpnZPaI9UmSHWMOTZKMWPFcOrra+/esWZUhU2fmv73ZPSkTB9t873ujfW/H3oH87IGX2zGvnpNNuwaSJAcsf7odZ8ac5Im7kyRf3HB8+je036HzLzqmzXvxmlwxfd8gBwfz2mBba9++vdU9dEjeN76t6cXVW3Pw/l1JkiM72nWkY0j2jn5PkmToxr78dk3bPyftey6fWrkhZ0xrz/iEt5cm2zclSRZ2HZEkWb5+W+ZObXvizU07M7vvsbbu1avakvbsys6PfT1JMvzhW5MkXSefl6fPuTBJMuc3N+aOt9v8Ls2iJMnA9q0Z2LQuSfLKXX/M4V9tr+178NG2zmGdmXT5V9rxh4/K1gfvSZKsW7w8SXLwcYdn6Jg2h5XHXpBtl89Lksz4/Nnt/LPOyM+POj9JMu2QUZkxr/3u2LNjV5LkJ9cvzO2DK/YN9f/nEy8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIWEFwAKCS8AFBJeACgkvABQSHgBoJDwAkAh4QWAQsILAIU6BgcH3+k1AMC7hk+8AFBIeAGgkPACQCHhBYBCwgsAhYQXAAoJLwAUEl4AKCS8AFBIeAGgkPACQCHhBYBCwgsAhYQXAAoJLwAUEl4AKCS8AFBIeAGgkPACQCHhBYBCwgsAhYQXAAoJLwAU+jf22J4A3pjEdgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x216 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "#Name of CSV file containing the data (one column only)\n",
    "\n",
    "savename = 'Stuttgart'\n",
    "\n",
    "lat = 48.7758\n",
    "lon = 9.1829\n",
    "\n",
    "results = berkeley_local_annual(berkeley_file, lat, lon, data)\n",
    "\n",
    "year_loc = np.where((results['years'] >= 1750) & (results['years'] <= 2018))[0]\n",
    "temps = results['temps'][year_loc] \n",
    "#np.genfromtxt(temp_data, delimiter=\",\")[1:]\n",
    "temps_normed = ((temps - temps.min(0)) / temps.ptp(0)) * (len(temps) - 1)\n",
    "elements = len(temps)\n",
    "\n",
    "x_lbls = np.arange(elements)\n",
    "y_vals = temps_normed / (len(temps) - 1)\n",
    "y_vals2 = np.full(elements, 1)\n",
    "bar_wd  = 1\n",
    "\n",
    "my_cmap = plt.cm.RdBu_r #choose colormap to use for bars\n",
    "norm = Normalize(vmin=0, vmax=elements - 1)\n",
    "\n",
    "def colorval(num):\n",
    "    return my_cmap(norm(num))\n",
    "\n",
    "fig=plt.figure(figsize=(6,3))\n",
    "plt.axis('off')\n",
    "plt.axis('tight')\n",
    "\n",
    "#Plot warming stripes. Change y_vals2 to y_vals to plot stripes under the line only.\n",
    "plt.bar(x_lbls, y_vals2, color = list(map(colorval, temps_normed)), width=1.0, edgecolor = \"none\")\n",
    "\n",
    "#Plot temperature timeseries. Comment out to only plot stripes\n",
    "#plt.plot((x_lbls + 0.5), y_vals - 0.002, color='black', linewidth=2)\n",
    "\n",
    "plt.xticks( x_lbls + bar_wd, x_lbls)\n",
    "plt.ylim(0, 1)\n",
    "fig.subplots_adjust(bottom = 0)\n",
    "fig.subplots_adjust(top = 1)\n",
    "fig.subplots_adjust(right = 1)\n",
    "fig.subplots_adjust(left = 0)\n",
    "fig.savefig(savename+'.png', dpi=300)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.7.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}