; 步骤1:云检测
function Cloud_Detection, band4, band5, TIR
return, (band4 + band5 > 0.70) and (TIR < 240)
end
; 步骤2:水体检测
function Water_Detection, band3, band5, TIR
NDWI = (band3 - band5) / (band3 + band5)
return, (band5 < 0.17) and (NDWI > 0.10) and (TIR < 305)
end
; 步骤3:高温可疑像元检测
function High_Temperature_Suspect, TIR; vegetation_mask
return, (TIR > 260) ;and vegetation_mask
end
; 步骤4:背景像元的亮温统计
function Background_Temperature_Stats, TIR, suspect_fire_mask, window_size
dim = size(TIR, /dimensions)
bg_temp_avg = fltarr(dim[0], dim[1])
bg_temp_stddev = fltarr(dim[0], dim[1])
for i=0, dim[0]-1 do begin
for j=0, dim[1]-1 do begin
if suspect_fire_mask[i, j] ne 0 then begin
x_start = max(i-window_size, 0)
x_end = min(i+window_size, dim[0]-1)
y_start = max(j-window_size, 0)
y_end = min(j+window_size, dim[1]-1)
window = TIR[x_start:x_end, y_start:y_end]
; 检查window是否为空
if n_elements(window) eq 0 or n_elements(where(window ne TIR[i,j])) eq 0 then continue
window = window[where(window ne TIR[i,j])]
bg_temp_avg[i,j] = mean(window)
bg_temp_stddev[i,j] = stddev(window)
endif
endfor
endfor
return, [bg_temp_avg, bg_temp_stddev]
end
; 步骤5:着火像元确认
function Fire_Pixel_Confirmation, TIR, bg_temp_avg, bg_temp_stddev, suspect_fire_mask
dim = size(TIR, /dimensions)
fire_mask = bytarr(dim[0], dim[1])
for i=0, dim[0]-1 do begin
for j=0, dim[1]-1 do begin
if suspect_fire_mask[i, j] and TIR[i,j] > (bg_temp_avg[i,j] + 3*bg_temp_stddev[i,j]) then begin
fire_mask[i,j] = 1
endif
endfor
endfor
return, fire_mask
end
; 主程序
pro Main
; 加载和预处理GF-4 PMI数据,假设已完成,各波段赋值给相应变量
; 加载植被掩膜,假设为vegetation_mask
file = 'E:\桌面1\大创\dataprocess\预处理后数据\贵州\08\MERGE.tif'
image = READ_TIFF(file)
; 注意:您需要根据实际情况调整下面的代码来匹配您数据的波段结构
band3 = image[*, *, 2]
band4 = image[*, *, 3]
band5 = image[*, *, 4]
TIR = image[*, *, 5] ; 假设TIR为第6波段
cloud_mask = Cloud_Detection(band4, band5, TIR)
water_mask = Water_Detection(band3, band5, TIR)
suspect_fire_mask = High_Temperature_Suspect(TIR)
result = Background_Temperature_Stats(TIR, suspect_fire_mask, 2)
bg_temp_avg = result[0]
bg_temp_stddev = result[1]
fire_mask = Fire_Pixel_Confirmation(TIR, bg_temp_avg, bg_temp_stddev, suspect_fire_mask)
; 输出结果
TVScl, fire_mask
end
代码应该是有问题 显示如图
请问应该怎么解决