Welcome to visit Anhui Agricultural Science Bulletin!

Research on estimation of rice leaf nitrogen content based on multi-source satellite data

  • JIA Tao 1 ,
  • TANG Lifang 1 ,
  • WANG Jianfei 2 ,
  • LU Lijiang 2
Expand
  • 1. Agriculture and Rural Bureau of Fengyang, Fengyang 233100, China
  • 2. Anhui Science and Technology University, Fengyang 233100, China

Received date: 2024-12-27

  Online published: 2025-08-28

Abstract

To investigate the accuracy of estimating rice leaf nitrogen concentration using multi-source satellite data, a rice LNC estimation experiment was conducted from 2020 to 2022 in Xiaogang Village, Fengyang County, Anhui Province. The HH2 spectroradiometer was used to obtain in situ hyperspectral data of rice leaves at different growth stages (tillering, booting, heading, flowering, and filling stages), while synchronous measurements of rice nitrogen content were collected. The band response functions of two satellite sensors, QuickBird and SPOT-6, were selected, and the reflectance data of these two sensors was simulated using ENVI/IDL resampling technology.Based on the simulated multispectral reflectance data of rice leaves, 12 spectral indices were constructed by combining any two bands. Spectral indices highly correlated with nitrogen content at each growth stage were screened through correlation fitting, and the established rice LNC estimation models were validated using the coefficient of determination (R2) and root mean square error (RMSE). The results showed that, across the entire growth period, rice LNC exhibited a gradually decreasing trend (except during the tillering stage) and was negatively correlated with spectral reflectance.For the QuickBird sensor, during the tillering and booting stages, the normalized difference spectral index (NDSI-BG) and ratio spectral index (RSI-BR) showed the highest correlations with LNC, with R2 values of 0.351 8, 0.617 0 and 0.352 9, 0.628 1, respectively. During the heading, flowering, and filling stages, NDSI-RN and RSI-RN achieved the highest R2 values, reaching 0.714 0, 0.553 2, 0.637 2 and 0.707 5, 0.542 9, 0.639 5, respectively. Similarly, for the SPOT-6 sensor, NDSI-BG and RSI-BR exhibited the highest R2 values during the tillering and booting stages (0.351 8, 0.612 9 and 0.355 4, 0.622 2, respectively), while NDSI-RN and RSI-RN performed best during the heading, flowering, and filling stages (0.702 9, 0.543 6, 0.584 6 and 0.697 8, 0.537 9, 0.584 8, respectively).Model validation results indicated that the QuickBird-simulated data, R2 = 0.56 and RMSE = 0.42, whereas for the SPOT-6 model, R2=0.54 and RMSE = 0.28. This suggests that the LNC estimation model based on QuickBird-simulated data achieved slightly higher accuracy than that of SPOT-6, though with a marginally larger RMSE. This method can be applied for relatively rough estimation of rice nitrogen content.

Cite this article

JIA Tao , TANG Lifang , WANG Jianfei , LU Lijiang . Research on estimation of rice leaf nitrogen content based on multi-source satellite data[J]. Anhui Agricultural Science Bulletin, 2025 , 31(16) : 1 -8 . DOI: 10.16377/j.cnki.issn1007-7731.2025.16.001

氮是水稻生长必需的营养元素之一,对其生长过程以及品质和产量具有重要的影响[1]。基于常规方法获取作物氮素含量数据的人力和时间成本较高,不利于大田作物观测。近年来,基于近地平台的作物氮素遥感监测得到越来越多地关注,但受限于数据尺度和成本,该方法的推广具有一定难度。为解决上述问题,相关研究人员基于近地面高光谱数据和卫星响应函数构建了作物氮浓度估测模型[2-4],这种方法具有较高的精度和效率,可在大范围作物监测研究中应用[5-7]。Huang等[4]研究指出,卫星遥感技术在大规模的作物生长监测和精准管理中具有广阔的应用前景。李粉玲等[1]以大田冬小麦为研究对象,基于GF-1号卫星的光谱响应函数结合实地采样与测量,构建了基于常见光谱指数的多个氮素模型,可用于冬小麦叶片的氮含量估测。Wu等[8]对比SPAD-502叶绿素仪测量数据与Quick Bird卫星遥感影像在评估不同氮素处理下马铃薯冠层氮素状况的表现,结果发现,Quick Bird卫星数据能够更为灵敏地反映土豆氮素变化。贾良良等[9]采用IKONOS卫星获取冬小麦拔节期图像,田间试验分析多个植被指数与氮素含量之间的相关性,并基于此构建氮素诊断模型,结果表明,该模型能准确反映冬小麦拔节期的氮素营养状况。王备战等[10]和武婕等[11]利用SPOT-5遥感影像数据,建立了冬小麦和玉米氮积累量的估测模型,并生成了研究区的遥感监测专题图。Eitel等[12]对Rapid Eye卫星的宽带数据中包含红边区域的多个植被指数进行了评估,结果表明,采用宽带Rapid Eye卫星数据的MCARI/MTVI2指数在预测小麦氮素营养状态上展现出一定的优势。陈拉等[13]利用不同卫星传感器的模拟数据,构建了多种光谱指数,并深入分析了这些植被指数在预测水稻叶面积指数时的估测精度和敏感性,结果表明,与红波段植被指数相比,红边比光谱指数(RVI)和绿波光谱指数(GRVI)与估测模型展现出更强的相关性。杨海波等[14]基于地面实测的马铃薯冠层高光谱数据,模拟了WorldView-2和VENμS卫星的不同波段数据,并基于这些波段构建了估测模型,结果表明,红边波段的融合光谱指数能够提高马铃薯氮素含量的监测效果。基于此,本研究以卫星响应函数与地面高光谱数据融合的方式,模拟反射率并建立氮含量估测模型,以实现大范围作物监测,为大田范围的水稻氮素含量监测提供参考。

1 材料与方法

1.1 研究区基本情况

试验于2020—2022年在安徽滁州凤阳小岗村进行,该区属于北亚热带向温带气候的过渡地带,气候温和宜人,四季变化明显,阳光照射充分,无霜期较长。水稻地块被划分为36个小区,每个小区长×宽为2 m×8 m,氮、磷肥一次性底施,结合大田管理促早期生长,确保作物健康高产。重点关注水稻生长的分蘖期、孕穗期、抽穗期、开花期和灌浆期。各个生育期内,均进行了水稻叶片光谱测定和植株采样。

1.2 数据来源

1.2.1 光谱数据

在11:00—12:00,采用HH2光谱仪(ASD公司),对水稻叶片进行高光谱反射率的数据采集,其参数配置如表1所示。数据采集前需使用白板校正,每个小区按“Z”字形测量3次。传感器对准水稻冠层,高度保持0.5 m,确保水稻覆盖视场角,避免阴影。将获取的光谱数据,通过ViewSpec Pro处理后保存为表格数据。在分蘖期、孕穗期、抽穗期、开花期和灌浆期均进行上述采样工作。
表1 HH2 光谱仪的相关参数
参数 性能参数
波长范围 325~1 075 nm
波长精度 ±1 nm
光谱分辨率 <3 nm @700 nm
等效噪声辐射(NEdL) 标准HandHeld 2:5×10-9 W/cm2/nm/sr @ 700 nm HandHeld 2 Pro:5×10-10 W/cm2/nm/sr @ 700 nm
采样时间间隔 最小8.5 ms(可选择设置)
视场角 25°(可选购其他可用前置光学系统)
温度范围

工作温度:0~40 ℃(32~104 ℉)

储存温度:0~45 ℃(32~113 ℉)

工作和保存湿度:90%非冷凝

电脑软件 RS3TM控制软件ViewSpec ProTM后续处理软件,HH2 Sync
采样间隔 1 nm
光谱文件大小 约30 KB
储存空间 >2 000个光谱文件

1.2.2 氮含量数据

首先获取不同时段内水稻的高光谱影像数据,在同一时间对每个小区3个均匀且有代表性的0.5 m双行区域进行采样,剪穗后带回实验室,对其样本进行预处理,包括清洗、烘干和研磨等步骤。使用凯氏定氮法进行全氮含量(以质量分数计)的测定。

1.3 数据处理及分析

1.3.1 高光谱重采样

Quick Bird卫星的传感器包括蓝光(B,450~520 nm)、绿光(G,520~600 nm)、红光(R,630~690 nm)和近红外(NIR,760~900nm)波段,空间分辨率为0.61 m。SPOT-6卫星的传感器包括蓝光(455~525 nm)、绿光(530~590 nm)、红光(625~695 nm)和近红外(760~890 nm)波段,空间分辨率为6 m。使用Quick Bird和SPOT-6官网获取的波段响应函数,利用ENVI 5.3软件制成SLI格式的光谱响应函数库。使用ENVI 5.3软件,对地面实测的高光谱数据进行重采样,模拟Quick Bird和SPOT-6卫星的波段,计算如式(1)
P b a n d = b a n d m i n b a n d m a x s b a n d ( λ ) p ( λ ) d λ
式(1)中,Pband 是模拟卫星波段的反射率;sband (λ)是卫星的光谱响应函数;bandmaxbandmin 分别代表波段的上限和下限波长;p(λ)是实测的冠层高光谱数据。

1.3.2 构建光谱指数

基于Quick Bird和SPOT-6卫星模拟数据,以任意波段两两组合的方式构建了由不同波段组合的归一化植被指数(NDSI)和比值植被指数(RSI),计算如式(2)~(3)。
N D S I ( R λ 1 ; R λ 2 ) =   ( R λ 1 R λ 2 ) /   ( R λ 1 + R λ 2 )
R S I ( R λ 1 ; R λ 2 )   =   R λ 1 / R λ 2
式中, R λ 1 R λ 2分别代表波长λ1λ2处的光谱反射率。

1.3.3 模型的构建与验证

通过计算NDSI和RSI与氮素含量的拟合决定系数R2,提取氮素含量估测的中心敏感波段,并匹配Quick Bird和SPOT-6卫星的敏感波段信息,以确定与氮素含量相关性最高的光谱指数,建立偏最小二乘回归(PLSR)模型。基于前期构建的具体模型,本研究主要对其进行验证。PLSR模型是通过决定系数(R2)和均方根误差(RMSE)进行精度评估。随机分割2/3为训练集,1/3为测试集,重复此过程10次以评估模型性能与稳定性。模型评价指标的计算如式(4)~(5)。
R 2 = i = 1 N Y i a c t - Y a c t ¯ Y i p r e - Y p r e ¯ i = 1 N Y i a c t - Y a c t ¯ 2 i = 1 N Y i p r e - Y p r e ¯ 2
R M S E = i = 1 N Y i a c t - Y i p r e 2 N
式中, Y i a c t   Y i p r e分别表示第某个样本的实测值和其预测值; Y a c t ¯ Y p r e ¯分别表示实测值的平均值与模型预测值的平均值;N表示样本的个数。

2 结果与分析

2.1 水稻氮素含量及高光谱特征的变化

2.1.1 水稻氮素含量变化特征

图1可知,分蘖期氮含量在1.34%~4.18%,平均含量为3.15%;孕穗期氮含量在2.12%~4.31%,平均含量为3.29%;抽穗期氮含量在1.72%~4.12%,平均含量为3.10%;开花期氮含量在1.35%~4.23%,平均含量为2.72%;灌浆期氮含量在1.36%~3.93%,平均含量为2.70%。除分蘖期外,氮浓度平均值从孕穗期到灌浆期逐步降低,表明生育期对水稻氮素含量有一定的影响。
图1 不同生育期水稻氮素含量

2.1.2 高光谱变化

通过高光谱数据处理,得到了分蘖期、孕穗期、抽穗期、开花期、灌浆期的反射率。由图2可知,曲线在480和670 nm处形成吸收谷,550 nm处有小峰值,680~750 nm波段反射率急剧上升,700~900 nm波段反射率较高且稳定。不同生育期的光谱反射率存在差异,分蘖期最低,孕穗期和抽穗期数值相近。整个生育期(除分蘖期外)水稻叶片氮素含量与光谱反射率呈负相关。生育期是影响反射率的重要因素,因此,在建模过程中充分考虑不同生育期的卫星模拟数据来估测氮含量。
图2 不同生育期的水稻冠层高光谱曲线

2.2 高光谱数据重采样

基于水稻叶片高光谱曲线数据和响应函数,重采样生成模拟数据。由图3可知,Quick Bird和SPOT-6传感器都包含了蓝光、绿光、红光和近红外4个波段。
图3 高光谱数据重采样模拟数据

(A)、(B)分别为Quick Bird和SPOT-6卫星模拟数据。

2.3 光谱指数与氮素含量的拟合分析

2.3.1 Quick Bird模拟数据

基于Quick Bird模拟数据,不同生育期NDSI与氮素含量的拟合结果见图4表2。分蘖期中,NDSI-BG的R2最高,为0.351 8;孕穗期的R2以NDSI-BG最高,为0.617 0;在抽穗期、开花期、灌浆期,NDSI-RN的R2均最高,分别为0.714 0、0.553 2、0.637 2。整个生育期中,NDSI-BG和NDSI-RN对水稻氮含量的线性拟合效果较好。
图4 Quick Bird模拟数据的NDSI与氮素含量的拟合图

(A)~(E)分别为分蘖期、孕穗期、抽穗期、开花期、灌浆期。

表2 不同生育期NDSI与氮素含量的线性拟合R2
生育期 NDSI-BG NDSI-BR NDSI-BN NDSI-GR NDSI-GN NDSI-RN
分蘖期 0.351 8 0.201 6 0.285 5 0.249 9 0.295 3 0.340 0
孕穗期 0.617 0 0.245 7 0.542 7 0.525 8 0.577 4 0.591 9
抽穗期 0.591 5 0.014 7 0.623 7 0.395 9 0.615 4 0.714 0
开花期 0.536 0 0.081 1 0.417 2 0.354 6 0.487 7 0.553 2
灌浆期 0.536 8 0.060 4 0.424 5 0.307 0 0.581 2 0.637 2
不同生育期RSI与氮素含量的拟合结果见图5表3。在分蘖期和孕穗期,RSI-BG的R2最高,分别为0.352 9和0.628 1;在抽穗期、开花期、灌浆期,RSI-NR组合的R2最高,分别为0.707 5、0.542 9、0.639 5。
图5 Quick Bird模拟数据的RSI与氮素含量的拟合图
表3 不同生育期RSI与氮素线性拟合的决定系数R2
生育期 RSI-BG RSI-BR RSI-BN RSI-GR RSI-GN RSI-RN
分蘖期 0.352 9 0.206 2 0.285 0 0.252 9 0.293 5 0.339 4
孕穗期 0.628 1 0.250 0 0.541 6 0.522 6 0.571 9 0.587 0
抽穗期 0.606 6 0.395 9 0.616 5 0.609 3 0.707 5 0.707 5
开花期 0.532 1 0.081 2 0.415 5 0.355 1 0.483 1 0.542 9
灌浆期 0.533 7 0.056 4 0.421 3 0.307 7 0.573 1 0.639 5
基于Quick Bird卫星数据的NDSI和RSI线性拟合R2在不同生育期的趋势相似。在分蘖期和孕穗期表现较好的植被指数是由蓝光(B)和绿光(G)波段组成的NDSI-BG和RSI-BG;在生育后期(抽穗期、开花期、灌浆期),表现较好的植被指数则是由红光(R)和近红外(NIR)波段组成的NDSI-NR和RSI-NR。说明Quick Bird卫星数据在生育前期对B和G波段较为敏感,而在生育中后期则对R和NIR波段较为敏感。

2.3.2 SPOT-6模拟数据

基于SPOT-6模拟数据,不同生育期NDSI与氮素含量的拟合结果见图6表4。在分蘖期和孕穗期,NDSI-BR的R2最高,分别为0.351 8和0.612 9;而抽穗期、开花期、灌浆期的R2均以NDSI-GN最高,分别为0.702 9、0.543 6和0.584 6。
图6 SPOT-6模拟数据的NDSI与氮素含量的拟合图
表 4 不同生育期NDSI与氮素线性拟合R2
生育期 NDSI-BG NDSI-BR NDSI-BN NDSI-GR NDSI-GN NDSI-RN
分蘖期 0.050 0 0.350 0 0.290 0 0.183 9 0.330 0 0.300 0
孕穗期 0.000 9 0.612 9 0.542 7 0.436 1 0.586 3 0.578 0
抽穗期 0.383 6 0.604 3 0.630 8 0.237 8 0.702 9 0.619 9
开花期 0.032 4 0.528 0 0.425 4 0.246 9 0.543 6 0.488 9
灌浆期 0.228 2 0.014 8 0.136 1 0.076 2 0.584 6 0.579 5
RSI方面(图7表5),在分蘖期和孕穗期,RSI-BR的R2最高,分别为0.355 4和0.622 2;在抽穗期、开花期、灌浆期,R2最高的是RSI-GN,分别为0.697 8、0.537 9和0.584 8。
图7 SPOT-6模拟数据的RSI与氮素含量的拟合图
表5 不同生育期比植被指数与氮素线性拟合R2
生育期 RSI-BG RSI-BR RSI-BN RSI-GR RSI-GN RSI-RN
分蘖期 0.055 4 0.355 4 0.285 1 0.165 8 0.330 3 0.293 8
孕穗期 0.001 0 0.622 2 0.541 6 0.424 5 0.583 0 0.572 6
抽穗期 0.382 1 0.606 8 0.628 2 0.223 4 0.697 8 0.613 7
开花期 0.031 8 0.524 5 0.425 4 0.245 3 0.537 9 0.484 3
灌浆期 0.045 1 0.057 2 0.128 2 0.040 0 0.584 8 0.004 0
基于SPOT-6模拟数据的NDSI和RSI线性拟合R2在不同生育期的趋势相似。在生育前期(分蘖期和孕穗期),表现较好的植被指数是由B和R波段组成的NDSI-BR和RSI-BR;在生育后期(抽穗期、开花期、灌浆期),表现较好的植被指数则是由G和NIR波段组成的NDSI-GN和RSI-GN。说明SPOT-6卫星数据在生育前期对B和R波段较为敏感,而在生育中后期则对G和NIR波段较为敏感。

2.4 模型验证

2.4.1 Quick Bird模拟数据

图8可知,分蘖期训练集(tr)和测试集(te)的R2分别为0.38和0.30,RMSE均为0.02;孕穗期R2分别为0.58和0.61,RMSE分别为0.32和0.30;抽穗期R2为0.72和0.63,RMSE分别为0.30和0.34;开花期R2为0.55和0.47,RMSE分别为0.56和0.64;灌浆期R2为0.68和0.53,RMSE分别为0.29和0.38。鉴于分蘖期水稻生理状态的快速变化可能影响模型稳定性,建议评估时排除该期数据。在训练集中,R2在0.55~0.72,RMSE在0.29~0.56;在测试集中,R2在0.47~0.63,RMSE在0.30~0.64。说明抽穗期模型的氮素含量估测精度较高,而在开花期,模型的估测精度则相对较低。
图8 基于 Quick Bird卫星模拟数据各生育期氮素含量估测模型精度评估

2.4.2 SPOT-6模拟数据

图9所示,分蘖期训练集和测试集的R2分别为0.36和0.26,RMSE为0.40和0.36;孕穗期R2分别为0.61和0.58,RMSE分别为0.30和0.32;抽穗期R2为0.70和0.69,RMSE分别为0.36和0.38;开花期R2为0.55和0.44,RMSE均为0.03(出现异常);灌浆期R2为0.61和0.44,RMSE分别为0.31和0.40。由于分蘖期水稻生理的急剧变化可能影响模型的准确性,因此在评估模型时,排除该阶段数据。在排除异常值后不同生育期的性能评估显示,训练集的R2在0.55~0.70,RMSE在0.30~0.40;在测试集中,R2值在0.44~0.69,RMSE在0.32~0.40。结果表明,基于SPOT-6数据的氮素估测模型在抽穗期表现出最高的预测精度,在开花期的预测精度则相对较低。
图9 基于SPOT-6卫星模拟数据的各生育期氮素估测模型精度评估

2.4.3 两种估测模型精度比较

Quick Bird和SPOT-6卫星数据被用于建立氮含量估测模型,在孕穗期、抽穗期、开花期和灌浆期评估中,在训练集中,Quick Bird模型的平均R2为0.63,平均RMSE为0.37,而SPOT-6模型的平均R2为0.61,平均RMSE为0.26。在测试集中,Quick Bird的验证R2平均值为0.56,RMSE平均值为0.42;SPOT-6的验证R2为0.54,RMSE为0.28,与训练集结果趋势一致(表6)。综合表明,Quick Bird模型精度较高但误差大,SPOT-6模型效果略差,只适用于粗略估测。
表6 基于两种卫星模拟数据的氮素估测模型精度
卫星模拟 tr-R2 tr-RMSE te-R2 te-RMSE
Quick Bird 0.63 0.37 0.56 0.42
SPOT-6 0.61 0.26 0.54 0.28

3 结论

本研究使用Quick Bird和SPOT-6卫星的波段响应函数对实测数据进行重采样,以模拟两种传感器的波段反射率数据,构建了不同光谱指数,与氮素含量进行线性拟合。结果表明,水稻在5个生育期中,使用不同的光谱指数组合进行氮含量的相关性分析会产生差异。Quick Bird卫星传感器在分蘖期、孕穗期NDSI-BG和RSI-BR与氮含量相关性最高,在抽穗期、开花期、灌浆期NDSI-RN和RSI-RN的R2最高,SPOT-6卫星传感器与Quick Bird的相关性表现一致。说明在生育前期(分蘖期和孕穗期)蓝光波段对于氮素含量估测具有重要意义,在生育中后期(抽穗期、开花期、灌浆期)近红外光是氮素含量估测不可缺少的波段。在考虑生育期影响的前提下(去除分蘖期),对比分析了Quick Bird和SPOT-6卫星传感器对于水稻孕穗期、抽穗期、开花期、灌浆期的氮素含量估测精度。结果表明,Quick Bird模型和SPOT-6模型在抽穗期的精度较高,R2分别为0.72和0.70。整体而言,Quick Bird的R2平均值为0.56,RMSE平均值为0.42;SPOT-6的R2为0.54,RMSE为0.28,说明Quick Bird的氮含量估测模型比SPOT-6更精确,但RMSE值大,适合粗略估测。
[1]
李粉玲,常庆瑞,申健,等. 模拟多光谱卫星宽波段反射率的冬小麦叶片氮含量估算[J]. 农业机械学报201647(2):302-308.

[2]
史舟,梁宗正,杨媛媛,等. 农业遥感研究现状与展望[J]. 农业机械学报201546(2):247-260.

[3]
WU J D WANG D BAUER M E. Assessing broadband vegetation indices and QuickBird data in estimating leaf area index of corn and potato canopies[J]. Field crops research2007102(1):33-42.

[4]
HUANG S Y MIAO Y X ZHAO G M,et al. Satellite remote sensing-based in-season diagnosis of rice nitrogen status in Northeast China[J]. Remote sensing20157(8):10646-10667.

[5]
FERNANDEZ S VIDAL D SIMON E,et al. Radiometric characteristics of Triticum aestivum cv,Astral under water and nitrogen stress[J]. International journal of remote sensing199415(9):1867-1884.

[6]
STONE M L SOLIE J B RAUN W R,et al. Use of spectral radiance for correcting in-season fertilizer nitrogen deficiencies in winter wheat[J]. Transactions of the asae199639(5):1623-1631.

[7]
ZHU Y TIAN Y C YAO X,et al. Analysis of common canopy reflectance spectra for indicating leaf nitrogen concentrations in wheat and rice[J]. Plant production science200710(4):400-411.

[8]
WU J D WANG D ROSEN C J,et al. Comparison of petiole nitrate concentrations,SPAD chlorophyll readings,and QuickBird satellite imagery in detecting nitrogen status of potato canopies[J]. Field crops research2007101(1):96-103.

[9]
贾良良,李斐,陈新平,等. 应用IKONOS卫星影像监测冬小麦氮营养状况[J]. 中国土壤与肥料2013(6):68-71.

[10]
王备战,冯晓,温暖,等. 基于SPOT-5影像的冬小麦拔节期生物量及氮积累量监测[J]. 中国农业科学201245(15):3049-3057.

[11]
武婕,李玉环,李增兵,等. 基于SPOT-5遥感影像估算玉米成熟期地上生物量及其碳氮累积量[J]. 植物营养与肥料学报201420(1):64-74.

[12]
EITEL J U H LONG D S GESSLER P E,et al. Using in-situ measurements to evaluate the new RapidEye satellite series for prediction of wheat nitrogen status[J]. International journal of remote sensing200728(18):4183-4190.

[13]
陈拉,黄敬峰,王秀珍. 不同传感器的模拟植被指数对水稻叶面积指数的估测精度和敏感性分析[J]. 遥感学报200812(1):143-151.

[14]
杨海波,高兴,黄绍福,等. 基于卫星波段的马铃薯植株氮素含量估测[J]. 光谱学与光谱分析201939(9):2686.

Outlines

/