MATLAB: Help with creating a function out of a working script

MATLAB and Simulink Student Suitematlab function

Hey everyone.
I have a homework and I've never worked in MATLAB before. My script works as it should but the only problem is that I want to put it in a function (by that I mean that I want to write function(Amax, Amin, Fk, Fkappa) and the outputs should be coefficients, magnitude and phase frequency response of the proposed filter. It is basically a script, I will paste the most important stuff here (please ignore the comments, I just needed to comment it in my language)
prompt = 'Zadajte hodnotu Amax [dB] :  ' ;
Amax = input(prompt) ;                      
prompt1 = 'Zadajte hodnotu Amin [dB] : ' ;
Amin = input(prompt1) ;                     
prompt2 = 'Zadajte hodnotu Fkappa [kHz] : ' ;
Fkappa = input(prompt2) ;                   
prompt3 = 'Zadajte hodnotu Fk [kHz] :  ' ;
Fk = input(prompt3) ;                                                   
omegaK = (Fk/Fkappa) ;                         
Epsilon = sqrt(10^(Amax/10)-1) ;               
n = ((log10((10^(Amin/10)-1)/(10^(Amax/10)-1))/(2*log10(omegaK)))) ; 
nz = ceil(n) ;                                 
syms s ;                                       
vzorec = 1 + Epsilon^2 * s^(2*nz) * (-1)^nz ;  % Vzorec na výpočet koreňov
vypocet = solve(vzorec, s) ;                   % Pomocou funkcie solve vyrátame korene
korene = vpa(vypocet, 4) ;                     % Numerické aproximačné číselne korene na 4 des. cisla
%display(korene) ;                             % Vypíše nám korene na command window
X = korene(real(korene)<0) ;                   % Takto si do premennej X dam tie, kt. majú záporné reálne časti
X1 = double(X) ;                               % Dám si ich vypísať do X1 ako double (čísla) 
syms P ;                                       % Definujeme si premennú P
Gp1 = Epsilon ;                                % Do Gp1 si priradím Epsilon, pretože mi to inak nešlo
w = P/(2 * pi * Fkappa) ;                      % w = vzorec ktorý máme v zošite
for i = 1:nz                                
    Gp1 = Gp1 * (w - X(i)) ;                   % vzorec na vypočet Gp (u mňa Gp1)    
end
Gp1 = vpa(Gp1, 4) ;                            % Znova vpa, korene zaokrúhli na 4 čísla
Gp = expand(Gp1) ;                             % Expand nám roznásobí Gp (malý lifehack čo som našiel)
Hp = 1/Gp ;                                    % Prenosová funkcia Hp = 1/Gp
Hp = vpa(Hp, 4) ;                              % Znova nech nám zaokrúhli na 4 čísla
Nasa_Prenosova_Funkcia = Hp                    % Aby nám vypísalo nasu prenosovú funkciu
koef2 = coeffs(Gp) ;                           % Beriem si koeficienty z Gp, je to ľahšie ako z Hp
koef1 = fliplr(koef2) ;                         % Toto mi však vezme aj imag. čast, musim ošetriť:
koef = real(koef1) ;                            % Do koef si dám iba reálne časti polynómu
Prenosova_Funkcia_Koeficienty = vpa(koef)       % Necháme vypísať pre istotu
cit = 1 ;                                       % Citatel nasej prenosovej funkcie je 1
men = double(Prenosova_Funkcia_Koeficienty) ;   % menovatel bude číselne z prenosovej funkcie
B = tf(cit, men) ;                              % tf = prenosová funkcia z čitatela a menovatela
bode(B), grid ;                                 % vykreslenie Bodeho.
It is fully working, I just need to put it in a function, so it can be used in the command window as function(Amax, Amin, Fk, Fkappa) and the outputs should be the Prenosova_Funkcia_Koeficienty and the magnitude and phase frequency response (which is Bode(B).
I tried almost everything but I cannot figure it out.
Thanks for a response.

Best Answer

Read here: https://www.mathworks.com/help/matlab/ref/function.html to see how to define a function in MATLAB. Replace function_name with any name
function Prenosova_Funkcia_Koeficienty = function_name(Amax, Amin, Fk, Fkappa)
omegaK = (Fk/Fkappa) ;                         
Epsilon = sqrt(10^(Amax/10)-1) ;               
n = ((log10((10^(Amin/10)-1)/(10^(Amax/10)-1))/(2*log10(omegaK)))) ; 
nz = ceil(n) ;                                 
syms s ;                                       
vzorec = 1 + Epsilon^2 * s^(2*nz) * (-1)^nz ;  % Vzorec na výpočet koreňov
vypocet = solve(vzorec, s) ;                   % Pomocou funkcie solve vyrátame korene
korene = vpa(vypocet, 4) ;                     % Numerické aproximačné číselne korene na 4 des. cisla
%display(korene) ;                             % Vypíše nám korene na command window
X = korene(real(korene)<0) ;                   % Takto si do premennej X dam tie, kt. majú záporné reálne časti
X1 = double(X) ;                               % Dám si ich vypísať do X1 ako double (čísla) 
syms P ;                                       % Definujeme si premennú P
Gp1 = Epsilon ;                                % Do Gp1 si priradím Epsilon, pretože mi to inak nešlo
w = P/(2 * pi * Fkappa) ;                      % w = vzorec ktorý máme v zošite
for i = 1:nz                                
    Gp1 = Gp1 * (w - X(i)) ;                   % vzorec na vypočet Gp (u mňa Gp1)    
end
Gp1 = vpa(Gp1, 4) ;                            % Znova vpa, korene zaokrúhli na 4 čísla
Gp = expand(Gp1) ;                             % Expand nám roznásobí Gp (malý lifehack čo som našiel)
Hp = 1/Gp ;                                    % Prenosová funkcia Hp = 1/Gp
Hp = vpa(Hp, 4) ;                              % Znova nech nám zaokrúhli na 4 čísla
Nasa_Prenosova_Funkcia = Hp                    % Aby nám vypísalo nasu prenosovú funkciu
koef2 = coeffs(Gp) ;                           % Beriem si koeficienty z Gp, je to ľahšie ako z Hp
koef1 = fliplr(koef2) ;                         % Toto mi však vezme aj imag. čast, musim ošetriť:
koef = real(koef1) ;                            % Do koef si dám iba reálne časti polynómu
Prenosova_Funkcia_Koeficienty = vpa(koef)       % Necháme vypísať pre istotu
cit = 1 ;                                       % Citatel nasej prenosovej funkcie je 1
men = double(Prenosova_Funkcia_Koeficienty) ;   % menovatel bude číselne z prenosovej funkcie
B = tf(cit, men) ;                              % tf = prenosová funkcia z čitatela a menovatela
bode(B), grid ;                                 % vykreslenie Bodeho.
end
Then you can call it like this
prompt = 'Zadajte hodnotu Amax [dB] :  ' ;
Amax = input(prompt) ;                      
prompt1 = 'Zadajte hodnotu Amin [dB] : ' ;
Amin = input(prompt1) ;                     
prompt2 = 'Zadajte hodnotu Fkappa [kHz] : ' ;
Fkappa = input(prompt2) ;                   
prompt3 = 'Zadajte hodnotu Fk [kHz] :  ' ;
Fk = input(prompt3) ;
coeffs = function_name(Amax, Amin, Fk, Fkappa)
Note that the function will only plot the bode diagram. Do you want the frequency and phase response as numeric output?
Related Question