Source code for schrodinger.test.stu.outcomes.compareDesmondSystems

"""script used to compare two desmond systems"""

from past.utils import old_div

from schrodinger.application.desmond.cms import get_box
from schrodinger.application.desmond.cms import get_boxvolume
from schrodinger.application.desmond.constants import CT_TYPE
from schrodinger.structure import MaestroReader as _MaestroReader
from schrodinger.structutils.rmsd import superimpose as _superimpose


class System(dict):
    """
    A very basic class to hold cms data.  Inherits the dict class to do most
    of the heavy lifting.
    """

    def __init__(self, fileName):
        """Initializes the class by loading the cts of the filename"""
        self.name = fileName
        for ct in _MaestroReader(fileName):
            self[ct.property[CT_TYPE]] = ct

    def getFFIOblock(self):
        """Gets a dictionary of ffio block values"""

    def _findFFIOdiffs(self, other):
        """Returns a list of diffs between the ffio blocks of self and other"""

    def _findCTdiffs(self, other, name):
        """
            Returns a list of diffs between the named ct in self and other.
            Only the solute CT has its positions compared.  The other CTs use
            comparisons of number of molecules.
            """
        diffs = []
        #Does this ct type exist in both?
        if name not in self:
            diffs.append(f"{name} missing from {self.name}")
        if name not in other:
            diffs.append(f"{name} missing from {other.name}")
        #Compare box sizes
        size1 = get_boxvolume(get_box(self[name]))
        size2 = get_boxvolume(get_box(other[name]))
        sizeDiff = abs(old_div((size1 - size2), size1))
        if sizeDiff > 0.04:
            diffs.append("Box size difference is %.1f%%" % (sizeDiff * 100))
        #Compare number of molecules
        nmol1 = self[name].mol_total
        nmol2 = other[name].mol_total
        if name == "solute":
            if nmol1 != nmol2:
                diffs.append(
                    "Different number of %s molecules: %s:%s and %s:%s" %
                    (name, self.name, nmol1, other.name, nmol2))
            else:
                diffs.extend(self._RMSDdiff(other, name))
        elif abs(old_div(float(nmol1 - nmol2), nmol1)) > 0.05:
            #Need to cast as float - int division always results in an int.
            diffs.append("Different number of %s molecules: %s:%s and %s:%s" %
                         (name, self.name, nmol1, other.name, nmol2))

        return diffs

    def _RMSDdiff(self, other, name):
        atoms = list(range(1, len(self[name].atom) + 1))
        rmsd = _superimpose(self[name], atoms, other[name], atoms)

        if rmsd > 1.04:
            return ["RMSD above tolerance (RMSD = %s)" % rmsd]
        else:
            return []

    def diff(self, other):
        """Returns a list of differences between self and other."""
        diffs = []
        for name in set(list(self.keys()) + list(other)):
            diffs.extend(self._findCTdiffs(other, name))
        return diffs


def compareDesmondSystems(filename1, filename2, log=None):
    """
    Function used to compare two desmond systems.  Takes as arguments
    the file names of the two systems and an optional argument
    -log=loglocation.  If the log file location is provided,
    differences are saved to file.  Otherwise, they are printed to screen
    """
    log = ""
    try:
        cms1 = System(filename1)
        cms2 = System(filename2)
    except:
        msg = f"Problem reading files {filename1} and {filename2}"
        raise AssertionError(msg)
    diffs = cms1.diff(cms2)

    if diffs:
        if log:
            file = open(log, "w")
            file.write("The cms files are not identical")
            for line in diffs:
                file.write(line + "\n")
            raise AssertionError("See diffs in %s" % log)
        else:
            msg = "The cms files are not identical:"
            for line in diffs:
                msg += "\n\t" + line
            raise AssertionError(msg)
    return True


if __name__ == "__main__":
    from sys import argv

    if compareDesmondSystems(argv[1:]):
        print("   Same (pass)")
    else:
        print("   Difference (fail)")