### nbdev Development Setup

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/index.ipynb

Set up a development environment using nbdev. This includes installing necessary tools like jupyterlab and nbdev, and configuring git hooks.

```zsh
mamba install jupyterlab nbdev -c fastai -c conda-forge 
nbdev_install_quarto  # To build docs
nbdev_install_hooks   # Make notebooks git-friendly
```

```zsh
nbdev_preview  # Render docs locally and inspect in browser
nbdev_clean    # NECESSARY BEFORE PUSHING
nbdev_test     # tests notebooks
nbdev_export   # builds package and builds docs
```

--------------------------------

### Load example CT data and set device

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/tutorials/trilinear.ipynb

Loads an example CT scan and sets the computation device to CUDA if available, otherwise CPU.

```python
subject = load_example_ct()
device = torch.device("cuda" if torch.cuda.is_available() else "cpu")
```

--------------------------------

### Install DiffDRR Development Tools

Source: https://github.com/eigenvivek/diffdrr/blob/main/README.md

Installs necessary tools for developing DiffDRR using nbdev. Ensure you have mamba installed.

```zsh
mamba install jupyterlab nbdev -c fastai -c conda-forge 
nbdev_install_quarto  # To build docs
nbdev_install_hooks  # Make notebooks git-friendly
```

--------------------------------

### Hello, World! DRR Generation in PyTorch

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/index.ipynb

A minimal example demonstrating DiffDRR's DRR module. It initializes the DRR with a CT volume, sets camera parameters, and generates a DRR image. Ensure CUDA is available for GPU acceleration.

```python
import matplotlib.pyplot as plt
import torch

from diffdrr.drr import DRR
from diffdrr.data import load_example_ct
from diffdrr.visualization import plot_drr

# Read in the volume and get its origin and spacing in world coordinates
subject = load_example_ct()

# Initialize the DRR module for generating synthetic X-rays
device = torch.device("cuda" if torch.cuda.is_available() else "cpu")
drr = DRR(
    subject,     # An object storing the CT volume, origin, and voxel spacing
    sdd=1020.0,  # Source-to-detector distance (i.e., focal length)
    height=200,  # Image height (if width is not provided, the generated DRR is square)
    delx=2.0,    # Pixel spacing (in mm).
).to(device)

# Set the camera pose with rotations (yaw, pitch, roll) and translations (x, y, z)
rotations = torch.tensor([[0.0, 0.0, 0.0]], device=device)
translations = torch.tensor([[0.0, 850.0, 0.0]], device=device)

# 📸 Also note that DiffDRR can take many representations of SO(3) 📸
# For example, quaternions, rotation matrix, axis-angle, etc...
img = drr(rotations, translations, parameterization="euler_angles", convention="ZXY")
plot_drr(img, ticks=False)
plt.show()
```

--------------------------------

### Install DiffDRR

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/index.ipynb

Install the latest stable release of DiffDRR using pip. For development, clone the repository and install in editable mode.

```zsh
pip install diffdrr
```

```zsh
git clone https://github.com/eigenvivek/DiffDRR.git --depth 1
pip install -e 'DiffDRR/[dev]'
```

--------------------------------

### Import necessary libraries

Source: https://github.com/eigenvivek/diffdrr/blob/main/experiments/animations.ipynb

Imports all required libraries for DRR visualization and computation. Ensure these are installed before running.

```python
%matplotlib inline
```

```python
from pathlib import Path

import matplotlib.pyplot as plt
import matplotlib.gridspec as gridspec
import numpy as np
import pandas as pd
import torch
from IPython.display import HTML

from diffdrr import DRR, load_example_ct
from diffdrr.visualization import animate, plot_drr
```

--------------------------------

### load_example_ct

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/api/03_data.ipynb

Loads a predefined example chest CT scan and its corresponding labelmap for demonstration purposes. This function simplifies the process of getting started with sample data.

```APIDOC
## load_example_ct

### Description
Loads an example chest CT for demonstration purposes.

### Parameters
- **labels** (any) - Optional - A subset of structures from the labelmap that you want to render.
- **orientation** (str) - Optional - A frame-of-reference change for the C-arm. Defaults to "AP".
- **bone_attenuation_multiplier** (float) - Optional - A constant multiplier to the estimated density of bone voxels. Defaults to 1.0.
- **kwargs** (any) - Optional - Any additional kwargs can be passed to the `torchio.Subject` and accessed as a dictionary.

### Returns
- **Subject** - A torchio.Subject object containing the example CT scan and labelmap.
```

--------------------------------

### Prepare CT Data and Pose

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/tutorials/timing.ipynb

Loads example CT data and sets up the computation device (GPU if available, otherwise CPU). It then defines the rotation and translation parameters for the X-ray source and converts them into a pose format compatible with the DRR computation.

```python
# Read in the volume
subject = load_example_ct()
device = "cuda" if torch.cuda.is_available() else "cpu"

# Get parameters for the detector
rotations = torch.tensor([[0.0, 0.0, 0.0]], device=device)
translations = torch.tensor([[0.0, 850.0, 0.0]], device=device)
pose = convert(rotations, translations, parameterization="euler_angles", convention="ZXY")
```

--------------------------------

### Hello, World! DRR Generation

Source: https://github.com/eigenvivek/diffdrr/blob/main/README.md

Minimal example to generate a DRR using DiffDRR. Initializes the DRR module with subject data, source-to-detector distance, and image dimensions. Sets camera pose and renders the DRR. Requires matplotlib for visualization.

```python
import matplotlib.pyplot as plt
import torch

from diffdrr.drr import DRR
from diffdrr.data import load_example_ct
from diffdrr.visualization import plot_drr

# Read in the volume and get its origin and spacing in world coordinates
subject = load_example_ct()

# Initialize the DRR module for generating synthetic X-rays
device = torch.device("cuda" if torch.cuda.is_available() else "cpu")
drr = DRR(
    subject,     # An object storing the CT volume, origin, and voxel spacing
    sdd=1020.0,  # Source-to-detector distance (i.e., focal length)
    height=200,  # Image height (if width is not provided, the generated DRR is square)
    delx=2.0,    # Pixel spacing (in mm)
).to(device)

# Set the camera pose with rotations (yaw, pitch, roll) and translations (x, y, z)
rotations = torch.tensor([[0.0, 0.0, 0.0]], device=device)
translations = torch.tensor([[0.0, 850.0, 0.0]], device=device)

# 📸 Also note that DiffDRR can take many representations of SO(3) 📸
# For example, quaternions, rotation matrix, axis-angle, etc...
img = drr(rotations, translations, parameterization="euler_angles", convention="ZXY")
plot_drr(img, ticks=False)
plt.show()
```

--------------------------------

### Example: Logarithmic Geodesic SE(3) Calculation

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/api/05_metrics.ipynb

Demonstrates how to use the LogGeodesicSE3 class to compute the distance between two poses. Ensure poses are converted to RigidTransform objects.

```python
# SE(3) distance
geodesic_se3 = LogGeodesicSE3()

pose_1 = convert(
    torch.tensor([[0.1, 1.0, torch.pi]]),
    torch.ones(1, 3),
    parameterization="euler_angles",
    convention="ZYX",
)
pose_2 = convert(
    torch.tensor([[0.1, 1.1, torch.pi]]),
    torch.zeros(1, 3),
    parameterization="euler_angles",
    convention="ZYX",
)

geodesic_se3(pose_1, pose_2)
```

--------------------------------

### Load Example CT Scan with DiffDRR

Source: https://context7.com/eigenvivek/diffdrr/llms.txt

Loads a pre-packaged chest CT scan for testing with DiffDRR. Supports custom orientation and bone attenuation scaling.

```python
from diffdrr.data import load_example_ct

# Load with default AP (anteroposterior) orientation
subject = load_example_ct()

# Load with PA orientation and boosted bone contrast
subject_pa = load_example_ct(orientation="PA", bone_attenuation_multiplier=2.0)

# Load with a specific label subset from the bundled labelmap
subject_bone = load_example_ct(labels=[1, 2, 3])  # e.g., specific structure indices

print(subject.volume.shape)     # torch.Size([1, D, H, W]) raw CT tensor
print(subject.density.shape)    # torch.Size([1, D, H, W]) normalized density
print(subject.reorient)         # 4x4 reorientation matrix (torch.Tensor)
```

--------------------------------

### Plot Nelder-Mead Optimization Results

Source: https://github.com/eigenvivek/diffdrr/blob/main/experiments/figure.ipynb

Plots the results of Nelder-Mead optimization. Requires scipy to be installed.

```python
plot("neldermead", scipy=True, max_itr=500, lossfn="loss", cutoff=-0.999)
```

--------------------------------

### Example: Double Geodesic SE(3) Calculation

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/api/05_metrics.ipynb

Illustrates the usage of DoubleGeodesicSE3 to obtain angular, translational, and combined geodesic distances. The source-to-detector distance (sdd) is a required parameter.

```python
# Angular distance and translational distance both in mm
double_geodesic = DoubleGeodesicSE3(1020.0)

pose_1 = convert(
    torch.tensor([[0.1, 1.0, torch.pi]]),
    torch.ones(1, 3),
    parameterization="euler_angles",
    convention="ZYX",
)
pose_2 = convert(
    torch.tensor([[0.1, 1.1, torch.pi]]),
    torch.zeros(1, 3),
    parameterization="euler_angles",
    convention="ZYX",
)

double_geodesic(pose_1, pose_2)  # Angular, translational, double geodesics
```

--------------------------------

### Run Registration Experiment

Source: https://github.com/eigenvivek/diffdrr/blob/main/experiments/figure.ipynb

Command to execute the registration experiment from the home directory. Specify the number of DRRs to generate.

```bash
experiments/registration/registration.py --n_drrs=1000
```

--------------------------------

### Optimization Methods Comparison

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/tutorials/registration.ipynb

Runs and compares various optimization methods including base SGD, SGD with momentum, SGD with momentum and dampening, Adam, L-BFGS, and L-BFGS with line search. Each method initializes a DRR and Registration object, performs optimization, and then cleans up the DRR object.

```python
#| code-fold: true
# Base SGD
drr = DRR(**kwargs).to(device)
reg = Registration(
    drr,
    rotations.clone(),
    translations.clone(),
    parameterization="euler_angles",
    convention="ZXY",
)
params_base = optimize(reg, ground_truth)
del drr

# SGD + momentum
drr = DRR(**kwargs).to(device)
reg = Registration(
    drr,
    rotations.clone(),
    translations.clone(),
    parameterization="euler_angles",
    convention="ZXY",
)
params_momentum = optimize(reg, ground_truth, momentum=5e-1)
del drr

# SGD + momentum + dampening
drr = DRR(**kwargs).to(device)
reg = Registration(
    drr,
    rotations.clone(),
    translations.clone(),
    parameterization="euler_angles",
    convention="ZXY",
)
params_momentum_dampen = optimize(reg, ground_truth, momentum=5e-1, dampening=1e-4)
del drr

# Adam
drr = DRR(**kwargs).to(device)
reg = Registration(
    drr,
    rotations.clone(),
    translations.clone(),
    parameterization="euler_angles",
    convention="ZXY",
)
params_adam = optimize(reg, ground_truth, 1e-1, 5e0, optimizer="adam")
del drr

# L-BFGS
drr = DRR(**kwargs).to(device)
reg = Registration(
    drr,
    rotations.clone(),
    translations.clone(),
    parameterization="euler_angles",
    convention="ZXY",
)
params_lbfgs = optimize_lbfgs(reg, ground_truth, lr=3e-1)
del drr

# L-BFGS + line search
drr = DRR(**kwargs).to(device)
reg = Registration(
    drr,
    rotations.clone(),
    translations.clone(),
    parameterization="euler_angles",
    convention="ZXY",
)
params_lbfgs_wolfe = optimize_lbfgs(
    reg, ground_truth, lr=1e0, line_search_fn="strong_wolfe"
)
del drr
```

--------------------------------

### Import necessary libraries for visualization

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/api/04_visualization.ipynb

Imports PyVista, VTK, and TorchIO for 3D visualization and data handling. Sets VTK logger verbosity.

```python
#| export
import pyvista
import vtk
from torchio import Subject

vtk.vtkLogger.SetStderrVerbosity(vtk.vtkLogger.ConvertToVerbosity(-1))
```

--------------------------------

### Generate Ground Truth DRR

Source: https://github.com/eigenvivek/diffdrr/blob/main/experiments/animations.ipynb

Sets up the DRR object and generates a ground truth DRR image from example CT data. This requires the `diffdrr` library and example data.

```python
# Make the ground truth DRR
volume, spacing = load_example_ct()
bx, by, bz = np.array(volume.shape) * np.array(spacing) / 2
true_params = {
    "sdr": 200.0,
    "theta": torch.pi,
    "phi": 0,
    "gamma": torch.pi / 2,
    "bx": bx,
    "by": by,
    "bz": bz,
}

drr = DRR(volume, spacing, height=100, delx=10.0, device="cuda")
ground_truth = drr(**true_params)
```

--------------------------------

### Initialize Registration Module for Pose Optimization

Source: https://context7.com/eigenvivek/diffdrr/llms.txt

Wraps a DRR module with learnable pose parameters for optimization. Minimizes image similarity loss to find the optimal camera pose. Requires a DRR module, initial pose parameters, and a similarity criterion.

```python
from diffdrr.registration import Registration
from diffdrr.metrics import NormalizedCrossCorrelation2d
import torch

# Initialize pose parameters (initial guess)
init_rot   = torch.tensor([[0.05, -0.03, 0.02]], device=device)
init_trans = torch.tensor([[0.0, 860.0, 5.0]], device=device)

registration = Registration(
    drr,
    rotation=init_rot,
    translation=init_trans,
    parameterization="euler_angles",
    convention="ZXY",
).to(device)

# Fixed reference image (ground truth X-ray or reference DRR)
rot_gt  = torch.tensor([[0.0, 0.0, 0.0]], device=device)
trans_gt = torch.tensor([[0.0, 850.0, 0.0]], device=device)
fixed = drr(rot_gt, trans_gt, parameterization="euler_angles", convention="ZXY").detach()

optimizer = torch.optim.Adam(registration.parameters(), lr=1e-2)
criterion = NormalizedCrossCorrelation2d()

for step in range(100):
    optimizer.zero_grad()
    moving = registration()            # Renders DRR at current pose estimate
    loss = -criterion(fixed, moving).mean()
    loss.backward()
    optimizer.step()
    if step % 20 == 0:
        print(f"Step {step:03d} | Loss: {loss.item():.4f}")

# Retrieve final estimated pose
final_pose = registration.pose
print(final_pose.matrix)
```

--------------------------------

### Initialize and Use MutualInformation

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/api/05_metrics.ipynb

Demonstrates initializing and using the MutualInformation metric for image similarity.

```python
mi = MutualInformation()
mi(x1, x2)
```

--------------------------------

### Quaternion Multiplication (Composition)

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/api/06_pose.ipynb

Multiplies two quaternions representing rotations to get the quaternion representing their composition. The result is standardized to have a non-negative real part.

```python
def quaternion_multiply(a: torch.Tensor, b: torch.Tensor) -> torch.Tensor:
    """
    Multiply two quaternions representing rotations, returning the quaternion
    representing their composition, i.e. the versor with nonnegative real part.
    Usual torch rules for broadcasting apply.

    Args:
        a: Quaternions as tensor of shape (..., 4), real part first.
        b: Quaternions as tensor of shape (..., 4), real part first.

    Returns:
        The product of a and b, a tensor of quaternions of shape (..., 4).
    """
    ab = quaternion_raw_multiply(a, b)
    return standardize_quaternion(ab)
```

--------------------------------

### Get Index from Axis Letter

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/api/06_pose.ipynb

Helper function to map axis letters ('X', 'Y', 'Z') to their corresponding integer indices (0, 1, 2).

```python
def _index_from_letter(letter: str) -> int:
    if letter == "X":
        return 0
    if letter == "Y":
        return 1
    if letter == "Z":
        return 2
    raise ValueError("letter must be either X, Y or Z.")
```

--------------------------------

### Initialize and Use NormalizedCrossCorrelation2d

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/api/05_metrics.ipynb

Demonstrates initializing and using the NormalizedCrossCorrelation2d metric with default and custom parameters like epsilon and patch size.

```python
x1 = torch.randn(8, 1, 128, 128)
x2 = torch.randn(8, 1, 128, 128)

cc = NormalizedCrossCorrelation2d()
cc(x1, x2)

cc = NormalizedCrossCorrelation2d(eps=1e-1)
cc(x1, x2)

cc = NormalizedCrossCorrelation2d(patch_size=9)
cc(x1, x2)
```

--------------------------------

### Initialize DRR and Generate Target DRR

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/tutorials/metrics.ipynb

Sets up the DRR module with subject data, source-to-detector distance, and image dimensions. It then generates a target DRR from a known pose.

```python
#| code-fold: true
import matplotlib.pyplot as plt
import torch
from mpl_toolkits.mplot3d import Axes3D
from tqdm import tqdm

from diffdrr.data import load_example_ct
from diffdrr.drr import DRR
from diffdrr.metrics import MultiscaleNormalizedCrossCorrelation2d
from diffdrr.visualization import plot_drr

# Read in the volume and get its origin and spacing in world coordinates
subject = load_example_ct()

# Initialize the DRR module for generating synthetic X-rays
device = torch.device("cuda" if torch.cuda.is_available() else "cpu")
drr = DRR(
    subject,  # A diffdrr.data.Subject object storing the CT volume, origin, and voxel spacing
    sdd=1020,  # Source-to-detector distance (i.e., the C-arm's focal length)
    height=200,  # Height of the DRR (if width is not seperately provided, the generated image is square)
    delx=2.0,  # Pixel spacing (in mm)
).to(device)

# Get parameters for the detector
alpha, beta, gamma = 0.0, 0.0, 0.0
bx, by, bz = 0.0, 850.0, 0.0
rotations = torch.tensor([[alpha, beta, gamma]], device=device)
translations = torch.tensor([[bx, by, bz]], device=device)

# MNake the DRR
target_drr = drr(
    rotations, translations, parameterization="euler_angles", convention="ZYX"
)
plot_drr(target_drr, ticks=False)
plt.show()
```

--------------------------------

### Get shape of batched DRR output

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/tutorials/introduction.ipynb

Prints the shape of the output tensor when rendering multiple DRRs, showing the batch size, channels, height, and width.

```python
img.shape
```

--------------------------------

### Visualize camera, detector, and principal ray in 3D

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/tutorials/geometry.ipynb

Creates a 3D visualization using PyVista, plotting the CT mesh, camera frustum, detector plane with DRR texture, and the principal ray. Exports the visualization to an HTML file.

```python
# Make a mesh from the camera and detector plane
camera, detector, texture, principal_ray = img_to_mesh(drr, pose)

# Make the plot
plotter = pyvista.Plotter()
plotter.add_mesh(ct)
plotter.add_mesh(camera, show_edges=True, line_width=1.5)
plotter.add_mesh(principal_ray, color="lime", line_width=3)
plotter.add_mesh(detector, texture=texture)

# Render the plot
plotter.add_axes()
plotter.add_bounding_box()
plotter.show_bounds(grid="front", location="outer", all_edges=True)

# plotter.show()  # If running Jupyter locally
# plotter.show(jupyter_backend="server")  # If running Jupyter remotely
plotter.export_html("render.html")
```

--------------------------------

### Initialize and Use GradientNormalizedCrossCorrelation2d

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/api/05_metrics.ipynb

Illustrates initializing and using the GradientNormalizedCrossCorrelation2d metric, with and without a specified patch size.

```python
gncc = GradientNormalizedCrossCorrelation2d()
gncc(x1, x2)

gncc = GradientNormalizedCrossCorrelation2d(patch_size=9)
gncc(x1, x2)
```

--------------------------------

### Convert: Unified Rotation Representation Conversion

Source: https://context7.com/eigenvivek/diffdrr/llms.txt

Utilizes the `convert` function for unified conversion of rotation and translation between various parameterizations. Examples include 6D rotation, axis-angle, and SE(3) log-map.

```python
from diffdrr.pose import convert
import torch

rot_6d = torch.randn(2, 6)   # 6D continuous rotation (Zhou et al.)
trans  = torch.randn(2, 3)
pose = convert(rot_6d, trans, parameterization="rotation_6d")
print(pose.matrix.shape)   # torch.Size([2, 4, 4])

# Axis-angle (3-vector, magnitude = angle in radians)
aa  = torch.tensor([[0.0, 0.0, 1.5708]])   # 90° about Z
t   = torch.tensor([[0.0, 850.0, 0.0]])
pose_aa = convert(aa, t, parameterization="axis_angle")

# SE(3) log-map (6-vector: [tx, ty, tz, rx, ry, rz])
log = torch.tensor([[0.0, 850.0, 0.0, 0.0, 0.0, 0.0]])
pose_se3 = convert(log[:, 3:], log[:, :3], parameterization="se3_log_map")
```

--------------------------------

### Generate Target X-ray with DiffDRR

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/tutorials/reconstruction.ipynb

Generates a target X-ray image from an example CT volume using DiffDRR. Requires matplotlib, torch, and DiffDRR libraries. The output is a normalized DRR image.

```python
#| code-fold: true
import matplotlib.pyplot as plt
import torch
from tqdm import tqdm

from diffdrr.data import load_example_ct
from diffdrr.drr import DRR
from diffdrr.visualization import plot_drr

device = torch.device("cuda" if torch.cuda.is_available() else "cpu")

subject = load_example_ct()
drr = DRR(subject, sdd=1020.0, height=200, delx=2.0).to(device=device)

rot = torch.tensor([[0.0, 0.0, 0.0]], device=device)
xyz = torch.tensor([[0.0, 850.0, 0.0]], device=device)
gt = drr(rot, xyz, parameterization="euler_angles", convention="ZXY")
gt = (gt - gt.min()) / (gt.max() - gt.min())
plot_drr(gt, ticks=False)
plt.show()
```

--------------------------------

### RigidTransform: Batched SE(3) Rigid Transform

Source: https://context7.com/eigenvivek/diffdrr/llms.txt

Demonstrates the creation and application of RigidTransform for batched SE(3) transformations. Includes point cloud transformation, composition, inversion, and conversion to/from Euler angles and quaternions.

```python
from diffdrr.pose import RigidTransform, convert
import torch

# Build from a 4x4 matrix
matrix = torch.eye(4).unsqueeze(0)   # batch of 1 identity transform
T = RigidTransform(matrix)

# Apply to a point cloud (b, n, 3)
pts = torch.randn(1, 100, 3)
pts_transformed = T(pts)
print(pts_transformed.shape)   # torch.Size([1, 100, 3])

# Compose two transforms: T_total = T2 ∘ T1
T2 = RigidTransform(torch.eye(4).unsqueeze(0))
T_composed = T.compose(T2)

# Closed-form inversion
T_inv = T.inverse()

# Convert to Euler angles and back
rot, trans = T.convert("euler_angles", convention="ZXY")
print(rot.shape)    # torch.Size([1, 3])
print(trans.shape)  # torch.Size([1, 3])

# Create from Euler angles
T_from_euler = convert(
    torch.tensor([[0.1, 0.2, 0.3]]),
    torch.tensor([[0.0, 850.0, 0.0]]),
    parameterization="euler_angles",
    convention="ZXY",
)

# Create from quaternion (real part first)
T_from_quat = convert(
    torch.tensor([[1.0, 0.0, 0.0, 0.0]]),
    torch.tensor([[0.0, 850.0, 0.0]]),
    parameterization="quaternion",
)

# SE(3) log map (6-vector)
log_vec = T.get_se3_log()
print(log_vec.shape)   # torch.Size([1, 6])
```

--------------------------------

### Initialize Detector Class

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/api/02_detector.ipynb

Constructs a 6 DoF X-ray detector system. Initializes intrinsic parameters, reorientation matrix, and default source/detector plane positions.

```python
class Detector(torch.nn.Module):
    """Construct a 6 DoF X-ray detector system. This model is based on a C-Arm."""

    def __init__(
        self,
        sdd: float,  # Source-to-detector distance (in units length)
        height: int,  # Y-direction length (in units pixels)
        width: int,  # X-direction length (in units pixels)
        delx: float,  # X-direction spacing (in units length / pixel)
        dely: float,  # Y-direction spacing (in units length / pixel)
        x0: float,  # Principal point x-coordinate (in units length)
        y0: float,  # Principal point y-coordinate (in units length)
        reorient: torch.Tensor,  # Frame-of-reference change matrix
        n_subsample: int | None = None,  # Number of target points to randomly sample
        reverse_x_axis: bool = False,  # If pose includes reflection (in E(3) not SE(3)), reverse x-axis
    ):
        super().__init__()
        self.height = height
        self.width = width
        self.n_subsample = n_subsample
        if self.n_subsample is not None:
            self.subsamples = []
        self.reverse_x_axis = reverse_x_axis

        # Initialize the source and detector plane in default positions (along the x-axis)
        source, target = self._initialize_carm()
        self.register_buffer("source", source)
        self.register_buffer("target", target)

        # Create a pose to reorient the scanner
        self.register_buffer("_reorient", reorient)

        # Create a calibration matrix that holds the detector's intrinsic parameters
        self.register_buffer(
            "_calibration",
            torch.tensor(
                [
                    [delx, 0, 0, x0],
                    [0, dely, 0, y0],
                    [0, 0, sdd, 0],
                    [0, 0, 0, 1],
                ]
            ),
        )

    @property
    def sdd(self):
        return self._calibration[2, 2].item()

    @property
    def delx(self):
        return self._calibration[0, 0].item()

    @property
    def dely(self):
        return self._calibration[1, 1].item()

    @property
    def x0(self):
        return -self._calibration[0, -1].item()

    @property
    def y0(self):
        return -self._calibration[1, -1].item()

    @property
    def reorient(self):
        return RigidTransform(self._reorient)

    @property
    def calibration(self):
        """A 4x4 matrix that rescales the detector plane to world coordinates."""
        return RigidTransform(self._calibration)

    @property
    def intrinsic(self):
        """The 3x3 intrinsic matrix."""
        return make_intrinsic_matrix(self).to(self.source)
```

--------------------------------

### Get SE(3) V Matrix Inputs

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/api/06_pose.ipynb

Helper function to compute the necessary inputs for the `_se3_V_matrix` function from the logarithm of the rotation. It calculates rotation angles, hat operator, and squared hat operator.

```python
def _get_se3_V_input(log_rotation: torch.Tensor, eps: float = 1e-4):
    """
    A helper function that computes the input variables to the `_se3_V_matrix`
    function.
    """
    # pyre-fixme[58]: `**` is not supported for operand types `Tensor` and `int`.
    nrms = (log_rotation**2).sum(-1)
    rotation_angles = torch.clamp(nrms, eps).sqrt()
    log_rotation_hat = hat(log_rotation)
    log_rotation_hat_square = torch.bmm(log_rotation_hat, log_rotation_hat)
    return log_rotation, log_rotation_hat, log_rotation_hat_square, rotation_angles
```

--------------------------------

### Initialize DRR with Trilinear Renderer

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/tutorials/trilinear.ipynb

Initializes the DRR object with the 'trilinear' rendering mode. Ensure the device is set correctly.

```python
drr = DRR(
    subject,
    sdd=1020.0,
    height=200,
    delx=2.0,
    renderer="trilinear",  # Switch the rendering mode
).to(device)

source, target = drr.detector(pose, calibration=None)
_ = drr(pose)  # Initialize drr.density
```

--------------------------------

### Utility Functions for DRR Plotting

Source: https://github.com/eigenvivek/diffdrr/blob/main/experiments/figure.ipynb

Defines utility functions to load example CT data, initialize a DRR object, and plot DRRs with associated loss information. Assumes 'cuda' is available for DRR computation.

```python
from diffdrr import load_example_ct, DRR, plot_drr


# Utility functions for plotting DRRs
volume, spacing = load_example_ct()
drr = DRR(volume, spacing, height=100, delx=10.0, device="cuda")


def plot_df(df, idx, ax, sdr=200):
    
    # Make the DRR
    params = df.iloc[idx][["theta", "phi", "gamma", "bx", "by", "bz"]].values
    img = drr(sdr, *params)
    
    # Make the title
    loss = df.iloc[idx]["loss"]
    if idx == 0:
        itr = "Initial Guess"
    elif idx == -1:
        itr = "Final Guess"
    else:
        itr = f"Iteration {idx}"
    title = f"{itr}\n-ZNCC={loss:.2f}"
    
    # Make the plot
    plot_drr(img, ax=ax, ticks=False)
    ax.set(xlabel=title)


def plot_a_few_drrs(df, int_1, int_2, axs):
    plot_df(df, 0, axs[0])
    plot_df(df, int_1, axs[1])
    plot_df(df, int_2, axs[2])
    plot_df(df, -1, axs[3])
```

--------------------------------

### Generate mesh from CT volume using SurfaceNets

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/tutorials/geometry.ipynb

Reads an example CT volume and generates a 3D mesh using the 'surface_nets' method. Requires DiffDRR and PyVista. The 'threshold' parameter controls the iso-surface extraction.

```python
# Read in the CT volume
subject = load_example_ct()

# Make a mesh from the CT volume
ct = drr_to_mesh(subject, "surface_nets", threshold=225, verbose=True)
```

--------------------------------

### Initialize and Use MultiscaleNormalizedCrossCorrelation2d

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/api/05_metrics.ipynb

Shows how to initialize and use the MultiscaleNormalizedCrossCorrelation2d metric with specified patch sizes and weights.

```python
msncc = MultiscaleNormalizedCrossCorrelation2d(
    patch_sizes=[9, None], patch_weights=[0.5, 0.5]
)
msncc(x1, x2)
```

--------------------------------

### Load Example CT Scan Data

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/api/03_data.ipynb

Loads a sample chest CT scan along with its corresponding labelmap and structure information. Supports specifying labels to render and C-arm orientation. Bone attenuation multiplier can also be adjusted.

```python
def load_example_ct(
    labels=None,
    orientation="AP",
    bone_attenuation_multiplier=1.0,
    **kwargs,
) -> Subject:
    """Load an example chest CT for demonstration purposes."""
    datadir = Path(__file__).resolve().parent / "data"
    volume = datadir / "cxr.nii.gz"
    labelmap = datadir / "mask.nii.gz"
    structures = pd.read_csv(datadir / "structures.csv")
    return read(
        volume,
        labelmap,
        labels,
        orientation=orientation,
        bone_attenuation_multiplier=bone_attenuation_multiplier,
        structures=structures,
        **kwargs,
    )
```

--------------------------------

### Initialize DRR and Estimate Parameters

Source: https://github.com/eigenvivek/diffdrr/blob/main/experiments/reg_test.ipynb

Initializes the DRR model for estimation with sparse subsampling and generates an initial estimate. It also prints the geodesic distance between the ground truth and the estimated DRR.

```python
# Initialize the DRR and optimization parameters
sdr, theta, phi, gamma, bx, by, bz = get_initial_parameters(true_params)
drr_est = DRR(volume, spacing, height=100, delx=10.0, p_subsample=0.1, device="cuda")
estimate = drr_est(sdr, theta, phi, gamma, bx, by, bz)

print(geodesic(drr_gt, drr_est))

plot_drr(estimate)
plt.show()
```

--------------------------------

### load_example_ct

Source: https://context7.com/eigenvivek/diffdrr/llms.txt

Loads a pre-packaged chest CT scan as a torchio.Subject. It applies HU-to-density conversion, reorients the volume, and centers the isocenter. Optional label selection and bone attenuation scaling are supported.

```APIDOC
## `load_example_ct` — Load a bundled chest CT for testing

Loads a pre-packaged chest CT scan (NIfTI format) as a `torchio.Subject` ready for use with DiffDRR. Applies HU-to-density conversion, reorients the volume to the chosen frame of reference, and centers the isocenter at the world origin. Accepts optional label selection and bone attenuation scaling.

### Usage

```python
from diffdrr.data import load_example_ct

# Load with default AP (anteroposterior) orientation
subject = load_example_ct()

# Load with PA orientation and boosted bone contrast
subject_pa = load_example_ct(orientation="PA", bone_attenuation_multiplier=2.0)

# Load with a specific label subset from the bundled labelmap
subject_bone = load_example_ct(labels=[1, 2, 3])  # e.g., specific structure indices

print(subject.volume.shape)     # torch.Size([1, D, H, W]) raw CT tensor
print(subject.density.shape)    # torch.Size([1, D, H, W]) normalized density
print(subject.reorient)         # 4x4 reorientation matrix (torch.Tensor)
```
```

--------------------------------

### Sample Volume at XYZ Coordinates

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/api/01_renderers.ipynb

Wraps `torch.nn.functional.grid_sample` to sample a 3D volume at specified XYZ coordinates. Handles optional integration with an existing image tensor.

```python
def _get_voxel(volume, xyzs, img, mode, align_corners):
    """Wraps torch.nn.functional.grid_sample to sample a volume at XYZ coordinates."""
    batch_size = len(xyzs)
    voxels = grid_sample(
        input=volume.permute(2, 1, 0)[None, None].expand(batch_size, -1, -1, -1, -1),
        grid=xyzs,
        mode=mode,
        align_corners=align_corners,
    )[:, 0, 0]
    if img is not None:
        img = torch.einsum("bcn, bnj -> bnj", img, voxels)
    else:
        img = voxels
    return img
```

--------------------------------

### Initialize DRR module and get detector geometry

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/tutorials/trilinear.ipynb

Initializes the DRR module with the subject CT data, source-to-detector distance, and detector parameters. It then calculates the source and target points for ray casting, transforming them into the volume's coordinate system.

```python
# Initialize the DRR module for generating synthetic X-rays
drr = DRR(
    subject,
    sdd=1020.0,
    height=200,
    delx=2.0,
).to(device)
_ = drr(pose)  # Initialize drr.density

source, target = drr.detector(pose, calibration=None)
source = drr.affine_inverse(source)
target = drr.affine_inverse(target)
```

--------------------------------

### Initialize Pose Regressor Backbone

Source: https://context7.com/eigenvivek/diffdrr/llms.txt

Builds a pose regression network using a timm image backbone with linear heads for rotation and translation. Allows direct prediction of camera pose from an X-ray image.

```python
from diffdrr.registration import PoseRegressor
import torch

regressor = PoseRegressor(
    model_name="resnet18",
    parameterization="rotation_6d",
    convention=None,
    pretrained=True,
    height=256,
)
regressor.eval()

```

--------------------------------

### Time DRR Computation with Patching (Height 750)

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/tutorials/timing.ipynb

Initializes a DRR object with a height of 750 pixels and a patch size of 150. This example further explores the impact of patch size on performance for larger DRR resolutions, highlighting memory management strategies.

```python
# |cuda
height = 750
patch_size = 150

drr = DRR(subject, sdd=1020, height=height, delx=2.0, patch_size=patch_size).to(device=device, dtype=torch.float32)
%timeit drr(pose)
del drr
```

--------------------------------

### Time DRR Computation with Patching (Height 1500)

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/tutorials/timing.ipynb

Initializes a DRR object with a height of 1500 pixels and a patch size of 250. This example pushes the limits of DRR generation size with patching, illustrating that computation time becomes the primary constraint rather than memory.

```python
# |cuda
height = 1500
patch_size = 250

drr = DRR(subject, sdd=1020, height=height, delx=2.0, patch_size=patch_size).to(device=device, dtype=torch.float32)
%timeit drr(pose)
del drr
```

--------------------------------

### Plotting DRR Optimization Trajectories and Convergence

Source: https://github.com/eigenvivek/diffdrr/blob/main/experiments/figure.ipynb

Generates plots visualizing DRR optimization loss curves and convergence timing histograms. It analyzes multiple runs, calculates convergence metrics, and displays example DRRs for successful and failed convergence cases. Requires pandas, seaborn, and matplotlib.

```python
fig = plt.figure(dpi=100, figsize=(9, 12), constrained_layout=True)
subfigs = fig.subfigures(4, 1, wspace=0.1, height_ratios=[1.15, 1, 1, 1])

charts = subfigs[0].subfigures(1, 2, width_ratios=[2, 1])
ax_plot = charts[0].subplots()
ax_hist = charts[1].subplots()

# Plot the training trajectories
dfs = []
kwargs = {"alpha": 0.1, "lw": 1.25}
for file in Path("results/momentum_dampen/runs/").glob("*.csv"):
    df = pd.read_csv(file)
    df["runtime"] = df["time"].cumsum()
    df["converged"] = converged(df, cutoff=-0.999, max_itr=250)
    color = colors[converged(df, cutoff=-0.999, max_itr=250)]
    ax_plot.plot(df["itr"], df["loss"], color, **kwargs)
    dfs.append(df.iloc[-1])
ax_plot.set(xlabel="Iteration", ylabel="Negative ZNCC")
ax_plot.set_title(f"(a) Loss curves for {len(dfs)} simulated DRRs", loc="left")

# Plot the convergence timing histogram
df = pd.concat(dfs, axis=1).T.reset_index(drop=True)
sns.histplot(data=df, x="runtime", hue="converged", hue_order=[True, False], multiple="stack", palette=colors, legend=True, ax=ax_hist)
sns.move_legend(ax_hist, loc="upper left", bbox_to_anchor=(-1, 1), title="Converged", frameon=False)
ax_hist.set(xlabel="Runtime (s)", xlim=(0, 8.5))
ax_hist.set_title("(b) Total optimization time", loc="left")

# Get optimization timing metrics
n_conv = df["converged"].astype(int).sum()
t_conv = df.query("converged == True")["runtime"].mean()
s_conv = df.query("converged == True")["runtime"].std()
nitr_conv = df.query("converged == True")["itr"].mean()
sitr_conv = df.query("converged == True")["itr"].std()
print(f"{n_conv}/{len(df)} ({n_conv/len(df):.2%}) runs achieved convergence in {t_conv:.2f}±{s_conv:.2f}s ({nitr_conv:.2f}±{sitr_conv:.2f} iterations) on average")

# Plot a few examples of convergence
subfigs[1].suptitle(
    "(c) Examples of successful convergence",
    x=0.0125,
    y=1,
    horizontalalignment='left',
    verticalalignment='top',
    fontsize=13,
)
axs = subfigs[1].subplots(ncols=4)
df = pd.read_csv("results/momentum_dampen/runs/34.csv")
plot_a_few_drrs(df, 20, 40, axs)

axs = subfigs[2].subplots(ncols=4)
df = pd.read_csv("results/momentum_dampen/runs/22.csv")
plot_a_few_drrs(df, 20, 40, axs)

# And one example of failure to converge
subfigs[3].suptitle(
    "(d) Examples of failure to converge",
    x=0.0125,
    y=1,
    horizontalalignment='left',
    verticalalignment='top',
    fontsize=13,
)
axs = subfigs[3].subplots(ncols=4)
df = pd.read_csv("results/momentum_dampen/runs/59.csv")
plot_a_few_drrs(df, 20, 40, axs)

# Save the final figure!
plt.savefig("../../figures/optimization.pdf", bbox_inches="tight")
plt.show()
```

--------------------------------

### Siddon Class Initialization

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/api/01_renderers.ipynb

Initializes the Siddon class with various parameters controlling the rendering process, such as voxel shift, interpolation mode, and gradient handling.

```APIDOC
## Siddon

### Description
Differentiable X-ray renderer implemented with Siddon's method for exact raytracing.

### Parameters
- **voxel_shift** (float) - Optional - 0 or 0.5, depending if the voxel is at the top left corner or the center. Defaults to 0.5.
- **mode** (str) - Optional - Interpolation mode for grid_sample. Defaults to "nearest".
- **stop_gradients_through_grid_sample** (bool) - Optional - Apply torch.no_grad when calling grid_sample. Defaults to False.
- **filter_intersections_outside_volume** (bool) - Optional - Use alphamin/max to filter the intersections. Defaults to False.
- **reducefn** (str) - Optional - Function for combining samples along each ray. Defaults to "sum".
- **eps** (float) - Optional - Small constant to avoid div by zero errors. Defaults to 1e-8.
```

--------------------------------

### Detector Class Initialization

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/api/02_detector.ipynb

Initializes a 6 DoF X-ray detector system. This model is based on a C-Arm.

```APIDOC
## Detector

### Description
Construct a 6 DoF X-ray detector system. This model is based on a C-Arm.

### Parameters

- **sdd** (float) - Source-to-detector distance (in units length)
- **height** (int) - Y-direction length (in units pixels)
- **width** (int) - X-direction length (in units pixels)
- **delx** (float) - X-direction spacing (in units length / pixel)
- **dely** (float) - Y-direction spacing (in units length / pixel)
- **x0** (float) - Principal point x-coordinate (in units length)
- **y0** (float) - Principal point y-coordinate (in units length)
- **reorient** (torch.Tensor) - Frame-of-reference change matrix
- **n_subsample** (int | None, optional) - Number of target points to randomly sample. Defaults to None.
- **reverse_x_axis** (bool, optional) - If pose includes reflection (in E(3) not SE(3)), reverse x-axis. Defaults to False.
```

--------------------------------

### Common nbdev Development Commands

Source: https://github.com/eigenvivek/diffdrr/blob/main/README.md

Essential commands for nbdev development, including previewing docs, cleaning, testing, and exporting the package.

```zsh
nbdev_preview  # Render docs locally and inspect in browser
nbdev_clean    # NECESSARY BEFORE PUSHING
nbdev_test     # tests notebooks
nbdev_export   # builds package and builds docs
```

--------------------------------

### Supported Parameterizations for SE(3) Conversion

Source: https://github.com/eigenvivek/diffdrr/blob/main/notebooks/api/06_pose.ipynb

Lists the supported parameterizations for converting SE(3) transformations.

```python
#| exporti
PARAMETERIZATIONS = [
    "axis_angle",
    "euler_angles",

```