diff --git a/src/org/broad/igv/PreferenceManager.java b/src/org/broad/igv/PreferenceManager.java index 32e7a53d20..b7ed8b8a49 100644 --- a/src/org/broad/igv/PreferenceManager.java +++ b/src/org/broad/igv/PreferenceManager.java @@ -26,6 +26,7 @@ import org.broad.igv.ui.color.PaletteColorTable; import org.broad.igv.ui.util.PropertyManager; import org.broad.igv.util.HttpUtils; +import org.broad.igv.sam.AlignmentTrack.ShadeBasesOption; import java.awt.*; import java.io.File; @@ -68,7 +69,7 @@ public class PreferenceManager implements PropertyManager { public static final String SAM_SHOW_CENTER_LINE = "SAM.SHOW_CENTER_LINE"; public static final String SAM_SHOW_REF_SEQ = "SAM.SHOW_REF_SEQ"; public static final String SAM_SHOW_COV_TRACK = "SAM.SHOW_COV_TRACK"; - public static final String SAM_SHADE_BASE_QUALITY = "SAM.SHADE_BASE_QUALITY"; + public static final String SAM_SHADE_BASES = "SAM.SHADE_BASE_QUALITY"; public static final String SAM_BASE_QUALITY_MIN = "SAM.BASE_QUALITY_MIN"; public static final String SAM_BASE_QUALITY_MAX = "SAM.BASE_QUALITY_MAX"; public static final String SAM_FILTER_ALIGNMENTS = "SAM.FILTER_ALIGNMENTS"; @@ -939,7 +940,7 @@ private void initDefaultValues() { defaultValues.put(SAM_SHOW_REF_SEQ, "false"); defaultValues.put(SAM_SHOW_CENTER_LINE, "true"); defaultValues.put(SAM_SHOW_COV_TRACK, "true"); - defaultValues.put(SAM_SHADE_BASE_QUALITY, "true"); + defaultValues.put(SAM_SHADE_BASES, ShadeBasesOption.QUALITY.toString()); defaultValues.put(SAM_FILTER_ALIGNMENTS, "false"); defaultValues.put(SAM_FILTER_FAILED_READS, "true"); defaultValues.put(SAM_DOWNSAMPLE_READS, "true"); diff --git a/src/org/broad/igv/sam/AbstractAlignment.java b/src/org/broad/igv/sam/AbstractAlignment.java index dd2c2d8e17..0028a34010 100644 --- a/src/org/broad/igv/sam/AbstractAlignment.java +++ b/src/org/broad/igv/sam/AbstractAlignment.java @@ -121,26 +121,28 @@ public byte getPhred(double position) { } return 0; } + private void bufAppendFlowSignals(AlignmentBlock block, StringBuffer buf, int offset) { if (block.hasFlowSignals()) { // flow signals int i, j, n = 0; - short[][] flowSignalContext = block.getFlowSignalContext(offset); - if (null != flowSignalContext) { + FlowSignalSubContext f = block.getFlowSignalSubContext(offset); + if (null != f && null != f.signals && null != f.bases) { buf.append("FZ = "); - for (i=0;i required return buf.toString(); diff --git a/src/org/broad/igv/sam/AlignmentBlock.java b/src/org/broad/igv/sam/AlignmentBlock.java index fcf793715c..bd511f9868 100644 --- a/src/org/broad/igv/sam/AlignmentBlock.java +++ b/src/org/broad/igv/sam/AlignmentBlock.java @@ -37,8 +37,8 @@ public static AlignmentBlock getInstance(int start, byte[] bases, byte[] qualiti return new AlignmentBlock(start, bases, qualities, baseAlignment); } - public static AlignmentBlock getInstance(int start, byte[] bases, byte[] qualities, short[][][] flowSignals, Alignment baseAlignment) { - return new AlignmentBlockFS(start, bases, qualities, flowSignals, baseAlignment); + public static AlignmentBlock getInstance(int start, byte[] bases, byte[] qualities, FlowSignalContext fContext, Alignment baseAlignment) { + return new AlignmentBlockFS(start, bases, qualities, fContext, baseAlignment); } protected AlignmentBlock(int start, byte[] bases, byte[] qualities, Alignment baseAlignment) { @@ -103,7 +103,7 @@ public boolean hasFlowSignals() { } // Default implementation -- to be overriden - public short[][] getFlowSignalContext(int offset) { + public FlowSignalSubContext getFlowSignalSubContext(int offset) { return null; } diff --git a/src/org/broad/igv/sam/AlignmentBlockFS.java b/src/org/broad/igv/sam/AlignmentBlockFS.java index 761b292a77..86f050ba6e 100644 --- a/src/org/broad/igv/sam/AlignmentBlockFS.java +++ b/src/org/broad/igv/sam/AlignmentBlockFS.java @@ -1,28 +1,30 @@ package org.broad.igv.sam; +import java.util.Arrays; + /** - * Represents an alignment block which contains flow signals. Added to suppor IonTorrent alignments. + * Represents an alignment block which contains flow signals. Added to support IonTorrent alignments. * * @author Jim Robinson * @date 3/19/12 */ public class AlignmentBlockFS extends AlignmentBlock { - public short[][][] flowSignals = null; + public FlowSignalContext fContext = null; - protected AlignmentBlockFS(int start, byte[] bases, byte[] qualities, short[][][] flowSignals, Alignment baseAlignment) { + protected AlignmentBlockFS(int start, byte[] bases, byte[] qualities, FlowSignalContext fContext, Alignment baseAlignment) { super(start, bases, qualities, baseAlignment); - if (flowSignals != null && flowSignals.length == bases.length) { - this.flowSignals = flowSignals; + if (fContext != null && fContext.signals.length == bases.length) { + this.fContext = fContext; } } - public short[][] getFlowSignalContext(int offset) { - return this.flowSignals[offset]; + public FlowSignalSubContext getFlowSignalSubContext(int offset) { + return new FlowSignalSubContext(this.fContext.signals[offset], this.fContext.bases[offset]); } public boolean hasFlowSignals() { - return (null != this.flowSignals); + return (null != this.fContext); } } diff --git a/src/org/broad/igv/sam/AlignmentInterval.java b/src/org/broad/igv/sam/AlignmentInterval.java index bb716ff06c..dfa8edf08c 100644 --- a/src/org/broad/igv/sam/AlignmentInterval.java +++ b/src/org/broad/igv/sam/AlignmentInterval.java @@ -477,4 +477,4 @@ public void remove() { // ignore } } -} \ No newline at end of file +} diff --git a/src/org/broad/igv/sam/AlignmentRenderer.java b/src/org/broad/igv/sam/AlignmentRenderer.java index ee8e6bb222..add348bc47 100644 --- a/src/org/broad/igv/sam/AlignmentRenderer.java +++ b/src/org/broad/igv/sam/AlignmentRenderer.java @@ -18,6 +18,7 @@ import org.broad.igv.feature.genome.GenomeManager; import org.broad.igv.renderer.ContinuousColorScale; import org.broad.igv.renderer.GraphicUtils; +import org.broad.igv.sam.AlignmentTrack.ShadeBasesOption; import org.broad.igv.sam.AlignmentTrack.ColorOption; import org.broad.igv.sam.AlignmentTrack.RenderOptions; import org.broad.igv.sam.BisulfiteBaseInfo.DisplayStatus; @@ -635,7 +636,7 @@ private void drawAlignment( /** * Draw bases for an alignment block. The bases are "overlaid" on the block with a transparency value (alpha) - * that is proportional to the base quality score. + * that is proportional to the base quality score, or flow signal deviation, whichever is selected. * * @param context * @param rect @@ -650,7 +651,7 @@ private void drawBases(RenderContext context, RenderOptions renderOptions) { - boolean shadeBases = renderOptions.shadeBases; + ShadeBasesOption shadeBasesOption = renderOptions.shadeBasesOption; ColorOption colorOption = renderOptions.getColorOption(); // Disable showAllBases in bisulfite mode @@ -731,10 +732,64 @@ private void drawBases(RenderContext context, color = Color.black; } - if (shadeBases) { + if (ShadeBasesOption.QUALITY == shadeBasesOption) { byte qual = block.qualities[loc - start]; color = getShadedColor(qual, color, alignmentColor, prefs); - } + } else if (ShadeBasesOption.FLOW_SIGNAL_DEVIATION_READ == shadeBasesOption || ShadeBasesOption.FLOW_SIGNAL_DEVIATION_REFERENCE == shadeBasesOption) { + if (block.hasFlowSignals()) { + int flowSignal = (int)block.getFlowSignalSubContext(loc - start).signals[1][0]; + int expectedFlowSignal; + if (ShadeBasesOption.FLOW_SIGNAL_DEVIATION_READ == shadeBasesOption) { + expectedFlowSignal = 100 * (short)((flowSignal + 50.0) / 100.0); + } else { + // NB: this may estimate the reference homopolymer length incorrect in some cases, especially when we have + // an overcall/undercall situation. Proper estimation of the reads observed versus expected homopolymer + // length should use flow signal alignment (SamToFlowspace): https://github.com/iontorrent/Ion-Variant-Hunter + if (!misMatch) { + final byte readbase = read[idx]; + byte refbase = reference[idx]; + int pos; // zero based + expectedFlowSignal = 100; + + // Count HP length + pos = start + idx - 1; + while (0 <= pos && genome.getReference(chr, pos) == refbase) { + pos--; + expectedFlowSignal += 100; + } + pos = start + idx + 1; + while (pos < genome.getChromosome(chr).getLength() && genome.getReference(chr, pos) == refbase) { + pos++; + expectedFlowSignal += 100; + } + } else { + expectedFlowSignal = 0; + } + } + int flowSignalDiff = (expectedFlowSignal < flowSignal) ? (flowSignal - expectedFlowSignal) : (expectedFlowSignal - flowSignal); + // NB: this next section is some mangling to use the base quality color preferences... + if (flowSignalDiff <= 0) { + flowSignalDiff = 0; + } else if (50 < flowSignalDiff) { + flowSignalDiff = 50; + } + flowSignalDiff = 50 - flowSignalDiff; // higher is better + int minQ = prefs.getAsInt(PreferenceManager.SAM_BASE_QUALITY_MIN); + int maxQ = prefs.getAsInt(PreferenceManager.SAM_BASE_QUALITY_MAX); + flowSignalDiff = flowSignalDiff * (maxQ - minQ) / 50; + byte qual; + int pos = start + idx; + if (Byte.MAX_VALUE < flowSignalDiff) { + qual = Byte.MAX_VALUE; + } else if (flowSignalDiff < Byte.MIN_VALUE) { + qual = Byte.MIN_VALUE; + } else { + qual = (byte)flowSignalDiff; + } + // Finally, get the color + color = getShadedColor(qual, color, alignmentColor, prefs); + } + } double bisulfiteXaxisShift = (bisulfiteMode) ? bisinfo.getXaxisShift(idx) : 0; diff --git a/src/org/broad/igv/sam/AlignmentTrack.java b/src/org/broad/igv/sam/AlignmentTrack.java index 65a77e317c..f2ff198c30 100755 --- a/src/org/broad/igv/sam/AlignmentTrack.java +++ b/src/org/broad/igv/sam/AlignmentTrack.java @@ -65,6 +65,10 @@ public class AlignmentTrack extends AbstractTrack implements AlignmentTrackEvent static final int DOWNAMPLED_ROW_HEIGHT = 3; static final int DS_MARGIN_2 = 5; + public enum ShadeBasesOption { + NONE, QUALITY, FLOW_SIGNAL_DEVIATION_READ, FLOW_SIGNAL_DEVIATION_REFERENCE; + } + public enum ExperimentType { RNA, BISULFITE, OTHER; @@ -153,6 +157,7 @@ static BisulfiteContext strToValue(String str) { bisulfiteContextToContextString.put(BisulfiteContext.WCG, new Pair(new byte[]{'W'}, new byte[]{'G'})); } + static final ShadeBasesOption DEFAULT_SHADE_BASES_OPTION = ShadeBasesOption.QUALITY; static final ColorOption DEFAULT_COLOR_OPTION = ColorOption.INSERT_SIZE; static final boolean DEFAULT_SHOWALLBASES = false; static final BisulfiteContext DEFAULT_BISULFITE_CONTEXT = BisulfiteContext.CG; @@ -480,6 +485,46 @@ public void copyToClipboard(final TrackClickEvent e, Alignment alignment, double } } + + /** + * Copy the contents of the popup text to the system clipboard. + */ + public void copyFlowSignalDistribution(final TrackClickEvent e, double location) { + TreeMap map = new TreeMap(); + for (AlignmentInterval interval : dataManager.getLoadedIntervals()) { + Iterator alignmentIterator = interval.getAlignmentIterator(); + while (alignmentIterator.hasNext()) { + Alignment alignment = alignmentIterator.next(); + if (!alignment.contains(location)) { + continue; + } + AlignmentBlock[] blocks = alignment.getAlignmentBlocks(); + for (int i = 0; i < blocks.length; i++) { + AlignmentBlock block = blocks[i]; + if (!block.contains((int)location) || !block.hasFlowSignals()) { + continue; + } + short flowSignal = block.getFlowSignalSubContext((int)location- block.getStart()).signals[1][0]; + if (map.containsKey(flowSignal)) { + // increment + map.put(flowSignal, map.get(flowSignal) + 1); + } else { + // insert + map.put(flowSignal, 1); + } + } + } + } + StringBuffer buf = new StringBuffer(); + buf.append("{\n"); + for (Short key : map.keySet()) { + buf.append(" \"" + key + "\" : \"" + map.get(key) + "\"\n"); + } + buf.append("}\n"); + StringSelection stringSelection = new StringSelection(buf.toString()); + Clipboard clipboard = Toolkit.getDefaultToolkit().getSystemClipboard(); + clipboard.setContents(stringSelection, null); + } /** * Jump to the mate region @@ -813,7 +858,7 @@ public void setPairedArcView(boolean option) { public static class RenderOptions { - boolean shadeBases; + ShadeBasesOption shadeBasesOption; boolean shadeCenters; boolean flagUnmappedPairs; boolean showAllBases; @@ -840,7 +885,7 @@ public static class RenderOptions { RenderOptions() { PreferenceManager prefs = PreferenceManager.getInstance(); - shadeBases = prefs.getAsBoolean(PreferenceManager.SAM_SHADE_BASE_QUALITY); + shadeBasesOption = ShadeBasesOption.valueOf(prefs.get(PreferenceManager.SAM_SHADE_BASES)); shadeCenters = prefs.getAsBoolean(PreferenceManager.SAM_SHADE_CENTER); flagUnmappedPairs = prefs.getAsBoolean(PreferenceManager.SAM_FLAG_UNMAPPED_PAIR); computeIsizes = prefs.getAsBoolean(PreferenceManager.SAM_COMPUTE_ISIZES); @@ -884,11 +929,11 @@ public static class RenderOptions { public Map getPersistentState() { Map attributes = new HashMap(); PreferenceManager prefs = PreferenceManager.getInstance(); - if (shadeBases != prefs.getAsBoolean(PreferenceManager.SAM_SHADE_BASE_QUALITY)) { - attributes.put("shadeBases", String.valueOf(shadeBases)); + if (shadeBasesOption != DEFAULT_SHADE_BASES_OPTION) { + attributes.put("shadeBasesOption", shadeBasesOption.toString()); } if (shadeCenters != prefs.getAsBoolean(PreferenceManager.SAM_SHADE_CENTER)) { - attributes.put("shadeCenters", String.valueOf(shadeBases)); + attributes.put("shadeCenters", String.valueOf(shadeCenters)); } if (flagUnmappedPairs != prefs.getAsBoolean(PreferenceManager.SAM_FLAG_UNMAPPED_PAIR)) { attributes.put("flagUnmappedPairs", String.valueOf(flagUnmappedPairs)); @@ -937,9 +982,9 @@ public void restorePersistentState(Map attributes) { if (value != null) { setMinInsertSize(Integer.parseInt(value)); } - value = attributes.get("shadeBases"); + value = attributes.get("shadeBasesOption"); if (value != null) { - shadeBases = Boolean.parseBoolean(value); + shadeBasesOption = ShadeBasesOption.valueOf(value); } value = attributes.get("shadeCenters"); if (value != null) { @@ -1091,6 +1136,7 @@ class PopupMenu extends IGVPopupMenu { addSeparator(); add(TrackMenuUtils.getTrackRenameItem(tracks)); addCopyToClipboardItem(e); + addCopyFlowSignalDistributionToClipboardItem(e); addSeparator(); addGroupMenuItem(); @@ -1098,7 +1144,7 @@ class PopupMenu extends IGVPopupMenu { addColorByMenuItem(); addSeparator(); - addShadeBaseMenuItem(); + addShadeBaseByMenuItem(); addShowMismatchesMenuItem(); addShowAllBasesMenuItem(); @@ -1441,6 +1487,26 @@ public void actionPerformed(ActionEvent aEvt) { add(item); } + public void addCopyFlowSignalDistributionToClipboardItem(final TrackClickEvent te) { + final MouseEvent me = te.getMouseEvent(); + JMenuItem item = new JMenuItem("Copy the flow signal distrubtion for this base to the clipboard"); + final ReferenceFrame frame = te.getFrame(); + if (frame == null) { + item.setEnabled(false); + } else { + final double location = frame.getChromosomePosition(me.getX()); + + // Change track height by attribute + item.addActionListener(new ActionListener() { + + public void actionPerformed(ActionEvent aEvt) { + copyFlowSignalDistribution(te, location); + } + }); + } + add(item); + } + public void addViewPairedArcsMenuItem() { final JMenuItem item = new JCheckBoxMenuItem("View paired arcs"); item.setSelected(isPairedArcView()); @@ -1617,28 +1683,53 @@ public void actionPerformed(ActionEvent aEvt) { item.setEnabled(dataManager.isPairedEnd()); add(item); } + + private JCheckBoxMenuItem getShadeBasesMenuItem(String label, final ShadeBasesOption option) { + final JCheckBoxMenuItem mi = new JCheckBoxMenuItem(label); + mi.setSelected(renderOptions.shadeBasesOption == option); - - public void addShadeBaseMenuItem() { - // Change track height by attribute - final JMenuItem item = new JCheckBoxMenuItem("Shade base by quality"); - item.setSelected(renderOptions.shadeBases); - item.addActionListener(new ActionListener() { + if (option == ShadeBasesOption.NONE) { + mi.setSelected(renderOptions.shadeBasesOption == null || renderOptions.shadeBasesOption == option); + } + mi.addActionListener(new ActionListener() { public void actionPerformed(ActionEvent aEvt) { UIUtilities.invokeOnEventThread(new Runnable() { public void run() { - renderOptions.shadeBases = item.isSelected(); + if (mi.isSelected()) { + renderOptions.shadeBasesOption = option; + } else { + renderOptions.shadeBasesOption = ShadeBasesOption.NONE; + } refresh(); } }); + } }); - add(item); + return mi; } + + public void addShadeBaseByMenuItem() { + JMenu groupMenu = new JMenu("Shade bases by..."); + ButtonGroup group = new ButtonGroup(); + + Map mappings = new LinkedHashMap(); + mappings.put("none", ShadeBasesOption.NONE); + mappings.put("quality", ShadeBasesOption.QUALITY); + mappings.put("read flow signal deviation", ShadeBasesOption.FLOW_SIGNAL_DEVIATION_READ); + mappings.put("reference flow signal deviation", ShadeBasesOption.FLOW_SIGNAL_DEVIATION_REFERENCE); + for (Map.Entry el : mappings.entrySet()) { + JCheckBoxMenuItem mi = getShadeBasesMenuItem(el.getKey(), el.getValue()); + groupMenu.add(mi); + group.add(mi); + } + + add(groupMenu); + } public void addShowCoverageItem() { // Change track height by attribute diff --git a/src/org/broad/igv/sam/FlowSignalContext.java b/src/org/broad/igv/sam/FlowSignalContext.java new file mode 100644 index 0000000000..457afbf7ab --- /dev/null +++ b/src/org/broad/igv/sam/FlowSignalContext.java @@ -0,0 +1,19 @@ +package org.broad.igv.sam; + +import java.util.Arrays; + +/** + * Represents a flow signals context in an alignment block. Added to support IonTorrent alignments. + * + * @author Nils Homer + * @date 4/11/12 + */ +public class FlowSignalContext { + public short[][][]signals = null; + public char[][][]bases = null; + + public FlowSignalContext(short[][][] signals, char[][][] bases) { + this.signals = signals; + this.bases = bases; + } +} diff --git a/src/org/broad/igv/sam/FlowSignalContextBuilder.java b/src/org/broad/igv/sam/FlowSignalContextBuilder.java new file mode 100644 index 0000000000..b81679d90e --- /dev/null +++ b/src/org/broad/igv/sam/FlowSignalContextBuilder.java @@ -0,0 +1,221 @@ +package org.broad.igv.sam; + +import java.util.Arrays; + +/** + * Builds a flow signals context in an alignment block. Added to support IonTorrent alignments. + * + * @author Nils Homer + * @date 4/11/12 + */ +public class FlowSignalContextBuilder { + + private short[] flowSignals = null; + private String flowOrder = null; + private int flowSignalsIndex = -1; + private int flowOrderIndex = -1; + private int prevFlowSignalsStart = -1; + private int prevFlowSignalsEnd = -1; + private int flowOrderStart = -1; + private boolean readNegativeStrandFlag; + private boolean[] incorporations = null; // required for the reverse strand + + public FlowSignalContextBuilder(short[] flowSignals, String flowOrder, int flowOrderStart, byte[] readBases, int fromIdx, boolean readNegativeStrandFlag) { + if (null == flowSignals || null == flowOrder || flowOrderStart < 0) { + return; + } + + this.flowSignals = flowSignals; + this.flowOrder = flowOrder; + this.flowOrderIndex = this.flowOrderStart = flowOrderStart; + this.flowSignalsIndex = 0; // NB: the key sequence/barcode sequence should have been removed for the signals already + this.readNegativeStrandFlag = readNegativeStrandFlag; + + // init + if (this.readNegativeStrandFlag) { + int i; + this.incorporations = new boolean[this.flowSignals.length]; + // go to the end of the signals to find the first sequenced base + for (i=readBases.length-1;0<=i;i--) { + while (this.flowOrder.charAt(this.flowOrderIndex) != SamAlignment.NT2COMP[readBases[i]]) { + this.flowOrderIndex++; + this.flowSignalsIndex++; + this.incorporations[this.flowSignalsIndex] = false; + if (this.flowOrder.length() <= this.flowOrderIndex) { + this.flowOrderIndex = 0; + } + } + this.incorporations[this.flowSignalsIndex] = true; + } + this.prevFlowSignalsStart = this.flowSignalsIndex + 1; + this.prevFlowSignalsEnd = this.flowSignals.length - 1; + } else { + this.prevFlowSignalsStart = this.prevFlowSignalsEnd = 0; + while (this.flowOrder.charAt(this.flowOrderIndex) != readBases[0]) { + this.flowOrderIndex++; + this.flowSignalsIndex++; + if (this.flowOrder.length() <= this.flowOrderIndex) { + this.flowOrderIndex = 0; + } + } + this.prevFlowSignalsEnd = this.flowSignalsIndex - 1; + } + /* + if (0 < fromIdx) { // skip over leading bases (ex. soft clipped bases) + int i = 0; + while (0 <= this.flowSignalsIndex && this.flowSignalsIndex < this.flowSignals.length && i < fromIdx) { + short s = this.flowSignals[this.flowSignalsIndex]; + int nextFlowSignalsStart = -1, nextFlowSignalsEnd = -1; + int j = i + 1; + if (j < readBases.length) { + if (this.readNegativeStrandFlag) { + nextFlowSignalsEnd = this.flowSignalsIndex - 1; + // NB: loop condition is not symmetric to the forward, as we must respect the directionality of sequencing. + // For example, if our flow order is TACAG, and our read bases are TAG, then the flow signal vector is + // approximately 100,100,0,0,100. Since we move in the reverse direction with respect to the flow signal + // vector we must pre-compute where the flows incorporations are expected to occur, instead of just looking + // for the next flow that matches our next read base (we would place the A incorporation flow in the fourth flow, + // which is wrong). + while (!this.incorporations[this.flowSignalsIndex] || + this.flowOrder.charAt(this.flowOrderIndex) != SamAlignment.NT2COMP[readBases[j]]) { // NB: malicious input can cause infinite loops here + this.flowOrderIndex--; + this.flowSignalsIndex--; + if (this.flowOrderIndex < 0) { + this.flowOrderIndex = this.flowOrder.length() - 1; + } + } + nextFlowSignalsStart = this.flowSignalsIndex + 1; + } else { + nextFlowSignalsStart = this.flowSignalsIndex + 1; + while (this.flowOrder.charAt(this.flowOrderIndex) != readBases[j]) { // NB: malicious input can cause infinite loops here + this.flowOrderIndex++; + this.flowSignalsIndex++; + if (this.flowOrder.length() <= this.flowOrderIndex) { + this.flowOrderIndex = 0; + } + } + nextFlowSignalsEnd = this.flowSignalsIndex - 1; + } + } + // update for the next iteration + this.prevFlowSignalsStart = nextFlowSignalsStart; + this.prevFlowSignalsEnd = nextFlowSignalsEnd; + i++; // next base + } + } + */ + } + + // TODO: + // - support IUPAC bases + // - support lower/upper cases (is this necessary)? + public FlowSignalContext getFlowSignalContext(byte[] readBases, int fromIdx, int nBases) { + int i, idx; + short[][][] blockFlowSignals = null; + char[][][] blockFlowOrder = null; + + if (null == this.flowSignals) { + return null; + } + + blockFlowSignals = new short[nBases][][]; + blockFlowOrder = new char[nBases][][]; + //Default value + Arrays.fill(blockFlowSignals, null); + Arrays.fill(blockFlowOrder, null); + + // NB: should be at the first base of a HP + // Go through the bases + i = fromIdx; + idx = 0; + while (0 <= this.flowSignalsIndex && this.flowSignalsIndex < this.flowSignals.length && i < fromIdx + nBases) { + short s = this.flowSignals[this.flowSignalsIndex]; + char f = this.flowOrder.charAt((this.flowSignalsIndex + this.flowOrderStart) % this.flowOrder.length()); + int nextFlowSignalsStart = -1, nextFlowSignalsEnd = -1; + int j = i + 1; + if (j < readBases.length) { + if (this.readNegativeStrandFlag) { + nextFlowSignalsEnd = this.flowSignalsIndex - 1; + // NB: loop condition is not symmetric to the forward, as we must respect the directionality of sequencing. + // For example, if our flow order is TACAG, and our read bases are TAG, then the flow signal vector is + // approximately 100,100,0,0,100. Since we move in the reverse direction with respect to the flow signal + // vector we must pre-compute where the flows incorporations are expected to occur, instead of just looking + // for the next flow that matches our next read base (we would place the A incorporation flow in the fourth flow, + // which is wrong). + while (!this.incorporations[this.flowSignalsIndex] || + this.flowOrder.charAt(this.flowOrderIndex) != SamAlignment.NT2COMP[readBases[j]]) { // NB: malicious input can cause infinite loops here + this.flowOrderIndex--; + this.flowSignalsIndex--; + if (this.flowOrderIndex < 0) { + this.flowOrderIndex = this.flowOrder.length() - 1; + } + } + nextFlowSignalsStart = this.flowSignalsIndex + 1; + } else { + nextFlowSignalsStart = this.flowSignalsIndex + 1; + while (this.flowOrder.charAt(this.flowOrderIndex) != readBases[j]) { // NB: malicious input can cause infinite loops here + this.flowOrderIndex++; + this.flowSignalsIndex++; + if (this.flowOrder.length() <= this.flowOrderIndex) { + this.flowOrderIndex = 0; + } + } + nextFlowSignalsEnd = this.flowSignalsIndex - 1; + } + } + // set-up block + blockFlowSignals[idx] = new short[3][]; + blockFlowOrder[idx] = new char[3][]; + // this.previous context + if (0 <= this.prevFlowSignalsStart && this.prevFlowSignalsStart <= this.prevFlowSignalsEnd && this.prevFlowSignalsEnd < this.flowSignals.length) { + blockFlowSignals[idx][0] = new short[this.prevFlowSignalsEnd - this.prevFlowSignalsStart + 1]; + blockFlowOrder[idx][0] = new char[this.prevFlowSignalsEnd - this.prevFlowSignalsStart + 1]; + if (this.readNegativeStrandFlag) { + for (j = this.prevFlowSignalsEnd; this.prevFlowSignalsStart <= j; j--) { + blockFlowSignals[idx][0][this.prevFlowSignalsEnd - j] = this.flowSignals[j]; + blockFlowOrder[idx][0][this.prevFlowSignalsEnd - j] = this.flowOrder.charAt((j + this.flowOrderStart) % this.flowOrder.length()); + } + } else { + for (j = this.prevFlowSignalsStart; j <= this.prevFlowSignalsEnd; j++) { + blockFlowSignals[idx][0][j-this.prevFlowSignalsStart] = this.flowSignals[j]; + blockFlowOrder[idx][0][j-this.prevFlowSignalsStart] = this.flowOrder.charAt((j + this.flowOrderStart) % this.flowOrder.length()); + } + } + } else { + blockFlowSignals[idx][0] = null; + blockFlowOrder[idx][0] = null; + } + // current context + blockFlowSignals[idx][1] = new short[1]; + blockFlowOrder[idx][1] = new char[1]; + blockFlowSignals[idx][1][0] = s; + blockFlowOrder[idx][1][0] = f; + // next context + if (0 <= nextFlowSignalsStart && nextFlowSignalsStart <= nextFlowSignalsEnd && nextFlowSignalsEnd < this.flowSignals.length) { + blockFlowSignals[idx][2] = new short[nextFlowSignalsEnd - nextFlowSignalsStart + 1]; + blockFlowOrder[idx][2] = new char[nextFlowSignalsEnd - nextFlowSignalsStart + 1]; + if (this.readNegativeStrandFlag) { + for (j = nextFlowSignalsEnd; nextFlowSignalsStart <= j; j--) { + blockFlowSignals[idx][2][nextFlowSignalsEnd - j] = this.flowSignals[j]; + blockFlowOrder[idx][2][nextFlowSignalsEnd - j] = this.flowOrder.charAt((j + this.flowOrderStart) % this.flowOrder.length()); + } + } else { + for (j = nextFlowSignalsStart; j <= nextFlowSignalsEnd; j++) { + blockFlowSignals[idx][2][j-nextFlowSignalsStart] = this.flowSignals[j]; + blockFlowOrder[idx][2][j-nextFlowSignalsStart] = this.flowOrder.charAt((j + this.flowOrderStart) % this.flowOrder.length()); + } + } + } else { + blockFlowSignals[idx][2] = null; + blockFlowOrder[idx][2] = null; + } + // update for the next iteration + this.prevFlowSignalsStart = nextFlowSignalsStart; + this.prevFlowSignalsEnd = nextFlowSignalsEnd; + i++; // next base + idx++; // next base + } + + return new FlowSignalContext(blockFlowSignals, blockFlowOrder); + } +} diff --git a/src/org/broad/igv/sam/FlowSignalSubContext.java b/src/org/broad/igv/sam/FlowSignalSubContext.java new file mode 100644 index 0000000000..4df3461278 --- /dev/null +++ b/src/org/broad/igv/sam/FlowSignalSubContext.java @@ -0,0 +1,19 @@ +package org.broad.igv.sam; + +import java.util.Arrays; + +/** + * Represents a flow signals context in an alignment block focused on a given base. Added to support IonTorrent alignments. + * + * @author Nils Homer + * @date 4/11/12 + */ +public class FlowSignalSubContext { + short[][] signals = null; + char[][] bases = null; + + public FlowSignalSubContext(short[][] signals, char[][] bases) { + this.signals = signals; + this.bases = bases; + } +} diff --git a/src/org/broad/igv/sam/SamAlignment.java b/src/org/broad/igv/sam/SamAlignment.java index 155fa4244b..22bce35243 100644 --- a/src/org/broad/igv/sam/SamAlignment.java +++ b/src/org/broad/igv/sam/SamAlignment.java @@ -254,158 +254,6 @@ private void setPairStrands() { } } - // Helper class to get the flow signals context - private class FlowOrderAndSignalsBlockHelper { - - private short[] flowSignals = null; - private String flowOrder = null; - private int flowSignalsIndex = -1; - private int flowOrderIndex = -1; - private int prevFlowSignalsStart = -1; - private int prevFlowSignalsEnd = -1; - private boolean readNegativeStrandFlag; - private boolean[] incorporations = null; // required for the reverse strand - - public FlowOrderAndSignalsBlockHelper(short[] flowSignals, String flowOrder, int flowOrderStart, byte[] readBases, boolean readNegativeStrandFlag) { - if (null == flowSignals || null == flowOrder || flowOrderStart < 0) { - return; - } - - this.flowSignals = flowSignals; - this.flowOrder = flowOrder; - this.flowOrderIndex = flowOrderStart; - this.flowSignalsIndex = 0; // NB: the key sequence/barcode sequence should have been remove for the signals already - this.readNegativeStrandFlag = readNegativeStrandFlag; - - // init - if (this.readNegativeStrandFlag) { - int i; - this.incorporations = new boolean[this.flowSignals.length]; - // go to the end of the signals - for (i = readBases.length - 1; 0 <= i; i--) { - while (this.flowOrder.charAt(this.flowOrderIndex) != NT2COMP[readBases[i]]) { - this.flowOrderIndex++; - this.flowSignalsIndex++; - this.incorporations[this.flowSignalsIndex] = false; - if (this.flowOrder.length() <= this.flowOrderIndex) { - this.flowOrderIndex = 0; - } - } - this.incorporations[this.flowSignalsIndex] = true; - } - this.prevFlowSignalsStart = this.flowSignalsIndex + 1; - this.prevFlowSignalsEnd = this.flowSignals.length - 1; - } else { - this.prevFlowSignalsStart = this.prevFlowSignalsEnd = 0; - while (this.flowOrder.charAt(this.flowOrderIndex) != readBases[0]) { - this.flowOrderIndex++; - this.flowSignalsIndex++; - if (this.flowOrder.length() <= this.flowOrderIndex) { - this.flowOrderIndex = 0; - } - } - this.prevFlowSignalsEnd = this.flowSignalsIndex - 1; - } - } - - // TODO: - // - support IUPAC bases - // - support lower/upper cases (is this necessary)? - public short[][][] createBlockFlowSignals(byte[] readBases, int fromIdx, int nBases) { - int i, idx; - short[][][] blockFlowSignals = null; - - if (null == this.flowSignals) { - return null; - } - - blockFlowSignals = new short[nBases][][]; - //Default value - Arrays.fill(blockFlowSignals, null); - - // NB: should be at the first base of a HP - // Go through the bases - i = fromIdx; - idx = 0; - while (0 <= this.flowSignalsIndex && this.flowSignalsIndex < this.flowSignals.length && i < fromIdx + nBases) { - short s = this.flowSignals[this.flowSignalsIndex]; - int nextFlowSignalsStart = -1, nextFlowSignalsEnd = -1; - int j = i + 1; - if (j < readBases.length) { - if (this.readNegativeStrandFlag) { - nextFlowSignalsEnd = this.flowSignalsIndex - 1; - // NB: loop condition is not symmetric to the forward, as we must respect the directionality of sequencing. - // For example, if our flow order is TACAG, and our read bases are TAG, then the flow signal vector is - // approximately 100,100,0,0,100. Since we move in the reverse direction with respect to the flow signal - // vector we must pre-compute where the flows incorporations are expected to occur, instead of just looking - // for the next flow that matches our next read base (we would place the A incorporation flow in the fourth flow, - // which is wrong). - while (!this.incorporations[this.flowSignalsIndex] || - this.flowOrder.charAt(this.flowOrderIndex) != NT2COMP[readBases[j]]) { // NB: malicious input can cause infinite loops here - this.flowOrderIndex--; - this.flowSignalsIndex--; - if (this.flowOrderIndex < 0) { - this.flowOrderIndex = this.flowOrder.length() - 1; - } - } - nextFlowSignalsStart = this.flowSignalsIndex + 1; - } else { - nextFlowSignalsStart = this.flowSignalsIndex + 1; - while (this.flowOrder.charAt(this.flowOrderIndex) != readBases[j]) { // NB: malicious input can cause infinite loops here - this.flowOrderIndex++; - this.flowSignalsIndex++; - if (this.flowOrder.length() <= this.flowOrderIndex) { - this.flowOrderIndex = 0; - } - } - nextFlowSignalsEnd = this.flowSignalsIndex - 1; - } - } - // set-up block - blockFlowSignals[idx] = new short[3][]; - // this.previous context - if (0 <= this.prevFlowSignalsStart && this.prevFlowSignalsStart <= this.prevFlowSignalsEnd && this.prevFlowSignalsEnd < this.flowSignals.length) { - blockFlowSignals[idx][0] = new short[this.prevFlowSignalsEnd - this.prevFlowSignalsStart + 1]; - if (this.readNegativeStrandFlag) { - for (j = this.prevFlowSignalsEnd; this.prevFlowSignalsStart <= j; j--) { - blockFlowSignals[idx][0][this.prevFlowSignalsEnd - j] = this.flowSignals[j]; - } - } else { - for (j = this.prevFlowSignalsStart; j <= this.prevFlowSignalsEnd; j++) { - blockFlowSignals[idx][0][j - this.prevFlowSignalsStart] = this.flowSignals[j]; - } - } - } else { - blockFlowSignals[idx][0] = null; - } - // current context - blockFlowSignals[idx][1] = new short[1]; - blockFlowSignals[idx][1][0] = s; - // next context - if (0 <= nextFlowSignalsStart && nextFlowSignalsStart <= nextFlowSignalsEnd && nextFlowSignalsEnd < this.flowSignals.length) { - blockFlowSignals[idx][2] = new short[nextFlowSignalsEnd - nextFlowSignalsStart + 1]; - if (this.readNegativeStrandFlag) { - for (j = nextFlowSignalsEnd; nextFlowSignalsStart <= j; j--) { - blockFlowSignals[idx][2][nextFlowSignalsEnd - j] = this.flowSignals[j]; - } - } else { - for (j = nextFlowSignalsStart; j <= nextFlowSignalsEnd; j++) { - blockFlowSignals[idx][2][j - nextFlowSignalsStart] = this.flowSignals[j]; - } - } - } else { - blockFlowSignals[idx][2] = null; - } - // update for the next iteration - this.prevFlowSignalsStart = nextFlowSignalsStart; - this.prevFlowSignalsEnd = nextFlowSignalsEnd; - i++; // next base - idx++; // next base - } - return blockFlowSignals; - } - } - /** * Create the alignment blocks from the read bases and alignment information in the CIGAR * string. The CIGAR string encodes insertions, deletions, skipped regions, and padding. @@ -503,10 +351,10 @@ private void createAlignmentBlocks(String cigarString, byte[] readBases, byte[] int blockIdx = 0; int insertionIdx = 0; int gapIdx = 0; - FlowOrderAndSignalsBlockHelper fBlockHelper = null; + FlowSignalContextBuilder fBlockBuilder = null; if (null != flowSignals) { if (0 < readBases.length) { - fBlockHelper = new FlowOrderAndSignalsBlockHelper(flowSignals, flowOrder, flowOrderStart, readBases, this.readNegativeStrandFlag); + fBlockBuilder = new FlowSignalContextBuilder(flowSignals, flowOrder, flowOrderStart, readBases, fromIdx, this.readNegativeStrandFlag); } } prevOp = 0; @@ -520,7 +368,6 @@ private void createAlignmentBlocks(String cigarString, byte[] readBases, byte[] byte[] blockBases = new byte[op.nBases]; byte[] blockQualities = new byte[op.nBases]; - short[][][] blockFlowSignals = null; AlignmentBlock block = null; //Default value @@ -543,10 +390,10 @@ private void createAlignmentBlocks(String cigarString, byte[] readBases, byte[] } else { System.arraycopy(readBaseQualities, fromIdx, blockQualities, 0, op.nBases); } - - if (null != fBlockHelper) { - blockFlowSignals = fBlockHelper.createBlockFlowSignals(readBases, fromIdx, op.nBases); - block = AlignmentBlock.getInstance(blockStart, blockBases, blockQualities, blockFlowSignals, this); + + if (null != fBlockBuilder) { + block = AlignmentBlock.getInstance(blockStart, blockBases, blockQualities, + fBlockBuilder.getFlowSignalContext(readBases, fromIdx, op.nBases), this); } else { block = AlignmentBlock.getInstance(blockStart, blockBases, blockQualities, this); } @@ -568,7 +415,6 @@ private void createAlignmentBlocks(String cigarString, byte[] readBases, byte[] gapTypes[gapIdx++] = op.operator; } else if (op.operator == INSERTION) { AlignmentBlock block = null; - short[][][] blockFlowSignals = null; // This gap is between blocks split by insertion. It is a zero // length gap but must be accounted for. @@ -588,10 +434,10 @@ private void createAlignmentBlocks(String cigarString, byte[] readBases, byte[] } else { System.arraycopy(readBaseQualities, fromIdx, blockQualities, 0, op.nBases); } - - if (null != fBlockHelper) { - blockFlowSignals = fBlockHelper.createBlockFlowSignals(readBases, fromIdx, op.nBases); - block = AlignmentBlock.getInstance(blockStart, blockBases, blockQualities, blockFlowSignals, this); + + if (null != fBlockBuilder) { + block = AlignmentBlock.getInstance(blockStart, blockBases, blockQualities, + fBlockBuilder.getFlowSignalContext(readBases, fromIdx, op.nBases), this); } else { block = AlignmentBlock.getInstance(blockStart, blockBases, blockQualities, this); } diff --git a/src/org/broad/igv/session/IGVSessionReader.java b/src/org/broad/igv/session/IGVSessionReader.java index 4b7c3254ec..92f18ba579 100644 --- a/src/org/broad/igv/session/IGVSessionReader.java +++ b/src/org/broad/igv/session/IGVSessionReader.java @@ -211,7 +211,7 @@ public static enum SessionAttribute { //TODO Add the following into the Attributes /* - boolean shadeBases; + ShadeBasesOption shadeBases; boolean shadeCenters; boolean flagUnmappedPairs; boolean showAllBases; diff --git a/src/org/broad/igv/track/TrackLoader.java b/src/org/broad/igv/track/TrackLoader.java index e742276c29..c47e39e61e 100644 --- a/src/org/broad/igv/track/TrackLoader.java +++ b/src/org/broad/igv/track/TrackLoader.java @@ -277,6 +277,7 @@ public List load(ResourceLocator locator, Genome genome) { return newTracks; } catch (IOException e) { log.error(e); + e.printStackTrace(); throw new DataLoadException(e.getMessage(), path); } diff --git a/src/org/broad/igv/ui/PreferencesEditor.java b/src/org/broad/igv/ui/PreferencesEditor.java index baf3031cba..3ee2188e21 100644 --- a/src/org/broad/igv/ui/PreferencesEditor.java +++ b/src/org/broad/igv/ui/PreferencesEditor.java @@ -26,6 +26,8 @@ import org.broad.igv.PreferenceManager; import org.broad.igv.data.expression.ProbeToLocusMap; import org.broad.igv.batch.CommandListener; +import org.broad.igv.sam.AlignmentTrack; +import org.broad.igv.sam.AlignmentTrack.ShadeBasesOption; import org.broad.igv.feature.genome.GenomeManager; import org.broad.igv.sam.CachingQueryReader; import org.broad.igv.ui.color.ColorUtilities; @@ -2581,12 +2583,22 @@ private void samDownsampleCountFieldActionPerformed(java.awt.event.ActionEvent e } private void samShadeMismatchedBaseCBActionPerformed(java.awt.event.ActionEvent evt) { - updatedPreferenceMap.put( - PreferenceManager.SAM_SHADE_BASE_QUALITY, - String.valueOf(samShadeMismatchedBaseCB.isSelected())); - samMinBaseQualityField.setEnabled(samShadeMismatchedBaseCB.isSelected()); - samMaxBaseQualityField.setEnabled(samShadeMismatchedBaseCB.isSelected()); - + if (samShadeMismatchedBaseCB.isSelected()) { + updatedPreferenceMap.put( + PreferenceManager.SAM_SHADE_BASES, + ShadeBasesOption.QUALITY.toString()); + samMinBaseQualityField.setEnabled(samShadeMismatchedBaseCB.isSelected()); + samMaxBaseQualityField.setEnabled(samShadeMismatchedBaseCB.isSelected()); + } else { + PreferenceManager prefMgr = PreferenceManager.getInstance(); + if (ShadeBasesOption.QUALITY == ShadeBasesOption.valueOf(prefMgr.get(PreferenceManager.SAM_SHADE_BASES))) { + updatedPreferenceMap.put( + PreferenceManager.SAM_SHADE_BASES, + ShadeBasesOption.NONE.toString()); + samMinBaseQualityField.setEnabled(false); + samMaxBaseQualityField.setEnabled(false); + } + } } private void showCenterLineCBActionPerformed(ActionEvent e) { @@ -3323,7 +3335,7 @@ private void initValues() { showSoftClippedCB.setSelected(prefMgr.getAsBoolean(PreferenceManager.SAM_SHOW_SOFT_CLIPPED)); samFlagUnmappedPairCB.setSelected(prefMgr.getAsBoolean(PreferenceManager.SAM_FLAG_UNMAPPED_PAIR)); showCenterLineCB.setSelected(prefMgr.getAsBoolean(PreferenceManager.SAM_SHOW_CENTER_LINE)); - samShadeMismatchedBaseCB.setSelected(prefMgr.getAsBoolean(PreferenceManager.SAM_SHADE_BASE_QUALITY)); + samShadeMismatchedBaseCB.setSelected(ShadeBasesOption.QUALITY == ShadeBasesOption.valueOf(prefMgr.get(PreferenceManager.SAM_SHADE_BASES))); samMinBaseQualityField.setText((String.valueOf(prefMgr.getAsInt(PreferenceManager.SAM_BASE_QUALITY_MIN)))); samMaxBaseQualityField.setText((String.valueOf(prefMgr.getAsInt(PreferenceManager.SAM_BASE_QUALITY_MAX)))); samMinBaseQualityField.setEnabled(samShadeMismatchedBaseCB.isSelected());