First, the paper is Bruedern and Wooley!
Next ... by Cauchy's inequality, one has
$$|\Phi(\alpha)|^2\le 6H\Psi(\alpha),$$
say, where
$$\Psi(\alpha)=\sum_{h\le 6H}\left| \sum_{P<x\le 2P}e(\alpha h(3x^2+3xh+h^2))\right|^2.$$
Thus the mean value in question is bounded above by $6H$ times
$$\int_0^1\Psi(\alpha)|\Phi(\alpha)|^2\, d\alpha ,$$
and by orthogonality that counts the number of solutions of the Diophantine equation
$$3h(x_1-x_2)(x_1+x_2+h)=g_1(3y_1^2+3y_1g_1+g_1^2)-g_2(3y_2^2+3y_2g_2+g_2^2),$$
with $h,g_1,g_2\le 6H$ and $P<x_1,x_2,y_1,y_2\le 2P$. When $x_1=x_2$ and $h$ is free, then given $g_2$ and $y_2$, one finds that $g_1$ and $y_1$ are determined by a divisor estimate. So the number of diagonal solutions is $\ll HP\cdot (HP)^{1+\epsilon}=(HP)^{2+\epsilon}$. Meanwhile, when $x_1\ne x_2$, the right hand side here is also non-zero. Then given a fixed choice for $g_1,g_2,y_1,y_2$, one finds that $x_1-x_2$ and $h$ are determined by a divisor function estimate. There are consequently $\ll (HP)^2\cdot (HP)^\epsilon$ such solutions. Combine these estimates, throw back in the factor $6H$ that arose from applying Cauchy's inequality, and one gets the estimate claimed.