Update:
As proposed by @dmckee, I added equation numbers and improved the display of some equations.
The answer by @Trimok inspired me to look at coordinate systems which are not specific to the detector locations, but more general, and this indeed made things clearer since $\delta(\cos \theta + 1)$ is a $\delta(\cos \theta + \cos 0)$ in disguise.
We have a relation between a Dirac delta of a ray in spherical coordinates, $\delta^2$, and a Dirac delta of a point, $\delta^3 = \delta(r-r_0)\delta(\theta-\theta_0)\delta(\phi-\phi_0)/(r^2 \sin \theta)$, such that
$\delta^2 = \int_0^\infty \delta^3 \operatorname{d}r_0 = \delta(\theta-\theta_0) \delta(\phi-\phi_0) / (r^2 \sin\theta) \tag{E1}$
(Using this definition, one can easily show that the sifting property of $\delta^2$ holds and the result is a line integral of some $f$.)
The major problem (of understanding) which I am having is now reduced to this: We have two expressions for $\delta^3$, which are $\delta(r-r_0)\delta(\theta-\theta_0)\delta(\phi-\phi_0)/(r^2 \sin \theta) \tag{E2}$ and $\delta(r-r_0)\delta(\theta-\theta_0)\delta(\phi-\phi_0)/(r_0^2 \sin \theta_0) \tag{E3}$ These are totally symmetric in $x$ and $x_0$ (note that $\delta$ is an even function).
However, for $\delta^2$ we have $\delta(\theta-\theta_0)\delta(\phi-\phi_0)/(r^2 \sin \theta) \tag{E4}$ and $\delta(\theta-\theta_0)\delta(\phi-\phi_0)/(r^2 \sin \theta_0) \tag{E5}$ which are not symmetric any longer. So what makes the $r$ stand out compared to $r_0$?
One explanation I worked out with the help of a colleague is directed towards the fact that there is no $r_0$ after integration, so having any $r_0$ would "overdetermine" the set of rays through the origin. Does anyone have a better (for example, more intuitive) explanation why the symmetry is lost here?
Thank you!
Original question:
This is a long one. I have been working on this for quite a while now, so I summarized most of what I found out. I am stuck at a certain point, however, and it would be great if someone could help me out.
Assume the following situation: Two perfectly collinear photons $p_1$ and $p_2$ are emitted in the origin of some Cartesian coordinate system. We consider two infinitesimal detectors $d_1$ and $d_2$, with surfaces $\operatorname{d}\!A_1$ and $\operatorname{d}\!A_2$ directed towards the origin so that we do not have to deal with angles. The solid angles of these detectors be $\operatorname{d}\!\Omega_1$ and $\operatorname{d}\!\Omega_2$.
I now express some basic detection probability densities, where $p_id_i$ means that "$p_i$ is emitted towards $d_i$".
$\frac{\operatorname{d}\!P(p_1d_1)}{\operatorname{d}\!\Omega_1}=\frac{1}{4\pi} \tag{1}$
$\frac{\operatorname{d}\!P(p_1d_1)}{\operatorname{d}\!A_1}=\frac{1}{4\pi (r_1)^2} \tag{2}$
$\frac{\operatorname{d}\!P(p_2d_2)}{\operatorname{d}\!\Omega_2}=\frac{1}{4\pi} \tag{3}$
$\frac{\operatorname{d}\!P(p_2d_2)}{\operatorname{d}\!A_2}=\frac{1}{4\pi (r_2)^2} \tag{4}$
I think there is no doubt that these probabilities are well-defined and the integrated detection probabilities in a spherical $4\pi$-detector equals 1, no matter what its radius is. So a photon is surely emitted towards a detector surrounding the origin.
Now, for some reason we want to calculate $P(p_2d_2 | p_1d_1)$, the conditional probability of $p_2$ emitted towards $d_2$, given that $p_1$ has been emitted towards $d_1$. Due to collinearity of the photons, the individual events are not independent: $P(p_2d_2 | p_1d_1) \neq P(p_2d_2)$.
I now propose that
$\frac{\operatorname{d}\!P(p_2d_2 | p_1d_1)}{\operatorname{d}\!\Omega_2} = \frac{\delta(\cos \theta_2(d_1)+1)}{2\pi}, \tag{5}$
where $\theta_2(d_1)$ is the polar angle of the position of $d_2$ expressed in spherical coordinates with polar axis from the origin through $d_1$. $\delta$ is the Dirac function, so obviously, this expression equals 0 if $\cos \theta_2(d_1) \neq -1$, that is, when $\theta_2(d_1) \neq \pi$.
We can verify that the integral of this probability over a spherical $4\pi$-detector $d_2$ equals 1: with $\operatorname{d}\!\Omega_2 = \operatorname{d}\!\cos \theta_2(d_1) \operatorname{d}\!\phi_2(d_1)$, we find
$\int \frac{\operatorname{d}\!P(p_2d_2 | p_1d_1)}{\operatorname{d}\!\Omega_2} \operatorname{d}\!\Omega_2 = \int \frac{\delta(\cos \theta_2(d_1)+1)}{2\pi} \operatorname{d}\!\Omega_2 \tag{6}$
$= \int_0^{2\pi} \frac{1}{2\pi} \operatorname{d}\!\phi_2(d_1) \cdot \int_0^\pi \delta(\cos \theta_2(d_1)+1) \operatorname{d}\!\cos \theta_2(d_1) = 1. \tag{6}$
So $p_2$ is still surely emitted towards a detector surrounding the origin, no matter where $p_1$ is emitted towards.
Further, we would like to verify that the law of total probability is fulfilled by integrating $P(p_1d_1 \wedge p_2d_2) = P(p_2d_2 | p_1d_1) \cdot P(p_1d_1)$ over a spherical $4\pi$-detector $d_1$. This works regardless of which formulation of $\operatorname{d}\!P(p_1d_1)$ (with respect to $\operatorname{d}\!\Omega_1$ or $\operatorname{d}\!A_1$) we use -- I chose $\operatorname{d}\!A_1 = (r_1)^2 \operatorname{d}\!\cos \theta_1(d_2) \operatorname{d}\!\phi_1(d_2)$ here.
Note, however, that integration over $d_1$ changes the polar axis used in the formulation of $\operatorname{d}\!P(p_2d_2 | p_1d_1)$ above, so we need a different, but equivalent formulation. This is more than a change of coordinate system (from one polar axis to another), because $\theta_2(d_2)$ would not help us a lot (it's 0). Instead, one can use
$\frac{\operatorname{d}\!P(p_2d_2 | p_1d_1)}{\operatorname{d}\!\Omega_2}=\frac{\delta(\cos \theta_1(d_2)+1)}{2\pi}, \tag{7}$
which uses both a different variable in the argument of the $\delta$ as well as a different coordinate system: $\theta_1(d_2)$ is the polar angle of detector $d_1$'s location in spherical coordinates with polar axis from the origin through $d_2$. Note that both expressions of $\frac{\operatorname{d}\!P(p_2d_2 | p_1d_1)}{\operatorname{d}\!\Omega_2}$ are symmetric in the indices 1 and 2.
Finally, we arrive at
$\int \frac{\operatorname{d}\!P(p_1d_1 \wedge p_2d_2)}{\operatorname{d}\!A_1 \operatorname{d}\!\Omega_2} \operatorname{d}\!A_1 = \int_0^{2\pi} \int_0^\pi \frac{1}{4\pi (r_1)^2} \cdot \frac{\delta(\cos \theta_1(d_2)+1)}{2\pi} (r_1)^2 \operatorname{d}\!\cos \theta_1(d_2) \operatorname{d}\!\phi_1(d_2) \tag{8}$
$= \frac{1}{4\pi} \int_0^{2\pi} \frac{1}{2\pi} \operatorname{d}\!\phi_1(d_2) \cdot \int_0^\pi \delta(\cos \theta_1(d_2)-1) \operatorname{d}\!\cos \theta_1(d_2) = \frac{1}{4\pi} = \frac{\operatorname{d}\!P(p_2d_2)}{\operatorname{d}\!\Omega_2} \tag{8}$
So, the law of total probability is fulfilled: the probability of $p_2$ being emitted towards $d_2$ can be partitioned into all conditional probabilities of $p_2$ being emitted towards $d_2$ given that $p_1$ has been emitted towards $d_1$, weighted by the probability that $p_1$ is actually emitted towards $d_1$.
So far, so good. Congratulations for reading up to this point :)
As expressed above (and verifiable relatively easily), it does not matter which formulation of $\operatorname{d}\!P(p_1d_1)$ I use -- both $\operatorname{d}\!\Omega_1$ and $\operatorname{d}\!A_1$ work fine. However, when I try to base my calculations on $\frac{\operatorname{d}\!P(p_2d_2)}{\operatorname{d}\!A_2}$ instead of $\frac{\operatorname{d}\!P(p_2d_2)}{\operatorname{d}\!\Omega_2}$, hence, I try to calculate $\frac{\operatorname{d}\!P(p_2d_2 | p_1d_1)}{\operatorname{d}\!A_2}$ instead of $\frac{\operatorname{d}\!P(p_2d_2 | p_1d_1)}{\operatorname{d}\!\Omega_2}$, I arrive at the following problem.
So I started with
$\frac{\operatorname{d}\!P(p_2d_2 | p_1d_1)}{\operatorname{d}\!A_2}=\frac{\delta(\cos \theta_2(d_1)+1)}{2\pi (r_2)^2}, \tag{9}$
and (6) is still easily verified (integrating over $\operatorname{d}\!A_2$) since the two additional $(r_2)^2$ simply cancel out.
However, working out (8) is now more complex when it comes to changing from one variable/coordinate system to another. I figured that
$\frac{\operatorname{d}\!P(p_2d_2 | p_1d_1)}{\operatorname{d}\!A_2}=\frac{\delta(\cos \theta_1(d_2)+1)}{2\pi (r_2)^2} \tag{10}$
is the expression which fulfills total probability (I can show this, if necessary), but now I wonder why this is not symmetric with the first expression in the indices 1 and 2 (as it was above). Specifically, I wonder where the connection to a Dirac of a ray along Oz in spherical coordinates is, which is $\frac{\delta(\cos \theta - 1)}{2\pi r^2}$ (compare http://www.physicsforums.com/showthread.php?t=215263)
More generally, I am working on expressing the probability density in 3D, so I expect to obtain something like $\operatorname{d}\!\Omega_2 \cdot \delta_{L(x_1, x_2)}(x_0)$ or $\operatorname{d}\!A_2 \cdot \delta_{L(x_1, x_2)}(x_0)$, where $L(x_1, x_2)$ is the line segment from $x_1$ and $x_2$, and $x_0$ is the point of emission.
So my specific questions are:
Why are the expressions involving $\operatorname{d}\!A_2$ not symmetric, while those involving $\operatorname{d}\!\Omega_2$ are?
Specifically, where is the connection to the Dirac of a ray in spherical coordinates, which involves an $r^2$ and which I would expect to be changed when converting from one coordinate system to another.
How (else) can I obtain a 3D probability density involving the three locations, which can be used to integrate over all three volumes in order to describe the number of coincidences detected from emissions from a volume in two finite detectors?
