Subject: Re: data resampling From: Sebastian Ewert <s.ewert@xxxxxxxx> Date: Fri, 15 Aug 2014 12:09:28 +0100 List-Archive:<http://lists.mcgill.ca/scripts/wa.exe?LIST=AUDITORY>Hi, to understand the difference between the three methods, looking at a swept sine could also be helpful. The main difference between decimate and resample in this case is the type of filtering applied before the downsampling. These filters can be adjusted. For example, an order 4 decimate can be more similar to the results of resample with default settings than to decimate with default settings. For my applications, resample with N set to 50-100 is usually fine. fs = 8000; len_sec = 3; t = [0:1/fs:len_sec]; y = sin(2*pi*( fs/(4*len_sec) * t.^2 )); lenwin = 2^nextpow2(fs * 0.04); window = hann(lenwin); [z,f,t] = spectrogram(y,window,round(lenwin*3/4),[],fs); figure; imagesc(t,f,abs(z)); title('input') fs = fs/2; lenwin = 2^nextpow2(fs * 0.04); window = hann(lenwin); [z,f,t] = spectrogram(decimate(y,2),window,round(lenwin*3/4),[],fs); figure; imagesc(t,f,abs(z)); title('decimate default') [z,f,t] = spectrogram(decimate(y,2,5),window,round(lenwin*3/4),[],fs); figure; imagesc(t,f,abs(z)); title('decimate order 5') [z,f,t] = spectrogram(resample(y,1,2),window,round(lenwin*3/4),[],fs); figure; imagesc(t,f,abs(z)); title('resample default') [z,f,t] = spectrogram(resample(y,1,2,50),window,round(lenwin*3/4),[],fs); figure; imagesc(t,f,abs(z)); title('resample length 50') [z,f,t] = spectrogram(downsample(y,2),window,round(lenwin*3/4),[],fs); figure; imagesc(t,f,abs(z)); title('downsample') Sebastian On 15/08/2014 08:26, Alain de Cheveigne wrote: > Hi Hari, > > Note that matlab's resample() has inadequate antialiasing. Try: > figure(1); clf; > x=randn(10000,1); > pwelch(decimate(x,2)); hold on > pwelch(resample(x,1,2)); > c=get(gca,'children'); set(c(1), 'color', [0 0 0]); > legend('decimate','resample'); > > decimate() does the right thing which is to apply a lowpass filter with a sufficiently low cutoff so that power is negligible at Nyquist frequency and beyond. resample()'s apparently applies a filter that attenuates by 3dB at Nyquist, so that there may be serious aliasing if the signal is wideband. > I filed a bug report to MathWorks, they answered that they think it's OK (!). Annoyingly decimate() works only for single channel data. > > It's an example of matlab's tendency to have lots of different functions for a similar purpose, each limited or buggy in a different way (see http://abandonmatlab.wordpress.com/2012/08/07/matlab-cant-read-plain-text-data-out-of-a-wet-paper-bag/ for an enlighening rant). > > Concerning the original post, the best form of downsampling depends on your requirements. If the goal is to preserve time-domain properties of the waveform, it may be best to simply replace sets of N samples by their average (where N is the downsampling factor). If the goal is to preserve frequency-domain properties (within the band below Nyquist), decimate is best. Both are similar in that each new sample is a function of a set (window) of samples of the original series. In the first case this set is a compact as makes sense (N). In the second it is much more extended, so temporal features may be blurred. In the example given, the window includes samples before and after the beginning and end of the signal (treated as zero) so the outcome is contaminated by the response to the onset and offset step. This effect is severe because of the large DC offset. > > downsample() is unlikely to be useful in any situation as it ignores all but 1/N samples of the original data. It's a rather silly function because its functionality is easy to get by writing x(1:N:end). > > > Alain > > > On 14 Aug 2014, at 15:41, Hari Bharadwaj <hari@xxxxxxxx> wrote: > >> Hi J, >> downsample(.) just picks out every nth sample (300th in your case). Though visually it may look like it worked well, it is quite likely that some of the high frequency portions of the signal are now aliased into the low frequency portions. >> >> resample(.) and decimate(.) guard against this aliasing by filtering first.. But filtering would cause artifacts at the edges of your signal as you don't have the input for infinite past (or in your FIR case, enough samples that you don't care about, for your filter transients to die out) . However the middle portions should be properly filtered and decimated. So if you throw out the first few samples (technically, the length of the filter used) and the last few sample, the rest should match reasonably with what you have in mind. This is not a problem with the software implementations but rather is a fundamental limitation of filtering and time-frequency uncertainty. If you care about those first few seconds, your input should include additional/extra samples at the beginning and end that you don't care about that would absorb the filter transients. >> >> I would recommend consulting a signal processing text book for more elaborations on these issues (I.e, aliasing, filter design, resampling etc.). >> Briefly, these might help: >> http://en.m.wikipedia.org/wiki/Anti-aliasing_filter >> https://ccrma.stanford.edu/~jos/fp/Transient_Response_Steady_State.html >> The help for resample(.) talks about what kind of filter is employed: http://www.mathworks.com/help/signal/ref/resample.html >> >> >> Hari >> >>> On Aug 13, 2014, at 7:05 PM, Juan Wu <wujuan22@xxxxxxxx> wrote: >>> >>> List experts, >>> >>> I am attempting to down sample my data from 300 Hz to 1 Hz. I just tried two functions of Matlab: 1) a = resample(Data,1,300); 2) b = downsample(Data,300). The results are quite different between each others. >>> >>> <image.png> >>> >>> Apparently, the result from "downsampled" is close to the the original data. However, quite a few persons suggested using the re-sampling method - I get this information from google searching, and very agree with this view. Personally I think the downsample method is too simple and very arbitrary. I also do not believe the "Nearest Neighbor" method. I assume that the method resampling takes both "FIR interpolator and decimator" implementation - this is what I expected. Am I right? But now my output from the re-sampling method is really very terrible. So I assume that the resample function from Matlab does not do a good job for this. I am not sure whether I need change to use other softwares or try other functions in the Matlab. >>> >>> Any opinions or references are much appreciated. >>> J >> >> The information in this e-mail is intended only for the person to whom it is >> addressed. If you believe this e-mail was sent to you in error and the e-mail >> contains patient information, please contact the Partners Compliance HelpLine at >> http://www.partners.org/complianceline . If the e-mail was sent to you in error >> but does not contain patient information, please contact the sender and properly >> dispose of the e-mail. -- Dr. Sebastian Ewert Centre for Digital Music, Queen Mary University of London s.ewert@xxxxxxxx Phone: +44 207 882 8287