1r.mapcalc(1) Grass User's Manual r.mapcalc(1)
2
6 r.mapcalc
7
9 r.mapcalc performs arithmetic on raster map layers. New raster map
10 layers can be created which are arithmetic expressions involving exist‐
11 ing raster map layers, integer or floating point constants, and func‐
12 tions.
13
15 If used without command line arguments, r.mapcalc will read its input,
16 one line at a time, from standard input (which is the keyboard, unless
17 redirected from a file or across a pipe). Otherwise, the expression on
18 the command line is evaluated. r.mapcalc expects its input to have the
19 form:
20 result=expression
22 where result is the name of a raster map layer to contain the result of
24 the calculation and expression is any legal arithmetic expression
25 involving existing raster map layers, integer or floating point con‐
26 stants, and functions known to the calculator. Parentheses are allowed
27 in the expression and may be nested to any depth. result will be cre‐
28 ated in the user's current mapset.
29 The formula entered to r.mapcalc by the user is recorded both in the
31 result map title (which appears in the category file for result) and in
32 the history file for result.
33 Some characters have special meaning to the command shell. If the user
35 is entering input to r.mapcalc on the command line, expressions should
36 be enclosed within single quotes. See NOTES, below.
37
39 The following operators are supported:
40 Operator Meaning Type Precedence
41 --------------------------------------------------------------
42 - negation Arithmetic 12
43 ~ one's complement Bitwise 12
44 ! not Logical 12
45 ^ exponentiation Arithmetic 11
46 % modulus Arithmetic 10
47 / division Arithmetic 10
48 * multiplication Arithmetic 10
49 + addition Arithmetic 9
50 - subtraction Arithmetic 9
51 << left shift Bitwise 8
52 >> right shift Bitwise 8
53 >>> right shift (unsigned) Bitwise 8
54 > greater than Logical 7
55 >= greater than or equal Logical 7
56 < less than Logical 7
57 <= less than or equal Logical 7
58 == equal Logical 6
59 != not equal Logical 6
60 & bitwise and Bitwise 5
61 | bitwise or Bitwise 4
62 && logical and Logical 3
63 &&& logical and[1] Logical 3
64 || logical or Logical 2
65 ||| logical or[1] Logical 2
66 ?: conditional Logical 1
67 (modulus is the remainder upon division)
68 [1] The &&& and ||| operators handle null values differently to
70 other operators. See the section entitled NULL support below for more
71 details.
72 The operators are applied from left to right, with those of higher
74 precedence applied before those with lower precedence. Division by 0
75 and modulus by 0 are acceptable and give a NULL result. The logical
76 operators give a 1 result if the comparison is true, 0 otherwise.
77
79 Anything in the expression which is not a number, operator, or function
80 name is taken to be a raster map layer name. Examples:
81 elevation
84 x3
85 3d.his
86 Most GRASS raster map layers meet this naming convention. However, if
89 a raster map layer has a name which conflicts with the above rule, it
90 should be quoted. For example, the expression
91 x = a-b
94 would be interpreted as: x equals a minus b, whereas
97 x = "a-b"
100 would be interpreted as: x equals the raster map layer named a-b
103 Also
105 x = 3107
108 would create x filled with the number 3107, while
111 x = "3107"
114 would copy the raster map layer 3107 to the raster map layer x.
117 Quotes are not required unless the raster map layer names look like
119 numbers or contain operators, OR unless the program is run non-interac‐
120 tively. Examples given here assume the program is run interactively.
121 See NOTES, below.
122 r.mapcalc will look for the raster map layers according to the user's
124 current mapset search path. It is possible to override the search path
125 and specify the mapset from which to select the raster map layer. This
126 is done by specifying the raster map layer name in the form:
127 name@mapset
130 For example, the following is a legal expression:
133 result = x@PERMANENT / y@SOILS
136 The mapset specified does not have to be in the mapset search path.
139 (This method of overriding the mapset search path is common to all
140 GRASS commands, not just r.mapcalc.)
141
143 Maps and images are data base files stored in raster format, i.e., two-
144 dimensional matrices of integer values. In r.mapcalc, maps may be fol‐
145 lowed by a neighborhood modifier that specifies a relative offset from
146 the current cell being evaluated. The format is map[r,c], where r is
147 the row offset and c is the column offset. For example, map[1,2]
148 refers to the cell one row below and two columns to the right of the
149 current cell, map[-2,-1] refers to the cell two rows above and one col‐
150 umn to the left of the current cell, and map[0,1] refers to the cell
151 one column to the right of the current cell. This syntax permits the
152 development of neighborhood-type filters within a single map or across
153 multiple maps.
154
156 Sometimes it is desirable to use a value associated with a category's
157 label instead of the category value itself. If a raster map layer name
158 is preceded by the @ operator, then the labels in the category file for
159 the raster map layer are used in the expression instead of the category
160 value.
161 For example, suppose that the raster map layer soil.ph (representing
163 soil pH values) has a category file with labels as follows:
164 cat label
167 ------------------
168 0 no data
169 1 1.4
170 2 2.4
171 3 3.5
172 4 5.8
173 5 7.2
174 6 8.8
175 7 9.4
176 Then the expression:
179 result = @soils.ph
182 would produce a result with category values 0, 1.4, 2.4, 3.5, 5.8, 7.2,
185 8.8 and 9.4.
186 Note that this operator may only be applied to raster map layers and
188 produces a floating point value in the expression. Therefore, the cat‐
189 egory label must start with a valid number. If the category label is
190 integer, it will be represented by a floating point number. I the cate‐
191 gory label does not start with a number or is missing, it will be rep‐
192 resented by NULL (no data) in the resulting raster map.
193
195 It is often helpful to manipulate the colors assigned to map cate‐
196 gories. This is particularly useful when the spectral properties of
197 cells have meaning (as with imagery data), or when the map category
198 values represent real quantities (as when category values reflect true
199 elevation values). Map color manipulation can also aid visual recogni‐
200 tion, and map printing.
201 The # operator can be used to either convert map category values to
203 their grey scale equivalents or to extract the red, green, or blue com‐
204 ponents of a raster map layer into separate raster map layers.
205 result = #map
208 converts each category value in map to a value in the range 0-255 which
211 represents the grey scale level implied by the color for the category.
212 If the map has a grey scale color table, then the grey level is what
213 #map evaluates to. Otherwise, it is computed as:
214 0.10 * red + 0.81 * green + 0.01 * blue
217 Alternatively, you can use:
220 result = y#map
223 to use the NTSC weightings:
226 0.30 * red + 0.59 * green + 0.11 * blue
229 Or, you can use:
232 result = i#map
235 to use equal weightings:
238 0.33 * red + 0.33 * green + 0.33 * blue
241 The # operator has three other forms: r#map, g#map, b#map. These
244 extract the red, green, or blue components in the named raster map,
245 respectively. The GRASS shell script r.blend extracts each of these
246 components from two raster map layers, and combines them by a user-
247 specified percentage. These forms allow color separates to be made.
248 For example, to extract the red component from map and store it in the
249 new 0-255 map layer red, the user could type:
250 red = r#map
253 To assign this map grey colors type:
256 r.colors map=red color=rules
259 black
260 white
261 To assign this map red colors type:
264 r.colors map=red color=rules
267 black
268 red
269
272 The functions currently supported are listed in the table below. The
273 type of the result is indicated in the last column. F means that the
274 functions always results in a floating point value, I means that the
275 function gives an integer result, and * indicates that the result is
276 float if any of the arguments to the function are floating point values
277 and integer if all arguments are integer.
278 function description type
281 ---------------------------------------------------------------------------
282 abs(x) return absolute value of x *
283 acos(x) inverse cosine of x (result is in degrees) F
284 asin(x) inverse sine of x (result is in degrees) F
285 atan(x) inverse tangent of x (result is in degrees) F
286 atan(x,y) inverse tangent of y/x (result is in degrees) F
287 cos(x) cosine of x (x is in degrees) F
288 double(x) convert x to double-precision floating point F
289 eval([x,y,...,]z) evaluate values of listed expr, pass results to z
290 exp(x) exponential function of x F
291 exp(x,y) x to the power y F
292 float(x) convert x to single-precision floating point F
293 graph(x,x1,y1[x2,y2..]) convert the x to a y based on points in a
294 graph F
295 if decision options: *
296 if(x) 1 if x not zero, 0 otherwise
297 if(x,a) a if x not zero, 0 otherwise
298 if(x,a,b) a if x not zero, b otherwise
299 if(x,a,b,c) a if x > 0, b if x is zero, c if x < 0
300 int(x) convert x to integer [ truncates ] I
301 isnull(x) check if x = NULL
302 log(x) natural log of x F
303 log(x,b) log of x base b F
304 max(x,y[,z...]) largest value of those listed *
305 median(x,y[,z...]) median value of those listed *
306 min(x,y[,z...]) smallest value of those listed *
307 mode(x,y[,z...]) mode value of those listed
308 *
309 not(x) 1 if x is zero, 0 otherwise
310 pow(x,y) x to the power y *
311 rand(a,b) random value x : a <= x < b
312 round(x) round x to nearest integer I
313 sin(x) sine of x (x is in degrees) F
314 sqrt(x) square root of x F
315 tan(x) tangent of x (x is in degrees) F
316 xor(x,y) exclusive-or (XOR) of x and y I
317 Internal variables:
319 row() current row of moving window
320 col() current col of moving window
321 x() current x-coordinate of moving window
322 y() current y-coordinate of moving window
323 ewres() current east-west resolution
324 nsres() current north-south resolution
325 null() NULL value
326 Note, that the row() and col() indexing starts with 1.
327
329 Floating point numbers are allowed in the expression. A floating point
330 number is a number which contains a decimal point:
331 2.3 12.0 12. .81
332 Floating point values in the expression are handled in a special way.
333 With arithmetic and logical operators, if either operand is float, the
334 other is converted to float and the result of the operation is float.
335 This means, in particular that division of integers results in a (trun‐
336 cated) integer, while division of floats results in an accurate float‐
337 ing point value. With functions of type * (see table above), the
338 result is float if any argument is float, integer otherwise.
339 Note: If you calculate with integer numbers, the resulting map will be
341 integer. If you want to get a float result, add the decimal point to
342 integer number(s).
343 If you want floating point division, at least one of the arguments has
345 to be a floating point value. Multiplying one of them by 1.0 will pro‐
346 duce a floating-point result, as will using float():
347 r.mapcalc "ndvi=float(lsat.4 - lsat.3) / (lsat.4 + lsat.3)"
348
351 Division by zero should result in NULL. Modulus by zero should result
352 in NULL. NULL-values in any arithmetic or logical operation should
353 result in NULL. (however, &&& and ||| are treated specially, as
354 described below). The &&& and ||| operators observe the following
355 axioms even when x is NULL:
356 x &&& false == false
357 false &&& x == false
358 x ||| true == true
359 true ||| x == true
360 NULL-values in function arguments should result in NULL (however,
361 if(), eval() and isnull() are treated specially, as described below).
362 The eval() function always returns its last argument The situation for
363 if() is:
364 if(x)
365 NULL if x is NULL; 0 if x is zero; 1 otherwise
366 if(x,a)
367 NULL if x is NULL; a if x is non-zero; 0 otherwise
368 if(x,a,b)
369 NULL if x is NULL; a if x is non-zero; b otherwise
370 if(x,n,z,p)
371 NULL if x is NULL; n if x is negative;
372 z if x is zero; p if x is positive
373 The (new) function isnull(x) returns: 1 if x is NULL; 0 otherwise. The
374 (new) function null() (which has no arguments) returns an integer NULL.
375 Non-NULL, but invalid, arguments to functions should result in NULL.
376 Examples:
377 log(-2)
378 sqrt(-2)
379 pow(a,b) where a is negative and b is not an integer
380 NULL support: Please note that any math performed with NULL cells
383 always results in a NULL value for these cells. If you want to replace
384 a NULL cell on-the-fly, use the isnull() test function in a if-state‐
385 ment.
386 Example: The users wants the NULL-valued cells to be treated like
388 zeros. To add maps A and B (where B contains NULLs) to get a map C the
389 user can use a construction like:
390 C=A + if(isnull(B),0,B)
393 NULL and conditions:
396 For the one argument form:
398 if(x) = NULL if x is NULL
399 if(x) = 0 if x = 0
400 if(x) = 1 otherwise (i.e. x is neither NULL nor 0).
401 For the two argument form:
404 if(x,a) = NULL if x is NULL
405 if(x,a) = 0 if x = 0
406 if(x,a) = a otherwise (i.e. x is neither NULL nor 0).
407 For the three argument form:
410 if(x,a,b) = NULL if x is NULL
411 if(x,a,b) = b if x = 0
412 if(x,a,b) = a otherwise (i.e. x is neither NULL nor 0).
413 For the four argument form:
416 if(x,a,b,c) = NULL if x is NULL
417 if(x,a,b,c) = a if x > 0
418 if(x,a,b,c) = b if x = 0
419 if(x,a,b,c) = c if x < 0
420 More generally, all operators and most functions return NULL if *any*
421 of their arguments are NULL.
422 The functions if(), isnull() and eval() are exceptions.
423 The function isnull() returns 1 if its argument is NULL and 0 other‐
424 wise. If the user wants the opposite, the ! operator, e.g.
425 "!isnull(x)" must be used.
426 All forms of if() return NULL if the first argument is NULL. The 2, 3
428 and 4 argument forms of if() return NULL if the "selected" argument is
429 NULL, e.g.:
430 if(0,a,b) = b regardless of whether a is NULL
431 if(1,a,b) = a regardless of whether b is NULL
432 eval() always returns its last argument, so it only returns NULL if
433 the last argument is NULL.
434 Note: The user cannot test for NULL using the == operator, as that
436 returns NULL if either or both arguments are NULL, i.e. if x and y are
437 both NULL, then "x == y" and "x != y" are both NULL rather than 1 and 0
438 respectively.
439 The behaviour makes sense if the user considers NULL as representing an
440 unknown quantity. E.g. if x and y are both unknown, then the values of
441 "x == y" and "x != y" are also unknown; if they both have unknown val‐
442 ues, the user doesn't know whether or not they both have the same
443 value.
444
446 To compute the average of two raster map layers a and b:
447 ave = (a + b)/2
448 To form a weighted average:
449 ave = (5*a + 3*b)/8.0
450 To produce a binary representation of the raster map layer a so that
451 category 0 remains 0 and all other categories become 1:
452 mask = a != 0
453 This could also be accomplished by:
454 mask = if(a)
455 To mask raster map layer b by raster map layer a:
456 result = if(a,b)
457 To change all values below 5 to NULL:
458 newmap = if(map<5, null(), 5)
459 The graph function allows users to specify a x-y conversion using
460 pairs of x,y coordinates. In some situations a transformation from one
461 value to another is not easily established mathematically, but can be
462 represented by a 2-D graph. The graph() function provides the opportu‐
463 nity to accomplish this. An x-axis value is provided to the graph
464 function along with the associated graph represented by a series of x,y
465 pairs. The x values must be monotonically increasing (each larger than
466 or equal to the previous). The graph function linearly interpolates
467 between pairs. Any x value lower the lowest x value (i.e. first) will
468 have the associated y value returned. Any x value higher than the last
469 will similarly have the associated y value returned. Consider the
470 request:
471 newmap = graph(map, 1,10, 2,25, 3,50)
472 X (map) values supplied and y (newmap) values returned:
473 0, 10
474 1, 10,
475 1.5, 16.5
476 2.9, 47.5
477 4, 50
478 100, 50
479
482 Extra care must be taken if the expression is given on the command
483 line. Some characters have special meaning to the UNIX shell. These
484 include, among others:
485 * ( ) > & |
487 It is advisable to put single quotes around the expression; e.g.:
489 result = 'elevation * 2'
490 Without the quotes, the *, which has special meaning to the UNIX
491 shell, would be altered and r.mapcalc would see something other than
492 the *.
493 If the input comes directly from the keyboard and the result raster map
495 layer exists, the user will be asked if it can be overwritten. Other‐
496 wise, the result raster map layer will automatically be overwritten if
497 it exists.
498 Quoting result is not allowed. However, it is never necessary to quote
500 result since it is always taken to be a raster map layer name.
501 For formulas that the user enters from standard input (rather than from
503 the command line), a line continuation feature now exists. If the user
504 adds \e to the end of an input line, r.mapcalc assumes that the formula
505 being entered by the user continues on to the next input line. There
506 is no limit to the possible number of input lines or to the length of a
507 formula.
508 If the r.mapcalc formula entered by the user is very long, the map
510 title will contain only some of it, but most (if not all) of the for‐
511 mula will be placed into the history file for the result map.
512 When the user enters input to r.mapcalc non-interactively on the com‐
514 mand line, the program will not warn the user not to overwrite existing
515 map layers. Users should therefore take care to assign program outputs
516 raster map names that do not yet exist in their current mapsets.
517 The environment variable GRASS_RND_SEED is read to initialise the ran‐
519 dom number generator.
520
522 Continuation lines must end with a \ and have NO trailing white space
523 (blanks or tabs). If the user does leave white space at the end of
524 continuation lines, the error messages produced by r.mapcalc will be
525 meaningless and the equation will not work as the user intended. This
526 is important for the eval() function.
527 Error messages produced by r.mapcalc are almost useless. In future,
529 r.mapcalc should make some attempt to point the user to the offending
530 section of the equation, e.g.:
531 x = a * b ++ c
532 ERROR: somewhere in line 1: ... b ++ c ...
533 Currently, there is no comment mechanism in r.mapcalc. Perhaps adding
536 a capability that would cause the entire line to be ignored when the
537 user inserted a # at the start of a line as if it were not present,
538 would do the trick.
539 The function should require the user to type "end" or "exit" instead of
541 simply a blank line. This would make separation of multiple scripts
542 separable by white space.
543 r.mapcalc does not print a warning in case of operations on NULL cells.
545 It is left to the user to utilize the isnull() function.
546
548 r.mapcalc: An Algebra for GIS and Image Processing, by Michael Shapiro
549 and Jim Westervelt, U.S. Army Construction Engineering Research Labora‐
550 tory (March/1991).
551 Performing Map Calculations on GRASS Data: r.mapcalc Program Tutorial,
553 by Marji Larson, Michael Shapiro and Scott Tweddale, U.S. Army Con‐
554 struction Engineering Research Laboratory (December 1991)
555 Grey scale conversion is based on the C.I.E. x,y,z system where y rep‐
557 resents luminance. See "Fundamentals of Digital Image Processing," by
558 Anil K. Jain (Prentice Hall, NJ, 1989; p 67).
559 g.region, r.bitpattern, r.blend, r.colors, r.fillnulls
561
563 Michael Shapiro, U.S.Army Construction Engineering Research Laboratory
564 Glynn Clements
566 Last changed: $Date: 2007-08-24 06:13:51 +0200 (Fri, 24 Aug 2007) $
568GRASS 6.3.0 r.mapcalc(1)