different attributes in a raster map

[Roger Bivand]

Francois:

On Mon, 7 Aug 2000, Francois Delclaux wrote:

> Hi, grass users,
> 
> May be someone has a solution for my problem :
> 
> I have a base raster map ( 475 watersheds)
> I have a cover raster map (dem)
> I would like to get an output  raster map giving, for each area of base map
> (watershed),
> the number of different category of cover map. (elevation)
> 
> I know that r.statistics is able to generate an output map including statistics
> of cover map
> for each area of  base map. In fact, I would like to get a similar result for
> the number
> of different category values.
> 

These kinds of things have been done using the R/GRASS interface, using
the R statistics environment to apply a chosen function (here a frequency
count) to a vector (here the DEM) by a category variable (a factor in R's
terminology, here the watersheds). You can find some details on:

http://www.geog.uni-hannover.de/grass/statsgrasslist.html

and a later draft paper on:

ftp://reclus.nhh.no/pub/gc00/gc009_e.tar.gz,

including the latest interface package. A paper on the interface is
forthcoming: Bivand, R. S., (2000) Using the R statistical data analysis
language on GRASS 5.0 GIS database files. Computers and Geosciences,
26, pp. 1043-1052. 

A practical question: how many cells are there in the DEM? 

As an example, below the variable under analysis is the spearfish
elevation.dem at its original resolution, classified by the fields
layer. First they are imported, then the counts of the numbers of cells by
field are reported, and finally tapply() used to apply the
summary() function to the elevation.dem data classified by field (in your
case it would be the table() function, and classification by watershed):

GRASS:~ > R --vsize=20M
> library(GRASS)
Running in /usr/local/grass-5.0b/data/spearfish/rsb 
> G <- gmeta()
> summary(G)
Data from GRASS 5.0 LOCATION spearfish with 633 columns and 466 rows;
UTM, zone: 13
The west-east range is: 590010, 609000, and the south-north: 4914020, 
4928000;
West-east cell sizes are 30 units, and south-north 30 units.
> spear <- rast.get(G, c("elevation.dem", "fields"), c(F,T))
> table(spear$fields)

                 no data               C. Smith#1               C. Smith#2 
                       0                      756                      700 
             P. Biggam#1              P. Biggam#2              P. Biggam#4 
                     350                     3379                      729 
             P. Biggam#5              P. Biggam#6              P. Biggam#7 
                     299                     2906                      621 
...


         B. Cartwright#1          B. Cartwright#2          B. Cartwright#3 
                    3123                      328                      312 
Black Hills Natl. Forest 
                  121379 
> tapply(spear$elevation.dem, spear$fields, summary)
$"no data"
NULL

$"C. Smith#1"
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
   1162    1175    1179    1182    1185    1224      14 

$"C. Smith#2"
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
   1135    1154    1161    1164    1172    1210 

...........

$"B. Cartwright#1"
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
   1199    1236    1260    1260    1282    1319 

$"B. Cartwright#2"
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
   1195    1208    1214    1214    1220    1243 

$"B. Cartwright#3"
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
   1193    1202    1207    1209    1216    1225 

$"Black Hills Natl. Forest"
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
   1158    1408    1494    1505    1600    1840      23 

>

(here using the summary() function rather than tabulation)

Hope this helps,

Roger

-- 
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Breiviksveien 40, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 93 93
e-mail: Roger.Bivand at nhh.no
and: Department of Geography and Regional Development, University of
Gdansk, al. Mar. J. Pilsudskiego 46, PL-81 378 Gdynia, Poland.