I have UTM x,y coordinate contour files representing probability
boundaries for animal home range utilization distributions. I can add
elevation data from digitized USGS quads. I need to be able to calculate
accurate topographic areas within these contours in GRASS. Any
information helpful to this project would be most appreciated.
Jim Taulman
Jim,
Here's a simplistic and rather laborious way of calculating the surface
area of 3-dimensional polygons. If anyone knows of an easier way of doing
this, I would appreciate hearing about it. First use your elevation map to
calculate out a slope map using percent slope as the output. Then,
assuming your elevation map cells are square instead of rectangular, the
3-d area of the cell will be run * sqrt(run^2 + rise^2), where run = the
length of elevation cell and rise = the slope value * run.
/@\
/@@@@@\
/@@@@@@@@@\
/@@@@@@@@@@@@/|
\@@@@@@@@@@@/ |
r \@@@@@@@/ | rise
u \@@@/ | (run * %_slope)
n /________|
run
In r.mapcalc the formula would be:
area = x * sqrt(x^2 + (x * slope)^2); substituting cell length for x
Next, determine which cells lie in the interior of your polygon and those
that are on the periphery. There are a number of ways to do this. I
typically import my polygon boundaries as line features, run v.to.rast,
reclass the output map assigning 0 to the polygon boundary and 1 to
everything else, run r.clump, then reclass the output map assigning 1 to
0, 2 to the interior clump, and 0 to everything else. This results in a
map where the periphery cells are coded 1 and the interior cells are
coded 2.
The final area of your 3-d polygon will be approximately the sum of the
areas of the interior cells plus 1/2 the total of the area of the
periphery cells. Since this methodology works remarkably well on flat
surfaces, I've made the assumption that it works equally well in the 3rd
dimension. I usually output the cell area values to do the final summation
in either a spreadsheet or database since the final area can be quite
large. To do so I run 2 commands:
r.stats -1zg in=polygon out=polygon.tmp
r.what area < polygon.tmp > polygon.3d
This outputs an ascii file with one row for each cell consisting of its
easting, northing, polygon code, and area value.
Hope this helps.
Don Ebert
Research Associate
National Biological Service
Cooperative Research Unit
University of Nevada - Las Vegas
On Mon, 1 Apr 1996, James Taulman wrote:
I have UTM x,y coordinate contour files representing probability
boundaries for animal home range utilization distributions. I can add
elevation data from digitized USGS quads. I need to be able to calculate
accurate topographic areas within these contours in GRASS. Any
information helpful to this project would be most appreciated.
Jim Taulman
This is really useful - I have been looking for a way of doing this for
quite some time. Just wanted one clarification however - why do you
divide the area occupied by the boundary cells by half?
------------------------------------------------------------------------
| Harini Nagendra E-MAIL : harini@ces.iisc.ernet.in |
| Center for Ecological PHONE : 91 80 309 2639 (Hostel: LR-94)|
| Sciences : 91 80 309 2506 (Department) |
| Indian Institute of Science FAX : 91 80 334 1683 |
| Bangalore 560 012 India TELEX : 845 8349 IISc IN |
------------------------------------------------------------------------
On Tue, 2 Apr 1996, DON EBERT wrote:
Jim,
Here's a simplistic and rather laborious way of calculating the surface
area of 3-dimensional polygons. If anyone knows of an easier way of doing
this, I would appreciate hearing about it. First use your elevation map to
calculate out a slope map using percent slope as the output. Then,
assuming your elevation map cells are square instead of rectangular, the
3-d area of the cell will be run * sqrt(run^2 + rise^2), where run = the
length of elevation cell and rise = the slope value * run.
/@\
/@@@@@\
/@@@@@@@@@\
/@@@@@@@@@@@@/|
\@@@@@@@@@@@/ |
r \@@@@@@@/ | rise
u \@@@/ | (run * %_slope)
n /________|
run
In r.mapcalc the formula would be:area = x * sqrt(x^2 + (x * slope)^2); substituting cell length for x
Next, determine which cells lie in the interior of your polygon and those
that are on the periphery. There are a number of ways to do this. I
typically import my polygon boundaries as line features, run v.to.rast,
reclass the output map assigning 0 to the polygon boundary and 1 to
everything else, run r.clump, then reclass the output map assigning 1 to
0, 2 to the interior clump, and 0 to everything else. This results in a
map where the periphery cells are coded 1 and the interior cells are
coded 2.The final area of your 3-d polygon will be approximately the sum of the
areas of the interior cells plus 1/2 the total of the area of the
periphery cells. Since this methodology works remarkably well on flat
surfaces, I've made the assumption that it works equally well in the 3rd
dimension. I usually output the cell area values to do the final summation
in either a spreadsheet or database since the final area can be quite
large. To do so I run 2 commands:r.stats -1zg in=polygon out=polygon.tmp
r.what area < polygon.tmp > polygon.3dThis outputs an ascii file with one row for each cell consisting of its
easting, northing, polygon code, and area value.Hope this helps.
Don Ebert
Research Associate
National Biological Service
Cooperative Research Unit
University of Nevada - Las VegasOn Mon, 1 Apr 1996, James Taulman wrote:
> I have UTM x,y coordinate contour files representing probability
> boundaries for animal home range utilization distributions. I can add
> elevation data from digitized USGS quads. I need to be able to calculate
> accurate topographic areas within these contours in GRASS. Any
> information helpful to this project would be most appreciated.
> Jim Taulman
>
If you overlay the vector boundary line on the raster polygon, you'll see
that the proportion of the boundary cells that are on the inside and
outside of the vector boundary varies from cell to cell. By dividing the
area of the boundary cells by half, you are saying that, on average, half
of the cell is inside the vector boundary and half is on the outside.
Another note about area calculations. Beware of small cell sizes. When
GRASS truncates r.mapcalc's area calculations to integer values, you can
lose a lot of area over a large polygon. You might consider exporting out
the percent slope values for each cell instead of its area value and do
the area calculations yourself using a spreadsheet or database.
Don Ebert
Research Associate
National Biological Service
Cooperative Research Unit
University of Nevada - Las Vegas
On Wed, 3 Apr 1996, Harini Nagendra wrote:
This is really useful - I have been looking for a way of doing this for
quite some time. Just wanted one clarification however - why do you
divide the area occupied by the boundary cells by half?------------------------------------------------------------------------
| Harini Nagendra E-MAIL : harini@ces.iisc.ernet.in |
| Center for Ecological PHONE : 91 80 309 2639 (Hostel: LR-94)|
| Sciences : 91 80 309 2506 (Department) |
| Indian Institute of Science FAX : 91 80 334 1683 |
| Bangalore 560 012 India TELEX : 845 8349 IISc IN |
------------------------------------------------------------------------