Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
def test_atomtyping(self, mol_name):
top_filename = '{}.top'.format(mol_name)
gro_filename = '{}.gro'.format(mol_name)
top_path = os.path.join(self.resource_dir, top_filename)
gro_path = os.path.join(self.resource_dir, gro_filename)
structure = pmd.gromacs.GromacsTopologyFile(top_path, xyz=gro_path,
parametrize=False)
known_opls_types = [atom.type for atom in structure.atoms]
print("Typing {}...".format(mol_name))
typed_structure = OPLSAA.apply(structure)
generated_opls_types = list()
for i, atom in enumerate(typed_structure.atoms):
message = ('Found multiple or no OPLS types for atom {} in {} ({}): {}\n'
'Should be atomtype: {}'.format(
i, mol_name, top_filename, atom.type, known_opls_types[i]))
assert atom.type, message
generated_opls_types.append(atom.type)
both = zip(generated_opls_types, known_opls_types)
def by_name(cls, basename, parameterized=False):
mol2_path, top_path, gro_path = cls.available()[basename]
if parameterized:
# structure = pmd.gromacs.GromacsTopologyFile(top_path, xyz=gro_path,
# parametrize=parameterized)
structure = pmd.gromacs.GromacsTopologyFile(top_path, xyz=gro_path,
parametrize=False)
structure.title = structure.title.replace(' GAS', '')
else:
# structure = pmd.gromacs.GromacsTopologyFile(pdb_path, xyz=pdb_path,
# parametrize=False)
structure = pmd.load_file(mol2_path, structure=True)
structure.title = structure.title.replace(' GAS', '')
return structure
#Import ParmEd
import parmed
#Require version 2.0.4 or later of ParmEd, otherwise ParmEd corrupts [ defaults ] section in GROMACS topologies with incorrect FudgeLJ/FudgeQQ
try:
ver = parmed.version
except:
oldParmEd = Exception('ERROR: ParmEd is too old, please upgrade to 2.0.4 or later')
raise oldParmEd
if ver < (2,0,4):
raise RuntimeError("ParmEd is too old, please upgrade to 2.0.4 or later")
#Read AMBER to ParmEd object
structure = parmed.amber.AmberParm( in_prmtop, in_crd )
#Make GROMACS topology
gromacs_topology = parmed.gromacs.GromacsTopologyFile.from_structure( structure )
#Write
parmed.gromacs.GromacsTopologyFile.write( gromacs_topology, out_top )
if precision == None:
parmed.gromacs.GromacsGroFile.write( gromacs_topology, out_gro )
else:
parmed.gromacs.GromacsGroFile.write( gromacs_topology, out_gro, precision = precision )
return out_top, out_gro
#Require version 2.0.4 or later of ParmEd, otherwise ParmEd corrupts [ defaults ] section in GROMACS topologies with incorrect FudgeLJ/FudgeQQ
try:
ver = parmed.version
except:
oldParmEd = Exception('ERROR: ParmEd is too old, please upgrade to 2.0.4 or later')
raise oldParmEd
if ver < (2,0,4):
raise RuntimeError("ParmEd is too old, please upgrade to 2.0.4 or later")
#Read AMBER to ParmEd object
structure = parmed.amber.AmberParm( in_prmtop, in_crd )
#Make GROMACS topology
gromacs_topology = parmed.gromacs.GromacsTopologyFile.from_structure( structure )
#Write
parmed.gromacs.GromacsTopologyFile.write( gromacs_topology, out_top )
if precision == None:
parmed.gromacs.GromacsGroFile.write( gromacs_topology, out_gro )
else:
parmed.gromacs.GromacsGroFile.write( gromacs_topology, out_gro, precision = precision )
return out_top, out_gro
Merging should work fine on general topology files. However, if molecule_names and/or molecule_numbers are provided the input topologies must contain single-molecule, single-residue topologies. Molecule_numbers are used to replicate the input topologies via ParmEd, and molecule_names are used to change the residue names from the input topologies.
"""
#PRELIMINARIES
N_tops = len( input_topologies )
#Check for obvious input problems - do we have the right number of everything, do all the input files exist
if molecule_numbers != None:
assert len( molecule_numbers ) == N_tops, "Must provide same number of molecule numbers as topology files."
for filenm in input_topologies:
assert( os.path.isfile( filenm )), "Error: Can't find input file %s provided to merge_topologies." % filenm
#WORK ON TOPOLOGIES
tops = []
for filenm in input_topologies:
top = parmed.gromacs.GromacsTopologyFile( filenm )
tops.append( top )
#List numbers of each molecule if not provided
if molecule_numbers == None:
molecule_numbers = [ 1] * N_tops
#Check that we've been provided with the correct number of molecule_names if any
if molecule_names != None:
total_molecules = 0
for topnr in range(N_tops):
total_molecules += len( tops[topnr].residues )
assert total_molecules == len( molecule_names ), "Must provide a number of molecule names equal to your total number of residues, but you have %s and %s, respectively." % ( len( molecule_names), total_molecules )
#Rename residues
ctr = 0
for topnr in range(N_tops):