|
From: | dastew |
Subject: | RE: OCST: Gain Margin |
Date: | Fri, 3 Jul 2009 09:49:40 +0000 |
From: address@hidden To: address@hidden Subject: Re: OCST: Gain Margin Date: Fri, 3 Jul 2009 09:48:16 +0200
Thanks for your help! If I understood your code correctly, the following code should do the job: Regards, Lukas I cleaned this up and added some comments. If it works the way you want, then we should put it in the controls tool box. Doug Stewart ## Copyright (C) 2009 Doug Stewart and Lukas ## ## This program is free software; you can redistribute it and/or modify ## it under the terms of the GNU General Public License as published by ## the Free Software Foundation; either version 2 of the License, or ## (at your option) any later version. ## ## This program is distributed in the hope that it will be useful, ## but WITHOUT ANY WARRANTY; without even the implied warranty of ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ## GNU General Public License for more details. ## ## You should have received a copy of the GNU General Public License ## along with this program; if not, write to the Free Software ## Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA % %function function [gamma, phi, w_gamma, w_phi] = margin(sys); % gamma = Gain Margin ( as gain not dbs)% w_gamma = radian frequency for the Gain Margin % phi = Phase Margin ( in degrees) % w_phi = radian frequency for the Phase margin % % This function calls bode - see help bode for details about bode. % This Function was written for continuous time systems function [gamma, phi, w_gamma, w_phi] = margin(sys); if(! isstruct(sys) ) error(" The input must be a system type variable. margin(sys)") % Get Frequency Response [mag, pha, w] = bode(sys); % fix the phase wrap around phw = unwrap(pha * pi/180); pha = phw * 180/pi; % find the Gain Margin and its location [x, ix] = min(abs(pha + 180)) if (x > 1) gamma = "INF" ix = length(ph); else gamma = 1/mag(ix); endif w_gamma = w(ix); % find the Phase Margin and its location [pmx, ipmx] = min(abs(mag - 1)); phi = 180 + pha(ipmx); w_phi = w(ipmx); endfunction |
[Prev in Thread] | Current Thread | [Next in Thread] |