[chimera-dev] implementing MMTK normal calculations in chimera

[Victor Muñoz]

Dear developpers,

I am trying to make an extension that allows normal mode calculations with
MMTK into chimera. I have several problems beginning with the following.

I have first used the MMTKinter.py and added a new function called NM_ffm
that calculates normal modes with a full atom force field. This works
perfectly (although I am not yet able to make an animation in chimera that
follows the vectors calculated...that will be a future question for you).
Then I added  another function (NM_mcm) that allows me to do NM calculations
with the mechanical constraint (see
http://dirac.cnrs-orleans.fr/MMTK/using-mmtk/mmtk-example-scripts/normal-modes/normal-modes-with-mechanical-constraints).
For this I need to define a subspace that requires a universe.protein object
(the functions are given at the end of the mail). Unfortunately, I can not
pass the molecular object from the chimera space as a universe.protein
object. We thought in putting self.universe.protein call into the
_makeUniverse function in MMTKinter.py but that did not work. We tried other
combinations but we now face the fact that we are blocked. May be we need
another function in order to define subspaces but yet we still don't know
how to pass chimera infos to MMTK.

Could you give me a hand so that I can define correctly subspaces on chimera
molecules?

All the best and thank in advance for any help,

Victor

*************
Computational Biotechnological Chemistry @transmet





####Function for normal modes calculations in full atom ff
    def NM_ffm(self, nsteps, stepsize=0.02,
            interval=None, action=None, **kw):
        from chimera import replyobj
        timestamp("_minimize")
        from MMTK import Units
        from MMTK.ForceFields.Amber import AmberData
        from MMTK.Minimization import SteepestDescentMinimizer
        from MMTK.Trajectory import LogOutput
        # from JD to make NM calculations
        from MMTK.NormalModes import NormalModes
        from chimera import runCommand as run
        from chimera import selection as sele
        from chimera import Point, Vector
        import os, sys
        # JD'end
        import sys
        if not interval:
            actions = []
        else:
            actions = [ LogOutput(sys.stdout, ["energy"],
                        interval, None, interval) ]
        kw["step_size"] = stepsize * Units.Ang
        minimizer = SteepestDescentMinimizer(self.universe,
                            actions=actions, **kw)
        if action is None or not interval:
            interval = None
        msg = "Initial energy: %f" % self.universe.energy()
        replyobj.status(msg)
        replyobj.info(msg)
        saveNormalizeName = AmberData._normalizeName
        AmberData._normalizeName = simpleNormalizeName
        remaining = nsteps
        while remaining > 0:
            timestamp(" minimize interval")
            if interval is None:
                realSteps = remaining
            else:
                realSteps = min(remaining, interval)
            minimizer(steps=realSteps)
            remaining -= realSteps
            if action is not None:
                action(self)
            timestamp(" finished %d steps" % realSteps)
            msg = "Finished %d of %d minimization steps" % (
                        nsteps - remaining, nsteps)
            replyobj.status(msg)
            replyobj.info(msg)
        replyobj.info("\n")
        #Calculate normal modes JD's trucho
        def despl(vec):
                run("sel #0")
                atoms=sele.currentAtoms()
                for a in atoms:
                       a.setCoord(a.xformCoord()-10*(vec))
        modes = NormalModes(self.universe)
        for mode in modes:
            print mode
            print mode.array
        L1=[]
        L2=[]
        L3=[]
        for mode in modes:
               y=6
                u=len(mode.array)
                while y<u:
                    L1.append(mode.array[y][0])
                    L2.append(mode.array[y][1])
                    L3.append(mode.array[y][2])
                    y+=1
        y=0
        run("sel:all")
        atoms=sele.currentAtoms()
        run("~sel")
        u=len(L1)
        a=0
        o=len(atoms)
        while y<u:
                   vec=Vector(float(L1[y]),float(L2[y]),float(L3[y]))
                   atoms[a].setCoord(atoms[a].coord()-(vec))
                   a+=1
            if a==o:
                a=0
                raw_input()
            else:
                pass
            y+=1
        #JD's end
        return modes
        AmberData._normalizeName = saveNormalizeName
        timestamp("end _minimize")

####Function for mechanical constraint
    def NM_mcm(self, nsteps, stepsize=0.02,
            interval=None, action=None, **kw):
        from MMTK import *
        from MMTK.Proteins import Protein
        from MMTK.ForceFields import Amber94ForceField
        from MMTK.NormalModes import NormalModes, SubspaceNormalModes
        from MMTK.Subspace import RigidMotionSubspace
        from MMTK.Minimization import ConjugateGradientMinimizer
        from MMTK.Trajectory import StandardLogOutput
        import numpy
        # Victor Modifications
        from chimera import replyobj
        timestamp("_minimize")
        from MMTK import Units
        from MMTK.ForceFields.Amber import AmberData
        from MMTK.Minimization import SteepestDescentMinimizer
        from MMTK.Trajectory import LogOutput
        from MMTK.NormalModes import NormalModes
        from chimera import runCommand as run
        from chimera import selection as sele
        from chimera import Point, Vector
        import os, sys
        #/end
        if not interval:
            actions = []
        else:
            actions = [ LogOutput(sys.stdout, ["energy"],
                        interval, None, interval) ]
        kw["step_size"] = stepsize * Units.Ang
        minimizer = SteepestDescentMinimizer(self.universe,
                            actions=actions, **kw)
        if action is None or not interval:
            interval = None
        msg = "Initial energy: %f" % self.universe.energy()
        replyobj.status(msg)
        replyobj.info(msg)
        saveNormalizeName = AmberData._normalizeName
        AmberData._normalizeName = simpleNormalizeName
        remaining = nsteps
        while remaining > 0:
            timestamp(" minimize interval")
            if interval is None:
                realSteps = remaining
            else:
                realSteps = min(remaining, interval)
            minimizer(steps=realSteps)
            remaining -= realSteps
            if action is not None:
                action(self)
            timestamp(" finished %d steps" % realSteps)
            msg = "Finished %d of %d minimization steps" % (
                        nsteps - remaining, nsteps)
            replyobj.status(msg)
            replyobj.info(msg)
        replyobj.info("\n")
        # Construct system
        #self.universe = InfiniteUniverse(Amber94ForceField())
        #self.universe.protein = self.universe
        # Set up the subspace: rigid-body translation and rotation for each
residue
        self.subspace = RigidMotionSubspace(self.universe,
self.universe.protein.residues())
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://plato.cgl.ucsf.edu/pipermail/chimera-dev/attachments/20091202/5593e7f5/attachment.html>