To design a high-performance indoor infrared wireless nondirected data communications system, an accurate estimation of the impulse response is required. This problem can be studied by means of the deterministic^{1} algorithm or the statistical approach.^{2} All surfaces in the channel are considered to be Lambertian diffuse reflectors. However, the Lambertian model cannot describe the reflection pattern of many polished surfaces common in the indoor space, e.g., the painted, varnished surfaces and glass, which presents a strong directional-diffuse (DD) reflection besides diffuse or just a simple mirror reflection (MR). Lopez-Hernandez et al.;^{3} considered MR and used a Monte Carlo ray-tracing algorithm to reduce the numerical burden. In their model only one new ray is generated for each reflection. Lomba^{4} and Yang,^{5} respectively, experimentally demonstrated that the empirical Phong’s scattering model can effectively approximate any nondiffuse reflection besides MR and as diffuse. Phong’s model consists of DD and diffuse reflection components. For an accurate and fast simulation of the impulse response in a real channel, we use a split Monte Carlo (SMC) algorithm to generate two split random rays for each reflection: one diffuse and another DD. Compared with the deterministic algorithm with exponential ray growth, the SMC algorithm with nearly geometrical ray growth can greatly reduce the numerical burden.

First, the surface reflection property, the relevant light source data, and the receiver in the channel are defined. The intensity profile of the light source may have a shape that is a generalized Lambertian, Guassian, or uniform, etc. The receiver typically employs silicon p-i-n detectors encapsulated in hemispherical or plano-cylindrical lenses, which provide proper optical gain as well as maintain a wide field of view (FOV). Based on the work of Lomba^{4} and Yang,^{5} the expression of Phong’s model to approximate the radiation intensity distribution can be given by

## (1)

$$I=\rho {P}_{s}/\pi \cdot [{r}_{d}\text{\hspace{0.17em}}\text{cos}({\theta}_{o})+C(1-{r}_{d}){\text{cos}}^{m}(\alpha )]$$_{s}is the incident optical power, r

_{d}is the percentage of the incident signal that is reflected diffusely, θ

_{o}is the angle between the observation direction vector

**v**and the surface normal vector

**n**,

*m*is the exponent of the DD reflection, C=(m+1)/2 is the normalization factor, and α is the angle between vector

**v**and the MR direction vector

**u**of the incident direction. The first and the second terms in Eq. (1) denote the reflection components of the homogeneous diffuse and DD, respectively. Two corrections are given in Refs. 4 and 5. First, the second term is weighted by factor

*C*to preserve energy conservation. Second, α is substituted for θ

_{o}−θ

_{i}(θ

_{i}is the angle between vector

**u**and

**n**) to extend the applicable reflection range from the primary incidence plane to three-dimensional space. Phong’s model can also denote Lambertian when r

_{d}=1 and MR when r

_{d}=0 and m→∞, respectively.

We use the light source with a generalized Lambertian radiation intensity pattern as an example. Many starting rays can be generated by using two uniform random parameters μ and υ (∈[0,1]). For a random starting ray, the components *x*, *y*, and *z* of its unit direction vector **s** relative to a coordinate system with the source normal as *z* axis can be chosen by

## (2)

$$z=\sqrt[{}^{n+1}]{\mu},\text{\hspace{0.17em}}\text{\hspace{0.17em}}x=\sqrt{1-{z}^{2}}\text{\hspace{0.17em}}\text{cos}(2\pi \upsilon ),\text{\hspace{0.17em}}\text{\hspace{0.17em}}y=xtg(2\pi \upsilon )$$*n*is the mode number of the generalized Lambertian radiation lobe, which specifies the directionality of the source. Convert the components of

**s**into the absolute coordinate of the room. Then P

_{los}, the line-of-sight (LOS) power contribution from the source to the receiver, can be calculated.

In the normal Monte Carlo raytracing algorithm, when a ray strikes an opaque surface, the intersection point is treated as a new point source. Thus a new random ray is generated, and the process continues until the ray reaches the receiver. However, the probability of the event that the random ray reaching the receiver before the maximum simulation time is very low, usually less than 10^{−6}. To improve the convergence, we sum P_{r}, the power directly contributed from each intersection point to the receiver. Substituting the corresponding surface reflection parameters approximated by Phong’s model into Eq. (1), we can get the radiation intensity *I* and further achieve P_{r}:

_{r}from the total reflected power P

_{s}weighted by the reflection coefficient ρ to get the value of the power carried by the next new ray. However, a complicated algorithm is required to generate only one new ray from the nonrotation symmetric intensity distribution of Phong’s model. Fortunately, the first and second terms are rotationally symmetric with the directions of the surface normal vector and MR, respectively, which can be considered as two generalized Lambertian patterns with n=1 and n=m, respectively. Therefore we can simplify the ray generation process of Phong’s model by means of the SMC technique. In detail, for each incident ray, two new random rays leaving the surface are generated simultaneously: one diffuse component with the power of r

_{d}(ρP

_{s}−P

_{r}) and another DD component with the power of (ρP

_{s}−P

_{r})(1−r

_{d}). The direction vectors of the two new rays are generated by Eq. (2) with n=1 and n=m, respectively. The direction vector of the diffuse ray is based on a coordinate system normal to the reflected surface, while the

*z*component of that of the DD ray is in the MR direction. Both of them need to be transformed into a unified room coordinate; Fig. 1 summarizes this procedure. The ray is reflected and split until the power it carries is less than the threshold. The power contribution from each intersection point is added to the impulse response curve. The delay spread is equal to the path length divided by the speed of light. For the special case of Phong’s model, diffuse or MR, only one new ray, diffuse or mirror component, will be generated. As a result, this approach can also be used to validate previous simulation schemes.

To show the effects of the non-Lambertian reflection on the impulse function, as an example, we investigate configurations *D* and D^{′}. The former is the same as that used in Ref. 1, in which all the surfaces are diffuse reflectors. And the latter has both diffuse and DD reflectors, in which, except for the glass ceiling with r_{d}=0.0001 and m=280 and the south Formica wall with r_{d}=0.9019 and m=112, the parameters are the same as those in configuration *D*. We use the SMC raytracing algorithm and Phong’s model to approximate the radiation patterns of all the indoor reflection surfaces.

The simulated impulse responses due to a 1-W source for configurations *D* and D^{′}
are shown in Fig. 2, where the time origin is the transmitting start time. The solid curve for configuration *D* obtained by means of the split Monte Carlo method with Phong’s model is similar to that obtained by the deterministic algorithm with the Lambertian model in Ref. 1, but with much less numerical time and memory requirements. Figure 2 also shows that the difference between the dashed curve and the solid curve can no longer be neglected, which indicates the DD reflection has significant effects on the simulation of the impulse response. For instance, the dashed curve exhibits a narrow peak of 1736 s^{−1} near 20 ns, which is the cause of the strong DD reflection of the glass ceiling. Hence, the distribution of the power is more focused than that of configuration *D*. It can be predicted that the channel bandwith will increase. In fact, the −3-dB bandwith of configuration D^{′}
is 55 MHz, the total received power is 0.66 μW, and the root mean square delay spread is 0.3 ns, which are respectively 23 MHz greater, 0.19 dB less, and 1.9 ns less than those of configuration *D*.

In this letter, the split Monte Carlo method based on Phong’s model to calculate the impulse response of a wireless indoor IR communications system is presented. It allows evaluation of not only Lambertian diffuse but also DD reflections. The simulation results show the important impact of the DD reflector. This approach can be used to validate previous simulation schemes with less computer complexity and offers increased flexibility for a complex channel with various reflecting patterns.