Mercurial > repos > chemteam > mdanalysis_cosine_analysis
changeset 2:e39bc1f90d8f draft
"planemo upload for repository https://github.com/galaxycomputationalchemistry/galaxy-tools-compchem/ commit 12823af06b7926cc56aaf9b59cfea9f16d342b8c"
author | chemteam |
---|---|
date | Thu, 06 Feb 2020 19:41:02 -0500 |
parents | d18e7db42633 |
children | 9a3a1f698fc6 |
files | distance.py distance_multiple.py distance_single.py test-data/list1.txt test-data/list2.txt |
diffstat | 5 files changed, 120 insertions(+), 59 deletions(-) [+] |
line wrap: on
line diff
--- a/distance.py Mon Oct 07 12:50:44 2019 -0400 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,59 +0,0 @@ -#!/usr/bin/env python - -import argparse -import sys - -import MDAnalysis as mda - -import matplotlib -matplotlib.use('Agg') # noqa -import matplotlib.pyplot as plt - -import numpy as np - - -def parse_command_line(argv): - parser = argparse.ArgumentParser() - parser.add_argument('--itraj', help='input traj') - parser.add_argument('--istr', help='input str') - parser.add_argument('--itrajext', help='input traj ext') - parser.add_argument('--istrext', help='input str ext') - parser.add_argument('--isegid1', help='segid 1') - parser.add_argument('--iresid1', help='resid 1') - parser.add_argument('--iname1', help='name 1') - parser.add_argument('--isegid2', help='segid 2') - parser.add_argument('--iresid2', help='resid 2') - parser.add_argument('--iname2', help='name 2') - parser.add_argument('--output', help='output') - parser.add_argument('--odistance_plot', help='odistance plot') - return parser.parse_args() - - -args = parse_command_line(sys.argv) - -atom1 = "(segid %s and resid %s and name %s)" % \ - (args.isegid1, args.iresid1, args.iname1) -atom2 = "(segid %s and resid %s and name %s)" % \ - (args.isegid2, args.iresid2, args.iname2) - -u = mda.Universe(args.istr, args.itraj, - topology_format=args.istrext, format=args.itrajext) -x = u.select_atoms(atom1) -y = u.select_atoms(atom2) - -with open(args.output, 'w') as f: - for t in u.trajectory: - r = x.positions - y.positions - d = np.linalg.norm(r) - f.write(str(t.frame) + '\t ') - f.write(str(d) + '\n') - -with open(args.output) as f: - g = [xtmp.strip() for xtmp in f] - data = [tuple(map(float, xtmp.split())) for xtmp in g[0:]] - time = [xtmp[0] for xtmp in data] - distance = [xtmp[1] for xtmp in data] - plt.plot(time, distance) - plt.xlabel('Frame No.') - plt.ylabel(r'Distance ($\AA$)') - plt.savefig(args.odistance_plot, format='png')
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/distance_multiple.py Thu Feb 06 19:41:02 2020 -0500 @@ -0,0 +1,56 @@ +import argparse +import sys + +import MDAnalysis as mda +from MDAnalysis.analysis import distances + +import numpy as np + + +def parse_command_line(argv): + parser = argparse.ArgumentParser() + parser.add_argument('--itraj', help='input traj') + parser.add_argument('--istr', help='input str') + parser.add_argument('--itrajext', help='input traj ext') + parser.add_argument('--istrext', help='input str ext') + parser.add_argument('--list1', help='list 2') + parser.add_argument('--list2', help='list 2') + parser.add_argument('--output', help='output') + parser.add_argument('--header', dest='header', action='store_true') + return parser.parse_args() + + +args = parse_command_line(sys.argv) + +u = mda.Universe(args.istr, args.itraj, + topology_format=args.istrext, format=args.itrajext) + +list1 = np.loadtxt(args.list1, dtype=str, delimiter="\t", ndmin=1) +list2 = np.loadtxt(args.list2, dtype=str, delimiter="\t", ndmin=1) + +sel1 = [u.select_atoms(selection) for selection in list1] +sel2 = [u.select_atoms(selection) for selection in list2] + +d = np.empty((u.trajectory.n_frames, list1.shape[0], list2.shape[0]),) + +for ts in u.trajectory: + c_o_m1 = np.array([selection.center_of_mass() for selection in sel1]) + c_o_m2 = np.array([selection.center_of_mass() for selection in sel2]) + distances.distance_array(c_o_m1, c_o_m2, result=d[ts.frame]) + +d = np.hstack(( + np.array(np.reshape(np.arange( + 0, d.shape[0]), (d.shape[0], 1)), dtype=int), # add column w frame + np.reshape(d, (d.shape[0], d.shape[1] * d.shape[2])) +)) + +if args.header: + header = 'Frame\t' + '\t'.join( + ['-'.join(pair) for pair in zip( + sum([[n, ] * len(list2) for n in list1], []), + list(list2) * len(list1),)]).replace(' ', '_') +else: + header = '' + +np.savetxt(args.output, d, header=header, comments='', + fmt=['%d'] + ['%f'] * (d.shape[1] - 1), delimiter='\t')
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/distance_single.py Thu Feb 06 19:41:02 2020 -0500 @@ -0,0 +1,62 @@ +#!/usr/bin/env python + +import argparse +import sys + +import MDAnalysis as mda + +import matplotlib +matplotlib.use('Agg') # noqa +import matplotlib.pyplot as plt + +import numpy as np + + +def parse_command_line(argv): + parser = argparse.ArgumentParser() + parser.add_argument('--itraj', help='input traj') + parser.add_argument('--istr', help='input str') + parser.add_argument('--itrajext', help='input traj ext') + parser.add_argument('--istrext', help='input str ext') + parser.add_argument('--isegid1', help='segid 1') + parser.add_argument('--iresid1', help='resid 1') + parser.add_argument('--iname1', help='name 1') + parser.add_argument('--isegid2', help='segid 2') + parser.add_argument('--iresid2', help='resid 2') + parser.add_argument('--iname2', help='name 2') + parser.add_argument('--output', help='output') + parser.add_argument('--odistance_plot', help='odistance plot') + parser.add_argument('--header', dest='header', action='store_true') + return parser.parse_args() + + +args = parse_command_line(sys.argv) + +atom1 = "(segid %s and resid %s and name %s)" % \ + (args.isegid1, args.iresid1, args.iname1) +atom2 = "(segid %s and resid %s and name %s)" % \ + (args.isegid2, args.iresid2, args.iname2) + +u = mda.Universe(args.istr, args.itraj, + topology_format=args.istrext, format=args.itrajext) +x = u.select_atoms(atom1) +y = u.select_atoms(atom2) + +with open(args.output, 'w') as f: + if args.header: + f.write('Frame\tDistance') + for t in u.trajectory: + r = x.positions - y.positions + d = np.linalg.norm(r) + f.write(str(t.frame) + '\t ') + f.write(str(d) + '\n') + +with open(args.output) as f: + g = [xtmp.strip() for xtmp in f] + data = [tuple(map(float, xtmp.split())) for xtmp in g[0:]] + time = [xtmp[0] for xtmp in data] + distance = [xtmp[1] for xtmp in data] + plt.plot(time, distance) + plt.xlabel('Frame No.') + plt.ylabel(r'Distance ($\AA$)') + plt.savefig(args.odistance_plot, format='png')