Merge pull request #1837 from ghutchis/add-conformer-properties

Add an initial conformer property window

avogadro/io/xyzformat.cpp

@@ -260,6 +260,10 @@ bool XyzFormat::read(std::istream& inStream, Core::Molecule& mol)
       }
 
       std::getline(inStream, buffer); // Skip the blank
+      // check for energies
+      if (findEnergy(buffer, energy)) {
+        energies.push_back(energy);
+      }
       positions.clear();
     }
   }
@@ -280,9 +284,8 @@ bool XyzFormat::read(std::istream& inStream, Core::Molecule& mol)
     mol.setPartialCharges("From File", chargesMatrix);
   }
 
-  if (energies.size() > 1) {
+  if (energies.size() > 1)
     mol.setData("energies", energies);
-  }
 
   return true;
 }

avogadro/qtplugins/openbabel/openbabel.cpp

@@ -706,9 +706,18 @@ void OpenBabel::onGenerateConformersFinished(const QByteArray& output)
     return;
   }
 
-  //@todo .. multiple coordinate sets
   m_molecule->undoMolecule()->setAtomPositions3d(mol.atomPositions3d(),
                                                  tr("Generate Conformers"));
+
+  // copy the coordinate sets
+  m_molecule->clearCoordinate3d();
+  for (size_t i = 0; i < mol.coordinate3dCount(); ++i)
+    m_molecule->setCoordinate3d(mol.coordinate3d(i), i);
+
+  // energy data too
+  // TODO: check if other properties are needed
+  m_molecule->setData("energies", mol.data("energies"));
+
   m_molecule->emitChanged(QtGui::Molecule::Atoms | QtGui::Molecule::Modified);
   m_progress->reset();
 }

avogadro/qtplugins/propertytables/propertymodel.cpp

@@ -6,6 +6,7 @@
 #include "propertymodel.h"
 
 #include <avogadro/calc/chargemanager.h>
+#include <avogadro/core/array.h>
 #include <avogadro/core/atom.h>
 #include <avogadro/core/bond.h>
 #include <avogadro/core/elements.h>
@@ -25,6 +26,7 @@
 
 namespace Avogadro {
 
+using Avogadro::Core::Array;
 using Avogadro::QtGui::Molecule;
 using QtGui::Molecule;
 using QtGui::RWAtom;
@@ -45,8 +47,23 @@ const int AngleColumns = 5;
 const int TorsionColumns = 6;
 // name, number, chain, secondary structure, heterogen, color
 const int ResidueColumns = 6;
-// number, energy, ??
-const int ConformerColumns = 2;
+// number, rmsd, energy or more depending on available properties
+const int ConformerColumns = 1;
+
+// compute the RMSD between the two sets of coordinates
+inline double calcRMSD(const Array<Vector3>& v1, const Array<Vector3>& v2)
+{
+  // if they're not the same length, it's an error
+  if (v1.size() != v2.size())
+    return numeric_limits<double>::quiet_NaN();
+
+  double sum = 0.0;
+  for (size_t i = 0; i < v1.size(); ++i) {
+    Vector3 diff = v1[i] - v2[i];
+    sum += diff.squaredNorm();
+  }
+  return sqrt(sum / v1.size());
+}
 
 inline double distance(Vector3 v1, Vector3 v2)
 {
@@ -119,8 +136,12 @@ int PropertyModel::columnCount(const QModelIndex& parent) const
       return TorsionColumns;
     case ResidueType:
       return ResidueColumns;
-    case ConformerType:
-      return ConformerColumns;
+    case ConformerType: {
+      if (m_molecule->hasData("energies"))
+        return ConformerColumns + 1;
+      else
+        return ConformerColumns;
+    }
     default:
       return 0;
   }
@@ -204,6 +225,9 @@ QVariant PropertyModel::data(const QModelIndex& index, int role) const
         return toVariant(Qt::AlignHCenter | Qt::AlignVCenter);
     } else if (m_type == ResidueType) {
       return toVariant(Qt::AlignHCenter | Qt::AlignVCenter);
+    } else if (m_type == ConformerType) {
+      return toVariant(Qt::AlignRight |
+                       Qt::AlignVCenter); // RMSD or other properties
     }
   }
 
@@ -227,9 +251,6 @@ QVariant PropertyModel::data(const QModelIndex& index, int role) const
   if (role != Qt::UserRole && role != Qt::DisplayRole && role != Qt::EditRole)
     return QVariant();
 
-  //  if (!m_validCache)
-  //    updateCache();
-
   if (m_type == AtomType) {
     auto column = static_cast<AtomColumn>(index.column());
 
@@ -288,7 +309,7 @@ QVariant PropertyModel::data(const QModelIndex& index, int role) const
       case BondDataOrder:
         return bond.order();
       default: // length, rounded to 4 decimals
-        return QString("%L1").arg(
+        return QString("%L1 Å").arg(
           distance(atom1.position3d(), atom2.position3d()), 0, 'f', 3);
     }
   } else if (m_type == ResidueType) {
@@ -342,7 +363,7 @@ QVariant PropertyModel::data(const QModelIndex& index, int role) const
       case AngleDataAtom3:
         return QVariant::fromValue(std::get<2>(angle) + 1);
       case AngleDataValue:
-        return QString("%L1").arg(calcAngle(a1, a2, a3), 0, 'f', 3);
+        return QString("%L1 °").arg(calcAngle(a1, a2, a3), 0, 'f', 3);
       default:
         return QVariant();
     }
@@ -378,10 +399,42 @@ QVariant PropertyModel::data(const QModelIndex& index, int role) const
       case TorsionDataAtom4:
         return QVariant::fromValue(std::get<3>(torsion) + 1);
       case TorsionDataValue:
-        return QString("%L1").arg(calcDihedral(a1, a2, a3, a4), 0, 'f', 3);
+        return QString("%L1 °").arg(calcDihedral(a1, a2, a3, a4), 0, 'f', 3);
       default:
         return QVariant();
     }
+  } else if (m_type == ConformerType) {
+    auto column = static_cast<ConformerColumn>(index.column());
+    if (row >= static_cast<int>(m_molecule->coordinate3dCount()) ||
+        column > ConformerColumns) {
+      return QVariant(); // invalid index
+    }
+
+    switch (column) {
+      case ConformerDataRMSD: { // rmsd
+        double rmsd = 0.0;
+        if (row > 0) {
+          rmsd = calcRMSD(m_molecule->coordinate3d(row),
+                          m_molecule->coordinate3d(0));
+        }
+
+        return QString("%L1 Å").arg(rmsd, 0, 'f', 3);
+      }
+      case ConformerDataEnergy: {
+        double energy = 0.0;
+        if (m_molecule->hasData("energies")) {
+          std::vector<double> energies = m_molecule->data("energies").toList();
+          // calculate the minimum
+          double minEnergy = std::numeric_limits<double>::max();
+          for (double e : energies) {
+            minEnergy = std::min(minEnergy, e);
+          }
+          if (row < static_cast<int>(energies.size()))
+            energy = energies[row] - minEnergy;
+        }
+        return QString("%L1").arg(energy, 0, 'f', 4);
+      }
+    }
   }
 
   return QVariant();
@@ -500,6 +553,20 @@ QVariant PropertyModel::headerData(int section, Qt::Orientation orientation,
       }
     } else // row headers
       return QString("%L1").arg(section + 1);
+  } else if (m_type == ConformerType) {
+    // check if we have energies
+    bool hasEnergies = (m_molecule->hasData("energies"));
+    if (orientation == Qt::Horizontal) {
+      unsigned int column = static_cast<ConformerColumn>(section);
+      switch (column) {
+        case ConformerDataRMSD:
+          return tr("RMSD (Å)", "root mean squared displacement in Angstrom");
+        case ConformerDataEnergy:
+          // should only hit this if we have energies anyway
+          return hasEnergies ? tr("Energy (kcal/mol)") : tr("Property");
+      }
+    } else // row headers
+      return QString("%L1").arg(section + 1);
   }
 
   return QVariant();

avogadro/qtplugins/propertytables/propertymodel.h

@@ -81,8 +81,7 @@ public slots:
   std::vector<int> m_fragment;
   Eigen::Affine3d m_transform;
   bool fragmentHasAtom(int uid) const;
-  void buildFragment(const QtGui::RWBond& bond,
-                      const QtGui::RWAtom& startAtom);
+  void buildFragment(const QtGui::RWBond& bond, const QtGui::RWAtom& startAtom);
   bool fragmentRecurse(const QtGui::RWBond& bond,
                        const QtGui::RWAtom& startAtom,
                        const QtGui::RWAtom& currentAtom);
@@ -152,7 +151,7 @@ public slots:
   // Conformer Data
   enum ConformerColumn
   {
-    ConformerDataType = 0,
+    ConformerDataRMSD = 0,
     ConformerDataEnergy
   };
 

avogadro/qtplugins/propertytables/propertytables.cpp

@@ -58,6 +58,14 @@ PropertyTables::PropertyTables(QObject* parent_)
   action->setEnabled(false); // changed by molecule
   connect(action, SIGNAL(triggered()), SLOT(showDialog()));
   m_actions.append(action);
+
+  action = new QAction(this);
+  action->setText(tr("Conformer Properties…"));
+  action->setData(PropertyType::ConformerType);
+  action->setProperty("menu priority", 830);
+  action->setEnabled(false); // changed by molecule
+  connect(action, SIGNAL(triggered()), SLOT(showDialog()));
+  m_actions.append(action);
 }
 
 PropertyTables::~PropertyTables() {}
@@ -85,12 +93,26 @@ void PropertyTables::setMolecule(QtGui::Molecule* mol)
 
   m_molecule = mol;
 
-  // check if there are residues
-  if (m_molecule->residueCount() > 0) {
-    for (const auto& action : m_actions) {
-      if (action->data().toInt() == PropertyType::ResidueType)
-        action->setEnabled(true);
-    }
+  updateActions();
+
+  // update if the molecule changes
+  connect(m_molecule, SIGNAL(changed(unsigned int)), SLOT(updateActions()));
+}
+
+void PropertyTables::updateActions()
+{
+  if (m_molecule == nullptr)
+    return;
+
+  // check if we enable / disable the residue and conformer actions
+  bool haveResidues = (m_molecule->residueCount() > 0);
+  // technically coordinate sets
+  bool haveConformers = (m_molecule->coordinate3dCount() > 1);
+  for (const auto& action : m_actions) {
+    if (action->data().toInt() == PropertyType::ResidueType)
+      action->setEnabled(haveResidues);
+    else if (action->data().toInt() == PropertyType::ConformerType)
+      action->setEnabled(haveConformers);
   }
 }
 
@@ -104,6 +126,10 @@ void PropertyTables::showDialog()
       m_molecule->residueCount() == 0)
     return;
 
+  if (action->data().toInt() == PropertyType::ConformerType &&
+      m_molecule->coordinate3dCount() < 2)
+    return;
+
   auto* dialog = new QDialog(qobject_cast<QWidget*>(parent()));
   auto* layout = new QVBoxLayout(dialog);
   dialog->setLayout(layout);
@@ -130,6 +156,7 @@ void PropertyTables::showDialog()
   view->setSourceModel(model);
 
   view->setHorizontalScrollBarPolicy(Qt::ScrollBarAlwaysOff);
+  view->horizontalHeader()->setStretchLastSection(true);
   view->resizeColumnsToContents();
   layout->addWidget(view);
   dialog->setWindowTitle(view->windowTitle());
@@ -149,4 +176,4 @@ void PropertyTables::showDialog()
   dialog->show();
 }
 
-} // namespace Avogadro
+} // namespace Avogadro::QtPlugins

avogadro/qtplugins/propertytables/propertytables.h

@@ -33,12 +33,13 @@ class PropertyTables : public Avogadro::QtGui::ExtensionPlugin
 
 public slots:
   void setMolecule(QtGui::Molecule* mol) override;
+  void updateActions();
 
 private slots:
   void showDialog();
 
 private:
-  QList<QAction *> m_actions;
+  QList<QAction*> m_actions;
   QtGui::Molecule* m_molecule;
 };
 

avogadro/qtplugins/propertytables/propertyview.cpp

@@ -5,6 +5,7 @@
 
 #include "propertyview.h"
 
+#include <avogadro/core/residue.h>
 #include <avogadro/qtgui/molecule.h>
 
 #include <QAction>
@@ -74,7 +75,10 @@ PropertyView::PropertyView(PropertyType type, QWidget* parent)
   // You can select everything (e.g., to copy, select all, etc.)
   setCornerButtonEnabled(true);
   setSelectionBehavior(QAbstractItemView::SelectRows);
-  setSelectionMode(QAbstractItemView::ExtendedSelection);
+  if (type == ConformerType)
+    setSelectionMode(QAbstractItemView::SingleSelection);
+  else
+    setSelectionMode(QAbstractItemView::ExtendedSelection);
   // Alternating row colors
   setAlternatingRowColors(true);
   // Allow sorting the table
@@ -137,6 +141,19 @@ void PropertyView::selectionChanged(const QItemSelection& selected,
         m_molecule->undoMolecule()->setAtomSelected(std::get<2>(torsion), true);
         m_molecule->undoMolecule()->setAtomSelected(std::get<3>(torsion), true);
       }
+    } else if (m_type == PropertyType::ResidueType) {
+      // select all the atoms in the residue
+      if (m_model != nullptr) {
+        const auto residue = m_molecule->residue(rowNum);
+        auto atoms = residue.residueAtoms();
+        for (Index i = 0; i < atoms.size(); ++i) {
+          const auto atom = atoms[i];
+          m_molecule->undoMolecule()->setAtomSelected(atom.index(), true);
+        }
+      }
+    } else if (m_type == PropertyType::ConformerType) {
+      // selecting a row means switching to that conformer
+      m_molecule->setCoordinate3d(rowNum);
     }
   } // end loop through selected