0
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
1 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
2 % TEST ndeigcsd
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
3 %
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
4 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
5 % HISTORY: 23-04-2009 L Ferraioli
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
6 % Creation
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
7 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
8
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
9
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
10
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
11 %% 2 Dim Test
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
12
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
13 %%% Define startinf TFs
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
14
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
15 fs = 10;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
16 f = [logspace(-6,log10(fs/2),2000)]';
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
17
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
18 % Model Stefano TF11
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
19 dRes11 = [2.44554138162509e-011 - 1.79482547894083e-011i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
20 2.44554138162509e-011 + 1.79482547894083e-011i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
21 2.66402334803101e-009 + 1.1025122049153e-009i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
22 2.66402334803101e-009 - 1.1025122049153e-009i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
23 -7.3560293387644e-009;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
24 -1.82811618589835e-009 - 1.21803627800855e-009i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
25 -1.82811618589835e-009 + 1.21803627800855e-009i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
26 1.16258677367555e-009;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
27 1.65216557639319e-016;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
28 -1.78092396888606e-016;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
29 -2.80420398962379e-017;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
30 9.21305973049041e-013 - 8.24686706827269e-014i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
31 9.21305973049041e-013 + 8.24686706827269e-014i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
32 5.10730060739905e-010 - 3.76571756625722e-011i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
33 5.10730060739905e-010 + 3.76571756625722e-011i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
34 3.45893698149735e-009;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
35 3.98139182134446e-014 - 8.25503935419059e-014i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
36 3.98139182134446e-014 + 8.25503935419059e-014i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
37 -1.40595719147164e-011];
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
38
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
39 dPoles11 = [0.843464045655194 - 0.0959986292915475i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
40 0.843464045655194 + 0.0959986292915475i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
41 0.953187595424927 - 0.0190043625473383i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
42 0.953187595424927 + 0.0190043625473383i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
43 0.967176277937188;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
44 0.995012027005247 - 0.00268322602801729i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
45 0.995012027005247 + 0.00268322602801729i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
46 0.996564761885673;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
47 0.999999366165445;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
48 0.999981722418555;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
49 0.999921882627659;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
50 0.999624431675213 - 0.000813407848742761i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
51 0.999624431675213 + 0.000813407848742761i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
52 0.997312006278751 - 0.00265611346834941i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
53 0.997312006278751 + 0.00265611346834941i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
54 0.990516544257531;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
55 0.477796923118318 - 0.311064085401834i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
56 0.477796923118318 + 0.311064085401834i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
57 0];
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
58
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
59 dDTerms11 = 0;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
60
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
61 % Model Stefano TF12
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
62 dRes12 = [1.44258422208796e-017 + 7.07359428613009e-019i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
63 1.44258422208796e-017 - 7.07359428613009e-019i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
64 -3.4918408053655e-021 - 1.05662874569329e-021i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
65 -3.4918408053655e-021 + 1.05662874569329e-021i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
66 -7.61773292876976e-021;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
67 4.84357724603939e-020 + 2.38824204294595e-019i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
68 4.84357724603939e-020 - 2.38824204294595e-019i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
69 -4.07088520945753e-020 - 2.31474543846105e-019i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
70 -4.07088520945753e-020 + 2.31474543846105e-019i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
71 8.73316588658882e-023;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
72 -5.21840635377469e-020;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
73 1.8461911504859e-023;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
74 5.20105247464461e-020;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
75 -4.68960092394415e-022;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
76 -1.44261407664171e-017 + 6.8922564526833e-019i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
77 -1.44261407664171e-017 - 6.8922564526833e-019i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
78 3.13688133935426e-022];
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
79
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
80 dPoles12 = [0.477546340377332 - 0.310830571032376i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
81 0.477546340377332 + 0.310830571032376i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
82 0.99790715414307 - 0.0028490561287024i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
83 0.99790715414307 + 0.0028490561287024i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
84 0.998014205354671 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
85 0.999585354543332 - 0.000780408757425194i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
86 0.999585354543332 + 0.000780408757425194i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
87 0.99966003029931 - 0.000830944038363768i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
88 0.99966003029931 + 0.000830944038363768i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
89 0.999962770401331 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
90 0.999981881865521 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
91 0.999999365763457 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
92 0.999981706320212 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
93 0.99992421574188 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
94 0.477898460791003 + 0.311001926610074i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
95 0.477898460791003 - 0.311001926610074i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
96 0];
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
97 dDTerms12 = 0;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
98
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
99 % Model Stefano Tf21
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
100 dRes21 = [-1.80035241582968e-016 + 1.99543917791863e-015i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
101 -1.80035241582968e-016 - 1.99543917791863e-015i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
102 -1.85590889333759e-013 - 1.23844418827409e-014i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
103 -1.85590889333759e-013 + 1.23844418827409e-014i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
104 5.03656596876842e-013 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
105 -2.62470963499904e-013 + 2.30024232938878e-012i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
106 -2.62470963499904e-013 - 2.30024232938878e-012i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
107 -9.83780507870955e-018 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
108 3.40426735130194e-021 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
109 9.78322351492755e-018 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
110 -1.65010934542937e-020 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
111 2.60918565203438e-015 + 1.0546609464659e-015i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
112 2.60918565203438e-015 - 1.0546609464659e-015i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
113 3.18105585405455e-014 + 2.48839990780042e-013i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
114 3.18105585405455e-014 - 2.48839990780042e-013i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
115 3.23021641947666e-013 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
116 4.81265000078114e-016 - 3.18269170053848e-017i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
117 4.81265000078114e-016 + 3.18269170053848e-017i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
118 5.16260024128201e-018];
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
119
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
120 dPoles21 = [0.872004077421604 - 0.110344282822693i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
121 0.872004077421604 + 0.110344282822693i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
122 0.956884129232757 - 0.0225532091775074i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
123 0.956884129232757 + 0.0225532091775074i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
124 0.966514825697177 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
125 0.995140550419744 - 0.000755127639524413i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
126 0.995140550419744 + 0.000755127639524413i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
127 0.999981802194393 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
128 0.99999936576546 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
129 0.999981722418555 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
130 0.999921882627659 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
131 0.999624431675213 + 0.000813407848742761i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
132 0.999624431675213 - 0.000813407848742761i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
133 0.997312006278751 + 0.00265611346834941i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
134 0.997312006278751 - 0.00265611346834941i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
135 0.990516544257531 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
136 0.477796923118318 + 0.311064085401834i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
137 0.477796923118318 - 0.311064085401834i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
138 0];
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
139
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
140 dDTerms21 = 0;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
141
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
142 % Model Stefano Tf22
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
143 dRes22 = [1.1284521501259e-014;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
144 -3.72133611555879e-014 - 2.08232683444075e-014i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
145 -3.72133611555879e-014 + 2.08232683444075e-014i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
146 9.84930639106637e-014 - 1.46640810672565e-013i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
147 9.84930639106637e-014 + 1.46640810672565e-013i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
148 2.51323684013671e-014 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
149 -5.64078525288305e-014 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
150 -1.62476406586366e-014 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
151 -1.08424815979566e-011 + 8.32328079357669e-012i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
152 -1.08424815979566e-011 - 8.32328079357669e-012i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
153 2.41831559776112e-011];
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
154
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
155 dPoles22 = [0.988511243978897;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
156 0.997305870640646 + 0.00211760900132725i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
157 0.997305870640646 - 0.00211760900132725i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
158 0.999626453270255 + 0.0008125673525946i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
159 0.999626453270255 - 0.0008125673525946i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
160 0.999999366366222 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
161 0.999981706320212 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
162 0.99992421574188 ;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
163 0.477898460791003 - 0.311001926610074i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
164 0.477898460791003 + 0.311001926610074i;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
165 0];
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
166
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
167 dDTerms22 = 0;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
168
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
169 % response calculation
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
170
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
171 pfparams.type = 'disc';
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
172 pfparams.freq = f;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
173 pfparams.fs = 10;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
174
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
175 % TF11
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
176 pfparams.res = dRes11;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
177 pfparams.pol = dPoles11;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
178 pfparams.dterm = dDTerms11;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
179 pfr = utils.math.pfresp(pfparams);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
180 mtf11 = pfr.resp;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
181
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
182 % TF12
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
183 pfparams.res = dRes12;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
184 pfparams.pol = dPoles12;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
185 pfparams.dterm = dDTerms12;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
186 pfr = utils.math.pfresp(pfparams);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
187 mtf12 = pfr.resp;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
188
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
189 % TF21
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
190 pfparams.res = dRes21;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
191 pfparams.pol = dPoles21;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
192 pfparams.dterm = dDTerms21;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
193 pfr = utils.math.pfresp(pfparams);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
194 mtf21 = pfr.resp;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
195
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
196 % TF22
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
197 pfparams.res = dRes22;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
198 pfparams.pol = dPoles22;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
199 pfparams.dterm = dDTerms22;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
200 pfr = utils.math.pfresp(pfparams);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
201 mtf22 = pfr.resp;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
202
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
203 % CSD calculation with Papoulis Method Style
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
204 % csd11 = mtf11.*conj(mtf11)+mtf12.*conj(mtf12);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
205 % csd12 = mtf11.*conj(mtf21)+mtf12.*conj(mtf22);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
206 % csd22 = mtf22.*conj(mtf22)+mtf21.*conj(mtf21);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
207 % csd21 = conj(csd12);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
208
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
209 % CSD = zeros(2,2,length(f));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
210 % for hh = 1:length(f)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
211 % CSD(:,:,hh) = [csd11(hh) csd12(hh);csd21(hh) csd22(hh)];
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
212 % end
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
213
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
214 % CSD = zeros(2,2,length(f));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
215 % for hh = 1:length(f)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
216 % CSD(:,:,hh) = [csd11(hh) 0;0 csd22(hh)];
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
217 % end
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
218
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
219
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
220 %% Get CSD Model
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
221
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
222 fs = 10;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
223 f = [logspace(-6,log10(fs/2),2000)]';
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
224
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
225 % currPath = cd;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
226 % cd MDC1 % assuming to start from ..\ltpda\software\m-toolbox\test
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
227 % [TF,CSD] = mdc1_tf_models(plist('f',f,'fs',fs));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
228 % cd(currPath)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
229
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
230 % get noise psd shapes
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
231 % parameters names
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
232 parnames = {'DH','S11','S1D','SD1','SDD','Dh1','Dh2','w1','w2',...
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
233 'TH','Th1','Th2'};
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
234
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
235 % nominal params values
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
236 parvalues = {1,1,0,-1e-4,1,1,1,-1.3e-6,-2.0e-6,0.35,0.25,0.28};
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
237
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
238 % Noise model
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
239 N = matrix(plist('built-in','mdc3_ifo2ifo_noise'));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
240
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
241 % evalueat model
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
242 [rw,cl]=size(N.objs);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
243 for ii=1:rw
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
244 for jj=1:cl
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
245 N.objs(ii,jj).setParams(parnames,parvalues);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
246 N.objs(ii,jj).setXvals(f);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
247 sp(ii,jj)=eval(N.objs(ii,jj));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
248 if ii==jj
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
249 sp(ii,jj)= abs(sp(ii,jj));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
250 end
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
251 end
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
252 end
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
253
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
254 csd11 = sp(1,1).y;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
255 csd12 = sp(1,2).y;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
256 csd21 = sp(2,1).y;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
257 csd22 = sp(2,2).y;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
258
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
259 CSD = zeros(2,2,length(f));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
260 for hh = 1:length(f)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
261 CSD(:,:,hh) = [csd11(hh) csd12(hh);csd21(hh) csd22(hh)];
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
262 end
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
263
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
264 %% Eigcsd 2-dim
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
265
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
266 h = utils.math.ndeigcsd(CSD,'MTD','PAP');
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
267
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
268 %% Test eigenshuffle
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
269
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
270 [l,m,npts] = size(CSD);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
271
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
272 % Finding suppression
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
273 k = min(sqrt(csd11./csd22));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
274 if k>=1
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
275 suppr = floor(k);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
276 else
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
277 n=0;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
278 while k<1
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
279 k=k*10;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
280 n=n+1;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
281 end
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
282 k = floor(k);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
283 suppr = k*10^(-n);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
284 end
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
285
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
286 supmat = [1 0;0 suppr];
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
287 % isupmat = [1 0;0 1/suppr];
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
288 isupmat = inv(supmat);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
289
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
290 PP = CSD;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
291 for ii=1:npts
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
292 PP(:,:,ii)=supmat*CSD(:,:,npts)*supmat;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
293 end
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
294
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
295 [V,D] = eigenshuffle(CSD);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
296 h = ones(l,m,npts);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
297
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
298 for ii=1:npts
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
299 HH = isupmat*V(:,:,ii)*sqrt(diag(D(:,ii)));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
300 h(:,:,ii) = HH;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
301 end
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
302
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
303 %% Plot CSD
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
304
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
305 % CSD calculation
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
306 ncsd11 = squeeze(h(1,1,:).*conj(h(1,1,:))+h(1,2,:).*conj(h(1,2,:)));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
307 ncsd12 = squeeze(h(1,1,:).*conj(h(2,1,:))+h(1,2,:).*conj(h(2,2,:)));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
308 ncsd22 = squeeze(h(2,2,:).*conj(h(2,2,:))+h(2,1,:).*conj(h(2,1,:)));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
309 ncsd21 = conj(csd12);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
310
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
311 figure()
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
312 loglog(f,abs(csd11),'k')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
313 hold on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
314 grid on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
315 loglog(f,abs(ncsd11),'r')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
316
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
317 figure()
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
318 loglog(f,abs(csd12),'k')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
319 hold on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
320 grid on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
321 loglog(f,abs(ncsd12),'r')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
322
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
323 figure()
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
324 semilogx(f,angle(csd12),'k')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
325 hold on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
326 grid on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
327 semilogx(f,angle(ncsd12),'r')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
328
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
329 figure()
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
330 loglog(f,abs(csd22),'k')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
331 hold on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
332 grid on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
333 loglog(f,abs(ncsd22),'r')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
334
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
335 %% plot TFs
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
336
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
337 figure()
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
338 loglog(f,abs(mtf11),'k')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
339 hold on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
340 grid on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
341 loglog(f,abs(squeeze(h(1,1,:))),'r')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
342
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
343 figure()
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
344 loglog(f,abs(mtf12),'k')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
345 hold on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
346 grid on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
347 loglog(f,abs(squeeze(h(1,2,:))),'r')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
348
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
349 figure()
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
350 loglog(f,abs(mtf21),'k')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
351 hold on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
352 grid on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
353 loglog(f,abs(squeeze(h(2,1,:))),'r')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
354
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
355 figure()
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
356 loglog(f,abs(mtf22),'k')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
357 hold on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
358 grid on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
359 loglog(f,abs(squeeze(h(2,2,:))),'r')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
360
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
361 %% plot TFs phase
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
362
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
363 figure()
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
364 semilogx(f,angle(mtf11),'k')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
365 hold on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
366 grid on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
367 semilogx(f,angle(squeeze(h(1,1,:))),'r')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
368
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
369 figure()
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
370 semilogx(f,angle(mtf12),'k')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
371 hold on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
372 grid on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
373 semilogx(f,angle(squeeze(h(1,2,:))),'r')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
374
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
375 figure()
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
376 semilogx(f,angle(mtf21),'k')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
377 hold on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
378 grid on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
379 semilogx(f,angle(squeeze(h(2,1,:))),'r')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
380
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
381 figure()
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
382 semilogx(f,angle(mtf22),'k')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
383 hold on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
384 grid on
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
385 semilogx(f,angle(squeeze(h(2,2,:))),'r')
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
386
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
387
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
388
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
389 %% 3 Dim Test
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
390
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
391
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
392 % sampling frequency
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
393 fs = 1; % Hz
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
394 % frequency vector
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
395 f = logspace(-6,log10(0.5),300);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
396 f = f.';
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
397
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
398 %%% Set the 3 channel model
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
399
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
400 % d = [1 -1.1 0.5 0.02];
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
401 %
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
402 % h11 = miir(1e-1.*[1 -0.5],d,fs);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
403 % h12 = miir(1e-4.*[0 -0.5],d,fs);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
404 % h13 = miir(1e-2.*[1 0.2],d,fs);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
405 %
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
406 % h21 = miir(1e-3.*[0 0.4],d,fs);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
407 % h22 = miir(1e-5.*[1 -0.6],d,fs);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
408 % h23 = miir(1e-6.*[1 0.3],d,fs);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
409 %
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
410 % h31 = miir(1e-4.*[0 0.1 -0.2],d,fs);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
411 % h32 = miir(1e-7.*[0 -0.6],d,fs);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
412 % h33 = miir(1e-5.*[1 0.05 0.3],d,fs);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
413
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
414
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
415 h11 = pzmodel(plist('GAIN', [0.01], 'POLES', [pz(9.9999999999999995e-007,NaN) pz(0.01,NaN) pz(0.02,NaN)], 'ZEROS', [pz(0.0001,NaN) pz(0.002,NaN)]));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
416 h12 = pzmodel(plist('GAIN', [0.0001], 'POLES', [pz(0.10000000000000001,NaN) pz(0.0030000000000000001,NaN) pz(0.0050000000000000001,NaN) pz(1.0000000000000001e-005,NaN)], 'ZEROS', [pz(0.00050000000000000001,NaN) pz(0.001,NaN)]));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
417 h13 = pzmodel(plist('GAIN', [0.00001], 'POLES', [pz(5.0000000000000004e-006,NaN) pz(0.02,NaN)], 'ZEROS', pz(0.00040000000000000002,NaN)));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
418
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
419 h21 = pzmodel(plist('GAIN', [0.00001], 'POLES', [pz(5.0000000000000004e-006,NaN) pz(0.10000000000000001,NaN) pz(0.050000000000000003,NaN)], 'ZEROS', [pz(0.002,NaN) pz(0.00020000000000000001,NaN)]));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
420 h22 = pzmodel(plist('GAIN', [0.001], 'POLES', [pz(0.01,NaN) pz(9.9999999999999995e-008,NaN) pz(0.002,NaN) pz(0.001,NaN)], 'ZEROS', [pz(1.0000000000000001e-005,NaN) pz(2.0000000000000002e-005,NaN) pz(5.0000000000000002e-005,NaN) pz(0.20000000000000001,NaN)]));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
421 h23 = pzmodel(plist('GAIN', [0.0001], 'POLES', [pz(0.01,NaN) pz(0.029999999999999999,NaN) pz(0.10000000000000001,NaN) pz(5.0000000000000002e-005,NaN)], 'ZEROS', [pz(0.0001,NaN) pz(0.0050000000000000001,NaN)]));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
422
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
423 h31 = pzmodel(plist('GAIN', [0.001], 'POLES', [pz(0.01,NaN) pz(0.029999999999999999,NaN) pz(0.10000000000000001,NaN) pz(0.00050000000000000001,NaN)], 'ZEROS', [pz(0.0001,NaN) pz(0.0050000000000000001,NaN)]));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
424 h32 = pzmodel(plist('GAIN', [0.00001], 'POLES', [pz(0.01,NaN) pz(0.029999999999999999,NaN) pz(0.10000000000000001,NaN) pz(0.00050000000000000001,NaN) pz(0.00029999999999999997,NaN)], 'ZEROS', [pz(0.002,NaN) pz(0.059999999999999998,NaN)]));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
425 h33 = pzmodel(plist('GAIN', [0.05], 'POLES', [pz(0.0001,NaN) pz(9.9999999999999995e-008,NaN) pz(0.002,NaN)], 'ZEROS', [pz(3.0000000000000001e-005,NaN) pz(5.0000000000000002e-005,NaN) pz(0.10000000000000001,NaN)]));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
426
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
427 %%% TFs resps
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
428
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
429 % filters response
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
430 plresp = plist('f',f);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
431
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
432 rh11 = resp(h11,plresp);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
433 rh12 = resp(h12,plresp);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
434 rh13 = resp(h13,plresp);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
435
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
436 rh21 = resp(h21,plresp);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
437 rh22 = resp(h22,plresp);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
438 rh23 = resp(h23,plresp);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
439
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
440 rh31 = resp(h31,plresp);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
441 rh32 = resp(h32,plresp);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
442 rh33 = resp(h33,plresp);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
443
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
444 % iplot(rh11,rh12,rh13)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
445 % iplot(rh21,rh22,rh23)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
446 % iplot(rh31,rh32,rh33)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
447
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
448 % iplot(rh11,rh22,rh33)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
449
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
450 %%% Spectra calculation with Kay style method
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
451
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
452 % Note that the definition of cross-spectral matrix follows that of s M
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
453 % Kay, Modern Spectral Estimation
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
454 fs = 1; % Hz
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
455 % csd11
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
456 G11 = rh11.y.*conj(rh11.y) + rh12.y.*conj(rh12.y) + rh13.y.*conj(rh13.y);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
457 G11 = ao(plist('xvals', f, 'yvals', G11, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1')));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
458 G11.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
459
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
460 % csd12
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
461 G12 = conj(rh11.y).*rh21.y + conj(rh12.y).*rh22.y + conj(rh13.y).*rh23.y;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
462 G12 = ao(plist('xvals', f, 'yvals', G12, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1')));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
463 G12.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
464
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
465 % csd13
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
466 G13 = conj(rh11.y).*rh31.y + conj(rh12.y).*rh32.y + conj(rh13.y).*rh33.y;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
467 G13 = ao(plist('xvals', f, 'yvals', G13, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1')));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
468 G13.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
469
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
470 % csd21
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
471 G21 = conj(G12);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
472 G21.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
473
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
474 % csd22
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
475 G22 = rh21.y.*conj(rh21.y) + rh22.y.*conj(rh22.y) + rh23.y.*conj(rh23.y);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
476 G22 = ao(plist('xvals', f, 'yvals', G22, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1')));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
477 G22.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
478
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
479 % csd23
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
480 G23 = conj(rh21.y).*rh31.y + conj(rh22.y).*rh32.y + conj(rh23.y).*rh33.y;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
481 G23 = ao(plist('xvals', f, 'yvals', G23, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1')));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
482 G23.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
483
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
484 % csd31
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
485 G31 = conj(G13);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
486 G31.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
487
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
488 % csd32
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
489 G32 = conj(G23);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
490 G32.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
491
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
492 % csd33
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
493 G33 = conj(rh31.y).*rh31.y + conj(rh32.y).*rh32.y + conj(rh33.y).*rh33.y;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
494 G33 = ao(plist('xvals', f, 'yvals', G33, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1')));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
495 G33.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
496
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
497 % iplot(G11,G12,G13,G22,G23,G33)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
498
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
499 % iplot(G11,G22,G33)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
500
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
501 %%% Build CSD matrix
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
502
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
503 CSD = zeros(3,3,length(f));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
504 for hh = 1:length(f)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
505 CSD(:,:,hh) = [G11.y(hh) G12.y(hh) G13.y(hh);G21.y(hh) G22.y(hh) G23.y(hh);G31.y(hh) G32.y(hh) G33.y(hh)];
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
506 end
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
507
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
508 % for hh = 1:length(f)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
509 % CSD(:,:,hh) = [G11.y(hh) 0 0;0 G22.y(hh) 0;0 0 G33.y(hh)];
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
510 % end
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
511
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
512 %%% 3-dim Eigcsd with kay method stule
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
513 h = utils.math.ndeigcsd(CSD,'MTD','KAY');
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
514
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
515 % %%% Eigenshuffle
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
516 %
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
517 % [l,m,kk] = size(CSD);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
518 % suppr = ones(l,l);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
519 % for ii = 2:l
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
520 % k = ones(l,1);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
521 %
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
522 % for jj = ii-1:-1:1
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
523 % k(jj) = min(sqrt(CSD(jj,jj,:)./CSD(ii,ii,:)));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
524 % if k(jj)>=1
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
525 % suppr(jj,ii) = floor(k(jj));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
526 % else
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
527 % n=0;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
528 % while k(jj)<1
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
529 % k(jj)=k(jj)*10;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
530 % n=n+1;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
531 % end
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
532 % k(jj) = floor(k(jj));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
533 % suppr(jj,ii) = k(jj)*10^(-n);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
534 % % suppr(ii) = suppr(ii)*suppr(ii-1);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
535 % end
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
536 % end
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
537 % % csuppr(ii) = prod(suppr(:,ii));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
538 % end
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
539 % csuppr = prod(suppr,2);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
540 % supmat = diag(csuppr);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
541 % % ssup = rot90(rot90(supmat));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
542 % ssup = supmat;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
543 % % ssup = eye(3);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
544 % % ssup = supmat*supmat.';
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
545 % issup = inv(ssup);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
546 %
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
547 % for jj = 1:kk
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
548 % nCSD(:,:,jj) = ssup*CSD(:,:,jj)*ssup;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
549 % end
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
550 %
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
551 % [Vseq,Dseq] = eigenshuffle(nCSD);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
552 %
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
553 % % get h = V*sqrt(D);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
554 % [nn,mm,kk] = size(Vseq);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
555 % h = zeros(nn,mm,kk);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
556 % for dd = 1:kk
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
557 % % h(:,:,dd) = rot90(issup*Vseq(:,:,dd)*sqrt(diag(Dseq(:,dd))));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
558 % h(:,:,dd) = issup*Vseq(:,:,dd)*sqrt(diag(Dseq(:,dd)));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
559 % end
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
560
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
561 %% Build TF AOs
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
562
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
563 nh11 = ao(plist('xvals', f, 'yvals', squeeze(h(1,1,:)), 'fs', fs, 'dtype', 'fsdata'));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
564 nh11.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
565
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
566 nh12 = ao(plist('xvals', f, 'yvals', squeeze(h(1,2,:)), 'fs', fs, 'dtype', 'fsdata'));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
567 nh12.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
568
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
569 nh13 = ao(plist('xvals', f, 'yvals', squeeze(h(1,3,:)), 'fs', fs, 'dtype', 'fsdata'));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
570 nh13.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
571
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
572 nh21 = ao(plist('xvals', f, 'yvals', squeeze(h(2,1,:)), 'fs', fs, 'dtype', 'fsdata'));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
573 nh21.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
574
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
575 nh22 = ao(plist('xvals', f, 'yvals', squeeze(h(2,2,:)), 'fs', fs, 'dtype', 'fsdata'));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
576 nh22.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
577
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
578 nh23 = ao(plist('xvals', f, 'yvals', squeeze(h(2,3,:)), 'fs', fs, 'dtype', 'fsdata'));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
579 nh23.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
580
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
581 nh31 = ao(plist('xvals', f, 'yvals', squeeze(h(3,1,:)), 'fs', fs, 'dtype', 'fsdata'));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
582 nh31.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
583
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
584 nh32 = ao(plist('xvals', f, 'yvals', squeeze(h(3,2,:)), 'fs', fs, 'dtype', 'fsdata'));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
585 nh32.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
586
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
587 nh33 = ao(plist('xvals', f, 'yvals', squeeze(h(3,3,:)), 'fs', fs, 'dtype', 'fsdata'));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
588 nh33.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
589
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
590 iplot(rh11,nh11)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
591 iplot(rh12,nh12)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
592 iplot(rh13,nh13)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
593 iplot(rh21,nh21)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
594 iplot(rh22,nh22)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
595 iplot(rh23,nh23)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
596 iplot(rh31,nh31)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
597 iplot(rh32,nh32)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
598 iplot(rh33,nh33)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
599
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
600 %% Build CSD AOs
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
601
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
602 % Note that the definition of cross-spectral matrix follows that of s M
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
603 % Kay, Modern Spectral Estimation
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
604
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
605 % csd11
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
606 H11 = nh11.y.*conj(nh11.y) + nh12.y.*conj(nh12.y) + nh13.y.*conj(nh13.y);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
607 H11 = ao(plist('xvals', f, 'yvals', H11, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1')));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
608 H11.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
609
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
610 % csd12
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
611 H12 = conj(nh11.y).*nh21.y + conj(nh12.y).*nh22.y + conj(nh13.y).*nh23.y;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
612 H12 = ao(plist('xvals', f, 'yvals', H12, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1')));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
613 H12.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
614
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
615 % csd13
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
616 H13 = conj(nh11.y).*nh31.y + conj(nh12.y).*nh32.y + conj(nh13.y).*nh33.y;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
617 H13 = ao(plist('xvals', f, 'yvals', H13, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1')));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
618 H13.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
619
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
620 % csd21
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
621 H21 = conj(H12);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
622 H21.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
623
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
624 % csd22
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
625 H22 = nh21.y.*conj(nh21.y) + nh22.y.*conj(nh22.y) + nh23.y.*conj(nh23.y);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
626 H22 = ao(plist('xvals', f, 'yvals', H22, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1')));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
627 H22.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
628
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
629 % csd23
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
630 H23 = conj(nh21.y).*nh31.y + conj(nh22.y).*nh32.y + conj(nh23.y).*nh33.y;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
631 H23 = ao(plist('xvals', f, 'yvals', H23, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1')));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
632 H23.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
633
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
634 % csd31
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
635 H31 = conj(H13);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
636 H31.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
637
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
638 % csd32
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
639 H32 = conj(H23);
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
640 H32.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
641
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
642 % csd33
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
643 H33 = conj(nh31.y).*nh31.y + conj(nh32.y).*nh32.y + conj(nh33.y).*nh33.y;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
644 H33 = ao(plist('xvals', f, 'yvals', H33, 'fs', fs, 'dtype', 'fsdata', 'yunits', unit('m^2').*unit('Hz^-1')));
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
645 H33.setName;
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
646
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
647 iplot(G11,H11)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
648 iplot(G12,H12)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
649 iplot(G13,H13)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
650 iplot(G22,H22)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
651 iplot(G23,H23)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
652 iplot(G33,H33)
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
653
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
654
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
655
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
656
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
657
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
658
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
659
|
Daniele Nicolodi <nicolodi@science.unitn.it>
parents:
diff
changeset
|
660
|