HWP_sys¶
This module implements HWP non-idealities both using Jones formalism (as described in Giardiello et al. 2021) and the Mueller one. In the Jones formalism, a non-ideal HWP is described by
where:
\(h_1\) and \(h_2\) are the efficiencies, describing the deviation from the unitary transmission of light components \(E_x\), \(E_y\). In the ideal case, \(h_1 = h_2 = 0\);
\(\beta=\phi-\pi\), where \(\phi\) is the phase shift between the two directions. It accounts for variations of the phase difference between \(E_x\) and \(E_y\) with respect to the nominal value of \(\pi\) for an ideal HWP. In the ideal case, \(\beta=0\);
\(\zeta_{1,2}\) and \(\chi_{1,2}\) are amplitudes and phases of the off-diagonal terms, coupling \(E_x\) and \(E_y\). In practice, if the incoming wave is fully polarized along x(y), a spurious y(x) component would show up in the outgoing wave. In the ideal case, \(\zeta_{1,2}=\chi_{1,2}=0\).
In the Mueller formalism, we have a general matrix
which, in the ideal case, would be
In the most general case, the Jones non-ideal parameters and the Mueller matrix elements can vary inside a bandpass.
The class hwp_sys.HwpSys
is a contanier for the parameters of the HWP systematics.
It defines three methods:
hwp_sys.HwpSys.set_parameters()
which sets the defaults and handles the interface with the parameter file of the simulation. The relevant section is tagged by[hwp_sys]
. The HWP parameters can be passed both in the Jones and Mueller formalism. This choice is regulated by the flagmueller_or_jones
, which has to be set to either"mueller"
or"jones"
. In case the Jones formalism is chosen, it is converted automatically into the Mueller one through the function.hwp_sys.JonesToMueller
. There is also the possibility to pass precomputed input maps (as a numpy array) through themaps
argument. Otherwise, input maps are internally computed through the litebird_sim map based simulator (mbs). The argumentintegrate_in_band
sets whether band integration is performed in the TOD computation; ifbuilt_map_on_the_fly = True
, the map-making can be performed internally (instead of using the litebird_sim binner);correct_in_solver
sets whether non-ideal parameters can be used in the map-making (map-making assuming a non-ideal HWP, generally using different HWP non-ideal parameters than the one used in the TOD, representing our estimate of their true value);integrate_in_band_solver
regulates whether band integration is performed in the map-making (to compute the \(B^T B\) and \(B^T d\) terms, see below).hwp_sys.HwpSys.fill_tod()
which fills the tod in a given Observation. Thepointings
angles passed have to include no rotating HWP, since the effect of the rotating HWP to the polarization angle is included in the TOD computation. The TOD is computed performing this operation:\(d_{obs}\left(t_{i}\right)\,=\,\frac{\int d\nu\,\frac{\partial BB(\nu,T)}{\partial T_{CMB}}\,\tau\left(\nu\right)\,M_{i}^{TX}(\nu)\left(m_{CMB}+m_{FG}\left(\nu\right)\right)}{\int d\nu \frac{\partial BB(\nu,T)}{\partial T_{CMB}}\,\tau \left(\nu\right)}\),
where \(\tau(\nu)\) is the bandpass, \(\frac{\partial BB(\nu,T)}{\partial T_{CMB}}\) converts from CMB thermodynamic temperature to differential source intensity (see eq.8 of https://arxiv.org/abs/1303.5070) and \(M_{i}^{TX}(\nu)\) is the Mueller matrix element including the non-ideal HWP.
If
built_map_on_the_fly = True
, the code computes also\(m_{out} = {\,\left(\sum_{i} B_{i}^{T} B_{i} \right)^{-1} \left( \sum_{i} B_{i}^{T} d_{obs}(t_{i}) \right)}\),
where the map-making matrix \(B^X = \left(\frac{\int d\nu \,\frac{\partial BB(\nu,T)}{\partial T_{CMB}}\,\tau_{s}\left(\nu\right)\,M_{i,s}^{TX}(\nu)}{\int d\nu \frac{\partial BB(\nu,T)}{\partial T_{CMB}}\,\tau_{s}\left(\nu\right)}\,\right)\).
\(\tau_s(\nu)\) and \(M_{i,s}^{TX}(\nu)\) are the estimate of the bandpass and Mueller matrix elements used in the map-making.
hwp_sys.HwpSys.make_map()
which can bin the observations in a map. This is available only ifbuilt_map_on_the_fly
variable is set toTrue
. With this method, it is possible to include non-ideal HWP knowledge in the map-making procedure, so use that instead of the generallitebird_sim
binner if you want to do so.
Defining a bandpass profile in hwp_sys
¶
It is possible to define more complex bandpass profiles than a top-hat when using hwp_sys
.
This can be done both for the TOD computation (\(\tau\)) and the map-making procedure
(\(\tau_s\)). All you have to do is create a dictionary with key “hwp_sys” in the parameter
file (a toml file) assigned to the simulation:
sim = lbs.Simulation(
parameter_file=toml_filename,
random_seed=0,
)
The dictionary under the key “hwp_sys” will also contain the paths to the files from which the HWP
parameters are read in the multifrequency case (under the keys “band_filename/band_filename_solver”), or their values in the single frequency one. See the notebook hwp_sys/examples/simple_scan
for more details.
To define the bandpasses to use, you need to have a dictionary with key “bandpass”
(for \(\tau\)) or “bandpass_solver” (for \(\tau_s\)) under “hwp_sys”:
paramdict = {...
"hwp_sys": {...
"band_filename": path_to_HWP_param_file,
"band_filename_solver": path_to_HWP_solver_param_file,
"bandpass": {"band_type": "cheby",
"band_low_edge": band_low_edge,
"band_high_edge": band_high_edge,
"bandcenter_ghz": bandcenter_ghz,
"band_ripple_dB": ripple_dB_tod,
"band_order": args.order_tod},
"bandpass_solver": {...},
...}}
The above example is for a bandpass with Chebyshev filter, but there are other parameters to define
different bandpass profile. It is important to define the “band_type”, which can be “top-hat”,
“top-hat-exp”, “top-hat-cosine” and “cheby” (see the bandpass
module for more details)
and the band edges, which define the frequency range over which the bandpass transmission
is close or equal to 1. If not assigned, the “band_type” is automatically set to “top-hat”
and the band edges will correspond to the limits of the frequency array used (which, in the
hwp_sys
module, is read from the HWP parameter files). There are default values also for
the parameters defining the specific bandpass profiles (see the
hwp_sys/hwp_sys/bandpass_template_module
code).
There is also the possibility to read the bandpass profile from an external file, which has to be a .txt file with two columns, the frequency and the bandpass transmission. It is important that the frequency array used for “bandpass/bandpass_solver” coincides with the ones passed in the “band_filename/band_filename_solver” file. Here is how to pass the bandpass file:
paramdict = {...
"hwp_sys": {...
"band_filename": path_to_HWP_param_file,
"band_filename_solver": path_to_HWP_solver_param_file,
"bandpass": {"bandpass_file": path_to_bandpass_file },
"bandpass_solver": {"bandpass_file": path_to_bandpass_solver_file},
...}}
You can find more examples for the bandpass construction in the hwp_sys/examples/simple_scan
notebook.
API reference¶
HWP_sys¶
- litebird_sim.hwp_sys.hwp_sys.JonesToMueller(jones)¶
Converts a Jones matrix to a Mueller matrix. The Jones matrix is assumed to be a 2x2 complex matrix (np.array). The Mueller matrix is a 4x4 real matrix. Credits to Yuki Sakurai for the function.
- litebird_sim.hwp_sys.hwp_sys.get_mueller_from_jones(h1, h2, z1, z2, beta)¶
Converts the (frequency-dependent) input Jones matrix to a Mueller matrix. Returns Mueller matrix (3x3xNfreq), V-mode related terms are discarded, given the assumption of vanishing circular polarization.
Inputs: \(h_1\), \(h_2\), \(\zeta_1\), \(\zeta_2\), \(\beta\) (i.e. systematics of the HWP, not the full Jones matrix) Returns: \(M^{II}\), \(M^{QI}\), \(M^{UI}\), \(M^{IQ}\), \(M^{IU}\), \(M^{QQ}\), \(M^{UU}\), \(M^{UQ}\), \(M^{QU}\) (single/multi-frequency Mueller matrix terms)
- litebird_sim.hwp_sys.hwp_sys.compute_signal_for_one_sample(T, Q, U, mII, mQI, mUI, mIQ, mIU, mQQ, mUU, mUQ, mQU, c2ThPs, s2ThPs, c2PsXi, s2PsXi, c2ThXi, s2ThXi, c4Th, s4Th)¶
Bolometric equation, tod filling for a single (time) sample
- litebird_sim.hwp_sys.hwp_sys.compute_signal_for_one_detector(tod_det, mII, mQI, mUI, mIQ, mIU, mQQ, mUU, mUQ, mQU, pixel_ind, theta, psi, xi, maps)¶
Single-frequency case: compute the signal for a single detector, looping over (time) samples.
- litebird_sim.hwp_sys.hwp_sys.integrate_inband_signal_for_one_sample(T, Q, U, freqs, band, mII, mQI, mUI, mIQ, mIU, mQQ, mUU, mUQ, mQU, c2ThPs, s2ThPs, c2PsXi, s2PsXi, c2ThXi, s2ThXi, c4Th, s4Th)¶
Multi-frequency case: band integration with trapezoidal rule, \(\sum (f(i) + f(i+1)) \cdot (\nu_(i+1) - \nu_i)/2\) for a single (time) sample.
- litebird_sim.hwp_sys.hwp_sys.integrate_inband_signal_for_one_detector(tod_det, freqs, band, mII, mQI, mUI, mIQ, mIU, mQQ, mUU, mUQ, mQU, pixel_ind, theta, psi, xi, maps)¶
Multi-frequency case: band integration of the signal for a single detector, looping over (time) samples.
- litebird_sim.hwp_sys.hwp_sys.compute_TQUsolver_for_one_sample(mIIs, mQIs, mUIs, mIQs, mIUs, mQQs, mUUs, mUQs, mQUs, c2ThPs, s2ThPs, c2PsXi, s2PsXi, c2ThXi, s2ThXi, c4Th, s4Th)¶
Single-frequency case: computes \(A^T A\) and \(A^T d\) for a single detector, for one (time) sample.
- litebird_sim.hwp_sys.hwp_sys.compute_atd_ata_for_one_detector(atd, ata, tod, mIIs, mQIs, mUIs, mIQs, mIUs, mQQs, mUUs, mUQs, mQUs, pixel_ind, theta, psi, xi)¶
Single-frequency case: compute \(A^T A\) and \(A^T d\) for a single detector, looping over (time) samples.
- litebird_sim.hwp_sys.hwp_sys.integrate_inband_TQUsolver_for_one_sample(freqs, band, mIIs, mQIs, mUIs, mIQs, mIUs, mQQs, mUUs, mUQs, mQUs, c2ThPs, s2ThPs, c2PsXi, s2PsXi, c2ThXi, s2ThXi, c4Th, s4Th)¶
Multi-frequency case: band integration with trapezoidal rule, \(\sum (f(i) + f(i+1)) \cdot (\nu_(i+1) - \nu_i)/2\) for a single (time) sample.
- litebird_sim.hwp_sys.hwp_sys.integrate_inband_atd_ata_for_one_detector(atd, ata, tod, freqs, band, mIIs, mQIs, mUIs, mIQs, mIUs, mQQs, mUUs, mUQs, mQUs, pixel_ind, theta, psi, xi)¶
Multi-frequency case: band integration of \(A^T A\) and \(A^T d\) for a single detector, looping over (time) samples.
- class litebird_sim.hwp_sys.hwp_sys.HwpSys(simulation)¶
Bases:
object
A container object for handling tod filling in presence of hwp non-idealities following the approach of Giardiello et al. 2021 https://arxiv.org/abs/2106.08031 :param simulation: an instance of the class
Simulation
:type simulation:Simulation
- set_parameters(nside: int | None = None, mueller_or_jones: str = 'mueller', Mbsparams: MbsParameters | None = None, integrate_in_band: bool | None = None, built_map_on_the_fly: bool | None = None, correct_in_solver: bool | None = None, integrate_in_band_solver: bool | None = None, Channel: FreqChannelInfo | None = None, maps: None | ndarray = None, parallel: bool | None = None)¶
It sets the input paramters reading a dictionary sim.parameters with key “hwp_sys” and the following input arguments
- Parameters:
nside (integer) – nside used in the analysis
mueller_or_jones (str) – “mueller” or “jones” (case insensitive) it is the kind of HWP matrix to be injected as a starting point if ‘jones’ is chosen, the parameters \(h_1\), \(h_2\), \(\beta\), \(\zeta_1\), \(\zeta_2\) are used to build the Jones matrix \(\begin{pmatrix} 1 + h_1 & \zeta_1 \\ \zeta_2 & - (1 + h_2) e^{i \beta} \\ \end{pmatrix}\) and then converted to Mueller. \(\zeta_1\), \(\zeta_2\) are assumed to be complex \(h_1\), \(h_2\), \(\beta\) are assumed to be real \(\beta\) is assumed to be in degrees (later converted to radians. To reproduce the ideal HWP case, set all Jones parameters to 0. If Mueller parameters are used, set \(M^{II/QQ} = 1\), \(M^{UU} = -1\) and all the others to 0.
Mbsparams (
Mbs
) – an instance of theMbs
class Input maps needs to be in galactic (mbs default)integrate_in_band (bool) – performs the band integration for tod generation
built_map_on_the_fly (bool) – fills \(A^T A\) and \(A^T d\)
correct_in_solver (bool) – if the map is computed on the fly, fills \(A^T A\) using map-making (solver) HWP parameters
integrate_in_band_solver (bool) – performs the band integration for the map-making solver
Channel (
FreqChannelInfo
) – an instance of theFreqChannelInfo
classmaps (float) – input maps (3, npix) coherent with nside provided, Input maps needs to be in galactic (mbs default) if maps is not None, Mbsparams is ignored (i.e. input maps are not generated)
parallel (bool) – uses parallelization if set to True
- fill_tod(obs: Observation | List[Observation], hwp_radpsec: float, pointings: ndarray | List[ndarray] | None = None, save_tod: bool = False)¶
It fills tod and/or \(A^T A\) and \(A^T d\) for the “on the fly” map production
Args:
- obs (
Observation
): container for tod. If the tod is not required, you can avoid allocating
obs.tod
i.e. inlbs.Observation
useallocate_tod=False
.- pointings (float): pointing for each sample and detector
generated by
lbs.get_pointings()
Optional if already allocated inobs
. When generating pointing information, set the variablehwp
to None since the hwp rotation angle is added to the orientation angle within thefill_tod
function.
hwp_radpsec (float): hwp rotation speed in radiants per second
- save_tod (bool): if True,
tod
is saved inobs.tod
; if False, tod
gets deleted
- obs (
- make_map(obss)¶
It generates “on the fly” map. This option is only availabe if built_map_on_the_fly is set to True.
- Parameters:
class (obss list of) – Observations: only necessary for the communicator
pointings (float) – pointing for each sample and detector generated by lbs.get_pointings
hwp_radpsec (float) – hwp rotation speed in radiants per second
- Returns:
rebinned T,Q,U maps
- Return type:
map (float)
Bandpass template¶
- litebird_sim.hwp_sys.bandpass_template_module.lowpass_chebyshev(freqs, f0, order=3, ripple_dB=0.2)¶
Define a lowpass with chebyshev prototype
- Parameters:
freqs (np.array) – frequency in GHz
f0 (float) – low-frequency edge of the band in GHz
order (float) – chebyshev filter order (default: 3)
ripple_dB (float) – maximum ripple amplitude in decibels (default: 0.2)
- Returns:
The bandpass transmission
- litebird_sim.hwp_sys.bandpass_template_module.find_central_frequency(freqs, bandpass)¶
- Find the effective central frequency of
a bandpass profile as defined in https://arxiv.org/abs/1303.5070
- Parameters:
freqs (np.array) – frequency in GHz
bandpass (np.array) – transmission profile
- Returns:
central frequency
- litebird_sim.hwp_sys.bandpass_template_module.add_high_frequency_transmission(freqs, bandpass, location=3, transmission=0.5)¶
Add high frequency leakage
- Parameters:
freqs (np.array) – frequency in GHz
bandpass (np.array) – transmission profile
location (float) – multiple of the central frequency of the bandpass profile where add the leakage (default: 3)
transmission (float) – relative amplitude of the high frequency leakage with respect to the nominal band (default: 0.5)
- Returns:
Frequency range and bandpass transmission
- litebird_sim.hwp_sys.bandpass_template_module.beam_throughtput(freqs)¶
Beam throughtput factor
- Parameters:
freqs (np.array) – frequency in GHz
- Returns:
\(1/\nu^2\)
- litebird_sim.hwp_sys.bandpass_template_module.bandpass_profile(freqs: array | None = None, bandpass: dict | None = None, include_beam_throughput: bool | None = None)¶
Returns the bandpass profile either reading it from a file or computing it through the class BandPassInfo
- Parameters:
freqs (np.array) – the array of frequencies considered
bandpass (dict) – dictionary with either the key “bandpass_file” pointing to a txt file with a column of freqs and a column of bandpass profile (attention: the frequencies read from the file have to coincide with the input array freqs) or with a key “band_type” (same naming convention as bandtype in BandPassInfo). In the second case, you should also define the edges of the bandpass, band_high_edge and band_low_edge, that would otherwise be set to the lowest and highest value of freqs.
include_beam_throughput (bool) – if True, divides by \(\nu^2\)
- Returns:
Array of frequency range (same as freqs) and bandpass profile