Bugfix for issue- CrossCorrelation Gives Strange Results (#592) #647
Conversation
Obtained time_lags in using scipy.signal.correlation_lags tool.
Changes in docstrings
|
Hello @mihirtripathi97! Thanks for updating this PR. We checked the lines you've touched for PEP 8 issues, and found:
Comment last updated at 2022-03-18 13:06:30 UTC |
|
@mihirtripathi97 thanks! |
Added test cases for finding stingray.crosscorrelation.time_shift for identical light curves with odd and even number of data points.
|
@matteobachetti As you had asked, I have added two new test cases in stingray/tests/test_crosscorrelation.py. They check if the timeshift corresponding to maximum cross-correlation between two identical light curves (having odd/even number of data points) is zero. |
There was a problem hiding this comment.
Thanks @mihirtripathi97! One remaining request, then it's good to go.
stingray/crosscorrelation.py
Outdated
| @@ -193,9 +194,15 @@ def _make_corr(self, lc1, lc2): | |||
| # Calculates cross-correlation of two lightcurves | |||
| self.corr = \ | |||
| signal.correlate(lc1_counts, lc2_counts, self.mode) | |||
| # Obtains correlation lags | |||
| # signal.correlation_lags() method uses SciPy versions >= 1.6.1 | |||
| x_lags = \ | |||
There was a problem hiding this comment.
Please note that we used the internal method cal_timeshift here, that remains buggy despite your fix. It would be better to move this operation to that method, so that it's solved everywhere.
There was a problem hiding this comment.
@matteobachetti
I have a little bit of confusion here. The function scipy.signal.correlation_lags() takes lc1_counts.size and lc2_counts.size as input parameters to calculate correlation lags (x_lags).
The method cal_timeshift() is also called internally by method _make_cross_corr(), which is for cross-spectrum and does not compute cross-correlation using scipy.signal.correlate(). So if I move the operations to cal_timeshift() method then the calculation of timeshift for _make_cross_corr() will also get affected. Should I still move the operations to cal_timeshift() ?
There was a problem hiding this comment.
The number you can input for both light curve sizes is self.n in that case
There was a problem hiding this comment.
Thanks, @matteobachetti, I have shifted the operations to the internal method cal_timeshift as you suggested. Let me know what would be the next step.
There was a problem hiding this comment.
Understanding why the test fails :).
Moved operations to find correlation lags to internal method cal_timeshift()
| dur = int(self.n / 2) | ||
| # Correlation against all possible lags, positive as well as negative lags are stored | ||
| x_lags = np.linspace(-dur, dur, self.n) | ||
|
|
There was a problem hiding this comment.
@matteobachetti As you suggested, I have moved these operations to cal_timeshift internal method.
There was a problem hiding this comment.
@matteobachetti
All tests are failing for one test case - test_cross_correlation_with_unequal_lc() in test_crosscorrelation.py. The reason is - time_lags obtained by scipy.signal.correlation_lag() are diffrent ([-1, 0, 1, 2, 3]) from of time_lags given in the test case ([-2, -1, 0, 1, 2]). Since the crosscorealtion is calculated using scipy.signal.correlate() it makes more sense to change the time_lags in the test case accordingly. Should I make the changes in the test case?
Codecov Report
@@ Coverage Diff @@
## main #647 +/- ##
==========================================
- Coverage 97.21% 97.12% -0.09%
==========================================
Files 41 41
Lines 7535 7546 +11
==========================================
+ Hits 7325 7329 +4
- Misses 210 217 +7
📣 Codecov can now indicate which changes are the most critical in Pull Requests. Learn more |
|
Ok. Thanks @mihirtripathi97! |
The reason why the previous code did not provide zero time_shift for cross-correlation between two identical light curves of even length is as follows:
The correlation lags (x_legs) required to calculate time_shift was obtained by creating a linear space in cal_timeshift method.
The following line stores possible correlation lags in x_lags.
x_lags = np.linspace(-dur, dur, self.n)The correlation lags obtained using this method depend on the evenness/oddity of self.n (i.e. number of time bins in the light curve). For an odd number of time bins, np.linspace, will provide an ndarray of integers centered at 0, but for an even number of time bins, it will generate an ndarray of floats which will not include 0. As a result, the time_shift will never be 0.
The correlation lags (x_legs) required to calculate time_shift can be obtained using signal.correlation_lags.
`
`
The signal.correlation_lags takes care of evenness/oddity of input light curves for all possible correlation modes (i.e 'same', 'valid' or 'full') . This method was introduced in SciPy version 1.6.1.
closes #592