GRASS GIS 8 Programmer's Manual 8.2.1(2023)-exported
secpar2d.c
Go to the documentation of this file.
1
2/*!
3 * \file secpar2d.c
4 *
5 * \author H. Mitasova, L. Mitas, I. Kosinovsky, D. Gerdes Fall 1994 (original authors)
6 * \author modified by McCauley in August 1995
7 * \author modified by Mitasova in August 1995
8 * \author H. Mitasova (University of Illinois)
9 * \author L. Mitas (University of Illinois)
10 * \author I. Kosinovsky, (USA-CERL)
11 * \author D.Gerdes (USA-CERL)
12 *
13 * \copyright
14 * (C) 1994-1995 by Helena Mitasova and the GRASS Development Team
15 *
16 * \copyright
17 * This program is free software under the
18 * GNU General Public License (>=v2).
19 * Read the file COPYING that comes with GRASS
20 * for details.
21 *
22 */
23
24
25#include <stdio.h>
26#include <math.h>
27#include <unistd.h>
28#include <grass/gis.h>
29#include <grass/bitmap.h>
30#include <grass/interpf.h>
31
32
33/*!
34 * Compute slope aspect and curvatures
35 *
36 * Computes slope, aspect and curvatures (depending on cond1, cond2) for
37 * derivative arrays adx,...,adxy between columns ngstc and nszc.
38 */
40 int ngstc, /*!< starting column */
41 int nszc, /*!< ending column */
42 int k, /*!< current row */
43 struct BM *bitmask,
44 double *gmin, double *gmax,
45 double *c1min, double *c1max,
46 double *c2min, double *c2max, /*!< min,max interp. values */
47 int cond1,
48 int cond2 /*!< determine if particular values need to be computed */
49 )
50{
51 double dnorm1, ro, /* rad to deg conv */
52 dx2 = 0, dy2 = 0, grad2 = 0, /* gradient squared */
53 slp = 0, grad, /* gradient */
54 oor = 0, /* aspect (orientation) */
55 curn = 0, /* profile curvature */
56 curh = 0, /* tangential curvature */
57 curm = 0, /* mean curvature */
58 temp, /* temp variable */
59 dxy2; /* temp variable square of part diriv. */
60
61 double gradmin;
62 int i, got, bmask = 1;
63 static int first_time_g = 1;
64
65 ro = M_R2D;
66 gradmin = 0.001;
67
68
69 for (i = ngstc; i <= nszc; i++) {
70 if (bitmask != NULL) {
71 bmask = BM_get(bitmask, i, k);
72 }
73 got = 0;
74 if (bmask == 1) {
75 while ((got == 0) && (cond1)) {
76 dx2 = (double)(params->adx[i] * params->adx[i]);
77 dy2 = (double)(params->ady[i] * params->ady[i]);
78 grad2 = dx2 + dy2;
79 grad = sqrt(grad2);
80 /* slope in % slp = 100. * grad; */
81 /* slope in degrees */
82 slp = ro * atan(grad);
83 if (grad <= gradmin) {
84 oor = 0.;
85 got = 3;
86 if (cond2) {
87 curn = 0.;
88 curh = 0.;
89 got = 3;
90 break;
91 }
92 }
93 if (got == 3)
94 break;
95
96 /***********aspect from r.slope.aspect, with adx, ady computed
97 from interpol. function RST **************************/
98
99 if (params->adx[i] == 0.) {
100 if (params->ady[i] > 0.)
101 oor = 90;
102 else
103 oor = 270;
104 }
105 else {
106 oor = ro * atan2(params->ady[i], params->adx[i]);
107 if (oor <= 0.)
108 oor = 360. + oor;
109 }
110
111 got = 1;
112 } /* while */
113 if ((got != 3) && (cond2)) {
114
115 dnorm1 = sqrt(grad2 + 1.);
116 dxy2 =
117 2. * (double)(params->adxy[i] * params->adx[i] *
118 params->ady[i]);
119
120
121 curn =
122 (double)(params->adxx[i] * dx2 + dxy2 +
123 params->adyy[i] * dy2) / (grad2 * dnorm1 *
124 dnorm1 * dnorm1);
125
126 curh =
127 (double)(params->adxx[i] * dy2 - dxy2 +
128 params->adyy[i] * dx2) / (grad2 * dnorm1);
129
130 temp = grad2 + 1.;
131 curm =
132 .5 * ((1. + dy2) * params->adxx[i] - dxy2 +
133 (1. + dx2) * params->adyy[i]) / (temp * dnorm1);
134 }
135 if (first_time_g) {
136 first_time_g = 0;
137 *gmin = *gmax = slp;
138 *c1min = *c1max = curn;
139 *c2min = *c2max = curh;
140 }
141 *gmin = amin1(*gmin, slp);
142 *gmax = amax1(*gmax, slp);
143 *c1min = amin1(*c1min, curn);
144 *c1max = amax1(*c1max, curn);
145 *c2min = amin1(*c2min, curh);
146 *c2max = amax1(*c2max, curh);
147 if (cond1) {
148 params->adx[i] = (FCELL) slp;
149 params->ady[i] = (FCELL) oor;
150 if (cond2) {
151 params->adxx[i] = (FCELL) curn;
152 params->adyy[i] = (FCELL) curh;
153 params->adxy[i] = (FCELL) curm;
154 }
155 }
156 } /* bmask == 1 */
157 }
158 return 1;
159}
int BM_get(struct BM *map, int x, int y)
Gets 'val' from the bitmap.
Definition: bitmap.c:223
#define NULL
Definition: ccmath.h:32
double amin1(double, double)
Definition: minmax.c:67
double amax1(double, double)
Definition: minmax.c:54
int IL_secpar_loop_2d(struct interp_params *params, int ngstc, int nszc, int k, struct BM *bitmask, double *gmin, double *gmax, double *c1min, double *c1max, double *c2min, double *c2max, int cond1, int cond2)
Definition: secpar2d.c:39
DCELL * adxy
Definition: interpf.h:79
DCELL * adyy
Definition: interpf.h:79
DCELL * adx
Definition: interpf.h:78
DCELL * ady
Definition: interpf.h:78
DCELL * adxx
Definition: interpf.h:79