# Non-Uniform 1D-Gaussian Filter
gaussfilt <- function(xdata, ydata, sigma, x = xdata, M = 3*sigma, extrap = 0)
{
  y <- rep(extrap, length(x))
  for (i in 1:length(y))
  {
    indices <- which(xdata > (x[i]-M) & xdata < (x[i]+M))
    if (length(indices)>0)
    {
      gaussWeights <- gausskernel(xdata[indices]-x[i], sigma)
      y[i] <- sum(gaussWeights*ydata[indices])/sum(gaussWeights)
    }
  }
  return(y)
}

gausskernel <- function(x, sigma)
{
  return(1/(sigma * sqrt(2*pi)) * exp(-0.5 * (x/sigma)^2))
}