Simulating observations with MUSTANG-2

MUSTANG-2 is a bolometric array on the Green Bank Telescope. In this notebook we simulate an observation of the Crab Nebula (M1).

import maria

input_map = maria.map.get("maps/M1.h5").to("uK_RJ")

input_map.plot(slices=dict(nu=[0, 1, 2]))
print(input_map)

2026-06-02 15:49:15.106 INFO: Fetching https://github.com/thomaswmorris/maria-data/raw/master/maps/M1.h5
Downloading: 100%|██████████| 21.8M/21.8M [00:00<00:00, 169MB/s]

ProjectionMap:
  data(1, 3, 1205, 1187):
    min: -4.145e+09
    max: 6.287e+08
    units: uK_RJ
    quantity: rayleigh_jeans_temperature
  stokes(1):
    components: ['I']
  nu(3):
    values: [149.8962 214.1375 272.5386] GHz
  eta(1205):
    height: 20.07’
    res: -1”
  xi(1187):
    width: 19.77’
    res: 1”
  frame: ra/dec
  center:
    ra: 05ʰ34ᵐ31.95ˢ
    dec: 22°00’52.16”
  beam(maj, min, psi): ragged
  memory: 34.33 MB

from maria import Planner

planner = Planner(target=input_map, site="green_bank", constraints={"el": (60, 90)})
plans = planner.generate_plans(total_duration=600, sample_rate=100)

plans[0].plot()
print(plans)

PlanList(1 plans, 600 s):
                           start_time duration    target(ra,dec)     center(az,el)
chunk
0      2026-06-02 16:13:43.462 +00:00    600 s  (83.63°, 22.01°)  (116.2°, 60.94°)

instrument = maria.get_instrument("MUSTANG-2")

print(instrument)
instrument.plot()

Instrument(1 array)
├ arrays:
│            n   FOV baseline      bands polarized
│  array1  217  4.2’      0 m  [m2/f093]     False
│
└ bands:
         name     center      width    η      NEP         NET_RJ         NET_CMB    FWHM
   0  m2/f093  86.21 GHz  20.98 GHz  0.1  15 aW√s  571.1 uK_RJ√s  690.5 uK_CMB√s  9.133”

sim = maria.Simulation(
    instrument,
    plans=plans,
    site="green_bank",
    map=input_map,
    atmosphere="2d",
)

print(sim)

2026-06-02 15:49:27.084 INFO: Fetching https://github.com/thomaswmorris/maria-data/raw/master/atmosphere/spectra/am/v3/green_bank.h5
Downloading: 100%|██████████| 22.0M/22.0M [00:00<00:00, 185MB/s]
2026-06-02 15:49:27.832 INFO: Fetching https://github.com/thomaswmorris/maria-data/raw/master/atmosphere/weather/era5/green_bank.h5
Downloading: 100%|██████████| 12.0M/12.0M [00:00<00:00, 138MB/s]

Simulation
├ Instrument(1 array)
│ ├ arrays:
│ │            n   FOV baseline      bands polarized
│ │  array1  217  4.2’      0 m  [m2/f093]     False
│ │
│ └ bands:
│          name     center      width    η      NEP         NET_RJ         NET_CMB    FWHM
│    0  m2/f093  86.21 GHz  20.98 GHz  0.1  15 aW√s  571.1 uK_RJ√s  690.5 uK_CMB√s  9.133”
├ Site:
│   region: green_bank
│   timezone: America/New_York
│   location:
│     longitude: 79°50’23.28” W
│     latitude:  38°25’59.16” N
│     altitude: 825 m
│   seasonal: True
│   diurnal: True
├ PlanList(1 plans, 600 s):
│                            start_time duration    target(ra,dec)     center(az,el)
│ chunk
│ 0      2026-06-02 16:13:43.462 +00:00    600 s  (83.63°, 22.01°)  (116.2°, 60.94°)
├ '2d'
└ ProjectionMap:
    data(1, 3, 1, 1205, 1187):
      min: -4.145e+09
      max: 6.287e+08
      units: uK_RJ
      quantity: rayleigh_jeans_temperature
    stokes(1):
      components: ['I']
    nu(3):
      values: [149.8962 214.1375 272.5386] GHz
    t(1):
      values: [1.78041535e+09] s
    eta(1205):
      height: 20.07’
      res: -1”
    xi(1187):
      width: 19.77’
      res: 1”
    frame: ra/dec
    center:
      ra: 05ʰ34ᵐ31.95ˢ
      dec: 22°00’52.16”
    beam(maj, min, psi): ragged
    memory: 34.33 MB

tods = sim.run()
tods[0].plot()

2026-06-02 15:49:28.482 INFO: Simulating observation 1 of 1
Constructing atmosphere: 100%|██████████| 8/8 [00:00<00:00, 13.37it/s]
Generating turbulence: 100%|██████████| 8/8 [00:00<00:00, 47.99it/s]
Sampling turbulence: 100%|██████████| 8/8 [00:03<00:00,  2.04it/s]
Computing atmospheric emission: 100%|██████████| 1/1 [00:00<00:00,  1.44it/s, band=m2/f093]
Sampling map: 100%|██████████| 1/1 [00:06<00:00,  6.49s/it, band=m2/f093, channel=(0 Hz, 182 GHz)]
Generating noise: 100%|██████████| 1/1 [00:01<00:00,  1.30s/it, band=m2/f093]
2026-06-02 15:49:50.395 INFO: Simulated observation 1 of 1 in 21.9 s

from maria.mappers import MaximumLikelihoodMapper

mapper = MaximumLikelihoodMapper(
    stokes="I",
    tod_preprocessing={
        "remove_spline": {"knot_spacing": 60, "remove_el_gradient_order": 1},
    },
    units="uK_RJ",
    tods=tods,
    resolution=2 * input_map.resolution,
)

mapper.map.plot()

2026-06-02 15:49:56.114 INFO: Inferring center {'ra': '05ʰ34ᵐ31.55ˢ', 'dec': '22°00’47.24”'} for mapper
2026-06-02 15:49:56.124 INFO: Inferring mapper width 23.66’ for mapper from observation patch
2026-06-02 15:49:56.124 INFO: Inferring mapper height 23.66’ to match supplied width
Preprocessing TODs: 100%|██████████| 1/1 [00:04<00:00,  4.87s/it]
Computing pointing matrices: 100%|██████████| 1/1 [00:02<00:00,  2.45s/it]

mapper.fit(epochs=2, steps_per_epoch=25)
mapper.map.plot()

Updating noise model: 100%|██████████| 1/1 [00:02<00:00,  2.95s/it, tod=1/1]
Fitting map (epoch 1/2): 100%|██████████| 25/25 [00:54<00:00,  2.17s/it, mll=-8.275e-12, pgrad=0.932, step=1.503e+21]
Updating noise model: 100%|██████████| 1/1 [00:02<00:00,  2.91s/it, tod=1/1]
Fitting map (epoch 2/2): 100%|██████████| 25/25 [00:56<00:00,  2.25s/it, mll=-6.398e-12, pgrad=0.549, step=1.326e+21]