Java Code Examples for htsjdk.samtools.reference.IndexedFastaSequenceFile#getSubsequenceAt()

The following examples show how to use htsjdk.samtools.reference.IndexedFastaSequenceFile#getSubsequenceAt() . You can vote up the ones you like or vote down the ones you don't like, and go to the original project or source file by following the links above each example. You may check out the related API usage on the sidebar.
Example 1
Source File: ReferenceResource.java    From VarDictJava with MIT License 6 votes vote down vote up 
/**
 * Method convert fasta data to array contains header and nucleotide bases.
 * @param fasta path to fasta file
 * @param chr chromosome name of region
 * @param start start position of region
 * @param end end position of region
 * @return array of nucleotide bases in the region of fasta
 */
public String[] retrieveSubSeq(String fasta, String chr, int start, int end) {
    try {
        IndexedFastaSequenceFile idx = fetchFasta(fasta);
        ReferenceSequence seq = idx.getSubsequenceAt(chr, start, end);
        byte[] bases = seq.getBases();
        return new String[]{">" + chr + ":" + start + "-" + end, bases != null ? new String(bases) : ""};
    } catch (SAMException e){
        if (UNABLE_FIND_CONTIG.matcher(e.getMessage()).find()){
            throw new WrongFastaOrBamException(chr, e);
        } else if (WRONG_START_OR_END.matcher(e.getMessage()).find()){
            throw new RegionBoundariesException(chr, start, end, e);
        } else {
            throw e;
        }
    }
}
Example 2
Source File: RefGenomeEnrichedSomaticVariantFactory.java    From hmftools with GNU General Public License v3.0 5 votes vote down vote up 
private static void addTrinucleotideContext(@NotNull final ImmutableEnrichedSomaticVariant.Builder builder,
        @NotNull final SomaticVariant variant, @NotNull IndexedFastaSequenceFile reference) {
    final int chromosomeLength = reference.getSequenceDictionary().getSequence(variant.chromosome()).getSequenceLength();
    if (variant.position() < chromosomeLength) {
        final ReferenceSequence sequence =
                reference.getSubsequenceAt(variant.chromosome(), Math.max(1, variant.position() - 1), variant.position() + 1);
        builder.trinucleotideContext(sequence.getBaseString());
    } else {
        LOGGER.warn("Requested ref sequence beyond contig length! variant = {}", variant);
    }
}