Feynman rules from interaction Lagrangian with electromagnetic tensor (vertex)

I am currently studying for my QFT exam and in particular learning the methods of reading the Feynman rules directly off the Lagrangian.

However, I'm still a bit uncertain how to deal with derivative terms in the interaction Lagrangian in complicated cases.

In this example, the interaction Lagrangian is coupling a real (pseudo)scalar field to the photon field, supposing to represent the effective coupling induced by the Loop of complex scalar particles.

The effective interaction Lagrangian is given by $$\mathcal{L}_{int,eff}=-\frac{g_{eff}}{4}\phi F^{\mu\nu}F_{\mu\nu}$$ in addition to the QED and scalar field lagrangian. The resulting vertex is: Now my general ansatz has been to decompose the part containing derivatives so basically $$F^{\mu\nu}F_{\mu\nu}=(\partial^\mu\,A^\nu-\partial^\nu\,A^\mu)(\partial_\mu\,A_\nu-\partial_\nu\,A_\mu)$$ and going for the usual Fourier decomposition to the replacement $$\partial_\mu \rightarrow ip_{\mu,i}$$ and multiplying with $2!$ for the interchangability of the outgoing photons. This gives naively the following terms:

$$p_1^\mu p_{\mu,1}+p_1^\mu p_{\nu,2}+p_2^\nu p_{\mu,1}+p_2^\nu p_{\nu,2} = 2p_1\cdot\,p_2 +p_1^\mu p_{\nu,2}+p_2^\nu p_{\mu,1}$$

Which only leads almost to the desired result. A main problem is that this somewhat causes an index mess that hardly looks promising.

What am I missing? Maybe I can not simply cast the derivative on the photon field in co- and contravariant Expression? Extracting a factor $g^{\mu\nu}$ does not help if I didn't make an algebraic error.

By recasting the Lagrangian vertex in Fourier space one has (ignoring $-g_{eff}/4$) $$-\phi(p_3) \left(p_1^\mu A^\nu(p_1) - p_1^\nu A^\mu(p_1)\right)\left(p_{2\mu} A_\nu(p_2) - p_{2\nu} A_\mu(p_2)\right)\,,$$ which by a simple relabeling can be recast as $$-2\,\phi(p_3) \,p_1^\mu A^\nu(p_1)\left(p_{2\mu} A_\nu(p_2) - p_{2\nu} A_\mu(p_2)\right) =$$ $$= -2\,\phi(p_3) \,\left((p_1\cdot p_2)\, A^\nu(p_1)A_\nu(p_2) - p_{1\mu} p_{2\nu} \, A^\nu(p_1)A^\mu(p_2)\right)$$ $$= -2\,\phi(p_3) \,\left((p_1\cdot p_2)\, g_{\rho\sigma} A^\sigma(p_1)A^\rho(p_2) - p_{1\rho} p_{2\sigma} \, A^\sigma(p_1)A^\rho(p_2)\right)\,.$$ At this point we just take the functional derivative $$\frac{(-\mathrm{i})^3\,\delta^3}{\delta\phi(p_3)\,\delta A^\mu(p_2)\,\delta A^\nu(p_1)}\,.$$ Each term, being quadratic in $A$, gets a factor two, resulting in $$= -4\mathrm{i}\,\left((p_1\cdot p_2)\, g_{\mu\nu} - p_{1\mu} p_{2\nu} \right)\,.$$ Which is the desired result.