Pixel scaling in healpix_interpolate

[aelanman]

If from_jansky is set in healpix_interpolate, then the convert_kelvin_jansky function is run to get the conversion scalar. https://github.com/EoRImaging/FHD/blob/61bca89ff15b0fa82347b0f15991cca08db69fa1/fhd_core/HEALPix/healpix_interpolate.pro#L51 https://github.com/EoRImaging/FHD/blob/61bca89ff15b0fa82347b0f15991cca08db69fa1/fhd_utils/format_conversion/convert_kelvin_jansky.pro#L5-L10

This gives a conversion factor from Kelvin to Jy/Ω, where Ω is the input pixel area = 4π/(Nside^2 * 12). The input map is interpolated to the instrument's sky coordinates, and then scaled by this factor to get map_interp.

The last steps of healpix_interpolate then apply a factor of pixel_area: https://github.com/EoRImaging/FHD/blob/61bca89ff15b0fa82347b0f15991cca08db69fa1/fhd_core/HEALPix/healpix_interpolate.pro#L117-L120

This factor is calculated from the pixel edges, and represents the projected area of each healpix pixel (call this ω): https://github.com/EoRImaging/FHD/blob/61bca89ff15b0fa82347b0f15991cca08db69fa1/fhd_utils/modified_astro/pixel_area.pro#L42-L49

This factor translates from Jy/ω to Jy/(degpix)^2. Prior to this, map_interp is in Jy/Ω. So this means that we're off by a factor of ω/Ω for each pixel.

I'm submitting this issue for discussion, but this might not be a problem. Maybe the interpolation step effectively scales the fluxes to make up the difference? However, I have seen improved agreement among maps of different resolutions by replacing the conversion factor in pixel_area.pro.

(Yes, I finally learned to embed code in an issue!)

[isullivan]

Is your correction that you replace IF Keyword_Set(relative) THEN area_map/=(obs.degpix*!DtoR)^2. with IF Keyword_Set(relative) THEN area_map/=4.*!Pi/nside2npix(nside)?

[aelanman]

My correction was:

area_map *= 0.0
area_map += (obs.degpix*!DtoR)^2 / 4.*!Pi/nside2npix(nside)

So I completely skipped the projected pixel areas and converted direction from Jy/(healpix_pixel_area) to Jy/(degpix * !DtoR)^2.

[aelanman]

For comparison: Here's the 1D power spectrum comparison without my correction:

and with:

(As usual – this is with MWA-128 layout, MWA beams, gaussian Healpix shells of mean 0 and variance 1, interpreted as Kelvin. Each curve is a different NSIDE for the input shells)

[aelanman]

Some more plots. For these, I put the pixel area scaling correction back into pixel_area.pro and repeated simulations for upgraded-resolution gaussian maps (ugrade) and smoothed gaussian maps (obtained by expanding in spherical harmonic basis up to lmax=3*nside-1 at each frequency, then transforming back, effectively convolving with the pixel window function). The second operation changes the variance of the maps a little:

NSIDE	σ^2
256	2.93692171573
512	2.93660915602
1024	2.93218540872

Comparison power spectra for smoothed (linear y scale):

Comparison power spectra for upgraded maps:

These input maps have identical variance (2).

[aelanman]

2D power spectra, first from eppsilon's plotter, then from a python plotter I've been working on:

Dimensional power spectrum (in mK^2 Mpc^3)

Dimensionless (mK^2). All are xx polarization.

[bhazelton]

@aelanman eppsilon will plot in linear bins if you'd like. You just need to set log_kperp=0 and log_kpar=0. It's really not reasonable to plot all the way down to kperp=0 because that's below the antenna size.

[bhazelton]

You can also adjust the color bar range by setting the data_range keyword. For simulations, I usually set /sim, which bypasses the standard color bar ranges which are data-motivated, not simulation motivated.

[aelanman]

I fixed the colorbar ranges (thanks Bryna), and got this.

I'm not sure... but it looks like the coarse band structure might still be there?

[bhazelton]

I don't think the coarse band structure is there. I could be wrong, but I think they don't fall where the vaguely linear things you're seeing are and I think it's just that your eye is good at seeing patterns. The lowest coarse band is generally at about kpar=.35. Also notice how small the range of the color bar is -- the fluctuations we're looking at here are pretty small, at least below 50 lambda.

[aelanman]

Ratio plots!

256/512: 512/1024: 256/1024:

[aelanman]

With a ridiculously compressed colorbar (0.9 to 1.1):

[aelanman]

Here are the ratio plots of the original results with the upgraded maps, without my pixel area correction:

256/512:

512/1024:

256/1024:

I've tried to adjust the color scale to best show the differences, but it's non-negligible. I'm going to look at through pixel_area.pro again and see how its area calculation compares with the healpix pixel areas.

[aelanman]

Linear color scales on the 2d power spectra ratios:

[aelanman]

For nside=1024, this is the variance vs baseline length for channel 192. The vertical lines indicate the resolution of the input map (red) and the orthoslant map (green). Note that these are for the "upgraded" maps with my pixel area fix in place (see eor_bubble_fixes branch).

The same plot, but for nside=512:

I've added a multipole taper in eor_bubble_sim which downweights higher order multipoles with a [1-tanh(l-l_0)] function, in the following way:

1. Calculate lmax for the orthoslant map (as 180/degpix) and lmax of the input map (as 3*nside-1).
2. If lmax_ortho < lmax_healpix, then this means the input map has higher resolution than the orthoslant map. Smooth the input map by downweighting l>lmax_ortho using the tanh function.
3. If lmax_ortho > lmax_healpix (this is usually the case), then the orthoslant map has higher resolution. Smooth the healpix map to its maximum resolution (lmax_healpix), weighting the highest lmodes with the tanh function for a smooth cutoff.
4. Rescale the map to preserve variance.

I ran simulations for nside=512 and nside=1024, and this seems to have removed the angular dependence for low kperp modes:

[aelanman]

Lastly, I tried varying the orthoslant resolution by changing the dimension parameter, keeping kbinsize fixed (0.5). The resolution is 1/(kbinsize * dimension), so higher dimension = higher resolution. Repeating simulations with an NSIDE=1024 shell, the 1D power spectra are:

And for nside=512: