diff --git a/development.md b/development.md
new file mode 100644
index 0000000..c15e200
--- /dev/null
+++ b/development.md
@@ -0,0 +1,110 @@
+# Coordinate Handling Analysis (2026-06-09)
+
+Scope: OpenDRIVE import, normal/aerial view editing, OpenDRIVE export, pixel↔geo translation corner cases. Findings ordered by importance. Each was verified directly in code.
+
+## Architecture summary
+
+Four coordinate frames: **image pixels** (y-down, origin top-left) ↔ **geo lon/lat (WGS84)** ↔ **local metric meters** (y-north; equirectangular around control-point mean, OR absolute projected coords when `export_proj_string` is set) ↔ **OpenDRIVE x/y** (meters, x-east y-north, hdg CCW from east). The pixel↔geo mapping is a fitted Affine / Homography / Hybrid / DroneAssisted transformer (`orbit/utils/coordinate_transform.py`). Geometry is dual-stored: `polyline.points` (pixels, current background image) + `polyline.geo_points` (lon/lat). The y-flip is absorbed by the fitted matrices (no explicit negation except synthetic import, `opendrive_coordinate_transform.py:229` — correct). Heading conversion pixel→metric uses a 1-px differential step (`coordinate_transform.py:449-474`) — correct.
+
+## HIGH priority
+
+### H1. Aerial-view transformer ignores Web-Mercator nonlinearity — ✅ IMPLEMENTED (2026-06-10)
+**Change:** new `WebMercatorTransformer` (`orbit/utils/coordinate_transform.py`, subclass of `CoordinateTransformer`) maps pixel_x linearly in longitude and pixel_y linearly in spherical-Mercator y = ln(tan(π/4+φ/2)) — exact for stitched slippy-map tiles, no fitting involved. `create_transformer_from_bounds()` now returns it instead of a corner-fitted `AffineTransformer` (the aerial view switch was its only production caller, both call sites). Supports adjustments, `geo_to_pixel_unadjusted`, and `get_scale_factor` (meters/px at center latitude; x and y scales agree for square-Mercator tile grids). Earth radius cancels in the pixel↔geo mapping, so no datum constant enters it. Invalid bounds (non-positive dims, inverted bounds, |lat| ≥ 90) make the factory return None as before. Tests updated: the old assertions baked in the linear-latitude assumption (image center == latitude midpoint) — replaced with Mercator-midpoint assertions and resolution-proportionality checks, plus new tests for isotropic scale and invalid bounds (`tests/unit/test_utils/test_transformer_from_bounds.py`). Effect: editing in aerial view no longer bakes the mid-image N/S offset (~0.05–0.2 m at 57°N for 1–2 km extents, quadratically worse for larger areas) into geo coordinates.
+
+Original finding:
+`create_transformer_from_bounds()` (`orbit/utils/coordinate_transform.py:1415-1463`) fits an affine transform **linear in latitude** from the 4 tile-image corners. ESRI tiles are EPSG:3857: pixel_y is linear in mercator-y = ln(tan(π/4+φ/2)), not in φ. Corners are exact; mid-image geometry is shifted north/south by ≈ R·tan(φ)·Δφ²/8 — at 57°N: ~0.05 m for 1 km extent, ~0.2 m at 2 km, ~1.2 m at 5 km. **Because editing happens in aerial view, this error is baked into geo_points** via `pixel_to_geo` on the aerial transformer and survives the switch back. For survey-grade GCP work (<0.2 m), this alone can consume the error budget.
+
+**Fix:** make the aerial transformer mercator-aware — convert pixel_y↔lat through the mercator formula (tile y is exactly linear in mercator-y), or synthesize a grid of control points (not just corners) so the affine/homography fit averages the error down, or define the aerial transformer analytically from tile (zoom, x, y) instead of corner bounds.
+
+### H2. View switching trusts `geo_points` unconditionally — unsynced pixel edits silently reverted — ✅ IMPLEMENTED (2026-06-09)
+**Change:** new shared module `orbit/utils/geo_sync.py` + a `Polyline.geo_stale` dirty flag (set automatically by `add_point`/`insert_point`/`remove_point`/`update_point`, serialized in .orbit files, default False so imported polylines keep their geo precision). `reproject_project_geometry()` now calls `refresh_stale_geo_points(project, old_transformer, edited_only=True)` before reprojecting, so pixel edits whose geo sync was skipped are reconciled instead of reverted. **Design decision:** the view-switch refresh is `edited_only` — it reconciles only flagged (edited) polylines and length mismatches; *unflagged* drift keeps geo authoritative, because such drift can also come from adjustment/transformer changes where the pixels are the stale side (this matches the documented adjustment round-trip behavior and its tests). `Project` endpoint snapping (`enforce_road_links`) was also routed through `update_point` so it sets the flag. Covered by `tests/unit/test_utils/test_geo_sync.py` incl. an end-to-end "unsynced edit survives reprojection" test.
+
+Original finding:
+`reproject_project_geometry()` (`orbit/utils/reproject.py:44-51` and same pattern for junctions/signals/objects/parking): if an entity has geo coords, pixel positions are recomputed **purely from geo** — pixel edits whose geo update was skipped are discarded on aerial switch. Export protects against exactly this with `_validate_and_refresh_geo_points()` (`orbit/export/opendrive_writer.py:159-229`), but the aerial switch only resyncs junction/CR geo (`main_window.py:3147`), not polylines. Geo updates can be skipped because:
+
+- `on_polyline_modified` (`main_window.py:2415-2430`) only rebuilds geo when geo is *absent*, and only if `_cached_transformer` is non-None — it is frequently invalidated (`_invalidate_cached_transformer`, 8+ call sites) and only lazily recreated.
+- "Smooth Road Curve" sets `geo_points = None` (`image_view.py:3487`) and relies on the modified-handler rebuild; if the cached transformer is None at that moment the polyline silently loses geo.
+
+**Fix:** run the same pixel-vs-geo staleness refresh (using the *old* transformer) at the start of `reproject_project_geometry`, or share `_validate_and_refresh_geo_points` between writer and reproject so both consumers see one consistency contract.
+
+### H3. Exports lack the aerial-view guard that save has — ✅ IMPLEMENTED (2026-06-09)
+**Change:** `_ensure_original_view_for_save()` is now called at the top of all four export functions in `main_window.py` (`export_to_opendrive`, `export_to_osm`, `export_georeferencing`, `export_layout_mask`), matching `save_project`. Exporting while aerial view is active first switches back to the original image so the transformer, image dimensions, drone metadata, and adjustment are all consistent.
+
+Original finding:
+`save_project`/`save_project_as` call `_ensure_original_view_for_save()` (`main_window.py:757-765, 770, 793`); **none of the four exports do** (`export_to_opendrive:901`, `export_to_osm:932`, `export_georeferencing:1020`, `export_layout_mask:1130`). Exporting while aerial is active builds a transformer from control points that were reprojected into aerial pixel space (`reproject.py:156-160`) but with the **original** image dimensions (hybrid blend margin wrong, `main_window.py:2091-2092`), the original-image **drone metadata** (drone-assisted fit becomes meaningless), and the original-space adjustment.
+
+**Fix:** call `_ensure_original_view_for_save()` at the top of each export (cheap, matches the save behaviour).
+
+### H4. Export uses `geo_points` without length check; refresh skips mismatched lengths — ✅ IMPLEMENTED (2026-06-09)
+**Change:** length mismatches are now *rebuilt* from pixels (with warning log) instead of skipped, in the shared `refresh_stale_geo_points()` (used by both writer init and reproject). New helper `polyline_to_metric_points(polyline, transformer)` encodes the rule "geo only when lengths match, else pixels (warn)" and replaced the hand-rolled geo/pixel selection in `opendrive_writer._calculate_bounds`/`_create_road`, `object_builder._get_road_metrics`, and `parking_builder._get_road_metrics`. The two single-point sites (`_get_contact_heading`, `_get_road_endpoint_meters` in opendrive_writer) got explicit length-match guards. The writer's old `_validate_and_refresh_geo_points` method was deleted in favor of the shared function.
+
+Original finding:
+`_create_road` (`opendrive_writer.py:693-699`) uses `geo_points` whenever non-empty; `_validate_and_refresh_geo_points` deliberately **skips** polylines where `len(geo_points) != len(points)` (`opendrive_writer.py:178`). A mismatched polyline is exported from possibly-stale geo, and `cumulative_metric_s` (`opendrive_writer.py:720-724`) is indexed by geo-point count while lane-section boundaries index pixel points — silent section misalignment. Same unchecked pattern at `opendrive_writer.py:656, 1344, 1381`, `object_builder.py:110, 568` (568 does check), `parking_builder.py:141`, `osm_writer.py` callers.
+
+**Fix:** one helper that resolves a polyline's metric points with an explicit rule: lengths match → geo; mismatch → pixels (and log). Use it everywhere.
+
+## MEDIUM priority
+
+### M1. Auto-georeference control points: falsy-zero bug — ✅ IMPLEMENTED (2026-06-09)
+**Change:** `_generate_suggested_control_points` now uses `is None` checks on all four bounds (`data_min_x/max_x/min_y/max_y`) instead of truthiness, so OpenDRIVE files whose data bounds include 0.0 get suggested control points.
+
+Original finding:
+`_generate_suggested_control_points` (`orbit/import/opendrive_coordinate_transform.py:322`): `if not self.data_min_x or not self.data_max_x:` returns `[]` when a bound is exactly `0.0`. OpenDRIVE files often start at x=0 → importing such a file with geoReference silently produces **no** suggested control points. Fix: `is None` checks (and include `data_min_y/max_y`).
+
+### M2. `header_offset_hdg` parsed but never applied — ✅ IMPLEMENTED (2026-06-09)
+**Change:** `_metric_to_latlon` now rotates (x, y) by `header_offset_hdg` before adding the x/y offsets (rigid transform, rotate-first convention: `absolute = R(hdg)·p + T`; documented as an assumption in the docstring). All geo-point storage and suggested control points get the rotation automatically since they go through `_metric_to_latlon`. The synthetic (non-georeferenced) path is unaffected — a global rotation has no absolute frame to be wrong against there.
+
+**Caveat (untested convention):** the ASAM spec wording is ambiguous about operation order. We assume rotate-then-translate. If a test file with nonzero header `hdg` imports with geometry rotated about the wrong pivot, the spec intends translate-first — swap to `absolute = R(hdg)·(p + T)`, i.e. in `_metric_to_latlon` move the `x_absolute/y_absolute` offset addition *before* the rotation block and rotate the summed coordinates instead. Two-line change in `orbit/import/opendrive_coordinate_transform.py::_metric_to_latlon`.
+
+Original finding:
+`opendrive_coordinate_transform.py:72` stores it; `_metric_to_latlon` (`:247-249`) applies only x/y offsets. A .xodr with a heading offset in `<header><offset hdg=…>` imports rotated wrongly. Fix: rotate (x,y) by hdg before translating (per spec the offset is applied as rotation+translation), or explicitly warn "unsupported".
+
+### M3. `latlon_to_meters` returns two different frames depending on mode — ✅ IMPLEMENTED (2026-06-10)
+**Change:** the two frames are now explicit methods on `CoordinateTransformer`: `latlon_to_local_meters`/`local_meters_to_latlon` (always equirectangular around the control-point mean) and `latlon_to_projected`/`projected_to_latlon` (require the export projection, raise otherwise). `latlon_to_meters`/`meters_to_latlon` remain as mode-dependent dispatchers for the existing call sites, but their docstrings (and `pixel_to_meters`/`meters_to_pixel`) now state the frame ambiguity loudly instead of claiming a local origin unconditionally. Export call sites keep going through the dispatcher with the export projection set, which resolves to the projected frame end-to-end. Migrating callers to the explicit variants can happen incrementally.
+
+Original finding:
+`coordinate_transform.py:321-337`: with `_export_proj` set → **absolute projected** coords (e.g. UTM eastings ~10⁵–10⁶); without → **local** equirectangular meters around the control-point mean. Docstrings (`pixel_to_meters`, `:399-418`) claim local origin unconditionally. Consequences: the GUI transformer (no proj) and export transformer (proj) are *fitted in different metric frames* (homography DLT absorbs it, but results differ slightly), so geo_points written by one and validated by the other can drift; equirectangular (sphere R=6371 km) vs UTM differ by up to ~0.1% in scale. The a13a6d7 hybrid fix addressed the worst symptom; the dual-frame design remains a trap. Fix direction: separate `latlon_to_local_meters` from `latlon_to_projected`, and make export always use the projected variant end-to-end.
+
+### M4. Inconsistent meters-per-degree constants — ✅ IMPLEMENTED (2026-06-10)
+**Change:** new stdlib-only module `orbit/utils/geo_constants.py` with `EARTH_RADIUS_M = 6371000.0` and `METERS_PER_DEGREE = R·π/180 ≈ 111195`. All four `111000` sites now use it (`coordinate_transform.get_scale_factor`, `uncertainty_estimator`, `osm_to_orbit`, the import fallback in `opendrive_coordinate_transform`), and the three hardcoded `6371000.0` literals in `coordinate_transform.py` reference `EARTH_RADIUS_M`. Scale displays, lane-width rendering, the aerial resize ratio, and uncertainty radii now agree with `latlon_to_meters` (~0.18% shift from before).
+
+Original finding:
+`111000` m/deg in `get_scale_factor` (`coordinate_transform.py:810-811`), `uncertainty_estimator.py:644-645`, `osm_to_orbit.py:95-96`, import fallback (`opendrive_coordinate_transform.py:273-274`) vs R·π/180 ≈ **111195** in `latlon_to_meters`. ~0.18% systematic error in m/px scale → lane-width rendering, the aerial resize ratio (`main_window.py:3118-3127`), and uncertainty radii are all slightly off. Fix: one shared constant derived from R.
+
+### M5. Export staleness threshold leaves ≤2 px edits unexported — ✅ IMPLEMENTED (2026-06-09)
+**Change:** per-polyline dirty tracking via `Polyline.geo_stale` (see H2). The refresh in `geo_sync.py` uses a tight 0.1 px round-trip tolerance (`EDIT_EPSILON_PX`) for flagged polylines — sub-2px edits now reach exported geometry — while keeping the loose 2.0 px threshold (`TRANSFORM_DRIFT_THRESHOLD_PX`) for unflagged polylines so imported geo precision isn't degraded by transformer noise. Note: connecting-road `inline_path` has no dirty flag yet (Road model); CRs keep threshold-based reconciliation at export and rely on the GUI's edit-time geo sync + `_resync_junction_geo_coords`.
+
+Original finding:
+`PIXEL_THRESHOLD = 2.0` (`opendrive_writer.py:172`): pixel nudges smaller than 2 px (≈0.1–0.3 m typically) never reach the exported geometry when geo_points exist, because geo wins below the threshold. With sub-meter accuracy goals this is in the noise-vs-signal gray zone. Better contract: mark geo stale per-point at edit time (dirty flag) instead of distance-detection at export.
+
+### M6. Hybrid transformer: blend factor evaluated in different frames forward vs inverse — ✅ IMPLEMENTED (2026-06-10, finding partially revised)
+**Change + revision:** deeper analysis showed the suspected inconsistency strip does not exist: when `geo_to_pixel` takes the t≥1.0 shortcut it returns the *same* pixel the blend factor was evaluated at, so `pixel_to_geo` at that pixel also resolves to pure homography — consistent. In the t<1.0 band the Newton refinement already guarantees round-trip consistency. What remained was the hardening item: the adjusted `geo_to_pixel` called `self._affine.geo_to_pixel` (adjusted variant) where `geo_to_pixel_unadjusted` called the unadjusted one — correct only as long as sub-transformers never carry an adjustment. Both affine fallback calls now use `geo_to_pixel_unadjusted`, making the two variants identical except for the final hybrid-level adjustment.
+
+Original finding:
+`HybridTransformer.geo_to_pixel` picks the blend factor at the *homography-predicted* pixel (`coordinate_transform.py:1106`), `pixel_to_geo` at the *actual* pixel (`:1207`). Newton refinement (`_refine_inverse`, `:1165-1195`) restores round-trip consistency in the blend band, but is skipped when t≥1.0 (`:1108-1109`) even though the result can land in the band where `pixel_to_geo` blends — a small inconsistency strip just inside the image edge. Also `geo_to_pixel` (adjusted) calls `self._affine.geo_to_pixel` (`:1111`) where the unadjusted variant calls `geo_to_pixel_unadjusted` (`:1142`) — harmless only while sub-transformer adjustments are guaranteed None (`set_adjustment` doesn't propagate, `:1225-1227`).
+
+### M7. Aerial round-trip restore paths are asymmetric for pixel-only entities — ✅ IMPLEMENTED (2026-06-10, test added; GUI compensations unchanged)
+**Change:** added `tests/unit/test_utils/test_aerial_edit_roundtrip.py` — the integration test this finding asked for: polyline edited in aerial space *without* geo sync (worst case, relies on the H2 `geo_stale` flag), CR edited *with* edit-time geo sync (as the GUI does), and a no-edit round-trip invariance check; all through `reproject_project_geometry` + `refresh_stale_geo_points` against a real affine original transformer and a `WebMercatorTransformer` aerial. The GUI restore compensations in `_switch_to_original` (CR/junction pixel restores for drone out-of-FOV cases) were deliberately left unchanged — they remain necessary for the drone-assisted transformer; the test pins the model-level contract underneath them.
+
+Original finding:
+Entities that enter aerial view *without* geo coords get geo assigned via the **adjusted** original transformer (`reproject.py:52-62`); the return trip reprojects via the **unadjusted** original transformer (`main_window.py:3182-3188`), and `update_all_from_geo_coords` with the adjusted transformer only runs when an adjustment is active (`:3221-3222`). The CR/junction pixel restores (`:3230-3246`) then overwrite reprojected values with saved pre-aerial pixels — correct for out-of-FOV drone cases, but it also discards CR edits made in aerial view unless `_regenerate_all_junction_crs` rebuilds them (`:3248-3254`). This block of compensations is where aerial-edit anomalies are most likely to hide; deserves an integration test: edit polyline+CR in aerial → switch back → export → compare geo.
+
+## LOW priority / corner cases
+
+- **L1. UTM zone unclamped**: `int((lon+180)/6)+1` → zone 61 at lon=180°, no [1,60] clamp (`coordinate_transform.py:630, 653`). — ✅ IMPLEMENTED (2026-06-10): `get_utm_zone()` clamps to [1, 60]; `get_utm_projection_string()` now calls it instead of duplicating the formula.
+- **L2. Pole/antimeridian**: `cos(ref_lat)` division with no guard (`:365`); no longitude wraparound in equirectangular or `tile_fetcher.deg2num`. Irrelevant for Sweden, cheap to guard. — ✅ PARTIALLY IMPLEMENTED (2026-06-10): the equirectangular pair now wraps the longitude difference to ±180° on the way in, wraps the output longitude on the way out, and raises `ValueError` when the reference latitude is too close to a pole. **Not done:** `tile_fetcher` still cannot fetch a bbox spanning the antimeridian (would need tile-range split logic; irrelevant for current use).
+- **L3. No-pyproj import fallback** returns meters **as degrees** as "last resort" (`opendrive_coordinate_transform.py:282`) — should raise instead of importing garbage. — ✅ IMPLEMENTED (2026-06-09): the last-resort branch now raises `ValueError` ("pyproj not installed and geoReference has no +lat_0/+lon_0 origin"). pyproj is a hard dependency in pyproject.toml, so this path only fires in broken environments — failing loudly beats importing garbage coordinates.
+- **L4. Aerial resize uses scale_x only** (`main_window.py:3121-3127`): assumes the original image has square effective pixels; for oblique homography/drone images scale varies across the image, so the "match pixels/meter" resize is approximate (display-only effect). — ✅ IMPLEMENTED (2026-06-10): the resize ratio now uses the geometric mean of the x/y scale factors of both transformers.
+- **L5. Synthetic import scale not persisted**: `scale_pixels_per_meter` lives in `ImportResult` (`opendrive_importer.py:292`) but not in the Project — a synthetic (non-georeferenced) import has no recorded pixel↔meter scale for later export; round-trip relies on georeferencing being added. — ✅ IMPLEMENTED (2026-06-10): new `Project.import_scale_pixels_per_meter` field (serialized in .orbit), set by the importer whenever the transform setup yields a synthetic scale. No consumer yet — it preserves the data so a future synthetic-export path can use it.
+- **L6. Insert-point geo fallback copies the neighbour's geo** (`image_view.py:4074-4080`) when no transformer and interpolation fails — produces a duplicate geo vertex; would be caught by the H2/M5 staleness contract. — ✅ IMPLEMENTED (2026-06-09): `insert_point` now sets `geo_stale`, so the duplicate placeholder is recomputed from the inserted pixel at the next refresh (view switch or export). No change to the fallback itself.
+- **L7. `_export_proj` init duplicated** in `HybridTransformer.__init__` (`:1028-1030`) vs base `_set_reference_point` (`:304-305`) — the a13a6d7 bug class; consider having Hybrid call a shared init helper so the next subclass can't repeat it. — ✅ IMPLEMENTED (2026-06-10): single `CoordinateTransformer._init_export_proj()` helper; base `__init__`, `HybridTransformer.__init__` use it, and the redundant re-init in `_set_reference_point` was removed.
+
+## Suggested fix order
+
+1. ~~H3 (one-line guard per export) and M1, M2, L3~~ — ✅ done 2026-06-09.
+2. ~~H2 + H4 + M5 + L6 together — single geo/pixel sync contract~~ — ✅ done 2026-06-09 (`orbit/utils/geo_sync.py`, `Polyline.geo_stale`, shared by writer + reproject; 2531 tests pass).
+3. ~~H1 — mercator-aware aerial transformer~~ — ✅ done 2026-06-10 (`WebMercatorTransformer`; 2533 tests pass).
+4. ~~M3 + M4 — metric-frame and constant unification~~ — ✅ done 2026-06-10 (explicit frame methods + `geo_constants.py`).
+5. ~~M6, M7, L-items~~ — ✅ done 2026-06-10 (M6 finding partially revised; M7 integration test in `test_aerial_edit_roundtrip.py`; L1/L2/L4/L5/L7 implemented, L2 tile-fetcher antimeridian left open; 2536 tests pass).
+
+**Remaining open items:** L2 tile-fetcher antimeridian bbox split (irrelevant for current use); L3 follow-up none; M5 note — connecting-road `inline_path` still has no dirty flag (relies on edit-time sync); M2 caveat — header-offset rotation convention untested against a real file.
diff --git a/orbit/export/object_builder.py b/orbit/export/object_builder.py
index ec806ef..19246b7 100644
--- a/orbit/export/object_builder.py
+++ b/orbit/export/object_builder.py
@@ -11,6 +11,7 @@
 
 from orbit.models import Road
 from orbit.models.object import ObjectType, RoadObject
+from orbit.utils.geo_sync import polyline_to_metric_points
 
 
 class ObjectBuilder:
@@ -106,14 +107,8 @@ def _get_meter_centerline(self, road: Road) -> tuple:
         if not centerline:
             return [], 0.0
 
-        # Use geo coords directly if available (more precise)
-        if centerline.geo_points:
-            all_points_meters = [
-                self.transformer.latlon_to_meters(lat, lon)
-                for lon, lat in centerline.geo_points
-            ]
-        else:
-            all_points_meters = self.transformer.pixels_to_meters_batch(centerline.points)
+        # Uses geo coords when consistent with pixels (more precise)
+        all_points_meters = polyline_to_metric_points(centerline, self.transformer)
         geometry_elements = self.curve_fitter.fit_polyline(all_points_meters)
         road_length = sum(elem.length for elem in geometry_elements)
         return all_points_meters, road_length
diff --git a/orbit/export/opendrive_writer.py b/orbit/export/opendrive_writer.py
index 7bcd54a..e259cd9 100644
--- a/orbit/export/opendrive_writer.py
+++ b/orbit/export/opendrive_writer.py
@@ -14,6 +14,7 @@
 
 from orbit.models import Junction, Project, Road
 from orbit.utils import CoordinateTransformer
+from orbit.utils.geo_sync import polyline_to_metric_points, refresh_stale_geo_points
 from orbit.utils.logging_config import get_logger
 
 from .curve_fitting import CurveFitter, GeometryElement, GeometryType
@@ -150,83 +151,7 @@ def __init__(
         self.reference_warnings: List[str] = []
 
         # Ensure geo_points are consistent with the current transformer
-        self._validate_and_refresh_geo_points()
-
-    # ------------------------------------------------------------------
-    # Geo-point consistency
-    # ------------------------------------------------------------------
-
-    def _validate_and_refresh_geo_points(self):
-        """Refresh stale geo_points so the export uses accurate geo coords.
-
-        Checks each polyline's geo_points against the current transformer
-        by comparing ``geo_to_pixel(geo_point)`` with the stored pixel
-        position. Points that diverge beyond a threshold are recomputed
-        from pixels via ``pixel_to_geo``.
-
-        Also validates ``inline_geo_path`` on connecting roads.
-        """
-        if self.transformer is None:
-            return
-
-        PIXEL_THRESHOLD = 2.0  # px — flags even ~0.2 m drift at typical resolution
-
-        refreshed_polylines = 0
-        refreshed_points = 0
-
-        for polyline in self.project.polylines:
-            if not polyline.geo_points or len(polyline.geo_points) != len(polyline.points):
-                continue
-
-            stale_indices = []
-            for i, (px, py) in enumerate(polyline.points):
-                lon, lat = polyline.geo_points[i]
-                try:
-                    rpx, rpy = self.transformer.geo_to_pixel(lon, lat)
-                except Exception:
-                    stale_indices.append(i)
-                    continue
-                if abs(rpx - px) > PIXEL_THRESHOLD or abs(rpy - py) > PIXEL_THRESHOLD:
-                    stale_indices.append(i)
-
-            if stale_indices:
-                refreshed_polylines += 1
-                refreshed_points += len(stale_indices)
-                for i in stale_indices:
-                    px, py = polyline.points[i]
-                    new_lon, new_lat = self.transformer.pixel_to_geo(px, py)
-                    polyline.geo_points[i] = (new_lon, new_lat)
-
-        # Validate inline_geo_path on connecting roads
-        for road in self.project.roads:
-            if not road.inline_geo_path or not road.inline_path:
-                continue
-            if len(road.inline_geo_path) != len(road.inline_path):
-                continue
-
-            stale_indices = []
-            for i, (px, py) in enumerate(road.inline_path):
-                lon, lat = road.inline_geo_path[i]
-                try:
-                    rpx, rpy = self.transformer.geo_to_pixel(lon, lat)
-                except Exception:
-                    stale_indices.append(i)
-                    continue
-                if abs(rpx - px) > PIXEL_THRESHOLD or abs(rpy - py) > PIXEL_THRESHOLD:
-                    stale_indices.append(i)
-
-            if stale_indices:
-                refreshed_points += len(stale_indices)
-                for i in stale_indices:
-                    px, py = road.inline_path[i]
-                    new_lon, new_lat = self.transformer.pixel_to_geo(px, py)
-                    road.inline_geo_path[i] = (new_lon, new_lat)
-
-        if refreshed_points > 0:
-            logger.info(
-                "Refreshed %d stale geo_points across %d polyline(s)",
-                refreshed_points, refreshed_polylines,
-            )
+        refresh_stale_geo_points(self.project, self.transformer)
 
     def write(self, output_path: str) -> bool:
         """
@@ -652,14 +577,8 @@ def _calculate_bounds(self) -> dict:
 
         # Collect all polyline points and convert to meters
         for polyline in self.project.polylines:
-            # Use geo coords directly if available (more precise)
-            if polyline.geo_points:
-                for lon, lat in polyline.geo_points:
-                    x_m, y_m = self.transformer.latlon_to_meters(lat, lon)
-                    all_points_meters.append((x_m, y_m))
-            else:
-                points_meters = self.transformer.pixels_to_meters_batch(polyline.points)
-                all_points_meters.extend(points_meters)
+            all_points_meters.extend(
+                polyline_to_metric_points(polyline, self.transformer))
 
         if not all_points_meters:
             return {'north': 0.0, 'south': 0.0, 'east': 0.0, 'west': 0.0}
@@ -689,14 +608,8 @@ def _create_road(self, road: Road) -> Optional[etree.Element]:
         centerline_points_pixel = centerline.points
 
         # Transform centerline points to metric coordinates (meters)
-        # Use geo coords directly if available (more precise, avoids double conversion)
-        if centerline.geo_points:
-            all_points_meters = [
-                self.transformer.latlon_to_meters(lat, lon)
-                for lon, lat in centerline.geo_points
-            ]
-        else:
-            all_points_meters = self.transformer.pixels_to_meters_batch(centerline_points_pixel)
+        # Uses geo coords when consistent with pixels (more precise)
+        all_points_meters = polyline_to_metric_points(centerline, self.transformer)
 
         # Fit curves to metric coordinates
         geometry_elements = self.curve_fitter.fit_polyline(all_points_meters)
@@ -1341,7 +1254,7 @@ def _get_road_heading_at_contact_meters(self, centerline_id, contact_point):
 
         use_end = (contact_point == "end")
 
-        if polyline.geo_points:
+        if polyline.geo_points and len(polyline.geo_points) == len(polyline.points):
             idx0 = -2 if use_end else 0
             idx1 = -1 if use_end else 1
             lon0, lat0 = polyline.geo_points[idx0]
@@ -1378,7 +1291,7 @@ def _get_road_endpoint_meters(
 
         use_end = (contact_point == "end")
 
-        if polyline.geo_points:
+        if polyline.geo_points and len(polyline.geo_points) == len(polyline.points):
             idx = -1 if use_end else 0
             lon, lat = polyline.geo_points[idx]
             return self.transformer.latlon_to_meters(lat, lon)
diff --git a/orbit/export/parking_builder.py b/orbit/export/parking_builder.py
index f4bfa02..90dd71d 100644
--- a/orbit/export/parking_builder.py
+++ b/orbit/export/parking_builder.py
@@ -12,6 +12,7 @@
 
 from orbit.models import Road
 from orbit.models.parking import ParkingSpace
+from orbit.utils.geo_sync import polyline_to_metric_points
 
 
 def _project_point_onto_polyline(px: float, py: float, pts: List[tuple]):
@@ -137,14 +138,8 @@ def _get_road_metrics(self, road: Road):
         if not centerline:
             return 0.0, [], []
 
-        # Use geo coords directly if available (more precise)
-        if centerline.geo_points:
-            all_points_meters = [
-                self.transformer.latlon_to_meters(lat, lon)
-                for lon, lat in centerline.geo_points
-            ]
-        else:
-            all_points_meters = self.transformer.pixels_to_meters_batch(centerline.points)
+        # Uses geo coords when consistent with pixels (more precise)
+        all_points_meters = polyline_to_metric_points(centerline, self.transformer)
         geometry_elements = self.curve_fitter.fit_polyline(all_points_meters)
         return sum(elem.length for elem in geometry_elements), all_points_meters, geometry_elements
 
diff --git a/orbit/gui/main_window.py b/orbit/gui/main_window.py
index 9c367c2..a3124c2 100644
--- a/orbit/gui/main_window.py
+++ b/orbit/gui/main_window.py
@@ -902,6 +902,8 @@ def export_to_opendrive(self):
         """Export project to OpenDrive format."""
         from .dialogs.export_dialog import ExportDialog
 
+        self._ensure_original_view_for_save()
+
         if not self._prompt_and_handle_unapplied_adjustment():
             return
 
@@ -935,6 +937,8 @@ def export_to_osm(self):
 
         from orbit.export.osm_writer import export_to_osm
 
+        self._ensure_original_view_for_save()
+
         if not self._prompt_and_handle_unapplied_adjustment():
             return
 
@@ -1021,6 +1025,8 @@ def export_georeferencing(self):
         """Export georeferencing parameters to JSON file."""
         from orbit.export import export_georeferencing
 
+        self._ensure_original_view_for_save()
+
         if not self._check_provenance_ready():
             return
 
@@ -1129,6 +1135,8 @@ def export_georeferencing(self):
 
     def export_layout_mask(self):
         """Export lane segmentation mask and metadata JSON."""
+        self._ensure_original_view_for_save()
+
         # Validate prerequisites
         if not self.image_view.image_item:
             show_warning(self, "Cannot export: No image loaded.", "No Image")
@@ -3114,12 +3122,20 @@ def _switch_to_aerial(self):
             return
 
         # Resize aerial image so its pixels/meter matches the original image.
+        # Geometric mean of x/y scales: oblique source imagery has
+        # anisotropic scale, so neither axis alone is representative.
+        import math as _math
+
         import cv2 as _cv2
         orig_scale_x, orig_scale_y = self._original_transformer.get_scale_factor()
         aerial_scale_x, aerial_scale_y = aerial_transformer_raw.get_scale_factor()
         aerial_image = result.image
-        if orig_scale_x > 0 and aerial_scale_x > 0:
-            resize_ratio = aerial_scale_x / orig_scale_x
+        orig_scale = (_math.sqrt(orig_scale_x * orig_scale_y)
+                      if orig_scale_x > 0 and orig_scale_y > 0 else 0.0)
+        aerial_scale = (_math.sqrt(aerial_scale_x * aerial_scale_y)
+                        if aerial_scale_x > 0 and aerial_scale_y > 0 else 0.0)
+        if orig_scale > 0 and aerial_scale > 0:
+            resize_ratio = aerial_scale / orig_scale
             if abs(resize_ratio - 1.0) > 0.01:
                 new_w = max(1, round(w * resize_ratio))
                 new_h = max(1, round(h * resize_ratio))
diff --git a/orbit/import/opendrive_coordinate_transform.py b/orbit/import/opendrive_coordinate_transform.py
index af79855..a9c75d1 100644
--- a/orbit/import/opendrive_coordinate_transform.py
+++ b/orbit/import/opendrive_coordinate_transform.py
@@ -12,6 +12,8 @@
 from dataclasses import dataclass
 from typing import List, Optional, Tuple
 
+from orbit.utils.geo_constants import METERS_PER_DEGREE
+
 
 @dataclass
 class TransformMode:
@@ -237,6 +239,9 @@ def _metric_to_latlon(self, x_meters: float, y_meters: float) -> Tuple[float, fl
         Applies header offset before conversion - OpenDRIVE coordinates are
         relative to the offset specified in the header.
 
+        Assumes the offset is a rigid transform applied rotate-first:
+        absolute = R(hdg) * (x, y) + (offset_x, offset_y).
+
         Args:
             x_meters: X coordinate in meters (relative to header offset)
             y_meters: Y coordinate in meters (relative to header offset)
@@ -245,6 +250,13 @@ def _metric_to_latlon(self, x_meters: float, y_meters: float) -> Tuple[float, fl
             Tuple of (longitude, latitude) in decimal degrees
         """
         # Apply header offset - coordinates in file are relative to this offset
+        if self.header_offset_hdg:
+            cos_h = math.cos(self.header_offset_hdg)
+            sin_h = math.sin(self.header_offset_hdg)
+            x_meters, y_meters = (
+                x_meters * cos_h - y_meters * sin_h,
+                x_meters * sin_h + y_meters * cos_h,
+            )
         x_absolute = x_meters + self.header_offset_x
         y_absolute = y_meters + self.header_offset_y
 
@@ -268,18 +280,19 @@ def _metric_to_latlon(self, x_meters: float, y_meters: float) -> Tuple[float, fl
             origin_lat, origin_lon = self._extract_origin_from_proj4()
 
             if origin_lat and origin_lon:
-                # Simple approximation: meters to degrees
-                # At mid-latitudes: ~111000 m per degree latitude, ~111000*cos(lat) m per degree longitude
-                meters_per_degree_lat = 111000.0
-                meters_per_degree_lon = 111000.0 * math.cos(math.radians(origin_lat))
+                # Simple equirectangular approximation: meters to degrees
+                meters_per_degree_lat = METERS_PER_DEGREE
+                meters_per_degree_lon = METERS_PER_DEGREE * math.cos(math.radians(origin_lat))
 
                 lat = origin_lat + (y_absolute / meters_per_degree_lat)
                 lon = origin_lon + (x_absolute / meters_per_degree_lon)
 
                 return (lon, lat)
 
-            # Last resort: return metric coords as if they were degrees (will be wrong)
-            return (x_absolute, y_absolute)
+            raise ValueError(
+                "Cannot convert OpenDRIVE coordinates to lat/lon: pyproj is not "
+                "installed and the geoReference string has no +lat_0/+lon_0 origin."
+            )
 
     def _extract_origin_from_proj4(self) -> Tuple[Optional[float], Optional[float]]:
         """
@@ -319,7 +332,8 @@ def _generate_suggested_control_points(
         Returns:
             List of (pixel_x, pixel_y, lon, lat) tuples
         """
-        if not self.data_min_x or not self.data_max_x:
+        if (self.data_min_x is None or self.data_max_x is None
+                or self.data_min_y is None or self.data_max_y is None):
             return []
 
         # Create 4 control points at corners of data bounds
diff --git a/orbit/import/opendrive_importer.py b/orbit/import/opendrive_importer.py
index d376172..60ba3ef 100644
--- a/orbit/import/opendrive_importer.py
+++ b/orbit/import/opendrive_importer.py
@@ -290,6 +290,10 @@ def _setup_coordinate_transform(self, options: ImportOptions, result: ImportResu
         transform_result = self.coord_transform.setup_transform(sample_points)
         result.transform_mode = transform_result.mode
         result.scale_used = transform_result.scale_pixels_per_meter
+        # Persist the synthetic pixel scale so it survives save/load
+        if transform_result.scale_pixels_per_meter is not None:
+            self.project.import_scale_pixels_per_meter = (
+                transform_result.scale_pixels_per_meter)
 
         if not transform_result.success:
             if transform_result.mode == TransformMode.AUTO_GEOREFERENCE:
diff --git a/orbit/import/osm_to_orbit.py b/orbit/import/osm_to_orbit.py
index 79b0146..71c96a2 100644
--- a/orbit/import/osm_to_orbit.py
+++ b/orbit/import/osm_to_orbit.py
@@ -18,6 +18,7 @@
 from orbit.models.signal import SignalType
 from orbit.utils import CoordinateTransformer
 from orbit.utils.enum_formatting import format_enum_name
+from orbit.utils.geo_constants import METERS_PER_DEGREE
 from orbit.utils.geometry import calculate_path_length, find_point_at_distance_along_path, shorten_geo_points
 from orbit.utils.logging_config import get_logger
 
@@ -92,8 +93,8 @@ def calculate_bbox_from_center(center_lat: float, center_lon: float,
     Returns:
         Tuple of (min_lat, min_lon, max_lat, max_lon)
     """
-    dlat = radius_m / 111_000
-    dlon = radius_m / (111_000 * math.cos(math.radians(center_lat)))
+    dlat = radius_m / METERS_PER_DEGREE
+    dlon = radius_m / (METERS_PER_DEGREE * math.cos(math.radians(center_lat)))
     return (center_lat - dlat, center_lon - dlon,
             center_lat + dlat, center_lon + dlon)
 
diff --git a/orbit/models/polyline.py b/orbit/models/polyline.py
index d13ba15..96b78e9 100644
--- a/orbit/models/polyline.py
+++ b/orbit/models/polyline.py
@@ -131,24 +131,29 @@ class Polyline:
     s_offsets: Optional[List[float]] = None  # S-coordinate for each point (if available)
     osm_node_ids: Optional[List[Optional[int]]] = None  # OSM node IDs for each point (from OSM import)
     geometry_segments: Optional[List[GeometrySegment]] = None  # Preserved geometry from OpenDRIVE import
+    geo_stale: bool = False  # Pixel points edited since geo_points last synced
 
     def add_point(self, x: float, y: float) -> None:
         """Add a point to the end of the polyline."""
         self.points.append((x, y))
+        self.geo_stale = True
 
     def insert_point(self, index: int, x: float, y: float) -> None:
         """Insert a point at the specified index."""
         self.points.insert(index, (x, y))
+        self.geo_stale = True
 
     def remove_point(self, index: int) -> None:
         """Remove a point at the specified index."""
         if 0 <= index < len(self.points):
             self.points.pop(index)
+            self.geo_stale = True
 
     def update_point(self, index: int, x: float, y: float) -> None:
         """Update the coordinates of a point at the specified index."""
         if 0 <= index < len(self.points):
             self.points[index] = (x, y)
+            self.geo_stale = True
 
     def get_point(self, index: int) -> Tuple[float, float]:
         """Get the coordinates of a point at the specified index."""
@@ -230,6 +235,8 @@ def to_dict(self) -> Dict[str, Any]:
             data['osm_node_ids'] = self.osm_node_ids
         if self.geometry_segments is not None:
             data['geometry_segments'] = [seg.to_dict() for seg in self.geometry_segments]
+        if self.geo_stale:
+            data['geo_stale'] = True
         return data
 
     @classmethod
@@ -270,7 +277,8 @@ def from_dict(cls, data: Dict[str, Any]) -> 'Polyline':
             elevations=data.get('elevations'),
             s_offsets=data.get('s_offsets'),
             osm_node_ids=data.get('osm_node_ids'),
-            geometry_segments=geometry_segments
+            geometry_segments=geometry_segments,
+            geo_stale=data.get('geo_stale', False)
         )
 
     def __repr__(self) -> str:
diff --git a/orbit/models/project.py b/orbit/models/project.py
index 2117483..8ed5710 100644
--- a/orbit/models/project.py
+++ b/orbit/models/project.py
@@ -219,6 +219,7 @@ class Project:
     imported_geo_reference: Optional[str] = None  # Preserved geoReference from OpenDRIVE import
     imported_origin_latitude: Optional[float] = None  # Back-projected origin lat from imported OpenDRIVE header offset
     imported_origin_longitude: Optional[float] = None  # Back-projected origin lon from imported OpenDRIVE header offset
+    import_scale_pixels_per_meter: Optional[float] = None  # Pixel scale of synthetic OpenDRIVE import
     enabled_sign_libraries: List[str] = field(default_factory=lambda: ['se'])  # Enabled sign library IDs
     synthetic_canvas_width: Optional[int] = None  # Synthetic canvas width in pixels (no real image)
     synthetic_canvas_height: Optional[int] = None  # Synthetic canvas height in pixels (no real image)
@@ -1509,7 +1510,7 @@ def enforce_road_link_coordinates(self, road_id: str) -> bool:
                                   else pred_cl.points[0])
                     # This road's start connects to predecessor
                     if centerline.points[0] != pred_point:
-                        centerline.points[0] = pred_point
+                        centerline.update_point(0, pred_point[0], pred_point[1])
                         changed = True
 
         # Enforce successor link
@@ -1524,7 +1525,9 @@ def enforce_road_link_coordinates(self, road_id: str) -> bool:
                                   else succ_cl.points[-1])
                     # This road's end connects to successor
                     if centerline.points[-1] != succ_point:
-                        centerline.points[-1] = succ_point
+                        centerline.update_point(
+                            len(centerline.points) - 1,
+                            succ_point[0], succ_point[1])
                         changed = True
 
         return changed
@@ -1732,6 +1735,7 @@ def to_dict(self) -> Dict[str, Any]:
             'imported_geo_reference': self.imported_geo_reference,
             'imported_origin_latitude': self.imported_origin_latitude,
             'imported_origin_longitude': self.imported_origin_longitude,
+            'import_scale_pixels_per_meter': self.import_scale_pixels_per_meter,
             'enabled_sign_libraries': self.enabled_sign_libraries,
             'synthetic_canvas_width': self.synthetic_canvas_width,
             'synthetic_canvas_height': self.synthetic_canvas_height,
@@ -1794,6 +1798,7 @@ def from_dict(cls, data: Dict[str, Any]) -> 'Project':
             imported_geo_reference=data.get('imported_geo_reference'),
             imported_origin_latitude=data.get('imported_origin_latitude'),
             imported_origin_longitude=data.get('imported_origin_longitude'),
+            import_scale_pixels_per_meter=data.get('import_scale_pixels_per_meter'),
             enabled_sign_libraries=data.get('enabled_sign_libraries', ['se']),
             synthetic_canvas_width=data.get('synthetic_canvas_width'),
             synthetic_canvas_height=data.get('synthetic_canvas_height'),
diff --git a/orbit/utils/coordinate_transform.py b/orbit/utils/coordinate_transform.py
index 3f89559..5ff21e8 100644
--- a/orbit/utils/coordinate_transform.py
+++ b/orbit/utils/coordinate_transform.py
@@ -25,6 +25,7 @@
 import numpy as np
 from pyproj import Proj
 
+from .geo_constants import EARTH_RADIUS_M, METERS_PER_DEGREE
 from .logging_config import get_logger
 
 if TYPE_CHECKING:
@@ -245,8 +246,7 @@ def __init__(self, control_points: List['ControlPoint'], use_validation: bool =
         self.use_validation = use_validation
 
         # Export projection (pyproj-based conversion instead of equirectangular)
-        self._export_proj_string: Optional[str] = export_proj_string
-        self._export_proj: Optional[Proj] = None
+        self._init_export_proj(export_proj_string)
 
         # Separate training (GCP) and validation (GVP) points
         if use_validation:
@@ -291,25 +291,29 @@ def has_adjustment(self) -> bool:
         """Check if an adjustment is currently applied."""
         return self.adjustment is not None and not self.adjustment.is_identity()
 
-    def _set_reference_point(self):
-        """Set reference point as mean of all control points.
+    def _init_export_proj(self, export_proj_string: Optional[str]):
+        """Initialize the export projection (single place — see a13a6d7 bug class)."""
+        self._export_proj_string: Optional[str] = export_proj_string
+        self._export_proj: Optional[Proj] = (
+            Proj(export_proj_string) if export_proj_string else None)
 
-        Also initializes the pyproj Proj if export_proj_string was provided.
-        """
+    def _set_reference_point(self):
+        """Set reference point as mean of all control points."""
         if not self.all_control_points:
             return
         self.reference_lat = np.mean([cp.latitude for cp in self.all_control_points])
         self.reference_lon = np.mean([cp.longitude for cp in self.all_control_points])
 
-        if self._export_proj_string:
-            self._export_proj = Proj(self._export_proj_string)
-
     def latlon_to_meters(self, lat: float, lon: float) -> Tuple[float, float]:
         """
         Convert lat/lon to metric coordinates.
 
-        When export_proj_string is set, uses pyproj for accurate projection.
-        Otherwise uses equirectangular approximation (suitable for areas <100km).
+        FRAME DEPENDS ON MODE: with export_proj_string set this returns
+        absolute projected coordinates (e.g. UTM eastings/northings) via
+        pyproj; otherwise local meters relative to the control-point mean
+        (equirectangular, suitable for areas <100 km). Use
+        latlon_to_local_meters / latlon_to_projected when a specific
+        frame is required.
 
         Args:
             lat: Latitude in decimal degrees
@@ -319,29 +323,15 @@ def latlon_to_meters(self, lat: float, lon: float) -> Tuple[float, float]:
             Tuple of (east, north) in meters
         """
         if self._export_proj is not None:
-            return self._export_proj(lon, lat)
-
-        if self.reference_lat is None or self.reference_lon is None:
-            raise ValueError("Reference point not set. Call compute transformation first.")
-
-        R = 6371000.0  # Earth radius in meters
-
-        lat_rad = math.radians(lat)
-        lon_rad = math.radians(lon)
-        ref_lat_rad = math.radians(self.reference_lat)
-        ref_lon_rad = math.radians(self.reference_lon)
-
-        east = R * (lon_rad - ref_lon_rad) * math.cos(ref_lat_rad)
-        north = R * (lat_rad - ref_lat_rad)
-
-        return east, north
+            return self.latlon_to_projected(lat, lon)
+        return self.latlon_to_local_meters(lat, lon)
 
     def meters_to_latlon(self, east: float, north: float) -> Tuple[float, float]:
         """
         Convert metric coordinates back to lat/lon.
 
-        When export_proj_string is set, uses pyproj inverse projection.
-        Otherwise uses equirectangular approximation.
+        FRAME DEPENDS ON MODE — exact inverse of latlon_to_meters
+        (projected when export_proj_string is set, local otherwise).
 
         Args:
             east: East coordinate in meters
@@ -351,20 +341,64 @@ def meters_to_latlon(self, east: float, north: float) -> Tuple[float, float]:
             Tuple of (latitude, longitude) in decimal degrees
         """
         if self._export_proj is not None:
-            lon, lat = self._export_proj(east, north, inverse=True)
-            return lat, lon
+            return self.projected_to_latlon(east, north)
+        return self.local_meters_to_latlon(east, north)
+
+    def latlon_to_projected(self, lat: float, lon: float) -> Tuple[float, float]:
+        """Convert lat/lon to absolute projected coords via the export projection."""
+        if self._export_proj is None:
+            raise ValueError("No export projection set (export_proj_string).")
+        return self._export_proj(lon, lat)
+
+    def projected_to_latlon(self, east: float, north: float) -> Tuple[float, float]:
+        """Convert absolute projected coords back to (lat, lon)."""
+        if self._export_proj is None:
+            raise ValueError("No export projection set (export_proj_string).")
+        lon, lat = self._export_proj(east, north, inverse=True)
+        return lat, lon
+
+    def latlon_to_local_meters(self, lat: float, lon: float) -> Tuple[float, float]:
+        """Convert lat/lon to local meters around the control-point mean.
+
+        Equirectangular approximation; longitude difference is wrapped so
+        areas near the antimeridian convert correctly.
+        """
+        if self.reference_lat is None or self.reference_lon is None:
+            raise ValueError("Reference point not set. Call compute transformation first.")
+
+        R = EARTH_RADIUS_M
 
+        lat_rad = math.radians(lat)
+        lon_rad = math.radians(lon)
+        ref_lat_rad = math.radians(self.reference_lat)
+        ref_lon_rad = math.radians(self.reference_lon)
+
+        dlon_rad = (lon_rad - ref_lon_rad + math.pi) % (2 * math.pi) - math.pi
+        east = R * dlon_rad * math.cos(ref_lat_rad)
+        north = R * (lat_rad - ref_lat_rad)
+
+        return east, north
+
+    def local_meters_to_latlon(self, east: float, north: float) -> Tuple[float, float]:
+        """Convert local meters back to (lat, lon)."""
         if self.reference_lat is None or self.reference_lon is None:
             raise ValueError("Reference point not set. Call compute transformation first.")
 
-        R = 6371000.0
+        R = EARTH_RADIUS_M
         ref_lat_rad = math.radians(self.reference_lat)
         ref_lon_rad = math.radians(self.reference_lon)
 
+        cos_ref = math.cos(ref_lat_rad)
+        if abs(cos_ref) < 1e-9:
+            raise ValueError(
+                "Reference latitude too close to a pole for equirectangular conversion")
+
         lat_rad = ref_lat_rad + (north / R)
-        lon_rad = ref_lon_rad + (east / (R * math.cos(ref_lat_rad)))
+        lon_rad = ref_lon_rad + (east / (R * cos_ref))
+        lon = math.degrees(lon_rad)
+        lon = (lon + 180.0) % 360.0 - 180.0
 
-        return math.degrees(lat_rad), math.degrees(lon_rad)
+        return math.degrees(lat_rad), lon
 
     def compute_transformation(self):
         """Compute transformation matrix. Must be implemented by subclass."""
@@ -398,16 +432,18 @@ def geo_to_pixel(self, longitude: float, latitude: float) -> Tuple[float, float]
 
     def pixel_to_meters(self, pixel_x: float, pixel_y: float) -> Tuple[float, float]:
         """
-        Convert pixel coordinates to local metric coordinates (meters).
+        Convert pixel coordinates to metric coordinates (meters).
 
-        Uses the center of the control points as the local origin (0, 0).
+        Frame follows latlon_to_meters: absolute projected coords when an
+        export projection is set, otherwise local meters around the
+        control-point mean.
 
         Args:
             pixel_x: X coordinate in pixels
             pixel_y: Y coordinate in pixels
 
         Returns:
-            Tuple of (x_meters, y_meters) in local metric coordinate system
+            Tuple of (x_meters, y_meters)
         """
         # First convert to geographic coordinates
         lon, lat = self.pixel_to_geo(pixel_x, pixel_y)
@@ -419,7 +455,10 @@ def pixel_to_meters(self, pixel_x: float, pixel_y: float) -> Tuple[float, float]
 
     def meters_to_pixel(self, x_meters: float, y_meters: float) -> Tuple[float, float]:
         """
-        Convert local metric coordinates (meters) to pixel coordinates.
+        Convert metric coordinates (meters) to pixel coordinates.
+
+        Frame follows meters_to_latlon (projected or local, see
+        latlon_to_meters).
 
         Args:
             x_meters: X coordinate in meters (east)
@@ -626,8 +665,7 @@ def get_utm_projection_string(self) -> str:
         Returns:
             PROJ4 UTM projection string (e.g., "+proj=utm +zone=33 ...")
         """
-        # Calculate UTM zone from longitude
-        zone = int((self.reference_lon + 180) / 6) + 1
+        zone = self.get_utm_zone()
 
         # Determine hemisphere
         hemisphere = "+north" if self.reference_lat >= 0 else "+south"
@@ -650,7 +688,9 @@ def get_utm_zone(self) -> int:
         Returns:
             UTM zone number (1-60)
         """
-        return int((self.reference_lon + 180) / 6) + 1
+        # Clamp: lon == 180 would yield zone 61, lon == -180 belongs to zone 1
+        zone = int((self.reference_lon + 180) / 6) + 1
+        return min(max(zone, 1), 60)
 
     def get_transformation_info(self) -> Dict:
         """
@@ -807,8 +847,8 @@ def get_scale_factor(self) -> Tuple[float, float]:
         a3, a4 = self.transform_matrix[1, 0], self.transform_matrix[1, 1]
 
         # Convert degrees per pixel to meters per pixel
-        meters_per_degree_lat = 111000
-        meters_per_degree_lon = 111000 * math.cos(math.radians(self.reference_lat))
+        meters_per_degree_lat = METERS_PER_DEGREE
+        meters_per_degree_lon = METERS_PER_DEGREE * math.cos(math.radians(self.reference_lat))
 
         # Scale in x direction
         scale_x = math.sqrt((a0 * meters_per_degree_lon)**2 + (a3 * meters_per_degree_lat)**2)
@@ -1024,10 +1064,7 @@ def __init__(self, control_points: List['ControlPoint'],
                  image_width: int = 0, image_height: int = 0):
         # Don't call super().__init__ — we delegate to sub-transformers
         self.all_control_points = control_points
-        self._export_proj_string = export_proj_string
-        self._export_proj = None
-        if export_proj_string:
-            self._export_proj = Proj(export_proj_string)
+        self._init_export_proj(export_proj_string)
 
         # Create both sub-transformers from the same control points
         self._homography = HomographyTransformer(
@@ -1100,7 +1137,7 @@ def geo_to_pixel(self, longitude: float, latitude: float) -> Tuple[float, float]
         w = pixel_homo[2]
 
         if abs(w) < 1e-10:
-            px, py = self._affine.geo_to_pixel(longitude, latitude)
+            px, py = self._affine.geo_to_pixel_unadjusted(longitude, latitude)
         else:
             hpx, hpy = pixel_homo[0] / w, pixel_homo[1] / w
             t = self._blend_factor(hpx, hpy, w)
@@ -1108,7 +1145,7 @@ def geo_to_pixel(self, longitude: float, latitude: float) -> Tuple[float, float]
             if t >= 1.0:
                 px, py = hpx, hpy
             else:
-                apx, apy = self._affine.geo_to_pixel(longitude, latitude)
+                apx, apy = self._affine.geo_to_pixel_unadjusted(longitude, latitude)
                 if t <= 0.0:
                     px, py = apx, apy
                 else:
@@ -1412,6 +1449,92 @@ def create_transformer(
         return None
 
 
+class WebMercatorTransformer(CoordinateTransformer):
+    """Exact pixel<->geo mapping for imagery on a Web-Mercator pixel grid.
+
+    Assumes pixel_x is linear in longitude and pixel_y is linear in
+    spherical-Mercator y = ln(tan(pi/4 + lat/2)), which holds exactly for
+    stitched slippy-map tiles (EPSG:3857). A latitude-linear affine fit on
+    the same corner bounds misplaces mid-image latitudes by
+    ~R*tan(lat)*dlat^2/8 (meters), growing quadratically with bbox height.
+    """
+
+    def __init__(self, image_width: int, image_height: int,
+                 min_lon: float, min_lat: float,
+                 max_lon: float, max_lat: float,
+                 export_proj_string: Optional[str] = None):
+        if image_width <= 0 or image_height <= 0:
+            raise ValueError("Image dimensions must be positive")
+        if not (max_lon > min_lon and max_lat > min_lat):
+            raise ValueError("Bounds must satisfy max > min")
+        if not (-90.0 < min_lat and max_lat < 90.0):
+            raise ValueError("Latitudes must be strictly within (-90, 90)")
+
+        from orbit.models.project import ControlPoint
+        corners = [
+            ControlPoint(pixel_x=0.0, pixel_y=0.0,
+                         longitude=min_lon, latitude=max_lat, name="NW"),
+            ControlPoint(pixel_x=float(image_width), pixel_y=0.0,
+                         longitude=max_lon, latitude=max_lat, name="NE"),
+            ControlPoint(pixel_x=float(image_width), pixel_y=float(image_height),
+                         longitude=max_lon, latitude=min_lat, name="SE"),
+            ControlPoint(pixel_x=0.0, pixel_y=float(image_height),
+                         longitude=min_lon, latitude=min_lat, name="SW"),
+        ]
+        super().__init__(corners, use_validation=False,
+                         export_proj_string=export_proj_string)
+
+        self.image_width = float(image_width)
+        self.image_height = float(image_height)
+        self.min_lon = min_lon
+        self.max_lon = max_lon
+        self._merc_top = self._lat_to_merc(max_lat)
+        self._merc_bottom = self._lat_to_merc(min_lat)
+
+        self._set_reference_point()
+        self.compute_reprojection_error()
+
+    @staticmethod
+    def _lat_to_merc(lat: float) -> float:
+        return math.log(math.tan(math.pi / 4 + math.radians(lat) / 2))
+
+    @staticmethod
+    def _merc_to_lat(merc: float) -> float:
+        return math.degrees(2 * math.atan(math.exp(merc)) - math.pi / 2)
+
+    def compute_transformation(self):
+        """Analytic transform — nothing to fit."""
+
+    def pixel_to_geo(self, pixel_x: float, pixel_y: float) -> Tuple[float, float]:
+        if self.adjustment is not None:
+            pixel_x, pixel_y = self.adjustment.apply_inverse_to_point(pixel_x, pixel_y)
+        lon = self.min_lon + (pixel_x / self.image_width) * (self.max_lon - self.min_lon)
+        merc = self._merc_top - (pixel_y / self.image_height) * (self._merc_top - self._merc_bottom)
+        return lon, self._merc_to_lat(merc)
+
+    def geo_to_pixel_unadjusted(self, longitude: float, latitude: float) -> Tuple[float, float]:
+        pixel_x = (longitude - self.min_lon) / (self.max_lon - self.min_lon) * self.image_width
+        merc = self._lat_to_merc(latitude)
+        pixel_y = (self._merc_top - merc) / (self._merc_top - self._merc_bottom) * self.image_height
+        return pixel_x, pixel_y
+
+    def geo_to_pixel(self, longitude: float, latitude: float) -> Tuple[float, float]:
+        pixel_x, pixel_y = self.geo_to_pixel_unadjusted(longitude, latitude)
+        if self.adjustment is not None:
+            pixel_x, pixel_y = self.adjustment.apply_to_point(pixel_x, pixel_y)
+        return pixel_x, pixel_y
+
+    def get_scale_factor(self) -> Tuple[float, float]:
+        """Meters per pixel at the reference (center) latitude."""
+        R = EARTH_RADIUS_M
+        cos_ref = math.cos(math.radians(self.reference_lat))
+        lon_span_rad = math.radians(self.max_lon - self.min_lon)
+        scale_x = R * cos_ref * lon_span_rad / self.image_width
+        # dNorth = R*cos(lat)*dMerc for spherical Mercator
+        scale_y = R * cos_ref * (self._merc_top - self._merc_bottom) / self.image_height
+        return scale_x, scale_y
+
+
 def create_transformer_from_bounds(
     image_width: int,
     image_height: int,
@@ -1420,13 +1543,12 @@ def create_transformer_from_bounds(
     max_lon: float,
     max_lat: float,
     export_proj_string: Optional[str] = None,
-) -> Optional[AffineTransformer]:
+) -> Optional[WebMercatorTransformer]:
     """
-    Create an AffineTransformer from image dimensions and geographic bounds.
+    Create a WebMercatorTransformer from image dimensions and geographic bounds.
 
-    Synthesizes corner control points mapping image corners to geographic
-    coordinates. Ideal for nadir (orthophoto/satellite) imagery where the
-    relationship between pixels and geography is a simple affine transform.
+    Intended for stitched web-map tile imagery (aerial view), where pixel
+    rows are spaced linearly in Web-Mercator y rather than in latitude.
 
     Image layout convention: pixel (0,0) = top-left = (min_lon, max_lat),
     pixel (W,H) = bottom-right = (max_lon, min_lat).
@@ -1438,27 +1560,15 @@ def create_transformer_from_bounds(
         export_proj_string: Optional pyproj string for metric conversions.
 
     Returns:
-        AffineTransformer if successful, None on error.
+        WebMercatorTransformer if successful, None on error.
     """
-    from orbit.models.project import ControlPoint
-
-    corners = [
-        ControlPoint(pixel_x=0.0, pixel_y=0.0,
-                     longitude=min_lon, latitude=max_lat, name="NW"),
-        ControlPoint(pixel_x=float(image_width), pixel_y=0.0,
-                     longitude=max_lon, latitude=max_lat, name="NE"),
-        ControlPoint(pixel_x=float(image_width), pixel_y=float(image_height),
-                     longitude=max_lon, latitude=min_lat, name="SE"),
-        ControlPoint(pixel_x=0.0, pixel_y=float(image_height),
-                     longitude=min_lon, latitude=min_lat, name="SW"),
-    ]
-
     try:
-        return AffineTransformer(
-            corners, use_validation=False,
+        return WebMercatorTransformer(
+            image_width, image_height,
+            min_lon, min_lat, max_lon, max_lat,
             export_proj_string=export_proj_string,
         )
-    except (ValueError, np.linalg.LinAlgError) as e:
+    except (ValueError, OverflowError) as e:
         logger.error(f"Error creating bounds transformer: {e}")
         return None
 
diff --git a/orbit/utils/geo_constants.py b/orbit/utils/geo_constants.py
new file mode 100644
index 0000000..5cf1abd
--- /dev/null
+++ b/orbit/utils/geo_constants.py
@@ -0,0 +1,13 @@
+"""
+Shared geodetic constants.
+
+Single source for the spherical-Earth approximation used across the
+codebase so scale computations agree with latlon_to_meters (M4 in
+development.md). Stdlib-only on purpose — importable from modules that
+must not pull in numpy/pyproj.
+"""
+
+import math
+
+EARTH_RADIUS_M = 6371000.0
+METERS_PER_DEGREE = EARTH_RADIUS_M * math.pi / 180.0  # ~111195 m/deg
diff --git a/orbit/utils/geo_sync.py b/orbit/utils/geo_sync.py
new file mode 100644
index 0000000..ab859f6
--- /dev/null
+++ b/orbit/utils/geo_sync.py
@@ -0,0 +1,124 @@
+"""
+Geo/pixel consistency for dual-stored geometry.
+
+Contract: geo coords are authoritative for output, pixel coords for
+editing. Point-mutating edits set Polyline.geo_stale; refresh recomputes
+only points whose geo no longer round-trips to the stored pixel, so
+imported geo precision is preserved for untouched points.
+"""
+
+import logging
+from typing import List, Tuple
+
+logger = logging.getLogger(__name__)
+
+# Drift below this is treated as transformer noise for un-edited polylines
+# (refreshing would needlessly degrade imported geo precision).
+TRANSFORM_DRIFT_THRESHOLD_PX = 2.0
+# Round-trip tolerance for edited (geo_stale) polylines: catches real
+# sub-2px edits while absorbing float/projection-frame noise.
+EDIT_EPSILON_PX = 0.1
+
+
+def _refresh_point_pairs(points, geo_points, transformer, threshold) -> int:
+    """Recompute geo entries that diverge from their pixel point (mutates geo_points)."""
+    refreshed = 0
+    for i, (px, py) in enumerate(points):
+        lon, lat = geo_points[i]
+        try:
+            rpx, rpy = transformer.geo_to_pixel(lon, lat)
+            stale = abs(rpx - px) > threshold or abs(rpy - py) > threshold
+        except Exception:
+            stale = True
+        if stale:
+            geo_points[i] = transformer.pixel_to_geo(px, py)
+            refreshed += 1
+    return refreshed
+
+
+def refresh_polyline_geo_points(polyline, transformer,
+                                edited_only: bool = False) -> int:
+    """Bring polyline.geo_points back in sync with polyline.points.
+
+    Length mismatch rebuilds geo entirely from pixels; an edited
+    (geo_stale) polyline gets a tight round-trip scan.  Unflagged drift
+    is reconciled with a loose threshold only when ``edited_only`` is
+    False — such drift can also come from transformer/adjustment changes
+    where geo is the authority and pixels are the stale side.
+    """
+    if not polyline.geo_points:
+        polyline.geo_stale = False
+        return 0
+    if len(polyline.geo_points) != len(polyline.points):
+        logger.warning(
+            "Polyline %s: geo_points length %d != points length %d — "
+            "rebuilding geo from pixels",
+            polyline.id, len(polyline.geo_points), len(polyline.points),
+        )
+        polyline.geo_points = [
+            transformer.pixel_to_geo(px, py) for px, py in polyline.points
+        ]
+        polyline.geo_stale = False
+        return len(polyline.points)
+    if polyline.geo_stale:
+        threshold = EDIT_EPSILON_PX
+    elif edited_only:
+        return 0
+    else:
+        threshold = TRANSFORM_DRIFT_THRESHOLD_PX
+    refreshed = _refresh_point_pairs(
+        polyline.points, polyline.geo_points, transformer, threshold)
+    polyline.geo_stale = False
+    return refreshed
+
+
+def refresh_stale_geo_points(project, transformer,
+                             edited_only: bool = False) -> int:
+    """Refresh stale geo coords on all polylines and connecting-road paths.
+
+    With ``edited_only`` (used before view-switch reprojection) only
+    explicitly edited polylines and length mismatches are reconciled;
+    unflagged drift keeps geo as the authority.
+    """
+    if transformer is None:
+        return 0
+    refreshed = 0
+    for polyline in project.polylines:
+        refreshed += refresh_polyline_geo_points(
+            polyline, transformer, edited_only=edited_only)
+    for road in project.roads:
+        if not road.inline_geo_path or not road.inline_path:
+            continue
+        if len(road.inline_geo_path) != len(road.inline_path):
+            logger.warning(
+                "Connecting road %s: inline_geo_path length %d != inline_path "
+                "length %d — rebuilding geo from pixels",
+                road.id, len(road.inline_geo_path), len(road.inline_path),
+            )
+            road.inline_geo_path = [
+                transformer.pixel_to_geo(px, py) for px, py in road.inline_path
+            ]
+            refreshed += len(road.inline_path)
+            continue
+        if not edited_only:
+            refreshed += _refresh_point_pairs(
+                road.inline_path, road.inline_geo_path, transformer,
+                TRANSFORM_DRIFT_THRESHOLD_PX)
+    if refreshed:
+        logger.info("Refreshed %d stale geo point(s)", refreshed)
+    return refreshed
+
+
+def polyline_to_metric_points(polyline, transformer) -> List[Tuple[float, float]]:
+    """Convert a polyline to metric points; geo only when consistent with pixels."""
+    if polyline.geo_points and len(polyline.geo_points) == len(polyline.points):
+        return [
+            transformer.latlon_to_meters(lat, lon)
+            for lon, lat in polyline.geo_points
+        ]
+    if polyline.geo_points:
+        logger.warning(
+            "Polyline %s: geo/pixel length mismatch — converting from pixel coords",
+            polyline.id,
+        )
+    return transformer.pixels_to_meters_batch(polyline.points)
diff --git a/orbit/utils/reproject.py b/orbit/utils/reproject.py
index 3add82f..4212d4d 100644
--- a/orbit/utils/reproject.py
+++ b/orbit/utils/reproject.py
@@ -10,6 +10,7 @@
 
 from orbit.models.project import Project
 from orbit.utils.coordinate_transform import CoordinateTransformer
+from orbit.utils.geo_sync import refresh_stale_geo_points
 
 logger = logging.getLogger(__name__)
 
@@ -38,6 +39,15 @@ def reproject_project_geometry(
     Returns:
         Number of entities re-projected.
     """
+    # Bring geo coords in sync with any pixel edits made in the current
+    # pixel space before geo_to_pixel overwrites the pixel positions —
+    # otherwise unsynced pixel edits would be silently reverted.
+    # edited_only: unflagged drift may come from adjustment/transformer
+    # changes where geo is the authority, so only explicit edits and
+    # length mismatches are reconciled here.
+    if old_transformer is not None:
+        refresh_stale_geo_points(project, old_transformer, edited_only=True)
+
     count = 0
 
     # --- Polylines ---
diff --git a/orbit/utils/uncertainty_estimator.py b/orbit/utils/uncertainty_estimator.py
index 05c1621..270871b 100644
--- a/orbit/utils/uncertainty_estimator.py
+++ b/orbit/utils/uncertainty_estimator.py
@@ -13,6 +13,7 @@
 from scipy.spatial import ConvexHull, distance
 
 from .coordinate_transform import CoordinateTransformer
+from .geo_constants import METERS_PER_DEGREE
 from .logging_config import get_logger
 
 logger = get_logger(__name__)
@@ -641,8 +642,8 @@ def _geo_to_meters(self, lon: float, lat: float) -> Tuple[float, float]:
         ref_lon = np.mean(ref_lons)
 
         # Simple equirectangular approximation
-        lat_m_per_deg = 111000.0
-        lon_m_per_deg = 111000.0 * math.cos(math.radians(ref_lat))
+        lat_m_per_deg = METERS_PER_DEGREE
+        lon_m_per_deg = METERS_PER_DEGREE * math.cos(math.radians(ref_lat))
 
         mx = (lon - ref_lon) * lon_m_per_deg
         my = (lat - ref_lat) * lat_m_per_deg
diff --git a/pyproject.toml b/pyproject.toml
index 4f18375..a580ba0 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -1,6 +1,6 @@
 [project]
 name = "orbit"
-version = "0.11.0"
+version = "0.11.1"
 description = "OpenDrive Road Builder from Imagery Tool"
 readme = "README.md"
 requires-python = ">=3.10"
diff --git a/tests/unit/test_utils/test_aerial_edit_roundtrip.py b/tests/unit/test_utils/test_aerial_edit_roundtrip.py
new file mode 100644
index 0000000..73cde86
--- /dev/null
+++ b/tests/unit/test_utils/test_aerial_edit_roundtrip.py
@@ -0,0 +1,117 @@
+"""Integration test: edits made in aerial view must survive the round-trip
+back to the original image and reach export-ready geo coordinates (M7)."""
+
+import pytest
+
+from orbit.models.junction import Junction
+from orbit.models.polyline import LineType, Polyline
+from orbit.models.project import ControlPoint, Project
+from orbit.models.road import Road
+from orbit.utils.coordinate_transform import (
+    create_transformer,
+    create_transformer_from_bounds,
+)
+from orbit.utils.geo_sync import refresh_stale_geo_points
+from orbit.utils.reproject import reproject_project_geometry
+
+
+def _make_project():
+    """Project with an affine 'original image' transformer, one road polyline
+    and one connecting road, all geo-synced."""
+    cps = [
+        ControlPoint(pixel_x=100.0, pixel_y=100.0,
+                     longitude=12.940, latitude=57.720, name="A"),
+        ControlPoint(pixel_x=900.0, pixel_y=100.0,
+                     longitude=12.950, latitude=57.720, name="B"),
+        ControlPoint(pixel_x=100.0, pixel_y=700.0,
+                     longitude=12.940, latitude=57.714, name="C"),
+        ControlPoint(pixel_x=900.0, pixel_y=700.0,
+                     longitude=12.950, latitude=57.714, name="D"),
+    ]
+    project = Project(control_points=cps)
+    orig_t = create_transformer(cps, "affine")
+
+    pl = Polyline(id="pl1", line_type=LineType.CENTERLINE)
+    pl.points = [(200.0, 200.0), (400.0, 300.0), (600.0, 400.0)]
+    pl.geo_points = [orig_t.pixel_to_geo(x, y) for x, y in pl.points]
+    project.polylines.append(pl)
+
+    cr = Road(name="CR1", junction_id="j1",
+              inline_path=[(600.0, 400.0), (650.0, 420.0), (700.0, 450.0)])
+    cr.inline_geo_path = [orig_t.pixel_to_geo(x, y) for x, y in cr.inline_path]
+    project.add_road(cr)
+
+    j = Junction(center_point=(650.0, 420.0))
+    j.geo_center_point = orig_t.pixel_to_geo(650.0, 420.0)
+    j.add_connecting_road(cr.id)
+    project.junctions.append(j)
+
+    return project, orig_t, pl, cr
+
+
+AERIAL_BOUNDS = dict(min_lon=12.938, min_lat=57.712, max_lon=12.952, max_lat=57.722)
+
+
+class TestAerialEditRoundTrip:
+
+    def test_polyline_edit_without_geo_sync_survives(self):
+        """Worst case: a point dragged in aerial view whose geo sync was
+        skipped — the geo_stale flag must carry the edit back (H2/M7)."""
+        project, orig_t, pl, _ = _make_project()
+        aerial_t = create_transformer_from_bounds(1400, 1000, **AERIAL_BOUNDS)
+
+        reproject_project_geometry(project, orig_t, aerial_t)
+
+        # Edit in aerial pixel space through the model mutator (sets geo_stale),
+        # deliberately without the GUI's edit-time geo sync.
+        ax, ay = pl.points[1]
+        pl.update_point(1, ax + 25.0, ay - 10.0)
+        edited_geo = aerial_t.pixel_to_geo(ax + 25.0, ay - 10.0)
+
+        reproject_project_geometry(project, aerial_t, orig_t)
+
+        assert pl.geo_points[1] == pytest.approx(edited_geo, abs=1e-9)
+        expected_px = orig_t.geo_to_pixel(*edited_geo)
+        assert pl.points[1] == pytest.approx(expected_px, abs=1e-6)
+
+        # Export-time refresh must not move anything further
+        assert refresh_stale_geo_points(project, orig_t) == 0
+
+    def test_cr_edit_with_geo_sync_survives(self):
+        """CR mid-point edit, geo synced at edit time (as the GUI does)."""
+        project, orig_t, _, cr = _make_project()
+        aerial_t = create_transformer_from_bounds(1400, 1000, **AERIAL_BOUNDS)
+
+        reproject_project_geometry(project, orig_t, aerial_t)
+
+        ax, ay = cr.inline_path[1]
+        cr.inline_path[1] = (ax + 15.0, ay + 20.0)
+        cr.inline_geo_path[1] = aerial_t.pixel_to_geo(ax + 15.0, ay + 20.0)
+        edited_geo = cr.inline_geo_path[1]
+
+        reproject_project_geometry(project, aerial_t, orig_t)
+
+        assert cr.inline_geo_path[1] == pytest.approx(edited_geo, abs=1e-9)
+        expected_px = orig_t.geo_to_pixel(*edited_geo)
+        assert cr.inline_path[1] == pytest.approx(expected_px, abs=1e-6)
+
+    def test_unedited_geometry_returns_exactly(self):
+        """Round-trip without edits must reproduce original pixels and
+        leave geo coordinates untouched."""
+        project, orig_t, pl, cr = _make_project()
+        original_pl_geo = list(pl.geo_points)
+        original_pl_px = list(pl.points)
+        original_cr_px = list(cr.inline_path)
+        aerial_t = create_transformer_from_bounds(1400, 1000, **AERIAL_BOUNDS)
+
+        reproject_project_geometry(project, orig_t, aerial_t)
+        reproject_project_geometry(project, aerial_t, orig_t)
+
+        # Geo untouched (it is the round-trip invariant)
+        for got, exp in zip(pl.geo_points, original_pl_geo):
+            assert got == pytest.approx(exp, abs=1e-12)
+        # Pixels reproduced from geo via the same transformer
+        for got, exp in zip(pl.points, original_pl_px):
+            assert got == pytest.approx(exp, abs=1e-6)
+        for got, exp in zip(cr.inline_path, original_cr_px):
+            assert got == pytest.approx(exp, abs=1e-6)
diff --git a/tests/unit/test_utils/test_geo_sync.py b/tests/unit/test_utils/test_geo_sync.py
new file mode 100644
index 0000000..27e274a
--- /dev/null
+++ b/tests/unit/test_utils/test_geo_sync.py
@@ -0,0 +1,210 @@
+"""
+Unit tests for geo/pixel consistency helpers (orbit.utils.geo_sync).
+"""
+
+import pytest
+
+from orbit.models import Project, Road
+from orbit.models.polyline import Polyline
+from orbit.utils.geo_sync import (
+    polyline_to_metric_points,
+    refresh_polyline_geo_points,
+    refresh_stale_geo_points,
+)
+
+# Degrees per pixel for the fake transformer (~1.1 m/px at the equator)
+DEG_PER_PX = 1e-5
+
+
+class FakeTransformer:
+    """Linear pixel<->geo mapping with a configurable pixel offset."""
+
+    def __init__(self, offset_px: float = 0.0):
+        self.offset_px = offset_px
+
+    def pixel_to_geo(self, px, py):
+        return ((px - self.offset_px) * DEG_PER_PX,
+                -(py - self.offset_px) * DEG_PER_PX)
+
+    def geo_to_pixel(self, lon, lat):
+        return (lon / DEG_PER_PX + self.offset_px,
+                -lat / DEG_PER_PX + self.offset_px)
+
+    def latlon_to_meters(self, lat, lon):
+        return (lon * 111320.0, lat * 111320.0)
+
+    def pixels_to_meters_batch(self, pixels):
+        result = []
+        for px, py in pixels:
+            lon, lat = self.pixel_to_geo(px, py)
+            result.append(self.latlon_to_meters(lat, lon))
+        return result
+
+
+def make_synced_polyline(transformer, points):
+    """Polyline whose geo_points match its pixel points exactly."""
+    return Polyline(
+        id='pl1',
+        points=list(points),
+        geo_points=[transformer.pixel_to_geo(px, py) for px, py in points],
+    )
+
+
+class TestGeoStaleFlag:
+    """Point mutators must mark geo as stale; constructor must not."""
+
+    def test_constructor_default_not_stale(self):
+        assert Polyline(id='p', points=[(0, 0)]).geo_stale is False
+
+    @pytest.mark.parametrize('mutate', [
+        lambda pl: pl.add_point(5.0, 5.0),
+        lambda pl: pl.insert_point(1, 5.0, 5.0),
+        lambda pl: pl.update_point(0, 5.0, 5.0),
+        lambda pl: pl.remove_point(0),
+    ])
+    def test_mutators_set_stale(self, mutate):
+        pl = Polyline(id='p', points=[(0.0, 0.0), (10.0, 10.0)])
+        mutate(pl)
+        assert pl.geo_stale is True
+
+    def test_serialization_roundtrip(self):
+        pl = Polyline(id='p', points=[(0.0, 0.0)])
+        pl.update_point(0, 1.0, 1.0)
+        restored = Polyline.from_dict(pl.to_dict())
+        assert restored.geo_stale is True
+        # Flag omitted from dict when not stale
+        assert 'geo_stale' not in Polyline(id='q', points=[]).to_dict()
+
+
+class TestRefreshPolyline:
+    """refresh_polyline_geo_points threshold behavior."""
+
+    def test_synced_polyline_untouched(self):
+        t = FakeTransformer()
+        pl = make_synced_polyline(t, [(0.0, 0.0), (100.0, 50.0)])
+        original_geo = list(pl.geo_points)
+        assert refresh_polyline_geo_points(pl, t) == 0
+        assert pl.geo_points == original_geo
+
+    def test_unflagged_subthreshold_drift_preserved(self):
+        """Sub-2px drift on un-edited polylines keeps imported geo precision."""
+        t = FakeTransformer()
+        pl = make_synced_polyline(t, [(0.0, 0.0), (100.0, 50.0)])
+        pl.points[1] = (101.0, 50.0)  # direct mutation: 1px, no flag
+        original_geo = list(pl.geo_points)
+        assert refresh_polyline_geo_points(pl, t) == 0
+        assert pl.geo_points == original_geo
+
+    def test_flagged_subthreshold_edit_refreshed(self):
+        """A sub-2px edit through the mutators reaches geo (M5)."""
+        t = FakeTransformer()
+        pl = make_synced_polyline(t, [(0.0, 0.0), (100.0, 50.0)])
+        pl.update_point(1, 101.0, 50.0)  # 1px edit, sets geo_stale
+        assert refresh_polyline_geo_points(pl, t) == 1
+        assert pl.geo_points[1] == pytest.approx(t.pixel_to_geo(101.0, 50.0))
+        assert pl.geo_stale is False
+
+    def test_unflagged_large_drift_refreshed(self):
+        t = FakeTransformer()
+        pl = make_synced_polyline(t, [(0.0, 0.0), (100.0, 50.0)])
+        pl.points[0] = (10.0, 0.0)  # 10px, beyond 2px threshold
+        assert refresh_polyline_geo_points(pl, t) == 1
+        assert pl.geo_points[0] == pytest.approx(t.pixel_to_geo(10.0, 0.0))
+
+    def test_length_mismatch_rebuilds_geo(self):
+        """H4: mismatched geo_points are rebuilt from pixels, not skipped."""
+        t = FakeTransformer()
+        pl = make_synced_polyline(t, [(0.0, 0.0), (100.0, 50.0)])
+        pl.points.append((200.0, 100.0))  # pixel added without geo entry
+        assert refresh_polyline_geo_points(pl, t) == 3
+        assert len(pl.geo_points) == len(pl.points)
+        assert pl.geo_points[2] == pytest.approx(t.pixel_to_geo(200.0, 100.0))
+
+    def test_edited_only_skips_unflagged_drift(self):
+        """In edited_only mode (view switch) geo stays authoritative for
+        unflagged drift, e.g. from adjustment/transformer changes."""
+        t = FakeTransformer()
+        pl = make_synced_polyline(t, [(0.0, 0.0), (100.0, 50.0)])
+        pl.points[0] = (10.0, 0.0)  # large drift, but no edit flag
+        original_geo = list(pl.geo_points)
+        assert refresh_polyline_geo_points(pl, t, edited_only=True) == 0
+        assert pl.geo_points == original_geo
+
+    def test_edited_only_still_refreshes_flagged(self):
+        t = FakeTransformer()
+        pl = make_synced_polyline(t, [(0.0, 0.0), (100.0, 50.0)])
+        pl.update_point(0, 10.0, 0.0)
+        assert refresh_polyline_geo_points(pl, t, edited_only=True) == 1
+        assert pl.geo_points[0] == pytest.approx(t.pixel_to_geo(10.0, 0.0))
+
+    def test_no_geo_points_noop(self):
+        pl = Polyline(id='p', points=[(0.0, 0.0)])
+        pl.geo_stale = True
+        assert refresh_polyline_geo_points(pl, FakeTransformer()) == 0
+        assert pl.geo_points is None
+        assert pl.geo_stale is False
+
+
+class TestRefreshProject:
+    """refresh_stale_geo_points over a whole project."""
+
+    def test_connecting_road_mismatch_rebuilt(self):
+        t = FakeTransformer()
+        project = Project()
+        cr = Road(id='cr1', junction_id='j1')
+        cr.inline_path = [(0.0, 0.0), (10.0, 10.0), (20.0, 20.0)]
+        cr.inline_geo_path = [t.pixel_to_geo(0.0, 0.0)]  # wrong length
+        project.roads.append(cr)
+        refreshed = refresh_stale_geo_points(project, t)
+        assert refreshed == 3
+        assert len(cr.inline_geo_path) == 3
+        assert cr.inline_geo_path[2] == pytest.approx(t.pixel_to_geo(20.0, 20.0))
+
+    def test_none_transformer_noop(self):
+        assert refresh_stale_geo_points(Project(), None) == 0
+
+    def test_endpoint_snap_marks_stale(self):
+        """Project.enforce_road_links endpoint snapping flags geo as stale."""
+        t = FakeTransformer()
+        pl = make_synced_polyline(t, [(0.0, 0.0), (100.0, 50.0)])
+        pl.update_point(0, 1.0, 1.0)
+        assert pl.geo_stale is True
+
+
+class TestPolylineToMetricPoints:
+    """Geo used only when consistent; pixel fallback on mismatch."""
+
+    def test_uses_geo_when_lengths_match(self):
+        t = FakeTransformer()
+        pl = make_synced_polyline(t, [(0.0, 0.0), (100.0, 50.0)])
+        expected = [t.latlon_to_meters(lat, lon) for lon, lat in pl.geo_points]
+        assert polyline_to_metric_points(pl, t) == pytest.approx(expected)
+
+    def test_falls_back_to_pixels_on_mismatch(self):
+        t = FakeTransformer()
+        pl = make_synced_polyline(t, [(0.0, 0.0), (100.0, 50.0)])
+        pl.points.append((200.0, 100.0))
+        result = polyline_to_metric_points(pl, t)
+        assert len(result) == 3
+        assert result == pytest.approx(t.pixels_to_meters_batch(pl.points))
+
+
+class TestReprojectPreservesEdits:
+    """H2: pixel edits survive a view switch even when geo was not synced."""
+
+    def test_unsynced_edit_survives_reprojection(self):
+        from orbit.utils.reproject import reproject_project_geometry
+
+        old_t = FakeTransformer(offset_px=0.0)
+        new_t = FakeTransformer(offset_px=500.0)
+        project = Project()
+        pl = make_synced_polyline(old_t, [(0.0, 0.0), (100.0, 50.0)])
+        pl.update_point(1, 110.0, 50.0)  # edit without geo sync
+        project.polylines.append(pl)
+
+        reproject_project_geometry(project, old_t, new_t)
+
+        # The edited point must land at the new-space equivalent of the
+        # EDITED pixel position, not the pre-edit one.
+        expected = new_t.geo_to_pixel(*old_t.pixel_to_geo(110.0, 50.0))
+        assert pl.points[1] == pytest.approx(expected)
diff --git a/tests/unit/test_utils/test_transformer_from_bounds.py b/tests/unit/test_utils/test_transformer_from_bounds.py
index 65cdae5..0e6ef9d 100644
--- a/tests/unit/test_utils/test_transformer_from_bounds.py
+++ b/tests/unit/test_utils/test_transformer_from_bounds.py
@@ -1,5 +1,7 @@
 """Tests for create_transformer_from_bounds and geometry re-projection."""
 
+import math
+
 import pytest
 
 from orbit.models.polyline import Polyline
@@ -7,8 +9,16 @@
 from orbit.utils.coordinate_transform import create_transformer_from_bounds
 
 
+def merc(lat_deg: float) -> float:
+    return math.log(math.tan(math.pi / 4 + math.radians(lat_deg) / 2))
+
+
+def inv_merc(m: float) -> float:
+    return math.degrees(2 * math.atan(math.exp(m)) - math.pi / 2)
+
+
 class TestCreateTransformerFromBounds:
-    """Tests for creating an affine transformer from geographic bounds."""
+    """Tests for the Web-Mercator transformer created from geographic bounds."""
 
     def test_corners_roundtrip(self):
         """Image corners should map to the expected geographic coordinates."""
@@ -21,24 +31,29 @@ def test_corners_roundtrip(self):
 
         # Top-left pixel → NW corner
         lon, lat = t.pixel_to_geo(0, 0)
-        assert lon == pytest.approx(11.0, abs=0.01)
-        assert lat == pytest.approx(58.0, abs=0.01)
+        assert lon == pytest.approx(11.0, abs=1e-9)
+        assert lat == pytest.approx(58.0, abs=1e-9)
 
         # Bottom-right pixel → SE corner
         lon, lat = t.pixel_to_geo(1000, 800)
-        assert lon == pytest.approx(12.0, abs=0.01)
-        assert lat == pytest.approx(57.0, abs=0.01)
+        assert lon == pytest.approx(12.0, abs=1e-9)
+        assert lat == pytest.approx(57.0, abs=1e-9)
 
-    def test_center_point(self):
-        """Image center should map to geographic center."""
+    def test_center_point_is_mercator_midpoint(self):
+        """Tile imagery is linear in Mercator y: the image center maps to the
+        Mercator midpoint, which lies slightly north of the latitude midpoint
+        in the northern hemisphere (H1)."""
         t = create_transformer_from_bounds(
             2000, 1000,
             min_lon=10.0, min_lat=55.0,
             max_lon=12.0, max_lat=57.0,
         )
         lon, lat = t.pixel_to_geo(1000, 500)
-        assert lon == pytest.approx(11.0, abs=0.01)
-        assert lat == pytest.approx(56.0, abs=0.01)
+        assert lon == pytest.approx(11.0, abs=1e-9)
+        expected_lat = inv_merc((merc(55.0) + merc(57.0)) / 2)
+        assert lat == pytest.approx(expected_lat, abs=1e-9)
+        # ~0.013° north of the linear-lat midpoint at this latitude/extent
+        assert lat > 56.0
 
     def test_geo_to_pixel_roundtrip(self):
         """pixel_to_geo → geo_to_pixel should round-trip."""
@@ -50,8 +65,27 @@ def test_geo_to_pixel_roundtrip(self):
         for px, py in [(0, 0), (400, 300), (800, 600), (200, 100)]:
             lon, lat = t.pixel_to_geo(px, py)
             px2, py2 = t.geo_to_pixel(lon, lat)
-            assert px2 == pytest.approx(px, abs=0.5)
-            assert py2 == pytest.approx(py, abs=0.5)
+            assert px2 == pytest.approx(px, abs=1e-6)
+            assert py2 == pytest.approx(py, abs=1e-6)
+
+    def test_scale_factors_isotropic_for_square_mercator_tiles(self):
+        """Slippy tiles are square in Mercator: when the pixel grid matches
+        (lon span / W == merc span / H), x and y scales must be equal."""
+        min_lat, max_lat = 57.0, 57.5
+        merc_span = merc(max_lat) - merc(min_lat)
+        lon_span = math.degrees(merc_span)  # square mercator pixels
+        width, height = 1000, 1000
+        t = create_transformer_from_bounds(
+            width, height, 11.0, min_lat, 11.0 + lon_span, max_lat)
+        scale_x, scale_y = t.get_scale_factor()
+        assert scale_x == pytest.approx(scale_y, rel=1e-6)
+        assert scale_x > 0
+
+    def test_invalid_bounds_return_none(self):
+        assert create_transformer_from_bounds(0, 100, 11.0, 57.0, 12.0, 58.0) is None
+        assert create_transformer_from_bounds(100, 100, 12.0, 57.0, 11.0, 58.0) is None
+        assert create_transformer_from_bounds(100, 100, 11.0, 58.0, 12.0, 57.0) is None
+        assert create_transformer_from_bounds(100, 100, 11.0, 57.0, 12.0, 90.0) is None
 
 
 class TestReprojectGeometry:
@@ -61,7 +95,7 @@ def test_polyline_with_geo_points(self):
         """Polyline with geo_points should re-project correctly."""
         from orbit.utils.reproject import reproject_project_geometry
 
-        # Create two transformers with different pixel spaces
+        # Same bounds at double resolution: positions scale exactly 2x
         t_old = create_transformer_from_bounds(
             1000, 800, 11.0, 57.0, 12.0, 58.0,
         )
@@ -69,19 +103,20 @@ def test_polyline_with_geo_points(self):
             2000, 1600, 11.0, 57.0, 12.0, 58.0,
         )
 
+        geo = (11.5, 57.5)
         project = Project()
         poly = Polyline(
             id="test",
-            points=[(500, 400)],
-            geo_points=[(11.5, 57.5)],
+            points=[t_old.geo_to_pixel(*geo)],
+            geo_points=[geo],
         )
         project.polylines.append(poly)
 
         reproject_project_geometry(project, t_old, t_new)
 
-        # In the new image (2x size), center should be at (1000, 800)
-        assert poly.points[0][0] == pytest.approx(1000.0, abs=2.0)
-        assert poly.points[0][1] == pytest.approx(800.0, abs=2.0)
+        old_x, old_y = t_old.geo_to_pixel(*geo)
+        assert poly.points[0][0] == pytest.approx(2 * old_x, abs=1e-6)
+        assert poly.points[0][1] == pytest.approx(2 * old_y, abs=1e-6)
 
     def test_polyline_pixel_only(self):
         """Polyline with only pixel coords gets geo_points created."""
@@ -103,9 +138,9 @@ def test_polyline_pixel_only(self):
         # Should now have geo_points
         assert poly.geo_points is not None
         assert len(poly.geo_points) == 1
-        # And pixel coords should be re-projected (halved image → halved coords)
-        assert poly.points[0][0] == pytest.approx(250.0, abs=2.0)
-        assert poly.points[0][1] == pytest.approx(200.0, abs=2.0)
+        # Same bounds at half resolution → exactly halved coords
+        assert poly.points[0][0] == pytest.approx(250.0, abs=1e-6)
+        assert poly.points[0][1] == pytest.approx(200.0, abs=1e-6)
 
     def test_control_points_reprojected(self):
         """Control point pixel positions should be updated."""
@@ -118,9 +153,10 @@ def test_control_points_reprojected(self):
             2000, 1600, 11.0, 57.0, 12.0, 58.0,
         )
 
+        old_px, old_py = t_old.geo_to_pixel(11.5, 57.5)
         project = Project()
         cp = ControlPoint(
-            pixel_x=500, pixel_y=400,
+            pixel_x=old_px, pixel_y=old_py,
             longitude=11.5, latitude=57.5,
             name="test",
         )
@@ -128,8 +164,8 @@ def test_control_points_reprojected(self):
 
         reproject_project_geometry(project, t_old, t_new)
 
-        assert cp.pixel_x == pytest.approx(1000.0, abs=2.0)
-        assert cp.pixel_y == pytest.approx(800.0, abs=2.0)
+        assert cp.pixel_x == pytest.approx(2 * old_px, abs=1e-6)
+        assert cp.pixel_y == pytest.approx(2 * old_py, abs=1e-6)
         # Geographic coords unchanged
         assert cp.longitude == 11.5
         assert cp.latitude == 57.5
diff --git a/uv.lock b/uv.lock
index d526eac..d24de0f 100644
--- a/uv.lock
+++ b/uv.lock
@@ -524,7 +524,7 @@ wheels = [
 
 [[package]]
 name = "orbit"
-version = "0.10.1"
+version = "0.11.0"
 source = { editable = "." }
 dependencies = [
     { name = "geomag" },