/tools/samtools/bam_to_sam.py
https://bitbucket.org/cistrome/cistrome-harvard/ · Python · 129 lines · 98 code · 9 blank · 22 comment · 31 complexity · 8bef67ec54e69d1a42ac689464c63890 MD5 · raw file
- #!/usr/bin/env python
- """
- Converts BAM data to sorted SAM data.
- usage: bam_to_sam.py [options]
- --input1: SAM file to be converted
- --output1: output dataset in bam format
- """
- import optparse, os, sys, subprocess, tempfile, shutil
- from galaxy import eggs
- import pkg_resources; pkg_resources.require( "bx-python" )
- from bx.cookbook import doc_optparse
- #from galaxy import util
- def stop_err( msg ):
- sys.stderr.write( '%s\n' % msg )
- sys.exit()
- def __main__():
- #Parse Command Line
- parser = optparse.OptionParser()
- parser.add_option( '', '--input1', dest='input1', help='The input SAM dataset' )
- parser.add_option( '', '--output1', dest='output1', help='The output BAM dataset' )
- parser.add_option( '', '--header', dest='header', action='store_true', default=False, help='Write SAM Header' )
- ( options, args ) = parser.parse_args()
- # output version # of tool
- try:
- tmp = tempfile.NamedTemporaryFile().name
- tmp_stdout = open( tmp, 'wb' )
- proc = subprocess.Popen( args='samtools 2>&1', shell=True, stdout=tmp_stdout )
- tmp_stdout.close()
- returncode = proc.wait()
- stdout = None
- for line in open( tmp_stdout.name, 'rb' ):
- if line.lower().find( 'version' ) >= 0:
- stdout = line.strip()
- break
- if stdout:
- sys.stdout.write( 'Samtools %s\n' % stdout )
- else:
- raise Exception
- except:
- sys.stdout.write( 'Could not determine Samtools version\n' )
- tmp_dir = tempfile.mkdtemp()
- try:
- # exit if input file empty
- if os.path.getsize( options.input1 ) == 0:
- raise Exception, 'Initial BAM file empty'
- # Sort alignments by leftmost coordinates. File <out.prefix>.bam will be created. This command
- # may also create temporary files <out.prefix>.%d.bam when the whole alignment cannot be fitted
- # into memory ( controlled by option -m ).
- tmp_sorted_aligns_file = tempfile.NamedTemporaryFile( dir=tmp_dir )
- tmp_sorted_aligns_file_base = tmp_sorted_aligns_file.name
- tmp_sorted_aligns_file_name = '%s.bam' % tmp_sorted_aligns_file.name
- tmp_sorted_aligns_file.close()
- command = 'samtools sort %s %s' % ( options.input1, tmp_sorted_aligns_file_base )
- tmp = tempfile.NamedTemporaryFile( dir=tmp_dir ).name
- tmp_stderr = open( tmp, 'wb' )
- proc = subprocess.Popen( args=command, shell=True, cwd=tmp_dir, stderr=tmp_stderr.fileno() )
- returncode = proc.wait()
- tmp_stderr.close()
- # get stderr, allowing for case where it's very large
- tmp_stderr = open( tmp, 'rb' )
- stderr = ''
- buffsize = 1048576
- try:
- while True:
- stderr += tmp_stderr.read( buffsize )
- if not stderr or len( stderr ) % buffsize != 0:
- break
- except OverflowError:
- pass
- tmp_stderr.close()
- if returncode != 0:
- raise Exception, stderr
- # exit if sorted BAM file empty
- if os.path.getsize( tmp_sorted_aligns_file_name) == 0:
- raise Exception, 'Intermediate sorted BAM file empty'
- except Exception, e:
- #clean up temp files
- if os.path.exists( tmp_dir ):
- shutil.rmtree( tmp_dir )
- stop_err( 'Error sorting alignments from (%s), %s' % ( options.input1, str( e ) ) )
- try:
- # Extract all alignments from the input BAM file to SAM format ( since no region is specified, all the alignments will be extracted ).
- if options.header:
- view_options = "-h"
- else:
- view_options = ""
- command = 'samtools view %s -o %s %s' % ( view_options, options.output1, tmp_sorted_aligns_file_name )
- tmp = tempfile.NamedTemporaryFile( dir=tmp_dir ).name
- tmp_stderr = open( tmp, 'wb' )
- proc = subprocess.Popen( args=command, shell=True, cwd=tmp_dir, stderr=tmp_stderr.fileno() )
- returncode = proc.wait()
- tmp_stderr.close()
- # get stderr, allowing for case where it's very large
- tmp_stderr = open( tmp, 'rb' )
- stderr = ''
- buffsize = 1048576
- try:
- while True:
- stderr += tmp_stderr.read( buffsize )
- if not stderr or len( stderr ) % buffsize != 0:
- break
- except OverflowError:
- pass
- tmp_stderr.close()
- if returncode != 0:
- raise Exception, stderr
- except Exception, e:
- #clean up temp files
- if os.path.exists( tmp_dir ):
- shutil.rmtree( tmp_dir )
- stop_err( 'Error extracting alignments from (%s), %s' % ( options.input1, str( e ) ) )
- #clean up temp files
- if os.path.exists( tmp_dir ):
- shutil.rmtree( tmp_dir )
- # check that there are results in the output file
- if os.path.getsize( options.output1 ) > 0:
- sys.stdout.write( 'BAM file converted to SAM' )
- else:
- stop_err( 'The output file is empty, there may be an error with your input file.' )
- if __name__=="__main__": __main__()