Here a general way how to compute Massey Products on nice spaces.
The first step would be to realize the space as the geometric realization of a simplicial set. I used the simplicial set given by the following picture. It has one 0-simplex $v$, its 1-simplices are $a,b,c,e,f$ and its 2-simplices are $A,B,C,D$ and their boundaries are indicated in the picture. Note that the arrowtips indicate the orientations on the edges. Note that each 2-simplex really is totally ordered, so we really have a simplicial set here.
Also note that there are identifications in this picture, for example the edge a appears four times and so on.
Now let us have a look at the dual of the simplicial chain complex. Let the dual bases be denoted by $v^*,a^*,...$. Let us read off the boundaries of those basis vectors and note that $\mathbb{F}_2$-coefficients lets us ignore signs. We have:
$$\begin{eqnarray*}
d(v^*)&=& 0,\\
d(a^*)&=& A^*+B^*+C^*+D^*,\\
d(b^*)&=& A^*+D^*\\
d(c^*)&=& A^*+B^*\\
d(e^*)&=& B^*+C^*\\
d(f^*)&=& C^*+D^*\\
\end{eqnarray*}$$
To find the DGA-structure on that cochain complex, we still have to write down the multiplication. In degree $0$ everything is a multiple of one ($v^*$) and everything above degree $2$ is zero. So it suffices to compute the products of all the basis vectors in degree one.
To determine this, we need to recall the definition of the cup product. To compute say $a^*\cdot b^*$,we have to look at the 2-simplices and add up the duals of those which start with $a$ and end with $b$. So we get
$$A^* = c^*a^*, \quad B^* = a^*e^*, \quad C^* = f^*a^*, \quad D^* =a^*b^*,$$
and all other products of those base vectors are zero (especially also the ones above after switching the order). Let us first compute the cohomology ring.
Since the cohomology in concentrated in degrees $\le 2$ and $H^0=\mathbb{F}_2$, we just have to compute products of classes from $H^1$.
The classes $[a^*+c^*+f^*]$ and $[b^*+c^*+e^*+f^*]$ form a basis of $H^1$ (Those are alled $a$ and $b$ in the question).
$$\begin{eqnarray*}
(a^*+c^*+f^*)^2&=& A^*+C^* = d(c^*+e^*),\\
(b^*+c^*+e^*+f^*)^2 &=& 0,\\
(a^*+c^*+f^*)(b^*+c^*+e^*+f^*)&=& B^*+D^* = d(e^*+f^*),\\
(b^*+c^*+e^*+f^*)(a^*+c^*+f^*)&=& A^*+C^* = d(c^*+e^*).\\
\end{eqnarray*}$$
So we see that the product in cohomology of all those classes is zero (They are boundary). Especially the Massey Product $\langle a,a,b \rangle$ really lives in $H^2$ and not in a quotient of it (We quotient out the zero module).
Following the definition of Massey product, we seee that it is represented by the cycle
$$ (c^*+e^*)(b^*+c^*+e^*+f^*)+(a^*+c^*+f^*)(e^*+f^*)= B^*$$
and $0\neq [B^*] \in H^2$.
One small final remark: That cochain complex also is an algebra over the sequence operad from McClure-Smith and thus we could also compute the Steenrod squares on $H^*$ in a similar fashion.
