smooth.c File Reference

Functions for smoothing data and removing spikes from data arrays. More...

#include "mdb.h"

Go to the source code of this file.

Functions
void	smoothData (double *data, long rows, long smoothPoints, long smoothPasses)
	Smooth a data array using a moving average.
long	despikeData (double *data, long rows, long neighbors, long passes, long averageOf, double threshold, long countLimit)
	Remove spikes from a data array by comparing each point to its neighbors.

Detailed Description

Functions for smoothing data and removing spikes from data arrays.

This file provides two main functions:

• smoothData(): Smooths a data set using a simple moving average over a defined number of points and passes.
• despikeData(): Attempts to remove spike values from a data set by comparing each point to its neighbors and replacing it if it exceeds a defined threshold.

Copyright

• (c) 2002 The University of Chicago, as Operator of Argonne National Laboratory.
• (c) 2002 The Regents of the University of California, as Operator of Los Alamos National Laboratory.

License

This file is distributed under the terms of the Software License Agreement found in the file LICENSE included with this distribution.

Author

M. Borland, R. Soliday, H. Shang

Definition in file smooth.c.

Function Documentation

despikeData()

long despikeData	(	double *	data,
		long	rows,
		long	neighbors,
		long	passes,
		long	averageOf,
		double	threshold,
		long	countLimit )

Remove spikes from a data array by comparing each point to its neighbors.

This function identifies and replaces "spikes" in a data array. A spike is defined as a value that differs significantly from its neighboring values. The function compares each point to its neighbors over multiple passes and replaces values exceeding a given threshold with an average of their neighbors. The process can be halted if too many spikes are found (based on countLimit).

Parameters

data	Pointer to the data array.
rows	The number of data points.
neighbors	The number of neighboring points to consider around each point.
passes	The maximum number of passes to attempt when removing spikes.
averageOf	The number of points to average when replacing a spiked value.
threshold	The threshold for determining if a value is a spike.
countLimit	The maximum number of spikes allowed before the process stops.

Returns

The number of spikes removed on the final pass.

Definition at line 86 of file smooth.c.

                                                    {
  long i0, i1, i2, i, j, i1a, i2a;
  int64_t imin, imax;
  double *deltaSum, sum, *tempdata;
  long despikeCount;
 
  neighbors = 2 * ((neighbors + 1) / 2);
  if (!(tempdata = (double *)malloc(sizeof(*tempdata) * rows)) ||
      !(deltaSum = (double *)malloc(sizeof(*deltaSum) * (neighbors + 1))))
    bomb("despikeData: memory allocation failure", NULL);
  memcpy(tempdata, data, sizeof(*data) * rows);
  despikeCount = 0;
  while (passes-- > 0) {
    despikeCount = 0;
    for (i0 = 0; i0 < rows; i0 += neighbors / 2) {
      i1 = i0 - neighbors / 2;
      i2 = i0 + neighbors / 2;
      if (i1 < 0)
        i1 = 0;
      if (i2 >= rows)
        i2 = rows - 1;
      if (i2 - i1 == 0)
        continue;
      for (i = i1; i <= i2; i++) {
        deltaSum[i - i1] = 0;
        for (j = i1; j <= i2; j++)
          deltaSum[i - i1] += fabs(tempdata[i] - tempdata[j]);
      }
      if (index_min_max(&imin, &imax, deltaSum, i2 - i1 + 1)) {
        if ((imax += i1) < 0 || imax > rows) {
          fprintf(stderr, "Error: index out of range in despikeData (sddssmooth)\n");
          fprintf(stderr, "  imax = %" PRId64 ", rows=%ld, i1=%ld, i2=%ld, neighbors=%ld\n",
                  imax - 1, rows, i1, i2, neighbors);
          exit(1);
        }
        if (threshold == 0 || threshold * neighbors < deltaSum[imax - i1]) {
          if ((i1a = imax - averageOf / 2) < 0)
            i1a = 0;
          if ((i2a = imax + averageOf / 2) >= rows)
            i2a = rows - 1;
          for (sum = 0, i = i1a; i <= i2a; i++) {
            if (i != imax)
              sum += tempdata[i];
          }
          despikeCount++;
          tempdata[imax] = sum / (i2a - i1a);
        }
      }
    }
    if (!despikeCount || (countLimit > 0 && despikeCount > countLimit))
      break;
  }
  if (countLimit <= 0 || despikeCount < countLimit)
    for (i = 0; i < rows; i++)
      data[i] = tempdata[i];
  else
    despikeCount = 0;
  free(tempdata);
  free(deltaSum);
  return despikeCount;
}

smoothData()

void smoothData	(	double *	data,
		long	rows,
		long	smoothPoints,
		long	smoothPasses )

Smooth a data array using a moving average.

This function applies a moving average smoothing operation to an array of data points a specified number of times. It uses a window defined by smoothPoints, averaging over this number of points and performing the operation for smoothPasses iterations.

Parameters

data	Pointer to the array of data to be smoothed.
rows	The number of data points.
smoothPoints	The number of points to include in the smoothing window.
smoothPasses	The number of smoothing passes to perform.

Definition at line 35 of file smooth.c.

                                                                               {
  long lower, upper, row, pass, smoothPoints2, terms;
  double sum;
  static double *smoothedData = NULL;
 
  smoothedData = trealloc(smoothedData, rows * sizeof(*smoothedData));
 
  smoothPoints2 = smoothPoints / 2;
 
  for (pass = 0; pass < smoothPasses; pass++) {
    for (row = sum = 0; row < smoothPoints2; row++)
      sum += data[row];
 
    terms = row;
    lower = -smoothPoints2;
    upper = smoothPoints2;
    for (row = 0; row < rows; row++, lower++, upper++) {
      if (upper < rows) {
        sum += data[upper];
        terms += 1;
      }
      smoothedData[row] = sum / terms;
      if (lower >= 0) {
        sum -= data[lower];
        terms -= 1;
      }
    }
 
    for (row = 0; row < rows; row++)
      data[row] = smoothedData[row];
  }
}