Mercurial > repos > miller-lab > genome_diversity
comparison draw_variants.py @ 27:8997f2ca8c7a
Update to Miller Lab devshed revision bae0d3306d3b
author | Richard Burhans <burhans@bx.psu.edu> |
---|---|
date | Mon, 15 Jul 2013 10:47:35 -0400 |
parents | 248b06e86022 |
children | a631c2f6d913 |
comparison
equal
deleted
inserted
replaced
26:91e835060ad2 | 27:8997f2ca8c7a |
---|---|
1 #!/usr/bin/env python | 1 #!/usr/bin/env python |
2 | 2 |
3 import gd_util | |
3 import sys | 4 import sys |
4 import subprocess | |
5 from Population import Population | 5 from Population import Population |
6 | |
7 def run_program(prog, args, stdout_file=None, space_to_tab=False): | |
8 #print "args: ", ' '.join(args) | |
9 p = subprocess.Popen(args, bufsize=-1, executable=prog, stdin=None, stdout=subprocess.PIPE, stderr=subprocess.PIPE) | |
10 (stdoutdata, stderrdata) = p.communicate() | |
11 rc = p.returncode | |
12 | |
13 if stdout_file is not None: | |
14 with open(stdout_file, 'w') as ofh: | |
15 lines = stdoutdata.split('\n') | |
16 for line in lines: | |
17 line = line.strip() | |
18 if line: | |
19 if space_to_tab: | |
20 line = line.replace(' ', '\t') | |
21 print >> ofh, line | |
22 | |
23 if rc != 0: | |
24 print >> sys.stderr, "FAILED: rc={0}: {1}".format(rc, ' '.join(args)) | |
25 print >> sys.stderr, stderrdata | |
26 sys.exit(1) | |
27 | 6 |
28 ################################################################################ | 7 ################################################################################ |
29 | 8 |
30 if len(sys.argv) < 10: | 9 if len(sys.argv) != 10: |
31 print >> sys.stderr, "Usage" | 10 gd_util.die('Usage') |
32 sys.exit(1) | |
33 | 11 |
34 snp_input, indel_input, coverage_input, annotation_input, indiv_input, ref_name, min_coverage, output = sys.argv[1:9] | 12 snp_input, indel_input, coverage_input, annotation_input, indiv_input, ref_name, min_coverage, output, ind_arg = sys.argv[1:] |
35 individual_metadata = sys.argv[9:] | |
36 | 13 |
37 p_total = Population() | 14 p_total = Population() |
38 p_total.from_tag_list(individual_metadata) | 15 p_total.from_wrapped_dict(ind_arg) |
39 | 16 |
40 p1 = Population() | 17 p1 = Population() |
41 p1.from_population_file(indiv_input) | 18 p1.from_population_file(indiv_input) |
42 if not p_total.is_superset(p1): | 19 if not p_total.is_superset(p1): |
43 print >> sys.stderr, 'There is an individual in the population individuals that is not in the SNP table' | 20 gd_util.die('There is an individual in the population individuals that is not in the SNP table') |
44 sys.exit(1) | |
45 | 21 |
46 ################################################################################ | 22 ################################################################################ |
47 | 23 |
48 prog = 'mk_Ji' | 24 prog = 'mk_Ji' |
49 | 25 |
50 args = [ prog ] | 26 args = [ prog ] |
51 | |
52 args.append(snp_input) | 27 args.append(snp_input) |
53 args.append(indel_input) | 28 args.append(indel_input) |
54 args.append(coverage_input) | 29 args.append(coverage_input) |
55 args.append(annotation_input) | 30 args.append(annotation_input) |
56 args.append(min_coverage) | 31 args.append(min_coverage) |
57 args.append(ref_name) | 32 args.append(ref_name) |
58 | 33 |
59 for tag in p1.tag_list(): | 34 for tag in p1.tag_list(): |
60 args.append(tag) | 35 args.append(tag) |
61 | 36 |
62 run_program(prog, args, stdout_file='mk_Ji.out') | 37 with open('mk_Ji.out', 'w') as fh: |
38 gd_util.run_program(prog, args, stdout=fh) | |
63 | 39 |
64 ################################################################################ | 40 ################################################################################ |
65 | 41 |
66 prog = 'varplot' | 42 prog = 'varplot' |
67 | 43 |
68 args = [ prog ] | 44 args = [ prog ] |
69 args.append('-w') | 45 args.append('-w') |
70 args.append('3') | 46 args.append(3) |
71 args.append('-s') | 47 args.append('-s') |
72 args.append('0.3') | 48 args.append(0.3) |
73 args.append('-g') | 49 args.append('-g') |
74 args.append('0.2') | 50 args.append(0.2) |
75 args.append('mk_Ji.out') | 51 args.append('mk_Ji.out') |
76 | 52 |
77 run_program(prog, args, stdout_file=output) | 53 with open(output, 'w') as fh: |
54 gd_util.run_program(prog, args, stdout=fh) | |
55 | |
78 sys.exit(0) | 56 sys.exit(0) |