Parse the Bismark M-Bias plot data
331022
def parse_bismark_mbias(self, f):
        
        s = f['s_name']
        self.bismark_mbias_data['meth']['CpG_R1'][s] = {}
        self.bismark_mbias_data['meth']['CHG_R1'][s] = {}
        self.bismark_mbias_data['meth']['CHH_R1'][s] = {}
        self.bismark_mbias_data['cov']['CpG_R1'][s] = {}
        self.bismark_mbias_data['cov']['CHG_R1'][s] = {}
        self.bismark_mbias_data['cov']['CHH_R1'][s] = {}
        self.bismark_mbias_data['meth']['CpG_R2'][s] = {}
        self.bismark_mbias_data['meth']['CHG_R2'][s] = {}
        self.bismark_mbias_data['meth']['CHH_R2'][s] = {}
        self.bismark_mbias_data['cov']['CpG_R2'][s] = {}
        self.bismark_mbias_data['cov']['CHG_R2'][s] = {}
        self.bismark_mbias_data['cov']['CHH_R2'][s] = {}
        key = None
        for l in f['f']:
            if 'context' in l:
                if 'CpG' in l:
                    key = 'CpG'
                elif 'CHG' in l:
                    key = 'CHG'
                elif 'CHH' in l:
                    key = 'CHH'
                if '(R1)' in l:
                    key += '_R1'
                elif '(R2)' in l:
                    key += '_R2'
                else:
                    key += '_R1'
            if key is not None:
                sections = l.split()
                try:
                    pos = int(sections[0])
                    self.bismark_mbias_data['meth'][key][s][pos] = float(sections[3])
                    self.bismark_mbias_data['cov'][key][s][pos] = int(sections[4])
                except (IndexError, ValueError):
                    continue

        # Remove empty dicts (eg. R2 for SE data)
        for t in self.bismark_mbias_data:
            for k in self.bismark_mbias_data[t]:
                self.bismark_mbias_data[t][k] = {
                    s_name: self.bismark_mbias_data[t][k][s_name]
                    for s_name in self.bismark_mbias_data[t][k]
                    if len(self.bismark_mbias_data[t][k][s_name]) > 0
                }
949979
def parse_samtools_idxstats(self):
        

        self.samtools_idxstats = dict()
        for f in self.find_log_files('samtools/idxstats'):
            parsed_data = parse_single_report(f['f'])
            if len(parsed_data) > 0:
                if f['s_name'] in self.samtools_idxstats:
                    log.debug("Duplicate sample name found! Overwriting: {}".format(f['s_name']))
                self.add_data_source(f, section='idxstats')
                self.samtools_idxstats[f['s_name']] = parsed_data

        # Filter to strip out ignored sample names
        self.samtools_idxstats = self.ignore_samples(self.samtools_idxstats)

        if len(self.samtools_idxstats) > 0:

            # Write parsed report data to a file (restructure first)
            self.write_data_file(self.samtools_idxstats, 'multiqc_samtools_idxstats')

            # Prep the data for the plots
            keys = list()
            pdata = dict()
            pdata_norm = dict()
            xy_counts = dict()
            # Count the total mapped reads for every chromosome
            chrs_mapped = defaultdict(lambda:0)
            sample_mapped = defaultdict(lambda:0)
            total_mapped = 0
            # Cutoff, can be customised in config
            cutoff = float(getattr(config, 'samtools_idxstats_fraction_cutoff', 0.001))
            if cutoff != 0.001:
                log.info('Setting idxstats cutoff to: {}%'.format(cutoff*100.0))
            for s_name in self.samtools_idxstats:
                for chrom in self.samtools_idxstats[s_name]:
                    chrs_mapped[chrom] += self.samtools_idxstats[s_name][chrom]
                    sample_mapped[s_name] += self.samtools_idxstats[s_name][chrom]
                    total_mapped += self.samtools_idxstats[s_name][chrom]
            req_reads = float(total_mapped)*cutoff
            chr_always = getattr(config, 'samtools_idxstats_always', [])
            if len(chr_always) > 0:
                log.info('Trying to include these chromosomes in idxstats: {}'.format(', '.join(chr_always)))
            chr_ignore = getattr(config, 'samtools_idxstats_ignore', [])
            if len(chr_ignore) > 0:
                log.info('Excluding these chromosomes from idxstats: {}'.format(', '.join(chr_ignore)))
            xchr = getattr(config, 'samtools_idxstats_xchr', False)
            if xchr:
                log.info('Using "{}" as X chromosome name'.format(xchr))
            ychr = getattr(config, 'samtools_idxstats_ychr', False)
            if ychr:
                log.info('Using "{}" as Y chromosome name'.format(ychr))
            # Go through again and collect all of the keys that have enough counts
            # Also get the X/Y counts if we find them
            for s_name in self.samtools_idxstats:
                x_count = False
                y_count = False
                for chrom in self.samtools_idxstats[s_name]:
                    if float(chrs_mapped[chrom]) > req_reads or chrom in chr_always:
                        if chrom not in chr_ignore and chrom not in keys:
                            keys.append(chrom)
                    # Collect X and Y counts if we have them
                    mapped = self.samtools_idxstats[s_name][chrom]
                    if xchr is not False :
                        if str(xchr) == str(chrom):
                            x_count = mapped
                    else:
                        if chrom.lower() == 'x' or chrom.lower() == 'chrx':
                            x_count = mapped
                    if ychr is not False :
                        if str(ychr) == str(chrom):
                            y_count = mapped
                    else:
                        if chrom.lower() == 'y' or chrom.lower() == 'chry':
                            y_count = mapped
                # Only save these counts if we have both x and y
                if x_count and y_count:
                    xy_counts[s_name] = {'x': x_count, 'y': y_count }
            # Ok, one last time. We have the chromosomes that we want to plot,
            # now collect the counts
            for s_name in self.samtools_idxstats:
                pdata[s_name] = OrderedDict()
                pdata_norm[s_name] = OrderedDict()
                for k in keys:
                    try:
                        pdata[s_name][k] = self.samtools_idxstats[s_name][k]
                        pdata_norm[s_name][k] = float(self.samtools_idxstats[s_name][k]) / sample_mapped[s_name]
                    except (KeyError, ZeroDivisionError):
                        pdata[s_name][k] = 0
                        pdata_norm[s_name][k] = 0

            # X/Y ratio plot
            if len(xy_counts) > 0:
                xy_keys = OrderedDict()
                xy_keys['x'] = { 'name': xchr if xchr else 'Chromosome X' }
                xy_keys['y'] = { 'name': ychr if ychr else 'Chromosome Y' }
                pconfig = {
                    'id': 'samtools-idxstats-xy-plot',
                    'title': 'Samtools idxstats: chrXY mapped reads',
                    'ylab': 'Percent of X+Y Reads',
                    'cpswitch_counts_label': 'Number of Reads',
                    'cpswitch_percent_label': 'Percent of X+Y Reads',
                    'cpswitch_c_active': False
                }
                self.add_section (
                    name = 'XY counts',
                    anchor = 'samtools-idxstats-xy-counts',
                    plot = bargraph.plot(xy_counts, xy_keys, pconfig)
                )


            # Mapped reads per chr line plot
            pconfig = {
                'id': 'samtools-idxstats-mapped-reads-plot',
                'title': 'Samtools idxstats: Mapped reads per contig',
                'ylab': '# mapped reads',
                'xlab': 'Chromosome Name',
                'categories': True,
                'tt_label': '<strong>{point.category}:</strong> {point.y:.2f}',
                'data_labels': [
                    {'name': 'Normalised Counts', 'ylab': 'Fraction of total count'},
                    {'name': 'Counts', 'ylab': '# mapped reads'}
                ]
            }
            self.add_section (
                name = 'Mapped reads per contig',
                anchor = 'samtools-idxstats',
                description = 'The <code>samtools idxstats</code> tool counts the number of mapped reads per chromosome / contig. ' +
                    'Chromosomes with &lt; {}% of the total aligned reads are omitted from this plot.'.format(cutoff*100),
                plot = linegraph.plot([pdata_norm, pdata], pconfig)
            )

        # Return the number of logs that were found
        return len(self.samtools_idxstats)
1119142
def parse_plotCorrelation(self):
        """"""
        self.deeptools_plotCorrelationData = dict()
        for f in self.find_log_files('deeptools/plotCorrelationData', filehandles=False):
            parsed_data, samples = self.parsePlotCorrelationData(f)
            for k, v in parsed_data.items():
                if k in self.deeptools_plotCorrelationData:
                    log.warning("Replacing duplicate sample {}.".format(k))
                self.deeptools_plotCorrelationData[k] = v
            if len(parsed_data) > 0:
                self.add_data_source(f, section='plotCorrelation')

        if len(self.deeptools_plotCorrelationData) > 0:
            config = {
                'id': 'deeptools_correlation_plot',
                'title': 'deeptools: Correlation Plot',
            }
            data = []
            for s_name in samples:
                try:
                    data.append(self.deeptools_plotCorrelationData[s_name])
                except KeyError:
                    pass
            if len(data) == 0:
                log.debug('No valid data for correlation plot')
                return None

            self.add_section(
                name="Correlation heatmap",
                anchor="deeptools_correlation",
                description="Pairwise correlations of samples based on distribution of sequence reads",
                plot=heatmap.plot(data, samples, samples, config)
            )

        return len(self.deeptools_plotCorrelationData)
894891
def _plot_x(self, iabscissa=1, x_opt=None, remark=None,
                annotations=None):
        
        if not hasattr(self, 'x'):
            _print_warning('no x-attributed found, use methods ' +
                           'plot_xrecent or plot_mean', 'plot_x',
                           'CMADataLogger')
            return
        from matplotlib.pyplot import plot, semilogy, hold, text, grid, axis, title
        dat = self  # for convenience and historical reasons
        # modify fake last entry in x for line extension-annotation
        if dat.x.shape[1] < 100:
            minxend = int(1.06 * dat.x[-2, iabscissa])
            # write y-values for individual annotation into dat.x
            dat.x[-1, iabscissa] = minxend  # TODO: should be ax[1]
            if x_opt is None:
                idx = np.argsort(dat.x[-2, 5:])
                idx2 = np.argsort(idx)
                dat.x[-1, 5 + idx] = np.linspace(np.min(dat.x[:, 5:]),
                            np.max(dat.x[:, 5:]), dat.x.shape[1] - 5)
            else: # y-axis is in log
                xdat = np.abs(dat.x[:, 5:] - np.array(x_opt, copy=False))
                idx = np.argsort(xdat[-2, :])
                idx2 = np.argsort(idx)
                xdat[-1, idx] = np.logspace(np.log10(np.min(abs(xdat[xdat!=0]))),
                            np.log10(np.max(np.abs(xdat))),
                            dat.x.shape[1] - 5)
        else:
            minxend = 0
        self._enter_plotting()
        if x_opt is not None:  # TODO: differentate neg and pos?
            semilogy(dat.x[:, iabscissa], abs(xdat), '-')
        else:
            plot(dat.x[:, iabscissa], dat.x[:, 5:], '-')
        hold(True)
        grid(True)
        ax = array(axis())
        # ax[1] = max(minxend, ax[1])
        axis(ax)
        ax[1] -= 1e-6  # to prevent last x-tick annotation, probably superfluous
        if dat.x.shape[1] < 100:
            yy = np.linspace(ax[2] + 1e-6, ax[3] - 1e-6, dat.x.shape[1] - 5)
            # yyl = np.sort(dat.x[-1,5:])
            if x_opt is not None:
                # semilogy([dat.x[-1, iabscissa], ax[1]], [abs(dat.x[-1, 5:]), yy[idx2]], 'k-')  # line from last data point
                semilogy(np.dot(dat.x[-2, iabscissa], [1, 1]),
                         array([ax[2] * (1+1e-6), ax[3] / (1+1e-6)]), 'k-')
            else:
                # plot([dat.x[-1, iabscissa], ax[1]], [dat.x[-1,5:], yy[idx2]], 'k-') # line from last data point
                plot(np.dot(dat.x[-2, iabscissa], [1, 1]),
                     array([ax[2] + 1e-6, ax[3] - 1e-6]), 'k-')
            # plot(array([dat.x[-1, iabscissa], ax[1]]),
            #      reshape(array([dat.x[-1,5:], yy[idx2]]).flatten(), (2,4)), '-k')
            for i in rglen(idx):
                # TODOqqq: annotate phenotypic value!?
                # text(ax[1], yy[i], 'x(' + str(idx[i]) + ')=' + str(dat.x[-2,5+idx[i]]))

                text(dat.x[-1, iabscissa], dat.x[-1, 5 + i]
                            if x_opt is None else np.abs(xdat[-1, i]),
                     ('x(' + str(i) + ')=' if annotations is None
                        else str(i) + ':' + annotations[i] + "=")
                     + str(dat.x[-2, 5 + i]))
        i = 2  # find smallest i where iteration count differs (in case the same row appears twice)
        while i < len(dat.f) and dat.f[-i][0] == dat.f[-1][0]:
            i += 1
        title('Object Variables (' +
                (remark + ', ' if remark is not None else '') +
                str(dat.x.shape[1] - 5) + '-D, popsize~' +
                (str(int((dat.f[-1][1] - dat.f[-i][1]) / (dat.f[-1][0] - dat.f[-i][0])))
                    if len(dat.f.T[0]) > 1 and dat.f[-1][0] > dat.f[-i][0] else 'NA')
                + ')')
        self._finalize_plotting()
505565
function parse_PtgMemArea(blob, length, opts) {
	var type = (blob.read_shift(1) >>> 5) & 0x03;
	blob.l += (opts && opts.biff == 2 ? 3 : 4);
	var cce = blob.read_shift(opts && opts.biff == 2 ? 1 : 2);
	return [type, cce];
}
1280603
def _plot_data_to_ax(
    data_all,
    ax1,
    e_unit=None,
    sed=True,
    ylabel=None,
    ulim_opts={},
    errorbar_opts={},
):
    

    if e_unit is None:
        e_unit = data_all["energy"].unit

    f_unit, sedf = sed_conversion(
        data_all["energy"], data_all["flux"].unit, sed
    )

    if "group" not in data_all.keys():
        data_all["group"] = np.zeros(len(data_all))

    groups = np.unique(data_all["group"])

    for g in groups:
        data = data_all[np.where(data_all["group"] == g)]
        _, sedfg = sed_conversion(data["energy"], data["flux"].unit, sed)

        # wrap around color and marker cycles
        color = color_cycle[int(g) % len(color_cycle)]
        marker = marker_cycle[int(g) % len(marker_cycle)]

        ul = data["ul"]
        notul = ~ul

        # Hack to show y errors compatible with 0 in loglog plot
        yerr_lo = data["flux_error_lo"][notul]
        y = data["flux"][notul].to(yerr_lo.unit)
        bad_err = np.where((y - yerr_lo) <= 0.0)
        yerr_lo[bad_err] = y[bad_err] * (1.0 - 1e-7)
        yerr = u.Quantity((yerr_lo, data["flux_error_hi"][notul]))
        xerr = u.Quantity((data["energy_error_lo"], data["energy_error_hi"]))

        opts = dict(
            zorder=100,
            marker=marker,
            ls="",
            elinewidth=2,
            capsize=0,
            mec=color,
            mew=0.1,
            ms=5,
            color=color,
        )
        opts.update(**errorbar_opts)

        ax1.errorbar(
            data["energy"][notul].to(e_unit).value,
            (data["flux"][notul] * sedfg[notul]).to(f_unit).value,
            yerr=(yerr * sedfg[notul]).to(f_unit).value,
            xerr=xerr[:, notul].to(e_unit).value,
            **opts
        )

        if np.any(ul):
            if "elinewidth" in errorbar_opts:
                ulim_opts["elinewidth"] = errorbar_opts["elinewidth"]

            _plot_ulims(
                ax1,
                data["energy"][ul].to(e_unit).value,
                (data["flux"][ul] * sedfg[ul]).to(f_unit).value,
                (xerr[:, ul]).to(e_unit).value,
                color,
                **ulim_opts
            )

    ax1.set_xscale("log")
    ax1.set_yscale("log")
    xmin = 10 ** np.floor(
        np.log10(
            np.min(data["energy"] - data["energy_error_lo"]).to(e_unit).value
        )
    )
    xmax = 10 ** np.ceil(
        np.log10(
            np.max(data["energy"] + data["energy_error_hi"]).to(e_unit).value
        )
    )
    ax1.set_xlim(xmin, xmax)
    # avoid autoscaling to errorbars to 0
    notul = ~data_all["ul"]
    if np.any(data_all["flux_error_lo"][notul] >= data_all["flux"][notul]):
        elo = (data_all["flux"][notul] * sedf[notul]).to(f_unit).value - (
            data_all["flux_error_lo"][notul] * sedf[notul]
        ).to(f_unit).value
        gooderr = np.where(
            data_all["flux_error_lo"][notul] < data_all["flux"][notul]
        )
        ymin = 10 ** np.floor(np.log10(np.min(elo[gooderr])))
        ax1.set_ylim(bottom=ymin)

    if ylabel is None:
        if sed:
            ax1.set_ylabel(
                r"$E^2\mathrm{{d}}N/\mathrm{{d}}E$"
                " [{0}]".format(u.Unit(f_unit).to_string("latex_inline"))
            )
        else:
            ax1.set_ylabel(
                r"$\mathrm{{d}}N/\mathrm{{d}}E$"
                " [{0}]".format(u.Unit(f_unit).to_string("latex_inline"))
            )
    else:
        ax1.set_ylabel(ylabel)
720042
def extract_ms(self, x, y):
        
        ms = np.zeros(self.conf.ndims)
        ms[self.expl_dims] = x
        ms[self.inf_dims] = y
        return ms[self.conf.m_dims], ms[self.conf.s_dims]
1309565
def parse_bismark_report(self, report, regexes):
        
        parsed_data = {}
        for k, r in regexes.items():
            r_search = re.search(r, report, re.MULTILINE)
            if r_search:
                try:
                    parsed_data[k] = float(r_search.group(1))
                except ValueError:
                    parsed_data[k] = r_search.group(1) # NaN
        if len(parsed_data) == 0: return None
        return parsed_data
643921
def parse_plotPCA(self):
        """"""
        self.deeptools_plotPCAData = dict()
        for f in self.find_log_files('deeptools/plotPCAData', filehandles=False):
            parsed_data = self.parsePlotPCAData(f)
            for k, v in parsed_data.items():
                if k in self.deeptools_plotPCAData:
                    log.warning("Replacing duplicate sample {}.".format(k))
                self.deeptools_plotPCAData[k] = v
            if len(parsed_data) > 0:
                self.add_data_source(f, section='plotPCA')

        if len(self.deeptools_plotPCAData) > 0:
            config = {
                'id': 'deeptools_pca_plot',
                'title': 'deeptools: PCA Plot',
                'xlab': 'PC1',
                'ylab': 'PC2',
                'tt_label': 'PC1 {point.x:.2f}: PC2 {point.y:.2f}',
            }
            data = dict()
            for s_name in self.deeptools_plotPCAData:
                try:
                    data[s_name] = {'x': self.deeptools_plotPCAData[s_name][1], 'y': self.deeptools_plotPCAData[s_name][2]}
                except KeyError:
                    pass
            if len(data) == 0:
                log.debug('No valid data for PCA plot')
                return None

            self.add_section(
                name="PCA plot",
                anchor="deeptools_pca",
                description="PCA plot with the top two principal components calculated based on genome-wide distribution of sequence reads",
                plot=scatter.plot(data, config)
            )

        return len(self.deeptools_plotPCAData)
1125416
def mstmap(args):
    
    from jcvi.assembly.geneticmap import MSTMatrix

    p = OptionParser(mstmap.__doc__)
    p.add_option("--population_type", default="RIL6",
                 help="Type of population, possible values are DH and RILd")
    p.add_option("--missing_threshold", default=.5,
                 help="Missing threshold, .25 excludes any marker with >25% missing")
    p.set_outfile()
    opts, args = p.parse_args(args)

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

    lmd, = args
    fp = open(lmd)
    next(fp)  # Header
    table = {"0": "-", "1": "A", "2": "B", "3": "X"}
    mh = ["locus_name"] + fp.next().split()[4:]
    genotypes = []
    for row in fp:
        atoms = row.split()
        chr, pos, ref, alt = atoms[:4]
        locus_name = ".".join((chr, pos))
        codes = [table[x] for x in atoms[4:]]
        genotypes.append([locus_name] + codes)

    mm = MSTMatrix(genotypes, mh, opts.population_type, opts.missing_threshold)
    mm.write(opts.outfile, header=True)
887017
def macs_filtered_reads_plot(self):
        
        data = dict()
        req_cats = ['control_fragments_total', 'control_fragments_after_filtering', 'treatment_fragments_total', 'treatment_fragments_after_filtering']
        for s_name, d in self.macs_data.items():
            if all([c in d for c in req_cats]):
                data['{}: Control'.format(s_name)] = dict()
                data['{}: Treatment'.format(s_name)] = dict()
                data['{}: Control'.format(s_name)]['fragments_filtered'] = d['control_fragments_total'] - d['control_fragments_after_filtering']
                data['{}: Control'.format(s_name)]['fragments_not_filtered'] = d['control_fragments_after_filtering']
                data['{}: Treatment'.format(s_name)]['fragments_filtered'] = d['treatment_fragments_total'] - d['treatment_fragments_after_filtering']
                data['{}: Treatment'.format(s_name)]['fragments_not_filtered'] = d['treatment_fragments_after_filtering']

        # Check that we have something to plot
        if len(data) == 0:
            return

        # Specify the order of the different possible categories
        keys = OrderedDict()
        keys['fragments_not_filtered'] = { 'color': '#437BB1', 'name': 'Remaining fragments' }
        keys['fragments_filtered'] =     { 'color': '#B1084C', 'name': 'Filtered fragments' }

        # Config for the plot
        pconfig = {
            'id': 'macs2_filtered',
            'title': 'MACS2: Filtered Fragments',
            'ylab': '# Fragments',
            'cpswitch_counts_label': 'Number of Fragments',
            'hide_zero_cats': False
        }

        self.add_section(
            plot = bargraph.plot(data, keys, pconfig)
        )
10602
def parse_tophat_log (self, raw_data):
        

        if 'Aligned pairs' in raw_data:
            # Paired end data
            regexes = {
                'overall_aligned_percent': r"([\d\.]+)% overall read mapping rate.",
                'concordant_aligned_percent': r"([\d\.]+)% concordant pair alignment rate.",
                'aligned_total': r"Aligned pairs:\s+(\d+)",
                'aligned_multimap': r"Aligned pairs:\s+\d+\n\s+of these:\s+(\d+)",
                'aligned_discordant': r"(\d+) \([\s\d\.]+%\) are discordant alignments",
                'total_reads': r"[Rr]eads:\n\s+Input\s*:\s+(\d+)",
            }
        else:
            # Single end data
            regexes = {
                'total_reads': r"[Rr]eads:\n\s+Input\s*:\s+(\d+)",
                'aligned_total': r"Mapped\s*:\s+(\d+)",
                'aligned_multimap': r"of these\s*:\s+(\d+)",
                'overall_aligned_percent': r"([\d\.]+)% overall read mapping rate.",
            }

        parsed_data = {}
        for k, r in regexes.items():
            r_search = re.search(r, raw_data, re.MULTILINE)
            if r_search:
                parsed_data[k] = float(r_search.group(1))
        if len(parsed_data) == 0: return None
        parsed_data['concordant_aligned_percent'] = parsed_data.get('concordant_aligned_percent', 0)
        parsed_data['aligned_total'] = parsed_data.get('aligned_total', 0)
        parsed_data['aligned_multimap'] = parsed_data.get('aligned_multimap', 0)
        parsed_data['aligned_discordant'] = parsed_data.get('aligned_discordant', 0)
        parsed_data['unaligned_total'] = parsed_data['total_reads'] - parsed_data['aligned_total']
        parsed_data['aligned_not_multimapped_discordant'] = parsed_data['aligned_total'] - parsed_data['aligned_multimap'] - parsed_data['aligned_discordant']
        return parsed_data
710262
def plot_irm(fignum, B, M, title):
    
    rpars = {}
    Mnorm = []
    backfield = 0
    X, Y = [], []
    for k in range(len(B)):
        if M[k] < 0:
            break
    if k <= 5:
        kmin = 0
    else:
        kmin = k - 5
    for k in range(kmin, k + 1):
        X.append(B[k])
        if B[k] < 0:
            backfield = 1
        Y.append(M[k])
    if backfield == 1:
        poly = np.polyfit(X, Y, 1)
        if poly[0] != 0:
            bcr = (old_div(-poly[1], poly[0]))
        else:
            bcr = 0
        rpars['remanence_mr_moment'] = '%8.3e' % (M[0])
        rpars['remanence_bcr'] = '%8.3e' % (-bcr)
        rpars['magic_method_codes'] = 'LP-BCR-BF'
        if M[0] != 0:
            for m in M:
                Mnorm.append(old_div(m, M[0]))  # normalize to unity Msat
            title = title + ':' + '%8.3e' % (M[0])
    else:
        if M[-1] != 0:
            for m in M:
                Mnorm.append(old_div(m, M[-1]))  # normalize to unity Msat
            title = title + ':' + '%8.3e' % (M[-1])
# do plots if desired
    if fignum != 0 and M[0] != 0:  # skip plot for fignum = 0
        plt.figure(num=fignum)
        plt.clf()
        if not isServer:
            plt.figtext(.02, .01, version_num)
        plt.plot(B, Mnorm)
        plt.axhline(0, color='k')
        plt.axvline(0, color='k')
        plt.xlabel('B (T)')
        plt.ylabel('M/Mr')
        plt.title(title)
        if backfield == 1:
            plt.scatter([bcr], [0], marker='s', c='b')
            bounds = plt.axis()
            n1 = 'Bcr: ' + '%8.2e' % (-bcr) + ' T'
            plt.figtext(.2, .5, n1)
            n2 = 'Mr: ' + '%8.2e' % (M[0]) + ' Am^2'
            plt.figtext(.2, .45, n2)
    elif fignum != 0:
        plt.figure(num=fignum)
        # plt.clf()
        if not isServer:
            plt.figtext(.02, .01, version_num)
        print('M[0]=0,  skipping specimen')
    return rpars
768221
def _parse_track_dim(self,d1,interp=True,phys=False):
        """"""
        if interp: interpStr= 'interpolated'
        else: interpStr= ''
        if d1.lower() == 'x':
            tx= self.__dict__['_%sObsTrackXY' % interpStr][:,0]
        elif d1.lower() == 'y':
            tx= self.__dict__['_%sObsTrackXY' % interpStr][:,1]
        elif d1.lower() == 'z':
            tx= self.__dict__['_%sObsTrackXY' % interpStr][:,2]
        elif d1.lower() == 'r':
            tx= self.__dict__['_%sObsTrack' % interpStr][:,0]
        elif d1.lower() == 'phi':
            tx= self.__dict__['_%sObsTrack' % interpStr][:,5]
        elif d1.lower() == 'vx': 
            tx= self.__dict__['_%sObsTrackXY' % interpStr][:,3]
        elif d1.lower() == 'vy':
            tx= self.__dict__['_%sObsTrackXY' % interpStr][:,4]
        elif d1.lower() == 'vz':
            tx= self.__dict__['_%sObsTrackXY' % interpStr][:,5]
        elif d1.lower() == 'vr':
            tx= self.__dict__['_%sObsTrack' % interpStr][:,1]
        elif d1.lower() == 'vt':
            tx= self.__dict__['_%sObsTrack' % interpStr][:,2]
        elif d1.lower() == 'll':
            tx= self.__dict__['_%sObsTrackLB' % interpStr][:,0]
        elif d1.lower() == 'bb':
            tx= self.__dict__['_%sObsTrackLB' % interpStr][:,1]
        elif d1.lower() == 'dist':
            tx= self.__dict__['_%sObsTrackLB' % interpStr][:,2]
        elif d1.lower() == 'pmll':
            tx= self.__dict__['_%sObsTrackLB' % interpStr][:,4]
        elif d1.lower() == 'pmbb':
            tx= self.__dict__['_%sObsTrackLB' % interpStr][:,5]
        elif d1.lower() == 'vlos':
            tx= self.__dict__['_%sObsTrackLB' % interpStr][:,3]
        if phys and (d1.lower() == 'x' or d1.lower() == 'y' \
                         or d1.lower() == 'z' or d1.lower() == 'r'):
            tx= copy.copy(tx)
            tx*= self._ro
        if phys and (d1.lower() == 'vx' or d1.lower() == 'vy' \
                         or d1.lower() == 'vz' or d1.lower() == 'vr' \
                         or d1.lower() == 'vt'):
            tx= copy.copy(tx)
            tx*= self._vo
        return tx
340413
def bismark_strand_chart (self):
        

        # Specify the order of the different possible categories
        keys = OrderedDict()
        keys['strand_ob']   = { 'name': 'Original bottom strand' }
        keys['strand_ctob'] = { 'name': 'Complementary to original bottom strand' }
        keys['strand_ctot'] = { 'name': 'Complementary to original top strand' }
        keys['strand_ot']   = { 'name': 'Original top strand' }

        # See if we have any directional samples
        directional = 0
        d_mode = ''
        for sn in self.bismark_data['alignment'].values():
            if 'strand_directional' in sn.keys():
                directional += 1
        if directional == len(self.bismark_data['alignment']):
            keys.pop('strand_ctob', None)
            keys.pop('strand_ctot', None)
            d_mode = 'All samples were run with <code>--directional</code> mode; alignments to complementary strands (CTOT, CTOB) were ignored.'
        elif directional > 0:
            d_mode = '{} samples were run with <code>--directional</code> mode; alignments to complementary strands (CTOT, CTOB) were ignored.'.format(directional)

        # Config for the plot
        config = {
            'id': 'bismark_strand_alignment',
            'title': 'Bismark: Alignment to Individual Bisulfite Strands',
            'ylab': '% Reads',
            'cpswitch_c_active': False,
            'cpswitch_counts_label': 'Number of Reads'
        }

        self.add_section (
            name = 'Strand Alignment',
            anchor = 'bismark-strands',
            description = d_mode,
            plot = bargraph.plot(self.bismark_data['alignment'], keys, config)
        )
1342055
def plot_probe(mea, ax=None, xlim=None, ylim=None, color_currents=False, top=None, bottom=None,
               cmap='viridis', type='shank', alpha_elec=0.7, alpha_prb=0.3):
    '''
    
    

    '''


    if ax is None:
        fig = plt.figure()
        ax = fig.add_subplot(111)

    n_elec = mea.positions.shape[0]
    elec_size = mea.size

    if mea.type == 'mea':
        mea_pos = np.array([np.dot(mea.positions, mea.main_axes[0]), np.dot(mea.positions, mea.main_axes[1])]).T

        min_x, max_x = [np.min(np.dot(mea.positions, mea.main_axes[0])),
                        np.max(np.dot(mea.positions, mea.main_axes[0]))]
        center_x = (min_x + max_x)/2.
        min_y, max_y = [np.min(np.dot(mea.positions, mea.main_axes[1])),
                        np.max(np.dot(mea.positions, mea.main_axes[1]))]
        center_y = (min_y + max_y)/2.

        if type == 'shank':
            if top is None:
                probe_height = 200
            else:
                probe_height = top
            probe_top = max_y + probe_height
            if bottom is None:
                probe_bottom = min_y - probe_height
            else:
                probe_bottom = min_y - bottom
            probe_corner = min_y - 0.1*probe_height
            probe_left = min_x - 0.1*probe_height
            probe_right = max_x + 0.1*probe_height

            verts = [
                (min_x - 2*elec_size, probe_top),  # left, bottom
                (min_x - 2*elec_size, probe_corner),  # left, top
                (center_x, probe_bottom),  # right, top
                (max_x + 2*elec_size, probe_corner),  # right, bottom
                (max_x + 2*elec_size, probe_top),
                (min_x - 2 * elec_size, max_y + 2 * elec_size) # ignored
            ]

            codes = [Path.MOVETO,
                     Path.LINETO,
                     Path.LINETO,
                     Path.LINETO,
                     Path.LINETO,
                     Path.CLOSEPOLY,
                     ]

        elif type == 'planar':
            probe_top = max_y + 2 * elec_size
            probe_bottom = min_y - 2 * elec_size
            probe_left = min_x - 2 * elec_size
            probe_right = max_x + 2 * elec_size

            verts = [
                (min_x - 2 * elec_size, max_y + 2 * elec_size),  # left, bottom
                (min_x - 2 * elec_size, min_y - 2 * elec_size),  # left, top
                (max_x + 2 * elec_size, min_y - 2 * elec_size),  # right, bottom
                (max_x + 2 * elec_size, max_y + 2 * elec_size),  # ignored
                (max_x + 2 * elec_size, max_y + 2 * elec_size)   # ignored
            ]

            codes = [Path.MOVETO,
                     Path.LINETO,
                     Path.LINETO,
                     Path.LINETO,
                     Path.CLOSEPOLY,
                     ]
        else:
            raise AttributeError("'type' can be 'shank' or 'planar'")

        path = Path(verts, codes)

        patch = patches.PathPatch(path, facecolor='green', edgecolor='k', lw=0.5, alpha=alpha_prb)
        ax.add_patch(patch)

        if color_currents:
            norm_curr = mea.currents / np.max(np.abs(mea.currents))
            colormap = plt.get_cmap(cmap)
            elec_colors = colormap(norm_curr)
        else:
            elec_colors = ['orange'] * mea.number_electrodes

        if mea.shape == 'square':
            for e in range(n_elec):
                elec = patches.Rectangle((mea_pos[e, 0] - elec_size, mea_pos[e, 1] - elec_size), 2*elec_size,
                                         2*elec_size, alpha=alpha_elec, facecolor=elec_colors[e],
                                         edgecolor=[0.3, 0.3, 0.3], lw=0.5)

                ax.add_patch(elec)
        elif mea.shape == 'circle':
            for e in range(n_elec):
                elec = patches.Circle((mea_pos[e, 0], mea_pos[e, 1]), elec_size, alpha=alpha_elec,
                                      facecolor=elec_colors[e], edgecolor=[0.3, 0.3, 0.3], lw=0.5)

                ax.add_patch(elec)
    else:
        raise NotImplementedError('Wire type plotting not implemented')

    ax.axis('equal')

    if xlim:
        ax.set_xlim(xlim)
    else:
        ax.set_xlim(probe_left - 5 * elec_size, probe_right + 5 * elec_size)
    if ylim:
        ax.set_ylim(ylim)
    else:
        ax.set_ylim(probe_bottom - 5 * elec_size, probe_top + 5 * elec_size)

    return ax
1072714
def bbt_simple_plot(self):
        

        # First, sum the different types of alignment counts
        data = OrderedDict()
        cats = OrderedDict()
        for s_name in self.bbt_data:
            data[s_name] = OrderedDict()
            for org in self.bbt_data[s_name]:
                data[s_name][org] = self.bbt_data[s_name][org]['hits'] - self.bbt_data[s_name][org]['shared']
                if org not in cats and org != 'multiMatch' and org != 'noMatch':
                    if org.lower().endswith('.fa'):
                        cname = org[:-3]
                    elif org.lower().endswith('.fasta'):
                        cname = org[:-6]
                    else:
                        cname = org
                    cats[org] = { 'name': cname }

        pconfig = {
            'id': 'biobloom_tools',
            'title': 'BioBloom Tools: Alignment counts per species',
            'ylab': 'Number of hits',
            'hide_zero_cats': False
        }
        cats['multiMatch'] = { 'name': 'Multiple Genomes', 'color': '#820000' }
        cats['noMatch'] = { 'name': 'No Match', 'color': '#cccccc' }

        return bargraph.plot(data, cats, pconfig)