Turbulence flux loss correction calculation method and system for improved vortex correlation observation system
Through the sensor geometric correction and path averaging algorithm of the non-orthogonal three-dimensional ultrasonic acemeter, combined with analytical solution and lookup table linear interpolation, the problem of insufficient correction accuracy of the turbulent energy loss of the non-orthogonal three-dimensional ultrasonic acemeter is solved, and flow field calculation with higher accuracy and speed is achieved.
Patent Information
- Application Number
- CN202510645969.2
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2025-05-20
- Publication Date
- 2025-08-15
AI Technical Summary
In the prior art, when the path averages of non-orthogonal three-dimensional ultrasonic acemeters, the turbulent energy loss correction calculation accuracy is insufficient and the speed is insufficient, so they cannot effectively take into account the accuracy and speed of the actual flow field calculation.
The sensor geometric shape correction method of a non-orthogonal three-dimensional ultrasonic acemeter is used, combined with the path averaging algorithm, the wind speed spectrum and temperature spectrum losses in the synthetic wind direction, side wind direction and vertical wind direction are corrected, and the turbulent energy loss is corrected by analytical solution and look-up table linear interpolation.
The turbulence energy loss correction accuracy of the non-orthogonal three-dimensional ultrasonic acemeter is improved, ensuring the accuracy of flow field calculation, and at the same time improving the calculation speed.
Smart Images

Figure CN120490529A_ABST
Abstract
Description
Technical Field
[0001] The present invention belongs to the technical field of turbulent flux correction, and in particular relates to an improved turbulent flux loss correction calculation method and system for an eddy correlation observation system. Background Art
[0002] The eddy covariance method uses high-frequency sampling (10Hz-100Hz) to collect the pulsating quantities of atmospheric variables, thereby obtaining the flux of turbulent transport of momentum, energy, and matter between the Earth's surface and the atmosphere. After more than 70 years of development, it is now widely used and the preferred high-precision direct observation method. Related observational research has achieved unprecedented development, and the international flux network (FLUXNET) focusing on ecosystem exchange has been established globally.
[0003] The eddy-correlation observation system mainly includes three-dimensional ultrasonic wind and temperature meters, open-circuit or closed-circuit water vapor / carbon dioxide analyzers or other trace gas analyzers, and related auxiliary observations. The high-frequency pulsation data of atmospheric variables observed by the eddy system also need to be combined with the instrument parameter settings and related records of the eddy system (such as the instrument firmware number or whether sensor blockage correction has been performed, instrument orientation, spacing between different instruments, observation altitude, dynamic records of surrounding underlying surface topography and vegetation cover, maintenance logs, etc.), and undergo outlier elimination, some necessary turbulent energy loss corrections (such as sensor blockage correction, tilt correction, frequency loss correction, ultrasonic virtual temperature correction, air density pulsation correction, etc.), quality evaluation and quality control (including stability and turbulence development sufficiency detection, footprint inspection, etc.). Among them, turbulent energy loss correction refers to the calculation of turbulent fluxes over a period of time, such as 30 minutes, based on the original observations of the vortex system. From the frequency characteristics, the statistical time is limited, which may result in the loss of low-frequency energy. In the high-frequency range, the sampling frequency is limited. The response time, acoustic path, or optical path average of the instrument makes it impossible to identify high-frequency vortices smaller than the path. The separation between the gas analyzer and the three-dimensional ultrasonic wind temperature meter along the transverse direction of the wind direction leads to phase differences and high-frequency loss. The active observation of closed-circuit observation leads to attenuation of high-frequency signals, etc., which may cause the observed turbulent flux to be less than the actual turbulent flux. Therefore, various corrections are made for this purpose. With the development and improvement of vortex instruments and related algorithms, the uncertainty of observations is constantly being recognized. In the existing technology, there are many software programs for post-processing turbulent fluxes.
[0004] Traditional correction algorithms for turbulent fluxes, including sensor obstruction correction and path averaging correction, are derived based on orthogonal three-dimensional ultrasonic anemometers. Currently, the industry mostly uses non-orthogonal three-dimensional ultrasonic anemometers to measure high-frequency pulsations of atmospheric wind speed and temperature. Due to the complexity of non-orthogonal ultrasonic anemometers, algorithms derived from orthogonal three-dimensional ultrasonic anemometers are generally used.
[0005] The traditional orthogonal three-dimensional ultrasonic anemometer is based on three pairs of independent sensors for observation. There is a pair of sensors in the vertical direction to observe the vertical wind speed u z , there are two pairs of sensors on the horizontal plane observing u x and u y In a non-orthogonal 3D ultrasonic anemometer, all three pairs of sensors are at a fixed angle to the vertical axis, and their projection angles on the horizontal plane are also fixed (the geometric shape of the three pairs of sensors in a 3D ultrasonic anemometer, such as IRGASON). The three pairs of sensors jointly determine the 3D wind speed.
[0006] Assuming that the acoustic path length of each pair of sensors is l, the vector directions of the three pairs of probes in the xyz coordinates of the ultrasound system are:
[0007]
[0008] Where, is the horizontal angle between the three pairs of sensors and the positive half axis of the x-axis, and is related to the type and model of the three-dimensional ultrasonic anemometer. For example, in IRGASON, The angle θ between the ultrasonic probe and the vertical direction (z axis) is 30°.
[0009] In a non-orthogonal three-dimensional ultrasonic anemometer, the three pairs of sensors do not directly detect the vertical velocity u z , what is actually detected is the velocity along the sound path direction of the three pairs of sensors (u A ,u B ,u C ), the wind speed and the output three-dimensional wind (u x ,u y ,u z ) needs to be transformed through a matrix.
[0010] The technical defects of the existing technology are: in the existing technology, the non-orthogonal three-dimensional ultrasonic anemometer generally uses the one-dimensional path averaging algorithm in the orthogonal three-dimensional ultrasonic anemometer to correct the turbulent energy loss caused by path averaging in the x1, x2 and x3 directions, and assumes that the x1 axis is parallel to the horizontal plane. This cannot guarantee the accuracy of the complex path averaging calculated by the non-orthogonal three-dimensional ultrasonic anemometer when the actual synthetic wind speed is non-horizontal, and cannot simultaneously take into account the calculation speed. Summary of the Invention
[0011] In order to overcome the problems existing in the related art, the disclosed embodiments of the present invention provide an improved turbulent flux loss correction calculation method and system for an eddy correlation observation system.
[0012] The technical solution is as follows: improving the turbulent flux loss correction calculation method of the eddy covariance observation system, including:
[0013] S1, based on the instantaneous wind speed and direction in the original turbulence observation, the flow field distortion effect is corrected by the sensor geometry of the non-orthogonal three-dimensional ultrasonic anemometer;
[0014] S2, using the path averaging algorithm of the non-orthogonal three-dimensional ultrasonic anemometer sensor, corrects the wind speed spectrum loss caused by path averaging in the synthetic wind direction x1, crosswind direction x2, and vertical wind direction x3 averaged over a period of time, such as 30 minutes;
[0015] S3, using the non-orthogonal three-dimensional ultrasonic wind temperature sensor path averaging correction algorithm to correct the ultrasonic temperature spectrum loss caused by path averaging;
[0016] S4, the numerical integration of wind speed spectrum loss correction and ultrasonic temperature spectrum loss correction at different dimensionless frequencies is made into a lookup table. In the frequency loss correction calculation of the turbulence spectrum, a method combining the analytical solution with the linear interpolation of the above lookup table is used to obtain the turbulent energy loss correction information caused by path averaging in the synthetic wind direction x1, crosswind direction x2, and vertical wind direction x3.
[0017] In step S1, the effect of the sensor geometry of the non-orthogonal 3D ultrasonic anemometer on flow field distortion is corrected, including:
[0018] Before path averaging correction, the effect of the sensor wake on the flow field distortion on the sensor acoustic path is calculated. The affected wind speed u detected by the first pair of sensors is Am The angle θ between the actual wind direction and the path of the first pair of sensors A The expression is:
[0019] u Am =u A (0.84+0.16sinθ A )
[0020] or
[0021]
[0022] Where u Am is the wind speed detected by the first pair of sensors after the impact, θ A is the angle between the actual wind direction and the path of the first pair of sensors, u A It is the wind speed projected from the actual wind speed to the direction of the first pair of sensors. The algorithms for the other two pairs of sensors are the same as those for the first pair of sensors.
[0023] After obtaining the original wind speed of the three-dimensional ultrasonic anemometer, when the flow field distortion correction is required, the three-dimensional wind speed (u x ,u y ,u z ) Refer to the following formula for conversion:
[0024]
[0025] Where u B ,u C are the wind speeds projected onto the paths of the second and third pairs of sensors, respectively. x ,u y ,u z are the wind speeds along the x-axis, y-axis, and z-axis of the three-dimensional coordinate system in the ultrasonic system, respectively; A is the wind field conversion matrix;
[0026] Converted to wind speed along the wind direction of three pairs of sensors (u A ,u B ,u C ), after correcting the flow field distortion, it is converted into a new (u x ,u y ,u z ):
[0027]
[0028] Where A -1 That is the inverse transformation matrix of the wind field.
[0029] In step S2, the non-orthogonal three-dimensional ultrasonic anemometer sensor wind speed path averaging algorithm includes:
[0030] Assume that the actual horizontal wind direction angle averaged over a period of time (such as 30 minutes) in the ultrasonic system is α. Its three-dimensional synthetic wind speed may not be located on the horizontal plane of the ultrasonic system, and there is an angle γ with the horizontal plane of the ultrasonic system. Here, the average vertical speed in the ultrasonic system is taken as the positive angle γ, and a new coordinate system is defined as the synthetic wind direction x1, side wind direction x2, and vertical wind direction x3 coordinate system. The actual synthetic wind speed u1 is along the synthetic wind direction x1 axis, u2 is along the side wind direction x2 axis, parallel to the actual horizontal plane, and u3 is along the vertical wind direction x3 axis.
[0031] If the ultrasound system is installed vertically on the ground and does not tilt over time, the xy coordinate system of the ultrasound system is parallel to the actual horizontal plane, then the x2 axis is the result of rotating the y axis of the ultrasound system counterclockwise about an angle α on its horizontal plane. If the ultrasonic system is tilted, for example, in the ultrasonic system, the ultrasonic instrument is at the γ0 position of the xy coordinate of the ultrasonic system, and the z axis is tilted at an angle of θ0 with the actual vertical direction. Then, after the above rotation, it can be rotated counterclockwise along the synthetic wind direction x1 axis. Angle so that the generated x2 axis is parallel to the actual horizontal plane.
[0032] (u1,u2,u3) and (u x ,uy ,u z ) exists
[0033]
[0034] T represents the transpose of the matrix. Therefore, in the (x1, x2, x3) coordinate system, the sound path vectors along the three pairs of sensor directions become:
[0035]
[0036] Where, l a ,l b ,l c are the sound path vectors of the first, second and third pairs of sensors in the (x1, x2, x3) coordinate system, respectively, and l is the sound path length (constant);
[0037] Without considering the influence of path averaging, the wind speed detected by the three pairs of ultrasonic sensors is (u A ,u B ,u C ), the transformation between this wind speed and the three-dimensional wind speed (u1, u2, u3) in the (x1, x2, x3) coordinate system satisfies the following matrix:
[0038]
[0039] Let matrix M = ABDE, N = E T D T B T A -1 After considering the influence of the path (or sound path) average, the wind speed (u A ,u B ,u C ), which is actually:
[0040]
[0041] Where M ij is the value of the i-th row and j-th column of matrix M, are the wind speed u1 along the direction of the first pair of sensors, the direction of the second pair of sensors, and the direction of the third pair of sensors after considering the average effect of their respective sound paths; are the wind speed u2 along the first pair of sensors, the second pair of sensors, and the third pair of sensors, respectively, after considering the average effect of their respective sound paths. are the wind speed u3 along the first pair of sensors, the second pair of sensors, and the third pair of sensors, respectively, after considering the average effect of their respective acoustic paths, and ~ is the path average;
[0042]
[0043] Where x0 is the common center position of the three pairs of sensors in the 3D ultrasound system;
[0044] In the (x1, x2, x3) coordinate system, the actual three-dimensional wind speed (u1 m ,u2 m ,u3 m ), the expression is:
[0045]
[0046] After substitution, we have:
[0047]
[0048] Considering that the observed wind speed is the contribution of eddies of different wave numbers, the Fourier transform yields:
[0049]
[0050] In the above three formulas, e is the exponent and i represents the imaginary number. is the three-dimensional wave number, and the corresponding wave numbers in the (x1, x2, x3) directions are (k1, k2, k3). Accordingly, the wind speed along the direction of the synthetic wind speed after path averaging is as follows: Three-dimensional wind speed components transformed in the frequency domain satisfy
[0051]
[0052] Where i represents an imaginary number; according to the above formula, we have:
[0053]
[0054] Where d is the component;
[0055]
[0056] j=1,2,3. The turbulent energy frequency spectrum along the three coordinate axes in the (x1,x2,x3) coordinate system is:
[0057] Φ 11 m =dU1 m dU1 m*
[0058] Φ 22 m =dU2 m dU2 m*
[0059] Φ 33m =dU3 m dU3 m*
[0060] Where, Φ 11 m ,Φ 22 m ,Φ 33 m They are the turbulent energy frequency spectrum along the synthetic wind direction x1 axis, the turbulent energy frequency spectrum along the crosswind direction x2 axis, and the turbulent energy frequency spectrum along the vertical wind direction x3 axis, U1 m ,U2 m ,U3 m are the path-averaged wind speed along the composite wind direction x1 axis, the path-averaged wind speed along the crosswind direction x2 axis, and the path-averaged wind speed along the perpendicular wind direction x3 axis, respectively. * indicates conjugate.
[0061] Calculate the transfer function of the wind speed spectrum along the x1, x2, and x3 axes, and the expression is:
[0062]
[0063] Where T1, T2, T3 are the transmission functions of the wind speed spectrum along the x1 axis, along the x2 axis, and along the x3 axis, respectively; k1, k2, k3 are the transmission wave numbers along the x1 axis, along the x2 axis, and along the x3 axis, respectively;
[0064] Assuming that the spectral density tensor is isotropic in the frequency band above the inertial region, the expression is:
[0065]
[0066] Where k represents the wave number vector The model, that is i and j represent coordinate axes, i=1,2,3, j=1,2,3; δ ij is the Kronecker symbol, which is 1 when i=j and 0 otherwise; is the spectral density under three-dimensional wave number conditions, k i ,k j are the transmission wave number of a certain i-axis and the transmission wave number of a certain j-axis respectively;
[0067]
[0068] Among them, A u is a constant, and ∈ is the turbulent energy dissipation rate. By integration, we can get
[0069]
[0070] If the turbulence is considered to be isotropic, Then we have:
[0071]
[0072] Among them, the wave vector (k1, k2, k3) is converted to the cylindrical coordinate system (k1, Ksinβ, Kcosβ), and then substituted into:
[0073]
[0074] Where β is the cylindrical coordinate conversion angle, and K is the number of transmitted waves in the cylindrical coordinate system. Therefore, the transmission functions T1, T2, and T3 are affected by the actual horizontal wind direction angle α and the inclination angle γ.
[0075] Numerical solution of the path average transfer function T1, T2, T3 with dimensionless frequency k1l or changes.
[0076] In step S3, the non-orthogonal three-dimensional ultrasonic anemometer sensor ultrasonic temperature path average correction algorithm includes:
[0077] The ultrasonic temperature of a non-orthogonal three-dimensional ultrasonic anemometer is the average value along the three acoustic paths, so the influence of the acoustic path average must be considered. When the temperature field variable X is isotropic, the expression is:
[0078]
[0079] Where X(s) is the three-dimensional temperature field variable, i is an imaginary number, s is the temperature field radius, is the three-dimensional wave vector, is a random complex function with orthogonality:
[0080]
[0081]
[0082] Where, is the three-dimensional spectrum distribution of the temperature scalar X, is the three-dimensional spectral density function of the temperature scalar X, Z * is three-dimensional conjugated, k i and k j Along a certain x i Axis and some x j The random transmission wave number of the axis;
[0083] The actual observed temperature is the average value of the acoustic path of three pairs of sensors in three-dimensional space, then:
[0084]
[0085] Where, is the actual temperature observed at radius s, x0 is the initial point of the temperature field (the common center position of the three pairs of sensors in the three-dimensional ultrasound system), l a ,l b ,l c are the sound path vectors along the direction of the first pair of sensors, the direction of the second pair of sensors, and the direction of the third pair of sensors respectively;
[0086] After simplification, we have:
[0087]
[0088] Where, H() is the function of the influence on the temperature field after considering the average of the sound path;
[0089] The correlation function of the temperature field along the synthetic wind direction x1 axis is:
[0090]
[0091] Where, is the correlation function of the temperature field along the synthetic wind direction, i is an imaginary number;
[0092] The one-dimensional spatial spectral density along the synthetic wind direction x1 is:
[0093]
[0094] Where i is an imaginary number, is the average value of the one-dimensional spatial spectrum density along the synthetic wind direction x1, ξ is the spatial position;
[0095] When the temperature field is a uniform free field, the frequency spectrum density above the inertial region satisfies the Kolmogorov law, which is similar to the turbulent velocity field. Then:
[0096]
[0097] Where A T is a constant;
[0098]
[0099] Where F1(k1) is the one-dimensional spatial spectral density value under the transmission wave number along the x1 axis;
[0100] Wave number vector Converted to the cylindrical coordinate system (k1, Ksinβ, Kcosβ), the final temperature transfer function is as follows:
[0101]
[0102] Where, T TT (k1) is the transfer function of the temperature spectrum.
[0103] In step S4, the numerical integration of the wind speed spectrum loss correction and the ultrasonic temperature spectrum loss correction at different dimensionless frequencies is used to create a lookup table, including:
[0104] The values of wind speed spectrum loss correction and ultrasonic temperature spectrum loss correction under different dimensionless frequency conditions are respectively integrated and calculated, and corresponding lookup tables are made; in the actual wind speed spectrum loss correction and ultrasonic temperature spectrum loss calculation, the corresponding dimensionless frequency is calculated according to the actual wind speed, sound path and frequency sequence. If the points at two adjacent positions in the start node and the end node of the dimensionless frequency sequence are within the dimensionless frequency range of the lookup table, linear interpolation is used between these positions; otherwise, the positions of other nodes are kept consistent with the maximum dimensionless frequency and minimum dimensionless frequency positions in the lookup table.
[0105] Another object of the present invention is to provide an improved eddy covariance observation system turbulent flux loss correction calculation system, which implements the improved eddy covariance observation system turbulent flux loss correction calculation method, and the system includes:
[0106] The flow field distortion correction module corrects the flow field distortion effect by using the sensor geometry of the non-orthogonal 3D ultrasonic anemometer based on the instantaneous wind speed and direction in the original turbulence observation;
[0107] The wind speed spectrum loss correction module is used to use the path averaging algorithm of the non-orthogonal three-dimensional ultrasonic anemometer sensor to correct the wind speed spectrum loss caused by path averaging in the synthetic wind direction x1, the crosswind direction x2, and the vertical wind direction x3;
[0108] Ultrasonic temperature spectrum loss correction module, used to use the non-orthogonal three-dimensional ultrasonic anemometer sensor path averaging correction algorithm to correct the ultrasonic temperature spectrum loss caused by path averaging;
[0109] The turbulent energy loss correction module is used to convert the numerical integration of wind speed spectrum loss correction and ultrasonic temperature spectrum loss correction at different dimensionless frequencies into a lookup table. In the frequency loss correction calculation of the turbulence spectrum, a method combining analytical solution with linear interpolation of the above lookup table is used to obtain the turbulent energy loss correction information caused by path averaging in the synthetic wind direction x1, crosswind direction x2, and vertical wind direction x3.
[0110] Furthermore, the system is mounted on a non-orthogonal three-dimensional ultrasonic anemometer to implement the functions of the above method.
[0111] Furthermore, the system is equipped with an open-circuit or closed-circuit water vapor / carbon dioxide analyzer or a trace gas analyzer to implement the functions of the method and obtain the flux of turbulent transport of momentum, energy and matter between the surface and the atmosphere.
[0112] Furthermore, the system is mounted on a computer-readable storage medium, which stores a computer program. When the computer program is executed by a processor, it can realize the functions for implementing the above method.
[0113] Furthermore, the system is mounted on an information data processing terminal, and when the information data processing terminal is used to be executed on an electronic device, it provides a user input interface to implement the functions of the above method.
[0114] Combining all the above technical solutions, the beneficial effects of the present invention are as follows: the present invention further improves the path average loss correction of the non-orthogonal three-dimensional ultrasonic anemometer, and adopts a combination of analytical and numerical integration lookup table linear interpolation in the actual flow field calculation to correct the turbulent energy loss caused by path averaging in the x1, x2 and x3 directions respectively, thereby ensuring the accuracy of the flow field calculation while taking into account the calculation speed. BRIEF DESCRIPTION OF THE DRAWINGS
[0115] The accompanying drawings, which are incorporated in and constitute a part of this specification, illustrate embodiments consistent with the present disclosure and, together with the description, serve to explain the principles of the present disclosure;
[0116] Figure 1 This is a flow chart of a method for calculating turbulent flux loss correction in an improved eddy covariance observation system provided by an embodiment of the present invention;
[0117] Figure 2 Schematic diagram of a turbulent flux loss correction calculation system for an improved eddy covariance observation system provided by an embodiment of the present invention;
[0118] In the figure: 1. Convection field distortion correction module; 2. Wind speed spectrum loss correction module; 3. Ultrasonic temperature spectrum loss correction module; 4. Turbulent energy loss correction module. DETAILED DESCRIPTION
[0119] To make the above-mentioned objects, features, and advantages of the present invention more readily apparent, specific embodiments of the present invention are described in detail below with reference to the accompanying drawings. The following description sets forth numerous specific details to facilitate a full understanding of the present invention. However, the present invention can be implemented in many other ways than those described herein, and those skilled in the art may make similar modifications without departing from the scope of the present invention. Therefore, the present invention is not limited to the specific embodiments disclosed below.
[0120] Example 1, as Figure 1As shown, the improved eddy covariance observation system turbulent flux loss correction calculation method provided by the embodiment of the present invention further improves the path average loss correction of non-orthogonal three-dimensional ultrasonic anemometers. During the actual flow field calculation, a combination of analytical and numerical integration lookup table linear interpolation is used to correct the turbulent energy loss caused by path averaging in the synthetic wind direction x1, crosswind direction x2, and perpendicular wind direction x3. This ensures the accuracy of the flow field calculation while taking into account the calculation speed. Specifically, the following steps are included:
[0121] S1, based on the instantaneous wind speed and direction in the original turbulence observation, the flow field distortion effect is corrected by the sensor geometry of the non-orthogonal three-dimensional ultrasonic anemometer;
[0122] S2, using the sound path averaging algorithm of the non-orthogonal three-dimensional ultrasonic anemometer sensor, the wind speed spectrum loss caused by the sound path averaging in the synthetic wind direction x1, the crosswind direction x2, and the vertical wind direction x3 is corrected respectively;
[0123] S3, using the acoustic path averaging correction algorithm of the non-orthogonal three-dimensional ultrasonic wind temperature meter sensor to correct the ultrasonic temperature spectrum loss caused by the acoustic path averaging;
[0124] S4. The numerical integration of the wind speed spectrum loss correction and the ultrasonic temperature spectrum loss correction at different dimensionless frequencies is made into a lookup table. In the frequency loss correction calculation of the turbulence spectrum, a method combining the analytical solution with the linear interpolation of the above lookup table is used to obtain the turbulent energy loss correction information caused by the average sound path in the synthetic wind direction x1, the crosswind direction x2, and the vertical wind direction x3.
[0125] In step S4, the numerical integration of wind speed spectrum loss correction and ultrasonic temperature spectrum loss correction at different dimensionless frequencies is made into a lookup table, including:
[0126] The values of wind speed spectrum loss correction and ultrasonic temperature spectrum loss correction under different dimensionless frequency conditions are integrated and calculated, and corresponding lookup tables are created. In the actual wind speed spectrum loss correction and ultrasonic temperature spectrum loss calculation, the corresponding dimensionless frequency is calculated based on the actual wind speed, sound path, and frequency sequence. If two adjacent points at the start and end nodes of the dimensionless frequency sequence are within the dimensionless frequency range of the lookup table, linear interpolation is used between these points. Otherwise, the positions of other nodes are consistent with the maximum dimensionless frequency (minimum dimensionless frequency) position of the lookup table.
[0127] Example 2: In practice, the correction of the effect of sensor geometry on flow field distortion includes:
[0128] Before path average correction, it is necessary to consider the effect of the sensor wake on the flow field distortion on the sensor sound path. The actual wind speed detected by the sensor is the affected wind speed u. Am, the true wind speed u along the sensor path A The ratio of the actual wind direction to the sensor path is θ A related:
[0129] u Am =u A (0.84+0.16sinθ A )
[0130] or
[0131]
[0132] Where u Am is the wind speed detected by the first pair of sensors after the impact, θ A is the angle between the actual wind direction and the path of the first pair of sensors, u A It is the wind speed projected from the actual wind speed to the direction of the first pair of sensors; the directions of the other two pairs of sensors are similar;
[0133] However, for non-orthogonal three-dimensional ultrasonic anemometers, the paths between the three pairs of sensors are three-dimensional geometric shapes, and there is cross-obstruction. The vertical velocity is obtained by matrix transformation after the obstruction of the sensors and supporting structures in the directions of the three pairs of sensors is reduced. The vertical velocity variance and sensible heat flux observed by non-orthogonal ultrasonic anemometers are 8-9% and more than 10% lower than those of orthogonal ones. Therefore, after obtaining the raw wind speed of the three-dimensional ultrasonic anemometer, it is also necessary to know whether the flow field distortion correction has been performed in the instrument settings or software settings, and then decide whether additional flow field distortion correction is needed for the wind speed data.
[0134] Example 3, wind speed path averaging algorithm for non-orthogonal three-dimensional ultrasonic anemometer sensor.
[0135] Assume that the actual horizontal wind direction angle averaged over a period of time (e.g., 30 minutes) in the ultrasonic system is α. Its three-dimensional synthetic wind speed may not be located on the horizontal plane of the ultrasonic system, but may be at an angle γ to the horizontal plane of the ultrasonic system. Take the positive angle γ when the average vertical speed in the ultrasonic system is upward. The new coordinate system is the synthetic wind direction x1, the crosswind direction x2, and the vertical wind direction x3 coordinate system. The actual synthetic wind speed u1 is along the synthetic wind direction (x1 axis), u2 is the crosswind direction (i.e., x2 axis), and u3 is the vertical wind direction (i.e., x3 axis). Then, in the (x1, x2, x3) coordinate system, the sound path vectors along the three pairs of sensor directions are:
[0136]
[0137] Where, l a ,l b ,l care the sound path vectors of the first, second and third pairs of sensors in the (x1, x2, x3) coordinate system, respectively. are the angles between the first, second, and third pairs of sensors when projected onto the horizontal plane and the positive x-axis in the ultrasound system. θ is the angle between the sensor and the positive z-axis of the ultrasound system. In fact, the angles between the three pairs of sensors and the positive z-axis are the same. l is the acoustic path length (constant). E, B, and D are rotation matrices. T represents the transpose of the matrix.
[0138]
[0139] The matrix E is related to whether the ultrasound system is tilted. If the ultrasound system is installed vertically on the ground and does not tilt over time, the xy coordinate system of the ultrasound system is parallel to the actual horizontal plane, and the x2 axis is the result of the y axis of the ultrasound system rotating counterclockwise by an angle α on its horizontal plane, that is, If the ultrasonic system is tilted, for example, in the ultrasonic system, the ultrasonic instrument is at the γ0 position of the xy coordinate of the ultrasonic system, and the z axis is tilted at an angle of θ0 with the actual vertical direction. Then, after the above rotation, it can always be rotated counterclockwise along the synthetic wind direction x1 axis. Angle so that the generated x2 axis is parallel to the actual horizontal plane.
[0140] If the effect of path averaging is not considered, the wind speed (u A ,u B ,u C ), the transformation between this wind speed and the three-dimensional wind speed (u1, u2, u3) in the (x1, x2, x3) coordinate system satisfies the following matrix:
[0141]
[0142] Among them, A is the wind field conversion matrix; A -1 is the inverse transformation matrix of the wind field.
[0143]
[0144] Where,
[0145] Let matrix M = ABDE, N = E T D T B T A -1 After considering the effect of path averaging, the wind speed (u A ,u B ,u C ), which is actually:
[0146]
[0147] Where M ij is the value of the i-th row and j-th column of matrix M, are the wind speed u1 along the direction of the first pair of sensors, the direction of the second pair of sensors, and the direction of the third pair of sensors after considering the average effect of their respective sound paths; are the wind speed u2 along the first pair of sensors, the second pair of sensors, and the third pair of sensors, respectively, after considering the average effect of their respective sound paths. are the wind speed u3 along the direction of the first pair of sensors, the direction of the second pair of sensors, and the direction of the third pair of sensors after considering the average effect of their respective sound paths, and ~ is the average sound path;
[0148]
[0149]
[0150] Where x0 is the common center position of the three pairs of sensors in the 3D ultrasound system;
[0151] In the (x1, x2, x3) coordinate system, the actual three-dimensional wind speed (u1 m ,u2 m ,u3 m ), the expression is:
[0152]
[0153] After substitution, we have:
[0154]
[0155] Considering that the observed wind speed is the contribution of eddies of different wave numbers, the Fourier transform yields:
[0156]
[0157] In the above three formulas, e is the exponent and i represents the imaginary number. is the three-dimensional wave number, and the corresponding wave numbers in the directions (x1, x2, x3) are (k1, k2, k3). Accordingly, the wind speed averaged along the direction of the synthetic wind speed is as follows: Three-dimensional wind speed components transformed in the frequency domain satisfy
[0158]
[0159] Where i represents an imaginary number; according to the above formula, we have:
[0160]
[0161] Where d is the component;
[0162]
[0163] j=1,2,3
[0164] The turbulent energy frequency spectrum along the three coordinate axes in the (x1, x2, x3) coordinate system:
[0165] Φ 11 m =dU1 m dU1 m*
[0166] Φ 22 m =dU2 m dU2 m*
[0167] Φ 33 m =dU3 m dU3 m*
[0168] Where, Φ 11 m ,Φ 22 m ,Φ 33 m They are the turbulent energy frequency spectrum along the synthetic wind direction x1 axis, the turbulent energy frequency spectrum along the crosswind direction x2 axis, and the turbulent energy frequency spectrum along the vertical wind direction x3 axis, U1 m ,U2 m ,U3 m are the sound path averaged wind speed along the synthetic wind direction x1 axis, the path averaged wind speed along the crosswind direction x2 axis, and the path averaged wind speed along the perpendicular wind direction x3 axis, * indicates conjugate;
[0169] Calculate the transfer function of the wind speed spectrum along the x1, x2, and x3 axes, and the expression is:
[0170]
[0171] Where T1, T2, T3 are the transmission functions of the wind speed spectrum along the x1 axis, along the x2 axis, and along the x3 axis, respectively; k1, k2, k3 are the transmission wave numbers along the x1 axis, along the x2 axis, and along the x3 axis, respectively;
[0172] Assuming that the spectral density tensor is isotropic in the frequency band above the inertial region, the expression is:
[0173]
[0174] Where k represents the wave number vector The model, that is i and j represent coordinate axes, i=1,2,3, j=1,2,3; δ ij is the Kronecker symbol, which is 1 when i=j and 0 otherwise; is the spectral density under three-dimensional wave number conditions, k i ,k j are the transmission wave number of a certain i-axis and the transmission wave number of a certain j-axis respectively;
[0175]
[0176] Among them, A u is a constant, and ∈ is the turbulent energy dissipation rate. By integration, we can get:
[0177]
[0178] If the turbulence is considered to be isotropic, Then we have:
[0179]
[0180] Among them, (k1, k2, k3) is converted to the cylindrical coordinate system (k1, Ksinβ, Kcosβ), and after being brought in,
[0181]
[0182]
[0183] Where β is the cylindrical coordinate conversion angle, and K is the transmitted wave number in the cylindrical coordinate system;
[0184] Numerical solution of the path average transfer function T1, T2, T3 with dimensionless frequency k1l or changes. The path average correction function of the non-orthogonal ultrasonic anemometer CSAT3 is compared with the one-dimensional correction. It can be seen that the velocity spectrum transfer function along the synthetic wind direction is affected by the angle α between the wind direction and the x-axis of the three-dimensional ultrasonic coordinate system and the angle of attack γ. Relatively speaking, the loss in the crosswind direction is very weak, and there are differences in the high-frequency part from the one-dimensional correction curve in the different directions of x1, x2, and x3. Therefore, there are obvious differences in the transfer function in the directions of different coordinate axes x1, x2, and x3 in the three-dimensional coordinate system, especially when the synthetic wind speed is weak. If the frequency-varying characteristics of the turbulence spectrum in all directions are considered, such as the anisotropic analysis of turbulence and the non-similar analysis of turbulent transport, it is very necessary to consider the x1, x2, and x3 directions separately. Generally, during observation, the direction of the synthetic wind is basically along the horizontal direction of the ultrasonic instrument, and the influence of the angle of attack γ is basically very weak, so the coordinate axis x1 can be approximately taken along the horizontal direction; if the three-dimensional ultrasonic instrument has a certain tilt or the angle of attack γ is large, the velocity spectrum transfer function along the synthetic wind direction is obviously affected by the angle α between the wind direction and the x-axis of the three-dimensional ultrasonic coordinate system and the angle of attack γ in the higher frequency range.
[0185] Example 4, an average correction algorithm for ultrasonic temperature paths of non-orthogonal three-dimensional ultrasonic anemometer sensors.
[0186] Correspondingly, the ultrasonic temperature of the non-orthogonal three-dimensional ultrasonic anemometer is also the average of three paths. It is different from the single path average of the gas analyzer and the single point observation of the thermocouple. The influence of the path average needs to be considered. Assume that the temperature field variable X is isotropic:
[0187]
[0188] Where X(s) is the three-dimensional temperature field variable, i is an imaginary number, s is the temperature field radius, is the three-dimensional wave vector, is a random complex function with orthogonality:
[0189]
[0190]
[0191] Where, is the three-dimensional spectrum distribution of the temperature scalar X, is the three-dimensional spectral density function of the temperature scalar X, Z * is three-dimensional conjugated, k i and k j Along a certain x i Axis and some x j The random transmission wave number of the axis;
[0192] The actual observed temperature is the average value of the acoustic path of three pairs of sensors in three-dimensional space, then:
[0193]
[0194] Where, is the actual temperature observed at radius s, r0 is the initial point of the temperature field (the common center position of the three pairs of sensors in the three-dimensional ultrasound system), l a ,l b ,l c are the sound path vectors along the direction of the first pair of sensors, the direction of the second pair of sensors, and the direction of the third pair of sensors respectively;
[0195] After simplification, we have:
[0196]
[0197] Where, H() is the function of the influence on the temperature field after considering the average of the sound path;
[0198] The correlation function of the temperature field along the synthetic wind direction x1 axis is:
[0199]
[0200] Where, is the correlation function of the temperature field along the synthetic wind direction, i is an imaginary number;
[0201] The one-dimensional spatial spectral density along the synthetic wind direction x1 is:
[0202]
[0203] Where, is the average value of the one-dimensional spatial spectrum density along the synthetic wind direction x1, ξ is the spatial position;
[0204] Assuming that the temperature field is a uniform free field, the frequency spectrum density above the inertial region satisfies the Kolmogorov law, which is similar to the turbulent velocity field. Then:
[0205]
[0206] Where A T is a constant;
[0207]
[0208] Where F1(k1) is the one-dimensional spatial spectral density value under the transmission wave number along the x1 axis;
[0209] Wave number vector Converted to the cylindrical coordinate system (k1, K sinβ, K cosβ), the final transfer function of the temperature spectrum is as follows:
[0210]
[0211] Where, T TT (k1) is the transfer function of the temperature spectrum.
[0212] Example 5, as Figure 2 As shown, an embodiment of the present invention provides an improved eddy covariance observation system turbulent flux loss correction calculation system, comprising:
[0213] Flow field distortion correction module 1, based on the instantaneous wind speed and direction in the original turbulence observation, corrects the flow field distortion effect through the sensor geometry of the non-orthogonal three-dimensional ultrasonic anemometer;
[0214] Wind speed spectrum loss correction module 2 is used to correct the wind speed spectrum loss caused by path averaging in the synthetic wind direction x1, the crosswind direction x2, and the vertical wind direction x3 respectively using the path averaging algorithm of the non-orthogonal three-dimensional ultrasonic anemometer sensor;
[0215] Ultrasonic temperature spectrum loss correction module 3, used to correct the ultrasonic temperature spectrum loss caused by path averaging using a non-orthogonal three-dimensional ultrasonic anemometer sensor path averaging correction algorithm;
[0216] Turbulent energy loss correction module 4 uses the numerical integration of wind speed spectrum loss correction and ultrasonic temperature spectrum loss correction at different dimensionless frequencies to create a lookup table. In the frequency loss correction calculation of the turbulence spectrum, a method combining an analytical solution with linear interpolation of the above lookup table is adopted to obtain the turbulent energy loss correction information caused by the average sound path in the synthetic wind direction x1, crosswind direction x2, and vertical wind direction x3.
[0217] It can be seen from the above embodiments that the path average loss correction of the non-orthogonal three-dimensional ultrasonic anemometer has been further improved. In actual calculations, a combination of analytical and numerical integration lookup table linear interpolation is adopted to correct the turbulent energy loss caused by path averaging in the x1, x2 and x3 directions respectively, thereby ensuring the accuracy of the calculation while taking into account the calculation speed.
[0218] In the above embodiments, the description of each embodiment has its own focus. For parts that are not described or recorded in detail in a certain embodiment, reference can be made to the relevant description of other embodiments.
[0219] The above description is only a preferred specific implementation method of the present invention, but the scope of protection of the present invention is not limited thereto. Any modifications, equivalent substitutions and improvements made by any technician familiar with this technical field within the technical scope disclosed by the present invention and within the spirit and principles of the present invention should be covered by the scope of protection of the present invention.
Claims
1. An improved turbulent flux loss correction calculation method for eddy covariance observation system, characterized by: The method comprises the following steps: S1, based on the instantaneous wind speed and direction in the original turbulence observation, the flow field distortion effect is corrected by the sensor geometry of the non-orthogonal three-dimensional ultrasonic anemometer; S2, using the path averaging algorithm of the non-orthogonal three-dimensional ultrasonic anemometer sensor, corrects the wind speed spectrum loss caused by path averaging in the synthetic wind direction x1, crosswind direction x2, and vertical wind direction x3 averaged over a period of time; S3, using the non-orthogonal three-dimensional ultrasonic anemometer sensor path averaging correction algorithm to correct the ultrasonic temperature spectrum loss caused by path averaging; S4, the numerical integration of wind speed spectrum loss correction and ultrasonic temperature spectrum loss correction at different dimensionless frequencies is made into a lookup table. In the frequency loss correction calculation of the turbulence spectrum, a method combining the analytical solution with the linear interpolation of the above lookup table is used to obtain the turbulent energy loss correction information caused by path averaging in the synthetic wind direction x1, crosswind direction x2, and vertical wind direction x3.
2. The improved eddy covariance observation system turbulent flux loss correction calculation method according to claim 1 is characterized in that: In step S1, the effect of the sensor geometry of the non-orthogonal 3D ultrasonic anemometer on flow field distortion is corrected, including: Before path averaging correction, calculate the effect of sensor wake on flow field distortion on sensor acoustic path, such as the wind speed u detected by the first pair of sensors after the effect. Am The angle θ between the actual wind direction and the path of the first pair of sensors A The expression is: u Am =u A (0.84+0.16sinθ A ) in, Where u Am is the wind speed detected by the first pair of sensors after the impact, θ A is the angle between the actual wind direction and the path of the first pair of sensors, u A The actual wind speed is projected onto the direction of the first pair of sensors. The algorithms for the other two pairs of sensors are the same as those for the first pair of sensors. After obtaining the original wind speed of the three-dimensional ultrasonic anemometer, when the flow field distortion correction is required, the three-dimensional wind speed (u x ,u y ,u z ) Refer to the following formula for conversion: Where u B ,u C are the wind speeds projected onto the paths of the second and third pairs of sensors, respectively. x ,u y ,u z are the wind speeds along the x-axis, y-axis, and z-axis of the three-dimensional coordinate system in the ultrasonic system, respectively. A is the wind field conversion matrix, which is expressed as follows: Where, are the angles between the first, second, and third pairs of sensors when projected onto the horizontal plane and the positive half-axis of the x-axis in the ultrasonic system. θ is the angle between the sensor and the positive half-axis of the z-axis in the ultrasonic system. The angles between the three pairs of sensors of the non-orthogonal ultrasonic anemometer and the positive half-axis of the z-axis in the ultrasonic system are the same. Converted to wind speed along the wind direction of three pairs of sensors (u A ,u B ,u C ), after correcting the flow field distortion, it is converted into the ultrasonic coordinate system (u x ,u y ,u z ): Where A -1 That is the inverse transformation matrix of the wind field, namely: in, 3. The improved eddy covariance observation system turbulent flux loss correction calculation method according to claim 2 is characterized in that: In step S2, the non-orthogonal three-dimensional ultrasonic anemometer sensor path averaging algorithm includes: Suppose the average horizontal wind direction angle over a period of time in the ultrasonic system is α, and its three-dimensional synthetic wind speed is not located on the horizontal plane of the ultrasonic system, but has an angle γ with the horizontal plane of the ultrasonic system. Take the average vertical speed in the ultrasonic system as the positive angle γ, and define the new coordinate system as the synthetic wind direction x1, side wind direction x2, and vertical wind direction x3 coordinate system. The actual synthetic wind speed u1 is along the synthetic wind direction x1 axis, u2 is along the side wind direction x2 axis, parallel to the actual horizontal plane, and u3 is along the vertical wind direction x3 axis. If the ultrasound system is installed vertically on the ground and does not tilt over time, the xy coordinate system of the ultrasound system is parallel to the actual horizontal plane, then the x2 axis is the result of rotating the y axis of the ultrasound system counterclockwise about an angle α on its horizontal plane. If the ultrasonic system is tilted, for example, in the ultrasonic system, the ultrasonic instrument is at the γ0 position of the xy coordinate of the ultrasonic system, and the z axis is tilted at an angle of θ0 with the actual vertical direction. Then, after the above rotation, it can be rotated counterclockwise along the synthetic wind direction x1 axis. Angle so that the generated x2 axis is parallel to the actual horizontal plane. The calculation of is as follows: Taking the ultrasonic system of the vertically mounted ultrasonic instrument as a reference, in this coordinate system, after the ultrasonic instrument is tilted, according to the Rodrigues matrix rotation algorithm, that is, rotating counterclockwise along the normal vector n1 by an angle θ0, the tilted ultrasonic system coordinate system is equivalent to the original vertical ultrasonic system multiplied by the rotation matrix. The rotation matrix is Then, the three-dimensional coordinate system is rotated counterclockwise by an angle α along its new z-axis (normal vector n2) (i.e., the new x-axis points to the horizontal synthetic wind direction of the ultrasonic system), and the rotation matrix is On this basis, the newly generated three-dimensional coordinate system is rotated counterclockwise along its new y-axis (normal vector n3) by an angle of γ (the new x-axis points to the direction of the three-dimensional synthetic wind, and the average wind speed in the new z-axis direction is 0). The rotation matrix is n3 is the second column of the matrix R2(n2,α)×R1(n1,θ0). Therefore, assuming n4 is the first column of the matrix R3(n3,γ)×R2(n2,α)×R1(n1,θ0), the newly generated three-dimensional coordinate system rotates counterclockwise along the new x-axis (normal vector n4). After the angle is adjusted, the new y-axis is parallel to the actual horizontal plane, then Then according to the vector n3 and After the cross product, the dot product value of vector n4 is determined to be positive or negative. If n4(1) and n4(2) are both 0, that is, the n4 normal vector is parallel to the actual vertical direction, and the yz plane of the new three-dimensional coordinate system is parallel to the actual horizontal plane, in this case, directly take is 0. Then, (u1,u2,u3) and (u x ,u y ,u z )The relationship is: Where T is the transpose of the matrix. In the (x1, u2, x3) coordinate system, the acoustic path vectors along the three pairs of sensor directions are: Where, l a ,l b ,l c are the sound path vectors of the first, second and third pairs of sensors in the (x1, x2, x3) coordinate system, respectively, and l is the sound path length; Without considering the influence of path averaging, the wind speed detected by the three pairs of ultrasonic sensors is (u A ,u B ,u C ), the transformation between this wind speed and the three-dimensional wind speed (u1, u2, u3) in the (x1, x2, x3) coordinate system satisfies the following matrix: Let matrix M = ABDE, N = E T D T B T A -1 After considering the influence of path or sound path average, the wind speed (u A ,u B ,u C ), which is actually: Where M ij is the value of the i-th row and j-th column of matrix M, are the wind speed u1 along the direction of the first pair of sensors, the direction of the second pair of sensors, and the direction of the third pair of sensors after considering the average effect of their respective sound paths; are the wind speed u2 along the direction of the first pair of sensors, the direction of the second pair of sensors, and the direction of the third pair of sensors after considering the average effect of their respective sound paths, are the wind speed u3 along the first pair of sensors, the second pair of sensors, and the third pair of sensors, respectively, after considering the average effect of their respective sound paths. ~ represents the path average; Where x0 is the common center position of the three pairs of sensors in the 3D ultrasound system; In the (x1, x2, x3) coordinate system, the actual three-dimensional wind speed (u1 m ,u2 m ,u3 m ), the expression is: After substitution, we have: Considering that the observed wind speed is the contribution of eddies of different wave numbers, the Fourier transform yields: In the above three formulas, e is the exponent and i represents the imaginary number. is the three-dimensional wave number, and the corresponding wave numbers in the (x1, x2, x3) directions are (k1, k2, k3) respectively; Correspondingly, the wind speed along the direction of the composite wind speed after path averaging is as follows: Three-dimensional wind speed components transformed in the frequency domain satisfy: According to the above formula, we have: Where d is the component; in, j=1,2,3, the turbulent energy frequency spectrum along the three coordinate axes in the (x1,x2,x3) coordinate system: Φ 11 m =dU1 m dU1 m* Φ 22 m =dU2 m dU2 m* Φ 33 m =dU3 m dU3 m* Where, Φ 11 m ,Φ 22 m ,Φ 33 m They are the turbulent energy frequency spectrum along the synthetic wind direction x1 axis, the turbulent energy frequency spectrum along the crosswind direction x2 axis, and the turbulent energy frequency spectrum along the vertical wind direction x3 axis, U1 m ,U2 m ,U3 m are the path-averaged wind speed along the composite wind direction x1 axis, the path-averaged wind speed along the crosswind direction x2 axis, and the path-averaged wind speed along the perpendicular wind direction x3 axis, respectively. * indicates conjugate. Calculate the transfer function of the wind speed spectrum along the x1, x2, and x3 axes, and the expression is: Where T1, T2, T3 are the transmission functions of the wind speed spectrum along the x1 axis, along the x2 axis, and along the x3 axis, respectively; k1, k2, k3 are the transmission wave numbers along the x1 axis, along the x2 axis, and along the x3 axis, respectively; In the frequency band above the inertial region, the spectral density tensor is isotropic and is expressed as: Where k represents the wave number vector The model, that is i and j represent coordinate axes, i=1,2,3, j=1,2,3; δ ij is the Kronecker symbol, which is 1 when i=j and 0 otherwise; is the spectral density under three-dimensional wave number conditions, k i ,k j are the transmission wave number of a certain i-axis and the transmission wave number of a certain j-axis respectively; Among them, A u is a constant, ∈ is the turbulent energy dissipation rate; From the integral we get: If the turbulence is considered to be isotropic, Then we have: Among them, the wave vector (k1, k2, k3) is converted to the cylindrical coordinate system (k1, Ksinβ, Kcosβ), and then it is: Where β is the cylindrical coordinate conversion angle, and K is the transmitted wave number in the cylindrical coordinate system; Numerical solution of the path average transfer function T1, T2, T3 with dimensionless frequency k1l or changes.
4. The improved eddy covariance observation system turbulent flux loss correction calculation method according to claim 1 is characterized in that: In step S3, the non-orthogonal three-dimensional ultrasonic anemometer sensor path average correction algorithm includes: The ultrasonic temperature of a non-orthogonal three-dimensional ultrasonic anemometer is the average value along the three acoustic paths, and the influence of the acoustic path average must be considered. When the temperature field variable X is isotropic, the expression is: Where X(s) is the three-dimensional temperature field variable, i is an imaginary number, s is the temperature field radius, is the three-dimensional wave vector, is a random complex function with orthogonality: Where, is the three-dimensional spectrum distribution of the temperature scalar X, is the three-dimensional spectral density function of the temperature scalar X, Z * is three-dimensional conjugated, k i and k j Along a certain x i Axis and some x j The random transmission wave number of the axis; The actual observed temperature is the average value of the acoustic path of three pairs of sensors in three-dimensional space, then: Where, is the actual temperature observed at radius s, x0 is the common center position of the three pairs of sensors of the three-dimensional ultrasonic system at the initial point of the temperature field, l a ,l b ,l c are the sound path vectors along the direction of the first pair of sensors, the direction of the second pair of sensors, and the direction of the third pair of sensors respectively; After simplification, we have: Where, H( ) is the function of the influence on the temperature field after considering the average of the sound path; The correlation function of the temperature field along the synthetic wind direction x1 axis is: Where, is the correlation function of the temperature field along the synthetic wind direction, i is an imaginary number; The one-dimensional spatial spectral density along the synthetic wind direction x1 is: Where i is an imaginary number, is the average value of the one-dimensional spatial spectrum density along the synthetic wind direction x1, ξ is the spatial position; When the temperature field is a uniform free field, the frequency spectrum density above the inertial region satisfies the Kolmogorov law, which is similar to the turbulent velocity field. Then: Where A T is a constant; Where F1(k1) is the one-dimensional spatial spectral density value under the transmission wave number along the x1 axis; Wave number vector Converted to the cylindrical coordinate system (k1, K sinβ, K cosβ), the final transfer function of the temperature spectrum is as follows: Where, T TT (k1) is the transfer function of the temperature spectrum.
5. The improved eddy covariance observation system turbulent flux loss correction calculation method according to claim 1 is characterized in that: In step S4, the numerical integration of the wind speed spectrum loss correction and the ultrasonic temperature spectrum loss correction at different dimensionless frequencies is used to create a lookup table, including: The values of wind speed spectrum loss correction and ultrasonic temperature spectrum loss correction under different dimensionless frequency conditions are respectively integrated and calculated, and corresponding lookup tables are made; in the actual wind speed spectrum loss correction and ultrasonic temperature spectrum loss calculation, the corresponding dimensionless frequency is calculated according to the actual wind speed, sound path and frequency sequence. If the points at two adjacent positions in the start node and the end node of the dimensionless frequency sequence are within the dimensionless frequency range of the lookup table, linear interpolation is used between these positions; otherwise, the positions of other nodes are kept consistent with the maximum dimensionless frequency and minimum dimensionless frequency positions in the lookup table.
6. An improved turbulent flux loss correction calculation system for eddy covariance observation system, characterized by: The system implements the improved turbulent flux loss correction calculation method of the eddy covariance observation system according to any one of claims 1 to 5, and the system comprises: The flow field distortion correction module (1) corrects the flow field distortion effect by using the sensor geometry of the non-orthogonal three-dimensional ultrasonic anemometer based on the instantaneous wind speed and direction in the original turbulence observation; A wind speed spectrum loss correction module (2) is used to respectively correct the wind speed spectrum loss caused by path averaging in the synthetic wind direction x1, the crosswind direction x2, and the vertical wind direction x3 using a non-orthogonal three-dimensional ultrasonic anemometer sensor path averaging algorithm; An ultrasonic temperature spectrum loss correction module (3) is used to correct the ultrasonic temperature spectrum loss caused by path averaging using a non-orthogonal three-dimensional ultrasonic anemometer sensor path averaging correction algorithm; The turbulent energy loss correction module (4) is used to convert the numerical integration of wind speed spectrum loss correction and ultrasonic temperature spectrum loss correction at different dimensionless frequencies into a lookup table. In the frequency loss correction calculation of the turbulence spectrum, a method combining an analytical solution with linear interpolation of the above lookup table is used to obtain the turbulent energy loss correction information caused by path averaging in the synthetic wind direction x1, the crosswind direction x2, and the vertical wind direction x3.
7. The improved eddy covariance observation system turbulent flux loss correction calculation system according to claim 6 is characterized in that: The system is equipped with a non-orthogonal three-dimensional ultrasonic anemometer and is used to implement the functions of the above method.
8. The improved eddy covariance observation system turbulent flux loss correction calculation system according to claim 6 is characterized in that: The system is equipped with an open-circuit or closed-circuit water vapor / carbon dioxide analyzer or a trace gas analyzer to implement the functions of the method and obtain the flux of turbulent transport of momentum, energy and matter between the surface and the atmosphere.
9. The improved eddy covariance observation system turbulent flux loss correction calculation system according to claim 6 is characterized in that: The system is mounted on a computer-readable storage medium, which stores a computer program. When the computer program is executed by a processor, it can realize the functions for implementing the above method.
10. The improved eddy covariance observation system turbulent flux loss correction calculation system according to claim 6, characterized in that: The system is mounted on an information data processing terminal, and when the information data processing terminal is used to be executed on an electronic device, it provides a user input interface to implement the functions of the above method.