User:Speleo3

From PyMOLWiki
Jump to: navigation, search

My name is Thomas Holder and I work for PyMOL at Schrödinger.

I was awarded the Warren L. DeLano Memorial PyMOL Open-Source Fellowship for 2011-2012.

Contact

  • speleo3/users.sourceforge.net
  • thomas.holder/schrodinger.com

Scripts written by me

Scripts Pastebin

Some random scripts with no dedicated PyMOLWiki page.

Launch interactive python terminal with PyMOL process: (see also Launching From a Script)

#!/usr/bin/python2.6 -i

import sys, os

# autocompletion
import readline
import rlcompleter
readline.parse_and_bind('tab: complete')

# pymol environment
moddir='/opt/pymol-svn/modules'
sys.path.insert(0, moddir)
os.putenv('PYMOL_PATH', os.path.join(moddir, 'pymol/pymol_path'))

# pymol launching
import pymol
pymol.pymol_argv = ['pymol','-qc'] + sys.argv[1:]
pymol.finish_launching()
cmd = pymol.cmd

Build FREEMOL (see also MovieSchool 6)

#!/bin/bash -e

src=/tmp
prefix=/opt/pymol-git
export FREEMOL=$prefix/freemol

freemoltrunk=$src/freemol-trunk
if [[ ! -e $freemoltrunk ]]; then
    svn co svn://bioinformatics.org/svnroot/freemol/trunk $freemoltrunk
fi

cd $freemoltrunk

sed -i 's/vdwtype\[11\]/vdwtype[14]/' src/mengine/src/field.h

for name in mpeg_encode mengine apbs pdb2pqr; do
    (cd src/$name && ./configure && make && make install)
done

cp -na freemol/libpy/freemol $prefix/modules/

ln -sfT $FREEMOL $prefix/modules/pymol/pymol_path/freemol

Download all PyMOL scripts from Robert L. Campbell's website:

wget -r -np -nd --level=1 -A .py \
    http://pldserver1.biochem.queensu.ca/~rlc/work/pymol/

Render movie from PNG files (save as png2mpeg1.sh):

#!/bin/bash

set -e

usage="usage: $(basename $0) [-w width] [-f fps] [-b vbitrate] <indir> <outfile.mpeg>"
width=""
fps=25
vbitrate=16000

args="$(getopt w:f:b:h "$@")" || args="-h"
set -- $args

while [[ $# > 0 ]]; do
    case "$1" in
        --) shift; break ;;
        -w) width=$2; shift 2 ;;
        -f) fps=$2; shift 2 ;;
        -b) vbitrate=$2; shift 2 ;;
        -h) echo $usage; exit 1 ;;
        *) echo "argument error: $1"; exit 1 ;;
    esac
done

if [[ $# > 2 ]]; then
    echo "too many arguments: $3 ..."
    echo $usage
    exit 1
fi

indir="$1"
outfile="$2"

if [[ -z "$indir" ]]; then
    echo "error: indir missing"
    echo $usage
    exit 2
fi

if [[ -z "$outfile" ]]; then
    echo "error: outfile missing"
    echo $usage
    exit 3
fi

MENCODER="mencoder -quiet"
MPEG1ARGS="-mf type=png:fps=$fps -ovc lavc -forceidx -noskip \
    -of rawvideo -mpegopts format=mpeg1 \
    -lavcopts vcodec=mpeg1video:vbitrate=$vbitrate:vhq:trell:keyint=25"

if [[ -n "$width" ]]; then
    MPEG1ARGS="-zoom -xy $width -sws 9 $MPEG1ARGS"
fi

pattern="mf://$indir/*.png"
$MENCODER "$pattern" $MPEG1ARGS:vpass=1 -o /dev/null
$MENCODER "$pattern" $MPEG1ARGS:vpass=2 -o "$outfile"

Export movie with transparent background

# make transparent pngs
set opaque_background, off
set ray_trace_frames
mpng foo

# create movie file (use codec "qtrle" or "png")
system ffmpeg -i foo%04d.png -vcodec qtrle foo.mov 

# clean up
system rm -f foo????.png

load_mtz_cctbx: Load MTZ files with a CCTBX wrapper (3 files)

1) ~/bin/mtz2ccp4.sh

#!/bin/bash
export PATH=/opt/ccp4/ccp4-6.5/bin:$PATH
exec cctbx.python ~/bin/mtz2ccp4.py "$@"

2) ~/bin/mtz2ccp4.py

#!/opt/ccp4/ccp4-6.5/bin/cctbx.python

import os
import sys
import tempfile

def mtz2ccp4maps(filename, prefix='map'):
    '''
Creates a temporary directory and dumps all maps from the given MTZ file
into this directory as CCP4 maps files. Returns the path of the temporary
directory.
    '''
    from iotbx.reflection_file_reader import any_reflection_file

    hkl_in = any_reflection_file(file_name=filename)

    temp_dir = tempfile.mkdtemp()

    for i_map, array in enumerate(hkl_in.as_miller_arrays()):
        if array.is_complex_array():
            fft_map = array.fft_map(resolution_factor=0.25).apply_sigma_scaling()
            map_filename = os.path.join(temp_dir,
                    prefix + '_' + '_'.join(array.info().labels) + '.ccp4')
            fft_map.as_ccp4_map(file_name=map_filename)

    return temp_dir

# print the name of the temporary directory to standard output
print mtz2ccp4maps(*sys.argv[1:])

3) ~/.pymolrc.py

@cmd.extend
def load_mtz_cctbx(filename, prefix=''):
    '''
DESCRIPTION

    Load all maps from an MTZ file, using the mtz2ccp4.sh wrapper which
    uses iotbx (cctbx).
    '''
    import subprocess
    import glob
    import shutil

    if not prefix:
        prefix = os.path.basename(filename).rpartition('.')[0]

    outdir = subprocess.Popen([os.path.expanduser('~/bin/mtz2ccp4.sh'),
        filename, prefix], stdout=subprocess.PIPE).stdout.readlines()[0].strip()

    for mapfilename in glob.glob(os.path.join(outdir, '*.ccp4')):
        cmd.load(mapfilename)

    shutil.rmtree(outdir)

ccmutate

@cmd.extend
def ccmutate(code, selection='??sele|?pk1', sculpt=1):
    '''
DESCRIPTION

    Mutate selected residue.

ARGUMENTS

    code = str: 3-letter PDBeChem chemical component identifier

    selection = str: single residue selection {default: pk1 or sele}

    sculpt = 0/1: try to adopt conformation of replaced sidechain, followed
    by relaxation using sculpting {default: 1}

EXAMPLE

    fetch 1ubq, async=0
    ccmutate 0HG, resi 24

SEE ALSO

    fetch ..., type=cc
    wizard mutagenesis
    '''
    code = code.upper()

    tmp_sele = cmd.get_unused_name('_sele')
    tmp_frag = cmd.get_unused_name('_frag')
    tmp_Nnbr = cmd.get_unused_name('_Nnbr')
    tmp_back = cmd.get_unused_name('_back')
    tmp_tmpl = cmd.get_unused_name('_tmpl')
    tmp_sc_o = cmd.get_unused_name('_sc_o')
    tmp_sc_n = cmd.get_unused_name('_sc_n')

    try:
        cmd.select(tmp_sele, 'byres (' + selection + ')', 0)

        # check input selection
        if cmd.count_atoms('name CA & ?' + tmp_sele) != 1:
            raise pymol.CmdException('selection must include exactly one residue')
        if cmd.count_atoms('name N+CA+C & ?' + tmp_sele) != 3:
            raise pymol.CmdException("selected residue doesn't have N+CA+C atoms")

        # PDBeChem fragment
        cmd.fetch(code, tmp_frag, type='cc', zoom=0)

        # check if fragment is amino acid
        if cmd.count_atoms('name N+CA+C & ?' + tmp_frag) != 3:
            raise pymol.CmdException("residue '%s' doesn't have N+CA+C atoms" % (code))

        # only keep hydrogens if target also has hydrogens
        if cmd.count_atoms('hydro & ?' + tmp_sele) == 0:
            cmd.remove('hydro & ?' + tmp_frag)

        # update residue name for old residue
        cmd.alter(tmp_sele, 'resn = ' + repr(code))

        # superpose fragment on backbone
        cmd.align(
                'name N+CA+C & ?' + tmp_frag,
                'name N+CA+C & ?' + tmp_sele)

        # extra N bonds, like in PRO
        cmd.select(tmp_Nnbr, 'neighbor (name N & ?' + tmp_frag + ')', 0)

        # backbone selection
        cmd.select(tmp_back, 'name CA+C+O+N+OXT', 0)
        cmd.select(tmp_back, 'hydro & neighbor ?' + tmp_back, 0, merge=1)

        # remove complementary atoms
        cmd.remove(           '?' + tmp_frag + ' & ' + tmp_back)
        cmd.extract(tmp_tmpl, '?' + tmp_sele + ' & !' + tmp_back, zoom=0)

        if cmd.count_atoms(tmp_frag):
            # attach new sidechain
            cmd.fuse('name CB & ?' + tmp_frag, 'name CA & ?' + tmp_sele, mode=1, move=0)
            cmd.unpick()

            # new atom selections
            cmd.select(tmp_sc_n, '(byres ?' + tmp_sele + ') & !?' + tmp_sele, 0)
            cmd.select(tmp_sc_o, '?' + tmp_sc_n + ' like ?' + tmp_tmpl, 0)

            # extra N bonds, like in PRO
            if cmd.count_atoms(tmp_Nnbr):
                cmd.bond('?' + tmp_sc_n + ' like ?' + tmp_Nnbr, 'name N & ?' + tmp_sele)

            # adopt old conformation, if possible
            if int(sculpt):
                model = cmd.get_object_list('?' + tmp_sele)[0]
                cmd.protect(model)
                cmd.deprotect('?' + tmp_sc_n + ' & !?' + tmp_sc_o)
                cmd.sculpt_activate(model)
                if cmd.count_atoms(tmp_sc_o):
                    cmd.update(tmp_sc_o, tmp_tmpl)
                    cmd.set('sculpt_field_mask', 63) # local geom + vdw
                    cmd.sculpt_iterate(model, cycles=100)
                    cmd.deprotect(tmp_sc_o)
                cmd.set('sculpt_field_mask', 0xff) # all
                cmd.sculpt_iterate(model, cycles=200)
                cmd.set('sculpt_field_mask', 31) # local geom
                cmd.sculpt_iterate(model, cycles=200)

    finally:
        cmd.delete(tmp_sele)
        cmd.delete(tmp_frag)
        cmd.delete(tmp_Nnbr)
        cmd.delete(tmp_back)
        cmd.delete(tmp_tmpl)
        cmd.delete(tmp_sc_o)
        cmd.delete(tmp_sc_n)