This is fundamentally no more difficult than understanding how quantum mechanics describes particle motion using plane waves. If you have a delocalized wavefunction $\exp(ipx)$ it describes a particle moving to the right with velocity p/m. But such a particle is already everywhere at once, and only superpositions of such states are actually moving in time.
Consider
$$\int \psi_k(p) e^{ipx - iE(p) t} dp$$
where $\psi_k(p)$ is a sharp bump at $p=k$, not a delta-function, but narrow. The superposition using this bump gives a wide spatial waveform centered about at x=0 at t=0. At large negative times, the fast phase oscillation kills the bump at x=0, but it creates a new bump at those x's where the phase is stationary, that is where
$${\partial\over\partial p}( p x - E(p)t ) = 0$$
or, since the superposition is sharp near k, where
$$ x = E'(k)t$$
which means that the bump is moving with a steady speed as determined by Hamilton's laws. The total probability is conserved, so that the integral of psi squared on the bump is conserved.
The actual time-dependent scattering event is a superposition of stationary states in the same way. Each stationary state describes a completely coherent process, where a particle in a perfect sinusoidal wave hits the target, and scatters outward, but because it is an energy eigenstate, the scattering is completely delocalized in time.
If you want a collision which is localized, you need to superpose, and the superposition produces a natural scattering event, where a wave-packet comes in, reflects and transmits, and goes out again. If the incoming wavepacked has an energy which is relatively sharply defined, all the properties of the scattering process can be extracted from the corresponding energy eigenstate.
Given the solutions to the stationary eigenstate problem $\psi_p(x)$ for each incoming momentum $p$, so that at large negative x, $\psi_p(x) = exp(ipx) + A \exp(-ipx)$ and $\psi_p(x) = B\exp(ipx)$ at large positive x, superpose these waves in the same way as for a free particle
$$\int dp \psi_k(p) \psi_p(x) e^{-iE(p)t}$$
At large negative times, the phase is stationary only for the incoming part, not for the outgoing or reflected part. This is because each of the three parts describes a free-particle motion, so if you understand where free particle with that momentum would classically be at that time, this is where the wavepacket is nonzero. So at negative times, the wavepacket is centered at
$$ x = E'(k)t$$
For large positive t, there are two places where the phase is stationary--- those x where
$$ x = - E'(k) t$$
$$ x = E_2'(k) t$$
Where $E_2'(k)$ is the change in phase of the transmitted k-wave in time (it can be different than the energy if the potential has an asymptotically different value at $+\infty$ than at $-\infty$). These two stationary phase regions are where the reflected and transmitted packet are located. The coefficient of the reflected and transmitted packets are A and B. If A and B were of unit magnitude, the superposition would conserve probability. So the actual transmission and reflection probability for a wavepacket is the square of the magnitude of A and of B, as expected.
The square modulus of the transmission coefficient ($|T|^2$) is the transmission probability, and you have the data to calculate that for your wave-packets (i.e. it's just $|\psi|^2$ as you surmised, if you normalized your wavefunction). You can compare that to the expected transmission coefficient for plane waves with the same average energy as your wave-packet. I don't know how well the comparison will be since you have uncertain energy, but you can use wider wave-packets and at least see if it seems to be converging to the plane wave solution.
Incidentally it is a good idea to check the norm of your wave function, it can help determine if anything is going wrong with your numerical simulation.
PS: This http://arxiv.org/abs/quant-ph/0301114 might be of some value. It is an explicit solution for a wave-packet on a square barrier (i.e. like your problem, but only one 'edge'.)
Best Answer
Sketched proof:
The first question the reader should ask him/herself is:
This is answered in e.g. this Phys.SE post. In particular, note that $e^{+ikx}$ and $e^{-ikx}$ are a right- and left-mover, respectively.
Secondly, to conserve probability over time, impose that the $S$-matrix should be unitary. This is e.g. done in my Phys.SE answer here. Unitarity implies with $G=0$ that $$|A|^2~=~|B|^2+|F|^2.$$
Thirdly, interpret $\frac{|B|^2}{|A|^2}$ and $\frac{|F|^2}{|A|^2}$ as probabilities for reflection and transmission, respectively, thereby adding up to 100 %.