东南大学matlab第三次大作业

Matlab Worksheet 3Part A1. Using function conv_m.m to make convolution between the following to functions (x and h):x=[3, 11, 7, 0, -1, 7, -5, 0, 2]; h=[11, 9, 0, -7, -3, 2, 0 -1]; nx=[-2:6]; nh=[0:7];Plot the functions and convolution results.x=[3, 11, 7, 0, -1, 7,5,0, 2]; nx=[-2:6];h=[11, 9, 0, -7, -3,2,0,-1]; nh=[0:7];[y, ny]=conv_m(x,nx,h,nh); subplot(3,1,1); stem(nx,x); ylabel('x[n]');axis([-6 10 -20 20]); subplot(3,1,2); stem(nh,h); ylabel('h[n]');axis([-4 10 -20 20]); subplot(3,1,3); stem(ny,y); xlabel('n'); ylabel('y[n]');axis([-6 15 -200 200]);2. Plot the frequency response over π≤Ω≤0for the following transfer function by lettingΩ=j e z , where Ωis the frequency (rad/sample)., with appropriate labels and title.9.06.1)(2++=z z zz H .delta=0.01;Omega=0:delta:pi;H= (exp(j .* Omega)) ./ ((exp(j .* Omega)).^2+1.6*exp(j .* Omega)+0.9); subplot(2,1,1);plot(Omega, abs(H)); xlabel('0<\Omega<\pi'); ylabel('|H(\Omega)|');axis([0 pi 0 max(abs(H))]); subplot(2,1,2);plot(Omega,atan2(imag(H),real(H))); xlabel('0<\Omega<\pi');ylabel(' -\pi < \Phi_H <\pi') axis([0 pi -pi pi]);3. Use fft to analyse following signal by plotting the original signal and its spectrum.⎪⎭⎫ ⎝⎛+⎪⎭⎫ ⎝⎛+⎪⎭⎫ ⎝⎛=10244672sin 10241372sin 1024322sin ][n n n n x πππ.% fft N=1024; dt=1/N;t=0:dt:1-dt;x=sin(2*pi*32 .*t) +sin(2*pi*137 .*t)+sin(2*pi*467 .*t) ; subplot(2,1,1);plot(t,x), xlabel('t sec'),ylabel('x'); title('Signal and its Fourier Transform'); axis([min(t) max(t) 1.5*min(x) 1.5*max(x)]); X=fft(x); f=0:N-1;subplot(2,1,2);stem(f,abs(X)), xlabel('Hz'), ylabel('|X|'); axis([ 0 N/2 0 1.5*max(abs(X))]);4. Use the fast Fourier transform function fft to analyse following signal. Plot the original signal,and the magnitude of its spectrum linearly and logarithmically. Apply Hamming window to reduce the leakage.⎪⎭⎫ ⎝⎛+⎪⎭⎫ ⎝⎛+⎪⎭⎫ ⎝⎛=10247.4672sin 10244.1372sin 10245.322sin ][n n n n x πππ.The hamming window can be coded in Matlab asfor n=1:Nhamming(n)=0.54+0.46*cos((2*n-N+1)*pi/N); end;where N is the data length in the FFT.%Examine signal components %sampling rate 1024Hz N=1024; dt=1/N;t=0:dt:1-dt;%****define the signal**x=0.1*sin(2*pi*32.5 .*t ) + ... 0.2*sin(2*pi*137.4 .*t )+ ... 0.15*sin(2*pi* 467.7 .*t ) ; %*********************** subplot(3,1,1); plot(t,x);axis([0 1 -1 1]);xlabel('nT (seconds)'); ylabel('x[n]');%hamming window for n=1:Nw(n)=0.54+0.46*cos((2*n-N+1)*pi/N); end ;x1=x .*w;subplot(3,1,2); plot(t,x1);axis([0 1 -5 5]);xlabel('nT (seconds)'); ylabel('x[n] Hamming win'); X=fft(x1);df=1;f=0:df:N-1;subplot(3,1,3);plot(f,20*log10(abs(X))); axis([0 N/2 -60 60]);xlabel('k (Hz)');ylabel('|X[k]| (db)');Part BSimulation Using SIMULINKINTRODUCTIONThe objective of this laboratory is to learn about various properties of signals and systems by doing simulations in SIMULINK.PART 1. BASICS OF SIMULINK1.1 UNDERSTANDING SIMPLE WAVEFORMS AND INTEGRATIONCreate a pulse of height 2 units from time 0 to 4 seconds by subtracting two unit steps and adding a gain. Connect this pulse to an integrator with a gain of 0.5 and a zero initial condition. Connect oscilloscopes to show the pulse and the output of the integrator. You may wish to name your simulation (block) diagram; to do so use the save as feature under edit. Your block diagram should be similar as below.Before simulating you need to pull down the simulation header and double click on parameters. Unless told otherwise always assume that you can use the ode 45 integration algorithm shown in this window and the other default parameters. Typically you will only alter the start and stop times but for this first simulation you can use the default values of 0 and 10 seconds. Double click on the oscilloscopes to get the windows in which the traces will appear, pull down the simulate menu and click on run. Plot below the integrator input and output waveforms.Repeat the experiment but with an initial condition of –2 on the integrator. Again draw the results.Repeat the experiment but with an initial condition of –2 on the integrator. Again draw the results.1.2 FIRST ORDER SYSTEMA single time constant may be simulated using the transfer function block in which you enter the coefficients of the numerator and denominator polynomials. Set up the configuration in a new SIMULINK window to realise the transfer function 1/(s+2) with the input unit step and an oscilloscope connected to the output of the transfer function block. Plot the block diagram in the space below. Simulate the system for 5 seconds and plot the response.1.3 SECOND ORDER SYSTEMFor the second order all pole transfer function )2/(200220ωωζω++s s you will recall that if a time scale of ω0t is used for plotting the step response, the response shape will only be affected by changes in the damping ratio ζ. This can also be shown if we normalise the transfer function by replacing (s/ω0) by s n togive )12/(12++n n s s ζ. To study the effect of varying ζ on the step response we will therefore use the transfer function 1/(s 2+2ζs+1). Set up the following configuration for a simulation study.i) A unit step input. ii) Connect the unit step to a transfer function of 1/(s 2+2ζs+1) with ζ = 1.0. iii) Take a summing block and connect the input step to a +ve input and the transfer functionoutput to a –ve input.iv) Connect the output of this summer to a square function. This is obtained by using f(u) in thenonlinear blocks.v)Connect the output of this squarer to an integrator.vi)Connect two oscilloscopes to the circuit, one to the transfer function output and the other to the output of the integrator.vii)Also connect simout blocks to the same signals as the oscilloscopes. Rename the one connected to the integrator output simout1.Simulate the system with the values of ζ listed in the table below and fill in the other figures from the simulation results. To get an accurate value at the end of each run you can type simout1 in the MATLAB window. You can also measure the overshoot by making use of the maximum command: simply type max(simout) in the MATLAB window.% overshoot 2.1076 1.2500 1.1778 1.0571 1.0000 1.1285 2.1431PART 2. SIMULATION OF AIRCRAFT PITCH ANGLE AND ALTITUTEThe purpose IS to use "SIMULINK" to simulate a (much simplified) model of the longitudinal motion of a fighter aircraft.The "angle of attack", is the angle between the direction a plane is pointing, and the direction in which it actually moves through the air. For a plane flying at approximately constant altitude, this is equivalent to the "pitch angle", φ, as illustrated in Fig.2.1. This angle is important because it produces a lift force perpendicular to the axis of the plane, and hence a "normal acceleration", a n , (also shown in the figure).vThe pilot wants to be able to control the pitch angle （倾斜角度）, and does so ultimately by rotating the front fins, and tail elevators of the aircraft, shown in Fig. 2.2. The first task is therefore to model theeffect of these movements on the "pitch rate" q ()=φ, and "normal acceleration" a n .Fig. 2.2 Illustration of control surfaces2.1: MODELLING THE AIRCRAFT AND ACTUATOR DYNAMICSi) Normal Acceleration, a nThe acceleration of the aircraft in a direction perpendicular to its axis, (the "normal accel.", a n ), is determined mainly by the angle θ1 of the tail elevators of the aircraft shown in Fig. 2.2. Indeed aerodynamic modelling shows that this relationship can be described by the differential equation,11116342.364.1982.81.3θθθ++=-+ n n n a a aConvert this relationship into a transfer function form:-ii) Pitch Rate, qThe rate at which the pitch angle changes, (the "pitch rate", q ), is determined mainly by the angle θ2 of the front fins of the aircraft shown in Fig. 2.2. Indeed aerodynamic modelling shows that this relationship can be described by the differential equation,2227.1125.782.81.3θθ++=-+ q q qConvert this relationship into a transfer function form:-iii) Actuator Dynamics + GearsThe tail elevators of the aircraft are driven directly by a hydraulic actuator which has a transfer function,G s s ()()=-+1414/ (3)To a unit step input U(s)=1/s , the response R(s) = G(s)⋅U(s)= 1411)14(14++-=+-s s s s . What form ofstep response r(t) would you expect from this actuator, and what is its time constant?The front fins and the tail elevators are both driven from the actuator by the same drive shaft through a gear box. The same gear on this drive shaft connects to the tail elevator gear wheel, which has 500 teeth, and the front fin gear wheel, which has 100 teeth, the relationship between θ1 and θ2 is52112==t t θθ(4) 2.2: SIMULATING THE AIRCRAFT AND ACTUATOR DYNAMICSi) F or the aircraft model you will need from the Continuous or Maths Library:-∙a Transfer Fcn block to represent the "tail elevator angle → normal accel" relationship (eqn (1))∙ a Transfer Fcn block to represent "front fin angle → pitch rate" relationship, (see eqn (2)) ∙ a Transfer Fcn block to represent the hydraulic "actuator dynamics", (see eqn (3)) ∙ a Gain block to represent the "gears", (see eqn (4))from the Sources Library:-∙ a Step Input block to act as a test input - (set Step time = 5sec, and Final value = 0.01 rad)from the Sinks Library:- ∙a ‘Scope block to monitor the pitch rate output - (set Horizontal Range = 20sec).ii) Now connect the blocks together, using the step input to drive the actuator, and the actuator to drive the tail elevators. The actuator is also used to drive the front fins via the gearbox. The 'Scope should be connected to the output of the "front fin angle → pitch rate" block. Note: you can take a branch from a signal line by dragging with the right mouse button from the point on the line at which you want the branch to start.Plot the block diagram you have created in the box below:-iii) From the pull-down menus of your simulation window, select Simulation/Parameters , in order to define how the simulation is to be performed. Set:-Stop Time = 20 sec; Initial Step Size = 0.001; Max Step Size = 0.01;iv) Finally double click on your ‘scope, so that you can watch the progress of the simulation, and then select Simulation/Start from the pull-down menu of your simulation window. You may want to adjustthe Vertical Range of your 'Scope, and re-run the simulation to get a good picture. Plot the response of the aircraft pitch rate below, and comment on these results. Note that it can be shown that if a transfer function has a denominator polynomial with a zero or negative coefficient then it must have a root with a zero or positive real part.2.3: SIMULATING THE PITCH RATE CONTROL SYSTEMThe designers are quite glad they simulated the response of the aircraft, before trying it for real. They now decide to improve the response by measuring the output (ie. what is actually happening to the pitch rate), and subtracting this from the input (ie. what they would like to happen), to produce an error signal. The error signal will be amplified, and then used to drive the actuator (see Fig. 2.3). In this way the actuator will automatically act so as to reduce any differences between the input demand, and the output response.Fig. 2.3 Feedback control of pitch ratei) Before you simulate the control system described above, you will probably find that you are short of space in your simulation window. Besides it would be nice to package the aircraft + actuator dynamics as a single block, so that they are distinct from anything added later as shown in Fig. 2.3.If you press the left mouse button in an empty space of the window, you get a "rubber band" box which you can drag to surround a group of elements. When you release the button, all elements inside are "selected". Use this technique to select all the aircraft + actuator dynamics, (but leave out the step input and ‘scope monitor). From the pull-down menu under Edit select Create a Subsystem, and your diagram will suddenly look as shown in Fig. 2.4.Fig. 2.4 New diagramii) Rename the "Subsystem" as "Actuator + Aircraft Dynamics" as shown in Fig. 2.3 by editing the text below the block. You can see what is inside by double clicking on the block itself, - try it! Notice that the input and output connections are labelled "in_1", "out_1" etc. Edit these to give them more meaningful names.iii) Now construct the planned control system, as illustrated in Fig. 2.3. Notice the pilot is now represented by a Signal Generator (from the Sources library), which should be set to produce a Square Wave, of frequency 0.1 Hz and peak amplitude 0.5.iv) The designers are unsure what value to use for the Gain. Open the Pitch Rate 'Scope, and run the simulation first for the Gain = -0.5 (as shown above), and then for the Gain = -5. Plot the two responses you obtain, and comment on the relative advantages and disadvantages of each:-Comments on responses (in terms of overshot and response time etc):PART 3: SIMILATION OF NONLINEAR SHIP ROLL DYNAMICSThe rolling motion of ships is of considerable interest to naval architects because even today around 50% of ships lost at sea, sink as a result of a capsize. The aim of this laboratory is to study such behaviour by:-∙constructing a nonlinear differential equation model for the system∙converting this model into a phase variable block diagram form suitable for simulation∙simulating the (nonlinear) response to sinusoidal excitation∙calculating the steady state response for a ship whose cargo has shifted dangerously to one sideFig. 3.1 Ship schematic3.1: MODELLING AND SIMULATION OF SHIP DYNAMICSi) Consider the ship sketch shown in Fig. 3.1. Given that:-∙ The effective inertia of the ship about its roll axis = J (1) ∙ The damping moment / torque due to friction between the body of the ship and the water, and due to turbulence round the "bilge keels" = c c 123 θθ+ (2) ∙ The restoring moment / torque due to the buoyancy and shape of the ship, as it is pushed to oneside = k k k 12335θθθ++ (3) ∙ The input / forcing moment due to the wave forces acting on the ship = u t ()- Using Newton’s second law of motion ∑=Torque Moment J / θ , the nonlinear differential equation which describes the rolling motion of the ship can be written as:-()()53321321)(θθθθθθk k k c c t u J ++-+-= (4)(ii) Re-arrange equation (4) to obtain an expression for the highest derivative of the output signal θ:-()()()53321321)(1θθθθθθk k k c c t u J++-+-= (5)iii) Hence complete the drawing of the SIMULINK block diagram for the ship roll system:-3.2: SIMULATING THE SHIP DYNAMICSi) For the ship model you will need:-from the Continuous and Maths Library:-∙ a few Integrator blocks - (assume zero initial conditions for now)∙ a few Gain blocks (not absolutely necessary)from the Nonlinear Library:-∙ a couple of Fcn blocks to implement the damping and stiffness terms (eqns (2), and (3)) from the Sources Library:- ∙ a Sine Wave block to simulate the sea waves u(t) - (set Ampl = 0.3 , and Freq = 0.6 rad/s)from the Sinks Library:-∙ a To Workspace block for the output ship roll angle, θ, of the simulation. This stores the results in a MATLAB variable which we can plot later, -(set the Variable name = theta, and the Maximum number of timesteps (ie points stored) = 5000)∙ an oscilloscope if you wish to see a waveform whilst the simulation is taking placeii) Drag the appropriate blocks over into your simulation window using the mouse, and double click on them to enter the appropriate parameter values. The numerical values for the coefficients of eqn (4) are:-;150;20;8.0;8.0;02.0;132121======k k k c c J(6)Note these values are given in consistent units so that θ is in radians. Finally by clicking on the text under each block, you can also enter a suitable description for each.iii) Now connect the blocks together, according to the phase variable block diagram you have plotted in section A. Note that some blocks will need changing so they "face the other way"; you can do this by selecting the block and "Flipping Horizontally " using the Options pulldown menu. Remember too that you can take a branch from a signal line by dragging with the right mouse button from the point on the line at which you want the branch to start.iv) From the pull-down menus of your simulation window, select Simulation/Parameters , in order to define how the simulation is to be performed. Set:-Stop Time = 10*pi (sec); Initial Step Size = 0.001; Max Step Size = 0.1;v) Now run the simulation by selecting Simulation/Start from the pull-down menu of your simulation window. If you do not use an oscilloscope you probably won't see much happening because the results are being stored in the MATLAB variable called 'theta', instead of being displayed immediately. After a few seconds however, you should hear a slight beep indicating the simulation is complete. Use of an oscilloscope is, however, recommended so that you can see if things ‘look all right’.vi)Now plot your results by typing plot(theta)at the MATLAB prompt in the MATLAB command window. Plot the response below, and comment on how it differs from what you would expect from a linear system if this were also excited with a sine wave input.。