Monte carlo integration in spherical coordinates

I was playing around with writing a code for Montecarlo integration of a function defined in spherical coordinates. As a first simple rapid test I decided to write a test code to obtain the solid angle under an angle $\theta_m$. For two random number $u$ and $v$ in $[0,1)$. I generate an homogeneous random sampling of the spherical angle using $\phi=2\pi u$ and $\theta =\arccos(1-2v)$.

For N generated points, I have M points for which $\theta < \theta_m$. My first idea was that since I have an homogeneous sampling I should have obtained the correct solid angle $\Omega=2 \pi (1-\cos (\theta_m))$ simply as $4\pi\times \frac{M}{N}$. Actually it looks like I the correct result comes out only if I use:

$$\Omega=\sum_{i=1}^M \frac{4\pi}{N} 2 cos(\theta_i)$$

I can not see the reason why this should be correct. The probability distribution function in $\theta$ is $PDF=\frac{1}{2} \sin(\theta)$ so I would rather expect I should normalize each point of the sum by this function but this doesn't works. What am I doing wrong and how could I justify the cosinus? Many thanks!

-
Hmm.. This was double-posted to the stats.SE site. Which site wants to keep it? – whuber May 22 '11 at 22:00
Can you explain what's going on with the definition of $\theta$? What is $a$, and what is the rationale behind the argument of cosine? (Right now it takes values in [-1,1], so $\theta \in [a\cos(1),a]$.) – Gerben May 23 '11 at 0:00
I'msorry it was a missprint, I had to write $\theta=arccos(1-2v)$. I have corrected the question now. – astolfo May 23 '11 at 6:02

You can considerably simplify this. The angle $\phi$ doesn't occur anywhere, so you're not actually using $u$ or $\phi$, so we can ignore them. The angle $\theta$ only occurs in the form $\cos\theta$, so it makes sense to switch to $z=\cos\theta$ instead.
Restated in this way, your question is whether generating (presumably uniformly distributed) random numbers $v$ in $[0,1)$ and calculating $z=1-2v$ will lead to $z>z_m$ in a fraction $M/N=2\pi(1-z_m)/(4\pi)=(1-z_m)/2$ of cases. This linear relationship is clearly correct, so your formula is OK and it seems there must be something wrong in your program instead if this is not working.