TS Map

Classes

class cosipy.ts_map.TSMap(*args, **kwargs)[source]

Compute the TS map of using threeML package.

Load the model and plugins

Parameters:
  • dr (str) – Path to full detector response.

  • data (histpy.Histogram) – Binned data. Note: Eventually this should be a cosipy data class.

  • bkg (histpy.Histogram) – Binned background model. Note: Eventually this should be a cosipy data class.

  • sc_orientation (cosipy.spacecraftfile.SpacecraftFile) – Contains the information of the orientation: timestamps (astropy.Time) and attitudes (scoord.Attitude) that describe the spacecraft for the duration of the data included in the analysis.

  • piv (float) – The pivotal energy of the spectrum.

  • index (float) – The index of the spectrum.

  • other_plugins (threeML.plugins, optional) – The plugins from other instruments.

  • norm (int, optional) – The norm of the spectrum model (the default is 1).

  • ra (float, optional) – The RA of the source model (the default is 0).

  • dec (float, optional) – The Dec of the source model (the default is 0).

instantiate_plugin()[source]

Instantiate the likelihood plugin.

gather_all_plugins()[source]

Gather all the plugins togather into a DataList.

create_model()[source]

Create the source model.

Returns:

The source model.

Return type:

astromodels.core.model.Model

fix_index()[source]

Return the index of the source spectrum.

ts_fitting()[source]

Peform the ts fitting.

ts_grid_data()[source]

Perform the ts fitting using the data on the different pixels.

ts_grid_bkg()[source]

Perform the ts fitting using the background on the different pixels.

calculate_ts()[source]

Calculate the TS by the TS of data and background.

print_best_fit()[source]

Print the best fit location.

save_ts(output_file_name)[source]

Save the TS map.

Parameters:

output_file_name (str) – The path to save the ts map.

load_ts(input_file_name)[source]

Load a ts map from file.

Parameters:

input_file_name (str) – The path to the saved TS map file.

refit_best_fit()[source]

Refit the best fit to check norm.

plot_ts_map()[source]

Plot the TS map.

class cosipy.ts_map.FastTSMap(data, bkg_model, response_path, orientation=None, cds_frame='local', scheme='RING')[source]
static slice_energy_channel(hist, channel_start, channel_stop)[source]

Slice one or more bins along first axis of the histogram.

Parameters:
  • hist (histpy.Histogram) – The histogram object to be sliced.

  • channel_start (int) – The start of the slice (inclusive).

  • channel_stop (int) – The stop of the slice (exclusive).

Returns:

sliced_hist – The sliced histogram.

Return type:

histpy.Histogram

static get_hypothesis_coords(nside, scheme='RING', coordsys='galactic')[source]

Get a list of hypothesis coordinates.

Parameters:
  • nside (int) – The nside of the map.

  • scheme (str, optional) – The scheme of the map where the hypothesis coordinates are generated (the default is “RING”, which implies the “RING” scheme is used to get the hypothesis coordinates).

  • coordsys (str, optional) – The coordinate system used in the map where the hypothesis coordinates are generated (the default is “galactic”, which implies the galactic coordinates system is used).

Returns:

hypothesis_coords – The list of the hypothesis coordinates at the center of each pixel.

Return type:

list

static get_cds_array(hist, energy_channel)[source]

Get the flattened cds array from input Histogram.

Parameters:
  • hist (histpy.Histogram) – The input Histogram.

  • energy_channel (list) – The format is [lower_channel, upper_chanel]. The lower_channel is inclusive while the upper_channel is exclusive.

Returns:

cds_array – The flattended Compton data space (CDS) array.

Return type:

numpy.ndarray

static get_psr_in_galactic(hypothesis_coord, response_path, spectrum)[source]

Get the point source response (psr) in galactic. Please be aware that you must use a galactic response! To do: to make the weight parameter not hardcoded

Parameters:
  • hypothesis_coord (astropy.coordinates.SkyCoord) – The hypothesis coordinate.

  • response_path (str or path.lib.Path) – The path to the response.

  • spectrum (astromodels.functions) – The spectrum of the source to be placed at the hypothesis coordinate.

Returns:

psr – The point source response of the spectrum at the hypothesis coordinate.

Return type:

histpy.Histogram

static get_ei_cds_array(hypothesis_coord, energy_channel, response_path, spectrum, cds_frame, orientation=None)[source]

Get the expected counts in CDS in local or galactic frame.

Parameters:
  • hypothesis_coord (astropy.coordinates.SkyCoord) – The hypothesis coordinate.

  • energy_channel (list) – The format is [lower_channel, upper_chanel]. The lower_channel is inclusive while the upper_channel is exclusive.

  • response_path (str or pathlib.Path) – The path to the response file.

  • spectrum (astromodels.functions) – The spectrum of the source.

  • cds_frame (str, optional) – “local” or “galactic”, it’s the Compton data space (CDS) frame of the data, bkg_model and the response. In other words, they should have the same cds frame.

  • orientation (cosipy.spacecraftfile.SpacecraftFile, optional) – The orientation of the spacecraft when data are collected (the default is None, which implies the orientation file is not needed).

Returns:

cds_array – The flattended Compton data space (CDS) array.

Return type:

numpy.ndarray

static fast_ts_fit(hypothesis_coord, energy_channel, data_cds_array, bkg_model_cds_array, orientation, response_path, spectrum, cds_frame, ts_nside, ts_scheme)[source]

Perform a TS fit on a single location at hypothesis_coord.

Parameters:
  • hypothesis_coord (astropy.coordinates.SkyCoord) – The hypothesis coordinate.

  • energy_channel (list) – The format is [lower_channel, upper_chanel]. The lower_channel is inclusive while the upper_channel is exclusive.

  • data_cds_array (numpy.ndarray) – The flattened Compton data space (CDS) array of the data.

  • bkg_model_cds_array (numpy.ndarray) – The flattened Compton data space (CDS) array of the background model.

  • orientation (cosipy.spacecraftfile.SpacecraftFile) – The orientation of the spacecraft when data are collected.

  • response_path (str or pathlib.Path) – The path to the response file.

  • spectrum (astromodels.functions) – The spectrum of the source.

  • cds_frame (str) – “local” or “galactic”, it’s the Compton data space (CDS) frame of the data, bkg_model and the response. In other words, they should have the same cds frame .

  • ts_nside (int) – The nside of the ts map.

  • ts_scheme (str) – The scheme of the Ts map.

Returns:

The list of the resulting TS fit: [pix number, ts value, norm, norm_err, failed, iterations, time_ei_cds_array, time_fit, time_fast_ts_fit]

Return type:

list

parallel_ts_fit(hypothesis_coords, energy_channel, spectrum, ts_scheme='RING', start_method='fork', cpu_cores=None, ts_nside=None)[source]

Perform parallel computation on all the hypothesis coordinates.

Parameters:
  • hypothesis_coords (list) – A list of the hypothesis coordinates to fit

  • energy_channel (list) – the energy channel you want to use: [lower_channel, upper_channel]. lower_channel is inclusive while upper_channel is exclusive.

  • spectrum (astromodels.functions) – The spectrum of the source.

  • ts_scheme (str, optional) – The scheme of the TS map (the default is “RING”, which implies a “RING” scheme of the TS map).

  • start_method (str, optional) – The starting method of the parallel computation (the default is “fork”, which implies using the fork method to start parallel computation).

  • cpu_cores (int, optional) – The number of cpu cores you wish to use for the parallel computation (the default is None, which implies using all the available number of cores -1 to perform the parallel computation).

  • ts_nside (int, optional) – The nside of the ts map. This must be given if the number of hypothesis_coords isn’t equal to the number of pixels of the total ts map, which means that you fit only a portion of the total ts map. (the default is None, which means that you fit the full ts map).

Returns:

results – The result of the ts fit over all the hypothesis coordinates.

Return type:

numpy.ndarray

plot_ts(ts_array=None, skycoord=None, containment=None, save_plot=False, save_dir='', save_name='ts_map.png', dpi=300)[source]

Plot the containment region of the TS map.

Parameters:
  • skyoord (astropy.coordinates.SkyCoord, optional) – The true location of the source (the default is None, which implies that there are no coordiantes to be printed on the TS map).

  • containment (float, optional) – The containment level of the source (the default is None, which will plot raw TS values).

  • save_plot (bool, optional) – Set True to save the plot (the default is False, which means it won’t save the plot.

  • save_dir (str or pathlib.Path, optional) – The directory to save the plot.

  • save_name (str, optional) – The file name of the plot to be save.

  • dpi (int, optional) – The dpi for plotting and saving.

static get_chi_critical_value(containment=0.9)[source]

Get the critical value of the chi^2 distribution based ob the confidence level.

Parameters:

containment (float, optional) – The confidence level of the chi^2 distribution (the default is 0.9, which implies that the 90% containment region).

Returns:

The critical value corresponds to the confidence level.

Return type:

float

static show_memory_info(hint)[source]
class cosipy.ts_map.FastNormFit(max_iter=1000, conv_frac_tol=0.001, zero_ts_tol=1e-05, allow_negative=False)[source]

Perform a fast poisson maximum likelihood ratio fit of the normalization of a source over background.

The likelihood ratio as a function of the norm is computed as follow

\[TS(N) = 2 \sum_i \left( \frac{\log P(d_i; b_i+N e_i)}{\log P(d_i; b_i)} \right)\]

where \(P(d; \lambda)\) is the Poisson probability of obtaining \(d\) count where \(\lambda\) is expected on average; \(b\) is the estimated number of background counts; \(N\) is the normalization; and \(e\) is the expected excess -i.e. signal- per normalization unit -i.e. the number of excess counts equals \(N\).

It can be shown that \(TS(N)\) has analytic derivative of arbitrary order and that the Newton’s method is guaranteed to converge if initialized at \(N=0\).

Note

The background is not allowed to float. It is assumed the error on the estimation of the background is small compared to the fluctuation of the background itself (i.e. \(N_{B}/N_{off} \lll 1\)).

Note

Because of the Poisson probability, \(TS(N)\) is only well-defined for \(N \geq 1\). By default, underfluctuations are set to \(TS(N=0) = 0\). For cases when there is benefit in letting the normalization float to negative values, you can use allow_negative, but in that case the results are only valid in the Gaussian regime.

Parameters:
  • max_iter (int) – Maximum number of iteration

  • conv_frac_tol (float) – Convergence stopping condition, expressed as the ratio between the size of the last step and the current norm value.

  • zero_ts_tol (float) – If zero_ts_tol < TS < 0, then TS is set to 0 without failed flag status (analytically, TS < 0 should never happen)

  • allow_negative (bool) – Allow the normalization to float toward negative values

static ts(data, bkg, unit_excess, norm)[source]

Get TS for a given normalization.

Parameters:
  • data (array) – Observed counts

  • bkg (array) – Background estimation. Same size as data. Every element should be >0

  • unit_excess (array) – Expected excess per normalization unit. Same size as data.

  • norm (float or array) – Normalization value

Returns:

TS, same shape as norm

Return type:

float or array

static dts(data, bkg, unit_excess, norm, order=1)[source]

Get the derivative of TS with respecto to the normalization, at given normalization.

Parameters:
  • data (array) – Observed counts

  • bkg (array) – Background estimation. Same size as data. Every element should be >0

  • unit_excess (array) – Expected excess per normalization unit. Same size as data.

  • norm (float or array) – Normalization value

  • order (int) – Derivative order

Returns:

d^n TS / dN^n, same shape as norm

Return type:

float or array

solve(data, bkg, unit_excess)[source]

Get the maximum TS, fitted normalization and normalization error (\(\Delta TS = 1\))

Note

The normalization error is obtained from approximating the TS function as as parabola (i.e. valid in the Gaussian regime). TS and norm are indeed valid in the Poisson regime.

Parameters:
  • data (array) – Observed counts

  • bkg (array) – Background estimation. Same size as data. Every element should be >0

  • unit_excess (array) – Expected excess per normalization unit. Same size as data.

Returns:

ts, norm, norm error and status (0 = good)

Return type:

(float, float, float, bool)