![Аватар пользователя [WYZ]](https://es.msx.org/sites/default/files/imagecache/avatar-comment/avatars/picture-775.png?v=1470481363 "Аватар пользователя [WYZ]")
The best sample I ever heard on MSX before. Very hight quality and noiseless. Congratulations Artrag and Dvik.
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.
but it works!!!
About CPU time, I think that there is spare time for many things,
maybe also for an encoding better than RLE.
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?
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
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 If 31 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 PS
ordered ones. This allows me to reconstruct a actual channel nmber to be stored in output. I will add this paer asap AR
Note that the remarks in the file about I are wrong!!
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...
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!!
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.