A method for calculating straight edge diffraction of acoustic waves based on path tracking
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2022-12-06
- Publication Date
- 2026-08-11
AI Technical Summary
但目前DLSM在实际声传播模拟器中的应用还没有出现
[0044]与传统方法相比,本文方法计算效率更高,且在处理大规模,遮挡关系复杂的场景时无需进行预求解波动方程,几何简化等预处理计算。本文方法使用的衍射模型无基尔霍夫近似,无限长直边假设等近似操作,因此结果精度更高,对于较简单几何体,其结果与理论衍射解完全符合。此外,基于路径跟踪的衍射算法是基于路径跟踪的传统声传播算法的自然扩展,因此本发明的方法也较易于添加至现有的声传播模拟器之中。
Smart Images

Figure CN117076827B_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the field of acoustic simulation technology, and in particular to a method for calculating the straight-edge diffraction of sound waves based on path tracking. Background Technology
[0002] Sound propagation simulation technology has wide applications in fields such as architectural design and virtual reality. Among sound propagation simulation algorithms, geometric acoustics-based methods are commonly used in practice due to their high computational efficiency. Geometric acoustics methods simplify the wave equation, assuming that sound waves propagate in straight lines. This simplification requires that phenomena such as diffraction and interference be handled independently. The lack of diffraction components results in noticeable discontinuities in the sound propagation simulation results when the listener enters the "shadow" area obscured by objects, which does not match the everyday acoustic experience and greatly reduces the realism of the simulation results.
[0003] Incorporating diffraction phenomena into the framework of geometric acoustics is a complex undertaking. In practical engineering, a diffraction-representing impulse response is typically superimposed on top of the impulse response calculated using ordinary geometric acoustics theory. This impulse response can be calculated using wave equations. However, in many cases, it's desirable to calculate diffraction in a manner similar to ordinary geometric methods. A natural approach is to first break down the scene into simpler, more computationally calculable diffraction structures, then piece together the results to form the final diffraction solution. Objects in the scene are typically represented using triangular meshes, and the shadow boundaries created by object occlusion are related to the edges of these meshes. Therefore, it's possible to attempt to calculate diffraction near each edge of the triangular mesh. This is the motivation behind geometric acoustics' interest in representing straight-edge diffraction.
[0004] Straight-edge diffraction, or wedge diffraction, refers to diffraction near two infinitely large half-planes connected by a straight edge and forming a certain angle. In the mid-20th century, Keller first proposed a method to describe convex straight-edge diffraction in geometrical optics, called the geometric diffraction theory. Keller's results are close to the exact solution in most directions, but the model fails near the "shadow boundary," that is, near the discontinuity of the geometrical optical solution. Later, others improved the GTD, solving its discontinuity problem near the shadow boundary. The improved new theory is called the uniform geometric diffraction theory, which can be found in MCNAMARA DA, PISTORIUS CWI, MALHERBE JA G. Introduction to the Uniform Geometrical Theory of Diffraction [M]. Artech House, 1990. Another attempt to express straight-edge diffraction began with the Biot-Tolstoy solution. It gives the analytical form of the echo when a spherical wave is incident on a convex dihedral angle. When applying the Biot-Tolstoy solution to acoustic research, Medwin, based on the intuition of Huygens' principle of the wave equation, decomposed the echo into a combination of a series of radial functions centered at points on the straight edge, and found that this decomposition helps to explain the diffraction problem of finite-length straight edges. Combining the Biot-Tolstoy solution and the Medwin decomposition, another descriptive model of straight-edge diffraction can be obtained, called the Biot-Tolstoy-Medwin (BTM) model, for details see SVENSSON UP, FRED RI, VANDERK00Y J. An analyticsecondary source model of edge diffraction impulse responses[J]. Journal of the Acoustical Society of America, 1999, 106: 2331-2344.
[0005] UTD (Unified Geometric Diffraction Theory) and BTM (Boundary Theory of Diffraction) are the two most common geometric diffraction models in acoustic propagation simulation. UTD, being an earlier model, has been extensively studied. Research articles can be found on plane and spherical incident waves, ideal half-planes or curved half-planes, and various boundary conditions. However, UTD uses Kirchhoff approximation and assumes infinite straight-side length. Therefore, in practical applications, UTD results will inevitably differ from the exact solution. Furthermore, the UTD solution is not represented in the time domain like a typical impulse response, but in the frequency domain. This makes it difficult to integrate with other parts of the geometric acoustic propagation simulator. Typical applications often only select results at a few frequency points of the UTD solution, discarding all phase information. This results in errors when applying UTD to practical acoustic propagation simulators that are much greater than theoretical errors. The Biot-Tolstoy solution upon which the UTD model relies is an analytical solution to the infinitely long straight-side diffraction problem, and the Medwin decomposition is not an approximation. Therefore, the theoretically achievable accuracy of the BTM method is higher than that of UTD. Furthermore, since the Medwin decomposition decomposes straight edges into a series of radial functions, the BTM method can also handle diffraction problems involving piecewise straight edges and even curved edges. The disadvantages of BTM compared to UTD are the limited amount of relevant research data and its less efficient implementation in existing acoustic propagation simulators.
[0006] Both UTD and BTM have wide applications in sound propagation simulation. At the beginning of this century, papers implemented single-shot diffraction based on UTD in acoustic virtual reality systems, and later researchers quickly extended this to higher-order diffraction cases. For more complex scenes, there are also methods for quickly finding feasible diffraction paths, but these require pre-simplification of the scene and direct use of the original BTM expression for sound propagation simulation. However, their computational efficiency is too low, making them only suitable for offline calculations in small-scale scenes. There are also methods based on BTM to calculate diffraction solutions for common geometries and use these solutions to approximate objects in the real environment, such as volume diffraction and propagation methods.
[0007] Besides UTD and BTM, there are some other diffraction models that have not yet gained widespread popularity. Some methods use Heisenberg's uncertainty principle to explain acoustic wave diffraction, using the diffraction angle probability density function to describe the non-linear propagation of diffracted waves. This method lacks a physical basis and is difficult to guarantee the correctness of the results; the recently proposed directional line source model is a promising model. It uses a few approximations to explain a large number of different cases of straight-edge diffraction (see MENOUNOU P, NIKOLAOUP. Analytical model for predicting edge diffraction in the time domain.[J].The Journal of the Acoustical Society of America, 2017, 1426:3580). However, the application of DLSM in actual acoustic propagation simulators has not yet appeared.
[0008] Therefore, we urgently need a solution for calculating high-precision straight-edge diffraction in large-scale scenarios within a geometric framework. Summary of the Invention
[0009] The purpose of this invention is to address the shortcomings of existing technologies by providing a method for simulating acoustic diffraction phenomena based on path tracking.
[0010] This invention is achieved through the following technical solution:
[0011] A method for calculating straight-edge diffraction of acoustic waves based on path tracking includes the following steps:
[0012] (1) Read the acoustic environment state related to diffraction, the acoustic environment state including the position of the sound source and the listener, the geometric description of the scene, and the interface features of the surface of the geometric body;
[0013] (2) A large number of diffraction paths are randomly generated in the scene. The path is composed of a series of nodes in the scene and line segments connecting the nodes. The starting node of the path is the sound source, the ending node is the listener, and the remaining intermediate nodes are located on the convex straight edge of the geometry in the scene. Except for the nodes that constitute the path, the other parts of the path do not intersect with the surface of the geometry in the scene.
[0014] (3) Calculate the final contribution of each path. Suppose a path contains n+1 nodes, and each node is numbered from x0 to x1 in sequence from the sound source to the listener. n Let x be the line segment between the nodes. i →x i+1 The complete path is x0→x1→…→x n Then the total contribution f(x0→x1→…→x) n The calculation formula is as follows:
[0015]
[0016] Where L0(x0→x)1 represents the sound source characteristics, and M(x i →x i+1 ) is x i →x i+1 The acoustic medium characteristics at ρ(x) i-1 →x i →x i+1 ) is x i The radial function representation of diffraction at the straight edge;
[0017] (4) Accumulate the contributions of all paths into the impulse response function: for the path x0→x1→…→x n It is known that its contribution is f(x0→x1→…→x) in step (3). n Let p(x0→x1→…→x) be the probability of its generation. n The propagation delay is t. S Then, the impulse response function is accumulated on top of the original impulse response function. Where δ(x) is the unit impulse function;
[0018] (5) Output the impulse response function generated in step (4) as the diffraction simulation result.
[0019] Specifically, step (2) includes the following sub-steps:
[0020] (2.1) Initialize the path so that the only node it contains is the starting point of the path;
[0021] (2.2) Let the current path state be x0→x1→…→x n From x n A ray is randomly emitted from the starting point; if x n If it is not the starting point of the path, the exit angle of the ray is determined by a sampling function that matches the radial function of the straight edge diffraction.
[0022] (2.3) Let the first intersection point between the ray and the scene geometry be denoted as . This is called a pseudo-intersection point; if the point does not exist, the intersection fails, and the process proceeds to step (2.8).
[0023] (2.4) Pseudo-intersections in step (2.3) Randomly select one of the three sides of the triangle; let the triangle be Δabc, the selected side be ab, and the vertex opposite the selected side be c. Then connect c with... Extend the line to intersect ab at a point, denoted as x. n+1 This is called the intersection point;
[0024] (2.5) If the side ab in step (2.4) is not a convex straight side, or the line segment x n →x n+1 If the intersection point x intersects with the scene geometry outside of the two endpoints, the intersection calculation fails, and proceed to step (2.8); otherwise, the intersection point x is set. n+1 Add as a node to the path x0→x1→…→x n At the end, a new path is formed: x0→x1→…→x n+1 ;
[0025] (2.6) Calculate the path x0→x1→…→x n+1 The generation probability: If the path x0→x1→…→x is known n The generation probability is p(x0→x1→…→x) n If the path is x0→x1→…→x n+1 The generation probability calculation method is as follows in, For x n The incident path direction is x n-1 →x n At that time, the direction of the emitted ray is The probability, if n=0, then this term must be replaced with the direction of the ray emitted from the ray source. The probability of; For rays Intersection with scene geometry The probability of p C (x n →x n+1 ) represents the probability correction term;
[0026] (2.7) If the path length does not meet the user's requirements, proceed to step (2.2);
[0027] (2.8) Use the current path as the diffraction path in the scene.
[0028] Further, in step (2), paths are generated simultaneously from both the sound source and the listener. The path originating from the sound source is called the forward sub-path, and the path originating from the listener is called the backward sub-path. Then, connections are randomly established between the nodes of the forward and backward sub-paths to form a complete path from the sound source to the listener. Let p be the probability of generating the portion of the forward sub-path from the sound source to the connected point. F The probability of generating the backward subpath from the listener to the connected point is p. B The probability of establishing a connection is p. L Then the probability of generating the complete path is p. F ·p B ·p L .
[0029] Furthermore, in step (5), the path contribution needs to be multiplied by a factor C used to suppress extreme values. S :
[0030]
[0031] in For a user-defined variable, C o C is the C corresponding to each randomly generated segment in the path. o The product of C for each path segment o A parameter used to measure the likelihood that a path containing this path segment will produce extreme values.
[0032] Furthermore, in step (3), calculating the contribution requires dividing by the [F] of the corresponding positions of the forward and reverse path connection nodes. j,i The product of elements in ]; the [F j,i Let [F] be an M x N array, where N is the maximum number of nodes in the forward and backward paths; the nodes of the paths are stored in an array with the same number of rows and columns, where each row corresponds to a single path. j,i The first column contains only 1s, and the first column of the node array represents the path starting point. Given that the nodes in the (i-1)th column have already been generated, the generation process of the node in the i-th column of the array includes the following sub-steps:
[0033] (5.1) Starting from the (i-1)th node x in each row of the node array i-2 Starting from the origin, a ray is randomly emitted; if i > 1, the emission angle of the ray is determined by a sampling function that matches the radial function of the straight edge diffraction.
[0034] (5.2) Let the first intersection point between the ray and the scene geometry be denoted as . This is called a pseudo-intersection point; if the point does not exist, the intersection fails, and the process proceeds to step (5.7).
[0035] (5.3) Pseudo-intersections in step (5.2) Randomly select one of the three sides of the triangle; let the triangle be Δabc, the selected side be ab, and the vertex opposite the selected side be c. Then connect c with... Extend the line to intersect ab at a point, denoted as x. i-1 This is called the intersection point;
[0036] (5.4) If the edge ab in step (5.3) is not a convex straight edge, or the line segment x i-2 →x i-1 If the intersection point x intersects with the scene geometry outside of the two endpoints, the intersection calculation fails, and proceed to step (5.7); otherwise, the intersection point x is set. i-1 The path x0→x1→…→x is added as a node to the current row.i-2 At the end, a new path is formed: x0→x1→…→x i-1 ;
[0037] (5.5) Calculate the path x0→x1→…→x i-1 The generation probability: If the path x0→x1→…→x is known i-2 The generation probability is p(x0→x1→…→x) i-2 If the path is x0→x1→…→x i-1 The generation probability is calculated as follows:
[0038]
[0039] in, For x i-2 The incident path direction is x i-3 →x i-2 At that time, the direction of the emitted ray is If i = 1, then this term must be replaced with the direction of the ray emitted from the ray source. The probability of; For rays Intersection with scene geometry The probability of p C (x i-2 →x i-1 ) represents the probability correction term;
[0040] (5.6) If the path length does not meet the user's requirements, proceed to step (5.1) to prepare to generate the next column of nodes;
[0041] (5.7) For each node x in the (i-1)th column i-2 The set S of rows where successor nodes were successfully generated is calculated. + (x i-2 and the set of failed rows S - (x i-2 );
[0042] (5.8) S + (x i-2 The node in the i-th column of each row of S is copied into S. - (x i-2 In the i-th column of each row of the array, if the node in the j-th row is copied k times, then let F... j,i =(k+1)F j,i-1 .
[0043] The beneficial effects of this invention are:
[0044] Compared to traditional methods, the proposed method offers higher computational efficiency and eliminates the need for preprocessing calculations such as solving wave equations and geometric simplification when handling large-scale scenes with complex occlusion relationships. The diffraction model used in this method avoids approximations such as Kirchhoff's approximation and the assumption of infinitely long straight edges, resulting in higher accuracy. For simpler geometries, the results perfectly match the theoretical diffraction solutions. Furthermore, the path-tracking-based diffraction algorithm is a natural extension of the traditional path-tracking sound propagation algorithm, making it easily integrateable into existing sound propagation simulators. Attached Figure Description
[0045] Figure 1 This is a schematic diagram of the straight-edge coordinate system in this invention;
[0046] Figure 2 This is the absolute value image of the straight-edge diffraction interface term function of this invention;
[0047] Figure 3 This is a schematic diagram of the sampling probability distribution of the straight-edge diffraction interface term in this invention;
[0048] Figure 4 This is a schematic diagram of the setup for a first-order diffraction test scenario;
[0049] Figure 5 This is a schematic diagram of the calculation results of the present invention and the control method under a first-order diffraction test scenario setup;
[0050] Figure 6 This is a schematic diagram of the setup for a second-order diffraction test scenario;
[0051] Figure 7 This is a schematic diagram of the calculation results of the present invention and the control method under a second-order diffraction test scenario. Detailed Implementation
[0052] This invention provides a method for calculating the straight-edge diffraction of acoustic waves based on path tracking, comprising the following steps:
[0053] (1) Read the acoustic environment state related to diffraction, the acoustic environment state including the position of the sound source and the listener, the geometric description of the scene, and the interface features of the surface of the geometric body;
[0054] (2) A large number of diffraction paths are randomly generated in the scene. The path is composed of a series of nodes in the scene and line segments connecting the nodes. The starting node of the path is the sound source, the ending node is the listener, and the remaining intermediate nodes are located on the convex straight edge of the geometry in the scene. Except for the nodes that constitute the path, the other parts of the path do not intersect with the surface of the geometry in the scene.
[0055] (3) Calculate the final contribution of each path. Suppose a path contains n+1 nodes, and each node is numbered from x0 to x1 in sequence from the sound source to the listener.n Let x be the line segment between the nodes. i →x i+1 The complete path is x0→x1→…→x n Then the total contribution f(x0→x1→…→x) n The calculation formula is as follows:
[0056]
[0057] Where L0(x0→x1) represents the sound source characteristics, and M(x i →x i+1 ) is x i →x i+1 The acoustic medium characteristics at ρ(x) i-1 →x i →x i+1 ) is x i The radial function representation of diffraction at the straight edge;
[0058] (4) Accumulate the contributions of all paths into the impulse response function: for the path x0→x1→…→x n It is known that its contribution is f(x0→x1→…→x) in step (3). n Let p(x0→x1→…→x) be the probability of its generation. n The propagation delay is t. S Then, the impulse response function is accumulated on top of the original impulse response function. Where δ(x) is the unit impulse function;
[0059] (5) Output the impulse response function generated in step (4) as the diffraction simulation result.
[0060] The method for calculating straight-edge diffraction of acoustic waves based on path tracking is as follows:
[0061] Path tracing: The path tracing technique in acoustic propagation simulation refers to the path tracing algorithm used in acoustic propagation simulation and the application of diffraction representation in the algorithm.
[0062] In sound propagation simulation, the path tracing algorithm solves the following rendering equation:
[0063]
[0064] In the formula, t represents time. The moment t = 0 is the moment when the sound source begins to emit sound; L(x′→x, t) is the energy response or impulse response excited by the sound source x′ at position x, which is the objective function to be solved in this equation; L0(x′→x, t) represents the direct component, that is, the response excited directly at x without interaction with the scene geometry; V(x″, x) is the visibility function, which is 1 when x″ and x are mutually visible, and 0 otherwise; ρ(x′→x″→x) is the interface term, which represents the intensity of the emitted sound near a point x″ on the interface in the direction x″→x when the incident sound direction is x′→x″; M t (x″, x′, t) is the medium term, representing the influence of the sound medium between x″ and x′, including energy absorption and propagation delay. The asterisk * is the convolution symbol; Ω represents the point set on the surface of the corresponding scene geometry, A Ω This is the area measure over the set. When calculating diffraction, for a homogeneous medium, the following medium term is used:
[0065]
[0066] Where a is the acoustic attenuation coefficient of the medium, ||x″-x′|| is the length of the vector x″-x′, δ(·) is the unit impulse response function, and c is the speed of sound in the medium.
[0067] Let the integral above be operator T, then the rendering equation can be simplified as L = L0 + TL; the equation can be expressed using Neumann series expansion:
[0068]
[0069] The operator T represents a single interaction between the sound wave and the scene geometry. Since T is essentially an integral, its value can be calculated using Monte Carlo integration.
[0070] In path tracing, this invention will use T n L0 represents the expected value of the contribution to the path; a path consists of a series of nodes in the scene and line segments connecting the nodes, with the starting node being the sound source, the ending node being the listener, and the intermediate nodes located on the surface of the scene geometry; a path containing n intermediate nodes (n+2 nodes) is called an n-order path, which is related to T in the series. n L0 corresponds to this item.
[0071] It is easy to see from the expression for the medium term that the total propagation delay of a path is the sum of the time delays of each segment of the path (length divided by the speed of sound). Therefore, it can be calculated separately. Ignoring the time term and only calculating the intensity of the path contribution, we use the medium term, which has no time parameter. Replace the original medium term M in the rendering equation. t At this point, an n-order path x0→x1→…→xn+1 The contribution expression is as follows:
[0072]
[0073] If the above path is randomly generated, its expected contribution needs to be divided by its generation probability. Considering the time delay, the expected contribution can be calculated as follows:
[0074]
[0075] Where t S This represents the total time delay for the path. When calculating the expectation using multiple paths, the expected value is the average of the results from each path.
[0076] There are many ways to generate random paths; either one-way path tracing or two-way path tracing can be used.
[0077] Interface term representation for straight-edge diffraction: The key to calculating diffraction using path tracing lies in constructing the corresponding interface term. This invention will present the form of the interface term here. For example... Figure 1 As shown, the coordinate system near the straight edge and the representation of the incident and outgoing path directions, along with related symbols, are given. In this invention, n is called the normal, b the binormal, and t the tangent (the direction of t coincides with the direction of the straight edge, and can be positive or negative). These three constitute a frame at a point on the straight edge. Half of the dihedral angle at the straight edge is denoted as φ. Since this paper only considers convex straight edges, the value of φ should be less than π / 2. The incident direction vector v... i and the outgoing direction vector v o Using a series of angles to represent:
[0078]
[0079]
[0080] In this coordinate system, write the interface term form for straight-edge diffraction; first, define the following auxiliary variables:
[0081]
[0082]
[0083] ω lu =θ i -π
[0084] ω lr =-θ i -π+2φ
[0085] ω ru =θ i +π
[0086] ω rr =-θ i +π-2φ
[0087]
[0088] If the geometric surface near a point on the straight edge satisfies the Dirichlet boundary conditions, then the interface term at that point is:
[0089]
[0090] If the geometric surface satisfies the Neumann boundary conditions, then the interface term is:
[0091]
[0092] The form of the above interface terms can be proven mathematically. Under other conditions, a form like aρ can be considered. D +bρ N The interface term, where a and b are two constants.
[0093] In path tracing, it's necessary to emit a ray from a point on a straight edge, ensuring that the probability distribution of its emission direction approximates the absolute value of the interface term function. This probability distribution requires specific design. To describe this distribution, some auxiliary variables need to be defined. Let...
[0094]
[0095]
[0096] When using Dirichlet boundary conditions, let D = D D When using Neumann boundary conditions, let D = D N The remaining variables are defined as follows:
[0097]
[0098]
[0099]
[0100]
[0101]
[0102] Then, the probability density function of the ray emission direction. This can be represented as follows:
[0103]
[0104] Where, p θ It is θo The probability density function, yes The probability density functions are as follows:
[0105]
[0106]
[0107] The function constructed above satisfies the requirements of nonnegativity and normalization of the probability density function.
[0108] Figure 2 This demonstrates the effect under Dirichlet conditions. θ i =0.4636, The absolute value of the interface term at that time, Figure 3 This is a schematic diagram of the sampling probability distribution of the straight-edge diffraction interface term in this invention, showing the sampling probability density function in the corresponding direction. It can be seen that the shapes of the functions in these two diagrams are very similar.
[0109] Diffraction Path Generation: In path tracing algorithms, a common path generation method is to generate each node in the path sequentially. Generally, a ray with a random direction is emitted from the end node of the current path; if it intersects with the scene geometry, the intersection point is taken as the successor node of the original end node. However, in simulating diffraction, the intermediate nodes of the path must be located on the edges of the geometry, and the general ray intersection method cannot guarantee that the intersection point meets this requirement. Therefore, it needs to be modified.
[0110] Assume the scene geometry is described by a triangular mesh. After a ray successfully intersects the geometry, the intersection point of the ray and the geometry is called a pseudo-intersection point. Examine the edges of the triangular face containing the pseudo-intersection point and determine if these edges are convex. Next, randomly select one of these convex edges. Let the triangular face be Δabc1, and the selected edge be ab. Then, denote the probability of selecting this edge as... After this, connect the pseudo-intersection point with point c1 using a straight line and extend it to intersect ab at a point. This point is called the (actual) intersection point corresponding to the pseudo-intersection point, and it is taken as the successor node of the original last node.
[0111] Now we calculate the generation probability of nodes under this method. Let the original path be x0→x1→…→x n The pseudo-intersection is The corresponding intersection point is x n+1 Let ab have two adjacent faces, and let Δabc2 be the other face that is different from Δabc1. Then, all corresponding x... n+1 The pseudo-intersections are all located on line segment x. n+1 c1 and x n+1 Above c2. Next, calculate the two auxiliary variables S(x)n x n+ 1c1) and S(x n x n+1 c2). S(x) n x n+1 The calculation method for c1) is as follows:
[0112] (3.1) Find x n+1 c1 for x n The visible (unobstructed) portion. Since the scene is composed of triangular facets, the visible portion...
[0113] The segment must be a set of non-intersecting line segments; by finding the triangle Δx in the scene geometry... n x n+1 Calculate all occluded line segments by using the intersecting triangular facets of c1, and then from x... n+1 After removing the occluded line segments in cell c1, the remaining set is the set of visible line segments. Let these visible line segments be denoted as _____. in compared to Closer to the c1 end.
[0114] (3.2) Calculate S(x) n x n+1 c1). The formula is
[0115]
[0116] S(x n x n+1 The calculation method for c2) is the same as above.
[0117] Finally, an estimate of the node generation probability is given:
[0118]
[0119] in, For x n The incident path direction is x n-1 →x n At that time, the emission direction is generated. The probability, specifically the probability distribution, is determined by the algorithm implementer. If n = 0, this term must be replaced with the direction of the ray emitted from the ray source. The probability of. For the emitted ray to intersect with the scene The probability of is calculated using the following formula:
[0120]
[0121] in For rays and The angle between the normals of the triangle containing the triangle. C (x n →x n+1 The probability correction term is calculated using the following formula:
[0122]
[0123] in Let Δabc1 be the area. Similarly. Because x n+1 corresponding There is more than one, and the probability mentioned above is not directly equivalent to the generation of node x. n+1 The generation probability, but its expectation is indeed related to the generation probability of x. n+1 The probabilities are equal.
[0124] The path generation probability after adding nodes is as follows:
[0125]
[0126] In the bidirectional path tracking, paths need to be generated simultaneously from both the sound source and the listener. Then, connections are established between the nodes of the two paths to generate a complete path. The probability of generating the complete path is p. F ·p B ·p L Where p F p represents the partial generation probability of the forward subpath from the sound source to the connected point. B p represents the generation probability of the portion of the backward subpath from the listener to the connected point. L The probability of establishing a connection.
[0127] Furthermore, this invention can also employ multiple importance sampling, a method to improve sample quality, which is often used in conjunction with bidirectional path tracing. For a path, multiple importance sampling requires providing the generation probability of this path under all generation methods. Assume that a path x0→x is generated through bidirectional path tracing. i →…→x n Its connecting segments can be any segment of the path, thus there are n-1 possible generation methods. The method for calculating the path generation probability under these methods is as follows: Based on the calculation of acoustic wave straight-edge diffraction and combined with the multiple importance sampling technique, to calculate the probability of various strategies in the multiple importance sampling technique, corresponding forward and backward pseudo-intersections must be generated for each non-starting or ending node. The forward pseudo-intersection is visible to the predecessor node of that node, and the backward pseudo-intersection is visible to the successor node of that node; specifically as follows:
[0128] Let the path be x0→x1→…→x nGenerated by bidirectional path tracing, its connection segment is x. i →x i+1 Then the probability of generating this path is p(x0→x1→…→x). i )p(x i+1 ←x i+2 ←…←x n )p L (x i →x i+l (x here) i+1 ←x i+2 ←…←x n The meaning is equivalent to x n →x n-1 →…→x i+1 The scene locations involved in calculating the path generation probability are listed below:
[0129]
[0130] This is similar to points and Similarly, these are pseudo-intersections, called "backward pseudo-intersections." They correspond to nodes on the reverse path. As can be seen from the process of generating the reverse path, For x i+1 It is evident, but not guaranteed for x. i-1 It can be seen that, with Such pseudo-intersections (called forward pseudo-intersections) are exactly the opposite.
[0131] By changing the value of 'i' in the diagram above, we can see that in order to calculate the generation probability of a connection segment being another path segment, each intermediate node needs corresponding forward and backward pseudo-intersections. Therefore, it is necessary to actively fill in the missing pseudo-intersections in the diagram. The method for filling in forward pseudo-intersections will be given here; the method for filling in backward pseudo-intersections is similar.
[0132] Assume x i+1 If there is no corresponding forward pseudo-intersection, first find x. i+1 Let's denote the edge containing this edge as ab. Let's denote the two adjacent faces of this edge as Δabc1 and Δabc2. Then find the line segment x. i+1 c1 and x i+1 c2 relative to x i The visible portion is represented as a set of line segments, each segment using... Representation. Define the function:
[0133]
[0134]
[0135] Next, for each line segment If it is x i+1 A subset of c1 is then passed through Map it to the set [0,1]; if it is x i+1 A subset of c2 is then passed through First, map it to the set [0, 1]; second, randomly select a point in the mapped set of line segments according to a uniform distribution; finally, use f lerp Map the selected point back to the original set of line segments and use it as the pseudo-intersection point.
[0136] After completing the pseudo-intersections, the path generation probability under all generation methods can be calculated.
[0137] Furthermore, other optimizations are also applicable to this invention:
[0138] When the probability of generating a path is too low, or the path's contribution is too high relative to the probability of generation, extreme values may appear in the results, potentially impacting the user experience. Therefore, it is necessary to identify and remove these extreme values.
[0139] Given a path x0→x1→…→x n+1 Define the auxiliary variables as follows:
[0140]
[0141] Where p V and The definition is the same as in Section 3. If i = 0, p V It also needs to be replaced using the method described in Chapter 3. This variable represents the probability that the path containing this segment will produce extreme values. The C of the entire path o The value is C for each segment of the path. o The product of values. If the path is composed of multiple sub-paths (such as bidirectional path tracing), then its C... o The value is C for each segment path. o The product of values. Next, multiply the final contribution of the path by a value C. s :
[0142]
[0143] Here This is a user-defined variable (usually between 100 and 1000). The multiplier above will suppress... This contributes to all paths, thereby reducing the generation of extreme values.
[0144] If the method for generating successor nodes for a path is not efficient enough, the following method can be used to improve node utilization efficiency in bidirectional path tracing: First, assume that the simulation algorithm generates M paths at a time, and the maximum order of the paths is N-2. Therefore, these paths have a total of M×N nodes. These nodes, regardless of whether they are actually successfully generated, are arranged in an M-N array. For each node in the array, define a variable F. j,i Let j = 1, 2, ..., M, i = 1, 2, ..., N, called the splitting degree, and let F j,1 =1. The nodes in the array need to be generated row by row. For the node in the i-th column of the array, the specific generation process is as follows: Try to generate all nodes in the i-th column using the method in Section 3; for each node x in the (i-1)-th column, count the set S of rows where its successor nodes were successfully generated. + (x) and the set of failed rows S - (x); S + The nodes in each row and i-th column of (x) are copied into S. - In (x), the i-th column of each row. If the node in the j-th row is copied k times, then let F j,i =(k+1)F j,i-1 Let L be the expected contribution of the path to this node. j,i It is recommended that the number of replications per node be related to L. j,i / F j,i-1 Proportional.
[0145] When establishing forward and reverse path connections and generating the final path, the expected contribution of the final path should be divided by the product of the two split degrees corresponding to the nodes connecting the forward and reverse paths. This method fully utilizes empty spaces in the array, saves storage space, and improves node connection efficiency.
[0146] Implementation Examples
[0147] The inventors implemented the algorithm described above on a computer equipped with an Intel i7 2.60GHz quad-core CPU and 12GB of RAM. Specifically, they chose bidirectional path tracing incorporating multiple importance sampling and utilized all the optimization techniques described in Section 5 of the implementation process. The invention compared its results with simulation results using UTD and BTM models.
[0148] Figure 4 and Figure 6 The first-order and second-order diffraction scenarios used for comparison are shown. Figure 4 The scene in the image contains a right dihedral with a straight side length of 2m. Figure 6 The scene in the image contains two parallel right dihedral angles, and the path from the sound source must undergo two diffractions before reaching the listener. Both scene geometries use Dirichlet boundary conditions. Figure 5 and Figure 7 The impulse response calculated using different methods in the scenario is shown. These methods include BTM, UTD, the method presented in this paper, and the results reconstructed from eight frequency points selected in UTD (located between 62.5-8000Hz, distributed at exponential intervals). The UTD method suffers from significant errors in calculating second-order and higher-order diffraction. Figure 7 Not shown in the figure. As can be seen, the computational quality of the algorithm in this paper for first and second order diffraction is not significantly different from existing methods, and is significantly better than the common simplified UTD scheme (which only calculates results at 8 frequencies). Since the BTM method is accurate for the solution of first-order diffraction, this comparison also verifies the correctness of the method in this invention. For large-scale complex scenes, this method exhibits near real-time computational efficiency. Using 12,000 paths, for a model containing 22,974 triangular faces, the time cost of calculating the diffraction impulse response in this invention is approximately 140 ms. For a model containing 279,133 faces, the time cost is approximately 205 ms. Therefore, the method of this invention is suitable for fast diffraction result calculation in complex scenes.
[0149] The above-described embodiments illustrate specific implementations of the present invention. Their detailed descriptions are intended to aid in understanding the method and core ideas of the present invention, but should not be construed as limiting the scope of the present invention. It should be noted that those skilled in the art can make various modifications and improvements without departing from the concept of the present invention, and these all fall within the scope of protection of the present invention. Therefore, the scope of protection of this patent should be determined by the appended claims.
Claims
1. A method for calculating straight-edge diffraction of acoustic waves based on path tracking, characterized in that, Includes the following steps: (1) Read the acoustic environment state related to diffraction, the acoustic environment state including the position of the sound source and the listener, the geometric description of the scene, and the interface features of the surface of the geometric body; (2) A large number of diffraction paths are randomly generated in the scene. The path is composed of a series of nodes in the scene and line segments connecting the nodes. The starting node of the path is the sound source, the ending node is the listener, and the remaining intermediate nodes are located on the convex straight edge of the geometry in the scene. Except for the nodes that constitute the path, the other parts of the path do not intersect with the surface of the geometry in the scene. (3) Calculate the final contribution of each path, assuming the path contains There are 10 nodes, each numbered sequentially from the sound source to the listener. arrive Let the line segments between nodes be denoted as . The complete path is Then the total contribution The calculation formula is as follows: ; in, For the characteristics of the sound source, for The acoustic medium characteristics at that location, for The radial function representation of diffraction at the straight edge; (4) Add the contributions of all paths to the impulse response function: for path Its contribution is known to be in step (3) And denote its generation probability as The propagation delay is Then, the impulse response function is accumulated on top of the original impulse response function. ,in The unit impulse function; (5) Output the impulse response function generated in step (4) as the diffraction simulation result.
2. The method for calculating acoustic wave straight-edge diffraction based on path tracking according to claim 1, characterized in that, Step (2) includes the following sub-steps: (2.1) Initialize the path so that the only node it contains is the starting point of the path; (2.2) Let the current path status be ,from A ray is randomly emitted from the starting point; if If it is not the starting point of the path, the exit angle of the ray is determined by a sampling function that matches the radial function of the straight edge diffraction. (2.3) Let the first intersection point between the ray and the geometric surface of the scene be denoted as . , and called it a pseudo-intersection; If the point does not exist, the intersection will fail, and the process will proceed to step (2.8). (2.4) Pseudo-intersections in step (2.3) Randomly select one of the three sides of the triangle; let the triangle be... The selected edge is The vertex opposite to this edge is Then connect with a straight line. and And extend, with They intersect at a point, denoted as This is called the intersection point; (2.5) If the edge in step (2.4) It is not a convex straight edge or a line segment. If the intersection with the scene geometry occurs outside of the two endpoints, the intersection calculation fails, and we proceed to step (2.8); otherwise, we find the intersection point. Add as a node to the path At the end, a new path is formed. ; (2.6) Calculation path The generation probability: if the path is known The generation probability is Then the path The generation probability calculation method is as follows in, for The incident path direction is At that time, the direction of the emitted ray is The probability, if Then this item must be replaced with the direction of the radiation emitted from the radiation source. The probability of; For rays Intersection with scene geometry The probability of; This is a probability correction term; (2.7) If the path length does not meet the user's requirements, proceed to step (2.2). (2.8) Use the current path as the diffraction path in the scene.
3. The method for calculating straight-edge diffraction of acoustic waves based on path tracking according to claim 1, characterized in that, In step (2), paths are generated simultaneously from both the sound source and the listener. The path originating from the sound source is called the forward sub-path, and the path originating from the listener is called the backward sub-path. Then, connections are randomly established between the nodes of the forward and backward sub-paths to form a complete path from the sound source to the listener. Let the probability of generating the portion of the forward sub-path from the sound source to the connected point be 1. The probability of generating the backward subpath from the listener to the connected point is: The probability of establishing a connection is Then the probability of generating the complete path is .
4. The method for calculating straight-edge diffraction of acoustic waves based on path tracking according to claim 1, characterized in that, In step (5), the path contribution needs to be multiplied by a factor used to suppress extreme values. : ; in For a user-defined variable, For each segment randomly generated in the path The product of each path segment A parameter used to measure the likelihood that a path containing this path segment will produce extreme values.
5. The method for calculating straight-edge diffraction of acoustic waves based on path tracking according to claim 3, characterized in that, In step (3), the contribution needs to be divided by the corresponding position of the forward and reverse path connection nodes. The product of elements; the stated For one OK An array of columns, where Let be the maximum number of nodes in the forward and backward paths; the nodes of the path are stored in an array with the same number of rows and columns, where each row corresponds to a single path, and let the array... The first column contains only 1s, and the first column of the node array represents the path's starting point; in the... If column nodes have already been generated, the first one in the array The process of generating column nodes includes the following sub-steps: (5.1) From the first row of each node array Nodes It departs, randomly shooting out a ray; like The exit angle of the ray is then determined using a sampling function that matches the radial function of the straight-edge diffraction. (5.2) Let the first intersection point between the ray and the scene geometry be denoted as . , and called it a pseudo-intersection; If the point does not exist, the intersection fails, and the process proceeds to step (5.7). (5.3) Pseudo-intersections in step (5.2) Randomly select one of the three sides of the triangle; let the triangle be... The selected edge is The vertex opposite to this edge is Then connect with a straight line. and And extend, with They intersect at a point, denoted as This is called the intersection point; (5.4) If the edge in step (5.3) Non-convex straight edge, or line segment If the intersection with the scene geometry occurs outside of the two endpoints, the intersection calculation fails, and proceed to step (5.7); otherwise, the intersection point is... The path added to the current line as a node. At the end, a new path is formed. ; (5.5) Calculation path The generation probability: if the path is known The generation probability is Then the path The generation probability is calculated as follows: in, for The incident path direction is At that time, the direction of the emitted ray is The probability, if Then this item must be replaced with the direction of the radiation emitted from the radiation source. The probability of; For rays Intersection with scene geometry The probability of; This is a probability correction term; (5.6) If the path length does not meet the user's requirements, proceed to step (5.1) to prepare to generate the next column of nodes; (5.7) For the first Each node of the column Calculate the set of rows where successor nodes were successfully generated. and the set of failed rows ; (5.8) will Each line in Copy the column nodes into The first of each bank Column; if the first Row nodes are copied Next, let .
Citation Information
Patent Citations
Real-time sound propagation simulation method based on bidirectional path tracking
CN105893719A
Method and system for virtual acoustic rendering by time-varying recursive filter structures
CN113348681A