GRASS 8 Programmer's Manual
8.6.0dev(2026)-f6f2c534ea
Loading...
Searching...
No Matches
gradient.c
Go to the documentation of this file.
1
/*!
2
\file gradient.c
3
4
\brief Gradient computation
5
6
(C) 2014 by the GRASS Development Team
7
8
This program is free software under the GNU General Public
9
License (>=v2). Read the file COPYING that comes with GRASS
10
for details.
11
12
\author Anna Petrasova
13
*/
14
15
/*!
16
\brief Gradient computation
17
18
Gradient computation (second order approximation)
19
using central differencing scheme (plus forward and backward
20
difference of second order approximation). When one or more of the cells,
21
from which the gradient for a particular cell is computed, is null,
22
gradient for that particular cell is set to 0.
23
24
\param array pointer to RASTER3D_Array with input values
25
\param step array of x, y, z steps for gradient (resolution values)
26
\param[out] grad_x pointer to RASTER3D_Array_double with gradient in x
27
direction \param[out] grad_y pointer to RASTER3D_Array_double with gradient
28
in y direction \param[out] grad_z pointer to RASTER3D_Array_double with
29
gradient in z direction
30
31
*/
32
#include <
grass/raster3d.h
>
33
34
void
Rast3d_gradient_double
(
RASTER3D_Array_double
*array,
double
*
step
,
35
RASTER3D_Array_double
*
grad_x
,
36
RASTER3D_Array_double
*
grad_y
,
37
RASTER3D_Array_double
*
grad_z
)
38
{
39
int
col
, row, depth;
40
double
val0
,
val1
,
val2
;
41
42
for
(depth = 0; depth < array->
sz
; depth++) {
43
for
(row = 0; row < array->
sy
; row++) {
44
/* row start */
45
val0
=
RASTER3D_ARRAY_ACCESS
(array, 0, row, depth);
46
val1
=
RASTER3D_ARRAY_ACCESS
(array, 1, row, depth);
47
val2
=
RASTER3D_ARRAY_ACCESS
(array, 2, row, depth);
48
if
(
Rast_is_d_null_value
(&
val0
))
49
Rast_set_null_value
(
50
&
RASTER3D_ARRAY_ACCESS
(
grad_x
, 0, row, depth), 1,
51
DCELL_TYPE
);
52
else
if
(
Rast_is_d_null_value
(&
val1
) ||
Rast_is_d_null_value
(&
val2
))
53
RASTER3D_ARRAY_ACCESS
(
grad_x
, 0, row, depth) = 0;
54
else
55
RASTER3D_ARRAY_ACCESS
(
grad_x
, 0, row, depth) =
56
(-3 *
val0
+ 4 *
val1
-
val2
) / (2 *
step
[0]);
57
58
/* row end */
59
val0
=
RASTER3D_ARRAY_ACCESS
(array, array->
sx
- 3, row, depth);
60
val1
=
RASTER3D_ARRAY_ACCESS
(array, array->
sx
- 2, row, depth);
61
val2
=
RASTER3D_ARRAY_ACCESS
(array, array->
sx
- 1, row, depth);
62
if
(
Rast_is_d_null_value
(&
val2
))
63
Rast_set_null_value
(
64
&
RASTER3D_ARRAY_ACCESS
(
grad_x
, array->
sx
- 1, row, depth),
65
1,
DCELL_TYPE
);
66
else
if
(
Rast_is_d_null_value
(&
val0
) ||
Rast_is_d_null_value
(&
val1
))
67
RASTER3D_ARRAY_ACCESS
(
grad_x
, array->
sx
- 1, row, depth) = 0;
68
else
69
RASTER3D_ARRAY_ACCESS
(
grad_x
, array->
sx
- 1, row, depth) =
70
(3 *
val2
- 4 *
val1
+
val0
) / (2 *
step
[0]);
71
72
/* row */
73
for
(
col
= 1;
col
< array->
sx
- 1;
col
++) {
74
val0
=
RASTER3D_ARRAY_ACCESS
(array,
col
- 1, row, depth);
75
val1
=
RASTER3D_ARRAY_ACCESS
(array,
col
, row, depth);
76
val2
=
RASTER3D_ARRAY_ACCESS
(array,
col
+ 1, row, depth);
77
if
(
Rast_is_d_null_value
(&
val1
))
78
Rast_set_null_value
(
79
&
RASTER3D_ARRAY_ACCESS
(
grad_x
,
col
, row, depth), 1,
80
DCELL_TYPE
);
81
else
if
(
Rast_is_d_null_value
(&
val0
) ||
82
Rast_is_d_null_value
(&
val2
))
83
RASTER3D_ARRAY_ACCESS
(
grad_x
,
col
, row, depth) = 0;
84
else
85
RASTER3D_ARRAY_ACCESS
(
grad_x
,
col
, row, depth) =
86
(
val2
-
val0
) / (2 *
step
[0]);
87
}
88
}
89
}
90
for
(depth = 0; depth < array->
sz
; depth++) {
91
for
(
col
= 0;
col
< array->
sx
;
col
++) {
92
/* col start */
93
val0
=
RASTER3D_ARRAY_ACCESS
(array,
col
, 0, depth);
94
val1
=
RASTER3D_ARRAY_ACCESS
(array,
col
, 1, depth);
95
val2
=
RASTER3D_ARRAY_ACCESS
(array,
col
, 2, depth);
96
if
(
Rast_is_d_null_value
(&
val0
))
97
Rast_set_null_value
(
98
&
RASTER3D_ARRAY_ACCESS
(
grad_y
,
col
, 0, depth), 1,
99
DCELL_TYPE
);
100
else
if
(
Rast_is_d_null_value
(&
val1
) ||
Rast_is_d_null_value
(&
val2
))
101
RASTER3D_ARRAY_ACCESS
(
grad_y
,
col
, 0, depth) = 0;
102
else
103
RASTER3D_ARRAY_ACCESS
(
grad_y
,
col
, 0, depth) =
104
-(-3 *
val0
+ 4 *
val1
-
val2
) / (2 *
step
[1]);
105
106
/* col end */
107
val0
=
RASTER3D_ARRAY_ACCESS
(array,
col
, array->
sy
- 3, depth);
108
val1
=
RASTER3D_ARRAY_ACCESS
(array,
col
, array->
sy
- 2, depth);
109
val2
=
RASTER3D_ARRAY_ACCESS
(array,
col
, array->
sy
- 1, depth);
110
if
(
Rast_is_d_null_value
(&
val2
))
111
Rast_set_null_value
(
112
&
RASTER3D_ARRAY_ACCESS
(
grad_y
,
col
, array->
sy
- 1, depth),
113
1,
DCELL_TYPE
);
114
else
if
(
Rast_is_d_null_value
(&
val0
) ||
Rast_is_d_null_value
(&
val1
))
115
RASTER3D_ARRAY_ACCESS
(
grad_y
,
col
, array->
sy
- 1, depth) = 0;
116
else
117
RASTER3D_ARRAY_ACCESS
(
grad_y
,
col
, array->
sy
- 1, depth) =
118
-(3 *
val2
- 4 *
val1
+
val0
) / (2 *
step
[1]);
119
120
/* col */
121
for
(row = 1; row < array->
sy
- 1; row++) {
122
val0
=
RASTER3D_ARRAY_ACCESS
(array,
col
, row - 1, depth);
123
val1
=
RASTER3D_ARRAY_ACCESS
(array,
col
, row, depth);
124
val2
=
RASTER3D_ARRAY_ACCESS
(array,
col
, row + 1, depth);
125
if
(
Rast_is_d_null_value
(&
val1
))
126
Rast_set_null_value
(
127
&
RASTER3D_ARRAY_ACCESS
(
grad_y
,
col
, row, depth), 1,
128
DCELL_TYPE
);
129
else
if
(
Rast_is_d_null_value
(&
val0
) ||
130
Rast_is_d_null_value
(&
val2
))
131
RASTER3D_ARRAY_ACCESS
(
grad_y
,
col
, row, depth) = 0;
132
else
133
RASTER3D_ARRAY_ACCESS
(
grad_y
,
col
, row, depth) =
134
-(
val2
-
val0
) / (2 *
step
[1]);
135
}
136
}
137
}
138
for
(row = 0; row < array->
sy
; row++) {
139
for
(
col
= 0;
col
< array->
sx
;
col
++) {
140
/* vertical col start */
141
val0
=
RASTER3D_ARRAY_ACCESS
(array,
col
, row, 0);
142
val1
=
RASTER3D_ARRAY_ACCESS
(array,
col
, row, 1);
143
val2
=
RASTER3D_ARRAY_ACCESS
(array,
col
, row, 2);
144
if
(
Rast_is_d_null_value
(&
val0
))
145
Rast_set_null_value
(&
RASTER3D_ARRAY_ACCESS
(
grad_z
,
col
, row, 0),
146
1,
DCELL_TYPE
);
147
else
if
(
Rast_is_d_null_value
(&
val1
) ||
Rast_is_d_null_value
(&
val2
))
148
RASTER3D_ARRAY_ACCESS
(
grad_z
,
col
, row, 0) = 0;
149
else
150
RASTER3D_ARRAY_ACCESS
(
grad_z
,
col
, row, 0) =
151
(-3 *
val0
+ 4 *
val1
-
val2
) / (2 *
step
[2]);
152
153
/* vertical col end */
154
val0
=
RASTER3D_ARRAY_ACCESS
(array,
col
, row, array->
sz
- 3);
155
val1
=
RASTER3D_ARRAY_ACCESS
(array,
col
, row, array->
sz
- 2);
156
val2
=
RASTER3D_ARRAY_ACCESS
(array,
col
, row, array->
sz
- 1);
157
if
(
Rast_is_d_null_value
(&
val2
))
158
Rast_set_null_value
(
159
&
RASTER3D_ARRAY_ACCESS
(
grad_z
,
col
, row, array->
sz
- 1), 1,
160
DCELL_TYPE
);
161
else
if
(
Rast_is_d_null_value
(&
val0
) ||
Rast_is_d_null_value
(&
val1
))
162
RASTER3D_ARRAY_ACCESS
(
grad_z
,
col
, row, array->
sz
- 1) = 0;
163
else
164
RASTER3D_ARRAY_ACCESS
(
grad_z
,
col
, row, array->
sz
- 1) =
165
(3 *
val2
- 4 *
val1
+
val0
) / (2 *
step
[2]);
166
/* vertical col */
167
for
(depth = 1; depth < array->
sz
- 1; depth++) {
168
val0
=
RASTER3D_ARRAY_ACCESS
(array,
col
, row, depth - 1);
169
val1
=
RASTER3D_ARRAY_ACCESS
(array,
col
, row, depth);
170
val2
=
RASTER3D_ARRAY_ACCESS
(array,
col
, row, depth + 1);
171
if
(
Rast_is_d_null_value
(&
val1
))
172
Rast_set_null_value
(
173
&
RASTER3D_ARRAY_ACCESS
(
grad_z
,
col
, row, depth), 1,
174
DCELL_TYPE
);
175
else
if
(
Rast_is_d_null_value
(&
val0
) ||
176
Rast_is_d_null_value
(&
val2
))
177
RASTER3D_ARRAY_ACCESS
(
grad_z
,
col
, row, depth) = 0;
178
else
179
RASTER3D_ARRAY_ACCESS
(
grad_z
,
col
, row, depth) =
180
(
val2
-
val0
) / (2 *
step
[2]);
181
}
182
}
183
}
184
}
AMI_STREAM
Definition
ami_stream.h:153
Rast_set_null_value
void Rast_set_null_value(void *, int, RASTER_MAP_TYPE)
To set one or more raster values to null.
Definition
null_val.c:98
Rast_is_d_null_value
#define Rast_is_d_null_value(dcellVal)
Definition
defs/raster.h:417
Rast3d_gradient_double
void Rast3d_gradient_double(RASTER3D_Array_double *array, double *step, RASTER3D_Array_double *grad_x, RASTER3D_Array_double *grad_y, RASTER3D_Array_double *grad_z)
Gradient computation.
Definition
gradient.c:34
raster3d.h
RASTER3D_ARRAY_ACCESS
#define RASTER3D_ARRAY_ACCESS(arr, x, y, z)
Definition
raster3d.h:267
DCELL_TYPE
#define DCELL_TYPE
Definition
raster.h:13
RASTER3D_Array_double
Definition
raster3d.h:259
RASTER3D_Array_double::sx
int sx
Definition
raster3d.h:261
RASTER3D_Array_double::sz
int sz
Definition
raster3d.h:263
RASTER3D_Array_double::sy
int sy
Definition
raster3d.h:262
lib
raster3d
gradient.c
Generated on Sun Apr 12 2026 06:59:17 for GRASS 8 Programmer's Manual by
1.9.8