Add flux-preserving mode to reproject

[cdeil]

Continuing the discussion from https://github.com/astrofrog/reproject/pull/67#issuecomment-98984750 and #3, I think it would be nice if the high-level reproject API contained a flag conserve_flux (or something like that) that's False by default and when set does this ("sb" = "surface brightness):

solid_angle_image = compute_solid_angle_image(flux_image)
sb_image = flux_image / solid_angle_image
sb_image2 = reproject(sb_image)
solid_angle_image2 = compute_solid_angle_image(sb_image2)
flux_image2 = solid_angle_image2 * sb_image2

@astrofrog – Do you think this is a common enough usecase to expose in the high-level API functions or is there some better way?

[cdeil]

I'm not sure this is the best way to implement this. In the reproject_exact method there might be a way to implement the flux-conserving mode without the separate extra computation of the solid angle images for the input and output image?

Since for HEALPIX there is no exact mode and even for normal images the interpolating reprojection is much faster, it might still be useful to have this method available?

[mwcraig]

@astrofrog @cdeil -- what is the status of flux conservation in reproject? We're getting ready to add a convenience function in ccdproc for combining images by reprojecting and I want to make sure I understand the current state of affairs.

[cdeil]

what is the status of flux conservation in reproject?

It's not available yet. Someone would have to implement a function to compute pixel area (i.e. solid angle) images first. There's different methods ranging from a few lines to complex code that result in difference levels of accuracy. Maybe we should just start by putting a function that assumes pixel areas are given by CDELT and are constant across the image, and leave the other methods for later, if / when someone has time? I think the API could remain unchanged: pixel_area(image, method='..').

[astrofrog]

@cdeil is referring to the fact we currently preserve surface brightness but not technically flux - just clarifying this since some people sometimes refer to flux-conserving more as an indication that the total flux of the image stays constant. But indeed someone still needs to implement the strictly flux-conserving mode.

[mwcraig]

I think I understand...so if I have an image in ADUs or photon counts I need to divide by pixel area, reproject, then multiply by pixel to get back to counts?

I know in some of the PRs @larrybradley has done for NDData he assumed the pixel area was constant across the image and used the area of the center pixel. I'll start there and see where it takes me.

[astrofrog]

@mwcraig - yes, that's correct. As @cdeil mentioned, we should have a function in reproject or probably actually in astropy core that can provide the exact area in each pixel.

[Cadair]

I think from a solar perspective flux is likely to be what people want, just ftr.