Archive for September, 2009
23 Sep 2009
Data volumes from multiple sky surveys have grown from gigabytes into terabytes during the past decade, and will grow from terabytes into tens (or hundreds) of petabytes in the next decade. ... For astronomy to effectively cope with and reap the maximum scientific return from existing and future large sky surveys, facilities, and data-producing projects, we need our own information science specialists. We therefore recommend the formal creation, recognition, and support of a major new discipline, which we call Astroinformatics. ... Now is the time for the recognition of Astroinformatics as an essential methodology of astronomical research. The future of astronomy depends on it.
4 Sep 2009
They're coming thick and fast now.
Here's a Python function to accompany the previous post. It's not maximally efficient, but should make sense...
from scipy import stats def gaussian_pixel(minxy, maxxy, sigma, meanxy=(0.,0.), norm=None): """Return the value of a pixel sampling a 2D Gaussian, normalized such that the area under the Gaussian is 1 (default) or such that the peak is given by norm.""" x1, y1 = minxy x2, y2 = maxxy x0, y0 = meanxy if norm is None: norm = 1. / 2 / math.pi / sigma ** 2 return norm * 2 * math.pi * sigma ** 2 / (x2 - x1) / (y2 - y1) * ( (1 - stats.erfc((x2 - x0) / math.sqrt(2) / sigma)) / 2. - (1 - stats.erfc((x1 - x0) / math.sqrt(2) / sigma)) / 2.) * ( (1 - stats.erfc((y2 - y0) / math.sqrt(2) / sigma)) / 2. - (1 - stats.erfc((y1 - y0) / math.sqrt(2) / sigma)) / 2.)
4 Sep 2009
Yesterday I said that the PRF for a map in Jy/beam (or similar) should be normalized so that that peak is 1. But this is true only for an idealised (not pixelated) PRF, or if the map has infinitesimally small pixels.
If the pixels are larger than infinitesimal, as is generally the case, then the maximum value of the pixelated PRF will be the average value over the pixel, which will be less than 1.
For example, if the PRF is a two-dimensional Gaussian, centred on , with standard deviation , then the value in a pixel with and will be
Ugh. Let's make that simpler. For a PRF centred on , and a pixel , this is
As an example, the fairly-Gaussian beam for the Herschel Space Observatory SPIRE instrument has an FWHM of around 18", which corresponds to a standard deviation of around . If we make a Jy/beam map with pixel size 6", then the peak value for a 1 Jy point source in the centre of a pixel will be
No big deal really...
3 Sep 2009
You have a map and you know what a point source looks like. How do you filter the map so that the value of each pixel is now the most likely flux of a point source centred on that pixel? (An isolated point source, to be more precise.)
First, find , which is the point response function (PRF), telling you what a point source of flux 1 will look like in the map. This may be normalized so that the peak is 1 (if your map is in Jy/beam or similar), or so that (if your map is in Jy/pixel or similar). If your map is in MJy/sr ... well, figure it out and add a comment below. Basically, if you normalize your PRF correctly, you won't need to worry about the map units in what follows. Phew.
Now the measured value of each pixel around the point source, , will be
where is the flux of the source and is the noise, drawn from a normal distribution with mean zero and standard deviation .
Now the badness of the fit is measured by the , which is given by
At the maximum likelihood value of the flux, , will be at a minimum, so
Solving this for , we find the maximum likelihood solution
Now just do this for each pixel in the map (corresponding to a point source centred on each pixel) and you're done.
Worked example. is 0.5, 1.0 and 0.5, for three adjacent pixels (you'll have realised that the map is in Jy/beam or similar), and is 1, 2 and 1 Jy/beam, for three adjacent pixels, with the same (tiny!) value of for each pixel (in this case, we can ignore the value of in what follows). So the flux at the central pixel is estimated to be
which is no surprise, since the maximum value of the map in Jy/beam is 2 at that position.
This is an example of a matched filter (I haven't read the page, but hopefully including the link will make me look clever). And, given that point sources are under no particular obligation to align themselves with the centres of the pixels of your map, can easily be re-estimated for a source with a certain offset from the pixel centre.