The einbein is just a Lagrange multiplier that enforces the mass-shell condition. This is most clearly seen in Hamiltonian formulation:
\begin{align}
S[x,p,e] &= \int d\sigma\,\, p_\mu \dot{x}^\mu - e\,\phi(x,p)
\\
0=\frac{\delta S}{\delta e} &\implies \phi=0,
\\
0=\frac{\delta S}{\delta x},\,\, 0=\frac{\delta S}{\delta p} &\implies \dot{f}=e\{f,\phi\},\,\,\{x^\mu,x^\nu\}=0,\,\,\{x^\mu,p_\nu\}=\delta^\mu_\nu,\,\,\{p_\mu,p_\nu\}=0.
\end{align}
$\phi(x,p)=\frac{1}{2}(g^{\mu\nu}(x) p_\mu p_\nu +m^2)$ is a Hamiltonian constraint that functions as both the time-evolution generator and a constraint on the phase space. It is also the generator of reparametrization gauge transformation, so we also call it by "gauge generator." When $m=0$, we get a massless relativistic particle.
Besides, phase space is just another name (physicists' name) for symplectic manifold. Constrained mechanics on phase space is a matter of describing a symplectic submanifold in a symplectic manifold. Hence we not only consider the constraint $\phi(x,p)$ but its pair constraint $\chi(x,p)$. For instance it is customary to let the gauge-fixing condition depend on the worldline time parameter and set $\chi(x,p,\sigma) = -u_\mu x^\mu - \sigma$ for a constant timelike one-form $u_\mu$. That is, we parameterize the worldline with the coordinate time in the co-moving frame of $u^\mu$. (In the "covariant phase space" setting, we set $\chi(x,p)=-u_\mu x^\mu$ and let $\chi(x,p)=0$ describe a time slice that the particle's initial position $x^\mu(0)$ takes its values.) By introducing this "gauge-fixing condition," the total number of constraints becomes an even number; recall that symplectic manifolds are even-dimensional.
In general, Lagrange multipliers can be determined by studying the following three conditions.
\begin{equation}
(1)\,\,\, \phi=0
\qquad
(2)\,\,\, \dot{\phi}=0
\qquad
(3)\,\,\, \dot{\chi}=0
\end{equation}
And in the case of $(2)\, \dot{\phi}=0$, it holds identically because $\dot{\phi} = e \{\phi,\phi\}=0$. So it does not give new information. Usually, $e(\sigma)$ is determined by enforcing (3). This gives the celebrated Dirac brackets. For example you may find it helpful to have a look at appendix A of arXiv:2102.07063.
The method (3) works for both massive and massless particles, so my answer to your question will be: use not only the mass-shell constraint but also a "gauge-fixing constraint" (chosen with your taste)!
The thing is, for massive particles, the einbein can also be determined by using (1) (which gives $e(\sigma) = \frac{1}{m} \sqrt{-\dot{x}^2(\sigma)}$), but for massless particles enforcing the mass-shell constraint gives nothing. So I think there is some kind of subtle difference between massive and massless particles, and this might be the confusing point addressed by your question.