#include <ruby.h>
#include <math.h>
#include <narray.h>
#include "ssl2.h"


scomplex num2scomplex(x)
     VALUE x;
{
  scomplex cx;

  if (rb_obj_is_kind_of(x,cComplex)!=Qtrue) raise(rb_eTypeError,"is not Complex");
  cx.r = NUM2FLT(rb_funcall(x, rb_intern("real"), 0));
  cx.i = NUM2FLT(rb_funcall(x, rb_intern("image"), 0));
  return cx;
}

dcomplex num2complex(x)
     VALUE x;
{
  dcomplex cx;

  if (rb_obj_is_kind_of(x,cComplex)!=Qtrue) raise(rb_eTypeError,"is not Complex");
  cx.r = NUM2DBL(rb_funcall(x, rb_intern("real"), 0));
  cx.i = NUM2DBL(rb_funcall(x, rb_intern("image"), 0));
  return cx;
}


static int ssl2_arysize(x)
     VALUE x;
{
  int n;
  int na;
  int i;
  VALUE *ptr;

  if ( TYPE(x)!=T_ARRAY )
    rb_raise(rb_eTypeError, "x is not array");

  n = RARRAY(x)->len;
  na = n;
  ptr = RARRAY(x)->ptr;
  for ( i=0; i<n; i++ ) {
    if ( TYPE(ptr[i])==T_ARRAY )
      na += ssl2_arysize(ptr[i])-1;
  }
  return na;
}

static int ssl2_aryrank(x)
     VALUE x;
{
  int rank;
  VALUE *ptr;
  int size;
  int n;

  if (TYPE(x)!=T_ARRAY)
    rb_raise(rb_eTypeError, "x in not array");

  n = RARRAY(x)->len;
  ptr = RARRAY(x)->ptr;
  if (TYPE(ptr[0])==T_ARRAY) {
    rank = ssl2_aryrank(ptr[0])+1;
  } else
    rank = 1;

  return rank;
}

static void ssl2_aryshape(x,shape)
     VALUE x; int* shape;
{
  int n;
  VALUE *ptr;

  n = RARRAY(x)->len;
  ptr = RARRAY(x)->ptr;
  shape[0] = n;
  if (TYPE(ptr[0])==T_ARRAY)
      ssl2_aryshape(ptr[0],shape[1]);
}  

static void ssl2_ary_check_form(x)
     VALUE x;
{
  int n;
  VALUE *ptr;
  int size, rank;

  int i;

  n = RARRAY(x)->len;
  ptr = RARRAY(x)->ptr;
  if (TYPE(ptr[0])==T_ARRAY) {
    size = ssl2_arysize(ptr[0]);
    rank = ssl2_aryrank(ptr[0]);
    for (i=1;i<n;i++){
      if (TYPE(ptr[i])!=T_ARRAY)
	rb_raise(rb_eTypeError, "x is not well formed array (type)");
      if (ssl2_arysize(ptr[i])!=size)
	rb_raise(rb_eRuntimeError, "x is not well formed array (size)");
      if (ssl2_aryrank(ptr[i])!=rank)
	rb_raise(rb_eRuntimeError, "x is not well formed array (rank)");
    }
    for (i=0;i<n;i++)
      ssl2_ary_check_form(ptr[i]);
  }
}

static float *ssl2_arysetvalue(x,cx)
     VALUE x; float *cx;
{
  VALUE *ptr;
  int n;
  int i;

  if ( TYPE(x)!=T_ARRAY )
    rb_raise(rb_eTypeError, "x is not array");

  n = RARRAY(x)->len;
  ptr = RARRAY(x)->ptr;
  for ( i=0; i<n; i++ ) {
    if ( TYPE(ptr[i])==T_ARRAY )
      cx = ssl2_arysetvalue(ptr[i],cx);
    else if ( TYPE(rb_Float(ptr[i]))==T_FLOAT ) {
      *cx = (float) RFLOAT(rb_Float(ptr[i]))->value;
      cx++;
    }
    else
      rb_raise(rb_eTypeError, "value is invalid type");
  }
  return cx;
}


static int *ssl2_arysetvaluei(x,cx)
     VALUE x; int *cx;
{
  VALUE *ptr;
  int n;
  int i;

  if ( TYPE(x)!=T_ARRAY )
    rb_raise(rb_eTypeError, "x is not array");

  n = RARRAY(x)->len;
  ptr = RARRAY(x)->ptr;
  for ( i=0; i<n; i++ ) {
    if ( TYPE(ptr[i])==T_ARRAY )
      cx = ssl2_arysetvaluei(ptr[i],cx);
    else if ( TYPE(rb_Integer(ptr[i]))==T_FIXNUM || TYPE(rb_Integer(ptr[i]))==T_BIGNUM ) {
      *cx = (int)NUM2INT(rb_Integer(ptr[i]));
      cx++;
    }
    else
      rb_raise(rb_eTypeError, "value is invalid type");
  }
  return cx;
}

static long *ssl2_arysetvaluel(x,cx)
     VALUE x; long *cx;
{
  VALUE *ptr;
  int n;
  int i;

  if ( TYPE(x)!=T_ARRAY )
    rb_raise(rb_eTypeError, "x is not array");

  n = RARRAY(x)->len;
  ptr = RARRAY(x)->ptr;
  for ( i=0; i<n; i++ ) {
    if ( TYPE(ptr[i])==T_ARRAY )
      cx = ssl2_arysetvaluel(ptr[i],cx);
    else if ( TYPE(rb_Integer(ptr[i]))==T_FIXNUM || TYPE(rb_Integer(ptr[i]))==T_BIGNUM ) {
      *cx = (long)NUM2INT(rb_Integer(ptr[i]));
      cx++;
    }
    else
      rb_raise(rb_eTypeError, "value is invalid type");
  }
  return cx;
}

static scomplex *ssl2_arysetvaluec(x,cx)
     VALUE x; scomplex *cx;
{
  VALUE *ptr,*child;
  scomplex c;
  int n;
  int i;

  if ( TYPE(x)!=T_ARRAY )
    rb_raise(rb_eTypeError, "x is not array");

  n = RARRAY(x)->len;
  ptr = RARRAY(x)->ptr;
  for ( i=0; i<n; i++ ) {
    child = RARRAY(ptr[i])->ptr;
    if ( TYPE(child[0])==T_ARRAY )
      cx = ssl2_arysetvaluec(ptr[i],cx);
    else if ( (RARRAY(ptr[i])->len)==2 && TYPE(rb_Float(child[0]))==T_FLOAT && TYPE(rb_Float(child[1]))==T_FLOAT ) {
      c.r = (float) RFLOAT(rb_Float(child[0]))->value;
      c.i = (float) RFLOAT(rb_Float(child[1]))->value;
      *cx = c;
      cx++;
    }
    else
      rb_raise(rb_eTypeError, "value is invalid type");
  }
  return cx;
}

static float *ssl2_ary2cary(x,pn)
     VALUE x; int *pn;
{
  float *cx;
  int n;
  int i,j;

  n = ssl2_arysize(x);
  *pn = n;
  cx = (float *) ALLOC_N(float,n);
  ssl2_arysetvalue(x,cx);
  return cx;
}

static int *ssl2_ary2caryi(x,pn)
     VALUE x; int *pn;
{
  int *cx;
  int n;
  int i,j;

  n = ssl2_arysize(x);
  *pn = n;
  cx = (int *) ALLOC_N(int,n);
  ssl2_arysetvaluei(x,cx);
  return cx;
}

static long *ssl2_ary2caryl(x,pn)
     VALUE x; int *pn;
{
  long *cx;
  int n;
  int i,j;

  n = ssl2_arysize(x);
  *pn = n;
  cx = (long *) ALLOC_N(long,n);
  ssl2_arysetvaluel(x,cx);
  return cx;
}

static scomplex *ssl2_ary2caryc(x,pn)
     VALUE x; int *pn;
{
  scomplex *cx;
  int n;
  int i,j;

  n = ssl2_arysize(x)/2;
  *pn = n;
  cx = (scomplex *) ALLOC_N(scomplex,n);
  ssl2_arysetvaluec(x,cx);
  return cx;
}


static float *ssl2_ary2cary_withshape(x,shape,rank)
     VALUE x; int **shape, *rank;
{
  float *cx;
  int pn;

  ssl2_ary_check_form(x);
  *rank = ssl2_aryrank(x);
  *shape = ALLOC_N(int,*rank);
  ssl2_aryshape(x,*shape);
  cx = ssl2_ary2cary(x,&pn);
  return cx;
}

static scomplex *ssl2_ary2caryc_withshape(x,shape,rank)
     VALUE x; int **shape, *rank;
{
  scomplex *cx;
  int pn;
  int *shape2;
  int i;

  ssl2_ary_check_form(x);
  *rank = ssl2_aryrank(x);
  shape2 = ALLOC_N(int,*rank);
  ssl2_aryshape(x,shape2);
  *rank = *rank-1;
  *shape = ALLOC_N(int,*rank);
  for (i=0;i<*rank;i++)
    *shape[i] = shape2[i+1];
  free(shape2);
  cx = ssl2_ary2caryc(x,&pn);
  return cx;
}


#ifdef HAVE_NARRAY_H
static float *ssl2_nary2cary(nx,pn)
     VALUE nx; int *pn;
{
  VALUE nxs;
  float *cx;
  struct NARRAY *nax;
  float *ptr;
  int n;
  int i;

  nxs = na_cast_object(nx, NA_SFLOAT);
  GetNArray(nxs,nax);
  n = nax->total;
  *pn = n;
  ptr = (float *)NA_PTR(nax,0);
  cx = ALLOC_N(float,n);
  for ( i=0; i<n; i++ ) {
    cx[i] = ptr[i];
  }
  return cx;
}


static int *ssl2_nary2caryi(nx,pn)
     VALUE nx; int *pn;
{
  VALUE nxi;
  int *cx;
  struct NARRAY *nax;
  int *ptr;
  int n;
  int i;

  nxi = na_cast_object(nx, NA_SINT);
  GetNArray(nxi,nax);
  n = nax->total; 
  *pn = n;
  ptr = (int *)NA_PTR(nax,0);
  cx = ALLOC_N(int,n);
  for ( i=0; i<n; i++ ) {
    cx[i] = ptr[i];
  }
  return cx;
}


static scomplex *ssl2_nary2caryc(nx,pn)
     VALUE nx; int *pn;
{
  VALUE nxc;
  scomplex *cx;
  struct NARRAY *nax;
  scomplex *ptr;
  int n;
  int i;

  nxc = na_cast_object(nx, NA_SCOMPLEX);
  GetNArray(nxc,nax);
  n = nax->total; 
  *pn = n;
  ptr = (scomplex *)NA_PTR(nax,0);
  cx = ALLOC_N(scomplex,n);
  for ( i=0; i<n; i++ ) {
    cx[i] = ptr[i];
  }
  return cx;
}


static float *ssl2_nary2cary_withshape(nx,shape,rank)
     VALUE nx; int **shape, *rank;
{
  VALUE nxs;
  float *cx;
  int pn;
  struct NARRAY *nax;
  int i;

  nxs = na_cast_object(nx, NA_SFLOAT);
  GetNArray(nxs,nax);
  *rank = nax->rank;
  *shape = ALLOC_N(int,*rank);
  for (i=0;i<*rank;i++)
    (*shape)[i] = nax->shape[i];
  cx = ssl2_nary2cary(nx,&pn);
  return cx;
}


static scomplex *ssl2_nary2caryc_withshape(nx,shape,rank)
     VALUE nx; int **shape, *rank;
{
  VALUE nxs;
  scomplex *cx;
  int pn;
  struct NARRAY *nax;
  int i;

  nxs = na_cast_object(nx, NA_SFLOAT);
  GetNArray(nxs,nax);
  *rank = nax->rank;
  *shape = ALLOC_N(int,*rank);
  for (i=0;i<*rank;i++)
    (*shape)[i] = nax->shape[i];
  cx = ssl2_nary2caryc(nx,&pn);
  return cx;
}
#endif


static VALUE ssl2_cary2ary(cx,dim,shape)
     float **cx; int dim; int *shape;
{
  VALUE ary;
  int i;

  if ( dim==0 ){
    ary = rb_float_new((double)*(*cx)++);
  }
  else {
    ary = rb_ary_new2(shape[dim-1]);
    for ( i=0;i<shape[dim-1];i++ )
      RARRAY(ary)->ptr[i] = ssl2_cary2ary(cx,dim-1,shape);
    RARRAY(ary)->len = shape[dim-1];
  }
  return ary;
}


static VALUE ssl2_cary2aryi(cx,dim,shape)
     int **cx; int dim; int *shape;
{
  VALUE ary;
  int i;

  if ( dim==0 ){
    ary = INT2NUM(*(*cx)++);
  }
  else {
    ary = rb_ary_new2(shape[dim-1]);
    for ( i=0;i<shape[dim-1];i++ )
      RARRAY(ary)->ptr[i] = ssl2_cary2aryi(cx,dim-1,shape);
    RARRAY(ary)->len = shape[dim-1];
  }
  return ary;
}


static VALUE ssl2_cary2aryl(cx,dim,shape)
     long **cx; int dim; int *shape;
{
  VALUE ary;
  int i;

  if ( dim==0 ){
    ary = LONG2NUM(*(*cx)++);
  }
  else {
    ary = rb_ary_new2(shape[dim-1]);
    for ( i=0;i<shape[dim-1];i++ )
      RARRAY(ary)->ptr[i] = ssl2_cary2aryl(cx,dim-1,shape);
    RARRAY(ary)->len = shape[dim-1];
  }
  return ary;
}


static VALUE ssl2_cary2aryc(cx,dim,shape)
     scomplex **cx; int dim; int *shape;
{
  VALUE ary;
  int i;

  if ( dim==0 ){
    ary = rb_ary_new3(2,(*(*cx)).r,(*(*cx)).i);
    (*cx)++;
  }
  else {
    ary = rb_ary_new2(shape[dim-1]);
    for ( i=0;i<shape[dim-1];i++ )
      RARRAY(ary)->ptr[i] = ssl2_cary2aryc(cx,dim-1,shape);
    RARRAY(ary)->len = shape[dim-1];
  }
  return ary;
}


#ifdef HAVE_NARRAY_H
static VALUE ssl2_cary2nary(cx,dim,shape)
     float *cx; int dim; int *shape;
{
  float *ptr;
  int total;
  struct NARRAY *nax;
  VALUE nx;
  int i;

  nx = na_make_object(NA_SFLOAT,dim,shape,cNArray);
  total = 1;
  for ( i=0; i<dim; i++ )
    total *= shape[i];
  GetNArray(nx,nax);
  ptr = (float *) NA_PTR(nax,0);
  for ( i=0; i<total; i++ ) {
    ptr[i] = cx[i];
  }
  return nx;
}


static VALUE ssl2_cary2naryi(cx,dim,shape)
     int *cx; int dim; int *shape;
{
  int *ptr;
  int total;
  struct NARRAY *nax;
  VALUE nx;
  int i;

  nx = na_make_object(NA_LINT,dim,shape,cNArray);
  total = 1;
  for ( i=0; i<dim; i++ )
    total *= shape[i];
  GetNArray(nx,nax);
  ptr = (int *) NA_PTR(nax,0);
  for ( i=0; i<total; i++ ) {
    ptr[i] = cx[i];
  }
  return nx;
}


static VALUE ssl2_cary2naryl(cx,dim,shape)
     long *cx; int dim; int *shape;
{
  long *ptr;
  int total;
  struct NARRAY *nax;
  VALUE nx;
  int i;

  nx = na_make_object(NA_SINT,dim,shape,cNArray);
  total = 1;
  for ( i=0; i<dim; i++ )
    total *= shape[i];
  GetNArray(nx,nax);
  ptr = (long *) NA_PTR(nax,0);
  for ( i=0; i<total; i++ ) {
    ptr[i] = cx[i];
  }
  return nx;
}

static VALUE ssl2_cary2naryc(cx,dim,shape)
     scomplex *cx; int dim; int *shape;
{
  scomplex *ptr;
  int total;
  struct NARRAY *nax;
  VALUE nx;
  int i;

  nx = na_make_object(NA_SCOMPLEX,dim,shape,cNArray);
  total = 1;
  for ( i=0; i<dim; i++ )
    total *= shape[i];
  GetNArray(nx,nax);
  ptr = (scomplex *) NA_PTR(nax,0);
  for ( i=0; i<total; i++ ) {
    ptr[i] = cx[i];
  }
  return nx;
}
#endif


float *ssl2_getcary(x,pn)
     VALUE x; int *pn;
{
  float *cx;

  if (TYPE(x) == T_FLOAT)
    x = rb_Array(x);

  if ( TYPE(x) == T_ARRAY )
    cx = ssl2_ary2cary(x,pn);
#ifdef HAVE_NARRAY_H
  else if ( rb_obj_is_kind_of(x,cNArray)== Qtrue)
    cx = ssl2_nary2cary(x,pn);
#endif
  else
    rb_raise(rb_eTypeError, "x is invalid type");
  return cx;
}


int *ssl2_getcaryi(x,pn)
     VALUE x; int *pn;
{
  int *cx;

  if ( TYPE(x)==T_FIXNUM || TYPE(x)==T_BIGNUM )
    x = rb_Array(x);

  if ( TYPE(x)==T_ARRAY )
    cx = ssl2_ary2caryi(x,pn);
#ifdef HAVE_NARRAY_H
  else if ( rb_obj_is_kind_of(x,cNArray)== Qtrue)
    cx = ssl2_nary2caryi(x,pn);
#endif
  else
    rb_raise(rb_eTypeError, "x is invalid type");
  return cx;
}


scomplex *ssl2_getcaryc(x,pn)
     VALUE x; int *pn;
{
  scomplex *cx;

  if (TYPE(x) == T_FLOAT)
    x = rb_Array(x);

  if ( TYPE(x) == T_ARRAY )
    cx = ssl2_ary2caryc(x,pn);
#ifdef HAVE_NARRAY_H
  else if ( rb_obj_is_kind_of(x,cNArray)== Qtrue)
    cx = ssl2_nary2caryc(x,pn);
#endif
  else
    rb_raise(rb_eTypeError, "x is invalid type");
  return cx;
}


float *ssl2_getcary_withshape(x,shape,rank)
     VALUE x; int **shape, *rank;
{
  float *cx;

  if (TYPE(x) == T_FLOAT)
    x = rb_Array(x);

  if ( TYPE(x) == T_ARRAY )
    cx = ssl2_ary2cary_withshape(x,shape,rank);
#ifdef HAVE_NARRAY_H
  else if ( rb_obj_is_kind_of(x,cNArray)== Qtrue)
    cx = ssl2_nary2cary_withshape(x,shape,rank);
#endif
  else
    rb_raise(rb_eTypeError, "x is invalid type");
  return cx;
}

scomplex *ssl2_getcaryc_withshape(x,shape,rank)
     VALUE x; int **shape, *rank;
{
  scomplex *cx;

  if (TYPE(x) == T_FLOAT)
    x = rb_Array(x);

  if ( TYPE(x) == T_ARRAY )
    cx = ssl2_ary2caryc_withshape(x,shape,rank);
#ifdef HAVE_NARRAY_H
  else if ( rb_obj_is_kind_of(x,cNArray)== Qtrue)
    cx = ssl2_nary2caryc_withshape(x,shape,rank);
#endif
  else
    rb_raise(rb_eTypeError, "x is invalid type");
  return cx;
}


VALUE ssl2_getary(cx,dim,shape)
     float *cx; int dim; int *shape;
{
#ifdef HAVE_NARRAY_H
  return ssl2_cary2nary(cx,dim,shape);
#else
  return ssl2_cary2ary(&cx,dim,shape);
#endif
}

VALUE ssl2_getaryi(cx,dim,shape)
     int *cx; int dim; int *shape;
{
#ifdef HAVE_NARRAY_H
  return ssl2_cary2naryi(cx,dim,shape);
#else
  return ssl2_cary2aryi(&cx,dim,shape);
#endif
}

VALUE ssl2_getaryl(cx,dim,shape)
     long *cx; int dim; int *shape;
{
#ifdef HAVE_NARRAY_H
  return ssl2_cary2naryl(cx,dim,shape);
#else
  return ssl2_cary2aryl(&cx,dim,shape);
#endif
}

VALUE ssl2_getaryc(cx,dim,shape)
     scomplex *cx; int dim; int *shape;
{
#ifdef HAVE_NARRAY_H
  return ssl2_cary2naryc(cx,dim,shape);
#else
  return ssl2_cary2aryc(&cx,dim,shape);
#endif
}


float ssl2_func(cx)
     float *cx;
{
  return (float)NUM2DBL(rb_yield(rb_float_new((double)*cx)));
}

void ssl2_sub(cx,cy)
     float *cx, *cy;
{
  VALUE y;
  int m;
  y = rb_yield(rb_float_new((double)*cx));
  cy = ssl2_getcary(y,&m);
}

void ssl2_error(icon)
     int icon;
{
  if ( icon!=0 ) raise(rb_eRuntimeError, "icon is %d\n", icon);
}

int ssl2_where(x,y,n)
     float x, *y; int n;
{
  int i,j;

  j=n/2;
  if (n>1)
    if (x<y[j])
      i = ssl2_where(x,y,j);
    else {
      i = ssl2_where(x,y+j,n-j);
      if (i!=-1) i += j;
    }
  else
    i=0;

  return i;
}