#pragma once /* * sketcherMinimizer.h * * Created by Nicola Zonta on 13/04/2010. * Copyright Schrodinger, LLC. All rights reserved. * */ #include #include #include #include "CoordgenConfig.hpp" #include "CoordgenFragmentBuilder.h" #include "CoordgenMinimizer.h" #include "sketcherMinimizerFragment.h" #include "sketcherMinimizerMarchingSquares.h" #include "sketcherMinimizerMolecule.h" #include "sketcherMinimizerResidue.h" #include "sketcherMinimizerResidueInteraction.h" #include "sketcherMinimizerRing.h" static const float SKETCHER_STANDARD_PRECISION = 1.f; static const float SKETCHER_QUICK_PRECISION = 0.2f; static const float SKETCHER_BEST_PRECISION = 3.f; class sketcherAtom; class sketcherBond; class sketcherMinimimizerInteraction; namespace schrodinger { namespace mae { class Block; } } // namespace schrodinger typedef struct { std::vector additionVectors; std::vector centers; std::vector counters; } proximityData; /* class to handle templates of common difficult ring structures */ class CoordgenTemplates { public: CoordgenTemplates() = default; ~CoordgenTemplates() { for (auto molecule : m_templates) { for (auto atom : molecule->_atoms) { delete atom; } for (auto bond : molecule->_bonds) { delete bond; } delete molecule; } m_templates.clear(); } std::vector& getTemplates() { return m_templates; } void setTemplateDir(std::string&& dir) { m_templateDir = std::move(dir); if (m_templateDir.back() != '/') { m_templateDir += "/"; } } std::string getTemplateDir() { return m_templateDir; } private: std::vector m_templates; std::string m_templateDir = ""; }; /* main class. Creates 2d coordinates for molecular inputs */ class EXPORT_COORDGEN sketcherMinimizer { public: sketcherMinimizer(float precision = SKETCHER_STANDARD_PRECISION); ~sketcherMinimizer(); CoordgenFragmentBuilder m_fragmentBuilder; CoordgenMinimizer m_minimizer; /* run coordinates generation and return true if the pose is considered * optimal */ bool runGenerateCoordinates(); /* return true if the molecules structure is reasonable (e.g. reasonable * amount of fused rings) */ bool structurePassSanityCheck() const; /* clear data and free memory */ void clear(); /* initialize data from given molecule */ void initialize(sketcherMinimizerMolecule* minMol); /* put atoms in a canonical order to reduce dependency from order in the * input vector */ static void canonicalOrdering(sketcherMinimizerMolecule* minMol); // void initializeFromMolecule(ChmMol& mol); /* if mol contains separate molecules, split them into a vector */ void splitIntoMolecules(sketcherMinimizerMolecule* mol, std::vector& mols); /* flag atoms that will be drawn with 90° angles (e.g. phosphate P) */ void flagCrossAtoms(); /* assign coordinates to all molecules and residues */ void minimizeAll(); /* assign coordinates to given molecule */ void minimizeMolecule(sketcherMinimizerMolecule* molecule); /* find the best angle to rotate each molecule */ void bestRotation(); /* add info to choose the best angle so that, if present, peptide chains are * horizontal */ void addBestRotationInfoForPeptides( std::vector>& angles, const std::vector& atoms); /* if a peptide chain is present make sure that the N term is on the left */ void maybeFlipPeptides(const std::vector& atoms, float& scoreX); /* consider flipping the molecule horizontally and/or vertically */ void maybeFlip(); /* mark atoms with wedges as above or below the plane to correctly draw * crossing bonds */ void assignPseudoZ(); /* write wedges and dashed bonds to mark stereochemistry */ void writeStereoChemistry(); /* arrange multiple molecules next to each other */ void arrangeMultipleMolecules(); /* arrange molecules that have parts that interact with each other so that * they are close */ void placeMoleculesWithProximityRelations( std::vector proximityMols); /* place molecules so that one is in the middle and other are around */ void placeMolResidueLigandStyle(sketcherMinimizerMolecule* mol, sketcherMinimizerMolecule* parent); /* if the molecule has more than one interaction and they cross mirror its coordinates so they don't cross anymore */ void flipIfCrossingInteractions(sketcherMinimizerMolecule* mol); /* build data vectors to place molecules with proximity relations */ std::vector buildProximityDataVector( std::vector& proximityMols, std::map& molMap); /* translate molecules with proximity relations */ void translateMoleculesWithProximityRelations( std::vector& proximityMols, std::map& molMap, std::map& templateCenters, std::vector& proximityDataVecto); /* rotate molecules with proximity relations */ void rotateMoleculesWithProximityRelations( std::vector& proximityMols, std::map& molMap, std::vector& proximityDataVector); /* place residues in concentric contours around the ligand */ void placeResiduesInCrowns(); /* place residues from the given strands of consecutive residues to fill the given path */ bool fillShape(std::vector>& SSEs, const std::vector& shape, int shapeN); /* place a single strand of consecutive residues */ void placeSSE(const std::vector& SSE, const std::vector& shape, int shapeN, std::vector& penalties, std::set& outliers, bool placeOnlyInteracting = false); /* score the position of the given strands */ float scoreSSEPosition(const std::vector& SSE, const std::vector& shape, int shapeN, std::vector& penalties, float f, float increment); /* score the distance between the two given points of connected residues */ float scoreSSEBondStretch(const sketcherMinimizerPointF& coordinates1, const sketcherMinimizerPointF& coordinates2); /* return the position of res, which is part of SSE, given that the first * residue of SSE is placed at startF and consecutive residues are placed * increment away from each other. All distances are expressed in floats, * where 0.f is an arbitrary starting point, 0.5 is the opposite side of the * curve and 1.0 is again the starting point */ float getResidueDistance(float startF, float increment, sketcherMinimizerResidue* res, const std::vector& SSE); /* return the vector index corresponding to floatPosition */ int getShapeIndex(const std::vector& shape, float floatPosition); /* solution represent the placement chosen for residues in SSE. Mark the * corresponding sections of the crown to prevent other residues to be * placed * there */ void markSolution(const std::pair& solution, const std::vector& SSE, const std::vector& shape, std::vector& penalties, std::set& outliers); /* return a concentric shape around the ligand. CrownN controls how far away from the ligand the shape is */ std::vector shapeAroundLigand(int crownN); /* group residues in strands of consecutive residues */ std::vector> groupResiduesInSSEs(const std::vector& residues); /* score the position of given residues */ float scoreResiduePosition(int index, const std::vector& shape, int shapeN, std::vector& penalties, sketcherMinimizerResidue* residue); /* assign coordinates to residues */ void placeResidues(const std::vector& atoms = std::vector(0)); /* assign coordinates to residues in the context of a protein-protein interaction diagram */ void placeResiduesProteinOnlyMode(); /* assign coordinates to residues in a protein-protein interaction diagram shaped as a circle */ void placeResiduesProteinOnlyModeCircleStyle( const std::map>& chains); /* assign coordinates to residues in a protein-protein interaction diagram shaped as a LID */ void placeResiduesProteinOnlyModeLIDStyle( const std::map>& chains); /* order residues for drawing, so that residues interacting together are * drawn one after the other and residues with more interactions are drawn * first */ std::vector orderResiduesOfChains( const std::map>& chains); /* find center position for each chain of residues using a meta-molecule * approach, building a molecule where each atom represents a chain and each * bond connects two interacting chains */ std::map computeChainsStartingPositionsMetaMol( const std::map>& chains); /* place interacting residues closer to each other, so they end up at the * periphery of the chain */ void shortenInteractions( const std::map>& chains); /* explore positions in a grid around the current one to ease clashes */ sketcherMinimizerPointF exploreGridAround( const sketcherMinimizerPointF& centerOfGrid, unsigned int levels, float gridD, float dx = 0.f, float dy = 0.f, float distanceFromAtoms = -1.f, bool watermap = false, sketcherMinimizerResidue* residueForInteractions = nullptr, const sketcherMinimizerPointF& direction = sketcherMinimizerPointF(0, 1)); sketcherMinimizerPointF exploreMolPosition(sketcherMinimizerMolecule* mol, unsigned int levels, float gridD, float distanceFromAtoms = -1.f); /* add given angle to a vector of angles, clustering them if a close angle * is already in the vector */ void addToVector(float weight, float angle, std::vector>& angles); /* return a score of the alignment between direction and templat.first, * weight on the angle between the two and templat.second */ static float testAlignment(const sketcherMinimizerPointF& direction, const std::pair& templat); /* find the best alignment of a fragment to its parent and set invert in * case the fragment needs to be flipped */ static sketcherMinimizerPointF scoreDirections( sketcherMinimizerFragment* fragment, float angle, const std::vector>& directions, bool& invert); /* align the fragment to its parent */ static void alignWithParentDirection(sketcherMinimizerFragment* f, const sketcherMinimizerPointF& position, float angle); /* align the fragment to its parent in the case of constrained coordinates */ static bool alignWithParentDirectionConstrained(sketcherMinimizerFragment* fragment, const sketcherMinimizerPointF& position, float angle); /* align the fragment to its parent in the case of unconstrained coordinates */ static bool alignWithParentDirectionUnconstrained(sketcherMinimizerFragment* fragment, float angle); /* get all bonds to a terminal atom */ static std::vector getAllTerminalBonds(sketcherMinimizerFragment* fragment); /* return a list of vectors the given fragment can be aligned with and a * score of the importance of each */ static std::vector> findDirectionsToAlignWith(sketcherMinimizerFragment* fragment); std::vector getAtoms() { return _atoms; } std::vector _atoms; std::vector _referenceAtoms; std::vector _residues; std::vector _residueInteractions; std::vector _fragments; std::vector _independentFragments; std::vector _bonds; std::vector _referenceBonds; std::vector m_proximityRelations; std::vector m_extraBonds; std::vector _molecules; void assignLongestChainFromHere(sketcherMinimizerFragment* f); void assignNumberOfChildrenAtomsFromHere(sketcherMinimizerFragment* f); // void exportCoordinates(ChmMol& molecule); /* split molecules into rigid fragments */ void findFragments(); /* initialize data and coordinates for each fragment */ void initializeFragments(); /* constrain coordinates on all atoms */ void constrainAllAtoms(); /* constrain coordinates on atoms corresponding to true */ void constrainAtoms(const std::vector& constrained); /* fix cooordinates (i.e. guarantee they will not change) on atoms marked as * true */ void fixAtoms(const std::vector& fixed); /* set a flag to enable/disable the scoring of interactions with residues */ void setScoreResidueInteractions(bool b); /* pick one atom out of the vector. Arbitrary criteria such as atomic number * and connectivity are used */ static sketcherMinimizerAtom* pickBestAtom(std::vector& atoms); /* if the three atoms share a ring, return it */ static sketcherMinimizerRing* sameRing(const sketcherMinimizerAtom* at1, const sketcherMinimizerAtom* at2, const sketcherMinimizerAtom* at3); /* if the two atoms share a ring, return it */ static sketcherMinimizerRing* sameRing(const sketcherMinimizerAtom* at1, const sketcherMinimizerAtom* at2); /* if the two atoms share a bond, return it */ static sketcherMinimizerBond* getBond(const sketcherMinimizerAtom* a1, const sketcherMinimizerAtom* a2); /* for each residue, find the closest atom among atoms, or all atoms if none are given */ void findClosestAtomToResidues(const std::vector& atoms = std::vector(0)); /* calculate root mean square deviation between templates and points */ static float RMSD(const std::vector& templates, const std::vector& points); /* singular value decomposition for 2x2 matrices. used for 2D alignment. */ static void svd(float* a, float* U, float* Sig, float* V); /* set m to a rotation matrix to align ref to points */ static void alignmentMatrix(const std::vector& ref, const std::vector& points, float* m); static void checkIdentity(std::vector solution, int newSol, std::vector& matrix, std::vector& templateCoordinates, std::vector>& molBonds, std::vector>& templateBonds, std::vector>& molCisTransChains, std::vector& molIsCis, size_t size, bool& found, std::vector& mapping); /* compare atoms and bonds to template and map which atom is which in case * of a positive match */ static bool compare(const std::vector& atoms, const std::vector& bonds, sketcherMinimizerMolecule* templ, std::vector& mapping); /* calculate morgan scores for the given input */ static int morganScores(const std::vector& atoms, const std::vector& bonds, std::vector& scores); std::string m_chainHint; /* load the templates from the template file */ static void setTemplateFileDir(std::string dir); static void loadTemplates(); static CoordgenTemplates m_templates; bool getTreatNonterminalBondsToMetalAsZOBs() {return m_treatNonterminalBondsToMetalAsZOBs;} void setTreatNonterminalBondsToMetalAsZOBs(bool b) {m_treatNonterminalBondsToMetalAsZOBs = b;} private: /*all non-terminal bonds to a metal atom are treated as if they were zero order bonds (this usually results in a longer bond*/ bool m_treatNonterminalBondsToMetalAsZOBs = true; };