Jammy2211 / PyAutoLens

PyAutoLens: Open Source Strong Gravitational Lensing
https://pyautolens.readthedocs.io/
MIT License
164 stars 32 forks source link

Galaxy class / RayTracingGalaxy / RayTracingPlane #21

Closed Jammy2211 closed 6 years ago

Jammy2211 commented 6 years ago

Okay - if you have a spare moment, it'd be good to get your opinion on this before I push ahead.

The galaxy class is a bit weird, if you look at the code, you'll see that we setup each Galaxy, and then feed them to the RayTracingPlane. The RayTracingPlane class then updates a lot of the properties of the galaxies, as they have a lot of properties that require knowledge of the entire ray tracing plane (e.g. the redshifts of the galaxies before / after each galaxy) in order to calculate.

This gets quite nasty, as it means the Galaxy class has a lot of methods which won't work until they've been passed through the RayTracingPlane class.

I'm going to split the Galaxy class in two, i.e. the Galaxy class as we use it now, and a RayTracingGalaxy, which the RayTracingPlane uses to set up its galaxies. The RayTracingGalaxy will inherite from the Galaxy class, and include all the methods which require knowledge of the RayTracingPlane to perform.

rhayes777 commented 6 years ago

Yeah that sounds really weird. I’ll try to take a look at it this evening.

On 11 Mar 2018, at 09:03, Jammy2211 notifications@github.com wrote:

Okay - if you have a spare moment, it'd be good to get your opinion on this before I push ahead.

The galaxy class is a bit weird, if you look at the code, you'll see that we setup each Galaxy, and then feed them to the RayTracingPlane. The RayTracingPlane class then updates a lot of the properties of the galaxies, as they have a lot of properties that require knowledge of the entire ray tracing plane (e.g. the redshifts of the galaxies before / after each galaxy) in order to calculate.

This gets quite nasty, as it means the Galaxy class has a lot of methods which won't work until they've been passed through the RayTracingPlane class.

I'm going to split the Galaxy class in two, i.e. the Galaxy class as we use it now, and a RayTracingGalaxy, which the RayTracingPlane uses to set up its galaxies. The RayTracingGalaxy will inherite from the Galaxy class, and include all the methods which require knowledge of the RayTracingPlane to perform.

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/Jammy2211/PyAutoLens/issues/21, or mute the thread https://github.com/notifications/unsubscribe-auth/AGETMp0zcwB9Yp-hcgR8fc-svLC2lJRQks5tdOhfgaJpZM4SlnTM.

Jammy2211 commented 6 years ago

I opted against making the changes now, as I figured we needed to have a clearer idea of the source analysis and what quantities we actually need. So I'm focusing on the source analysis.

But it is really weird. I think, to do a full lensing calculation, we need to know the angular diameter distance between a galaxy and all other galaxies, for every galaxy. From what, we could derive all the other physical quantities we want, but its whether we want to set those up in the RayTracingPlane class or generate them on the fly.

You'd probably want to set them up, but if you are varying your cosmology then they'd need refreshing for every new set of cosmological parameters 0_0.

rhayes777 commented 6 years ago

It sounds safer to generate them on the fly. So these properties could be calculated for the RayTracingPlane using the galaxies that it contains?

On 11 Mar 2018, at 10:12, Jammy2211 notifications@github.com wrote:

I opted against making the changes now, as I figured we needed to have a clearer idea of the source analysis and what quantities we actually need. So I'm focusing on the source analysis.

But it is really weird. I think, to do a full lensing calculation, we need to know the angular diameter distance between a galaxy and all other galaxies, for every galaxy. From what, we could derive all the other physical quantities we want, but its whether we want to set those up in the RayTracingPlane class or generate them on the fly.

You'd probably want to set them up, but if you are varying your cosmology then they'd need refreshing for every new set of cosmological parameters 0_0.

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/Jammy2211/PyAutoLens/issues/21#issuecomment-372103848, or mute the thread https://github.com/notifications/unsubscribe-auth/AGETMu-kNi9djheZRPrn3OwflrRUP9g9ks5tdPh8gaJpZM4SlnTM.

Jammy2211 commented 6 years ago

I think so. Basically, a number of the properties of a galaxy depend on the knowledge of the other galaxies. So like, a galaxy's critical surface density - its a quantity that 'belongs' to the galaxy, but requires knowledge of the angular diameter distances of it to the other galaxys. It's used to calculate that galaxy's mass, in physical units.

So, we could make these properties of the RayTracingPlane, but its a headache, as they are tied to galaxies at some level.

Like, if we have 20 galaxies in the RayTracingPlane, and we want the mass of 1, do we want to really compute all these quantities via the RayTracingPlane, or just have it stored as a variable for that 1 galaxy? Probably the latter, hence why I suggested having two Galaxy classes. But it always gets a bit confusing.

rhayes777 commented 6 years ago

Galaxies shouldn’t need to know about each other. How do angular diameter distances factor into critical surface density?

On 11 Mar 2018, at 10:23, Jammy2211 notifications@github.com wrote:

I think so. Basically, a number of the properties of a galaxy depend on the knowledge of the other galaxies. So like, a galaxy's critical surface density - its a quantity that 'belongs' to the galaxy, but requires knowledge of the angular diameter distances of it to the other galaxys. It's used to calculate that galaxy's mass, in physical units.

So, we could make these properties of the RayTracingPlane, but its a headache, as they are tied to galaxies at some level.

Like, if we have 20 galaxies in the RayTracingPlane, and we want the mass of 1, do we want to really compute all these quantities via the RayTracingPlane, or just have it stored as a variable for that 1 galaxy? Probably the latter, hence why I suggested having two Galaxy classes. But it always gets a bit confusing.

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/Jammy2211/PyAutoLens/issues/21#issuecomment-372104449, or mute the thread https://github.com/notifications/unsubscribe-auth/AGETMnQVoo6XN1kFkNG7PalN2eUY8LKJks5tdPsKgaJpZM4SlnTM.

Jammy2211 commented 6 years ago

Critical Surface density is (c^2 / 4piG) (D_l / D_s D_ls)

Where: c = speed of light G = Newton Constant D_l = Angular diameter distance of the lower redshift galaxy (i.e. the lens) to Earth D_s = Angular diameter distance of the higher redshift galaxy (i.e. the source) to Earth D_ls = Angular diameter distance between the two galaxies

Bare in mind, we might have situations with more than 2 galaxies, where we'll need to know D for different combinations or pairs. I need to determine if, for a large number of galaxies, we need every possible D, or just a reduced list of some variety.

So, as you can see above, for some lensing calculations the galaxies have to know each others redshifts. For example, when we compute deflection angles and trace coordinates from one galaxy's image-plane to the next, their redshifts (and therefore angular diameter distances) dictate the values of those deflection calculations. I guess this is why I made it a 'RayTracingPlane'.

If the galaxy redshifts and cosmology are fixed throughout the analysis, we only need to compute these quantities once, as they won't change. So, I guess we should compute them when we set up a new instance of RayTracingPlane - when we change galaxy redshift or cosmology we generate a new RayTracingPlane.

In terms of lens modeling, there are two common scenarios:

1) We have two planes (e..g the lens plane and source plane) - their redshifts might be known or might not. In this situation, we can perform all calculations in arcseconds and don't actually need any of the RayTracingPlane information, other than to turn arcsecond coordinates to physical units.

2) When we have 3 or more planes (e.g. two lenses at two different redshifts and a source). When this occurs, the redshifts of the lens galaxies and source (more specifically their ratio of angular diameter distances) determines key quantities like the deflection angles. So, we need the RayTracingPlane to perform any lens modeling.

rhayes777 commented 6 years ago

So is this necessary for the simple source - lens system? That is, if we focus on case 1 can we implement a solution without having to worry about this issue?

On 11 Mar 2018, at 10:52, Jammy2211 notifications@github.com wrote:

Critical Surface density is (c^2 / 4piG) (D_l / D_s D_ls)

Where: c = speed of light G = Newton Constant D_l = Angular diameter distance of the lower redshift galaxy (i.e. the lens) to Earth D_s = Angular diameter distance of the higher redshift galaxy (i.e. the source) to Earth D_ls = Angular diameter distance between the two galaxies

Bare in mind, we might have situations with more than 2 galaxies, where we'll need to know D for different combinations or pairs. I need to determine if, for a large number of galaxies, we need every possible D, or just a reduced list of some variety.

So, as you can see above, for some lensing calculations the galaxies have to know each others redshifts. For example, when we compute deflection angles and trace coordinates from one galaxy's image-plane to the next, their redshifts (and therefore angular diameter distances) dictate the values of those deflection calculations. I guess this is why I made it a 'RayTracingPlane'.

If the galaxy redshifts and cosmology are fixed throughout the analysis, we only need to compute these quantities once, as they won't change. So, I guess we should compute them when we set up a new instance of RayTracingPlane - when we change galaxy redshift or cosmology we generate a new RayTracingPlane.

In terms of lens modeling, there are two common scenarios:

We have two planes (e..g the lens plane and source plane) - their redshifts might be known or might not. In this situation, we can perform all calculations in arcseconds and don't actually need any of the RayTracingPlane information, other than to turn arcsecond coordinates to physical units.

When we have 3 or more planes (e.g. two lenses at two different redshifts and a source). When this occurs, the redshifts of the lens galaxies and source (more specifically their ratio of angular diameter distances) determines key quantities like the deflection angles. So, we need the RayTracingPlane to perform any lens modeling.

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/Jammy2211/PyAutoLens/issues/21#issuecomment-372106103, or mute the thread https://github.com/notifications/unsubscribe-auth/AGETMv1oexxfVUu2tCXHDUBiZcTdCoUvks5tdQH5gaJpZM4SlnTM.

Jammy2211 commented 6 years ago

For a source - lens system, we can perform all lens modeling without knowledge of the lens / source redshifts and other quantities. However, we use those quantities to translate coordinates from arcseconds to physical distances (parsecs) as well as compute physical quantities for things like mass. However, that's something we'd typically do after the actual lens modeling and analysis.

The 'mass_within_circles' routine is an example of one currently in Galaxy, where I'm using the critical density to convert the dimensionless mass to a physical mass. I guess we could move stuff like this to the RayTracingPlane class, e.g. 'mass_within_circle_of_lens_galaxies'. This makes me think having LensGalaxy, SourceGalaxy and LensAndSourceGalaxy classes could be worthwhile, as this dictates what quantities you can or can't measure for galaxy. Lets not do that yet though!

My immediate plan is to get the ray tracing calculation set up and working for a lens-source system, to give us a feel for how it all shapes up and looks. It gets pretty confusing when you're interfacing it with the image data, light/mass profiles and source-planes, so I think this'll help us make progress.

However, as said above, if we have >2 galaxies, the ray-tracing code is necessary for deflection angle calculations and therefore lens modeling. So well need to figure all this out sooner or later, but I suspect once we have the 2 system case sorted itll sorta design itself.

Jammy2211 commented 6 years ago

Okay, so we've set up a set of galaxies, with light profiles, mass profiles and redshifts. The mistake I've been making is trying to jump from the galaxies to an all encompassing 'RayTracingPlane' class.

What we need, are two intermediate classes:

1) ImagePlane - this is the set of coordinates that begins our ray tracing calculation (i.e. the image coordinates) and the galaxies that are in the image-plane (i.e. the main lens galaxy). By assigning galaxies to the image plane, we know how to compute deflection angles or make a light model in the image plane at its coordinates.

2) SourcePlane - this is the set of coordinates that we hit at each galaxy redshift (i.e. the image coordinates after deflection via foreground galaxies). Thus, for multiple galaxies with different redshifts, we will have multiple source-planes. A SourcePlane will be assigned Galaxys that define how a light model is fitted at its coordinates and how deflection angles are computed.

In the ray_tracing module, we can now have multiple classes for different ray_tracing setups (e.g. LensAndSource for just one image and source plane, LensAndSourceNoCosmology if we want to bypass cosmology, MultiplePlanes for multiple planes, etc.). If there are MultiplePlanes, the ImagePlane and SourcePlane can be set up automatically by the galaxy redshifts.The 'ray-tracing layer' therefore now takes care of all deflection angle computations.

This layer has no knowledge of the image data (only its coordinates). The next layer up ('data_fitting' or something) will then take our ray tracing, image and source-planes, and use them to generate model image data and fit it.

In https://github.com/Jammy2211/PyAutoLens/issues/11, I discussed how on a uniform coordinate grid we can make computation of things like light profiles efficiently, as the sub-gridding can be exploit the uniformity. This grid only occurs in the ImagePlane, thus we can directly tie these calculations to the plane a galaxy is assigned - its a lot easier.

Jammy2211 commented 6 years ago

Okay, so this is beginning to take shape, although there is lots more work to do!

When we call a RayTracing class, we pass it its image-plane coordinates using the PlaneCoordinates class (this groups the different coordinate sets, e.g. the image-grid, sub-grid, sparse-grid, blurring-region pixels etc..) We then use these PlaneCoordinates to set up an ImagePlane, with which we attach galaxies. Those ImagePlane galaxies are used in its constructor to compute the deflection angles at the image PlaneCoordinates.

The RayTracing class defines the number of lenses and sources, e.g. the most basic class is TraceLensSource, providing a single len-source system. It has a function, 'trace_to_next_plane', which uses the ImagePlane coordinates and deflection angles to trace the to the next pane, which for a lens-source system forms a SourcePlane. For a lens-source system the SourcePlane is the last plane, so we don't need to set up its deflection angles, only its LightProfile and Pixelization.

There will next be a TraceLensSourceCosmology class, which the above class inherites from. This class pretty much does the exact same as above, but sets up cosomology related information using an input set of cosmological parameters.

If our RayTracing class is a TraceMultiPlane, this means we have multiple lens galaxies (as many as we pass in). In this case, we set up the ImagePlane as before, but when we call 'trace_to_next_plane', there are two differences:

1) If the next galaxies is not the last one, a LensPlane instead of SourcePlane is set up. a LensPlane acts the same as an ImagePlane, in that it computes deflection angles, for the next ray-trace. However, unlike the Image-Plane, it is not on the uniform grid of image-coordinates.

2) The function 'trace_to_next_plane' now incorporates galaxy redshift and cosmological information to perform ray-tracing, as for multi-plane lens systems you must factor in the ratios of angulaar diameter distances between all the galaxies.

Once you know every set of Image, Lens and SourcePlanes, it is trivial to calculate either their image-plane model source images (via the Galaxy.LightProfiles) or their pixel-grid mapping matrices (via the Galaxy.Pixelization). These quantities will thus be generated in the DataModel layer, where things like the PSF will be available for convolution and the ImageData will be available for fitting the model.

Outstanding Issues

I think I will add a bool to the Galaxy class, e.g. 'fit_light_profile_in_image_plane', which when True means that the galaxy's light profile is passed to the ImagePlane class in a new attribute (e.g. 'higher_redshift_light_profiles'). Its light profile is then NOT passed to its actual LensPlane or SourcePlane. This still feels very messy though!

Jammy2211 commented 6 years ago

Following the last bullet point above, I envision each Plane class having their own routines for coordinate exaction, e.g.:

Jammy2211 commented 6 years ago

Currently, the RayTracing class has a class PlaneCoordinates, that represents the coordinates used for ray-tracing. This is messy, as the different sets of coordinates we use (image, sub, sparse, etc.) are being explicitly passed around the RayTracing class, when the class itself shouldn't care what coordinates we are actually ray-tracing.

The Planes that we set up however, do care, as to generate a model image from a light profile or source-plane pixelization we will be using multiple grids.

Therefore, I am going to do away with the PlaneCoordinates class, and in the image module have a class called AnalysisGrids. AnalysisGrids stores all the 2D coordinate grids (after translation to 1D) that we use for lens modeling. Thus, we pass RayTracing a set of AnalysisGrids, and it automatically performs ray-tracing on them, but the specifics of the grid types are irrelevent to RayTracing. The Plane classes are however set up with the complete set of AnalysisGrids, which have suitably been ray-tracing beforehand.

Jammy2211 commented 6 years ago

I'm going to call the class Grids and it will also store all the information on how different grids map to one another. E.g. grids.image, grids.sub, grids.sub_to_image, etc.

This is most likely building up to some all encompassing AnalysisData and Image class, which has the grids but also the psf, noise, etc.

Jammy2211 commented 6 years ago

The ray_tracing module and grids module have now fully taken shape, and look a lot like how I've described here. The key is that the grids module keeps the coordinates / data structure separate from the galaxy profiles.

Next, we need a module to analyse data using the output of the ray_tracing module.