Differences

This shows you the differences between two versions of the page.

--- calibration_simulations [2018/05/21 08:00] – mauriziotomasi
+++ calibration_simulations [2019/02/20 04:56] (current) – Put an "outdated" warning message at the top of the page. mauriziotomasi
@@ Line 1: / Line 1: @@
+**This page contains outdated information**. Refer to the report «Simulating dipole calibration for PICO», by Maurizio Tomasi, for the most up-to-date information {{::tomasi-calibration-note.pdf|}}.
 ====== Calibration simulations ======
-In this page we collect information and results about the simulation of the photometric calibration procedure.
+In this page we collect information and results about the simulation
+of the photometric calibration procedure.
 ===== Hypothesis =====
-We assume to calibrate the instrument using the solar dipole of the CMB, whose peak-to-peak amplitude is ~7 mK. We use [[https://github.com/hpc4cmb/toast|TOAST]] to produce the timelines, and [[https://github.com/ziotom78/dacapo_calibration|DaCapo]] to compute the calibration factors. The scripts used to run the simulation have been prepared by Andrea Zonca and are available at [[https://github.com/zonca/pico-simulations|pico-simulations]].
+I assume to calibrate the instrument using the solar dipole of the
+CMB, whose peak-to-peak amplitude is ~7 mK. I use
+[[https://github.com/hpc4cmb/toast|TOAST]] to produce the timelines,
+and [[https://github.com/ziotom78/dacapo_calibration|DaCapo]] to
+compute the calibration factors. The scripts used to run the
+simulation have been prepared by Andrea Zonca and are available at
+[[https://github.com/zonca/pico-simulations|pico-simulations]].
-As a reference, the dipole in Ecliptic coordinates has the following shape:
+As a reference, the dipole in Ecliptic coordinates has the following
+shape:
 {{ ::dipole-ecliptic.png?600 |}}
@@ Line 13: / Line 23: @@
 ===== Observation of the dipole =====
-I used the scanning strategy parameters listed at the page [[optimizingscanstrategy:start|Optimizing Scan Strategy]], with **no** kinematic dipole and the Galactic signal simulated by PySM for a W-band detector, which should be the best case in terms of Galactic contamination.
+I used the scanning strategy parameters listed at the page
+[[optimizingscanstrategy:start|Optimizing Scan Strategy]], with **no**
+kinematic dipole and the Galactic signal simulated by PySM for a
+W-band detector, which should be the best case in terms of Galactic
+contamination.
-The sky coverage as a function of time shows a sharp rise in the first few hours, reaching ~50% of the sky in one precession period. Then, the value increases slowly till it reaches 100% after nearly 6 months. (For comparison, fsky varies more linearly in the case of Planck.) The following figure shows fsky as a function of time:
+The sky coverage as a function of time shows a sharp rise in the first
+few hours, reaching ~50% of the sky in one precession period. Then,
+the value increases slowly till it reaches 100% after nearly 6
+months. (For comparison, fsky varies more linearly in the case of
+Planck.) The following figure shows fsky as a function of time:
 {{ ::pico-fsky-6months.png?600 |}}
-The precession allows PICO to observe wider portions of the sky in the same time. This is good for dipole calibration as well, as it allows the detectors to potentially sample larger variations of the dipole, thus increasing the S/N of the dipole measurement.
+The precession allows PICO to observe wider portions of the sky in the
+same time. This is good for dipole calibration as well, as it allows
+the detectors to potentially sample larger variations of the dipole,
+thus increasing the S/N of the dipole measurement.
-A simple parameter that encapsulates much of the details about photometric calibration is the peak-to-peak variation in the dipole signal during some fixed time span. As PICO detectors are not meant to perform absolute measurements, the //modulation// of the dipole signal consititutes our calibration source. The peak-to-peak amplitude of the dipole is therefore a good tracer of the S/N for the calibration signal.
+A simple parameter that encapsulates much of the details about
+photometric calibration is the peak-to-peak variation in the dipole
+signal during some fixed time span. As PICO detectors are not meant to
+perform absolute measurements, the //modulation// of the dipole signal
+consititutes our calibration source. The peak-to-peak amplitude of the
+dipole is therefore a good tracer of the S/N for the calibration
+signal.
-Consider this image, which shows the simulated timeline produced by TOAST for one of the W-band detectors (temperature in K as a function of time in seconds, no Galaxy is present):
+Consider this image, which shows the simulated timeline produced by
+TOAST for one of the W-band detectors (temperature in K as a function
+of time in seconds, no Galaxy is present):
 {{ :pico-nofg-totaltod.png?600 |}}
-The low-frequency modulation is due to the dipole signal, while the high-frequency noise is detector noise. As the plot shows ~100 seconds of data, the peak-to-peak variation of the dipole signal within this time frame is less than 3 mK. The following image shows how the peak-to-peak dipole amplitude changes as a function of the time window used to estimate the peak-to-peak value:
+The low-frequency modulation is due to the dipole signal, while the
+high-frequency noise is detector noise. As the plot shows ~100 seconds
+of data, the peak-to-peak variation of the dipole signal within this
+time frame is less than 3 mK. The following image shows how the
+peak-to-peak dipole amplitude changes as a function of the time window
+used to estimate the peak-to-peak value:
 {{ ::dipole-amplitude.gif?600 |}}
@@ Line 33: / Line 67: @@
 (Go to https://zzz.physics.umn.edu/ipsig/_media/dipole-amplitude.gif to see the animation.)
-There are several interesting facts that become apparent from the image above. First of all, there are high-frequency fluctuations which disappear for longer time windows. Specifically, as the length of the time window approaches the precession period, high-frequency fluctuations are present no longer. A second interesting fact is that the plot shows low-frequency fluctuations as well. These are unavoidable and are due to the yearly variation in the angle between the spacecraft spin axis and the CMB dipole axis. As a matter of fact, these fluctuations are present even for much longer time windows:
+There are several interesting facts that become apparent from the
+image above. First of all, there are high-frequency fluctuations which
+disappear for longer time windows. Specifically, as the length of the
+time window approaches the precession period, high-frequency
+fluctuations are present no longer. A second interesting fact is that
+the plot shows low-frequency fluctuations as well. These are
+unavoidable and are due to the yearly variation in the angle between
+the spacecraft spin axis and the CMB dipole axis. As a matter of fact,
+these fluctuations are present even for much longer time windows:
 {{ ::dipole-amplitude-overview.png?600 |}}
@@ Line 40: / Line 82: @@
 ===== Running DaCapo =====
-We assume that the calibration code to be used in a real mission would include a component-separation step able to reduce the bias due to the contamination of the dipole signal by Galactic large-scale structures. Therefore, in our simulations we assume that the only signal in the sky apart from the dipole is the CMB itself.
+I assume that the calibration code to be used in a real mission would
+include a component-separation step able to reduce the bias due to the
+contamination of the dipole signal by Galactic large-scale
+structures. Therefore, in our simulations I assume that the only
+signal in the sky apart from the dipole is the CMB itself. From the
+point of view of calibration, CMB acts as a noise signal which is
+correlated among detectors.
-Our simulations assume that each detector has a calibration factor equal to 1.0. The calibration code is therefore expected to retrieve the same value from the simulated data, plus some statistical noise.
+The DaCapo algorithm is explained in
+[[https://www.aanda.org/articles/aa/abs/2016/10/aa26632-15/aa26632-15.html|Planck
+results. V. LFI calibration]] and in the source code of the
+Python version of
+[[https://github.com/ziotom78/dacapo_calibration/raw/master/dacapo_calibration.pdf|DaCapo]],
+which has been written with the aim to be readable as an article. The
+algorithm works as follows:
-I have run a few simulations using DaCapo on the simulated timelines produced by Andrea Zonca. Here are the details:
+  - Apply an extended version of the destriping equation, which reconstructs the baselines of 1/f noise as well as the gain factors. This uses the following assumptions:
+    * 1/f baselines have zero mean
+    * Gains are constrained by a model of the dipole, which must be provided as input
+  - The equation is solved iteratively using the Conjugate Gradient algorithm. At the end, the following estimates are available:
+    * 1/f baselines
+    * Timeline of gain factors
+    * Sky map (temperature only) without the dipole
+  - Iterate the algorithm again, subtracting the sky map obtained in the previous step from the timelines, in order to remove any systematic due to spurious signals (CMB, foregrounds).
+  - When the 1/f baselines, the gain factors, and the sky map no longer change, save all the results and quit.
+  - DaCapo must be run once per each detector, so that the sky map produced as output is typically noisier than the maps expected from a typical configuration of detectors.
-  * 2 years of observation
+I employed the DaCapo algorithm in the following pipeline:
+  - Use TOAST to produce timelines using two detectors (0A and 0B) along the boresight direction, the nominal scanning strategy and realistic noise, and save them to FITS files. The gain used in this step is 1.0 (constant) for both detectors. TOAST produces a sky map (including the dipole) using the Madam map-maker (a destriper) and data from both detectors 0A and 0B.
+  - Run DaCapo on the FITS files produced during the previous stage. As the algorithm works on single detectors, the code must be run twice (once for detector 0A and once for 0B).
+  - If DaCapo converges, re-run TOAST as in the previous step, but this time use the gains calculated by DaCapo instead of assuming a constant (1.0) gain.
+  - At the end of the simulation, take the difference between the maps produced by TOAST in steps 1. and 3. and compute the power spectrum. These spectra represent the error caused by an imperfect estimate of gain drifts during the mission.
+Our simulations assume that each detector has a calibration factor
+equal to 1.0. The calibration code is therefore expected to retrieve
+the same value from the simulated data, plus some statistical
+noise. Of course, the code does not know that the gain is constant, so
+it always tries to estimate the dependency of gain on time.
+I have run a few simulations using DaCapo on the simulated timelines
+produced by Andrea Zonca. Here are the details:
+  * 1 year of observation
   * 2 boresight W-band detectors
   * No Galaxy
@@ Line 54: / Line 133: @@
   * Constant gain (1.0)
-As DaCapo implements a destriper as well, let's check that the 1/f is correctly estimated. The baseline I used was 10 s:
+As DaCapo implements a destriper as well, let's check that the 1/f is
+correctly estimated. The baseline I used was 10 s:
-{{ ::dacapo-cleaned-tod-dsp.png?600 |}}
+{{ :dacapo-cleaned-tod-dsp.png?600 |}}
-This plot is complicate to understand due to the presence of the sky signal (the dipole peak and its harmonics). Let's compute the PSD of the noise component alone: this is easy, since TOAST saves the TOD of the foreground signal in a separate column in the FITS files.
+(“Original” is the TOD produced by TOAST and fed to DaCapo, while
+“Cleaned” is the same TOD minus the 1/f baselines calculated by
+DaCapo.) This plot is complicate to understand due to the presence of
+the sky signal (the dipole peak and its harmonics). Let's compute the
+PSD of the noise component alone: this is easy, since TOAST saves the
+TOD of the foreground signal in a separate column in the FITS files.
 {{ ::dacapo-cleaned-tod-dsp-nofg.png?600 |}}
-The peak at 167 mHz corresponds to the spin period (60 s): note that the height of this residual peak matches the level of the high-frequency WN tail. Apart from that, 1/f is correctly suppressed by the code.
+The peak at 167 mHz corresponds to the spin period (60 s): note that
+the height of this residual peak matches the level of the
+high-frequency WN tail. Apart from that, 1/f is correctly suppressed
+by the code.
-The estimated gains for a time window of 10 hours (the same as the precession period) are shown here:
+The estimated gains for a time window of 10 hours (the same as the
+precession period) are shown here:
 {{ ::dacapo-gains-nofg-10h.png?600 |}}
-You can see that the RMS of the gain factors varies with time. This is anticorrelated with the peak-to-peak amplitude of the dipole (see the plot above). The average value of the gains shows some residual bias, as it is not exactly 1.0. The RMS of the timeline is 0.04%, which can be considered as an estimate of the relative calibration error.
+You can see that the RMS of the gain factors varies with time. This is
+anticorrelated with the peak-to-peak amplitude of the dipole (see the
+plot above). The average value of the gains shows some residual bias,
+as it is not exactly 1.0: this bias is due to the fact that TOAST
+convolves the sky and the dipole signals with a beam of some finite
+size, but DaCapo considers a pencil beam (i.e., a Dirac delta). The
+usual way to deal with this is to use a 4π convolver, like the one
+used in Appendix A of
+[[https://www.aanda.org/articles/aa/abs/2016/10/aa26632-15/aa26632-15.html|Planck
+results. V. LFI calibration]], but for our purpose it's easier to
+subtract an offset so that the average value of the gain is 1.0. The
+RMS of the timeline is 0.04%, which can be considered as an estimate
+of the relative calibration error.
+===== Output maps =====
+As sketched in the previous section, I produced full-sky maps with
+TOAST, calibrate them with DaCapo and re-run TOAST with the gains
+estimated by DaCapo to get a new full-sky map which contains the
+effect of imperfect gain drift reconstruction. The difference between
+the two maps show the projected effect of gain drifts on the sky.
+I have run two sets of simulations, using the following assumptions:
+  - Two noiseless receivers (0A and 0B);
+  - Two W-band receivers with nominal white noise and 1/f noise with knee frequency 10 mHz.
+The differenced maps produced in the two cases are quite different:
+{{ :systematicswg:cal_error_maps.png?600 |}}
+Of course, the scale of the effect considers the presence of two
+detectors on the whole focal plane, which is not the case for PICO.
+Computing the BB spectrum of the map in the two cases and confronting
+it with the input map shows the scale of the effect (as well as two
+power spectra assuming only tensor modes, provided as a reference):
+{{ :systematicswg:calibration_spectrum_bb.png?600 |}}
+In order to produce this plot, I rescaled the spectrum by 2/N, where
+N is the number of W-band bolometers in the focal plane of PICO, and
+then I performed a NET-weighted average over all the cosmology bands
+from 60 GHz (band 7) to 220 GHz (band 22): this
+assumes that the gain fluctuations are negligibly correlated among
+detectors. This is the case for the case with realistic noise, but not
+for the noiseless simulation; this should be expected, as the only
+noise here is caused by the CMB, which is the same for detectors 0A
+and 0B:
+{{ :systematicswg:gain_timeline.png?600 |}}

User Tools

Differences

Page Tools