PURPOSE: sample value of a raster at any coordinates sx/sy
-------------------------------------------------------------------------
USAGE: value = sample_center(sx,sy,xcoords,ycoords,raster);
where: [sy],[sx] are the sample coordinates
[xcoords],[ycoords] are coordinate vectors of grid to be sampled
[grid] is the raster to be sampled
-------------------------------------------------------------------------
OUTPUTS:
[value] is the returned value at the specified coords
-------------------------------------------------------------------------
NOTES: xcoords(1) and ycoords(1) should have the llc coords of raster.
if target coord lies exactely on the border of two cells,
lower left cell is sampled.
If coordinates are outside the raster domain, NaN is returned.
Felix Hebeler, Geography Dept., University Zurich, March 2006.0001 function value = sample_center(sx,sy,xcoords,ycoords,raster) 0002 % PURPOSE: sample value of a raster at any coordinates sx/sy 0003 % ------------------------------------------------------------------------- 0004 % USAGE: value = sample_center(sx,sy,xcoords,ycoords,raster); 0005 % where: [sy],[sx] are the sample coordinates 0006 % [xcoords],[ycoords] are coordinate vectors of grid to be sampled 0007 % [grid] is the raster to be sampled 0008 % ------------------------------------------------------------------------- 0009 % OUTPUTS: 0010 % [value] is the returned value at the specified coords 0011 % ------------------------------------------------------------------------- 0012 % NOTES: xcoords(1) and ycoords(1) should have the llc coords of raster. 0013 % if target coord lies exactely on the border of two cells, 0014 % lower left cell is sampled. 0015 % If coordinates are outside the raster domain, NaN is returned. 0016 % 0017 % Felix Hebeler, Geography Dept., University Zurich, March 2006. 0018 0019 % calculate source cellsize 0020 xcs= abs(xcoords(2)-xcoords(1)); ycs= abs(ycoords(2)-ycoords(1)); 0021 0022 % get distance on x and y axis 0023 dx=abs(xcoords(1)- sx); 0024 dy=abs(ycoords(1)- sy); 0025 if sx>max(xcoords) || sy>max(ycoords) || sx<min(xcoords) || sy<min(ycoords) 0026 value = NaN; 0027 return 0028 end 0029 value = raster(size(raster,1)-(floor(dy/ycs)+1)+1, (floor(dx/xcs)+1)); 0030