com.actelion.research.chem.sar.SARScaffold 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.sar;
import com.actelion.research.chem.*;
import com.actelion.research.chem.coords.CoordinateInventor;
import java.util.ArrayList;
import java.util.Arrays;
import java.util.TreeMap;
public class SARScaffold {
private StereoMolecule mQuery,mCore,mScaffold;
private int[] mCoreToQueryAtom,mQueryToCoreAtom;
private String mIDCodeWithRGroups, mIDCoordsWithRGroups;
private TreeMap mOldToNewMap;
private int mBridgeAtomRGroupCount;
private SARScaffoldGroup mScaffoldGroup;
private ExitVector[] mBridgeAtomExitVector;
private boolean[] mHasSeenSubstituentOnScaffold;
private int[] mSeenBondOrdersOnScaffold;
protected SARScaffold(StereoMolecule query, StereoMolecule core, int[] coreToQueryAtom, int[] queryToCoreAtom, SARScaffoldGroup scaffoldGroup) {
mQuery = query;
mCore = core;
mCoreToQueryAtom = coreToQueryAtom;
mQueryToCoreAtom = queryToCoreAtom;
mScaffoldGroup = scaffoldGroup;
mOldToNewMap = new TreeMap<>();
mBridgeAtomRGroupCount = -1;
analyzeAtomBridgeExitVectors(coreToQueryAtom);
mHasSeenSubstituentOnScaffold = new boolean[getExitVectorCount()];
mSeenBondOrdersOnScaffold = new int[getExitVectorCount()];
}
private void analyzeAtomBridgeExitVectors(int[] coreToQueryAtom) {
ArrayList evList = new ArrayList<>();
for (int atom=0; atom 1) {
hasExitPiBond = true;
break;
}
}
if (hasExitPiBond
&& !mol.isAtomStereoCenter(rootAtom)) {
for (int i=0; i 1)
|| (exitVector.getIndex() == 1 && mol.getConnBondOrder(rootAtom, i) == 1)))
return connAtom;
}
}
int count = 0;
for (int i=0; i getOldToNewMap() {
return mOldToNewMap;
}
protected int assignRGroupsToBridgeAtoms(int firstBridgeAtomRGroup) {
if (mBridgeAtomRGroupCount == -1) {
mBridgeAtomRGroupCount = 0;
for (ExitVector exitVector:mBridgeAtomExitVector)
if (exitVector.substituentVaries())
exitVector.setRGroupNo(++mBridgeAtomRGroupCount + firstBridgeAtomRGroup - 1);
}
return mBridgeAtomRGroupCount;
}
protected void addRGroupsToCoreStructure() {
mScaffold = new StereoMolecule(mCore);
mScaffold.ensureHelperArrays(Molecule.cHelperNeighbours);
double coreAVBL = mScaffold.getAverageBondLength();
int exitVectorCount = getExitVectorCount();
boolean[] closureCovered = new boolean[exitVectorCount];
for (int exitVectorIndex=0; exitVectorIndex attach an R group
ExitVector exitVector = getExitVector(exitVectorIndex);
if (exitVector.substituentVaries()) {
// But don't attach an R-group to one scaffold of a scaffold group,
// if that particular scaffold has never substituents at that position.
if (mHasSeenSubstituentOnScaffold[exitVectorIndex]) {
int rGroupNo = exitVector.getRGroupNo();
int newAtom = mScaffold.addAtom((rGroupNo<=3) ? 141 + rGroupNo : 125 + rGroupNo);
int bondType = calculateExitVectorCoordsAndBondType(exitVectorIndex, mScaffold.getAtomCoordinates(newAtom));
mScaffold.addBond(exitVector.getCoreAtom(mQueryToCoreAtom), newAtom, bondType);
}
}
else { // else => attach the non-varying substituent (if it is not null = 'unsubstituted')
if (!closureCovered[exitVectorIndex] && exitVector.getConstantSubstituent() != null) {
StereoMolecule substituent = new IDCodeParser(true).getCompactMolecule(exitVector.getConstantSubstituent());
// Substitutions, which connect back to the core fragment are decorated with atomicNo=0 atoms that
// carry a label with the respective exit vector index. Here we just copy those connection atoms,
// but mark the exit vector indexes as already attached (closureCovered) to avoid processing from the
// other end again. When all substituents are attached, we convert those labelled atoms into
// proper closure connections.
for (int atom=0; atom 0 && neighbourRank[index-1] > rank) {
neighbourRank[index] = neighbourRank[index-1];
neighbourBond[index] = neighbourBond[index-1];
neighbourAngle[index] = neighbourAngle[index-1];
index--;
}
neighbourRank[index] = rank;
neighbourBond[index] = connBond;
neighbourAngle[index] = mol.getBondAngle(rootAtom, connAtom);
totalNeighbourCount++;
}
if (totalNeighbourCount < 3 || totalNeighbourCount > 4)
return -1;
int stereoType = (mol.getBondType(stereoBond) == Molecule.cBondTypeUp) ? 2 : 1;
// Here we have one of the following neighbour counts in addition to the defined exitAtom:
// A: 1 neighbour that exists in core structure; 1 additional exit atom
// B: 2 neighbours that exists in core structure; no additional exit atom
// C: 2 neighbours that exists in core structure; 1 additional exit atom
// D: 3 neighbours that exists in core structure; no additional exit atom
// Case A and B:
// 3 non-H neighbours at stereo center.
// Thus, we don't need to change up/down bond type when shifting stereo bond to other neighbour.
// Cases C and D:
// 4 non-H neighbours at stereo center.
// Thus, we need to invert up/down bond type when shifting stereo bond to direct neighbour bond.
if (mol.getConnAtoms(rootAtom) == 4) {
if (otherExitAtomFound && stereoBond == neighbourBond[3]) {
if (areDirectNeighbours(neighbourAngle, 2, 3))
stereoType = 3 - stereoType;
}
else if (stereoBond != exitBond) {
// if the stereo bond is one of the core bonds and if the exit bond
int stereoBondIndex = -1;
for (int i=0; i angle1) && (angle[i] < angle2))
count++;
return count != 1;
}
/**
* If a substituent can be attached to a core structure in two distinguishable ways regarding
* stereo configuration, then this method calculates in a reproducible way a topicity (0 or 1)
* reflecting a given up/down bond stereo type and the bond angles from the stereo center to
* all neighbour atoms. The bond angle array is expected to contain sorted core neighbours first,
* followed by one or two exit vector angles. The last exit vector is considered to carry the
* stereo bond. If in reality the stereo bond is a different one, then the calling method is
* responsible to compensate for a stereo bond shift and/or to compensate for the second exit
* vector by potentially inverting stereoType.
* Topicity is defined as follows:
* - 1 neighbour atoms in core structure and 2 exit atoms:
* If walking from first atom (lowest relevant index) via root atom to first exit atom
* making a left turn, then an up-bond connecting second exit atom gives topicity=1
* - 2 neighbour atoms in core structure and 1 or 2 exit atoms:
* If walking from first atom (lowest relevant index) via root atom to second atom
* making a left turn, then an up-bond connecting first/only exist atom gives topicity=1
* - 3 neighbour atoms in core structure and one exit atom:
* If neighbours 1,2,3 are in counter-clockwise order and forth neighbour is connected
* with an up-bond, then topicity=1
* @param angle bond angles at stereo center sorted by relevant atom index
* @param coreNeighbourCount number of bonds at stereo center that have a counterpart in the core structure
* @param totalNeighbourCount number of bonds at stereo center including a second exit atom
* @param stereoType 1 (down) or 2 (up); corrected if original stereo bond is not the exit bond or/and second exit atom exists
* @return topicity 0 or 1
*/
private int calculateTHTopicity(double[] angle, int coreNeighbourCount, int totalNeighbourCount, int stereoType) {
for (int i=1; i Math.PI);
return leftTurn ^ (stereoType == 2) ? 1 : 0;
}
boolean clockwise = (angle[1] > angle[2]);
return clockwise ^ (stereoType == 2) ? 1 : 0;
}
private int calculateEZTopicity(StereoMolecule mol, int rootAtom, int exitAtom, int[] molToCoreAtom) {
if (mol.getAtomPi(rootAtom) != 1)
return -1;
if (mol.getBondOrder(mol.getBond(rootAtom, exitAtom)) != 1)
return -1;
int doubleBond = -1;
int rearDBAtom = -1;
for (int i=0; i candidate) {
oppositeAtom = candidate;
oppositeAngle = mol.getBondAngle(connAtom, rearDBAtom);
}
}
}
if (oppositeAtom == Integer.MAX_VALUE)
return -1;
double angleDif2 = Molecule.getAngleDif(oppositeAngle, dbAngle);
return (angleDif1 < 0) ^ (angleDif2 < 0) ? 0 : 1; // E:0 Z:1
}
/**
* If coreRootAtom matches a stereo center in the molecule and if coreConnAtom is one of coreRootAtom's
* neighbours in the core structure, then this method returns for this neighbour that atom index, which
* is used to determine the topicity for any exit vector (neighbours in mol, which are not part of the
* core structure). Typically, we use query structure atom indexes for this, i.e. the relevant atom index
* of a coreConnAtom is the index of its correscponding atom in the query structure. If, however,
* coreConnAtom is an atom of a matching bridge bond, then there is no corresponding query structure atom.
* In that case we walk along the bridge bond atoms in the core structure until we hit an atom that exists
* in the query, which is the remote bridge bond atom, whose index is then returned.
* @param coreRootAtom
* @param coreConnAtom
* @return
*/
private int getTopicityRelevantAtomIndex(int coreRootAtom, int coreConnAtom) {
int queryRoot = mCoreToQueryAtom[coreRootAtom];
// If the stereo center itself is not part of the query, then it is within a bridge bond
// and does not exist in all scaffolds of the scaffold group.
// In this case we use atom index of the core rather than the query as reference.
if (queryRoot == -1)
return coreConnAtom;
int queryAtom = mCoreToQueryAtom[coreConnAtom];
if (queryAtom != -1)
return queryAtom;
// If the stereo center neighbour in the core does not exist in the query, then it is part of
// a bridge bond. In this case we have to find that stereo center neighbour in the query
// that is connected with that bridge bond in the core, which copntains coreConnAtom.
// For that we build a graph from the core root atom adding only atoms that don't exist in
// the query until we hit a query core neighbour, which we return.
int[] bridgeNeighbour = new int[mQuery.getConnAtoms(queryRoot)];
int bridgeNeighbourCount = 0;
for (int i=0; i
* We define: If we have increasing atom indexes of query bonds in clockwise order,
* then topicity=0 is associated with an UP-bond and topicity=1 is associated with a DOWN-bond.
* @param exitVectorIndex
* @param coords receives suggested coordinates for first exit atom
* @return
*/
private int calculateExitVectorCoordsAndBondType(int exitVectorIndex, Coordinates coords) {
if ((mSeenBondOrdersOnScaffold[exitVectorIndex] & 2) == 0)
return ((mSeenBondOrdersOnScaffold[exitVectorIndex] & 4) == 0) ? 3 : 2;
ExitVector exitVector = getExitVector(exitVectorIndex);
int rootAtom = exitVector.getCoreAtom(mQueryToCoreAtom);
int[] neighbour = new int[3];
double[] angle = new double[3];
int coreNeighbourCount = 0;
int piBondSum = 0;
for (int i=0; i 0 && neighbour[index-1] > neighbourAtom) {
neighbour[index] = neighbour[index-1];
angle[index] = angle[index-1];
index--;
}
neighbour[index] = neighbourAtom;
angle[index] = mScaffold.getBondAngle(rootAtom, connAtom);
piBondSum += mScaffold.getConnBondOrder(rootAtom, i) - 1;
coreNeighbourCount++;
}
}
if (piBondSum != 0) {
if (coreNeighbourCount == 1)
calculateSP2ExitVectorCoords(rootAtom, piBondSum, exitVector.getTopicity(), angle, coords);
else
calculateSP3ExitVectorCoords(rootAtom, Arrays.copyOf(angle, coreNeighbourCount), coords);
// we assume that we don't have a stereo center with double bonds, e.g. at S or P
return Molecule.cBondTypeSingle;
}
calculateSP3ExitVectorCoords(rootAtom, Arrays.copyOf(angle, coreNeighbourCount), coords);
if (exitVector.getTopicity() == -1)
return Molecule.cBondTypeSingle;
angle[coreNeighbourCount] = Molecule.getAngle(mScaffold.getAtomX(rootAtom), mScaffold.getAtomY(rootAtom), coords.x, coords.y);
int totalNeighbourCount = coreNeighbourCount + 1;
if (coreNeighbourCount == 1) {
angle[coreNeighbourCount+1] = angle[coreNeighbourCount] + Math.PI * 2 / 3;
totalNeighbourCount++;
}
int topicity = calculateTHTopicity(angle, coreNeighbourCount, totalNeighbourCount, 1);
return (topicity == -1) ? Molecule.cBondTypeSingle : (topicity == exitVector.getTopicity()) ? Molecule.cBondTypeDown : Molecule.cBondTypeUp;
}
/**
* If there are at least two existing neighbours (angle.length >= 2), this method
* places the new neighbour at a position furthest away from any existing neighbour
* using the scaffolds average bond length. If there is only one neighbour, the new
* neighbour will be placed with a bond angle 120 degrees larger.
* @param atom
* @param angle
* @param coords
*/
private void calculateSP3ExitVectorCoords(int atom, double[] angle, Coordinates coords) {
double exitAngle = Math.PI * 2 / 3;
if (angle.length >= 2) {
Arrays.sort(angle);
double largestDiff = -1.0;
for (int i=0; i candidate) {
oppositeAtom = candidate;
oppositeAngle = mScaffold.getBondAngle(rearDBAtom, connAtom);
}
}
}
double angleDif = Molecule.getAngleDif(oppositeAngle, dbAngle);
exitAngle = angle[0] + ((angleDif < 0) ^ (topicity == 1) ? 0.6667 : 1.3333) * Math.PI;
}
double avbl = mScaffold.getAverageBondLength();
coords.x = mScaffold.getAtomX(atom) + avbl * Math.sin(exitAngle);
coords.y = mScaffold.getAtomY(atom) + avbl * Math.cos(exitAngle);
}
}