[med-svn] [Git][med-team/last-align][master] 4 commits: New upstream version 1282

Nilesh Patra (@nilesh) gitlab at salsa.debian.org
Wed May 4 11:25:32 BST 2022



Nilesh Patra pushed to branch master at Debian Med / last-align


Commits:
f4a6eafe by Nilesh Patra at 2022-05-04T14:40:05+05:30
New upstream version 1282
- - - - -
38cb2130 by Nilesh Patra at 2022-05-04T14:40:12+05:30
Update upstream source from tag 'upstream/1282'

Update to upstream version '1282'
with Debian dir d4c311804b8c5154b1935fe27cc2f91e9abd76dd
- - - - -
88d8ded5 by Nilesh Patra at 2022-05-04T15:08:15+05:30
Refresh patches

- - - - -
04dde6a2 by Nilesh Patra at 2022-05-04T15:10:05+05:30
Upload to unstable

- - - - -


19 changed files:

- debian/changelog
- debian/patches/simde
- doc/FAQ.rst
- doc/bisulfite.rst
- doc/last-papers.rst
- doc/last-split.rst
- doc/lastal.rst
- examples/last-bisulfite.sh
- src/lastal.cc
- src/makefile
- src/mcf_simd.hh
- src/split/last-split.cc
- src/tantan.cc
- test/SRR359290-1k.maf
- test/bs100.maf
- test/last-split-test.out
- test/maf-convert-test.out
- test/maf-swap-test.out
- test/maf-swap-test.sh


Changes:

=====================================
debian/changelog
=====================================
@@ -1,3 +1,10 @@
+last-align (1282-1) unstable; urgency=medium
+
+  * New upstream version 1282
+  * Refresh patches
+
+ -- Nilesh Patra <nilesh at debian.org>  Wed, 04 May 2022 15:09:54 +0530
+
 last-align (1268-1) unstable; urgency=medium
 
   * New upstream version 1268


=====================================
debian/patches/simde
=====================================
@@ -26,15 +26,17 @@ equivalents automatically
 -
  typedef __m256i SimdInt;
  typedef __m256i SimdUint1;
- 
-@@ -146,279 +141,6 @@ static inline SimdInt simdChoose1(SimdIn
+ typedef __m256d SimdDbl;
+@@ -181,358 +176,6 @@
+ static inline SimdInt simdChoose1(SimdInt items, SimdInt choices) {
    return _mm256_shuffle_epi8(items, choices);
  }
- 
+-
 -#elif defined __SSE4_1__
 -
 -typedef __m128i SimdInt;
 -typedef __m128i SimdUint1;
+-typedef __m128d SimdDbl;
 -
 -const int simdBytes = 16;
 -
@@ -46,6 +48,10 @@ equivalents automatically
 -  return _mm_setzero_si128();
 -}
 -
+-static inline SimdDbl simdZeroDbl() {
+-  return _mm_setzero_pd();
+-}
+-
 -static inline SimdInt simdOnes1() {
 -  return _mm_set1_epi32(-1);
 -}
@@ -58,6 +64,10 @@ equivalents automatically
 -  return _mm_loadu_si128((const SimdInt *)p);
 -}
 -
+-static inline SimdDbl simdLoadDbl(const double *p) {
+-  return _mm_loadu_pd(p);
+-}
+-
 -static inline void simdStore(void *p, SimdInt x) {
 -  _mm_storeu_si128((SimdInt *)p, x);
 -}
@@ -66,6 +76,10 @@ equivalents automatically
 -  _mm_storeu_si128((SimdInt *)p, x);
 -}
 -
+-static inline void simdStoreDbl(double *p, SimdDbl x) {
+-  _mm_storeu_pd(p, x);
+-}
+-
 -static inline SimdInt simdOr1(SimdInt x, SimdInt y) {
 -  return _mm_or_si128(x, y);
 -}
@@ -75,6 +89,7 @@ equivalents automatically
 -}
 -
 -const int simdLen = 4;
+-const int simdDblLen = 2;
 -
 -static inline SimdInt simdSet(int i3, int i2, int i1, int i0) {
 -  return _mm_set_epi32(i3, i2, i1, i0);
@@ -88,6 +103,10 @@ equivalents automatically
 -		      i7, i6, i5, i4, i3, i2, i1, i0);
 -}
 -
+-static inline SimdDbl simdSetDbl(double i1, double i0) {
+-  return _mm_set_pd(i1, i0);
+-}
+-
 -static inline SimdInt simdFill(int x) {
 -  return _mm_set1_epi32(x);
 -}
@@ -96,6 +115,10 @@ equivalents automatically
 -  return _mm_set1_epi8(x);
 -}
 -
+-static inline SimdDbl simdFillDbl(double x) {
+-  return _mm_set1_pd(x);
+-}
+-
 -static inline SimdInt simdGt(SimdInt x, SimdInt y) {
 -  return _mm_cmpgt_epi32(x, y);
 -}
@@ -116,6 +139,10 @@ equivalents automatically
 -  return _mm_adds_epu8(x, y);
 -}
 -
+-static inline SimdDbl simdAddDbl(SimdDbl x, SimdDbl y) {
+-  return _mm_add_pd(x, y);
+-}
+-
 -static inline SimdInt simdSub(SimdInt x, SimdInt y) {
 -  return _mm_sub_epi32(x, y);
 -}
@@ -124,6 +151,10 @@ equivalents automatically
 -  return _mm_sub_epi8(x, y);
 -}
 -
+-static inline SimdDbl simdMulDbl(SimdDbl x, SimdDbl y) {
+-  return _mm_mul_pd(x, y);
+-}
+-
 -static inline SimdInt simdQuadruple1(SimdInt x) {
 -  return _mm_slli_epi32(x, 2);
 -}
@@ -148,6 +179,10 @@ equivalents automatically
 -  return _mm_extract_epi16(x, 0);
 -}
 -
+-static inline double simdHorizontalAddDbl(SimdDbl x) {
+-  return _mm_cvtsd_f64(_mm_hadd_pd(x, x));
+-}
+-
 -static inline SimdInt simdChoose1(SimdInt items, SimdInt choices) {
 -  return _mm_shuffle_epi8(items, choices);  // SSSE3
 -}
@@ -157,6 +192,7 @@ equivalents automatically
 -typedef int32x4_t SimdInt;
 -typedef uint32x4_t SimdUint;
 -typedef uint8x16_t SimdUint1;
+-typedef float64x2_t SimdDbl;
 -
 -const int simdBytes = 16;
 -
@@ -168,6 +204,10 @@ equivalents automatically
 -  return vdupq_n_u8(0);
 -}
 -
+-static inline SimdDbl simdZeroDbl() {
+-  return vdupq_n_f64(0);
+-}
+-
 -static inline SimdUint1 simdOnes1() {
 -  return vdupq_n_u8(-1);
 -}
@@ -180,6 +220,10 @@ equivalents automatically
 -  return vld1q_u8(p);
 -}
 -
+-static inline SimdDbl simdLoadDbl(const double *p) {
+-  return vld1q_f64(p);
+-}
+-
 -static inline void simdStore(int *p, SimdInt x) {
 -  vst1q_s32(p, x);
 -}
@@ -188,6 +232,10 @@ equivalents automatically
 -  vst1q_u8(p, x);
 -}
 -
+-static inline void simdStoreDbl(double *p, SimdDbl x) {
+-  vst1q_f64(p, x);
+-}
+-
 -static inline SimdUint1 simdOr1(SimdUint1 x, SimdUint1 y) {
 -  return vorrq_u8(x, y);
 -}
@@ -197,6 +245,7 @@ equivalents automatically
 -}
 -
 -const int simdLen = 4;
+-const int simdDblLen = 2;
 -
 -static inline SimdInt simdSet(unsigned i3, unsigned i2,
 -                              unsigned i1, unsigned i0) {
@@ -225,6 +274,10 @@ equivalents automatically
 -  return vcombine_u8(vcreate_u8(lo), vcreate_u8(hi));
 -}
 -
+-static inline SimdDbl simdSetDbl(double i1, double i0) {
+-  return vcombine_f64(vdup_n_f64(i0), vdup_n_f64(i1));
+-}
+-
 -static inline SimdInt simdFill(int x) {
 -  return vdupq_n_s32(x);
 -}
@@ -233,6 +286,10 @@ equivalents automatically
 -  return vdupq_n_u8(x);
 -}
 -
+-static inline SimdDbl simdFillDbl(double x) {
+-  return vdupq_n_f64(x);
+-}
+-
 -static inline SimdUint simdGt(SimdInt x, SimdInt y) {
 -  return vcgtq_s32(x, y);
 -}
@@ -253,6 +310,10 @@ equivalents automatically
 -  return vqaddq_u8(x, y);
 -}
 -
+-static inline SimdDbl simdAddDbl(SimdDbl x, SimdDbl y) {
+-  return vaddq_f64(x, y);
+-}
+-
 -static inline SimdInt simdSub(SimdInt x, SimdInt y) {
 -  return vsubq_s32(x, y);
 -}
@@ -261,6 +322,10 @@ equivalents automatically
 -  return vsubq_u8(x, y);
 -}
 -
+-static inline SimdDbl simdMulDbl(SimdDbl x, SimdDbl y) {
+-  return vmulq_f64(x, y);
+-}
+-
 -static inline SimdUint1 simdQuadruple1(SimdUint1 x) {
 -  return vshlq_n_u8(x, 2);
 -}
@@ -281,6 +346,10 @@ equivalents automatically
 -  return vminvq_u8(x);
 -}
 -
+-static inline double simdHorizontalAddDbl(SimdDbl x) {
+-  return vaddvq_f64(x);
+-}
+-
 -static inline SimdUint1 simdChoose1(SimdUint1 items, SimdUint1 choices) {
 -  return vqtbl1q_u8(items, choices);
 -}
@@ -288,19 +357,29 @@ equivalents automatically
 -#else
 -
 -typedef int SimdInt;
+-typedef double SimdDbl;
 -const int simdBytes = 1;
 -const int simdLen = 1;
+-const int simdDblLen = 1;
 -static inline int simdZero() { return 0; }
+-static inline double simdZeroDbl() { return 0; }
 -static inline int simdSet(int x) { return x; }
+-static inline double simdSetDbl(double x) { return x; }
 -static inline int simdFill(int x) { return x; }
 -static inline int simdLoad(const int *p) { return *p; }
+-static inline double simdLoadDbl(const double *p) { return *p; }
 -static inline void simdStore(int *p, int x) { *p = x; }
+-static inline void simdStoreDbl(double *p, double x) { *p = x; }
+-static inline double simdFillDbl(double x) { return x; }
 -static inline int simdGt(int x, int y) { return x > y; }
 -static inline int simdAdd(int x, int y) { return x + y; }
+-static inline double simdAddDbl(double x, double y) { return x + y; }
 -static inline int simdSub(int x, int y) { return x - y; }
+-static inline double simdMulDbl(double x, double y) { return x * y; }
 -static inline int simdMax(int x, int y) { return x > y ? x : y; }
 -static inline int simdBlend(int x, int y, int mask) { return mask ? y : x; }
 -static inline int simdHorizontalMax(int a) { return a; }
+-static inline double simdHorizontalAddDbl(double x) { return x; }
 -
 -#endif
 -
@@ -309,7 +388,7 @@ equivalents automatically
  #endif
 --- a/src/GappedXdropAligner.cc
 +++ b/src/GappedXdropAligner.cc
-@@ -140,17 +140,13 @@ int GappedXdropAligner::align(const ucha
+@@ -140,17 +140,13 @@
      if (isAffine) {
        for (int i = 0; i < numCells; i += simdLen) {
  	SimdInt s = simdSet(
@@ -329,7 +408,7 @@ equivalents automatically
  	SimdInt y = simdSub(simdLoad(y1+i), mDelGrowCost);
 --- a/src/GappedXdropAlignerPssm.cc
 +++ b/src/GappedXdropAlignerPssm.cc
-@@ -91,17 +91,13 @@ int GappedXdropAligner::alignPssm(const
+@@ -91,17 +91,13 @@
      if (isAffine) {
        for (int i = 0; i < numCells; i += simdLen) {
  	SimdInt s = simdSet(
@@ -363,7 +442,7 @@ equivalents automatically
  CXXFLAGS += -std=c++11
  CXXFLAGS += -pthread -DHAS_CXX_THREAD
  # -fomit-frame-pointer ?
-@@ -56,11 +55,12 @@ split/last-split.o
+@@ -56,11 +55,12 @@
  PPOBJ = last-pair-probs.o last-pair-probs-main.o
  
  MBOBJ = last-merge-batches.o
@@ -380,7 +459,7 @@ equivalents automatically
  
  indexObj5 = $(indexObj4:.o=.o5)
  alignObj5 = $(alignObj4:.o=.o5)
-@@ -74,45 +74,45 @@ all: $(ALL)
+@@ -74,45 +74,45 @@
  last8: $(LAST8)
  
  indexAllObj4 = $(indexObj0) $(indexObj4)
@@ -448,7 +527,7 @@ equivalents automatically
  //#include <iostream>  // for debugging
  
  namespace cbrc {
-@@ -43,12 +41,10 @@ int GappedXdropAligner::alignDna(const u
+@@ -43,12 +41,10 @@
  
    const SimdUint1 scorer4x4 =
      simdSet1(
@@ -461,7 +540,7 @@ equivalents automatically
  		 scorer[3][3], scorer[3][2], scorer[3][1], scorer[3][0],
  		 scorer[2][3], scorer[2][2], scorer[2][1], scorer[2][0],
  		 scorer[1][3], scorer[1][2], scorer[1][1], scorer[1][0],
-@@ -126,7 +122,6 @@ int GappedXdropAligner::alignDna(const u
+@@ -126,7 +122,6 @@
  
        for (int i = 0; i < numCells; i += simdBytes) {
  	SimdUint1 s = simdSet1(
@@ -469,7 +548,7 @@ equivalents automatically
  			     scorer[s1[31]][s2[31]],
  			     scorer[s1[30]][s2[30]],
  			     scorer[s1[29]][s2[29]],
-@@ -143,7 +138,6 @@ int GappedXdropAligner::alignDna(const u
+@@ -143,7 +138,6 @@
  			     scorer[s1[18]][s2[18]],
  			     scorer[s1[17]][s2[17]],
  			     scorer[s1[16]][s2[16]],
@@ -477,7 +556,7 @@ equivalents automatically
  			     scorer[s1[15]][s2[15]],
  			     scorer[s1[14]][s2[14]],
  			     scorer[s1[13]][s2[13]],
-@@ -275,5 +269,3 @@ bool GappedXdropAligner::getNextChunkDna
+@@ -275,5 +269,3 @@
  }
  
  }
@@ -485,7 +564,7 @@ equivalents automatically
 -#endif
 --- a/src/Alignment.cc
 +++ b/src/Alignment.cc
-@@ -365,12 +365,10 @@ void Alignment::extend( std::vector< Seg
+@@ -365,12 +365,10 @@
  				  del.openCost, del.growCost,
  				  ins.openCost, ins.growCost,
  				  gap.pairCost, gap.isAffine, maxDrop, smMax)
@@ -498,7 +577,7 @@ equivalents automatically
      :           aligner.align(seq1, seq2, isForward, globality, sm,
  			      del.openCost, del.growCost,
  			      ins.openCost, ins.growCost,
-@@ -387,14 +385,12 @@ void Alignment::extend( std::vector< Seg
+@@ -387,14 +385,12 @@
        while( greedyAligner.getNextChunk( end1, end2, size ) )
  	chunks.push_back( SegmentPair( end1 - size, end2 - size, size ) );
      }
@@ -515,7 +594,7 @@ equivalents automatically
  				   del.openCost, del.growCost,
 --- a/src/GappedXdropAligner.hh
 +++ b/src/GappedXdropAligner.hh
-@@ -352,7 +352,6 @@ class GappedXdropAligner {
+@@ -352,7 +352,6 @@
    void initFrame();
  
    // Everything below here is for alignDna & getNextChunkDna
@@ -523,7 +602,7 @@ equivalents automatically
    std::vector<TinyScore> xTinyScores;
    std::vector<TinyScore> yTinyScores;
    std::vector<TinyScore> zTinyScores;
-@@ -402,7 +401,6 @@ class GappedXdropAligner {
+@@ -402,7 +401,6 @@
      while (*x2 != target) ++x2;
      bestSeq1position = x2 - x2beg + seq1beg;
    }
@@ -531,3 +610,33 @@ equivalents automatically
  };
  
  }
+--- a/src/tantan.cc
++++ b/src/tantan.cc
+@@ -324,13 +324,9 @@
+     int i = 0;
+     for (; i <= maxOffset - simdDblLen; i += simdDblLen) {
+       SimdDbl rV = simdSetDbl(
+-#if defined __SSE4_1__ || defined __ARM_NEON
+-#ifdef __AVX2__
+ 			      lrRow[sp[-i-4]],
+ 			      lrRow[sp[-i-3]],
+-#endif
+ 			      lrRow[sp[-i-2]],
+-#endif
+ 			      lrRow[sp[-i-1]]);
+       SimdDbl fV = simdLoadDbl(fp+i);
+       sV = simdAddDbl(sV, fV);
+@@ -368,13 +364,9 @@
+     int i = 0;
+     for (; i <= maxOffset - simdDblLen; i += simdDblLen) {
+       SimdDbl rV = simdSetDbl(
+-#if defined __SSE4_1__ || defined __ARM_NEON
+-#ifdef __AVX2__
+ 			      lrRow[sp[-i-4]],
+ 			      lrRow[sp[-i-3]],
+-#endif
+ 			      lrRow[sp[-i-2]],
+-#endif
+ 			      lrRow[sp[-i-1]]);
+       SimdDbl fV = simdMulDbl(simdLoadDbl(fp+i), rV);
+       sV = simdAddDbl(sV, simdMulDbl(simdLoadDbl(b2f+i), fV));


=====================================
doc/FAQ.rst
=====================================
@@ -1,6 +1,26 @@
 LAST FAQ
 ========
 
+:Q: Why does LAST sometimes miss strong alignments?
+
+:A: This is because it avoids getting huge numbers of repetitive
+    alignments, by just finding a few strongest matches for each part
+    of each query sequence.  To find more repetitive matches, increase
+    the ``lastal -m`` parameter, e.g. ``-m100``.
+
+..
+
+:Q: Why does LAST sometimes show two separate alignments, and
+    sometimes one merged alignment, between the same sequences?
+
+:A: This is because LAST forbids highly dissimilar regions inside
+    alignments, but the maximum-allowed dissimilarity depends on the
+    alignment score / E-value threshold (which depends on the database
+    size).  For more details, please see the discussion of "X-drop" in
+    `Parameters for accurate genome alignment`_.
+
+..
+
 :Q: Does it matter which sequence is used as the "reference" (given to
     lastdb) and which is used as the "query" (given to lastal)?
 
@@ -43,15 +63,6 @@ LAST FAQ
 
 ..
 
-:Q: I'd like to compare my queries to a database of known proteins.
-    Where can I get a database of known proteins?
-
-:A: You could try UniRef90 or UniRef50
-    (http://www.uniprot.org/help/uniref), which have reduced
-    redundancy.
-
-..
-
 :Q: How does LAST get the sequence names?  How can I get nice, short,
     unique names?
 
@@ -88,3 +99,5 @@ RAQ (Rarely Asked Questions)
     command with::
 
       numactl --interleave=all
+
+.. _Parameters for accurate genome alignment: https://doi.org/10.1186/1471-2105-11-80


=====================================
doc/bisulfite.rst
=====================================
@@ -4,8 +4,7 @@ Aligning bisulfite-converted DNA reads to a genome
 The easy way
 ------------
 
-Use `Bisulfighter <http://epigenome.cbrc.jp/bisulfighter>`_, which
-wraps LAST.
+Use Bisulfighter_, which wraps LAST.
 
 The hard way
 ------------
@@ -81,3 +80,5 @@ This assumes that reads1.fastq are all from the converted strand
 reverse-complement (i.e. they have G->A conversions).  It requires GNU
 parallel to be installed.  You are encouraged to customize this
 script.
+
+.. _Bisulfighter: https://github.com/yutaka-saito/Bisulfighter


=====================================
doc/last-papers.rst
=====================================
@@ -127,7 +127,7 @@ research to society.
 15. `Improved DNA-versus-protein homology search for protein fossils`__.
     Yao Y, Frith MC.
 
-    __ https://doi.org/10.1101/2021.01.25.428050
+    __ https://doi.org/10.1007/978-3-030-74432-8_11
 
     Describes "new-style" DNA-versus-protein search with
     ``last-train --codon``.


=====================================
doc/last-split.rst
=====================================
@@ -1,14 +1,22 @@
 last-split
 ==========
 
-This program finds "split alignments" (typically for DNA) or "spliced
+last-split finds "split alignments" (typically for DNA) or "spliced
 alignments" (typically for RNA).
 
-It reads candidate alignments of query sequences to a genome, and
-looks for a unique best alignment for each part of each query.  It
-allows different parts of one query to match different parts of the
-genome.  This is useful for DNA queries that cross rearrangement
-breakpoints, or RNA queries that cross splice junctions.
+It reads candidate alignments of query sequences to a genome, and cuts
+from them a unique best alignment for each part of each query.  This
+is useful for DNA queries that cross rearrangement breakpoints, or RNA
+queries that cross splice junctions.
+
+For a detailed explanation of what last-split does, please see `A
+survey of localized sequence rearrangements in human DNA
+<https://doi.org/10.1093/nar/gkx1266>`_ (update in `A pipeline for
+complete characterization of complex germline rearrangements from long
+DNA reads <https://doi.org/10.1186/s13073-020-00762-1>`_).  The
+algorithm, and application to whole genomes, is in `Split-alignment of
+genomes finds orthologies more accurately
+<https://doi.org/10.1186/s13059-015-0670-9>`_.
 
 Examples
 --------
@@ -34,12 +42,13 @@ tells it where the splice signals are (GT, AG, etc)::
 
   lastal -p train.out -D10 db q.fastq | last-split -g db > out.maf
 
-This will favour splices starting at GT (and to a lesser extent GC and
-AT), and ending at AG (and to a lesser extent AC).  However, it allows
-splices starting and ending anywhere.  It also favours splices with
-introns of typical length, specified by a log-normal distribution
-(i.e. cis-splices).  However, it allows arbitrary trans-splices
-between any two places in the genome.
+This will favor splices starting at GT (and to a lesser extent GC and
+AT), and ending at AG (and to a lesser extent AC).  The output shows
+the donor (``don``) and acceptor (``acc``) dinucleotides.
+
+It also favors splices with introns of typical length, specified by a
+log-normal distribution (cis-splices).  However, it allows arbitrary
+trans-splices between any two places in the genome.
 
 ``-D10`` sets a very loose significance threshold, so that we can find
 very short parts of a spliced alignment (e.g. short exons).  Note that
@@ -148,19 +157,20 @@ Faster spliced alignment
 ~~~~~~~~~~~~~~~~~~~~~~~~
 
 Spliced alignment can be slow.  It can be sped up, at a small cost in
-accuracy, by not favouring cis-splices::
+accuracy, by not favoring cis-splices::
 
   lastal -p train.out -D10 db q.fastq | last-split -c0 -t0.004 -g db > out.maf
 
 The ``-c0`` turns off cis-splicing, and the ``-t0.004`` specifies a
-higher probability of trans-splicing.
+higher probability of "trans-splicing" (which includes cis-splicing,
+just doesn't favor it).
 
 "Spliced" alignment of DNA reads to a genome
 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 
 If we do not wish to allow arbitrarily large unaligned parts in the
 middle of the query, we can do "spliced" alignment without considering
-splice signals or favouring cis-splices::
+splice signals or favoring cis-splices::
 
   lastal -p train.out db q.fastq | last-split -c0 > out.maf
 
@@ -193,12 +203,33 @@ Options
        determines whether forward and/or reverse-complement splice
        signals are used.
 
-       If you use ``-d2``, the output will have an extra "sense"
+       If you use ``-d2``, the output will have an extra ``sense``
        field, indicating the log-odds that the query is
        sense-stranded::
 
 	   log2[ prob(sense) / prob(antisense) ]
 
+       The donor and acceptor annotations also indicate strandedness.
+       This order means that the query sequence is from the forward
+       strand of a spliced RNA::
+
+           don
+           acc don
+           acc don
+           acc
+
+       And this order means it is from the reverse strand::
+
+           acc
+           don acc
+           don acc
+           don
+
+       It's possible for ``don`` and ``acc`` to conflict with
+       ``sense``.  That happens if one orientation is overall more
+       likely (indicated by ``sense``), but any one alignment with
+       that orientation is not the most likely.
+
 -c, --cis=PROB
        Do spliced alignment, and set the average probability per
        base of cis-splicing.  The default value roughly fits human
@@ -260,10 +291,10 @@ Details
   (of the kind produced by lastal) describing the alignment score
   parameters.
 
-* The program reads one batch of alignments at a time (by looking for
-  lines starting with "# batch").  If the batches are huge
-  (e.g. because there are no lines starting with "# batch"), it might
-  need too much memory.
+* The input must not mix alignments of different query sequences.  In
+  other words, all the alignments of one query must be next to each
+  other.  If you use ``-r``/``--reverse``, however, there is no such
+  requirement, and the whole input gets read into memory.
 
 * lastal can optionally write "p" lines, indicating the probability
   that each base is misaligned due to wrong gap placement.
@@ -285,8 +316,8 @@ last-split5
 -----------
 
 last-split5 is almost identical to last-split.  The only difference is
-the -g option: last-split can only read the output of lastdb, whereas
-last-split5 can only read the output of lastdb5_.
+the ``-g`` option: last-split can only read the output of lastdb,
+whereas last-split5 can only read the output of lastdb5_.
 
 Limitations
 -----------


=====================================
doc/lastal.rst
=====================================
@@ -318,7 +318,7 @@ Miscellaneous options
     single sequence exceeds this amount, however, it is not split.
     You can use suffixes K, M, and G to specify KibiBytes,
     MebiBytes, and GibiBytes.  This option has no effect on the
-    results (apart from their order).
+    results.
 
     If the reference was split into volumes by lastdb_, then each
     volume will be read into memory once per query batch.


=====================================
examples/last-bisulfite.sh
=====================================
@@ -30,8 +30,8 @@ trap 'rm -f $tmp.*' EXIT
 # Convert C to t, and all other letters to uppercase:
 perl -pe 'y/Cca-z/ttA-Z/ if $. % 4 == 2' "$@" > "$tmp".q
 
-lastal -pBISF -s1 -Q1 -D1000 "$my_f" "$tmp".q > "$tmp".f
-lastal -pBISR -s0 -Q1 -D1000 "$my_r" "$tmp".q > "$tmp".r
+lastal -pBISF -s1 -Q1 -D1000 -i1 "$my_f" "$tmp".q > "$tmp".f
+lastal -pBISR -s0 -Q1 -D1000 -i1 "$my_r" "$tmp".q > "$tmp".r
 
 last-merge-batches "$tmp".f "$tmp".r | last-split -m0.1 -fMAF+ |
 perl -F'(\s+)' -ane '$F[12] =~ y/ta/CG/ if /^s/ and $s++ % 2; print @F'


=====================================
src/lastal.cc
=====================================
@@ -93,6 +93,7 @@ namespace {
   sequenceFormat::Enum referenceFormat = sequenceFormat::fasta;
   int minScoreGapless;
   int isCaseSensitiveSeeds = -1;  // initialize it to an "error" value
+  unsigned numOfVolumes = -1;
   unsigned numOfIndexes = 1;  // assume this value, if unspecified
 }
 
@@ -341,11 +342,10 @@ static void calculateScoreStatistics(const std::string& matrixName,
 }
 
 // Read the .prj file for the whole database
-void readOuterPrj( const std::string& fileName, unsigned& volumes,
-                   size_t& refMinimizerWindow, size_t& minSeedLimit,
-		   bool& isKeepRefLowercase, int& refTantanSetting,
-                   countT& refSequences, countT& refLetters,
-		   countT& refMaxSeqLen ){
+void readOuterPrj(const std::string &fileName, size_t &refMinimizerWindow,
+		  size_t &minSeedLimit, bool &isKeepRefLowercase,
+		  int &refTantanSetting, countT &refSequences,
+		  countT &refLetters, countT &refMaxSeqLen) {
   std::ifstream f( fileName.c_str() );
   if( !f ) ERR( "can't open file: " + fileName );
   unsigned version = 0;
@@ -371,7 +371,7 @@ void readOuterPrj( const std::string& fileName, unsigned& volumes,
     if( word == "masklowercase" ) iss >> isCaseSensitiveSeeds;
     if( word == "sequenceformat" ) iss >> referenceFormat;
     if( word == "minimizerwindow" ) iss >> refMinimizerWindow;
-    if( word == "volumes" ) iss >> volumes;
+    if( word == "volumes" ) iss >> numOfVolumes;
     if( word == "numofindexes" ) iss >> numOfIndexes;
     if( word == "integersize" ) iss >> fileBitsPerInt;
   }
@@ -419,15 +419,15 @@ void readInnerPrj( const std::string& fileName,
 }
 
 // Write match counts for each query sequence
-void writeCounts( std::ostream& out ){
+void writeCounts() {
   for( indexT i = 0; i < matchCounts.size(); ++i ){
-    out << query.seqName(i) << '\n';
+    std::cout << query.seqName(i) << '\n';
 
     for( size_t j = args.minHitDepth; j < matchCounts[i].size(); ++j ){
-      out << j << '\t' << matchCounts[i][j] << '\n';
+      std::cout << j << '\t' << matchCounts[i][j] << '\n';
     }
 
-    out << '\n';  // blank line afterwards
+    std::cout << '\n';  // blank line afterwards
   }
 }
 
@@ -544,7 +544,7 @@ struct Dispatcher{
 
 static bool isCollatedAlignments() {
   return args.outputFormat == 'b' || args.outputFormat == 'B' ||
-    args.cullingLimitForFinalAlignments + 1;
+    args.cullingLimitForFinalAlignments + 1 || numOfVolumes > 1;
 }
 
 static void printAndDelete(char *text) {
@@ -913,11 +913,6 @@ static void printAndClear(std::vector<AlignmentText> &textAlns) {
   textAlns.clear();
 }
 
-static void printAndClearAll() {
-  for (size_t i = 0; i < aligners.size(); ++i)
-    printAndClear(aligners[i].textAlns);
-}
-
 void makeQualityPssm( LastAligner& aligner,
 		      size_t queryNum, const SubstitutionMatrices &matrices,
 		      const uchar* querySeq, bool isMask ){
@@ -1041,7 +1036,7 @@ void scan(LastAligner& aligner, size_t queryNum,
     }
   }
 
-  if( !isCollatedAlignments() ) gappedAlns.sort();  // sort by score
+  if (!isCollatedAlignments()) gappedAlns.sort();  // sort by score
   alignFinish(aligner, gappedAlns, queryNum, matrices, frameSize, dis3);
 }
 
@@ -1130,20 +1125,17 @@ static void alignOneQuery(LastAligner &aligner, size_t finalCullingLimit,
 		     args.isQueryStrandMatrix ? revMatrices : fwdMatrices);
 }
 
-static void alignSomeQueries(size_t chunkNum,
-			     unsigned volume, unsigned volumeCount) {
+static void alignSomeQueries(size_t chunkNum, unsigned volume) {
   size_t numOfChunks = aligners.size();
   LastAligner &aligner = aligners[chunkNum];
   std::vector<AlignmentText> &textAlns = aligner.textAlns;
   size_t beg = firstSequenceInChunk(query, numOfChunks, chunkNum);
   size_t end = firstSequenceInChunk(query, numOfChunks, chunkNum + 1);
-  bool isMultiVolume = (volumeCount > 1);
+  bool isMultiVolume = (numOfVolumes > 1);
   bool isFirstVolume = (volume == 0);
-  bool isFinalVolume = (volume + 1 == volumeCount);
   bool isFirstThread = (chunkNum == 0);
-  bool isSort = isCollatedAlignments();
-  bool isSortPerQuery = (isSort && !isMultiVolume);
-  bool isPrintPerQuery = (isFirstThread && !(isSort && isMultiVolume));
+  bool isSortPerQuery = (isCollatedAlignments() && !isMultiVolume);
+  bool isPrintPerQuery = (isFirstThread && !isMultiVolume);
   size_t finalCullingLimit = args.cullingLimitForFinalAlignments ?
     args.cullingLimitForFinalAlignments : isMultiVolume;
   for (size_t i = beg; i < end; ++i) {
@@ -1152,27 +1144,27 @@ static void alignSomeQueries(size_t chunkNum,
     if (isSortPerQuery) sort(textAlns.begin() + oldNumOfAlns, textAlns.end());
     if (isPrintPerQuery) printAndClear(textAlns);
   }
-  if (isFinalVolume && isMultiVolume) {
+  if (isMultiVolume && volume + 1 == numOfVolumes) {
     cullFinalAlignments(textAlns, 0, args.cullingLimitForFinalAlignments);
-    if (isSort) sort(textAlns.begin(), textAlns.end());
-    if (isFirstThread) printAndClear(textAlns);
+    sort(textAlns.begin(), textAlns.end());
   }
 }
 
-static void scanOneVolume(unsigned volume, unsigned volumeCount) {
-#ifdef HAS_CXX_THREADS
-  size_t numOfChunks = aligners.size();
-  std::vector<std::thread> threads(numOfChunks - 1);
-  for (size_t i = 1; i < numOfChunks; ++i)
-    threads[i - 1] = std::thread(alignSomeQueries, i, volume, volumeCount);
-  // Exceptions from threads are not handled nicely, but I don't
-  // think it matters much.
-#endif
-  alignSomeQueries(0, volume, volumeCount);
+static void scanOneVolume(unsigned volume, unsigned numOfThreadsLeft) {
+  if (numOfThreadsLeft > 1) {
 #ifdef HAS_CXX_THREADS
-  for (size_t i = 1; i < numOfChunks; ++i)
-    threads[i - 1].join();
+    std::thread t(scanOneVolume, volume, numOfThreadsLeft - 1);
+    // Exceptions from threads are not handled nicely, but I don't
+    // think it matters much.
+    alignSomeQueries(numOfThreadsLeft - 1, volume);
+    t.join();
 #endif
+  } else {
+    alignSomeQueries(0, volume);
+  }
+  if (volume + 1 == numOfVolumes) {
+    printAndClear(aligners[numOfThreadsLeft - 1].textAlns);
+  }
 }
 
 void readIndex( const std::string& baseName, indexT seqCount ) {
@@ -1242,23 +1234,18 @@ void readVolume( unsigned volumeNumber ){
 }
 
 // Scan one batch of query sequences against all database volumes
-void scanAllVolumes( unsigned volumes, std::ostream& out ){
+void scanAllVolumes() {
   if( args.outputType == 0 ){
     matchCounts.clear();
     matchCounts.resize( query.finishedSequences() );
   }
 
-  if( volumes+1 == 0 ) volumes = 1;
-  bool isMultiVolume = (volumes > 1);
-
-  for( unsigned i = 0; i < volumes; ++i ){
-    if( text.unfinishedSize() == 0 || isMultiVolume ) readVolume( i );
-    scanOneVolume( i, volumes );
-    if( !isCollatedAlignments() ) printAndClearAll();
+  for (unsigned i = 0; i < numOfVolumes; ++i) {
+    if (text.unfinishedSize() == 0 || numOfVolumes > 1) readVolume(i);
+    scanOneVolume(i, aligners.size());
   }
 
-  if( args.outputType == 0 ) writeCounts( out );
-  printAndClearAll();
+  if (args.outputType == 0) writeCounts();
 }
 
 void writeHeader( countT refSequences, countT refLetters, std::ostream& out ){
@@ -1311,7 +1298,6 @@ void lastal( int argc, char** argv ){
   args.fromArgs( argc, argv );
   args.resetCumulativeOptions();  // because we will do fromArgs again
 
-  unsigned volumes = unsigned(-1);
   size_t refMinimizerWindow = 1;  // assume this value, if not specified
   size_t minSeedLimit = 0;
   countT refSequences = -1;
@@ -1319,10 +1305,9 @@ void lastal( int argc, char** argv ){
   countT refMaxSeqLen = -1;
   bool isKeepRefLowercase = true;
   int refTantanSetting = 0;
-  readOuterPrj( args.lastdbName + ".prj", volumes,
-		refMinimizerWindow, minSeedLimit,
-		isKeepRefLowercase, refTantanSetting,
-		refSequences, refLetters, refMaxSeqLen );
+  readOuterPrj(args.lastdbName + ".prj",
+	       refMinimizerWindow, minSeedLimit, isKeepRefLowercase,
+	       refTantanSetting, refSequences, refLetters, refMaxSeqLen);
   bool isDna = (alph.letters == alph.dna);
   bool isProtein = alph.isProtein();
 
@@ -1354,7 +1339,7 @@ void lastal( int argc, char** argv ){
 
   aligners.resize( decideNumberOfThreads( args.numOfThreads,
 					  args.programName, args.verbosity ) );
-  bool isMultiVolume = (volumes+1 > 0 && volumes > 1);
+  bool isMultiVolume = (numOfVolumes + 1 > 0 && numOfVolumes > 1);
   args.setDefaultsFromAlphabet( isDna, isProtein,
 				isKeepRefLowercase, refTantanSetting,
                                 isCaseSensitiveSeeds, isMultiVolume,
@@ -1412,11 +1397,12 @@ void lastal( int argc, char** argv ){
 
   queryAlph.tr(query.seqWriter(), query.seqWriter() + query.seqBeg(0));
 
-  if( volumes+1 == 0 ) readIndex( args.lastdbName, refSequences );
+  if (numOfVolumes + 1 == 0) {
+    readIndex(args.lastdbName, refSequences);
+    numOfVolumes = 1;
+  }
 
-  std::ostream& out = std::cout;
-  writeHeader( refSequences, refLetters, out );
-  out.precision(3);  // print non-integers more compactly
+  writeHeader(refSequences, refLetters, std::cout);
   countT queryBatchCount = 0;
   countT sequenceCount = 0;
   indexT maxSeqLen = args.batchSize;
@@ -1437,24 +1423,24 @@ void lastal( int argc, char** argv ){
 	++sequenceCount;
       } else {
         // this enables downstream parsers to read one batch at a time:
-        out << "# batch " << queryBatchCount++ << "\n";
-	scanAllVolumes( volumes, out );
+	std::cout << "# batch " << queryBatchCount++ << "\n";
+	scanAllVolumes();
 	query.reinitForAppending();
       }
     }
   }
 
   if( query.finishedSequences() > 0 ){
-    out << "# batch " << queryBatchCount << "\n";
-    scanAllVolumes( volumes, out );
+    std::cout << "# batch " << queryBatchCount << "\n";
+    scanAllVolumes();
   }
 
   countT numOfNormalLetters = 0;
   for (size_t i = 0; i < aligners.size(); ++i) {
     numOfNormalLetters += aligners[i].numOfNormalLetters;
   }
-  out << "# Query sequences=" << sequenceCount
-      << " normal letters=" << numOfNormalLetters << "\n";
+  std::cout << "# Query sequences=" << sequenceCount
+	    << " normal letters=" << numOfNormalLetters << "\n";
 }
 
 int main( int argc, char** argv )


=====================================
src/makefile
=====================================
@@ -143,7 +143,7 @@ ScoreMatrixData.hh: ../data/*.mat
 	../build/mat-inc.sh ../data/*.mat > $@
 
 VERSION1 = git describe --dirty
-VERSION2 = echo ' (HEAD -> main, tag: 1268) ' | sed -e 's/.*tag: *//' -e 's/[,) ].*//'
+VERSION2 = echo ' (HEAD -> main, tag: 1282) ' | sed -e 's/.*tag: *//' -e 's/[,) ].*//'
 
 VERSION = \"`test -e ../.git && $(VERSION1) || $(VERSION2)`\"
 


=====================================
src/mcf_simd.hh
=====================================
@@ -18,6 +18,7 @@ namespace mcf {
 
 typedef __m256i SimdInt;
 typedef __m256i SimdUint1;
+typedef __m256d SimdDbl;
 
 const int simdBytes = 32;
 
@@ -29,6 +30,10 @@ static inline SimdInt simdZero1() {
   return _mm256_setzero_si256();
 }
 
+static inline SimdDbl simdZeroDbl() {
+  return _mm256_setzero_pd();
+}
+
 static inline SimdInt simdOnes1() {
   return _mm256_set1_epi32(-1);
 }
@@ -41,6 +46,10 @@ static inline SimdInt simdLoad1(const void *p) {
   return _mm256_loadu_si256((const SimdInt *)p);
 }
 
+static inline SimdDbl simdLoadDbl(const double *p) {
+  return _mm256_loadu_pd(p);
+}
+
 static inline void simdStore(void *p, SimdInt x) {
   _mm256_storeu_si256((SimdInt *)p, x);
 }
@@ -49,6 +58,10 @@ static inline void simdStore1(void *p, SimdInt x) {
   _mm256_storeu_si256((SimdInt *)p, x);
 }
 
+static inline void simdStoreDbl(double *p, SimdDbl x) {
+  _mm256_storeu_pd(p, x);
+}
+
 static inline SimdInt simdOr1(SimdInt x, SimdInt y) {
   return _mm256_or_si256(x, y);
 }
@@ -58,6 +71,7 @@ static inline SimdInt simdBlend(SimdInt x, SimdInt y, SimdInt mask) {
 }
 
 const int simdLen = 8;
+const int simdDblLen = 4;
 
 static inline SimdInt simdSet(int i7, int i6, int i5, int i4,
 			      int i3, int i2, int i1, int i0) {
@@ -78,6 +92,10 @@ static inline SimdInt simdSet1(char jF, char jE, char jD, char jC,
 			 i7, i6, i5, i4, i3, i2, i1, i0);
 }
 
+static inline SimdDbl simdSetDbl(double i3, double i2, double i1, double i0) {
+  return _mm256_set_pd(i3, i2, i1, i0);
+}
+
 static inline SimdInt simdFill(int x) {
   return _mm256_set1_epi32(x);
 }
@@ -86,6 +104,10 @@ static inline SimdInt simdFill1(char x) {
   return _mm256_set1_epi8(x);
 }
 
+static inline SimdDbl simdFillDbl(double x) {
+  return _mm256_set1_pd(x);
+}
+
 static inline SimdInt simdGt(SimdInt x, SimdInt y) {
   return _mm256_cmpgt_epi32(x, y);
 }
@@ -106,6 +128,10 @@ static inline SimdInt simdAdds1(SimdInt x, SimdInt y) {
   return _mm256_adds_epu8(x, y);
 }
 
+static inline SimdDbl simdAddDbl(SimdDbl x, SimdDbl y) {
+  return _mm256_add_pd(x, y);
+}
+
 static inline SimdInt simdSub(SimdInt x, SimdInt y) {
   return _mm256_sub_epi32(x, y);
 }
@@ -114,6 +140,10 @@ static inline SimdInt simdSub1(SimdInt x, SimdInt y) {
   return _mm256_sub_epi8(x, y);
 }
 
+static inline SimdDbl simdMulDbl(SimdDbl x, SimdDbl y) {
+  return _mm256_mul_pd(x, y);
+}
+
 static inline SimdInt simdQuadruple1(SimdInt x) {
   return _mm256_slli_epi32(x, 2);
 }
@@ -142,6 +172,12 @@ static inline int simdHorizontalMin1(SimdInt x) {
   return _mm_extract_epi16(z, 0);
 }
 
+static inline double simdHorizontalAddDbl(SimdDbl x) {
+  __m128d z = _mm256_castpd256_pd128(x);
+  z = _mm_add_pd(z, _mm256_extractf128_pd(x, 1));
+  return _mm_cvtsd_f64(_mm_hadd_pd(z, z));
+}
+
 static inline SimdInt simdChoose1(SimdInt items, SimdInt choices) {
   return _mm256_shuffle_epi8(items, choices);
 }
@@ -150,6 +186,7 @@ static inline SimdInt simdChoose1(SimdInt items, SimdInt choices) {
 
 typedef __m128i SimdInt;
 typedef __m128i SimdUint1;
+typedef __m128d SimdDbl;
 
 const int simdBytes = 16;
 
@@ -161,6 +198,10 @@ static inline SimdInt simdZero1() {
   return _mm_setzero_si128();
 }
 
+static inline SimdDbl simdZeroDbl() {
+  return _mm_setzero_pd();
+}
+
 static inline SimdInt simdOnes1() {
   return _mm_set1_epi32(-1);
 }
@@ -173,6 +214,10 @@ static inline SimdInt simdLoad1(const void *p) {
   return _mm_loadu_si128((const SimdInt *)p);
 }
 
+static inline SimdDbl simdLoadDbl(const double *p) {
+  return _mm_loadu_pd(p);
+}
+
 static inline void simdStore(void *p, SimdInt x) {
   _mm_storeu_si128((SimdInt *)p, x);
 }
@@ -181,6 +226,10 @@ static inline void simdStore1(void *p, SimdInt x) {
   _mm_storeu_si128((SimdInt *)p, x);
 }
 
+static inline void simdStoreDbl(double *p, SimdDbl x) {
+  _mm_storeu_pd(p, x);
+}
+
 static inline SimdInt simdOr1(SimdInt x, SimdInt y) {
   return _mm_or_si128(x, y);
 }
@@ -190,6 +239,7 @@ static inline SimdInt simdBlend(SimdInt x, SimdInt y, SimdInt mask) {
 }
 
 const int simdLen = 4;
+const int simdDblLen = 2;
 
 static inline SimdInt simdSet(int i3, int i2, int i1, int i0) {
   return _mm_set_epi32(i3, i2, i1, i0);
@@ -203,6 +253,10 @@ static inline SimdInt simdSet1(char iF, char iE, char iD, char iC,
 		      i7, i6, i5, i4, i3, i2, i1, i0);
 }
 
+static inline SimdDbl simdSetDbl(double i1, double i0) {
+  return _mm_set_pd(i1, i0);
+}
+
 static inline SimdInt simdFill(int x) {
   return _mm_set1_epi32(x);
 }
@@ -211,6 +265,10 @@ static inline SimdInt simdFill1(char x) {
   return _mm_set1_epi8(x);
 }
 
+static inline SimdDbl simdFillDbl(double x) {
+  return _mm_set1_pd(x);
+}
+
 static inline SimdInt simdGt(SimdInt x, SimdInt y) {
   return _mm_cmpgt_epi32(x, y);
 }
@@ -231,6 +289,10 @@ static inline SimdInt simdAdds1(SimdInt x, SimdInt y) {
   return _mm_adds_epu8(x, y);
 }
 
+static inline SimdDbl simdAddDbl(SimdDbl x, SimdDbl y) {
+  return _mm_add_pd(x, y);
+}
+
 static inline SimdInt simdSub(SimdInt x, SimdInt y) {
   return _mm_sub_epi32(x, y);
 }
@@ -239,6 +301,10 @@ static inline SimdInt simdSub1(SimdInt x, SimdInt y) {
   return _mm_sub_epi8(x, y);
 }
 
+static inline SimdDbl simdMulDbl(SimdDbl x, SimdDbl y) {
+  return _mm_mul_pd(x, y);
+}
+
 static inline SimdInt simdQuadruple1(SimdInt x) {
   return _mm_slli_epi32(x, 2);
 }
@@ -263,6 +329,10 @@ static inline int simdHorizontalMin1(SimdInt x) {
   return _mm_extract_epi16(x, 0);
 }
 
+static inline double simdHorizontalAddDbl(SimdDbl x) {
+  return _mm_cvtsd_f64(_mm_hadd_pd(x, x));
+}
+
 static inline SimdInt simdChoose1(SimdInt items, SimdInt choices) {
   return _mm_shuffle_epi8(items, choices);  // SSSE3
 }
@@ -272,6 +342,7 @@ static inline SimdInt simdChoose1(SimdInt items, SimdInt choices) {
 typedef int32x4_t SimdInt;
 typedef uint32x4_t SimdUint;
 typedef uint8x16_t SimdUint1;
+typedef float64x2_t SimdDbl;
 
 const int simdBytes = 16;
 
@@ -283,6 +354,10 @@ static inline SimdUint1 simdZero1() {
   return vdupq_n_u8(0);
 }
 
+static inline SimdDbl simdZeroDbl() {
+  return vdupq_n_f64(0);
+}
+
 static inline SimdUint1 simdOnes1() {
   return vdupq_n_u8(-1);
 }
@@ -295,6 +370,10 @@ static inline SimdUint1 simdLoad1(const unsigned char *p) {
   return vld1q_u8(p);
 }
 
+static inline SimdDbl simdLoadDbl(const double *p) {
+  return vld1q_f64(p);
+}
+
 static inline void simdStore(int *p, SimdInt x) {
   vst1q_s32(p, x);
 }
@@ -303,6 +382,10 @@ static inline void simdStore1(unsigned char *p, SimdUint1 x) {
   vst1q_u8(p, x);
 }
 
+static inline void simdStoreDbl(double *p, SimdDbl x) {
+  vst1q_f64(p, x);
+}
+
 static inline SimdUint1 simdOr1(SimdUint1 x, SimdUint1 y) {
   return vorrq_u8(x, y);
 }
@@ -312,6 +395,7 @@ static inline SimdInt simdBlend(SimdInt x, SimdInt y, SimdUint mask) {
 }
 
 const int simdLen = 4;
+const int simdDblLen = 2;
 
 static inline SimdInt simdSet(unsigned i3, unsigned i2,
                               unsigned i1, unsigned i0) {
@@ -340,6 +424,10 @@ static inline SimdUint1 simdSet1(unsigned char iF, unsigned char iE,
   return vcombine_u8(vcreate_u8(lo), vcreate_u8(hi));
 }
 
+static inline SimdDbl simdSetDbl(double i1, double i0) {
+  return vcombine_f64(vdup_n_f64(i0), vdup_n_f64(i1));
+}
+
 static inline SimdInt simdFill(int x) {
   return vdupq_n_s32(x);
 }
@@ -348,6 +436,10 @@ static inline SimdUint1 simdFill1(unsigned char x) {
   return vdupq_n_u8(x);
 }
 
+static inline SimdDbl simdFillDbl(double x) {
+  return vdupq_n_f64(x);
+}
+
 static inline SimdUint simdGt(SimdInt x, SimdInt y) {
   return vcgtq_s32(x, y);
 }
@@ -368,6 +460,10 @@ static inline SimdUint1 simdAdds1(SimdUint1 x, SimdUint1 y) {
   return vqaddq_u8(x, y);
 }
 
+static inline SimdDbl simdAddDbl(SimdDbl x, SimdDbl y) {
+  return vaddq_f64(x, y);
+}
+
 static inline SimdInt simdSub(SimdInt x, SimdInt y) {
   return vsubq_s32(x, y);
 }
@@ -376,6 +472,10 @@ static inline SimdUint1 simdSub1(SimdUint1 x, SimdUint1 y) {
   return vsubq_u8(x, y);
 }
 
+static inline SimdDbl simdMulDbl(SimdDbl x, SimdDbl y) {
+  return vmulq_f64(x, y);
+}
+
 static inline SimdUint1 simdQuadruple1(SimdUint1 x) {
   return vshlq_n_u8(x, 2);
 }
@@ -396,6 +496,10 @@ static inline int simdHorizontalMin1(SimdUint1 x) {
   return vminvq_u8(x);
 }
 
+static inline double simdHorizontalAddDbl(SimdDbl x) {
+  return vaddvq_f64(x);
+}
+
 static inline SimdUint1 simdChoose1(SimdUint1 items, SimdUint1 choices) {
   return vqtbl1q_u8(items, choices);
 }
@@ -403,19 +507,29 @@ static inline SimdUint1 simdChoose1(SimdUint1 items, SimdUint1 choices) {
 #else
 
 typedef int SimdInt;
+typedef double SimdDbl;
 const int simdBytes = 1;
 const int simdLen = 1;
+const int simdDblLen = 1;
 static inline int simdZero() { return 0; }
+static inline double simdZeroDbl() { return 0; }
 static inline int simdSet(int x) { return x; }
+static inline double simdSetDbl(double x) { return x; }
 static inline int simdFill(int x) { return x; }
 static inline int simdLoad(const int *p) { return *p; }
+static inline double simdLoadDbl(const double *p) { return *p; }
 static inline void simdStore(int *p, int x) { *p = x; }
+static inline void simdStoreDbl(double *p, double x) { *p = x; }
+static inline double simdFillDbl(double x) { return x; }
 static inline int simdGt(int x, int y) { return x > y; }
 static inline int simdAdd(int x, int y) { return x + y; }
+static inline double simdAddDbl(double x, double y) { return x + y; }
 static inline int simdSub(int x, int y) { return x - y; }
+static inline double simdMulDbl(double x, double y) { return x * y; }
 static inline int simdMax(int x, int y) { return x > y ? x : y; }
 static inline int simdBlend(int x, int y, int mask) { return mask ? y : x; }
 static inline int simdHorizontalMax(int a) { return a; }
+static inline double simdHorizontalAddDbl(double x) { return x; }
 
 #endif
 


=====================================
src/split/last-split.cc
=====================================
@@ -27,23 +27,39 @@ static std::istream& openIn(const std::string& fileName, std::ifstream& ifs) {
 }
 
 // Does the string start with the prefix?
-static bool startsWith(const std::string& s, const char* prefix) {
-  const char* t = s.c_str();
+static bool startsWith(const char *s, const char *prefix) {
   for (;;) {
     if (*prefix == 0) return true;
-    if (*prefix != *t) return false;
-    ++t;
+    if (*prefix != *s) return false;
+    ++s;
     ++prefix;
   }
 }
 
 // Does the string have no non-space characters?
-static bool isSpace(const std::string& s) {
-  const char* t = s.c_str();
+static bool isBlankLine(const char *s) {
   for (;;) {
-    if (*t == 0) return true;
-    if (!std::isspace(*t)) return false;
-    ++t;
+    if (*s == 0) return true;
+    if (!std::isspace(*s)) return false;
+    ++s;
+  }
+}
+
+static bool isSpace(char c) {
+  return c > 0 && c <= ' ';
+}
+
+static bool isSameName(const char *sLine1, const char *sLine2) {
+  do { ++sLine1; } while (isSpace(*sLine1));
+  do { ++sLine2; } while (isSpace(*sLine2));
+  for (;;) {
+    if (*sLine1 > ' ') {
+      if (*sLine2 != *sLine1) return false;
+    } else {
+      return *sLine2 <= ' ';
+    }
+    ++sLine1;
+    ++sLine2;
   }
 }
 
@@ -321,12 +337,15 @@ void lastSplit(LastSplitOptions& opts) {
   double genomeSize = 0;
   std::vector<std::string> mafLines;  // lines of multiple MAF blocks
   std::vector<unsigned> mafEnds(1);  // which lines are in which MAF block
+  unsigned sLineCount = 0;
+  unsigned qNameLineNum = 0;
   bool isAlreadySplit = false;  // has the input already undergone last-split?
 
   for (unsigned i = 0; i < opts.inputFileNames.size(); ++i) {
     std::ifstream inFileStream;
     std::istream& input = openIn(opts.inputFileNames[i], inFileStream);
     while (getline(input, line)) {
+      const char *linePtr = line.c_str();
       if (state == -1) {  // we are reading the score matrix within the header
 	std::istringstream ls(line);
 	std::vector<int> row;
@@ -352,7 +371,7 @@ void lastSplit(LastSplitOptions& opts) {
 	if (word == "#" && !names.empty() && !ls && scoreMatrix.empty()) {
 	  colNames = names;
 	  state = -1;
-	} else if (startsWith(line, "#")) {
+	} else if (linePtr[0] == '#') {
 	  std::istringstream ls(line);
 	  while (ls >> word) {
 	    std::istringstream ws(word);
@@ -367,8 +386,8 @@ void lastSplit(LastSplitOptions& opts) {
 	    if (key == "letters") ws >> genomeSize;
 	  }
 	  // try to determine if last-split was already run (fragile):
-	  if (startsWith(line, "# m=")) isAlreadySplit = true;
-	} else if (!isSpace(line)) {
+	  if (startsWith(linePtr, "# m=")) isAlreadySplit = true;
+	} else if (!isBlankLine(linePtr)) {
 	  if (scoreMatrix.empty())
 	    err("I need a header with score parameters");
 	  if (gapExistenceCost < 0 || gapExtensionCost < 0 ||
@@ -410,18 +429,20 @@ void lastSplit(LastSplitOptions& opts) {
 	}
       }
       if (state == 1) {  // we are reading alignments
-	if (startsWith(line, "# batch") && !opts.isTopSeqQuery) {
+	if (isBlankLine(linePtr)) {
 	  addMaf(mafEnds, mafLines);
-	  doOneBatch(mafLines, mafEnds, sa, opts, isAlreadySplit);
-	  mafLines.clear();
-	  mafEnds.resize(1);
-	} else if (isSpace(line)) {
-	  addMaf(mafEnds, mafLines);
-	} else if (std::strchr(opts.no_split ? "asqpc" : "sqpc", line[0])) {
+	} else if (std::strchr(opts.no_split ? "asqpc" : "sqpc", linePtr[0])) {
+	  if (!opts.isTopSeqQuery && linePtr[0] == 's' && sLineCount++ % 2 &&
+	      !isSameName(mafLines[qNameLineNum].c_str(), linePtr)) {
+	    doOneBatch(mafLines, mafEnds, sa, opts, isAlreadySplit);
+	    mafLines.erase(mafLines.begin(), mafLines.begin()+mafEnds.back());
+	    mafEnds.resize(1);
+	    qNameLineNum = mafLines.size();
+	  }
 	  mafLines.push_back(line);
 	}
       }
-      if (startsWith(line, "#") && !startsWith(line, "# batch"))
+      if (linePtr[0] == '#' && !startsWith(linePtr, "# batch"))
 	std::cout << line << "\n";
     }
   }


=====================================
src/tantan.cc
=====================================
@@ -1,6 +1,7 @@
 // Copyright 2010 Martin C. Frith
 
 #include "tantan.hh"
+#include "mcf_simd.hh"
 
 #include <algorithm>  // fill, max
 #include <cassert>
@@ -14,6 +15,8 @@
 
 namespace tantan {
 
+using namespace mcf;
+
 void multiplyAll(std::vector<double> &v, double factor) {
   for (std::vector<double>::iterator i = v.begin(); i < v.end(); ++i)
     *i *= factor;
@@ -308,15 +311,37 @@ struct Tantan {
     }
 
     double b = backgroundProb;
-    double fromForeground = 0;
-    double *foregroundBeg = BEG(foregroundProbs);
+    const double *b2f = BEG(b2fProbs);
+    double *fp = BEG(foregroundProbs);
     const double *lrRow = likelihoodRatioMatrix[*seqPtr];
     int maxOffset = maxOffsetInTheSequence();
-
-    for (int i = 0; i < maxOffset; ++i) {
-      double f = foregroundBeg[i];
+    const uchar *sp = seqPtr;
+
+    SimdDbl bV = simdFillDbl(b);
+    SimdDbl tV = simdFillDbl(f2f0);
+    SimdDbl sV = simdZeroDbl();
+
+    int i = 0;
+    for (; i <= maxOffset - simdDblLen; i += simdDblLen) {
+      SimdDbl rV = simdSetDbl(
+#if defined __SSE4_1__ || defined __ARM_NEON
+#ifdef __AVX2__
+			      lrRow[sp[-i-4]],
+			      lrRow[sp[-i-3]],
+#endif
+			      lrRow[sp[-i-2]],
+#endif
+			      lrRow[sp[-i-1]]);
+      SimdDbl fV = simdLoadDbl(fp+i);
+      sV = simdAddDbl(sV, fV);
+      SimdDbl xV = simdMulDbl(bV, simdLoadDbl(b2f+i));
+      simdStoreDbl(fp+i, simdMulDbl(simdAddDbl(xV, simdMulDbl(fV, tV)), rV));
+    }
+    double fromForeground = simdHorizontalAddDbl(sV);
+    for (; i < maxOffset; ++i) {
+      double f = fp[i];
       fromForeground += f;
-      foregroundBeg[i] = (b * b2fProbs[i] + f * f2f0) * lrRow[seqPtr[-i-1]];
+      fp[i] = (b * b2f[i] + f * f2f0) * lrRow[sp[-i-1]];
     }
 
     backgroundProb = b * b2b + fromForeground * f2b;
@@ -330,15 +355,36 @@ struct Tantan {
     }
 
     double toBackground = f2b * backgroundProb;
-    double toForeground = 0;
-    double *foregroundBeg = BEG(foregroundProbs);
+    const double *b2f = BEG(b2fProbs);
+    double *fp = BEG(foregroundProbs);
     const double *lrRow = likelihoodRatioMatrix[*seqPtr];
     int maxOffset = maxOffsetInTheSequence();
-
-    for (int i = 0; i < maxOffset; ++i) {
-      double f = foregroundBeg[i] * lrRow[seqPtr[-i-1]];
-      toForeground += b2fProbs[i] * f;
-      foregroundBeg[i] = toBackground + f2f0 * f;
+    const uchar *sp = seqPtr;
+
+    SimdDbl bV = simdFillDbl(toBackground);
+    SimdDbl tV = simdFillDbl(f2f0);
+    SimdDbl sV = simdZeroDbl();
+
+    int i = 0;
+    for (; i <= maxOffset - simdDblLen; i += simdDblLen) {
+      SimdDbl rV = simdSetDbl(
+#if defined __SSE4_1__ || defined __ARM_NEON
+#ifdef __AVX2__
+			      lrRow[sp[-i-4]],
+			      lrRow[sp[-i-3]],
+#endif
+			      lrRow[sp[-i-2]],
+#endif
+			      lrRow[sp[-i-1]]);
+      SimdDbl fV = simdMulDbl(simdLoadDbl(fp+i), rV);
+      sV = simdAddDbl(sV, simdMulDbl(simdLoadDbl(b2f+i), fV));
+      simdStoreDbl(fp+i, simdAddDbl(bV, simdMulDbl(tV, fV)));
+    }
+    double toForeground = simdHorizontalAddDbl(sV);
+    for (; i < maxOffset; ++i) {
+      double f = fp[i] * lrRow[sp[-i-1]];
+      toForeground += b2f[i] * f;
+      fp[i] = toBackground + f2f0 * f;
     }
 
     backgroundProb = b2b * backgroundProb + toForeground;


=====================================
test/SRR359290-1k.maf
=====================================
The diff for this file was not included because it is too large.

=====================================
test/bs100.maf
=====================================
@@ -10,6 +10,54 @@ s 1       0 87 +    87 TTATATAGAGGAGATAAGTCGTAATATGGTAAGTGTATTGGAAAGTGTATTTGGATG
 q 1                    BB at 1;?@B at BBBC@BBBCACBBA>B??=?A at B?B?B=7 at 4BB;B>@@@A?>=B at A@>@7>BA>BA6BA?@@B@@?B?AABBAB=<>B
 p                      $(-16:@FLRW]bddfgfccedffdedffdfgfdedeccfggggfdedecbcfgfdfgfddfggfdedfecba`]WTOKF@:40+'$
 
+a score=364 mismap=8.24e-08
+s chrM 12300 84 + 16571 GGCCCCAAAAATTTTGGTGCAACTCCAAATAAAAGTAATAACCATGCACACTACTATAACCACCCTAACCCTGACTTCCCTAAT
+s 2        3 84 +    87 GGTTTTAAAAATTTTGGTGTAATTTTAAATAAAAGTAATAATTATGTATATTATTATAATTATTCTAACTTTGATTTTTTTAAT
+q 2                     B at B>BB?B;=B@@=BBABA@<A5??'8@@>=3 at 47>@+AA9?;@@>@9?@(A>-9 at A=(=>?1;=7=96 at .4<(49=.2.1(==
+p                       &+.158>DJPUWZ]_dededfecbbbegfdfgggfdffdffccedededdbcdccedfeccdbbbbb_ZWTQNID@=:741.*$
+
+a score=258 mismap=8.71e-08
+s chrM 12453 69 + 16571 TGTCGCATCCACCTTTATTATCAGTCTCTTCCCCACAACAATATTCATGTGCCTAGACCAAGAAGTTAT
+s 3       18 69 -    87 AATTGCATCCACCATTACTATCAATCTCTTCCCCCCAACAATATTCATATACCCAAACCAAAAAATGAT
+q 3                     #52%7;3,5>BA4)9?+%+@<11>AA;*4:)383)CB=)+7*;8BB4>@@966)(BBAA<?@='?()>B
+p                       #%()-26<BGJOQQVZ[Z_adeccfggfggffd`Z^`acbbddfgfcc_\WTOIHHFEA<852/,)%%#
+
+a score=408 mismap=8.02e-08
+s chrM 5737 87 + 16571 CGCCGCCGGGAAAAAAGGCGGGAGAAGCCCCGGCAGGTTTGAAGCTGCTTCTTCGAATTTGCAATTCAATATGAAAATCACCTCGGA
+s 4       0 87 +    87 CGTCGTCGGGAAAAAAGGCGGGAGAAGTTTCGGTAGGTTTGAAGCTGTTTTTTCGAATTTGTAATTTAATATGAAAATTATTTCGGA
+q 4                    ?BCC@?CCBCC at CBCACCCCCCCC@ABCCAABBCCBC at A=B at CBCACCC@BCB@?CCACBCBBBBB>C>@>BBB?BBB:=<B at B>BB
+p                      $),059<BHNTZ_dfggfdgggggggfcbbcffdfgfcccfggfccdbbaaabcfffcbcedffcbcffdecda\VPJFC>:741,&
+
+a score=280 mismap=9.01e-08
+s chrM 6086 87 + 16571 TGCATTTGTAATAATCTTCTTCATAGTAATACCCATCATAATCGGAGGCTTTGGCAACTGACTAGTTCCCCTAATAATCGGTGCCCC
+s 5       0 87 -    87 TGGGTTCTTCTTAGCAACTGTTATAATAATACCCATCATAATCGAAAACTTTAACAACTAACAAATTCCCCTAATAATCGATACCCC
+q 5                    ######################A8</<B::6>?@>@@?5>7/?B68@/B.A4BC<BB;B3B?BB>C5BB=>@>CBBABAABCCC?CB
+p                      !!!"""####$$$%%&''(*+-169=BEINRX]abefdeccfecbbbcfggfccdcbb_ZVSNQU[`dfggea_]XURMGC@;72,&
+
+a score=368 mismap=8.65e-08
+s chrM 11106 87 + 16571 TCATATTTTATATCTTCTTCGAAACCACACTTATCCCCACCTTGGCTATCATCACCCGATGAGGCAACCAGCCAGAACGCCTGAACG
+s 6        0 87 +    87 TTATATTTTATATTTTTTTCGAAATTATATTTATTTTTATTTTGGTTATTATTATTCGATGAGGTAATTAGTTAGAATGTCTGAACG
+q 6                     BBB?C<@BA;5AB:AB at B?B<AB>A<AB:@AB9;=@AA43A>@;,5=@@A@/@>@=:?:@>><>-=6 at 746?=05?*?<<==9775=
+p                       $(-058;?BGKORTVWYZ\^begfccedecbcdbbbbbdbbbcefccecceccdcbcffdfggfdffccfecba]XRMID@=:50*&
+
+a score=340 mismap=9.13e-08
+s chrM 143 85 + 16571 CCTCATTCTATTATTTATCGCACCTACGTTCAATATTACAGGCGAACATACCTACTAAAGTGTGTTAATTAATTAATGCTTGTAG
+s 7      2 85 -    87 CCTAATTCTCTTATCTATCGCTCCTACATTCAATATTACAAGCAAACATACCTACTAAAATATATTAATTAATTAATACTTATAA
+q 7                   ####;672++4<=.060?>+/+22(9<;9;;94@@1=5=B:)?BA at BBBA98>4+:?<BB?@>:,6?B;BBCBCBBCBCC>@CCB
+p                     #%')-38;<=ADEFFKOTXYYW]accedfgfccedffdecbbdbbcededfgfdffcbbcededffcced`^\WQNJEB<71-($
+
+a score=368 mismap=8.66e-08
+s chrM 5006 87 + 16571 AGCATACTCCTCAATTACCCACATAGGATGAATAATAGCAGTTCTACCGTACAACCCTAACATAACCATTCTTAATTTAACTATTTA
+s 8       0 87 +    87 AGTATATTCTTTAATTATTTATATAGGATGAATAATAGTAGTTTTATCGTATAATTTTAATATATCTATTTTTAATTTAATTATTTA
+q 8                    B at C=BBCCBC?CCCC>BCC<B;@B@@BA2>C:<ABCC<A at A=@BB- at B?:B<B7>?AB:A@@@>(:AB@=B;B at BB>@<@>@A at ABB
+p                      &+/48=@CFILPUY\^`aabddedfggfdfgfdffdffdffcbbcdccededfecbbcefddb`\^`ba``^]ZVPMIFA<84/,(%
+
+a score=375 mismap=9.32e-08
+s chrM 7901 87 + 16571 TGAACCTACGAGTACACCGACTACGGCGGACTAATCTTCAACTCCTACATACTTCCCCCATTATTCCTAGAACCAGGCGACCTGCGA
+s 9       0 87 +    87 TGAATTTACGAGTATATCGATTACGGTGGATTAATTTTTAATTTTTATATATTTTTTTTATTATTTTTAGAATTAGGTGATTTGCGA
+q 9                    B=9=B+BBBB>BBB=:BB@;<@@B:<A4@?A:?=BABAAAA>>@<A>@AA>@@@?@@6878:6588?<:>=?9<::=<51<99>;48
+p                      $*058;?DHNTY\`bdccffccedffdfgfccfecbbbceecbbbcededebbaaaabbdccdbbbbcfggea_]YSMID>:74/+&
+
 a score=207 mismap=1.15e-07
 s chrM 3290 60 + 16571 ttcctcttcttaacaacaTACCCATGGCCAACCTCCTACTCCTCATTGTACCCATTCTAA
 s 10      0 60 +    87 TTTTTTTTTTTAATCATATATTTATGGTTAATTTTTTATTTTTTATTGTATTTATTCTAA
@@ -64,12 +112,6 @@ s 19       0 87 +    87 ACAATATATATTTTAATAAATAATGTTTAATTGGTAATTATTATTAATTAACGTTT
 q 19                    =@9 at B>C<';>)*8?'>67 at 1;(?<;A<@:@@(9A@@)(38>9+@<8?6 at 4:=?B976/5<@>@@/@<?1:0'9178--9=:426>@
 p                       &*/58=AEGJKKNRW\_cffdfedecbcec_]Z_addbbdcceccffccffdecbcedffccededfggeb`_^\ZURMIFA;72,&
 
-a score=364 mismap=8.24e-08
-s chrM 12300 84 + 16571 GGCCCCAAAAATTTTGGTGCAACTCCAAATAAAAGTAATAACCATGCACACTACTATAACCACCCTAACCCTGACTTCCCTAAT
-s 2        3 84 +    87 GGTTTTAAAAATTTTGGTGTAATTTTAAATAAAAGTAATAATTATGTATATTATTATAATTATTCTAACTTTGATTTTTTTAAT
-q 2                     B at B>BB?B;=B@@=BBABA@<A5??'8@@>=3 at 47>@+AA9?;@@>@9?@(A>-9 at A=(=>?1;=7=96 at .4<(49=.2.1(==
-p                       &+.158>DJPUWZ]_dededfecbbbegfdfgggfdffdffccedededdbcdccedfeccdbbbbb_ZWTQNID@=:741.*$
-
 a score=329 mismap=1.29e-07
 s chrM 5859 87 + 16571 TTACAGTCCAATGCTTCACTCAGCCATTTTACCTCACCCCCACTGATGTTCGCCGACCGTTGACTATTCTCTACAAACCACAAAGAC
 s 20      0 87 +    87 TTATAGTTTAATGTTTTATTTAGTTATTTTATTTTATTTTTACTGATGTTCGTCGATCGTTGATTATTATTTATAGATTATAAAGAT
@@ -130,12 +172,6 @@ s 29      0 81 +    87 AGTAATAAAATTAGGAATAGTTTTTTTTTATTTTTGAGTTTTAGAGGTTATTTAAGG
 q 29                   ?6 at C@5<(<CCB.CA?4BB=4?AA;=@CA=A<%=5>BB=15=B=>>6<B??):=>-<BA>(<@2;+%,?0805.'A7>(=;
 p                      &+/5:>DIOTW[`dfgfdfecbaaaaabbdba`abeffcbbcfgggfccdbacfggfddb```_^\ZWROIFB?:62.+'$
 
-a score=258 mismap=8.71e-08
-s chrM 12453 69 + 16571 TGTCGCATCCACCTTTATTATCAGTCTCTTCCCCACAACAATATTCATGTGCCTAGACCAAGAAGTTAT
-s 3       18 69 -    87 AATTGCATCCACCATTACTATCAATCTCTTCCCCCCAACAATATTCATATACCCAAACCAAAAAATGAT
-q 3                     #52%7;3,5>BA4)9?+%+@<11>AA;*4:)383)CB=)+7*;8BB4>@@966)(BBAA<?@='?()>B
-p                       #%()-26<BGJOQQVZ[Z_adeccfggfggffd`Z^`acbbddfgfcc_\WTOIHHFEA<852/,)%%#
-
 a score=349 mismap=1.08e-07
 s chrM 6055 87 + 16571 ACCACATCTACAACGTTATCGTCACAGCCCATGCATTTGTAATAATCTTCTTCATAGTAATACCCATCATAATCGGAGGCTTTGGCA
 s 30      0 87 +    87 ATTATATTTATAACGTTATCGTTATAGTTTATGTATTTGTAATAATTTTTTTTATAGTAATATTTATTAGAATCGGAGGTTTTGGTC
@@ -196,12 +232,6 @@ s 39      0 87 +    87 AATGGTATTGAATTTACGAGTATATCGATTACGGTGGATTAATTTTTAATTTTTATA
 q 39                   B?=BBBBB?B at BABBB=@BABBAAA at BBAA>;0===8>=;7.5:<;@@@A@<8A>351:@A2(;>@9-;73:7;<34.962+:=;.6
 p                      &+/5:>CFJPVZ]_addfgfdedeccffccedffdfgfccfecbbbceecbbbcededebaa`aaabdba`\YVSPMHB=7321.)&
 
-a score=408 mismap=8.02e-08
-s chrM 5737 87 + 16571 CGCCGCCGGGAAAAAAGGCGGGAGAAGCCCCGGCAGGTTTGAAGCTGCTTCTTCGAATTTGCAATTCAATATGAAAATCACCTCGGA
-s 4       0 87 +    87 CGTCGTCGGGAAAAAAGGCGGGAGAAGTTTCGGTAGGTTTGAAGCTGTTTTTTCGAATTTGTAATTTAATATGAAAATTATTTCGGA
-q 4                    ?BCC@?CCBCC at CBCACCCCCCCC@ABCCAABBCCBC at A=B at CBCACCC@BCB@?CCACBCBBBBB>C>@>BBB?BBB:=<B at B>BB
-p                      $),059<BHNTZ_dfggfdgggggggfcbbcffdfgfcccfggfccdbbaaabcfffcbcedffcbcffdecda\VPJFC>:741,&
-
 a score=354 mismap=8.27e-08
 s chrM 5186 87 + 16571 ACTAACACCCTTAATTCCATCCACCCTCCTCTCCCTAGGAGGCCTGCCCCCGCTAACCGGCTTTTTGCCCAAATGGGCCATTATCGA
 s 40      0 87 +    87 ATTAATATTTTTAATTTTATTTATTTTTTTTTTTTTAGGAGGTTTGTTTTCGTTAATCGGTTTTTTGTCTAAATGGGTTATTATTGA
@@ -262,12 +292,6 @@ s 49       0 87 +    87 TATGATTTATTTTATATTTTAATTTATGAGATTTATAATAAATAGTTTTTTTAAAT
 q 49                    BCBCBBCCCBB?BBB6CBCCB@@@BAA at C?BB at B@B@?BBC?B at B@@BBB>BBA5 at BB<>B=?@A4B>?>?BA==>?B?=@B:<=@B
 p                       $*-38;>BGJMPTY\``abceecbcedfggfcbcedffdfgfdfecbaaabcfgfdeccfecccfffcba``\YVSPKGD?<83.($
 
-a score=280 mismap=9.01e-08
-s chrM 6086 87 + 16571 TGCATTTGTAATAATCTTCTTCATAGTAATACCCATCATAATCGGAGGCTTTGGCAACTGACTAGTTCCCCTAATAATCGGTGCCCC
-s 5       0 87 -    87 TGGGTTCTTCTTAGCAACTGTTATAATAATACCCATCATAATCGAAAACTTTAACAACTAACAAATTCCCCTAATAATCGATACCCC
-q 5                    ######################A8</<B::6>?@>@@?5>7/?B68@/B.A4BC<BB;B3B?BB>C5BB=>@>CBBABAABCCC?CB
-p                      !!!"""####$$$%%&''(*+-169=BEINRX]abefdeccfecbbbcfggfccdcbb_ZVSNQU[`dfggea_]XURMGC@;72,&
-
 a score=151 mismap=1.69e-07
 s chrM 6025 36 + 16571 TGGGCCAGCCAGGCAACCTTCTAGGTAACGACCACA
 s 50      0 36 +    87 TGGGTTAGTTAGGTAATTTTTTAGGTAACGATGATA
@@ -328,12 +352,6 @@ s 59      5 82 -    87 CTACGCCTAATACTAACATTCTATAAATATAATTTAGCTATTTCTATATATCTCCAT
 q 59                   ###?1@<AA>/@)+;6=3,=4)7:<9A@@>:@A>B?)420>99A5 at AB<?<>6=?B?BBAB:*0A?B-B9BA>B:>@8=ACB
 p                      "#&+/5;@CGLPUY\^`^\WRVZ^_abdceccfgfcbefdfgggfdededfggggdfgfdfebcdbbaaaa^YSMID>82,&
 
-a score=368 mismap=8.65e-08
-s chrM 11106 87 + 16571 TCATATTTTATATCTTCTTCGAAACCACACTTATCCCCACCTTGGCTATCATCACCCGATGAGGCAACCAGCCAGAACGCCTGAACG
-s 6        0 87 +    87 TTATATTTTATATTTTTTTCGAAATTATATTTATTTTTATTTTGGTTATTATTATTCGATGAGGTAATTAGTTAGAATGTCTGAACG
-q 6                     BBB?C<@BA;5AB:AB at B?B<AB>A<AB:@AB9;=@AA43A>@;,5=@@A@/@>@=:?:@>><>-=6 at 746?=05?*?<<==9775=
-p                       $(-058;?BGKORTVWYZ\^begfccedecbcdbbbbbdbbbcefccecceccdcbcffdfggfdffccfecba]XRMID@=:50*&
-
 a score=383 mismap=8.82e-08
 s chrM 12614 87 + 16571 CATTGTTCGTTACATGGTCCATCATAGAATTCTCACTGTGATATATAAACTCAGACCCAAACATTAATCAGTTCTTCAAATATCTAC
 s 60       0 87 -    87 CATTATTCGTTACATAATCCATCATAAAATTCTCACTATAATATATAAACTCAAACCCAAACATTAATCAATTCTTCAAATCTCTAC
@@ -382,12 +400,6 @@ s 69      0 87 +    87 GATTAATGGGATTTAAATTTAAAAATATTTAGTTAATAGTTAAGTATTTTAATTAAT
 q 69                   B>BB;1ABB?@9ABAABAB at 4)?@A>(@>>*>>53A<>0@>@BAA@:7;)@<?AA at AB2',(@((0A@>?=7:>?A at B:@'?.<###
 p                      &+/28=AGMSWZ]`dfeba`^Z_cecdbbcdeccffdffccfgfdecbbbeeccffccedbaaacda_][VROLIFC@=852-*'$"
 
-a score=340 mismap=9.13e-08
-s chrM 143 85 + 16571 CCTCATTCTATTATTTATCGCACCTACGTTCAATATTACAGGCGAACATACCTACTAAAGTGTGTTAATTAATTAATGCTTGTAG
-s 7      2 85 -    87 CCTAATTCTCTTATCTATCGCTCCTACATTCAATATTACAAGCAAACATACCTACTAAAATATATTAATTAATTAATACTTATAA
-q 7                   ####;672++4<=.060?>+/+22(9<;9;;94@@1=5=B:)?BA at BBBA98>4+:?<BB?@>:,6?B;BBCBCBBCBCC>@CCB
-p                     #%')-38;<=ADEFFKOTXYYW]accedfgfccedffdecbbdbbcededfgfdffcbbcededffcced`^\WQNJEB<71-($
-
 a score=340 mismap=8.71e-08
 s chrM 7675 84 + 16571 CATTTTCCTTATCTGCTTCCTAGTCCTGTATGCCCTTTTCCTAACACTCACAACAAAACTAACTAATACTAACATCTCAGACGC
 s 70      0 84 +    87 TATTTTTTTTATTTGTTTTTTAGTTTTGTATGTTTTTTTTTTAATATTTATAATAAAATTAATTAATATTAATATTTTAGACGT
@@ -448,12 +460,6 @@ s 79       0 87 +    87 ATTTGAATTAATATAATTATTTATAGTTTAATTATTAGTATTATTTTTTTATTATT
 q 79                    BCC?BCBBBACB?B:B>4B:9B>A<@@B at BBBA=B?9 at AB@AA9BA9?@@A14=BB@>B>@B8@@?@:4:0=<A>@<A>><85:;22
 p                       %),06;@DHMRV[^bdccdcbcedffcbcefcceccffdeccdbbaaabbdccdbbaabcefccfgfcceeba\WSNHDA>952/*$
 
-a score=368 mismap=8.66e-08
-s chrM 5006 87 + 16571 AGCATACTCCTCAATTACCCACATAGGATGAATAATAGCAGTTCTACCGTACAACCCTAACATAACCATTCTTAATTTAACTATTTA
-s 8       0 87 +    87 AGTATATTCTTTAATTATTTATATAGGATGAATAATAGTAGTTTTATCGTATAATTTTAATATATCTATTTTTAATTTAATTATTTA
-q 8                    B at C=BBCCBC?CCCC>BCC<B;@B@@BA2>C:<ABCC<A at A=@BB- at B?:B<B7>?AB:A@@@>(:AB@=B;B at BB>@<@>@A at ABB
-p                      &+/48=@CFILPUY\^`aabddedfggfdfgfdffdffdffcbbcdccededfecbbcefddb`\^`ba``^]ZVPMIFA<84/,(%
-
 a score=359 mismap=9.42e-08
 s chrM 8356 87 + 16571 TTACAGTGAAATGCCCCAACTAAATACTACCGTATGGCCCACCATAATTACCCCCATACTCCTTACACTATTCCTCATCACCCAACT
 s 81      0 87 +    87 TTATAGTGAAATGTTTTAATTAAATATTATCGTATGGTTTATTATAATTATTTTTATATTTTTTATATTATTTTTTATTATTTAATT
@@ -502,12 +508,6 @@ s 89      0 87 +    87 TTTTATATTTTTTATATTATTTTTTATTATTTAATTAAAAATATTAAATATAAATTA
 q 89                   9B:):<CCCAC>C:CCAC==B?CCCCCCCCB at BBCBCB<CBCCCB9?<B6BBA;B>CCC at C@<74C=CA/(;8@;-<3B at C=)%51C
 p                      $'*-26:=@CGJNSV[]`bbaaabbdccdcbcefccfgggfdeccfgfdfdfgfcceccdcbcdbbaaa```a]ZXSMGA=:730,&
 
-a score=375 mismap=9.32e-08
-s chrM 7901 87 + 16571 TGAACCTACGAGTACACCGACTACGGCGGACTAATCTTCAACTCCTACATACTTCCCCCATTATTCCTAGAACCAGGCGACCTGCGA
-s 9       0 87 +    87 TGAATTTACGAGTATATCGATTACGGTGGATTAATTTTTAATTTTTATATATTTTTTTTATTATTTTTAGAATTAGGTGATTTGCGA
-q 9                    B=9=B+BBBB>BBB=:BB@;<@@B:<A4@?A:?=BABAAAA>>@<A>@AA>@@@?@@6878:6588?<:>=?9<::=<51<99>;48
-p                      $*058;?DHNTY\`bdccffccedffdfgfccfecbbbceecbbbcededebbaaaabbdccdbbbbcfggea_]YSMID>:74/+&
-
 a score=278 mismap=1.03e-07
 s chrM 1598 66 + 16571 CGAACCAGAGTGTAGCTTAACACAAAGCACCCAACTTACACTTAGGAGATTTCAACTTAACTTGAC
 s 90     21 66 -    87 CTAACCAAAATATAACTTAACACAAAACACCCAACTTACACTTAAAAAATTTCAACTTAACTTAAC


=====================================
test/last-split-test.out
=====================================
@@ -4816,13 +4816,6 @@ s SRR359290.9962       0 75 +        75 AGCGAGGCTCCCCAATGGTGCTTTTGGCTTTAGTGTACGA
 q SRR359290.9962                        AABA>=C?DBEEEEEEEEEEBEEEBAEBEE?EDEAEBDEDEEB?D=?B:@:ABA@@C@@A??B@@??BEB:??B2
 p                                       &,28>DJPV\bhntz~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ztnhb\VPJD>82,&
 
-# CPU time: 2.52 seconds
-a score=334 mismap=1e-10
-s chr12           110785491 75 + 133851895 AATGTTTGCGCATGTTCGAGATGAGTCTCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
-s SRR359290.10000         0 75 -        75 AATGTTTGCGCATGTTCGAGATGAGTCGCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
-q SRR359290.10000                          ############################B13BB:B1;*6 at CCC=CACCC3A9A>12:;+9A;>2+A6AA=C?:5?
-p                                          "$&(*,.02468:<>@BDFHJLNPRTVX^djpv|~~~~~~~~~~~~~~~~~~~~~~~~~~ztnhb\VPJD>82,&
-
 a score=450 mismap=1e-10
 s chr17          72765010 75 + 81195210 CCTCCTCCTCACTGAGTGCCTCATCCTACCGGGTGTCCCTTTGCCACCCTGCCTGGGACATCGCTGGAACCTGCA
 s SRR359290.9964        0 75 -       75 CCTCCTCCTCACTGAGTGCCTCATCCTACCGGGTGTCCCTTTGCCACCCTGCCTGGGACATCGCTGGAACCTGCA
@@ -5027,6 +5020,12 @@ s SRR359290.9998        0 75 -        75 GAAATGGAGATGGTGGCAGTACCAGTGAGACCCCTCAGC
 q SRR359290.9998                         ECCE?DECDD==A=CDEBAECDE?FFG?EEFFGFFGEEDF=EEGGGGEGFGFGGGGGGGGGGGFGFGGFGAFGDG
 p                                        &,15788888888888888888888888888888888888888888888888888888888888888888751,&
 
+a score=334 mismap=1e-10
+s chr12           110785491 75 + 133851895 AATGTTTGCGCATGTTCGAGATGAGTCTCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
+s SRR359290.10000         0 75 -        75 AATGTTTGCGCATGTTCGAGATGAGTCGCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
+q SRR359290.10000                          ############################B13BB:B1;*6 at CCC=CACCC3A9A>12:;+9A;>2+A6AA=C?:5?
+p                                          "$&(*,.02468:<>@BDFHJLNPRTVX^djpv|~~~~~~~~~~~~~~~~~~~~~~~~~~ztnhb\VPJD>82,&
+
 # LAST version 356
 #
 # a=21 b=9 A=21 B=9 c=100000 F=-1 e=120 d=115 x=119 y=44 z=119
@@ -8938,12 +8937,6 @@ s chr7           2702653 75 + 159138663 AGCGAGGCTCCCCAATGGTGCTTTTGGCTTTAGTGTACGA
 s SRR359290.9962       0 75 +        75 AGCGAGGCTCCCCAATGGTGCTTTTGGCTTTAGTGTACGATGTTTGCTGTGCTTCCCGCCGTGGAGGGCAGAGCC
 q SRR359290.9962                        AABA>=C?DBEEEEEEEEEEBEEEBAEBEE?EDEAEBDEDEEB?D=?B:@:ABA@@C@@A??B@@??BEB:??B2
 
-# CPU time: 2.52 seconds
-a score=334 mismap=1e-10
-s chr12           110785491 75 + 133851895 AATGTTTGCGCATGTTCGAGATGAGTCTCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
-s SRR359290.10000         0 75 -        75 AATGTTTGCGCATGTTCGAGATGAGTCGCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
-q SRR359290.10000                          ############################B13BB:B1;*6 at CCC=CACCC3A9A>12:;+9A;>2+A6AA=C?:5?
-
 a score=450 mismap=1e-10
 s chr17          72765010 75 + 81195210 CCTCCTCCTCACTGAGTGCCTCATCCTACCGGGTGTCCCTTTGCCACCCTGCCTGGGACATCGCTGGAACCTGCA
 s SRR359290.9964        0 75 -       75 CCTCCTCCTCACTGAGTGCCTCATCCTACCGGGTGTCCCTTTGCCACCCTGCCTGGGACATCGCTGGAACCTGCA
@@ -9099,6 +9092,11 @@ s chr15          79183831 75 + 102531392 GAAATGGAGATGGTGGCAGTACCAGTGAGACCCCTCAGC
 s SRR359290.9998        0 75 -        75 GAAATGGAGATGGTGGCAGTACCAGTGAGACCCCTCAGCCTCCTCGGAAGAAAAGGGCCCGGGTAGATCCTACTG
 q SRR359290.9998                         ECCE?DECDD==A=CDEBAECDE?FFG?EEFFGFFGEEDF=EEGGGGEGFGFGGGGGGGGGGGFGFGGFGAFGDG
 
+a score=334 mismap=1e-10
+s chr12           110785491 75 + 133851895 AATGTTTGCGCATGTTCGAGATGAGTCTCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
+s SRR359290.10000         0 75 -        75 AATGTTTGCGCATGTTCGAGATGAGTCGCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
+q SRR359290.10000                          ############################B13BB:B1;*6 at CCC=CACCC3A9A>12:;+9A;>2+A6AA=C?:5?
+
 # LAST version 356
 #
 # a=21 b=9 A=21 B=9 c=100000 F=-1 e=120 d=115 x=119 y=44 z=119
@@ -13059,12 +13057,6 @@ s chr7           2702653 75 + 159138663 AGCGAGGCTCCCCAATGGTGCTTTTGGCTTTAGTGTACGA
 s SRR359290.9962       0 75 +        75 AGCGAGGCTCCCCAATGGTGCTTTTGGCTTTAGTGTACGATGTTTGCTGTGCTTCCCGCCGTGGAGGGCAGAGCC
 q SRR359290.9962                        AABA>=C?DBEEEEEEEEEEBEEEBAEBEE?EDEAEBDEDEEB?D=?B:@:ABA@@C@@A??B@@??BEB:??B2
 
-# CPU time: 2.52 seconds
-a score=334 mismap=1e-10
-s chr12           110785491 75 + 133851895 AATGTTTGCGCATGTTCGAGATGAGTCTCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
-s SRR359290.10000         0 75 -        75 AATGTTTGCGCATGTTCGAGATGAGTCGCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
-q SRR359290.10000                          ############################B13BB:B1;*6 at CCC=CACCC3A9A>12:;+9A;>2+A6AA=C?:5?
-
 a score=450 mismap=1e-10
 s chr17          72765010 75 + 81195210 CCTCCTCCTCACTGAGTGCCTCATCCTACCGGGTGTCCCTTTGCCACCCTGCCTGGGACATCGCTGGAACCTGCA
 s SRR359290.9964        0 75 -       75 CCTCCTCCTCACTGAGTGCCTCATCCTACCGGGTGTCCCTTTGCCACCCTGCCTGGGACATCGCTGGAACCTGCA
@@ -13230,6 +13222,11 @@ s chr15          79183831 75 + 102531392 GAAATGGAGATGGTGGCAGTACCAGTGAGACCCCTCAGC
 s SRR359290.9998        0 75 -        75 GAAATGGAGATGGTGGCAGTACCAGTGAGACCCCTCAGCCTCCTCGGAAGAAAAGGGCCCGGGTAGATCCTACTG
 q SRR359290.9998                         ECCE?DECDD==A=CDEBAECDE?FFG?EEFFGFFGEEDF=EEGGGGEGFGFGGGGGGGGGGGFGFGGFGAFGDG
 
+a score=334 mismap=1e-10
+s chr12           110785491 75 + 133851895 AATGTTTGCGCATGTTCGAGATGAGTCTCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
+s SRR359290.10000         0 75 -        75 AATGTTTGCGCATGTTCGAGATGAGTCGCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
+q SRR359290.10000                          ############################B13BB:B1;*6 at CCC=CACCC3A9A>12:;+9A;>2+A6AA=C?:5?
+
 # LAST version 356
 #
 # a=21 b=9 A=21 B=9 c=100000 F=-1 e=120 d=115 x=119 y=44 z=119
@@ -17217,12 +17214,6 @@ s chr7           2702653 75 + 159138663 AGCGAGGCTCCCCAATGGTGCTTTTGGCTTTAGTGTACGA
 s SRR359290.9962       0 75 +        75 AGCGAGGCTCCCCAATGGTGCTTTTGGCTTTAGTGTACGATGTTTGCTGTGCTTCCCGCCGTGGAGGGCAGAGCC
 q SRR359290.9962                        AABA>=C?DBEEEEEEEEEEBEEEBAEBEE?EDEAEBDEDEEB?D=?B:@:ABA@@C@@A??B@@??BEB:??B2
 
-# CPU time: 2.52 seconds
-a score=334 mismap=1e-10
-s chr12           110785491 75 + 133851895 AATGTTTGCGCATGTTCGAGATGAGTCTCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
-s SRR359290.10000         0 75 -        75 AATGTTTGCGCATGTTCGAGATGAGTCGCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
-q SRR359290.10000                          ############################B13BB:B1;*6 at CCC=CACCC3A9A>12:;+9A;>2+A6AA=C?:5?
-
 a score=450 mismap=1e-10
 s chr17          72765010 75 + 81195210 CCTCCTCCTCACTGAGTGCCTCATCCTACCGGGTGTCCCTTTGCCACCCTGCCTGGGACATCGCTGGAACCTGCA
 s SRR359290.9964        0 75 -       75 CCTCCTCCTCACTGAGTGCCTCATCCTACCGGGTGTCCCTTTGCCACCCTGCCTGGGACATCGCTGGAACCTGCA
@@ -17388,6 +17379,11 @@ s chr15          79183831 75 + 102531392 GAAATGGAGATGGTGGCAGTACCAGTGAGACCCCTCAGC
 s SRR359290.9998        0 75 -        75 GAAATGGAGATGGTGGCAGTACCAGTGAGACCCCTCAGCCTCCTCGGAAGAAAAGGGCCCGGGTAGATCCTACTG
 q SRR359290.9998                         ECCE?DECDD==A=CDEBAECDE?FFG?EEFFGFFGEEDF=EEGGGGEGFGFGGGGGGGGGGGFGFGGFGAFGDG
 
+a score=334 mismap=1e-10
+s chr12           110785491 75 + 133851895 AATGTTTGCGCATGTTCGAGATGAGTCTCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
+s SRR359290.10000         0 75 -        75 AATGTTTGCGCATGTTCGAGATGAGTCGCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
+q SRR359290.10000                          ############################B13BB:B1;*6 at CCC=CACCC3A9A>12:;+9A;>2+A6AA=C?:5?
+
 # LAST version 356
 #
 # a=21 b=9 A=21 B=9 c=100000 F=-1 e=120 d=115 x=119 y=44 z=119
@@ -20918,12 +20914,6 @@ s chr7           2702653 75 + 159138663 AGCGAGGCTCCCCAATGGTGCTTTTGGCTTTAGTGTACGA
 s SRR359290.9962       0 75 +        75 AGCGAGGCTCCCCAATGGTGCTTTTGGCTTTAGTGTACGATGTTTGCTGTGCTTCCCGCCGTGGAGGGCAGAGCC
 q SRR359290.9962                        AABA>=C?DBEEEEEEEEEEBEEEBAEBEE?EDEAEBDEDEEB?D=?B:@:ABA@@C@@A??B@@??BEB:??B2
 
-# CPU time: 2.52 seconds
-a score=334 mismap=1e-10
-s chr12           110785491 75 + 133851895 AATGTTTGCGCATGTTCGAGATGAGTCTCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
-s SRR359290.10000         0 75 -        75 AATGTTTGCGCATGTTCGAGATGAGTCGCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
-q SRR359290.10000                          ############################B13BB:B1;*6 at CCC=CACCC3A9A>12:;+9A;>2+A6AA=C?:5?
-
 a score=450 mismap=1e-10
 s chr17          72765010 75 + 81195210 CCTCCTCCTCACTGAGTGCCTCATCCTACCGGGTGTCCCTTTGCCACCCTGCCTGGGACATCGCTGGAACCTGCA
 s SRR359290.9964        0 75 -       75 CCTCCTCCTCACTGAGTGCCTCATCCTACCGGGTGTCCCTTTGCCACCCTGCCTGGGACATCGCTGGAACCTGCA
@@ -21054,6 +21044,11 @@ s chrM           2277 75 + 16571 CTATCACCCTATAGAAGAACTAATGTTAGTATAAGTAACATGAAAAC
 s SRR359290.9997    0 75 +    75 CTATCACCCTATAGAAGAACTAATGTTAGTATAAGTAACATGAAAACATTCTCCTCCGCATAAGCCTGCGTCAGA
 q SRR359290.9997                 C55BCCCCCBC=CCBCB4CC8C@=CBBCCBBB=CC=CCCB=A?CBCCC4BC?BB?BCBBBABBBCC??BB?BBBB
 
+a score=334 mismap=1e-10
+s chr12           110785491 75 + 133851895 AATGTTTGCGCATGTTCGAGATGAGTCTCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
+s SRR359290.10000         0 75 -        75 AATGTTTGCGCATGTTCGAGATGAGTCGCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
+q SRR359290.10000                          ############################B13BB:B1;*6 at CCC=CACCC3A9A>12:;+9A;>2+A6AA=C?:5?
+
 # LAST version 356
 #
 # a=21 b=9 A=21 B=9 c=100000 F=-1 e=120 d=115 x=119 y=44 z=119
@@ -25850,13 +25845,6 @@ s SRR359290.9962       0 75 +        75 AGCGAGGCTCCCCAATGGTGCTTTTGGCTTTAGTGTACGA
 q SRR359290.9962                        AABA>=C?DBEEEEEEEEEEBEEEBAEBEE?EDEAEBDEDEEB?D=?B:@:ABA@@C@@A??B@@??BEB:??B2
 p                                       &,28>DJPV\bhntz~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ztnhb\VPJD>82,&
 
-# CPU time: 2.52 seconds
-a score=334 mismap=1e-10
-s chr12           110785491 75 + 133851895 AATGTTTGCGCATGTTCGAGATGAGTCTCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
-s SRR359290.10000         0 75 -        75 AATGTTTGCGCATGTTCGAGATGAGTCGCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
-q SRR359290.10000                          ############################B13BB:B1;*6 at CCC=CACCC3A9A>12:;+9A;>2+A6AA=C?:5?
-p                                          "$&(*,.02468:<>@BDFHJLNPRTVX^djpv|~~~~~~~~~~~~~~~~~~~~~~~~~~ztnhb\VPJD>82,&
-
 a score=450 mismap=1e-10
 s chr17          72765010 75 + 81195210 CCTCCTCCTCACTGAGTGCCTCATCCTACCGGGTGTCCCTTTGCCACCCTGCCTGGGACATCGCTGGAACCTGCA
 s SRR359290.9964        0 75 -       75 CCTCCTCCTCACTGAGTGCCTCATCCTACCGGGTGTCCCTTTGCCACCCTGCCTGGGACATCGCTGGAACCTGCA
@@ -26061,6 +26049,12 @@ s SRR359290.9998        0 75 -        75 GAAATGGAGATGGTGGCAGTACCAGTGAGACCCCTCAGC
 q SRR359290.9998                         ECCE?DECDD==A=CDEBAECDE?FFG?EEFFGFFGEEDF=EEGGGGEGFGFGGGGGGGGGGGFGFGGFGAFGDG
 p                                        &,15788888888888888888888888888888888888888888888888888888888888888888751,&
 
+a score=334 mismap=1e-10
+s chr12           110785491 75 + 133851895 AATGTTTGCGCATGTTCGAGATGAGTCTCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
+s SRR359290.10000         0 75 -        75 AATGTTTGCGCATGTTCGAGATGAGTCGCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
+q SRR359290.10000                          ############################B13BB:B1;*6 at CCC=CACCC3A9A>12:;+9A;>2+A6AA=C?:5?
+p                                          "$&(*,.02468:<>@BDFHJLNPRTVX^djpv|~~~~~~~~~~~~~~~~~~~~~~~~~~ztnhb\VPJD>82,&
+
 # LAST version 356
 #
 # a=21 b=9 A=21 B=9 c=100000 F=-1 e=120 d=115 x=119 y=44 z=119
@@ -30457,13 +30451,6 @@ s SRR359290.9962       0 75 +        75 AGCGAGGCTCCCCAATGGTGCTTTTGGCTTTAGTGTACGA
 q SRR359290.9962                        AABA>=C?DBEEEEEEEEEEBEEEBAEBEE?EDEAEBDEDEEB?D=?B:@:ABA@@C@@A??B@@??BEB:??B2
 p                                       ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 
-# CPU time: 2.52 seconds
-a score=334 mismap=1e-10
-s chr12           110785491 75 + 133851895 AATGTTTGCGCATGTTCGAGATGAGTCTCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
-s SRR359290.10000         0 75 -        75 AATGTTTGCGCATGTTCGAGATGAGTCGCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
-q SRR359290.10000                          ############################B13BB:B1;*6 at CCC=CACCC3A9A>12:;+9A;>2+A6AA=C?:5?
-p                                          ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-
 a score=450 mismap=1e-10
 s chr17          72765010 75 + 81195210 CCTCCTCCTCACTGAGTGCCTCATCCTACCGGGTGTCCCTTTGCCACCCTGCCTGGGACATCGCTGGAACCTGCA
 s SRR359290.9964        0 75 -       75 CCTCCTCCTCACTGAGTGCCTCATCCTACCGGGTGTCCCTTTGCCACCCTGCCTGGGACATCGCTGGAACCTGCA
@@ -30632,6 +30619,12 @@ s SRR359290.9997    0 75 +    75 CTATCACCCTATAGAAGAACTAATGTTAGTATAAGTAACATGAAAAC
 q SRR359290.9997                 C55BCCCCCBC=CCBCB4CC8C@=CBBCCBBB=CC=CCCB=A?CBCCC4BC?BB?BCBBBABBBCC??BB?BBBB
 p                                }}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}
 
+a score=334 mismap=1e-10
+s chr12           110785491 75 + 133851895 AATGTTTGCGCATGTTCGAGATGAGTCTCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
+s SRR359290.10000         0 75 -        75 AATGTTTGCGCATGTTCGAGATGAGTCGCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
+q SRR359290.10000                          ############################B13BB:B1;*6 at CCC=CACCC3A9A>12:;+9A;>2+A6AA=C?:5?
+p                                          ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
 # LAST version 356
 #
 # a=21 b=9 A=21 B=9 c=100000 F=-1 e=120 d=115 x=119 y=44 z=119
@@ -35028,13 +35021,6 @@ s SRR359290.9962       0 75 +        75 AGCGAGGCTCCCCAATGGTGCTTTTGGCTTTAGTGTACGA
 q SRR359290.9962                        AABA>=C?DBEEEEEEEEEEBEEEBAEBEE?EDEAEBDEDEEB?D=?B:@:ABA@@C@@A??B@@??BEB:??B2
 p                                       ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 
-# CPU time: 2.52 seconds
-a score=334 mismap=1e-10
-s chr12           110785491 75 + 133851895 AATGTTTGCGCATGTTCGAGATGAGTCTCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
-s SRR359290.10000         0 75 -        75 AATGTTTGCGCATGTTCGAGATGAGTCGCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
-q SRR359290.10000                          ############################B13BB:B1;*6 at CCC=CACCC3A9A>12:;+9A;>2+A6AA=C?:5?
-p                                          ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-
 a score=450 mismap=1e-10
 s chr17          72765010 75 + 81195210 CCTCCTCCTCACTGAGTGCCTCATCCTACCGGGTGTCCCTTTGCCACCCTGCCTGGGACATCGCTGGAACCTGCA
 s SRR359290.9964        0 75 -       75 CCTCCTCCTCACTGAGTGCCTCATCCTACCGGGTGTCCCTTTGCCACCCTGCCTGGGACATCGCTGGAACCTGCA
@@ -35203,6 +35189,12 @@ s SRR359290.9997    0 75 +    75 CTATCACCCTATAGAAGAACTAATGTTAGTATAAGTAACATGAAAAC
 q SRR359290.9997                 C55BCCCCCBC=CCBCB4CC8C@=CBBCCBBB=CC=CCCB=A?CBCCC4BC?BB?BCBBBABBBCC??BB?BBBB
 p                                }}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}
 
+a score=334 mismap=1e-10
+s chr12           110785491 75 + 133851895 AATGTTTGCGCATGTTCGAGATGAGTCTCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
+s SRR359290.10000         0 75 -        75 AATGTTTGCGCATGTTCGAGATGAGTCGCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
+q SRR359290.10000                          ############################B13BB:B1;*6 at CCC=CACCC3A9A>12:;+9A;>2+A6AA=C?:5?
+p                                          ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
 # LAST version 356
 #
 # a=21 b=9 A=21 B=9 c=100000 F=-1 e=120 d=115 x=119 y=44 z=119
@@ -39599,13 +39591,6 @@ s SRR359290.9962       0 75 +        75 AGCGAGGCTCCCCAATGGTGCTTTTGGCTTTAGTGTACGA
 q SRR359290.9962                        AABA>=C?DBEEEEEEEEEEBEEEBAEBEE?EDEAEBDEDEEB?D=?B:@:ABA@@C@@A??B@@??BEB:??B2
 p                                       ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 
-# CPU time: 2.52 seconds
-a score=334 mismap=1e-10 sense=0
-s chr12           110785491 75 + 133851895 AATGTTTGCGCATGTTCGAGATGAGTCTCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
-s SRR359290.10000         0 75 -        75 AATGTTTGCGCATGTTCGAGATGAGTCGCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
-q SRR359290.10000                          ############################B13BB:B1;*6 at CCC=CACCC3A9A>12:;+9A;>2+A6AA=C?:5?
-p                                          ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-
 a score=450 mismap=1e-10 sense=0
 s chr17          72765010 75 + 81195210 CCTCCTCCTCACTGAGTGCCTCATCCTACCGGGTGTCCCTTTGCCACCCTGCCTGGGACATCGCTGGAACCTGCA
 s SRR359290.9964        0 75 -       75 CCTCCTCCTCACTGAGTGCCTCATCCTACCGGGTGTCCCTTTGCCACCCTGCCTGGGACATCGCTGGAACCTGCA
@@ -39774,6 +39759,12 @@ s SRR359290.9997    0 75 +    75 CTATCACCCTATAGAAGAACTAATGTTAGTATAAGTAACATGAAAAC
 q SRR359290.9997                 C55BCCCCCBC=CCBCB4CC8C@=CBBCCBBB=CC=CCCB=A?CBCCC4BC?BB?BCBBBABBBCC??BB?BBBB
 p                                }}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}
 
+a score=334 mismap=1e-10 sense=0
+s chr12           110785491 75 + 133851895 AATGTTTGCGCATGTTCGAGATGAGTCTCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
+s SRR359290.10000         0 75 -        75 AATGTTTGCGCATGTTCGAGATGAGTCGCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
+q SRR359290.10000                          ############################B13BB:B1;*6 at CCC=CACCC3A9A>12:;+9A;>2+A6AA=C?:5?
+p                                          ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
 # LAST version 356
 #
 # a=21 b=9 A=21 B=9 c=100000 F=-1 e=120 d=115 x=119 y=44 z=119
@@ -43479,12 +43470,6 @@ s chr7           2702653 75 + 159138663 AGCGAGGCTCCCCAATGGTGCTTTTGGCTTTAGTGTACGA
 s SRR359290.9962       0 75 +        75 AGCGAGGCTCCCCAATGGTGCTTTTGGCTTTAGTGTACGATGTTTGCTGTGCTTCCCGCCGTGGAGGGCAGAGCC
 p                                       &,28>DJPV\bhntz~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ztnhb\VPJD>82,&
 
-# CPU time: 2.52 seconds
-a score=426 mismap=1e-10
-s chr12           110785491 75 + 133851895 AATGTTTGCGCATGTTCGAGATGAGTCTCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
-s SRR359290.10000         0 75 -        75 AATGTTTGCGCATGTTCGAGATGAGTCGCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
-p                                          &,28>DJPV\bhntz~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ztnhb\VPJD>82,&
-
 a score=450 mismap=1e-10
 s chr17          72765010 75 + 81195210 CCTCCTCCTCACTGAGTGCCTCATCCTACCGGGTGTCCCTTTGCCACCCTGCCTGGGACATCGCTGGAACCTGCA
 s SRR359290.9964        0 75 -       75 CCTCCTCCTCACTGAGTGCCTCATCCTACCGGGTGTCCCTTTGCCACCCTGCCTGGGACATCGCTGGAACCTGCA
@@ -43615,6 +43600,11 @@ s chrM           2277 75 + 16571 CTATCACCCTATAGAAGAACTAATGTTAGTATAAGTAACATGAAAAC
 s SRR359290.9997    0 75 +    75 CTATCACCCTATAGAAGAACTAATGTTAGTATAAGTAACATGAAAACATTCTCCTCCGCATAAGCCTGCGTCAGA
 p                                &,28>DJPV\bhnsxz{{|||||||||||||||||||||||||||||||||||||||{{zxsnhb\VPJD>82,&
 
+a score=426 mismap=1e-10
+s chr12           110785491 75 + 133851895 AATGTTTGCGCATGTTCGAGATGAGTCTCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
+s SRR359290.10000         0 75 -        75 AATGTTTGCGCATGTTCGAGATGAGTCGCACCAACAGTGTGTAAGTCATTAACAGTCCTAACTGTGGTGTTTTCC
+p                                          &,28>DJPV\bhntz~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ztnhb\VPJD>82,&
+
 # LAST version 761
 #
 # a=14 b=3 A=19 B=3 e=110 d=90 x=109 y=44 z=109 D=10 E=1.70365e+07


=====================================
test/maf-convert-test.out
=====================================
The diff for this file was not included because it is too large.

=====================================
test/maf-swap-test.out
=====================================
@@ -20,6 +20,54 @@ q 1                    BB at 1;?@B at BBBC@BBBCACBBA>B??=?A at B?B?B=7 at 4BB;B>@@@A?>=B at A@>
 s chrM 1543 87 + 16571 TTATATAGAGGAGACAAGTCGTAACATGGTAAGTGTACTGGAAAGTGCACTTGGACGAACCAGAGTGTAGCTTAACACAAAGCACCC
 p                      $(-16:@FLRW]bddfgfccedffdedffdfgfdedeccfggggfdedecbcfgfdfgfddfggfdedfecba`]WTOKF@:40+'$
 
+a score=364 mismap=8.24e-08
+s 2        3 84 +    87 GGTTTTAAAAATTTTGGTGTAATTTTAAATAAAAGTAATAATTATGTATATTATTATAATTATTCTAACTTTGATTTTTTTAAT
+q 2                     B at B>BB?B;=B@@=BBABA@<A5??'8@@>=3 at 47>@+AA9?;@@>@9?@(A>-9 at A=(=>?1;=7=96 at .4<(49=.2.1(==
+s chrM 12300 84 + 16571 GGCCCCAAAAATTTTGGTGCAACTCCAAATAAAAGTAATAACCATGCACACTACTATAACCACCCTAACCCTGACTTCCCTAAT
+p                       &+.158>DJPUWZ]_dededfecbbbegfdfgggfdffdffccedededdbcdccedfeccdbbbbb_ZWTQNID@=:741.*$
+
+a score=258 mismap=8.71e-08
+s 3       0 69 +    87 ATCATTTTTTGGTTTGGGTATATGAATATTGTTGGGGGGAAGAGATTGATAGTAATGGTGGATGCAATT
+q 3                    B>)(?'=@?<AABB()669@@>4BB8;*7+)=BC)383):4*;AA>11<@+%+?9)4AB>5,3;7%25#
+s chrM 4049 69 - 16571 ATAACTTCTTGGTCTAGGCACATGAATATTGTTGTGGGGAAGAGACTGATAATAAAGGTGGATGCGACA
+p                      #%%),/258<AEFHHIOTW\_ccfgfddbbca`^Z`dffggfggfcceda_Z[ZVQQOJGB<62-)(%#
+
+a score=408 mismap=8.02e-08
+s 4       0 87 +    87 CGTCGTCGGGAAAAAAGGCGGGAGAAGTTTCGGTAGGTTTGAAGCTGTTTTTTCGAATTTGTAATTTAATATGAAAATTATTTCGGA
+q 4                    ?BCC@?CCBCC at CBCACCCCCCCC@ABCCAABBCCBC at A=B at CBCACCC@BCB@?CCACBCBBBBB>C>@>BBB?BBB:=<B at B>BB
+s chrM 5737 87 + 16571 CGCCGCCGGGAAAAAAGGCGGGAGAAGCCCCGGCAGGTTTGAAGCTGCTTCTTCGAATTTGCAATTCAATATGAAAATCACCTCGGA
+p                      $),059<BHNTZ_dfggfdgggggggfcbbcffdfgfcccfggfccdbbaaabcfffcbcedffcbcffdecda\VPJFC>:741,&
+
+a score=280 mismap=9.01e-08
+s 5        0 87 +    87 GGGGTATCGATTATTAGGGGAATTTGTTAGTTGTTAAAGTTTTCGATTATGATGGGTATTATTATAACAGTTGCTAAGAAGAACCCA
+q 5                     BC?CCCBAABABBC>@>=BB5C>BB?B3B;BB<CB4A.B/@86B?/7>5?@@>@?>6::B</<8A######################
+s chrM 10398 87 - 16571 GGGGCACCGATTATTAGGGGAACTAGTCAGTTGCCAAAGCCTCCGATTATGATGGGTATTACTATGAAGAAGATTATTACAAATGCA
+p                       &,27;@CGMRUX]_aeggfd`[UQNSVZ_bbcdccfggfcbbbcefccedfeba]XRNIEB=961-+*(''&%%$$$####"""!!!
+
+a score=368 mismap=8.65e-08
+s 6        0 87 +    87 TTATATTTTATATTTTTTTCGAAATTATATTTATTTTTATTTTGGTTATTATTATTCGATGAGGTAATTAGTTAGAATGTCTGAACG
+q 6                     BBB?C<@BA;5AB:AB at B?B<AB>A<AB:@AB9;=@AA43A>@;,5=@@A@/@>@=:?:@>><>-=6 at 746?=05?*?<<==9775=
+s chrM 11106 87 + 16571 TCATATTTTATATCTTCTTCGAAACCACACTTATCCCCACCTTGGCTATCATCACCCGATGAGGCAACCAGCCAGAACGCCTGAACG
+p                       $(-058;?BGKORTVWYZ\^begfccedecbcdbbbbbdbbbcefccecceccdcbcffdfggfdffccfecba]XRMID@=:50*&
+
+a score=340 mismap=9.13e-08
+s 7        0 85 +    87 TTATAAGTATTAATTAATTAATATATTTTAGTAGGTATGTTTGCTTGTAATATTGAATGTAGGAGCGATAGATAAGAGAATTAGG
+q 7                     BCC@>CCBCBBCBCBB;B?6,:>@?BB<?:+4>89ABBB at AB?):B=5=1@@49;;9;<9(22+/+>?060.=<4++276;####
+s chrM 16343 85 - 16571 CTACAAGCATTAATTAATTAACACACTTTAGTAGGTATGTTCGCCTGTAATATTGAACGTAGGTGCGATAAATAATAGAATGAGG
+p                       $(-17<BEJNQW\^`deccffdedecbbcffdfgfdedecbbdbbcedffdeccfgfdecca]WYYXTOKFFEDA=<;83-)'%#
+
+a score=368 mismap=8.66e-08
+s 8       0 87 +    87 AGTATATTCTTTAATTATTTATATAGGATGAATAATAGTAGTTTTATCGTATAATTTTAATATATCTATTTTTAATTTAATTATTTA
+q 8                    B at C=BBCCBC?CCCC>BCC<B;@B@@BA2>C:<ABCC<A at A=@BB- at B?:B<B7>?AB:A@@@>(:AB@=B;B at BB>@<@>@A at ABB
+s chrM 5006 87 + 16571 AGCATACTCCTCAATTACCCACATAGGATGAATAATAGCAGTTCTACCGTACAACCCTAACATAACCATTCTTAATTTAACTATTTA
+p                      &+/48=@CFILPUY\^`aabddedfggfdfgfdffdffdffcbbcdccededfecbbcefddb`\^`ba``^]ZVPMIFA<84/,(%
+
+a score=375 mismap=9.32e-08
+s 9       0 87 +    87 TGAATTTACGAGTATATCGATTACGGTGGATTAATTTTTAATTTTTATATATTTTTTTTATTATTTTTAGAATTAGGTGATTTGCGA
+q 9                    B=9=B+BBBB>BBB=:BB@;<@@B:<A4@?A:?=BABAAAA>>@<A>@AA>@@@?@@6878:6588?<:>=?9<::=<51<99>;48
+s chrM 7901 87 + 16571 TGAACCTACGAGTACACCGACTACGGCGGACTAATCTTCAACTCCTACATACTTCCCCCATTATTCCTAGAACCAGGCGACCTGCGA
+p                      $*058;?DHNTY\`bdccffccedffdfgfccfecbbbceecbbbcededebbaaaabbdccdbbbbcfggea_]YSMID>:74/+&
+
 a score=207 mismap=1.15e-07
 s 10      0 60 +    87 TTTTTTTTTTTAATCATATATTTATGGTTAATTTTTTATTTTTTATTGTATTTATTCTAA
 q 10                   @6)>B+A9<)?>>A(?A=-(48/<5?=1>4>:B@@<-:@@5;/9==ABA47*5':B<###
@@ -74,12 +122,6 @@ q 19                    =@9 at B>C<';>)*8?'>67 at 1;(?<;A<@:@@(9A@@)(38>9+@<8?6 at 4:=?B9
 s chrM 14174 87 + 16571 ACAATATATACACCAACAAACAATGTTCAACCAGTAACCACTACTAATCAACGCCCATAATCATACAAAGCCCCCGCACCAATAGGA
 p                       &*/58=AEGJKKNRW\_cffdfedecbcec_]Z_addbbdcceccffccffdecbcedffccededfggeb`_^\ZURMIFA;72,&
 
-a score=364 mismap=8.24e-08
-s 2        3 84 +    87 GGTTTTAAAAATTTTGGTGTAATTTTAAATAAAAGTAATAATTATGTATATTATTATAATTATTCTAACTTTGATTTTTTTAAT
-q 2                     B at B>BB?B;=B@@=BBABA@<A5??'8@@>=3 at 47>@+AA9?;@@>@9?@(A>-9 at A=(=>?1;=7=96 at .4<(49=.2.1(==
-s chrM 12300 84 + 16571 GGCCCCAAAAATTTTGGTGCAACTCCAAATAAAAGTAATAACCATGCACACTACTATAACCACCCTAACCCTGACTTCCCTAAT
-p                       &+.158>DJPUWZ]_dededfecbbbegfdfgggfdffdffccedededdbcdccedfeccdbbbbb_ZWTQNID@=:741.*$
-
 a score=329 mismap=1.29e-07
 s 20      0 87 +    87 TTATAGTTTAATGTTTTATTTAGTTATTTTATTTTATTTTTACTGATGTTCGTCGATCGTTGATTATTATTTATAGATTATAAAGAT
 q 20                   )@,27 at CBC:1?>C?-2AA?59'<)?B7A46)@>:*?8,5)<24?)5>?;?(2@@2(1>-:<3 at B5,<,3>7)?2%81,;@@BAA?B
@@ -140,12 +182,6 @@ q 29                   ?6 at C@5<(<CCB.CA?4BB=4?AA;=@CA=A<%=5>BB=15=B=>>6<B??):=>-<
 s chrM 4775 81 + 16571 AGCAATAAAACTAGGAATAGCCCCCTTTCACTTCTGAGTCCCAGAGGTTACCCAAGGCACCCCTCTGACATCCGGCCTGCT
 p                      &+/5:>DIOTW[`dfgfdfecbaaaaabbdba`abeffcbbcfgggfccdbacfggfddb```_^\ZWROIFB?:62.+'$
 
-a score=258 mismap=8.71e-08
-s 3       0 69 +    87 ATCATTTTTTGGTTTGGGTATATGAATATTGTTGGGGGGAAGAGATTGATAGTAATGGTGGATGCAATT
-q 3                    B>)(?'=@?<AABB()669@@>4BB8;*7+)=BC)383):4*;AA>11<@+%+?9)4AB>5,3;7%25#
-s chrM 4049 69 - 16571 ATAACTTCTTGGTCTAGGCACATGAATATTGTTGTGGGGAAGAGACTGATAATAAAGGTGGATGCGACA
-p                      #%%),/258<AEFHHIOTW\_ccfgfddbbca`^Z`dffggfggfcceda_Z[ZVQQOJGB<62-)(%#
-
 a score=349 mismap=1.08e-07
 s 30      0 87 +    87 ATTATATTTATAACGTTATCGTTATAGTTTATGTATTTGTAATAATTTTTTTTATAGTAATATTTATTAGAATCGGAGGTTTTGGTC
 q 30                   BB@<BBCBBBB at ABBB>ABB at B>BBA3>-@>BBCABAB>B6:+4BB;B8C?@?B=@BB?A;@?A?B@@B at BB@>-@:@?:>@A81>#
@@ -206,12 +242,6 @@ q 39                   B?=BBBBB?B at BABBB=@BABBAAA at BBAA>;0===8>=;7.5:<;@@@A@<8A>35
 s chrM 7893 87 + 16571 AATGGTACTGAACCTACGAGTACACCGACTACGGCGGACTAATCTTCAACTCCTACATACTTCCCCCATTATTCCTAGAACCAGGCG
 p                      &+/5:>CFJPVZ]_addfgfdedeccffccedffdfgfccfecbbbceecbbbcededebaa`aaabdba`\YVSPMHB=7321.)&
 
-a score=408 mismap=8.02e-08
-s 4       0 87 +    87 CGTCGTCGGGAAAAAAGGCGGGAGAAGTTTCGGTAGGTTTGAAGCTGTTTTTTCGAATTTGTAATTTAATATGAAAATTATTTCGGA
-q 4                    ?BCC@?CCBCC at CBCACCCCCCCC@ABCCAABBCCBC at A=B at CBCACCC@BCB@?CCACBCBBBBB>C>@>BBB?BBB:=<B at B>BB
-s chrM 5737 87 + 16571 CGCCGCCGGGAAAAAAGGCGGGAGAAGCCCCGGCAGGTTTGAAGCTGCTTCTTCGAATTTGCAATTCAATATGAAAATCACCTCGGA
-p                      $),059<BHNTZ_dfggfdgggggggfcbbcffdfgfcccfggfccdbbaaabcfffcbcedffcbcffdecda\VPJFC>:741,&
-
 a score=354 mismap=8.27e-08
 s 40      0 87 +    87 ATTAATATTTTTAATTTTATTTATTTTTTTTTTTTTAGGAGGTTTGTTTTCGTTAATCGGTTTTTTGTCTAAATGGGTTATTATTGA
 q 40                   BB?BBBBBBBBCBBABBAB?B@@B<9BAB>BA at BBB<@BA at ABAAAB=:@BB:BAAB=<>@?4<@B at A@?9=?@=.4?9:+778;32
@@ -272,12 +302,6 @@ q 49                    BCBCBBCCCBB?BBB6CBCCB@@@BAA at C?BB at B@B@?BBC?B at B@@BBB>BBA5@
 s chrM 12896 87 + 16571 CATGATTTATCCTACACTCCAACTCATGAGACCCACAACAAATAGCCCTTCTAAACGCTAATCCAAGCCTCACCCCACTACTAGGCC
 p                       $*-38;>BGJMPTY\``abceecbcedfggfcbcedffdfgfdfecbaaabcfgfdeccfecccfffcba``\YVSPKGD?<83.($
 
-a score=280 mismap=9.01e-08
-s 5        0 87 +    87 GGGGTATCGATTATTAGGGGAATTTGTTAGTTGTTAAAGTTTTCGATTATGATGGGTATTATTATAACAGTTGCTAAGAAGAACCCA
-q 5                     BC?CCCBAABABBC>@>=BB5C>BB?B3B;BB<CB4A.B/@86B?/7>5?@@>@?>6::B</<8A######################
-s chrM 10398 87 - 16571 GGGGCACCGATTATTAGGGGAACTAGTCAGTTGCCAAAGCCTCCGATTATGATGGGTATTACTATGAAGAAGATTATTACAAATGCA
-p                       &,27;@CGMRUX]_aeggfd`[UQNSVZ_bbcdccfggfcbbbcefccedfeba]XRNIEB=961-+*(''&%%$$$####"""!!!
-
 a score=151 mismap=1.69e-07
 s 50      0 36 +    87 TGGGTTAGTTAGGTAATTTTTTAGGTAACGATGATA
 q 50                   ?)3<2A>:>*?6<B8>@7BB@<@B@>B0*(26)8B<
@@ -338,12 +362,6 @@ q 59                   BCA=8@>:B>AB9B-B?A0*:BABB?B?=6><?<BA at 5A99>024)?B>A@:>@@A9
 s chrM 6572 82 - 16571 AAAAGAGTAAGACCCTCATCAATAGATGGAGACATACAGAAATAGTCAAACCACATCTACAAAATGCCAGTATCAGGCGGCG
 p                      &,28>DIMSY^aaaabbdcbefdfgfdggggfdededfgggfdfebcfgfccecdba_^ZVRW\^`^\YUPLGC@;5/+&#"
 
-a score=368 mismap=8.65e-08
-s 6        0 87 +    87 TTATATTTTATATTTTTTTCGAAATTATATTTATTTTTATTTTGGTTATTATTATTCGATGAGGTAATTAGTTAGAATGTCTGAACG
-q 6                     BBB?C<@BA;5AB:AB at B?B<AB>A<AB:@AB9;=@AA43A>@;,5=@@A@/@>@=:?:@>><>-=6 at 746?=05?*?<<==9775=
-s chrM 11106 87 + 16571 TCATATTTTATATCTTCTTCGAAACCACACTTATCCCCACCTTGGCTATCATCACCCGATGAGGCAACCAGCCAGAACGCCTGAACG
-p                       $(-058;?BGKORTVWYZ\^begfccedecbcdbbbbbdbbbcefccecceccdcbcffdfggfdffccfecba]XRMID@=:50*&
-
 a score=383 mismap=8.82e-08
 s 60      0 87 +    87 GTAGAGATTTGAAGAATTGATTAATGTTTGGGTTTGAGTTTATATATTATAGTGAGAATTTTATGATGGATTATGTAACGAATAATG
 q 60                   B>B>=)=(@3:87/;7;(;:B8((:+@<)@7<B@:27(<ABAABAAB@=5()-@,1)>7.2A?@9=@?60((5=15A=31<@:?>@#
@@ -392,12 +410,6 @@ q 69                   B>BB;1ABB?@9ABAABAB at 4)?@A>(@>>*>>53A<>0@>@BAA@:7;)@<?AA at A
 s chrM 5660 87 + 16571 GACCAATGGGACTTAAACCCACAAACACTTAGTTAACAGCTAAGCACCCTAATCAACTGGCTTCAATCTACTTCTCCCGCCGCCGGG
 p                      &+/28=AGMSWZ]`dfeba`^Z_cecdbbcdeccffdffccfgfdecbbbeeccffccedbaaacda_][VROLIFC@=852-*'$"
 
-a score=340 mismap=9.13e-08
-s 7        0 85 +    87 TTATAAGTATTAATTAATTAATATATTTTAGTAGGTATGTTTGCTTGTAATATTGAATGTAGGAGCGATAGATAAGAGAATTAGG
-q 7                     BCC@>CCBCBBCBCBB;B?6,:>@?BB<?:+4>89ABBB at AB?):B=5=1@@49;;9;<9(22+/+>?060.=<4++276;####
-s chrM 16343 85 - 16571 CTACAAGCATTAATTAATTAACACACTTTAGTAGGTATGTTCGCCTGTAATATTGAACGTAGGTGCGATAAATAATAGAATGAGG
-p                       $(-17<BEJNQW\^`deccffdedecbbcffdfgfdedecbbdbbcedffdeccfgfdecca]WYYXTOKFFEDA=<;83-)'%#
-
 a score=340 mismap=8.71e-08
 s 70      0 84 +    87 TATTTTTTTTATTTGTTTTTTAGTTTTGTATGTTTTTTTTTTAATATTTATAATAAAATTAATTAATATTAATATTTTAGACGT
 q 70                   ?ACBCCCCCACBCBCCC=CCC>@7BCAACCBCBC?CBCCACCCC at ACC?AC at BA<C>CACC?@CBBBC67;@)=@>@A*=:B at 8
@@ -458,12 +470,6 @@ q 79                    BCC?BCBBBACB?B:B>4B:9B>A<@@B at BBBA=B?9 at AB@AA9BA9?@@A14=BB
 s chrM 10826 87 + 16571 ATTTGAATCAACACAACCACCCACAGCCTAATTATTAGCATCATCCCCCTACTATTTTTTAACCAAATCAACAACAACCTATTTAGC
 p                       %),06;@DHMRV[^bdccdcbcedffcbcefcceccffdeccdbbaaabbdccdbbaabcefccfgfcceeba\WSNHDA>952/*$
 
-a score=368 mismap=8.66e-08
-s 8       0 87 +    87 AGTATATTCTTTAATTATTTATATAGGATGAATAATAGTAGTTTTATCGTATAATTTTAATATATCTATTTTTAATTTAATTATTTA
-q 8                    B at C=BBCCBC?CCCC>BCC<B;@B@@BA2>C:<ABCC<A at A=@BB- at B?:B<B7>?AB:A@@@>(:AB@=B;B at BB>@<@>@A at ABB
-s chrM 5006 87 + 16571 AGCATACTCCTCAATTACCCACATAGGATGAATAATAGCAGTTCTACCGTACAACCCTAACATAACCATTCTTAATTTAACTATTTA
-p                      &+/48=@CFILPUY\^`aabddedfggfdfgfdffdffdffcbbcdccededfecbbcefddb`\^`ba``^]ZVPMIFA<84/,(%
-
 a score=359 mismap=9.42e-08
 s 81      0 87 +    87 TTATAGTGAAATGTTTTAATTAAATATTATCGTATGGTTTATTATAATTATTTTTATATTTTTTATATTATTTTTTATTATTTAATT
 q 81                   @C<?BCCBCCB at CB@;BB=4@@>BCAAB at .0;B@@?C at BBA3>AC>8/B=A=C at B<C226A1>8B at B;A<98AB>>)>A?<@@8@<5
@@ -512,12 +518,6 @@ q 89                   9B:):<CCCAC>C:CCAC==B?CCCCCCCCB at BBCBCB<CBCCCB9?<B6BBA;B>C
 s chrM 8407 87 + 16571 CCCCATACTCCTTACACTATTCCTCATCACCCAACTAAAAATATTAAACACAAACTACCACCTACCTCCCTCACCAAAGCCCATAAA
 p                      $'*-26:=@CGJNSV[]`bbaaabbdccdcbcefccfgggfdeccfgfdfdfgfcceccdcbcdbbaaa```a]ZXSMGA=:730,&
 
-a score=375 mismap=9.32e-08
-s 9       0 87 +    87 TGAATTTACGAGTATATCGATTACGGTGGATTAATTTTTAATTTTTATATATTTTTTTTATTATTTTTAGAATTAGGTGATTTGCGA
-q 9                    B=9=B+BBBB>BBB=:BB@;<@@B:<A4@?A:?=BABAAAA>>@<A>@AA>@@@?@@6878:6588?<:>=?9<::=<51<99>;48
-s chrM 7901 87 + 16571 TGAACCTACGAGTACACCGACTACGGCGGACTAATCTTCAACTCCTACATACTTCCCCCATTATTCCTAGAACCAGGCGACCTGCGA
-p                      $*058;?DHNTY\`bdccffccedffdfgfccfecbbbceecbbbcededebbaaaabbdccdbbbbcfggea_]YSMID>:74/+&
-
 a score=278 mismap=1.03e-07
 s 90       0 66 +    87 GTTAAGTTAAGTTGAAATTTTTTAAGTGTAAGTTGGGTGTTTTGTGTTAAGTTATATTTTGGTTAG
 q 90                    8;'6=ABA;1BA.789AA9A>'))4BA??7)?=:(?<1)3 at AA>6<@:'.B@?-?:=328124:##


=====================================
test/maf-swap-test.sh
=====================================
@@ -16,4 +16,4 @@ PATH=../bin:$PATH
     try maf-swap -n1 90089.maf
     try maf-swap -n3 ../examples/multiMito.maf
     try maf-swap frameshift-new.maf
-} 2>&1 | diff -u $(basename $0 .sh).out -
+} 2>&1 | diff -u maf-swap-test.out -



View it on GitLab: https://salsa.debian.org/med-team/last-align/-/compare/d66071ad838061ef3de19947217335b96bf2e96e...04dde6a254a5dacbed8eb4dc37426e605fecc37a

-- 
View it on GitLab: https://salsa.debian.org/med-team/last-align/-/compare/d66071ad838061ef3de19947217335b96bf2e96e...04dde6a254a5dacbed8eb4dc37426e605fecc37a
You're receiving this email because of your account on salsa.debian.org.


-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://alioth-lists.debian.net/pipermail/debian-med-commit/attachments/20220504/a4029625/attachment-0001.htm>


More information about the debian-med-commit mailing list