Scripts: placem.py

# -----------------------------------------------------------------------------
# Read a file listing models and placements (x,y,z translation and optionally
# quaternion for rotation) and open the models at the specified positions.
# Run this script from the Chimera command-line (Favorites / Command Line).
#
#       runscript /Users/smith/placem.py /mydata/placements.txt
#
# Example file placements.txt:
#
# map1.mrc   10.25  -15.44   24.5
# map2.mrc    9.24  12.24   23.1
# proteinX.pdb  3.24  43.243  11.11   .7071  0   0  .7071
#
# The final 4 fields specify a rotation as a quaternion:
#
#     cos(a/2)  vx*sin(a/2)  vy*sin(a/2)  vz*sin(a/2)
#
# where a is rotation angle and v = (vx,vy,vz) is normalized rotation axis.
# If not specified, not rotation is applied.
#
def place_models(path):

  from chimera import openModels, Xform, replyobj, viewer

  # Make placecements relative to an already opened model.
  mlist = openModels.list()
  rel = mlist[0].openState.xform if mlist else Xform()

  pxlist = read_placement_file(path)
  for path, xform in pxlist:
    xform.premultiply(rel)
    for m in openModels.open(path):
      m.openState.xform = xform

  replyobj.status('placed %d models' % len(pxlist))
  
  viewer.viewAll()
  
# -----------------------------------------------------------------------------
#
def read_placement_file(path):

  try:
    f = open(path, 'r')
  except IOError:
    from Midas import MidasError as err
    raise err, 'Could not open "%s"' % path
  lines = f.readlines()
  f.close()

  pxlist = []
  from os.path import exists, dirname, join
  dir = dirname(path)
  from chimera import Xform
  warnings = []
  for line in lines:
    fields = line.split()
    if line[0] == '#':
      continue
    mpath = join(dir, fields[0])
    if not exists(mpath):
      warnings.append('File does not exist "%s"' % mpath)
      continue
    if len(fields) != 4 and len(fields) != 8:
      warnings.append('Line must have 4 or 8 fields "%s"' % line.strip())
      continue
    try:
      x,y,z = [float(s) for s in fields[1:4]]
    except ValueError:
      warnings.append('Could not parse shift values "%s"' % line.strip())
      continue
    if len(fields) == 8:
      try:
        q = [float(s) for s in fields[4:]]
      except ValueError:
        warnings.append('Could not parse quaternion values "%s"' % line.strip())
        continue
    else:
      q = (1,0,0,0) # Identity rotation

    a,b,c,d = q
    xf = Xform.xform(a*a+b*b-c*c-d*d, 2*b*c-2*a*d, 2*b*d+2*a*c, x,
                     2*b*c+2*a*d, a*a-b*b+c*c-d*d, 2*c*d-2*a*b, y,
                     2*b*d-2*a*c, 2*c*d+2*a*b, a*a-b*b-c*c+d*d, z,
                     orthogonalize = True)
    pxlist.append((mpath, xf))

  if warnings:
    from chimera import replyobj
    replyobj.warning('\n'.join(warnings))

  return pxlist

# -----------------------------------------------------------------------------
#
from Midas import MidasError as err
if not 'arguments' in globals() or len(arguments) != 1:
  raise err, 'Syntax: runscript placem.py <placement-file>'

place_models(arguments[0])