My problem is to take the $d$-dimensional Ising Hamiltonian, $$H = -\sum_{i,j}\sigma_i J_{i,j} \sigma_j - \sum_{i} \tilde{h}_i \sigma_i$$ where $J_{ij}$ is a matrix describing the couplings between sites $i$ and $j$. Applying a Hubbard-Stratonovich transformation, rewrite the partition function as $$Z = N_0 \int d^N \psi \exp\left\{-\left[\frac{1}{4}\sum_{i,j} \psi_i K_{ij} \psi_j - \sum_{i} \ln[\cosh(h_i+\psi_i)]\right]\right\}$$ where $N_0$ is an overall normalization constant, $K_{ij} = (\beta J_{ij})^{-1}$, and $h_i = \beta\tilde{h}_i$. This much is relatively straightforward. We write the field as $\psi_i = \phi_i - h_i$, and we can show that $\left<\phi_i\right> \propto J_{ij} \left<\sigma_j\right>$, i.e. it can be interpreted as a "mean field" at site $i$ due to the interaction with all other sites.
Next we assume that the variation in the field is small, $\left|\phi_i\right|<<1$, we set $h_i = 0$, and expand $\ln \cosh(x) \approx \frac{1}{2}x^2 - \frac{1}{12}x^4$ to get $$Z \approx N_0\int d^N\psi \exp\left\{-\left[\frac{1}{4}\sum_{i,j}\phi_i K_{ij} \phi_j - \sum_i \left[\frac{\phi_i^2}{2} - \frac{\phi_i^4}{12}\right]\right]\right\}$$
Now we take the continuum limit, in units where the lattice spacing is unity, labeling each site by its position $\mathbf{r}$, which gives $$Z\rightarrow \mathcal{N} \int \mathcal{D}\phi\, \exp\left\{-\frac{1}{2}\left[\frac{1}{2}\int d\mathbf{r}\,d\mathbf{r}'\,\phi(\mathbf{r}) K(\mathbf{r}-\mathbf{r}') \phi(\mathbf{r}') - \int d\mathbf{r}\,\left[\phi(\mathbf{r})^2 - \frac{\phi(\mathbf{r})^4}{6}\right]\right]\right\}$$
This is where I am not sure how to proceed. I am told to expand $\phi(\mathbf{r}')$ as a small variation from the value at $\mathbf{r}$, i.e. $$\phi(\mathbf{r}') \approx \phi(\mathbf{r}) + (x_\mu'-x_\mu)\partial_\mu \phi(\mathbf{r}) + \frac{1}{2}(x_\mu' - x_\mu)(x_{\nu}'-x_\nu)\partial_\mu \partial_\nu \phi(\mathbf{r}) + \cdots$$ and introduce the Fourier transform $\tilde{K}(\mathbf{q}) = \int d\mathbf{r} K(\mathbf{r}) e^{-i\mathbf{q}\cdot\mathbf{r}}$ and write the continuum action as $$S = \int d^d\mathbf{r} \left[c_1 \left(\partial \phi\right)^2 + c_2 \phi^2 + c_4 \phi^4\right]$$ and find the coefficients in terms of $\tilde{K}(0)$ and $\tilde{K}''(0)$.
I believe that I can argue that $K$ is only a function of $\left|\mathbf{r}-\mathbf{r}'\right|$, in which case $K(\mathbf{r}-\mathbf{r}')(x_\mu'-x_\mu)$ is odd about the point $\mathbf{r}$, and so integrating over $d\mathbf{r}'$ (treating $\mathbf{r}$ as constant) will kill any term except those which depend on the square of the difference, leaving me with $$\int d\mathbf{r}\,d\mathbf{r}'\,\phi(\mathbf{r}) K(\mathbf{r}-\mathbf{r}') \phi(\mathbf{r}') = \int d\mathbf{r}\,d\mathbf{r}' K(\mathbf{r}-\mathbf{r}')\left(\phi(\mathbf{r})^2 + \frac{1}{2}(x_\mu'-x_\mu)^2\phi(\mathbf{r})\partial_\mu^2 \phi(\mathbf{r})\right) $$
The first term I can deal with, but it's the second term that I don't know how to deal with.