The point is that electromagnetism is a gauge theory. This is the most important feature to understand about it. Everything else in e+m is an application of this principle ;-).
Observable is only the electromagnetic field with the electric components \vec{E} and magnetic components \vec{B}. The homogeneous Maxwell equations (in Heaviside-Lorentz units with c=1),
\vec{\nabla} \times \vec{E}+\partial_t \vec{B}=0, \quad \vec{\nabla} \cdot \vec{B}=0
imply the existence of the four-vector potential with components (\phi, \vec{A}) such that
\vec{E}=-\partial_t \vec{A}-\vec{\nabla} \phi, \quad \vec{B}=\vec{\nabla} \times \vec{A}.
If for given (\vec{E},\vec{B}) you have found scalar and vector potential then any other set of potentials, connected with the old one via a gauge transformation of the form
\phi'=\phi+\partial_t \chi, \quad \vec{A}'=\vec{A}-\vec{\nabla} \chi
with an arbitrary scalar field \chi also leads to the same electromagnetic field, i.e., the one set of potentials is as good as any other connected to it by a gauge transformation.
You've seen this in your example of a homogeneous electrostatic field. Usually, of course, you have to solve the Maxwell equations for some given physical arrangement of charges and current densities. The inhomogeneous Maxwell equations read (in the vacuum!)
\vec{\nabla} \cdot \vec{E}=\rho, \quad \vec{\nabla} \times \vec{B}-\partial_t \vec{E}=\vec{j}.
Expressing the field components in terms of the potentials this gives
-\vec{\nabla} \cdot (\partial_t \vec{A}+\vec{\nabla} \phi)=\rho, \quad \vec{\nabla} (\vec{\nabla} \cdot \vec{A})-\Delta \vec{A}+\partial_t (\partial_t \vec{A}+\vec{\nabla} \phi)=j.
This is a pretty complicated set of coupled partial differential equations, but the point is that you have some freedom to choose the potentials since they are only defined modulo a gauge transformation. Thus you can simplify your work by constraining the potentials to fix the corresponding arbitrary gauge function \chi.
A very convenient choice is the Lorenz gauge condition
\partial_t \phi+\vec{\nabla} \cdot \vec{A}=0.
It's easy to show that the pretty complicated equations above decouple to inhomogeneous wave equations for each component of the potential separately:
\Box \phi=\rho, \quad \Box \vec{A}=\vec{j},
where \Box=\partial_t^2-\Delta is the D'Alembert operator.