Define attributes of nodes

Define attribute `coordinates` of nodes instead of separated
dictionary pos. Use function `report_graph.read_eddy_graph` to define
the attributes of nodes. So we do not need the functions
`compute_ishape` and `get_coord` any longer. Note that, for visible
eddies, `plot_traj.py` used to loop on nodes of the graph and access
corresponding shapes, while `report_graph.read_eddy_graph` loops on
the shapes and accesses corresponding nodes. This is probably quicker.

Analysis/plot_traj.py

@@ -5,40 +5,9 @@ and SHP-triplet directory for interpolated eddies.
 
 """
 
-import sys
-
-def compute_ishape(n, k1, ishape_last):
-    """n: node, tuple (date_index, eddy_index). k1: first date in triplet
-    of shapefiles.
-
-    """
-    
-    if n[0] == k1:
-        ishape_first = 0
-    elif n[0] > k1:
-        ishape_first = ishape_last[n[0] - k1 - 1] + 1
-    else:
-        print("n =", n)
-        print("k1 =", k1)
-        sys.exit("compute_ishape: bad dates")
-
-    ishape = ishape_first + n[1] - 1
-    
-    if ishape > ishape_last[n[0] - k1]:
-        # bad input n
-        return None
-    else:
-        return ishape
-
-def get_coord(n, k1, ishape_last, s_read):
-    ishape = compute_ishape(n, k1, ishape_last)
-    return s_read.shape(ishape).points[0]
-
 if __name__ == "__main__":
+    import sys
     import itertools
-    import shapefile
-    import numpy as np
-    from os import path
     import report_graph
     import matplotlib.pyplot as plt
     import networkx as nx
@@ -46,26 +15,8 @@ if __name__ == "__main__":
     if len(sys.argv) != 2:
         sys.exit("Required argument: directory containing triplet of "
                  "shapefiles, result of extracting the eddies at all dates")
-        
-    shp_tr_dir = sys.argv[1]
-    G = report_graph.read_eddy_graph("edgelist.csv")
-    ishape_last = np.loadtxt(path.join(shp_tr_dir, "ishape_last.txt"),
-                             dtype = int)
-    pos = {}
-    
-    # Interpolated eddies:
-    with shapefile.Reader(path.join("SHP_triplet", "extremum")) as s_read:
-        for shape_rec in s_read:
-            n = (shape_rec.record.date_index, shape_rec.record.eddy_index)
-            pos[n] = shape_rec.shape.points[0]
-            G.nodes[n]["interpolat"] = True
-
-    with shapefile.Reader(path.join(shp_tr_dir, "extremum")) as s_read:
-        k1 = s_read.record(0)["date_index"]
-        for n in G:
-            if "interpolat" not in G.nodes[n]:
-                pos[n] = get_coord(n, k1, ishape_last, s_read)
 
+    G = report_graph.read_eddy_graph("edgelist.csv", sys.argv[1])
     color_iter = itertools.cycle(('#1f77b4', '#aec7e8', '#ff7f0e', '#ffbb78',
                                   '#2ca02c', '#98df8a', '#d62728', '#ff9896',
                                   '#9467bd', '#c5b0d5', '#8c564b', '#c49c94',
@@ -74,14 +25,16 @@ if __name__ == "__main__":
     plt.figure()
 
     for component, color in zip(nx.weakly_connected_components(G), color_iter):
-        nx.draw_networkx(G, pos, edgelist = G.edges(component),
-                         nodelist = component, edge_color=color, node_size=10,
-                         node_color=color, with_labels=False)
+        nx.draw_networkx(G, G.nodes.data("coordinates"),
+                         edgelist = G.edges(component), nodelist = component,
+                         edge_color=color, node_size=10, node_color=color,
+                         with_labels=False)
         H = G.subgraph(component)
         t = nx.topological_sort(H)
         ancestor = next(t)
-        plt.annotate(str(ancestor), pos[ancestor], color = color,
-                     xytext = (2, 2), textcoords = 'offset points')
+        plt.annotate(str(ancestor), G.nodes[ancestor]["coordinates"],
+                     color = color, xytext = (2, 2),
+                     textcoords = 'offset points')
 
     plt.tick_params(reset = True)
     plt.show()