PageRenderTime 54ms CodeModel.GetById 43ms app.highlight 8ms RepoModel.GetById 2ms app.codeStats 0ms

/tools/metag_tools/megablast_xml_parser.py

https://bitbucket.org/cistrome/cistrome-harvard/
Python | 78 lines | 60 code | 8 blank | 10 comment | 13 complexity | ce64d58eceab6a8d45762afabb901377 MD5 | raw file
 1#!/usr/bin/env python
 2    
 3import sys, os, re
 4
 5if sys.version_info[:2] >= ( 2, 5 ):
 6    import xml.etree.cElementTree as ElementTree
 7else:
 8    from galaxy import eggs
 9    import pkg_resources; pkg_resources.require( "elementtree" )
10    from elementtree import ElementTree
11
12def stop_err( msg ):
13    sys.stderr.write( "%s\n" % msg )
14    sys.exit()
15
16def __main__():
17    source  = sys.argv[1]
18    hspTags = [
19           "Hsp_bit-score",
20           "Hsp_evalue",
21           "Hsp_query-from",
22           "Hsp_query-to",
23           "Hsp_hit-from",
24           "Hsp_hit-to",
25           "Hsp_query-frame",
26           "Hsp_hit-frame",
27           "Hsp_identity",
28           "Hsp_align-len",
29           "Hsp_qseq",
30           "Hsp_hseq",
31           "Hsp_midline"
32          ]
33    hspData = []
34
35    # get an iterable
36    try: 
37        context = ElementTree.iterparse( source, events=( "start", "end" ) )
38    except:
39        stop_err( "Invalid data format." )
40    # turn it into an iterator
41    context = iter( context )
42    # get the root element
43    try:
44        event, root = context.next()
45    except:
46        stop_err( "Invalid data format." )
47
48    outfile = open( sys.argv[2], 'w' )
49    try:
50        for event, elem in context:
51           # for every <Iteration> tag
52           if event == "end" and elem.tag == "Iteration":
53               query = elem.findtext( "Iteration_query-def" )
54               qLen = elem.findtext( "Iteration_query-len" )
55               # for every <Hit> within <Iteration>
56               for hit in elem.findall( "Iteration_hits/Hit" ):
57                   subject = hit.findtext( "Hit_id" )
58                   if re.search( '^gi', subject ):
59                       subject = subject.split('|')[1]
60                   sLen = hit.findtext( "Hit_len" )
61                   # for every <Hsp> within <Hit>
62                   for hsp in hit.findall( "Hit_hsps/Hsp" ):
63                        outfile.write( "%s\t%s\t%s\t%s" % ( query, qLen, subject, sLen ) )
64                        for tag in hspTags:
65                            outfile.write("\t%s" %(hsp.findtext( tag )))
66                            #hspData.append( hsp.findtext( tag ) )
67                        #hspData = []
68                        outfile.write('\n')
69               # prevents ElementTree from growing large datastructure
70               root.clear()
71               elem.clear()
72    except:
73        outfile.close()
74        stop_err( "The input data is malformed, or there is more than one dataset in the input file. Error: %s" % sys.exc_info()[1] )
75
76    outfile.close()
77
78if __name__ == "__main__": __main__()