feat: add border corridor analysis, gap-filler detection, PTZ camera geolocation, and 3D terrain/satellite visualization controls

2026-09-26 00:06:52 +03:00
parent cb91e5f0de
commit 5028d20881
5 changed files with 2713 additions and 252 deletions
@@ -50,6 +50,11 @@ export class ArtilleryMissionRequestDto {
  @IsOptional()
  @IsString()
  caliber?: string;
+
+  @ApiPropertyOptional({ description: 'Trajectory mode: auto, low, or high', default: 'auto' })
+  @IsOptional()
+  @IsString()
+  trajectoryMode?: 'auto' | 'low' | 'high';
 }

 export class TacticalSymbolDto {
@@ -300,7 +300,7 @@ export class TacticalService {
    observerHeight: number = 2,
    radiusMeters: number = 5000,
    rayCount: number = 72,
-    samplesPerRay: number = 25,
+    samplesPerRay: number = 60,
  ) {
    const R_earth = 6371000;
    const k_refraction = 0.13;
@@ -338,42 +338,66 @@ export class TacticalService {
    const centerElev = (await DemTileService.getElevation(centerLat, centerLng, 13)) ?? this.estimateElevation(centerLat, centerLng);
    const observerTotal = centerElev + observerHeight;

-    const polygonCoordinates: [number, number][] = [];
-    let totalVisibleSum = 0;
+    const multiPolygonCoords: [number, number][][][] = [];
+    const invisibleMultiPolygonCoords: [number, number][][][] = [];
+    let visibleAreaSum = 0;
+    let totalAreaSum = 0;
+
+    const raysVis: boolean[][] = [];
+    const raysPoints: {lat: number, lng: number}[][] = [];

    for (let r = 0; r < rayCount; r++) {
      const raySteps = rays[r];
+      const vis: boolean[] = [true];
+      const pts = [{lat: centerLat, lng: centerLng}];
      let maxTheta = -Infinity;
-      let visibleDist = radiusMeters;
-      let visibleLat = raySteps[raySteps.length - 1].lat;
-      let visibleLng = raySteps[raySteps.length - 1].lng;

      for (const step of raySteps) {
+        pts.push({lat: step.lat, lng: step.lng});
        const sElev = (await DemTileService.getElevation(step.lat, step.lng, 13)) ?? this.estimateElevation(step.lat, step.lng);
        const sagitta = (step.dist * step.dist) / (2 * effectiveRadius);
-        const apparentElev = sElev + sagitta;
+        // Earth curvature drops the apparent terrain elevation below observer tangent plane
+        const apparentElev = sElev - sagitta;
        const theta = (apparentElev - observerTotal) / step.dist;

        if (theta >= maxTheta) {
          maxTheta = theta;
-          visibleDist = step.dist;
-          visibleLat = step.lat;
-          visibleLng = step.lng;
+          vis.push(true);
+        } else {
+          vis.push(false);
+        }
+      }
+      raysVis.push(vis);
+      raysPoints.push(pts);
+    }
+
+    for (let r = 0; r < rayCount; r++) {
+      const nextR = (r + 1) % rayCount;
+      const vis1 = raysVis[r];
+      const vis2 = raysVis[nextR];
+      const pts1 = raysPoints[r];
+      const pts2 = raysPoints[nextR];
+
+      for (let s = 1; s <= samplesPerRay; s++) {
+        totalAreaSum += s;
+        const p1 = [pts1[s-1].lng, pts1[s-1].lat] as [number, number];
+        const p2 = [pts1[s].lng, pts1[s].lat] as [number, number];
+        const p3 = [pts2[s].lng, pts2[s].lat] as [number, number];
+        const p4 = [pts2[s-1].lng, pts2[s-1].lat] as [number, number];
+
+        // Cell is visible only if BOTH bounding radial rays have direct line-of-sight
+        if (vis1[s] && vis2[s]) {
+          visibleAreaSum += s;
+          multiPolygonCoords.push([[p1, p2, p3, p4, p1]]);
+        } else {
+          invisibleMultiPolygonCoords.push([[p1, p2, p3, p4, p1]]);
+        }
      }
    }

-      totalVisibleSum += visibleDist;
-      polygonCoordinates.push([Number(visibleLng.toFixed(6)), Number(visibleLat.toFixed(6))]);
-    }
-
-    if (polygonCoordinates.length > 0) {
-      polygonCoordinates.push(polygonCoordinates[0]);
-    }
-
-    const avgRadius = totalVisibleSum / rayCount;
-    const coveredAreaKm2 = Math.PI * Math.pow(avgRadius / 1000, 2);
+    const coveragePercent = Math.min(100, Math.round((visibleAreaSum / totalAreaSum) * 100));
    const maxAreaKm2 = Math.PI * Math.pow(radiusMeters / 1000, 2);
-    const coveragePercent = Math.min(100, Math.round((coveredAreaKm2 / maxAreaKm2) * 100));
+    const coveredAreaKm2 = maxAreaKm2 * (coveragePercent / 100);

    return {
      center: {
@@ -390,10 +414,18 @@ export class TacticalService {
        type: 'Feature',
        properties: { centerLat, centerLng, radiusMeters, coveragePercent, coveredAreaKm2: Math.round(coveredAreaKm2 * 10) / 10 },
        geometry: {
-          type: 'Polygon',
-          coordinates: [polygonCoordinates],
+          type: 'MultiPolygon',
+          coordinates: multiPolygonCoords,
        },
      },
+      invisiblePolygon: {
+        type: 'Feature',
+        properties: { centerLat, centerLng, radiusMeters },
+        geometry: {
+          type: 'MultiPolygon',
+          coordinates: invisibleMultiPolygonCoords,
+        },
+      }
    };
  }

@@ -428,15 +460,19 @@ export class TacticalService {
    const distanceMeters = this.haversineDistance(gunLat, gunLng, targetLat, targetLng);
    const azimuthDegrees = this.calculateBearing(gunLat, gunLng, targetLat, targetLng);

-    const gunGroundElev = this.estimateElevation(gunLat, gunLng);
-    const targetGroundElev = this.estimateElevation(targetLat, targetLng);
+    const gunGroundElev = (await DemTileService.getElevation(gunLat, gunLng, 13)) ?? this.estimateElevation(gunLat, gunLng);
+    const targetGroundElev = (await DemTileService.getElevation(targetLat, targetLng, 13)) ?? this.estimateElevation(targetLat, targetLng);

    const gunTotalElev = gunGroundElev + gunElevationOffset;
    const targetTotalElev = targetGroundElev + targetElevationOffset;
    const heightDiff = targetTotalElev - gunTotalElev;

    const g = 9.80665;
-    const v0 = muzzleVelocity;
+    const isMortar = caliber.toLowerCase().includes('mortar') || caliber.includes('هاون');
+    let v0 = muzzleVelocity;
+    if (isMortar && (v0 > 450 || !v0 || v0 === 827)) {
+      v0 = 240; // Standard 120mm mortar velocity (Charge 3 ~240 m/s)
+    }

    const term = Math.pow(v0, 4) - g * (g * Math.pow(distanceMeters, 2) + 2 * heightDiff * Math.pow(v0, 2));

@@ -451,7 +487,7 @@ export class TacticalService {
      highAngleRad = Math.atan((Math.pow(v0, 2) + sqrtTerm) / (g * distanceMeters));
    } else {
      lowAngleRad = (45 * Math.PI) / 180;
-      highAngleRad = (60 * Math.PI) / 180;
+      highAngleRad = (65 * Math.PI) / 180;
    }

    const lowAngleDeg = (lowAngleRad * 180) / Math.PI;
@@ -460,13 +496,16 @@ export class TacticalService {
    const highAngleDeg = (highAngleRad * 180) / Math.PI;
    const highAngleMils = highAngleDeg * (6400 / 360);

-    const timeOfFlightSeconds = distanceMeters / (v0 * Math.cos(lowAngleRad));
-    const apexHeightMeters = gunTotalElev + Math.pow(v0 * Math.sin(lowAngleRad), 2) / (2 * g);
+    const samples = 60;

-    const trajectoryPoints: any[] = [];
-    const samples = 50;
-    let hasCrestClearance = true;
-    let criticalObstacle: any = null;
+    // Helper to evaluate trajectory points, apex, and obstacle crest clearance
+    const evaluateTrajectory = async (angleRad: number) => {
+      const tof = distanceMeters / (v0 * Math.cos(angleRad));
+      const apex = gunTotalElev + Math.pow(v0 * Math.sin(angleRad), 2) / (2 * g);
+      const points: any[] = [];
+      let isClear = true;
+      let minClearance = Infinity;
+      let obstacle: any = null;

      for (let i = 0; i <= samples; i++) {
        const frac = i / samples;
@@ -474,17 +513,21 @@ export class TacticalService {
        const lat = gunLat + (targetLat - gunLat) * frac;
        const lng = gunLng + (targetLng - gunLng) * frac;

-      const t = frac * timeOfFlightSeconds;
-      const y = v0 * Math.sin(lowAngleRad) * t - 0.5 * g * Math.pow(t, 2);
+        const t = frac * tof;
+        const y = v0 * Math.sin(angleRad) * t - 0.5 * g * Math.pow(t, 2);
        const projectileAlt = gunTotalElev + y;

-      const terrainElev = this.estimateElevation(lat, lng);
+        const terrainElev = (await DemTileService.getElevation(lat, lng, 13)) ?? this.estimateElevation(lat, lng);
        const clearance = projectileAlt - terrainElev;

+        if (clearance < minClearance && i > 1 && i < samples) {
+          minClearance = clearance;
+        }
+
        if (clearance <= 0 && i > 1 && i < samples) {
-        hasCrestClearance = false;
-        if (!criticalObstacle || clearance < criticalObstacle.clearance) {
-          criticalObstacle = {
+          isClear = false;
+          if (!obstacle || clearance < obstacle.clearance) {
+            obstacle = {
              distanceMeters: Math.round(d),
              terrainElevMeters: Math.round(terrainElev),
              projectileAltMeters: Math.round(projectileAlt),
@@ -495,7 +538,7 @@ export class TacticalService {
          }
        }

-      trajectoryPoints.push({
+        points.push({
          distanceMeters: Math.round(d),
          lat,
          lng,
@@ -505,29 +548,67 @@ export class TacticalService {
        });
      }

+      return { tof, apex, points, isClear, minClearance, obstacle };
+    };
+
+    const lowEval = await evaluateTrajectory(lowAngleRad);
+    const highEval = await evaluateTrajectory(highAngleRad);
+
+    // Tactical trajectory selection: Mortars are always high-angle; howitzers auto-switch if low is blocked
+    let activeTrajectory: 'low' | 'high' = 'low';
+    let chosenEval = lowEval;
+
+    if (isMortar) {
+      activeTrajectory = 'high';
+      chosenEval = highEval;
+    } else if (dto.trajectoryMode === 'high') {
+      activeTrajectory = 'high';
+      chosenEval = highEval;
+    } else if (dto.trajectoryMode === 'low') {
+      activeTrajectory = 'low';
+      chosenEval = lowEval;
+    } else {
+      if (!lowEval.isClear && highEval.isClear) {
+        activeTrajectory = 'high';
+        chosenEval = highEval;
+      } else {
+        activeTrajectory = 'low';
+        chosenEval = lowEval;
+      }
+    }
+
    return {
      fireMissionId: `FM-${Date.now().toString().slice(-6)}`,
      caliber,
-      muzzleVelocityMs: muzzleVelocity,
+      muzzleVelocityMs: v0,
      distanceMeters: Math.round(distanceMeters),
      distanceKm: Math.round((distanceMeters / 1000) * 100) / 100,
      azimuthDegrees: Math.round(azimuthDegrees * 10) / 10,
      azimuthMils: Math.round((azimuthDegrees * (6400 / 360)) * 10) / 10,
      gunElevationMeters: Math.round(gunTotalElev),
      targetElevationMeters: Math.round(targetTotalElev),
-      apexAltitudeMeters: Math.round(apexHeightMeters),
-      timeOfFlightSeconds: Math.round(timeOfFlightSeconds * 10) / 10,
+      apexAltitudeMeters: Math.round(chosenEval.apex),
+      timeOfFlightSeconds: Math.round(chosenEval.tof * 10) / 10,
+      activeTrajectory,
+      isMortar,
+      hasCrestClearance: chosenEval.isClear,
+      clearanceMarginMeters: Math.round(chosenEval.minClearance),
+      criticalObstacle: chosenEval.obstacle,
+      trajectoryPoints: chosenEval.points,
      lowAngle: {
        degrees: Math.round(lowAngleDeg * 100) / 100,
        mils: Math.round(lowAngleMils * 10) / 10,
+        hasClearance: lowEval.isClear,
+        apexAltitudeMeters: Math.round(lowEval.apex),
+        timeOfFlightSeconds: Math.round(lowEval.tof * 10) / 10,
      },
      highAngle: {
        degrees: Math.round(highAngleDeg * 100) / 100,
        mils: Math.round(highAngleMils * 10) / 10,
+        hasClearance: highEval.isClear,
+        apexAltitudeMeters: Math.round(highEval.apex),
+        timeOfFlightSeconds: Math.round(highEval.tof * 10) / 10,
      },
-      hasCrestClearance,
-      criticalObstacle,
-      trajectoryPoints,
    };
  }

@@ -84,6 +84,25 @@
        "https://tiles.intaleqapp.com/places_iraq/{z}/{x}/{y}"
      ],
      "maxzoom": 14
+    },
+    "terrain-dem": {
+      "type": "raster-dem",
+      "tiles": [
+        "https://s3.amazonaws.com/elevation-tiles-prod/terrarium/{z}/{x}/{y}.png",
+        "https://tiles.intaleqapp.com/raster_dem/{z}/{x}/{y}.png"
+      ],
+      "encoding": "terrarium",
+      "tileSize": 256,
+      "maxzoom": 15
+    },
+    "esri-satellite": {
+      "type": "raster",
+      "tiles": [
+        "https://server.arcgisonline.com/ArcGIS/rest/services/World_Imagery/MapServer/tile/{z}/{y}/{x}"
+      ],
+      "tileSize": 256,
+      "maxzoom": 19,
+      "attribution": "© Esri"
    }
  },
  "layers": [
@@ -94,6 +113,34 @@
        "background-color": "#F6F4F0"
      }
    },
+    {
+      "id": "satellite-base-layer",
+      "type": "raster",
+      "source": "esri-satellite",
+      "layout": {
+        "visibility": "none"
+      },
+      "paint": {
+        "raster-opacity": 0.95,
+        "raster-saturation": 0.1
+      }
+    },
+    {
+      "id": "terrain-3d-hillshading",
+      "type": "hillshade",
+      "source": "terrain-dem",
+      "layout": {
+        "visibility": "visible"
+      },
+      "paint": {
+        "hillshade-illumination-direction": 315,
+        "hillshade-illumination-anchor": "viewport",
+        "hillshade-shadow-color": "#261d15",
+        "hillshade-highlight-color": "#fffbf2",
+        "hillshade-accent-color": "#784a28",
+        "hillshade-exaggeration": 0.75
+      }
+    },
    {
      "id": "admin-boundary-national",
      "type": "line",
@@ -1362,7 +1409,7 @@
      "type": "fill-extrusion",
      "source": "local-osm-polygons",
      "source-layer": "planet_osm_polygon",
-      "minzoom": 13,
+      "minzoom": 12,
      "filter": [
        "has",
        "building"
@@ -1410,7 +1457,7 @@
      "type": "fill-extrusion",
      "source": "overture_buildings",
      "source-layer": "overture_building",
-      "minzoom": 13,
+      "minzoom": 12,
      "layout": {
        "visibility": "visible"
      },
@@ -470,6 +470,7 @@ export interface Viewshed360Result {
  centerElevation: number;
  radiusMeters: number;
  polygonGeoJson: any;
+  invisiblePolygonGeoJson?: any;
  totalRays: number;
  visibleAreaKm2: number;
  visiblePercentage: number;
@@ -503,8 +504,10 @@ export async function calculateRadialViewshed(
      return {
        centerLat,
        centerLng,
+        centerElevation: data.center?.totalElevation || obsHeight,
        radiusMeters,
        polygonGeoJson: data.polygon,
+        invisiblePolygonGeoJson: data.invisiblePolygon,
        totalRays: numRays,
        visibleAreaKm2: data.coveredAreaKm2,
        visiblePercentage: data.coveragePercent
@@ -521,7 +524,7 @@ export async function calculateRadialViewshed(
  const k_refraction = 0.13;
  const effectiveEarthRadius = R_earth / (1 - k_refraction);

-  const samplesPerRay = 20;
+  const samplesPerRay = 60;
  const rayCoords: Array<{ rayIdx: number; step: number; sDist: number; sLat: number; sLng: number }> = [];

  for (let rayIdx = 0; rayIdx < numRays; rayIdx++) {
@@ -546,17 +549,17 @@ export async function calculateRadialViewshed(
  // Batch sample all ray points
  const elevations = await sampleElevationsBatch(rayCoords.map(rc => ({ lat: rc.sLat, lng: rc.sLng })), 13);

-  const polygonCoordinates: [number, number][] = [];
-  let totalVisibleDistanceSum = 0;
+  const raysVis: boolean[][] = [];
+  const raysPoints: {sLat: number, sLng: number, sDist: number}[][] = [];

  for (let rayIdx = 0; rayIdx < numRays; rayIdx++) {
+    const vis: boolean[] = [true];
+    const pts = [{sLat: centerLat, sLng: centerLng, sDist: 0}];
    let maxAngleSoFar = -Infinity;
-    let visibleHorizonDist = radiusMeters;
-    let visibleHorizonLat = centerLat;
-    let visibleHorizonLng = centerLng;

    const rayPoints = rayCoords.filter(rc => rc.rayIdx === rayIdx);
    rayPoints.forEach(rc => {
+      pts.push({sLat: rc.sLat, sLng: rc.sLng, sDist: rc.sDist});
      const sElev = elevations[rc.rayIdx * samplesPerRay + (rc.step - 1)];
      const earthCurvatureDrop = (rc.sDist * rc.sDist) / (2 * effectiveEarthRadius);
      const apparentElev = sElev - earthCurvatureDrop;
@@ -564,18 +567,43 @@ export async function calculateRadialViewshed(

      if (angle >= maxAngleSoFar) {
        maxAngleSoFar = angle;
-        visibleHorizonDist = rc.sDist;
-        visibleHorizonLat = rc.sLat;
-        visibleHorizonLng = rc.sLng;
+        vis.push(true);
+      } else {
+        vis.push(false);
      }
    });
-
-    totalVisibleDistanceSum += visibleHorizonDist;
-    polygonCoordinates.push([visibleHorizonLng, visibleHorizonLat]);
+    raysVis.push(vis);
+    raysPoints.push(pts);
  }

-  if (polygonCoordinates.length > 0) {
-    polygonCoordinates.push(polygonCoordinates[0]);
+  const multiPolygonCoordinates: [number, number][][][] = [];
+  const invisibleMultiPolygonCoordinates: [number, number][][][] = [];
+  let visibleAreaSum = 0;
+  let totalAreaSum = 0;
+
+  for (let rayIdx = 0; rayIdx < numRays; rayIdx++) {
+    const nextRayIdx = (rayIdx + 1) % numRays;
+    const vis1 = raysVis[rayIdx];
+    const vis2 = raysVis[nextRayIdx];
+    const pts1 = raysPoints[rayIdx];
+    const pts2 = raysPoints[nextRayIdx];
+
+    for (let step = 1; step <= samplesPerRay; step++) {
+      totalAreaSum += step;
+      
+      const p1 = [pts1[step-1].sLng, pts1[step-1].sLat] as [number, number];
+      const p2 = [pts1[step].sLng, pts1[step].sLat] as [number, number];
+      const p3 = [pts2[step].sLng, pts2[step].sLat] as [number, number];
+      const p4 = [pts2[step-1].sLng, pts2[step-1].sLat] as [number, number];
+      
+      // Strict conjunction: cell is visible only if BOTH bounding rays have direct line of sight
+      if (vis1[step] && vis2[step]) {
+        visibleAreaSum += step;
+        multiPolygonCoordinates.push([[p1, p2, p3, p4, p1]]);
+      } else {
+        invisibleMultiPolygonCoordinates.push([[p1, p2, p3, p4, p1]]);
+      }
+    }
  }

  const polygonGeoJson = {
@@ -583,18 +611,32 @@ export async function calculateRadialViewshed(
    properties: {
      centerLat,
      centerLng,
-      radiusMeters
+      radiusMeters,
+      visibility: 'visible'
    },
    geometry: {
-      type: 'Polygon',
-      coordinates: [polygonCoordinates]
+      type: 'MultiPolygon',
+      coordinates: multiPolygonCoordinates
    }
  };

-  const avgVisibleDist = totalVisibleDistanceSum / numRays;
+  const invisiblePolygonGeoJson = {
+    type: 'Feature',
+    properties: {
+      centerLat,
+      centerLng,
+      radiusMeters,
+      visibility: 'invisible'
+    },
+    geometry: {
+      type: 'MultiPolygon',
+      coordinates: invisibleMultiPolygonCoordinates
+    }
+  };
+
+  const visiblePercentage = Math.min(100, Math.round((visibleAreaSum / totalAreaSum) * 100));
  const theoreticalMaxArea = Math.PI * Math.pow(radiusMeters / 1000, 2);
-  const actualVisibleArea = Math.PI * Math.pow(avgVisibleDist / 1000, 2);
-  const visiblePercentage = Math.min(100, Math.round((actualVisibleArea / theoreticalMaxArea) * 100));
+  const actualVisibleArea = theoreticalMaxArea * (visiblePercentage / 100);

  return {
    centerLat,
@@ -602,12 +644,161 @@ export async function calculateRadialViewshed(
    centerElevation: Math.round(observerElevation),
    radiusMeters,
    polygonGeoJson,
+    invisiblePolygonGeoJson,
    totalRays: numRays,
    visibleAreaKm2: Math.round(actualVisibleArea * 10) / 10,
    visiblePercentage
  };
 }

+export interface BorderCorridorResult {
+  postA: {
+    lat: number;
+    lng: number;
+    height: number;
+    elevation: number;
+    viewshed: Viewshed360Result;
+  };
+  postB: {
+    lat: number;
+    lng: number;
+    height: number;
+    elevation: number;
+    viewshed: Viewshed360Result;
+  };
+  distanceKm: number;
+  azimuthDeg: number;
+  combinedVisibleKm2: number;
+  combinedBlindKm2: number;
+  combinedCoveragePercent: number;
+  gapFiller?: {
+    lat: number;
+    lng: number;
+    groundElevation: number;
+    recommendedHeight: number;
+    totalElevation: number;
+    prominenceAboveValley: number;
+    mitigationPercent: number;
+    viewshed: Viewshed360Result;
+    tacticalRationale: string;
+  };
+}
+
+/**
+ * Dual-Post Border Corridor Viewshed & Topographic Gap-Filler Optimization Engine
+ */
+export async function calculateBorderCorridorAnalysis(
+  postALat: number,
+  postALng: number,
+  postBLat: number,
+  postBLng: number,
+  mastHeight: number = 20,
+  radiusMeters: number = 12000,
+  findGapFiller: boolean = true
+): Promise<BorderCorridorResult> {
+  const distMeters = calculateDistance(postALat, postALng, postBLat, postBLng);
+  const distanceKm = Math.round((distMeters / 1000) * 10) / 10;
+  const azimuthDeg = Math.round(calculateAzimuth(postALat, postALng, postBLat, postBLng));
+
+  // Run viewsheds concurrently for Post A and Post B
+  const [viewshedA, viewshedB] = await Promise.all([
+    calculateRadialViewshed(postALat, postALng, mastHeight, radiusMeters, 72),
+    calculateRadialViewshed(postBLat, postBLng, mastHeight, radiusMeters, 72)
+  ]);
+
+  const combinedVisibleKm2 = Math.round((viewshedA.visibleAreaKm2 + viewshedB.visibleAreaKm2) * 0.85 * 10) / 10;
+  const theoreticalTotalArea = Math.PI * Math.pow(radiusMeters / 1000, 2) * 1.7;
+  const combinedBlindKm2 = Math.max(0, Math.round((theoreticalTotalArea - combinedVisibleKm2) * 10) / 10);
+  const combinedCoveragePercent = Math.min(100, Math.round((combinedVisibleKm2 / theoreticalTotalArea) * 100));
+
+  let gapFillerData: BorderCorridorResult['gapFiller'] | undefined;
+
+  if (findGapFiller) {
+    const fractions = [0.25, 0.38, 0.5, 0.62, 0.75];
+    const lateralOffsetsMeters = [-3000, -1500, 0, 1500, 3000];
+    const radAz = (azimuthDeg * Math.PI) / 180;
+    const perpRad = radAz + Math.PI / 2;
+
+    const candidates: Array<{ lat: number; lng: number; frac: number; offset: number }> = [];
+    const R = 6371000;
+
+    for (const f of fractions) {
+      const midLat = postALat + (postBLat - postALat) * f;
+      const midLng = postALng + (postBLng - postALng) * f;
+
+      for (const off of lateralOffsetsMeters) {
+        const dLat = (off / R) * (180 / Math.PI) * Math.cos(perpRad);
+        const dLng = (off / (R * Math.cos((midLat * Math.PI) / 180))) * (180 / Math.PI) * Math.sin(perpRad);
+        candidates.push({
+          lat: Number((midLat + dLat).toFixed(5)),
+          lng: Number((midLng + dLng).toFixed(5)),
+          frac: f,
+          offset: off
+        });
+      }
+    }
+
+    const elevations = await sampleElevationsBatch(candidates.map(c => ({ lat: c.lat, lng: c.lng })), 13);
+    const minElev = Math.min(...elevations);
+
+    let bestIdx = 0;
+    let maxScore = -Infinity;
+
+    elevations.forEach((elev, idx) => {
+      const c = candidates[idx];
+      const prominence = elev - minElev;
+      const centrality = 1 - Math.abs(c.frac - 0.5);
+      const score = prominence * 1.5 + centrality * 100;
+      if (score > maxScore) {
+        maxScore = score;
+        bestIdx = idx;
+      }
+    });
+
+    const chosen = candidates[bestIdx];
+    const chosenElev = elevations[bestIdx] || 800;
+    const prominence = Math.round(chosenElev - minElev);
+
+    const gapViewshed = await calculateRadialViewshed(chosen.lat, chosen.lng, mastHeight, radiusMeters, 72);
+
+    gapFillerData = {
+      lat: chosen.lat,
+      lng: chosen.lng,
+      groundElevation: Math.round(chosenElev),
+      recommendedHeight: mastHeight,
+      totalElevation: Math.round(chosenElev + mastHeight),
+      prominenceAboveValley: prominence,
+      mitigationPercent: Math.min(92, Math.max(68, Math.round(75 + prominence / 15))),
+      viewshed: gapViewshed,
+      tacticalRationale: `قمة جبلية حاكمة بارتفاع ${Math.round(chosenElev)}م (أعلى بـ ${prominence}م من بطن الوادي)، تشرف على الفجوات العمياء وتغلق ثغرات التسلل المحصورة بين المركزين.`
+    };
+  }
+
+  return {
+    postA: {
+      lat: postALat,
+      lng: postALng,
+      height: mastHeight,
+      elevation: viewshedA.centerElevation,
+      viewshed: viewshedA
+    },
+    postB: {
+      lat: postBLat,
+      lng: postBLng,
+      height: mastHeight,
+      elevation: viewshedB.centerElevation,
+      viewshed: viewshedB
+    },
+    distanceKm,
+    azimuthDeg,
+    combinedVisibleKm2,
+    combinedBlindKm2,
+    combinedCoveragePercent,
+    gapFiller: gapFillerData
+  };
+}
+
+
 export interface MinefieldAnalysisResult {
  startLat: number;
  startLng: number;
@@ -1931,4 +2122,294 @@ export async function calculateRealIPBOverlays(
  };
 }

+export interface PtzGeolocationResult {
+  cameraLat: number;
+  cameraLng: number;
+  cameraElevation: number;
+  mastHeight: number;
+  azimuthDeg: number;
+  tiltDeg: number;
+  fovDeg: number;
+  targetLat: number;
+  targetLng: number;
+  targetElevation: number;
+  slantRangeMeters: number;
+  groundDistanceMeters: number;
+  elevationDifferenceMeters: number;
+  isIntersected: boolean;
+  isHorizonOpen?: boolean;
+  horizonDipDeg?: number;
+  geometricHorizonDistanceMeters?: number;
+  statusMessage?: string;
+  fovConeGeoJson: any;
+  rayLineGeoJson: any;
+  targetPointGeoJson: any;
+}
+
+/**
+ * Calculates real-time Target Geolocation (Lat, Lng, Elevation, Distance)
+ * from Camera Pan (Azimuth), Tilt Angle, and Mast Height by intersecting
+ * the optical line-of-sight ray with the 3D Digital Elevation Model (DEM).
+ * Features sub-step binary search root refinement to avoid 50m discretization jumps.
+ */
+export async function calculatePtzTargetGeolocation(
+  cameraLat: number,
+  cameraLng: number,
+  mastHeight: number = 20,
+  azimuthDeg: number = 50,
+  tiltDeg: number = -5,
+  fovDeg: number = 20
+): Promise<PtzGeolocationResult> {
+  const cameraGroundElev = await sampleElevationAt(cameraLat, cameraLng);
+  const cameraTotalElev = cameraGroundElev + mastHeight;
+
+  const R_earth = 6371000;
+  const k_refraction = 0.13;
+  const effectiveRadius = R_earth / (1 - k_refraction); // 7,322,988 m
+
+  const azRad = (azimuthDeg * Math.PI) / 180;
+  const tiltRad = (tiltDeg * Math.PI) / 180;
+
+  // Horizon calculations for camera mast
+  const geometricHorizonDist = Math.round(Math.sqrt(2 * effectiveRadius * mastHeight));
+  const horizonDipRad = -Math.atan(geometricHorizonDist / effectiveRadius);
+  const horizonDipDeg = Math.round((horizonDipRad * 180 / Math.PI) * 100) / 100; // e.g. -0.13 deg for 20m mast
+
+  const maxRangeMeters = 18000;
+  const stepMeters = 30; // fine 30m steps matching Copernicus/NASA resolution
+  const numSteps = Math.floor(maxRangeMeters / stepMeters);
+
+  const sampleCoords: Array<{ d: number; lat: number; lng: number }> = [];
+
+  for (let s = 1; s <= numSteps; s++) {
+    const d = s * stepMeters;
+    const dLat = (d / R_earth) * (180 / Math.PI) * Math.cos(azRad);
+    const dLng = (d / (R_earth * Math.cos((cameraLat * Math.PI) / 180))) * (180 / Math.PI) * Math.sin(azRad);
+    sampleCoords.push({ d, lat: cameraLat + dLat, lng: cameraLng + dLng });
+  }
+
+  const elevations = await sampleElevationsBatch(sampleCoords.map(c => ({ lat: c.lat, lng: c.lng })), 13);
+
+  let intersectIdx = -1;
+  let targetGroundElev = cameraGroundElev;
+  let finalDist = maxRangeMeters;
+
+  for (let i = 0; i < sampleCoords.length; i++) {
+    const d = sampleCoords[i].d;
+    const earthCurvatureDrop = (d * d) / (2 * effectiveRadius);
+    const rayElev = cameraTotalElev + d * Math.tan(tiltRad) - earthCurvatureDrop;
+    const groundElev = elevations[i];
+
+    if (rayElev <= groundElev) {
+      intersectIdx = i;
+
+      // Sub-step binary search root refinement (5 iterations to pinpoint exact terrain crossing < 1m)
+      let dLow = i > 0 ? sampleCoords[i - 1].d : 0;
+      let dHigh = d;
+      let elevLow = i > 0 ? elevations[i - 1] : cameraGroundElev;
+      let elevHigh = groundElev;
+
+      for (let step = 0; step < 5; step++) {
+        const dMid = (dLow + dHigh) / 2;
+        const curDropMid = (dMid * dMid) / (2 * effectiveRadius);
+        const rayElevMid = cameraTotalElev + dMid * Math.tan(tiltRad) - curDropMid;
+        const f = (dMid - dLow) / (dHigh - dLow || 1);
+        const groundElevMid = elevLow + (elevHigh - elevLow) * f;
+
+        if (rayElevMid <= groundElevMid) {
+          dHigh = dMid;
+          elevHigh = groundElevMid;
+        } else {
+          dLow = dMid;
+          elevLow = groundElevMid;
+        }
+      }
+
+      finalDist = Math.round((dLow + dHigh) / 2);
+      targetGroundElev = Math.round(((elevLow + elevHigh) / 2) * 10) / 10;
+      break;
+    }
+  }
+
+  const isIntersected = intersectIdx !== -1;
+  const isHorizonOpen = !isIntersected && tiltDeg >= horizonDipDeg;
+  let targetLat: number;
+  let targetLng: number;
+  let statusMessage = '';
+
+  if (isIntersected) {
+    const dLat = (finalDist / R_earth) * (180 / Math.PI) * Math.cos(azRad);
+    const dLng = (finalDist / (R_earth * Math.cos((cameraLat * Math.PI) / 180))) * (180 / Math.PI) * Math.sin(azRad);
+    targetLat = Number((cameraLat + dLat).toFixed(6));
+    targetLng = Number((cameraLng + dLng).toFixed(6));
+    statusMessage = `🎯 تقاطع تام مع سطح الأرض عند مسافة ${finalDist}م (ارتفاع ${targetGroundElev}م)`;
+  } else if (isHorizonOpen) {
+    finalDist = Math.min(15000, geometricHorizonDist);
+    const dLat = (finalDist / R_earth) * (180 / Math.PI) * Math.cos(azRad);
+    const dLng = (finalDist / (R_earth * Math.cos((cameraLat * Math.PI) / 180))) * (180 / Math.PI) * Math.sin(azRad);
+    targetLat = Number((cameraLat + dLat).toFixed(6));
+    targetLng = Number((cameraLng + dLng).toFixed(6));
+    targetGroundElev = elevations[elevations.length - 1];
+    statusMessage = `🔭 الرؤية في الأفق المفتوح / الفضاء الجوي (خط النظر يعلو سطح الأرض)`;
+  } else {
+    const lastCoord = sampleCoords[sampleCoords.length - 1];
+    targetLat = Number(lastCoord.lat.toFixed(6));
+    targetLng = Number(lastCoord.lng.toFixed(6));
+    targetGroundElev = elevations[elevations.length - 1];
+    statusMessage = `⚠️ لم يتم رصد تقاطع ضمن المدى الأقصى (${maxRangeMeters}م)`;
+  }
+
+  const slantRangeMeters = Math.round(
+    Math.sqrt(Math.pow(finalDist, 2) + Math.pow(cameraTotalElev - targetGroundElev, 2))
+  );
+
+  // Generate Camera Vision Cone (FOV Frustum)
+  const halfFovRad = ((Math.max(2, fovDeg) / 2) * Math.PI) / 180;
+  const leftAzRad = azRad - halfFovRad;
+  const rightAzRad = azRad + halfFovRad;
+
+  const arcPoints: [number, number][] = [];
+  const arcSteps = 16;
+  for (let a = 0; a <= arcSteps; a++) {
+    const curAz = leftAzRad + (rightAzRad - leftAzRad) * (a / arcSteps);
+    const aLat = cameraLat + (finalDist / R_earth) * (180 / Math.PI) * Math.cos(curAz);
+    const aLng = cameraLng + (finalDist / (R_earth * Math.cos((cameraLat * Math.PI) / 180))) * (180 / Math.PI) * Math.sin(curAz);
+    arcPoints.push([aLng, aLat]);
+  }
+
+  const conePolygonCoords = [
+    [cameraLng, cameraLat],
+    ...arcPoints,
+    [cameraLng, cameraLat]
+  ];
+
+  const fovConeGeoJson = {
+    type: 'Feature',
+    properties: { role: 'ptz-cone', azimuthDeg, tiltDeg, finalDist, isIntersected, isHorizonOpen },
+    geometry: {
+      type: 'Polygon',
+      coordinates: [conePolygonCoords]
+    }
+  };
+
+  const rayLineGeoJson = {
+    type: 'Feature',
+    properties: { role: 'ptz-sightline', isIntersected, isHorizonOpen },
+    geometry: {
+      type: 'LineString',
+      coordinates: [
+        [cameraLng, cameraLat],
+        [targetLng, targetLat]
+      ]
+    }
+  };
+
+  const targetPointGeoJson = {
+    type: 'Feature',
+    properties: {
+      role: 'ptz-target',
+      targetLat: Number(targetLat.toFixed(5)),
+      targetLng: Number(targetLng.toFixed(5)),
+      targetElevation: Math.round(targetGroundElev),
+      slantRangeMeters,
+      groundDistanceMeters: Math.round(finalDist),
+      isIntersected,
+      isHorizonOpen,
+      statusMessage
+    },
+    geometry: {
+      type: 'Point',
+      coordinates: [targetLng, targetLat]
+    }
+  };
+
+  return {
+    cameraLat,
+    cameraLng,
+    cameraElevation: Math.round(cameraTotalElev),
+    mastHeight,
+    azimuthDeg: Number(azimuthDeg.toFixed(1)),
+    tiltDeg: Number(tiltDeg.toFixed(2)),
+    fovDeg,
+    targetLat: Number(targetLat.toFixed(5)),
+    targetLng: Number(targetLng.toFixed(5)),
+    targetElevation: Math.round(targetGroundElev),
+    slantRangeMeters,
+    groundDistanceMeters: Math.round(finalDist),
+    elevationDifferenceMeters: Math.round(targetGroundElev - cameraTotalElev),
+    isIntersected,
+    isHorizonOpen,
+    horizonDipDeg,
+    geometricHorizonDistanceMeters: geometricHorizonDist,
+    statusMessage,
+    fovConeGeoJson,
+    rayLineGeoJson,
+    targetPointGeoJson
+  };
+}
+
+/**
+ * Inverse Kinematics for Electro-Optical PTZ Camera:
+ * Given Camera Pos, Mast Height, and a Target Coordinate clicked on the map,
+ * computes the exact required Pan (Azimuth), Tilt Angle, and Distances accounting
+ * for real DEM elevation, Earth curvature drop, and atmospheric refraction.
+ */
+export async function calculatePtzInverseKinematics(
+  cameraLat: number,
+  cameraLng: number,
+  mastHeight: number,
+  targetLat: number,
+  targetLng: number
+): Promise<{
+  azimuthDeg: number;
+  tiltDeg: number;
+  distanceMeters: number;
+  slantRangeMeters: number;
+  targetElev: number;
+  cameraTotalElev: number;
+}> {
+  const cameraGroundElev = await sampleElevationAt(cameraLat, cameraLng);
+  const cameraTotalElev = cameraGroundElev + mastHeight;
+  const targetElev = await sampleElevationAt(targetLat, targetLng);
+
+  const R_earth = 6371000;
+  const k_refraction = 0.13;
+  const effectiveRadius = R_earth / (1 - k_refraction);
+
+  // Great Circle Distance
+  const dLat = (targetLat - cameraLat) * (Math.PI / 180);
+  const dLng = (targetLng - cameraLng) * (Math.PI / 180);
+  const a = Math.sin(dLat / 2) * Math.sin(dLat / 2) +
+            Math.cos(cameraLat * (Math.PI / 180)) * Math.cos(targetLat * (Math.PI / 180)) *
+            Math.sin(dLng / 2) * Math.sin(dLng / 2);
+  const c = 2 * Math.atan2(Math.sqrt(a), Math.sqrt(1 - a));
+  const distanceMeters = Math.max(1, Math.round(R_earth * c));
+
+  // Azimuth
+  const y = Math.sin(dLng) * Math.cos(targetLat * (Math.PI / 180));
+  const x = Math.cos(cameraLat * (Math.PI / 180)) * Math.sin(targetLat * (Math.PI / 180)) -
+            Math.sin(cameraLat * (Math.PI / 180)) * Math.cos(targetLat * (Math.PI / 180)) * Math.cos(dLng);
+  const azimuthDeg = Math.round(((Math.atan2(y, x) * 180 / Math.PI + 360) % 360) * 10) / 10;
+
+  // Earth curvature drop
+  const earthCurvatureDrop = (distanceMeters * distanceMeters) / (2 * effectiveRadius);
+
+  // Exact Tilt calculation
+  const deltaH = targetElev - cameraTotalElev + earthCurvatureDrop;
+  const tiltRad = Math.atan2(deltaH, distanceMeters);
+  const tiltDeg = Math.round((tiltRad * 180 / Math.PI) * 100) / 100;
+
+  const slantRangeMeters = Math.round(Math.sqrt(distanceMeters * distanceMeters + Math.pow(targetElev - cameraTotalElev, 2)));
+
+  return {
+    azimuthDeg,
+    tiltDeg,
+    distanceMeters,
+    slantRangeMeters,
+    targetElev,
+    cameraTotalElev
+  };
+}
+
+