Implementation and simulation of the pointing reconstruction model for the LSPE/Strip telescope.

This repository, together with PRMaps, contains the code used for the simulations in my thesis. Basically, it is a Julia port of some functions of the following Python libraries:

• PySM: which allows to simulate a map of the sky to observe.
• FGBuster: which allows you to perform a component separation analysis using frequencies associated with the following surveys: LSPE/Strip, LSPE/SWIPE, Planck, Quijote.

System requirements

Once the Julia environment has been successfully activated, it is necessary to configure a Python environment containing PySM (v.3.3.2) and FGBuster (v2.0.0). It will then be necessary to set PyCall.jl correctly by indicating the path to the correct Python environment (see the documentation here).

Generate a sky signal

This pipeline can be used to generate Healpix map in equatorial coordinates like the one produced by PySM:

# Choose the sky model to simulate
sky_model = "d0s0"

# Choose the frequencies by using one or more instruments:
strip = read_pickle("../instruments/lspe_strip_instrument.pkl")

# Choose the maps resolution
nside = 512

# Generate the maps
maps = get_foreground_maps(
    instruments
    sky_model, 
    nside
)

Simple component separation example

Below I report a simple example containing a component separation analysis implemented by introducing a systematic pointing error in the 43 GHz LSPE/Strip channel.

First it is necessary to choose which instruments you have to simulate; all of the aviable instruments are stored into .pkl file into the instruments directory:

using Pandas

swipe = read_pickle("../instruments/lspe_swipe_instrument.pkl")
strip = read_pickle("../instruments/lspe_strip_instrument.pkl")

instruments = concat([strip, swipe])

Then you have to decide the simulation setup for the LSPE/Strip telescope:

using PRMaps
using Dates
import Stripeline as Sl

# Choose a Healpix map resolution
nside = 128

# Choose the sky model to observe
sky_model = "c1s0d0"

# Choose the starting day of the simulation and the number of days to simulate 
t_start = DateTime(2022, 1, 1, 12, 0, 0)
days = 10
    
setup = PRMaps.Setup(
    sampling_freq_Hz = 50.0,
    total_time_s = 24. * 3600. * days
)

# Choose the detector to simulate (in this case H0)
cam_ang = Sl.CameraAngles()

# Choose the uncertainties on the configuration angles
tel_ang = Sl.TelescopeAngles(forkang_rad = deg2rad(1.0/60))

Then you can run the simulation that reconstruct the observed signal and perform a component separation analisys:

(result_dict, result_map) = run_fgbuster_with_error(
    instruments, 
    cam_ang, 
    tel_ang, 
    setup, 
    sky_model, 
    nside, 
    t_start
    )

At this point result_dict contain a dictionary with the results generated by FGBUster and result_map contains the maps of the observed signal in the sky.

Additional usage examples can be found in the notebooks folder.