A ray-tracing channel modeling method considering antenna and material depolarization effects

Abstract

This invention discloses a ray-tracing channel modeling method considering antenna and material depolarization effects, relating to the field of wireless communication technology. The method includes: 1) acquiring a discrete antenna pattern database and reconstructing a continuous antenna radiation field through scalar spherical harmonic expansion; 2) running ray-tracing simulations to extract the propagation path and material interaction information of the rays; 3) calculating the ray departure angle and arrival angle to obtain the radiated electric field at the transmitter and receiver at the corresponding angles; 4) identifying the materials interacting with the rays and calculating the propagation coefficient matrix based on the simplified generalized propagation matrix method; and 5) obtaining the total received electric field at the receiver and calculating and outputting channel characteristic information. This invention significantly improves the computational efficiency of ray-tracing channel modeling methods under complex material conditions and has high simulation accuracy for channel depolarization effects.

Description

Technical Field

[0001] This invention relates to the field of wireless communication technology, and in particular to a ray tracing channel modeling method that considers antenna and material depolarization effects. Background Technology

[0002] In sixth-generation (6G) mobile communication systems, polarization, as a critical physical resource, is crucial for meeting the demands of high reliability, low latency, and massive connectivity. Polarization multiplexing technology can significantly enhance the capacity of communication systems. However, in real-world wireless propagation environments, factors such as non-ideal antenna design, polarization mismatch caused by antenna rotation, and obstruction by building structures can significantly alter the polarization state of electromagnetic waves, leading to severe channel depolarization effects. Depolarization disrupts channel orthogonality and introduces cross-polarization interference, resulting in a sharp decline in communication system performance. To efficiently and accurately design and optimize communication systems, it is essential to establish a high-precision channel model capable of jointly characterizing the antenna and wireless propagation channel depolarization effects. Compared to traditional statistical channel models, ray tracing (RT), as a deterministic method, can provide high-precision polarization information at specific locations, thereby enabling accurate modeling of channel depolarization effects.

[0003] In ray tracing channel models, the antenna patterns at the transceiver end are typically derived from full-wave simulations or antenna pattern measurements and presented as discrete data. To characterize the depolarization effect caused by antenna rotation, the traditional approach is to transform the discrete antenna pattern from the antenna coordinate system to the global coordinate system and resample it to align the antenna pattern with the angular information from the ray tracing. However, during pattern resampling, traditional sampling methods such as weighted summation or bilinear interpolation cannot effectively handle the singularity problem at the poles of the spherical coordinate system and are prone to over-smoothing of features, thus reducing the accuracy of pattern reconstruction. Furthermore, matching the discrete antenna pattern with rays at arbitrary angles inevitably introduces quantization errors. Therefore, a precise and efficient method is urgently needed to address the depolarization modeling problem caused by antenna rotation in ray tracing.

[0004] On the other hand, the influence of the propagation environment on the polarization state of the received signal is equally significant. Existing ray-tracing channel models that consider channel depolarization effects are mostly based on the isotropic material assumption to analyze the depolarization effects generated by propagation mechanisms such as reflection, transmission, scattering, and diffraction. However, the isotropic assumption cannot accurately characterize the electromagnetic property direction dependence that may exist in complex materials (such as multilayer anisotropic materials) in real-world environments, leading to significant deviations in channel depolarization modeling. Although multimode ray tracing and generalized propagation matrix methods have been proposed for calculating the depolarization effects of multilayer anisotropic materials, both face the technical bottleneck of high computational complexity. Faced with the severe challenges posed by complex materials, existing ray-tracing methods struggle to balance modeling accuracy with computational efficiency. There is an urgent need to develop novel ray-tracing channel models that can efficiently and accurately characterize antenna rotation and the depolarization effects of complex materials to meet the high-precision modeling requirements of channel depolarization effects in complex environments. Summary of the Invention

[0005] The technical problem to be solved by this invention is to overcome the shortcomings of existing ray tracing channel modeling methods that cannot simultaneously characterize antenna rotation and the depolarization effect of complex materials. Instead, it provides a ray tracing channel modeling method that considers the depolarization effects of antennas and materials, accurately characterizes the mapping relationship between antenna rotation and multilayer anisotropic materials and the depolarization effect of wireless channels, and achieves a breakthrough in channel modeling accuracy and computational efficiency.

[0006] To solve the above-mentioned technical problems, the present invention adopts the following technical solution:

[0007] A ray tracing channel modeling method considering antenna and material depolarization effects proposed according to the present invention includes:

[0008] Step S1: Import the simulation scenario and set the transceiver antenna configuration and simulation parameters.

[0009] Step S2: Based on the antenna configuration and simulation parameters of the transceiver end, obtain the radiation pattern of the discrete antenna at the transceiver end, and reconstruct the radiation field of the continuous antenna at the transceiver end through scalar spherical harmonic expansion;

[0010] Step S3: Run the ray tracing channel simulation according to the simulation scene, the location of the transceiver and the simulation parameters, and extract the propagation path of each ray and the material interaction information between the ray and the scene;

[0011] Step S4: Based on the propagation path of the rays, calculate the departure angle at the transmitting end and the arrival angle at the receiving end for each ray, and obtain the radiated electric field at the transmitting end and the radiated electric field at the receiving end based on the continuous antenna radiation field at the transmitting and receiving ends according to the departure angle and the arrival angle.

[0012] Step S5: Identify the ray interaction materials corresponding to each ray based on the material interaction information, and calculate the electromagnetic propagation characteristics of the material based on the generalized propagation matrix method for the ray interaction materials to obtain the reflection coefficient matrix and transmission coefficient matrix that match the ray tracing model.

[0013] Step S6: Calculate the receiving electric field corresponding to each ray based on the radiated electric field at the transmitting end and the radiated electric field at the receiving end, as well as the reflection coefficient matrix and the transmission coefficient matrix. Coherently superimpose the receiving electric fields corresponding to all rays to obtain the total receiving electric field at the receiving end. Calculate and output the channel characteristic information based on the ray propagation path and the total receiving electric field at the receiving end.

[0014] As a further optimization scheme for the ray tracing channel modeling method that considers the depolarization effect of antenna and materials described in this invention, the simulation parameters include the transmitter and receiver positions, simulation frequency, propagation mechanism, material parameters, antenna attitude, and spherical harmonic expansion order.

[0015] As a further optimization scheme of the ray tracing channel modeling method considering antenna and material depolarization effects described in this invention, step S2 includes:

[0016] Step S201: Obtain the sampling data of the discrete antenna pattern of the transceiver end, extract the zenith angle component and azimuth angle component of the electric field in the spherical coordinate system for each sampling direction, and construct the coordinate transformation relationship according to the zenith angle and azimuth angle corresponding to each sampling direction to convert the zenith angle component and azimuth angle component of the electric field into the electric field sampling vector in the Cartesian coordinate system.

[0017] Step S202: Based on the sampling directions corresponding to the discrete antenna pattern sampling data and the maximum order of the spherical harmonic expansion, construct the spherical harmonic basis function vector corresponding to each sampling direction, and construct the spherical harmonic matrix based on the spherical harmonic basis function vector;

[0018] Step S203: Based on the electric field sampling vector and spherical harmonic matrix in the Cartesian coordinate system, establish a linear least squares problem, and solve the vector expansion coefficient matrix by calculating the pseudo-inverse of the spherical harmonic matrix, thereby reconstructing the continuous antenna radiation field in any direction.

[0019] As a further optimization scheme for the ray tracing channel modeling method considering antenna and material depolarization effects described in this invention, step S201 includes:

[0020] Let the discrete antenna pattern include The sampling direction, where the first... Each sampling direction is represented as ,and , Indicates the first The zenith angle of each sampling direction in spherical coordinates Indicates the first The azimuth angle of each sampling direction in spherical coordinates;

[0021] No. The electric field sampling data in spherical coordinates corresponding to each sampling direction is represented as follows:

[0022] ;

[0023] in, Indicates the first The electric field sampling vector in spherical coordinates corresponding to each sampling direction Indicates the first The zenith angle component of the electric field corresponding to each sampling direction Indicates the first The electric field azimuth component corresponding to each sampling direction;

[0024] Will The electric field sampling data in spherical coordinates corresponding to each sampling direction are arranged as follows:

[0025] ;

[0026] in, Indicates by The electric field sampling matrix is ​​composed of the electric field sampling vectors in the spherical coordinate system corresponding to each sampling direction;

[0027] Will The electric field sampling data corresponding to each sampling direction in the spherical coordinate system is transformed to the Cartesian coordinate system, and the specific expression is as follows:

[0028] ;

[0029] in, Indicates by The electric field sampling matrix is ​​composed of the electric field sampling vectors in the Cartesian coordinate system corresponding to each sampling direction. This represents the transformation matrix between spherical coordinates and Cartesian coordinates.

[0030] As a further optimization scheme of the ray tracing channel modeling method considering antenna and material depolarization effects described in this invention, step S202 includes:

[0031] Let the maximum order of the pre-defined spherical harmonic expansion be... The total number of basis functions contained in the spherical harmonic expansion For the first There are spherical harmonic basis functions, whose corresponding spherical harmonic order and mode are respectively expressed as: and ,in , , ;

[0032] according to Constructing the spherical harmonic matrix from the spherical harmonic basis function vectors corresponding to each sampling direction. spherical harmonic matrix Represented as:

[0033] ;

[0034] in, Representing the spherical harmonic matrix The The sampling direction of the first A spherical harmonic basis function For the reason Arrangement OK Column matrix; Represented as:

[0035] ;

[0036] in, Represents the imaginary unit. Indicates the first Each sampling direction Next Rank, number scalar spherical harmonic functions of the mode; The expression is:

[0037] ;

[0038] in, Let e ​​represent the Legendre function of the first kind of correlation, where e is the natural base;

[0039] Define the vector expansion coefficient matrix in the Cartesian coordinate system. The vector expansion coefficient matrix is ​​represented as:

[0040] ;

[0041] in, , , and These represent the x-axis, y-axis, and z-axis in the Cartesian coordinate system, respectively. In Cartesian coordinate system Directional electric field component Vector expansion coefficients on spherical harmonic basis functions For the reason Arrangement OK Column matrix.

[0042] As a further optimization scheme for the ray tracing channel modeling method considering antenna and material depolarization effects described in this invention, step S203 includes:

[0043] Based on the electric field sampling matrix in the Cartesian coordinate system obtained in step S201 The spherical harmonic matrix obtained in step S202 The linear least squares problem is established, specifically represented as:

[0044] ;

[0045] in, Represents the coefficient matrix of the vector expansion The least squares solution, Representing the Frobenius norm; by calculating the spherical harmonic matrix. The pseudo-inverse matrix is ​​used to solve for it. Represented as:

[0046] ;

[0047] in, Representing the spherical harmonic matrix The pseudo-inverse matrix;

[0048] For any direction to be reconstructed , This represents the zenith angle of the direction to be reconstructed in spherical coordinates. The azimuth angle of the direction to be reconstructed in spherical coordinates is defined as the direction to be reconstructed. The corresponding column vector of spherical harmonic basis functions The column vector of spherical harmonic basis functions is represented as:

[0049] ;

[0050] in, Represents the column vector of spherical harmonic basis functions The One spherical harmonic basis function;

[0051] Represented as:

[0052] ;

[0053] in, Indicates the first The order of the spherical harmonics corresponding to each spherical harmonic basis function. Indicates the first The spherical harmonic function modes corresponding to each spherical harmonic basis function. , , , Indicates any direction to be reconstructed Next Rank, number scalar spherical harmonic functions of the mode;

[0054] The least squares solution based on the vector expansion coefficient matrix Column vector of spherical harmonic basis functions The continuous antenna radiation field vector in Cartesian coordinates under the direction to be reconstructed is obtained. Specifically, it is expressed as:

[0055] .

[0056] As a further optimization scheme for the ray tracing channel modeling method considering antenna and material depolarization effects described in this invention, step S4 includes:

[0057] Step S401: Based on each ray propagation path, obtain the departure angle at the transmitting end and the arrival angle at the receiving end for each propagating ray in the global coordinate system, and determine the radial unit vector at the transmitting end corresponding to the departure angle and the radial unit vector at the receiving end corresponding to the arrival angle.

[0058] Step S402: Based on the transmit and receive antenna attitude determined in step S1 and the actual rotation angle of the transmit and receive antenna in three-dimensional space, calculate the unit rotation vector of each coordinate axis of the transmit antenna coordinate system and the receive antenna coordinate system relative to the global coordinate system.

[0059] Step S403: Using the unit rotation vector obtained in step S402, transform the radial unit vector leaving the transmitter and the radial unit vector arriving at the receiver obtained in step S401 to the transmitter antenna coordinate system and the receiver antenna coordinate system, respectively. Solve for the departure zenith angle and departure azimuth angle in the transmitter antenna coordinate system, and the arrival zenith angle and arrival azimuth angle in the receiver antenna coordinate system. Based on the vector expansion coefficient matrix obtained in step S2, calculate the continuous antenna radiation field of the transmitter corresponding to the departure zenith angle and departure azimuth angle, and the radiation electric field of the transmitter and the radiation electric field of the receiver corresponding to the arrival zenith angle and arrival azimuth angle.

[0060] Step S401 includes:

[0061] Based on the propagation rays between the transceivers, obtain the departure angle of each ray in the global coordinate system. and angle of arrival ,in and These represent the zenith angle and azimuth angle of the propagating ray in the global coordinate system at the transmitting end, respectively. and These represent the zenith angle and azimuth angle of the global coordinate system at the receiving end, respectively;

[0062] The radial unit vector corresponding to each ray is calculated based on its departure angle and arrival angle, expressed as follows:

[0063] ;

[0064] ;

[0065] in This represents the radial unit vector corresponding to the ray's departure angle in the global coordinate system. This represents the radial unit vector corresponding to the angle of arrival of a ray in the global coordinate system.

[0066] Step S402 includes:

[0067] The departure angle and arrival angle in the global coordinate system are transformed to the antenna coordinate system and matched with the antenna pattern. The antenna coordinate system is defined as follows:

[0068] ;

[0069] ;

[0070] in, In the antenna coordinate system Axial unit vector, Representing the antenna coordinate system Axial unit vector, , and These represent the antenna relative to the global coordinate system. axis, shaft and Rotation angle of the axis; calculate the departure angle and arrival angle of each ray in the antenna coordinate system based on the unit vector of the antenna coordinate system. :

[0071] ;

[0072] in, and These represent the zenith angle and azimuth angle of the propagating ray in the global coordinate system at the transmitting end, respectively. Represents the inverse cosine function. The radial unit vector representing the ray departure angle in the global coordinate system in step S401. Or the radial unit vector of the angle of arrival in the global coordinate system The superscript T indicates transpose;

[0073] Step S403 includes:

[0074] Determine the departure angle and arrival angle of each ray in the antenna coordinate system. Then, the antenna radiation field is calculated based on the vector expansion coefficient matrix and spherical harmonic matrix from step S2:

[0075] ;

[0076] in, This represents the antenna radiation field vector corresponding to the departure angle and arrival angle of each ray in the antenna coordinate system. This represents the least-squares solution of the vector expansion coefficient matrix obtained in step S2. This represents a column vector consisting of the spherical harmonic basis function values ​​corresponding to the departure angle and arrival angle of each ray in the antenna coordinate system.

[0077] As a further optimization scheme for the ray tracing channel modeling method considering antenna and material depolarization effects described in this invention, step S5 includes:

[0078] Step S501: Extract the incident angle at the interaction point between the ray and the multilayer anisotropic material and the electromagnetic parameters of the multilayer material. Construct a coordinate system rotation matrix according to the direction of incident wave propagation, rotate and align the dielectric constant and permeability of the anisotropic material, and obtain the dielectric constant and permeability of the rotated anisotropic material.

[0079] Step S502: Construct the coupling matrix of each layer of anisotropic material based on the dielectric constant and magnetic permeability of the rotated anisotropic material, and directly solve the eigenvalues ​​and eigenvectors of the corresponding coupling matrix according to the preset anisotropic material category.

[0080] Step S503: Using the boundary continuity conditions of the anisotropic materials in each layer, the eigenvalues ​​and eigenvectors of the anisotropic materials in each layer are combined into an eigenmatrix and concatenated to obtain the generalized propagation coefficient matrix. The generalized propagation coefficient matrix is ​​then converted into a reflection coefficient matrix and a transmission coefficient matrix that match the ray tracing channel model.

[0081] As a further optimization scheme for the ray tracing channel modeling method considering antenna and material depolarization effects described in this invention,

[0082] Step S501 includes:

[0083] Based on the propagation ray between the transceiver, information on the interaction points between the ray and obstacles is obtained, including the coordinates of the interaction points and the material parameters of the interaction. The surface material of the obstacle is assumed to be a uniformly thick plate. The direction dependence of the electromagnetic parameters of the anisotropic material is considered, along with the dielectric constant of the anisotropic material. and permeability Represented in tensor form, as shown below:

[0084] ;

[0085] in, , and Defined respectively as anisotropic materials Optical axis Optical axis and Optical axis The dielectric constant of anisotropic materials is represented by... Components along the optical axis The dielectric constant of anisotropic materials is represented by... Components along the optical axis The dielectric constant of anisotropic materials is represented by... Components along the optical axis This indicates the electrical conductivity of anisotropic materials at... Components along the optical axis This indicates the electrical conductivity of anisotropic materials at... Components along the optical axis This indicates the electrical conductivity of anisotropic materials at... Components along the optical axis;

[0086] Coordinate rotation is used to align the optical axis of the anisotropic material with the incident coordinate system, where the incident coordinate system is defined by the propagation direction of the incident wave. The transmitted wave in the anisotropic material is decomposed into two types, referred to as Type I and Type II waves. In the incident coordinate system, the wave vector... Simplified to ,in Indicates the wave vector in Components along the optical axis Indicates the wave vector in The components along the optical axis, and , Represents the angle between the incident ray and the normal vector of the medium surface; rotation matrix Represented as:

[0087] ;

[0088] in, This represents the counterclockwise rotation angle of the z-axis in the coordinate system.

[0089] Dielectric constant of the rotated anisotropic material and permeability Represented as:

[0090] ;

[0091] ;

[0092] in, and They represent and The line, number Column parameters, , ;

[0093] Step S502 includes:

[0094] For each layer in a multilayer anisotropic material, the dielectric constant of the anisotropic material after rotation of that layer is used. and permeability Directly construct the coupling matrix of this layer :

[0095] ;

[0096] in, Indicates angular frequency;

[0097] Based on coupling matrix Anisotropic materials are classified, and the coupling matrix corresponding to each type of anisotropic material is calculated. The eigenvalues ​​and eigenvectors; based on the anisotropic characteristics of materials, the electromagnetic properties of anisotropic materials are divided into five categories: isotropic materials, The axis is a uniaxial material under the condition of optical axis. The axis is a uniaxial material under the condition of optical axis. The axis is uniaxial and biaxial materials under optical axis conditions;

[0098] When anisotropic materials are isotropic materials or When the axis is the optical axis and the material is uniaxial, the coupling matrix is... The eigenvalues ​​and eigenvectors are represented as follows:

[0099] ;

[0100] ;

[0101] in, Represents the coupling matrix under the current material conditions. The first and second eigenvalues; Represents the coupling matrix under the current material conditions. The third and fourth eigenvalues; Indicates that it is made of isotropic materials or The eigenvector matrix of a uniaxial material under the condition that the axis is the optical axis;

[0102] When anisotropic materials are When the axis is the optical axis, for uniaxial or biaxial materials, the coupling matrix eigenvalues for:

[0103] ;

[0104] Coupling matrix eigenvectors for:

[0105] ;

[0106] ;

[0107] in, Represents the coupling matrix The 1 eigenvector Represents the coupling matrix under the current material conditions. The 1 eigenvalue, The material is an isotropic material or The analytical solution for uniaxial materials with the optical axis as the axis only exists in Established at that time, Let z be the counterclockwise rotation angle of the optical axis in the coordinate system; when At that time, isotropic materials or The eigenvalues ​​and eigenvectors of a uniaxial material with the axis being the optical axis are derived from the characteristics of isotropic materials or The analytical solution for uniaxial materials is determined when the axis is the optical axis.

[0108] Step S503 includes:

[0109] The electromagnetic wave propagation characteristics of single-layer anisotropic materials are extended to those of multilayer anisotropic materials; for the first... The tangential component of the electromagnetic field in a layered anisotropic material. Represented as:

[0110] ;

[0111] , ;

[0112] in, This represents the total thickness of a multilayer anisotropic material. Indicates the first The distance between the anisotropic material layer and the incident plane. Indicates the first Tangential components of the electromagnetic field in anisotropic materials. Indicates the first The eigenvector matrix of the anisotropic material is obtained by solving step S502. express The Middle Line 1 Column parameters, , , This indicates the position of each eigenwave in anisotropic materials. The magnitude column vector at that point, express The Middle The parameters of the row, , The column vector order corresponds to the eigenvalues The order in which They represent the first The first, second, third, and fourth eigenvalues ​​of the layer; The eigenvalues ​​of the layer are arranged according to the direction of electromagnetic wave propagation: and These correspond to the ascending waves of type I and type II waves, respectively. and These correspond to the descending waves of Type I and Type II waves, respectively; the ascending and descending waves are determined by the direction of the Poynting vector. The definition is as follows:

[0113] ;

[0114] in, Represents electric field, Indicates magnetic field, This represents the conjugate transpose of the magnetic field. The eigenvalue represents the real part; when the Poynting vector corresponding to a certain eigenvalue is greater than 0, the eigenvalue represents an upward wave; when the Poynting vector corresponding to a certain eigenvalue is less than 0, the eigenvalue represents a downward wave.

[0115] Based on the material information obtained after the ray interacts with the obstacle, a multi-layer anisotropic material model is constructed for the current interaction point. Layer 0 and layer (n+1) represent air, and the region between these two layers represents n layers of anisotropic material. The number of layers n is determined by the actual number of material layers in the environment. The propagation equation between layer 0 and layer (n+1) is shown below:

[0116] ;

[0117] in,

[0118] ,

[0119] This indicates the position of the anisotropic material in layer 0. The eigenwave amplitude column vector at that location, , Indicates the first Layered anisotropic materials and the first The transformation matrix between layers of anisotropic materials. , Indicates the first The propagation matrix within the anisotropic material. Indicates the first Layered anisotropic materials and the first The transformation matrix between layers of anisotropic materials. This indicates the position of each eigenwave in anisotropic materials. The magnitude column vector at that point, Indicates the first The thickness of the anisotropic material layer, This represents the forward propagation matrix of the anisotropic material from layer 0 to layer (n+1). That is The Middle Line 1 The corresponding parameters are listed below. , , Indicates the first The eigenvector matrix of anisotropic layered materials, Indicates the first The eigenvector matrix of anisotropic layered materials;

[0120] Since there are no reflected waves in the (n+1)th layer, that is... , and They represent the first The amplitude coefficients of the first and second eigenwaves in the layer; and They represent the first The eigenvalues ​​corresponding to the first and second eigenwaves in the layer are determined; the forward propagation matrix is ​​decomposed, and the generalized reflection coefficient matrix is ​​calculated. and generalized transmission coefficient matrix :

[0121] ;

[0122] ;

[0123] in, and They represent the first Line 1 Generalized reflection and transmission coefficients of the array; converting the generalized propagation coefficient matrix into a usable reflection and transmission coefficient matrix for RT; usable reflection coefficient matrix for RT. and transmission coefficient matrix Represented as:

[0124] ;

[0125] ;

[0126] , , and These represent the generalized reflection coefficient matrix. The elements in the first row and first column, the first row and second column, the second row and first column, and the second row and second column; , , and These represent the generalized reflection coefficient matrix. The elements in the first row and first column, the first row and second column, the second row and first column, and the second row and second column; , , , , and These represent the eigenvector matrices corresponding to the 0th layer of the medium. The element at the corresponding row and column position; and They represent the first eigenvector matrix corresponding to the layer medium The elements at the corresponding row and column positions; where the 0th layer is an air layer and located on the incident side of the multilayer anisotropic material, the 1st layer... The medium is an air layer and is located on the transmission side of a multilayer anisotropic material.

[0127] As a further optimization scheme for the ray tracing channel modeling method considering antenna and material depolarization effects described in this invention, step S6 includes:

[0128] Step S601: Distinguish between line-of-sight and non-line-of-sight paths based on the ray propagation path information, and calculate the channel propagation matrix for each ray for different propagation paths;

[0129] Step S602: Combine the continuous radiated electric field vectors of the transceiver at the corresponding angle obtained in step S4, perform polarization matching on the channel propagation matrix, and calculate the received electric field of a single propagation ray.

[0130] Step S603: Coherently superimpose the received electric fields corresponding to all propagation paths and calculate the channel characteristic information;

[0131] Step S601 includes:

[0132] Based on whether the ray interacts with the propagation environment, the propagation path is divided into line-of-sight (LoS) paths and non-line-of-sight (NLoS) paths. For LoS paths, the propagation electric field reaching the receiver is calculated using the free-space propagation attenuation of the radiated electric field at the transmitter, and the channel propagation matrix is ​​used. It is considered an identity matrix; in contrast, electric field calculations for non-line-of-sight paths involve a series of cascaded matrix transformations, including... Step reflection and Channel propagation matrix of first-order transmission Represented as:

[0133] ;

[0134] in, and They represent the first The reflection coefficient matrix at the nth reflection point and the nth reflection coefficient matrix The transmission coefficient matrices at each transmission point are obtained from step S5; angle Indicates the first The rotation angle between adjacent material coordinate systems at each reflection point, angle Indicates the first Rotation angle between adjacent material coordinate systems at each point of transmission action; channel propagation matrix The order of the matrix cascades is determined by the interaction between the rays and the propagation environment;

[0135] Step S602 includes:

[0136] For each propagation path, at the known departure angle of the ray and angle of arrival Under these conditions, the received electric field The expression is:

[0137] ;

[0138] in, and These represent the zenith angle and azimuth angle in the antenna coordinate system corresponding to the propagation path at the transmitting end, respectively. and These represent the zenith angle and azimuth angle in the antenna coordinate system corresponding to the propagation path at the receiving end, respectively. Indicates wavelength. and These represent the gains of the transmitting and receiving antennas, respectively. Indicates the transmission power. It is the receiving antenna impedance, reflecting the conjugate matching characteristics of the load impedance; and The radiated electric fields of the transmitting antenna and the receiving antenna at the departure angle and the arrival angle of the ray, respectively, are calculated by step S4;

[0139] Step S603 includes:

[0140] The received electric fields corresponding to all propagation paths are coherently superimposed to calculate the total received power at the receiver. Based on the common polarization and cross-polarization received electric field components, the channel's cross-polarization ratio, delay spread, and angle spread channel characteristic information are calculated and output.

[0141] Compared with the prior art, the present invention, employing the above technical solution, has the following technical effects:

[0142] This invention addresses the problem of insufficient accuracy caused by oversimplification in traditional ray tracing channel modeling methods when dealing with depolarization effects in real-world environments. It proposes a ray tracing channel modeling method that considers antenna and material depolarization effects, which significantly improves the computational efficiency of ray tracing channel modeling under complex material conditions and has high accuracy in modeling channel depolarization effects, thus meeting the requirements for wireless propagation channel modeling that considers polarization. Attached Figure Description

[0143] Figure 1 This is a flowchart of the method of the present invention;

[0144] Figure 2 This is a schematic diagram of a simulation scene;

[0145] Figure 3 A comparison of the received power calculated by different algorithms when the transmitter is rotated at a 45° angle;

[0146] Figure 4 The received power calculated by full-wave simulation in an anisotropic material scenario is compared with the received power calculated by the method of the present invention; wherein, (a) the received power calculated by full-wave simulation, and (b) the received power calculated by the method of the present invention.

[0147] Figure 5 The diagrams show a comparison of the channel cross-polarization ratio calculated by the method of the present invention in isotropic and anisotropic scenarios; where (a) represents the isotropic scenario and (b) represents the comparison of the channel cross-polarization ratio calculated by the method of the present invention in the anisotropic scenario. Detailed Implementation

[0148] To make the objectives, technical solutions, and advantages of this invention clearer, the invention will be described in detail below with reference to the accompanying drawings and specific embodiments.

[0149] like Figure 1 As shown, a ray tracing channel modeling method considering antenna and material depolarization effects includes:

[0150] Step S1: Import the simulation scenario, set the transceiver positions, and determine the simulation parameters;

[0151] Step S2: Obtain the discrete antenna pattern and reconstruct the continuous antenna radiation field through scalar spherical harmonic expansion;

[0152] Step S3: Run ray tracing simulation to extract the propagation path of the rays and the information on material interaction.

[0153] Step S4: Calculate the departure angle and arrival angle of the ray, and obtain the radiation electric field at the transmitting and receiving ends at the corresponding angles;

[0154] Step S5: Identify the ray-interacting material and calculate the propagation coefficient matrix based on the simplified generalized propagation matrix method;

[0155] Step S6: Obtain the total received electric field at the receiver and calculate and output the channel characteristic information.

[0156] Step S2 includes:

[0157] Step S201: Obtain the sampling data of the discrete antenna pattern of the transceiver end, extract the zenith angle component and azimuth component of the electric field of each sampling point in the spherical coordinate system, and use the transformation matrix to uniformly transform them to the Cartesian coordinate system.

[0158] Step S202: Extract each sampling angle of the discrete antenna pattern, and construct the spherical harmonic vector matrix corresponding to the sampling angle by combining the preset spherical harmonic function expansion order;

[0159] Step S203: Construct a cost function based on the least squares method. By calculating the pseudo-inverse of the spherical harmonic vector matrix constructed in step S202, the optimal vector expansion coefficient matrix is ​​obtained, thereby reconstructing the antenna radiation field at any angle.

[0160] Step S201 includes:

[0161] consider Point antenna radiation pattern measurement data, the data is represented in spherical coordinates as follows:

[0162]

[0163] in, and Let denot zenith angular component and azimuth angular component of the electric field, respectively. The first component in the antenna pattern... Use each angle and direction To facilitate solving the vector coefficient matrix, the antenna pattern measurement data in spherical coordinates needs to be transformed into Cartesian coordinates, expressed as follows:

[0164]

[0165] in This represents the transformation matrix between spherical coordinates and Cartesian coordinates.

[0166] Step S202 includes:

[0167] No. From one angle With a finite number energy level spherical harmonic vector They can be arranged as follows:

[0168]

[0169]

[0170] in, , Represents an imaginary number. The scalar spherical harmonic function is expressed as follows:

[0171]

[0172] in, This represents the Legendre function of the first type of correlation. This represents the order of the spherical harmonic function. The mode representing the spherical harmonic function, The zenith angle in spherical coordinates represents the direction of electromagnetic wave propagation. It represents the azimuth angle of the electromagnetic wave propagation direction in spherical coordinates.

[0173] According to the From one angle spherical harmonic vector Spherical harmonic matrices can be constructed Specifically, it can be expressed as:

[0174]

[0175] In the Cartesian coordinate system, the vector expansion coefficient matrix is ​​defined as:

[0176]

[0177] in, , and These represent the vector expansion coefficient matrices respectively. exist axis, shaft and The components of the axis.

[0178] Step S203 includes:

[0179] To solve for the optimal vector expansion coefficient matrix to reconstruct the antenna radiation field at arbitrary angles, a cost function is constructed using the least squares method:

[0180]

[0181] in, Let Frobenius norm be denoted. For the linear least squares problem described above, it can usually be solved by calculating the spherical harmonic matrix. The solution to this problem is obtained using the pseudo-inverse matrix. It can be represented as:

[0182]

[0183] in, Representing the spherical harmonic matrix The pseudo-inverse. Given the vector expansion coefficient matrix solution... Based on this, arbitrary zenith angles and azimuth angles can be calculated. The antenna radiation field under the given conditions is expressed as follows:

[0184]

[0185] Step S4 includes:

[0186] Step S401: Obtain the departure angle and arrival angle of each propagation ray in the global coordinate system, and determine the radial unit vector of the corresponding ray;

[0187] Step S402: Based on the actual rotation angle of the transceiver antenna in three-dimensional space, calculate the unit rotation vector of each coordinate axis of the antenna coordinate system relative to the global coordinate system;

[0188] Step S403: Using the unit rotation vector between the antenna coordinate system and the global coordinate system and the radial unit vector of the ray, solve for the zenith angle and azimuth angle of the ray in the antenna coordinate system. Based on the vector expansion coefficient and spherical harmonic matrix solved in step S2, calculate the antenna radiation field at the corresponding angle.

[0189] Step S401 includes:

[0190] Based on the propagation rays between the transceivers, the departure angle corresponding to each ray in the global coordinate system can be obtained. and angle of arrival The corresponding radial unit vector can be calculated based on the departure angle and arrival angle of each ray, expressed as follows:

[0191]

[0192]

[0193] in This represents the radial unit vector corresponding to the ray's departure angle in the global coordinate system. This represents the radial unit vector corresponding to the angle of arrival of a ray in the global coordinate system.

[0194] Step S402 includes:

[0195] In practice, antenna patterns are defined based on their own coordinate system. When the antenna rotates in three-dimensional space, the antenna coordinate system will no longer be aligned with the global coordinate system. Therefore, to calculate the antenna radiation field corresponding to the ray departure angle and arrival angle, the departure angle and arrival angle in the global coordinate system need to be transformed to the antenna coordinate system and matched with the antenna pattern. The antenna coordinate system is defined as follows:

[0196]

[0197]

[0198] in, In the antenna coordinate system Axial unit vector, Representing the antenna coordinate system Axial unit vector, , and These represent the antenna relative to the global coordinate system. axis, shaft and The rotation angle of the axis. Based on the unit vector of the antenna coordinate system, the ray departure angle and arrival angle in the antenna coordinate system can be calculated. :

[0199]

[0200] in, Represents the inverse cosine function. The radial unit vector representing the ray departure angle in the global coordinate system as described in step S401. Radial unit vector of the angle of arrival .

[0201] Step S403 includes:

[0202] Determine the zenith angle and azimuth angle in the antenna coordinate system. Then, the antenna radiation field at the corresponding angle can be calculated based on the vector expansion coefficients and spherical harmonic matrix described in step S2:

[0203]

[0204] Step S5 includes:

[0205] Step S501: Extract the incident angle at the interaction point between the ray and the multilayer anisotropic material and the electromagnetic parameters of the multilayer material. Construct a coordinate system rotation matrix according to the propagation direction of the incident wave, rotate and align the dielectric constant tensor and magnetic permeability tensor of the anisotropic material, and obtain the rotated electromagnetic characteristic tensor.

[0206] Step S502: Construct the coupling matrix of each layer of material based on the rotated electromagnetic parameters, and directly solve the eigenvalues ​​and eigenvectors of the corresponding coupling matrix according to the preset anisotropic material category;

[0207] Step S503: Using the continuity condition of each layer boundary, cascade the characteristic parameters of all material layers to obtain the generalized propagation coefficient matrix, and convert the generalized propagation coefficient matrix into a reflection coefficient matrix and a transmission coefficient matrix that match the ray tracing model.

[0208] Step S501 includes:

[0209] Based on the propagation ray between the transceiver, information on the interaction points between the ray and obstacles is obtained, including the coordinates of the interaction points and the material parameters of the interaction. The surface material of the obstacle is assumed to be a uniformly thick plate. The direction dependence of the electromagnetic parameters of the anisotropic material is considered, along with the dielectric constant of the anisotropic material. and permeability It can be represented in tensor form, as shown below:

[0210]

[0211] in, , and Defined respectively as anisotropic materials Optical axis Optical axis and Optical axis. Since anisotropic materials are treated as flat plates of uniform thickness, therefore... The axis is defined as the normal direction. The plane is defined as the tangential plane.

[0212] To reduce the complexity of subsequent matrix operations, coordinate rotation is used to reposition the optical axis of the anisotropic material. , and ) and the incident coordinate system ( , and Alignment, where the incident coordinate system is defined by the direction of propagation of the incident wave. Transmitted waves in anisotropic materials can be decomposed into two types, called Type I waves and Type II waves. In the incident coordinate system, the wave vector... Simplified to ,and , This represents the angle between the incident ray and the normal vector of the medium surface. Rotation matrix. Represented as:

[0213]

[0214] in, The z-axis represents the counterclockwise rotation angle in the coordinate system. The dielectric constant of the anisotropic material after rotation. and permeability It can be represented as:

[0215]

[0216]

[0217] in, and They represent and The Okay, number The parameters of the column.

[0218] Step S502 includes:

[0219] For each layer in a multilayer material, the dielectric constant of the anisotropic material after rotation of that layer is used. and permeability Directly construct the coupling matrix of this layer :

[0220]

[0221] in, It represents angular frequency.

[0222] Based on coupling matrix It can classify anisotropic materials and calculate the coupling matrix corresponding to each type of material. The eigenvalues ​​and eigenvectors. Based on the anisotropic characteristics of materials, the electromagnetic properties of anisotropic materials can be divided into five categories: isotropic materials, uniaxial materials (…). (axis is optical axis), uniaxial material ( (axis is optical axis), uniaxial material ( The electromagnetic properties of materials are categorized into uniaxial (optical axis) and biaxial materials. For isotropic materials, their electromagnetic properties remain consistent in all directions; in contrast, uniaxial materials are characterized by the existence of only a single optical axis, along which their electromagnetic properties differ from those in the other two orthogonal directions; and for biaxial materials, their electromagnetic properties are inconsistent in all three principal directions.

[0223] When the material is an isotropic material or a uniaxial material ( When the axis is the optical axis, the coupling matrix The eigenvalues ​​and eigenvectors can be represented as:

[0224]

[0225]

[0226] When anisotropic materials are When the axis is the optical axis, for uniaxial or biaxial materials, the coupling matrix eigenvalues for:

[0227] ;

[0228] Coupling matrix eigenvectors for:

[0229]

[0230]

[0231] in, Represents the coupling matrix The 1 eigenvector Represents the coupling matrix under the current material conditions. The 1 eigenvalue, The material is an isotropic material or The analytical solution for uniaxial materials with the optical axis as the axis only exists in Established at that time, Let z be the counterclockwise rotation angle of the optical axis in the coordinate system; when At that time, isotropic materials or The eigenvalues ​​and eigenvectors of a uniaxial material with the axis being the optical axis are derived from the characteristics of isotropic materials or The analytical solution for uniaxial materials is determined when the axis is the optical axis.

[0232] Step S503 includes:

[0233] To calculate the propagation coefficient of multilayer anisotropic materials, it is necessary to extend the electromagnetic wave propagation characteristics of single-layer anisotropic materials to those of multilayer anisotropic materials. For the th For anisotropic materials, the tangential components of the electric and magnetic fields can be expressed as:

[0234]

[0235] ,

[0236] in, Indicates the first The eigenvector matrix of the anisotropic material can be obtained by solving step S502. express The Middle Line 1 Column parameters, express The Middle The parameters of the row, The column vector order corresponds to the eigenvalues The order. To ensure the convergence of the solution, the first... The eigenvalues ​​of the layer are arranged according to the direction of electromagnetic wave propagation: and These correspond to the ascending waves of type I and type II waves, respectively. and These correspond to the descending waves of Type I and Type II waves, respectively. The ascending and descending waves are determined by the direction of the Poynting vector, which is defined as follows:

[0237]

[0238] in, Represents electric field, Indicates magnetic field, The eigenvalue represents the real part. When the Poynting vector corresponding to a certain eigenvalue is greater than 0, the eigenvalue represents an upward wave; when the Poynting vector corresponding to a certain eigenvalue is less than 0, the eigenvalue represents a downward wave.

[0239] Based on the material information obtained after the ray interacts with the obstacle, a multi-layer anisotropic material model is constructed for the current interaction point. Layer 0 and layer (n+1) represent air, and the region between these two layers represents n layers of anisotropic material. The number of layers n is determined by the actual number of material layers in the environment. The propagation equation between layer 0 and layer (n+1) is shown below:

[0240]

[0241] in,

[0242] ,

[0243] Indicates the first The thickness of the anisotropic material layer, , This represents the forward propagation matrix from layer 0 to layer (n+1). That is The Middle Line 1 Column corresponding parameters.

[0244] Since there are no reflected waves in the (n+1)th layer, that is... The forward propagation matrix can be decomposed, and the generalized reflection coefficient matrix can be calculated. and generalized transmission coefficient matrix Specifically, it can be expressed as:

[0245]

[0246]

[0247] in, and They represent from Type of incident wave to The generalized reflection coefficient and generalized transmission coefficient of the emitted wave. Generalized propagation coefficient matrix. Describes multilayer anisotropic materials and The relationship between them. However, the reflection coefficient and transmission coefficient in RT are based on the interface. and The relationship between them is defined. Therefore, the generalized propagation coefficient matrix must be converted into a reflection coefficient and transmission coefficient matrix usable by the RT. The reflection coefficient matrix usable by the RT and transmission coefficient matrix It can be represented as:

[0248]

[0249]

[0250] Step S6 includes:

[0251] Step S601: Distinguish between line-of-sight and non-line-of-sight paths based on the ray propagation path information, and calculate the channel propagation matrix for each ray for different propagation paths;

[0252] Step S602: Combine the continuous radiated electric field vector of the transceiver at the corresponding angle obtained in step S4, perform polarization matching on the channel propagation matrix, and calculate the received electric field of a single propagation ray;

[0253] Step S603: Coherently superimpose the received electric fields corresponding to all propagation paths and calculate the channel characteristics.

[0254] Step S601 includes:

[0255] Based on whether the ray interacts with the propagation environment, propagation paths can be divided into line-of-sight (LoS) paths and non-line-of-sight (non-LoS, NLoS) paths. For LoS paths, the propagation electric field reaching the receiver is calculated using the free-space propagation attenuation of the radiated electric field at the transmitter, and the channel propagation matrix... It is considered an identity matrix. In contrast, electric field calculations for non-line-of-sight paths involve a series of cascaded matrix transformations, including... Step reflection and Channel propagation matrix of first-order transmission It can be represented as:

[0256]

[0257]

[0258] in, and They represent the first The reflection coefficient matrix at the nth interaction point and the nth interaction point The transmission coefficient matrices at each interaction point are obtained from step S5. (Angle) Indicates the first Rotation angle between adjacent material coordinate systems at each interaction point. Channel propagation matrix. The order of the matrix cascades is determined by the interaction between the rays and the propagation environment.

[0259] Step S602 includes:

[0260] For each propagation path, at the known departure angle of the ray and angle of arrival In this case, the expression for the received electric field is:

[0261]

[0262] in, Indicates wavelength. and These represent the gains of the transmitting and receiving antennas, respectively. Indicates the transmission power. It is the receiving antenna impedance, reflecting the conjugate matching characteristics of the load impedance. and The antenna patterns representing the transmitting and receiving antennas, respectively, can be calculated in step S4. (Matrix) The channel propagation matrix is ​​calculated in step S601.

[0263] Step S603 includes:

[0264] The received electric fields corresponding to all propagation paths are coherently superimposed to calculate the total received power at the receiver. Based on the common polarization and cross-polarization received electric field components, the channel characteristics such as cross-polarization ratio, delay spread, and angle spread are calculated and output.

[0265] Example 1:

[0266] Step S1, Import as follows Figure 2 The indoor scene shown depicts the transceiver positions. The simulation frequency is set to 5 GHz. The transmitter uses a horn antenna equipped with a WR-187 rectangular waveguide, and the receiver uses a half-wave dipole antenna. The transmitter coordinates are (4.4, 1.4, 2.1) m, and the receiver positions range from (0.4, 2.4, 1.1) to (5.4, 2.4, 1.1) m, with a spacing of 0.02 m. The transmit power is set to 0 dBm. The main lobe direction of the horn antenna is set to... The transceiver rotates 45° counterclockwise around the main lobe. The propagation mechanism is set to 3rd-order reflection and 1st-order transmission. The spherical harmonic expansion orders of the transceiver are 15th and 6th, respectively. The materials of the multi-layered walls in the indoor scene are sequentially set along the positive x-axis as 0.1 m plasterboard, 0.2 m cement, and 0.1 m wood. The electromagnetic parameters of the isotropic materials are all from the International Telecommunication Union Recommendation P.2040. The simulation platform is a high-performance server with 512 GB RAM and an Intel Xeon Platinum 8369B CPU.

[0267] Step S2: Obtain the discrete antenna pattern based on the antenna configuration imported in Step S1, and reconstruct the continuous antenna radiation field through scalar spherical harmonic expansion.

[0268] Step S3: Run ray tracing simulation to extract the propagation path of the rays and the information on material interaction.

[0269] Step S4: Calculate the departure angle and arrival angle of the ray, and obtain the radiation electric field at the transmitting and receiving ends at the corresponding angles;

[0270] Step S5: Identify the ray-interacting material and calculate the propagation coefficient matrix based on the simplified generalized propagation matrix method. Since there are no multilayer anisotropic materials in this scenario, the electric field calculation is degenerated into the traditional multilayer isotropic material propagation coefficient calculation.

[0271] Step S6: Obtain the total received electric field at the receiver and calculate and output the channel characteristic information.

[0272] Figure 3A comparison chart showing the received power calculated by the theoretical solution, the conventional method, and the method of this invention is presented. The mean square error between the received power calculated by the method of this invention and the theoretical solution is approximately 2.91 dB, while the mean square error between the received power calculated by the conventional method and the theoretical solution is approximately 4.17 dB. The calculation results of the method of this invention are in high agreement with the theoretical solution, demonstrating high simulation accuracy.

[0273] Example 2:

[0274] Step S1, Import as follows Figure 2 The indoor scene shown has the transceiver positions set. The simulation frequency is set to 5 GHz. The transmitter uses a half-wave dipole antenna, and the receiver uses an ideal omnidirectional antenna. The transmitter coordinates are (4.4, 1.4, 2.1) m, and the receiver positions range from (0.4, 2.4, 1.1) to (5.4, 2.4, 1.1) m, with a spacing of 0.02 m. The transmit power is set to 0 dBm. The propagation mechanism is set to 3rd-order reflection and 1st-order transmission. The spherical harmonic expansion order of both the transceiver and receiver is 6th. The influence of multilayer anisotropic materials on depolarization is analyzed by considering isotropic and anisotropic materials in the scene: in the isotropic scene, the multilayer walls are composed of 0.02 m gypsum board; while in the anisotropic scene, the multilayer walls are successively set as 0.01 m gypsum board and 0.01 m biaxial material along the positive x-axis, with all other materials consistent with the isotropic scene. The electromagnetic parameters of isotropic materials are derived from the International Telecommunication Union Recommendation P.2040. The dielectric constant of biaxial materials is... Magnetic permeability is The simulation platform is a high-performance server with 512 GB of RAM and an Intel Xeon Platinum 8369B CPU.

[0275] Step S2: Obtain the discrete antenna pattern based on the antenna configuration imported in Step S1, and reconstruct the continuous antenna radiation field through scalar spherical harmonic expansion.

[0276] Step S3: Run ray tracing simulation to extract the propagation path of the rays and the information on material interaction.

[0277] Step S4: Calculate the departure angle and arrival angle of the ray, and obtain the radiation electric field at the transmitting and receiving ends at the corresponding angles;

[0278] Step S5: Identify the ray-interacting material and calculate the propagation coefficient matrix based on the simplified generalized propagation matrix method. Since there are no multilayer anisotropic materials in this scenario, the electric field calculation is degenerated into the traditional multilayer isotropic material propagation coefficient calculation.

[0279] Step S6: Obtain the total received electric field at the receiver and calculate and output the channel characteristic information.

[0280] Figure 4 The paper presents a comparison between the received power calculated by full-wave simulation and the received power calculated by the method of this invention in an anisotropic material scenario. Figure 4 In the diagram, (a) represents the received power calculated through full-wave simulation. Figure 4 (b) represents the received power calculated using the method of this invention. The mean square error between the calculation result of the method of this invention and the calculation result of the full-wave simulation is approximately 3.52 dB. Furthermore, the calculation time of the method of this invention is approximately 467.17 min, while the calculation time of the full-wave simulation is approximately 2370.22 min, representing a decrease in calculation time of 80.29%. This demonstrates that the method of this invention has good computational efficiency while maintaining high computational accuracy.

[0281] Figure 5 The method of this invention is shown to compare the channel cross-polarization ratio calculated in isotropic and anisotropic scenarios, so as to characterize the influence of multilayer anisotropic materials on the channel depolarization effect. Figure 5 In the example (a), the scene is isotropic. Figure 5 (b) shows a comparison of the channel cross-polarization ratio calculated by the method of this invention in an anisotropic scenario. When there is a Loss path between the transceivers, the difference between the two scenarios is small: the cross-polarization ratio in the anisotropic scenario is only 0.88 dB lower than that in the isotropic scenario. This is because when there is a Loss path between the transceivers, the received power is mainly dominated by the Loss path, and the NLoS path contributes less power than the Loss path. However, when there is no Loss path between the transceivers, the cross-polarization ratio in the anisotropic scenario is 2.76 dB lower than that in the isotropic scenario, reflecting that multilayer anisotropic materials lead to a more significant depolarization effect. This decrease is mainly due to the polarization conversion caused by multilayer anisotropic materials, which in turn enhances the coupling of the electric field co-polarization component to the cross-polarization component.

[0282] The above description is merely a specific embodiment of the present invention, but the scope of protection of the present invention is not limited thereto. Any variations or substitutions that can be easily conceived by those skilled in the art within the scope of the technology disclosed in the present invention should be included within the scope of protection of the present invention.

Claims

1. A ray tracing channel modeling method considering antenna and material depolarization effects, characterized in that, include: Step S1: Import the simulation scenario and set the transceiver antenna configuration and simulation parameters. Step S2: Based on the antenna configuration and simulation parameters of the transceiver end, obtain the radiation pattern of the discrete antenna at the transceiver end, and reconstruct the radiation field of the continuous antenna at the transceiver end through scalar spherical harmonic expansion; Step S3: Run the ray tracing channel simulation according to the simulation scene, the location of the transceiver and the simulation parameters, and extract the propagation path of each ray and the material interaction information between the ray and the scene; Step S4: Based on the propagation path of the rays, calculate the departure angle at the transmitting end and the arrival angle at the receiving end for each ray, and obtain the radiated electric field at the transmitting end and the radiated electric field at the receiving end based on the continuous antenna radiation field at the transmitting and receiving ends according to the departure angle and the arrival angle. Step S5: Identify the ray interaction materials corresponding to each ray based on the material interaction information, and calculate the electromagnetic propagation characteristics of the material based on the generalized propagation matrix method for the ray interaction materials to obtain the reflection coefficient matrix and transmission coefficient matrix that match the ray tracing model. Step S6: Calculate the receiving electric field corresponding to each ray based on the radiated electric field at the transmitting end and the radiated electric field at the receiving end, as well as the reflection coefficient matrix and the transmission coefficient matrix. Coherently superimpose the receiving electric fields corresponding to all rays to obtain the total receiving electric field at the receiving end. Calculate and output the channel characteristic information based on the ray propagation path and the total receiving electric field at the receiving end.

2. The ray tracing channel modeling method considering antenna and material depolarization effects according to claim 1, characterized in that, Simulation parameters include transmitter and receiver positions, simulation frequency, propagation mechanism, material parameters, antenna attitude, and spherical harmonic expansion order.

3. The ray tracing channel modeling method considering antenna and material depolarization effects according to claim 1, characterized in that, Step S2 includes: Step S201: Obtain the sampling data of the discrete antenna pattern of the transceiver end, extract the zenith angle component and azimuth angle component of the electric field in the spherical coordinate system for each sampling direction, and construct the coordinate transformation relationship according to the zenith angle and azimuth angle corresponding to each sampling direction to convert the zenith angle component and azimuth angle component of the electric field into the electric field sampling vector in the Cartesian coordinate system. Step S202: Based on the sampling directions corresponding to the discrete antenna pattern sampling data and the maximum order of the spherical harmonic expansion, construct the spherical harmonic basis function vector corresponding to each sampling direction, and construct the spherical harmonic matrix based on the spherical harmonic basis function vector; Step S203: Based on the electric field sampling vector and spherical harmonic matrix in the Cartesian coordinate system, establish a linear least squares problem, and solve the vector expansion coefficient matrix by calculating the pseudo-inverse of the spherical harmonic matrix, thereby reconstructing the continuous antenna radiation field in any direction.

4. The ray tracing channel modeling method considering antenna and material depolarization effects according to claim 3, characterized in that, Step S201 includes: Let the discrete antenna pattern include The sampling direction, where the first... Each sampling direction is represented as ,and , Indicates the first The zenith angle of each sampling direction in spherical coordinates Indicates the first The azimuth angle of each sampling direction in spherical coordinates; No. The electric field sampling data in spherical coordinates corresponding to each sampling direction is represented as follows: ； in, Indicates the first The electric field sampling vector in spherical coordinates corresponding to each sampling direction Indicates the first The zenith angle component of the electric field corresponding to each sampling direction Indicates the first The electric field azimuth component corresponding to each sampling direction; Will The electric field sampling data in spherical coordinates corresponding to each sampling direction are arranged as follows: ； in, Indicates by The electric field sampling matrix is ​​composed of the electric field sampling vectors in the spherical coordinate system corresponding to each sampling direction; Will The electric field sampling data corresponding to each sampling direction in the spherical coordinate system is transformed to the Cartesian coordinate system, and the specific expression is as follows: ； in, Indicates by The electric field sampling matrix is ​​composed of the electric field sampling vectors in the Cartesian coordinate system corresponding to each sampling direction. This represents the transformation matrix between spherical coordinates and Cartesian coordinates.

5. A ray tracing channel modeling method considering antenna and material depolarization effects according to claim 3, characterized in that, Step S202 includes: Let the maximum order of the pre-defined spherical harmonic expansion be... The total number of basis functions contained in the spherical harmonic expansion For the first There are spherical harmonic basis functions, whose corresponding spherical harmonic order and mode are respectively expressed as: and ,in , , ; according to Constructing the spherical harmonic matrix from the spherical harmonic basis function vectors corresponding to each sampling direction. spherical harmonic matrix Represented as: ； in, Representing the spherical harmonic matrix The The sampling direction of the first A spherical harmonic basis function For the reason Arrangement OK Column matrix; Represented as: ； in, Represents the imaginary unit. Indicates the first Each sampling direction Next Rank, number scalar spherical harmonic functions of the mode; The expression is: ； in, Let e ​​represent the Legendre function of the first kind of correlation, where e is the natural base; Define the vector expansion coefficient matrix in the Cartesian coordinate system. The vector expansion coefficient matrix is ​​represented as: ； in, , , and These represent the x-axis, y-axis, and z-axis in the Cartesian coordinate system, respectively. In Cartesian coordinate system Directional electric field component Vector expansion coefficients on spherical harmonic basis functions For the reason Arrangement OK Column matrix.

6. The ray tracing channel modeling method considering antenna and material depolarization effects according to claim 3, characterized in that, Step S203 includes: Based on the electric field sampling matrix in the Cartesian coordinate system obtained in step S201 The spherical harmonic matrix obtained in step S202 The linear least squares problem is established, specifically represented as: ； in, Represents the coefficient matrix of the vector expansion The least squares solution, Representing the Frobenius norm; by calculating the spherical harmonic matrix. The pseudo-inverse matrix is ​​used to solve for it. Represented as: ； in, Representing the spherical harmonic matrix The pseudo-inverse matrix; For any direction to be reconstructed , This represents the zenith angle of the direction to be reconstructed in spherical coordinates. The azimuth angle of the direction to be reconstructed in spherical coordinates is defined as the direction to be reconstructed. The corresponding column vector of spherical harmonic basis functions The column vector of spherical harmonic basis functions is represented as: ； in, Represents the column vector of spherical harmonic basis functions The One spherical harmonic basis function; Represented as: ； in, Indicates the first The order of the spherical harmonics corresponding to each spherical harmonic basis function. Indicates the first The spherical harmonic function modes corresponding to each spherical harmonic basis function. , , , Indicates any direction to be reconstructed Next Rank, number scalar spherical harmonic functions of the mode; The least squares solution based on the vector expansion coefficient matrix Column vector of spherical harmonic basis functions The continuous antenna radiation field vector in Cartesian coordinates under the direction to be reconstructed is obtained. Specifically, it is expressed as: 。 7. A ray tracing channel modeling method considering antenna and material depolarization effects according to claim 3, characterized in that, Step S4 includes: Step S401: Based on each ray propagation path, obtain the departure angle at the transmitting end and the arrival angle at the receiving end for each propagating ray in the global coordinate system, and determine the radial unit vector at the transmitting end corresponding to the departure angle and the radial unit vector at the receiving end corresponding to the arrival angle. Step S402: Based on the transmit and receive antenna attitude determined in step S1 and the actual rotation angle of the transmit and receive antenna in three-dimensional space, calculate the unit rotation vector of each coordinate axis of the transmit antenna coordinate system and the receive antenna coordinate system relative to the global coordinate system. Step S403: Using the unit rotation vector obtained in step S402, transform the radial unit vector leaving the transmitter and the radial unit vector arriving at the receiver obtained in step S401 to the transmitter antenna coordinate system and the receiver antenna coordinate system, respectively. Solve for the departure zenith angle and departure azimuth angle in the transmitter antenna coordinate system, and the arrival zenith angle and arrival azimuth angle in the receiver antenna coordinate system. Based on the vector expansion coefficient matrix obtained in step S2, calculate the continuous antenna radiation field of the transmitter corresponding to the departure zenith angle and departure azimuth angle, and the radiation electric field of the transmitter and the radiation electric field of the receiver corresponding to the arrival zenith angle and arrival azimuth angle. Step S401 includes: Based on the propagation rays between the transceivers, obtain the departure angle of each ray in the global coordinate system. and angle of arrival ,in and These represent the zenith angle and azimuth angle of the propagating ray in the global coordinate system at the transmitting end, respectively. and These represent the zenith angle and azimuth angle of the global coordinate system at the receiving end, respectively; The radial unit vector corresponding to each ray is calculated based on its departure angle and arrival angle, expressed as follows: ； ； in This represents the radial unit vector corresponding to the ray's departure angle in the global coordinate system. This represents the radial unit vector corresponding to the angle of arrival of a ray in the global coordinate system. Step S402 includes: The departure angle and arrival angle in the global coordinate system are transformed to the antenna coordinate system and matched with the antenna pattern. The antenna coordinate system is defined as follows: ； ； in, In the antenna coordinate system Axial unit vector, Representing the antenna coordinate system Axial unit vector, , and These represent the antenna relative to the global coordinate system. axis, shaft and Rotation angle of the axis; calculate the departure angle and arrival angle of each ray in the antenna coordinate system based on the unit vector of the antenna coordinate system. : ； in, and These represent the zenith angle and azimuth angle of the propagating ray in the global coordinate system at the transmitting end, respectively. Represents the inverse cosine function. The radial unit vector representing the ray departure angle in the global coordinate system in step S401. Or the radial unit vector of the angle of arrival in the global coordinate system The superscript T indicates transpose; Step S403 includes: Determine the departure angle and arrival angle of each ray in the antenna coordinate system. Then, the antenna radiation field is calculated based on the vector expansion coefficient matrix and spherical harmonic matrix from step S2: ； in, This represents the antenna radiation field vector corresponding to the departure angle and arrival angle of each ray in the antenna coordinate system. This represents the least-squares solution of the vector expansion coefficient matrix obtained in step S2. This represents a column vector consisting of the spherical harmonic basis function values ​​corresponding to the departure angle and arrival angle of each ray in the antenna coordinate system.

8. A ray tracing channel modeling method considering antenna and material depolarization effects according to claim 3, characterized in that, Step S5 includes: Step S501: Extract the incident angle at the interaction point between the ray and the multilayer anisotropic material and the electromagnetic parameters of the multilayer material. Construct a coordinate system rotation matrix according to the direction of incident wave propagation, rotate and align the dielectric constant and permeability of the anisotropic material, and obtain the dielectric constant and permeability of the rotated anisotropic material. Step S502: Construct the coupling matrix of each layer of anisotropic material based on the dielectric constant and magnetic permeability of the rotated anisotropic material, and directly solve the eigenvalues ​​and eigenvectors of the corresponding coupling matrix according to the preset anisotropic material category. Step S503: Using the boundary continuity conditions of the anisotropic materials in each layer, the eigenvalues ​​and eigenvectors of the anisotropic materials in each layer are combined into an eigenmatrix and concatenated to obtain the generalized propagation coefficient matrix. The generalized propagation coefficient matrix is ​​then converted into a reflection coefficient matrix and a transmission coefficient matrix that match the ray tracing channel model.

9. A ray tracing channel modeling method considering antenna and material depolarization effects according to claim 8, characterized in that, Step S501 includes: Based on the propagation ray between the transceiver, information on the interaction points between the ray and obstacles is obtained, including the coordinates of the interaction points and the material parameters of the interaction. The surface material of the obstacle is assumed to be a uniformly thick plate. The direction dependence of the electromagnetic parameters of the anisotropic material is considered, along with the dielectric constant of the anisotropic material. and permeability Represented in tensor form, as shown below: ； in, , and Defined respectively as anisotropic materials Optical axis Optical axis and Optical axis The dielectric constant of anisotropic materials is represented by... Components along the optical axis The dielectric constant of anisotropic materials is represented by... Components along the optical axis The dielectric constant of anisotropic materials is represented by... Components along the optical axis This indicates the electrical conductivity of anisotropic materials at... Components along the optical axis This indicates the electrical conductivity of anisotropic materials at... Components along the optical axis This indicates the electrical conductivity of anisotropic materials at... Components along the optical axis; Coordinate rotation is used to align the optical axis of the anisotropic material with the incident coordinate system, where the incident coordinate system is defined by the propagation direction of the incident wave. The transmitted wave in the anisotropic material is decomposed into two types, referred to as Type I and Type II waves. In the incident coordinate system, the wave vector... Simplified to ,in Indicates the wave vector in Components along the optical axis Indicates the wave vector in The components along the optical axis, and , Represents the angle between the incident ray and the normal vector of the medium surface; rotation matrix Represented as: ； in, This represents the counterclockwise rotation angle of the z-axis in the coordinate system. Dielectric constant of the rotated anisotropic material and permeability Represented as: ； ； in, and They represent and The line, number Column parameters, , ; Step S502 includes: For each layer in a multilayer anisotropic material, the dielectric constant of the anisotropic material after rotation of that layer is used. and permeability Directly construct the coupling matrix of this layer : ； in, Indicates angular frequency; Based on coupling matrix Anisotropic materials are classified, and the coupling matrix corresponding to each type of anisotropic material is calculated. The eigenvalues ​​and eigenvectors; based on the anisotropic characteristics of materials, the electromagnetic properties of anisotropic materials are divided into five categories: isotropic materials, The axis is a uniaxial material under the condition of optical axis. The axis is a uniaxial material under the condition of optical axis. The axis is uniaxial and biaxial materials under optical axis conditions; When anisotropic materials are isotropic materials or When the axis is the optical axis and the material is uniaxial, the coupling matrix is... The eigenvalues ​​and eigenvectors are represented as follows: ； ； in, Represents the coupling matrix under the current material conditions. The first and second eigenvalues; Represents the coupling matrix under the current material conditions. The third and fourth eigenvalues; Indicates that it is made of isotropic materials or The eigenvector matrix of a uniaxial material under the condition that the axis is the optical axis; When anisotropic materials are When the axis is the optical axis, for uniaxial or biaxial materials, the coupling matrix eigenvalues for: ； Coupling matrix eigenvectors for: ； ； in, Represents the coupling matrix The 1 eigenvector Represents the coupling matrix under the current material conditions. The 1 eigenvalue, The material is an isotropic material or The analytical solution for uniaxial materials with the optical axis as the axis only exists in Established at that time, Let z be the counterclockwise rotation angle of the optical axis in the coordinate system; when At that time, isotropic materials or The eigenvalues ​​and eigenvectors of a uniaxial material with the axis being the optical axis are derived from the characteristics of isotropic materials or The analytical solution for uniaxial materials is determined when the axis is the optical axis. Step S503 includes: The electromagnetic wave propagation characteristics of single-layer anisotropic materials are extended to those of multilayer anisotropic materials; for the first... The tangential component of the electromagnetic field in a layered anisotropic material. Represented as: ； , ； in, This represents the total thickness of a multilayer anisotropic material. Indicates the first The distance between the anisotropic material layer and the incident plane. Indicates the first Tangential components of the electromagnetic field in anisotropic materials. Indicates the first The eigenvector matrix of the anisotropic material is obtained by solving step S502. express The Middle Line 1 Column parameters, , , This indicates the position of each eigenwave in anisotropic materials. The magnitude column vector at that point, express The Middle The parameters of the row, , The column vector order corresponds to the eigenvalues The order in which They represent the first The first, second, third, and fourth eigenvalues ​​of the layer; The eigenvalues ​​of the layer are arranged according to the direction of electromagnetic wave propagation: and These correspond to the ascending waves of type I and type II waves, respectively. and These correspond to the descending waves of Type I and Type II waves, respectively; the ascending and descending waves are determined by the direction of the Poynting vector. The definition is as follows: ； in, Represents electric field, Indicates magnetic field, This represents the conjugate transpose of the magnetic field. The eigenvalue represents the real part; when the Poynting vector corresponding to a certain eigenvalue is greater than 0, the eigenvalue represents an upward wave; when the Poynting vector corresponding to a certain eigenvalue is less than 0, the eigenvalue represents a downward wave. Based on the material information obtained after the ray interacts with the obstacle, a multi-layer anisotropic material model is constructed for the current interaction point. Layer 0 and layer (n+1) represent air, and the region between these two layers represents n layers of anisotropic material. The number of layers n is determined by the actual number of material layers in the environment. The propagation equation between layer 0 and layer (n+1) is shown below: ； in, ， This indicates the position of the anisotropic material in layer 0. The eigenwave amplitude column vector at that location, , Indicates the first Layered anisotropic materials and the first The transformation matrix between layers of anisotropic materials. , Indicates the first The propagation matrix within the anisotropic material. Indicates the first Layered anisotropic materials and the first The transformation matrix between layers of anisotropic materials. This indicates the position of each eigenwave in anisotropic materials. The magnitude column vector at that point, Indicates the first The thickness of the anisotropic material layer, This represents the forward propagation matrix of the anisotropic material from layer 0 to layer (n+1). That is The Middle Line 1 The corresponding parameters are listed below. , , Indicates the first The eigenvector matrix of anisotropic layered materials, Indicates the first The eigenvector matrix of anisotropic layered materials; Since there are no reflected waves in the (n+1)th layer, that is... , and They represent the first The amplitude coefficients of the first and second eigenwaves in the layer; and They represent the first The eigenvalues ​​corresponding to the first and second eigenwaves in the layer are determined; the forward propagation matrix is ​​decomposed, and the generalized reflection coefficient matrix is ​​calculated. and generalized transmission coefficient matrix : ； ； in, and They represent the first Line 1 Generalized reflection and transmission coefficients of the array; converting the generalized propagation coefficient matrix into a usable reflection and transmission coefficient matrix for RT; usable reflection coefficient matrix for RT. and transmission coefficient matrix Represented as: ； ； , , and These represent the generalized reflection coefficient matrix. The elements in the first row and first column, the first row and second column, the second row and first column, and the second row and second column; , , and These represent the generalized reflection coefficient matrix. The elements in the first row and first column, the first row and second column, the second row and first column, and the second row and second column; , , , , and These represent the eigenvector matrices corresponding to the 0th layer of the medium. The element at the corresponding row and column position; and They represent the first eigenvector matrix corresponding to the layer medium The elements at the corresponding row and column positions; where the 0th layer is an air layer and located on the incident side of the multilayer anisotropic material, the 1st layer... The medium is an air layer and is located on the transmission side of a multilayer anisotropic material.

10. A ray tracing channel modeling method considering antenna and material depolarization effects according to claim 3, characterized in that, Step S6 includes: Step S601: Distinguish between line-of-sight and non-line-of-sight paths based on the ray propagation path information, and calculate the channel propagation matrix for each ray for different propagation paths; Step S602: Combine the continuous radiated electric field vectors of the transceiver at the corresponding angle obtained in step S4, perform polarization matching on the channel propagation matrix, and calculate the received electric field of a single propagation ray. Step S603: Coherently superimpose the received electric fields corresponding to all propagation paths and calculate the channel characteristic information; Step S601 includes: Based on whether the ray interacts with the propagation environment, the propagation path is divided into line-of-sight (LoS) paths and non-line-of-sight (NLoS) paths. For LoS paths, the propagation electric field reaching the receiver is calculated using the free-space propagation attenuation of the radiated electric field at the transmitter, and the channel propagation matrix is ​​used. It is considered an identity matrix; in contrast, electric field calculations for non-line-of-sight paths involve a series of cascaded matrix transformations, including... Step reflection and Channel propagation matrix of first-order transmission Represented as: ； in, and They represent the first The reflection coefficient matrix at the nth reflection point and the nth reflection coefficient matrix The transmission coefficient matrices at each transmission point are obtained from step S5; angle Indicates the first The rotation angle between adjacent material coordinate systems at each reflection point, angle Indicates the first Rotation angle between adjacent material coordinate systems at each point of transmission action; channel propagation matrix The order of the matrix cascades is determined by the interaction between the rays and the propagation environment; Step S602 includes: For each propagation path, at the known departure angle of the ray and angle of arrival Under these conditions, the received electric field The expression is: ； in, and These represent the zenith angle and azimuth angle in the antenna coordinate system corresponding to the propagation path at the transmitting end, respectively. and These represent the zenith angle and azimuth angle in the antenna coordinate system corresponding to the propagation path at the receiving end, respectively. Indicates wavelength. and These represent the gains of the transmitting and receiving antennas, respectively. Indicates the transmission power. It is the receiving antenna impedance, reflecting the conjugate matching characteristics of the load impedance; and The radiated electric fields of the transmitting antenna and the receiving antenna at the departure angle and the arrival angle of the ray, respectively, are calculated by step S4; Step S603 includes: The received electric fields corresponding to all propagation paths are coherently superimposed to calculate the total received power at the receiver. Based on the common polarization and cross-polarization received electric field components, the channel's cross-polarization ratio, delay spread, and angle spread channel characteristic information are calculated and output.

Images