There are at least 3 ways of seeing this - an elementary way, via the determinants of $\ell$-adic sheaves, and via Katz's calculus of finite field hypergeometric functions.
Elementary:
Subtraction inside a multipliciative character is too difficult. Let's add a new variable
$$=\chi(n) \sum_{r,t \in \mathbb F_p} \chi(t) \delta_{r^2-g,t} \zeta_p^{nr}$$
and detect using additive characters
$$=\chi(n) \frac{1}{p} \sum_{r,t,s \in \mathbb F_p} \chi(t) \zeta_p^{nr+s (r^2-g-t)}$$
and separate variables
$$ = \chi(n) \frac{1}{p} \sum_{s \in \mathbb F_p}\zeta_p^{-sg} \left( \sum_{t \in \mathbb F_p} \chi(t) \zeta_p^{-st} \right) \left(\sum_{r \in \mathbb F_p} \zeta_p^{n r + sr^2} \right)$$
Both sums vanish for $s=0$ and for other $s$ are standard Gauss sums
$$ = \chi(n) \frac{1}{p} \sum_{s \in \mathbb F_p} \zeta_p^{-sg} \left( \chi^{-1}(s) G(\chi) \right) \left(\chi_2(s)
\zeta_p^{- \frac{n^2}{4s} } G(\chi_2) \right) $$
$$ = \frac{ G(\chi) G(\chi_2) }{p} \chi(n) \sum_{s \in \mathbb F_p} \chi_2(s) \chi^{-1}(s) \zeta_p^{ -sg - \frac{n^2}{4s}}$$
We can now calculate the argument with the substitution $s \mapsto - \frac{n^2}{4gs}$, which sends $\zeta_p^{ -sg - \frac{n^2}{4s}}$ to its complex conjugate, $\chi^{-1}(s)$ to its complex conjugate times $\chi^{-1}(\frac{n^2}{4g})$, and $\chi_2(s)$ to its complex conjugate times $\chi_2(g)$, so it sends the whole sum to its complex conjugate times $\chi^{-1}(n^2) \chi(4g)/\chi_2(g)$, and hence the argument of the inner sum is half the argument of $\chi^{-1}(n^2) \chi(4g)/\chi_2(g)$, which when multiplied by $\chi(n)$ is half the argument of $\chi(4g)/\chi_2(g)$ and hence is independent of $n$.
Determinant:
This family of exponential sum is associated to a rank 2 sheaf (a slight variant of a hypergeometric sheaf). This means the exponential sum corresponds to the trace of Frobenius, a unitary matrix, on a rank two representation. The argument of the trace of a unitary matrix is half the argument of the determinant.
So your statement is equivalent to the statement that the determinant of the family of exponential sums is constant. The determinant of an individual exponential sum can be calculated using local epsilon factors. For a one-parameter exponential sum like this, if we only want to calculate the ratios between the different determinants, it can be easier. Because your family of sums can be expressed as a Fourier transform, it is easy to apply Laumon's theory of the $\ell$-adic Fourier transform to calculate the local monodromy of the sheaf, take its determinant, and verify that it is trivial. Because the determinant sheaf has trivial local monodromy it has trivial global monodromy and must be constant.
So there are many different determinant formulas like this that can be proven, especially for simple Fourier transforms. Only for rank $2$ sheaves will they give argument formulas. For example, the graph of $\chi(n) \sum_r \chi(r^d-g) \zeta_p^{nr}$ in the complex plane should form a $d$-pointed star (where a $2$-pointed star is a line segment, i.e. a set of bounded numbers with fixed argument modulo $\pi$.)
Hypergeometric calculus - TBC