Source code for moldesign.interfaces.qtrfit

# Copyright 2016 Autodesk Inc.
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
#     http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.

"""Interface to quaternion-based atom alignment routines.

This code has been translated from Fortran to Python
"""

import numpy as np

[docs]def qtrfit(mol, pos, ref, massweight=True): r""" Transform the atomic positions in ``pos`` to minimize the overall distance to those in ``ref`` This function minimizes the distance between pos and ref. If `massweight`=True, we minimize the mass-weighted distance: ..math:: D = \sqrt{ \sum_{i=1}^{\text{N}} m_i \left| \mathbf{pos}_i - \mathbf{ref}_i \right|^2 } Otherwise, the root-squared distance is used: ..math:: D = \sqrt{ \sum_{i=1}^{\text{N}} \left| \mathbf{pos}_i - \mathbf{ref}_i \right|^2 } Args: mol (moldesign.Molecule): molecule that these arrays describe pos (u.Array): The Nx3 array to be rotated ref (u.Array): The Nx3 array to match massweight (bool): whether to use mass-weighted coordinates (the default) Returns: u.Array: translated and rotated version of pos with minimum RMSD References: TODO: reference to original qtrfit code (and license, if it exists) """ N = 1.0 * len(pos) if massweight: weights = mol.masses else: weights = np.ones(mol.num_atoms) totweight = weights.sum() pcenter = (weights[:, None] * pos).sum(axis=0) / totweight rcenter = (weights[:, None] * ref).sum(axis=0) / totweight p = pos - pcenter r = ref - rcenter # yadda yadda yadda return p + rcenter