卫星姿态基于matlab卫星姿态控制系统【含Matlab源码 2751期】.zip

clear,clc, %% 转动惯量 global Ix Iy Iz If Ix=264;Iy=264;Iz=28;If=0.5; %% 轨道高度与质心系转速 global miu Ro Wo miu = 398600.44;Ro = 7000;Wo = sqrt(miu/Ro)/Ro; %% 期望姿态 global r2d Etheta1 Etheta2 Etheta3 r2d=57.3; Etheta1=50/r2d; Etheta2=10/r2d; Etheta3=-10/r2d; %% PID参数 global kp1 kd1 kp2 kd2 kp3 kd3 kp1=4000; kd1=4000; kp2=9500; kd2=9500; kp3=8000; kd3=8000; %% 初值 theta1=10/r2d; theta2=20/r2d; theta3=-30/r2d; w1=10/r2d; w2=20/r2d; w3=-30/r2d; stateInit0=[theta1;theta2;theta3;w1;w2;w3]; stateInit=[theta1;theta2;theta3;w1;w2;w3;0;0;0]; tspan = [0.01:0.01:20]; [m,n] = size(tspan); %% 无控制 [t,stateout]=ode45('staoutcontl',tspan,stateInit0); stateout = transout(stateout,n); th1out=stateout(:,1); th2out=stateout(:,2); th3out=stateout(:,3); w1out=stateout(:,4); w2out=stateout(:,5); w3out=stateout(:,6); figure(1),subplot(2,1,1), plot(t,th1out*r2d,'r:',t,th2out*r2d,'k--',t,th3out*r2d,'linewidth',2), xlabel('time(s)'),ylabel('姿态角(deg)'),title('质心坐标系'); legend('\theta_1','\theta_2','\theta_3'),grid on subplot(2,1,2), plot(t,w1out*r2d,t,w2out*r2d,t,w3out*r2d,'linewidth',2), legend('\omega_1','\omega_2','\omega_3'),grid on xlabel('time(s)'),ylabel('角速率(deg/s)') %% PID控制 [t,stateout]=ode45('stacontrol_feilun',tspan,stateInit); th1out=stateout(:,1); th2out=stateout(:,2); th3out=stateout(:,3); w1out=stateout(:,4); w2out=stateout(:,5); w3out=stateout(:,6); omig1=stateout(:,7); omig2=stateout(:,8); omig3=stateout(:,9); figure(2),subplot(2,1,1), plot(t,th1out*r2d,'r:',t,th2out*r2d,'k--',t,th3out*r2d,'linewidth',2), xlabel('time(s)'),ylabel('姿态角(deg)'),title('质心坐标系（未控制）'); legend('\theta_1','\theta_2','\theta_3'),grid on subplot(2,1,2), plot(t,w1out*r2d,t,w2out*r2d,t,w3out*r2d,'linewidth',2), legend('\omega_1','\omega_2','\omega_3'),grid on xlabel('time(s)'),ylabel('角速率(deg/s)') figure(3) plot(t,omig1*r2d/360,'r:',t,omig2*r2d/360,'k--',t,omig3*r2d/360,'linewidth',2), xlabel('time(s)'),ylabel('角速度(n/s)'),title('飞轮角速度'); legend('\theta_1','\theta_2','\theta_3'),grid on