How to use the mdtraj.rmsd function in mdtraj

t1_a = md.load(get_fn('traj.h5'), stride=10)
    t2_a = md.load(get_fn('traj.h5'), stride=10)
    t1_b = md.load(get_fn('traj.h5'), stride=10)
    t2_b = md.load(get_fn('traj.h5'), stride=10)
    # don't center t1, and use it without precentered
    # explicitly center t2, and use *with* precentered

    t2_a.center_coordinates()
    t2_b.center_coordinates()

    for parallel in [True, False]:
        for i in range(t1_b.n_frames):
            ref = np.zeros(t1_a.n_frames)
            for j in range(t1_a.n_frames):
                ref[j] = rmsd_qcp(t1_a.xyz[j], t1_b.xyz[i])
            val1 = md.rmsd(t1_a, t1_b, i, parallel=parallel, precentered=False)
            val2 = md.rmsd(t2_a, t2_b, i, parallel=parallel, precentered=True)

            eq(ref, val1, decimal=3)
            eq(val1, val2, decimal=4)

def test_rmsd_atom_indices(get_fn):
    native = md.load(get_fn('native.pdb'))
    t1 = md.load(get_fn('traj.h5'))

    atom_indices = np.arange(10)
    dist1 = md.rmsd(t1, native, atom_indices=atom_indices)

    t2 = md.load(get_fn('traj.h5'))
    t2.restrict_atoms(atom_indices)
    native.restrict_atoms(atom_indices)
    dist2 = md.rmsd(t2, native)

    eq(dist1, dist2)

raise exception.ImproperlyConfigured(
                "When --features is specified, --topology is unneccessary.")
        if args.atoms:
            raise exception.ImproperlyConfigured(
                "Option --atoms is only meaningful when clustering "
                "trajectories.")
        if not args.cluster_distance:
            raise exception.ImproperlyConfigured(
                "Option --cluster-distance is required when clustering "
                "features.")

    elif args.trajectories and args.topologies:
        args.trajectories = expand_files(args.trajectories)

        if not args.cluster_distance or args.cluster_distance == 'rmsd':
            args.cluster_distance = md.rmsd
        else:
            raise exception.ImproperlyConfigured(
                "Option --cluster-distance must be rmsd when clustering "
                "trajectories.")

        if not args.atoms:
            raise exception.ImproperlyConfigured(
                "Option --atoms is required when clustering trajectories.")
        elif len(args.atoms) == 1:
            args.atoms = args.atoms * len(args.trajectories)
        elif len(args.atoms) != len(args.trajectories):
            raise exception.ImproperlyConfigured(
                "Flag --atoms must be provided either once (selection is "
                "applied to all trajectories) or the same number of times "
                "--trajectories is supplied.")

def rmsd_ref_state(
            self,
            ref_frame: int
    ):
        self._rmsd_ref_state = ref_frame
        self.rmsd = mdtraj.rmsd(self, self, ref_frame)

xyz = xyz.reshape((1, xyz.shape[0], 3)) / 10.0
        # write to trajectory file
        mode = 'a' if os.path.isfile(self.filename) else 'w'
        t = mdtraj.formats.hdf5.HDF5TrajectoryFile(self.filename, mode=mode)
        t.write(xyz, time=len(t))
        t.close()

        if update_rmsd:
            self._xyz = np.append(self._xyz, xyz, axis=0)
            self._time = np.arange(len(self._xyz))

            self.mdtraj._xyz = self._xyz
            self.mdtraj._time = self.time
            new = self.mdtraj[-1]
            previous = self.mdtraj[-2]
            next_drmsd = mdtraj.rmsd(new, previous) * 10.0
            next_rmsd = mdtraj.rmsd(new, self.mdtraj[self.rmsd_ref_state]) * 10.0
        else:
            next_drmsd = [0.0]
            next_rmsd = [0.0]
        self.drmsd.append(float(next_drmsd[0]))
        self.rmsd.append(float(next_rmsd[0]))
        self.energy.append(energy)
        self.chi2r.append(energy_fret)
        if verbose:
            print("%.3f\t%.3f\t%.4f\t%.4f" % (energy, energy_fret, next_rmsd[0], next_drmsd[0]))

distance is minimized
       
    """
    path_out = []
    now = start
    assert now.n_frames == 1

    if selection is None:
       atom_indices = _np.arange(now.n_atoms)
    else:
       atom_indices = selection

    # For the list of candidates, extract the closest one
    history = now
    for ii, cands in enumerate(path_of_candidates):
        closest_to_now = _np.argmin(_md.rmsd(cands, now, atom_indices=atom_indices))
        path_out.append(closest_to_now)
        #print("choose frame %u from %u cands"%(path_out[-1], len(cands)))
        now = cands[closest_to_now]
        history = history.join(now)
        if history_aware:
           history.superpose(history, atom_indices=atom_indices)
           xyz = history.xyz.mean(0)
           now = _md.Trajectory(xyz, history.top)
    return path_out

def find_best(target, reference, atom_indices=None):
    """
    target and reference are both of type mdtraj.Trajectory
    reference is of length 1, target of arbitrary length

    returns a Structure object and the index within the trajectory

    Examples
    --------

    >>> import chisurf.settings as mfm
    >>> times = times = mfm.TrajectoryFile('./test/data/structure/2807_8_9_b.h5', reading_routine='r', stride=1)
    >>> find_best(times.mdtraj, times.mdtraj[2])
    (2, )
    """
    rmsds = md.rmsd(target, reference, atom_indices=atom_indices)
    iMin = np.argmin(rmsds)
    return iMin, target[iMin]

def _get_distance_method(metric):
    if metric == 'rmsd':
        return md.rmsd
    elif isinstance(metric, str):
        try:
            import msmbuilder.libdistance as libdistance
        except ImportError:
            raise ImproperlyConfigured(
                "To use '{}' as a clustering metric, STAG ".format(metric) +
                "uses MSMBuilder3's libdistance, but we weren't able to " +
                "import msmbuilder.libdistance.")

        def f(X, Y):
            return libdistance.dist(X, Y, metric)
        return f
    elif callable(metric):
        return metric
    else:
        raise ImproperlyConfigured(

def minimize_rmsd2ref_in_sample(sample, ref):
    # Candidate selection

    out_geoms = None
    for cand_geoms in sample:
        igeom = cand_geoms[(_np.argmin(_md.rmsd(cand_geoms, ref)))]
        if out_geoms is None:
            out_geoms = igeom
        else:
            out_geoms = out_geoms.join(igeom)

    return out_geoms.superpose(ref)