I'm trying to compute derivative of end effector positions $\mathbf{p}\in \mathbb{R}^{3}$ with respect to quaterion $\mathbf{q}\in \mathbb{R}^{4}$. When I use euler angles for joint rotations, the derivatives can be fairly easily computed using cross product (https://cseweb.ucsd.edu/classes/wi17/cse169-a/slides/CSE169_08.pdf). However, I want to use quaternion, and I'm struggling to compute Jacobian $\mathbf{J} =\frac{\partial \mathbf{p}}{\partial \mathbf{q}} $. Do I have to compute something like below for quaternion based Inverse Kinematics ?:
$$ \mathbf{J} = \begin{bmatrix} \frac{\partial p_{x}}{\partial q_{x}} & \frac{\partial p_{x}}{\partial q_{y}} &\frac{\partial p_{x}}{\partial q_{z}} & \frac{\partial p_{x}}{\partial q_{w}} \\ \frac{\partial p_{y}}{\partial q_{x}} & \frac{\partial p_{y}}{\partial q_{y}} &\frac{\partial p_{y}}{\partial q_{z}} & \frac{\partial p_{x}}{\partial q_{w}} \\ \frac{\partial p_{z}}{\partial q_{x}} & \frac{\partial p_{z}}{\partial q_{y}} &\frac{\partial p_{z}}{\partial q_{z}} & \frac{\partial p_{x}}{\partial q_{w}} \end{bmatrix}\in \mathbb{R}^{3 \times 4} $$
If so, please help me out to compute those values.
Thanks in advance!