harmonic_osc.h

Go to the documentation of this file.

/*
Copyright (c) 2020 Electrosmith, Corp, Emilie Gillet
 
Use of this source code is governed by an MIT-style
license that can be found in the LICENSE file or at
https://opensource.org/licenses/MIT.
*/
 
#pragma once
#ifndef DSY_HARMONIC_H
#define DSY_HARMONIC_H
 
#include <stdint.h>
#include "Utility/dsp.h"
#ifdef __cplusplus
 
 
namespace daisysp
{
template <int num_harmonics = 16>
class HarmonicOscillator
{
  public:
    HarmonicOscillator() {}
    ~HarmonicOscillator() {}
 
    void Init(float sample_rate)
    {
        sample_rate_ = sample_rate;
        phase_       = 0.0f;
 
        for(int i = 0; i < num_harmonics; ++i)
        {
            amplitude_[i]    = 0.0f;
            newamplitude_[i] = 0.f;
        }
        amplitude_[0]    = 1.f;
        newamplitude_[0] = 1.f;
 
        SetFirstHarmIdx(1);
        SetFreq(440.f);
 
        recalc_ = false;
    }
 
    float Process()
    {
        if(recalc_)
        {
            recalc_ = false;
            for(int i = 0; i < num_harmonics; ++i)
            {
                float f = frequency_
                          * static_cast<float>(first_harmonic_index_ + i);
                if(f >= 0.5f)
                {
                    f = 0.5f;
                }
                amplitude_[i] = newamplitude_[i] * (1.0f - f * 2.0f);
            }
        }
 
        phase_ += frequency_;
        if(phase_ >= 1.0f)
        {
            phase_ -= 1.0f;
        }
        const float two_x = 2.0f * sinf(phase_ * TWOPI_F);
        float       previous, current;
        if(first_harmonic_index_ == 1)
        {
            previous = 1.0f;
            current  = two_x * 0.5f;
        }
        else
        {
            const float k = first_harmonic_index_;
            previous      = sinf((phase_ * (k - 1.0f) + 0.25f) * TWOPI_F);
            current       = sinf((phase_ * k) * TWOPI_F);
        }
 
        float sum = 0.0f;
        for(int i = 0; i < num_harmonics; ++i)
        {
            sum += amplitude_[i] * current;
            float temp = current;
            current    = two_x * current - previous;
            previous   = temp;
        }
 
        return sum;
    }
 
    void SetFreq(float freq)
    {
        //convert from Hz to phase_inc / sample
        freq       = freq / sample_rate_;
        freq       = freq >= .5f ? .5f : freq;
        freq       = freq <= -.5f ? -.5f : freq;
        recalc_    = cmp(freq, frequency_) || recalc_;
        frequency_ = freq;
    }
 
    void SetFirstHarmIdx(int idx)
    {
        idx                   = idx < 1 ? 1 : idx;
        recalc_               = cmp(idx, first_harmonic_index_) || recalc_;
        first_harmonic_index_ = idx;
    }
 
    void SetAmplitudes(const float* amplitudes)
    {
        for(int i = 0; i < num_harmonics; i++)
        {
            recalc_          = cmp(newamplitude_[i], amplitudes[i]) || recalc_;
            newamplitude_[i] = amplitudes[i];
        }
    }
 
    void SetSingleAmp(const float amp, int idx)
    {
        if(idx < 0 || idx >= num_harmonics)
        {
            return;
        }
        recalc_            = cmp(amplitude_[idx], amp) || recalc_;
        newamplitude_[idx] = amp;
    }
 
 
  private:
    bool cmp(float a, float b) { return fabsf(a - b) > .000001f; }
 
    float sample_rate_;
    float phase_;
    float frequency_;
    float amplitude_[num_harmonics];
    float newamplitude_[num_harmonics];
    bool  recalc_;
 
    int first_harmonic_index_;
};
} // namespace daisysp
#endif
#endif