FPGA Central

Hello,

I want to change the sample rate of an audio signal. Where can I find
some information about it? Some mathematical model and some source code
would be perfect.



What I want to do is:



I have a device that can play audio signal with the sampling rate of 8KSPS.
If I want to play a signal that is samples at say 44KSPS I need to change
the sampling rate. How can I do it efficiently and mathematically correct
( no aliasing)



In other hand this device can generate audio signal with the sampling rate
of 8KBPS, How can change the sampling rate to say 44KSPS?



I am working on windows so any help from windows OS can be used.



Any suggestion on where I may start will be very appreciating?



Best regards

ma wrote:
>
> Hello,
>
> I want to change the sample rate of an audio signal. Where can I find
> some information about it? Some mathematical model and some source code
> would be perfect.

Secret Rabbit Code:

http://www.mega-nerd.com/SRC/

is released under the terms of the GNU General Public License.

> I am working on windows so any help from windows OS can be used.

It compiles for Linux, MacOSX and Win32.

Erik
--
+-----------------------------------------------------------+
Erik de Castro Lopo [email protected] (Yes it's valid)
+-----------------------------------------------------------+
Linux: the only OS that makes you feel guilty when you reboot
-- Kenneth Crudup in comp.os.linux.misc

ma wrote:
> Hello,
>
> I want to change the sample rate of an audio signal. Where can I find
> some information about it? Some mathematical model and some source code
> would be perfect.

There is an excellent tutorial at http://www.dspdimension.com/ that
supplements Eric's, as well as some example (working) code.

Jerry
--
Engineering is the art of making what you want from things you can get.
¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯ ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯

Thanks,

The tutorials in this page are very interesting and informative. I
enjoyed reading them but I couldn't find any tutorial on sampling rate
changing. Did I miss something?



Best regards



"Jerry Avins" <[email protected]> wrote in message
news:[email protected]..
> ma wrote:
>> Hello,
>>
>> I want to change the sample rate of an audio signal. Where can I
>> find some information about it? Some mathematical model and some source
>> code would be perfect.
>
> There is an excellent tutorial at http://www.dspdimension.com/ that
> supplements Eric's, as well as some example (working) code.
>
> Jerry
> --
> Engineering is the art of making what you want from things you can get.
> ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯ ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯

Thanks,

Interesting code. But I couldn't find any mathematical model on how they are
doing this. Why they need a whole file to do this? Why small chunk of data
generate some effect?

Where can I find some mathematical modeling? I need to know the theory and
then use the code as a learning practice.



Best regards





"Erik de Castro Lopo" <[email protected]> wrote in message
news:[email protected]..
> ma wrote:
>>
>> Hello,
>>
>> I want to change the sample rate of an audio signal. Where can I
>> find
>> some information about it? Some mathematical model and some source code
>> would be perfect.
>
> Secret Rabbit Code:
>
> http://www.mega-nerd.com/SRC/
>
> is released under the terms of the GNU General Public License.
>
>> I am working on windows so any help from windows OS can be used.
>
> It compiles for Linux, MacOSX and Win32.
>
> Erik
> --
> +-----------------------------------------------------------+
> Erik de Castro Lopo [email protected] (Yes it's valid)
> +-----------------------------------------------------------+
> Linux: the only OS that makes you feel guilty when you reboot
> -- Kenneth Crudup in comp.os.linux.misc

ma wrote:

(top-posting fixed)
> "Erik de Castro Lopo" <[email protected]> wrote in message
> news:[email protected]..
>
>>ma wrote:
>>
>>>Hello,
>>>
>>> I want to change the sample rate of an audio signal. Where can I
>>>find
>>>some information about it? Some mathematical model and some source code
>>>would be perfect.
>>
>>Secret Rabbit Code:
>>
>> http://www.mega-nerd.com/SRC/
>>
>>is released under the terms of the GNU General Public License.
>>
>>
>>>I am working on windows so any help from windows OS can be used.
>>
>>It compiles for Linux, MacOSX and Win32.
>>
>>Erik
>>--
>>+-----------------------------------------------------------+
>> Erik de Castro Lopo [email protected] (Yes it's valid)
>>+-----------------------------------------------------------+
>>Linux: the only OS that makes you feel guilty when you reboot
>> -- Kenneth Crudup in comp.os.linux.misc
>
>
>
> Thanks,
>
> Interesting code. But I couldn't find any mathematical model on how
> they are doing this. Why they need a whole file to do this? Why small
> chunk of data generate some effect?
>
> Where can I find some mathematical modeling? I need to know the theory
> and then use the code as a learning practice.
>

The underlying technique is called polyphase filtering. I just tried to
type up a synopsis and failed -- search for it on the web & you'll find
some good description, no doubt.

--

Tim Wescott
Wescott Design Services
http://www.wescottdesign.com

For some theory, start by reading the notes on analog devices on their hardware
sample rate converters. Here are a few links:
http://www.analog.com/UploadedFiles/...52667AN394.pdf
http://www.analog.com/en/content/0,2...5F9690,00.html

This focuses on real-time hardware implementations, but I find it to be a very
readable introduction.

--
Jon Harris
SPAM blocker in place:
Remove 99 (but leave 7) to reply

"ma" <[email protected]> wrote in message
news:VVLIe.108880$[email protected] .uk...
> Thanks,
>
> Interesting code. But I couldn't find any mathematical model on how they are
> doing this. Why they need a whole file to do this? Why small chunk of data
> generate some effect?
>
> Where can I find some mathematical modeling? I need to know the theory and
> then use the code as a learning practice.
>
>
>
> Best regards
>
>
>
> "Erik de Castro Lopo" <[email protected]> wrote in message
> news:[email protected]..
>> ma wrote:
>>>
>>> Hello,
>>>
>>> I want to change the sample rate of an audio signal. Where can I find
>>> some information about it? Some mathematical model and some source code
>>> would be perfect.
>>
>> Secret Rabbit Code:
>>
>> http://www.mega-nerd.com/SRC/
>>
>> is released under the terms of the GNU General Public License.
>>
>>> I am working on windows so any help from windows OS can be used.
>>
>> It compiles for Linux, MacOSX and Win32.
>>
>> Erik

in article VVLIe.108880$[email protected], ma at
[email protected] wrote on 08/05/2005 11:51:

> Where can I find some mathematical modeling? I need to know the theory and
> then use the code as a learning practice.

you want theory? i'll give you theory.

below is about as fundamental for the theory as you can get. it really is
just a result or application of the Nyquist/Shannon Sampling and
Reconstruction Theorem. imagine outputting your signal, at the original
sampling rate, to an ideal D/A converter and then resampling that analog
output signal with an ideal A/D converter (running at the new sampling
rate). now convert that imagination to mathematics (see below).

i could have just referred to the post (responding to the infamous
"Airy-head") but Google is messing up the ASCII spacing and you don't wanna
read what Airy had to say, anyway. so it's reposted below.

--

r b-j [email protected]

"Imagination is more important than knowledge."





(a monospaced font needed to view the following correctly)


__________________________________________________ __________________________

* * * * Bandlimited Interpolation, Sample Rate Conversion, or
* * * * * * Fractional-Sample Delay from the P.O.V. of the
* * * * * Nyquist/Shannon Sampling and Reconstruction Theorem

__________________________________________________ __________________________

Here is the mathematical expression of the sampling and reconstruction
theorem:



* * * * * * * * * x(t)*q(t) = T*SUM{x[k]*d(t-k*T)} * .------.
* * x(t)--->(*)------------------------------------->| H(f) |---> x(t)
* * * * * * *^ * * * * * * * * * * * * * * * * * * *'------'
* * * * * * *|
* * * * * * *| * * * * * * * *+inf
* * * * * * *'------- q(t) = T * SUM{ d(t - k*T) }
* * * * * * * * * * * * * * *k=-inf


where: d(t) = 'dirac' impulse function
and T = 1/Fs = sampling period
* * * * Fs = sampling frequency


* * * * * * * *+inf
* * q(t) = T * SUM{ d(t - k*T) } *is the "sampling function", is periodic
* * * * * * * k=-inf * * * * * * *with period T, and can be expressed as a
* * * * * * * * * * * * * * * * * Fourier series. *It turns out that ALL of
* * * * * * * * * * * * * * * * * the Fourier coefficients are equal to 1.


* * * * * +inf
* * q(t) = SUM{ c[n]*exp(j*2*pi*n/T*t) }
* * * * * n=-inf


where c[n] is the Fourier series coefficient:


* * * * * * * * * * t0+T
* * c[n] = *1/T * integral{q(t) * exp(-j*2*pi*n/T*t) dt}
* * * * * * * * * * *t0


and t0 can be any time. *We choose t0 to be -T/2.


* * * * * * * * * * *T/2 * * *+inf
* * c[n] = *1/T * integral{T * SUM{ d(t - k*T) } * exp(-j*2*pi*n/T*t) dt}
* * * * * * * * * * -T/2 * * k=-inf


* * * * * * * T/2
* * * * *= integral{ d(t - 0*T) * exp(-j*2*pi*n/T*t) dt}
* * * * * * *-T/2


* * * * * * * T/2
* * * * *= integral{ d(t) * exp(-j*2*pi*n/T*t) dt}
* * * * * * *-T/2


* * * * *= exp(-j*2*pi*n/T*0)


= 1


So all Fourier series coefficients for the periodic q(t) are 1.
This results in another valid expression for q(t):


* * * * * +inf
* * q(t) = SUM{ exp(j*2*pi*n*Fs*t) } .
* * * * *n=-inf


The sampled signal looks like:


* * * * * * * * * * * * * +inf
* * x(t)*q(t) = x(t) * T * SUM{ d(t - k*T) }
* * * * * * * * * * * * *k=-inf


* * * * * * * * * *+inf
* * * * * * * = T * SUM{ x(t) * d(t - k*T) }
* * * * * * * * * k=-inf


* * * * * * * * * *+inf
* * * * * * * = T * SUM{ x(k*T) * d(t - k*T) }
* * * * * * * * * k=-inf


* * * * * * * * * *+inf
* * * * * * * = T * SUM{ x[k] * d(t - k*T) }
* * * * * * * * * k=-inf


where x[k] =def x(k*T) by definition. (We are using the notation
that x[] is a discrete sequence and x() is a continuous function.)
BTW, if you where to compute the Laplace Transform of the sampled
signal x(t)*q(t), you would get the Z Transform of x[k] scaled by
the constant T.


An alternative (and equally important representation) of
the sampled signal is:


* * * * * * * * * * * +inf
* * x(t)*q(t) = x(t) * SUM{ exp(j*2*pi*n*Fs*t) }
* * * * * * * * * * *n=-inf


* * * * * * * *+inf
* * * * * * * = SUM{ x(t) * exp(j*2*pi*n*Fs*t) }
* * * * * * * n=-inf

Defining the spectrum of x(t) as its Fourier Transform:


* * * * * * * * * * * * *+inf
* * X(f) =def FT{ x(t) } =def integral{ x(t)*exp(-j*2*pi*f*t) dt}
* * * * * * * * * * * * * -inf


Then using the frequency shifting property of the Fourier Transform,


* * FT{ x(t)*exp(j*2*pi*f0*t) } = X(f - f0)


results in


* * * * * * * * * * *+inf
* * FT{ x(t)*q(t) } = SUM{ X(f - n*Fs) }
* * * * * * * * * * n=-inf


This says, what we all know, that the spectrum of our signal
being sampled is shifted and repeated forever at multiples
of the sampling frequency, Fs. *If x(t) or X(f) is bandlimited
to B (i.e. X(f) = 0 for all |f| > B) AND if there is no
overlap of the tails of adjacent images X(f), that is:


* right tail of nth image of X(f) < left tail of (n+1)th image


* n*Fs + B < (n+1)*Fs - B = n*Fs + Fs - B


* * * * * *B < Fs - B


* * * * *2*B < Fs = 1/T


Then we ought to be able to reconstruct X(f) (and also x(t)) by
low pass filtering out all of the images of X(f). *To do that,
Fs > 2*B (to prevent overlap) and H(f) must be:


* * * * * *{ *1 *for *|f| < *Fs/2 = 1/(2*T)
* * H(f) = {
* * * * * *{ *0 *for *|f| >= Fs/2 = 1/(2*T)


This is the "sampling" half of the Nyquist/Shannon sampling and
reconstruction theorem. *It says that the sampling frequency, Fs,
must be strictly greater than twice the bandwidth, B, of the
continuous-time signal, x(t), for no information to be lost (or
"aliased"). *The "reconstruction" half follows:


The impulse response of the reconstruction LPF, H(f), is the
inverse Fourier Transform of H(f), called h(t):


* * * * * * * * * * * * * *+inf
* * h(t) = iFT{ H(f) } = integral{ H(f)*exp(+j*2*pi*f*t) df}
* * * * * * * * * * * * * *-inf


* * * * * *1/(2*T)
* * * * = integral{ exp(+j*2*pi*f*t) df}
* * * * * -1/(2*T)


* * * * = (exp(+j*2*pi/(2*T)*t) - exp(+j*2*pi/(-2*T)*t))/(j*2*pi *t)


* * * * = ( exp(+j*pi/T*t) - exp(-j*pi/T*t) )/(j*2*pi*t)


* * * * = sin((pi/T)*t)/(pi*t)


* * * * = (1/T)*sin(pi*t/T)/(pi*t/T)


* * * * = (1/T)*sinc(t/T)


where sinc(u) =def sin(pi*u)/(pi*u)


The input to the LPF is x(t)*q(t) = T*SUM{x[k]*d(t - k*T)} .
Each d(t - k*T) impulse generates its own impulse response and
since the LPF is linear and time invariant, all we have to do
is add up the time-shifted impulse responses weighted by their
coefficients, x[k]. *The T and 1/T kill each other off.

The output of the LPF is


* * * * * +inf * * * * * * * * * * * * * +inf
* * x(t) = SUM{ x[k]*sinc((t-k*T)/T) } = SUM{ x[k]*sinc(t/T - k) }
* * * * *k=-inf * * * * * * * * * * * * k=-inf


This equation tells us explicitly how to reconstruct our
sampled, bandlimited input signal from the samples. *This is
the "reconstruction" half of the Nyquist/Shannon sampling and
reconstruction theorem.



Interpolation, Fractional Delay, Sample Rate Conversion, etc:

When interpolating (whether for fractional sample delay or
for doing sample rate conversion), you must evaluate this
equation for times that might not be integer multiples of your
sampling period, T.

The sinc(t/T) function is one for t = 0 and zero for t = k*T
for k = nonzero integer. *This means that if your new
sample time, t, happens to land exactly on an old sample
time, k*T, only that sample (and none of the neighbors)
contributes to the output sample and the output sample
is equal to the input sample. *Only in the case where the
output sample is in between input samples, do the neighbors
contribute to the calculation.

Since the sinc() function goes on forever to +/- infinity,
it must be truncated somewhere to be of practical use.
Truncating is actually applying the rectangular window (the
worst kind) so it is advantageous to window down the sinc()
function gradually using something like a Hamming or Kaiser
window. *In my experience, with a good window, you'll need to
keep the domain of the windowed sinc() function from -16 to +16
(a 32 tap FIR) and sample it 16384 times in that region (that
would be 16384/32 = 512 polyphases or "fractional delays" if
linear interpolation is subsequently used to get in between
neighboring polyphases).

Now we simply split the time, t, into an integer part, m, and
fractional part, a, and substitute:


* * * * * * * * +inf
* * x((m+a)*T) = SUM{ x[k]*sinc((m+a)*T/T - k) }
* * * * * * * * k=-inf


* * * * * * * * +inf
* * * * * * * *= SUM{ x[k]*sinc(a + m - k) }
* * * * * * * * k=-inf


where: t = (m + a)*T
* and m = floor(t/T)
* a = fract(t/T) = t/T - floor(t/T)

t-T < m*T <= t *which means a is always nonnegative and less
than 1: * 0 <= a < 1


The limits of the sum can be bumped around a bit without changing
it and by substituting k <- k+m :


* * * * * * * * +inf
* * x((m+a)*T) = SUM{ x[m+k] * sinc(a - k) }
* * * * * * * * k=-inf


Now, an approximation is made by truncating the summation to a
finite amount (which is effectively windowing):


* * * * * * * * *+16
* * x((m+a)*T) = SUM{ x[m+k] * sinc(a-k) * w((a-k)/16) }
* * * * * * * * k=-15


w(u) is the window function. *A Hamming window would be:

* * * * * *{ 0.54 + 0.46*cos(pi*u) * for |u| < 1
* * w(u) = {
* * * * * *{ 0 * * * * * * * * * * * for |u| >= 1

A Kaiser window would be better.

This requires two 32 point FIR computations to calculate two
adjacent interpolated samples (these are linearly interpolated
to get the final output sample). *Since the sinc(u)*w(u)
function is symmetrical, that means 8193 numbers stored somewhere
in memory. *When interpolating, the integer part of t/T (or 'm')
determines which 32 adjacent samples of x[k] to use, the
fractional part of t/T (or 'a') determines which set of 32 windowed
sinc() coefficients to be used to combine the 32 samples.

__________________________________________________ __________________________

__________________________________________________ __________________________

ma wrote:
> Thanks,
>
> The tutorials in this page are very interesting and informative. I
> enjoyed reading them but I couldn't find any tutorial on sampling rate
> changing. Did I miss something?

...

There's a close connection between pitch changing ans sample-rate
changing. I hoped you would see how they fit.

Jerry
--
Engineering is the art of making what you want from things you can get.
¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯ ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯

ma wrote:
> I want to change the sample rate of an audio signal. Where can I find
> some information about it? Some mathematical model and some source code
> would be perfect.

It's all about filtered interpolation. Most modern methods are based
on windowed-Sinc interpolation. One of the fundamental papers on this
subject is here:

http://ccrma-www.stanford.edu/~jos/resample/

The various polyphase methods are just an optimization. On a fast
PC, you can calculate exact Sinc filter coefficients on the fly for
some simpler types of windows.


IMHO. YMMV.
--
rhn A.T nicholson d.O.t C-o-M

ma wrote:
>
> Thanks,
>
> Interesting code. But I couldn't find any mathematical model on how they are
> doing this. Why they need a whole file to do this? Why small chunk of data
> generate some effect?
>
> Where can I find some mathematical modeling?

The Secret Rabbit Code web pages do state that the rabbit is a re-implementation
of the Julis O. Smith algorithm. Is is described here:

http://ccrma-www.stanford.edu/~jos/resample/


Erik
--
+-----------------------------------------------------------+
Erik de Castro Lopo [email protected] (Yes it's valid)
+-----------------------------------------------------------+
" ... new TV ad for Microsoft's Internet Explorer e-mail program which
uses the musical theme of the "Confutatis Maledictis" from Mozart's
Requiem. "Where do you want to go today?" is the cheery line on the
screen, while the chorus sings "Confutatis maledictis, flammis acribus
addictis,". This translates to "The damned and accursed are convicted
to the flames of hell."

The most understandable conceptual explanation I've come across on the
web is here:

http://dspguru.com/info/faqs/mrfaq.htm

In a nutshell, changing the sampling rate by an integer factor of N:

Interpolation:
1. Upsample: Insert N-1 0's between samples
2. Low-pass filter below original fs/2

Decimation:
1. Low-pass filter below desired fs/2
2. Downsample: Skip samples (keep only every Nth sample, and disregard
N-1 samples in between)

ma wrote:
> Hello,
>
> I want to change the sample rate of an audio signal. Where can I find
> some information about it? Some mathematical model and some source code
> would be perfect.

Tim Wescott wrote:
>
> The underlying technique is called polyphase filtering.

Actually the polyphase filtering usually means a fixed number of
phases. The Julius O. SMith technique can acutally handle irrational
ratios between input and output sample rate and therefore has an
infinite number of phases.

Erik
--
+-----------------------------------------------------------+
Erik de Castro Lopo [email protected] (Yes it's valid)
+-----------------------------------------------------------+
Traditional capital was stuck in a company's bank account or investments.
It could not walk away in disgust. Human capital has free will. It can
walk out the door; traditional capital cannot.

in article [email protected], Erik de Castro Lopo at
[email protected] wrote on 08/05/2005 19:12:

> Tim Wescott wrote:
>>
>> The underlying technique is called polyphase filtering.
>
> Actually the polyphase filtering usually means a fixed number of
> phases. The Julius O. Smith technique can acutally handle irrational
> ratios between input and output sample rate and therefore has an
> infinite number of phases.

but it still is not inaccurate to use the term "polyphase filtering" for
naming the underlying technique. so Smith and Gossett linearly interpolate
(or use some other low-order spline) between these polyphases to resample at
any arbitrary time - it's essentially a way to get a more accurate
representation of that windowed sinc or whatever the LPF impulse response
is.

i *really* have a lot of respect for JOS (and consider him to be a sorta
friend, but i only see him at the AES conventions), but i think the idea to
interpolate between the discrete and finite polyphases is incremental. i am
surprized at how many folks refer to Smith & Gossett as the seminal paper on
this. Rabiner, in the 70's, might be the seminal author.

--

r b-j [email protected]

"Imagination is more important than knowledge."

Thanks to everyone that replied.

It really helped me to refresh my signal processing knowledge and also add a
lot to it.

After reading most of the documents that was suggested as a good reading, I
grasp some information about the technique but I find new questions:



1- In Rabbit source code, it is mentioned that if I use small chunk of data
for transformation then I will have odd sound effects. Why is this
happening? I need to change the sampling rate of an online system which
means I only have access to some limited samples (small chunk of data). How
can I make the system not to generate odd sound effects?



2- What is the delay in this process? Since I want to change the sampling
rate online, if the delay is very high, I may get some echo effect or
similar problem. Any idea of how to solve this problem?



3- Are these algorithms fast enough to do processing on real-time on a PC?



Best regards



"ma" <[email protected]> wrote in message
news:URGIe.104715$[email protected] .uk...
> Hello,
>
> I want to change the sample rate of an audio signal. Where can I find
> some information about it? Some mathematical model and some source code
> would be perfect.
>
>
>
> What I want to do is:
>
>
>
> I have a device that can play audio signal with the sampling rate of
> 8KSPS. If I want to play a signal that is samples at say 44KSPS I need to
> change the sampling rate. How can I do it efficiently and mathematically
> correct ( no aliasing)
>
>
>
> In other hand this device can generate audio signal with the sampling rate
> of 8KBPS, How can change the sampling rate to say 44KSPS?
>
>
>
> I am working on windows so any help from windows OS can be used.
>
>
>
> Any suggestion on where I may start will be very appreciating?
>
>
>
> Best regards
>
>
>
>

"ma" <[email protected]> schrieb im Newsbeitrag
news:URGIe.104715$[email protected] .uk...
> Hello,
>
> I want to change the sample rate of an audio signal. Where can I find
> some information about it? Some mathematical model and some source code
> would be perfect.
>
>
>
> What I want to do is:
>
>
>
> I have a device that can play audio signal with the sampling rate of
8KSPS.
> If I want to play a signal that is samples at say 44KSPS I need to change
> the sampling rate. How can I do it efficiently and mathematically correct
> ( no aliasing)
>
>
>
> In other hand this device can generate audio signal with the sampling rate
> of 8KBPS, How can change the sampling rate to say 44KSPS?
>
>
>
> I am working on windows so any help from windows OS can be used.
>
>
>
> Any suggestion on where I may start will be very appreciating?
>
>
>
> Best regards
>
>
>
>

In a DSP Book there was the following receipe (i think)

Do a DFT of 44100 input Samples, store them in an array of 48000(x2)
(Complex)Samples, fill the remaining (44101 to 48000) bins with zeroes, make
an inverse DFT of this 48000 bins and you get 48000 output samples.

In article <gh%Ie.118455$[email protected] >,
ma <[email protected]> wrote:
>1- ... it is mentioned that if I use small chunk of data
>for transformation then I will have odd sound effects. Why is this
>happening? I need to change the sampling rate of an online system which
>means I only have access to some limited samples (small chunk of data). How
>can I make the system not to generate odd sound effects?

Typical poly-phase resampling techniques only mention linear phase and
symmetric interpolation filters. You can also try minimum phase filters
or even one-sided interpolators near the sample buffer boundries.

>2- What is the delay in this process? Since I want to change the sampling
>rate online, if the delay is very high, I may get some echo effect or
>similar problem. Any idea of how to solve this problem?

The delay depends on the width and symmetry of your interpolator, you
can always use a shorter or non-symmetric (minimum phase or one-sided)
interpolator, trading off against your various quality metrics.

>3- Are these algorithms fast enough to do processing on real-time on a PC?

Depends on the algorithm and the PC. I've wrote code which recalculated
each raised-cosine windowed sinc coefficient for dozens of taps for
real-time audio non-rational resampling just using standard sin/cos math
library calls on a current PowerMac. With various precalculated table
look-up optimizations, you can of course get more performance if needed.


IMHO. YMMV.
--
Ron Nicholson rhn AT nicholson DOT com http://www.nicholson.com/rhn/
#include <canonical.disclaimer> // only my own opinions, etc.

1. Why do you represent yourself as an escapee from the
school playground?

Airy has not been around for months and yet you are
repeatedly making rather silly and infantile remarks
in his direction.

If you don't understand the points that Airy makes, then
admit it, or better to keep quiet. As it is, by dealing with
technical points in an international public forum such as this
one by just bad-mouthing your opponents makes you appear to
be a charlatan peddling some snake-oil elixir.

2. Talking of snake-oil, where have your Diracian impulses come from?

The Diracian has some interesting properties - zero width,
area of unity, infinite sum of all possible cosines, a height
which is not discussed but which appears to be greater in magnitude
than a Bristow-Johnson ego (if that were indeed possible). In the
systems that you deal with, what experimental evidence do you
have that there are pulses with the attributes of Diracians upon
which to base your theory?

3. Talking of snake-oil, where did the factor of "T" come from in your
opening lines?

Consider a 16-bit ADC capable of 100 M Samples per sec. In the first
instance we'll use it to sample a geophysical signal of bandwidth
limited
to 300 HZ and sample at 1 kHz, with suitable analogue instrumentation
to match the input signal to the full range of the ADC. If we now
keep the circuit the same, but now sample at 65.536 MHZ, your claimed
factor of "T" will result in the 16 bit range being compressed down
to one bit. We know this doesn't happen -there will be more samples,
but they'll still be of the same magnitude and 16-bit range.

There

robert bristow-johnson wrote:
> in article VVLIe.108880$[email protected], ma at
> [email protected] wrote on 08/05/2005 11:51:
>
> > Where can I find some mathematical modeling? I need to know the theory and
> > then use the code as a learning practice.
>
> you want theory? i'll give you theory.
>
> below is about as fundamental for the theory as you can get. it really is
> just a result or application of the Nyquist/Shannon Sampling and
> Reconstruction Theorem. imagine outputting your signal, at the original
> sampling rate, to an ideal D/A converter and then resampling that analog
> output signal with an ideal A/D converter (running at the new sampling
> rate). now convert that imagination to mathematics (see below).
>
> i could have just referred to the post (responding to the infamous
> "Airy-head") but Google is messing up the ASCII spacing and you don't wanna
> read what Airy had to say, anyway. so it's reposted below.
>
> --
>
> r b-j [email protected]
>
> "Imagination is more important than knowledge."
>
>
>
>
>
> (a monospaced font needed to view the following correctly)
>
>
> __________________________________________________ __________________________
>
> Bandlimited Interpolation, Sample Rate Conversion, or
> Fractional-Sample Delay from the P.O.V. of the
> Nyquist/Shannon Sampling and Reconstruction Theorem
>
> __________________________________________________ __________________________
>
> Here is the mathematical expression of the sampling and reconstruction
> theorem:
>
>
>
> x(t)*q(t) = T*SUM{x[k]*d(t-k*T)} .------.
> x(t)--->(*)------------------------------------->| H(f) |---> x(t)
> ^ '------'
> |
> | +inf
> '------- q(t) = T * SUM{ d(t - k*T) }
> k=-inf
>
>
> where: d(t) = 'dirac' impulse function
> and T = 1/Fs = sampling period
> Fs = sampling frequency
>
>
> +inf
> q(t) = T * SUM{ d(t - k*T) } is the "sampling function", is periodic
> k=-inf with period T, and can be expressed as a
> Fourier series. It turns out that ALL of
> the Fourier coefficients are equal to 1.
>
>
> +inf
> q(t) = SUM{ c[n]*exp(j*2*pi*n/T*t) }
> n=-inf
>
>
> where c[n] is the Fourier series coefficient:
>
>
> t0+T
> c[n] = 1/T * integral{q(t) * exp(-j*2*pi*n/T*t) dt}
> t0
>
>
> and t0 can be any time. We choose t0 to be -T/2.
>
>
> T/2 +inf
> c[n] = 1/T * integral{T * SUM{ d(t - k*T) } * exp(-j*2*pi*n/T*t) dt}
> -T/2 k=-inf
>
>
> T/2
> = integral{ d(t - 0*T) * exp(-j*2*pi*n/T*t) dt}
> -T/2
>
>
> T/2
> = integral{ d(t) * exp(-j*2*pi*n/T*t) dt}
> -T/2
>
>
> = exp(-j*2*pi*n/T*0)
>
>
> = 1
>
>
> So all Fourier series coefficients for the periodic q(t) are 1.
> This results in another valid expression for q(t):
>
>
> +inf
> q(t) = SUM{ exp(j*2*pi*n*Fs*t) } .
> n=-inf
>
>
> The sampled signal looks like:
>
>
> +inf
> x(t)*q(t) = x(t) * T * SUM{ d(t - k*T) }
> k=-inf
>
>
> +inf
> = T * SUM{ x(t) * d(t - k*T) }
> k=-inf
>
>
> +inf
> = T * SUM{ x(k*T) * d(t - k*T) }
> k=-inf
>
>
> +inf
> = T * SUM{ x[k] * d(t - k*T) }
> k=-inf
>
>
> where x[k] =def x(k*T) by definition. (We are using the notation
> that x[] is a discrete sequence and x() is a continuous function.)
> BTW, if you where to compute the Laplace Transform of the sampled
> signal x(t)*q(t), you would get the Z Transform of x[k] scaled by
> the constant T.
>
>
> An alternative (and equally important representation) of
> the sampled signal is:
>
>
> +inf
> x(t)*q(t) = x(t) * SUM{ exp(j*2*pi*n*Fs*t) }
> n=-inf
>
>
> +inf
> = SUM{ x(t) * exp(j*2*pi*n*Fs*t) }
> n=-inf
>
> Defining the spectrum of x(t) as its Fourier Transform:
>
>
> +inf
> X(f) =def FT{ x(t) } =def integral{ x(t)*exp(-j*2*pi*f*t) dt}
> -inf
>
>
> Then using the frequency shifting property of the Fourier Transform,
>
>
> FT{ x(t)*exp(j*2*pi*f0*t) } = X(f - f0)
>
>
> results in
>
>
> +inf
> FT{ x(t)*q(t) } = SUM{ X(f - n*Fs) }
> n=-inf
>
>
> This says, what we all know, that the spectrum of our signal
> being sampled is shifted and repeated forever at multiples
> of the sampling frequency, Fs. If x(t) or X(f) is bandlimited
> to B (i.e. X(f) = 0 for all |f| > B) AND if there is no
> overlap of the tails of adjacent images X(f), that is:
>
>
> right tail of nth image of X(f) < left tail of (n+1)th image
>
>
> n*Fs + B < (n+1)*Fs - B = n*Fs + Fs - B
>
>
> B < Fs - B
>
>
> 2*B < Fs = 1/T
>
>
> Then we ought to be able to reconstruct X(f) (and also x(t)) by
> low pass filtering out all of the images of X(f). To do that,
> Fs > 2*B (to prevent overlap) and H(f) must be:
>
>
> { 1 for |f| < Fs/2 = 1/(2*T)
> H(f) = {
> { 0 for |f| >= Fs/2 = 1/(2*T)
>
>
> This is the "sampling" half of the Nyquist/Shannon sampling and
> reconstruction theorem. It says that the sampling frequency, Fs,
> must be strictly greater than twice the bandwidth, B, of the
> continuous-time signal, x(t), for no information to be lost (or
> "aliased"). The "reconstruction" half follows:
>
>
> The impulse response of the reconstruction LPF, H(f), is the
> inverse Fourier Transform of H(f), called h(t):
>
>
> +inf
> h(t) = iFT{ H(f) } = integral{ H(f)*exp(+j*2*pi*f*t) df}
> -inf
>
>
> 1/(2*T)
> = integral{ exp(+j*2*pi*f*t) df}
> -1/(2*T)
>
>
> = (exp(+j*2*pi/(2*T)*t) - exp(+j*2*pi/(-2*T)*t))/(j*2*pi *t)
>
>
> = ( exp(+j*pi/T*t) - exp(-j*pi/T*t) )/(j*2*pi*t)
>
>
> = sin((pi/T)*t)/(pi*t)
>
>
> = (1/T)*sin(pi*t/T)/(pi*t/T)
>
>
> = (1/T)*sinc(t/T)
>
>
> where sinc(u) =def sin(pi*u)/(pi*u)
>
>
> The input to the LPF is x(t)*q(t) = T*SUM{x[k]*d(t - k*T)} .
> Each d(t - k*T) impulse generates its own impulse response and
> since the LPF is linear and time invariant, all we have to do
> is add up the time-shifted impulse responses weighted by their
> coefficients, x[k]. The T and 1/T kill each other off.
>
> The output of the LPF is
>
>
> +inf +inf
> x(t) = SUM{ x[k]*sinc((t-k*T)/T) } = SUM{ x[k]*sinc(t/T - k) }
> k=-inf k=-inf
>
>
> This equation tells us explicitly how to reconstruct our
> sampled, bandlimited input signal from the samples. This is
> the "reconstruction" half of the Nyquist/Shannon sampling and
> reconstruction theorem.
>
>
>
> Interpolation, Fractional Delay, Sample Rate Conversion, etc:
>
> When interpolating (whether for fractional sample delay or
> for doing sample rate conversion), you must evaluate this
> equation for times that might not be integer multiples of your
> sampling period, T.
>
> The sinc(t/T) function is one for t = 0 and zero for t = k*T
> for k = nonzero integer. This means that if your new
> sample time, t, happens to land exactly on an old sample
> time, k*T, only that sample (and none of the neighbors)
> contributes to the output sample and the output sample
> is equal to the input sample. Only in the case where the
> output sample is in between input samples, do the neighbors
> contribute to the calculation.
>
> Since the sinc() function goes on forever to +/- infinity,
> it must be truncated somewhere to be of practical use.
> Truncating is actually applying the rectangular window (the
> worst kind) so it is advantageous to window down the sinc()
> function gradually using something like a Hamming or Kaiser
> window. In my experience, with a good window, you'll need to
> keep the domain of the windowed sinc() function from -16 to +16
> (a 32 tap FIR) and sample it 16384 times in that region (that
> would be 16384/32 = 512 polyphases or "fractional delays" if
> linear interpolation is subsequently used to get in between
> neighboring polyphases).
>
> Now we simply split the time, t, into an integer part, m, and
> fractional part, a, and substitute:
>
>
> +inf
> x((m+a)*T) = SUM{ x[k]*sinc((m+a)*T/T - k) }
> k=-inf
>
>
> +inf
> = SUM{ x[k]*sinc(a + m - k) }
> k=-inf
>
>
> where: t = (m + a)*T
> and m = floor(t/T)
> a = fract(t/T) = t/T - floor(t/T)
>
> t-T < m*T <= t which means a is always nonnegative and less
> than 1: 0 <= a < 1
>
>
> The limits of the sum can be bumped around a bit without changing
> it and by substituting k <- k+m :
>
>
> +inf
> x((m+a)*T) = SUM{ x[m+k] * sinc(a - k) }
> k=-inf
>
>
> Now, an approximation is made by truncating the summation to a
> finite amount (which is effectively windowing):
>
>
> +16
> x((m+a)*T) = SUM{ x[m+k] * sinc(a-k) * w((a-k)/16) }
> k=-15
>
>
> w(u) is the window function. A Hamming window would be:
>
> { 0.54 + 0.46*cos(pi*u) for |u| < 1
> w(u) = {
> { 0 for |u| >= 1
>
> A Kaiser window would be better.
>
> This requires two 32 point FIR computations to calculate two
> adjacent interpolated samples (these are linearly interpolated
> to get the final output sample). Since the sinc(u)*w(u)
> function is symmetrical, that means 8193 numbers stored somewhere
> in memory. When interpolating, the integer part of t/T (or 'm')
> determines which 32 adjacent samples of x[k] to use, the
> fractional part of t/T (or 'a') determines which set of 32 windowed
> sinc() coefficients to be used to combine the 32 samples.
>
> __________________________________________________ __________________________
>
> __________________________________________________ __________________________

A good question indeed!

And a question that the so-called experts and authors in this
field fail to answer time and time again.

The Diracian, or Unit Impulse is a very good mathematical tool
to analyse the response of systems once a mathmatical model of
those systems had been produced.

It is, however, a poor mathematical claim to make that such
impulses are found to be part of a system when neither the
area nor the magnitude of such impulses are found anywhere in such
systems.

To those who ask, "Who cares? I get good results." , I suggest that
their approach is unscientific and compares to the religious
loonies who sacrifice goats and virgins to stop the Sun falling
out of the sky and justify the continuing practice by the Sun
remaining in the sky.

So.....is the world of DSP a world of scientific men, or is
it a world of snake-oil charlatans and of religious loonies?

Where do these Diracian impulses come from?

Airy R.Bean wrote:

> 2. Talking of snake-oil, where have your Diracian impulses come from?
>
> The Diracian has some interesting properties - zero width,
> area of unity, infinite sum of all possible cosines, a height
> which is not discussed but which appears to be greater in magnitude
> than a Bristow-Johnson ego (if that were indeed possible). In the
> systems that you deal with, what experimental evidence do you
> have that there are pulses with the attributes of Diracians upon
> which to base your theory?
>

Polymath wrote:
> A good question indeed!
>
> And a question that the so-called experts and authors in this
> field fail to answer time and time again.
>
> The Diracian, or Unit Impulse is a very good mathematical tool
> to analyse the response of systems once a mathmatical model of
> those systems had been produced.
>
> It is, however, a poor mathematical claim to make that such
> impulses are found to be part of a system when neither the
> area nor the magnitude of such impulses are found anywhere in such
> systems.
>
> To those who ask, "Who cares? I get good results." , I suggest that
> their approach is unscientific and compares to the religious
> loonies who sacrifice goats and virgins to stop the Sun falling
> out of the sky and justify the continuing practice by the Sun
> remaining in the sky.
>
> So.....is the world of DSP a world of scientific men, or is
> it a world of snake-oil charlatans and of religious loonies?
>
> Where do these Diracian impulses come from?

They're all the same, so they come right off dirac.
--
john

Why aren't you both guys dealing with something more earthbound?I never
bothered to deal with diracian functions.You should better persuade people
*not* to use airconditions but use fans instead, preferably with some
arguments (the chilling effect the fans work with) etc.

--
Tzortzakakis Dimitrios
major in electrical engineering, freelance electrician
FH von Iraklion-Kreta, freiberuflicher Elektriker
dimtzort AT otenet DOT gr
? <[email protected]> ?????? ??? ??????
news:[email protected] oups.com...
>
> Polymath wrote:
> > A good question indeed!
> >
> > And a question that the so-called experts and authors in this
> > field fail to answer time and time again.
> >
> > The Diracian, or Unit Impulse is a very good mathematical tool
> > to analyse the response of systems once a mathmatical model of
> > those systems had been produced.
> >
> > It is, however, a poor mathematical claim to make that such
> > impulses are found to be part of a system when neither the
> > area nor the magnitude of such impulses are found anywhere in such
> > systems.
> >
> > To those who ask, "Who cares? I get good results." , I suggest that
> > their approach is unscientific and compares to the religious
> > loonies who sacrifice goats and virgins to stop the Sun falling
> > out of the sky and justify the continuing practice by the Sun
> > remaining in the sky.
> >
> > So.....is the world of DSP a world of scientific men, or is
> > it a world of snake-oil charlatans and of religious loonies?
> >
> > Where do these Diracian impulses come from?
>
> They're all the same, so they come right off dirac.
> --
> john
>

Dimitrios Tzortzakakis wrote:
> Why aren't you both guys dealing with something more earthbound?I never
> bothered to deal with diracian functions.You should better persuade people
> *not* to use airconditions but use fans instead, preferably with some
> arguments (the chilling effect the fans work with) etc.

It's no good. The arguments just add heat.
--
john

Talking of snake-oil, where did the factor of "T" come from in your
opening lines?

Consider a 16-bit ADC capable of 100 M Samples per sec. In the first
instance we'll use it to sample a geophysical signal of bandwidth
limited
to 300 HZ and sample at 1 kHz, with suitable analogue instrumentation
to match the input signal to the full range of the ADC. If we now
keep the circuit the same, but now sample at 65.536 MHZ, your claimed
factor of "T" will result in the 16 bit range being compressed down
to one bit. We know this doesn't happen -there will be more samples,
but they'll still be of the same magnitude and 16-bit range.

Assuming that we were able to generate a Diracian, and then produce
a comb of them by delays and auperposition, there wouldn't be a factor
of "T" in such superposition, so where does yours come from?

>
> robert bristow-johnson wrote:
> >
> > x(t)*q(t) = T*SUM{x[k]*d(t-k*T)} .------.
> > x(t)--->(*)------------------------------------->| H(f) |---> x(t)
> > ^ '------'
> > |
> > | +inf
> > '------- q(t) = T * SUM{ d(t - k*T) }
> > k=-inf
> >
> >
> > where: d(t) = 'dirac' impulse function
> > and T = 1/Fs = sampling period
> > Fs = sampling frequency
> >
> >
> > +inf
> > q(t) = T * SUM{ d(t - k*T) } is the "sampling function", is periodic
> > k=-inf with period T, and can be expressed as a
> > Fourier series. It turns out that ALL of
> > the Fourier coefficients are equal to 1.