1r3.mapcalc(1) GRASS GIS User's Manual r3.mapcalc(1)
2
6 r3.mapcalc - Raster map calculator.
7
9 raster, algebra
10
12 r3.mapcalc
13 r3.mapcalc --help
14 r3.mapcalc [-sl] [expression=string] [region=string] [file=name]
15 [seed=integer] [--overwrite] [--help] [--verbose] [--quiet]
16 [--ui]
17 Flags:
19 -s
20 Generate random seed (result is non-deterministic)
21 -l
23 List input and output maps
24 --overwrite
26 Allow output files to overwrite existing files
27 --help
29 Print usage summary
30 --verbose
32 Verbose module output
33 --quiet
35 Quiet module output
36 --ui
38 Force launching GUI dialog
39 Parameters:
41 expression=string
42 Expression to evaluate
43 region=string
45 The computational region that should be used.
46 - current uses the current region of the mapset.
47 - intersect computes the intersection region between
48 all input maps and uses the smallest resolution
49 - union computes the union extent of all map regions
50 and uses the smallest resolution
51 Options: current, intersect, union
52 Default: current
53 file=name
55 File containing expression(s) to evaluate
56 seed=integer
58 Seed for rand() function
59
61 r3.mapcalc performs arithmetic on raster map layers. New raster map
62 layers can be created which are arithmetic expressions involving exist‐
63 ing raster map layers, integer or floating point constants, and func‐
64 tions.
65 Program use
67 r3.mapcalc expression have the form:
68 result = expression
70 where result is the name of a raster map layer to contain the result of
72 the calculation and expression is any legal arithmetic expression
73 involving existing raster map layers (except result itself), integer or
74 floating point constants, and functions known to the calculator.
75 Parentheses are allowed in the expression and may be nested to any
76 depth. result will be created in the user’s current mapset.
77 As expression= is the first option, it is the default. This means that
79 passing an expression on the command line is possible as long as the
80 expression is quoted and a space is included before the first = sign.
81 Example (’foo’ is the resulting map):
82 r3.mapcalc "foo = 1"
83 or:
84 r3.mapcalc ’foo = 1’
85 An unquoted expression (i.e. split over multiple arguments) won’t work,
86 nor will omitting the space before the = sign:
87 r3.mapcalc ’foo=1’
88 Sorry, <foo> is not a valid parameter
89 If no options are given, it manufactures "file=-" (which reads from
90 stdin), so you can continue to use e.g.:
91 r3.mapcalc < file
92 or:
93 r3.mapcalc <<EOF
94 foo = 1
95 EOF
96 But unless you need compatibility with previous versions, use file=
97 explicitly, e.g.:
98 r3.mapcalc file=file
99 or:
100 r3.mapcalc file=- <<EOF
101 foo = 1
102 EOF
103 The formula entered to r3.mapcalc by the user is recorded both in the
105 result map title (which appears in the category file for result) and in
106 the history file for result.
107 Some characters have special meaning to the command shell. If the user
109 is entering input to r3.mapcalc on the command line, expressions should
110 be enclosed within single quotes. See NOTES, below.
111 Operators and order of precedence
113 The following operators are supported:
114 Operator Meaning Type Precedence
115 --------------------------------------------------------------
116 - negation Arithmetic 12
117 ~ one’s complement Bitwise 12
118 ! not Logical 12
119 ^ exponentiation Arithmetic 11
120 % modulus Arithmetic 10
121 / division Arithmetic 10
122 * multiplication Arithmetic 10
123 + addition Arithmetic 9
124 - subtraction Arithmetic 9
125 << left shift Bitwise 8
126 >> right shift Bitwise 8
127 >>> right shift (unsigned) Bitwise 8
128 > greater than Logical 7
129 >= greater than or equal Logical 7
130 < less than Logical 7
131 <= less than or equal Logical 7
132 == equal Logical 6
133 != not equal Logical 6
134 & bitwise and Bitwise 5
135 | bitwise or Bitwise 4
136 && logical and Logical 3
137 &&& logical and[1] Logical 3
138 || logical or Logical 2
139 ||| logical or[1] Logical 2
140 ?: conditional Logical 1
141 (modulus is the remainder upon division)
142 [1] The &&& and ||| operators handle null values differently to other
144 operators. See the section entitled NULL support below for more
145 details.
146 The operators are applied from left to right, with those of higher
148 precedence applied before those with lower precedence. Division by 0
149 and modulus by 0 are acceptable and give a NULL result. The logical
150 operators give a 1 result if the comparison is true, 0 otherwise.
151 3D Grid names
153 Anything in the expression which is not a number, operator, or function
154 name is taken to be a 3D grid name. Examples:
155 volume
157 x3
158 3d.his
159 Most GRASS raster map layers and 3D grids meet this naming convention.
161 However, if a 3D grid has a name which conflicts with the above rule,
162 it should be quoted. For example, the expression
163 x = a-b
165 would be interpreted as: x equals a minus b, whereas
167 x = "a-b"
169 would be interpreted as: x equals the 3D grid named a-b
171 Also
173 x = 3107
175 would create x filled with the number 3107, while
177 x = "3107"
179 would copy the 3D grid 3107 to the 3D grid x.
181 Quotes are not required unless the 3D grid names look like numbers or
183 contain operators, OR unless the program is run non-interactively.
184 Examples given here assume the program is run interactively. See NOTES,
185 below.
186 r3.mapcalc will look for the 3D grids according to the user’s current
188 mapset search path. It is possible to override the search path and
189 specify the mapset from which to select the 3D grid. This is done by
190 specifying the 3D grid name in the form:
191 name@mapset
193 For example, the following is a legal expression:
195 result = x@PERMANENT / y@SOILS
197 The mapset specified does not have to be in the mapset search path.
199 (This method of overriding the mapset search path is common to all
200 GRASS commands, not just r3.mapcalc.)
201 The neighborhood modifier
203 3D grids are data base files stored in voxel format, i.e., three-dimen‐
204 sional matrices of float/double values. In r3.mapcalc, 3D grids may be
205 followed by a neighborhood modifier that specifies a relative offset
206 from the current cell being evaluated. The format is map[r,c,d], where
207 r is the row offset, c is the column offset and d is the depth offset.
208 For example, map[1,2,3] refers to the cell one row below, two columns
209 to the right and 3 levels below of the current cell, map[-3,-2,-1]
210 refers to the cell three rows above, two columns to the left and one
211 level below of the current cell, and map[0,1,0] refers to the cell one
212 column to the right of the current cell. This syntax permits the devel‐
213 opment of neighborhood-type filters within a single 3D grid or across
214 multiple 3D grids.
215 Functions
217 The functions currently supported are listed in the table below. The
218 type of the result is indicated in the last column. F means that the
219 functions always results in a floating point value, I means that the
220 function gives an integer result, and * indicates that the result is
221 float if any of the arguments to the function are floating point values
222 and integer if all arguments are integer.
223 function description type
225 ---------------------------------------------------------------------------
226 abs(x) return absolute value of x *
227 acos(x) inverse cosine of x (result is in degrees) F
228 asin(x) inverse sine of x (result is in degrees) F
229 atan(x) inverse tangent of x (result is in degrees) F
230 atan(x,y) inverse tangent of y/x (result is in degrees) F
231 ceil(x) the smallest integral value not less than x *
232 cos(x) cosine of x (x is in degrees) F
233 double(x) convert x to double-precision floating point F
234 eval([x,y,...,]z) evaluate values of listed expr, pass results to z
235 exp(x) exponential function of x F
236 exp(x,y) x to the power y F
237 float(x) convert x to single-precision floating point F
238 floor(x) the largest integral value not greater than x *
239 graph(x,x1,y1[x2,y2..]) convert the x to a y based on points in a graph F
240 graph2(x,x1[,x2,..],y1[,y2..])
241 alternative form of graph() F
242 if decision options: *
243 if(x) 1 if x not zero, 0 otherwise
244 if(x,a) a if x not zero, 0 otherwise
245 if(x,a,b) a if x not zero, b otherwise
246 if(x,a,b,c) a if x > 0, b if x is zero, c if x < 0
247 int(x) convert x to integer [ truncates ] I
248 isnull(x) check if x = NULL
249 log(x) natural log of x F
250 log(x,b) log of x base b F
251 max(x,y[,z...]) largest value of those listed *
252 median(x,y[,z...]) median value of those listed *
253 min(x,y[,z...]) smallest value of those listed *
254 mode(x,y[,z...]) mode value of those listed *
255 nmax(x,y[,z...]) largest value of those listed, excluding NULLs *
256 nmedian(x,y[,z...]) median value of those listed, excluding NULLs *
257 nmin(x,y[,z...]) smallest value of those listed, excluding NULLs *
258 nmode(x,y[,z...]) mode value of those listed, excluding NULLs *
259 not(x) 1 if x is zero, 0 otherwise
260 pow(x,y) x to the power y *
261 rand(a,b) random value x : a <= x < b *
262 round(x) round x to nearest integer I
263 round(x,y) round x to nearest multiple of y
264 round(x,y,z) round x to nearest y*i+z for some integer i
265 sin(x) sine of x (x is in degrees) F
266 sqrt(x) square root of x F
267 tan(x) tangent of x (x is in degrees) F
268 xor(x,y) exclusive-or (XOR) of x and y I
269 Internal variables:
270 row() current row of moving window I
271 col() current col of moving window I
272 depth() return current depth I
273 nrows() number of rows in computation region I
274 ncols() number of columns in computation region I
275 ndepths() number of depth levels in computation region I
276 x() current x-coordinate of moving window F
277 y() current y-coordinate of moving window F
278 z() return current z value F
279 ewres() current east-west resolution F
280 nsres() current north-south resolution F
281 tbres() current top-bottom resolution F
282 area() area of current cell in square meters F
283 null() NULL value
284 Note, that the row(), col() and depth() indexing starts with 1.
285 Floating point values in the expression
287 Floating point numbers are allowed in the expression. A floating point
288 number is a number which contains a decimal point:
289 2.3 12.0 12. .81
290 Floating point values in the expression are handled in a special way.
291 With arithmetic and logical operators, if either operand is float, the
292 other is converted to float and the result of the operation is float.
293 This means, in particular that division of integers results in a (trun‐
294 cated) integer, while division of floats results in an accurate float‐
295 ing point value. With functions of type * (see table above), the
296 result is float if any argument is float, integer otherwise.
297 Note: If you calculate with integer numbers, the resulting map will be
299 integer. If you want to get a float result, add the decimal point to
300 integer number(s).
301 If you want floating point division, at least one of the arguments has
303 to be a floating point value. Multiplying one of them by 1.0 will pro‐
304 duce a floating-point result, as will using float():
305 r3.mapcalc "ratio = float(soil.4 - soil.3) / soil.3)"
306 NULL support
308 · Division by zero should result in NULL.
309 · Modulus by zero should result in NULL.
311 · NULL-values in any arithmetic or logical operation should
313 result in NULL. (however, &&& and ||| are treated specially, as
314 described below).
315 · The &&& and ||| operators observe the following axioms even
317 when x is NULL:
318 x &&& false == false
319 false &&& x == false
320 x ||| true == true
321 true ||| x == true
322 · NULL-values in function arguments should result in NULL (how‐
324 ever, if(), eval() and isnull() are treated specially, as
325 described below).
326 · The eval() function always returns its last argument
328 · The situation for if() is:
330 if(x)
331 NULL if x is NULL; 0 if x is zero; 1 otherwise
332 if(x,a)
333 NULL if x is NULL; a if x is non-zero; 0 otherwise
334 if(x,a,b)
335 NULL if x is NULL; a if x is non-zero; b otherwise
336 if(x,n,z,p)
337 NULL if x is NULL; n if x is negative;
338 z if x is zero; p if x is positive
339 · The (new) function isnull(x) returns: 1 if x is NULL; 0 other‐
341 wise. The (new) function null() (which has no arguments)
342 returns an integer NULL.
343 · Non-NULL, but invalid, arguments to functions should result in
345 NULL.
346 Examples:
347 log(-2)
348 sqrt(-2)
349 pow(a,b) where a is negative and b is not an integer
350 NULL support: Please note that any math performed with NULL cells
352 always results in a NULL value for these cells. If you want to replace
353 a NULL cell on-the-fly, use the isnull() test function in a if-state‐
354 ment.
355 Example: The users wants the NULL-valued cells to be treated like
357 zeros. To add maps A and B (where B contains NULLs) to get a map C the
358 user can use a construction like:
359 C = A + if(isnull(B),0,B)
361 NULL and conditions:
363 For the one argument form:
365 if(x) = NULL if x is NULL
366 if(x) = 0 if x = 0
367 if(x) = 1 otherwise (i.e. x is neither NULL nor 0).
368 For the two argument form:
370 if(x,a) = NULL if x is NULL
371 if(x,a) = 0 if x = 0
372 if(x,a) = a otherwise (i.e. x is neither NULL nor 0).
373 For the three argument form:
375 if(x,a,b) = NULL if x is NULL
376 if(x,a,b) = b if x = 0
377 if(x,a,b) = a otherwise (i.e. x is neither NULL nor 0).
378 For the four argument form:
380 if(x,a,b,c) = NULL if x is NULL
381 if(x,a,b,c) = a if x > 0
382 if(x,a,b,c) = b if x = 0
383 if(x,a,b,c) = c if x < 0
384 More generally, all operators and most functions return NULL if *any*
385 of their arguments are NULL.
386 The functions if(), isnull() and eval() are exceptions.
387 The function isnull() returns 1 if its argument is NULL and 0 other‐
388 wise. If the user wants the opposite, the ! operator, e.g.
389 "!isnull(x)" must be used.
390 All forms of if() return NULL if the first argument is NULL. The 2, 3
392 and 4 argument forms of if() return NULL if the "selected" argument is
393 NULL, e.g.:
394 if(0,a,b) = b regardless of whether a is NULL
395 if(1,a,b) = a regardless of whether b is NULL
396 eval() always returns its last argument, so it only returns NULL if the
397 last argument is NULL.
398 Note: The user cannot test for NULL using the == operator, as that
400 returns NULL if either or both arguments are NULL, i.e. if x and y are
401 both NULL, then "x == y" and "x != y" are both NULL rather than 1 and 0
402 respectively.
403 The behaviour makes sense if the user considers NULL as representing an
404 unknown quantity. E.g. if x and y are both unknown, then the values of
405 "x == y" and "x != y" are also unknown; if they both have unknown val‐
406 ues, the user doesn’t know whether or not they both have the same
407 value.
408
410 Usage from command line
411 Extra care must be taken if the expression is given on the command
412 line. Some characters have special meaning to the UNIX shell. These
413 include, among others:
414 * ( ) > & |
415 It is advisable to put single quotes around the expression; e.g.:
417 ’result = volume * 2’
418 Without the quotes, the *, which has special meaning to the UNIX shell,
419 would be altered and r3.mapcalc would see something other than the *.
420 Multiple computations
422 In general, it’s preferable to do as much as possible in each r3.map‐
423 calc command using multi-line input.
424 Backwards compatibility
426 For the backwards compatibility with GRASS 6, if no options are given,
427 it manufactures file=- (which reads from stdin), so you can continue to
428 use e.g.:
429 r3.mapcalc < file
430 or:
431 r3.mapcalc <<EOF
432 foo = 1
433 EOF
434 But unless you need compatibility with previous GRASS GIS versions, use
435 file= explicitly, as stated above.
436 When the map name contains uppercase letter(s) or a dot which are not
438 allowed to be in module option names, the r3.mapcalc command will be
439 valid also without quotes:
440 r3.mapcalc volume_A=1
441 r3.mapcalc volume.1=1
442 However, this syntax is not recommended as quotes as stated above more
443 safe. Using quotes is both backwards compatible and valid in future.
444 Interactive input in command line
446 For formulas that the user enters from standard input (rather than from
447 the command line), a line continuation feature now exists. If the user
448 adds a backslash to the end of an input line, r3.mapcalc assumes that
449 the formula being entered by the user continues on to the next input
450 line. There is no limit to the possible number of input lines or to
451 the length of a formula.
452 If the r3.mapcalc formula entered by the user is very long, the map
454 title will contain only some of it, but most (if not all) of the for‐
455 mula will be placed into the history file for the result map.
456 When the user enters input to r3.mapcalc non-interactively on the com‐
458 mand line, the program will not warn the user not to overwrite existing
459 map layers. Users should therefore take care to assign program outputs
460 raster map names that do not yet exist in their current mapsets.
461 3D GRID MASK handling
463 r3.mapcalc follows the common GRASS behavior of raster MASK handling,
464 so the MASK is only applied when reading an existing GRASS raster map.
465 This implies that, for example, the command:
466 r3.mapcalc "volume_amplified = volume * 3"
467 create a map respecting the masked pixels if MASK is active.
468 However, when creating a map which is not based on any map, e.g. a map
470 from a constant:
471 r3.mapcalc "volume_const = 200.0"
472 the created raster map is limited only by a computation region but it
473 is not affected by an active MASK. This is expected because, as men‐
474 tioned above, MASK is only applied when reading, not when writing a
475 raster map.
476 If also in this case the MASK should be applied, an if() statement
478 including the MASK should be used, e.g.:
479 r3.mapcalc "volume_const = if(MASK, 200.0, null())"
480 When testing MASK related expressions keep in mind that when MASK is
481 active you don’t see data in masked areas even if they are not NULL.
482 See r.mask for details.
483 Random number generator initialization
485 The pseudo-random number generator used by the rand() function can be
486 initialised to a specific value using the seed option. This can be
487 used to replicate a previous calculation.
488 Alternatively, it can be initialised from the system time and the PID
490 using the -r flag. This should result in a different seed being used
491 each time.
492 In either case, the seed will be written to the map’s history, and can
494 be seen using r.info.
495 If you want other people to be able to verify your results, it’s
497 preferable to use the seed option to supply a seed which is either
498 specified in the script or generated from a determenistic process such
499 as a pseudo-random number generator given an explicit seed.
500 Note that the rand() function will generate a fatal error if neither
502 the seed option nor the -s flag are given.
503
505 To compute the average of two 3D grids a and b:
506 ave = (a + b)/2
507 To form a weighted average:
508 ave = (5*a + 3*b)/8.0
509 To produce a binary representation of 3D grid a so that category 0
510 remains 0 and all other categories become 1:
511 mask = a != 0
512 This could also be accomplished by:
513 mask = if(a)
514 To mask 3D grid b by 3D grid a:
515 result = if(a,b)
516 To change all values below 5 to NULL, keep otherwise:
517 newmap = if(map < 5, null(), map)
518 The graph() function allows users to specify a x-y conversion using
519 pairs of x,y coordinates. In some situations a transformation from one
520 value to another is not easily established mathematically, but can be
521 represented by a 2-D graph and then linearly interpolated. The graph()
522 function provides the opportunity to accomplish this. An x-axis value
523 is provided to the graph function along with the associated graph rep‐
524 resented by a series of x,y pairs. The x values must be monotonically
525 increasing (each larger than or equal to the previous). The graph
526 function linearly interpolates between pairs. Any x value lower the
527 lowest x value (i.e. first) will have the associated y value returned.
528 Any x value higher than the last will similarly have the associated y
529 value returned. Consider the request:
530 newmap = graph(map, 1,10, 2,25, 3,50)
531 X (map) values supplied and y (newmap) values returned:
532 0, 10
533 1, 10
534 1.5, 17.5
535 2.9, 47.5
536 4, 50
537 100, 50
538
540 The result variable on the left hand side of the equation should not
541 appear in the expression on the right hand side.
542 mymap = if( mymap > 0, mymap, 0)
543 Any maps generated by a r3.mapcalc command only exist after the entire
545 command has completed. All maps are generated concurrently, row-by-row
546 (i.e. there is an implicit "for row in rows {...}" around the entire
547 expression). Thus the #, @, and [ ] operators cannot be used on a map
548 generated within same r3.mapcalc command run.
549 newmap = oldmap * 3.14
550 othermap = newmap[-1, 0] / newmap[1, 0]
551 Continuation lines must end with a \ and have no trailing white space
553 (blanks or tabs). If the user does leave white space at the end of con‐
554 tinuation lines, the error messages produced by r3.mapcalc will be
555 meaningless and the equation will not work as the user intended. This
556 is particularly important for the eval() function.
557 Currently, there is no comment mechanism in r3.mapcalc. Perhaps adding
559 a capability that would cause the entire line to be ignored when the
560 user inserted a # at the start of a line as if it were not present,
561 would do the trick.
562 The function should require the user to type "end" or "exit" instead of
564 simply a blank line. This would make separation of multiple scripts
565 separable by white space.
566 r3.mapcalc does not print a warning in case of operations on NULL
568 cells. It is left to the user to utilize the isnull() function.
569
571 g.region, r3.colors, r.mapcalc
572
574 r.mapcalc: An Algebra for GIS and Image Processing, by Michael Shapiro
575 and Jim Westervelt, U.S. Army Construction Engineering Research Labora‐
576 tory (March/1991).
577 Performing Map Calculations on GRASS Data: r.mapcalc Program Tutorial,
579 by Marji Larson, Michael Shapiro and Scott Tweddale, U.S. Army Con‐
580 struction Engineering Research Laboratory (December 1991)
581 Neteler, M. (2001): Volume modelling of soils using GRASS GIS 3D tools.
583 - in: Brovelli, M. (ed.)(2001): The Geomatics Workbook N. 2. Politec‐
584 nico di Milano, Italy (ISSN 1591-092X) (PDF)
585
587 Tomas Paudits & Jaro Hofierka, funded by GeoModel s.r.o., Slovakia
588 tpaudits@mailbox.sk, hofierka@geomodel.sk
589 Glynn Clements
591
593 Available at: r3.mapcalc source code (history)
594 Main index | 3D raster index | Topics index | Keywords index | Graphi‐
596 cal index | Full index
597 © 2003-2020 GRASS Development Team, GRASS GIS 7.8.5 Reference Manual
599GRASS 7.8.5 r3.mapcalc(1)