sambit-giri / tools21cm

MIT License
22 stars 18 forks source link

angular resolution when computing uv coverage #3

Closed dprelogo closed 4 years ago

dprelogo commented 4 years ago

In function telescope_functions.get_uv_coverage, could you explain why is Nb = Nbase*theta_max/2 and not only Nbase*theta_max?

I'm using your code for applying thermal noise to the 21cmfast lightcones, found this as final-checking, not sure if it is a mistake or I don't see something...

Thanks :)

sambit-giri commented 4 years ago

Sorry for such a late reply. I must have missed the notification for this issue.

Nbase*theta_max/2 has no physical meaning. This 1/2 factor cancels out in the condition statements (e.g. Nb[:,0]<ncells/2) following that part of the code. In the current form of the code to calculate the uv coverage, the baselines were counted twice when this 1/2 factor is removed from both Nbase*theta_max/2 and Nb[:,0]<ncells/2.

dprelogo commented 4 years ago

Thank you for a reply!

In the meantime I made a simple equivalent function, you can look at it here. I'm just creating uv frequency map with np.fft.fftfreq and then creating a histogram of baselines with np.histogram2d. The overall effect is that I get zoomed-in uv coverage by the factor of two in each dimension, which should be consistent with this 1/2.

As you said, factors cancel out, so if multiply by theta_max/2 and take conditions in range <-ncells/2, ncells/2>, that's effectively <-1/theta_max, 1/theta_max> so 2/theta_max interval in total in fourier space. If I multiply with theta_max, and take the same range, that will sum to 1/theta_max.

For the double counting, I think the final uv_map can be divided by two, but I can't see where this comes from.

sambit-giri commented 4 years ago

This factor of 2 comes from the way one takes care of mirror baselines.

steven-murray commented 4 years ago

Just to weigh in on this:

The size of a cell in fourier space (in inverse radians) should simply be D_c(z)/boxsize. In inverse Mpc, its 2pi/boxsize. Placing a cut-off at B (metres) means placing a cutoff at u = B/lambda. This should be at the n=u*boxsize / D_c(z) pixel from the centre.

So, the lines

theta_max = boxsize/cm.z_to_cdist(z)
Nb  = np.round(Nbase*theta_max/2)

are getting the baselines, Nb, almost in units of their cell-number, except for that divide-by-2.

Say we have a simulation box of size 100 Mpc, with 100 cells. Then the size of a Fourier-cell (in inverse radians) is D_c / 100, and we have 100 cells, from -D_c/2 to D_c/2.

We should expect a baseline with u=D_c/2 to be included in the last cell. In the code, we'd get Nb = u * 50 / D_c = 25 for this baseline. This passes all the cuts (which are that it has to be <50 and >-50. But it really should be in the 100th cell, not the 25th cell!

This doesn't seem correct, and I'm not sure it has much to do with the mirror baselines?

sambit-giri commented 4 years ago

Thanks @steven-murray and @dprelogo for pointing out the mistake. I have modified the code. Now the map should match with the one @dprelogo gets with their version of telescope_functions.get_uv_coverage (when mirror baselines are not included).

steven-murray commented 4 years ago

Awesome!