|
| 1 | +// |
| 2 | +// EPAAirQuality.swift |
| 3 | +// Meshtastic |
| 4 | +// |
| 5 | +// US EPA air quality math for PM2.5 (issue #2040 / design#54): |
| 6 | +// • NowCast — a 12-hour weighted rolling average that favors recent hours. |
| 7 | +// • AQI — the 0–500 index derived from a PM2.5 concentration via EPA Equation 1. |
| 8 | +// |
| 9 | +// References: |
| 10 | +// • EPA/AirNow "Technical Assistance Document for the Reporting of Daily Air Quality" |
| 11 | +// (Equation 1, truncation + rounding rules, NowCast data-sufficiency rule). |
| 12 | +// • EPA AQS breakpoint code table (parameter 88101, PM2.5 24-hour), current |
| 13 | +// post-2024 concentration breakpoints. |
| 14 | +// |
| 15 | + |
| 16 | +import Foundation |
| 17 | + |
| 18 | +enum EPAAirQuality { |
| 19 | + |
| 20 | + /// One row of the EPA AQI breakpoint table: a concentration range (µg/m³) mapped to an AQI range. |
| 21 | + struct PM25Breakpoint { |
| 22 | + let concentrationLow: Double |
| 23 | + let concentrationHigh: Double |
| 24 | + let aqiLow: Int |
| 25 | + let aqiHigh: Int |
| 26 | + } |
| 27 | + |
| 28 | + /// Current (2024) EPA PM2.5 24-hour AQI breakpoints (µg/m³ → AQI index). |
| 29 | + /// Contiguous at 0.1 µg/m³ resolution; concentrations are truncated to one decimal |
| 30 | + /// place before lookup so the 0.1 gaps between rows (e.g. 9.0 / 9.1) never fall through. |
| 31 | + /// The app's `Aqi` scale tops out at 500, so concentrations above the final breakpoint |
| 32 | + /// clamp to 500 rather than extending into the 501–999 "beyond the AQI" range. |
| 33 | + static let pm25Breakpoints: [PM25Breakpoint] = [ |
| 34 | + PM25Breakpoint(concentrationLow: 0.0, concentrationHigh: 9.0, aqiLow: 0, aqiHigh: 50), |
| 35 | + PM25Breakpoint(concentrationLow: 9.1, concentrationHigh: 35.4, aqiLow: 51, aqiHigh: 100), |
| 36 | + PM25Breakpoint(concentrationLow: 35.5, concentrationHigh: 55.4, aqiLow: 101, aqiHigh: 150), |
| 37 | + PM25Breakpoint(concentrationLow: 55.5, concentrationHigh: 125.4, aqiLow: 151, aqiHigh: 200), |
| 38 | + PM25Breakpoint(concentrationLow: 125.5, concentrationHigh: 225.4, aqiLow: 201, aqiHigh: 300), |
| 39 | + PM25Breakpoint(concentrationLow: 225.5, concentrationHigh: 325.4, aqiLow: 301, aqiHigh: 500) |
| 40 | + ] |
| 41 | + |
| 42 | + /// Minimum weight factor for the PM NowCast; approximates a 3-hour average when |
| 43 | + /// pollutant levels are changing rapidly. |
| 44 | + static let nowCastMinimumWeight = 0.5 |
| 45 | + |
| 46 | + /// AQI (0–500) from a PM2.5 concentration (µg/m³) using EPA Equation 1: |
| 47 | + /// I = (Ihi − Ilo) / (BPhi − BPlo) · (C − BPlo) + Ilo |
| 48 | + /// The concentration is truncated to one decimal place; the result is rounded to the |
| 49 | + /// nearest integer. Returns `nil` for negative input. |
| 50 | + static func aqi(fromPM25 concentration: Double) -> Int? { |
| 51 | + guard concentration >= 0 else { return nil } |
| 52 | + // Truncate to one decimal place (EPA rule for PM2.5). Direct sensor readings are integer |
| 53 | + // µg/m³ (Stage 1 stores pm25Standard as UInt32) so they truncate exactly; NowCast averages |
| 54 | + // are fractional, where an exact-boundary hit is measure-zero and at most ±1 AQI. `.towardZero` |
| 55 | + // == floor here because the value is already guarded non-negative. |
| 56 | + let c = (concentration * 10).rounded(.towardZero) / 10.0 |
| 57 | + for bp in pm25Breakpoints where c <= bp.concentrationHigh { |
| 58 | + let slope = Double(bp.aqiHigh - bp.aqiLow) / (bp.concentrationHigh - bp.concentrationLow) |
| 59 | + let index = slope * (c - bp.concentrationLow) + Double(bp.aqiLow) |
| 60 | + return min(500, max(0, Int(index.rounded()))) |
| 61 | + } |
| 62 | + // Above the top breakpoint: clamp to the ceiling of the app's AQI scale. |
| 63 | + return 500 |
| 64 | + } |
| 65 | + |
| 66 | + /// EPA NowCast for PM2.5 from up to 12 hourly averages. |
| 67 | + /// `hourly[0]` is the most recent clock hour, `hourly[11]` is 11 hours earlier; `nil` |
| 68 | + /// marks a missing hour. Returns `nil` when there isn't enough recent data — EPA requires |
| 69 | + /// at least two of the three most recent hours to be valid. |
| 70 | + static func nowCastPM25(hourly: [Double?]) -> Double? { |
| 71 | + let hours = Array(hourly.prefix(12)) |
| 72 | + // Data sufficiency: 2 of the 3 most recent hours must be present. |
| 73 | + guard hours.prefix(3).compactMap({ $0 }).count >= 2 else { return nil } |
| 74 | + |
| 75 | + let valid = hours.compactMap { $0 } |
| 76 | + guard let cMax = valid.max(), let cMin = valid.min() else { return nil } |
| 77 | + // If the max is 0 every reading is 0 → NowCast is 0 (avoids divide-by-zero on the weight). |
| 78 | + guard cMax > 0 else { return 0 } |
| 79 | + |
| 80 | + // Weight factor = 1 − (scaled rate of change), floored at the minimum weight. |
| 81 | + let weightFactor = max(cMin / cMax, nowCastMinimumWeight) |
| 82 | + |
| 83 | + var numerator = 0.0 |
| 84 | + var denominator = 0.0 |
| 85 | + for (hoursAgo, value) in hours.enumerated() { |
| 86 | + guard let concentration = value else { continue } |
| 87 | + let weight = pow(weightFactor, Double(hoursAgo)) |
| 88 | + numerator += weight * concentration |
| 89 | + denominator += weight |
| 90 | + } |
| 91 | + guard denominator > 0 else { return nil } |
| 92 | + return numerator / denominator |
| 93 | + } |
| 94 | + |
| 95 | + /// Buckets timestamped PM2.5 readings into the most recent 12 clock hours, computes the |
| 96 | + /// NowCast, and maps it to an AQI (0–500). Returns `nil` when NowCast can't be computed |
| 97 | + /// (insufficient recent history) — callers should fall back to showing the raw reading. |
| 98 | + static func nowCastAQI( |
| 99 | + from readings: [(date: Date, pm25: Double)], |
| 100 | + now: Date = Date(), |
| 101 | + calendar: Calendar = .current |
| 102 | + ) -> Int? { |
| 103 | + guard !readings.isEmpty else { return nil } |
| 104 | + let currentHourStart = calendar.dateInterval(of: .hour, for: now)?.start ?? now |
| 105 | + |
| 106 | + var buckets: [Int: (sum: Double, count: Int)] = [:] |
| 107 | + for reading in readings { |
| 108 | + guard let hourStart = calendar.dateInterval(of: .hour, for: reading.date)?.start else { continue } |
| 109 | + guard let hoursAgo = calendar.dateComponents([.hour], from: hourStart, to: currentHourStart).hour else { continue } |
| 110 | + guard hoursAgo >= 0, hoursAgo < 12 else { continue } |
| 111 | + var bucket = buckets[hoursAgo] ?? (sum: 0, count: 0) |
| 112 | + bucket.sum += reading.pm25 |
| 113 | + bucket.count += 1 |
| 114 | + buckets[hoursAgo] = bucket |
| 115 | + } |
| 116 | + |
| 117 | + let hourly: [Double?] = (0..<12).map { hoursAgo in |
| 118 | + guard let bucket = buckets[hoursAgo], bucket.count > 0 else { return nil } |
| 119 | + return bucket.sum / Double(bucket.count) |
| 120 | + } |
| 121 | + |
| 122 | + guard let nowCast = nowCastPM25(hourly: hourly) else { return nil } |
| 123 | + return aqi(fromPM25: nowCast) |
| 124 | + } |
| 125 | +} |
0 commit comments