function MKS_06a_v13_Gain
% MKS_06a – SUBTRACTIVE SYNTHESISER
% Dva nezavisle filtry FIR/IIR: LP, HP, BP, BS, RES (all-pole); SINGLE, CASCADE, PARALLEL.
% Bez Audio Toolboxu a Signal Processing Toolboxu. MATLAB audioplayer.
% Klavesnice QWERTZ: A W S E D F T G Z H U J K (13 pul tonu).
% Tato verze prehrava hotove tony; nejde o blokovy real-time engine.

fs=8000; zaklad=220; delka=1.2; klavesy='awsedftgzhujk';
barva=[.96 .97 .98];
S.fs=fs; S.klavesy=klavesy; S.hrac=[]; S.zaznam=[]; S.nahrava=false;
S.posledni=[]; S.hraje=false; S.cache=cell(1,13); S.pressed=false(1,13);
fig=figure('Name','MKS 06 - SUBTRACTIVE SYNTHESISER','NumberTitle','off', ...
 'MenuBar','none','ToolBar','none','Color',barva,'Position',[70 65 1400 810], ...
 'WindowKeyPressFcn',@stiskKlavesy,'WindowKeyReleaseFcn',@pustKlavesy, ...
 'CloseRequestFcn',@zavri);

uicontrol(fig,'Style','text','String','MKS 06  |  SUBTRACTIVE SYNTHESIS', ...
 'FontSize',19,'FontWeight','bold','BackgroundColor',barva, ...
 'Position',[35 757 1300 38]);

pZdroj=uipanel(fig,'Title','SOURCE / ZDROJ','Units','pixels','Position',[25 548 300 192]);
pF1=uipanel(fig,'Title','FILTER F1','Units','pixels','Position',[342 548 500 192]);
pF2=uipanel(fig,'Title','FILTER F2','Units','pixels','Position',[860 548 500 192]);
pEnv=uipanel(fig,'Title','AMPLITUDE ENVELOPE (ADSR)','Units','pixels','Position',[25 375 460 155]);
pMode=uipanel(fig,'Title','FILTER CONNECTION','Units','pixels','Position',[505 375 355 155]);
pTransport=uipanel(fig,'Title','AUDIO / RECORDING','Units','pixels','Position',[880 375 480 155]);

textAt(pZdroj,'Waveform',12,137,90);
S.src=popup(pZdroj,{'SAW','SINE','TRIANGLE','PULSE','NOISE'},110,134,160);
textAt(pZdroj,'Pulse width (%)',12,99,130);
S.pw=editAt(pZdroj,'50',150,97,110);
textAt(pZdroj,'f0 (Hz)',12,61,100);
S.f0=editAt(pZdroj,'220',150,59,110);
textAt(pZdroj,'fs (Hz)',12,23,100);
S.fsChoice=popup(pZdroj,{'8000','16000','22050','44100','48000'},150,19,110);
textAt(pTransport,'Tone duration (s)',10,14,130);
S.duration=editAt(pTransport,num2str(delka),140,10,90);
S.F1=filterControls(pF1,1); S.F2=filterControls(pF2,2);
textAt(pEnv,'Attack (s)',10,99,90);S.A=editAt(pEnv,'0.03',110,97,85);
textAt(pEnv,'Decay (s)',220,99,90);S.D=editAt(pEnv,'0.20',320,97,85);
textAt(pEnv,'Sustain (0-1)',10,56,100);S.S=editAt(pEnv,'0.65',110,54,85);
textAt(pEnv,'Release (s)',220,56,95);S.R=editAt(pEnv,'0.25',320,54,85);
textAt(pEnv,'ADSR applies after filtering',10,16,360);

textAt(pMode,'Connection',15,98,110);
S.mode=popup(pMode,{'SINGLE (F1)','CASCADE (F1 -> F2)','PARALLEL (G1*F1 + G2*F2)/2'},120,94,210);
set(S.mode,'Value',3); % Start with both filters contributing to the output
textAt(pMode,'The plots show F1, F2 and their combined response.',15,46,330);
textAt(pMode,'RES = all-pole IIR; other IIR = biquads.',15,19,330);

button(pTransport,'PLAY f0',10,93,105,@prehratA3);
button(pTransport,'REC',127,93,90,@zapniRec);
button(pTransport,'STOP',230,93,90,@zastav);
button(pTransport,'PLAY REC',333,93,110,@prehratRec);
button(pTransport,'SAVE WAV',10,43,180,@uloz);
S.status=uicontrol(pTransport,'Style','text','String','Ready','Position',[245 10 210 30], ...
 'HorizontalAlignment','left');

S.ax1=axes('Parent',fig,'Units','pixels','Position',[65 192 390 158]);
S.ax2=axes('Parent',fig,'Units','pixels','Position',[545 192 390 158]);
S.ax3=axes('Parent',fig,'Units','pixels','Position',[1025 192 315 158]);
S.piano=uipanel(fig,'Title','KLAVIATURA  |  a - a1  |  A W S E D F T G Z H U J K', ...
 'Units','pixels','Position',[75 22 1250 145]);
S.keyButtons=gobjects(1,13);
S.keyColors=zeros(13,3);
bile=[1 3 5 6 8 10 12 13];
cerne=[2 4 7 9 11];
sirka=149; vlevo=18;
for j=1:8
    k=bile(j);
    S.keyButtons(k)=uicontrol(S.piano,'Style','pushbutton', ...
        'String',upper(klavesy(k)),'Position',[vlevo+(j-1)*sirka 10 sirka-3 100], ...
        'BackgroundColor',[1 1 1],'ForegroundColor',[0 0 0], ...
        'FontSize',13,'Callback',{@klikNota,k});
    S.keyColors(k,:)=[1 1 1];
end
zaBilou=[1 2 4 5 6];
for j=1:5
    k=cerne(j);
    S.keyButtons(k)=uicontrol(S.piano,'Style','pushbutton', ...
        'String',upper(klavesy(k)), ...
        'Position',[vlevo+zaBilou(j)*sirka-43 51 86 59], ...
        'BackgroundColor',[.13 .15 .19],'ForegroundColor',[1 1 1], ...
        'FontSize',12,'Callback',{@klikNota,k});
    S.keyColors(k,:)=[.13 .15 .19];
end
set([S.src S.pw S.duration S.A S.D S.S S.R S.mode], 'Callback',@obnovGrafy);
for k=1:2
    q=S.F1;if k==2,q=S.F2;end
    set(q.typ,'Callback',{@zmenaTyp,k});
    set([q.metoda q.rad q.fc q.bw q.gain],'Callback',@obnovGrafy);
end
set(S.fsChoice,'Callback',@zmenaFs);
set(S.f0,'Callback',@zmenaF0);
guidata(fig,S);obnovGrafy();

    function q=filterControls(p,k)
        textAt(p,'Type',12,130,65);q.typ=popup(p,{'LP','HP','BP','BS','RES (2-pole)'},75,127,105);
        textAt(p,'Design',207,130,65);q.metoda=popup(p,{'IIR','FIR'},275,127,80);
        textAt(p,'Order',365,130,52);q.rad=editAt(p,'4',419,127,55);
        textAt(p,'Cutoff / center (Hz)',12,80,160);
        if k==1,fc='1200';else,fc='2400';end
        q.fc=editAt(p,fc,177,78,100);
        textAt(p,'Bandwidth (Hz)',285,80,115);q.bw=editAt(p,'600',399,78,75);
        textAt(p,'Gain (dB)',12,29,110);q.gain=editAt(p,'0',122,27,75);
        textAt(p,'RES: all-pole IIR, order 2',215,29,260);
    end
    function zmenaTyp(~,~,k)
        s=guidata(fig);
        if k==1,q=s.F1;else,q=s.F2;end
        if get(q.typ,'Value')==5
            set(q.rad,'String','2');
            set(q.metoda,'Value',1);
        end
        obnovGrafy();
    end

    function h=popup(parent,items,x,y,w)
        h=uicontrol(parent,'Style','popupmenu','String',items,'Position',[x y w 26]);
    end
    function h=editAt(parent,val,x,y,w)
        h=uicontrol(parent,'Style','edit','String',val,'Position',[x y w 27]);
    end
    function textAt(parent,str,x,y,w)
        uicontrol(parent,'Style','text','String',str,'Position',[x y w 22], ...
          'HorizontalAlignment','left');
    end
    function button(parent,str,x,y,w,cb)
        uicontrol(parent,'Style','pushbutton','String',str,'Position',[x y w 34], ...
          'FontWeight','bold','Callback',cb);
    end
    function zmenaFs(varargin)
        s=guidata(fig);
        moznosti=[8000 16000 22050 44100 48000];
        noveFs=moznosti(get(s.fsChoice,'Value'));
        if noveFs==fs,return;end
        if ~isempty(s.hrac),stop(s.hrac);end
        fs=noveFs;
        s.fs=fs;
        zaklad=cislo(s.f0,zaklad,20,fs/2-100);
        s.cache=cell(1,13);
        s.zaznam=[];
        s.posledni=[];
        s.nahrava=false;
        guidata(fig,s);
        obnovGrafy();
        s=guidata(fig);
        set(s.status,'String',sprintf('fs = %d Hz; cache a záznam vymazány',fs));
    end
    function zmenaF0(varargin)
        s=guidata(fig);
        zaklad=cislo(s.f0,zaklad,20,fs/2-100);
        s.cache=cell(1,13);
        guidata(fig,s);
        set(s.status,'String',sprintf('f0 = %.1f Hz | klaviatura: f0 az 2*f0',zaklad));
    end
    function v=cislo(h,def,lo,hi)
        v=str2double(strrep(get(h,'String'),',','.'));
        if ~isfinite(v),v=def;end
        v=max(lo,min(hi,v));
        set(h,'String',num2str(v));
    end
    function [b,a]=navrh(q)
        typ=get(q.typ,'Value');met=get(q.metoda,'Value');
        rad=round(cislo(q.rad,4,1,100));
        fc=cislo(q.fc,1200,30,fs/2-100);
        bw=cislo(q.bw,600,30,fs/2-200);
        if typ==5
            % All-pole formant resonator, no finite-frequency zeros.
            % Pole radius controls bandwidth; normalize peak gain to 1.
            r=exp(-pi*bw/fs); w0=2*pi*fc/fs;
            a=[1 -2*r*cos(w0) r^2];
            b=1;
            h0=abs(odpoved(b,a,w0));
            b=1/max(h0,realmin);
            return;
        end
        if typ>=3
            lo=max(30,fc-bw/2);hi=min(fs/2-100,fc+bw/2);
            if hi<=lo+20,hi=lo+20;end
        end
        if met==2
            % FIR s oknem Hamming; bez fir1 a bez toolboxu.
            % Pro HP a BS musi byt rad sudy (lichy pocet koeficientu).
            if (typ==2 || typ==4) && mod(rad,2)==1,rad=rad+1;end
            n=(-rad/2:rad/2)';
            if typ==1
                b=2*fc/fs*sinc0(2*fc/fs*n);
            elseif typ==2
                b=-2*fc/fs*sinc0(2*fc/fs*n);b(rad/2+1)=b(rad/2+1)+1;
            elseif typ==3
                b=2*hi/fs*sinc0(2*hi/fs*n)-2*lo/fs*sinc0(2*lo/fs*n);
            else
                b=2*lo/fs*sinc0(2*lo/fs*n)-2*hi/fs*sinc0(2*hi/fs*n);
                b(rad/2+1)=b(rad/2+1)+1;
            end
            win=.54-.46*cos(2*pi*(0:rad)'/rad);
            b=(b.*win)';a=1;
            if typ==1,ref=0;elseif typ==2,ref=pi;elseif typ==3,ref=2*pi*fc/fs;else,ref=0;end
            H=sum(b.*exp(-1i*ref*(0:rad)));
            if abs(H)>1e-8,b=b/abs(H);end
        else
            % IIR Butterworth-like cascade of RBJ biquads; order controls
            % number of second-order sections (order 2,4,..., max 20).
            nsec=max(1,min(10,ceil(rad/2)));
            b=1;a=1;
            for j=1:nsec
                if typ==1 || typ==2
                    % Butterworth pole Q for nsec second-order sections.
                    Q=1/(2*sin((2*j-1)*pi/(4*nsec)));
                else
                    Q=max(.15,fc/max(bw,1));
                end
                w0=2*pi*fc/fs; alpha=sin(w0)/(2*Q);
                if typ==1
                    bb=[(1-cos(w0))/2,1-cos(w0),(1-cos(w0))/2];
                elseif typ==2
                    bb=[(1+cos(w0))/2,-(1+cos(w0)),(1+cos(w0))/2];
                elseif typ==3
                    bb=[alpha,0,-alpha];
                else
                    bb=[1,-2*cos(w0),1];
                end
                aa=[1+alpha,-2*cos(w0),1-alpha];
                bb=bb/aa(1);aa=aa/aa(1);
                b=conv(b,bb);a=conv(a,aa);
            end
        end
    end
    function z=sinc0(x)
        z=ones(size(x));ii=abs(x)>1e-12;z(ii)=sin(pi*x(ii))./(pi*x(ii));
    end
    function [b,a,b1,a1,b2,a2]=koeficienty()
        s=guidata(fig);[b1,a1]=navrh(s.F1);[b2,a2]=navrh(s.F2);
        g1=10^(cislo(s.F1.gain,0,-30,30)/20);
        g2=10^(cislo(s.F2.gain,0,-30,30)/20);
        b1=g1*b1;b2=g2*b2;
        switch get(s.mode,'Value')
            case 1,b=b1;a=a1;
            case 2,b=conv(b1,b2);a=conv(a1,a2);
            case 3
                t1=conv(b1,a2);t2=conv(b2,a1);
                L=max(numel(t1),numel(t2));
                b=([t1 zeros(1,L-numel(t1))]+[t2 zeros(1,L-numel(t2))])/2;
                a=conv(a1,a2);
        end
    end
    function obnovGrafy(varargin)
        if ~isgraphics(fig),return;end
        s=guidata(fig);
        s.cache=cell(1,13);guidata(fig,s);
        try
            [b,a,b1,a1,b2,a2]=koeficienty();
            f=linspace(0,fs/2,1200);w=2*pi*f/fs;
            h1=odpoved(b1,a1,w);h2=odpoved(b2,a2,w);
            switch get(s.mode,'Value')
                case 1,h=h1;
                case 2,h=h1.*h2;
                case 3,h=(h1+h2)/2;
            end
            cla(s.ax1);plot(s.ax1,f,20*log10(max(abs(h1),1e-5)),'b', ...
                f,20*log10(max(abs(h2),1e-5)),'r');
            title(s.ax1,'F1 (blue) / F2 (red)');xlabel(s.ax1,'Hz');ylabel(s.ax1,'dB');
            grid(s.ax1,'on');ylim(s.ax1,[-85 25]);
            xlim(s.ax1,[0 fs/2]);
            gain=abs(h);
            gain=gain/max(max(gain),realmin); % Normalized combined magnitude: peak = 0 dB
            cla(s.ax2);plot(s.ax2,f,20*log10(max(gain,1e-5)),'k','LineWidth',1.4);
            title(s.ax2,'COMBINED RESPONSE');xlabel(s.ax2,'Hz');ylabel(s.ax2,'dB');
            grid(s.ax2,'on');ylim(s.ax2,[-85 25]);
            xlim(s.ax2,[0 fs/2]);
            % Z-plane: exact transfer after cancellation, with branch pole colors.
            % Equal-degree padding is essential: roots([1 0 0]) has two
            % zeros at the origin, whereas roots(1) has none.
            cla(s.ax3);hold(s.ax3,'on');
            theta=linspace(0,2*pi,361);
            plot(s.ax3,cos(theta),sin(theta),'k:','LineWidth',1.2);
            plot(s.ax3,[-1.4 1.4],[0 0],'-','Color',[.75 .75 .75]);
            plot(s.ax3,[0 0],[-1.4 1.4],'-','Color',[.75 .75 .75]);
            [zz,pp,valid,reason]=skutecneKoreny(b,a);
            if valid
                % The zeros are those of the RESULT, not the branches.
                if ~isempty(zz)
                    plot(s.ax3,real(zz),imag(zz),'ok','MarkerSize',7,'LineWidth',1.5);
                end
                % Color resulting poles by the branch from which they come.
                p1=roots(a1);p2=roots(a2);
                rezim=get(s.mode,'Value');
                for ip=1:numel(pp)
                    pole=pp(ip);
                    d1=inf;d2=inf;
                    if ~isempty(p1),d1=min(abs(p1-pole));end
                    if rezim~=1 && ~isempty(p2),d2=min(abs(p2-pole));end
                    if d1<=d2,col=[0 .3 .95];else,col=[.9 .15 .1];end
                    plot(s.ax3,real(pole),imag(pole),'x','Color',col, ...
                        'MarkerSize',9,'LineWidth',1.7);
                end
                title(s.ax3,'RESULT: o zeros black | x poles F1 blue, F2 red');
            else
                title(s.ax3,'z-plane: not displayed (numerical reliability)');
                text(s.ax3,-1.25,0,reason,'FontSize',9);
            end
            axis(s.ax3,'equal');
            xlim(s.ax3,[-1.4 1.4]);ylim(s.ax3,[-1.4 1.4]);
            grid(s.ax3,'on');hold(s.ax3,'off');
            set(s.status,'String','Filters updated');
        catch ME
            set(s.status,'String',['ERROR: ' ME.message]);
        end
    end
    function [zz,pp,valid,reason]=skutecneKoreny(b,a)
        % H(z) = B(z^-1)/A(z^-1). Multiply both by z^N,
        % N=max(deg B,deg A), to obtain polynomials in z.
        % This retains genuine zeros/poles at z=0 from order differences.
        zz=[];pp=[];valid=false;reason='';
        b=b(:).';a=a(:).';
        if any(~isfinite(b)) || any(~isfinite(a)) || ~any(abs(b)>0) || ~any(abs(a)>0)
            reason='Invalid transfer coefficients';return;
        end
        % Remove trailing coefficients only if EXACTLY zero; a tiny
        % coefficient can be a genuine high-order term.
        while numel(b)>1 && b(end)==0,b(end)=[];end
        while numel(a)>1 && a(end)==0,a(end)=[];end
        N=max(numel(b),numel(a))-1;
        if N>12,reason='Polynomial degree > 12';return;end
        bz=[b zeros(1,N+1-numel(b))];
        az=[a zeros(1,N+1-numel(a))];
        zz=roots(bz);pp=roots(az);
        % Cancel only numerically coincident roots, including common
        % origin roots introduced by the equal-degree representation.
        used=false(size(pp));keep=true(size(zz));
        for ii=1:numel(zz)
            if isempty(pp),break;end
            dist=abs(pp-zz(ii));dist(used)=inf;
            [d,jj]=min(dist);
            if d<1e-7*max([1,abs(zz(ii)),abs(pp(jj))])
                keep(ii)=false;used(jj)=true;
            end
        end
        zz=zz(keep);pp=pp(~used);
        valid=true;
    end
    function H=odpoved(b,a,w)
        zb=exp(-1i*(0:numel(b)-1)'*w);
        za=exp(-1i*(0:numel(a)-1)'*w);
        H=(b(:).'*zb)./(a(:).'*za);
    end
    function [x,y]=ton(f)
        s=guidata(fig);N=round(fs*cislo(s.duration,1.2,.1,10));
        t=(0:N-1)'/fs;phase=2*pi*f*t;typ=get(s.src,'Value');
        switch typ
            case 1
                x=2*(mod(f*t,1)-.5);
            case 2,x=sin(phase);
            case 3
                x=4*abs(mod(f*t,1)-.5)-1;
            case 4
                pw=cislo(s.pw,50,5,95)/100;x=2*(mod(f*t,1)<pw)-1;
            otherwise,x=.4*randn(N,1);
        end
        [~,~,b1,a1,b2,a2]=koeficienty();
        y1=filter(b1,a1,x);y2=filter(b2,a2,x);
        switch get(s.mode,'Value')
            case 1,y=y1;
            case 2,y=filter(b2,a2,y1);
            case 3,y=(y1+y2)/2;
        end
        % ADSR na pevne delce tonu, posledni usek je release.
        A=cislo(s.A,.03,.001,3);D=cislo(s.D,.2,.001,3);
        Su=cislo(s.S,.65,0,1);R=cislo(s.R,.25,.001,3);
        nA=round(A*fs);nD=round(D*fs);nR=round(R*fs);
        if nA+nD+nR>N
            factor=N/(nA+nD+nR);nA=floor(nA*factor);
            nD=floor(nD*factor);nR=N-nA-nD;
        end
        nS=N-nA-nD-nR;
        env=[linspace(0,1,nA)';linspace(1,Su,nD)'; ...
            Su*ones(nS,1);linspace(Su,0,nR)'];
        y=y.*env;
        % Peak-normalize every generated note, including strongly attenuated filters.
        peak=max(abs(y));
        if isfinite(peak) && peak>1e-12,y=0.90*y/peak;end
    end
    function zahraj(f,k)
        s=guidata(fig);
        try
            if nargin>1 && ~isempty(s.cache{k})
                y=s.cache{k};
            else
                [~,y]=ton(f);
                if nargin>1,s.cache{k}=y;end
            end
            if ~isempty(s.hrac),stop(s.hrac);end
            s.hrac=audioplayer(y,fs);s.posledni=y;
            if s.nahrava,s.zaznam=[s.zaznam;y;zeros(round(.06*fs),1)];end
            if nargin>1 && s.pressed(k)
                set(s.keyButtons(k),'BackgroundColor',[.35 .70 .95]);
            end
            guidata(fig,s);
            play(s.hrac);
            set(s.status,'String',sprintf('Tón %.1f Hz | délka %.2f s',f,numel(y)/fs));
        catch ME
            set(s.status,'String',['ERROR: ' ME.message]);
        end
    end
    function stiskKlavesy(~,ev)
        s=guidata(fig);key=lower(ev.Key);
        if numel(key)~=1,return;end
        k=find(s.klavesy==key,1);
        if isempty(k) || s.pressed(k),return;end
        s.pressed(k)=true;guidata(fig,s);
        zahraj(zaklad*2^((k-1)/12),k);
    end
    function pustKlavesy(~,ev)
        s=guidata(fig);key=lower(ev.Key);
        if numel(key)~=1,return;end
        k=find(s.klavesy==key,1);
        if isempty(k),return;end
        s.pressed(k)=false;
        set(s.keyButtons(k),'BackgroundColor',s.keyColors(k,:));
        guidata(fig,s);
    end
    function klikNota(~,~,k)
        zahraj(zaklad*2^((k-1)/12),k);
        % Mouse click has no WindowKeyRelease event: restore key color.
        if isgraphics(fig)
            s=guidata(fig);
            if isgraphics(s.keyButtons(k))
                set(s.keyButtons(k),'BackgroundColor',s.keyColors(k,:));
            end
        end
    end
    function prehratA3(varargin),zahraj(zaklad,1);end
    function zapniRec(varargin)
        s=guidata(fig);s.zaznam=[];s.nahrava=true;guidata(fig,s);
        set(s.status,'String','REC ON: play notes');
    end
    function zastav(varargin)
        s=guidata(fig);s.nahrava=false;
        if ~isempty(s.hrac),stop(s.hrac);end
        guidata(fig,s);set(s.status,'String','Stopped / REC OFF');
    end
    function prehratRec(varargin)
        s=guidata(fig);if isempty(s.zaznam),set(s.status,'String','No recording');return;end
        if ~isempty(s.hrac),stop(s.hrac);end
        rec=s.zaznam;pk=max(abs(rec));
        if isfinite(pk) && pk>1e-12,rec=0.90*rec/pk;end
        s.hrac=audioplayer(rec,fs);guidata(fig,s);play(s.hrac);
    end
    function uloz(varargin)
        s=guidata(fig);if ~isempty(s.zaznam),y=s.zaznam;else,y=s.posledni;end
        if isempty(y),set(s.status,'String','Nothing to save');return;end
        [name,path]=uiputfile('*.wav','Save synthesized sound','MKS_06a.wav');
        if isequal(name,0),return;end
        peak=max(abs(y));if isfinite(peak) && peak>1e-12,y=0.90*y/peak;end
        audiowrite(fullfile(path,name),y,fs);
        set(s.status,'String',['Saved: ' name]);
    end
    function zavri(varargin)
        s=guidata(fig);if ~isempty(s.hrac),stop(s.hrac);end
        delete(fig);
    end
end