generate_meanburnsize.pro

generate_meanburnsize.pro

posted 05-09-2001 02:14 PM
A script that inputs a “Percent of pixel burned” image and a “number of burnscar centroids in pixel” image, and divides area/count to get average burnscar area per pixel. If burn count pixel is 0 (or close to it), then the output pixel is also zero.
Also generates metadata in a .txt file for the run.

Thanks, Nick

code:

; PRO GENERATE_MEANBURNSIZE
;
; Initial hopes: This program will input two files:
;
; 1 km raster with % of pixel burned (From Landsat mapped burnscars)
; 1 km raster with # of distinct burnscars mapped.
;
; Both images must have the same dimensions, etc, and probably must be floating point
;
; well as ENVI_SELECT, and DIALOG_PICKFILE to allow the user to select input/output
; files
;
;
; Started May, 2001.
; written by Nick Matzke
pro generate_meanburnsize

; Get first image:
; use the envi_select widget to select the image that contains the data,
; fid will be file id of the open file
; pos will be the number of bands
envi_select, title=’Input Filename of burn area image’, fid=fid, pos=pos, dims=dims, /NO_SPEC;, /NO_DIMS

;use envi_file_query to find out number of rows, samples, bands,
;byte order, and note interleave: 0= bsq, 1=bil, 2=bip
;note data type: 1=byt, 2=integer, 4=floating point
; (see progguid.pdf page 330 for more)
envi_file_query, fid, data_type=data_type, nb=nb, interleave=interleave, sname=image1name, fname=fullpath;, h_map=h_map
;handle_value, h_map, map_info

; Get second image:
envi_select, title=’Input Filename of burn area image’, fid=fid2, pos=pos2, dims=dims2, /NO_SPEC;, /NO_DIMS
envi_file_query, fid2, data_type=data_type, nb=nb, interleave=interleave, sname=image2name, fname=fullpath;, h_map=h_map
;handle_value, h_map, map_info

;calculate the number of lines and samples of the input image (this deals with subsetting)
ns = dims(2)-dims(1)+1
nl = dims(4)-dims(3)+1

;use pickfile to select the output image (starting path will be same as input ascii rulesfile)
outimage = dialog_pickfile(title=’Enter Filename of Output Mean Burn Area Image (add *.img if needed)’, get_path=gp, path=filename)

nbands=nb

;set up and write the envi header for output image
bandname=strarr(1)
bandname[0]=’Mean area (‘+image1name+’)’
envi_setup_head, data_type=1, $
descrip=’Output image from GENERATE_MEANBURNAREA.PRO using 1 km Burn Percent and Burn Count maps’, $
file_type=0, fname=outimage+’.hdr’, interleave=0, nl=nl, ns=ns, nb=nbands, map_info=map_info, $
xstart=dims(1), ystart=dims(3), bnames=bandname, /WRITE
;nb = nb. Original nb=1

;stop
;open a metadata file that has the same name as the output image
;with a .meta at the end and start filli ng it
openw,u1,outimage+”, /GET_LUN
openw,u2,outimage+’.txt’, /GET_LUN
;BE CAREFUL WITH openw NOT TO BLAST A GOOD FILE

date_time = systime()
;print to file (= printf)
printf,u2,’Metadata file for BINARY_CLASS run ‘,date_time
printf,u2, ‘ ‘
printf,u2, ‘Input % 1km pixel burned image : ‘,image1name
printf,u2, ‘Input % 1km burn count image : ‘,image2name
printf,u2, ‘Input image pathname: ‘,fullpath
printf,u2, ‘Number of lines : ‘, nl
printf,u2, ‘Number of samples : ‘, ns
printf,u2, ‘Number of bands : ‘, nb
printf,u2, ‘Data type : ‘, data_type
printf,u2, ‘ ‘
printf,u2, ‘Output classified image : ‘,outimage
printf,u2, ‘Number of lines : ‘, nl
printf,u2, ‘Number of samples : ‘, ns
printf,u2, ‘Number of bands : ‘, 1
printf,u2, ‘Data type : byte’

;print to screen
print,’Metadata file for BINARY_CLASS run ‘,date_time
print, ‘ ‘
printf,u2, ‘Input % 1km pixel burned image : ‘,image1name
printf,u2, ‘Input % 1km burn count image : ‘,image2name
print, ‘Input image pathname: ‘,fullpath
print, ‘Number of lines : ‘, nl
print, ‘Number of samples : ‘, ns
print, ‘Number of bands : ‘, nb
print, ‘Data type : ‘, data_type

print, ‘ ‘
print, ‘Output classified image : ‘,outimage
print, ‘Number of lines : ‘, nl
print, ‘Number of samples : ‘, ns
print, ‘Number of bands : ‘, 1
print, ‘Data type : byte’

printf,u2, ‘ ‘
printf,u2, ‘ENVI Image bands used:
;printf,u2, FORMAT='(I2)’, envivars
printf,u2, ‘ ‘

print, ‘ ‘
print, ‘ENVI Image bands used:
;print, FORMAT='(I2)’, envivars
print, ‘ ‘

;make an array with the same dimensions as original image
inputimage1 = findgen(ns,nl,nbands)
inputimage2 = findgen(ns,nl,nbands)
S = SIZE(arr)

;help,inputimage
;retrieve image data for only those bands used in the tree
;this reduces the data volume read into memory to the bands that will actually be used
;for band=0,1 do begin
; ;****double check that this is getting the right bands****
inputimage1[*,*] = ENVI_GET_DATA(dims=dims,fid=fid,pos=pos) ;idlvars[i])
inputimage2[*,*] = ENVI_GET_DATA(dims=dims,fid=fid2,pos=pos2) ;idlvars[i])
; print,’importing envi image band number: ‘,band+1 ;envivars[i]
;endfor
help,inputimage1
help,inputimage2
;stop

;create an array to hold the float image (future output)
floatimage = findgen(ns,nl,nb)
;thresh20image = make_array(ns,nl,nb,/FLOAT)
;thresh30image = make_array(ns,nl,nb,/FLOAT)

print,’all parameters set, beginning calculation’

;get time that classification begins
starttime = systime(1)

;help,
help,floatimage
;help,inputimage

samp=0
line=0
for samp=0,ns-1 do begin
for line=0,nl-1 do begin
if inputimage2[samp,line] lt 0.9 then begin
floatimage[samp,line] = 0
endif else begin
floatimage[samp,line]=(inputimage1[samp,line]/inputimage2[samp,line])
endelse
endfor
endfor

help,floatimage

;fireimagesx = bindgen(ns,nl,nbands)
;fireimage20 = make_array(ns,nl,nb,/byte)
;fireimage30 = make_array(ns,nl,nb,/byte)
;help,fireimagesx

;for band=0,nb-1 do begin
; fireimagesx[*,*,0:nb-1] = bytscl(threshximage[*,*,0:nb-1],max=1.0000,min=0.9999,/NAN,TOP=1)
;endfor

;help,fireimagesx
;totalfiresx = bindgen(ns,nl)
;totalfiresx = TOTAL(fireimagesx,3,/NAN)

;bytscl

;****double check that this is getting the right bands****
;image[*,*,i] = ENVI_GET_DATA(dims=dims,fid=fid,pos=idlvars[i])
;print,’importing envi image band number: ‘,envivars[i]
;endfor

;create an array to hold the output mean image
;fltmeanarray = make_array(ns,nl, /float)
;fltmeanarray = AVERAGE(image)

print,’all calcs done set, beginning file writing’

;get time that classification begins
starttime = systime(1)

;begin a loop that processes the ENVI image pixel by pixel
;for i=0,nl-1 do begin

; show percent completed to screen
; percent = float(float(i)/float(nl)*100)
; print,FORMAT='($(1A))’,’percent completed: ‘
; print,FORMAT='(1F4.0)’,percent
;endfor

;write the array to the output image
print,’writing output image file’
;for
writeu,u1,floatimage
;calculate time that classification took
endtime = systime(1)
time = endtime – starttime
minutes = time/60
printf,u2,’Total time for classification: ‘,minutes, ‘ minutes’
print,’done, total time for classification: ‘,minutes, ‘ minutes’

close,u2,u1
end

; AVERAGE function sucked from help files & modified.
;FUNCTION AVERAGE, arr
; RETURN, TOTAL(arr)/N_ELEMENTS(arr)
;END

_________________

matzke
Geo Member

Leave a Reply