Author
| Crystal clean PCM 8bit samples on the poor PSG
| POISONIC
msx professional
Posts: 883 |  Posted: January 20 2006, 17:10 | WOW but i guess its very CPU time consuming but the quality overall is fantastic especialy that royksopp sample 
| | [WYZ]
msx lover
Posts: 94 | Posted: January 20 2006, 17:22 | The best sample I ever heard on MSX before. Very hight quality and noiseless. Congratulations Artrag and Dvik.
| | dvik
msx master
Posts: 1302 | Posted: January 20 2006, 17:28 | Quote:
| Note that the encoder program doesn't "just" compile on my LInux box;
st.h:25: error: 'uint32_t' does not name a type
for example. A bit strange, in st.h the types are defined only for non-GNU compilers??
If I reverse that ifdef, it compiles fine.
|
Maybe change line 20 in the same file from
#ifndef __GNUC__ )
to
# if 1
solves the problem.
I'm not sure why these types are defined when compiling with cygwin gcc.
| | ARTRAG
msx master
Posts: 1578 | Posted: January 20 2006, 17:35 | but it works!!!
| | ARTRAG
msx master
Posts: 1578 | Posted: January 20 2006, 18:03 | About CPU time, I think that there is spare time for many things,
maybe also for an encoding better than RLE.
| | dvik
msx master
Posts: 1302 | Posted: January 20 2006, 18:14 | Yes, the 8kHz replayer is quite similar to the one I used in Waves, so its definately possible to do something else in parallel. For demos I think that the replayer_core8000v2 is the best one since it leaves most unused CPU cycles to other things.
And actually when looking at it again that replayer plays 16kHz samples without any modifications (just add '-rto 2' to the command line options to the encoder). So you'll get a 16kHz replayer that leaves about 50% of the CPU for other tasks. Not bad huh?
| | ARTRAG
msx master
Posts: 1578 | Posted: January 24 2006, 00:45 | Actually it is!! the PSG can go even better if the replayer can change tha channel at each step!!
This viterbi core allows the change of any channel at each transition and gets 39,9-40dB on almost any input!!!!
u=1;
S=(zeros(16*16*16,3));
for i=0:15
for j=0:15
for k=0:15
S(u,:)=[ i j k ];
u=u+1;
end
end
end
for u=1:4096
S(u,:)=sort(S(u,:)')';
end
I=S(:,1)*256+S(:,2)*16+S(:,3);
[l i]=sort(I);
S=S(i,:);
Sv(1,:)=[0 0 0];
v=1;
for u=1:4095
if not(all(S(u,:)==S(u+1,:)))
Sv(v,:)=S(u,:);
v=v+1;
end
end
Sv(v,:)=S(4096,:);
%[Sv sum(vol(Sv+1)')']
I = Sv(:,1)*256+Sv(:,2)*16+Sv(:,3);
Ns=v;
nxtS = zeros(Ns,3*16);
for i=0:Ns-1
for in=0:15
ns=sort([in Sv(i+1,2) Sv(i+1,3)]')';
ins = ns(1)*256+ns(2)*16+ns(3);
nxtS(i+1,in+1)=find(ins==I)-1;
end
for in=0:15
ns=sort([Sv(i+1,1) in Sv(i+1,3)]')';
ins = ns(1)*256+ns(2)*16+ns(3);
nxtS(i+1,in+16+1)=find(ins==I)-1;
end
for in=0:15
ns=sort([Sv(i+1,1) Sv(i+1,2) in]')';
ins = ns(1)*256+ns(2)*16+ns(3);
nxtS(i+1,in+32+1)=find(ins==I)-1;
end
end
n=[0:15]';
vol=2.^(n/2)/2^7.5;
vol(1)=0;
curV = zeros(Ns,3*16);
for i=0:Ns-1
for in=0:15
curV(i+1,in+1) = vol(in+1)+vol(bitand(fix(i/16),15)+1)+vol(bitand(i,15)+1);
end
for in=0:15
curV(i+1,in+16+1) = vol(fix(i/256)+1)+vol(in+1)+vol(bitand(i,15)+1);
end
for in=0:15
curV(i+1,in+32+1) = vol(fix(i/256)+1)+vol(bitand(fix(i/16),15)+1)+vol(in+1);
end
end
f = 3579545;
inputfile = 'a1.wav';
inputfile = lower(inputfile);
if isempty(findstr(inputfile,'.wav'))
inputfile = strcat(inputfile,'.wav');
end
[y,Fs,nbit,opt]=wavread(inputfile);
y = (y-min(y)) /(max(y)-min(y)) * 1.15; % input in 0 e 1.15
f = 3579545;
fpcm = Fs;
tp = 447; % 3 PSG transitions per sample: put here the tp of your replayer
dt1 = (12+5+14)/tp;
dt2 = (12+14)/tp;
dt3 = 1-dt1-dt2;
dt = [dt1 dt2 dt3];
N = 3*length(y);
x = zeros(1,N);
x(1:3:end) = y;
x(2:3:end-3) = y(1:end-1)*dt1+y(2:end)*(1-dt1);
x(3:3:end-3) = y(1:end-1)*(dt1+dt2)+y(2:end)*(1-dt1-dt2);
x(end-2:end)= y(end);
disp(sprintf(' This corresponds to an input with samplig frequency = %3.1f Hz',f/tp))
Stt = uint16(zeros(Ns,N));
Itt = uint8(zeros(Ns,N));
L = zeros(Ns,1);
St = uint16(ones(Ns,1));
It = uint8(ones(Ns,1));
for t =1:N
Ln = inf*ones(Ns,1);
if mod(t-1,1000)==0
disp(sprintf(' Processing complete at %3.2f %%',(t/N*100)));
end
for cs=0:Ns-1
for in=0:47
cv = curV(cs+1,in+1);
ns = double(nxtS(cs+1,in+1));
Ltst = L(cs+1)+dt(mod(t-1,3)+1)*abs(x(t)-cv)^2;
if Ln(ns+1) >= Ltst
Ln(ns+1) = Ltst;
St(ns+1) = cs;
It(ns+1) = in;
end
end
end
L = Ln;
Stt(:,t) = St;
Itt(:,t) = It;
end
[l,i] = min(L);
P = uint16(zeros(1,N));
I = uint8(zeros(1,N));
P(N) = Stt(i,N);
I(N) = Itt(i,N);
for t = (N-1):-1:1
P(t) = Stt(double(P(t+1))+1,t);
I(t) = Itt(double(P(t+1))+1,t);
end
V = zeros(1,N);
for t = 1:N
V(t) = curV(double(P(t))+1,double(I(t))+1);
end
%er = l;
en = sum(abs(x(1:3:end)).^2)*dt1+...
sum(abs(x(2:3:end)).^2)*dt2+...
sum(abs(x(3:3:end)).^2)*dt3+...
-mean(x)^2;
er = sum(abs(x(1:3:end)-V(1:3:end)).^2)*dt1+...
sum(abs(x(2:3:end)-V(2:3:end)).^2)*dt2+...
sum(abs(x(3:3:end)-V(3:3:end)).^2)*dt3;
SNR = 10*log10(en/er);
disp(sprintf(' SNR = %2.1f dB',SNR))
plot(x,'r')
hold on
plot(V)
% meaning of I(t)
% If I(t) in [0 15] change to I(t) the volume in the triplet with minimum amplitude
% If I(t) in [16 31] change to I(t)-16 the volume in the triplet with middle amplitude
% I(t) in [32 47] change to I(t)-32 the volume in the triplet with max amplitude
| | ARTRAG
msx master
Posts: 1578 | Posted: January 24 2006, 08:32 | Actually I have found that the I(t) vector has to be further processed In order to get usable info about the channels.
The correct meaning of I(t) is
If I(t)<16 I have changed in I(t) the valume in the triplet with minimum amplitude
If 15<I(t)<32 I have changed in I(t)-16 the valume in the triplet with middle amplitude
If 31<I(t)<48 I have changed in I(t)-32 the valume in the triplet with max amplitude
I think that the easyest think to do to get a true channl number to be used in the ASM replayer is to run the whole sequence of inputs and build true triplets from the
ordered ones. This allows me to reconstruct a actual channel nmber to be stored in output. I will add this paer asap AR
PS
Note that the remarks in the file about I are wrong!!
| | ARTRAG
msx master
Posts: 1578 | Posted: January 24 2006, 16:11 | This "optimized" version produces a sequence ( Out ) that
encode in b7b6 which PSG channel should cange (0=A, 1=B, 2=C)
and the level in b3b2b1b0.
Actually the memory usage is huge, so the final version probably
will need some sub-optimal serch strategy (like truncated viterbi).
In any case the test wav performs very well, it gets 39,6dB on a1.wav
u = 1;
S = zeros(16*16*16,3);
for i=0:15
for j=0:15
for k=0:15
S(u,:)=[ i j k ];
u=u+1;
end
end
end
S = sort(S')';
I = S(:,1)*256+S(:,2)*16+S(:,3);
[l i] = sort(I);
S = S(i,:);
I = S(:,1)*256+S(:,2)*16+S(:,3);
i = [ find(diff(I)) ; 4096];
S = S(i,:);
I = I(i);
Ns = length(I);
nxtS = zeros(Ns,3*16);
for i=0:Ns-1
for in=0:15
ns = sort([in S(i+1,2) S(i+1,3)]')';
ins = ns(1)*256+ns(2)*16+ns(3);
nxtS(i+1,in+1)=find(ins==I)-1;
ns = sort([S(i+1,1) in S(i+1,3)]')';
ins = ns(1)*256+ns(2)*16+ns(3);
nxtS(i+1,in+16+1)=find(ins==I)-1;
ns = sort([S(i+1,1) S(i+1,2) in]')';
ins = ns(1)*256+ns(2)*16+ns(3);
nxtS(i+1,in+32+1)=find(ins==I)-1;
end
end
n=[0:15]';
vol=2.^(n/2)/2^7.5;
vol(1)=0;
curV = zeros(Ns,3*16);
for i=0:Ns-1
for in=0:15
curV(i+1,in+ +1) = vol(in+1) + vol(bitand(fix(i/16),15)+1) + vol(bitand(i,15)+1);
curV(i+1,in+16+1) = vol(fix(i/256)+1) + vol(in+1) + vol(bitand(i,15)+1);
curV(i+1,in+32+1) = vol(fix(i/256)+1) + vol(bitand(fix(i/16),15)+1) + vol(in+1);
end
end
f = 3579545;
inputfile = 'a1.wav';
inputfile = lower(inputfile);
if isempty(findstr(inputfile,'.wav'))
inputfile = strcat(inputfile,'.wav');
end
[y,Fs,nbit,opt]=wavread(inputfile);
%y = y(1:50);
y = (y-min(y))/(max(y)-min(y)) * 1.15; % input in 0 e 1.15
f = 3579545;
fpcm = Fs;
tp = 447; % 3 PSG transitions per sample: put here the tp of your replayer
dt1 = (12+5+14)/tp;
dt2 = (12+14)/tp;
dt3 = 1-dt1-dt2;
dt = [dt1 dt2 dt3];
N = 3*length(y);
x = zeros(1,N);
x(1:3:end) = y;
x(2:3:end-3) = y(1:end-1)*dt1+y(2:end)*(1-dt1);
x(3:3:end-3) = y(1:end-1)*(dt1+dt2)+y(2:end)*(1-dt1-dt2);
x(end-2:end) = y(end);
disp(sprintf(' This corresponds to an input with samplig frequency = %3.1f Hz',f/tp))
Stt = uint16(zeros(Ns,N));
Itt = uint8(zeros(Ns,N));
L = zeros(Ns,1);
St = uint16(ones(Ns,1));
It = uint8(ones(Ns,1));
for t =1:N
Ln = inf*ones(Ns,1);
if mod(t-1,1000)==0
disp(sprintf(' Processing complete at %3.2f %%',(t/N*100)));
end
for cs=0:Ns-1
for in=0:47
cv = curV(cs+1,in+1);
ns = double(nxtS(cs+1,in+1));
Ltst = L(cs+1)+dt(mod(t-1,3)+1)*abs(x(t)-cv)^2;
if Ln(ns+1) >= Ltst
Ln(ns+1) = Ltst;
St(ns+1) = cs;
It(ns+1) = in;
end
end
end
L = Ln;
Stt(:,t) = St;
Itt(:,t) = It;
end
[l,i] = min(L);
P = uint16(zeros(1,N));
I = uint8(zeros(1,N));
P(N) = Stt(i,N);
I(N) = Itt(i,N);
for t = (N-1):-1:1
P(t) = Stt(double(P(t+1))+1,t);
I(t) = Itt(double(P(t+1))+1,t);
end
V = zeros(1,N);
for t = 1:N
V(t) = curV(double(P(t))+1,double(I(t))+1);
end
%er = l;
en = sum(abs(x(1:3:end)).^2)*dt1+...
sum(abs(x(2:3:end)).^2)*dt2+...
sum(abs(x(3:3:end)).^2)*dt3+...
-mean(x)^2;
er = sum(abs(x(1:3:end)-V(1:3:end)).^2)*dt1+...
sum(abs(x(2:3:end)-V(2:3:end)).^2)*dt2+...
sum(abs(x(3:3:end)-V(3:3:end)).^2)*dt3;
SNR = 10*log10(en/er);
disp(sprintf(' SNR = %2.1f dB',SNR))
plot(x,'r')
hold on
plot(V)
% meaning of I
% I(t) in [0 15] change the min value of the triplet in I(t)
% I(t) in [16 31] change the middle value of the triplet in I(t)-16
% I(t) in [32 47] change the max value of the triplet in I(t)-32
Out = zeros(size(I));
S = [bitand(fix(double(P(1))/256),15) bitand(fix(double(P(1))/16),15) bitand(fix(double(P(1))),15)];
for t = 1:N
if (I(t)<=15)
[u,i] = min(S);
S(i) = double(I(t));
Out(t) = (i-1)*64 + S(i);
elseif (I(t)<=31) & (I(t)>=16)
[u,i] = min(S);
SS = S;
SS(i) = inf;
[u,i] = min(SS); % i is the second smallest
S(i) = double(I(t))-16 ;
Out(t) = (i-1)*64 + S(i);
elseif (I(t)<=47) & (I(t)>=32)
[u,i] = max(S);
S(i) = double(I(t))-32 ;
Out(t) = (i-1)*64 + S(i);
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Out has the channel # in b7b8 and the level in b3-b0
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
I am a bit not confident about the initial state of S in the loop
that compute Out. Probably I should reconsider that value...
| | ARTRAG
msx master
Posts: 1578 | Posted: January 25 2006, 09:27 | OK, after a small bug correction, here it is the encoder
complete of simulated ASM decoder.
THE SNR IS 43,1dB !!!!!
(correctly computed without the power of mean!!)
u = 1;
S = zeros(16*16*16,3);
for i=0:15
for j=0:15
for k=0:15
S(u,:)=[ i j k ];
u=u+1;
end
end
end
S = sort(S')';
I = S(:,1)*256+S(:,2)*16+S(:,3);
[l i] = sort(I);
S = S(i,:);
I = S(:,1)*256+S(:,2)*16+S(:,3);
i = [ find(diff(I)) ; 4096];
S = S(i,:);
I = I(i);
Ns = length(I);
nxtS = zeros(Ns,3*16);
for i=0:Ns-1
for in=0:15
ns = sort([in S(i+1,2:3)]);
ins = ns*[256 ; 16 ;1 ];
nxtS(i+1,in+1)=find(ins==I)-1;
ns = sort([S(i+1,1) in S(i+1,3)]);
ins = ns*[256 ; 16 ;1 ];
nxtS(i+1,in+16+1)=find(ins==I)-1;
ns = sort([S(i+1,1:2) in]);
ins = ns*[256 ; 16 ;1 ];
nxtS(i+1,in+32+1)=find(ins==I)-1;
end
end
n=[0:15]';
vol=2.^(n/2)/2^7.5;
vol(1)=0;
curV = zeros(Ns,3*16);
for i=0:Ns-1
for in=0:15
curV(i+1,in+ +1) = sum(vol(S(nxtS(i+1,in+ +1)+1,:)+1));
curV(i+1,in+16+1) = sum(vol(S(nxtS(i+1,in+16+1)+1,:)+1));
curV(i+1,in+32+1) = sum(vol(S(nxtS(i+1,in+32+1)+1,:)+1));
end
end
f = 3579545;
inputfile = 'a1.wav';
inputfile = lower(inputfile);
if isempty(findstr(inputfile,'.wav'))
inputfile = strcat(inputfile,'.wav');
end
[y,Fs,nbit,opt]=wavread(inputfile);
y = (y-min(y))/(max(y)-min(y)) * 1.15; % input in 0 e 1.15
f = 3579545;
fpcm = Fs;
tp = 447; % 3 PSG transitions per sample: put here the tp of your replayer
dt1 = (12+5+14)/tp;
dt2 = (12+14)/tp;
dt3 = 1-dt1-dt2;
dt = [dt1 dt2 dt3];
N = 3*length(y);
x = zeros(1,N);
x(1:3:end) = y;
x(2:3:end-3) = y(1:end-1)*dt1+y(2:end)*(1-dt1);
x(3:3:end-3) = y(1:end-1)*(dt1+dt2)+y(2:end)*(1-dt1-dt2);
x(end-2:end) = y(end);
disp(sprintf(' This corresponds to an input with samplig frequency = %3.1f Hz',f/tp))
Stt = uint16(zeros(Ns,N));
Itt = uint8(zeros(Ns,N));
L = zeros(Ns,1);
St = uint16(ones(Ns,1));
It = uint8(ones(Ns,1));
for t =1:N
Ln = inf*ones(Ns,1);
if mod(t-1,1000)==0
disp(sprintf(' Processing complete at %3.2f %%',(t/N*100)));
end
for cs=0:Ns-1
for in=0:47
cv = curV(cs+1,in+1);
ns = double(nxtS(cs+1,in+1));
Ltst = L(cs+1)+dt(mod(t-1,3)+1)*abs(x(t)-cv)^2;
if Ln(ns+1) >= Ltst
Ln(ns+1) = Ltst;
St(ns+1) = cs;
It(ns+1) = in;
end
end
end
L = Ln;
Stt(:,t) = St;
Itt(:,t) = It;
end
disp(sprintf(' Processing complete at %3.2f %%',100));
[l,i] = min(L);
P = uint16(zeros(1,N));
I = uint8(zeros(1,N));
P(N) = Stt(i,N);
I(N) = Itt(i,N);
for t = (N-1):-1:1
P(t) = Stt(double(P(t+1))+1,t);
I(t) = Itt(double(P(t+1))+1,t);
end
V = zeros(1,N);
for t = 1:N
V(t) = curV(double(P(t))+1,double(I(t))+1);
end
%er = l;
en = sum(abs(x(1:3:end)).^2)*dt1+...
sum(abs(x(2:3:end)).^2)*dt2+...
sum(abs(x(3:3:end)).^2)*dt3+...
-mean(x(1:3:end))^2*dt1+...
-mean(x(2:3:end))^2*dt2+...
-mean(x(3:3:end))^2*dt3;
er = sum(abs(x(1:3:end)-V(1:3:end)).^2)*dt1+...
sum(abs(x(2:3:end)-V(2:3:end)).^2)*dt2+...
sum(abs(x(3:3:end)-V(3:3:end)).^2)*dt3;
SNR = 10*log10(en/er);
disp(sprintf(' SNR = %2.1f dB',SNR))
plot(x,'r')
hold on
plot(V)
% meaning of I
% I(t) in [0 15] change the min value of the triplet in I(t)
% I(t) in [16 31] change the middle value of the triplet in I(t)-16
% I(t) in [32 47] change the max value of the triplet in I(t)-32
Out = zeros(size(I));
S = [9 12 15];
for t = 1:N
[u,i] = sort(S);
if (I(t)<=15)
i = i(1);
S(i) = double(I(t));
Out(t) = (i-1)*64 + S(i);
elseif (I(t)<=31) & (I(t)>=16)
i = i(2);
S(i) = double(I(t))-16 ;
Out(t) = (i-1)*64 + S(i);
elseif (I(t)<=47) & (I(t)>=32)
i = i(3);
S(i) = double(I(t))-32 ;
Out(t) = (i-1)*64 + S(i);
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Out has the channel # in b7b8 and the level in b3-b0
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% ASM replayer simulation
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Q = zeros(size(Out));
S = [9 12 15];
S = vol(S+1)';
for t = 1:N
c = fix(Out(t)/64);
S(c+1) = vol(bitand(Out(t),15)+1);
Q(t)= sum(S);
end
plot (Q,'.k')
The only unsolved problem is to find the initial state of the PSG
that in any case becomes irrelevant in 6-8 steps!! | | Arjan
msx addict
Posts: 453 | Posted: January 25 2006, 12:20 | Since you still have two empty bits ber byte, you could use them to indicate if there is more data needed to complete the sample transition. I don't know if it's worthwhile though, but I suppose you don't always need a triplet to change volume from one sample to another.
Or even better:
first byte: bit 6-7 is indicates first channel, bit 4-5 indicates second channel. bit 0-3 is volume data for first channel change
next byte: contains volume data for the other channels.
You could also use the value of bit4-5 to indicate you only need to change the volume of one channel (probably doesn't happen alot though...). If it's value is 11, no more data is needed to complete the sample transition.
| | ARTRAG
msx master
Posts: 1578 | Posted: January 25 2006, 17:45 | Quote:
| Since you still have two empty bits ber byte, you could use them to indicate if there is more data needed to complete the sample transition. I don't know if it's worthwhile though, but I suppose you don't always need a triplet to change volume from one sample to another.
|
More PSG transitions => better sample accuracy
Quote:
|
Or even better:
first byte: bit 6-7 is indicates first channel, bit 4-5 indicates second channel. bit 0-3 is volume data for first channel change
next byte: contains volume data for the other channels.
|
I should have a more aggressive compression : even assumig perfect bit packing you have 4+2=6 bit per PSG level, 3*6=18bit per input sample, 793800bit/sec or 99225
byte/sec for inputs a 44,1KHz.
It is 100kbyte per second!!! No MSX has enough ram for this replayer!!
Quote:
|
You could also use the value of bit4-5 to indicate you only need to change the volume of one channel (probably doesn't happen alot though...). If it's value is 11, no more data is needed to complete the sample transition.
|
IMHO variable number of transitions per sample increases complexity and doesn't help
| | dvik
msx master
Posts: 1302 | Posted: January 25 2006, 18:49 | Yeah, playing a 44.1kHz sample that requires 3 PSG channel updates per sample is impossible on an MSX. A 44.1kHz replayer with one channel update per sample for the above encoding uses almost all available CPU time:
;-------------------------------------
; Plays one sample
; IN HL - Encoded sample start address
; DE - Sample length (#pcm samples)
;-------------------------------------
PLAY_SAMPLE:
ld b,e
inc d
ld c,$a1
PsgLoop:
; Output one channel
ld a,(hl) ; 8
rlca ; 5
rlca ; 5
and 3 ; 8
add 8 ; 8
out ($a0),a ; 12
nop ; 5
outi ; 18
jr nz,PsgLoop ; 8/13 Total: 82 (81 is best fit) 3579545/82 => 43653Hz
dec d ; 5
jp nz,PsgLoop ; 11
ret
| | msd
msx professional
Posts: 600 | Posted: January 25 2006, 18:57 | He my MSX could play 42 seconds on 44.1Khz | | dvik
msx master
Posts: 1302 | Posted: January 25 2006, 19:05 | But it would sound very good  Almost like an Ipod. Only a bit heavier to carry arround. | |
| |
| |
