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

#e 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. #e vibration responses of the target points are obtained. #e datadriven 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.


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 largescale applications (e.g., civil structures) where input generation is impractical.is 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.e 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].
is 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][5][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].
e 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.e modal similarity is calculated and differentiated by a new index based on tolerance and modal confidence.e modal parameters are clustered and automatically identified with statistical analysis of results from dynamic vibration sequences.

Principle of Binocular Photogrammetry
Figure 1 shows the imaging model of binocular vision.Subscripts 1 and 2 represent the left and right cameras, respectively.A three-dimensional world coordinate system (O W − X W Y W Z W ) coincides with the right camera coordinate system.e coordinates of a point P in the three-dimensional world coordinate system are P(X W , Y W , Z W ).
Taking the left camera as an example, the camera coordinate system is (O C1 − X C1 Y C1 Z C1 ). e origin O C1 is defined as the optical center of the camera, and Z C1 axis coincides with the optical axis.u 1 and v 1 axes of the computer pixel coordinate system (o 1 − u 1 v 1 ) are parallel to the X C1 and Y C1 axes.
e abscissa and ordinate values represent the columns and rows of the image, respectively.e image plane coordinate system is ). e ideal principal point O 1 lies at the center of the image plane.e ideal image point of the point P in the left image plane is P 1 (x p1 , y p1 ).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 P 1 ′ (x p1 ′ , y p1 ′ ), and the actual principal point is expressed as o 1 ′ (u 1 ′ , v 1 ′ ), respectively.e right camera coordinate system is (O C2 − X C2 Y C2 Z C2 ).If it rotates around X C2 , Y C2 , and Z C2 axes with angles ε, ψ, and κ in turn, the new formed coordinate system (O − XYZ) will be parallel to the left camera coordinate system (O C1 − X C1 Y C1 Z C1 ), where R is the rotation matrix as in equation (1).e origin O C1 can be translated to Transformation of the 3D world coordinate system to the left camera coordinate system can be described using 3 × 3 orthogonal matrix R and 3 × 1 vector T as follows: Equation ( 2) can be rewritten in the matrix form by using homogeneous styles: According to the mapping relationship between computer pixel coordinates and three-dimensional world coordinates, the world coordinates of space point P can be calculated as 2 Shock and Vibration where f 1 and f 2 are the focal lengths of the left and right cameras and (X 1 , Y 1 ) and (X 2 , Y 2 ) are the coordinates of the target points in the plane coordinates of the left and right images.e object movements are captured by using a binocular vision system.And the dynamic images are treated with segmentation and identi cation.
e authors developed a kind of tracking and matching technology for coded targets pasted on the blade [17].e cameras are calibrated with high accuracy by bundle alignment.An adaptive Wiener lter and threshold segmentation were explored to e ectively distinguish independent mark points.
ree-dimensional coordinates with subpixel precision are obtained by homogeneous correspondence.
e 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.

Automatic Identification of Modal Parameters
3.1.Data-Driven Stochastic Subspace Identi cation.After the measurement of binocular vision, the vibration data are used to identify the operational modal parameters.e 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 identi cation methods and indicted that parametric approach outperforms the nonparametric approaches.Here, we use the parametric approach, data-driven stochastic subspace identi cation (SSI-Data), to extract the modal parameters from the timedomain data.
e discrete SSI model of the structure can be described as follows [12]: where t represents the index of sampled data and x t ∈ R 2n×1 and y t ∈ R k×1 are the state vector and output vector corresponding to t, respectively.w t and v t are white noises with mean 0, A is the state space matrix, and C is the output matrix.e vibration response data of multiple locations are gathered in a block Hankel matrix H. Hankel matrix H can

Shock and Vibration
be divided into a past data matrix Y p ∈ R ka×b and future part matrix Y f ∈ R ka×b , which are defined as where y i � [y i1 , y i2 , . . ., y ik ] T ∈ R k×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 of Y p and Y f , we can have where R ∈ R 2ka×b is a lower triangular matrix; Q ∈ R b×b is an orthogonal matrix; R 11 � R(1 : ka, 1 : ka) includes the first ka rows and ka columns of R; and R 21 � R((ka + 1 : ka + k), 1 : ka) includes from ka + 1 to ka + k rows of R and 1 to ka columns of R; R 22 , R 31 , R 32 , and R 33 can be determined similarly.Q 1 � Q(1 : ka, :) includes the first ka rows of Q; Q 2 and Q 3 can be determined similarly.
From equation ( 8), we have e projection of row space of Y p into row space of Y f is defined as χ a ∈ R ka×b , which is By applying singular value decomposition (SVD) to the projection matrix, where U 1 ∈ R ka×ka , V 1 ∈ R b×b are the unitary matrices composed of left and right singular vectors corresponding to nonzero singular values, and S 1 is a diagonal matrix with singular values arranged in the descending order.
e projection matrix can be factorized as the product of the observability matrix and the Kalman filter sequence: Combining equations ( 11) and (12) gives e matrixes Y − f and Y + p can be obtained by moving the first k rows of Y f matrix to Y p .According to equations ( 8) and ( 9), If the row space of Y + p is projected to Y − f , the projection matrix χ a−1 can be obtained: And we have Defining matrix Y a|a as the rows from ka + 1 to k(a + 1) of H matrix, it is expressed in the form QR-factorization: System matrixes A and C can be obtained by solving the set of equations in least squares sense.

A
System state matrix A is factorized to produce the eigen values: where Λ � diag(λ s ), s � 1, . . ., n.If λ c s is a continuous-time eigen value, then there is where ξ s is the damping ratio and ω s is the circular frequency.ey can be computed from the following equation: where W f and W Φ are the weights of natural frequencies and mode shapes, df and dΦ are the tolerances of natural frequencies and mode shapes, f ki and σ(f ki ) 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 r k ≤ 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.e average value of physical modal parameters for each group is set as the final identified modal parameters.

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 (Figure 2). Figure 3 shows the blade image of labeled targets on the surface taken by using the cameras.ere are 10 pairs of homogeneous points.e 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.e sampling frequency of the cameras is 120 Hz, and 1000 continuous images are recorded.
e coordinates of the targets on the image plane can be extracted from the image.
e three-dimensional coordinates of the object can be obtained by solving the correspondence problem with the point tracking (PT).en, the spatial coordinate changes of the same target point can be calculated through 1000 continuous images.e vibration response of the blade can be obtained by normalizing displacement changes.Figures 4-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.

Analysis of Modal Identification Results.
e SSI-Data is performed to identify the modal parameters of the blade structure from photogrammetric response data.A Hankel matrix H 1 with 20 rows and 900 columns is constructed.
e maximum model order is set to 100.Stabilization diagrams are drawn with frequency (abscissa coordinate) and the model order (ordinate coordinate).e results in X, Y, and Z directions under random excitation are shown in Figures 9(a ere 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.e Shock and Vibration 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 H 2 with row number 30 and column number 850 and drew stabilization diagrams in the Z direction under random excitation (Figure 10).e 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.
After eliminating spurious poles, we can further extract more accurate modal parameters by spectral clustering.e frequency is taken as abscissa coordinate, and the damping ratio is the ordinate coordinate; we have the clustering results as shown in Figure 12. e final modal parameters are set as the mean values of these clustered samples.e first       Shock and Vibration 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.
e 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.e vibration of the blade tip is detected, and the modal parameters are identified with the presented algorithm.
e 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.

4
Shock and Vibration where Redenotes the real part and the corresponding modal shape Φ is Φ � Cψ. (22) 3.2.Spurious Pole Discrimination Method.At the process of the modal parameters identification, the maximal model order N 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.e 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 H i .A set of poles P ij will be obtained from every H i by using SSI-Data process.en the N poles P ij (j � 1, . . ., N) are clustered to produce N clusters, and the center of each cluster gives a possible pole.e possible poles which have nearest natural frequency will be regrouped and each group G k (k � 1, . . ., N) has m poles g ki (i � 1, . . ., m).By introducing tolerance and modal confidence, the poles are judged and discriminated.e similarity between these m poles is defined as )-9(c), while diagrams under impact excitation in Figures9(d)-9(f ).

Figure 3 :Figure 4 :
Figure 3: Blade images taken by using the cameras: (a) left camera and (b) right camera.

Figure 5 :
Figure 5: Vibration displacement curves of point 10 under random excitation: (a) X direction, (b) Y direction, and (c) Z direction.

Figure 6 :
Figure 6: Vibration displacement curves of point 1 under impact excitation: (a) X direction, (b) Y direction, and (c) Z direction.

Figure 7 :
Figure 7: Vibration displacement curves of point 10 under impact excitation: (a) X direction, (b) Y direction, and (c) Z direction.
Table 1 lists the lowest three frequencies of the blade identified with our method and hammer test.ey are quite consistent in Table 1, and the relative errors of the three frequencies are less than 3.2%.e results show that the vibration response data measured by the