{
  "cells": [
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "%matplotlib inline"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "\n# Plotting in Different Projections\n\nThis example shows how to overlay data and graphics in different projections,\ndemonstrating various features of Iris, Cartopy and matplotlib.\n\nWe wish to overlay two datasets, defined on different rotated-pole grids.\nTo display both together, we make a pseudocoloured plot of the first, overlaid\nwith contour lines from the second.\nWe also add some lines and text annotations drawn in various projections.\n\nWe plot these over a specified region, in two different map projections.\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "import cartopy.crs as ccrs\nimport matplotlib.pyplot as plt\nimport numpy as np\n\nimport iris\nimport iris.plot as iplt\n\n# Define a Cartopy 'ordinary' lat-lon coordinate reference system.\ncrs_latlon = ccrs.PlateCarree()\n\n\ndef make_plot(projection_name, projection_crs):\n\n    # Create a matplotlib Figure.\n    plt.figure()\n\n    # Add a matplotlib Axes, specifying the required display projection.\n    # NOTE: specifying 'projection' (a \"cartopy.crs.Projection\") makes the\n    # resulting Axes a \"cartopy.mpl.geoaxes.GeoAxes\", which supports plotting\n    # in different coordinate systems.\n    ax = plt.axes(projection=projection_crs)\n\n    # Set display limits to include a set region of latitude * longitude.\n    # (Note: Cartopy-specific).\n    ax.set_extent((-80.0, 20.0, 10.0, 80.0), crs=crs_latlon)\n\n    # Add coastlines and meridians/parallels (Cartopy-specific).\n    ax.coastlines(linewidth=0.75, color=\"navy\")\n    ax.gridlines(crs=crs_latlon, linestyle=\"-\")\n\n    # Plot the first dataset as a pseudocolour filled plot.\n    maindata_filepath = iris.sample_data_path(\"rotated_pole.nc\")\n    main_data = iris.load_cube(maindata_filepath)\n    # NOTE: iplt.pcolormesh calls \"pyplot.pcolormesh\", passing in a coordinate\n    # system with the 'transform' keyword:  This enables the Axes (a cartopy\n    # GeoAxes) to reproject the plot into the display projection.\n    iplt.pcolormesh(main_data, cmap=\"RdBu_r\")\n\n    # Overplot the other dataset (which has a different grid), as contours.\n    overlay_filepath = iris.sample_data_path(\"space_weather.nc\")\n    overlay_data = iris.load_cube(overlay_filepath, \"total electron content\")\n    # NOTE: as above, \"iris.plot.contour\" calls \"pyplot.contour\" with a\n    # 'transform' keyword, enabling Cartopy reprojection.\n    iplt.contour(\n        overlay_data, 20, linewidths=2.0, colors=\"darkgreen\", linestyles=\"-\"\n    )\n\n    # Draw a high resolution margin line, inset from the pcolormesh border.\n    # First calculate rectangle corners, 7% in from each corner of the data.\n    x_coord, y_coord = main_data.coord(axis=\"x\"), main_data.coord(axis=\"y\")\n    x_start, x_end = np.min(x_coord.points), np.max(x_coord.points)\n    y_start, y_end = np.min(y_coord.points), np.max(y_coord.points)\n    margin = 0.07\n    margin_fractions = np.array([margin, 1.0 - margin])\n    x_lower, x_upper = x_start + (x_end - x_start) * margin_fractions\n    y_lower, y_upper = y_start + (y_end - y_start) * margin_fractions\n    steps = np.linspace(0, 1)\n    zeros, ones = np.zeros(steps.size), np.ones(steps.size)\n    x_delta, y_delta = (x_upper - x_lower), (y_upper - y_lower)\n    x_points = x_lower + x_delta * np.concatenate(\n        (steps, ones, steps[::-1], zeros)\n    )\n    y_points = y_lower + y_delta * np.concatenate(\n        (zeros, steps, ones, steps[::-1])\n    )\n    # Get the Iris coordinate sytem of the X coordinate (Y should be the same).\n    cs_data1 = x_coord.coord_system\n    # Construct an equivalent Cartopy coordinate reference system (\"crs\").\n    crs_data1 = cs_data1.as_cartopy_crs()\n    # Draw the rectangle in this crs, with matplotlib \"pyplot.plot\".\n    # NOTE: the 'transform' keyword specifies a non-display coordinate system\n    # for the plot points (as used by the \"iris.plot\" functions).\n    plt.plot(\n        x_points,\n        y_points,\n        transform=crs_data1,\n        linewidth=2.0,\n        color=\"white\",\n        linestyle=\"--\",\n    )\n\n    # Mark some particular places with a small circle and a name label...\n    # Define some test points with latitude and longitude coordinates.\n    city_data = [\n        (\"London\", 51.5072, 0.1275),\n        (\"Halifax, NS\", 44.67, -63.61),\n        (\"Reykjavik\", 64.1333, -21.9333),\n    ]\n    # Place a single marker point and a text annotation at each place.\n    for name, lat, lon in city_data:\n        plt.plot(\n            lon,\n            lat,\n            marker=\"o\",\n            markersize=7.0,\n            markeredgewidth=2.5,\n            markerfacecolor=\"black\",\n            markeredgecolor=\"white\",\n            transform=crs_latlon,\n        )\n        # NOTE: the \"plt.annotate call\" does not have a \"transform=\" keyword,\n        # so for this one we transform the coordinates with a Cartopy call.\n        at_x, at_y = ax.projection.transform_point(\n            lon, lat, src_crs=crs_latlon\n        )\n        plt.annotate(\n            name,\n            xy=(at_x, at_y),\n            xytext=(30, 20),\n            textcoords=\"offset points\",\n            color=\"black\",\n            backgroundcolor=\"white\",\n            size=\"large\",\n            arrowprops=dict(arrowstyle=\"->\", color=\"white\", linewidth=2.5),\n        )\n\n    # Add a title, and display.\n    plt.title(\n        \"A pseudocolour plot on the {} projection,\\n\"\n        \"with overlaid contours.\".format(projection_name)\n    )\n    iplt.show()\n\n\ndef main():\n    # Demonstrate with two different display projections.\n    make_plot(\"Equidistant Cylindrical\", ccrs.PlateCarree())\n    make_plot(\"North Polar Stereographic\", ccrs.NorthPolarStereo())\n\n\nif __name__ == \"__main__\":\n    main()"
      ]
    }
  ],
  "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.10.6"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}