Created
January 20, 2015 22:17
-
-
Save teddokano/17c53e6da8c4b5cb8e70 to your computer and use it in GitHub Desktop.
Noise shaper (1st and 2nd order) simulation as writing practice of GNU Octave
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
% 1st and 2nd order noise shaper time/frequency domain plot by GNU Octave | |
clear -all | |
% | |
% parameters | |
% | |
calc_length = 2**16; % number of samples | |
t_domain_start = 280; | |
t_domain_stop = 850; | |
signal_freq = 10000; % Hz | |
sampling_freq = 44100; % Hz | |
oversampling_ratio = 256; % Fs | |
attenuation = 1; % dB to full-scale | |
signal_duration = calc_length / sampling_freq; | |
oversampling_freq = sampling_freq * oversampling_ratio; | |
amplitude = 10 ** (-attenuation/20); | |
% | |
% input signal | |
% | |
in = sinetone( signal_freq, oversampling_freq, signal_duration, amplitude ); | |
% | |
% 1st order noise shaper | |
% | |
noise = 0; | |
for i = 1 : calc_length | |
qin = in(i) - noise; | |
if ( qin < 0 ) | |
out1(i) = -1; | |
else | |
out1(i) = 1; | |
endif | |
noise = out1(i) - qin; | |
end | |
% | |
% 2nd order noise shaper (SAA7320) | |
% | |
s2 = 0; | |
s1 = 0; | |
for i = 1 : calc_length | |
if ( s2 < 0 ) | |
out2(i) = -1; | |
else | |
out2(i) = 1; | |
endif | |
noise = out2(i) - s2; | |
s2 = s1 - noise * 2; | |
s1 = in(i) + noise; | |
end | |
% | |
% window function (blackman-harris) | |
% | |
N = calc_length-1; | |
a0 = 0.35875; | |
a1 = 0.48829; | |
a2 = 0.14128; | |
a3 = 0.01168; | |
for i = 0 : N | |
wndw(i+1) =0.35875-a1*cos(2*pi*i/N)+a2*cos(4*pi*i/N)-a3*cos(6*pi*i/N); | |
end | |
% | |
% applying window function to 1 bit output data | |
% | |
w_mean = mean( wndw ); | |
wndw_in1 = (out1.') .* (wndw.'); | |
wndw_in2 = (out2.') .* (wndw.'); | |
% | |
% get power spectrum by FFT | |
% | |
spectrum1 = 20 * log10( abs( fft( wndw_in1,calc_length ) ) ./ ((calc_length / 2) * w_mean) ); | |
spectrum2 = 20 * log10( abs( fft( wndw_in2,calc_length ) ) ./ ((calc_length / 2) * w_mean) ); | |
spectrum1 = spectrum1( 1:calc_length/2 ); | |
spectrum2 = spectrum2( 1:calc_length/2 ); | |
% | |
% plot timedomain (waveform) | |
% | |
subplot( 2, 1, 1 ); | |
t = t_domain_start:1:t_domain_stop; | |
ph_t0 = plot( t, in(t_domain_start:t_domain_stop) ); | |
hold on | |
ph_t1 = plot( t, out1(t_domain_start:t_domain_stop) ); | |
ph_t2 = plot( t, out2(t_domain_start:t_domain_stop) ); | |
hold off | |
set( ph_t1, "color", "blue" ); | |
set( ph_t2, "color", "red" ); | |
legend( ph_t1, "1st order" ); | |
legend( ph_t2, "2nd order" ); | |
axis( [ t_domain_start, t_domain_stop, -1.1, 1.1 ] ); | |
xlabel( "time [sample]" ); | |
ylabel( "amplitude [fullscale = -1..+1]" ); | |
% | |
% plot frequency domain (spectrum) | |
% | |
subplot( 2, 1, 2 ); | |
f = 0 : oversampling_freq / (calc_length - 1) : oversampling_freq / 2; | |
ph_f1 = semilogx( (f.'), spectrum1 ); | |
hold on; | |
ph_f2 = semilogx( (f.'), spectrum2 ); | |
hold off | |
set( ph_f1, "color", "blue" ); | |
set( ph_f2, "color", "red" ); | |
legend( ph_f1, "1st order" ); | |
legend( ph_f2, "2nd order" ); | |
axis([ 100, oversampling_freq / 2, -150, 0 ] ); | |
xlabel( "frequency [Hz]" ); | |
ylabel( "power [dBFS]" ); | |
grid on; |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Info about this code (Japanese): http://goo.gl/Tb7Pud