Skip to content

Commit

Permalink
Merge pull request #1777 from ghutchis/place-template
Browse files Browse the repository at this point in the history 
Update the template tool to place ligands or functional groups
  • Loading branch information
@ghutchis
ghutchis authored Nov 12, 2024
2 parents 4941a21 + fa12cf2 commit 5d8cf11
Show file tree
Hide file tree
Showing 3 changed files with 188 additions and 87 deletions.
268 changes: 181 additions & 87 deletions avogadro/qtplugins/templatetool/templatetool.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -269,36 +269,113 @@ void TemplateTool::reset()

void TemplateTool::emptyLeftClick(QMouseEvent* e)
{
QFile templ(":/templates/centers/" + m_toolWidget->coordinationString() +
".cjson");
if (!templ.open(QFile::ReadOnly | QFile::Text))
// Get the coordinates of the clicked position
if (m_renderer == nullptr)
return;
QTextStream templateStream(&templ);

m_toolWidget->selectedUIDs().clear();
Vector2f windowPos(e->localPos().x(), e->localPos().y());
Vector3f atomPos = m_renderer->camera().unProject(windowPos);
// center of inserted template
Vector3 center(0.0f, 0.0f, 0.0f);

CjsonFormat ff;
Molecule templateMolecule;
if (!ff.readString(templateStream.readAll().toStdString(), templateMolecule))
return;

m_toolWidget->selectedUIDs().clear();
// before we do anything, check if it's a metal or a ligand
// in the dialog
int currentTab = m_toolWidget->currentTab();
if (currentTab == 0) { // metal center
QFile templ(":/templates/centers/" + m_toolWidget->coordinationString() +
".cjson");
if (!templ.open(QFile::ReadOnly | QFile::Text))
return;
QTextStream templateStream(&templ);

// Add an atom at the clicked position
Vector2f windowPos(e->localPos().x(), e->localPos().y());
Vector3f atomPos = m_renderer->camera().unProject(windowPos);
if (!ff.readString(templateStream.readAll().toStdString(),
templateMolecule))
return;

// Add hydrogens around it following template
Vector3 center(0.0f, 0.0f, 0.0f);
size_t centerIndex = 0;
for (size_t i = 0; i < templateMolecule.atomCount(); i++) {
if (templateMolecule.atomicNumber(i) != 1) {
center = templateMolecule.atomPosition3d(i);
centerIndex = i;
templateMolecule.setAtomicNumber(i, m_toolWidget->atomicNumber());
templateMolecule.setFormalCharge(i, m_toolWidget->formalCharge());
continue;
// Add the atom and hydrogens around it following template
size_t centerIndex = 0;
for (size_t i = 0; i < templateMolecule.atomCount(); i++) {
if (templateMolecule.atomicNumber(i) != 1) {
center = templateMolecule.atomPosition3d(i);
centerIndex = i;
templateMolecule.setAtomicNumber(i, m_toolWidget->atomicNumber());
templateMolecule.setFormalCharge(i, m_toolWidget->formalCharge());
continue;
}
}
// done with metal center and coordination
} else {
// ligand
// check if it's clipboard first
if (m_toolWidget->ligandString() == tr("Clipboard")) {
const QMimeData* mimeData(QApplication::clipboard()->mimeData());

if (!mimeData) {
return;
}

// Try to find a reader that can handle the available mime-types.
Io::FileFormatManager& mgr = Io::FileFormatManager::instance();
QStringList mimeTypes(mimeData->formats());
Io::FileFormat* pastedFormat = nullptr;
QByteArray pastedData;
Io::FileFormat::Operations ops(Io::FileFormat::Read |
Io::FileFormat::String);
foreach (const QString& mimeType, mimeTypes) {
if ((pastedFormat =
mgr.newFormatFromMimeType(mimeType.toStdString(), ops))) {
pastedData = mimeData->data(mimeType);
break;
}
}

// No mime-type match, default to cjson.
if (!pastedFormat && mimeData->hasText()) {
pastedFormat = new Io::CjsonFormat;
pastedData = mimeData->text().toLatin1();
}

if (pastedFormat == nullptr)
return;

// we have a format, so try to insert the new bits into the molecule
bool success = pastedFormat->readString(
std::string(pastedData.constData(), pastedData.size()),
templateMolecule);

if (!success)
return;

center = templateMolecule.centerOfGeometry();
// done with clipboard ligands
} else { // a ligand file
QString path;
if (m_toolWidget->ligandString().endsWith(".cjson")) {
// we already have the full path .. from the insert browser
path = m_toolWidget->ligandString();
} else {
path = ":/templates/ligands/" + m_toolWidget->ligandString() + ".cjson";
}

QFile templ(path);
if (!templ.open(QFile::ReadOnly | QFile::Text))
return;
QTextStream templateStream(&templ);

if (!ff.readString(templateStream.readAll().toStdString(),
templateMolecule))
return;

center = templateMolecule.centerOfGeometry();
// done with ligand
}
}

// move the template to the clicked position
for (size_t i = 0; i < templateMolecule.atomCount(); i++) {
Vector3 pos =
templateMolecule.atomPosition3d(i) - center + atomPos.cast<double>();
Expand All @@ -317,7 +394,7 @@ void TemplateTool::emptyLeftClick(QMouseEvent* e)
// Update the clicked object
m_clickedObject.type = Rendering::AtomType;
m_clickedObject.molecule = m_molecule;
m_clickedObject.index = firstIndex + centerIndex;
m_clickedObject.index = firstIndex;

// Emit changed signal
m_molecule->emitChanged(changes);
Expand Down Expand Up @@ -371,8 +448,10 @@ Matrix3 applyKabsch(std::vector<Vector3> templatePoints,
void TemplateTool::atomLeftClick(QMouseEvent*)
{
size_t selectedIndex = m_clickedObject.index;
// if it's a valid selected atom and a hydrogen or dummy atom
if (m_molecule->atom(selectedIndex).isValid() &&
m_molecule->atomicNumber(selectedIndex) == 1) {
(m_molecule->atomicNumber(selectedIndex) == 1 ||
m_molecule->atomicNumber(selectedIndex) == 0)) {
m_toolWidget->selectedUIDs().push_back(
m_molecule->atomUniqueId(selectedIndex));
if (static_cast<int>(m_toolWidget->selectedUIDs().size()) !=
Expand All @@ -384,7 +463,7 @@ void TemplateTool::atomLeftClick(QMouseEvent*)
// - otherwise use the template
Molecule templateMolecule;

if (m_toolWidget->ligandString() == "Clipboard") {
if (m_toolWidget->ligandString() == tr("Clipboard")) {
const QMimeData* mimeData(QApplication::clipboard()->mimeData());

if (!mimeData) {
Expand Down Expand Up @@ -461,36 +540,46 @@ void TemplateTool::atomLeftClick(QMouseEvent*)
}

// Find center atom in molecule and get all necessary info
size_t moleculeCenterIndex =
m_molecule->bonds(selectedIndex)[0].getOtherAtom(selectedIndex).index();
size_t moleculeCenterUID = m_molecule->atomUniqueId(moleculeCenterIndex);
// - first check to see if there is a bond
Vector3 moleculeLigandOutVector(0.0, 0.0, 0.0);
for (size_t UID : m_toolWidget->selectedUIDs()) {
size_t index = m_molecule->atomByUniqueId(UID).index();
Vector3 newPos = m_molecule->atomPosition3d(index);
moleculeLigandOutVector +=
newPos - m_molecule->atomPosition3d(moleculeCenterIndex);
}

// Estimate and try to realize bond distances
Vector3 displacement(0.0, 0.0, 0.0);
for (size_t i = 0; i < templateLigandIndices.size(); i++) {
unsigned char ligandAtomicNumber =
templateMolecule.atomicNumber(templateLigandIndices[i]);
ligandAtomicNumber = (ligandAtomicNumber == 0) ? 6 : ligandAtomicNumber;
// Estimate as the sum of covalent radii
double bondDistance =
Elements::radiusCovalent(ligandAtomicNumber) +
Elements::radiusCovalent(m_molecule->atomicNumber(moleculeCenterIndex));
Vector3 inVector =
templateMolecule.atomPosition3d(templateDummyIndex) -
templateMolecule.atomPosition3d(templateLigandIndices[i]);
Vector3 correctionVector = inVector;
correctionVector.normalize();
correctionVector *= bondDistance - inVector.norm();
displacement += correctionVector;
Vector3 centerPosition = m_molecule->atomPosition3d(selectedIndex);
size_t moleculeCenterIndex = selectedIndex;
size_t moleculeCenterUID = m_molecule->atomUniqueId(moleculeCenterIndex);

if (m_molecule->bonds(selectedIndex).size() != 0) {
moleculeCenterIndex =
m_molecule->bonds(selectedIndex)[0].getOtherAtom(selectedIndex).index();
moleculeCenterUID = m_molecule->atomUniqueId(moleculeCenterIndex);
for (size_t UID : m_toolWidget->selectedUIDs()) {
size_t index = m_molecule->atomByUniqueId(UID).index();
Vector3 newPos = m_molecule->atomPosition3d(index);
moleculeLigandOutVector +=
newPos - m_molecule->atomPosition3d(moleculeCenterIndex);
}

// Estimate and try to realize bond distances
for (size_t i = 0; i < templateLigandIndices.size(); i++) {
unsigned char ligandAtomicNumber =
templateMolecule.atomicNumber(templateLigandIndices[i]);
ligandAtomicNumber = (ligandAtomicNumber == 0) ? 6 : ligandAtomicNumber;
// Estimate as the sum of covalent radii
double bondDistance = Elements::radiusCovalent(ligandAtomicNumber) +
Elements::radiusCovalent(
m_molecule->atomicNumber(moleculeCenterIndex));
Vector3 inVector =
templateMolecule.atomPosition3d(templateDummyIndex) -
templateMolecule.atomPosition3d(templateLigandIndices[i]);
Vector3 correctionVector = inVector;
correctionVector.normalize();
correctionVector *= bondDistance - inVector.norm();
displacement += correctionVector;
}
displacement *= 1.0 / templateLigandIndices.size();
} else {
// direction can be random
displacement = Eigen::Vector3d::Random();
}
displacement *= 1.0 / templateLigandIndices.size();
Vector3 newPos =
templateMolecule.atomPosition3d(templateDummyIndex) + displacement;
templateMolecule.setAtomPosition3d(templateDummyIndex, newPos);
Expand All @@ -505,43 +594,45 @@ void TemplateTool::atomLeftClick(QMouseEvent*)
}
}

// Create arrays with the points to align and apply Kabsch algorithm
std::vector<Vector3> templateLigandPositions;
for (size_t index : templateLigandIndices)
templateLigandPositions.push_back(
templateMolecule.atomPosition3d(index) -
m_molecule->atomPosition3d(moleculeCenterIndex));
std::vector<Vector3> moleculeLigandPositions;
for (size_t UID : m_toolWidget->selectedUIDs())
moleculeLigandPositions.push_back(
m_molecule->atomPosition3d(m_molecule->atomByUniqueId(UID).index()) -
m_molecule->atomPosition3d(moleculeCenterIndex));
Matrix3 rotation =
applyKabsch(templateLigandPositions, moleculeLigandPositions);
for (size_t i = 0; i < templateMolecule.atomCount(); i++) {
if (i != templateDummyIndex) {
templateMolecule.setAtomPosition3d(
i, rotation * (templateMolecule.atomPosition3d(i) -
m_molecule->atomPosition3d(moleculeCenterIndex)) +
m_molecule->atomPosition3d(moleculeCenterIndex));
if (m_molecule->bonds(selectedIndex).size() != 0) {
// Create arrays with the points to align and apply Kabsch algorithm
std::vector<Vector3> templateLigandPositions;
for (size_t index : templateLigandIndices)
templateLigandPositions.push_back(
templateMolecule.atomPosition3d(index) -
m_molecule->atomPosition3d(moleculeCenterIndex));
std::vector<Vector3> moleculeLigandPositions;
for (size_t UID : m_toolWidget->selectedUIDs())
moleculeLigandPositions.push_back(
m_molecule->atomPosition3d(m_molecule->atomByUniqueId(UID).index()) -
m_molecule->atomPosition3d(moleculeCenterIndex));
Matrix3 rotation =
applyKabsch(templateLigandPositions, moleculeLigandPositions);
for (size_t i = 0; i < templateMolecule.atomCount(); i++) {
if (i != templateDummyIndex) {
templateMolecule.setAtomPosition3d(
i, rotation * (templateMolecule.atomPosition3d(i) -
m_molecule->atomPosition3d(moleculeCenterIndex)) +
m_molecule->atomPosition3d(moleculeCenterIndex));
}
}
}

// Rotate partially aligned template to align "out" vectors
Vector3 templateLigandOutVector(0.0, 0.0, 0.0);
for (size_t index : templateLigandIndices) {
Vector3 pos = templateMolecule.atomPosition3d(index);
templateLigandOutVector +=
pos - m_molecule->atomPosition3d(moleculeCenterIndex);
}
for (size_t i = 0; i < templateMolecule.atomCount(); i++) {
if (templateMolecule.atomicNumber(i) != 0) {
templateMolecule.setAtomPosition3d(
i,
rotateLigandCoords(templateMolecule.atomPosition3d(i) -
m_molecule->atomPosition3d(moleculeCenterIndex),
templateLigandOutVector, moleculeLigandOutVector) +
m_molecule->atomPosition3d(moleculeCenterIndex));
// Rotate partially aligned template to align "out" vectors
Vector3 templateLigandOutVector(0.0, 0.0, 0.0);
for (size_t index : templateLigandIndices) {
Vector3 pos = templateMolecule.atomPosition3d(index);
templateLigandOutVector +=
pos - m_molecule->atomPosition3d(moleculeCenterIndex);
}
for (size_t i = 0; i < templateMolecule.atomCount(); i++) {
if (templateMolecule.atomicNumber(i) != 0) {
templateMolecule.setAtomPosition3d(
i, rotateLigandCoords(
templateMolecule.atomPosition3d(i) -
m_molecule->atomPosition3d(moleculeCenterIndex),
templateLigandOutVector, moleculeLigandOutVector) +
m_molecule->atomPosition3d(moleculeCenterIndex));
}
}
}

Expand All @@ -561,8 +652,11 @@ void TemplateTool::atomLeftClick(QMouseEvent*)
}

// Remove selected atoms and insert ligand
for (size_t UID : m_toolWidget->selectedUIDs())
m_molecule->removeAtom(m_molecule->atomByUniqueId(UID).index());
// (unless there wasn't a bond to begin with)
if (m_molecule->bonds(selectedIndex).size() != 0) {
for (size_t UID : m_toolWidget->selectedUIDs())
m_molecule->removeAtom(m_molecule->atomByUniqueId(UID).index());
}
size_t moleculeBaseIndex = m_molecule->atomCount();
m_molecule->appendMolecule(templateMolecule, tr("Insert Ligand"));

Expand Down
5 changes: 5 additions & 0 deletions avogadro/qtplugins/templatetool/templatetoolwidget.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -156,6 +156,11 @@ QString TemplateToolWidget::coordinationString() const
return m_centers.at(m_ui->coordinationComboBox->currentIndex());
}

int TemplateToolWidget::currentTab() const
{
return m_ui->tabWidget->currentIndex();
}

unsigned char TemplateToolWidget::ligand() const
{
return static_cast<unsigned char>(m_ui->ligandComboBox->currentIndex());
Expand Down
2 changes: 2 additions & 0 deletions avogadro/qtplugins/templatetool/templatetoolwidget.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -43,6 +43,8 @@ class TemplateToolWidget : public QWidget
int denticity() const;
std::vector<size_t>& selectedUIDs();

int currentTab() const;

private slots:
void elementChanged(int index);
void updateElementCombo();
Expand Down

0 comments on commit 5d8cf11

Please sign in to comment.