Przeglądaj źródła

Replace r.mfilter with r.mfilter.fp

git-svn-id: https://svn.osgeo.org/grass/grass/trunk@33131 15284696-431f-4ddb-bdfa-cd5b030d7da7
Glynn Clements 16 lat temu
rodzic
commit
a652f9429e

+ 0 - 10
raster/r.mfilter.fp/Makefile

@@ -1,10 +0,0 @@
-MODULE_TOPDIR = ../..
-
-PGM = r.mfilter.fp
-
-LIBES = $(ROWIOLIB) $(GISLIB)
-DEPENDENCIES = $(ROWIODEP) $(GISDEP)
-
-include $(MODULE_TOPDIR)/include/Make/Module.make
-
-default: cmd

+ 0 - 50
raster/r.mfilter.fp/apply.c

@@ -1,50 +0,0 @@
-#include "filter.h"
-
-/**************************************************************
- * apply_filter: apply the filter to a single neighborhood
- *
- *  filter:    filter to be applied
- *  input:     input buffers
- **************************************************************/
-DCELL apply_filter(FILTER * filter, DCELL ** input)
-{
-    int size = filter->size;
-    double **matrix = filter->matrix;
-    double divisor = filter->divisor;
-    int r, c;
-    DCELL v;
-
-    v = 0;
-
-    if (divisor == 0) {
-	int have_result = 0;
-
-	for (r = 0; r < size; r++)
-	    for (c = 0; c < size; c++) {
-		if (G_is_d_null_value(&input[r][c]))
-		    continue;
-		v += input[r][c] * matrix[r][c];
-		divisor += filter->dmatrix[r][c];
-		have_result = 1;
-	    }
-
-	if (have_result)
-	    v /= divisor;
-	else
-	    G_set_d_null_value(&v, 1);
-    }
-    else {
-	for (r = 0; r < size; r++)
-	    for (c = 0; c < size; c++) {
-		if (G_is_d_null_value(&input[r][c])) {
-		    G_set_d_null_value(&v, 1);
-		    return v;
-		}
-		v += input[r][c] * matrix[r][c];
-	    }
-
-	v /= divisor;
-    }
-
-    return v;
-}

+ 0 - 120
raster/r.mfilter.fp/execute.c

@@ -1,120 +0,0 @@
-#include <unistd.h>
-#include <grass/rowio.h>
-#include "glob.h"
-#include "filter.h"
-
-int execute_filter(ROWIO * r, int out, FILTER * filter, DCELL * cell)
-{
-    int i;
-    int count;
-    int size;
-    int row, rcount;
-    int col, ccount;
-    int startx, starty;
-    int dx, dy;
-    int mid;
-    DCELL **bufs, **box, *cp;
-
-    size = filter->size;
-    mid = size / 2;
-    bufs = (DCELL **) G_malloc(size * sizeof(DCELL *));
-    box = (DCELL **) G_malloc(size * sizeof(DCELL *));
-
-    switch (filter->start) {
-    case UR:
-	startx = ncols - size;
-	starty = 0;
-	dx = -1;
-	dy = 1;
-	break;
-    case LL:
-	startx = 0;
-	starty = nrows - size;
-	dx = 1;
-	dy = -1;
-	break;
-    case LR:
-	startx = ncols - size;
-	starty = nrows - size;
-	dx = -1;
-	dy = -1;
-	break;
-    case UL:
-    default:
-	startx = 0;
-	starty = 0;
-	dx = 1;
-	dy = 1;
-	break;
-    }
-    direction = dy;
-
-    G_debug(3, "direction %d, dx=%d, dy=%d", direction, dx, dy);
-
-    rcount = nrows - (size - 1);
-    ccount = ncols - (size - 1);
-
-    /* rewind output */
-    lseek(out, 0L, 0);
-
-    /* copy border rows to output */
-    row = starty;
-    for (i = 0; i < mid; i++) {
-	cp = (DCELL *) rowio_get(r, row);
-	write(out, cp, buflen);
-	row += dy;
-    }
-
-    /* for each row */
-    for (count = 0; count < rcount; count++) {
-	G_percent(count, rcount, 2);
-	row = starty;
-	starty += dy;
-	/* get "size" rows */
-	for (i = 0; i < size; i++) {
-	    bufs[i] = (DCELL *) rowio_get(r, row);
-	    box[i] = bufs[i] + startx;
-	    row += dy;
-	}
-	if (filter->type == SEQUENTIAL)
-	    cell = bufs[mid];
-	/* copy border */
-	cp = cell;
-	for (i = 0; i < mid; i++)
-	    *cp++ = bufs[mid][i];
-
-	/* filter row */
-	col = ccount;
-	while (col--) {
-	    if (null_only) {
-		if (G_is_d_null_value(&box[mid][mid]))
-		    *cp++ = apply_filter(filter, box);
-		else
-		    *cp++ = box[mid][mid];
-	    }
-	    else {
-		*cp++ = apply_filter(filter, box);
-	    }
-	    for (i = 0; i < size; i++)
-		box[i] += dx;
-	}
-
-	/* copy border */
-	for (i = ncols - mid; i < ncols; i++)
-	    *cp++ = bufs[mid][i];
-
-	/* write row */
-	write(out, cell, buflen);
-    }
-    G_percent(count, rcount, 2);
-
-    /* copy border rows to output */
-    row = starty + mid * dy;
-    for (i = 0; i < mid; i++) {
-	cp = (DCELL *) rowio_get(r, row);
-	write(out, cp, buflen);
-	row += dy;
-    }
-
-    return 0;
-}

+ 0 - 9
raster/r.mfilter.fp/filter

@@ -1,9 +0,0 @@
-MATRIX 5
-1 1 1 1 1
-1 1 1 1 1
-1 1 1 1 1
-1 1 1 1 1
-1 1 1 1 1
-DIVISOR  0
-TYPE     P
-

+ 0 - 6
raster/r.mfilter.fp/filter.ave

@@ -1,6 +0,0 @@
-MATRIX 3
-1 1 1
-1 1 1
-1 1 1
-DIVISOR 9
-START LR

+ 0 - 30
raster/r.mfilter.fp/filter.h

@@ -1,30 +0,0 @@
-#include <grass/gis.h>
-#include <grass/rowio.h>
-typedef struct
-{
-    int size;			/* size of filter matrix */
-    double **matrix;		/* filter coefficient matrix */
-    double **dmatrix;		/* divisor coefficient matrix */
-    double divisor;		/* filter scale factor */
-    int type;			/* sequential or parallel */
-    int start;			/* starting corner */
-} FILTER;
-
-#define PARALLEL 1
-#define SEQUENTIAL 2
-#define UL 1
-#define UR 2
-#define LL 3
-#define LR 4
-
-/* apply.c */
-DCELL apply_filter(FILTER *, DCELL **);
-
-/* getfilt.c */
-FILTER *get_filter(char *, int *, char *);
-
-/* perform.c */
-int perform_filter(char *, char *, char *, FILTER *, int, int);
-
-/* execute.c */
-int execute_filter(ROWIO *, int, FILTER *, DCELL *);

+ 0 - 162
raster/r.mfilter.fp/getfilt.c

@@ -1,162 +0,0 @@
-#include <string.h>
-#include <grass/glocale.h>
-#include "filter.h"
-#include "local_proto.h"
-
-FILTER *get_filter(char *name, int *nfilters, char *title)
-{
-    FILE *fd;
-    FILTER *filter, *f;
-    char buf[300];
-    char temp[100];
-    char label[50];
-    int row, col;
-    int n;
-    int have_divisor;
-    int have_type;
-    int have_start;
-    int count;
-    double div;
-
-    f = filter = 0;
-    count = *nfilters = 0;
-    *title = 0;
-
-    fd = fopen(name, "r");
-    if (!fd) {
-	G_fatal_error(_("Cannot open filter file '%s'"), name);
-    }
-
-    while (fgets(buf, sizeof buf, fd)) {
-	G_strip(buf);
-	if (*buf == '#' || *buf == 0)
-	    continue;
-
-	if (sscanf(buf, "%s %[^\n]", label, temp) == 2) {
-	    uppercase(label);
-	    if (strcmp(label, "TITLE") == 0) {
-		G_strip(temp);
-		strcpy(title, temp);
-		continue;
-	    }
-	}
-
-	uppercase(buf);
-	if (sscanf(buf, "MATRIX %d", &n) == 1) {
-	    if (n < 3) {
-		G_fatal_error(_("Illegal filter matrix size specified"));
-	    }
-	    if (n % 2 == 0) {
-		G_fatal_error(_("Even filter matrix size specified"));
-	    }
-
-	    count++;
-	    filter = (FILTER *) G_realloc(filter, count * sizeof(FILTER));
-	    f = &filter[count - 1];
-	    f->size = n;
-	    f->divisor = 1;
-	    f->dmatrix = NULL;
-	    f->type = PARALLEL;
-	    f->start = UL;
-	    have_divisor = 0;
-	    have_type = 0;
-	    have_start = 0;
-
-	    f->matrix = (DCELL **) G_malloc(n * sizeof(DCELL *));
-	    for (row = 0; row < n; row++)
-		f->matrix[row] = (DCELL *) G_malloc(n * sizeof(DCELL));
-
-	    for (row = 0; row < n; row++)
-		for (col = 0; col < n; col++)
-		    if (fscanf(fd, "%lf", &f->matrix[row][col]) != 1) {
-			G_fatal_error(_("Illegal filter matrix"));
-		    }
-	    continue;
-	}
-	if (sscanf(buf, "DIVISOR %lf", &div) == 1)
-	    if (sscanf(buf, "%s", label) == 1 &&
-		strcmp(label, "DIVISOR") == 0) {
-		if (!filter) {
-		    G_fatal_error(_("Filter file format error"));
-		}
-		if (have_divisor) {
-		    G_fatal_error(_("Duplicate filter divisor specified"));
-		}
-		have_divisor = 1;
-		if (sscanf(buf, "DIVISOR %lf", &div) == 1) {
-		    f->divisor = div;
-		    if (n == 0)
-			f->dmatrix = f->matrix;
-		    continue;
-		}
-		f->divisor = 0;
-		f->dmatrix = (DCELL **) G_malloc(f->size * sizeof(DCELL *));
-		for (row = 0; row < f->size; row++)
-		    f->dmatrix[row] =
-			(DCELL *) G_malloc(f->size * sizeof(DCELL));
-
-		for (row = 0; row < f->size; row++)
-		    for (col = 0; col < f->size; col++)
-			if (fscanf(fd, "%lf", &f->dmatrix[row][col]) != 1) {
-			    G_fatal_error(_("Illegal divisor matrix"));
-			}
-		continue;
-	    }
-	if (sscanf(buf, "TYPE %s", temp) == 1) {
-	    if (!filter) {
-		G_fatal_error(_("Filter file format error"));
-	    }
-	    if (have_type) {
-		G_fatal_error(_("Duplicate filter type specified"));
-	    }
-	    if (strcmp(temp, "P") == 0)
-		f->type = PARALLEL;
-	    else if (strcmp(temp, "S") == 0)
-		f->type = SEQUENTIAL;
-	    else {
-		G_fatal_error(_("Illegal filter type specified"));
-	    }
-	    have_type = 1;
-	    continue;
-	}
-	if (sscanf(buf, "START %s", temp) == 1) {
-	    if (!filter) {
-		G_fatal_error(_("Filter file format error"));
-	    }
-	    if (have_start) {
-		G_fatal_error(_("Duplicate filter start specified"));
-	    }
-	    if (strcmp(temp, "UL") == 0)
-		f->start = UL;
-	    /* disable any other starting corner until the rest of
-	     * this program handles them properly
-	     */
-	    else
-		G_warning(_("Filter start %s ignored, using UL"), temp);
-
-	    /* disable these others
-	       else if (strcmp (temp, "UR") == 0)
-	       f->start  = UR ;
-	       else if (strcmp (temp, "LL") == 0)
-	       f->start  = LL ;
-	       else if (strcmp (temp, "LR") == 0)
-	       f->start  = LR ;
-	       else
-	       {
-	       ERROR ("illegal filter type specified");
-	       }
-	     */
-	    have_start = 1;
-	    continue;
-	}
-
-	/* other lines are ignored */
-    }
-
-    if (!filter) {
-	G_fatal_error(_("Illegal filter file format"));
-    }
-
-    *nfilters = count;
-    return filter;
-}

+ 0 - 25
raster/r.mfilter.fp/getrow.c

@@ -1,25 +0,0 @@
-#include <stdlib.h>
-#include <unistd.h>
-#include <sys/types.h>
-#include <grass/gis.h>
-#include <grass/glocale.h>
-#include "glob.h"
-#include "local_proto.h"
-
-int getmaprow(int fd, void *buf, int row, int len)
-{
-    if (G_get_d_raster_row(fd, (DCELL *) buf, row) < 0)
-	G_fatal_error(_("Cannot read raster row %d"), row);
-    return 1;
-}
-
-int getrow(int fd, void *buf, int row, int len)
-{
-    if (direction > 0)
-	lseek(fd, (off_t) row * len, 0);
-    else
-	lseek(fd, (off_t) (nrows - row - 1) * len, 0);
-    if (read(fd, (DCELL *) buf, len) != len)
-	G_fatal_error(_("Error reading temporary file"));
-    return 1;
-}

+ 0 - 6
raster/r.mfilter.fp/glob.h

@@ -1,6 +0,0 @@
-
-extern int nrows, ncols;
-extern int buflen;
-extern int direction;
-extern int null_only;
-extern int preserve_edges;

+ 0 - 6
raster/r.mfilter.fp/local_proto.h

@@ -1,6 +0,0 @@
-/* getrow.c */
-int getmaprow(int, void *, int, int);
-int getrow(int, void *, int, int);
-
-/* uppercase.c */
-int uppercase(char *);

+ 0 - 141
raster/r.mfilter.fp/main.c

@@ -1,141 +0,0 @@
-
-/****************************************************************************
- *
- * MODULE:       r.mfilter
- * AUTHOR(S):    Michael Shapiro, CERL (original contributor)
- *               Roberto Flor <flor itc.it>, Markus Neteler <neteler itc.it>
- *               Glynn Clements <glynn gclements.plus.com>, Jachym Cepicky <jachym les-ejk.cz>,
- *               Jan-Oliver Wagner <jan intevation.de>
- * PURPOSE:      
- * COPYRIGHT:    (C) 1999-2006 by the GRASS Development Team
- *
- *               This program is free software under the GNU General Public
- *               License (>=v2). Read the file COPYING that comes with GRASS
- *               for details.
- *
- *****************************************************************************/
-
-#include <stdlib.h>
-#include <string.h>
-#include <stdio.h>
-#include <grass/gis.h>
-#include <grass/glocale.h>
-
-#include "filter.h"
-#include "glob.h"
-
-int nrows, ncols;
-int buflen;
-int direction;
-int null_only;
-int preserve_edges;
-
-int main(int argc, char **argv)
-{
-    FILTER *filter;
-    int nfilters;
-    int repeat;
-    char *in_name, *in_mapset;
-    char *filt_name;
-    char *out_name;
-    char title[1024];
-    char temp[300];
-    int i;
-    struct GModule *module;
-    struct Flag *flag2;
-    struct Option *opt1;
-    struct Option *opt2;
-    struct Option *opt3;
-    struct Option *opt4;
-    struct Option *opt5;
-
-    G_gisinit(argv[0]);
-
-    module = G_define_module();
-    module->keywords = _("raster, map algebra");
-    module->description = _("Raster map matrix filter.");
-
-    /* Define the different options */
-
-    opt1 = G_define_standard_option(G_OPT_R_INPUT);
-
-    opt2 = G_define_standard_option(G_OPT_R_OUTPUT);
-
-    opt3 = G_define_standard_option(G_OPT_F_INPUT);
-    opt3->key = "filter";
-    opt3->required = YES;
-    opt3->description = _("Name of filter file");
-
-    opt4 = G_define_option();
-    opt4->key = "repeat";
-    opt4->type = TYPE_INTEGER;
-    opt4->multiple = NO;
-    opt4->required = NO;
-    opt4->answer = "1";
-    opt4->description = _("Number of times to repeat the filter");
-
-    opt5 = G_define_option();
-    opt5->key = "title";
-    opt5->type = TYPE_STRING;
-    opt5->required = NO;
-    opt5->description = _("Output raster map title");
-
-    /* Define the different flags */
-
-    /* this isn't implemented at all 
-       flag3 = G_define_flag() ;
-       flag3->key         = 'p' ;
-       flag3->description = _("Preserved edge") ;
-     */
-
-    flag2 = G_define_flag();
-    flag2->key = 'z';
-    flag2->description = _("Apply filter only to null data values");
-
-    if (G_parser(argc, argv))
-	exit(EXIT_FAILURE);
-
-    /*
-       preserve_edges = flag3->answer;
-     */
-    null_only = flag2->answer;
-
-    sscanf(opt4->answer, "%d", &repeat);
-    out_name = opt2->answer;
-    filt_name = opt3->answer;
-
-    in_name = opt1->answer;
-
-    in_mapset = G_find_cell2(in_name, "");
-    if (in_mapset == NULL)
-	G_fatal_error(_("Raster map <%s> not found"), in_name);
-
-    nrows = G_window_rows();
-    ncols = G_window_cols();
-    buflen = ncols * sizeof(DCELL);
-
-    /* get the filter */
-    filter = get_filter(filt_name, &nfilters, temp);
-
-    /* make sure filter matrix won't extend outside the raster map */
-    for (i = 0; i < nfilters; i++) {
-	if (filter[i].size > ncols || filter[i].size > nrows)
-	    G_fatal_error(_("Raster map too small for the size of the filter"));
-    }
-
-
-    /* make a title for result */
-    if (opt5->answer)
-	strcpy(title, opt5->answer);
-    else {
-	if (*temp == 0)
-	    strcpy(temp, "unknown filter");
-	sprintf(title, "%s filtered using %s", in_name, temp);
-    }
-
-    perform_filter(in_name, in_mapset, out_name, filter, nfilters, repeat);
-
-    G_put_cell_title(out_name, title);
-
-    exit(EXIT_SUCCESS);
-}

+ 0 - 103
raster/r.mfilter.fp/perform.c

@@ -1,103 +0,0 @@
-#include <unistd.h>
-#include <fcntl.h>
-#include <grass/rowio.h>
-#include <grass/glocale.h>
-#include "glob.h"
-#include "filter.h"
-#include "local_proto.h"
-
-int perform_filter(char *in_name, char *in_mapset, char *out_name,
-		   FILTER * filter, int nfilters, int repeat)
-{
-    int in;
-    int out;
-    int n;
-    int pass;
-    ROWIO r;
-    char *tmp1, *tmp2;
-    int count;
-    int row;
-    DCELL *cell;
-
-
-    cell = G_allocate_d_raster_buf();
-
-    count = 0;
-    for (pass = 0; pass < repeat; pass++) {
-	G_debug(1, "Pass %d", pass + 1);
-	for (n = 0; n < nfilters; n++, count++) {
-	    G_debug(1, "Filter %d", n + 1);
-
-	    if (count == 0) {
-		in = G_open_cell_old(in_name, in_mapset);
-
-		G_debug(1, "Open raster map %s in %s = %d", in_name,
-			in_mapset, in);
-
-		if (in < 0) {
-		    G_fatal_error(_("Cannot open raster map <%s>"), in_name);
-		}
-		close(creat(tmp1 = G_tempfile(), 0666));
-		out = open(tmp1, 2);
-		if (out < 0)
-		    G_fatal_error(_("Unable to create temporary file"));
-	    }
-	    else if (count == 1) {
-
-		G_debug(1, "Closing raster map");
-
-		G_close_cell(in);
-		in = out;
-		close(creat(tmp2 = G_tempfile(), 0666));
-		out = open(tmp2, 2);
-		if (out < 0)
-		    G_fatal_error(_("Unable to create temporary file"));
-	    }
-	    else {
-		int fd;
-
-		G_debug(1, "Swap temp files");
-
-		fd = in;
-		in = out;
-		out = fd;
-	    }
-
-	    rowio_setup(&r, in, filter[n].size, buflen,
-			count ? getrow : getmaprow, NULL);
-
-	    execute_filter(&r, out, &filter[n], cell);
-
-	    rowio_release(&r);
-	}
-    }
-
-    if (count == 1)
-	G_close_cell(in);
-    else if (count > 1)
-	close(in);
-
-    /* copy final result to output raster map */
-    in = out;
-    out = G_open_fp_cell_new(out_name);
-    if (out < 0) {
-	G_fatal_error(_("Cannot create raster map <%s>"), out_name);
-    }
-
-    G_message(_("Writing raster map <%s>"), out_name);
-    for (row = 0; row < nrows; row++) {
-	getrow(in, cell, row, buflen);
-	G_put_d_raster_row(out, cell);
-    }
-
-    /* remove the temporary files before closing so that the G_close_cell()
-       has more disk to work with
-     */
-    if (count > 0)
-	unlink(tmp1);
-    if (count > 1)
-	unlink(tmp2);
-    G_close_cell(out);
-
-    return 0;
-}

+ 0 - 157
raster/r.mfilter.fp/r.mfilter.fp.html

@@ -1,157 +0,0 @@
-<h2>DESCRIPTION</h2>
-
-
-<em>r.mfilter.fp</em> filters the raster <em>input</em> to produce the
-raster <em>output</em> according to the matrix <em>filter</em> designed
-by the user (see <em>FILTERS</em> below).
-The filter is applied <em>repeat</em> times (default <em>value</em> is 1).
-The <em>output</em> raster map layer can be given a <em>TITLE</em> if desired.
-(This TITLE should be put in quotes if it contains more than one word.)
-
-With <b>-z</b> flag the filter is applied only to null values in
-the input raster map layer.  The non-null category values are not changed.
-Note that if there is more than one filter step, this rule is applied to the
-intermediate raster map layer -- only null category values which result from
-the first filter will be changed.  In most cases this will NOT be the
-desired result. Hence -z should be used only with single step filters.
-<p>
-
-The <b>filter</b> parameter defines the name of an existing, user-created
-UNIX ASCII file whose contents is a matrix defining the way in which the
-<em>input</em> file will be filtered. The format of this file is described
-below, under FILTERS.
-<p>
-
-The <b>repeat</b> parameter defines the number of times the <em>filter</em>
-is to be applied to the <em>input</em> data.
-
-<h2>FILTERS</h2>
-
-The <em>filter</em> file is a normal UNIX ASCII file designed by the user.
-It has the following format:
-<pre>
-     TITLE      TITLE
-     MATRIX     n
-                  .
-     n lines of n values
-                  .
-     DIVISOR    d
-     TYPE        S/P
-</pre>
-
-<dl>
-<dt>TITLE 
-
-<dd>A one-line TITLE for the filter.
-If a TITLE was not specified on the command line, it can be specified here.
-This TITLE would be used to construct a TITLE for the resulting raster map
-layer.  It should be a one-line description of the filter.
-
-<dt>MATRIX 
-
-<dd>The matrix (n x n) follows on the next n lines.  <em>n</em> must be
-an odd integer greater than or equal to 3.
-The matrix itself consists of n rows of n values.
-The values must be separated from each other by at least 1 blank.
-
-<dt>DIVISOR 
-
-<dd>The filter divisor is <em>d</em>.  If not specified, the default is 1.
-If the divisor is zero (0), then the divisor is dependent on the
-category values in the neighborhood
-(see HOW THE FILTER WORKS below).
-
-<dt>TYPE 
-
-<dd>The filter type.  <em>S</em> means sequential, while <em>P</em> mean parallel.
-If not specified, the default is S.
-
-
-<p>
-
-Sequential filtering happens in place.  As the filter is applied to the
-raster map layer, the category values that were changed in neighboring
-cells affect the resulting category value of the current
-cell being filtered.
-
-
-<p>
-Parallel filtering happens in such a way that the original raster
-map layer category values are used to produce the new category value.
-
-
-<p>
-More than one filter may be specified in the filter file.
-The additional filter(s) are described just like the first.
-For example, the following describes two filters:
-</p>
-
-</dl>
-
-
-<h2>EXAMPLE FILTER FILE</h2>
-
-<pre>
-      TITLE     3x3 average, non-null data only, followed by 5x5 average
-     MATRIX    3
-     1 1 1
-     1 1 1
-     1 1 1
-     DIVISOR   0
-     TYPE      P
-
-     MATRIX    5
-     1 1 1 1 1
-     1 1 1 1 1
-     1 1 1 1 1
-     1 1 1 1 1
-     1 1 1 1 1
-     DIVISOR   25
-     TYPE      P
-</pre>
-
-<h2>HOW THE FILTER WORKS</h2>
-
-The filter process produces a new category value for each cell
-in the input raster map layer by multiplying the category values of the
-cells in the n x n neighborhood around the center cell
-by the corresponding matrix value and adding them together.
-If a divisor is specified, the sum is divided by this divisor.
-(If a zero divisor was specified, then
-the divisor is computed for each cell as the sum of the MATRIX
-values where the corresponding input cell is non-null.)
-
-
-<p>
-
-If more than one filter step is specified, either because the
-repeat value was greater than one or because the filter file
-contained more than one matrix, these steps are performed
-sequentially. This means that first one filter is applied to
-the entire input raster map layer to produce an intermediate result;
-then the next filter is applied to the intermediate result to
-produce another intermediate result;  and so on, until the
-final filter is applied.  Then the output cell is written.
-
-<h2>NOTES</h2>
-
-If the resolution of the geographic region does not agree with the
-resolution of the raster map layer, unintended resampling of the original
-data may occur.  The user should be sure that the geographic region
-is set properly.
-
-<h2>SEE ALSO</h2>
-
-<em><a href="g.region.html">g.region</a></em>,
-<em><a href="r.clump.html">r.clump</a></em>,
-<em><a href="r.neighbors.html">r.neighbors</a></em>
-<em><a href="r.mfilter.html">r.mfilter</a></em>
-
-<h2>AUTHOR</h2>
-
-Glynn Clements.
-Based upon r.mfilter, by Michael Shapiro,
-U.S.Army Construction Engineering Research Laboratory
-
-<p>
-<i>Last changed: $Date$</i>

+ 0 - 8
raster/r.mfilter.fp/uppercase.c

@@ -1,8 +0,0 @@
-int uppercase(char *s)
-{
-    for (; *s; s++)
-	if (*s >= 'a' && *s <= 'z')
-	    *s += 'A' - 'a';
-
-    return 0;
-}

+ 0 - 5
raster/r.mfilter/TODO

@@ -1,5 +0,0 @@
-TODO:
-
-should accept FP filter rules
-
-6/2001

+ 30 - 34
raster/r.mfilter/apply.c

@@ -6,49 +6,45 @@
  *  filter:    filter to be applied
  *  input:     input buffers
  **************************************************************/
-CELL apply_filter(FILTER * filter, CELL ** input)
+DCELL apply_filter(FILTER * filter, DCELL ** input)
 {
-    int **matrix;
-    int size;
-    int divisor;
-    int round;
-    register int r, c;
-    register CELL v;
-
-    size = filter->size;
-    divisor = filter->divisor;
-    matrix = filter->matrix;
+    int size = filter->size;
+    double **matrix = filter->matrix;
+    double divisor = filter->divisor;
+    int r, c;
+    DCELL v;
 
     v = 0;
-    for (r = 0; r < size; r++)
-	for (c = 0; c < size; c++)
-	    v += input[r][c] * matrix[r][c];
-    /* zero divisor means compute a divisor. this is done by adding the
-       numbers is the divisor matrix where the corresponding cell value
-       is not zero.
-     */
+
     if (divisor == 0) {
-	matrix = filter->dmatrix;
+	int have_result = 0;
+
 	for (r = 0; r < size; r++)
-	    for (c = 0; c < size; c++)
-		if (input[r][c])
-		    divisor += matrix[r][c];
+	    for (c = 0; c < size; c++) {
+		if (G_is_d_null_value(&input[r][c]))
+		    continue;
+		v += input[r][c] * matrix[r][c];
+		divisor += filter->dmatrix[r][c];
+		have_result = 1;
+	    }
+
+	if (have_result)
+	    v /= divisor;
+	else
+	    G_set_d_null_value(&v, 1);
     }
+    else {
+	for (r = 0; r < size; r++)
+	    for (c = 0; c < size; c++) {
+		if (G_is_d_null_value(&input[r][c])) {
+		    G_set_d_null_value(&v, 1);
+		    return v;
+		}
+		v += input[r][c] * matrix[r][c];
+	    }
 
-    /* now round the result to nearest integer. negative numbers are rounded
-       a little differently than non-negative numbers
-     */
-    if (divisor) {
-	if (round = divisor / 2) {
-	    if ((round > 0 && v > 0) || (round < 0 && v < 0))
-		v += round;
-	    else
-		v -= round;
-	}
 	v /= divisor;
     }
-    else
-	v = 0;
 
     return v;
 }

+ 9 - 10
raster/r.mfilter/execute.c

@@ -3,7 +3,7 @@
 #include "glob.h"
 #include "filter.h"
 
-int execute_filter(ROWIO * r, int out, FILTER * filter, CELL * cell)
+int execute_filter(ROWIO * r, int out, FILTER * filter, DCELL * cell)
 {
     int i;
     int count;
@@ -13,13 +13,12 @@ int execute_filter(ROWIO * r, int out, FILTER * filter, CELL * cell)
     int startx, starty;
     int dx, dy;
     int mid;
-    CELL **bufs, **box, *cp;
-    CELL apply_filter();
+    DCELL **bufs, **box, *cp;
 
     size = filter->size;
     mid = size / 2;
-    bufs = (CELL **) G_malloc(size * sizeof(CELL *));
-    box = (CELL **) G_malloc(size * sizeof(CELL *));
+    bufs = (DCELL **) G_malloc(size * sizeof(DCELL *));
+    box = (DCELL **) G_malloc(size * sizeof(DCELL *));
 
     switch (filter->start) {
     case UR:
@@ -61,7 +60,7 @@ int execute_filter(ROWIO * r, int out, FILTER * filter, CELL * cell)
     /* copy border rows to output */
     row = starty;
     for (i = 0; i < mid; i++) {
-	cp = (CELL *) rowio_get(r, row);
+	cp = (DCELL *) rowio_get(r, row);
 	write(out, cp, buflen);
 	row += dy;
     }
@@ -73,7 +72,7 @@ int execute_filter(ROWIO * r, int out, FILTER * filter, CELL * cell)
 	starty += dy;
 	/* get "size" rows */
 	for (i = 0; i < size; i++) {
-	    bufs[i] = (CELL *) rowio_get(r, row);
+	    bufs[i] = (DCELL *) rowio_get(r, row);
 	    box[i] = bufs[i] + startx;
 	    row += dy;
 	}
@@ -87,8 +86,8 @@ int execute_filter(ROWIO * r, int out, FILTER * filter, CELL * cell)
 	/* filter row */
 	col = ccount;
 	while (col--) {
-	    if (zero_only) {
-		if (!box[mid][mid])
+	    if (null_only) {
+		if (G_is_d_null_value(&box[mid][mid]))
 		    *cp++ = apply_filter(filter, box);
 		else
 		    *cp++ = box[mid][mid];
@@ -112,7 +111,7 @@ int execute_filter(ROWIO * r, int out, FILTER * filter, CELL * cell)
     /* copy border rows to output */
     row = starty + mid * dy;
     for (i = 0; i < mid; i++) {
-	cp = (CELL *) rowio_get(r, row);
+	cp = (DCELL *) rowio_get(r, row);
 	write(out, cp, buflen);
 	row += dy;
     }

+ 7 - 8
raster/r.mfilter/filter.h

@@ -1,10 +1,11 @@
 #include <grass/gis.h>
+#include <grass/rowio.h>
 typedef struct
 {
     int size;			/* size of filter matrix */
-    int **matrix;		/* filter coefficient matrix */
-    int **dmatrix;		/* divisor coefficient matrix */
-    int divisor;		/* filter scale factor */
+    double **matrix;		/* filter coefficient matrix */
+    double **dmatrix;		/* divisor coefficient matrix */
+    double divisor;		/* filter scale factor */
     int type;			/* sequential or parallel */
     int start;			/* starting corner */
 } FILTER;
@@ -17,15 +18,13 @@ typedef struct
 #define LR 4
 
 /* apply.c */
-CELL apply_filter(FILTER *, CELL **);
+DCELL apply_filter(FILTER *, DCELL **);
 
 /* getfilt.c */
 FILTER *get_filter(char *, int *, char *);
 
 /* perform.c */
-int perform_filter(char *, char *, char *, FILTER *, int, int);
+int perform_filter(const char *, const char *, FILTER *, int, int);
 
-#ifdef GRASS_ROWIO_H
 /* execute.c */
-int execute_filter(ROWIO *, int, FILTER *, CELL *);
-#endif
+int execute_filter(ROWIO *, int, FILTER *, DCELL *);

+ 11 - 9
raster/r.mfilter/getfilt.c

@@ -16,6 +16,7 @@ FILTER *get_filter(char *name, int *nfilters, char *title)
     int have_type;
     int have_start;
     int count;
+    double div;
 
     f = filter = 0;
     count = *nfilters = 0;
@@ -61,18 +62,18 @@ FILTER *get_filter(char *name, int *nfilters, char *title)
 	    have_type = 0;
 	    have_start = 0;
 
-	    f->matrix = (int **)G_malloc(n * sizeof(int *));
+	    f->matrix = (DCELL **) G_malloc(n * sizeof(DCELL *));
 	    for (row = 0; row < n; row++)
-		f->matrix[row] = (int *)G_malloc(n * sizeof(int));
+		f->matrix[row] = (DCELL *) G_malloc(n * sizeof(DCELL));
 
 	    for (row = 0; row < n; row++)
 		for (col = 0; col < n; col++)
-		    if (fscanf(fd, "%d", &f->matrix[row][col]) != 1) {
+		    if (fscanf(fd, "%lf", &f->matrix[row][col]) != 1) {
 			G_fatal_error(_("Illegal filter matrix"));
 		    }
 	    continue;
 	}
-	if (sscanf(buf, "DIVISOR %d", &n) == 1)
+	if (sscanf(buf, "DIVISOR %lf", &div) == 1)
 	    if (sscanf(buf, "%s", label) == 1 &&
 		strcmp(label, "DIVISOR") == 0) {
 		if (!filter) {
@@ -82,20 +83,21 @@ FILTER *get_filter(char *name, int *nfilters, char *title)
 		    G_fatal_error(_("Duplicate filter divisor specified"));
 		}
 		have_divisor = 1;
-		if (sscanf(buf, "DIVISOR %d", &n) == 1) {
-		    f->divisor = n;
+		if (sscanf(buf, "DIVISOR %lf", &div) == 1) {
+		    f->divisor = div;
 		    if (n == 0)
 			f->dmatrix = f->matrix;
 		    continue;
 		}
 		f->divisor = 0;
-		f->dmatrix = (int **)G_malloc(f->size * sizeof(int *));
+		f->dmatrix = (DCELL **) G_malloc(f->size * sizeof(DCELL *));
 		for (row = 0; row < f->size; row++)
-		    f->dmatrix[row] = (int *)G_malloc(f->size * sizeof(int));
+		    f->dmatrix[row] =
+			(DCELL *) G_malloc(f->size * sizeof(DCELL));
 
 		for (row = 0; row < f->size; row++)
 		    for (col = 0; col < f->size; col++)
-			if (fscanf(fd, "%d", &f->dmatrix[row][col]) != 1) {
+			if (fscanf(fd, "%lf", &f->dmatrix[row][col]) != 1) {
 			    G_fatal_error(_("Illegal divisor matrix"));
 			}
 		continue;

+ 2 - 2
raster/r.mfilter/getrow.c

@@ -8,7 +8,7 @@
 
 int getmaprow(int fd, void *buf, int row, int len)
 {
-    if (G_get_map_row(fd, (CELL *) buf, row) < 0)
+    if (G_get_d_raster_row(fd, (DCELL *) buf, row) < 0)
 	G_fatal_error(_("Cannot read raster row %d"), row);
     return 1;
 }
@@ -19,7 +19,7 @@ int getrow(int fd, void *buf, int row, int len)
 	lseek(fd, (off_t) row * len, 0);
     else
 	lseek(fd, (off_t) (nrows - row - 1) * len, 0);
-    if (read(fd, (CELL *) buf, len) != len)
+    if (read(fd, (DCELL *) buf, len) != len)
 	G_fatal_error(_("Error reading temporary file"));
     return 1;
 }

+ 2 - 1
raster/r.mfilter/glob.h

@@ -1,5 +1,6 @@
+
 extern int nrows, ncols;
 extern int buflen;
 extern int direction;
-extern int zero_only;
+extern int null_only;
 extern int preserve_edges;

+ 8 - 11
raster/r.mfilter/main.c

@@ -20,21 +20,22 @@
 #include <stdio.h>
 #include <grass/gis.h>
 #include <grass/glocale.h>
+
 #include "filter.h"
 #include "glob.h"
 
 int nrows, ncols;
 int buflen;
 int direction;
-int zero_only;
+int null_only;
 int preserve_edges;
 
-int main(int argc, char *argv[])
+int main(int argc, char **argv)
 {
     FILTER *filter;
     int nfilters;
     int repeat;
-    char *in_name, *in_mapset;
+    char *in_name;
     char *filt_name;
     char *out_name;
     char title[1024];
@@ -89,7 +90,7 @@ int main(int argc, char *argv[])
 
     flag2 = G_define_flag();
     flag2->key = 'z';
-    flag2->description = _("Apply filter only to zero data values");
+    flag2->description = _("Apply filter only to null data values");
 
     if (G_parser(argc, argv))
 	exit(EXIT_FAILURE);
@@ -97,7 +98,7 @@ int main(int argc, char *argv[])
     /*
        preserve_edges = flag3->answer;
      */
-    zero_only = flag2->answer;
+    null_only = flag2->answer;
 
     sscanf(opt4->answer, "%d", &repeat);
     out_name = opt2->answer;
@@ -105,13 +106,9 @@ int main(int argc, char *argv[])
 
     in_name = opt1->answer;
 
-    in_mapset = G_find_cell2(in_name, "");
-    if (in_mapset == NULL)
-	G_fatal_error(_("Raster map <%s> not found"), in_name);
-
     nrows = G_window_rows();
     ncols = G_window_cols();
-    buflen = ncols * sizeof(CELL);
+    buflen = ncols * sizeof(DCELL);
 
     /* get the filter */
     filter = get_filter(filt_name, &nfilters, temp);
@@ -132,7 +129,7 @@ int main(int argc, char *argv[])
 	sprintf(title, "%s filtered using %s", in_name, temp);
     }
 
-    perform_filter(in_name, in_mapset, out_name, filter, nfilters, repeat);
+    perform_filter(in_name, out_name, filter, nfilters, repeat);
 
     G_put_cell_title(out_name, title);
 

+ 7 - 8
raster/r.mfilter/perform.c

@@ -6,7 +6,7 @@
 #include "filter.h"
 #include "local_proto.h"
 
-int perform_filter(char *in_name, char *in_mapset, char *out_name,
+int perform_filter(const char *in_name, const char *out_name,
 		   FILTER * filter, int nfilters, int repeat)
 {
     int in;
@@ -17,10 +17,10 @@ int perform_filter(char *in_name, char *in_mapset, char *out_name,
     char *tmp1, *tmp2;
     int count;
     int row;
-    CELL *cell;
+    DCELL *cell;
 
 
-    cell = G_allocate_cell_buf();
+    cell = G_allocate_d_raster_buf();
 
     count = 0;
     for (pass = 0; pass < repeat; pass++) {
@@ -29,10 +29,9 @@ int perform_filter(char *in_name, char *in_mapset, char *out_name,
 	    G_debug(1, "Filter %d", n + 1);
 
 	    if (count == 0) {
-		in = G_open_cell_old(in_name, in_mapset);
+		in = G_open_cell_old(in_name, "");
 
-		G_debug(1, "Open raster map %s in %s = %d", in_name,
-			in_mapset, in);
+		G_debug(1, "Open raster map %s = %d", in_name, in);
 
 		if (in < 0) {
 		    G_fatal_error(_("Cannot open raster map <%s>"), in_name);
@@ -79,7 +78,7 @@ int perform_filter(char *in_name, char *in_mapset, char *out_name,
 
     /* copy final result to output raster map */
     in = out;
-    out = G_open_cell_new(out_name);
+    out = G_open_fp_cell_new(out_name);
     if (out < 0) {
 	G_fatal_error(_("Cannot create raster map <%s>"), out_name);
     }
@@ -87,7 +86,7 @@ int perform_filter(char *in_name, char *in_mapset, char *out_name,
     G_message(_("Writing raster map <%s>"), out_name);
     for (row = 0; row < nrows; row++) {
 	getrow(in, cell, row, buflen);
-	G_put_raster_row(out, cell, CELL_TYPE);
+	G_put_d_raster_row(out, cell);
     }
 
     /* remove the temporary files before closing so that the G_close_cell()

+ 13 - 12
raster/r.mfilter/r.mfilter.html

@@ -1,17 +1,17 @@
 <h2>DESCRIPTION</h2>
 
 
-<em>r.mfilter</em> filters the raster <em>input</em> to produce the
+<em>r.mfilter.fp</em> filters the raster <em>input</em> to produce the
 raster <em>output</em> according to the matrix <em>filter</em> designed
 by the user (see <em>FILTERS</em> below).
 The filter is applied <em>repeat</em> times (default <em>value</em> is 1).
 The <em>output</em> raster map layer can be given a <em>TITLE</em> if desired.
 (This TITLE should be put in quotes if it contains more than one word.)
 
-With <b>-z</b> flag the filter is applied only to zero category values in
-the input raster map layer.  The non-zero category values are not changed.
+With <b>-z</b> flag the filter is applied only to null values in
+the input raster map layer.  The non-null category values are not changed.
 Note that if there is more than one filter step, this rule is applied to the
-intermediate raster map layer -- only zero category values which result from
+intermediate raster map layer -- only null category values which result from
 the first filter will be changed.  In most cases this will NOT be the
 desired result. Hence -z should be used only with single step filters.
 <p>
@@ -33,7 +33,7 @@ It has the following format:
      TITLE      TITLE
      MATRIX     n
                   .
-     n lines of n integers
+     n lines of n values
                   .
      DIVISOR    d
      TYPE        S/P
@@ -51,8 +51,8 @@ layer.  It should be a one-line description of the filter.
 
 <dd>The matrix (n x n) follows on the next n lines.  <em>n</em> must be
 an odd integer greater than or equal to 3.
-The matrix itself consists of n rows of n integers.
-The integers must be separated from each other by at least 1 blank.
+The matrix itself consists of n rows of n values.
+The values must be separated from each other by at least 1 blank.
 
 <dt>DIVISOR 
 
@@ -92,7 +92,7 @@ For example, the following describes two filters:
 <h2>EXAMPLE FILTER FILE</h2>
 
 <pre>
-      TITLE     3x3 average, non-zero data only, followed by 5x5 average
+      TITLE     3x3 average, non-null data only, followed by 5x5 average
      MATRIX    3
      1 1 1
      1 1 1
@@ -116,11 +116,10 @@ The filter process produces a new category value for each cell
 in the input raster map layer by multiplying the category values of the
 cells in the n x n neighborhood around the center cell
 by the corresponding matrix value and adding them together.
-If a divisor is specified, the sum is divided by this divisor,
-rounding to the nearest integer.
+If a divisor is specified, the sum is divided by this divisor.
 (If a zero divisor was specified, then
 the divisor is computed for each cell as the sum of the MATRIX
-values where the corresponding input cell is non-zero.)
+values where the corresponding input cell is non-null.)
 
 
 <p>
@@ -146,10 +145,12 @@ is set properly.
 <em><a href="g.region.html">g.region</a></em>,
 <em><a href="r.clump.html">r.clump</a></em>,
 <em><a href="r.neighbors.html">r.neighbors</a></em>
+<em><a href="r.mfilter.html">r.mfilter</a></em>
 
 <h2>AUTHOR</h2>
 
-Michael Shapiro,
+Glynn Clements.
+Based upon r.mfilter, by Michael Shapiro,
 U.S.Army Construction Engineering Research Laboratory
 
 <p>