use of htsjdk.samtools.fastq.FastqRecord in project gridss by PapenfussLab.
the class StructuralVariationCallBuilderTest method breakpoint_assembly_should_be_written.
@Test
@Ignore("Issue #17: Output breakpoint assembly sequences")
public void breakpoint_assembly_should_be_written() {
ProcessingContext pc = getContext();
String seq = "CATTAATCGCAAGAGCGGGTTGTATTCGcCGCCAAGTCAGCTGAAGCACCATTACCCGAtCAAAACATATCAGAAATGATTGACGTATCACAAGCCGGATTTTGTTTACAGCCTGTCTTATATCCTGAATAACGCACCGCCTATTCG";
int anchor = 78;
SAMRecord ass = AssemblyFactory.createAnchoredBreakend(getContext(), AES(), new SequentialIdGenerator("asm"), FWD, null, 6, 78, anchor, B(seq), B(40, seq.length()));
SingleReadEvidence ae = incorporateRealignment(AES(), ass, ImmutableList.of(withQual(B(40, seq.length() - anchor), Read(2, anchor + 1, String.format("%dM", seq.length() - anchor)))[0]));
StructuralVariationCallBuilder builder = new StructuralVariationCallBuilder(pc, (VariantContextDirectedEvidence) new IdsvVariantContextBuilder(getContext()) {
{
breakpoint(new BreakpointSummary(6, FWD, 78, 2, BWD, 79), "");
phredScore(50);
}
}.make());
builder.addEvidence(ae);
VariantContextDirectedEvidence e = builder.make();
FastqRecord fr = e.getBreakendAssemblyFastq();
assertEquals(e.getEvidenceID() + "#" + Integer.toString(anchor), fr.getReadName());
assertEquals(seq, fr.getReadString());
assertEquals(fr.getReadName(), e.getAttribute(VcfSvConstants.BREAKPOINT_ID_KEY));
}
use of htsjdk.samtools.fastq.FastqRecord in project gridss by PapenfussLab.
the class SplitReadRealigner method processAlignmentRecord.
private void processAlignmentRecord(StreamingAligner aligner, SplitReadFastqExtractor recursiveExtractor, Map<String, SplitReadRealignmentInfo> realignments, SAMFileWriter writer) throws IOException {
SAMRecord supp = aligner.getAlignment();
String lookupkey = SplitReadIdentificationHelper.getOriginatingAlignmentUniqueName(supp);
SplitReadRealignmentInfo info = realignments.get(lookupkey);
if (supp.getSupplementaryAlignmentFlag() || supp.getNotPrimaryAlignmentFlag()) {
// only consider the best mapping location reported by the aligner
// log.debug(String.format("%s: ignoring supp alignment", supp.getReadName()));
} else {
assert (info.outstandingRealignments > 0);
info.outstandingRealignments--;
if (!supp.getReadUnmappedFlag()) {
info.realignments.add(supp);
List<FastqRecord> nestedRealignments = recursiveExtractor.extract(supp);
for (FastqRecord fq : nestedRealignments) {
aligner.asyncAlign(fq);
info.outstandingRealignments++;
// log.debug(String.format("%s: performing nested realignment. %d realignments now outstanding", info.originatingRecord.getReadName(), info.outstandingRealignments));
}
}
// all splits identified
if (info.outstandingRealignments == 0) {
if (info.realignments.size() > 0) {
SplitReadIdentificationHelper.convertToSplitRead(info.originatingRecord, info.realignments);
}
writer.addAlignment(info.originatingRecord);
for (SAMRecord sar : info.realignments) {
writer.addAlignment(sar);
}
// log.debug(String.format("%s: %d supp alignments found.", info.originatingRecord.getReadName(), info.realignments.size()));
realignments.remove(lookupkey);
} else {
// log.debug(String.format("%s: %d outstanding alignments", info.originatingRecord.getReadName(), info.outstandingRealignments));
}
}
}
use of htsjdk.samtools.fastq.FastqRecord in project gridss by PapenfussLab.
the class SplitReadIdentificationHelper method getSplitReadRealignments.
/**
* Extract the unaligned portions of the read requiring realignment to identify split reads
*
* @param r alignment record
* @param recordIsPartialAlignment true if the record is the result of aligning FastqRecords
* from a previous call to getSplitReadRealignments()
* @return bases requiring alignment to identify split reads
*/
public static List<FastqRecord> getSplitReadRealignments(SAMRecord r, boolean recordIsPartialAlignment, EvidenceIdentifierGenerator eidgen) {
int startClipLength = SAMRecordUtil.getStartSoftClipLength(r);
int endClipLength = SAMRecordUtil.getEndSoftClipLength(r);
if (startClipLength + endClipLength == 0 || r.getReadUnmappedFlag()) {
return Collections.emptyList();
}
int offset;
String name;
if (recordIsPartialAlignment) {
offset = getRealignmentFirstAlignedBaseReadOffset(r);
name = getOriginatingAlignmentUniqueName(r);
} else {
offset = 0;
name = eidgen.getAlignmentUniqueName(r);
}
List<FastqRecord> list = new ArrayList<>(2);
if (startClipLength > 0) {
byte[] startbases = Arrays.copyOfRange(r.getReadBases(), 0, startClipLength);
byte[] startqual = Arrays.copyOfRange(r.getBaseQualities(), 0, startClipLength);
int startOffset = offset;
if (r.getReadNegativeStrandFlag()) {
SequenceUtil.reverseComplement(startbases);
ArrayUtils.reverse(startqual);
// 5432
// SS
startOffset = startOffset + r.getReadLength() - startClipLength;
}
FastqRecord start = new FastqRecord(getSplitAlignmentFastqName(name, startOffset), new String(startbases, StandardCharsets.US_ASCII), "", SAMUtils.phredToFastq(startqual));
list.add(start);
}
if (endClipLength > 0) {
byte[] endbases = Arrays.copyOfRange(r.getReadBases(), r.getReadLength() - endClipLength, r.getReadLength());
byte[] endqual = Arrays.copyOfRange(r.getBaseQualities(), r.getReadLength() - endClipLength, r.getReadLength());
int endOffset = offset + r.getReadLength() - endClipLength;
if (r.getReadNegativeStrandFlag()) {
SequenceUtil.reverseComplement(endbases);
ArrayUtils.reverse(endqual);
// 5432
// SS
endOffset = offset;
}
FastqRecord end = new FastqRecord(getSplitAlignmentFastqName(name, endOffset), new String(endbases, StandardCharsets.US_ASCII), "", SAMUtils.phredToFastq(endqual));
list.add(end);
}
return list;
}
use of htsjdk.samtools.fastq.FastqRecord in project jvarkit by lindenb.
the class FastqRevComp method run.
private void run(FastqReader r, PrintStream out) {
String s;
long nRec = 0L;
r.setValidationStringency(ValidationStringency.LENIENT);
while (r.hasNext()) {
if (++nRec % 1E6 == 0) {
LOG.info("N-Reads:" + nRec);
}
FastqRecord fastq = r.next();
out.print(FastqConstants.SEQUENCE_HEADER);
out.println(fastq.getReadName());
s = fastq.getReadString();
if (// interleaced
(this.only_R2 && nRec % 2 == 1) || (this.only_R1 && nRec % 2 == 0)) {
out.print(s);
} else {
for (int i = s.length() - 1; i >= 0; i--) {
out.print(AcidNucleics.complement(s.charAt(i)));
}
}
out.println();
out.print(FastqConstants.QUALITY_HEADER);
s = fastq.getBaseQualityHeader();
if (s != null)
out.print(s);
out.println();
s = fastq.getBaseQualityString();
if (// interleaced
(this.only_R2 && nRec % 2 == 1) || (this.only_R1 && nRec % 2 == 0)) {
out.print(s);
} else {
for (int i = s.length() - 1; i >= 0; i--) {
out.print(s.charAt(i));
}
}
out.println();
if (out.checkError())
break;
}
out.flush();
LOG.info("Done. N-Reads:" + nRec);
}
use of htsjdk.samtools.fastq.FastqRecord in project jvarkit by lindenb.
the class FastqToFasta method run.
private void run(FastqReader r, PrintStream out) {
int wsp = 0;
long nRec = 0L;
r.setValidationStringency(ValidationStringency.LENIENT);
while (r.hasNext()) {
if (++nRec % 1E6 == 0) {
LOG.info("N-Reads:" + nRec);
}
FastqRecord fastq = r.next();
out.print(">");
if (!trim_after_space || (wsp = fastq.getReadName().indexOf(' ')) == -1) {
out.println(fastq.getReadName());
} else {
out.println(fastq.getReadName().substring(0, wsp));
}
int readLen = fastq.getReadString().length();
int i = 0;
while (i < readLen) {
int end = Math.min(i + fastaLineLen, readLen);
out.println(fastq.getReadString().substring(i, end));
i = end;
}
if (out.checkError())
break;
}
out.flush();
LOG.info("Done. N-Reads:" + nRec);
}
Aggregations