Let me say a couple of things which might help here:
(1) The information you give here does not uniquely determine a joint distribution, and in fact information of this sort will never uniquely determine a joint distribution. For instance, if I say that I have two variables X ~ N(0, 1) and Y ~ N(0, 1) with a given correlation between X and Y, it does not follow that (X, Y) follows a multivariate normal distribution -- there are infinitely many joint distributions with this property. (This is in some part a consequence of the fact that covariance only gives us information up to the second order.) You will have to make up a joint density with the desired properties; it will not be uniquely determined.
(2) There are different ways that covariance could be introduced into this model, and they will give the model different theoretical properties. For instance, there could be heavy correlation between the Bernoulli variables that decide which distribution to sample from and no correlation between the actual samples. Conversely the Bernoulli variables for X and Y could be totally uncorrelated, and some or all of the normal/binomial distributions could be correlated to one another. One could also envision a mixture of both scenarios. If you have any more fundamental information about your processes (which I'm assuming you do -- how else would you know that it's a mixture of normal and binomial samples?) then you should use that to decide here.
The main problem here is coming up with a joint PDF -- you are in reality building a mathematical model here, and that's as much art as science. Once you have the joint distribution, there are well-established ways to sample from it. The most straightforward conceptually is to sample x_1, then take a conditional probability distribution for (x_2, x_3, ... x_n), then repeat. There are more sophisticated ways to deal with this problem as well, and in any case most of these algorithms are built into MATLAB or whatever you're going to use to do your sampling already.