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
|
//Author: Parthasarathi Panda
//parthasarathipanda314@gmail.com
function [phi, varargout]=phasedelay(varargin)
//cas variable is 2 if sos form is involved and 1 if direct rational form is given
//(sos,n) or (sos,w) or (sos,'whole')or (b,a) is the input
//cas variable is 2 if sos form is involved and 1 if direct rational form is given
//cas1 variable is 1 if f is to be given as output, 2 other wise
[nargout,nargin]=argn();
//do not forget to execute 'phaseInputParseAs_sos' and 'phaseInputParseAs_ab' before running
v=size(varargin(1));
if size(v)>2 then
error ('invalid input dimension');
end
[n,k]=size(varargin(1));
if type(varargin(1))~=1 then
error ('check the input type');
end
if (n==1 & k==6) then //not clear if sos or (a,b)
v=size(varargin(2));
if (nargin==1) //(sos) is the input
cas=2;
[sos,w,cas1,fs]=phaseInputParseAs_sos(varargin,nargin);
elseif (varargin(2)=='whole') //(sos,'whole')is the input
cas=2;
[sos,w,cas1,fs]=phaseInputParseAs_sos(varargin,nargin);
else //taking it as (a,b)
cas=1;
[a,b,w,cas1,fs]=phaseInputParseAs_ab(varargin,nargin);
end
elseif (n==1 | k==1) then
cas=1;
[a,b,w,cas1,fs]=phaseInputParseAs_ab(varargin,nargin);
elseif k==6 then //first variable is sos
cas=2;
[sos,w,cas1,fs]=phaseInputParseAs_sos(varargin,nargin);
end
//cas,cas1,fs,w,[(a,b),sos]
if cas==1 then
[m,n]=size(a);
N=[0:n-1];
M=N'*w;//computing matrix Mij=(i-1)*wj
num=(a*exp(%i*M));//the operation computes phase of sum(ak*exp(i*w*k))
num_dot=(%i*(a.*N)*exp(%i*M));//computing the derivative on those points
phdel_num=(imag(num_dot).*real(num)-imag(num).*real(num_dot))./(abs(num).*abs(num));
[m,n]=size(b);
N=[0:n-1];
M=N'*w;
den=(b*exp(%i*M));//the operation computes phase of sum(ak*exp(i*w*k))
den_dot=(%i*(b.*N)*exp(%i*M));//computing the derivative on those points
phdel_den=(imag(den_dot).*real(den)-imag(den).*real(den_dot))./(abs(den).*abs(den));
phi=(phdel_num-phdel_den);
0;
else
[n,k]=size(sos)
N=[0,1,2];
M=N'*w;
num=(sos(:,4:6)*exp(%i*M));//the operation computes phase of sum(ak*exp(i*w*k))
num_dot=(%i*(sos(:,4:6).*(ones(n,1)*N))*exp(%i*M));//computing the derivative on those points
phdel_num=(imag(num_dot).*real(num)-imag(num).*real(num_dot))./(abs(num).*abs(num));
den=(sos(:,1:3)*exp(%i*M));//the operation computes phase of sum(ak*exp(i*w*k))
den_dot=(%i*(sos(:,1:3).*(ones(n,1)*N))*exp(%i*M));//computing the derivative on those points
phdel_den=(imag(den_dot).*real(den)-imag(den).*real(den_dot))./(abs(den).*abs(den));
phi_mat=phdel_num-phdel_den;
phi=sum(phi_mat,1);//summing each of the componenet second order system phases
end
if cas1==1 then
varargout(2)=w*fs/(2*%pi);
end
endfunction
|
