# A Method of Modal Parameter Identification for Wind Turbine Blade Based on Binocular Dynamic Photogrammetry

Authors: Wenyun Wang; Xuejun Li; Anhua Chen
Journal: Shock and Vibration
 (2019)

---

## Abstract

The identification of operational modal parameters of a wind turbine blade is fundamental for online damage detection. In this paper, we use binocular photogrammetry technology instead of traditional contact sensors to measure the vibration of blade and apply the advanced stochastic system identification technique to identify the blade modal frequencies automatically when only output data are available. Image feature extraction and target point tracking (PT) are carried out to acquire the displacement of labeled targets on the wind turbine blade. The vibration responses of the target points are obtained. The data-driven stochastic subspace identification (SSI-Data) method based on the Kalman filter prediction sequence is explored to extract modal parameters from vibration response under unknown excitation. Hankel matrixes are reconstructed with different dimensions, so different modal parameters are produced. Similarity of these modal parameters is compared and used to cluster modes into groups. Under appropriate tolerance thresholds, spurious modes can be eliminated. Experiment results show that good effects and stable accuracy can also be achieved with the presented photogrammetry vibration measurement and automatic modal identification algorithm.


---

## Body

## 1. Introduction

Modal analysis of a wind turbine blade is a key step for analysis and testing of the dynamic performance of blade structure and is also fundamental for blade design and manufacturing [1]. Estimating modal properties, i.e., natural frequencies, mode shapes, and damping ratios, using vibration measurements have been extensively studied. Operational modal analysis means to analyze the modal parameters only based on structural dynamic response, without knowing excitation signals. New methods which relaxed the requirement of input measurement started to be developed and were found to be useful especially for large-scale applications (e.g., civil structures) where input generation is impractical. This kind of method can reveal dynamic performance of the structure in real conditions without interruption and easy to be implemented [2]. Vibration measurement is essential in providing data for modal analysis. With the gradual increase of the power of a wind turbine, the blade size is getting larger and larger. The traditional contact method is limited by the number of sensors, the acquisition channels, and the difficulty of installation and cannot meet the requirements of measuring for large-scale blades. Photogrammetry can achieve dynamic positioning and measurement for spatial objects by taking photos of the measured object through multiple cameras [3]. This technology is noncontact, full field measuring and especially applicable for light, flexible structures with inconvenient sensor installation, which has no adding weights and no affection to the dynamics of measured objects [4–6]. Yang et al. [7] pasted artificial signs on the blades of large wind turbines and used photogrammetric technology to analyze their deformation. Jurjo et al. [8] marked the film surface and used photogrammetry to detect the deformation of the film structure under different tensions. Kim and Kim [9] achieved the dynamic deformation of the cable of Gwangan Bridge through image processing with videos captured by using a portable digital camera. In order to acquire dynamic displacements of structures in 3D space, it is important to effectively identify features and track points through continuous object matching.On the contrary, the stochastic subspace identification (SSI) is a classical time-domain method to identify modal parameters for a linear system directly from structural vibration response data [10]. The calculated poles due to overestimation of the model order and the noise poles introduced by the measurement noise will greatly affect the determination of the final modal parameters [11]. Higher requirements on the relevant professional skills of operators were asked in the method of artificial determination of modal order and real modes. Reynders et al. [12] and Bakir [13] obtained accurate stabilization diagram by introducing criterions such as modal consistency indicator, modal phase collinearity, and modal contribution to distinguish physical poles from spurious poles. Li et al. [14] introduced a method based on the eigen-decomposition of the state matrix constructed in SSI-Data or autoregressive model. Xu et al. [15] used a digital camera to measure the displacement vibration of an antenna and applied the eigensystem realization algorithm to identify the modal parameters, which were susceptible to background noises. Based on video measurement, Wang et al. [16] realized geometric modeling of a light fan blade and some leaves, and then modal analysis with PolyMax method provided by commercial software was carried out.From the analysis of pole discrimination above, we find that it is complex and time-consuming for the state decomposition method or inaccurate and limited scenarios with certain criterions in the stabilization diagram method and still need high skilled manipulations. In this paper, a new method of dynamics detection and automatic modal parameter identification is proposed. Photogrammetry is used to measure marked points vibration, and modal analysis is carried out with data-driven stochastic subspace identification (SSI-Data). Considering the impact of noise on vibration data, spurious modes are discriminated through reconstructing Hankel matrix multiple times and fuzzy clustering method. The modal similarity is calculated and differentiated by a new index based on tolerance and modal confidence. The modal parameters are clustered and automatically identified with statistical analysis of results from dynamic vibration sequences.


## 2. Principle of Binocular Photogrammetry

Figure1 shows the imaging model of binocular vision. Subscripts 1 and 2 represent the left and right cameras, respectively. A three-dimensional world coordinate system OW−XWYWZW coincides with the right camera coordinate system. The coordinates of a point P in the three-dimensional world coordinate system are PXW,YW,ZW.Figure 1
Stereo imaging model of binocular vision.Taking the left camera as an example, the camera coordinate system isOC1−XC1YC1ZC1. The origin OC1 is defined as the optical center of the camera, and ZC1 axis coincides with the optical axis. u1 and v1 axes of the computer pixel coordinate system o1−u1v1 are parallel to the XC1andYC1 axes. The abscissa and ordinate values represent the columns and rows of the image, respectively. The image plane coordinate system is O1−X1Y1. The ideal principal point O1 lies at the center of the image plane. The ideal image point of the point P in the left image plane is P1xp1,yp1. For high accuracy measuring, the processing and assembly errors of the camera and optical lens cannot be omitted, then the actual image point is expressed as P1′xp1′,yp1′, and the actual principal point is expressed as o1′u1′,v1′, respectively.The right camera coordinate system isOC2−XC2YC2ZC2. If it rotates around XC2,YC2,andZC2 axes with angles ε,ψ,andκ in turn, the new formed coordinate system O−XYZ will be parallel to the left camera coordinate system OC1−XC1YC1ZC1, where R is the rotation matrix as in equation (1). The origin OC1 can be translated to OC2 with TtX,tY,tZ:(1)R=r1r2r3r4r5r6r7r8r9=cosψcosk−cosψsinksinψsinεsinψcosk+cosεsink−sinεsinψsink+coskcosε−cosψsinε−sinψcosεcosk+sinεsinksinψsinkcosε+cosksinεcosψcosε.Transformation of the 3D world coordinate system to the left camera coordinate system can be described using 3 × 3 orthogonal matrixR and 3 × 1 vector T as follows:(2)XC1YC1ZC1=R⋅XC2YC2ZC2+T.Equation (2) can be rewritten in the matrix form by using homogeneous styles:(3)XC1YC1ZC1=r1r2r3txr4r5r6tyr7r8r9tz⋅XC2YC2ZC21.According to the mapping relationship between computer pixel coordinates and three-dimensional world coordinates, the world coordinates of space pointP can be calculated as(4)XW=ZWX1f1,YW=ZWY1f1,ZW=f1f2tx−X2tzX2r7X1+r8Y1+f1r9−f2r1X1+r2Y1+f1r3=f1f2ty−X2tzY2r7X1+r8Y1+f1r9−f2r4X1+r5Y1+f1r6,where f1 and f2 are the focal lengths of the left and right cameras and X1,Y1 and X2,Y2 are the coordinates of the target points in the plane coordinates of the left and right images.The object movements are captured by using a binocular vision system. And the dynamic images are treated with segmentation and identification. The authors developed a kind of tracking and matching technology for coded targets pasted on the blade [17]. The cameras are calibrated with high accuracy by bundle alignment. An adaptive Wiener filter and threshold segmentation were explored to effectively distinguish independent mark points. Three-dimensional coordinates with subpixel precision are obtained by homogeneous correspondence. The results are normalized to achieve continuous displacement. So the vibration data of measured object are acquired which will be used for modal analysis in Section 3.


## 3. Automatic Identification of Modal Parameters

### 3.1. Data-Driven Stochastic Subspace Identification

After the measurement of binocular vision, the vibration data are used to identify the operational modal parameters. The measured vibration is discrete time series data. Here, we assume that the excitation is white noise with zero mean value (if the input contains dominant frequency components, these frequencies will appear as poles of the equation (5)). Lardiès [18] discussed modal identification methods and indicted that parametric approach outperforms the nonparametric approaches. Here, we use the parametric approach, data-driven stochastic subspace identification (SSI-Data), to extract the modal parameters from the time-domain data.The discrete SSI model of the structure can be described as follows [12]:(5)xt+1=Axt+wt,yt=Cxt+vt,where t represents the index of sampled data and xt∈R2n×1 and yt∈Rk×1 are the state vector and output vector corresponding to t, respectively. wt and vt are white noises with mean 0, A is the state space matrix, and C is the output matrix.The vibration response data of multiple locations are gathered in a block Hankel matrixH. Hankel matrix H can be divided into a past data matrix Yp∈Rka×b and future part matrix Yf∈Rka×b, which are defined as(6)Yp=1by0y1y2⋯yb−1y1y2y3⋯yby2y3y4⋯yb+1⋮⋮⋮⋮⋮ya−1yaya+1⋯ya+b−1,(7)Yf=1byaya+1ya+2⋯ya+bya+1ya+2ya+3⋯ya+b+1ya+2ya+3ya+4⋯ya+b+2⋮⋮⋮⋮⋮y2a−1y2ay2a+1⋯y2a+b−1,where yi=yi1,yi2,…,yikT∈Rk×1, k is the number of measuring points, and a and b represent the number of rows and columns of the Hankel matrix constructed by single point response. For statistical reasons, it is assumed that b⟶∞. Typically, b>20a is generally chosen [13].By taking the QR-factorization of the Hankel matrix which consisted ofYp and Yf, we can have(8)H=YpYf=RQT=R1100R21R220R31R32R33Q1TQ2TQ3T,where R∈R2ka×b is a lower triangular matrix; Q∈Rb×b is an orthogonal matrix; R11=R1:ka,1:ka includes the first ka rows and ka columns of R; and R21=Rka+1:ka+k,1:ka includes from ka+1 to ka+k rows of R and 1 to ka columns of R; R22,R31,R32,andR33 can be determined similarly. Q1=Q1:ka,: includes the first ka rows of Q; Q2andQ3 can be determined similarly.From equation (8), we have(9)Yp=R11Q1T,Yf=R21R31Q1T+R220R32R33Q2TQ3T.The projection of row space ofYp into row space of Yf is defined as χa∈Rka×b, which is(10)χa=YfYp=YfYpTYpYpT+Yp=R21R31Q1T.By applying singular value decomposition (SVD) to the projection matrix,(11)χa=U1S1V1T,where U1∈Rka×ka, V1∈Rb×b are the unitary matrices composed of left and right singular vectors corresponding to nonzero singular values, and S1 is a diagonal matrix with singular values arranged in the descending order.The projection matrix can be factorized as the product of the observability matrix and the Kalman filter sequence:(12)χa=OaX^a.Combining equations (11) and (12) gives(13)Oa=U1S11/2,X^a=Oa−1χa.The matrixesYf− and Yp+ can be obtained by moving the first k rows of Yf matrix to Yp. According to equations (8) and (9),(14)Yp+=R110R21R22Q1TQ2T,Yf−=R31R32R33Q1TQ2TQ3T.If the row space ofYp+ is projected to Yf−, the projection matrix χa−1 can be obtained:(15)χa−1=Yf−Yp+=Yf−Yp+TYp+Yp+T+Yp+=R31R32Q1TQ2T.And we have(16)X^a+1=Oa−1−1χa−1.Defining matrixYaa as the rows from ka+1 to ka+1 of H matrix, it is expressed in the form QR-factorization:(17)Yaa=R21R22Q1TQ2T.System matrixes A and C can be obtained by solving the set of equations in least squares sense.(18)AC=X^a+1YaaX^a−1.System state matrixA is factorized to produce the eigen values:(19)A=ψΛψ−1,where Λ=diagλs,s=1,…,n. If λsc is a continuous-time eigen value, then there is(20)λsc=−ξsωs+jωs1−ξs2,where ξs is the damping ratio and ωs is the circular frequency. They can be computed from the following equation:(21)fs=ωs2π=λsc2π,ξs=−Reλscωs,where Redenotes the real part and the corresponding modal shape Φ is(22)Φ=Cψ.


### 3.2. Spurious Pole Discrimination Method

At the process of the modal parameters identification, the maximal model orderN of the state space model should be large enough to conserve all possible practical poles, and thus a lot of spurious poles may be introduced during the calculation. Because of the measurement noise exists, if we choose different dimensions (a, b) for the Hankel matrix, the poles obtained will also be quite different. The modal parameters are inherent properties of the structure. Ideally, they do not change with the external environment, so the physical poles remain stable, while the spurious poles are scattered from time to time. Based on this, we present a spurious pole discrimination method by reconstructing Hankel matrix many times followed with response sequences. Each time (i = 1, …, m), we choose different (a, b) and construct the Hankel matrix Hi. A set of poles Pij will be obtained from every Hi by using SSI-Data process. Then the N poles Pijj=1,…,N are clustered to produce N clusters, and the center of each cluster gives a possible pole. The possible poles which have nearest natural frequency will be regrouped and each group Gkk=1,…,N has m poles gkii=1,…,m. By introducing tolerance and modal confidence, the poles are judged and discriminated. The similarity between these m poles is defined as(23)rk=Wfdfσfkimaxfki+WΦdΦ1−∑MACΦki,Φ¯kim,where Wf and WΦ are the weights of natural frequencies and mode shapes, df and dΦ are the tolerances of natural frequencies and mode shapes, fki and σfki are the natural frequency and its variance value, and MACΦki,Φ¯ki is the modal confidence values of modes Φki and the average value Φ¯ki. If rk≤1, then the pole is a stable pole. On the contrary, it is eliminated. Here, we do not consider the similarity of damping ratio because the large flexible structure uses a lot of composite materials and glass fiber and has a great uncertainty of the damping ratio. According to the criterion of pole similarity, the results after eliminating spurious poles are analyzed by spectral clustering. The average value of physical modal parameters for each group is set as the final identified modal parameters.


## 3.1. Data-Driven Stochastic Subspace Identification

After the measurement of binocular vision, the vibration data are used to identify the operational modal parameters. The measured vibration is discrete time series data. Here, we assume that the excitation is white noise with zero mean value (if the input contains dominant frequency components, these frequencies will appear as poles of the equation (5)). Lardiès [18] discussed modal identification methods and indicted that parametric approach outperforms the nonparametric approaches. Here, we use the parametric approach, data-driven stochastic subspace identification (SSI-Data), to extract the modal parameters from the time-domain data.The discrete SSI model of the structure can be described as follows [12]:(5)xt+1=Axt+wt,yt=Cxt+vt,where t represents the index of sampled data and xt∈R2n×1 and yt∈Rk×1 are the state vector and output vector corresponding to t, respectively. wt and vt are white noises with mean 0, A is the state space matrix, and C is the output matrix.The vibration response data of multiple locations are gathered in a block Hankel matrixH. Hankel matrix H can be divided into a past data matrix Yp∈Rka×b and future part matrix Yf∈Rka×b, which are defined as(6)Yp=1by0y1y2⋯yb−1y1y2y3⋯yby2y3y4⋯yb+1⋮⋮⋮⋮⋮ya−1yaya+1⋯ya+b−1,(7)Yf=1byaya+1ya+2⋯ya+bya+1ya+2ya+3⋯ya+b+1ya+2ya+3ya+4⋯ya+b+2⋮⋮⋮⋮⋮y2a−1y2ay2a+1⋯y2a+b−1,where yi=yi1,yi2,…,yikT∈Rk×1, k is the number of measuring points, and a and b represent the number of rows and columns of the Hankel matrix constructed by single point response. For statistical reasons, it is assumed that b⟶∞. Typically, b>20a is generally chosen [13].By taking the QR-factorization of the Hankel matrix which consisted ofYp and Yf, we can have(8)H=YpYf=RQT=R1100R21R220R31R32R33Q1TQ2TQ3T,where R∈R2ka×b is a lower triangular matrix; Q∈Rb×b is an orthogonal matrix; R11=R1:ka,1:ka includes the first ka rows and ka columns of R; and R21=Rka+1:ka+k,1:ka includes from ka+1 to ka+k rows of R and 1 to ka columns of R; R22,R31,R32,andR33 can be determined similarly. Q1=Q1:ka,: includes the first ka rows of Q; Q2andQ3 can be determined similarly.From equation (8), we have(9)Yp=R11Q1T,Yf=R21R31Q1T+R220R32R33Q2TQ3T.The projection of row space ofYp into row space of Yf is defined as χa∈Rka×b, which is(10)χa=YfYp=YfYpTYpYpT+Yp=R21R31Q1T.By applying singular value decomposition (SVD) to the projection matrix,(11)χa=U1S1V1T,where U1∈Rka×ka, V1∈Rb×b are the unitary matrices composed of left and right singular vectors corresponding to nonzero singular values, and S1 is a diagonal matrix with singular values arranged in the descending order.The projection matrix can be factorized as the product of the observability matrix and the Kalman filter sequence:(12)χa=OaX^a.Combining equations (11) and (12) gives(13)Oa=U1S11/2,X^a=Oa−1χa.The matrixesYf− and Yp+ can be obtained by moving the first k rows of Yf matrix to Yp. According to equations (8) and (9),(14)Yp+=R110R21R22Q1TQ2T,Yf−=R31R32R33Q1TQ2TQ3T.If the row space ofYp+ is projected to Yf−, the projection matrix χa−1 can be obtained:(15)χa−1=Yf−Yp+=Yf−Yp+TYp+Yp+T+Yp+=R31R32Q1TQ2T.And we have(16)X^a+1=Oa−1−1χa−1.Defining matrixYaa as the rows from ka+1 to ka+1 of H matrix, it is expressed in the form QR-factorization:(17)Yaa=R21R22Q1TQ2T.System matrixes A and C can be obtained by solving the set of equations in least squares sense.(18)AC=X^a+1YaaX^a−1.System state matrixA is factorized to produce the eigen values:(19)A=ψΛψ−1,where Λ=diagλs,s=1,…,n. If λsc is a continuous-time eigen value, then there is(20)λsc=−ξsωs+jωs1−ξs2,where ξs is the damping ratio and ωs is the circular frequency. They can be computed from the following equation:(21)fs=ωs2π=λsc2π,ξs=−Reλscωs,where Redenotes the real part and the corresponding modal shape Φ is(22)Φ=Cψ.


## 3.2. Spurious Pole Discrimination Method

At the process of the modal parameters identification, the maximal model orderN of the state space model should be large enough to conserve all possible practical poles, and thus a lot of spurious poles may be introduced during the calculation. Because of the measurement noise exists, if we choose different dimensions (a, b) for the Hankel matrix, the poles obtained will also be quite different. The modal parameters are inherent properties of the structure. Ideally, they do not change with the external environment, so the physical poles remain stable, while the spurious poles are scattered from time to time. Based on this, we present a spurious pole discrimination method by reconstructing Hankel matrix many times followed with response sequences. Each time (i = 1, …, m), we choose different (a, b) and construct the Hankel matrix Hi. A set of poles Pij will be obtained from every Hi by using SSI-Data process. Then the N poles Pijj=1,…,N are clustered to produce N clusters, and the center of each cluster gives a possible pole. The possible poles which have nearest natural frequency will be regrouped and each group Gkk=1,…,N has m poles gkii=1,…,m. By introducing tolerance and modal confidence, the poles are judged and discriminated. The similarity between these m poles is defined as(23)rk=Wfdfσfkimaxfki+WΦdΦ1−∑MACΦki,Φ¯kim,where Wf and WΦ are the weights of natural frequencies and mode shapes, df and dΦ are the tolerances of natural frequencies and mode shapes, fki and σfki are the natural frequency and its variance value, and MACΦki,Φ¯ki is the modal confidence values of modes Φki and the average value Φ¯ki. If rk≤1, then the pole is a stable pole. On the contrary, it is eliminated. Here, we do not consider the similarity of damping ratio because the large flexible structure uses a lot of composite materials and glass fiber and has a great uncertainty of the damping ratio. According to the criterion of pole similarity, the results after eliminating spurious poles are analyzed by spectral clustering. The average value of physical modal parameters for each group is set as the final identified modal parameters.


## 4. Measurement Experiments and Results

### 4.1. Vibration Response Measurement Analysis

In order to verify the validity and reliability of our method, the vibration measurement experiment of a 3 kW wind turbine blade model was carried out, and the modal parameters were identified based on the measured vibration response data. Two Dahua A5131MU210 industrial cameras were selected to form a stereo camera station (Figure2). Figure 3 shows the blade image of labeled targets on the surface taken by using the cameras. There are 10 pairs of homogeneous points. The blade plane is defined as the XY plane, and the direction perpendicular to the blade plane is the Z direction. Two kinds of inputs, continuous random excitation, and instantaneous impact are exerted to excite the blade along the Z direction. The sampling frequency of the cameras is 120 Hz, and 1000 continuous images are recorded.Figure 2
Experiment configuration.Figure 3
Blade images taken by using the cameras: (a) left camera and (b) right camera.
(a)
(b)The coordinates of the targets on the image plane can be extracted from the image. The three-dimensional coordinates of the object can be obtained by solving the correspondence problem with the point tracking (PT). Then, the spatial coordinate changes of the same target point can be calculated through 1000 continuous images. The vibration response of the blade can be obtained by normalizing displacement changes. Figures4–7 show the measurement results of vibration response of points 1 and 10 in three directions of X, Y, and Z under random excitation and impact excitation, respectively. Trajectories of 10 homogeneous points are shown in Figures 8(a) and 8(b). From these figures, we can find that all the measured points exhibit damped vibration in the three directions of X, Y, and Z under impact excitation, and the amplitude in the Z direction is much larger than that in the X and Y directions under both excitations.Figure 4
Vibration displacement curves of point 1 under random excitation: (a)X direction, (b) Y direction, and (c) Z direction.
(a)
(b)
(c)Figure 5
Vibration displacement curves of point 10 under random excitation: (a)X direction, (b) Y direction, and (c) Z direction.
(a)
(b)
(c)Figure 6
Vibration displacement curves of point 1 under impact excitation: (a)X direction, (b) Y direction, and (c) Z direction.
(a)
(b)
(c)Figure 7
Vibration displacement curves of point 10 under impact excitation: (a)X direction, (b) Y direction, and (c) Z direction.
(a)
(b)
(c)Figure 8
3D space motion trajectories: (a) random excitation and (b) impact excitation.
(a)
(b)


### 4.2. Analysis of Modal Identification Results

The SSI-Data is performed to identify the modal parameters of the blade structure from photogrammetric response data. A Hankel matrixH1 with 20 rows and 900 columns is constructed. The maximum model order is set to 100. Stabilization diagrams are drawn with frequency (abscissa coordinate) and the model order (ordinate coordinate). The results in X, Y, and Z directions under random excitation are shown in Figures 9(a)–9(c), while diagrams under impact excitation in Figures 9(d)–9(f). There are a lot of spurious poles scattered in the figures. By comparing diagrams in Figure 9, we can see that random excitation is better than impact excitation and can excite more order modes. The vibration energy in the Z direction is larger than the other directions, and the first three modes under the two excitation types are clearly visible. Keeping other conditions unchanged, we reconstructed Hankel matrix H2 with row number 30 and column number 850 and drew stabilization diagrams in the Z direction under random excitation (Figure 10). The row and column numbers are changed many times to produce groups of possible poles. By using the spurious pole discrimination method in Section 3.2, the spurious poles are eliminated. Figure 11 shows the diagrams after elimination of spurious poles. It can be seen that a stabilization diagram with better recognition effect can be obtained after applying the proposed spurious poles discrimination method.Figure 9
Stabilization diagram based on SSI-Data:X, Y, and Z directions under (a–c) random excitation and (d–f) impact excitation.
(a)
(b)
(c)
(d)
(e)
(f)Figure 10
Stabilization diagram (H2).Figure 11
Stabilization diagram after eliminating spurious poles.After eliminating spurious poles, we can further extract more accurate modal parameters by spectral clustering. The frequency is taken as abscissa coordinate, and the damping ratio is the ordinate coordinate; we have the clustering results as shown in Figure12. The final modal parameters are set as the mean values of these clustered samples. The first three clusters with most samples of the poles are coincidently complying with the first three physical modes of the blade from simulation. In order to verify the accuracy of this method, a modal test of the same blade with the force hammer is carried out, and Figure 13 shows the result. A hammer (B&K 8206-002), an acceleration sensor (B&K 4514B-001), and an acquisition card (B&K 3053-B-120) are used. The modal parameters are identified by the MTC hammer module and reflex postprocessing module of PULSE software. Table 1 lists the lowest three frequencies of the blade identified with our method and hammer test. They are quite consistent in Table 1, and the relative errors of the three frequencies are less than 3.2%. The results show that the vibration response data measured by the proposed method based on photogrammetric point tracking are reliable; the modes can be identified by the improved data-driven SSI method, and the accuracy are improved by spectral clustering.Figure 12
Spectral clustering.Figure 13
Stabilization diagram by the force hammer method.Table 1
The recognition results of blade’s modal frequencies.
Mode
Measurement method
Frequency (Hz)
Proposed method
Impact test
Relative error (%)
First order
4.42
4.28
3.2
Second order
11.4
11.76
−3.1
Third order
28.23
29.01
−2.7


## 4.1. Vibration Response Measurement Analysis

In order to verify the validity and reliability of our method, the vibration measurement experiment of a 3 kW wind turbine blade model was carried out, and the modal parameters were identified based on the measured vibration response data. Two Dahua A5131MU210 industrial cameras were selected to form a stereo camera station (Figure2). Figure 3 shows the blade image of labeled targets on the surface taken by using the cameras. There are 10 pairs of homogeneous points. The blade plane is defined as the XY plane, and the direction perpendicular to the blade plane is the Z direction. Two kinds of inputs, continuous random excitation, and instantaneous impact are exerted to excite the blade along the Z direction. The sampling frequency of the cameras is 120 Hz, and 1000 continuous images are recorded.Figure 2
Experiment configuration.Figure 3
Blade images taken by using the cameras: (a) left camera and (b) right camera.
(a)
(b)The coordinates of the targets on the image plane can be extracted from the image. The three-dimensional coordinates of the object can be obtained by solving the correspondence problem with the point tracking (PT). Then, the spatial coordinate changes of the same target point can be calculated through 1000 continuous images. The vibration response of the blade can be obtained by normalizing displacement changes. Figures4–7 show the measurement results of vibration response of points 1 and 10 in three directions of X, Y, and Z under random excitation and impact excitation, respectively. Trajectories of 10 homogeneous points are shown in Figures 8(a) and 8(b). From these figures, we can find that all the measured points exhibit damped vibration in the three directions of X, Y, and Z under impact excitation, and the amplitude in the Z direction is much larger than that in the X and Y directions under both excitations.Figure 4
Vibration displacement curves of point 1 under random excitation: (a)X direction, (b) Y direction, and (c) Z direction.
(a)
(b)
(c)Figure 5
Vibration displacement curves of point 10 under random excitation: (a)X direction, (b) Y direction, and (c) Z direction.
(a)
(b)
(c)Figure 6
Vibration displacement curves of point 1 under impact excitation: (a)X direction, (b) Y direction, and (c) Z direction.
(a)
(b)
(c)Figure 7
Vibration displacement curves of point 10 under impact excitation: (a)X direction, (b) Y direction, and (c) Z direction.
(a)
(b)
(c)Figure 8
3D space motion trajectories: (a) random excitation and (b) impact excitation.
(a)
(b)


## 4.2. Analysis of Modal Identification Results

The SSI-Data is performed to identify the modal parameters of the blade structure from photogrammetric response data. A Hankel matrixH1 with 20 rows and 900 columns is constructed. The maximum model order is set to 100. Stabilization diagrams are drawn with frequency (abscissa coordinate) and the model order (ordinate coordinate). The results in X, Y, and Z directions under random excitation are shown in Figures 9(a)–9(c), while diagrams under impact excitation in Figures 9(d)–9(f). There are a lot of spurious poles scattered in the figures. By comparing diagrams in Figure 9, we can see that random excitation is better than impact excitation and can excite more order modes. The vibration energy in the Z direction is larger than the other directions, and the first three modes under the two excitation types are clearly visible. Keeping other conditions unchanged, we reconstructed Hankel matrix H2 with row number 30 and column number 850 and drew stabilization diagrams in the Z direction under random excitation (Figure 10). The row and column numbers are changed many times to produce groups of possible poles. By using the spurious pole discrimination method in Section 3.2, the spurious poles are eliminated. Figure 11 shows the diagrams after elimination of spurious poles. It can be seen that a stabilization diagram with better recognition effect can be obtained after applying the proposed spurious poles discrimination method.Figure 9
Stabilization diagram based on SSI-Data:X, Y, and Z directions under (a–c) random excitation and (d–f) impact excitation.
(a)
(b)
(c)
(d)
(e)
(f)Figure 10
Stabilization diagram (H2).Figure 11
Stabilization diagram after eliminating spurious poles.After eliminating spurious poles, we can further extract more accurate modal parameters by spectral clustering. The frequency is taken as abscissa coordinate, and the damping ratio is the ordinate coordinate; we have the clustering results as shown in Figure12. The final modal parameters are set as the mean values of these clustered samples. The first three clusters with most samples of the poles are coincidently complying with the first three physical modes of the blade from simulation. In order to verify the accuracy of this method, a modal test of the same blade with the force hammer is carried out, and Figure 13 shows the result. A hammer (B&K 8206-002), an acceleration sensor (B&K 4514B-001), and an acquisition card (B&K 3053-B-120) are used. The modal parameters are identified by the MTC hammer module and reflex postprocessing module of PULSE software. Table 1 lists the lowest three frequencies of the blade identified with our method and hammer test. They are quite consistent in Table 1, and the relative errors of the three frequencies are less than 3.2%. The results show that the vibration response data measured by the proposed method based on photogrammetric point tracking are reliable; the modes can be identified by the improved data-driven SSI method, and the accuracy are improved by spectral clustering.Figure 12
Spectral clustering.Figure 13
Stabilization diagram by the force hammer method.Table 1
The recognition results of blade’s modal frequencies.
Mode
Measurement method
Frequency (Hz)
Proposed method
Impact test
Relative error (%)
First order
4.42
4.28
3.2
Second order
11.4
11.76
−3.1
Third order
28.23
29.01
−2.7


## 5. Conclusions

Traditional contact measurement has many demerits for large-scale structures. A vibration detection method based on photogrammetry and modal parameters identification with improved SSI-Data is presented in this paper. The images of marked object points are continuously taken with binocular vision system and processed. Point tracking and matching are carried out to obtain the vibration response of marked points. The Hankel matrix is constructed with different dimensions for many times to produce different results of modal parameters. Similarity of these modal parameters is computed and clustered to form several stabilized groups. Experiments of a wind turbine blade are performed. The vibration of the blade tip is detected, and the modal parameters are identified with the presented algorithm. The modal frequencies are quite consistent with results of hammer pulse impact testing. It illustrates that photogrammetric vibration measurement and the new SSI-Data method are simple, efficient, and feasible for analyzing the dynamic characteristics of a large-scale irregular structure in applications.



---

*Source: 7610930-2019-04-15.xml*