[GRASS-SVN] r60494 - grass-addons/grass7/raster/r.findtheriver
svn_grass at osgeo.org
svn_grass at osgeo.org
Mon May 26 11:40:30 PDT 2014
Author: hcho
Date: 2014-05-26 11:40:29 -0700 (Mon, 26 May 2014)
New Revision: 60494
Modified:
grass-addons/grass7/raster/r.findtheriver/main.c
grass-addons/grass7/raster/r.findtheriver/point_list.c
grass-addons/grass7/raster/r.findtheriver/point_list.h
Log:
r.findtheriver: DEBUG => G_debug, use the standard coordinates option
Modified: grass-addons/grass7/raster/r.findtheriver/main.c
===================================================================
--- grass-addons/grass7/raster/r.findtheriver/main.c 2014-05-26 17:50:19 UTC (rev 60493)
+++ grass-addons/grass7/raster/r.findtheriver/main.c 2014-05-26 18:40:29 UTC (rev 60494)
@@ -20,14 +20,12 @@
#include <stdlib.h>
#include <string.h>
-#include <assert.h>
#include <math.h>
#include <grass/raster.h>
#include <grass/glocale.h>
#include "point_list.h"
#define THRESHOLD_DISTANCE 100
-#define DEFAULT_COORD_SEP ' '
/*
* global function declaration
@@ -59,36 +57,30 @@
int main(int argc, char *argv[])
{
+ struct GModule *module; /* GRASS module for parsing arguments */
+ struct
+ {
+ struct Option *input;
+ struct Option *window;
+ struct Option *threshold;
+ struct Option *coords;
+ struct Option *separator;
+ } opt;
struct Cell_head cellhd; /* it stores region information,
and header information of rasters */
char name[GNAME_MAX]; /* input raster name */
-
const char *mapset; /* mapset name */
-
int nrows, ncols;
-
int rowIdx, colIdx, nrows_less_one, ncols_less_one;
-
int infd; /* file descriptor */
-
RASTER_MAP_TYPE data_type; /* type of the map (CELL/DCELL/...) */
-
- struct GModule *module; /* GRASS module for parsing arguments */
-
- struct Option *optInput, *optWindow, *optThreshold, *optE, *optN, *optSep;
-
double E, N;
-
struct Cell_head window;
-
int windowSize;
-
double threshold;
+ char sep;
+ int debug;
- size_t SEP_SIZE = 2;
-
- char *sep = (char *)G_malloc(SEP_SIZE * sizeof(char));
-
/* initialize GIS environment */
G_gisinit(argv[0]); /* reads grass env, stores program name to G_program_name() */
@@ -100,82 +92,66 @@
_("Find the stream pixel nearest the input coordinate");
/* Define command options */
- optInput = G_define_option();
- optInput->key = "accumulation";
- optInput->type = TYPE_STRING;
- optInput->required = YES;
- optInput->gisprompt = "old,cell,raster";
- optInput->description =
+ opt.input = G_define_option();
+ opt.input->key = "accumulation";
+ opt.input->type = TYPE_STRING;
+ opt.input->required = YES;
+ opt.input->gisprompt = "old,cell,raster";
+ opt.input->description =
_("Name of input upstream accumulation area raster map");
- optWindow = G_define_option();
- optWindow->key = "window";
- optWindow->type = TYPE_INTEGER;
- optWindow->key_desc = "x";
- optWindow->multiple = NO;
- optWindow->required = NO;
- optWindow->description =
+ opt.window = G_define_option();
+ opt.window->key = "window";
+ opt.window->type = TYPE_INTEGER;
+ opt.window->key_desc = "x";
+ opt.window->multiple = NO;
+ opt.window->required = NO;
+ opt.window->description =
_
- ("The size of the window, in pixels, to search in for stream pixels. Must be an odd integer. If not supplied, window will be inferred based on raster resolution.");
+ ("The size of the window in pixels to search in for stream pixels. Must be an odd integer. If not supplied, window will be inferred based on raster resolution.");
- optThreshold = G_define_option();
- optThreshold->key = "threshold";
- optThreshold->type = TYPE_DOUBLE;
- optThreshold->key_desc = "x";
- optThreshold->multiple = NO;
- optThreshold->required = NO;
- optThreshold->description =
+ opt.threshold = G_define_option();
+ opt.threshold->key = "threshold";
+ opt.threshold->type = TYPE_DOUBLE;
+ opt.threshold->key_desc = "x";
+ opt.threshold->multiple = NO;
+ opt.threshold->required = NO;
+ opt.threshold->description =
_
("The threshold for distinguishing log(UAA) values of stream and non-stream pixels. If not supplied, threshold will be inferred from minimum and maximum raster values.");
- optE = G_define_option();
- optE->key = "easting";
- optE->type = TYPE_STRING;
- optE->key_desc = "x";
- optE->multiple = NO;
- optE->required = YES;
- optE->description = _("The map E grid coordinates");
+ opt.coords = G_define_standard_option(G_OPT_M_COORDS);
+ opt.coords->description = _("Coordinates of outlet point");
+ opt.coords->required = YES;
- optN = G_define_option();
- optN->key = "northing";
- optN->type = TYPE_STRING;
- optN->key_desc = "y";
- optN->multiple = NO;
- optN->required = YES;
- optN->description = _("The map N grid coordinates");
+ opt.separator = G_define_standard_option(G_OPT_F_SEP);
+ opt.separator->answer = "space";
- optSep = G_define_option();
- optSep->key = "separator";
- optSep->type = TYPE_STRING;
- optSep->key_desc = "y";
- optSep->multiple = NO;
- optSep->required = NO;
- optSep->description =
- _
- ("Coordinate separator. Defaults to ' '. Must be 1 character in length.");
-
/* options and flags parser */
if (G_parser(argc, argv))
exit(EXIT_FAILURE);
+ if (G__getenv("DEBUG"))
+ debug = atoi(G__getenv("DEBUG"));
+ else
+ debug = 0;
+
G_get_window(&window);
/* stores options and flags to variables */
- strncpy(name, optInput->answer, GNAME_MAX);
+ strncpy(name, opt.input->answer, GNAME_MAX);
/* returns NULL if the map was not found in any mapset,
* mapset name otherwise */
mapset = G_find_raster(name, "");
- if (mapset == NULL) {
+ if (mapset == NULL)
G_fatal_error(_("Raster map <%s> not found"), name);
- }
/* Get raster metadata */
Rast_get_cellhd(name, mapset, &cellhd);
- if (NULL != optWindow->answer) {
- windowSize = atoi(optWindow->answer);
- }
+ if (NULL != opt.window->answer)
+ windowSize = atoi(opt.window->answer);
else {
/* Determine window size */
double cellRes = (cellhd.ew_res + cellhd.ns_res) / 2;
@@ -184,55 +160,44 @@
if (!(windowSize & 1))
windowSize++;
}
- if ((windowSize < 2) || !(windowSize & 1)) {
- G_warning(_
- ("Invalid window size %s. Window size must be an odd integer >= 3\n"),
- optWindow->answer);
- G_usage();
- exit(EXIT_FAILURE);
- }
+ if (windowSize < 2 || !(windowSize & 1))
+ G_fatal_error(_
+ ("Invalid window size %s. Window size must be an odd integer >= 3"),
+ opt.window->answer);
G_verbose_message(_("Stream search window size %d\n"), windowSize);
- if (NULL != optThreshold->answer) {
- threshold = atof(optThreshold->answer);
- }
- else {
+ if (NULL != opt.threshold->answer)
+ threshold = atof(opt.threshold->answer);
+ else
/* Automatically determine the threshold */
threshold = -1.0;
- }
- if (threshold != -1.0 && threshold <= 0.0) {
- G_warning(_("Invalid threshold %s. Threshold must be > 0.0\n"),
- optThreshold->answer);
- G_usage();
- exit(EXIT_FAILURE);
- }
- if (!G_scan_easting(*optE->answers, &E, G_projection())) {
- G_warning(_("Illegal east coordinate <%s>\n"), optE->answer);
- G_usage();
- exit(EXIT_FAILURE);
- }
- if (!G_scan_northing(*optN->answers, &N, G_projection())) {
- G_warning(_("Illegal north coordinate <%s>\n"), optN->answer);
- G_usage();
- exit(EXIT_FAILURE);
- }
+ if (threshold != -1.0 && threshold <= 0.0)
+ G_fatal_error(_("Invalid threshold %s. Threshold must be > 0.0."),
+ opt.threshold->answer);
+
+ if (!G_scan_easting(opt.coords->answers[0], &E, G_projection()))
+ G_fatal_error(_("Illegal east coordinate '%s'"), opt.coords->answers[0]);
+
+ if (!G_scan_northing(opt.coords->answers[1], &N, G_projection()))
+ G_fatal_error(_("Illegal north coordinate '%s'"), opt.coords->answers[1]);
+
G_verbose_message(_("Input coordinates, easting %f, northing %f\n"), E, N);
- if (NULL == optSep->answer) {
- sep[0] = DEFAULT_COORD_SEP;
- }
- else {
- if (strlen(optSep->answer) != 1) {
- G_warning(_("Separator must be 1 character in length\n"));
- G_usage();
- exit(EXIT_FAILURE);
- }
- strncpy(sep, optSep->answer, SEP_SIZE);
- }
+ if (strcmp(opt.separator->answer, "newline") == 0)
+ sep = '\n';
+ else if (strcmp(opt.separator->answer, "comma") == 0)
+ sep = ',';
+ else if (strcmp(opt.separator->answer, "space") == 0)
+ sep = ' ';
+ else if (strcmp(opt.separator->answer, "tab") == 0)
+ sep = '\t';
+ else
+ sep = opt.separator->answer[0];
/* Determine the inputmap type (CELL/FCELL/DCELL) */
- data_type = Rast_map_type(name, mapset);
+ if ((data_type = Rast_map_type(name, mapset)) < 0)
+ G_fatal_error(_("Unable to determine the type of raster map <%s>"), name);
/* Open the raster - returns file descriptor (>0) */
if ((infd = Rast_open_old(name, mapset)) < 0)
@@ -254,49 +219,41 @@
nrows_less_one, ncols_less_one,
rowIdx, colIdx);
-#ifdef DEBUG
- fprintf(stderr, "Stream pixels: ");
- print_list(stderr, streamPixels, " ");
- fprintf(stderr, "\n");
- fflush(stderr);
-#endif
+ G_debug(1, "Stream pixels: ");
+ if (debug)
+ print_list(streamPixels, " ");
+
PointList_t *nearestStreamPixel =
find_nearest_point(streamPixels, colIdx, rowIdx);
if (NULL != nearestStreamPixel) {
+ if (debug) {
+ double nearestValue;
+ void *tmpRow = Rast_allocate_buf(data_type);
+ int currCol = nearestStreamPixel->col;
+ int currRow = nearestStreamPixel->row;
-#ifdef DEBUG
- double nearestValue;
+ G_debug(1, "Nearest pixel col: %d, row: %d", currCol, currRow);
- void *tmpRow = Rast_allocate_buf(data_type);
+ /* Get value of central cell */
+ Rast_get_row(infd, tmpRow, currRow, data_type);
- int currCol = nearestStreamPixel->col;
+ switch (data_type) {
+ case FCELL_TYPE:
+ nearestValue = (double)((FCELL *) tmpRow)[currCol];
+ break;
+ case DCELL_TYPE:
+ nearestValue = (double)((DCELL *) tmpRow)[currCol];
+ break;
+ default:
+ nearestValue = (double)((CELL *) tmpRow)[currCol];
+ break;
+ }
- int currRow = nearestStreamPixel->row;
-
- fprintf(stderr, "Nearest pixel col: %d, row: %d\n", currCol, currRow);
- fflush(stderr);
-
- /* Get value of central cell */
- Rast_get_row(infd, tmpRow, currRow, data_type);
-
- switch (data_type) {
- case FCELL_TYPE:
- nearestValue = (double)((FCELL *) tmpRow)[currCol];
- break;
- case DCELL_TYPE:
- nearestValue = (double)((DCELL *) tmpRow)[currCol];
- break;
- default:
- nearestValue = (double)((CELL *) tmpRow)[currCol];
- break;
+ G_debug(1, "Nearest stream pixel UAA value: %f", nearestValue);
+ G_free(tmpRow);
}
- fprintf(stderr, "Nearest stream pixel UAA value: %f\n", nearestValue);
- fflush(stderr);
- G_free(tmpRow);
-#endif
-
/* Get center of each column */
double nearestEasting =
Rast_col_to_easting(nearestStreamPixel->col + 0.5, &window);
@@ -304,17 +261,16 @@
Rast_row_to_northing(nearestStreamPixel->row + 0.5, &window);
/* Print snapped coordinates */
- fprintf(stdout, _("%f%s%f\n"), nearestEasting, sep, nearestNorthing);
+ fprintf(stdout, "%f%c%f\n", nearestEasting, sep, nearestNorthing);
fflush(stdout);
-
}
- /* Clean up */
- destroy_list(streamPixels);
-
/* closing raster maps */
Rast_close(infd);
+ /* Clean up */
+ destroy_list(streamPixels);
+
exit(EXIT_SUCCESS);
}
@@ -329,39 +285,31 @@
int ncols_less_one, int currRow,
int currCol)
{
- assert(DCELL_TYPE == dataType);
- assert(windowSize & 1);
-
PointList_t *streamPixels = NULL;
-
DCELL centralValue, tmpValue;
-
double logCentralValue, logTmpValue;
-
void *tmpRow = Rast_allocate_buf(dataType);
/* Get value of central cell */
Rast_get_row(fd, tmpRow, currRow, dataType);
switch (dataType) {
- case CELL_TYPE:
- centralValue = (double)((CELL *) tmpRow)[currCol];
- break;
case FCELL_TYPE:
centralValue = (double)((FCELL *) tmpRow)[currCol];
break;
case DCELL_TYPE:
centralValue = (double)((DCELL *) tmpRow)[currCol];
break;
+ default:
+ centralValue = (double)((CELL *) tmpRow)[currCol];
+ break;
}
if (centralValue <= 0)
centralValue = 1;
logCentralValue = log10(centralValue);
-#ifdef DEBUG
- fprintf(stderr, "logCentralValue: %f\n", logCentralValue);
- fflush(stderr);
-#endif
+ G_debug(1, "logCentralValue: %f", logCentralValue);
+
/* Determine threshold if need be */
if (-1.0 == threshold) {
struct FPRange *range =
@@ -386,18 +334,13 @@
threshold = 2.0;
}
-#ifdef DEBUG
- fprintf(stderr, "logMax: %f\n", logMax);
- fflush(stderr);
-#endif
-
+ G_debug(1, "logMax: %f", logMax);
}
G_verbose_message(_("Threshold: %f\n"), threshold);
/* Define window bounds */
int windowOffset = (windowSize - 1) / 2;
-
int minCol = currCol - windowOffset;
if (minCol < 0)
@@ -414,67 +357,57 @@
if (maxRow > nrows_less_one)
maxRow = nrows_less_one;
-#ifdef DEBUG
- fprintf(stderr, "currCol: %d, currRow: %d\n", currCol, currRow);
- fprintf(stderr, "min. col: %d, max. col: %d\n", minCol, maxCol);
- fprintf(stderr, "min. row: %d, max. row: %d\n", minRow, maxRow);
- fflush(stderr);
-#endif
+ G_debug(1, "currCol: %d, currRow: %d", currCol, currRow);
+ G_debug(1, "min. col: %d, max. col: %d", minCol, maxCol);
+ G_debug(1, "min. row: %d, max. row: %d", minRow, maxRow);
+
/* Search for stream pixels within the window */
int row, col;
for (row = minRow; row <= maxRow; row++) {
-#ifdef DEBUG
- fprintf(stderr, "row: %d\n", row);
- fflush(stderr);
-#endif
+ G_debug(1, "row: %d", row);
/* Get the current row */
Rast_get_row(fd, tmpRow, row, dataType);
for (col = minCol; col <= maxCol; col++) {
-#ifdef DEBUG
- fprintf(stderr, "\tcol: %d\n", col);
- fflush(stderr);
-#endif
+ G_debug(1, " col: %d", col);
+
switch (dataType) {
- case CELL_TYPE:
- tmpValue = (double)((CELL *) tmpRow)[col];
- break;
case FCELL_TYPE:
tmpValue = (double)((FCELL *) tmpRow)[col];
break;
case DCELL_TYPE:
tmpValue = (double)((DCELL *) tmpRow)[col];
break;
+ default:
+ tmpValue = (double)((CELL *) tmpRow)[col];
+ break;
}
logTmpValue = log10(tmpValue);
/* Test for nearby pixels that are stream pixels when compared to the central pixel */
-#ifdef DEBUG
- fprintf(stderr, "\t\tlogTmpValue: %f, logCentralValue: %f\n",
+ G_debug(1, " tmpValue: %f, centralValue: %f",
+ tmpValue, centralValue);
+ G_debug(1, " logTmpValue: %f, logCentralValue: %f",
logTmpValue, logCentralValue);
- fflush(stderr);
-#endif
- if ((logTmpValue - logCentralValue) > threshold) {
+
+ if (logTmpValue - logCentralValue > threshold) {
/* Add to list of stream pixels */
- if (NULL == streamPixels) {
+ if (NULL == streamPixels)
streamPixels = create_list(col, row);
- }
- else {
+ else
append_point(streamPixels, col, row);
- }
}
- else if ((logCentralValue - logTmpValue) > threshold) {
+ else if (logCentralValue - logTmpValue > threshold) {
/* Add to list of stream pixels */
- if (NULL == streamPixels) {
+ if (NULL == streamPixels)
streamPixels = create_list(currCol, currRow);
- }
- else {
+ else
append_point(streamPixels, currCol, currRow);
- }
}
}
}
G_free(tmpRow);
+
return streamPixels;
}
Modified: grass-addons/grass7/raster/r.findtheriver/point_list.c
===================================================================
--- grass-addons/grass7/raster/r.findtheriver/point_list.c 2014-05-26 17:50:19 UTC (rev 60493)
+++ grass-addons/grass7/raster/r.findtheriver/point_list.c 2014-05-26 18:40:29 UTC (rev 60494)
@@ -15,14 +15,9 @@
*
*****************************************************************************/
-#include <stdio.h>
#include <stdlib.h>
-#include <assert.h>
#include <math.h>
-
-#include <grass/config.h>
#include <grass/gis.h>
-
#include "point_list.h"
PointList_t *create_list(int col, int row)
@@ -37,8 +32,6 @@
PointList_t *append_point(PointList_t * const head, int col, int row)
{
- assert(NULL != head);
-
PointList_t *tmp = head;
while (NULL != tmp->next)
@@ -70,9 +63,7 @@
{
PointList_t *nearest = NULL;
-
double tmpDistance, minDistance = HUGE_VAL;
-
PointList_t *curr = head;
while (NULL != curr) {
@@ -86,12 +77,12 @@
return nearest;
}
-void print_list(FILE * fp, PointList_t * const head, const char *const sep)
+void print_list(PointList_t * const head, const char *const sep)
{
PointList_t *curr = head;
while (NULL != curr) {
- fprintf(fp, "%d%s%d\n", curr->col, sep, curr->row);
+ G_debug(1, "%d%s%d", curr->col, sep, curr->row);
curr = curr->next;
}
}
Modified: grass-addons/grass7/raster/r.findtheriver/point_list.h
===================================================================
--- grass-addons/grass7/raster/r.findtheriver/point_list.h 2014-05-26 17:50:19 UTC (rev 60493)
+++ grass-addons/grass7/raster/r.findtheriver/point_list.h 2014-05-26 18:40:29 UTC (rev 60494)
@@ -18,8 +18,6 @@
#ifndef __POINT_LIST_H__
#define __POINT_LIST_H__
-#include <stdio.h>
-
typedef struct PointNode
{
int col;
@@ -28,12 +26,9 @@
} PointList_t;
PointList_t *create_list(int col, int row);
-
PointList_t *append_point(PointList_t * const head, int col, int row);
-
void destroy_list(PointList_t * head);
-
PointList_t *find_nearest_point(PointList_t * const head, int col, int row);
+void print_list(PointList_t * const head, const char *const sep);
-void print_list(FILE * fp, PointList_t * const head, const char *const sep);
#endif
More information about the grass-commit
mailing list