2

I have generated a 3D graph network using the following code and Mayavi has been used for visualization.

import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
import networkx as nx
from mayavi import mlab


pos = [[0.1, 2, 0.3], [40, 0.5, -10],
       [0.1, -40, 0.3], [-49, 0.1, 2],
       [10.3, 0.3, 0.4], [-109, 0.3, 0.4]]
pos = pd.DataFrame(pos, columns=['x', 'y', 'z'])

ed_ls = [(x, y) for x, y in zip(range(0, 5), range(1, 6))]

G = nx.Graph()
G.add_edges_from(ed_ls)

nx.draw(G)
plt.show()


# plot 3D in mayavi
edge_size = 0.2
edge_color = (0.8, 0.8, 0.8)
bgcolor = (0, 0, 0)


mlab.figure(1, bgcolor=bgcolor)
mlab.clf()

for i, e in enumerate(G.edges()):
    # ----------------------------------------------------------------------------
    # the x,y, and z co-ordinates are here
    pts = mlab.points3d(pos['x'], pos['y'], pos['z'],
                        scale_mode='none',
                        scale_factor=1)
    # ----------------------------------------------------------------------------
    pts.mlab_source.dataset.lines = np.array(G.edges())
    tube = mlab.pipeline.tube(pts, tube_radius=edge_size)

    mlab.pipeline.surface(tube, color=edge_color)

mlab.show()

I would like to ask for suggestions on how to save this 3D graph in VTK format/ how to convert Networkx graph object to a VTK file for visualization in Paraview.

EDIT: I've tried to adapt the code available here for the input Networkx graph shared above. However, I am not able to obtain the output. I just get an empty window and the vtkpolyData is not plotted in the window.

"""
This code converts netwrokx graph to vtk polyData
ref: https://networkx.github.io/documentation/networkx-0.37/networkx.drawing.nx_vtk-pysrc.html
"""

import vtk
import pandas as pd
import networkx as nx
import matplotlib.pyplot as plt

from vtk.util.colors import banana, plum


def draw_nxvtk(G, node_pos):
    """
    Draw networkx graph in 3d with nodes at node_pos.

    See layout.py for functions that compute node positions.

    node_pos is a dictionary keyed by vertex with a three-tuple
    of x-y positions as the value.

    The node color is plum.
    The edge color is banana.

    All the nodes are the same size.

    """
    # set node positions
    np={}
    for n in G.nodes():
       try:
           np[n]=node_pos[n]
       except nx.NetworkXError:
           print("node %s doesn't have position"%n)

    nodePoints = vtk.vtkPoints()

    i=0
    for (x,y,z) in np.values():
       nodePoints.InsertPoint(i, x, y, z)
       i=i+1

    # Create a polydata to be glyphed.
    inputData = vtk.vtkPolyData()
    inputData.SetPoints(nodePoints)

    # Use sphere as glyph source.
    balls = vtk.vtkSphereSource()
    balls.SetRadius(.05)
    balls.SetPhiResolution(20)
    balls.SetThetaResolution(20)

    glyphPoints = vtk.vtkGlyph3D()
    glyphPoints.SetInputData(inputData)
    glyphPoints.SetSourceData(balls.GetOutput())

    glyphMapper = vtk.vtkPolyDataMapper()
    glyphMapper.SetInputData(glyphPoints.GetOutput())

    glyph = vtk.vtkActor()
    glyph.SetMapper(glyphMapper)
    glyph.GetProperty().SetDiffuseColor(plum)
    glyph.GetProperty().SetSpecular(.3)
    glyph.GetProperty().SetSpecularPower(30)

    # Generate the polyline for the spline.
    points = vtk.vtkPoints()
    edgeData = vtk.vtkPolyData()

    # Edges

    lines = vtk.vtkCellArray()
    i = 0
    for e in G.edges():
        # The edge e can be a 2-tuple (Graph) or a 3-tuple (Xgraph)
        u = e[0]
        v = e[1]
        if v in node_pos and u in node_pos:
            lines.InsertNextCell(2)
            for n in (u, v):
                (x, y, z) = node_pos[n]
                points.InsertPoint(i, x, y, z)
                lines.InsertCellPoint(i)
                i = i+1

    edgeData.SetPoints(points)
    edgeData.SetLines(lines)

    # Add thickness to the resulting line.
    Tubes = vtk.vtkTubeFilter()
    Tubes.SetNumberOfSides(16)
    Tubes.SetInputData(edgeData)
    Tubes.SetRadius(.01)
    #
    profileMapper = vtk.vtkPolyDataMapper()
    profileMapper.SetInputData(Tubes.GetOutput())

    #
    profile = vtk.vtkActor()
    profile.SetMapper(profileMapper)
    profile.GetProperty().SetDiffuseColor(banana)
    profile.GetProperty().SetSpecular(.3)
    profile.GetProperty().SetSpecularPower(30)

    # Now create the RenderWindow, Renderer and Interactor
    ren = vtk.vtkRenderer()
    renWin = vtk.vtkRenderWindow()
    renWin.AddRenderer(ren)

    iren = vtk.vtkRenderWindowInteractor()
    iren.SetRenderWindow(renWin)

    # Add the actors
    ren.AddActor(glyph)
    ren.AddActor(profile)

    renWin.SetSize(640, 640)

    iren.Initialize()
    renWin.Render()
    iren.Start()


if __name__ == "__main__":

    pos = [[0.1, 2, 0.3], [40, 0.5, -10],
           [0.1, -40, 0.3], [-49, 0.1, 2],
           [10.3, 0.3, 0.4], [-109, 0.3, 0.4]]
    pos = pd.DataFrame(pos, columns=['x', 'y', 'z'])
    pos_d = pos.T.to_dict(orient='list')
    
    ed_ls = [(x, y) for x, y in zip(range(0, 5), range(1, 6))]

    G = nx.Graph()
    G.add_edges_from(ed_ls)
    # nx.draw(G, with_labels=True, pos=nx.spring_layout(G))
    # plt.show()
    draw_nxvtk(G=G, node_pos=pos_d)

Suggestions on how to view the output with the polyData displayed while running the above code and how to save the vtkPolyData for importing in Paraview will be really helpful.

Natasha
  • 1,111
  • 5
  • 28
  • 66

1 Answers1

2

If you're ok with using vedo, which is built on top of vtk, this becomes easy:

import networkx as nx

pos = [[0.1, 2, 0.3],    [40, 0.5, -10],
       [0.1, -40, 0.3],  [-49, 0.1, 2],
       [10.3, 0.3, 0.4], [-109, 0.3, 0.4]]

ed_ls = [(x, y) for x, y in zip(range(0, 5), range(1, 6))]

G = nx.Graph()
G.add_edges_from(ed_ls)
nxpos = nx.spring_layout(G)
nxpts = [nxpos[pt] for pt in sorted(nxpos)]
# nx.draw(G, with_labels=True, pos=nxpos)
# plt.show()

raw_lines = [(pos[x],pos[y]) for x, y in ed_ls]
nx_lines = []
for x, y in ed_ls:
    p1 = nxpos[x].tolist() + [0] # add z-coord
    p2 = nxpos[y].tolist() + [0]
    nx_lines.append([p1,p2])

from vedo import *
raw_pts = Points(pos, r=12)
raw_edg = Lines(raw_lines).lw(2)
show(raw_pts, raw_edg, raw_pts.labels('id'),
     at=0, N=2, axes=True, sharecam=False)

nx_pts = Points(nxpts, r=12)
nx_edg = Lines(nx_lines).lw(2)
show(nx_pts, nx_edg, nx_pts.labels('id'),
     at=1, interactive=True)

write(nx_edg, 'afile.vtk') # save the lines

enter image description here

The package also supports DirectedGraphs, so a second options is:

from vedo import *
from vedo.pyplot import DirectedGraph

# Layouts: [2d, fast2d, clustering2d, circular, circular3d, cone, force, tree]
g = DirectedGraph(layout='fast2d')
g.arrow_scale =0.1
for i in range(6): g.add_child(i)
g.build()
show(g, axes=1)

write(g.unpack(0), 'afile.vtk')

EDIT: Following up the request,

How to include color mapping of the lines based on a scalar :

# ... from the first example

from vedo import *
raw_pts = Points(pos, r=12)
raw_edg = Lines(raw_lines).lw(3)

nx_pts = Points(nxpts, r=12).c('red').alpha(0.5)
nx_edg = Lines(nx_lines).lw(2)

v1 = [sin(x)  for x in range(6)]
v2 = [sqrt(x) for x in range(6)]
vc = [x for x in range(nx_edg.NCells())]
labs1 = nx_pts.labels(v1, scale=.05).c('green').shift(0.02,.04,0)
labs2 = nx_pts.labels(v2, scale=.05).c('red').shift(0.02,-.04,0)
labsc = nx_edg.labels(vc, on="cells", scale=.04, precision=1, rotZ=-45)
labsc.c('black')

nx_edg.cmap('viridis', vc, on="cells").add_scalarbar3d(c='k').shift(.2,0,0)
# nx_edg.cmap('jet', vc).add_scalarbar() # this is a 2D scalarbar

show(nx_pts, nx_edg, labs1, labs2, labsc, axes=1)

enter image description here

How to interpolate the line color to the node values:

(NB: here clean() removes duplicate points so please double check with possible mismatches with the initial array)

from vedo import *

nx_pts = Points(nxpts, r=12).c('grey').alpha(0.5)
nx_edg = Lines(nx_lines).lw(5)

v1 = [sin(x) for x in range(6)]
labs1 = nx_pts.labels(v1, scale=.05).c('green').shift(0.02,.04,0)

nx_edg.clean().cmap('viridis', v1).add_scalarbar()

show(nx_pts, nx_edg, labs1, axes=1)

enter image description here

mmusy
  • 1,292
  • 9
  • 10
  • Thanks a lot. This helps! May I ask an additional question? I have values corresponding to each node in the above-mentioned graph (time-series data) and I would like to use an interpolator to interpolate the value between two nodes and animate it over time. I would like to know if this can be done and visualized in vedo (I saw a couple of examples on github repo, but I couldn't locate an example of this kind). I can load the time series data corresponding to each node in an np array. I'm not sure how to take it from here. – Natasha Jul 19 '20 at 03:22
  • As a followup to my original question, I would like to know how to display node attributes. For instance, ` raw_pts.labels('id'),` displays node ids. I'd like to know if additional attributes of the nodes can be displayed? – Natasha Jul 19 '20 at 03:32
  • Hi @Natasha - yes you can add anything you like: change 'id' with a python list of same length as points, you can also add edge values, scale them, change color etc.. – mmusy Jul 19 '20 at 09:28
  • Hi @mmusy Thank you. So if I have another pointt attribute in a list `value` = [2.1,2.2,2.3,2.4,2.5,2.6] .When I try, `show(raw_pts, raw_edg, raw_pts.labels('id'), raw_pts.labels([2.0,2.1,2.2,2.3,2.4,2.5]), at=0, N=2, axes=True, sharecam=False)` graphics window closes abruptly for some reason. Could you please suggest what's the right way to do this? – Natasha Jul 19 '20 at 10:24
  • 1
    oh.. try adding `interactive=True` in show(). Note that you can modify the position of the labels by adding an offset: `mylabels = raw_pts.labels([2.0,2.1,2.2,2.3,2.4,2.5]).addPos(0,-1,0)` – mmusy Jul 19 '20 at 10:28
  • Thank you. Now I could display both the labels. Is there an option to display the label only while hovering over the point and how do we change the font size of the labels?. Also, is there a way to interpolate the `value` between two points (this would color the edge)? Please let me know if I have to post this as a separate question. – Natasha Jul 19 '20 at 11:07
  • 1
    Yes you can do something like that - please check the updated answer. – mmusy Jul 19 '20 at 11:50
  • Thank you very much for the update. Could you please explain how the edges are colored? If my understanding is right, the scalarbar assigns 1 color to each edge based on the edge id in `vc`. But I'd like to know how to color the edge based on the values assigned for two nodes that form the edge. For instance, `v1` has the values of nodes : node 0 = 0 and node1 = 0.841; the color of `edge 0` will show a gradient scale that varies from 0 to 0.841 (the scalarbar will range from 0 to max(v1)). I believe this has to be done through interpolation (similar to finite element mesh interpolation). – Natasha Jul 19 '20 at 12:13
  • Edges are colored based on the scalar value (can be any float) assigned to it in `vc`, in your case the simplest option is to assign the average value of the two nodes of the edge. You can also color points with: `nx_pts.pointColors(v1, cmap='viridis')` You may also try to color the points of the segments (to have a gradient) but then you have to pay attention to the fact that segments of lines have duplicate points! PS: in general you'll get a faster support if you post to vedo/issues ;) – mmusy Jul 19 '20 at 12:29
  • Thanks a lot for the suggestion. But I am really looking for options to color the edge based on a gradient (using values obtained by interpolation e.g [ref](https://blog.kitware.com/point-and-smoothed-particle-hydrodynamics-sph-interpolation-in-paraview/)) and don't want to assign one scalar value. – Natasha Jul 19 '20 at 12:41
  • 1
    In that case check again the updated answer and the caveat. – mmusy Jul 19 '20 at 12:55
  • Thanks a ton for the update! Could you please explain how `nx_edg.clean().pointColors(v1, cmap='viridis').addScalarBar() ` this works and also elaborate a bit on the note `duplicate points so please double check with possible mismatches with the initial array`? i'm afraid I don't understand why duplicates occur. If I may ask for 1 more clarification, say `v1` is a transient parameter observed over time from t=0 to t=10 sec; stored as [v1_0, ., v1_10]. I'd like to know if it is possible to use `show` to animate over time. I can use for loop and plot it over time, but is there an alternate way? – Natasha Jul 19 '20 at 14:02
  • Let us [continue this discussion in chat](https://chat.stackoverflow.com/rooms/218150/discussion-between-mmusy-and-natasha). – mmusy Jul 19 '20 at 14:39