diff --git a/java/src/org/broadinstitute/sting/utils/BasicPileup.java b/java/src/org/broadinstitute/sting/utils/BasicPileup.java index a86cedf3a..763ca3bf4 100755 --- a/java/src/org/broadinstitute/sting/utils/BasicPileup.java +++ b/java/src/org/broadinstitute/sting/utils/BasicPileup.java @@ -140,6 +140,56 @@ abstract public class BasicPileup implements Pileup { return quals2.toString(); } + public static double[][] probDistPileup( List reads, List offsets ) { + double[][] dist = new double[reads.size()][4]; + + for (int readIndex = 0; readIndex < dist.length; readIndex++) { + SAMRecord read = reads.get(readIndex); + + String bases = read.getReadString(); + int offset = offsets.get(readIndex); + + int bestBaseIndex = BaseUtils.simpleBaseToBaseIndex(bases.charAt(offset)); + + if (bestBaseIndex >= 0 && bestBaseIndex < 4) { + dist[readIndex][bestBaseIndex] = QualityUtils.qualToProb(read.getBaseQualities()[offset]); + + byte[] sqs = (byte[]) read.getAttribute("SQ"); + if (sqs != null) { + int secondBestBaseIndex = QualityUtils.compressedQualityToBaseIndex(sqs[offset]); + dist[readIndex][secondBestBaseIndex] = (1.0 - dist[readIndex][bestBaseIndex]); + } else { + for (int baseIndex = 0; baseIndex < 4; baseIndex++) { + if (baseIndex != bestBaseIndex) { + dist[readIndex][baseIndex] = (dist[readIndex][bestBaseIndex]/3.0); + } + } + } + } else { + for (int baseIndex = 0; baseIndex < 4; baseIndex++) { + dist[readIndex][baseIndex] = 0.25; + } + } + } + + return dist; + } + + public static String probDistPileupAsString( List reads, List offsets ) { + double[][] dist = probDistPileup(reads, offsets); + + String distString = ""; + for (int readIndex = 0; readIndex < dist.length; readIndex++) { + distString += "[ "; + for (int baseIndex = 0; baseIndex < 4; baseIndex++) { + distString += String.format("%3.3f ", dist[readIndex][baseIndex]); + } + distString += "]\n"; + } + + return distString; + } + public static String pileupDiff(final Pileup a, final Pileup b) { return pileupDiff(a,b,true); diff --git a/java/src/org/broadinstitute/sting/utils/ReadBackedPileup.java b/java/src/org/broadinstitute/sting/utils/ReadBackedPileup.java index a34831b29..05cd0c723 100755 --- a/java/src/org/broadinstitute/sting/utils/ReadBackedPileup.java +++ b/java/src/org/broadinstitute/sting/utils/ReadBackedPileup.java @@ -69,7 +69,11 @@ public class ReadBackedPileup extends BasicPileup { public String getSecondaryQualPileup() { return secondaryQualPileupAsString(reads, offsets); } - + + public String getProbDistPileup() { + return probDistPileupAsString(reads, offsets); + } + public String getPileupString() { return String.format("%s: %s %s %s", getLocation(), getRef(), getBases(), getQuals());