// ----------------------------------------------------------------------------
//
// Copyright (C) Columbia University, 1994. All Rights Reserved.
// Sameer A. Nene, Shree K. Nayar, Hiroshi Murase
//
// See file LICENSE for details of software license agreement.
//
// ----------------------------------------------------------------------------
//
// binarysearch.C
//
// Author:               Sameer Nene
// Date:                 05/26/94
// Version:              1.0
// Modification History:
//   06/26/94: Changes necessary because of change in interface.
//   05/02/95: Major changes in search(). Added code to select smallest range
//             first.
//
// Bugs:
//
// Classes:
//   BinarySearch
//
// Notes:
//   This module contains implementation of classes declared in binarysearch.h
//
// ----------------------------------------------------------------------------

#include <stdio.h>
#include <stdlib.h>
#include <values.h>
#include "errorscope.h"
#include "registry.h"
#include "persistent.h"
#include "vector.h"
#include "dataset.h"
#include "persscheme.h"
#include "interpolation.h"
#include "binarysearch.h"

void* BinarySearch::d_protocol = (void*)(&BinarySearch::d_protocol);
int BinarySearch::d_pos;
double *BinarySearch::d_sa_p;
char *BinarySearch::d_name_p = "BinarySearch";
Registrar BinarySearch::d_registrar(BinarySearch::d_name_p,
				    &(BinarySearch::createFunc));
Vector BinarySearch::d_pv;				       
				       
void BinarySearch::sample(const Interpolation &ip, int level, int d, const Vector &f)
{
  int i, s = int(f.element(level)) - 1;
  double t, inc;
  Vector samp;

  if(level == f.length() - 1) {
    for(i = 0, t = 0., inc = 1. / s; i < s; ++i, t += inc) {
      d_pv.element(level, t);
      d_sa_p[d_pos++] = (ip.sample(d_pv)).element(d);
    }
    d_pv.element(level, 1.);
    d_sa_p[d_pos++] = (ip.sample(d_pv)).element(d);
  }
  else {
    if(level == 0) {
      d_pv.length(f.length());
      d_pos = 0;
    }
    for(i = 0, t = 0., inc = 1. / s; i < s; ++i, t += inc) {
      d_pv.element(level, t);
      sample(ip, level + 1, d, f);
    }
    d_pv.element(level, 1.);
    sample(ip, level + 1, d, f);
  }
}

Persistent* BinarySearch::createFunc()
{
  return new BinarySearch;
}

int BinarySearch::hasProtocol(void* protocol) const
{
  return d_protocol == protocol || PersistentSearchScheme::hasProtocol(protocol);
}

int BinarySearch::compare(const void *e1, const void *e2)
{
  return d_sa_p[*((int *)e1)] - d_sa_p[*((int *)e2)] < 0 ? -1 : 1;
}

int BinarySearch::bsearchL(int col, double d) const
{
  Vector &v = d_data[col];
  int b = 0, c, t = v.length() - 1;
  double q;
  
  while(t - b > 1) {
    c = (b + t) >> 1;
    q = v.element(c);
    if(d < q)
      t = c;
    else if(d > q)
      b = c;
    else
      return c;
  }

  return ((d <= v.element(b)) ? b : t);
}

int BinarySearch::bsearchR(int col, double d) const
{
  Vector &v = d_data[col];
  int b = 0, c, t = v.length() - 1;
  double q;
  
  while(t - b > 1) {
    c = (b + t) >> 1;
    q = v.element(c);
    if(d < q)
      t = c;
    else if(d > q)
      b = c;
    else
      return c;
  }

  return ((d >= v.element(t)) ? t : b);
}

BinarySearch::BinarySearch() : d_bmaps_p(0), d_fmaps_p(0), d_thresh(.1)
{
}

BinarySearch::BinarySearch(const DataSet &ds, const Vector &rangemin,
			   const Vector &rangemax) : d_data(ds()[0].length()),
			   d_n(ds.getSize()),
			   d_rangemin(rangemin), d_rangemax(rangemax),
			   d_thresh(.1)
{
  int i, j, k, d_i, *p1, *p2;
  VectAry &va = ds();

  for(i = 0, j = ds.getSize(); i < j; ++i)
    d_n.element(i, ds.getSize(i));
  
  d_data[0].length(j = va.getSize());
  d_bmaps_p = new int*[k = va[0].length()];
  d_fmaps_p = new int*[k];
  d_sa_p = new double[j];

  for(d_i = 0; d_i < k; ++d_i) {
    p1 = d_bmaps_p[d_i] = new int[j];
    d_data[d_i].length(j);
    for(i = 0; i < j; ++i) {
      p1[i] = i;
      d_sa_p[i] = va[i].element(d_i);
    }
    qsort(p1, j, sizeof(int), compare);
    for(i = 0; i < j; ++i)
      d_data[d_i].element(i, d_sa_p[p1[i]]);
    p2 = d_fmaps_p[d_i] = new int[j];
    for(i = 0; i < j; ++i)
      p2[p1[i]] = i;
  }
}

BinarySearch::BinarySearch(const Interpolation &ip, const Vector &f,
			   const Vector &n, const Vector &rangemin,
			   const Vector &rangemax) :
			   d_data(ip.spaceDimension()), d_n(n),
			   d_rangemin(rangemin), d_rangemax(rangemax),
			   d_thresh(.1)
{
  int i, j, k, d_i, *p1, *p2;

  for(i = 0, k = f.length(), j = 1; i < k; j *= int(f.element(i++)));
  
  d_data[0].length(j);
  d_bmaps_p = new int*[k = ip.spaceDimension()];
  d_fmaps_p = new int*[k];
  d_sa_p = new double[j];

  for(d_i = 0; d_i < k; ++d_i) {
    p1 = d_bmaps_p[d_i] = new int[j];
    
    d_data[d_i].length(j);

    sample(ip, 0, d_i, f);

    for(i = 0; i < j; ++i)
      p1[i] = i;

    qsort(p1, j, sizeof(int), compare);

    for(i = 0; i < j; ++i)
      d_data[d_i].element(i, d_sa_p[p1[i]]);
    p2 = d_fmaps_p[d_i] = new int[j];
    for(i = 0; i < j; ++i)
      p2[p1[i]] = i;
  }
}

BinarySearch::~BinarySearch()
{
  for(int i = 0, j = d_data.getSize(); i < j; ++i) {
    delete[] d_bmaps_p[i];
    delete[] d_fmaps_p[i];
  }
  delete[] d_bmaps_p;
  delete[] d_fmaps_p;
  delete[] d_sa_p;
}

BinarySearch* BinarySearch::safeCast(Persistent *p)
{
  return p -> hasProtocol(d_protocol) ? (BinarySearch*)p : 0;
}

Vector BinarySearch::search(const Vector &v) const
{
  int *b, *d, i, j, k, l, n, o, *t, *p, u, bot, top, *bmap, *list;
  double m, q, tsq;
  Vector r(d_n.length());

  b = new int[j = d_data.getSize()];
  t = new int[j];
  for(i = 0, k = d_data[0].length() - 1, l = 0; i < j; ++i) {
    b[i] = bsearchL(i, v.element(i) - d_thresh);
    t[i] = bsearchR(i, v.element(i) + d_thresh);
  }

  d = new int[j];
  for(i = 0; i < j; ++i)
    d[i] = i;
  for(i = 0, n = 1; i < j - 1 && n != 0; ++i)
    for(k = n = 0; k < j - i - 1; ++k)
      if(t[d[k + 1]] - b[d[k + 1]] < t[d[k]] - b[d[k]]) {
	l = d[k];
	d[k] = d[k + 1];
	d[k + 1] = l;
	n = 1;
      }

  bmap = d_bmaps_p[d[0]];
  p = d_fmaps_p[d[1]];
  list = new int[t[d[0]] - b[d[0]] + 1];
  
  for(i = b[d[0]], k = t[d[0]], n = 0, bot = b[d[1]], top = t[d[1]]; i <= k; ++i) {
    o = p[l = bmap[i]]; 
    if(o >= bot && o <= top)
      list[n++] = l;
  }

  for(i = 2, j = d_data.getSize(); i < j; ++i) {
    p = d_fmaps_p[d[i]];
    bot = b[d[i]];
    top = t[d[i]];
    for(k = 0, l = n, n = 0; k < l; ++k) {
      o = p[u = list[k]]; 
      if(o >= bot && o <= top)
	list[n++] = u;
    }
  }

  delete[] b;
  delete[] t;
  delete[] d;

  Vector s(j);
  
  m = MAXDOUBLE;
  for(k = 0, l = n, tsq = d_thresh * d_thresh; k < l; ++k) {
    o = list[k];
    for(i = 0; i < j; ++i)
      s.element(i, d_data[i].element(d_fmaps_p[i][o]));
    if((q = s.distance(v)) < m)
      if(q <= tsq) {
	m = q;
	n = list[k];
      }
  }
  
  delete[] list;

  for(i = 0, j = d_n.length(), k = 1; i < j; k *= int(d_n.element(i++)));

  for(i = 0; i < j; ++i) {
    k /= int(d_n.element(i));
    r.element(i, (n / k) * (d_rangemax.element(i) - d_rangemin.element(i)) /
	      (d_n.element(i) - 1.) + d_rangemin.element(i));
    n %= k;
  }
  
  return r;
}

void BinarySearch::setSearchParameters(const Vector &v)
{
  d_thresh = v.element(0);
}

ErrorScope::Error BinarySearch::get(FILE *file)
{
  int i, j, k;
  ErrorScope::Error error;

  if((error = d_n.get(file)) != OK)
    return error;
  
  if((error = d_rangemin.get(file)) != OK)
    return error;
  
  if((error = d_rangemax.get(file)) != OK)
    return error;

  j = d_data.getSize();

  if((error = d_data.get(file)) != OK)
    return error;

  for(i = 0; i < j; ++i) {
    delete[] d_bmaps_p[i];
    delete[] d_fmaps_p[i];
  }
  delete[] d_bmaps_p;
  delete[] d_fmaps_p;
  d_bmaps_p = new int*[j = d_data.getSize()];
  d_fmaps_p = new int*[j = d_data.getSize()];
  k = d_data[0].length();

  for(i = 0; i < j; ++i) {
    d_bmaps_p[i] = new int[k];
    if(fread((char *)d_bmaps_p[i], sizeof(int), k, file) != k)
      return FILE_READ_ERROR;
  }
  
  for(i = 0; i < j; ++i) {
    d_fmaps_p[i] = new int[k];
    if(fread((char *)d_fmaps_p[i], sizeof(int), k, file) != k)
      return FILE_READ_ERROR;
  }
  
  return OK;
}

ErrorScope::Error BinarySearch::put(FILE *file) const
{
  int i, j, k;
  ErrorScope::Error error;
  
  if((error = d_n.put(file)) != OK)
    return error;

  if((error = d_rangemin.put(file)) != OK)
    return error;
  
  if((error = d_rangemax.put(file)) != OK)
    return error;

  if((error = d_data.put(file)) != OK)
    return error;

  k = d_data[0].length();

  for(i = 0, j = d_data.getSize(); i < j; ++i)
    if(fwrite((char *)d_bmaps_p[i], sizeof(int), k, file) != k)
      return FILE_WRITE_ERROR;
  
  for(i = 0, j = d_data.getSize(); i < j; ++i)
    if(fwrite((char *)d_fmaps_p[i], sizeof(int), k, file) != k)
      return FILE_WRITE_ERROR;
  
  return OK;
}

const char* BinarySearch::name() const
{
  return d_name_p;
}

const VectAry& BinarySearch::getData() const
{
  return d_data;
}
