root/galaxy-central/eggs/bx_python-0.5.0_dev_f74aec067563-py2.6-macosx-10.6-universal-ucs2.egg/EGG-INFO/scripts/maf_to_axt.py

#!/usr/bin/python2.6
"""
Application to convert MAF file to AXT file, projecting to any two species. 
Reads a MAF file from standard input and writes an AXT file to standard out;  
some statistics are written to standard error.  The user must specify the 
two species of interest.

usage: %prog primary_species secondary_species < maf_file > axt_file
"""

__author__ = "Bob Harris (rsharris@bx.psu.edu)"

import sys
import copy
import bx.align.maf
import bx.align.axt

def usage(s=None):
        message = """
maf_to_axt primary_species secondary_species < maf_file > axt_file
"""
        if (s == None): sys.exit (message)
        else:           sys.exit ("%s\n%s" % (s,message))


def main():

        # parse the command line

        primary   = None
        secondary = None

        args = sys.argv[1:]
        while (len(args) > 0):
                arg = args.pop(0)
                val = None
                fields = arg.split("=",1)
                if (len(fields) == 2):
                        arg = fields[0]
                        val = fields[1]
                        if (val == ""):
                                usage("missing a value in %s=" % arg)

                if (primary == None) and (val == None):
                        primary = arg
                elif (secondary == None) and (val == None):
                        secondary = arg
                else:
                        usage("unknown argument: %s" % arg)

        if (primary == None):
                usage("missing primary species")

        if (secondary == None):
                usage("missing secondary species")

        # read the alignments and other info

        out = bx.align.axt.Writer(sys.stdout)

        axtsRead = 0
        mafsWritten = 0
        for mafBlock in bx.align.maf.Reader(sys.stdin):
                axtsRead += 1

                p = mafBlock.get_component_by_src_start(primary)
                if (p == None): continue
                s = mafBlock.get_component_by_src_start(secondary)
                if (s == None): continue

                axtBlock = bx.align.Alignment (mafBlock.score, mafBlock.attributes)
                axtBlock.add_component (clone_component(p))
                axtBlock.add_component (clone_component(s))

                remove_mutual_gaps (axtBlock)
                if (axtBlock.text_size == 0):
                        continue

                out.write (axtBlock)
                mafsWritten += 1

        sys.stderr.write ("%d blocks read, %d written\n" % (axtsRead,mafsWritten))


def clone_component(c):
        return bx.align.Component (c.src, c.start, c.size, c.strand, c.src_size, \
                                   copy.copy(c.text))


def remove_mutual_gaps (block):

        if (len(block.components) == 0): return

        nonGaps = []

        for c in block.components:
                for ix in range(0,block.text_size):
                        if (ix not in nonGaps) and (c.text[ix] != "-"):
                                nonGaps.append(ix)

        nonGaps.sort()

        for c in block.components:
                c.text = "".join([c.text[ix] for ix in nonGaps])

        block.text_size = len(nonGaps)


if __name__ == "__main__": main()