render the sequences from the mussa analysis
authorDiane Trout <diane@caltech.edu>
Tue, 28 Feb 2006 07:03:09 +0000 (07:03 +0000)
committerDiane Trout <diane@caltech.edu>
Tue, 28 Feb 2006 07:03:09 +0000 (07:03 +0000)
I implemented a GlSequence which attempts to handle the opengl rendering
for a sequence object. Also to make Sequence a bit more consistent
with STL I changed len() to length() and added size() (and changed
the code using it to size()

I added some unit tests for GlSequence as I was not remotely confident
that const_cast was actually the right way to implement the assignment
operator= when one of the member variables is a const reference.

The updateAnalysis code currently is what is handling the layout of the
various sequence objects.

16 files changed:
alg/mussa_class.cxx
alg/mussa_class.hh
alg/sequence.cxx
alg/sequence.hh
gui/ConnView.cxx
mussagl.pro
qui/GlSequence.cxx [new file with mode: 0644]
qui/GlSequence.h [new file with mode: 0644]
qui/PathScene.cxx
qui/PathScene.h
qui/PathWindow.cxx
qui/PathWindow.h
test/module.mk
test/test_glsequence.cxx [new file with mode: 0644]
test/test_mussa.cxx
test/test_sequence.cxx

index 2efdcddc80fb234916ef6b52da1f831adb1fed5d..8992446af61b8f8aa2359efbeb0a924d351ff7fe 100644 (file)
@@ -153,6 +153,12 @@ Mussa::set_seq_info(string seq_file, string annot_file, int fa_i, int a_start,
   sub_seq_ends.push_back(the_end);
 }
 
+const vector<Sequence>& 
+Mussa::sequences() const
+{
+  return the_seqs;
+}
+
 string
 Mussa::load_mupa_file(string para_file_path)
 {
@@ -414,7 +420,7 @@ Mussa::seqcomp()
       all_comps[i].push_back(dummy_comp);
   }
   for(i = 0; i < the_seqs.size(); i++)
-    seq_lens.push_back(the_seqs[i].len());
+    seq_lens.push_back(the_seqs[i].size());
 
   for(i = 0; i < the_seqs.size(); i++)
     for(i2 = i+1; i2 < the_seqs.size(); i2++)
@@ -584,7 +590,7 @@ Mussa::load(string ana_file)
     for (i = 1; i <= the_paths.seq_num(); i++)
     {
       tmp_seq.clear();
-      cout << "mussa_class: loading the fucking museq frag...\n";
+      cout << "mussa_class: loading museq frag...\n";
       tmp_seq.load_museq(a_file_path, i);
       the_seqs.push_back(tmp_seq);
     }
index c85ceed04e134b91129008b523f389cd339d323c..dcbda559abd39022ac9d182ca02870c2ae3841ff 100644 (file)
@@ -87,6 +87,8 @@ class Mussa
     // sets info to load a seq and annotations from a fasta file 
     void set_seq_info(std::string seq_file, std::string annot_file, 
                       int fa_i, int a_start, int the_end);
+    //! allow examining the sequences we have loaded
+    const std::vector<Sequence>& sequences() const;
 
     // deprecated - support bridge for python version of mussa
     // these save & load from the old file format
index fb68dbed7e7ae59d2c53dd7ac9976d40d82bf01c..5d150e3eb22d07144fd04fca4c4549924fa427dd 100644 (file)
@@ -244,11 +244,19 @@ Sequence::load_annot(string file_path, int start_index, int end_index)
   }
 }
 
+const std::list<annot> Sequence::annotations() const
+{
+  return annots;
+}
+
+string::size_type Sequence::length() const
+{
+  return size();
+}
 
-int 
-Sequence::len() const
+string::size_type Sequence::size() const
 {
-  return sequence.length();
+  return sequence.size();
 }
 
 
index a48584372253fd04a458f3bede191de4ae1f2a95..5c6329617b02a7f45cfda172dafd35795c40a442 100644 (file)
@@ -48,6 +48,8 @@ class Sequence
   public:
     Sequence();
     Sequence(std::string seq);
+    //! assignment to constant sequences
+    //Sequence &operator=(const Sequence&);
     //! load sequence AGCT from fasta file
     //! throws sequence_load_error if something goes wrong
     void load_fasta(const std::string file_path, int seq_num=1, 
@@ -55,12 +57,16 @@ class Sequence
     //! load sequence annotations
     //! throws sequence_load_error if something goes wrong
     void load_annot(const std::string file_path, int start_index, int end_index);
+    const std::list<annot> annotations() const;
     void set_seq(const std::string& a_seq);
     const std::string& get_seq() const;
     const char * c_seq() const;
     std::string subseq(int start, int end) const;
     std::string rev_comp() const;
-    int len() const;
+    //! the number of base pairs in this sequence
+    std::string::size_type size() const;
+    //! alias for size, (included as this is much like a string)
+    std::string::size_type length() const;
     const std::string& get_header() const;
     //std::string sp_name() const;
     std::vector<int> find_motif(std::string a_motif);
index e0390a31586f9befc258e8f56668e46605875630..c11db9932b0e380c4a8b21b7d1f95593d8b1916e 100644 (file)
@@ -44,8 +44,8 @@ ConnView::setup(string name, int sq_num, int win_len,
   for(i = 0; i < seq_num; ++i)
   {
     a_seq = (*S)[i];
-    cout << a_seq.len() << endl;
-    seq_length = a_seq.len();
+    cout << a_seq.size() << endl;
+    seq_length = a_seq.size();
     seq_lens.push_back(seq_length);
     if (seq_length > max_seq_len)
       max_seq_len = seq_length;
index 62efcaed1cad6f84cecb1dc3b8aada0025f22e65..3fe90529553854b2b0f43df7d6b4f6ec1d27c1e7 100644 (file)
@@ -12,17 +12,19 @@ INCLUDEPATH += . alg qui
 
 # Input
 HEADERS += mussa_exceptions.hh \
+           qui/GlSequence.h \
            qui/PathWindow.h \
            qui/PathScene.h \
-          qui/ThresholdWidget.h \
+           qui/ThresholdWidget.h \
            alg/flp.hh \
            alg/mussa_class.hh \
            alg/nway_paths.hh \
            alg/sequence.hh
 SOURCES += mussagl.cxx \
+           qui/GlSequence.cxx \
            qui/PathWindow.cxx \
            qui/PathScene.cxx \
-          qui/ThresholdWidget.cxx \
+           qui/ThresholdWidget.cxx \
            alg/flp.cxx \
            alg/flp_seqcomp.cxx \
            alg/mussa_class.cxx \
diff --git a/qui/GlSequence.cxx b/qui/GlSequence.cxx
new file mode 100644 (file)
index 0000000..c206bd4
--- /dev/null
@@ -0,0 +1,98 @@
+#include "qui/GlSequence.h"
+
+#include <iostream>
+using namespace std;
+
+GlSequence::GlSequence(const Sequence &s) 
+  : seq(s),
+    seq_x(0), 
+    seq_y(0), 
+    seq_z(0), 
+    seq_width(s.size()),
+    seq_height(gl_track_height) 
+{
+}
+
+GlSequence::GlSequence(const GlSequence &s) 
+  : seq(s.seq),
+    seq_x(s.seq_x),
+    seq_y(s.seq_y),
+    seq_z(s.seq_z),
+    seq_width(s.seq_width),
+    seq_height(s.seq_height)
+{
+}
+
+GlSequence &GlSequence::operator=(const GlSequence & s)
+{
+  if (this != &s) {
+    const_cast<Sequence &>(seq) = s.seq;
+    seq_x = s.seq_x;
+    seq_y = s.seq_y;
+    seq_z = s.seq_z;
+    seq_width = s.seq_width;
+    seq_height = s.seq_height;
+  }
+  return *this;
+}
+
+const Sequence &GlSequence::sequence()
+{
+  return seq;
+}
+
+void GlSequence::setX(GLfloat value)
+{
+  seq_x = value;
+}
+
+GLfloat GlSequence::x()
+{
+  return seq_x;
+}
+
+void GlSequence::setY(GLfloat value)
+{
+  seq_y = value;
+}
+
+GLfloat GlSequence::y()
+{
+  return seq_y;
+}
+
+void GlSequence::setWidth(GLfloat value)
+{
+  seq_width = value;
+}
+
+GLfloat GlSequence::width()
+{
+  return seq_width;
+}
+
+
+void GlSequence::draw() const
+{
+  glLineWidth(seq_height);
+  glColor3f(0.0, 0.0, 0.0);
+  // draw main sequence track
+  glBegin(GL_LINES);
+    glVertex3f(seq_x,           seq_y, seq_z);
+    glVertex3f(seq_x+seq_width, seq_y, seq_z);
+    clog << "drawing " << seq_x << " " << seq_y << " " << seq_width
+         << std::endl;
+  glEnd();
+  // draw annotations
+  std::list<annot> annots = seq.annotations();
+  for (std::list<annot>::const_iterator annot_itor = annots.begin();
+       annot_itor != annots.end();
+       ++annot_itor)
+  {
+    glColor3f(0.0, 0.5, 0.0);
+    glBegin(GL_LINES);
+      glVertex3f(annot_itor->start, seq_y, seq_z);
+      glVertex3f(annot_itor->end  , seq_y, seq_z);
+    glEnd();
+  }
+}
diff --git a/qui/GlSequence.h b/qui/GlSequence.h
new file mode 100644 (file)
index 0000000..519465b
--- /dev/null
@@ -0,0 +1,39 @@
+#ifndef _GL_SEQUENCE_H_
+#define _GL_SEQUENCE_H_
+
+#include "alg/sequence.hh"
+#include <GL/gl.h>
+
+//! Manage rendering a mussa sequence track
+/*! The idea is this will keep track of the location of where the sequence
+ *  is being rendered, and handle displaying annotations on that track
+ */
+class GlSequence
+{
+public: 
+  GlSequence(const Sequence & s);
+  GlSequence(const GlSequence & s);
+  GlSequence &operator=(const GlSequence &s);
+
+  void draw() const;
+
+  const Sequence &sequence();
+  void setX(GLfloat);
+  GLfloat x();
+  void setY(GLfloat);
+  GLfloat y();
+  void setWidth(GLfloat);
+  GLfloat width();
+
+private:
+  const Sequence& seq;
+  GLfloat seq_x;
+  GLfloat seq_y;
+  GLfloat seq_z;
+  GLfloat seq_width;
+  GLfloat seq_height;
+};
+
+const float gl_track_height = 5.0;
+
+#endif
index 8e78c61294fade580c3ed51fb51241c995689a2c..b221b7b0502b0181612d260143fccda077c1d2f3 100644 (file)
@@ -8,6 +8,10 @@
 #include <GL/gl.h>
 #include <math.h>
 
+#include "qui/GlSequence.h"
+
+using namespace std;
+
 PathScene::PathScene(int frags, int len, QWidget *parent) : 
   QGLWidget(parent),
   X(0),
@@ -21,6 +25,11 @@ PathScene::PathScene(int frags, int len, QWidget *parent) :
   debugBand(true),
   drawingBand(false)
 { 
+  // Hack in loading an analysis 
+  //mussaAnalysis.load_mupa_file( "examples/mck3test.mupa" );
+  //mussaAnalysis.analyze();
+  mussaAnalysis.load("mck3test_w30_t20");
+  updateAnalysis();
 }
 
 QSize PathScene::sizeHint() const
@@ -87,7 +96,7 @@ void PathScene::setClipPlane(int newZ)
 void PathScene::initializeGL()
 {
   glEnable(GL_DEPTH_TEST|GL_BLEND);
-  glClearColor(0.0, 0.0, 0.0, 0.0);
+  glClearColor(1.0, 1.0, 1.0, 0.0);
   glBlendFunc(GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA);
   glShadeModel(GL_FLAT);
   glLoadIdentity();
@@ -130,6 +139,31 @@ void PathScene::paintGL()
   glFlush();
 }
 
+void PathScene::updateAnalysis()
+{
+  GlSequence *gl_seq;
+  float y = mussaAnalysis.sequences().size() * 100 ;
+  float max_base_pairs = 0;
+  // Delete old glsequences
+  // FIXME: does this actually free the memory from the new'ed GlSequences?
+  tracks.clear();
+  // save a reference to our GlSequences
+  for(vector<Sequence>::const_iterator seq_itor = mussaAnalysis.sequences().begin();
+      seq_itor != mussaAnalysis.sequences().end();
+      ++seq_itor)
+  {
+    y = y - 100;
+    gl_seq = new GlSequence(*seq_itor);
+    gl_seq->setX(0);
+    gl_seq->setY(y);
+    if (gl_seq->width() > max_base_pairs ) max_base_pairs = gl_seq->width();
+    tracks.push_back( *gl_seq );
+  }
+  maxOrtho2d.setWidth(max_base_pairs + 100.0);
+  maxOrtho2d.setHeight((mussaAnalysis.sequences().size()) * 100 );
+  curOrtho2d = maxOrtho2d;
+}
+
 static void processHits(GLint hits, GLuint buffer[])
 {
   GLuint names, *ptr;
@@ -240,15 +274,14 @@ static void seq_track(long length, long height)
   glEnd();
 }
 
-static void tracks()
+void PathScene::draw_tracks()
 {
   glPushMatrix();
-  glLoadName(1);
-  seq_track(2000000, 0);
-  glLoadName(2);
-  seq_track(2500000, 100);
-  glLoadName(3);
-  seq_track(1900000, 200);
+  for (vector<GlSequence>::size_type i = 0; i != tracks.size(); ++i )
+  {
+    glLoadName(i);
+    tracks[i].draw();
+  }
   glPopMatrix();
 }
 
@@ -262,8 +295,6 @@ static void draw_lines(int fragsize, int length)
   GLfloat z;
   const long short_sequence_len = 1800000;
   long blocks = short_sequence_len / fragsize;
-  std::cout << "blocks " << blocks << std::endl;
-  std::cout << "line count " << blocks * length << std::endl;
   for (long base_block=0; base_block < blocks; ++base_block)
   {
     for (long base_offset=0; base_offset < length; ++base_offset)
@@ -304,9 +335,9 @@ void PathScene::mussaesque(int fragments, int length)
 
   glInitNames();
   glPushName(0);
-  tracks();
+  draw_tracks();
   glLoadName(10);
-  glCallList(theLines);
+  //glCallList(theLines);
   glPopName();
   //draw_lines();
 }
index 77737abcf26aad455e323c82b5532a944d1342ed..7f6b014007a764b570f831bac8869b0ab192a6c3 100644 (file)
@@ -4,6 +4,10 @@
 #include <QGLWidget>
 #include <QRectF>
 #include <QPoint>
+#include <vector>
+
+#include "alg/mussa_class.hh"
+#include "qui/GlSequence.h"
 
 class QMouseEvent;
 class QRubberBand;
@@ -19,10 +23,19 @@ public:
 
   QSize sizeHint() const;
 
+  Mussa mussaAnalysis;
+  std::vector<GlSequence> tracks;
+
 public slots:
   void setX(int x);
   void setClipPlane(int z);
   void setZoom(int);
+  //! called when we've changed the analysis
+  void updateAnalysis();
+
+signals:
+  //! emitted when our analysis has changed
+  void analysisUpdated();
 
 protected:
   int X;
@@ -37,6 +50,12 @@ protected:
   void resizeGL(int width, int height);
   void paintGL();
 
+  void mussaesque(int, int);
+  //! draw all of our sequence tracks
+  void draw_tracks();
+
+
+  //! \defgroup Selection
   QRubberBand *rubberBand;
   QPoint bandOrigin;
   QRectF previousBand;
@@ -46,6 +65,5 @@ protected:
   void mouseMoveEvent(QMouseEvent *);
   void mouseReleaseEvent(QMouseEvent *);
 
-  void mussaesque(int, int);
 };
 #endif
index 488a3f9db03740921624ff16da9bf379ae076578..92beaa67a08715cff8dbc5fb5ef2757909e25fa4 100644 (file)
@@ -162,7 +162,7 @@ void PathWindow::loadMupa()
                                                    caption, 
                                                    QDir::currentPath(),
                                                    filter);
-  mussaAnalysis.load_mupa_file(mupa_path.toStdString());
+  scene->mussaAnalysis.load_mupa_file(mupa_path.toStdString());
 }
 
 void PathWindow::loadSavedAnalysis()
index 86d73525749b7d83328d85a96814f3d1255c5371..1aefe0cda3ce9fc54204fd1fb601e5c6856eea66 100644 (file)
@@ -3,8 +3,6 @@
 
 #include <QMainWindow>
 
-#include "alg/mussa_class.hh"
-
 class QAction;
 class PathScene;
 
@@ -13,7 +11,6 @@ class PathWindow : public QMainWindow
   Q_OBJECT
 
 public: 
-  Mussa mussaAnalysis;
   PathWindow(QWidget *parent=0);
 
 public slots:
index 5fb4ede6d9d8cfd2e06b42bb0774626222e39257..d418c41dfc9b3c88ba3040fa71cc0faaf423d968 100644 (file)
@@ -1,10 +1,12 @@
 CURDIR := $(BASEDIR)test/
 
-SOURCES.cxx := test_main.cxx \
-               test_flp.cxx \
+SOURCES.cxx := test_flp.cxx \
+               test_glsequence.cxx \
+               test_main.cxx \
                test_mussa.cxx \
                test_nway.cxx \
-               test_sequence.cxx
+               test_sequence.cxx \
+               ../qui/GlSequence.cxx
 
 TESTSRC := $(addprefix $(CURDIR), $(SOURCES.cxx))
 
@@ -15,5 +17,5 @@ TEST := $(BASEDIR)/unittests$(BINEXT)
 TARGETBINS += $(TEST)
 
 $(TEST): $(TESTSRC:.cxx=$(OBJEXT))  $(MUSSA_ALG_LIB)
-       g++ $(CXXFLAGS) -lboost_unit_test_framework -lboost_filesystem -o $@ $^
+       g++ $(CXXFLAGS) -lGL -lboost_unit_test_framework -lboost_filesystem -o $@ $^
 
diff --git a/test/test_glsequence.cxx b/test/test_glsequence.cxx
new file mode 100644 (file)
index 0000000..0c46812
--- /dev/null
@@ -0,0 +1,31 @@
+#include <boost/test/auto_unit_test.hpp>
+
+#include <string>
+
+#include "qui/GlSequence.h"
+#include "alg/sequence.hh"
+
+using namespace std;
+
+BOOST_AUTO_TEST_CASE ( glsequence_operator_equal )
+{
+  // I don't trust my operator = hack so lets make sure it works.
+  string s0("AAGGCCTT");
+  string s1("TTGGCCAA");
+  Sequence seq0(s0);
+  Sequence seq1(s1);
+
+  GlSequence glseq0(seq0);
+  BOOST_CHECK (glseq0.sequence().get_seq() == s0);
+  // width of a sequence should be number of base pairs (aka chars)
+  BOOST_CHECK (glseq0.width() == s0.size());
+  GlSequence glseq1(seq1);
+  GlSequence glseq_copy0(glseq0);
+
+  BOOST_CHECK(glseq_copy0.sequence().get_seq() == glseq0.sequence().get_seq());
+  BOOST_CHECK( &(glseq_copy0.sequence()) == &(glseq0.sequence()));
+
+  glseq0 = glseq1;
+
+  BOOST_CHECK( glseq0.sequence().get_seq() == s1 );
+}
index 608c6a297eaf464d311513aaab43d2ffca314686..645fe3ed2b8dcc759a7a4c56480f0f47eca1ec3f 100644 (file)
@@ -1,5 +1,7 @@
 #include <boost/test/auto_unit_test.hpp>
 
+#include <string>
+
 #include "alg/mussa_class.hh"
 
 //! can we initialize a mussa object?
@@ -39,6 +41,23 @@ BOOST_AUTO_TEST_CASE( mussa_analysis_name )
   BOOST_CHECK_EQUAL( m.get_analysis_mode_name(), "[deprecated] Recursive" );
 }
 
+BOOST_AUTO_TEST_CASE( mussa_sequences )
+{
+  std::string s0("AAAANNNN");
+  std::string s1("GGGGNNNN");
+  std::string s2("TTTTNNNN");
+
+  Mussa analysis;
+  analysis.add_a_seq(s0);
+  analysis.add_a_seq(s1);
+  analysis.add_a_seq(s2);
+
+  BOOST_CHECK_EQUAL( analysis.sequences().size(), 3 );
+  BOOST_CHECK_EQUAL( analysis.sequences()[0].get_seq(), s0);
+  BOOST_CHECK_EQUAL( analysis.sequences()[1].get_seq(), s1);
+  BOOST_CHECK_EQUAL( analysis.sequences()[2].get_seq(), s2);
+}
+
 #include <unistd.h>
 BOOST_AUTO_TEST_CASE( mussa_load_mupa )
 {
index ded5fe797590b5e8be30a535ad1bd678bd02d0a2..cc9a9a3256d4cc51517be2b218f026891bceb131 100644 (file)
@@ -21,7 +21,7 @@ BOOST_AUTO_TEST_CASE( sequence_filter )
   Sequence s2("aattggcc");
   BOOST_CHECK_EQUAL(s2.get_seq(), "AATTGGCC");
   BOOST_CHECK_EQUAL(s2.rev_comp(), "GGCCAATT");
-  BOOST_CHECK_EQUAL(s2.len(), s2.get_seq().size());
+  BOOST_CHECK_EQUAL(s2.size(), s2.get_seq().size());
   BOOST_CHECK_EQUAL(s2.c_seq(), s2.get_seq().c_str());
 
   Sequence s3("asdfg");