ENVI/IDL——获取AVIRIS数据波段值和经纬度

来源：互联网发布：mac 邮箱软件编辑：程序博客网时间：2024/05/16 06:19

以前，这一步，我都是手动导出的，详细请看另一篇博文：http://blog.chinaunix.net/uid-28490468-id-3462577.html

最近要处理较多的遥感数据，所以要把这一步通过代码自动化。

看了看help文档，发现ENVI支持一种叫IDL的语言，完全没有听说过，可能是我孤陋寡闻了，呵呵，有同感的举手。

有好玩的东西，那么就开始学呗。

打开ENVI+IDLX.X，编写IDL的IDE会自动和ENVI一起打开，这IDE很面熟，像eclipse。

IDL有自己的help文档，通过IDE可以打开。

另外网上也有些资料：

简洁明了： http://www2.geog.ucl.ac.uk/~mdisney/teaching/unix/idl/idl.html

一个国内的论坛，要注册，审核时间要几天（最好快点，呵呵）：http://bbs.esrichina-bj.cn/ESRI/forumdisplay.php?fid=28

各种资料胡乱看了一通，稍微总结几点：

IDL

1）IDL的数组是列优先的array[cols, rows]

2）数组是0开始的，和matlab不一样

3）procedure和function都叫routine，其中文件命名

File naming—a program file must be named the same as the main routine.

4）Primary routine must be the last in the file，这个关系到自动编译

5）对routine的参数传递有两种，argument和keyword

具体大家看help文档吧。

但是光有IDL这些基本的知识，根本无法处理像AVIRIS这种超光谱图片，我甚至尝试着去手动解析hdr头文件。

不过在愚蠢之余，试着google大法了一把，结果发现原来ENVI用IDL另外开发了一些API，特意用来处理硕大的遥感图片，

太激动了，呵呵。http://www.360doc.com/content/10/1101/20/472115_65790354.shtml

这些API在IDL的help文档里没有，要在ENVI的help文档里有，真是坑爹。

ENVI+IDL

怎么使用这些API来完成我的任务呢？我是参考：http://www.360doc.com/content/10/0405/21/472115_21741788.shtml

用第一题的做法就能完成我的任务了，不过读取波段值可以用另一个函数，更方便。

具体API

获取波段值

envi_open_file

filename, r_fid=fid, /no_realize, /no_interactive_query

if (fid eq - 1) then return

envi_file_query

fid, ns=ns, nl=nl, nb=nb, wl=wl

envi_get_slice

获取经纬度

envi_convert_file_coordinates

envi_convert_projection_coordinates

envi_get_projection

envi_proj_create

下面是我的实验代码：（要把图片和对应的hdr头文件放在一个文件目录下）

PRO file_test  data_dir = "E:\Data\OilLeakData\data\f100524t01p00r11rdn_b\"  filename = data_dir + "f100524t01p00r11rdn_b_sc01_ort_img"  envi_open_file, filename, r_fid=fid, /no_realize, /no_interactive_query  if (fid eq -1) then return ;;check fid  envi_file_query, fid, ns=ns, nl=nl, nb=nb  pos = lindgen(nb);  xs = 225  xe = 344  ys = 6183  ye = 6193  iproj = envi_get_projection(fid=fid)  oproj = envi_proj_create(/geographic)  print, 'content'  for y = ys, ye do begin    for x = xs, xe do begin      data = envi_get_slice(fid=fid, line=y, pos=pos, $        xs=x, xe=x, /bip)   ;; get bands' value      envi_convert_file_coordinates, fid, x, y, xmap, ymap, /to_map      envi_convert_projection_coordinates, xmap, ymap, iproj, oxmap, oymap, oproj ;;Lon = oxmap, Lat = oymap      help, x, xmap, oxmap, data      print, x, y, xmap, ymap, oxmap, oymap, data    endfor  endfor  envi_file_mng, id=fid, /remove ;;Don't forget to close fileEND