Where does the really nice '8-dimensional' description of the $E_7$ root system come from? The Wikipedia page on $E_7$ tells me:

Even though the roots span a 7-dimensional space, it is more symmetric and convenient to represent them as vectors lying in a 7-dimensional subspace of an 8-dimensional vector space.
  The roots are all the 8×7 permutations of (1,−1,0,0,0,0,0,0) and all the $\begin{pmatrix}8\\4\end{pmatrix}$ permutations of (½,½,½,½,−½,−½,−½,−½)
  Note that the 7-dimensional subspace is the subspace where the sum of all the eight coordinates is zero. There are 126 roots.

This presentation is indeed more symmetrical than anything 7-dimensional I could have come up with myself. But more importantly it makes computations really easy because, magically, the Killing form is just the restriction of the ordinary inner product on the ambient 8-dimensional space to the 7-dimensional subspace where the roots live. (Wikipedia does not state this explicitly but it can be easily checked from the list of simple roots corresponding to the nodes of the Dynkin diagram that Wikipedia is kind enough to list.)
The question in the title ('where does this come from?') is really two questions:
1) A reference request: if I use this presentation of the root system to simplify my computations, who do I give credit to?
2) A more 'philosophical' question: where does the nice presentation in terms of a bigger space 'come from', mathematically?
The idea of realizing a degree $n$ root system in the $n$-dimensional 'sum of coordinates equals zero'-subspace of an $(n+1)$-dimensional space is of course very familiar: it is how we normally describe the root systems of type $A_n$. 
But in the $A_n$-case the appearance of the extra dimension seems very natural. Thinking about the root system as coming from the Lie algebra $\mathfrak{g} = \mathfrak{sl}_{n+1}$, we can either argue that the $(n+1)$-dimensional space 'surrounding' the $n$-dimensional Cartan subalgebra is `really' the Cartan of the central extension $\mathfrak{gl}_{n+1}$ of $\mathfrak{g}$ or simply accept that everything there is to understand about $\mathfrak{g}$ can be seen inside its 'defining' representation, which happens to be $(n+1)$-dimensional. (I say it a bit sloppy but you hopefully get what I mean.)
It seems that neither of these explanations is available in the $\mathfrak{g} = \mathfrak{e}_7$ case. It definitely does not have a non-trivial 8-dimensional representation (defining or otherwise) and I also never heard about any interesting central extensions.
So is there another explanation of this '$A_n$-like' behavior of $E_7$? The best I could come up with is that the nice 8-dimensional presentation of the $E_7$-root system is 'inherited' from the $A_7$-root system sitting inside of it. But as far as explanations go this feels like cheating since I only found out that there is an $A_7$ sitting inside $E_7$ by looking at the very description of the root system I am trying to explain!
 A: [edited to include more examples of the construction and fix a typo]
This can also be seen purely in terms of lattices or quadratic forms.
As you note, for any $n$ the vectors in the slice $\sum_{i=0}^n x_i = 0$ of
${\bf Z}^{n+1}$ constitute the root lattice $A_n$; so we'll use
the case $n=7$ of the following general picture.
The dual lattice $A_n^*$ consists of the intersection of this slice with the
sublattice of $(n+1)^{-1} {\bf Z}^{n+1}$ with all $x_i \equiv x_j \bmod\bf Z$;
the map taking $(x_0,\ldots,x_n) \in A_n^*$ to the common value of
$((n+1) x_i) \bmod n+1$ descends to an isomorphism
$A_n^*/A_n \cong {\bf Z}/(n+1){\bf Z}$.  Any intermediate lattice
then consists of the intersection of $A_n^*$ with $d^{-1} {\bf Z}^{n+1}$
for some factor $d$ of $n+1$; this lattice is integral iff $d^2 | n+1$,
in which case it is even unless $n$ and $(n+1)/d^2$ are both odd.
In particular, taking $n=7$ and $d=2$ we obtain an even lattice
of discriminant $8/2^2 = 2$ that contains $A_7$ with index $2$;
this must be the lattice $E_7$ (for instance because we can count
${8 \choose 4} = 70$ new roots which together with the $8 \cdot 7 = 56$
roots of $A_7$ give a total of $126$).
The root lattice $E_8$ can be constructed similarly with $(n,d) = (8,3)$
as a lattice intermediate between $A_8^{\phantom*}$ and $A_8^*$,
with index $3$ in both directions.
Other examples: $(n,d) = (2,3)$ gives a lattice isometric with ${\bf Z}^3$;
for $(n,d) = (15,4)$ we get the unique unimodular lattice of rank $15$
with no vectors of norm $1$; for $(17,3)$, one of the two indecomposable
even lattices of rank $17$ and discriminant $2$ (the other one has the
same theta series, and contains $D_{10} \oplus E_7$ with index $2$);
and for $(24,5)$, the Niemeier lattice with root system $A_{24}$.
(If you like coding theory or the Fano plane, you can also use a
symmetrical picture of $E_7$ in ${\bf R}^7$ that consists of
all vectors $2^{-1/2} a$ with $a \in {\bf Z}^7$ such that $a \bmod 2$ is
an element of the $[7,3,4]$ code, or equivalently either zero or the
complement of a line in the Fano plane.  NB $[7,3,4]$ is the dual Hamming code;
using the Hamming code
$[7,4,3]$ itself yields $E_7^*$, while the extended [8,4,4] Hamming code
yields $E_8$.)
A: The span of the roots of E7 is the pentacontihexapentacosiheptacontihexaexon, aka Gosset polytope 2_13 = x3o3o3o *c3o3o3o. That one allows for a complemental compound decomposition into the small petated hexadecaexon = expanded octaexon = x3o3o3o3o3o3x and the hexadecaexon = trirectified octaexon = o3o3o3x3o3o3o. 
In fact, the first one has 56 vertices, the latter one has 70 vertices. If these are compounded in the given orientation wrt. A7 (which most symmetrically can be represented within cartesian coordinates of one dimension plus), then the convex hull of these 126 vertices is nothing than that Gosset polytope 2_13.
You'll find that also on https://bendwavy.org/klitzing/incmats/laq.htm
A: I don't know who found this presentation first but I can imagine that already Cartan knew it since it comes from a symmetric space. More precisely, $\mathfrak g=E_7$ has an involution $\theta$ whose fixed point set is $\mathfrak k=\mathfrak{sl}(8)$. The $-1$-eigenspace of $\theta$, customarily denoted by $\mathfrak p$, is a representation of $\mathfrak k$. In fact, it is the fourth fundamental representation $\wedge^4\mathbb C^8$. So we have a decomposition $E_7=\mathfrak{sl}(8)\oplus\wedge^4\mathbb C^8$.
Now your presentation is clear: the permutations of $(1,-1,0,0,0,0,0,0)$ are the root vectors of $\mathfrak{sl}(8)$ while the vectors $(\pm1/2,\ldots\pm1/2)$ are the weights of $\wedge^4\mathbb C^8$.
Similar games can be played with all symmetric spaces. For example $E_8=\mathfrak{so}(16)\oplus\mathbb C^{128}$ (the latter is the spin representation) is a very popular choice. This leads to the root vectors all permutations of $(\pm1,\pm1,0,0,0,0,0,0)$ and all $(\pm1/2,\ldots\pm1/2)$ with an even number of minus signs.
