$f$ is a function such that $f:\mathbb{R}^n\to\mathbb{R}$, and function $\boldsymbol{g}:\mathbb{R}^n\to\mathbb{R}^n:\delta\mapsto\boldsymbol{q}\cos\|\boldsymbol{\delta}\| + \frac{\boldsymbol{\delta}}{\|\boldsymbol{\delta}\|}\sin\|\boldsymbol{\delta}\|$, where the constant vector $\boldsymbol{q}\in\mathbb{R}^n$. Here obviously $\boldsymbol{g}(\boldsymbol{0}) = \boldsymbol{q}$.
Assume all involved derivatives exist. In a small neighbourhood of $\boldsymbol{0}$, consider the Taylor expansion to quadratic term of the function $$f\big(\boldsymbol{g}(\boldsymbol{\delta})\big) = f\bigg(\boldsymbol{q}\cos\|\boldsymbol{\delta}\| + \frac{\boldsymbol{\delta}}{\|\boldsymbol{\delta}\|}\sin\|\boldsymbol{\delta}\| \bigg)$$.
We use the notation $\nabla f|_{\boldsymbol{q}}$ to mean gradient of $f$ evaluated at the point $\boldsymbol{q}$, and similar for Hessian $\nabla^2 f|_\boldsymbol{q}$.
How do I see the following: $$f\big(\boldsymbol{g}(\boldsymbol{0}+\boldsymbol{\delta})\big) \approx f(\boldsymbol{q}) + \langle \nabla\boldsymbol{f}|_{\boldsymbol{q}}\,,\boldsymbol{\delta}\rangle + \frac{1}{2}\boldsymbol{\delta}^T\bigg( \nabla^2f|_{\boldsymbol{q}} - \langle\nabla f|_{\boldsymbol{q}}\,,q\rangle \boldsymbol{I} \bigg)\boldsymbol{\delta} $$
First, observe we have \begin{align} f(g(0+\delta)) = f\left(q+ Dg(0)\delta+\frac{1}{2!}\delta^TD^2g(t^\ast\delta)\delta\right)=: f(q + \Delta) \end{align} for some $0<t^\ast<1$. Here where I have used the first order Taylor expansion with Lagrange remainder \begin{align} g(x+h) = g(x) + Dg(x)h+\frac{1}{2}h^TD^2g(x+t^\ast h) h. \end{align} Next, again using first order Taylor expansion with Lagrange remainder we see \begin{align} f(q+\Delta)= f(q) + Df(q)\Delta + \frac{1}{2!}\Delta^TD^2f(q+\theta^\ast \Delta)\Delta. \end{align} where $0<\theta^\ast<1$. Note that $D^2f$ is a bilinear form which takes in two vectors which in our case is $\Delta$ and $\Delta$ (in this case, it's just a matrix). Hence it follows \begin{align} f(g(0+\delta)) =&\ f(q) + Df(q)\left[Dg(0)\delta+\frac{1}{2}\delta^TD^2g(t^\ast\delta)\delta\right]\\ &+\frac{1}{2!}\left[Dg(0)\delta+\frac{1}{2!}\delta^TD^2g(t^\ast\delta)\delta\right]^T[ D^2f(q+\theta^\ast \Delta)]\left[Dg(0)\delta+\frac{1}{2!}\delta^TD^2g(t^\ast\delta)\delta\right]. \end{align} Since everything is as smooth as a like and $\|\delta\|$ is a small as a want, then it follows that \begin{align} D^2g(t^\ast\delta)\approx D^2g(0), \ \ \text{ and } \ \ D^2f(q+\theta^\ast \Delta)\approx D^2f(q) \end{align} which means \begin{align} f(g(0+\delta)) \approx&\ f(q) + Df(q)\left[Dg(0)\delta+\frac{1}{2!}\delta^TD^2g(0)\delta\right] \\ &+ \frac{1}{2!} \left[Dg(0)\delta+\frac{1}{2!}\delta^TD^2g(0)\delta\right]^T[D^2f(q)]\left[Dg(0)\delta+\frac{1}{2!}\delta^TD^2g(0)\delta\right]\\ \approx&\ f(q)+Df(q)[Dg(0)\delta]+\frac{1}{2!}Df(q)[\delta^T D^2g(0)\delta]+\frac{1}{2!}[Dg(0)\delta]^T [D^2f(q)][Dg(0)\delta]\\ &+ \mathcal{O}(\|\delta\|^3) \end{align} Finally, note that \begin{align} Dg(0)\delta=\frac{d}{dt}g(t\delta)\bigg|_{t=0} =&\ \frac{d}{dt}\left[q \cos(t\|\delta\|)+ \frac{\delta}{\|\delta\|} \sin(t\|\delta\|)\right]\bigg|_{t=0}\\ =&\ -q\|\delta\|\sin(t\|\delta\|)+ \delta\cos(t\|\delta\|)\bigg|_{t=0} =\delta \end{align} and \begin{align} \delta^TD^2g(0)\delta = \frac{d^2}{dt^2}g(t\delta)\bigg|_{t=0} = -\|\delta\|^2q = - (\delta^T\delta)\ q \end{align} Combining everything gives us the desired result \begin{align} f(g(\delta)) \approx&\ f(q) + Df(q)\delta-\frac{1}{2!}(\delta^T\delta)Df(q)q + \frac{1}{2!} \delta^T D^2f(q)\delta\\ &=\ f(q) + Df(q)\delta +\frac{1}{2!}\delta^T\left(D^2f(q) - Df(q)[q] I\right)\delta \end{align}