blog (fast cosine osc)

[xoxos]

http://xoxosvst.wordpress.com/

first entry, i revised my mass-spring algo yesterday, blogging is more effective than waiting on a pdf that my never be finished

as such, it's also the fastest sine (well, initialised as cosine) oscillator i'm aware of, one multiply, two adds, and includes amplification.

my thanks to you

[duncanparsons]

good stuff!

[asomers]

good stuff!

+1. Here's a lil matlab script to demonstrate. Two questions:

1. Know how to calculate decay time as a function of gain?
2. Is there a way to start it off from 0 (i.e. so that it's a sine rather than cos?

Code: Select all

clear all;

fs = 44100;
f = 440;
w = 2*pi*f/fs;
k = -w*w;
mv = 0;
mp = 1;
t = 1;
gain = .9995; % controls decay time.  1 = sine wave
output = zeros(1, fs * t);
for (i=1:fs * t - 1)
    mv = k*mp + gain*mv;
    mp = mp + mv;
    output(i) = mp * 10^(-3/20);
end
plot(output);
sound(output, fs);

[duncanparsons]

1. given that gain^samples=value,
it follows that:
log[base:gain](value) = samples, and
gain = samples(th) root of value.
So using samplerate & the ms you want to decay over, you get number of samples that takes, you can decide what value you think is quiet enough to decay to (without bit truncation you'll never hit zero).. -60dB is commonly used for reverb, but for this you'd probably want to go a bit quieter.
I'll let you do what you need to with all that

2. Starting from zero - allow the thing to run for 3/4 of a cycle & pick it up from there. The thing needs an impulse to get it going, and zero doesn't fill the bill too well, whereas 1 is a superb choice.. for convenience, you could pre-calculate the values of mv & mp for 3/4 of a cycle and start there, but bear in mind that you'd need different values for each frequency/samplerate combination that is used, but that could all be worked out & tabulated during initialisation.

I'm curious as to why the output in the Matlab script is converted to amp; to my mind, it's already in that scale not in dB as suggested - without trying it, I think this would generate a slightly lumpy transcendental, but I could be wrong (I just don't remember seeing it like that in these filters before)

[asomers]

'm curious as to why the output in the Matlab script is converted to amp; to my mind, it's already in that scale not in dB as suggested - without trying it, I think this would generate a slightly lumpy transcendental, but I could be wrong (I just don't remember seeing it like that in these filters before)

It's not converting, just scaling by -3dB

[duncanparsons]

oh, so it is! sometimes it would help me no end if I actually read what was in front of me properly! [tho' I'd probably premult that as a constant in code..]

'slightly lumpy transcendental' would be a great track name..

[asomers]

1. given that gain^samples=value,
it follows that:
log[base:gain](value) = samples, and
gain = samples(th) root of value.
So using samplerate & the ms you want to decay over, you get number of samples that takes, you can decide what value you think is quiet enough to decay to (without bit truncation you'll never hit zero).. -60dB is commonly used for reverb, but for this you'd probably want to go a bit quieter.
I'll let you do what you need to with all that

Well, I had tried the classic way calculating t60 assuming exponential decay (as you would do in a reverb or filter pole), as follows:

Code: Select all

T = 1/fs;
t60 = 1;
tau = t60/6.91;
gain = exp(-T/tau);

Is the mass-spring exponential decay? By the time we hit fs*t60, the peak value is at -33 db. The result is the same if we do it your way:

Code: Select all

t = 1;
samps = t * fs;
minLevel = -60;
gain = nthroot(10^(minLevel/20), samps);

2. Starting from zero - allow the thing to run for 3/4 of a cycle & pick it up from there. The thing needs an impulse to get it going, and zero doesn't fill the bill too well, whereas 1 is a superb choice.. for convenience, you could pre-calculate the values of mv & mp for 3/4 of a cycle and start there, but bear in mind that you'd need different values for each frequency/samplerate combination that is used, but that could all be worked out & tabulated during initialisation.

I was hoping there was some trick using an initial velocity and 0 start position. In other words the impulse is in the form of velocity rather than position.

[duncanparsons]

the second part of the assigment to mv makes it exponential, but with a variable DC provided by the first bit - that's what probably throws it off. As I recall most slew rate filtering of this kind use the time constant as a marker for when you get to ~37% of the target, which I've always found kinda irritating!
The variation to the formula so that it works fully could probably be worked out empirically, either that or search the mathworld forums..

the sine thing.. you possibly could do it like that, but the maths is beyond me for the little intake of tea I've had so far today It's certainly viable, in that all you have to do is prime mv & mp with the correct state values for that moment at 270degrees - an iterative method is by far the simplest way

[xoxos]

*edit*

thank you, had a moment to consider what was being said myself. eg. condition #1 is an expression that had escaped my attention, having never used logarithms before dsp. i'll have a go at these calibrations before making more things

[xoxos]

my uninformed observations on gain. i'll have to think more before these gel mathematically for me.

i'm using .01 as the base, which works well for the time for the peak to reach 1/10th (guessing this is known). it seems multiplying the # of samples seems accurate (this by holding a piece of paper up to the screen on my frequency analyser and setting 90 second decays ) by 1 + the 10s of decibels.. eg. -60dB is times 7.. -20dB is times 3.. and by 0.5 is the point where the peak decays to .316227

(of course, if it is exponential, then # of samples would be a viable target.. i'll have to figure out later why i need two)

knowing the way i do things, i wouldn't be surprised if i can devise a velocity based initialiser. arbitrary phase initialisation would be versatile

[xoxos]

a sine can be derived from a cursory observation of angular frequency and peak amplitude. thereafter, arbitrary phase is trivial. good questions there is a minor discrepancy in amplitude produced by non-90º angles.. 45º slightly over, 30º and 60º slightly under. that's new territory for me and i don't think i'm going to fix it today.

[mahaya]

Really interesting article, thanks for that!

I was hoping there was some trick using an initial velocity and 0 start position. In other words the impulse is in the form of velocity rather than position.

You'll get that with the following starting position/velocity:

Code: Select all

mp0 = 0; 
mv0 = sqrt(-k);

For an arbitrary angle use:

Code: Select all

mp0 = cos(angle);
pv0 = sin(angle)*sqrt(-k);

Chris

[aciddose]

this is just a state-variable filter, or a allpass with feedback. it isn't one multiplication, it's two, and two subtractions. you can obviously change the signs of the multiplications to change them to addition or subtraction as required to move parts of the equation around.

fundamentally however it's very, very well known. it's also two multiplies.

[asomers]

[quote="mahaya"]For an arbitrary angle use:

Code: Select all

mp0 = cos(angle);
pv0 = sin(angle)*sqrt(-k);

Thanks For the record, sqrt(-k) = w

[xoxos]

it's also two multiplies.

the counting was in reference to oscillation.

..any mavens interested in contributing a correction to the change in gain??

45° offsets slightly increase gain, and the intermittent multiples of 30° i tried slightly decreased gain. in both cases by less than 1/10th, so suitability for audio applications may not be strongly affected.