dmsky package

Module contents

Dark matter skymaps.

Priors on J-factor

class dmsky.priors.PriorFunctor(funcname, scale=1.0)

Bases: object

A functor class that wraps simple functions we use to make priors on parameters.

funcname

A string identifying the function.

log_value(x)

Return the log of the function value

Parameters:x (numpy.ndarray) – Input values Returns:y – Output values, same shape as x Return type:numpy.ndarray

marginalization_bins()

Binning to use to do the marginalization integrals

Default is to marginalize over two decades, centered on mean, using 1000 bins

mean()

Mean value of the function.

normalization()

Normalization, i.e. the integral of the function over the normalization_range.

profile_bins()

The binning to use to do the profile fitting

Default is to profile over +-5 sigma, Centered on mean, using 100 bins

scale

The scale factor applied to input values

sigma()

The ‘width’ of the function. What this means depend on the function being used.

class dmsky.priors.FunctionPrior(funcname, mu, sigma, fn, lnfn=None, scale=1.0)

Bases: dmsky.priors.PriorFunctor

Implementation of a prior that simply wraps an existing function

log_value(x)

“Return the log of the function value

Parameters:x (numpy.ndarray) – Input values Returns:y – Output values, same shape as x Return type:numpy.ndarray

mean()

Return the mean value of the function.

normalization()

The normalization i.e., the intergral of the function over the normalization_range

sigma()

Return the ‘width’ of the function. What this means depend on the function being used.

class dmsky.priors.GaussPrior(mu, sigma, scale=1.0)

Bases: dmsky.priors.FunctionPrior

Implemenation of a Prior that wraps a Gaussian

class dmsky.priors.LGaussPrior(mu, sigma, scale=1.0)

Bases: dmsky.priors.FunctionPrior

Implemenation of a Prior that wraps a log Gaussian

class dmsky.priors.LGaussLogPrior(mu, sigma, scale=1.0)

Bases: dmsky.priors.FunctionPrior

Implemenation of a Prior that wraps the inverse of the log of a Gaussian (i.e., x and y axes are swapped) The prior is implemented in log-space.

class dmsky.priors.LGaussLikePrior(mu, sigma, scale=1.0)

Bases: dmsky.priors.FunctionPrior

Implemenation of a Prior that wraps the inverse of the log of a Gaussian (i.e., x and y axes are swapped)

class dmsky.priors.LognormPrior(mu, sigma, scale=1.0)

Bases: dmsky.priors.PriorFunctor

A wrapper around the lognormal function.

A note on the highly confusing scipy.stats.lognorm function… The three inputs to this function are:

s : This is the variance of the underlying gaussian distribution

scale = 1.0 : This is the mean of the linear-space lognormal distribution. The mean of the underlying normal distribution occurs at ln(scale)

loc = 0 : This linearly shifts the distribution in x (DO NOT USE)

The convention is different for numpy.random.lognormal

mean : This is the mean of the underlying normal distribution (so mean = log(scale))

sigma : This is the standard deviation of the underlying normal distribution (so sigma = s)

For random sampling: numpy.random.lognormal(mean, sigma, size)

mean : This is the mean of the underlying normal distribution (so mean = exp(scale))

sigma : This is the standard deviation of the underlying normal distribution (so sigma = s)

scipy.stats.lognorm.rvs(s, scale, loc, size)

s : This is the standard deviation of the underlying normal distribution

scale : This is the mean of the generated random sample scale = exp(mean)

Remember, pdf in log space is plot( log(x), stats.lognorm(sigma,scale=exp(mean)).pdf(x)*x )

Parameters:

• mu (float) – Mean value of the function
• sigma (float) – Variance of the underlying gaussian distribution

mean()

Mean value of the function.

normalization()

Normalization, i.e. the integral of the function over the normalization_range.

sigma()

The ‘width’ of the function. What this means depend on the function being used.

class dmsky.priors.FileFuncPrior(filename, scale=1.0)

Bases: dmsky.priors.PriorFunctor

A wrapper around the interpolated function.

Parameters:filename (string) – File with the function parameters

mean()

Mean value of the function.

sigma()

The ‘width’ of the function. What this means depend on the function being used.

Dark Matter Denisty Porfiles

class dmsky.density.DensityProfile(**kwargs)

Bases: pymodeler.model.Model

Am abstract base class for DM density profiles

At a minimum sub-classes need to implement the self._rho(r) method to compute the density as a function of radius from the center of the halo

deriv_params

Return the list of paramters we can take derivatives w.r.t.

mass(r=None)

Compute the mass of the object out to a particular radius.

Parameters:r (numpy.array or float) – The radii Returns:values – Return values, same shape as the input radii Return type:numpy.array

rho(r)

Return the density for given radii.

Parameters:r (numpy.array or float) – The radii Returns:values – Return values, same shape as the input radii Return type:numpy.array

rho_deriv(r, paramNames)

Return the derivatives of the density as a function of radius,

w.r.t. a list of parameters

Parameters:

• r (numpy.array or float) – The radii
• paramNames (list) – The names of the parameters to differentiation w.r.t.

Returns:

matrix – An n x m array, where: n is the number of radii m is the number of parameters

Return type:

numpy.array

rho_uncertainty(r)

Calculate the uncertainty of the density at given radii

Parameters:r (numpy.array or float) – The radii Returns:values – Return values, same shape as the input radii Return type:numpy.array

set_mvir_c(mvir, c)

Fix the mass inside the virial radius.

Parameters:

• mvir (float) – The virial radius
• c (float) – Scale factor

set_rho_r(rho, r)

Fix the density normalization at a given radius.

Parameters:

• rho (float) – The normalization density
• r (float) – The corresponding radius

class dmsky.density.UniformProfile(**kwargs)

Bases: dmsky.density.DensityProfile

Uniform spherical profile rho(r) = rhos for r <= rs rho(r) = 0 otherwise

class dmsky.density.IsothermalProfile(**kwargs)

Bases: dmsky.density.DensityProfile

Non-Singular Isothermal Profile: Begeman et al. MNRAS 249, 523 (1991) http://adsabs.harvard.edu/full/1991MNRAS.249..523B rho(r) = rhos/(1+(r/rs))**2

class dmsky.density.BurkertProfile(**kwargs)

Bases: dmsky.density.DensityProfile

Burkert ApJ 447, L25 (1995) [Eqn. 2] http://arxiv.org/abs/astro-ph/9504041 rho(r) = rho0 * r0**3 / ( (r + r0)*(r**2+r0**2) ) ==> rho(r) = rhos / ( (1+r/rs)*(1+(r/rs)**2) )

class dmsky.density.NFWProfile(**kwargs)

Bases: dmsky.density.DensityProfile

Navarro, Frenk, and White, ApJ 462, 563 (1996) http://arxiv.org/abs/astro-ph/9508025 rho(r) = rhos / ((r/rs) * (1+r/rs)**2)

jvalue_fast(r=None)

Fast integrated J-factor computation

class dmsky.density.GNFWProfile(**kwargs)

Bases: dmsky.density.DensityProfile

Generalized NFW Profile Strigari et al. ApJ 678, 614 (2008) [Eqn. 3] http://arxiv.org/abs/0709.1510 rho(r) = rhos / ( (r/rs)**gamma * (1+r/rs)**(3-gamma))

deriv_params

Return the list of paramters we can take derivatives w.r.t.

class dmsky.density.EinastoProfile(**kwargs)

Bases: dmsky.density.DensityProfile

Einasto profile Einasto Trudy Inst. Astrofiz. Alma-Ata 5, 87 (1965) (Russian) [Eqn. 4] http://adsabs.harvard.edu/abs/1965TrAlm…5…87E rho(r) = rhos*exp(-2*((r/rs)**alpha-1)/alpha) ==>

deriv_params

Return the list of paramters we can take derivatives w.r.t.

class dmsky.density.ZhouProfile(**kwargs)

Bases: dmsky.density.DensityProfile

Generalized double power-law models Zhou MNRAS 278, 488 (1996) [Eqn. 1] http://arxiv.org/abs/astro-ph/9509122 rho(r) = C * (r/rs)**-gamma * (1 + (r/rs)**1/alpha))**-(beta-gamma)*alpha C = 4 * rhos

also see… Zhou MNRAS 287, 525 (1997) [Eqn. 2] http://arxiv.org/abs/astro-ph/9605029 Strigari et al., Nature 454 (2008) [Eqn. 8] http://arxiv.org/abs/0808.3772

deriv_params

Return the list of paramters we can take derivatives w.r.t.

Line of sight integration

class dmsky.jcalc.LoSFn(dp, d, xi, alpha=3.0)

Bases: object

Integrand function (luminosity density) for LoS integration. The parameter alpha introduces a change of variables:

x’ = x^(1/alpha).

A value of alpha > 1 samples the integrand closer to x = 0 (distance of closest approach). The change of variables requires the substitution:

dx = alpha * (x’)^(alpha-1) dx’

func(r)

Function to compute the integrand

class dmsky.jcalc.LoSAnnihilate(dp, d, xi, alpha=3.0)

Bases: dmsky.jcalc.LoSFn

Integrand function for LoS annihilation (J-factor).

func(r)

Function to compute the line-of-sight integrand

class dmsky.jcalc.LoSAnnihilate_Deriv(dp, d, xi, paramNames, alpha=3.0)

Bases: dmsky.jcalc.LoSFn

Integrand function for LoS annihilation (J-factor).

func(r)

Function to compute the line-of-sight integrand

class dmsky.jcalc.LoSDecay(dp, d, xi, alpha=3.0)

Bases: dmsky.jcalc.LoSFn

Integrand function for LoS decay (D-factor).

func(r)

Function to compute the line-of-sight integrand

class dmsky.jcalc.LoSDecay_Deriv(dp, d, xi, paramNames, alpha=1.0)

Bases: dmsky.jcalc.LoSFn

Integrand function for LoS decay (D-factor).

func(r)

Function to compute the line-of-sight integrand

class dmsky.jcalc.LoSIntegral(density, dhalo, alpha=3.0, ann=True, derivPar=None)

Bases: object

Slowest (and most accurate?) LoS integral. Uses scipy.integrate.quad with a change of variables to better sample the LoS close to the halo center.

angularIntegral(angle=None)

Compute the solid-angle integrated j-value within a given radius

Parameters:angle (numpy.ndarray or None) – Maximum integration angle (in degrees) If None, use the ‘rmax’ and ‘dhalo’ parameters. Returns:values – Return values, same shape as the input xp Return type:numpy.array

name

Return the name of this profile

rmax

Return the maximum integration radius

class dmsky.jcalc.LoSIntegralFast(density, dhalo, alpha=3.0, ann=True, nsteps=400, derivPar=None)

Bases: dmsky.jcalc.LoSIntegral

Vectorized version of LoSIntegral that performs midpoint integration with a fixed number of steps.

rmax

Return the maximum integration radius

class dmsky.jcalc.LoSIntegralInterp(density, dhalo, alpha=3.0, ann=True, nsteps=400, derivPar=None)

Bases: dmsky.jcalc.LoSIntegralFast

Interpolate fast integral a for even faster look-up.

create_func(dhalo)

Create the spline function

Parameters:dhalo (numpy.ndarray) – Array of halo distances Returns:func – A function that return J-factor as a function of psi and dhalo Return type:function

create_profile(dhalo, nsteps=None)

Create a spatial J-factor profile

Parameters:

• dhalo (numpy.ndarray) – Array of halo distances
• nsteps (int) – Number of steps for vectorization

Returns:

• dhalo, psi (numpy.meshgrid) – Array of halo distances and angular offsets
• jval (numpy.array) – Corresponding J-factors

class dmsky.jcalc.LoSIntegralFile(dp, dist, filename, ann=True)

Bases: dmsky.jcalc.LoSIntegralInterp

Interpolate over a pre-generated file.

NOT IMPLEMENTED YET

create_profile(dhalo, nsteps=300)

Build the profile values

class dmsky.jcalc.ROIIntegrator(jspline, lat_cut, lon_cut, source_list=None)

Bases: object

Class to integrate a J-factor over a region of interest

compute()

Integrate the ROI

eval(rgc, decay=False)

Evaluate the J-factor

print_profile(decay=False)

Print the profile

Dark matter targets

class dmsky.targets.Target(**kwargs)

Bases: pymodeler.model.Model

A DM search target

create_dmap(npix=150, subsample=4, coordsys='CEL', projection='AIT')

Create a D-factor map

Parameters:

• npix (int) – Number of pixels along one axis of output map
• subsample (int) – Number of subsamples to take per pixel
• coordsys (str) – Coordniate system: ‘GAL’ or ‘CEL’
• projection (str) – Map projection

Returns:

• image (numpy.ndarray) – Image data
• pix (numpy.ndarray) – Pixel coordinatates
• wcs (WCS.wcs) – WCS object for map

create_jmap(npix=150, subsample=4, coordsys='CEL', projection='AIT')

Create a J-factor map

Parameters:

• npix (int) – Number of pixels along one axis of output map
• subsample (int) – Number of subsamples to take per pixel
• coordsys (str) – Coordniate system: ‘GAL’ or ‘CEL’
• projection (str) – Map projection

Returns:

• image (numpy.ndarray) – Image data
• pix (numpy.ndarray) – Pixel coordinatates
• wcs (WCS.wcs) – WCS object for map

dsigma(ra, dec)

Return the uncertainty of J in any direction

Parameters:

• ra (numpy.ndarray) – Right-accension (in degrees)
• dec (numpy.ndarray) – Declination (in degrees)

Returns:

values – Return values, same shape as the input ra, dec

Return type:

numpy.array

dvalue(ra, dec)

Return the D-factor in any direction

Parameters:

• ra (numpy.ndarray) – Right-accension (in degrees)
• dec (numpy.ndarray) – Declination (in degrees)

Returns:

values – Return values, same shape as the input ra, dec

Return type:

numpy.array

jsigma(ra, dec)

Return the uncertainty of J in any direction

Parameters:

• ra (numpy.ndarray) – Right-accension (in degrees)
• dec (numpy.ndarray) – Declination (in degrees)

Returns:

values – Return values, same shape as the input ra, dec

Return type:

numpy.array

jvalue(ra, dec)

Return the J-factor in any direction

Parameters:

• ra (numpy.ndarray) – Right-accension (in degrees)
• dec (numpy.ndarray) – Declination (in degrees)

Returns:

values – Return values, same shape as the input ra, dec

Return type:

numpy.array

write_d_rad_file(d_rad_file=None, npts=50, minpsi=0.0001)

Write a text file with the D-fractor radial profile

Parameters:

• d_rad_file (str) – Filename to write to
• npts (int) – Number of angles to write
• minpsi (float) – Value for smallest angle to write

write_dmap(filename, npix=150, clobber=False, map_kwargs=None, file_kwargs=None)

Write the D-factor to a template map.

write_dmap_hpx(filename)

Write the D-factor to a template map.

write_dmap_wcs(filename, npix=150, clobber=False, map_kwargs=None, file_kwargs=None)

Write the D-factor to a template map.

write_j_rad_file(j_rad_file=None, npts=50, minpsi=0.0001)

Write a text file with the J-fractor radial profile

Parameters:

• j_rad_file (str) – Filename to write to
• npts (int) – Number of angles to write
• minpsi (float) – Value for smallest angle to write

write_jmap(filename, npix=150, clobber=False, map_kwargs=None, file_kwargs=None)

Write the J-factor to a template map.

write_jmap_hpx(filename)

Write the J-factor to a template map.

write_jmap_wcs(filename, npix=150, clobber=False, map_kwargs=None, file_kwargs=None)

Write the J-factor to a template map.

class dmsky.targets.Galactic(**kwargs)

Bases: dmsky.targets.Target

Class to add specifics for Galactic DM targets

class dmsky.targets.Dwarf(**kwargs)

Bases: dmsky.targets.Target

Class to add specifics for Dwarf Galaxy DM targets

class dmsky.targets.Galaxy(**kwargs)

Bases: dmsky.targets.Target

Class to add specifics for Galaxy DM targets

class dmsky.targets.Cluster(**kwargs)

Bases: dmsky.targets.Target

Class to add specifics for Galaxy Cluster DM targets

class dmsky.targets.Isotropic(**kwargs)

Bases: dmsky.targets.Target

Class to add specifics for Isotropic DM targets

Plotting fucntions

Module for plotting stuff in dmsky.

dmsky.plotting.plot_density(name, targets, xlims=None, nstep=100)

Make a plot of the density as a function of distance from the target center.

Parameters:

• name (str) – Name for the plot
• targets (list) – List of targets to include in the plot
• xlims (tuple) – Range for the x-axis of the plot
• nteps (int) – Number of points to include in the plot

Returns:

• fig (matplotlib.Figure)
• ax (matplotlib.Axes)
• leg (matplotlib.Legend)

dmsky.plotting.plot_j_integ(name, targets, xlims=None, nstep=100, ylims=None)

Make a plot of the integrated J-factor as a function of the angle from the target center.

Parameters:

• name (str) – Name for the plot
• targets (list) – List of targets to include in the plot
• xlims (tuple) – Range for the x-axis of the plot
• nteps (int) – Number of points to include in the plot
• ylims (tuple) – Range for the y-axis of the plot

Returns:

• fig (matplotlib.Figure)
• ax (matplotlib.Axes)
• leg (matplotlib.Legend)

dmsky.plotting.plot_j_profile(name, targets, xlims=None, nstep=100, ylims=None)

Make a plot of the J-factor as a function of the angle from the target center.

Parameters:

• name (str) – Name for the plot
• targets (list) – List of targets to include in the plot
• xlims (tuple) – Range for the x-axis of the plot
• nteps (int) – Number of points to include in the plot
• ylims (tuple) – Range for the y-axis of the plot

Returns:

• fig (matplotlib.Figure)
• ax (matplotlib.Axes)
• leg (matplotlib.Legend)

Skymap class

Roster class

class dmsky.roster.Roster(*targets)

Bases: collections.OrderedDict

A self-consistent set of search targets, typically with exactly one version for any given target.

Top-level bookkeeping

Factory for generating instances of classes

dmsky.factory.factory(cls, module=None, **kwargs)

Factory for creating objects. Arguments are passed directly to the constructor of the chosen class.

Parameters:

• cls (str) – Class name of the object to create
• module (str) – python module defining the class in question

Returns:

object – Newly created object

Return type:

object

class dmsky.library.ObjectLibrary(path=None)

Bases: object

Library that keeps track of object we are building

classmethod load_library(paths)

Build the library by walking through paths

class dmsky.targets.TargetLibrary(path=None)

Bases: dmsky.library.ObjectLibrary

A top-level object, keeping track of all the Target objects that we have created

create_target(name, version=None, **kwargs)

Create a Target

Parameters:

• name (str) – A name for the Target
• version (str) – Key that species which set of parameters to used for this target

Returns:

target – The newly created Target

Return type:

Target

get_target_dict(name, version=None, **kwargs)

Step through the various levels of dependencies to get the full dictionary for a target.

target: version -> … -> target: default -> default: type

class dmsky.roster.RosterLibrary(path=None)

Bases: dmsky.library.ObjectLibrary

A top-level object, keeping track of all the Roster objects that we have created

create_roster(name, *targets)

Create a roster

Parameters:

• name (str) – A name for the Roster
• targets (list) – Objects that will go on the Roster

Returns:

roster – The newly created Roster

Return type:

Roster

get_roster_list(name, *targets)

Get the list of objects on a roster, creating them in needed.

Parameters:

• name (str) – A name for the list
• targets (list) – Objects that will go on the list

Returns:

items – List of items needed to create a Roster

Return type:

list