nanograv / enterprise

ENTERPRISE (Enhanced Numerical Toolbox Enabling a Robust PulsaR Inference SuitE) is a pulsar timing analysis code, aimed at noise analysis, gravitational-wave searches, and timing model analysis.
https://enterprise.readthedocs.io
MIT License
64 stars 65 forks source link

ECORR incorrectly only applied when `if slc.stop - slc.start > 1` in an epoch #355

Open vhaasteren opened 11 months ago

vhaasteren commented 11 months ago

In ShermanMorrison (and BlockMatrix?), the ECORR signal is only applied to epochs where there is more than 1 observation. Although the reasoning seems to make sense that this would otherwise just be an EQUAD at that epoch, this is not a correct way to restrict application. The ECORR signal in our generative model affects an observation whether or not there are more observations in that epoch. It always affects the observation, whether or not it is distinguishable from an EQUAD or not.

So, for instance, in the below code we should change the if condition to > 0, so allow for single observations.

    def _get_logdet(self):
        """Returns log determinant of :math:`N+UJU^{T}` where :math:`U`
        is a quantization matrix.
        """
        logdet = np.einsum("i->", np.log(self._nvec))
        for slc, jv in zip(self._slices, self._jvec):
            if slc.stop - slc.start > 1:
                niblock = 1 / self._nvec[slc]
                beta = 1.0 / (np.einsum("i->", niblock) + 1.0 / jv)
                logdet += np.log(jv) - np.log(beta)
        return logdet