重点还是掌握envitask,此过程我自己写的是按照之前lai反演的模式,并没有引入新的方法来完成,但可以引入envi——file_query可以查询头文件中的bias和offset.两种算法算出的差距不大,甚至在第四位小数的时候还是一样的结果。
本篇利用lai反演检测的过程 利用function进行计算vfc,在计算对应像元的比辐射率,得到不同landcover下的表面亮温,在根据反射亮度l10(第十通道)结合水汽进行反演。
水汽和比辐射率还用于算一些连带参数,加强计算的可靠性。
温度反演
pro try
;fn=dialog_pickfile()
;读取数据
fn_TIRS='TIRS’
envi_open_file,fn_TIRS,r_fid=fid_TIRS
envi_file_query,fid_TIRS,dims=dims, data_gains=gains, data_offsets=offsets
map_info=envi_get_map_info(fidfid_TIRS)
fn_Landcover='Landcover'
envi_open_file,fn_Landcover,r_fid=fid_Landcover
envi_file_query, fid_Landcover, dims=dims
fn_Fv='VFC'
envi_open_file,fn_Fv,r_fid=fid_Fv
envi_file_query, fid_Fv, dims=dims
;辐射定标
DN=envi_get_data(fid=fid_TIRS, dims=dims, pos=0)
L=gains[0]*DN+offsets[O]
Landcover=envi_get_data(fid=fid_landcover,dims=dims, pos-0)
Fv=envi_get_data(fid=fid_Fv,dims=dims, pos=0)
Emiss=cal_Emiss(Fv,Landcover)
;计算地表温度
Lup=3.09
Ldown=4.81
Trans=0.66
Ts=cal_Ts(L, Emiss, Lup, Ldown, Trans)
;保存结果
envi_write_envi_file,Emiss,out_name='emiss',map_info=map_info
envi_write_envi_file,ts,out_name='ts',map_info=map_info
end
function cal_Ts,L, Emiss, Lup, Ldown, Trans
K1=774.8853
K2=1321.0789
Ls=(L-(1-Emiss)*Ldown*Trans-Lup)/Emiss*Trans
Ts=K2/alog(1+K1/Ls)
return, Ts
end
function cal_Emiss,Fv,Landcover
dE=0.0038*Fy
w=where(Fv gt 0.5 ,count)
dE[w]=0.0038*(1-Fv[w])
Rv=0.9332+0.0565*Fv
Rs = 0.9902+0.1068*Fv
Rm=0.9886+0.1287*Fv
Emiss=Fv*0+0.9968
w=where(Landcover eq 1)
Emiss[w]=0.9867*Fv[w]*Rv[w]+0.9649*(1-Fv[w])*Rm[w]+dE
w=where(Landcover eq 2)
Emiss[w]=0.9867*Fv[w]*Rv[w]+0.9677*(1-Fv[w])*Rm[w]+dE
return, Emiss
end