My ultimate goal is to calculate how much mass contains a 3D (dx,dy) cut of a spherically symmetric nucleus for which I assume a radius-dependent mass distribution. I believe that solving a simplified problem - calculate surface of 2D (dx) strip of a disc - would allow me to generalize to the true one.
For the case of a uniform mass distribution on a disc of radius R, it is easy to solve the problem via MC integration:
- generate phi from uniform[0,2pi] distribution and r from f® = r distribution for r[0…R]
If I am correct points (x(phi,r), y(phi,r)) would be uniformly distributed on a disc in x,y space
- calculate x(phi, r) values and count what fraction falls into dx interval
Now I am wondering if it is possible to solve it via numerical, non-MC integration, in a faster (in terms of CPU usage) and more stable way than the MC solution. For that I wanted to define a RooFit PDF(x,y) representing uniform distribution on a disc, that could be integrated in rectangular range (dx, dy). I tried the following:
RooFormulaVar phi(“phi”, “TMath::ATan2(y, x)”, RooArgSet(x, y));
RooFormulaVar r(“r”, “TMath::Sqrt(x * x + y * y)”, RooArgSet(x, y));
RooUniform phiPdf(“phiPdf”, “phiPdf”, phi);
RooGenericPdf rPdf(“rPdf”, “r”, RooArgList®);
RooProdPdf model(“model”, “model”, RooArgSet(phiPdf, rPdf));
But I don’t know how to limit “r” to “< 5”. I tried setting limits on y:
RooFormulaVar ylo(“ylo”,"-sqrt(25 - xx)",x) ;
RooFormulaVar yhi(“yhi”,"sqrt(25 - xx)",x) ;
In both cases model.createHistogram(“x,y”)->Draw(“colz”) shows something that is not a uniform disc.
Another idea is to define the model in (phi,r) space - that shows no problem - and integrate in parametrized ranges. However here defining range limits is a bit complex and while (most probably) doable also for the original case of volume integration over sphere, it is not so elegant as if it was possible to define the model in Cartesian coordinates.
Can anyone advise me on how to solve the integration issue - via MC as other solutions are bound to be less efficient in CPU, via some clever definition of the PDF, or via parametrized ranges?