%prog mask fastafile

    Mask the contaminants. By default, this will compare against UniVec_Core and
    Ecoli.fasta. Merge the contaminant results, and use `maskFastaFromBed`. Can
    perform FASTA tidy if requested.
17902
def mask(args):
    
    p = OptionParser(mask.__doc__)
    p.add_option("--db",
                 help="Contaminant db other than Ecoli K12 [default: %default]")
    opts, args = p.parse_args(args)

    if len(args) != 1:
        sys.exit(not p.print_help())

    fastafile, = args
    assert op.exists(fastafile)

    outfastafile = fastafile.rsplit(".", 1)[0] + ".masked.fasta"
    vecbedfile = blast([fastafile])
    ecoliurl = \
    "ftp://ftp.ncbi.nih.gov/genomes/Bacteria/Escherichia_coli_K_12_substr__DH10B_uid58979/NC_010473.fna"
    ecolifile = opts.db or download(ecoliurl, filename="Ecoli.fasta")
    assert op.exists(ecolifile)
    ecolibedfile = blast([fastafile, "--db={0}".format(ecolifile)])

    cmd = "cat {0} {1}".format(vecbedfile, ecolibedfile)
    cmd += " | mergeBed -nms -d 100 -i stdin"
    cmd += " | maskFastaFromBed -fi {0} -bed stdin -fo {1}".\
            format(fastafile, outfastafile)
    sh(cmd)

    return tidy([outfastafile])
173778
def execute(self, eopatch):
        
        feature_type, feature_name, new_feature_name = next(self.feature(eopatch))
        mask_feature_type, mask_feature_name = next(self.mask_feature(eopatch))

        data = np.copy(eopatch[feature_type][feature_name])
        mask = eopatch[mask_feature_type][mask_feature_name]

        if not isinstance(self.mask_values, list):
            raise ValueError('Incorrect format or values of argument `mask_values`')

        for value in self.mask_values:
            data[mask.squeeze() == value] = self.no_data_value

        eopatch.add_feature(feature_type, new_feature_name, data)

        return eopatch
642252
def fasta2bed(fastafile):
    
    dustfasta = fastafile.rsplit(".", 1)[0] + ".dust.fasta"
    for name, seq in parse_fasta(dustfasta):
        for islower, ss in groupby(enumerate(seq), key=lambda x: x[-1].islower()):
            if not islower:
                continue
            ss = list(ss)
            ms, mn = min(ss)
            xs, xn = max(ss)
            print("\t".join(str(x) for x in (name, ms, xs)))
442345
def patch(args):
    
    from jcvi.formats.base import write_file
    from jcvi.formats.fasta import format

    p = OptionParser(patch.__doc__)
    p.add_option("--cleanfasta", default=False, action="store_true",
                 help="Clean FASTA to remove description [default: %default]")
    p.add_option("--highqual", default=False, action="store_true",
                 help="Reads are of high quality [default: %default]")
    p.set_home("pbjelly")
    p.set_cpus()
    opts, args = p.parse_args(args)

    if len(args) != 2:
        sys.exit(not p.print_help())

    ref, reads = args
    cpus = opts.cpus
    cmd = op.join(opts.pbjelly_home, "setup.sh")
    setup = "source {0}".format(cmd)
    if not which("fakeQuals.py"):
        sh(setup)

    pf = ref.rsplit(".", 1)[0]
    pr, px = reads.rsplit(".", 1)
    # Remove description line
    if opts.cleanfasta:
        oref = pf + ".f.fasta"
        oreads = pr + ".f.fasta"
        format([ref, oref])
        format([reads, oreads])
        ref, reads = oref, oreads

    # Check if the FASTA has qual
    ref, refq = fake_quals(ref)
    convert_reads = not px in ("fq", "fastq", "txt")
    if convert_reads:
        reads, readsq = fake_quals(reads)
        readsfiles = " ".join((reads, readsq))
    else:
        readsfiles = reads

    # Make directory structure
    dref, dreads = "data/reference", "data/reads"
    cwd = os.getcwd()
    reference = op.join(cwd, "{0}/{1}".format(dref, ref))
    reads = op.join(cwd, "{0}/{1}".format(dreads, reads))
    if not op.exists(reference):
        sh("mkdir -p {0}".format(dref))
        sh("cp {0} {1}/".format(" ".join((ref, refq)), dref))
    if not op.exists(reads):
        sh("mkdir -p {0}".format(dreads))
        sh("cp {0} {1}/".format(readsfiles, dreads))

    outputDir = cwd
    p = Protocol(outputDir, reference, reads, highqual=opts.highqual)
    p.write_xml()

    # Build the pipeline
    runsh = [setup]
    for action in "setup|mapping|support|extraction".split("|"):
        runsh.append("Jelly.py {0} Protocol.xml".format(action))

    runsh.append('Jelly.py assembly Protocol.xml -x "--nproc={0}"'.format(cpus))
    runsh.append("Jelly.py output Protocol.xml")

    runfile = "run.sh"
    contents = "\n".join(runsh)
    write_file(runfile, contents)
843161
def replaceMaskedValue(data, mask, replacementValue, copyOnReplace=True):
    
    if mask is False:
        result = data
    elif mask is True:
        result = np.copy(data) if copyOnReplace else data
        result[:] = replacementValue
    else:
        #logger.debug("############ count_nonzero: {}".format(np.count_nonzero(mask)))
        if copyOnReplace and np.any(mask):
            #logger.debug("Making copy")
            result = np.copy(data)
        else:
            result = data

        result[mask] = replacementValue

    return result
62368
def fastq(args):
    
    from jcvi.formats.fastq import FastqLite

    p = OptionParser(fastq.__doc__)
    p.add_option("--qv", type="int",
                 help="Use generic qv value [dafault: %default]")

    opts, args = p.parse_args(args)

    if len(args) != 1:
        sys.exit(not p.print_help())

    fastafile, = args
    fastqfile = fastafile.rsplit(".", 1)[0] + ".fastq"
    fastqhandle = open(fastqfile, "w")
    num_records = 0

    if opts.qv is not None:
        qv = chr(ord('!') + opts.qv)
        logging.debug("QV char '{0}' ({1})".format(qv, opts.qv))
    else:
        qv = None

    if qv:
        f = Fasta(fastafile, lazy=True)
        for name, rec in f.iteritems_ordered():
            r = FastqLite("@" + name, str(rec.seq).upper(), qv * len(rec.seq))
            print(r, file=fastqhandle)
            num_records += 1

    else:
        qualfile = get_qual(fastafile)
        for rec in iter_fasta_qual(fastafile, qualfile):
            SeqIO.write([rec], fastqhandle, "fastq")
            num_records += 1

    fastqhandle.close()
    logging.debug("A total of %d records written to `%s`" % \
            (num_records, fastqfile))
1293481
def do_rmfilter(self, arg):
        
        if arg in self.curargs["tfilter"]:
            if arg == "*":
                msg.warn("The default filter cannot be removed.")
            else:
                self.curargs["tfilter"].remove(arg)
                self.do_filter("list")
1046951
def _match_mask(mask, ctype):
    

    # Handle the simple cases first
    if '*' not in mask:
        return ctype == mask
    elif mask == '*/*':
        return True
    elif not mask.endswith('/*'):
        return False

    mask_major = mask[:-2]
    ctype_major = ctype.split('/', 1)[0]
    return ctype_major == mask_major
666804
def fastq(args):
    
    p = OptionParser(fastq.__doc__)
    opts, args = p.parse_args(args)

    if len(args) != 2:
        sys.exit(not p.print_help())

    bamfile, pf = args
    singletons = pf + ".se.fastq"
    a = pf + ".read1.fastq"
    b = pf + ".read2.fastq"

    cmd  = "samtools collate -uOn 128 {} tmp-prefix".format(bamfile)
    cmd += " | samtools fastq -s {} -1 {} -2 {} -"\
                .format(singletons, a, b)
    sh(cmd)

    if os.stat(singletons).st_size == 0:  # singleton file is empty
        os.remove(singletons)
    return a, b
1235091
def remove_option(self, mask):
        
        if not isinstance(mask, int):
            raise TypeError("mask must be an int")
        self.__check_okay_to_chain()

        if mask & _QUERY_OPTIONS["exhaust"]:
            self.__exhaust = False

        self.__query_flags &= ~mask
        return self
1024609
def mask(args):
    
    p = OptionParser(mask.__doc__)
    p.add_option("--hard", dest="hard", default=False, action="store_true",
            help="Hard mask the low-complexity bases [default: %default]")
    opts, args = p.parse_args(args)

    if len(args) != 1:
        sys.exit(not p.print_help())

    genomefile, = args

    # entire pipeline
    countsfile = genomefile + ".counts"
    wm_mk_counts(infile=genomefile, outfile=countsfile)

    maskedfastafile = "%s.masked%s" % op.splitext(genomefile)
    wm_mk_masks(infile=countsfile, outfile=maskedfastafile, genomefile=genomefile)

    if opts.hard:
        hardmask(maskedfastafile)
888641
def pasteprepare(args):
    
    p = OptionParser(pasteprepare.__doc__)
    p.add_option("--flank", default=5000, type="int",
                 help="Get the seq of size on two ends [default: %default]")
    opts, args = p.parse_args(args)

    if len(args) != 1:
        sys.exit(not p.print_help())

    goodfasta, = args
    flank = opts.flank
    pf = goodfasta.rsplit(".", 1)[0]
    extbed = pf + ".ext.bed"

    sizes = Sizes(goodfasta)
    fw = open(extbed, "w")
    for bac, size in sizes.iter_sizes():
        print("\t".join(str(x) for x in \
                               (bac, 0, min(flank, size), bac + "L")), file=fw)
        print("\t".join(str(x) for x in \
                               (bac, max(size - flank, 0), size, bac + "R")), file=fw)
    fw.close()

    fastaFromBed(extbed, goodfasta, name=True)
1114815
def strip_inserts(fasta):
    
    for seq in parse_fasta(fasta):
        seq[1] = ''.join([b for b in seq[1] if b == '-' or b.isupper()])
        yield seq
547261
static int mask(Class<? extends ChannelHandler> clazz) {
        // Try to obtain the mask from the cache first. If this fails calculate it and put it in the cache for fast
        // lookup in the future.
        Map<Class<? extends ChannelHandler>, Integer> cache = MASKS.get();
        Integer mask = cache.get(clazz);
        if (mask == null) {
            mask = mask0(clazz);
            cache.put(clazz, mask);
        }
        return mask;
    }
193900
def _try_mask_first_value(value, row, all_close):
    '''
    
    '''
    # Compare value to row
    for i, value2 in enumerate(row):
        if _value_equals(value, value2, all_close):
            row[i] = ma.masked
            return True
    return False
261136
def get_submask(self, *args, **kwargs) -> masktools.CustomMask:
        
        if args or kwargs:
            masks = self.availablemasks
            mask = masktools.CustomMask(numpy.full(self.shape, False))
            for arg in args:
                mask = mask + self._prepare_mask(arg, masks)
            for key, value in kwargs.items():
                mask = mask + self._prepare_mask(key, masks, **value)
            if mask not in self.mask:
                raise ValueError(
                    f'Based on the arguments `{args}` and `{kwargs}` '
                    f'the mask `{repr(mask)}` has been determined, '
                    f'which is not a submask of `{repr(self.mask)}`.')
        else:
            mask = self.mask
        return mask
886330
def fasta(args):
    
    from jcvi.formats.sizes import Sizes

    p = OptionParser(fasta.__doc__)
    p.add_option("--extend", default=1000, type="int",
                 help="Extend seq flanking the gaps [default: %default]")
    opts, args = p.parse_args(args)

    if len(args) != 2:
        sys.exit(not p.print_help())

    mapout, sfasta = args
    Flank = opts.extend
    pf = mapout.split(".")[0]
    mapbed = pf + ".bed"
    bm = BinMap(mapout)
    bm.print_to_bed(mapbed)

    bed = Bed(mapbed, sorted=False)
    markersbed = pf + ".markers.bed"
    fw = open(markersbed, "w")
    sizes = Sizes(sfasta).mapping
    for b in bed:
        accn = b.accn
        scf, pos = accn.split(".")
        pos = int(pos)
        start = max(0, pos - Flank)
        end = min(pos + Flank, sizes[scf])
        print("\t".join(str(x) for x in \
                    (scf, start, end, accn)), file=fw)

    fw.close()

    fastaFromBed(markersbed, sfasta, name=True)