I believe there is another problem in the code, LLtwoCIR, in the calculation of the Q matrix, the (1,1) element of which is,
(theta1*sigma1*sigma1*(1-exp(-kappa1*dt))^2/(2*kappa1)+sigma1*sigma1/kappa1*(exp(-kappa1*dt)-exp(-2*kappa1*dt)))*AdjS(1)
I think the correct code is,
theta1*sigma1*sigma1*(1-exp(-kappa1*dt))^2/(2*kappa1)+sigma1*sigma1/kappa1*(exp(-kappa1*dt)-exp(-2*kappa1*dt))*AdjS(1)
The AdjS only multiplies the second term, not both the first and second, see Chen and Scott (2003), page 147, second formula.
Thanks
I think the correct code is,
theta1*sigma1*sigma1*(1-exp(-kappa1*dt))^2/(2*kappa1)+sigma1*sigma1/kappa1*(exp(-kappa1*dt)-exp(-2*kappa1*dt))*AdjS(1)
The AdjS only multiplies the second term, not both the first and second, see Chen and Scott (2003), page 147, second formula.
I like the simple code, however I donot understand the following code segment:
VarS=VarS*(1-KalmanGain*H);
Given the KalmanGain definition:
KalmanGain=VarS*H'*InvVarY;
The above VarS is not symmetric, and (1-KalmanGain*H) is not correct either. I believe the corrrect calculation should be:
VarS=(eye(2)-KalmanGain*H)*VarS;
Very clean example of KF, but not general enough to deal with state vectoc. For example, s.x = inv(s.H)*s.z; and s.P = inv(s.H)*s.R*inv(s.H') would not work if number of state and number of measurement are not the same.
