Next, reverse this bijection. Take any orthonormal matrix M with top row r and multiply on the right by S inverse, (also known as S transpose), and the result is orthonormal, with 1 in the upper left. Therefore the number of matrices in O with top row r is equal to the number of orthonormal matrices of dimension n-1. Of course this assumes a standard matrix S exists.
Look at the case n = 3 over a finite field F of size w. We would like to produce a standard matrix S for each top row r. But how many rows are there? Remember that the length of r is 1. When the upper left entry is ±1, the other two have squares that sum to 0. This contributes 4w-2, or 2, as w is 1 or 3 mod 4. For anything else in the upper left, the remaining two squares sum to something nonzero. We showed earlier that sums of squares are distributed evenly across the nonzero entries of F. This yields w-2 times (w-1 or w+1), as w is 1 or 3 mod 4. Add these up and get w2+w or w2-w, as w is 1 or 3 mod 4. Finally, multiply this expression by the size of O for n = 2, which was determined earlier. The result is w(w-1)(w+1), regardless of w mod 4. But we still have to demonstrate a standard orthonormal matrix S for each row r.
Start with r and extend this to an n dimensional basis. Scale one of the rows, so that the determinant is 1. Then apply the gram schmidt process to create an orthogonal basis. The determinant is still 1. Dot each row with itself, and multiply the results together, and get the square of the determinant, which is 1. An even number of dot products are nonsquares.
Scale each row so that its dot product is 1 or q, where q is your favorite nonsquare in F. If all dot products are 1 then S is orthonormal, and we're done. Otherwise consider a pair of rows with dot products equal to q. Call these two vectors u and v. Replace u with au+bv, and v with -bu+av, where a and b are arbitrary scalars in F. Dot either of these linear combinations with any other row and the result is still 0. Verify that (au+bv).(-bu+av) = 0. Now evaluate (au+bv).(au+bv) and get (a2+b2)q. Similarly, (-bu+av).(-bu+av) = (a2+b2)q. Since sums of squares are evenly distributed across the nonzero entries of F, there is an a and b such that a2+b2 = 1/q. The two rows now have length 1. Do this for every pair of nonsquare rows, and S becomes orthonormal. That completes the proof.
With the size of O3 established, you might think it is isomorphic to G2. After all, both groups have the same order, namely w(w-1)(w+1). But G2 has a center, namely ±1, while O3 has only 1 as its center. So the groups are not isomorphic.
Just to review:
z1 = 1, y1 = 1, o1 = 1.
z2 = 2w-1 | 1, y2 = w-1 | w+1, o2 = w-1 | w+1.
z3 = w2, y3 = w2+w | w2-w, q3 = w2-w | w2+w, o3 = w3-w
When n is even, build zn in two ways. Let the first n-2 entries have squares that sum to 0, and then bring in the last two entries. This gives zn-2z2. Then the first n-2 entries could produce something nonzero, and the last two must take up the slack. Similar reasoning holds for yn.
zn = zn-2z2 + (w-1)yn-2y2
yn = zn-2y2 + yn-2z2 + (w-2)yn-2y2
Even though z2 and y2 depend on the value of w mod 4, z4 and y4 do not.
z4 = w3+w2-w, y4 = w3-w, o4 = w2(w-1)2(w+1)2.
Verify that our formula for y4 comes out the same for a nonsquare length. Like y2, y4 is the same across the nonzero lengths. In fact, as you step through n = 6, 8, 10, and so on, qn always equals yn. The distribution of vectors is even across the nonzero lengths.
zn = zn-1 + (w-1)yn-1
yn = (w-2)yn-1 + 2zn-1
qn = wyn-1
z5 = w4, y5 = w4+w2, q5 = w4-w2, o5 = w4(w-1)2(w+1)2(w2+1).
For dimension 7, separate the first entry and use the formulas that were applied to dimension 5.
Continue this process as far as you like. In 6, as with 2, the formula for o depends on w mod 4.
z6 = w5+w3-w2 | w5-w3+w2, y6 = w5-w2 | w5+w2, o6 = z6o5
z7 = w6, y7 = w6+w3 | w6-w3, q7 = w6-w3 | w6+w3, o7 = z7o6
z8 = w7+w4-w3, y8 = w7-w3, o8 = z8o7