Using map algebra in GRASS

Using Map Algebra and the command line in GRASS

Introduction

GRASS is geographic information system (GIS) software that is used for geospatial data management and analysis. If you find yourself working on raster data specifically in GRASS it may be very useful to know how to use map algebra in order to manipulate the data so that you can properly analyze it. The advantages to understanding map algebra in GRASS is that it allows for manipulation of raster layers in a program that is focused on data processing and analysis. Also it is important to know how to use this tool using the command line instead of the graphical user interface (GUI) and so to that end for this whole page all of the examples will be done using the command line commands.

Map algebra

What is Map Algebra?

Map algebra is a way in GIS to perform arithmetic on raster layers in maps. The way that it works is that it will take two raster layers and it performs an operation on them based on what you tell it to do. The basics to know before starting map algebra is to know that map algebra is like a language and in order to communicate the statement to GRASS, the program follows a rules to effectively read the command given. Using r.mapcalc, the expression can either be typed out in the command line or use the preset buttons on top to assemble the command. One way access the r.mapcalc is through going to the Raster menu, then to Raster map calculator, and then the dropdown menu will have the option for Raster map calculator (r.mapcalc). For example, if you wanted to cut a raster image down to the study area of another you could do that by using the following command:

r.mapcalc --overwrite rasterAnew=”if(rasterB, rasterA)” 

this command basically tells the program that for the new raster that you are creating, if raster B(the study area) is in that location then put raster A(the raster you are cutting) there, meaning that the end result should be whatever you are measuring in raster A in just the area of raster B.

What can you do with Map Algebra?

There is a lot that you actually can do with map algebra; it is not just limited to changing the shape of raster files, it can also be used to shift raster files geographically, or to edit the values of each pixel of raster files. You could use it, for example, to subtract or add raster layers to Landsat images. This is one way to create colour photos from Landsat imagery - if you add the red green and blue layers together you get true colour photos. You can also use this technique to make different indexes from satellite images which can be used to tell things like plant life, from transforming the infrared layer and the red layer to create an NDVI index.

this table shows all of the different operations that can be done on raster files using r.mapcalc

The method for using these operators works like this:

1. Start with the command which in this case is: r.mapcalc
2. Then you want to add the file you are creating: newFile=
3. Then you add in your operator and variables for example: !fileA which means not file A

r.mapcalc newFile= !fileA

this will give you a raster that is the inverse of file A

The other way to do mapcalc is without the command line using the map calculator built into the GUI, which is accessed using the button shown below.

location of map calculator

map calculator

Examples

Subtracting raster to fit into study area

For this first example the DEMs (Digital Elevation Model) that are used can be obtained from the Government of Canada website. However, this example will work with any raster file in the same situation. To start, I have two DEMs one that is smaller than the other as you can see in the picture.

In order to subtract the excess from the large DEM to make it the same size as the as the smaller DEM which for this example will be our study area I used the command

"r.mapcalc --overwrite expression=DEM = DEM@elmin if( DEM@elmin, DEM__2_@elmin)" (elmin is the mapset name just happens that is my username on my computer).

the result as you can see in this picture is that both rasters fit in the study area.

Creating NDVI index using Landsat data

For this example I have obtained some data from USGS which is free to get you just need to create account with them. We will be calculating NDVI which stands for normalized difference vegetation index from Landsat images. This index is used to help tell the health of plants from satellite images and as you can see from the image below you need to start by uploading your Landsat imagery into GRASS. For this particular index you need bands four and five I have bands one through five here. You can use this method to create any index but for this example I have chosen NDVI because it is very common.

This is what the GRASS work station will look like before doing the mapcalc

Next, you have to open the command line so that you can enter the following command:

NDVI = float(band5 - band4)/float(band5 + band4)

In the example you can see that I forgot to add the float for each calculation the first time I tried it. It is important to add that because without that your result won't be correct, this is because the actual number for each pixel in an index is not an integer it is usually a decimal number so making sure the type of number is a float is important.

My command line

If you were to use the mapcalc calculator in the GUI it would like this:

then after you run the command it will look like this:

final result of the map calc

That is how you make an index using r.mapcalc.

Conclusion

r.mapcalc is a very useful skill to know in GRASS as GRASS is used to do a lot of analysis of raster data and in order to do that you need to know how to manipulate raster data. I made this tutorial because I had trouble with GRASS and hopefully this will help someone else have not as much trouble.

References

“Mapcalc Manual.” GRASS GIS Manual: R.mapcalc, GRASS Development Team, grass.osgeo.org/grass76/manuals/r.mapcalc.html.