EDITTTT: two research level articles on this, one by Kontorovich and one by Fuchs, can be downloaded for free at AMS BULLETIN APRIL 2013. There is also a short survey article by Peter Sarnak, in the April 2011 MAA Monthly. I have a pdf of that, if anyone is interested.
Alright, I see part of what is going on. The standard recipe, stereographic projection around a given integral solution, is guaranteed to give all rational solutions; this method parametrizes solutions with four components by three parameters, or generally $k$ components by $k-1$ parameters. And, some choices of central point for projection give more efficient recipes than others; the one individ chose was very good.
However, it appears that such a recipe is not going to give all integral solutions, nor will a finite number of such recipes. So, here is my version, with rational solutions. That is, I took coprime $(x,y,z)$ and calculated
$$ a = x^2 + y^2 + 2 z^2 - y z - z x - 2 x y, $$
$$ b = x^2 - y z + z x,$$
$$ c = y^2 + y z - z x,$$
$$ d = 2 z^2 + y z + z x.$$
Then I took $$ g = \gcd(a,b,c,d), $$ and divided all four of $(a,b,c,d)$ by that. As a result, I found all the "root" solutions given at TABLE. I am encouraged that $g$ was always the sum of two squares. Hmmm; actually, it is easy to show that $\gcd(x,y,z) = 1$ implies that $ g = \gcd(a,b,c,d) $ is not divisible by $4$ or by any prime $q \equiv 3 \pmod 4,$ so we can always write $g = u^2 + v^2$ with integers and $\gcd(u,v)=1.$
Oh, root solutions have
$$ a \leq 0 \leq b \leq c \leq d, $$
$$ a+b+c+d > 0, $$
$$ a + b + c \geq d. $$ These were defined in Graham, Lagarias, Mallows, Wilks, Yan, Apollonian circle packings: Number Theory, Journal of Number Theory, volume 100 (2003) pages 1-45. It was shown that every integral solution can be connected to such a root solution by Vieta jumps, thus dividing the solutions into a forest of countably many rooted trees. An unusual feature is that these Vieta jumps are (they must be) elements in a group of $4$ by $4$ invertible integral matrices, called the Apollonian Group, which is just the (orthogonal or rotation or automorphism) group for the quadratic form/indefinite lattice given by $$ a^2 + b^2 + c^2 + d^2 -2ab-2ac-2ad-2bc-2bd-2cd, $$ with Gram matrix
$$ H \; = \;
\left( \begin{array}{rrrr}
1 & -1 & -1 & -1 \\
-1 & 1 & -1 & -1 \\
-1 & -1 & 1 & -1 \\
-1 & -1 & -1 & 1
\end{array}
\right).
$$
As a matrix, eigenvalues are $-2,2,2,2$ with orthogonal (but not orthonormal) eigenvectors as columns of
$$ W \; = \;
\left( \begin{array}{rrrr}
1 & -1 & -1 & -1 \\
1 & 1 & -1 & -1 \\
1 & 0 & 2 & -1 \\
1 & 0 & 0 & 3
\end{array}
\right).
$$
a b c d g x y z
-1 2 2 3 gcd was 1 1 0 1
-2 3 6 7 gcd was 5 -3 -6 2
-3 4 12 13 gcd was 10 4 12 -3
-3 5 8 8 gcd was 1 2 3 -1
-4 5 20 21 gcd was 5 -2 1 -7
-4 8 9 9 gcd was 29 -1 12 9
-5 6 30 31 gcd was 26 -6 -30 5
-5 7 18 18 gcd was 10 -4 2 -9
-6 7 42 43 gcd was 37 -7 -42 6
-6 10 15 19 gcd was 13 -10 -15 6
-6 11 14 15 gcd was 2 0 -4 -3
-7 8 56 57 gcd was 13 -3 2 -19
-7 9 32 32 gcd was 17 5 -3 16
-7 12 17 20 gcd was 5 -8 -9 -4
-8 9 72 73 gcd was 65 9 72 -8
-8 12 25 25 gcd was 25 -14 -27 8
-8 13 21 24 gcd was 5 7 11 -4
-9 10 90 91 gcd was 82 -10 -90 9
-9 11 50 50 gcd was 5 -4 -17 3
-9 14 26 27 gcd was 26 -16 -28 9
-9 18 19 22 gcd was 17 -1 13 11
-10 11 110 111 gcd was 25 -4 3 -37
-10 14 35 39 gcd was 10 -6 2 -13
-10 18 23 27 gcd was 1 -4 -5 2
-11 12 132 133 gcd was 122 -12 -132 11
-11 13 72 72 gcd was 37 -7 5 -36
-11 16 36 37 gcd was 34 -18 -38 11
-11 21 24 28 gcd was 5 -7 1 -7
-12 13 156 157 gcd was 145 -13 -156 12
-12 16 49 49 gcd was 5 -6 -17 4
-12 17 41 44 gcd was 1 -3 -7 2
-12 21 28 37 gcd was 25 21 28 -12
-12 21 29 32 gcd was 2 -6 -8 3
-12 25 25 28 gcd was 1 -5 -5 -2
-13 14 182 183 gcd was 41 -5 4 -61
-13 15 98 98 gcd was 50 -8 6 -49
-13 18 47 50 gcd was 2 -2 -8 -5
-13 23 30 38 gcd was 26 -16 2 -19
-14 15 210 211 gcd was 197 -15 -210 14
-14 18 63 67 gcd was 53 -18 -63 14
-14 19 54 55 gcd was 1 -2 1 -5
-14 22 39 43 gcd was 37 -24 -41 14
-14 27 31 34 gcd was 10 -16 -18 7
-15 16 240 241 gcd was 226 -16 -240 15
-15 17 128 128 gcd was 13 -6 -43 5
-15 24 40 49 gcd was 34 -24 -40 15
-15 24 41 44 gcd was 10 2 16 11
-15 28 33 40 gcd was 2 0 -6 -5
-15 32 32 33 gcd was 2 8 8 -3
a b c d g x y z