# include "data.h"

int vsmooth (int argc, char **argv) {
  
  int i, n, N, Nx, Ns, Ngauss;
  float *vi, *vo, *gauss, *gaussnorm;
  float g, s, sigma, Nsigma;
  Vector *in;

  Nsigma = 3;
  if ((N = get_argument (argc, argv, "-Nsigma"))) {
    remove_argument (N, &argc, argv);
    Nsigma = atof (argv[N]);
    remove_argument (N, &argc, argv);
  }

  if (argc != 3) {
    gprint (GP_ERR, "USAGE: vsmooth (input) sigma\n");
    return (FALSE);
  }
  
  if ((in  = SelectVector (argv[1], OLDVECTOR, TRUE)) == NULL) return (FALSE);
  sigma = atof (argv[2]);
  Nx = in[0].Nelements;
  vi = in[0].elements;

  /* build a 1D gaussian */
  Ns = (int) (Nsigma*sigma + 0.5);
  Ngauss = 2*Ns + 1;
  ALLOCATE (gaussnorm, float, Ngauss);
  gauss = &gaussnorm[Ns];
  for (i = -Ns; i < Ns + 1; i++) {
    gauss[i] = exp ((i*i)/(-2*sigma*sigma));
  }

  ALLOCATE (vo, float, Nx);

  for (i = 0; i < Nx; i++) {
    g = s = 0;
    for (n = -Ns; n < Ns + 1; n++) {
      if (i+n < 0) continue;
      if (i+n >= Nx) continue;
      s += gauss[n]*vi[i+n];
      g += gauss[n];
    }
    vo[i] = s / g;
  }

  free (in[0].elements);
  free (gaussnorm);

  in[0].elements = vo;
  return (TRUE);
}

