Re: [OctDev] signal package: possible bug on fir2 with patch (please review)

[Carnë D. <car...@gm...>]

On 14 August 2012 22:19, Mike Miller <mtm...@ie...> wrote:
> On Tue, Aug 14, 2012 at 4:54 PM, Carnë Draug wrote:
>> On 14 August 2012 21:33, Mike Miller <mtm...@ie...> wrote:
>>> On Tue, Aug 14, 2012 at 3:59 PM, Carnë Draug wrote:
>>>> On 14 August 2012 20:02, Mike Miller <mtm...@ie...> wrote:
>>>>> On Tue, Aug 14, 2012 at 1:53 PM, Carnë Draug wrote:
>>>>>> On 9 August 2012 03:11, Mike Miller <mtm...@ie...> wrote:
>>>>>>>
>>>>>>> Finally got back to this, it's now closer to what you showed here.
>>>>>>> It's clear that Matlab is forcing the 4th argument to always be a
>>>>>>> power of 2, even though that's not in the help text. I implemented
>>>>>>> that, added an error message for the case of the argument being too
>>>>>>> small (after conversion to a power of 2), and added a test case based
>>>>>>> on your examples above. Thanks for those!
>>>>>>
>>>>>> Hi Mike
>>>>>>
>>>>>> I have just noticed that this introduced another bug (which I think
>>>>>> was already there before, this just made it more apparent. For
>>>>>> example, if the minimum value for 4th argument is more than the
>>>>>> default, it also automatically (and silently) increases though I'm not
>>>>>> sure for how much. For example:
>>>>>
>>>>> Thanks for checking, I'll take a look at this again today when I get
>>>>> back to my dev box.
>>>>>
>>>>>>>> f = [0 0.1 0.1 0.2 0.2 1]; m = [0 0 1 1 0 0];
>>>>>>>> bdef = fir2(1024, f, m);
>>>>>>
>>>>>> should make the 4th argument default to 512 which would be too low and
>>>>>> return an error (current behaviour with Octave). However, this does
>>>>>> not happen in Matlab. But if the default value is specified, then it
>>>>>> gives an error. And it doesn't default to the next possible value
>>>>>> (1024) as I first expected, it defaults to 2048.
>>>>>>
>>>>>> [...]
>>>>>>
>>>>>> Let me know if you need more tests.
>>>>>
>>>>> Yes, can you repeat your tests using 1023 and 2047 instead? The limit
>>>>> is calculated based on the length of the filter which is order+1, I
>>>>> think that's why you're seeing it go up to the next power of two.
>>>>>
>>>>> bdef = fir2(1023, f, m);
>>>>> b512 = fir2(1023, f, m, 512); % should not error now
>>>>
>>>> Yes, it doesn't give an error.
>>>>
>>>>> b1024 = fir2(1023, f, m, 1024);
>>>>> isequal(bdef, b512) % should be true
>>>>
>>>> No, it's false
>>>>
>>>>> bdef = fir2(2047, f, m);
>>>>> b1024 = fir2(2047, f, m, 1024);
>>>>> isequal(bdef, b1024) % should be true
>>>>
>>>> No. Also false.
>>>
>>> Hmm, ok, so it allows it to be as low as 1/2 but defaults to a larger
>>> value, how about these, let's try and find the cutoff:
>>
>> isequal (fir2 (511, f, m), fir2 (511, f, m, 512))   # true
>> isequal (fir2 (512, f, m), fir2 (512, f, m, 512))   # true
>> isequal (fir2 (513, f, m), fir2 (513, f, m, 512))   # true
>> isequal (fir2 (512, f, m), fir2 (512, f, m, 1024))  # false
>> isequal (fir2 (513, f, m), fir2 (513, f, m, 1024))   #false
>> isequal (fir2 (1022, f, m), fir2 (1022, f, m, 512))  #true
>> isequal (fir2 (1022, f, m), fir2 (1022, f, m, 1024))  # false
>> isequal (fir2 (1023, f, m), fir2 (1023, f, m, 1024))  #true
>> isequal (fir2 (2046, f, m), fir2 (2046, f, m, 2048))  #true
>> isequal (fir2 (2046, f, m), fir2 (2046, f, m, 2048))  #true
>> isequal (fir2 (2047, f, m), fir2 (2047, f, m, 2048))  #true
>
> Great, I think I have enough to fix that one now.
>
>>>>> If you have patience for any more, I can come up with some tests to
>>>>> verify the 5th argument too.
>>>>
>>>> Sure, no problem. Give me whatever tests you need and I'll run them.
>>>> Same for other functions of the package or if you want to take a
>>>> closer look at differences between results.
>>>
>>> Based on how our fir2 is working now, the 5th arg defaults to the 4th
>>> arg / 20 (and ML's help says the defaults are 512 and 25), so:
>>> isequal (fir2 (63, f, m, 512), fir2 (63, f, m, 512, 25))
>>> isequal (fir2 (63, f, m, 1024), fir2 (63, f, m, 1024, 51))
>>> isequal (fir2 (63, f, m, 2048), fir2 (63, f, m, 2048, 102))
>>>
>>> I don't like my chances that these are right out of the box :)
>>
>> You are right to not make any bets. All 3 are false in Matlab 2011b.
>
> Ok, the docs are wrong yet again, brute force it:
>
> for i = [512, 1024, 2048, 4096]
>   for j = 1:fix(i/10)
>     if isequal (fir2 (63, f, m, i), fir2 (63, f, m, i, j))
>       sprintf('default for %d is %d\n', i, j)
>     end
>   end
> end

There you go:

default for 512 is 20
default for 1024 is 40
default for 2048 is 81
default for 4096 is 163

>From this, it seems to me that the default of j is "floor (i/25)".
Where on documentation did you saw that it was /20?

Carnë