Rendering equation for spherical lights I am currently implementing spherical lights in my DirectX game engine. I decided to build the shading formula from the classical rendering equation. Here is what I get:


The simplifications show that the cosine term should be removed for spherical lights since power is a Flux.
Here the problem. While looking for recourses online I found that engines always multiply by the cosine term, even for spherical lights.

Is it wrong or am I missing something?
If you just want to explicitly sample an area light, then here's the general procedure you should follow. Pick light $i$ out of $L$ lights with some probability $p_i$ (the other probabilities being $p_1,...,p_L$, a light may be picked through inverse transform sampling). Pick a point $\pmb{y}$ on the surface of the light with some probability $q_i(\pmb{y})$. Then your estimator is:

$$I_N = \frac{1}{N}\sum_{k=1}^{N}\frac{1}{p_{i_k}q_{i_k}(\pmb{y}_k)}f(\pmb{z} \leftarrow \pmb{x} \leftarrow \pmb{y}_k)L_e(\pmb{y}_k\rightarrow\pmb{x})\frac{\cos\theta_x\cos\theta_{y_k}}{\|\pmb{x}-\pmb{y}_k\|^2_2}V(\pmb{x},\pmb{y}_k)$$

Note that this is only for direct illumination at point $\pmb{x}$. I have generalized it to a secondary estimator using $N$ samples, so essentially you pick $N$ points $\pmb{y}_1,...,\pmb{y}_N$ on possibly different lights. This is an estimator that can directly be derived from the area formulation of the rendering equation:

$$L(\pmb{z} \leftarrow \pmb{x}) = L_e(\pmb{x} \leftarrow \pmb{z}) + \int_{\Omega}f(\pmb{z} \leftarrow \pmb{x} \leftarrow \pmb{y})L(\pmb{y}\rightarrow\pmb{x})\frac{\cos\theta_x\cos\theta_y}{\|\pmb{x}-\pmb{y}_k\|^2_2}V(\pmb{x},\pmb{y})\,dA(\pmb{y})$$

The area formulation can be derived from the solid angle one, by using the identity:

$$d\omega = \frac{\cos\theta_y}{r^2}dA$$,

as well as:

$$L_i(\pmb{x}, \omega) = L(r(\pmb{x},\omega) \rightarrow \pmb{x}) = L(\pmb{y} \rightarrow \pmb{x})V(\pmb{x},\pmb{y})$$