com.actelion.research.chem.contrib.HydrogenHandler Maven / Gradle / Ivy
Go to download
Show more of this group Show more artifacts with this name
Show all versions of openchemlib Show documentation
Show all versions of openchemlib Show documentation
Open Source Chemistry Library
package com.actelion.research.chem.contrib;
import com.actelion.research.chem.Molecule;
import com.actelion.research.chem.StereoMolecule;
public class HydrogenHandler {
final static double cos30 = Math.cos(30.0 / 360.0 * Math.PI * 2);
final static double sin30 = Math.sin(30.0 / 360.0 * Math.PI * 2);
final static double cos330 = Math.cos(-30.0 / 360.0 * Math.PI * 2);
final static double sin330 = Math.sin(-30.0 / 360.0 * Math.PI * 2);
final static double cos60 = Math.cos(60.0 / 360.0 * Math.PI * 2);
final static double sin60 = Math.sin(60.0 / 360.0 * Math.PI * 2);
final static double cos300 = Math.cos(-60.0 / 360.0 * Math.PI * 2);
final static double sin300 = Math.sin(-60.0 / 360.0 * Math.PI * 2);
final static double cos45 = Math.cos(45.0 / 360.0 * Math.PI * 2);
final static double sin45 = Math.sin(45.0 / 360.0 * Math.PI * 2);
final static double cos315 = Math.cos(-45.0 / 360.0 * Math.PI * 2);
final static double sin315 = Math.sin(-45.0 / 360.0 * Math.PI * 2);
final static double H_BOND_RATIO=0.7;
// changed from method name from expandHydrogensAtAtom(); TLS 9.Nov.2015
public static void addImplicitHydrogens(StereoMolecule molecule, int iAtom) {
final int atomicNo = 1;
// Following line let crash the "cleanup"
// molecule.setAtomMarker(iAtom, true);
int nHydrogens = molecule.getImplicitHydrogens(iAtom);
double x = molecule.getAtomX(iAtom);
double y = molecule.getAtomY(iAtom);
switch (nHydrogens) {
case 1: {
int nNeighbours = molecule.getConnAtoms(iAtom);
// this is for salts like HCL
double dx, dy;
int iNew;
if (nNeighbours == 0) {
double distMin = Double.MAX_VALUE;
int iMin = -1;
for (int atm = 0; atm < molecule.getAllAtoms(); atm++) {
if (atm == iAtom)
continue;
double deltaX = x - molecule.getAtomX(atm);
double deltaY = y - molecule.getAtomY(atm);
double currDist = Math.sqrt((deltaX * deltaX) + (deltaY * deltaY));
if (distMin > currDist) {
distMin = currDist;
iMin = atm;
}
}
dx = x - molecule.getAtomX(iMin);
dy = y - molecule.getAtomY(iMin);
} else {
dx = x - molecule.getAtomX(molecule.getConnAtom(iAtom, 0));
dy = y - molecule.getAtomY(molecule.getConnAtom(iAtom, 0));
}
if (nNeighbours == 1) {
// rotate by 45
iNew = molecule.addAtom((float)(x + cos45 * dx + sin45 * dy),(float)(y - sin45 * dx + cos45 * dy));
} else if (nNeighbours == 2) {
dx = x - 0.5 * (molecule.getAtomX(molecule.getConnAtom(iAtom, 0)) + molecule.getAtomX(molecule.getConnAtom(iAtom, 1)));
dy = y - 0.5 * (molecule.getAtomY(molecule.getConnAtom(iAtom, 0)) + molecule.getAtomY(molecule.getConnAtom(iAtom, 1)));
iNew = molecule.addAtom((float)(x + dx),(float)(y + dy));
} else if (nNeighbours == 3) { // case of existing 3 atoms. We need to put the hydrogen next to the chiral bond !
// looking for a stereobond ...
int stereoAtom1=molecule.getConnAtom(iAtom, 0);
for (int i=1; i<3; i++) {
int bond=molecule.getConnBond(iAtom, i);
if ((molecule.getBondType(bond)== Molecule.cBondTypeDown) ||
(molecule.getBondType(bond)==Molecule.cBondTypeUp)) {
stereoAtom1=molecule.getConnAtom(iAtom, i);
}
}
// should we put the hydrogen next to the chiral bond ? It depends ot the sum of the 3 vectors
double angle1=Math.abs(Molecule.getAngleDif(
Molecule.getAngle(molecule.getAtomX(iAtom),
molecule.getAtomY(iAtom),
molecule.getAtomX(molecule.getConnAtom(iAtom, 0)),
molecule.getAtomY(molecule.getConnAtom(iAtom, 0))
),
Molecule.getAngle(molecule.getAtomX(iAtom),
molecule.getAtomY(iAtom),
molecule.getAtomX(molecule.getConnAtom(iAtom, 1)),
molecule.getAtomY(molecule.getConnAtom(iAtom, 1))
)
));
double angle2=Math.abs(Molecule.getAngleDif(
Molecule.getAngle(molecule.getAtomX(iAtom),
molecule.getAtomY(iAtom),
molecule.getAtomX(molecule.getConnAtom(iAtom, 0)),
molecule.getAtomY(molecule.getConnAtom(iAtom, 0))
),
Molecule.getAngle(molecule.getAtomX(iAtom),
molecule.getAtomY(iAtom),
molecule.getAtomX(molecule.getConnAtom(iAtom, 2)),
molecule.getAtomY(molecule.getConnAtom(iAtom, 2))
)
));
double angle3=Math.abs(Molecule.getAngleDif(
Molecule.getAngle(molecule.getAtomX(iAtom),
molecule.getAtomY(iAtom),
molecule.getAtomX(molecule.getConnAtom(iAtom, 1)),
molecule.getAtomY(molecule.getConnAtom(iAtom, 1))
),
Molecule.getAngle(molecule.getAtomX(iAtom),
molecule.getAtomY(iAtom),
molecule.getAtomX(molecule.getConnAtom(iAtom, 2)),
molecule.getAtomY(molecule.getConnAtom(iAtom, 2))
)
));
// System.out.println(angle1+" - "+angle2+" - "+angle3);
boolean normal=true;
if ((angle1>angle2) && (angle1>angle3)) {
if ((angle2+angle3) < Math.PI) {
normal=false;
dx = x - 0.5 * (molecule.getAtomX(molecule.getConnAtom(iAtom, 0)) + molecule.getAtomX(molecule.getConnAtom(iAtom, 1)));
dy = y - 0.5 * (molecule.getAtomY(molecule.getConnAtom(iAtom, 0)) + molecule.getAtomY(molecule.getConnAtom(iAtom, 1)));
}
} else if ((angle2>angle1) && (angle2>angle3)) {
if ((angle1+angle3) < Math.PI) {
normal=false;
dx = x - 0.5 * (molecule.getAtomX(molecule.getConnAtom(iAtom, 0)) + molecule.getAtomX(molecule.getConnAtom(iAtom, 2)));
dy = y - 0.5 * (molecule.getAtomY(molecule.getConnAtom(iAtom, 0)) + molecule.getAtomY(molecule.getConnAtom(iAtom, 2)));
}
} else {
if ((angle1+angle2) < Math.PI) {
normal=false;
dx = x - 0.5 * (molecule.getAtomX(molecule.getConnAtom(iAtom, 1)) + molecule.getAtomX(molecule.getConnAtom(iAtom, 2)));
dy = y - 0.5 * (molecule.getAtomY(molecule.getConnAtom(iAtom, 1)) + molecule.getAtomY(molecule.getConnAtom(iAtom, 2)));
}
}
if (normal) {
// what is the closest atom to the atom 1 ?
int stereoAtom2=molecule.getConnAtom(iAtom, 0);
double distance=Double.MAX_VALUE;
for (int i=0; i<3; i++) {
int atom=molecule.getConnAtom(iAtom, i);
if (atom!=stereoAtom1) {
double currentDistance=Math.pow(molecule.getAtomX(iAtom)-molecule.getAtomX(atom),2)+Math.pow(molecule.getAtomY(iAtom)-molecule.getAtomY(atom),2);
if (currentDistance-1) {
molecule.addBond(iAtom, iNew, Molecule.cBondTypeSingle);
} else {
molecule.addBond(iAtom, iNew, Molecule.cBondTypeUp);
}
iNew = molecule.addAtom((float)(x + cos330 * dx - sin330 * dy),(float)(y + sin330 * dx + cos330 * dy));
molecule.setAtomicNo(iNew, atomicNo);
molecule.addBond(iAtom, iNew, Molecule.cBondTypeSingle);
/*
iNew = molecule.addAtom(x + dx1, y + dy1);
molecule.setAtomicNo(iNew, atomicNo);
molecule.addBond(iAtom, iNew, Molecule.cBondTypeUp);
iNew = molecule.addAtom(x + dx2, y + dy2);
molecule.setAtomicNo(iNew, atomicNo);
molecule.addBond(iAtom, iNew, Molecule.cBondTypeSingle);
*/
} else {
for (int iHydrogen = 0; iHydrogen < nHydrogens; iHydrogen++) {
int iNew = molecule.addAtom((float)x, (float)y);
molecule.setAtomicNo(iNew, atomicNo);
molecule.addBond(iAtom, iNew, Molecule.cBondTypeSingle);
}
}
break;
case 3:
// CH3, NH3 and the like
{
int iNew;
double dx;
double dy;
if (molecule.getConnAtoms(iAtom)>0) {
dx = (x - molecule.getAtomX(molecule.getConnAtom(iAtom, 0)))*H_BOND_RATIO;
dy = (y - molecule.getAtomY(molecule.getConnAtom(iAtom, 0)))*H_BOND_RATIO;
// backwards
iNew = molecule.addAtom((float)(x + dx),(float)(y + dy));
molecule.setAtomicNo(iNew, atomicNo);
molecule.addBond(iAtom, iNew, Molecule.cBondTypeSingle);
// to the left
iNew = molecule.addAtom((float)(x - dy),(float)(y + dx));
molecule.setAtomicNo(iNew, atomicNo);
molecule.addBond(iAtom, iNew, Molecule.cBondTypeSingle);
// to the right
iNew = molecule.addAtom((float)(x + dy),(float)(y - dx));
molecule.setAtomicNo(iNew, atomicNo);
molecule.addBond(iAtom, iNew, Molecule.cBondTypeSingle);
} else {
dx = molecule.getAverageBondLength(true);
dy = molecule.getAverageBondLength(true);
// forward
iNew = molecule.addAtom((float)(x + dx),(float)(y + dy));
molecule.setAtomicNo(iNew, atomicNo);
molecule.addBond(iAtom, iNew, Molecule.cBondTypeSingle);
// backwards top
iNew = molecule.addAtom((float)(x - dy*cos60),(float)(y + dx*sin60));
molecule.setAtomicNo(iNew, atomicNo);
molecule.addBond(iAtom, iNew, Molecule.cBondTypeSingle);
// backwards bottom
iNew = molecule.addAtom((float)(x - dy*cos60),(float)(y - dx*sin60));
molecule.setAtomicNo(iNew, atomicNo);
molecule.addBond(iAtom, iNew, Molecule.cBondTypeSingle);
}
}
break;
default:
//
{
for (int iHydrogen = 0; iHydrogen < nHydrogens; iHydrogen++) {
int iNew = molecule.addAtom((float)x,(float)y);
molecule.setAtomicNo(iNew, atomicNo);
molecule.addBond(iAtom, iNew, Molecule.cBondTypeSingle);
}
break;
}
}
}
// changed from method name from expandHydrogens(); TLS 9.Nov.2015
public static void addImplicitHydrogens(StereoMolecule molecule) {
molecule.ensureHelperArrays(Molecule.cHelperNeighbours);
int nAtomsBefore = molecule.getAtoms();
for (int iAtom = 0; iAtom < nAtomsBefore; iAtom++) {
addImplicitHydrogens(molecule, iAtom);
}
}
/* This seems strange: explicit hydrogens are counted twice.
* I suggest using one of the following two functions: TLS 9.Nov.2015
*
public static int getNumberOfExplicitHydrogens(StereoMolecule molecule) {
molecule.ensureHelperArrays(Molecule.cHelperNeighbours);
return molecule.getAllAtoms() - molecule.getAtoms();
}
public static int getNumberOfAllHydrogens(StereoMolecule molecule) {
molecule.ensureHelperArrays(Molecule.cHelperNeighbours);
int nbHydrogens = 0;
for (int iAtom = 0; iAtom < molecule.getAtoms(); iAtom++)
nbHydrogens += molecule.getAllHydrogens(iAtom);
return nbHydrogens;
}
*/
public static int getNumberOfHydrogens(StereoMolecule molecule) {
molecule.ensureHelperArrays(Molecule.cHelperNeighbours);
int nbHydrogens = 0;
for (int iAtom = 0; iAtom < molecule.getAllAtoms(); iAtom++) {
if (molecule.getAtomicNo(iAtom) == 1)
nbHydrogens++;
else
nbHydrogens += molecule.getPlainHydrogens(iAtom);
}
return nbHydrogens;
}
}