First a little ugly math that you're already using in your code.
Define x and y are bits with probability of being 1 of X = p(x=1), Y = p(y=1) respectively.
Then we have that
p( x & y = 1) = X Y
p( x | y = 1) = 1 - (1-X) (1-Y)
p( x ^ y = 1) = X (1 - Y) + Y (1 - X)
Now if we let Y = 1/2 we get
P( x & y ) = X/2
P( x | y ) = (X+1)/2
Now set the RHS to the probability we want and we have two cases that we can solve for X
X = 2 p // if we use &
X = 2 p - 1 // if we use |
Next we assume we can use this again to obtain X in terms of another variable Z...
And then we keep iterating until we've done "enough".
Thats a bit unclear but consider p = 0.375
0.375 * 2 = 0.75 < 1.0 so our first operation is &
0.75 * 2 = 1.5 > 1.0 so our second operation is |
0.5 is something we know so we stop.
Thus we can get a variable with p=0.375 by X1 & (X2 | X3 )
The problem is that for most variables this will not terminate. e.g.
0.333 *2 = 0.666 < 1.0 so our first operation is &
0.666 *2 = 1.333 > 1.0 so our second operation is |
0.333 *2 = 0.666 < 1.0 so our third operation is &
etc...
so p=0.333 can be generated by
X1 & ( X2 | (X3 & (X4 | ( ... ) ) ) )
Now I suspect that taking enough terms in the series will give you enough accuracy, and this can be written as a recursive function. However there might be a better way that that too... I think the order of the operations is related to the binary representation of p, I'm just not sure exactly how... and dont have time to think about it deeper.
Anyway heres some untested C++ code that does this. You should be able to javaify it easily.
uint bitsWithProbability( float p )
{
return bitsWithProbabilityHelper( p, 0.001, 0, 10 );
}
uint bitsWithProbabilityHelper( float p, float tol, int cur_depth, int max_depth )
{
uint X = randbits();
if( cur_depth >= max_depth) return X;
if( p<0.5-tol)
{
return X & bitsWithProbabilityHelper( 2*p, 0.001, cur_depth+1, max_depth );
}
if(p>0.5+tol)
{
return X | bitsWithProbabilityHelper( 2*p-1, 0.001, cur_depth+1, max_depth );
}
return X;
}