Goal is to be able to both generate sky maps and to source subtract them across multiple frequency bands. Included code does just that. Run help_yubo for help with arguments of various functions.
Documentation/coding style will probably be revised as merging into general kSZ simulation
The flow goes roughly as follows:
- Generate a PSD using
spdensgen
- Generate a noise realisation using
gennoise
- Compute a fit using
submatwrap
in the single band case and eithersubmatwrap_multi
in the multi band fixed-SED case oremissivity_multi
/bbody_multi
in the beta-T_dust case (I never implementedsubmatwrap
for either of those two since I focused on trying to break the degenerate chi2 first)
Analysis techniques include:
- Run best-estimator routine multiple times and create multi-dimensional histogram, take a Gaussian fit for covariance matrix
- For a fixed noise realization, compute the chi2 surface as a function of the various parameters, interpret as a likelihood function to obtain a covariance matrix
- Analytically compute, given a signal map, spectral density etc., what the expected covariance matrix should look like.
The covariance matricies for all three of these should agree assuming a quadratic chi2 surface, which does not hold for beta-T_dust simultaneously unconstrained
Some examples of these routines in action can be found under my Bitbucket
help_yubo.pro
helper
to2d.pro
sfit_p.pro
invert_matuncert.pro
getcovars.pro
hist_wrap.pro
contour_wrap.pro
single_band
helper
subtractmat.pro
meansub.pro
estopt.pro
iterimp.pro
subtractmax.pro
addfunc.pro
submatwrap.pro
spdensgen.pro
fitquad.pro
addgauss.pro
gennoise.pro
multi_band
helper
iterimp_multi.pro
estopt_multi.pro
tnmin_wrap.pro
subtractmat_multi.pro
amps_multi.pro
sigm_multi.pro
spdensgen_multi.pro
addfunc_multi.pro
getsedpts.pro
subtractknown_multi.pro
bbody_multi.pro
subtractmax_multi.pro
gennoise_multi.pro
emissivity_multi.pro
addgauss_multi.pro
submatwrap_multi.pro
plot_bbody_snu.pro
help_yubo.pro
(procedure)- Prints out a short documentation on most of the routines that I use. Slightly outdated now, this readme might be more useful.
helper
to2d.pro
- Short little routine to convert a 1D array index to a 2D array index
sfit_p.pro
- Canned IDL
sfit.pro
routine. Copied tosfit_p
so would not override built-insfit
, but do not remember why I made a copy. I think it was anIDL_PATH
conflict
- Canned IDL
invert_matuncert.pro
- Inverts a 2x2 matrix that has uncertainties in its elements. The return also has uncertainties attached.
getcovars.pro
- For a fit to the chi2 surface, computes the covariance matrix obtained by treating the chi2 as a likelihood function
hist_wrap.pro
- A wrapper around the build in idl
hist
routine that has some nice fonting and colors and a useful options concerning what values to write where. Automatically computes placement of text for both logarithmic and non-logarithmic scales.
- A wrapper around the build in idl
contour_wrap.pro
- A wrapper around the build in idl
contour
routine that has some nice fonting and colors and a useful options concerning where to draw the contours.
- A wrapper around the build in idl
single_band
Routines for the single frequency band case. Most of these are described in 2D.pdf.helper
subtractmat.pro
- Given the X, Y positions of some Gaussians, find the set of amplitudes of each Gaussian that minimize the chi2. Also known as simultaneously fitting for Gaussian amplitudes given fixed positions.
meansub.pro
- This is an artifact of an old hypothesis. With
MINSEP=0
, it is possible for two nearly-coincident sources to look like a single source.meansub
computes a list of sources by runningsubmatwrap
(the best subtraction routine I have, discussed later), then examines each individual source to see whether the mean of the residual map decreases if we add another source. - This is a faulty hypothesis, as is discussed in my thesis, since the mean of the residual map is a meaningless number.
- This is an artifact of an old hypothesis. With
estopt.pro
- This is in charge of the "iterative improvement" procedure described above. Given a set of flux estimates, it iteratively adds and re-subtracts them until positions converge
iterimp.pro
- Calls
estopt
but also computes initial flux estimates upon which to improve
- Calls
subtractmax.pro
- For a given 1D signal, spectral density and binwidth, find the maximum signal for a Gaussian with width
sigm
. Can constrainreal_pos
if the position of the Gaussian is to be fixed (useful for certain exercises)
- For a given 1D signal, spectral density and binwidth, find the maximum signal for a Gaussian with width
addfunc.pro
- Add a bunch of Gaussians of a fixed height to a
noise
array. Can specify aMINSEP
to make sure the Gaussians do not overlap too much
- Add a bunch of Gaussians of a fixed height to a
submatwrap.pro
- Feeds iterimp (which iterates until position estimators converge) into subtractmat (which computes optimal flux estimators given positions.
spdensgen.pro
- Given some parameters, computes a 1D spectral power density
fitquad.pro
- Fits a 4x4 patch to a quadratic. Fits the log of a Gaussian.
addgauss.pro
- Add a gaussian of width
sigm
and heighta
to a given signal array
- Add a gaussian of width
gennoise.pro
- Given a spectral density and a binwidth, generate a random noise realisation
multi_band
multi-band versions of the above routines. I will only elaborate on those routines that don't exist in the single-band versions. The formalism is developed mostly only in my thesis, as by this time I had grown very bad at maintaining 2D.pdfhelper
iterimp_multi.pro
estopt_multi.pro
tnmin_wrap.pro
- Since with multiple frequency bands we must also fit for beta, T_dust, we must use a gradient descent fitter. IDL has a built in
tnmin
, andtnmin_wrap
is a wrapper that specializestnmin
for our specific application. Supports both beta-only fitting or beta+T_dust (/bbody
flag)
- Since with multiple frequency bands we must also fit for beta, T_dust, we must use a gradient descent fitter. IDL has a built in
subtractmat_multi.pro
amps_multi.pro
- Given a beta or a beta+T_dust, compute amplitudes relative to lowest frequency band.
sigm_multi.pro
- Compute beam widening relative to lowest frequency band due to diffraction and aperture size.
spdensgen_multi.pro
addfunc_multi.pro
getsedpts.pro
- computes the single-band estimators + uncertainties for a multi-band signal. Used to compare the multi-band flux estimators to the single-band estimators.
subtractknown_multi.pro
- Given a known SED, subtract it from the map. Used to subtract the true parameters and compare the chi2 to that obtained from subtracting the best fit parameters
bbody_multi.pro
- calls
subtractmax_multi
and performs gradient descent to find the best-fit beta/T_dust. Can fix either parameter or let both float (bad degeneracies...)
- calls
subtractmax_multi.pro
- Should be noted this only works for a fixed beta/T_dust
gennoise_multi.pro
emissivity_multi.pro
- Same as bbody_multi but only performs gradient descent on beta
addgauss_multi.pro
submatwrap_multi.pro
plot_bbody_snu.pro
- Plots SEDs for a list of emissitivities and T_dusts. Useful to compare two SEDs that have similar chi2 values
knose
b/c used to generate noise in k-space, k-noisesigm
b/csigma
is a reserved word in IDL