Source code for xga.sourcetools.deproj

#  This code is a part of X-ray: Generate and Analyse (XGA), a module designed for the XMM Cluster Survey (XCS).
#  Last modified by David J Turner ( 20/02/2023, 14:04. Copyright (c) The Contributors

from typing import Union

import numpy as np
from astropy.units import Quantity, UnitConversionError

[docs]def shell_ann_vol_intersect(shell_radii: Union[np.ndarray, Quantity], ann_radii: Union[np.ndarray, Quantity]) \ -> Union[np.ndarray, Quantity]: """ This function calculates the volume intersection matrix of a set of circular annuli and a set of spherical shells. It is assumed that the annuli and shells have the same x and y origin. The intersection is derived using simple geometric considerations, have a look in the appendix of DOI 10.1086/300836. :param ndarray/Quantity shell_radii: The radii of the spherical shells. :param ndarray/Quantity ann_radii: The radii of the circular annuli (DOES NOT need to be the same length as shell_radii). :return: A 2D array containing the volumes of intersections between the circular annuli defined by i_ann and o_ann, and the spherical shells defined by i_sph and o_sph. Annular radii are along the 'x' axis and shell radii are along the 'y' axis. :rtype: Union[np.ndarray, Quantity] """ if all([type(shell_radii) == Quantity, type(ann_radii) == Quantity]) and shell_radii.unit != ann_radii.unit: raise UnitConversionError("If quantities are passed, they must be in the same units.") elif all([type(shell_radii) == Quantity, type(ann_radii) == Quantity]): pass elif all([type(shell_radii) == np.ndarray, type(ann_radii) == np.ndarray]): pass else: raise TypeError("shell_radii and ann_radii must either both be astropy quantities or numpy arrays, " "you cannot mix the two") i_ann, i_sph = np.meshgrid(ann_radii[0:-1], shell_radii[0:-1]) o_ann, o_sph = np.meshgrid(ann_radii[1:], shell_radii[1:]) # The main term which makes use of the radii of the shells and annuli. # The use of clip enforces that none of the terms can be less than 0, as we don't care # about those intersections, you can't have a negative volume. The None passed to clip is just to tell it # that there is no upper limit that we wish to enforce. main_term = np.power(np.clip((o_sph ** 2 - i_ann ** 2), 0, None), 3 / 2) - \ np.power(np.clip((o_sph ** 2 - o_ann ** 2), 0, None), 3 / 2) + \ np.power(np.clip((i_sph ** 2 - o_ann ** 2), 0, None), 3 / 2) - \ np.power(np.clip((i_sph ** 2 - i_ann ** 2), 0, None), 3 / 2) # Multiply by the necessary constants and return. return (4 / 3) * np.pi * main_term
[docs]def shell_volume(inn_radius: Quantity, out_radius: Quantity) -> Union[Quantity, np.ndarray]: """ Silly little function that calculates the volume of a spherical shell with inner radius inn_radius and outer radius out_radius. :param Quantity/np.ndarray inn_radius: The inner radius of the spherical shell. :param Quantity/np.ndarray out_radius: The outer radius of the spherical shell. :return: The volume of the specified shell :rtype: Union[Quantity, np.ndarray] """ if all([type(inn_radius) == Quantity, type(out_radius) == Quantity]) and inn_radius.unit != out_radius.unit: raise UnitConversionError("If quantities are passed, they must be in the same units.") elif all([type(inn_radius) == Quantity, type(out_radius) == Quantity]): pass elif all([type(inn_radius) == np.ndarray, type(out_radius) == np.ndarray]): pass else: raise TypeError("inn_radius and out_radius must either both be astropy quantities or numpy arrays, " "you cannot mix the two") outer_vol = (4/3) * np.pi * out_radius**3 inner_vol = (4/3) * np.pi * inn_radius**3 return outer_vol - inner_vol
# def temp_onion(proj_prof: ProjectedGasTemperature1D) -> GasTemperature3D: # """ # This function will generate deprojected, three-dimensional, gas temperature profile from a projected profile using # the 'onion peeling' deprojection method. The function is an implementation of a fairly old technique, though it # has been used recently in For a more in depth discussion of this # technique and its uses I would currently recommend # # :param ProjectedGasTemperature1D proj_prof: A projected cluster temperature profile, which the # user wants to use to infer the 3D temperature profile. # :return: The deprojected temperature profile. # :rtype: GasTemperature3D # """ # # pass