The subject line says it all. I need to resample some data in Scilab -- I'm getting data sampled at 5MHz that I need to resample to 4.9 some-odd MHz. When I was using 20ms long data sets it worked very well to do an FFT, trim bins out of the middle to reduce the frequency by the right factor, then do an IFFT. Now we're doing the "real" algorithm with 1/2 second long data sets, and it's taking way too long to do the resampling. I'm thinking of doing the resampling in chunks of 3000 or so (to hit the least common denominator of the two frequencies), either using my magic FFT method, or using a polyphase filtering method based on a multiply by a 3000^2 sparse matrix. Is this the right way to go? Do the Scilab experts have any suggestions from Scilab space? Do the DSP experts have any suggestions from DSP space? I'd like to do this in a "Scilab kosher" way initially, although I'm getting close to porting the whole mess over to C++. Thanks a bundle. -- Tim Wescott Control systems and communications consulting http://www.wescottdesign.com Need to learn how to apply control theory in your embedded system? "Applied Control Theory for Embedded Systems" by Tim Wescott Elsevier/Newnes, http://www.wescottdesign.com/actfes/actfes.html
Scilab, sparse matricies, and polyphase filtering
Started by ●December 15, 2008
Reply by ●December 15, 20082008-12-15
On Dec 15, 11:18 am, Tim Wescott <t...@seemywebsite.com> wrote:> The subject line says it all. > > I need to resample some data in Scilab -- I'm getting data sampled at > 5MHz that I need to resample to 4.9 some-odd MHz. When I was using 20ms > long data sets it worked very well to do an FFT, trim bins out of the > middle to reduce the frequency by the right factor, then do an IFFT. > > Now we're doing the "real" algorithm with 1/2 second long data sets, and > it's taking way too long to do the resampling. > > I'm thinking of doing the resampling in chunks of 3000 or so (to hit the > least common denominator of the two frequencies), either using my magic > FFT method, or using a polyphase filtering method based on a multiply by > a 3000^2 sparse matrix. > > Is this the right way to go? Do the Scilab experts have any suggestions > from Scilab space? Do the DSP experts have any suggestions from DSP > space? I'd like to do this in a "Scilab kosher" way initially, although > I'm getting close to porting the whole mess over to C++. > > Thanks a bundle. > > -- > Tim Wescott > Control systems and communications consultinghttp://www.wescottdesign.com > > Need to learn how to apply control theory in your embedded system? > "Applied Control Theory for Embedded Systems" by Tim Wescott > Elsevier/Newnes,http://www.wescottdesign.com/actfes/actfes.htmlCombine your FFT size change with the overlap-add or overlap-save filter methods. Then your smaller FFTs may fit in cache and run faster than a 2.5 megasample FFT. Assuming your data are well anti-aliased, the filter function you perform can be an easy one to require little overlap. Dale B. Dalrymple http://dbdimages.com
Reply by ●December 15, 20082008-12-15
On Dec 15, 2:18�pm, Tim Wescott <t...@seemywebsite.com> wrote:> The subject line says it all. > > I need to resample some data in Scilab -- I'm getting data sampled at > 5MHz that I need to resample to 4.9 some-odd MHz. �When I was using 20ms > long data sets it worked very well to do an FFT, trim bins out of the > middle to reduce the frequency by the right factor, then do an IFFT. > > Now we're doing the "real" algorithm with 1/2 second long data sets, and > it's taking way too long to do the resampling. > > I'm thinking of doing the resampling in chunks of 3000 or so (to hit the > least common denominator of the two frequencies), either using my magic > FFT method, or using a polyphase filtering method based on a multiply by > a 3000^2 sparse matrix. > > Is this the right way to go? �Do the Scilab experts have any suggestions > from Scilab space? �Do the DSP experts have any suggestions from DSP > space? �I'd like to do this in a "Scilab kosher" way initially, although > I'm getting close to porting the whole mess over to C++. > > Thanks a bundle. > > -- > Tim Wescott > Control systems and communications consultinghttp://www.wescottdesign.com > > Need to learn how to apply control theory in your embedded system? > "Applied Control Theory for Embedded Systems" by Tim Wescott > Elsevier/Newnes,http://www.wescottdesign.com/actfes/actfes.htmlHello Tim, A couple of basic questions. You mention 5.0 to approx 4.9 MHz - is this ratio known (I assume it is) and does it work out to be something irrational or nearly so. You mention 3000 as a common chuncking size. Of course 3000 doesn't divide evenly into 5MHz. If 5.0 to 4.9 MHz then it is simply 49/50 resampling which you can factor into a cascade of two 7/5 interpolator/decimator stages followed by a 1/2 stage. But you hint that your ratio is not that nice. Also how is the spectral occupancy. Is everything filled right up to 2.5MHz? If you have some excess bandwidth, then you may use a raised cosine interpolator to design all of your coefs for your various phases. The advantage of a RC interpolator is if you have many phases (caused by sampling rates linked through large prime numbers) is you can design the coefs without using a remez method which is not very good if you end up with 100s of phases. So more details may be helpful. Clay
Reply by ●December 15, 20082008-12-15
Tim Wescott wrote:> The subject line says it all. > > I need to resample some data in Scilab -- I'm getting data sampled at > 5MHz that I need to resample to 4.9 some-odd MHz. �When I was using 20ms > long data sets it worked very well to do an FFT, trim bins out of the > middle to reduce the frequency by the right factor, then do an IFFT. > > Now we're doing the "real" algorithm with 1/2 second long data sets, and > it's taking way too long to do the resampling. > > I'm thinking of doing the resampling in chunks of 3000 or so (to hit the > least common denominator of the two frequencies), either using my magic > FFT method, or using a polyphase filtering method based on a multiply by > a 3000^2 sparse matrix. > > Is this the right way to go?I doubt it.>�Do the Scilab experts have any suggestions > from Scilab space? �Do the DSP experts have any suggestions from DSP > space? �I'd like to do this in a "Scilab kosher" way initially, although > I'm getting close to porting the whole mess over to C++Polyphase filtering is nothing more than leaving out multiplications with zero (I know you know). Since these are rational sampling rate changes, a program for such sampling rate changes has a period, so everything can be encoded in a huge "for" loop. Matlab has "upfirdn", I can't believe somebody hasn't already done the same thing in Scilab. However, there are lot of audio sample rate converters available for free - there must also be one that does arbitrary ratios (try r8brain). Just store your 5 MHz data in a wav file (with nominal 50 kHz sample rate, for example) and convert to the appropriate frequency (49 something kHz). That's probably also a lot faster than rolling your own in Scilab (does Scilab have a compiler, similar like the MEX C/C++ interface for Matlab?). Regards, Andor
Reply by ●December 15, 20082008-12-15
On Mon, 15 Dec 2008 12:55:03 -0800, clay wrote:> On Dec 15, 2:18 pm, Tim Wescott <t...@seemywebsite.com> wrote: >> The subject line says it all. >> >> I need to resample some data in Scilab -- I'm getting data sampled at >> 5MHz that I need to resample to 4.9 some-odd MHz. When I was using >> 20ms long data sets it worked very well to do an FFT, trim bins out of >> the middle to reduce the frequency by the right factor, then do an >> IFFT. >> >> Now we're doing the "real" algorithm with 1/2 second long data sets, >> and it's taking way too long to do the resampling. >> >> I'm thinking of doing the resampling in chunks of 3000 or so (to hit >> the least common denominator of the two frequencies), either using my >> magic FFT method, or using a polyphase filtering method based on a >> multiply by a 3000^2 sparse matrix. >> >> Is this the right way to go? Do the Scilab experts have any >> suggestions from Scilab space? Do the DSP experts have any suggestions >> from DSP space? I'd like to do this in a "Scilab kosher" way >> initially, although I'm getting close to porting the whole mess over to >> C++. >> >> Thanks a bundle. >> >> -- >> Tim Wescott >> Control systems and communications >> consultinghttp://www.wescottdesign.com >> >> Need to learn how to apply control theory in your embedded system? >> "Applied Control Theory for Embedded Systems" by Tim Wescott >> Elsevier/Newnes,http://www.wescottdesign.com/actfes/actfes.html > > Hello Tim, > > A couple of basic questions. You mention 5.0 to approx 4.9 MHz - is this > ratio known (I assume it is) and does it work out to be something > irrational or nearly so. You mention 3000 as a common chuncking size. Of > course 3000 doesn't divide evenly into 5MHz. If 5.0 to 4.9 MHz then it > is simply 49/50 resampling which you can factor into a cascade of two > 7/5 interpolator/decimator stages followed by a 1/2 stage. But you hint > that your ratio is not that nice. > > Also how is the spectral occupancy. Is everything filled right up to > 2.5MHz? If you have some excess bandwidth, then you may use a raised > cosine interpolator to design all of your coefs for your various phases. > The advantage of a RC interpolator is if you have many phases (caused by > sampling rates linked through large prime numbers) is you can design the > coefs without using a remez method which is not very good if you end up > with 100s of phases. So more details may be helpful. > > ClayThe bandwidth is about 2.5MHz, so 1.2 some-odd on either side of zero. The data is right at 5.0MHz, I need it at a rate around 4.9, and the greatest common divisor is around 3000. -- Tim Wescott Control systems and communications consulting http://www.wescottdesign.com Need to learn how to apply control theory in your embedded system? "Applied Control Theory for Embedded Systems" by Tim Wescott Elsevier/Newnes, http://www.wescottdesign.com/actfes/actfes.html
Reply by ●December 16, 20082008-12-16
On Mon, 15 Dec 2008 13:18:37 -0600, Tim Wescott wrote:> The subject line says it all. > > I need to resample some data in Scilab -- I'm getting data sampled at > 5MHz that I need to resample to 4.9 some-odd MHz. When I was using 20ms > long data sets it worked very well to do an FFT, trim bins out of the > middle to reduce the frequency by the right factor, then do an IFFT. > > Now we're doing the "real" algorithm with 1/2 second long data sets, and > it's taking way too long to do the resampling. > > I'm thinking of doing the resampling in chunks of 3000 or so (to hit the > least common denominator of the two frequencies), either using my magic > FFT method, or using a polyphase filtering method based on a multiply by > a 3000^2 sparse matrix. > > Is this the right way to go? Do the Scilab experts have any suggestions > from Scilab space? Do the DSP experts have any suggestions from DSP > space? I'd like to do this in a "Scilab kosher" way initially, although > I'm getting close to porting the whole mess over to C++. > > Thanks a bundle.FYI, for this application (spread spectrum, and bandlimited), I can accept any end effects from doing an unwindowed FFT. So I'm just concatenating rather sizable chunks, and I think than the noise at the end of each chunk will get lost in my coding gain. If not, I'll come whine some more until I get help... -- Tim Wescott Control systems and communications consulting http://www.wescottdesign.com Need to learn how to apply control theory in your embedded system? "Applied Control Theory for Embedded Systems" by Tim Wescott Elsevier/Newnes, http://www.wescottdesign.com/actfes/actfes.html
Reply by ●December 16, 20082008-12-16
On Dec 15, 8:04 pm, Tim Wescott <t...@seemywebsite.com> wrote:> On Mon, 15 Dec 2008 13:18:37 -0600, Tim Wescott wrote: > > The subject line says it all. > > > I need to resample some data in Scilab -- I'm getting data sampled at > > 5MHz that I need to resample to 4.9 some-odd MHz. When I was using 20ms > > long data sets it worked very well to do an FFT, trim bins out of the > > middle to reduce the frequency by the right factor, then do an IFFT. > > > Now we're doing the "real" algorithm with 1/2 second long data sets, and > > it's taking way too long to do the resampling. > > > I'm thinking of doing the resampling in chunks of 3000 or so (to hit the > > least common denominator of the two frequencies), either using my magic > > FFT method, or using a polyphase filtering method based on a multiply by > > a 3000^2 sparse matrix. > > > Is this the right way to go? Do the Scilab experts have any suggestions > > from Scilab space? Do the DSP experts have any suggestions from DSP > > space? I'd like to do this in a "Scilab kosher" way initially, although > > I'm getting close to porting the whole mess over to C++. > > > Thanks a bundle. > > FYI, for this application (spread spectrum, and bandlimited), I can > accept any end effects from doing an unwindowed FFT. > > So I'm just concatenating rather sizable chunks, and I think than the > noise at the end of each chunk will get lost in my coding gain. > > If not, I'll come whine some more until I get help... > > -- > Tim Wescott > Control systems and communications consultinghttp://www.wescottdesign.com > > Need to learn how to apply control theory in your embedded system? > "Applied Control Theory for Embedded Systems" by Tim Wescott > Elsevier/Newnes,http://www.wescottdesign.com/actfes/actfes.htmlA more detailed description than mine and some Matlab code is with figure 3 in: Turning overlap-save into a multiband mixing, downsampling filter bank Borgerding, M. Signal Processing Magazine, IEEE March 2006, Vol. 23, Issue: 2 On page(s): 158-161 Dale B. Dalrymple
