<p>Let <span class="math-container">$W_{\text e}$</span> denote the camera sensor responsivity (i.e. importance) and <span class="math-container">$h_j$</span> and <span class="math-container">$W_{\text e}^{(j)}$</span> denote the  image reconstruction filter and sensor response associated with the <span class="math-container">$j$</span>th pixel, respectively, i.e. <span class="math-container">$$W_{\text e}^{(j)}=h_jW_{\text e}.\tag1$$</span> The measurement of the <span class="math-container">$j$</span>th pixel is given by <span class="math-container">\begin{equation}\begin{split}\Phi_j&amp;:=\sigma_M({\rm d}x_0)\int\sigma_M({\rm d}x_1)g(x_0\leftrightarrow x_1)W_{\text e}^{(j)}(x_1\to x_0)L_{\text o}(x_1\to x_0)\\&amp;=\int\lambda({\rm d}x)h_j(x_1\to x_0)f(x),\end{split}\tag2\end{equation}</span> where <span class="math-container">$g$</span> denotes the geometry term, <span class="math-container">$\lambda$</span> is the infinite-product area measure on the path space and <span class="math-container">$f$</span> is the measurement contribution function, i.e. <span class="math-container">$$f(x)=g(x_0\leftrightarrow x_1)W_{\text e}(x_1\to x_0)\prod_{i=2}^kg(x_{i-1}\leftrightarrow x_i)f_{\text s}(x_i\to x_{i-1}\to x_{i-2})L_{\text e}(x_k\to x_{k-1}),\tag3$$</span> if <span class="math-container">$x$</span> is a path of length <span class="math-container">$k\in\mathbb N$</span>.</p>

<blockquote>
  <p>Usually, in the context of Metropolis Light Transport, the scalar contribution function <span class="math-container">$p$</span> is defined to be the luminance of <span class="math-container">$f$</span>. So, by definition, <span class="math-container">$p(x)=0$</span> whenever <span class="math-container">$(x_0,x_1)$</span> doesn't belong to <span class="math-container">$$B:=\left\{(x_0,x_1):W_{\text e}(x_1\to x_0)\ne 0\right\}.$$</span> Now let <span class="math-container">$$B_j:=\left\{(x_0,x_1):W_{\text e}^{(j)}(x_1\to x_0)\ne 0\right\}.$$</span> We clearly have <span class="math-container">$$B=\bigcup_jB_j.\tag4$$</span> However, in light of the way samples are splat to image (see <a href="http://www.pbr-book.org/3ed-2018/Light_Transport_III_Bidirectional_Methods/Metropolis_Light_Transport.html#fragment-Splatbothcurrentandproposedsamplestomonofilm-0" rel="nofollow noreferrer">http://www.pbr-book.org/3ed-2018/Light_Transport_III_Bidirectional_Methods/Metropolis_Light_Transport.html#fragment-Splatbothcurrentandproposedsamplestomonofilm-0</a> or <a href="https://pdfs.semanticscholar.org/a386/55316ee7f438ba115b38e8d0b1410c691a26.pdf" rel="nofollow noreferrer">https://pdfs.semanticscholar.org/a386/55316ee7f438ba115b38e8d0b1410c691a26.pdf</a> and <a href="https://computergraphics.stackexchange.com/q/9589/9254">my other question</a>), it seems like it is assumed that the union in <span class="math-container">$(4)$</span> is disjoint. Is this really guaranteed to be the case?</p>
</blockquote>

<p>Taking a look in Section 2.6.4, on page 51, of <a href="https://backend.orbit.dtu.dk/ws/portalfiles/portal/51112612/phd240_awk_vers2.pdf" rel="nofollow noreferrer">this PhD thesis</a>, this seems at least to be the case when the camera is a pinhole camera. The authors write: <span class="math-container">$W_{\text e}^{(j)}$</span> <em>is non-zero only for the set of directions within the pyramid formed by the pinhole and the four corners of the corresponding pixel</em> <span class="math-container">$j$</span>.</p>

<p>This seems to indicate that (at least for the pinhole camera) the "support" of <span class="math-container">$W_{\text e}^{(j)}$</span> is limited to the image rectangle region occupied by the <span class="math-container">$j$</span>th pixel. But is this really true for <em>any</em> filter <span class="math-container">$h_j$</span> or only for a special one?</p>

CC BY-SA 4.0
2020-02-21T12:25:16.040Z
The $B_j$ are not necessarily disjoint - consider having an aperture with some area. Furthermore, a reconstruction filter $h$ in theory (I am considering what PBRT does here), can have an effect over the whole screen if desired. Consider a Laplace or Poisson reconstruction problem from sparse samples (the Dirichlet data). The latter is even the case if you consider gradient domain path tracing.
lightxbulb
CC BY-SA 4.0
2020-02-21T14:28:29.257Z
@lightxbulb You're right. Meanwhile, I've figured the source of my understandings out: https://computergraphics.stackexchange.com/a/9598/9254.
0xbadf00d
CC BY-SA 4.0
2020-02-22T08:34:43.493Z
On second thought - the $B_j$ are disjoint. I was thinking in terms of directions: (x_1 - x_0), in which case they would not have been. But as a tuple, they are, since one of the points (depending on your convention) is on the film, and those sets (of the points corresponding to a specific pixel) are disjoint. Basically pixels do not overlap. Now I am curious what the result will be of having overlapping pixels.
lightxbulb
CC BY-SA 4.0
2020-02-22T13:26:16.420Z
@lightxbulb I need to think about your comment, but let me ask you a simple question. Say each $h_j$ is a box filter with radius $\frac12$ centered at the pixel center. Won't this yield that $h_j$ (considered as a function on the raster space) is simply the indicator function $1_{S_j}$ of the region $S_j$ of the image plane (in raster space) occupied by the $j$th pixel?.
0xbadf00d
CC BY-SA 4.0
2020-02-22T14:11:06.100Z
That is correct if you assume your pixels are of size 1. In the PBRT implementation this disjointedness probably also holds - you'll have to look into the internals to see if this is the case. They have a convolution of the normalized filter $\frac{f}{\sum_i f_i}$ with $wL$. If they make sure that the filter cutoff radii cover only a single pixel then the sets are disjoint. One can however pick a global filter - then the above simply does not hold. Ultimately, one wants to compute the integral within each pixel - so you do not even need to consider a filter in the continuous case.
lightxbulb
CC BY-SA 4.0
2020-02-22T14:13:23.107Z
In the discrete case the filter becomes useful as a means to make up for the lack of infinitely many estimates over the image. I believe the two things should be considered separately - think of the filtering as a numerical integration rule and you won't be that far off the mark.
lightxbulb
CC BY-SA 4.0
2020-02-22T14:24:07.357Z
@lightxbulb In raster space, the "pixels" are always of size $1$. My point is that I think that PBRT is using the Box filter with radius $\frac12$ implicitly in their implementation of (M)MLT without stating this at any point. Another point which is worrying me in the path tracing implemention: If we use a primary state space sampler which doesn't produce uniformly distributed samples (like the Halton sampler, which is the default one used by PBRT), doesn't this introduce bias in the final estimate?
0xbadf00d
CC BY-SA 4.0
2020-02-22T14:38:47.293Z
Depends on your convention/implementation - pixels can have any size you would like, as long as you account for this properly (I am not sure whether pbrt uses size 1). You don't have to think too hard - just check PBRT's film functions. As for Halton - it's deterministic - then taking the expectation of the error obviously does nothing and you get that it's biased, sure. However, low-discrepancy samplers have a different property that guarantees bounds on the error and the equivalent of consistency also. Look into Koksma-Hlawka's inequality.
lightxbulb
CC BY-SA 4.0
2020-02-23T16:01:13.210Z
@lightxbulb I'm sorry to bother you again (it's the last time on this, I promise), but the whole scaling thing is continuing to confuse me. Can you tell me what they mean by the "image radiance function" $L$ here http://www.pbr-book.org/3ed-2018/Monte_Carlo_Integration/Bias.html#?
0xbadf00d
CC BY-SA 4.0
2020-02-23T16:01:27.350Z
@lightxbulb I really got problems to understand how precisely the image function $\Phi$ on the (continuous) raster space $S$ is defined, since the radiance arriving at an point on $S$ (which is infinitely small) should be $0$. There is clearly a map $\psi$ from the path space $E$ to $S$ (such that $\psi(x)$ only depends on $x_0,x_1$) and we should clearly have $\int_B\Phi=\int_{\psi^{-1}(B)}f\:{\rm d}\mu$, where $f$ is the measurement contribution function and the first integral is integrated wrt the Lebesgue measure on $S$.
0xbadf00d
CC BY-SA 4.0
2020-02-23T17:24:45.447Z
That's the radiance function assuming a pinhole camera, then $L(x,y) = L((c_x,c_y) \leftarrow (x,y))$, where $c$ is the point of the aperture. Then there is a unique correspondence between a ray direction and a point on the film. At the end of the day, you have discrete pixels, so you want to integrate the radiance coming from all directions at a pixel point, and integrate over the area of the pixel:
$$I_{i,j} = \int_{x_i}^{x_{i+1}}\int_{y_j}^{y_{j+1}}\int_{\Omega}W(x,y,\omega)L(x,y,\omega)\,d\omega\,dx\,dy$$
As you will notice there is no reconstruction filter - only sensor function $W$.
lightxbulb
CC BY-SA 4.0
2020-02-23T17:29:16.600Z
@lightxbulb But the reconstruction filter is "built into" your $W$, isn't it? So $W$ depends on $i,j$ or am I missing something?
0xbadf00d
CC BY-SA 4.0
2020-02-23T17:30:02.127Z
It is not built into it, and has no reason to be. I suggest you check out the CG discord so we could actually chat in real time.
lightxbulb
CC BY-SA 4.0
2020-02-23T17:41:43.073Z
@lightxbulb Everything I'm asking is related to my scaling problem in my other question nhttps://computergraphics.stackexchange.com/q/9589/9254. There's obviously a crucial thing I'm missing. My understanding (which seems to be wrong) is that the measurement of the $j$th pixel is $\Phi_j=\int(h_j\circ\psi)f\:{\rm d}\lambda$, $h_j:S\to[0,\infty)$ being the filter for pixel $j$, $\psi$ the map from path space $E$ to raster space $S$, $f:E\to[0,\infty)^3$ being the measurement contribution function, $\lambda$ the infinite-area measure.
0xbadf00d
CC BY-SA 4.0
2020-02-23T17:41:59.770Z
@lightxbulb Let's only consider paths of length $k$ (so $E$ and $\lambda$ replaced accordingly) and say we generate samples $X_n\sim q\lambda$ ($q$ being any bidirectional path tracing strategy; for simplicity, assume path tracing if you like). Then the estimate for $\Phi_j$ should simply be $\frac1n\sum_{i=1}^n\frac{(h_j\circ\psi)f}q(X_i)$. Now PBRT is implicitly assuming a box filter with radius $r=1/2$. So, $h_j$ is simply the indicator function of the $j$th pixel region $S_j$ in raster space. That's fine.
0xbadf00d
CC BY-SA 4.0
2020-02-23T17:42:10.100Z
@lightxbulb Now PBRT is summing $\sum_{j\in J}\sum_{i=1}^n1_{S_j}(\psi\circ X_i)\frac fq(X_i)$ (so each summand contributes to precisely one pixel, since the $S_j$ are disjoint). Now we should be left with multiplying each pixel value by $\frac1n$, but PBRT multiplies by $\frac{|J|}n$ ($|J|$ being the number of pixels which is equal to the area of the raster space). It cost me a while, but since PBRT takes $n$ to be of the form $n=\text{samplesPerPixel}\cdot|J|$, $\frac{|J|}n$ happens to be equal to $\text{samplesPerPixel}$.
0xbadf00d
CC BY-SA 4.0
2020-02-23T17:42:15.447Z
@lightxbulb So, where does the need for the multiplication by $|J|$ came from? Sure, $h_j$ needs to be normalized, but since the measure of $S_j$ is $1$, there is nothing to normalize here ... (Sorry for extreme long comment. I'll delete them later).
0xbadf00d
CC BY-SA 4.0
2020-02-23T18:11:38.683Z
This seems to be purely PBRT specific, and as I mentioned I haven't looked into the details of its code in a while. However you mentioned that $n = spp \cdot |J|$, and that the final result is multiplied by $\frac{|J|}{n} = \frac{1}{spp}$. Then everything should be fine, right? What is the issue? Did I maybe misunderstand some part? I highly recommend joining the computer graphics discord so it will be easier to clarify things, rather than producing many comments that are somewhat unrelated to the question.
lightxbulb
CC BY-SA 4.0
2020-02-23T18:49:00.567Z
@lightxbulb I will take a look at the discord thing. No, it's not PBRT specific. When I divide by $n$ then the produced image is black (since obviously $n$ is the wrong factor to divide by). But when I multiply by $|J|/n$ (as PBRT does), the output is correct. I simply don't get why we need to multiply by $|J|/n$ (and not by $1/n$). So, I guess my formula for $\Phi_j$ or $h_j$ is wrong.
0xbadf00d
CC BY-SA 4.0
2020-02-23T18:53:46.533Z
But you just said that $n = spp \cdot |J|$, and multiplying by $\frac{|J|}{n}$ is the same as multiplying by $\frac{1}{spp}$. Which is exactly what you want.
lightxbulb
CC BY-SA 4.0
2020-02-23T19:00:23.117Z
@lightxbulb I said that I figured out (from testing) that multiplying by $|J|/n$ is the correct thing to do, but I have no theoretical justification where the factor $|J|$ comes from. The theoretical Monte Carlo estimate of $\Phi_j=\int(h_j\circ\psi)f\:{\rm d}\lambda$ is $\frac1n\sum_{i=1}^n\frac{(h_j\circ\psi)f}q(X_i)$. So, *it seems like* either my formula for $\Phi_j$ is wrong or I need to "normalize" $h_j$ by multiplying it with $\frac1{|J|}$ (but since $h_j$ should be equal to $1_{S_j}$ and the region $S_j$ of the $j$th pixel has measure $1$, that doesn't make sense to me).
0xbadf00d
CC BY-SA 4.0
2020-02-23T19:17:09.620Z
But your formula has $n$ be $spp$, and in the implementation $n$ was $spp \cdot |J|$, no?
lightxbulb
CC BY-SA 4.0
2020-02-23T19:20:53.727Z
@lightxbulb My formula works for every $n$. So, in the implementation I *theoretically* need to divide by $n=spp\cdot|J|$. Either my $h_j$ or my $\Phi_j$ is wrong.
0xbadf00d
CC BY-SA 4.0
2020-02-23T19:26:09.820Z
The $n$ in your formula needs to be the number of samples, it cannot be arbitrary. That is, it is exactly $spp$. And since in pbrt$ n= spp \cdot |J|$ it obviously requires a multiplication by $|J|$ in the numerator.
lightxbulb
CC BY-SA 4.0
2020-02-23T19:26:12.437Z
@lightxbulb This is the essential part of the code I'm working with: https://coliru.stacked-crooked.com/view?id=7f55a67c5e1dab10. The code uses the suggested number of samples $n=spp\cdot|J|$ and divides by $spp$ at the end of the `render` function. It works, but I don't understand why. From my understanding it should be `m_camera->film->WriteImage(pbrt::Float{ 1 } / m_n);` instead.
0xbadf00d
CC BY-SA 4.0
2020-02-23T19:27:11.750Z
@lightxbulb Sure, $n$ is the number of samples. What I meant was that my $n$ can be any number of samples; it doesn't need to be of the specific form $n=spp\cdot|J|$.
0xbadf00d
CC BY-SA 4.0
2020-02-23T19:30:18.903Z
It seems we're misunderstanding each other. PBRT can also use any number of samples, it's just that in your formula your $n$ variable is not their $n$ variable, it's their $spp$ variable. Why do you believe you should divide by the film area?
lightxbulb
CC BY-SA 4.0
2020-02-23T19:32:21.953Z
@lightxbulb No, that's not correct. When we talk about the number of generated samples, we mean the number of generated paths $X_i$ and they generate $n = spp\cdot|J|$ many; not only $spp$. Please take a look at the code.
0xbadf00d
CC BY-SA 4.0
2020-02-23T19:34:48.453Z
But this is over the whole screen, since for every pixel you get $spp$ samples. So the weight per pixel is $\frac{1}{spp}$.
lightxbulb
CC BY-SA 4.0
2020-02-23T19:39:56.630Z
@lightxbulb What you are doing now is arguing with intuition. And I agree to you; that makes sense. However, the reason why every pixel is getting $spp$ samples is the multiplication by the filter $(h_j\circ\psi)(X_i)$. But, once again, the correct Monte Carlo estimate is $\frac1n\sum_{i=1}^n\frac{(h_j\circ\psi)f}q(X_i)$. From my understanding, $h_j=1_{S_j}$, but it seems like it should be $|J|1_{S_j}$ instead.
0xbadf00d
CC BY-SA 4.0
2020-02-23T19:45:26.457Z
I am not arguing from intuition. I am arguing based on the math that I wrote several replies back. You reconstruct each pixel separately - you do not have an estimator over the whole screen - that would make no sense unless you have a single pixel. So if you want to argue that $|J|$ has to be in there you need to motivate that. This is precisely the reason why I told you you should separate the filter from the estimate. In the original problem there is no filter - a filter is required only after you discretise the problem to have finitely many samples.
lightxbulb
CC BY-SA 4.0
2020-02-23T19:57:45.230Z
@lightxbulb Okay, let me continue to denote the region of the $j$th pixel in raster space by $S_j$. Assume a pinhole camera with eye $x_0$. You agree to me that the measurement $\Phi_j$ of the $j$th pixel is $\Phi_j=\int_{D_j}\sigma^\perp(x_0,{\rm d}\omega)W_e(x_0,\omega)L_i(x_0,\omega)$, where $\sigma_M^\perp$ is the projected solid angle kernel on the scene surfaces $M$ and $D_j$ the set of directions corresponding to the $j$th pixel, right? (My $\Phi_j$ is your $I_{x,y}$ for $j=(x,y)$ and my $S_j$ is your $[x_i,x_{i+1}]\times[y_i,y_{i+1}]$.)
0xbadf00d
CC BY-SA 4.0
2020-02-23T20:03:20.207Z
@lightxbulb Using the map $\psi:E\to S$, this can be written as $\Phi_j=\int\sigma_M({\rm d}x_1)(1_{S_j}\circ\psi)(x)g(x_0\leftrightarrow x_1)W_e(x_1\to x_0)L_o(x_1\to x_0)=\int\lambda({\rm d}x)(1_{S_j}\circ\psi)(x)f(x)$. The estimator for this is $\frac1n\sum_{i=1}^n\frac{(1_{S_j}\circ\psi)f}q(X_i)$ which is the same as before, since we assumed the filter $h_j$ is equal to $1_{S_j}$.
0xbadf00d
CC BY-SA 4.0
2020-02-23T20:05:35.963Z
Here $n$ is $spp$ in PBRT, that is the only thing you keep believing is not the case for some reason. I'll make it simple for you: here's the problem you want to solve: $I_{i,j} = \int_{x_i}^{x_{i+1}}\int_{y_j}^{y_{j+1}}\int_{\Omega}W(x,y,\omega)L(x,y,\omega)\,d\omega\,dx\,dy$
An estimator for it is:
$$\frac{1}{spp}\sum_{k}^{spp}W(X_k, Y_k, Z_k) L(X_k, Y_k, Z^k)$$
lightxbulb
CC BY-SA 4.0
2020-02-23T20:06:21.857Z
@lightxbulb Let me stress that $U_i:=\psi(X_i)$ are the projections of the sample paths $X_i$ to the raster space and the sample sequence may vary over the whole screen. What you write corresponds to what one is usually doing ordinary path tracing (more generally, in what's been done by the `pbrt::SampleIntegrator`). In path tracing we generate a sample sequence $(X_n)_{n\in\mathbb N}$ for each pixel separately. It seems like you're confusing these approaches.
0xbadf00d
CC BY-SA 4.0
2020-02-23T20:11:20.563Z
@lightxbulb Your estimator is only correct if you fix $j$ and use it for a sequence $(X^{(j)}_n)_{n\in\mathbb N}$ of paths such that $\psi(X^{(j)}_n)\in S_j$ for all $n$. This is what `pbrt::SampleIntegrator` is doing. But in MLT only *one* path sequence $(X_n)_{n\in\mathbb N}$ is generated and $\psi(X_n)$ varies over $S=\biguplus_jS_j$.
0xbadf00d
CC BY-SA 4.0
2020-02-23T20:12:03.133Z
Your reconstruction is still performed per pixel. I believe your whole confusion arises from the fact that you have merged the reconstruction process with the actual problem you are trying to solve - you're basically modelling the problem wrong, and confusing yourself further with PBRT's implementation. The reconstruction process is a pure 2D problem. It's signal reconstruction from sparse data. In pbrt they have decided to perform the reconstruction per pixel. Example - signal reconstruction through radial basis functions, interpolation, approximation, PDEs etc. Anyways I have to go, try to
lightxbulb
CC BY-SA 4.0
2020-02-23T20:12:57.080Z
figure out what I mean. If you still haven't figured it out I can try to make some time tomorrow to chat with you.
lightxbulb
CC BY-SA 4.0
2020-02-23T20:14:34.460Z
@lightxbulb Thanks for trying to help. Ignore the reconstruction. If you agree that $\Phi_j=\int\sigma_M({\rm d}x_1)(1_{S_j}\circ\psi)(x)g(x_0\leftrightarrow x_1)W_e(x_1\to x_0)L_o(x_1\to x_0)=\int\lambda({\rm d}x)(1_{S_j}\circ\psi)(x)f(x)$, then the correct estimator for this is $\frac1n\sum_{i=1}^n\frac{(1_{S_j}\circ\psi)f}q(X_i)$. I'm sure you're confusing this with sampling per pixel as it is done by the `SamplerIntegrator`: http://www.pbr-book.org/3ed-2018/Introduction/pbrt_System_Overview.html#fragment-SamplerIntegratorMethodDefinitions-0.
0xbadf00d
CC BY-SA 4.0
2020-02-23T20:17:06.830Z
You can't ignore the reconstruction when you've literally plugged it into your formula. What you just wrote is a box filter that was multiplied with $\frac{1}{|J|}$ for unknown reasons. Anyways, drop me a line on the CG discord tomorrow and I'll try to clarify it.
lightxbulb
CC BY-SA 4.0
2020-02-23T20:18:57.193Z
@lightxbulb What do you mean? Which thing is a box filter multiplied by $|J|^{-1}$? There is no reconstruction in my latter formula for $\Phi_j$. It should be the same quantity as you denoted by $I_{x,y}$ (for $j=(x,y)$). Are you really sure you don't confuse this with generating path samples per pixel individually (as done by the `SamplerIntegrator`)?
0xbadf00d
CC BY-SA 4.0
2020-02-24T06:13:56.600Z
@lightxbulb Might the missing $|J|$ come from the Jacobian of the transformation $\psi$ from path to raster space?
0xbadf00d
CC BY-SA 4.0
2020-02-24T15:14:38.050Z
Your multiplication with the characteristic function of $S_j$ and the subsequent averaging over $n$ is exactly that - box filtering. If $n$ is not $spp$, but rather $spp \cdot |J|$, then it's also mixed with brightness decrease based on the number of pixels in the image. Additionally you keep mentioning MLT, but the code you linked is not for MLT. MLT has a brightness correction factor $b$ anyways, since it works over the whole sample space, and relies on the Metropolis algorithm. I highly suggest separating this from the implementation and figuring what is the question that you have.
lightxbulb
CC BY-SA 4.0
2020-02-24T15:19:58.687Z
On the other hand, if you're not considering MLT, and instead want to apply your formula for vanilla PT (by considering the sampling sets for ALL pixels), you get something akin to rejection sampling though the characteristic function of $S_j$ that you introduced. Then as you are aware - you need to account for this probability (which happens to be exactly $\frac{1}{|J|}$). Can't say I understand why you would do that (considering samples that have obviously 0 contribution to the pixel) - but you are free to do so, sure - you just have to account for the probability introduced by this.
lightxbulb
CC BY-SA 4.0
2020-02-24T16:04:50.917Z
@lightxbulb The probability introduced by this is already accounted for by the division by $q(X_i)$ in $\frac1n\sum_{i=1}^n\frac{(1_{S_j}\circ\psi)f}q(X_i)$. I'll try to clarify things in a separate thread, once I've figured out what's going on. Trust me, the basic problem is the transformation $\psi$. Let's simplify things: Assume a pinhole camera model. In camera space, there is a natural projection $\pi$ of the hemisphere $H^2:=\{\omega\in S^2:\omega_3>0\}$ onto $\mathbb R^2\times\{1\}$. The directions $\omega$ in $H^2$ correspond to the "forward-facing" directions and
0xbadf00d
CC BY-SA 4.0
2020-02-24T16:05:04.637Z
@lightxbulb $\pi(\omega)$ is the mapped raster position (it lies on the infinite extension of the image plane). Now, our filter functions $h_j$ are defined on the raster space $S\subseteq\mathbb R^2$. If $\Phi:S\to[0,\infty)^3$ is the image function, then the value of the $j$th pixel is defined to be $\Phi_j:=\int_Sh_j\Phi\:{\rm d}\lambda^2$ ($\lambda^2$ denoting the Lebesgue measure on $\mathbb R^2$).
0xbadf00d
CC BY-SA 4.0
2020-02-24T16:05:13.493Z
@lightxbulb If we want to express this quantity in terms of the measurement contribution function $f$ and as an integral over the product area measure $\mu$, we need to transform from path space to $S$ by a transformation $\psi$ and account for that transformation in the integral. My guess is that $\Phi_j=\int_Sh_j\Phi\:{\rm d}\lambda^2=\frac1{\lambda^2(S)}\int(h_j\circ\psi)f\:{\rm d}\mu$, noting that $\lambda^2(S)=|J|$. So, the remaining question is: How precisely is $\psi$ defined.
0xbadf00d
CC BY-SA 4.0
2020-02-24T16:37:20.203Z
I can not know how $\psi$ is defined when this is something you came up with. But I think I finally got what you're trying to do (or at least I hope I did). You want to map $x_0, x_1$ to $[0,1]^2$ I guess? But in your notation they are in $[x_0, x_N] \times [y_0, y_M]$. Then you apply the change of variables where your Jacobian pops out from. Assuming that each pixel is of size $1$, then $(x_N-x_0)(y_M-y_0) = N\cdot M = |J|$. Note that this is equivalent to what I said about rejection sampling, but without requiring the change of variables. I am unsure what the point of this is though.
lightxbulb
CC BY-SA 4.0
2020-02-24T16:41:26.640Z
@lightxbulb No, I want to map $x_0,x_1$ to the raster space $S=[0,a)\times[0,b)$, where $a,b\in\mathbb N$ correspond to the horizontal and vertical image resoulution (as described here http://www.pbr-book.org/3ed-2018/Camera_Models/Projective_Camera_Models.html#). Assuming a pinhole camera, $x_0$ is fixed though.
0xbadf00d
CC BY-SA 4.0
2020-02-24T16:43:04.287Z
Nevermind, your $x_0, x_1$ are something different from what I have here. It's probably just best if you rewrite what exactly it is that you want, because it is currently unclear. I suggest rewriting your integral/estimator so that all $X_i$ are in $[0,1]^D$.
lightxbulb
CC BY-SA 4.0
2020-02-24T16:44:55.607Z
@lightxbulb $x_0$ is the camera position and $x_1$ is a surface point. As long as the direction from $x_0$ to $x_1$ is inside the viewing frustum, there is a map mapping $x_1$ to its raster space coordinates.
0xbadf00d
CC BY-SA 4.0
2020-02-24T19:25:21.623Z
Something you might find interesting btw: https://agraphicsguy.wordpress.com/2016/02/04/the-missing-primary-ray-pdf-in-path-tracing/
lightxbulb
CC BY-SA 4.0
2020-03-05T07:56:51.510Z
@lightxbulb Okay, I think I'm close to figuring out what the source of my scaling confusion is. I've completely rewritten the question: https://computergraphics.stackexchange.com/q/9589/9254. Maybe you can take a look.
0xbadf00d
<raytracing><physically-based><global-illumination><monte-carlo><filtering>
Relation between camera sensor responsivity and image reconstruction filter associated with a pixel
0xbadf00d
# Relation between camera sensor responsivity and image reconstruction filter associated with a pixel

(Asked by: username_0 on 2020-02-20)

Let $W\_{\text e}$ denote the camera sensor responsivity (i.e. importance) and $h\_j$ and $W\_{\text e}^{(j)}$ denote the image reconstruction filter and sensor response associated with the $j$th pixel, respectively, i.e. $$W\_{\text e}^{(j)}=h\_jW\_{\text e}.\tag1$$ The measurement of the $j$th pixel is given by \begin{equation}\begin{split}\Phi\_j&:=\sigma\_M({\rm d}x\_0)\int\sigma\_M({\rm d}x\_1)g(x\_0\leftrightarrow x\_1)W\_{\text e}^{(j)}(x\_1\to x\_0)L\_{\text o}(x\_1\to x\_0)\\&=\int\lambda({\rm d}x)h\_j(x\_1\to x\_0)f(x),\end{split}\tag2\end{equation} where $g$ denotes the geometry term, $\lambda$ is the infinite-product area measure on the path space and $f$ is the measurement contribution function, i.e. $$f(x)=g(x\_0\leftrightarrow x\_1)W\_{\text e}(x\_1\to x\_0)\prod\_{i=2}^kg(x\_{i-1}\leftrightarrow x\_i)f\_{\text s}(x\_i\to x\_{i-1}\to x\_{i-2})L\_{\text e}(x\_k\to x\_{k-1}),\tag3$$ if $x$ is a path of length $k\in\mathbb N$.

> Usually, in the context of Metropolis Light Transport, the scalar contribution function $p$ is defined to be the luminance of $f$. So, by definition, $p(x)=0$ whenever $(x\_0,x\_1)$ doesn't belong to $$B:=\left\{(x\_0,x\_1):W\_{\text e}(x\_1\to x\_0)\ne 0\right\}.$$ Now let $$B\_j:=\left\{(x\_0,x\_1):W\_{\text e}^{(j)}(x\_1\to x\_0)\ne 0\right\}.$$ We clearly have $$B=\bigcup\_jB\_j.\tag4$$ However, in light of the way samples are splat to image (see <http://www.pbr-book.org/3ed-2018/Light_Transport_III_Bidirectional_Methods/Metropolis_Light_Transport.html#fragment-Splatbothcurrentandproposedsamplestomonofilm-0> or <https://pdfs.semanticscholar.org/a386/55316ee7f438ba115b38e8d0b1410c691a26.pdf> and [my other question](https://computergraphics.stackexchange.com/q/9589/9254)), it seems like it is assumed that the union in $(4)$ is disjoint. Is this really guaranteed to be the case?

Taking a look in Section 2.6.4, on page 51, of [this PhD thesis](https://backend.orbit.dtu.dk/ws/portalfiles/portal/51112612/phd240_awk_vers2.pdf), this seems at least to be the case when the camera is a pinhole camera. The authors write: $W\_{\text e}^{(j)}$ *is non-zero only for the set of directions within the pyramid formed by the pinhole and the four corners of the corresponding pixel* $j$.

This seems to indicate that (at least for the pinhole camera) the "support" of $W\_{\text e}^{(j)}$ is limited to the image rectangle region occupied by the $j$th pixel. But is this really true for *any* filter $h\_j$ or only for a special one?

**Tags:** raytracing, physically-based, global-illumination, monte-carlo, filtering

> **Comments:**
> **username_1**: The $B_j$ are not necessarily disjoint - consider having an aperture with some area. Furthermore, a reconstruction filter $h$ in theory (I am considering what PBRT does here), can have an effect over the whole screen if desired. Consider a Laplace or Poisson reconstruction problem from sparse samples (the Dirichlet data). The latter is even the case if you consider gradient domain path tracing.  
> **username_0**: @username_1 You're right. Meanwhile, I've figured the source of my understandings out: https://computergraphics.stackexchange.com/a/9598/9254.  
> **username_1**: On second thought - the $B_j$ are disjoint. I was thinking in terms of directions: (x_1 - x_0), in which case they would not have been. But as a tuple, they are, since one of the points (depending on your convention) is on the film, and those sets (of the points corresponding to a specific pixel) are disjoint. Basically pixels do not overlap. Now I am curious what the result will be of having overlapping pixels.  
> **username_0**: @username_1 I need to think about your comment, but let me ask you a simple question. Say each $h_j$ is a box filter with radius $\frac12$ centered at the pixel center. Won't this yield that $h_j$ (considered as a function on the raster space) is simply the indicator function $1_{S_j}$ of the region $S_j$ of the image plane (in raster space) occupied by the $j$th pixel?.  
> **username_1**: That is correct if you assume your pixels are of size 1. In the PBRT implementation this disjointedness probably also holds - you'll have to look into the internals to see if this is the case. They have a convolution of the normalized filter $\frac{f}{\sum_i f_i}$ with $wL$. If they make sure that the filter cutoff radii cover only a single pixel then the sets are disjoint. One can however pick a global filter - then the above simply does not hold. Ultimately, one wants to compute the integral within each pixel - so you do not even need to consider a filter in the continuous case.  
> **username_1**: In the discrete case the filter becomes useful as a means to make up for the lack of infinitely many estimates over the image. I believe the two things should be considered separately - think of the filtering as a numerical integration rule and you won't be that far off the mark.  
> **username_0**: @username_1 In raster space, the "pixels" are always of size $1$. My point is that I think that PBRT is using the Box filter with radius $\frac12$ implicitly in their implementation of (M)MLT without stating this at any point. Another point which is worrying me in the path tracing implemention: If we use a primary state space sampler which doesn't produce uniformly distributed samples (like the Halton sampler, which is the default one used by PBRT), doesn't this introduce bias in the final estimate?  
> **username_1**: Depends on your convention/implementation - pixels can have any size you would like, as long as you account for this properly (I am not sure whether pbrt uses size 1). You don't have to think too hard - just check PBRT's film functions. As for Halton - it's deterministic - then taking the expectation of the error obviously does nothing and you get that it's biased, sure. However, low-discrepancy samplers have a different property that guarantees bounds on the error and the equivalent of consistency also. Look into Koksma-Hlawka's inequality.  
> **username_0**: @username_1 I'm sorry to bother you again (it's the last time on this, I promise), but the whole scaling thing is continuing to confuse me. Can you tell me what they mean by the "image radiance function" $L$ here http://www.pbr-book.org/3ed-2018/Monte_Carlo_Integration/Bias.html#?  
> **username_0**: @username_1 I really got problems to understand how precisely the image function $\Phi$ on the (continuous) raster space $S$ is defined, since the radiance arriving at an point on $S$ (which is infinitely small) should be $0$. There is clearly a map $\psi$ from the path space $E$ to $S$ (such that $\psi(x)$ only depends on $x_0,x_1$) and we should clearly have $\int_B\Phi=\int_{\psi^{-1}(B)}f\:{\rm d}\mu$, where $f$ is the measurement contribution function and the first integral is integrated wrt the Lebesgue measure on $S$.  
> **username_1**: That's the radiance function assuming a pinhole camera, then $L(x,y) = L((c_x,c_y) \leftarrow (x,y))$, where $c$ is the point of the aperture. Then there is a unique correspondence between a ray direction and a point on the film. At the end of the day, you have discrete pixels, so you want to integrate the radiance coming from all directions at a pixel point, and integrate over the area of the pixel: $$I_{i,j} = \int_{x_i}^{x_{i+1}}\int_{y_j}^{y_{j+1}}\int_{\Omega}W(x,y,\omega)L(x,y,\omega)\,d\omega\,dx\,dy$$ As you will notice there is no reconstruction filter - only sensor function $W$.  
> **username_0**: @username_1 But the reconstruction filter is "built into" your $W$, isn't it? So $W$ depends on $i,j$ or am I missing something?  
> **username_1**: It is not built into it, and has no reason to be. I suggest you check out the CG discord so we could actually chat in real time.  
> **username_0**: @username_1 Everything I'm asking is related to my scaling problem in my other question nhttps://computergraphics.stackexchange.com/q/9589/9254. There's obviously a crucial thing I'm missing. My understanding (which seems to be wrong) is that the measurement of the $j$th pixel is $\Phi_j=\int(h_j\circ\psi)f\:{\rm d}\lambda$, $h_j:S\to[0,\infty)$ being the filter for pixel $j$, $\psi$ the map from path space $E$ to raster space $S$, $f:E\to[0,\infty)^3$ being the measurement contribution function, $\lambda$ the infinite-area measure.  
> **username_0**: @username_1 Let's only consider paths of length $k$ (so $E$ and $\lambda$ replaced accordingly) and say we generate samples $X_n\sim q\lambda$ ($q$ being any bidirectional path tracing strategy; for simplicity, assume path tracing if you like). Then the estimate for $\Phi_j$ should simply be $\frac1n\sum_{i=1}^n\frac{(h_j\circ\psi)f}q(X_i)$. Now PBRT is implicitly assuming a box filter with radius $r=1/2$. So, $h_j$ is simply the indicator function of the $j$th pixel region $S_j$ in raster space. That's fine.  
> **username_0**: @username_1 Now PBRT is summing $\sum_{j\in J}\sum_{i=1}^n1_{S_j}(\psi\circ X_i)\frac fq(X_i)$ (so each summand contributes to precisely one pixel, since the $S_j$ are disjoint). Now we should be left with multiplying each pixel value by $\frac1n$, but PBRT multiplies by $\frac{|J|}n$ ($|J|$ being the number of pixels which is equal to the area of the raster space). It cost me a while, but since PBRT takes $n$ to be of the form $n=\text{samplesPerPixel}\cdot|J|$, $\frac{|J|}n$ happens to be equal to $\text{samplesPerPixel}$.  
> **username_0**: @username_1 So, where does the need for the multiplication by $|J|$ came from? Sure, $h_j$ needs to be normalized, but since the measure of $S_j$ is $1$, there is nothing to normalize here ... (Sorry for extreme long comment. I'll delete them later).  
> **username_1**: This seems to be purely PBRT specific, and as I mentioned I haven't looked into the details of its code in a while. However you mentioned that $n = spp \cdot |J|$, and that the final result is multiplied by $\frac{|J|}{n} = \frac{1}{spp}$. Then everything should be fine, right? What is the issue? Did I maybe misunderstand some part? I highly recommend joining the computer graphics discord so it will be easier to clarify things, rather than producing many comments that are somewhat unrelated to the question.  
> **username_0**: @username_1 I will take a look at the discord thing. No, it's not PBRT specific. When I divide by $n$ then the produced image is black (since obviously $n$ is the wrong factor to divide by). But when I multiply by $|J|/n$ (as PBRT does), the output is correct. I simply don't get why we need to multiply by $|J|/n$ (and not by $1/n$). So, I guess my formula for $\Phi_j$ or $h_j$ is wrong.  
> **username_1**: But you just said that $n = spp \cdot |J|$, and multiplying by $\frac{|J|}{n}$ is the same as multiplying by $\frac{1}{spp}$. Which is exactly what you want.  
> **username_0**: @username_1 I said that I figured out (from testing) that multiplying by $|J|/n$ is the correct thing to do, but I have no theoretical justification where the factor $|J|$ comes from. The theoretical Monte Carlo estimate of $\Phi_j=\int(h_j\circ\psi)f\:{\rm d}\lambda$ is $\frac1n\sum_{i=1}^n\frac{(h_j\circ\psi)f}q(X_i)$. So, *it seems like* either my formula for $\Phi_j$ is wrong or I need to "normalize" $h_j$ by multiplying it with $\frac1{|J|}$ (but since $h_j$ should be equal to $1_{S_j}$ and the region $S_j$ of the $j$th pixel has measure $1$, that doesn't make sense to me).  
> **username_1**: But your formula has $n$ be $spp$, and in the implementation $n$ was $spp \cdot |J|$, no?  
> **username_0**: @username_1 My formula works for every $n$. So, in the implementation I *theoretically* need to divide by $n=spp\cdot|J|$. Either my $h_j$ or my $\Phi_j$ is wrong.  
> **username_1**: The $n$ in your formula needs to be the number of samples, it cannot be arbitrary. That is, it is exactly $spp$. And since in pbrt$ n= spp \cdot |J|$ it obviously requires a multiplication by $|J|$ in the numerator.  
> **username_0**: @username_1 This is the essential part of the code I'm working with: https://coliru.stacked-crooked.com/view?id=7f55a67c5e1dab10. The code uses the suggested number of samples $n=spp\cdot|J|$ and divides by $spp$ at the end of the `render` function. It works, but I don't understand why. From my understanding it should be `m_camera->film->WriteImage(pbrt::Float{ 1 } / m_n);` instead.  
> **username_0**: @username_1 Sure, $n$ is the number of samples. What I meant was that my $n$ can be any number of samples; it doesn't need to be of the specific form $n=spp\cdot|J|$.  
> **username_1**: It seems we're misunderstanding each other. PBRT can also use any number of samples, it's just that in your formula your $n$ variable is not their $n$ variable, it's their $spp$ variable. Why do you believe you should divide by the film area?  
> **username_0**: @username_1 No, that's not correct. When we talk about the number of generated samples, we mean the number of generated paths $X_i$ and they generate $n = spp\cdot|J|$ many; not only $spp$. Please take a look at the code.  
> **username_1**: But this is over the whole screen, since for every pixel you get $spp$ samples. So the weight per pixel is $\frac{1}{spp}$.  
> **username_0**: @username_1 What you are doing now is arguing with intuition. And I agree to you; that makes sense. However, the reason why every pixel is getting $spp$ samples is the multiplication by the filter $(h_j\circ\psi)(X_i)$. But, once again, the correct Monte Carlo estimate is $\frac1n\sum_{i=1}^n\frac{(h_j\circ\psi)f}q(X_i)$. From my understanding, $h_j=1_{S_j}$, but it seems like it should be $|J|1_{S_j}$ instead.  
> **username_1**: I am not arguing from intuition. I am arguing based on the math that I wrote several replies back. You reconstruct each pixel separately - you do not have an estimator over the whole screen - that would make no sense unless you have a single pixel. So if you want to argue that $|J|$ has to be in there you need to motivate that. This is precisely the reason why I told you you should separate the filter from the estimate. In the original problem there is no filter - a filter is required only after you discretise the problem to have finitely many samples.  
> **username_0**: @username_1 Okay, let me continue to denote the region of the $j$th pixel in raster space by $S_j$. Assume a pinhole camera with eye $x_0$. You agree to me that the measurement $\Phi_j$ of the $j$th pixel is $\Phi_j=\int_{D_j}\sigma^\perp(x_0,{\rm d}\omega)W_e(x_0,\omega)L_i(x_0,\omega)$, where $\sigma_M^\perp$ is the projected solid angle kernel on the scene surfaces $M$ and $D_j$ the set of directions corresponding to the $j$th pixel, right? (My $\Phi_j$ is your $I_{x,y}$ for $j=(x,y)$ and my $S_j$ is your $[x_i,x_{i+1}]\times[y_i,y_{i+1}]$.)  
> **username_0**: @username_1 Using the map $\psi:E\to S$, this can be written as $\Phi_j=\int\sigma_M({\rm d}x_1)(1_{S_j}\circ\psi)(x)g(x_0\leftrightarrow x_1)W_e(x_1\to x_0)L_o(x_1\to x_0)=\int\lambda({\rm d}x)(1_{S_j}\circ\psi)(x)f(x)$. The estimator for this is $\frac1n\sum_{i=1}^n\frac{(1_{S_j}\circ\psi)f}q(X_i)$ which is the same as before, since we assumed the filter $h_j$ is equal to $1_{S_j}$.  
> **username_1**: Here $n$ is $spp$ in PBRT, that is the only thing you keep believing is not the case for some reason. I'll make it simple for you: here's the problem you want to solve: $I_{i,j} = \int_{x_i}^{x_{i+1}}\int_{y_j}^{y_{j+1}}\int_{\Omega}W(x,y,\omega)L(x,y,\omega)\,d\omega\,dx\,dy$ An estimator for it is: $$\frac{1}{spp}\sum_{k}^{spp}W(X_k, Y_k, Z_k) L(X_k, Y_k, Z^k)$$  
> **username_0**: @username_1 Let me stress that $U_i:=\psi(X_i)$ are the projections of the sample paths $X_i$ to the raster space and the sample sequence may vary over the whole screen. What you write corresponds to what one is usually doing ordinary path tracing (more generally, in what's been done by the `pbrt::SampleIntegrator`). In path tracing we generate a sample sequence $(X_n)_{n\in\mathbb N}$ for each pixel separately. It seems like you're confusing these approaches.  
> **username_0**: @username_1 Your estimator is only correct if you fix $j$ and use it for a sequence $(X^{(j)}_n)_{n\in\mathbb N}$ of paths such that $\psi(X^{(j)}_n)\in S_j$ for all $n$. This is what `pbrt::SampleIntegrator` is doing. But in MLT only *one* path sequence $(X_n)_{n\in\mathbb N}$ is generated and $\psi(X_n)$ varies over $S=\biguplus_jS_j$.  
> **username_1**: Your reconstruction is still performed per pixel. I believe your whole confusion arises from the fact that you have merged the reconstruction process with the actual problem you are trying to solve - you're basically modelling the problem wrong, and confusing yourself further with PBRT's implementation. The reconstruction process is a pure 2D problem. It's signal reconstruction from sparse data. In pbrt they have decided to perform the reconstruction per pixel. Example - signal reconstruction through radial basis functions, interpolation, approximation, PDEs etc. Anyways I have to go, try to  
> **username_1**: figure out what I mean. If you still haven't figured it out I can try to make some time tomorrow to chat with you.  
> **username_0**: @username_1 Thanks for trying to help. Ignore the reconstruction. If you agree that $\Phi_j=\int\sigma_M({\rm d}x_1)(1_{S_j}\circ\psi)(x)g(x_0\leftrightarrow x_1)W_e(x_1\to x_0)L_o(x_1\to x_0)=\int\lambda({\rm d}x)(1_{S_j}\circ\psi)(x)f(x)$, then the correct estimator for this is $\frac1n\sum_{i=1}^n\frac{(1_{S_j}\circ\psi)f}q(X_i)$. I'm sure you're confusing this with sampling per pixel as it is done by the `SamplerIntegrator`: http://www.pbr-book.org/3ed-2018/Introduction/pbrt_System_Overview.html#fragment-SamplerIntegratorMethodDefinitions-0.  
> **username_1**: You can't ignore the reconstruction when you've literally plugged it into your formula. What you just wrote is a box filter that was multiplied with $\frac{1}{|J|}$ for unknown reasons. Anyways, drop me a line on the CG discord tomorrow and I'll try to clarify it.  
> **username_0**: @username_1 What do you mean? Which thing is a box filter multiplied by $|J|^{-1}$? There is no reconstruction in my latter formula for $\Phi_j$. It should be the same quantity as you denoted by $I_{x,y}$ (for $j=(x,y)$). Are you really sure you don't confuse this with generating path samples per pixel individually (as done by the `SamplerIntegrator`)?  
> **username_0**: @username_1 Might the missing $|J|$ come from the Jacobian of the transformation $\psi$ from path to raster space?  
> **username_1**: Your multiplication with the characteristic function of $S_j$ and the subsequent averaging over $n$ is exactly that - box filtering. If $n$ is not $spp$, but rather $spp \cdot |J|$, then it's also mixed with brightness decrease based on the number of pixels in the image. Additionally you keep mentioning MLT, but the code you linked is not for MLT. MLT has a brightness correction factor $b$ anyways, since it works over the whole sample space, and relies on the Metropolis algorithm. I highly suggest separating this from the implementation and figuring what is the question that you have.  
> **username_1**: On the other hand, if you're not considering MLT, and instead want to apply your formula for vanilla PT (by considering the sampling sets for ALL pixels), you get something akin to rejection sampling though the characteristic function of $S_j$ that you introduced. Then as you are aware - you need to account for this probability (which happens to be exactly $\frac{1}{|J|}$). Can't say I understand why you would do that (considering samples that have obviously 0 contribution to the pixel) - but you are free to do so, sure - you just have to account for the probability introduced by this.  
> **username_0**: @username_1 The probability introduced by this is already accounted for by the division by $q(X_i)$ in $\frac1n\sum_{i=1}^n\frac{(1_{S_j}\circ\psi)f}q(X_i)$. I'll try to clarify things in a separate thread, once I've figured out what's going on. Trust me, the basic problem is the transformation $\psi$. Let's simplify things: Assume a pinhole camera model. In camera space, there is a natural projection $\pi$ of the hemisphere $H^2:=\{\omega\in S^2:\omega_3>0\}$ onto $\mathbb R^2\times\{1\}$. The directions $\omega$ in $H^2$ correspond to the "forward-facing" directions and  
> **username_0**: @username_1 $\pi(\omega)$ is the mapped raster position (it lies on the infinite extension of the image plane). Now, our filter functions $h_j$ are defined on the raster space $S\subseteq\mathbb R^2$. If $\Phi:S\to[0,\infty)^3$ is the image function, then the value of the $j$th pixel is defined to be $\Phi_j:=\int_Sh_j\Phi\:{\rm d}\lambda^2$ ($\lambda^2$ denoting the Lebesgue measure on $\mathbb R^2$).  
> **username_0**: @username_1 If we want to express this quantity in terms of the measurement contribution function $f$ and as an integral over the product area measure $\mu$, we need to transform from path space to $S$ by a transformation $\psi$ and account for that transformation in the integral. My guess is that $\Phi_j=\int_Sh_j\Phi\:{\rm d}\lambda^2=\frac1{\lambda^2(S)}\int(h_j\circ\psi)f\:{\rm d}\mu$, noting that $\lambda^2(S)=|J|$. So, the remaining question is: How precisely is $\psi$ defined.  
> **username_1**: I can not know how $\psi$ is defined when this is something you came up with. But I think I finally got what you're trying to do (or at least I hope I did). You want to map $x_0, x_1$ to $[0,1]^2$ I guess? But in your notation they are in $[x_0, x_N] \times [y_0, y_M]$. Then you apply the change of variables where your Jacobian pops out from. Assuming that each pixel is of size $1$, then $(x_N-x_0)(y_M-y_0) = N\cdot M = |J|$. Note that this is equivalent to what I said about rejection sampling, but without requiring the change of variables. I am unsure what the point of this is though.  
> **username_0**: @username_1 No, I want to map $x_0,x_1$ to the raster space $S=[0,a)\times[0,b)$, where $a,b\in\mathbb N$ correspond to the horizontal and vertical image resoulution (as described here http://www.pbr-book.org/3ed-2018/Camera_Models/Projective_Camera_Models.html#). Assuming a pinhole camera, $x_0$ is fixed though.  
> **username_1**: Nevermind, your $x_0, x_1$ are something different from what I have here. It's probably just best if you rewrite what exactly it is that you want, because it is currently unclear. I suggest rewriting your integral/estimator so that all $X_i$ are in $[0,1]^D$.  
> **username_0**: @username_1 $x_0$ is the camera position and $x_1$ is a surface point. As long as the direction from $x_0$ to $x_1$ is inside the viewing frustum, there is a map mapping $x_1$ to its raster space coordinates.  
> **username_1**: Something you might find interesting btw: https://agraphicsguy.wordpress.com/2016/02/04/the-missing-primary-ray-pdf-in-path-tracing/  
> **username_0**: @username_1 Okay, I think I'm close to figuring out what the source of my scaling confusion is. I've completely rewritten the question: https://computergraphics.stackexchange.com/q/9589/9254. Maybe you can take a look.  

---

## 0 answers

