I wouldn't call it GL, it is the orthogonal group of the lattice we are discussing. References, as i said, include Lattices and Codes by W. Ebeling, Rational Quadratic Forms by Cassels, these two being available and inexpensive.
We double the quadratic form to get an integral Gram matrix:
$$ G \; = \;
\left( \begin{array}{rrr}
2 & 1 & 0 \\\
1 & 2 & 0 \\\
0 & 0 & -2
\end{array}
\right).
$$
Your solutions, the light cone, are column vectors $v$ such that $ v^T G v = 0.$ A root, since we have doubled everything to get an "even" lattice, is a vector $r$ with norm 2, $ r^T G r = 2.$ The general definition of reflection (Cassels calls this a symmetry, page 19) in any vector $w$ is that
$$ x \mapsto \; \; x \; - \; \frac{2 \, x^T G w}{w^T G w} \; w.$$
As a result, when we take $w$ to be a root, the factors of 2 cancel and we are taking lattice points to other lattice points. In your original form, a root $(x,y,z)$ solves $x^2 + x y + y^2 = 1 + z^2.$ Then the reflection in the root is just a linear map, determinant $-1,$ and is therefore given by a square matrix with respect to the original basis. Finally, the reflection is an isometry, part of the orthogonal group of the quadratic form, and if we call the matix $A,$ it solves $$A^T G A = G.$$
I think it wise to include reflection in the root $(1,0,0)^T$ to get some negative values taken care of,
$$ A_0 \; = \;
\left( \begin{array}{rrr}
-1 & -1 & 0 \\\
0 & 1 & 0 \\\
0 & 0 & 1
\end{array}
\right).
$$
We also include reflection in the root $(0,1,0)^T,$
$$ A_{00} \; = \;
\left( \begin{array}{rrr}
1 & 0 & 0 \\\
-1 & -1 & 0 \\\
0 & 0 & 1
\end{array}
\right).
$$
Note that we have already included the automorph that interchanges the first two items, $$ (x,y,z) \mapsto (y,x,z). $$
$$ A_0 A_{00} A_0 \; = \;
\left( \begin{array}{rrr}
0 & 1 & 0 \\\
1 & 0 & 0 \\\
0 & 0 & 1
\end{array}
\right).
$$
Although it may not be called a root, having norm $-2,$ we include reflection in $(0,0,1)^T,$ or
$$ A_{000} \; = \;
\left( \begin{array}{rrr}
1 & 0 & 0 \\\
0 & 1 & 0 \\\
0 & 0 & -1
\end{array}
\right).
$$
Also, as is usual, we explicitly include negation,
$$ (x,y,z) \mapsto (-x,-y,-z), $$ so
$$ A_{0000} \; = \;
\left( \begin{array}{rrr}
-1 & 0 & 0 \\\
0 & -1 & 0 \\\
0 & 0 & -1
\end{array}
\right).
$$
It occurred to me that we would get some reflections with quite small matrix entries by taking other vectors of norm $-2,$ for example $(1,1,2)^T.$
$$ A_{-1} \; = \;
\left( \begin{array}{rrr}
4 & 3 & -4 \\\
3 & 4 & -4 \\\
6 & 6 & -7
\end{array}
\right).
$$ This still solves, as with all the $A'$s, both $A^T G A = G$ and $A^2 = I.$
Taking the root $(3,4,6)^T$ gives the reflection
$$ A_1 \; = \;
\left( \begin{array}{rrr}
-29 & -33 & 36 \\\
-40 & -43 & 48 \\\
-60 & -66 & 73
\end{array}
\right).
$$
If we take the initial triple in the light cone to be $c = (3,5,7)^T,$ we get $A_1 c = (0,1,1)^T$ which is indeed another solution, and shows that care must be used in constructing the "tree."
Taking the root $(5,15,18)^T$ gives the reflection
$$ A_2 \; = \;
\left( \begin{array}{rrr}
-124 & -175 & 180 \\\
-375 & -524 & 540 \\\
-450 & -630 & 649
\end{array}
\right).
$$
we get $A_2 c = (13,35,43)^T$
Taking the root $(8,19,24)^T$ gives the reflection
$$ A_3 \; = \;
\left( \begin{array}{rrr}
-279 & -368 & 384 \\\
-665 & -873 & 912 \\\
-840 & -1104 & 1153
\end{array}
\right).
$$
we get $A_3 c = (11,24,31)^T$
Taking the root $(7,32,36)^T$ gives the reflection
$$ A_4 \; = \;
\left( \begin{array}{rrr}
-321 & -497 & 504 \\\
-1472 & -2271 & 2304 \\\
-1656 & -2556 & 2593
\end{array}
\right).
$$
we get $A_4 c = (80,357,403)^T$
Taking the root $(12,47,54)^T$ gives the reflection
$$ A_5 \; = \;
\left( \begin{array}{rrr}
-851 & -1272 & 1296 \\\
-3337 & -4981 & 5076 \\\
-3834 & -5724 & 5833
\end{array}
\right).
$$
we get $A_5 c = (159,616,709)^T$
One typically includes $\pm 1$ anyway. Ian Agol would know how many reflections are enough, but I suspect this will do. Nothing really wrong with finding too many.
Note that the squarefree parts in $1+z^2$ in the five nontrivial roots I chose are $37, 13, 577, 1297, 2917.$
