Mercurial > hg > ltpda

comparison m-toolbox/test/test_mdc1_conv2acc_utn.m @ 0:f0afece42f48

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Import.
author Daniele Nicolodi <nicolodi@science.unitn.it>
date Wed, 23 Nov 2011 19:22:13 +0100
parents
children
comparison
equal deleted inserted replaced
-1:000000000000 0:f0afece42f48
1 % MDC1 Conversion to acceleration
2 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3 % Test script for the MDC1 conversion to acceleration of IFO data
4 %
5 % $Id: test_mdc1_conv2acc_utn.m,v 1.2 2009/03/02 15:53:46 luigi Exp $
6 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
7 %% Clear
8
9 clear all
10
11 %% Loading data
12
13 currPath = cd;
14 cd 'C:\Dati\mock_data\interferometer\';
15
16 j=1;
17 fs = 10;
18
19 disp(sprintf('*** Reading data %d', j));
20 in = load(sprintf('mockdata_16_48_17_11_2007_%d.dat', j));
21 d1 = in(:,2);
22 d12 = in(:,3);
23
24 % Build O1
25 eval('ts = tsdata(d1, fs);');
26 eval(sprintf('o1_%d = ao(ts);', j));
27 eval(sprintf('o1_%d = setName(o1_%d, ''o1'');', j, j));
28 eval(sprintf('o1_%d = setXunits(o1_%d, ''s'');', j, j));
29 eval(sprintf('o1_%d = setYunits(o1_%d, ''m'');', j, j));
30
31 % Build O12
32 eval('ts = tsdata(d12, fs);');
33 eval(sprintf('o12_%d = ao(ts);', j));
34 eval(sprintf('o12_%d = setName(o12_%d, ''o12'');', j, j));
35 eval(sprintf('o12_%d = setXunits(o12_%d, ''s'');', j, j));
36 eval(sprintf('o12_%d = setYunits(o12_%d, ''m'');', j, j));
37
38 global LTPDAinvar
39
40 LTPDAinvar = { o1_1 , 1 , 'From file';...
41 o12_1 , 1 , 'From file'};
42
43 cd(currPath)
44
45
46 %% =============================================================
47
48 channel1 = LTPDAinvar{1,1};
49 channel2 = LTPDAinvar{2,1};
50
51 %%
52
53 CH1 = channel1.lpsd;
54 CH2 = channel2.lpsd;
55 iplot(CH1)
56 iplot(CH2)
57
58 %% DRAG FREE
59
60 rat1pl = plist('NUM', [0.0004659;0.1349;4.37;0.8304;0.07449;0.002978;4.403e-005], 'DEN', [1;5.046;9.609;11.05;0.01221;3.401e-006;0]);
61 control1 = mdc1_ifo2cont_utn(channel1,rat1pl);
62
63 %% retarded actuaction - thrusters:
64
65 retard1pl = plist('TAU', 0.1, 'DELTAT', 0.315, 'FS', 10, 'TOL', 1e-007, 'AMP', 1);
66 actuation1 = mdc1_cont2act_utn(control1,retard1pl);
67
68 %% LOW FREQ SUSP
69
70 rat2pl = plist('NUM', [-2.726e-007;1.665e-005;1.303e-007;8.381e-010], 'DEN', [1;0.2189;0.01922;0.0007803;0]);
71 control2 = mdc1_ifo2cont_utn(channel2,rat2pl);
72
73 %% retarded actuaction - electrical actuators:
74
75 retard2pl = plist('TAU', 0.01, 'DELTAT', 0.305, 'FS', 10, 'TOL', 1e-007, 'AMP', 1);
76 actuation2 = mdc1_cont2act_utn(control2,retard2pl);
77
78 %% Free dynamics:
79
80 plfd = plist('pstiff1', -1.3e-6,'pstiff2', -2e-6,'cross_talk', -1e-4,'METHOD', 'PARFIT');
81 freeDynamics = mdc1_ifo2acc_fd_utn(channel1,channel2,plfd);
82
83 %% Output
84
85 output1 = freeDynamics.index(1) + actuation1;
86 output2 = freeDynamics.index(2) + actuation2;
87 % output1 = freeDynamics(1) + control1;
88 % output2 = freeDynamics(2) + control2;
89
90 %% ******************************************
91 % plotting output psd
92
93 psd1 = output1.lpsd;
94 psd2 = output2.lpsd;
95 iplot(psd1)
96 iplot(psd2)
97
98 %% ====================== More checks =====================================
99
100 %% Extracting filters
101
102 miir1Obj = find(control1.procinfo,'CONT_FILTER');
103 retard1 = find(actuation1.procinfo,'ACT_FILTER');
104 miir2Obj = find(control2.procinfo,'CONT_FILTER');
105 retard2 = find(actuation2.procinfo,'ACT_FILTER');
106
107 %% Checking controllers
108
109 % frequencies
110 f = logspace(-6,log10(5),300);
111 % f = logspace(-6,0,300);
112 f = f.';
113
114 % response of continuous drag free
115 Rcdf = resp(rat1Obj,f);
116 Rcdf.setName('Continuous DF');
117
118 % response of discrete drag free
119 tRddf = resp(miir1Obj,plist('f',f));
120 Rddf = tRddf(1);
121 for ii = 2:numel(tRddf)
122 Rddf = Rddf + tRddf(ii);
123 end
124 Rddf.setName('Discrete DF');
125
126 % response of continuous low freq susp
127 Rclfs = resp(rat2Obj,f);
128 Rclfs.setName('Continuous LFS');
129
130 % response of discrete low freq susp
131 tRdlfs = resp(miir2Obj,plist('f',f));
132 Rdlfs = tRdlfs(1);
133 for ii = 2:numel(tRdlfs)
134 Rdlfs = Rdlfs + tRdlfs(ii);
135 end
136 Rdlfs.setName('Discrete LFS');
137
138 % response of actuation filter 1
139 Ract1 = resp(retard1,plist('f',f));
140 Ract1.setName('Thrusters');
141
142 % response of actuation filter 2
143 Ract2 = resp(retard2,plist('f',f));
144 Ract1.setName('Actuators');
145
146 %% Plotting reposnse
147
148 % Continuous and discrete controllers
149 pl = plist('Legends', {'cH_{df}','dH_{df}','cH_{lfs}','dH_{lfs}'},...
150 'LineStyles', {'', '--', '', '--'});
151 iplot(Rcdf,Rddf,Rclfs,Rdlfs,pl)
152
153 % Ratio between discrete and continuous
154 pl = plist('Legends', {'dH_{df}/cH_{df}','dH_{lfs}/cH_{lfs}'},...
155 'YScales',{'All','lin'});
156 iplot(Rddf./Rddf,Rdlfs./Rclfs,pl)
157
158 % Actuators
159 pl = plist('Legends', {'Thrusters','El. Act.'},...
160 'YScales',{'All','lin'});
161 iplot(Ract1,Ract2,pl)
162
163 % discrete controllers + actuation
164 pl = plist('Legends', {'dH_{df}','dH_{df} + Tht','dH_{lfs}','dH_{lfs} + Act'},...
165 'LineStyles', {'', '--', '', '--'});
166 iplot(Rddf,Rddf.*Ract1,Rdlfs,Rdlfs.*Ract2,pl)
167
168 %% Calculating transfer functions
169
170 % controller ***
171 pl = plist('Nfft', fs*1000);
172 tfdf = tfe(channel1,control1,pl);
173 tflfs = tfe(channel2,control2,pl);
174 plpl = plist('Legends', {'H_{df}','H_{df} data','H_{lfs}','H_{lfs} data'},...
175 'LineStyles', {'', '--', '', '--'});
176 iplot(Rddf,tfdf(1,2),Rdlfs,tflfs(1,2),plpl)
177
178 % Actuation ***
179 pl = plist('Nfft', fs*1000);
180 tfact1 = tfe(control1,actuation1,pl);
181 tfact2 = tfe(control2,actuation2,pl);
182 plpl = plist('Legends', {'Tht','Tht data','El Act','El Act data'},...
183 'LineStyles', {'', '--', '', '--'});
184 iplot(Ract1,tfact1(1,2),Ract2,tfact2(1,2),plpl)
185
186 % Free Dynamics TFs ***
187 pl = plist('Nfft', fs*10000);
188 tffd1 = tfe(channel1,freeDynamics(1),pl);
189 tffd2 = tfe(channel2,freeDynamics(2),pl);
190
191 % Theorethical
192 % w1 = -1.3e-6;
193 % w2 = -2e-6;
194 % D1 = -1e-4;
195 % s = 2.*pi.*1i.*f;
196 % FD1 = s.^2 + w1;
197 % FD2 = -D1.*s.^2 + w2 - w1 - D1*w2 + s.^2 + w2;
198 FD1 = ao(plist('fsfcn', '(2.*pi.*1i.*f).^2 + -1.3e-6', 'f', f));
199 FD2 = ao(plist('fsfcn', '-(-1e-4).*(2.*pi.*1i.*f).^2 + (-2e-6) - (-1.3e-6) - (-1e-4)*(-2e-6) + (2.*pi.*1i.*f).^2 + (-2e-6)', 'f', f));
200
201 plpl = plist('Legends', {'FD1','FD1 data','FD2','FD2 data'},...
202 'LineStyles', {'', '--', '', '--'});
203 iplot(abs(FD1),abs(tffd1(1,2)),abs(FD2),abs(tffd2(1,2)),plpl)
204