# Copyright (c) 2009 Greg Pintilie - pintilie@mit.edu # Permission is hereby granted, free of charge, to any person obtaining a copy # of this software and associated documentation files (the "Software"), to deal # in the Software without restriction, including without limitation the rights # to use, copy, modify, merge, publish, distribute, sublicense, and/or sell # copies of the Software, and to permit persons to whom the Software is # furnished to do so, subject to the following conditions: # # The above copyright notice and this permission notice shall be included in # all copies or substantial portions of the Software. # # THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR # IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, # FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE # AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER # LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, # OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN # THE SOFTWARE. import chimera import os import os.path import Tkinter from CGLtk import Hybrid import VolumeData import _multiscale import MultiScale.surface import _surface import numpy import _contour import Matrix import VolumeViewer from sys import stderr from time import clock from axes import prAxes import regions import graph from Segger import dev_menus, timing, seggerVersion OML = chimera.openModels.list REG_OPACITY = 0.45 from segment_dialog import current_segmentation, segmentation_map def umsg ( txt ) : print txt status ( txt ) def status ( txt ) : txt = txt.rstrip('\n') msg.configure(text = txt) msg.update_idletasks() class FlexFit_Dialog ( chimera.baseDialog.ModelessDialog ): title = "Flexible Fitting (Segger v" + seggerVersion + ")" name = "flexfit" buttons = ( "Close") #help = '' def fillInUI(self, parent): self.group_mouse_mode = None tw = parent.winfo_toplevel() self.toplevel_widget = tw tw.withdraw() parent.columnconfigure(0, weight = 1) row = 0 menubar = Tkinter.Menu(parent, type = 'menubar', tearoff = False) tw.config(menu = menubar) f = Tkinter.Frame(parent) f.grid(column=0, row=row, sticky='ew') l = Tkinter.Label(f, text=' ') l.grid(column=0, row=row, sticky='w') if 1 : row += 1 ff = Tkinter.Frame(f) ff.grid(column=0, row=row, sticky='w') l = Tkinter.Label(ff, text=' Map:') l.grid(column=0, row=0, sticky='w') self.dmap = Tkinter.StringVar(parent) self.mb = Tkinter.Menubutton ( ff, textvariable=self.dmap, relief=Tkinter.RAISED ) self.mb.grid (column=1, row=0, sticky='we', padx=5) self.mb.menu = Tkinter.Menu ( self.mb, tearoff=0, postcommand=self.MapMenu ) self.mb["menu"] = self.mb.menu # set first visible map by default from VolumeViewer import Volume mlist = OML(modelTypes = [Volume]) for m in mlist : if m.display : self.dmap.set ( m.name ) break if 1 : row += 1 ff = Tkinter.Frame(f) ff.grid(column=0, row=row, sticky='w') l = Tkinter.Label(ff, text=' Structure:') l.grid(column=0, row=0, sticky='w') self.struc = Tkinter.StringVar(parent) self.strucMB = Tkinter.Menubutton ( ff, textvariable=self.struc, relief=Tkinter.RAISED ) self.strucMB.grid (column=1, row=0, sticky='we', padx=5) self.strucMB.menu = Tkinter.Menu ( self.strucMB, tearoff=0, postcommand=self.StrucMenu ) self.strucMB["menu"] = self.strucMB.menu from chimera import Molecule mlist = OML(modelTypes = [Molecule]) for m in mlist : if m.display : self.struc.set ( m.name ) self.cur_mol = m break row += 1 dummyFrame = Tkinter.Frame(f, relief='groove', borderwidth=1) Tkinter.Frame(dummyFrame).pack() dummyFrame.grid(row=row,column=0,columnspan=1, pady=5, sticky='we') if 1 : row += 1 ff = Tkinter.Frame(f) ff.grid(column=0, row=row, sticky='w') b = Tkinter.Label(ff, text="Visible chains: ") b.grid (column=0, row=0, sticky='w', padx=10) b = Tkinter.Button(ff, text="Join", command=self.JoinMols) b.grid (column=1, row=0, sticky='w', padx=10) b = Tkinter.Button(ff, text="Split", command=self.Split) b.grid (column=2, row=0, sticky='w', padx=10) row += 1 dummyFrame = Tkinter.Frame(f, relief='groove', borderwidth=1) Tkinter.Frame(dummyFrame).pack() dummyFrame.grid(row=row,column=0,columnspan=1, pady=5, sticky='we') if 1 : row += 1 ff = Tkinter.Frame(f) ff.grid(column=0, row=row, sticky='w') b = Tkinter.Label(ff, text="Molecular Dynamics Flexible Fitting (MDFF)") b.grid (column=0, row=0, sticky='w', padx=10) row += 1 ff = Tkinter.Frame(f) ff.grid(column=0, row=row, sticky='w') b = Tkinter.Label(ff, text=" density weight: ") b.grid (column=0, row=0, sticky='w', padx=0, pady=5) self.mdffDensityWeight = Tkinter.StringVar(f) self.mdffDensityWeight.set ( "0.3" ) e = Tkinter.Entry(ff, width=5, textvariable=self.mdffDensityWeight) e.grid(column=1, row=0, sticky='w', padx=5, pady=5) b = Tkinter.Label(ff, text=" min steps: ") b.grid (column=2, row=0, sticky='w', padx=0, pady=5) self.mdffMinSteps = Tkinter.StringVar(f) self.mdffMinSteps.set ( "10000" ) e = Tkinter.Entry(ff, width=6, textvariable=self.mdffMinSteps) e.grid(column=3, row=0, sticky='w', padx=5, pady=5) b = Tkinter.Label(ff, text=" # md steps: ") b.grid (column=4, row=0, sticky='w', padx=0, pady=5) self.mdffMDSteps = Tkinter.StringVar(f) self.mdffMDSteps.set ( "100000" ) e = Tkinter.Entry(ff, width=6, textvariable=self.mdffMDSteps) e.grid(column=5, row=0, sticky='w', padx=5, pady=5) row += 1 ff = Tkinter.Frame(f) ff.grid(column=0, row=row, sticky='w') b = Tkinter.Button(ff, text="Setup", command=self.ToMDFF) b.grid (column=0, row=0, sticky='w', padx=10, pady=5) b = Tkinter.Button(ff, text="Fix Sel", command=self.FixSel) b.grid (column=1, row=0, sticky='w', padx=5, pady=5) b = Tkinter.Button(ff, text="Target Sel", command=self.TargetSel) b.grid (column=2, row=0, sticky='w', padx=5, pady=5) b = Tkinter.Button(ff, text="Minimize", command=self.Minimize) b.grid (column=3, row=0, sticky='w', padx=5, pady=5) b = Tkinter.Button(ff, text="SSEs", command=self.SSEs) b.grid (column=4, row=0, sticky='w', padx=5, pady=5) b = Tkinter.Button(ff, text="Syms", command=self.PlaceSym) b.grid (column=5, row=0, sticky='w', padx=5, pady=5) row += 1 dummyFrame = Tkinter.Frame(f, relief='groove', borderwidth=1) Tkinter.Frame(dummyFrame).pack() dummyFrame.grid(row=row,column=0,columnspan=1, pady=5, sticky='we') if 1 : row += 1 ff = Tkinter.Frame(f) ff.grid(column=0, row=row, sticky='w') e = Tkinter.Label(ff, text='DCD File: ') e.grid(column=0, row=0, sticky='w', pady=5, padx=10) self.dcdFile = Tkinter.StringVar(ff) self.dcdFile.set ("_-step1.dcd") e = Tkinter.Entry(ff, width=40, textvariable=self.dcdFile) e.grid(column=1, row=0, sticky='we', pady=0) e = Tkinter.Button(ff, text="Load", command=self.LoadDCD) e.grid(column=2, row=0, sticky='we', pady=0) if 1 : row += 1 ff = Tkinter.Frame(f) ff.grid(column=0, row=row, sticky='w') l = Tkinter.Label(ff, text='Frames: ') l.grid(column=0, row=0, pady=5, sticky='w', padx=10) self.dcdFramesPerStep = Tkinter.StringVar(ff) self.dcdFramesPerStep.set ('1') e = Tkinter.Label(ff, text='Frames Per Step') e.grid(column=0, row=0, sticky='w', pady=0, padx=10) e = Tkinter.Entry(ff, width=5, textvariable=self.dcdFramesPerStep) e.grid(column=1, row=0, sticky='e', pady=0) self.movieOn = Tkinter.IntVar() c = Tkinter.Checkbutton(ff, text="Record frames", variable=self.movieOn) c.grid (column=2, columnspan=1, row=0, pady=5, sticky='w') row += 1 ff = Tkinter.Frame(f) ff.grid(column=0, row=row, sticky='w') e = Tkinter.Button(ff, text="Start", command=self.StartMovie) e.grid(column=0, row=0, sticky='e', pady=0, padx=10) b = Tkinter.Button(ff, text="<", command=self.Reverse) b.grid(column=1, row=0, padx=5, pady=5, sticky='we') b = Tkinter.Button(ff, text="|<", command=self.StepReverse) b.grid(column=2, row=0, padx=5, pady=5, sticky='we') b = Tkinter.Button(ff, text="||", command=self.Stop) b.grid(column=3, row=0, padx=5, pady=5, sticky='we') b = Tkinter.Button(ff, text=">|", command=self.StepForward) b.grid(column=4, row=0, padx=5, pady=5, sticky='we') b = Tkinter.Button(ff, text=">", command=self.Forward) b.grid(column=5, row=0, padx=5, pady=5, sticky='we') b = Tkinter.Button(ff, text="Last", command=self.ForwardEnd) b.grid(column=6, row=0, padx=5, pady=5, sticky='we') row += 1 dummyFrame = Tkinter.Frame(f, relief='groove', borderwidth=1) Tkinter.Frame(dummyFrame).pack() dummyFrame.grid(row=row,column=0,columnspan=7, pady=5, sticky='we') global msg msg = Tkinter.Label(parent, width = 60, anchor = 'w', justify = 'left', fg="red") msg.grid(column=0, row=row, sticky='ew') self.msg = msg row += 1 def MapMenu ( self ) : self.mb.menu.delete ( 0, 'end' ) # Clear menu from VolumeViewer import Volume mlist = OML(modelTypes = [Volume]) for m in mlist : self.mb.menu.add_radiobutton ( label=m.name, variable=self.dmap, command=lambda m=m: self.MapSelected(m) ) def MapSelected ( self, dmap ) : self.cur_dmap = dmap if dmap: dmap.display = True def StrucSelected ( self, mol ) : print "Selected " + mol.name self.cur_mol = mol def StrucMenu ( self ) : self.strucMB.menu.delete ( 0, 'end' ) # Clear menu from chimera import Molecule mlist = OML(modelTypes = [Molecule]) for m in mlist : self.strucMB.menu.add_radiobutton ( label=m.name, variable=self.struc, command=lambda m=m: self.StrucSelected(m) ) def JoinMols ( self ) : nmol = None chainIDs = "ABCDEFGHIJKLMNOPQRSTUVWXYZ0123456789abcdefghijklmnopqrstuvwxyz" #chainIDs = "BCDEFGHIJKLMNOPQRSTUVWXYZ0123456789abcdefghijklmnopqrstuvwxyz" #chainIDs = "CDEFGHIJKLMNOPQRSTUVWXYZ0123456789abcdefghijklmnopqrstuvwxyz" chainAt = 0 for mod in chimera.openModels.list() : if type(mod) == chimera.Molecule and mod.display == True : chains = Chains (mod) for chain in chains : try : newChainId = chainIDs[chainAt] except : continue chainAt += 1 print mod.name, chain, " --> ", newChainId from random import random clr = (random(), random(), random(), 1) nmol = copyMolChain ( nmol, mod, chain, newChainId, mod.openState.xform, clr ) nmol.name = "visible_new" chimera.openModels.add ( [nmol], noprefs = True ) nmol.openState.xform = chimera.Xform() #path = os.path.dirname ( refMod.data.path ) + os.path.sep #chimera.PDBio().writePDBfile ( [nmol], path + newName ) #print ( "Saved to " + path + newName ) #chimera.openModels.close ( [nmol] ) #chimera.openModels.open ( path + newName ) def Split ( self ) : nmol = None chain_at = 0 for mod in chimera.openModels.list() : if type(mod) == chimera.Molecule and mod.display == True : chains = Chains (mod) for chain in chains : print mod.name, chain, " --> ", chain_at from random import random clr = (random(), random(), random(), 1) nmol = copyMolChain ( None, mod, chain, "%d"%chain_at, mod.openState.xform, clr ) nmol.name = "_%d" % chain_at chimera.openModels.add ( [nmol], noprefs = True ) nmol.openState.xform = chimera.Xform() chain_at += 1 #path = os.path.dirname ( refMod.data.path ) + os.path.sep #chimera.PDBio().writePDBfile ( [nmol], path + newName ) #print ( "Saved to " + path + newName ) #chimera.openModels.close ( [nmol] ) #chimera.openModels.open ( path + newName ) def ToMDFF ( self ) : mol_name = self.struc.get() map_name = self.dmap.get() mol = getMod ( mol_name ) dmap = getMod ( map_name ) print "Fitting " + mol_name + " in " + map_name + " path: " + dmap.data.path import os mapName = os.path.splitext(map_name)[0] molName = os.path.splitext(mol_name)[0] mdir, mfile = os.path.split(dmap.data.path) dpath = mdir + "/__mdff__" + mapName + "__" + molName print " - path " + mdir print " - dpath " + dpath try : os.mkdir ( dpath ) print " - made path: " + dpath except : print " - path exists: " + dpath import shutil print dmap.data.path + " -> " + dpath + "/" + mfile shutil.copyfile ( dmap.data.path, dpath + "/" + mfile ) molpath = mol.openedAs[0] mdir, mfile = os.path.split(molpath) print molpath + " -> " + dpath + "/" + mfile shutil.copyfile ( molpath, dpath + "/" + mfile ) # these commands have to be executed in VMD to generate the NAMD files f = open ( dpath + "/vmd.txt", "w" ) f.write ("package require mdff\n") f.write ("mdff griddx -i "+mapName+".situs -o "+mapName+"-grid.dx\n") f.write ("mdff gridpdb -psf "+molName+"_autopsf.psf -pdb "+molName+"_autopsf.pdb -o "+molName+"_autopsf-grid.pdb\n") f.write ("package require ssrestraints\n") f.write ("mol new "+molName+"_autopsf.psf\n") f.write ("mol addfile "+molName+"_autopsf.pdb\n") f.write ("ssrestraints -psf "+molName+"_autopsf.psf -pdb "+molName+"_autopsf.pdb -o "+molName+"-eb.txt -hbonds\n") f.write ("cispeptide restrain -o "+molName+"-cis.txt\n") f.write ("chirality restrain -o "+molName+"-chi.txt\n") f.write ("mdff setup -o _ -psf "+molName+"_autopsf.psf -pdb "+molName+"_autopsf.pdb -griddx "+mapName+"-grid.dx -gridpdb "+molName+"_autopsf-grid.pdb -extrab {"+molName+"-eb.txt "+molName+"-cis.txt "+molName+"-chi.txt} -gscale " + self.mdffDensityWeight.get() + " -minsteps "+self.mdffMinSteps.get()+" -numsteps "+self.mdffMDSteps.get() +"\n") f.close() self.mdStep = 1 def FixSel ( self ) : if not hasattr(self, 'mdStep'): self.mdStep = 1 #umsg ( "Create vmd file first" ) #return mol_name = self.struc.get() map_name = self.dmap.get() mol = getMod ( mol_name ) dmap = getMod ( map_name ) print "Fitting " + mol_name + " in " + map_name + " path: " + dmap.data.path import os mapName = os.path.splitext(map_name)[0] molName = os.path.splitext(mol_name)[0] mdir, mfile = os.path.split(dmap.data.path) dpath = mdir + "/__mdff__" + mapName + "__" + molName print " - path " + mdir print " - dpath " + dpath print dmap.data.path + " -> " + dpath + "/" + mfile molpath = mol.openedAs[0] mdir, mfile = os.path.split(molpath) print molpath + " -> " + dpath + "/" + mfile namdFile = dpath + "/_-step%d.namd" % self.mdStep print "Namd file: " + namdFile newNamdFile = dpath + "/_-step%d_fx.namd" % (self.mdStep) print "New namd file: " + newNamdFile coorPath = dpath + "/_-step%d.coor" % self.mdStep molName = molName.replace ( "_autopsf", "" ) pdbPath = dpath + "/" + molName+"_autopsf.pdb" newPdbPath = dpath + "/" + molName+"_autopsf_fx.pdb" #oldPdbPath = dpath + "/" + ("_-step%d.b.coor" % self.mdStep) #newPdbPath = dpath + "/" + ("_-step%d.coor" % self.mdStep) #import shutil #try : # shutil.move (newPdbPath, ) #except : # print " - did not move original coor file" fp = open ( namdFile ) fpn = open ( newNamdFile, "w" ) for line in fp : if "set FIXPDB" in line : fpn.write ( "set FIXPDB " + molName+"_autopsf_fx.pdb" + "\n" ) fpn.write ( "set FIXCOL O\n" ) ## X, Y, Z, O, or B; 0 = not fixed elif "set FIXCOL" in line : pass elif "set INPUTNAME" in line : #fpn.write ( "set INPUTNAME step%d\n" % (self.mdStep) ) pass elif "set OUTPUTNAME" in line : #fpn.write ( "set INPUTNAME step%d\n" % (self.mdStep) ) fpn.write ( "set OUTPUTNAME _-step%d_fx\n" % (self.mdStep) ) pass else : fpn.write ( line ) # http://www.wwpdb.org/documentation/format33/sect9.html#ATOM selatoms = chimera.selection.currentAtoms() if len(selatoms) == 0 : umsg ("No selected atoms found") return print "Fixing %d atoms" % len (selatoms) fixAtoms = {} for at in selatoms : #print " res %d - %d %s bf:%.3f" % (at.residue.id.position, at.serialNumber, at.name, at.bfactor) fixAtoms [ "%d_%d" % (at.residue.id.position, at.serialNumber) ] = at.name #molid = selatoms[0].molecule.id #aats = chimera.selection.OSLSelection ( "#%d" % molid ).atoms() #allAtoms = {} #for at in aats : # allAtoms [ "%d_%d" % (at.residue.id.position, at.serialNumber) ] = at fp = open ( pdbPath ) fpn = open ( newPdbPath, "w" ) for line in fp : newLine = line if line[0:4] == "ATOM" : atSerial = int ( line[6:11] ) atResIdPos = int ( line[22:26] ) atName = line[12:16].replace(" ","") fkey = "%d_%d" % (atResIdPos, atSerial) # occupancy: 54:60 newLine = newLine[0:54] + "%6.2f" % 0.0 + newLine[60:] if fixAtoms.has_key ( fkey ) : fatName = fixAtoms [fkey] if fatName == atName : fixAtoms [fkey] = fatName + " OK" newLine = newLine[0:54] + "%6.2f" % 1.0 + newLine[60:] else : fixAtoms [fkey] = fatName + " name: " + atName fpn.write ( newLine ) for k,n in fixAtoms.iteritems () : if not "OK" in n : print k, " - ", n print "Wrote " + newPdbPath def TargetSel ( self ) : if not hasattr(self, 'mdStep'): self.mdStep = 1 #umsg ( "Create vmd file first" ) #return mol_name = self.struc.get() map_name = self.dmap.get() mol = getMod ( mol_name ) dmap = getMod ( map_name ) print "Fitting " + mol_name + " in " + map_name + " path: " + dmap.data.path import os mapName = os.path.splitext(map_name)[0] molName = os.path.splitext(mol_name)[0] mdir, mfile = os.path.split(dmap.data.path) dpath = mdir + "/__mdff__" + mapName + "__" + molName print " - path " + mdir print " - dpath " + dpath print dmap.data.path + " -> " + dpath + "/" + mfile molpath = mol.openedAs[0] mdir, mfile = os.path.split(molpath) print molpath + " -> " + dpath + "/" + mfile namdFile = dpath + "/_-step%d.namd" % self.mdStep print "Namd file: " + namdFile newNamdFile = dpath + "/_-step%d_tmd.namd" % (self.mdStep) print "New namd file: " + newNamdFile coorPath = dpath + "/_-step%d.coor" % self.mdStep molName = molName.replace ( "_autopsf", "" ) pdbPath = dpath + "/" + molName+"_target_autopsf.pdb" newPdbPath = dpath + "/" + molName+"_target_autopsf_conspdb.pdb" #oldPdbPath = dpath + "/" + ("_-step%d.b.coor" % self.mdStep) #newPdbPath = dpath + "/" + ("_-step%d.coor" % self.mdStep) #import shutil #try : # shutil.move (newPdbPath, ) #except : # print " - did not move original coor file" fp = open ( namdFile ) fpn = open ( newNamdFile, "w" ) for line in fp : if "set CONSPDB" in line : fpn.write ( "set CONSPDB " + molName+"_target_autopsf_conspdb.pdb" + "\n" ) fpn.write ( "set CONSCOL O\n" ) ## X, Y, Z, O, or B; 0 = not fixed elif "set CONSCOL" in line : pass elif "set INPUTNAME" in line : #fpn.write ( "set INPUTNAME step%d\n" % (self.mdStep) ) pass elif "set OUTPUTNAME" in line : #fpn.write ( "set INPUTNAME step%d\n" % (self.mdStep) ) fpn.write ( "set OUTPUTNAME _-step%d_tmd\n" % (self.mdStep) ) pass else : fpn.write ( line ) # http://www.wwpdb.org/documentation/format33/sect9.html#ATOM selatoms = chimera.selection.currentAtoms() if len(selatoms) == 0 : umsg ("No selected atoms found") return print "Fixing %d atoms" % len (selatoms) fixAtoms = {} for at in selatoms : keystr = "%d_%s" % (at.residue.id.position, at.name) print " res %d - %d %s bf:%.3f - [%s]" % (at.residue.id.position, at.serialNumber, at.name, at.bfactor, keystr) fixAtoms [ keystr ] = at.name #molid = selatoms[0].molecule.id #aats = chimera.selection.OSLSelection ( "#%d" % molid ).atoms() #allAtoms = {} #for at in aats : # allAtoms [ "%d_%d" % (at.residue.id.position, at.serialNumber) ] = at fp = open ( pdbPath ) fpn = open ( newPdbPath, "w" ) for line in fp : newLine = line if line[0:4] == "ATOM" : #atSerial = int ( line[6:11] ) atResIdPos = int ( line[22:26] ) atName = line[12:16].replace(" ","") fkey = "%d_%s" % (atResIdPos, atName) # occupancy: 54:60 newLine = newLine[0:54] + "%6.2f" % 0.0 + newLine[60:] if fixAtoms.has_key ( fkey ) : fatName = fixAtoms [fkey] fixAtoms [fkey] = fatName + " OK" newLine = newLine[0:54] + "%6.2f" % 1.0 + newLine[60:] fpn.write ( newLine ) print "selected atoms that were not used:" for k,n in fixAtoms.iteritems () : if not "OK" in n : print k, " - ", n print "Wrote " + newPdbPath def Minimize ( self ) : mol_name = self.struc.get() map_name = self.dmap.get() mol = getMod ( mol_name ) dmap = getMod ( map_name ) print "Fitting " + mol_name + " in " + map_name + " path: " + dmap.data.path import os mapName = os.path.splitext(map_name)[0] molName = os.path.splitext(mol_name)[0] mdir, mfile = os.path.split(dmap.data.path) dpath = mdir + "/__mdff__" + mapName + "__" + molName print " - path " + mdir print " - dpath " + dpath print dmap.data.path + " -> " + dpath + "/" + mfile molpath = mol.openedAs[0] mdir, mfile = os.path.split(molpath) print molpath + " -> " + dpath + "/" + mfile def SSEs ( self ) : mol_name = self.struc.get() mol = getMod ( mol_name ) start = None parts = [] print mol.name for r in mol.residues : #print "%d - " % r.id.position, #if r.isSheet or r.isStrand : #print "S", #elif r.isHelix : #print "H", if r.isSheet or r.isStrand or r.isHelix : if start == None : start = r.id.position #print " - start", else : if start != None : #print " - end", parts.append ( [start, r.id.position-1] ) start = None #print "" for start, end in parts : print start, end def LoadDCD(self) : global mdtMol, DCD from Trajectory.DCD.MDToolsMarch97 import md dcd = self.dcdFile.get() #psf = self.psfFile.get() mol_name = self.struc.get() mol = getMod ( mol_name ) map_name = self.dmap.get() dmap = getMod ( map_name ) import os mapName = os.path.splitext(map_name)[0] molName = os.path.splitext(mol_name)[0] mdir, mfile = os.path.split(dmap.data.path) mdffPath = mdir + "/__mdff__" + mapName + "__" + molName + "/" dpath = mdir + "/" molpath = mol.openedAs[0] moldir, molfile = os.path.split(molpath) print molpath + " -> " + dpath + "/" + mfile pdbFile = mdffPath + molName + "_autopsf.pdb" psfFile = mdffPath + molName + "_autopsf.psf" dcdFile = mdffPath + self.dcdFile.get() print ("Loading structure from %s" % pdbFile) opened = chimera.openModels.open( pdbFile ) self.dcdMol = opened[0] print ( "Opened %d models from %s:" % (len(opened), self.dcdMol.name) ) FindChains ( self.dcdMol ) #self.phipsi0 = self.GetPhiPsi () print "Reading PSF file" try: mdtMol = md.Molecule(psf = psfFile) finally: print "Done reading PSF file" print "Processing PSF file" try: mdtMol.buildstructure() finally: umsg ("Processed PSF file") print "Reading DCD header" try: DCD = md.DCD( dcdFile ) finally: print "Done reading DCD header" if DCD.numatoms != len(mdtMol.atoms): raise ValueError, "PSF has different number of atoms (%d) than DCD (%d)!" % (len(mdtMol.atoms), self.dcd.numatoms) umsg ( "DCD: %d atoms, %d frames" %(DCD.numatoms, DCD.numframes) ) ## for at in self.dcdMol.atoms[0:1] : ## print "at ", at.name ## api = 0 ## alist = DCD[0] ## while api < DCD.numatoms : ## va = alist[api] ## v = at.coord() - chimera.Point(va[0], va[1], va[2]) ## if v.length < .1 : ## print "found in dcd, index ", api ## ## api = api + 1 print "done" self.fat = 0; def StepReverse (self) : self.Step ( - int ( self.dcdFramesPerStep.get() ) ); def StepForward (self) : self.Step ( int ( self.dcdFramesPerStep.get() ) ); def ForwardEnd (self) : self.AxesSet = False self.Step ( DCD.numframes - self.fat - 2 ); def StartMovie (self) : print "Start" self.iNumFrames = 0 self.fat = 0 self.Step ( 0 ) self.fat = -1 return fname = self.psfFile.get() + ".pdb" print ("Getting structure from %s" % fname) tmol = chimera._openPDBModel ( fname )[0] for ri, res in enumerate( self.dcdMol.residues ) : #if self.buildSurfs.get() and res.id.chainId != 'X' : # res.ribbonDisplay = False #else : # res.ribbonDisplay = True dres = tmol.residues[ri] for dat in dres.atoms : try : at = res.atomsMap[ dat.name ][0] at.setCoord ( chimera.Point(dat.x, dat.y, dat.z) ) except : # print "Atom not found: ", dat.name pass if self.bAlignDomains.get() : self.AlignDomains () self.dcdMol.updateRibbonData() chimera.viewer.postRedisplay() self.iNumFrames = 0 if self.movieOn.get() : chimera.printer.saveImage ( "frames/%05d.png" % self.iNumFrames ) self.iNumFrames = self.iNumFrames + 1 def Step (self, incr, goto=None) : global mdtMol, DCD if goto != None : self.fat = goto else : self.fat = self.fat + incr if 0 : [phis_l, psis_l] = self.GetPhiPsi () if self.fat < DCD.numframes-1 and self.fat >= 0 : print "%d, " % self.fat, mdtMol.getmolframe( DCD[self.fat] ) for ri, res in enumerate( self.dcdMol.residues ) : #if self.buildSurfs.get() and res.id.chainId != 'X' : # res.ribbonDisplay = False #else : # res.ribbonDisplay = True dres = mdtMol.residues[ri] for dat in dres.atoms : try : at = res.atomsMap[ dat.name ][0] at.setCoord ( chimera.Point(dat.x, dat.y, dat.z) ) except : # print "Atom not found: ", dat.name pass else : print "---end---" self.Stop = 1 self.dcdMol.updateRibbonData() chimera.viewer.postRedisplay() try : i = self.iNumFrames except : self.iNumFrames = 0 if self.movieOn.get() : chimera.printer.saveImage ( "frames/%05d.png" % self.iNumFrames ) self.iNumFrames = self.iNumFrames + 1 def Forward (self) : self.StepForward () if self.Stop : self.Stop = 0 else : self._master.after ( 33, self.Forward ) def Reverse (self) : self.StepReverse () if self.Stop : self.Stop = 0 else : self._master.after ( 33, self.Reverse ) def Stop (self) : print "Stop" self.Stop = 1 def PlaceSym ( self ) : nmol = None for mod in chimera.openModels.list() : if type(mod) == chimera.Molecule and mod.display == True : nmol = mod break; nmol = None if hasattr ( self, "cur_mol" ) : nmol = self.cur_mol if nmol == None : umsg ( "Please select a Molecule model in the Structure: field" ) return if nmol != None : esym = "C12" print "Sym: ", esym, " to mol: ", nmol.name syms = [] #syms.append ( Matrix.identity_matrix () ) if ( esym == "C12" ) : ax = chimera.Vector ( 0, 0, 1 ) #ax = dmap.openState.xform.inverse().apply ( ax ) for i in range ( 1, 12 ) : print "%d:" % i rm1 = Matrix.rotation_transform ( (ax.x,ax.y,ax.z), float(i)*360.0/12.0, [0,0,0] ) print rm1 syms.append ( rm1 ) #syms.append ( Matrix.rotation_transform ( (1.0,0.0,0.0), 2.0*360.0/3.0 ) ) #centers, xyz, w = centers_and_points(dmap) #print " - center:", centers #ctf = Matrix.translation_matrix([-x for x in COM[0]]) #syms = Matrix.coordinate_transform_list(syms, ctf) smols = [] for si, sym in enumerate ( syms [0 : ] ) : T = numpy.array ( sym ) #print "\nSym %d\n" % si, T xf = chimera.Xform.xform ( T[0,0], T[0,1], T[0,2], T[0,3], T[1,0], T[1,1], T[1,2], T[1,3], T[2,0], T[2,1], T[2,2], T[2,3] ) mols = self.PlaceCopy ( [nmol], xf, si+1 ) for m in mols : m.openState.xform = nmol.openState.xform smols = smols + mols return smols def PlaceCopy(self, molecules, xf, sind, clr=None): new_mols = [] for molecule in molecules : mol = CopyMol ( molecule ) mol.name = os.path.splitext ( mol.name )[0] + "_sym_%d.pdb" % sind if clr : r, g, b, a = clr mclr = chimera.MaterialColor ( r, g, b, a ) else : mclr = molecule.residues[0].ribbonColor if xf != None : for at in mol.atoms : c = xf.apply ( at.coord() ) at.setCoord ( c ) new_mols.append ( mol ) for res in mol.residues : res.ribbonDisplay = True res.ribbonDrawMode = 2 res.ribbonColor = mclr for at in res.atoms : at.display = False mol.display = True chimera.openModels.add ( new_mols, noprefs = True ) return new_mols def show_dialog ( closeOld = True ): from chimera import dialogs d = dialogs.find ( "flexfit", create=False ) if d : if closeOld : d.toplevel_widget.update_idletasks () d.Close() d.toplevel_widget.update_idletasks () else : # is there a way to bring it to front? return d dialogs.register ( "flexfit", FlexFit_Dialog, replace = True) d = dialogs.find ( "flexfit", create=True ) d.toplevel_widget.update_idletasks () d.enter() return d def getMod ( name ) : import chimera mlist = chimera.openModels.list () for mol in mlist : if mol.name == name : return mol return None def Chains ( m ) : ct = {} for r in m.residues : ct[r.id.chainId] = 1 clist = ct.keys() clist.sort() return clist def CopyMol ( mol ) : nmol = chimera.Molecule() nmol.name = mol.name aMap = dict() from random import random as rand clr = ( rand(), rand(), rand() ) for res in mol.residues : nres = nmol.newResidue (res.type, chimera.MolResId(res.id.chainId, res.id.position)) # print "New res: %s %d" % (nres.id.chainId, nres.id.position) for at in res.atoms : nat = nmol.newAtom (at.name, chimera.Element(at.element.number)) aMap[at] = nat nres.addAtom( nat ) nat.setCoord ( at.coord() ) nat.drawMode = nat.Sphere nat.color = chimera.MaterialColor( clr[0], clr[1], clr[2], 1.0 ) nat.display = True nres.isHelix = res.isHelix nres.isHet = res.isHet nres.isSheet = res.isSheet nres.isStrand = res.isStrand nres.ribbonDisplay = True nres.ribbonDrawMode = 2 nres.ribbonColor = chimera.MaterialColor( clr[0], clr[1], clr[2], 1.0 ); for bond in mol.bonds : nb = nmol.newBond ( aMap[bond.atoms[0]], aMap[bond.atoms[1]] ) nb.display = nb.Smart return nmol def copyMolChain ( nmol, mol, cid, new_cid, xf, clr ) : if nmol == None : nmol = chimera.Molecule() nmol.name = "complex" print "Copying chain %s (%s) to %s" % (cid, mol.name, new_cid) aMap = dict() for res in mol.residues : if res.id.chainId == cid : nres = nmol.newResidue(res.type, chimera.MolResId(new_cid, res.id.position)) # print "New res: %s %d" % (nres.id.chainId, nres.id.position) for at in res.atoms : nat = nmol.newAtom (at.name, chimera.Element(at.element.number)) aMap[at] = nat nres.addAtom( nat ) if xf : nat.setCoord ( xf.apply( at.coord() ) ) else : nat.setCoord ( at.xformCoord() ) nat.drawMode = nat.Sphere nat.color = chimera.MaterialColor( clr[0], clr[1], clr[2], 1.0 ) nat.display = True nres.isHelix = res.isHelix nres.isHet = res.isHet nres.isSheet = res.isSheet nres.isStrand = res.isStrand nres.ribbonDisplay = True nres.ribbonDrawMode = 2 nres.ribbonColor = chimera.MaterialColor( clr[0], clr[1], clr[2], 1.0 ); for bond in mol.bonds : if (bond.atoms[0].residue.id.chainId == cid and bond.atoms[1].residue.id.chainId == cid) : nb = nmol.newBond ( aMap[bond.atoms[0]], aMap[bond.atoms[1]] ) nb.display = True return nmol #chain_colors = MakeDiffRandColors ( len(chain_ids) ) #color_indexes = range ( len(chain_ids) ) class Chain : def __init__ (self, mol) : self.molecule = mol #self.bbox = MyBBox () def FindMols () : nmols = len ( chimera.openModels.list() ) print "%d Molecules" % nmols colors = MakeColors ( nmols ) for i, m in enumerate ( chimera.openModels.list() ) : cid = chain_ids[i] nm = CopyChain ( m, ' ', cid, chimera.Xform(), colors[i] ) chains[cid] = Chain ( nm ) chimera.openModels.add ( [nm] ) chimera.openModels.close ( m ) chains = {} def FindChains ( mol=None ) : ct = {} if mol == None : m = chimera.openModels.list()[0] else : m = mol for r in m.residues: ct[r.id.chainId] = 1 clist = ct.keys() clist.sort() from random import random as rand print "%d chains:" % len(clist) for ci, cid in enumerate ( clist ) : #clr = chain_colors [ color_indexes[ci] clr = ( rand(), rand(), rand(), 1.0 ) print "- %s: clr(%.2f, %.2f, %.2f)" % (cid, clr[0], clr[1], clr[2]) new_chain = Chain (m) new_chain.color = clr chains[ cid ] = new_chain # color atoms rtypes = {} for r in m.residues : rtypes[r.type] = 1 clr = chains[ r.id.chainId ].color mclr = chimera.MaterialColor( clr[0], clr[1], clr[2], 1.0 ) r.ribbonDisplay = True r.ribbonDrawMode = 2 clr = chains[ r.id.chainId ].color r.ribbonColor = chimera.MaterialColor( *clr ); for at in r.atoms : at.display = False at.color = r.ribbonColor print "\nFound residues: ", for rt in rtypes.keys() : print "r.type=='%s' or " % rt, print "\n" def DispRibbon ( mol, clr ) : for r in mol.residues : r.ribbonDisplay = True r.ribbonDrawMode = 2 r.ribbonColor = chimera.MaterialColor( clr[0], clr[1], clr[2], 1.0 ); for at in r.atoms : at.color = r.ribbonColor at.display = False