#include <math.h>
#include <grass/gis.h>
#include <grass/Vect.h>
#include "global.h"


/* probability for gaussian distribution */
double gaussian2dBySigma(double d, double sigma)
{
  double res;

  res=1./(2.*M_PI*sigma*sigma)*exp(-d*d/(2.*sigma*sigma));

  return(res);
}


double gaussianFunction(double x, double sigma, double dimension)
{
  return((1./(pow(2.*M_PI,dimension/2.)*pow(sigma,dimension)))*exp(-0.5*pow(x/sigma,2.)));
}


/* probability for gaussian distribution */
double gaussianKernel(double x, double term)
{
  return(term * exp(-(x*x)/2.));
}


/*
    term1 = 1./(2.*M_PI*sigma*sigma)
    term2 = 2.*sigma*sigma;
*/
double gaussian2dByTerms(double d, double term1, double term2)
{
  double res;

  res=term1*exp(-d*d/term2);

  return(res);
}


double segno(double x)
{
  double y; 

  y = (x > 0 ? 1. : 0.) + (x < 0 ? -1. : 0.);

  return y;
}


double kernel1(double d,double rs, double lambda)
{
  double res;
  double a = lambda - 1.;

  if(lambda==1.){
    res=1./(M_PI*(d*d+rs*rs));
  }else{
    res=segno(a)*(a/M_PI)*(pow(rs,2.*a))*(1/pow(d*d+rs*rs,lambda));
  }

    /*  res=1./(M_PI*(d*d+rs*rs));*/
  return(res);
}


double invGaussian2d(double sigma, double prob)
{
  double d;

  d = sqrt(-2*sigma*sigma*log(prob*M_PI*2*sigma*sigma));
  return(d);
}     


/* euclidean distance between vectors x and y of length n */
double euclidean_distance(double *x, double *y, int n)
{
  int j;
  double out = 0.0;
  double tmp;

  for(j=0;j<n;j++){
    tmp = x[j] - y[j];
    out += tmp * tmp;
  }

  return sqrt(out);
}