Skip to content

Commit

Permalink
Debugging prepare_mols and create_bonds
Browse files Browse the repository at this point in the history
  • Loading branch information
choglass committed Apr 9, 2024
1 parent faf1d25 commit dfc66d9
Show file tree
Hide file tree
Showing 2 changed files with 3,690 additions and 484 deletions.
35 changes: 19 additions & 16 deletions cell2mol/classes.py
Original file line number Diff line number Diff line change
Expand Up @@ -146,7 +146,7 @@ def reset_charge(self):
if hasattr(self,"totcharge"): delattr(self,"totcharge")
if hasattr(self,"atomic_charges"): delattr(self,"atomic")
if hasattr(self,"smiles"): delattr(self,"smiles")
if hasattr(self,"rdkit_mol"): delattr(self,"rdkit_mol")
if hasattr(self,"rdkit_obj"): delattr(self,"rdkit_obj")
if hasattr(self,"poscharges"): delattr(self,"poscharges")
for a in self.atoms:
a.reset_charge()
Expand Down Expand Up @@ -296,7 +296,7 @@ def get_possible_cs(self, debug: int=0):

############
def create_bonds(self, debug: int=0):
if not hasattr(self,"rdkit_obj"): self.get_parent("cell").assign_charges()
# if not hasattr(self,"rdkit_obj"): self.get_parent("cell").assign_charges()
for idx, atom in enumerate(self.atoms):
# Security Check. Confirms that the labels are the same
if debug >= 2: print("BUILD BONDS: atom", idx, atom.label)
Expand All @@ -309,14 +309,16 @@ def create_bonds(self, debug: int=0):
bond_startatom = b.GetBeginAtomIdx()
bond_endatom = b.GetEndAtomIdx()
bond_order = b.GetBondTypeAsDouble()
if debug >= 2: print("BUILD BONDS: bond", bond_startatom, bond_endatom, bond_order, self.atoms[bond_startatom].label, self.atoms[bond_endatom].label, self.rdkit_obj.GetAtomWithIdx(bond_endatom).GetSymbol())
if debug >= 2: print("BUILD BONDS: bond", bond_startatom, bond_endatom, bond_order,
self.atoms[bond_startatom].label, self.atoms[bond_endatom].label,
self.rdkit_obj.GetAtomWithIdx(bond_endatom).GetSymbol())
if (self.subtype == "ligand") and (bond_startatom >= self.natoms or bond_endatom >= self.natoms):
continue
else:
if self.atoms[bond_endatom].label != self.rdkit_obj.GetAtomWithIdx(bond_endatom).GetSymbol():
pass
if debug >= 1:
print("Error with Bond EndAtom", self.atoms[bond_endatom].label, self.rdkit_obj.GetAtomWithIdx(bond_endatom).GetSymbol())
print("Error with Bond EndAtom", self.atoms[bond_endatom].label,
self.rdkit_obj.GetAtomWithIdx(bond_endatom).GetSymbol())
else:
if bond_endatom == idx:
start = bond_endatom
Expand Down Expand Up @@ -480,7 +482,7 @@ def correct_smiles(self):
if not hasattr(self,"smiles"): self.smiles = []
if not self.iscomplex: return self.smiles
for lig in self.ligands:
lig.smiles, lig.rdkit_mol = correct_smiles_ligand(lig)
lig.smiles, lig.rdkit_obj = correct_smiles_ligand(lig)
self.smiles.append(lig.smiles)

###############
Expand Down Expand Up @@ -780,16 +782,21 @@ def __init__(self, atom1: object, atom2: object, bond_order: int=1):
self.atom1 = atom1
self.atom2 = atom2
self.order = bond_order
self.distance = np.linalg.norm(np.array(atom1.coord) - np.array(atom2.coord))

def __repr__(self):
to_print = ""
to_print += f'------------- Cell2mol BOND Object --------------\n'
to_print += f' Version = {self.version}\n'
to_print += f' Type = {self.type}\n'
to_print += f' Version = {self.version}\n'
to_print += f' Type = {self.type}\n'
idx1 = self.atom1.get_parent_index("molecule")
idx2 = self.atom2.get_parent_index("molecule")
to_print += f' Molecule Atom 1 = {idx1}\n'
to_print += f' Molecule Atom 2 = {idx2}\n'
to_print += f' Bond Order = {self.order}\n'
to_print += f' Molecule Atom 1 label = {self.atom1.label}\n'
to_print += f' Molecule Atom 2 label = {self.atom2.label}\n'
to_print += f' Molecule Atom 1 index = {idx1}\n'
to_print += f' Molecule Atom 2 index = {idx2}\n'
to_print += f' Bond Order = {self.order}\n'
to_print += f' Distance = {round(self.distance,3)}\n'
to_print += '----------------------------------------------------\n'
return to_print

Expand Down Expand Up @@ -1440,14 +1447,13 @@ def create_bonds(self, debug: int=0):
count = 0
for met in mol.metals:
isconnected = at.check_connectivity(met, debug=debug)
# isconnected = check_connectivity(at, met)
if isconnected:
newbond = bond(at, met, 0.5)
at.add_bond(newbond)
met.add_bond(newbond)
count += 1
if count != at.mconnec:
print(f"CELL.CREATE_BONDS: error creating bonds for atom: \n{atom}\n of ligand: \n{lig}\n")
print(f"CELL.CREATE_BONDS: error creating bonds for atom: \n{at}\n of ligand: \n{lig}\n")
print(f"CELL.CREATE_BONDS: count differs from atom.mconnec: {count}, {at.mconnec}")

# Adds Metal-Metal Bonds, with an arbitrary 0.5 order:
Expand All @@ -1459,13 +1465,10 @@ def create_bonds(self, debug: int=0):
for jdx, met2 in enumerate(mol.metals):
if idx <= jdx: continue
isconnected = met1.check_connectivity(met2, debug=debug)
# isconnected = check_connectivity(met1, met2)
if isconnected:
newbond = bond(met1, met2, 0.5)
met1.add_bond(newbond)
met2.add_bond(newbond)
else :
pass

# Fourth part : correction smiles of ligands
mol.smiles_with_H = [lig.smiles for lig in mol.ligands]
Expand Down
Loading

0 comments on commit dfc66d9

Please sign in to comment.