Are you trying to calculate interaction coupling constants as well? 

In simple cubic lattices, which are the simplest, you can take one of the directions, say x, calculate the Bloch functions, and extract the band structure epsilon_n(x). At low energy, it has the form epsilon_0(x) = E_0 -2Jcos(ka) -2J_2 cos(2ka)+...

So you should fit the numerically obtained band structure to that series expansion and you'll get J, J_2, etc.  Then the tight-binding model comes out straightforward, without interactions.

If you need to fit the interactions too, then it depends on the density of your system. Say, again, in a sc lattice, with V(x,y,z) = V(x) + V(y) +V(z) (so that all tunnelings are identical), the low-energy limit of your scattering amplitude is pure s-wave. You can fit that to the "Hubbard" U interaction by calculating the scattering amplitude in the tight binding model and then fitting one to the other. This works for gas-like (metallic) phases at low densities. At high densities you are better off doing strong-coupling perturbation theory.

Another way of fitting interactions, especially for high-density cases, or more than one band (in the tight-binding limit), is to Taylor expand a single lattice well, 

V(x,y,z) = V_0 + V1x^2 + V2y^2+V3y^2 +...

And calculate energy shifts due to interactions. Then you have a straight-away mapping to the on-site energy shifts. 

Hope this helps