Ok, you can, after some calculations see that
the absolute value of your sum is less than or equal to
$$\sum_{n=-\infty}^\infty e^{\pi(-tn^2 + 2n|z|)}$$
where $t = im(\tau) >0.$
The part of the sum where $n<0$ will be less than some finite constant
depending only on $\tau,$ so we only need to worry about $n>0.$
(This you need to prove.)
Ok, so, divide that sum into two parts:
Let $N$ be the least integer greater than $4|z|/t$
(where $t = im(\tau) >0.$)
If $n\geq N$ we have $nt/4 > |z|.$
This will imply that $-tn^2 + 2n|z| \leq -n^2t/2.$
Use this on the second sum, $N...\infty$ and use some other estimations on the first part.
Edit: Continuation:
First, you should see that
$$\sum_{n=-\infty}^0 e^{\pi(-tn^2 + 2n|z|)} \leq
\sum_{n=-\infty}^0 e^{-\pi t n^2}
$$
and this clearly converges to something that only depends on $\tau.$
(The exponent will diverge to $-\infty$ rapidly, so the sum must converge.)
Now, if $n\geq N$ then
$$\sum_{n=N}^\infty e^{\pi(-tn^2 + 2n|z|)}\leq \sum_{n=N}^\infty e^{\pi(-n^2 t/2)}
\leq \sum_{n=0}^\infty e^{\pi(-n^2 t/2)}
$$
which converges to a finite number independent of $z$.
Now, you only need to do something about the last part:
$$\sum_{n=0}^{N-1} e^{\pi(-tn^2 + 2n|z|)}$$
(Hint: use $N-1< 4|z|/t$ which follows from the conditions above).
Edit: Continuation:
Ok, choose $r \in \mathbb{R}$ so that $\theta(r,\tau) \neq 0.$
By quasi-periodicity, we have for integers $m$ that
$$\theta(r+m\tau,\tau) = e^{-i\pi m^2\tau - 2im \pi r}\theta(r,\tau)$$
so by taking absolute values on each side, we see that
$$|\theta(r+m\tau,\tau)| = e^{t\pi m^2}|\theta(r,\tau)|$$
where $t = im(\tau)>0.$
Hence,
$$\lim_{m \rightarrow \infty} \frac{|\theta(r+m\tau,\tau)|}{e^{t\pi m^2}}
= |\theta(r,\tau)| > 0.$$
Thus, the order is at least 2.
$$\wp_{\tau}(z) = \frac1{z^2} + \sum_{(n,m) \ne (0,0)}\frac1{(z+n\tau+m)^2}-\frac1{(n\tau+m)^2}$$ is $1$ periodic and analytic on $\Im(z) \in (0,|\Im(\tau)|)$ so it has a Fourier series valid for $\Im(z) \in (0,|\Im(\tau)|)$ which is a Laurent series for $\wp_{\tau}(\frac{\log u}{2i\pi})$ valid for $|e^{2i \pi u}| \in (1,e^{|\Im(\tau)|})$
The Fourier coefficients are found from the pole expansion of inverse of trigonometric functions. For $\Im(z) > 0$ $$\sum_{m}\frac{1}{(z+m)^2}= \frac{d}{dz} \frac1{1-e^{2i \pi z}}= \sum_{k \ge 0} (2i\pi k) e^{2 i \pi kz}$$
For $\Im(z) < 0$ $$\sum_{m}\frac{1}{(z+m)^2}= \frac{d}{dz} \frac1{1-e^{2i \pi z}}= \sum_{k \ge 1} (-2i\pi k) e^{-2 i \pi kz}$$
Best Answer
First of all, there is an error in the edition of Stein and Shakarchi that you are using. In later editions, the first two derivatives of $\Theta(z|\tau)$ has been replaced by the first three derivatives of $\Theta(z|\tau)$.
Here's my attempt to answer your question. I am not sure whether it is satisfactory because the expression of $c_\tau$ also includes three constants, $e_1$, $e_2$, and $e_3$, defined as
$\wp(1/2) = e_1$, $\wp(\tau/2) = e_2$, $\wp(1/2+\tau/2) = e_3$.
For completeness, I will also briefly go through the proof of equality itself.
From Corollary 1.5, we know the LHS, hereafter denoted as $L(z)$, is an elliptic function of order 2 with periods 1 and $\tau$, and with a double pole at $z=1/2+\tau/2$. (From here on I will denote $1/2+\tau/2$ as $z_0$.)
Furthermore, multiplying $L(z)$ by $(z-z_0)^2$ and letting $z\rightarrow z_0$, we can see the coefficient of the double pole $\frac{1}{(z-z_0)^2}$ is exactly 1.
We also know that $\wp(z-z_0)$ is an elliptic function of order 2 with periods 1 and $\tau$, and with a double pole at $z=z_0$. Its coefficient of the double pole is also 1. (This can be easily seen from the definition of $\wp(z)$.)
Therefore, $L(z)-\wp(z)$ is an elliptic function that is entire. It must be a constant. This establishes the desired equality.
To get $c_\tau$, we take the derivative of both sides with respect to $z$, and square both sides. We use Theorem 1.7 from the Chapter 9.
$(L')^2 = (\wp')^2 = 4(\wp - e_1)(\wp - e_2)(\wp - e_3) = 4(L-c_\tau-e_1)(L-c_\tau-e_2)(L-c_\tau-e_3)$.
Setting $z=z_0$, all $\Theta(z|\tau)$ vanish. What is left is an equation of the first three derivatives of $\Theta(z|\tau)$ at $z=z_0$, $c_\tau$, and the three constants $e_1$, $e_2$, and $e_3$.