Mercurial > repos > peterjc > mira_assembler

comparison tools/sr_assembly/mira.py @ 4:117cce3296af

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Uploaded wrapper v0.0.3, which is for MIRA v3.4.x which includes support for Ion Torrent. The Galaxy wrapper will no longer work with MIRA v3.2.x - if you are still using the old version of MIRA, please continue to use v0.0.2 of the wrapper.
author peterjc
date Wed, 21 Dec 2011 11:33:19 -0500
parents 298f5c1d9521
children 3e7eca1f5d04
comparison
equal deleted inserted replaced
3:298f5c1d9521 4:117cce3296af
9 9
10 def stop_err(msg, err=1): 10 def stop_err(msg, err=1):
11 sys.stderr.write(msg+"\n") 11 sys.stderr.write(msg+"\n")
12 sys.exit(err) 12 sys.exit(err)
13 13
14 def tcs_to_tabular(old, new):
15 in_handle = open(old, "rU")
16 out_handle = open(new, "w")
17 assert in_handle.readline() == "#TCS V1.0\n"
18 assert in_handle.readline() == "#\n"
19 assert in_handle.readline() == "# contig name padPos upadPos | B Q | tcov covA covC covG covT cov* | qA qC qG qT q* | S | Tags\n"
20 assert in_handle.readline() == "#\n"
21 out_handle.write("#%s\n" % "\t".join(["contig", "pasPos", "upadPos", "B", "Q",
22 "tcov", "covA", "covC", "covG", "covT", "cov*",
23 "qA", "qC", "qG", "qT", "q*", "S", "Tags"]))
24 for line in in_handle:
25 parts = line.rstrip("\n").split(None,22)
26 assert parts[3] == parts[6] == parts[13] == parts[19] == parts[21] == "|"
27 wanted = parts[:3] + parts[4:6]+parts[7:13]+parts[14:19]+parts[20:21]+parts[22:]
28 out_handle.write("%s\n" % "\t".join(wanted))
29 out_handle.close()
30 in_handle.close()
31
32 def collect_output(temp, name): 14 def collect_output(temp, name):
33 n3 = (temp, name, name, name) 15 n3 = (temp, name, name, name)
34 f = "%s/%s_assembly/%s_d_results" % (temp, name, name) 16 f = "%s/%s_assembly/%s_d_results" % (temp, name, name)
35 if not os.path.isdir(f): 17 if not os.path.isdir(f):
36 stop_err("Missing output folder") 18 stop_err("Missing output folder")
37 if not os.listdir(f): 19 if not os.listdir(f):
38 stop_err("Empty output folder") 20 stop_err("Empty output folder")
21 missing = []
39 for old, new in [("%s/%s_out.unpadded.fasta" % (f, name), out_fasta), 22 for old, new in [("%s/%s_out.unpadded.fasta" % (f, name), out_fasta),
40 ("%s/%s_out.unpadded.fasta.qual" % (f, name), out_qual), 23 ("%s/%s_out.unpadded.fasta.qual" % (f, name), out_qual),
41 ("%s/%s_out.wig" % (f, name), out_wig), 24 ("%s/%s_out.wig" % (f, name), out_wig),
42 ("%s/%s_out.caf" % (f, name), out_caf), 25 ("%s/%s_out.caf" % (f, name), out_caf),
43 ("%s/%s_out.ace" % (f, name), out_ace)]: 26 ("%s/%s_out.ace" % (f, name), out_ace)]:
44 if not os.path.isfile(old): 27 if not os.path.isfile(old):
45 stop_err("Missing %s output file" % os.path.splitext(old)[-1]) 28 missing.append(os.path.splitext(old)[-1])
46 else: 29 else:
47 shutil.move(old, new) 30 shutil.move(old, new)
48 tcs_to_tabular("%s/%s_assembly/%s_d_results/%s_out.tcs" % n3, out_tcs) 31 if missing:
32 stop_err("Missing output files: %s" % ", ".join(missing))
49 33
50 def clean_up(temp, name): 34 def clean_up(temp, name):
51 folder = "%s/%s_assembly" % (temp, name) 35 folder = "%s/%s_assembly" % (temp, name)
52 if os.path.isdir(folder): 36 if os.path.isdir(folder):
53 shutil.rmtree(folder) 37 shutil.rmtree(folder)
54 38
55 #TODO - Run MIRA in /tmp or a configurable directory? 39 #TODO - Run MIRA in /tmp or a configurable directory?
56 #Currently Galaxy puts us somewhere safe like: 40 #Currently Galaxy puts us somewhere safe like:
57 #/opt/galaxy-dist/database/job_working_directory/846/ 41 #/opt/galaxy-dist/database/job_working_directory/846/
58 temp = "." 42 temp = "."
59 name, out_fasta, out_qual, out_tcs, out_ace, out_caf, out_wig, out_log = sys.argv[1:9] 43 name, out_fasta, out_qual, out_ace, out_caf, out_wig, out_log = sys.argv[1:8]
60 44
61 start_time = time.time() 45 start_time = time.time()
62 cmd = " ".join(sys.argv[9:]) 46 cmd_list =sys.argv[8:]
47 cmd = " ".join(cmd_list)
63 48
64 assert os.path.isdir(temp) 49 assert os.path.isdir(temp)
65 d = "%s_assembly" % name 50 d = "%s_assembly" % name
66 assert not os.path.isdir(d) 51 assert not os.path.isdir(d), "Path %s already exists" % d
67 try: 52 try:
68 #Check path access 53 #Check path access
69 os.mkdir(d) 54 os.mkdir(d)
70 except Exception, err: 55 except Exception, err:
71 sys.stderr.write("Error making directory %s\n%s" % (d, err)) 56 sys.stderr.write("Error making directory %s\n%s" % (d, err))
75 #print cmd 60 #print cmd
76 61
77 handle = open(out_log, "w") 62 handle = open(out_log, "w")
78 try: 63 try:
79 #Run MIRA 64 #Run MIRA
80 child = subprocess.Popen(sys.argv[9:], 65 child = subprocess.Popen(cmd_list,
81 stdout=handle, 66 stdout=handle,
82 stderr=subprocess.STDOUT) 67 stderr=subprocess.STDOUT)
83 except Exception, err: 68 except Exception, err:
84 sys.stderr.write("Error invoking command:\n%s\n\n%s\n" % (cmd, err)) 69 sys.stderr.write("Error invoking command:\n%s\n\n%s\n" % (cmd, err))
85 #TODO - call clean up? 70 #TODO - call clean up?
100 clean_up(temp, name) 85 clean_up(temp, name)
101 stop_err("Return error code %i from command:\n%s" % (return_code, cmd), 86 stop_err("Return error code %i from command:\n%s" % (return_code, cmd),
102 return_code) 87 return_code)
103 handle.close() 88 handle.close()
104 89
90 #print "Collecting output..."
105 collect_output(temp, name) 91 collect_output(temp, name)
92
93 #print "Cleaning up..."
106 clean_up(temp, name) 94 clean_up(temp, name)
95
107 print "Done" 96 print "Done"