Define $b_{n,m}(x)=\sum_{i=1}^m \frac{b_{n-1,i}(1)}{4i^2-1}x^i$, and let $b_{1,m}(x)=\sum_{i=1}^m \frac{x^i}{4i^2-1}$. Does we have a closed form for $b_{n,m}(1)$? I just know this from these functions: $$f(x)=\prod_{n=1}^\infty(1+\frac{x}{4n^2-1})=1+b_{1,\infty}(1)x+(b_{1,\infty}(1)^2-b_{2,\infty}(1))x^2+(b_{1,\infty}(1)^3-2b_{1,\infty}(1)b_{2,\infty}(1)+b_{3,\infty}(1))x^3+...$$ $$f(1)=\frac{2\times 2 \times 4 \times 4 \times 6 \times 6 ...}{1\times 3\times 3 \times 5 \times 5 \times 7 ...}=\frac{\pi}{2}$$
EDIT: Since if we have a function that satisfies the conditions $f(0)=1, f(a_n)=0$, for $n \in [1,\infty)$, then we have that the function f is equal to a infinite product and a infinite sum: $$f(x)=\prod_{n=1}^\infty (1-\frac{x}{a_n})=1-x\sum_{n=1}^\infty \frac{1}{a_n}+x^2\sum_{n=1}^\infty \sum_{k=n+1}^\infty \frac{1}{a_n a_k}-x^3\sum_{n=1}^\infty \sum_{k=n+1}^\infty \sum_{j=k+1}^\infty\frac{1}{a_n a_k a_j}+...$$ So if we have that $a_n=1-4n^2$, then: $$f(x)=\prod_{n=1}^\infty \frac{4n^2-1+x}{4n^2-1}=1+x\sum_{n=1}^\infty \frac{1}{4n^2-1}+...=1+xb_{1,\infty}(1)+...$$
I will only concern myself with the values of the requested quantities at $x=1$, so I will identify henceforth $b_{n,m}(1)\equiv b_{n,m}$. First, note that all quantities $b_{n,\infty}$ have a closed form in terms of even powers of $\pi$. This can be seen by rewriting the function $f$ as
$$f(x)=\prod_{n=0}^\infty \frac{\left(n+1-\sqrt{\frac{1-x}{4}}\right)\left(n+1+\sqrt{\frac{1-x}{4}}\right)}{(n+1/2)(n+3/2)}$$
which by a very general theorem is equal to
$$f(x)=\frac{\Gamma(1/2)\Gamma(3/2)}{\Gamma\left(1-\sqrt{\frac{1-x}{4}}\right)\Gamma\left(1+\sqrt{\frac{1-x}{4}}\right)}=\frac{\sin\frac{\pi}{2}\sqrt{1-x}}{\sqrt{1-x}}$$
(the reflection formula for the gamma function has been used to reach the second equality). Hence, one can find expressions order by order for $b_{n,\infty}$. The expansion to the first few orders is
$$f(x)=1 + \frac{x}{2} + (\frac{3}{8} - \frac{\pi^2}{32}) x^2 + \frac{1}{32} (10 - \pi^2) x^3+\mathcal{O}(x^4)$$
From this we conclude that for example
$$b_{1,\infty} = 1/2~~,~~ b_{2,\infty} = \frac{\pi^2-4}{32}~~,~~ b_{3,\infty} =1/16$$
Also, some of the more general quantities $b_{n,m}(1)$ seem to have closed forms. For example the series for $b_{1,m}$ telescopes:
$$b_{1,m}=\frac{1}{2}-\frac{1}{4m+2}$$
On the other hand, after some simplifications
$$b_{2,m}=\frac{1}{8}\frac{2m}{2m+1}+\frac{1}{4}\sum_{k=1}^m\frac{1}{(2k+1)^2}=\frac{1}{8}\frac{2m}{2m+1}+\frac{1}{4}\left(H^{(2)}_{2m+1}-\frac{H^{(2)}_m}{4}-1\right)$$
has a closed form in terms of 2nd order harmonic numbers. Given the final expressions, it is natural to conjecture that all above sums can be reduced to linear combinations of the sums $$\sum {1/(4k^2-1)},\sum {1/(2k+1)^{2j}}, j\leq \lfloor n/2\rfloor $$ and hence are linear combinations of $\{1,\zeta(2j)\}$. It would be interesting to see if there is some explicit formula that generates the above coefficients explicitly.
Update: I have managed to find explicit expressions for the series coefficients of $f(x)$ around $x=0$. The expansion reads
$$\frac{\sin \frac{\pi}{2}\sqrt{1-x}}{\frac{\pi}{2}\sqrt{1-x}}=\sum_{n=0}^\infty \left(\frac{\pi x}{2}\right)^{n}\frac{j_n(\pi/2)}{2^n n!}$$
where $j_n$ is the n-th order spherical Bessel function. As a matter of fact, these functions are expressible in terms of elementary functions and their values at $\pi/2$ are polynomials in even powers of $\pi$ as we will show below. Using the simple formula for the spherical Hankel function
$$h_n^{(1)}(x)=e^{ix}\sum_{k=0}^n \frac{i^{k-n-1}}{z^{k+1}}\frac{(n+k)!}{2^k k!(n-k)!}$$
Taking the real part of this, we can write
$$[x^n]f(x)=\frac{1}{4^n n!}\sum_{k=0}^{\lfloor n/2 \rfloor}(-1)^k\pi^{2k}\frac{(2n-2k)!}{(2k)!(n-2k)!}$$
which proves part of the conjecture. This however isn't of great help in calculating exact expressions for $b_{n,\infty}$, as the structure of the coefficients of $f$ in terms of those is rather complicated.
What remains is work to be done directly on the $b_{n,m}$'s. In fact, it can be shown with some algebra that the genearating function
$$F(x,y)=\sum_{n,m\geq 1}b_{n,m}x^n y^m$$
obeys a second order differential equation in the $y-$variable:
$$\left[\left(y\frac{\partial}{\partial y}\right)^2-\frac{1}{4}\right]\Bigg[(1-y)(F(x,y)-xH(y))\Bigg]=xF(x,y)$$
where $H(y)=\frac{1}{2(1-y)}+\frac{1}{4\sqrt{y}}\ln\left(\frac{1-\sqrt{y}}{1+\sqrt{y}}\right)$
and the initial conditions are
$$F(x,0)=0~~,~ \partial_y F(x,0)=\frac{1}{1-x/3}$$