scitools / cartopy Goto Github PK

View Code? Open in Web Editor NEW

1.4K 55.0 360.0 82.47 MB

Cartopy - a cartographic python library with matplotlib support

Home Page: https://scitools.org.uk/cartopy/docs/latest

License: BSD 3-Clause "New" or "Revised" License

Python 97.70% Cython 2.30%

cartopy matplotlib python geometry maps spatial projections

cartopy's Introduction

Cartopy is a Python package designed to make drawing maps for data analysis and visualisation easy.

Overview
Get in touch
License and copyright

Overview

Cartopy is a Python package designed to make drawing maps for data analysis and visualisation easy.

It features:

object oriented projection definitions
point, line, polygon and image transformations between projections
integration to expose advanced mapping in Matplotlib with a simple and intuitive interface
powerful vector data handling by integrating shapefile reading with Shapely capabilities

Documentation can be found at https://scitools.org.uk/cartopy/docs/latest/.

Get in touch

Ask usage questions on StackOverflow.
For less well defined questions, ideas, general discussion or announcements of related projects use the Cartopy category on Matplotlib's Discourse.
Report bugs, suggest features or view the source code on GitHub.
To chat with developers and other users you can use the Gitter Chat.

Credits, copyright and license

Cartopy is developed collaboratively under the SciTools umberella.

A full list of codecontributors ("Cartopy contributors") can be found at https://github.com/SciTools/cartopy/graphs/contributors.

Code is just one of many ways of positively contributing to Cartopy, please see our contributing guide for more details on how you can get involved.

Cartopy is released under the 3-Clause BSD license with a shared copyright model. See LICENSE for full terms.

The Met Office has made a significant contribution to the development, maintenance and support of this library. All Met Office contributions are copyright on behalf of the British Crown.

cartopy's People

Contributors

Stargazers

Watchers

Forkers

pelson bblay esc24 rhattersley dmcdougall bjlittle pp-mo marqh cpelley ajdawson rockdoc dpeterk munslowa sgillies susanjl-zz ahill818 octeufer mikewalker92 rachelnorth 5n1p imclab seilis sigmamonster krischer zak-k testvidya11 rsignell-usgs fandres70 njwilson23 perrette lbdreyer qulogic om-henners micolous snorfalorpagus nltyrrell philip-riggs scollis grabkiak ckmo expertanalytics kesariverma odermatt gepcel dreamalligator akrherz williemaddox ocefpaf ozayn jeitan josephhogg scmc72 gavin971 timmie dvalters duncanwp ztessler rabernat djkirkham jonathf philippjfr ljwolf corinnebosley belteshassar ltiao volcanohawk lidi100 siyh ahaberlie cequencer dopplershift bssrdf solter skjaeve aeroevan dennissergeev melissaebrooks swatchai 8828 olivier-membrive heyifei cvelascof russsty xigrug mhearne-usgs rouvenb matthewstidham lpk-py choldgraf xclu iamasam has2k1 ghisvail jpcanright r5r3 lukec92 aluhamaa kaedonkers stefmiles huayue21

cartopy's Issues

Write about the image interpolation used in cartopy

Discuss:

The process
Possible alternatives (http://www.ncl.ucar.edu/Document/Functions/Built-in/natgrid.shtml?)
How users can make use of this for data gridding and warping.

Fix external boundary calculation

The external boundary calculation is broken since a major re-write of the path transformation code.

import matplotlib.pyplot as plt
import cartopy.crs as ccrs

rob = ccrs.Robinson()
ax = plt.axes(projection=rob)
ax.set_global()
print ax.boundary_poly()
print ax.ll_boundary_poly()
plt.show()

Investigate and fix, and ensure that this resolves the gshhs dataset plotting issue (#3).

Ability to do a repr on a CRS.

Currently there is not enought infrastructure in place to be able to do a code generating repr on a CRS. This would be a very useful feature and is something that we should certainly strive for.

Low-level handling of degenerate polygons

Projection._project_polygon() is rather unhelpfully converting degenerate polygons (such as those produced by matplotlib.contourf) into global extent.
e.g.
Polygon([(20, 20), (20, 20), (20, 20), (20, 20)])
Polygon([(20, 20), (20, 20), (20, 20.1), (20, 20)])

Whilst it's easy to filter these out in the MPL integration layer, ideally the low level re-projection code should be able to deal with them.

Ingest GeoTiff files

Currently, the only tool I have seen to read the full header of a geotiff file in python is GDAL. Implement a simple interface to load GeoTiffs easily.

axes range bug

The code below sets the dataLim as expected.
However, when using the PlateCarree object instead (commented out here),
the upper limits remain at 1,1.

import matplotlib.pyplot as plt
import cartopy.crs as ccrs

#ax = plt.axes(projection=ccrs.PlateCarree())
ax = plt.axes()

print ax.dataLim

ax.update_datalim([(-0.42041015714311403, -0.12780000269415132), (-0.33129882901803986, -0.0386999994516301)])

print ax.dataLim

pcolormesh can only be drawn in native coordinates

more specifically, in coordinates that don't wrap.

Consider whether this can be worked around by drawing more than one pcolormesh.

Non native extents of a plot

Allow users to define the region of their plot in non-native coordinates.

Currently, users must use the standard interface of ax.set_xlim/ax.set_ylim with native values. Accept bounds in another coordinate system, producing a polygon, warping the polygon to the native projection, and then subsequently applying the bounds of the polygon as the native extents.

Use of image_comparison decorator in v0.4.x causes test error.

The use of the image_comparison decorator in test_mpl_integration.py causes:

======================================================================
ERROR: Failure: NameError (name 'image_comparison' is not defined)
----------------------------------------------------------------------
Traceback (most recent call last):
  File "/usr/local/sci/lib/python2.7/site-packages/nose-1.1.2-py2.7.egg/nose/loader.py", line 390, in loadTestsFromName
    addr.filename, addr.module)
  File "/usr/local/sci/lib/python2.7/site-packages/nose-1.1.2-py2.7.egg/nose/importer.py", line 39, in importFromPath
    return self.importFromDir(dir_path, fqname)
  File "/usr/local/sci/lib/python2.7/site-packages/nose-1.1.2-py2.7.egg/nose/importer.py", line 86, in importFromDir
    mod = load_module(part_fqname, fh, filename, desc)
  File "/net/home/h06/ecamp/cartopy/lib/cartopy/tests/mpl/test_mpl_integration.py", line 189, in <module>
    @image_comparison(baseline_images=['pcolormesh_global_wrap1'])
NameError: name 'image_comparison' is not defined

as the decorator is not imported. I suggest modifying test_mpl_integration.py to use the ImageTesting decorator.

In addition, the image_comparison decorator is imported but not used in the following:

test_image_tiles.py
test_img_nest.py
test_shapely_to_mpl.py
test_caching.py
test_gridliner.py

Global extent failure

ax = plt.axes(projection=PlateCarree())

# Works OK
ax.set_extent((-40, 40, -30, 30))

# Crashes
ax.set_extent((-180, 180, -90, 90))
...
File ".../cartopy/mpl_integration/geoaxes.py", line 244, in set_extent
    x1, y1, x2, y2 = r.bounds
ValueError: need more than 0 values to unpack

Tile merging misses interiors

To reproduce:


import matplotlib.pyplot as plt

import cartopy.crs as ccrs

ax = plt.axes(projection=ccrs.Mercator())
ax.set_extent([-80, 80, 0, 60], ccrs.PlateCarree())
import cartopy.io.img_tiles as cit
ax.add_image(cit.MapQuestOSM(), 3)

ax.coastlines()

plt.show()

At the same time, try investigating the failing (probably all related...):

import matplotlib.pyplot as plt

import cartopy.crs as ccrs

ax = plt.axes(projection=ccrs.Mercator())
ax.set_extent([3e5, 3.7e5, 1.44e5, 2.14e5], ccrs.OSGB())
import cartopy.io.img_tiles as cit
ax.add_image(cit.MapQuestOSM(), 5)

ax.coastlines(resolution='50m')

plt.show()

and

import matplotlib.pyplot as plt
import cartopy.crs as ccrs

ax = plt.axes(projection=ccrs.Mercator())
ax.set_global()
import cartopy.io.img_tiles as cit
ax.add_image(cit.MapQuestOSM(), 1)
ax.coastlines()
plt.show()

Contour labelling

Implement contour labelling on the mpl interface. If it all works, provide an example of doing contour labelling.

Travis testing

Implement the appropriate hooks to get travis-ci running for cartopy.

Patch collections are transformed twice

    ax = plt.axes(projection=ccrs.PlateCarree())
    ax.set_global()

    import matplotlib.path as mpath 
    pth = mpath.Path([[0, 45], [60, 45], [60, -45], [0, -45], [0, -45], 
                [10, 20], [10, -20], [40, -20], [40, 20], [10, -20]], 
                [1, 2, 2, 2, 79, 1, 2, 2 , 2, 79])

    import matplotlib.collections as mcollections
    import matplotlib.patches as mpatches

#    ax.add_patch(mpatches.PathPatch(pth, transform=ccrs.PlateCarree()))
    ax.add_collection(mcollections.PathCollection([pth], transform=ccrs.PlateCarree()))

    plt.show()

The transform_path_non_affine method is called twice for a simple show, where a Patch only calls it once.

This might be a matplotlib fix, or may be doable in cartopy.

This finding came from a message from @rhattersley that collection drawing was slow - so there is some hope that this can be sped up...

Allow user to control transform resolution

Currently, the great circle transform resolution has been set to an arbitrary value. Tweak the value to some heuristic, and make it possible for users to define the value themselves.

Vector handling

Consider the ability to plot vector fields in Cartopy. The initial usecase should:

Take data on a rotated pole grid (for example Met Office NAE output) of wind data (u and v component) which is co-located and where the u and v component are in the direction of the grid rather than east/north.
Plot the data on its native grid (no transformation necessary)
Plot the data on a PlateCarree map (the vector magnitude is preserved, but its direction needs to be transformed).
Plot the data on an arbitrary map (may be no more difficult than the PlateCarree case)

To others: Please feel free to update this issue with clarifications/questions.

Implement visual tests

Implement visual unit tests which follow a similar pattern to those in the matplotlib project.

coastlines_land doesn't play with Orthographic

>>> import cartopy.crs as ccrs
>>> import matplotlib.pyplot as plt
>>> ax = plt.axes(projection=ccrs.Orthographic())
>>> ax.coastlines_land()
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File ".../cartopy/lib/cartopy/mpl_integration/geoaxes.py", line 177, in coastlines_land
    paths.extend(patch.geos_to_path(self.projection.project_geometry(geom)))
  File ".../cartopy/lib/cartopy/crs.py", line 141, in project_geometry
    return getattr(self, method_name)(geometry, src_crs)
  File ".../cartopy/lib/cartopy/crs.py", line 202, in _project_multipolygon
    return sgeom.MultiPolygon(geoms)
  File ".../lib/python2.7/site-packages/Shapely-1.2.14-py2.7.egg/shapely/geometry/multipolygon.py", line 57, in __init__
    self._geom, self._ndim = geos_multipolygon_from_polygons(polygons)
  File ".../lib/python2.7/site-packages/Shapely-1.2.14-py2.7.egg/shapely/geometry/multipolygon.py", line 127, in geos_multipolygon_from_polygons
    exemplar = obs[0]
IndexError: list index out of range

Stuck polygon transform

from matplotlib.path import Path
import numpy

import cartopy.mpl_integration.patch as cpatch
from cartopy.crs import Geodetic, PlateCarree

x = [260.625, 360.0, 360.0, 0.0, 210.0, 135.0, 260.625]
y = [68.90383337092122, 88.790680473372788, 90.0, 88.790680473372788,
     77.768481754584982, 79.855609199690093, 68.90383337092122]
codes = [1, 2, 2, 2, 2, 2, 2]

vertices = numpy.vstack([numpy.array(x), numpy.array(y)]).T

target_projection = PlateCarree()
source_crs = Geodetic()

path = Path(vertices, codes)

for geom in cpatch.path_to_geos(path):
    print 'Reprojecting...'
    print target_projection.project_geometry(geom, source_crs)

Implement/source a geodetic contouring algorithm

Currently mpl cartopy contouring is very 2d. Find or write a 3d geodetic contour algorithm.

Another stuck polygon

import cartopy.crs as ccrs
import matplotlib.pyplot as plt

p = ccrs.RotatedPole(pole_latitude=37.5, pole_longitude=357.5-180)
ax = plt.axes(projection=p)
ax.natural_earth_shp('ocean',
                         resolution='110m',
                         facecolor=(192/256., 232/256., 254/256.),
                         edgecolor='none')

ax.coastlines()
plt.show()

Takes an age (I've never seen it finish).

Probably another case of #77 and #60...

Improve the examples bundled with cartopy

For 0.4 the examples which came with cartopy were a mish-mash of code, most of which were badly coded and didn't work for that release (I wrote them myself so am comfortable saying this).

Format and make good the current suite of examples, and add each to the test_examples test to maintain the quality of the ouput in further releases.

Cursor is not displaying the correct coordinate values

The GenericProjectionAxes.format_coord method is not producing the correct results. To reproduce:

import matplotlib.pyplot as plt
import cartopy.crs as ccrs

rob = ccrs.Robinson()
ax = plt.axes(projection=rob)
ax.coastlines()
ax.set_global()
plt.show()

Placing the mouse on the map should show the correct coordinates in the status bar.

Provide geodetic polygon "contains" functionality

It is a common requirement to test containment of two polygons on the globe.

Provide geodetic polygon "contains" functionality. Make it possible to extend this to other geodetic polygon operations in the future.

Infinite (project?) loop

The following code never terminates:

import matplotlib.pyplot as plt
import numpy as np

import cartopy.crs as ccrs


def sample_data(shape=(73, 145)):
    """Returns lons, lats and data of some interesting data on a regular grid."""
    nlats, nlons = shape 
    lats = np.linspace(-np.pi/2, np.pi/2, nlats)
    lons = np.linspace(0, 2*np.pi, nlons)
    lons, lats = np.meshgrid(lons, lats)
    wave = 0.75*(np.sin(2*lats)**8) * np.cos(4*lons)
    mean = 0.5*np.cos(2*lats) * ((np.sin(2*lats))**2 + 2)

    lats = np.rad2deg(lats)
    lons = np.rad2deg(lons)
    data = wave + mean

    return lons, lats, data

def main():
    ax = plt.axes(projection=ccrs.NorthPolarStereo())

    ax.contourf(*sample_data(), nlev=2, transform=ccrs.PlateCarree())
    ax.coastlines()
    plt.show()


if __name__ == '__main__':
    main()

Tiled image accumulation and plotting

Implement a mechanism to retrieve and join tiled images (such as those from a web tile service (WTS)), and return a single image for inclusion on the geoaxes.

Consider reversing the tranform API

Currently:

target_projection.transform_points(source_projection, x, y)

Consider changing this to:

source_projection.transform_points(target_projection, x, y)

custom limits

I believe we'll need to specify custom limits for plotting osgb, rotated and regular lat lon data.
Is this possible? Perhaps this information is provided alongside the projection, rather than inside it?

Missing big polygons

Once a polygon is large enough to entirely contain a map domain (not just the extents currently visible) then it disappears. ~~This can only be a problem for non-global maps, such as OSGB.~~

import matplotlib.pyplot as plt
import shapely.geometry as sgeom

import cartopy.crs as ccrs

pc = ccrs.PlateCarree()
osgb = ccrs.OSGB()

small = [sgeom.box(-3, 51, 1, 55)]
big = [sgeom.box(-20, 10, 30, 80)]

ax = plt.subplot(2, 2, 1, projection=pc)
ax.add_geometries(small, ccrs.PlateCarree())
ax.scatter(0, 0, c='g', marker='+', transform=osgb)
ax.set_global()

ax = plt.subplot(2, 2, 2, projection=osgb)
ax.add_geometries(small, ccrs.PlateCarree())

ax = plt.subplot(2, 2, 3, projection=pc)
ax.add_geometries(big, ccrs.PlateCarree())
ax.scatter(0, 0, c='g', marker='+', transform=osgb)
ax.set_global()

ax = plt.subplot(2, 2, 4, projection=osgb)
ax.add_geometries(big, ccrs.PlateCarree())

plt.show()

Be cleverer with pcolormesh

Plotting data with pcolormesh can have really nasty effects if not plotting in the native coordinate system (180 degree shifts included). Consider what can be done to improve this situation. Be it rolling the data, issuing a warning, or plotting the pcolormesh more than once to get the expected overall effect.

Handle interior polygon transforming

cartopy failed to install on OpenSUSE 12.2

Received a question regarding installing cartopy on OpenSUSE 12.2:

trying to install this on my home PC (OpenSUSE 12.2), including all of the stated dependencies.... 
But I still get exceptions running any example code relating to the projections. 
The ccrs.<projection> method returns an object, but matplotlib knows nothing
about these projections.

Implement gridlines/axis artists for arbitrary projections

It would be nice to unify the interface to draw graticules for all projections under the same API.
e.g. make it possible to use the same interface to add autoscaling axis artists for both osgb and plate carree (extension: both on the same plot).

This is an extension of #5.

Convert the deleted main tests from #69 into unit tests

Original PR: #69

Implement axis/gridlines which adjust according to zoom level/pan position

Currently, the only means of adding graticules is to give a specific degree interval. Implement latitude longitude graticules which auto adjust with mpl integration.

Provide data attribution capabilities.

Many of the datasets used by cartopy require attribution. Provide an interface which makes it easy to attribute a dataset which is automatically displayed at the bottom of a map.

Gshhs plotting is broken

The primary cause is that the external boundary polygon is being incorrectly calculated.

Expand the install instructions

Cartopy installation is quite simple, but it would be nice to have clear instructions for fresh Ubuntu 12.04 and OSX 10.8 systems.

As part of this work, try deploying cartopy to PyPI for easy pip installation.

Writeup Ubuntu 12.04 install instructions

They are:

wget http://python-distribute.org/distribute_setup.py
sudo python distribute_setup.py
wget https://raw.github.com/pypa/pip/master/contrib/get-pip.py
sudo python get-pip.py

# to install from git source
sudo apt-get update
sudo apt-get install git

sudo apt-get install python-dev
sudo pip install numpy

sudo pip install cython
sudo apt-get install libproj-dev

# install geos 
sudo apt-add-repository ppa:sharpie/for-science
sudo apt-get update
sudo apt-get install libgeos-dev

sudo pip install shapely
sudo apt-get install python-scipy
sudo pip install pyshp

sudo apt-get install python-imaging

# for matplotlib 1.2 (no pre-build package available yet)
sudo apt-get build-dep matplotlib
git clone https://github.com/matplotlib/matplotlib.git
cd matplotlib
git checkout v1.2.0rc2
sudo python setup.py install

Optimise plot size

Because cartopy uses a fixed aspect ratio in mpl, the actual size of the axes in pixels depends on the viewLim values. Instead, make it possible to define the ideal size of the plot, constrained by just one of the bounds in viewLim.

This will involve close mpl work.

Negative X bug

Uncomment the x flip and the code fails.

import iris

import matplotlib
matplotlib.use('tkagg')
import matplotlib.pyplot as plt

import iris.tests as tests
import iris.plot as iplt


cube = iris.load(tests.get_data_path(('GRIB', 'ij_directions', "ipos_jpos.grib2")))[0]
#cube = cube[:,::-1]

iplt.contourf(cube)
plt.gca().coastlines()
plt.title("ipos_jpos cube")
iplt.show()

clockwise/anti-clockwise polygons

In order to make a polygon unambiguous on a sphere, cartopy must take the direction (clockwise or anti-clockwise) of a polygon into account. Currently that direction is fixed. Make it easy for a user to change that on a polygon by polygon basis.

mpl interface: make it possible to add filled polygons (coastlines)

Writeup 0.4 install on OSX 10.8 (mountain lion) using homebrew

ruby -e "$(curl -fsSkL raw.github.com/mxcl/homebrew/go)"
brew doctor
brew install python
export PATH=/usr/local/bin:$PATH
curl -O https://raw.github.com/pypa/virtualenv/master/virtualenv.py
python virtualenv.py cartopy-python


brew install gfortran
cartopy-python/bin/pip install -e git+https://github.com/scipy/scipy#egg=scipy-dev



cartopy-python/bin/pip install nose cython numpy pyshp
brew install geos
brew install proj
cartopy-python/bin/pip install shapely


# install PIL's dependencies (and then remove the actual pil build. There must be a better way...)
brew install --devel pil
brew remove pil
cartopy-python/bin/pip install pil


# for basic matplotlib
git clone https://github.com/matplotlib/matplotlib.git
cd matplotlib
cartopy-python/bin/python setup.py install


git clone https://github.com/SciTools/cartopy.git
cd cartopy
cartopy-python/bin/python setup.py install

Add CONTRIBUTING[.md]

https://github.com/blog/1184-contributing-guidelines

Path <-> LineString vs Polygon conversion

Consider the following:

p = Path([[0, 0], [0, 2], [0, 0]])
print cpatch.path_to_geos(p)

In this case, the result should be a LineString, but the following should be a polygon:

p = Path([[0, 0], [0, 2], [2, 2], [2, 0], [0, 0]])
print cpatch.path_to_geos(p)

Update cartopy.mpl_integration.patch (possibly renaming it) to handle these cases in the best possible way (which way is, as yet, unclear).

Add tests for these cases.

Projection documentation

Add the cartopy documentation for each projection.

It would be great to:

have groups of projections by features (conformal, equal area etc.) (perhaps as a graphical table which links to the primary documentation)
to have visual representations of each of the projections
to explain how one would go about adding a new projection and/or projection keyword

Tidy up the cartopy mpl interface.

There are currently a lot of methods which should not be there. Rationalise and simplify the geoaxes API.

ax.coastlines() produces empty plot

In [1]: import matplotlib

In [2]: matplotlib.use('gtkagg')

In [3]: import cartopy.crs as ccrs

In [4]: import matplotlib.pyplot as plt

In [5]: ax = plt.axes(projection=ccrs.Mercator())

In [6]: ax.coastlines()

In [7]: plt.show()

Here's sample output:

I have tried both the Mercator and Orthographic projections, the result is the same in both cases.

Mollweide global not producing a global plot

    import cartopy.crs as ccrs
    import matplotlib.pyplot as plt

    ax = plt.axes(projection=ccrs.Mollweide())
    ax.stock_img()

    # toggle this line
    ax.set_global()

    plt.show()

scitools / cartopy Goto Github PK

cartopy's Introduction

Cartopy is a Python package designed to make drawing maps for data analysis and visualisation easy.

Table of contents

Overview

Get in touch

Credits, copyright and license

cartopy's People

Contributors

Stargazers

Watchers

Forkers

cartopy's Issues

Recommend Projects

Recommend Topics

Recommend Org