fix: backfill missing sidewalk corners

scripts/build-osm2streets-qgis.js

@@ -458,6 +458,7 @@ function splitLayers(dir, arrowScaleValue, maxCornerDimensionMeters, osm) {
   const lanePolygons = JSON.parse(fs.readFileSync(path.join(dir, "lane_polygons.geojson"), "utf8"));
   const markings = JSON.parse(fs.readFileSync(path.join(dir, "lane_markings.geojson"), "utf8"));
   const intersections = JSON.parse(fs.readFileSync(path.join(dir, "intersection_markings.geojson"), "utf8"));
+  const network = JSON.parse(fs.readFileSync(path.join(dir, "network.json"), "utf8"));
   const crosswalkData = buildCrosswalks(osm);
   const serviceWayIds = new Set([...osm.ways.values()]
     .filter((way) => way.tags.highway === "service")
@@ -470,7 +471,7 @@ function splitLayers(dir, arrowScaleValue, maxCornerDimensionMeters, osm) {
     centerLines: emptyCollection(),
     vehicleStopLines: crosswalkData.stopLines,
     laneArrows: emptyCollection(),
-    sidewalkCorners: filteredSidewalkCorners(intersections, maxCornerDimensionMeters),
+    sidewalkCorners: buildSidewalkCorners(intersections, plain, network, maxCornerDimensionMeters),
     crosswalks: crosswalkData.stripes,
   };
   const serviceDrivingPolygons = [];
@@ -527,6 +528,389 @@ function filteredSidewalkCorners(intersections, maxDimensionMeters) {
   return out;
 }
 
+function buildSidewalkCorners(intersections, plain, network, maxDimensionMeters) {
+  const out = filteredSidewalkCorners(intersections, maxDimensionMeters);
+  const fallback = synthesizeMissingSidewalkCorners(out, plain, network, maxDimensionMeters);
+  out.features.push(...fallback);
+  return out;
+}
+
+function synthesizeMissingSidewalkCorners(existing, plain, network, maxDimensionMeters) {
+  const roadFeatures = new Map((plain.features || [])
+    .filter((feature) => feature.properties?.type === "road")
+    .map((feature) => [Number(feature.properties.id), feature]));
+  const intersectionFeatures = new Map((plain.features || [])
+    .filter((feature) => feature.properties?.type === "intersection")
+    .map((feature) => [Number(feature.properties.id), feature]));
+  const roads = new Map((network.roads || []).map(([id, road]) => [Number(id), road]));
+  const intersections = new Map((network.intersections || []).map(([id, intersection]) => [Number(id), intersection]));
+  const existingByIntersection = assignCornersToIntersections(existing.features || [], intersectionFeatures);
+  const synthesized = [];
+
+  for (const [intersectionId, intersection] of intersections.entries()) {
+    const intersectionFeature = intersectionFeatures.get(intersectionId);
+    if (!intersectionFeature || intersection.roads.length < 2) continue;
+
+    const edges = buildIntersectionEdges(intersection, roads, roadFeatures, intersectionFeature);
+    if (!edges.length) continue;
+    const qualifyingPairs = qualifyingCornerPairs(edges);
+    if (qualifyingPairs.length < 3) continue;
+
+    const current = existingByIntersection.get(intersectionId) || [];
+    if (qualifyingPairs.length - current.length !== 1) continue;
+    const currentCenters = current.map((entry) => entry.center);
+    const targetDimension = median(current.map((entry) => maxFeatureDimensionMeters(entry.feature)).filter(Number.isFinite));
+    const candidates = [];
+
+    for (const [one, two] of qualifyingPairs) {
+      const candidate = synthesizeCornerFeature(one, two, intersectionFeature, maxDimensionMeters);
+      if (!candidate) continue;
+
+      const candidateCenter = featureCenter(candidate);
+      if (!candidateCenter || !pointInPolygon(candidateCenter, intersectionFeature.geometry.coordinates)) continue;
+      const dimension = maxFeatureDimensionMeters(candidate);
+      if (dimension === null || dimension > maxDimensionMeters) continue;
+      if (polygonAreaMeters2(candidate) < 0.4) continue;
+      if (currentCenters.some((point) => pointDistance(point, candidateCenter) <= 0.6)) continue;
+      candidates.push({
+        feature: candidate,
+        center: candidateCenter,
+        dimension,
+        score: Math.abs(dimension - targetDimension),
+      });
+    }
+
+    if (!candidates.length) continue;
+    candidates.sort((a, b) => a.score - b.score || a.dimension - b.dimension);
+    synthesized.push(candidates[0].feature);
+  }
+
+  return synthesized;
+}
+
+function assignCornersToIntersections(features, intersectionFeatures) {
+  const out = new Map();
+  for (const feature of features) {
+    const point = featureCenter(feature);
+    if (!point) continue;
+    for (const [intersectionId, intersectionFeature] of intersectionFeatures.entries()) {
+      if (!pointInPolygon(point, intersectionFeature.geometry.coordinates)) continue;
+      const bucket = out.get(intersectionId) || [];
+      bucket.push({
+        feature,
+        center: point,
+      });
+      out.set(intersectionId, bucket);
+      break;
+    }
+  }
+  return out;
+}
+
+function buildIntersectionEdges(intersection, roads, roadFeatures, intersectionFeature) {
+  const edges = [];
+  for (const roadId of intersection.roads || []) {
+    const road = roads.get(Number(roadId));
+    const roadFeature = roadFeatures.get(Number(roadId));
+    if (!road || !roadFeature) return [];
+    const geometry = roadEndpointGeometry(road, roadFeature, intersectionFeature, intersection.id);
+    if (!geometry) return [];
+    const first = road.dst_i === intersection.id
+      ? makeRoadEdge(road, geometry, "right")
+      : makeRoadEdge(road, geometry, "left");
+    const second = road.dst_i === intersection.id
+      ? makeRoadEdge(road, geometry, "left")
+      : makeRoadEdge(road, geometry, "right");
+    if (!first || !second) return [];
+    edges.push(first, second);
+  }
+  return edges;
+}
+
+function qualifyingCornerPairs(edges) {
+  if (!edges.length) return [];
+  const loop = [...edges, edges[0]];
+  const pairs = [];
+  for (let i = 0; i < loop.length - 1; i += 1) {
+    const one = loop[i];
+    const two = loop[i + 1];
+    if (one.roadId === two.roadId) continue;
+    if (!isWalkableOuterLane(one.laneType) || !isWalkableOuterLane(two.laneType)) continue;
+    if (one.laneCount === 1 || two.laneCount === 1) continue;
+    pairs.push([one, two]);
+  }
+  return pairs;
+}
+
+function isWalkableOuterLane(type) {
+  return type === "Sidewalk" || type === "Shoulder";
+}
+
+function roadEndpointGeometry(road, roadFeature, intersectionFeature, intersectionId) {
+  const ring = normalizedRing(roadFeature.geometry.coordinates?.[0]);
+  if (ring.length !== 4) return null;
+
+  const shortEdges = shortEdgePairs(ring);
+  if (!shortEdges) return null;
+  const intersectionCenter = ringCenter(intersectionFeature.geometry.coordinates[0]);
+  const candidates = shortEdges.map(([a, b]) => {
+    const near = [ring[a], ring[b]];
+    return {
+      pair: [a, b],
+      center: midpoint(near[0], near[1]),
+      distance: pointDistanceMeters(midpoint(near[0], near[1]), intersectionCenter, metersForLat(intersectionCenter[1])),
+    };
+  });
+  candidates.sort((a, b) => a.distance - b.distance);
+  const nearPair = candidates[0].pair;
+  const farPair = candidates[1].pair;
+
+  const nearPoints = nearPair.map((idx) => ring[idx]);
+  const farPoints = farPair.map((idx) => ring[idx]);
+  const nearCenter = midpoint(nearPoints[0], nearPoints[1]);
+  const farCenter = midpoint(farPoints[0], farPoints[1]);
+  const roadDirection = road.src_i === intersectionId
+    ? normalizeLonLatVector([farCenter[0] - nearCenter[0], farCenter[1] - nearCenter[1]], nearCenter[1])
+    : normalizeLonLatVector([nearCenter[0] - farCenter[0], nearCenter[1] - farCenter[1]], nearCenter[1]);
+  if (!roadDirection) return null;
+
+  const correspondences = nearPair.map((idx) => {
+    const farIdx = farPair.find((candidate) => circularIndexDistance(idx, candidate, ring.length) === 1);
+    return farIdx === undefined ? null : [ring[idx], ring[farIdx]];
+  });
+  if (correspondences.some((pair) => !pair)) return null;
+
+  const classified = correspondences.map(([nearPoint, farPoint]) => ({
+    near: nearPoint,
+    far: farPoint,
+    cross: signedSide(roadDirection, nearCenter, nearPoint, nearCenter[1]),
+  })).sort((a, b) => a.cross - b.cross);
+
+  return {
+    nearCenter,
+    nearLeft: classified[1].near,
+    farLeft: classified[1].far,
+    nearRight: classified[0].near,
+    farRight: classified[0].far,
+  };
+}
+
+function shortEdgePairs(ring) {
+  const lengths = ring.map((point, index) => lineLengthMeters(point, ring[(index + 1) % ring.length]));
+  const optionA = lengths[0] + lengths[2];
+  const optionB = lengths[1] + lengths[3];
+  if (!Number.isFinite(optionA) || !Number.isFinite(optionB)) return null;
+  return optionA <= optionB
+    ? [[0, 1], [2, 3]]
+    : [[1, 2], [3, 0]];
+}
+
+function circularIndexDistance(a, b, size) {
+  const distance = Math.abs(a - b);
+  return Math.min(distance, size - distance);
+}
+
+function makeRoadEdge(road, geometry, side) {
+  const lane = side === "left"
+    ? road.lane_specs_ltr?.[0]
+    : road.lane_specs_ltr?.[road.lane_specs_ltr.length - 1];
+  if (!lane) return null;
+  const outerNear = side === "left" ? geometry.nearLeft : geometry.nearRight;
+  const outerFar = side === "left" ? geometry.farLeft : geometry.farRight;
+  const oppositeNear = side === "left" ? geometry.nearRight : geometry.nearLeft;
+  const oppositeFar = side === "left" ? geometry.farRight : geometry.farLeft;
+  const widthMeters = Number(lane.width) / 10000;
+  const innerNear = moveTowards(outerNear, oppositeNear, widthMeters);
+  const innerFar = moveTowards(outerFar, oppositeFar, widthMeters);
+  return {
+    roadId: road.id,
+    laneType: lane.lt,
+    laneCount: road.lane_specs_ltr?.length || 0,
+    outerNear,
+    innerNear,
+    innerFar,
+  };
+}
+
+function synthesizeCornerFeature(one, two, intersectionFeature, maxDimensionMeters) {
+  const ring = normalizedRing(intersectionFeature.geometry.coordinates?.[0]);
+  const slice = shorterRingSliceBetween(ring, one.outerNear, two.outerNear);
+  if (!slice || slice.length < 2) return null;
+
+  const meetPoint = lineIntersection(one.innerFar, one.innerNear, two.innerFar, two.innerNear);
+  const points = dedupeSequentialPoints([
+    ...slice,
+    two.innerNear,
+    ...(meetPoint && pointInPolygon(meetPoint, intersectionFeature.geometry.coordinates) ? [meetPoint] : []),
+    one.innerNear,
+    slice[0],
+  ]);
+  if (points.length < 4) return null;
+
+  const feature = {
+    type: "Feature",
+    properties: {
+      type: "sidewalk corner",
+      source: "fallback",
+    },
+    geometry: {
+      type: "Polygon",
+      coordinates: [points],
+    },
+  };
+  const dimension = maxFeatureDimensionMeters(feature);
+  if (dimension === null || dimension > maxDimensionMeters) return null;
+  return feature;
+}
+
+function shorterRingSliceBetween(ring, start, end) {
+  if (!ring.length) return null;
+  const startIndex = nearestRingPointIndex(ring, start, 0.8);
+  const endIndex = nearestRingPointIndex(ring, end, 0.8);
+  if (startIndex === null || endIndex === null) return null;
+  if (startIndex === endIndex) return [ring[startIndex]];
+  const forward = walkRing(ring, startIndex, endIndex, 1);
+  const backward = walkRing(ring, startIndex, endIndex, -1);
+  return pathLengthMeters(forward) <= pathLengthMeters(backward) ? forward : backward;
+}
+
+function walkRing(ring, startIndex, endIndex, direction) {
+  const out = [ring[startIndex]];
+  let index = startIndex;
+  while (index !== endIndex) {
+    index = (index + direction + ring.length) % ring.length;
+    out.push(ring[index]);
+  }
+  return out;
+}
+
+function pathLengthMeters(points) {
+  let total = 0;
+  for (let i = 1; i < points.length; i += 1) total += lineLengthMeters(points[i - 1], points[i]);
+  return total;
+}
+
+function lineIntersection(a1, a2, b1, b2) {
+  const originLat = (a1[1] + a2[1] + b1[1] + b2[1]) / 4;
+  const meters = metersForLat(originLat);
+  const ax1 = 0;
+  const ay1 = 0;
+  const ax2 = (a2[0] - a1[0]) * meters.lon;
+  const ay2 = (a2[1] - a1[1]) * meters.lat;
+  const bx1 = (b1[0] - a1[0]) * meters.lon;
+  const by1 = (b1[1] - a1[1]) * meters.lat;
+  const bx2 = (b2[0] - a1[0]) * meters.lon;
+  const by2 = (b2[1] - a1[1]) * meters.lat;
+  const denominator = (ax2 - ax1) * (by2 - by1) - (ay2 - ay1) * (bx2 - bx1);
+  if (Math.abs(denominator) < 1e-9) return null;
+  const ua = ((bx2 - bx1) * (ay1 - by1) - (by2 - by1) * (ax1 - bx1)) / denominator;
+  return [
+    a1[0] + ((ax1 + ua * (ax2 - ax1)) / meters.lon),
+    a1[1] + ((ay1 + ua * (ay2 - ay1)) / meters.lat),
+  ];
+}
+
+function normalizedRing(ring) {
+  if (!Array.isArray(ring) || ring.length < 4) return [];
+  const out = ring.map((point) => [point[0], point[1]]);
+  if (pointDistance(out[0], out[out.length - 1]) <= 0.02) out.pop();
+  return out;
+}
+
+function nearestRingPointIndex(ring, point, maxDistanceMeters) {
+  let bestIndex = null;
+  let bestDistance = Infinity;
+  for (let i = 0; i < ring.length; i += 1) {
+    const distance = pointDistance(ring[i], point);
+    if (distance < bestDistance) {
+      bestDistance = distance;
+      bestIndex = i;
+    }
+  }
+  return bestDistance <= maxDistanceMeters ? bestIndex : null;
+}
+
+function nearestRingPoint(ring, point, maxDistanceMeters) {
+  const normalized = normalizedRing(ring);
+  const index = nearestRingPointIndex(normalized, point, maxDistanceMeters);
+  return index === null ? null : normalized[index];
+}
+
+function dedupeSequentialPoints(points, toleranceMeters = 0.02) {
+  const out = [];
+  for (const point of points) {
+    if (!out.length || pointDistance(out[out.length - 1], point) > toleranceMeters) out.push(point);
+  }
+  if (out.length >= 2 && pointDistance(out[0], out[out.length - 1]) > toleranceMeters) out.push(out[0]);
+  return out;
+}
+
+function dedupePointList(points, toleranceMeters) {
+  const out = [];
+  for (const point of points) {
+    if (out.some((other) => pointDistance(point, other) <= toleranceMeters)) continue;
+    out.push(point);
+  }
+  return out;
+}
+
+function midpoint(a, b) {
+  return [(a[0] + b[0]) / 2, (a[1] + b[1]) / 2];
+}
+
+function moveTowards(from, to, distanceMeters) {
+  const meters = metersForLat((from[1] + to[1]) / 2);
+  const unit = normalizeMetersVector([to[0] - from[0], to[1] - from[1]], meters);
+  if (!unit) return from;
+  return addMeters(from, unit, distanceMeters, meters);
+}
+
+function normalizeLonLatVector([dxLon, dyLat], lat) {
+  return normalizeMetersVector([dxLon, dyLat], metersForLat(lat));
+}
+
+function signedSide(direction, origin, point, lat) {
+  const meters = metersForLat(lat);
+  const dx = (point[0] - origin[0]) * meters.lon;
+  const dy = (point[1] - origin[1]) * meters.lat;
+  return direction[0] * dy - direction[1] * dx;
+}
+
+function pointDistance(a, b) {
+  return lineLengthMeters(a, b);
+}
+
+function lineLengthMeters(a, b) {
+  const meters = metersForLat((a[1] + b[1]) / 2);
+  return Math.hypot((a[0] - b[0]) * meters.lon, (a[1] - b[1]) * meters.lat);
+}
+
+function ringCenter(ring) {
+  const points = normalizedRing(ring);
+  const xs = points.map((point) => point[0]);
+  const ys = points.map((point) => point[1]);
+  return [(Math.min(...xs) + Math.max(...xs)) / 2, (Math.min(...ys) + Math.max(...ys)) / 2];
+}
+
+function polygonAreaMeters2(feature) {
+  const ring = normalizedRing(feature.geometry?.coordinates?.[0]);
+  if (ring.length < 3) return 0;
+  const meters = metersForLat(ring.reduce((sum, point) => sum + point[1], 0) / ring.length);
+  let area = 0;
+  for (let i = 0; i < ring.length; i += 1) {
+    const a = ring[i];
+    const b = ring[(i + 1) % ring.length];
+    area += (a[0] * meters.lon) * (b[1] * meters.lat) - (b[0] * meters.lon) * (a[1] * meters.lat);
+  }
+  return Math.abs(area) / 2;
+}
+
+function median(values) {
+  if (!values.length) return 0;
+  const sorted = [...values].sort((a, b) => a - b);
+  return sorted[Math.floor(sorted.length / 2)];
+}
+
 function maxFeatureDimensionMeters(feature) {
   const points = [];
   collectCoords(feature.geometry?.coordinates, points);
@@ -599,6 +983,15 @@ function representativePoint(feature) {
   return coords[Math.floor(coords.length / 2)];
 }
 
+function featureCenter(feature) {
+  const coords = [];
+  collectCoords(feature.geometry?.coordinates, coords);
+  if (!coords.length) return null;
+  const xs = coords.map((point) => point[0]);
+  const ys = coords.map((point) => point[1]);
+  return [(Math.min(...xs) + Math.max(...xs)) / 2, (Math.min(...ys) + Math.max(...ys)) / 2];
+}
+
 function pointInPolygon(point, rings) {
   if (!rings?.length || !pointInRing(point, rings[0])) return false;
   return !rings.slice(1).some((ring) => pointInRing(point, ring));