On 08/05/12 13:26, Daniel J Sebald wrote:
> On 05/07/2012 09:29 PM, James Steward wrote:
>> On 08/05/12 11:54, Daniel J Sebald wrote:
>>
>>>
>>> octave:33> N=7
>>> N = 7
>>> octave:34> n = -ceil(N/2):N/2
>>> n =
>>>
>>> -4 -3 -2 -1 0 1 2 3
>>>
>>> Instead of defining 'n' the way the code does, it might be better to
>>> use the description from the Wiki page that says "w(n) = w_0(n -
>>> (N-1)/2)", this N not being the same as the N in the lines above. Note
>>> that blackmanharris.m uses an input L which is the length of window
>>> and N = L - 1.
>>
>> The definition of "n" in the wiki page says "n is an integer, with
>> values *0 ≤ n ≤ N-1*."
>
> The N in the Wiki is the length, equivalent to L in blackmanharris.m.
> So the N in blackmanharras.m is the same as N-1 in the Wiki.
> Confusing, yes.
>
It is the definition of "n" from the wiki that does not match the way it
is defined in the blackmanharris.m file or your proposal (just below)
that I was pointing out.
>>>
>>> I think the definition of n inside blackmanharris.m should be:
>>>
>>> n = (0:N) - N/2;
>>>
>>> Does that produce the result you are expecting James?
>>
>> I haven't tested this inside the signal toolbox blackmanharris.m file,
>> but it yields non integer values of "n" for L = even numbers, N being
>> then odd, and N/2 resulting in a non integer value.
>>
>> My version of blackmanharris.m, that I simply overloaded the system
>> default with, looks like this;
>>
>> function w = blackmanharris (N)
>>
>> w = 1;
>>
>> a0 = 0.35875;
>> a1 = 0.48829;
>> a2 = 0.14128;
>> a3 = 0.01168;
>>
>> for n=0:N-1
>> w(n+1) = a0 - a1 * cos(2*pi()*n/(N-1)) + a2 * cos(4*pi()*n/(N-1)) - a3 *
>> cos(6*pi()*n/(N-1));
>> endfor
>> endfunction;
>
> OK, so you are using N as the length. That matches with the Wiki
> then. However, loops are discouraged. Did you use a loop to avoid a
> conditional test that N be greater than zero?
I did say my code was not beautifully efficient.
> Otherwise, I'd suggest the following:
>
> function w = blackmanharris (N)
>
> if (N < 1)
> error('Window length must be 1 or greater');
> elseif (N == 1)
> w = 1;
> else
> a0 = 0.35875;
> a1 = 0.48829;
> a2 = 0.14128;
> a3 = 0.01168;
> n=0:N-1;
> w = a0 - a1 * cos(2*pi*n./(N-1)) + a2 * cos(4*pi*n./(N-1)) -
> a3 * cos(6*pi*n./(N-1));
> endif
>
> endfunction;
>
> The result looks similar to what I was getting from the change I
> suggested in the previous post. Note one peculiar artifact, which is
> that changing the phasing of the inputs (theoretically correct)
> results in a small round off error as reported by Octave:
>
> octave:181> for n=0:N-1
> * cos(6*pi()*n/(N-1)); a1 * cos(2*pi()*n/(N-1)) + a2 *
> cos(4*pi()*n/(N-1)) - a3
>> endfor
> octave:182> w
> w =
>
> Columns 1 through 6:
>
> 6.0000e-05 5.5645e-02 5.2057e-01 1.0000e+00 5.2058e-01
> 5.5645e-02
>
> Column 7:
>
> 6.0000e-05
>
> The above is indicating 5.2057e-01 and 5.2058e-01 for coefficients,
> but when I list them in long format, why the rounding occurs isn't
> evident. The rounding effect will likely show up somewhere no matter
> the formula.
Nice. Perhaps this can make it into the next signal processing package?
Regards,
James.
------------------------------------------------------------------------------
Live Security Virtual Conference
Exclusive live event will cover all the ways today's security and
threat landscape has changed and how IT managers can respond. Discussions
will include endpoint security, mobile security and the latest in malware
threats. http://www.accelacomm.com/jaw/sfrnl04242012/114/50122263/
_______________________________________________
Octave-dev mailing list
Octave-dev@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/octave-dev
