{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Introduction to GeoPandas\n",
    "\n",
    "This quick tutorial provides an introduction to the key concepts of GeoPandas. In a few minutes, we'll describe the basics which allow you to start your projects.\n",
    "\n",
    "## Concepts\n",
    "\n",
    "GeoPandas, as the name suggests, extends popular data science library [pandas](https://pandas.pydata.org) by adding support for geospatial data. If you are not familiar with `pandas`, we recommend taking a quick look at its [Getting started documentation](https://pandas.pydata.org/docs/getting_started/index.html#getting-started) before proceeding.\n",
    "\n",
    "The core data structure in GeoPandas is `geopandas.GeoDataFrame`, a subclass of `pandas.DataFrame` able to store geometry columns and perform spatial operations. Geometries are handled by `geopandas.GeoSeries`, a subclass of `pandas.Series`. Therefore, your `GeoDataFrame` is a combination of `Series` with your data (numerical, boolean, text etc.) and `GeoSeries` with geometries (points, polygons etc.). You can have as many columns with geometries as you wish, there's no limit typical for desktop GIS software.\n",
    "\n",
    "![geodataframe schema](../_static/dataframe.svg)\n",
    "\n",
    "Each `GeoSeries` can contain any geometry type (we can even mix them within a single array) and has a `GeoSeries.crs` attribute, which stores information on the projection (CRS stands for Coordinate Reference System). Therefore, each `GeoSeries` in a `GeoDataFrame` can be in a different projection, allowing you to have, for example, multiple versions of the same geometry, just in a different CRS.\n",
    "\n",
    "One `GeoSeries` within a `GeoDataFrame` is seen as the _active_ geometry, which means that all geometric operations applied to a `GeoDataFrame` use the specified column.\n",
    "\n",
    "\n",
    "<div class=\"alert alert-info\">\n",
    "User Guide\n",
    "    \n",
    "See more on [data structures in the User Guide](../docs/user_guide/data_structures.rst).\n",
    "</div>\n",
    "\n",
    "\n",
    "Let's see how this works in practice.\n",
    "\n",
    "## Reading and writing files\n",
    "\n",
    "First, we need to read some data.\n",
    "\n",
    "### Read files\n",
    "\n",
    "Assuming we have a file containing both data and geometry (e.g. GeoPackage, GeoJSON, Shapefile), we can easily read it using `geopandas.read_file` function, which automatically detects filetype and creates a `GeoDataFrame`. In this example, we'll use the `\"nybb\"` dataset, a map of New York boroughs which is part of GeoPandas installation. Therefore we need to get the path to the actual file. With your file, you specify a path as a string (`\"my_data/my_file.geojson\"`)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import geopandas\n",
    "\n",
    "path_to_data = geopandas.datasets.get_path(\"nybb\")\n",
    "gdf = geopandas.read_file(path_to_data)\n",
    "\n",
    "gdf"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Write files\n",
    "\n",
    "Writing a `GeoDataFrame` back to file is similarly simple, using `GeoDataFrame.to_file`. The default file format is Shapefile, but you can specify your own using `driver` keyword."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdf.to_file(\"my_file.geojson\", driver=\"GeoJSON\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-info\">\n",
    "User Guide\n",
    "    \n",
    "See more on [reading and writing data in the User Guide](../docs/user_guide/io.rst).\n",
    "</div>\n",
    "\n",
    "\n",
    "\n",
    "## Simple methods\n",
    "\n",
    "Now we have our `GeoDataFrame` and can start working with its geometry. \n",
    "\n",
    "Since we have only one geometry column read from the file, it is automatically seen as the active geometry and methods used on `GeoDataFrame` will be applied to the `\"geometry\"` column.\n",
    "\n",
    "### Measuring area\n",
    "\n",
    "To measure the area of each polygon (or MultiPolygon in this specific case), we can use `GeoDataFrame.area` attribute, which returns a `pandas.Series`. Note that `GeoDataFrame.area` is just `GeoSeries.area` applied to an active geometry column.\n",
    "\n",
    "But first, we set the names of boroughs as an index, to make the results easier to read."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdf = gdf.set_index(\"BoroName\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdf[\"area\"] = gdf.area\n",
    "gdf[\"area\"]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Getting polygon boundary and centroid\n",
    "\n",
    "To get just the boundary of each polygon (LineString), we can call `GeoDataFrame.boundary`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdf['boundary'] = gdf.boundary\n",
    "gdf['boundary']"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since we have saved boundary as a new column, we now have two geometry columns in the same `GeoDataFrame`.\n",
    "\n",
    "We can also create new geometries, which could be, for example, a buffered version of the original one (i.e., `GeoDataFrame.buffer(10)`) or its centroid:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdf['centroid'] = gdf.centroid\n",
    "gdf['centroid']"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Measuring distance\n",
    "\n",
    "We can also measure how far is each centroid from the first one."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "first_point = gdf['centroid'].iloc[0]\n",
    "gdf['distance'] = gdf['centroid'].distance(first_point)\n",
    "gdf['distance']"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "It's still a DataFrame, so we have all the pandas functionality available to use on the geospatial dataset, and to do data manipulations with the attributes and geometry information together.\n",
    "\n",
    "For example, we can calculate average of the distance measured above (by accessing the `'distance'` column, and calling the `mean()` method on it):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdf['distance'].mean()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Making maps\n",
    "\n",
    "GeoPandas can also plot maps, so we can check how our geometries look like in space. The key method here is `GeoDataFrame.plot()`. In the example below, we plot the `\"area\"` we measured earlier using the active geometry column. We also want to show a legend (`legend=True`)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdf.plot(\"area\", legend=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Switching the active geometry (`GeoDataFrame.set_geometry`) to centroids, we can plot the same data using point geometry."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdf = gdf.set_geometry(\"centroid\")\n",
    "gdf.plot(\"area\", legend=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And we can also layer both `GeoSeries` on top of each other. We just need to use one plot as an axis for the other."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ax = gdf[\"geometry\"].plot()\n",
    "gdf[\"centroid\"].plot(ax=ax, color=\"black\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we set the active geometry back to the original `GeoSeries`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdf = gdf.set_geometry(\"geometry\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-info\">\n",
    "User Guide\n",
    "    \n",
    "See more on [mapping in the User Guide](../docs/user_guide/mapping.rst).\n",
    "</div>\n",
    "\n",
    "## Geometry creation\n",
    "\n",
    "We can further work with the geometry and create new shapes based on those we already have. \n",
    "\n",
    "### Convex hull\n",
    "\n",
    "If we are interested in the convex hull of our polygons, we can call `GeoDataFrame.convex_hull`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdf[\"convex_hull\"] = gdf.convex_hull"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ax = gdf[\"convex_hull\"].plot(alpha=.5)  # saving the first plot as an axis and setting alpha (transparency) to 0.5\n",
    "gdf[\"boundary\"].plot(ax=ax, color=\"white\", linewidth=.5)  # passing the first plot and setting linewitdth to 0.5"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Buffer\n",
    "\n",
    "In other cases, we may need to buffer the geometry using `GeoDataFrame.buffer()`. Geometry methods are automatically applied to the active geometry, but we can apply them directly to any `GeoSeries` as well. Let's buffer the boroughs and their centroids and plot both on top of each other."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# buffering the active geometry by 10 000 feet (geometry is already in feet)\n",
    "gdf[\"buffered\"] = gdf.buffer(10000)\n",
    "\n",
    "# buffering the centroid geometry by 10 000 feet (geometry is already in feet)\n",
    "gdf[\"buffered_centroid\"] = gdf[\"centroid\"].buffer(10000)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ax = gdf[\"buffered\"].plot(alpha=.5)  # saving the first plot as an axis and setting alpha (transparency) to 0.5\n",
    "gdf[\"buffered_centroid\"].plot(ax=ax, color=\"red\", alpha=.5)  # passing the first plot as an axis to the second\n",
    "gdf[\"boundary\"].plot(ax=ax, color=\"white\", linewidth=.5)  # passing the first plot and setting linewitdth to 0.5"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div class=\"alert alert-info\">\n",
    "User Guide\n",
    "    \n",
    "See more on [geometry creation and manipulation in the User Guide](../docs/user_guide/geometric_manipulations.rst).\n",
    "</div>\n",
    "\n",
    "## Geometry relations\n",
    "\n",
    "We can also ask about the spatial relations of different geometries. Using the geometries above, we can check which of the buffered boroughs intersect the original geometry of Brooklyn, i.e., is within 10 000 feet from Brooklyn.\n",
    "\n",
    "First, we get a polygon of Brooklyn."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "brooklyn = gdf.loc[\"Brooklyn\", \"geometry\"]\n",
    "brooklyn"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The polygon is a [shapely geometry object](https://shapely.readthedocs.io/en/stable/manual.html#geometric-objects), as any other geometry used in GeoPandas."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "type(brooklyn)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then we can check which of the geometries in `gdf[\"buffered\"]` intersects it."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdf[\"buffered\"].intersects(brooklyn)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Only Bronx (on the north) is more than 10 000 feet away from Brooklyn. All the others are closer and intersect our polygon.\n",
    "\n",
    "Alternatively, we can check which buffered centroids are entirely within the original boroughs polygons. In this case, both `GeoSeries` are aligned, and the check is performed for each row."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdf[\"within\"] = gdf[\"buffered_centroid\"].within(gdf)\n",
    "gdf[\"within\"]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can plot the results on the map to confirm the finding."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdf = gdf.set_geometry(\"buffered_centroid\")\n",
    "ax = gdf.plot(\"within\", legend=True, categorical=True, legend_kwds={'loc': \"upper left\"})  # using categorical plot and setting the position of the legend\n",
    "gdf[\"boundary\"].plot(ax=ax, color=\"black\", linewidth=.5)  # passing the first plot and setting linewitdth to 0.5"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Projections\n",
    "\n",
    "Each `GeoSeries` has the Coordinate Reference System (CRS) accessible as `GeoSeries.crs`. CRS tells GeoPandas where the coordinates of geometries are located on the Earth. In some cases, CRS is geographic, which means that coordinates are in latitude and longitude. In those cases, its CRS is WGS84, with the authority code `EPSG:4326`. Let's see the projection of our NY boroughs `GeoDataFrame`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdf.crs"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Geometries are in `EPSG:2263` with coordinates in feet. We can easily re-project a `GeoSeries` to another CRS, like `EPSG:4326` using `GeoSeries.to_crs()`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdf = gdf.set_geometry(\"geometry\")\n",
    "boroughs_4326 = gdf.to_crs(\"EPSG:4326\")\n",
    "boroughs_4326.plot()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "boroughs_4326.crs"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Notice the difference in coordinates along the axes of the plot. Where we had 120 000 - 280 000 (feet) before, we have 40.5 - 40.9 (degrees) now. In this case, `boroughs_4326` has a `\"geometry\"` column in WGS84 but all the other (with centroids etc.) remains in the original CRS.\n",
    "\n",
    "<div class=\"alert alert-warning\">\n",
    "Warning\n",
    "    \n",
    "For operations that rely on distance or area, you always need to use projected CRS (in meters, feet, kilometers etc.) not a geographic one. GeoPandas operations are planar, and degrees reflect the position on a sphere. Therefore the results may not be correct. For example, the result of `gdf.area.sum()` (projected CRS) is 8 429 911 572 ft<sup>2</sup> but the result of `boroughs_4326.area.sum()` (geographic CRS) is 0.083.\n",
    "</div>\n",
    "\n",
    "<div class=\"alert alert-info\">\n",
    "User Guide\n",
    "    \n",
    "See more on [projections in the User Guide](../docs/user_guide/projections.rst).\n",
    "</div>\n",
    "\n",
    "## What next?\n",
    "\n",
    "With GeoPandas we can do much more that this, from [aggregations](../docs/user_guide/aggregation_with_dissolve.rst), to [spatial joins](../docs/user_guide/mergingdata.rst), [geocoding](../docs/user_guide/geocoding.rst) and [much more](../gallery/index.rst).\n",
    "\n",
    "Head to the [User Guide](../docs/user_guide.rst) for to learn more about different functionality of GeoPandas, to the [Examples](../gallery/index.rst) to see how it can be used or the the [API reference](../docs/reference.rst) for the details."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "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.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}