Retrieve the reference genome file location from galaxy configuration.
483341
def get_genome_ref(genome_build, aligner, galaxy_base):
    
    ref_files = dict(
            bowtie = "bowtie_indices.loc",
            bwa = "bwa_index.loc",
            samtools = "sam_fa_indices.loc",
            maq = "bowtie_indices.loc")
    remap_fns = dict(
            maq = _remap_to_maq
            )
    out_info = []
    for ref_get in [aligner, "samtools"]:
        ref_file = os.path.join(galaxy_base, "tool-data", ref_files[ref_get])
        with open(ref_file) as in_handle:
            for line in in_handle:
                if not line.startswith("#"):
                    parts = line.strip().split()
                    if parts[0] == "index":
                        parts = parts[1:]
                    if parts[0] == genome_build:
                        out_info.append(parts[-1])
                        break
        try:
            out_info[-1] = remap_fns[ref_get](out_info[-1])
        except KeyError:
            pass
        except IndexError:
            raise IndexError("Genome %s not found in %s" % (genome_build,
                ref_file))

    if len(out_info) != 2:
        raise ValueError("Did not find genome reference for %s %s" %
                (genome_build, aligner))
    else:
        return tuple(out_info)
1372479
def _get_galaxy_tool_info(galaxy_base):
    
    ini_file = os.path.join(galaxy_base, "universe_wsgi.ini")
    info = {"tool_data_table_config_path": os.path.join(galaxy_base, "tool_data_table_conf.xml"),
            "tool_data_path": os.path.join(galaxy_base, "tool-data")}
    config = configparser.ConfigParser()
    config.read(ini_file)
    if "app:main" in config.sections():
        for option in config.options("app:main"):
            if option in info:
                info[option] = os.path.join(galaxy_base, config.get("app:main", option))
    return info
373653
def list_installed_genomes(genome_dir=None):
    
    if not genome_dir:
        genome_dir = config.get("genome_dir", None)
    if not genome_dir:
        raise norns.exceptions.ConfigError("Please provide or configure a genome_dir")

    return [f for f in os.listdir(genome_dir) if 
            _is_genome_dir(genome_dir + "/" + f)]
1146875
def _file_to_folder(gi, fname, sample_info, libitems, library, folder):
    
    full_name = os.path.join(folder["name"], os.path.basename(fname))

    # Handle VCF: Galaxy reports VCF files without the gzip extension
    file_type = "vcf_bgzip" if full_name.endswith(".vcf.gz") else "auto"
    if full_name.endswith(".vcf.gz"):
        full_name = full_name.replace(".vcf.gz", ".vcf")

    for item in libitems:
        if item["name"] == full_name:
            return item
    logger.info("Uploading to Galaxy library '%s': %s" % (library.name, full_name))
    return gi.libraries.upload_from_galaxy_filesystem(str(library.id), fname, folder_id=str(folder["id"]),
                                                      link_data_only="link_to_files",
                                                      dbkey=sample_info["genome_build"],
                                                      file_type=file_type,
                                                      roles=str(library.roles) if library.roles else None)
446938
def _findLocation(self, reference_name, start, end):
        
        try:
            # TODO - sequence_annotations does not have build?
            return self._locationMap['hg19'][reference_name][start][end]
        except:
            return None
717408
def location(self):
        
        try:
            lat, lng = self["metadata"]["latitude"], self["metadata"]["longitude"]
        except KeyError:
            return None

        if not lat or not lng:
            return None

        return lat, lng
146830
def _get_G2P(data):
    
    G2P_file = os.path.realpath(tz.get_in(("genome_resources", "variation", "genotype2phenotype"), data))
    args = ["--plugin", "G2P,file:%s" % (G2P_file)]
    if G2P_file:
        return args
    else:
        return []
914834
def get_ref_bedtool(ref_file, config, chrom=None):
    
    broad_runner = broad.runner_from_path("picard", config)
    ref_dict = broad_runner.run_fn("picard_index_ref", ref_file)
    ref_lines = []
    with pysam.Samfile(ref_dict, "r") as ref_sam:
        for sq in ref_sam.header["SQ"]:
            if not chrom or sq["SN"] == chrom:
                ref_lines.append("%s\t%s\t%s" % (sq["SN"], 0, sq["LN"]))
    return pybedtools.BedTool("\n".join(ref_lines), from_string=True)
1008378
def get_runinfo(galaxy_url, galaxy_apikey, run_folder, storedir):
    
    galaxy_api = GalaxyApiAccess(galaxy_url, galaxy_apikey)
    fc_name, fc_date = flowcell.parse_dirname(run_folder)
    galaxy_info = galaxy_api.run_details(fc_name, fc_date)
    if "error" in galaxy_info:
        return galaxy_info
    if not galaxy_info["run_name"].startswith(fc_date) and not galaxy_info["run_name"].endswith(fc_name):
        raise ValueError("Galaxy NGLIMS information %s does not match flowcell %s %s" %
                         (galaxy_info["run_name"], fc_date, fc_name))
    ldetails = _flatten_lane_details(galaxy_info)
    out = []
    for item in ldetails:
        # Do uploads for all non-controls
        if item["description"] != "control" or item["project_name"] != "control":
            item["upload"] = {"method": "galaxy", "run_id": galaxy_info["run_id"],
                              "fc_name": fc_name, "fc_date": fc_date,
                              "dir": storedir,
                              "galaxy_url": galaxy_url, "galaxy_api_key": galaxy_apikey}
            for k in ["lab_association", "private_libs", "researcher", "researcher_id", "sample_id",
                      "galaxy_library", "galaxy_role"]:
                item["upload"][k] = item.pop(k, "")
        out.append(item)
    return out
406143
def _get_chrbands(self, limit, taxon):
        
        model = Model(self.graph)
        line_counter = 0
        myfile = '/'.join((self.rawdir, self.files[taxon]['file']))
        LOG.info("Processing Chr bands from FILE: %s", myfile)
        geno = Genotype(self.graph)

        # build the organism's genome from the taxon
        genome_label = self.files[taxon]['genome_label']
        taxon_id = 'NCBITaxon:' + taxon

        # add the taxon as a class.  adding the class label elsewhere
        model.addClassToGraph(taxon_id, None)
        model.addSynonym(taxon_id, genome_label)

        genome_id = geno.makeGenomeID(taxon_id)
        geno.addGenome(taxon_id, genome_label)
        model.addOWLPropertyClassRestriction(
            genome_id, self.globaltt['in taxon'], taxon_id)

        placed_scaffold_pattern = r'chr(\d+|X|Y|Z|W|MT|M)'
        # currently unused patterns
        # unlocalized_scaffold_pattern = placed_scaffold_pattern + r'_(\w+)_random'
        # unplaced_scaffold_pattern = r'chrUn_(\w+)'

        col = ['chrom', 'start', 'stop', 'band', 'rtype']
        with gzip.open(myfile, 'rb') as reader:
            for line in reader:
                line_counter += 1
                # skip comments
                line = line.decode().strip()
                if line[0] == '#':
                    continue
                # chr13	4500000	10000000	p12	stalk
                row = line.split('\t')
                chrom = row[col.index('chrom')]
                band = row[col.index('band')]
                rtype = row[col.index('rtype')]
                # NOTE
                # some less-finished genomes have placed and unplaced scaffolds
                # * Placed scaffolds:
                #    Scaffold has an oriented location within a chromosome.
                # * Unlocalized scaffolds:
                #     scaffold 's chromosome  is known,
                #     scaffold's position, orientation or both is not known.
                # *Unplaced scaffolds:
                #   it is not known which chromosome the scaffold belongs to.
                # find out if the thing is a full on chromosome, or a scaffold:
                # ex: unlocalized scaffold: chr10_KL568008v1_random
                # ex: unplaced scaffold: chrUn_AABR07022428v1

                mch = re.match(placed_scaffold_pattern+r'$', chrom)
                if mch is not None and len(mch.groups()) == 1:
                    # the chromosome is the first match of the pattern
                    # chrom = m.group(1)  # TODO unused
                    pass
                else:
                    # let's skip over anything that isn't a placed_scaffold
                    LOG.info("Skipping non-placed chromosome %s", chrom)
                    continue
                # the chrom class, taxon as the reference
                cclassid = makeChromID(chrom, taxon, 'CHR')

                # add the chromosome as a class
                geno.addChromosomeClass(chrom, taxon_id, genome_label)
                model.addOWLPropertyClassRestriction(
                    cclassid, self.globaltt['member of'], genome_id)

                # add the band(region) as a class
                maplocclass_id = cclassid+band
                maplocclass_label = makeChromLabel(chrom+band, genome_label)
                if band is not None and band.strip() != '':
                    region_type_id = self.map_type_of_region(rtype)
                    model.addClassToGraph(
                        maplocclass_id, maplocclass_label,
                        region_type_id)
                else:
                    region_type_id = self.globaltt['chromosome']
                # add the staining intensity of the band
                if re.match(r'g(neg|pos|var)', rtype):
                    if region_type_id in [
                            self.globaltt['chromosome_band'],
                            self.globaltt['chromosome_subband']]:
                        stain_type = self.resolve(rtype)
                        if stain_type is not None:
                            model.addOWLPropertyClassRestriction(
                                maplocclass_id,
                                self.globaltt['has_sequence_attribute'],
                                self.resolve(rtype))
                    else:
                        # usually happens if it's a chromosome because
                        # they don't actually have banding info
                        LOG.info("feature type %s != chr band", region_type_id)
                else:
                    LOG.warning('staining type not found: %s', rtype)

                # get the parent bands, and make them unique
                parents = list(self.make_parent_bands(band, set()))
                # alphabetical sort will put them in smallest to biggest
                parents.sort(reverse=True)

                # print("PARENTS of", maplocclass_id, "=", parents)
                # add the parents to the graph, in hierarchical order
                # TODO this is somewhat inefficient due to
                # re-adding upper-level nodes when iterating over the file
                for prnt in parents:
                    parent = prnt.strip()
                    if parent is None or parent == "":
                        continue
                    pclassid = cclassid + parent  # class chr parts
                    pclass_label = makeChromLabel(chrom + parent, genome_label)
                    rti = getChrPartTypeByNotation(parent, self.graph)
                    model.addClassToGraph(pclassid, pclass_label, rti)

                    # for canonical chromosomes,
                    # then the subbands are subsequences of the full band
                    # add the subsequence stuff as restrictions

                    if prnt != parents[-1]:
                        grandparent = 1 + parents.index(prnt)
                        pid = cclassid + parents[grandparent]   # the instance
                        model.addOWLPropertyClassRestriction(
                            pclassid, self.globaltt['is subsequence of'], pid)
                        model.addOWLPropertyClassRestriction(
                            pid, self.globaltt['has subsequence'], pclassid)
                    else:
                        # add the last one (p or q usually)
                        # as attached to the chromosome
                        model.addOWLPropertyClassRestriction(
                            pclassid, self.globaltt['is subsequence of'], cclassid)
                        model.addOWLPropertyClassRestriction(
                            cclassid, self.globaltt['has subsequence'], pclassid)

                # connect the band here to the first one in the parent list
                if len(parents) > 0:
                    model.addOWLPropertyClassRestriction(
                        maplocclass_id, self.globaltt['is subsequence of'],
                        cclassid + parents[0])
                    model.addOWLPropertyClassRestriction(
                        cclassid + parents[0], self.globaltt['has subsequence'],
                        maplocclass_id)

                if limit is not None and line_counter > limit:
                    break

        # TODO figure out the staining intensities for the encompassing bands

        return
254543
def get_coords(data):
    
    for category, vtypes in [("LOH", {"LOSS", "HETEROZYGOSITY"}),
                             ("amplification", {"AMPLIFICATION"})]:
        out = tz.get_in([category, dd.get_genome_build(data)], _COORDS, {})
        priority_file = dd.get_svprioritize(data)
        if priority_file:
            if os.path.basename(priority_file).find("civic") >= 0:
                for chrom, start, end, gene in _civic_regions(priority_file, vtypes, dd.get_disease(data)):
                    out[gene] = (chrom, start, end)
            elif os.path.basename(priority_file).find(".bed") >= 0:
                for line in utils.open_gzipsafe(priority_file):
                    parts = line.strip().split("\t")
                    if len(parts) >= 4:
                        chrom, start, end, gene = parts[:4]
                        out[gene] = (chrom, int(start), int(end))
        yield category, out
467246
def write_versions(dirs, items):
    
    genomes = {}
    for d in items:
        genomes[d["genome_build"]] = d.get("reference", {}).get("versions")
    out_file = _get_out_file(dirs)
    found_versions = False
    if genomes and out_file:
        with open(out_file, "w") as out_handle:
            writer = csv.writer(out_handle)
            writer.writerow(["genome", "resource", "version"])
            for genome, version_file in genomes.items():
                if not version_file:
                    genome_dir = install.get_genome_dir(genome, dirs.get("galaxy"), items[0])
                    if genome_dir:
                        version_file = os.path.join(genome_dir, "versions.csv")
                if version_file and os.path.exists(version_file):
                    found_versions = True
                    with open(version_file) as in_handle:
                        reader = csv.reader(in_handle)
                        for parts in reader:
                            if len(parts) >= 2:
                                resource, version = parts[:2]
                                writer.writerow([genome, resource, version])
    if found_versions:
        return out_file
1109902
def ref_file_from_bam(bam_file, data):
    
    new_ref = os.path.join(utils.safe_makedir(os.path.join(dd.get_work_dir(data), "inputs", "ref")),
                           "%s-subset.fa" % dd.get_genome_build(data))
    if not utils.file_exists(new_ref):
        with file_transaction(data, new_ref) as tx_out_file:
            contig_file = "%s-contigs.txt" % utils.splitext_plus(new_ref)[0]
            with open(contig_file, "w") as out_handle:
                for contig in [x.contig for x in idxstats(bam_file, data) if x.contig != "*"]:
                    out_handle.write("%s\n" % contig)
            cmd = "seqtk subseq -l 100 %s %s > %s" % (dd.get_ref_file(data), contig_file, tx_out_file)
            do.run(cmd, "Subset %s to BAM file contigs" % dd.get_genome_build(data))
    ref.fasta_idx(new_ref, data["config"])
    runner = broad.runner_from_path("picard", data["config"])
    runner.run_fn("picard_index_ref", new_ref)
    return {"base": new_ref}
737519
def reference_index(self):
        """"""
        if self._db_location:
            ref_indices = glob.glob(os.path.join(self._db_location, "*", self._REF_INDEX))
            if ref_indices:
                return ref_indices[0]
1213816
def package_data_location():
    '''

    '''
    pkg_dir = os.path.split(__file__)[0]
    pkg_data_dir = os.path.join(pkg_dir,'data')
    return pkg_data_dir
98860
def _get_data_versions(data):
    
    genome_dir = install.get_genome_dir(data["genome_build"], data["dirs"].get("galaxy"), data)
    if genome_dir:
        version_file = os.path.join(genome_dir, "versions.csv")
        if version_file and os.path.exists(version_file):
            return version_file
    return None
559135
def parse_gc_content(self):
        """"""
        # Find and parse GC content:
        for f in self.find_log_files('homer/GCcontent', filehandles=True):
            # Get the s_name from the parent directory
            s_name = os.path.basename(f['root'])
            s_name = self.clean_s_name(s_name, f['root'])
            parsed_data = self.parse_twoCol_file(f)
            if parsed_data is not None:
                if s_name in self.tagdir_data['GCcontent']:
                    log.debug("Duplicate GCcontent sample log found! Overwriting: {}".format(s_name))

                self.add_data_source(f, s_name, section='GCcontent')
                self.tagdir_data['GCcontent'][s_name] = parsed_data

        ## get esimated genome content distribution:
        for f in self.find_log_files('homer/genomeGCcontent', filehandles=True):
            parsed_data = self.parse_twoCol_file(f)
            if parsed_data is not None:
                if s_name + "_genome" in self.tagdir_data['GCcontent']:
                    log.debug("Duplicate genome GCcontent sample log found! Overwriting: {}".format(s_name+ "_genome"))

                self.add_data_source(f, s_name + "_genome", section='GCcontent')
                self.tagdir_data['GCcontent'][s_name + "_genome"] = parsed_data

        self.tagdir_data['GCcontent'] = self.ignore_samples(self.tagdir_data['GCcontent'])

        if len(self.tagdir_data['GCcontent']) > 0:
            self.add_section (
                name = 'Per Sequence GC Content',
                anchor = 'homer_per_sequence_gc_content',
                description = 'This plot shows the distribution of GC content.',
                plot = self.GCcontent_plot()
            )