[GRASS-SVN] r68024 - in grass/branches/releasebranch_7_0: include/defs lib/vector/neta
svn_grass at osgeo.org
svn_grass at osgeo.org
Tue Mar 8 12:56:07 PST 2016
Author: mmetz
Date: 2016-03-08 12:56:07 -0800 (Tue, 08 Mar 2016)
New Revision: 68024
Modified:
grass/branches/releasebranch_7_0/include/defs/neta.h
grass/branches/releasebranch_7_0/lib/vector/neta/components.c
grass/branches/releasebranch_7_0/lib/vector/neta/path.c
grass/branches/releasebranch_7_0/lib/vector/neta/spanningtree.c
grass/branches/releasebranch_7_0/lib/vector/neta/utils.c
Log:
netalib: backport bugfixes from trunk
Modified: grass/branches/releasebranch_7_0/include/defs/neta.h
===================================================================
--- grass/branches/releasebranch_7_0/include/defs/neta.h 2016-03-08 20:53:29 UTC (rev 68023)
+++ grass/branches/releasebranch_7_0/include/defs/neta.h 2016-03-08 20:56:07 UTC (rev 68024)
@@ -40,6 +40,8 @@
/*path.c */
int NetA_distance_from_points(dglGraph_s * graph, struct ilist *from, int *dst,
dglInt32_t ** prev);
+int NetA_distance_to_points(dglGraph_s * graph, struct ilist *from, int *dst,
+ dglInt32_t ** prev);
int NetA_find_path(dglGraph_s * graph, int from, int to, int *edges,
struct ilist *list);
Modified: grass/branches/releasebranch_7_0/lib/vector/neta/components.c
===================================================================
--- grass/branches/releasebranch_7_0/lib/vector/neta/components.c 2016-03-08 20:53:29 UTC (rev 68023)
+++ grass/branches/releasebranch_7_0/lib/vector/neta/components.c 2016-03-08 20:56:07 UTC (rev 68024)
@@ -11,8 +11,28 @@
(>=v2). Read the file COPYING that comes with GRASS for details.
\author Daniel Bundala (Google Summer of Code 2009)
+ \author Markus Metz
*/
+/* example:
+ *
+ * X -->-- X ---- X --<-- X ---- X
+ * N1 N2 N3 N4 N5
+ *
+ * -->--, --<-- one-way
+ * ---- both ways
+ *
+ * weakly connected:
+ * all 5 nodes, even though there is no direct path from N1 to N4, 5
+ * but N1 connects to N2, 3, and N4, 5 also connect to N2, 3
+ *
+ * strongly connected:
+ * no path from N2 to N1, no path from N3 to N4
+ * component 1: N1
+ * component 2: N2, 3
+ * Component3: N4, 5
+ */
+
#include <stdio.h>
#include <stdlib.h>
#include <grass/gis.h>
@@ -21,43 +41,53 @@
#include <grass/dgl/graph.h>
/*!
- \brief Computes weekly connected components
+ \brief Computes weakly connected components
\param graph input graph
- \param[out] component list of components
+ \param[out] component array of component ids
\return number of components
\return -1 on failure
*/
int NetA_weakly_connected_components(dglGraph_s * graph, int *component)
{
- int nnodes;
+ int nnodes, i;
dglInt32_t *stack;
- int *visited;
int stack_size, components;
dglInt32_t *cur_node;
dglNodeTraverser_s nt;
+ int have_node_costs;
+ dglInt32_t ncost;
+ if (graph->Version < 2) {
+ G_warning("Directed graph must be version 2 or 3 for NetA_weakly_connected_components()");
+ return -1;
+ }
+
components = 0;
nnodes = dglGet_NodeCount(graph);
stack = (dglInt32_t *) G_calloc(nnodes + 1, sizeof(dglInt32_t));
- visited = (int *)G_calloc(nnodes + 1, sizeof(int));
- if (!stack || !visited) {
+ if (!stack) {
G_fatal_error(_("Out of memory"));
return -1;
}
+ for (i = 1; i <= nnodes; i++)
+ component[i] = 0;
+
+ ncost = 0;
+ have_node_costs = dglGet_NodeAttrSize(graph);
+
dglNode_T_Initialize(&nt, graph);
for (cur_node = dglNode_T_First(&nt); cur_node;
cur_node = dglNode_T_Next(&nt)) {
- dglInt32_t node_id = dglNodeGet_Id(graph, cur_node);
+ dglInt32_t cur_node_id = dglNodeGet_Id(graph, cur_node);
- if (!visited[node_id]) {
- visited[node_id] = 1;
- stack[0] = node_id;
+ if (!component[cur_node_id]) {
+ stack[0] = cur_node_id;
stack_size = 1;
- component[node_id] = ++components;
+ component[cur_node_id] = ++components;
while (stack_size) {
dglInt32_t *node, *edgeset, *edge;
dglEdgesetTraverser_s et;
@@ -70,84 +100,131 @@
dglInt32_t to;
to = dglNodeGet_Id(graph, dglEdgeGet_Tail(graph, edge));
- if (!visited[to]) {
- visited[to] = 1;
+ if (!component[to]) {
component[to] = components;
- stack[stack_size++] = to;
+ /* do not go through closed nodes */
+ if (have_node_costs) {
+ memcpy(&ncost, dglNodeGet_Attr(graph, dglEdgeGet_Tail(graph, edge)),
+ sizeof(ncost));
+ }
+ if (ncost >= 0)
+ stack[stack_size++] = to;
}
}
dglEdgeset_T_Release(&et);
+
+ edgeset = dglNodeGet_InEdgeset(graph, node);
+ dglEdgeset_T_Initialize(&et, graph, edgeset);
+ for (edge = dglEdgeset_T_First(&et); edge;
+ edge = dglEdgeset_T_Next(&et)) {
+ dglInt32_t to;
+
+ to = dglNodeGet_Id(graph, dglEdgeGet_Head(graph, edge));
+ if (!component[to]) {
+ component[to] = components;
+ /* do not go through closed nodes */
+ if (have_node_costs) {
+ memcpy(&ncost, dglNodeGet_Attr(graph, dglEdgeGet_Tail(graph, edge)),
+ sizeof(ncost));
+ }
+ if (ncost >= 0)
+ stack[stack_size++] = to;
+ }
+ }
+ dglEdgeset_T_Release(&et);
}
}
}
dglNode_T_Release(&nt);
- G_free(visited);
+
+ G_free(stack);
return components;
}
/*!
- \brief Computes strongly connected components
+ \brief Computes strongly connected components with Kosaraju's
+ two-pass algorithm
\param graph input graph
- \param[out] component list of components
+ \param[out] component array of component ids
\return number of components
\return -1 on failure
*/
int NetA_strongly_connected_components(dglGraph_s * graph, int *component)
{
- int nnodes;
+ int nnodes, i;
dglInt32_t *stack, *order;
- int *visited, *processed;
+ int *processed;
int stack_size, order_size, components;
- dglInt32_t *node;
+ dglInt32_t *cur_node;
dglNodeTraverser_s nt;
+ int have_node_costs;
+ dglInt32_t ncost;
+ if (graph->Version < 2) {
+ G_warning("Directed graph must be version 2 or 3 for NetA_strongly_connected_components()");
+ return -1;
+ }
+
components = 0;
nnodes = dglGet_NodeCount(graph);
stack = (dglInt32_t *) G_calloc(nnodes + 1, sizeof(dglInt32_t));
order = (dglInt32_t *) G_calloc(nnodes + 1, sizeof(dglInt32_t));
- visited = (int *)G_calloc(nnodes + 1, sizeof(int));
processed = (int *)G_calloc(nnodes + 1, sizeof(int));
- if (!stack || !visited || !order || !processed) {
+ if (!stack || !order || !processed) {
G_fatal_error(_("Out of memory"));
return -1;
}
+ for (i = 1; i <= nnodes; i++) {
+ component[i] = 0;
+ }
+
+ ncost = 0;
+ have_node_costs = dglGet_NodeAttrSize(graph);
+
order_size = 0;
dglNode_T_Initialize(&nt, graph);
- for (node = dglNode_T_First(&nt); node; node = dglNode_T_Next(&nt)) {
- dglInt32_t node_id = dglNodeGet_Id(graph, node);
+ for (cur_node = dglNode_T_First(&nt); cur_node;
+ cur_node = dglNode_T_Next(&nt)) {
+ dglInt32_t cur_node_id = dglNodeGet_Id(graph, cur_node);
- component[node_id] = 0;
- if (!visited[node_id]) {
- visited[node_id] = 1;
- stack[0] = node_id;
+ if (!component[cur_node_id]) {
+ component[cur_node_id] = --components;
+ stack[0] = cur_node_id;
stack_size = 1;
while (stack_size) {
dglInt32_t *node, *edgeset, *edge;
dglEdgesetTraverser_s et;
- dglInt32_t cur_node_id = stack[stack_size - 1];
+ dglInt32_t node_id = stack[stack_size - 1];
- if (processed[cur_node_id]) {
+ if (processed[node_id]) {
stack_size--;
- order[order_size++] = cur_node_id;
+ order[order_size++] = node_id;
continue;
}
- processed[cur_node_id] = 1;
- node = dglGetNode(graph, cur_node_id);
+ processed[node_id] = 1;
+
+ node = dglGetNode(graph, node_id);
edgeset = dglNodeGet_OutEdgeset(graph, node);
dglEdgeset_T_Initialize(&et, graph, edgeset);
for (edge = dglEdgeset_T_First(&et); edge;
edge = dglEdgeset_T_Next(&et)) {
dglInt32_t to;
- if (dglEdgeGet_Id(graph, edge) < 0)
- continue; /*ignore backward edges */
to = dglNodeGet_Id(graph, dglEdgeGet_Tail(graph, edge));
- if (!visited[to]) {
- visited[to] = 1;
+ if (!component[to]) {
+ component[to] = components;
+ /* do not go through closed nodes */
+ if (have_node_costs) {
+ memcpy(&ncost, dglNodeGet_Attr(graph, dglEdgeGet_Tail(graph, edge)),
+ sizeof(ncost));
+ }
+ if (ncost < 0)
+ processed[to] = 1;
+
stack[stack_size++] = to;
}
}
@@ -158,41 +235,48 @@
dglNode_T_Release(&nt);
+ components = 0;
+ dglNode_T_Initialize(&nt, graph);
+
while (order_size) {
- dglInt32_t node_id = order[--order_size];
+ dglInt32_t cur_node_id = order[--order_size];
+ int cur_comp = component[cur_node_id];
- if (component[node_id])
- continue;
- components++;
- component[node_id] = components;
- stack[0] = node_id;
- stack_size = 1;
- while (stack_size) {
- dglInt32_t *node, *edgeset, *edge;
- dglEdgesetTraverser_s et;
- dglInt32_t cur_node_id = stack[--stack_size];
+ if (cur_comp < 1) {
+ component[cur_node_id] = ++components;
+ stack[0] = cur_node_id;
+ stack_size = 1;
+ while (stack_size) {
+ dglInt32_t *node, *edgeset, *edge;
+ dglEdgesetTraverser_s et;
+ dglInt32_t node_id = stack[--stack_size];
- node = dglGetNode(graph, cur_node_id);
- edgeset = dglNodeGet_OutEdgeset(graph, node);
- dglEdgeset_T_Initialize(&et, graph, edgeset);
- for (edge = dglEdgeset_T_First(&et); edge;
- edge = dglEdgeset_T_Next(&et)) {
- dglInt32_t to;
+ node = dglGetNode(graph, node_id);
+ edgeset = dglNodeGet_InEdgeset(graph, node);
+ dglEdgeset_T_Initialize(&et, graph, edgeset);
+ for (edge = dglEdgeset_T_First(&et); edge;
+ edge = dglEdgeset_T_Next(&et)) {
+ dglInt32_t to;
- if (dglEdgeGet_Id(graph, edge) > 0)
- continue; /*ignore forward edges */
- to = dglNodeGet_Id(graph, dglEdgeGet_Tail(graph, edge));
- if (!component[to]) {
- component[to] = components;
- stack[stack_size++] = to;
+ to = dglNodeGet_Id(graph, dglEdgeGet_Head(graph, edge));
+ if (component[to] == cur_comp) {
+ component[to] = components;
+ /* do not go through closed nodes */
+ if (have_node_costs) {
+ memcpy(&ncost, dglNodeGet_Attr(graph, dglEdgeGet_Head(graph, edge)),
+ sizeof(ncost));
+ }
+ if (ncost >= 0)
+ stack[stack_size++] = to;
+ }
}
+ dglEdgeset_T_Release(&et);
}
- dglEdgeset_T_Release(&et);
}
}
+ dglNode_T_Release(&nt);
G_free(stack);
- G_free(visited);
G_free(order);
G_free(processed);
return components;
Modified: grass/branches/releasebranch_7_0/lib/vector/neta/path.c
===================================================================
--- grass/branches/releasebranch_7_0/lib/vector/neta/path.c 2016-03-08 20:53:29 UTC (rev 68023)
+++ grass/branches/releasebranch_7_0/lib/vector/neta/path.c 2016-03-08 20:56:07 UTC (rev 68024)
@@ -11,6 +11,7 @@
(>=v2). Read the file COPYING that comes with GRASS for details.
\author Daniel Bundala (Google Summer of Code 2009)
+ \author Markus Metz
*/
#include <stdio.h>
@@ -22,16 +23,16 @@
#include <grass/neta.h>
/*!
- \brief Computes shortests paths to every node from nodes in "from".
+ \brief Computes shortest paths to every node from nodes in "from".
- Array "dst" contains the length of the path or -1 if the node is not
+ Array "dst" contains the cost of the path or -1 if the node is not
reachable. Prev contains edges from predecessor along the shortest
path.
\param graph input graph
\param from list of 'from' positions
- \param dst list of 'to' positions
- \param[out] prev list of edges from predecessor along the shortest path
+ \param[out] dst array of costs to reach nodes
+ \param[out] prev array of edges from predecessor along the shortest path
\return 0 on success
\return -1 on failure
@@ -41,6 +42,8 @@
{
int i, nnodes;
dglHeap_s heap;
+ int have_node_costs;
+ dglInt32_t ncost;
nnodes = dglGet_NodeCount(graph);
dglEdgesetTraverser_s et;
@@ -51,13 +54,16 @@
prev[i] = NULL;
}
+ ncost = 0;
+ have_node_costs = dglGet_NodeAttrSize(graph);
+
dglHeapInit(&heap);
for (i = 0; i < from->n_values; i++) {
int v = from->value[i];
if (dst[v] == 0)
- continue; /*ingore duplicates */
+ continue; /* ignore duplicates */
dst[v] = 0; /* make sure all from nodes are processed first */
dglHeapData_u heap_data;
@@ -68,6 +74,8 @@
dglInt32_t v, dist;
dglHeapNode_s heap_node;
dglHeapData_u heap_data;
+ dglInt32_t *edge;
+ dglInt32_t *node;
if (!dglHeapExtractMin(&heap, &heap_node))
break;
@@ -76,11 +84,20 @@
if (dst[v] < dist)
continue;
- dglInt32_t *edge;
+ node = dglGetNode(graph, v);
+ if (have_node_costs && prev[v]) {
+ memcpy(&ncost, dglNodeGet_Attr(graph, node),
+ sizeof(ncost));
+ if (ncost > 0)
+ dist += ncost;
+ /* do not go through closed nodes */
+ if (ncost < 0)
+ continue;
+ }
+
dglEdgeset_T_Initialize(&et, graph,
- dglNodeGet_OutEdgeset(graph,
- dglGetNode(graph, v)));
+ dglNodeGet_OutEdgeset(graph, node));
for (edge = dglEdgeset_T_First(&et); edge;
edge = dglEdgeset_T_Next(&et)) {
@@ -105,16 +122,123 @@
}
/*!
- \brief Find a path (minimum number of edges) from 'from' to 'to' using only edges in 'edges'.
+ \brief Computes shortest paths from every node to nodes in "to".
+ Array "dst" contains the cost of the path or -1 if the node is not
+ reachable. Nxt contains edges from successor along the shortest
+ path. This method does reverse search starting with "to" nodes and
+ going backward.
+
+ \param graph input graph
+ \param to list of 'to' positions
+ \param[out] dst array of costs to reach nodes
+ \param[out] nxt array of edges from successor along the shortest path
+
+ \return 0 on success
+ \return -1 on failure
+ */
+int NetA_distance_to_points(dglGraph_s *graph, struct ilist *to,
+ int *dst, dglInt32_t **nxt)
+{
+ int i, nnodes;
+ dglHeap_s heap;
+ dglEdgesetTraverser_s et;
+ int have_node_costs;
+ dglInt32_t ncost;
+
+ nnodes = dglGet_NodeCount(graph);
+
+ /* initialize costs and edge list */
+ for (i = 1; i <= nnodes; i++) {
+ dst[i] = -1;
+ nxt[i] = NULL;
+ }
+
+ if (graph->Version < 2) {
+ G_warning("Directed graph must be version 2 or 3 for NetA_distance_to_points()");
+ return -1;
+ }
+
+ ncost = 0;
+ have_node_costs = dglGet_NodeAttrSize(graph);
+
+ dglHeapInit(&heap);
+
+ for (i = 0; i < to->n_values; i++) {
+ int v = to->value[i];
+
+ if (dst[v] == 0)
+ continue; /* ignore duplicates */
+ dst[v] = 0; /* make sure all to nodes are processed first */
+ dglHeapData_u heap_data;
+
+ heap_data.ul = v;
+ dglHeapInsertMin(&heap, 0, ' ', heap_data);
+ }
+ while (1) {
+ dglInt32_t v, dist;
+ dglHeapNode_s heap_node;
+ dglHeapData_u heap_data;
+ dglInt32_t *edge;
+ dglInt32_t *node;
+
+ if (!dglHeapExtractMin(&heap, &heap_node))
+ break;
+ v = heap_node.value.ul;
+ dist = heap_node.key;
+ if (dst[v] < dist)
+ continue;
+
+ node = dglGetNode(graph, v);
+
+ if (have_node_costs && nxt[v]) {
+ memcpy(&ncost, dglNodeGet_Attr(graph, node),
+ sizeof(ncost));
+ if (ncost > 0)
+ dist += ncost;
+ /* do not go through closed nodes */
+ if (ncost < 0)
+ continue;
+ }
+
+ dglEdgeset_T_Initialize(&et, graph,
+ dglNodeGet_InEdgeset(graph, node));
+
+ for (edge = dglEdgeset_T_First(&et); edge;
+ edge = dglEdgeset_T_Next(&et)) {
+ dglInt32_t *from = dglEdgeGet_Head(graph, edge);
+ dglInt32_t from_id = dglNodeGet_Id(graph, from);
+ dglInt32_t d = dglEdgeGet_Cost(graph, edge);
+
+ if (dst[from_id] < 0 || dst[from_id] > dist + d) {
+ dst[from_id] = dist + d;
+ nxt[from_id] = edge;
+ heap_data.ul = from_id;
+ dglHeapInsertMin(&heap, dist + d, ' ', heap_data);
+ }
+ }
+
+ dglEdgeset_T_Release(&et);
+ }
+
+ dglHeapFree(&heap, NULL);
+
+ return 0;
+}
+
+/*!
+ \brief Find a path (minimum number of edges) from 'from' to 'to'
+ using only edges flagged as valid in 'edges'. Edge costs are not
+ considered. Closed nodes are not traversed.
+
Precisely, edge with id I is used if edges[abs(i)] == 1. List
- stores the indices of lines on the path. Method return number of
- edges or -1 if no path exist.
+ stores the indices of lines on the path. The method returns the
+ number of edges or -1 if no path exists.
\param graph input graph
\param from 'from' position
\param to 'to' position
- \param edges list of available edges
+ \param edges array of edges indicating wether an edge should be used
\param[out] list list of edges
\return number of edges
@@ -127,6 +251,8 @@
dglEdgesetTraverser_s et;
char *vis;
int begin, end, cur, nnodes;
+ int have_node_costs;
+ dglInt32_t ncost;
nnodes = dglGet_NodeCount(graph);
prev = (dglInt32_t **) G_calloc(nnodes + 1, sizeof(dglInt32_t *));
@@ -138,6 +264,9 @@
}
Vect_reset_list(list);
+ ncost = 0;
+ have_node_costs = dglGet_NodeAttrSize(graph);
+
begin = 0;
end = 1;
vis[from] = 'y';
@@ -145,22 +274,32 @@
prev[from] = NULL;
while (begin != end) {
dglInt32_t vertex = queue[begin++];
+ dglInt32_t *edge, *node;
if (vertex == to)
break;
- dglInt32_t *edge, *node = dglGetNode(graph, vertex);
+ /* do not go through closed nodes */
+ if (have_node_costs && prev[vertex]) {
+ memcpy(&ncost, dglNodeGet_Attr(graph, dglEdgeGet_Tail(graph, edge)),
+ sizeof(ncost));
+ if (ncost < 0)
+ continue;
+ }
+
+ node = dglGetNode(graph, vertex);
+
dglEdgeset_T_Initialize(&et, graph,
dglNodeGet_OutEdgeset(graph, node));
for (edge = dglEdgeset_T_First(&et); edge;
edge = dglEdgeset_T_Next(&et)) {
- dglInt32_t id = abs(dglEdgeGet_Id(graph, edge));
- dglInt32_t to =
+ dglInt32_t edge_id = abs(dglEdgeGet_Id(graph, edge));
+ dglInt32_t node_id =
dglNodeGet_Id(graph, dglEdgeGet_Tail(graph, edge));
- if (edges[id] && !vis[to]) {
- vis[to] = 'y';
- prev[to] = edge;
- queue[end++] = to;
+ if (edges[edge_id] && !vis[node_id]) {
+ vis[node_id] = 'y';
+ prev[node_id] = edge;
+ queue[end++] = node_id;
}
}
dglEdgeset_T_Release(&et);
Modified: grass/branches/releasebranch_7_0/lib/vector/neta/spanningtree.c
===================================================================
--- grass/branches/releasebranch_7_0/lib/vector/neta/spanningtree.c 2016-03-08 20:53:29 UTC (rev 68023)
+++ grass/branches/releasebranch_7_0/lib/vector/neta/spanningtree.c 2016-03-08 20:56:07 UTC (rev 68024)
@@ -136,6 +136,7 @@
Vect_list_append(tree_list, dglEdgeGet_Id(graph, perm[i].edge));
}
}
+ G_percent(index, index, 1);
G_free(perm);
uf_release(&uf);
return edges;
Modified: grass/branches/releasebranch_7_0/lib/vector/neta/utils.c
===================================================================
--- grass/branches/releasebranch_7_0/lib/vector/neta/utils.c 2016-03-08 20:53:29 UTC (rev 68023)
+++ grass/branches/releasebranch_7_0/lib/vector/neta/utils.c 2016-03-08 20:56:07 UTC (rev 68024)
@@ -1,5 +1,5 @@
/*!
- \file vector/neta/timetables.c
+ \file vector/neta/utils.c
\brief Network Analysis library - utils
@@ -91,8 +91,8 @@
the array node_costs. If there is no point with a category,
node_costs=0.
- node_costs are multiplied by 1000000 and truncated to integers (as
- is done in Vect_net_build_graph)
+ node_costs are multiplied by the graph's cost multiplier and
+ truncated to integers (as is done in Vect_net_build_graph)
\param In pointer to Map_info structure
\param layer layer number
@@ -141,8 +141,12 @@
if (!Vect_cat_get(Cats, layer, &cat))
continue;
Vect_get_line_nodes(In, i, &node, NULL);
- if (db_CatValArray_get_value_double(&vals, cat, &value) == DB_OK)
- node_costs[node] = value * 1000000.0;
+ if (db_CatValArray_get_value_double(&vals, cat, &value) == DB_OK) {
+ if (value < 0)
+ node_costs[node] = -1;
+ else
+ node_costs[node] = value * In->dgraph.cost_multip;
+ }
}
}
@@ -159,12 +163,13 @@
nodes_to_features conains the index of a feature adjecent to each
node or -1 if no such feature specified by varray
exists. Nodes_to_features might be NULL, in which case it is left
- unitialised.
+ unitialised. Nodes_to_features will be wrong if several lines connect
+ to the same node.
\param map pointer to Map_info structure
\param varray pointer to varray structure
\param[out] nodes list of node ids
- \param node_to_features ?
+ \param[out] nodes_to_features maps nodes to features
*/
void NetA_varray_to_nodes(struct Map_info *map, struct varray *varray,
struct ilist *nodes, int *nodes_to_features)
More information about the grass-commit
mailing list