Mercurial > repos > peterjc > mira4_assembler
comparison tools/mira4_0/mira4_bait.py @ 4:1713289d9908 draft default tip
v0.0.11 tweak for use with bioconda dependencies
author | peterjc |
---|---|
date | Thu, 10 Aug 2017 11:09:10 -0400 |
parents | 4eb32a3d67d1 |
children |
comparison
equal
deleted
inserted
replaced
3:a4f602cc3aa9 | 4:1713289d9908 |
---|---|
1 #!/usr/bin/env python | 1 #!/usr/bin/env python |
2 """A simple wrapper script to call MIRA4's mirabait and collect its output. | 2 """A simple wrapper script to call MIRA4's mirabait and collect its output. |
3 """ | 3 """ |
4 | |
5 from __future__ import print_function | |
6 | |
4 import os | 7 import os |
8 import shutil | |
9 import subprocess | |
5 import sys | 10 import sys |
6 import subprocess | |
7 import shutil | |
8 import time | 11 import time |
9 | 12 |
10 WRAPPER_VER = "0.0.5" #Keep in sync with the XML file | 13 WRAPPER_VER = "0.0.11" # Keep in sync with the XML file |
11 | |
12 def sys_exit(msg, err=1): | |
13 sys.stderr.write(msg+"\n") | |
14 sys.exit(err) | |
15 | 14 |
16 | 15 |
17 def get_version(mira_binary): | 16 def get_version(mira_binary): |
18 """Run MIRA to find its version number""" | 17 """Run MIRA to find its version number""" |
19 # At the commend line I would use: mira -v | head -n 1 | 18 # At the commend line I would use: mira -v | head -n 1 |
20 # however there is some pipe error when doing that here. | 19 # however there is some pipe error when doing that here. |
21 cmd = [mira_binary, "-v"] | 20 cmd = [mira_binary, "-v"] |
22 try: | 21 try: |
23 child = subprocess.Popen(cmd, | 22 child = subprocess.Popen(cmd, universal_newlines=True, |
24 stdout=subprocess.PIPE, | 23 stdout=subprocess.PIPE, |
25 stderr=subprocess.STDOUT) | 24 stderr=subprocess.STDOUT) |
26 except Exception, err: | 25 except Exception as err: |
27 sys.stderr.write("Error invoking command:\n%s\n\n%s\n" % (" ".join(cmd), err)) | 26 sys.stderr.write("Error invoking command:\n%s\n\n%s\n" % (" ".join(cmd), err)) |
28 sys.exit(1) | 27 sys.exit(1) |
29 ver, tmp = child.communicate() | 28 ver, tmp = child.communicate() |
30 del child | 29 del child |
31 #Workaround for -v not working in mirabait 4.0RC4 | 30 # Workaround for -v not working in mirabait 4.0RC4 |
32 if "invalid option" in ver.split("\n", 1)[0]: | 31 if "invalid option" in ver.split("\n", 1)[0]: |
33 for line in ver.split("\n", 1): | 32 for line in ver.split("\n", 1): |
34 if " version " in line: | 33 if " version " in line: |
35 line = line.split() | 34 line = line.split() |
36 return line[line.index("version")+1].rstrip(")") | 35 return line[line.index("version") + 1].rstrip(")") |
37 sys_exit("Could not determine MIRA version:\n%s" % ver) | 36 sys.exit("Could not determine MIRA version:\n%s" % ver) |
38 return ver.split("\n", 1)[0] | 37 return ver.split("\n", 1)[0] |
39 | 38 |
40 try: | 39 |
40 if "MIRA4" in os.environ: | |
41 mira_path = os.environ["MIRA4"] | 41 mira_path = os.environ["MIRA4"] |
42 except KeyError: | 42 mira_binary = os.path.join(mira_path, "mirabait") |
43 sys_exit("Environment variable $MIRA4 not set") | 43 if not os.path.isfile(mira_binary): |
44 mira_binary = os.path.join(mira_path, "mirabait") | 44 sys.exit("Missing mirabait under $MIRA4, %r\nFolder contained: %s" |
45 if not os.path.isfile(mira_binary): | 45 % (mira_binary, ", ".join(os.listdir(mira_path)))) |
46 sys_exit("Missing mirabait under $MIRA4, %r\nFolder contained: %s" | 46 else: |
47 % (mira_binary, ", ".join(os.listdir(mira_path)))) | 47 sys.stderr.write("DEBUG: Since $MIRA4 is not set, assuming mira binaries are on $PATH.\n") |
48 mira_path = None | |
49 mira_binary = "mirabait" | |
50 | |
48 mira_ver = get_version(mira_binary) | 51 mira_ver = get_version(mira_binary) |
49 if not mira_ver.strip().startswith("4.0"): | 52 if not mira_ver.strip().startswith("4.0"): |
50 sys_exit("This wrapper is for MIRA V4.0, not:\n%s" % mira_ver) | 53 sys.exit("This wrapper is for MIRA V4.0, not:\n%s" % mira_ver) |
51 if "-v" in sys.argv or "--version" in sys.argv: | 54 if "-v" in sys.argv or "--version" in sys.argv: |
52 print "%s, MIRA wrapper version %s" % (mira_ver, WRAPPER_VER) | 55 print("%s, MIRA wrapper version %s" % (mira_ver, WRAPPER_VER)) |
53 sys.exit(0) | 56 sys.exit(0) |
54 | 57 |
55 | 58 |
56 format, output_choice, strand_choice, kmer_length, min_occurance, bait_file, in_file, out_file = sys.argv[1:] | 59 format, output_choice, strand_choice, kmer_length, min_occurance, bait_file, in_file, out_file = sys.argv[1:] |
57 | 60 |
58 if format.startswith("fastq"): | 61 if format.startswith("fastq"): |
59 format = "fastq" | 62 format = "fastq" |
60 elif format == "mira": | 63 elif format == "mira": |
61 format = "maf" | 64 format = "maf" |
62 elif format != "fasta": | 65 elif format != "fasta": |
63 sys_exit("Was not expected format %r" % format) | 66 sys.exit("Was not expected format %r" % format) |
64 | 67 |
65 assert out_file.endswith(".dat") | 68 assert out_file.endswith(".dat") |
66 out_file_stem = out_file[:-4] | 69 out_file_stem = out_file[:-4] |
67 | 70 |
68 cmd_list = [mira_binary, "-f", format, "-t", format, | 71 cmd_list = [mira_binary, "-f", format, "-t", format, |
69 "-k", kmer_length, "-n", min_occurance, | 72 "-k", kmer_length, "-n", min_occurance, |
70 bait_file, in_file, out_file_stem] | 73 bait_file, in_file, out_file_stem] |
71 if output_choice == "pos": | 74 if output_choice == "pos": |
72 pass | 75 pass |
73 elif output_choice == "neg": | 76 elif output_choice == "neg": |
74 #Invert the selection... | 77 # Invert the selection... |
75 cmd_list.insert(1, "-i") | 78 cmd_list.insert(1, "-i") |
76 else: | 79 else: |
77 sys_exit("Output choice should be 'pos' or 'neg', not %r" % output_choice) | 80 sys.exit("Output choice should be 'pos' or 'neg', not %r" % output_choice) |
78 if strand_choice == "both": | 81 if strand_choice == "both": |
79 pass | 82 pass |
80 elif strand_choice == "fwd": | 83 elif strand_choice == "fwd": |
81 #Ingore reverse strand... | 84 # Ingore reverse strand... |
82 cmd_list.insert(1, "-r") | 85 cmd_list.insert(1, "-r") |
83 else: | 86 else: |
84 sys_exit("Strand choice should be 'both' or 'fwd', not %r" % strand_choice) | 87 sys.exit("Strand choice should be 'both' or 'fwd', not %r" % strand_choice) |
85 | 88 |
86 cmd = " ".join(cmd_list) | 89 cmd = " ".join(cmd_list) |
87 #print cmd | 90 # print cmd |
88 start_time = time.time() | 91 start_time = time.time() |
89 try: | 92 try: |
90 #Run MIRA | 93 # Run MIRA |
91 child = subprocess.Popen(cmd_list, | 94 child = subprocess.Popen(cmd_list, universal_newlines=True, |
92 stdout=subprocess.PIPE, | 95 stdout=subprocess.PIPE, |
93 stderr=subprocess.STDOUT) | 96 stderr=subprocess.STDOUT) |
94 except Exception, err: | 97 except Exception as err: |
95 sys.stderr.write("Error invoking command:\n%s\n\n%s\n" % (cmd, err)) | 98 sys.stderr.write("Error invoking command:\n%s\n\n%s\n" % (cmd, err)) |
96 sys.exit(1) | 99 sys.exit(1) |
97 #Use .communicate as can get deadlocks with .wait(), | 100 # Use .communicate as can get deadlocks with .wait(), |
98 stdout, stderr = child.communicate() | 101 stdout, stderr = child.communicate() |
99 assert stderr is None # Due to way we ran with subprocess | 102 assert stderr is None # Due to way we ran with subprocess |
100 run_time = time.time() - start_time | 103 run_time = time.time() - start_time |
101 return_code = child.returncode | 104 return_code = child.returncode |
102 print "mirabait took %0.2f minutes" % (run_time / 60.0) | 105 print("mirabait took %0.2f minutes" % (run_time / 60.0)) |
103 | 106 |
104 if return_code: | 107 if return_code: |
105 sys.stderr.write(stdout) | 108 sys.stderr.write(stdout) |
106 sys_exit("Return error code %i from command:\n%s" % (return_code, cmd), | 109 sys.exit("Return error code %i from command:\n%s" % (return_code, cmd), |
107 return_code) | 110 return_code) |
108 | 111 |
109 #Capture output | 112 # Capture output |
110 out_tmp = out_file_stem + "." + format | 113 out_tmp = out_file_stem + "." + format |
111 if not os.path.isfile(out_tmp): | 114 if not os.path.isfile(out_tmp): |
112 sys.stderr.write(stdout) | 115 sys.stderr.write(stdout) |
113 sys_exit("Missing output file from mirabait: %s" % out_tmp) | 116 sys.exit("Missing output file from mirabait: %s" % out_tmp) |
114 shutil.move(out_tmp, out_file) | 117 shutil.move(out_tmp, out_file) |