X-Git-Url: https://git.sesse.net/?p=c64tapwav;a=blobdiff_plain;f=synth.cpp;h=90c16e05d956d055e2fead4f1d01b78d1626a470;hp=f40a76b4837d61158582defe3c965f00305ddfed;hb=866751f2909657ee7d62f936e3d2cb2ecb7806fb;hpb=aec567fb778f2a3dc479dc7bd47551b96b361ff0;ds=sidebyside

diff --git a/synth.cpp b/synth.cpp
index f40a76b..90c16e0 100644
--- a/synth.cpp
+++ b/synth.cpp
@@ -4,17 +4,20 @@
 #include <vector>
 #include <assert.h>
 
+#include "common.h"
+
 using namespace std;
 
 #define LANCZOS_RADIUS 10
 #define RESOLUTION 256
-#define WAVE_FREQ 44100
-#define C64_FREQ 985248
 
 // Cutoff frequency of low-pass filter (must be max. WAVE_FREQ / 2 or you
 // will get aliasing)
 #define LPFILTER_FREQ 22050
 
+// Cutoff frequency of high-pass filter
+#define HPFILTER_FREQ 40
+
 #define NORMALIZED_LPFILTER_FREQ (float(LPFILTER_FREQ) / float(WAVE_FREQ))
 #define LANCZOS_EFFECTIVE_RADIUS (LANCZOS_RADIUS / (NORMALIZED_LPFILTER_FREQ * 2.0))
 
@@ -80,45 +83,62 @@ struct pulse {
 };
 vector<pulse> pulses;
 
-int main(int argc, char **argv)
+static float a1, a2, b0, b1, b2;
+static float d0, d1;
+
+static void filter_init(float cutoff_radians)
 {
-	make_table();
+	float resonance = 1.0f / sqrt(2.0f);
+	float sn, cs;
+	sincosf(cutoff_radians, &sn, &cs);
+
+	float alpha = float(sn / (2 * resonance));
+
+	// coefficients for highpass filter
+        float a0 = 1 + alpha;
+	b0 = (1 + cs) * 0.5f;
+	b1 = -(1 + cs);
+	b2 = b0;
+        a1 = -2 * cs;
+        a2 = 1 - alpha;
+
+	float invA0 = 1.0f / a0;
+	b0 *= invA0;
+	b1 *= invA0;
+	b2 *= invA0;
+	a1 *= invA0;
+	a2 *= invA0;
+
+	// reset filter delays
+	d0 = d1 = 0.0f;
+}
 
-	FILE *fp = fopen(argv[1], "rb");
-	fseek(fp, 14, SEEK_SET);
-
-	int x = 0;
-
-	while (!feof(fp)) {
-		int len = getc(fp);
-		int cycles;
-		if (len == 0) {
-			int a = getc(fp);
-			int b = getc(fp);
-			int c = getc(fp);
-			cycles = a | (b << 8) | (c << 16);
-		} else {
-			cycles = len * 8;
-		}
-		pulse p;
-		p.start = float(x) * WAVE_FREQ / C64_FREQ;
-		p.end = (float(x) + cycles * 0.5) * WAVE_FREQ / C64_FREQ;
-		pulses.push_back(p);
-		x += cycles;
-	}
+static float filter_update(float in)
+{
+	float out = b0*in + d0;
+	d0 = b1 * in - a1 * out + d1;
+	d1 = b2 * in - a2 * out;
+	return out;
+}
 
-	int len_cycles = x;
-	int len_samples = int(ceil(float(len_cycles) * WAVE_FREQ / C64_FREQ));
+vector<float> synth(const vector<pulse> &pulses)
+{
+	make_table();
+
+	int len_samples = int(ceil(pulses.back().start + (pulses.back().end - pulses.back().start) * 2));
 	
-	fprintf(stderr, "%d pulses, total %.2f seconds (%d samples)\n", pulses.size(), len_cycles / float(C64_FREQ), len_samples);
+	fprintf(stderr, "%d pulses, total %.2f seconds (%d samples)\n", int(pulses.size()), len_samples / float(WAVE_FREQ), len_samples);
 
 	int pulse_begin = 0;
 
+	vector<float> samples;
+	samples.reserve(len_samples);
+
 	for (int i = 0; i < len_samples; ++i) {
 //	for (int i = 50000000; i < 51000000; ++i) {
 //		double t = i * 0.01;
 		double t = i;
-		double sample = 0.0;
+		double sample = -0.5;
 		for (unsigned j = pulse_begin; j < pulses.size(); ++j) {
 			if (t - pulses[j].end > LANCZOS_EFFECTIVE_RADIUS) {
 				++pulse_begin;
@@ -130,8 +150,22 @@ int main(int argc, char **argv)
 			float contribution = window_integral(pulses[j].start - t, pulses[j].end - t);
 			sample += contribution;
 		}
-		printf("%f\n", sample);
-//		short s = lrintf((sample - 0.5f) * 16384.0f);
-//		fwrite(&s, 2, 1, stdout);
+		samples.push_back(sample);
 	}
+
+	// filter forwards, then backwards (perfect phase filtering)
+	vector<float> filtered_samples, refiltered_samples;
+	filtered_samples.resize(len_samples);
+	refiltered_samples.resize(len_samples);
+
+	filter_init(M_PI * HPFILTER_FREQ / WAVE_FREQ);
+	for (int i = 0; i < len_samples; ++i) {
+		filtered_samples[i] = filter_update(samples[i]);
+	}
+	filter_init(M_PI * HPFILTER_FREQ / WAVE_FREQ);
+	for (int i = len_samples; i --> 0; ) {
+		refiltered_samples[i] = filter_update(filtered_samples[i]);
+	}
+
+	return refiltered_samples;
 }