Fix a bug in the PSF initialization

[aguinot]

In the current implementation there is an assumption that only one image will be made per file. If one wants to simulate multiple SCAs and save them into a single file the PSF will only be initialized with the settings of the first image (SCA, WCS, bandpass). This is because the __init__ of the class in the InputLoader is only called at the start of a new file in GalSim. This hasn't been a problem so far because the simulation was always run with one image per file. With this update we can now safely make multiple images of different SCA within a file and the PSF will be updated accordingly by doing the initialization in the BuildImage function using the method setupImage (see where this happen here).
The update also include a small change regarding how the initialization is done. Now that we do it from the image builder we have access to some image characteristic which include its size. Instead of hard coding the image size to roman.n_pix for the corners, we now use the actual image size. This gives more flexibility to the code and allow to use images with a different sizes and have a more "appropriate" PSF interpolation.

Extra note: the __init__ of RomanPSF could almost be removed but I prefer to keep it because I have in mind features that would require it.

Validation: This implementation has been tested against an images generated using the config from this repo before modifications.

[aguinot]

The commit 9a0f746 propose a new implementation of the Roman PSF Builder.
This implementation adds a single register object: RomanPSF but with a lot of flexibility. It can call directly the galsim.roman.getPSF function with arguments passed in the config file. Or call the PSF builder through an interpolator (the current implementation that is used). The interpolator as been changed to allow for new interpolation scheme. I am planing to add a "constant" interpollation in the future where the PSF would be built only at the center of the image and reused everywhere for the shear simulation on small fields.
Bellow is a diagram of how the implementation works:

RomanPSF
    |
    |------ No interpolation ---- Call the galsim.roman.getPSF function and parameters of this function
    |                             can be set from the config file
    |
    |
    |------ Interpolation ------- It requires the `interpolator` keyword in the config file and you need
                |                 to specify the type of interpolator, at the moment there is only
                |                 `RomanPSFInterpolator` available.
                |
                |------- Interpolator loader ---- To use an interpolator it needs to be initialized and
                                |                 for that we need to declare it in the `input` section of
                                |                 the config file as: `RomanPSFInterpolator`.
                                |
                                |----- Interpolation kind ---- For the `RomanPSFInterpolator` in the input
                                                                you need to specify which `kind` of
                                                                interpolation to use. At the moment only
                                                                `corners` is availalbe.

(It can be found at the top of the psf.py file)

To reproduce the current behavior with this modification, the config file would look like this:

psf:
    type: RomanPSF
    interpolator:
        type: RomanPSFInterpolator

input:
    RomanPSFInterpolator:
        kind: corners
        n_waves: 5

On top of allowing for new interpolator I think this implementation also open the door to more parametrization of the interpolation. Now you can also use any PSF already existing in GalSim and also making fast simulation by not providing an interpolator but setting a wavelength in the RomanPSF, such as:

psf:
    type: RomanPSF
    wavelength: 1048

Technical note:
The implementation might look a bit convoluted by declaring a new ValueType to handle the interpolator, but I did that to be as future and safe proof as possible. A simpler alternative that could work rely on the fact that GalSim does not do a strict check of all the inputs. For example one can use a registered InputType even if the ObjectType is not registered with an input_type. Similarly, you can register an ObjectType with an input_type and not use it. If the input is not used in the function call, it will not create an error. This could simplify the implementation because we could remove the intermediate new ValueType. I don't like this idea because if GalSim implement more safty check in the future it could break this implementation. The way it is implemented here goes around this problem by doing 2 things:

• The RomanPSF ObjectType does not use an input_type directly
• The RomanPSFInterpolator InputType is pass through a ValueType RomanPSFInterpolator which you don't have to use. So the InputType is only needed if this variable is declared in the config file for the RomanPSF.

This is very technical and I am not sure my explanation is very clear here. I'll be happy to discuss it more. Please ask any questions about this and I would be very interested in any feedbacks.
(edit: I did most of the tests during the implementation using the python interface so I might be bypassing some checks that GalSim does when called through the CLI. So I am not sure that the alternative approach proposed above would work 🤔. The validation is done using the CLI so at least the implementation I did for this PR works.)

Final note:
I am really bad at giving names to things, feel free to suggest better ideas!

Validation: This implementation has been tested against an images generated using the config from this repo before modifications.

[arunkannawadi]

Some quick questions before I dig into the details - from the PR title, it appears that you're fixing a bug. But it's not necessarily a bug that you're fixing but you're extending it to support a new feature of having multiple SCAs in one file, right? Also, one SCA per file makes sense to me. This is also how romanisim seems to be doing it. What's the use case for storing multiple SCAs in one file?

[aguinot]

What's the use case for storing multiple SCAs in one file?

The direct interest I have in doing that is because of the way my simulation works. I simulate a small field with multiple observations from multiple SCAs and keep them in memory. This is how it looks in python:

# Read config
config = ReadYaml(config_path)[0]

# Import the required modules
galsim.config.ImportModules(config)

# Process the InputType (such as the PSF or SkyCatalog)
galsim.config.ProcessInput(config)

# Make the images
image = galsim.config.BuildImages(3, config)

This is consider as "1 file". To do this with the current implementation I would have to change the config file everytime I need to change the SCA.

But I want to emphasize on the fact that, in my opinion, the current implementation could lead to some issues. There is nothing that says that one should process only one SCA per file and per run. Even changing the nfiles in the config would not work. And the note in the config default.yaml is also misleading:

# We're just doing one SCA here.
# If you wanted to do all of them in each of three filters (given below), you could use:
#
# SCA:
#     type: Sequence
#     first: 1
#     last: 18
#     repeat: 3  # repeat each SCA num 3 times before moving on, for the 3 filters.

that would not work. You would only get the PSF from the first SCA and it would be very difficult to figure it out in the output.

There are other things in the implementation that could cause issues, such as querying the WCS. This is fine given how the simulation is used right now but with my simulation the WCS need the image size, which is only determined when the image is built, which happen after the inputs are set.

From a pure implementation perspective I think it is a bit dangerous to have methods that do more than what they are supose to, such as the getKwargs in the PSFLoader that read the WCS and bandpass. Or the PSF interpolation that is done in the init. Finally, galsim provide the perfect method for this case within the InputLoader called setupImage which would be perfectly adapted for this case (the SkyCatalog loader suffer from the same problem, I am working on a similar fix).