How quickly will billiard trajectories cluster? Suppose you launch $n$ point-particles on 
distinct reflecting nonperiodic billiard trajectories 
inside a convex polygon. Assume they all have the same speed.
Define an $\epsilon$-cluster as a configuration of the particles
in which they all simultaneously lie within a disk of radius $\epsilon$.
It is my understanding that Poincaré's Recurrence Theorem
implies that at some time after launch, the particles
will form an $\epsilon$-cluster somewhere.
(Please correct me if I am wrong here, in which case the remainder
is moot.)
Picturesquely, if I sit in my office long enough, all the air molecules
will cluster into a corner of the room. :-)
The reason I specify that the trajectories be distinct is to
exclude
the particles being shot in a stream all on the same trajectory.
The reason I specify nonperiodic is to exclude sending
the particles on parallel periodic trajectories whose length
ratios are rational, in which case no clustering need occur.
My question is:

How long must one wait for an $\epsilon$-cluster to occur?

Essentially I am seeking a quantitative version of
Poincaré's Recurrence Theorem, quantitative enough
to actually make a calculation.
I would like to put a number of years to the air-molecule example
(air molecules move perhaps 700 mph or 300 m/s).
It could serve as a useful pedagogical anecdote.
I found a beautiful paper that should help
me answer this question:
Benoit Saussol,
"An Introduction to Quantitative Poincaré Recurrence in Dynamical Systems,"
Reviews in Mathematical Physics, Volume 21, Issue 08, pp. 949-979 (2009).
But I am having difficulty making the leap from the abstract
theorems to an explicit calculation.
Any help or additional pointers would be appreciated!
Addendum. Vaughn's analysis, although leaving a few loose ends (as he notes), largely
answers my question.  Thanks to all for the astute comments and responses!
 A: I'll offer a few partial answers, which may eventually lead to a complete answer.
As observed in Saussol's paper (Theorem 3, Kac's lemma), if you have an ergodic invariant measure $\mu$ then the mean return time to a set $A$ is equal to $1/\mu(A)$.  (Note, however, that this is the mean return time for trajectories that begin in $A$, rather than for an arbitrary initial condition, for which you'd need the hitting time).
So we need to know what the ergodic invariant measure $\mu$ is.  Many topological dynamical systems have lots of invariant measures lying around, in which case it's a non-trivial task to pick the one you're really interested in.  However, if we assume that your convex polygon satisfies the Veech dichotomy (see, for example, this article by John Smillie and Barak Weiss), then it's true that the flow in every direction is either periodic or uniquely ergodic.  Regular polygons satisfy the Veech dichotomy, so in that setting your non-periodicity condition would be enough to guarantee that Lebesgue measure is the only invariant measure.
Of course, Lebesgue measure is always invariant for the billiard flow, so the main thing we need this for is the fact that Lebesgue measure is ergodic under these assumptions.  (There may be an easier way to get ergodicity, but this is the first thing that came to my mind.)
If we write $X_\theta$ for the phase space of the billiard flow with angle $\theta$, and let $\theta_1,\dots,\theta_n$ be the angles that your $n$ particles are launched at, then the above discussion implies that the billiard flow on each $(X_{\theta_n},\text{Leb})$ is ergodic.  The phase space for your entire system is $\prod_n X_{\theta_n}$, and the direct product of the $n$ different Lebesgue measures is an invariant measure for your overall system.
What we'd like to do next is say that this measure is actually ergodic.  Unfortunately, it's not quite that simple, since the direct product of two ergodic systems need not be ergodic (just consider $R_\alpha\times R_\alpha$, the direct product of two circle rotations by irrational multiples of $\pi$).  So to say that Lebesgue measure is ergodic for your whole system requires something more, which is where your condition that the trajectories be distinct ought to come in.  I'm not sure exactly how this step should go, but you should be able to get something along the lines of "for almost every set of angles $(\theta_1,\theta_2,\dots,\theta_n)$, Lebesgue measure is ergodic for the whole system".  Billiard flows have something to do with interval exchange transformations (IETs), and so this paper by Jon Chaika may well have the result that's needed at this point.
Once you know that (product) Lebesgue measure is invariant for the whole system, you're in the clear:  given a ball of radius $\epsilon$, the (normalised) Lebesgue measure of that ball in each $X_\theta$ is $\pi \epsilon^2 / C$, where $C$ is the area of the polygon, and so the set in $\prod_n X_{\theta_n}$ corresponding to those configurations for which all particles lie in this ball is just $(\pi \epsilon^2 / C)^n$.  The inverse of this is your expected return time, provided we have appropriate hypotheses to justify all the above steps.  (Expected hitting time might be a different story, I'd need to think a bit.)
A: In general, we cannot make an explicit estimate. For example, you specified that the trajectories should not follow periodic orbits which do not intersect each other. However, by choosing the trajectories sufficiently close to periodic trajectories we can make them shadow those trajectories to within $\varepsilon$ until time $T$, and hence before time $T$ they will not come closer to each other than the separation of the periodic orbits minus $2\epsilon$. So no universal upper bound exists.
