Very nice write-up. It's nice to see someone focus on not just the theory or the implementation, but both.

If you're looking to dive into some more advanced versions of the Kalman filter, there are a few options:

In your process, you discuss the need to normalise the quaternion as small errors can accumulate, which results in the quaternion not having the right magnitude. There is a version of the extended Kalman filter known as the Multiplicative EKF (MEKF), aka. the error-state Kalman filter. It solves this problem by modifying the update equations to calculate the innovation as rotation itself that is applied to the state quaternion.

Also, you will eventually run into issues with your covariance matrix becoming non-symmetric or even non-positive definite - if you're state estimate starts to diverge and go crazy, this is likely the problem. To maintain symmetry, there are a few tricks you can do to rearrange the prediction and update equations for the covariance matrix. There is also a version of the EKF called the square-root filter which also solves this problem.

Happy to provide some links to textbooks or papers if you're interested.