Examples with TwoBitBufferedReferenceSequenceFile au.edu.wehi.idsv.picard.TwoBitBufferedReferenceSequenceFile used on opensource projects

Example 1 with TwoBitBufferedReferenceSequenceFile

use of au.edu.wehi.idsv.picard.TwoBitBufferedReferenceSequenceFile in project gridss by PapenfussLab.

the class ByReadNameSinglePassSamProgram method makeItSo.

public static void makeItSo(final File input, final File referenceSequence, final boolean assumeSorted, final long stopAfter, final Collection<ByReadNameSinglePassSamProgram> programs) throws FileNotFoundException {
    // Setup the standard inputs
    IOUtil.assertFileIsReadable(input);
    SamReader in = SamReaderFactory.makeDefault().referenceSequence(referenceSequence).open(input);
    // Optionally load up the reference sequence and double check sequence dictionaries
    final ReferenceLookup lookup;
    if (referenceSequence == null) {
        lookup = null;
    } else {
        IOUtil.assertFileIsReadable(referenceSequence);
        lookup = new TwoBitBufferedReferenceSequenceFile(new IndexedFastaSequenceFile(referenceSequence));
        if (!in.getFileHeader().getSequenceDictionary().isEmpty()) {
            SequenceUtil.assertSequenceDictionariesEqual(in.getFileHeader().getSequenceDictionary(), lookup.getSequenceDictionary());
        }
    }
    // Check on the sort order of the BAM file
    final SortOrder sort = in.getFileHeader().getSortOrder();
    if (sort != SortOrder.queryname) {
        if (assumeSorted) {
            log.warn("File reports sort order '" + sort + "', assuming it's queryname sorted anyway.");
        } else {
            throw new PicardException("File " + input.getAbsolutePath() + " should be queryname sorted but " + "the header says the sort order is " + sort + ". If you believe the file " + "to be queryname sorted you may pass ASSUME_SORTED=true");
        }
    }
    for (final ByReadNameSinglePassSamProgram program : programs) {
        program.setReference(lookup);
        program.setup(in.getFileHeader(), input);
    }
    final ProgressLogger progress = new ProgressLogger(log);
    final SAMRecordIterator rawit = in.iterator();
    final CloseableIterator<SAMRecord> it = new AsyncBufferedIterator<SAMRecord>(rawit, "ByReadNameSinglePassSamProgram " + input.getName());
    try {
        List<SAMRecord> currentRecords = new ArrayList<>();
        String currentReadName = null;
        while (it.hasNext()) {
            SAMRecord r = it.next();
            String readname = r.getReadName();
            // if read name we have to just treat it as a single read
            if (readname == null || !readname.equals(currentReadName)) {
                if (currentRecords.size() > 0) {
                    for (final ByReadNameSinglePassSamProgram program : programs) {
                        program.acceptFragment(currentRecords, lookup);
                    }
                }
                currentRecords.clear();
                currentReadName = readname;
                if (stopAfter > 0 && progress.getCount() >= stopAfter) {
                    break;
                }
            }
            currentRecords.add(r);
            progress.record(r);
        }
        if (currentRecords.size() > 0) {
            for (final ByReadNameSinglePassSamProgram program : programs) {
                program.acceptFragment(currentRecords, lookup);
            }
        }
    } finally {
        CloserUtil.close(it);
        CloserUtil.close(rawit);
        CloserUtil.close(in);
    }
    for (final ByReadNameSinglePassSamProgram program : programs) {
        program.finish();
    }
}

Example 2 with TwoBitBufferedReferenceSequenceFile

use of au.edu.wehi.idsv.picard.TwoBitBufferedReferenceSequenceFile in project gridss by PapenfussLab.

the class BreakpointHomologyTest method should_not_exceed_contig_bounds.

@Test
public void should_not_exceed_contig_bounds() {
    InMemoryReferenceSequenceFile underlying = new InMemoryReferenceSequenceFile(new String[] { "0", "1" }, new byte[][] { B("CCCCCCCCCCC"), B("AAAAAAAAAAA") });
    // BreakpointHomology bh0 = BreakpointHomology.calculate(ref, new BreakpointSummary(0, FWD, 1, 1, BWD, 20), "", 50, 50);
    TwoBitBufferedReferenceSequenceFile ref = new TwoBitBufferedReferenceSequenceFile(underlying);
    for (int b1pos = 1; b1pos <= ref.getSequenceDictionary().getSequences().get(0).getSequenceLength(); b1pos++) {
        for (int b2pos = 1; b2pos <= ref.getSequenceDictionary().getSequences().get(0).getSequenceLength(); b2pos++) {
            for (BreakendDirection b1dir : BreakendDirection.values()) {
                for (BreakendDirection b2dir : BreakendDirection.values()) {
                    for (int maxBreakendLength = 1; maxBreakendLength < ref.getSequenceDictionary().getSequences().get(0).getSequenceLength() + 2; maxBreakendLength++) {
                        for (int margin = 0; margin < ref.getSequenceDictionary().getSequences().get(0).getSequenceLength() + 2; margin++) {
                            BreakpointHomology bh = BreakpointHomology.calculate(ref, new BreakpointSummary(0, b1dir, b1pos, 1, b2dir, b2pos), "", maxBreakendLength, margin);
                            assertEquals(0, bh.getLocalHomologyLength());
                            assertEquals(0, bh.getRemoteHomologyLength());
                        }
                    }
                }
            }
        }
    }
}