Re: A question regarding unwrap

[Joe Koski]

on 9/24/05 7:44 PM, Doug Stewart at address@hidden wrote:

> Ron Crummett wrote:
> 
>> Hi -
>> 
>> I have a recent problem with unwrap; I tried plotting the phase of a
>> transfer function and it jumped halfway through the plot, from -180 to
>> 180.  I tried to unwrap the phase using unwrap and plotted the phase
>> again, but still had the same jump.  Reading the help on unwrap, it
>> sounds like it will only really work if the phase is positive (it adds
>> 2*pi).  For me, I have a negative phase, so I need to subtract 2*pi.
>> Any idea of what the solution could be?  I tried looking at the code
>> but couldn't make too much sense out of it from my few glances.
>> 
>> Thanks,
>> 
>> Ron Crummett
>> 
>> 
> This is not about unwrap bur it might do what you want  and more  :-)
> Let me know if you want more help.
> Doug
> 
> 
> ## Copyright (C) 2004 doug
> ##
> ## 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     bode1 (n,d)
> %n= numerator
> % d=denominator
> %
> 
> ## bode1w
> 
> ## Author: doug
> ##
> ## 2004-09-27 doug
> ## * Initial revision
> ## Rev #1  Nov 29 2004 Doug Stewart
> 
> function     bode1w (n,d);
> 
> sys=tf2sys(n,d);
> [mag,ph,w]=bode(sys);
> 
> # fix the phase wrap arownd
> ll=length(ph);
> [qm iqm]=max(abs(diff(ph)));
> if(qm>300)
>  ph(iqm+1:ll)=ph(iqm+1:ll)-360;
> endif
> 
> 
> # make a 180 deg line
> v=ones(ll,1)*(-180);
> db=20*log10(mag);
> 
> 
Ron, Doug,

My experience with an alternative (Stearns and David, Signal Processing
Routines in Matlab, which is compatible with Octave) is that the phase
unwrap routines are good for 2*pi+, but then get confused. Maybe a more
robust version of unwrap is in order with an "if" (or some such) to check
the curve continuation.

Joe




-------------------------------------------------------------
Octave is freely available under the terms of the GNU GPL.

Octave's home on the web:  http://www.octave.org
How to fund new projects:  http://www.octave.org/funding.html
Subscription information:  http://www.octave.org/archive.html
-------------------------------------------------------------