ipyrad/assemble/util.py
def merge_after_pysam(data, clust):
    """
    This is for pysam post-flight merging. The input is a cluster
    for an individual locus. We have to split the clusters, write
    R1 and R2 to files then call merge_pairs(). This is not ideal,
    it's slow, but it works. This is the absolute worst way to do this,
    it bounces all the files for each locus off the disk. I/O _hog_.
    """
    try:
        r1file = tempfile.NamedTemporaryFile(mode='wb', delete=False,
                                            dir=data.dirs.edits,
                                            suffix="_R1_.fastq")
        r2file = tempfile.NamedTemporaryFile(mode='wb', delete=False,
                                            dir=data.dirs.edits,
                                            suffix="_R2_.fastq")

        r1dat = []
        r2dat = []
        for locus in clust:
            sname, seq = locus.split("\n")
            ## Have to munge the sname to make it look like fastq format
            sname = "@" + sname[1:]
            r1, r2 = seq.split("nnnn")
            r1dat.append("{}\n{}\n{}\n{}".format(sname, r1, "+", "B"*(len(r1))))
            r2dat.append("{}\n{}\n{}\n{}".format(sname, r2, "+", "B"*(len(r2))))

        r1file.write("\n".join(r1dat))
        r2file.write("\n".join(r2dat))
        r1file.close()
        r2file.close()

        ## Read in the merged data and format to return as a clust
        merged_file = tempfile.NamedTemporaryFile(mode='wb',
                                            dir=data.dirs.edits,
                                            suffix="_merged.fastq").name

        clust = []
        merge_pairs(data, [(r1file.name, r2file.name)], merged_file, 0, 1)

        with open(merged_file) as infile:
            quarts = itertools.izip(*[iter(infile)]*4)

            while 1:
                try:
                    sname, seq, _, _ = quarts.next()
                    ## Vsearch expects R2 oriented how it would be in a raw data file
                    ## i.e. revcomp, and that's also how it returns it
                    ## but we want to maintain the genomic orientation so R1 and
                    ## R2 are both on the + strand and both in ascending positional order
                    ## so here if the reads don't merge we have to revcomp R2
                    if not "_m1" in sname.rsplit(";", 1)[1]:
                        try:
                            R1, R2 = seq.split("nnnn")
                            seq = R1 + "nnnn" + revcomp(R2)
                        except ValueError as inst:
                            LOGGER.error("Failed merge_after_pysam: {} {}".format(sname, seq))
                            raise

                except StopIteration:
                    break
                ## put sname back
                sname = ">" + sname[1:]
                clust.extend([sname.strip(), seq.strip()])
    except:
        LOGGER.info("Error in merge_pairs post-refmap.")
        raise
    finally:
        for i in [r1file.name, r2file.name, merged_file]:
            if os.path.exists(i):
                log_level = logging.getLevelName(LOGGER.getEffectiveLevel())
                ## uncomment this to shit ALL over the filesystem. 
                ## if not log_level == "DEBUG":
                os.remove(i)

    return clust
def merge_after_pysam(data, clust):
    """
    This is for pysam post-flight merging. The input is a cluster
    for an individual locus. We have to split the clusters, write
    R1 and R2 to files then call merge_pairs(). This is not ideal,
    it's slow, but it works. This is the absolute worst way to do this,
    it bounces all the files for each locus off the disk. I/O _hog_.
    """
    try:
        r1file = tempfile.NamedTemporaryFile(mode='wb', delete=False,
                                            dir=data.dirs.edits,
                                            suffix="_R1_.fastq")
        r2file = tempfile.NamedTemporaryFile(mode='wb', delete=False,
                                            dir=data.dirs.edits,
                                            suffix="_R2_.fastq")

        r1dat = []
        r2dat = []
        for locus in clust:
            sname, seq = locus.split("\n")
            ## Have to munge the sname to make it look like fastq format
            sname = "@" + sname[1:]
            r1, r2 = seq.split("nnnn")
            r1dat.append("{}\n{}\n{}\n{}".format(sname, r1, "+", "B"*(len(r1))))
            r2dat.append("{}\n{}\n{}\n{}".format(sname, r2, "+", "B"*(len(r2))))

        r1file.write("\n".join(r1dat))
        r2file.write("\n".join(r2dat))
        r1file.close()
        r2file.close()

        ## Read in the merged data and format to return as a clust
        merged_file = tempfile.NamedTemporaryFile(mode='wb',
                                            dir=data.dirs.edits,
                                            suffix="_merged.fastq").name

        clust = []
        merge_pairs(data, [(r1file.name, r2file.name)], merged_file, 0, 1)

        with open(merged_file) as infile:
            quarts = itertools.izip(*[iter(infile)]*4)

            while 1:
                try:
                    sname, seq, _, _ = quarts.next()
                    ## Vsearch expects R2 oriented how it would be in a raw data file
                    ## i.e. revcomp, and that's also how it returns it
                    ## but we want to maintain the genomic orientation so R1 and
                    ## R2 are both on the + strand and both in ascending positional order
                    ## so here if the reads don't merge we have to revcomp R2
                    if not "_m1" in sname.rsplit(";", 1)[1]:
                        try:
                            R1, R2 = seq.split("nnnn")
                            seq = R1 + "nnnn" + revcomp(R2)
                        except ValueError as inst:
                            LOGGER.error("Failed merge_after_pysam: {} {}".format(sname, seq))
                            raise

                except StopIteration:
                    break
                ## put sname back
                sname = ">" + sname[1:]
                clust.extend([sname.strip(), seq.strip()])
    except:
        LOGGER.info("Error in merge_pairs post-refmap.")
        raise
    finally:
        for i in [r1file.name, r2file.name, merged_file]:
            if os.path.exists(i):
                log_level = logging.getLevelName(LOGGER.getEffectiveLevel())
                ## uncomment this to shit ALL over the filesystem. 
                ## if not log_level == "DEBUG":
                os.remove(i)

    return clust
def merge_after_pysam ( data , clust ) : try : r1file = tempfile . NamedTemporaryFile ( mode = 'wb' , delete = False , dir = data . dirs . edits , suffix = "_R1_.fastq" ) r2file = tempfile . NamedTemporaryFile ( mode = 'wb' , delete = False , dir = data . dirs . edits , suffix = "_R2_.fastq" ) r1dat = [ ] r2dat = [ ] for locus in clust : sname , seq = locus . split ( "\n" ) ## Have to munge the sname to make it look like fastq format sname = "@" + sname [ 1 : ] r1 , r2 = seq . split ( "nnnn" ) r1dat . append ( "{}\n{}\n{}\n{}" . format ( sname , r1 , "+" , "B" * ( len ( r1 ) ) ) ) r2dat . append ( "{}\n{}\n{}\n{}" . format ( sname , r2 , "+" , "B" * ( len ( r2 ) ) ) ) r1file . write ( "\n" . join ( r1dat ) ) r2file . write ( "\n" . join ( r2dat ) ) r1file . close ( ) r2file . close ( ) ## Read in the merged data and format to return as a clust merged_file = tempfile . NamedTemporaryFile ( mode = 'wb' , dir = data . dirs . edits , suffix = "_merged.fastq" ) . name clust = [ ] merge_pairs ( data , [ ( r1file . name , r2file . name ) ] , merged_file , 0 , 1 ) with open ( merged_file ) as infile : quarts = itertools . izip ( * [ iter ( infile ) ] * 4 ) while 1 : try : sname , seq , _ , _ = quarts . next ( ) ## Vsearch expects R2 oriented how it would be in a raw data file ## i.e. revcomp, and that's also how it returns it ## but we want to maintain the genomic orientation so R1 and ## R2 are both on the + strand and both in ascending positional order ## so here if the reads don't merge we have to revcomp R2 if not "_m1" in sname . rsplit ( ";" , 1 ) [ 1 ] : try : R1 , R2 = seq . split ( "nnnn" ) seq = R1 + "nnnn" + revcomp ( R2 ) except ValueError as inst : LOGGER . error ( "Failed merge_after_pysam: {} {}" . format ( sname , seq ) ) raise except StopIteration : break ## put sname back sname = ">" + sname [ 1 : ] clust . extend ( [ sname . strip ( ) , seq . strip ( ) ] ) except : LOGGER . info ( "Error in merge_pairs post-refmap." ) raise finally : for i in [ r1file . name , r2file . name , merged_file ] : if os . path . exists ( i ) : log_level = logging . getLevelName ( LOGGER . getEffectiveLevel ( ) ) ## uncomment this to shit ALL over the filesystem.  ## if not log_level == "DEBUG": os . remove ( i ) return clust
This is for pysam post-flight merging. The input is a cluster
    for an individual locus. We have to split the clusters, write
    R1 and R2 to files then call merge_pairs(). This is not ideal,
    it's slow, but it works. This is the absolute worst way to do this,
    it bounces all the files for each locus off the disk. I/O _hog_.
This is for pysam post - flight merging . The input is a cluster for an individual locus . We have to split the clusters write R1 and R2 to files then call merge_pairs () . This is not ideal it s slow but it works . This is the absolute worst way to do this it bounces all the files for each locus off the disk . I / O _hog_ .