All Things RF Pulse: From Single Wire to Birdcage Coil Simulation¶
The Magnetic Field Around a Wire: Biot–Savart’s Discovery¶
In 1820, a simple classroom demonstration changed science: Hans Christian Ørsted noticed a compass needle jump when he switched on an electric current in a wire. This surprise revelation – that electricity can create magnetism – set the stage for a new era in physics. Shortly after, French scientists Jean-Baptiste Biot and Félix Savart took up the mystery. Through careful experiments, they quantified the magnetic field around a current-carrying wire, formulating what we now call the Biot–Savart law. In plain terms, their law says that each tiny segment of current produces a magnetic field that circles around the wire. The strength of this field falls off with distance, and its direction wraps around the wire according to the right-hand rule (point your thumb along the current, and your fingers curl in the field’s direction).
Mathematically, the Biot–Savart law is expressed as an integral:
$$ \mathbf{B}(\mathbf{r}) = \frac{\mu_0}{4\pi} \int_C \frac{I\,d\boldsymbol{\ell} \times \hat{\mathbf{r}}'}{|\mathbf{r}'|^2}\,, $$
where $I,d\boldsymbol{\ell}$ is a current element on the wire path $C$, and $\hat{\mathbf{r}}'$ is a unit vector pointing from that current element to the point $\mathbf{r}$ where we calculate the field. This looks complicated, but for a very long straight wire it simplifies nicely. If a steady current $I$ runs through an infinitely long wire, the magnetic field at a distance $r$ from the wire is:
$$ B(r) = \frac{\mu_0\,I}{2\pi\,r}\,, $$
pointing in circles around the wire (tangent to a circle of radius $r$ around the wire). The $\mu_0$ here is the permeability of free space – a constant that sets the scale of magnetic effects. This formula tells us a crucial fact: electric currents create magnetic fields, and those fields form closed loops around the current.
And if the wire has a finite length then the magnetic field at a point in space is given by
Analytic finite-wire formula
$$ \mathbf B(\mathbf r)=\frac{\mu_0 I}{4\pi R_\perp}\, (\sin\theta_2+\sin\theta_1)\,\hat{\boldsymbol\phi}, $$
where $\theta_1,\theta_2$ are the end-angles subtended at the
for full derivation of this formula, see this video.
Given this background, we can now explore how to compute the magnetic field around a wire in practice, especially in the context of MRI coil simulations.
- Analytic finite-wire formula – exact, fast.
- Numerical line-integral – slower but works for any bent or segmented wire if you refine the sampling.
Each snippet quotes the corresponding Biot–Savart equation inside the comments so you can see how code ↔ math line up.
1 Analytic field of a finite straight wire¶
For a straight segment from p₁ to p₂ (end-points in 3-D) carrying steady complex current I:
$$ \boxed{\; \mathbf B(\mathbf r)=\frac{\mu_0 I}{4\pi R}\; \bigl(\sin\theta_2+\sin\theta_1\bigr)\, \hat{\boldsymbol\phi}\;} $$
Definitions
- $\mathbf R = \mathbf r - \mathbf r_\text{mid}$ is the vector from the segment centre to the field point; $R = |\mathbf R|$.
- $\theta_1,\theta_2$ are the angles subtended by the extensions of p₁-r and p₂-r to the observation point.
- $\hat{\boldsymbol\phi} = \dfrac{\hat z \times \mathbf R}{R}$ is the azimuthal unit-vector (right-hand rule).
Comments
- Matches the text-book “finite straight conductor” formula.
- Works everywhere except exactly on the wire axis (
R_perp → 0).
import numpy as np
import pyvista as pv
pv.set_jupyter_backend("html") # <- no trame server
from scipy.constants import mu_0
def biot_savart_segment_exact(p1, p2, r_obs, I):
"""
Exact B-field (phasor) of a finite straight wire.
p1, p2 : ndarray (3,), endpoints [m]
r_obs : ndarray (N,3), observation points [m]
I : complex scalar current [A]
returns : ndarray (N,3), complex B-field [T]
"""
# Vector from p1 to p2 (wire direction)
dl = p2 - p1
L = np.linalg.norm(dl) # segment length
dl_hat = dl / L # unit vector along wire
z_hat = dl_hat # choose local z along the wire
# Vectors from ends to observation points
r1 = r_obs - p1 # (N,3)
r2 = r_obs - p2
# Distances
R1 = np.linalg.norm(r1, axis=1)
R2 = np.linalg.norm(r2, axis=1)
# Perpendicular (radial) component magnitude
# R_perp = |r × dl_hat| (distance from wire axis)
R_perp = np.linalg.norm(np.cross(r1, dl_hat), axis=1)
# Angles θ1, θ2 between wire and lines to obs-point
# sinθ = R_perp / R
sin_th1 = R_perp / (R1 + 1e-12)
sin_th2 = R_perp / (R2 + 1e-12)
# Azimuthal unit-vector phî = (dl_hat × r̂)/|dl_hat × r̂| the right hand rule
phi_hat = np.cross(dl_hat, r1) # (N,3) still not unit
phi_hat /= (np.linalg.norm(phi_hat, axis=1)[:, None] + 1e-12)
prefac = mu_0 * I / (4*np.pi*R_perp + 1e-12) # μ0 I / (4π R⊥)
coeff = (sin_th2 + sin_th1) # (N,)
return (prefac * coeff)[:, None] * phi_hat # (N,1)*(N,3)
2 Numerical Biot–Savart line integral¶
When your wire is bent, segmented, or you just want brute-force fidelity, sample it into small pieces and sum:
$$ \mathbf B(\mathbf r)=\frac{\mu_0}{4\pi}\sum_{j} I\,\frac{\Delta\boldsymbol\ell_j \times (\mathbf r-\mathbf r_j)} {|\,\mathbf r-\mathbf r_j\,|^{3}} . $$
Comments
- Converges to the analytic result as n_steps → ∞.
- Automatically handles wires of arbitrary shape if you break them into straight segments and call this repeatedly.
def biot_savart_segment_numeric(p1, p2, r_obs, I, n_steps=200):
"""
Numerical Biot–Savart integral along a straight segment.
Subdivides the wire into n_steps tiny elements Δℓ.
"""
# Parametric points along the wire
s = np.linspace(0, 1, n_steps, endpoint=False) + 0.5/n_steps
dl = (p2 - p1) / n_steps # tiny segment vector (3,)
# Observation minus mid-point of each sub-segment
r_mid = p1 + np.outer(s, (p2 - p1)) # (n_steps,3)
r_vec = r_obs[:, None, :] - r_mid[None, :, :] # (N, n_steps, 3)
R = np.linalg.norm(r_vec, axis=2)[..., None] # (N, n_steps, 1)
cross = np.cross(dl, r_vec) # (N, n_steps, 3)
dB = mu_0 * I / (4*np.pi) * cross / (R**3 + 1e-18)
return dB.sum(axis=1) # sum over segments
When to use which?¶
Use the analytic version whenever the segment-length is not comparable to voxel distance (most MRI coil grids). Switch to the numeric version if you need curved wires, end effects, or simply want a ground-truth benchmark.
Additionally, if you want to compute the field around a small segment of wire, you could also use: The integral form of the Biot–Savart law is basically a summation of the contributions from each infinitesimal segment of the wire. For a small segment, we can approximate the field as:
# Biot-Savart function for a straight wire small segment
def biot_savart_segment(p1, p2, r_obs, I):
dl = p2 - p1
r0 = r_obs - (p1 + p2) / 2
norm_r0 = np.linalg.norm(r0, axis=1)[:, np.newaxis] + 1e-12
cross = np.cross(dl, r0)
return mu_0 / (4 * np.pi) * I * cross / (norm_r0**3)
# This script is designed to visualize the magnetic field of a coil using PyVista.
gamma = 42.58e6 # Hz/T for proton
omega = 3 * gamma # Larmor frequency ~128 MHz for 3T for our MAGNETOM Vida Fit scanner
# Compute field from all wires at all observational positions
def compute_field_at_t(wire_actors, obs_positions,L, ti, g_norm=None):
B_total = np.zeros((obs_positions.shape[0], 3), dtype=np.complex128)
for (_, phi, (x, y)) in wire_actors:
I = np.exp(1j * (omega * ti + phi))
p1 = np.array([x, y, -L/2])
p2 = np.array([x, y, L/2])
B = biot_savart_segment_exact(p1, p2, obs_positions, I)
B_total += B
if g_norm is not None:
B_total *= g_norm
return B_total
# Field evaluation grid at center plane
grid_span = 0.2
Ngrid = 50
gx = np.linspace(-grid_span, grid_span, Ngrid)
gy = np.linspace(-grid_span, grid_span, Ngrid)
gx, gy = np.meshgrid(gx, gy)
gz = np.zeros_like(gx)
positions = np.stack([gx.ravel(), gy.ravel(), gz.ravel()], axis=1)
# Parameters
Nr = 1 # Number of wires
L = 0.03 # Length of wires
wire_radius = 0.002 # Thickness of wires
# omega = 0
Nt = 600 # Number of animation frames
t = np.linspace(0, 1e-6, Nt) # Time from 0 to 1 µs
# Spin parameters
M0 = np.array([0, 0, 1]) # Initial magnetization
# Create wire geometry
angles = np.linspace(0, 2*np.pi, Nr, endpoint=False)
wire_positions = np.array([[0,0]]) # Center wire at origin
current_phases = [0] # Phase shift per wire
anim_name = "b1_1_wire.gif"
html = pv_field_animation(
wire_positions,
current_phases,
t,
positions,
anim_name=anim_name,
)
EmbeddableWidget(value='<iframe srcdoc="<!DOCTYPE html>\n<html>\n <head>\n <meta http-equiv="Content-…
Rotating B1 Field Visualization
↳ Source: b1_1_wire.gif.html
One Wire is Not Enough: From Loops to the Birdcage Coil¶
A single straight wire is neat, but early engineers realized it’s not very efficient at delivering magnetic fields to an entire sample (like a human body). One improvement is to bend the wire into a loop. A loop of current produces a more useful field: inside the loop (along its central axis), the magnetic field lines bunch up and point mostly in one direction (similar to a bar magnet). In fact, the on-axis field of a loop of radius $R$ carrying current $I$ is given by:
$$ B_{\text{on-axis}}(x) = \frac{\mu_0\,I\,R^2}{2(R^2 + x^2)^{3/2}}\,, $$
directed along the loop’s axis. At the center of the loop ($x=0$), this simplifies to $B = \frac{\mu\_0 I}{2R}$ along the axis. So a loop concentrates the field in a specific region – a handy trait for an RF coil in MRI, which needs to bathe the tissue in a strong magnetic RF field. Early NMR experiments in the 1940s (by pioneers like Felix Bloch and Edward Purcell) often used loop-like coils to transmit and receive signals. They knew that sending a radiofrequency pulse through such a coil could tip the tiny “magnetization” of nuclei and detect it – the birth of NMR and eventually MRI.
Magnetic field created through a loop is shown in the figure below. This is a classic textbook example of a magnetic field around a current-carrying loop. The field lines are denser inside the loop, indicating a stronger magnetic field there.
Magnetic field around a current-carrying loop
</img>Spin Excitation and Transverse¶
To stimulate the NMR spin system, an RF-coil must produce a time-varying excitation field B1(t) with the following characteristics:
✔ B1(t) must have components that rotate near the resonant frequency (ωo), and
✔ B1(t) must have components perpendicular to the static magnetic field (Bo)
The simplest form of such an RF-transmit coil is a single loop oriented at right angles to the main magnetic field (see figure right). By driving a sinusoidal alternating current through this loop at the Larmor frequency, an oscillating magnetic field perpendicular to Bo is produced. Somewhat more sophisticated variations of this coil can be easily imagined, such as 2-loop (Helmholz) or multi-loop (solenoid) configurations.
Magnetic field around a current-carrying loop
Linearly Polarized Fields and Counter-Rotating Decomposition¶
What kind of field excites NMR?¶
To excite nuclear spins efficiently, the transmit RF field B₁(t) must satisfy two physical constraints:
$$ \begin{aligned} &\textbf{(i)}\quad B_1(t)\;\perp\;B_0\quad &\text{(must be transverse)}\\ &\textbf{(ii)}\quad B_1(t)\;\text{rotates at } \omega_0 \quad &\text{(must match Larmor frequency)} \end{aligned} $$
Consider the nuclear magnetisation vector $\vec{M}$ in a rotating frame. The Bloch equation tells us that only the component of the RF field that co-rotates with $\vec{M}$ at the Larmor frequency contributes to tipping it:
$$ \frac{d\vec{M}}{dt} = \gamma\,\vec{M} \times \big[B_0\,\hat{z} + 2\,\Re\{B_1^{+}(t)\,\hat{x}_\text{rot}\} \big] $$
Simulating the LP field and decomposition using two wires π apart¶
Let’s simulate a minimal LP configuration: two long wires placed π radians apart on a circle, both carrying in-phase current. This is like a simplified Helmholtz pair, producing a transverse oscillating B₁ field along a fixed axis.
This setup approximates a linear-polarized coil with a B₁ field along x.
leak = np.linalg.norm(B1_minus) / np.linalg.norm(B1_plus)
print(f"Fractional B1- leak (LP drive): {leak:.2%}")
You should see:
Fractional B1- leak (LP drive): 100.00%
This confirms the equal-magnitude counter-rotating decomposition of a linearly oscillating field. It's a real physical effect — the field itself is oscillating, but only half of it can rotate with the spins and drive NMR transitions.
When “just two wires” won’t cut it¶
Two wires driven in-phase → linearly polarized field → 50% wasted energy
Decomposition shows equal $B_1^+$ and $B_1^-$ components
Efficient MRI excitation demands a rotating field — achieved by
- adding a second coil orthogonal and 90° shifted, or
- building a proper birdcage or quadrature array
# Parameters
N_rungs = 2
R = 0.1
Nt = 600
t = np.linspace(0, 0.1e-6, Nt) # Time from 0 to 200 µs
# Geometry: positions of rungs
angles = np.linspace(0, 2*np.pi, N_rungs, endpoint=False)
x_rungs = R * np.cos(angles) # x-coordinates of rungs
y_rungs = R * np.sin(angles) # y-coordinates of rungs: shape (N_rungs,)
# Phase shifts for each rung
current_phases = angles # mode-1 in-phase current i.e no delay in current
wire_positions = np.stack([x_rungs, y_rungs], axis=1) # Center wire at origin
anim_name = "b1_2_rungs_in_phase.gif"
pv_field_animation(
wire_positions,
current_phases,
t*10, # Cscale the time to see animation faster
positions,
anim_name=anim_name,
)
EmbeddableWidget(value='<iframe srcdoc="<!DOCTYPE html>\n<html>\n <head>\n <meta http-equiv="Content-…
# Usage:
pretty_iframe(f"{anim_name}.html", title="2 Wires Field Visualization")
2 Wires Field Visualization
↳ Source: b1_2_rungs_in_phase.gif.html
pretty_gif(anim_name)
Scaled Magnetic Field Simulation¶
anim_name = "b1_2_rungs_in_phase_scaled.gif"
pv_field_animation(
wire_positions,
current_phases,
t*10, # Cscale the time to see animation faster
positions,
scaled=True,
anim_name=anim_name,
)
EmbeddableWidget(value='<iframe srcdoc="<!DOCTYPE html>\n<html>\n <head>\n <meta http-equiv="Content-…
pretty_gif(anim_name)
A linearly polarized field is inherently inefficient¶
Suppose we drive a single loop (or equivalently two opposite wires) with a sinusoidal current:
$$ I(t) = I_0 \cos(\omega_0 t) $$
The resulting magnetic field at the center is also linearly polarized:
$$ \vec{B}_1(t) = B_0 \cos(\omega_0 t)\, \hat{x} $$
But this can be algebraically decomposed into two counter-rotating circular fields:
$$ \boxed{ \vec{B}_1(t) = \frac{B_0}{2}\left[\hat{x} + i\hat{y}\right] e^{-i\omega_0 t} + \frac{B_0}{2}\left[\hat{x} - i\hat{y}\right] e^{+i\omega_0 t} = B_1^+(t) + B_1^-(t) } $$
Here:
- $B_1^+$ rotates with the spin precession (resonant)
- $B_1^-$ rotates against it (off-resonant → ignored)
Because only $B_1^+$ contributes to excitation, 50% of the transmitted power is wasted in a linear coil. Yet, this is the default field generated by many simple configurations.
B_total = np.zeros((2500, 3), dtype=np.complex128)
for φ, (x, y) in zip(current_phases, wire_positions):
I = np.exp(1j * φ) # current phasor
p1 = np.array([x, y, -L/2])
p2 = np.array([x, y, +L/2])
B_total += biot_savart_segment_exact(p1, p2, positions, I)
# Compute B1+ and B1- from complex transverse field
Bx, By = B_total[:, 0], B_total[:, 1]
B1_plus = 0.5 * (Bx + 1j * By)
B1_minus = 0.5 * (Bx - 1j * By)
# Now compute |B1+| map
B1_magnitude = np.abs(B1_plus) # shape = (2500,)
leak = np.linalg.norm(B1_minus) / np.linalg.norm(B1_plus)
print(f"Fractional B1- leak (LP drive): {leak:.2%}")
Fractional B1- leak (LP drive): 100.00%
# Assume you have (x, y) coordinates corresponding to each point
# obs_positions is (2500, 3)
x = positions[:, 0]
y = positions[:, 1]
# Create a 2D scatter or interpolate to a grid
import matplotlib.pyplot as plt
plt.figure(figsize=(5, 5))
plt.scatter(x, y, c=B1_magnitude, cmap='viridis', s=10)
plt.colorbar(label='|B1+| (T)')
plt.title('B1+ Field Map')
plt.xlabel('x (m)')
plt.ylabel('y (m)')
plt.axis('equal')
plt.show()
B1_grid = B1_magnitude.reshape(50, 50)
plt.imshow(B1_grid, extent=[x.min(), x.max(), y.min(), y.max()], origin='lower', cmap='viridis')
plt.colorbar(label='|B1+|')
plt.title('B1+ Field Map')
plt.xlabel('x (m)')
plt.ylabel('y (m)')
Text(0, 0.5, 'y (m)')
This confirms the equal-magnitude counter-rotating decomposition of a linearly oscillating field. It's a real physical effect — the field itself is oscillating, but only half of it can rotate with the spins and drive NMR transitions.
Transverse Magnetic Field: The Birdcage Coil¶
As we saw, a single loop or linearly polarized coil still has limitations. The field it produces might be strong at the center, but it falls off near the edges. Also, a simple loop’s field can be non-uniform across a large object. And it wastes 50% of the energy as well. To image larger volumes (like a human torso) and to get a uniform excitation, engineers needed a more sophisticated design. This is where the birdcage coil comes in – a design that looks like its namesake and revolutionized MRI in the 1980s.
In 1985, Cecil E. Hayes and colleagues at GE introduced the birdcage resonator. Picture a birdcage: it has two circular end-rings connected by a number of straight bars (the “rungs”) going down the sides radioeng.cz. The MRI birdcage coil is exactly that: two conductive loops (end-rings) and an even number of straight wires (rungs) connecting them. Typically there might be 8, 12, or 16 rungs spaced evenly around a cylinder. Each rung carries current, and capacitors (strategically placed either in the rungs or end-rings) help tune the coil to the desired frequency. The result is a resonant structure that can create a very uniform magnetic field in its center. The name “birdcage” really comes from its appearance.
Birdcage coil design
</img>Visualizing the Current in the Rungs¶
Let's consider a birdcage coil with N rungs, each carrying a steady current I. First, let's just focus on the Current in the rungs, ignoring the resultant magnetic field for now. If we assume all rungs carry the same current I in the same direction.
# Parameters
N_rungs = 8
R = 0.1
Nt = 600
t = np.linspace(0, 0.1e-6, Nt) # Time from 0 to 200 µs
Angular Position of Each Rung¶
we are going to place the N_rungs rungs evenly around a circle. The angular position of each rung can be calculated as:
$$\theta_j = \frac{2\pi j}{N_{\text{rungs}}}\,,
$$
where $j$ is the rung index (from 0 to $N_{\text{rungs}}-1$).
Starting the phase $\phi_j$ of the AC current in each rung is set to zero, meaning all rungs are in phase. The current in each rung can be expressed as: $$ I_j = I \cdot e^{i\omega t + \phi_j} = I \cdot e^{i\omega t + 0} = I \cdot e^{i\omega t}, $$ where $I$ is the magnitude of the current in each rung.
# Geometry: positions of rungs
angles = np.linspace(0, 2*np.pi, N_rungs, endpoint=False)
x_rungs = R * np.cos(angles) # x-coordinates of rungs
y_rungs = R * np.sin(angles) # y-coordinates of rungs: shape (N_rungs,)
# Phase shifts for each rung
phases = np.zeros_like(angles) # mode in-phase current i.e no delay in current
from matplotlib.animation import FuncAnimation
import matplotlib.pyplot as plt
from IPython.display import HTML
# Setup figure
fig = plt.figure(figsize=(12,6))
ax_current = fig.add_subplot(1, 1,1)
ax_current.set_xlim(-1.5*R, 1.5*R)
ax_current.set_ylim(-1.5*R, 1.5*R)
ax_current.set_aspect('equal')
lines = [ax_current.plot([], [], lw=3)[0] for _ in range(N_rungs)]
ax_current.set_title("Rotating Current Pattern (Quadrature Drive)")
ax_current.set_xticks([])
ax_current.set_yticks([])
# Animation update function
def update(i):
for j, line in enumerate(lines):
I = 1*np.cos(omega * t[i] + phases[j]) # Current at time t[i] for rung j
# Represent current as color intensity and vector length
x = x_rungs[j]
y = y_rungs[j]
dx = 0.05 * I * np.cos(angles[j])
dy = 0.05 * I * np.sin(angles[j])
line.set_data([x - dx, x + dx], [y - dy, y + dy])
line.set_color(plt.cm.plasma((I + 1)/2)) # map current to color
return lines
ani = FuncAnimation(fig, update, frames=Nt, interval=5,
blit=True,
repeat=True)
plt.close() # Close the figure to avoid displaying it in Jupyter
HTML(ani.to_jshtml())
# ani.save("anim1.gif", writer="ffmpeg", fps=30)
# pretty_gif("anim1.gif", caption="")