在面对大批量遥感影像数据重复操作的时候,我们会想到批处理的方式。尽管遥感软件提供了一些批处理的方式,就小部分需求而言,单一的批处理方式往往是不够的,这时候程序化处理就派上用场了。
(当然,也可以使用建模的方式做这个事情)
使用程序化处理的好处在于可以将一个复杂的地理处理工作流一次执行完成,灵活,可以按照自己需求来定,而无需中间过多的数据拷贝、转换处理过程。可是另一个问题也接踵而至,当我们借助于程序的时候,ENVI+IDL那些恼人、庞杂的程序知识也迫使大部分人放弃这一手段。如果使用ENVI+IDL的话,需要学习IDL语言、ENVI软件,万一使用的人跑偏了呢,那就不得不在IDL编辑器、面向过程和面向对象程序设计、界面编程、可视化编程与ArcGIS集成的沼泽里苦苦挣扎了。
(其实我们的目的很简单,只要把这批数据处理好了就行啊)
对于一个批处理地理任务问题而言,无外乎就是数据的输入、数据处理和数据导出三个小任务(这一思路适用于所有的数据程序化处理过程),至于那些界面、可视化的编码实在是根本用不到有木有。多年前,IDL使用的是过程编码,后续面向对象兴起,马上转投面向对象编码,然后又突然告诉使用的这帮人,面向对象比面向过程方便,也是让使用ENVI+IDL的人又心累了一把。
(带别人跳坑,然后又把人从坑里捞出来,你不嫌累啊)
本文使用面向对象的程序设计,大致讲述ENVI+IDL编码过程,尽量省去那些界面编程、可视化编程的繁琐细节,内容有:IDL工作台、数据输入、数据处理、数据导出四个方面。至于IDL编码规范也会部分涉及,但此项不是重点。
1、IDL编辑器使用
在安装ENVI的时候会自动安装IDL编程软件。在软件安装完毕后,在开始菜单上会新增ENVI和IDL的应用程序,点击IDL8.5(ENVI版本不一样,IDL版本也会不一样),就可以看到IDL的工作台了,工作台主要有资源管理器、过程编辑界面和控制台三个。可以在菜单栏的新建文件按钮上,点击新建文件,写入如下语句:
pro hellow
e=envi()
;核心语句,获取一个ENVI对象,加载到内存中。
end
(注:IDL过程语句语法结构为:pro xxx end )
点击右上角的编译、运行,一个ENVI对象就启动运行了。
在程序执行完毕的时候,是需要回收ENVI内存资源的,应加上ENVI对象关闭语句:e.Close
pro hellow
e=envi()
e.Close
end
e=envi()和e.Close之间就是编写具体的代码了。
2、数据的输入(数据访问)
(1)单个数据输入
面向对象ENVI+IDL本身也是有数据访问函数的,查看ENVI程序设计帮助,如下:
单个数据导入语法为:Result = ENVI.OpenRaster(URI [, Keywords=value]),该语句返回一个ENVI栅格数据对象。
注:ENVI函数的形参有必选参数和非必选参数,必选参数和非必选参数用逗号隔开,位于中括号外面的为必选参数(教程一般说是位置参数),中括号里面的为非必选参数(教程说是关键词),部分非必选参数是不能同时存在的,只能选一个。
URL:位置参数,对于本地数据而言,URL是遥感影像数据在本地磁盘的详细路径,对于远程数据而言,URL是一个远程访问地址,如果带有访问权限的远程数据,还需要输入用户名和密码。
Keywords :关键词,指在创建这个栅格数据对象时的一些初始化设置。如数据类型、需要忽略的值,数据集名字、扩展类型、元数据、远程数据用户和密码等等。
假如数据路径为:D:\Program Files\Exelis\ENVI53\data\qb_boulder_msi.hdr
那么该数据的导入为:
e = ENVI()
;如果在参数里加/headless 关键字,为不显示ENVI软件
file='D:\Program Files\Exelis\ENVI53\data\qb_boulder_msi.hdr'
raster=e.OpenRaster(file)
如果是使用单个数据集,写好这一步就算完成数据的导入了。
(可以再加两条代码测试,如:view1 = e.GetView()
;获取视图对象
layer1 = view1.CreateLayer(raster)
;在视图中显示栅格数据 )
(注:与C/C++、JAVA、python不同的是,反斜杆并不是转义符,IDL路径分割符斜杆与反斜杠都是一样的,也可以写成这样:file='D:/Program Files/Exelis/ENVI53/data/qb_boulder_msi.hdr',但是IDL函数关键词引用为斜杆,IDL不区分字母大小写)
(2)批量数据导入
面向对象ENVI+IDL并不支持批量数据的导入,如果需要批量导入数据,需要使用IDL文件操作函数配合ENVI::OpenRaster一起使用。其意义在于遍历文件夹下所有数据,并且使用通配符进行查询筛选,将所需要的数据路径保存至一个字符串数组里,然后再使用循环语句遍历该字符串数组,通过ENVI::OpenRaster依次打开。
(IDL文件操作函数:
比较常用的文件函数有:FILE_COPY(拷贝文件),FILE_SEARCH (文件搜索)、FILEPATH(文件完全路径)。)
批量导入数据使用的核心文件操作函数是FILE_SEARCH (文件搜索)。
语法为:Result = FILE_SEARCH(Path_Specification),该函数返回一个字符串数组,参数Path_Specification文件夹路径,并且可以使用通配符,使得字符串数组包含所有Path_Specification文件夹匹配的文件路径。IDL通配符与shell通配符一致,有:* ?{}等。
比如D:\Program Files\Exelis\ENVI53\data文件夹下:
搜索所有带后缀名为.hdr的文件:
filesearch='D:\Program Files\Exelis\ENVI53\data'
files=FILE_SEARCH(filesearch,'*.hdr')
查询files字符串结果如下:
(补充通配符搜索:
如果是搜索所有hdr文件中含有set的文件:
files=FILE_SEARCH(filesearch,'*set*.hdr')
如果是搜索所有包含.hdr和.enp后缀的文件:
files=FILE_SEARCH(filesearch,'*.{hdr,enp}')
)
以上也可以使用界面DIALOG_PICKFILE函数来实现多个文件的选择,按住ctrl键手动选择多个文件,返回一个字符串数组。
files = DIALOG_PICKFILE(/MULTIPLE_FILES)
由于是先在遍历该文件夹下文件路径的时候,未获取文件数目多少,可以使用数组处理函数来获取整个的文件数目。
1)size()函数,返回数据信息的一个整型数组,整型数组第二个数字存储了该整型数组的大小:
stringsize=size(files)
count1=stringsize(1)
或者使用/N_elements关键字参数获取字符串数组大小
count2=size(files,/N_elements)
2)N_elements()函数,返回数组大小
count3=N_elements(files)
3)在使用FILE_SEARCH()函数时,直接将关键字遍历文件数目count=xx,传出。
files=FILE_SEARCH(filesearch,'*.hdr',count = count4)
3、数据处理
(1)循环语句
获取到符合条件的影像路径后就可以使用循环语句进行影像数据的遍历与任务执行了。
循环语句语法:
1)for语句 for i=m,n(,step) do begin 语句 endfor ;
当i从m以步长inc(若不设置,默认为1)逐步变化到n时,循环执行语句。
2)while语句 while 条件 do begin 语句 endwhile ;
当条件表达式的值为真时,while语句循环体执行语句序列,先判断循环条件,如果条件满足,再执行循环体,否则跳出。
3)repeat语句 repeat begin 语句 endrep until 条件
先执行循环体,然后在判断条件表达式。哪怕是条件表达式为假,循环体也要执行一次。
在这里使用第一种循环语句,加入参数m,是为了查看是否循环执行完毕。
pro hello
e=envi(/headless)
filesearch='D:\Program Files\Exelis\ENVI53\data'
files=FILE_SEARCH(filesearch,'*.hdr',count=countall)
m=0
for i=0,countall-1 do begin
tmp=e.OpenRaster(files[i])
m = m+1
endfor
print,m
end
m参数最后结果是16,说明执行完毕。
(2)任务处理(TASK)
ENVI的任务处理包含了众多的TASK,如:影像平滑、增强、几何纠正、光谱纠正、波段运算、变化检测、融合、配准、拼接、掩模、裁剪,分类、分类后处理(矢量栅格转化、矢量小面积处理)、格式转换等等,ENVI5.3有将近150个TASK模块。
由于各模块需要的参数不一致,TASK的参数装配千差万别,为说明TASK执行任务过程,以下以一个栅格数据为例说明TASK装配过程:
1)引用TASK任务,如:百分比线性拉伸
Task=ENVITask('LinearPercentStretchRaster')
2)装配TASK任务参数,TASK任务参数既包括输入、输出参数,也包含该任务的特殊参数,比如百分比线性拉伸按照多少百分比的像素值分布对原始图像进行拉伸。
pro hello
compile_opt idl2
e=envi()
;导入栅格数据
file='D:\Program Files\Exelis\ENVI53\data\qb_boulder_msi.hdr'
raster=e.OpenRaster(file)
;引入LinearPercentStretchRaster TASK任务
Task=envitask('LinearPercentStretchRaster')
;设置TASK参数
Task.INPUT_RASTER=raster
Task.PERCENT = 10.0
Task.OUTPUT_RASTER_URI='f:\file\qb_boulder_LinearPercentStretch'
;执行TASK任务
Task.Execute
end
3)执行TASK
Task.Execute
10%线性拉伸效果如下:
(3)数据导出
在执行TASK任务之前,只要设置好Task.OUTPUT_RASTER_URI的输出路径,即可将数据输出。由于执行的是批量处理任务,因此,在输出路径时,还需要自动的获取输出路径。由于路径是一个个的字符串,可以使用字符串处理函数来截取数据。假定其中一个原始栅格数据的路径为:file='D:\Program Files\Exelis\ENVI53\data\qb_boulder_msi.hdr',截取文件名方法如下:
1)直接使用文件操作函数FILE_BASENAME()截取,返回文件名
str=FILE_BASENAME(file,'.hdr')
2)我们需要截取该栅格数据的文件名qb_boulder_msi,由于qb_boulder_msi位于反斜杠和点号之间,因此可以通过这三者之间的关系得到。
获取最右边斜杠的位置:
pos=strpos(file,'\',/reverse_search)
截取反斜杠右边的字符串:
str1=strmid(file,pos+1)
获取点号的位置:
pos=strpos(str1,'.',/reverse_search)
截取文件名:
str=strmid(str1,0,pos)
上面两种方法用到了文件操作函数或者字符串处理函数,这是在IDL处理路径字符串中最常用的方法。
假如指定批量输出的文件夹路径为:files1='F:\file' 那么输出的完整路径应为:
file=files1+'\'+str
最后再装配到TASK处理路径中
Task.OUTPUT_RASTER_URI=files2
因此,配合循环,整个执行过程为:
pro hello
compile_opt idl2
e=envi(/headless)
filesearch='D:\Program Files\Exelis\ENVI53\data'
files=FILE_SEARCH(filesearch,'*.hdr',count=countall)
files1='F:\file'
for i=0,countall-1 do begin
str=FILE_BASENAME(files[i],'.hdr')
file=files1+'\'+str
raster=e.OpenRaster(files[i])
Task=envitask('LinearPercentStretchRaster')
Task.INPUT_RASTER=raster
Task.PERCENT = 10.0
Task.OUTPUT_RASTER_URI=file
Task.Execute
endfor
e.close
end
处理结果如下:
上述过程中都是一些代码的搬运过程,整个任务流程也不需要什么深层干预,因为ENVI+IDL已经将相关函数封装好了。如果真正涉及到内部算法,没有深厚的数据结构功底,纯粹做底层算法研发,也实在是痴人说梦。从某个角度说,这或许是如今python如此流行的原因,基本上在网上弄一些函数库,然后几行代码就搞定了,至于说研究,不是搞笑嘛,各种COPY,各种搬运。。。
波段运算对应数组处理
文件路径对应字符串处理
控制语句对应批处理
以上三个需要与ENVI的函数灵活搭配。