1 // ==========================================================================
2 //                 SeqAn - The Library for Sequence Analysis
3 // ==========================================================================
4 // Copyright (c) 2006-2010, Knut Reinert, FU Berlin
5 // All rights reserved.
6 //
7 // Redistribution and use in source and binary forms, with or without
8 // modification, are permitted provided that the following conditions are met:
9 //
10 //     * Redistributions of source code must retain the above copyright
11 //       notice, this list of conditions and the following disclaimer.
12 //     * Redistributions in binary form must reproduce the above copyright
13 //       notice, this list of conditions and the following disclaimer in the
14 //       documentation and/or other materials provided with the distribution.
15 //     * Neither the name of Knut Reinert or the FU Berlin nor the names of
16 //       its contributors may be used to endorse or promote products derived
17 //       from this software without specific prior written permission.
18 //
19 // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
20 // AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
21 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
22 // ARE DISCLAIMED. IN NO EVENT SHALL KNUT REINERT OR THE FU BERLIN BE LIABLE
23 // FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
24 // DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
25 // SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
26 // CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
27 // LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
28 // OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
29 // DAMAGE.
30 //
31 // ==========================================================================
32 //  Author: Andreas Gogol-Doering <andreas.doering@mdc-berlin.de>
33 // ==========================================================================
34 
35 #ifndef SEQAN_HEADER_GAPS_ARRAY_H
36 #define SEQAN_HEADER_GAPS_ARRAY_H
37 
38 namespace SEQAN_NAMESPACE_MAIN
39 {
40 
41 //////////////////////////////////////////////////////////////////////////////
42 // Tag
43 
44 struct ArrayGaps;
45 
46 
47 //////////////////////////////////////////////////////////////////////////////
48 // Gaps - ArrayGaps Spec
49 //////////////////////////////////////////////////////////////////////////////
50 
51 /**
52 .Spec.ArrayGaps:
53 ..cat:Alignments
54 ..general:Class.Gaps
55 ..summary:Stores gaps sizes in an array.
56 ..signature:Gaps<TSource, ArrayGaps >
57 ..param.TSource:Type of the ungapped sequence.
58 ...metafunction:Metafunction.Source
59 ..include:seqan/align.h
60 */
61 template <typename TSource>
62 class Gaps<TSource, ArrayGaps>
63 {
64 public:
65 	typedef typename Size<Gaps>::Type TSize;
66 	typedef String<TSize> TArr;
67 	typedef typename Position<TSource>::Type TSourcePosition;
68 	typedef typename Position<Gaps>::Type TViewPosition;
69 
70 public:
71 	String<TSize> data_arr; //a list of gap and non-gap region lengths
72 	TViewPosition data_end_position;
73 	TViewPosition data_unclipped_end_position;
74 
75 	Holder<TSource> data_source;
76 	TSourcePosition clipped_source_begin;
77 	TSourcePosition clipped_source_end;
78 
79 public:
Gaps()80 	Gaps():
81 		data_unclipped_end_position(0),
82 		clipped_source_begin(0),
83 		clipped_source_end(0)
84 	{
85 SEQAN_CHECKPOINT
86 	}
Gaps(TSize _size)87 	Gaps(TSize _size):
88 		data_unclipped_end_position(0),
89 		clipped_source_begin(0),
90 		clipped_source_end(0)
91 	{
92 		_initToResize(*this, _size);
93 	}
Gaps(TSource & source_)94 	Gaps(TSource & source_):
95 		data_unclipped_end_position(0),
96 		data_source(source_),
97 		clipped_source_begin(beginPosition(source_)),
98 		clipped_source_end(endPosition(source_))
99 	{
100 SEQAN_CHECKPOINT
101 		_initToResize(*this, length(source_));
102 	}
103 
104 	template <typename TSource2>
Gaps(TSource2 const & source_)105 	Gaps(TSource2 const & source_):
106 		data_unclipped_end_position(0),
107 		clipped_source_begin(0),
108 		clipped_source_end(length(source_))
109 		//clipped_source_begin(beginPosition(source_)),
110 		//clipped_source_end(endPosition(source_))
111 	{
112 SEQAN_CHECKPOINT
113 		data_source = source_;
114 		_initToResize(*this, length(source_));
115 	}
116 
Gaps(Gaps const & other_)117 	Gaps(Gaps const & other_):
118 		data_arr(other_.data_arr),
119 		data_end_position(other_.data_end_position),
120 		data_unclipped_end_position(other_.data_unclipped_end_position),
121 //		data_source(value(other_.data_source)), //variant: setValue => Align benutzen gemeinsame Sources
122 		data_source(other_.data_source),		//variant: assignValue => Align kopieren Sources
123 		clipped_source_begin(other_.clipped_source_begin),
124 		clipped_source_end(other_.clipped_source_end)
125 	{
126 SEQAN_CHECKPOINT
127 	}
128 	Gaps & operator = (Gaps const & other_)
129 	{
130 SEQAN_CHECKPOINT
131 		data_arr = other_.data_arr;
132 		data_end_position = other_.data_end_position;
133 		data_unclipped_end_position = other_.data_unclipped_end_position,
134 		setValue(data_source, source(other_));
135 		clipped_source_begin = other_.clipped_source_begin;
136 		clipped_source_end = other_.clipped_source_end;
137 		return *this;
138 	}
139 
~Gaps()140 	~Gaps()
141 	{
142 SEQAN_CHECKPOINT
143 	}
144 
145 	inline typename Reference<Gaps>::Type
146 	operator[](TViewPosition view_pos)
147 	{
148 SEQAN_CHECKPOINT
149 		return value(*this, view_pos);
150 	}
151 	inline typename Reference<Gaps const>::Type
152 	operator[](TViewPosition view_pos) const
153 	{
154 SEQAN_CHECKPOINT
155 		return value(*this, view_pos);
156 	}
157 //____________________________________________________________________________
158 };
159 
160 //////////////////////////////////////////////////////////////////////////////
161 
162 template <typename TSource>
163 inline String<typename Size<Gaps<TSource, ArrayGaps> >::Type> &
_dataArr(Gaps<TSource,ArrayGaps> & me)164 _dataArr(Gaps<TSource, ArrayGaps> & me)
165 {
166 SEQAN_CHECKPOINT
167 	return me.data_arr;
168 }
169 template <typename TSource>
170 inline String<typename Size<Gaps<TSource, ArrayGaps> >::Type> const &
_dataArr(Gaps<TSource,ArrayGaps> const & me)171 _dataArr(Gaps<TSource, ArrayGaps> const & me)
172 {
173 SEQAN_CHECKPOINT
174 	return me.data_arr;
175 }
176 
177 //////////////////////////////////////////////////////////////////////////////
178 
179 template <typename TSource>
180 inline Holder<TSource> &
_dataSource(Gaps<TSource,ArrayGaps> & me)181 _dataSource(Gaps<TSource, ArrayGaps> & me)
182 {
183 SEQAN_CHECKPOINT
184 	return me.data_source;
185 }
186 template <typename TSource>
187 inline Holder<TSource> const &
_dataSource(Gaps<TSource,ArrayGaps> const & me)188 _dataSource(Gaps<TSource, ArrayGaps> const & me)
189 {
190 SEQAN_CHECKPOINT
191 	return me.data_source;
192 }
193 
194 //////////////////////////////////////////////////////////////////////////////
195 
196 template <typename TSource>
197 inline typename Position<Gaps<TSource, ArrayGaps> >::Type
endPosition(Gaps<TSource,ArrayGaps> & me)198 endPosition(Gaps<TSource, ArrayGaps> & me)
199 {
200 SEQAN_CHECKPOINT
201 	return me.data_end_position;
202 }
203 template <typename TSource>
204 inline typename Position<Gaps<TSource, ArrayGaps> >::Type
endPosition(Gaps<TSource,ArrayGaps> const & me)205 endPosition(Gaps<TSource, ArrayGaps> const & me)
206 {
207 SEQAN_CHECKPOINT
208 	return me.data_end_position;
209 }
210 
211 //////////////////////////////////////////////////////////////////////////////
212 
213 template <typename TSource>
214 inline typename Position<Gaps<TSource, ArrayGaps> >::Type
_getTrailingGaps(Gaps<TSource,ArrayGaps> const & me)215 _getTrailingGaps(Gaps<TSource, ArrayGaps> const & me)
216 {
217 SEQAN_CHECKPOINT
218 	return me.data_unclipped_end_position;
219 }
220 
221 template <typename TSource, typename TSize>
222 inline void
_setTrailingGaps(Gaps<TSource,ArrayGaps> & me,TSize const & size)223 _setTrailingGaps(Gaps<TSource, ArrayGaps> & me, TSize const & size)
224 {
225 	TSize zero = 0;
226 	if (size >= zero)
227 	{
228 		me.data_unclipped_end_position = size;
229 	}
230 	else
231 	{
232 		me.data_unclipped_end_position = 0;
233 	}
234 }
235 
236 //////////////////////////////////////////////////////////////////////////////
237 /**
238  * Returns the length of the gaps sequence including the trailing gaps
239  */
240 template <typename TSource>
241 inline typename Position<Gaps<TSource, ArrayGaps> >::Type
_unclippedLength(Gaps<TSource,ArrayGaps> & me)242 _unclippedLength(Gaps<TSource, ArrayGaps> & me)
243 {
244 SEQAN_CHECKPOINT
245 	return me.data_end_position + me.data_unclipped_end_position;
246 }
247 
248 //////////////////////////////////////////////////////////////////////////////
249 
250 template <typename TSource, typename TPosition>
251 inline void
_setEndPosition(Gaps<TSource,ArrayGaps> & me,TPosition _pos)252 _setEndPosition(Gaps<TSource, ArrayGaps> & me, TPosition _pos)
253 {
254 SEQAN_CHECKPOINT
255 	me.data_end_position = _pos;
256 }
257 
258 //////////////////////////////////////////////////////////////////////////////
259 
260 /**
261 .Function.clippedBeginPosition:
262 ..summary:Begin position of the source segment.
263 ..cat:Alignments
264 ..signature:clippedBeginPosition(object)
265 ..param.object:An object that has a source
266 ...type:Class.Gaps
267 ..returns:The position of the first item in @Function.source.source(object)@ that is used in object.
268 ..see:Function.begin
269 ..see:Function.source
270 ..see:Function.clippedEndPosition
271 ..include:seqan/align.h
272 */
273 template <typename TSource>
274 inline typename Position<TSource>::Type
clippedBeginPosition(Gaps<TSource,ArrayGaps> const & me)275 clippedBeginPosition(Gaps<TSource, ArrayGaps> const & me)
276 {
277 SEQAN_CHECKPOINT
278 	return me.clipped_source_begin;
279 }
280 
281 //////////////////////////////////////////////////////////////////////////////
282 
283 template <typename TSource, typename TSourcePosition>
284 inline void
_setClippedBeginPosition(Gaps<TSource,ArrayGaps> & me,TSourcePosition _pos)285 _setClippedBeginPosition(Gaps<TSource, ArrayGaps> & me, TSourcePosition _pos)
286 {
287 SEQAN_CHECKPOINT
288 	me.clipped_source_begin = _pos;
289 }
290 
291 //////////////////////////////////////////////////////////////////////////////
292 
293 /**
294 .Function.clippedEndPosition:
295 ..summary:Position of the end of the source segment.
296 ..cat:Alignments
297 ..signature:clippedEndPosition(object)
298 ..param.object:An object that has a source
299 ...type:Class.Gaps
300 ..returns:The position behind the last element of the source segment.
301 ..see:Function.end
302 ..see:Function.sourceEnd
303 ..see:Function.clippedBeginPosition
304 ..include:seqan/align.h
305 */
306 template <typename TSource>
307 inline typename Position<TSource>::Type
clippedEndPosition(Gaps<TSource,ArrayGaps> const & me)308 clippedEndPosition(Gaps<TSource, ArrayGaps> const & me)
309 {
310 SEQAN_CHECKPOINT
311 	return me.clipped_source_end;
312 }
313 //////////////////////////////////////////////////////////////////////////////
314 
315 template <typename TSource, typename TSourcePosition>
316 inline void
_setClippedEndPosition(Gaps<TSource,ArrayGaps> & me,TSourcePosition _pos)317 _setClippedEndPosition(Gaps<TSource, ArrayGaps> & me, TSourcePosition _pos)
318 {
319 SEQAN_CHECKPOINT
320 	me.clipped_source_end = _pos;
321 }
322 
323 //////////////////////////////////////////////////////////////////////////////
324 
325 template <typename TSource, typename TSize>
326 inline void
_initToResize(Gaps<TSource,ArrayGaps> & me,TSize _size)327 _initToResize(Gaps<TSource, ArrayGaps> & me,
328 				TSize _size)
329 {
330 SEQAN_CHECKPOINT
331 	resize(_dataArr(me), 2);
332 	_dataArr(me)[0] = 0;
333 	_dataArr(me)[1] = _size;
334 	_setEndPosition(me, _size);
335 	_setTrailingGaps(me, 0);
336 }
337 
338 //////////////////////////////////////////////////////////////////////////////
339 
340 template <typename TSource>
341 inline void
clearGaps(Gaps<TSource,ArrayGaps> & me)342 clearGaps(Gaps<TSource, ArrayGaps> & me)
343 {
344 SEQAN_CHECKPOINT
345 	_initToResize(me, clippedEndPosition(me) - clippedBeginPosition(me));
346 }
347 
348 //////////////////////////////////////////////////////////////////////////////
349 
350 template <typename TSource>
351 inline void
clear(Gaps<TSource,ArrayGaps> & me)352 clear(Gaps<TSource, ArrayGaps> & me)
353 {
354 SEQAN_CHECKPOINT
355 	_initToResize(me, 0);
356 	_setClippedBeginPosition(me, 0);
357 	_setClippedEndPosition(me, 0);
358 }
359 
360 //////////////////////////////////////////////////////////////////////////////
361 
362 template <typename TSource>
363 inline void
clearClipping(Gaps<TSource,ArrayGaps> & me)364 clearClipping(Gaps<TSource, ArrayGaps> & me)
365 {
366 SEQAN_CHECKPOINT
367 	_setEndPosition(me, length(source(me)));
368 	_setTrailingGaps(me, 0);
369 	_setClippedBeginPosition(me, 0);
370 	_setClippedEndPosition(me, length(source(me)));
371 }
372 
373 //////////////////////////////////////////////////////////////////////////////
374 
375 /**
376 .Function.toViewPosition:
377 ..summary:Transforms source to view position.
378 ..cat:Alignments
379 ..signature:toViewPosition(gaps, pos)
380 ..param.gap:A Gaps object, e.g. a row in the alignment.
381 ...type:Class.Gaps
382 ..param.pos:Position in the original sequence to get the view position for.
383 ..returns:The position in the view/gaps position.
384 ..remarks:If $gap$ is a clipped alignment row, gaps in the clipped part will influence the result. The position $pos$ is counted from the unclipped begin position and must be greater or equal the clipped begin position of $gap$.
385 ..see:Function.toSourcePosition
386 ..include:seqan/align.h
387 */
388 template <typename TSource>
389 inline typename Position< Gaps<TSource, ArrayGaps> >::Type
toViewPosition(Gaps<TSource,ArrayGaps> const & gaps,typename Position<TSource>::Type pos)390 toViewPosition(Gaps<TSource, ArrayGaps> const & gaps,
391 			   typename Position<TSource>::Type pos)
392 {
393 SEQAN_CHECKPOINT
394 
395 	SEQAN_ASSERT_GEQ(pos, clippedBeginPosition(gaps));
396 	pos -= clippedBeginPosition(gaps);
397 
398 	typedef Gaps<TSource, ArrayGaps> const TGaps;
399 	typedef typename Size<TGaps>::Type TSize;
400 	typedef String<TSize> const TArr;
401 
402 	typename Iterator<TArr, Standard>::Type arr_begin = begin(_dataArr(gaps));
403 	typename Iterator<TArr, Standard>::Type arr_end = end(_dataArr(gaps));
404 	typename Position<TGaps>::Type view_pos = pos;
405 	typename Position<TSource>::Type source_pos = pos;
406 
407 	while (true)
408 	{
409 		if (arr_begin == arr_end) return view_pos;
410 		TSize step = *arr_begin;
411 		view_pos += step;
412 
413 		++arr_begin;
414 
415 		step = *arr_begin;
416 		if (source_pos < step) return view_pos;
417 		source_pos -= step;
418 
419 		++arr_begin;
420 	}
421 }
422 
423 //____________________________________________________________________________
424 
425 /**
426 .Function.toSourcePosition:
427 ..summary:Transforms view to source position, if the view position is a gap, the original position of the next non-gap entry is returned.
428 ..cat:Alignments
429 ..signature:toSourcePosition(gaps, pos)
430 ..param.gap:A Gaps object, e.g. a row in the alignment.
431 ...type:Class.Gaps
432 ..param.pos:Position in the view sequence (this includes gaps) to get the original position for.
433 ..returns:The position in the source sequence.
434 ..remarks:If $gap$ is a clipped alignment row, gaps in the clipped part will influence the result. The position $pos$ is counted from the unclipped begin position.
435 ..see:Function.toViewPosition
436 ..include:seqan/align.h
437 */
438 template <typename TSource>
439 inline typename Position<TSource>::Type
toSourcePosition(Gaps<TSource,ArrayGaps> const & gaps,typename Position<Gaps<TSource,ArrayGaps>>::Type pos)440 toSourcePosition(Gaps<TSource, ArrayGaps> const & gaps,
441 				 typename Position< Gaps<TSource, ArrayGaps> >::Type pos)
442 {
443 SEQAN_CHECKPOINT
444 	typedef Gaps<TSource, ArrayGaps> const TGaps;
445 	typedef typename Size<TGaps>::Type TSize;
446 	typedef String<TSize> const TArr;
447 
448 	typename Iterator<TArr, Standard>::Type arr_begin = begin(_dataArr(gaps));
449 	typename Iterator<TArr, Standard>::Type arr_end = end(_dataArr(gaps));
450 	typename Position<TSource>::Type source_pos = clippedBeginPosition(gaps);
451 	typename Position<TGaps>::Type view_pos = pos;
452 
453 	while (true)
454 	{
455 		if (arr_begin == arr_end) return source_pos;
456 
457 		TSize step = *arr_begin;
458 		if (view_pos <= step) return source_pos;
459 		view_pos -= step;
460 
461 		++arr_begin;
462 		step = *arr_begin;
463 		if (view_pos <= step ) return source_pos + view_pos;
464 		source_pos += step;
465 		view_pos -= step;
466 
467 		++arr_begin;
468 	}
469 }
470 
471 //////////////////////////////////////////////////////////////////////////////
472 // returns an iterator to view-Position
473 
474 template <typename TGaps, typename TPosition>
475 inline typename Iterator<TGaps, Standard>::Type
_iteratorGapsArray(TGaps & gaps,TPosition view_position)476 _iteratorGapsArray(TGaps & gaps,
477 					 TPosition view_position)
478 {
479 	typedef typename Size<TGaps>::Type TSize;
480 	typedef typename TGaps::TArr const TArr;
481 	typedef typename Source<TGaps>::Type TSource;
482 	typedef typename Iterator<TArr, Standard>::Type TArrIterator;
483 	typedef typename Value<TArrIterator>::Type TArrIteratorType;
484 
485  	TArrIterator arr_begin = begin(_dataArr(gaps), Standard());
486 	TArrIterator arr_end = end(_dataArr(gaps), Standard());
487 	typename Position<TGaps>::Type block = 0;
488 
489 	TArrIteratorType view_pos = view_position;
490 
491 	if (emptySource(gaps))
492 	{
493 		for (; arr_begin != arr_end; ++arr_begin)
494 		{
495 			if (view_pos < *arr_begin) break;
496 			++block;
497 			view_pos -= *arr_begin;
498 		}
499 		return typename Iterator<TGaps, Standard>::Type(
500 			gaps,
501 			block,
502 			view_pos);
503 	}
504 	else
505 	{
506 		typename Position<TSource>::Type source_pos = clippedBeginPosition(gaps);
507 
508 		while (true)
509 		{
510 			if ((arr_begin == arr_end) || (view_pos < *arr_begin))
511 			{
512 				break;
513 			}
514 			++block;
515 			view_pos -= *arr_begin;
516 			++arr_begin;
517 
518 			if (view_pos < *arr_begin)
519 			{
520 				source_pos += view_pos;
521 				break;
522 			}
523 			++block;
524 			view_pos -= *arr_begin;
525 			source_pos += *arr_begin;
526 			++arr_begin;
527 		}
528 
529 		return typename Iterator<TGaps, Standard>::Type(
530 			gaps,
531             iter(source(gaps), source_pos),
532 			block,
533 			view_pos);
534 	}
535 
536 }
537 
538 template <typename TSource, typename TPosition, typename TTag>
539 inline typename Iterator<Gaps<TSource, ArrayGaps>, Tag<TTag> const>::Type
iter(Gaps<TSource,ArrayGaps> & gaps,TPosition view_pos,Tag<TTag> const)540 iter(Gaps<TSource, ArrayGaps> & gaps,
541 	 TPosition view_pos,
542 	 Tag<TTag> const)
543 {
544 SEQAN_CHECKPOINT
545 	return _iteratorGapsArray(gaps, view_pos);
546 }
547 template <typename TSource, typename TPosition, typename TTag>
548 inline typename Iterator<Gaps<TSource, ArrayGaps> const, Tag<TTag> const>::Type
iter(Gaps<TSource,ArrayGaps> const & gaps,TPosition view_pos,Tag<TTag> const)549 iter(Gaps<TSource, ArrayGaps> const & gaps,
550 	 TPosition view_pos,
551 	 Tag<TTag> const)
552 {
553 SEQAN_CHECKPOINT
554 	return _iteratorGapsArray<Gaps<TSource, ArrayGaps> const>(gaps, view_pos);
555 }
556 
557 //////////////////////////////////////////////////////////////////////////////
558 
559 template <typename TSource>
560 inline typename Position< Gaps<TSource, ArrayGaps> >::Type
beginPosition(Gaps<TSource,ArrayGaps> & gaps)561 beginPosition(Gaps<TSource, ArrayGaps> & gaps)
562 {
563 SEQAN_CHECKPOINT
564 	return value(_dataArr(gaps), 0);
565 }
566 template <typename TSource>
567 inline typename Position< Gaps<TSource, ArrayGaps> const>::Type
beginPosition(Gaps<TSource,ArrayGaps> const & gaps)568 beginPosition(Gaps<TSource, ArrayGaps> const & gaps)
569 {
570 SEQAN_CHECKPOINT
571 	return value(_dataArr(gaps), 0);
572 }
573 
574 //////////////////////////////////////////////////////////////////////////////
575 //////////////////////////////////////////////////////////////////////////////
576 
577 
578 template <typename TSource, typename TPosition>
579 inline void
setBeginPosition(Gaps<TSource,ArrayGaps> & me,TPosition view_position)580 setBeginPosition(Gaps<TSource, ArrayGaps> & me,
581 				 TPosition view_position)
582 {
583 SEQAN_CHECKPOINT
584 	TPosition old_pos = beginPosition(me);
585 	if (length(_dataArr(me)) == 0)
586 	{
587 		_initToResize(me, length(source(me)));
588 	}
589 	_dataArr(me)[0] = view_position;
590 	_setEndPosition(me, endPosition(me) + view_position - old_pos);
591 }
592 
593 //////////////////////////////////////////////////////////////////////////////
594 
595 template <typename TSource, typename TPosition>
596 inline void
setClippedBeginPosition(Gaps<TSource,ArrayGaps> & me,TPosition source_position)597 setClippedBeginPosition(Gaps<TSource, ArrayGaps> & me,
598 					   TPosition source_position)
599 {
600 	SEQAN_ASSERT(length(_dataArr(me)))
601 
602 	typedef Gaps<TSource, ArrayGaps> TGaps;
603 	typedef typename Position<TGaps>::Type TViewPosition;
604 	typedef typename TGaps::TArr TArr;
605 
606 	TPosition old_clipped_begin_pos = clippedBeginPosition(me);
607 	if (old_clipped_begin_pos == source_position) return;
608 	else if (source_position < old_clipped_begin_pos)
609 	{
610 SEQAN_CHECKPOINT
611 		TViewPosition delta = old_clipped_begin_pos - source_position;
612 		TViewPosition view_begin_pos = beginPosition(me);
613 		if (view_begin_pos <= delta)
614 		{
615 			setBeginPosition(me, 0);
616 		}
617 		else
618 		{
619 			setBeginPosition(me, view_begin_pos - delta);
620 		}
621 		_setEndPosition(me, endPosition(me) + delta);
622 		_dataArr(me)[1] += delta;
623 		_setClippedBeginPosition(me, source_position);
624 	}
625 	else //(source_position > old_clipped_begin_pos)
626 	{
627 SEQAN_CHECKPOINT
628 		TViewPosition view_pos = toViewPosition(me, source_position);
629 		TViewPosition source_pos_left = source_position - old_clipped_begin_pos;
630 		TViewPosition gaps_count = 0;
631 
632 		typename Iterator<TArr, Standard>::Type it_arr_begin = begin(_dataArr(me));
633 		typename Iterator<TArr, Standard>::Type it_arr_end = end(_dataArr(me));
634 		typename Iterator<TArr, Standard>::Type it_arr = it_arr_begin;
635 
636 		while (it_arr != it_arr_end)
637 		{
638 			gaps_count += *it_arr;
639 			++it_arr;
640 			if (*it_arr > source_pos_left)
641 			{
642 				*it_arr_begin = view_pos;
643 				*it_arr -= source_pos_left;
644 				if (it_arr - it_arr_begin > 1)
645 				{
646 					replace(_dataArr(me), 1, it_arr - it_arr_begin, "");
647 				}
648 				_setClippedBeginPosition(me, source_position);
649 				return;
650 			}
651 			gaps_count += *it_arr;
652 			source_pos_left -= *it_arr;
653 			++it_arr;
654 		}
655 		//alignment is empty
656 		clear(me);
657 	}
658 }
659 
660 //////////////////////////////////////////////////////////////////////////////
661 
662 template <typename TSource, typename TPosition>
663 inline void
setClippedEndPosition(Gaps<TSource,ArrayGaps> & me,TPosition source_position)664 setClippedEndPosition(Gaps<TSource, ArrayGaps> & me,
665 					 TPosition source_position)
666 {
667 	typedef Gaps<TSource, ArrayGaps> TGaps;
668 	typedef typename Position<TGaps>::Type TViewPosition;
669 	typedef typename TGaps::TArr TArr;
670 
671 	TArr arr = _dataArr(me);
672 	SEQAN_ASSERT(length(arr));
673 
674 	TPosition old_end_begin_pos = clippedEndPosition(me);
675 	if (old_end_begin_pos == source_position) return;
676 	else if (source_position < old_end_begin_pos)
677 	{
678 SEQAN_CHECKPOINT
679 		typename Iterator<TArr, Standard>::Type it_arr_begin = begin(_dataArr(me));
680 		typename Iterator<TArr, Standard>::Type it_arr_end = end(_dataArr(me));
681 		typename Iterator<TArr, Standard>::Type it_arr = it_arr_begin;
682 		TViewPosition end_pos = 0;
683 		TViewPosition chars_to_scan = source_position - clippedBeginPosition(me);
684 
685 		while (it_arr != it_arr_end)
686 		{
687 			end_pos += *it_arr;
688 			++it_arr;
689 			if (*it_arr >= chars_to_scan)
690 			{
691 				resize(_dataArr(me), it_arr - it_arr_begin + 1);
692 				*it_arr = chars_to_scan;
693 				_setEndPosition(me, end_pos + chars_to_scan);
694 				break;
695 			}
696 			end_pos += *it_arr;
697 			chars_to_scan -= *it_arr;
698 			++it_arr;
699 		}
700 	}
701 	else //(source_position > old_end_begin_pos)
702 	{
703 SEQAN_CHECKPOINT
704 		*(end(_dataArr(me)) - 1) += (source_position - old_end_begin_pos);
705 		_setEndPosition(me, endPosition(me) + source_position - old_end_begin_pos);
706 	}
707 	_setClippedEndPosition(me, source_position);
708 }
709 
710 //////////////////////////////////////////////////////////////////////////////
711 
712 template <typename TSource>
713 inline typename Size<TSource>::Type
sourceLength(Gaps<TSource,ArrayGaps> & me)714 sourceLength(Gaps<TSource, ArrayGaps> & me)
715 {
716 SEQAN_CHECKPOINT
717 	return clippedEndPosition(me) - clippedBeginPosition(me);
718 }
719 
720 //////////////////////////////////////////////////////////////////////////////
721 //////////////////////////////////////////////////////////////////////////////
722 //////////////////////////////////////////////////////////////////////////////
723 // Gaps Iterator
724 //////////////////////////////////////////////////////////////////////////////
725 
726 
727 template <typename TGaps>
728 class Iter<TGaps, GapsIterator<ArrayGaps> >
729 {
730 public:
731 	typedef typename Size<TGaps>::Type TGapsSize;
732 	typedef String<TGapsSize> TArr;
733 	typedef typename Position<TArr>::Type TArrPosition;
734 	typedef typename Source<Iter>::Type TSourceIterator;
735 
736 	TGaps * data_container;								//the gaps object
737 	mutable TSourceIterator data_source;
738 	mutable TArrPosition data_block;					//block in array of container
739 	mutable TGapsSize data_sub;							//position block
740 
741 public:
Iter()742 	Iter()
743 	{
744 SEQAN_CHECKPOINT
745 	}
Iter(Iter const & other_)746 	Iter(Iter const & other_):
747 		data_container(other_.data_container),
748 		data_source(other_.data_source),
749 		data_block(other_.data_block),
750 		data_sub(other_.data_sub)
751 	{
752 SEQAN_CHECKPOINT
753 	}
754 	Iter(TGaps & container_,
755 			TSourceIterator source_iterator,
756 			TArrPosition block_ = 0,
757 			TGapsSize sub_ = 0):
758 		data_container(&container_),
759 		data_source(source_iterator),
760 		data_block(block_),
761 		data_sub(sub_)
762 	{
763 SEQAN_CHECKPOINT
764 	}
765 
766 	Iter(TGaps & container_,
767 			TArrPosition block_,
768 			TGapsSize sub_ = 0):
769 		data_container(&container_),
770 		data_block(block_),
771 		data_sub(sub_)
772 	{
773 SEQAN_CHECKPOINT
774 	}
775 
~Iter()776 	~Iter()
777 	{
778 SEQAN_CHECKPOINT
779 	}
780 
781 	Iter const & operator = (Iter const & other_)
782 	{
783 SEQAN_CHECKPOINT
784 		data_container = other_.data_container;
785 		data_block = other_.data_block;
786 		data_sub = other_.data_sub;
787 		data_source = other_.data_source;
788 		return *this;
789 	}
790 };
791 
792 //////////////////////////////////////////////////////////////////////////////
793 
794 template <typename TGaps>
795 inline typename GetValue< Iter<TGaps, GapsIterator<ArrayGaps> > >::Type
getValue(Iter<TGaps,GapsIterator<ArrayGaps>> & me)796 getValue(Iter<TGaps, GapsIterator<ArrayGaps> > & me)
797 {
798 SEQAN_CHECKPOINT
799 	typedef typename Value<Iter<TGaps, GapsIterator<ArrayGaps> > >::Type TValue;
800 	if (isGap(me)) return gapValue<TValue>();
801 	else return getValue(source(me));
802 }
803 template <typename TGaps>
804 inline typename GetValue< Iter<TGaps, GapsIterator<ArrayGaps> > const>::Type
getValue(Iter<TGaps,GapsIterator<ArrayGaps>> const & me)805 getValue(Iter<TGaps, GapsIterator<ArrayGaps> > const & me)
806 {
807 SEQAN_CHECKPOINT
808 	typedef typename Value<Iter<TGaps, GapsIterator<ArrayGaps> > const>::Type TValue;
809 	if (isGap(me)) return gapValue<TValue>();
810 	else return getValue(source(me));
811 }
812 
813 //____________________________________________________________________________
814 
815 template <typename TGaps>
816 inline bool
isGap(Iter<TGaps,GapsIterator<ArrayGaps>> const & me)817 isGap(Iter<TGaps, GapsIterator<ArrayGaps> > const & me)
818 {
819 SEQAN_CHECKPOINT
820 	return !(me.data_block & 1);
821 }
822 
823 //____________________________________________________________________________
824 
825 template <typename TGaps>
826 inline typename Size<TGaps>::Type
countGaps(Iter<TGaps,GapsIterator<ArrayGaps>> const & me)827 countGaps(Iter<TGaps, GapsIterator<ArrayGaps> > const & me)
828 {
829 	if (isGap(me))
830 	{
831 		typedef typename Size<TGaps>::Type TGapsSize;
832 		typedef String<TGapsSize> const TArr;
833 
834 		TArr & arr = _dataArr(container(me));
835 		typename Position<TArr>::Type pos = me.data_block;
836 		typename Position<TArr>::Type block_pos = me.data_sub;
837 
838 		if (length(arr) == pos)
839 		{//counting trailing gaps
840 			if (block_pos <= _getTrailingGaps(container(me)))
841 			{//is in range of trailing gaps
842 SEQAN_CHECKPOINT
843 				return _getTrailingGaps(container(me)) - block_pos;
844 			}
845 			else
846 			{//no trailing gaps here
847 SEQAN_CHECKPOINT
848 				return 0;
849 			}
850 		}
851 		else
852 		{
853 SEQAN_CHECKPOINT
854 			return arr[pos] - me.data_sub;
855 		}
856 	}
857 	else
858 	{
859 SEQAN_CHECKPOINT
860 		return 0;
861 	}
862 }
863 
864 template <typename TGaps>
865 inline typename Size<TGaps>::Type
countCharacters(Iter<TGaps,GapsIterator<ArrayGaps>> const & me)866 countCharacters(Iter<TGaps, GapsIterator<ArrayGaps> > const & me)
867 {
868 	if (isGap(me))
869 	{// only count chars
870 SEQAN_CHECKPOINT
871 		return 0;
872 	}
873 	else
874 	{//count chars
875 SEQAN_CHECKPOINT
876 		typedef typename Size<TGaps>::Type TGapsSize;
877 		typedef String<TGapsSize> const TArr;
878 
879 		TArr & arr = _dataArr(container(me));
880 		typename Position<TArr>::Type pos = me.data_block;
881 		typename Position<TArr>::Type block_pos = me.data_sub;
882 
883 		return arr[pos] - block_pos;
884 	}
885 }
886 
887 //____________________________________________________________________________
888 
889 template <typename T>
890 inline void
_goNextArrayGapsIterator(T const & me)891 _goNextArrayGapsIterator(T const & me)
892 {
893 	if (!isGap(me))
894 	{
895 SEQAN_CHECKPOINT
896 		goNext(me.data_source);
897 	}
898 
899 	++me.data_sub;
900 	if ((me.data_block < length(_dataArr(container(me)))) && (me.data_sub >= _dataArr(container(me))[me.data_block]))
901 	{
902 SEQAN_CHECKPOINT
903 		++me.data_block;
904 		me.data_sub = 0;
905 	}
906 }
907 template <typename TGaps>
908 inline void
goNext(Iter<TGaps,GapsIterator<ArrayGaps>> & me)909 goNext(Iter<TGaps, GapsIterator<ArrayGaps> > & me)
910 {
911 	_goNextArrayGapsIterator(me);
912 }
913 template <typename TGaps>
914 inline void
goNext(Iter<TGaps,GapsIterator<ArrayGaps>> const & me)915 goNext(Iter<TGaps, GapsIterator<ArrayGaps> > const & me)
916 {
917 	_goNextArrayGapsIterator(me);
918 }
919 
920 //____________________________________________________________________________
921 
922 template <typename T>
923 inline void
_goPreviousArrayGapsIterator(T const & me)924 _goPreviousArrayGapsIterator(T const & me)
925 {
926 	if (me.data_sub > 0)
927 	{
928 SEQAN_CHECKPOINT
929 		--me.data_sub;
930 	}
931 	else
932 	{
933 		if (me.data_block > 0)
934 		{
935 SEQAN_CHECKPOINT
936 			--me.data_block;
937 			me.data_sub = _dataArr(container(me))[me.data_block] - 1;
938 		}
939 	}
940 	if (!isGap(me))
941 	{
942 SEQAN_CHECKPOINT
943 		goPrevious(me.data_source);
944 	}
945 }
946 template <typename TGaps>
947 inline void
goPrevious(Iter<TGaps,GapsIterator<ArrayGaps>> & me)948 goPrevious(Iter<TGaps, GapsIterator<ArrayGaps> > & me)
949 {
950 	_goPreviousArrayGapsIterator(me);
951 }
952 template <typename TGaps>
953 inline void
goPrevious(Iter<TGaps,GapsIterator<ArrayGaps>> const & me)954 goPrevious(Iter<TGaps, GapsIterator<ArrayGaps> > const & me)
955 {
956 	_goPreviousArrayGapsIterator(me);
957 }
958 //____________________________________________________________________________
959 
960 template <typename TGaps>
961 inline bool
atBegin(Iter<TGaps,GapsIterator<ArrayGaps>> const & me)962 atBegin(Iter<TGaps, GapsIterator<ArrayGaps> > const & me)
963 {
964 SEQAN_CHECKPOINT
965 	return ((me.data_block == 1) && (me.data_sub == 0));
966 }
967 
968 //____________________________________________________________________________
969 
970 template <typename TGaps>
971 inline bool
atEnd(Iter<TGaps,GapsIterator<ArrayGaps>> const & me)972 atEnd(Iter<TGaps, GapsIterator<ArrayGaps> > const & me)
973 {
974 SEQAN_CHECKPOINT
975 	return ((me.data_block == length(_dataArr(container(me)))) && (me.data_sub == 0));
976 }
977 
978 //____________________________________________________________________________
979 
980 template <typename TGaps, typename TCount>
981 inline void
insertGaps(Iter<TGaps,GapsIterator<ArrayGaps>> const & me,TCount size)982 insertGaps(Iter<TGaps, GapsIterator<ArrayGaps> > const & me,
983 		   TCount size)
984 {
985 	typedef typename Size<TGaps>::Type TGapsSize;
986 	typedef String<TGapsSize> TArr;
987 
988 	TArr & arr = _dataArr(container(me));
989 	typename Position<TArr>::Type pos = me.data_block;
990 	typename Position<TArr>::Type block_pos = me.data_sub;
991 	typename Iterator<TArr, Rooted>::Type it = iter(arr, pos);
992 
993 	if (isGap(me))
994 	{
995 		if (pos < length(arr))
996 		{//here is a gap: expand it
997 SEQAN_CHECKPOINT
998 			value(it) += size;
999 		}
1000 		else if (pos == length(arr))
1001 		{//gap at end: add trailing gaps
1002 			if (block_pos <= _getTrailingGaps(container(me)))
1003 			{
1004 SEQAN_CHECKPOINT
1005 				_setTrailingGaps(container(me), _getTrailingGaps(container(me)) + size);
1006 			}
1007 			return; // do not adjust right position
1008 		}
1009 	}
1010 	else
1011 	{
1012 		if (me.data_sub)
1013 		{//here is no gap: insert one
1014 SEQAN_CHECKPOINT
1015 			_clearSpace(arr, 2, pos + 1, pos + 1, Generous());
1016 
1017 			it = iter(arr, pos); //reload, since iterator could be invalid
1018 			it[2] = it[0] - me.data_sub;
1019 			it[0] = me.data_sub;
1020 			it[1] = size;
1021 
1022 			++me.data_block;
1023 			me.data_sub = 0;
1024 		}
1025 		else
1026 		{//insert gaps at begin of a char block
1027 SEQAN_CHECKPOINT
1028 			me.data_sub = arr[pos - 1];
1029 			arr[pos - 1] += size; //note: pos > 0 because this is a char block
1030 			--me.data_block;
1031 		}
1032 	}
1033 
1034 	//adjust right position
1035 	_setEndPosition(container(me), endPosition(container(me)) + size);
1036 }
1037 
1038 //____________________________________________________________________________
1039 
1040 //delete up to size gaps
1041 
1042 template <typename TGaps, typename TCount>
1043 inline void
removeGaps(Iter<TGaps,GapsIterator<ArrayGaps>> const & me,TCount _size)1044 removeGaps(Iter<TGaps, GapsIterator<ArrayGaps> > const & me,
1045 		   TCount _size)
1046 {
1047 	typedef typename Size<TGaps>::Type TGapsSize;
1048 	typedef String<TGapsSize> TArr;
1049 
1050 	TGapsSize size = _size;
1051 
1052 	if (isGap(me))
1053 	{//here is a gap
1054 		TArr & arr = _dataArr(container(me));
1055 		if (me.data_block < length(arr))
1056 		{//not trailing gaps
1057 			typename Iterator<TArr, Standard>::Type it = iter(arr, me.data_block);
1058 
1059 			TGapsSize rest_size = *it - me.data_sub;
1060 
1061 			if (rest_size <= size)
1062 			{//remove rest of gap block
1063 				if (me.data_sub || !me.data_block)
1064 				{//keep this gap block
1065 SEQAN_CHECKPOINT
1066 					*it = me.data_sub;
1067 					++me.data_block;
1068 					me.data_sub = 0;
1069 				}
1070 				else
1071 				{//remove complete gap block, merge two char blocks
1072 SEQAN_CHECKPOINT
1073 					TGapsSize data_sub = *(it-1); //(note: this well defined, since "it" was not the first block)
1074 
1075 					*(it-1) += *(it+1);
1076 					replace(_dataArr(container(me)), me.data_block, me.data_block + 2, "");
1077 
1078 					--me.data_block;
1079 					me.data_sub = data_sub;
1080 				}
1081 				_setEndPosition(container(me), endPosition(container(me)) - rest_size);
1082 			}
1083 			else
1084 			{//remove a part of this block
1085 SEQAN_CHECKPOINT
1086 				*it -= size;
1087 				_setEndPosition(container(me), endPosition(container(me)) - size);
1088 			}
1089 		}
1090 		else if (me.data_block == length(arr))
1091 		{
1092 			if (me.data_sub <= _getTrailingGaps(container(me)))
1093 			{//remove trailing gaps
1094 				if (size > (_getTrailingGaps(container(me)) - me.data_sub))
1095 				{//remove more than exists
1096 SEQAN_CHECKPOINT
1097 					_setTrailingGaps(container(me), me.data_sub);
1098 				}
1099 				else
1100 				{//remove less or equal than exists
1101 SEQAN_CHECKPOINT
1102 					_setTrailingGaps(container(me), _getTrailingGaps(container(me)) - size);
1103 				}
1104 			}
1105 		}//else: trailing gaps: infinite gaps here - nothing to do
1106 	}
1107 	//else: here is no gap - nothing to do
1108 }
1109 
1110 //____________________________________________________________________________
1111 
1112 template <typename TGaps1, typename TGaps2>
1113 inline bool
1114 operator == (Iter<TGaps1, GapsIterator<ArrayGaps> > & _left,
1115 			 Iter<TGaps2, GapsIterator<ArrayGaps> > & _right)
1116 {
1117 SEQAN_CHECKPOINT
1118 	return (_left.data_block == _right.data_block) && (_left.data_sub == _right.data_sub);
1119 }
1120 template <typename TGaps1, typename TGaps2>
1121 inline bool
1122 operator == (Iter<TGaps1, GapsIterator<ArrayGaps> > const & _left,
1123 			 Iter<TGaps2, GapsIterator<ArrayGaps> > & _right)
1124 {
1125 SEQAN_CHECKPOINT
1126 	return (_left.data_block == _right.data_block) && (_left.data_sub == _right.data_sub);
1127 }
1128 template <typename TGaps1, typename TGaps2>
1129 inline bool
1130 operator == (Iter<TGaps1, GapsIterator<ArrayGaps> > & _left,
1131 			 Iter<TGaps2, GapsIterator<ArrayGaps> > const & _right)
1132 {
1133 SEQAN_CHECKPOINT
1134 	return (_left.data_block == _right.data_block) && (_left.data_sub == _right.data_sub);
1135 }
1136 template <typename TGaps1, typename TGaps2>
1137 inline bool
1138 operator == (Iter<TGaps1, GapsIterator<ArrayGaps> > const & _left,
1139 			 Iter<TGaps2, GapsIterator<ArrayGaps> > const & _right)
1140 {
1141 SEQAN_CHECKPOINT
1142 	return (_left.data_block == _right.data_block) && (_left.data_sub == _right.data_sub);
1143 }
1144 //____________________________________________________________________________
1145 
1146 template <typename TGaps1, typename TGaps2>
1147 inline bool
1148 operator != (Iter<TGaps1, GapsIterator<ArrayGaps> > & _left,
1149 			 Iter<TGaps2, GapsIterator<ArrayGaps> > & _right)
1150 {
1151 SEQAN_CHECKPOINT
1152 	return (_left.data_block != _right.data_block) || (_left.data_sub != _right.data_sub);
1153 }
1154 template <typename TGaps1, typename TGaps2>
1155 inline bool
1156 operator != (Iter<TGaps1, GapsIterator<ArrayGaps> > const & _left,
1157 			 Iter<TGaps2, GapsIterator<ArrayGaps> > & _right)
1158 {
1159 SEQAN_CHECKPOINT
1160 	return (_left.data_block != _right.data_block) || (_left.data_sub != _right.data_sub);
1161 }
1162 template <typename TGaps1, typename TGaps2>
1163 inline bool
1164 operator != (Iter<TGaps1, GapsIterator<ArrayGaps> > & _left,
1165 			 Iter<TGaps2, GapsIterator<ArrayGaps> > const & _right)
1166 {
1167 SEQAN_CHECKPOINT
1168 	return (_left.data_block != _right.data_block) || (_left.data_sub != _right.data_sub);
1169 }
1170 template <typename TGaps1, typename TGaps2>
1171 inline bool
1172 operator != (Iter<TGaps1, GapsIterator<ArrayGaps> > const & _left,
1173 			 Iter<TGaps2, GapsIterator<ArrayGaps> > const & _right)
1174 {
1175 SEQAN_CHECKPOINT
1176 	return (_left.data_block != _right.data_block) || (_left.data_sub != _right.data_sub);
1177 }
1178 
1179 
1180 //////////////////////////////////////////////////////////////////////////////
1181 }// namespace SEQAN_NAMESPACE_MAIN
1182 
1183 //////////////////////////////////////////////////////////////////////////////
1184 #endif //#ifndef SEQAN_HEADER_...
1185