Commit edd39ec2 by Peter Eastman

Added option to add a water box

parent 907e7708
<script> <script>
function validateForm() { function validateForm() {
if (document.getElementById("addCheckbox").checked) { if (document.getElementById("addHydrogensCheckbox").checked) {
var ph = document.getElementById("phfield").value; var ph = document.getElementById("phfield").value
if (ph > 0 && ph < 14) if (!(ph > 0 && ph < 14)) {
return true;
alert("pH must be a number between 0 and 14.") alert("pH must be a number between 0 and 14.")
return false; return false
} }
}
if (document.getElementById("addWaterCheckbox").checked) {
var xsize = document.getElementById("boxxfield").value
var ysize = document.getElementById("boxyfield").value
var zsize = document.getElementById("boxzfield").value
var strength = document.getElementById("ionicstrengthfield").value
if (!(xsize > 0 && ysize > 0 && zsize > 0)) {
alert("Box dimensions must be positive numbers.")
return false
}
if (!(strength >= 0)) {
alert("Ionic strength must be a nonnegative number.")
return false
}
}
return true
} }
</script> </script>
<form id="mainform" method="post" action="/" onsubmit="return validateForm()"> <form id="mainform" method="post" action="/">
<h1>Delete Heterogens</h1> <h1>Delete Heterogens</h1>
A heterogen is any residue other than a standard amino acid or nucleotide. Do you want to delete heterogens? A heterogen is any residue other than a standard amino acid or nucleotide. Do you want to delete heterogens?
<p> <p>
...@@ -21,12 +36,31 @@ A heterogen is any residue other than a standard amino acid or nucleotide. Do y ...@@ -21,12 +36,31 @@ A heterogen is any residue other than a standard amino acid or nucleotide. Do y
<h1>Add Missing Hydrogens</h1> <h1>Add Missing Hydrogens</h1>
Add missing hydrogen atoms? Add missing hydrogen atoms?
<p> <p>
<input type="checkbox" id="addCheckbox" name="add" checked> Add hydrogens appropriate for pH <input type="text" id="phfield" name="ph" value="7.0" size="5"> <input type="checkbox" id="addHydrogensCheckbox" name="addhydrogens" onchange="document.getElementById('phfield').disabled=!document.getElementById('addHydrogensCheckbox').checked" checked> Add hydrogens appropriate for pH <input type="text" id="phfield" name="ph" value="7.0" size="5">
<h1>Add Water</h1>
Add a water box surrounding the model?
<p>
<input type="checkbox" id="addWaterCheckbox" name="addwater" onchange="enableControls(document.getElementById('waterInputs'), document.getElementById('addWaterCheckbox').checked)"> Add water
<div id="waterInputs" style="margin-left:50px">
<table style="text-align:right">
<tr><td>Box dimensions (in nm):</td><td><input type="text" id="boxxfield" name="boxx" size="8"></td><td><input type="text" id="boxyfield" name="boxy" size="8"></td><td><input type="text" id="boxzfield" name="boxz" size="8"></td></tr>
%s
</table>
<p>
Ions will be added to neutralize the model. You can optionally add more ions based on a desired bulk ionic strength.
<p>
<table>
<tr><td style="text-align:right">Ionic strength (molar):</td><td><input type="text" id="ionicstrengthfield" name="ionicstrength" size="8" value="0.0"></td></tr>
<tr><td style="text-align:right">Positive ion:</td><td><select name="positiveion"><option value="Cs">Cs+</option><option value="K">K+</option><option value="Li">Li+</option><option value="Na" selected>Na+</option><option value="Rb">Rb+</option></select></td></tr>
<tr><td style="text-align:right">Negative ion:</td><td><select name="negativeion"><option value="Cl" selected>Cl-</option><option value="Br">Br-</option><option value="F">F-</option><option value="I">I-</option></select></td></tr>
</table>
</div>
<p> <p>
<input type="button" value="Continue" onclick="submitWithSpinner()"/> <input type="button" value="Continue" onclick="if (validateForm()) submitWithSpinner()"/>
</form> </form>
<script> <script>
setCurrentStep(6) setCurrentStep(6)
enableControls(document.getElementById('waterInputs'), false)
</script> </script>
</body> </body>
<html> <html>
...@@ -16,6 +16,14 @@ function setCheckboxes(selected) { ...@@ -16,6 +16,14 @@ function setCheckboxes(selected) {
for (var i = 0; i < checkboxes.length; i++) for (var i = 0; i < checkboxes.length; i++)
checkboxes[i].checked = selected checkboxes[i].checked = selected
} }
function enableControls(container, enabled) {
tags = ["input", "select"]
for (var i = 0; i < tags.length; i++) {
controls = container.getElementsByTagName(tags[i])
for (var j = 0; j < controls.length; j++)
controls[j].disabled = !enabled
}
}
function submitWithSpinner() { function submitWithSpinner() {
animateSpinner.animateIndex = 0 animateSpinner.animateIndex = 0
setTimeout(function() {animateSpinner();document.getElementById('overlay').style.visibility='visible';}, 1500) setTimeout(function() {animateSpinner();document.getElementById('overlay').style.visibility='visible';}, 1500)
...@@ -83,7 +91,7 @@ body {font-family:sans-serif} ...@@ -83,7 +91,7 @@ body {font-family:sans-serif}
<p id="step3">Add Residues</p> <p id="step3">Add Residues</p>
<p id="step4">Convert Residues</p> <p id="step4">Convert Residues</p>
<p id="step5">Add Heavy Atoms</p> <p id="step5">Add Heavy Atoms</p>
<p id="step6">Add Hydrogens</p> <p id="step6">Add Water/Hydrogens</p>
<p id="step7">Save File</p> <p id="step7">Save File</p>
</span> </span>
<p/> <p/>
\ No newline at end of file
...@@ -36,6 +36,7 @@ import simtk.openmm.app as app ...@@ -36,6 +36,7 @@ import simtk.openmm.app as app
import simtk.unit as unit import simtk.unit as unit
from simtk.openmm.app.internal.pdbstructure import PdbStructure from simtk.openmm.app.internal.pdbstructure import PdbStructure
from simtk.openmm.app.element import hydrogen, oxygen from simtk.openmm.app.element import hydrogen, oxygen
from simtk.openmm.app.forcefield import NonbondedGenerator
import numpy as np import numpy as np
import numpy.linalg as lin import numpy.linalg as lin
import sys import sys
...@@ -401,7 +402,7 @@ class PDBFixer(object): ...@@ -401,7 +402,7 @@ class PDBFixer(object):
# Create a System for energy minimizing it. # Create a System for energy minimizing it.
res = list(newTopology.residues()) res = list(newTopology.residues())
forcefield = self._createForceField(newTopology) forcefield = self._createForceField(newTopology, False)
system = forcefield.createSystem(newTopology) system = forcefield.createSystem(newTopology)
# Set any previously existing atoms to be massless, they so won't move. # Set any previously existing atoms to be massless, they so won't move.
...@@ -479,7 +480,20 @@ class PDBFixer(object): ...@@ -479,7 +480,20 @@ class PDBFixer(object):
self.topology = modeller.topology self.topology = modeller.topology
self.positions = modeller.positions self.positions = modeller.positions
def _createForceField(self, newTopology): def addSolvent(self, boxSize, positiveIon='Na+', negativeIon='Cl-', ionicStrength=0*unit.molar):
modeller = app.Modeller(self.topology, self.positions)
forcefield = self._createForceField(self.topology, True)
modeller.addSolvent(forcefield, boxSize=boxSize, positiveIon=positiveIon, negativeIon=negativeIon, ionicStrength=ionicStrength)
self.topology = modeller.topology
self.positions = modeller.positions
def _createForceField(self, newTopology, water):
if water:
forcefield = app.ForceField('amber10.xml', 'tip3p.xml')
nonbonded = [f for f in forcefield._forces if isinstance(f, NonbondedGenerator)][0]
radii = {'H':0.198, 'Li':0.203, 'C':0.340, 'N':0.325, 'O':0.299, 'F':0.312, 'Na':0.333, 'Mg':0.141,
'P':0.374, 'S':0.356, 'Cl':0.347, 'K':0.474, 'Br':0.396, 'Rb':0.527, 'I':0.419, 'Cs':0.605}
else:
forcefield = app.ForceField(os.path.join(os.path.dirname(__file__), 'soft.xml')) forcefield = app.ForceField(os.path.join(os.path.dirname(__file__), 'soft.xml'))
# The Topology may contain residues for which the ForceField does not have a template. # The Topology may contain residues for which the ForceField does not have a template.
...@@ -513,6 +527,14 @@ class PDBFixer(object): ...@@ -513,6 +527,14 @@ class PDBFixer(object):
if element not in atomTypes: if element not in atomTypes:
atomTypes[element] = (typeName, 0.0, element) atomTypes[element] = (typeName, 0.0, element)
forcefield._atomTypes[typeName] = atomTypes[element] forcefield._atomTypes[typeName] = atomTypes[element]
if water:
# Select a reasonable vdW radius for this atom type.
if element.symbol in radii:
sigma = radii[element.symbol]
else:
sigma = 0.5
nonbonded.typeMap[typeName] = (0.0, sigma, 0.0)
indexInResidue[atom.index] = len(template.atoms) indexInResidue[atom.index] = len(template.atoms)
template.atoms.append(app.ForceField._TemplateAtomData(atom.name, typeName, element)) template.atoms.append(app.ForceField._TemplateAtomData(atom.name, typeName, element))
for atom in residue.atoms(): for atom in residue.atoms():
......
import simtk.openmm.app as app import simtk.openmm.app as app
from simtk.openmm.app.internal.pdbstructure import PdbStructure from simtk.openmm.app.internal.pdbstructure import PdbStructure
import simtk.unit as unit
from pdbfixer import PDBFixer, substitutions, proteinResidues, dnaResidues, rnaResidues from pdbfixer import PDBFixer, substitutions, proteinResidues, dnaResidues, rnaResidues
import uiserver import uiserver
import webbrowser import webbrowser
...@@ -69,9 +70,15 @@ def addHydrogensPageCallback(parameters, handler): ...@@ -69,9 +70,15 @@ def addHydrogensPageCallback(parameters, handler):
fixer.removeHeterogens(False) fixer.removeHeterogens(False)
elif heterogens == 'water': elif heterogens == 'water':
fixer.removeHeterogens(True) fixer.removeHeterogens(True)
if 'add' in parameters: if 'addhydrogens' in parameters:
pH = float(parameters.getfirst('ph')) pH = float(parameters.getfirst('ph'))
fixer.addMissingHydrogens(pH) fixer.addMissingHydrogens(pH)
if 'addwater' in parameters:
boxSize = (float(parameters.getfirst('boxx')), float(parameters.getfirst('boxy')), float(parameters.getfirst('boxz')))*unit.nanometer
ionicStrength = float(parameters.getfirst('ionicstrength'))*unit.molar
positiveIon = parameters.getfirst('positiveion')+'+'
negativeIon = parameters.getfirst('negativeion')+'-'
fixer.addSolvent(boxSize, positiveIon, negativeIon, ionicStrength)
displaySaveFilePage() displaySaveFilePage()
def saveFilePageCallback(parameters, handler): def saveFilePageCallback(parameters, handler):
...@@ -171,11 +178,16 @@ def displayMissingAtomsPage(): ...@@ -171,11 +178,16 @@ def displayMissingAtomsPage():
if residue in fixer.missingTerminals: if residue in fixer.missingTerminals:
atoms.extend(atom for atom in fixer.missingTerminals[residue]) atoms.extend(atom for atom in fixer.missingTerminals[residue])
table += ' <tr><td>%d</td><td>%s %d</td><td>%s</td></tr>\n' % (residue.chain.index+1, residue.name, indexInChain[residue], ', '.join(atoms)) table += ' <tr><td>%d</td><td>%s %d</td><td>%s</td></tr>\n' % (residue.chain.index+1, residue.name, indexInChain[residue], ', '.join(atoms))
uiserver.setContent(header+loadHtmlFile("addHeavyAtoms.html")% table) uiserver.setContent(header+loadHtmlFile("addHeavyAtoms.html") % table)
def displayAddHydrogensPage(): def displayAddHydrogensPage():
uiserver.setCallback(addHydrogensPageCallback) uiserver.setCallback(addHydrogensPageCallback)
uiserver.setContent(header+loadHtmlFile("addHydrogens.html")) dimensions = ""
if fixer.topology.getUnitCellDimensions() is not None:
dimensions = "<tr><td>Crystallographic unit cell:</td><td>%.3f</td><td>%.3f</td><td>%.3f</td></tr>" % fixer.topology.getUnitCellDimensions().value_in_unit(unit.nanometer)
sizeRange = tuple(max((pos[i] for pos in fixer.positions))-min((pos[i] for pos in fixer.positions)) for i in range(3))
dimensions += "<tr><td>Box containing all atoms:</td><td>%.3f</td><td>%.3f</td><td>%.3f</td></tr>" % tuple(x.value_in_unit(unit.nanometer) for x in sizeRange)
uiserver.setContent(header+loadHtmlFile("addHydrogens.html") % dimensions)
def displaySaveFilePage(): def displaySaveFilePage():
uiserver.setCallback(saveFilePageCallback) uiserver.setCallback(saveFilePageCallback)
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment