5SSG2060_Practical_2_Processing spatial data.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# This Week’s Overview\n",
    "\n",
    "\n",
    "You will further developing you skills of data handling and processing in this practical by finishing 10 small tasks about:\n",
    "\n",
    "- Working with streaming data. You will load weather data through API. It is harder conceptually because API data is harder to understand -- we've simplified it quite a bit but it's still got some parts that are going to be hard going.\n",
    "- Making an interactive map using skills learned last week.\n",
    "- Converting queried API data into a well-formatted DataFrame. Some DataFrame operation such as 'join','append',’merge’ you learned in last term will be used.\n",
    "- Create a ShapeFile using Shapely and GeoPandas - get a sense of projected and geographic coordinate system.\n",
    "- Calculating distance between points."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Learning Outcomes\n",
    "\n",
    "By the end of this practical you should have:\n",
    "- Used API to request streaming data\n",
    "- Enhanced your skills of manipulating data frame \n",
    "- Known calculate geographical distance in various ways"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "(always remember) The first thing we need to do is setup our working environment. Run the scripts below to import the relevant libraries."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "import folium  \n",
    "import os\n",
    "\n",
    "print (folium.__version__)\n",
    "\n",
    "%matplotlib inline\n",
    "\n",
    "pd.set_option('display.max_rows', 300) # specifies number of rows to show\n",
    "pd.options.display.float_format = '{:5,.4f}'.format # specifies default number format to 4 decimal places\n",
    "\n",
    "import warnings\n",
    "warnings.simplefilter('ignore')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Weather Data  - background knowledge you should know\n",
    "\n",
    "The UK's Met Office is a world-leading weather and climate research centre, and even if it doesn't always seem like their forecasts are very accurate that's because Britain's weather is inherently _unpredictable_. They've also done a lot of work to make their weather data widely available to people like us.\n",
    "\n",
    "I probably don't need to say a _lot_ about weather data because you've probably been making use of forecasts for much of your life! But it's _still_ worth understanding something about how weather data is gathered and reported: many organisations operate weather stations where data on wind speed, temperature, rain, and amount of sun are collected and then transmitted to a server to be integrated into a larger data set of weather _observations_ at a national or global scale. Of course, any _one_ station might be in the 'wrong' place (somewhere shady or protected from the rain) or it might even break down, but the idea is that if you have enough of them you can collect a pretty good range of data for the country and begin to look for patterns and, potentially, make predictions.\n",
    "\n",
    "We will be accessing data from the MetOffice from weather stations where observations, such as the ones below, are collected:\n",
    "* <Param name=\"F\" units=\"C\">Feels Like Temperature (units: degrees Celsius)\n",
    "* <Param name=\"G\" units=\"mph\">Wind Gust (units: mph)</Param>\n",
    "* <Param name=\"H\" units=\"%\">Screen Relative Humidity (units: percent)</Param> \n",
    "* <Param name=\"T\" units=\"C\">Temperature (units: degrees Celsius)</Param> \n",
    "* <Param name=\"V\" units=\"\">Visibility (units: km?)</Param> \n",
    "* <Param name=\"D\" units=\"compass\">Wind Direction (units: compass degrees)</Param>  \n",
    "* <Param name=\"S\" units=\"mph\">Wind Speed (units: mph)</Param> \n",
    "* <Param name=\"U\" units=\"\">Max UV Index (units: index value)</Param> \n",
    "* <Param name=\"W\" units=\"\">Weather Type (units: categorical)</Param> \n",
    "* <Param name=\"Pp\" units=\"%\">Precipitation Probability (units: percent)</Param>\n",
    "\n",
    "These observations are associated with a particular station (where did we see these values/where _will_ we see these values?), they will also be associated with _either_ a particular time in the past (when were they collected?) or, if they're forecasts, with a particular time in the future (when do we expect to see them?). \n",
    "\n",
    "So although weather data might seem more 'objective' than data on social class (though for obvious reasons it turns out that both are just attempts to capture data about reality, not reality itself), it may also turn out to be very complex to store and manage beccause of the temporal element _and_ the fact that it's not just a count of one thing, each of these observations uses a very different set of units.\n",
    "\n",
    "To really get to grips with the MetOffice API you will need to RTM (Read The Manual): http://www.metoffice.gov.uk/media/pdf/3/0/DataPoint_API_reference.pdf. The ruder version of that, which you will sometimes see on StackOverflow and elsewhere, is: RTFM."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Getting Weather Data via an API\n",
    "\n",
    "Because the weather is changing all the time, so is the data! And, 'worse', it's becoming obsolete: the forecast from 2 years ago isn't particularly useful to us now. *And* asking for \"yesterday's weather\" depends on the day that we're asking! When you have data that is always changing from minute to minute or day to day then you use an API (Application Programming Interface) to access it: the API knows that \"yesterday's weather\" means \"work out what day it is right now and then get the weather from the day before\", and it also knows that \"give me the current weather from station X\" means \"look up station X and find the latest weather report that I've received\". In other words, an API is  designed with programmatic, dynamic interaction in mind right from the start.\n",
    "\n",
    "#### So What _is_ an API?\n",
    "\n",
    "There's a nice, friendly introduction to APIs over at [Free Code Camp](https://medium.freecodecamp.com/what-is-an-api-in-english-please-b880a3214a82#.rmjnli2nn). \n",
    "\n",
    "Helpfully, the MetOffice provides a lot of documentation about their API (I'd suggest bookmarking it): http://www.metoffice.gov.uk/datapoint/support/api-reference\n",
    "\n",
    "This type of data requires a lot more research up front to work with, but it's very flexible once you know how to 'speak API' because you can _customise_ the API request so that the server responds with _only_ the data we're interested in instead of being 'stuck' with what the provider wants to give us.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "### Task 1: Obtaining an API Key\n",
    "\n",
    "The first step to working with the API from the MetOffice is to obtain an API key: [Click to register and obtain the API Key here](http://www.metoffice.gov.uk/datapoint/API). \n",
    "\n",
    "Make a note of this key in your notebook. Right here:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "api_key   = \"???\" # your API key"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "That way your API key is saved somewhere easy to access.\n",
    "\n",
    "We _always_ have to use the key as part of an API request: the process by which we _ask_ for data. Think of the key as being _your_ unique identifier: no two people share the same key and that way the MetOffice can cut off people who abuse the system or look at which APIs are popular with lots of users... **Twitter and Facebook do the same thing.**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Task 2: Obtaining a List of Sites from the API\n",
    "\n",
    "How to start? Well, the first thing that we need to know is: for what locations can I get weather data? For this to work, we need to know how to ask the API nicely for a list of available sites... "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "First we import two new modules: one that makes requests to a web server, and one that will parse JSON* responses from the server in order to turn them into something that we can work with more easily."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import json, requests # Libraries we need"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then we set up some default variables (api_url, site_url) that will help us to build our request to the MetOffice's server. The comments help us to remember what each of these variables holds."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "api_url  = \"http://datapoint.metoffice.gov.uk/public/data/\" # base URL\n",
    "site_url = \"val/wxobs/all/json/sitelist\" # sites API URL\n",
    "payload = {'key': api_key} # Dictionary to hold request parameters"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You'll notice that the payload is just a dictionary and that this dictionary is then passed to the requests library (the get function). All it does it convert this dictionary to a key-value pair in the format expected by the API. Think of it as a kind of translation between languages: from the language of Python to the language of the web (HTTP, to be precise).\n",
    "We issue our request and it returns a response that we store in s (short for sites data)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s = requests.get(api_url + site_url, params=payload) # Do the request"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Lastly, we ask the response object to convert the reply into a JSON data structure... more on JSON in a second, but first let's look at what we got from our request!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sites = s.json() # Capture the output\n",
    "print(\"Done!\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Show the requested URL\n",
    "# your code here\n",
    "???\n",
    "# Click on the link below to see it nicely formatted automatically!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Capture the returned data\n",
    "# your code here\n",
    "???"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Task 3: Parse JSON files\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We have played around GeoJSON file last week. JSON is just like GeoJSON, but without embedded geographic data structures. Now that you have a better idea of what data the reply contains, let's see if we can convert the JSON reply into something useful for Python; if you scroll back up to where we printed out the reply you'll notice that it all starts with a '{', meaning that it's a dictionary. \n",
    "\n",
    "Yes, look just like the GeoJSON file. (So you should known how to parse it)\n",
    "\n",
    "Let's start by printing out the keys in the dictionary and the _type_ of data associated as a value to that key:\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for k in sites.keys():\n",
    "    print(\"Key: \" + ???(k))\n",
    "    print(\"Value: \" + str(???(sites[k])))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Not the most exciting answer, but at least we know that the value is a dictionary. Let's try moving down a level:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for k in sites.keys():\n",
    "    print(\"Key: \" + ???(k))\n",
    "    print(\"Value: \" + str(???(sites[k])))\n",
    "    for k2 in sites[k].keys():\n",
    "        print(\"\\tKey: \" + ???(k2))\n",
    "        print(\"\\tValue: \" + str(???(sites???)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**This isn't great design in the reply: we have a dictionary that contains only one key/value pair, and _that_ dictionary in turn contains only one key/value pair. But after that we get to a long, long list containing the data...**\n",
    "\n",
    "The MetOffice is not making life easy for us here: there's a _lot_ of extra 'baggage' in this API response. But we at least know that the next level down is a list and _that_ suggests that things are about to get a bit more interesting... Let's simplify our code at the same time:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "apiList = ???['Locations']['Location']\n",
    "print(\"List in API response is \" + ???(???(apiList)) + \" long\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now _that_ is a rather more interesting response. What it means is that our JSON reply has this structure:\n",
    "```\n",
    "{\n",
    "    'Locations': {\n",
    "        'Location': [\n",
    "            ... lots of data here ...\n",
    "        ]\n",
    "    }\n",
    "}\n",
    "```\n",
    "\n",
    "If you scroll back up to the JSON reply you should now be able to read a little bit more of the response... and this should give you a clue as to how to print out the `name`, `id`, and `longitude` of the first five sites. I'll get you started:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for i in range(5):\n",
    "    location = apiList[???]\n",
    "    print(\"Location: \" + location[???] + \" (id: \" + location[???] + \") is at longitude: \" + location[???])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Your answer should look like this:\n",
    "\n",
    "`\n",
    "Location: Kinloss (id: 3066) is at longitude: -3.5606\n",
    "Location: Lossiemouth (id: 3068) is at longitude: -3.322\n",
    "Location: Wick John O Groats Airport (id: 3075) is at longitude: -3.089\n",
    "Location: Baltasound (id: 3002) is at longitude: -0.854\n",
    "Location: Lerwick (S. Screen) (id: 3005) is at longitude: -1.183\n",
    "`\n",
    "\n",
    "And _now_ that we know where all the data was 'hidden', we can convert this to a proper data structure in which it is possible to _interact_ with it. To do that, we'll put the site data into a pandas data frame..."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Task 4: Turning API data into a Pandas DataFrame\n",
    "\n",
    "Pandas is remarkably intelligent and will _often_ -- though not always -- work out the sensible thing to do from many kinds of data structures (list-of-lists, dictionary-of-lists, list-of-dictionaries...). So let's see what happens when we simply pass `apiList` (a LoD) directly to the `DataFrame` 'constructor' instead of passing the data through, for instance, the `read_csv` function as we did above with a CSV file."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "site_df = pd.???(???)\n",
    "print(site_df.head())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Wow, that's... almost scarily easy. You can see that pandas worked out the structure of our LoD and then automatically converted that to columns in a data frame. So it got the hardest part of the process exactly right and has saved us a lot of work. _That_ is the point of functions and of code: to be constructively lazy.\n",
    "\n",
    "Of course, we could have predicted that pandas would cope since there is a whole section in the documentaiton [devoted to creating data frames from different structures](http://pandas.pydata.org/pandas-docs/stable/dsintro.html#dataframe).\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The problem is that pandas didn't know what we expected the columns to be, so it's treated them all as 'objects' (basically: strings) and not as numeric data types. To fix that you need to know that there's a function called `'astype'` that allows pandas to convert between different data types where it's fairly easy for pandas to figure out what we want to do:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for c in ['region','unitaryAuthArea']:\n",
    "    site_df[c] = site_df[c].???('str')\n",
    "for c in ['elevation','latitude','longitude']:\n",
    "    site_df[c] = site_df[c].???(???)\n",
    "for c in ['id']:\n",
    "    site_df[c] = site_df[c].???(???)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "site_df.describe() # only the numeric paratmers shows here"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Task 5: Making an interactive map using Folium"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Before we go further, lets make use of your mapping skills learned last week to visually identify where are the locations collects weather data."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "With the map to be created, we are able to answer the first question: for what locations can I get weather data? - Let's add all sites (points) to the map."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "MAP_COORDINATES = (51.5113, -0.1160) \n",
    "\n",
    "m = folium.???(location=???, zoom_start=8, tiles=\"Stamen Toner\")\n",
    "\n",
    "for index, loc in site_df.iterrows():\n",
    "    # And now add a marker\n",
    "    folium.???((???['latitude'], ???['longitude']), icon=folium.Icon(color='green',icon='cloud'), \n",
    "              popup = 'weather station'\n",
    "             ).add_to(???)\n",
    "m.save('Weather Station Map.html') # Print the map"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Zoom out, and what did you see there?**\n",
    "\n",
    "The weather station are located across the whole UK, not just in London."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can use spatial join to find out all stations in GLA, but not there yet. Lets use the key word 'unitaryAuthArea' to query out stations in 'Greater London'."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Lets find out the weather station in the Greater London area by querying a key word\n",
    "site_df_london = site_df.???[(site_df.??? =='Greater London')]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "MAP_COORDINATES = (51.5113, -0.1160) \n",
    "\n",
    "m = folium.???(location=???, zoom_start=10, tiles=\"Stamen Toner\")\n",
    "\n",
    "for index, loc in site_df_london.iterrows():\n",
    "    # And now add a marker\n",
    "    folium.???((???['latitude'], ???['longitude']), icon=folium.Icon(color='green',icon='cloud'), \n",
    "              popup= loc['unitaryAuthArea']\n",
    "             ).add_to(???)\n",
    "m.save('Weather Stations in Greater London.html') # Print the map"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Task 6: Obtaining Weather Data in London\n",
    "\n",
    "The next step in this process is a bit more complicated because weather data is a bit more complicated than a list of locations...\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "Here, the MetOffice has _not_ made our lives very easy because the data is packaged in a way that doesn't allow us to easily load it into pandas. If you search online, you'll find plenty of people complaining about how the MetOffice API works. Or doesn't work, if you prefer.\n",
    "\n",
    "So we're not going to ask you to sort this out for yourselves. Instead, we're going to provide you with a function (!) to take the observation data and convert it into a data frame.\n",
    "\n",
    "Well, you should at least know how to use a function."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from datetime import datetime, timedelta\n",
    "def processMetOfficeObservations(loc): \n",
    "    \"\"\"\n",
    "    Process a series of 'reports' for a single\n",
    "    location using the datetime object as the \n",
    "    reference time against which to build the \n",
    "    timedelta (i.e. we start from midnight and \n",
    "    the timedelta is the number of minutes past \n",
    "    midnight)\n",
    "    \"\"\"\n",
    "    observations = {} # Stores results\n",
    "    \n",
    "    for d in loc['Period']: # d for day\n",
    "        if 'value' not in d or 'Rep' not in d:\n",
    "            # print 'no \"value\" or \"Rep\" key, skip'\n",
    "            continue\n",
    "            \n",
    "        dt = datetime.strptime(d['value'],'%Y-%m-%dZ')# Convert date to datetime object \n",
    "        \n",
    "        # Now deal with the actual observations (i.e. 'Reports')\n",
    "        the_type = type(d['Rep'])\n",
    "        if the_type is dict:\n",
    "            reports = [d['Rep']]\n",
    "        elif the_type is list:\n",
    "            reports = d['Rep']\n",
    "        else:\n",
    "            print (\"***warning: d['Rep'] type: \", the_type, \" not support\")\n",
    "            continue\n",
    "        \n",
    "        for report in d['Rep']:\n",
    "            \n",
    "            # Find the timestampe and add it to the date\n",
    "            minutes_after_midnight = int(report['$'])\n",
    "            ts = dt + timedelta(minutes=minutes_after_midnight)\n",
    "                \n",
    "            # For each of the possible values, set a default value\n",
    "            # if the weather station doesn't actually collect that\n",
    "            # parameter... can you see a problem with our defaults?\n",
    "            if 'ts' not in observations:\n",
    "                observations['ts'] = []\n",
    "            observations['ts'].append(str(ts))\n",
    "            for key in ['D','Pt']:\n",
    "                if key not in report:\n",
    "                    report[key] = u\"\"\n",
    "                if key not in observations:\n",
    "                    observations[key] = []\n",
    "                observations[key].append(report[key])\n",
    "            for key in ['W','V','S','G']:\n",
    "                if key not in report or report[key] == \"\":\n",
    "                    report[key] = 0\n",
    "                if key not in observations:\n",
    "                    observations[key] = []\n",
    "                observations[key].append(report[key])\n",
    "            for key in ['T','Dp','H']:\n",
    "                if key not in report or report[key] == \"\":\n",
    "                    report[key] = 0.0\n",
    "                if key not in observations:\n",
    "                    observations[key] = []\n",
    "                observations[key].append(report[key])\n",
    "    \n",
    "    return observations"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "First, just in case you want to only run this section again (and not revisit the content above), I'd suggest saving a copy of your API key here as well:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "api_key = \"???\" # Here you need to replace this with your unique API key"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's start from a simple task. - querying weather data at one weather station. And we choose one located in London."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "london_ids = (site_df_london.id).unique()\n",
    "london_ids  \n",
    "# how many stations with hourly updated data in London? You have got the answer from your map."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Just like what we did before for querying site locations. Here, we just added one small parameter - location id.\n",
    "\n",
    "**if you are unsure where to add it, go to check data point [API](https://www.metoffice.gov.uk/datapoint/support/api-reference)**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import json, requests # Libraries we need\n",
    "\n",
    "api_url  = \"http://datapoint.metoffice.gov.uk/public/data/\" # Base URL\n",
    "obs_json= \"val/wxobs/all/json/\" # Observations URL\n",
    "\n",
    "station = str(3772)  #3772 # This is heathrow airport weather station\n",
    "\n",
    "payload = {'res': 'hourly', 'key': api_key} # Dictionary to hold request parameters\n",
    "\n",
    "# pay attention here, station id added\n",
    "r = requests.get(??? + ??? + ???, ???=payload)\n",
    "\n",
    "#print(r.url)\n",
    "\n",
    "weather = r.json() # Capture the reply\n",
    "\n",
    "print(\"Done!\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "print(weather)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#using the function\n",
    "data = processMetOfficeObservations(weather['SiteRep']['DV']['Location'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df3 = pd.???.from_dict( data )\n",
    "# add one more parameter to indicate the station id - say, where the data from\n",
    "df3['id'] = station\n",
    "df3.head(4)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Task 7: Obtaining Weather Data in the UK ( its your turn)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here you are going to repeat what has been done in Task 6. But instead of querying only 1 station, you need to write a loop to query all weather station in the UK and put the data into a DataFrame called weather_df.\n",
    "\n",
    "Hints: first you get all unique ids of weather stations; then you write a loop to iterate all station ids; Each time you request data of one station, you convert the JSON into a new DataFrame; then you join the new dataframe to the weather_df (not sure what operation to use, then check [FROM HERE](https://pandas.pydata.org/pandas-docs/stable/merging.html)); After you gather all data from all station, rename the columns to make it easy understood and reset index to make it ready for query. That's it!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "weather_df = ???.iloc[0:0] # so the big container is prepared for you, you just have to fill it with data.\n",
    "weather_df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# get station ids here\n",
    "uk_ids = site_df.???.unique()\n",
    "# get the total number of records\n",
    "???(uk_ids)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "the answer should be something like 148 (data updates all the time as we are using streaming data)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# loop here\n",
    "for stn in uk_ids:\n",
    "    station = ???(stn)\n",
    "    r = requests.get(??? + ??? + ???, ???=payload)\n",
    "    weather = ??? # Capture the reply    \n",
    "    if 'Location' not in weather['SiteRep']['DV']: # if there is no keyword - location, means no data coming in\n",
    "        continue\n",
    "    else:\n",
    "        print('weather station id is %s' %(???) )\n",
    "        test = processMetOfficeObservations(weather['SiteRep']['DV']['Location'])\n",
    "        df3_new = pd.???.from_dict(???)\n",
    "        df3_new['id'] = station \n",
    "        weather_df = weather_df.???(df3_new)  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#change column names here. \n",
    "cnames = {'D' : 'WindDirection', \n",
    "          'G' : 'WindGust', \n",
    "          'Dp': 'DewPoint', \n",
    "          'H' : 'Humidity', \n",
    "          'P' : 'PressureTendency', \n",
    "          'S' : 'WindSpeed', \n",
    "          'T' : 'Temperature', \n",
    "          'W' : 'WeatherType', \n",
    "          'V' : 'Visibility', \n",
    "          'Pt': 'PressureTendency',\n",
    "          'id': 'id',\n",
    "          'ts': 'DateTime'\n",
    "         }\n",
    "\n",
    "weather_df.???(???=cnames, inplace=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "#check you data, see if it makes sense\n",
    "???"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Did you notice that the row index is not with unique numbers? Why? This is a common issue when combining DataFrame together and may cause problems when we index data. We just need to `reset the index`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "weather_df = weather_df.???(drop=True)\n",
    "weather_df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#how many records you got? \n",
    "???(weather_df)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "you should get an answer around 3471. Again number changes slightly as we are using streaming data."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Task 8: Join Weather Data (weather_df) with Site Data (df2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The weather data we get does not have location tag. You can just use the weather station id and always go back the site location table to find the coordinates. It is a good strategy for data storage, especially for large data sets. And it is probably a bad idea to combine them into one, as in such case both station id and timetag cannot be used as index any more.  \n",
    "\n",
    "But (just) for practice, let's combine them into one to create a spatial file for later use."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "weather_df.head(3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "site_df.head(3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "**What is the common key word for joining two tables?**\n",
    "\n",
    "**What sort of join should we perform? inner? outter? left? or right?**\n",
    "\n",
    "**You can also do it use merge, if you forgot, go to geocomputation week 6 - practical**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For your convinience, read more about merge [here](http://pandas.pydata.org/pandas-docs/stable/merging.html). \n",
    "\n",
    "![Illustration of the Pandas merge function](http://pandas.pydata.org/pandas-docs/stable/_images/merging_merge_on_key.png)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df_in_one = pd.???(weather_df, site_df, on = '???', how = 'left')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You may get a ValueError:You are trying to merge on object and int64 columns. If you wish to proceed you should use pd.concat.\n",
    "\n",
    "Why??? Think really hard about that before you raise your hand to get an answer from us."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "site_df.??? # check the data types for site_df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "weather_df.??? # check the data types for weather_df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for c in ['WindDirection','PressureTendency']:\n",
    "    weather_df[c] = weather_df[c].astype('str')\n",
    "for c in ['DewPoint','WindGust','Humidity','WindSpeed','Temperature','Visibility','WeatherType']:\n",
    "    weather_df[c] = ???.to_numeric(weather_df[c], errors='coerce') # change numeric data into integar and float types."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "weather_df.??? # check the changes of data types for weather_df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "weather_df['id'] = weather_df['id'].astype('int')\n",
    "df_in_one = pd.???(weather_df, site_df, how = 'left', on = '???')\n",
    "df_in_one.head(5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Task 9: Generate Shapefile/Geojson File from DataFrame"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The examples we gave here probably not a great one. But it is true that quite often we are in the situation that no spatial file is available and we have to generate one by yourself. \n",
    "\n",
    "You have learned how to use QGIS to import CSV and generate a vector map from there. Python can do that as well, with simple steps. One more spatial data package you will learn here is Shapely. We will first show you an example and you will create a shapefile of site station following the example.\n",
    "\n",
    "It should be installed already, if not, do it by yourself."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "step 1: you use the package shapely to generate points from latitute and longitude. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from shapely.geometry import Point\n",
    "\n",
    "# combine lat and lon column to a shapely Point() object\n",
    "df_in_one['geometry'] = df_in_one.apply(lambda x: Point((float(x.longitude), float(x.latitude))), axis=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# there is a column added as geometry\n",
    "df_in_one.head(5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "step 2: Then you use Geopandas to generate 'Geopandas dataframe' from 'Pandas dataframe'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import geopandas\n",
    "df_in_one = geopandas.GeoDataFrame(df_in_one, geometry=???)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "step 3: Then you set a projected and geographic coordination system here to convert the file into a shapefile"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df_in_one.crs= \"+proj=longlat +ellps= EPSG:27700 +datum=WGS84 +no_defs\"\n",
    "df_in_one.to_file('weather.shp', driver='ESRI Shapefile')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "step 4: The generated shapefile should be in the same folder as your ipython notebook. and the shapefile can be openned by other GIS tools, such as QGIS. Try it. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "???(df_in_one) # how many records this time?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You have learned how to use QGIS to import CSV and generate a vector map from there. Python can do that as well, with simple steps. One more spatial data package you will learn here is Shapely. \n",
    "\n",
    "It should be installed already, if not, do it by yourself.\n",
    "\n",
    "\n",
    "Though there are c.3531 records, there are only a few points shown in the map. why? because may records have the exact same coordinates. They are weather records from one weather station collected at different time. So there should be only 148 points (== number of weather stations) shown in the map.\n",
    "\n",
    "There are two directions for your to further explore **after the practical**. \n",
    "\n",
    "1 - go for time series analysis, take one weather station as example - time series analysis and modeling is a challenge task to do.\n",
    "2 - generate a series of maps by time series - which sounds more spatial data related, let's try this using what you have learned in Geocomputation week 8 - making maps."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Could you prepare a shapefile of weather station locations using just 'site_df' for task 10, following what we did to generete weather shapefile?**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "site_df['???'] = site_df.apply(lambda x: Point((float(x.longitude), float(x.latitude))), axis=1)\n",
    "site_df = ???.GeoDataFrame(site_df, geometry=???)\n",
    "site_df.crs= \"+proj=longlat +ellps=EPSG:27700 +datum=WGS84 +no_defs\"\n",
    "site_df.to_file('weather_station.shp', driver='ESRI Shapefile')\n",
    "???(site_df)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "### Task 10: Generating a Distance Matrix "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A distance matrix is a square matrix (two-dimensional array) containing the distances, taken pairwise, between the elements of a set. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We will calculate a distance matrix for all weather stations using methods we mentioned in lecture this week."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "site_df.head(5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Could you calculate the distance between station id - 3066 to station id 3068?**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#x - longitude y - latitude\n",
    "\n",
    "ori = [{'y': site_df.iloc[???]['latitude'], 'x':site_df.???[???]['longitude']}]\n",
    "des = [{'y': site_df.iloc[???]['latitude'], 'x':site_df.???[???]['longitude']}]\n",
    "\n",
    "ori = pd.DataFrame(ori)\n",
    "des = pd.DataFrame(des)\n",
    "\n",
    "# to make your life easier - ori is the station 3066, des is the station 3068"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To verify what we taught in the class, could you complete the function below? Calculate the distance based on a Euclidean method. the squre_root method.\n",
    "\n",
    "**implement Euclidean method**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import math\n",
    "def distance_euclidean(ori,des):\n",
    "    \"\"\"\n",
    "    calculate the euclidean distant between two points. -  \n",
    "    \"\"\"\n",
    "    distance = np.sqrt((ori[???]-des[???])**2+(ori[???]-des[???])**2)\n",
    "    return distance"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "haversine method, using package `pyproj`, and cosines method are provided, for you the compare the results. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import math\n",
    "from math import radians, cos, sin, asin, sqrt\n",
    "\n",
    "def distance_haversine(ori,des):\n",
    "    \"\"\"\n",
    "    Calculate the great circle distance between two points \n",
    "    on the earth (specified in decimal degrees)\n",
    "    \"\"\"\n",
    "    # convert decimal degrees to radians \n",
    "    lon1, lat1, lon2, lat2 = map(radians, [ori.x, ori.y, des.x, des.y])\n",
    "\n",
    "    # haversine formula \n",
    "    dlon = lon2 - lon1 \n",
    "    dlat = lat2 - lat1 \n",
    "    a = sin(dlat/2)**2 + cos(lat1) * cos(lat2) * sin(dlon/2)**2\n",
    "    c = 2 * asin(sqrt(a)) \n",
    "    r = 6371 # Radius of earth in kilometers. Use 3956 for miles. Check the value by yourself.\n",
    "    return c * r"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pyproj\n",
    "\n",
    "def distance_pyproj(ori, des):\n",
    "    geod = pyproj.Geod(ellps='WGS84')\n",
    "    angle1,angle2,distance = geod.inv(ori.x[0], ori.y[0], des.x[0], des.y[0])\n",
    "    print (\"{0:8.4f}\".format(distance/1000))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def distance_cosines(ori, des):\n",
    "    # law of cosines\n",
    "    distance = math.acos(math.sin(math.radians(ori.y))*math.sin(math.radians(des.y))+math.cos(math.radians(ori.y))*math.cos(math.radians(des.y))*math.cos(math.radians(des.x)-math.radians(ori.x)))*6371\n",
    "    print (\"{0:8.4f}\".format(distance))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "???(ori,des) # function you defined above"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "???(ori, des) # function you defined above"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "???(ori, des) # function you defined above"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**they are slightly different! but very close**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now, we calculate distance matrix using Scipy package. 'euclidean' method is used among [many other methods](https://docs.scipy.org/doc/scipy/reference/generated/scipy.spatial.distance.cdist.html) provided. You already tested the Euclidean method using lat and long, the answer seems quite wrong. (because it is in degree, not in meters). \n",
    "\n",
    "To get the distance in meters, we will first project the data into two-dimensional space (from long, lat to x, y). Then calculate the distance matrix."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pyproj\n",
    "from scipy.spatial import distance"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(site_df.???) # what is the projected and geographic coordinate system?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define some common projections using EPSG codes\n",
    "wgs84=pyproj.Proj(\"+init=EPSG:4326\") # LatLon with WGS84 datum used by GPS units and Google Earth\n",
    "osgb36=pyproj.Proj(\"+init=EPSG:27700\") # UK Ordnance Survey, 1936 datum - this is the one used in the UK"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "site_df['x'], site_df['y'] = pyproj.transform(wgs84, osgb36, list(site_df['longitude']), list(site_df['latitude']))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "coords = site_df[['x','y']] "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "coords.???  #check it"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "distance.cdist(coords, coords, 'euclidean')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Verifying it using measure tools in QGIS!!!** The shapefile you created can be opened there."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Credits"
   ]
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "Geocomp",
   "language": "python",
   "name": "gsa2018"
  },
  "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.6.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}