TIMESAT提取物候信息操作流程
软件环境:Matlab R2014a+TIMESAT3.2
数据介绍:MODIS A3或Q1的NVI(NDVI)均测试过这个流程,可行(大拇指)。
TIMESAT输入n年数据,提取n-1年的物候参数。通常用三年的数据,取中间一年的物候影像。因为软件无论提取的是像元的前两年物候,还是后两年,均有中间的年份,像元的物候更完整;还能保证是完整的物候周期,结果更准确。
如果是一年的数据,倒也是可以用一年的数据复制成三年,骗过软件。
本文介绍:操作过程中的小记录,害怕自己忘记,所以是“傻瓜式”教程。不涉及软件安装与配置,不涉及理论原理和软件原理,只是从准备TIMESAT可兼容的数据,到生成物候影像的操作流程。
操作流程
1数据准备
①研究区影像提取。进行反演的遥感影像最好是矩形的。因为不规则的裁剪,边缘像元的缺损容易使反演结果产生错误。在TIMESAT中,不规则裁剪可能会因为nodata值太多,生成时序曲线失败。所以解决方案就是使用包含研究区域的矩形影像。(如果是矩形区域还失败的话,可能是因为研究区沿海。扩大数据范围,降低水体在影像中的占比,可以解决。)
下图影像的范围就是进行作业的影像,矢量是研究区(河北省)。
②波段提取。影像一定是只是植被指数的单波段影像。这个就没什么好说的了,就八仙过海,各显神通吧。
③转化为dat或img格式。TIMESAT兼容这两种数据格式。使用ArcGIS中【栅格转其他格式(批量)】工具,将提取波段后的数据(.tif)转化为可兼容的格式(.dat)。
④建立数据列表。TIMESAT用数据列表的txt读取时间序列影像,所以需要建立数据列表txt,可以用excel实现。路径要完全正确,必须带有后缀,首行是影像景数。
2生成物候二进制文件
①timesat界面介绍。
②TSM_GUI生成时序曲线,保存设置文件(* .set)。
横坐标为影像期数,纵坐标为像元值。提取的点,不同专业有不同叫法。选择拟合函数,大部分用的是S-G,我用的是逻辑斯蒂;根据需要调整阈值提取点位,操作手册上建议的阈值是0.2,我用的是0.14。
·第一个点【返青期、生长季开始期等】就是植被指数曲线开始上升时,叶面积开始增大。
·第二个点【抽穗期、生殖生长转折期等】就是曲线最大值,叶面积登顶后开始减小。
·第三个点【成熟期、收获期】就是曲线下降到最小值,叶面积降至最小。
在这里我有个疑问,不过以后再解决吧。
如果从原理出发,NVI影像像元值的值域该是[-1,1],可我处理的影像像元值的值域是(-500,4000),查了文献都没有强调,看了修改值域的博文,嗯,并不是很明白原理,逻辑没有走通,所以先认为可能对提取物候结果没有影响
A.调整参数,生成时序曲线
B.保存设置文件(* .set)。
③TSF_process拟合物候参数文件( * _TS.tpa)。单击TSF_process选择刚刚保存的*.set文件,软件开始按行拟合。出现文件路径后,拟合完成
④TSF_seas2img生成物候二进制文件。单击TSF_seas2img选择拟合好的文件*_TS.tpa,接下来直接走图。
3生成物候影像
①二进制文件另存为TIFF文件。打开ENVI,打开一景同区域影像(需要它的头文件),以图示方式打开生成的二进制文件“*_s1”文件。
弹出Header Info对话框,输入头文件,设置参数。
设置完成后,小OK鼠标一点,影像就显示了。此时生成的影像是过程文件,需要另存为TIFF文件。(就不用多suo了吧)
②定义投影。此时物候影像没有坐标系,需要在ArcGIS中对TIFF影像【定义投影】。
③影像期数转化为天数。主要使用ArcGIS的栅格计算器
-剔除负值(生成二进制文件时设置的Nodata的值)。
-此时像元值表示为影像期数,所以值域应该在中间年份的影像期数之间。我的就应该是[24,46]。但实际像元值的值域是[0,43]。
-所以再剔除前23景影像的参数值,保证留下的是中间年份的。我用的方法是计算【影像-23】后,再次剔除负值。
-再以计算公式【影像*时间分辨率+1】,将期数转化为天数。
经过统计,返青期主要集中在81-162天之间。有异常值也是正常的,毕竟数据有云、有水、有积雪、有……,【但】像元数很少(也就是在误差范围内的意思)。
④按研究区裁剪。