Add a routine to compute the solid angle of all pixels in a FITS image

[astrofrog]

This would need to take a HDU object since WCS doesn't technically contain the dimensions of the image, only the projection information. I think astropy.wcs.WCS does contain the naxis info, but I've been told that this could be removed at any point since WCS doesn't say anything about the dimensions of the image.

[astrofrog]

Actually I think it should take a Header object.

[cdeil]

I started working on this (doesn't run yet and I hadn't seen your comments that it should take a Header as input): https://github.com/cdeil/python-reprojection/commit/5aeb598479410414a267489277b9eb8f0abfa679 https://github.com/cdeil/python-reprojection/commit/a23d1d3a5059796534dc4914721e8e61802d4b23

What do you think about the C function and Cython wrapper API though? I think having 2D lon and lat arrays of shape (naxis1 + 1, naxis2 + 2) of pixel corner positions like I implemented it now is most memory-efficient? Creating arrays of shape 4 * ((naxis1 + 1) and 4 * (naxis2 + 1)) has a memory overhead factor of 4 and would probably be slower?

[astrofrog]

Right, so you mean the difference is between having the four vertices listed for each pixel, which would mean any pixel is listed 4 times (except at the edge) or being smart and passing all the vertices only once to the C code, which will then have to do a 2-d loop? The second would indeed be more memory efficient, so I agree that's best.

[cdeil]

Yes, a 2D loop in C is what I had in mind, probably in combination with an inlined function to compute the solid angle that takes the four pixel corner coordinates as 8 input arguments.

I'll probably have some time to try this tonight.

[astrofrog]

Great, thanks for working on this!