r.mapcalc/r3.mapcalc: graph2() function added (contributed by glynnc, trac https://trac.osgeo.org/grass/ticket/1802)

git-svn-id: https://svn.osgeo.org/grass/grass/trunk@59362 15284696-431f-4ddb-bdfa-cd5b030d7da7

raster/r.mapcalc/func_proto.h

@@ -74,6 +74,7 @@ extern func_t f_isnull;
 extern args_t c_isnull;
 
 extern func_t f_graph;
+extern func_t f_graph2;
 extern args_t c_graph;
 
 extern func_t f_min;

raster/r.mapcalc/function.c

@@ -74,6 +74,7 @@ func_desc func_descs[] = {
     {"nmode", c_varop, f_nmode},
 
     {"graph", c_graph, f_graph},
+    {"graph2", c_graph, f_graph2},
 
     {"rand", c_binop, f_rand},
 

raster/r.mapcalc/r.mapcalc.html

@@ -301,6 +301,8 @@ exp(x)                  exponential function of x                       F
 exp(x,y)                x to the power y                                F
 float(x)                convert x to single-precision floating point    F
 graph(x,x1,y1[x2,y2..]) convert the x to a y based on points in a graph F
+graph2(x,x1[,x2,..],y1[,y2..])
+                        alternative form of graph()                     F
 if                      decision options:                               *
 if(x)                   1 if x not zero, 0 otherwise
 if(x,a)                 a if x not zero, 0 otherwise

raster/r.mapcalc/r3.mapcalc.html

@@ -203,6 +203,8 @@ exp(x)                  exponential function of x                       F
 exp(x,y)                x to the power y                                F
 float(x)                convert x to single-precision floating point    F
 graph(x,x1,y1[x2,y2..]) convert the x to a y based on points in a graph F
+graph2(x,x1[,x2,..],y1[,y2..])
+                        alternative form of graph()                     F
 if                      decision options:                               *
 if(x)                   1 if x not zero, 0 otherwise
 if(x,a)                 a if x not zero, 0 otherwise

raster/r.mapcalc/xgraph.c

@@ -90,6 +90,9 @@ int f_graph(int argc, const int *argt, void **args)
 
 	    break;
 	}
+#undef X
+#undef Y
+#undef x
 
 	continue;
 
@@ -99,3 +102,79 @@ int f_graph(int argc, const int *argt, void **args)
 
     return 0;
 }
+
+int f_graph2(int argc, const int *argt, void **args)
+{
+    DCELL **argz = (DCELL **) args;
+    DCELL *res = argz[0];
+    int n = (argc - 1) / 2;
+    int i, j;
+
+    if (argc < 3)
+	return E_ARG_LO;
+
+    if (argc % 2 == 0)
+	return E_ARG_NUM;
+
+    if (argt[0] != DCELL_TYPE)
+	return E_RES_TYPE;
+
+    for (i = 1; i <= argc; i++)
+	if (argt[i] != DCELL_TYPE)
+	    return E_ARG_TYPE;
+
+    for (i = 0; i < columns; i++) {
+#define X(j) (argz[2 + (j) + 0][i])
+#define Y(j) (argz[2 + (j) + n][i])
+#define x (argz[1][i])
+
+	if (IS_NULL_D(&x))
+	    goto null;
+
+	for (j = 0; j < n; j++)
+	    if (IS_NULL_D(&X(j)))
+		goto null;
+
+	for (j = 0; j < n - 1; j++)
+	    if (X(j + 1) <= X(j))
+		goto null;
+
+	if (x <= X(0)) {
+	    if (IS_NULL_D(&Y(0)))
+		goto null;
+	    res[i] = Y(0);
+	    continue;
+	}
+
+	if (x >= X(n - 1)) {
+	    if (IS_NULL_D(&Y(n - 1)))
+		goto null;
+	    res[i] = Y(n - 1);
+	    continue;
+	}
+
+	for (j = 0; j < n - 1; j++) {
+	    if (x > X(j + 1))
+		continue;
+
+	    if (IS_NULL_D(&Y(j)) || IS_NULL_D(&Y(j + 1)))
+		goto null;
+
+	    res[i] =
+		Y(j) + (x - X(j)) * (Y(j + 1) - Y(j)) / (X(j + 1) - X(j));
+
+	    break;
+	}
+#undef X
+#undef Y
+#undef x
+
+	continue;
+
+      null:
+	SET_NULL_D(&res[i]);
+    }
+
+    return 0;
+}
+