Our calculation of Q and P matrices of size 2 x 2 is based on a Hamiltonian function corresponding to traveltime as the variable along the ray (see, e.g., Cerveny, 2001). For isotropic media, the differential equations that we use are the following: 

dQ/dt = v^2 P 
dP/dt = -v^{-1}VQ 

ds/dt = v 

with the initializations
Q(0)=0 
P(0)=I . 

The matrices Q and P refer to a local ray-centered system which is continuously rotated along the ray (Popov and Psencik, 1979). This means that we also need an additional differential equation describing the change of rotation of the basis vectors of the system. For details, see Cerveny's books from 1985, 2001 or his earlier work on the subject. 

In the equations above, v is the velocity in the current point on the ray, and V is a 2 x 2 matrix containing second derivatives of velocity with respect to local ray-centered coordinates. The "0" and "I" denote, respectively, the zero matrix and identy matrix of size 2 x 2.