| 181 | | % %hlavni iteracni smycka |
| 182 | | % for iterace = 1:It, |
| 183 | | % %generovani stavu |
| 184 | | % for n = 1:N, |
| 185 | | % xn(:, 1, n) = x0' + neznalost*sqrtm(P)*randn(4,1); |
| 186 | | % for k = 1:Kt-1, |
| 187 | | % tu = L(:, :, k)*(xn(:, k, n) - [0 0 OMEGAt 0]'); |
| 188 | | % xn(1, k+1, n) = a*xn(1, k, n) + b*xn(3, k, n)*sin(xn(4, k, n)) + c*tu(1) + sumsim*sqrt(Q(1, 1))*randn(); |
| 189 | | % xn(2, k+1, n) = a*xn(2, k, n) - b*xn(3, k, n)*cos(xn(4, k, n)) + c*tu(2) + sumsim*sqrt(Q(2, 2))*randn(); |
| 190 | | % xn(3, k+1, n) = d*xn(3, k, n) + e*(xn(2, k, n)*cos(xn(4, k, n)) - xn(1, k, n)*sin(xn(4, k, n))) + sumsim*sqrt(Q(3, 3))*randn(); |
| 191 | | % xn(4, k+1, n) = xn(4, k, n) + xn(3, k, n)*DELTAt + sumsim*sqrt(Q(4, 4))*randn(); |
| 192 | | % % xn(1, k+1, n) = a*xn(1, k, n) + b*(xn(3, k, n)+OMEGAt)*sin(xn(4, k, n)) + c*tu(1) + sum*sqrt(Q(1, 1))*randn(); |
| 193 | | % % xn(2, k+1, n) = a*xn(2, k, n) - b*(xn(3, k, n)+OMEGAt)*cos(xn(4, k, n)) + c*tu(2) + sum*sqrt(Q(2, 2))*randn(); |
| 194 | | % % xn(3, k+1, n) = -OMEGAt + d*(xn(3, k, n)+OMEGAt) + e*(xn(2, k, n)*cos(xn(4, k, n)) - xn(1, k, n)*sin(xn(4, k, n))) + sum*sqrt(Q(3, 3))*randn(); |
| 195 | | % % xn(4, k+1, n) = xn(4, k, n) + (xn(3, k, n)+OMEGAt)*DELTAt + sum*sqrt(Q(4, 4))*randn(); |
| 196 | | % end |
| 197 | | % end |
| 198 | | % |
| 199 | | % %prumerny stav |
| 200 | | % xs = mean(xn, 3); |
| 201 | | % |
| 202 | | % %napocteni S a L |
| 203 | | % for kt = 1:Kt-1, |
| 204 | | % % receding horizon |
| 205 | | % if((K-1+kt)<Kt) |
| 206 | | % S(:, :, K) = X; |
| 207 | | % for k = K-1+kt-1:-1:1+kt-1, |
| 208 | | % A(3, 1) = -e*sin(xs(4, k)); |
| 209 | | % A(3, 2) = e*cos(xs(4, k)); |
| 210 | | % A(1, 3) = b*sin(xs(4, k)); |
| 211 | | % A(2, 3) = -b*cos(xs(4, k)); |
| 212 | | % A(1, 4) = b*(xs(3, k))*cos(xs(4, k)); |
| 213 | | % A(2, 4) = b*(xs(3, k))*sin(xs(4, k)); |
| 214 | | % % A(1, 4) = b*(xs(3, k)+OMEGAt)*cos(xs(4, k)); |
| 215 | | % % A(2, 4) = b*(xs(3, k)+OMEGAt)*sin(xs(4, k)); |
| 216 | | % A(3, 4) = -e*(xs(2, k)*sin(xs(4, k) + xs(1,k)*cos(xs(4, k)))); |
| 217 | | % S(:, :, k-kt+1) = A'*(S(:, :, k-kt+2) - S(:, :, k-kt+2)*B*inv(B'*S(:, :, k-kt+2)*B + Y)*B'*S(:, :, k-kt+2))*A + X; |
| 218 | | % end |
| 219 | | % L(:, :, kt) = -inv(B'*S(:, :, 1)*B + Y)*B'*S(:, :, 1)*A; |
| 220 | | % else |
| 221 | | % L(:, :, kt) = L(:, :, kt-1); %kopiruje poslednich K kroku z Kt kde to nejde na K predpocitat |
| 222 | | % end |
| 223 | | % |
| 224 | | % end |
| 225 | | % end |
| 226 | | % toc |
| 227 | | % |
| 228 | | % %vysledky |
| 229 | | % clf |
| 230 | | % subplot(2, 3, 3); |
| 231 | | % hold all |
| 232 | | % plot(1:Kt, OMEGAt*ones(1,Kt)); |
| 233 | | % % L(:,:,1) |
| 234 | | % for n = 1:N, |
| 235 | | % xn(:, 1, n) = x0' + neznalost*sqrtm(P)*randn(4,1); |
| 236 | | % for k = 1:Kt-1, |
| 237 | | % tu = L(:, :, 1)*(xn(:, k, n) - [0 0 OMEGAt 0]'); |
| 238 | | % xn(1, k+1, n) = a*xn(1, k, n) + b*xn(3, k, n)*sin(xn(4, k, n)) + c*tu(1) + sumsim*sqrt(Q(1, 1))*randn(); |
| 239 | | % xn(2, k+1, n) = a*xn(2, k, n) - b*xn(3, k, n)*cos(xn(4, k, n)) + c*tu(2) + sumsim*sqrt(Q(2, 2))*randn(); |
| 240 | | % xn(3, k+1, n) = d*xn(3, k, n) + e*(xn(2, k, n)*cos(xn(4, k, n)) - xn(1, k, n)*sin(xn(4, k, n))) + sumsim*sqrt(Q(3, 3))*randn(); |
| 241 | | % xn(4, k+1, n) = xn(4, k, n) + xn(3, k, n)*DELTAt + sumsim*sqrt(Q(4, 4))*randn(); |
| 242 | | % % xn(1, k+1, n) = a*xn(1, k, n) + b*(xn(3, k, n)+OMEGAt)*sin(xn(4, k, n)) + c*tu(1) + sum*sqrt(Q(1, 1))*randn(); |
| 243 | | % % xn(2, k+1, n) = a*xn(2, k, n) - b*(xn(3, k, n)+OMEGAt)*cos(xn(4, k, n)) + c*tu(2) + sum*sqrt(Q(2, 2))*randn(); |
| 244 | | % % xn(3, k+1, n) = -OMEGAt + d*(xn(3, k, n)+OMEGAt) + e*(xn(2, k, n)*cos(xn(4, k, n)) - xn(1, k, n)*sin(xn(4, k, n))) + sum*sqrt(Q(3, 3))*randn(); |
| 245 | | % % xn(4, k+1, n) = xn(4, k, n) + (xn(3, k, n)+OMEGAt)*DELTAt + sum*sqrt(Q(4, 4))*randn(); |
| 246 | | % |
| 247 | | % u(:, k) = tu; |
| 248 | | % end |
| 249 | | % |
| 250 | | % % xn(3, :, n) = xn(3, :, n) + OMEGAt*ones(1, Kt); |
| 251 | | % |
| 252 | | % subplot(2, 3, 1); |
| 253 | | % hold all |
| 254 | | % plot(1:Kt, xn(1, :, n)); |
| 255 | | % title('i_{\alpha}'); |
| 256 | | % subplot(2, 3, 2); |
| 257 | | % hold all |
| 258 | | % plot(1:Kt, xn(2, :, n)); |
| 259 | | % title('i_{\beta}'); |
| 260 | | % subplot(2, 3, 3); |
| 261 | | % hold all |
| 262 | | % plot(1:Kt, xn(3, :, n)); |
| 263 | | % title('\omega'); |
| 264 | | % subplot(2, 3, 4); |
| 265 | | % hold all |
| 266 | | % plot(1:Kt, xn(4, :, n)); |
| 267 | | % title('\theta'); |
| 268 | | % subplot(2, 3, 5); |
| 269 | | % hold all |
| 270 | | % plot(1:Kt, u(1, :)); |
| 271 | | % title('u_{\alpha}'); |
| 272 | | % subplot(2, 3, 6); |
| 273 | | % hold all |
| 274 | | % plot(1:Kt, u(2, :)); |
| 275 | | % title('u_{\beta}'); |
| 276 | | % end |
| | 186 | if(errnans > 0) |
| | 187 | lstore |
| | 188 | disp('L is NaN') |
| | 189 | errnans |
| | 190 | end |
| | 191 | |
| | 192 | |