Mercurial > repos > peterjc > mira4_assembler

comparison tools/mira4_0/mira4_convert.py @ 31:fd95aaef8818 draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
planemo upload for repository https://github.com/peterjc/galaxy_mira/tree/master/tools/mira4_0 commit bc3d484c5cd68ddcf456db2fff489d584aa2034c
author peterjc
date Wed, 10 Feb 2016 09:07:39 -0500
parents 55ae131c5862
children 56b421d59805
comparison
equal deleted inserted replaced
30:b506e3b779fa 31:fd95aaef8818
5 """ 5 """
6 import os 6 import os
7 import sys 7 import sys
8 import subprocess 8 import subprocess
9 import shutil 9 import shutil
10 import time
11 import tempfile
12 from optparse import OptionParser 10 from optparse import OptionParser
13 try: 11 try:
14 from io import BytesIO 12 from io import BytesIO
15 except ImportError: 13 except ImportError:
16 #Should we worry about Python 2.5 or older? 14 #Should we worry about Python 2.5 or older?
19 #Do we need any PYTHONPATH magic? 17 #Do we need any PYTHONPATH magic?
20 from mira4_make_bam import depad 18 from mira4_make_bam import depad
21 19
22 WRAPPER_VER = "0.0.7" # Keep in sync with the XML file 20 WRAPPER_VER = "0.0.7" # Keep in sync with the XML file
23 21
24 def sys_exit(msg, err=1):
25 sys.stderr.write(msg+"\n")
26 sys.exit(err)
27 22
28 def run(cmd): 23 def run(cmd):
29 #Avoid using shell=True when we call subprocess to ensure if the Python 24 #Avoid using shell=True when we call subprocess to ensure if the Python
30 #script is killed, so too is the child process. 25 #script is killed, so too is the child process.
31 try: 26 try:
32 child = subprocess.Popen(cmd, stdout=subprocess.PIPE, stderr=subprocess.PIPE) 27 child = subprocess.Popen(cmd, stdout=subprocess.PIPE, stderr=subprocess.PIPE)
33 except Exception, err: 28 except Exception, err:
34 sys_exit("Error invoking command:\n%s\n\n%s\n" % (" ".join(cmd), err)) 29 sys.exit("Error invoking command:\n%s\n\n%s\n" % (" ".join(cmd), err))
35 #Use .communicate as can get deadlocks with .wait(), 30 #Use .communicate as can get deadlocks with .wait(),
36 stdout, stderr = child.communicate() 31 stdout, stderr = child.communicate()
37 return_code = child.returncode 32 return_code = child.returncode
38 if return_code: 33 if return_code:
39 cmd_str = " ".join(cmd) # doesn't quote spaces etc 34 cmd_str = " ".join(cmd) # doesn't quote spaces etc
40 if stderr and stdout: 35 if stderr and stdout:
41 sys_exit("Return code %i from command:\n%s\n\n%s\n\n%s" % (return_code, cmd_str, stdout, stderr)) 36 sys.exit("Return code %i from command:\n%s\n\n%s\n\n%s" % (return_code, cmd_str, stdout, stderr))
42 else: 37 else:
43 sys_exit("Return code %i from command:\n%s\n%s" % (return_code, cmd_str, stderr)) 38 sys.exit("Return code %i from command:\n%s\n%s" % (return_code, cmd_str, stderr))
44 39
45 def get_version(mira_binary): 40 def get_version(mira_binary):
46 """Run MIRA to find its version number""" 41 """Run MIRA to find its version number"""
47 # At the commend line I would use: mira -v | head -n 1 42 # At the commend line I would use: mira -v | head -n 1
48 # however there is some pipe error when doing that here. 43 # however there is some pipe error when doing that here.
96 parser.add_option("-v", "--version", dest="version", 91 parser.add_option("-v", "--version", dest="version",
97 default=False, action="store_true", 92 default=False, action="store_true",
98 help="Show version and quit") 93 help="Show version and quit")
99 options, args = parser.parse_args() 94 options, args = parser.parse_args()
100 if args: 95 if args:
101 sys_exit("Expected options (e.g. --input example.maf), not arguments") 96 sys.exit("Expected options (e.g. --input example.maf), not arguments")
102 97
103 input_maf = options.input 98 input_maf = options.input
104 out_maf = options.maf 99 out_maf = options.maf
105 out_bam = options.bam 100 out_bam = options.bam
106 out_fasta = options.fasta 101 out_fasta = options.fasta
108 out_cstats = options.cstats 103 out_cstats = options.cstats
109 104
110 try: 105 try:
111 mira_path = os.environ["MIRA4"] 106 mira_path = os.environ["MIRA4"]
112 except KeyError: 107 except KeyError:
113 sys_exit("Environment variable $MIRA4 not set") 108 sys.exit("Environment variable $MIRA4 not set")
114 mira_convert = os.path.join(mira_path, "miraconvert") 109 mira_convert = os.path.join(mira_path, "miraconvert")
115 if not os.path.isfile(mira_convert): 110 if not os.path.isfile(mira_convert):
116 sys_exit("Missing miraconvert under $MIRA4, %r\nFolder contained: %s" 111 sys.exit("Missing miraconvert under $MIRA4, %r\nFolder contained: %s"
117 % (mira_convert, ", ".join(os.listdir(mira_path)))) 112 % (mira_convert, ", ".join(os.listdir(mira_path))))
118 113
119 mira_convert_ver = get_version(mira_convert) 114 mira_convert_ver = get_version(mira_convert)
120 if not mira_convert_ver.strip().startswith("4.0"): 115 if not mira_convert_ver.strip().startswith("4.0"):
121 sys_exit("This wrapper is for MIRA V4.0, not:\n%s\n%s" % (mira_convert_ver, mira_convert)) 116 sys.exit("This wrapper is for MIRA V4.0, not:\n%s\n%s" % (mira_convert_ver, mira_convert))
122 if options.version: 117 if options.version:
123 print("%s, MIRA wrapper version %s" % (mira_convert_ver, WRAPPER_VER)) 118 print("%s, MIRA wrapper version %s" % (mira_convert_ver, WRAPPER_VER))
124 sys.exit(0) 119 sys.exit(0)
125 120
126 if not input_maf: 121 if not input_maf:
127 sys_exit("Input MIRA file is required") 122 sys.exit("Input MIRA file is required")
128 elif not os.path.isfile(input_maf): 123 elif not os.path.isfile(input_maf):
129 sys_exit("Missing input MIRA file: %r" % input_maf) 124 sys.exit("Missing input MIRA file: %r" % input_maf)
130 125
131 if not (out_maf or out_bam or out_fasta or out_ace or out_cstats): 126 if not (out_maf or out_bam or out_fasta or out_ace or out_cstats):
132 sys_exit("No output requested") 127 sys.exit("No output requested")
133 128
134 129
135 def check_min_int(value, name): 130 def check_min_int(value, name):
136 try: 131 try:
137 i = int(value) 132 i = int(value)
138 except: 133 except ValueError:
139 sys_exit("Bad %s setting, %r" % (name, value)) 134 sys.exit("Bad %s setting, %r" % (name, value))
140 if i < 0: 135 if i < 0:
141 sys_exit("Negative %s setting, %r" % (name, value)) 136 sys.exit("Negative %s setting, %r" % (name, value))
142 return i 137 return i
143 138
144 min_length = check_min_int(options.min_length, "minimum length") 139 min_length = check_min_int(options.min_length, "minimum length")
145 min_cover = check_min_int(options.min_cover, "minimum cover") 140 min_cover = check_min_int(options.min_cover, "minimum cover")
146 min_reads = check_min_int(options.min_reads, "minimum reads") 141 min_reads = check_min_int(options.min_reads, "minimum reads")
174 cmd_list.append("cstats") 169 cmd_list.append("cstats")
175 run(cmd_list) 170 run(cmd_list)
176 171
177 def collect(old, new): 172 def collect(old, new):
178 if not os.path.isfile(old): 173 if not os.path.isfile(old):
179 sys_exit("Missing expected output file %s" % old) 174 sys.exit("Missing expected output file %s" % old)
180 shutil.move(old, new) 175 shutil.move(old, new)
181 176
182 if out_maf: 177 if out_maf:
183 collect(os.path.join(temp, "converted.maf"), out_maf) 178 collect(os.path.join(temp, "converted.maf"), out_maf)
184 if out_fasta: 179 if out_fasta:
188 collect(old, out_fasta) 183 collect(old, out_fasta)
189 else: 184 else:
190 #Might the output be filtered down to zero contigs? 185 #Might the output be filtered down to zero contigs?
191 old = os.path.join(temp, "converted.fasta") 186 old = os.path.join(temp, "converted.fasta")
192 if not os.path.isfile(old): 187 if not os.path.isfile(old):
193 sys_exit("Missing expected output FASTA file") 188 sys.exit("Missing expected output FASTA file")
194 elif os.path.getsize(old) == 0: 189 elif os.path.getsize(old) == 0:
195 print("Warning - no contigs (harsh filters?)") 190 print("Warning - no contigs (harsh filters?)")
196 collect(old, out_fasta) 191 collect(old, out_fasta)
197 else: 192 else:
198 sys_exit("Missing expected output FASTA file (only generic file present)") 193 sys.exit("Missing expected output FASTA file (only generic file present)")
199 if out_ace: 194 if out_ace:
200 collect(os.path.join(temp, "converted.maf"), out_ace) 195 collect(os.path.join(temp, "converted.maf"), out_ace)
201 if out_cstats: 196 if out_cstats:
202 collect(os.path.join(temp, "converted_info_contigstats.txt"), out_cstats) 197 collect(os.path.join(temp, "converted_info_contigstats.txt"), out_cstats)
203 198
205 assert os.path.isfile(out_fasta) 200 assert os.path.isfile(out_fasta)
206 old = os.path.join(temp, "converted.samnbb") 201 old = os.path.join(temp, "converted.samnbb")
207 if not os.path.isfile(old): 202 if not os.path.isfile(old):
208 old = os.path.join(temp, "converted.sam") 203 old = os.path.join(temp, "converted.sam")
209 if not os.path.isfile(old): 204 if not os.path.isfile(old):
210 sys_exit("Missing expected intermediate file %s" % old) 205 sys.exit("Missing expected intermediate file %s" % old)
211 h = BytesIO() 206 h = BytesIO()
212 msg = depad(out_fasta, old, out_bam, h) 207 msg = depad(out_fasta, old, out_bam, h)
213 if msg: 208 if msg:
214 print(msg) 209 print(msg)
215 print(h.getvalue()) 210 print(h.getvalue())