  Globular clusters (GCs) are found in all types of galaxies and harbor some of
the most extreme stellar systems, including black holes that may dynamically
assemble into merging binaries (BBHs). Uncertain GC properties, including when
they formed, their initial masses and sizes, affect their production rate of
BBH mergers. Using the gravitational-wave catalog GWTC-3, we measure that
dynamically-assembled BBHs -- those that are consistent with isotropic spin
directions -- make up ${61^{+29}_{-44}\%}$ of the total merger rate, with a
local merger rate of ${10.9^{+16.8}_{-9.3}}$ Gpc$^{-3}$ yr$^{-1}$ rising to
${58.9^{+149.4}_{-46.0}}$ Gpc$^{-3}$ yr$^{-1}$ at $z = 1$. We assume this
inferred rate describes the contribution from GCs and compare it against the
Cluster Monte Carlo (CMC) simulation catalog to directly fit for the GC initial
mass function, virial radius distribution, and formation history. We find that
GC initial masses are consistent with a Schechter function with slope ${\beta_m
= -1.9^{+0.8}_{-0.8}}$. Assuming a mass function slope of $\beta_m = -2$ and a
mass range between $10^4$--$10^8\,M_\odot$, we infer a GC formation rate at $z
= 2$ of ${5.0^{+9.4}_{-4.0}}$ Gpc$^{-3}$ yr$^{-1}$, or
${2.1^{+3.9}_{-1.7}}\times 10^6\,M_\odot$ Gpc$^{-3}$ yr$^{-1}$ in terms of mass
density. We find that the GC formation rate probably rises more steeply than
the global star formation rate between $z = 0$ and $z = 3$ ({82\%} credibility)
and implies a local number density that is ${f_\mathrm{ev} =
22.6^{+29.9}_{-16.2}}$ times higher than the observed density of survived GCs.
This is consistent with expectations for cluster evaporation, but may suggest
that other environments contribute to the rate of BBH mergers with
significantly tilted spins.
