Figure 5.1. Moving window method for neighborhood operations in raster map algebra. The raster cell value X of the new map is calculated from a 3x3 matrix of the old map

The resulting map newmapl will be stored as integer, while newmap2 will be stored as floating point. To transform an integer map into a floating point map, simply multiply it by 1.0 or use the floatO or doubleO functions:

inapcalc> newmap = 1.0 * old integer map mapcalc> newmap = float(old integer map)

The calculation of the normalized difference vegetation index (NDVI from LANDSAT-TM5) is a good example of an application where the function of integer maps needs to be stored as a floating point map:

mapcalo ndvi = 1.0 * {tm4 - tm3) / (tm4 + tmS) mapcalo ndvi = float ( (tm4 - tm3) / (tm4 + tm3) )

There is an alternative NDVI algorithm which uses a different function: mapcalo ndvi2 = 255.0 / 90.0 atan((tm4 - tm3) / (tm4 + tm3))

The maps ndvi and ndvi2 are the new floating point raster output maps, tm3 and tm4 are LANDSAT-TM5 channels used as integer input maps. Without the multiplication by 1.0, the result would be saved as integer and important information would be lost.

Examples of basic calculations. Cell-wise addition of two or more raster maps is one of the common map algebra tasks. For example, we can add the height of the buildings stored in a map buildings to an elevation model stored in a map dem:

r.mapcalc dem with buildings = buildings + dem

Or we can calculate a weighted average of two maps (here decimal points are used to ensure that the resulting newmap is stored in floating point format):

Handling of NULL values in r.mapcalc. The basic rule to remember when working with NULL data in map algebra is that operations on NULL cells lead to NULL cells. For example, if one of the maps included in the r.mapcalc expression has NULLs in the given area, the resulting map will have NULLs in this area too (both for addition and multiplication functions). In this way, NULL behaves differently from zero, which will have in this area zeroes for multiplication but not necessarily for addition. Therefore, if we want to do operations with NULL data, we need to use a special function isnull (). For example, if we want to fill the NULLs in mapl with values from map2 (in other words, when cell value in mapl is NULL, then write corresponding value of map2, otherwise use value in mapl) we run:

r.mapcalc new=if(isnull(mapl), map2, mapl)

We can also apply a function of map2 to the replaced NULL values (in our example we just add a constant 1000.0) as follows:

r.mapcalc new=if(isnull(mapl), map2 + 1000., mapl)

If you dont use the isnull () function, the NULL values will remain in the output map.

In another example we want to add the maps mapl and map2, where the map2 contains NULLs. To get a new map with NULLs replaced by 0 (zero) you have to enter:

r.mapcalc new=mapl + if(isnull(map2), 0, map2)

The above examples show that it is important to carefully evaluate the use of the function isnull() when applying map algebra to raster maps containing

NULLs.

Working with if conditions. Various logical operations can be performed with raster data by combining the operators with the if () functions. For example, we can create a new raster newmap by applying the if () function to an existing raster map and a set of other raster maps or values a, b, c:

if map = a then b else c is coded:

newmap = if((map == a),b,c)

if map is not equal a then b else c is coded:

newmap = if((map != a),b,c)

if map >= a then b else c is coded:

newmap = if((map >= a),b,c)

newmap = if ( (map==l),1,0)

Conditions with NULL values:

if((map==2) , 1,0)

- Change NULL values into a new value (if map=NULL then 3 else map):

newmap = if(isnull(map),3,map)

- Change all cell values smaller than 5 into NULL value, all other values to 5 (if map < 5 then NULL else 5):

newmap = if(map<5,null (),5)

Neighborhood operations with relative coordinates. The usage of the relative coordinates of the moving window is another useful option provided by r.mapcalc (see Figure 5.1). Neighboring cells can be used in calculations, and larger, possibly asymmetrical moving windows beyond the common 3x3 matrix can be defined. An example:

newmap = oldmap[1,1] + oldmap[-1,-1]

will read only the offset cells bottom right [1,1) and top left f-1,-1] in an 3x3 matrix for calculating the new map. This option is applicable to different input maps as well. Current row and column values can be integrated into expressions using row () and col (), the current coordinates of the moving window with x () and y ().

Note that the offset format is map c], where r is the row (y) offset and c is the column (x) offset. For example, map [1,2] refers to the cell one row below and two columns to the right of the current cell, map [-2, -1] refers to the cell two rows above and one column to the left of the current cell, and map [0,1] refers to the cell one column to the right of the current cell.

if map > = a and map <= b then c else d is coded:

newmap = if(((map >= a) && (map <= b)),c,d)

The if() functions can be combined to define more complex logical operations:

Select the values 1 and 2 from map and save them in newmap while setting all other values to 0:

newmap = if((map==l),1,0) + if((map==2),2,0)

Select the values 1 and 2 from map and save them in newmap as a binary map (only the values 0 and 1, representing yes and no ):