2

For a project I need to create a visualization that draws a circle around some locations on a map. The visualization used Cartopy v.0.18.0 to render the map. It uses the GoogleTiles class to fetch and display the tiles in the relevant region, and the add_patch(Patch.Circle(..., transform=ccrs.PlateCarree())) method to draw the circle.

tiles = GoogleTiles()
fig = plt.figure(figsize=(15,15))
ax = fig.add_subplot(1, 1, 1, projection=tiles.crs)

ax.set_extent((-121.8,-122.55,37.25,37.85))

ax.add_image(tiles, 11)

ax.add_patch(Patch.Circle(xy=[-122.4015173428571, 37.78774634285715], radius = 0.021709041989311614 + 0.005, alpha=0.3, zorder=30, transform=ccrs.PlateCarree()))

plt.show()

However, although I tried several transform objects I either got a ellipse instead of a circle (e.g. using ccrs.PlateCarree()) or no circle at all (e.g. using ccrs.Mercator()).

I found several different solutions online (e.g. Drawing Circles with cartopy in orthographic projection), however, these were not for the Mercator projection and I sadly lack the projection/transformation knowledge to adapt these to my problem.

Circle on Mercator projection

The only way I was able to produce a circular patch, was when I set the projection parameter on fig.add_subplot to ccrs.PlateCarree(). This, however, distorts the map and the labels become blured, so this is sadly not an acceptable solution.

As the project is due soon, a speedy reply would be much appreciated.

AxelEckenberger
  • 16,628
  • 3
  • 48
  • 70
  • Use a ScaledTranslation as detailed in https://matplotlib.org/stable/tutorials/advanced/transforms_tutorial.html – Jody Klymak Feb 27 '21 at 15:27
  • Thanks, gona try this, but I think this is not a very elegant solution, as the values for the transform will have to be determined empirically from what I see at the moment. Do you know whether these the necessary x and y scaling values could be obtaind from somewhere? The difference in length of one degree in latitudal/logitudal direction does not accuount for the distortion of the circle. – AxelEckenberger Feb 27 '21 at 15:38
  • @Jody Klymak sadly does not work with cartopy - or better to say I did not get it to work. – AxelEckenberger Feb 27 '21 at 15:47
  • 1
    Maybe `tissot's indicatrix` is better, try `ax.tissot(rad_km=15, lons=[-122.4015], lats=[37.7877])`. It will be a circle on Mercator projection. – swatchai Feb 27 '21 at 16:33

1 Answers1

2

Thanks @swatchai this was the missing hint, so for those intested the code looks like this right now, and it does work! Hooray!

tiles = GoogleTiles()
fig = plt.figure(figsize=(15,15))
ax = fig.add_subplot(1, 1, 1, projection=tiles.crs)

ax.set_extent((-121.8,-122.55,37.25,37.85))

ax.add_image(tiles, 11)

# The diameter is in degrees in EPSG:4326 coordinates therefore, the degrees have 
# to be converted to km. At 37N the degree latitude is 11.0977 km.
ax.tissot(rad_km=(0.021709041989311614 + 0.005) * 11.0977, lons=[-122.4015], lats=[37.7877], alpha=0.3)

plt.show()

When executing the above code the following warning is thrown but it has visible effect on the result:

/opt/conda/lib/python3.8/site-packages/cartopy/mpl/geoaxes.py:761: UserWarning: Approximating coordinate system <cartopy._crs.Geodetic object at 0x7fa4c7529770> with the PlateCarree projection.
  warnings.warn('Approximating coordinate system {!r} with the '

Circle on Mercator map

So thanks again @swatchai you saved my day!

AxelEckenberger
  • 16,628
  • 3
  • 48
  • 70