VHDL sine wave oscillator
Published on 20180922 22:00
11 min read
In category
FPGA
In my last article a VHDL I2S transmitter was presented which allows one to playback arbitrary wave data. Eventually the goal is to develop an additive synthesis engine. In this article the goal is being able to generate a sine wave.
As a starting point I use this VHDL sine generator sub component. This generator is nicely parametrized. You can specify the number of bits in the output resolution as well as the phase input.
The additive synthesis engine which is to be developed needs to be able to control the frequency of the oscillator. I am not sure yet on the required frequency resolution. I will pick a reasonable number and keep the formulas general so that the value can be easily changed later. The challenge in the development of this oscillator is to implement the following features:
 Configurable design parameter 'frequency resolution' which allows specifying how precisely the frequency can be controlled.
 Real time configurable parameter 'frequency', value between 1 Hz and half the sample rate.
 No floating point
 No general divisions
 Where possible use bit shifts to do multiplication / division.
Phase step
Theory
The used sine wave generator sub component has two design parameters:
 Amplitude resolution: 24 bit
 Phase space size. This value has to be a power of two and determines the phase resolution as well as frequency resolution. Given a target frequency resolution the minimum required phase space bit width can be determined.
 Each sample a small step will be made to advance the phase of the sine. This phase step is a rational number and can be represented by a numerator and divisor. The algorithm will calculate the divisor and numerator.
Step three has to be performed whenever the frequency is changed while the design is 'running'. Step two can be performed designtime. The next steps calculate the static designtime parameters. In the formulas below I use '>' to show its value given a sample rate and target frequency resolution.
Design time constants
Frequency resolution:
target_frequency_resolution = 0.01 Hz
Sample Rate:
sample_rate = 48000 Hz
max_frequency = sample_rate / 2
nax_frequency > 24000 Hz
Phase space size:
phase_space_size = sample_rate / target_frequency_resolution
phase_space_size > 4800000
The sine generator sub component requires the phase space to be a power of two. As phase_space_size is the minimum required phase space size it is required to round up to the nearest power of two.
Power of 2 phase space size:
power2_phase_space_bits = ceiling(log(phase_space_size) / log(2))
power2_phase_space_bits > 23 bits
power2_phase_space_size > 2^23 > 8388608
This results in a slightly better frequency resolution:
Desired frequency resolution:
target_frequency_resolution = 0.01 Hz
Quantized frequency resolution:
quantized_frequency_resolution = sample_rate / power2_phase_space_size
quantized_frequency_resolution > 0.0057220459
The reciprocal of the frequency resolution is the phase step. At each played back sample at 48 kHz sample rate the phase is increased by the phase step.
Phase step:
phase_step = 1 / quantized_frequency_resolution
phase_step > 174.7626666667
phase_step_bits = log2(phase_step)
phase_step_bits = 7.449 bits
assert: sample_rate * phase_step == power2_phase_space_size
The phase step needs to be rounded to a power of two so that it can be easily used in integer division and multiplication. It is possible to round down if the resulting frequency resolution is still above the target. If not round up:
power2_phase_step_bits = floor(log(phase_step)/log(2))
power2_phase_step_bits > 7 bits
power2_phase_step > 128
frequency resolution would be: 1/128 > 0.0078125
power2_phase_step_bits = ceil(log(phase_step)/log(2))
power2_phase_step_bits > 8 bits
power2_phase_step > 256
frequency resolution would be: 1/256 > 0.00390625
With 7 bits the frequency resolution is still above the target. At runtime the frequency of the oscillator can be set as an integer value by multiplying the target frequency with the scaling factor. For example to get 440 Hz, the frequency parameter of the oscillator will be set to 440 * 128 = 56320.
The combinations of the sample_rate and power2_phase_space_size vales imply a noninteger phase_step value (174.76...). The phase_step value needs to be scaled up to keep accuracy when using it in division. To determine the number of bits in the scaling factor first the determine the maximum error introduced by quantisation.
scaled_phase_step = trunc(phase_step * phase_step_scaling_factor)
quantised_phase_step_error = phase_step  scaled_phase_step / phase_step_scaling_factor
maximum_phase_error = maximum_frequency * quantised_phase_step_error
This maximum_error must be below 1. So the quantised_phase_step_error must be lower than:
max_quantised_phase_step_error = 1 / maximum_frequency
max_quantised_phase_step_error > 0.00004166...
max_quantised_phase_step_error_bits > 14.5507... bits
This implies the required number of bits in the scaled phase step:
scaled_phase_step_bits = ceiling( phase_step_bits + max_quantised_phase_step_error_bits)
scaled_phase_step_bits > 22
For the phase step scaling factor there are two options. Choose the lowest number of bits that still results in an error lower than one.
phase_step_scaling_factor_bits = ceiling ( max_quantised_phase_step_error_bits )
phase_step_scaling_factor_bits > 15
phase_step_scaling_factor > 2^15 > 32768
scaled_phase_step = trunc(phase_step * phase_step_scaling_factor)
scaled_phase_step > 5726623
phase_step_error = phase_step  scaled_phase_step / phase_step_scaling_factor
phase_step_error > 0.000001871744786
max_error = maximum_frequency * phase_step_error
max_error > 0.04492187486 must be <1
Or
phase_step_scaling_factor_bits = floor ( max_quantised_phase_step_error_bits )
phase_step_scaling_factor_bits > 14
phase_step_scaling_factor > 2^14 > 16384
scaled_phase_step = trunc(phase_step * phase_step_scaling_factor)
scaled_phase_step > 2863311
phase_step_error = phase_step  scaled_phase_step / phase_step_scaling_factor
phase_step_error > 0.00003238932291
max_error = maximum_frequency * phase_step_error
max_error > 0.7773437499 must be <1
So 14 bits can be used for the phase step scaling factor.
Run time parameters
Let's say the frequency to generate is 440.0078125 Hz. Then the input scaled frequency would be:
frequency_scaled = frequency * power2_phase_step
frequency_scaled > 56321
If it would be possible to use floating point arithmetic then the phase step would be:
phase_step_fp = ( power2_phase_space_size / sample_rate ) * frequency
or rewritten:
phase_step_fp > ( power2_phase_space_size * frequency ) / sample_rate
phase_step_fp > 76896.93867
The integer version of phase_step_fp consists of two parts: phase_step_decimal and phase_step_numerator. Phase_step_decimal will give the decimal part (in the example 76895) while the fraction (0.57333...). The decimal part is calculated as follows:
scaled_phase = frequency_scaled * scaled_phase_step
scaled_phase > 161264538831
decimal_divider_bits = power2_phase_step_bits + phase_step_scaling_bits
decimal_divider_bits > 21
phase_step_decimal = shift_right ( scaled_phase, decimal_deivider_bits)
phase_step_decimal > 76896
Generally the phase_step_decimal will equal the decimal part of the phase_step_fp. In corner cases when the fractional part of phase_step_fp is near zero the decimal part will be off by one. For example with a scaled frequency of 440573 (3441.97 Hz) the actual phase step will be 601529.0027. However the quantised phase step will be floor(601528.8912) = 601528 which is off by one.
Now the fractional part is calculated as a rational value. The value consists of a numerator and divisor. Recall the calculation of phase_step_fp:
phase_step_fp > ( power2_phase_space_size * frequency ) / sample_rate
This value can be converted to a rational number:
phase_step_divisor = sample_rate
phase_step_divisor > 48000
using frequency_scaled (an integer value) instead of the floating point frequency:
phase_step_numerator_incl_decimal = ( power2_phase_space_size * frequency_scaled ) / power2_phase_step
phase_step_numerator_incl_decimal = shift_right ( power2_phase_space_size * frequency_scaled, power2_phase_step_bits)
phase_step_numerator_incl_decimal > 3691053056
The fact that phase_step_fp is larger than one implies that the numerator is larger than the divisor. To get the fractional part without the decimal part the decimal value is subtracted:
phase_step_numerator = phase_step_numerator_incl_decimal  phase_step_decimal * sample_rate
phase_step_numerator > 45056
Assert that the rational number is indeed the fractional part of phase_step_fp:
phase_step_fp > 76896.93867
phase_step_numerator / phase_step_divisor > 0.93867
In the earlier described corner cases it could be that the numerator is still larger than the divider after subtraction of the decimal part. It will however always be smaller than 2 * divider. So when the numerator is larger than the divider it is possible to simply subtract the divider once to get the numerator corresponding to the fractional part. In that case the phase_step_decimal value is increased by one.
The sine phase step calculations in this section can also be found on Google sheets.
VHDL Design
In this article I focus on generating a single sine wave. The following will be tested:
 Setting the sine generator to a certain frequency.
 Updating the frequency : This should not introduce noticeable audio 'glitches'.
The top level design simply connects an oscillator (the sine_wave component) to the I2S_sender. Two buttons enable toggling sine wave frequency. The sine_wave component takes a frequency as its input and presents the sine wave samples at its output.
Interface
The sine wave generator has the following interface:
entity sine_wave is
Port ( resetn : std_logic;
MCLK_in : in std_logic;  Master clock of the I2S Sender
freq_in_ce: in std_logic;  if '1' freq_in will be sampled on rising edge
freq_in: in frequency_t;  the scaled frequency of the sine
wave_out : out sample_t  the generated output samples
);
end sine_wave;
Architecture
The architecture consists of the sincos generator sub component and three processes:

Calculating the sample clock. Each tick the sine phase will need to be updated.

Accepting new frequency updates and update the phase step value. Calculating the phase step value is accomplished through a procedure
procedure Calculate_Phase_Step( constant frequency_scaled: in frequency_t; variable phase_step: out phase_step_t)

Generate sine samples by advancing the phase with the phase step. The phase value is presented to the sincos generator sub component. Updating the phase is again accomplished with a procedure. The implementation will use the decimal phase step and fractional phase step values. Each sample the phase will be increased with the decimal part: phase_step_decimal. Each sample the phase_step_numerator will be added to a counter. When the counter value is above the phase_step_divisor (the sample rate) value the phase will be advanced by one and the counter is decreased by the divisor.
type phase_step_t is record decimal : phase_step_decimal_t; fraction : phase_step_fraction_t;  fractional / sample_rate end record; type phase_state_t is record step : phase_step_t;  the decimal part, added each step current: phase_t;  the fractional part, if it overflows (above sample_rate)  then the it is reset and current is increased by one current_fraction: phase_fraction_t; end record; procedure Advance_Phase( variable phase: inout phase_state_t)
Testing
Two testbenches have been setup.
Testing phase step calculation.
The 'sine_sim' tests the calculation of the design time constants and calculating the phase step decimal and fractional values. This is fairly simple:
 calculate the real phase step using standard floating point arithmetic.
 calculate the same value using the phase step calculation algorithm.
 compare the values: if the difference is too large report an error.
The phase step will be added to the phase each time a sample is generated. A test is performed to check that the phase is updated correctly:
 Advance the phase a number of times.
 Update a version of the phase given the phase step as a floating point value.
 Update a version of the phase given the phase step as a decimal and fractional part.
 Each iteration compare if the difference between both both version of the phase is small.
Testing the sine wave generator
The simulation sets up a sine wave generator and configures its frequency to 440 Hz. The output samples are shown in Vivado's wave view. You can configure the way a signal is presented in the wave view. By showing the wave output as an analog signal with radix 'signed decimal' a nice sine wave is readily visible. Two markers were added to check the frequency of the sine wave. The check indeed shows that the sine wave has a frequency of 440 Hz.
Final thoughts

It would be interesting to check if the fixed point value types in VHDL 2008 could simplify things. Getting the fixed point package working in Vivado seemed a hassle.

Using an embedded bit scope (Xilinx ila) helped in finding a bug in the i2s_sender. It was simulating correctly, but synthesized with an bug. One of the steps setting up an ila is marking signals to be debugged with an attribute 'mark_debug'. I was able to enable and disable debugging using a generic variable:
attribute mark_debug of signal_name : signal is boolean'image(debug_generic_var)

In VHDL setting a constant to a value given a condition is not straight forward. You can use a function to return a value based on a condition and use that function to set the constant:
function sel(Cond: BOOLEAN; If_True, If_False: real) return real is begin if (Cond = TRUE) then return(If_True); else return(If_False); end if; end function sel; constant X : real := sel( TRUE, 1.0, 2.0);  will set X to 1.0
The VHDL code can be found in the i2s_sender repository.