I'm attempting to use convol() to implement a simple filter. Scilab's Help entry is terse to say the least. I suspect I should be using the overlap and add variant, but Scilab gives no example. I've experimented and wrote the following code. Calculating the convolution by both methods. The results can only be described as "similar". What am I missing? TIA t=[0:.00001:1]; s=size(t,2); h= exp(-t(1:101)/.001); // impulse response f1=123; / arbitrary test frequencies f2=5678; f3=76767; signal=sin(2*%pi*f1*t) +sin(2*%pi*f2*t) +sin(2*%pi*f3*t); // test signal rco=convol(h,signal); sizerco=size(rco); rco=rco(1:s)/max(rco); N=int(s/100); for i=[1:N] k=i-1; x(i,:)=signal((1+100*k):(101+100*k)); // I was guessing that h and x(i) // should be same size end xx=x(1,:); [yy,ee] =convol(h,xx); // calc y(1) y=yy; for k=[2:N-1] xx=x(k,:); [yy,ee] =convol(h,xx,ee); // calc y(2)->y(N-1) y=[y,yy]; end k=k+1; xx=x(k,:); yy=convol(h,xx,ee); // calc y(N) y=[y,yy]; sizey=size(y) sizerco yin=signal(1:s)/max(signal(1:s)); yout=y(1:s)/max(y(1:s)); // plot2d(t(1:s)',[yin(1:s)' yout(1:s)']) plot2d(t(1:s)',[rco(1:s)' yout(1:s)'])

# Attempting DSP in Scilab -- filtering using convol()

Started by ●July 9, 2008

Reply by ●July 9, 20082008-07-09

On Jul 9, 12:52�pm, Richard Owlett <rowl...@atlascomm.net> wrote:> I'm attempting to use convol() to implement a simple filter. > Scilab's Help entry is terse to say the least. > > I suspect I should be using the overlap and add variant, but Scilab > gives no example. > > I've experimented and wrote the following code. Calculating the > convolution by both methods. The results can only be described as > "similar". What am I missing? > > TIA > > t=[0:.00001:1]; > s=size(t,2); > h= exp(-t(1:101)/.001); � � � � � � � � � � // impulse response > > f1=123; � / arbitrary test frequencies > f2=5678; > f3=76767; > signal=sin(2*%pi*f1*t) +sin(2*%pi*f2*t) +sin(2*%pi*f3*t); // test signal > > rco=convol(h,signal); > sizerco=size(rco); > rco=rco(1:s)/max(rco); > > N=int(s/100); > for i=[1:N] > � �k=i-1; > � �x(i,:)=signal((1+100*k):(101+100*k)); � // I was guessing that h and > x(i) > � � � � � � � � � � � � � � � � � � � � � �// should be same size > end > > xx=x(1,:); > [yy,ee] =convol(h,xx); � � � � // calc y(1) > y=yy; > > for k=[2:N-1] > � �xx=x(k,:); > � �[yy,ee] =convol(h,xx,ee); � // calc y(2)->y(N-1) > � �y=[y,yy]; > end > > k=k+1; > xx=x(k,:); > yy=convol(h,xx,ee); � � � � � // calc y(N) > y=[y,yy]; > > sizey=size(y) > sizerco > > yin=signal(1:s)/max(signal(1:s)); > yout=y(1:s)/max(y(1:s)); > // plot2d(t(1:s)',[yin(1:s)' yout(1:s)']) > plot2d(t(1:s)',[rco(1:s)' yout(1:s)'])Richard, I only looked at the first few lines of what you have. It looks like you have a sample rate of 100Ksps and one input frequency of 76.767KHz, which is > 1/2 the sample rate, so it will alias, possibly into the passband of the filter (which you didn't specify). That might account for some conflict between your expectations and results. Dirk.

Reply by ●July 9, 20082008-07-09

On Jul 9, 2:40�pm, dbell <bellda2...@cox.net> wrote:> On Jul 9, 12:52�pm, Richard Owlett <rowl...@atlascomm.net> wrote: > > > > > > > I'm attempting to use convol() to implement a simple filter. > > Scilab's Help entry is terse to say the least. > > > I suspect I should be using the overlap and add variant, but Scilab > > gives no example. > > > I've experimented and wrote the following code. Calculating the > > convolution by both methods. The results can only be described as > > "similar". What am I missing? > > > TIA > > > t=[0:.00001:1]; > > s=size(t,2); > > h= exp(-t(1:101)/.001); � � � � � � � � � � // impulse response > > > f1=123; � / arbitrary test frequencies > > f2=5678; > > f3=76767; > > signal=sin(2*%pi*f1*t) +sin(2*%pi*f2*t) +sin(2*%pi*f3*t); // test signal > > > rco=convol(h,signal); > > sizerco=size(rco); > > rco=rco(1:s)/max(rco); > > > N=int(s/100); > > for i=[1:N] > > � �k=i-1; > > � �x(i,:)=signal((1+100*k):(101+100*k)); � // I was guessing that h and > > x(i) > > � � � � � � � � � � � � � � � � � � � � � �// should be same size > > end > > > xx=x(1,:); > > [yy,ee] =convol(h,xx); � � � � // calc y(1) > > y=yy; > > > for k=[2:N-1] > > � �xx=x(k,:); > > � �[yy,ee] =convol(h,xx,ee); � // calc y(2)->y(N-1) > > � �y=[y,yy]; > > end > > > k=k+1; > > xx=x(k,:); > > yy=convol(h,xx,ee); � � � � � // calc y(N) > > y=[y,yy]; > > > sizey=size(y) > > sizerco > > > yin=signal(1:s)/max(signal(1:s)); > > yout=y(1:s)/max(y(1:s)); > > // plot2d(t(1:s)',[yin(1:s)' yout(1:s)']) > > plot2d(t(1:s)',[rco(1:s)' yout(1:s)']) > > Richard, > > I only looked at the first few lines of what you have. It looks like > you have a sample rate of 100Ksps and one input frequency of > 76.767KHz, which is > 1/2 the sample rate, so it will alias, possibly > into the passband of the filter (which you didn't specify). �That > might account for some conflict between your expectations and results. > > Dirk.- Hide quoted text - > > - Show quoted text -Richard, I see that you actually did specify your filter, so I will correct myself. Aliasing comment still correct. Dirk

Reply by ●July 9, 20082008-07-09

Richard Owlett wrote:> I'm attempting to use convol() to implement a simple filter. > Scilab's Help entry is terse to say the least. > > I suspect I should be using the overlap and add variant, but Scilab > gives no example. > > I've experimented and wrote the following code. Calculating the > convolution by both methods. The results can only be described as > "similar". What am I missing? > > TIA >(code snipped) What are you trying to accomplish? I tested it out on a few simple cases; it appears to do a nicely padded-out convolution, so if you want to do 'straight' convolution you just do the call without the overlap-add nonsense. I really don't know what the overlap-add is about, unless it's to more efficiently convolve a pair of sequences with widely disparate lengths; I'd be tempted to just do it all in one step for starters. You could always examine the code to figure out it's functionality, and file a bug against the help entry with a suggested wording. -- Tim Wescott Wescott Design Services http://www.wescottdesign.com Do you need to implement control loops in software? "Applied Control Theory for Embedded Systems" gives you just what it says. See details at http://www.wescottdesign.com/actfes/actfes.html

Reply by ●July 9, 20082008-07-09

dbell wrote:> On Jul 9, 2:40 pm, dbell <bellda2...@cox.net> wrote: > >>On Jul 9, 12:52 pm, Richard Owlett <rowl...@atlascomm.net> wrote: >> >> >> >> >> >> >>>I'm attempting to use convol() to implement a simple filter. >>>Scilab's Help entry is terse to say the least. >> >>>I suspect I should be using the overlap and add variant, but Scilab >>>gives no example. >> >>>I've experimented and wrote the following code. Calculating the >>>convolution by both methods. The results can only be described as >>>"similar". What am I missing? >> >>>TIA >>[snip code] >>Richard, >> >>I only looked at the first few lines of what you have. It looks like >>you have a sample rate of 100Ksps and one input frequency of >>76.767KHz, which is > 1/2 the sample rate, so it will alias, possibly >>into the passband of the filter (which you didn't specify). That >>might account for some conflict between your expectations and results. >> >>Dirk.- Hide quoted text - >> >>- Show quoted text - > > > Richard, > > I see that you actually did specify your filter, so I will correct > myself. Aliasing comment still correct. > > DirkAlias folded into "harmless" location. Just gave "another" noise source. Moved it down to ~36 kHz. Thanks.

Reply by ●July 9, 20082008-07-09

Tim Wescott wrote:> Richard Owlett wrote: > >> I'm attempting to use convol() to implement a simple filter. >> Scilab's Help entry is terse to say the least. >> >> I suspect I should be using the overlap and add variant, but Scilab >> gives no example. >> >> I've experimented and wrote the following code. Calculating the >> convolution by both methods. The results can only be described as >> "similar". What am I missing? >> >> TIA >> > (code snipped) > > What are you trying to accomplish?Immediate goal was to figure out how to use Scilab's convol() routine. Motivation was need to filter a signal with code speed taking a distant second to code correctness. Not being expert in programming NOR DSP that indicated using a canned routine.> I tested it out on a few simple > cases; it appears to do a nicely padded-out convolution, so if you want > to do 'straight' convolution you just do the call without the > overlap-add nonsense. > > I really don't know what the overlap-add is about, unless it's to more > efficiently convolve a pair of sequences with widely disparate lengths; > I'd be tempted to just do it all in one step for starters.If it's mathematically correct to convolve sequences with widely differing lengths, I'll do it.> > You could always examine the code to figure out it's functionality, and > file a bug against the help entry with a suggested wording. >First I'd have to figure out DSP ;) Thanks