comparison distance_multiple.py @ 2:b348dfa55e0a draft

"planemo upload for repository https://github.com/galaxycomputationalchemistry/galaxy-tools-compchem/ commit 12823af06b7926cc56aaf9b59cfea9f16d342b8c"
author chemteam
date Thu, 06 Feb 2020 19:40:37 -0500
parents
children
comparison
equal deleted inserted replaced
1:ce0728b92289 2:b348dfa55e0a
1 import argparse
2 import sys
3
4 import MDAnalysis as mda
5 from MDAnalysis.analysis import distances
6
7 import numpy as np
8
9
10 def parse_command_line(argv):
11 parser = argparse.ArgumentParser()
12 parser.add_argument('--itraj', help='input traj')
13 parser.add_argument('--istr', help='input str')
14 parser.add_argument('--itrajext', help='input traj ext')
15 parser.add_argument('--istrext', help='input str ext')
16 parser.add_argument('--list1', help='list 2')
17 parser.add_argument('--list2', help='list 2')
18 parser.add_argument('--output', help='output')
19 parser.add_argument('--header', dest='header', action='store_true')
20 return parser.parse_args()
21
22
23 args = parse_command_line(sys.argv)
24
25 u = mda.Universe(args.istr, args.itraj,
26 topology_format=args.istrext, format=args.itrajext)
27
28 list1 = np.loadtxt(args.list1, dtype=str, delimiter="\t", ndmin=1)
29 list2 = np.loadtxt(args.list2, dtype=str, delimiter="\t", ndmin=1)
30
31 sel1 = [u.select_atoms(selection) for selection in list1]
32 sel2 = [u.select_atoms(selection) for selection in list2]
33
34 d = np.empty((u.trajectory.n_frames, list1.shape[0], list2.shape[0]),)
35
36 for ts in u.trajectory:
37 c_o_m1 = np.array([selection.center_of_mass() for selection in sel1])
38 c_o_m2 = np.array([selection.center_of_mass() for selection in sel2])
39 distances.distance_array(c_o_m1, c_o_m2, result=d[ts.frame])
40
41 d = np.hstack((
42 np.array(np.reshape(np.arange(
43 0, d.shape[0]), (d.shape[0], 1)), dtype=int), # add column w frame
44 np.reshape(d, (d.shape[0], d.shape[1] * d.shape[2]))
45 ))
46
47 if args.header:
48 header = 'Frame\t' + '\t'.join(
49 ['-'.join(pair) for pair in zip(
50 sum([[n, ] * len(list2) for n in list1], []),
51 list(list2) * len(list1),)]).replace(' ', '_')
52 else:
53 header = ''
54
55 np.savetxt(args.output, d, header=header, comments='',
56 fmt=['%d'] + ['%f'] * (d.shape[1] - 1), delimiter='\t')