resample_grid

 resample high resolution grid to lower target resolution
 Usage: resample_grid(sourcename,targetname,asciiout)
        where [sourcename] is the filename (incl path) of the source ascii
        grid to resample, [targetname] is the name of the lower resolution
        target grid to snap the resampled grid to, and [asciiout] is the
        filename to write the resampled output grid to.

0001 function resample_grid(sourcename,targetname,asciiout)
0002 % resample high resolution grid to lower target resolution
0003 % Usage: resample_grid(sourcename,targetname,asciiout)
0004 %        where [sourcename] is the filename (incl path) of the source ascii
0005 %        grid to resample, [targetname] is the name of the lower resolution
0006 %        target grid to snap the resampled grid to, and [asciiout] is the
0007 %        filename to write the resampled output grid to.
0008 
0009 %% Preallocate memory for large arrays
0010 % % note: this backfires when the allocated memory is not big enough!
0011 % disp('Preallocatinng memory for input DEMs')
0012 % [rows, cols]=get_gridsize([sourcename]);
0013 % srtm=zeros(cols, rows);
0014 % [rows, cols]=get_gridsize([targetname]);
0015 % globe=zeros(cols,rows);
0016 
0017 
0018 %% read in source
0019 disp('Reading in Source')
0020 [sxcoords,sycoords,srtm,sheader]=asciigrid2mat([sourcename]);
0021 disp('Done...')
0022 %% read in globe
0023 disp('Reading in Target')
0024 [gxcoords,gycoords,globe,gheader]=asciigrid2mat([targetname]);
0025 disp('Done...')
0026 
0027 
0028 %% resample high resolution grid to lower target resolution
0029 
0030 disp('Resampling Source to Target resolution using mean of cells and snapping to Target grid.')
0031 newsrtm=zeros(size(gycoords,2),size(gxcoords,2));
0032 for rows = 1:size(globe,1);
0033     for cols = 1:size(globe,2);
0034         newsrtm(rows,cols) = resample_mean(sheader.cellsize,gheader.cellsize, gxcoords(cols), gycoords(rows),min(sxcoords),min(sycoords), srtm);
0035     end;
0036 end; 
0037 disp('Done...')
0038 %clear srtm;
0039 
0040 %% crop
0041 % determine larger raster and crop NaNs
0042 % Note: if GLOBE originally had a larger extent then SRTM, border cells
0043 % are filled with NaNs. If SRTM was larger, everything is ok
0044 if (max(gxcoords)-min(gxcoords))*(max(gycoords)-min(gycoords))  >  (max(sxcoords)-min(sxcoords))*(max(sycoords)-min(sycoords))
0045     disp('Cropping NaN borders on grids.')
0046     [gxcoords,gycoords,globe,newsrtm] = crop_compare_nan2(gxcoords,gycoords,globe,newsrtm);  
0047 end
0048 %% write out
0049 mat2grid([asciiout],gxcoords,gycoords,newsrtm)