Martin's notes

Extending axsim to wavelength dependent template direct images

axesim creates from all images in the "image template list" a so-called Model Image File.
This is done in the axesim task 'prepimages'.
Currently there is only one image given for all wavelengths in the image template list, and 'prepimages' just stacks these images together into the Model Image File, which is a a multi-extension fits file.

The modelling of the beams is done in the executable 'aXe_PETCONT' (source code in 'aXe_PETCONT.c'), there the Model Image File is read into the structure 'object_models' ('specmodel_utils.h' l.60).
This structure consists of a list of 'dirim_emission' objects ('specmodel_utils.h' l.52), with one dirim_emission-object for each extension in the Model Image File.
As you can easily see, a dirim_emission-object consists mainly of one 'gsl_matrix' ('specmodel_utils.h' l.50), which is a simple matrix that stores the image values.

The routines to read in an entire Model Image File is "load_object_models" ('specmodel_utils.h' l.63).
The code is in the file 'specmodel_utils.c'.

To evaluate the brightness of a 'dirim_emission'-object at a given position (x,y) the routine 'get_diremission_value()' ('specmodel_utils.h' l.90) is used.
In the modelling of slitless beams, this function is used in the subroutine 'make_gauss_spectra2' (see "spc_model.c, l.391").
In the modelling direct images, this function is used in the subroutine "make_dirimage()" ('dirimage_model.c' l.192.).

Extending axsim for wavelength dependent template images has two aspects.
The first aspect is computational, which is how the wavelength dependency is take into account when modelling sliltess spectra and direct images.
The second aspect is organisational, and covers the way the user can specify in the input that there are now several images per direct image model, how they are organized and transported, how the user can address a specific direct image model in the Model Object Table, how each model image know which wavelength it corresponds to and so on.
For the organisation, there are many solutions possible and I will not cover this here (the programmer has to find his/her own way), and I will cover only the computational aspect here.

So I could do the following:
  • Extend the structure 'dirim_emission' such that it accepts several images ('gsl_matrix * modimage' --> 'gsl_matrix * * modimage').
    I would add an integer storing the number of matrixes in ' * * modimage' and a 'gsl_vector' that carries the wavelength for each matrix in ' * * modimage'.
    The matrixes should be sorted in wavelength to easily pick out the two (or one) relevant for any specific wavelength.
    This requires changing the subroutines 'load_object_models', 'load_diremission_list' and 'load_diremission' according to your organisational structure (don't forget also the corresponding 'free_ * ' routines.
  • Change the subroutine 'get_diremission_value()' such that it spits out the value for a wavelength given as input.
    This is not that compliacated, you just have to select the two matrixes sandwiching the wavelength, intepolate each to get the pixel value at the right (x,y) position and then interpolate between these values to get the value at the correct wavelength (however see open questions below).
  • Change the subroutin 'make_gauss_spectra2' around line 391 such that it uses the new, wavelength dependent form of 'get_diremission_value()'.
    Then it should iterate over the wavelength given in the structure 'tracedata' (parameter 'acttrace' and fill into 'acttrace->gvalue' NOT a constant value but the value corresponding to the wavelength.
  • For modelling direct images, change the subroutine 'make_dirimage()' (file 'dimimage_model.c') such that it uses the new, wavelength dependent form of 'get_diremission_value()' (around line 192).
    I am not totally sure how exactly this should be done (see open questions below).

Most of the code that has to be touched is comparably well written and documented, thus it should not be too complicted to get through.

Open questions:

  • The matrixes used for representing a direct image model are normalized (the sum is 1.0), which is done in 'prepimages'.
    The normalization is important for the flux conservation.
    However if you interpolate for one (x,y) position between two images specified at different wavelength.
    I am not sure whether the interpolated 'meta'-image is still normalized.
    My stomach says no, it is no longer normalized, but the differences are small and can be neglected.
    But my stomach can be wrong.
  • When modelling a direct image, the total intensity value (in e/s==cps) for an objcect is determined in the routine 'get_cps_for_dirobject()' ('dirimage_model.c', l.155) by integrating the SED over the system throughput and then distributing this total value according to the spatial distribution (double for-loop starting in l.171), which is either a Gaussian distribution or a direct image.
    I am not sure how to insert a wavelength dependency of a direct image there.
    Either it is sufficient to also integrate over the system throughput and thus generating an 'effective image' for the given passband or the wavelength dependency of the direct images have to be taken into account when computing the total intensity value in 'get_cps_for_dirobject()'.
    One should write down the equations and discuss....