signature:
def main ( ) :
    """main function"""

implementation:
# # parse params file input ( returns to stdout if - - help or - - version )
args = parse_command_line ( )
print ( HEADER . format ( ip . __version__ ) )
# # set random seed
np . random . seed ( args . rseed )
# # debugger - - - - -
if os . path . exists ( ip . __debugflag__ ) :
    os . remove ( ip . __debugflag__ )
if args . debug :
    print ( "\n  ** Enabling debug mode ** " )
    ip . _debug_on ( )
# # if JSON , load existing Tetrad analysis - - - - -
if args . json :
    data = ipa . tetrad ( name = args . name , workdir = args . workdir , load = True )
    # # if force then remove all results
    if args . force :
        data . _refresh ( )
# # else create a new tmp assembly for the seqarray - - - - -
else : # # create new Tetrad class Object if it doesn ' t exist
    newjson = os . path . join ( args . workdir , args . name + '.tet.json' )
    # # if not quiet . . .
    print ( "tetrad instance: {}" . format ( args . name ) )
    if ( not os . path . exists ( newjson ) ) or args . force : # # purge any files associated with this name if forced
        if args . force : # # init an object in the correct location just to refresh
            ipa . tetrad ( name = args . name , workdir = args . workdir , data = args . seq , initarr = False , save_invariants = args . invariants , cli = True , quiet = True ) . _refresh ( )
        # # create new tetrad object
        data = ipa . tetrad ( name = args . name , workdir = args . workdir , method = args . method , data = args . seq , resolve = args . resolve , mapfile = args . map , guidetree = args . tree , nboots = args . boots , nquartets = args . nquartets , cli = True , save_invariants = args . invariants , )
    else :
        raise SystemExit ( QUARTET_EXISTS . format ( args . name , args . workdir , args . workdir , args . name , args . name ) )
# # boots can be set either for a new object or loaded JSON to continue it
if args . boots :
    data . params . nboots = int ( args . boots )
# # if ipyclient is running ( and matched profile ) then use that one
if args . ipcluster :
    ipyclient = ipp . Client ( profile = args . ipcluster )
    data . _ipcluster [ "cores" ] = len ( ipyclient )
# # if not then we need to register and launch an ipcluster instance
else : # # set CLI ipcluster terms
    ipyclient = None
    data . _ipcluster [ "cores" ] = args . cores if args . cores else detect_cpus ( )
    data . _ipcluster [ "engines" ] = "Local"
    if args . MPI :
        data . _ipcluster [ "engines" ] = "MPI"
        if not args . cores :
            raise IPyradWarningExit ( "must provide -c argument with --MPI" )
    # # register to have a cluster - id with " ip - name "
    data = register_ipcluster ( data )
# # message about whether we are continuing from existing
if data . checkpoint . boots :
    print ( LOADING_MESSAGE . format ( data . name , data . params . method , data . checkpoint . boots ) )
# # run tetrad main function within a wrapper . The wrapper creates an
# # ipyclient view and appends to the list of arguments to run ' run ' .
data . run ( force = args . force , ipyclient = ipyclient )


signature:
def _command_list ( self ) :
    """build the command list"""

implementation:
# # base args
cmd = [ self . params . binary , "-i" , OPJ ( self . workdir , self . name + ".treemix.in.gz" ) , "-o" , OPJ ( self . workdir , self . name ) , ]
# # addon params
args = [ ]
for key , val in self . params :
    if key not in [ "minmap" , "binary" ] :
        if key == "g" :
            if val [ 0 ] :
                args += [ "-" + key , str ( val [ 0 ] ) , str ( val [ 1 ] ) ]
        elif key == "global_" :
            if val :
                args += [ "-" + key [ : - 1 ] ]
        elif key in [ "se" , "global" , "noss" ] :
            if val :
                args += [ "-" + key ]
        else :
            if val :
                args += [ "-" + key , str ( val ) ]
return cmd + args


signature:
def _subsample ( self ) :
    """returns a subsample of unlinked snp sites"""

implementation:
spans = self . maparr
samp = np . zeros ( spans . shape [ 0 ] , dtype = np . uint64 )
for i in xrange ( spans . shape [ 0 ] ) :
    samp [ i ] = np . random . randint ( spans [ i , 0 ] , spans [ i , 1 ] , 1 )
return samp


signature:
def copy ( self , name ) :
    """Returns a copy of the treemix object with the same parameter settings
 but with the files attributes cleared , and with a new ' name ' attribute .
 Parameters
 name ( str ) :
 A name for the new copied treemix bject that will be used for the
 output files created by the object ."""

implementation:
# # make deepcopy of self . _ _ dict _ _ but do not copy async objects
subdict = { i : j for i , j in self . __dict__ . iteritems ( ) if i != "asyncs" }
newdict = copy . deepcopy ( subdict )
# # make back into a Treemix object
# if name = = self . name :
# raise Exception ( " new object name must be different from its parent " )
newobj = Treemix ( data = newdict [ "data" ] , name = name , workdir = newdict [ "workdir" ] , imap = { i : j for i , j in newdict [ "imap" ] . items ( ) } , mapfile = newdict [ 'mapfile' ] , minmap = { i : j for i , j in newdict [ "minmap" ] . items ( ) } , seed = np . random . randint ( 0 , int ( 1e9 ) ) , )
# # update special dict attributes but not files
for key , val in newobj . params . __dict__ . iteritems ( ) :
    newobj . params . __setattr__ ( key , self . params . __getattribute__ ( key ) )
# for key , val in newobj . filters . _ _ dict _ _ . iteritems ( ) :
# newobj . filters . _ _ setattr _ _ ( key , self . filters . _ _ getattribute _ _ ( key ) )
# # new object must have a different name than it ' s parent
return newobj


signature:
def draw ( self , axes ) :
    """Returns a treemix plot on a toyplot . axes object ."""

implementation:
# # create a toytree object from the treemix tree result
tre = toytree . tree ( newick = self . results . tree )
tre . draw ( axes = axes , use_edge_lengths = True , tree_style = 'c' , tip_labels_align = True , edge_align_style = { "stroke-width" : 1 } ) ;
# # get coords
for admix in self . results . admixture : # # parse admix event
    pidx , pdist , cidx , cdist , weight = admix
    a = _get_admix_point ( tre , pidx , pdist )
    b = _get_admix_point ( tre , cidx , cdist )
    # # add line for admixture edge
    mark = axes . plot ( a = ( a [ 0 ] , b [ 0 ] ) , b = ( a [ 1 ] , b [ 1 ] ) , style = { "stroke-width" : 10 * weight , "stroke-opacity" : 0.95 , "stroke-linecap" : "round" } )
    # # add points at admixture sink
    axes . scatterplot ( a = ( b [ 0 ] ) , b = ( b [ 1 ] ) , size = 8 , title = "weight: {}" . format ( weight ) , )
# # add scale bar for edge lengths
axes . y . show = False
axes . x . ticks . show = True
axes . x . label . text = "Drift parameter"
return axes


signature:
def _resolveambig ( subseq ) :
    """Randomly resolves iupac hetero codes . This is a shortcut
 for now , we could instead use the phased alleles in RAD loci ."""

implementation:
N = [ ]
for col in subseq :
    rand = np . random . binomial ( 1 , 0.5 )
    N . append ( [ _AMBIGS [ i ] [ rand ] for i in col ] )
return np . array ( N )


signature:
def _count_PIS ( seqsamp , N ) :
    """filters for loci with > = N PIS"""

implementation:
counts = [ Counter ( col ) for col in seqsamp . T if not ( "-" in col or "N" in col ) ]
pis = [ i . most_common ( 2 ) [ 1 ] [ 1 ] > 1 for i in counts if len ( i . most_common ( 2 ) ) > 1 ]
if sum ( pis ) >= N :
    return sum ( pis )
else :
    return 0


signature:
def write_nexus_files ( self , force = False , quiet = False ) :
    """Write nexus files to { workdir } / { name } / [ 0 - N ] . nex , If the directory already
 exists an exception will be raised unless you use the force flag which
 will remove all files in the directory .
 Parameters :
 force ( bool ) :
 If True then all files in { workdir } / { name } / * . nex * will be removed ."""

implementation:
# # clear existing files
existing = glob . glob ( os . path . join ( self . workdir , self . name , "*.nex" ) )
if any ( existing ) :
    if force :
        for rfile in existing :
            os . remove ( rfile )
    else :
        path = os . path . join ( self . workdir , self . name )
        raise IPyradWarningExit ( EXISTING_NEX_FILES . format ( path ) )
# # parse the loci or alleles file
with open ( self . files . data ) as infile :
    loci = iter ( infile . read ( ) . strip ( ) . split ( "|\n" ) )
# # use entered samples or parse them from the file
if not self . samples :
    with open ( self . files . data ) as infile :
        samples = set ( ( i . split ( ) [ 0 ] for i in infile . readlines ( ) if "//" not in i ) )
else :
    samples = set ( self . samples )
# # keep track of how many loci pass filtering
totn = len ( samples )
nloci = 0
# # this set is just used for matching , then we randomly
# # subsample for real within the locus so it varies
if self . _alleles :
    msamples = { i + rbin ( ) for i in samples }
else :
    msamples = samples
# # write subsampled set of loci
for loc in loci : # # get names and seqs from locus
    dat = loc . split ( "\n" ) [ : - 1 ]
    try :
        names = [ i . split ( ) [ 0 ] for i in dat ]
        snames = set ( names )
        seqs = np . array ( [ list ( i . split ( ) [ 1 ] ) for i in dat ] )
    except IndexError :
        print ( ALLELESBUGFIXED )
        continue
    # # check name matches
    if len ( snames . intersection ( msamples ) ) == totn : # # prune sample names if alleles . Done here so it is randomly
    # # different in every locus which allele is selected from
    # # each sample ( e . g . , 0 or 1)
        if self . _alleles :
            _samples = [ i + rbin ( ) for i in samples ]
        else :
            _samples = samples
        # # re - order seqs to be in set order
        seqsamp = seqs [ [ names . index ( tax ) for tax in _samples ] ]
        # # resolve ambiguities randomly if . loci file otherwise
        # # sample one of the alleles if . alleles file .
        if not self . _alleles :
            seqsamp = _resolveambig ( seqsamp )
        # # find parsimony informative sites
        if _count_PIS ( seqsamp , self . params . minsnps ) : # # keep the locus
            nloci += 1
            # # remove empty columns given this sampling
            copied = seqsamp . copy ( )
            copied [ copied == "-" ] == "N"
            rmcol = np . all ( copied == "N" , axis = 0 )
            seqsamp = seqsamp [ : , ~ rmcol ]
            # # write nexus file
            if self . _alleles : # # trim off the allele number
                samps = [ i . rsplit ( "_" , 1 ) [ 0 ] for i in _samples ]
                mdict = dict ( zip ( samps , [ i . tostring ( ) for i in seqsamp ] ) )
            else :
                mdict = dict ( zip ( _samples , [ i . tostring ( ) for i in seqsamp ] ) )
            self . _write_nex ( mdict , nloci )
            # # quit early if using maxloci
            if nloci == self . params . maxloci :
                break
# # print data size
if not quiet :
    path = os . path . join ( self . workdir , self . name )
    path = path . replace ( os . path . expanduser ( "~" ) , "~" )
    print ( "wrote {} nexus files to {}" . format ( nloci , path ) )


signature:
def run ( self , steps = None , ipyclient = None , force = False , quiet = False ) :
    """Submits an ordered list of jobs to a load - balancer to complete
 the following tasks , and reports a progress bar :
 (1 ) Write nexus files for each locus
 (2 ) Run mrBayes on each locus to get a posterior of gene trees
 (3 ) Run mbsum ( a bucky tool ) on the posterior set of trees
 (4 ) Run Bucky on the summarized set of trees for all alpha values .
 Parameters :
 ipyclient ( ipyparallel . Client ( ) )
 A connected ipyparallel Client object used to distribute jobs
 force ( bool ) :
 Whether to overwrite existing files with the same name and workdir
 if they exist . Default is False .
 quiet ( bool ) :
 Whether to suppress progress information . Default is False .
 steps ( list ) :
 A list of integers of steps to perform . This is useful if a
 job was interrupted , or you created a new bucky object copy ,
 or you wish to run an analysis under a new set of parameters ,
 after having run it once . For example , if you finished running
 steps 1 and 2 ( write nexus files and infer mrbayes posteriors ) ,
 but you want to rerun steps 3 and 4 with new settings , then you
 could enter ` steps = [ 3,4 ] ` and also ` force = True ` to run steps 3
 and 4 with a new set of parameters . Default argument is None
 which means run all steps ."""

implementation:
# # require ipyclient
if not ipyclient :
    raise IPyradWarningExit ( "an ipyclient object is required" )
# # check the steps argument
if not steps :
    steps = [ 1 , 2 , 3 , 4 ]
if isinstance ( steps , ( int , str ) ) :
    steps = [ int ( i ) for i in [ steps ] ]
if isinstance ( steps , list ) :
    if not all ( isinstance ( i , int ) for i in steps ) :
        raise IPyradWarningExit ( "steps must be a list of integers" )
# # run steps - - - - -
# # TODO : wrap this function so it plays nice when interrupted .
if 1 in steps :
    self . write_nexus_files ( force = force , quiet = quiet )
if 2 in steps :
    self . run_mrbayes ( force = force , quiet = quiet , ipyclient = ipyclient )
if 3 in steps :
    self . run_mbsum ( force = force , quiet = quiet , ipyclient = ipyclient )
if 4 in steps :
    self . run_bucky ( force = force , quiet = quiet , ipyclient = ipyclient )
# # make sure jobs are done if waiting ( TODO : maybe make this optional )
ipyclient . wait ( )


signature:
def _write_nex ( self , mdict , nlocus ) :
    """function that takes a dictionary mapping names to sequences ,
 and a locus number , and writes it as a NEXUS file with a mrbayes
 analysis block given a set of mcmc arguments ."""

implementation:
# # create matrix as a string
max_name_len = max ( [ len ( i ) for i in mdict ] )
namestring = "{:<" + str ( max_name_len + 1 ) + "} {}\n"
matrix = ""
for i in mdict . items ( ) :
    matrix += namestring . format ( i [ 0 ] , i [ 1 ] )
# # ensure dir
minidir = os . path . realpath ( os . path . join ( self . workdir , self . name ) )
if not os . path . exists ( minidir ) :
    os . makedirs ( minidir )
# # write nexus block
handle = os . path . join ( minidir , "{}.nex" . format ( nlocus ) )
with open ( handle , 'w' ) as outnex :
    outnex . write ( NEXBLOCK . format ( ** { "ntax" : len ( mdict ) , "nchar" : len ( mdict . values ( ) [ 0 ] ) , "matrix" : matrix , "ngen" : self . params . mb_mcmc_ngen , "sfreq" : self . params . mb_mcmc_sample_freq , "burnin" : self . params . mb_mcmc_burnin , } ) )


signature:
def run_mbsum ( self , ipyclient , force = False , quiet = False ) :
    """Sums two replicate mrbayes runs for each locus"""

implementation:
minidir = os . path . realpath ( os . path . join ( self . workdir , self . name ) )
trees1 = glob . glob ( os . path . join ( minidir , "*.run1.t" ) )
trees2 = glob . glob ( os . path . join ( minidir , "*.run2.t" ) )
# # clear existing files
existing = glob . glob ( os . path . join ( self . workdir , self . name , "*.sumt" ) )
if any ( existing ) :
    if force :
        for rfile in existing :
            os . remove ( rfile )
    else :
        path = os . path . join ( self . workdir , self . name )
        raise IPyradWarningExit ( EXISTING_SUMT_FILES . format ( path ) )
# # load balancer
lbview = ipyclient . load_balanced_view ( )
# # submit each to be processed
asyncs = [ ]
for tidx in xrange ( len ( trees1 ) ) :
    rep1 = trees1 [ tidx ]
    rep2 = trees2 [ tidx ]
    outname = os . path . join ( minidir , str ( tidx ) + ".sumt" )
    async = lbview . apply ( _call_mbsum , * ( rep1 , rep2 , outname ) )
    asyncs . append ( async )
# # track progress
start = time . time ( )
printstr = "[mbsum] sum replicate runs      | {} | "
while 1 :
    ready = [ i . ready ( ) for i in asyncs ]
    elapsed = datetime . timedelta ( seconds = int ( time . time ( ) - start ) )
    if not quiet :
        progressbar ( len ( ready ) , sum ( ready ) , printstr . format ( elapsed ) , spacer = "" )
    if len ( ready ) == sum ( ready ) :
        if not quiet :
            print ( "" )
        break
    else :
        time . sleep ( 0.1 )
# # check success
for async in asyncs :
    if not async . successful ( ) :
        raise IPyradWarningExit ( async . result ( ) )


signature:
def run_mrbayes ( self , ipyclient , force = False , quiet = False ) :
    """calls the mrbayes block in each nexus file ."""

implementation:
# # get all the nexus files for this object
minidir = os . path . realpath ( os . path . join ( self . workdir , self . name ) )
nexus_files = glob . glob ( os . path . join ( minidir , "*.nex" ) )
# # clear existing files
# existing = glob . glob ( os . path . join ( self . workdir , self . name , " * . nex " ) )
existing = glob . glob ( os . path . join ( minidir , "*.nex.*" ) )
if any ( existing ) :
    if force :
        for rfile in existing :
            os . remove ( rfile )
    else :
        raise IPyradWarningExit ( EXISTING_NEXdot_FILES . format ( minidir ) )
# # write new nexus files , or should users do that before this ?
# self . write _ nexus _ files ( force = True )
# # load balancer
lbview = ipyclient . load_balanced_view ( )
# # submit each to be processed
asyncs = [ ]
for nex in nexus_files :
    async = lbview . apply ( _call_mb , nex )
    asyncs . append ( async )
# # track progress
start = time . time ( )
printstr = "[mb] infer gene-tree posteriors | {} | "
while 1 :
    ready = [ i . ready ( ) for i in asyncs ]
    elapsed = datetime . timedelta ( seconds = int ( time . time ( ) - start ) )
    if not quiet :
        progressbar ( len ( ready ) , sum ( ready ) , printstr . format ( elapsed ) , spacer = "" )
    if len ( ready ) == sum ( ready ) :
        if not quiet :
            print ( "" )
        break
    else :
        time . sleep ( 0.1 )
# # check success
for async in asyncs :
    if not async . successful ( ) :
        raise IPyradWarningExit ( async . result ( ) )


signature:
def run_bucky ( self , ipyclient , force = False , quiet = False , subname = False ) :
    """Runs bucky for a given set of parameters and stores the result
 to the ipa . bucky object . The results will be stored by default
 with the name ' { name } - { alpha } ' unless a argument is passed for
 ' subname ' to customize the output name .
 Parameters :
 subname ( str ) :
 A custom name prefix for the output files produced by the bucky
 analysis and output into the { workdir } / { name } directory .
 force ( bool ) :
 If True then existing result files with the same name prefix
 will be overwritten .
 quiet ( bool ) :
 If True the progress bars will be suppressed .
 ipyclient ( ipyparallel . Client )
 An active ipyparallel client to distribute jobs to ."""

implementation:
# # check for existing results files
minidir = os . path . realpath ( os . path . join ( self . workdir , self . name ) )
infiles = glob . glob ( os . path . join ( minidir , "*.sumt" ) )
outroot = os . path . realpath ( os . path . join ( self . workdir , self . name ) )
# # build alpha list
if isinstance ( self . params . bucky_alpha , list ) :
    alphas = self . params . bucky_alpha
else :
    alphas = [ self . params . bucky_alpha ]
# # load balancer
lbview = ipyclient . load_balanced_view ( )
# # submit each to be processed
asyncs = [ ]
for alpha in alphas :
    pathname = os . path . join ( outroot , "CF-a" + str ( alpha ) )
    if ( os . path . exists ( pathname ) ) and ( force != True ) :
        print ( "BUCKy results already exist for this object at alpha={}\n" . format ( alpha ) + "use force=True to overwrite existing results" )
    else :
        args = [ alpha , self . params . bucky_nchains , self . params . bucky_nreps , self . params . bucky_niter , pathname , infiles ]
        async = lbview . apply ( _call_bucky , * args )
        asyncs . append ( async )
# # track progress
start = time . time ( )
printstr = "[bucky] infer CF posteriors     | {} | "
while 1 :
    ready = [ i . ready ( ) for i in asyncs ]
    elapsed = datetime . timedelta ( seconds = int ( time . time ( ) - start ) )
    if not quiet :
        progressbar ( len ( ready ) , sum ( ready ) , printstr . format ( elapsed ) , spacer = "" )
    if len ( ready ) == sum ( ready ) :
        if not quiet :
            print ( "" )
        break
    else :
        time . sleep ( 0.1 )
# # check success
for async in asyncs :
    if not async . successful ( ) :
        raise IPyradWarningExit ( async . result ( ) )


signature:
def _get_samples ( self , samples ) :
    """Internal function . Prelude for each step ( ) to read in perhaps
 non empty list of samples to process . Input is a list of sample names ,
 output is a list of sample objects ."""

implementation:
# # if samples not entered use all samples
if not samples :
    samples = self . samples . keys ( )
# # Be nice and allow user to pass in only one sample as a string ,
# # rather than a one element list . When you make the string into a list
# # you have to wrap it in square braces or else list makes a list of
# # each character individually .
if isinstance ( samples , str ) :
    samples = list ( [ samples ] )
# # if sample keys , replace with sample obj
assert isinstance ( samples , list ) , "to subselect samples enter as a list, e.g., [A, B]."
newsamples = [ self . samples . get ( key ) for key in samples if self . samples . get ( key ) ]
strnewsamples = [ i . name for i in newsamples ]
# # are there any samples that did not make it into the dict ?
badsamples = set ( samples ) . difference ( set ( strnewsamples ) )
if badsamples :
    outstring = ", " . join ( badsamples )
    raise IPyradError ( "Unrecognized Sample name(s) not linked to {}: {}" . format ( self . name , outstring ) )
# # require Samples
assert newsamples , "No Samples passed in and none in assembly {}" . format ( self . name )
return newsamples


signature:
def _name_from_file ( fname , splitnames , fields ) :
    """internal func : get the sample name from any pyrad file"""

implementation:
# # allowed extensions
file_extensions = [ ".gz" , ".fastq" , ".fq" , ".fasta" , ".clustS" , ".consens" ]
base , _ = os . path . splitext ( os . path . basename ( fname ) )
# # remove read number from name
base = base . replace ( "_R1_." , "." ) . replace ( "_R1_" , "" ) . replace ( "_R1." , "." )
# # remove extensions , retains ' . ' in file names .
while 1 :
    tmpb , tmpext = os . path . splitext ( base )
    if tmpext in file_extensions :
        base = tmpb
    else :
        break
if fields :
    namebits = base . split ( splitnames )
    base = [ ]
    for field in fields :
        try :
            base . append ( namebits [ field ] )
        except IndexError :
            pass
    base = splitnames . join ( base )
if not base :
    raise IPyradError ( """
    Found invalid/empty filename in link_fastqs. Check splitnames argument.
    """ )
return base


signature:
def _read_sample_names ( fname ) :
    """Read in sample names from a plain text file . This is a convenience
 function for branching so if you have tons of sample names you can
 pass in a file rather than having to set all the names at the command
 line ."""

implementation:
try :
    with open ( fname , 'r' ) as infile :
        subsamples = [ x . split ( ) [ 0 ] for x in infile . readlines ( ) if x . strip ( ) ]
except Exception as inst :
    print ( "Failed to read input file with sample names.\n{}" . format ( inst ) )
    raise inst
return subsamples


signature:
def _expander ( namepath ) :
    """expand . / ~ and . . / designators in location names"""

implementation:
if "~" in namepath :
    namepath = os . path . expanduser ( namepath )
else :
    namepath = os . path . abspath ( namepath )
return namepath


signature:
def merge ( name , assemblies ) :
    """Creates and returns a new Assembly object in which samples from two or more
 Assembly objects with matching names are ' merged ' . Merging does not affect
 the actual files written on disk , but rather creates new Samples that are
 linked to multiple data files , and with stats summed ."""

implementation:
# # checks
assemblies = list ( assemblies )
# # create new Assembly as a branch ( deepcopy )
merged = assemblies [ 0 ] . branch ( name )
# # get all sample names from all Assemblies
allsamples = set ( merged . samples . keys ( ) )
for iterass in assemblies [ 1 : ] :
    allsamples . update ( set ( iterass . samples . keys ( ) ) )
# # Make sure we have the max of all values for max frag length
# # from all merging assemblies .
merged . _hackersonly [ "max_fragment_length" ] = max ( [ x . _hackersonly [ "max_fragment_length" ] for x in assemblies ] )
# # warning message ?
warning = 0
# # iterate over assembly objects , skip first already copied
for iterass in assemblies [ 1 : ] : # # iterate over allsamples , add if not in merged
    for sample in iterass . samples : # # iterate over stats , skip ' state '
        if sample not in merged . samples :
            merged . samples [ sample ] = copy . deepcopy ( iterass . samples [ sample ] )
            # # if barcodes data present then keep it
            if iterass . barcodes . get ( sample ) :
                merged . barcodes [ sample ] = iterass . barcodes [ sample ]
        else : # # merge stats and files of the sample
            for stat in merged . stats . keys ( ) [ 1 : ] :
                merged . samples [ sample ] . stats [ stat ] += iterass . samples [ sample ] . stats [ stat ]
            # # merge file references into a list
            for filetype in [ 'fastqs' , 'edits' ] :
                merged . samples [ sample ] . files [ filetype ] += iterass . samples [ sample ] . files [ filetype ]
            if iterass . samples [ sample ] . files [ "clusters" ] :
                warning += 1
# # print warning if clusters or later was present in merged assembly
if warning :
    print ( """\
    Warning: the merged Assemblies contained Samples that are identically named,
    and so ipyrad has attempted to merge these Samples. This is perfectly fine to
    do up until step 3, but not after, because at step 3 all reads for a Sample
    should be included during clustering/mapping. Take note, you can merge Assemblies
    at any step *if they do not contain the same Samples*, however, here that is not
    the case. If you wish to proceed with this merged Assembly you will have to
    start from step 3, therefore the 'state' of the Samples in this new merged
    Assembly ({}) have been set to 2.
    """ . format ( name ) )
    for sample in merged . samples :
        merged . samples [ sample ] . stats . state = 2
        # # clear stats
        for stat in [ "refseq_mapped_reads" , "refseq_unmapped_reads" , "clusters_total" , "clusters_hidepth" , "hetero_est" , "error_est" , "reads_consens" ] :
            merged . samples [ sample ] . stats [ stat ] = 0
        # # clear files
        for ftype in [ "mapped_reads" , "unmapped_reads" , "clusters" , "consens" , "database" ] :
            merged . samples [ sample ] . files [ ftype ] = [ ]
# # Set the values for some params that don ' t make sense inside
# # merged assemblies
merged_names = ", " . join ( [ x . name for x in assemblies ] )
merged . paramsdict [ "raw_fastq_path" ] = "Merged: " + merged_names
merged . paramsdict [ "barcodes_path" ] = "Merged: " + merged_names
merged . paramsdict [ "sorted_fastq_path" ] = "Merged: " + merged_names
# # return the new Assembly object
merged . save ( )
return merged


signature:
def _bufcountlines ( filename , gzipped ) :
    """fast line counter . Used to quickly sum number of input reads when running
 link _ fastqs to append files ."""

implementation:
if gzipped :
    fin = gzip . open ( filename )
else :
    fin = open ( filename )
nlines = 0
buf_size = 1024 * 1024
read_f = fin . read
# loop optimization
buf = read_f ( buf_size )
while buf :
    nlines += buf . count ( '\n' )
    buf = read_f ( buf_size )
fin . close ( )
return nlines


signature:
def _zbufcountlines ( filename , gzipped ) :
    """faster line counter"""

implementation:
if gzipped :
    cmd1 = [ "gunzip" , "-c" , filename ]
else :
    cmd1 = [ "cat" , filename ]
cmd2 = [ "wc" ]
proc1 = sps . Popen ( cmd1 , stdout = sps . PIPE , stderr = sps . PIPE )
proc2 = sps . Popen ( cmd2 , stdin = proc1 . stdout , stdout = sps . PIPE , stderr = sps . PIPE )
res = proc2 . communicate ( ) [ 0 ]
if proc2 . returncode :
    raise IPyradWarningExit ( "error zbufcountlines {}:" . format ( res ) )
LOGGER . info ( res )
nlines = int ( res . split ( ) [ 0 ] )
return nlines


signature:
def _tuplecheck ( newvalue , dtype = str ) :
    """Takes a string argument and returns value as a tuple .
 Needed for paramfile conversion from CLI to set _ params args"""

implementation:
if isinstance ( newvalue , list ) :
    newvalue = tuple ( newvalue )
if isinstance ( newvalue , str ) :
    newvalue = newvalue . rstrip ( ")" ) . strip ( "(" )
    try :
        newvalue = tuple ( [ dtype ( i . strip ( ) ) for i in newvalue . split ( "," ) ] )
    # # Type error is thrown by tuple if it ' s applied to a non - iterable .
    except TypeError :
        newvalue = tuple ( dtype ( newvalue ) )
    # # If dtype fails to cast any element of newvalue
    except ValueError :
        LOGGER . info ( "Assembly.tuplecheck() failed to cast to {} - {}" . format ( dtype , newvalue ) )
        raise
    except Exception as inst :
        LOGGER . info ( inst )
        raise SystemExit ( "\nError: Param`{}` is not formatted correctly.\n({})\n" . format ( newvalue , inst ) )
return newvalue


signature:
def _paramschecker ( self , param , newvalue ) :
    """Raises exceptions when params are set to values they should not be"""

implementation:
if param == 'assembly_name' : # # Make sure somebody doesn ' t try to change their assembly _ name , bad
# # things would happen . Calling set _ params on assembly _ name only raises
# # an informative error . Assembly _ name is set at Assembly creation time
# # and is immutable .
    raise IPyradWarningExit ( CANNOT_CHANGE_ASSEMBLY_NAME )
elif param == 'project_dir' :
    expandpath = _expander ( newvalue )
    if not expandpath . startswith ( "/" ) :
        if os . path . exists ( expandpath ) :
            expandpath = _expander ( expandpath )
    # # Forbid spaces in path names
    if " " in expandpath :
        raise IPyradWarningExit ( BAD_PROJDIR_NAME . format ( expandpath ) )
    self . paramsdict [ "project_dir" ] = expandpath
    self . dirs [ "project" ] = expandpath
# # ` Merged : ` in newvalue for raw _ fastq _ path indicates that this
# # assembly is a merge of several others , so this param has no
# # value for this assembly
elif param == 'raw_fastq_path' :
    if newvalue and not "Merged:" in newvalue :
        fullrawpath = _expander ( newvalue )
        if os . path . isdir ( fullrawpath ) :
            raise IPyradWarningExit ( RAW_PATH_ISDIR . format ( fullrawpath ) )
        # # if something is found in path
        elif glob . glob ( fullrawpath ) :
            self . paramsdict [ 'raw_fastq_path' ] = fullrawpath
        # # else allow empty , tho it can still raise an error in step1
        else :
            raise IPyradWarningExit ( NO_RAW_FILE . format ( fullrawpath ) )
    else :
        self . paramsdict [ 'raw_fastq_path' ] = ""
# # ` Merged : ` in newvalue for barcodes _ path indicates that this
# # assembly is a merge of several others , so this param has no
# # value for this assembly
elif param == 'barcodes_path' : # # if a value was entered check that it exists
    if newvalue and not "Merged:" in newvalue : # # also allow for fuzzy match in names using glob
        fullbarpath = glob . glob ( _expander ( newvalue ) ) [ 0 ]
        # # raise error if file is not found
        if not os . path . exists ( fullbarpath ) :
            raise IPyradWarningExit ( BARCODE_NOT_FOUND . format ( fullbarpath ) )
        else :
            self . paramsdict [ 'barcodes_path' ] = fullbarpath
            self . _link_barcodes ( )
    # # if no path was entered then set barcodes path to empty .
    # # this is checked again during step 1 and will raise an error
    # # if you try demultiplexing without a barcodes file
    else :
        self . paramsdict [ 'barcodes_path' ] = newvalue
# # ` Merged : ` in newvalue for sorted _ fastq _ path indicates that this
# # assembly is a merge of several others , so this param has no
# # value for this assembly
elif param == 'sorted_fastq_path' :
    if newvalue and not "Merged:" in newvalue :
        fullsortedpath = _expander ( newvalue )
        if os . path . isdir ( fullsortedpath ) :
            raise IPyradWarningExit ( SORTED_ISDIR . format ( fullsortedpath ) )
        elif glob . glob ( fullsortedpath ) :
            self . paramsdict [ 'sorted_fastq_path' ] = fullsortedpath
        else :
            raise IPyradWarningExit ( SORTED_NOT_FOUND . format ( fullsortedpath ) )
    # # if no value was entered then set to " " .
    else :
        self . paramsdict [ 'sorted_fastq_path' ] = ""
elif param == 'assembly_method' : # # TEMPORARY BLOCK ON DENOVO + REFERENCE METHOD
# if newvalue = = " denovo + reference " :
# raise IPyradWarningExit ( " " "
# Error : The ' denovo + reference ' method is temporarily blocked while we
# refactor it to greatly improve the speed . You can either revert to an
# older version ( pre v . 0.7.0 ) or wait for the next update to resume using
# this method .
    methods = [ "denovo" , "reference" , "denovo+reference" , "denovo-reference" ]
    assert newvalue in methods , BAD_ASSEMBLY_METHOD . format ( newvalue )
    self . paramsdict [ 'assembly_method' ] = newvalue
elif param == 'reference_sequence' :
    if newvalue :
        fullrawpath = _expander ( newvalue )
        if not os . path . isfile ( fullrawpath ) :
            LOGGER . info ( "reference sequence file not found." )
            raise IPyradWarningExit ( REF_NOT_FOUND . format ( fullrawpath ) )
        self . paramsdict [ 'reference_sequence' ] = fullrawpath
    # # if no value was entered the set to " " . Will be checked again
    # # at step3 if user tries to do refseq and raise error
    else :
        self . paramsdict [ 'reference_sequence' ] = ""
elif param == 'datatype' : # # list of allowed datatypes
    datatypes = [ 'rad' , 'gbs' , 'ddrad' , 'pairddrad' , 'pairgbs' , 'merged' , '2brad' , 'pair3rad' ]
    # # raise error if something else
    if str ( newvalue ) not in datatypes :
        raise IPyradError ( """
    datatype {} not recognized, must be one of: {}
    """ . format ( newvalue , datatypes ) )
    else :
        self . paramsdict [ 'datatype' ] = str ( newvalue )
        # # link _ barcodes is called before datatypes is set
        # # we need to know the datatype so we can read in
        # # the multiplexed barcodes for 3rad . This seems
        # # a little annoying , but it was better than any
        # # alternatives I could think of .
        if "3rad" in self . paramsdict [ 'datatype' ] and not self . paramsdict [ 'sorted_fastq_path' ] . strip ( ) :
            if not "Merged:" in self . paramsdict [ 'barcodes_path' ] :
                self . _link_barcodes ( )
elif param == 'restriction_overhang' :
    newvalue = _tuplecheck ( newvalue , str )
    assert isinstance ( newvalue , tuple ) , """
    cut site must be a tuple, e.g., (TGCAG, '') or (TGCAG, CCGG)"""
    # # Handle the special case where the user has 1
    # # restriction overhang and does not include the trailing comma
    if len ( newvalue ) == 1 : # # for gbs users might not know to enter the second cut site
    # # so we do it for them .
        if self . paramsdict [ "datatype" ] == "gbs" :
            newvalue += newvalue
        else :
            newvalue += ( "" , )
    # newvalue = ( newvalue [ 0 ] , " " )
    # > > > > > d40a5d5086a0d0aace04dd08338ec4ba5341d1f2
    # # Handle 3rad datatype with only 3 cutters
    if len ( newvalue ) == 3 :
        newvalue = ( newvalue [ 0 ] , newvalue [ 1 ] , newvalue [ 2 ] , "" )
    assert len ( newvalue ) <= 4 , """
    most datasets require 1 or 2 cut sites, e.g., (TGCAG, '') or (TGCAG, CCGG).
    For 3rad/seqcap may be up to 4 cut sites."""
    self . paramsdict [ 'restriction_overhang' ] = newvalue
elif param == 'max_low_qual_bases' :
    assert isinstance ( int ( newvalue ) , int ) , """
    max_low_qual_bases must be an integer."""
    self . paramsdict [ 'max_low_qual_bases' ] = int ( newvalue )
elif param == 'phred_Qscore_offset' :
    assert isinstance ( int ( newvalue ) , int ) , "phred_Qscore_offset must be an integer."
    self . paramsdict [ 'phred_Qscore_offset' ] = int ( newvalue )
elif param == 'mindepth_statistical' :
    assert isinstance ( int ( newvalue ) , int ) , "mindepth_statistical must be an integer."
    # # do not allow values below 5
    if int ( newvalue ) < 5 :
        raise IPyradError ( """
    mindepth statistical cannot be set < 5. Use mindepth_majrule.""" )
    else :
        self . paramsdict [ 'mindepth_statistical' ] = int ( newvalue )
elif param == 'mindepth_majrule' :
    assert isinstance ( int ( newvalue ) , int ) , "mindepth_majrule must be an integer."
    self . paramsdict [ 'mindepth_majrule' ] = int ( newvalue )
elif param == 'maxdepth' :
    self . paramsdict [ 'maxdepth' ] = int ( newvalue )
elif param == 'clust_threshold' :
    newvalue = float ( newvalue )
    assert ( newvalue < 1 ) & ( newvalue > 0 ) , "clust_threshold must be a decimal value between 0 and 1."
    self . paramsdict [ 'clust_threshold' ] = newvalue
elif param == 'max_barcode_mismatch' :
    self . paramsdict [ 'max_barcode_mismatch' ] = int ( newvalue )
elif param == 'filter_adapters' :
    self . paramsdict [ 'filter_adapters' ] = int ( newvalue )
elif param == 'filter_min_trim_len' :
    self . paramsdict [ "filter_min_trim_len" ] = int ( newvalue )
elif param == 'max_alleles_consens' :
    self . paramsdict [ 'max_alleles_consens' ] = int ( newvalue )
elif param == 'max_Ns_consens' :
    newvalue = _tuplecheck ( newvalue , int )
    assert isinstance ( newvalue , tuple ) , "max_Ns_consens should be a tuple e.g., (8, 8)"
    self . paramsdict [ 'max_Ns_consens' ] = newvalue
elif param == 'max_Hs_consens' :
    newvalue = _tuplecheck ( newvalue , int )
    assert isinstance ( newvalue , tuple ) , "max_Hs_consens should be a tuple e.g., (5, 5)"
    self . paramsdict [ 'max_Hs_consens' ] = newvalue
elif param == 'min_samples_locus' :
    self . paramsdict [ 'min_samples_locus' ] = int ( newvalue )
elif param == 'max_shared_Hs_locus' :
    if isinstance ( newvalue , str ) :
        if newvalue . isdigit ( ) :
            newvalue = int ( newvalue )
        else :
            try :
                newvalue = float ( newvalue )
            except Exception as inst :
                raise IPyradParamsError ( """
    max_shared_Hs_locus must be int or float, you put: {}""" . format ( newvalue ) )
    self . paramsdict [ 'max_shared_Hs_locus' ] = newvalue
elif param == 'max_SNPs_locus' :
    newvalue = _tuplecheck ( newvalue , int )
    assert isinstance ( newvalue , tuple ) , "max_SNPs_locus should be a tuple e.g., (20, 20)"
    self . paramsdict [ 'max_SNPs_locus' ] = newvalue
elif param == 'max_Indels_locus' :
    newvalue = _tuplecheck ( newvalue , int )
    assert isinstance ( newvalue , tuple ) , "max_Indels_locus should be a tuple e.g., (5, 100)"
    self . paramsdict [ 'max_Indels_locus' ] = newvalue
# # deprecated but retained for legacy , now uses trim _ reads ( below )
elif param == 'edit_cutsites' : # # Force into a string tuple
    newvalue = _tuplecheck ( newvalue )
    # # try converting each tup element to ints
    newvalue = list ( newvalue )
    for i in range ( 2 ) :
        try :
            newvalue [ i ] = int ( newvalue [ i ] )
        except ( ValueError , IndexError ) :
            newvalue . append ( 0 )
            pass
    newvalue = tuple ( newvalue )
    # # make sure we have a nice tuple
    if not isinstance ( newvalue , tuple ) :
        raise IPyradWarningExit ( """
    Error: edit_cutsites should be a tuple e.g., (0, 5) or ('TGCAG', 6),
    you entered {}
    """ . format ( newvalue ) )
    self . paramsdict [ 'edit_cutsites' ] = newvalue
elif param == 'trim_reads' : # # Force into a string tuple
    newvalue = _tuplecheck ( newvalue )
    # # try converting each tup element to ints
    newvalue = list ( newvalue )
    for i in range ( 4 ) :
        try :
            newvalue [ i ] = int ( newvalue [ i ] )
        except ( ValueError , IndexError ) :
            newvalue . append ( 0 )
            pass
    newvalue = tuple ( newvalue )
    # # make sure we have a nice tuple
    if not isinstance ( newvalue , tuple ) :
        raise IPyradWarningExit ( """
    Error: trim_reads should be a tuple e.g., (0, -5, -5, 0)
    or (0, 90, 0, 90), or (0, 0, 0, 0).
    You entered {}\n""" . format ( newvalue ) )
    self . paramsdict [ 'trim_reads' ] = newvalue
# # deprecated but retained for legacy , now named trim _ loci
elif param == 'trim_overhang' :
    newvalue = _tuplecheck ( newvalue , str )
    assert isinstance ( newvalue , tuple ) , "trim_overhang should be a tuple e.g., (4, *, *, 4)"
    self . paramsdict [ 'trim_overhang' ] = tuple ( [ int ( i ) for i in newvalue ] )
elif param == 'trim_loci' :
    newvalue = _tuplecheck ( newvalue , str )
    assert isinstance ( newvalue , tuple ) , "trim_overhang should be a tuple e.g., (0, -5, -5, 0)"
    self . paramsdict [ 'trim_loci' ] = tuple ( [ int ( i ) for i in newvalue ] )
elif param == 'output_formats' : # # let ' s get whatever the user entered as a tuple of letters
    allowed = assemble . write_outfiles . OUTPUT_FORMATS . keys ( )
    # < < < < < HEAD
    # # Handle the case where output formats is an empty string
    if isinstance ( newvalue , str ) : # # strip commas and spaces from string so we have only letters
        newvalue = newvalue . replace ( "," , "" ) . replace ( " " , "" )
        newvalue = list ( newvalue )
        if not newvalue :
            newvalue = [ "*" ]
    if isinstance ( newvalue , tuple ) :
        newvalue = list ( newvalue )
    # if isinstance ( newvalue , tuple ) :
    # newvalue = list ( newvalue )
    # if isinstance ( newvalue , str ) :
    # newvalue = [ i . strip ( ) for i in newvalue . split ( " , " ) ]
    # # # Handle the case where output formats is empty
    # if not any ( newvalue ) :
    # newvalue = " * "
    # > > > > > 488144d1d97240b8b6f6caf9cfb6c023bb6ebb36
    if isinstance ( newvalue , list ) : # # if more than letters , raise an warning
        if any ( [ len ( i ) > 1 for i in newvalue ] ) :
            LOGGER . warning ( """
    'output_formats' params entry is malformed. Setting to * to avoid errors.""" )
            newvalue = allowed
        newvalue = tuple ( newvalue )
        # newvalue = tuple ( [ i for i in newvalue if i in allowed ] )
    if "*" in newvalue :
        newvalue = allowed
    # # set the param
    self . paramsdict [ 'output_formats' ] = newvalue
elif param == 'pop_assign_file' :
    fullpoppath = _expander ( newvalue )
    # # if a path is entered , raise exception if not found
    if newvalue :
        if not os . path . isfile ( fullpoppath ) :
            LOGGER . warn ( "Population assignment file not found." )
            raise IPyradWarningExit ( """
    Warning: Population assignment file not found. This must be an
    absolute path (/home/wat/ipyrad/data/my_popfile.txt) or relative to
    the directory where you're running ipyrad (./data/my_popfile.txt)
    You entered: {}\n""" . format ( fullpoppath ) )
        # # should we add a check here that all pop samples are in samples ?
        self . paramsdict [ 'pop_assign_file' ] = fullpoppath
        self . _link_populations ( )
    else :
        self . paramsdict [ 'pop_assign_file' ] = ""
        # # Don ' t forget to possibly blank the populations dictionary
        self . populations = { }
return self


signature:
def stats ( self ) :
    """Returns a data frame with Sample data and state ."""

implementation:
nameordered = self . samples . keys ( )
nameordered . sort ( )
# # Set pandas to display all samples instead of truncating
pd . options . display . max_rows = len ( self . samples )
statdat = pd . DataFrame ( [ self . samples [ i ] . stats for i in nameordered ] , index = nameordered ) . dropna ( axis = 1 , how = 'all' )
# ensure non h , e columns print as ints
for column in statdat :
    if column not in [ "hetero_est" , "error_est" ] :
        statdat [ column ] = np . nan_to_num ( statdat [ column ] ) . astype ( int )
return statdat


signature:
def files ( self ) :
    """Returns a data frame with Sample files . Not very readable . . ."""

implementation:
nameordered = self . samples . keys ( )
nameordered . sort ( )
# # replace curdir with . for shorter printing
# fullcurdir = os . path . realpath ( os . path . curdir )
return pd . DataFrame ( [ self . samples [ i ] . files for i in nameordered ] , index = nameordered ) . dropna ( axis = 1 , how = 'all' )


signature:
def _build_stat ( self , idx ) :
    """Returns a data frame with Sample stats for each step"""

implementation:
nameordered = self . samples . keys ( )
nameordered . sort ( )
newdat = pd . DataFrame ( [ self . samples [ i ] . stats_dfs [ idx ] for i in nameordered ] , index = nameordered ) . dropna ( axis = 1 , how = 'all' )
return newdat


signature:
def _link_fastqs ( self , path = None , force = False , append = False , splitnames = "_" , fields = None , ipyclient = None ) :
    """Create Sample objects from demultiplexed fastq files in sorted _ fastq _ path ,
 or append additional fastq files to existing Samples . This provides
 more flexible file input through the API than available in step1 of the
 command line interface . If passed ipyclient it will run in parallel .
 Note
 This function is called during step 1 if files are specified in
 ' sorted _ fastq _ path ' .
 Parameters
 path : str
 Path to the fastq files to be linked to Sample objects . The default
 location is to select all files in the ' sorted _ fastq _ path ' .
 Alternatively a different path can be entered here .
 append : bool
 The default action is to overwrite fastq files linked to Samples if
 they already have linked files . Use append = True to instead append
 additional fastq files to a Sample ( file names should be formatted
 the same as usual , e . g . , [ name ] _ R1 _ [ optional ] . fastq . gz ) .
 splitnames : str
 A string character used to file names . In combination with the
 fields argument can be used to subselect filename fields names .
 fields : list
 A list of indices for the fields to be included in names after
 filnames are split on the splitnames character . Useful for appending
 sequence names which must match existing names . If the largest index
 is greater than the number of split strings in the name the index
 if ignored . e . g . , [ 2,3,4 ] # # excludes 0 , 1 and > 4
 force : bool
 Overwrites existing Sample data and statistics .
 Returns
 str
 Prints the number of new Sample objects created and the number of
 fastq files linked to Sample objects in the Assembly object ."""

implementation:
# # cannot both force and append at once
if force and append :
    raise IPyradError ( "Cannot use force and append at the same time." )
if self . samples and not ( force or append ) :
    raise IPyradError ( "Files already linked to `{}`." . format ( self . name ) + " Use force=True to replace all files, or append=True to add" + " additional files to existing Samples." )
# # make sure there is a workdir and workdir / fastqdir
self . dirs . fastqs = os . path . join ( self . paramsdict [ "project_dir" ] , self . name + "_fastqs" )
if not os . path . exists ( self . paramsdict [ "project_dir" ] ) :
    os . mkdir ( self . paramsdict [ "project_dir" ] )
# # get path to data files
if not path :
    path = self . paramsdict [ "sorted_fastq_path" ]
# # but grab fastq / fq / gz , and then sort
fastqs = glob . glob ( path )
# # Assert files are not . bz2 format
if any ( [ i for i in fastqs if i . endswith ( ".bz2" ) ] ) :
    raise IPyradError ( NO_SUPPORT_FOR_BZ2 . format ( path ) )
fastqs = [ i for i in fastqs if i . endswith ( ".gz" ) or i . endswith ( ".fastq" ) or i . endswith ( ".fq" ) ]
fastqs . sort ( )
LOGGER . debug ( "Linking these fastq files:\n{}" . format ( fastqs ) )
# # raise error if no files are found
if not fastqs :
    raise IPyradError ( NO_FILES_FOUND_PAIRS . format ( self . paramsdict [ "sorted_fastq_path" ] ) )
# # link pairs into tuples
if 'pair' in self . paramsdict [ "datatype" ] : # # check that names fit the paired naming convention
# # trying to support flexible types ( _ R2 _ , _ 2 . fastq )
    r1_try1 = [ i for i in fastqs if "_R1_" in i ]
    r1_try2 = [ i for i in fastqs if i . endswith ( "_1.fastq.gz" ) ]
    r1_try3 = [ i for i in fastqs if i . endswith ( "_R1.fastq.gz" ) ]
    r2_try1 = [ i for i in fastqs if "_R2_" in i ]
    r2_try2 = [ i for i in fastqs if i . endswith ( "_2.fastq.gz" ) ]
    r2_try3 = [ i for i in fastqs if i . endswith ( "_R2.fastq.gz" ) ]
    r1s = [ r1_try1 , r1_try2 , r1_try3 ]
    r2s = [ r2_try1 , r2_try2 , r2_try3 ]
    # # check that something was found
    if not r1_try1 + r1_try2 + r1_try3 :
        raise IPyradWarningExit ( "Paired filenames are improperly formatted. See Documentation" )
    if not r2_try1 + r2_try2 + r2_try3 :
        raise IPyradWarningExit ( "Paired filenames are improperly formatted. See Documentation" )
    # # find the one with the right number of R1s
    for idx , tri in enumerate ( r1s ) :
        if len ( tri ) == len ( fastqs ) / 2 :
            break
    r1_files = r1s [ idx ]
    r2_files = r2s [ idx ]
    if len ( r1_files ) != len ( r2_files ) :
        raise IPyradWarningExit ( R1_R2_name_error . format ( len ( r1_files ) , len ( r2_files ) ) )
    fastqs = [ ( i , j ) for i , j in zip ( r1_files , r2_files ) ]
# # data are not paired , create empty tuple pair
else : # # print warning if _ R2 _ is in names when not paired
    idx = 0
    if any ( [ "_R2_" in i for i in fastqs ] ) :
        print ( NAMES_LOOK_PAIRED_WARNING )
    fastqs = [ ( i , "" ) for i in fastqs ]
# # counters for the printed output
linked = 0
appended = 0
# # clear samples if force
if force :
    self . samples = { }
# # track parallel jobs
linkjobs = { }
if ipyclient :
    lbview = ipyclient . load_balanced_view ( )
# # iterate over input files
for fastqtuple in list ( fastqs ) :
    assert isinstance ( fastqtuple , tuple ) , "fastqs not a tuple."
    # # local counters
    createdinc = 0
    linkedinc = 0
    appendinc = 0
    # # remove file extension from name
    if idx == 0 :
        sname = _name_from_file ( fastqtuple [ 0 ] , splitnames , fields )
    elif idx == 1 :
        sname = os . path . basename ( fastqtuple [ 0 ] . rsplit ( "_1.fastq.gz" , 1 ) [ 0 ] )
    elif idx == 2 :
        sname = os . path . basename ( fastqtuple [ 0 ] . rsplit ( "_R1.fastq.gz" , 1 ) [ 0 ] )
    LOGGER . debug ( "New Sample name {}" . format ( sname ) )
    if sname not in self . samples : # # create new Sample
        LOGGER . debug ( "Creating new sample - " . format ( sname ) )
        self . samples [ sname ] = Sample ( sname )
        self . samples [ sname ] . stats . state = 1
        self . samples [ sname ] . barcode = None
        self . samples [ sname ] . files . fastqs . append ( fastqtuple )
        createdinc += 1
        linkedinc += 1
    else : # # if not forcing , shouldn ' t be here with existing Samples
        if append : # if fastqtuple not in self . samples [ sname ] . files . fastqs :
            self . samples [ sname ] . files . fastqs . append ( fastqtuple )
            appendinc += 1
        elif force : # # overwrite / create new Sample
            LOGGER . debug ( "Overwriting sample - " . format ( sname ) )
            self . samples [ sname ] = Sample ( sname )
            self . samples [ sname ] . stats . state = 1
            self . samples [ sname ] . barcode = None
            self . samples [ sname ] . files . fastqs . append ( fastqtuple )
            createdinc += 1
            linkedinc += 1
        else :
            print ( """
        The files {} are already in Sample. Use append=True to append additional
        files to a Sample or force=True to replace all existing Samples.
        """ . format ( sname ) )
    # # support serial execution w / o ipyclient
    if not ipyclient :
        if any ( [ linkedinc , createdinc , appendinc ] ) :
            gzipped = bool ( fastqtuple [ 0 ] . endswith ( ".gz" ) )
            nreads = 0
            for alltuples in self . samples [ sname ] . files . fastqs :
                nreads += _zbufcountlines ( alltuples [ 0 ] , gzipped )
            self . samples [ sname ] . stats . reads_raw = nreads / 4
            self . samples [ sname ] . stats_dfs . s1 [ "reads_raw" ] = nreads / 4
            self . samples [ sname ] . state = 1
            LOGGER . debug ( "Got reads for sample - {} {}" . format ( sname , self . samples [ sname ] . stats . reads_raw ) )
            # created + = createdinc
            linked += linkedinc
            appended += appendinc
    # # do counting in parallel
    else :
        if any ( [ linkedinc , createdinc , appendinc ] ) :
            gzipped = bool ( fastqtuple [ 0 ] . endswith ( ".gz" ) )
            for sidx , tup in enumerate ( self . samples [ sname ] . files . fastqs ) :
                key = sname + "_{}" . format ( sidx )
                linkjobs [ key ] = lbview . apply ( _bufcountlines , * ( tup [ 0 ] , gzipped ) )
            LOGGER . debug ( "sent count job for {}" . format ( sname ) )
            # created + = createdinc
            linked += linkedinc
            appended += appendinc
# # wait for link jobs to finish if parallel
if ipyclient :
    start = time . time ( )
    printstr = ' loading reads         | {} | s1 |'
    while 1 :
        fin = [ i . ready ( ) for i in linkjobs . values ( ) ]
        elapsed = datetime . timedelta ( seconds = int ( time . time ( ) - start ) )
        progressbar ( len ( fin ) , sum ( fin ) , printstr . format ( elapsed ) , spacer = self . _spacer )
        time . sleep ( 0.1 )
        if len ( fin ) == sum ( fin ) :
            print ( "" )
            break
    # # collect link job results
    sampdict = { i : 0 for i in self . samples }
    for result in linkjobs :
        sname = result . rsplit ( "_" , 1 ) [ 0 ]
        nreads = linkjobs [ result ] . result ( )
        sampdict [ sname ] += nreads
    for sname in sampdict :
        self . samples [ sname ] . stats . reads_raw = sampdict [ sname ] / 4
        self . samples [ sname ] . stats_dfs . s1 [ "reads_raw" ] = sampdict [ sname ] / 4
        self . samples [ sname ] . state = 1
# # print if data were linked
# print ( " { } new Samples created in ' { } ' . " . format ( created , self . name ) )
if linked : # # double for paired data
    if 'pair' in self . paramsdict [ "datatype" ] :
        linked = linked * 2
    if self . _headers :
        print ( "{}{} fastq files loaded to {} Samples." . format ( self . _spacer , linked , len ( self . samples ) ) )
    # # save the location where these files are located
    self . dirs . fastqs = os . path . realpath ( os . path . dirname ( path ) )
if appended :
    if self . _headers :
        print ( "{}{} fastq files appended to {} existing Samples." . format ( self . _spacer , appended , len ( self . samples ) ) )
# # save step - 1 stats . We don ' t want to write this to the fastq dir , b / c
# # it is not necessarily inside our project dir . Instead , we ' ll write
# # this file into our project dir in the case of linked _ fastqs .
self . stats_dfs . s1 = self . _build_stat ( "s1" )
self . stats_files . s1 = os . path . join ( self . paramsdict [ "project_dir" ] , self . name + '_s1_demultiplex_stats.txt' )
with open ( self . stats_files . s1 , 'w' ) as outfile :
    self . stats_dfs . s1 . fillna ( value = 0 ) . astype ( np . int ) . to_string ( outfile )


signature:
def _link_barcodes ( self ) :
    """Private function . Links Sample barcodes in a dictionary as
 [ Assembly ] . barcodes , with barcodes parsed from the ' barcodes _ path '
 parameter . This function is called during set _ params ( ) when setting
 the barcodes _ path ."""

implementation:
# # parse barcodefile
try : # # allows fuzzy match to barcodefile name
    barcodefile = glob . glob ( self . paramsdict [ "barcodes_path" ] ) [ 0 ]
    # # read in the file
    bdf = pd . read_csv ( barcodefile , header = None , delim_whitespace = 1 , dtype = str )
    bdf = bdf . dropna ( )
    # # make sure bars are upper case
    bdf [ 1 ] = bdf [ 1 ] . str . upper ( )
    # # if replicates are present then print a warning
    reps = bdf [ 0 ] . unique ( ) . shape [ 0 ] != bdf [ 0 ] . shape [ 0 ]
    if reps :
        print ( "{spacer}Warning: technical replicates (same name) will be combined." . format ( ** { 'spacer' : self . _spacer } ) )
        # # add - technical - replicate - N to replicate names
        reps = [ i for i in bdf [ 0 ] if list ( bdf [ 0 ] ) . count ( i ) > 1 ]
        ureps = list ( set ( reps ) )
        for name in ureps :
            idxs = bdf [ bdf [ 0 ] == ureps [ 0 ] ] . index . tolist ( )
            for num , idx in enumerate ( idxs ) :
                bdf . ix [ idx ] [ 0 ] = bdf . ix [ idx ] [ 0 ] + "-technical-replicate-" + str ( num + 1 )
    # # make sure chars are all proper
    if not all ( bdf [ 1 ] . apply ( set ( "RKSYWMCATG" ) . issuperset ) ) :
        LOGGER . warn ( BAD_BARCODE )
        raise IPyradError ( BAD_BARCODE )
    # # 3rad / seqcap use multiplexed barcodes
    # # We ' ll concatenate them with a plus and split them later
    if "3rad" in self . paramsdict [ "datatype" ] :
        try :
            bdf [ 2 ] = bdf [ 2 ] . str . upper ( )
            self . barcodes = dict ( zip ( bdf [ 0 ] , bdf [ 1 ] + "+" + bdf [ 2 ] ) )
        except KeyError as inst :
            msg = "    3rad assumes multiplexed barcodes. Doublecheck your barcodes file."
            LOGGER . error ( msg )
            raise IPyradError ( msg )
    else : # # set attribute on Assembly object
        self . barcodes = dict ( zip ( bdf [ 0 ] , bdf [ 1 ] ) )
except ( IOError , IndexError ) :
    raise IPyradWarningExit ( "    Barcodes file not found. You entered: {}" . format ( self . paramsdict [ "barcodes_path" ] ) )
except ValueError as inst :
    msg = "    Barcodes file format error."
    LOGGER . warn ( msg )
    raise IPyradError ( inst )


signature:
def _link_populations ( self , popdict = None , popmins = None ) :
    """Creates self . populations dictionary to save mappings of individuals to
 populations / sites , and checks that individual names match with Samples .
 The self . populations dict keys are pop names and the values are lists
 of length 2 . The first element is the min number of samples per pop
 for final filtering of loci , and the second element is the list of
 samples per pop .
 Population assigments are used for heirarchical clustering , for
 generating summary stats , and for outputing some file types ( . treemix
 for example ) . Internally stored as a dictionary .
 Note
 By default a File is read in from ` pop _ assign _ file ` with one individual
 per line and space separated pairs of ind pop :
 ind1 pop1
 ind2 pop2
 ind3 pop3
 etc . . .
 Parameters
 TODO : NB : Using API and passing in popdict and popmins is currently
 unimplemented , or at least looks like it doesn ' t work . Leaving
 these docs cuz Deren might have ideas about it being useful .
 popdict : dict
 When using the API it may be easier to simply create a dictionary
 to pass in as an argument instead of reading from an input file .
 This can be done with the ` popdict ` argument like below :
 pops = { ' pop1 ' : [ ' ind1 ' , ' ind2 ' , ' ind3 ' ] , ' pop2 ' : [ ' ind4 ' , ' ind5 ' ] }
 [ Assembly ] . _ link _ populations ( popdict = pops ) .
 popmins : dict
 If you want to apply a minsamples filter based on populations
 you can add a popmins dictionary . This indicates the number of
 samples in each population that must be present in a locus for
 the locus to be retained . Example :
 popmins = { ' pop1 ' : 3 , ' pop2 ' : 2}"""

implementation:
if not popdict : # # glob it in case of fuzzy matching
    popfile = glob . glob ( self . paramsdict [ "pop_assign_file" ] ) [ 0 ]
    if not os . path . exists ( popfile ) :
        raise IPyradError ( "Population assignment file not found: {}" . format ( self . paramsdict [ "pop_assign_file" ] ) )
    try : # # parse populations file
        popdat = pd . read_csv ( popfile , header = None , delim_whitespace = 1 , names = [ "inds" , "pops" ] , comment = "#" )
        popdict = { key : group . inds . values . tolist ( ) for key , group in popdat . groupby ( "pops" ) }
        # # parse minsamples per population if present ( line with # )
        mindat = [ i . lstrip ( "#" ) . lstrip ( ) . rstrip ( ) for i in open ( popfile , 'r' ) . readlines ( ) if i . startswith ( "#" ) ]
        if mindat :
            popmins = { }
            for i in range ( len ( mindat ) ) :
                minlist = mindat [ i ] . replace ( "," , "" ) . split ( )
                popmins . update ( { i . split ( ':' ) [ 0 ] : int ( i . split ( ':' ) [ 1 ] ) for i in minlist } )
        else :
            raise IPyradError ( MIN_SAMPLES_PER_POP_MALFORMED )
    except ( ValueError , IOError ) :
        LOGGER . warn ( "Populations file may be malformed." )
        raise IPyradError ( MIN_SAMPLES_PER_POP_MALFORMED )
else : # # pop dict is provided by user
    pass
# # check popdict . Filter for bad samples
# # Warn user but don ' t bail out , could be setting the pops file
# # on a new assembly w / o any linked samples .
badsamples = [ i for i in itertools . chain ( * popdict . values ( ) ) if i not in self . samples . keys ( ) ]
if any ( badsamples ) :
    LOGGER . warn ( "Some names from population input do not match Sample " + "names: " . format ( ", " . join ( badsamples ) ) )
    LOGGER . warn ( "If this is a new assembly this is normal." )
# # If popmins not set , just assume all mins are zero
if not popmins :
    popmins = { i : 0 for i in popdict . keys ( ) }
# # check popmins
# # cannot have higher min for a pop than there are samples in the pop
popmax = { i : len ( popdict [ i ] ) for i in popdict }
if not all ( [ popmax [ i ] >= popmins [ i ] for i in popdict ] ) :
    raise IPyradWarningExit ( " minsample per pop value cannot be greater than the " + " number of samples in the pop. Modify the populations file." )
# # return dict
self . populations = { i : ( popmins [ i ] , popdict [ i ] ) for i in popdict }


signature:
def get_params ( self , param = "" ) :
    """pretty prints params if called as a function"""

implementation:
fullcurdir = os . path . realpath ( os . path . curdir )
if not param :
    for index , ( key , value ) in enumerate ( self . paramsdict . items ( ) ) :
        if isinstance ( value , str ) :
            value = value . replace ( fullcurdir + "/" , "./" )
        sys . stdout . write ( "{}{:<4}{:<28}{:<45}\n" . format ( self . _spacer , index , key , value ) )
else :
    try :
        if int ( param ) : # sys . stdout . write ( self . paramsdict . values ( ) [ int ( param ) - 1 ] )
            return self . paramsdict . values ( ) [ int ( param ) ]
    except ( ValueError , TypeError , NameError , IndexError ) :
        try :
            return self . paramsdict [ param ]
        except KeyError :
            return 'key not recognized'


signature:
def set_params ( self , param , newvalue ) :
    """Set a parameter to a new value . Raises error if newvalue is wrong type .
 Note
 Use [ Assembly ] . get _ params ( ) to see the parameter values currently
 linked to the Assembly object .
 Parameters
 param : int or str
 The index ( e . g . , 1 ) or string name ( e . g . , " project _ dir " )
 for the parameter that will be changed .
 newvalue : int , str , or tuple
 The new value for the parameter selected for ` param ` . Use
 ` ipyrad . get _ params _ info ( ) ` to get further information about
 a given parameter . If the wrong type is entered for newvalue
 ( e . g . , a str when it should be an int ) , an error will be raised .
 Further information about each parameter is also available
 in the documentation .
 Examples
 # # param ' project _ dir ' takes only a str as input
 [ Assembly ] . set _ params ( ' project _ dir ' , ' new _ directory ' )
 # # param ' restriction _ overhang ' must be a tuple or str , if str it is
 # # converted to a tuple with the second entry empty .
 [ Assembly ] . set _ params ( ' restriction _ overhang ' , ( ' CTGCAG ' , ' CCGG ' )
 # # param ' max _ shared _ Hs _ locus ' can be an int or a float :
 [ Assembly ] . set _ params ( ' max _ shared _ Hs _ locus ' , 0.25)"""

implementation:
# # this includes current params and some legacy params for conversion
legacy_params = [ "edit_cutsites" , "trim_overhang" ]
current_params = self . paramsdict . keys ( )
allowed_params = current_params + legacy_params
# # require parameter recognition
# if not ( ( param in range ( 50 ) ) or \
# ( param in [ str ( i ) for i in range ( 50 ) ] ) or \
# ( param in allowed _ params ) ) :
if not param in allowed_params :
    raise IPyradParamsError ( "Parameter key not recognized: {}" . format ( param ) )
# # make string
param = str ( param )
# # get index if param is keyword arg ( this index is now zero based ! )
if len ( param ) < 3 :
    param = self . paramsdict . keys ( ) [ int ( param ) ]
# # run assertions on new param
try :
    self = _paramschecker ( self , param , newvalue )
except Exception as inst :
    raise IPyradWarningExit ( BAD_PARAMETER . format ( param , inst , newvalue ) )


signature:
def write_params ( self , outfile = None , force = False ) :
    """Write out the parameters of this assembly to a file properly
 formatted as input for ` ipyrad - p < params . txt > ` . A good and
 simple way to share / archive parameter settings for assemblies .
 This is also the function that ' s used by _ _ main _ _ to
 generate default params . txt files for ` ipyrad - n `"""

implementation:
if outfile is None :
    outfile = "params-" + self . name + ".txt"
# # Test if params file already exists ?
# # If not forcing , test for file and bail out if it exists
if not force :
    if os . path . isfile ( outfile ) :
        raise IPyradWarningExit ( PARAMS_EXISTS . format ( outfile ) )
with open ( outfile , 'w' ) as paramsfile : # # Write the header . Format to 80 columns
    header = "------- ipyrad params file (v.{})" . format ( ip . __version__ )
    header += ( "-" * ( 80 - len ( header ) ) )
    paramsfile . write ( header )
    # # Whip through the current paramsdict and write out the current
    # # param value , the ordered dict index number . Also ,
    # # get the short description from paramsinfo . Make it look pretty ,
    # # pad nicely if at all possible .
    for key , val in self . paramsdict . iteritems ( ) : # # If multiple elements , write them out comma separated
        if isinstance ( val , list ) or isinstance ( val , tuple ) :
            paramvalue = ", " . join ( [ str ( i ) for i in val ] )
        else :
            paramvalue = str ( val )
        # # skip deprecated params
        if key in [ "edit_cutsites" , "trim_overhang" ] :
            continue
        padding = ( " " * ( 30 - len ( paramvalue ) ) )
        paramkey = self . paramsdict . keys ( ) . index ( key )
        paramindex = " ## [{}] " . format ( paramkey )
        LOGGER . debug ( key , val , paramindex )
        name = "[{}]: " . format ( paramname ( paramkey ) )
        description = paraminfo ( paramkey , short = True )
        paramsfile . write ( "\n" + paramvalue + padding + paramindex + name + description )


signature:
def branch ( self , newname , subsamples = None , infile = None ) :
    """Returns a copy of the Assembly object . Does not allow Assembly
 object names to be replicated in namespace or path ."""

implementation:
# # subsample by removal or keeping .
remove = 0
# # is there a better way to ask if it already exists ?
if ( newname == self . name or os . path . exists ( os . path . join ( self . paramsdict [ "project_dir" ] , newname + ".assembly" ) ) ) :
    print ( "{}Assembly object named {} already exists" . format ( self . _spacer , newname ) )
else : # # Make sure the new name doesn ' t have any wacky characters
    self . _check_name ( newname )
    # # Bozo - check . Carve off ' params - ' if it ' s in the new name .
    if newname . startswith ( "params-" ) :
        newname = newname . split ( "params-" ) [ 1 ]
    # # create a copy of the Assembly obj
    newobj = copy . deepcopy ( self )
    newobj . name = newname
    newobj . paramsdict [ "assembly_name" ] = newname
    if subsamples and infile :
        print ( BRANCH_NAMES_AND_INPUT )
    if infile :
        if infile [ 0 ] == "-" :
            remove = 1
            infile = infile [ 1 : ]
        if os . path . exists ( infile ) :
            subsamples = _read_sample_names ( infile )
    # # if remove then swap the samples
    if remove :
        subsamples = list ( set ( self . samples . keys ( ) ) - set ( subsamples ) )
    # # create copies of each subsampled Sample obj
    if subsamples :
        for sname in subsamples :
            if sname in self . samples :
                newobj . samples [ sname ] = copy . deepcopy ( self . samples [ sname ] )
            else :
                print ( "Sample name not found: {}" . format ( sname ) )
        # # reload sample dict w / o non subsamples
        newobj . samples = { name : sample for name , sample in newobj . samples . items ( ) if name in subsamples }
    # # create copies of each subsampled Sample obj
    else :
        for sample in self . samples :
            newobj . samples [ sample ] = copy . deepcopy ( self . samples [ sample ] )
    # # save json of new obj and return object
    newobj . save ( )
    return newobj


signature:
def _step1func ( self , force , ipyclient ) :
    """hidden wrapped function to start step 1"""

implementation:
# # check input data files
sfiles = self . paramsdict [ "sorted_fastq_path" ]
rfiles = self . paramsdict [ "raw_fastq_path" ]
# # do not allow both a sorted _ fastq _ path and a raw _ fastq
if sfiles and rfiles :
    raise IPyradWarningExit ( NOT_TWO_PATHS )
# # but also require that at least one exists
if not ( sfiles or rfiles ) :
    raise IPyradWarningExit ( NO_SEQ_PATH_FOUND )
# # print headers
if self . _headers :
    if sfiles :
        print ( "\n{}Step 1: Loading sorted fastq data to Samples" . format ( self . _spacer ) )
    else :
        print ( "\n{}Step 1: Demultiplexing fastq data to Samples" . format ( self . _spacer ) )
# # if Samples already exist then no demultiplexing
if self . samples :
    if not force :
        print ( SAMPLES_EXIST . format ( len ( self . samples ) , self . name ) )
    else : # # overwrite existing data else do demux
        if glob . glob ( sfiles ) :
            self . _link_fastqs ( ipyclient = ipyclient , force = force )
        else :
            assemble . demultiplex . run2 ( self , ipyclient , force )
# # Creating new Samples
else : # # first check if demultiplexed files exist in sorted path
    if glob . glob ( sfiles ) :
        self . _link_fastqs ( ipyclient = ipyclient )
    # # otherwise do the demultiplexing
    else :
        assemble . demultiplex . run2 ( self , ipyclient , force )


signature:
def _step2func ( self , samples , force , ipyclient ) :
    """hidden wrapped function to start step 2"""

implementation:
# # print header
if self . _headers :
    print ( "\n  Step 2: Filtering reads " )
# # If no samples in this assembly then it means you skipped step1,
if not self . samples . keys ( ) :
    raise IPyradWarningExit ( FIRST_RUN_1 )
# # Get sample objects from list of strings , if API .
samples = _get_samples ( self , samples )
if not force : # # print warning and skip if all are finished
    if all ( [ i . stats . state >= 2 for i in samples ] ) :
        print ( EDITS_EXIST . format ( len ( samples ) ) )
        return
# # Run samples through rawedit
assemble . rawedit . run2 ( self , samples , force , ipyclient )


signature:
def _step3func ( self , samples , noreverse , maxindels , force , ipyclient ) :
    """hidden wrapped function to start step 3"""

implementation:
# # print headers
if self . _headers :
    print ( "\n  Step 3: Clustering/Mapping reads" )
# # Require reference seq for reference - based methods
if self . paramsdict [ 'assembly_method' ] != "denovo" :
    if not self . paramsdict [ 'reference_sequence' ] :
        raise IPyradError ( REQUIRE_REFERENCE_PATH . format ( self . paramsdict [ "assembly_method" ] ) )
    else : # # index the reference sequence
    # # Allow force to reindex the reference sequence
    # # send to run on the cluster .
        lbview = ipyclient . load_balanced_view ( )
        async = lbview . apply ( index_reference_sequence , * ( self , force ) )
        # # print a progress bar for the indexing
        start = time . time ( )
        while 1 :
            elapsed = datetime . timedelta ( seconds = int ( time . time ( ) - start ) )
            printstr = " {}    | {} | s3 |" . format ( "indexing reference" , elapsed )
            finished = int ( async . ready ( ) )
            progressbar ( 1 , finished , printstr , spacer = self . _spacer )
            if finished :
                print ( "" )
                break
            time . sleep ( 0.9 )
        # # error check
        if not async . successful ( ) :
            raise IPyradWarningExit ( async . result ( ) )
# # Get sample objects from list of strings
samples = _get_samples ( self , samples )
# # Check if all / none in the right state
if not self . _samples_precheck ( samples , 3 , force ) :
    raise IPyradError ( FIRST_RUN_2 )
elif not force : # # skip if all are finished
    if all ( [ i . stats . state >= 3 for i in samples ] ) :
        print ( CLUSTERS_EXIST . format ( len ( samples ) ) )
        return
# # run the step function
assemble . cluster_within . run ( self , samples , noreverse , maxindels , force , ipyclient )


signature:
def _step4func ( self , samples , force , ipyclient ) :
    """hidden wrapped function to start step 4"""

implementation:
if self . _headers :
    print ( "\n  Step 4: Joint estimation of error rate and heterozygosity" )
# # Get sample objects from list of strings
samples = _get_samples ( self , samples )
# # Check if all / none in the right state
if not self . _samples_precheck ( samples , 4 , force ) :
    raise IPyradError ( FIRST_RUN_3 )
elif not force : # # skip if all are finished
    if all ( [ i . stats . state >= 4 for i in samples ] ) :
        print ( JOINTS_EXIST . format ( len ( samples ) ) )
        return
# # send to function
assemble . jointestimate . run ( self , samples , force , ipyclient )


signature:
def _step5func ( self , samples , force , ipyclient ) :
    """hidden wrapped function to start step 5"""

implementation:
# # print header
if self . _headers :
    print ( "\n  Step 5: Consensus base calling " )
# # Get sample objects from list of strings
samples = _get_samples ( self , samples )
# # Check if all / none in the right state
if not self . _samples_precheck ( samples , 5 , force ) :
    raise IPyradError ( FIRST_RUN_4 )
elif not force : # # skip if all are finished
    if all ( [ i . stats . state >= 5 for i in samples ] ) :
        print ( CONSENS_EXIST . format ( len ( samples ) ) )
        return
# # pass samples to rawedit
assemble . consens_se . run ( self , samples , force , ipyclient )


signature:
def _step6func ( self , samples , noreverse , force , randomseed , ipyclient , ** kwargs ) :
    """Hidden function to start Step 6."""

implementation:
# # Get sample objects from list of strings
samples = _get_samples ( self , samples )
# # remove samples that aren ' t ready
csamples = self . _samples_precheck ( samples , 6 , force )
# # print CLI header
if self . _headers :
    print ( "\n  Step 6: Clustering at {} similarity across {} samples" . format ( self . paramsdict [ "clust_threshold" ] , len ( csamples ) ) )
# # Check if all / none in the right state
if not csamples :
    raise IPyradError ( FIRST_RUN_5 )
elif not force : # # skip if all are finished
    if all ( [ i . stats . state >= 6 for i in csamples ] ) :
        print ( DATABASE_EXISTS . format ( len ( samples ) ) )
        return
# # run if this point is reached . We no longer check for existing
# # h5 file , since checking Sample states should suffice .
assemble . cluster_across . run ( self , csamples , noreverse , force , randomseed , ipyclient , ** kwargs )


signature:
def _step7func ( self , samples , force , ipyclient ) :
    """Step 7 : Filter and write output files"""

implementation:
# # Get sample objects from list of strings
samples = _get_samples ( self , samples )
if self . _headers :
    print ( "\n  Step 7: Filter and write output files for {} Samples" . format ( len ( samples ) ) )
# # Check if all / none of the samples are in the self . database
try :
    with h5py . File ( self . clust_database , 'r' ) as io5 :
        dbset = set ( io5 [ "seqs" ] . attrs [ 'samples' ] )
        iset = set ( [ i . name for i in samples ] )
        # # TODO : Handle the case where dbdiff is not empty ?
        # # This might arise if someone tries to branch and remove
        # # samples at step 7.
        dbdiff = dbset . difference ( iset )
        idiff = iset . difference ( dbset )
        if idiff :
            print ( NOT_CLUSTERED_YET . format ( self . database , ", " . join ( list ( idiff ) ) ) )
            # # The the old way that failed unless all samples were
            # # clustered successfully in step 6 . Adding some flexibility
            # # to allow writing output even if some samples failed .
            # # raise IPyradWarningExit ( msg )
            # # Remove the samples that aren ' t ready for writing out
            # # i . e . only proceed with the samples that are actually
            # # present in the db
            samples = [ x for x in samples if x . name not in idiff ]
except ( IOError , ValueError ) :
    raise IPyradError ( FIRST_RUN_6 . format ( self . clust_database ) )
if not force :
    outdir = os . path . join ( self . dirs . project , self . name + "_outfiles" )
    if os . path . exists ( outdir ) :
        raise IPyradWarningExit ( OUTPUT_EXISTS . format ( outdir ) )
# # Run step7
assemble . write_outfiles . run ( self , samples , force , ipyclient )


signature:
def _samples_precheck ( self , samples , mystep , force ) :
    """Return a list of samples that are actually ready for the next step .
 Each step runs this prior to calling run , makes it easier to
 centralize and normalize how each step is checking sample states .
 mystep is the state produced by the current step ."""

implementation:
subsample = [ ]
# # filter by state
for sample in samples :
    if sample . stats . state < mystep - 1 :
        LOGGER . debug ( "Sample {} not in proper state." . format ( sample . name ) )
    else :
        subsample . append ( sample )
return subsample


signature:
def _compatible_params_check ( self ) :
    """check for mindepths after all params are set , b / c doing it while each
 is being set becomes complicated"""

implementation:
# # do not allow statistical < majrule
val1 = self . paramsdict [ "mindepth_statistical" ]
val2 = self . paramsdict [ 'mindepth_majrule' ]
if val1 < val2 :
    msg = """
    Warning: mindepth_statistical cannot not be < mindepth_majrule.
    Forcing mindepth_majrule = mindepth_statistical = {}
    """ . format ( val1 )
    LOGGER . warning ( msg )
    print ( msg )
    self . paramsdict [ "mindepth_majrule" ] = val1


signature:
def run ( self , steps = 0 , force = False , ipyclient = None , show_cluster = 0 , ** kwargs ) :
    """Run assembly steps of an ipyrad analysis . Enter steps as a string ,
 e . g . , " 1 " , " 123 " , " 12345 " . This step checks for an existing
 ipcluster instance otherwise it raises an exception . The ipyparallel
 connection is made using information from the _ ipcluster dict of the
 Assembly class object ."""

implementation:
# # check that mindepth params are compatible , fix and report warning .
self . _compatible_params_check ( )
# # wrap everything in a try statement to ensure that we save the
# # Assembly object if it is interrupted at any point , and also
# # to ensure proper cleanup of the ipyclient .
inst = None
try : # # use an existing ipcluster instance
    if not ipyclient :
        args = self . _ipcluster . items ( ) + [ ( "spacer" , self . _spacer ) ]
        ipyclient = ip . core . parallel . get_client ( ** dict ( args ) )
    # # print a message about the cluster status
    # # if MPI setup then we are going to wait until all engines are
    # # ready so that we can print how many cores started on each
    # # host machine exactly .
    if ( self . _cli ) or show_cluster :
        ip . cluster_info ( ipyclient = ipyclient , spacer = self . _spacer )
    # # get the list of steps to run
    if isinstance ( steps , int ) :
        steps = str ( steps )
    steps = sorted ( list ( steps ) )
    # # print an Assembly name header if inside API
    if not self . _cli :
        print ( "Assembly: {}" . format ( self . name ) )
    # # store ipyclient engine pids to the Assembly so we can
    # # hard - interrupt them later if assembly is interrupted .
    # # Only stores pids of engines that aren ' t busy at this moment ,
    # # otherwise it would block here while waiting to find their pids .
    self . _ipcluster [ "pids" ] = { }
    for eid in ipyclient . ids :
        engine = ipyclient [ eid ]
        if not engine . outstanding :
            pid = engine . apply ( os . getpid ) . get ( )
            self . _ipcluster [ "pids" ] [ eid ] = pid
    # ipyclient [ : ] . apply ( os . getpid ) . get _ dict ( )
    # # has many fixed arguments right now , but we may add these to
    # # hackerz _ only , or they may be accessed in the API .
    if '1' in steps :
        self . _step1func ( force , ipyclient )
        self . save ( )
        ipyclient . purge_everything ( )
    if '2' in steps :
        self . _step2func ( samples = None , force = force , ipyclient = ipyclient )
        self . save ( )
        ipyclient . purge_everything ( )
    if '3' in steps :
        self . _step3func ( samples = None , noreverse = 0 , force = force , maxindels = 8 , ipyclient = ipyclient )
        self . save ( )
        ipyclient . purge_everything ( )
    if '4' in steps :
        self . _step4func ( samples = None , force = force , ipyclient = ipyclient )
        self . save ( )
        ipyclient . purge_everything ( )
    if '5' in steps :
        self . _step5func ( samples = None , force = force , ipyclient = ipyclient )
        self . save ( )
        ipyclient . purge_everything ( )
    if '6' in steps :
        self . _step6func ( samples = None , noreverse = 0 , randomseed = 12345 , force = force , ipyclient = ipyclient , ** kwargs )
        self . save ( )
        ipyclient . purge_everything ( )
    if '7' in steps :
        self . _step7func ( samples = None , force = force , ipyclient = ipyclient )
        self . save ( )
        ipyclient . purge_everything ( )
# # handle exceptions so they will be raised after we clean up below
except KeyboardInterrupt as inst :
    print ( "\n  Keyboard Interrupt by user" )
    LOGGER . info ( "assembly interrupted by user." )
except IPyradWarningExit as inst :
    LOGGER . error ( "IPyradWarningExit: %s" , inst )
    print ( "\n  Encountered an error (see details in ./ipyrad_log.txt)" + "\n  Error summary is below -------------------------------" + "\n{}" . format ( inst ) )
except Exception as inst :
    LOGGER . error ( inst )
    print ( "\n  Encountered an unexpected error (see ./ipyrad_log.txt)" + "\n  Error message is below -------------------------------" + "\n{}" . format ( inst ) )
# # close client when done or interrupted
finally :
    try : # # save the Assembly
        self . save ( )
        # # can ' t close client if it was never open
        if ipyclient : # # send SIGINT ( 2 ) to all engines
            try :
                ipyclient . abort ( )
                time . sleep ( 1 )
                for engine_id , pid in self . _ipcluster [ "pids" ] . items ( ) :
                    if ipyclient . queue_status ( ) [ engine_id ] [ "tasks" ] :
                        os . kill ( pid , 2 )
                        LOGGER . info ( 'interrupted engine {} w/ SIGINT to {}' . format ( engine_id , pid ) )
                time . sleep ( 1 )
            except ipp . NoEnginesRegistered :
                pass
            # # if CLI , stop jobs and shutdown . Don ' t use _ cli here
            # # because you can have a CLI object but use the - - ipcluster
            # # flag , in which case we don ' t want to kill ipcluster .
            if 'ipyrad-cli' in self . _ipcluster [ "cluster_id" ] :
                LOGGER . info ( "  shutting down engines" )
                ipyclient . shutdown ( hub = True , block = False )
                ipyclient . close ( )
                LOGGER . info ( "  finished shutdown" )
            else :
                if not ipyclient . outstanding :
                    ipyclient . purge_everything ( )
                else : # # nanny : kill everything , something bad happened
                    ipyclient . shutdown ( hub = True , block = False )
                    ipyclient . close ( )
                    print ( "\nwarning: ipcluster shutdown and must be restarted" )
    # # if exception is close and save , print and ignore
    except Exception as inst2 :
        print ( "warning: error during shutdown:\n{}" . format ( inst2 ) )
        LOGGER . error ( "shutdown warning: %s" , inst2 )


signature:
def _to_fulldict ( self ) :
    """Write to dict including data frames . All sample dicts
 are combined in save ( ) to dump JSON output"""

implementation:
returndict = OrderedDict ( [ ( "name" , self . name ) , ( "barcode" , self . barcode ) , ( "files" , self . files ) , ( "stats_dfs" , { "s1" : self . stats_dfs . s1 . to_dict ( ) , "s2" : self . stats_dfs . s2 . to_dict ( ) , "s3" : self . stats_dfs . s3 . to_dict ( ) , "s4" : self . stats_dfs . s4 . to_dict ( ) , "s5" : self . stats_dfs . s5 . to_dict ( ) , } ) , ( "stats" , self . stats . to_dict ( ) ) , ( "depths" , self . depths ) ] )
return returndict


signature:
def combinefiles ( filepath ) :
    """Joins first and second read file names"""

implementation:
# # unpack seq files in filepath
fastqs = glob . glob ( filepath )
firsts = [ i for i in fastqs if "_R1_" in i ]
# # check names
if not firsts :
    raise IPyradWarningExit ( "First read files names must contain '_R1_'." )
# # get paired reads
seconds = [ ff . replace ( "_R1_" , "_R2_" ) for ff in firsts ]
return zip ( firsts , seconds )


signature:
def findbcode ( cutters , longbar , read1 ) :
    """find barcode sequence in the beginning of read"""

implementation:
# # default barcode string
for cutter in cutters [ 0 ] : # # If the cutter is unambiguous there will only be one .
    if not cutter :
        continue
    search = read1 [ 1 ] [ : int ( longbar [ 0 ] + len ( cutter ) + 1 ) ]
    barcode = search . rsplit ( cutter , 1 )
    if len ( barcode ) > 1 :
        return barcode [ 0 ]
# # No cutter found
return barcode [ 0 ]


signature:
def find3radbcode ( cutters , longbar , read1 ) :
    """find barcode sequence in the beginning of read"""

implementation:
# # default barcode string
for ambigcuts in cutters :
    for cutter in ambigcuts : # # If the cutter is unambiguous there will only be one .
        if not cutter :
            continue
        search = read1 [ 1 ] [ : int ( longbar [ 0 ] + len ( cutter ) + 1 ) ]
        splitsearch = search . rsplit ( cutter , 1 )
        if len ( splitsearch ) > 1 :
            return splitsearch [ 0 ]
# # No cutter found
return splitsearch [ 0 ]


signature:
def make_stats ( data , perfile , fsamplehits , fbarhits , fmisses , fdbars ) :
    """Write stats and stores to Assembly object ."""

implementation:
# # out file
outhandle = os . path . join ( data . dirs . fastqs , 's1_demultiplex_stats.txt' )
outfile = open ( outhandle , 'w' )
# # write the header for file stats - - - - -
outfile . write ( '{:<35}  {:>13}{:>13}{:>13}\n' . format ( "raw_file" , "total_reads" , "cut_found" , "bar_matched" ) )
# # write the file stats
r1names = sorted ( perfile )
for fname in r1names :
    dat = perfile [ fname ]
    # dat = [ perfile [ fname ] [ i ] for i in [ " ftotal " , " fcutfound " , " fmatched " ] ]
    outfile . write ( '{:<35}  {:>13}{:>13}{:>13}\n' . format ( fname , dat [ 0 ] , dat [ 1 ] , dat [ 2 ] ) )
    # # repeat for pairfile
    if 'pair' in data . paramsdict [ "datatype" ] :
        fname = fname . replace ( "_R1_" , "_R2_" )
        outfile . write ( '{:<35}  {:>13}{:>13}{:>13}\n' . format ( fname , dat [ 0 ] , dat [ 1 ] , dat [ 2 ] ) )
# # spacer , how many records for each sample - - - - -
outfile . write ( '\n{:<35}  {:>13}\n' . format ( "sample_name" , "total_reads" ) )
# # names alphabetical . Write to file . Will save again below to Samples .
snames = set ( )
for sname in data . barcodes :
    if "-technical-replicate-" in sname :
        sname = sname . rsplit ( "-technical-replicate" , 1 ) [ 0 ]
    snames . add ( sname )
for sname in sorted ( list ( snames ) ) :
    outfile . write ( "{:<35}  {:>13}\n" . format ( sname , fsamplehits [ sname ] ) )
# # spacer , which barcodes were found - - - - -
outfile . write ( '\n{:<35}  {:>13} {:>13} {:>13}\n' . format ( "sample_name" , "true_bar" , "obs_bar" , "N_records" ) )
# # write sample results
for sname in sorted ( data . barcodes ) :
    if "-technical-replicate-" in sname :
        fname = sname . rsplit ( "-technical-replicate" , 1 ) [ 0 ]
    else :
        fname = sname
    # # write perfect hit
    hit = data . barcodes [ sname ]
    offhitstring = ""
    # # write off - n hits
    # # sort list of off - n hits
    if fname in fdbars :
        offkeys = list ( fdbars . get ( fname ) )
        for offhit in offkeys [ : : - 1 ] : # # exclude perfect hit
            if offhit not in data . barcodes . values ( ) :
                offhitstring += '{:<35}  {:>13} {:>13} {:>13}\n' . format ( sname , hit , offhit , fbarhits [ offhit ] / 2 )
                # sumoffhits + = fbarhits [ offhit ]
        # # write string to file
        outfile . write ( '{:<35}  {:>13} {:>13} {:>13}\n' . # format ( sname , hit , hit , fsamplehits [ fname ] - sumoffhits ) )
        format ( sname , hit , hit , fbarhits [ hit ] / 2 ) )
        outfile . write ( offhitstring )
# # write misses
misskeys = list ( fmisses . keys ( ) )
misskeys . sort ( key = fmisses . get )
for key in misskeys [ : : - 1 ] :
    outfile . write ( '{:<35}  {:>13}{:>13}{:>13}\n' . format ( "no_match" , "_" , key , fmisses [ key ] ) )
outfile . close ( )
# # Link Sample with this data file to the Assembly object
for sname in snames : # # make the sample
    sample = Sample ( )
    sample . name = sname
    # # allow multiple barcodes if its a replicate .
    barcodes = [ ]
    for n in xrange ( 500 ) :
        fname = sname + "-technical-replicate-{}" . format ( n )
        fbar = data . barcodes . get ( fname )
        if fbar :
            barcodes . append ( fbar )
    if barcodes :
        sample . barcode = barcodes
    else :
        sample . barcode = data . barcodes [ sname ]
    # # file names
    if 'pair' in data . paramsdict [ "datatype" ] :
        sample . files . fastqs = [ ( os . path . join ( data . dirs . fastqs , sname + "_R1_.fastq.gz" ) , os . path . join ( data . dirs . fastqs , sname + "_R2_.fastq.gz" ) ) ]
    else :
        sample . files . fastqs = [ ( os . path . join ( data . dirs . fastqs , sname + "_R1_.fastq.gz" ) , "" ) ]
    # # fill in the summary stats
    sample . stats [ "reads_raw" ] = int ( fsamplehits [ sname ] )
    # # fill in the full df stats value
    sample . stats_dfs . s1 [ "reads_raw" ] = int ( fsamplehits [ sname ] )
    # # Only link Sample if it has data
    if sample . stats [ "reads_raw" ] :
        sample . stats . state = 1
        data . samples [ sample . name ] = sample
    else :
        print ( "Excluded sample: no data found for" , sname )
# # initiate s1 key for data object
data . stats_dfs . s1 = data . _build_stat ( "s1" )
data . stats_files . s1 = outhandle


signature:
def barmatch2 ( data , tups , cutters , longbar , matchdict , fnum ) :
    """cleaner barmatch func . . ."""

implementation:
# # how many reads to store before writing to disk
waitchunk = int ( 1e6 )
# # pid name for this engine
epid = os . getpid ( )
# # counters for total reads , those with cutsite , and those that matched
filestat = np . zeros ( 3 , dtype = np . int )
# # store observed sample matches
samplehits = { }
# # dictionaries to store first and second reads until writing to file
dsort1 = { }
dsort2 = { }
# # dictionary for all bars matched in sample
dbars = { }
# # fill for sample names
for sname in data . barcodes :
    if "-technical-replicate-" in sname :
        sname = sname . rsplit ( "-technical-replicate" , 1 ) [ 0 ]
    samplehits [ sname ] = 0
    dsort1 [ sname ] = [ ]
    dsort2 [ sname ] = [ ]
    dbars [ sname ] = set ( )
# # store observed bars
barhits = { }
for barc in matchdict :
    barhits [ barc ] = 0
# # store others
misses = { }
misses [ '_' ] = 0
# # build func for finding barcode
getbarcode = get_barcode_func ( data , longbar )
# # get quart iterator of reads
if tups [ 0 ] . endswith ( ".gz" ) :
    ofunc = gzip . open
else :
    ofunc = open
# # create iterators
ofile1 = ofunc ( tups [ 0 ] , 'r' )
fr1 = iter ( ofile1 )
quart1 = itertools . izip ( fr1 , fr1 , fr1 , fr1 )
if tups [ 1 ] :
    ofile2 = ofunc ( tups [ 1 ] , 'r' )
    fr2 = iter ( ofile2 )
    quart2 = itertools . izip ( fr2 , fr2 , fr2 , fr2 )
    quarts = itertools . izip ( quart1 , quart2 )
else :
    quarts = itertools . izip ( quart1 , iter ( int , 1 ) )
# # go until end of the file
while 1 :
    try :
        read1 , read2 = quarts . next ( )
        read1 = list ( read1 )
        filestat [ 0 ] += 1
    except StopIteration :
        break
    barcode = ""
    # # Get barcode _ R2 and check for matching sample name
    if '3rad' in data . paramsdict [ "datatype" ] : # # Here we ' re just reusing the findbcode function
    # # for R2 , and reconfiguring the longbar tuple to have the
    # # maxlen for the R2 barcode
    # # Parse barcode . Use the parsing function selected above .
        barcode1 = find3radbcode ( cutters = cutters , longbar = longbar , read1 = read1 )
        barcode2 = find3radbcode ( cutters = cutters , longbar = ( longbar [ 2 ] , longbar [ 1 ] ) , read1 = read2 )
        barcode = barcode1 + "+" + barcode2
    else : # # Parse barcode . Uses the parsing function selected above .
        barcode = getbarcode ( cutters , read1 , longbar )
    # # find if it matches
    sname_match = matchdict . get ( barcode )
    if sname_match : # sample _ index [ filestat [ 0 ] - 1 ] = snames . index ( sname _ match ) + 1
    # # record who matched
        dbars [ sname_match ] . add ( barcode )
        filestat [ 1 ] += 1
        filestat [ 2 ] += 1
        samplehits [ sname_match ] += 1
        barhits [ barcode ] += 1
        if barcode in barhits :
            barhits [ barcode ] += 1
        else :
            barhits [ barcode ] = 1
        # # trim off barcode
        lenbar = len ( barcode )
        if '3rad' in data . paramsdict [ "datatype" ] : # # Iff 3rad trim the len of the first barcode
            lenbar = len ( barcode1 )
        if data . paramsdict [ "datatype" ] == '2brad' :
            overlen = len ( cutters [ 0 ] [ 0 ] ) + lenbar + 1
            read1 [ 1 ] = read1 [ 1 ] [ : - overlen ] + "\n"
            read1 [ 3 ] = read1 [ 3 ] [ : - overlen ] + "\n"
        else :
            read1 [ 1 ] = read1 [ 1 ] [ lenbar : ]
            read1 [ 3 ] = read1 [ 3 ] [ lenbar : ]
        # # Trim barcode off R2 and append . Only 3rad datatype
        # # pays the cpu cost of splitting R2
        if '3rad' in data . paramsdict [ "datatype" ] :
            read2 = list ( read2 )
            read2 [ 1 ] = read2 [ 1 ] [ len ( barcode2 ) : ]
            read2 [ 3 ] = read2 [ 3 ] [ len ( barcode2 ) : ]
        # # append to dsort
        dsort1 [ sname_match ] . append ( "" . join ( read1 ) )
        if 'pair' in data . paramsdict [ "datatype" ] :
            dsort2 [ sname_match ] . append ( "" . join ( read2 ) )
    else :
        misses [ "_" ] += 1
        if barcode :
            filestat [ 1 ] += 1
    # # how can we make it so all of the engines aren ' t trying to write to
    # # ~ 100-200 files all at the same time ? This is the I / O limit we hit . .
    # # write out at 100K to keep memory low . It is fine on HPC which can
    # # write parallel , but regular systems might crash
    if not filestat [ 0 ] % waitchunk : # # write the remaining reads to file "
        writetofile ( data , dsort1 , 1 , epid )
        if 'pair' in data . paramsdict [ "datatype" ] :
            writetofile ( data , dsort2 , 2 , epid )
        # # clear out dsorts
        for sample in data . barcodes :
            if "-technical-replicate-" in sname :
                sname = sname . rsplit ( "-technical-replicate" , 1 ) [ 0 ]
            dsort1 [ sname ] = [ ]
            dsort2 [ sname ] = [ ]
        # # reset longlist
        # longlist = np . zeros ( waitchunk , dtype = np . uint32)
# # close open files
ofile1 . close ( )
if tups [ 1 ] :
    ofile2 . close ( )
# # write the remaining reads to file
writetofile ( data , dsort1 , 1 , epid )
if 'pair' in data . paramsdict [ "datatype" ] :
    writetofile ( data , dsort2 , 2 , epid )
# # return stats in saved pickle b / c return _ queue is too small
# # and the size of the match dictionary can become quite large
samplestats = [ samplehits , barhits , misses , dbars ]
outname = os . path . join ( data . dirs . fastqs , "tmp_{}_{}.p" . format ( epid , fnum ) )
with open ( outname , 'w' ) as wout :
    pickle . dump ( [ filestat , samplestats ] , wout )
return outname


signature:
def get_quart_iter ( tups ) :
    """returns an iterator to grab four lines at a time"""

implementation:
if tups [ 0 ] . endswith ( ".gz" ) :
    ofunc = gzip . open
else :
    ofunc = open
# # create iterators
ofile1 = ofunc ( tups [ 0 ] , 'r' )
fr1 = iter ( ofile1 )
quart1 = itertools . izip ( fr1 , fr1 , fr1 , fr1 )
if tups [ 1 ] :
    ofile2 = ofunc ( tups [ 1 ] , 'r' )
    fr2 = iter ( ofile2 )
    quart2 = itertools . izip ( fr2 , fr2 , fr2 , fr2 )
    quarts = itertools . izip ( quart1 , quart2 )
else :
    ofile2 = 0
    quarts = itertools . izip ( quart1 , iter ( int , 1 ) )
# # make a generator
def feedme ( quarts ) :
    for quart in quarts :
        yield quart
genquarts = feedme ( quarts )
# # return generator and handles
return genquarts , ofile1 , ofile2


signature:
def writetofastq ( data , dsort , read ) :
    """Writes sorted data ' dsort dict ' to a tmp files"""

implementation:
if read == 1 :
    rrr = "R1"
else :
    rrr = "R2"
for sname in dsort : # # skip writing if empty . Write to tmpname
    handle = os . path . join ( data . dirs . fastqs , "{}_{}_.fastq" . format ( sname , rrr ) )
    with open ( handle , 'a' ) as out :
        out . write ( "" . join ( dsort [ sname ] ) )