Botteldooren, Dick

de Greve, Bram

De Muer, Tom

WILLEMS, T

English

Ghent University Academic Bibliography

THE EFFECT OF DIFFUSE REFLECTIONS ON NOISE MAPSUSING PHONON MAPPINGBram de Greve∗, Tom Willems, Tom De Muer, and Dick BotteldoorenDepartment of Information Technology, Ghent UniversitySint-Pietersnieuwstraat 41, B-9000 Ghent, Belgiumbram.degreve@intec.ugent.beAbstractConstructing noise maps for large urban areas is a CPU intensive task. Engineering approxi-mations to outdoor sound propagation are commonly introduced to reduce the computationalcomplexity. In particular, diffuse reflections are often neglected because of their great com-putational expense. This paper reports on an effort to approximate the diffusely reflected fieldby using a technique called phonon mapping.In noise mapping, many techniques can be and are adapted from the field of computer graph-ics where various techniques to approximate indirect lighting exist. A popular technique incomputer graphics is photon mapping. It was introduced in 1995 by Jensen and is known to bean efficient algorithm for simulating indirect lighting. In this paper, it will be shown how thisphoton mapping algorithm can be adapted to the field of noise mapping to simulate diffusereflections.The phonon mapping algorithm consists of two passes. In the first pass, a phonon map isbuilt by bouncing phonons, abstract sound particles, through the environment. On each in-teraction with a (partially) diffusely reflecting surface, it is stored in the phonon map. TheRussian roulette technique is used to determine if the phonon is absorbed, specularly or dif-fusely reflected. Upon reflection, the phonon is traced through the environment again. Thisprocess is continued until the phonon map is sufficiently populated. The second pass is a moretraditional backward ray tracer. Rays are emitted from the receivers. Specular reflections arehandled as usual but for diffuse reflections, no secondary, reflected rays are cast. Instead, thephonon map is queried for nearby phonons. Their presence is used to estimate the energy thatis available from this diffuse reflector.DEFINITIONSDirectional sound power D: sound power per unit solid angle dω. Its units are Wsr.B. de Greve, T. Willems, T. De Muer, D. BotteldoorenAcoustic radiance L: sound power per unit projected area dA⊥ and per unit solid angledω, where dA⊥ = dA cos θ is the projected area of dA on the hypothecial surfaceorthogonal to dω. Its units are Wm2sr.INTRODUCTIONPresent methods for constructing noise maps of large areas are complex and computational in-tensive, even when engineering methods are employed. Classic ray/beam tracers [3, 4, 6] oftenimplement only specular reflections based on image sources. Diffuse reflections are neglectedbecause of their great computational cost, despite the fact that the diffusely reflected field cansignificantly contribute to urban sound fields. This paper presents an effort to approximate thediffusely reflected field using a technique called phonon mapping.Phonon mapping belongs to the family of geometric acoustics, together with methodslike ray tracing and beam tracing, all making the same high frequency approximation. Phononmapping is a particle tracing technique in which acoustic power quanta, called phonons,are bounced through the environment. At diffuse surfaces, incident phonons are stored ina phonon map. Later, the phonon map is queried to estimate the exitant acoustic radiance Lat the diffuse reflector.This technique is adapted from the field of computer graphics where it is known asphoton mapping. It was introduced in 1995 by Jensen [7, 8] as an efficient alternative toradiosity [5] and Monte Carlo ray tracing [9] for simulating global illumination.The goal of noise mapping is to determine how much acoustic power W reaches areceiver located at r. It is given by the integral of the incoming directional power D over alldirections:W (r) =∫ΩD (r, ω) dω (1)In geometric acoustics, this integral is written as the sum over all sound paths that terminateat r:W (r) =∑kDk (r, ωk) (2)This sum is taken over the set of all sound paths from a source to the receiver, via any numberof specular or diffuse reflections. It contains an infinite number of sound paths, therefore forpractical applications a finite approximation to this sum must be made.A traditional ray or beam tracer handles only direct paths, or paths with only specularreflections. The sum of their contributions is Ws. Such a ray tracer disregards all sound pathswith at least one diffuse reflection. The goal of the phonon mapper presented in this paper isprecisely to estimate the sum Wd of the contributions due to these omitted paths. Finally, bothresults must be summed to get W :W (r) = Ws (r) +Wd (r) (3)ICSV13, July 2-6, 2006, Vienna, AustriaFIRST PASS – PHONON TRACINGEmitting Phonons From A Single SourceA number of phonons are emitted from the source with power P until the phonon map ispopulated withNmap phonons.Nemit is the total number of phonons emitted and is not knowna priori since some phonons do not hit a surface at all and thus are not stored in the phononmap, while others may be stored several times.All emitted phonons should have the same power. This is important for two reasons:firstly computational resources are not spent on low power phonons which make little con-tribution to the result, and secondly it reduces the variance of the radiance estimate in thegathering pass. Consequently, phonons must be emitted with directional densities propor-tional to power output. For omnidirectional sources, the phonon emission directions ω mustbe uniformly distributed over the unit sphere. This is accomplished by transforming a uniformvariate (ξ1, ξ2) over [0, 1]× [0, 1) according to:ω (ξ1, ξ2) = 2√ξ1 (1− ξ1) cos 2piξ22√ξ1 (1− ξ1) sin 2piξ21− 2ξ1 (4)The sum of all phonon powers must equal P , so each phonon should have power Pph =PNemit. Since Nemit is initially unknown, this can only be computed after the phonon tracingpass, in a normalisation pass over all stored phonons.Multiple SourcesIn case of multiple sources having powersPi, an additional random variable ξ3 is used to selectan emitting source. Each source is selected with probability pi = PiPtot with Ptot =∑Pi, suchthat a source will emit proportionally many phonons to its (relative) power output. The powerof each phonon will be Pph = PtotNemit .Surface InteractionIncident sound power at a surface point r having normal n is partially reflected back intothe scene. The coefficients Rs and Rd correspond to the fraction of specularly and diffuselyreflected sound power, respectively1.The sound power P carried by phonons hitting a surface must thus be reflected indifferent directions. One way to accomplish this is path splitting: several new phonons areemitted carrying the reflected power: one with the specularly reflected power RsP and nphonons with diffusively reflected power 1nRdP . Although this technique is straightforward,it has two major drawbacks: it requires an exponentially growing number of paths to be traced,and the phonons in the phonon map will no longer have similar powers, which increases thevariance of the radiance estimate in the second pass.1In this paper, the spectral and spatial dependency of both Rs and Rd, and the angular dependency of Rs arenot explicitly mentioned to simplify the notation. The implementation takes this into account, however.B. de Greve, T. Willems, T. De Muer, D. BotteldoorenInstead, another technique called Russian roulette [1, 10] is used to determine whetherthe phonon should be absorbed (discarded) or reflected. A uniform random number over [0, 1)is used to determine which mode of reflection, if any, is chosen:ξ ∈ [0, Rs) → specular reflectionξ ∈ [Rs, Rs +Rd) → diffuse reflectionξ ∈ [Rs +Rd, 1) → absorption(5)In case of diffuse reflection, the new direction ω is sampled using the diffuse BRDF2component ρd. If ρd is Lambertian, this function is the constant Rdpi and the new direction isuniformly sampled from the entire sphere using (4), hemiflipping it to be in the hemisphereabout n if necessary (ω · n < 0).If ρd is sufficiently similar to a Lambertian BRDF, uniform sampling can still be used ifthe power of the outgoing phonon is modified according to the real ρd. However, the phononswill no longer have equal powers, increasing the variance of the radiance estimate in thegathering pass.Multiple FrequenciesIn this paper, the spectral dependency of source powers and reflection coefficients are notexplicitly taken into account. Its effect is however trivial to include using (1/3) octavebands.Consider a source whose spectral power distribution is given by a number of bands havingcenter frequencies fk and powers Pi,k, with Pi =∑k Pi,k. A straightforward approach wouldbe to emit phonons with full spectral information. However, if Rs and Rd are frequencydependent, this complicates surface interaction: Russian roulette requires scalar thresholds in(5), and reflection will change the phonon’s power spectrum.A better approach, used in the implementation, is to emit single-frequency phonons.Analogous to selecting one from many sources, f is selected from the set fk with probabil-ities pk =Pi,kPi. This approach leaves the Russian roulette computation unchanged, and thephonon’s power is also unchanged upon reflection.PHONON MAP CONSTRUCTIONIn the first pass, each time a phonon hits a surface with Rd > 0, it is stored in a list forbuilding the phonon map. It is important to do this independently of the Russian rouletteoutcome, since the phonon map represents an estimation of the incident, not exitant, powerflux. Each phonon’s power Pph, its position rph and its incoming direction ωph are stored.Once sufficiently many (Nmap) phonon interactions have been stored, the first passis terminated and the phonon map must be built. A balanced, adaptive k-d tree [2, 10] isconstructed from the stored phonons. This is a compact and efficient data structure that enableslocating the M nearest phonons around a point r in O (M logNmap) time on average. The2BRDF (Bidirectional Reflectance Distribution Function) ρ (r, ωi, ωo): reflection coefficient as a function ofposition r, incident and exitant angles ωi and ωo.ICSV13, July 2-6, 2006, Vienna, Austriacost of storing the tree is O (Nmap), and constructing it (which is performed only once) isO (Nmap logNmap).SECOND PASS – GATHERING PASSIn the second pass, Wd (r) is evaluated at each receiver. The integral gets approximated as asum over M random direction ωk, chosen from a distribution with probability density p (ω):Wd (r) =∫ΩDd (r, ω) dω ≈1MM∑k=1Dd (r, ωk)p (ωk)(6)If the directions are uniformly distributed over the unit sphere using (4), then p (ω) = 14pi, and(6) simplifies to:Wd (r) ≈4piMM∑k=1Dd (r, ωk) (7)Dd (r, ωk) is computed by casting a ray (r, ωk) into the scene. It intersects a surface at r′having normal n′. Dd (r, ωk) equals the acoustic radiance Lo (r′,−ωk) leaving the surface atr′ towards the receiver. It is evaluated as an integral over the unit hemisphere ℧ about n′ ofthe incoming radiance Li (r′, ω′′) and the BRDF ρ (r′,−ωk, ω′′):Dd (r, ωk) = Lo(r′,−ωk)=∫℧ρ(r′,−ωk, ω′′)Li(r′, ω′′) (ω′′ · n)dω′′ (8)The BRDF is split into specular and diffuse components ρs and ρd. The specular partis non-zero only for the perfectly reflected direction, and is evaluated by tracing the ray(r′, ωk − 2ωk · n′) and evaluating (9) again at the new intersection.Lo(r′,−ωk)= RsLi(r′, ωk − 2ωk · n′)+∫℧ρd(r′,−ωk, ω′′)Li(r′, ω′′) (ω′′ · n′)dω′′(9)The diffuse term is estimated by querying the phonon map. In the neighbourhood of r′, then nearest phonons (Pph, rph, ωph)j are found with ‖rph − r′‖ 6 rmax and n 6 nmax. Thephonons provide information about the incident power flux Pi, which must be converted toincident radiance Li:Li(r′, ω′′)=d2Pi (r′, ω′′)(ω′′ · n′) dAdω′′(10)∫℧ρd(r′,−ωk, ω′′)Li(r′, ω′′) (ω′′ · n′)dω′′ ≈1AM∑j=1ρd(r′,−ωk,−ωph,j)Pph,j (11)where A is the area of the circular neighbourhood in which the n phonons are found. If r isthe largest distance between r′ and the n phonons, then A = pir2. In case of a Lambertiansurface, this simplifies to:∫℧ρd(r′,−ωk, ω′′)Li(r′, ω′′) (ω′′ · n′)dω′′ ≈RdApiM∑j=1Pph,j (12)B. de Greve, T. Willems, T. De Muer, D. Botteldooren-40 -20 0 20 40position [m]-48-46-44-42-40-38-36-34-32SPL [dBA]referencepam 00001=Npam 000005=N-40 -20 0 20 40position [m]-1.5-1-0.500.511.522.5error [dBA]referencepam 00001=Npam 000005=NFigure 1: SPL and error in receivers r (x, 0, z0) with z0 = 10 m; plane z = 0 with Rd = 0.95,Rs = 0; source s (0, 0, z0) with P = 0 dB; M = 225, nmax = 20, rmax = 10 m.Stratified SamplingWhen evaluating (6) with directions ωk generated using (4) with (ξ1, ξ2) directly sampledfrom Λ = [0, 1]× [0, 1), the distribution of directions over the full sphere isn’t as good as onemay think. The samples will tend to clump, causing some areas to have too many samples,while others are undersampled.It is possible to counteract this, thereby also reducing the variance of the result, byusing stratified sampling [10]. The sampling space Λ is split into M1 ×M2 equal strata Λjk(13), and each sample (ξ1, ξ2)jk is taken from a uniform distribution over Λjk. This way,samples are much more evenly spread over the sampling space and no directions are severelyundersampled.Λjk =[jM1,j + 1M1)×[kM2,k + 1M2)(13)APPLICATIONS AND RESULTSValidationA source with power P = 0 dB is positioned at s (0, 0, z0) with z0 = 10 m, above a very largeLambertian reflector with Rd = 0.95 and Rs = 0. Using the phonon mapper, (6) is evaluatedfor receivers at positions r (x, 0, z0). As a reference solution, (7) is evaluated using a directsolution for Dd (r, ωk) with M = 80000:Dd (r, ωk) = Lo(r′,−ωk)=RspiP‖s− r′‖2(s− r′) · n′‖s− r′‖(14)Figure 1 shows the sound level and error for two phonon mapping simulations against thereference. It can be seen that the photon mapper slightly overestimates the levels. This is dueto an understimation of A in (11). A is a too aggressive a bound for the phonons.ICSV13, July 2-6, 2006, Vienna, Austria(a) diffuse component Wd of noise map (b) local expected error on WdFigure 2: Diffuse component Wd of noise map and local expected error for one source s (30, 30, 4)with P = 0 dB; receivers r (x, y, 4) with ∆x = ∆y = 1 m; ground & buildings Rd = 0.05, Rs = 0;Nmap = 40000, M = 100, nmax = 20, rmax = 10 m.Variance StatisticsMonte Carlo simulations often sufferFigure 3: percentiles of expected error distributionfor variousM (number of rays cast in gather pass),stratified vs. unstratified sampling of rays. Otherparameters as in figure 2.from noise or variance in the estimate. Forthe phonon mapper as described in this pa-per, this is mainly influenced by the num-ber M of rays cast in the gathering pass.The local expected error at each point of themap is estimated using the standard devia-tion of the results of 20 different simulationruns with different random seeds, since theaverage will converge to the exact solution.The simulation consists of one source, 4 mabove a ground surface, and in the vicinityof some buildings. Figure 2 shows the resultof one run and the local expected error forM = 100. It can clearly be seen that thelargest errors are found where the source occluded by an obstacle, since fewer phonons canbe sampled. Fortunately, the sound level significantly drops in such regions, reducing the im-portance of the error.To get an idea of the overall expected error, the distribution of the local error over themap is examined. The 90th and 95th percentiles measure how well the largest differencesbehave, while the 50th percentile is a measure of the average spatial error. Figure 3 showsthese percentiles. It can clearly be seen that the expected error decreases with M , however ata linearly increasing cost. It can also be seen that stratified sampling decreases the expectedB. de Greve, T. Willems, T. De Muer, D. Botteldoorenerror for the same M , without incurring significant additional computational cost.POSSIBLE IMPROVEMENTSIn this paper, the diffraction of sound paths is not taken into account for simplicity of pre-sentation. A polygonal beam tracer [3, 4] can easily accommodate this effect in Ws for pathswithout diffuse reflection. However, the phonon mapper must be extended to handle pathswith both diffraction and diffuse reflection. Based on Fresnel zones, Russian roulette may beapplied to subject phonons to diffraction.CONCLUSIONSIn this paper, it is shown how the photon mapping technique is adapted to the field of noisemapping for efficiently simulating the effect of diffuse reflections. The new phonon mappingtechnique is orthogonal to traditional ray/beam tracing, which handles the strictly specularpaths.The phonon mapper consists of two passes: a phonon tracing pass, where phonons arepropagated in an environment to form the phonon map, and a gathering pass where, at eachreceiver, the incident diffusive sound field is reconstructed from the phonon map.The variance of the result is mainly influenced by the number and distribution of gath-ering rays, the latter of which can be significantly improved using stratified sampling.REFERENCES[1] J. Arvo, D. Kirk, “Particle Transport and Image Synthesis”, SIGGRAPH ‘90, 63–66 (1990)[2] J. L. Bentley, “Multidimensional Binary Search Trees Used for Associative Searching”. Comm.of the ACM 18 (9), 509–517 (1975)[3] B. de Greve, “3D polygonal beam tracer for acoustic simulations” (in Dutch), Master thesis,Ghent University, Belgium (2003)[4] T. De Muer, “Policy Supporting Tools for Urban Noise Assessment”, PhD thesis, Ghent Univer-sity, Belgium (2005)[5] C. M. Goral, K. E. Torrance, D. P. Greenberg, B. Battaile, “Modeling the interaction of lightbetween diffuse surfaces”, SIGGRAPH ‘84, 213–232 (1984)[6] P. S. Heckbert, P. Hanrahan. “Beam tracing polygonal objects”. SIGGRAPH ‘84, 119–127(1984)[7] H. W. Jensen, N. J. Christensen, “Photon Maps in Bidirectional Monte Carlo Ray Tracing ofComplex Objects”, Computers & Graphics, 19 (2), 215–224 (1995)[8] H. W. Jensen, “Global Illumination using Photon Maps”, Rendering Techniques ‘96, Eds. X.Pueyo and P. Schrder, Springer-Verlag, 21–30 (1996)[9] J. T. Kajiya, “The Rendering Equation”, SIGGRAPH ‘86, 143–150 (1986)[10] M. Pharr, G. Humphreys, Physically Based Rendering: From Theory to Implementation. (Mor-gan Kaufmann, San Fransisco, ISBN 0–1255–3180–X, 2004) 

The effect of diffuse reflections on noise maps using phonon mapping

https://biblio.ugent.be/publication/361073/file/573970.pdf

The effect of diffuse reflections on noise maps using phonon mapping

Abstract

Similar works

Full text

Available Versions

Ghent University Academic Bibliography