What you call finite field contour sums are more usually called exponential sums. Like complex analytic contour integrals, they do not always (or even usually) have a nice exact formula. Instead, sometimes they have a completely satisfactory exact formula, sometimes they have an expression in terms of special functions, and sometimes they are are completely new functions with no simplified expression. Your sum is of the second type - it can be expressed in terms of hypergeometric functions.
This is the trace of Frobenius on the cohomology of a certain sheaf on the affine curve $a^2 - b^2 \delta =\delta/4 ( 4-r^2)$, where I am taking $z=a+ \sqrt{\delta} b$. As $r$ varies this curve varies in a family of curves, and you can take the cohomology of the family of sheaves, giving a sheaf on $\mathbb A^1$. How the sum varies with $r$ depends on the monodromy of this sheaf.
The sheaf is:
$$\mathcal L_\chi \left(a + b \sqrt{\delta} + (r-1) \sqrt{\delta} \right) \otimes \mathcal L_{\overline{\chi}} \left(a +b \sqrt{\delta} + (r+1) \sqrt{\delta} \right)$$
To study this we can look at it over a larger finite field. In particular, we can look at this sum over $\mathbb F_p ( \sqrt{\delta})$. Let's call the square root of $\delta$ that we are adjoining $\sigma$. Over this field, $a^2 - b^2 \delta$ factors as $(a+b \sigma)$ $(a-b\sigma)$. The sheaf $\mathcal L_\chi (a + b \sqrt{\delta} + (r-1) \sqrt{\delta})$ factors under the norm map as
$$\mathcal L_\chi (a + b \sqrt{\delta} + (r-1) \sqrt{\delta}) \otimes \mathcal L_\chi (a^* + b^* \sqrt{\delta} + (r-1) \sqrt{\delta})$$
where $*$ is conjugation by the field automorphism. This may also be written:
$$\mathcal L_\chi (a + b \sqrt{\delta} + (r-1) \sqrt{\delta}) \otimes \mathcal L_{\chi^*} (a- b\sqrt{\delta} - (r-1) \sqrt{\delta})$$
Over the field $F_{p^2}$, your sum:
$$ \sum_{a^2 - b^2 \delta =\frac{\delta}{4} ( 4-r^2)} \chi \left(a + b \sqrt{\delta} + (r-1) \sqrt{\delta} \right) \overline{\chi} \left(a +b \sqrt{\delta} + (r+1) \sqrt{\delta} \right)$$
becomes:
$$ \sum_{(a+b\sigma)(a-b\sigma) =\frac{\delta}{4} ( 4-r^2)} \chi \left(a + b \sigma+ (r-1) \sigma\right) \chi^* \left(a - b \sigma - (r-1) \sigma\right) \overline{\chi} \left(a +b \sigma + (r+1) \sigma \right)\overline{\chi}^* \left(a -b \sigma - (r+1) \sigma \right)$$
Now we can simplify this by setting $\alpha = a + \sigma b$ so $a-\sigma b = \frac{\frac{\delta}{4} (4-r^2)}{\alpha}$
$$ \sum_{\alpha} \chi \left(\alpha + (r-1) \sigma\right) \chi^* \left(\frac{\frac{\delta}{4} (4-r^2)}{\alpha} - (r-1) \sigma\right) \overline{\chi} \left(\alpha + (r+1) \sigma \right)\overline{\chi}^* \left(\frac{\frac{\delta}{4} (4-r^2)}{\alpha} - (r+1) \sigma \right)$$
or:
$$ \sum_{\alpha} \chi \left(\frac{ \alpha + (r-1) \sigma} {\alpha + (r+1) \sigma} \right) \chi^* \left(\frac{\frac{\frac{\delta}{4} (4-r^2)}{\alpha} - (r-1) \sigma}{\frac{\frac{\delta}{4} (4-r^2)}{\alpha} - (r+1) \sigma} \right) $$
This looks ugly but it's just $\chi$ of a rational function with one zero and one pole times $\chi^*$ of a rational function with one zero and one pole. In sheaf-theoretic terms, we have a sheaf of rank $1$ on $\mathbb P^1$ with four singularities with tame monodromy at each, for an Euler characteristic of $-2$, so their are two eigenvalues of Frobenius of size $\sqrt{p}$. (Here I am including $\alpha=0$ and $\alpha=\infty$ to make the curve $\mathbb P^1$ to simplify things, which does not affect the original sum because those points are not defined over $\mathbb F_p$.)
So the exponential sum is a complex number of size at most $2\sqrt{p}$.
How this sum depends on $r$ is determined by the monodromy of this cohomology group as a sheaf over the $r$ line. We can compute this by noting that the sum can be expressed as:
$$\sum_{\alpha} \chi \left( \frac{ \alpha-b}{\alpha-a} \right) \chi^* \left( \frac{ \alpha-d}{\alpha-c} \right)$$
for some $a,b,c,d$, up to $\chi$ of something times $\chi^*$ of something else. We may change variables by an automorphism of $\mathbb P^1$ to send $a$ to $\infty$, $b$ to $1$, $c$ to $0$, and $d$ to a parameter $t$, so we obtain
$$ \sum_{\alpha} \chi \left( \alpha-1 \right) \chi^* \left( 1- \frac{t}{\alpha} \right)$$
We have written this sum as a multiplicative convolution of the function $\chi(\alpha-1)$ and $\chi^* (1-\beta)$. Thus it is expressed as a multiplicative convolution sheaf. Because this formula is analogous to the integral expression of a hypergeometric sum as a multiplicative convolution, this is called a hypergeometric sheaf. In Katz's notation for hypergeometric sheaves, it is $\mathcal H( \chi, \chi^* ; 1, 1 )$. All hypergeometric sheaves are irreducible, so the monodromy group acts irreducibly on the $2$-dimensional sheaf. Because there is a repeated character at $\infty$, the local monodromy at $\infty$ contains a nontrivial unipotent element, so the monodromy group must contain $SL_2$, because it acts irreducibly and contains a unipotent.
So the monodromy group of the hypergeometric sum contains $SL_2$. As long as the cross-ratio map that expresses the cross-ratio of $a$, $b$, $c$, $d$ in terms of $r$ is nonconstant, your sheaf will have monodromy containing $SL_2$ as well. In Mathematica I got that the cross ratio of the roots and poles of the polynomials is
$$\frac{9r^4-16 r^2 + 64}{9 r^4}$$
That is nonconstant, so the monodromy is indeed large, and the sheaf is nonsingular, meaning the exponential form is expressible by the sum of two numbers of norm $\sqrt{p}$ instead of a different way, when this function is not $0$, $1$, or $\infty$.
Because the monodromy group contains $SL_2$, the sum cannot be expressed in terms of simpler sums that have monodromy groups different from $SL_2$, and it will behave in a seemingly random way, with the absolute value of the exponential sum distributed according to the Sato-Tate measure.
There is a lot of theory of this type of sum that I have only barely scratched the surface of, and have definitely not really given a satisfactory explanation here. If you are interested in this sort of thing you should look into it.
