DSPRelated.com
Forums

Restorative Upsampling

Started by Boyboy 8 years ago13 replieslatest reply 8 years ago280 views

Hi.

A new website shares a Matlab/Octave routine that upsamples time-delimited signals while eliminating aliasing errors and restoring lost high-frequency components. Applies to time-delimited signals for which (Fsamp > Bandlimit). Potential applications are legion. http://www.contextualsignalprocessing.com

[ - ]
Reply by drmikeJuly 3, 2016

Physics makes the statement "restoring lost high-frequency components" impossible.  Lost is lost.  It's kind of like free energy and perpetual motion.  Entropy wins.

Patience, persistence, truth,

Dr. mike

[ - ]
Reply by SlartibartfastJuly 3, 2016

As has already been pointed out, there are some unsubstantiated grandiose claims on the linked site and the example code doesn't appear to do anything that would support the claims.  The solicitation for donations on the main page makes me suspicious that this is essentially click bait to get people to donate.


[ - ]
Reply by Rick LyonsJuly 3, 2016

Hello Boyboy,

Your post piqued my interest. I have two questions:

[1] With regard to a signal, can you give us a precise definition for your word "Bandlimit"?

[2] In your code, the commands 

'n = rows(y0);' and 'wc = 2*[wb;zeros((N-1),1)];'

are illegal commands in MATLAB. What are the purposes of those commands?

[ - ]
Reply by BoyboyJuly 3, 2016

Rick:

"Bandlimit" would be the highest frequency component in the original signal. If the original signal were sampled at 44.1kHz and properly low-pass filtered for non-aliasing sampling, the "Bandlimit" would be 22.05kHz. Down-sampling this signal by x2 without filtration introduces aliasing and "throws-away" high-frequency components which this algorithm seems to recover, while erasing aliasing errors, at least in large measure (as evidenced by the percentage variance [error] from the original signal cited on the Results page of the website.)

The code was written in Octave, and there may be some adaptations to be made for execution in the Matlab platform. I appreciate you taking the time to try out the code.

[ - ]
Reply by Rick LyonsJuly 3, 2016

Boyboy, I sense a possible "learning moment" here.

I'm trying to figure out EXACTLY what your code does. To do that I need to replace your illegal commands with legal MATLAB commands.

Can you tell me, what are the purposes of the commands:

'n = rows(y0);', 'wc = 2*[wb;zeros((N-1),1)];' and 'w3 = [w2a;w2c];'?

If I know what those commands are supposed to do then I can replace them with appropriate MATLAB commands.

[ - ]
Reply by BoyboyJuly 3, 2016

'n=rows(y0);' assigns to 'n' the number of rows in the row vector 'y0'

'wc=2*[wb;zeros((N-1),1)];' takes 'wb', concats it with (N-1) zeros, multiplies by scalar '2' and assigns the result to 'wc'

'w3=[w2a;w2c];' concats 'w2a' and 'w2c' and assigns the result to 'w3'

As an aside, if you cut and pasted the code from the webpage, there may be a hidden character at the end of the last line of code that needs to be deleted before the code will run in Octave, or Matlab (with your translations).

Thanks again for the time you have invested. Have a brewski on me and our sponsors.

[ - ]
Reply by BoyboyJuly 3, 2016

Rick: 

I thought I should clarify my response. In my Octave code, all non-scalars are row vectors with 1-column width. The concatenate commands result in an end to end concatenation with the result still being only 1-column in width. Octave performs its fft and ifft commands on row vectors.

Dave

[ - ]
Reply by Rick LyonsJuly 3, 2016

Hi Boyboy,

After modifying the MATLAB-illegal commands (including the illegal 

  'w2c = conj((fliplr(w2b'))');'

command) I as was able to run your code.

Your code generates a dozen, or so, different time and frequency vectors. I can't say that I understand the "overall concept" of your processing at this point.

In any case, I used a low-frequency single sinusoid, riding on a DC bias, as an 'Input' signal vector to your code. Then I generated my output 'C' vector that appears to merely be (roughly) an interpolated-by-two (with the DC bias removed) version of the 'Input' vector as shown below.

boyboy_i_2079.jpg

Boyboy, is my 'C' output correct?

[ - ]
Reply by BoyboyJuly 3, 2016

I think we have a error translating into MATLAB, and a potential misunderstanding on my part of what a row and column vector are. My output was a single 1-column vector (C) that could be readily stored as a .wav file using the Octave wavwrite command. The variance with the original wavfile was checked with the equation:

                                   rms(Signal R - Signal A)

           0.0168% =  % Error = ------------------------------- x 100%

                                        rms(Signal A)

Where Signal A is the original, and Signal R is the restorative upsampling output. This, in comparison with the Error for a conventionally Fourier upsampled (x2) Signal C:

                                rms(Signal C - Signal A)

         7.18% = % Error = ----------------------------------- x 100%

                                    rms(Signal A)

Further, I don't think a generated sinusoid will work as it is not naturally-windowed (time-delimited).

I appreciate the time you spent on this, Rick, but I am going to check out of this forum to contemplate my navel for a while.

[ - ]
Reply by Rick LyonsJuly 3, 2016

Boyboy, in MATLAB the following

x = 1 2 3 4 5

is called a "row vector" because it has one row and multiple columns. And in MATLAB the following

x = 

1

2

3

4

5

is called a "column vector" because it has one column and multiple rows.

[ - ]
Reply by Tim WescottJuly 3, 2016

I see very slick web-design skills, I see grammatically-correct use of buzzwords, I see some very questionable use of the FFT, and I see absolutely no math to back things up.

You take fft's, and you chop up the results, and you "downsample" -- but nowhere do I see you taking the result of an IFFT and taking it's real part.  Which means that you're retaining the complex-valued signal, which actually has two values for each sample.

If you pay attention to the Nyquist/Shannon sampling theorem, it says that you need to sample at an average rate of Fs >= 2B.  It does NOT say that the samples have to be identical.  You could, in theory, take some audio with a bandwidth not greater than 22kHz and grab one sample each second -- as long as the one sample is of the audio and it's first 44000 derivatives.

So I think that you may have discovered some mildly interesting property of the FFT, but not some really deep truth about the limitations of the Nyquist/Shannon sampling theorem.

[ - ]
Reply by kazJuly 3, 2016

The terminology of rows/columns of a matrix can be confusing.  I want to add these notes:

Matlab(and hence Octave) views x(r,c) as r = number of rows, c =  number of columns.

however,visually r then indicates number of samples per column and c number of samples per row.

Matlab(and Octave) functions work columnwise i.e. all samples in one column are viewed as one vector.


Kaz


[ - ]
Reply by jyaJuly 3, 2016

I learned from Einstein to value gedanken experiments. Suppose that a pure 15 KHz signal is sampled at 20 KHz. (I like to keep to round numbers.) Normally, the sampler would include an antialias filter to remove components below 10 KHz; an effective one would leave nothing at all. How could that be upsampled and reconstructed? Without an antialias filter, there would be components at 15 and 25 KHz. What would reconstruction do with that?

Jerry