Custom map simulations¶
In this tutorial we will build a simulation from scratch.
We start by defining a Band that will determine our array’s sensitivity to different spectra. We then generate an array by specifying a field of view, which will be populated by evenly-spaced beams of the given band.
[1]:
import maria
from maria.instrument import Band
f090 = Band(
center=90e9, # in Hz
width=20e9, # in Hz
NET_RJ=40e-6, # in K sqrt(s)
knee=1e0, # in Hz
gain_error=5e-2)
f150 = Band(
center=150e9,
width=30e9,
NET_RJ=60e-6,
knee=1e0,
gain_error=5e-2)
[2]:
array = {"field_of_view": 0.5,
"shape": "circle",
"beam_spacing": 1.5,
"primary_size": 25,
"bands": [f090, f150]}
instrument = maria.get_instrument(array=array)
print(instrument)
instrument.plot()
Instrument(1 array)
├ arrays:
│ n FOV baseline bands polarized
│ array1 2116 29.85’ 0 m [f090,f150] False
│
└ bands:
name center width η NEP NET_RJ NET_CMB FWHM
0 f090 90 GHz 20 GHz 0.5 5.445 aW√s 40 uK√s 49.13 uK√s 0.5832’
1 f150 150 GHz 30 GHz 0.5 12.25 aW√s 60 uK√s 104 uK√s 21”
Here, the fetch function downloads a map to the path map_filename, but map_filename can be any .h5 or .fits file of an image that corresponds to the maria map convention (see Maps).
[3]:
from maria.io import fetch
map_filename = fetch("maps/cluster1.fits")
input_map = maria.map.load(
filename=map_filename,
nu=150e9,
center=(291.156, -31.23))
input_map.data *= 1e1
print(input_map)
input_map.to("K_RJ").plot()
Downloading https://github.com/thomaswmorris/maria-data/raw/master/maps/cluster1.fits: 100%|██████████| 4.20M/4.20M [00:00<00:00, 162MB/s]
ProjectedMap:
shape(nu, y, x): (1, 1024, 1024)
stokes: naive
nu: [150.] GHz
t: naive
z: naive
quantity: spectral_flux_density_per_pixel
units: Jy/pixel
min: -3.845e-05
max: -5.005e-08
center:
ra: 19ʰ24ᵐ37.44ˢ
dec: -31°13’48.00”
size(y, x): (1°, 1°)
resolution(y, x): (3.516”, 3.516”)
beam(maj, min, rot): (0°, 0°, 0°)
memory: 16.78 MB
[4]:
site = maria.get_site("llano_de_chajnantor", altitude=5065)
print(site)
site.plot()
Site:
region: chajnantor
location: 23°01’45.84”S 67°45’17.28”W
altitude: 5065 m
seasonal: True
diurnal: True
Downloading https://github.com/thomaswmorris/maria-data/raw/master/world_heightmap.h5: 100%|██████████| 7.34M/7.34M [00:00<00:00, 191MB/s]
[5]:
plan = maria.Plan(
start_time="2024-08-06T03:00:00",
scan_pattern="daisy",
scan_options={"radius": 0.5, "speed": 0.1}, # in degrees
duration=900, # in seconds
sample_rate=50, # in Hz
scan_center=(291.156, -31.23),
frame="ra_dec")
print(plan)
plan.plot()
Plan:
start_time: 2024-08-06 03:00:00.000 +00:00
duration: 900 s
sample_rate: 50 Hz
center:
ra: 19ʰ24ᵐ37.44ˢ
dec: -31°13’48.00”
scan_pattern: daisy
scan_radius: 0.9989°
scan_kwargs: {'radius': 0.5, 'speed': 0.1}
[6]:
sim = maria.Simulation(
instrument,
plan=plan,
site=site,
atmosphere="2d",
atmosphere_kwargs={"weather": {"pwv": 0.5}},
map=input_map)
print(sim)
Downloading https://github.com/thomaswmorris/maria-data/raw/master/atmosphere/spectra/am/v2/chajnantor.h5: 100%|██████████| 14.3M/14.3M [00:00<00:00, 185MB/s]
Downloading https://github.com/thomaswmorris/maria-data/raw/master/atmosphere/weather/era5/chajnantor.h5: 100%|██████████| 8.00M/8.00M [00:00<00:00, 164MB/s]
Constructing atmosphere: 100%|██████████| 10/10 [00:58<00:00, 5.81s/it]
Simulation
├ Instrument(1 array)
│ ├ arrays:
│ │ n FOV baseline bands polarized
│ │ array1 2116 29.85’ 0 m [f090,f150] False
│ │
│ └ bands:
│ name center width η NEP NET_RJ NET_CMB FWHM
│ 0 f090 90 GHz 20 GHz 0.5 5.445 aW√s 40 uK√s 49.13 uK√s 0.5832’
│ 1 f150 150 GHz 30 GHz 0.5 12.25 aW√s 60 uK√s 104 uK√s 21”
├ Site:
│ region: chajnantor
│ location: 23°01’45.84”S 67°45’17.28”W
│ altitude: 5065 m
│ seasonal: True
│ diurnal: True
├ Plan:
│ start_time: 2024-08-06 03:00:00.000 +00:00
│ duration: 900 s
│ sample_rate: 50 Hz
│ center:
│ ra: 19ʰ24ᵐ37.44ˢ
│ dec: -31°13’48.00”
│ scan_pattern: daisy
│ scan_radius: 0.9989°
│ scan_kwargs: {'radius': 0.5, 'speed': 0.1}
├ Atmosphere(10 processes with 10 layers):
│ ├ spectrum:
│ │ region: chajnantor
│ └ weather:
│ region: chajnantor
│ altitude: 5065 m
│ time: Aug 5 23:07:29 -04:00
│ pwv[mean, rms]: (0.5 mm, 15 um)
└ ProjectedMap:
shape(stokes, nu, y, x): (1, 1, 1024, 1024)
stokes: I
nu: [150.] GHz
t: naive
z: naive
quantity: rayleigh_jeans_temperature
units: K_RJ
min: -1.915e-04
max: -2.492e-07
center:
ra: 19ʰ24ᵐ37.44ˢ
dec: -31°13’48.00”
size(y, x): (1°, 1°)
resolution(y, x): (3.516”, 3.516”)
beam(maj, min, rot): (0°, 0°, 0°)
memory: 16.78 MB
[7]:
tod = sim.run()
print(tod)
tod.plot()
Generating turbulence: 100%|██████████| 10/10 [00:00<00:00, 13.20it/s]
Sampling turbulence: 100%|██████████| 10/10 [00:21<00:00, 2.13s/it]
Computing atmospheric emission: 100%|██████████| 2/2 [00:15<00:00, 7.84s/it, band=f150]
Sampling map: 100%|██████████| 2/2 [00:28<00:00, 14.40s/it, channel=[ 0. inf] Hz]
Generating noise: 100%|██████████| 2/2 [00:03<00:00, 1.53s/it, band=f150]
TOD(shape=(2116, 45000), fields=['atmosphere', 'map', 'noise'], units='pW', start=2024-08-06 03:14:59.979 +00:00, duration=900.0s, sample_rate=50.0Hz, metadata={'atmosphere': True, 'sim_time': <Arrow [2025-06-20T14:46:53.364164+00:00]>, 'altitude': 5065.0, 'region': 'chajnantor', 'pwv': 0.5, 'base_temperature': 272.523})
[8]:
from maria.mappers import BinMapper
mapper = BinMapper(
center=(291.156, -31.23),
frame="ra_dec",
width=1.,
height=1.,
resolution=1. / 256,
tod_preprocessing={
"window": {"name": "tukey", "kwargs": {"alpha": 0.1}},
"remove_spline": {"knot_spacing": 30, "remove_el_gradient": True},
"remove_modes": {"modes_to_remove": [0]},
},
map_postprocessing={
"gaussian_filter": {"sigma": 1},
},
units="mK_RJ",
)
mapper.add_tods(tod)
output_map = mapper.run()
Mapping band f090: 100%|██████████| 1/1 [00:00<00:00, 2.46it/s, band=f090, stokes=I]
2025-06-20 14:47:19.140 INFO: Ran mapper for band f090 in 20.25 s.
Mapping band f150: 100%|██████████| 1/1 [00:00<00:00, 2.43it/s, band=f150, stokes=I]
2025-06-20 14:47:38.143 INFO: Ran mapper for band f150 in 19 s.
We can see the recovered map with
[9]:
output_map.plot(nu_index=[0, 1])
