All Downloads are FREE. Search and download functionalities are using the official Maven repository.

com.actelion.research.chem.sar.ScaffoldData Maven / Gradle / Ivy

There is a newer version: 2024.12.1
Show newest version
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 ScaffoldData {
	private StereoMolecule mQuery,mCore,mScaffold;
	private int[] mCoreToQueryAtom,mQueryToCoreAtom;
	private String mIDCodeWithRGroups, mIDCoordsWithRGroups;
	private TreeMap mOldToNewMap;
	private int mBridgeAtomRGroupCount;
	private ScaffoldGroup mScaffoldGroup;
	private ExitVector[] mBridgeAtomExitVector;
	private boolean[] mHasSeenSubstituentOnScaffold;
	private int[] mSeenBondOrdersOnScaffold;

	protected ScaffoldData(StereoMolecule query, StereoMolecule core, int[] coreToQueryAtom, int[] queryToCoreAtom, ScaffoldGroup 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 && neighbour[index-1] > neighbourAtom) {
					neighbour[index] = neighbour[index-1];
					angle[index] = angle[index-1];
					index--;
				}

				if (stereoBondIndex >= index)
					stereoBondIndex++;

				neighbour[index] = neighbourAtom;
				angle[index] = mol.getBondAngle(rootAtom, connAtom);

				if (connBond == stereoBond)
					stereoBondIndex = index;

				count++;
			}
		}

		if (count < 3 || count > 4)
			return -1;

		int stereoType = (mol.getBondType(stereoBond) == Molecule.cBondTypeUp) ? 2 : 1;

		if (stereoBond == otherExitBond) {
			stereoBondIndex = count-1;   // exit vector bond
			stereoType = 3 - stereoType;
		}

		return calculateTHTopicity(angle, count, stereoBondIndex, stereoType);
	}

	/**
	 * Use Canonizer's logic to calculate parities of 2D stereo center.
	 * The difference is that we assign the highest atom index to the exit atom,
	 * and neglect an optional second exit vector. If a second exit vector exists and
	 * if that is connected with a stereo bond, then we adapt accordingly.
	 * @param angle
	 * @return
	 */
	private int calculateTHTopicity(double[] angle, int count, int stereoBondIndex, int stereoType) {
		final int[][] up_down = { { 2,1,2,1 },	// direction of stereobond
				{ 1,2,2,1 },	// for topicity = 1
				{ 1,1,2,2 },	// first dimension: order of
				{ 2,1,1,2 },	// angles to connected atoms
				{ 2,2,1,1 },	// second dimension: number of
				{ 1,2,1,2 } };  // mMol.getConnAtom that has stereobond

// No support for Fisher projections here!
//		byte parity = (byte)mMol.getFisherProjectionParity(atom, remappedConn, angle, null);
//		if (parity != Molecule.cAtomParityUnknown)
//			return parity;

		for (int i=1; i angle[2]) && (angle[1] - angle[2] > Math.PI)))
						stereoType = 3 - stereoType;
					break;
				case 1:
					if (angle[2] - angle[0] > Math.PI)
						stereoType = 3 - stereoType;
					break;
				case 2:
					if (angle[1] - angle[0] < Math.PI)
						stereoType = 3 - stereoType;
					break;
			}

			return (stereoType == 1) ? 0 : 1;
		}

		int order = 0;
		if		(angle[1] <= angle[2] && angle[2] <= angle[3]) order = 0;
		else if (angle[1] <= angle[3] && angle[3] <= angle[2]) order = 1;
		else if (angle[2] <= angle[1] && angle[1] <= angle[3]) order = 2;
		else if (angle[2] <= angle[3] && angle[3] <= angle[1]) order = 3;
		else if (angle[3] <= angle[1] && angle[1] <= angle[2]) order = 4;
		else if (angle[3] <= angle[2] && angle[2] <= angle[1]) order = 5;
		return (up_down[order][stereoBondIndex] == stereoType) ? 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;
		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 macthing 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 count = 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;
				count++;
			}
		}

		if (piBondSum != 0) {
			if (count == 1)
				calculateEZExitVectorCoords(rootAtom, piBondSum, exitVector.getTopicity(), angle, coords);
			else
				calculateFurthestAwayExitVectorCoords(rootAtom, Arrays.copyOf(angle, count), coords);

			// we assume that we don't have a stereo center with double bonds, e.g. at S or P
			return Molecule.cBondTypeSingle;
		}

		calculateFurthestAwayExitVectorCoords(rootAtom, Arrays.copyOf(angle, count), coords);

		if (exitVector.getTopicity() == -1)
			return Molecule.cBondTypeSingle;

		angle[count] = Molecule.getAngle(mScaffold.getAtomX(rootAtom), mScaffold.getAtomY(rootAtom), coords.x, coords.y);

		int topicity = calculateTHTopicity(angle, count+1, count, 1);
		return (topicity == -1) ? Molecule.cBondTypeSingle : (topicity == exitVector.getTopicity()) ? Molecule.cBondTypeDown : Molecule.cBondTypeUp;
	}

	/**
	 * Assuming at least two existing neighbours, and using the scaffolds average bond length, this method
	 * places the new neighbour at a position furthest away from any existing neighbour.
	 * @param atom
	 * @param angle
	 * @param coords
	 */
	private void calculateFurthestAwayExitVectorCoords(int atom, double[] angle, Coordinates coords) {
		double exitAngle = 0.0;

		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);
	}
}




© 2015 - 2025 Weber Informatics LLC | Privacy Policy