1 /****************************************************************************
2  *
3  * MODULE:	r.li.ritchness
4  *
5  * PURPOSE:	calculates richness diversity index
6  *
7  * AUTHOR(S):	Serena Pallecchi student of Computer Science University of Pisa (Italy)
8  *              Commission from Faunalia Pontedera (PI) www.faunalia.it
9  *              Rewrite: Markus Metz
10  *
11  * COPYRIGHT:
12  *		This program is free software under the GPL (>=v2)
13  *		Read the COPYING file that comes with GRASS for details.
14  *
15  ***************************************************************************/
16 
17 #include <stdlib.h>
18 #include <fcntl.h>
19 #include <math.h>
20 
21 #include <grass/gis.h>
22 #include <grass/raster.h>
23 #include <grass/glocale.h>
24 
25 #include "../r.li.daemon/daemon.h"
26 #include "../r.li.daemon/avlDefs.h"
27 #include "../r.li.daemon/avl.h"
28 
29 rli_func richness;
30 int calculate(int fd, struct area_entry *ad, double *result);
31 int calculateD(int fd, struct area_entry *ad, double *result);
32 int calculateF(int fd, struct area_entry *ad, double *result);
33 
main(int argc,char * argv[])34 int main(int argc, char *argv[])
35 {
36     struct Option *raster, *conf, *output;
37     struct GModule *module;
38 
39     G_gisinit(argv[0]);
40     module = G_define_module();
41     module->description =
42 	_("Calculates richness index on a raster map");
43     G_add_keyword(_("raster"));
44     G_add_keyword(_("landscape structure analysis"));
45     G_add_keyword(_("diversity index"));
46 
47     /* define options */
48 
49     raster = G_define_standard_option(G_OPT_R_INPUT);
50 
51     conf = G_define_standard_option(G_OPT_F_INPUT);
52     conf->key = "config";
53     conf->description = _("Configuration file");
54     conf->required = YES;
55 
56     output = G_define_standard_option(G_OPT_R_OUTPUT);
57 
58     if (G_parser(argc, argv))
59 	exit(EXIT_FAILURE);
60 
61     return calculateIndex(conf->answer, richness, NULL, raster->answer,
62 			  output->answer);
63 }
64 
richness(int fd,char ** par,struct area_entry * ad,double * result)65 int richness(int fd, char **par, struct area_entry *ad, double *result)
66 {
67     int ris = RLI_OK;
68     double indice = 0;
69 
70     switch (ad->data_type) {
71     case CELL_TYPE:
72 	{
73 	    ris = calculate(fd, ad, &indice);
74 	    break;
75 	}
76     case DCELL_TYPE:
77 	{
78 	    ris = calculateD(fd, ad, &indice);
79 	    break;
80 	}
81     case FCELL_TYPE:
82 	{
83 	    ris = calculateF(fd, ad, &indice);
84 	    break;
85 	}
86     default:
87 	{
88 	    G_fatal_error("data type unknown");
89 	    return RLI_ERRORE;
90 	}
91     }
92 
93     if (ris != RLI_OK) {
94 	return RLI_ERRORE;
95     }
96 
97     *result = indice;
98 
99     return RLI_OK;
100 }
101 
102 
calculate(int fd,struct area_entry * ad,double * result)103 int calculate(int fd, struct area_entry *ad, double *result)
104 {
105     CELL *buf;
106     CELL corrCell;
107     CELL precCell;
108 
109     int i, j;
110     int mask_fd = -1, *mask_buf = NULL;
111     int ris = 0;
112     int masked = FALSE;
113     long m = 0;
114     long area = 0;
115     avl_tree albero = NULL;
116     generic_cell uc;
117 
118     uc.t = CELL_TYPE;
119 
120     /* open mask if needed */
121     if (ad->mask == 1) {
122 	if ((mask_fd = open(ad->mask_name, O_RDONLY, 0755)) < 0)
123 	    return RLI_ERRORE;
124 	mask_buf = G_malloc(ad->cl * sizeof(int));
125 	if (mask_buf == NULL) {
126 	    G_fatal_error("malloc mask_buf failed");
127 	    return RLI_ERRORE;
128 	}
129 	masked = TRUE;
130     }
131 
132     Rast_set_c_null_value(&precCell, 1);
133     for (j = 0; j < ad->rl; j++) {	/* for each row */
134 	if (masked) {
135 	    if (read(mask_fd, mask_buf, (ad->cl * sizeof(int))) < 0) {
136 		G_fatal_error("mask read failed");
137 		return RLI_ERRORE;
138 	    }
139 	}
140 
141 	buf = RLI_get_cell_raster_row(fd, j + ad->y, ad);
142 
143 	for (i = 0; i < ad->cl; i++) {	/* for each cell in the row */
144 	    corrCell = buf[i + ad->x];
145 
146 	    if ((masked) && (mask_buf[i] == 0)) {
147 		Rast_set_c_null_value(&corrCell, 1);
148 	    }
149 	    else {
150 		/* total sample area */
151 		area++;
152 	    }
153 
154 	    if (!(Rast_is_null_value(&precCell, uc.t)) &&
155 		 corrCell != precCell) {
156 
157 		if (albero == NULL) {
158 		    uc.val.c = precCell;
159 		    albero = avl_make(uc, (long)1);
160 		    if (albero == NULL) {
161 			G_fatal_error("avl_make error");
162 			return RLI_ERRORE;
163 		    }
164 		    m++;
165 		}
166 		else {
167 		    uc.val.c = precCell;
168 		    ris = avl_add(&albero, uc, (long)1);
169 		    switch (ris) {
170 		    case AVL_ERR:
171 			{
172 			    G_fatal_error("avl_add error");
173 			    return RLI_ERRORE;
174 			}
175 		    case AVL_ADD:
176 			{
177 			    m++;
178 			    break;
179 			}
180 		    case AVL_PRES:
181 			{
182 			    break;
183 			}
184 		    default:
185 			{
186 			    G_fatal_error("avl_make unknown error");
187 			    return RLI_ERRORE;
188 			}
189 		    }
190 		}
191 	    }		/* endif not equal cells */
192 	    precCell = corrCell;
193 	}
194     }
195 
196     /* last closing */
197     if (!(Rast_is_null_value(&precCell, uc.t))) {
198 	if (albero == NULL) {
199 	    uc.val.c = precCell;
200 	    albero = avl_make(uc, (long)1);
201 	    if (albero == NULL) {
202 		G_fatal_error("avl_make error");
203 		return RLI_ERRORE;
204 	    }
205 	    m++;
206 	}
207 	else {
208 	    uc.val.c = precCell;
209 	    ris = avl_add(&albero, uc, (long)1);
210 	    switch (ris) {
211 	    case AVL_ERR:
212 		{
213 		    G_fatal_error("avl_add error");
214 		    return RLI_ERRORE;
215 		}
216 	    case AVL_ADD:
217 		{
218 		    m++;
219 		    break;
220 		}
221 	    case AVL_PRES:
222 		{
223 		    break;
224 		}
225 	    default:
226 		{
227 		    G_fatal_error("avl_add unknown error");
228 		    return RLI_ERRORE;
229 		}
230 	    }
231 	}
232     }
233 
234     if (area)
235 	*result = m;
236     else
237 	Rast_set_d_null_value(result, 1);
238 
239     avl_destroy(albero);
240     if (masked) {
241 	close(mask_fd);
242 	G_free(mask_buf);
243     }
244 
245     return RLI_OK;
246 }
247 
248 
calculateD(int fd,struct area_entry * ad,double * result)249 int calculateD(int fd, struct area_entry *ad, double *result)
250 {
251     DCELL *buf;
252     DCELL corrCell;
253     DCELL precCell;
254 
255     int i, j;
256     int mask_fd = -1, *mask_buf = NULL;
257     int ris = 0;
258     int masked = FALSE;
259     long m = 0;
260     long area = 0;
261     avl_tree albero = NULL;
262     generic_cell uc;
263 
264     uc.t = DCELL_TYPE;
265 
266     /* open mask if needed */
267     if (ad->mask == 1) {
268 	if ((mask_fd = open(ad->mask_name, O_RDONLY, 0755)) < 0)
269 	    return RLI_ERRORE;
270 	mask_buf = G_malloc(ad->cl * sizeof(int));
271 	if (mask_buf == NULL) {
272 	    G_fatal_error("malloc mask_buf failed");
273 	    return RLI_ERRORE;
274 	}
275 	masked = TRUE;
276     }
277 
278     Rast_set_d_null_value(&precCell, 1);
279     for (j = 0; j < ad->rl; j++) {	/* for each row */
280 	if (masked) {
281 	    if (read(mask_fd, mask_buf, (ad->cl * sizeof(int))) < 0) {
282 		G_fatal_error("mask read failed");
283 		return RLI_ERRORE;
284 	    }
285 	}
286 
287 	buf = RLI_get_dcell_raster_row(fd, j + ad->y, ad);
288 
289 	for (i = 0; i < ad->cl; i++) {	/* for each dcell in the row */
290 	    corrCell = buf[i + ad->x];
291 
292 	    if ((masked) && (mask_buf[i] == 0)) {
293 		Rast_set_d_null_value(&corrCell, 1);
294 	    }
295 	    else {
296 		/* total sample area */
297 		area++;
298 	    }
299 
300 	    if (!(Rast_is_null_value(&precCell, uc.t)) &&
301 		 corrCell != precCell) {
302 
303 		if (albero == NULL) {
304 		    uc.val.dc = precCell;
305 		    albero = avl_make(uc, (long)1);
306 		    if (albero == NULL) {
307 			G_fatal_error("avl_make error");
308 			return RLI_ERRORE;
309 		    }
310 		    m++;
311 		}
312 		else {
313 		    uc.val.dc = precCell;
314 		    ris = avl_add(&albero, uc, (long)1);
315 		    switch (ris) {
316 		    case AVL_ERR:
317 			{
318 			    G_fatal_error("avl_add error");
319 			    return RLI_ERRORE;
320 			}
321 		    case AVL_ADD:
322 			{
323 			    m++;
324 			    break;
325 			}
326 		    case AVL_PRES:
327 			{
328 			    break;
329 			}
330 		    default:
331 			{
332 			    G_fatal_error("avl_make unknown error");
333 			    return RLI_ERRORE;
334 			}
335 		    }
336 		}
337 	    }		/* endif not equal dcells */
338 	    precCell = corrCell;
339 	}
340     }
341 
342     /* last closing */
343     if (!(Rast_is_null_value(&precCell, uc.t))) {
344 	if (albero == NULL) {
345 	    uc.val.dc = precCell;
346 	    albero = avl_make(uc, (long)1);
347 	    if (albero == NULL) {
348 		G_fatal_error("avl_make error");
349 		return RLI_ERRORE;
350 	    }
351 	    m++;
352 	}
353 	else {
354 	    uc.val.dc = precCell;
355 	    ris = avl_add(&albero, uc, (long)1);
356 	    switch (ris) {
357 	    case AVL_ERR:
358 		{
359 		    G_fatal_error("avl_add error");
360 		    return RLI_ERRORE;
361 		}
362 	    case AVL_ADD:
363 		{
364 		    m++;
365 		    break;
366 		}
367 	    case AVL_PRES:
368 		{
369 		    break;
370 		}
371 	    default:
372 		{
373 		    G_fatal_error("avl_add unknown error");
374 		    return RLI_ERRORE;
375 		}
376 	    }
377 	}
378     }
379 
380     if (area)
381 	*result = m;
382     else
383 	Rast_set_d_null_value(result, 1);
384 
385     avl_destroy(albero);
386     if (masked) {
387 	close(mask_fd);
388 	G_free(mask_buf);
389     }
390 
391     return RLI_OK;
392 }
393 
394 
calculateF(int fd,struct area_entry * ad,double * result)395 int calculateF(int fd, struct area_entry *ad, double *result)
396 {
397     FCELL *buf;
398     FCELL corrCell;
399     FCELL precCell;
400 
401     int i, j;
402     int mask_fd = -1, *mask_buf = NULL;
403     int ris = 0;
404     int masked = FALSE;
405     long m = 0;
406     long area = 0;
407     avl_tree albero = NULL;
408     generic_cell uc;
409 
410     uc.t = FCELL_TYPE;
411 
412     /* open mask if needed */
413     if (ad->mask == 1) {
414 	if ((mask_fd = open(ad->mask_name, O_RDONLY, 0755)) < 0)
415 	    return RLI_ERRORE;
416 	mask_buf = G_malloc(ad->cl * sizeof(int));
417 	if (mask_buf == NULL) {
418 	    G_fatal_error("malloc mask_buf failed");
419 	    return RLI_ERRORE;
420 	}
421 	masked = TRUE;
422     }
423 
424     Rast_set_f_null_value(&precCell, 1);
425     for (j = 0; j < ad->rl; j++) {	/* for each row */
426 	if (masked) {
427 	    if (read(mask_fd, mask_buf, (ad->cl * sizeof(int))) < 0) {
428 		G_fatal_error("mask read failed");
429 		return RLI_ERRORE;
430 	    }
431 	}
432 
433 	buf = RLI_get_fcell_raster_row(fd, j + ad->y, ad);
434 
435 	for (i = 0; i < ad->cl; i++) {	/* for each fcell in the row */
436 
437 	    corrCell = buf[i + ad->x];
438 
439 	    if ((masked) && (mask_buf[i] == 0)) {
440 		Rast_set_f_null_value(&corrCell, 1);
441 	    }
442 	    else {
443 		/* total sample area */
444 		area++;
445 	    }
446 
447 	    if (!(Rast_is_null_value(&precCell, uc.t)) &&
448 		 corrCell != precCell) {
449 
450 		if (albero == NULL) {
451 		    uc.val.fc = precCell;
452 		    albero = avl_make(uc, (long)1);
453 		    if (albero == NULL) {
454 			G_fatal_error("avl_make error");
455 			return RLI_ERRORE;
456 		    }
457 		    m++;
458 		}
459 		else {
460 		    uc.val.fc = precCell;
461 		    ris = avl_add(&albero, uc, (long)1);
462 		    switch (ris) {
463 		    case AVL_ERR:
464 			{
465 			    G_fatal_error("avl_add error");
466 			    return RLI_ERRORE;
467 			}
468 		    case AVL_ADD:
469 			{
470 			    m++;
471 			    break;
472 			}
473 		    case AVL_PRES:
474 			{
475 			    break;
476 			}
477 		    default:
478 			{
479 			    G_fatal_error("avl_make unknown error");
480 			    return RLI_ERRORE;
481 			}
482 		    }
483 		}
484 
485 	    }		/* endif not equal fcells */
486 	    precCell = corrCell;
487 	}
488     }
489 
490     /* last closing */
491     if (!(Rast_is_null_value(&precCell, uc.t))) {
492 	if (albero == NULL) {
493 	    uc.val.fc = precCell;
494 	    albero = avl_make(uc, (long)1);
495 	    if (albero == NULL) {
496 		G_fatal_error("avl_make error");
497 		return RLI_ERRORE;
498 	    }
499 	    m++;
500 	}
501 	else {
502 	    uc.val.fc = precCell;
503 	    ris = avl_add(&albero, uc, (long)1);
504 	    switch (ris) {
505 	    case AVL_ERR:
506 		{
507 		    G_fatal_error("avl_add error");
508 		    return RLI_ERRORE;
509 		}
510 	    case AVL_ADD:
511 		{
512 		    m++;
513 		    break;
514 		}
515 	    case AVL_PRES:
516 		{
517 		    break;
518 		}
519 	    default:
520 		{
521 		    G_fatal_error("avl_add unknown error");
522 		    return RLI_ERRORE;
523 		}
524 	    }
525 	}
526     }
527 
528     if (area)
529 	*result = m;
530     else
531 	Rast_set_d_null_value(result, 1);
532 
533     avl_destroy(albero);
534     if (masked) {
535 	close(mask_fd);
536 	G_free(mask_buf);
537     }
538 
539     return RLI_OK;
540 }
541