bugGNU Octave - Bugs: bug #57359, [octave forge] (signal) residuez...

 
 

bug #57359: [octave forge] (signal) residuez gives an incorrect result for [1 1 1.5 .5],[1 1.5 .5]

Submitter:  None
Submitted:  Thu 05 Dec 2019 03:01:35 AM UTC
   
 
Category:  Octave Package Severity:  3 - Normal
Priority:  5 - Normal Item Group:  Incorrect Result
Status:  Fixed Assigned to:  lostbard
Originator Name:  Andrew Stone Originator Email:  -email is unavailable-
Open/Closed:  * Closed Release:  * 5.1.0
Operating System:  * Any Fixed Release:  None
Planned Release:  None
* Mandatory Fields

Add a New Comment Rich Markup
   

Jump to the original submission

Fri 28 Oct 2022 12:06:05 PM UTC, comment #8: 


Signal 1.4.3 released - closing as fixed

John Donoghue <lostbard>
Group Member
Mon 19 Sep 2022 12:20:31 PM UTC, comment #7: 

ok - that makes it clear to me.

I pushed up the change and added a test for it.

https://sourceforge.net/p/octave/signal/ci/c48ea62659c49782b18d8eb7244fc5f6ade5bdc1/

John Donoghue <lostbard>
Group Member
Mon 19 Sep 2022 12:23:54 AM UTC, comment #6: 

Coming back to the original example from the OP. It's solution is


H(z) = (1+z^-1+1.5*z^-2+0.5*z^-3)/(1+1.5*z^-1+0.5*z^-2)
     = -1/(1-(-0.5)*z^-1) + 2/(1-(-1)*z^-1) + 0+1*z^-1


Therefore, the result from the residuez.m function must be either


r=[-1,2] and p=[-0.5,-1]

or

r=[2,-1] and p=[-1,-0.5]


Both expressions are correct because you can swap the fractions while the result remains the same:


[-1/(1-(-0.5)*z^-1)] + [2/(1-(-1)*z^-1)] = [2/(1-(-1)*z^-1)] + [-1/(1-(-0.5)*z^-1)]


What's wrong is the k (or f) argument, because the polynomial part in the PFE is


0+1*z^-1

and not

1+0*z^-1


In Matlab the k part is expressed as k(1)+k(2)*z^-1+.... Therefore the correct answer is [0,1] and not [1,0]. The reason why this is happening is explained in comment #2.

In the patch i have removed the if statement because fliplr([]) is IMO already a valid Octave statement.

Hg200 <hg200>
Sun 18 Sep 2022 08:49:55 PM UTC, comment #5: 

Same output on 8.0.0.
I see that all output vectors are flipped w.r.t. Matlab's answer.

Philip Nienhuis <philipnienhuis>
Group Member
Sun 18 Sep 2022 12:25:27 PM UTC, comment #4: 

With octave 7.2.1 for me, without the patch now returns:


[r,p,k] = residuez([1 1 1.5 .5],[1 1.5 .5])
r =

  -1
   2

p =

  -0.5000
  -1.0000

k =

   1   0


John Donoghue <lostbard>
Group Member
Sat 22 Feb 2020 07:14:13 PM UTC, comment #3: 

i will never learn. patch (untested) attached.

(file #48475)

Hg200 <hg200>
Sat 22 Feb 2020 06:56:41 PM UTC, comment #2: 

The if statement does not work with arrays. In the mentioned example the reverse flip of the direct term in residuez.m is not done because


if ([1,0])
  flip()
endif


is treated as false. A fix would be


if (! isempty (f))
  flip()
endif


or better: Removing the if clause completely, presuming flipr() can also treat empty variables.

I cannot install "signal" package for testing but an offline test (see attachment) looks fine.

Some reverse engineering of residuez.m as far as i can understand:

i.) Presume a Z-Transform

H(z)=[b0+b1*z^(-1)...+bm*z^(-m)]/[a0+a1*z^(-1)..+an*z^(-n)]
Parameters to call residuez (B,A) are:
B=[b0,b1,..,bm] and A=[a0,a1,..,an]


ii.)Substitute s:=z^-1 and do a PFE using the residue () function

H'(s)=H(1/s)=[b0+b1*s^(1)...+bm*s^(m)]/[a0+a1*s^(1)..+an*s^(n)]
Parameters to call residue (BF,AF) are:
BF=[bm,..,b1,b0] and AF=[an,..,a1,a0]


=> Conclusion is that the inputs BF=flip(B) and AF=flip(A) must be flipped before calling residue ()

iii.) Return values from residue () are poles(s) and direct terms

F(s)=f(m-n)*s^(m-n)+...+f0
F=[fm-n,..,f1,f0]
If we resubstitue s=z^(-1) we get
F(1/z)=f(m-n)*z^-(m-n)+...f1*z^-1+f0

and therefore, because residuez () returns the result as [f0,f1, ...] instead as [fm-n,...,f1,f0] a "back-flip" must be done.

This is not happening because the mentioned if [1,0] does not work.

Hg200 <hg200>
Thu 05 Dec 2019 10:50:36 AM UTC, comment #1: 

Confirmed with Matlab 2014a:

>> [r,p,k] = residuez([1 1 1.5 .5],[1 1.5 .5])
r =
     2
    -1
p =
   -1.0000
   -0.5000
k =
     0     1
>> ver
----------------------------------------------------------------------------------------------------
MATLAB Version: 8.3.0.532 (R2014a)
MATLAB License Number: ######
Operating System: Microsoft Windows 7 Enterprise  Version 6.1 (Build 7601: Service Pack 1)
Java Version: Java 1.7.0_11-b21 with Oracle Corporation Java HotSpot(TM) 64-Bit Server VM mixed mode
----------------------------------------------------------------------------------------------------
MATLAB                                                Version 8.3        (R2014a)
Signal Processing Toolbox                             Version 6.21       (R2014a)


Philip Nienhuis <philipnienhuis>
Group Member
Thu 05 Dec 2019 03:01:35 AM UTC, original submission:  

the octave residuez function when run on the transfer function [1 1 1.5 .5],[1 1.5 .5] gives results consistent with the matlab command, except for the K polynomial (or the FIR part, as octave calls it). In Octave the K polynomial is [1 0] whereas in matlab it is [0 1]. This is of course, assuming that the F row vector Octave returns is in the format: f(0)+f(1)z^-1... Working it out on paper shows the MATLAB answer to be correct, and the Octave answer to be incorrect.
minimal example:

pkg load signal
[r,p,k] = residuez([1 1 1.5 .5],[1 1.5 .5])

output:

r =

   2
  -1

p =

  -1.00000
  -0.50000

k =

   1   0

expected output (from MATLAB):

r =

   2
  -1

p =

  -1.00000
  -0.50000

k =

   0   1


Anonymous

 

(Note: upload size limit is set to 16384 kB, after insertion of the required escape characters.)

Attach Files:
   
   
Comment:
   

Attached Files
file #48475:  patch.diff added by hg200 (617B - text/x-patch)

 

Depends on the following items: None found

Items that depend on this one: None found

 

Carbon-Copy List
  • -email is unavailable- added by lostbard (Posted a comment)
  • -email is unavailable- added by hg200 (Posted a comment)
  • -email is unavailable- added by siko1056 (Updated the item)
  • -email is unavailable- added by None (Submitted the item)
  •  

    There are 0 votes so far. Votes easily highlight which items people would like to see resolved in priority, independently of the priority of the item set by tracker managers.

    Only group members can vote.

     

    Follow 9 latest changes.

    Date Changed by Updated Field Previous Value => Replaced by
    2022-10-28 lostbard StatusReady For Test Fixed
        Open/ClosedOpen Closed
    2022-09-19 lostbard StatusPatch Submitted Ready For Test
        Assigned toNone lostbard
    2020-03-27 mtmiller StatusConfirmed Patch Submitted
    2020-02-22 hg200 Attached File- Added patch.diff, #48475
    2019-12-05 philipnienhuis StatusNone Confirmed
        Operating SystemGNU/Linux Any
    2019-12-05 siko1056 Summaryresiduez gives an incorrect result for [1 1 1.5 .5],[1 1.5 .5] [octave forge] (signal) residuez gives an incorrect result for [1 1 1.5 .5],[1 1.5 .5]

    Back to the top

    Powered by Savane 3.13-758e.
    Corresponding source code