Determinant of a specific $4 \times 4$ symmetric matrix

Here's a method for calculating the determinant, explaining at least why it ends up as a product. I don't know if there's any significance to your determinant being a square.

Define $$H= \frac{1}{\sqrt{2}} \begin{pmatrix} 1 & 1 & 0 & 0 \\ 1 & -1 & 0 & 0 \\ 0 & 0 & 1 & 1 \\ 0 & 0 & 1 & -1 \\ \end{pmatrix}. $$ (The tensor product of a one-dimensional Hadamard matrix with a two-by-two identity matrix.)

Then $\det H=1$ and for any $a,b,c,d,e,f,g,h$, $$ H \begin{pmatrix} a & b & c & d \\ b & a & d & c \\ e & f & g & h \\ f & e & h & g \end{pmatrix} H\\ =\begin{pmatrix} a+b & 0 & c+d & 0 \\ 0 & a-b & 0 & c-d \\ e+f & 0 & g+h & 0 \\ 0 & e-f & 0 & g-h \end{pmatrix}$$ which is (similar to) $$\begin{pmatrix} a+b & c+d \\ e+f & g+h \\ \end{pmatrix} \oplus \begin{pmatrix} a-b & c-d \\ e-f & g-h \\ \end{pmatrix}. $$

Plugging in a rotated version of your matrix gives $$\begin{vmatrix} 2 & x & b & a \\ x & 2 & a & b \\ b & a & 2 & ab-x \\ a & b & ab-x & 2 \end{vmatrix} \\= \begin{vmatrix} 2+x & a+b \\ a+b & 2+ab-x \\ \end{vmatrix} \cdot \begin{vmatrix} 2-x & b-a \\ b-a & 2-ab+x \\ \end{vmatrix} \\ =(4-x^2+abx-a^2-b^2)(4-x^2+abx-a^2-b^2). $$

The following answer is inspired by Colin's.

As noted by Colin MacQuillan, the matrix under consideration turns out to be similar to
the block matrix $\begin{bmatrix} A & C \\ C & B \end{bmatrix}$ where $A=2I_2+x K$, $C=a I_2+ b K$, $B=2 I_2+(ab-x) K$, with $K:=\begin{bmatrix} 0 & 1 \\ 1 & 0 \end{bmatrix}$. Since $I_2$ and $K$ commute, it is then standard that $$\begin{vmatrix} A & C \\ C & B \end{vmatrix}=\det(AB-C^2).$$ Here, we see that $AB=(4+x(ab-x))I_2+2ab K$ and $C^2=(a^2+b^2) I_2+2ab K$. Hence, $AB-C^2=(4-a^2-b^2+x(ab-x))I_2$, which yields the claimed result.

Generalizing the diagonal entry of the initial matrix to a $c$ gives, with a similar method, that the determinant equals that of $(c^2-a^2-b^2+x(ab-x))I_2+(c-2)ab K$, yielding Andreas's formula.

In short: determinant of every symmetric matrix is a square!

Consider arbitrary three matrices $a, b, c \in Sl_2.$ One can wonder, what relations do the traces of their products satisfy. The answer is given by the famous Jimbo-Fricke cubic: \begin{equation} \begin{split} &tr(ab)\ tr(bc) \ tr(ac)+ tr(ab)^2+tr(bc)^2+tr(ac)^2\\ &+tr(a)^2+tr(b)^2+tr(c)^2+tr(abc)^2\\ &-(tr(a)tr(b)+tr(c)tr(abc))tr(ab)\\ &-(tr(b)tr(c)+tr(a)tr(abc))tr(bc)\\ &-(tr(a)tr(c)+tr(b)tr(abc))tr(ac)\\ &+tr(a)tr(b)tr(c)tr(abc)-4=0.\\ \end{split} \end{equation}

Every determinant of a symmetric matrix can be written in the following form for some matrices $a, b, c$: $$ G=\begin{vmatrix} 2 & -tr(a) & -tr(b) & -tr(bc)\\ -tr(a) & 2 & -tr(ab)& -tr(abc)\\ -tr(b) & -tr(ab) & 2 & -tr(c)\\ -tr(bc) & -tr(abc) & -tr(c)& 2\\ \end{vmatrix}. $$ The relation above is equivalent to the following: $$(2tr(ac)+tr(ab)tr(bc)-tr(a)tr(c)-tr(b)tr(abc))^2=G.$$

Usually a symmetric determinant is not a square, because $tr(ac)$ is not a polynomial in the entries of $G.$ The case of the matrix in the question corresponds to $c=a^{-1},$ because $tr(c)=tr(c^{-1}),$ $tr(aba^{-1})=tr(a)$ and $tr(ab)+tr(ba^{-1})=tr(a)+tr(b).$ The square root of $G$ is algebraic, because $tr(ac)=2.$

I have seen this presentation of the Jimbo-Fricke cubic only in one place: https://arxiv.org/pdf/1308.4092.pdf, formula (3.9), and I will be really grateful for any references.