Identiﬁcation of Wind Loads on Structures Based on Modal Kalman Filter with Unknown Inputs

: Wind loads on structures are difﬁcult to directly measure, so it is practical to identify structural wind loads based on the measurements of structural responses. However, this inversed problem is challenging compared with conventional load identiﬁcation as wind loads are time-space coupled and spatially distributed dynamic loads on structures. An improved method is proposed for identifying wind loads on structures using only partial measurements of structural acceleration responses in this paper. First, the wind loads on a structure are decomposed by proper orthogonal decomposition as a series of time-space decoupled sub-distributed dynamic loads with independent basic spatial distribution functions and time history functions. Herein, structural modes are adopted as the basic spatial distribution functions and structural modes of discretized and continuous structural systems are investigated. Then, a history function of the decomposed wind load is identiﬁed in the modal domain based on modal Kalman ﬁlter with unknown inputs, which is proposed by the authors. Finally, the distributed wind loads are reconstructed for discrete or continuous structural systems. The feasibility of the proposed algorithm is veriﬁed by two numerical examples of identiﬁcation of wind loads on a discrete shear frame and a wind turbine tower, respectively.


Introduction
Wind load information is important for the design and assessment of engineering structures, especially for high-rise buildings and large-span bridges. However, wind loads on the structures are usually unknown and are quite complicated because they are affected by many factors, such as topographic and geomorphic conditions, surrounding environment and structural morphology. In addition, the direct measurements of wind loads are not feasible in practical engineering. However, it is easier to measure the windinduced responses of structures. Therefore, it is more realistic to identify/reconstruct wind loads using partially measured wind-induced responses as the inversed problem.
In recent years, some wind load identification methods have been developed [1][2][3]. Law et al. [4] proposed the problem of time-varying wind load identification in state space and solved it with a regularization method. They simulated the fluctuating wind load as a multi-variable random function and identified the wind load at the full height of the structure by observing the displacement and strain response of the structure. Luigi et al. [5] and Federica et al. [6] used the orthogonal decomposition theory to establish a turbulent field of wind loads and considered the coherence of turbulence. Hwang et al. [7] proposed a method of modal wind load identification based on the Kalman filter and verified it with the data obtained from wind tunnel experiments. Subsequently, the Kalman filterbased method was used to estimate the modal wind load of structures with tuned mass dampers under typhoon action [8]. Zhi et al. [9][10][11] regarded the fluctuating wind load as zero-mean Gaussian white noise, and successfully identified the unknown fluctuating wind load of the structure by using continuous Kalman, discrete Kalman and Taylor polynomial expansion methods. It is a good method to assume the fluctuating wind load as zero-mean white noise or random walk process, but in practical engineering, it is difficult to distinguish the fluctuating wind response from the observed sensor response for wind load calculation. Niu et al. [12] measured the acceleration data of the structure through the long-term health monitoring system and identified the equivalent wind load acting on the Guangzhou Tower. Xu et al. [13] applied the direct inverse method to solve the generalized wind load by fully measuring the non-stationary wind vibration acceleration responses. Li et al. [14] studied the wind load characteristics of oval high-rise structures. Wu et al. [15] proposed a new wind load inversion method, which decomposed wind load into the product of orthogonal basis of position and time history function and took fully observed structural displacements, derived velocity and acceleration responses to obtain the wind load of the structure by least square estimation. Xue et al. [16,17] used the method of unbiased estimation to identify wind load under both known and unknown structural parameters. In a case study of the 1310 m long Hardanger Bridge, Petersen et al. [18,19] studied wind loads in time domain and frequency domain and introduced the Potential Gaussian Process model (GP-LFMS) to characterize the evolution of wind loads. Nasrabad et al. [20] presented a new wind speed estimation method by observing the deflection response of the tower and used the neural network to help estimate. Charisi et al. [21] determined the variation of wind pressure coefficient at different locations along the wind direction by on-site measurement of a double medium-sized building complex, which is of great significance to the subsequent calculation of wind pressure and wind load. For flexible slender structures, the interaction of the wind with the structure also causes aeroelastic effects. In this regard, many scholars [22][23][24] have also carried out corresponding research on aerodynamic damping and aeroelastic loads. Currently, in the study of identifying wind loads, some studies [9][10][11] assume that the wind load was zero mean white noise or random walk process, and the influence of the average wind on the structure was ignored. Therefore, it is necessary to identify both fluctuating and mean winds. In addition, modal truncation is also a very common method to reduce the number of observations. When the modal wind load is transformed into the physical coordinate wind load, the modal shape is no longer a square array after truncation. In some studies [16,17] the pseudo-inverse of the modal matrix was calculated and the modal wind loads were converted to wind loads in physical coordinates, but the use of the pseudo-inverse resulted in non-unique results. How to avoid the pseudo-inverse of the matrix is also content to be studied. Some studies [18,19,25] only identify modal wind loads and do not convert them into loads in physical coordinates. Although the modal wind load also has certain practical significance, its effect is not as large as that of the wind load in physical coordinates. Moreover, most of the current studies are aimed at shear frames or high-rise buildings, which are generally regarded as discrete structures, and the identification of wind loads acting on continuous system structures is not very involved.
An improved method for identifying wind loads is proposed in this paper to realize real-time identification of unknown wind loads of structures based on the partial windinduced responses of the structures. The wind loads are assumed to be composed of two mutually independent functions about position and time, and the spatial distribution function about position is expanded by orthogonal polynomial about vibration mode. In this paper, formulas are derived and calculated for both continuous and discrete systems. After that, the time function obtained by decomposition is combined with the modal Kalman filter with unknown inputs (MKF-UI) developed by authors [26] to identify unknown forces, so as to achieve the purpose of identifying unknown wind loads.
The contents of this paper are organized as follows. First, the improved identification method of wind loads is proposed in Section 2. The decomposition of wind loads decoupled sub-distributed dynamic loads with structural modes and time history functions for discrete system and continuous system, and modal Kalman filter with unknown inputs (MKF-UI) are presented. In Section 3, two numerical identifications of wind load on a 60 m continuous wind tower system and a 20-story discrete shear frame system are carried out respectively to verify the effectiveness of the proposed method.

Simulation of Wind Loads on Structures
In this paper, only one-dimensional wind loads along the vertical space field are considered. Then, the wind speed at level h denoted as v(h, t) at time t can be expressed as: where v(h) is the mean wind speed at level h; v(h, t) is the fluctuating wind speed at level h, which varies over time.
The mean wind speed can be calculated according to the power law [27] as follows: where h and v(h ) represent the reference height and the mean wind speed at the corresponding height according to the specification. In this study, h is set to 10 m above the ground; h and v(h) are the arbitrary height and its corresponding average wind speed; α is the power law exponent. In this paper, the power spectrum of fluctuating wind speed v(h, t) is davenport spectrum, then the self-power spectrum S ii (ω) can be expressed as follows [28]: where ω is the circular frequency, k is ground roughness, v 10 is the basic wind speed at a standard reference height of 10 m, taking into account the local ground roughness.
The cross spectrum of longitudinal fluctuating wind loads at different heights can be expressed as [29]: where S ii (ω) and S jj (ω) are the wind speed spectrum at h i and h j height respectively, and γ ij is the coherence function. The expression of the coherence function is as follows [29]: where C z is the exponential decay coefficient, and ∆h is the height difference between the i-th point and the j-th point. The spectral density function matrix can be expressed as: Perform a Cholesky decomposition of this matrix, that is: According to Shinozuka's theory [28,30], the pulsating wind speed can be simulated by: in which ∆ω = (ω u − ω k )/N, ω l = ω k + (l − 1 + j/N)∆ω, ω u and ω k are the upper and lower limits of the interception frequency, respectively; N is a sufficiently large positive integer and θ jl is a random number whose value range is (0, 2π). Based on Equation (1), the total wind pressure w(h, t) can be obtained: in which ρ is the air density, and w(h, t) can be divided into average wind pressure w(h) and the part of fluctuating wind pressure w(h, t) as follows: Then at level h, the wind load can be calculated as where µ s (h) is the structure drag coefficient at level h; A s (h) is the windward area at level h.

Wind Loads Decomposition and Modal Transformation
Since wind load changes with height, full observation is required to identify all unknown wind loads. However, it is almost impossible to deploy a lot of sensors in practical engineering. In this paper, the motion equation is transformed into modal coordinates, the wind loads are decomposed into a series of independent functions of space and time, and modal truncation is carried out to reduce the identification numbers and achieve a partial observation effect. In this paper, the following formulas are mainly derived for two cases: continuous system and discrete system.

Identification of Wind Loads on Continuous Structural Systems
A vertical cantilever structure can be used to verify the unknown wind loads identification method. It is assumed that the aerodynamic damping effect is not considered in this paper. The equation of motion is given as: where ρ and A(h) are the mass density of the cantilever structure and the structure crosssectional area at level h. Symbols EI(h) and c are the flexural rigidity and the damping of the cantilever at level h, respectively. f(h, t) are the distributed wind loads calculated in Equation (12). w(h, t) is the displacement of the cantilever at level h and time t.
Due to the limited number of measurement sensors, structural modal reduction is adopted according to the modal coordinate transformation theory [31]. The structural acceleration response vector ..
x(t) can be obtained as follows: ..

q(t)
is the p × 1 modal acceleration vector, Φ is the mode shape of continuous system structure after p order truncation and p is the number of modal truncations. Truncated modal number has an impact on the result, and many scholars [11,17] have done relevant research on it. It is pointed out that for structural vibration, the first several modes are mainly involved in the vibration of the structure. The more the number of intercepted modes, the higher the identification accuracy will be, but the more structural responses need to be measured. In order to intercept appropriate modal numbers and obtain relatively accurate identification results, this paper mainly uses the proper orthogonal decomposition (POD) technique to determine the number of modal truncations of a multi-degrees of freedom (MDOF) system. The following Equation can be used to calculate the energy contribution η: where p means the number of truncated modes and n is the number of structural degrees of freedom. λ i is the eigenvalues of acceleration response covariance matrix, which represent the proportion of energy in the response. In general, η should exceed 99% of the whole energy in the structural responses and its corresponding modal order is the number of truncations needed. After transforming physical coordinates into modal coordinates, Equation (13) can be rewritten as: ..
For Rayleigh damped system, some matrix parameters in the above formula can be expressed as: where ω i and ξ i mean the natural frequency and the modal damping ratio of the ith mode, and m i is the ith modal mass.
It is well known that the distributed wind loads are the function of position and time. This article is proposed to decompose the distributed wind loads f(h, t) by proper orthogonal decomposition (POD) as a series of time-space decoupled sub-distributed dynamic loads with independent basic spatial distribution functions T j (h) and time history functions d j (t). i.e., Then, the ith modal force of the distributed wind loads can be expressed as follows: where ϕ i (h) is the ith mode shape function of the vertical cantilever structure, L is the length of the continuous structure, and f m i (t) is the ith modal force for the distributed wind loads. If T j (h) and ϕ i (h) are orthogonality, Equation (18) can be easily solved. Herein, the orthogonal basis is defined as: Due to the orthogonality of the mode of vibration: where δ ij is the Kronecker delta and m i means the ith modal mass. Therefore, the expression of the ith modal force of the wind loads can be written as: After transforming Equation (22) into modal space, the motion equation is rewritten as: ..
Herein, d(t) can be obtained by the method [26] in Section 2.3. Therefore, after the estimation of d(t), the distributed loads can be reconstructed as:

Identification of Wind Loads on Discrete Structural Systems
In many cases, the modal shape function of complex continuous systems cannot be analytically analyzed, so structures need to be discretized into finite element elements.
The motion equation of n degrees of freedom discrete structure can be expressed as follows: M ..
where M is the mass matrices, C is the damping matrices, and K means the stiffness matrices; F(t) means the n × 1 external excitation vector; the structural acceleration, velocity, and displacement vector can be expressed as ..
x(t), and x(t). Similar to the continuous system, the motion equation in modal coordinates can be expressed as: where ..

q(t),
. q(t), q(t) are the p × 1 modal acceleration, velocity and displacement vector, Φ is the n × p mode shape, p is the number of modal truncation, M * = Φ T MΦ = I is the modal mass matrix, C * = Φ T CΦ = Γ is the modal damping matrix, K * = Φ T KΦ = Λ is the modal stiffness matrix, Equation (20) can be rewritten as follows: ..
Similar to Equation (17), the following can be obtained: Wind loads in discrete system can also be decomposed into a series of two functions of time and space: The jth modal force f m j (t) of the distributed wind load F(t) can be obtained as follows: where ϕ j is the jth modal shape. The basis function can be defined by using the orthogonality of the modes [15]: Then the Equation (30) can be converted to as follows: (32) in which m j is the jth modal mass. Substitute Equation (26) into Equation (21), which can be rewritten as: ..
Equation (33) is similar to Equation (23) in the continuous system. When d(t) is identified by MKF-UI, the distributed wind loads can be reconstructed as:

Wind Loads Identification Based on MKF-UI
When structural state vector X(t) is defined as X(t) = q(t) .

q(t)
, the state equation of the modal Equations (23) and (33) can be expressed as .

X(t) = AX(t) + Bd(t)
where: If samples are taken at equal intervals at each t = [t 0 , t 1 , t 2 , · · · t k , · · ·] and the excitation is assumed to remain constant within the time interval of ∆t = t k+1 − t k , the continuous state motion equation can be rewritten as the following discrete state equation: where X k+1 and X k mean the structural state vectors at time t = (k + 1)∆t and t = k∆t respectively. d k is the modal force of the excitation acting on the structure at time t = k∆t. w k is the zero mean model error matrix and its covariance matrix is Q k = E(w k w T k ). The correlation coefficient matrix of state and excitation from t = (k + 1)∆t to t = k∆t can be discretized as: In this paper, only part of the acceleration responses are observed, so Y k+1 denotes r × 1 acceleration responses under the observed physical coordinates, r is the number of response observation points. The measurement equation at time t = (k + 1)∆t can be given as: where v k+1 denotes zero mean Gaussian measurement noise vector with and covariance matrix is R k+1 = E(v k+1 v T k+1 ), the relevant observation matrix can be expressed as: where D a is the r × n position matrix based on the observed acceleration. Then, at time t = (k + 1)∆t, the predicted state estimate X k+1|k can be expressed as follows: whereX k|k andd k|k are the state estimates and the modal forces at time t = k∆t, respectively. Then the state estimation equation can be expressed as: where Y k+1 represents the observation vector in Equation (38), that is, the structural response data collected by the sensor installed on the structure. K k+1 is Kalman gain matrix, which is consistent with traditional Kalman filter and can be expressed as: In the above equation, P X k+1|k andP X k+1|k+1 are the covariance matrices of the prediction error and estimation error of the state vector respectively.
If the number of observed responses is greater than the unknown number of modal forces, then, using the least square method, the unknown modal forces can be identified as: The error covariance matrix of state vector and modal force P X k+1|k ,P X k+1|k+1 , P Xd k+1|k+1P dX k+1|k+1 ,P d k+1|k+1 can be expressed as: According to Equation (43),d k+1|k+1 can be identified, then through Equations (24) and (34), unknown wind loads at physical coordinates in continuous and discrete systems can be successfully identified respectively.

Identification of Wind Loads on Continuous Structure of Wind Tower
To verify the feasibility of this method in a continuous system, a 60 m-high wind turbine tower structure is established for numerical calculation. In practical engineering, similar high-rise structures have a nonlinear nature and the dynamic equation describing the vibration process of structures should be a nonlinear differential equation. However, in consideration of the low wind speed in this paper, such buildings basically only show the linear elastic response [32] and non-linear response is basically not reflected. Therefore, wind power towers are regarded as linear structures for wind load identification in this paper. At the same time, on the premise of not affecting the dynamic characteristics of the wind power tower structure, the auxiliary components that have little influence on it are ignored, and the complex parts such as connections are simplified. Figure 1 shows the schematic diagram of its structure.
The diameter of the bottom and top of the structure is 4 m and 3 m respectively, and the thickness is 25 mm and 10 mm. The diameter and thickness of the wind turbine tower linearly change from the bottom to the top. The structure is made of Q345 steel and its Young's modulus is E = 2.06 × 10 11 N/m 2 . The mass density is ρ = 7850 kg/m 3 . The first-order frequency of the wind tower structure is 1.16 Hz. It is assumed that aerodynamic damping is not considered in this numerical simulation.  Young's modulus is E = 2.06 × 10 11 N/m 2 . The mass density is ρ = 7850 kg/m order frequency of the wind tower structure is 1.16 Hz. It is assumed that a damping is not considered in this numerical simulation.
The structure is a cantilever beam with variable section, and its modal sh solved by Rayleigh-Ritz method. Based on the boundary conditions confor structure, the i th modal test function is defined as the following expression where L is the total height of the structure. Substitute the trial function into ing equation to obtain the corresponding stiffness and mass matrices s The natural frequency i ω and eigenvector i a of the structure can be s ing the eigen-equation , and the modal shape function The structure is a cantilever beam with variable section, and its modal shape curve is solved by Rayleigh-Ritz method. Based on the boundary conditions conforming to the structure, the ith modal test function is defined as the following expression [15]: where L is the total height of the structure. Substitute the trial function into the following equation to obtain the corresponding stiffness and mass matrices K s , M s : The natural frequency ω i and eigenvector a i of the structure can be solved by using the eigen-equation K s − ω 2 M s · a = 0, and the modal shape function of the structure can be obtained from the following expression: where p is the order of modal truncation. Since the harmonic superposition method cannot produce a continuous wind loads action curve, in this numerical example, based on the spatial correlation of wind loads height direction, it is assumed that the wind loads of adjacent 2 m is correlated, and the wind loads within the interval changes linearly. In this numerical simulation, the fluctuating wind speed is simulated by harmonic superposition method using the spectral density of Davenport spectrum. During the simulation, the time step is 0.2 s. The mean wind speed is calculated by exponential formula. Based on the Chinese National Load Code [33], the parameters α, h and v(h ) in Equation (2) are set to 0.22, 10 m and 10 m/s. The ground roughness is set to 0.0045. The parameter µ s in Equation (12) is set to be 1.6 and the air density is set to be 1.23 kg/m 3 . Figure 2 shows the pulsating wind speed at 30 m and 60 m obtained by harmonic superposition. action curve, in this numerical example, based on the spatial correlation of wind loads height direction, it is assumed that the wind loads of adjacent 2 m is correlated, and the wind loads within the interval changes linearly. In this numerical simulation, the fluctuating wind speed is simulated by harmonic superposition method using the spectral density of Davenport spectrum. During the simulation, the time step is 0.2 s. The mean wind speed is calculated by exponential formula. Based on the Chinese National Load Code [33], the parameters α, ' h and ( ') v h in Equation (2) are set to 0.22, 10 m and 10 m/s. The ground roughness is set to 0.0045. The parameter s μ in Equation (12) is set to be 1.6 and the air density is set to be 1.23 kg/m 3 . Figure 2 shows the pulsating wind speed at 30 m and 60 m obtained by harmonic superposition. The simulated wind loads are applied to the structure and the acceleration responses can be obtained through calculation. According to the calculation, the energy contribution rate of the first 15 orders of acceleration response of shear frame in this section has reached 99%. Therefore, the modal truncation order of wind tower structure is calculated as 15 order.
Because the number of sensors measuring points must be larger than the number of modal truncations, 15 acceleration observation points are arranged along the structure every four meters. The wind load results obtained by identification are compared with the real wind load obtained by original simulation. Figure 3 shows the comparison of timevarying wind load at different height. Figure 4 shows the comparison of wind load identification in the frequency domain. As can be seen from Figure 3, in the continuous system structure, unknown wind loads can be effectively identified by using the method in this paper only through partial observation. Table 1 shows the RMS error of wind load at different heights. Although it is not completely consistent, the identification error is in an acceptable range for wind loads identification. As can be seen from the table and figure, there is a large error in identifying the wind load at the low height. The analysis may be caused by a small average wind and a large proportion of fluctuating wind at the low height leading to a decrease in the identification accuracy. However, since the significance of identifying the wind load at the low height is relatively small, the proposed method is still effective in identifying the wind loads on the structure. The simulated wind loads are applied to the structure and the acceleration responses can be obtained through calculation. According to the calculation, the energy contribution rate of the first 15 orders of acceleration response of shear frame in this section has reached 99%. Therefore, the modal truncation order of wind tower structure is calculated as 15 order.
Because the number of sensors measuring points must be larger than the number of modal truncations, 15 acceleration observation points are arranged along the structure every four meters. The wind load results obtained by identification are compared with the real wind load obtained by original simulation. Figure 3 shows the comparison of time-varying wind load at different height. Figure 4 shows the comparison of wind load identification in the frequency domain. As can be seen from Figure 3, in the continuous system structure, unknown wind loads can be effectively identified by using the method in this paper only through partial observation. Table 1 shows the RMS error of wind load at different heights. Although it is not completely consistent, the identification error is in an acceptable range for wind loads identification. As can be seen from the table and figure, there is a large error in identifying the wind load at the low height. The analysis may be caused by a small average wind and a large proportion of fluctuating wind at the low height leading to a decrease in the identification accuracy. However, since the significance of identifying the wind load at the low height is relatively small, the proposed method is still effective in identifying the wind loads on the structure.

Identification of Wind Loads on a Discrete 20-Story Shear Frame
To verify the accuracy and availability of the above method on the identification of wind loads on discrete structural system, a 20-story shear frame is adopted as a demonstration example. At each floor, the mass coefficient m and the stiffness coefficient k are 50 t and 7.5 × 10 6 N/m, respectively. The damping is assumed to be Rayleigh damping and the damping ratio of the building is set to be 0.47%. The story height is 4.8 m for each floor.
In this numerical simulation, the wind load simulation parameters are consistent with the first calculation example. Figure 5 shows the simulated time history diagrams of fluctuating wind speed at different height. The parameters µ s and A s in Equation (12) are set to be 1.6 and 24 m 2 of each floor. The air density is set the same as the previous example.
The simulated wind loads are applied to the structure and the acceleration responses can be obtained through calculation. According to the calculation, the energy contribution rate of the first seven orders of acceleration response of shear frame in this section has reached 99%, so the first seven orders of modes are selected for identification. Since the number of intercepted modes is seven, the number of observed responses should be greater than or equal to seven. Therefore, acceleration observation on the 2nd, 4th, 7th, 10th, 13th and 16th floors are selected herein.

Identification of Wind Loads on a Discrete 20-Story Shear Frame
To verify the accuracy and availability of the above method on the identification of wind loads on discrete structural system, a 20-story shear frame is adopted as a demonstration example. At each floor, the mass coefficient m and the stiffness coefficient k are 50 t and 7.5 × 10 6 N/m, respectively. The damping is assumed to be Rayleigh damping and the damping ratio of the building is set to be 0.47%. The story height is 4.8 m for each floor.
In this numerical simulation, the wind load simulation parameters are consistent with the first calculation example. Figure 5 shows the simulated time history diagrams of fluctuating wind speed at different height. The parameters s μ and s A in Equation (12) are set to be 1.6 and 24 m 2 of each floor. The air density is set the same as the previous example. The simulated wind loads are applied to the structure and the acceleration responses can be obtained through calculation. According to the calculation, the energy contribution rate of the first seven orders of acceleration response of shear frame in this section has reached 99%, so the first seven orders of modes are selected for identification. Since the Then, the unknown force can be successfully identified by the proposed method introduced in Section 2 using only seven acceleration responses. In Figures 6 and 7, it shows the comparison of simulation with inversion results the of inverse problem in the time domain and frequency domain respectively. It can be seen from the figure that the wind load identification results of other layers are in good agreement with the wind load obtained by harmonic superposition method, except for the large error of the bottom layer. Table 2 shows the RMS error of wind load of each floor.
The simulated wind loads are applied to the structure and the acceleration responses can be obtained through calculation. According to the calculation, the energy contribution rate of the first seven orders of acceleration response of shear frame in this section has reached 99%, so the first seven orders of modes are selected for identification. Since the number of intercepted modes is seven, the number of observed responses should be greater than or equal to seven. Therefore, acceleration observation on the 2nd, 4th, 7th, 10th, 13th and 16th floors are selected herein.
Then, the unknown force can be successfully identified by the proposed method introduced in Section 2 using only seven acceleration responses. In Figures 6 and 7, it shows the comparison of simulation with inversion results the of inverse problem in the time domain and frequency domain respectively. It can be seen from the figure that the wind load identification results of other layers are in good agreement with the wind load obtained by harmonic superposition method, except for the large error of the bottom layer. Table 2 shows the RMS error of wind load of each floor.

Conclusions
In this paper, a method is presented for the identification of unknown wind loads on structures only using partial observations of structural acceleration responses. The method is based on the decomposition of wind load as a series of time-space decoupled spatial distribution functions and time history functions. Structural modes are adopted as the spatial distribution functions and structural modal truncation can be used to reduce the numbers of unknown time history functions, which can be identified by the method of MKF-UI. The proposed method can be used for the identification of wind loads on both continuous and discrete system structures, and only partial structural acceleration responses which are more than the modal truncation modes are required, which is applicable for actual engineering situations. Besides, it is necessary to conduct more numerical and experimental validations of the proposed method for the identification of wind loads on complex structures. In addition, the aerodynamic damping is not considered in this paper, which is very important. Wind load identification related to aerodynamic damping is also being investigated.