A great big script — generate_burnimages4.pro

A great big script — generate_burnimages4.pro

posted 06-04-2001 11:30 AM
Hi,
Here’s a great big script I used to generate a whole bunch of coarsened versions of my classified Landsat burnscar maps.

It’s not pretty, but it works. It took about 2.5 hours to run on a 6334×6700 floating-point Landsat image.

Thanks, Nick

code:

; PRO GENERATE_BURNIMAGES4
;
; by Nick Matzke, May 2001
;
; PURPOSE: To convert a high-resolution (30-m) resolution image to multiple coarse-
; resolution products (1-km, 3-km, and 5-km), for comparison with coarse-resolution
; DMSP (1-km, or resized to 3 and 5 km). The coarse-res output will, for example,
; represent “number of fire scars centered in each 1-km pixel with mapped area over’
; 500,000 square meters,” or “total area of fires above 500,000”, or “average area of
; fires above 500,000 square meters.”
;
; As you can see, if we’ve got 10 different thresholds, 3 different resolution outputs,
; and three types of outputs, that makes for 90 (10x3x3) output images — more than can be
; practically done by hand with ENVI. And if the initial classification has to be
; changed, then one would have to go back and do the whole thing over again, which would
; be a drag.
;
; TO REPLICATE (hopefully!) BY HAND:
; 1) Build mask, DN > threshold output to temp-mask.img file
; 2) Apply mask, output to temp-masked.img file
; 3) Via Band Math, take images temp-mask (for counts) and temp-masked (for area),
; and multiply each by the number of ETM pixels in coarse-resolution power. E.g.:
; a) for 1-km, “b1*1111.111111”
; b) for 3-km, “b1*10,000” (or, for (a), “b1 * 9”)
; c) for 5-km, “b1*27,777.777777” (or, for (a), “b1 * 25”)
; 4) ENVI Resize Data, for 30-m –> 1-km, set resize dimensions to 0.03, and resize
; with pixel aggregate resampling. This produces the average value of the input
; pixels in the output pixel, but since we multiplied by the number of pixels,
; this output represents the total.
; …This process isn’t difficult, but it is tedious to repeat by hand, and naming and
; keeping track of all the files, not making mistakes, etc., is difficult.
;
; This script can currently replicate these steps, but there are slight differences
; compared with the Pixel Aggregrate function. Whether these are mistakes/differences
; in the script or RSI’s aggregate function has yet to be determined.
;
; POSSIBLE FUTURE ADDITIONS: This will probably go in another script, but it would be
; worthwhile to have a program that compares all of these products to the co-extensive
; 1-km, 3-km, and 5-km products (e.g., extracts the points and performs a regression
; against each of them, and pops out the best combinations
;
; INPUTS: The input is derived thus:
; 1) Take a Landsat scene, process it (unmixing via reference endmembers in this case),
; and classify it via binary decision-tree (method of Dar Roberts).
; 2) Take the classified image, and ENVI export it to ArcView BIL.
; 3) Import to ArcView
; 4) Convert to Grid
; 5) Query to extract all points of the “fresh burn scar” class. Eliminate other
; points.
; 6) Convert burn scar points to polygon shapefile.
; 7) Run ArcView Avenue script to attribute each polygon with AREA and PERIMETER
; attributes.
; 8) Run another ArcView Avenue script to output the centroid of each polygon into a
; point shapefile, each centroid retaining the attributes of its polygon.
; 9) Convert centroid point shapefile to a grid with the same dimensions as the original
; image.
; 10) Save as floating-point binary image.
; 11) Import to ENVI. You now have an ENVI raster with a nonzero pixel at the centroid
; of every mapped burnscar. This can be inputted into this script and output
; a multitude of images comparable to the DMSP count images (or coarser-res area
; images)
;
; THE GENERAL PLAN FOR THIS SCRIPT:
; 1) IF necessary, subset input image to match maximum dimensions used in generating
; coarse-resolution images.
;
; 2) a) Input raster with pixels either containing burn area or zero (= no fire)
; b) Input fire threshold size (this will be an array of increasing thresholds).
;
; [BEGIN LOOP HERE, for every one of 10 (or whatever) thresholds]
;
; 3) Generate mask of centroid points image (output: temporary)
; Output image: If input pixel < threshold, output = 0; if >, output = 1.0000.
;
; 4) Apply mask to centroid burn area image, temporary output image only contains pixels
; gt set area threshold
;
; 5) Resample to 1 km, 3 km and 5 km images (SKIP THIS STEP as we’re using TOTAL, not
; AVERAGE.
; a) multiply by 1111.1111, 10000, and 27777.77777
; b) resize (PIXEL AGGREGATE) by factors 0.03, 0.01, and 0.006, respectively
; c) this generates images of total count of burnscars above area threshold
;
; 4) Use mask to generate 30-m image with only centroid pixels above threshold area
;
; 5) Resample masked image to 1 km via resampling function written in this script, and
; 3- and 5-km outputs via REBINning of 1-km output.
; a) this generates total area above the given threshold
;
; 6) Perform band math on coarse-resolution images
; a) Divide b1(total masked burn area) by b2 (count of burns over threshold)
; b) generates 1, 3, and 5 km images with average average area of masked burns
;
pro generate_burnimages4
; Program title (must match filename)
; Started May, 2001.
; Written by Nick Matzke

; use envi_check_save to load the save files with the registration and band math routines
ENVI_CHECK_SAVE, /registration
ENVI_CHECK_SAVE, /envi_math
ENVI_CHECK_SAVE, /utility
ENVI_CHECK_SAVE, /spectral

; OPEN INPUT HIGH-RESOLUTION 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 30-m pixelated burnscar centroid image’, fid=fid1, 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, fid1, data_type=data_type, nb=nb, interleave=interleave, sname=sname1, fname=fname1;, h_map=h_map
;handle_value, h_map, map_info
map_info = envi_get_map_info(fid=fid1)

;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 name (starting path will be same as input ascii rulesfile)
outimagename=dialog_pickfile(title=’Enter Path/row/version (e.g. p158r73v2) of Output Image (DO NOT add *.img)’, get_path=gp, path=filename)
; don’t put apostrophes inside IDL quotes or they will truncate the quote!!!)
nbands=nb

; CHECK TO MAKE SURE IMAGE IS CORRECT SIZE (ALREADY SUBSETTED)
; To minimize processing and issues, the high-resolution image should be subset to the
; minimum number of pixels needed to fill up the coarse-resolution output array.

;figure out max fine dimensions that can be fit into ouput coarse image
resfine=30 ; high-resolution, in meters
rescoarse=1000 ; low-resolution, in meters
finefract=(float(resfine)/rescoarse)
ns_coarse_max=(floor(ns*finefract))
ns_fine_max=CEIL(ns_coarse_max/(float(finefract)))
nl_coarse_max=(floor(nl*resfine/rescoarse))
nl_fine_max=CEIL(nl_coarse_max/(float(finefract)))

print,’ns_coarse_max = ‘,ns_coarse_max
print,’nl_coarse_max = ‘,nl_coarse_max
print,’ ‘
print,’ns_fine_max = ‘,ns_fine_max
print,’nl_fine_max = ‘,nl_fine_max
print,’ ‘

; Run test to see if subsetted:
print,’Running test to see if subsetting needed…’
if ( (ns ne ns_fine_max) OR (nl ne nl_fine_max) ) then begin ; ‘ne’ = not equal
print,’High-resolution image subsetting is necessary…’
print,’Actually, do this manually to save time (ENVI Tools –> Resize data).’
print,’The high-resolution image dimensions should be:’
print,’ Number of samples = ‘,ns_fine_max
print,’ Number of lines = ‘,nl_fine_max
print,’Ending program; re-run when ready.’

;print,’In-script subsetting enabled.’
;print,’Subsetting commencing…’
;subsetimage(0:ns_fine_max-1,0:nl_fine_max-1) = $
;envi_get_data(fid=fid1, dims=[-1,0,ns_fine_max-1,0,nl_fine_max-1], pos=0)
;print,’Subsetting done.’
;print,’Entering ENVI data…”
;ENVI_ENTER_DATA, subsetimage, descrip=’Max coarsenable subset of 30-m image’, $
;file_type=0, map_info=map_info, xstart=0, ystart=0, bnames=bandname, r_fid=fid2
;fname=outimagename+’-30m_subset.hdr’,
;print,’ENVI_ENTER_DATA done.’

stop
;ending program early (saves time if subsetting done manually)
endif else begin
print,’Test passed.’
print,’Entering ENVI data…’
inputimage = envi_get_data(fid=fid1, dims=[-1,0,ns_fine_max-1,0,nl_fine_max-1], pos=0)
ENVI_ENTER_DATA, inputimage, descrip=’Max coarsenable subset of 30-m image’, $
file_type=0, map_info=map_info, xstart=0, ystart=0, bnames=bandname, r_fid=fid2;
;ENVI_SETUP_HEAD,
;, fname=outimagename+’-30m_input.hdr’
print,’ENVI_ENTER_DATA done.’
endelse

;FROM HELP FILE
;envi_open_file, ‘/data/test.image’, r_fid=fid
;envi_file_query, fid, ns=ns, nl=nl

;FROM HALLIGAN’S RANGECALC_part2.pro
;select the images
;envi_file_query, sr1_fid, nl=nl, ns=ns, h_map=h_map, fname=fname1, sname=sname1
;map_info = envi_get_map_info(fid=first_fid)

; SET THE MINIMUM AREA THRESHOLD FOR BURN SCARS
; Create an array of thresholds, to produce a group of output images
numberthresholds=byte(11)
burnareathresh = findgen(numberthresholds)
burnareathresh[0] = 1000000
burnareathresh[1] = 500000
burnareathresh[2] = 250000
burnareathresh[3] = 100000
burnareathresh[4] = 50000
burnareathresh[5] = 25000
burnareathresh[6] = 10000
burnareathresh[7] = 5000
burnareathresh[8] = 2500
burnareathresh[9] = 1000
burnareathresh[10] = 50 ;Using 0 screws up count images

;burnareathresh_max = float(10^37) ; set in program code; not currently used
;maximum in floating point (single precision) = +/- 10^38

; WRITING METADATA FILES
; open a metadata file that has the same name as the output image
; with a *-metadata.txt at the end and start filling it
; One metadata file for the whole shebang.

openw,u1,outimagename+’-metadata.txt’, /GET_LUN

;BE CAREFUL WITH openw NOT TO BLAST A GOOD FILE

;get time/date that program run begins
starttime = systime(1)

;print to file (= printf)
printf,u1, ‘Metadata file for GENERATE_BURNIMAGES.PRO run’
printf,u1, ‘Time program run performed: ‘,SYSTIME()
printf,u1, ‘ ‘
printf,u1, ‘Input pixelated burn centroid (with area attribute) image: ‘,sname1
printf,u1, ‘ fname = ‘,fname1
printf,u1, ‘ ‘
printf,u1, ‘Number of lines : ‘, nl
printf,u1, ‘Number of samples : ‘, ns
printf,u1, ‘Number of bands : ‘, nb
printf,u1, ‘Data type : ‘, data_type
printf,u1, ‘ ‘
printf,u1, ‘OUTPUT IMAGES’
printf,u1, ‘Thresholds:’
for a=0,numberthresholds-1 do begin
printf,u1, ‘Area threshold Output image names’
printf,u1, burnareathresh[a],’ ‘,outimagename,’-gt’,burnareathresh[a],’.img’
endfor
printf,u1, ‘Output mask image name : ‘,outimagename,’-tempmask.img’
printf,u1, ‘Output masked image name : ‘,outimagename,’-tempmasked.img’
printf,u1, ‘Number of lines : ‘, nl
printf,u1, ‘Number of samples : ‘, ns
printf,u1, ‘Number of bands : ‘, 1
printf,u1, ‘Data type : byte’

;print to screen
print,’ ‘
print,’===========================================================================’
print,’Metadata file for GENERATE_BURNIMAGES.PRO run’
print, ‘Time program run performed: ‘,SYSTIME()
print, ‘ ‘
print, ‘Input pixelated burn centroid (with area attribute) image: ‘,sname1
print, ‘ Full path name = fname = ‘,fname1
print, ‘ ‘
print, ‘Number of lines : ‘, nl
print, ‘Number of samples : ‘, ns
print, ‘Number of bands : ‘, nb
print, ‘Data type : ‘, data_type
print, ‘ ‘
print, ‘ ‘
print, ‘OUTPUT IMAGES’
print, ‘Thresholds:’
for a=0,numberthresholds-1 do begin
print, ‘Area threshold Output image names’
print, burnareathresh[a],’ ‘,outimagename,’-gt’,burnareathresh[a],’.img’
endfor
print, ‘Output mask image name : ‘,outimagename,’-tempmask.img’
print, ‘Output masked image name : ‘,outimagename,’-tempmasked.img’
print, ‘Number of lines : ‘, nl
print, ‘Number of samples : ‘, ns
print, ‘Number of bands : ‘, 1
print, ‘Data type : byte’
print,’ ‘
print,’ACTUAL SUBSETTED IMAGE DIMS’
print,’ns_coarse_max = ‘,ns_coarse_max
print,’nl_coarse_max = ‘,nl_coarse_max
print,’ ‘
print,’ns_fine_max = ‘,ns_fine_max
print,’nl_fine_max = ‘,nl_fine_max
print,’ ‘
print,’===========================================================================’
print,’ ‘

;This shouldn’t be necessary if above test was passed
;ns=ns_fine_max
;nl=nl_fine_max

;FOR EACH THRESHOLD
for a=0,10 do begin

; fnum=1+a
threshstring=string(burnareathresh[a],FORMAT='(I7.7)’)

openw,u2,outimagename+’-ttl_area-th’+threshstring+’.img’, /GET_LUN
openw,u3,outimagename+’-ttl_count-th’+threshstring+’.img’, /GET_LUN
openw,u4,outimagename+’-avg_area-th’+threshstring+’.img’, /GET_LUN

; run1string=string(‘openw,u’+fnum+’,’+outimagename+’-thresh’+threshstring+’.img’)
;print,run1string
;un1=’openw,u’,fnum,’,’+outimagename+’-thresh’+threshstring+’.img’ ; /GET_LUN
;dashes = ‘print,FORMAT=”(1A,””,’ + string(nclasses) + ‘A10,””,A)”,dash,dasharray

,dash’
;=execute(run1)

;MASKING IMAGE (start function here)

;FROM HALLIGAN:
; make a mask from the slant range image using the nearslant and farslant ranges and the math_doit routine
; CREATE MASK IMAGE — PIXELS > burnareathresh go into image “outimagename+’-mask.img'”, where DN=1 means TRUE

;r_fid = file id of output file
;m_fid = file id of mask file (or perhaps other input as well)
;

;length=strlen(lengthtest)
;print,burnareathresh[a]
;print,lengthtest
;print,length
;stop

;zo=’string(burnareathresh[a],/PRINT,FORMAT = ‘(I1)’))’
;astring
;change back to fid2 if subsetting done:
maskexp='(b1 ge ‘+threshstring+’)’ ; and (b1 le burnareathresh_max)’ (left out ’cause nothing gets that high)
math_doit, fid=fid1, exp=maskexp, out_name=outimagename+’-mask.img’, $
dims=[-1,0,ns_fine_max-1,0,nl_fine_max-1], pos=0, r_fid=fid3

;MASK BURN AREAS IMAGE
; where previous mask = 1, transfer value to new image…output “outimagename+’-masked.img'”
envi_mask_apply_doit, fid=fid1, dims=[-1,0,ns_fine_max-1,0,nl_fine_max-1],$
m_fid=fid3, m_pos=0, out_name=outimagename+’-masked.img’, pos=0, r_fid=fid4, value=0
;gr = mask * envi_get_data(fid=gr1_fid, dims=dims, pos=0)

;MULTIPLY BURN AREAS BY 1111.1111 TO COUNTERACT EFFECTS OF AVERAGING
;mathexp='(b1 * 1111.1111)’ ; 1111.111 = Number of TM pixels 1 km pixel
;math_doit, fid=fid3, dims=dims, exp=mathexp, out_name=outimagename+’-maskedx1111.img’, pos=0, r_fid=fid4

;RESIZE FILE: 30 m –> 10 m
;Set the resize factors (opposite of ENVI, e.g. 30m–>1 km in ENVI = (0.03,0.03); here, = (33.33,33.33)
;rfact4=[0.33333333,0.33333333]

; Call the doit
;resize_doit, fid=fid5, pos=pos, out_name=outimagename+’-areapixels-10m.img’, dims=dims, rfact=rfact4, interp=1, r_fid=fid6
;interp=1 is, I think, nearest neighbor resampling, which will not modify values…

;Import array from fid6
;outimagename2 = ENVI_GET_DATA(dims=dims,fid=fid6,pos=pos)
;outimagename3 = findgen(ns_coarse_max,nl_coarse_max)
;outimagename3 = rebin(outimagename2,ns_coarse_max,nl_coarse_max,/SAMPLE)

;ENVI_ENTER_DATA, outimagename3, data_type=4, $
;descrip=’1-km avg burn area, thresh>100,000 — from GENERATE_BURNIMAGES using 30-m burn scar area input’, $
;file_type=0, fname=outimagename+’-area_gt_100000-1km.hdr’, interleave=0, nl=nl, ns=ns, nb=nbands, map_info=map_info, $
;xstart=dims(1), ystart=dims(3), bnames=bandname, r_fid=fid7 /WRITE
;nb = nb. Original nb=1

;rfact5=[1,1]
;resize_doit, fid=fid7, pos=pos, out_name=outimagename+’-area_gt_100000-1km.img’, dims=dims, rfact=rfact5, interp=1, r_fid=fid6

;Create array for coarse-res output image
lowresimage=findgen(ns_coarse_max,nl_coarse_max)
outcoarseimage=findgen(ns_coarse_max,nl_coarse_max)
outcoarsecountimage=findgen(ns_coarse_max,nl_coarse_max)

;Create array for fine-res subsetted, thresholded input image
highresareaimage = envi_get_data(fid=fid4, dims=[-1,0,ns_fine_max-1,0,nl_fine_max-1], pos=0)
threshcountimage = envi_get_data(fid=fid3, dims=[-1,0,ns_fine_max-1,0,nl_fine_max-1], pos=0)

;INVOKING RESAMPLE FUNCTION
asdf=bytarr(1)
;asdf =call_function(
;CALLING FUNCTION
asdf=call_function(‘RESAMP30_1000’, rescoarse,resfine,nl_coarse_max,ns_coarse_max,highresareaimage,$
outcoarseimage)

asdf=call_function(‘RESAMP30_1000′, rescoarse,resfine,nl_coarse_max,ns_coarse_max,threshcountimage,$
outcoarsecountimage)

;ENVI_ENTER_DATA, outcoarseimage, descrip=’output image’+threshstring+’-thresh’, $
;file_type=0, map_info=map_info, r_fid=fid8
;nl=nl_coarse_max, ns=ns_coarse_max, $ ;Standard stuff…
;nb=1, xstart=1, ystart=1, ;bnames=bandname,

;file_type=0, $ ;0 = ENVI Standard File
;map_info=map_info, $ ;Map_info retrieved from input image, at top of script

print,’1st ENVI_ENTER DATA DONE’

;envi_setup_head, fname=out_name, ns=230, nl=400, nb=3, $
;data_type=2, interleave=0, bnames=[‘Munsell Hue’, $
;’Munsell Sat’, ‘Munsell Val’], /write, /open, $
;descrip=’Munsell Transform’

;set up and write the envi header for output image

;bandname will hold the name of the band — threshold, etc.

burnareathreshstring=string(burnareathresh[10])

bandname=strarr(1)
description=strarr(1)
bandname[0]=’Mean area > ‘+burnareathreshstring+’ square meters’
description[0]=’Output ‘+burnareathreshstring+’ — (‘+outimagename+’.img)’

print,”+bandname[0]+”
print,”+description[0]+”

;B = {dstretch_struct, type:2, vals:fltarr(3)}

;begin writing header
; envi_setup_head, data_type=4,

;4 = floating point
;descrip=description[0], $ ;description is string w/ description of file
; file_type=0, $ ;0 = ENVI Standard File
; fname=fid8, $ ;Name of input file
; fname=outimagename+’.hdr’,$ ;Name of input file
; interleave=0, $ ;0 = Band Sequential format (BSQ)
; nl=nl_coarse_max, ns=ns_coarse_max, $ ;Standard stuff…
; nb=1, $
; map_info=map_info, $ ;Map_info retrieved from input image, at top of script
; xstart=1, ystart=1;, $ ;xstart & ystart not derived from high-res “dims”
;bnames=bandname[0], $ ;bandname is string with the name of the band
;/WRITE, $ ;Presumably this writes the header…
;/OPEN ;$ ;Puts the file in the ENVI list
;DEF_STRETCH=dstretch_struct
;nb = nb. Original nb=1

print,”
print,’ENVI_ENTER_DATA done.’
print,”
;fname=outimagename+’-30m_subset.hdr’,

envi_file_mng, id=fid3, /remove, /delete
envi_file_mng, id=fid4, /remove, /delete

;write the array to the output image
print,’writing output image file, threshold = ‘,burnareathresh[a]

;run1=’writeu,u’+fnum+’,’+outcoarseimage+”
;outimagename+’-thresh’+threshstring+’.img’ ; /GET_LUN
;r=execute(run1)

; writearray=outcoarseimage
writeu,u2,outcoarseimage
writeu,u3,outcoarsecountimage

;Create average image…

print,’Creating average image…’
outcoarseavgareaimage=findgen(ns_coarse_max,nl_coarse_max)
;Check for zeros in the denominator
print, ‘Checking for zeros in the denominator’
samp=0
line=0
for samp=0,ns_coarse_max-1 do begin
for line=0,nl_coarse_max-1 do begin
if outcoarsecountimage[samp,line] eq 0 then begin
outcoarseavgareaimage[samp,line] = 0
endif else begin
outcoarseavgareaimage[samp,line]=(outcoarseimage[samp,line]/outcoarsecountimage[samp,line])
endelse
endfor
endfor

; finitethresh=FINITE(outcoarseavgareaimage)
; outcoarseavgareaimage=finitethresh*outcoarseavgareaimage

;ENVI_ENTER_DATA, outcoarseavgareaimage, descrip=’avg_area image’+threshstring+’-thresh’, $
;file_type=0, map_info=map_info, r_fid=fid6, $
;nl=nl_coarse_max, ns=ns_coarse_max $ ;Standard stuff…
;nb=1, xstart=1, ystart=1, ;bnames=bandname,

;file_type=0, $ ;0 = ENVI Standard File
;map_info=map_info, $ ;Map_info retrieved from input image, at top of script

;DO I NEED TO REMOVE -NaN?

;maskexp='(b1 ne ‘+nanthresh+’)’ ; and (b1 le burnareathresh_max)’ (left out ’cause nothing gets that high)
;maskexp='(b1 ge ‘+threshstring+’)’
;math_doit, fid=fid6, exp=maskexp, dims=[-1,0,ns_coarse_max-1,0,nl_coarse_max-1], pos=0, r_fid=fid7

;threshcountimage = envi_get_data(fid=fid3, dims=[-1,0,ns_fine_max-1,0,nl_fine_max-1], pos=0)

;MASK BURN AREAS IMAGE
; where previous mask = 1, transfer value to new image…output “outimagename+’-masked.img'”
;envi_mask_apply_doit, fid=fid6, dims=[-1,0,ns_coarse_max-1,0,nl_coarse_max-1],$
;m_fid=fid7, m_pos=0, out_name=outimagename+’-avg_area-th’+threshstring+’.img’, pos=0, r_fid=fid8, value=0

;envi_file_mng, id=fid5, /remove, /delete
;envi_file_mng, id=fid6, /remove, /delete
;envi_file_mng, id=fid7, /remove, /delete
;envi_file_mng, id=fid8, /remove, /delete

writeu,u4,outcoarseavgareaimage

;RESAMPLE EVERYTHING TO 3 & 5 KM
;
; Work out ns & nl for coarser images…
; CHECK TO MAKE SURE IMAGE IS CORRECT SIZE (ALREADY SUBSETTED)
; To minimize processing and issues, the high-resolution image should be subset to the
; minimum number of pixels needed to fill up the coarse-resolution output array.

;figure out max fine dimensions that can be fit into ouput coarse image
res3km=3000 ; high-resolution, in meters
rescoarse=1000 ; low-resolution, in meters
km3fract=(float(rescoarse)/res3km)
ns_coarse3km_max=(floor(ns_coarse_max*km3fract))
ns_coarse1km_max=CEIL(ns_coarse3km_max/(float(km3fract)))
nl_coarse3km_max=(floor(nl_coarse_max*km3fract))
nl_coarse1km_max=CEIL(nl_coarse3km_max/(float(km3fract)))

;1 km->3km, AVERAGE AREA
outcoarseavgareaimage1kmmax=findgen(ns_coarse1km_max,nl_coarse1km_max)
;Subset by max 1km dims
outcoarseavgareaimage1kmmax[0:ns_coarse1km_max-1,0:nl_coarse1km_max-1] = $
outcoarseavgareaimage[0:ns_coarse1km_max-1,0:nl_coarse1km_max-1]

outcoarseavgareaimage3km=findgen(ns_coarse3km_max,nl_coarse3km_max)
outcoarseavgareaimage3km=REBIN(outcoarseavgareaimage1kmmax,ns_coarse3km_max,nl_coarse3km_max)
help,outcoarseavgareaimage3km
;Straight avg, no zero check…
openw,u5,outimagename+’-avg_area-th’+threshstring+’-3km.img’, /GET_LUN
writeu,u5,outcoarseavgareaimage3km
close,u5
FREE_LUN,u5

;1 km->3km, TOTAL AREA
outcoarseimage1kmmax=findgen(ns_coarse1km_max,nl_coarse1km_max)
;Subset by max 1km dims
outcoarseimage1kmmax[0:ns_coarse1km_max-1,0:nl_coarse1km_max-1] = $
outcoarseimage[0:ns_coarse1km_max-1,0:nl_coarse1km_max-1]
outcoarseimage1kmmax=outcoarseimage1kmmax*9
outcoarseimage3km=REBIN(outcoarseimage1kmmax,ns_coarse3km_max,nl_coarse3km_max)
;Straight avg, no zero check…
openw,u6,outimagename+’-ttl_area-th’+threshstring+’-3km.img’, /GET_LUN
writeu,u6,outcoarseimage3km
close,u6
FREE_LUN,u6

;1 km->3km, TOTAL COUNT
outcoarsecountimage1kmmax=findgen(ns_coarse1km_max,nl_coarse1km_max)
;Subset by max 1km dims
outcoarsecountimage1kmmax[0:ns_coarse1km_max-1,0:nl_coarse1km_max-1] = $
outcoarsecountimage[0:ns_coarse1km_max-1,0:nl_coarse1km_max-1]
outcoarsecountimage1kmmax=outcoarsecountimage1kmmax*9
outcoarsecountimage3km=REBIN(outcoarsecountimage1kmmax,ns_coarse3km_max,nl_coarse3km_max)
;Straight avg, no zero check…
openw,u7,outimagename+’-ttl_count-th’+threshstring+’-3km.img’, /GET_LUN
writeu,u7,outcoarsecountimage3km
close,u7
FREE_LUN,u7

;1KM->5KM
;figure out max fine dimensions that can be fit into ouput coarse image
res5km=5000 ; very-low-resolution, in meters
rescoarse=1000 ; low-resolution, in meters
km5fract=(float(rescoarse)/res5km)
ns_coarse5km_max=(floor(ns_coarse_max*km5fract))
ns_coarse1km_max=CEIL(ns_coarse5km_max/(float(km5fract)))
nl_coarse5km_max=(floor(nl_coarse_max*km5fract))
nl_coarse1km_max=CEIL(nl_coarse5km_max/(float(km5fract)))

;1 km->5km, AVERAGE AREA
outcoarseavgareaimage1kmmax=findgen(ns_coarse1km_max,nl_coarse1km_max)
;Subset by max 1km dims
outcoarseavgareaimage1kmmax[0:ns_coarse1km_max-1,0:nl_coarse1km_max-1] = $
outcoarseavgareaimage[0:ns_coarse1km_max-1,0:nl_coarse1km_max-1]

outcoarseavgareaimage5km=findgen(ns_coarse5km_max,nl_coarse5km_max)
outcoarseavgareaimage5km=REBIN(outcoarseavgareaimage1kmmax,ns_coarse5km_max,nl_coarse5km_max)
help,outcoarseavgareaimage5km
;Straight avg, no zero check…
openw,u8,outimagename+’-avg_area-th’+threshstring+’-5km.img’, /GET_LUN
writeu,u8,outcoarseavgareaimage5km
close,u8
FREE_LUN,u8

;1 km->5km, TOTAL AREA
outcoarseimage1kmmax=findgen(ns_coarse1km_max,nl_coarse1km_max)
;Subset by max 1km dims
outcoarseimage1kmmax[0:ns_coarse1km_max-1,0:nl_coarse1km_max-1] = $
outcoarseimage[0:ns_coarse1km_max-1,0:nl_coarse1km_max-1]
outcoarseimage1kmmax=outcoarseimage1kmmax*25
outcoarseimage5km=REBIN(outcoarseimage1kmmax,ns_coarse5km_max,nl_coarse5km_max)
;Straight avg, no zero check…
openw,u9,outimagename+’-ttl_area-th’+threshstring+’-5km.img’, /GET_LUN
writeu,u9,outcoarseimage5km
close,u9
FREE_LUN,u9

;1 km->5km, TOTAL COUNT
outcoarsecountimage1kmmax=findgen(ns_coarse1km_max,nl_coarse1km_max)
;Subset by max 1km dims
outcoarsecountimage1kmmax[0:ns_coarse1km_max-1,0:nl_coarse1km_max-1] = $
outcoarsecountimage[0:ns_coarse1km_max-1,0:nl_coarse1km_max-1]
outcoarsecountimage1kmmax=outcoarsecountimage1kmmax*25
outcoarsecountimage5km=REBIN(outcoarsecountimage1kmmax,ns_coarse5km_max,nl_coarse5km_max)
;Straight avg, no zero check…
openw,u10,outimagename+’-ttl_count-th’+threshstring+’-5km.img’, /GET_LUN
writeu,u10,outcoarsecountimage5km
close,u10
FREE_LUN,u10

;outimagename3=findgen(191,201)
;outimagename3=congrid(outimagename2,33.33333333,33.33333333,/INTERP)
;fid6=outimagename2
;LOAD outimagename2 with ENVI RESIZE_DOIT
;Set the resize factors (opposite of ENVI, e.g. 30m–>1 km in ENVI = (0.03,0.03); here, = (33.33,33.33)
;rfact6=[1,1]
; Call the doit
;resize_doit, fid=fid6, pos=pos, out_name=outimagename+’-ge100000-1kmb.img’, dims=dims, rfact=rfact6, r_fid=fid7

close,u2,u3,u4
;close,/ALL
FREE_LUN,u2,u3,u4,fid3,fid4;,fid5;,fid6,fid7,fid8;u4,

endfor

;stop

;END FOR EACH THRESHOLD
;endfor

;if ((k/coarsefract)-floor(i/coarsefract)) eq 0.99 then begin

; endif
; if ((k/coarsefract)-floor(i/coarsefract)) eq 0.99 then begin

; endif
; if ((k/coarsefract)-floor(i/coarsefract)) eq 0.99 then begin

; endif

;Now we have to multiply the mask and the masked image by 1111.1111 and 10,000 to make up for averaging in a moment…

;maskexp='(b1 ge ‘+string(burnareathresh)+’)’ ; and (b1 le burnareathresh_max)’ (left out ’cause nothing gets that high)
;math_doit, fid=fid1, dims=dims, exp=maskexp, out_name=outimagename+’-mask.img’, pos=0, r_fid=fid3

;make an array with the same dimensions as original image
;inputimage = 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,nb-1 do begin
; print,’Importing image…’
; inputimage = ENVI_GET_DATA(dims=dims,fid=fid,pos=pos) ;idlvars[i])
; print,’Image imported successfully’ ;envivars[i]
;endfor
;help,inputimage

;set up and write the envi header for output image
;bandname=strarr(1)
;bandname[0]=’Temp Threshold Mask image (‘+outimagename+’.img)’
;envi_setup_head, data_type=1, $
;descrip=’Output binary mask image from GENERATE_BURNIMAGES.PRO using threshold’, $
;file_type=1, fname=outimagename+’.hdr’, interleave=0, nl=nl, ns=ns, nb=1, map_info=map_info, $
;xstart=dims(1), ystart=dims(3), bnames=bandname, /WRITE
;nb = nb. Original nb=1

;threshximage= bindgen(ns,nl)

;i=float(0)
;j=float(0)

;print,’ ‘
;print,’Generating THRESHXIMAGE (temp file)’
;print,’Line number = ‘
;for i=0,nl-1 do begin
; for j=0,ns-1 do begin
;help,i
;print,i
;help,j
;print,j
; threshximage=float((inputimage)/burnareathresh) ;)[*,*,*] [(0:ns-1),(0:nl-1),band]
; endfor
; print,i
;endfor
;print,’Threshold high-resolution image completed’
;help,threshximage

;print,’ ‘
;print,’Generating outimagename:’
;for i=0,nl-1 do begin
; for j=0,ns-1 do begin
; outimagename = bytscl(threshximage,max=1.0000,min=0.9999,/NAN,TOP=1)
;
;ENVI_APPLY_MASK

; endfor
;endfor
;print,’outimagename complete’
;help,outimagename

;inputimage1 = findgen(ns,nl,nbands)

;make an array with the same dimensions as original image

;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’

;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’

;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

;calculate time that classification took
endtime = systime(1)
time = endtime – starttime
minutes = time/60
printf,u1,’Total time for classification: ‘,minutes, ‘ minutes’
print,’done, total time for classification: ‘,minutes, ‘ minutes’

close,u1;,u2
end

;START REBIN FUNCTION
;function

;START RESAMPLE FUNCTION
function RESAMP30_1000,rescoarse,resfine,nl_coarse_max,ns_coarse_max,highresimage,$
lowresimage,threshcountimage

;SUBSAMPLING NON-EASILY OVERLAYABLE GRIDS
;
;Go through 30-m image pixel-by-pixel
k=float(0) ; = samps for fine (=30m) image
l=float(0) ; = lines for fine (=30m) image
i=0 ; = samps for coarse (=1 km) image
j=0 ; = lines for coarse (=1 km) image
x=float(0) ; TOGGLE: if x=0, subtract 1 from k (= fine samp)
y=float(0) ; TOGGLE: if y=0, subtract 1 from l (= fine line)
samps_processed=float(0) ; Tally of # samps processed in a line
coarsefract=(float(float(rescoarse)/resfine)) ; Number of fine pixels in coarse pixel
print,’Beginning IF statements’
print,’Progress displayed by coarse and fine image line number…’

;FOR EACH LINE IN COARSE IMAGE:
for j=0,nl_coarse_max-1 do begin

;FOR EACH SAMPLE IN COARSE IMAGE:
for i=0,ns_coarse_max-1 do begin

;UPPER THIRD OF 1-km PIXEL
if ((round(100*((l/coarsefract)-floor(l/coarsefract))))) eq 0 then begin
if ((round(100* ((k/coarsefract)-floor(k/coarsefract)) ))) eq 0 then begin
;UPPER LEFT
;BTW, “l” = L, NOT 1 (one)…wasted an hour or two there!
;print,1 ;option #1 chosen
if (TOTAL(highresimage[k:k+33,l:l+33]) eq 0) then begin
lowresimage[i,j]=0
endif else begin
lowresimage[i,j]=( (TOTAL(highresimage[k:k+32,l:l+32]))$
+ (TOTAL(highresimage[k+33,l:l+32])/3) $ ;right edge
+ (TOTAL(highresimage[k:k+32,l+33])/3) $ ;bottom edge
+ (highresimage[k+33,l+33]/9) ) ;lower right corner
endelse
endif
if ((round(100* ((k/coarsefract)-floor(k/coarsefract))))) eq 99 then begin
;UPPER MIDDLE
;print,2 ;option #2 chosen
if (TOTAL(highresimage[k:k+33,l:l+33]) eq 0) then begin
lowresimage[i,j]=0
endif else begin
lowresimage[i,j]=( (TOTAL(highresimage[k+1:k+32,l:l+32]))$
+ (TOTAL(highresimage[k,l:l+32])*2/3) $ ;left edge
+ (TOTAL(highresimage[k+33,l:l+32])*2/3) $ ;right edge
+ (TOTAL(highresimage[k+1:k+32,l+33])/3) $ ;bottom edge
+ (highresimage[k,l+33]*2/9) + (highresimage[k+33,l+33]*2/9) ) ;bottom corners
endelse
endif
if ((round(100* ((k/coarsefract)-floor(k/coarsefract))))) eq 98 then begin
;UPPER RIGHT
;print,3 ;option #3 chosen
if (TOTAL(highresimage[k:k+33,l:l+33]) eq 0) then begin
lowresimage[i,j]=0
endif else begin
lowresimage[i,j]=( (TOTAL(highresimage[k+1:k+33,l:l+32]))$
+ (TOTAL(highresimage[k,l:l+32])/3) $ ;left edge
+ (TOTAL(highresimage[k+1:k+33,l+33])/3) $ ;bottom edge
+ (highresimage[k,l+33]/9) ) ;lower left corner
endelse
endif

endif

;MIDDLE THIRD OF 1-km PIXEL
if ((round(100*((l/coarsefract)-floor(l/coarsefract))))) eq 99 then begin
if ((round(100* ((k/coarsefract)-floor(k/coarsefract)) ))) eq 0 then begin
;MIDDLE LEFT
;print,4 ;option #4 chosen
if (TOTAL(highresimage[k:k+33,l:l+33]) eq 0) then begin
lowresimage[i,j]=0
endif else begin
lowresimage[i,j]=( (TOTAL(highresimage[k:k+32,l+1:l+32]))$
+ (TOTAL(highresimage[k+33,l+1:l+32])*1/3) $;right 1/3 edge
+ (TOTAL(highresimage[k:k+32,l+0])*2/3) $ ;top 2/3 edge
+ (TOTAL(highresimage[k:k+32,l+33])*2/3) $ ;bottom 2/3 edge
+ (highresimage[k+33,l]*2/9) $ ;top right corner
+ (highresimage[k+33,l+33]*2/9) ) ;bottom right corner
endelse
endif
if ((round(100* ((k/coarsefract)-floor(k/coarsefract))))) eq 99 then begin
;MIDDLE MIDDLE
;print,5 ;option #5 chosen
if (TOTAL(highresimage[k:k+33,l:l+33]) eq 0) then begin
lowresimage[i,j]=0
endif else begin
lowresimage[i,j]=( (TOTAL(highresimage[k+1:k+32,l+1:l+32]))$
+ (TOTAL(highresimage[k+0,l+1:l+32])*2/3) $ ;left edge
+ (TOTAL(highresimage[k+33,l+1:l+32])*2/3) $ ;right edge
+ (TOTAL(highresimage[k+1:k+32,l+0])*2/3) $ ;top edge
+ (TOTAL(highresimage[k+1:k+32,l+33])*2/3) $ ;bottom edge
+ ((highresimage[k,l]+highresimage[k+33,l] $;four corners
+ highresimage[k+33,l+33]+highresimage[k,l+33])*4/9));2/3*2/3=4/9
endelse
endif
if ((round(100* ((k/coarsefract)-floor(k/coarsefract))))) eq 98 then begin
;MIDDLE RIGHT
;print,6 ;option #6 chosen
if (TOTAL(highresimage[k:k+33,l:l+33]) eq 0) then begin
lowresimage[i,j]=0
endif else begin
lowresimage[i,j]=( (TOTAL(highresimage[k+1:k+33,l+1:l+32]))$
+ (TOTAL(highresimage[k,l+1:l+32])*1/3) $ ;left 1/3 edge
+ (TOTAL(highresimage[k+1:k+33,l+0])*2/3) $ ;top 2/3 edge
+ (TOTAL(highresimage[k+1:k+33,l+33])*2/3) $ ;bottom 2/3 edge
+ (highresimage[k,l]*2/9) $ ;top left corner
+ (highresimage[k,l+33]*2/9) ) ;bottom left corner

endelse
endif
endif

;LOWER THIRD OF 1-km PIXEL
if ((round(100*((l/coarsefract)-floor(l/coarsefract))))) eq 98 then begin
if ((round(100* ((k/coarsefract)-floor(k/coarsefract)) ))) eq 0 then begin
;LOWER LEFT
;print,7 ;option #7 chosen
if (TOTAL(highresimage[k:k+33,l:l+33]) eq 0) then begin
lowresimage[i,j]=0
endif else begin
lowresimage[i,j]=( (TOTAL(highresimage[k:k+32,l+1:l+33]))$
+ (TOTAL(highresimage[k+33,l+1:l+33])/3) $ ;right edge
+ (TOTAL(highresimage[k:k+32,l+0])/3) $ ;top edge
+ (highresimage[k+33,l]/9) ) ;top right corner
endelse
endif
if ((round(100* ((k/coarsefract)-floor(k/coarsefract))))) eq 99 then begin
;LOWER MIDDLE
;print,8 ;option #8 chosen
if (TOTAL(highresimage[k:k+33,l:l+33]) eq 0) then begin
lowresimage[i,j]=0
endif else begin
lowresimage[i,j]=( (TOTAL(highresimage[k+1:k+32,l+1:l+33]))$
+ (TOTAL(highresimage[k+0,l+1:l+33])*2/3) $ ;left edge
+ (TOTAL(highresimage[k+33,l+1:l+33])*2/3) $ ;right edge
+ (TOTAL(highresimage[k+1:k+32,l+0])/3) $ ;top edge
+ (highresimage[k,l]*2/9) + (highresimage[k+33,l]*2/9) );top corners
endelse
endif
if ((round(100* ((k/coarsefract)-floor(k/coarsefract))))) eq 98 then begin
;LOWER RIGHT
;print,9 ;option #9 chosen
if (TOTAL(highresimage[k:k+33,l:l+33]) eq 0) then begin
lowresimage[i,j]=0
endif else begin
lowresimage[i,j]=( (TOTAL(highresimage[k+1:k+33,l+1:l+33]))$
+ (TOTAL(highresimage[k,l+1:l+33])/3) $ ;left edge
+ (TOTAL(highresimage[k+1:k+33,l])/3) $ ;top edge
+ (highresimage[k,l]/9) ) ;top left corner
endelse
endif
endif

;COUNT SAMPLES
samps_processed = samps_processed+1

if ( (float(i+1)/3)-(floor((i+1)/3)) ) eq 0 then begin ;1st +34 at 3rd i (i=2)
; print, (float(i+1)/3)
; print, (floor((i+1)/3))
k=k+34
endif else begin
k=k+33
endelse

;if x eq 0 then begin
; k=k-1
; x=1
;endif

;END SAMPLES (coarse) LOOP
endfor

k=0
x=0
if ( (float(j+1)/3)-(floor((j+1)/3)) ) eq 0 then begin
l=l+34
endif else begin
l=l+33
endelse

;if y eq 0 then begin
; l=l-1
; y=1
;endif

;END LINES (coarse) LOOP
endfor

print,’coarse line # =’,j,’ fine line # =’,l,’ +y =’,l+y,’ coarse samp# =’,i,’ fine samp# =’,k,’ +x=’,k+x
print,’resampling complete’
;print
;print,’ACCURACIES (%)’
;print,’—————————————————————————————-‘
;print,”Class Producer’s User’s”
;print,’—————————————————————————————-‘
;for row=0,nclasses-1 do begin
;print,alphanamearray[row],” “,producers[row],” “,users[row]
; print,FORMAT='(A,” “,F6.2,” “,F6.2)’,alphanamearray[row],producers[row],users[row]
;endfor
print,’—————————————————————————————-‘

print
print
print,”Congradulations! You’re finished. Have a nice day.”
print
print

return,1
end

; $Id: import_create_varname.pro,v 1.2 2000/07/14 16:37:05 chris Exp $
;
; Copyright (c) 2000, Research Systems, Inc. All rights reserved.
; Unauthorized reproduction prohibited.
;+
; NAME:
; IMPORT_CREATE_VARNAME
;
; PURPOSE:
; This routine takes a string file name and constructs a valid
; variable name. For use with the IMPORT_ routines.
;
; CATEGORY:
; Input/Output
;
; CALLING SEQUENCE:
;
; IMPORT_CREATE_VARNAME, Name [, Path [, Suffix]]
;
; OUTPUTS:
; A string with a valid variable name.
;
; ARGUMENTS:
;
; Name = A string containing the file name, possibly including a file path.
;
; Path = A string containing the file path to be removed from Name.
; The default is ” (Null string).
;
; Suffix = A string containing a suffix to be appended onto the variable
; name. The default is ” (Null string).
;
; MODIFICATION HISTORY:
; Written by: CT, RSI, July, 2000
;-

FUNCTION import_create_varname, file_nameIn, file_path, suffix

; COMPILE_OPT hidden, strictarr

ON_ERROR, 2 ; return to caller

; strip off file_path if necessary
IF (SIZE(file_path,/TNAME) NE ‘STRING’) THEN file_path = ”
varName = (file_name = STRMID(file_nameIn, STRLEN(file_path)))

; strip off filetype suffix if necessary
period = STRPOS(varName, ‘.’, /REVERSE_SEARCH) ; is there a filetype?
IF (period GT 0) THEN $ ; strip off
varName = STRMID(varName, 0, period) $
ELSE $ ; or, for a “dot” file, strip off the period
IF (period EQ 0) THEN varName = STRMID(varName, 1)

; remove illegal variable name characters
; first character must be a A-Z,a-z letter
firstLetter = STREGEX(varName, ‘[a-z]+’, /FOLD_CASE)
varName = (firstLetter LT 0) ? ‘var’+varName : STRMID(varName, firstLetter)

IF (STRLEN(varName) GT 0) THEN BEGIN
varName = STRMID(varName,LINDGEN(STRLEN(varName)),1) ; split into chars

; Replace any spaces with underscores
spaces = WHERE(varName EQ ‘ ‘,nspace)
IF (nspace GT 0) THEN varName[spaces] = ‘_’

; Remove illegal characters
legalChars = ‘[a-z_$0-9]’
good = WHERE(STRMATCH(varName, legalChars, /FOLD_CASE), nmatch)
varName = (nmatch EQ 0) ? ” : STRJOIN(varName[good],/SINGLE)
ENDIF

; If null name, throw an error…
IF (STRLEN(varName) EQ 0) THEN MESSAGE, /NONAME, $
‘Unable to construct variable name from file name “‘ + file_name + ‘”‘

; append the suffix (note that the suffix does not get checked!)
IF (SIZE(suffix,/TNAME) EQ ‘STRING’) THEN varName = varName + suffix

RETURN, varName
END

Nick

______________

matzke
Geo Member

Leave a Reply