burst_model can fail in non-sorted observations

[Kamuish]

I was trying to load a GP model from disk, using an unordered dataset

>>> self.data.obs
array([2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
       2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
       2, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
       1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
       1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
       1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3,
       3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3,
       3, 3, 3])

And I am running into the following error when calling the following code:

>>> p = res.maximum_likelihood_sample()
>>>    res.phase_plot(p)

Traceback (most recent call last):
  File "/home/amiguel/RV_analysis/EpsIndi/kima_runs/plot_kima.py", line 69, in <module>
  File "/home/amiguel/RV_analysis/tools/kima/src/kima/pykima/display.py", line 2611, in phase_plot
    no_planets_model = res.burst_model(sample, tt, no_planets_model)
  File "/home/amiguel/RV_analysis/tools/kima/src/kima/pykima/results.py", line 2709, in burst_model
  File "/home/amiguel/RV_analysis/venv/lib64/python3.13/site-packages/numpy/_core/fromnumeric.py", line 207, in take
    return _wrapfunc(a, 'take', indices, axis=axis, out=out, mode=mode)
  File "/home/amiguel/RV_analysis/venv/lib64/python3.13/site-packages/numpy/_core/fromnumeric.py", line 57, in _wrapfunc
    return bound(*args, **kwds)
IndexError: index 4 is out of bounds for axis 0 with size 3

I have managed to trace back the source of the issue to the following lines of code (which you dubbed as an "Hack"):

https://github.com/kima-org/kima/blob/f5af4467775c187837020931db09e36a0fb88001/src/kima/pykima/results.py#L2703C1-L2706C25

I have no idea how to fix this, so I will leave it up to you (and ensure that I use ordered datasets). However, I noticed that this does not occur for the RV_model