EoRImaging / eppsilon

eppsilon - error propagated power spectrum with interleaved observed noise
BSD 2-Clause "Simplified" License
5 stars 4 forks source link

rotate with lomb scargle, then fold #110

Closed nicholebarry closed 4 years ago

nicholebarry commented 5 years ago

Currently, the code calculates the Lomb-Scargle rotation phase for the frequency to folded kz transform.

If the distribution of the errors in the cos-sin basis depend on frequency, then perhaps the rotation phase should be calculated for the frequency to unfolded kz transform. Doing so assumes that the error distribution is independent in frequency. This doubles the number of rotation phases calculated.

I changed the rotation phase to calculated on the freq to unfolded kz. This required a bit of reordering; I had to transform, rotate, then fold rather than transform, fold, then rotate. The code in this pull request could be better...I tried to implement this with minimal changes, but that required overwriting variables in a clunky manner. Code also assumes the true z dft is used if even_freq is set (default).

Here is the difference between the original Lomb-Scargle and the updated Lomb-Scargle. There is improvement in the window. ps__master_compare_minus_rotate_then_fold

If I force the rotation phase to be 0, I have the standard Fourier Transform. When I compared the original Lomb-Scargle to the unrotated Fourier Transform, the original Lomb-Scargle did worse (can provide PS if wanted). Now, the updated Lomb-Scargle does better. Here is the difference between the unrotated Fourier Transform and the updated Lomb-Scargle. ps__temp_minus_rotate_then_fold

bhazelton commented 5 years ago

Actually, the current code calculated sine & cosine modes and then calculates the LS on them.

The sine and cosine modes are constructed using the postive and negative exponential modes to take advantage of the speed of the FFT, which calculates exponential modes (which go positive & negative). The cosine and sine modes are only positive, resulting in the same number of numbers (N complex modes go to N/2 cosine modes and N/2 sine modes). The positive/negative complex modes are only a computational convenience.

nicholebarry commented 5 years ago

Everything about the Lomb-Scargle is the same, except for whether or not the rotation is calculated for x number of frequencies -> x number of kz (positive and negative) or x number of frequencies -> x/2 number of kz (just positive).

I think I have proof that it should be x number -> x number. But the reason why is kinda surprising!

Here is the cosine of the rotation phase as a function of kz for an example pixel ~roughly within the EoR window. This is using the original Lomb-Scargle, so there are only x/2 rotations (originally 118 frequency channels, so 59 kz rotations).

Screen Shot 2019-07-02 at 1 10 29 pm

Here is the same data for the same pixel for the updated Lomb-Scargle. There are x rotations. However, even though I treated every index as an independent rotation in the code, there only seems to be x/2 independent rotations...it's symmetric about the center.

Screen Shot 2019-07-02 at 12 17 42 pm

Extra note: x title should be kz index, not freq index

So, do these give roughly the same PS? Not at all. The updated Lomb-Scargle does significantly better in the EoR window (see PS above). In addition, the original Lomb-Scargle did worse than a regular FT! Here is the original LS minus the regular FT. Blue indicates the regular FT did better... ps_master__compare_minus_stdft

So why are they different? I think it's because of the sin(theta) term. It's mirrored and flipped about the center.

This is using the original Lomb-Scargle, so there are only x/2 rotations (originally 118 frequency channels, so 59 kz rotations).

Screen Shot 2019-07-02 at 5 53 29 pm

Here is the updated Lomb-Scargle, so there are x rotations.

Screen Shot 2019-07-02 at 5 50 53 pm

This isn't a symmetric function, so I don't think we can apply just the positive half to folded data and expect the same answer. I would be fine if this the sin(theta) term was applied to data that we flipped + folded, but that's not the case.

This is my going theory. What do you all think? Does it also matter that the FT of the data isn't perfectly symmetric as well?

miguelfmorales commented 5 years ago

I think to sort this out, we actually need to go back to the error covariances. There is a very specific reason behind the LS approach: diagonalizing the error covariance between the two terms on each k mode of a complex FT. (Bryna thinks sine and cosine, I think real and imaginary.)

So the question is whether we are correctly calculating the diagonalization, and whether we have enough rotation terms. I think this is a good question that we'll have to think about. But I'm not swayed by the PS, this is really a statistics question about the covariance of the error bars and I think it is more fruitful to think about it that way.

bhazelton commented 4 years ago

closing in favor of #115