How to add a legend for a GeoAxes that adds a Cartopy shapely feature?

Question

I copied the code for adding legend via proxy artists from matplotlib's documentation but it doesn't work. I also tried the rest in matplotlib's legends guide but nothing works. I guess it's because the element is a shapely feature which ax.legend() somehow doesn't recognize.

Code

bounds = [116.9283371, 126.90534668, 4.58693981, 21.07014084]
stamen_terrain = cimgt.Stamen('terrain-background')

fault_line = ShapelyFeature(Reader('faultLines.shp').geometries(), ccrs.epsg(32651), 
                             linewidth=1, edgecolor='black', facecolor='none') # geometry is multilinestring
fig = plt.figure(figsize=(15,10))   
ax  = fig.add_subplot(1, 1, 1, projection=ccrs.PlateCarree())
ax.set_extent(bounds)
ax.add_image(stamen_terrain, 8)
a = ax.add_feature(fault_line, zorder=1, label='test')
ax.legend([a], loc='lower left', fancybox=True) #plt.legend() has the same result
plt.show()

Result

Answer 1

When copying the matplotlib example, you omitted the actual "proxy" artist line!

red_patch = mpatches.Patch(color='red', label='The red data')
plt.legend(handles=[red_patch])

That red_patch is the proxy artist. You have to create a dummy artist to pass to legend(). Your code as written is still passing the unrecognized Shapely feature.

It's tedious, but the relevant code would be something like:

fault_line = ShapelyFeature(Reader('faultLines.shp').geometries(), ccrs.epsg(32651), linewidth=1, edgecolor='black', facecolor='none')
ax.add_feature(fault_line, zorder=1)

# Now make a dummy object that looks as similar as possible
import matplotlib.patches as mpatches
proxy_artist = mpatches.Rectangle((0, 0), 1, 0.1, linewidth=1, edgecolor='black', facecolor='none')

# And manually add the labels here
ax.legend([proxy_artist], ['test'], loc='lower left', fancybox=True)

Here I just used a Rectangle, but depending on the feature, you can use various supported matplotlib "artists".