That's a trivial simplification of the [Johnson-Lindenstrauss lemma][1] . 
The matrix $X$ can be seen as a set of $q$ points in $p$ dimensions. Let us first prove the following result, than your question will follow directly with a union bound: (by the way, I'm going to assume that your matrix $G$ is normalized)

> **Theorem**
> Let $x \in \mathbb{R}^p$ and assume that the entries of the matrix $\tilde{G} \in \mathbb{R}^{s \times p}$ are drawn at random from a Gaussian distribution. Let $G = \frac{1}{\sqrt{s}}\tilde{G}$. Then $$ \mathbb{P}\left( |\|Gx\|_2^2 - \|x\|_2^2| \leq \varepsilon \|x\|_2^2 \right) \geq 1-2e^{-(\varepsilon^2-\varepsilon^3)s/4} $$


The proof uses a small lemma about the $\chi^2$ distribution which can be proved via Markov's inequality, after applying an exponential. 

> **Proof**
> First of all, remember that $\mathbb{E}[\|Gx\|_2^2] = \|x\|_2^2$.
It is also worth mentioning that the individual components of the image vector $(Gx)_j/\|x\|_2$ for $1 \leq j \leq s$ are independent Normal variable with variance $1/\sqrt{s}$. As a consequence $\|Gx\|_2^2/\|x\|_2^2$ behaves as a $\chi^2/s$ random variable with $s$ degrees of freedom. 
Thus it follows, using classical results on the $\chi^2$ distribution, that $$ \mathbb{P}\left(\|Gx\|_2^2 \geq (1+\varepsilon) \|x\|_2^2 \right) = \mathbb{P}(\sum_{j=1}^s \frac{(Gx)_j}{\|x\|_2} \geq 1 + \varepsilon) = \mathbb{P}(\frac{1}{s}\chi^2 \geq 1 + \varepsilon) = \mathbb{P}(\chi^2 \geq (1+\varepsilon)s) \leq e^{-(\varepsilon^2-\varepsilon^3)s/4}. $$
Similarly, the following bounds hold
$$ \mathbb{P}\left(\|Gx\|_2^2 \leq (1-\varepsilon) \|x\|_2^2 \right) =  \mathbb{P}(\chi^2 \leq (1-\varepsilon)s) \leq e^{-(\varepsilon^2-\varepsilon^3)s/4}. $$
Combining both bounds together yields the expected result. 

The result you asked for can be obtained by considering the norm $\|X\|_2^2 = \sum_{j = 1}^s\|x_j\|_2^2$ and using a union bound, playing with the constance


  [1]: https://en.wikipedia.org/wiki/Johnson%E2%80%93Lindenstrauss_lemma