Very simple atmosphere module

[barnabytprowe]

Assigned: Joerg, Barney, Rachel, Mike

Create an image of a very simple PSF that approximates the effects due to the atmosphere (e.g. a simple parametric Gaussian, Kolmogorov - I think Kaiser 2000 has some good material for this). This could be circular to begin with, or elliptical if we want to be ambitious. This PSF will convolved together with the PSF from the optics module to provide the overall PSF (pixelation not included) for the GalSim General Milestone 1 images.

[rmandelb]

Shall we discuss our approach here on github so we have a record of the conversation that others can see if they so choose, without spamming the entire great3-code list?

Given that SBProfile has analytic treatment of Gaussians, it seems like that's the obvious place to start for a placeholder. If we want atmosphere profiles that are more complex / extended than that while still relying only on existing SBProfile code, we could make it a double Gaussian.

If we want Kolmogorov for this milestone, we might still find it useful to initially code it up in terms of a Gaussian.

[joergdietrich]

I think this is a sensible approach. I'll start a new branch now that I have a little more time (although I still have to write a VLT proposal until Thursday), and have started looking into SBProfile and its Python bindings. This should also give me a safe playground to learn more about git.

[rmandelb]

Hi Joerg - thanks for this. Feel free to e-mail, or post here if you have questions. (I'm going to pull @rmjarvis into this too, since I think he also volunteered for this module, but if I'm wrong he can always turn off notifications.) Barney and I were talking about this earlier, and we think that you should have what you need to make a Gaussian / double Gaussian / Moffat as SBProfiles. If we want to implement Kolmogorov (which would be nice, but not required for this milestone if we're lacking time / energy) then we have a choice to make:

• do we want to implement it in C++ and wrap it, similar to the many other parametric models that are in SBProfile already, or
• do we want to simply implement this in the python layer for the atmosphere module, which should still be reasonably fast / efficient / etc.?

I think I would prefer the former, but it's probably worth discussing. I'm also willing to put in some work to make that happen, e.g. if you want to take responsibility for an initial version that takes advantage of existing SBProfile functionality, and I could plan to extend it to Kolmogorov.

[joergdietrich]

I don't speak C++, which settles the question if I write the code. I'm not sure I'll manage to implement Kolmogorov before the milestone. Why would you prefer C++ if python is fast enough? In the end it shouldn't matter to the user whether something is implemented in C++ or python.

[barnabytprowe]

So the reason I could see for adding a C++ generation only makes sense if we were to create an SBKolmogorov class as a derived class of SBProfile, but then I think it does make sense.

Currently, all objects we wish to convolve must be made into some sort of SBProfile object in order for us to use the convolution class that's already built in to SBProfiles. At convolution, these classes are able to use some information about the Surface Brightness Profiles they contain to find efficiency savings, such as whether or not the profile is Axisymmetric, or whether it's analytic in real space, in k space, etc.

The options for making atmospheric PSFs based on multiple Gaussians, or Moffats, are super easy: we should be able simply to generate instances of the wrapped Python classes SBGaussian & SBMoffat. For the Kolmogorov we have two options:

i) Create an image of the Kolmogorov PSF as a numpy array in Python, initialize an Image class object using this array and then use this Image class object to initialize an SBPixel instance, which then contains the Kolmogorov PSF profile expressed as a 2D lookup table plus interpolant.

ii) Write a new C++ class SBKolmogorov, which would doubtless be something somewhat similar to SBMoffat or SBExponential. This would contain just the radial profile, but could be made elliptical etc. using the standard SBDistort class. When SBConvolve comes to use SBKolmogorov in convolution with other objects, and then you ask it to draw the image, Gary's code is smart enough to take relevant shortcuts when dealing with certain objects that can be expressed simply as analytic profiles. This may mean that there are some efficiencies over using an SBPixel object to store the Kolmogorov PSF, as well as the more hard-to-quantify benefit that this way you have a full description of the profile in its functional form, rather than just a lookup table for it.

I think ii) has some arguments in favour of it, for consistency and accuracy of representation of the object as well as possibly efficiency. The problem with ii) is that it is significantly more time consuming than i) [hence why I went for option i) when writing the optics module, where it would have been worse still], and moreover we are not expecting these simple parametric forms of atmospheric PSFs to be the final story for the atmosphere module. We're expecting to at least try to have some sort of turbulent propagation model (maybe ray-tracing) to generate atmospheric PSF patterns, we're just not sure where from yet. So the Kolmogorov profile is only one of the many types of PSF we'll be generating in the atmosphere module.

I think whatever you go for we can upgrade (or not) later anyway if required. I agree with you that the end user won't care if something is in C++ or Python, and I think that most of the "plumbing" for GalSim should happen in Python. This is just one of those occasions where a dedicated SBKolmogorov derived class might make a little sense. But it's clearly not a crucial issue, so I am very happy if you go with option i) for the time being.

On Mar 26, 2012, at 1:05 PM, Jrg Dietrich wrote:

I don't speak C++, which settles the question if I write the code. I'm not sure I'll manage to implement Kolmogorov before the milestone. Why would you prefer C++ if python is fast enough? In the end it shouldn't matter to the user whether something is implemented in C++ or python.

---

Reply to this email directly or view it on GitHub: https://github.com/GalSim-developers/GalSim/issues/25#issuecomment-4702887

Barnaby Rowe Postdoctoral Research Associate

Jet Propulsion Laboratory California Institute of Technology MS 169-237 4800 Oak Grove Drive Pasadena CA 91109 United States of America

Department of Physics & Astronomy University College London Gower Street London WC1E 6BT United Kingdom

[rmjarvis]

If Jörg is more comfortable writing it in Python, I think that would be fine. Then if it turns out that this is too slow, I don't mind converting it to C++, which is more my native language, so that wouldn't be so bad for me. I suspect there will be lots of things like this where we start out in Python and then just convert the things that are too slow into C++.

[rmandelb]

Thank you, Barney, for explaining all that more clearly than I did! And I certainly agree that given time constraints, the initial goal should be something that uses existing SB, then if there's time, a Kolmogorov implemented in the way that works best for whoever has time to implement it. We can always decide later to change the details of that implementation... this is our initial atmosphere module, not the final one!

[rmandelb]

Hi Joerg - just wanted to touch base about this. I know you had a proposal deadline just now, but would you like to chat about the atmosphere module more at the beginning of next week?

[joergdietrich]

Hi Rachel,

Sure. I don't think I'll manage to have Kolmogorov for the first milestone, but I hope not too long after that.

As Barney said, multiple Gaussians or Moffats are indeed super easy, we can just return the SBProfile instances. Very simple code to do that for single Gaussian/Moffat is in #25_simple_atmosphere. The question in my mind is now, how do we want to encapsulate that, and that is, I think, closely related to the still open discussion on issue #59. I don't think it would be good practice to define a GSAtmo object, which just duplicates code that is in GSObject and/or GSGaussian/GSMoffat. On the other hand Kolmogorov (or more detailed raytracing) is not something that will appear in any other part of the code.

So in my mind the question is not so much, what will go into the first milestone (that's easily coded up in an afternoon), but rather what should the relationship between the atmosphere part of the code and the other python GSObjects be?

Cheers, Joerg

On Fri, Mar 30, 2012 at 02:06:30PM -0700, Rachel Mandelbaum wrote:

Hi Joerg - just wanted to touch base about this. I know you had a proposal deadline just now, but would you like to chat about the atmosphere module more at the beginning of next week?

---

Reply to this email directly or view it on GitHub: https://github.com/GalSim-developers/GalSim/issues/25#issuecomment-4850914

Jörg Dietrich jorgd@umich.edu http://www-personal.umich.edu/~jorgd/ Physics Department, University of Michigan, Tel. +1 734 615 4256 450 Church Street Fax +1 734 763 9694 Ann Arbor, MI 48109-1040, USA

[rmandelb]

No worries about Kolmogorov, we can put that on the list for the future.

I agree about duplication of GSGaussian/GSMoffat code with a bunch of case statements being undesirable. Perhaps we should see how the rest of the discussion of #59 goes today (I believe we will wrap it up) and that will clear up the way forward on the atmosphere.

[barnabytprowe]

Joerg, you raise a very good point. The way the discussion progressed regarding base classes, and the proposal I made (after some thought) has ended up making it unclear what parts of a separate atmospheric module are needed to meet the milestone requirements of a simple profile (a Gaussian or Moffat).

Clearly we are going to need a good bit of code in the atmos module when we attempt some form of raytracing/photon-shooting treatment of the atmosphere (or some approximation to that!). But right now, something along the lines of Mike's suggested user calls in the discussion regarding pull request #60 would achieve most of what we want (paraphrasing):

if psfType == 'Gaussian'
    atmospsf = galsim.Gaussian()
elif psfType == 'Moffat'
    atmospsf = galsim.Moffat()
else
    raise 'Unknown atmospheric PSF Type'

...and this would not need to take place in a separate module, but could happen in a top level script for generating the postage stamp.

This makes me feel bad, as I have unintentionally duplicated work that you have already done in this branch - unfortunately, the class hierarchy I fixed upon in issues #59/#60 seems to be the best way of preserving the common methods of the SBProfile objects by inheritance (which is desirable) without inheriting the C++ classes directly (which seems not to be). I'm sure it is not the last word on the low-level design, and I fully expect it to be extended, but it does provide a very direct and consistent way of instantiating these simple objects.

I don't want to waffle, and I think what we do now is totally up to you Joerg! Things we could maybe consider might be... (all optional)

i) Pausing the atmosphere.py module until we start thinking about real atmospheres - which will be possibly as early as next week I think.

ii) Coding up a double Gaussian class (Rachel says that in Sloan this was often very useful) and leaving it there!

iii) Going ahead and allowing an Atmos class (I've dropped the GS thanks to Mike) to contain SBMoffat, SBGaussian and/or SBAdd (for double Gaussian) objects as attributes if initialized that way, BUT also having additional freedom to also contain more complex profiles in future depending on the init call.

Item iii) would result in a little duplication, although this is not necessarily the end of the world. Gary's code has SBSersic, SBExponential and SBDeVaucouleurs classes, obviously with the latter two as special cases of the former. Given that the definition of new classes is so lightweight given the base class hierearchy outlined in #59/#60, we are not likely to bring in huge potential for inconsistencies by having potential overlap between Moffat/Gaussian and some instances of the Atmos class. On the other hand: "There should be one-- and preferably only one --obvious way to do it.", from the Zen of Python, is a good principle...

I have totally waffled! Joerg, your input on this has already been brilliant and I hope you haven't felt your efforts were wasted - it's a steep learning curve for all of us, frankly. What we do next is totally up to you, and depends on your schedule and wishes. Currently we now don't 100% need an atmosphere module to meet the terms of the milestone, and it is my fault that this was not clear at the beginning: sorry for that, and for making a set of python classes that cover moffats and Gaussians. However, what I am not sorry about is your involvement - I'm delighted you're engaged with this problem. Please feel free to make up your mind about the next steps! As I said, if you wish, we can wait until next week... But a double Gaussian (or sum of arbitrary Gaussians) would be cool too :)

P.S. The unit tests you have written are also very useful and could be added to SBProfile unit tests immediately!

[joergdietrich]

Barney,

No worries about the way we ended up with you duplicating some of the things I did. This is a learning experience for me at least as much as it is for you. To get a little bit more learning out of this, I decided to go with your option (ii) for the time being, and while we think a little more how we want to integrate the atmospheric PSF into the overall structure of GalSim.

What I did now is to create a GSAdd class as a python class that has an SBAdd attribute, very similar to what you did for GSObject/SBProfile. Then I derive a DoubleGaussian class from that. I pushed this to #25_simple_atmosphere but will sleep once over it before I open a pull request. I mostly see this as a training branch for myself, with its usefulness to the entire project TBD.

On Mon, Apr 02, 2012 at 05:21:06PM -0700, Barnaby Rowe wrote:

Joerg, you raise a very good point. The way the discussion progressed regarding base classes, and the proposal I made (after some thought) has ended up making it unclear what parts of a separate atmospheric module are needed to meet the milestone requirements of a simple profile (a Gaussian or Moffat).

Clearly we are going to need a good bit of code in the atmos module when we attempt some form of raytracing/photon-shooting treatment of the atmosphere (or some approximation to that!). But right now, something along the lines of Mike's suggested user calls in the discussion regarding pull request #60 would achieve most of what we want (paraphrasing):

if psfType == 'Gaussian'
    atmospsf = galsim.Gaussian()
elif psfType == 'Moffat'
    atmospsf = galsim.Moffat()
else
    raise 'Unknown atmospheric PSF Type'

...and this would not need to take place in a separate module, but could happen in a top level script for generating the postage stamp.

This makes me feel bad, as I have unintentionally duplicated work that you have already done in this branch - unfortunately, the class hierarchy I fixed upon in issues #59/#60 seems to be the best way of preserving the common methods of the SBProfile objects by inheritance (which is desirable) without inheriting the C++ classes directly (which seems not to be). I'm sure it is not the last word on the low-level design, and I fully expect it to be extended, but it does provide a very direct and consistent way of instantiating these simple objects.

I don't want to waffle, and I think what we do now is totally up to you Joerg! Things we could maybe consider might be... (all optional)

i) Pausing the atmosphere.py module until we start thinking about real atmospheres - which will be possibly as early as next week I think.

ii) Coding up a double Gaussian class (Rachel says that in Sloan this was often very useful) and leaving it there!

iii) Going ahead and allowing an Atmos class (I've dropped the GS thanks to Mike) to contain SBMoffat, SBGaussian and/or SBAdd (for double Gaussian) objects as attributes if initialized that way, BUT also having additional freedom to also contain more complex profiles in future depending on the init call.

Item iii) would result in a little duplication, although this is not necessarily the end of the world. Gary's code has SBSersic, SBExponential and SBDeVaucouleurs classes, obviously with the latter two as special cases of the former. Given that the definition of new classes is so lightweight given the base class hierearchy outlined in #59/#60, we are not likely to bring in huge potential for inconsistencies by having potential overlap between Moffat/Gaussian and some instances of the Atmos class. On the other hand: "There should be one-- and preferably only one --obvious way to do it.", from the Zen of Python, is a good principle...

I have totally waffled! Joerg, your input on this has already been brilliant and I hope you haven't felt your efforts were wasted - it's a steep learning curve for all of us, frankly. What we do next is totally up to you, and depends on your schedule and wishes. Currently we now don't 100% need an atmosphere module to meet the terms of the milestone, and it is my fault that this was not clear at the beginning: sorry for that, and for making a set of python classes that cover moffats and Gaussians. However, what I am not sorry about is your involvement - I'm delighted you're engaged with this problem. Please feel free to make up your mind about the next steps! As I said, if you wish, we can wait until next week... But a double Gaussian (or sum of arbitrary Gaussians) would be cool too :)

P.S. The unit tests you have written are also very useful and could be added to SBProfile unit tests immediately!

---

Reply to this email directly or view it on GitHub: https://github.com/GalSim-developers/GalSim/issues/25#issuecomment-4891904

Jörg Dietrich jorgd@umich.edu http://www-personal.umich.edu/~jorgd/ Physics Department, University of Michigan, Tel. +1 734 615 4256 450 Church Street Fax +1 734 763 9694 Ann Arbor, MI 48109-1040, USA

[barnabytprowe]

Merged, thanks Joerg!