Iwahori versus Bruhat decompositions
Your Iwahori-decomposition computation is a bit too free with quotient computations. It would work fine for vector spaces, which is, in some sense, why you get the correct leading term in your count; but there are additional subtleties on the group level that it does not take into account—among other things, that some of the entries denoted by $\mathcal O$ can be $0$, but some cannot!
The same computation would suggest that the quotient of $\operatorname{GL}_2(k)$ by its Borel subgroup of upper-triangular matrices should have as coset representatives the lower uni-triangular matrices $u_-(c) = \begin{pmatrix} 1 & 0 \\ c & 1 \end{pmatrix}$, whereas actually there is another coset, the one containing $w_0 = \begin{pmatrix} 0 & 1 \\ -1 & 0 \end{pmatrix}$. Actually a better way to think about it is that most coset representatives are of the form $u_+(c)w_0$, where $u_+(c)$ is the transpose of $u_-(c)$, and only the trivial coset is not of this form; it just happens that $u_+(c)w_0$ and $u_-(c^{-1})$ lie in the same coset when $c \ne 0$.
Your count via the Bruhat decomposition shows what’s going on here: you get $(1 + N(p))(1 + N(p)^2) = 1 + N(p) + N(p)^2 + N(p)^3$, which suggests, correctly, that we have stratified $\mathbb P^3k$ by affine spaces of the obvious dimensions, corresponding to the Weyl-group elements $1$, $x$, $x y$, and $x y x$ (which are minimal-length representatives in the quotient by an appropriate parabolic subgroup, which is doubtless what led you to choose them!).