[GRASS-SVN] r42949 - in grass-addons/raster/r.pi: . r.pi.index
svn_grass at osgeo.org
svn_grass at osgeo.org
Thu Jul 29 12:21:45 EDT 2010
Author: wegmann
Date: 2010-07-29 16:21:45 +0000 (Thu, 29 Jul 2010)
New Revision: 42949
Added:
grass-addons/raster/r.pi/r.pi.index/
grass-addons/raster/r.pi/r.pi.index/Makefile
grass-addons/raster/r.pi/r.pi.index/description.html
grass-addons/raster/r.pi/r.pi.index/frag.c
grass-addons/raster/r.pi/r.pi.index/func.c
grass-addons/raster/r.pi/r.pi.index/local_proto.h
grass-addons/raster/r.pi/r.pi.index/main.c
Log:
initial commit: r.pi.index provides basic patch based spatial indices like area, shape, euclidean distance to nearest neighbour etc.
Added: grass-addons/raster/r.pi/r.pi.index/Makefile
===================================================================
--- grass-addons/raster/r.pi/r.pi.index/Makefile (rev 0)
+++ grass-addons/raster/r.pi/r.pi.index/Makefile 2010-07-29 16:21:45 UTC (rev 42949)
@@ -0,0 +1,10 @@
+MODULE_TOPDIR = ../..
+
+PGM = r.pi.index
+
+LIBES = $(STATSLIB) $(GISLIB)
+DEPENDENCIES = $(STATSDEP) $(GISDEP)
+
+include $(MODULE_TOPDIR)/include/Make/Module.make
+
+default: cmd
Property changes on: grass-addons/raster/r.pi/r.pi.index/Makefile
___________________________________________________________________
Added: svn:executable
+ *
Added: grass-addons/raster/r.pi/r.pi.index/description.html
===================================================================
--- grass-addons/raster/r.pi/r.pi.index/description.html (rev 0)
+++ grass-addons/raster/r.pi/r.pi.index/description.html 2010-07-29 16:21:45 UTC (rev 42949)
@@ -0,0 +1,40 @@
+<H2>DESCRIPTION</H2>
+
+Available options for the index to be computed for patches within a certain class are: area (area), perimeter (perim), SHAPE (shape), Border-Index (bor), Compactness (comp), Asymmetry (asym), area-perimeter ratio (apr), fractal dimension (fract), distance to euclidean nearest neighbour (ENN)</DD>
+
+<P>
+
+The program will be run non-interactively if the user specifies program arguments (see OPTIONS) on the command
+line. Alternately, the user can simply type <B>r.pi.index</B> on the command line, without program
+arguments. In this case, the user will be prompted for flag settings and parameter values.
+
+
+<H2>NOTES</H2>
+
+The <EM>Nearest Neighbour Index</EM> (ENN) analyse the Euclidean Nearest
+Neighbour to the first neighbouring patch. The output value is in pixel and can
+be converted to a distance values using g.region resolution information.
+<EM>r.pi.ENN</EM> and <EM>r.pi.FNN</EM> provide the same analysis concerning the
+ first nearest neighbour (NN), but are extended to the n-th NN. However due to
+code construction does the <EM>r.pi.index</EM> distance analysis to first ENN perform faster.
+
+
+<P>
+
+<H2>BUGS</H2>
+Landscapes with more than 10 000 patches might cause a memory
+allocation error depending on your system.
+
+
+<H2>SEE ALSO</H2>
+
+<EM><A HREF="r.pi.ENN.html">r.pi.ENN</A></EM><br>
+<EM><A HREF="r.pi.FNN.html">r.pi.FNN</A></EM><br>
+<EM><A HREF="r.pi.dist.html">r.pi.dist</A></EM><br>
+<EM><A HREF="r.li.setup.html">r.li.setup</A></EM><br>
+
+<H2>AUTHOR</H2>
+Programming: Elshad Shirinov<br>
+Scientific concept: Dr. Martin Wegmann <br>
+Department of Remote Sensing <br>Remote Sensing and Biodiversity Unit<br> University of Wuerzburg, Germany
+
Property changes on: grass-addons/raster/r.pi/r.pi.index/description.html
___________________________________________________________________
Added: svn:executable
+ *
Added: grass-addons/raster/r.pi/r.pi.index/frag.c
===================================================================
--- grass-addons/raster/r.pi/r.pi.index/frag.c (rev 0)
+++ grass-addons/raster/r.pi/r.pi.index/frag.c 2010-07-29 16:21:45 UTC (rev 42949)
@@ -0,0 +1,108 @@
+#include <string.h>
+#include <stdlib.h>
+#include <unistd.h>
+#include <grass/gis.h>
+#include <grass/glocale.h>
+#include <grass/stats.h>
+#include <math.h>
+#include "local_proto.h"
+
+typedef struct {
+ int x, y;
+ } Position;
+
+int getNeighbors(Position* res, int x, int y, int nx, int ny, int nbr_cnt);
+
+int getNeighbors(Position* res, int x, int y, int nx, int ny, int nbr_cnt) {
+ int left, right, top, bottom;
+ int i, j;
+ int cnt = 0;
+
+ switch(nbr_cnt) {
+ case 4: /* von Neumann neighborhood */
+ if (x > 0 && flagbuf[y * nx + x - 1] == 1) {
+ res[cnt].x = x - 1; res[cnt].y = y; cnt++;
+ }
+ if (y > 0 && flagbuf[(y - 1) * nx + x] == 1) {
+ res[cnt].x = x; res[cnt].y = y - 1; cnt++;
+ }
+ if (x < nx - 1 && flagbuf[y * nx + x + 1] == 1) {
+ res[cnt].x = x + 1; res[cnt].y = y; cnt++;
+ }
+ if (y < ny - 1 && flagbuf[(y + 1) * nx + x] == 1) {
+ res[cnt].x = x; res[cnt].y = y + 1; cnt++;
+ }
+ break;
+ case 8: /* Moore neighborhood */
+ left = x > 0 ? x - 1 : 0;
+ top = y > 0 ? y - 1 : 0;
+ right = x < nx - 1 ? x + 1 : nx - 1;
+ bottom = y < ny - 1 ? y + 1 : ny - 1;
+ for (i = left; i <= right; i++) {
+ for (j = top; j <= bottom; j++) {
+ if(!(i == x && j == y) && flagbuf[j * nx + i] == 1) {
+ res[cnt].x = i; res[cnt].y = j; cnt++;
+ }
+ }
+ }
+ break;
+ }
+
+ return cnt;
+}
+
+void writeFrag(int row, int col, int nbr_cnt) {
+ int x, y, i;
+ Position* list = (Position *) G_malloc(nrows * ncols * sizeof (Position));
+ Position* first = list;
+ Position* last = list;
+ Position* nbr_list = (Position *) G_malloc(8 * sizeof (Position));
+
+ /* count neighbors */
+ int neighbors = 0;
+ if (col > 0 && flagbuf[row * ncols + col - 1] != 0) neighbors++;
+ if (row > 0 && flagbuf[(row - 1) * ncols + col] != 0) neighbors++;
+ if (col < ncols - 1 && flagbuf[row * ncols + col + 1] != 0) neighbors++;
+ if (row < nrows - 1 && flagbuf[(row + 1) * ncols + col] != 0) neighbors++;
+
+ /* write first cell */
+ actpos->x = col; actpos->y = row; actpos->neighbors = neighbors; actpos++;
+ flagbuf[row * ncols + col] = -1;
+
+ /* push position on fifo-list */
+ last->x = col; last->y = row; last++;
+
+ while (first < last) {
+ /* get position from fifo-list */
+ int r = first->y; int c = first->x; first++;
+
+ int left = c > 0 ? c - 1 : 0;
+ int top = r > 0 ? r - 1 : 0;
+ int right = c < ncols - 1 ? c + 1 : ncols - 1;
+ int bottom = r < nrows - 1 ? r + 1 : nrows - 1;
+
+ /* add neighbors to fifo-list */
+ int cnt = getNeighbors(nbr_list, c, r, ncols, nrows, nbr_cnt);
+
+ for(i = 0; i < cnt; i++) {
+ x = nbr_list[i].x; y = nbr_list[i].y;
+
+ /* add position to fifo-list */
+ last->x = x; last->y = y; last++;
+
+ /* count neighbors */
+ neighbors = 0;
+ if (x > 0 && flagbuf[y * ncols + x - 1] != 0) neighbors++;
+ if (y > 0 && flagbuf[(y - 1) * ncols + x] != 0) neighbors++;
+ if (x < ncols - 1 && flagbuf[y * ncols + x + 1] != 0) neighbors++;
+ if (y < nrows - 1 && flagbuf[(y + 1) * ncols + x] != 0) neighbors++;
+
+ /* set values */
+ actpos->x = x; actpos->y = y; actpos->neighbors = neighbors; actpos++;
+ flagbuf[y * ncols + x] = -1;
+ }
+ }
+
+ G_free(list);
+ G_free(nbr_list);
+}
Added: grass-addons/raster/r.pi/r.pi.index/func.c
===================================================================
--- grass-addons/raster/r.pi/r.pi.index/func.c (rev 0)
+++ grass-addons/raster/r.pi/r.pi.index/func.c 2010-07-29 16:21:45 UTC (rev 42949)
@@ -0,0 +1,264 @@
+#include <string.h>
+#include <stdlib.h>
+#include <unistd.h>
+#include <grass/gis.h>
+#include <grass/glocale.h>
+#include <grass/stats.h>
+#include <math.h>
+#include "local_proto.h"
+
+int f_area(DCELL *vals, Coords **frags, int count) {
+ int i;
+ /* for all patches */
+ for(i = 0; i < count; i++) {
+ // display progress
+ if (verbose)
+ G_percent (i, count, 2);
+
+ vals[i] = (DCELL) (frags[i+1] - frags[i]);
+ }
+
+ return 0;
+}
+
+int f_perim(DCELL *vals, Coords **frags, int count) {
+ Coords *p;
+ int i;
+ /* for all patches */
+ for(i = 0; i < count; i++) {
+ // display progress
+ if (verbose)
+ G_percent (i, count, 2);
+
+ int border = 0;
+ /* for all cells in a patch */
+ for(p = frags[i]; p < frags[i+1]; p++) {
+ border+= 4 - p->neighbors;
+ }
+ vals[i] = (DCELL) border;
+ }
+
+ return 0;
+}
+
+int f_shapeindex(DCELL *vals, Coords **frags, int count) {
+ Coords *p;
+ int i;
+ /* for all patches */
+ for(i = 0; i < count; i++) {
+ // display progress
+ if (verbose)
+ G_percent (i, count, 2);
+
+ int border = 0;
+ int area = (DCELL) (frags[i+1] - frags[i]);
+ /* for all cells in a patch */
+ for(p = frags[i]; p < frags[i+1]; p++) {
+ border+= 4 - p->neighbors;
+ }
+ vals[i] = (DCELL) border / ( 4 * sqrt((DCELL) area) );
+ }
+
+ return 0;
+}
+
+int f_borderindex(DCELL *vals, Coords **frags, int count) {
+ Coords *p;
+ int i;
+ /* for all patches */
+ for(i = 0; i < count; i++) {
+ // display progress
+ if (verbose)
+ G_percent (i, count, 2);
+
+ int border = 0;
+ int maxx, maxy, minx, miny;
+ int l = 0;
+ int w = 0;
+ maxx = minx = frags[i]->x;
+ maxy = miny = frags[i]->y;
+ /* for all cells in a patch */
+ for(p = frags[i]; p < frags[i+1]; p++) {
+ border+= 4 - p->neighbors;
+ maxx = p->x > maxx ? p->x : maxx;
+ minx = p->x < minx ? p->x : minx;
+ maxy = p->y > maxy ? p->y : maxy;
+ miny = p->y < miny ? p->y : miny;
+ }
+ l = maxx - minx + 1; w = maxy - miny + 1;
+ vals[i] = (DCELL) border / ( 2 * ( (DCELL)l + (DCELL)w ) );
+ }
+
+ return 0;
+}
+
+int f_compactness(DCELL *vals, Coords **frags, int count) {
+ Coords *p;
+ int i;
+ /* for all patches */
+ for(i = 0; i < count; i++) {
+ // display progress
+ if (verbose)
+ G_percent (i, count, 2);
+
+ int area = 0;
+ int maxx, maxy, minx, miny;
+ int l = 0;
+ int w = 0;
+ maxx = minx = frags[i]->x;
+ maxy = miny = frags[i]->y;
+ /* for all cells in a patch */
+ for(p = frags[i]; p < frags[i+1]; p++) {
+ area++;
+ maxx = p->x > maxx ? p->x : maxx;
+ minx = p->x < minx ? p->x : minx;
+ maxy = p->y > maxy ? p->y : maxy;
+ miny = p->y < miny ? p->y : miny;
+ }
+ l = maxx - minx + 1; w = maxy - miny + 1;
+ vals[i] = (DCELL)l * (DCELL)w / (DCELL) area;
+ }
+
+ return 0;
+}
+
+int f_asymmetry(DCELL *vals, Coords **frags, int count) {
+ Coords *p;
+ int i;
+ /* for all patches */
+ for(i = 0; i < count; i++) {
+ // display progress
+ if (verbose)
+ G_percent (i, count, 2);
+
+ /* compute variance for x,y and xy in the patch */
+ /* formula: a(x) = sum(x_i), b(x) = sum(x_i²), var(x) = (b(x) - a(x)² / n) / n */
+ /* covar(x*y) = (a(x * y) - a(x) * a(y) / n) /n */
+ int ax, ay, axy;
+ int bx, by;
+ DCELL vx, vy, vxy, vsum, invn;
+ int n = 0;
+ ax = ay = axy = 0; bx = by = 0;
+ /* for all cells in a patch */
+// fprintf(stderr, "\npatch %d: ", i);
+ for(p = frags[i]; p < frags[i+1]; p++, n++) {
+ int x = p->x; int y = p->y; int xy = p->x * p->y;
+ ax += x; ay += y; axy += xy;
+ bx += x * x; by += y * y;
+// fprintf(stderr, "x_%d = %d, y_%d = %d; ", n, x, n, y);
+ }
+ invn = 1.0 / (DCELL)n;
+ vx = ( (DCELL)bx - (DCELL)ax * (DCELL)ax * invn ) * invn;
+ vy = ( (DCELL)by - (DCELL)ay * (DCELL)ay * invn ) * invn;
+ vxy = ( (DCELL)axy - (DCELL)ax * (DCELL)ay * invn ) * invn;
+// fprintf(stderr, " axy = %d, ax = %d, ay = %d, n = %d", axy, ax, ay, n);
+ vsum = vx + vy;
+ vals[i] = 2 * sqrt(0.25 * vsum * vsum + vxy * vxy - vx * vy) / vsum;
+ }
+}
+
+int f_area_perim_ratio(DCELL *vals, Coords **frags, int count) {
+ Coords *p;
+ int i;
+ /* for all patches */
+ for(i = 0; i < count; i++) {
+ // display progress
+ if (verbose)
+ G_percent (i, count, 2);
+
+ int border = 0;
+ int area = (DCELL) (frags[i+1] - frags[i]);
+ /* for all cells in a patch */
+ for(p = frags[i]; p < frags[i+1]; p++) {
+ border+= 4 - p->neighbors;
+ }
+ vals[i] = (DCELL) area / (DCELL) border;
+ }
+
+ return 0;
+}
+
+int f_frac_dim(DCELL *vals, Coords **frags, int count) {
+ Coords *p;
+ int i;
+ /* for all patches */
+ for(i = 0; i < count; i++) {
+ // display progress
+ if (verbose)
+ G_percent (i, count, 2);
+
+ int border = 0;
+ int area = (DCELL) (frags[i+1] - frags[i]);
+ /* for all cells in a patch */
+ for(p = frags[i]; p < frags[i+1]; p++) {
+ border+= 4 - p->neighbors;
+ }
+ vals[i] = 2 * log(0.25 * (DCELL) border) / log((DCELL) area);
+ }
+
+ return 0;
+}
+
+DCELL dist(Coords *p1, Coords *p2) {
+ int x1 = p1->x; int y1 = p1->y;
+ int x2 = p2->x; int y2 = p2->y;
+ int dx = x2 - x1; int dy = y2 - y1;
+ return sqrt(dx * dx + dy * dy);
+}
+
+DCELL min_dist(Coords **frags, int n1, int n2) {
+ Coords *p1, *p2;
+ DCELL min = 1000000.0;
+ // for all cells in the first patch
+ for (p1 = frags[n1]; p1 < frags[n1+1]; p1++) {
+ // if cell at the border
+ if (p1->neighbors < 4) {
+ // for all cells in the second patch
+ for (p2 = frags[n2]; p2 < frags[n2+1]; p2++) {
+ // if cell at the border
+ if (p2->neighbors < 4) {
+ DCELL d = dist(p1, p2);
+ if (d < min) {
+ min = d;
+ }
+ }
+ }
+ }
+ }
+ return min;
+}
+
+int f_nearest_dist(DCELL *vals, Coords **frags, int count) {
+ int i, j;
+ /* for all patches */
+ for (i = 0; i < count; i++) {
+ // display progress
+ if (verbose)
+ G_percent (i, count, 2);
+
+ DCELL min = 1000000.0;
+ for (j = 0; j < count; j++) {
+ if (i != j) {
+ DCELL d = min_dist(frags, i, j);
+ if (d < min) {
+ min = d;
+ }
+ }
+ }
+ vals[i] = min;
+ }
+
+ return 0;
+}
+
+int f_proximity(DCELL *vals, Coords **frags, int count, DCELL param) {
+ int i, j;
+ for (i = 0; i < count; i++) {
+ for (j = 0; j < count; j++) {
+ if (i != j) {
+ }
+ }
+ }
+
+ return 0;
+}
Added: grass-addons/raster/r.pi/r.pi.index/local_proto.h
===================================================================
--- grass-addons/raster/r.pi/r.pi.index/local_proto.h (rev 0)
+++ grass-addons/raster/r.pi/r.pi.index/local_proto.h 2010-07-29 16:21:45 UTC (rev 42949)
@@ -0,0 +1,30 @@
+#ifndef LOCAL_PROTO_H
+#define LOCAL_PROTO_H
+
+#define MAX_DOUBLE 1000000.0
+
+typedef struct {
+ int x, y;
+ int neighbors;
+ } Coords;
+
+extern void writeFrag(int row, int col, int nbr_cnt);
+
+extern int f_area(DCELL *vals, Coords **frags, int);
+extern int f_perim(DCELL *vals, Coords **frags, int);
+extern int f_shapeindex(DCELL*, Coords**, int);
+extern int f_borderindex(DCELL*, Coords**, int);
+extern int f_compactness(DCELL*, Coords**, int);
+extern int f_asymmetry(DCELL*, Coords**, int);
+extern int f_area_perim_ratio(DCELL*, Coords**, int);
+extern int f_frac_dim(DCELL*, Coords**, int);
+extern int f_nearest_dist(DCELL*, Coords**, int);
+
+/* global variables */
+int nrows, ncols;
+Coords ** fragments;
+int *flagbuf;
+Coords * actpos;
+int verbose;
+
+#endif /* LOCAL_PROTO_H */
Added: grass-addons/raster/r.pi/r.pi.index/main.c
===================================================================
--- grass-addons/raster/r.pi/r.pi.index/main.c (rev 0)
+++ grass-addons/raster/r.pi/r.pi.index/main.c 2010-07-29 16:21:45 UTC (rev 42949)
@@ -0,0 +1,306 @@
+#include <string.h>
+#include <stdlib.h>
+#include <unistd.h>
+#include <grass/gis.h>
+#include <grass/glocale.h>
+#include <grass/stats.h>
+#include <math.h>
+#include "local_proto.h"
+
+/*
+ Fragmentation analysis - basic indices in r.pi.index
+
+ by Elshad Shirinov.
+*/
+
+typedef int (*f_func)(DCELL*, Coords**, int);
+
+struct menu
+{
+ f_func method; /* routine to compute new value */
+ char *name; /* method name */
+ char *text; /* menu display - full description */
+};
+
+static struct menu menu[] =
+{
+ {f_area, "area", "area of specified patches"},
+ {f_perim, "perimeter", "perimeter of specified patches"},
+ {f_shapeindex, "shape", "shape index of specified patches"},
+ {f_borderindex, "border", "border index of specified patches"},
+ {f_compactness, "compactness", "compactness of specified patches"},
+ {f_asymmetry, "asymmetry", "asymmetry of specified patches"},
+ {f_area_perim_ratio, "area-perimeter", "area-perimeter ratio of specified patches"},
+ {f_frac_dim, "fractal", "fractal dimension of specified patches"},
+ {f_nearest_dist, "ENN", "euclidean distance to nearest patch"},
+ {0,0,0}
+};
+
+int main (int argc, char *argv[])
+{
+ /* input */
+ char *newname, *oldname, *newmapset, *oldmapset;
+ char title[1024];
+
+ /* in and out file pointers */
+ int in_fd;
+ int out_fd;
+ DCELL *result, res[30];
+
+ /* map_type and categories */
+ RASTER_MAP_TYPE map_type;
+ struct Categories cats;
+
+ int method;
+ f_func compute_values;
+
+ /* neighbors count */
+ int neighb_count;
+
+ char *p;
+
+ int row, col, i, j;
+ int readrow;
+ int keyval;
+ DCELL range = MAX_DOUBLE;
+
+ int n;
+ int copycolr;
+ struct Colors colr;
+ struct GModule *module;
+ struct
+ {
+ struct Option *input, *output;
+ struct Option *keyval, *method;
+ struct Option *title;
+ } parm;
+ struct
+ {
+ struct Flag *adjacent, *quiet;
+ } flag;
+
+ DCELL *values;
+ Coords *cells;
+ int fragcount = 0;
+
+ struct Cell_head ch, window;
+
+ G_gisinit (argv[0]);
+
+ module = G_define_module();
+ module->keywords = _("raster");
+ module->description =
+ _("Provides basic patch based indices, like area, SHAPE or Nearest Neighbour distance.");
+
+ parm.input = G_define_option() ;
+ parm.input->key = "input" ;
+ parm.input->type = TYPE_STRING ;
+ parm.input->required = YES ;
+ parm.input->gisprompt = "old,cell,raster" ;
+ parm.input->description= _("Raster map containing categories") ;
+
+ parm.output = G_define_option() ;
+ parm.output->key = "output" ;
+ parm.output->type = TYPE_STRING ;
+ parm.output->required = YES ;
+ parm.output->gisprompt = "new,cell,raster" ;
+ parm.output->description= _("Output patch-based result as raster map") ;
+
+ parm.keyval = G_define_option() ;
+ parm.keyval->key = "keyval" ;
+ parm.keyval->type = TYPE_INTEGER ;
+ parm.keyval->required = YES ;
+ parm.keyval->description= _("The value of the class to be analysed") ;
+
+ parm.method = G_define_option() ;
+ parm.method->key = "method" ;
+ parm.method->type = TYPE_STRING ;
+ parm.method->required = YES ;
+ p = parm.method->options = G_malloc(1024);
+ for (n = 0; menu[n].name; n++)
+ {
+ if (n)
+ strcat (p, ",");
+ else
+ *p = 0;
+ strcat (p, menu[n].name);
+ }
+ parm.method->description= _("Operation to perform on fragments") ;
+
+ parm.title = G_define_option() ;
+ parm.title->key = "title" ;
+ parm.title->key_desc = "\"phrase\"" ;
+ parm.title->type = TYPE_STRING ;
+ parm.title->required = NO ;
+ parm.title->description= _("Title of the output raster file") ;
+
+ flag.adjacent = G_define_flag() ;
+ flag.adjacent->key = 'a' ;
+ flag.adjacent->description= _("Set for 8 cell-neighbors. 4 cell-neighbors are default.") ;
+
+ flag.quiet = G_define_flag();
+ flag.quiet->key = 'q';
+ flag.quiet->description = _("Run quietly");
+
+ if (G_parser(argc,argv))
+ exit(1);
+
+ /* get names of input files */
+ oldname = parm.input->answer;
+
+ /* test input files existance */
+ if(NULL == (oldmapset = G_find_cell2(oldname,"")))
+ {
+ fprintf (stderr, "%s: <%s> raster file not found\n",
+ G_program_name(), oldname);
+ exit(1);
+ }
+
+ /* check if the new file name is correct */
+ newname = parm.output->answer;
+ if (G_legal_filename(newname) < 0)
+ {
+ fprintf (stderr, "%s: <%s> illegal file name\n",
+ G_program_name(), newname);
+ exit(1);
+ }
+ newmapset = G_mapset();
+
+ /* get size */
+ nrows = G_window_rows();
+ ncols = G_window_cols();
+
+ /* open cell files */
+ if ((in_fd = G_open_cell_old (oldname, oldmapset)) < 0)
+ {
+ char msg[200];
+ sprintf(msg,"can't open cell file <%s> in mapset %s\n",
+ oldname, oldmapset);
+ G_fatal_error (msg);
+ exit(-1);
+ }
+
+ /* get map type */
+ map_type = DCELL_TYPE;//G_raster_map_type(oldname, oldmapset);
+
+ /* copy color table */
+ copycolr = (G_read_colors (oldname, oldmapset, &colr) > 0);
+
+ /* get key value */
+ sscanf(parm.keyval->answer, "%d", &keyval);
+
+ /* get the method */
+ for (method = 0; (p = menu[method].name); method++)
+ if ((strcmp(p, parm.method->answer) == 0))
+ break;
+ if (!p)
+ {
+ G_warning (_("<%s=%s> unknown %s"),
+ parm.method->key, parm.method->answer, parm.method->key);
+ G_usage();
+ exit(EXIT_FAILURE);
+ }
+
+ /* establish the newvalue routine */
+ compute_values = menu[method].method;
+
+ /* get number of cell-neighbors */
+ neighb_count = flag.adjacent->answer ? 8 : 4;
+
+ /* allocate the cell buffers */
+ cells = (Coords *) G_malloc ( nrows * ncols * sizeof(Coords));
+ actpos = cells;
+ fragments = (Coords **) G_malloc (nrows * ncols * sizeof (Coords*));
+ fragments[0] = cells;
+ flagbuf = (int *) G_malloc (nrows * ncols * sizeof (int));
+ result = G_allocate_d_raster_buf();
+
+ /* get title, initialize the category and stat info */
+ if (parm.title->answer)
+ strcpy (title, parm.title->answer);
+ else
+ sprintf (title,"Fragmentation of file: %s", oldname);
+
+ /* open the new cellfile */
+ out_fd = G_open_raster_new (newname, map_type);
+ if (out_fd < 0) {
+ char msg[200];
+ sprintf(msg,"can't create new cell file <%s> in mapset %s\n",
+ newname, newmapset);
+ G_fatal_error (msg);
+ exit(1);
+ }
+
+ if (verbose = !flag.quiet->answer)
+ fprintf (stderr, "Loading Patches ... ");
+
+ /* find fragments */
+ for (row = 0; row < nrows; row++) {
+ G_get_d_raster_row (in_fd, result, row);
+ for (col = 0; col < ncols; col++) {
+ if(result[col] == keyval)
+ flagbuf[row * ncols + col] = 1;
+ }
+
+ if (verbose)
+ G_percent (row, nrows, 2);
+ }
+
+ for(row = 0; row < nrows; row++) {
+ for(col = 0; col < ncols; col++) {
+ if(flagbuf[row * ncols + col] == 1) {
+ fragcount++;
+ writeFrag(row, col, neighb_count);
+ fragments[fragcount] = actpos;
+ }
+ }
+ }
+ if (verbose)
+ G_percent (nrows, nrows, 2);
+
+ /* perform actual function on the patches */
+ if (verbose = !flag.quiet->answer)
+ fprintf (stderr, "Performing operation ... ");
+ values = (DCELL *) G_malloc ( fragcount * sizeof(DCELL));
+ compute_values(values, fragments, fragcount);
+ if (verbose)
+ G_percent (fragcount, fragcount, 2);
+
+ /* write the output file */
+ if (verbose = !flag.quiet->answer)
+ fprintf (stderr, "Writing output ... ");
+ for(row = 0; row < nrows; row++) {
+ G_set_d_null_value(result, ncols);
+
+ for(i = 0; i < fragcount; i++) {
+ for(actpos = fragments[i]; actpos < fragments[i+1]; actpos++) {
+ if(actpos->y == row) {
+ result[actpos->x] = values[i];
+ }
+ }
+ }
+
+ G_put_d_raster_row(out_fd, result);
+
+ if (verbose)
+ G_percent (row, nrows, 2);
+ }
+
+ if (verbose)
+ G_percent (nrows, nrows, 2);
+
+ G_close_cell (out_fd);
+ G_close_cell (in_fd);
+
+ G_free(cells);
+ G_free(fragments);
+ G_free(flagbuf);
+
+ G_init_cats (0, title, &cats);
+ G_write_cats (newname, &cats);
+
+ if(copycolr)
+ G_write_colors (newname, newmapset, &colr);
+
+ exit(0);
+}
Property changes on: grass-addons/raster/r.pi/r.pi.index/main.c
___________________________________________________________________
Added: svn:executable
+ *
More information about the grass-commit
mailing list