plurigrid Map Projection Skill
Category theory of map projections: functors between manifolds with distortion analysis.
Trigger
- Map projection selection and analysis
- Distortion metrics (Tissot's indicatrix)
- Coordinate system transformations
- Cartographic design decisions
GF(3) Trit: +1 (Generator)
Generates projections from sphere to plane, creating new coordinate representations.
Category Theory of Projections
A map projection is a functor:
P: Sphere → Plane
S² → ℝ²
Different projections preserve different properties:
- Conformal (angle-preserving): Mercator, Stereographic
- Equal-area: Albers, Lambert, Mollweide
- Equidistant: Azimuthal equidistant
- Compromise: Robinson, Winkel Tripel
Projection Functors
import math
class Projection:
"""Base projection functor."""
def forward(self, lat, lon):
"""S² → ℝ²"""
raise NotImplementedError
def inverse(self, x, y):
"""ℝ² → S²"""
raise NotImplementedError
@property
def distortion_type(self):
raise NotImplementedError
class Mercator(Projection):
"""Conformal cylindrical projection."""
def forward(self, lat, lon):
x = math.radians(lon)
y = math.log(math.tan(math.pi/4 + math.radians(lat)/2))
return x, y
def inverse(self, x, y):
lon = math.degrees(x)
lat = math.degrees(2 * math.atan(math.exp(y)) - math.pi/2)
return lat, lon
@property
def distortion_type(self):
return "conformal" # Preserves angles
class LambertAzimuthal(Projection):
"""Equal-area azimuthal projection."""
def __init__(self, lat0=0, lon0=0):
self.lat0 = math.radians(lat0)
self.lon0 = math.radians(lon0)
def forward(self, lat, lon):
phi = math.radians(lat)
lam = math.radians(lon)
k = math.sqrt(2 / (1 + math.sin(self.lat0)*math.sin(phi) +
math.cos(self.lat0)*math.cos(phi)*math.cos(lam - self.lon0)))
x = k * math.cos(phi) * math.sin(lam - self.lon0)
y = k * (math.cos(self.lat0)*math.sin(phi) -
math.sin(self.lat0)*math.cos(phi)*math.cos(lam - self.lon0))
return x, y
@property
def distortion_type(self):
return "equal-area" # Preserves area
class Stereographic(Projection):
"""Conformal azimuthal projection."""
def __init__(self, lat0=90, lon0=0):
self.lat0 = math.radians(lat0)
self.lon0 = math.radians(lon0)
def forward(self, lat, lon):
phi = math.radians(lat)
lam = math.radians(lon)
k = 2 / (1 + math.sin(self.lat0)*math.sin(phi) +
math.cos(self.lat0)*math.cos(phi)*math.cos(lam - self.lon0))
x = k * math.cos(phi) * math.sin(lam - self.lon0)
y = k * (math.cos(self.lat0)*math.sin(phi) -
math.sin(self.lat0)*math.cos(phi)*math.cos(lam - self.lon0))
return x, y
@property
def distortion_type(self):
return "conformal"
Tissot's Indicatrix (Distortion Analysis)
def tissot_indicatrix(projection, lat, lon, delta=0.01):
"""
Compute Tissot's indicatrix at a point.
Returns semi-major axis a, semi-minor axis b, and rotation theta.
"""
# Jacobian via finite differences
x0, y0 = projection.forward(lat, lon)
x1, y1 = projection.forward(lat + delta, lon)
x2, y2 = projection.forward(lat, lon + delta)
# Partial derivatives
dx_dlat = (x1 - x0) / delta
dy_dlat = (y1 - y0) / delta
dx_dlon = (x2 - x0) / delta
dy_dlon = (y2 - y0) / delta
# Scale factors
h = math.sqrt(dx_dlat**2 + dy_dlat**2) # meridian scale
k = math.sqrt(dx_dlon**2 + dy_dlon**2) / math.cos(math.radians(lat)) # parallel scale
# Angular distortion
sin_theta = (dx_dlat * dy_dlon - dy_dlat * dx_dlon) / (h * k * math.cos(math.radians(lat)))
# Area distortion
area_factor = h * k * sin_theta
return {
'h': h, # meridian scale
'k': k, # parallel scale
'area_factor': area_factor,
'angular_distortion': math.degrees(math.asin(1 - abs(sin_theta)))
}
Projection with GF(3) Coloring
def project_with_color(projection, lat, lon, seed):
"""Project point and assign GF(3) color."""
x, y = projection.forward(lat, lon)
# Derive color from seed + projected coords
point_seed = int((seed + hash((x, y))) & 0x7FFFFFFFFFFFFFFF)
hue = point_seed % 360
if hue < 60 or hue >= 300:
trit = 1 # Red → Generator
elif hue < 180:
trit = 0 # Green → Ergodic
else:
trit = -1 # Blue → Validator
return {
'lat': lat,
'lon': lon,
'x': x,
'y': y,
'projection': projection.__class__.__name__,
'distortion_type': projection.distortion_type,
'seed': point_seed,
'trit': trit
}
Natural Transformations Between Projections
def projection_morphism(P1, P2, lat, lon):
"""
Natural transformation between projections.
P1 → P2 via S² (the universal object).
"""
# Forward through P1
x1, y1 = P1.forward(lat, lon)
# Inverse to sphere
lat_s, lon_s = P1.inverse(x1, y1)
# Forward through P2
x2, y2 = P2.forward(lat_s, lon_s)
return {
'source': (x1, y1),
'target': (x2, y2),
'sphere_point': (lat_s, lon_s),
'transformation': f"{P1.__class__.__name__} → {P2.__class__.__name__}"
}
# All projections form a category with S² as terminal object
# The "best" projection is context-dependent (no universal winner)
DuckDB Integration
-- Create projection lookup table
CREATE TABLE projections (
projection_id VARCHAR PRIMARY KEY,
name VARCHAR,
type VARCHAR, -- conformal, equal-area, equidistant, compromise
suitable_for VARCHAR[],
gf3_trit INTEGER DEFAULT 1 -- Generator
);
INSERT INTO projections VALUES
('mercator', 'Mercator', 'conformal', ['navigation', 'web_maps'], 1),
('albers', 'Albers Equal-Area', 'equal-area', ['thematic_maps', 'usa'], 1),
('robinson', 'Robinson', 'compromise', ['world_maps', 'education'], 1),
('utm', 'UTM', 'conformal', ['surveying', 'military'], 1);
-- Select projection based on use case
SELECT * FROM projections
WHERE 'navigation' = ANY(suitable_for);
Triads
map-projection (+1) ⊗ duckdb-spatial (0) ⊗ osm-topology (-1) = 0 ✓
map-projection (+1) ⊗ geodesic-manifold (0) ⊗ geohash-coloring (-1) = 0 ✓
References
- Snyder, "Map Projections: A Working Manual"
- PROJ library documentation
- Tissot's Indicatrix theory