#include "tjcstd.h"
#include "tjutils.h"
#include "tjstring.h"

#ifdef USING_WIN32
#include <direct.h>
#include <windows.h>
#include <io.h>
#endif


#ifndef HAVE_CTYPE_H
int isdigit(int c) {
  if(c>='0' && c<='9') return true;
  else return false;
}

int isalpha(int c) {
  return islower(c) || isupper(c);
}

int isalnum(int c) {
  return isdigit(c) || isalpha(c);
}

int isspace(int c) {
  if(c==' ' || c=='\n' || c=='\r' || c=='\t') return true;
  else return false;
}

int islower(int c) {
  return (c>='a' && c<='z');
}

int tolower(int c) {
  if(isupper(c)) c+=32;
  return c;
}

int isupper(int c) {
  return (c>='A' && c<='Z');
}

int toupper(int c) {
  if(islower(c)) c-=32;
  return c;
}

#endif

#ifndef HAVE_J1
double j1(double x) {
  double xa=fabs(x);

  double result=0.0;

  static double a1[] = {
       0.1171875,
      -0.1441955566406250,
       0.6765925884246826,
      -6.883914268109947,
       1.215978918765359e2,
      -3.302272294480852e3,
       1.276412726461746e5,
      -6.656367718817688e6,
       4.502786003050393e8,
      -3.833857520742790e10,
       4.011838599133198e12,
      -5.060568503314727e14,
       7.572616461117958e16,
      -1.326257285320556e19};
  static double b1[] = {
      -0.1025390625,
       0.2775764465332031,
      -1.993531733751297,
       2.724882731126854e1,
      -6.038440767050702e2,
       1.971837591223663e4,
      -8.902978767070678e5,
       5.310411010968522e7,
      -4.043620325107754e9,
       3.827011346598605e11,
      -4.406481417852278e13,
       6.065091351222699e15,
      -9.833883876590679e17,
       1.855045211579828e20};

  if(xa==0.0) return result;

  double x2=xa*xa;
  if(xa<=12.0) { // simple series expansion for small numbers
    result=1.0;
    double r=1.0;
    for(int i=1;i<=30;i++) {
      r *= -0.25*x2/(i*(i+1));
      result += r;
      if (fabs(r) < fabs(result)*1e-15) break;
    }
    result *= 0.5*xa;
  } else {
    int n;
    if (xa >= 50.0) n = 8;
    else if (xa >= 35.0) n = 10;
    else n = 12;
    double t2 = xa-0.75*PII;
    double p1 = 1.0;
    double q1 = 0.375/xa;
    for (int i=0; i<n; i++) {
      p1 += a1[i]*pow(xa,-2*i-2);
      q1 += b1[i]*pow(xa,-2*i-3);
    }
    result = sqrt(2.0/PII/xa)*(p1*cos(t2)-q1*sin(t2));
  }
  if(x<0.0) return -result;
  return result;
}
#endif





#ifndef HAVE_ACOSH

// copy & pasted from the GSL library
#define GSL_SQRT_DBL_EPSILON   1.4901161193847656e-08
#define GSL_M_LN2      0.69314718055994530941723212146

double gsl_log1p (const double x)
{
  volatile double y;
  y = 1 + x;
  return log(y) - ((y-1)-x)/y ;  /* cancels errors with IEEE arithmetic */
}


double acosh(double x) {

  if (x > 1.0 / GSL_SQRT_DBL_EPSILON)
    {
      return log (x) + GSL_M_LN2;
    }
  else if (x > 2)
    {
      return log (2 * x - 1 / (sqrt (x * x - 1) + x));
    }
  else if (x > 1)
    {
      double t = x - 1;
      return gsl_log1p (t + sqrt (2 * t + t * t));
    }
  else if (x == 1)
    {
      return 0;
    }
  else
    {
      return 0; // GSL_NAN;
    }
}
#endif



#ifdef VXWORKS
double exp4vxworks(double x) {
  if(x<-50.0) return 0.0; // otherwise the C-library on VxWorks throws an exception (no comment ;-)
#undef exp
  else return exp(x);
#define exp exp4vxworks
}
#endif




////////////////////////////////////////////////////////////////////
#ifndef HAVE_DL

void *dlopen (const char *filename, int flag) {
#ifdef USING_WIN32
  HINSTANCE inst=LoadLibraryA(filename);
  return (void*)inst;
#else
  return 0;
#endif
}

const char *dlerror(void) {
  return lasterr();
}

void *dlsym(void *handle, char *symbol) {
#ifdef USING_WIN32
  FARPROC fp=GetProcAddress((HMODULE)handle,symbol);
  return (void*)fp;
#else
  return 0;
#endif
}


int dlclose (void *handle) {
#ifdef USING_WIN32
  BOOL b=FreeLibrary((HMODULE)handle);
  return (int)!b;
#else
  return 0;
#endif
}

#endif