Underlying Maths
Underlying Maths
Note: This article needs a number of diagrams to finish it off. If after reading this you think you can create appropriate diagrams (.jpg or .gif files) then please send them to me.
This article follows from the introduction to understanding our DSP module at
http://humanise.org/demining/dsp
Gentle Maths Version:
Now lets try describing the maths behind the Phase Lock Loop and the Two Channel Lock-In Amplifier using nothing more than simple high school maths.
It has been stated that the heart of all Fourier analysis is, amazingly, all based on the following single high school trigonometry formula for the product of two sines:
sin(X1) * sin(X2) = ½cos(X1-X2) - ½cos(X1+X2)
The * is used here to indicate multiplication. In the maths below we leave out the * unless we want to emphasize the multiplication. For sine waveforms, the above X1 and X2 are based on time and frequency.
For simplicity, let us assume that Frequency is given to us in units of degrees per second where a complete cycle of the sine wave is 360 degrees. This is sometimes called Angular Velocity. We will use F to represent this Frequency. Let us assume that the sine table or sine lookup we use (sin for short) is based on degrees. Let us also assume that time represented by t is given in seconds. For these sine waves, the X1 and X2 above are as follows:
(F t) + θ
Where θ = Initial Phase Angle In Degrees
As F is an Angle per time (degrees per second), at each time point we multiply the Angle per time by the current time (in seconds) we get an angle (in degrees) which is what we need for looking up the the sin lookup.
Unfortunately, frequencies are not normally specified as degrees per second. Rather, they are specified in cycles per second. This creates a mishmash between the units used for angle with the the sin lookup and the units used for angle per time in the frequency specification. Let us represent these cycles per second frequency by f. Then
F = 360 f
With the sin lookup still being based on angles that have 360 units per waveform unit, X1 and X2 now needs to be as follows:
(360 f t) + θ
We have had to place in the factor of 360 to take account of the mishmash with 360 degrees in each cycle in the sin lookup but frequency being specified with whole cycles per second.
Note that with each of the above definitions of X1 and X2, once we have a fixed input frequency and phase, the only thing that is varying is time. For frequencies, the initial equation relates the Amplitude of these waves to time.
Now lets start using the initial equation with waveforms.
sin(X1) sin(X2) = ½cos(X1-X2) - ½cos(X1+X2)
When this equation is used on waveforms, it tells us what happens when we have one sin(X1) waveform input plus another another sin(X2) waveform input and we multiply them together as they are being input. We are just multiplying the same point in time on one waveform to the same point in time of the other waveform. The resultant waveform is ½cos(X1-X2) - ½cos(X1+X2)
Lets see what happens when X1 and X2 are different frequencies. Lets assume that these frequencies are f1 and f2 and that both of these are of an amplitude of 1. If we look at the first resultant term
½cos(X1-X2)
With waveform X1 being a frequency f1 and X2 being a frequency f2, this will be equal to
½cos((360 f1 t) - (360 f2 t))
This is simplified to
½cos(360 (f1-f2) t)
This is a waveform with a frequency that is the difference between the two input frequencies. As we know these frequencies are different (we made that assumption at the beginning) then this must produce a wave form that has a non zero frequency that is the difference between the two frequencies. As this is a standard sine waveform, it is symmetric about the X axis. This is, its an AC only waveform with an amplitude of one half (due to the ½ previous to the cos(X1-X2)).
Similarly with the ½cos(X1+X2) term, the result
½cos((360 (f1+f2) t) + (2 θ))
is a waveform that is the addition of the two frequencies. This produces a second waveform of half the amplitude of the input waveforms. In this case it has a frequency that is the sum of the two frequencies. As this is also a standard sine waveform, it is also symmetric about the X axis.
The addition of two pure AC complete waveforms always produces a resultant waveform that is pure AC only. The word 'complete' is here used to indicate that we are looking at complete cycles, not just part of a cycle of the waveform. As can be seen from the above, so long as the frequencies are different, there is no DC offset in the result.
Lets assume that X1 and X2 above are exactly equal. Start again with the initial equation:
sin(X1) sin(X2) = ½cos(X1-X2) - ½cos(X1+X2)
As cos(X1-X1) = cos(0) = 1
sin(X1) sin(X1) = ½ - ½cos(2 X1)
In this case, the same input sine waves multiplied by each other produce a DC offset (equal to half the amplitude of the input sine waves) plus a wave at double the frequency i.e. (2 * X1) and half the amplitude of the input sine waves.
One strange aspect of sine waves is that whether a wave is added or subtracted doesn't change its position with respect to the Y axis. That is, it doesn't give or change the wave's DC offset. As the wave is minus, it is equivalent to the wave being 180 degrees out of phase. That is, the following equation is equivalent to the last one:
sin(X1) sin(X1) = ½ + ½cos(2 X1 + 180 degrees)
Lets assume that X1 and X2 above are equal except that X2 is phase shifted forward by 90 degrees.
sin(X1) sin(X1 + 90) = ½cos(0 - 90) - ½cos(2 X1 + 90)
As cos(-90) = 0
sin(X1) sin(X1 + 90) = - ½cos(2 X1 + 90)
= ½cos(2 X1 - 90)
= ½sin(2 X1)
In the case of waveforms this means that when input sine waves of the same frequency but with a phase offset of 90 degrees are multiplied by each other, a wave at double the frequency i.e. (2 * X1) and half the amplitude of the input sine waves is produced. As the wave is minus, it is equivalent to the wave being 180 degrees out of phase but the extra 90 degrees in (2 X1 + 90 degrees) then shifts it back 90 degrees. Note that this wave is symmetrical about the X axis. There is no DC offset.
If we did the same as above but with X2 phase shifted back by 90 degrees (the same as phase shifted forward by 270 degrees) then we would get a similar result in that the resultant wave has no DC offset.
sin(X1) sin(X1 - 90) = - ½cos(2 X1 - 90)
= - ½sin(2 X1)
Lets assume that X1 and X2 above are equal except that X2 is phase shifted by 180 degrees. This is the equivalent to being inverted. Starting at our initial equation we see that:
sin(X1) * (-1 * sin(X1)) = -1 * (½cos(X1-X1) - ½cos(X1+X1))
= ½cos(X1+X1) - ½cos(X1-X1)
Again, as cos(X1-X1) = cos(0) = 1
= ½cos(2 X1) - ½
In the case of waveforms this means that when input sine waves that are equal except one is inverted from the other, are multiplied by each other, they produce a minus DC offset (equal to half the amplitude of the input sine waves) plus a wave at double the frequency i.e. (2 * X1) and half the amplitude of the input sine waves.
From these equations we can see that we only obtain a DC offset when we multiply two waveforms if they are exactly at the same frequency. Even then, we only obtain the DC offset if the waves are not at 90 or 270 degrees phase difference from each other. The DC offset reaches a maximum positive value when the waves are exactly in phase with each other and reach a maximum negative value when the waves are 180 degrees out of phase.
Also note that we have simplified the basic equation by leaving off the amplitudes. We have assumed that the amplitudes of the two input waves are both 1. Assume the amplitude of the sin(X1) input wave is A and the amplitude of the sin(X2) input wave is B then the original equation would be
A sin(X1) * B sin(X2)) = (A B) (½cos(X1-X2) - ½cos(X1+X2))
If the amplitudes of the input sin(X1) and sin(X2) waves is not 1 then all the previous formulas can be corrected by simply multiply the results by a factor of (A * B).
The previous equations should have explained the workings of the Phase Lock Loop. The following expands this to show the workings of the Two Channel Lock-In Amplifier.
Lets see what happens if we create an X2 waveform that is the same as an input X1 except that it is phase shifted an angle θ from X1.
sin(X1) sin(X1 + θ) = ½cos(X1-X1 - θ) - ½cos(2 X1 + θ)
= ½cos(- θ) - ½cos(2 X1 + θ)
The second term ½cos(2 X1 + θ) by itself is a pure sine wave that is symmetrical about the X axis. It has no DC component. If we are only interested in the resultant DC component then we can throw away this term. The first term (½cos(- θ)) gives us the value of the resultant waveform if we averaged it, filtering out the AC component. The DC component would be
½cos(- θ) = ½cos(θ)
½cos(θ) is a constant so we have a constant DC offset of ½cos(θ). If we perform this with an input wave amplitude of A and an internally created wave of amplitude B, and measure the resultant DC offset as Y1 then
Y1 = ½ A B cos(θ)
If we now create an X3 waveform that is the same as X2 except that it is phase shifted by (θ + 90) degrees from X1 and use that to multiply against another copy of the input X1 waveform.
sin(X1) sin(X1 + θ + 90) = ½cos(X1-X1 - (θ + 90)) - ½cos(2 X1 + θ + 90)
Again we are only interested in the DC component so we throw away the second term. The first term gives us
½cos( - (θ + 90)) = ½cos(θ + 90) = ½sin(θ)
If we perform this with an input wave amplitude of A and an internally created wave of amplitude of B, (assume that X3 is the same amplitude as X2) and measure the resultant DC offset as Y2 then
Y2 = ½ A B sin(θ).
If we now square Y1 and Y2, add these squares together and then take the square root
√ (Y12 + Y22)
= √ ((½ A B cos(θ))2 + (½ A B sin(θ))2)
= √ (¼ A A B B (cos(θ)2) + ¼ A A B B sin(θ)2)
= √ (¼ A A B B ((sin(θ)2) + (cos(θ)2)))
Note: sin(X)2 + cos(X)2 = 1
This is a trigonometrical identity that is easily derived from the main trigonometrical identity we have been using. Using this reduces the above to
= √ (¼ A A B B)
= ½ A B
Consequently, if we know B then the Amplitude A can be calculated from
A = 2 √ (Y12 + Y22) / B
Using the above, if we have an input waveform X1 that contained a known frequency but we didn't know the phase or amplitude of this frequency, we can find the amplitude by creating two reference waveforms, both of Amplitude B and frequency the same as X1. The second internally created reference waveform would be created with a 90 degree phase shift from the first reference waveform. We also create a second copy of the input waveform.
We now multiply the input waveform by the first reference waveform and average the result so that we measure the resultant Y1. We also multiply the copy of the input waveform by the second reference waveform and average the result so that we measure the resultant Y2. The amplitude A of the known frequency within the input waveform can now be calculated as
A = 2 √ (Y12 + Y22) / B
Lets denote by Y3 what would have been measured had the internally created reference waveform been exactly in phase with the input frequency. This gets rid of the ½ that's needed when comparing to the input waveform. The measured Y1 and Y2 are the sine and cosine components of the Y3 amplitude. They have the same relationship as the clock hand (of Y3 length) has to the X and Y axis or the same relationship as the hypotenuse (of Y3 length) of a right angle triangle.
Consequently, we can find the true amplitude or "magnitude" as the hypotenuse of this right triangle by taking the square root of the sum of the squares of the sine and cosine components. This is called the "vector sum" of the components. Similarly, we can find the phase as the angle whose tangent is equal to the sine component divided by the cosine component.
The sine and cosine components are often referred to as the "in-phase" and "quadrature" components. Some engineering and mathematicians prefer to use "real" and "imaginary" for the sine and cosine components. From this they create some very intimidating expressions using an imaginary number (√ -1) in a style of maths called "complex notation".
Also note that the first internally created reference frequency can be at any phase difference to the actual signal. Given any particular pure sinewave, without a DC component, at a particular frequency f1, it can be fully described by specifying its amplitude and phase. Alternatively, the sine wave can be described by specifying two amplitudes, these being the equivalent in-phase and quadrature components of that frequency. As we can change where (which phase angle) our first reference frequency is, there is effectively an infinite number of pairs of amplitudes that can be used as an equivalent to the phase plus amplitude pair. As long as we keep in mind where the the reference is in each case, all these pairs of amplitudes are equivalents.
With VLF metal detectors, where we are both transmitting the frequency on one coil and receiving the resultant frequency on another coil, the sine wave being output is normally used as our reference (sine component) frequency for the primary received frequency. Consequently, in this case it clearly makes sense to call the sine wave component the in-phase component. Changes in the received amplitude of this component are normally caused by resistive property effects in the soil and anomalies near the coil. Note that resistive is here used to indicated properties that change the received voltage in phase with the signal like resistance does. It doesn't necessarily mean resistance itself. Changes in the received amplitude of the cosine component are normally caused by reactive (such as capacitance and inductance) property effects of the soil and anomalies near the coil.
When programming using any of the above equations, you need to recognise the order in which the sine or cosine waves will be processed. The sine graphs with time increasing on the X axis are a representation of waves that are travelling from right to left. Generally, we expect things on paper to be represented from left to right so this works out to be opposite to what many people expect.
To confirm that it is this direction note what happens if we place the sine or cosine graph onto a long tape with time in seconds marked as normal along the X axis. We cover over the tape except that we leave a small view port that allows us to see the tape at the zero time mark (normally the Y axis). We now fire a starters gun to mark our zero time point. Following the gun firing we pull the tape so that the time marks of 1 second, 2 seconds etc are viewable at our view port at 1 second, 2 seconds etc after our zero time. To do this, we have to pull the tape to the left so that the graphed sine wave is moving right to left.
The conclusion here is that if you have a sine wave and a cosine wave being input to the processor, the cosine wave will reach the processor at ¼ wave cycle ahead of the sine wave or ¾ wave cycle after the sine wave. The equation cos(X) = sin(X + 90) indicates that cosine is before sine by 90 degrees.
We should also point out the usefulness of the earlier discussed fact that multiplying two frequencies produces a sine wave that is the difference of the two frequencies plus a sine wave that is the sum of the two frequencies.
As an example, Let us assume that we are interested in using an audio card input to look at an input that has a frequency of 300Khz. This is a higher frequency than current audio cards can handle correctly. To make it suitable we could create an oscillator that is for example a fixed 295khz. We then multiply this oscillator frequency to the 300Khz input. The result is a waveform at (300 + 295)Khz plus a wave form at (300 - 295)Khz. We can then use a low pass filter to filter out the 595Khz waveform leaving us with just a 5Khz waveform which is very easy for the audio card to handle. As the original 300Khz either changes in frequency or changes in amplitude, the 5Khz input to the audio card will similarly change. This principle is used in many applications such as the super-heterodyne receiver
It also explains a common problem in many wave systems called Intermodulation Distortion. This distortion occurs when two sine waves of frequencies X1 and X2 in the input result in the creation of other unwanted frequency components. These frequencies include (X1+X2), (X1-X2), (2X1-X2), (2X2-X1), and generally (mX1 ± nX2) for integer m and n. Generally, the size of the unwanted frequency components falls rapidly as m and n increase.
A similar explanation for some of the above is given at the following link:
http://www.daqarta.com/eex01.htm
A follow on article on understanding our DSP module is at
http://humanise.org/demining/averaging
This analyses a number of requirements associated with the averaging part of the DSP module.