@ James, lidnariq, tepples: Thanks for all your help. I think I finally got it.
Using GNU Octave, I generated
this 13th-order Elliptic filter consisting of 7 biquad segments, with the following code:
Code:
clear
fc = 20000; % Cut-off frequency (Hz)
fs = 19687500.0 / 11.0; % NTSC Sampling rate (Hz)
%fs = 53203425.0 / 32.0; % PAL Sampling rate (Hz)
order = 13; % Filter order
[z, p, k] = ellip(order, 0.1, 100, 2*fc/fs);
[sos, g] = zp2sos(z, p, k);
Bs = sos(:,1:3)
num2hex(Bs)
As = sos(:,4:6)
num2hex(As)
[nsec,temp]=size(sos);
nsamps = 100000;
x = g*[1,zeros(1,nsamps-1)];
for i=1:nsec
x = filter(Bs(i,:),As(i,:),x);
end
X = fft(x);
f=[0:nsamps-1]*fs/nsamps;
figure(1);
grid('on');
axis([0 fc*2 -100 5]);
legend('off');
pscale = 48
plot(f(1:nsamps/pscale),20*log10(X(1:nsamps/pscale)));
The cutoff frequency is set to 20000 Hz and it attenuates down to -100 dB by 22000 Hz, satisfying the Nyquist criterion for the traditional audio sampling rates of 44100+ Hz.
lidnariq wrote:
James was using a sample rate of 1786830 intentionally instead of the NTSC 1789772⁸/₁₁: the NTSC NES produces video at 60.1 Hz, not 60 Hz as the end-user's monitor does. If one wants to avoid visible tearing, you need to slow down the emulated NES by 1646ppm:
39375000÷11 Hz × 1½ ÷ (341×262-0.5) = 60.0988…
39375000÷11 Hz ÷ 2 × 60 ÷ 60.0988… = 1786830 Hz exactly
In my emulator, the CPU, PPU and APU use an independent timer to run at the correct rates. The monitor is updated at its refresh rate on a separate thread and they coordinate using the thread disruptor pattern. This is a topic that can be discussed in a separate thread if anyone is interesting in further details.
Below is a version in Java that uses 64-bit doubles. I tested out 32-bit floats and it remains stable and mostly accurate. But, I found that the 64-bit version does not consume a lot of CPU on my box; so, I'll stick with the extra accuracy for now.
Like the Butterworth version that I posted earlier, this one provides an addInputSample() method that is invoked once per CPU cycle and the listener will be called back at the output sampling frequency supplied to the constructor.
Code:
public class Decimator {
private static final double I2 = 1.0 / 2.0;
private static final double I3 = 1.0 / 3.0;
private static final double I6 = 1.0 / 6.0;
private static final double NTSC_SAMPLING_FREQUENCY = 19687500.0 / 11.0;
private static final double PAL_SAMPLING_FREQUENCY = 53203425.0 / 32.0;
private static final long[] LONG_NTSC_A1 = {
0xbfff67701c84a3deL,
0xbfff8b62242e56e4L,
0xbfffaf25694786a4L,
0xbfffc959eca889e2L,
0xbfffda7652351a06L,
0xbfffe6495d699ef1L,
0xbfef576a39439bdaL,
};
private static final long[] LONG_NTSC_A2 = {
0x3feed6d89d961d3bL,
0x3fef28cd096bf76aL,
0x3fef7a3a53cd4f67L,
0x3fefb5a2213e3abbL,
0x3fefdbe39ed2f06eL,
0x3feff554ab2f0006L,
};
private static final long[] LONG_NTSC_B1 = {
0xbfff2a2df4fd2ad2L,
0xbfffbe3738d86ac8L,
0xbfffd94c2258d55cL,
0xbfffe22ba0cf6805L,
0xbfffe5b14c11db72L,
0xbfffe6ffa75a14ebL,
};
private static final long LONG_NTSC_G = 0x3ecc208a2d079678L;
private static final double[] NTSC_A1;
private static final double[] NTSC_A2;
private static final double[] NTSC_B1;
private static final double NTSC_G;
static {
NTSC_A1 = toDoubleArray(LONG_NTSC_A1);
NTSC_A2 = toDoubleArray(LONG_NTSC_A2);
NTSC_B1 = toDoubleArray(LONG_NTSC_B1);
NTSC_G = Double.longBitsToDouble(LONG_NTSC_G);
}
private static final long[] LONG_PAL_A1 = {
0xbfff5baa1b81a309L,
0xbfff81d9b2eb15feL,
0xbfffa7df44761c9dL,
0xbfffc3c2600f5556L,
0xbfffd5ff6fbab417L,
0xbfffe2a551b417d0L,
0xbfef4aa6ceb6c246L,
};
private static final long[] LONG_PAL_A2 = {
0x3feec08d9ac7eb51L,
0x3fef18943f900329L,
0x3fef7018a2a55ed4L,
0x3fefaffb4199e44aL,
0x3fefd9237788b71bL,
0x3feff4844bf56456L,
};
private static final long[] LONG_PAL_B1 = {
0xbfff08b370e2ff76L,
0xbfffb3ce724731b7L,
0xbfffd3295d0112b5L,
0xbfffdd7036a20286L,
0xbfffe184a77ef413L,
0xbfffe307f8bb95fdL,
};
private static final long LONG_PAL_G = 0x3ece3d98b822795aL;
private static final double[] PAL_A1;
private static final double[] PAL_A2;
private static final double[] PAL_B1;
private static final double PAL_G;
static {
PAL_A1 = toDoubleArray(LONG_PAL_A1);
PAL_A2 = toDoubleArray(LONG_PAL_A2);
PAL_B1 = toDoubleArray(LONG_PAL_B1);
PAL_G = Double.longBitsToDouble(LONG_PAL_G);
}
private final double g;
private final double[] a1;
private final double[] a2;
private final double[] b1;
private final double[] d1 = new double[7];
private final double[] d2 = new double[6];
private final double[] ys = new double[4];
private final double inputSamplingFrequency;
private final double inputSamplingPeriod;
private final double outputSamplingPeriod;
private final DecimatorListener listener;
private double time;
private int index;
public Decimator(boolean ntsc, double outputSamplingFrequency,
DecimatorListener listener) {
if (ntsc) {
a1 = NTSC_A1;
a2 = NTSC_A2;
b1 = NTSC_B1;
g = NTSC_G;
inputSamplingFrequency = NTSC_SAMPLING_FREQUENCY;
} else {
a1 = PAL_A1;
a2 = PAL_A2;
b1 = PAL_B1;
g = PAL_G;
inputSamplingFrequency = PAL_SAMPLING_FREQUENCY;
}
this.listener = listener;
inputSamplingPeriod = 1.0 / inputSamplingFrequency;
outputSamplingPeriod = 1.0 / outputSamplingFrequency;
}
private static double[] toDoubleArray(long[] values) {
double[] ds = new double[values.length];
for(int i = values.length - 1; i >= 0; i--) {
ds[i] = Double.longBitsToDouble(values[i]);
}
return ds;
}
public void addInputSample(double x) {
double y;
for(int i = 0; i < 6; i++) {
y = x + d1[i];
d1[i] = b1[i] * x - a1[i] * y + d2[i];
d2[i] = x - a2[i] * y;
x = y;
}
y = x + d1[6];
d1[6] = x - a1[6] * y;
ys[index] = g * y;
time += inputSamplingPeriod;
if (time >= outputSamplingPeriod) {
listener.processOutputSample(computeSpline(
1.0 - (time - outputSamplingPeriod) * inputSamplingFrequency,
ys[(index - 3) & 3],
ys[(index - 2) & 3],
ys[(index - 1) & 3],
ys[index]));
time -= outputSamplingPeriod;
}
index = (index + 1) & 3;
}
private double computeSpline(double t, double y0, double y1, double y2,
double y3) {
double c0 = y1;
double c1 = y2 - I3 * y0 - I2 * y1 - I6 * y3;
double c2 = I2 * (y0 + y2) - y1;
double c3 = I6 * (y3 - y0) + I2 * (y1 - y2);
return ((c3 * t + c2) * t + c1) * t + c0;
}
}
public interface DecimatorListener {
void processOutputSample(double sample);
}