mmagnus / rna-tools

    if method == 'fit':

        return (pymol.cmd.fit('s1', 's2'), 'fit')

def calc_rmsd(a, b, target_selection, target_ignore_selection, model_selection, model_ignore_selection, way, verbose):

    """

    Calculate RMSD between two XYZ files

    by: Jimmy Charnley Kromann <[email protected]> and Lars Andersen Bratholm <[email protected]>

    project: https://github.com/charnley/rmsd

    license: https://github.com/charnley/rmsd/blob/master/LICENSE

    a is model

    b is target

    :params: a = filename of structure a

    :params: b = filename of structure b

    :return: rmsd, number of atoms

    """

    if verbose: print('in:', a)

    atomsP, P = get_coordinates(a, model_selection, model_ignore_selection, 'pdb', True, way)

    atomsQ, Q = get_coordinates(b, target_selection,target_ignore_selection,  'pdb', True, way)

    if verbose:

        print(atomsP, P)

        print(atomsQ, Q)

    if atomsQ != atomsP:

        print('Error: number of atoms is not equal target (' + b + '):' + str(atomsQ) + ' vs model (' + a + '):' + str(atomsP))

        return (-1,0) # skip this RNA

    # Calculate 'dumb' RMSD

    normal_rmsd = rmsd(P, Q)

    # Create the centroid of P and Q which is the geometric center of a

    # N-dimensional region and translate P and Q onto that center.

    # http://en.wikipedia.org/wiki/Centroid

    Pc = centroid(P)

    Qc = centroid(Q)

    P -= Pc

    Q -= Qc

    if False:

        V = rotate(P, Q)

        V += Qc

        write_coordinates(atomsP, V)

        quit()

    return round(kabsch_rmsd(P, Q),2), atomsP

def get_parser():

    import argparse

    if method == 'fit':

        return (pymol.cmd.fit('s1', 's2'), 'fit')

def calc_rmsd(a, b, target_selection, target_ignore_selection, model_selection, model_ignore_selection, way, verbose):

    """

    Calculate RMSD between two XYZ files

    by: Jimmy Charnley Kromann <[email protected]> and Lars Andersen Bratholm <[email protected]>

    project: https://github.com/charnley/rmsd

    license: https://github.com/charnley/rmsd/blob/master/LICENSE

    a is model

    b is target

    :params: a = filename of structure a

    :params: b = filename of structure b

    :return: rmsd, number of atoms

    """

    if verbose: print('in:', a)

    atomsP, P, atoms = get_coordinates(a, model_selection, model_ignore_selection, 'pdb', True, way)

    atomsQ, Q, atoms = get_coordinates(b, target_selection,target_ignore_selection,  'pdb', True, way)

    if verbose:

        print(atomsP, P)

        print(atomsQ, Q)

    if atomsQ != atomsP:

        print('Error: number of atoms is not equal target (' + b + '):' + str(atomsQ) + ' vs model (' + a + '):' + str(atomsP))

        return (-1,0) # skip this RNA

    # Calculate 'dumb' RMSD

    normal_rmsd = rmsd(P, Q, atoms)

    # Create the centroid of P and Q which is the geometric center of a

    # N-dimensional region and translate P and Q onto that center.

    # http://en.wikipedia.org/wiki/Centroid

    Pc = centroid(P)

    Qc = centroid(Q)

    P -= Pc

    Q -= Qc

    if False:

        V = rotate(P, Q)

        V += Qc

        write_coordinates(atomsP, V)

        quit()

    return round(kabsch_rmsd(P, Q, atoms),2), atomsP

def get_parser():

    import argparse