# Compute RMSD between two molecules#

Let us compute the root-mean-square deviation (RMSD) of atomic positions between the atoms of two superimposed molecules.

```from math import sqrt

def compute_rmsd(mol1_indexer, mol2_indexer):
"""Computes RMSD between two superimposed molecules and returns RMSD in angstrom"""

# check if indexers are of the same size
if len(mol1_indexer) != len(mol2_indexer):
print('Molecules are not the same')
return -1

sqrd_sum = SBQuantity.squareLength(0.0)
# go through all atoms in the indexer
for i in range(len(mol1_indexer)):
if mol1_indexer[i].elementType != mol2_indexer[i].elementType:
print('Atoms are not the same')
return -1
# get positions of atoms from indexers
pos1 = mol1_indexer[i].getPosition()
pos2 = mol2_indexer[i].getPosition()
# compute distance between atoms
# and accumulate squared norm for distances
sqrd_sum += (pos1 - pos2).norm2()

# take a square root of the sum of squared distances divided by the number of atoms
rmsd = SBQuantity.sqrt(sqrd_sum / len(mol1_indexer))

# return RMSD in angstrom
return rmsd.angstrom.value

def compute_rmsd_2mol(mol1, mol2):
"""Computes several RMSD between two superimposed molecules:
1) between all atoms,
2) between backbone atoms,
3) between alpha-Carbon atoms.
"""
# get indexers of all atoms
mol1_all_indexer      = mol1.getNodes('node.type atom')
mol2_all_indexer      = mol2.getNodes('n.t a')
# get indexers of all atoms in the backbone
mol1_backbone_indexer = mol1.getNodes('node.type atom in node.type backbone')
mol2_backbone_indexer = mol2.getNodes('n.t a in n.t bb')
# get indexers of all alpha Carbons in the backbone
mol1_CA_indexer       = mol1.getNodes('"CA" in node.type backbone')
mol2_CA_indexer       = mol2.getNodes('"CA" in n.t bb')

rmsd_all      = compute_rmsd(mol1_all_indexer, mol2_all_indexer);
rmsd_backbone = compute_rmsd(mol1_backbone_indexer, mol2_backbone_indexer);
rmsd_CA       = compute_rmsd(mol1_CA_indexer, mol2_CA_indexer);

if (rmsd_all != -1):      print("RMSD (all atoms): %e Angstrom" % rmsd_all)
if (rmsd_backbone != -1): print("RMSD (backbone) : %e Angstrom" % rmsd_backbone)
if (rmsd_CA != -1):       print("RMSD (C_alpha)  : %e Angstrom" % rmsd_CA)

# return a list of RMSDs
return [rmsd_all, rmsd_backbone, rmsd_CA]

# get all molecules in the active document
indexer = SAMSON.getNodes('node.type molecule')

if len(indexer) < 2:
print("Less than two molecules are present in the data graph")
else:
# we assume here that we need to compute RMSD between the first two molecules in the indexer
mol1 = indexer[0]    # get the first molecule
mol2 = indexer[1]    # get the second molecule
# compute RMSD between the two superimposed molecules
rmsd = compute_rmsd_2mol(mol1, mol2)
```