GRASS Programmer's Manual  6.5.svn(2014)-r66266
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
c_reg.c
Go to the documentation of this file.
1 #include <grass/gis.h>
2 
3 #define REGRESSION_SLOPE 0
4 #define REGRESSION_OFFSET 1
5 #define REGRESSION_COEFF_DET 2
6 
7 static void regression(DCELL * result, DCELL * values, int n, int which)
8 {
9  DCELL xsum, ysum;
10  DCELL xbar, ybar;
11  DCELL numer, denom, denom2;
12  int count;
13  int i;
14 
15  xsum = ysum = 0.0;
16  count = 0;
17 
18  for (i = 0; i < n; i++) {
19  if (G_is_d_null_value(&values[i]))
20  continue;
21 
22  xsum += i;
23  ysum += values[i];
24  count++;
25  }
26 
27  if (count < 2) {
28  G_set_d_null_value(result, 1);
29  return;
30  }
31 
32  xbar = xsum / count;
33  ybar = ysum / count;
34 
35  numer = 0.0;
36  for (i = 0; i < n; i++)
37  if (!G_is_d_null_value(&values[i]))
38  numer += i * values[i];
39  numer -= count * xbar * ybar;
40 
41  denom = 0.0;
42  for (i = 0; i < n; i++)
43  if (!G_is_d_null_value(&values[i]))
44  denom += (DCELL) i *i;
45 
46  denom -= count * xbar * xbar;
47 
48  if (which == REGRESSION_COEFF_DET) {
49  denom2 = 0.0;
50  for (i = 0; i < n; i++)
51  if (!G_is_d_null_value(&values[i]))
52  denom2 += values[i] * values[i];
53  denom2 -= count * ybar * ybar;
54  }
55 
56  switch (which) {
57  case REGRESSION_SLOPE:
58  *result = numer / denom;
59  break;
60  case REGRESSION_OFFSET:
61  *result = ybar - xbar * numer / denom;
62  break;
63  case REGRESSION_COEFF_DET:
64  *result = (numer * numer) / (denom * denom2);
65  break;
66  default:
67  G_set_d_null_value(result, 1);
68  break;
69  }
70 
71  /* Check for NaN */
72  if (*result != *result)
73  G_set_d_null_value(result, 1);
74 }
75 
76 void c_reg_m(DCELL * result, DCELL * values, int n, const void *closure)
77 {
78  regression(result, values, n, REGRESSION_SLOPE);
79 }
80 
81 void c_reg_c(DCELL * result, DCELL * values, int n, const void *closure)
82 {
83  regression(result, values, n, REGRESSION_OFFSET);
84 }
85 
86 void c_reg_r2(DCELL * result, DCELL * values, int n, const void *closure)
87 {
88  regression(result, values, n, REGRESSION_COEFF_DET);
89 }
90 
91 static void regression_w(DCELL * result, DCELL(*values)[2], int n, int which)
92 {
93  DCELL xsum, ysum;
94  DCELL xbar, ybar;
95  DCELL numer, denom, denom2;
96  int count;
97  int i;
98 
99  xsum = ysum = 0.0;
100  count = 0;
101 
102  for (i = 0; i < n; i++) {
103  if (G_is_d_null_value(&values[i][0]))
104  continue;
105 
106  xsum += i * values[i][1];
107  ysum += values[i][0] * values[i][1];
108  count += values[i][1];
109  }
110 
111  if (count < 2) {
112  G_set_d_null_value(result, 1);
113  return;
114  }
115 
116  xbar = xsum / count;
117  ybar = ysum / count;
118 
119  numer = 0.0;
120  for (i = 0; i < n; i++)
121  if (!G_is_d_null_value(&values[i][0]))
122  numer += i * values[i][0] * values[i][1];
123  numer -= count * xbar * ybar;
124 
125  denom = 0.0;
126  for (i = 0; i < n; i++)
127  if (!G_is_d_null_value(&values[i][0]))
128  denom += (DCELL) i *i * values[i][1];
129 
130  denom -= count * xbar * xbar;
131 
132  if (which == REGRESSION_COEFF_DET) {
133  denom2 = 0.0;
134  for (i = 0; i < n; i++)
135  if (!G_is_d_null_value(&values[i][0]))
136  denom2 += values[i][0] * values[i][0] * values[i][1];
137  denom2 -= count * ybar * ybar;
138  }
139 
140  switch (which) {
141  case REGRESSION_SLOPE:
142  *result = numer / denom;
143  break;
144  case REGRESSION_OFFSET:
145  *result = ybar - xbar * numer / denom;
146  break;
147  case REGRESSION_COEFF_DET:
148  *result = (numer * numer) / (denom * denom2);
149  break;
150  default:
151  G_set_d_null_value(result, 1);
152  break;
153  }
154 
155  /* Check for NaN */
156  if (*result != *result)
157  G_set_d_null_value(result, 1);
158 }
159 
160 void w_reg_m(DCELL * result, DCELL(*values)[2], int n, const void *closure)
161 {
162  regression_w(result, values, n, REGRESSION_SLOPE);
163 }
164 
165 void w_reg_c(DCELL * result, DCELL(*values)[2], int n, const void *closure)
166 {
167  regression_w(result, values, n, REGRESSION_OFFSET);
168 }
169 
170 void w_reg_r2(DCELL * result, DCELL(*values)[2], int n, const void *closure)
171 {
172  regression_w(result, values, n, REGRESSION_COEFF_DET);
173 }
w_reg_c
void w_reg_c(DCELL *result, DCELL(*values)[2], int n, const void *closure)
Definition: c_reg.c:165
c_reg_r2
void c_reg_r2(DCELL *result, DCELL *values, int n, const void *closure)
Definition: c_reg.c:86
G_set_d_null_value
void G_set_d_null_value(DCELL *dcellVals, int numVals)
Definition: null_val.c:176
count
int count
w_reg_m
void w_reg_m(DCELL *result, DCELL(*values)[2], int n, const void *closure)
Definition: c_reg.c:160
c_reg_c
void c_reg_c(DCELL *result, DCELL *values, int n, const void *closure)
Definition: c_reg.c:81
w_reg_r2
void w_reg_r2(DCELL *result, DCELL(*values)[2], int n, const void *closure)
Definition: c_reg.c:170
G_is_d_null_value
int G_is_d_null_value(const DCELL *dcellVal)
Returns 1 if dcell is NULL, 0 otherwise. This will test if the value dcell is a NaN. Same test as in G_is_f_null_value().
Definition: null_val.c:306
REGRESSION_OFFSET
#define REGRESSION_OFFSET
Definition: c_reg.c:4
c_reg_m
void c_reg_m(DCELL *result, DCELL *values, int n, const void *closure)
Definition: c_reg.c:76
n
int n
Definition: dataquad.c:291
REGRESSION_COEFF_DET
#define REGRESSION_COEFF_DET
Definition: c_reg.c:5
REGRESSION_SLOPE
#define REGRESSION_SLOPE
Definition: c_reg.c:3