30 lines
750 B
Python
30 lines
750 B
Python
|
#!/usr/bin/env python3
|
||
|
|
# -*- coding: UTF-8 -*-
|
||
|
|
import sys
|
||
|
|
|
||
|
|
import pysam
|
||
|
|
|
||
|
|
# 提取fr+rr>=2的点
|
||
|
|
if len(sys.argv) != 3:
|
||
|
|
print(" ".join(['usage:python', sys.argv[0], 'outputDir', 'tumor']))
|
||
|
|
sys.exit()
|
||
|
|
|
||
|
|
# output_dir = sys.argv[1]
|
||
|
|
# tumor = sys.argv[2]
|
||
|
|
|
||
|
|
# infile = "".join([output_dir, '/alignment/', tumor, '.rmdup.bam'])
|
||
|
|
# outfile = "".join([output_dir, '/alignment/', tumor, '.2r.rmdup.bam'])
|
||
|
|
samfile = pysam.AlignmentFile(sys.argv[1], "rb")
|
||
|
|
bamfile = pysam.AlignmentFile(sys.argv[2], "wb", template=samfile)
|
||
|
|
for s in samfile:
|
||
|
|
fr = 0
|
||
|
|
rr = 0
|
||
|
|
if s.has_tag('FR'):
|
||
|
|
fr = s.get_tag('FR')
|
||
|
|
if s.has_tag('RR'):
|
||
|
|
rr = s.get_tag('RR')
|
||
|
|
if fr + rr >= 2:
|
||
|
|
bamfile.write(s)
|
||
|
|
|
||
|
|
# os.system('samtools index ' + outfile)
|