用GEE和Sentinel-5P数据，5分钟搞定城市空气质量变化趋势分析（以NO2、O3为例）

城市空气质量动态监测基于GEE与Sentinel-5P的高效分析实战清晨打开天气预报APP时那些跳动的空气质量指数背后其实隐藏着卫星每天扫描地球大气层产生的海量数据。作为环境研究者我们完全可以通过Google Earth EngineGEE平台直接调用Sentinel-5P卫星的实时观测数据用代码代替传统手工处理快速捕捉城市上空NO2和O3的微妙变化。本文将带你体验这个高效工作流——无需下载TB级原始数据不用搭建本地计算环境5分钟生成专业级趋势图表的技术魔法。1. 环境准备与数据基础1.1 GEE平台的核心优势与传统遥感数据处理相比GEE的云端计算模式彻底改变了游戏规则。当我们需要分析某城市2020-2023年的空气质量变化时存储零负担所有Sentinel-5P数据已预处理并存储在Google云服务器算力无上限调用ee.ImageCollection时的并行计算速度是本地MATLAB的百倍即时可视化内置的ui.Chart模块可直接生成交互式时间序列图提示注册GEE账号时建议使用机构邮箱教育用户通常能获得更高频的API调用权限1.2 Sentinel-5P数据特性解析这颗专门监测大气成分的卫星其TROPOMI传感器提供了这些关键参数气体类型数据层名称单位典型城市浓度范围NO2NO2_column_number_densitymol/m²0.0001-0.0005O3O3_column_number_densitymol/m²0.12-0.15SO2SO2_column_number_densitymol/m²0.00001-0.0001// 验证数据可用性的基础代码 var testCollection ee.ImageCollection(COPERNICUS/S5P/OFFL/L3_NO2) .filterDate(2023-01-01, 2023-01-31); print(一月数据量, testCollection.size());2. 城市边界数据融合技巧2.1 矢量数据精准匹配分析北京五环内空气质量时粗糙的矩形研究区会导致数据失真。推荐两种获取城市边界的方法GEE内置数据集加载FAO/GAUL/2015/level2行政边界自定义GeoJSON从OpenStreetMap导出后通过ee.FeatureCollection上传// 加载北京市辖区边界示例 var beijing ee.FeatureCollection(FAO/GAUL/2015/level2) .filter(ee.Filter.eq(ADM1_NAME, Beijing)); Map.addLayer(beijing, {color: red}, Beijing Boundary);2.2 时空范围智能筛选针对不同分析目的时间粒度的选择至关重要突发事件响应使用NRTI近实时数据3小时延迟长期趋势分析选择OFFL再处理数据经过质量校正// 最佳实践创建日期序列过滤器 var dateRange ee.DateRange(2020-01-01, 2023-12-31); var yearFilter ee.Filter.calendarRange(2020, 2023, year); var monthFilter ee.Filter.calendarRange(6, 8, month); // 夏季分析3. 空气质量指标深度分析3.1 NO2时空变化解码作为交通排放的指示剂NO2浓度在工作日早高峰会出现典型峰值。这段代码可提取二环内NO2的昼夜差异var no2Daily ee.ImageCollection(COPERNICUS/S5P/OFFL/L3_NO2) .filterBounds(secondRingRoad) .select(NO2_column_number_density); var timeSeries no2Daily.map(function(image) { var date image.date(); var hour date.get(hour); var isWeekday date.getDay().gte(1).and(date.getDay().lte(5)); return image.set(hour, hour).set(isWeekday, isWeekday); }); // 生成小时均值统计图表 var chart ui.Chart.image.seriesByRegion({ imageCollection: timeSeries, regions: secondRingRoad, reducer: ee.Reducer.mean(), band: NO2_column_number_density, scale: 1000, xProperty: hour, seriesProperty: isWeekday }).setOptions({ title: 二环内NO2浓度小时变化, hAxis: {title: 小时}, vAxis: {title: NO2浓度 (mol/m²)}, lineWidth: 2, series: { 0: {color: red, lineDash: [4, 4]}, // 周末 1: {color: blue} // 工作日 } }); print(chart);3.2 O3与气象因子的耦合分析不同于NO2臭氧浓度往往与温度正相关。这个工作流整合了ERA5气象数据var ozone ee.ImageCollection(COPERNICUS/S5P/OFFL/L3_O3); var temperature ee.ImageCollection(ECMWF/ERA5/DAILY); var joinedCollection ee.Join.saveAll(matches).apply({ primary: ozone.filterDate(2022-06-01, 2022-08-31), secondary: temperature.filterDate(2022-06-01, 2022-08-31), condition: ee.Filter.equals({ leftField: system:time_start, rightField: system:time_start }) }); var correlationChart ui.Chart.image.doySeriesByYear({ imageCollection: joinedCollection.map(function(image) { var matches ee.ImageCollection.fromImages(image.get(matches)); var temp matches.first().select(mean_2m_air_temperature); return image.addBands(temp); }), bandName: [O3_column_number_density, mean_2m_air_temperature], region: studyArea, regionReducer: ee.Reducer.mean(), scale: 1000 }).setOptions({ title: 夏季O3浓度与气温相关性, series: { 0: {targetAxisIndex: 0, color: purple}, 1: {targetAxisIndex: 1, color: orange} }, vAxes: { 0: {title: O3浓度 (mol/m²)}, 1: {title: 温度 (°C)} } });4. 工业事件的可视化追踪4.1 疫情期间排放变化2020年春季的封锁措施在卫星数据中留下清晰信号。这段代码对比了疫情前后同期数据var preCovid ee.ImageCollection(COPERNICUS/S5P/OFFL/L3_NO2) .filter(ee.Filter.date(2019-02-15, 2019-04-15)); var duringCovid ee.ImageCollection(COPERNICUS/S5P/OFFL/L3_NO2) .filter(ee.Filter.date(2020-02-15, 2020-04-15)); var diff duringCovid.mean() .subtract(preCovid.mean()) .clip(cityBoundary); var diffViz { min: -0.0002, max: 0.0002, palette: [red, white, blue] }; Map.addLayer(diff, diffViz, NO2浓度变化);4.2 节假日效应分析春节期间的工业活动停滞会显著改变空气成分构成。使用reduceRegion进行多年度对比var springFestival ee.Filter.or( ee.Filter.date(2019-02-04, 2019-02-10), ee.Filter.date(2020-01-24, 2020-01-30), ee.Filter.date(2021-02-11, 2021-02-17) ); var festivalStats ee.ImageCollection(COPERNICUS/S5P/OFFL/L3_SO2) .filter(springFestival) .map(function(image) { var stats image.reduceRegion({ reducer: ee.Reducer.mean(), geometry: industrialZone, scale: 1000 }); return image.set(stats); }); print(ui.Chart.feature.byFeature(festivalStats, system:time_start, [SO2_column_number_density]));5. 自动化报告生成策略5.1 动态仪表板构建将上述分析整合为交互式面板使用ui.Panel实现一键刷新var controlPanel ui.Panel([ ui.Label(时间范围选择), ui.DatePicker({id: startDate, value: 2020-01-01}), ui.DatePicker({id: endDate, value: 2020-12-31}), ui.Button(更新图表, updateCharts) ]); function updateCharts() { var start controlPanel.widgets().get(1).getValue(); var end controlPanel.widgets().get(2).getValue(); // 动态更新所有图表 }5.2 数据导出最佳实践当需要本地进一步分析时这些方法能保持数据完整性CSV导出使用Export.table.toDrive输出统计结果GeoTIF导出对空间分布数据用Export.image.toDrive自动命名技巧在导出任务中添加时间戳变量// 典型导出配置 Export.table.toDrive({ collection: timeSeriesStats, description: NO2_Monthly_Stats, fileFormat: CSV, selectors: [date, mean, max, min] });在最近为某环保组织做的咨询项目中发现将GEE分析结果与地面监测站数据结合时需要注意卫星过境时间当地时间上午10:30左右与地面监测时段的匹配问题。建议在制作对比图表时统一使用世界协调时(UTC)以避免时区混淆。