IDL计算NDVI(landsat8 OLI)。。。

其实用来用去,还是觉得python最上手,R语言我怎么都爱不起来,虽然有时候R比python简便很多,IDL只是一开始的时候无意中接触的,使用频率不高,使用也不顺手。

以下在计算NDVI的过程中使用的是反射率,不以原始的像元DN值作为计算方式,实际使用过程中也可以增加大气校正,该过程比较适用于进行批处理计算的基础模板和附件。裁剪处理并未附带。。。行列可以通过获取的方式,从而改变一开始手工输入

pro aaabbb_44

work_dir=dialog_pickfile(title='D:\idltemp\nb1',/directory);设置工作空间

cd,work_dir

envi,/restore_base_save_files;导入envi函数

compile_opt idl2;强制执行function

envi_batch_init;ENVI_GET_DATA(fid=fid, dims=dims, pos=0)就不会出错

;----

fna=dialog_pickfile(title='D:\idltemp\nb1\2010nb1.dat');;;;;输入dat影像,不知道为什么dat无法输入

envi_open_file,fna,r_fid=fid

envi_file_query,fid, ns=ns, nl=nl, nb=nb, dims=dims, data_type=data_type, interleave=interleave, offset=offset

map_info=envi_get_map_info(fid=fid)

b3= ENVI_GET_DATA(fid=fid, dims=dims, pos=3)

b4= ENVI_GET_DATA(fid=fid, dims=dims, pos=4)

Date=[2017,4,2]

zenith=23.24*!dtor

gain=[0.012869,0.0099997,0.0061193,0.00051293];;;;;;;;

bias=[-64.66,-59.2,-50.14,-2.49];;;;;;;;;

esun=[0.9995528,0.9995528,0.9995528,0.9995528];;;;;;;;;;;

se=56.5759*!dtor

apr_r=cal_apr(b3,gain[1],bias[1],se)

apr_nir=cal_apr(b4,gain[2],bias[2],se)

;ndvi=fltarr(ns,nl)

ndvi=(apr_nir-apr_r)/(apr_nir+apr_r)

o_fn='D:\idl\526.tif'

envi_write_envi_file,ndvi,out_name=o_fn,ns=ns,nl=nl,nb=1,data_type=4,interleave=0,offest=offset,map_info=map_info

end

function cal_apr,data,gain,bias,se

result=(gain*data+bias)/(sin(se))

return,result;计算表观反射率

end