def smoother_step(F, filt, next_pred, next_smth):
    """Smoothing step of Kalman filter/smoother.

    Parameters
    ----------
    F:  (dx, dx) numpy array
        Mean of X_t | X_{t-1} is F * X_{t-1} 
    filt: MeanAndCov object
        filtering distribution at time t
    next_pred: MeanAndCov object
        predictive distribution at time t+1
    next_smth: MeanAndCov object
        smoothing distribution at time t+1

    Returns
    -------
    smth: MeanAndCov object
        smoothing distribution at time t 
    """
    J = dotdot(filt.cov, F.T, inv(next_pred.cov))
    smth_cov = filt.cov + dotdot(J, next_smth.cov - next_pred.cov, J.T)
    smth_mean = filt.mean + np.matmul(next_smth.mean - next_pred.mean, J.T)
    return MeanAndCov(mean=smth_mean, cov=smth_cov)