/**
 * Tactical Elevation & Line of Sight (LOS) Calculation Service
 * خدمة حساب مقطع الارتفاع التضاريسي وتبادل الرؤية العسكري (Intervisibility)
 */

export interface ElevationPoint {
  distance: number;       // Distance from observer (meters)
  lat: number;
  lng: number;
  elevation: number;      // Terrain elevation AMSL (meters)
  rayHeight: number;      // Line of Sight ray elevation at this distance (meters)
  isVisible: boolean;     // Can observer see this terrain point?
  isTargetRayBlocked: boolean; // Does this terrain point block the ray to the final target?
  clearance: number;      // Clearance distance (rayHeight - elevation) in meters
}

export interface ObstacleInfo {
  distance: number;
  elevation: number;
  lat: number;
  lng: number;
  excessHeight: number;  // How much the obstacle penetrates above the ray (meters)
}

export interface LineOfSightResult {
  points: ElevationPoint[];
  totalDistance: number;         // Total distance in meters
  isDirectlyVisible: boolean;    // Is target visible from observer?
  observerElevation: number;     // Ground elevation + observer height
  targetElevation: number;       // Ground elevation + target height
  observerGroundElev: number;    // Raw ground elevation
  targetGroundElev: number;      // Raw ground elevation
  minElevation: number;
  maxElevation: number;
  highestObstacle: ObstacleInfo | null;
  deadGroundPercentage: number;  // % of line hidden behind crests
  angleDegrees: number;          // Vertical angle (degrees)
  angleMils: number;             // Military Artillery Mils (6400 mils = 360 deg)
  azimuthDegrees: number;        // Compass Bearing (0-360 deg)
}

// In-memory cache for DEM tile image data to avoid re-fetching
const tileCache = new Map<string, ImageData>();

/**
 * Calculates Great-Circle Haversine distance in meters
 */
export function calculateDistance(lat1: number, lon1: number, lat2: number, lon2: number): number {
  const R = 6371000; // Earth radius in meters
  const dLat = (lat2 - lat1) * (Math.PI / 180);
  const dLon = (lon2 - lon1) * (Math.PI / 180);
  const a =
    Math.sin(dLat / 2) * Math.sin(dLat / 2) +
    Math.cos(lat1 * (Math.PI / 180)) * Math.cos(lat2 * (Math.PI / 180)) *
    Math.sin(dLon / 2) * Math.sin(dLon / 2);
  const c = 2 * Math.atan2(Math.sqrt(a), Math.sqrt(1 - a));
  return R * c;
}

/**
 * Calculates Forward Azimuth / Bearing (0-360 degrees)
 */
export function calculateAzimuth(lat1: number, lon1: number, lat2: number, lon2: number): number {
  const phi1 = lat1 * (Math.PI / 180);
  const phi2 = lat2 * (Math.PI / 180);
  const deltaLambda = (lon2 - lon1) * (Math.PI / 180);
  const y = Math.sin(deltaLambda) * Math.cos(phi2);
  const x = Math.cos(phi1) * Math.sin(phi2) - Math.sin(phi1) * Math.cos(phi2) * Math.cos(deltaLambda);
  const theta = Math.atan2(y, x);
  return (theta * (180 / Math.PI) + 360) % 360;
}

/**
 * Samples elevation from Terrarium DEM tile or fallback topographic model
 */
/**
 * Fetches or retrieves cached DEM tile ImageData
 */
async function fetchTileImageData(zoom: number, x: number, y: number): Promise<ImageData | null> {
  const tileKey = `${zoom}/${x}/${y}`;
  const cached = tileCache.get(tileKey);
  if (cached) return cached;

  const apiUrl = (import.meta as any).env.VITE_API_URL || '/api';
  const urlsToTry = [
    `${apiUrl}/tactical/dem/${zoom}/${x}/${y}.png`,
    `https://s3.amazonaws.com/elevation-tiles-prod/terrarium/${zoom}/${x}/${y}.png`
  ];

  for (const tileUrl of urlsToTry) {
    try {
      const img = new Image();
      img.crossOrigin = 'anonymous';

      const loadPromise = new Promise<HTMLImageElement>((resolve, reject) => {
        img.onload = () => resolve(img);
        img.onerror = (e) => reject(e);
        img.src = tileUrl;
      });

      const loadedImg = await Promise.race([
        loadPromise,
        new Promise<never>((_, reject) => setTimeout(() => reject(new Error('DEM Timeout')), 4000))
      ]);

      const canvas = document.createElement('canvas');
      canvas.width = 256;
      canvas.height = 256;
      const ctx = canvas.getContext('2d');
      if (ctx) {
        ctx.drawImage(loadedImg, 0, 0);
        const imgData = ctx.getImageData(0, 0, 256, 256);
        tileCache.set(tileKey, imgData);
        return imgData;
      }
    } catch {
      // Continue to next URL
    }
  }

  return null;
}

/**
 * Samples sub-pixel elevation from ImageData using Bilinear Interpolation
 */
function interpolateElevation(imgData: ImageData, subX: number, subY: number): number {
  const clampedX = Math.max(0, Math.min(254.99, subX));
  const clampedY = Math.max(0, Math.min(254.99, subY));

  const x0 = Math.floor(clampedX);
  const x1 = x0 + 1;
  const y0 = Math.floor(clampedY);
  const y1 = y0 + 1;

  const fx = clampedX - x0;
  const fy = clampedY - y0;

  const decodePixel = (px: number, py: number): number => {
    const idx = (py * 256 + px) * 4;
    const r = imgData.data[idx];
    const g = imgData.data[idx + 1];
    const b = imgData.data[idx + 2];
    return (r * 256 + g + b / 256) - 32768;
  };

  const z00 = decodePixel(x0, y0);
  const z10 = decodePixel(x1, y0);
  const z01 = decodePixel(x0, y1);
  const z11 = decodePixel(x1, y1);

  const zTop = z00 * (1 - fx) + z10 * fx;
  const zBottom = z01 * (1 - fx) + z11 * fx;
  return Math.round((zTop * (1 - fy) + zBottom * fy) * 10) / 10;
}

/**
 * Samples elevation for multiple coordinates in batch with tile grouping & cache
 */
export async function sampleElevationsBatch(
  coords: Array<{ lat: number; lng: number }>,
  zoom: number = 13
): Promise<number[]> {
  // 1. Try sovereign backend API for zero-CORS 100% accurate satellite DEM
  try {
    const apiUrl = (import.meta as any).env.VITE_API_URL || '/api';
    const controller = new AbortController();
    const timeoutId = setTimeout(() => controller.abort(), 2500);

    const apiKey = localStorage.getItem('map_admin_key') || localStorage.getItem('intaleq_api_key') || (import.meta as any).env.VITE_ADMIN_API_KEY || (import.meta as any).env.VITE_API_KEY || 'zP9vL5mK2nQ8xR7jT4wS1yB6hG3fV0cX';
    const res = await fetch(`${apiUrl}/tactical/elevations?key=${encodeURIComponent(apiKey)}`, {
      method: 'POST',
      headers: { 'Content-Type': 'application/json', 'x-api-key': apiKey },
      body: JSON.stringify({ coordinates: coords }),
      signal: controller.signal
    });
    clearTimeout(timeoutId);

    if (res.ok) {
      const data = await res.json();
      if (Array.isArray(data) && data.length === coords.length) {
        return data;
      }
    }
  } catch {
    // Fall back to client tile decoding
  }

  const n = Math.pow(2, zoom);

  // Group coordinates by tile
  const tileMap = new Map<string, { zoom: number; x: number; y: number; points: Array<{ index: number; subX: number; subY: number; lat: number; lng: number }> }>();

  coords.forEach((coord, index) => {
    const x = Math.floor(((coord.lng + 180) / 360) * n);
    const latRad = (coord.lat * Math.PI) / 180;
    const y = Math.floor((1 - Math.log(Math.tan(latRad) + 1 / Math.cos(latRad)) / Math.PI) / 2 * n);

    const subX = (((coord.lng + 180) / 360) * n - x) * 256;
    const subY = ((1 - Math.log(Math.tan(latRad) + 1 / Math.cos(latRad)) / Math.PI) / 2 * n - y) * 256;

    const tileKey = `${zoom}/${x}/${y}`;
    if (!tileMap.has(tileKey)) {
      tileMap.set(tileKey, { zoom, x, y, points: [] });
    }
    tileMap.get(tileKey)!.points.push({ index, subX, subY, lat: coord.lat, lng: coord.lng });
  });

  const results: number[] = new Array(coords.length).fill(0);

  // Load all unique tiles in parallel
  await Promise.all(
    Array.from(tileMap.entries()).map(async ([_key, entry]) => {
      const imgData = await fetchTileImageData(entry.zoom, entry.x, entry.y);
      entry.points.forEach((pt) => {
        if (imgData) {
          results[pt.index] = interpolateElevation(imgData, pt.subX, pt.subY);
        } else {
          results[pt.index] = getApproximateElevation(pt.lat, pt.lng);
        }
      });
    })
  );

  return results;
}

/**
 * Samples elevation from Terrarium DEM tile or fallback topographic model
 */
async function sampleElevationAt(lat: number, lng: number): Promise<number> {
  const zoom = 13;
  const n = Math.pow(2, zoom);
  const x = Math.floor(((lng + 180) / 360) * n);
  const latRad = (lat * Math.PI) / 180;
  const y = Math.floor((1 - Math.log(Math.tan(latRad) + 1 / Math.cos(latRad)) / Math.PI) / 2 * n);

  const subX = (((lng + 180) / 360) * n - x) * 256;
  const subY = ((1 - Math.log(Math.tan(latRad) + 1 / Math.cos(latRad)) / Math.PI) / 2 * n - y) * 256;

  const imgData = await fetchTileImageData(zoom, x, y);
  if (imgData) {
    return Math.round(interpolateElevation(imgData, subX, subY));
  }
  return getApproximateElevation(lat, lng);
}

/**
 * Topographic estimation model for Jordan terrain when tiles are offline
 */
function getApproximateElevation(lat: number, lng: number): number {
  if (lng < 35.65 && lat < 32.5 && lat > 30.8) {
    const riftCenter = 35.50;
    const distFromRift = Math.abs(lng - riftCenter);
    const riftElev = -420.0 + distFromRift * 3200.0;
    if (distFromRift < 0.15) {
      return Math.round(Math.max(-430.0, Math.min(1000.0, riftElev)));
    }
  }

  const DEM_BENCHMARKS: [number, number, number][] = [
    [32.0720, 36.0880, 610.0],  // الزرقاء
    [32.1150, 35.9550, 750.0],  // بيرين
    [32.1250, 36.1150, 580.0],  // الهاشمية
    [32.0250, 36.0350, 670.0],  // الرصيفة
    [31.9615, 35.9130, 740.0],  // وسط البلد
    [32.0220, 35.8450, 1060.0], // صويلح
    [32.2780, 35.8950, 600.0],  // جرش
    [32.3250, 35.7350, 1150.0], // عجلون
    [32.5450, 35.8550, 620.0],  // إربد
    [32.3560, 36.2590, 700.0],  // المفرق
  ];

  let num = 0.0;
  let den = 0.0;
  const latRad = lat * (Math.PI / 180.0);

  for (const c of DEM_BENCHMARKS) {
    const dLat = (lat - c[0]) * 111.0;
    const dLng = (lng - c[1]) * 111.0 * Math.cos(latRad);
    const distKm = Math.max(Math.sqrt(dLat * dLat + dLng * dLng), 0.15);
    const w = 1.0 / Math.pow(distKm, 2.0);
    num += w * c[2];
    den += w;
  }
  return Math.round(den > 0 ? num / den : 620.0);
}

/**
 * Calculates Line of Sight and Elevation Profile between two coordinates
 */
export async function calculateLineOfSight(
  startLat: number,
  startLng: number,
  endLat: number,
  endLng: number,
  obsHeight: number = 2,
  tgtHeight: number = 2,
  samples: number = 60
): Promise<LineOfSightResult> {
  const totalDistance = calculateDistance(startLat, startLng, endLat, endLng);
  const azimuthDegrees = calculateAzimuth(startLat, startLng, endLat, endLng);

  // Try fetching high-precision result from Backend Tactical API
  try {
    const apiUrl = (import.meta as any).env.VITE_API_URL || '/api';
    const apiKey = (import.meta as any).env.VITE_ADMIN_API_KEY || (import.meta as any).env.VITE_API_KEY || localStorage.getItem('map_admin_key') || localStorage.getItem('intaleq_api_key') || 'zP9vL5mK2nQ8xR7jT4wS1yB6hG3fV0cX';

    const controller = new AbortController();
    const timeoutId = setTimeout(() => controller.abort(), 8000);

    const res = await fetch(`${apiUrl}/tactical/line-of-sight?observerLat=${startLat}&observerLng=${startLng}&targetLat=${endLat}&targetLng=${endLng}&observerHeight=${obsHeight}&targetHeight=${tgtHeight}&samples=${samples}&key=${encodeURIComponent(apiKey)}`, {
      headers: { 'x-api-key': apiKey },
      signal: controller.signal
    });
    clearTimeout(timeoutId);

    if (res.ok) {
      const data = await res.json();
      return {
        points: (data.profile || []).map((p: any) => ({
          distance: p.distanceMeters,
          lat: p.lat,
          lng: p.lng,
          elevation: p.groundElevationMeters,
          rayHeight: p.sightRayElevationMeters,
          isVisible: p.isVisible,
          isTargetRayBlocked: !p.isVisible && p.marginMeters < 0,
          clearance: p.marginMeters
        })),
        totalDistance: data.summary.totalDistanceMeters,
        isDirectlyVisible: data.isDirectlyVisible,
        observerElevation: data.summary.observerTotalElevationMeters,
        targetElevation: data.summary.targetTotalElevationMeters,
        observerGroundElev: data.summary.observerGroundElevationMeters,
        targetGroundElev: data.summary.targetGroundElevationMeters,
        minElevation: data.summary.minElevationMeters,
        maxElevation: data.summary.maxElevationMeters,
        highestObstacle: data.highestObstacle ? {
          distance: data.highestObstacle.distanceMeters ?? data.highestObstacle.distance,
          elevation: data.highestObstacle.elevationMeters ?? data.highestObstacle.groundElevationMeters ?? data.highestObstacle.elevation,
          lat: data.highestObstacle.lat,
          lng: data.highestObstacle.lng,
          excessHeight: data.highestObstacle.penetrationMeters ?? data.highestObstacle.excessHeightMeters ?? data.highestObstacle.excessHeight
        } : null,
        deadGroundPercentage: data.summary.deadGroundPercentage,
        angleDegrees: data.summary.verticalAngleDegrees,
        angleMils: data.summary.verticalAngleMilsNato,
        azimuthDegrees: data.summary.azimuthDegrees
      };
    }
  } catch {
    // Fall back to local DEM processing
  }

  // Generate sample coordinates along the geodesic path
  const sampleCoords: { lat: number; lng: number; dist: number }[] = [];
  for (let i = 0; i <= samples; i++) {
    const fraction = i / samples;
    const lat = startLat + (endLat - startLat) * fraction;
    const lng = startLng + (endLng - startLng) * fraction;
    const dist = totalDistance * fraction;
    sampleCoords.push({ lat, lng, dist });
  }

  // Fetch elevations for all sample points in fast batch
  const elevations = await sampleElevationsBatch(sampleCoords, 13);

  const observerGroundElev = elevations[0];
  const targetGroundElev = elevations[elevations.length - 1];
  const observerElevation = observerGroundElev + obsHeight;
  const targetElevation = targetGroundElev + tgtHeight;

  const R_earth = 6371000;
  const k_refraction = 0.13;
  const effectiveEarthRadius = R_earth / (1 - k_refraction);

  let isDirectlyVisible = true;
  let highestObstacle: ObstacleInfo | null = null;
  let maxObstacleExcess = 0;
  let deadGroundCount = 0;
  let maxAngleSoFar = -Infinity;

  const points: ElevationPoint[] = [];
  let minElev = Infinity;
  let maxElev = -Infinity;

  for (let i = 0; i <= samples; i++) {
    const d = sampleCoords[i].dist;
    const elev = elevations[i];
    minElev = Math.min(minElev, elev);
    maxElev = Math.max(maxElev, elev);

    const earthCurvatureDrop = (d * (totalDistance - d)) / (2 * effectiveEarthRadius);
    const rayHeight = observerElevation + ((targetElevation - observerElevation) * (d / totalDistance)) - earthCurvatureDrop;
    const clearance = rayHeight - elev;

    let isTargetRayBlocked = false;
    if (i > 1 && i < samples) {
      if (elev > rayHeight) {
        isDirectlyVisible = false;
        isTargetRayBlocked = true;
        const excess = elev - rayHeight;
        if (excess > maxObstacleExcess) {
          maxObstacleExcess = excess;
          highestObstacle = {
            distance: Math.round(d),
            elevation: Math.round(elev),
            lat: sampleCoords[i].lat,
            lng: sampleCoords[i].lng,
            excessHeight: Math.round(excess * 10) / 10,
          };
        }
      }
    }

    let isVisible = true;
    if (i === 0) {
      isVisible = true;
    } else {
      const angleFromObs = (elev - observerElevation) / d;
      if (angleFromObs >= maxAngleSoFar) {
        maxAngleSoFar = angleFromObs;
        isVisible = true;
      } else {
        isVisible = false;
        deadGroundCount++;
      }
    }

    points.push({
      distance: Math.round(d),
      lat: sampleCoords[i].lat,
      lng: sampleCoords[i].lng,
      elevation: Math.round(elev),
      rayHeight: Math.round(rayHeight * 10) / 10,
      isVisible,
      isTargetRayBlocked,
      clearance: Math.round(clearance * 10) / 10,
    });
  }

  const verticalDiff = targetElevation - observerElevation;
  const angleRad = Math.atan2(verticalDiff, totalDistance);
  const angleDegrees = Math.round((angleRad * (180 / Math.PI)) * 100) / 100;
  const angleMils = Math.round((angleDegrees * (6400 / 360)) * 10) / 10;
  const deadGroundPercentage = Math.round((deadGroundCount / samples) * 100);

  return {
    points,
    totalDistance: Math.round(totalDistance),
    isDirectlyVisible,
    observerElevation: Math.round(observerElevation),
    targetElevation: Math.round(targetElevation),
    observerGroundElev: Math.round(observerGroundElev),
    targetGroundElev: Math.round(targetGroundElev),
    minElevation: Math.round(minElev),
    maxElevation: Math.round(maxElev),
    highestObstacle,
    deadGroundPercentage,
    angleDegrees,
    angleMils,
    azimuthDegrees: Math.round(azimuthDegrees * 10) / 10,
  };
}

export interface Viewshed360Result {
  centerLat: number;
  centerLng: number;
  centerElevation: number;
  radiusMeters: number;
  polygonGeoJson: any;
  totalRays: number;
  visibleAreaKm2: number;
  visiblePercentage: number;
}

/**
 * Calculates 360-degree Radial Viewshed around observer
 */
export async function calculateRadialViewshed(
  centerLat: number,
  centerLng: number,
  obsHeight: number = 2,
  radiusMeters: number = 5000,
  numRays: number = 72
): Promise<Viewshed360Result> {
  // 1. Try sovereign backend API for zero-CORS 100% accurate satellite DEM Viewshed
  try {
    const apiUrl = (import.meta as any).env.VITE_API_URL || '/api';
    const controller = new AbortController();
    const timeoutId = setTimeout(() => controller.abort(), 3500);

    const apiKey = (import.meta as any).env.VITE_ADMIN_API_KEY || (import.meta as any).env.VITE_API_KEY || localStorage.getItem('map_admin_key') || localStorage.getItem('intaleq_api_key') || 'zP9vL5mK2nQ8xR7jT4wS1yB6hG3fV0cX';
    const res = await fetch(`${apiUrl}/tactical/viewshed?lat=${centerLat}&lng=${centerLng}&height=${obsHeight}&radius=${radiusMeters}&rays=${numRays}&key=${encodeURIComponent(apiKey)}`, {
      headers: { 'x-api-key': apiKey },
      signal: controller.signal
    });
    clearTimeout(timeoutId);

    if (res.ok) {
      const data = await res.json();
      return {
        centerLat,
        centerLng,
        radiusMeters,
        polygonGeoJson: data.polygon,
        totalRays: numRays,
        visibleAreaKm2: data.coveredAreaKm2,
        visiblePercentage: data.coveragePercent
      };
    }
  } catch {
    // Fall back to client calculation
  }

  const centerGroundElev = await sampleElevationAt(centerLat, centerLng);
  const observerElevation = centerGroundElev + obsHeight;

  const R_earth = 6371000;
  const k_refraction = 0.13;
  const effectiveEarthRadius = R_earth / (1 - k_refraction);

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

  for (let rayIdx = 0; rayIdx < numRays; rayIdx++) {
    const angleDeg = (rayIdx * 360) / numRays;
    const angleRad = (angleDeg * Math.PI) / 180;

    const dLat = (radiusMeters / R_earth) * (180 / Math.PI) * Math.cos(angleRad);
    const dLng = (radiusMeters / (R_earth * Math.cos((centerLat * Math.PI) / 180))) * (180 / Math.PI) * Math.sin(angleRad);

    const endLat = centerLat + dLat;
    const endLng = centerLng + dLng;

    for (let step = 1; step <= samplesPerRay; step++) {
      const frac = step / samplesPerRay;
      const sLat = centerLat + (endLat - centerLat) * frac;
      const sLng = centerLng + (endLng - centerLng) * frac;
      const sDist = radiusMeters * frac;
      rayCoords.push({ rayIdx, step, sDist, sLat, sLng });
    }
  }

  // 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;

  for (let rayIdx = 0; rayIdx < numRays; rayIdx++) {
    let maxAngleSoFar = -Infinity;
    let visibleHorizonDist = radiusMeters;
    let visibleHorizonLat = centerLat;
    let visibleHorizonLng = centerLng;

    const rayPoints = rayCoords.filter(rc => rc.rayIdx === rayIdx);
    rayPoints.forEach(rc => {
      const sElev = elevations[rc.rayIdx * samplesPerRay + (rc.step - 1)];
      const earthCurvatureDrop = (rc.sDist * rc.sDist) / (2 * effectiveEarthRadius);
      const apparentElev = sElev - earthCurvatureDrop;
      const angle = (apparentElev - observerElevation) / rc.sDist;

      if (angle >= maxAngleSoFar) {
        maxAngleSoFar = angle;
        visibleHorizonDist = rc.sDist;
        visibleHorizonLat = rc.sLat;
        visibleHorizonLng = rc.sLng;
      }
    });

    totalVisibleDistanceSum += visibleHorizonDist;
    polygonCoordinates.push([visibleHorizonLng, visibleHorizonLat]);
  }

  if (polygonCoordinates.length > 0) {
    polygonCoordinates.push(polygonCoordinates[0]);
  }

  const polygonGeoJson = {
    type: 'Feature',
    properties: {
      centerLat,
      centerLng,
      radiusMeters
    },
    geometry: {
      type: 'Polygon',
      coordinates: [polygonCoordinates]
    }
  };

  const avgVisibleDist = totalVisibleDistanceSum / numRays;
  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));

  return {
    centerLat,
    centerLng,
    centerElevation: Math.round(observerElevation),
    radiusMeters,
    polygonGeoJson,
    totalRays: numRays,
    visibleAreaKm2: Math.round(actualVisibleArea * 10) / 10,
    visiblePercentage
  };
}

export interface MinefieldAnalysisResult {
  startLat: number;
  startLng: number;
  endLat: number;
  endLng: number;
  frontageMeters: number;
  azimuthDegrees: number;
  startElevation: number;
  endElevation: number;
  avgSlopeDegrees: number;
  suitabilityScore: number; // 0-100%
  suitabilityVerdict: string;
  suitabilityVerdictAr: string;
  isChokePoint: boolean;
  quantities: {
    antiTankMines: number;
    antiPersonnelMines: number;
    fuseSets: number;
    warningSigns: number;
    rowsCount: number;
    depthMeters: number;
    estimatedDeploymentHours: number;
  };
}

/**
 * Tactical Minefield Barrier & Sapper Engineering Analysis
 */
export async function calculateMinefieldAnalysis(
  startLat: number,
  startLng: number,
  endLat: number,
  endLng: number,
  densityLevel: 'standard' | 'dense' | 'light' = 'standard'
): Promise<MinefieldAnalysisResult> {
  const frontageMeters = calculateDistance(startLat, startLng, endLat, endLng);
  const azimuthDegrees = calculateAzimuth(startLat, startLng, endLat, endLng);

  const samples = 12;
  const sampleCoords: Array<{ lat: number; lng: number }> = [];
  for (let i = 0; i <= samples; i++) {
    const frac = i / samples;
    const lat = startLat + (endLat - startLat) * frac;
    const lng = startLng + (endLng - startLng) * frac;
    sampleCoords.push({ lat, lng });
  }

  const elevs = await sampleElevationsBatch(sampleCoords, 13);
  const startElev = elevs[0];
  const endElev = elevs[elevs.length - 1];

  const elevDiff = Math.abs(endElev - startElev);
  const avgSlopeRad = frontageMeters > 0 ? Math.atan(elevDiff / frontageMeters) : 0;
  const avgSlopeDegrees = Math.round((avgSlopeRad * (180 / Math.PI)) * 10) / 10;

  const midElev = elevs[Math.floor(samples / 2)];
  const isChokePoint = midElev < Math.min(startElev, endElev) - 10 || avgSlopeDegrees <= 10;

  const densityFactor = densityLevel === 'dense' ? 1.5 : densityLevel === 'light' ? 0.6 : 1.0;
  const atMinesPerMeter = 1.0 * densityFactor;
  const apMinesPerMeter = 2.0 * densityFactor;

  const antiTankMines = Math.round(frontageMeters * atMinesPerMeter);
  const antiPersonnelMines = Math.round(frontageMeters * apMinesPerMeter);
  const warningSigns = Math.max(4, Math.round(frontageMeters / 30));
  const deploymentHours = Math.round(((frontageMeters / 150) * 2 * densityFactor) * 10) / 10;

  let suitabilityScore = 85;
  let suitabilityVerdict = 'Highly Effective Defensive Barrier';
  let suitabilityVerdictAr = 'موضع دفاعي نموذجي يسد محور اقتراب الدروع بفعالية';

  if (avgSlopeDegrees > 20) {
    suitabilityScore = 55;
    suitabilityVerdict = 'Steep Terrain: Armor movement already restricted by natural slope';
    suitabilityVerdictAr = 'انحدار جبلي شديد: التضاريس تشكل مانعاً طبيعياً لحركة الدروع';
  } else if (isChokePoint) {
    suitabilityScore = 96;
    suitabilityVerdict = 'Optimal Choke Point Defile Barrier';
    suitabilityVerdictAr = 'ممر ومضيق إجباري مثالي (Choke Point) يمنع التفاف آليات العدو';
  }

  return {
    startLat,
    startLng,
    endLat,
    endLng,
    frontageMeters: Math.round(frontageMeters),
    azimuthDegrees: Math.round(azimuthDegrees * 10) / 10,
    startElevation: Math.round(startElev),
    endElevation: Math.round(endElev),
    avgSlopeDegrees,
    suitabilityScore,
    suitabilityVerdict,
    suitabilityVerdictAr,
    isChokePoint,
    quantities: {
      antiTankMines,
      antiPersonnelMines,
      fuseSets: Math.round((antiTankMines + antiPersonnelMines) * 1.1),
      warningSigns,
      rowsCount: densityLevel === 'dense' ? 4 : 3,
      depthMeters: densityLevel === 'dense' ? 75 : 50,
      estimatedDeploymentHours: Math.max(1, deploymentHours)
    }
  };
}

export interface CliffFeature {
  lat: number;
  lng: number;
  dropMeters: number;
  slopeDegrees: number;
  label: string;
  tacticalImpact: string;
}

export interface TerrainStudyResult {
  centerLat: number;
  centerLng: number;
  radiusMeters: number;
  areaKm2: number;
  sectorName: string;
  centerElevation: number;
  minElevation: number;
  maxElevation: number;
  reliefMeters: number;
  avgSlopeDegrees: number;
  maxSlopeDegrees: number;
  slopeDistribution: {
    flatPercentage: number;     // < 7 deg (Go)
    moderatePercentage: number; // 7-15 deg (Slow-Go / Rolling)
    steepPercentage: number;    // 15-25 deg (Severe Slow-Go / Rugged)
    cliffPercentage: number;    // > 25 deg (No-Go / Cliffs & Rock Barriers)
  };
  cliffsCount: number;
  cliffs: CliffFeature[];
  terrainClassification: string;
  terrainClassificationAr: string;
  mobilityStatus: 'GO' | 'SLOW-GO' | 'NO-GO';
  mobilityStatusAr: string;
  highestPoint: {
    lat: number;
    lng: number;
    elevation: number;
    label: string;
  };
  lowestPoint: {
    lat: number;
    lng: number;
    elevation: number;
    label: string;
  };
  spatialGeoJson: {
    type: 'FeatureCollection';
    features: any[];
  };
  naturalObstacles: Array<{
    name: string;
    type: string;
    description: string;
    impact: string;
  }>;
  manMadeObstacles: Array<{
    name: string;
    type: string;
    description: string;
    tacticalNote: string;
  }>;
  urbanAndDemographics: {
    density: 'HIGH' | 'MEDIUM' | 'LOW' | 'SPARSE';
    densityAr: string;
    settlements: string[];
    moutComplexity: string;
    collateralDamageRisk: string;
  };
  landCoverAndCover: {
    soilType: string;
    trafficability: string;
    concealmentRating: string;
    airObservationExposure: string;
  };
  oakocAssessment: {
    obstacles: string;
    avenuesOfApproach: string;
    keyTerrain: string;
    observationAndFieldsOfFire: string;
    coverAndConcealment: string;
  };
  tacticalRecommendations: string[];
}

/**
 * Military Tactical Terrain Intelligence & Environmental Study (OAKOC Doctrine)
 * دراسة الأرض الشاملة عالية الدقة: حساب الميول الدقيقة، كشف القواطع والجروف الصخرية، وممرات الحركة
 */
export async function calculateTerrainStudy(
  centerLat: number,
  centerLng: number,
  radiusMeters: number = 3000
): Promise<TerrainStudyResult> {
  const centerElev = await sampleElevationAt(centerLat, centerLng);
  const areaKm2 = Math.round(Math.PI * Math.pow(radiusMeters / 1000, 2) * 10) / 10;

  // High-Density Sampling Matrix: 25x25 = 625 nodes (576 terrain cells)
  const gridSize = 25;
  const R_earth = 6371000;
  const dLatTotal = (radiusMeters / R_earth) * (180 / Math.PI);
  const dLngTotal = (radiusMeters / (R_earth * Math.cos((centerLat * Math.PI) / 180))) * (180 / Math.PI);

  const minLat = centerLat - dLatTotal;
  const maxLat = centerLat + dLatTotal;
  const minLng = centerLng - dLngTotal;
  const maxLng = centerLng + dLngTotal;

  // 1. Generate all grid coordinates
  const gridCoords: Array<{ r: number; c: number; lat: number; lng: number }> = [];
  for (let r = 0; r < gridSize; r++) {
    const lat = Number((minLat + (r / (gridSize - 1)) * (maxLat - minLat)).toFixed(6));
    for (let c = 0; c < gridSize; c++) {
      const lng = Number((minLng + (c / (gridSize - 1)) * (maxLng - minLng)).toFixed(6));
      gridCoords.push({ r, c, lat, lng });
    }
  }

  // 2. High-speed parallel batch sampling from DEM tiles
  const sampleResults = await sampleElevationsBatch(gridCoords, 13);
  const elevationGrid: number[][] = Array(gridSize).fill(0).map(() => Array(gridSize).fill(0));

  gridCoords.forEach((pt, idx) => {
    elevationGrid[pt.r][pt.c] = sampleResults[idx];
  });

  // 3. Precise Peak and Valley Identification
  let maxElev = -9999;
  let minElev = 9999;
  let highestPoint = { lat: centerLat, lng: centerLng, elevation: centerElev, label: 'أعلى قمة' };
  let lowestPoint = { lat: centerLat, lng: centerLng, elevation: centerElev, label: 'أخفض نقطة' };

  gridCoords.forEach((p, idx) => {
    const el = sampleResults[idx];
    const dist = calculateDistance(centerLat, centerLng, p.lat, p.lng);
    if (dist <= radiusMeters * 1.05) {
      if (el > maxElev) {
        maxElev = el;
        highestPoint = { lat: p.lat, lng: p.lng, elevation: Math.round(el), label: `🔺 أعلى قمة: ${Math.round(el)}م` };
      }
      if (el < minElev) {
        minElev = el;
        lowestPoint = { lat: p.lat, lng: p.lng, elevation: Math.round(el), label: `🔻 أخفض نقطة: ${Math.round(el)}م` };
      }
    }
  });

  const relief = Math.round(maxElev - minElev);

  // 4. Compute Local Gradients, Slopes, and Detect Natural Rock Barriers / Cliffs
  const spatialFeatures: any[] = [];
  const slopeValues: number[] = [];
  const detectedCliffs: CliffFeature[] = [];

  const dxMeters = (2 * radiusMeters) / (gridSize - 1);
  const dyMeters = (2 * radiusMeters) / (gridSize - 1);

  let flatCount = 0;      // < 7 deg (GO)
  let rollingCount = 0;   // 7-15 deg (SLOW-GO)
  let ruggedCount = 0;    // 15-25 deg (SEVERE SLOW-GO)
  let cliffCount = 0;     // > 25 deg (NO-GO / Cliffs & Rock Barriers)

  for (let r = 0; r < gridSize - 1; r++) {
    for (let c = 0; c < gridSize - 1; c++) {
      const p1Coord = gridCoords[r * gridSize + c];
      const p2Coord = gridCoords[r * gridSize + (c + 1)];
      const p3Coord = gridCoords[(r + 1) * gridSize + (c + 1)];
      const p4Coord = gridCoords[(r + 1) * gridSize + c];

      const cellCenterLat = (p1Coord.lat + p3Coord.lat) / 2;
      const cellCenterLng = (p1Coord.lng + p3Coord.lng) / 2;
      const distFromCenter = calculateDistance(centerLat, centerLng, cellCenterLat, cellCenterLng);

      if (distFromCenter > radiusMeters * 1.03) continue;

      const z00 = elevationGrid[r][c];
      const z01 = elevationGrid[r][c + 1];
      const z10 = elevationGrid[r + 1][c];
      const z11 = elevationGrid[r + 1][c + 1];

      const cellElevAvg = Math.round((z00 + z01 + z10 + z11) / 4);

      // Central difference gradient
      const dz_dx = ((z01 + z11) - (z00 + z10)) / (2 * dxMeters);
      const dz_dy = ((z10 + z11) - (z00 + z01)) / (2 * dyMeters);

      const slopeRad = Math.atan(Math.sqrt(dz_dx * dz_dx + dz_dy * dz_dy));
      const slopeDeg = Math.round((slopeRad * (180 / Math.PI)) * 10) / 10;
      slopeValues.push(slopeDeg);

      // Local maximum vertical step drop across cell corners
      const localMaxZ = Math.max(z00, z01, z10, z11);
      const localMinZ = Math.min(z00, z01, z10, z11);
      const localDrop = Math.round(localMaxZ - localMinZ);

      const p1 = [p1Coord.lng, p1Coord.lat];
      const p2 = [p2Coord.lng, p2Coord.lat];
      const p3 = [p3Coord.lng, p3Coord.lat];
      const p4 = [p4Coord.lng, p4Coord.lat];

      // Military Grade Slope & Barrier Classification:
      // - Flat (GO): < 7° (< 12%)
      // - Rolling (SLOW-GO): 7° - 15° (12% - 27%)
      // - Rugged / Mountainous (SEVERE SLOW-GO): 15° - 28° (27% - 53%)
      // - Natural Cliffs & Rock Escarpments (NO-GO): >= 28° (Physical vertical barrier)
      const isCliff = slopeDeg >= 28;

      if (isCliff) {
        cliffCount++;
        detectedCliffs.push({
          lat: Number(cellCenterLat.toFixed(5)),
          lng: Number(cellCenterLng.toFixed(5)),
          dropMeters: localDrop,
          slopeDegrees: slopeDeg,
          label: `🧗‍♂️ جرف صخري حاد (${localDrop}م / ${slopeDeg}°)`,
          tacticalImpact: 'مانع طبيعي قطعي - غير قابل لاجتياز الدروع والآليات، يفرض ممرات إجبارية'
        });

        // Add cliff barrier polygon without title label spam
        spatialFeatures.push({
          type: 'Feature',
          properties: {
            layerType: 'cliff',
            category: 'natural-obstacle',
            color: '#ef4444',
            slope: slopeDeg,
            drop: localDrop,
            elevation: cellElevAvg
          },
          geometry: {
            type: 'Polygon',
            coordinates: [[p1, p2, p3, p4, p1]]
          }
        });
      } else if (slopeDeg >= 15) {
        ruggedCount++;
        spatialFeatures.push({
          type: 'Feature',
          properties: {
            layerType: 'slope-sector',
            category: 'severe-slow-go',
            color: '#f97316',
            slope: slopeDeg,
            elevation: cellElevAvg
          },
          geometry: {
            type: 'Polygon',
            coordinates: [[p1, p2, p3, p4, p1]]
          }
        });
      } else if (slopeDeg >= 7) {
        rollingCount++;
        spatialFeatures.push({
          type: 'Feature',
          properties: {
            layerType: 'slope-sector',
            category: 'slow-go',
            color: '#eab308',
            slope: slopeDeg,
            elevation: cellElevAvg
          },
          geometry: {
            type: 'Polygon',
            coordinates: [[p1, p2, p3, p4, p1]]
          }
        });
      } else {
        flatCount++;
        spatialFeatures.push({
          type: 'Feature',
          properties: {
            layerType: 'slope-sector',
            category: 'go',
            color: '#22c55e',
            slope: slopeDeg,
            elevation: cellElevAvg
          },
          geometry: {
            type: 'Polygon',
            coordinates: [[p1, p2, p3, p4, p1]]
          }
        });
      }
    }
  }

  // Add Key Terrain POI Points (Only these will display text badges on the map)
  spatialFeatures.push({
    type: 'Feature',
    properties: {
      category: 'peak',
      layerType: 'poi',
      title: `▲ أعلى قمة: ${highestPoint.elevation}م`,
      elevation: highestPoint.elevation,
      color: '#22c55e'
    },
    geometry: {
      type: 'Point',
      coordinates: [highestPoint.lng, highestPoint.lat]
    }
  });

  spatialFeatures.push({
    type: 'Feature',
    properties: {
      category: 'valley',
      layerType: 'poi',
      title: `▼ أخفض نقطة: ${lowestPoint.elevation}م`,
      elevation: lowestPoint.elevation,
      color: '#06b6d4'
    },
    geometry: {
      type: 'Point',
      coordinates: [lowestPoint.lng, lowestPoint.lat]
    }
  });

  // 5. Statistics of Slope Distribution
  const totalCells = Math.max(1, slopeValues.length);
  const flatPct = Math.round((flatCount / totalCells) * 100);
  const modPct = Math.round((rollingCount / totalCells) * 100);
  const ruggedPct = Math.round((ruggedCount / totalCells) * 100);
  const cliffPct = Math.max(0, 100 - flatPct - modPct - ruggedPct);

  const avgSlope = slopeValues.length > 0 ? Math.round((slopeValues.reduce((a, b) => a + b, 0) / slopeValues.length) * 10) / 10 : 2.5;
  const maxSlope = slopeValues.length > 0 ? Math.round(Math.max(...slopeValues) * 10) / 10 : 5.0;

  // 6. Scientific Hydrological D8 Flow Routing & Natural Wadi (Thalweg) Extraction
  // Compute steepest descent neighbor for every grid cell
  interface FlowTarget {
    r: number;
    c: number;
    slope: number;
  }

  const flowTarget: (FlowTarget | null)[][] = Array(gridSize).fill(null).map(() => Array(gridSize).fill(null));

  for (let r = 0; r < gridSize; r++) {
    for (let c = 0; c < gridSize; c++) {
      const zCurr = elevationGrid[r][c];
      let maxDropSlope = 0;
      let bestTarget: FlowTarget | null = null;

      for (let dr = -1; dr <= 1; dr++) {
        for (let dc = -1; dc <= 1; dc++) {
          if (dr === 0 && dc === 0) continue;
          const nr = r + dr;
          const nc = c + dc;
          if (nr < 0 || nr >= gridSize || nc < 0 || nc >= gridSize) continue;

          const zNeighbor = elevationGrid[nr][nc];
          const distM = Math.hypot(dr * dyMeters, dc * dxMeters);
          const drop = zCurr - zNeighbor;
          const slope = drop / distM;

          if (slope > maxDropSlope) {
            maxDropSlope = slope;
            bestTarget = { r: nr, c: nc, slope };
          }
        }
      }

      flowTarget[r][c] = bestTarget;
    }
  }

  // Calculate Flow Accumulation (Upslope contributing area)
  const flowAcc: number[][] = Array(gridSize).fill(1).map(() => Array(gridSize).fill(1));

  // Sort cells from highest elevation to lowest
  const sortedCells: Array<{ r: number; c: number; elev: number }> = [];
  for (let r = 0; r < gridSize; r++) {
    for (let c = 0; c < gridSize; c++) {
      sortedCells.push({ r, c, elev: elevationGrid[r][c] });
    }
  }
  sortedCells.sort((a, b) => b.elev - a.elev);

  for (const cell of sortedCells) {
    const target = flowTarget[cell.r][cell.c];
    if (target) {
      flowAcc[target.r][target.c] += flowAcc[cell.r][cell.c];
    }
  }

  // Extract Thalweg Valley & Stream Channels
  // Threshold for channel initiation: cells where accumulated flow exceeds minimum threshold
  const minAccThreshold = Math.max(5, Math.floor(gridSize * 0.3));
  const visitedStream = Array(gridSize).fill(false).map(() => Array(gridSize).fill(false));

  // Find all channel heads (cells where flowAcc >= threshold, but none of upstream neighbors >= threshold)
  const channelHeads: Array<{ r: number; c: number; acc: number }> = [];
  for (let r = 1; r < gridSize - 1; r++) {
    for (let c = 1; c < gridSize - 1; c++) {
      if (flowAcc[r][c] >= minAccThreshold) {
        // Check if any upstream neighbor that flows to (r, c) already had flowAcc >= minAccThreshold
        let hasUpstreamThreshold = false;
        for (let dr = -1; dr <= 1; dr++) {
          for (let dc = -1; dc <= 1; dc++) {
            if (dr === 0 && dc === 0) continue;
            const nr = r + dr;
            const nc = c + dc;
            if (nr >= 0 && nr < gridSize && nc >= 0 && nc < gridSize) {
              const target = flowTarget[nr][nc];
              if (target && target.r === r && target.c === c && flowAcc[nr][nc] >= minAccThreshold) {
                hasUpstreamThreshold = true;
                break;
              }
            }
          }
          if (hasUpstreamThreshold) break;
        }

        if (!hasUpstreamThreshold) {
          channelHeads.push({ r, c, acc: flowAcc[r][c] });
        }
      }
    }
  }

  // Trace streams downstream from each channel head
  for (const head of channelHeads) {
    const rawPath: [number, number][] = [];
    const streamElevations: number[] = [];
    let currR = head.r;
    let currC = head.c;
    let maxStreamAcc = head.acc;
    let steps = 0;

    while (currR >= 0 && currR < gridSize && currC >= 0 && currC < gridSize && steps < gridSize * 2) {
      steps++;
      const pt = gridCoords[currR * gridSize + currC];
      const distFromCenter = calculateDistance(centerLat, centerLng, pt.lat, pt.lng);
      
      if (distFromCenter <= radiusMeters * 1.05) {
        rawPath.push([pt.lng, pt.lat]);
        streamElevations.push(elevationGrid[currR][currC]);
        if (flowAcc[currR][currC] > maxStreamAcc) {
          maxStreamAcc = flowAcc[currR][currC];
        }
      }

      visitedStream[currR][currC] = true;
      const next = flowTarget[currR][currC];
      if (!next) break; // Reached sink or local depression
      if (visitedStream[next.r][next.c]) {
        // Confluence with existing river channel
        const confluencePt = gridCoords[next.r * gridSize + next.c];
        rawPath.push([confluencePt.lng, confluencePt.lat]);
        streamElevations.push(elevationGrid[next.r][next.c]);
        break;
      }
      currR = next.r;
      currC = next.c;
    }

    // Only keep streams that have at least 3 points and actual length
    if (rawPath.length >= 3) {
      // Smooth the stream line using moving average
      const smoothedPath: [number, number][] = [];
      for (let i = 0; i < rawPath.length; i++) {
        if (i === 0 || i === rawPath.length - 1) {
          smoothedPath.push(rawPath[i]);
        } else {
          const pPrev = rawPath[i - 1];
          const pCurr = rawPath[i];
          const pNext = rawPath[i + 1];
          smoothedPath.push([
            Number((pPrev[0] * 0.25 + pCurr[0] * 0.5 + pNext[0] * 0.25).toFixed(6)),
            Number((pPrev[1] * 0.25 + pCurr[1] * 0.5 + pNext[1] * 0.25).toFixed(6))
          ]);
        }
      }

      const streamMinZ = Math.min(...streamElevations);
      const streamMaxZ = Math.max(...streamElevations);

      spatialFeatures.push({
        type: 'Feature',
        properties: {
          layerType: 'wadi',
          category: 'natural-obstacle',
          title: `🌊 مجرى وادٍ ومصرف سيل طبيعي (${Math.round(streamMinZ)}م - ${Math.round(streamMaxZ)}م)`,
          flowAcc: maxStreamAcc,
          color: '#0284c7',
          width: 3.5
        },
        geometry: {
          type: 'LineString',
          coordinates: smoothedPath
        }
      });
    }
  }



  // 8. Regional Sector Intelligence Context
  let sectorName = 'قطاع عمليات ميداني';
  let terrainClassificationAr = 'هضاب وتلال متوسطة الوعورة';
  let terrainClassification = 'Rolling Hills & Ridges';
  let urbanDensity: 'HIGH' | 'MEDIUM' | 'LOW' | 'SPARSE' = 'MEDIUM';
  let urbanDensityAr = 'كثافة متوسطة';
  let settlements: string[] = ['تجمعات سكنية وضواحي متصلة'];
  let soilType = 'تربة طينية صخرية جافة';
  let trafficability = 'جيدة للمركبات المدولبة والمجنزرة مع الحذر في مجاري السيول';

  if (centerLat >= 32.05 && centerLat <= 32.20 && centerLng >= 36.00 && centerLng <= 36.18) {
    sectorName = 'قاطع الزرقاء والسخنة / حوض سيل الزرقاء';
    terrainClassificationAr = 'حوض وادي صدعي واسع تحيط به تلال وهضاب حاكمة ومجرى سيل مائي';
    terrainClassification = 'River Valley Basin flanked by Ridge Plateaus';
    urbanDensity = 'HIGH';
    urbanDensityAr = 'كثافة عالية (بلدات وضواحي ممتدة)';
    settlements = ['بلدة السخنة', 'خربة قنتول', 'ضاحية البستان', 'أحياء الهاشمية', 'مدينة الزرقاء'];
    soilType = 'طمي رسوبي في بطن الوادي وصخور كلسية على الأكتاف';
    trafficability = 'ممتازة على شبكة الطرق وبطون الأودية المنبسطة (GO) مع وجود قواطع صخرية على الأكتاف';
  } else if (centerLat > 32.3 && centerLng < 36.0) {
    sectorName = 'قاطع الشمال / حوض اليرموك وإربد';
    terrainClassificationAr = 'مرتفعات جبلية تتخللها أودية سحيقة وخوانق مائية وجروف صخرية';
    terrainClassification = 'Highland Plateaus with Deep Gorges & Rock Cliffs';
    urbanDensity = 'MEDIUM';
    urbanDensityAr = 'كثافة متوسطة (قرى وضواحي ممتدة)';
    settlements = ['بلدات وضواحي إربد', 'قرى حوض اليرموك', 'مزارع ومحميات وادي الشلالة'];
    soilType = 'تربة زراعية خصبة وصخور كلسية صلبة';
    trafficability = 'مقيدة بالجروف الصخرية الحادة - تتطلب محاور طرق معبدة';
  } else if (centerLat >= 31.85 && centerLat <= 32.05 && centerLng >= 35.8 && centerLng <= 36.02) {
    sectorName = 'قاطع الوسط / إقليم العاصمة عمان';
    terrainClassificationAr = 'هضاب وكتل جبلية حضرية مكتظة بالعمران';
    terrainClassification = 'Urbanized Highland Terrain';
    urbanDensity = 'HIGH';
    urbanDensityAr = 'كثافة عالية (بيئة قتال مباني MOUT)';
    settlements = ['أحياء ومجمعات العاصمة عمان', 'المناطق التجارية والتقاطعات المرورية', 'مناطق صناعية متقدمة'];
    soilType = 'أرض معبدة وأرصفة إسمنتية مع تربة صخرية';
    trafficability = 'ممتازة على شبكة الطرق ولكنها مقيدة بالمباني والاختناقات';
  } else if (centerLng > 36.18) {
    sectorName = 'قاطع البادية والشرق / المفرق والصفاوي';
    terrainClassificationAr = 'سهول صحراوية مفتوحة وأراضي حرات بازلتية';
    terrainClassification = 'Open Desert Steppe & Basalt Harrah';
    urbanDensity = 'SPARSE';
    urbanDensityAr = 'شبه خالية / تجمعات بدوية ومزارع متفرقة';
    settlements = ['تجمعات ريفية بدوية', 'مواقع عسكرية ونقاط حدودية', 'مزارع صحراوية'];
    soilType = 'رمال وحصى بازلتية صلبة';
    trafficability = 'حركة حرة مفتوحة للدروع والآليات الصحراوية (Cross-Country GO)';
  } else if (centerLat < 31.8 && centerLng < 35.7) {
    sectorName = 'قاطع الغور / المنخفض الصدعي والبحر الميت';
    terrainClassificationAr = 'أخفـض نقطة على سطح الأرض (جروف صخرية وقواطع حادة)';
    terrainClassification = 'Rift Valley Escarpment (-400m MSL) with Rock Cliffs';
    urbanDensity = 'LOW';
    urbanDensityAr = 'كثافة منخفضة (مزارع الغور والمنتجعات)';
    settlements = ['مزارع الأغوار الجنوبية', 'مجمعات الفنادق ومصانع البوتاس', 'قرى غور الصافي'];
    soilType = 'طمي نهري وسبخات ملحية وصخور رسوبية قاسية';
    trafficability = 'جروف صخرية مانعة (NO-GO) على الحواف الغربية والشرقية';
  }

  // 9. Mobility Assessment
  let mobilityStatus: 'GO' | 'SLOW-GO' | 'NO-GO' = 'GO';
  let mobilityStatusAr = 'حركة حرة ممتازة للآليات والدروع (Cross-Country GO)';
  if (cliffPct > 12 || ruggedPct > 25 || avgSlope > 14) {
    mobilityStatus = 'NO-GO';
    mobilityStatusAr = 'مانع تضاريسي - قواطع صخرية وانحدارات حادة تعيق المناورة (NO-GO)';
  } else if (cliffPct > 3 || ruggedPct > 10 || modPct > 30 || avgSlope > 6) {
    mobilityStatus = 'SLOW-GO';
    mobilityStatusAr = 'تضاريس متموجة مقيدة تتطلب مسارب وتجهيزاً هندسياً (SLOW-GO)';
  }

  // Sort and pick top distinct cliffs
  const topCliffs = detectedCliffs
    .sort((a, b) => b.dropMeters - a.dropMeters || b.slopeDegrees - a.slopeDegrees)
    .slice(0, 8);

  const cliffSummaryText = topCliffs.length > 0
    ? `تم رصد ${detectedCliffs.length} قاطعاً وجرفاً صخرياً حاداً (أقصى سقوط ${topCliffs[0].dropMeters}م بميل ${topCliffs[0].slopeDegrees}°). تشكل هذه القواطع موانع طبيعية مانعة لحركة الآليات والدروع.`
    : `القطاع خالٍ من الجروف الصخرية العمودية، والميول العامة تتراوح بين منبسطة ومتموجة معتدلة (أقصى انحدار ${maxSlope}°).`;

  return {
    centerLat,
    centerLng,
    radiusMeters,
    areaKm2,
    sectorName,
    centerElevation: Math.round(centerElev),
    minElevation: Math.round(minElev),
    maxElevation: Math.round(maxElev),
    reliefMeters: relief,
    avgSlopeDegrees: avgSlope,
    maxSlopeDegrees: maxSlope,
    slopeDistribution: {
      flatPercentage: flatPct,
      moderatePercentage: modPct,
      steepPercentage: ruggedPct,
      cliffPercentage: cliffPct
    },
    cliffsCount: detectedCliffs.length,
    cliffs: topCliffs,
    terrainClassification,
    terrainClassificationAr,
    mobilityStatus,
    mobilityStatusAr,
    highestPoint,
    lowestPoint,
    spatialGeoJson: {
      type: 'FeatureCollection',
      features: spatialFeatures
    },
    naturalObstacles: [
      {
        name: 'القواطع والجروف الصخرية الطبيعية (Cliffs & Rock Escarpments)',
        type: 'موانع وقواطع صخرية حادة',
        description: cliffSummaryText,
        impact: topCliffs.length > 0 
          ? 'إغلاق محاور الالتفاف الجانبي وإجبار القوات المهاجمة على استخدام ممرات إجبارية (Choke Points)' 
          : 'حرية مناورة واسعة بدون عوائق صخرية عمودية'
      },
      {
        name: 'مجاري الوديان والسيول الطبيعية (Drainage Defiles)',
        type: 'مصارف مائية وأودية',
        description: `بطن الوادي عند منسوب ${minElev}م يشكل مصرفاً مائياً ومحور حركة منخفضاً، يوفر ممراً محجوباً ولكن قد يعيق المناورة عند هطول الأمطار.`,
        impact: 'محور تسلل مناسب للمشاة ومقيد للآليات الثقيلة عند تشكل الأوحال والسيول'
      },
      {
        name: 'التلال والقمم الحاكمة (Dominant Ridges & Crests)',
        type: 'أرض مسيطرة (Key Terrain)',
        description: `المرتفع الحاكم بارتفاع ${maxElev}م يشرف بشكل كامل على بطن الوادي ومحاور التحرك على مسافة ${radiusMeters / 1000} كم بفارق تضاريسي ${relief}م.`,
        impact: 'موقع استراتيجي لتمركز أسلحة الإسناد ونقاط الملاحظة والاستطلاع والسيطرة النارية'
      }
    ],
    manMadeObstacles: [
      {
        name: 'محاور الطرق الرئيسية وخطوط الإمداد (MSR)',
        type: 'بنية تحتية للمواصلات',
        description: 'الطرق المعبدة المارة بالقطاع تتيح تدفقاً سريعاً للأرتال وإعادة التزود اللوجستي.',
        tacticalNote: 'وجوب إعداد خطط قطع وغلق هندسي سريع للتحكم بحركة العدو'
      },
      {
        name: 'المنشآت والكتل العمرانية',
        type: 'موانع مبنية وسواتر خرسانية',
        description: `البلدات والتجمعات السكانية (${settlements.slice(0, 3).join('، ')}) توفر سواتر ممتازة ضد النيران والشظايا.`,
        tacticalNote: 'تتطلب عمليات تطهير متقدمة وحماية للمدنيين'
      },
      {
        name: 'أبراج وشبكات نقل الطاقة والاتصالات',
        type: 'عوائق صناعية بارزة',
        description: 'خطوط الضغط العالي والاتصالات تشكل خطراً على الطيران المنخفض ومواقع محتملة للتشويش.',
        tacticalNote: 'تجنب تحديد مهابط الطيران بالقرب من مسارات خطوط الطاقة'
      }
    ],
    urbanAndDemographics: {
      density: urbanDensity,
      densityAr: urbanDensityAr,
      settlements,
      moutComplexity: urbanDensity === 'HIGH' ? 'بيئة قتال مباني حضرية متقدمة (MOUT)' : 'بيئة قتال ريفية مفتوحة تتخللها مزارع متفرقة',
      collateralDamageRisk: urbanDensity === 'HIGH' ? 'عالية وتتطلب مراعاة قواعد الاشتباك الدقيقة' : 'منخفضة مع حماية الممتلكات المدنية'
    },
    landCoverAndCover: {
      soilType,
      trafficability,
      concealmentRating: urbanDensity === 'HIGH' || relief > 100 ? 'ممتازة خلف التلال والكتل الخرسانية' : 'متوسطة إلى مكشوفة في الأراضي المفتوحة',
      airObservationExposure: urbanDensity === 'HIGH' ? 'مختلطة (المباني توفر إخفاء جزئياً)' : 'عالية للرصد الجوي والمسيرات نهاراً'
    },
    oakocAssessment: {
      obstacles: `القطاع يضم ${detectedCliffs.length > 0 ? `${detectedCliffs.length} قاطعاً وجرفاً صخرياً و` : ''}مجرى وادٍ رئيسي عند منسوب ${minElev}م، مع انحدارات تصل إلى ${maxSlope}°، مما يشكل محاور حركة طبيعية محددة.`,
      avenuesOfApproach: 'محور الاقتراب الرئيسي يتبع بطن الوادي وشبكة الطرق الرئيسية المنبسطة، مع تفادي الجروف الصخرية الحادة.',
      keyTerrain: `القمة الحاكمة بارتفاع ${maxElev}م توفر ميزة كشف وهيمنة نارية حاسمة على كامل القطاع.`,
      observationAndFieldsOfFire: `حقول الرماية مكشوفة ومثالية بنسبة ${flatPct}% في الأراضي المنبسطة، ومحجوبة داخل بطون الأودية وخلف الجروف.`,
      coverAndConcealment: 'تتوفر السواتر التضاريسية في المنحدرات الخلفية (Reverse Slopes) وأسفل الجروف الصخرية وداخل الكتل العمرانية.'
    },
    tacticalRecommendations: [
      `نشر مرصد استطلاع رئيسي ورادار كشف أرضي على المرتفع الحاكم (${maxElev}م) لتأمين إنذار مبكر ومراقبة القطاع.`,
      detectedCliffs.length > 0 
        ? `استغلال القواطع والجروف الصخرية المكتشفة كحواجز طبيعية موجهة لحصر تقدم العدو داخل المخانق المستهدفة بنيران الأسلحة الثقيلة.` 
        : `تجهيز نقاط كمائن وموانع هندسية وحقول ألغام موجهة في مخانق بطن الوادي (${minElev}م) لصد أي تقدم مدرع.`,
      `تأمين الطرق الرئيسية والمحاور الحيوية بنقاط سيطرة وحراسة ثابتة لمنع عمليات الالتفاف السريعة.`,
      `استغلال المناطق المنبسطة (${flatPct}%) لإنشاء مهابط المروحيات (HLZ) ومناطق الإخلاء الطبي والتجمع اللوجستي.`
    ]
  };
}

export interface IPBOverlayFeatureSet {
  aoBoundary: any;
  layers: {
    layer_1: any[]; // Mountain Scarps & Steep Slopes (>20°)
    layer_2: any[]; // Basalt / Rocky Boulder Fields & Cliffs
    layer_3: any[]; // Urban Built-up Centers & Demographics
    layer_4: any[]; // Wadis, Drainage Corridors & Waterways
    layer_5: any[]; // Minefields & Man-made Choke Point Barriers
    layer_6: any[]; // MCOO Combined Obstacle Grid (Green, Yellow, Red)
    layer_7: any[]; // Avenues of Approach / Armor Mobility Corridors
    layer_8: any[]; // Key Terrain Peaks (K) & Command Summits
    layer_9: any[]; // Dead Ground & Viewshed Shadowed Zones
    layer_10: any[]; // Threat SITTEMP, Kill Zones & Artillery Arcs
  };
  metrics: {
    maxElev: number;
    minElev: number;
    relief: number;
    maxSlope: number;
    steepSlopePct: number;
    ruggedPct: number;
    unrestrictedPct: number;
    restrictedPct: number;
    severelyRestrictedPct: number;
    cliffCount: number;
    keyTerrainSummit: { lat: number; lng: number; elevation: number };
    killZoneCenter: { lat: number; lng: number };
  };
  allFeatures: any[];
}

/**
 * Real Mathematical IPB Digital Overlay Engine
 * محرك استخبارات إعداد ساحة المعركة الرياضي عالي الدقة:
 * يشتق الشفافات الـ 10 الحقيقية مباشرة من شبكة الارتفاعات الرقمية DEM، خطوط الكنتور، وحسابات الميول الدقيقة
 */
export async function calculateRealIPBOverlays(
  centerLat: number,
  centerLng: number,
  radiusMeters: number = 6000
): Promise<IPBOverlayFeatureSet> {
  const centerElev = await sampleElevationAt(centerLat, centerLng);
  
  // Matrix sampling resolution (25x25 grid = 625 nodes)
  const gridSize = 25;
  const R_earth = 6371000;
  const dLatTotal = (radiusMeters / R_earth) * (180 / Math.PI);
  const dLngTotal = (radiusMeters / (R_earth * Math.cos((centerLat * Math.PI) / 180))) * (180 / Math.PI);

  const minLat = centerLat - dLatTotal;
  const maxLat = centerLat + dLatTotal;
  const minLng = centerLng - dLngTotal;
  const maxLng = centerLng + dLngTotal;

  const gridCoords: Array<{ r: number; c: number; lat: number; lng: number }> = [];
  for (let r = 0; r < gridSize; r++) {
    const lat = Number((minLat + (r / (gridSize - 1)) * (maxLat - minLat)).toFixed(6));
    for (let c = 0; c < gridSize; c++) {
      const lng = Number((minLng + (c / (gridSize - 1)) * (maxLng - minLng)).toFixed(6));
      gridCoords.push({ r, c, lat, lng });
    }
  }

  const sampleResults = await sampleElevationsBatch(gridCoords, 13);
  const elevationGrid: number[][] = Array(gridSize).fill(0).map(() => Array(gridSize).fill(0));
  gridCoords.forEach((pt, idx) => {
    elevationGrid[pt.r][pt.c] = sampleResults[idx];
  });

  // Calculate Extremes (Key Terrain Peak & Lowest Valley Point)
  let maxElev = -9999;
  let minElev = 9999;
  let highestPoint = { lat: centerLat, lng: centerLng, elevation: centerElev };
  let lowestPoint = { lat: centerLat, lng: centerLng, elevation: centerElev };

  gridCoords.forEach((p, idx) => {
    const el = sampleResults[idx];
    const dist = calculateDistance(centerLat, centerLng, p.lat, p.lng);
    if (dist <= radiusMeters * 1.05) {
      if (el > maxElev) {
        maxElev = el;
        highestPoint = { lat: p.lat, lng: p.lng, elevation: Math.round(el) };
      }
      if (el < minElev) {
        minElev = el;
        lowestPoint = { lat: p.lat, lng: p.lng, elevation: Math.round(el) };
      }
    }
  });

  const relief = Math.round(maxElev - minElev);

  // Initialize Layer Feature Buckets
  const layer_1: any[] = []; // Mountains & Steep Slopes (>20°)
  const layer_2: any[] = []; // Basalt / Rocky Escarpments
  const layer_3: any[] = []; // Urban Demographics
  const layer_4: any[] = []; // Wadis & Waterways
  const layer_5: any[] = []; // Minefield Barriers
  const layer_6: any[] = []; // MCOO Trafficability Grid
  const layer_7: any[] = []; // Avenues of Approach
  const layer_8: any[] = []; // Key Terrain
  const layer_9: any[] = []; // Dead Space / Viewshed
  const layer_10: any[] = []; // SITTEMP & Kill Zones

  const dxMeters = (2 * radiusMeters) / (gridSize - 1);
  const dyMeters = (2 * radiusMeters) / (gridSize - 1);

  let flatCells = 0;
  let slowGoCells = 0;
  let severeCells = 0;
  let cliffCount = 0;
  let maxCalculatedSlope = 0;

  // Track valley drainage nodes for Layer 4
  const valleyNodes: Array<{ lat: number; lng: number; elevation: number; r: number; c: number }> = [];
  // Track choke point cells for Layer 5
  const chokePoints: Array<{ lat: number; lng: number; elevation: number }> = [];

  for (let r = 0; r < gridSize - 1; r++) {
    // Find local elevation minimum in this row for valley tracing
    let rowMinElev = 9999;
    let rowMinNode: any = null;

    for (let c = 0; c < gridSize - 1; c++) {
      const p1Coord = gridCoords[r * gridSize + c];
      const p2Coord = gridCoords[r * gridSize + (c + 1)];
      const p3Coord = gridCoords[(r + 1) * gridSize + (c + 1)];
      const p4Coord = gridCoords[(r + 1) * gridSize + c];

      const cellCenterLat = (p1Coord.lat + p3Coord.lat) / 2;
      const cellCenterLng = (p1Coord.lng + p3Coord.lng) / 2;
      const distFromCenter = calculateDistance(centerLat, centerLng, cellCenterLat, cellCenterLng);

      if (distFromCenter > radiusMeters * 1.02) continue;

      const z00 = elevationGrid[r][c];
      const z01 = elevationGrid[r][c + 1];
      const z10 = elevationGrid[r + 1][c];
      const z11 = elevationGrid[r + 1][c + 1];

      const cellElevAvg = Math.round((z00 + z01 + z10 + z11) / 4);

      if (cellElevAvg < rowMinElev) {
        rowMinElev = cellElevAvg;
        rowMinNode = { lat: cellCenterLat, lng: cellCenterLng, elevation: cellElevAvg, r, c };
      }

      // Central difference slope gradient
      const dz_dx = ((z01 + z11) - (z00 + z10)) / (2 * dxMeters);
      const dz_dy = ((z10 + z11) - (z00 + z01)) / (2 * dyMeters);
      const slopeRad = Math.atan(Math.sqrt(dz_dx * dz_dx + dz_dy * dz_dy));
      const slopeDeg = Math.round((slopeRad * (180 / Math.PI)) * 10) / 10;

      if (slopeDeg > maxCalculatedSlope) maxCalculatedSlope = slopeDeg;

      const localDrop = Math.round(Math.max(z00, z01, z10, z11) - Math.min(z00, z01, z10, z11));

      const p1 = [p1Coord.lng, p1Coord.lat];
      const p2 = [p2Coord.lng, p2Coord.lat];
      const p3 = [p3Coord.lng, p3Coord.lat];
      const p4 = [p4Coord.lng, p4Coord.lat];
      const polyCoords = [[p1, p2, p3, p4, p1]];

      // Classification
      if (slopeDeg >= 20) {
        severeCells++;
        if (slopeDeg >= 28) cliffCount++;

        // 1. Layer 1: Real Mountain Scarps & Steep Slope Polygon
        layer_1.push({
          type: 'Feature',
          properties: {
            layer: 'layer_1',
            color: '#dc2626',
            strokeColor: '#991b1b',
            opacity: 0.78,
            slope: slopeDeg,
            elevation: cellElevAvg,
            label: slopeDeg >= 28 ? `▲ جرف صخري حاد (${slopeDeg}° / ${cellElevAvg}م)` : `▲ انحدار شديد (${slopeDeg}°)`
          },
          geometry: { type: 'Polygon', coordinates: polyCoords }
        });

        // 2. Layer 2: Basalt / Rocky Escarpments & Mountain Crest Cliffs
        if (slopeDeg >= 18 || localDrop >= 18) {
          layer_2.push({
            type: 'Feature',
            properties: {
              layer: 'layer_2',
              color: '#b45309',
              strokeColor: '#78350f',
              opacity: 0.78,
              drop: localDrop,
              label: `■ مقطع صخري وعر (${localDrop}م هبوط / ${slopeDeg}° ميل)`
            },
            geometry: { type: 'Polygon', coordinates: polyCoords }
          });
        }

        // 6. Layer 6: MCOO - Severely Restricted Red (Only on steep slopes)
        layer_6.push({
          type: 'Feature',
          properties: {
            layer: 'layer_6',
            color: '#dc2626',
            strokeColor: '#b91c1c',
            opacity: 0.28
          },
          geometry: { type: 'Polygon', coordinates: polyCoords }
        });
      } else if (slopeDeg >= 9) {
        slowGoCells++;
        // 6. Layer 6: MCOO - Restricted Amber (Only on rugged slopes)
        layer_6.push({
          type: 'Feature',
          properties: {
            layer: 'layer_6',
            color: '#d97706',
            strokeColor: '#b45309',
            opacity: 0.20
          },
          geometry: { type: 'Polygon', coordinates: polyCoords }
        });

        // Detect Choke Points where passable ground meets steep walls
        if (c > 2 && c < gridSize - 3) {
          const leftSlope = Math.abs(elevationGrid[r][c - 2] - z00) / (2 * dxMeters);
          const rightSlope = Math.abs(elevationGrid[r][c + 2] - z01) / (2 * dxMeters);
          if (leftSlope > 0.35 && rightSlope > 0.35) {
            chokePoints.push({ lat: cellCenterLat, lng: cellCenterLng, elevation: cellElevAvg });
          }
        }
      } else {
        flatCells++;
        // Unrestricted flat land remains unshaded and clean so the map shows naturally
      }
    }

    if (rowMinNode) {
      valleyNodes.push(rowMinNode);
    }
  }

  const totalCells = Math.max(1, flatCells + slowGoCells + severeCells);
  const unrestrictedPct = Math.round((flatCells / totalCells) * 100);
  const restrictedPct = Math.round((slowGoCells / totalCells) * 100);
  const severelyRestrictedPct = Math.round((severeCells / totalCells) * 100);

  // 3. Layer 3: Urban Built-up Centers (Rendered directly by OSM/Overture style layers)

  // 4. Layer 4: Real Valley Drainage Corridors & Wadis (مجاري السيول والأودية الحقيقية)
  if (valleyNodes.length >= 3) {
    const wadiLineCoords = valleyNodes.map(n => [n.lng, n.lat]);
    layer_4.push({
      type: 'Feature',
      properties: {
        layer: 'layer_4',
        color: '#0284c7',
        strokeColor: '#0369a1',
        opacity: 0.85,
        label: `≈≈ مجرى وادٍ وسيل رئيسي (${lowestPoint.elevation}م)`
      },
      geometry: { type: 'LineString', coordinates: wadiLineCoords }
    });
  }

  // 5. Layer 5: Tactical Obstacles & Retaining Barriers (Populated if obstacles exist)

  // 7. Layer 7: Real Approach Corridor & Mobility (المقترب التعبوي الرئيسي AA-1)
  let approachLengthKm = Math.round((radiusMeters * 1.4 / 1000) * 10) / 10;
  let approachDirection = 'الغرب نحو الشرق';
  let approachAzimuth = 85;
  let approachWidthM = Math.round(dxMeters * 3);
  let approachCapacity = 'كتيبة مدرعة بنسق رتل';

  if (valleyNodes.length >= 3) {
    const approachCoords = valleyNodes.slice(0, Math.floor(valleyNodes.length * 0.8)).map(n => [n.lng, n.lat]);
    
    // Calculate real length
    let lenKm = 0;
    for (let i = 0; i < approachCoords.length - 1; i++) {
      const p1 = approachCoords[i];
      const p2 = approachCoords[i + 1];
      lenKm += Math.hypot((p2[1] - p1[1]) * 111.139, (p2[0] - p1[0]) * 111.139 * Math.cos(centerLat * Math.PI / 180));
    }
    if (lenKm > 0.5) approachLengthKm = Math.round(lenKm * 10) / 10;

    // Calculate azimuth from start to end
    const startP = approachCoords[0];
    const endP = approachCoords[approachCoords.length - 1];
    const dLngRad = (endP[0] - startP[0]) * (Math.PI / 180);
    const lat1Rad = startP[1] * (Math.PI / 180);
    const lat2Rad = endP[1] * (Math.PI / 180);
    const y = Math.sin(dLngRad) * Math.cos(lat2Rad);
    const x = Math.cos(lat1Rad) * Math.sin(lat2Rad) - Math.sin(lat1Rad) * Math.cos(lat2Rad) * Math.cos(dLngRad);
    const azDeg = Math.round((Math.atan2(y, x) * 180 / Math.PI + 360) % 360);
    approachAzimuth = azDeg;
    
    if (azDeg >= 337.5 || azDeg < 22.5) approachDirection = 'الجنوب نحو الشمال';
    else if (azDeg >= 22.5 && azDeg < 67.5) approachDirection = 'الجنوب الغربي نحو الشمال الشرقي';
    else if (azDeg >= 67.5 && azDeg < 112.5) approachDirection = 'الغرب نحو الشرق';
    else if (azDeg >= 112.5 && azDeg < 157.5) approachDirection = 'الشمال الغربي نحو الجنوب الشرقي';
    else if (azDeg >= 157.5 && azDeg < 202.5) approachDirection = 'الشمال نحو الجنوب';
    else if (azDeg >= 202.5 && azDeg < 247.5) approachDirection = 'الشمال الشرقي نحو الجنوب الغربي';
    else if (azDeg >= 247.5 && azDeg < 292.5) approachDirection = 'الشرق نحو الغرب';
    else approachDirection = 'الجنوب الشرقي نحو الشمال الغربي';

    if (approachWidthM > 1200) approachCapacity = 'لواء مدرع بنسق كتائب';
    else if (approachWidthM > 600) approachCapacity = 'كتيبة مدرعة بنسق رتل سرية';
    else approachCapacity = 'سرية دبابات بنسق رتل واحد';

    layer_7.push({
      type: 'Feature',
      properties: {
        layer: 'layer_7',
        color: '#2563eb',
        strokeColor: '#1d4ed8',
        opacity: 0.90,
        label: `➔➔ المقترب التعبوي الرئيسي (AA-1: بطول ${approachLengthKm} كم / عرض ${approachWidthM}م / سعة ${approachCapacity})`
      },
      geometry: { type: 'LineString', coordinates: approachCoords }
    });
  }

  // 8. Layer 8: Real Key Terrain (الأرض الحيوية والقمم الحاكمة الحقيقية)
  layer_8.push({
    type: 'Feature',
    properties: {
      layer: 'layer_8',
      color: '#9333ea',
      strokeColor: '#6b21a8',
      opacity: 0.95,
      label: `★ قمة حاكمة: ${highestPoint.elevation}م [Key Terrain]`
    },
    geometry: {
      type: 'Point',
      coordinates: [highestPoint.lng, highestPoint.lat]
    }
  });

  // 9. Layer 9: Real DEM Line of Sight & Viewshed Shadows from Key Terrain Peak
  // Soft, continuous terrain shadow masking without harsh black borders
  const summitR = Math.min(gridSize - 1, Math.max(0, Math.round(((centerLat + dLatTotal - highestPoint.lat) / (2 * dLatTotal)) * (gridSize - 1))));
  const summitC = Math.min(gridSize - 1, Math.max(0, Math.round(((highestPoint.lng - (centerLng - dLngTotal)) / (2 * dLngTotal)) * (gridSize - 1))));
  const summitElev = highestPoint.elevation;

  for (let r = 0; r < gridSize; r += 2) {
    for (let c = 0; c < gridSize; c += 2) {
      if (Math.abs(r - summitR) <= 1 && Math.abs(c - summitC) <= 1) continue;
      const distCells = Math.hypot(r - summitR, c - summitC);
      if (distCells === 0) continue;

      const targetElev = elevationGrid[r][c];
      const angleToTarget = (targetElev - summitElev) / distCells;

      let isShadowed = false;
      const steps = Math.floor(distCells);
      for (let s = 1; s < steps; s++) {
        const interR = Math.round(summitR + (r - summitR) * (s / distCells));
        const interC = Math.round(summitC + (c - summitC) * (s / distCells));
        if (interR >= 0 && interR < gridSize && interC >= 0 && interC < gridSize) {
          const interElev = elevationGrid[interR][interC];
          const angleToInter = (interElev - summitElev) / s;
          if (angleToInter > angleToTarget + 0.8) {
            isShadowed = true;
            break;
          }
        }
      }

      if (isShadowed) {
        const cellLat = centerLat + dLatTotal - (r / (gridSize - 1)) * (2 * dLatTotal);
        const cellLng = centerLng - dLngTotal + (c / (gridSize - 1)) * (2 * dLngTotal);
        const distM = Math.hypot(
          (cellLat - centerLat) * (R_earth * Math.PI / 180),
          (cellLng - centerLng) * (R_earth * Math.cos(centerLat * Math.PI / 180) * Math.PI / 180)
        );

        if (distM <= radiusMeters * 0.95) {
          const hLat = (dLatTotal / (gridSize - 1)) * 1.05;
          const hLng = (dLngTotal / (gridSize - 1)) * 1.05;
          layer_9.push({
            type: 'Feature',
            properties: {
              layer: 'layer_9',
              color: '#475569',
              strokeColor: '#334155',
              opacity: 0.18
            },
            geometry: {
              type: 'Polygon',
              coordinates: [[
                [cellLng - hLng, cellLat - hLat],
                [cellLng + hLng, cellLat - hLat],
                [cellLng + hLng, cellLat + hLat],
                [cellLng - hLng, cellLat + hLat],
                [cellLng - hLng, cellLat - hLat]
              ]]
            }
          });
        }
      }
    }
  }

  // 10. Layer 10: SITTEMP & Engagement Area (منطقة التقتيل الرئيسية Kill Zone Alpha)
  // Placed at the natural choke point near lowest valley within line of sight of key terrain
  const killZoneLat = (highestPoint.lat + lowestPoint.lat * 2) / 3;
  const killZoneLng = (highestPoint.lng + lowestPoint.lng * 2) / 3;
  const kzR = (radiusMeters * 0.14 / R_earth) * (180 / Math.PI);
  const kzRLng = kzR / Math.cos((killZoneLat * Math.PI) / 180);

  const kzCoords: [number, number][] = [];
  for (let i = 0; i <= 24; i++) {
    const a = (i / 24) * 2 * Math.PI;
    kzCoords.push([killZoneLng + kzRLng * Math.sin(a), killZoneLat + kzR * Math.cos(a)]);
  }

  // Direct Fire Line from Key Terrain Peak to Kill Zone Center
  const fireLineCoords = [
    [highestPoint.lng, highestPoint.lat],
    [killZoneLng, killZoneLat]
  ];
  
  const fireDistKm = Math.round(Math.hypot(
    (killZoneLat - highestPoint.lat) * 111.139,
    (killZoneLng - highestPoint.lng) * 111.139 * Math.cos(centerLat * Math.PI / 180)
  ) * 10) / 10;
  
  const fireAzimuth = Math.round((Math.atan2(
    (killZoneLng - highestPoint.lng) * Math.cos(centerLat * Math.PI / 180),
    killZoneLat - highestPoint.lat
  ) * 180 / Math.PI + 360) % 360);

  // Add Direct Fire Line
  layer_10.push({
    type: 'Feature',
    properties: {
      layer: 'layer_10',
      color: '#dc2626',
      strokeColor: '#991b1b',
      opacity: 0.95,
      label: `⚡ خط نار ورصد مباشر مسيطر (أزيموث ${fireAzimuth}° / مدى ${fireDistKm} كم)`
    },
    geometry: {
      type: 'LineString',
      coordinates: fireLineCoords
    }
  });

  // Add Kill Zone Area
  layer_10.push({
    type: 'Feature',
    properties: {
      layer: 'layer_10',
      color: '#ef4444',
      strokeColor: '#b91c1c',
      opacity: 0.75,
      label: '⚔️ منطقة التقتيل الرئيسية (KILL ZONE ALPHA)'
    },
    geometry: {
      type: 'Polygon',
      coordinates: [kzCoords]
    }
  });

  // AO Boundary Ring Feature
  const aoCoords: [number, number][] = [];
  for (let i = 0; i <= 48; i++) {
    const a = (i / 48) * 2 * Math.PI;
    aoCoords.push([centerLng + dLngTotal * Math.sin(a), centerLat + dLatTotal * Math.cos(a)]);
  }

  const aoBoundary = {
    type: 'Feature',
    properties: {
      color: '#38bdf8',
      strokeColor: '#0284c7',
      opacity: 0.08,
      label: `قاطع العمليات (${Math.round(radiusMeters / 1000)} كم / منسوب ${centerElev}م)`
    },
    geometry: {
      type: 'Polygon',
      coordinates: [aoCoords]
    }
  };

  const allFeatures = [
    aoBoundary,
    ...layer_1,
    ...layer_2,
    ...layer_3,
    ...layer_4,
    ...layer_5,
    ...layer_6,
    ...layer_7,
    ...layer_8,
    ...layer_9,
    ...layer_10
  ];

  return {
    aoBoundary,
    layers: {
      layer_1,
      layer_2,
      layer_3,
      layer_4,
      layer_5,
      layer_6,
      layer_7,
      layer_8,
      layer_9,
      layer_10
    },
    metrics: {
      maxElev,
      minElev,
      relief,
      maxSlope: maxCalculatedSlope,
      steepSlopePct: severelyRestrictedPct,
      ruggedPct: restrictedPct,
      unrestrictedPct,
      restrictedPct,
      severelyRestrictedPct,
      cliffCount,
      keyTerrainSummit: highestPoint,
      killZoneCenter: { lat: killZoneLat, lng: killZoneLng },
      approachLengthKm,
      approachWidthM,
      approachDirection,
      approachAzimuth,
      approachCapacity,
      fireDistKm,
      fireAzimuth
    },
    allFeatures
  };
}