Finding local orthonormal frame on a Pseudo-Riemannian Manifold

Yes, this is basically fine. One small correction: You wrote

Because of this $|X_1| = \sqrt{|g(X_1,X_1)|}$ is a non-vanishing smooth vector field on $U_1$, ...

But $|X_1|$ is not a vector field. It should say

Because of this $|X_1| = \sqrt{|g(X_1,X_1)|}$ is a non-vanishing smooth function on $U_1$, ...

Let $(M,g)$ be a pseudo-Riemannian $n$-manifold and $p\in M$ be an arbitrary point. Consider a local frame $(\partial_1,\dots,\partial_n)$ in some coordinate neighbourhood $U\ni p$ and an arbitrary orthonormal basis $(V_1,\dots,V_n)$ in $T_pM$. We have $(V_1,\dots,V_n)=((\partial_1)_p,\dots,(\partial_n)_p) A$ for the transition matrix $A\in\mathbb{R}^{n\times n}$, so we can set a new local frame $(X_1,\dots,X_n)$ with $((X_1)_q,\dots,(X_n)_q)=((\partial_1)_q,\dots,(\partial_n)_q) A$ for all $q\in U$. Since $g_p((X_i)_p,(X_j)_p)=g_p(V_i,V_j)$ and $g(X_i,X_j)\in\mathfrak{F}(U)$, then by continuity there exists an open neighbourhood $W\subseteq U$ of $p$ such that \begin{equation*} \lvert{g(X_i,X_i)}\rvert >a=\frac{3}{4}+\frac{1}{3n}<1 \mbox{ and } \lvert{g(X_i,X_j)}\rvert <b=\frac{1}{6n} \mbox{ for all } 1\leq i\neq j\leq n. \end{equation*} Let us apply the Gram-Schmidt process to create an orthogonal local frame $(Y_1,\dots,Y_n)$ on $W$ such that \begin{equation*} \lvert{g(Y_j,Y_j)}\rvert >d=\frac{3}{4} \mbox{ and } \lvert{g(Y_j,X_k)}\rvert <c=\frac{1}{2n} \mbox{ for all } 1\leq j < k\leq n. \end{equation*} We start with $Y_1=X_1$. Suppose that we inductively set $(Y_1,\dots,Y_{j-1})$ for some $2\leq j\leq n$ and define \begin{equation*} Y_j=X_j-\sum_{i=1}^{j-1}\frac{g(X_j,Y_i)}{g(Y_i,Y_i)}Y_i. \end{equation*} It is clear that $g(Y_j,Y_i)=0$ for $1\leq i\leq j-1$. Additionally we have, \begin{align*} \lvert{g(Y_j,Y_j)}\rvert&=\left\lvert{g(X_j,X_j)-\sum_{i=1}^{j-1}\frac{(g(X_j,Y_i))^2}{g(Y_i,Y_i)}}\right\rvert\\ &\geq \lvert{g(X_j,X_j)}\rvert-\sum_{i=1}^{j-1}\frac{\lvert{g(X_j,Y_i)}\rvert^2}{\lvert{g(Y_i,Y_i)}\rvert}>a-n\frac{c^2}{d}=d,\\ \lvert{g(Y_j,X_k)}\rvert&=\left\lvert{g(X_j,X_k)-\sum_{i=1}^{j-1}\frac{g(X_j,Y_i)g(Y_i,X_k)}{g(Y_i,Y_i)}}\right\rvert\\ &\leq \lvert{g(X_j,X_k)}\rvert+\sum_{i=1}^{j-1}\frac{\lvert{g(X_j,Y_i)\rvert} \lvert{g(Y_i,X_k)}\rvert}{\lvert{g(Y_i,Y_i)}\rvert}<b+n\frac{c^2}{d}=c. \end{align*} Since $\lVert{Y_i}\rVert=\sqrt{\lvert{g(Y_i,Y_i)}\rvert}\in\mathfrak{F}(W)$ is non-vanishing we can set $E_i=Y_i/\lVert{Y_i}\rVert\in\mathfrak{X}(W)$ for each $1\leq i\leq n$, and get a desired orthonormal local frame $(E_1,\dots,E_n)$ on $W\ni p$.