Skip to content
Snippets Groups Projects
BamView.java 136 KiB
Newer Older
tjc's avatar
tjc committed
/* BamView
tjc's avatar
tjc committed
 *
 * created: 2009
 *
 * This file is part of Artemis
 *
 * Copyright(C) 2009  Genome Research Limited
 *
 * This program is free software; you can redistribute it and/or
 * modify it under the terms of the GNU General Public License
 * as published by the Free Software Foundation; either version 2
 * of the License, or(at your option) any later version.
 *
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this program; if not, write to the Free Software
 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
 *
 */
tjc's avatar
tjc committed
package uk.ac.sanger.artemis.components.alignment;

import java.awt.AlphaComposite;
tjc's avatar
tjc committed
import java.awt.BasicStroke;
import java.awt.BorderLayout;
import java.awt.Color;
import java.awt.Composite;
tjc's avatar
tjc committed
import java.awt.Dimension;
tjc's avatar
tjc committed
import java.awt.FlowLayout;
tjc's avatar
tjc committed
import java.awt.FontMetrics;
tjc's avatar
tjc committed
import java.awt.Graphics;
import java.awt.Graphics2D;
import java.awt.GridBagConstraints;
import java.awt.GridBagLayout;
import java.awt.Insets;
tjc's avatar
tjc committed
import java.awt.Point;
import java.awt.Rectangle;
tjc's avatar
tjc committed
import java.awt.Stroke;
tjc's avatar
tjc committed
import java.awt.event.ActionEvent;
import java.awt.event.ActionListener;
tjc's avatar
tjc committed
import java.awt.event.AdjustmentEvent;
import java.awt.event.AdjustmentListener;
tjc's avatar
tjc committed
import java.awt.event.ItemEvent;
import java.awt.event.ItemListener;
import java.awt.event.KeyAdapter;
import java.awt.event.KeyEvent;
tjc's avatar
tjc committed
import java.awt.event.MouseAdapter;
import java.awt.event.MouseEvent;
import java.awt.event.MouseMotionListener;
import java.awt.image.BufferedImage;
tjc's avatar
tjc committed
import java.io.BufferedReader;
tjc's avatar
tjc committed
import java.io.File;
tjc's avatar
tjc committed
import java.io.FileOutputStream;
import java.io.IOException;
tjc's avatar
tjc committed
import java.io.InputStream;
import java.io.InputStreamReader;
import java.lang.management.ManagementFactory;
import java.lang.management.MemoryMXBean;
import java.lang.reflect.Field;
tjc's avatar
tjc committed
import java.net.URL;
tjc's avatar
tjc committed
import java.util.Collections;
tjc's avatar
tjc committed
import java.util.Enumeration;
import java.util.HashMap;
tjc's avatar
tjc committed
import java.util.Hashtable;
import java.util.List;
import java.util.Map;
tjc's avatar
tjc committed
import java.util.Vector;
import java.util.concurrent.CountDownLatch;
import java.util.concurrent.ExecutorService;
import java.util.concurrent.Executors;
tjc's avatar
tjc committed

import javax.swing.BorderFactory;
import javax.swing.Box;
tjc's avatar
tjc committed
import javax.swing.BoxLayout;
import javax.swing.ButtonGroup;
import javax.swing.ImageIcon;
tjc's avatar
tjc committed
import javax.swing.JButton;
tjc's avatar
tjc committed
import javax.swing.JCheckBox;
tjc's avatar
tjc committed
import javax.swing.JCheckBoxMenuItem;
import javax.swing.JComponent;
tjc's avatar
tjc committed
import javax.swing.JFrame;
import javax.swing.JMenu;
tjc's avatar
tjc committed
import javax.swing.JMenuBar;
tjc's avatar
tjc committed
import javax.swing.JOptionPane;
tjc's avatar
tjc committed
import javax.swing.JPanel;
import javax.swing.JPopupMenu;
import javax.swing.JRadioButton;
tjc's avatar
tjc committed
import javax.swing.JScrollBar;
tjc's avatar
tjc committed
import javax.swing.JScrollPane;
import javax.swing.JSeparator;
tcarver's avatar
tcarver committed
import javax.swing.JSlider;
tjc's avatar
tjc committed
import javax.swing.JTextField;
import javax.swing.SwingUtilities;
tjc's avatar
tjc committed
import javax.swing.UIManager;
import javax.swing.border.Border;
import javax.swing.border.EmptyBorder;
tcarver's avatar
tcarver committed
import javax.swing.event.ChangeEvent;
import javax.swing.event.ChangeListener;
tjc's avatar
tjc committed

import net.sf.picard.reference.ReferenceSequenceFile;
import net.sf.picard.sam.BuildBamIndex;
import net.sf.samtools.AlignmentBlock;
import net.sf.samtools.SAMException;
tjc's avatar
tjc committed
import net.sf.samtools.SAMFileHeader;
import net.sf.samtools.SAMFileReader;
tcarver's avatar
tcarver committed
import net.sf.samtools.SAMReadGroupRecord;
tjc's avatar
tjc committed
import net.sf.samtools.SAMRecord;
import net.sf.samtools.SAMSequenceRecord;
import net.sf.samtools.SAMFileReader.ValidationStringency;
import net.sf.samtools.util.CloseableIterator;

tjc's avatar
tjc committed
import uk.ac.sanger.artemis.Entry;
import uk.ac.sanger.artemis.EntryGroup;
import uk.ac.sanger.artemis.FeatureVector;
tjc's avatar
tjc committed
import uk.ac.sanger.artemis.Options;
import uk.ac.sanger.artemis.Selection;
import uk.ac.sanger.artemis.SelectionChangeEvent;
import uk.ac.sanger.artemis.SelectionChangeListener;
tjc's avatar
tjc committed
import uk.ac.sanger.artemis.SimpleEntryGroup;
import uk.ac.sanger.artemis.circular.TextFieldInt;
import uk.ac.sanger.artemis.components.DisplayAdjustmentEvent;
import uk.ac.sanger.artemis.components.DisplayAdjustmentListener;
import uk.ac.sanger.artemis.components.EntryEdit;
tjc's avatar
tjc committed
import uk.ac.sanger.artemis.components.EntryFileDialog;
import uk.ac.sanger.artemis.components.FeatureDisplay;
import uk.ac.sanger.artemis.components.FileViewer;
import uk.ac.sanger.artemis.components.IndexReferenceEvent;
tjc's avatar
tjc committed
import uk.ac.sanger.artemis.components.MessageDialog;
import uk.ac.sanger.artemis.components.NonModalDialog;
import uk.ac.sanger.artemis.components.SequenceComboBox;
import uk.ac.sanger.artemis.components.SwingWorker;
import uk.ac.sanger.artemis.editor.MultiLineToolTipUI;
tjc's avatar
tjc committed
import uk.ac.sanger.artemis.io.EntryInformation;
tjc's avatar
tjc committed
import uk.ac.sanger.artemis.io.Range;
tjc's avatar
tjc committed
import uk.ac.sanger.artemis.sequence.Bases;
import uk.ac.sanger.artemis.sequence.MarkerRange;
tjc's avatar
tjc committed
import uk.ac.sanger.artemis.sequence.NoSequenceException;
import uk.ac.sanger.artemis.util.Document;
import uk.ac.sanger.artemis.util.DocumentFactory;
import uk.ac.sanger.artemis.util.FTPSeekableStream;
tjc's avatar
tjc committed
import uk.ac.sanger.artemis.util.OutOfRangeException;

public class BamView extends JPanel
tjc's avatar
tjc committed
                     implements DisplayAdjustmentListener, SelectionChangeListener
tjc's avatar
tjc committed
{
  private static final long serialVersionUID = 1L;
tjc's avatar
tjc committed

  private List<BamViewRecord> readsInView;
tjc's avatar
tjc committed
  private Hashtable<String, SAMFileReader> samFileReaderHash = new Hashtable<String, SAMFileReader>();
tcarver's avatar
tcarver committed
  private List<SAMReadGroupRecord> readGroups = new Vector<SAMReadGroupRecord>();
tjc's avatar
tjc committed

  private HashMap<String, Integer> seqLengths = new HashMap<String, Integer>();
  private HashMap<String, Integer> offsetLengths;
tjc's avatar
tjc committed
  private Vector<String> seqNames = new Vector<String>();
tcarver's avatar
tcarver committed
  protected List<Short> hideBamList = new Vector<Short>();
tjc's avatar
tjc committed

tjc's avatar
tjc committed
  private SAMRecordPredicate samRecordFlagPredicate;
tjc's avatar
tjc committed
  private SAMRecordMapQPredicate samRecordMapQPredicate;
tjc's avatar
tjc committed

  private SAMRecordFilter filterFrame;

tjc's avatar
tjc committed
  private Bases bases;
tjc's avatar
tjc committed
  private JScrollPane jspView;
tjc's avatar
tjc committed
  private JScrollBar scrollBar;
  
  private SequenceComboBox combo;
  private boolean isOrientation = false;
tjc's avatar
tjc committed
  private boolean isSingle = false;
  private boolean isSNPs = false;
  private boolean isCoverage = false;
tjc's avatar
tjc committed
  private boolean isSNPplot = false;
tjc's avatar
tjc committed
  
  private EntryEdit entry_edit;
  private FeatureDisplay feature_display;
  private Selection selection;
tjc's avatar
tjc committed
  private JPanel mainPanel = new JPanel();
  private CoveragePanel coveragePanel;
tjc's avatar
tjc committed
  private SnpPanel snpPanel;
  private boolean logScale = false;
tjc's avatar
tjc committed
  private int nbasesInView;
  
  private int startBase = -1;
  private int endBase   = -1;
tjc's avatar
tjc committed
  private int laststart;
  private int lastend;
  private boolean asynchronous = true;
  private boolean showBaseAlignment = false;
  private JMenu bamFilesMenu = new JMenu("BAM files");
  private JCheckBoxMenuItem logMenuItem = new JCheckBoxMenuItem("Use Log Scale", logScale);
  
  private JCheckBoxMenuItem cbStackView = new JCheckBoxMenuItem("Stack", true);
  private JCheckBoxMenuItem cbPairedStackView = new JCheckBoxMenuItem("Paired Stack");
  private JCheckBoxMenuItem cbStrandStackView = new JCheckBoxMenuItem("Strand Stack");
  private JCheckBoxMenuItem cbIsizeStackView = new JCheckBoxMenuItem("Inferred Size", false);
  private JCheckBoxMenuItem cbCoverageView = new JCheckBoxMenuItem("Coverage", false);
  private JCheckBoxMenuItem cbCoverageStrandView = new JCheckBoxMenuItem("Coverage by Strand", false);
tcarver's avatar
tcarver committed
  private JCheckBoxMenuItem cbCoverageHeatMap = new JCheckBoxMenuItem("Coverage Heat Map", false);
  private JCheckBoxMenuItem cbLastSelected;
  
  private ButtonGroup buttonGroup = new ButtonGroup();
  
tcarver's avatar
tcarver committed
  private JCheckBoxMenuItem colourByReadGrp = new JCheckBoxMenuItem("Read Group");
  private JCheckBoxMenuItem colourByStrandTag = new JCheckBoxMenuItem("RNASeq Strand Specific Tag (XS)");
  private JCheckBoxMenuItem colourByCoverageColour = new JCheckBoxMenuItem("Coverage Plot Colours");
  private JCheckBoxMenuItem baseQualityColour = new JCheckBoxMenuItem("Base Quality");
tjc's avatar
tjc committed
  private JCheckBoxMenuItem markInsertions = new JCheckBoxMenuItem("Mark Insertions", true);
    AlphaComposite.getInstance(AlphaComposite.SRC_OVER, 0.6f);
  
tcarver's avatar
tcarver committed
  private ReadGroupsFrame readGrpFrame;
tcarver's avatar
tcarver committed
  private GroupBamFrame groupsFrame = new GroupBamFrame(this, bamFilesMenu);
  private CoveragePanel coverageView = new CoveragePanel();
  
tcarver's avatar
tcarver committed
  protected static String BAM_SUFFIX = ".*\\.(bam|cram)$";
  /** Used to colour the frames. */
  private static Color LIGHT_GREY = new Color(200, 200, 200);
  private static Color DARK_GREEN = new Color(0, 150, 0);
  private static Color DARK_ORANGE = new Color(255,140,0);
  private static Color DEEP_PINK   = new Color(139,10,80);
  private static Color NON_SELECTED_READ_HIGHLIGHT_COLOUR = new Color(189,103,107);
  
  private Point lastMousePoint = null;
  private volatile BamViewRecord mouseOverSAMRecord = null;
  private volatile BamViewRecord highlightSAMRecord = null;
  private String mouseOverInsertion;
  // record of where a mouse drag starts
  protected int dragStart = -1;
  private static int MAX_BASES = 26000;
  private int maxHeight = 800;
tjc's avatar
tjc committed
  private boolean concatSequences = false;
tjc's avatar
tjc committed
  private int ALIGNMENT_PIX_PER_BASE;
  private JPopupMenu popup;
  private PopupMessageFrame popFrame = new PopupMessageFrame();
  private PopupMessageFrame waitingFrame = new PopupMessageFrame("waiting...");
  private ExecutorService bamReadTaskExecutor;
  private int MAX_COVERAGE = Integer.MAX_VALUE;
tcarver's avatar
tcarver committed
  private float readLnHgt = 2.0f;
  
tjc's avatar
tjc committed
  public static org.apache.log4j.Logger logger4j = 
    org.apache.log4j.Logger.getLogger(BamView.class);
  
  public BamView(List<String> bamList, 
                String reference,
                int nbasesInView,
                final EntryEdit entry_edit,
                final FeatureDisplay feature_display,
                final Bases bases,
                final JPanel containerPanel,
                final JFrame frame)
  {
    this(bamList, reference, nbasesInView, feature_display, bases, containerPanel, frame);
    this.entry_edit = entry_edit;
  }
  
  public BamView(List<String> bamList, 
tjc's avatar
tjc committed
                 String reference,
tjc's avatar
tjc committed
                 int nbasesInView,
                 final FeatureDisplay feature_display,
                 final Bases bases,
tjc's avatar
tjc committed
                 final JPanel containerPanel,
tjc's avatar
tjc committed
                 final JFrame frame)
tjc's avatar
tjc committed
  {
    super();
    setBackground(Color.white);
    this.bamList = bamList;
tjc's avatar
tjc committed
    this.nbasesInView = nbasesInView;
tjc's avatar
tjc committed
    this.feature_display = feature_display;
    this.bases = bases;
tjc's avatar
tjc committed

tjc's avatar
tjc committed
    containerPanel.setLayout(new BoxLayout(containerPanel, BoxLayout.Y_AXIS));
    containerPanel.add(mainPanel);
    // filter out unmapped reads by default
    setSamRecordFlagPredicate(
        new SAMRecordFlagPredicate(SAMRecordFlagPredicate.READ_UNMAPPED_FLAG));
    
tjc's avatar
tjc committed
    if(reference != null)
    {
      System.setProperty("reference", reference); // for CRAM
tjc's avatar
tjc committed
      EntryGroup entryGroup = new SimpleEntryGroup();
tjc's avatar
tjc committed
      try
      {
        getEntry(reference,entryGroup);
      }
      catch (NoSequenceException e)
      {
        e.printStackTrace();
      }
    }
tjc's avatar
tjc committed
    
    if(Options.getOptions().getIntegerProperty("bam_read_thread") != null)
    { 
      logger4j.debug("BAM READ THREADS="+Options.getOptions().getIntegerProperty("bam_read_thread"));
      bamReadTaskExecutor = Executors.newFixedThreadPool(
          Options.getOptions().getIntegerProperty("bam_read_thread"));
    }
    else
      bamReadTaskExecutor = Executors.newFixedThreadPool(1);
    
    
    if(Options.getOptions().getIntegerProperty("bam_max_coverage") != null)
    { 
      logger4j.debug("BAM MAX COVERAGE="+Options.getOptions().getIntegerProperty("bam_max_coverage"));
      MAX_COVERAGE = Options.getOptions().getIntegerProperty("bam_max_coverage");
    }  

    try
    {
      readHeaderPicard();
    }
    catch(java.lang.UnsupportedClassVersionError err)
    {
      JOptionPane.showMessageDialog(null, 
          "This requires Java 1.6 or higher.", 
          "Check Java Version", JOptionPane.WARNING_MESSAGE);
    }
tjc's avatar
tjc committed
    catch (IOException e)
    {
      e.printStackTrace();
    }
tjc's avatar
tjc committed

    final javax.swing.plaf.FontUIResource font_ui_resource =
      Options.getOptions().getFontUIResource();
    
tjc's avatar
tjc committed
    Enumeration<Object> keys = UIManager.getDefaults().keys();
tjc's avatar
tjc committed
    while(keys.hasMoreElements()) 
    {
      Object key = keys.nextElement();
      Object value = UIManager.get(key);
      if(value instanceof javax.swing.plaf.FontUIResource) 
        UIManager.put(key, font_ui_resource);
    }
tjc's avatar
tjc committed

    setFont(Options.getOptions().getFont());
tjc's avatar
tjc committed
    FontMetrics fm  = getFontMetrics(getFont());
    ALIGNMENT_PIX_PER_BASE = fm.charWidth('M');
    BASE_HEIGHT = fm.getMaxAscent();
    selection = new Selection(null);
    
    MultiLineToolTipUI.initialize();
    setToolTipText("");
    
    buttonGroup.add(cbStackView);
    buttonGroup.add(cbPairedStackView);
    buttonGroup.add(cbStrandStackView);
    buttonGroup.add(cbIsizeStackView);
    buttonGroup.add(cbCoverageView);
    buttonGroup.add(cbCoverageStrandView);
tcarver's avatar
tcarver committed
    buttonGroup.add(cbCoverageHeatMap);
tjc's avatar
tjc committed
    addMouseListener(new PopupListener());
tjc's avatar
tjc committed
    jspView = new JScrollPane(this, 
        JScrollPane.VERTICAL_SCROLLBAR_AS_NEEDED,
        JScrollPane.HORIZONTAL_SCROLLBAR_NEVER);
    
tjc's avatar
tjc committed
    jspView.setViewportBorder(BorderFactory.createMatteBorder(0, 0, 1, 0, Color.DARK_GRAY));
tjc's avatar
tjc committed
    Border empty = new EmptyBorder(0,0,0,0);
    jspView.setBorder(empty);
    jspView.getVerticalScrollBar().setUnitIncrement(8);
tjc's avatar
tjc committed
    addBamToPanel(frame);

    // apply command line options
    if(System.getProperty("show_snps") != null)
      isSNPs = true;
    if(System.getProperty("show_snp_plot") != null)
    {
      isSNPplot = true;
      snpPanel.setVisible(isSNPplot);
    }
    if(System.getProperty("show_cov_plot") != null)
    {
      isCoverage = true;
      coveragePanel.setVisible(isCoverage);
    }
  }
  
  public String getToolTipText()
  {	
		if(isCoverageView(getPixPerBaseByWidth()) && lastMousePoint != null)
	      return coverageView.getToolTipText(
	          lastMousePoint.y-getJspView().getViewport().getViewPosition().y);
	    
	    if(mouseOverSAMRecord == null)
	      return null;
	    
	    String msg = 
	        mouseOverSAMRecord.sam.getReadName() + "\n" + 
	        mouseOverSAMRecord.sam.getAlignmentStart() + ".." +
	        mouseOverSAMRecord.sam.getAlignmentEnd() + 
	       (mouseOverSAMRecord.sam.getReadGroup() != null ? "\nRG="+mouseOverSAMRecord.sam.getReadGroup().getId() : "") +
	        "\nisize=" + mouseOverSAMRecord.sam.getInferredInsertSize() + 
	        "\nmapq=" + mouseOverSAMRecord.sam.getMappingQuality()+
	        "\nrname="+ mouseOverSAMRecord.sam.getReferenceName();

	    if( mouseOverSAMRecord.sam.getReadPairedFlag() && 
	        mouseOverSAMRecord.sam.getProperPairFlag() && 
	       !mouseOverSAMRecord.sam.getMateUnmappedFlag())
	    {
	      msg = msg +
	        "\nstrand (read/mate): "+
	       (mouseOverSAMRecord.sam.getReadNegativeStrandFlag() ? "-" : "+")+" / "+
	       (mouseOverSAMRecord.sam.getMateNegativeStrandFlag() ? "-" : "+");
	    }
	    else
	      msg = msg +
	        "\nstrand (read/mate): "+
	       (mouseOverSAMRecord.sam.getReadNegativeStrandFlag() ? "-" : "+");
	    
	    if(msg != null && mouseOverInsertion != null)
	      msg = msg + "\nInsertion at:" +mouseOverInsertion;
	    
	    return msg;
tjc's avatar
tjc committed
  }
tjc's avatar
tjc committed
  
  /**
   * Get the BAM index file from the list
   * @param bam
   * @return
   * @throws IOException
   */
  private File getBamIndexFile(String bam) throws IOException
  {
    File bamIndexFile = null;
    if (bam.startsWith("http") || bam.startsWith("ftp"))
tjc's avatar
tjc committed
    {
      final URL urlBamIndexFile = new URL(bam + ".bai");
      InputStream is = urlBamIndexFile.openStream();
tjc's avatar
tjc committed

tjc's avatar
tjc committed
      // Create temp file.
      bamIndexFile = File.createTempFile(urlBamIndexFile.getFile().replaceAll(
          "[\\/\\s]", "_"), ".bai");
      bamIndexFile.deleteOnExit();

      FileOutputStream out = new FileOutputStream(bamIndexFile);
      int c;
      while ((c = is.read()) != -1)
        out.write(c);
      out.flush();
      out.close();
      is.close();

tjc's avatar
tjc committed
      logger4j.debug("create... " + bamIndexFile.getAbsolutePath());
tjc's avatar
tjc committed
    }
    else
tjc's avatar
tjc committed
      bamIndexFile = new File(bam + ".bai");
      if(!bamIndexFile.exists())
tcarver's avatar
tcarver committed
      {
        final File cramIndexFile = new File(bam + ".crai");
        if(cramIndexFile.exists())
tcarver's avatar
tcarver committed
        {
          logger4j.debug(
                 "ERROR: CRAM INDEX FILE ("+cramIndexFile.getName()+ 
                 ") EXPECTING A BAM INDEX FILE (USE THIS OPTION --bam-style-index) ");
tcarver's avatar
tcarver committed
          return cramIndexFile;
tcarver's avatar
tcarver committed
        }
tcarver's avatar
tcarver committed
      }
tjc's avatar
tjc committed

    return bamIndexFile;
  }
tjc's avatar
tjc committed
    
tjc's avatar
tjc committed
  /**
   * Get the SAM file reader.
   * @param bam
   * @return
   * @throws IOException
   */
  private SAMFileReader getSAMFileReader(final String bam) throws IOException
  {
    // parsing of the header happens during SAMFileReader construction, 
    // so need to set the default stringency
    SAMFileReader.setDefaultValidationStringency(ValidationStringency.LENIENT);
    
tjc's avatar
tjc committed
    if(samFileReaderHash.containsKey(bam))
      return samFileReaderHash.get(bam);
tjc's avatar
tjc committed
    File bamIndexFile = getBamIndexFile(bam);
      try
      {
        logger4j.warn("Index file not found so using picard to index the BAM.");
        // Use Picard to index the file
        // requires reads to be sorted by coordinate
        new BuildBamIndex().instanceMain(
tcarver's avatar
tcarver committed
          new String[]{ "I="+bam, "O="+bamIndexFile.getAbsolutePath(), "MAX_RECORDS_IN_RAM=50000", "VALIDATION_STRINGENCY=SILENT" });
      }
      catch(SAMException e)
      {
        String ls = System.getProperty("line.separator");
        String msg = 
            "BAM index file is missing. The BAM file needs to be sorted and indexed"+ls+
            "This can be done using samtools (http://samtools.sf.net/):"+ls+ls+
            "samtools sort <in.bam> <out.prefix>"+ls+
            "samtools index <sorted.bam>";
        
        throw new SAMException(msg);
      }
tjc's avatar
tjc committed
    final SAMFileReader samFileReader;
    if(feature_display != null && bam.endsWith("cram"))
    {
tcarver's avatar
tcarver committed
      // set log level
      net.sf.picard.util.Log.setGlobalLogLevel( 
          net.sf.picard.util.Log.LogLevel.ERROR);
      final CRAMReferenceSequenceFile ref = new CRAMReferenceSequenceFile(
        feature_display.getEntryGroup().getSequenceEntry(), this);
      
      final Map<Object, ReferenceSequenceFile> referenceFactory = 
          new HashMap<Object, ReferenceSequenceFile>();
      referenceFactory.put(bamIndexFile, ref);

      try
      {
tcarver's avatar
tcarver committed
        Class<?> cls = getClass().getClassLoader().loadClass("net.sf.samtools.ReferenceDiscovery");
        Field f = cls.getDeclaredField("referenceFactory");
        f.set(null, referenceFactory);
      }
      catch (ClassNotFoundException e)
      {
        System.err.println("Check cramtools.jar is in the CLASSPATH. "+e.getMessage());
      }
      catch (SecurityException e)
      {
        e.printStackTrace();
      }
      catch (NoSuchFieldException e)
      {
        e.printStackTrace();
      }
      catch (IllegalArgumentException e)
      {
        e.printStackTrace();
      }
      catch (IllegalAccessException e)
      {
        e.printStackTrace();
      }

      
      //net.sf.samtools.ReferenceDiscovery.referenceFactory.put(bamIndexFile, ref);
    }
    
    if(bam.startsWith("ftp"))
    {
      FTPSeekableStream fss = new FTPSeekableStream(new URL(bam));
      samFileReader = new SAMFileReader(fss, bamIndexFile, false);
    }
    else if(!bam.startsWith("http"))
tjc's avatar
tjc committed
    {
      File bamFile = new File(bam);
tjc's avatar
tjc committed
      samFileReader = new SAMFileReader(bamFile, bamIndexFile);
tjc's avatar
tjc committed
    }
    else
    {
      final URL urlBamFile = new URL(bam);
      samFileReader = new SAMFileReader(urlBamFile, bamIndexFile, false);
tjc's avatar
tjc committed
    }
tjc's avatar
tjc committed
    samFileReader.setValidationStringency(ValidationStringency.SILENT);
    samFileReaderHash.put(bam, samFileReader);
tcarver's avatar
tcarver committed

    readGroups.addAll(samFileReader.getFileHeader().getReadGroups());
tjc's avatar
tjc committed
    return samFileReader;
tjc's avatar
tjc committed
  }
tjc's avatar
tjc committed

tjc's avatar
tjc committed
  private void readHeaderPicard() throws IOException
tjc's avatar
tjc committed
  {
tcarver's avatar
tcarver committed
    final SAMFileReader inputSam = getSAMFileReader(bamList.get(0));
    final SAMFileHeader header = inputSam.getFileHeader();

    for(SAMSequenceRecord seq: header.getSequenceDictionary().getSequences())
tjc's avatar
tjc committed
    {
tcarver's avatar
tcarver committed
      seqLengths.put(seq.getSequenceName(),
                     seq.getSequenceLength());
      seqNames.add(seq.getSequenceName());
tjc's avatar
tjc committed
    }
tjc's avatar
tjc committed
  }
  
  class BamReadTask implements Runnable 
  {
    private int start; 
    private int end; 
    private short bamIndex; 
    private float pixPerBase;
    private CountDownLatch latch;
    BamReadTask(int start, int end, short bamIndex, float pixPerBase, CountDownLatch latch)  
    {
      this.start = start;
      this.end = end;
      this.bamIndex = bamIndex;
      this.pixPerBase = pixPerBase;
      this.latch = latch;
    }

    public void run() 
    {
      try
      {
        readFromBamPicard(start, end, bamIndex, pixPerBase) ;
      }
      catch (OutOfMemoryError ome)
      {
        throw ome;
      }
      catch(IOException me)
      {
        me.printStackTrace();
      }
      finally
      {
        latch.countDown();
      }
    }
  }
tjc's avatar
tjc committed

tjc's avatar
tjc committed
  /**
   * Read a SAM or BAM file.
tjc's avatar
tjc committed
   * @throws IOException 
tjc's avatar
tjc committed
   */
  private void readFromBamPicard(int start, int end, short bamIndex, float pixPerBase) 
tjc's avatar
tjc committed
          throws IOException
tjc's avatar
tjc committed
  {
    // Open the input file.  Automatically detects whether input is SAM or BAM
    // and delegates to a reader implementation for the appropriate format.
tcarver's avatar
tcarver committed
    final String bam = bamList.get(bamIndex);
tjc's avatar
tjc committed
    final SAMFileReader inputSam = getSAMFileReader(bam);
    
    //final SAMFileReader inputSam = new SAMFileReader(bamFile, indexFile);
    if(isConcatSequences())
tjc's avatar
tjc committed
    {
tcarver's avatar
tcarver committed
      for(String seq: seqNames)
tjc's avatar
tjc committed
      {
tcarver's avatar
tcarver committed
        int sLen = seqLengths.get(seq);
        int offset = getSequenceOffset(seq); 
        int sBeg = offset+1;
        int sEnd = sBeg+sLen-1;
tjc's avatar
tjc committed

        if( (sBeg >= start && sBeg <= end) ||
            (sBeg <= start && sEnd >= start) )
tjc's avatar
tjc committed
        {
          int thisStart = start - offset;
          if(thisStart < 1)
            thisStart = 1;
          int thisEnd   = end - offset;
tcarver's avatar
tcarver committed
          iterateOverBam(inputSam, seq, thisStart, thisEnd, bamIndex, pixPerBase, bam);
          //System.out.println("READ "+seq+"  "+thisStart+".."+thisEnd+" "+start+" --- "+offset);
tjc's avatar
tjc committed
        }
      }
    }
    else
    {
      String refName = (String) combo.getSelectedItem();
      iterateOverBam(inputSam, refName, start, end, bamIndex, pixPerBase, bam);
tjc's avatar
tjc committed
    }
tjc's avatar
tjc committed
    //inputSam.close();
tjc's avatar
tjc committed
  }
  
  /**
   * Iterate over BAM file and load into the <code>List</code> of
   * <code>SAMRecord</code>.
   * @param inputSam
   * @param refName
   * @param start
   * @param end
   */
  private void iterateOverBam(final SAMFileReader inputSam, 
                             final String refName, final int start, final int end,
                             final short bamIndex, final float pixPerBase,
                             final String bam)
  {
    final MemoryMXBean memory = ManagementFactory.getMemoryMXBean();
tcarver's avatar
tcarver committed
    final int checkMemAfter = 8000;
    final int seqOffset = getSequenceOffset(refName);
    final int offset = seqOffset- getBaseAtStartOfView();
    final boolean isCoverageView = isCoverageView(pixPerBase);
    int binSize = ((end-start)/nbins)+1;
    if(binSize < 1)
    {
      binSize = 1;
      nbins = end-start+1;
    }
    int max = MAX_COVERAGE*binSize;
    if(max < 1)
      max = Integer.MAX_VALUE;
    final int cov[] = new int[nbins];
    for(int i=0; i<nbins; i++)
      cov[i] = 0;
tcarver's avatar
tcarver committed

    final CloseableIterator<SAMRecord> it = inputSam.queryOverlapping(refName, start, end);
    try
tjc's avatar
tjc committed
    {
      while ( it.hasNext() )
tjc's avatar
tjc committed
      {
          cnt++;
          SAMRecord samRecord = it.next();

tcarver's avatar
tcarver committed
          if(readGrpFrame != null && !readGrpFrame.isReadGroupVisible(samRecord.getReadGroup()))
            continue;
          
          if( samRecordFlagPredicate == null ||
             !samRecordFlagPredicate.testPredicate(samRecord))
            if(samRecordMapQPredicate == null ||
               samRecordMapQPredicate.testPredicate(samRecord))
            {
              int abeg = samRecord.getAlignmentStart();
              int aend = samRecord.getAlignmentEnd();
              boolean over = false;
              
              for(int i=abeg; i<aend; i++)
              {
                int bin = ((i-start)/binSize)-1;
                if(bin < 0)
tcarver's avatar
tcarver committed
                  continue;
                else if(bin > nbins-1)
                  bin = nbins-1;
                cov[bin]++;
                if(cov[bin] > max)
                {
                  over = true;
                  break;
                }
              }
             
              if(over)
                continue;

              if(isCoverageView)
                coverageView.addRecord(samRecord, offset, bam, colourByStrandTag.isSelected());
                coveragePanel.addRecord(samRecord, offset, bam, colourByStrandTag.isSelected());
              if(isSNPplot)
                snpPanel.addRecord(samRecord, seqOffset);
              if(!isCoverageView)
                readsInView.add(new BamViewRecord(samRecord, bamIndex));
            }
          if(cnt > checkMemAfter)
tjc's avatar
tjc committed
          {
            cnt = 0;
            float heapFraction =
              (float)((float)memory.getHeapMemoryUsage().getUsed()/
                      (float)memory.getHeapMemoryUsage().getMax());
            logger4j.debug("Heap memory usage (used/max): "+heapFraction);
          
            if(readsInView.size() > checkMemAfter*2 && !waitingFrame.isVisible())
              waitingFrame.showWaiting("loading...", mainPanel);

            if(heapFraction > 0.90) 
            {
              popFrame.show(
                "Using > 90 % of the maximum memory limit:"+
                (memory.getHeapMemoryUsage().getMax()/1000000.f)+" Mb.\n"+
                "Not all reads in this range have been read in. Zoom in or\n"+
                "consider increasing the memory for this application.",
                mainPanel,
                15000);
              break;
            }
tjc's avatar
tjc committed
          }
        catch(Exception e)
        {
          System.err.println(e.getMessage());
        }
tjc's avatar
tjc committed
      }
    }
tjc's avatar
tjc committed
  }

  private int getSequenceLength()
  {
    if(isConcatSequences())
tjc's avatar
tjc committed
    {
      int len = 0;
tcarver's avatar
tcarver committed
      for(String seq: seqNames)
        len += seqLengths.get(seq);
tjc's avatar
tjc committed
      return len;
    }
    else
      return seqLengths.get((String) combo.getSelectedItem());
  }
  
tjc's avatar
tjc committed
  /**
   * For BAM files with multiple references sequences, calculate
   * the offset from the start of the concatenated sequence for 
   * a given reference.
   * @param refName
   * @return
   */
tjc's avatar
tjc committed
  protected int getSequenceOffset(String refName)
  {
    if(!isConcatSequences())
tjc's avatar
tjc committed
      return 0;
    {
      if(feature_display == null)
        offsetLengths = new HashMap<String, Integer>(combo.getItemCount());
        int offset = 0;
        for(int i=0; i<combo.getItemCount(); i++)
        {
          String thisSeqName = (String) combo.getItemAt(i);
          offsetLengths.put(thisSeqName, offset);
          offset += seqLengths.get(combo.getItemAt(i));
        }
        return offsetLengths.get(refName);
      }
      final FeatureVector features = feature_display.getEntryGroup().getAllFeatures();
      final HashMap<String, Integer> lookup = new HashMap<String, Integer>();
      for(int i=0; i<features.size(); i++)
        lookup.put(features.elementAt(i).getIDString(), features.elementAt(i).getFirstBase());
        
      offsetLengths = new HashMap<String, Integer>(seqNames.size());
      for(int i=0; i<seqNames.size(); i++)
      {
        final Integer pos = lookup.get(seqNames.get(i));
        if(pos != null)
          offsetLengths.put(seqNames.get(i), pos-1);
       /*final FeatureContigPredicate predicate = new FeatureContigPredicate(seqNames.get(i).trim());
        for(int j=0; j<features.size(); j++)
        {
          if(predicate.testPredicate(features.elementAt(j)))
          {
            offsetLengths.put(seqNames.get(i), features.elementAt(j).getFirstBase()-1);
            break;
          }
      
      if(offsetLengths.size() != seqNames.size())
tjc's avatar
tjc committed
      {
tjc's avatar
tjc committed
        System.err.println("Found: "+offsetLengths.size() +" of "+ seqNames.size());
        SwingUtilities.invokeLater(new Runnable() 
        {
          public void run() 
          {
            JOptionPane.showMessageDialog(BamView.this, 
                "There is a problem matching the reference sequences\n"+
                "to the names in the BAM file. This may mean the labels\n"+
                "on the reference features do not match those in the in\n"+
                "the BAM file.", 
                "Problem Found", JOptionPane.WARNING_MESSAGE);
          }
        });
        //concatSequences = false;
        int offset = 0;
        for(int i=0; i<combo.getItemCount(); i++)
        {
          String thisSeqName = (String) combo.getItemAt(i);
          offsetLengths.put(thisSeqName, offset);
          offset += seqLengths.get(combo.getItemAt(i));
        }
        //return 0;
tjc's avatar
tjc committed
      }
tjc's avatar
tjc committed
    }
    return offsetLengths.get(refName);
tjc's avatar
tjc committed
  }
tjc's avatar
tjc committed
  /**
   * Override
   */
  protected void paintComponent(Graphics g)
tjc's avatar
tjc committed
  {
	super.paintComponent(g);
	Graphics2D g2 = (Graphics2D)g;

	mouseOverSAMRecord = null;
tcarver's avatar
tcarver committed
    final int seqLength = getSequenceLength();
    int start;
tjc's avatar
tjc committed
    int end;
    
    if(startBase > 0)
      start = startBase;
    else
tjc's avatar
tjc committed
      start = getBaseAtStartOfView();
    
    if(endBase > 0)
      end = endBase;
    else
tjc's avatar
tjc committed
    {
tjc's avatar
tjc committed
      end   = start + nbasesInView - 1;
      if(end > seqLength)
        end = seqLength;
      
      if(feature_display != null && nbasesInView < feature_display.getMaxVisibleBases())
        nbasesInView = feature_display.getMaxVisibleBases();
tjc's avatar
tjc committed
    }
tcarver's avatar
tcarver committed
    final float pixPerBase = getPixPerBaseByWidth();
    boolean changeToStackView = false;
    MemoryMXBean memory = ManagementFactory.getMemoryMXBean();
tjc's avatar
tjc committed
    if(laststart != start ||
       lastend   != end ||
       CoveragePanel.isRedraw())
tjc's avatar
tjc committed
    {
      if(isCoverageView(pixPerBase))
        coverageView.init(this, pixPerBase, start, end);
      if(isCoverage)
        coveragePanel.init(this, pixPerBase, start, end);
      if(isSNPplot)
        snpPanel.init(this, pixPerBase, start, end);
tjc's avatar
tjc committed
      {
        try
        {
          float heapFractionUsedBefore = (float) ((float) memory.getHeapMemoryUsage().getUsed() / 
                                                  (float) memory.getHeapMemoryUsage().getMax());
          if(readsInView == null)
            readsInView = new Vector<BamViewRecord>();
          else
            readsInView.clear();
tcarver's avatar
tcarver committed
          final CountDownLatch latch = new CountDownLatch(bamList.size()-hideBamList.size());
          //long ms = System.currentTimeMillis();
          for(short i=0; i<bamList.size(); i++)
          {
            if(!hideBamList.contains(i))
              bamReadTaskExecutor.execute(
                  new BamReadTask(start, end, i, pixPerBase, latch));

          try 
          {
            latch.await();
          }
          catch (InterruptedException e) {} // TODO 

          //System.out.println("===== NO. THREADS="+
          //     ((java.util.concurrent.ThreadPoolExecutor)bamReadTaskExecutor).getPoolSize()+" TIME="+(System.currentTimeMillis()-ms));

          float heapFractionUsedAfter = (float) ((float) memory.getHeapMemoryUsage().getUsed() / 
                                                 (float) memory.getHeapMemoryUsage().getMax());

          // System.out.println("Heap Max  : "+memory.getHeapMemoryUsage().getMax());
          // System.out.println("Heap Used : "+memory.getHeapMemoryUsage().getUsed());
          // System.out.println("Heap memory used "+heapFractionUsedAfter);

          if ((heapFractionUsedAfter - heapFractionUsedBefore) > 0.06
              && !isStackView() && heapFractionUsedAfter > 0.8)
            cbStackView.setSelected(true);
          if((!isStackView() && !isStrandStackView()) || isBaseAlignmentView(pixPerBase))
          {
            Collections.sort(readsInView, new SAMRecordComparator());
          }
          else if( (isStackView() || isStrandStackView()) &&
tjc's avatar
tjc committed
              bamList.size() > 1)
          {
            // merge multiple BAM files
            Collections.sort(readsInView, new SAMRecordPositionComparator(BamView.this));
tjc's avatar
tjc committed
          }