1 #include <stdio.h>
2 #include <stdlib.h>
3 #include <string.h>
4 #include <assert.h>
5 
6 #include "growing.h"
7 
8 double **P, **cvxHull, **punti_bordo;
9 
10 /*----------------------------------------------------------------------------------------------------*/
regGrow8(struct Cell_head Elaboration,struct element_grow ** mat,double ** punti,int * lung,int r,int c,int v,double Th_j,int maxP)11 void regGrow8(struct Cell_head Elaboration, struct element_grow **mat,
12 	      double **punti, int *lung, int r, int c, int v, double Th_j,
13 	      int maxP)
14 {
15     extern int count_obj;
16 
17     mat[r][c].clas = v;
18     mat[r][c].obj = count_obj;
19 
20     punti[*lung][0] = c;
21     punti[*lung][1] = r;
22     punti[*lung][2] = mat[r][c].interp;
23 
24     assert((*lung)++ < maxP - 1);	/* Condition to finish regGrow8 */
25 
26     if (r - 1 >= 0) {
27 	if ((mat[r - 1][c].clas > Th_j) && (mat[r - 1][c].clas < v) &&
28 	    (mat[r - 1][c].fi != 0))
29 	    regGrow8(Elaboration, mat, punti, lung, r - 1, c, v, Th_j, maxP);
30     }
31 
32     if (c - 1 >= 0) {
33 	if ((mat[r][c - 1].clas > Th_j) && (mat[r][c - 1].clas < v) &&
34 	    (mat[r][c - 1].fi != 0))
35 	    regGrow8(Elaboration, mat, punti, lung, r, c - 1, v, Th_j, maxP);
36     }
37 
38     if (c + 1 < Elaboration.cols) {
39 	if ((mat[r][c + 1].clas > Th_j) && (mat[r][c + 1].clas < v) &&
40 	    (mat[r][c + 1].fi != 0))
41 	    regGrow8(Elaboration, mat, punti, lung, r, c + 1, v, Th_j, maxP);
42     }
43 
44     if (r + 1 < Elaboration.rows) {
45 	if ((mat[r + 1][c].clas > Th_j) && (mat[r + 1][c].clas < v) &&
46 	    (mat[r + 1][c].fi != 0))
47 	    regGrow8(Elaboration, mat, punti, lung, r + 1, c, v, Th_j, maxP);
48     }
49 
50     if ((r - 1 >= 0) && (c - 1 >= 0)) {
51 	if ((mat[r - 1][c - 1].clas > Th_j) && (mat[r - 1][c - 1].clas < v) &&
52 	    (mat[r - 1][c - 1].fi != 0))
53 	    regGrow8(Elaboration, mat, punti, lung, r - 1, c - 1, v, Th_j,
54 		     maxP);
55     }
56 
57     if ((r - 1 >= 0) && (c + 1 < Elaboration.cols)) {
58 	if ((mat[r - 1][c + 1].clas > Th_j) && (mat[r - 1][c + 1].clas < v) &&
59 	    (mat[r - 1][c + 1].fi != 0))
60 	    regGrow8(Elaboration, mat, punti, lung, r - 1, c + 1, v, Th_j,
61 		     maxP);
62     }
63 
64     if ((r + 1 < Elaboration.rows) && (c - 1 >= 0)) {
65 	if ((mat[r + 1][c - 1].clas > Th_j) && (mat[r + 1][c - 1].clas < v) &&
66 	    (mat[r + 1][c - 1].fi != 0))
67 	    regGrow8(Elaboration, mat, punti, lung, r + 1, c - 1, v, Th_j,
68 		     maxP);
69     }
70 
71     if ((r + 1 < Elaboration.rows) && (c + 1 < Elaboration.cols)) {
72 	if ((mat[r + 1][c + 1].clas > Th_j) && (mat[r + 1][c + 1].clas < v) &&
73 	    (mat[r + 1][c + 1].fi != 0))
74 	    regGrow8(Elaboration, mat, punti, lung, r + 1, c + 1, v, Th_j,
75 		     maxP);
76     }
77 }
78 
ccw(double ** P,int i,int j,int k)79 int ccw(double **P, int i, int j, int k)
80 /* true if points i, j, k counterclockwise */
81 {
82     double a, b, c, d;
83 
84     /* It compares coord differences */
85     a = P[i][0] - P[j][0];
86     b = P[i][1] - P[j][1];
87     c = P[k][0] - P[j][0];
88     d = P[k][1] - P[j][1];
89     return a * d - b * c <= 0;
90 }
91 
92 
cmpl(const void * a,const void * b)93 int cmpl(const void *a, const void *b)
94 {
95     double v;
96 
97     CMPM(0, a, b);
98     CMPM(1, b, a);
99     return 0;
100 }
101 
cmph(const void * a,const void * b)102 int cmph(const void *a, const void *b)
103 {
104     return cmpl(b, a);
105 }
106 
make_chain(double ** V,int n,int (* cmp)(const void *,const void *))107 int make_chain(double **V, int n, int (*cmp) (const void *, const void *))
108 {
109     int i, j, s = 1;
110     double *t;
111 
112     qsort(V, (size_t) n, sizeof(double *), cmp);	/* It sorts the V's index */
113 
114     /* */
115     for (i = 2; i < n; i++) {
116 	for (j = s; j >= 1 && ccw(V, i, j, j - 1); j--) ;
117 	s = j + 1;
118 	t = V[s];
119 	V[s] = V[i];
120 	V[i] = t;
121     }
122     return s;
123 }
124 
ch2d(double ** P,int n)125 int ch2d(double **P, int n)
126 {
127     int u = make_chain(P, n, cmpl);	/* make lower hull */
128 
129     if (!n)
130 	return 0;
131     P[n] = P[0];
132 
133     return u + make_chain(P + u, n - u + 1, cmph);	/* make upper hull */
134 }
135 
print_hull(double ** P,double ** pti,int m,double ** h)136 void print_hull(double **P, double **pti, int m, double **h)
137 {
138     int i;
139 
140     for (i = 0; i < m; i++) {
141 	h[i][0] = pti[(P[i] - pti[0]) / 2][0];
142 	h[i][1] = pti[(P[i] - pti[0]) / 2][1];
143 	h[i][2] = pti[(P[i] - pti[0]) / 2][2];
144     }
145 }
146 
checkHull(int cR,int cC,double ** oldHull,int lungOld)147 int checkHull(int cR, int cC, double **oldHull, int lungOld)
148 {
149     double **newP;
150     double **newPoint;
151     int count, lungHullNew;
152 
153     newP = Pvector(0, lungOld + 1);
154     newPoint = G_alloc_matrix(lungOld + 1, 2);
155 
156     for (count = 0; count < lungOld; count++) {
157 	newPoint[count][0] = oldHull[count][0];
158 	newPoint[count][1] = oldHull[count][1];
159 	newP[count] = newPoint[count];
160     }
161 
162     newPoint[lungOld][0] = cC;
163     newPoint[lungOld][1] = cR;
164 
165     newP[lungOld] = newPoint[lungOld];
166 
167     lungHullNew = ch2d(newP, lungOld + 1);
168 
169     if (lungOld != lungHullNew) {
170 	G_free_matrix(newPoint);
171 	free_Pvector(newP, 0, lungOld + 1);
172 	return 0;
173     }
174     else {
175 	for (count = 0; count < lungOld; count++) {
176 	    if ((oldHull[count][0] != newP[count][0]) ||
177 		(oldHull[count][1] != newP[count][1])) {
178 		G_free_matrix(newPoint);
179 		free_Pvector(newP, 0, lungOld + 1);
180 		return 0;
181 	    }
182 	}
183     }
184     G_free_matrix(newPoint);
185     free_Pvector(newP, 0, lungOld + 1);
186     return 1;
187 }
188 
189 /*---------------------------------------------------------------------------------------------------*/
pianOriz(double ** punti,int obsNum,double * minNS,double * minEW,double * maxNS,double * maxEW,struct element_grow ** mat,int CBordo)190 double pianOriz(double **punti, int obsNum, double *minNS, double *minEW,
191 		double *maxNS, double *maxEW, struct element_grow **mat,
192 		int CBordo)
193 {
194     int c1;
195     double minBordo, medioBordo;	/*, minBordo1; */
196 
197     /*Calcola coordinate min e max della zona e media delle righe e delle colonne */
198     *minNS = punti[0][1];
199     *maxNS = punti[0][1];
200     *minEW = punti[0][0];
201     *maxEW = punti[0][0];
202 
203     medioBordo = 0;
204 
205     minBordo = punti[0][2];
206     /*minBordo1 = punti[0][2]; */
207 
208     for (c1 = 0; c1 < obsNum; c1++) {
209 	if (punti[c1][0] > *maxEW)
210 	    *maxEW = punti[c1][0];
211 
212 	if (punti[c1][0] < *minEW)
213 	    *minEW = punti[c1][0];
214 
215 	if (punti[c1][1] > *maxNS)
216 	    *maxNS = punti[c1][1];
217 
218 	if (punti[c1][1] < *minNS)
219 	    *minNS = punti[c1][1];
220 	/*
221 	   if ((punti[c1][2] < minBordo1) && (mat[(int)(punti[c1][1])][(int)(punti[c1][0])].clas >= 1)
222 	   && (mat[(int)(punti[c1][1])][(int)(punti[c1][0])].clas < CBordo)) {
223 	   minBordo1 = punti[c1][2];
224 	   }
225 	 */
226 	if (punti[c1][2] < minBordo)
227 	    minBordo = punti[c1][2];
228 
229 	medioBordo += punti[c1][2];
230     }
231     medioBordo /= obsNum;
232     return medioBordo;
233 }
234 
235 /*----------------------------------------------------------------------------------------------*/
Pvector(long nl,long nh)236 double **Pvector(long nl, long nh)
237 {
238     double **v;
239 
240     v = (double **)calloc((size_t) (nh - nl + 1 + NR_END), sizeof(double *));
241     if (!v)
242 	nrerror("allocation failure in dvector()");
243 
244     return v - nl + NR_END;
245 }
246 
P_alloc_element(int rows,int cols)247 struct element_grow **P_alloc_element(int rows, int cols)
248 {
249     struct element_grow **m;
250     int i;
251 
252     m = (struct element_grow **)G_calloc((rows + 1),
253 					 sizeof(struct element_grow *));
254     m[0] =
255 	(struct element_grow *)G_calloc(rows * cols,
256 					sizeof(struct element_grow));
257 
258     for (i = 1; i <= rows; i++)
259 	m[i] = m[i - 1] + cols;
260 
261     return m;
262 }
263 
nrerror(char error_text[])264 void nrerror(char error_text[])
265 /* standard error handler */
266 {
267     G_debug(1, "run-time error...");
268     G_debug(1, "%s", error_text);
269     G_fatal_error(_("...now exiting to system..."));
270     exit(EXIT_FAILURE);
271 }
272 
structMatrix(long nrl,long nrh,long ncl,long nch)273 struct element_grow **structMatrix(long nrl, long nrh, long ncl, long nch)
274 {
275     long i, nrow = nrh - nrl + 1, ncol = nch - ncl + 1;
276     struct element_grow **m;
277 
278     /* allocate pointers to rows */
279     m = (struct element_grow **)calloc((size_t) (nrow + NR_END),
280 				       sizeof(struct element_grow *));
281     if (!m)
282 	nrerror("allocation failure 1 in matrix()");
283 
284     m += NR_END;
285     m -= nrl;
286 
287     /* allocate rows and set pointers to them */
288     m[nrl] =
289 	(struct element_grow *)calloc((size_t) (nrow * ncol + NR_END),
290 				      sizeof(struct element_grow));
291     if (!m[nrl])
292 	nrerror("allocation failure 2 in matrix()");
293 
294     m[nrl] += NR_END;
295     m[nrl] -= ncl;
296 
297     for (i = nrl + 1; i <= nrh; i++)
298 	m[i] = m[i - 1] + ncol;
299 
300     /* return pointer to array of pointers to rows */
301     return m;
302 }
303 
free_Pvector(double ** v,long nl,long nh)304 void free_Pvector(double **v, long nl, long nh)
305 {
306 
307     free((FREE_ARG) (v + nl - NR_END));
308 }
309 
free_structmatrix(struct element_grow ** m,long nrl,long nrh,long ncl,long nch)310 void free_structmatrix(struct element_grow **m, long nrl, long nrh, long ncl,
311 		       long nch)
312 {
313     free((FREE_ARG) (m[nrl] + ncl - NR_END));
314     free((FREE_ARG) (m + nrl - NR_END));
315 }
316