diff --git a/autogalaxy/__init__.py b/autogalaxy/__init__.py index 5d1295df2..c19d6c983 100644 --- a/autogalaxy/__init__.py +++ b/autogalaxy/__init__.py @@ -60,6 +60,7 @@ from .ellipse.dataset_interp import DatasetInterp from .ellipse.ellipse.ellipse import Ellipse from .ellipse.ellipse.ellipse_multipole import EllipseMultipole +from .ellipse.ellipse.ellipse_multipole import EllipseMultipoleScaled from .ellipse.fit_ellipse import FitEllipse from .ellipse.model.analysis import AnalysisEllipse from .operate.image import OperateImage diff --git a/autogalaxy/config/notation.yaml b/autogalaxy/config/notation.yaml index 593c44555..e4481b452 100644 --- a/autogalaxy/config/notation.yaml +++ b/autogalaxy/config/notation.yaml @@ -19,6 +19,8 @@ label: ell_comps_1: \epsilon_{\rm 2} multipole_comps_0: M_{\rm 1} multipole_comps_1: M_{\rm 2} + scaled_multipole_comps_0: M_{\rm 1} + scaled_multipole_comps_1: M_{\rm 2} flux: F gamma: \gamma gamma_1: \gamma @@ -94,6 +96,8 @@ label_format: ell_comps_1: '{:.4f}' multipole_comps_0: '{:.4f}' multipole_comps_1: '{:.4f}' + input_multipole_comps_0: '{:.4f}' + input_multipole_comps_1: '{:.4f}' flux: '{:.4e}' gamma: '{:.4f}' inner_coefficient: '{:.4f}' diff --git a/autogalaxy/config/priors/ellipse/ellipse_multipole.yaml b/autogalaxy/config/priors/ellipse/ellipse_multipole.yaml index 3f5d34739..ca5dff9b7 100644 --- a/autogalaxy/config/priors/ellipse/ellipse_multipole.yaml +++ b/autogalaxy/config/priors/ellipse/ellipse_multipole.yaml @@ -19,3 +19,25 @@ EllipseMultipole: limits: lower: -inf upper: inf + +EllipseMultipoleRelative: + scaled_multipole_comps_0: + type: Uniform + lower_limit: -0.1 + upper_limit: 0.1 + width_modifier: + type: Absolute + value: 0.05 + gaussian_limits: + lower: -inf + upper: inf + scaled_multipole_comps_1: + type: Uniform + lower_limit: -0.1 + upper_limit: 0.1 + width_modifier: + type: Absolute + value: 0.05 + gaussian_limits: + lower: -inf + upper: inf diff --git a/autogalaxy/ellipse/ellipse/ellipse.py b/autogalaxy/ellipse/ellipse/ellipse.py index 8986d085f..a1eb347e4 100644 --- a/autogalaxy/ellipse/ellipse/ellipse.py +++ b/autogalaxy/ellipse/ellipse/ellipse.py @@ -41,7 +41,7 @@ class representing an ellispe, which is used to perform ellipse fitting to 2D da @property def circular_radius(self) -> float: """ - The radius of the circle that bounds the ellipse, assuming that the `major_axis` is the radius of the circle. + The circumference of the circle that bounds the ellipse, assuming that the `major_axis` is the radius of the circle. """ return 2.0 * np.pi * np.sqrt((2.0 * self.major_axis**2.0) / 2.0) diff --git a/autogalaxy/ellipse/ellipse/ellipse_multipole.py b/autogalaxy/ellipse/ellipse/ellipse_multipole.py index 89ac08517..f9cfbb690 100644 --- a/autogalaxy/ellipse/ellipse/ellipse_multipole.py +++ b/autogalaxy/ellipse/ellipse/ellipse_multipole.py @@ -1,5 +1,6 @@ import numpy as np from typing import Tuple +from autogalaxy.convert import multipole_comps_from, multipole_k_m_and_phi_m_from from autogalaxy.ellipse.ellipse.ellipse import Ellipse @@ -12,7 +13,7 @@ def __init__( multipole_comps: Tuple[float, float] = (0.0, 0.0), ): """ - class representing the multipole of an ellispe with, which is used to perform ellipse fitting to + class representing the multipole of an ellipse, which is used to perform ellipse fitting to 2D data (e.g. an image). The multipole is added to the (y,x) coordinates of an ellipse that are already computed via the `Ellipse` class. @@ -37,6 +38,35 @@ class representing the multipole of an ellispe with, which is used to perform el self.m = m self.multipole_comps = multipole_comps + def get_shape_angle( + self, + ellipse: Ellipse, + ) -> float: + """ + The shape angle is the offset between the angle of the ellipse and the angle of the multipole, + this defines the shape that the multipole takes. + + In the case of the m=4 multipole, angles of 0 indicate pure diskiness, angles +- 45 + indicate pure boxiness. + + Parameters + ---------- + ellipse + The base ellipse profile that is perturbed by the multipole. + + Returns + ------- + The angle between the ellipse and the multipole, in degrees between +- 180/m. + """ + + angle = ellipse.angle-multipole_k_m_and_phi_m_from(self.multipole_comps, self.m)[1] + if angle < -180/self.m: + angle += 360/self.m + elif angle > 180/self.m: + angle -= 360/self.m + + return angle + def points_perturbed_from( self, pixel_scale, points, ellipse: Ellipse, n_i: int = 0 ) -> np.ndarray: @@ -56,15 +86,116 @@ def points_perturbed_from( ------- The (y,x) coordinates of the input points, which are perturbed by the multipole. """ + symmetry = 360 / self.m + k_orig, phi_orig = multipole_k_m_and_phi_m_from(self.multipole_comps, self.m) + comps_adjusted = multipole_comps_from(k_orig, + symmetry-2*phi_orig+(symmetry-(ellipse.angle-phi_orig)), #Re-align light to match mass + self.m) + + # 1) compute cartesian (polar) angle + theta = np.arctan2(points[:,0], points[:,1]) # <- true polar angle + + # 2) multipole in that same frame + delta_theta = self.m * (theta - ellipse.angle_radians) + radial = ( + comps_adjusted[1] * np.cos(delta_theta) + + comps_adjusted[0] * np.sin(delta_theta) + ) - angles = ellipse.angles_from_x0_from(pixel_scale=pixel_scale, n_i=n_i) + # 3) perturb along the true radial direction + x = points[:, 1] + radial * np.cos(theta) + y = points[:, 0] + radial * np.sin(theta) - radial = np.add( - self.multipole_comps[1] * np.cos(self.m * (angles - ellipse.angle_radians)), - self.multipole_comps[0] * np.sin(self.m * (angles - ellipse.angle_radians)), + return np.stack((y, x), axis=-1) + +class EllipseMultipoleScaled(EllipseMultipole): + def __init__( + self, + m=4, + scaled_multipole_comps: Tuple[float, float] = (0.0, 0.0), + major_axis=1., + ): + """ + class representing the multipole of an ellipse, which is used to perform ellipse fitting to + 2D data (e.g. an image). This multipole is fit with its strength held relative to an ellipse with a + major_axis of 1, allowing for a set of ellipse multipoles to be fit at different major axes but with + the same scaled strength k/a. + + The scaled_multipole_comps (for all ellipses) are converted to a k value, which is then reset to + its `true' value for a multipole at the given major axis value, which is then used to perturb an ellipse + as per the normal `EllipseMultipole' class and below. + + The multipole is added to the (y,x) coordinates of an ellipse that are already computed via the `Ellipse` class. + + The addition of the multipole is performed as follows: + + :math: r_m = \sum_{i=1}^{m} \left( a_i \cos(i(\theta - \phi)) + b_i \sin(i(\theta - \phi)) \right) + :math: y_m = r_m \sin(\theta) + :math: x_m = r_m \cos(\theta) + + Where: + + m = The order of the multipole. + r = The radial coordinate of the ellipse perturbed by the multipole. + \phi = The angle of the ellipse. + a = The amplitude of the cosine term of the multipole. + b = The amplitude of the sine term of the multipole. + y = The y-coordinate of the ellipse perturbed by the multipole. + x = The x-coordinate of the ellipse perturbed by the multipole. + """ + + + self.scaled_multipole_comps = scaled_multipole_comps + k, phi = multipole_k_m_and_phi_m_from(multipole_comps=scaled_multipole_comps, m=m) + k_adjusted = k*major_axis + + specific_multipole_comps = multipole_comps_from(k_adjusted, phi, m) + + super().__init__(m, specific_multipole_comps) + + self.specific_multipole_comps = specific_multipole_comps + self.m = m + + def points_perturbed_from( + self, pixel_scale, points, ellipse: Ellipse, n_i: int = 0 + ) -> np.ndarray: + """ + Returns the (y,x) coordinates of the input points, which are perturbed by the multipole of the ellipse. + + Parameters + ---------- + pixel_scale + The pixel scale of the data that the ellipse is fitted to and interpolated over. + points + The (y,x) coordinates of the ellipse that are perturbed by the multipole. + ellipse + The ellipse that is perturbed by the multipole, which is used to compute the angles of the ellipse. + + Returns + ------- + The (y,x) coordinates of the input points, which are perturbed by the multipole. + """ + symmetry = 360/self.m + k_orig, phi_orig = multipole_k_m_and_phi_m_from(self.specific_multipole_comps, self.m) + comps_adjusted = multipole_comps_from(k_orig, + symmetry-2*phi_orig+(symmetry-(ellipse.angle-phi_orig)), #Re-align light to match mass + self.m) + + # 1) compute cartesian (polar) angle + theta = np.arctan2(points[:, 0], points[:, 1]) # <- true polar angle + + # 2) multipole in that same frame + delta_theta = self.m * (theta - ellipse.angle_radians) + radial = ( + comps_adjusted[1] * np.cos(delta_theta) + + comps_adjusted[0] * np.sin(delta_theta) ) - x = points[:, 1] + (radial * np.cos(angles)) - y = points[:, 0] + (radial * np.sin(angles)) + # 3) perturb along the true radial direction + x = points[:, 1] + radial * np.cos(theta) + y = points[:, 0] + radial * np.sin(theta) + + return np.stack((y, x), axis=-1) + + - return np.stack(arrays=(y, x), axis=-1) diff --git a/autogalaxy/ellipse/fit_ellipse.py b/autogalaxy/ellipse/fit_ellipse.py index 764bd9dbc..1a0f3edaf 100644 --- a/autogalaxy/ellipse/fit_ellipse.py +++ b/autogalaxy/ellipse/fit_ellipse.py @@ -109,7 +109,7 @@ def points_from_major_axis_from(self) -> np.ndarray: raise ValueError( """ - The code has attempted to add over 1000 points to the ellipse and still not found a set of points that + The code has attempted to add over 300 extra points to the ellipse and still not found a set of points that do not hit the mask with the expected number of points. This is likely due to the mask being too large or a strange geometry, and the code is unable to find a diff --git a/autogalaxy/ellipse/plot/fit_ellipse_plotters.py b/autogalaxy/ellipse/plot/fit_ellipse_plotters.py index d4dbbae7a..39bd1c0b0 100644 --- a/autogalaxy/ellipse/plot/fit_ellipse_plotters.py +++ b/autogalaxy/ellipse/plot/fit_ellipse_plotters.py @@ -116,14 +116,15 @@ def figures_2d( for_subplot=for_subplot, ) - def subplot_fit_ellipse(self): + def subplot_fit_ellipse(self, disable_data_contours: bool = False): """ Standard subplot of the attributes of the plotter's `FitEllipse` object. """ self.open_subplot_figure(number_subplots=2) - self.figures_2d(data=True) + self.mat_plot_2d.use_log10 = True + self.figures_2d(data=True, disable_data_contours=disable_data_contours) self.figures_2d(ellipse_residuals=True, for_subplot=True) self.mat_plot_2d.output.subplot_to_figure(auto_filename="subplot_fit_ellipse") @@ -227,7 +228,7 @@ def subplot_ellipse_errors(self): visuals_2d=visuals_2d, auto_labels=aplt.AutoLabels( title=f"Ellipse Fit", - filename=f"subhplot_ellipse_errors", + filename=f"subplot_ellipse_errors", ), ) diff --git a/autogalaxy/plot/mat_plot/two_d.py b/autogalaxy/plot/mat_plot/two_d.py index e96af4758..62682ca8d 100644 --- a/autogalaxy/plot/mat_plot/two_d.py +++ b/autogalaxy/plot/mat_plot/two_d.py @@ -25,7 +25,6 @@ def __init__( output: Optional[aplt.Output] = None, array_overlay: Optional[aplt.ArrayOverlay] = None, contour: Optional[aplt.Contour] = None, - fill: Optional[aplt.Fill] = None, grid_scatter: Optional[aplt.GridScatter] = None, grid_plot: Optional[aplt.GridPlot] = None, vector_yx_quiver: Optional[aplt.VectorYXQuiver] = None, @@ -197,7 +196,6 @@ def __init__( xlabel=xlabel, text=text, annotate=annotate, - fill=fill, output=output, origin_scatter=origin_scatter, mask_scatter=mask_scatter, diff --git a/autogalaxy/plot/visuals/two_d.py b/autogalaxy/plot/visuals/two_d.py index 1c33e86bf..1a44e245e 100644 --- a/autogalaxy/plot/visuals/two_d.py +++ b/autogalaxy/plot/visuals/two_d.py @@ -1,230 +1,230 @@ -from typing import List, Union, Optional - -import autoarray as aa -import autoarray.plot as aplt - - -class Visuals2D(aplt.Visuals2D): - def __init__( - self, - origin: aa.Grid2D = None, - border: aa.Grid2D = None, - mask: aa.Mask2D = None, - lines: Optional[Union[List[aa.Array1D], aa.Grid2DIrregular]] = None, - positions: Optional[Union[aa.Grid2DIrregular, List[aa.Grid2DIrregular]]] = None, - grid: Union[aa.Grid2D] = None, - mesh_grid: aa.Grid2D = None, - vectors: aa.VectorYX2DIrregular = None, - patches: "Union[ptch.Patch]" = None, - fill_region: Optional[List] = None, - array_overlay: aa.Array2D = None, - light_profile_centres: aa.Grid2DIrregular = None, - mass_profile_centres: aa.Grid2DIrregular = None, - multiple_images: aa.Grid2DIrregular = None, - tangential_critical_curves: Optional[ - Union[aa.Grid2DIrregular, List[aa.Grid2DIrregular]] - ] = None, - radial_critical_curves: Optional[ - Union[aa.Grid2DIrregular, List[aa.Grid2DIrregular]] - ] = None, - tangential_caustics: Optional[ - Union[aa.Grid2DIrregular, List[aa.Grid2DIrregular]] - ] = None, - radial_caustics: Optional[ - Union[aa.Grid2DIrregular, List[aa.Grid2DIrregular]] - ] = None, - parallel_overscan=None, - serial_prescan=None, - serial_overscan=None, - indexes: Union[List[int], List[List[int]]] = None, - ): - super().__init__( - mask=mask, - positions=positions, - grid=grid, - lines=lines, - mesh_grid=mesh_grid, - vectors=vectors, - patches=patches, - fill_region=fill_region, - array_overlay=array_overlay, - origin=origin, - border=border, - parallel_overscan=parallel_overscan, - serial_prescan=serial_prescan, - serial_overscan=serial_overscan, - indexes=indexes, - ) - - self.light_profile_centres = light_profile_centres - self.mass_profile_centres = mass_profile_centres - self.multiple_images = multiple_images - self.tangential_critical_curves = tangential_critical_curves - self.radial_critical_curves = radial_critical_curves - self.tangential_caustics = tangential_caustics - self.radial_caustics = radial_caustics - - def plot_via_plotter(self, plotter, grid_indexes=None): - super().plot_via_plotter( - plotter=plotter, - grid_indexes=grid_indexes, - ) - - if self.light_profile_centres is not None: - plotter.light_profile_centres_scatter.scatter_grid( - grid=self.light_profile_centres - ) - - if self.mass_profile_centres is not None: - plotter.mass_profile_centres_scatter.scatter_grid( - grid=self.mass_profile_centres - ) - - if self.multiple_images is not None: - plotter.multiple_images_scatter.scatter_grid( - grid=self.multiple_images.array - ) - - if self.tangential_critical_curves is not None: - try: - plotter.tangential_critical_curves_plot.plot_grid( - grid=self.tangential_critical_curves - ) - except TypeError: - pass - - if self.radial_critical_curves is not None: - try: - plotter.radial_critical_curves_plot.plot_grid( - grid=self.radial_critical_curves - ) - except TypeError: - pass - - if self.tangential_caustics is not None: - try: - try: - plotter.tangential_caustics_plot.plot_grid( - grid=self.tangential_caustics - ) - except (AttributeError, ValueError): - plotter.tangential_caustics_plot.plot_grid( - grid=self.tangential_caustics.array - ) - except TypeError: - pass - - if self.radial_caustics is not None: - try: - plotter.radial_caustics_plot.plot_grid(grid=self.radial_caustics) - except TypeError: - pass - - def add_critical_curves_or_caustics( - self, mass_obj, grid: aa.type.Grid2DLike, plane_index: int - ): - """ - From a object with mass profiles (e.g. mass profile, galaxy) extract the critical curves or caustics and - returns them in a `Visuals2D` object. - - This includes support for a `plane_index`, which specifies the index of the plane in the tracer, which is - an object used in PyAutoLens to represent a lensing system with multiple planes (e.g. an image plane and a - source plane). The `plane_index` allows for the extraction of quantities from a specific plane in the tracer. - - When plotting a `Tracer` it is common for plots to only display quantities corresponding to one plane at a time - (e.g. the convergence in the image plane, the source in the source plane). Therefore, quantities are only - extracted from one plane, specified by the input `plane_index`. - - Parameters - ---------- - mass_obj - The mass object (e.g. mass profile, galaxy, tracer) object which has attributes extracted for plotting. - grid - The 2D grid of (y,x) coordinates used to plot the tracer's quantities in 2D. - plane_index - The index of the plane in the tracer which is used to extract quantities, as only one plane is plotted - at a time. - - Returns - ------- - vis.Visuals2D - A collection of attributes that can be plotted by a `Plotter` object. - """ - if plane_index == 0: - return self.add_critical_curves(mass_obj=mass_obj, grid=grid) - return self.add_caustics(mass_obj=mass_obj, grid=grid) - - def add_critical_curves(self, mass_obj, grid: aa.type.Grid2DLike): - """ - From a object with mass profiles (e.g. mass profile, galaxy) extract the critical curves and - returns them in a `Visuals2D` object. - - When plotting a `Tracer` it is common for plots to only display quantities corresponding to one plane at a time - (e.g. the convergence in the image plane, the source in the source plane). Therefore, quantities are only - extracted from one plane, specified by the input `plane_index`. - - Parameters - ---------- - mass_obj - The mass object (e.g. mass profile, galaxy, tracer) object which has attributes extracted for plotting. - grid - The 2D grid of (y,x) coordinates used to plot the tracer's quantities in 2D. - plane_index - The index of the plane in the tracer which is used to extract quantities, as only one plane is plotted - at a time. - - Returns - ------- - vis.Visuals2D - A collection of attributes that can be plotted by a `Plotter` object. - """ - - tangential_critical_curves = mass_obj.tangential_critical_curve_list_from( - grid=grid - ) - - radial_critical_curves = None - radial_critical_curve_area_list = mass_obj.radial_critical_curve_area_list_from( - grid=grid - ) - - if any([area > grid.pixel_scale for area in radial_critical_curve_area_list]): - radial_critical_curves = mass_obj.radial_critical_curve_list_from(grid=grid) - - return self + self.__class__( - tangential_critical_curves=tangential_critical_curves, - radial_critical_curves=radial_critical_curves, - ) - - def add_caustics(self, mass_obj, grid: aa.type.Grid2DLike): - """ - From a object with mass profiles (e.g. mass profile, galaxy) extract the caustics and - returns them in a `Visuals2D` object. - - When plotting a `Tracer` it is common for plots to only display quantities corresponding to one plane at a time - (e.g. the convergence in the image plane, the source in the source plane). Therefore, quantities are only - extracted from one plane, specified by the input `plane_index`. - - Parameters - ---------- - mass_obj - The mass object (e.g. mass profile, galaxy, tracer) object which has attributes extracted for plotting. - grid - The 2D grid of (y,x) coordinates used to plot the tracer's quantities in 2D. - plane_index - The index of the plane in the tracer which is used to extract quantities, as only one plane is plotted - at a time. - - Returns - ------- - vis.Visuals2D - A collection of attributes that can be plotted by a `Plotter` object. - """ - - tangential_caustics = mass_obj.tangential_caustic_list_from(grid=grid) - radial_caustics = mass_obj.radial_caustic_list_from(grid=grid) - - return self + self.__class__( - tangential_caustics=tangential_caustics, - radial_caustics=radial_caustics, - ) +from typing import List, Union, Optional + +import autoarray as aa +import autoarray.plot as aplt + + +class Visuals2D(aplt.Visuals2D): + def __init__( + self, + origin: aa.Grid2D = None, + border: aa.Grid2D = None, + mask: aa.Mask2D = None, + lines: Optional[Union[List[aa.Array1D], aa.Grid2DIrregular]] = None, + positions: Optional[Union[aa.Grid2DIrregular, List[aa.Grid2DIrregular]]] = None, + grid: Union[aa.Grid2D] = None, + mesh_grid: aa.Grid2D = None, + vectors: aa.VectorYX2DIrregular = None, + patches: "Union[ptch.Patch]" = None, + fill_region: Optional[List] = None, + array_overlay: aa.Array2D = None, + light_profile_centres: aa.Grid2DIrregular = None, + mass_profile_centres: aa.Grid2DIrregular = None, + multiple_images: aa.Grid2DIrregular = None, + tangential_critical_curves: Optional[ + Union[aa.Grid2DIrregular, List[aa.Grid2DIrregular]] + ] = None, + radial_critical_curves: Optional[ + Union[aa.Grid2DIrregular, List[aa.Grid2DIrregular]] + ] = None, + tangential_caustics: Optional[ + Union[aa.Grid2DIrregular, List[aa.Grid2DIrregular]] + ] = None, + radial_caustics: Optional[ + Union[aa.Grid2DIrregular, List[aa.Grid2DIrregular]] + ] = None, + parallel_overscan=None, + serial_prescan=None, + serial_overscan=None, + indexes: Union[List[int], List[List[int]]] = None, + ): + super().__init__( + mask=mask, + positions=positions, + grid=grid, + lines=lines, + mesh_grid=mesh_grid, + vectors=vectors, + patches=patches, + fill_region=fill_region, + array_overlay=array_overlay, + origin=origin, + border=border, + parallel_overscan=parallel_overscan, + serial_prescan=serial_prescan, + serial_overscan=serial_overscan, + indexes=indexes, + ) + + self.light_profile_centres = light_profile_centres + self.mass_profile_centres = mass_profile_centres + self.multiple_images = multiple_images + self.tangential_critical_curves = tangential_critical_curves + self.radial_critical_curves = radial_critical_curves + self.tangential_caustics = tangential_caustics + self.radial_caustics = radial_caustics + + def plot_via_plotter(self, plotter, grid_indexes=None): + super().plot_via_plotter( + plotter=plotter, + grid_indexes=grid_indexes, + ) + + if self.light_profile_centres is not None: + plotter.light_profile_centres_scatter.scatter_grid( + grid=self.light_profile_centres + ) + + if self.mass_profile_centres is not None: + plotter.mass_profile_centres_scatter.scatter_grid( + grid=self.mass_profile_centres + ) + + if self.multiple_images is not None: + plotter.multiple_images_scatter.scatter_grid( + grid=self.multiple_images.array + ) + + if self.tangential_critical_curves is not None: + try: + plotter.tangential_critical_curves_plot.plot_grid( + grid=self.tangential_critical_curves + ) + except TypeError: + pass + + if self.radial_critical_curves is not None: + try: + plotter.radial_critical_curves_plot.plot_grid( + grid=self.radial_critical_curves + ) + except TypeError: + pass + + if self.tangential_caustics is not None: + try: + try: + plotter.tangential_caustics_plot.plot_grid( + grid=self.tangential_caustics + ) + except (AttributeError, ValueError): + plotter.tangential_caustics_plot.plot_grid( + grid=self.tangential_caustics.array + ) + except TypeError: + pass + + if self.radial_caustics is not None: + try: + plotter.radial_caustics_plot.plot_grid(grid=self.radial_caustics) + except TypeError: + pass + + def add_critical_curves_or_caustics( + self, mass_obj, grid: aa.type.Grid2DLike, plane_index: int + ): + """ + From a object with mass profiles (e.g. mass profile, galaxy) extract the critical curves or caustics and + returns them in a `Visuals2D` object. + + This includes support for a `plane_index`, which specifies the index of the plane in the tracer, which is + an object used in PyAutoLens to represent a lensing system with multiple planes (e.g. an image plane and a + source plane). The `plane_index` allows for the extraction of quantities from a specific plane in the tracer. + + When plotting a `Tracer` it is common for plots to only display quantities corresponding to one plane at a time + (e.g. the convergence in the image plane, the source in the source plane). Therefore, quantities are only + extracted from one plane, specified by the input `plane_index`. + + Parameters + ---------- + mass_obj + The mass object (e.g. mass profile, galaxy, tracer) object which has attributes extracted for plotting. + grid + The 2D grid of (y,x) coordinates used to plot the tracer's quantities in 2D. + plane_index + The index of the plane in the tracer which is used to extract quantities, as only one plane is plotted + at a time. + + Returns + ------- + vis.Visuals2D + A collection of attributes that can be plotted by a `Plotter` object. + """ + if plane_index == 0: + return self.add_critical_curves(mass_obj=mass_obj, grid=grid) + return self.add_caustics(mass_obj=mass_obj, grid=grid) + + def add_critical_curves(self, mass_obj, grid: aa.type.Grid2DLike): + """ + From a object with mass profiles (e.g. mass profile, galaxy) extract the critical curves and + returns them in a `Visuals2D` object. + + When plotting a `Tracer` it is common for plots to only display quantities corresponding to one plane at a time + (e.g. the convergence in the image plane, the source in the source plane). Therefore, quantities are only + extracted from one plane, specified by the input `plane_index`. + + Parameters + ---------- + mass_obj + The mass object (e.g. mass profile, galaxy, tracer) object which has attributes extracted for plotting. + grid + The 2D grid of (y,x) coordinates used to plot the tracer's quantities in 2D. + plane_index + The index of the plane in the tracer which is used to extract quantities, as only one plane is plotted + at a time. + + Returns + ------- + vis.Visuals2D + A collection of attributes that can be plotted by a `Plotter` object. + """ + + tangential_critical_curves = mass_obj.tangential_critical_curve_list_from( + grid=grid + ) + + radial_critical_curves = None + radial_critical_curve_area_list = mass_obj.radial_critical_curve_area_list_from( + grid=grid + ) + + if any([area > grid.pixel_scale for area in radial_critical_curve_area_list]): + radial_critical_curves = mass_obj.radial_critical_curve_list_from(grid=grid) + + return self + self.__class__( + tangential_critical_curves=tangential_critical_curves, + radial_critical_curves=radial_critical_curves, + ) + + def add_caustics(self, mass_obj, grid: aa.type.Grid2DLike): + """ + From a object with mass profiles (e.g. mass profile, galaxy) extract the caustics and + returns them in a `Visuals2D` object. + + When plotting a `Tracer` it is common for plots to only display quantities corresponding to one plane at a time + (e.g. the convergence in the image plane, the source in the source plane). Therefore, quantities are only + extracted from one plane, specified by the input `plane_index`. + + Parameters + ---------- + mass_obj + The mass object (e.g. mass profile, galaxy, tracer) object which has attributes extracted for plotting. + grid + The 2D grid of (y,x) coordinates used to plot the tracer's quantities in 2D. + plane_index + The index of the plane in the tracer which is used to extract quantities, as only one plane is plotted + at a time. + + Returns + ------- + vis.Visuals2D + A collection of attributes that can be plotted by a `Plotter` object. + """ + + tangential_caustics = mass_obj.tangential_caustic_list_from(grid=grid) + radial_caustics = mass_obj.radial_caustic_list_from(grid=grid) + + return self + self.__class__( + tangential_caustics=tangential_caustics, + radial_caustics=radial_caustics, + ) diff --git a/autogalaxy/profiles/mass/total/power_law_multipole.py b/autogalaxy/profiles/mass/total/power_law_multipole.py index 16e5a2a88..780fe0b44 100644 --- a/autogalaxy/profiles/mass/total/power_law_multipole.py +++ b/autogalaxy/profiles/mass/total/power_law_multipole.py @@ -7,6 +7,7 @@ from autogalaxy import convert from autogalaxy.profiles.mass.abstract.abstract import MassProfile +from autogalaxy.profiles.mass.total import PowerLaw def radial_and_angle_grid_from( @@ -128,6 +129,36 @@ def __init__( ) self.angle_m *= units.deg.to(units.rad) + def get_shape_angle( + self, + base_profile: PowerLaw, + ) -> float: + """ + The shape angle is the offset between the angle of the ellipse and the angle of the multipole, + this defines the shape that the multipole takes. + + In the case of the m=4 multipole, angles of 0 indicate pure diskiness, angles +- 45 + indicate pure boxiness. + + Parameters + ---------- + base_profile + The base power-law mass profile that is perturbed by the multipole. + + Returns + ------- + The angle between the ellipse and the multipole, in degrees, between +- 180/m. + """ + + angle = convert.angle_from(base_profile.ell_comps) - convert.multipole_k_m_and_phi_m_from(self.multipole_comps, self.m)[1] + if angle < -180/self.m: + angle += 360/self.m + elif angle > 180/self.m: + angle -= 360/self.m + + return angle + + def jacobian( self, a_r: np.ndarray, a_angle: np.ndarray, polar_angle_grid: np.ndarray ) -> Tuple[np.ndarray, Tuple]: @@ -226,3 +257,5 @@ def potential_2d_from(self, grid: aa.type.Grid2DLike, **kwargs) -> np.ndarray: The grid of (y,x) arc-second coordinates the deflection angles are computed on. """ return jnp.zeros(shape=grid.shape[0]) + + diff --git a/test_autogalaxy/ellipse/ellipse/test_multipole.py b/test_autogalaxy/ellipse/ellipse/test_multipole.py index b3d075a11..aa30c9899 100644 --- a/test_autogalaxy/ellipse/ellipse/test_multipole.py +++ b/test_autogalaxy/ellipse/ellipse/test_multipole.py @@ -24,5 +24,20 @@ def test__points_perturbed_from(): pixel_scale=pixel_scale, points=points, ellipse=ellipse ) - assert points_perturbed[1, 0] == pytest.approx(-0.982728, 1.0e-4) + assert points_perturbed[1, 0] == pytest.approx(-0.919384, 1.0e-4) assert points_perturbed[1, 1] == pytest.approx(0.298726, 1.0e-4) + + + ellipse = ag.Ellipse(major_axis=1.5) + + points = ellipse.points_from_major_axis_from(pixel_scale=pixel_scale) + + multipole_scaled = ag.EllipseMultipoleScaled(m=4, scaled_multipole_comps=(0.1, 0.2), major_axis=1.5) + + points_perturbed = multipole_scaled.points_perturbed_from( + pixel_scale=pixel_scale, points=points, ellipse=ellipse + ) + + assert points_perturbed[1, 0] == pytest.approx(-0.848528, 1.0e-4) + assert points_perturbed[1, 1] == pytest.approx(0.848528, 1.0e-4) + diff --git a/test_autogalaxy/ellipse/test_fit_ellipse.py b/test_autogalaxy/ellipse/test_fit_ellipse.py index 110a4e4a7..245b64b2f 100644 --- a/test_autogalaxy/ellipse/test_fit_ellipse.py +++ b/test_autogalaxy/ellipse/test_fit_ellipse.py @@ -63,8 +63,8 @@ def test___points_from_major_axis__multipole(imaging_lh): dataset=imaging_lh, ellipse=ellipse_0, multipole_list=[multipole] ) - assert fit._points_from_major_axis[1, 0] == pytest.approx(-0.542453, 1.0e-4) - assert fit._points_from_major_axis[1, 1] == pytest.approx(-0.038278334, 1.0e-4) + assert fit._points_from_major_axis[1, 0] == pytest.approx(-0.119588, 1.0e-4) + assert fit._points_from_major_axis[1, 1] == pytest.approx(0.038856679, 1.0e-4) # def test__mask_interp(imaging_lh, imaging_lh_masked): @@ -197,3 +197,5 @@ def test__log_likelihood(imaging_lh, imaging_lh_masked): fit = ag.FitEllipse(dataset=imaging_lh_masked, ellipse=ellipse_0) assert fit.log_likelihood == pytest.approx(-0.169821080058, 1.0e-4) + +