wing-ops/frontend/src/common/utils/geo.ts

398 lines
14 KiB
TypeScript
Executable File

import type { BoomLine, BoomLineCoord, AlgorithmSettings, ContainmentResult } from '../types/boomLine'
const DEG2RAD = Math.PI / 180
const RAD2DEG = 180 / Math.PI
const EARTH_RADIUS = 6371000 // meters
// ============================================================
// Convex Hull + 면적 계산
// ============================================================
interface LatLon { lat: number; lon: number }
/** Convex Hull (Graham Scan) — 입자 좌표 배열 → 외곽 다각형 좌표 반환 */
export function convexHull(points: LatLon[]): LatLon[] {
if (points.length < 3) return [...points]
// 가장 아래(lat 최소) 점 찾기 (동일하면 lon 최소)
const sorted = [...points].sort((a, b) => a.lat - b.lat || a.lon - b.lon)
const pivot = sorted[0]
// pivot 기준 극각으로 정렬
const rest = sorted.slice(1).sort((a, b) => {
const angleA = Math.atan2(a.lon - pivot.lon, a.lat - pivot.lat)
const angleB = Math.atan2(b.lon - pivot.lon, b.lat - pivot.lat)
if (angleA !== angleB) return angleA - angleB
// 같은 각도면 거리 순
const dA = (a.lat - pivot.lat) ** 2 + (a.lon - pivot.lon) ** 2
const dB = (b.lat - pivot.lat) ** 2 + (b.lon - pivot.lon) ** 2
return dA - dB
})
const hull: LatLon[] = [pivot]
for (const p of rest) {
while (hull.length >= 2) {
const a = hull[hull.length - 2]
const b = hull[hull.length - 1]
const cross = (b.lon - a.lon) * (p.lat - a.lat) - (b.lat - a.lat) * (p.lon - a.lon)
if (cross <= 0) hull.pop()
else break
}
hull.push(p)
}
return hull
}
/** Shoelace 공식으로 다각형 면적 계산 (km²) — 위경도 좌표를 미터 변환 후 계산 */
// export function polygonAreaKm2(polygon: LatLon[]): number {
// if (polygon.length < 3) return 0
// // 중심 기준 위경도 → 미터 변환
// const centerLat = polygon.reduce((s, p) => s + p.lat, 0) / polygon.length
// const mPerDegLat = 111320
// const mPerDegLon = 111320 * Math.cos(centerLat * DEG2RAD)
// const pts = polygon.map(p => ({
// x: (p.lon - polygon[0].lon) * mPerDegLon,
// y: (p.lat - polygon[0].lat) * mPerDegLat,
// }))
// // Shoelace
// let area = 0
// for (let i = 0; i < pts.length; i++) {
// const j = (i + 1) % pts.length
// area += pts[i].x * pts[j].y
// area -= pts[j].x * pts[i].y
// }
// return Math.abs(area) / 2 / 1_000_000 // m² → km²
// }
/** 오일 궤적 입자 → Convex Hull 외곽 다각형 + 면적 + 둘레 계산 */
export function analyzeSpillPolygon(trajectory: LatLon[]): {
hull: LatLon[]
areaKm2: number
perimeterKm: number
particleCount: number
} {
if (trajectory.length < 3) {
return { hull: [], areaKm2: 0, perimeterKm: 0, particleCount: trajectory.length }
}
const hull = convexHull(trajectory)
const areaKm2 = polygonAreaKm2(hull)
// 둘레 계산
let perimeter = 0
for (let i = 0; i < hull.length; i++) {
const j = (i + 1) % hull.length
perimeter += haversineDistance(
{ lat: hull[i].lat, lon: hull[i].lon },
{ lat: hull[j].lat, lon: hull[j].lon },
)
}
return {
hull,
areaKm2,
perimeterKm: perimeter / 1000,
particleCount: trajectory.length,
}
}
/** 두 좌표 간 Haversine 거리 (m) */
export function haversineDistance(p1: BoomLineCoord, p2: BoomLineCoord): number {
const dLat = (p2.lat - p1.lat) * DEG2RAD
const dLon = (p2.lon - p1.lon) * DEG2RAD
const a =
Math.sin(dLat / 2) ** 2 +
Math.cos(p1.lat * DEG2RAD) * Math.cos(p2.lat * DEG2RAD) * Math.sin(dLon / 2) ** 2
return EARTH_RADIUS * 2 * Math.atan2(Math.sqrt(a), Math.sqrt(1 - a))
}
/** 두 좌표 간 방위각 (degrees, 0=북, 시계방향) */
export function computeBearing(from: BoomLineCoord, to: BoomLineCoord): number {
const dLon = (to.lon - from.lon) * DEG2RAD
const fromLat = from.lat * DEG2RAD
const toLat = to.lat * DEG2RAD
const y = Math.sin(dLon) * Math.cos(toLat)
const x = Math.cos(fromLat) * Math.sin(toLat) - Math.sin(fromLat) * Math.cos(toLat) * Math.cos(dLon)
return ((Math.atan2(y, x) * RAD2DEG) + 360) % 360
}
/** 기준점에서 방위각+거리로 목적점 계산 */
export function destinationPoint(origin: BoomLineCoord, bearingDeg: number, distanceM: number): BoomLineCoord {
const brng = bearingDeg * DEG2RAD
const lat1 = origin.lat * DEG2RAD
const lon1 = origin.lon * DEG2RAD
const d = distanceM / EARTH_RADIUS
const lat2 = Math.asin(Math.sin(lat1) * Math.cos(d) + Math.cos(lat1) * Math.sin(d) * Math.cos(brng))
const lon2 = lon1 + Math.atan2(Math.sin(brng) * Math.sin(d) * Math.cos(lat1), Math.cos(d) - Math.sin(lat1) * Math.sin(lat2))
return { lat: lat2 * RAD2DEG, lon: lon2 * RAD2DEG }
}
/** 폴리라인 총 길이 (m) */
export function computePolylineLength(coords: BoomLineCoord[]): number {
let total = 0
for (let i = 1; i < coords.length; i++) {
total += haversineDistance(coords[i - 1], coords[i])
}
return total
}
/** 두 선분 교차 판정 (2D) */
export function segmentsIntersect(
a1: BoomLineCoord, a2: BoomLineCoord,
b1: BoomLineCoord, b2: BoomLineCoord
): boolean {
const cross = (o: BoomLineCoord, a: BoomLineCoord, b: BoomLineCoord) =>
(a.lon - o.lon) * (b.lat - o.lat) - (a.lat - o.lat) * (b.lon - o.lon)
const d1 = cross(b1, b2, a1)
const d2 = cross(b1, b2, a2)
const d3 = cross(a1, a2, b1)
const d4 = cross(a1, a2, b2)
if (((d1 > 0 && d2 < 0) || (d1 < 0 && d2 > 0)) &&
((d3 > 0 && d4 < 0) || (d3 < 0 && d4 > 0))) {
return true
}
// Collinear cases
const onSegment = (p: BoomLineCoord, q: BoomLineCoord, r: BoomLineCoord) =>
Math.min(p.lon, r.lon) <= q.lon && q.lon <= Math.max(p.lon, r.lon) &&
Math.min(p.lat, r.lat) <= q.lat && q.lat <= Math.max(p.lat, r.lat)
if (d1 === 0 && onSegment(b1, a1, b2)) return true
if (d2 === 0 && onSegment(b1, a2, b2)) return true
if (d3 === 0 && onSegment(a1, b1, a2)) return true
if (d4 === 0 && onSegment(a1, b2, a2)) return true
return false
}
/** AI 자동 배치 생성 */
export function generateAIBoomLines(
trajectory: Array<{ lat: number; lon: number; time: number; particle?: number }>,
incident: BoomLineCoord,
settings: AlgorithmSettings
): BoomLine[] {
if (trajectory.length === 0) return []
// 1. 최종 시간 입자들의 중심점 → 주요 확산 방향
const maxTime = Math.max(...trajectory.map(p => p.time))
const finalPoints = trajectory.filter(p => p.time === maxTime)
if (finalPoints.length === 0) return []
const centroid: BoomLineCoord = {
lat: finalPoints.reduce((s, p) => s + p.lat, 0) / finalPoints.length,
lon: finalPoints.reduce((s, p) => s + p.lon, 0) / finalPoints.length,
}
const mainBearing = computeBearing(incident, centroid)
const totalDist = haversineDistance(incident, centroid)
// 입자 분산 폭 계산 (최종 시간 기준)
// eslint-disable-next-line @typescript-eslint/no-unused-vars
const perpBearing = (mainBearing + 90) % 360
let maxSpread = 0
for (const p of finalPoints) {
const bearing = computeBearing(incident, p)
const dist = haversineDistance(incident, p)
const perpDist = dist * Math.abs(Math.sin((bearing - mainBearing) * DEG2RAD))
if (perpDist > maxSpread) maxSpread = perpDist
}
const correction = settings.currentOrthogonalCorrection
const waveFactor = settings.waveHeightCorrectionFactor
// 2. 1차 방어선 (긴급) - 30% 지점, 해류 직교
const dist1 = totalDist * 0.3
const center1 = destinationPoint(incident, mainBearing, dist1)
const halfLen1 = Math.max(600, maxSpread * 0.4)
const line1Left = destinationPoint(center1, (mainBearing + 90 + correction) % 360, halfLen1)
const line1Right = destinationPoint(center1, (mainBearing - 90 + correction + 360) % 360, halfLen1)
const coords1 = [line1Left, line1Right]
const eff1 = Math.min(95, Math.max(60, 92 / waveFactor))
// 3. 2차 방어선 (중요) - 50% 지점, U형 3포인트
const dist2 = totalDist * 0.5
const center2 = destinationPoint(incident, mainBearing, dist2)
const halfLen2 = Math.max(450, maxSpread * 0.5)
const line2Left = destinationPoint(center2, (mainBearing + 90 + correction) % 360, halfLen2)
const line2Right = destinationPoint(center2, (mainBearing - 90 + correction + 360) % 360, halfLen2)
const line2Front = destinationPoint(center2, mainBearing, halfLen2 * 0.5)
const coords2 = [line2Left, line2Front, line2Right]
const eff2 = Math.min(90, Math.max(55, 85 / waveFactor))
// 4. 3차 방어선 (보통) - 80% 지점, 해안선 평행
const dist3 = totalDist * 0.8
const center3 = destinationPoint(incident, mainBearing, dist3)
const halfLen3 = Math.max(350, maxSpread * 0.6)
const line3Left = destinationPoint(center3, (mainBearing + 90) % 360, halfLen3)
const line3Right = destinationPoint(center3, (mainBearing - 90 + 360) % 360, halfLen3)
const coords3 = [line3Left, line3Right]
const eff3 = Math.min(85, Math.max(50, 78 / waveFactor))
const boomLines: BoomLine[] = [
{
id: `boom-ai-1-${Date.now()}`,
name: '1차 방어선 (고강도 차단형)',
priority: 'CRITICAL',
type: '고강도 차단형',
coords: coords1,
length: computePolylineLength(coords1),
angle: (mainBearing + 90 + correction) % 360,
efficiency: Math.round(eff1),
status: 'PLANNED',
},
{
id: `boom-ai-2-${Date.now()}`,
name: '2차 방어선 (외해용 포위망)',
priority: 'HIGH',
type: '외해용 중형 포위망',
coords: coords2,
length: computePolylineLength(coords2),
angle: (mainBearing + 90 + correction) % 360,
efficiency: Math.round(eff2),
status: 'PLANNED',
},
{
id: `boom-ai-3-${Date.now()}`,
name: '3차 방어선 (연안 경량형)',
priority: 'MEDIUM',
type: '연안 경량형',
coords: coords3,
length: computePolylineLength(coords3),
angle: (mainBearing + 90) % 360,
efficiency: Math.round(eff3),
status: 'PLANNED',
},
]
// 최소 효율 필터 경고 (라인은 유지하되 표시용)
for (const line of boomLines) {
if (line.efficiency < settings.minContainmentEfficiency) {
line.name += ' ⚠'
}
}
return boomLines
}
/** Ray casting — 점이 다각형 내부인지 판정 */
export function pointInPolygon(
point: { lat: number; lon: number },
polygon: { lat: number; lon: number }[]
): boolean {
if (polygon.length < 3) return false
let inside = false
const x = point.lon
const y = point.lat
for (let i = 0, j = polygon.length - 1; i < polygon.length; j = i++) {
const xi = polygon[i].lon, yi = polygon[i].lat
const xj = polygon[j].lon, yj = polygon[j].lat
const intersect = ((yi > y) !== (yj > y)) && (x < (xj - xi) * (y - yi) / (yj - yi) + xi)
if (intersect) inside = !inside
}
return inside
}
/** 다각형 면적 (km²) — Shoelace formula, 구면 보정 포함 */
export function polygonAreaKm2(polygon: { lat: number; lon: number }[]): number {
if (polygon.length < 3) return 0
const n = polygon.length
const latCenter = polygon.reduce((s, p) => s + p.lat, 0) / n
const cosLat = Math.cos(latCenter * DEG2RAD)
let area = 0
for (let i = 0; i < n; i++) {
const j = (i + 1) % n
const x1 = polygon[i].lon * cosLat * EARTH_RADIUS * DEG2RAD / 1000
const y1 = polygon[i].lat * EARTH_RADIUS * DEG2RAD / 1000
const x2 = polygon[j].lon * cosLat * EARTH_RADIUS * DEG2RAD / 1000
const y2 = polygon[j].lat * EARTH_RADIUS * DEG2RAD / 1000
area += x1 * y2 - x2 * y1
}
return Math.abs(area) / 2
}
/** 원 면적 (km²) */
export function circleAreaKm2(radiusM: number): number {
return Math.PI * (radiusM / 1000) ** 2
}
/** 차단 시뮬레이션 실행 */
export function runContainmentAnalysis(
trajectory: Array<{ lat: number; lon: number; time: number; particle?: number }>,
boomLines: BoomLine[]
): ContainmentResult {
if (trajectory.length === 0 || boomLines.length === 0) {
return {
totalParticles: 0,
blockedParticles: 0,
passedParticles: 0,
overallEfficiency: 0,
perLineResults: [],
}
}
// 입자별 그룹핑
const particleMap = new Map<number, Array<{ lat: number; lon: number; time: number }>>()
for (const pt of trajectory) {
const pid = pt.particle ?? 0
if (!particleMap.has(pid)) particleMap.set(pid, [])
particleMap.get(pid)!.push(pt)
}
// 입자별 시간 정렬
for (const points of particleMap.values()) {
points.sort((a, b) => a.time - b.time)
}
const perLineBlocked = new Map<string, Set<number>>()
for (const line of boomLines) {
perLineBlocked.set(line.id, new Set())
}
const blockedParticleIds = new Set<number>()
// 각 입자의 이동 경로와 오일펜스 라인 교차 판정
for (const [pid, points] of particleMap) {
for (let i = 0; i < points.length - 1; i++) {
const segA1: BoomLineCoord = { lat: points[i].lat, lon: points[i].lon }
const segA2: BoomLineCoord = { lat: points[i + 1].lat, lon: points[i + 1].lon }
for (const line of boomLines) {
for (let j = 0; j < line.coords.length - 1; j++) {
if (segmentsIntersect(segA1, segA2, line.coords[j], line.coords[j + 1])) {
blockedParticleIds.add(pid)
perLineBlocked.get(line.id)!.add(pid)
}
}
}
}
}
const totalParticles = particleMap.size
const blocked = blockedParticleIds.size
const passed = totalParticles - blocked
return {
totalParticles,
blockedParticles: blocked,
passedParticles: passed,
overallEfficiency: totalParticles > 0 ? Math.round((blocked / totalParticles) * 100) : 0,
perLineResults: boomLines.map(line => {
const lineBlocked = perLineBlocked.get(line.id)!.size
return {
boomLineId: line.id,
boomLineName: line.name,
blocked: lineBlocked,
passed: totalParticles - lineBlocked,
efficiency: totalParticles > 0 ? Math.round((lineBlocked / totalParticles) * 100) : 0,
}
}),
}
}