From 55c95cdb772969bde9a1952745778b607b6864ed Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Hannes=20P=C3=A9tur=20Eggertsson?= <hannespetur@gmail.com>
Date: Wed, 18 Oct 2023 13:23:18 +0000
Subject: [PATCH] modified skyr such that indels are represented together

---
 include/paw/align/alignment_options.hpp | 10 ++++-----
 include/paw/align/event.hpp             |  5 +++++
 include/paw/align/skyr.hpp              | 28 ++++++++++++-------------
 include/paw/align/variant.hpp           |  7 +++++++
 4 files changed, 30 insertions(+), 20 deletions(-)

diff --git a/include/paw/align/alignment_options.hpp b/include/paw/align/alignment_options.hpp
index 934a3c8..f7f4c05 100644
--- a/include/paw/align/alignment_options.hpp
+++ b/include/paw/align/alignment_options.hpp
@@ -280,11 +280,11 @@ set_query(AlignmentOptions<Tuint> & opt, AlignmentCache<Tuint> & aln_cache, Tseq
 
   // TODO add more tests
   //// set free snps
-  //for (Event2 const & e : opt.free_edits)
-  //{
-  //  assert(e.is_snp());
-  //  aln_cache.set_free_snp(e.pos, e.alt[0]);
-  //}
+  for (Event2 const & e : opt.free_edits)
+  {
+    assert(e.is_snp());
+    aln_cache.set_free_snp(e.pos, e.alt[0]);
+  }
 
   aln_cache.vH_up = Tvec_pack(static_cast<std::size_t>(aln_cache.num_vectors),
                               static_cast<Tpack>(simdpp::make_int(2 * aln_cache.gap_open_val +
diff --git a/include/paw/align/event.hpp b/include/paw/align/event.hpp
index f1de4c4..8fc6156 100644
--- a/include/paw/align/event.hpp
+++ b/include/paw/align/event.hpp
@@ -45,6 +45,11 @@ class Event2
     return ref.size() == 0;
   }
 
+  inline bool
+  is_indel() const
+  {
+    return is_deletion() || is_insertion();
+  }
 
 };
 
diff --git a/include/paw/align/skyr.hpp b/include/paw/align/skyr.hpp
index 5711afe..eae4fb3 100644
--- a/include/paw/align/skyr.hpp
+++ b/include/paw/align/skyr.hpp
@@ -335,33 +335,31 @@ Skyr::find_variants_from_edits()
 
     if (e.pos == new_var.pos)
     {
-      if (e.is_insertion() && new_var.is_insertion())
+      if (e.is_indel() && new_var.is_indel())
       {
-        // Case 1: insertions at the same position.
-        new_var.add_event(e);
-      }
-      else if (e.is_deletion() && new_var.is_deletion())
-      {
-        // Case 2: Deletions at the same position.
+        // Case 1: indel at the same position.
         assert(e.ref.size() > new_var.seqs[0].size());
-        /// Extend all sequences.
-        auto begin_extend = e.ref.cbegin() + new_var.seqs[0].size();
 
-        for (auto & s : new_var.seqs)
-          s.append(begin_extend, e.ref.cend());
-        ///
+        // Extend all variant sequences if needed
+        {
+          auto const begin_extend = std::next(e.ref.cbegin(), new_var.seqs[0].size());
+          auto const end_extend = e.ref.cend();
+
+          for (auto & s : new_var.seqs)
+            s.append(begin_extend, end_extend);
+        }
 
         new_var.add_event(e);
       }
       else if (e.is_snp() && new_var.is_snp())
       {
-        // Case 3: both are SNPs at the same position.
+        // Case 2: both are SNPs at the same position.
         new_var.add_event(e);
       }
       else
       {
         // Otherwise, they are two different events.
-        vars.push_back(std::move(new_var));
+        vars.push_back(std::move(new_var)); // push old variant
         new_var = {static_cast<uint32_t>(e.pos), {e.ref} };
         new_var.add_event(e);
       }
@@ -369,7 +367,7 @@ Skyr::find_variants_from_edits()
     else
     {
       // Otherwise, they are two different events.
-      vars.push_back(std::move(new_var));
+      vars.push_back(std::move(new_var)); // push old variant
       new_var = {static_cast<uint32_t>(e.pos), {e.ref} };
       new_var.add_event(e);
     }
diff --git a/include/paw/align/variant.hpp b/include/paw/align/variant.hpp
index e02268c..e72c72c 100644
--- a/include/paw/align/variant.hpp
+++ b/include/paw/align/variant.hpp
@@ -31,6 +31,7 @@ class Variant
   bool has_sequences() const;
   bool is_deletion() const;
   bool is_insertion() const;
+  bool is_indel() const;
   bool is_snp() const;
   long get_max_del_reach() const;
   long get_max_reach() const;
@@ -92,6 +93,12 @@ Variant::is_insertion() const
   return has_sequences() && seqs[0].size() == 0;
 }
 
+bool
+Variant::is_indel() const
+{
+  return is_deletion() || is_insertion();
+}
+
 
 bool
 Variant::is_snp() const