contourf in Cartopy with ccrs.Mercator()

Question

fig = pl.figure(figsize=(8, 8))
ax1 = fig.add_subplot(1, 1, 1, projection=ccrs.Mercator())

cf = ax1.contourf(lon, lat, chl ,levels=np.linspace(0,5,50),cmap=cmaps.MPL_rainbow,
                transform=ccrs.Mercator(),extend='both')
ax1.coastlines()
ax1.add_feature(cfeature.LAND)
ax1.set_extent([9,17,67.5,69.75])

lat_formatter = LatitudeFormatter(degree_symbol='')
lon_formatter = LongitudeFormatter(degree_symbol='')
ax1.gridlines(draw_labels={"left":"y","top":'x'},xformatter=lon_formatter,yformatter=lat_formatter)

cax = fig.add_axes([0.92,0.2,0.03,0.6])
fig.colorbar(cf, cax=cax, ticks=np.arange(0,6,5),shrink=0.7,orientation='vertical')

the figure turns out blank, but if I replace Mercator() with PlateCarree(), the figure is right.

So what goes wrong with contourf using ccrs.Mercator()?

Answer 1

To get correct result, the coordinate transformation must be correct. In your case I suspect that the location data is in latitude/longitude. If I am correct, your error is "transform=ccrs.Mercator()" specified in your code.

Here is a demo code that shows a proper approach:

import matplotlib.pyplot as plt
import numpy as np
import cartopy.crs as ccrs
import cartopy.feature as cfeature

def sample_data(shape=(73, 145)):
    """Return ``lons``, ``lats`` and ``data`` of some fake data."""
    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

# Create a data set ...
# location data units: decimal degrees ==> ccrs.PlateCarree()
lons, lats, data = sample_data()

# Plot specs: Mercator projection
fig = plt.figure(figsize=(8, 8))
ax1 = fig.add_subplot(1, 1, 1, projection=ccrs.Mercator())

# Data specs: ccrs.PlateCarree()
# Location data is in decimal degrees of latitude/longitude ...
# which requires correct option ...
# transform=ccrs.PlateCarree()
cf = ax1.contourf(lons, lats, data,
                  levels=np.linspace(-1,1,50),
                  cmap='gist_rainbow',
                  transform=ccrs.PlateCarree(), 
                  extend='both', alpha=0.65, zorder=21)
ax1.coastlines(lw=0.3, zorder=22)
ax1.add_feature(cfeature.LAND, zorder=10)

# Use wider extent for clearer rendering
ax1.set_extent([-4+9, 4+17, -2+67.5, 2+69.75])
ax1.gridlines(draw_labels=True)

Output: