Transformations

def warp_array(array, target_proj, source_proj=None, target_res=(400, 200),
               source_extent=None, target_extent=None, mask_extrapolated=False):
    """
    Warp numpy array to target projection.
    
    Parameters:
    - array: numpy.ndarray, source data array
    - target_proj: cartopy.crs.CRS, target projection
    - source_proj: cartopy.crs.CRS, source projection (default PlateCarree)
    - target_res: tuple, target resolution (width, height)
    - source_extent: tuple, source extent (x0, x1, y0, y1)  
    - target_extent: tuple, target extent (x0, x1, y0, y1)
    - mask_extrapolated: bool, mask extrapolated values
    
    Returns:
    Tuple of (warped_array, target_extent)
    """

def warp_img(fname, target_proj, source_proj=None, target_res=(400, 200)):
    """
    Warp image file to target projection.
    
    Parameters:
    - fname: str, path to image file
    - target_proj: cartopy.crs.CRS, target projection
    - source_proj: cartopy.crs.CRS, source projection (default PlateCarree)
    - target_res: tuple, target resolution (width, height)
    
    Returns:
    Tuple of (warped_array, target_extent)
    """

def regrid(array, source_x_coords, source_y_coords, source_proj, target_proj,
           target_res, mask_extrapolated=False):
    """
    Regrid array between coordinate systems.
    
    Parameters:
    - array: numpy.ndarray, source data
    - source_x_coords: array-like, source x coordinates
    - source_y_coords: array-like, source y coordinates  
    - source_proj: cartopy.crs.CRS, source projection
    - target_proj: cartopy.crs.CRS, target projection
    - target_res: tuple, target resolution (width, height)
    - mask_extrapolated: bool, mask extrapolated values
    
    Returns:
    Tuple of (regridded_array, target_x_coords, target_y_coords)
    """

def mesh_projection(projection, nx, ny, x_extents=(None, None), y_extents=(None, None)):
    """
    Generate coordinate mesh for projection.
    
    Parameters:
    - projection: cartopy.crs.CRS, target projection
    - nx: int, number of x grid points
    - ny: int, number of y grid points
    - x_extents: tuple, x extent limits (min, max)
    - y_extents: tuple, y extent limits (min, max)
    
    Returns:
    Tuple of (x_mesh, y_mesh) coordinate arrays
    """

def vector_scalar_to_grid(src_crs, target_proj, regrid_shape, x, y, u, v, *scalars, **kwargs):
    """
    Transform vector and scalar fields to regular grid.
    
    Parameters:
    - src_crs: cartopy.crs.CRS, source coordinate system
    - target_proj: cartopy.crs.CRS, target projection
    - regrid_shape: tuple, output grid shape (nx, ny)
    - x: array-like, source x coordinates
    - y: array-like, source y coordinates
    - u: array-like, x-component of vector field
    - v: array-like, y-component of vector field  
    - *scalars: additional scalar fields to transform
    - **kwargs: additional interpolation parameters
    
    Returns:
    Tuple of (target_x, target_y, target_u, target_v, *target_scalars)
    """

def add_cyclic_point(data, coord=None, axis=-1):
    """
    Add cyclic point to data array.
    
    Parameters:
    - data: array-like, input data array
    - coord: array-like, coordinate array (optional)
    - axis: int, axis along which to add cyclic point
    
    Returns:
    If coord provided: tuple of (extended_data, extended_coord)
    Otherwise: extended_data
    """

def add_cyclic(data, x=None, y=None, axis=-1, cyclic=360, precision=1e-4):
    """
    Add cyclic point with coordinate handling.
    
    Parameters:
    - data: array-like, input data
    - x: array-like, x coordinates (optional)
    - y: array-like, y coordinates (optional) 
    - axis: int, cyclic axis
    - cyclic: float, cyclic interval (default 360 for longitude)
    - precision: float, precision for cyclic detection
    
    Returns:
    Tuple of (extended_data, extended_x, extended_y) or subsets based on inputs
    """

def has_cyclic(x, axis=-1, cyclic=360, precision=1e-4):
    """
    Check if coordinates already have cyclic point.
    
    Parameters:
    - x: array-like, coordinate array
    - axis: int, axis to check
    - cyclic: float, expected cyclic interval
    - precision: float, precision for comparison
    
    Returns:
    bool, True if cyclic point exists
    """

def add_cyclic(data, x=None, y=None, axis=-1, cyclic=360, precision=1e-4):
    """
    Add cyclic point with coordinate handling (newer version of add_cyclic_point).
    
    Parameters:
    - data: array-like, input data
    - x: array-like, x coordinates (optional)
    - y: array-like, y coordinates (optional) 
    - axis: int, cyclic axis
    - cyclic: float, cyclic interval (default 360 for longitude)
    - precision: float, precision for cyclic detection
    
    Returns:
    Tuple of (extended_data, extended_x, extended_y) or subsets based on inputs
    """

class Geodesic:
    """Geodesic calculations on ellipsoid."""
    def __init__(self, a, f):
        """
        Parameters:
        - a: float, semi-major axis
        - f: float, flattening
        """
    
    def inverse(self, lat1, lon1, lat2, lon2):
        """
        Compute inverse geodesic problem.
        
        Parameters:
        - lat1, lon1: float, first point coordinates
        - lat2, lon2: float, second point coordinates
        
        Returns:
        Dict with distance, azimuth angles
        """
    
    def direct(self, lat1, lon1, azi1, s12):
        """
        Compute direct geodesic problem.
        
        Parameters:
        - lat1, lon1: float, starting point
        - azi1: float, initial azimuth
        - s12: float, distance
        
        Returns:
        Dict with end point coordinates and final azimuth
        """

import numpy as np
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
from cartopy.img_transform import warp_array

# Create sample data on regular lat/lon grid
lons = np.linspace(-180, 180, 360)
lats = np.linspace(-90, 90, 180)
lon_grid, lat_grid = np.meshgrid(lons, lats)
data = np.sin(np.radians(lat_grid)) * np.cos(np.radians(lon_grid))

# Warp to different projections
source_proj = ccrs.PlateCarree()
target_proj = ccrs.Orthographic(central_longitude=-90, central_latitude=45)

warped_data, target_extent = warp_array(
    data, 
    target_proj, 
    source_proj=source_proj,
    target_res=(400, 400)
)

# Plot results
fig, axes = plt.subplots(1, 2, figsize=(15, 6))

# Original data
ax1 = plt.subplot(1, 2, 1, projection=source_proj)
ax1.contourf(lon_grid, lat_grid, data, transform=source_proj)
ax1.coastlines()
ax1.set_global()
ax1.set_title('Original Data (PlateCarree)')

# Warped data  
ax2 = plt.subplot(1, 2, 2, projection=target_proj)
ax2.imshow(warped_data, extent=target_extent, transform=target_proj)
ax2.coastlines()
ax2.set_global()
ax2.set_title('Warped Data (Orthographic)')

plt.tight_layout()
plt.show()

import numpy as np
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
from cartopy.util import add_cyclic_point

# Create longitude data that doesn't include 360°
lons = np.arange(0, 360, 5)  # 0 to 355
lats = np.arange(-90, 91, 5)
lon_grid, lat_grid = np.meshgrid(lons, lats)

# Sample data
data = np.sin(np.radians(lat_grid)) * np.cos(np.radians(lon_grid))

# Add cyclic point to avoid gap at 180°/-180°
cyclic_data, cyclic_lons = add_cyclic_point(data, coord=lons)

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

# Without cyclic point
lon_mesh, lat_mesh = np.meshgrid(lons, lats)
axes[0].pcolormesh(lon_mesh, lat_mesh, data, transform=ccrs.PlateCarree())
axes[0].coastlines()
axes[0].set_global()
axes[0].set_title('Without Cyclic Point (Gap at Dateline)')

# With cyclic point
cyclic_lon_mesh, cyclic_lat_mesh = np.meshgrid(cyclic_lons, lats)
axes[1].pcolormesh(cyclic_lon_mesh, cyclic_lat_mesh, cyclic_data, 
                   transform=ccrs.PlateCarree())
axes[1].coastlines()
axes[1].set_global()
axes[1].set_title('With Cyclic Point (No Gap)')

plt.tight_layout()
plt.show()

import numpy as np
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
from cartopy.vector_transform import vector_scalar_to_grid

# Create irregular vector field data
np.random.seed(42)
x_pts = np.random.uniform(-180, 180, 100)
y_pts = np.random.uniform(-90, 90, 100)
u_pts = 10 * np.sin(np.radians(y_pts))
v_pts = 5 * np.cos(np.radians(x_pts))

# Transform to regular grid
source_crs = ccrs.PlateCarree()
target_proj = ccrs.PlateCarree()

# Regular grid output
x_reg, y_reg, u_reg, v_reg = vector_scalar_to_grid(
    source_crs, target_proj, 
    regrid_shape=(20, 15),  # 20x15 regular grid
    x=x_pts, y=y_pts, u=u_pts, v=v_pts
)

# Plot results
fig = plt.figure(figsize=(15, 10))

# Original irregular data
ax1 = plt.subplot(2, 1, 1, projection=ccrs.PlateCarree())
ax1.quiver(x_pts, y_pts, u_pts, v_pts, transform=ccrs.PlateCarree(),
           alpha=0.7, color='blue', scale=200)
ax1.coastlines()
ax1.set_global()
ax1.set_title('Original Irregular Vector Field')

# Interpolated regular grid
ax2 = plt.subplot(2, 1, 2, projection=ccrs.PlateCarree())
ax2.quiver(x_reg, y_reg, u_reg, v_reg, transform=ccrs.PlateCarree(),
           color='red', scale=200)
ax2.coastlines()
ax2.set_global()
ax2.set_title('Interpolated Regular Grid')

plt.tight_layout()
plt.show()

import numpy as np
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
from cartopy.img_transform import mesh_projection

# Generate coordinate mesh for Lambert Conformal projection
proj = ccrs.LambertConformal(central_longitude=-95, central_latitude=39)

# Create mesh over continental US region
x_mesh, y_mesh = mesh_projection(proj, nx=50, ny=40,
                                x_extents=(-2.5e6, 2.5e6),
                                y_extents=(-1.5e6, 1.5e6))

# Plot the mesh
fig = plt.figure(figsize=(12, 8))
ax = plt.axes(projection=proj)

# Plot mesh points
ax.scatter(x_mesh.flatten(), y_mesh.flatten(), s=1, 
           transform=proj, alpha=0.5)

# Add geographic features
ax.coastlines()
ax.add_feature(cfeature.BORDERS)
ax.set_extent([-125, -65, 25, 55], ccrs.PlateCarree())

plt.title('Lambert Conformal Projection Mesh')
plt.show()

from cartopy.geodesic import Geodesic
import matplotlib.pyplot as plt
import cartopy.crs as ccrs

# WGS84 ellipsoid parameters
geod = Geodesic(6378137.0, 1/298.257223563)

# Calculate great circle between two cities
# New York to London
ny_lat, ny_lon = 40.7128, -74.0060
london_lat, london_lon = 51.5074, -0.1278

result = geod.inverse(ny_lat, ny_lon, london_lat, london_lon)
distance_km = result['distance'] / 1000
azimuth = result['azi1']

print(f"Distance: {distance_km:.0f} km")
print(f"Initial bearing: {azimuth:.1f}°")

# Plot great circle path
fig = plt.figure(figsize=(12, 8))
ax = plt.axes(projection=ccrs.PlateCarree())

# Plot cities
ax.scatter([ny_lon, london_lon], [ny_lat, london_lat], 
           s=100, c=['red', 'blue'], transform=ccrs.PlateCarree())

# Add great circle (simplified - would need more points for accuracy)
ax.plot([ny_lon, london_lon], [ny_lat, london_lat], 
        'r--', transform=ccrs.Geodetic(), linewidth=2)

ax.coastlines()
ax.set_global()
ax.gridlines(draw_labels=True)
plt.title(f'Great Circle: New York to London ({distance_km:.0f} km)')
plt.show()

from cartopy.vector_transform import vector_scalar_to_grid

# Custom interpolation with different methods
x_reg, y_reg, u_reg, v_reg, temp_reg = vector_scalar_to_grid(
    ccrs.PlateCarree(), ccrs.LambertConformal(),
    regrid_shape=(30, 25),
    x=lon_data, y=lat_data, u=u_wind, v=v_wind, temperature_data,
    method='linear',  # or 'nearest', 'cubic'
    fill_value=np.nan
)

# Warp with extrapolation masking
warped_data, extent = warp_array(
    original_data, target_projection,
    source_proj=ccrs.PlateCarree(),
    mask_extrapolated=True  # Mask values outside original domain
)

# Masked values will be np.nan
valid_mask = ~np.isnan(warped_data)