#!/usr/bin/env python # $Id: kin2vmd.py 1159 2007-09-27 17:21:38Z oliver $ # (c) 2006 Oliver Beckstein # This software is made available under the terms of the # GNU Public License --- http://www.gnu.org/licenses/gpl.txt # This program is free software; you can redistribute it and/or modify # it under the terms of the GNU General Public License as published by # the Free Software Foundation; either version 2 of the License, or (at # your option) any later version. # # This program is distributed in the hope that it will be useful, # but WITHOUT ANY WARRANTY; without even the implied warranty of # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the # GNU General Public License for more details. # # You should have received a copy of the GNU General Public License # along with this program; if not, write to the Free Software # Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA __doc__ = """ Simple (=quick hack) converter for networks in kinemage to VMD. Vertices are built as C atoms, edges are graphics objects. Labels are drawn as VMD text objects. You might have to prune your kinemage input file; for instance, remove parts of the graphs that are not displayed. PDB format (fixed column) http://www.wwpdb.org/documentation/format23/sect9.html ATOM 159 N GLU E 19 26.999 15.659 21.956 1.00 25.64 2IFB 218 ....,....1....,....2....,....3....,....4....,....5....,....6....,....7....,....8 6----| 4--| 2| 3-| 1 3--|4___8------|8------|8------| """ import sys,os.path,string,getopt import re import networkx as NX # https://networkx.lanl.gov/wiki #------------------------------------------------------------ # Global variables #------------------------------------------------------------ verbose = 1 #------------------------------------------------------------ # Classes #------------------------------------------------------------ class Network: """Translate a kinemage [1] graphic of a network produced in Pajek [2] to VMD [3], using a pdb file [4] as an intermediate. It reconstructs the network graph based on the positions and provides it as a networkX graph [5]. This fails when nodes are located on the same position because then the edge information is ambivalent. This is not important for visual output, though. The only sections this script parses in a kinemage file are: @balls @vectors @arrows [1] Pajek http://vlado.fmf.uni-lj.si/pub/networks/pajek/ [2] kinemage: http://kinemage.biochem.duke.edu/ [3] VMD: http://www.ks.uiuc.edu/Research/vmd/ [4] pdb format: http://www.wwpdb.org/documentation/format23/sect9.html [5] https://networkx.lanl.gov/wiki """ def __init__(self,kinemagefile,vmdfile='n.tcl',pdbfile='n.pdb', write_conect=True, write_spheres=False, line_mode = 'cylinder' ): """ Constructor for Network class: >>> n = Network('h108.kin',vmdfile='h108.tcl',pdbfile='h108.pdb') write_conect True use pdb CONECT records for edges (fast and flexible but is buggy; VMD also only displays up to 12 bonds per node) False draw edges as VMD graphics objects write_spheres spheres as graphic objects (True) or treat them as atoms (False), which is more flexible line_mode cylinder, line, arrow (how to draw edges) The kinemage file is mandatory. The vmd file contains tcl commands for VMD, and the pdb file the coordinates of the vertices as pseudo atoms. Methods: Attributes: files filename dict opts options dict graph networkX graph (https://networkx.lanl.gov/wiki) """ self.file = dict(kinemage=kinemagefile, vmd=vmdfile, pdb=pdbfile) self.opts = dict(write_conect=write_conect, write_spheres=write_spheres) # translate colour names to numbers in VMD colour2id = { 'blue': 0, 'red': 1, 'gray': 2, 'orange':3, 'yellow':4, \ 'tan':5, 'silver':6, 'green':7, 'white':8, 'pink':9, \ 'cyan':10, 'purple':11, 'lime':12, 'mauve':13, 'ochre':14, \ 'iceblue':15, 'black':16, 'violet': 25 } material = { 'balls': 'HardPlastic', 'arrows': 'opaque' } radius = 1.0 # default radius colour = colour2id['red'] # default colour filebonds = {True:1, False:0} # index with write_conect # regular expressions used for parsing the kinemage file # (note use of P labels for the re.group method) number = '([-+]?[0-9]+\.[0-9]+)' # name label -- more than one letter, CAPS, '-' included # arrow{end|point} only work if mode == 'arrows' xyz = '(?P'+number+'),\s+(?P'+number+'),\s+(?P'+number+')' # /x, y, z/ r = { 'name' : re.compile('^{(?P[A-Z][A-Z-]+)}'), 'positions' : re.compile(number+"s\+"+number+"\s+"+number), 'balllist' : re.compile(\ '@balllist\s+({.*}|[0-9]+)\s+color=\s*(?P\w+)\s+'\ 'radius=\s+(?P\d+\.\d+)' ), 'vectorlist': re.compile('@vectorlist\s+{.*}\s+color=\s*(?P\w+)'), 'vectorpoint' : re.compile('^P\s+'+xyz), # anchored 'vectorend' : re.compile('^(?\w+)\s+'\ 'radius=\s*(?P\d+\.\d+)\s+' \ 'angle=\s*(?P\d+)'), 'arrowpoint' : re.compile('{[a-zA-Z"]*}P\s'), 'arrowend' : re.compile('{[a-zA-Z"]*}'), } try: kin = file(kinemagefile,'r') except IOError: msg(0,"Failed to read kinemage file %s." % kinemagefile) return 1 try: vmd = file(vmdfile,"w") except IOError: msg(0,"Failed to open vmd file %s." % vmdfile) return 1 # pseudo pdb: # residue number: changes with balllist # connections and text: use draw (tcl) file try: pdb = file(pdbfile,"w") except: msg(0,"Failed to open pdb file %s." % pdbfile) return 1 # reconstruct graph (a node is defined by the cartesian coordinates) # (not great but it's a quick hack...) self.graph = g = NX.Graph() self.nodeprops = nodeprops = dict() g.name = kinemagefile vmd.write('# Network\n') vmd.write('mol new %(file)s type pdb first 0 last -1 step 1 '\ 'autobonds 0 filebonds %(bonds)d\n' % \ {'file':pdbfile,'bonds':filebonds[write_conect]}) vmd.write('mol delrep 0 top\n') vmd.write('draw color %d\n' % colour) pdb.write('HEADER pseudo pdb file to display network nodes\n') group = inode = iarrow = ivector = nedges = 0 mode = None # state of the kinemage list of objects have_vector_point = False # state of the vector 'parser' for l in kin: m = r['balllist'].search(l) if m is not None: mode = 'balls' group = group + 1 representation = group - 1 (colour,radius) = m.group('colour', 'radius') radius = float(radius) vmd.write("draw color %s\n" % colour) vmd.write('material %s\n' % material['balls']) if write_conect: vmd.write('mol representation CPK %f 0.2 8.0 6.0\n' %\ radius) else: vmd.write('mol representation VDW %f 8\n' % (radius)) vmd.write('mol selection {chain %s}\n' % (chr(64+group))) vmd.write('mol color ColorID %s\n' % (colour2id[colour])) vmd.write('mol addrep top\n') msg(4,'Processing "balllist"') continue m = r['arrowlist'].search(l) if m is not None: mode = 'arrows' (colour,radius,angle) = m.group('colour', 'radius', 'angle') radius = float(radius) angle = float(angle) vmd.write("draw color %s\n" % colour) vmd.write('material %s\n' % material['arrows']) msg(4,'Processing "arrowlist"') continue m = r['vectorlist'].search(l) if m is not None: mode = 'vectors' colour = m.group('colour') vmd.write("draw color %s\n" % colour) vmd.write('material %s\n' % material['arrows']) msg(4,'Processing "vectorlist"') continue if mode == 'balls': m = r['name'].search(l) if m is None: continue name = m.group('name') inode = inode+1 msg(6,"[match %d] ball name = %s" % (inode,name)) l = r['name'].sub('',l) # remove label for split [x, y, z] = l.split() (x,y,z) = map(float, (x,y,z)) # add to graph node = (x,y,z) if node not in g: g.add_node(node) nodeprops[node] = {'name':name,'inode':inode} msg(3, "Added node '"+name+"' to the graph.") msg(10," with coordinates"+str(node)) else: msg(1,"WARNING: node '%s' is already in the graph "\ "as node '%s' so we skip it for the graph g." % \ (name, nodeprops[node]['name'])) if write_spheres: vmd.write("draw sphere {%f %f %f} radius %f\n" % \ (x,y,z,radius)) vmd.write("draw text {%f %f %f} %s\n" %\ (x,y,z,name)) # write a balllist as separate residues and chains # in a pdb file (see top doc string for PDB format) # * kinemage representations become chains A,B,C...Z # * each node becomes a separate residue consisting of a single atom # * the atom is taken as a carbon 'C' pdb.write("ATOM %4d %-2s %3s %1s %3d %8.3f%8.3f%8.3f\n" %\ (inode,'C',name[:3], chr(64+group), inode, x,y,z)) continue if mode == 'arrows': # NOTE: arrow mode is 'under developed'. If in doubt, look # at 'vectors' below. -- Oliver 2007-09-27 # # must be in arrow mode and must check point before end # assume that kinemage always has Point -> End m = r['arrowpoint'].search(l) if m is not None: l = r['arrowpoint'].sub('',l) # remove label for split arrow = [] # start new arrow x, y, z = l.split() p = tuple(map(float, (x,y,z))) arrow.append(p) msg(10,"[arrow %d] start %f %f %f" % (iarrow+1, p[0],p[1],p[2])) m = r['arrowend'].search(l) if m is not None: l = r['arrowend'].sub('',l) # remove label for split x, y, z = l.split() p = tuple(map(float, (x,y,z))) arrow.append(p) iarrow += 1 if not write_conect: vmd.write("draw line {%f %f %f} {%f %f %f}\n" %\ (arrow[0][0],arrow[0][1],arrow[0][2],\ arrow[1][0],arrow[1][1],arrow[1][2])) #vmd.write("draw text {%f %f %f} '%d'\n" % \ # (arrow[1][0],arrow[1][1],arrow[1][2],iarrow)) msg(10,"[arrow %d] end %f %f %f" % (iarrow, p[0],p[1],p[2])) continue # crappy redundant code: unify with arrows above!! if mode == 'vectors': # must be in vector mode and must check point before end m = r['vectorpoint'].search(l) if not have_vector_point and m is not None: have_vector_point = True vector = [] # start new vector x, y, z = m.group('x','y','z') p = tuple(map(float, (x,y,z))) vector.append(p) ivector += 1 msg(10,"[vector %d] start %f %f %f" % (ivector, p[0],p[1],p[2])) m = r['vectorend'].search(l) if have_vector_point and m is not None: have_vector_point = False x, y, z = m.group('x','y','z') p = tuple(map(float, (x,y,z))) vector.append(p) n1,n2 = vector if n1 in g and n2 in g: g.add_edge(n1,n2) nedges += 1 msg(3,"Added edge/vector %d/%d (%s,%s) to graph." % \ (nedges,ivector, nodeprops[n1]['name'],nodeprops[n2]['name'])) else: msg(1,"Failed to add edge %d to graph." % (ivector)) for n in vector: if n not in g: msg(3," node "+str(n1)+" not in graph.") else: msg(3," node "+nodeprops[n]['name']+" is in graph.") if not write_conect: vmd.write('# edge %d (%s,%s)\n' % \ (ivector,nodeprops[n1]['name'],nodeprops[n2]['name'])) if line_mode is 'cylinder': vmd.write("draw cylinder {%f %f %f} {%f %f %f} "\ "radius 0.1 resolution 6 filled no\n" %\ (vector[0][0],vector[0][1],vector[0][2],\ vector[1][0],vector[1][1],vector[1][2])) elif line_mode is 'arrow': pass else: # line_mode == 'line' vmd.write("draw line {%f %f %f} {%f %f %f}\n" %\ (vector[0][0],vector[0][1],vector[0][2],\ vector[1][0],vector[1][1],vector[1][2])) msg(10,"[vector %d] end %f %f %f" % (ivector, p[0],p[1],p[2])) continue vmd.write('set graph [atomselect top all]\n$graph set radius 1.0\n') msg(1,"Wrote '%s' with VMD tcl graphics commands." % (vmdfile)) if write_conect: if verbose > 2: msg(3,"Information about the graph g:") g.info() iedges = map( lambda (n1,n2): (nodeprops[n1]['inode'],nodeprops[n2]['inode']), g.edges()) # translate nodes to node numbers eback = map(lambda (i1,i2):(i2,i1), iedges) iedges += eback # CONECT wants forward & backward iedges.sort() # CONECT must be sorted by atom number for i1,i2 in iedges: pdb.write('CONECT%5d%5d\n' % (i1, i2)) msg(1,"Wrote pseudo pdb '%s' for nodes, including bonds." % (pdbfile)) else: msg(1,"Wrote pseudo pdb '%s' for nodes." % (pdbfile)) def usage(rc): str = """usage: kin2vmd [OPTIONS] FILE.kin Simple (=quick hack) converter for networks from Pajek [1], output as kinemage files [1], to VMD [2]. Vertices are built as C atoms, edges are bonds in a pdb file [4]. Labels are drawn as VMD text objects. You might have to prune your kinemage input file; for instance, remove parts of the graphs that are not displayed. Display in VMD with vmd -e vmdfile or load the vmd file from the menu File -> Load State. -h,--help this help -v,--verbose verbose messages --debug=N set verbosity level (1 default, verbose == 3, 10 full debugging) -f,--file=FILE.kin kinemage input file (instead of commandline) -V,--vmd=FILE.tcl name for VMD command file -P,--pdb=FILE.pdb name for pdb file [3] with vertices Limitations: The only sections this script parses in a kinemage file are: @balllist @vectorlist @arrowlist References: [1] Pajek http://vlado.fmf.uni-lj.si/pub/networks/pajek/ [2] kinemage: http://kinemage.biochem.duke.edu/ [3] VMD: http://www.ks.uiuc.edu/Research/vmd/ [4] pdb format: http://www.wwpdb.org/documentation/format23/sect9.html """ sys.stderr.write(str) sys.exit(rc) def msg(level,m): """ Print a message if the verbosity level >= the level of the message. The verbosity level 'verbose' must be declared globally """ if (verbose >= level): if (verbose > 3): """If in debug mode then show debug level with message""" m = "[dbglvl%3d] %s" % (level,m) print("%s" % m) def mainCommand(): """ Main driver for running program from the command line. (Note: python's getopt stops processing options at the FIRST non-opt argument, unlike GNU getopt) """ global verbose shortOptlist = "hvD:f:V:P:" longOptlist = ["help","verbose","debug=", "file=","vmd=","pdb="] try: opts, args = getopt.getopt(sys.argv[1:], shortOptlist, longOptlist) except getopt.GetoptError, details: msg(0,"GetoptError: %s\n" % details) usage(2) kinemagefile = pdbfile = vmdfile = None try: kinemagefile = args[0] except IndexError: pass # hopefully someone has set -f ... for o,a in opts: if o in ("-v","--verbose"): verbose = 3 elif o in ("-D","--debug"): try: verbose = int(a) except ValueError: msg(0,"Option error: Debug level must be a number\n") sys.exit(2) elif o in ("-h","--help"): usage(0) sys.exit() elif o in ("-f","--file"): kinemagefile = a elif o in ("-V","--vmd"): vmdfile = a elif o in ("-P","--pdb"): pdbfile = a if not kinemagefile: raise IOError, "Missing inputfile. See '--help' for usage." base = os.path.basename(kinemagefile) base = base.replace('.kin','') if not pdbfile: pdbfile = base + '.pdb' if not vmdfile: vmdfile = base + '.tcl' n = Network(kinemagefile,vmdfile=vmdfile,pdbfile=pdbfile) msg(3,"Now start up VMD. On the command line use\n"+ " vmd -e "+vmdfile+"\n"+ "From within VMD, load '"+vmdfile+"' with File -> Load State\n") #------------------------------------------------------------ # main #------------------------------------------------------------ if __name__ == "__main__": mainCommand()