// ==========================================================================
// SeqAn - The Library for Sequence Analysis
// ==========================================================================
// Copyright (c) 2006-2018, Knut Reinert, FU Berlin
// All rights reserved.
//
// Redistribution and use in source and binary forms, with or without
// modification, are permitted provided that the following conditions are met:
//
// * Redistributions of source code must retain the above copyright
// notice, this list of conditions and the following disclaimer.
// * Redistributions in binary form must reproduce the above copyright
// notice, this list of conditions and the following disclaimer in the
// documentation and/or other materials provided with the distribution.
// * Neither the name of Knut Reinert or the FU Berlin nor the names of
// its contributors may be used to endorse or promote products derived
// from this software without specific prior written permission.
//
// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
// AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
// ARE DISCLAIMED. IN NO EVENT SHALL KNUT REINERT OR THE FU BERLIN BE LIABLE
// FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
// DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
// SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
// CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
// LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
// OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
// DAMAGE.
//
// ==========================================================================
// Author: Andreas Gogol-Doering <andreas.doering@mdc-berlin.de>
// Author: Manuel Holtgrewe <manuel.holtgrewe@fu-berlin.de>
// ==========================================================================
// SEQAN_NO_GENERATED_FORWARDS
// TODO(holtgrew): Currently, operations are a function of the whole gap count, could be of clipped region only.
// TODO(holtgrew): Problem with the gap value, getValue(), value().
#ifndef SEQAN_INCLUDE_SEQAN_ALIGN_GAPS_ARRAY_H_
#define SEQAN_INCLUDE_SEQAN_ALIGN_GAPS_ARRAY_H_
namespace seqan {
// ============================================================================
// Forwards
// ============================================================================
// Internally used tag for creating iterators at the begin of containers.
struct Begin__;
typedef Tag<Begin__> Begin_;
// Internally used tag for creating iterators at the end of containers.
struct End__;
typedef Tag<End__> End_;
// Internally used tag for creating iterators inside of containers.
struct Position__;
typedef Tag<Position__> Position_;
struct ArrayGaps_;
typedef Tag<ArrayGaps_> ArrayGaps;
template <typename TSequence> class Gaps<TSequence, ArrayGaps>;
template <typename TSequence>
inline void _reinitArrayGaps(Gaps<TSequence, ArrayGaps> & gaps);
// ============================================================================
// Tags, Classes, Enums
// ============================================================================
struct ArrayGaps_;
typedef Tag<ArrayGaps_> ArrayGaps;
/*!
* @class ArrayGaps
* @headerfile <seqan/align.h>
* @extends Gaps
* @brief Stores length of gap- and non-gap runs in an array.
*
* @signature template <typename TSequence>
* class Gaps<TSequence, ArrayGaps>
*
* @tparam TSequence The type of the underling sequence.
*/
/*!
* @fn ArrayGaps::Gaps
* @headerfile <seqan/align.h>
* @brief Constructor.
*
* @signature Gaps::Gaps([other]);
* @signature Gaps::Gaps(seq);
*
* @param[in] other Other Gaps object to copy from.
* @param[in] seq Sequence concept to construct the gaps for.
*/
template <typename TSequence>
class Gaps<TSequence, ArrayGaps>
{
public:
// -----------------------------------------------------------------------
// Internal Typedefs
// -----------------------------------------------------------------------
typedef typename Size<Gaps>::Type TSize_;
typedef typename Size<TSequence>::Type TSequenceSize_;
typedef typename Position<Gaps>::Type TPosition_;
typedef typename Position<TSequence>::Type TSequencePosition_;
typedef typename Value<Gaps>::Type TValue_;
typedef String<TSequenceSize_> TArray_;
typedef typename Position<TArray_>::Type TArrayPos_;
// -----------------------------------------------------------------------
// Member Variables
// -----------------------------------------------------------------------
// Holder of the underlying sequence.
Holder<TSequence> _source;
// The array with the alternating gap/source char counts.
TArray_ _array;
// Begin and end position in the source.
TSequencePosition_ _sourceBeginPos, _sourceEndPos;
// Begin and end position in the view.
TPosition_ _clippingBeginPos, _clippingEndPos;
// TODO(holtgrew): The following is a possible optimization.
// // Index of clipping begin and end in the _array seq/gap char count array.
// // This identifies a slice of the view.
// TArrayPos_ _clippingBeginIdx, _clippingEndIdx;
// // Offset within the slice.
// TSequenceSize_ _clippingBeginOffset, _clippingEndOffset;
// -----------------------------------------------------------------------
// Constructors
// -----------------------------------------------------------------------
Gaps() : _sourceBeginPos(0), _sourceEndPos(0), _clippingBeginPos(0), _clippingEndPos(0)//,
// _clippingBeginIdx(0), _clippingEndIdx(0), _clippingBeginOffset(0), _clippingEndOffset(0)
{}
explicit
Gaps(TSequence & seq) :
_source(seq), _sourceBeginPos(0), _sourceEndPos(length(seq)),
_clippingBeginPos(0), _clippingEndPos(length(seq))//,
// _clippingBeginIdx(0), _clippingEndIdx(0), _clippingBeginOffset(0),
// _clippingEndOffset(0)
{
// Initialize array gaps object for ungapped sequence.
_reinitArrayGaps(*this);
}
Gaps(Gaps const & other) :
_source(other._source), _array(other._array), _sourceBeginPos(other._sourceBeginPos),
_sourceEndPos(other._sourceEndPos), _clippingBeginPos(other._clippingBeginPos),
_clippingEndPos(other._clippingEndPos)
{}
// -----------------------------------------------------------------------
// Array Subscript Operator
// -----------------------------------------------------------------------
inline Gaps &
operator=(Gaps const & other)
{
setValue(_source, source(other));
_array = other._array;
_sourceBeginPos = other._sourceBeginPos;
_sourceEndPos = other._sourceEndPos;
_clippingBeginPos = other._clippingBeginPos;
_clippingEndPos = other._clippingEndPos;
return *this;
}
// -----------------------------------------------------------------------
// Array Subscript Operator
// -----------------------------------------------------------------------
inline typename Reference<Gaps>::Type
operator[](TPosition_ const clippedViewPos)
{
return value(*this, clippedViewPos);
}
inline typename Reference<Gaps const>::Type
operator[](TPosition_ const clippedViewPos) const
{
return value(*this, clippedViewPos);
}
};
// ----------------------------------------------------------------------------
// Function swap()
// ----------------------------------------------------------------------------
template <typename TSequence>
void swap(Gaps<TSequence, ArrayGaps> & lhs, Gaps<TSequence, ArrayGaps> & rhs)
{
swap(lhs._source, rhs._source);
swap(lhs._array, rhs._array);
std::swap(lhs._sourceBeginPos, rhs._sourceBeginPos);
std::swap(lhs._sourceEndPos, rhs._sourceEndPos);
std::swap(lhs._clippingBeginPos, rhs._clippingBeginPos);
std::swap(lhs._clippingEndPos, rhs._clippingEndPos);
}
// ============================================================================
// Metafunctions
// ============================================================================
// ============================================================================
// Functions
// ============================================================================
// ----------------------------------------------------------------------------
// Function detach()
// ----------------------------------------------------------------------------
// TODO(holtgrew): Remove? Only used by module blast.
template <typename TSequence>
void detach(Gaps<TSequence, ArrayGaps> & gaps)
{
detach(gaps._source);
}
// ----------------------------------------------------------------------------
// Function _setLength()
// ----------------------------------------------------------------------------
// Set the length, only use if TSequence is Nothing.
template <typename TSequence, typename TSize>
inline void _setLength(Gaps<TSequence, ArrayGaps> & gaps, TSize newLen)
{
// Reset array.
resize(gaps._array, 3);
gaps._array[0] = 0;
gaps._array[1] = newLen;
gaps._array[2] = 0;
// Reset clipping information.
gaps._clippingBeginPos = 0;
gaps._clippingEndPos = newLen;
gaps._sourceBeginPos = 0;
gaps._sourceEndPos = gaps._clippingEndPos;
// gaps._clippingBeginIdx = 1;
// gaps._clippingBeginOffset = 0;
// gaps._clippingEndIdx = 1;
// gaps._clippingEndOffset = value(gaps._source)[1];
}
// ----------------------------------------------------------------------------
// Helper Function _reinitArrayGaps()
// ----------------------------------------------------------------------------
// Reset the array gaps DS such that represents the ungapped sequence.
template <typename TSequence>
inline void _reinitArrayGaps(Gaps<TSequence, ArrayGaps> & gaps)
{
_setLength(gaps, length(value(gaps._source)));
}
// ----------------------------------------------------------------------------
// Function begin()
// ----------------------------------------------------------------------------
// TODO(holtgrew): We'd rather have "TTag const &" here.
template <typename TSequence, typename TTag>
inline typename Iterator<Gaps<TSequence, ArrayGaps> >::Type
begin(Gaps<TSequence, ArrayGaps> & gaps, Tag<TTag> const /*tag*/)
{
typedef typename Iterator<Gaps<TSequence, ArrayGaps> >::Type TIter;
return TIter(gaps, Begin_());
}
// TODO(holtgrew): We'd rather have "TTag const &" here.
template <typename TSequence, typename TTag>
inline typename Iterator<Gaps<TSequence, ArrayGaps> const>::Type
begin(Gaps<TSequence, ArrayGaps> const & gaps, Tag<TTag> const /*tag*/)
{
typedef typename Iterator<Gaps<TSequence, ArrayGaps> const>::Type TIter;
return TIter(gaps, Begin_());
}
// ----------------------------------------------------------------------------
// Function end()
// ----------------------------------------------------------------------------
// TODO(holtgrew): We'd rather have "TTag const &" here.
template <typename TSequence, typename TTag>
inline typename Iterator<Gaps<TSequence, ArrayGaps> >::Type
end(Gaps<TSequence, ArrayGaps> & gaps, Tag<TTag> const /*tag*/)
{
typedef typename Iterator<Gaps<TSequence, ArrayGaps> >::Type TIter;
return TIter(gaps, End_());
}
// TODO(holtgrew): We'd rather have "TTag const &" here.
template <typename TSequence, typename TTag>
inline typename Iterator<Gaps<TSequence, ArrayGaps> const>::Type
end(Gaps<TSequence, ArrayGaps> const & gaps, Tag<TTag> const /*tag*/)
{
typedef typename Iterator<Gaps<TSequence, ArrayGaps> const>::Type TIter;
return TIter(gaps, End_());
}
// ----------------------------------------------------------------------------
// Function iter()
// ----------------------------------------------------------------------------
// TODO(holtgrew): We'd rather have "TTag const &" here.
template <typename TSequence, typename TTag, typename TPosition>
inline typename Iterator<Gaps<TSequence, ArrayGaps> >::Type
iter(Gaps<TSequence, ArrayGaps> & gaps, TPosition pos, Tag<TTag> const /*tag*/)
{
typedef typename Iterator<Gaps<TSequence, ArrayGaps> >::Type TIter;
return TIter(gaps, pos, Position_());
}
// TODO(holtgrew): We'd rather have "TTag const &" here.
template <typename TSequence, typename TTag, typename TPosition>
inline typename Iterator<Gaps<TSequence, ArrayGaps> const>::Type
iter(Gaps<TSequence, ArrayGaps> const & gaps, TPosition pos, Tag<TTag> const /*tag*/)
{
typedef typename Iterator<Gaps<TSequence, ArrayGaps> const>::Type TIter;
return TIter(gaps, pos, Position_());
}
// ----------------------------------------------------------------------------
// Function length()
// ----------------------------------------------------------------------------
template <typename TSequence>
inline typename Size<Gaps<TSequence, ArrayGaps> >::Type
length(Gaps<TSequence, ArrayGaps> const & gaps)
{
SEQAN_ASSERT_GEQ(gaps._clippingEndPos, gaps._clippingBeginPos);
return gaps._clippingEndPos - gaps._clippingBeginPos;
}
// ----------------------------------------------------------------------------
// Function unclippedLength()
// ----------------------------------------------------------------------------
template <typename TSequence>
inline typename Size<Gaps<TSequence, ArrayGaps> >::Type
unclippedLength(Gaps<TSequence, ArrayGaps> const & gaps)
{
typedef typename Size<Gaps<TSequence, ArrayGaps> >::Type TSize;
TSize result = 0;
for (unsigned i = 0; i < length(gaps._array); ++i)
result += gaps._array[i];
return result;
}
// ----------------------------------------------------------------------------
// Function createSource()
// ----------------------------------------------------------------------------
// TODO(holtgrew): Remove? Switch to Hosted Type Interface?
template <typename TSequence>
inline void createSource(Gaps<TSequence, ArrayGaps> & gaps)
{
create(gaps._source);
}
// ----------------------------------------------------------------------------
// Function source()
// ----------------------------------------------------------------------------
template <typename TSequence>
inline typename Source<Gaps<TSequence, ArrayGaps> const>::Type &
source(Gaps<TSequence, ArrayGaps> const & gaps)
{
return value(gaps._source);
}
template <typename TSequence>
inline typename Source<Gaps<TSequence, ArrayGaps> >::Type &
source(Gaps<TSequence, ArrayGaps> & gaps)
{
return value(gaps._source);
}
// ----------------------------------------------------------------------------
// Function setSource()
// ----------------------------------------------------------------------------
// TODO(holtgrew): Test with clippings, also for AnchorGaps.
template <typename TSequence>
inline void
setSource(Gaps<TSequence, ArrayGaps> & gaps, TSequence & source)
{
setValue(gaps._source, source);
_reinitArrayGaps(gaps);
}
template <typename TSequence>
inline void
setSource(Gaps<TSequence const, ArrayGaps> & gaps, TSequence & source)
{
setValue(gaps._source, source);
_reinitArrayGaps(gaps);
}
// ----------------------------------------------------------------------------
// Function assignSource()
// ----------------------------------------------------------------------------
template <typename TSequence, typename TSequence2>
inline void
assignSource(Gaps<TSequence, ArrayGaps> & gaps, TSequence2 const & source)
{
value(gaps._source) = source;
_reinitArrayGaps(gaps);
}
// ----------------------------------------------------------------------------
// Function toSourcePosition()
// ----------------------------------------------------------------------------
template <typename TSequence, typename TPosition>
inline typename Position<TSequence>::Type
toSourcePosition(Gaps<TSequence, ArrayGaps> const & gaps, TPosition const clippedViewPos, LeftOfViewPos const & /*tag*/)
{
typedef Gaps<TSequence, ArrayGaps> TGaps;
typedef typename Position<TGaps>::Type TGapsPos;
typedef typename TGaps::TArrayPos_ TArrayPos;
typedef typename Position<TSequence>::Type TSourcePos;
// Translate from clipped view position to unclipped view position.
TGapsPos unclippedViewPos = clippedViewPos + clippedBeginPosition(gaps);
// Get index i of the according bucket and offset within bucket.
TSourcePos result = 0;
TArrayPos i = 1;
TSourcePos const iEnd = length(gaps._array);
if (unclippedViewPos < static_cast<TGapsPos>(gaps._array[0]))
return 0;
TSourcePos counter = unclippedViewPos - gaps._array[0];
for (; counter > TGapsPos(0) && i < iEnd; ++i)
{
if (counter < gaps._array[i])
break;
if (i % 2) // character bucket
result += gaps._array[i];
counter -= gaps._array[i];
}
if (i % 2) // character bucket.
{
if (i == length(gaps._array)) // We pointing behind the last bucket.
return --result;
result += counter; // If maps into char bucket.
if (counter == gaps._array[i])
--result;
return result;
}
return --result;
}
template <typename TSequence, typename TPosition>
inline typename Position<TSequence>::Type
toSourcePosition(Gaps<TSequence, ArrayGaps> const & gaps, TPosition const clippedViewPos, RightOfViewPos const & /*tag*/)
{
typedef Gaps<TSequence, ArrayGaps> TGaps;
typedef typename Position<TGaps>::Type TGapsPos;
typedef typename TGaps::TArrayPos_ TArrayPos;
typedef typename Position<TSequence>::Type TSourcePos;
// Translate from clipped view position to unclipped view position.
TGapsPos unclippedViewPos = clippedViewPos + clippedBeginPosition(gaps);
// Get index i of the according bucket and offset within bucket.
TSourcePos result = 0;
TArrayPos i = 0;
TSourcePos const iEnd = length(gaps._array);
for (TSourcePos counter = unclippedViewPos; counter > TGapsPos(0) && i < iEnd;)
{
if (counter > gaps._array[i])
{
if (i % 2) // character bucket
result += gaps._array[i];
counter -= gaps._array[i];
i += 1;
}
else if (counter <= gaps._array[i])
{
if (i % 2) // character bucket
result += counter;
counter = 0;
}
}
return result;
}
// ----------------------------------------------------------------------------
// Function toViewPosition()
// ----------------------------------------------------------------------------
// Parameter rightOfGaps moves to the right end of gaps if the character at sourcePosition is followed by a gap in the
// view.
template <typename TSequence, typename TPosition>
inline typename Position<Gaps<TSequence, ArrayGaps> >::Type
toViewPosition(Gaps<TSequence, ArrayGaps> const & gaps, TPosition sourcePosition, bool rightOfGaps = true)
{
typedef Gaps<TSequence, ArrayGaps> TGaps;
typedef typename Position<TGaps>::Type TGapsPosition;
typedef typename TGaps::TArray_ TArray;
typedef typename TGaps::TArrayPos_ TArrayPos;
typedef typename Value<TArray>::Type TArrayValue;
if (sourcePosition == TPosition(0))
return gaps._array[0] - clippedBeginPosition(gaps);
// First, convert to unclipped source position.
TGapsPosition unclippedViewPosition = 0;
TArrayPos i = 0;
for (TArrayValue counter = sourcePosition; counter > TArrayValue(0); ++i)
{
if (i % 2 /*== 1*/) // sequence bucket
{
if (counter > gaps._array[i])
{
unclippedViewPosition += gaps._array[i];
counter -= gaps._array[i];
}
else if (counter < gaps._array[i])
{
unclippedViewPosition += counter;
counter = 0;
}
else // counter == gaps._array[i]
{
unclippedViewPosition += counter;
if (rightOfGaps && i + 2 < length(gaps._array))
unclippedViewPosition += gaps._array[i + 1];
counter = 0;
}
}
else // gaps bucket
{
unclippedViewPosition += gaps._array[i];
}
}
// Return after clipping.
return unclippedViewPosition - clippedBeginPosition(gaps);
}
// ----------------------------------------------------------------------------
// Function insertGaps()
// ----------------------------------------------------------------------------
template <typename TSequence, typename TPosition, typename TCount>
inline void
insertGaps(Gaps<TSequence, ArrayGaps> & gaps, TPosition clippedViewPos, TCount count)
{
typedef Gaps<TSequence, ArrayGaps> TGaps;
typedef typename Position<TGaps>::Type TGapsPosition;
typedef typename TGaps::TArray_ TArray;
typedef typename TGaps::TArrayPos_ TArrayPos;
typedef typename Position<TSequence>::Type TSeqPos;
// Translate from clipped view position to unclipped view position.
TGapsPosition unclippedViewPos = clippedViewPos + clippedBeginPosition(gaps);
// Get index i of the according bucket and offset within bucket.
TArrayPos i = 0;
TSeqPos offset = 0;
for (TSeqPos counter = unclippedViewPos; counter > 0;)
{
SEQAN_ASSERT_LT(i, length(gaps._array));
if (counter > gaps._array[i])
{
counter -= gaps._array[i];
i += 1;
}
else
{
offset = counter;
counter = 0;
}
}
SEQAN_ASSERT_GEQ(gaps._array[i], offset);
// Insert gaps, simple and fast if we are in a gaps bucket, a bit harder
// otherwise.
if (i % 2) // character bucket
{
if (gaps._array[i] > offset) // In the middle of the bucket.
{
TArray arr;
resize(arr, 2, 0);
arr[0] = count;
arr[1] = gaps._array[i] - offset;
gaps._array[i] = offset;
insert(gaps._array, i + 1, arr);
}
else // At the end of the bucket.
{
if (i + 1 < length(gaps._array)) // Not at end of array.
{
gaps._array[i + 1] += count;
}
else // At end of array.
{
resize(gaps._array, length(gaps._array) + 2, 0);
gaps._array[i + 1] = count;
gaps._array[i + 2] = 0;
}
}
}
else // gap bucket
{
gaps._array[i] += count;
}
// Adjust clipping information.
gaps._clippingEndPos += count;
}
// ----------------------------------------------------------------------------
// Function removeGaps()
// ----------------------------------------------------------------------------
template <typename TSequence, typename TPosition, typename TCount>
inline typename Size<Gaps<TSequence, ArrayGaps> >::Type
removeGaps(Gaps<TSequence, ArrayGaps> & gaps, TPosition clippedViewPos, TCount count)
{
typedef Gaps<TSequence, ArrayGaps> TGaps;
typedef typename Position<TGaps>::Type TGapsPosition;
typedef typename TGaps::TArray_ TArray;
typedef typename TGaps::TArrayPos_ TArrayPos;
typedef typename Value<TArray>::Type TArrayValue;
typedef typename Position<TSequence>::Type TSeqPos;
// Translate from clipped view position to unclipped view position.
TGapsPosition pos = clippedViewPos + clippedBeginPosition(gaps);
// Get index i of the according bucket and offset within bucket.
SEQAN_ASSERT_GEQ(length(gaps._array), 2u);
// Start at position 1 if there are no leading gaps.
TArrayPos i = (gaps._array[0] == 0);
TSeqPos offset = 0;
for (TSeqPos counter = pos; counter > 0;)
{
SEQAN_ASSERT_LT(i, length(gaps._array));
if (counter > gaps._array[i])
{
counter -= gaps._array[i];
i += 1;
}
else
{
offset = counter;
counter = 0;
}
}
// Advance into next bucket if at end of current.
if (offset > 0 && offset == gaps._array[i])
{
i += 1;
offset = 0;
}
// If we are inside a non-gap bucket then we cannot remove any gaps.
if (i % 2)
return 0;
// Otherwise, we can remove gaps right of the current position but not
// more than there are.
TSeqPos toRemove = count;
if (toRemove > gaps._array[i] - offset)
toRemove = gaps._array[i] - offset;
gaps._array[i] -= toRemove;
// In some cases, we remove the whole gap and merge the character buckets.
if (gaps._array[i] == TArrayValue(0))
{
// No merging for leading and trailing gap.
if (i != TArrayPos(0) && i != TArrayPos(length(gaps._array) - 1))
{
gaps._array[i - 1] += gaps._array[i + 1];
erase(gaps._array, i, i + 2);
}
}
// Also update the right clipping position.
gaps._clippingEndPos -= toRemove;
// Finally, return number of removed gaps.
return toRemove;
}
// ----------------------------------------------------------------------------
// Function clearGaps()
// ----------------------------------------------------------------------------
template <typename TSequence>
inline void
clearGaps(Gaps<TSequence, ArrayGaps> & gaps)
{
_reinitArrayGaps(gaps);
}
// ----------------------------------------------------------------------------
// Function clear()
// ----------------------------------------------------------------------------
template <typename TSequence>
inline void
clear(Gaps<TSequence, ArrayGaps> & gaps)
{
clear(gaps._source);
clear(gaps._array);
gaps._sourceBeginPos = 0;
gaps._sourceEndPos = 0;
gaps._clippingBeginPos = 0;
gaps._clippingEndPos = 0;
// cannot use clearGaps() here, since that calls value() on _source
// which instates the Holder to Owner; we want it to be empty
}
// ----------------------------------------------------------------------------
// Function isGap()
// ----------------------------------------------------------------------------
template <typename TSequence, typename TPosition>
inline bool
isGap(Gaps<TSequence, ArrayGaps> const & gaps, TPosition clippedViewPos)
{
typedef Gaps<TSequence, ArrayGaps> TGaps;
typedef typename Position<TGaps>::Type TGapsPosition;
typedef typename TGaps::TArrayPos_ TArrayPos;
typedef typename Position<TSequence>::Type TSeqPos;
// Translate from clipped view position to unclipped view position.
TGapsPosition pos = clippedViewPos + clippedBeginPosition(gaps);
// Get index i of the according bucket and offset within bucket.
SEQAN_ASSERT_GEQ(length(gaps._array), 2u);
// Start at position 1 if there are no leading gaps.
TArrayPos i = (gaps._array[0] == 0);
TSeqPos offset = 0;
for (TSeqPos counter = pos; counter > TSeqPos(0);)
{
SEQAN_ASSERT_LT(i, length(gaps._array));
if (counter > gaps._array[i])
{
counter -= gaps._array[i];
i += 1;
}
else
{
offset = counter;
counter = 0;
}
}
// Advance into next bucket if at end of current.
if (offset > TSeqPos(0) && offset == gaps._array[i])
i += 1;
return !(i % 2);
}
// ----------------------------------------------------------------------------
// Function clearClipping()
// ----------------------------------------------------------------------------
template <typename TSequence>
inline void
clearClipping(Gaps<TSequence, ArrayGaps> & gaps)
{
typedef Gaps<TSequence, ArrayGaps> TGaps;
typedef typename TGaps::TArrayPos_ TArrayPos;
gaps._sourceBeginPos = 0;
gaps._sourceEndPos = length(value(gaps._source));
gaps._clippingBeginPos = 0;
gaps._clippingEndPos = 0;
for (TArrayPos i = 0; i < length(gaps._array); ++i)
gaps._clippingEndPos += gaps._array[i];
SEQAN_ASSERT_LEQ(static_cast<int>(gaps._sourceEndPos), static_cast<int>(gaps._clippingEndPos));
}
// ----------------------------------------------------------------------------
// Function setClippedBeginPosition()
// ----------------------------------------------------------------------------
template <typename TSequence, typename TPosition>
inline void
setClippedBeginPosition(Gaps<TSequence, ArrayGaps> & gaps, TPosition unclippedViewPosition)
{
gaps._sourceBeginPos = toSourcePosition(gaps, unclippedViewPosition - clippedBeginPosition(gaps));
gaps._clippingBeginPos = unclippedViewPosition;
}
// ----------------------------------------------------------------------------
// Function setClippedEndPosition()
// ----------------------------------------------------------------------------
template <typename TSequence, typename TPosition>
inline void
setClippedEndPosition(Gaps<TSequence, ArrayGaps> & gaps, TPosition unclippedViewPosition)
{
gaps._sourceEndPos = toSourcePosition(gaps, unclippedViewPosition - clippedBeginPosition(gaps));
//if (isGap(gaps, unclippedViewPosition - clippedBeginPosition(gaps)))
// gaps._sourceEndPos += 1;
gaps._clippingEndPos = unclippedViewPosition;
}
// ----------------------------------------------------------------------------
// Function clippedBeginPosition()
// ----------------------------------------------------------------------------
template <typename TSequence>
inline typename Position<Gaps<TSequence, ArrayGaps> >::Type
clippedBeginPosition(Gaps<TSequence, ArrayGaps> const & gaps)
{
return gaps._clippingBeginPos;
}
// ----------------------------------------------------------------------------
// Function clippedEndPosition()
// ----------------------------------------------------------------------------
template <typename TSequence>
inline typename Position<Gaps<TSequence, ArrayGaps> >::Type
clippedEndPosition(Gaps<TSequence, ArrayGaps> const & gaps)
{
return gaps._clippingEndPos;
}
// ----------------------------------------------------------------------------
// Function setBeginPosition()
// ----------------------------------------------------------------------------
template <typename TSequence, typename TPosition>
inline void
setBeginPosition(Gaps<TSequence, ArrayGaps> & gaps, TPosition sourcePosition)
{
setClippedBeginPosition(gaps, toViewPosition(gaps, sourcePosition) + clippedBeginPosition(gaps));
}
// ----------------------------------------------------------------------------
// Function setEndPosition()
// ----------------------------------------------------------------------------
template <typename TSequence, typename TPosition>
inline void
setEndPosition(Gaps<TSequence, ArrayGaps> & gaps, TPosition sourcePosition)
{
setClippedEndPosition(gaps, toViewPosition(gaps, sourcePosition) + clippedBeginPosition(gaps));
}
// ----------------------------------------------------------------------------
// Function beginPosition()
// ----------------------------------------------------------------------------
// TODO(holtgrew): We would rather like to have the const version only :(
template <typename TSequence>
inline typename Position<TSequence>::Type
beginPosition(Gaps<TSequence, ArrayGaps> const & gaps)
{
return gaps._sourceBeginPos;
}
template <typename TSequence>
inline typename Position<TSequence>::Type
beginPosition(Gaps<TSequence, ArrayGaps> & gaps)
{
return gaps._sourceBeginPos;
}
// ----------------------------------------------------------------------------
// Function endPosition()
// ----------------------------------------------------------------------------
// TODO(holtgrew): We would rather like to have the const version only :(
template <typename TSequence>
inline typename Position<TSequence>::Type
endPosition(Gaps<TSequence, ArrayGaps> const & gaps)
{
return gaps._sourceEndPos;
}
template <typename TSequence>
inline typename Position<ArrayGaps>::Type
endPosition(Gaps<TSequence, ArrayGaps> & gaps)
{
return gaps._sourceEndPos;
}
} // namespace seqan
#endif // #ifndef SEQAN_INCLUDE_SEQAN_ALIGN_GAPS_ARRAY_H_
