For giving proper weights to the embedded events one needs to know the pT and y distribution of Upsilons.
For 200 GeV pp, dAu, AuAu a sigma=1.0 estimate was used which is supported by simulations. However, the 500 GeV case can be quite different.

I have ran pythia8176 in order to get the y distribution for 500 GeV p+p collisions.
The tails of the distribution are falling faster than Gaussian, so the fit depends on the range. Since the acceptance dies out with higher rapidity one should mostly concentrate on the middle range. However, fitting wider part of the distribution can give a handle on the error.

As a comparison, the same method was applied in the 193 GeV U+U case (there I chose the range -1..1 since the distribution is somewhat narrower) and sigma=-1.16 was obtained. Fitting the range -2.5..2.5 (omitting only the very tail of the distribution), one gets sigma~1.0. This was used for an error estimation, causing less than 1% error generally (but 2.1% in the lowest pT bin).

In the 500 GeV case:
* fit -2..+2 : sigma = 1.6
* fit -3..+3 : sigma = 1.4
The limits are somewhat arbitrary but the "real" distribution is apparently contained within these two fits. One way for error estimation on the efficiency would be to take 1.5 and then vary it with +-0.1.