summaryrefslogtreecommitdiff
path: root/macros/tf2zp.sci
blob: f607ba0acbdc27271add49e5a1e4c97c519f7f6b (plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
// Copyright (C) 2018 - IIT Bombay - FOSSEE
// This file must be used under the terms of the CeCILL.
// This source file is licensed as described in the file COPYING, which
// you should have received as part of this distribution.  The terms
// are also available at
// http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
// Original Source : https://octave.sourceforge.io/signal/
// Modifieded by: Abinash Singh Under FOSSEE Internship
// Date of Modification: 3 Feb 2024
// Organization: FOSSEE, IIT Bombay
// Email: toolbox@scilab.in
function [z,p,k]=tf2zp(num,den)

    // Transfer function to zero pole conversion

    //Calling Sequence
    //[z,p,k]= tf2zp(b,a);

    //Parameters
    //z=zeros of the corrsponding tf
    //p=poles of the corresponding tf
    //k=gain of the tf
    //b=vector containing the numerator coefficients of the transfer function in descending powers of s/z
    //a=vector containing the denominator coefficients of the transfer function in descending powers of s/z

    //For discrete-time transfer functions, it is highly recommended to
    //make the length of the numerator and denominator equal to ensure
    //correct results.  You can do this using the function EQTFLENGTH in
    //the Signal Processing Toolbox.

    //Author : Debdeep Dey

    //Example :
    //b = [1 2 3];
    //a = [4 5 6];
    //[z p k] = tf2zp(b,a)
    //Output:
    //k  =
    //
    //    0.25
    // p  =
    //
    //  - 0.625 + 1.0532687i
    //  - 0.625 - 1.0532687i
    // z  =
    //
    //  - 1. + 1.4142136i
    //  - 1. - 1.4142136i

    numop=argn(1);
    //take only the first row of numerator into consideration
    num=num(1,:);
    //remove leading columns of zeros from numerator
    if ~isempty(num) then
        while(num(:,1)==0 & length(num)>1)
            num(:,1)=[];
        end
    end
    [rd,cod]=size(den);
    [ny,np]=size(num);

    if(rd>1) then
        error("Denominator must be row vector");
        
    // elseif np>cod then
    //      error("Improper transfer function"); it is a bug it should be removed
    end
    if (~isempty(den)) then
        coef=den(1);
    else
        coef=1;
    end
    if coef ==0 then
        error("Denominator must have non zero leading coefficient");
    end
    //remove leading columns of zeros from numerator
    if ~isempty(num) then
        while(num(:,1)==0 & length(num)>1)
            num(:,1)=[];
        end
    end

    if (find(den==%inf) ~= [] | find(num==%inf) ~= []) then
        error("Input to ROOTS must not contain NaN or Inf")
    end
    //poles

    p=roots(den);
    //zeros and gain

    k=zeros(ny,1);
    linf=%inf;

    z=linf(ones(np-1,1),ones(ny,1));
    for i=1:ny
        zz=roots(num(i,:));
        if ~isempty(zz), z(1:length(zz), i) = zz; end
        ndx = find(num(i,:)~=0);
        if ~isempty(ndx), k(i,1) = num(i,ndx(1))./coef; end
    end
    z=roots(num);
endfunction
/*
Note : This function is tested with Octave's outputs as a reference.
[z p k] = tf2zp([9.6500e+02  -1.1773e+05   3.8648e+06  -4.5127e+071.5581e+08],[ 1  -3   2]) //pass
[z p k] = tf2zp([1 2 3 4],[ 1  -3   2]) //pass
[z p k] = tf2zp([4 5 6],[1]) // pass
*/