From e080a99087dd18d84a2d26119546197ba0e9e9a6 Mon Sep 17 00:00:00 2001
From: Christian Faber <chris-faber@gmx.de>
Date: Mon, 21 Oct 2024 09:47:05 +0200
Subject: [PATCH] added rmsd per residue

---
 src/BioHelpers_FABER/rmsd.py | 12 ++++++------
 test_rmsd.py                 |  1 -
 2 files changed, 6 insertions(+), 7 deletions(-)

diff --git a/src/BioHelpers_FABER/rmsd.py b/src/BioHelpers_FABER/rmsd.py
index 2a066e3..5f4281d 100644
--- a/src/BioHelpers_FABER/rmsd.py
+++ b/src/BioHelpers_FABER/rmsd.py
@@ -117,10 +117,10 @@ def rmsd_from_rot_tran(
     :return: RMSD like value (:math:`\sqrt{E}`)
     :rtype: float
     """
-    print(trans.shape)
-    print(np.matmul(rot, x).shape, y.shape)
-
-    return 0.0
+    num_vec = x.shape[1]
+    intermediate = np.matmul(rot, x) + trans[:, np.newaxis] + y
+    squared = np.trace(np.matmul(np.transpose(intermediate), intermediate)) / num_vec
+    return np.sqrt(squared)
 
 
 def rmsd(x: np.ndarray, y: np.ndarray) -> float:
@@ -220,9 +220,9 @@ def rmsd_per_residue(file1: str, file2: str) -> list[float]:
             res2_AA = [list(a.get_vector()) for a in res2.get_atoms()]
             res1_AA = np.transpose(np.array(res1_AA))
             res2_AA = np.transpose(np.array(res2_AA))
-            rmsd_from_rot_tran(res1_AA, res2_AA, rot, trans)
+            values.append(rmsd_from_rot_tran(res1_AA, res2_AA, rot, trans))
 
-    return [0.0]
+    return values
 
 
 if __name__ == "__main__":
diff --git a/test_rmsd.py b/test_rmsd.py
index d42a0a5..32a194c 100644
--- a/test_rmsd.py
+++ b/test_rmsd.py
@@ -32,7 +32,6 @@ def main():
     rot, tran = sup.get_rotran()
     # rmsd per residue
     print(rmsd.rmsd_per_residue(file1, file2))
-    print("Hallo Welt", val, val2)
 
 
 if __name__ == "__main__":
-- 
GitLab