Patch fixes bug in the original linear interpolation,

along with generating a proper "streaming" resampled output.

src/spicelib/devices/vsrc/vsjack.c

@@ -7,7 +7,10 @@
 #include <math.h>
 #include <sndfile.h>
 #include <inttypes.h>
-#define VS_BUFSIZ 1024
+
+// Resampling can be rather slow. Don't resample
+// the whole audio file, do it in smaller chunks
+#define VS_RESAMPLING_CHUNK 1024
 
 #include "ngspice/ngspice.h"
 
@@ -21,18 +24,18 @@ static int m_channels[MAX_D]; //< number of channles in src-file
 static uint32_t m_samplerate[MAX_D]; //< samplerate of source
 static uint32_t m_frames[MAX_D]; //< duration of source in frames
 static float* (interleaved[MAX_D]); //< internal soundfile buffer
-static uint32_t ilb_start[MAX_D]; //< first sample in buffer 
-static uint32_t ilb_end[MAX_D]; //< last sample in buffer
 
 #define HAVE_SRC
 
 #ifdef HAVE_SRC
 #include <samplerate.h>
 static double src_ratio = 64.0;
-#define SRC_RATIO 64
+#define SRC_RATIO (src_ratio)
 static SRC_STATE* rabbit[MAX_D];
 static int rabbit_err[MAX_D];
 static float* (resampled[MAX_D]); //< internal soundfile buffer
+static uint32_t input_frames_used[MAX_D];
+static uint32_t output_frames_generated[MAX_D];
 #endif
 
 void vsjack_initialize(void) {
@@ -73,6 +76,7 @@ void closefile_sf(int d) {
 #endif
 
 int openfile_sf(int d, char* filename) {
+    uint32_t nframes;
     SF_INFO sfinfo;
     if (!m_sndfile[d])
         sf_close(m_sndfile[d]);
@@ -88,7 +92,6 @@ int openfile_sf(int d, char* filename) {
 
     m_sndfile[d] = sf_open(path, SFM_READ, &sfinfo);
     txfree(path);
-    ilb_end[d] = ilb_start[d] = 0;
 
     if (SF_ERR_NO_ERROR != sf_error(m_sndfile[d])) {
         fprintf(stderr, "Error: This is not a sndfile supported audio file format\n");
@@ -98,100 +101,96 @@ int openfile_sf(int d, char* filename) {
         fprintf(stderr, "Error: This is an empty audio file\n");
         return (-1);
     }
+    nframes = sfinfo.frames;
+
     m_channels[d] = sfinfo.channels;
     m_samplerate[d] = sfinfo.samplerate;
-    m_frames[d] = (uint32_t)sfinfo.frames;
+    m_frames[d] = nframes;
-    realloc_sf(d, VS_BUFSIZ);
+    realloc_sf(d, nframes);
 #ifdef HAVE_SRC
-    realloc_src(d, VS_BUFSIZ * SRC_RATIO);
+    SRC_DATA src_data;
+
+    realloc_src(d, nframes * SRC_RATIO);
     rabbit[d] = src_new(SRC_SINC_BEST_QUALITY, m_channels[d], &(rabbit_err[d]));
     src_set_ratio(rabbit[d], SRC_RATIO);
     src_reset(rabbit[d]);
+
 #endif
+    nframes = sf_readf_float(m_sndfile[d], (interleaved[d]), nframes);
+    if (nframes < 0) {
+	fprintf(stderr, "Error: Failed to read audio frames\n");
+	return (-1);
+    }
+    m_frames[d] = nframes;
     return (0);
 }
 
-void load_buffer(int d, uint32_t sample) {
-    sf_seek(m_sndfile[d], sample, SEEK_SET);
-    ilb_start[d] = sample;
-    uint32_t nframes;
-    if ((nframes = (uint32_t)sf_readf_float(m_sndfile[d], (interleaved[d]), VS_BUFSIZ)) > 0) {
-        ilb_end[d] = ilb_start[d] + nframes;
-    }
-    else {
-        ilb_end[d] = ilb_start[d];
-        printf("Decoder error.\n");
-    }
-#ifdef HAVE_SRC
-    SRC_DATA src_data;
-    src_data.data_in = interleaved[d];
-    src_data.data_out = resampled[d];
-    src_data.input_frames = VS_BUFSIZ;
-    src_data.output_frames = VS_BUFSIZ * SRC_RATIO;
-    src_data.end_of_input = ((ilb_end[d] - ilb_start[d]) < VS_BUFSIZ);
-    src_data.src_ratio = SRC_RATIO;
-    src_data.input_frames_used = 0;
-    src_data.output_frames_gen = 0;
-
-    src_process(rabbit[d], &src_data);
-#endif
-}
-
 double get_value(int d, double time, int channel) {
-    uint32_t sample = (uint32_t)floor(time * ((double)m_samplerate[d]));
+    uint32_t nframes = m_frames[d];
+    double sample_fp = time * ((double)m_samplerate[d]);
+    uint32_t sample = (uint32_t)floor(sample_fp);
 
-    // TODO: print EOF warning (once). FIXME move to load_buffer
+    if (sample >= nframes) return (0.0);
-    if (sample > m_frames[d]) return (0.0);
-
-    if (sample < ilb_start[d] || sample >= ilb_end[d])
-        load_buffer(d, sample);
-
-    if (sample < ilb_start[d] || sample >= ilb_end[d]) {
-        printf("no such value buffered for file:%i.\n", d);
-        return (0.0); // nan ?
-    }
 
 #ifdef HAVE_SRC
-    int offset = (int)floor((sample - ilb_start[d]) * SRC_RATIO);
+    sample_fp *= SRC_RATIO;
-    if (offset > VS_BUFSIZ * SRC_RATIO || offset < 0) {
+    sample = (uint32_t)floor(sample_fp);
-        printf("value not in buffer:%i.\n", d);
+
-        return (0.0); // nan ?
+    // Do we need to generate more output frames?
+    while (sample >= output_frames_generated[d]) {
+	SRC_DATA src_data;
+	uint32_t output_generated = output_frames_generated[d];
+	uint32_t input_used = input_frames_used[d];
+	uint32_t input_frames_left = nframes - input_used;
+
+	// Not enough output frames, and nothing more to input?
+	// Give up.
+	if (!input_frames_left)
+	    return (0.0);
+
+	// Do the resampling in smaller chunks
+	src_data.end_of_input = 1;
+	if (input_frames_left > VS_RESAMPLING_CHUNK) {
+	    input_frames_left = VS_RESAMPLING_CHUNK;
+	    src_data.end_of_input = 0;
 	}
-    float val = ((float*)(resampled[d]))[m_channels[d] * offset + channel];
+
-#   if 0 // DEBUG 
+	src_data.data_in = interleaved[d] + m_channels[d] * input_used;
-#   define SQUARE(A) ((A)*(A))
+	src_data.data_out = resampled[d] + m_channels[d] * output_generated;
-    static double stride = 0;
+	src_data.input_frames = input_frames_left;
-    static double last = 0;
+	src_data.output_frames = nframes * SRC_RATIO - output_generated;
-    static double deviation = 0;
+	src_data.src_ratio = SRC_RATIO;
-    static int dev_cnt = 0;
+	src_data.output_frames_gen = 0;
-    if (channel == 0) {
+	src_data.input_frames_used = 0;
-        stride += (SRC_RATIO * time * ((double)m_samplerate[d])) - last;
+
-        last = (SRC_RATIO * time * ((double)m_samplerate[d]));
+	if (src_process(rabbit[d], &src_data)) {
-        deviation += SQUARE((SRC_RATIO * time * ((double)m_samplerate[d])) - floor(SRC_RATIO * time * ((double)m_samplerate[d])));
+	    fprintf(stderr, "src_process() failed on sound file");
-        dev_cnt++;
+	    return -1;
-        if ((dev_cnt % (12000)) == 0)
-            printf("read time dev= %f - stride= %f\n", sqrt(deviation / (double)dev_cnt), stride / (double)dev_cnt);
 	}
-#   endif 
+
-# if 0 // zero order hold.
+	output_frames_generated[d] += src_data.output_frames_gen;
-    return((double)val);
+	input_frames_used[d] += src_data.input_frames_used;
-# else
+	if (src_data.end_of_input)
-    // linear interpolation
+		break;
-    float val1 = ((float*)(resampled[d]))[(m_channels[d] * (offset + 1)) + channel];
+    }
-    double diff = (SRC_RATIO * time * ((double)m_samplerate[d])) -
+
-        floor(SRC_RATIO * time * ((double)m_samplerate[d]));
+    // Are we past all the generated samples?
+    if (sample >= output_frames_generated[d])
+	return (0.0);
+    float val = ((float*)(resampled[d]))[m_channels[d] * sample + channel];
+    // Are we the last sample?
+    if (sample + 1 == output_frames_generated[d])
+	return val;
+
+    // linear interpolation between samples
+    double diff = sample_fp - sample;
+    float val1 = ((float*)(resampled[d]))[(m_channels[d] * (sample + 1)) + channel];
     double rv = ((double)val) * (1.0 - diff) + ((double)val1) * diff;
     return(rv);
-# endif
 
 #else // no upsampling.
 
-    int offset = sample - ilb_start[d];
+    return((double)(((float*)(interleaved[d]))[m_channels[d] * sample + channel]));
-    if (offset > VS_BUFSIZ || offset < 0) {
-        printf("value not in buffer:%i.\n", d);
-        return (0.0); // nan ?
-    }
-    return((double)(((float*)(interleaved[d]))[m_channels[d] * offset + channel]));
 #endif
     }