Java SAMSequenceDictionary.getSequences方法代码示例

import htsjdk.samtools.SAMSequenceDictionary; //导入方法依赖的package包/类
private void writeHeader(final SAMFileHeader header) {
	binaryCodec.writeBytes("BAM\001".getBytes(Charset.forName("UTF8")));

	final Writer sw = new StringWriter();
	new SAMTextHeaderCodec().encode(sw, header);

	binaryCodec.writeString(sw.toString(), true, false);

	final SAMSequenceDictionary dict = header.getSequenceDictionary();

	binaryCodec.writeInt(dict.size());
	for (final SAMSequenceRecord rec : dict.getSequences()) {
		binaryCodec.writeString(rec.getSequenceName(), true, true);
		binaryCodec.writeInt   (rec.getSequenceLength());
	}
}
 

import htsjdk.samtools.SAMSequenceDictionary; //导入方法依赖的package包/类
/**
 * Returns a compact String representation of the sequence dictionary it's passed
 *
 * The format of the returned String is:
 * [ contig1Name(length: contig1Length) contig2Name(length: contig2Length) ... ]
 *
 * @param dict a non-null SAMSequenceDictionary
 * @return A String containing all of the contig names and lengths from the sequence dictionary it's passed
 */
public static String getDictionaryAsString( final SAMSequenceDictionary dict ) {
    Utils.nonNull(dict, "Sequence dictionary must be non-null");

    StringBuilder s = new StringBuilder("[ ");

    for ( SAMSequenceRecord dictionaryEntry : dict.getSequences() ) {
        s.append(dictionaryEntry.getSequenceName());
        s.append("(length:");
        s.append(dictionaryEntry.getSequenceLength());
        s.append(") ");
    }

    s.append("]");

    return s.toString();
}
 

import htsjdk.samtools.SAMSequenceDictionary; //导入方法依赖的package包/类
/**
 * Create a genome loc parser based on seqDict with the specified level of validation
 * @param seqDict the sequence dictionary to use when creating genome locs
 * @param validationLevel how much validation should we do of the genome locs at runtime? Purely for testing purposes
 */
protected GenomeLocParser(SAMSequenceDictionary seqDict, final ValidationLevel validationLevel) {
    Utils.nonNull(validationLevel, "validation level cannot be null");
    if (seqDict == null) { // we couldn't load the reference dictionary
        //logger.info("Failed to load reference dictionary, falling back to lexicographic order for contigs");
        throw new CommandLineException("Failed to load reference dictionary");
    }

    this.validationLevel = validationLevel;
    this.contigInfo = new MRUCachingSAMSequenceDictionary(seqDict);
    if ( logger.isDebugEnabled() ) {
        logger.debug(String.format("Prepared reference sequence contig dictionary"));
        for (SAMSequenceRecord contig : seqDict.getSequences()) {
            logger.debug(String.format(" %s (%d bp)", contig.getSequenceName(), contig.getSequenceLength()));
        }
    }
}
 

import htsjdk.samtools.SAMSequenceDictionary; //导入方法依赖的package包/类
/**
 * Compares two non-null sequence dictionaries using sequence index, name, and length only.
 * Less stringent than {@link SAMSequenceDictionary#isSameDictionary}.
 */
public static boolean isSameDictionary(final SAMSequenceDictionary dictionary1,
                                       final SAMSequenceDictionary dictionary2) {
    Utils.nonNull(dictionary1);
    Utils.nonNull(dictionary2);
    if (dictionary1 == dictionary2) {
        return true;
    }

    final Iterator<SAMSequenceRecord> dictionary1Sequences = dictionary1.getSequences().iterator();
    for (final SAMSequenceRecord dictionary2Sequence : dictionary2.getSequences()) {
        if (!dictionary1Sequences.hasNext()) {
            return false;
        } else {
            final SAMSequenceRecord dictionary1Sequence = dictionary1Sequences.next();
            if (!isSameSequence(dictionary1Sequence, dictionary2Sequence)) {
                return false;
            }
        }
    }
    return !dictionary1Sequences.hasNext();
}
 

import htsjdk.samtools.SAMSequenceDictionary; //导入方法依赖的package包/类
/**
 * Creates a random reference and writes it in FASTA format into a {@link Writer}.
 * @param out the output writer.
 * @param dict the dictionary indicating the number of contigs and their lengths.
 * @param basesPerLine number of base to print in each line of the output FASTA file.
 *
 * @throws IOException if such an exception was thrown while accessing and writing into the temporal file.
 * @throws IllegalArgumentException if {@code dict} is {@code null}, or {@code out } is {@code null}
 *    or {@code basesPerLine} is 0 or negative.
 */
public void nextFasta(final Writer out, final SAMSequenceDictionary dict, final int basesPerLine)
        throws IOException {
    Utils.nonNull(out);
    Utils.nonNull(dict);
    ParamUtils.isPositive(basesPerLine, "number of base per line must be strictly positive: " + basesPerLine);
    final byte[] buffer = new byte[basesPerLine];
    final String lineSeparator = System.lineSeparator();
    for (final SAMSequenceRecord sequence : dict.getSequences()) {
        int pendingBases = sequence.getSequenceLength();
        out.append(">").append(sequence.getSequenceName()).append(lineSeparator);
        while (pendingBases > 0) {
            final int lineLength = pendingBases < basesPerLine ? pendingBases : basesPerLine;
            nextBases(buffer, 0, lineLength);
            out.append(new String(buffer, 0, lineLength)).append(lineSeparator);
            pendingBases -= lineLength;
        }
    }
}
 

import htsjdk.samtools.SAMSequenceDictionary; //导入方法依赖的package包/类
@Test(groups = "sv")
public void miniRefTest() throws IOException {
    final JavaSparkContext ctx = SparkContextFactory.getTestSparkContext();
    final ReferenceMultiSource ref = new ReferenceMultiSource(
            REFERENCE_FILE_NAME, ReferenceWindowFunctions.IDENTITY_FUNCTION);
    final SAMSequenceDictionary dict = ref.getReferenceSequenceDictionary(null);
    if ( dict == null ) throw new GATKException("No reference dictionary available.");

    final Map<SVKmer, Long> kmerMap = new LinkedHashMap<>();
    for ( final SAMSequenceRecord rec : dict.getSequences() ) {
        final SimpleInterval interval = new SimpleInterval(rec.getSequenceName(), 1, rec.getSequenceLength());
        final byte[] bases = ref.getReferenceBases(interval).getBases();
        SVKmerizer.canonicalStream(bases, KMER_SIZE, new SVKmerLong())
                .forEach(kmer -> kmerMap.put(kmer, kmerMap.getOrDefault(kmer, 0L) + 1));
    }
    kmerMap.entrySet().removeIf( x -> x.getValue() <= FindBadGenomicKmersSpark.MAX_KMER_FREQ);

    final List<SVKmer> badKmers =
            FindBadGenomicKmersSpark.findBadGenomicKmers(ctx, KMER_SIZE, Integer.MAX_VALUE, ref, null);
    final Set<SVKmer> badKmerSet = new HashSet<>(badKmers);
    Assert.assertEquals(badKmers.size(), badKmerSet.size());
    Assert.assertEquals(badKmerSet, kmerMap.keySet());
}
 

import htsjdk.samtools.SAMSequenceDictionary; //导入方法依赖的package包/类
private void assertFastaContent(final Path path, final boolean withDescriptions, final SAMSequenceDictionary dictionary, final int defaultBpl,
                                final Map<String, byte[]> bases, final Map<String, Integer> basesPerLine)
        throws IOException {
    try (final BufferedReader reader = new BufferedReader(new InputStreamReader(path.getFileSystem().provider().newInputStream(path)))) {
        for (final SAMSequenceRecord sequence : dictionary.getSequences()) {
            final String description = String.format("index=%d\tlength=%d",
                    dictionary.getSequenceIndex(sequence.getSequenceName()), sequence.getSequenceLength());
            final String expectedHeader =
                    FastaReferenceWriter.HEADER_START_CHAR + sequence.getSequenceName()
                            + ((withDescriptions) ? FastaReferenceWriter.HEADER_NAME_AND_DESCRIPTION_SEPARATOR + description : "");
            Assert.assertEquals(reader.readLine(), expectedHeader);
            final byte[] expectedBases = bases.get(sequence.getSequenceName());
            final int bpl_ = basesPerLine.get(sequence.getSequenceName());
            final int bpl = bpl_ < 0 ? (defaultBpl < 0 ? FastaReferenceWriter.DEFAULT_BASES_PER_LINE : defaultBpl) : bpl_;
            int offset = 0;
            while (offset < expectedBases.length) {
                final int expectedLength = Math.min(expectedBases.length - offset, bpl);
                final byte[] expectedBaseLine = SequenceUtil.upperCase(Arrays.copyOfRange(expectedBases, offset, offset + expectedLength));
                final byte[] actualBaseLine = SequenceUtil.upperCase(reader.readLine().getBytes());
                Assert.assertEquals(actualBaseLine, expectedBaseLine);
                offset += expectedLength;
            }
        }
    }
}
 

import htsjdk.samtools.SAMSequenceDictionary; //导入方法依赖的package包/类
private void assertFastaIndexContent(final Path path, final Path indexPath, final SAMSequenceDictionary dictionary,
                                     final Map<String, byte[]> bases)
        throws IOException {
    final FastaSequenceIndex index = new FastaSequenceIndex(indexPath);
    final IndexedFastaSequenceFile indexedFasta = new IndexedFastaSequenceFile(path, index);
    for (final SAMSequenceRecord sequence : dictionary.getSequences()) {
        final String name = sequence.getSequenceName();
        final int length = sequence.getSequenceLength();
        final ReferenceSequence start = indexedFasta.getSubsequenceAt(name, 1, Math.min(length, 30));
        final ReferenceSequence end = indexedFasta.getSubsequenceAt(name, Math.max(1, length - 29), length);
        final int middlePos = Math.max(1, Math.min(length, length / 2));
        final ReferenceSequence middle = indexedFasta.getSubsequenceAt(name, middlePos, Math.min(middlePos + 29, length));
        Assert.assertEquals(start.getBases(), Arrays.copyOfRange(bases.get(name), 0, start.length()));
        Assert.assertEquals(end.getBases(), Arrays.copyOfRange(bases.get(name), Math.max(0, length - 30), length));
        Assert.assertEquals(middle.getBases(), Arrays.copyOfRange(bases.get(name), middlePos - 1, middlePos - 1 + middle.length()));
    }
}
 

import htsjdk.samtools.SAMSequenceDictionary; //导入方法依赖的package包/类
/**
 * create a list of genomic locations, given a reference sequence
 *
 * @param dict the sequence dictionary to create a collection from
 *
 * @return the GenomeLocSet of all references sequences as GenomeLoc's
 */
public static GenomeLocSortedSet createSetFromSequenceDictionary(final SAMSequenceDictionary dict) {
    final GenomeLocParser parser = new GenomeLocParser(dict);
    final GenomeLocSortedSet returnSortedSet = new GenomeLocSortedSet(parser);
    for ( final SAMSequenceRecord sequence : dict.getSequences() ) {
        returnSortedSet.add(parser.createOverEntireContig(sequence.getSequenceName()));
    }
    return returnSortedSet;
}
 

import htsjdk.samtools.SAMSequenceDictionary; //导入方法依赖的package包/类
public static String prettyPrintSequenceRecords ( SAMSequenceDictionary sequenceDictionary ) {
    String[] sequenceRecordNames = new String[sequenceDictionary.size()];
    int sequenceRecordIndex = 0;
    for (SAMSequenceRecord sequenceRecord : sequenceDictionary.getSequences())
        sequenceRecordNames[sequenceRecordIndex++] = sequenceRecord.getSequenceName();
    return Arrays.deepToString(sequenceRecordNames);
}
 

import htsjdk.samtools.SAMSequenceDictionary; //导入方法依赖的package包/类
public static void setSequenceDictionary(Configuration conf, SAMSequenceDictionary dict) throws IOException, URISyntaxException {
    int counter = 0;
    for(SAMSequenceRecord seq : dict.getSequences()) {
        conf.set(dictionarySequenceName + counter, seq.getSequenceName());
        conf.setInt(dictionarySequenceLength + counter, seq.getSequenceLength());
        counter++;
    }
    conf.setInt(dictionaryCount, counter);
}
 

import htsjdk.samtools.SAMSequenceDictionary; //导入方法依赖的package包/类
/**
 * Helper function to print out a sequence dictionary
 */
private void printDictionary(String name, SAMSequenceDictionary dict) {
    log.info(name);
    for (final SAMSequenceRecord contig : dict.getSequences()) {
        log.info("  SN=%s LN=%d%n", contig.getSequenceName(), contig.getSequenceLength());
    }
}
 

import htsjdk.samtools.SAMSequenceDictionary; //导入方法依赖的package包/类
private static List<SAMSequenceRecord> readSequences(final File vcf) {
    final VCFFileReader reader = new VCFFileReader(vcf);
    final VCFHeader header = reader.getFileHeader();
    final SAMSequenceDictionary dict = header.getSequenceDictionary();
    reader.close();
    return dict.getSequences();
}
 

import htsjdk.samtools.SAMSequenceDictionary; //导入方法依赖的package包/类
/**
 * Utility function that tests whether dict1's set of contigs is a superset of dict2's
 *
 * @param dict1 first sequence dictionary
 * @param dict2 second sequence dictionary
 * @return true if dict1's set of contigs supersets dict2's
 */
private static boolean supersets( SAMSequenceDictionary dict1, SAMSequenceDictionary dict2 ) {
    // Cannot rely on SAMSequenceRecord.equals() as it's too strict (takes extended attributes into account).
    for ( final SAMSequenceRecord dict2Record : dict2.getSequences() ) {
        final SAMSequenceRecord dict1Record = dict1.getSequence(dict2Record.getSequenceName());
        if ( dict1Record == null || ! sequenceRecordsAreEquivalent(dict2Record, dict1Record) ) {
            return false;
        }
    }

    return true;
}
 

import htsjdk.samtools.SAMSequenceDictionary; //导入方法依赖的package包/类
/**
 * A very simple (and naive) algorithm to determine (1) if the dict is a human reference (hg18, hg19, b36, or b37) and if it's
 * lexicographically sorted.  Works by matching lengths of the static chr1, chr10, and chr2, and then if these
 * are all matched, requiring that the order be chr1, chr2, chr10.
 *
 * @param dict
 * @return
 */
private static boolean nonCanonicalHumanContigOrder(SAMSequenceDictionary dict) {
    SAMSequenceRecord chr1 = null, chr2 = null, chr10 = null;
    for ( SAMSequenceRecord elt : dict.getSequences() ) {
        if ( isHumanSeqRecord(elt, CHR1_HG18, CHR1_HG19, CHR1_B36, CHR1_B37) ) chr1 = elt;
        if ( isHumanSeqRecord(elt, CHR2_HG18, CHR2_HG19, CHR2_B36, CHR2_B37) ) chr2 = elt;
        if ( isHumanSeqRecord(elt, CHR10_HG18, CHR10_HG19, CHR10_B36, CHR10_B37) ) chr10 = elt;
    }
    if ( chr1 != null  && chr2 != null && chr10 != null) {
        return ! ( chr1.getSequenceIndex() < chr2.getSequenceIndex() && chr2.getSequenceIndex() < chr10.getSequenceIndex() );
    }

    return false;
}