<s> '''<STR_LIT>''' <EOL> import sys <EOL> import os <EOL> import pysam <EOL> from ngsutils . bam import bam_iter <EOL> def bam_removeclipping ( infile , outfile ) : <EOL> bam = pysam . Samfile ( infile , "<STR_LIT:rb>" ) <EOL> out = pysam . Samfile ( outfile , "<STR_LIT:wb>" , template = bam ) <EOL> total = <NUM_LIT:0> <EOL> count = <NUM_LIT:0> <EOL> unmapped = <NUM_LIT:0> <EOL> for read in bam_iter ( bam ) : <EOL> code = read_removeclipping ( read ) <EOL> if code == <NUM_LIT:1> : <EOL> unmapped += <NUM_LIT:1> <EOL> elif code == <NUM_LIT:2> : <EOL> count += <NUM_LIT:1> <EOL> total += <NUM_LIT:1> <EOL> out . write ( read ) <EOL> bam . close ( ) <EOL> out . close ( ) <EOL> sys . stderr . write ( '<STR_LIT>' % ( total , count , unmapped ) ) <EOL> def read_removeclipping ( read ) : <EOL> if read . is_unmapped : <EOL> return <NUM_LIT:1> <EOL> newcigar = [ ] <EOL> clip_5 = <NUM_LIT:0> <EOL> clip_3 = <NUM_LIT:0> <EOL> changed = False <EOL> inseq = False <EOL> for op , length in read . cigar : <EOL> if op == <NUM_LIT:5> : <EOL> changed = True <EOL> elif op == <NUM_LIT:4> : <EOL> changed = True <EOL> if not inseq : <EOL> clip_5 = length <EOL> else : <EOL> clip_3 = length <EOL> else : <EOL> inseq = True <EOL> newcigar . append ( ( op , length ) ) <EOL> if not changed : <EOL> return <NUM_LIT:0> <EOL> read . cigar = newcigar <EOL> orig_length = len ( read . seq ) <EOL> s = read . seq <EOL> q = read . qual <EOL> if clip_3 : <EOL> read . seq = s [ clip_5 : - clip_3 ] <EOL> if q : <EOL> read . qual = q [ clip_5 : - clip_3 ] <EOL> else : <EOL> read . seq = s [ clip_5 : ] <EOL> if q : <EOL> read . qual = q [ clip_5 : ] <EOL> newtags = [ ] <EOL> if clip_5 : <EOL> newtags . append ( ( '<STR_LIT>' , clip_5 ) ) <EOL> if clip_3 : <EOL> newtags . append ( ( '<STR_LIT>' , clip_3 ) ) <EOL> newtags . append ( ( '<STR_LIT>' , float ( clip_5 + clip_3 ) / orig_length ) ) <EOL> read . tags = read . tags + newtags <EOL> return <NUM_LIT:2> <EOL> def usage ( ) : <EOL> print __doc__ <EOL> print """<STR_LIT>""" <EOL> sys . exit ( - <NUM_LIT:1> ) <EOL> if __name__ == "<STR_LIT:__main__>" : <EOL> infile = None <EOL> outfile = None <EOL> force = False <EOL> for arg in sys . argv [ <NUM_LIT:1> : ] : <EOL> if arg == "<STR_LIT>" : <EOL> usage ( ) <EOL> elif arg == "<STR_LIT>" : <EOL> force = True <EOL> elif not infile : <EOL> if os . path . exists ( os . path . expanduser ( arg ) ) : <EOL> infile = os . path . expanduser ( arg ) <EOL> else : <EOL> sys . stderr . write ( "<STR_LIT>" % arg ) <EOL> usage ( ) <EOL> elif not outfile : <EOL> if force or not os . path . exists ( os . path . expanduser ( arg ) ) : <EOL> outfile = arg <EOL> else : <EOL> sys . stderr . write ( <EOL> "<STR_LIT>" % arg <EOL> ) <EOL> usage ( ) <EOL> else : <EOL> usage ( ) <EOL> if not infile or not outfile : <EOL> usage ( ) <EOL> bam_removeclipping ( infile , outfile ) </s>
