Foreground transfer function

docs/notebooks/XFaster_Tutorial.ipynb

@@ -404,7 +404,7 @@
     "\n",
     "Another option you might wish to use is `weighted_bins`, which changes the default $\\chi_b(\\ell)$ binning operator from a tophat to one that weights by $\\ell(\\ell+1)$.\n",
     "\n",
-    "To enable fitting for a foreground component in the harmonic domain, set `foreground_fit=True`, and use the corresponding `bin_width_fg` to set the bin width for the foreground amplitude bins, and set `lmin_fg` and `lmax_fg` to optionally limit the bandwidth over which foregrounds are to be fit.  Optionally, also set `beta_fit=True` to enable fitting for a frequency spectral index component for the foreground amplitude; otherwise, the frequency dependence is assumed to follow a single component dust model as defined in [scale_dust()](../api.rst#xfaster.spec_tools.scale_dust).  If `foreground_fit` is True, the foreground model is assumed to be a simple power law model, as defined in [dust_model()](../api.rst#xfaster.spec_tools.dust_model), with the same transfer function as the CMB component."
+    "To enable fitting for a foreground component in the harmonic domain, set `foreground_fit=True`, and use the corresponding `bin_width_fg` to set the bin width for the foreground amplitude bins, and set `lmin_fg` and `lmax_fg` to optionally limit the bandwidth over which foregrounds are to be fit.  Optionally, also set `beta_fit=True` to enable fitting for a frequency spectral index component for the foreground amplitude; otherwise, the frequency dependence is assumed to follow a single component dust model as defined in [scale_dust()](../api.rst#xfaster.spec_tools.scale_dust).  If `foreground_fit` is True, use the `foreground_type` and `foreground_transfer_type` file options to enable the use of a particular foreground signal model and computation of a separate transfer function for the foreground component.  If the file options are not set, the foreground model is assumed to be a simple power law model, as defined in [dust_model()](../api.rst#xfaster.spec_tools.dust_model), with the same transfer function as the CMB component."
    ]
   },
   {
@@ -541,7 +541,7 @@
    "cell_type": "markdown",
    "metadata": {},
    "source": [
-    "Now we want to get the ensemble average of all of our signal and noise simulations, which we're going to use to calculate the filter transfer function and the noise shape, respectively. This is done with [get_masked_sims()](../api.rst#xfaster.xfaster_class.XFaster.get_masked_sims). The method will also compute the signal cross noise terms, which are used for null tests, where they can contribute significantly to the expected residuals that are subtracted from the data."
+    "Now we want to get the ensemble average of all of our signal (and/or optionally foreground) and noise simulations, which we're going to use to calculate the filter transfer function and the noise shape, respectively. This is done with [get_masked_sims()](../api.rst#xfaster.xfaster_class.XFaster.get_masked_sims). The method will also compute the signal cross noise terms, which are used for null tests, where they can contribute significantly to the expected residuals that are subtracted from the data."
    ]
   },
   {
@@ -574,7 +574,8 @@
    "source": [
     "The resulting outputs are:\n",
     "\n",
-    "* `cls_signal`: the average of the signal-only cross spectra, used to compute the transfer function\n",
+    "* `cls_signal`: the average of the signal-only cross spectra, used to compute the signal transfer function\n",
+    "* `cls_fg`: the average of the foreground-only cross spectra, used to compute the foreground transfer function\n",
     "* `cls_noise`: the average of the noise-only cross spectra, used as the noise model, $N_\\ell$\n",
     "* `cls_sim`: the average of the signal+noise spectra, where signal and noise maps are added in $a_{\\ell m}$s and thus the spectra include SxN terms.\n",
     "* `cls_med`: the median of the signal+noise spectra-- this is mainly used for debugging potential biases seen in the pipeline\n",
@@ -719,9 +720,9 @@
     "\n",
     "We have $K_{\\ell, \\ell'}$ and $B_{\\ell}$. For the transfer function calculation, we're going to set $F_\\ell$ to 1 so that we measure $q_b$s as the deviation from a uniform transfer function for our simulations. Binning, $\\chi_b$ has been chosen. All that's left is the full sky signal shape, $\\mathcal{C}_{\\ell'}^{XY (S)}$, loaded with [get_signal_shape()](../api.rst#xfaster.xfaster_class.XFaster.get_signal_shape). \n",
     "\n",
-    "For calculating the transfer function, this is just the shape spectrum that went into making our simulations. This spectrum can be specified by setting the [xfaster_run()](../api.rst#xfaster.xfaster_exec.xfaster_run) argument `signal_transfer_spec` to a file containing the spectrum for the signal component. If not provided, the code will look in the maps directory for the signal sims for a file labeled `spec_signal_<signal_type>.dat`. The file is expected to look like a CAMB output file, as demonstrated in `make_example_maps.py`, which writes such a file to the proper location in the signal sims directory.\n",
+    "For calculating the transfer function, this is just the shape spectrum that went into making our simulations. This spectrum can be specified by setting the [xfaster_run()](../api.rst#xfaster.xfaster_exec.xfaster_run) argument `signal_transfer_spec` or `foreground_transfer_spec` to a file containing the spectrum for the signal or foreground component, respectively. If not provided, the code will look in the maps directory for the signal sims for a file labeled `spec_signal_<signal_type>.dat`, and similarly for foreground sims. The file is expected to look like a CAMB output file, as demonstrated in `make_example_maps.py`, which writes such a file to the proper location in the signal and foreground sims directories.\n",
     "\n",
-    "For foreground fitting, this function also applies a frequency-dependent scaling to the input signal shape, using the [scale_dust()](../api.rst#xfaster.spec_tools.scale_dust) function.  The arguments `freq_ref` and `beta_ref` are used to set the reference frequency and spectral index for the scaling function."
+    "For foreground fitting, this function also applies a frequency-dependent scaling to the input signal shape, using the [scale_dust()](../api.rst#xfaster.spec_tools.scale_dust) function.  The arguments `freq_ref` and `beta_ref` are used to set the reference frequency and spectral index for the scaling function.  Different reference points may be used for the transfer function computation, by setting the arguments `freq_ref_transfer` and `beta_ref_transfer`."
    ]
   },
   {
@@ -781,13 +782,13 @@
     "The expression for the Fisher matrix does not change, other than the fact that its constituents are the same as detailed above for the transfer function.\n",
     "\n",
     "\n",
-    "Within the code, `get_transfer` basically has two steps within the function itself, which it performs per map. \n",
+    "Within the code, `get_transfer` basically has two steps within the function itself, which it performs per map and per signal component (either `\"cmb\"` or `\"fg\"`). \n",
     "\n",
     "1. Load up the $\\tilde{\\mathcal{C}}_{b\\ell}$: \n",
     "```python \n",
     "cbl = self.bin_cl_template(map_tag=m0, transfer=True)\n",
     "```\n",
-    "This uses the signal_shape internally that we calculated earlier, `m0` is the map, which is used to select the beam and kernel, and `transfer=True` sets the $F_\\ell$ term to 1.\n",
+    "This uses the signal_shape internally that we calculated earlier, `m0` is the map, which is used to select the beam and kernel, and `transfer=True` sets the $F_\\ell$ term for the CMB component to 1.\n",
     "\n",
     "2. Run [fisher_iterate()](../api.rst#xfaster.xfaster_class.XFaster.fisher_iterate).\n",
     "```python\n",
@@ -804,7 +805,9 @@
     "```python\n",
     "qb['cmb_eb'] = np.sqrt(np.abs(qb['cmb_ee'] * qb['cmb_bb']))\n",
     "qb['cmb_tb'] = np.sqrt(np.abs(qb['cmb_tt'] * qb['cmb_bb']))\n",
-    "```"
+    "```\n",
+    "\n",
+    "When computing the transfer function for the foreground component, the TE transfer function is also handled similarly to the EB and TB terms."
    ]
   },
   {

example/xfaster_example_fg.py

@@ -34,6 +34,7 @@
     "noise_type": "gaussian",
     "mask_type": "rectangle",
     "signal_type": "synfast",
+    "foreground_type": "gaussian",
     "data_root2": None,
     "data_subset2": None,
     # residual fitting
@@ -50,6 +51,7 @@
     # spectrum
     "ensemble_mean": False,
     "tbeb": True,
+    "fix_fg_transfer": True,
     "converge_criteria": 0.005,
     "iter_max": 200,
     "save_iters": False,

xfaster/parse_tools.py

@@ -315,6 +315,9 @@ def load_and_parse(filename, check_version=True):
             ]:
                 data.pop(k, None)
 
+        # In version 1, these parameters referred to ensembles used for
+        # sim_index runs.  In version 2+, these parameters refer to ensembles
+        # used for the foreground component of the covariance model
         if "foreground_type" in data:
             data["foreground_type_sim"] = data.pop("foreground_type")
         if "foreground_root" in data: