As you are only working with 2D physics (no Z velocity), this problem can be greatly simplified. The easy way to do this is to stop thinking about both source and target moving relative to world co-ordinates and to just think of the target moving relative to the source (and keep the source stationary).
Vector TargetInitialPosition = new Vector ( target.X - source.X ,
target.Y - source.Y );
Vector TargetApparentVelocity = new Vector( target.velocityX - source.velocityX ,
target.velocityY - source.velocityY );
Normally, a bullet's velocity would be much higher than the shooter's velocity so it is usually assumed that the bullet is independent but there are occasions where this is not true, such as firing out of a helicopter or fighter jet.
Then we need to work out the bullet velocity:
// Your directional vector MUST be normalized...
Vector BulletVelocity = new Vector( source.directionX * Bullet::StaticSpeed + source.velocityX ,
source.directionY * Bullet::StaticSpeed + source.velocityY );
The problem you're having is that the target has moved by the time the bullet reaches them.
TargetPosition = TargetInitialPosition + TargetApparentVelocity * t
BulletPosition = BulletInitialPosition + BulletVelocity * t
= BulletVelocity * t
and solve for TargetPosition == BulletPosition because then the bullet would have hit the target. Now you have three unknowns and only two equations. We can remove 't' by taking the first order derivative:
TargetInitialPosition + ( TargetApparentVelocity - BulletVelocity ) * t == 0
dV / dt = TargetApparentVelocity - BulletVelocity
Now to hit the target, you'd want
dV/dt == -TargetInitialPosition * k. The constant has to be the same in the X and Y coordinates and is the number of seconds the bullet will take to hit the target.
TargetApparentVelocity.X - BulletVelocity.X == k * -TargetInitialPosition.X
k = ( BulletVelocity.X - TargetApparentVelocity.X ) / TargetInitialPosition.X
TargetApparentVelocity.Y - BulletVelocity.Y == k * -TargetInitialPosition.Y
k = ( BulletVelocity.Y - TargetApparentVelocity.Y ) / TargetInitialPosition.Y
make them equal:
( BulletVelocity.X - TargetApparentVelocity.X ) / TargetInitialPosition.X
= ( BulletVelocity.Y - TargetApparentVelocity.Y ) / TargetInitialPosition.Y
or to expand the variables:
( source.directionX * Bullet::StaticSpeed + source.velocityX - target.velocityX + source.velocityX ) / ( target.X - source.X )
== ( source.directionY * Bullet::StaticSpeed + source.velocityY - target.velocityY + source.velocityY ) / ( target.Y - source.Y )
Then algebra gives you your final equation:
source.directionY = ( target.velocityY * ( source.X - target.X ) - 2 * source.velocityY * ( source.X - target.X ) + ( Bullet::Speed * source.directionX + 2 * source.velocityX - target.velocityX ) * ( source.Y - target.Y ) ) / ( Bullet::Speed * ( source.X - target.X ) )
The next part is messy and it's up to you how you want to implement it in your code, but we just substitute this in and normalize the vector.
sqrt( source.directionX ^ 2 + source.directionY ^ 2 ) == 1
You end up with an equation with just one unknown (source.directionX), and you can solve it for directionX then substitute back in to get directionY.
I haven't tested any of this code and feel free to point out any methematical misstakes I've made, but the theory should be sound :).