Reference: Gaussian Modeling with MIRIAD maths

A note to self. Sometimes I want to use maths to generate an image with a large Gaussian component. Let’s say that the total flux I want to capture is F (in Janskys). Using what we know about normalized Gaussians the expression I want is something of the form

$\frac{F}{2\pi\sigma_x\sigma_y} \exp\left(-\frac{1}{2}\left[\frac{x^2}{\sigma_x^2}+\frac{y^2}{\sigma_y^2}\right]\right)$

. Assuming circularity and simplifying, we get

$A\exp\left(B\left[x^2+y^2\right]\right)$

, where

$A = \frac{F}{2\pi\sigma^2}$

and

$B=-\frac{1}{2\sigma^2}$

.

For now I’ve thought through the implementation using arcsecond units. Say we have a 2047×2047 image with a pixel scale of 10 arcsec that’s uniform across the image. If I want σ to measure the ATA primary beam, with a FWHM of 3.5°/GHz, I get σ = 5350.7 arcsec/GHz. To generate a model-type image, we need to work in Jy/pixel, so A must be scaled by the arcsec²/pixel conversion ratio, 100 in this case.

Implementing this in maths is straightforward. You use the xrange and yrange keywords to implement the Gaussian expression, ranging them from -10230 to +10230 in the example here.