[GRASS-SVN] r68791 - sandbox/moritz/r.object.geometry

[svn_grass at osgeo.org]

Author: mlennert
Date: 2016-06-28 05:56:33 -0700 (Tue, 28 Jun 2016)
New Revision: 68791

Modified:
   sandbox/moritz/r.object.geometry/main.c
Log:
r.object.geometry: add circular compactness and fractal dimension; add header to output


Modified: sandbox/moritz/r.object.geometry/main.c
===================================================================
--- sandbox/moritz/r.object.geometry/main.c	2016-06-28 09:36:38 UTC (rev 68790)
+++ sandbox/moritz/r.object.geometry/main.c	2016-06-28 12:56:33 UTC (rev 68791)
@@ -120,6 +120,7 @@
     if (Rast_is_c_null_value(&min) || Rast_is_c_null_value(&max))
 	G_fatal_error(_("Empty input map <%s>"), opt_in->answer);
     
+    /* REMARK: The following is only true if object ids are continuously numbered */
     n_objects = max - min + 1;
     obj_geos = G_malloc(n_objects * sizeof(struct obj_geo));
     
@@ -148,6 +149,7 @@
     compute_areas = flag_m->answer && !planimetric;
     G_begin_distance_calculations();
 
+    G_message(_("Calculating statistics"));
     for (row = 0; row < nrows; row++) {
 	G_percent(row, nrows + 1, 2);
 
@@ -279,13 +281,21 @@
     G_free(cur_in);
     G_free(prev_in);
 
+    G_message(_("Writing output"));
     /* print table */
     /* TODO: use separator */
     /* print table header (column names) */
+    fprintf(out_fp, "cat%s", sep);
+    fprintf(out_fp, "area%s", sep);
+    fprintf(out_fp, "perimeter%s", sep);
+    fprintf(out_fp, "compact_square%s", sep);
+    fprintf(out_fp, "compact_circle%s", sep);
+    fprintf(out_fp, "fd", sep);
     fprintf(out_fp, "\n");
 
     /* print table body */
     for (i = 0; i < n_objects; i++) {
+	 G_percent(i, n_objects - 1, 1);
 	/* skip emtpy objects */
 	if (obj_geos[i].area == 0)
 	    continue;
@@ -293,7 +303,12 @@
 	fprintf(out_fp, "%d%s", min+i, sep);
 	fprintf(out_fp, "%f%s", obj_geos[i].area, sep);
 	fprintf(out_fp, "%f%s", obj_geos[i].perimeter, sep);
-	fprintf(out_fp, "%f", 4 * sqrt(obj_geos[i].area) / obj_geos[i].perimeter);
+	fprintf(out_fp, "%f%s", 4 * sqrt(obj_geos[i].area) / obj_geos[i].perimeter, sep);
+	fprintf(out_fp, "%f%s", obj_geos[i].perimeter / (2 * sqrt(M_PI * obj_geos[i].area)), sep);
+	/* log 1 = 0, so avoid that: */
+	if (obj_geos[i].area == 1)
+		obj_geos[i].area += 0.000001;
+	fprintf(out_fp, "%f", 2 * log(obj_geos[i].perimeter) / log (obj_geos[i].area));
 	/* TODO: convert bounding box area and perimeter to (square) meters for -m flag */