This is srate_convert.c in view mode; [Download] [Up]
#include <stdio.h>
#include <sound/sound.h>
#include <math.h>
/* This is a command line c function call to implement sampling rate conversion on a sound.
The conversion is performed by sinc interpolation (ideal low-pass filter response).
One method uses sin() function calls from the floating point library, and the other method
uses linearly interpolated table lookup as described in:
J. O. Smith and P. Gossett, "A Flexible Sampling-Rate Conversion Method,"
Proc. of the IEEE Conference on Acoustics, Speech, and Signal Processing,
San Diego, CA, March, 1984.
The basic parameters select how wide the sinc function for interpolation is (how
many neighboring samples are used), whether to use a table, and how finely the sinc
is sampled in forming the table. Usage is:
srate_convert factor filein fileout
The factor which specifies the conversion is a floating point number expressing the ratio
of the old sampling rate to the new one. The sampling rate in the soundfile header is
left alone, so this operation causes pitch-shifting as implemented, and the factor is
the ratio of the processed pitch to the original. The factor can also be viewed as the
new sampling period, relative to the original sampling period.
The sound is assumed to be band limited to 1/2 the highest sampling rate, so no extra
filtering is done in the interpolation. Thus sampling rate conversion factors of less than
one are fine, and the signal must be pre-filtered to use sampling rate conversion factors
which are greater than one. To implement low-pass filtering as part of the conversion
operation, the low-pass filter coefficients can be integrated into the sinc table by
convolution. This is described in the Smith/Gossett article, and requires a generic
filter design library such as the IEEE filter design routines. In general, a low-pass
filter with cutoff below 1/2 the new sampling rate is designed, and convolved with the
sinc table. This yield the new sinc table to be used for interpolation.
The real 'inner loop', where the work gets done is enclosed in asterisk lines. If the
sound file format and number of channels is known, the unneeded code could be removed,
yielding about 80 lines for the entire sampling rate conversion routine, the function
which fills the sinc table, the linear interpolation function, and the sinc table access
function. The file codex_convert.c is an example of a pared-down example.
*/
#define WIDTH 5 /* this controls the number of neighboring samples
which are used to interpolate the new samples. The
processing time is linearly related to this width */
#define USE_TABLE TRUE /* this controls whether a linearly interpolated lookup
table is used for sinc function calculation, or the
sinc is calculated by floating point trig function calls. */
#define SAMPLES_PER_ZERO_CROSSING 5 /* this defines how finely the sinc function is sampled
for storage in the table */
float sinc_table[WIDTH * SAMPLES_PER_ZERO_CROSSING] = { 0.0 };
main(int argc,char *argv[]) /* sampling rate convert filename by factor.
factor should be less than 1.0 or file must
be pre-filtered to band-limit to NewSrate/2. */
{
extern double sinc(double x),make_sinc();
char *file_name_in,*file_name_out;
float fact;
double temp1=0.0,temp2=0.0,time=0.0,factor;
unsigned char *codex_data_in,*codex_data_out;
short *linear_data_in,*linear_data_out;
int i=0,j,left_limit,right_limit,num_samples,num_channels;
SNDSoundStruct *sound_in,*sound_out;
if (argc!=4) {
fprintf(stderr,"usage: srate_convert factor filein fileout\n");
exit(0);
}
file_name_in = argv[2];
file_name_out = argv[3];
sscanf(argv[1],"%f",&fact);
factor = fact;
if (factor<0.01 || factor>100) {
fprintf(stderr,"Your factor is a little extreme!!!\n");
exit(0);
}
printf("Converting by: %f file: %s to file: %s\n",factor,file_name_in,file_name_out);
SNDReadSoundfile(file_name_in,&sound_in);
num_channels = sound_in->channelCount;
num_samples = sound_in->dataSize / num_channels;
/* allocate an output sound and set data pointers */
if (sound_in->dataFormat==1) {
SNDAlloc(&sound_out, (int) ((double) num_samples / factor),
sound_in->dataFormat, sound_in->samplingRate, num_channels, 4);
codex_data_in = (unsigned char *) sound_in + sound_in->dataLocation;
codex_data_out = (unsigned char *) sound_out + sound_out->dataLocation;
}
else {
num_samples /= 2;
SNDAlloc(&sound_out, (int) (num_samples * 2 * num_channels / factor),
sound_in->dataFormat, sound_in->samplingRate, num_channels, 4);
linear_data_in = (short *) sound_in + sound_in->dataLocation;
linear_data_out = (short *) sound_out + sound_out->dataLocation;
}
if(USE_TABLE) make_sinc();
/* Everything up to here was setup
Entry at this point assumes an input sound band-limited to 1/2 the highest sampling rate,
an output sound allocated to the correct size (input sound size * conversion factor)
pointers to the sounds,
and a factor which is the new sampling period relative to the original */
/************************* INNER LOOP HERE ********************************************/
while (time < num_samples) {
temp1 = 0.0;
temp2 = 0.0;
left_limit = time - WIDTH + 1; /* leftmost neighboring sample used for interp.*/
right_limit = time + WIDTH; /* rightmost leftmost neighboring sample used for interp.*/
if (left_limit<0) left_limit = 0;
if (right_limit>num_samples) right_limit = num_samples;
/* Here comes the real work. True interpolation of sampled data
signals involves ideal low-pass filtering in the frequency
domain, which corresponds to convolution with a sinc function
in the time domain (consult any text on Fourier transforms). */
if (num_channels==1) {
if (sound_in->dataFormat==1) {
for (j=left_limit;j<right_limit;j++)
temp1 += SNDiMulaw(codex_data_in[j]) * sinc(time - j);
codex_data_out[i] = SNDMulaw((short) temp1);
}
else {
for (j=left_limit;j<right_limit;j++)
temp1 += linear_data_in[j] * sinc(time - j);
linear_data_out[i] = (short) temp1;
}
}
else {
for (j=left_limit;j<right_limit;j++) {
temp1 += linear_data_in[j * 2] * sinc(time - j);
temp2 += linear_data_in[j * 2 + 1] * sinc(time - j);
}
linear_data_out[i * 2] = (short) temp1;
linear_data_out[i * 2 + 1] = (short) temp2;
}
i++;
time += factor;
}
/***********************************************************************************************/
SNDWait(SNDWriteSoundfile(file_name_out,sound_out));
SNDFree(sound_in);
SNDFree(sound_out);
return;
}
#define PI 3.14159263
make_sinc()
{
int i;
double temp,win_freq,win;
win_freq = PI / WIDTH / SAMPLES_PER_ZERO_CROSSING;
sinc_table[0] = 1.0;
for (i=1;i<WIDTH * SAMPLES_PER_ZERO_CROSSING;i++) {
temp = (double) i * PI / SAMPLES_PER_ZERO_CROSSING;
sinc_table[i] = sin(temp) / temp;
win = 0.5 + 0.5 * cos(win_freq * i);
sinc_table[i] *= win;
}
return;
}
double linear_interp(double first_number,double second_number,double fraction)
{
return (first_number + ((second_number - first_number) * fraction));
}
double t_sinc(double x)
{
int low;
double temp,delta;
if (fabs(x)>WIDTH-1)
return 0.0;
else {
temp = fabs(x) * (double) SAMPLES_PER_ZERO_CROSSING;
low = temp; /* these are interpolation steps */
delta = temp - low; /* and can be ommited if desired */
return linear_interp(sinc_table[low],sinc_table[low + 1],delta);
/* for no interpolation case, return sinc_table[low] */
}
}
double sinc(double x)
{
double temp;
if(USE_TABLE) return t_sinc(x);
else {
if (x==0.0) return 1.0;
else {
temp = PI * x;
return sin(temp) / (temp);
}
}
}
These are the contents of the former NiCE NeXT User Group NeXTSTEP/OpenStep software archive, currently hosted by Netfuture.ch.