How to Apply Instrumental Broadening Effects to Simulated Spectra in ExoJAX

To apply instrumental broadening in ExoJAX, convert your target spectral resolution R to a Gaussian velocity dispersion using resolution_to_gaussian_std, instantiate SopInstProfile with your wavenumber grid, and convolve your high-resolution spectrum using ipgauss.

ExoJAX treats instrumental broadening as a Gaussian convolution applied after the radiative-transfer step. When you need to compare synthetic exoplanet spectra against real observations, knowing how to apply instrumental broadening effects to simulated spectra in ExoJAX ensures your models match the resolving power of instruments like IGRINS, CRIRES+, or JWST NIRSpec.

Understanding the Instrumental Broadening Workflow

The pipeline follows a discrete mathematical progression from physical resolution to convolved flux. First, define your target instrument's resolving power (R = \lambda/\Delta\lambda). Next, translate this resolving power into a velocity-space standard deviation (km s⁻¹) using the built-in utility. Finally, instantiate a spectral operator that constructs the convolution kernel and applies it to your high-resolution model.

The architecture relies on two primary modules:

• exojax/utils/instfunc.py: Contains the resolution-to-dispersion conversion logic
• exojax/postproc/specop.py: Implements the SopInstProfile class that performs the convolution

Converting Spectral Resolution to Gaussian Sigma

Before convolving, you must convert the dimensionless resolving power (R) into a physical velocity width. The function resolution_to_gaussian_std in exojax/utils/instfunc.py (lines 15-26) implements the standard relation:

[ \sigma_v = \frac{c}{2\sqrt{2\ln2},R} ]

This returns the Gaussian standard deviation in km s⁻¹ required by the convolution kernel.

from exojax.utils.instfunc import resolution_to_gaussian_std

# Convert R = 100,000 to velocity dispersion

R_inst = 100_000.0
sigma_ip = resolution_to_gaussian_std(R_inst)  # Returns ~0.127 km/s

The SopInstProfile Class

The SopInstProfile class, defined in exojax/postproc/specop.py (lines 44-77), manages the velocity grid and convolution infrastructure. When instantiated, it automatically determines the spectral resolution of your input nu_grid and constructs a velocity array spanning (\pm) vrmax km s⁻¹.

Key parameters:

• nu_grid: The high-resolution wavenumber grid from your radiative transfer calculation
• vrmax: Maximum velocity range in km s⁻¹ (should exceed several times the Gaussian sigma)

The class provides two essential methods:

• ipgauss(spectrum, sigma): Performs the Gaussian convolution using the internally built velocity grid
• sampling(broadened_spectrum, radial_velocity, nu_grid_sampling): Resamples the convolved spectrum onto the exact instrumental wavenumber grid, accounting for Doppler shifts

Convolution Methods and Performance

SopInstProfile inherits from SopCommonConv and supports two convolution algorithms selected via the convolution_method parameter:

• "exojax.signal.convolve": The default FFT-based JAX convolution (optimal for most spectra)
• "exojax.signal.ola": An overlap-add implementation for very long spectra where memory constraints prevent standard FFT-based approaches

The ipgauss method (lines 65-74 in specop.py) automatically handles the dispatch between these backends based on your initialization settings.

Minimal Working Example

This example demonstrates the complete workflow from high-resolution synthesis to instrumentally broadened output:

import numpy as np
import jax.numpy as jnp
from exojax.utils.grids import wavenumber_grid
from exojax.utils.instfunc import resolution_to_gaussian_std
from exojax.postproc.specop import SopInstProfile

# 1. Define a high-resolution wavenumber grid

nu_min, nu_max, ngrid = 2000.0, 2500.0, 50000
nu_grid, _, _ = wavenumber_grid(nu_min, nu_max, ngrid, xsmode="premodit")

# 2. Generate a dummy high-resolution spectrum

Tplanet = 1200.0
sigma = 5.670374e-5
F_high = sigma * Tplanet**4 * jnp.ones_like(nu_grid)

# 3. Convert instrumental resolution to Gaussian sigma

R_inst = 100_000.0
sigma_ip = resolution_to_gaussian_std(R_inst)

# 4. Create the instrument operator

sop_inst = SopInstProfile(nu_grid, vrmax=5.0)

# 5. Apply Gaussian instrumental broadening

F_broadened = sop_inst.ipgauss(F_high, sigma_ip)

# 6. Resample onto observational grid

nu_inst, _, _ = wavenumber_grid(nu_min, nu_max, 3000, xsmode="eslog")
F_sampled = sop_inst.sampling(F_broadened, radial_velocity=0.0, nu_grid_sampling=nu_inst)

Full Retrieval Workflow

In practice, instrumental broadening is often combined with rotational broadening. This excerpt from examples/spec_with_photo.py shows the canonical ordering of operations:

from exojax.postproc.specop import SopInstProfile, SopRotation
from exojax.utils.instfunc import resolution_to_gaussian_std

# Instrumental settings

Rinst = 100_000.0
beta_inst = resolution_to_gaussian_std(Rinst)

# Build spectral operators

sop_rot = SopRotation(nu_grid_spec, vsini_max=100.0)
sop_inst = SopInstProfile(nu_grid_spec, vrmax=100.0)

# Apply physical effects in sequence

flux_rot = sop_rot.rigid_rotation(flux, vsini=5.0, u1=0.0, u2=0.0)
flux_ip = sop_inst.ipgauss(flux_rot, beta_inst)

# Sample to observational grid with radial velocity shift

flux_obs = sop_inst.sampling(flux_ip, radial_velocity=0.0, nu_grid_sampling=nu_grid_obs)

Note that SopRotation is applied before SopInstProfile to ensure the rotational broadening (which has extended wings) is properly convolved with the instrumental profile.

Summary

• Convert resolution to velocity: Use resolution_to_gaussian_std in exojax/utils/instfunc.py to translate resolving power (R) into km s⁻¹
• Initialize the operator: Instantiate SopInstProfile with your model nu_grid and a vrmax value exceeding the Gaussian sigma by a factor of 3-5
• Convolve: Call ipgauss(spectrum, sigma) to apply the Gaussian instrumental profile using JAX-optimized FFT or overlap-add methods
• Resample: Use sampling() to interpolate onto the exact instrumental grid and account for radial velocity shifts
• Chain effects: Apply rotational broadening (SopRotation) before instrumental broadening for physical accuracy

Frequently Asked Questions

What is the mathematical relationship between spectral resolution and Gaussian sigma in ExoJAX?

ExoJAX uses the standard Gaussian line-shape definition where (\sigma_v = c / (2\sqrt{2\ln2},R)). The function resolution_to_gaussian_std in exojax/utils/instfunc.py implements this exactly, returning values in km s⁻¹. For a resolution of (R = 100,000), this yields approximately 1.27 km s⁻¹ FWHM, or 0.54 km s⁻¹ standard deviation.

How do I handle very long spectra without memory issues?

Set convolution_method="exojax.signal.ola" when instantiating SopInstProfile. This invokes the overlap-add algorithm instead of the standard FFT-based convolution, reducing memory consumption for extremely long wavenumber grids at a slight computational cost.

Can I combine instrumental broadening with rotational broadening?

Yes. Instantiate both SopRotation and SopInstProfile with the same nu_grid, then chain the operations: first apply sop_rot.rigid_rotation() to your high-resolution flux, then pass the result to sop_inst.ipgauss(). This ordering ensures the rotational profile (which has wider wings than typical instrument profiles) is properly convolved.

What value should I use for the vrmax parameter?

Set vrmax to at least 3-5 times the Gaussian sigma returned by resolution_to_gaussian_std. For (R = 100,000) (where (\sigma \approx 0.54) km s⁻¹), a vrmax of 5.0 km s⁻¹ is sufficient. If vrmax is too small, the convolution kernel will truncate prematurely, causing artifacts in the wings of spectral lines.

{
  "mcpServers": {
    "instagit": {
      "command": "npx",
      "args": ["-y", "instagit@latest"]
    }
  }
}