Python molecule.Molecule类代码示例

def read_molecule(r):
    """ (reader) -> Molecule

    Read a single molecule from r and return it,
    or return None to signal end of file.
    """

    # If there isn't another line, we're at the end of the file.
    line = r.readline()
    if not line:
        return None

    # Name of the molecule: "COMPND   name"
    key, name = line.split()

    # Other lines are either "END" or "ATOM num kind x y z"
    molecule = Molecule(name)
    reading = True

    while reading:
        line = r.readline()
        if line.startswith('END'):
            reading = False
        else:
            key, num, kind, x, y, z = line.split()
            molecule.add(Atom(num, kind, float(x), float(y), float(z)))

    return molecule

def parse_lines():
    #f = open(filename,'r')
    isParsing = True
    mol = None
    line_counter = 0
    n_atoms = 0
    mol = Molecule("")
    title = ""
    to_angs = False
    datas = []
    for line in fileinput.input():
        data = line.split()
        if isParsing:
            line_counter += 1
            if line_counter == 1: n_atoms = int(data[0])
            if line_counter == 2:
                uline = line.upper()
                to_angs = "AU" in uline
            if line_counter > 2:
                toangs = 1.0
                if to_angs: toangs = 0.529177249
                char = data[0]
                data = map(float, data[1:])
                atom = Atom(char, data[0]*toangs, data[1]*toangs, data[2]*toangs)
                if mol is not None:
                    mol.addAtom(atom)
    #f.close()
    return mol

class LeastSquaresCharges(LeastSquaresBasic):

    def __init__(self, data):
        self.molecule = Molecule(atoms=data['atoms'])
        self.grid = Grid(data)
        self.setA()
        reference = self.grid.get_properties(property_name=data['property'],\
                theory='%s_%s' % (data['theory'], data['basis']))
        LeastSquaresBasic.__init__(self, A=self.A, b=reference)

    def setA(self):
        n_p = len(self.grid.points)
        n_s = len(self.molecule.sites_noneq)
        self.A =  np.zeros((n_p, n_s))
        for i in xrange(n_p):
            proton = Site(coordinates=self.grid.points[i].coordinates, name='H+',index=1)
            for j, name in enumerate(self.molecule.sites_names_noneq):
                for site in self.molecule.get_sites(name=name):
                    self.A[i,j] += 1/site.distance_to(proton)

    def setA_fast(self):
        grid_coordinates = self.grid.get_coordinates()
        sites_coordinates = self.molecule.get_coordinates()
        self.A = fast.set_inversed(grid_coordinates, sites_coordinates, \
                len(self.molecule.sites_names_eq), self.molecule.sym_sites)

    @property
    def charges(self):
        charges = {}
        for q, name in zip(self.solution, self.molecule.sites_names):
            charges[name] = q
        return charges

def createMolecules(arr,params):
    retarr=[]
    for elem in arr:
        mol=Molecule()
        mol.RDMol=elem["RDMol"]
        retarr.append(mol)
    return retarr

 def test_unit_conversion_symmetry(self):
     """ 
     Does converting back and forth between bohrs and angstroms introduce and compound rounding errors? Each 
     operation should exactly reverse its counterpart.
     """
             
     with open('../../extra-files/molecule.xyz', 'r') as file1:
         mol1 = Molecule(file1.read())
     mol2 = mol1.copy()
     for count in range(10000) :
         mol2.to_bohr()
         mol2.to_angstrom()
     self.assertEqual(mol1.geom, mol2.geom)

    def test_unit_conversion_accuracy(self):
        """ 
        1.0 Angstrom is approximately 1.889725989 Bohr. Is this conversion (and its reverse) carried out correctly?
        """

        with open('../../extra-files/molecule.xyz', 'r') as file1:
            mol1 = Molecule(file1.read())
        mol2 = mol1.copy()
        mol2.to_bohr()
        for i in range(mol1.natom) :
            self.assertAlmostEqual(mol1.geom[i][0] * 1.889725989, mol2.geom[i][0])
            self.assertAlmostEqual(mol1.geom[i][1] * 1.889725989, mol2.geom[i][1])
            self.assertAlmostEqual(mol1.geom[i][2] * 1.889725989, mol2.geom[i][2])

 def read_cartesian(self,fname):
     molc=Molecule()
     f = open(fname,'r')
     while (True):
         str = f.readline()
         if (str==''):
             break
         sp = str.split()
         if (len(sp) != 4):
             continue
         molc.add_atom(Atom(sp[0].lower(),0,float(sp[1]),float(sp[2]),float(sp[3])))
     self.mol = molc
     return  molc

def solve(mol_path, hess_path):

    with open(mol_path, 'r') as f:
        molecule = Molecule(f.read())
        molecule.to_angstrom()

    with open(hess_path, 'r') as f:
        str = (f.read()).replace("\n",";")
        while str[-1] == ';' :
            str = str[:-1]
        mat = matrix(str)

    output_frequencies(molecule, mat)

 def __init__(self, data):
     self.molecule = Molecule(atoms=data['atoms'])
     self.grid = Grid(data)
     self.setA()
     reference = self.grid.get_properties(property_name=data['property'],\
             theory='%s_%s' % (data['theory'], data['basis']))
     LeastSquaresBasic.__init__(self, A=self.A, b=reference)

 def __init__(self, name=None, ff=None, sym=False, g_settings=None, m_settings=None, multipoles=None, ls_jobs=['w', 'non-w']):
     GM_logger.info("Creating MoleculeOnGrid instance")
     Molecule.__init__(self, name=name, ff=ff, sym=sym, settings=m_settings, multipoles=multipoles)
     self.grid = grid.Grid(molecule_name=name, settings=g_settings)
     self.grid_coordinates = self.grid.get_coordinates()
     reference = {}
     reference['non-w'] = self.grid.get_properties(property_name='esp', theory='reference')
     try:
         reference['w']= reference['non-w']*np.sqrt(self.grid.weights)
     except AttributeError:
         pass
     self.A = {}
     self.set_A(ls_jobs)
     self.LS = {}
     for key, A in self.A.items():
         self.LS[key] = LeastSquares(name=key, A=A, b=reference[key])

    def test_random_molecule(self):
        m = Molecule.random_molecule(12, 16)

        assert m.__class__ == Molecule

        assert m.r.x <= 12
        assert m.r.y <= 16
        assert m.v.x <= 3
        assert m.v.y <= 4

    def test_copy_function(self):
        """ 
        Does the copy function correctly initialize all variables of a new molecule? Also, is the new molecule
        truly a different object? (ie: changing properties of the original does not affect the copy and vice versa)?
        """
                
        with open('../../extra-files/molecule.xyz', 'r') as file1:
            mol1 = Molecule(file1.read())
        mol2 = mol1.copy()
        self.assertEqual(mol1.units, mol2.units, 'checking units')
        self.assertEqual(mol1.natom, mol2.natom, 'checking natom')
        self.assertEqual(mol1.labels, mol2.labels, 'checking labels')
        self.assertEqual(mol1.masses, mol2.masses, 'checking masses')
        self.assertEqual(mol1.charges, mol2.charges, 'checking charges')
        self.assertEqual(mol1.geom, mol2.geom, 'checking geometry')

        mol2.to_bohr()
        self.assertNotEqual(mol1.units, mol2.units)
        self.assertNotEqual(mol1.geom, mol2.geom)

    def test_random_molecules(self):
        ms = Molecule.random_molecules(100, 12, 16)

        assert len(ms) == 100
        for m in ms:
            assert m.__class__ == Molecule

            assert m.r.x <= 12
            assert m.r.y <= 16
            assert m.v.x <= 3
            assert m.v.y <= 4

def filterByActivity(arr,params):
    retarr=[]
    discard=[0,0,0,0,0,0,0,0]
    mols={}
    for line in arr:
        #print(line)
        if not line[u'bioactivity_type'] == u'IC50':
            discard[0]+=1
            continue
        if line[u"target_confidence"]<7:
            discard[1]+=1
            continue
        if line[u"units"]!=u"nM":
            discard[2]+=1
            continue
        value=0.0
        try:
            value=float(line[u"value"])
        except ValueError:
            discard[3]+=1
            continue
        if value > 1000:
            discard[4]+=1
            continue
        mol=Molecule(line[u'parent_cmpd_chemblid'])

        #mol.id=line[u'ingredient_cmpd_chemblid']
        if(mol.id in mols):
            discard[5]+=1
            continue
        else:
            mols[mol.id]=True           
        mol.pIC=-math.log(value)
        retarr.append(mol)
        discard[6]+=1
    print(discard)
    return retarr

def harvest_zmat(zmat):
    """Parses the contents of the Cfour ZMAT file into array and
    coordinate information. The coordinate info is converted into a
    rather dinky Molecule (no fragment, but does read charge, mult,
    unit). Return qcdb.Molecule. Written for findif zmat* where
    geometry always Cartesian and Bohr.

    """
    zmat = zmat.splitlines()[1:]  # skip comment line
    Nat = 0
    readCoord = True
    isBohr = ''
    charge = 0
    mult = 1
    molxyz = ''
    cgeom = []
    for line in zmat:
        if line.strip() == '':
            readCoord = False
        elif readCoord:
            lline = line.split()
            molxyz += line + '\n'
            Nat += 1
        else:
            if line.find('CHARGE') > -1:
                idx = line.find('CHARGE')
                charge = line[idx + 7:]
                idxc = charge.find(',')
                if idxc > -1:
                    charge = charge[:idxc]
                charge = int(charge)
            if line.find('MULTIPLICITY') > -1:
                idx = line.find('MULTIPLICITY')
                mult = line[idx + 13:]
                idxc = mult.find(',')
                if idxc > -1:
                    mult = mult[:idxc]
                mult = int(mult)
            if line.find('UNITS=BOHR') > -1:
                isBohr = ' bohr'

    molxyz = '%d%s\n%d %d\n' % (Nat, isBohr, charge, mult) + molxyz
    mol = Molecule.init_with_xyz(molxyz, no_com=True, no_reorient=True, contentsNotFilename=True)

    return mol