Complex numbers, phasors and the meaning of j
Every IQ sample, every FFT bin and every phase-based delay estimate is an arrow in a plane. This module rebuilds that arrow from nothing: j as a quarter turn, multiplication as rotation, Euler's formula three ways, the rotating phasor, why real signals carry negative frequencies, and what the "imaginary" part of a signal physically is.
You will be able to
- Do complex arithmetic by hand and draw each operation as an arrow move
- Prove \(e^{j\theta} = \cos\theta + j\sin\theta\) by series, by a differential equation and by a limit
- Explain negative frequencies and the quadrature (Q) component in physical terms
- Demodulate four samples to I/Q by hand, convert phase to delay, and unwrap phase safely
- Predict the effect of adding phasors, of beats and of differentiation (\(\times j\omega\))
01.1Numbers as arrows; j is a quarter turn
Start with the real line. A real number is an arrow from 0: its length is the magnitude and its direction (left or right) is the sign. Multiplying by a positive number stretches the arrow. Multiplying by a negative number also turns it around. That second effect, a rotation, is the whole seed of complex numbers.
Step 1 · Multiplying by −1 is a half turn
Multiplication by −1 keeps the length and reverses the direction, a rotation by 180°. Two half turns make a full turn, which is why \((-1)(-1) = +1\). Nothing about this is a rule to memorize; it is geometry.
(−2) × (−1) = +2 (another 180°: back home)
A sign is a direction, and multiplying by −1 is a rotation by 180°.
Step 2 · j is the operation that, done twice, is a half turn
Ask for an operation that, applied twice, gives −1. Geometrically: what done twice is a half turn? A quarter turn. Call it j. The equation \(j^2 = -1\) is the sentence "two quarter turns are a half turn." No arrow that stays on the real line can turn 90°, so the arrow must leave the line: the second dimension is forced on us. Engineers write j because i is current; mathematicians write i.
- \(j\)
- imaginary unit, \(j^2=-1\)
\(1/j = -j\), because \(j\cdot(-j) = -j^2 = +1\) (a quarter turn backwards)
j is not a quantity you imagine; it is the instruction "turn 90° counterclockwise."
Step 3 · The complex plane and the rectangular form
Once arrows can point anywhere in the plane, any arrow is "a units along 1, then b units along j": \(z = a + jb\). The horizontal coordinate a is the real part, the vertical b the imaginary part. The word "imaginary" goes back to Descartes (1637) and was meant dismissively; the b coordinate is exactly as real as a, it is just perpendicular to it.
\(\operatorname{Re} z = 3,\ \operatorname{Im} z = 4\) (Im z is the number 4, not 4j)
A complex number is a 2-D arrow written as (along 1) + j(along j).
Step 4 · Polar form: length and angle
The same arrow is described by its length \(r = \lvert z\rvert\) (Pythagoras) and its angle \(\theta = \angle z\) from the positive real axis. Going back uses the definitions of cos and sin. Section 01.3 shows that \(\cos\theta + j\sin\theta\) is \(e^{j\theta}\), which gives the compact polar form.
- \(a,\ b\)
- real and imaginary parts (rectangular coordinates)
- \(r = \lvert z\rvert\)
- length (magnitude, modulus)
- \(\theta = \angle z\)
- angle from the positive real axis, rad, in \((-\pi, \pi]\)
- \(z=a+jb\)
- complex number
- \(\operatorname{atan2}(b,a)\)
- four-quadrant arctangent, rad
- \(j\)
- imaginary unit, \(j^2=-1\)
back: 5 cos 53.13° = 3.000, 5 sin 53.13° = 4.000
−3 − 4j: atan(b/a) = atan(4/3) = 53.13° (wrong), atan2 = −126.87°
Rectangular form is for adding, polar form is for multiplying; use atan2, never atan(b/a).
Read "−" as "turn around" and "j" as "turn left a quarter." Then \((-1)(-1) = 1\) and \(j\cdot j = -1\) are the same kind of statement: turning twice. The algebra of complex numbers is the algebra of turning and stretching arrows in a plane.
01.2Arithmetic is geometry: add, multiply, mirror, divide
Each arithmetic operation is one arrow move. Addition slides arrows tip to tail. Multiplication rotates and stretches. Conjugation mirrors. Division undoes a rotation and a stretch. Work each case once by hand, then drag the arrows in Demo 01.1.
Step 5 · Addition: tip to tail
Add the real parts and add the imaginary parts. Geometrically, put the tail of w at the tip of z; the sum runs from the origin to the new tip. This is vector addition, nothing more.
- \(a,\ b,\ c,\ d\)
- real and imaginary parts of the two numbers
- \(j\)
- imaginary unit, \(j^2=-1\)
\(\lvert z + w\rvert = \sqrt{16 + 36} = 7.211\)
\(\le \lvert z\rvert + \lvert w\rvert = 5 + 2.236 = 7.236\) (triangle inequality)
Adding complex numbers is adding displacement arrows; echoes from several scatterers add this way (01.8).
Step 6 · Multiplication: rotate and scale
Expand with \(j^2 = -1\). Then look at what multiplying by \(w = c + jd\) does to z: it makes c copies of z plus d copies of z turned 90° (that is what jz is). A stretched copy plus a stretched perpendicular copy is exactly a rotation by \(\angle w\) combined with a scaling by \(\lvert w\rvert\). The 2×2 matrix shows it.
- \(a,\ b,\ c,\ d\)
- real and imaginary parts of the two factors, as in (01.3)
- \(j\)
- imaginary unit, \(j^2=-1\)
- \(w = c + jd\)
- the multiplier; \(\lvert w\rvert = \sqrt{c^2 + d^2}\)
- \(\theta_w = \angle w\)
- its angle, rad
- \(r_1, r_2\)
- lengths of the two factors
- \(\theta_1, \theta_2\)
- angles of the two factors, rad
- \(z\)
- complex number (arrow) being multiplied by \(w\)
- \(c,\ d\)
- real and imaginary parts of \(w\)
- \(j\)
- imaginary unit, \(j^2=-1\)
- \(e^{j\theta}\)
- unit arrow at angle \(\theta\) (01.8)
- \(\Leftrightarrow\)
- is equivalent to
polar check: \(\lvert zw\rvert = 5 \times 2.236 = 11.18 = \sqrt{25 + 100}\) ✓
∠zw = 53.13° + 63.43° = 116.57° = atan2(10, −5) ✓
Multiply the lengths, add the angles. Multiplying by \(e^{j\alpha}\) (length 1) is a pure rotation: a phase shift.
Step 7 · Conjugate: the mirror, and |z|² = z z*
The conjugate \(z^* = a - jb\) flips the sign of the imaginary part: the same arrow mirrored in the real axis, same length, negated angle. Multiplying z by z* rotates by θ and then by −θ, so the result lands on the positive real axis with length \(r\cdot r\).
- \(z^*\)
- complex conjugate, \(a - jb = r\,e^{-j\theta}\)
- \(z=a+jb\)
- complex number; \(a,\ b\) real and imaginary parts
- \(\lvert z\rvert^2\)
- squared magnitude (power)
- \(\operatorname{Re},\ \operatorname{Im}\)
- real part, imaginary part
- \(j\)
- imaginary unit, \(j^2=-1\)
but \(z^2 = (3 + 4j)^2 = 9 + 24j - 16 = -7 + 24j\) (not the power)
Power is the arrow times its mirror; the conjugate undoes a rotation.
Step 8 · Division: unrotate and unscale
To divide by w, multiply top and bottom by w*: the denominator becomes the real number \(\lvert w\rvert^2\), and the top is a single multiplication. Geometrically z/w has length \(\lvert z\rvert/\lvert w\rvert\) and angle \(\angle z - \angle w\).
- \(z,\ w\)
- complex numbers, \(w\neq 0\)
- \(w^*\)
- complex conjugate of \(w\)
- \(\lvert w\rvert\)
- magnitude of \(w\)
\(= \dfrac{3 - 6j + 4j - 8j^2}{5} = \dfrac{11 - 2j}{5} = 2.2 - 0.4j\)
polar check: 5/2.236 = 2.236, 53.13° − 63.43° = −10.30° ✓
Division is multiplication by the mirror, divided by a length squared.
The four moves on one arrow, as a reference drawing:
Drag the arrows: add, multiply, mirror, divide
Exact formulaEvery IQ sample is one of these arrows. A phase rotation in an IQ beamformer is a multiplication by \(e^{j\alpha}\) (Step 6 with \(\lvert w\rvert = 1\)). Power and B-mode intensity are \(z z^*\) (Step 7). The Kasai lag-one autocorrelation \(R(1) = \sum_n \mathrm{IQ}_{n+1}\,\mathrm{IQ}_n^*\) uses the conjugate to subtract angles: \(\angle(z_2 z_1^*) = \theta_2 - \theta_1\) is the phase advanced between two pulses (Module 11).
01.3Euler's formula, three ways
The claim is (01.8): raising e to an imaginary power gives a unit-length arrow at angle θ radians. Three independent arguments lead to the same picture. Each one is useful later: the series explains numerical evaluation, the differential equation explains why phasors are the natural signals of linear systems (Module 03), and the limit explains why naive time-stepping of oscillators spirals outward.
- \(\theta\)
- angle in radians (not degrees)
- \(e^{j\theta}\)
- unit-length arrow at angle θ: real part \(\cos\theta\), imaginary part \(\sin\theta\)
- \(r\)
- length of a general arrow z
- \(j\)
- imaginary unit, \(j^2=-1\)
- \(z\)
- complex number, \(z=r\,e^{j\theta}\)
- \(\pi\)
- 3.14159, rad
Step 9 · The Taylor series, side by side
Write the three known series and substitute \(x = j\theta\) into the exponential. Because the powers of j cycle 1, j, −1, −j, the even terms are real with alternating signs (the cosine series) and the odd terms are imaginary with alternating signs (the sine series). Each term is the previous one turned 90° and multiplied by θ/k, so the partial sums spiral in on the answer.
- \(x\)
- any real (or complex) number
- \(\theta\)
- angle, rad
- \(j\)
- imaginary unit, \(j^2=-1\)
- \(n!\)
- factorial, \(n!=1\cdot2\cdots n\)
- \(\cdots\)
- continuing series
partial sums: 1 → 1 + j → 0.5 + j → 0.5 + 0.8333j → 0.5417 + 0.8333j → 0.5417 + 0.8417j → 0.5403 + 0.8417j
cos 1 = 0.5403, sin 1 = 0.8415; 8 terms leave an error of \(2.5\times10^{-5}\)
Plug jθ into the exponential series and the real and imaginary parts sort themselves into cos and sin.
Step 10 · The differential equation: velocity is perpendicular to position
The exponential is the function whose rate of change is proportional to itself: \(\frac{d}{d\theta} e^{a\theta} = a\,e^{a\theta}\). With a = j, the velocity \(f' = j f\) is the position turned 90°. A point whose velocity is always perpendicular to its position cannot get closer to or farther from the origin, so it moves on a circle. Its speed \(\lvert f'\rvert = \lvert f\rvert = 1\), so after "time" θ it has travelled an arc of length θ, which on a unit circle is an angle of θ radians.
- \(f(\theta)\)
- test function \(e^{j\theta}\), complex
- \(f'\)
- derivative \(df/d\theta\)
- \(f^*\)
- complex conjugate of \(f\)
- \(\theta\)
- angle, rad
- \(j\)
- imaginary unit, \(j^2=-1\)
- \(\Rightarrow\)
- implies
θ = π: f = −1 after an arc of length π, half the circumference 2π
\(e^{j\theta}\) is the point that moves around the unit circle at unit speed; this is why θ must be in radians.
Step 11 · The limit: n small rotations
Compound interest at an imaginary rate. Each factor \(1 + j\theta/n\) is a tiny rotation by \(\operatorname{atan}(\theta/n) \approx \theta/n\) with a tiny stretch \(\sqrt{1 + \theta^2/n^2}\). Multiply n of them: the angles add to \(n\operatorname{atan}(\theta/n) \to \theta\) and the length \((1 + \theta^2/n^2)^{n/2} \approx e^{\theta^2/(2n)} \to 1\).
- \(n\)
- number of small rotation steps, integer
- \(\theta\)
- total angle, rad
- \(j\)
- imaginary unit, \(j^2=-1\)
- \(\lim\)
- limit as \(n\to\infty\)
- \(\lvert\cdot\rvert\)
- magnitude
n = 4: 0.6289 + 0.9375j; exact 0.5403 + 0.8415j
θ = π, n = 4: length 2.614 at 152.6°; n = 100: length 1.0506 (\(e^{\pi^2/200} = 1.0506\)) at 179.94°
Many tiny rotate-and-stretch steps converge to a pure rotation; the length error shrinks like θ²/(2n).
All three routes give \(e^{j\pi} = -1\): a half turn. And they explain the polar form of Step 4: \(z = r(\cos\theta + j\sin\theta) = r\,e^{j\theta}\), so multiplication adds exponents and therefore angles.
Euler three ways: series, velocity, limit
Exact formulaEuler's formula is in radians. \(e^{j\cdot 90}\) is a rotation by 90 rad (≈ 5157°, i.e. 116.6° after removing whole turns), not by 90°. In code, every phase term is \(2\pi f t\) or \(2\pi f\tau\) with f in Hz and t in s; a missing 2π is the most common phase bug, and a degree/radian mix-up is the second.
01.4The rotating phasor \(e^{j\omega t}\)
Now let the angle grow with time, \(\theta = \omega t + \varphi\). The arrow spins. This single object, a phasor, is the building block of every signal in the rest of the course.
Step 12 · Frequency is turns per second
\(e^{j(\omega t + \varphi)}\) rotates at ω radians per second, i.e. \(f = \omega/(2\pi)\) full turns per second. One turn takes the period \(T = 1/f\). The starting angle φ is the phase.
- \(A\)
- amplitude (arrow length)
- \(f\)
- frequency, turns per second (Hz); the sign gives the direction
- \(\omega\)
- angular frequency, rad/s
- \(\varphi\)
- phase (angle at t = 0), rad
- \(T\)
- period, s
- \(z(t)\)
- complex phasor, rotating arrow
- \(t\)
- time, s
- \(j\)
- imaginary unit, \(j^2=-1\)
t = 50 ns: angle = 2π × 5 × 10⁶ × 50 × 10⁻⁹ = π/2 → z = j
t = 30 ns: angle = 54° → z = 0.588 + 0.809j
A phasor is a clock hand: its length is the amplitude, its speed is the frequency, its starting angle is the phase.
Step 13 · The sign of ω is the direction
\(e^{-j\omega t}\) turns clockwise. Its real part is identical to that of \(e^{+j\omega t}\) because cos is even, \(\cos(-x) = \cos x\); only the imaginary part changes sign, because sin is odd. So a negative frequency is not time running backwards; it is clockwise rotation.
both have real part 0; t = 30 ns: 0.588 ± 0.809j
The real part alone cannot tell which way the arrow turns; the imaginary part can.
Step 14 · cos and sin are shadows
Shine a light from above: the shadow of the arrow tip on the real axis moves as \(A\cos(\omega t + \varphi)\). Shine it from the side: the shadow on the imaginary axis moves as \(A\sin(\omega t + \varphi)\). The sine is the cosine a quarter period later: \(\sin\omega t = \cos(\omega t - \pi/2)\).
- \(A,\ \omega,\ \varphi\)
- amplitude, angular frequency (rad/s), phase (rad), as in (01.12)
- \(t\)
- time, s
- \(\operatorname{Re},\ \operatorname{Im}\)
- real part, imaginary part
- \(j\)
- imaginary unit, \(j^2=-1\)
real shadow cos 45° = 0.707, imaginary shadow sin 45° = 0.707
A real sinusoid is the shadow of a rotating arrow; the shadow loses the direction of rotation.
A rotating phasor and its two shadows
Exact formula01.5Why real signals need negative frequencies, and what "imaginary" means physically
One rotating arrow never stays on the real axis. A measured RF voltage is real. So a real sinusoid must be built from something whose imaginary part is zero at every instant. The answer is two arrows turning in opposite directions.
Step 15 · cos is two counter-rotating halves
Add Euler's formula for +θ and for −θ: the sines cancel and the cosines double. Subtract and the cosines cancel.
- \(e^{j\theta}\)
- counterclockwise unit arrow (positive frequency when \(\theta = \omega t\))
- \(e^{-j\theta}\)
- its clockwise mirror (negative frequency)
- \(\theta\)
- angle, rad
- \(j\)
- imaginary unit, \(j^2=-1\)
half-sum = (1.0 + 0j)/2 = 0.5 = cos 60°: the ±0.433j parts cancel
A real cosine is a pair of half-length arrows, one at +f and one at −f, whose vertical parts cancel.
Step 16 · So the spectrum of every real signal has two mirror halves
For \(A\cos(2\pi f_0 t + \varphi)\) the two arrows are \((A/2)e^{j\varphi}\) at \(+f_0\) and \((A/2)e^{-j\varphi}\) at \(-f_0\): a conjugate pair. Any real signal is a sum of such pairs, so its spectrum obeys \(X(-f) = X^*(f)\). The negative-frequency half carries no new information; it is the mirror image needed to cancel the imaginary parts. If the two halves are unequal, the sum leaves the real axis and the signal is complex.
- \(f_0\)
- frequency of the sinusoid, Hz
- \(X(f)\)
- spectrum (Fourier transform) of x(t)
- \(A,\ \varphi\)
- amplitude and phase of the cosine, as in (01.12)
- \(t\)
- time, s
- \(x(t)\)
- real signal
- \(X(-f),\ X^*(f)\)
- spectrum at \(-f\); complex conjugate of the spectrum at \(f\)
- \(j\)
- imaginary unit, \(j^2=-1\)
- \(\Leftrightarrow\)
- if and only if
halves 0.5 and 0.25: imaginary part (0.5 − 0.25) sin θ, at most 0.25 (an ellipse)
halves 0.5 and 0: a circle, the complex phasor \(0.5\,e^{j\theta}\)
Negative frequencies exist because real signals need a clockwise partner for every counterclockwise arrow.
Step 17 · I and Q: the imaginary part is the quadrature component
Expand a sinusoid of arbitrary phase with the angle-sum identity. It is always a mix of two fixed references: the cosine itself and the cosine shifted 90° ahead. The phasor \(A e^{j\varphi} = I + jQ\) stores the two mixing amounts. The real part says how much in-phase reference is present; the imaginary part says how much quadrature (90°-shifted) reference is present. Hardware measures Q with a second mixer whose local oscillator is delayed by a quarter period. Nothing about it is imaginary.
- \(A,\ \varphi\)
- amplitude and phase of the sinusoid
- \(\omega\)
- angular frequency of the reference, rad/s
- \(I\)
- in-phase amount: weight of \(\cos\omega t\)
- \(Q\)
- quadrature amount: weight of \(\cos(\omega t + 90^\circ)\)
- \(t\)
- time, s
- \(90^\circ\)
- quarter-cycle shift, \(\pi/2\) rad
- \(\operatorname{Re},\ \operatorname{Im}\)
- real part, imaginary part
- \(j\)
- imaginary unit, \(j^2=-1\)
check at ωt = 0: 0.866·1 + 0.5·cos 90° = 0.866 = cos 30° ✓
The imaginary part of a phasor is the amount of the 90°-shifted copy in the signal.
Step 18 · Q is what the signal was doing a quarter period earlier
Go back a quarter period from the reference instant \(t_0\) (where the reference cosine peaks): the phasor turns back by 90°, which is multiplication by −j, and the real part of \(-j(I + jQ)\) is Q. So for a narrowband signal the quadrature value is simply the real signal sampled T/4 earlier. This is the basis of quadrature sampling at \(f_s = 4f_0\).
- \(t_0\)
- time at which \(x\) is read, s
- \(T\)
- carrier period, \(1/f\), s
- \(x(t)\)
- real signal \(I\cos\omega t+Q\cos(\omega t+90^\circ)\) of (01.16)
- \(I,\ Q\)
- in-phase and quadrature amounts, as in (01.16)
- \(\operatorname{Re}\)
- real part
- \(j\)
- imaginary unit, \(j^2=-1\)
x(−T/4) = cos(−90° + 30°) = cos(−60°) = 0.500 = Q ✓
For a narrowband echo, I is "now" and Q is "a quarter period ago"; both are ordinary real measurements.
The same reading applies to a spectrum. \(X(f) = \int x(t)\cos(2\pi f t)\,dt - j\int x(t)\sin(2\pi f t)\,dt\), so \(\operatorname{Re} X(f)\) is the correlation of the signal with a cosine at f and \(\operatorname{Im} X(f)\) is minus its correlation with a sine at f: the quadrature content at that frequency. For a real signal, the even part of x(t) produces Re X and the odd part produces Im X. Module 02 builds the DFT on exactly this.
Two counter-rotating phasors make a real cosine
Exact formulaDoppler direction needs Q. A scatterer moving toward the probe shifts the echo phasor to a positive Doppler frequency, away from the probe to a negative one. A single real channel sees \(\cos(\pm 2\pi f_D t)\), identical for both signs (Step 13), so the spectral Doppler display cannot tell flow toward from flow away without the quadrature channel. A gain or phase mismatch between the I and Q channels makes the \(A_+ \ne A_-\) situation in reverse: some of the \(+f_D\) arrow leaks into \(-f_D\) (a 1 % gain and 1° phase error leave a mirror image at −40 dB).
01.6I/Q: two projections of one arrow, and demodulation as stopping the rotation
An RF echo is a phasor spinning at the carrier \(f_0\) whose length and angle change slowly. I/Q demodulation removes the fast spin so that only the slow arrow is left. Module 08 covers the filters, decimation and beamforming on IQ; here the goal is to see the mechanism with four numbers.
Step 19 · Multiply by \(e^{-j\omega_d t}\): spin backwards at the carrier rate
Write the echo as its two counter-rotating halves and multiply by a clockwise phasor at the demodulation frequency \(f_d = f_0\). The counterclockwise half stops dead; the clockwise half now spins twice as fast. A low-pass filter removes the fast one, and a factor 2 restores the amplitude. In hardware this is two mixers: \(I = 2\cdot\operatorname{LPF}\{x\cos\omega_d t\}\), \(Q = -2\cdot\operatorname{LPF}\{x\sin\omega_d t\}\).
- \(\omega_d = 2\pi f_d\)
- demodulation (mixer) frequency; here \(f_d = f_0\)
- \(\operatorname{LPF}\)
- low-pass filter that keeps the band near 0 Hz and rejects \(-2f_0\)
- \(x(t)\)
- real RF signal, \(A\cos(\omega_0 t+\varphi)\)
- \(t\)
- time, s
- \(\omega_0\)
- carrier angular frequency, rad/s
- \(A,\ \varphi\)
- amplitude and phase
- \(I,\ Q\)
- in-phase and quadrature parts of the IQ sample, as in (01.16)
- \(\mathrm{IQ}\)
- complex baseband sample, \(I+jQ\)
- \(j\)
- imaginary unit, \(j^2=-1\)
- \(\times 2\)
- gain restoring the factor \(A/2\) to \(A\)
at f₀ = 5 MHz the unwanted line sits at −10 MHz, far outside a ±5 MHz low-pass
Demodulation is a frequency shift by −f₀ followed by throwing away the shifted mirror.
Step 20 · I/Q by hand from four samples (\(f_s = 4f_0\))
Sample the echo once per quarter period. The clockwise phasor \(e^{-j\omega_0 t}\) at those instants is 1, −j, −1, j: no multiplications needed. Averaging four products is a crude low-pass that removes the \(2f_0\) term exactly, because that term alternates +, −, +, − at quarter-period sampling and sums to zero.
× (1, −j, −1, j): 0.866, +0.5j, 0.866, +0.5j
mean = 0.433 + 0.25j; × 2 = 0.866 + 0.500j = \(1\cdot e^{j30^\circ}\) ✓
Four quarter-period samples and sign flips give the phasor: I and Q are averages, nothing more mysterious.
Step 21 · A delay is a rotation of the IQ arrow
Delay the echo by τ and demodulate at \(f_0\): the slow arrow turns by \(-2\pi f_0\tau\) while its length (the envelope) is unchanged. One full turn is one carrier period, so the phase measures delay only modulo T.
- \(\tau\)
- extra round-trip delay of the echo, s
- \(A(t)\)
- envelope of the echo
- \(T_e\)
- reference echo time; choosing \(f_0 T_e\) an integer makes the reference phase zero
- \(\Delta\angle\)
- change of the IQ angle, rad (wrapped into \((-\pi, \pi]\))
- \(\Delta z\)
- depth shift for round-trip delay τ; c = 1540 m/s
- \(x(t)\)
- real RF echo
- \(t\)
- time, s
- \(f_0\)
- demodulation frequency, Hz
- \(\mathrm{IQ}(t)\)
- complex baseband echo
- \(j\)
- imaginary unit, \(j^2=-1\)
- \(c\)
- sound speed, 1540 m/s
- \(\to\)
- after demodulation gives
τ = 50 ns → −90°; τ = 200 ns → −360° ≡ 0° (the phase cannot tell)
depth: Δz = cτ/2 = 1540 × 13 ns / 2 = 10 µm for −23.4°
Phase is delay in disguise: −360° per carrier period, i.e. per λ/2 of round-trip depth.
Step 22 · Demodulating at the wrong frequency leaves a slow spin
If \(f_d \ne f_0\) the arrow is not quite stopped: it keeps turning at \(f_0 - f_d\). The envelope \(\lvert\mathrm{IQ}\rvert\) is untouched; the phase gains a ramp. That ramp is real information (a frequency offset), or a bias if it is mistaken for motion.
- \(f_0\)
- carrier frequency of the echo, Hz
- \(f_d\)
- demodulation frequency, Hz
- \(t\)
- time, s
- \(A(t-T_e)\)
- echo envelope, delayed by \(T_e\)
- \(T_e\)
- echo arrival time, s
- \(\mathrm{IQ}(t)\)
- complex baseband signal
- \(j\)
- imaginary unit, \(j^2=-1\)
m01_phasors_and_iq.m measures 0.501 MHz inside the echo
The envelope ignores the demodulation frequency; the phase does not.
I/Q of an RF echo: delay, demodulation frequency, bandwidth
SimulationJS: Gaussian-envelope 5 MHz echo sampled at 160 MHz, mixed with \(e^{-j2\pi f_d t}\), ideal FFT low-pass at 5 MHz, ×2. Dots show every 4th sample (40 MHz).01.7Phase: wrapping and unwrapping
An arrow at angle θ and one at θ + 2π are the same arrow. So any measured angle is known only modulo 2π, and turning measured angles back into a continuous phase history (unwrapping) needs an assumption.
Step 23 · Measured angles live in (−π, π]
angle() and atan2 return the principal value. A phase that grows steadily therefore comes back as a sawtooth: every time it passes +π it reappears near −π.
wrapped: 0, 2, 4 − 2π = −2.283, 6 − 2π = −0.283
A single phase value cannot tell how many whole turns have happened.
Step 24 · Unwrap by assuming small steps
Walk along the sequence. Whenever two neighbours differ by more than π, assume a wrap happened and add the multiple of 2π that brings the difference back into (−π, π]. MATLAB unwrap does exactly this.
−4.283 < −π → add 2π → 2.000
cumulative: 0, 2, 4, 6 ✓
Unwrapping is a guess that the true phase never jumps by more than half a turn between samples.
Step 25 · The guess fails when the true step exceeds π
A step of more than half a turn looks exactly like a step the other way, so unwrap picks the wrong one and every later sample inherits the error. This is the sampling theorem applied to phase: at least two samples per turn. In echo tracking a round-trip displacement Δz rotates the phase by \(-4\pi f_0\Delta z/c\), so the limit is λ/4 per sample; in Kasai Doppler it is the aliasing velocity.
- \(\Delta\varphi\)
- true phase step between consecutive samples (frames or pulses), rad
- \(\lambda = c/f_0\)
- wavelength, 0.308 mm at 5 MHz
- \(\mathrm{PRF}\)
- pulse repetition frequency, Hz
- \(v_{\text{alias}}\)
- largest unambiguous axial velocity, m/s
- \(\Delta z\)
- axial displacement between samples, m
- \(c\)
- sound speed, 1540 m/s
- \(f_0\)
- center frequency, Hz
- \(\lvert\cdot\rvert\)
- magnitude
PRF = 5 kHz: \(v_{\text{alias}}\) = 1540 × 5000 / (4 × 5 × 10⁶) = 0.385 m/s
step growing 0.1 rad per sample: first failure at sample 32 (true step 3.15 rad)
Phase is unambiguous only below half a turn per sample; above that it aliases like any undersampled signal.
Wrapped, unwrapped and true phase of an accelerating target
SimulationJS: phase ½·a·n² plus Gaussian phase noise (fixed seed), wrapped with atan2, then unwrapped with the same rule as MATLAB unwrap.01.8Adding sinusoids: phasor sums and beats
Real sinusoids of the same frequency are awkward to add with trigonometric identities and trivial to add as arrows. Sinusoids of different frequencies add as arrows that turn relative to each other, which produces beats.
Step 26 · Same frequency: add the phasors
Both signals share the factor \(e^{j\omega t}\), so it factors out. Add the stationary arrows, then let the result spin.
- \(A_1,\ A_2\)
- amplitudes of the two echoes
- \(\varphi_1,\ \varphi_2\)
- their phases, rad
- \(\omega\)
- common angular frequency, rad/s
- \(t\)
- time, s
- \(\operatorname{Re}\)
- real part
- \(j\)
- imaginary unit, \(j^2=-1\)
1∠0° + 0.5∠120° = 1 + (−0.25 + 0.433j) = 0.75 + 0.433j = 0.866∠30°
1∠0° + 1∠120° + 1∠240° = 0 (three-phase balance)
Amplitudes of equal-frequency sinusoids do not add; their arrows do.
Step 27 · Two scatterers in one resolution cell
Two equal echoes whose depths differ by Δz arrive with a round-trip phase difference \(4\pi\Delta z/\lambda\). Their sum has length \(2\lvert\cos(2\pi\Delta z/\lambda)\rvert\): full addition when the separation is a multiple of λ/2, complete cancellation at λ/4. Many random arrows of this kind are speckle (Module 07).
- \(\Delta z\)
- axial separation of the two scatterers, m
- \(\lambda\)
- wavelength at the centre frequency, m
- \(j\)
- imaginary unit
- \(\lvert\cdot\rvert\)
- magnitude of the summed phasor; the second echo has unit relative amplitude
A quarter wavelength in depth is half a turn of phase: two equal scatterers can erase each other.
Step 28 · Different frequencies: beats
In a frame that turns with the first arrow, the second arrow turns at the difference frequency Δf. The resultant length, which is the envelope, swings between \(A_1 + A_2\) and \(\lvert A_1 - A_2\rvert\) once per 1/Δf. The product identity says the same in real form.
- \(\Delta\omega = 2\pi\Delta f\)
- difference angular frequency; the envelope repeats every \(1/\lvert\Delta f\rvert\)
- \(\bar\omega\)
- mean (carrier) angular frequency
- \(\omega_1,\ \omega_2\)
- the two angular frequencies, rad/s
- \(t\)
- time, s
envelope period 1/Δf = 2 µs (the cosine factor has period 4 µs; its magnitude 2 µs)
nulls at 1, 3, 5 µs; m01_phase_and_derivative.m measures 2.000 µs
Two frequencies beat at their difference because one arrow laps the other Δf times per second.
Phasor addition and beats
Exact formula01.9Derivative = multiply by jω; preview of the analytic signal
The phasor makes calculus into arithmetic. Differentiating a rotating arrow gives its velocity, which (Step 10) is the arrow turned 90° and scaled by its angular speed.
Step 29 · d/dt of a phasor is jω times the phasor
The chain rule on \(e^{j\omega t}\) gives \(j\omega\,e^{j\omega t}\): multiply by ω (longer for faster arrows) and by j (90° ahead). Taking the real part gives the familiar result for cosines.
- \(\omega = 2\pi f\)
- angular frequency of the phasor, rad/s
- \(j\omega\)
- the derivative as a multiplier: ×ω in length, +90° in angle
- \(t\)
- time, s
- \(j\)
- imaginary unit, \(j^2=-1\)
it leads by 90°, i.e. by T/4 = 50 ns
Differentiation rotates every phasor by +90° and stretches it by its own angular frequency.
Step 30 · So a differentiator is a filter with H(f) = j2πf
Every frequency component is scaled by \(2\pi\lvert f\rvert\) and rotated by +90° (−90° for negative f, since j times a negative number points down). Doubling f doubles the gain: +6.02 dB per octave, +20 dB per decade. An integrator, \(1/(j2\pi f)\), does the opposite: −6 dB per octave, −90°. Differentiation therefore shifts a pulse's spectrum upward and amplifies high-frequency noise.
- \(H(f)\)
- frequency response of the derivative operator
- \(f\)
- frequency, Hz
- \(\angle H\)
- phase of \(H\)
- \(j\)
- imaginary unit, \(j^2=-1\)
- \(\lvert H\rvert\)
- magnitude of \(H\)
Gaussian pulse, σ = 100 ns, f₀ = 5 MHz: spectral peak 5.00 → 5.46 MHz after d/dt
(peak of \(f\,e^{-(f - f_0)^2/(2\sigma_f^2)}\) at \(\tfrac12\big(f_0 + \sqrt{f_0^2 + 4\sigma_f^2}\big)\) = 5.46 MHz, \(\sigma_f = 1/(2\pi\sigma)\) = 1.59 MHz)
d/dt is a high-pass filter with a constant +90° phase.
Step 31 · Preview: the analytic signal keeps one arrow
A real signal is pairs of counter-rotating arrows (Step 16). Delete every clockwise arrow and double every counterclockwise one: the result \(x_a = x + j\mathcal{H}\{x\}\) has the original signal as its real part, the quadrature copy as its imaginary part, the envelope as its length and the instantaneous phase as its angle. The Hilbert transform \(\mathcal{H}\) multiplies positive frequencies by −j (−90°) and negative ones by +j. IQ demodulation (01.6) is this analytic signal turned backwards at \(f_d\). Module 08 develops it.
- \(x_a\)
- analytic signal of the real signal x
- \(\mathcal{H}\)
- Hilbert transform, frequency response \(-j\operatorname{sgn} f\)
- \(u(f)\)
- unit step: 1 for f > 0, 0 for f < 0
- \(f_d\)
- demodulation frequency, Hz
- \(x(t)\)
- real signal
- \(t\)
- time, s
- \(X(f),\ X_a(f)\)
- spectra of \(x\) and \(x_a\)
- \(f\)
- frequency, Hz
- \(\mathrm{IQ}(t)\)
- complex baseband signal
- \(j\)
- imaginary unit, \(j^2=-1\)
at t = 0: \(x_a\) = 0.866 + 0.5j, the I/Q of Step 17 before demodulation
The analytic signal is the single counterclockwise arrow whose shadow is your real signal.
On axis in the far field, the pressure radiated by a small baffled source is proportional to the time derivative of its surface velocity (the Rayleigh integral gives \(p = \rho_0 S\,\dot v(t - r/c)/(2\pi r)\)), so the radiated pulse is a differentiated, upward-shifted version of the face motion. Strain imaging differentiates displacement in space, with the same noise penalty (Module 11); estimators use least-squares slopes over a kernel instead of point differences for that reason.
01.10Visual summary
Fifteen pictures that carry the module, each computed from its formula: the question you should answer in five seconds, the answer with its number, and the figure to keep in mind.
01.11How experts think about it
An expert never sees a complex number as a pair of numbers with a strange rule; they see an arrow, and every operation as a move of that arrow: add slides, multiply turns and stretches, conjugate mirrors, d/dt turns by 90° and stretches by ω. A real signal is the shadow of arrows turning both ways; IQ is the arrow with the carrier spin removed; phase is delay modulo one period. Most bugs in phase-based processing are a missing 2π, a sign convention, or an unwrap that assumed small steps.
Multiplication by −1 rotates an arrow by 180°. An operation that applied twice gives −1 is a 90° rotation, which forces arrows out of the line into a plane (01.1). The powers cycle with period 4, which is why Taylor terms of \(e^{j\theta}\) sort into cos and sin, and why the DFT twiddle factors at N = 4 are just 1, −j, −1, j.
\(w z = c\,z + d\,(jz)\): a scaled copy plus a scaled perpendicular copy, i.e. the rotation matrix \(\begin{bmatrix} c & -d \\ d & c\end{bmatrix} = \lvert w\rvert R(\angle w)\) (01.5). Worked: 5 × 2.236 = 11.18 and 53.13° + 63.43° = 116.57°, matching −5 + 10j. A unit-length multiplier \(e^{j\alpha}\) only rotates: every phase shift, steering phase or demodulation is this.
b/a loses the individual signs: 3 + 4j and −3 − 4j have the same ratio. atan2(b, a) uses both signs and returns the true angle in (−π, π] (01.2). A phase-based delay estimate built on atan(Q/I) is wrong by π in half the plane, which at 5 MHz is a 100 ns (77 µm round-trip) error.
The conjugate has the negated angle, so \(z z^*\) turns by θ and back by −θ and lands on the positive real axis with length r² (01.6). Intensity, \(\lvert X(f)\rvert^2\), coherence magnitudes and power Doppler are all \(z z^*\). Cross products \(z_2 z_1^*\) carry the angle difference, which is how lag-one autocorrelation measures Doppler phase.
jf is f turned 90°, so the velocity is tangent; then \(d\lvert f\rvert^2/d\theta = j\lvert f\rvert^2 - j\lvert f\rvert^2 = 0\) and the length stays 1 (01.10). Unit speed means arc length θ after "time" θ, hence angle θ in radians. The same argument is why phasors are eigenfunctions of linear time-invariant systems: d/dt (and any LTI operation) maps \(e^{j\omega t}\) to a multiple of itself (Module 03).
Each term is the previous one turned 90° and scaled by θ/k (01.9), so the terms grow while k < θ and then shrink factorially. At θ = π the partial sums swing out to −3.93 + 3.14j before spiralling in on −1. Libraries therefore reduce the argument modulo 2π (and to a small interval by symmetry) before using a short polynomial.
Each factor turns by atan(θ/n) and stretches by \(\sqrt{1 + \theta^2/n^2}\), so the product has angle \(n\operatorname{atan}(\theta/n) \to \theta\) and length \((1 + \theta^2/n^2)^{n/2} \approx e^{\theta^2/(2n)}\) (01.11). This is exactly forward-Euler integration of \(f' = jf\): every step leaves the circle along the tangent, so amplitude grows by \(\sqrt{1 + \omega^2\Delta t^2}\) per step for any step size. Wave-equation solvers use leapfrog or other energy-conserving schemes for this reason.
\(e^{j2\pi f t}\) turns counterclockwise for f > 0 and clockwise for f < 0 (01.12). Their real parts are identical (cos is even) and their imaginary parts opposite (sin is odd). Hence a real channel cannot measure the sign of a frequency, while an IQ channel can: Doppler direction and the side of the spectrum in complex baseband both rely on it.
One arrow never stays on the real axis; two arrows of equal length turning opposite ways keep their vertical parts cancelled at every instant (01.14). For \(A\cos(2\pi f_0 t + \varphi)\) the halves are \((A/2)e^{\pm j\varphi}\) at \(\pm f_0\) (01.15). Unequal halves trace an ellipse with semi-axes \(A_+ + A_-\) and \(\lvert A_+ - A_-\rvert\), i.e. a complex signal. The FFT of a real cosine of amplitude A therefore shows two peaks of A/2 (A·N/2 in unnormalized bins).
\(A\cos(\omega t + \varphi) = I\cos\omega t + Q\cos(\omega t + 90^\circ)\) (01.16). Going back T/4 multiplies the phasor by −j, and \(\operatorname{Re}\{-j(I + jQ)\} = Q\) (01.17). Worked: cos(ωt + 30°) gives I = x(0) = 0.866 and Q = x(−T/4) = 0.5. Hardware gets Q from a mixer with a quarter-period-delayed oscillator; nothing about it is imaginary.
Mixing shifts the \(+f_0\) band to 0 and the \(-f_0\) band to \(-2f_0\); the low-pass keeps the first (01.18). At four samples per period the mixer values are 1, −j, −1, j, and averaging one period cancels the \(2f_0\) term exactly: samples (0.866, −0.5, −0.866, 0.5) give 2 × mean = 0.866 + 0.5j = 1∠30°. With \(f_d \ne f_0\) the arrow keeps turning at \(f_0 - f_d\); the envelope is unchanged (01.20).
A delayed narrowband echo is the same slow arrow multiplied by \(e^{-j2\pi f_0\tau}\) (01.19). τ = 13 ns at 5 MHz is −23.4°, a 10 µm depth shift. Phase gives sub-nanosecond delay resolution from a single sample pair, which is why Kasai Doppler, phase-based speckle tracking and fine aberration estimates use it, and why coarse delays come from cross-correlation first.
A step of π + ε is indistinguishable from −π + ε, so unwrap chooses the smaller one. This is the sampling theorem for phase: at least two samples per turn (01.21). For round-trip tracking \(\Delta\varphi = -4\pi f_0\Delta z/c\), so the limit is λ/4 per frame; at PRF 5 kHz and 5 MHz the Doppler aliasing velocity is 0.385 m/s. Errors propagate: one wrong decision offsets every later sample by 2π.
The common \(e^{j\omega t}\) factors out, so the sum is a single sinusoid with phasor \(A_1 e^{j\varphi_1} + A_2 e^{j\varphi_2}\) (01.22). Two equal scatterers Δz apart in depth give \(\lvert\text{sum}\rvert = 2\lvert\cos(2\pi\Delta z/\lambda)\rvert\) (01.23): 2 at 0 or λ/2, zero at λ/4 = 77 µm at 5 MHz. Speckle is the same sum with many random arrows, and coherent compounding works because aligned arrows add while misaligned ones do not.
In the frame of the first arrow the second turns at Δf, so the resultant length swings between \(A_1 + A_2\) and \(\lvert A_1 - A_2\rvert\) once per 1/Δf. The product form \(2\cos(\pi\Delta f\,t)\cos(2\pi\bar f t)\) (01.24) has a 4 µs cosine factor whose magnitude repeats every 2 µs. The same mechanism, with time and frequency swapped, gives the spectral ripple (scalloping) of two echoes Δt apart, with period 1/Δt, used to size thin layers; close tones in a short DFT window interfere the same way (Module 02).
\(\frac{d}{dt}e^{j\omega t} = j\omega\,e^{j\omega t}\) (01.25). Gain doubles when f doubles (20 log₁₀ 2 = 6.02 dB), so differentiation shifts a 5 MHz, σ = 100 ns Gaussian pulse's spectral peak from 5.00 to 5.46 MHz and amplifies broadband noise: for white noise band-limited to B, the derivative variance grows as B³. Velocity from displacement, strain from displacement gradient and instantaneous frequency from phase all pay this price.
Delete the negative-frequency half and double the positive half (01.27): the real part is unchanged, the imaginary part is the 90°-shifted copy \(\mathcal{H}\{x\}\), the magnitude is the envelope and the angle is the instantaneous phase. cos(ωt + 30°) becomes \(e^{j(\omega t + 30^\circ)}\). IQ data is this arrow multiplied by \(e^{-j2\pi f_d t}\), so B-mode envelopes from analytic signals and from IQ agree exactly.
01.12Cheat sheet
| Concept | Formula | Ultrasound use |
|---|---|---|
| j | \(j^2 = -1\); ×j = +90° | Quadrature mixer, 90° phase shifts |
| Rectangular ↔ polar | \(z = a + jb = r e^{j\theta}\), \(r = \sqrt{a^2 + b^2}\), \(\theta = \operatorname{atan2}(b, a)\) | I/Q ↔ envelope and phase |
| Product | \((ac - bd) + j(ad + bc)\); \(r_1 r_2\,e^{j(\theta_1 + \theta_2)}\) | Phase rotation in IQ beamforming, steering |
| Conjugate, power | \(z^* = a - jb\); \(\lvert z\rvert^2 = z z^*\) | Intensity, \(\lvert X(f)\rvert^2\), power Doppler |
| Phase difference | \(\angle(z_2 z_1^*) = \theta_2 - \theta_1\) | Kasai lag-one autocorrelation, phase tracking |
| Division | \(z/w = z w^*/\lvert w\rvert^2\) | Deconvolution, transfer functions |
| Euler | \(e^{j\theta} = \cos\theta + j\sin\theta\); \(e^{j\pi} = -1\) | Every phasor, DFT kernel, steering vector |
| Rotating phasor | \(A e^{j(2\pi f t + \varphi)}\), f in turns/s | Carrier, Doppler, plane-wave phase |
| Real sinusoid | \(\cos\theta = \frac{e^{j\theta} + e^{-j\theta}}{2}\); \(X(-f) = X^*(f)\) | Why RF spectra are two-sided |
| Quadrature | \(A\cos(\omega t + \varphi) = I\cos\omega t + Q\cos(\omega t + 90^\circ)\) | Meaning of the Q channel |
| Quarter-period sample | \(x(t_0 - T/4) = Q\) (narrowband) | Quadrature sampling at \(f_s = 4f_0\) |
| Demodulation | \(\mathrm{IQ} = 2\operatorname{LPF}\{x\,e^{-j2\pi f_d t}\}\) | Baseband IQ data, decimation |
| Delay ↔ phase | \(\Delta\angle = -2\pi f_0\tau\); \(\Delta z = c\tau/2\) | Displacement, Doppler, fine aberration delays |
| Unwrap limit | \(\lvert\Delta\varphi\rvert < \pi \Leftrightarrow \Delta z < \lambda/4\) per sample | Aliasing velocity \(c\cdot\mathrm{PRF}/(4f_0)\) |
| Phasor sum | \(\sum_k A_k e^{j\varphi_k}\) | Speckle, coherent compounding, interference |
| Beats | \(\cos\omega_1 t + \cos\omega_2 t = 2\cos(\Delta\omega t/2)\cos(\bar\omega t)\) | Thin-layer scalloping, close-tone leakage |
| Derivative | \(d/dt \leftrightarrow \times j2\pi f\) (+90°, +6 dB/oct) | Radiated pressure ∝ dv/dt, strain noise |
| Analytic signal | \(x_a = x + j\mathcal{H}\{x\}\), \(X_a = 2u(f)X\) | Envelope detection, instantaneous phase |
01.13Interview and qualifying-exam questions
I1What does j mean geometrically, and why is j² = −1 not mysterious? · screen
j is the operator "rotate by +90°". Multiplying by −1 is a 180° rotation, so \(j^2 = -1\) says two quarter turns make a half turn. The plane is forced because no arrow on the real line can turn 90°. Follow-up: "So what is 1/j?" A quarter turn back, −j, since \(j(-j) = 1\).
I2Compute (3 + 4j)(1 + 2j) and explain the result in polar form. · screen
\(3 + 6j + 4j + 8j^2 = -5 + 10j\). Polar: \(\lvert z\rvert = 5\), \(\lvert w\rvert = 2.236\), product 11.18; angles 53.13° + 63.43° = 116.57°, and atan2(10, −5) = 116.57°. Multiplication multiplies lengths and adds angles because \(wz = c\,z + d\,(jz)\) is a rotation matrix times \(\lvert w\rvert\). Follow-up: "And the quotient?" \((3 + 4j)(1 - 2j)/5 = 2.2 - 0.4j\), length 2.236, angle −10.30°.
I3Why does the FFT of a real cosine show two peaks? What is a negative frequency? · screen
\(\cos\theta = (e^{j\theta} + e^{-j\theta})/2\): two half-amplitude phasors turning in opposite directions, whose imaginary parts cancel at every instant. A negative frequency is clockwise rotation, not time reversal. For any real signal \(X(-f) = X^*(f)\), so the negative half is redundant. Follow-up: "Amplitude of each peak for A = 1, N points, no normalization?" N/2 each, when the tone is exactly on a bin.
I4What is the physical meaning of the Q channel, and how is it measured? · onsite
\(A\cos(\omega t + \varphi) = I\cos\omega t + Q\cos(\omega t + 90^\circ)\), with \(I + jQ = A e^{j\varphi}\). Q is the amount of the 90°-shifted reference in the signal. It is measured with a second mixer whose oscillator is shifted by a quarter period (\(Q = -2\operatorname{LPF}\{x\sin\omega_d t\}\)), or, for a narrowband signal, by sampling the RF a quarter period earlier: \(x(t_0 - T/4) = Q\). Follow-up: "When does the quarter-period trick fail?" When the envelope changes within T/4: about 20 % rms error at 60 % bandwidth for a 5 MHz Gaussian echo.
I5Derive Euler's formula without the Taylor series. · onsite
Let \(f(\theta) = e^{j\theta}\), so \(f' = jf\), \(f(0) = 1\). Then \(d\lvert f\rvert^2/d\theta = f'f^* + f f'^* = j\lvert f\rvert^2 - j\lvert f\rvert^2 = 0\), so \(\lvert f\rvert = 1\); the speed \(\lvert f'\rvert = 1\), so the point covers arc length θ on the unit circle starting at 1, counterclockwise because the velocity at 1 is j. Its coordinates are (cos θ, sin θ). Alternatively, \((1 + j\theta/n)^n\) is n rotations by atan(θ/n) with total stretch \(e^{\theta^2/(2n)} \to 1\). Follow-up: "What does the limit tell you about forward Euler on an oscillator?" It spirals outward by \(\sqrt{1 + \omega^2\Delta t^2}\) per step for any Δt.
I6You demodulate at fd = 4.5 MHz but the echo is centred at 5 MHz. What changes in I, Q and the envelope? · onsite
\(\mathrm{IQ}(t) = A(t - T_e)\,e^{j[2\pi(f_0 - f_d)t - 2\pi f_0 T_e]}\): the envelope \(\lvert\mathrm{IQ}\rvert\) is unchanged (as long as the band stays inside the low-pass), but the arrow keeps turning at +0.5 MHz, one turn per 2 µs. I and Q oscillate under the envelope; phase-based estimators see a frequency offset. Follow-up: "Why does depth matter?" Attenuation lowers the echo centre frequency with depth, so a fixed \(f_d\) creates a depth-dependent residual rotation that biases phase-to-delay conversions.
I7Two equal point scatterers lie in one resolution cell, separated axially by Δz. When do they cancel at 5 MHz? · onsite
Round-trip phase difference \(4\pi\Delta z/\lambda\); \(\lvert\text{sum}\rvert = 2\lvert\cos(2\pi\Delta z/\lambda)\rvert\), zero at Δz = λ/4 = 77 µm (λ = 0.308 mm), maximum at multiples of λ/2 = 154 µm. Follow-up: "So what is speckle?" The sum of many such arrows with random phases: a 2-D random walk whose length is Rayleigh distributed (Module 07).
I8In phase-based displacement estimation, what limits the displacement per frame, and how do you extend it? · research
The phase step must stay below π, i.e. \(\lvert\Delta z\rvert < \lambda/4\) per frame (77 µm at 5 MHz), otherwise unwrapping chooses the wrong turn and the error propagates. Extensions: higher frame rate; lower \(f_0\) for the coarse estimate; coarse delay from cross-correlation of envelopes or RF and fine delay from phase; spatial or temporal continuity priors for 2-D/3-D unwrapping; or using the phase of a lower difference frequency (two-frequency methods). Follow-up: "What does noise do?" A phase noise σ near 1 rad produces steps beyond π by itself, so unwrap fails even for slow motion.
I9Your I and Q channels have a 1 % gain mismatch and a 1° phase error. What artefact appears in spectral Doppler, and how large is it? · research
\(z = \cos\theta + j g\sin(\theta + \delta)\) \(= e^{j\theta}\frac{1 + g e^{j\delta}}{2} + e^{-j\theta}\frac{1 - g e^{-j\delta}}{2}\). The second term is a mirror image at \(-f_D\) with relative level \(\lvert 1 - g e^{-j\delta}\rvert / \lvert 1 + g e^{j\delta}\rvert\) = −40.0 dB for g = 1.01, δ = 1° (−26 dB for 5 % and 5°). Geometrically the counter-rotating halves are unequal in the wrong direction, so the IQ circle becomes an ellipse. Follow-up: "How do you correct it?" Estimate g and δ from a calibration tone and apply the inverse 2×2 transform to (I, Q).
I10Why do differentiation-based estimates (velocity from displacement, strain, instantaneous frequency) amplify noise, and by how much? · research
d/dt multiplies each component by \(j2\pi f\), so the noise power spectrum is weighted by \((2\pi f)^2\). For white noise of density \(N_0\) band-limited to ±B, the derivative variance is \(\int_{-B}^{B} N_0 (2\pi f)^2\,df = N_0 (2\pi)^2\,\frac{2B^3}{3}\), growing as B³ while the signal's derivative grows only with its own bandwidth. Strain estimators therefore fit a least-squares slope over a kernel (a low-pass differentiator), trading resolution for SNR. Follow-up: "Discrete first difference?" \(\lvert H\rvert = 2 f_s\lvert\sin(\pi f/f_s)\rvert\) with a half-sample delay.
I11How does a time delay appear on IQ data, and why must you rotate after interpolating IQ samples? · onsite
RF delayed by τ has analytic signal \(A(t - \tau)e^{j2\pi f_0(t - \tau)}\). After demodulation at \(f_d\), \(\mathrm{IQ}_{\text{delayed}}(t) = A(t - \tau)e^{j2\pi f_0(t - \tau)}e^{-j2\pi f_d t}\) \(= \mathrm{IQ}(t - \tau)\,e^{-j2\pi f_d\tau}\). A beamformer undoes the delay: interpolating a channel at \(t + \tau\) realigns the envelope but not the carrier phase, which still carries the factor \(e^{-j2\pi f_d\tau}\). So the interpolated sample must be multiplied by \(e^{+j2\pi f_d\tau}\) (the sign is for demodulation with \(e^{-j2\pi f_d t}\); it flips with the opposite convention). Without it, channels add with phase errors of \(2\pi f_d\tau\) and the beam is destroyed. Follow-up: "Which frequency, \(f_0\) or \(f_d\)?" \(f_d\), the frequency actually removed (Module 08).
01.14MATLAB lab
Base MATLAB only (no toolboxes); each script runs as-is in MATLAB R2019b+ and GNU Octave 8. Every hand calculation on this page is printed by one of these scripts.
m01_complex_basics.mRectangular/polar, add, multiply, conjugate, divide, powers of j, atan vs atan2, drawn as arrows · change z and w · runtime ≈ 1 s
%% m01_complex_basics.m -- complex numbers as arrows
% Purpose : check by computer every hand calculation of Module 01, sections 01.1-01.2:
% rectangular <-> polar, addition (tip to tail), multiplication
% (rotate and scale), conjugate (mirror), |z|^2 = z z*, division,
% powers of j, and the atan vs atan2 pitfall.
% Change : z and w below; try w = 1j (pure quarter turn), w = -1 (half turn),
% w = 2 (pure scale), w = exp(1j*pi/6) (pure 30 deg rotation).
% Expect : |z| = 5, angle 53.13 deg; z*w = -5 + 10j, |zw| = 11.180 = 5 x 2.236,
% angle 116.57 = 53.13 + 63.43 deg; z/w = 2.2 - 0.4j; z z* = 25;
% atan(4/3) = 53.13 deg but the true angle of -3-4j is -126.87 deg.
% Runtime : < 1 s. Base MATLAB / GNU Octave, no toolboxes.
clear; close all; rng(1);
%% 1. One number, two descriptions
z = 3 + 4j; % rectangular: 3 units right, 4 units up
w = 1 + 2j;
r_z = abs(z); % length sqrt(3^2 + 4^2)
th_z = angle(z); % angle atan2(4, 3), in rad, range (-pi, pi]
fprintf('z = %g %+gj -> |z| = %.4f, angle = %.4f rad = %.2f deg\n', real(z), imag(z), r_z, th_z, th_z*180/pi);
fprintf('w = %g %+gj -> |w| = %.4f, angle = %.2f deg\n', real(w), imag(w), abs(w), angle(w)*180/pi);
z_back = r_z*cos(th_z) + 1j*r_z*sin(th_z); % polar -> rectangular
fprintf('back to rectangular: %.4f %+.4fj (Euler: r e^{j theta})\n\n', real(z_back), imag(z_back));
%% 2. Arithmetic by hand vs built-in
a = real(z); b = imag(z); c = real(w); d = imag(w);
s_hand = (a + c) + 1j*(b + d); % add components (tip to tail)
p_hand = (a*c - b*d) + 1j*(a*d + b*c); % FOIL with j^2 = -1
q_hand = ((a*c + b*d) + 1j*(b*c - a*d)) / (c^2 + d^2); % z w* / |w|^2
fprintf('z + w = %g %+gj (built-in %g %+gj)\n', real(s_hand), imag(s_hand), real(z+w), imag(z+w));
fprintf('z * w = %g %+gj (built-in %g %+gj)\n', real(p_hand), imag(p_hand), real(z*w), imag(z*w));
fprintf(' |zw| = %.4f = |z||w| = %.4f ; angle(zw) = %.2f = %.2f + %.2f deg\n', ...
abs(z*w), abs(z)*abs(w), angle(z*w)*180/pi, angle(z)*180/pi, angle(w)*180/pi);
fprintf('z / w = %g %+gj (built-in %g %+gj)\n', real(q_hand), imag(q_hand), real(z/w), imag(z/w));
fprintf(' |z/w| = %.4f = |z|/|w| ; angle = %.2f deg = 53.13 - 63.43\n', abs(z/w), angle(z/w)*180/pi);
fprintf('conj(z) = %g %+gj (mirror in the real axis); z*conj(z) = %g %+gj = |z|^2\n', ...
real(conj(z)), imag(conj(z)), real(z*conj(z)), imag(z*conj(z)));
fprintf('but z^2 = %g %+gj (squaring is not |z|^2)\n\n', real(z^2), imag(z^2));
%% 3. j is the quarter-turn operator
for k = 0:4
v = [round(real(1j^k)), round(imag(1j^k))];
v(v == 0) = 0; % print +0, not -0
fprintf('j^%d = %+g %+gj (angle %4.0f deg)\n', k, v(1), v(2), k*90);
end
fprintf('1/j = %g %+gj (a quarter turn backwards)\n', real(1/1j), imag(1/1j));
fprintf('j*z = %g %+gj : z turned by +90 deg, same length %.1f\n\n', real(1j*z), imag(1j*z), abs(1j*z));
%% 4. atan vs atan2
zz = -3 - 4j;
fprintf('z = -3-4j: atan(b/a) = %.2f deg (wrong quadrant), atan2(b,a) = %.2f deg, angle() = %.2f deg\n', ...
atan(imag(zz)/real(zz))*180/pi, atan2(imag(zz), real(zz))*180/pi, angle(zz)*180/pi);
%% 5. Figures: the same algebra drawn as arrows
hd = @(p0, p1, s, sg) p1 - s*(p1 - p0)/abs(p1 - p0)*exp(sg*0.45j); % arrowhead corner
arrow = @(p0, p1, col) plot(real([p0 p1 NaN hd(p0,p1,0.45,1) p1 hd(p0,p1,0.45,-1)]), ...
imag([p0 p1 NaN hd(p0,p1,0.45,1) p1 hd(p0,p1,0.45,-1)]), 'Color', col, 'LineWidth', 2);
uc = exp(1j*linspace(0, 2*pi, 200));
figure('Position', [80 80 1100 360]);
subplot(1,3,1); hold on; grid on;
arrow(0, z, [0.12 0.37 0.82]); arrow(z, z+w, [0.76 0.21 0.44]); arrow(0, z+w, [0.06 0.48 0.39]);
axis equal; axis([-1 8 -1 8]); xlabel('real part'); ylabel('imaginary part');
title('z + w: tip to tail (4 + 6j)');
subplot(1,3,2); hold on; grid on;
plot(real(uc), imag(uc), ':', 'Color', [0.5 0.5 0.5]);
arrow(0, z, [0.12 0.37 0.82]); arrow(0, w, [0.76 0.21 0.44]); arrow(0, z*w, [0.06 0.48 0.39]);
axis equal; axis([-7 6 -1 11]); xlabel('real part'); ylabel('imaginary part');
title('z w: lengths multiply, angles add (-5 + 10j)');
subplot(1,3,3); hold on; grid on;
arrow(0, z, [0.12 0.37 0.82]); arrow(0, conj(z), [0.66 0.43 0.02]); arrow(0, z/w, [0.06 0.48 0.39]);
plot([-1 6], [0 0], 'k-');
axis equal; axis([-1 6 -5 5]); xlabel('real part'); ylabel('imaginary part');
title('z* = mirror (3 - 4j); z / w = 2.2 - 0.4j');
What to look for: |z| = 5, 53.13°; z·w = −5 + 10j with |zw| = 11.1803 = |z||w| and 116.57° = 53.13° + 63.43°; z/w = 2.2 − 0.4j; z z* = 25 but z² = −7 + 24j; j⁰…j⁴ cycle +1, j, −1, −j, +1; atan gives 53.13° for −3 − 4j while atan2 gives −126.87°.
m01_euler_three_ways.mTaylor partial sums, the ODE f′ = j f, and the limit (1 + jθ/n)ⁿ = forward Euler · change θ and n · runtime ≈ 1 s
%% m01_euler_three_ways.m -- e^{j theta} = cos(theta) + j sin(theta), three ways
% Purpose : verify Euler's formula numerically by the three routes of section 01.3:
% (1) Taylor series: partial sums of sum (j theta)^k / k!
% (2) ODE: f' = j f means velocity is position turned by 90 deg, so |f| stays 1
% (3) limit: (1 + j theta/n)^n = n small rotate-and-scale steps
% (this is exactly forward-Euler integration of f' = j f with step theta/n)
% Change : theta (try pi, 2*pi, 10); the list of n.
% Expect : theta = 1: 8 terms give error 2.5e-5; theta = pi: 13 terms for error < 1e-3.
% Limit at theta = pi: n = 4 -> |.| = 2.614, angle 152.6 deg;
% n = 100 -> |.| = 1.0506 (close to exp(theta^2/(2n)) = 1.0506), angle 179.94 deg.
% Runtime : < 1 s. Base MATLAB / GNU Octave, no toolboxes.
clear; close all; rng(1);
%% 1. Taylor series, term by term
for theta = [1, pi]
exact = exp(1j*theta);
S = 0;
fprintf('Taylor, theta = %.4f rad (exact %.4f %+.4fj)\n', theta, real(exact), imag(exact));
fprintf(' terms term (j theta)^k/k! partial sum |error|\n');
for k = 0:13
term = (1j*theta)^k / factorial(k); % powers of j cycle 1, j, -1, -j
S = S + term;
fprintf(' %4d %+9.4f %+9.4fj %+8.4f %+8.4fj %9.2e\n', k+1, real(term), imag(term), real(S), imag(S), abs(S - exact));
end
end
fprintf('\n');
%% 2. ODE: velocity is perpendicular to position
theta = linspace(0, 2*pi, 2001);
f = exp(1j*theta);
h = theta(2) - theta(1);
df = (f(3:end) - f(1:end-2)) / (2*h); % central difference d f / d theta
fm = f(2:end-1);
fprintf('ODE check: max |df/dtheta - j f| = %.2e (finite-difference error ~ h^2/6)\n', max(abs(df - 1j*fm)));
fprintf(' max |Re(conj(f) df)| = %.2e (velocity _|_ position, so |f| is constant)\n', max(abs(real(conj(fm).*df))));
fprintf(' max ||f| - 1| = %.2e ; speed |df| = %.6f rad per rad\n\n', max(abs(abs(f) - 1)), mean(abs(df)));
%% 3. Limit (1 + j theta / n)^n and forward Euler
theta = pi;
ns = [1 2 4 8 16 32 100 1000 10000];
lim_err = zeros(size(ns));
fprintf('Limit, theta = pi\n n |(1+j pi/n)^n| exp(pi^2/2n) angle (deg) |error|\n');
for i = 1:numel(ns)
n = ns(i);
v = (1 + 1j*theta/n)^n;
y = 1; for k = 1:n, y = y + (theta/n)*1j*y; end % forward Euler on f' = j f
lim_err(i) = abs(v - exp(1j*theta));
fprintf('%6d %12.4f %12.4f %10.2f %9.2e (Euler-step diff %.1e)\n', n, abs(v), exp(theta^2/(2*n)), angle(v)*180/pi, lim_err(i), abs(y - v));
end
fprintf('(1 + j/2)^2 by hand = 1 + j - 1/4 = %.4f %+.4fj\n', real((1 + 0.5j)^2), imag((1 + 0.5j)^2));
%% 4. Figures
uc = exp(1j*linspace(0, 2*pi, 400));
figure('Position', [80 80 1100 380]);
subplot(1,3,1); hold on; grid on;
plot(real(uc), imag(uc), ':', 'Color', [0.5 0.5 0.5]);
S = cumsum((1j*pi).^(0:14) ./ factorial(0:14));
plot(real([0 S]), imag([0 S]), '-o', 'Color', [0.12 0.37 0.82], 'MarkerSize', 4, 'LineWidth', 1.5);
plot(-1, 0, 'p', 'MarkerSize', 12, 'Color', [0.76 0.21 0.44]);
axis equal; axis([-4.5 1.5 -2.5 3.5]); xlabel('real part'); ylabel('imaginary part');
title('Taylor partial sums, \theta = \pi -> -1');
subplot(1,3,2); hold on; grid on;
plot(real(uc), imag(uc), ':', 'Color', [0.5 0.5 0.5]);
cols = [0.76 0.21 0.44; 0.66 0.43 0.02; 0.12 0.37 0.82];
nlist = [4 8 32];
for i = 1:3
n = nlist(i); pts = (1 + 1j*pi/n).^(0:n);
plot(real(pts), imag(pts), '-o', 'Color', cols(i,:), 'MarkerSize', 3, 'LineWidth', 1.3);
end
axis equal; axis([-2.6 1.4 -0.6 3]); xlabel('real part'); ylabel('imaginary part');
legend('unit circle', 'n = 4', 'n = 8', 'n = 32', 'Location', 'northwest');
title('(1 + j\pi/n)^k, k = 0..n');
subplot(1,3,3);
Nt = 1:20; terr = zeros(size(Nt));
for i = Nt, terr(i) = abs(sum((1j*pi).^(0:i-1) ./ factorial(0:i-1)) - (-1)); end
loglog(Nt, terr, 'o-', 'Color', [0.12 0.37 0.82], 'LineWidth', 1.5); hold on; grid on;
loglog(ns, lim_err, 's-', 'Color', [0.76 0.21 0.44], 'LineWidth', 1.5);
loglog(ns, pi^2./(2*ns), 'k--');
xlabel('number of terms N or steps n'); ylabel('|error| at \theta = \pi');
legend('Taylor, N terms', 'limit, n steps', '\pi^2/(2n)', 'Location', 'southwest');
title('Taylor converges factorially, the limit like 1/n');
What to look for: at θ = 1 the partial sums reproduce the hand table (0.5 + 0.8333j after 4 terms) and 8 terms leave \(2.47\times10^{-5}\); at θ = π 13 terms leave \(4.56\times10^{-4}\). The ODE check gives \(\operatorname{Re}(f^* f') \approx 10^{-14}\) (velocity ⟂ position). The limit at θ = π has length 2.6142 at n = 4 and 1.0506 at n = 100, matching \(e^{\pi^2/(2n)}\), and is identical to forward Euler to machine precision.
m01_phasors_and_iq.mPhasor direction, cos as two halves, Q = quarter-period-earlier sample, I/Q by hand at 4f₀, demodulation of a 5 MHz echo · change τ, fd, bandwidth · runtime ≈ 1 s
%% m01_phasors_and_iq.m -- rotating phasors, the imaginary part, and I/Q of an RF echo
% Purpose : sections 01.4-01.6 by computer.
% (1) e^{j 2 pi f t}: angle at a given time, the sign of f flips sin but not cos
% (2) cos = (e^{+j} + e^{-j})/2: the imaginary parts cancel to machine precision
% (3) the imaginary part (Q) of a narrowband signal = the real signal a quarter period earlier
% (4) I/Q by hand at fs = 4 f0 (multiply by 1, -j, -1, j, average, x2)
% (5) I/Q demodulation of a Gaussian RF echo: mix with e^{-j 2 pi fd t}, low-pass, x2;
% phase at the peak = -2 pi f0 tau, envelope = |IQ|; residual rotation if fd ~= f0
% Change : tau (echo delay), fd (demodulation frequency), bw (fractional -6 dB bandwidth).
% Expect : angle at 50 ns = 90 deg; hand I/Q = 0.8660 + 0.5000j = 1 at 30 deg;
% tau = 13 ns -> phase -23.4 deg; fd = 4.5 MHz -> IQ spins at +0.5 MHz.
% Runtime : < 2 s. Base MATLAB / GNU Octave, no toolboxes.
clear; close all; rng(1);
%% 1. A rotating arrow
f0 = 5e6; % carrier (Hz) = turns per second
t1 = 50e-9; % 50 ns = a quarter of the 200 ns period
p_pos = exp(1j*2*pi*f0*t1);
p_neg = exp(-1j*2*pi*f0*t1);
fprintf('f0 = 5 MHz, period %.0f ns. At t = 50 ns the angle is %.1f deg\n', 1e9/f0, 2*pi*f0*t1*180/pi);
fprintf(' +f: %+.3f %+.3fj ; -f: %+.3f %+.3fj (same cos, opposite sin)\n\n', real(p_pos), imag(p_pos), real(p_neg), imag(p_neg));
%% 2. Two counter-rotating halves make a real cosine
t = (0:999)*1e-10; % 0 .. 100 ns, 0.1 ns steps
x = 0.5*exp(1j*2*pi*f0*t) + 0.5*exp(-1j*2*pi*f0*t);
fprintf('max |imag( e^{+j}/2 + e^{-j}/2 )| = %.1e ; max |real - cos| = %.1e\n', max(abs(imag(x))), max(abs(real(x) - cos(2*pi*f0*t))));
x_unbal = 0.5*exp(1j*2*pi*f0*t) + 0.25*exp(-1j*2*pi*f0*t); % unbalanced: no longer real
fprintf('unbalanced (0.5, 0.25): max |imag| = %.3f = 0.5 - 0.25\n\n', max(abs(imag(x_unbal))));
%% 3. The imaginary part is the quarter-period-earlier copy
A = 1; phi = 30*pi/180; T0 = 1/f0;
xr = @(tt) A*cos(2*pi*f0*tt + phi); % the real RF
I = A*cos(phi); Q = A*sin(phi); % phasor A e^{j phi} = I + jQ
fprintf('A cos(wt + 30 deg): I = A cos(phi) = %.4f, Q = A sin(phi) = %.4f\n', I, Q);
fprintf(' x(0) = %.4f (= I); x(-T/4) = %.4f (= Q)\n', xr(0), xr(-T0/4));
fprintf(' and x(t) = I cos(wt) + Q cos(wt + 90 deg): error %.1e\n\n', max(abs(xr(t) - (I*cos(2*pi*f0*t) + Q*cos(2*pi*f0*t + pi/2)))));
%% 4. I/Q by hand at fs = 4 f0
n = 0:3;
s = cos(n*pi/2 + phi); % four samples, one per quarter period
lo = exp(-1j*n*pi/2); % e^{-j 2 pi f0 t_n} = 1, -j, -1, j
prod_ = s .* lo;
iq = 2*mean(prod_);
fprintf('fs = 4 f0 samples : %s\n', sprintf('%+.4f ', s));
fprintf('times 1, -j, -1, j : %s\n', sprintf('%+.4f%+.4fj ', [real(prod_); imag(prod_)]));
fprintf('2 x mean = %.4f %+.4fj -> |IQ| = %.4f, angle %.2f deg\n\n', real(iq), imag(iq), abs(iq), angle(iq)*180/pi);
%% 5. I/Q demodulation of an RF echo
fs = 40e6; % sampling rate (Hz)
N = 256; % 6.4 us record
tt = (0:N-1)/fs;
t0 = 3.2e-6; % reference echo time: 16 whole periods of f0
bw = 0.6; % -6 dB fractional bandwidth (60 %)
sig = sqrt(log(2)/2)/(pi*bw*f0/2);% Gaussian envelope std (s): 125 ns for 60 %
tau = 13e-9; % extra delay (s)
rf = exp(-(tt - t0 - tau).^2/(2*sig^2)) .* cos(2*pi*f0*(tt - t0 - tau));
fk = [0:N/2-1, -N/2:-1]*fs/N; % FFT frequency axis (Hz)
lpf = double(abs(fk) < 5e6); % ideal low-pass, cutoff 5 MHz (halfway to the 2f0 image)
demod = @(fd) 2*ifft(fft(rf .* exp(-1j*2*pi*fd*tt)) .* lpf);
iq0 = demod(f0);
[~, ipk] = max(abs(iq0));
fprintf('Gaussian echo: f0 = 5 MHz, 60 %% bandwidth (sigma = %.0f ns), fs = 40 MHz, tau = %.0f ns\n', sig*1e9, tau*1e9);
fprintf(' fd = f0 : phase at the peak = %.2f deg ; predicted -360 f0 tau = %.2f deg\n', angle(iq0(ipk))*180/pi, -360*f0*tau);
fprintf(' peak |IQ| = %.4f (envelope peak 1)\n', abs(iq0(ipk)));
fd2 = 4.5e6; iq2 = demod(fd2);
win = abs(tt - t0) < 1.5*sig; % inside the echo
ph2 = unwrap(angle(iq2(win)));
rate = (ph2(end) - ph2(1)) / (2*pi*(sum(win) - 1)/fs); % turns per second
fprintf(' fd = 4.5 MHz: IQ keeps spinning at %.3f MHz (= f0 - fd); max envelope change %.1e\n', rate/1e6, max(abs(abs(iq2) - abs(iq0))));
% quarter-period trick on the broadband echo (no Hilbert filter), error relative to the analytic signal
for bwq = [0.1 0.6] % fractional -6 dB bandwidth
sq = sqrt(log(2)/2)/(pi*bwq*f0/2); % envelope std (s)
rq = exp(-(tt - t0).^2/(2*sq^2)) .* cos(2*pi*f0*(tt - t0));
xaq = ifft(fft(rq) .* (1 + sign(fk))); % analytic signal (keep +f, x2)
xqq = interp1(tt, rq, tt - T0/4, 'spline', 0); % real signal a quarter period earlier
fprintf(' quarter-period copy vs Hilbert at %2.0f %% bandwidth: rms error %.1f %% of the analytic-signal rms\n', 100*bwq, 100*sqrt(sum(abs(rq + 1j*xqq - xaq).^2)/sum(abs(xaq).^2)));
end
%% 6. Figures
figure('Position', [80 80 1100 380]);
m = tt > 2.6e-6 & tt < 3.8e-6;
subplot(1,3,1);
plot(tt(m)*1e6, rf(m), 'o-', 'Color', [0.12 0.37 0.82], 'MarkerSize', 3); hold on; grid on;
plot(tt(m)*1e6, abs(iq0(m)), 'k--', tt(m)*1e6, -abs(iq0(m)), 'k--');
xlabel('time (\mus)'); ylabel('amplitude'); title('RF echo at 40 MHz and envelope |IQ|');
subplot(1,3,2);
plot(tt(m)*1e6, real(iq0(m)), '-', 'Color', [0.12 0.37 0.82], 'LineWidth', 1.5); hold on; grid on;
plot(tt(m)*1e6, imag(iq0(m)), '-', 'Color', [0.76 0.21 0.44], 'LineWidth', 1.5);
plot(tt(m)*1e6, real(iq2(m)), ':', 'Color', [0.12 0.37 0.82]);
legend('I (fd = f0)', 'Q (fd = f0)', 'I (fd = 4.5 MHz)', 'Location', 'southwest');
xlabel('time (\mus)'); ylabel('amplitude'); title('I and Q: the slow arrow');
subplot(1,3,3);
plot(real(iq0(m)), imag(iq0(m)), '-', 'Color', [0.12 0.37 0.82], 'LineWidth', 3); hold on; grid on;
plot(real(iq2(m)), imag(iq2(m)), '-', 'Color', [0.66 0.43 0.02]);
plot(real(iq0(ipk)), imag(iq0(ipk)), 'ko', 'MarkerSize', 7);
axis equal; axis([-1.1 1.1 -1.1 1.1]); xlabel('I'); ylabel('Q');
legend('fd = f0: one fixed direction', 'fd = 4.5 MHz: loop (residual turn)', 'Location', 'southoutside');
title('IQ plane: phase -23.4 deg at the peak');
What to look for: 90.0° at 50 ns; imaginary residue 0 for balanced halves and 0.25 for (0.5, 0.25); x(0) = 0.8660 = I and x(−T/4) = 0.5000 = Q; the hand I/Q 0.8660 + 0.5000j; the echo phase −23.40° for τ = 13 ns (predicted −23.40°); residual rotation 0.501 MHz for \(f_d\) = 4.5 MHz; quarter-period copy vs Hilbert 3.3 % rms at 10 % bandwidth and 19.8 % at 60 %.
m01_phase_and_derivative.mPhasor sums, two-scatterer interference, beats, phase unwrap failure, d/dt = j2πf by FFT · change the phasors, Δf, the step growth a · runtime ≈ 1 s
%% m01_phase_and_derivative.m -- adding phasors, beats, phase unwrapping, d/dt = j omega
% Purpose : sections 01.7-01.9 by computer.
% (1) equal-frequency sinusoids add as arrows (phasor sum vs brute-force time sum)
% (2) two scatterers dz apart in depth: |1 + e^{-j 4 pi dz / lambda}| = 2|cos(2 pi dz/lambda)|
% (3) beats: 5.0 + 5.5 MHz, envelope period 1/df = 2 us
% (4) unwrap works only while the true phase step per sample is below pi
% (5) differentiation multiplies each phasor by j omega: +90 deg and +6.02 dB per octave
% Change : the phasors in (1), df in (3), the acceleration a in (4).
% Expect : cos + cos(+90) = 1.4142 cos(wt + 45 deg); 1 + 0.5 at 120 deg = 0.8660 at 30 deg;
% cancellation at dz = lambda/4 = 77 um (5 MHz); beat period 2.00 us;
% unwrap first fails at sample 32 for a = 0.1 rad/sample^2; derivative gain 6.02 dB per octave;
% Gaussian pulse spectral peak 5.00 -> 5.46 MHz after d/dt.
% Runtime : < 2 s. Base MATLAB / GNU Octave, no toolboxes.
clear; close all; rng(1);
%% 1. Same frequency: add the arrows
f0 = 5e6; c = 1540; lambda = c/f0; % 0.308 mm
t = (0:3999)*1e-10; % 0 .. 400 ns
pairs = {[1, 0; 1, 90], [1, 0; 0.5, 120], [1, 0; 1, 120; 1, 240], [1, 0; 1, 180]};
for k = 1:numel(pairs)
P = pairs{k}; % rows: [amplitude, phase in deg]
S = sum(P(:,1) .* exp(1j*P(:,2)*pi/180)); % phasor sum
xsum = zeros(size(t));
for r = 1:size(P,1), xsum = xsum + P(r,1)*cos(2*pi*f0*t + P(r,2)*pi/180); end
err = max(abs(xsum - abs(S)*cos(2*pi*f0*t + angle(S))));
if abs(S) < 1e-12, S = 0; end % exact cancellation: the angle is undefined
fprintf('phasors %-28s -> sum %.4f at %7.2f deg (time-domain check %.1e)\n', mat2str(P), abs(S), angle(S)*180/pi, err);
end
%% 2. Two scatterers in one resolution cell
dz = [0 0.25 0.5 1] * lambda/2; % axial separations
amp = abs(1 + exp(-1j*4*pi*dz/lambda)); % round trip: phase 4 pi dz / lambda
fprintf('\nlambda = %.3f mm. Two equal echoes, separation dz:\n', lambda*1e3);
for i = 1:numel(dz)
fprintf(' dz = %6.1f um (%4.2f lambda) -> |sum| = %.3f\n', dz(i)*1e6, dz(i)/lambda, amp(i));
end
%% 3. Beats
fs = 200e6; tb = (0:1999)/fs; % 10 us: whole periods of both tones
df = 0.5e6;
xb = cos(2*pi*f0*tb) + cos(2*pi*(f0 + df)*tb);
fk = [0:999, -1000:-1]*fs/2000;
xa = ifft(fft(xb) .* (1 + sign(fk))); % analytic signal
env = abs(xa);
mins = find(env(2:end-1) < env(1:end-2) & env(2:end-1) < env(3:end) & env(2:end-1) < 0.2) + 1;
fprintf('\nbeats 5.0 + 5.5 MHz: envelope minima spacing %.3f us (1/df = %.3f us); max %.3f, min %.3f\n', ...
mean(diff(tb(mins)))*1e6, 1e6/df, max(env(200:end-200)), min(env(200:end-200)));
%% 4. Phase wrapping and unwrapping
a = 0.1; % rad per sample^2 : phase step grows by a each sample
n = 0:59;
ph_true = 0.5*a*n.^2; % step between samples n-1 and n = a(n - 0.5)
ph_wrap = angle(exp(1j*ph_true)); % what angle() returns, in (-pi, pi]
ph_unw = unwrap(ph_wrap);
bad = find(abs(ph_unw - ph_true) > 1e-9, 1);
fprintf('\nunwrap: a = %.2f rad/sample^2, first wrong sample n = %d (true step there %.3f rad > pi)\n', a, n(bad), ph_true(bad) - ph_true(bad-1));
fprintf(' at 5 MHz a step of pi is a round-trip displacement of lambda/4 = %.1f um per sample\n', lambda/4*1e6);
fprintf(' wrapped values of 0, 2, 4, 6 rad: %s\n', sprintf('%.4f ', angle(exp(1j*[0 2 4 6]))));
%% 5. d/dt is multiplication by j omega
fs2 = 400e6; N = 4096; td = ((0:N-1) - N/2)/fs2;
sig = 100e-9;
g = exp(-td.^2/(2*sig^2)) .* cos(2*pi*f0*td); % Gaussian RF pulse
dg_exact = (-td/sig^2).*exp(-td.^2/(2*sig^2)).*cos(2*pi*f0*td) - 2*pi*f0*exp(-td.^2/(2*sig^2)).*sin(2*pi*f0*td);
fk2 = [0:N/2-1, -N/2:-1]*fs2/N;
dg_fft = real(ifft(1j*2*pi*fk2 .* fft(g))); % multiply each phasor by j 2 pi f
fprintf('\nFFT derivative vs exact: max error %.1e of max %.3e /s\n', max(abs(dg_fft - dg_exact)), max(abs(dg_exact)));
fprintf('|H| = 2 pi f: 5 MHz -> %.3e /s, 10 MHz -> %.3e /s, ratio %.2f dB (one octave)\n', 2*pi*5e6, 2*pi*10e6, 20*log10(2));
tc = (0:799)/fs2; % exactly 10 periods of f0
fc = [0:399, -400:-1]*fs2/800;
d1 = real(ifft(1j*2*pi*fc .* fft(cos(2*pi*f0*tc)))); % derivative of a pure cosine
fprintf('d/dt cos: phase of H at f0 = %.1f deg ; |d/dt cos(wt) - w cos(wt + 90 deg)| / w = %.1e\n', ...
angle(1j*2*pi*f0)*180/pi, max(abs(d1 - 2*pi*f0*cos(2*pi*f0*tc + pi/2)))/(2*pi*f0));
G = abs(fft(g)); D = abs(fft(dg_exact));
Np = 2^17; fkp = (0:Np/2-1)*fs2/Np; % zero-padded grid (3 kHz bins): peaks not bin-quantized
Gp = abs(fft(g, Np)); Dp = abs(fft(dg_exact, Np));
[~, ig] = max(Gp(1:Np/2)); [~, id] = max(Dp(1:Np/2));
sigf = 1/(2*pi*sig); % spectral std of the Gaussian envelope, Hz
fprintf('spectral peak: pulse %.2f MHz, its derivative %.2f MHz (the jw weighting raises it)\n', fkp(ig)/1e6, fkp(id)/1e6);
fprintf(' predicted (f0 + sqrt(f0^2 + 4 sigf^2))/2 = %.2f MHz, sigf = %.2f MHz\n', (f0 + sqrt(f0^2 + 4*sigf^2))/2/1e6, sigf/1e6);
%% 6. Figures
figure('Position', [80 80 1100 380]);
subplot(1,3,1);
plot(tb*1e6, xb, 'Color', [0.12 0.37 0.82]); hold on; grid on;
plot(tb*1e6, env, 'k', 'LineWidth', 1.5);
xlim([0 6]); xlabel('time (\mus)'); ylabel('amplitude');
title('beats: 5.0 + 5.5 MHz, envelope period 2 \mus');
subplot(1,3,2);
plot(n, ph_true, 'k-', 'LineWidth', 1.5); hold on; grid on;
plot(n, ph_wrap, 'o', 'Color', [0.76 0.21 0.44], 'MarkerSize', 4);
plot(n, ph_unw, '-', 'Color', [0.12 0.37 0.82], 'LineWidth', 1.5);
plot([n(bad) n(bad)], [-5 180], ':', 'Color', [0.66 0.43 0.02]);
xlabel('sample n'); ylabel('phase (rad)'); legend('true', 'wrapped', 'unwrapped', 'step > \pi', 'Location', 'northwest');
title('unwrap fails once the step exceeds \pi');
subplot(1,3,3);
ff = fk2(1:N/2)/1e6;
plot(ff, G(1:N/2)/max(G), 'Color', [0.12 0.37 0.82], 'LineWidth', 1.5); hold on; grid on;
plot(ff, D(1:N/2)/max(D), 'Color', [0.76 0.21 0.44], 'LineWidth', 1.5);
xlim([0 12]); xlabel('frequency (MHz)'); ylabel('normalized |spectrum|');
legend('pulse', 'd/dt pulse', 'Location', 'northeast'); title('d/dt weights each line by 2\pi f');
What to look for: 1.4142 at 45°, 0.8660 at 30°, and 0 for the three-phase set; |sum| = 1.414, 0, 2 at Δz = 38.5, 77, 154 µm; envelope minima every 2.000 µs; unwrap first wrong at sample 32 (true step 3.150 rad); the FFT derivative matches the analytic one to \(1.8\times10^{-7}\ \mathrm{s^{-1}}\) out of \(2.9\times10^{7}\ \mathrm{s^{-1}}\) (\(\approx 6\times10^{-15}\) relative); 6.02 dB per octave; spectral peak 5.00 → 5.46 MHz on a zero-padded FFT grid, matching the predicted 5.46 MHz.
01.15Self-check
Q1Convert −2 + 2j to polar form.
\(r = \sqrt{4 + 4} = 2.828\); θ = atan2(2, −2) = 135° (second quadrant; atan(2/−2) would give −45°). So \(-2 + 2j = 2.828\,e^{j3\pi/4}\).
Q2Compute (1 + j)⁸ without expanding.
\(1 + j = \sqrt2\,\angle 45^\circ\). Eighth power: \((\sqrt2)^8 = 16\) and 8 × 45° = 360° ≡ 0°. Result: 16.
Q3Show that jj is a real number and find its principal value.
\(j = e^{j\pi/2}\), so \(j^j = e^{j\cdot j\pi/2} = e^{-\pi/2} = 0.2079\). (Other branches, \(j = e^{j(\pi/2 + 2\pi k)}\), give \(e^{-\pi/2 - 2\pi k}\), all real.)
Q4A 5 MHz phasor starts at angle 0. Where is it at t = 30 ns, and what changes if its frequency is −5 MHz?
Angle = 360° × 5 × 10⁶ × 30 × 10⁻⁹ = 54°: 0.588 + 0.809j. At −5 MHz: 0.588 − 0.809j. Same real part (cos even), opposite imaginary part (sin odd).
Q5Four samples taken at fs = 4f₀ are 0.5, 0.5, −0.5, −0.5. Find I, Q, amplitude and phase.
Multiply by 1, −j, −1, j: 0.5, −0.5j, 0.5, −0.5j. Sum 1 − j; mean 0.25 − 0.25j; ×2 → I + jQ = 0.5 − 0.5j: amplitude 0.707, phase −45°. Check: 0.707 cos(ωt − 45°) at 0°, 90°, 180°, 270° gives 0.5, 0.5, −0.5, −0.5 ✓.
Q6The phase of an echo at 5 MHz goes from 170° in one frame to −175° in the next. What is the true step and the axial displacement?
Raw difference −345°, beyond ±180°, so add 360°: +15° (assuming |step| < 180°). Phase = \(-2\pi f_0\tau\), so τ decreased by (15/360) × 200 ns = 8.33 ns, and the scatterer moved \(\Delta z = c\Delta\tau/2\) = 6.4 µm toward the probe.
Q7Simplify 2cos(ωt) + 2cos(ωt + 120°) + 2cos(ωt + 240°), and 3cos(ωt) + 4sin(ωt).
First: phasors 2∠0°, 2∠120°, 2∠240° form a closed triangle, sum 0. Second: sin ωt = cos(ωt − 90°), so the phasor is 3 − 4j = 5∠−53.13°: 5cos(ωt − 53.13°).
01.16Figures from the literature
For a fundamentals module the best published figures are textbook ones: the same arrows drawn by the people who taught them best, and one ultrasound paper where the phase rotator of (01.19) sits inside a real beamformer. Figure and section numbers were checked against the online editions; images are not reproduced, follow the links.
Feynman builds eiθ by repeated multiplication
R. P. Feynman, R. B. Leighton and M. Sands, The Feynman Lectures on Physics, Vol. I, Ch. 22 "Algebra", Sec. 22-5 and 22-6, Table 22-3, Table 22-4 and Fig. 22-1 (free online edition, Caltech). feynmanlectures.caltech.edu/I_22.html
Table 22-3 starts from \(10^{i/1024} \approx 1 + 0.0022486i\) and squares it repeatedly; Table 22-4 lists the powers of \(10^{i/8}\), and Fig. 22-1 plots them as dots on a circle. This is the limit route of Step 11 done by hand with real numbers: a tiny rotate-and-stretch factor, multiplied many times, walks around the unit circle. Eq. (22.9) then states \(e^{i\theta} = \cos\theta + i\sin\theta\), and Fig. 22-2 draws \(x + iy = re^{i\theta}\), the polar form of (01.2).
A complex sinusoid and its two projections
J. O. Smith III, Mathematics of the Discrete Fourier Transform (DFT), with Audio Applications, 2nd ed., W3K Publishing, 2007, "Projection of Circular Motion", Fig. 4.9 (free online edition, CCRMA Stanford). ccrma.stanford.edu/~jos/mdft/Projection_Circular_Motion.html
Fig. 4.9 plots \(e^{j\omega t}\) as a helix in (Re, Im, t) space with its three shadows: a circle in the complex plane, a cosine on the real-versus-time wall and a sine on the imaginary-versus-time wall. It is Step 14 and Demo 01.3 in one 3-D picture. The neighbouring page "Positive and Negative Frequencies" gives (01.14) and the symmetry \(\lvert X(-\omega)\rvert = \lvert X(\omega)\rvert\) for real signals.
Building the analytic signal in the frequency domain
J. O. Smith III, Mathematics of the DFT, 2nd ed., 2007, "Analytic Signals and Hilbert Transform Filters", Fig. 4.16. ccrma.stanford.edu/~jos/mdft/Analytic_Signals_Hilbert_Transform.html
Four panels: the spectrum of \(\cos\omega_0 t\) (lines at \(\pm\omega_0\)), of the quadrature copy \(\sin\omega_0 t\) (imaginary lines of opposite sign), of \(j\sin\omega_0 t\), and of their sum, where the negative-frequency lines cancel and the positive ones add to a single line at \(+\omega_0\). That is (01.27), \(X_a = 2u(f)X\), and the bottom panel of the Step 16 figure, derived line by line.
Quadrature signals with pictures: counter-rotating phasors and the I/Q mixer
R. G. Lyons, "A Quadrature Signals Tutorial: Complex, But Not Complicated," DSPRelated.com, 2013, Figs. 5–7 and 12–15. dsprelated.com/showarticle/192.php
Fig. 5 shows \(e^{j2\pi f_0 t}\) and \(e^{-j2\pi f_0 t}\) orbiting in opposite directions, Fig. 6 the corkscrew of a single phasor along the time axis with its cosine and sine projections, and Fig. 7 a cosine built from two half-magnitude counter-rotating phasors (Steps 15 and 16). Fig. 12 shows quadrature mixing shifting a spectrum by \(f_0\), and Figs. 14–15 the I/Q down-conversion block diagram with the spectra in each path: compare with the demodulator drawing in 01.6. The same material appears in Lyons, Understanding Digital Signal Processing, 3rd ed., Ch. 8.
The phase rotator inside an ultrasound IQ beamformer
V. Perrot, M. Polichetti, F. Varray and D. Garcia, "So you think you can DAS? A viewpoint on delay-and-sum beamforming," Ultrasonics 111:106309, 2021, Eq. (10) and Figs. 6, 10 and 14. doi:10.1016/j.ultras.2020.106309 · preprint arXiv:2007.11960
Eq. (10) beamforms I/Q signals as \(\sum_i \mathrm{IQ}_i(\tau_i)\,e^{2i\pi f_c\tau_i}\): interpolate each channel at its delay, then rotate by the carrier phase removed by the delay, exactly the I11 answer (their sign follows their demodulation convention). Fig. 14 shows the sparse DAS matrix that holds those interpolation weights and phase rotators. Fig. 6 shows how the assumed speed of sound changes the DAS output and introduces a phase-based quality metric; Fig. 10 applies that metric to choose the speed of sound in eight carotids in vivo: the IQ phase used as a sound-speed sensor, the step before aberration correction.
01.17Go deeper
- T. Needham, Visual Complex Analysis, Oxford University Press, 1997. Complex numbers as rotations and stretches throughout; chapter 1 is this module done properly.
- P. J. Nahin, An Imaginary Tale: The Story of √−1, Princeton University Press, 1998. History and the geometric meaning of i, including Euler's formula.
- R. P. Feynman, R. B. Leighton, M. Sands, The Feynman Lectures on Physics, Vol. I, Ch. 22 "Algebra", Addison-Wesley, 1963. Builds eiθ numerically from repeated multiplication, close to the limit route of Step 11.
- A. V. Oppenheim, A. S. Willsky, S. H. Nawab, Signals and Systems, 2nd ed., Prentice Hall, 1997. Complex exponentials as eigenfunctions of LTI systems (Module 03).
- R. G. Lyons, Understanding Digital Signal Processing, 3rd ed., Prentice Hall, 2010. Quadrature signals and complex down-conversion explained with pictures.
- J. O. Smith III, Mathematics of the Discrete Fourier Transform (DFT), W3K Publishing, 2007. Complex numbers, Euler's identity and sinusoids as a lead-in to the DFT (Module 02).
- T. L. Szabo, Diagnostic Ultrasound Imaging: Inside Out, 2nd ed., Academic Press, 2014. Analytic signals, IQ and Doppler in the ultrasound context.