Suppose $f$ is an entire function and $f(0)\ne 0$. Let $z_1, z_2,\ldots$ be the zeros of $f$ and $p_1, p_2,\ldots$ be any sequence of nonnegative integers such that $$r\gt 0\implies \sum_{n=1}^\infty \left(\frac{r}{|z_n|}\right)^{p_n+1}\lt\infty.$$ Then there exists an entire function $g$ such that $$f(z)=e^{g(z)}\prod_{n=1}^\infty E_{p_n}\left(\frac{z}{z_n}\right)$$ where $$E_{p}(z)=(1-z)\exp\sum_{v=1}^p \frac{z^v}{v}.$$
If $f$ has an $k$-order zero at $z=0$, the theorem is applied to $f(z)/z^k$.
This can be easily applied to the sine function which has zeros at integer multiples of $\pi$: $$\sin z=z\prod_{n\ne 0}\left(1-\frac{z}{n\pi}\right)e^{\frac{z}{n\pi}}=z\prod_{n=1}^\infty \left(1-\frac{z^2}{n^2\pi ^2}\right).$$
Question
How should this theorem be used when the zeros are "double-indexed"? Let $$f(z)=\theta_1(z|\tau)=2\sum_{n=0}^\infty (-1)^nq^{\left(n+\frac{1}{2}\right)^2}\sin ((2n+1)z).$$ where $q=e^{\pi i \tau}$ with $\Im\tau\gt 0$. Then $f(z)$ is entire and has zeros precisely at $$z=(m+n\tau)\pi$$ where $m,n\in\mathbb{Z}$. What would the Weierstrass factorization of $f(z)$ be?