Source

bx-python / scripts / maf_limit_to_species.py

#!/usr/bin/env python2.3

"""
Read a maf file from stdin and write out a new maf with only blocks having all
of the required in species, after dropping any other species and removing
columns containing only gaps.

usage: %prog species,species2,... < maf
"""

import psyco_full

import bx.align.maf
import copy
import sys

from itertools import *

def main():

    species = sys.argv[1].split( ',' )

    maf_reader = bx.align.maf.Reader( sys.stdin )
    maf_writer = bx.align.maf.Writer( sys.stdout )

    for m in maf_reader:        
        new_components = []    
        for comp in m.components:
            if comp.src.split( '.' )[0] in species:
                new_components.append( comp )
        m.components = new_components
        m.remove_all_gap_columns()
        if len( m.components ) > 1:
            maf_writer.write( m )
        
    maf_reader.close()
    maf_writer.close()

if __name__ == "__main__": 
    main()