Mercurial > repos > devteam > sam2interval

comparison sam2interval.py @ 1:75557c0908a9 draft default tip

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
"planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/sam2interval commit 206cd8245e7619b0e924c5066d0172129222993d"
author devteam
date Wed, 05 Feb 2020 06:58:53 -0500
parents 8c737b8ddc45
children
comparison
equal deleted inserted replaced
0:8c737b8ddc45 1:75557c0908a9
2 2
3 import sys 3 import sys
4 import optparse 4 import optparse
5 import re 5 import re
6 6
7 def stop_err( msg ): 7
8 sys.stderr.write( msg ) 8 def stop_err(msg):
9 sys.exit() 9 sys.exit(msg)
10
10 11
11 def main(): 12 def main():
12 usage = """%prog [options] 13 usage = """%prog [options]
13 14
14 options (listed below) default to 'None' if omitted 15 options (listed below) default to 'None' if omitted
15 """ 16 """
16 parser = optparse.OptionParser(usage=usage) 17 parser = optparse.OptionParser(usage=usage)
17 18
18 parser.add_option( 19 parser.add_option(
19 '-f','--input_sam_file', 20 '-f', '--input_sam_file',
20 metavar="INPUT_SAM_FILE", 21 metavar="INPUT_SAM_FILE",
21 dest='input_sam', 22 dest='input_sam',
22 default = False, 23 default=False,
23 help='Name of the SAM file to be filtered. STDIN is default') 24 help='Name of the SAM file to be filtered. STDIN is default')
24 25
25 parser.add_option( 26 parser.add_option(
26 '-c','--flag_column', 27 '-c', '--flag_column',
27 dest='flag_col', 28 dest='flag_col',
28 default = '2', 29 default='2',
29 help='Column containing SAM bitwise flag. 1-based') 30 help='Column containing SAM bitwise flag. 1-based')
30 31
31 parser.add_option( 32 parser.add_option(
32 '-s','--start_column', 33 '-s', '--start_column',
33 dest='start_col', 34 dest='start_col',
34 default = '4', 35 default='4',
35 help='Column containing position. 1-based') 36 help='Column containing position. 1-based')
36 37
37 parser.add_option( 38 parser.add_option(
38 '-g','--cigar_column', 39 '-g', '--cigar_column',
39 dest='cigar_col', 40 dest='cigar_col',
40 default = '6', 41 default='6',
41 help='Column containing CIGAR or extended CIGAR string') 42 help='Column containing CIGAR or extended CIGAR string')
42 43
43 parser.add_option( 44 parser.add_option(
44 '-r','--ref_column', 45 '-r', '--ref_column',
45 dest='ref_col', 46 dest='ref_col',
46 default = '3', 47 default='3',
47 help='Column containing name of the reference sequence coordinate. 1-based') 48 help='Column containing name of the reference sequence coordinate. 1-based')
48 49
49 parser.add_option( 50 parser.add_option(
50 '-e','--read_column', 51 '-e', '--read_column',
51 dest='read_col', 52 dest='read_col',
52 default = '1', 53 default='1',
53 help='Column containing read name. 1-based') 54 help='Column containing read name. 1-based')
54 55
55 parser.add_option( 56 parser.add_option(
56 '-p','--print_all', 57 '-p', '--print_all',
57 dest='prt_all', 58 dest='prt_all',
58 action='store_true', 59 action='store_true',
59 default = False, 60 default=False,
60 help='Print coordinates and original SAM?') 61 help='Print coordinates and original SAM?')
61 62
62 options, args = parser.parse_args() 63 options, args = parser.parse_args()
63 64
64 if options.input_sam: 65 if options.input_sam:
65 infile = open ( options.input_sam, 'r') 66 infile = open(options.input_sam, 'r')
66 else: 67 else:
67 infile = sys.stdin 68 infile = sys.stdin
68 69
69 cigar = re.compile( '\d+M|\d+N|\d+D|\d+P' ) 70 cigar = re.compile('\d+M|\d+N|\d+D|\d+P')
70 71
71 print '#chrom\tstart\tend\tstrand\tread_name' # provide a (partial) header so that strand is automatically set in metadata 72 print('#chrom\tstart\tend\tstrand\tread_name') # provide a (partial) header so that strand is automatically set in metadata
72 73
73 for line in infile: 74 for line in infile:
74 line = line.rstrip( '\r\n' ) 75 line = line.rstrip('\r\n')
75 if line and not line.startswith( '#' ) and not line.startswith( '@' ) : 76 if line and not line.startswith('#') and not line.startswith('@'):
76 fields = line.split( '\t' ) 77 fields = line.split('\t')
77 start = int( fields[ int( options.start_col ) - 1 ] ) - 1 78 start = int(fields[int(options.start_col) - 1]) - 1
78 end = 0 79 end = 0
79 for op in cigar.findall( fields[ int( options.cigar_col) - 1 ] ): 80 for op in cigar.findall(fields[int(options.cigar_col) - 1]):
80 end += int( op[ 0:len( op ) - 1 ] ) 81 end += int(op[0:len(op) - 1])
81 82
82 strand = '+' 83 strand = '+'
83 if bool( int( fields[ int( options.flag_col ) - 1 ] ) & 0x0010 ): 84 if bool(int(fields[int(options.flag_col) - 1]) & 0x0010):
84 strand = '-' 85 strand = '-'
85 read_name = fields[ int( options.read_col ) - 1 ] 86 read_name = fields[int(options.read_col) - 1]
86 ref_name = fields[ int( options.ref_col ) - 1 ] 87 ref_name = fields[int(options.ref_col) - 1]
87 88
88 if ref_name != '*': 89 if ref_name != '*':
89 # Do not print lines with unmapped reads that contain '*' instead of chromosome name 90 # Do not print lines with unmapped reads that contain '*' instead of chromosome name
90 if options.prt_all: 91 if options.prt_all:
91 print '%s\t%s\t%s\t%s\t%s' % (ref_name, str(start), str(end+start), strand, line) 92 print('%s\t%s\t%s\t%s\t%s' % (ref_name, str(start), str(end + start), strand, line))
92 else: 93 else:
93 print '%s\t%s\t%s\t%s\t%s' % (ref_name, str(start), str(end+start), strand, read_name) 94 print('%s\t%s\t%s\t%s\t%s' % (ref_name, str(start), str(end + start), strand, read_name))
94 95
95 if __name__ == "__main__": main()
96 96
97 if __name__ == "__main__":
98 main()