专利名称：一种大面积作物种植面积及其产量的遥感检测估算方法
技术领域：
本发明涉及对大面积单一作物种植面积及其产量的遥感检测和估算方法。
背景技术：
目前，对单一作物纯像元识别，除了小范围的实地考察方法以外，多数都是利用遥感影象的监督或非监督分类方法，这些方法要么涉及大量的野外工作，要么对混合像元的界定十分不准确，而且选择监督分类时训练区也受人为因素影响很大，对训练样区的特征不能很好地把握，致使在遥感应用中不同尺度的影象转换很难满足应用的要求，这使在大范围的作物面积和产量估算中不得不使用精度不高的传统分类数据。
发明内容
本发明的目的是提供一种通过检测单一作物纯像元，对大面积作物种植面积及其产量进行估算的遥感检测新方法。
实现本发明的技术方案是对于大面积单一作物种植面积的遥感估算方法，其估算步骤如下(1)从卫星地面接收站获取估算区的TM卫星遥感图像；(2)参考波谱获取在估算区的地面预先选定若干同一作物的探测样区，获取各探测样区的实地冠层光谱，即参考波谱，所说的参考波谱通过对所述的各探测样区的作物实地测量获得(由同步ASD获取的实地被检测作物的冠层光谱，但需要通过TM的波段响应函数转换为TM对应波段的反射率值)或由波谱数据库根据先验知识和配套参数模拟获得；(3)影象波谱获取从遥感图像中找出与各探测样区相对应的作物像元样区，通过对遥感图像进行预处理，包括常规辐射、大气校正和几何处理，获得作物像元样区的影像均值波谱，即反射率矢量Ri，所说的反射率矢量Ri为宽波段TM除热红外波段之外的各波段反射率组成的向量；(4)计算距离阈值将获得的每个作物样区像元的影像波谱及其对应的参考波谱所形成的波谱对，通过下面的公式计算该波谱对的光谱距离阈值AV＝SQRT(∑((Ri-Rmi)^2)/N)式中，AV为距离阈值，Ri为作物像元样区的影像波谱反射率矢量(这里就是宽波段TM的各波段反射率组成的向量)，Rmi为参考波谱的反射率值均值矢量，N为波段数(这里针对TM为可见光、近红外6个波段即TM1-TM5和TM7)；(5)影象像元与距离阈值的比较将所计算的距离阈值AV与不同影像像元的影像波谱与参考波谱的距离AV’比较，如果某像元的影像波谱与参考波谱的距离计算值小于距离阈值AV则定义为纯作物TM像元，在计算不同影像像元的影像波谱与参考波谱的距离AV’时，其计算公式与计算光谱距离阈值AV的公式相同，不过公式中Ri是影像的像元反射率矢量；(6)从同期同地区的中分辨率成像光谱数据MODIS的每像元中统计相应TM纯作物像元的数量，确定MODIS影象中作物的种植成数；(7)根据种植成数和MODIS中每像元的面积确定MODIS中作物的种植面积。
在上述方案的基础上进一步估算被检测作物产量的方法如下通过影像-参考光谱距离阈值方法识别并统计MODIS每像元中TM纯作物像元数、然后结合每TM像元土壤调整的植被指数，确定每MODIS像元累积的植被指数SAVI，用以确定每MODIS像元作物的长势，再结合前面的面积估算，计算作物的产量，SAVI=Σi=1n[1.5*(NIRi-Redi)/(NIRi+Redi+0.5)]]]>上式中，n为每个MODIS像元的TM纯像元个数，NIR和Red分别表示卫星的近红外和红光波段反射率。
图1为本发明对作物纯像元检测方法的实施流程图。
对该流程的实施步骤说明如下(1)获取遥感图象；(2)利用常规方法进行遥感图象预处理，得到反射率影象；(3)在影象上寻找纯作物像元样区，并得到相应的影象波谱和参考波谱对(参考波谱可野外同步实测或从波谱库得到)；(4)通过距离阈值公式计算波谱对距离阈值；(5)影象上逐一像元与参考波谱对比，并计算距离，小于等于上述阈值则定义为纯像元。
本发明的特征在于由于反应地物特征的实测端元均值矢量取代传统分类法中训练样区的均值向量，因而更好地代表了目标地物光谱特征本质，这样遥感影象波谱和波谱库提取或实地获取的波谱对比才更具代表性。
上述所说的参考波谱还可以通过我国典型地物波谱库(spl.bnu.edu.cn)直接模拟得到；影象端元波谱可以通过地面先验知识(如土地利用现状图)采集。
本发明所采用的技术方案，因为参考波谱和对应影象端元波谱的获取都很确定，计算机程序自动实现，比野外方法省时省力效益高，也避免了传统监督分类通过概率等方法确定混合像元的缺点，能方便可靠地实现农作物纯像元的自动高精度检测，可以为较低分辨率影象的作物长势和产量估算监测决策所需要的尺度转换和其他遥感应用提供快速、准确的纯像元信息。
图1为本发明对作物纯像元检测方法的实施流程图；
图2为玉米实验采样点布局示意图；图3为地面测量的12个地块的玉米冠层光谱变异；图4为在124/34TM卫星轨道上实验区位置略图；图5为实验区TM影像453波段假彩色合成；图6为六个地点的玉米冠层均值波谱；图7为TM影像上提取的主要地物影像端元；图8为基于影像和测量光谱距离阈值的分类结果；图9为基于六个光谱测量地点的最大似然分类法的结果；图10为三个玉米纯作物像元检验地块在1∶10000地形图上的位置；图11为地块1两种分类结果对比(左A方法，右B方法)；图12为地块2两种分类结果对比(左A方法，右B方法)；图13为地块3两种分类结果对比(左A方法，右B方法)。
具体实施例方式
1.实验地点在中国科学院栾城农业生态系统试验站进行(如图2-1)，该站位于北纬37°53’，东经114°40’，海拔高度50.1m，地处石家庄市东南部的栾城县城东3KM处，种植制度为冬小麦—夏玉米1年两熟轮作高产区。试区面积为5km×5km，并与周围农田相联成片，为大面积均匀一致的夏玉米，有利于光谱数据准同步采集，并保证瞬时资料的精确可靠性。供试夏玉米品种为“郑单958”和“农大108”，按当地常规栽培措施管理。
2.实测光谱数据获取为使实验观测的数据样本具有一定的代表性，实验区落实面积为5*5平方公里，地面选择了六个测量点。
观测项目包括冠层光谱测量、农田小气候观测、作物结构参数、作物生化参数、土壤物理化学参数、背景光谱的测量、组分光谱、其他典型地物光谱测量、探空数据。
3.遥感图象获取使用的遥感数据是2003年9月13日成像的Landsat TM影像数据，空间分辨率是30m，使用的是除热红外波段外的可见光、近红外六个波段(TM1，TM2，TM3，TM4，TM5和TM7)。本研究选取了其中在地面有6个光谱测量实验地块的5KM×5KM区域为试验区(图示的范围略有扩大)。测试波谱和假彩色遥感影像见(图3，图4，图5)。
4.参考波谱在影像成像日期，我们获取了6个地点12个地块的玉米冠层和环境背景光谱参数，以及大气探空数据。经过处理，这12个地块-6个地点的平均光谱如图6所示。从图6可见，尤其是在近红外波段，其反射率波谱值的区间范围非常明显，其绝对值最大可达到0.2以上。
5.影像波谱的获取基于预处理后的影像，影像相关地面先验知识(如地形图、土地利用图等)和参考波谱的地面坐标，可以在TM影像上找到相应目标的位置，从而提取到六条目标影像波谱(图7)，玉米、居民地、苗圃、草地、果园、水体)，并得到玉米的均值波谱如表1。
表1实验区影像上提取的玉米影像波谱
6.计算距离阈值将获得的每个作物样区像元的影像波谱及其对应的参考波谱所形成的波谱对，通过下面的公式计算该波谱对的距离阈值AV＝SQRT(∑((Ri-Rmi)^2)/N)……………………(1)式中，Ri为作物像元样区的影像波谱反射率矢量，Rmi为参考波谱的反射率值均值矢量；对试验数据，我们得到该距离阈值为0.010726。
7.影象像元与距离阈值的比较将所计算的距离阈值AV与不同像元影象-参考波谱对的距离AV’比较，如果某像元的影像波谱与参考波谱的距离AV’的计算值小于距离阈值AV则定义为纯作物TM像元，这样就得到了实验区的纯作物像元图；把这一结果(图8)与传统最大似然监督分类结果(图9)的比较，表明距离阈值方法更符合地面情况。
8.距离阈值方法分类结果的验证通过获取当地实验区(约90平方公里)土地利用图，我们得到其耕地面积比例为87.3％。基于6个样点的训练区最大似然监督分类传统方法(以下简称A方法)计算的玉米地面积为37.8％；基于6个样点玉米波谱库方法(以下简称B方法)提取的纯玉米作物像元面积为55.6％，遗憾的是不能获取到相应的作物种植分布图进行检验。另外考虑到目前还没有纯像元的验证方法，为此，我们采用地面调查的真实地块进行效验，结果表明总体分类精度达到92％。
验证地块位于实验区东牛村南一个324亩地块(实地测量光谱样地1附近)和中国科学院农业生态实验站附近两个(样地6附近)97亩和94亩的两个纯作物像元区域((图10，三个长方形方框区)。
从两种分类方法的结果图上提取三个地块的分类结果区，我们得到如下的分类图和精度对比表(图11、图12、图13)，分类精度对比表，见表2。
表2玉米单一作物像元分类精度的对比 另外，对没有地面波谱对应的地区，可以借助从波谱库获取近似生长期的作物波谱或从波谱库通过模拟的方法来解决。
该方法虽然是针对玉米TM遥感影像纯作物像元识别的，但对其他作物和影像(只要光谱波段足够多，起码和TM相近)，只要能获取气候条件相似的生长期作物波谱，我们认为也是适用的，只要具备完备的波谱库支持就可以。
9.从同期同地区的中分辨率成像光谱数据MODIS的每像元中统计相应TM纯作物像元的数量。
10.根据MODIS中每像元的TM纯作物像元个数统计，获取MODIS影象中作物的种植成数 Mi为每个MODIS像元中的TM纯像元的个数，areapixel为每个像元的面积，Aimage为整个影象区的面积。
11.统计MODIS中每个行政区的种植成数和相应的区域面积就可以确定MODIS中各行政区的作物种植面积。
12.作物产量的估算每个TM像元作物的长势可以利用前人发明的土壤调整的植被指数SAVI的计算来定量，这样MODIS像元中作物的累积SAVI也可以得到，因为SAVI提供了在卫星影像中每个像元地块上绿色植被覆盖率，在某段作物生长阶段这是重要的作物长势信息，和作物产量具有重要关联；从而结合前面的面积估算，可以计算作物的产量SAVI=Σi=1n[1.5*(NIRi-Redi)/(NIRi+Redi+0.5)].........(2)]]>
Yield＝-208830*SAVI2+136104*SAVI-21667……………………(3)上式中，n为每个MODIS像元的TM纯像元个数，NIR和Red分别表示卫星的近红外和红光波段反射率。
根据我们对河北栾城实验区TM影象6个地点产量和SAVI的多项式非线性回归分析，通过上式(3)可以模拟9月中旬玉米的产量TM影象的SAVI与产量(Yield)的相关关系达到0.68。
权利要求
1.一种大面积作物种植面积的遥感估算方法，其估算步骤如下(1)从卫星地面接收站获取估算区的TM卫星遥感图像；(2)参考波谱获取在估算区的地面预先选定若干同一作物的探测样区，获取各探测样区的实地冠层光谱，即参考波谱，所说的参考波谱通过对所述的各探测样区的作物实地测量获得或由波谱数据库根据先验知识和配套参数模拟获得；(3)影象波谱获取从遥感图像中找出与各探测样区相对应的作物像元样区，通过对遥感图像进行预处理，包括常规辐射、大气校正和几何处理，获得作物像元样区的影像均值波谱，即反射率矢量Ri，所说的反射率矢量Ri为宽波段TM除热红外波段之外的各波段反射率组成的向量；(4)计算距离阈值将获得的每个作物样区像元的影像波谱及其对应的参考波谱所形成的波谱对，通过下面的公式计算该波谱对的光谱距离阈值AV＝SQRT(∑((Ri-Rmi)^2)/N)式中，Ri为作物像元样区的影像波谱反射率矢量，Rmi为参考波谱的反射率值均值矢量，N为波段数；(5)影象像元与距离阈值的比较将所计算的距离阈值AV与不同影像像元的影像波谱与参考波谱的距离AV’比较，如果某像元的影像波谱与参考波谱的距离AV’的计算值小于距离阈值AV，则定义为纯作物TM像元，在计算不同影像像元的影像波谱与参考波谱的距离AV’时，其计算公式与计算光谱距离阈值AV的公式相同，不过公式中Ri是影像的像元反射率矢量；(6)从同期同地区的中分辨率成像光谱数据MODIS的每像元中统计相应TM纯作物像元的数量，确定MODIS影象中作物的种植成数；(7)根据种植成数和MODIS中每像元的面积确定MODIS中作物的种植面积。
2.根据权利要求1所述的方法，其特征在于通过影像-参考光谱距离阈值方法识别并统计MODIS每像元中TM纯作物像元数、然后结合每TM像元土壤调整的植被指数，确定每MODIS像元累积的植被指数SAVI，用以确定每MODIS像元作物的长势，再结合前面的面积估算，计算作物的产量，SAVI=Σi=1n[1.5*(NIRi-Redi)/(NIRi+Redi+0.5)]]]>上式中，n为每个MODIS像元的TM纯像元个数，NIR和Red分别表示卫星的近红外和红光波段反射率。
全文摘要
本发明提供一种大面积作物种植面积及其产量的遥感检测估算方法，该方法是在对遥感影像进行常规辐射和大气校正的基础上，建立一个包括与遥感影像作物生长期一致的实地冠层光谱对应的纯作物影像光谱在内的光谱对，计算该光谱对的距离阈值，将影象像元波谱与参考波谱的距离与距离阈值比较，获得TM纯作物像元，从同期同地区的中分辨率成像光谱数据MODIS的每像元中统计相应TM纯作物像元的数量，确定MODIS影象中作物的种植成数，再根据种植成数和MODIS中每像元的面积确定MODIS中作物的种植面积。在此基础上结合每TM像元土壤调整的植被指数，确定每MODIS像元累积的植被指数SAVI，用以确定每MODIS像元作物的长势，计算作物的产量。
文档编号G01S17/00GK1731216SQ20051003656
公开日2006年2月8日 申请日期2005年8月19日 优先权日2005年8月19日
发明者陈水森, 柳钦火, 陈良富, 谭启宇, 方立刚 申请人:广州地理研究所