Mercurial > repos > miller-lab > genome_diversity

comparison draw_variants.py @ 27:8997f2ca8c7a

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
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)