matlab-estimate-sar

# SAR Estimation Skill

You are an expert RF and antenna engineer assisting with Specific Absorption Rate (SAR) estimation. Use MATLAB Antenna Toolbox to model antennas near biological tissue, compute internal E-fields, and calculate SAR for regulatory compliance assessment.

## When to Use

- User wants to compute SAR from an antenna near or inside biological tissue
- User wants to assess RF exposure compliance (FCC/ICNIRP limits)
- User wants to model a birdcage coil for MRI SAR analysis
- User asks about tissue absorption, power deposition, or RF safety
- User wants to compute point SAR or mass-averaged SAR (1g/10g)
- User is designing an implantable antenna and needs SAR estimation

## When NOT to Use

- User wants to design the antenna itself — use `matlab-design-antenna`
- User wants RF propagation/coverage analysis — use `matlab-analyze-rf-propagation`
- User wants RCS or scattering — use `matlab-analyze-rcs`

## SAR Formula

$$\mathrm{SAR} = \frac{\sigma |E|^2}{2\rho} \quad [\mathrm{W/kg}]$$

where:
- $\sigma$ = tissue conductivity (S/m)
- $|E|$ = electric field magnitude inside tissue (V/m)
- $\rho$ = tissue mass density (kg/m^3)

The conductivity relates to loss tangent via: $\sigma = \omega \varepsilon_0 \varepsilon_r \tan\delta$

## Three Approaches

| Approach | Antenna | Tissue Model | Tissue Properties | Best For |
|----------|---------|--------------|-------------------|----------|
| `birdcage` + `Phantom` | Birdcage MRI coil only | Volumetric tetrahedral mesh (struct) | Arbitrary εr and LossTangent (no cap) | MRI SAR with realistic tissue |
| `conformalArray` + `shape.Custom3D` | Any antenna (dipole, PIFA, etc.) | Surface triangulation with catalog material | Limited to catalog (LossTangent ≤ 0.03) | Phone/device SAR (antenna outside tissue) |
| Direct `EHfields` + SAR formula | Any antenna (pcbStack, catalog) | Post-processing only (no tissue in solver) | Real values in formula (no cap) | Implantable antenna SAR |

**Key limitation:** The `dielectric` class enforces `LossTangent <= 0.03`. Real tissue (brain at 900 MHz: tanδ ≈ 0.36) exceeds this cap. Only `birdcage.Phantom` bypasses this limit using a custom struct format. The direct `EHfields` approach avoids this cap by applying tissue conductivity in the SAR formula rather than in the solver.

## Approach 1: birdcage + Phantom (Full Tissue Properties)

The `birdcage` antenna is the ONLY object that supports volumetric dielectric bodies via its `Phantom` property. The phantom is included directly in the MoM solver.

### Phantom Data Format

```matlab
phantom = struct( ...
    Points=vertices, ...      % N x 3 vertex coordinates (meters)
    Tetrahedra=elements, ...  % M x 4 tetrahedral element connectivity
    EpsilonR=epsR, ...        % Relative permittivity (scalar)
    LossTangent=tanD);        % Dielectric loss tangent (scalar, no cap)
```

### Shipped Phantoms

- `humanheadcoarse.mat` -- 584 vertices, 2818 tetrahedra (fast)
- `humanheadfine.mat` -- finer resolution (more accurate, slower)

Both provide variables `P` (vertices) and `T` (tetrahedra). Apply `scaleFactor = 0.003` to get physical dimensions.

### Workflow

```matlab
freq = 128e6;  % 3T MRI Larmor frequency

% Tissue properties (gray matter at 128 MHz)
tissueEpsR = 77;
tissueSigma = 0.51;  % S/m
tissueRho = 1040;    % kg/m^3
omega = 2 * pi * freq;
eps0 = 8.854e-12;
tanD = tissueSigma / (omega * eps0 * tissueEpsR);

% Load phantom
load humanheadcoarse.mat
scaleFactor = 0.003;
phantom = struct( ...
    Points=scaleFactor * P, ...
    Tetrahedra=T, ...
    EpsilonR=tissueEpsR, ...
    LossTangent=tanD);

% Create birdcage with phantom
bc = birdcage(Phantom=phantom);
figure;
show(bc);

% Compute E-fields inside head
gridSpacing = 0.02;
[obsPoints, insideMask] = createObservationGrid(phantom.Points, gridSpacing);
[E, ~] = EHfields(bc, freq, obsPoints');

% Calculate SAR
E_mag_sq = abs(E(1,:)).^2 + abs(E(2,:)).^2 + abs(E(3,:)).^2;
SAR_point = tissueSigma * E_mag_sq / (2 * tissueRho);
```

## Approach 2: conformalArray + shape.Custom3D (Any Antenna)

Use `conformalArray` to place any antenna alongside a dielectric head phantom. The `shape.Custom3D` object acts as a passive dielectric scatterer in the MoM problem -- no feed is required.

### Key Insight

`shape` objects can be placed directly as elements in `conformalArray.Element` without wrapping in `customAntenna`. They do not need a feed and participate as passive dielectric bodies in the full-wave solver.

### Workflow

```matlab
freq = 2.4e9;
c = physconst("LightSpeed");
lambda = c / freq;

% Load and scale head phantom
load humanheadcoarse.mat
scaleFactor = 0.003;
pts = scaleFactor * P;
pts = pts * (0.18 / (max(pts(:,1)) - min(pts(:,1))));  % 180mm width

% Extract surface triangulation
TR = triangulation(T, pts);
[surfFaces, surfVertices] = freeBoundary(TR);
headTri = triangulation(surfFaces, surfVertices);

% Create dielectric shape
headShape = shape.Custom3D(headTri);
headShape.Dielectric = "TMM10";  % highest-permittivity catalog material

% Design antenna
d = design(dipole, freq);

% Build conformalArray
arr = conformalArray;
arr.Element = {d, headShape};
arr.ElementPosition = [0.10 0 0; 0 0 0];
arr.Reference = "origin";

figure;
show(arr);
figure;
pattern(arr, freq);

% Compute E-fields inside head
gridSpacing = 0.015;
[obsPoints, insideMask] = createObservationGrid(pts, gridSpacing);
[E, ~] = EHfields(arr, freq, obsPoints');

% SAR computation (using TMM10 properties)
epsR = 9.8; tanD_mat = 0.0022; rho = 2270;
omega = 2 * pi * freq; eps0 = 8.854e-12;
sigma = omega * eps0 * epsR * tanD_mat;
E_mag_sq = abs(E(1,:)).^2 + abs(E(2,:)).^2 + abs(E(3,:)).^2;
SAR_point = sigma * E_mag_sq / (2 * rho);
```

## Observation Grid Generation

Filter a 3D grid to only include points inside the head volume:

```matlab
function [obsPoints, insideMask] = createObservationGrid(pts, gridSpacing)
    margin = gridSpacing;
    xVec = (min(pts(:,1))+margin) : gridSpacing : (max(pts(:,1))-margin);
    yVec = (min(pts(:,2)