diff SNV/annotate.py @ 0:74f5ea818cea

Uploaded

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/SNV/annotate.py	Wed Oct 12 19:50:38 2011 -0400
@@ -0,0 +1,80 @@
+#!/usr/bin/env python
+
+"""
+Annotates filtered SNVMix output that has been attached to codon information (outputs the same data with additional columns describing the predicted effect of the SNV)
+
+usage: %prog [options]
+   -i, --input1=i: bam file
+   -o, --output1=o: Output pileup
+"""
+
+import os, shutil, subprocess, sys, tempfile
+from galaxy import eggs
+import pkg_resources; pkg_resources.require( "bx-python" )
+from bx.cookbook import doc_optparse
+
+
+
+def stop_err( msg ):
+    sys.stderr.write( '%s\n' % msg )
+    sys.exit()
+
+def check_seq_file( dbkey, GALAXY_DATA_INDEX_DIR ):
+    seqFile = '%s/sam_fa_indices.loc' % GALAXY_DATA_INDEX_DIR
+    seqPath = ''
+    for line in open( seqFile ):
+        line = line.rstrip( '\r\n' )
+        if line and not line.startswith( '#' ) and line.startswith( 'index' ):
+            fields = line.split( '\t' )
+            if len( fields ) < 3:
+                continue
+            if fields[1] == dbkey:
+                seqPath = fields[2].strip()
+                break
+    return seqPath
+
+def __main__():
+    #Parse Command Line
+    options, args = doc_optparse.parse( __doc__ )
+    #make temp dir
+    tmpDir = tempfile.mkdtemp()
+
+    cmd = 'identify_nonsynonymous_mutations.pl < %s > %s'
+    try:
+        # have to nest try-except in try-finally to handle 2.4
+        try:
+            cmd = cmd % ( options.input1, options.output1 )
+            print(cmd, "\n")
+            tmp = tempfile.NamedTemporaryFile( dir=tmpDir ).name
+            tmp_stderr = open( tmp, 'wb' )
+
+            proc = subprocess.Popen( args=cmd, shell=True, cwd=tmpDir, stderr=tmp_stderr.fileno() )
+            returncode = proc.wait()
+            tmp_stderr.close()
+            #did it succeed?
+            # get stderr, allowing for case where it's very large
+            tmp_stderr = open( tmp, 'rb' )
+            stderr = ''
+            try:
+                while True:
+                    stderr += tmp_stderr.read( )
+                    if not stderr:
+                        break
+            except OverflowError:
+                pass
+            tmp_stderr.close()
+            if returncode != 0:
+                raise Exception, stderr
+        except Exception, e:
+            stop_err( 'Error running annotator tool\n' + str( e ) )
+    finally:
+        #clean up temp files
+        if os.path.exists( tmpDir ):
+            shutil.rmtree( tmpDir )
+    # check that there are results in the output file
+    if os.path.getsize( options.output1 ) > 0:
+        sys.stdout.write( 'wrote annotated output' )
+    else:
+        stop_err( 'The output file is empty' )
+
+if __name__ == "__main__" : __main__()