www.pudn.com > batch_reproject_resize_v03.rar > batch_reproject_resize_v03.pro
; batch_reproject_resize v 0.3 and batch_reproject_resize_reverse v 0.3
; Author: Jonathan A. Greenberg (greenberg@ucdavis.edu)
; batch_reproject_resize is designed to batch reproject and crop (or expand) a set
; of images to a single base image's projection and size. This is
; helpful with preparing large sets of data for the purposes of change detection.
; batch_reproject_resize_reverse is designed to reproject and crop (or expand) a single
; image based on a set of input base images. This can be helpful if, for instance,
; you want a set of DEMs associated with a large number of flightline images, all pulled
; from the same source DEM. batch_reproject_resize_reverse requires batch_reproject_resize.
; Dependencies: COMPARE_STUCT (included),
; from http://idlastro.gsfc.nasa.gov/ftp/pro/structure/compare_struct.pro
; Known issues: I haven't done a lot of error checking yet.
; - if an image has an arbitrary projection it will throw an error.
; - hemispheres may cause problems (can someone check this?)
; - no output directory parameter yet, all output files go where the input files are.
; - possibly do a header check and set background value to the data ignore value.
; Usage:
; Open in ENVI/IDL, go to run-> Compile batch_reproject_resize_vXX.pro or
; place this in your save_add directory. At the ENVI command line:
;
; batch_reproject_resize[,INFILES=string array,BASEFILE=string,INSUFFIX=string,OUTSUFFIX=string,$
; WARP_METHOD={0 | 1 | 2 | 3},RESAMPLING={0 | 1 | 2},background=double,outbasename=string array
;
; INFILES (optional): a string array containing the full path to each of the input files. If
; no INFILES are given, a GUI window to select the files will be opened. All infiles must
; have a proper projection set up in the header.
; BASEFILE (optional): a string containing the full path to the base file that all of the INFILES
; will be warped and cropped (or expanded) to. If no BASEFILE is given, a GUI window to select
; the file will be opened. The basefile must have a proper projection set up in the header.
; INSUFFIX (optional): a string containing the input files' suffix (default is '.img'). The code
; uses this to set the output filename by replacing the INSUFFIX portion of the filename with OUTSUFFIX.
; OUTSUFFIX (optional): a string containing the output files' suffix (default is '_reproject.img'). The code
; uses this to set the output filename by replacing the INSUFFIX portion of the filename with OUTSUFFIX.
; WARP_METHOD (optional): keyword to set warp method as described in ENVI_CONVERT_FILE_MAP_PROJECTION. We
; HIGHLY recommend leaving this set to its default (3, Rigorous (pixel-by-pixel)), as other settings
; can cause poor results.
; RESAMPLING (optional): keyword to set resampling method as described in ENVI_CONVERT_FILE_MAP_PROJECTION.
; Default is 0 (nearest neighbor). Other settings can be 1: bilinear or 2: cubic convolution.
; BACKGROUND (optional): value to set the background values to during the reproject and any expansion
; of the input files. The default is 0.
; OUTBASENAME (optional): don't use this, its there for compatibility with batch_reproject_resize_reverse.
; batch_reproject_resize_reverse,[,INFILES=string,BASEFILE=string,INSUFFIX=string,OUTSUFFIX=string,$
; WARP_METHOD={0 | 1 | 2 | 3},RESAMPLING={0 | 1 | 2},background=double
; INFILES (optional): a string array containing the full path to each of the input files. If
; no INFILES are given, a GUI window to select the files will be opened. All infiles must
; have a proper projection set up in the header.
; BASEFILE (optional): a string containing the full path to the base file that will be warped and
; cropped/expanded to each of the INFILES. If no BASEFILE is given, a GUI window to select the
; file will be opened. The basefile must have a proper projection set up in the header.
; INSUFFIX (optional): a string containing the input files' suffix (default is '.img'). The code
; uses this to set the output filename by replacing the INSUFFIX portion of the filename with OUTSUFFIX.
; OUTSUFFIX (optional): a string containing the output files' suffix (default is '_reproject.img'). The code
; uses this to set the output filename by replacing the INSUFFIX portion of the filename with OUTSUFFIX.
; WARP_METHOD (optional): keyword to set warp method as described in ENVI_CONVERT_FILE_MAP_PROJECTION. We
; HIGHLY recommend leaving this set to its default (3, Rigorous (pixel-by-pixel)), as other settings
; can cause poor results.
; RESAMPLING (optional): keyword to set resampling method as described in ENVI_CONVERT_FILE_MAP_PROJECTION.
; Default is 0 (nearest neighbor). Other settings can be 1: bilinear or 2: cubic convolution.
; BACKGROUND (optional): value to set the background values to during the reproject and any expansion
; of the input files. The default is 0.
; batch_reproject_resize and batch_reproject_resize_reverse is Copyright (C) 2008 Jonathan A. Greenberg
; This program is free software: you can redistribute it and/or modify
; it under the terms of the GNU General Public License as published by
; the Free Software Foundation, either version 3 of the License, or
; (at your option) any later version.
; This program is distributed in the hope that it will be useful,
; but WITHOUT ANY WARRANTY; without even the implied warranty of
; MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
; GNU General Public License for more details.
; You should have received a copy of the GNU General Public License
; along with this program. If not, see .
;+
; NAME:
; COMPARE_STRUCT
; PURPOSE:
; Compare all matching tag names and return differences
;
; EXPLANATION:
; Compare all matching Tags names (except for "except_Tags")
; between two structure arrays (may have different struct.definitions),
; and return a structured List of fields found different.
;
; CALLING SEQUENCE:
; diff_List = compare_struct( struct_A, struct_B [ EXCEPT=, /BRIEF,
; /FULL, /NaN, /RECUR_A, /RECUR_B )
; INPUTS:
; struct_A, struct_B : the two structure arrays to compare.
; Struct_Name : for internal recursion use only.
; OPTIONAL INPUT KEYWORDS:
; EXCEPT = string array of Tag names to ignore (NOT to compare).
; /BRIEF = number of differences found for each matching field
; of two structures is printed.
; /FULL = option to print even if zero differences found.
; /NaN = if set, then tag values are considered equal if they
; are both set to NaN
; /RECUR_A = option to search for Tag names
; in sub-structures of struct_A,
; and then call compare_struct recursively
; for those nested sub-structures.
; /RECUR_B = search for sub-structures of struct_B,
; and then call compare_struct recursively
; for those nested sub-structures.
; Note:
; compare_struct is automatically called recursively
; for those nested sub-structures in both struct_A and struct_B
; (otherwise cannot take difference)
; OUTPUT:
; Returns a structure array describing differences found,
; which can be examined using print,diff_List or help,/st,diff_List.
; PROCEDURE:
; Match Tag names and then use where function on tags.
; EXAMPLE:
; Find the tags in the !X system variable which are changed after a
; simple plot.
; IDL> x = !X ;Save original values
; IDL> plot, indgen(25) ;Make a simple plot
; IDL> help,/str,compare_struct(x,!X) ;See how structure has changed
;
; and one will see that the tags !X.crange and !X.S are changed
; by the plot.
; MODIFICATION HISTORY:
; written 1990 Frank Varosi STX @ NASA/GSFC (using copy_struct)
; modif Aug.90 by F.V. to check and compare same # of elements only.
; Converted to IDL V5.0 W. Landsman September 1997
; Added /NaN keyword W. Landsman March 2004
;-
function compare_struct, struct_A, struct_B, EXCEPT=except_Tags, Struct_Name, $
FULL=full, BRIEF=brief, NaN = NaN, $
RECUR_A = recur_A, RECUR_B = recur_B
common compare_struct, defined
if N_params() LT 2 then begin
print,'Syntax - diff_List = compare_struct(struct_A, struct_B '
print,' [EXCEPT=, /BRIEF, /FULL, /NaN, /RECUR_A, /RECUR_B ]'
if N_elements(diff_List) GT 0 then return, diff_List else return, -1
endif
if N_elements( defined ) NE 1 then begin
diff_List = { DIFF_LIST, Tag_Num_A:0, Tag_Num_B:0, $
Field:"", Ndiff:0L }
defined = N_tags( diff_List )
endif else diff_List = replicate( {DIFF_LIST}, 1 )
Ntag_A = N_tags( struct_A )
if (Ntag_A LE 0) then begin
message," 1st argument must be a structure variable",/CONTIN
return,diff_List
endif
Ntag_B = N_tags( struct_B )
if (Ntag_B LE 0) then begin
message," 2nd argument must be a structure variable",/CONTIN
return,diff_List
endif
N_A = N_elements( struct_A )
N_B = N_elements( struct_B )
if (N_A LT N_B) then begin
message,"comparing "+strtrim(N_A,2)+" of first structure",/CON
message,"to first "+strtrim(N_A,2)+" of "+strtrim(N_B,2)+ $
" in second structure",/CONTIN
diff_List = compare_struct( struct_A, struct_B[0:N_A-1], $
EXCEPT=except_Tags, $
RECUR_A = recur_A, $
RECUR_B = recur_B, $
FULL=full, BRIEF=brief )
return,diff_List
endif else if (N_A GT N_B) then begin
message,"comparing first "+strtrim(N_B,2)+" of "+ $
strtrim(N_A,2)+" in first structure",/CON
message,"to "+strtrim(N_B,2)+" in second structure",/CONTIN
diff_List = compare_struct( struct_A[0:N_B-1], struct_B, $
EXCEPT=except_Tags, $
RECUR_A = recur_A, $
RECUR_B = recur_B, $
FULL=full, BRIEF=brief )
return,diff_List
endif
Tags_A = tag_names( struct_A )
Tags_B = tag_names( struct_B )
wB = indgen( N_elements( Tags_B ) )
Nextag = N_elements( except_Tags )
if (Nextag GT 0) then begin
except_Tags = [strupcase( except_Tags )]
for t=0,Nextag-1 do begin
w = where( Tags_B NE except_Tags[t], Ntag_B )
Tags_B = Tags_B[w]
wB = wB[w]
endfor
endif
if N_elements( struct_name ) NE 1 then sname = "." $
else sname = struct_name + "."
for t = 0, Ntag_B-1 do begin
wA = where( Tags_A EQ Tags_B[t] , nf )
if (nf GT 0) then begin
tA = wA[0]
tB = wB[t]
NtA = N_tags( struct_A.(tA) )
NtB = N_tags( struct_B.(tB) )
if (NtA GT 0 ) AND (NtB GT 0) then begin
if keyword_set( full ) OR keyword_set( brief ) then $
print, sname + Tags_A[tA], " :"
diffs = compare_struct( struct_A.(tA), struct_B.(tB), $
sname + Tags_A[tA], $
EXCEPT=except_Tags, $
FULL=full, BRIEF=brief )
diff_List = [ diff_List, diffs ]
endif else if (NtA LE 0) AND (NtB LE 0) then begin
if keyword_set(NaN) then begin
x1 = struct_b.(tB)
x2 = struct_a.(tA)
g = where( finite(x1) or finite(x2), Ndiff )
if Ndiff GT 0 then $
w = where( x1[g] NE x2[g], Ndiff )
endif else $
w = where( struct_B.(tB) NE struct_A.(tA) , Ndiff )
if (Ndiff GT 0) then begin
diff = replicate( {DIFF_LIST}, 1 )
diff.Tag_Num_A = tA
diff.Tag_Num_B = tB
diff.Field = sname + Tags_A[tA]
diff.Ndiff = Ndiff
diff_List = [ diff_List, diff ]
endif
if keyword_set( full ) OR $
(keyword_set( brief ) AND (Ndiff GT 0)) then $
print, Tags_A[tA], Ndiff, FORM="(15X,A15,I9)"
endif else print, Tags_A[ta], " not compared"
endif
endfor
if keyword_set( recur_A ) then begin
for tA = 0, Ntag_A-1 do begin
if N_tags( struct_A.(tA) ) GT 0 then begin
diffs = compare_struct( struct_A.(tA), struct_B, $
sname + Tags_A[tA], $
EXCEPT=except_Tags, $
RECUR_A = recur_A, $
RECUR_B = recur_B, $
FULL=full, BRIEF=brief )
diff_List = [ diff_List, diffs ]
endif
endfor
endif
if keyword_set( recur_B ) then begin
for tB = 0, Ntag_B-1 do begin
if N_tags( struct_B.(tB) ) GT 0 then begin
diffs = compare_struct( struct_A, struct_B.(tB), $
sname + Tags_B[tB], $
EXCEPT=except_Tags, $
RECUR_A = recur_A, $
RECUR_B = recur_B, $
FULL=full, BRIEF=brief )
diff_List = [ diff_List, diffs ]
endif
endfor
endif
w = where( [diff_List.Ndiff] GT 0, np )
if (np LE 0) then w = [0]
return, diff_List[w]
end
pro batch_reproject_resize,infiles=raster,basefile=basefilename,insuffix=insuffix,outsuffix=outsuffix,$
warp_method=warp_method,resampling=resampling,background=background,outbasename=outbasename
if n_elements(raster) eq 0 then begin
rasterfilenames=envi_pickfile(title='Please select input raster(s):',/multiple_files)
if rasterfilenames[0] eq '' then return
endif else rasterfilenames=raster
numfiles=n_elements(rasterfilenames)
rasterfids=lonarr(numfiles)
for i=0,numfiles-1 do begin
envi_open_file, rasterfilenames[i], r_fid=tempfid
rasterfids[i]=tempfid
endfor
if n_elements(basefilename) eq 0 then begin
basefilename=envi_pickfile(title='Please select base image to crop to...')
if basefilename eq '' then return
endif
if n_elements(insuffix) eq 0 then insuffix='.img'
if n_elements(outsuffix) eq 0 then outsuffix='_reproject.img'
if insuffix eq outsuffix then begin
print,'insuffix must be different from outsuffix, exiting.'
return
end
tempsuffix='_temp.img'
if n_elements(warp_method) eq 0 then warp_method=3
if n_elements(resampling) eq 0 then resampling=0
if n_elements(background) eq 0 then background=0
rasterfilebasenames=strarr(numfiles)
rasterfiletempnames=strarr(numfiles)
if n_elements(outbasename) lt numfiles then begin
for i=0,numfiles-1 do begin
rasterfilebasenames[i]=$
file_dirname(rasterfilenames[i],/mark_directory)+file_basename(rasterfilenames[i],insuffix)+outsuffix
rasterfiletempnames[i]=$
file_dirname(rasterfilenames[i],/mark_directory)+file_basename(rasterfilenames[i],insuffix)+tempsuffix
endfor
endif else begin
for i=0,numfiles - 1 do begin
rasterfilebasenames[i]=$
outbasename+outsuffix
rasterfiletempnames[i]=$
outbasename+tempsuffix
endfor
endelse
envi_open_file,basefilename,r_fid=basefid
envi_file_query,basefid,ns=basens,nl=basenl,nb=basenb,xstart=xstart,ystart=ystart,dims=basedims
o_proj = envi_get_projection(fid=basefid)
map_info = envi_get_map_info(fid=basefid)
o_pixel_size=map_info.ps[0:1]
inherit = envi_set_inheritance(basefid, basedims,[0],/full)
hemispherical_correction=[1,-1]
base_corners=fltarr(2,4)
; Right now we assume N. hemisphere -- this should be changed
base_corners[*,0]=map_info.mc[2:3]-(map_info.mc[0:1]*o_pixel_size)*hemispherical_correction
base_corners[*,1]=base_corners[*,0]+(o_pixel_size*[basens,0])*hemispherical_correction
base_corners[*,2]=base_corners[*,0]+(o_pixel_size*[0,basenl])*hemispherical_correction
base_corners[*,3]=base_corners[*,0]+(o_pixel_size*[basens,basenl])*hemispherical_correction
for i=0,numfiles-1 do begin
; First, we reproject the entire image.
; If the image is already in the right projection, skip it...
tempfid=rasterfids[i]
tempoutname=rasterfiletempnames[i]
envi_file_query,tempfid,ns=ns,nl=nl,nb=nb,data_ignore_value=data_ignore_value,data_type=data_type
temp_proj=envi_get_projection(fid=tempfid)
if n_elements(compare_struct(temp_proj,o_proj)) gt 1 then begin
dims = [-1l, 0, ns-1, 0, nl-1]
pos=lindgen(nb)
if n_elements(data_ignore_value) eq 0 then data_ignore_value=background
envi_convert_file_map_projection, fid=tempfid, $
pos=pos, dims=dims, o_proj=o_proj, $
o_pixel_size=o_pixel_size, $
out_name=tempoutname, warp_method=warp_method, $
resampling=resampling, background=data_ignore_value, $
r_fid=tempoutfid
endif else begin
tempoutfid=tempfid
endelse
; Now we crop the image to the base image.
envi_file_query,tempoutfid,ns=outns,nl=outnl,nb=outnb
tempout_map_info=envi_get_map_info(fid=tempoutfid)
temp_o_proj = envi_get_projection(fid=tempoutfid)
tempout_corners=fltarr(2,4)
; Right now we assume N. hemisphere -- this should be changed
tempout_corners[*,0]=tempout_map_info.mc[2:3]-(tempout_map_info.mc[0:1]*o_pixel_size)*hemispherical_correction
tempout_corners[*,1]=tempout_corners[*,0]+(o_pixel_size*[outns,0])*hemispherical_correction
tempout_corners[*,2]=tempout_corners[*,0]+(o_pixel_size*[0,outnl])*hemispherical_correction
tempout_corners[*,3]=tempout_corners[*,0]+(o_pixel_size*[outns,outnl])*hemispherical_correction
if tempout_corners[0,0] ge base_corners[0,0] then begin
tempout_xs = 0
tempout_xspos = long(abs(tempout_corners[0,0] - base_corners[0,0])/o_pixel_size[0])+1
endif else begin
tempout_xs = long(abs(tempout_corners[0,0] - base_corners[0,0])/o_pixel_size[0])
tempout_xspos = 0
endelse
if tempout_corners[0,1] lt base_corners[0,1] then begin
tempout_xe = outns-1
endif else begin
tempout_xe = outns-long(abs(tempout_corners[0,1] - base_corners[0,1])/o_pixel_size[0])-2
endelse
tempout_xepos=tempout_xspos+(tempout_xe-tempout_xs)
if tempout_corners[1,0] lt base_corners[1,0] then begin
tempout_ys = 0
tempout_yspos = long(abs(tempout_corners[1,0] - base_corners[1,0])/o_pixel_size[1])
endif else begin
tempout_ys = long(abs(tempout_corners[1,0] - base_corners[1,0])/o_pixel_size[1])
tempout_yspos = 0
endelse
if tempout_corners[1,3] ge base_corners[1,3] then begin
tempout_ye = outnl-1
tempout_yepos = long(abs(tempout_corners[1,3] - base_corners[1,3])/o_pixel_size[1])-1
endif else begin
tempout_ye = outnl-long(abs(tempout_corners[1,3] - base_corners[1,3])/o_pixel_size[1])-1
tempout_yepos = 0
endelse
; Write out file
temppos=lindgen(outnb)
tempoutname=rasterfilebasenames[i]
openw,tempoutunit,tempoutname,/get_lun
tempoutarr=make_array([basens,outnb],type=data_type,value=background)
; fltarr(basens,nb)+float(background)
for j=0,tempout_yspos-1 do writeu,tempoutunit,tempoutarr
for j=tempout_ys,tempout_ye do begin
tempoutarr[tempout_xspos:tempout_xepos,*]=$
envi_get_slice(fid=tempoutfid, line=j, pos=temppos, xs=tempout_xs, xe=tempout_xe)
writeu,tempoutunit,tempoutarr
endfor
tempoutarr=fltarr(basens,nb)+float(background)
for j=0,tempout_yepos-1 do writeu,tempoutunit,tempoutarr
close,tempoutunit
free_lun,tempoutunit
envi_setup_head,fname=tempoutname,nb=nb,/open,/write,ns=basens,nl=basenl,$
data_type=data_type,interleave=1,inherit=inherit,data_ignore_value=data_ignore_value
endfor
end
pro batch_reproject_resize_reverse,infiles=raster,basefile=basefile,insuffix=insuffix,outsuffix=outsuffix,$
warp_method=warp_method,resampling=resampling,background=background
if n_elements(raster) eq 0 then begin
rasterfilenames=envi_pickfile(title='Please select base raster(s) to crop to...:',/multiple_files)
if rasterfilenames[0] eq '' then return
endif else rasterfilenames=raster
numfiles=n_elements(rasterfilenames)
rasterfids=lonarr(numfiles)
; for i=0,numfiles-1 do begin
; envi_open_file, rasterfilenames[i], r_fid=tempfid
; rasterfids[i]=tempfid
; endfor
if n_elements(basefile) eq 0 then begin
basefilename=envi_pickfile(title='Please select input raster...')
if basefilename eq '' then return
endif
if n_elements(outdirectory) eq 0 then begin
outdirectory=dialog_pickfile(title='Please select output directory...',/directory)
if outdirectory eq '' then return
endif
if n_elements(insuffix) eq 0 then insuffix='.img'
if n_elements(outsuffix) eq 0 then outsuffix='_reproject.img'
tempsuffix='_temp.img'
if n_elements(warp_method) eq 0 then warp_method=3
if n_elements(resampling) eq 0 then resampling=0
if n_elements(background) eq 0 then background=0
for i=0,numfiles-1 do begin
tempbasename=outdirectory+file_basename(rasterfilenames[i],insuffix)
batch_reproject_resize,infiles=basefilename,basefile=rasterfilenames[i],insuffix=insuffix,outsuffix='_crop.img',$
warp_method=warp_method,resampling=resampling,background=background,outbasename=tempbasename
endfor
end