ftp.nice.ch/pub/next/unix/audio/cmix.s.tar.gz#/cmix/track/median.c

This is median.c in view mode; [Download] [Up]

#include <math.h>
#include <stdio.h>
/*
#include <carl/carl.h>
#include <carl/sndio.h>
#include <carl/defaults.h>
*/


/*-------------------------------------------------------

			median.c

This program implements a median filter.  The median filter
	is a nonlinear filter; its output at each point in
	time is the median of the N most recent input values.
	The attraction of such a filter is that - unlike a linear
	filter - it can remove glitches while still following
	step changes perfectly.  On the other hand, the median
	filter does introduce a signal-dependent distortion
	which depends also on the value of N.  

	cc -O median.c -lcarl -lm -o median

-------------------------------------------------------*/


main(argc, argv)
	int argc; char **argv;
{

	int	i,		/* pointer to current sample in buffer */
		j,
		N = 9,		/* number of sample to mediate over */
		N2;		/* N / 2 */

	float	*buf,		/* buffer for N most recent samples */
		min,		/* temporary minimum */
		max,		/* temporary maximum */
		med = 0.,	/* median */
		new,		/* current input sample */
		old,		/* sample being replaced */
		srate = 1.;	/* sample rate */

	char	ch;

	char	*dbuf; 		/* buffer for strings read from header */

	PROP	*proplist;	/* from header on stdin */

	float	dummy;		/* pl, to throw away rmsamps */

/* Read header from stdin. */

	if (isatty(0))
		usage(1);

	if ((proplist =  getheader(stdin)) != NULL) {	/* there is a header */

		if ((dbuf = getprop(stdin,H_SRATE)) != NULL) {
			srate = atof(dbuf);		/* get input srate */
		}
	}

/* Call crack to interpret commandline. */

	while ((ch = crack(argc, argv, "R|", 1)) != NULL){
		switch (ch){
			case 'R': srate = sfexpr(arg_option, 1.0); break;
		}
	}

	arg_index = 0;	/* reset for second call to crack() */

	while ((ch = crack(argc, argv, "R|w|m|h", 0)) != NULL){
		switch (ch){ 
			case 'w': 	N = sfexpr(arg_option, srate); break;
			case 'R': 	break;
			case 'h':	usage(0);
			default: 	usage(1);
		}
	}

	if (N%2 == 0)
		N--;
	N2 = (N-1) / 2;
printf("%d %d\n",N,N2);
/* Initialization. */

	if ((buf = (float *) calloc(N,sizeof(float))) == NULL)
		malerr("median: insufficient memory",1);

	for (i = 0; i < N2; i++){
		if (getfloat(&new) <= 0) new = 0;
		getfloat(&dummy);
		buf[i] = new;
		buf[N-i-1] = new;
		if (new > med){
			min = 1e6;
			for (j = 0; j <= i; j++)
				if ((buf[j] >= med) && (buf[j] < min))
					min = buf[j];
			med = min;
		}
		if (new < med){
			max = -1e6;
			for (j = 0; j <= i; j++)
				if ((buf[j] <= med) && (buf[j] > max))
					max = buf[j];
			med = max;
		}
	}

	buf[N2] = med;
	i = N2;

/* Main loop:  buf contains the N most recent inputs.  If input(n) > median
	and input(n-N) > median, then median is unchanged.  The same applies
	if both are less than the median.  Hence the median is only re-
	calculated when input(n) > median and input(n-N) < median or (vice
	versa).  The new median is simply the next largest (or smallest) of
	the N most recent inputs. */

	while (getfloat(&new) > 0){
		getfloat(&dummy);

		old = buf[i];
		buf[i] = new;
		if (++i >= N) i = 0;

		if ((new > med) && (old <= med)){
			min = 1e6;
			for (j = 0; j < N; j++)
				if ((buf[j] >= med) && (buf[j] < min))
					min = buf[j];
			med = min;
		}

		if ((new < med) && (old >= med)){
			max = -1e6;
			for (j = 0; j < N; j++)
				if ((buf[j] <= med) && (buf[j] > max))
					max = buf[j];
			med = max;
		}

		putfloat(&med);
		putfloat(&dummy);

	}

	flushfloat();

}

usage(exitcode)
	int exitcode;
{
	fprintf(stderr, "%s%s%s%s",
	"median [flag] < floatsams > output\n",
	" flags:\n",
	" -wN\t set window size to N seconds (default: 9S).\n",
	" -RN\t set sample rate to N (default = read from stdin).\n"
	);
	exit(exitcode);
}

malerr(str, ex)
	char *str;
	int ex;
{
	fprintf(stderr, "%s\n", str);
	exit(ex);
}

These are the contents of the former NiCE NeXT User Group NeXTSTEP/OpenStep software archive, currently hosted by Netfuture.ch.