Urban-Air-Quality-Digital-Twin

一座城市的呼吸剖面

大多數城市只靠十幾座官方空氣品質監測站來代表數百萬人的呼吸環境。這些監測站按照法規設置在特定高度和位置——通常是遠離直接排放源的「背景站」——但它們完全無法反映街道峽谷中行人實際吸入的空氣品質。同一條街，上午八點車流密集的一側 PM2.5 可以比背風側高出三倍；中午陽光催化下，地面臭氧在通風不良的巷弄中累積到超過標準值。

城市空污數位孿生的核心任務，是建立一個米級解析度、分鐘級更新頻率的虛擬空氣實驗室。輸入端來自固定監測站、移動感測器（公車、計程車、無人機搭載）、衛星氣溶膠光學厚度反演、交通車流即時數據、氣象場（風場、混合層高度、溫度逆層），以及工業排放申報資料。這些多源異質資料經過同化演算法，餵入計算流體力學（CFD）微尺度模型後，產出每一條街道、每一個小時的六種污染物濃度場。

排放溯源的兩難

知道空氣差是一回事，知道誰排放了什麼才是政策干預的關鍵。但排放溯源在技術和法律上都極其敏感——一組 PM2.5 濃度異常可能來自三公里外的工廠、五條街外的工地揚塵、或是鄰近夜市集中油炸產生的油煙。如果模型將污染歸因到特定排放源並觸發罰則，必須承受非常高的證據標準。

實務上採用伴隨法（adjoint method）進行反向追蹤：從監測到的高濃度網格點出發，沿著風場和擴散參數逆向積分，計算每個潛在排放源對該濃度的貢獻敏感度。但伴隨模擬的計算成本極高，不適合即時溯源。折衷方案是預先計算全年度不同氣象條件下的「足跡資料庫」，即時監測到異常時只需查表插值，將溯源延遲從數小時壓縮到數分鐘。目前一線城市的系統可以在檢測到 PM2.5 異常峰值後 8 分鐘內，產出前三大嫌疑排放源的排序清單與各自的貢獻機率。

從數據到行動的轉譯

空污數位孿生最容易失敗的環節不是模型不準，而是產出的資訊沒有人使用。當系統以每 15 分鐘一次的頻率更新數百萬個網格點的污染物濃度時，如果沒有針對不同使用者的需求進行轉譯，這些資料就只是螢幕上漂亮的熱力圖。

對一般市民，轉譯後的資訊應該是——「今天下午三點接小孩放學，建議走路線 B 而不是路線 A，暴露量可減少 40%」。對學校，應該是——「明日清晨逆溫層穩定，懸浮微粒不易擴散，建議取消上午的戶外體育課」。對環保單位，則是——「過去一週東工業區夜間排放超標三次，鄰近居民區暴露增加 28%，建議啟動稽查」。同一個孿生模型，不同使用者看到的是完全不同顆粒度的決策介面。

這層轉譯的困難在於：健康風險溝通不能製造恐慌，但也不能淡化風險。實務上越來越傾向使用「空氣品質健康指數」（AQHI）而非傳統的 AQI——AQHI 將污染物濃度直接轉換為超額死亡或急診就醫風險的增加百分比，這對一般民眾來說比「127 μg/m³」更有意義。

import numpy as np
from scipy.sparse.linalg import lsqr

class AdjointSourceTracer:
    # Adjoint-based pollution source attribution.
    # Given a concentration spike at a receptor, estimates source contributions.

    def __init__(self, footprint_library, wind_fields):
        self.footprints = footprint_library  # Precomputed adjoint footprints
        self.wind = wind_fields

    def trace(self, receptor_coords, concentration, timestamp):
        # Match wind condition to select appropriate footprint.
        wind_key = self._discretize_wind(self.wind[timestamp])
        footprint = self.footprints[(receptor_coords, wind_key)]

        # Solve: concentration = footprint @ emissions + noise
        sources = list(footprint.keys())
        A = np.array([footprint[s] for s in sources]).T
        emissions, _, _, _ = lsqr(A, concentration)
        ranking = sorted(zip(sources, emissions), key=lambda x: -x[1])
        return [(src, round(contrib, 2)) for src, contrib in ranking[:3]]

    def _discretize_wind(self, w):
        return (round(w[0] / 45) * 45, round(w[1], 1))