/**
 * Environmental model: air temperature, humidity, wind, cloud, dust and solar
 * irradiance for a region at an instant.
 *
 * Physically motivated rather than random: irradiance comes from actual solar geometry
 * (declination + hour angle at the region's latitude), temperature from a seasonal +
 * diurnal harmonic, and everything carries slow fractal drift so successive samples look
 * like weather rather than noise.
 */

import { clamp, fbm } from './rng'
import type { RegionSeed } from '../data/regions'

/** Saudi Arabia observes UTC+3 year-round, with no daylight saving. */
export const KSA_UTC_OFFSET_HOURS = 3
const MS_PER_HOUR = 3_600_000
const MS_PER_DAY = 86_400_000

export interface EnvironmentSample {
  ts: number
  tempC: number
  humidityPct: number
  windSpeedMs: number
  windDirDeg: number
  irradianceWm2: number
  cloudCoverPct: number
  dustIndex: number
  visibilityKm: number
}

/** Local (Asia/Riyadh) clock hour as a float, e.g. 14.5 for 14:30. */
export function localHour(ts: number): number {
  return (((ts / MS_PER_HOUR + KSA_UTC_OFFSET_HOURS) % 24) + 24) % 24
}

/** Local day-of-year in [1, 366]. */
export function dayOfYear(ts: number): number {
  const local = new Date(ts + KSA_UTC_OFFSET_HOURS * MS_PER_HOUR)
  const start = Date.UTC(local.getUTCFullYear(), 0, 0)
  return Math.floor((local.getTime() - start) / MS_PER_DAY)
}

/** Continuous days since epoch — the time axis for slow weather drift. */
function dayIndex(ts: number): number {
  return ts / MS_PER_DAY
}

const DEG = Math.PI / 180

/**
 * Clear-sky global horizontal irradiance (W/m²) from solar geometry.
 * Uses the Cooper declination approximation and a simple air-mass attenuation.
 */
export function clearSkyIrradiance(region: RegionSeed, ts: number): number {
  const doy = dayOfYear(ts)
  const declination = 23.45 * Math.sin(DEG * ((360 * (284 + doy)) / 365))
  // Correct clock time to local solar time using the longitude offset from the
  // UTC+3 standard meridian (45°E).
  const solarTime = localHour(ts) + (region.lon - 45) / 15
  const hourAngle = 15 * (solarTime - 12)
  const sinElevation =
    Math.sin(DEG * region.lat) * Math.sin(DEG * declination) +
    Math.cos(DEG * region.lat) * Math.cos(DEG * declination) * Math.cos(DEG * hourAngle)
  if (sinElevation <= 0) return 0
  // 1090 W/m² peak with a mild air-mass exponent gives realistic desert clear-sky GHI.
  return 1090 * Math.pow(sinElevation, 1.15)
}

/** Solar elevation in degrees (negative at night). Used for dusk/dawn UI cues. */
export function solarElevation(region: RegionSeed, ts: number): number {
  const doy = dayOfYear(ts)
  const declination = 23.45 * Math.sin(DEG * ((360 * (284 + doy)) / 365))
  const solarTime = localHour(ts) + (region.lon - 45) / 15
  const hourAngle = 15 * (solarTime - 12)
  const sinEl =
    Math.sin(DEG * region.lat) * Math.sin(DEG * declination) +
    Math.cos(DEG * region.lat) * Math.cos(DEG * declination) * Math.cos(DEG * hourAngle)
  return Math.asin(clamp(sinEl, -1, 1)) / DEG
}

function regionSeed(region: RegionSeed): number {
  let h = 0
  for (let i = 0; i < region.code.length; i += 1) h = (h * 31 + region.code.charCodeAt(i)) | 0
  return h
}

export interface EnvironmentOverrides {
  tempDeltaC?: number
  cloudDeltaPct?: number
  dustDelta?: number
  windDeltaPct?: number
}

/** Full environmental sample for a region at `ts`, with optional scenario overrides. */
export function sampleEnvironment(
  region: RegionSeed,
  ts: number,
  overrides: EnvironmentOverrides = {},
): EnvironmentSample {
  const seed = regionSeed(region)
  const d = dayIndex(ts)
  const hour = localHour(ts)
  const doy = dayOfYear(ts)

  // Seasonal harmonic peaks around day 200 (mid-July).
  const seasonal = region.seasonalAmpC * Math.cos((2 * Math.PI * (doy - 200)) / 365)
  // Diurnal harmonic peaks at 15:00 local.
  const diurnal = region.dailyAmpC * Math.cos((2 * Math.PI * (hour - 15)) / 24)
  // Weather fronts drift over ~2.5 days.
  const front = fbm(seed + 11, d / 2.5, 3) * 2.6

  const tempC = region.annualMeanC + seasonal + diurnal + front + (overrides.tempDeltaC ?? 0)

  const cloudCoverPct = clamp(
    region.baseCloudPct + fbm(seed + 23, d / 1.6, 3) * 26 + (overrides.cloudDeltaPct ?? 0),
    0,
    100,
  )
  const dustIndex = clamp(
    region.baseDustIndex + fbm(seed + 37, d / 3.1, 2) * 0.22 + (overrides.dustDelta ?? 0),
    0,
    1,
  )

  // Humidity falls as air warms; coastal regions stay high.
  const humidityPct = clamp(
    region.baseHumidityPct - diurnal * 1.6 + fbm(seed + 41, d / 2, 2) * 9,
    3,
    98,
  )

  // Wind: daily sea-breeze component plus fractal gusting.
  const breeze = 1 + 0.35 * Math.cos((2 * Math.PI * (hour - 16)) / 24)
  const windSpeedMs = clamp(
    (region.baseWindMs * breeze + fbm(seed + 53, d * 3.4, 3) * 2.4) *
      (1 + (overrides.windDeltaPct ?? 0) / 100),
    0.2,
    28,
  )
  const windDirDeg = (((fbm(seed + 67, d / 1.2, 2) + 1) / 2) * 360 + hour * 3) % 360

  const irradianceWm2 =
    clearSkyIrradiance(region, ts) *
    (1 - 0.78 * (cloudCoverPct / 100)) *
    (1 - 0.32 * dustIndex)

  const visibilityKm = clamp(28 - dustIndex * 22 - (humidityPct / 100) * 4, 1.5, 30)

  return {
    ts,
    tempC,
    humidityPct,
    windSpeedMs,
    windDirDeg,
    irradianceWm2: Math.max(0, irradianceWm2),
    cloudCoverPct,
    dustIndex,
    visibilityKm,
  }
}

/**
 * Photovoltaic output as a fraction of nameplate capacity.
 * Includes the temperature derate real PV modules suffer in desert heat
 * (≈ -0.38 %/°C above the 25 °C reference), which is exactly why NABDH treats a
 * heatwave as *both* a demand and a supply problem.
 */
export function solarCapacityFactor(env: EnvironmentSample): number {
  if (env.irradianceWm2 <= 1) return 0
  const cellTempC = env.tempC + (env.irradianceWm2 / 800) * 25
  const tempDerate = clamp(1 - 0.0038 * (cellTempC - 25), 0.6, 1)
  const soiling = clamp(1 - 0.12 * env.dustIndex, 0.85, 1)
  return clamp((env.irradianceWm2 / 1000) * tempDerate * soiling, 0, 1)
}

/**
 * Wind turbine output as a fraction of nameplate capacity, using a standard
 * cut-in / rated / cut-out power curve (3 / 12 / 25 m/s).
 */
export function windCapacityFactor(env: EnvironmentSample): number {
  const v = env.windSpeedMs
  if (v < 3 || v >= 25) return 0
  if (v >= 12) return 0.98
  return clamp(Math.pow((v - 3) / 9, 3) * 0.98, 0, 0.98)
}
