Delete nodes if applicable.

Author: SebKrantz

src/main/apply_geography.jl

@@ -1,6 +1,6 @@
 
 """
-    apply_geography(graph, geography; kwargs...) -> updated_graph
+    apply_geography(graph, geography; kwargs...) -> Dict
 
 Update the network building costs of a graph based on geographical features and remove edges impeded by geographical barriers. 
 Aversion to altitude changes rescales building infrastructure costs `delta_i` by (see also user manual to MATLAB toolbox):\n
@@ -90,6 +90,7 @@ function apply_geography(graph, geography; kwargs...)
     end
 
     # Remove edges where geographical barriers are (rivers)
+    remove_edge = false
     if obstacles !== nothing
 
         # Store initial delta matrics (avoid double counting)
@@ -101,6 +102,9 @@ function apply_geography(graph, geography; kwargs...)
         across_obstacle = falses(graph.J, graph.J)
         along_obstacle = falses(graph.J, graph.J)
         remove_edge = isinf(op.across_obstacle_delta_i) || isinf(op.across_obstacle_delta_tau) # remove the edge
+        if remove_edge
+            rm_nodes = falses(graph.J)
+        end
 
         for i in 1:graph.J
             neighbors = graph.nodes[i]
@@ -131,16 +135,16 @@ function apply_geography(graph, geography; kwargs...)
                                 rmi = findfirst(==(i), nodes_new[j])
                                 if rmi !== nothing
                                     deleteat!(nodes_new[j], rmi)
-                                    # if isempty(nodes_new[j])
-                                    #     deleteat!(nodes_new, j)
-                                    # end
+                                    if isempty(nodes_new[j])
+                                        rm_nodes[j] = true # deleteat!(nodes_new, j)
+                                    end
                                 end
                                 rmj = findfirst(==(j), nodes_new[i])
                                 if rmj !== nothing
                                     deleteat!(nodes_new[i], rmj)
-                                    # if isempty(nodes_new[i])
-                                    #     deleteat!(nodes_new, i)
-                                    # end
+                                    if isempty(nodes_new[i])
+                                        rm_nodes[i] = true # deleteat!(nodes_new, i)
+                                    end
                                 end
 
                                 adjacency_new[i, j] = false
@@ -165,20 +169,21 @@ function apply_geography(graph, geography; kwargs...)
         # This is for edges along obstable: allowing different delta (e.g., water transport on river)
         if isinf(op.along_obstacle_delta_i) || isinf(op.along_obstacle_delta_tau)
             # if infinite, remove edge
+            remove_edge = true
             for (io, jo) in zip(obstacles[:, 1], obstacles[:, 2])
                 rmio = findfirst(==(io), nodes_new[jo])
                 if rmio !== nothing
                     deleteat!(nodes_new[jo], rmio)
-                    # if isempty(nodes_new[jo])
-                    #     deleteat!(nodes_new, jo)
-                    # end
+                    if isempty(nodes_new[jo])
+                        rm_nodes[jo] = true # deleteat!(nodes_new, jo)
+                    end
                 end
                 rmjo = findfirst(==(jo), nodes_new[io])
                 if rmjo !== nothing
                     deleteat!(nodes_new[io], rmjo)
-                    # if isempty(nodes_new[io])
-                    #     deleteat!(nodes_new, io)
-                    # end
+                    if isempty(nodes_new[io])
+                        rm_nodes[io] = true # deleteat!(nodes_new, io)
+                    end
                 end
                 adjacency_new[io, jo] = false
                 adjacency_new[jo, io] = false
@@ -200,16 +205,47 @@ function apply_geography(graph, geography; kwargs...)
 
     # Creating new object
     graph_new = namedtuple_to_dict(graph)
-    graph_new[:delta_i] = delta_i_new
-    graph_new[:delta_tau] = delta_tau_new
 
-    if obstacles !== nothing
-        graph_new[:nodes] = nodes_new
-        graph_new[:adjacency] = adjacency_new
+    if remove_edge && any(rm_nodes)
+        keep = findall(.!rm_nodes)
+        graph_new[:delta_i] = delta_i_new[keep, keep]
+        graph_new[:delta_tau] = delta_tau_new[keep, keep]
+        graph_new[:nodes] = nodes_new[keep]
+        graph_new[:adjacency] = adjacency_new[keep, keep]
         # make sure that the degrees of freedom of the updated graph match the # of links
-        graph_new[:ndeg] = sum(tril(adjacency_new))
-        graph_new[:across_obstacle] = across_obstacle
-        graph_new[:along_obstacle] = along_obstacle
+        graph_new[:ndeg] = sum(tril(graph_new[:adjacency]))
+        graph_new[:across_obstacle] = across_obstacle[keep, keep]
+        graph_new[:along_obstacle] = along_obstacle[keep, keep]
+
+        if haskey(graph, :Lj)
+            graph_new[:Lj] = graph[:Lj][keep]
+        end
+        if haskey(graph, :Hj)
+            graph_new[:Hj] = graph[:Hj][keep]
+        end
+        if haskey(graph, :hj)
+            graph_new[:hj] = graph_new[:Hj] ./ graph_new[:Lj]
+        end
+        if haskey(graph, :omegaj)
+            graph_new[:omegaj] = graph[:omegaj][keep]
+        end
+        if haskey(graph, :Zjn)
+            graph_new[:Zjn] = graph[:Zjn][keep, :]
+        end
+        if haskey(graph, :region)
+            graph_new[:region] = graph[:region][keep]
+        end
+    else
+        graph_new[:delta_i] = delta_i_new
+        graph_new[:delta_tau] = delta_tau_new
+        if obstacles !== nothing
+            graph_new[:nodes] = nodes_new
+            graph_new[:adjacency] = adjacency_new
+            # make sure that the degrees of freedom of the updated graph match the # of links
+            graph_new[:ndeg] = sum(tril(adjacency_new))
+            graph_new[:across_obstacle] = across_obstacle
+            graph_new[:along_obstacle] = along_obstacle
+        end
     end
 
     return graph_new