$$\phantom{******** If you see this on the webpage then the browser could not locate *********}$$
$$\phantom{******** jsMath/easy/load.js or the variable root is set wrong in this file *********}$$
$$\newcommand{\vector}[1]{\left[\begin{array}{c} #1 \end{array}\right]}$$ $$\newenvironment{matrix}{\left[\begin{array}{cccccccccc}} {\end{array}\right]}$$ $$\newcommand{\A}{{\cal A}}$$ $$\newcommand{\W}{{\cal W}}$$

/* Copyright 2002-2014 Research Foundation State University of New York   */

/* This file is part of QUB Express.                                      */

/* QUB Express is free software; you can redistribute it and/or modify    */
/* it under the terms of the GNU General Public License as published by   */
/* the Free Software Foundation, either version 3 of the License, or      */
/* (at your option) any later version.                                    */

/* QUB Express is distributed in the hope that it will be useful,         */
/* but WITHOUT ANY WARRANTY; without even the implied warranty of         */
/* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the          */
/* GNU General Public License for more details.                           */

/* You should have received a copy of the GNU General Public License,     */
/* named LICENSE.txt, in the QUB Express program directory.  If not, see  */
/* <http://www.gnu.org/licenses/>.                                        */

#include "qubfast.h"
#include "qubidl.h"

/*

*/

#include <iostream>

#include <map>
#include <math.h>
#include <vector>

#include "QUB_Freezer.h"
using namespace std;

#define QUB_IDLSEG_CHUNKSHIFT 14
#define QUB_IDLSEG_CHUNKSIZE (1 << QUB_IDLSEG_CHUNKSHIFT)


typedef struct
{
	int cls;
	int first;
	int last;
} QUB_IdlDwell;


typedef map< int, QUB_IdlDwell >  QUB_IdlMap;


class QUB_IdlChunk
{
protected:
	void setComplete(int count, QUB_IdlDwell *buf, int offset, int pointsInSeg);
	void setComplete(int count, int *ff, int *ll, int *cc, int offset, int pointsInSeg);

public:
	int first;

	int storage;
	QUB_IdlMap *asMap;
        vector<QUB_IdlDwell> asArr;

	QUB_IdlChunk(int afirst);
	virtual ~QUB_IdlChunk();

	int size();
	void setupStorage(int lvl);
	void clearStorage(int lvl);
	int  getStorage();

	// the rest of the methods do nothing unless you have storage level 0
	void clear(int f, int l);
	void set(int f, int l, int c, bool atBegin, bool atEnd); // assumes area is cleared

	// clears affected, returns offset in buf of first dwell not completely consumed
	int set(int count, QUB_IdlDwell *buf, int offset, int pointsInSeg); // offset == seg.first
	int set(int count, int *ff, int *ll, int *cc, int offset, int pointsInSeg);
};

typedef QUB_Freezer<QUB_IdlChunk> QUB_IdlFreezer;
typedef QUB_FreezerSeries<QUB_IdlChunk> QUB_IdlFreezerSeries;


class QUB_IdlMaps;

class QUB_IdlMapsIter
{
private:
	QUB_IdlMaps* maps;

public:
	QUB_IdlChunk *chunk;
	int          ix;
	
	QUB_IdlMap::iterator iter;

	QUB_IdlDwell *arr;
	int          i, count;

	QUB_IdlMapsIter()
		: maps(0), chunk(0), ix(-1), arr(0)
	{}
	QUB_IdlMapsIter(const QUB_IdlMapsIter& other)
		: maps(0), chunk(0), ix(-1), arr(0)
	{
		*this = other;
	}

	QUB_IdlMapsIter(QUB_IdlMaps& amaps, int point );
	QUB_IdlMapsIter(QUB_IdlMaps& amaps, bool atBegin);
	virtual ~QUB_IdlMapsIter();

	QUB_IdlMapsIter& operator=(const QUB_IdlMapsIter& other);

	bool operator==(const QUB_IdlMapsIter& other) const { return arr ? (arr == other.arr && i == other.i) : iter == other.iter; }
	bool operator!=(const QUB_IdlMapsIter& other) const { return ! (*this == other); }
	operator bool () const { return chunk && (arr ? (i != count) : (iter != chunk->asMap->end())); }
	QUB_IdlDwell& operator* () const { return arr ? arr[i] : iter->second; }
	QUB_IdlDwell* operator->() const { return &(arr ? arr[i] : iter->second); }

	QUB_IdlMapsIter& operator++ ();
	QUB_IdlMapsIter& operator-- ();
};



class QUB_IdlMaps
{
private:
	void setSize( int n );
public:

	QUB_IdlFreezer& freezer;
	QUB_IdlFreezerSeries chunks;
	int maxmap;
	int points;

	QUB_IdlMaps(QUB_IdlFreezer& freezer, int points);
	virtual ~QUB_IdlMaps();

	void insert(int f, int l);
	void remove(int f, int l);
	
	void set(int f, int l, int c);
	void join(int f, int l);

	// return offset in buf of first dwell not completely consumed
	int set(int count, QUB_IdlDwell *buf, int offset); // offset == seg.first
	int set(int count, int *ff, int *ll, int *cc, int offset);

	typedef QUB_IdlMapsIter iterator;

	iterator begin()         { return QUB_IdlMapsIter( *this, true ); }  // { return maps.size() ? iterator(maps, maxmap, 0, maps[0].begin()) : end(); }
	iterator end()           { return QUB_IdlMapsIter( *this, false ); } // { return iterator(maps, maxmap, maxmap, (maxmap+1) ? maps[maxmap].end() : QUB_IdlMap::iterator()); }
	iterator find(int point) { return QUB_IdlMapsIter( *this, point ); }
};

class QUB_IdlSeg
{
public:
	QUB_IdlFreezer&      freezer;
	double               sampling;
	int                  first;
	int                  last;
	QUB_IdlMaps          maps;


	QUB_IdlSeg(QUB_IdlFreezer& freezer, int first, int last, double samp);

	void      setFirst(int f);
	void      setLast(int l);
	void      offset(int points);

	void      insert(int f, int l);
	void      remove(int f, int l);
	
	void      set(int f, int l, int c);
	void      join(int f, int l);

        int       getDwellCount(int f, int l, int fragments);
        int       getDwellCountAndGaps(int f, int l);
        int       getDwells(int f, int l, int fragments, int *firsts, int *lasts, int *classes, float *durations);
        int       getDwellsAndGaps(int f, int l, int *firsts, int *lasts, int *classes, float *durations);
	void      setDwells(int count, QUB_IdlDwell *buf);
	void      setDwells(int count, int *ff, int *ll, int *cc);

	QUB_IdlMaps::iterator begin()         { return maps.begin(); }
	QUB_IdlMaps::iterator end()           { return maps.end(); }
	QUB_IdlMaps::iterator find(int point) { return maps.find(point - first); }
private:
  void defrag(int& f, int& l);
};




class QUB_IdlSegmentsIter;

class QUB_IdlSegments
{
private:
	QUB_IdlFreezer      freezer;

	void       getSegRange(int f, int l, int& sf, int& sl);
	void       resetSeg(int s);

public:
	std::vector<QUB_IdlSeg*> segs;
	double              sampling;

	QUB_IdlSegments(double sampling);
	~QUB_IdlSegments();

	double     getSampling();
	void       setSampling(double samp);
	int        getSegCount();
        int        getSegFirst(int i);
        int        getSegLast(int i);

        int        getDwellCount(int f, int l, int fragments);
        int        getDwellCountAndGaps(int f, int l);
        int        getDwells(int f, int l, int fragments, int *firsts, int *lasts, int *classes, float *durations);
        int        getDwellsAndGaps(int f, int l, int *firsts, int *lasts, int *classes, float *durations);
        int        setDwells(int count, int *ff, int *ll, int *cc);
        void       sampleDwells(int f, int l, double *amp, int Nsample, float *los, float *his, unsigned int *class_bits);
        void       sampleDwells(int f, int l, float *amp, float *std, int Nsample, float *los, float *his);
  
	int        clear();
	int        addSeg(int f, int l);
	int        removePoints(int f, int l);
	int        insertPoints(int f, int l);
	int        breakSeg(int at);
	int        joinSegs(int at);

	int        erase(int f, int l);
	int        join(int f, int l);

	typedef QUB_IdlSegmentsIter iterator;
	iterator           begin();
	iterator           end();
	iterator           find(int point);
};


class QUB_IdlSegmentsIter
{
private:
	QUB_IdlSegments&       segs;
	int                    seg_ix;
	QUB_IdlMaps::iterator  iter;
	QUB_IdlSeg*            seg;

public:
	QUB_IdlSegmentsIter(QUB_IdlSegments& asegs, int aseg_ix, const QUB_IdlMaps::iterator& aiter )
		: segs( asegs ), seg_ix( aseg_ix ), iter( aiter ), seg( (seg_ix>=0) ? segs.segs[seg_ix] : NULL )
	{
	}

	QUB_IdlSegmentsIter(const QUB_IdlSegmentsIter &other)
		: segs( other.segs ), seg_ix( other.seg_ix ), iter( other.iter ), seg( other.seg )
	{
	}

	virtual ~QUB_IdlSegmentsIter()
	{
	}

	QUB_IdlSegmentsIter& operator= (QUBFAST_VAR_NOT_USED const QUB_IdlSegmentsIter& other)
	{
		throw "undefined";
	}

	bool operator==(const QUB_IdlSegmentsIter& other) const { return iter == other.iter; }
	operator bool () const { return seg && iter; } // (iter != seg->end()); }
	QUB_IdlDwell& operator* () const { return *iter; }
	QUB_IdlDwell* operator->() const { return &(*iter); }
	// iter->first, last are in seg.local terms -- add segfirst() to find location in whole data file
	int segfirst() const { return seg->first; }
	int segindex() const { return seg_ix; }

	QUB_IdlSegmentsIter& operator++ ();
	QUB_IdlSegmentsIter& operator-- ();

	int valid();
	void read( int *f, int *l, int *c );
	void next();
	void prev();
};


extern "C" QUBFAST_API void*  qubidl_create(double sampling)
{
  return new QUB_IdlSegments(sampling);
}

extern "C" QUBFAST_API void   qubidl_destroy(void *idl)
{
  delete (QUB_IdlSegments*) idl;
}

extern "C" QUBFAST_API double qubidl_get_sampling(void *idl)
{
  return ((QUB_IdlSegments*)idl)->sampling;
}

extern "C" QUBFAST_API void   qubidl_set_sampling(void *idl, double sampling)
{
  ((QUB_IdlSegments*)idl)->setSampling(sampling);
}

extern "C" QUBFAST_API int    qubidl_get_seg_count(void *idl)
{
  return ((QUB_IdlSegments*)idl)->getSegCount();
}

extern "C" QUBFAST_API int    qubidl_get_seg_first(void *idl, int i)
{
  return ((QUB_IdlSegments*)idl)->getSegFirst(i);
}

extern "C" QUBFAST_API int    qubidl_get_seg_last(void *idl, int i)
{
  return ((QUB_IdlSegments*)idl)->getSegLast(i);
}

extern "C" QUBFAST_API int    qubidl_count_dwells(void *idl, int f, int l, int fragments)
{
  return ((QUB_IdlSegments*)idl)->getDwellCount(f, l, fragments);
}

extern "C" QUBFAST_API int    qubidl_count_dwells_and_gaps(void *idl, int f, int l)
{
  return ((QUB_IdlSegments*)idl)->getDwellCountAndGaps(f, l);
}

extern "C" QUBFAST_API int    qubidl_get_dwells(void *idl, int f, int l, int fragments, int *firsts, int *lasts, int *classes, float *durations)
{
  return ((QUB_IdlSegments*)idl)->getDwells(f, l, fragments, firsts, lasts, classes, durations);
}

extern "C" QUBFAST_API int    qubidl_get_dwells_and_gaps(void *idl, int f, int l, int *firsts, int *lasts, int *classes, float *durations)
{
  return ((QUB_IdlSegments*)idl)->getDwellsAndGaps(f, l, firsts, lasts, classes, durations);
}

extern "C" QUBFAST_API void   qubidl_set_dwells(void *idl, int count, int *firsts, int *lasts, int *classes)
{
  ((QUB_IdlSegments*)idl)->setDwells(count, firsts, lasts, classes);
}

extern "C" QUBFAST_API void   qubidl_sample_dwells(void *idl, int f, int l, double *amp, int Nsample, float *lo, float *hi)
{
  ((QUB_IdlSegments*)idl)->sampleDwells(f, l, amp, Nsample, lo, hi, NULL);
}

extern "C" QUBFAST_API void   qubidl_sample_dwells_classes(void *idl, int f, int l, double *amp, int Nsample, float *lo, float *hi, unsigned int *class_bits)
{
  ((QUB_IdlSegments*)idl)->sampleDwells(f, l, amp, Nsample, lo, hi, class_bits);
}

extern "C" QUBFAST_API void   qubidl_sample_dwells_std(void *idl, int f, int l, float *amp, float *std, int Nsample, float *lo, float *hi)
{
  ((QUB_IdlSegments*)idl)->sampleDwells(f, l, amp, std, Nsample, lo, hi);
}

extern "C" QUBFAST_API void   qubidl_clear(void *idl)
{
  ((QUB_IdlSegments*)idl)->clear();
}

extern "C" QUBFAST_API void   qubidl_add_seg(void *idl, int f, int l)
{
  ((QUB_IdlSegments*)idl)->addSeg(f, l);
}

extern "C" QUBFAST_API void   qubidl_erase(void *idl, int f, int l)
{
  ((QUB_IdlSegments*)idl)->erase(f, l);
}

extern "C" QUBFAST_API void   qubidl_join(void *idl, int f, int l)
{
  ((QUB_IdlSegments*)idl)->join(f, l);
}

extern "C" QUBFAST_API void*  qubidl_begin(void *idl)
{
  return new QUB_IdlSegments::iterator( ((QUB_IdlSegments*)idl)->begin() );
}

extern "C" QUBFAST_API void*  qubidl_end(void *idl)
{
  return new QUB_IdlSegments::iterator( ((QUB_IdlSegments*)idl)->end() );
}

extern "C" QUBFAST_API void*  qubidl_find(void *idl, int f)
{
  return new QUB_IdlSegments::iterator( ((QUB_IdlSegments*)idl)->find(f) );
}

extern "C" QUBFAST_API void*  qubidl_iter_clone(void *iter)
{
  return new QUB_IdlSegments::iterator( *((QUB_IdlSegments::iterator *) iter) );
}

extern "C" QUBFAST_API void   qubidl_iter_destroy(void *iter)
{
  delete (QUB_IdlSegments::iterator *) iter;
}

extern "C" QUBFAST_API int    qubidl_iter_valid(void *iter)
{
  return ((QUB_IdlSegments::iterator *) iter)->valid();
}

extern "C" QUBFAST_API int    qubidl_iter_equals(void *iter, void *iter2)
{
  return *((QUB_IdlSegments::iterator *) iter) == *((QUB_IdlSegments::iterator *) iter2);
}

extern "C" QUBFAST_API int    qubidl_iter_segindex(void *iter)
{
  return ((QUB_IdlSegments::iterator *) iter)->segindex();
}

extern "C" QUBFAST_API int    qubidl_iter_segfirst(void *iter)
{
  return ((QUB_IdlSegments::iterator *) iter)->segfirst();
}

extern "C" QUBFAST_API void   qubidl_iter_read(void *iter, int *f, int *l, int *c)
{
  ((QUB_IdlSegments::iterator *) iter)->read(f, l, c);
}

extern "C" QUBFAST_API void   qubidl_iter_next(void *iter)
{
  ((QUB_IdlSegments::iterator *) iter)->next();
}

extern "C" QUBFAST_API void   qubidl_iter_prev(void *iter)
{
  ((QUB_IdlSegments::iterator *) iter)->prev();
}





inline int WhichMap(int first)
{
	return first >> QUB_IDLSEG_CHUNKSHIFT;// / QUB_IDLSEG_CHUNKSIZE;
}

inline int HowManyMaps(int count)
{
	if ( count <= 0 )
		return 0;
	else
		return (count >> QUB_IDLSEG_CHUNKSHIFT) + ((count % QUB_IDLSEG_CHUNKSIZE) ? 1 : 0);
}



QUB_IdlChunk::QUB_IdlChunk(int afirst)
	: first( afirst ), storage( QUB_FS_0 )
{
	asMap = new QUB_IdlMap;
}

QUB_IdlChunk::~QUB_IdlChunk()
{
	if ( asMap )
		delete asMap;
}

int QUB_IdlChunk::size()
{
	if ( asMap )
	  return (int) asMap->size();
	return (int) asArr.size();
}

void QUB_IdlChunk::setupStorage(int lvl)
{
  switch ( lvl )
    {
    case 0:
      if ( ! asMap ) {
        asMap = new QUB_IdlMap;
        for ( int i = ((int) asArr.size()) - 1; i >= 0; --i )
          (*asMap)[ asArr[i].first ] = asArr[i];
      }
      storage |= QUB_FS_0;
      break;
    case 1:
      if ( asArr.size() != asMap->size() ) {
        asArr.resize(asMap->size());
        if ( asArr.size() ) {
          QUB_IdlDwell *dwells = & asArr[0];
          for ( QUB_IdlMap::iterator di = asMap->begin(); di != asMap->end(); ++di )
            *(dwells++) = di->second;
        }
	  }
      storage |= QUB_FS_1;
      break;
    }
}

void QUB_IdlChunk::clearStorage(int lvl)
{
	switch ( lvl )
	{
	case 0:
		if ( asMap ) {
			delete asMap;
			asMap = NULL;
		}
		storage &= ~ QUB_FS_0;
		break;
	case 1:
		asArr.clear();
		storage &= ~ QUB_FS_1;
		break;
	}
}

int QUB_IdlChunk::getStorage()
{
	return storage;
}

void QUB_IdlChunk::clear(int f, int l)
{
	if ( ! asMap ) {
		if ( (f <= first) && (l >= (first + QUB_IDLSEG_CHUNKSIZE - 1)) ) {
			asArr.clear();
			return;
		}
		else {
			setupStorage(0);
		}
	}

	if ( ! asMap->size() )
		return;

	f = max(f, first);
	l = min(l, first + QUB_IDLSEG_CHUNKSIZE - 1);

	QUB_IdlMap::iterator di = asMap->upper_bound(f);
	--di;
	if ( di->second.last < f )
		return;

	if ( di->second.first < f ) {
		if ( l < di->second.last ) {
			QUB_IdlDwell& dwell = (*asMap)[l+1];
			dwell.first = l+1;
			dwell.last = di->second.last;
			dwell.cls = di->second.cls;
		}
		di->second.last = f-1;
		++di;
	}
	while ( di != asMap->end() && di->second.last <= l ) {
		QUB_IdlMap::iterator toErase = di;
		++di;
		asMap->erase(toErase);
	}
	if ( di != asMap->end() && di->second.first <= l ) {
		QUB_IdlDwell& dwell = (*asMap)[l+1];
		dwell.first = l+1;
		dwell.last = di->second.last;
		dwell.cls = di->second.cls;
		asMap->erase(di);
	}
}

void QUB_IdlChunk::set(int f, int l, int c, QUBFAST_VAR_NOT_USED bool atBegin, QUBFAST_VAR_NOT_USED bool atEnd) // assuming already cleared
{
	if ( ! asMap )
		return;

	f = max(f, first);
	l = min(l, first + QUB_IDLSEG_CHUNKSIZE - 1 );
	c = max(-1, c);
	
	QUB_IdlDwell *dwell = & (*asMap)[f];
	dwell->first = f;
	dwell->last = l;
	dwell->cls = c;

	if ( /*atBegin &&*/ f > first ) {
		QUB_IdlMap::iterator di = asMap->upper_bound(f-1);
		if ( di != asMap->begin() ) {
			--di;
			if ( di->second.cls == c && (di->second.last == (f-1)) ) {
				di->second.last = l;
				asMap->erase( asMap->find(f) );
				dwell = & di->second;
			}
		}
	}
	if ( /*atEnd &&*/ l < (first + QUB_IDLSEG_CHUNKSIZE - 1) ) {
		QUB_IdlMap::iterator di = asMap->upper_bound(l+1);
		--di;
		if ( di->second.cls == c && (di->first == (l+1)) ) {
			dwell->last = di->second.last;
			asMap->erase(di);
		}
	}
}

void QUB_IdlChunk::setComplete(int count, QUB_IdlDwell *buf, int offset, int pointsInSeg)
{
	if ( asMap )
		asMap->clear();
	setupStorage(1);
	clearStorage(0);

	asArr.resize(count);
	for ( int i=0; i<count; ++i ) {
		asArr[i].first = buf[i].first - offset;
		asArr[i].last = buf[i].last - offset;
		asArr[i].cls = max(-1, buf[i].cls);
	}
	asArr[0].first = max(first, asArr[0].first);
	asArr[count-1].last = min(asArr[count-1].last, min(first + QUB_IDLSEG_CHUNKSIZE, pointsInSeg) - 1);
}

void QUB_IdlChunk::setComplete(int count, int *ff, int *ll, int *cc, int offset, int pointsInSeg)
{
	if ( asMap )
		asMap->clear();
	setupStorage(1);
	clearStorage(0);

	asArr.resize(count);
	for ( int i=0; i<count; ++i ) {
		asArr[i].first = ff[i] - offset;
		asArr[i].last = ll[i] - offset;
		asArr[i].cls = max(-1, cc[i]);
	}
	asArr[0].first = max(first, asArr[0].first);
	asArr[count-1].last = min(asArr[count-1].last, min(first + QUB_IDLSEG_CHUNKSIZE, pointsInSeg) - 1);
}

int QUB_IdlChunk::set(int count, QUB_IdlDwell *buf, int offset, int pointsInSeg)
{
	int i0 = 0;
	while ( i0 < count && (buf[i0].last - offset) < first )
		++i0;
	if ( i0 == count )
		return count;

	int n = i0;
	while ( (n+1) < count && (buf[n].last - offset - first) < QUB_IDLSEG_CHUNKSIZE )
		++n; // n: first event ending past this chunk
	if ( (n+1) < count && (buf[n].first - offset - first) >= QUB_IDLSEG_CHUNKSIZE )
		--n; // n: last event fully/partially in chunk

	if ( (buf[i0].first - offset <= first) && (buf[n].last - offset >= (first + QUB_IDLSEG_CHUNKSIZE - 1)) ) {
		setComplete(n - i0 + 1, buf+i0, offset, pointsInSeg);
	}
	else {
		setupStorage(0);
		clear( buf[i0].first - offset, buf[n].last - offset );
		set( buf[i0].first - offset, min(buf[i0].last - offset, pointsInSeg-1), buf[i0].cls, true, i0==n );
		int i = i0 + 1;
		while ( i < n ) {
		  set( buf[i].first - offset, buf[i].last - offset, buf[i].cls, false, false );
		  ++i;
		}
		if ( i0 != n )
		  set( buf[n].first - offset,
		       min(buf[n].last - offset, pointsInSeg-1), buf[n].cls, false, true );
	}
	return n + ((buf[n].last - offset - first) < QUB_IDLSEG_CHUNKSIZE);
}

int QUB_IdlChunk::set(int count, int *ff, int *ll, int *cc, int offset, int pointsInSeg)
{
	int i0 = 0;
	while ( i0 < count && (ll[i0] - offset) < first )
		++i0;
	if ( i0 == count )
		return count;

	int n = i0;
	while ( (n+1) < count && (ll[n] - offset - first) < QUB_IDLSEG_CHUNKSIZE )
		++n;
	if ( (n+1) < count && (ff[n] - offset - first) >= QUB_IDLSEG_CHUNKSIZE )
		--n;

	if ( (ff[i0] - offset - first <= 0) && (ll[n] - offset - first >= (QUB_IDLSEG_CHUNKSIZE - 1)) ) {
		setComplete(n - i0 + 1, ff+i0, ll+i0, cc+i0, offset, pointsInSeg);
	}
	else {
		setupStorage(0);
		clear( ff[i0] - offset, ll[n] - offset );
		set( ff[i0] - offset, min(ll[i0] - offset, pointsInSeg-1), cc[i0], true, i0==n );
		int i = i0 + 1;
		while ( i < n ) {
			set( ff[i] - offset, ll[i] - offset, cc[i], false, false );
			++i;
		}
		if ( i0 != n )
			set( ff[n] - offset, 
			min(ll[n] - offset, pointsInSeg-1), cc[n], false, true );
	}
	return n + ((ll[n] - offset - first) < QUB_IDLSEG_CHUNKSIZE);
}


QUB_IdlMaps::QUB_IdlMaps(QUB_IdlFreezer& afreezer, int apoints)
: freezer( afreezer ), chunks( &afreezer ), points( 0 )
{
	if ( apoints ) {
		setSize( apoints );
		for ( int i=0; i<=maxmap; ++i ) {
			QUB_IdlChunk *chunk = chunks.checkOut( i, 0x1 );
			chunk->set( i*QUB_IDLSEG_CHUNKSIZE, min(points, (i+1) * QUB_IDLSEG_CHUNKSIZE) - 1, -1, false, false );
			chunks.checkIn( i, true );
		}
	}
	else
		maxmap = -1;
}

QUB_IdlMaps::~QUB_IdlMaps()
{
}

void QUB_IdlMaps::setSize( int n )
{
	if ( n > points ) {
		maxmap = HowManyMaps(n) - 1;
		for ( int i=chunks.size(); i<=maxmap; ++i )
			chunks.insert( chunks.size(), new QUB_IdlChunk( i * QUB_IDLSEG_CHUNKSIZE ) );
		points = n;
	}
	else if ( n < points ) {
		chunks.truncate( HowManyMaps(n) );
		maxmap = chunks.size() - 1;
		points = n;

		if ( points % QUB_IDLSEG_CHUNKSIZE ) {
			QUB_IdlChunk *chunk = chunks.checkOut(maxmap, 0x1);
			chunk->clear(points, points + QUB_IDLSEG_CHUNKSIZE);
			chunks.checkIn(maxmap, true);
		}
	}
}

void QUB_IdlMaps::insert(int f, int l)
{
	vector<QUB_IdlDwell> dwells;
	if ( true ) {
		iterator di = find(f);
		if ( di && di->first < f ) {
			dwells.resize(1);
			dwells[0].first = f;
			dwells[0].last = di->last;
			dwells[0].cls = di->cls;
			++di;
		}
		while ( di ) {
			dwells.push_back( *di );
			++di;
		}
	}

	int oldLen = points;
	int offset = l - f + 1;
	if ( f < points )
		remove(f, points-1);

	setSize( oldLen + offset );
	set(f, l, -1);
	if ( dwells.size() )
	  set((int) dwells.size(), &(dwells[0]), - offset);
}

void QUB_IdlMaps::remove(int f, int l)
{
	vector<QUB_IdlDwell> dwells;
	if ( true ) { // make di deleted before making changes
		iterator di = find(l+1);
		if ( di && di->first <= l ) {
			dwells.resize(1);
			dwells[0].first = l+1;
			dwells[0].last = di->last;
			dwells[0].cls = di->cls;
			++di;
		}
		while ( di ) {
			dwells.push_back( *di );
			++di;
		}
	}
	
	int origPoints = points;
	int offset = (l - f + 1);
	setSize( f );
	setSize( origPoints - offset );
	if ( dwells.size() )
	  set((int) dwells.size(), &(dwells[0]), offset);
}

void QUB_IdlMaps::set(int f, int l, int c)
{
	set(1, &f, &l, &c, 0);
}

void QUB_IdlMaps::join(int f, int l)
{
	iterator di = find(f);
	if ( di ) {
		if ( f && (di->first == f) )
			--di;
		set(f, l, di->cls);
	}
}

int QUB_IdlMaps::set(int count, QUB_IdlDwell *buf, int offset)
{
	if ( ! count )
		return 0;

	int kf = max( WhichMap(buf[0].first - offset), 0 );
	int kl = min( WhichMap(buf[count-1].last - offset), maxmap );
	int totalConsumed = 0;

	for ( int ki=kf; totalConsumed<count && ki<=kl; ++ki ) {
		QUB_IdlChunk *chunk = chunks.checkOut(ki, 0x1);
		int consumed = chunk->set(count - totalConsumed, buf + totalConsumed, offset, points);
		totalConsumed += consumed;
		chunks.checkIn( ki, true );
	}
	return totalConsumed;
}

int QUB_IdlMaps::set(int count, int *ff, int *ll, int *cc, int offset)
{
	if ( ! count )
		return 0;

	int kf = max( WhichMap(ff[0] - offset), 0 );
	int kl = min( WhichMap(ll[count-1] - offset), maxmap );
	int totalConsumed = 0;

	for ( int ki=kf; totalConsumed<count && ki<=kl; ++ki ) {
		QUB_IdlChunk *chunk = chunks.checkOut(ki, 0x1);
		int consumed = chunk->set(count - totalConsumed, ff + totalConsumed, ll + totalConsumed, cc + totalConsumed, offset, points);
		totalConsumed += consumed;
		chunks.checkIn( ki, true );
	}
	return totalConsumed;
}

QUB_IdlMapsIter::QUB_IdlMapsIter(QUB_IdlMaps& amaps, int point)
: maps( &amaps )
{
	if ( point < 0 || point >= maps->points ) {
		ix = maps->maxmap;
		if ( ix >= 0 ) {
			chunk = maps->chunks.checkOut( ix, 0x1 );
			if ( chunk->getStorage() & 0x1 ) {
				iter = chunk->asMap->end();
				arr = NULL;
			}
			else {
			        arr = (QUB_IdlDwell *) & chunk->asArr[0];
				count = chunk->size();
				i = count;
			}
		}
		else {
			chunk = NULL;
		}
		return;
	}

	ix = WhichMap(point);
	chunk = maps->chunks.checkOut( ix, 0x1 );
	if ( chunk->getStorage() & QUB_FS_0 ) {
		iter = chunk->asMap->upper_bound(point);
		if ( chunk->asMap->size() )
			--iter;
		arr = NULL;
	}
	else {
		arr = (QUB_IdlDwell *) & chunk->asArr[0];
		count = chunk->size();
		for ( i=0; i < count && arr[i].last < point; ++i )
			;
	}
}

QUB_IdlMapsIter::QUB_IdlMapsIter(QUB_IdlMaps& amaps, bool atBegin)
: maps( &amaps )
{
	if ( maps->maxmap < 0 ) {
		ix = maps->maxmap;
		chunk = NULL;
	}
	else {
		if ( atBegin ) {
			ix = 0;
			chunk = maps->chunks.checkOut( ix, 0x1 );
			if ( chunk->getStorage() & QUB_FS_0 ) {
				iter = chunk->asMap->begin();
				arr = NULL;
			}
			else {
				arr = (QUB_IdlDwell *) & chunk->asArr[0];
				count = chunk->size();
				i = 0;
			}
		}
		else {
			ix = maps->maxmap;
			chunk = maps->chunks.checkOut( ix, 0x1 );
			if ( chunk->getStorage() & QUB_FS_0 ) {
				iter = chunk->asMap->end();
				arr = NULL;
			}
			else {
				arr = (QUB_IdlDwell *) & chunk->asArr[0];
				count = chunk->size();
				i = count;
			}
		}
	}
}

QUB_IdlMapsIter::~QUB_IdlMapsIter()
{
	if ( ix >= 0 )
		maps->chunks.checkIn( ix, false );
}

QUB_IdlMapsIter& QUB_IdlMapsIter::operator=(const QUB_IdlMapsIter& other)
{
	if ( maps != other.maps || ix != other.ix ) {
		if ( ix >= 0 )
			maps->chunks.checkIn( ix, false );
		maps = other.maps;
		ix = other.ix;
		if ( ix >= 0 )
			chunk = maps->chunks.checkOut( ix, 0x1 );
		else
			chunk = NULL;
	}
	iter = other.iter;
	arr = other.arr;
	count = other.count;
	i = other.i;

	return *this;
}

QUB_IdlMapsIter& QUB_IdlMapsIter::operator++ ()
{
	if ( ix < 0 || (arr ? (i == count) : (iter == chunk->asMap->end())) )
		return *this;

	bool atEnd;
	if ( arr ) {
		++i;
		atEnd = i == count;
	}
	else {
		++iter;
		atEnd = iter == chunk->asMap->end();
	}
	if ( atEnd && ix < maps->maxmap ) {
		maps->chunks.checkIn( ix, false );

		++ix;
		chunk = maps->chunks.checkOut( ix, 0x1 );
		if ( chunk->getStorage() & QUB_FS_0 ) {
			iter = chunk->asMap->begin();
			arr = NULL;
		}
		else {
			arr = (QUB_IdlDwell *) & chunk->asArr[0];
			count = chunk->size();
			i = 0;
		}
	}
	return *this;
}

QUB_IdlMapsIter& QUB_IdlMapsIter::operator-- ()
{
	if ( arr ? i : iter != chunk->asMap->begin() ) {
		if ( arr )
			--i;
		else
			--iter;
	}
	else if ( ix ) {
		maps->chunks.checkIn( ix, false );

		--ix;
		chunk = maps->chunks.checkOut( ix, 0x1 );
		if ( chunk->getStorage() & QUB_FS_0 ) {
			iter = chunk->asMap->end();
			--iter;
			arr = NULL;
		}
		else {
			arr = (QUB_IdlDwell *) & chunk->asArr[0];
			count = chunk->size();
			i = count - 1;
		}
	}
	return *this;
}


QUB_IdlSeg::QUB_IdlSeg(QUB_IdlFreezer& afreezer, int afirst, int alast, double samp)
: freezer( afreezer ), sampling( samp ), first( afirst ), last( alast ), maps( afreezer, alast - afirst + 1 )
{
}

void      QUB_IdlSeg::setFirst(int f)
{
	if ( f > first ) {
		if ( f > last )
			f = last + 1;
		maps.remove(0, (f-first)-1);
	}
	else if ( f < first ) {
		maps.insert(0, (first - f)-1);
	}
	first = f;
}

void      QUB_IdlSeg::setLast(int l)
{
	if ( l < last ) {
		if ( l < first )
			l = first - 1;
		maps.remove(l - first + 1, last - first);
	}
	else if ( l > last ) {
		maps.insert(last - first + 1, l - first);
	}
	last = l;
}

void      QUB_IdlSeg::offset(int points)
{
	first += points;
	last += points;
}

void      QUB_IdlSeg::insert(int f, int l)
{
	maps.insert(f-first, l-first);
	last += (l - f + 1);
}

void      QUB_IdlSeg::remove(int f, int l)
{
	maps.remove(f-first, l-first);
	last -= (l - f + 1);
}

void      QUB_IdlSeg::set(int f, int l, int c)
{
	maps.set(f-first, l-first, c);
}

void      QUB_IdlSeg::join(int f, int l)
{
	maps.join(f-first, l-first);
}

void     QUB_IdlSeg::defrag(int &f, int &l)
  // fragments: start in middle of dwell, or
  //            start anywhere in repeat-class dwell
  //            end in middle of dwell, or
  //            end anywhere in prepeat-class dwell
{
  QUB_IdlMaps::iterator counter = find(f+first);
  int lastCls = counter->cls;
  bool skip = counter->first < f; // in middle of dwell
  if ( (! skip) && counter != begin() ) {     
    --counter;
    skip = (lastCls == counter->cls); // or prev dwell has same class
  }
  if ( skip ) {
    do
      ++counter;
    while ( counter && (counter->cls == lastCls) );
    f = counter ? counter->first : (l+1);
  }
  
  counter = find(l+first);
  lastCls = counter->cls;
  if ( counter->last > l )
    skip = true; // in middle of dwell
  else {
    ++counter;
    skip = (counter && (counter->cls == lastCls)); // or next dwell has same class
  }
  if ( skip ) {
    QUB_IdlMaps::iterator thefirst = begin();
    do
      --counter;
    while ( (counter->cls == lastCls) && (counter != thefirst) );
    l = counter ? counter->last : (f-1);
  }
}

//    new: class -1 is also a fragment
int       QUB_IdlSeg::getDwellCount(int f, int l, int fragments)
{
  f = max(f, first) - first;
  l = min(l, last)  - first;
  if ( ! fragments )
    defrag(f, l);
  
  QUB_IdlMaps::iterator counter = find(f+first);
  int count = 0;
  int lastCls = counter ? (counter->cls - 1) : -1;
  
  // count upon entering a new class
  while ( counter && (counter->first <= l) ) {
    if ( fragments || (counter->cls >= 0) ) {
      if ( counter->cls != lastCls && counter->cls >= 0 )
	++count;
    }
    lastCls = counter->cls;
    ++counter;
  }
  return count;
}
  
int       QUB_IdlSeg::getDwellCountAndGaps(int f, int l)
{
  f = max(f, first) - first;
  l = min(l, last)  - first;
  
  int count = 0;
  int lastCls = -382; // shouldn't be any negative classes anyway
  
  // count upon entering a new class
  for (QUB_IdlMaps::iterator counter = find(f+first); counter && (counter->first <= l); ++counter ) {
    if ( counter->cls != lastCls ) {
      ++count;
      lastCls = counter->cls;
    }
  }
  return count;
}
  
int       QUB_IdlSeg::getDwells(int f, int l, int fragments, int *firsts, int *lasts, int *classes, float *durations)
{
  f = max(f, first) - first;
  l = min(l, last)  - first;
  if ( ! fragments )
    defrag(f, l);
  
  QUB_IdlMaps::iterator counter = find(f+first);
  int lastCls = counter ? (counter->cls - 1) : -1;
  int i = -1;
  
  // count upon entering a new class
  while ( counter && (counter->first <= l) ) {
    if ( fragments || (counter->cls >= 0) ) {
      if ( counter->cls >= 0 ) {
	if ( counter->cls != lastCls ) {
	  ++i;
	  firsts[i] = counter->first + first;
	  classes[i] = counter->cls;
	}
	lasts[i] = counter->last + first;
      }
      lastCls = counter->cls;
    }
    ++counter;
  }
  
  if ( durations ) {
    for (int j=i; j>=0; --j)
      durations[j] = (float) ((lasts[j] - firsts[j] + 1) * sampling);
  }
  
  return i+1;
}
  
int       QUB_IdlSeg::getDwellsAndGaps(int f, int l, int *firsts, int *lasts, int *classes, float *durations)
{
  f = max(f, first) - first;
  l = min(l, last)  - first;
  
  int lastCls = -382;
  int i = -1;
  
  for (QUB_IdlMaps::iterator counter = find(f+first); counter && (counter->first <= l); ++counter ) {
    if ( counter->cls != lastCls ) {
      ++i;
      firsts[i] = counter->first + first;
      classes[i] = counter->cls;
    }
    lastCls = counter->cls;
    lasts[i] = counter->last + first;
  }
  
  if ( durations ) {
    for (int j=i; j>=0; --j)
      durations[j] = (float) ((lasts[j] - firsts[j] + 1) * sampling);
  }
  
  return ( i + 1 );
}
  
void      QUB_IdlSeg::setDwells(int count, QUB_IdlDwell *buf)
{
	int a = 0;
	while ( a<count && buf[a].first < first )
		++a;
	int b = count-1;
	while ( b>=a && buf[b].last > last )
		--b;

	maps.set(b - a + 1, buf+a, first);
}

void      QUB_IdlSeg::setDwells(int count, int *ff, int *ll, int *cc)
{
	if ( ! (ff && ll && cc) )
		return;

	int a = 0;
	while ( a<count && ff[a] < first )
		++a;
	int b = count-1;
	while ( b>=a && ll[b] > last )
		--b;

	maps.set(b - a + 1, ff+a, ll+a, cc+a, first);
}


QUB_IdlSegmentsIter& QUB_IdlSegmentsIter::operator++ ()
{
	if ( *this ) {
		++iter;
		if ( (seg_ix+1) < int(segs.segs.size()) && ! iter ) {
			++seg_ix;
			seg = segs.segs[seg_ix];
			iter = seg->begin();
		}
	}
	return *this;
}

QUB_IdlSegmentsIter& QUB_IdlSegmentsIter::operator-- ()
{
	if ( seg ) {
		if ( seg_ix > 0 && iter == seg->begin() ) {
			--seg_ix;
			seg = segs.segs[ seg_ix ];
			iter = seg->end();
		}
		--iter;
	}
	return *this;
}


int QUB_IdlSegmentsIter::valid()
{
	return *this ? 1 : 0;
}

void QUB_IdlSegmentsIter::read( int *f, int *l, int *c )
{
	QUB_IdlDwell& dwell = **this;
	*f = dwell.first + seg->first;
	*l = dwell.last + seg->first;
	*c = dwell.cls;
}

void QUB_IdlSegmentsIter::next()
{
	++(*this);
}

void QUB_IdlSegmentsIter::prev()
{
	--(*this);
}


QUB_IdlSegments::QUB_IdlSegments(double samp)
: sampling(samp)
{
	freezer.setThreshold(0, 100000);
}

QUB_IdlSegments::~QUB_IdlSegments()
{
	for ( int i=0; i<int(segs.size()); ++i )
		delete segs[i];
}

void       QUB_IdlSegments::getSegRange(int f, int l, int& sf, int& sl)
{
  int nseg = (int) segs.size();
	if ( nseg == 0 ) {
		sf = 0; sl = -1; return;
	}
	if ( segs[0]->last < segs[0]->first ) {
		sf = 0; sl = 0; return;
	}

	sf = 0;
	while ( (sf < nseg) && (segs[sf]->last < f) )
		++sf;
	sl = sf;
	while ( (sl < nseg) && (segs[sl]->last < l) )
		++sl;
	if ( (sf < nseg) && (segs[sf]->last < f) )
		++sf;

	// Chris 12-2005 this is making weird stuff happen
	//if ( (sf == nseg) && (segs[nseg-1]->last == (f-1)) )
	//	sf = sl = nseg-1;
	
	if ( sl == nseg )
	  --sl;
}

void       QUB_IdlSegments::resetSeg(int s)
{
	int f = segs[s]->first;
	int l = segs[s]->last;
	delete segs[s];
	segs[s] = new QUB_IdlSeg(freezer, f, l, sampling);
}

double     QUB_IdlSegments::getSampling()
{
	return sampling;
}

void       QUB_IdlSegments::setSampling(double samp)
{
	sampling = samp;
	for ( int i=0; i<int(segs.size()); ++i )
		segs[i]->sampling = samp;
}


int        QUB_IdlSegments::getSegCount()
{
  return (int) segs.size();
}

int  QUB_IdlSegments::getSegFirst(int i)
{
  return segs[i]->first;
}

int  QUB_IdlSegments::getSegLast(int i)
{
  return segs[i]->last;
}

inline int iround( double x )
{
	return (int) floor( x + 0.5 );
}

int        QUB_IdlSegments::getDwellCount(int f, int l, int fragments)
{
  int sf, sl;
  getSegRange(f, l, sf, sl);
  if ( sf == (int) segs.size() )
    sf = sl = ((int) segs.size()) - 1;
  if ( sf < 0 )
    return 0;
  
  int count = 0;
  for (int i=sf; i<=sl; ++i)
    count += segs[i]->getDwellCount(f, l, fragments);

  return count;
}

int        QUB_IdlSegments::getDwellCountAndGaps(int f, int l)
{
  int sf, sl;
  getSegRange(f, l, sf, sl);
  if ( sf == (int) segs.size() )
    sf = sl = ((int) segs.size()) - 1;
  if ( sf < 0 )
    return 0;
  
  int count = 0;
  for (int i=sf; i<=sl; ++i)
    count += segs[i]->getDwellCountAndGaps(f, l);

  return count;
}

int        QUB_IdlSegments::getDwells(int f, int l, int fragments, int *firsts, int *lasts, int *classes, float *durations)
{
  int sf, sl;
  getSegRange(f, l, sf, sl);
  if ( sf == (int) segs.size() )
    sf = sl = ((int) segs.size()) - 1;
  if ( sf < 0 )
    return 0;
  
  int count = 0;
  for (int i=sf; i<=sl; ++i)
    count += segs[i]->getDwells(f, l, fragments, firsts+count, lasts+count, classes+count, durations?(durations+count):NULL);
  
  if ( count && (firsts[0] < f) ) {
    if ( durations )
      durations[0] -= (float) ((f - firsts[0]) * sampling);
    firsts[0] = f;
  }
  if ( count && (lasts[count-1] > l) ) {
    if ( durations )
      durations[count-1] -= (float) ((lasts[count-1] - l) * sampling);
    lasts[count-1] = l;
  }

  return count;
}

int        QUB_IdlSegments::getDwellsAndGaps(int f, int l, int *firsts, int *lasts, int *classes, float *durations)
{
  int sf, sl;
  getSegRange(f, l, sf, sl);
  if ( sf == (int) segs.size() )
    sf = sl = ((int) segs.size()) - 1;
  if ( sf < 0 )
    return 0;
  
  int count = 0;
  for (int i=sf; i<=sl; ++i)
    count += segs[i]->getDwellsAndGaps(f, l, firsts+count, lasts+count, classes+count, durations?(durations+count):NULL);
  
  if ( count && (firsts[0] < f) ) {
    if ( durations )
      durations[0] -= (float) ((f - firsts[0]) * sampling);
    firsts[0] = f;
  }
  if ( count && (lasts[count-1] > l) ) {
    if ( durations )
      durations[count-1] -= (float) ((lasts[count-1] - l) * sampling);
    lasts[count-1] = l;
  }

  return count;
}

int        QUB_IdlSegments::setDwells(int count, int *ff, int *ll, int *cc)
{
	int s = 0;

	int f = ff[0];
	int l = ll[count-1];

	while ( (s < int(segs.size())) && (segs[s]->last < f) )
		++s;
	while ( (s < int(segs.size())) && (segs[s]->first <= l) ) {
		if ( segs[s]->first >= f && segs[s]->last <= l )
			resetSeg(s);
		segs[s]->setDwells(count, ff, ll, cc);
		++s;
	}

	return 0;
}

// (l - f + 1) == Np :: Nsample;       pps = Np/Nsample;            spp = Nsample/Np;
// s = (p-f) * spp
// p = s*pps + f

# define UNSET_F 262.144f

void QUB_IdlSegments::sampleDwells(int f, int l, double *amp, int Nsample, float *los, float *his, unsigned int *class_bits)
{
  int Np = (l - f + 1);
  double spp = Nsample*1.0/Np;
  iterator ii = find(f);
  int sam = (int) ((ii->first + ii.segfirst() - f) * spp);
  if ( sam < 0 )
    sam = 0;
  int p_next = (int) (sam/spp + f);
  float last_a = UNSET_F;
  for (; sam<Nsample; ++sam) {
    int p = p_next;
    p_next = (int) ((sam+1) / spp + f);
    float lo = last_a;
    float hi = last_a;
    unsigned int clbit = 0;
    int ef, el, ec;
    while ( ii ) {
      ii.read(&ef, &el, &ec);
      // cout << "sam: " << sam << "\t" << ef << ", " << el << ", " << ec << "\tp: " << p << ", " << p_next << endl;
      if ( (ef > p) && (ef >= p_next) ) {
	last_a = UNSET_F;
	break; // starts after this sample
      } else if ( ec >= 0 ) {
	clbit |= (1 << ec);
	float a = (float) amp[ec];
	if ( lo == UNSET_F ) {
	  lo = hi = a;
	} else {
	  if ( a < lo )
	    lo = a;
	  if ( a > hi )
	    hi = a;
	}
	last_a = a;
      }
      else {
	last_a = UNSET_F;
      }
      if ( el <= (p_next-1) ) {
	++ii;
      }
      if ( el >= (p_next-1) ) {
	break;
      }
    }
    if ( class_bits )
      class_bits[sam] = clbit;
    if ( lo != UNSET_F ) {
      los[sam] = lo;
      his[sam] = hi;
    }
  }
}

void QUB_IdlSegments::sampleDwells(int f, int l, float *amp, float *std, int Nsample, float *los, float *his)
{
  int Np = (l - f + 1);
  double spp = Nsample*1.0/Np;
  iterator ii = find(f);
  int sam = (int) ((ii->first + ii.segfirst() - f) * spp);
  if ( sam < 0 )
    sam = 0;
  int p_next = (int) (sam/spp + f);
  float last_a_lo = UNSET_F;
  float last_a_hi = UNSET_F;
  for (; sam<Nsample; ++sam) {
    int p = p_next;
    p_next = (int) ((sam+1) / spp + f);
    float lo = last_a_lo;
    float hi = last_a_hi;
    int ef, el, ec;
    while ( ii ) {
      ii.read(&ef, &el, &ec);
      // cout << "sam: " << sam << "\t" << ef << ", " << el << ", " << ec << "\tp: " << p << ", " << p_next << endl;
      if ( (ef > p) && (ef >= p_next) ) {
	last_a_lo = last_a_hi = UNSET_F;
	break; // starts after this sample
      } else if ( ec >= 0 ) {
	float a = amp[ec];
	float s = std[ec];
	last_a_lo = a - s;
	last_a_hi = a + s;
	if ( lo == UNSET_F ) {
	  lo = last_a_lo;
	  hi = last_a_hi;
	} else {
	  if ( a < lo )
	    lo = last_a_lo;
	  if ( a > hi )
	    hi = last_a_hi;
	}
      }
      else {
	last_a_lo = last_a_hi = UNSET_F;
      }
      if ( el <= (p_next-1) ) {
	++ii;
      }
      if ( el >= (p_next-1) ) {
	break;
      }
    }
    if ( lo != UNSET_F ) {
      los[sam] = lo;
      his[sam] = hi;
    }
  }
}

int        QUB_IdlSegments::clear()
{
	int rslt = 0;
	if ( segs.size() )
		rslt = erase(segs[0]->first, segs[segs.size()-1]->last);
	return rslt;
}

int        QUB_IdlSegments::addSeg(int f, int l)
{
	int sf, sl, s;
	getSegRange(f, l, sf, sl);

	if ( sf >= int(segs.size()) ) {
	        segs.push_back( new QUB_IdlSeg(freezer, f, l, sampling) );
		return 0;
	}

	if ( segs[sf]->first < f ) {
		breakSeg(f);
		++sf;
	}
	
	vector<QUB_IdlSeg*>::iterator segi = segs.begin();
	for ( s=0; s<sf; ++s )
		++segi;
	segs.insert( segi, new QUB_IdlSeg(freezer, f, l, sampling) );

	for ( s=sf+1; s<int(segs.size()); ++s ) {
		segs[s]->offset(l-f+1);
	}

	return 0;
}

int        QUB_IdlSegments::removePoints(int f, int l)
{
	int sf, sl, s;
	getSegRange(f, l, sf, sl);
	int count = l - f + 1;

	vector<QUB_IdlSeg*>::iterator segi = segs.begin();
	for ( s=0; s<sf; ++s )
		++ segi;

	for ( s=sf; s<=sl; ++s, ++segi ) {
		QUB_IdlSeg& seg = * segs[s];
		bool hasBefore = seg.first < f;
		bool hasAfter = l < seg.last;

		if ( hasBefore && hasAfter ) {
			breakSeg(f);
			segs[s+1]->setFirst(l+1);
		}
		else if ( hasBefore ) {
			seg.setLast(f-1);
		}
		else if ( hasAfter ) {
			seg.setFirst(l+1);
			--sl;
		}
		else {
			if ( segs.size() > 1 ) {
				delete segs[s];
				segs.erase(segi);
				if ( s )
					--segi;
				--s;
				--sl;
			}
			else {
				segs[s]->setLast( segs[s]->first - 1 );
			}
		}
	}

	for ( s=sl+1; s<int(segs.size()); ++s ) {
		segs[s]->offset( - count );
	}

	return 0;
}

int        QUB_IdlSegments::insertPoints(int f, int l)
{
	int sf, s;
	getSegRange(f, f, sf, sf);

	int nseg = (int) segs.size();
	if ( (sf == nseg) && (segs[nseg-1]->last == (f-1)) )
		sf = nseg-1;

	int rslt = 0;
	if ( sf == int(segs.size()) )
		rslt = addSeg(f, l);
	else {
		segs[sf]->insert(f, l);
		for ( s=sf+1; s<int(segs.size()); ++s )
			segs[s]->offset(l - f + 1);
	}
	return rslt;
}

int        QUB_IdlSegments::breakSeg(int at)
{
	int sf, s;
	getSegRange(at, at, sf, sf);
	int rslt = 0;
	
	if ( sf == int(segs.size()) )
		rslt = -1;
	else if ( segs[sf]->first == at )
		rslt = 1;
	else {
		vector<QUB_IdlSeg*>::iterator segi = segs.begin();
		for ( s=0; s<sf; ++s )
			++segi;
		
		int segfirst = segs[s]->first;
		segs.insert(segi, new QUB_IdlSeg(freezer, segfirst, at-1, sampling));
		vector<QUB_IdlDwell> dwells;
		if ( true ) {
			QUB_IdlMaps::iterator di = segs[sf+1]->begin();
			while ( di && (di->first + segfirst) <= at ) {
				dwells.push_back( *di );
				++di;
			}
		}
		segs[sf+1]->setFirst( at );
	}
	return rslt;
}

int        QUB_IdlSegments::joinSegs(int at)
{
	int s, sf, sl;
	getSegRange(at-1, at, sf, sl);
	int rslt = 0;

	if ( sf == int(segs.size()) )
		rslt = -1;
	if ( sf == sl )
		rslt = -2;
	else {
		vector<QUB_IdlSeg*>::iterator segi = segs.begin();
		for ( s=0; s<sl; ++s )
			++segi;

		QUB_IdlSeg& segA = * segs[sf];
		int oldApoints = segA.last - segA.first + 1;
		vector<QUB_IdlDwell> dwells;
		if ( true ) {
			QUB_IdlMaps::iterator di = segs[sl]->begin();
			while ( di ) {
				dwells.push_back( *di );
				++di;
			}
		}

		segA.setLast( segs[sl]->last );
		segA.maps.set((int) dwells.size(), &(dwells[0]), - oldApoints);

		delete segs[sl];
		segs.erase(segi);
	}
	return rslt;
}

int        QUB_IdlSegments::erase(int f, int l)
{
	if ( f > l ) return 0;
	int s, sf, sl, ff, ll;
	getSegRange(f, l, sf, sl);
	for ( s=sf; s<=sl; ++s ) {
		QUB_IdlSeg& seg = * segs[s];
		ff = max(f, seg.first);
		ll = min(l, seg.last);
		if ( seg.first == ff && seg.last == ll )
			resetSeg(s);
		else if ( ff <= ll )
			seg.set(ff, ll, -1);
	}
	return 0;
}

int        QUB_IdlSegments::join(int f, int l)
{
	int s, sf, sl, ff, ll;
	getSegRange(f, l, sf, sl);
	for ( s=sf; s<=sl; ++s ) {
		QUB_IdlSeg& seg = * segs[s];
		ff = max(f, seg.first);
		ll = min(l, seg.last);
		if ( ff <= ll )
			seg.join(ff, ll);
	}
	return 0;
}

QUB_IdlSegments::iterator QUB_IdlSegments::begin()
{
	if ( segs.size() )
		return iterator(*this, 0, segs[0]->begin());
	return iterator(*this, -1, QUB_IdlMaps::iterator());
}

QUB_IdlSegments::iterator QUB_IdlSegments::end()
{
	if ( segs.size() )
	  return iterator(*this, ((int) segs.size())-1, segs[segs.size()-1]->end());
	return iterator(*this, -1, QUB_IdlMaps::iterator());
}

QUB_IdlSegments::iterator QUB_IdlSegments::find(int point)
{
	for ( int i=0; i<int(segs.size()); ++i )
		if ( segs[i]->first <= point && point <= segs[i]->last )
			return iterator(*this, i, segs[i]->find(point));
	return iterator(*this, -1, QUB_IdlMaps::iterator());
}