def match_pairs(reads, out_fhand, orphan_out_fhand, out_format, ordered=True,
                check_order_buffer_size=0, max_reads_memory=None,
                temp_dir=None):
    '''It matches the seq pairs in an iterator and splits the orphan seqs.'''
    counts = 0
    check_order_buffer = KeyedSet()
    for pair in _get_paired_and_orphan(reads, ordered, max_reads_memory,
                                       temp_dir):
        if len(pair) == 1:
            write_seqs(pair, orphan_out_fhand, out_format)
            try:
                name = _parse_pair_direction_and_name(pair[0])[0]
            except PairDirectionError:
                name = get_name(pair[0])
            if ordered and counts < check_order_buffer_size:
                counts += 1
                if not check_order_buffer.check_add(name):
                    msg = 'Reads are not ordered by pairs.Use unordered option'
                    raise ItemsNotSortedError(msg)
            elif ordered and counts >= check_order_buffer_size:
                if name in check_order_buffer:
                    msg = 'Reads are not ordered by pairs.Use unordered option'
                    raise ItemsNotSortedError(msg)
        elif len(pair) == 2:
            write_seqs(pair, out_fhand, out_format)
    flush_fhand(orphan_out_fhand)
    flush_fhand(out_fhand) def interleave_pairs(seqs1, seqs2, skip_checks=False):
    '''A generator that interleaves the paired reads found in two iterators.

    It will fail if forward and reverse reads do not match in both sequence
    iterators.
    '''
    for seq1, seq2 in izip_longest(seqs1, seqs2, fillvalue=None):
        if not skip_checks:
            if seq1 is None or seq2 is None:
                msg = 'The files had a different number of sequences'
                raise InterleaveError(msg)
            _check_name_and_direction_match(seq1, seq2)
        if seq1 is not None:
            yield seq1
        if seq2 is not None:
            yield seq2 def deinterleave_pairs(seqs, out_fhand1, out_fhand2, out_format):
    '''It splits a sequence iterator with alternating paired reads in two.

    It will fail if forward and reverse reads are not alternating.
    '''
    for pair in group_pairs(seqs, n_seqs_in_pair=2):
        write_seqs((pair[0],), out_fhand1, out_format)
        write_seqs((pair[1],), out_fhand2, out_format)
    out_fhand1.flush()
    out_fhand2.flush() def group_pairs(seqs, n_seqs_in_pair=None, check_all_same_n_seqs=True,
                check_name_matches=True):

    seqs = iter(seqs)
    if n_seqs_in_pair is None:
        first_pair, next_read = _get_first_pair_by_name(seqs)
        if first_pair is None:
            n_seqs_in_pair = None
        else:
            yield first_pair
            n_seqs_in_pair = len(first_pair)
            seqs = chain([next_read], seqs)

    if n_seqs_in_pair == 1:
        # No need to check anything, a pair cannot have less than one read
        # or more than one name
        check_all_same_n_seqs = False
        check_name_matches = False

    if n_seqs_in_pair:
        pairs = group_in_packets_fill_last(seqs, packet_size=n_seqs_in_pair)
        for pair in pairs:
            pair = filter(lambda seq: seq is not None, pair)
            if check_all_same_n_seqs and n_seqs_in_pair != len(pair):
                msg = 'The last pair has fewer reads'
                raise InterleaveError(msg)
            if check_name_matches:
                _check_name_and_direction_match(*pair)
            yield pair def group_pairs_by_name(seqs, all_pairs_same_n_seqs=False):
    paired_seqs = []
    prev_name = None
    n_seqs_per_pair = None
    for seq in iter(seqs):
        try:
            name = _parse_pair_direction_and_name(seq)[0]
        except PairDirectionError:
            name = None
        if name is None or (paired_seqs and name != prev_name):
            if all_pairs_same_n_seqs:
                if n_seqs_per_pair is None:
                    n_seqs_per_pair = len(paired_seqs)
                elif n_seqs_per_pair != len(paired_seqs):
                    msg = 'Pair had different number of reads: '
                    msg += prev_name
                    raise InterleaveError(msg)
            if paired_seqs:
                yield paired_seqs
                paired_seqs = []
        paired_seqs.append(seq)
        prev_name = name
    if paired_seqs:
        if all_pairs_same_n_seqs and n_seqs_per_pair is not None:
            if n_seqs_per_pair != len(paired_seqs):
                msg = 'Pair had different number of reads: '
                msg += prev_name
                raise InterleaveError(msg)
        yield paired_seqs def _parse_pair_direction_and_name_from_title(title):
    'It guesses the direction from the title line'
    reg_exps = ['(.+)[/|\\\\](\d+)', '(.+)\s(\d+):.+', '(.+)\.(\w)\s?']
    for reg_exp in reg_exps:
        match = re.match(reg_exp, title)
        if match:
            name, direction = match.groups()
            if direction in ('1', 'f'):
                direction = FWD
            elif direction in ('2', 'r'):
                direction = REV
            else:
                raise PairDirectionError('unknown direction descriptor')
            return name, direction
    raise PairDirectionError('Unable to detect the direction of the seq') def _parse_pair_direction_and_name(seq):
    'It parses the description field to get the name and the pair direction'
    return _parse_pair_direction_and_name_from_title(get_title(seq)) def flat_zip_longest(iter1, iter2, fillvalue=None):
    '''It yields items alternatively from both iterators.

    It won't return the items equal to the fillvalue.
    '''
    for item1, item2 in izip_longest(iter1, iter2, fillvalue=fillvalue):
        if item1 != fillvalue:
            yield item1
        if item2 != fillvalue:
            yield item2 FWD = 'fwd' SEQRECORD = 'seqrecord' SEQITEM = 'seqitem' SEQ_BIN_DIR = os.path.join(BIN_DIR, 'seq') TEST_DATA_DIR = os.path.abspath(os.path.join(os.path.dirname(__file__), '..',
                                             '..', 'test', 'test_data')) def get_str_seq(seq):
    seq_class = seq.kind
    if seq_class == SEQITEM:
        seq = seq.object.lines[1].strip()
    elif seq_class == SEQRECORD:
        seq = str(seq.object.seq)
    return seq.strip() def _clean_seqrecord_stream(seqs):
def _write_seqrecords(seqs, fhand=None, file_format='fastq'):
def _write_seqrecord_packets(fhand, seq_packets, file_format='fastq',
                            workers=None):
def write_seq_packets(fhand, seq_packets, file_format='fastq', workers=None):
def _write_filter_trim_packets(passed_fhand, diverted_fhand, packets,
                               file_format='fastq', workers=None,
                               seqs_diverted=SEQS_FILTERED_OUT):
def write_filter_packets(passed_fhand, filtered_fhand, filter_packets,
                         file_format='fastq', workers=None):
def write_trim_packets(passed_fhand, orphan_fhand, trim_packets,
                       file_format='fastq', workers=None):
def title2ids(title):
def read_seq_packets(fhands, size=get_setting('PACKET_SIZE'), out_format=None,
                     prefered_seq_classes=None):
def _read_seqrecord_packets(fhands, size=get_setting('PACKET_SIZE')):
def _read_seqrecords(fhands):
def seqio(in_fhands, out_fhand, out_format, copy_if_same_format=True):
def fastaqual_to_fasta(seq_fhand, qual_fhand, out_fhand):
def guess_seq_type(fhand):
def _get_name_from_lines(lines):
def _get_name_from_chunk_fastq(lines):
def _line_is_not_empty(line):
def _itemize_fastq_singleline(fhand):
def _itemize_fastx(fhand):
def _itemize_fastx_multiline(fhand):  # this is a generator function
def _read_seqitems(fhands):
def _write_seqitems(items, fhand, file_format):
def write_seqs(seqs, fhand=None, file_format=None):
def read_seqs(fhands, out_format=None, prefered_seq_classes=None): class InterleaveError(Exception):
    'Error in the alternation of forward and reverse reads'
    pass class PairDirectionError(Exception):
    'Error in parser trying to get direction from description'
    pass class ItemsNotSortedError(Exception):
    'A list or generator of items is not sorted as expected'
    pass class SeqItem(_SeqItem):
    def __new__(cls, name, lines, annotations=None):
def get_title(seq):
def get_description(seq):
def get_name(seq):
def get_file_format(seq):
def _break():
def _is_fastq_plus_line(line, seq_name):
def _get_seqitem_quals(seq):
def get_str_seq(seq):
def get_length(seq):
def _get_seqitem_qualities(seqwrap):
def get_int_qualities(seq):
def _int_quals_to_str_quals(int_quals, out_format):
def get_str_qualities(seq, out_format=None):
def get_annotations(seq):
def _copy_seqrecord(seqrec, seq=None, name=None, id_=None):
def _copy_seqitem(seqwrapper, seq=None, name=None):
def copy_seq(seqwrapper, seq=None, name=None):
def _slice_seqitem(seqwrap, start, stop):
def slice_seq(seq, start=None, stop=None):
def assing_kind_to_seqs(kind, seqs, file_format):
SANGER_QUALS = {chr(i): i - 33 for i in range(33, 127)}
ILLUMINA_QUALS = {chr(i): i - 64 for i in range(64, 127)}
SANGER_STRS = {i - 33: chr(i)  for i in range(33, 127)}
ILLUMINA_STRS = {i - 64: chr(i) for i in range(64, 127)} def set_format(fhand, file_format):
    'It sets the file format in the global inventory variable'
    id_ = _get_fhand_id(fhand)
    if id_ in FILEFORMAT_INVENTORY:
        msg = 'The given instance already setted its file format'
        raise RuntimeError(msg)
    FILEFORMAT_INVENTORY[id_] = file_format
import unittest
import os
from cStringIO import StringIO
from subprocess import check_output, call
from tempfile import NamedTemporaryFile
from Bio.bgzf import BgzfReader
from Bio.SeqRecord import SeqRecord
from Bio.Seq import Seq
from crumbs.seq.pairs import (match_pairs, interleave_pairs,
                              deinterleave_pairs,
                              group_pairs, group_pairs_by_name,
                              _parse_pair_direction_and_name_from_title,
                              _parse_pair_direction_and_name)
from crumbs.iterutils import flat_zip_longest
from crumbs.utils.tags import FWD, SEQRECORD, SEQITEM
from crumbs.utils.bin_utils import SEQ_BIN_DIR
from crumbs.utils.test_utils import TEST_DATA_DIR
from crumbs.seq.seq import get_str_seq
from crumbs.seq.seqio import read_seqs, assing_kind_to_seqs
from crumbs.exceptions import (InterleaveError, PairDirectionError,
                               ItemsNotSortedError)
from crumbs.seq.seq import SeqWrapper, SeqItem
from crumbs.seq.utils.file_formats import set_format

        # with reads with no direcction
        file1 = os.path.join(TEST_DATA_DIR, 'pairend7.sfastq')
        file2 = os.path.join(TEST_DATA_DIR, 'pairend2.sfastq')
        fhand = NamedTemporaryFile()
        fhand.write(open(file1).read())
        fhand.write(open(file2).read())
        fhand.flush()
        seqs = read_seqs([fhand])

        out_fhand = StringIO()
        orphan_out_fhand = StringIO()
        out_format = 'fastq'

        match_pairs(seqs, out_fhand, orphan_out_fhand, out_format,
                    ordered=False)
        output = out_fhand.getvalue()
        assert '@seq8:136:FC706VJ:2:2104:15343:197393 1:Y:18:ATCACG' in output
        assert '@seq8:136:FC706VJ:2:2104:15343:197393 2:Y:18:ATCACG' in output
        assert '@seq1:136:FC706VJ:2:2104:15343:197393 2:Y:18:ATCACG' in output
        assert '@seq1:136:FC706VJ:2:2104:15343:197393 2:Y:18:ATCACG' in output

        orp = orphan_out_fhand.getvalue()
        assert '@seq6:136:FC706VJ:2:2104:15343:197393.mpl_1' in orp
        assert '@seq7:136:FC706VJ:2:2104:15343:197393.hhhh' in orp
        assert '@seq2:136:FC706VJ:2:2104:15343:197393 2:Y:18:ATCAC' in orp

    def test_pair_direction_and_name(self):
        'it test the pair_name parser'
        title = 'seq8:136:FC706VJ:2:2104:15343:197393 1:Y:18:ATCACG'
        name, dir_ = _parse_pair_direction_and_name_from_title(title)
        assert name == 'seq8:136:FC706VJ:2:2104:15343:197393'
        assert dir_ == FWD

        title = 'seq8:136:FC706VJ:2:2104:15343:197393/1'
        name, dir_ = _parse_pair_direction_and_name_from_title(title)
        assert name == 'seq8:136:FC706VJ:2:2104:15343:197393'
        assert dir_ == FWD

        title = 'seq8:136:FC706VJ:2:2104:15343:197393.f'
        name, dir_ = _parse_pair_direction_and_name_from_title(title)
        assert name == 'seq8:136:FC706VJ:2:2104:15343:197393'
        assert dir_ == FWD

        title = 'seq8:136:FC706VJ:2:2104:15343:197393.mp12'
        try:
            name, dir_ = _parse_pair_direction_and_name_from_title(title)
            self.fail()
        except PairDirectionError:
            pass

        title = r'seq8:136:FC706VJ:2:2104:15343:197393\1'
        name, dir_ = _parse_pair_direction_and_name_from_title(title)
        assert name == 'seq8:136:FC706VJ:2:2104:15343:197393'
        assert dir_ == FWD

        # With SeqRecord
        seq = SeqRecord(id=r'seq8:136:FC706VJ:2:2104:15343:197393\1',
                        seq=Seq('ACT'))
        name, dir_ = _parse_pair_direction_and_name(SeqWrapper(SEQRECORD, seq,
                                                               None))
        assert name == 'seq8:136:FC706VJ:2:2104:15343:197393'
        assert dir_ == FWD


class PairMatcherbinTest(unittest.TestCase):
    'It test the matepair binary'
    def test_pair_matcher_bin(self):
        'It test the pair matcher binary'
        pair_matcher_bin = os.path.join(SEQ_BIN_DIR, 'pair_matcher')
        assert 'usage' in check_output([pair_matcher_bin, '-h'])

        in_fpath = os.path.join(TEST_DATA_DIR, 'pairend5.sfastq')
        out_fhand = NamedTemporaryFile()
        orphan_fhand = NamedTemporaryFile()
        check_output([pair_matcher_bin, '-o', out_fhand.name,
                      '-p', orphan_fhand.name, in_fpath])

        result = open(out_fhand.name).read()
        assert '@seq1:136:FC706VJ:2:2104:15343:197393 1:Y:18:ATCACG' in result
        assert '@seq1:136:FC706VJ:2:2104:15343:197393 2:Y:18:ATCACG' in result

        orp = open(orphan_fhand.name).read()
        assert '@seq8:136:FC706VJ:2:2104:15343:197393 2:Y:18:ATCACG' in orp

        # compressed output
        in_fpath = os.path.join(TEST_DATA_DIR, 'pairend5.sfastq')
        out_fhand = NamedTemporaryFile()
        orphan_fhand = NamedTemporaryFile()
        check_output([pair_matcher_bin, '-o', out_fhand.name,
                      '-p', orphan_fhand.name, in_fpath, '-Z'])
        result = BgzfReader(out_fhand.name).read(2000)
        assert '@seq1:136:FC706VJ:2:2104:15343:197393 1:Y:18:ATCACG' in result
        assert '@seq1:136:FC706VJ:2:2104:15343:197393 2:Y:18:ATCACG' in result

        orp = BgzfReader(orphan_fhand.name).read(2000)
        assert '@seq8:136:FC706VJ:2:2104:15343:197393 2:Y:18:ATCACG' in orp

        # unordered file
        in_fpath = os.path.join(TEST_DATA_DIR, 'pairend6.sfastq')
        out_fhand = NamedTemporaryFile()
        orphan_fhand = NamedTemporaryFile()
        check_output([pair_matcher_bin, '-o', out_fhand.name,
                      '-p', orphan_fhand.name, in_fpath, '-u'])

        result = open(out_fhand.name).read()
        assert '@seq1:136:FC706VJ:2:2104:15343:197393 1:Y:18:ATCACG' in result
        assert '@seq1:136:FC706VJ:2:2104:15343:197393 2:Y:18:ATCACG' in result

        orp = open(orphan_fhand.name).read()
        assert '@seq8:136:FC706VJ:2:2104:15343:197393 2:Y:18:ATCACG' in orp


class InterleavePairsTest(unittest.TestCase):
    'It tests the interleaving and de-interleaving of pairs'
    def test_interleave(self):
        'It interleaves two iterators with paired reads'
        file1 = os.path.join(TEST_DATA_DIR, 'pairend1.sfastq')
        file2 = os.path.join(TEST_DATA_DIR, 'pairend2.sfastq')
        fwd_seqs = list(read_seqs([open(file1)], 'fastq'))

        rev_seqs = list(read_seqs([open(file2)], 'fastq'))