Legacy geospatial plotting workflows frequently depend on Basemap, which encounters severe compilation failures on modern macOS architectures and has been official deprecated. Switching to Cartopy eliminates these dependency conflicts while delivering a projection-aware API that integrates natively with Matplotlib axes.

The migration requires shifting from Basemap's coordinate transformation methods to Cartopy's axis-level projection handling. Key architectural adjustments include:

• Initializing subplots using the projection keyword argument instead of instantiating a standalone map object.
• Replacing manual parallel and meridian drawing with ax.gridlines(), which automatically handles label placement and projection distortion.
• Masking invalid data points with np.nan rather than arbitrary sentinel values like -9999. This allows Matplotlib's contouring engine to graceful skip missing cells without rendering artifacts.
• Passing the transform=ccrs.PlateCarree() argument too plotting functions so Cartopy correctly reprojects the raw longitude/latitude coordinates to the target axis projection.

Legacy Implementation (Basemap)

import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap
import numpy as np

def render_spatial_comparison(longitudes, latitudes, forecast_data, ground_truth, day_idx, region_mask, settings, plot_flag):
    if plot_flag != 'Fig.2':
        return

    lon_grid, lat_grid = np.meshgrid(longitudes, latitudes)
    bounds = {
        'llcrnrlon': lon_grid.min(), 'llcrnrlat': lat_grid.min(),
        'urcrnrlon': lon_grid.max(), 'urcrnrlat': lat_grid.max()
    }

    fig = plt.figure(figsize=(14, 6))

    for idx, (data_array, title_str) in enumerate([(forecast_data, 'Model Forecast'), (ground_truth, 'Ground Truth')], start=1):
        ax = fig.add_subplot(1, 2, idx)
        geo_map = Basemap(**bounds)
        geo_map.drawcoastlines(linewidth=0.5)
        geo_map.drawcountries(linewidth=0.5)

        lat_ticks = np.arange(-90, 91, 18)
        lon_ticks = np.arange(-180, 180, 36)
        geo_map.drawparallels(lat_ticks, labels=[0, 1, 1, 0], dashes=[2, 2])
        geo_map.drawmeridians(lon_ticks, labels=[1, 0, 0, 1], dashes=[2, 2])

        x_proj, y_proj = geo_map(lon_grid, lat_grid)
        daily_slice = data_array[day_idx, :, :].copy()
        daily_slice[region_mask == 0] = -9999.0

        var_name = settings.get('target_variable')
        if var_name in ['soil_moisture_l1', 'soil_moisture_l2', 'soil_moisture_20cm']:
            levels = np.arange(0, 0.65, 0.05)
            cmap_type = 'YlGnBu'
            unit_label = 'm³/m³'
        elif var_name == 'sensible_heat_flux':
            levels = np.arange(-140, 150, 20)
            cmap_type = 'jet'
            unit_label = 'W/m²'
        else:
            continue

        contour_plot = geo_map.contourf(x_proj, y_proj, daily_slice, levels=levels, cmap=cmap_type)
        cbar = geo_map.colorbar(contour_plot, location='bottom', pad='8%')
        cbar.set_label(unit_label)
        ax.set_title(title_str)

    plt.tight_layout()
    plt.show()

Modern Implementation (Cartopy)

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

def render_spatial_comparison(longitudes, latitudes, forecast_data, ground_truth, day_idx, region_mask, settings, plot_flag):
    if plot_flag != 'Fig.2':
        return

    lon_grid, lat_grid = np.meshgrid(longitudes, latitudes)
    data_bounds = [longitudes.min(), longitudes.max(), latitudes.min(), latitudes.max()]

    fig, axes = plt.subplots(1, 2, figsize=(14, 6), subplot_kw={'projection': ccrs.PlateCarree()})

    for ax, (data_array, title_str) in zip(axes, [(forecast_data, 'Model Forecast'), (ground_truth, 'Ground Truth')]):
        ax.add_feature(cfeature.COASTLINE, linewidth=0.5)
        ax.add_feature(cfeature.BORDERS, linewidth=0.5)
        ax.set_extent(data_bounds, crs=ccrs.PlateCarree())

        gl = ax.gridlines(draw_labels=True, dms=True, x_inline=False, y_inline=False, linestyle='--', alpha=0.6)
        gl.top_labels = False
        gl.right_labels = False

        daily_slice = data_array[day_idx, :, :].copy()
        daily_slice[region_mask == 0] = np.nan

        var_name = settings.get('target_variable')
        if var_name in ['soil_moisture_l1', 'soil_moisture_l2', 'soil_moisture_20cm']:
            levels = np.arange(0, 0.65, 0.05)
            cmap_type = 'YlGnBu'
            unit_label = 'm³/m³'
        elif var_name == 'sensible_heat_flux':
            levels = np.arange(-140, 150, 20)
            cmap_type = 'jet'
            unit_label = 'W/m²'
        else:
            continue

        contour_plot = ax.contourf(lon_grid, lat_grid, daily_slice, levels=levels, cmap=cmap_type, transform=ccrs.PlateCarree())
        cbar = fig.colorbar(contour_plot, ax=ax, orientation='horizontal', pad=0.05, shrink=0.8)
        cbar.set_label(unit_label)
        ax.set_title(title_str)

    plt.tight_layout()
    plt.show()