idl-programming-language - IDL 中不同分辨率的掩码/CLIP 光栅图像
问题描述
在剪辑/蒙版过程图像的分辨率应为 500 米之后,我想使用 IDL 编程将光栅图像(500 米分辨率)剪辑/蒙版为另一个光栅图像(10 公里分辨率)。我有 365 个图像对,我想通过 IDL 编程进行处理。任何人都可以在 IDL 中编写此代码吗?提前致谢。
解决方案
@MrFuppes
最初,我在 IDL 中编写了一个代码,它能够屏蔽/剪切相同分辨率的图像,因为它们的尺寸相同。但是我不能为不同的分辨率编写掩码函数,这在哪里是尺寸和分辨率的问题。
在我的代码中:
波段 1(500 米)是我的基本图像,我想用云遮罩图像(500 米)遮罩这个图像(效果很好),然后我想用 AOD_MASK 图像(10 公里)分辨率做另一个遮罩。但是如何为不同维度的图像编写 MASK 函数我不知道。如果您对此问题有任何建议,请建议我。
在这里我附上了我的测试数据,请检查它。
以下是我的 IDL 代码:
`Pro MASK
COMPILE_OPT idl2
ENVI, /RESTORE_BASE_SAVE_FILES
ENVI_BATCH_INIT, /NO_STATUS_WINDOW
print, ' ' + 'Batch Processing Start'
print, ' ' + systime()
;********************************************************************************
;Input location of different file
InputAOD = 'F:\DB_test_data\VAR2\TEST_DATA\Clip\'
Input_cloud = 'F:\DB_test_data\VAR2\TEST_DATA\Clip\'
Input_TOA_B1 = 'F:\DB_test_data\VAR2\TEST_DATA\Clip\'
;Output location
Output = 'F:\DB_test_data\OUT_PUT\'
; Search for inpiut data
FileCloud = File_search(Input_cloud + "*Aerosol_Cldmask_Land_Ocean-
Aerosol_Cldmask_Land_Ocean.tif",COUNT=Count_in1)
FileAOD = File_Search(InputAOD + '*Corrected_Optical_Depth_Land_2-
Corrected_Optical_Depth_Land.tif',COUNT=count_in)
FileBAND1 = File_Search(Input_TOA_B1 + '*pssgrp_000501330352.EV_500_RefSB_1-
EV_500_RefSB.tif',COUNT=COUNT)
nfiles = n_elements(count_in1)
;$$$$$$$$$$$$$$$$$$$$$$$$$$$$ $$$$$$$$$$$$$$$$$
FOR k=0, nfiles-1 DO BEGIN
Print, 'Processing File =', k +1
Print, 'Total Files=', count_in1
string_count = strtrim(nfiles,2)
;$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
;cloud(500_METER)
ENVI_OPEN_FILE, FileCloud[k], r_FiD=FiD_cld
IF (FiD_cld EQ -1) then return
ENVI_FILE_QUERY, FiD_cld, dims=dims, nb=nb, ns=cld_ncols, nl=cld_nrows
pos = [0]
map_info = envi_get_map_info(FiD=FiD_cld)
Cloud = Double(ENVI_GET_DATA(FiD=FiD_cld, dims=dims, pos=[0]))
;$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
;BAND_1(500_METER)
envi_open_file,FileBAND1[k],r_fid=fid_T_B1
if (fid_T_B1 eq -1)then return
envi_file_query,fid_T_B1,dims=dims,nb=nb, ns=toa_B1_ncols, nl=toa_B1_nrows
pos=[0]
map_info1=envi_get_map_info(fid=fid_T_B1)
layer1 = fltarr(toa_B1_ncols,toa_B1_nrows,nb)
FOR i=0,nb-1 DO BEGIN
layer1[*,*,i]=ENVI_GET_DATA(fid=fid_T_B1,dims=dims,pos=pos[i])
ENDFOR
;To replace 65535 with NAN
i=WHERE(layer1 eq 65535,count)
if (count gt 0)then layer1[i]=!VALUES.F_NAN
TOA_B1=layer1*0.000053811826d
;$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
;AOD(10KM_RESOLUTION)
ENVI_OPEN_FILE, FileAOD[k], r_FiD=FiD_AOD
IF (FiD_AOD EQ -1) then return
ENVI_FILE_QUERY, FiD_AOD, dims=dims, nb=nb, ns=aod_ncols, nl=aod_nrows
pos = [0]
map_info = envi_get_map_info(FiD=FiD_AOD)
layer2 = fltarr(aod_ncols,aod_nrows,nb)
FOR i=0,nb-1 DO BEGIN
layer2[*,*,i]=ENVI_GET_DATA(fid=FiD_AOD,dims=dims,pos=pos[i])
ENDFOR
;To replace -9999 with NAN
i=WHERE(layer2 eq -9999,count)
if (count gt 0)then layer2[i]=!VALUES.F_NAN
AOD=layer2*0.0010000000474974513d
;set the specific range for mask
AOD_MASK = 1.0*(AOD GE 0.0 and AOD LE 0.1) + (AOD GT 0.1)*0.0
;$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
; to get same dimensions for B1TOA and CLD data
if (cld_ncols le toa_B1_ncols) then begin
xsize = cld_ncols
endif else begin
xsize = toa_B1_ncols
endelse
if (cld_nrows le toa_B1_nrows) then begin
ysize = cld_nrows
endif else begin
ysize = toa_B1_nrows
endelse
;$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
; to call parameters again with same dims
Cloud = FSC_Resize_Image(Cloud,xsize, ysize)
Band1 = FSC_Resize_Image(TOA_B1,xsize, ysize)
;$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
;First masking processes with cloud mask images(500 meters)
Cloud_masked = Cloud * Band1
;To replace 0 with NAN
i=WHERE(Cloud_masked eq 0,count)
if (count gt 0)then Cloud_masked[i]=!VALUES.F_NAN
Cloud_masked_Band1=Cloud_masked
;Second masking processes with AOD_MASK images(10KM)
;$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
; Define output FileName
FileName = STRJOIN(STRSPLIT(STRMID(FILE_BASENAME(FileBAND1[k], '.tif'), 0, 75), '.', /EXTRACT), '_')
print, FileName
; write output File
out_name_out = Output + FileName +'_masked_image.tif'
map_info = envi_get_map_info(FiD=fid_T_B1)
ENVI_WRITE_ENVI_FILE, Cloud_masked, out_name=out_name_out, map_info=map_info
;FileName = STRJOIN(STRSPLIT(STRMID(FILE_BASENAME(FileAOD[k], '.tif'), 0, 75), '.', /EXTRACT), '_')
;print, FileName
; write output File
;out_name_out = Output + FileName +'_AOD.tif'
;map_info = envi_get_map_info(FiD=FiD_AOD)
;ENVI_WRITE_ENVI_FILE, AOD,out_name=out_name_out, map_info=map_info
;$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
Print, FileName + ' ' + strtrim(k+1,2) + ' of ' + string_count + ' Processed'
Endfor
Print, ' ' + systime()
Print, ' ' + 'Batch Processing End'
End
;$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
FUNCTION FSC_Resize_Image, image, cols, rows, $
INTERPOLATE=interp, $
MINUS_ONE=minus_one
Compile_Opt idl2
; Check parameters.
IF N_Params() EQ 0 THEN BEGIN
Print, 'USE SYNTAX: FSC_Resize_Image, image, cols, rows'
RETURN, "FSC_Resize_Image Syntax Error"
ENDIF
; Return to the caller on an error.
On_Error, 2
; Check image parameter for size. Only 2D and 3D images allowed.
ndim = SIZE(image, /N_DIMENSIONS)
dims = SIZE(image, /DIMENSIONS)
IF ((ndim LT 2) OR (ndim GT 3)) THEN $
Message, 'Input image must have two or three dimensions.'
; Check for keywords. Default for minus_one is 0.
interp = Keyword_Set(interp)
IF N_Elements(minus_one) EQ 0 THEN minus_one = 0
m1 = Keyword_Set(minus_one)
; 2D images are immediately passed to IDL's Congrid command, with
; the CENTER keyword set.
IF ndim EQ 2 THEN BEGIN
RETURN, Congrid(image, cols, rows, CENTER=1, MINUS_ONE=minus_one, $
INTERP=interp)
ENDIF
; 24-bit images are handled differently. The "3" dimension of a 24-bit
; image should not be interpolated.
offset = 0.5 ; To center the pixel location (i.e., CENTER=1 for Congrid)
index3 = Where(dims EQ 3)
CASE index3 OF
0: BEGIN
srx = [0,1,2]
sry = Float(dims[1] - m1) / (cols - m1) *(Findgen(cols) + offset) - offset
srz = Float(dims[2] - m1) / (rows - m1) *(Findgen(rows) + offset) - offset
END
1: BEGIN
srx = Float(dims[0] - m1) / (cols - m1) *(Findgen(cols) + offset) - offset
sry = [0,1,2]
srz = Float(dims[2] - m1) / (rows - m1) *(Findgen(rows) + offset) - offset
END
2: BEGIN
srx = Float(dims[0] - m1) / (cols - m1) *(Findgen(cols) + offset) - offset
sry = Float(dims[1] - m1) / (rows - m1) *(Findgen(rows) + offset) - offset
srz = [0,1,2]
END
ENDCASE
; Do the interpolation. Preserve nearest neighbor sampling, if required.
IF interp THEN BEGIN
RETURN, Interpolate(image, srx, sry, srz, /GRID)
ENDIF ELSE BEGIN
RETURN, Interpolate(image, Round(srx), Round(sry), Round(srz), /GRID)
ENDELSE
END`
推荐阅读
- yarnpkg - Yarn 如何指定包的安装路径
- java - MongoDB 通过 $in 命令获取匹配数
- qt - 从资源文件夹或 url 加载图像?
- javascript - 选中各自的行复选框时如何启用和禁用两个输入
- hadoop - 使用 Hive JDBC 驱动程序通过 Squirrel 连接到 Knox
- perl - c shell脚本不退出
- sql-server - 包括只返回一行而不是关闭 json
- jquery - JQuery 不适用于 Node.js
- c# - JSON.NET 没有突然转换
- python - Bees with Machine Guns 涉及 paramiko 的语法错误。(self._close(async=True))