您的位置:首页 > 编程语言

IDL波段分解与合成源代码

2017-06-30 10:00 106 查看
转载:http://www.cnblogs.com/zhzhx/p/3141031.html

一、波段的分解:

pro band_decompose
compile_opt IDL2

envi,/Restore_Base_Save_Files
envi_batch_init,log_file='batch.txt'

inputfile='E:\L5_12332_20040908\tm'
envi_open_file,inputfile,r_fid=fid
if(fid eq -1) then begin
envi_batch_exit
return
endif

envi_file_query,fid,dims=dims,nl=nl,ns=ns,data_type=data_type

;get the map projection information from the resource file
mapinfo=envi_get_map_info(fid=fid)

;output to memory and to other format.
;data=envi_get_data(fid=fid,dims=dims,pos=[0])
;envi_enter_data,data,r_fid=fid
;out_name='E:\L5_12332_20040908\Band\band1_otherformat.tif'
;envi_output_to_external_format,dims=dims,fid=fid,out_name=out_name,pos=[0],/tiff

;output to ENVI format
data=envi_get_data(fid=fid,dims=dims,pos=[7])

outputfile='E:\L5_12332_20040908\Band\band2_1'
;outputfile=data

OpenW, unit, outputfile, /Get_LUN
WriteU, unit, data
Free_LUN, unit

;create map projection information
;ps=[30,30]
;mc=[342198.231,4489929.046]
;proj=envi_proj_create(/geographic)
;mapinfo=envi_map_info_create($
;name='Geographic',$
;mc=mc,$
;ps=ps,$
;proj=proj,$
;/GEOGRAPHIC)

;create the header of file
ENVI_SETUP_HEAD, fname=outputfile, $
ns=ns, nl=nl, nb=1, $
interleave=0, data_type=data_type, $
map_info=mapinfo,$
offset=0, /write, /open

end


可以保持为ENVI格式或者其他格式,代码中有注释。

波段分解批处理:

pro band_decompose_batch
compile_opt idl2

envi,/Resotre_Base_Save_Files
envi_batch_init,log_file='batch.txt'

inputfilefolder=dialog_pickfile(/directory,title='Choose the input directory!')
;inputfile=dialog_pickfile(path=EXPAND_PATH('<IDL_DIR>')+PATH_SEP()+'examples\data',title='打开文件')
outputfilefolder=dialog_pickfile(/directory,title='Choose the output directory!')

inputfilefolderfiles=file_search(inputfilefolder,'*.img',count=num)
for j=0,num-1 do begin
inputfile=inputfilefolderfiles[j]

envi_open_file,inputfile,r_fid=fid
if(fid eq -1) then begin
envi_batch_exit
return
endif

envi_file_query,fid,bnames=bnames,dims=dims,nl=nl,ns=ns,nb=nb,data_type=data_type

mapinfo=envi_get_map_info(fid=fid)

for i=0,nb-1 do begin
data=envi_get_data(fid=fid,dims=dims,pos=[i])

;outputfile=outputfile1+'band_'+strmid(string(i+1),strlen(string(i+1))-1,1)
outputfilename1=file_basename(inputfile)
outputfilename2=strmid(outputfilename1,0,strlen(outputfilename1)-4)
outputfile=outputfilefolder+outputfilename2+'_'+bnames[i]

OpenW, unit, outputfile, /Get_LUN
WriteU, unit, data
Free_LUN, unit

ENVI_SETUP_HEAD, fname=outputfile, $
ns=ns, nl=nl, nb=1, $
interleave=0, data_type=data_type, $
map_info=mapinfo,$
offset=0, /write, /open

endfor

endfor

end


二、波段合成

IDL帮助文档里的一种方法:

PRO Band_ENVI_LAYER_STACKING_DOIT

compile_opt IDL2

;

; First restore all the base save files.

;

envi, /restore_base_save_files

;

; Initialize ENVI Classic and send all errors

; and warnings to the file batch.txt

;

envi_batch_init, log_file='batch.txt'

;

; Open the first input file. Also open the

; single-band DEM file to layer stack with

; this file.

;

envi_open_file, 'E:\L5_12332_20040908\Band\tm_band1', r_fid=t_fid

if (t_fid eq -1) then begin

envi_batch_exit

return

endif

;

; Open the second input file.

;

envi_open_file, 'E:\L5_12332_20040908\Band\tm_band2', r_fid=d_fid

if (d_fid eq -1) then begin

envi_batch_exit

return

endif

;

; Use all the bands from both files

; and all spatial pixels. First build the

; array of FID, POS and DIMS for both

; files.

;

envi_file_query, t_fid, $

ns=t_ns, nl=t_nl, nb=t_nb

envi_file_query, d_fid, $

ns=d_ns, nl=d_nl, nb=d_nb

nb = t_nb + d_nb

fid = lonarr(nb)

pos = lonarr(nb)

dims = lonarr(5,nb)

for i=0L,t_nb-1 do begin

fid[i] = t_fid

pos[i] = i

dims[0,i] = [-1,0,t_ns-1,0,t_nl-1]

endfor

;

for i=t_nb,nb-1 do begin

fid[i] = d_fid

pos[i] = i-t_nb

dims[0,i] = [-1,0,d_ns-1,0,d_nl-1]

endfor

;

; Set the output projection and

; pixel size from the TM file. Save

; the result to disk and use floating

; point output data.

;

out_proj = envi_get_projection(fid=t_fid, $

pixel_size=out_ps)

out_name = 'C:\Band1_2'

out_dt = 4

;

; Call the layer stacking routine. Do not

; set the exclusive keyword allow for an

; inclusive result. Use cubic convolution

; for the interpolation method.

;

envi_doit, 'envi_layer_stacking_doit', $

fid=fid, pos=pos, dims=dims, $

out_dt=out_dt, $

out_name=out_name, $

interp=2, $

out_ps=out_ps, $

out_proj=out_proj, r_fid=r_fid

;

; Exit ENVI Classic

;

;envi_batch_exit

END


另一种方法:

pro band_envi_layer_stacking_doit_test
compile_opt idl2

envi,/restore_base_save_files
envi_batch_init,log_file='batch.txt'

inputfiles=['E:\L5_12332_20040908\Band\tm_band1','E:\L5_12332_20040908\Band\tm_band2']
outputfile='c:\tm_band_12'

fids=lonarr(n_elements(inputfiles))
dimses=lonarr(5,n_elements(inputfiles))
poses=lonarr(n_elements(inputfiles))

for i=0,n_elements(inputfiles)-1 do begin
envi_open_file,inputfiles[i],r_fid=fids1
envi_file_query,fids1,ns=ns,nl=nl,nb=nb

fids[i]=fids1
dimses[0,i]=[-1,0,ns-1,0,nl-1]
proj=envi_get_projection(fid=fids,pixel_size=out_ps)

endfor

envi_doit,'envi_layer_stacking_doit',$
fid=fids,pos=poses,dims=dimses,$
out_dt=4,out_name=outputfile,$
interp=2,out_ps=out_ps,$
out_proj=proj,r_fid=r_fid

end
内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息