copy motifs to a new subanalysis

ticket:199
Also the second problem was that the motif_scan wasn't using the right
offsets to compute the motif start/stop locations.

diff --git a/alg/sequence.cpp b/alg/sequence.cpp
index 270cc4e9ed5900db1eed1f7c4e37978fa47b9409..7538e4649cb2852d468e1549a53469d60a5c94b1 100644 (file)
--- a/alg/sequence.cpp
+++ b/alg/sequence.cpp
@@ -921,12 +921,12 @@ Sequence::motif_scan(const Sequence& a_motif, std::vector<int> * motif_match_sta
   Sequence::value_type motif_char;
   Sequence::value_type seq_char;
 
-  while (seq_i < seq->length())
+  while (seq_i < size())
   {
     // this is pretty much a straight translation of Nora's python code
     // to match iupac letter codes
     motif_char = toupper(a_motif[motif_i]);
-    seq_char = toupper(seq->at(seq_i));
+    seq_char = toupper(seq->at(seq_start+seq_i));
     if (motif_char =='N')
       motif_i++;
     else if (motif_char == seq_char)

diff --git a/alg/test/test_sequence.cpp b/alg/test/test_sequence.cpp
index 6317144c026510016c8e7f1fe376f1105aab6d65..a3c63b6d87d8c4710a4cb76777a0e9b6775b711f 100644 (file)
--- a/alg/test/test_sequence.cpp
+++ b/alg/test/test_sequence.cpp
@@ -612,6 +612,37 @@ BOOST_AUTO_TEST_CASE( sequence_motifs )
   */
 }
 
+BOOST_AUTO_TEST_CASE( sequence_motif_subseq)
+{
+  // when searching for a motif on a subsequence we should 
+  // only search the subsequence ticket:199
+  string aaaa("AAAA");
+  string cccc("CCCC");
+  Sequence s1("AAAANCCCC", Sequence::reduced_dna_alphabet);
+
+  // this shouldn't show up
+  s1.add_motif(cccc);
+  BOOST_CHECK_EQUAL( s1.motifs().size(), 1 );
+
+  s1.add_motif(aaaa);
+  BOOST_CHECK_EQUAL( s1.motifs().size(), 2 );
+
+  Sequence subseq1 = s1.subseq(4,5);
+  BOOST_CHECK_EQUAL(subseq1.motifs().size(), 2);
+  subseq1.clear_motifs();
+  BOOST_CHECK_EQUAL(subseq1.motifs().size(), 0);
+  // this is outside of our subsequence, and so shouldn't be found    
+  subseq1.add_motif(aaaa);
+  BOOST_CHECK_EQUAL( subseq1.motifs().size(), 0 );
+  
+  subseq1.add_motif(cccc);
+  BOOST_CHECK_EQUAL( subseq1.motifs().size(), 1);
+  std::list<motif>::const_iterator motif_i = subseq1.motifs().begin();
+  BOOST_REQUIRE(motif_i != subseq1.motifs().end());
+  BOOST_CHECK_EQUAL(motif_i->begin, 1);
+  BOOST_CHECK_EQUAL(motif_i->end, 5);
+}
+
 BOOST_AUTO_TEST_CASE( annot_test )
 {
   annot a(0, 10, "test", "thing");

diff --git a/qui/SubanalysisWindow.cpp b/qui/SubanalysisWindow.cpp
index 6d5f3ce99e80277181d02150c7492f22fcb3dba8..57dafc9792d4e494599ec1bb9c84e82d0cfb991e 100644 (file)
--- a/qui/SubanalysisWindow.cpp
+++ b/qui/SubanalysisWindow.cpp
@@ -84,7 +84,7 @@ void SubanalysisWindow::run()
                              "empty model now.");
   }
 
-  MussaRef m(new Mussa);
+  MussaRef new_m(new Mussa);
   
   for(SequenceLocationModel::iterator itor = model.begin();
       itor != model.end();
@@ -92,14 +92,30 @@ void SubanalysisWindow::run()
   {
     // append_sequence from a const Sequence & will make a copy 
     // for the shared pointer.
-    m->append_sequence(itor->getSelectedSequence());
+    Sequence s = itor->getSelectedSequence();
+    s.clear_motifs();
+    new_m->append_sequence(s);
   }
 
   try {
-    m->set_window(window->value());
-    m->set_threshold(threshold->value());
-    m->analyze();
-    MussaWindow *mw = new MussaWindow(m);
+    new_m->set_window(window->value());
+    new_m->set_threshold(threshold->value());
+    new_m->analyze();
+    
+    // copy over motifs
+    std::vector<Sequence> motifs_copy;
+    std::vector<Color> color_copy;
+    const Mussa::motif_set motifs = analysis->motifs();
+    boost::shared_ptr<AnnotationColors> mapper = analysis->colorMapper();
+    for(Mussa::motif_set::const_iterator motif_i = motifs.begin();
+        motif_i != motifs.end();
+        ++motif_i)
+    {
+      motifs_copy.push_back(*motif_i);
+      color_copy.push_back(mapper->lookup("motif", motif_i->get_sequence()));
+    }
+    new_m->set_motifs(motifs_copy, color_copy);
+    MussaWindow *mw = new MussaWindow(new_m);
     mw->show();
     model.clear();
     hide();