Research Article  Open Access
Gang Yu, "An Underdetermined Blind Source Separation Method with Application to Modal Identification", Shock and Vibration, vol. 2019, Article ID 1637163, 15 pages, 2019. https://doi.org/10.1155/2019/1637163
An Underdetermined Blind Source Separation Method with Application to Modal Identification
Abstract
In structural dynamic analysis, the blind source separation (BSS) technique has been accepted as one of the most effective ways for modal identification, in which how to extract the modal parameters using very limited sensors is a highly challenging task in this field. In this paper, we first review the drawbacks of the conventional BSS methods and then propose a novel underdetermined BSS method for addressing the modal identification with limited sensors. The proposed method is established on the clustering features of timefrequency (TF) transform of modal response signals. This study finds that the TF energy belonging to different monotone modals can cluster into distinct straight lines. Meanwhile, we provide the detailed theorem to explain the clustering features. Moreover, the TF coefficients of each modal are employed to reconstruct all monotone signals, which can benefit to individually identify the modal parameters. In experimental validations, two experimental validations demonstrate the effectiveness of the proposed method.
1. Introduction
One of the main issues in structural dynamic analysis is to identify the modal parameters, e.g., frequency, damping, and modal shape. Starting from [1, 2], the blind source separation (BSS) technique has become more and more popular in modal identification, due to its straightforward, computationally fairly efficient, nonparametric, and requiring no prior information of the dynamic system. From early literatures, it can be known that the applications of conventional BSS methods, e.g., independent component analysis, secondorder blind identification, and their improved versions, mainly focus on fundamental research and theorem analysis [3–7]. As research interest increased, some limitations of conventional BSS techniques are gradually recognized. An essential assumption that guarantees successful application of conventional BSS methods is that the number of sensors should be not less than that of sources, namely, the issue of the determined BSS. However, considering various limitations in practice, e.g., costs, and difficulties in access, the installation of multiple sensors may not be feasible. Furthermore, without the prior knowledge of active modals in most cases, the requirement of adequate sensors is too strict to meet [8, 9]. Therefore, there is a strong requirement to develop the underdetermined BSS techniques for structural dynamic analysis, which can address more sources with limited sensors.
Recently, an underdetermined BSS method called sparse component analysis (SCA) draws many attentions in structural dynamics [10, 11]. By transforming the original signal from time domain into sparse domain, the mixing matrix and monotone modal sources can be precisely obtained. Then, we can estimate the modal parameters by the monotone modal identification method. Sadhu utilized wavelet transform (WT) to transform raw signals into wavelet domain to achieve the sparsity, and then principal component analysis and parallel factor decomposition method were employed for determining the mode shape vectors, natural frequency, and damping, respectively [12, 13]. Yang proposed a scheme using the clustering algorithm to estimate the mode shape matrix in frequency domain and then separated monotone modal responses via linear programming techniques [14]. Yu estimated the mode shape combining with single source point method and timefrequency ratio of mixtures (TIFROM) after executing shorttime Fourier transform (STFT), and then natural frequency and damping can be identified by detecting TF ridges [15]. Qin proposed a hybrid method combining [14, 15] to improve the accuracy of matrix estimation and source separation [16]. With more papers being carried out, the SCA method shows more and more powerful ability in modal identification than conventional determined BSS techniques [17, 18].
However, the SCA method only works well in the situation of instantaneous mixture, which means that the modal response signals should reach the sensors at the same time [10, 11]. In practical applications, this requirement is also too strict to meet. Therefore, there are two main purposes in this paper: first to provide detailed analysis on the theorem of SCA and point out its limitations in practical applications and then to establish a novel BSS method that can overcome the limitations. With these prime objectives, this paper is organized in the following manner. Section 2 illustrates the detailed theorem of SCA. Section 3 introduces a novel BSS method. Two experimental validations are carried out in Section 4. Then, the conclusion is followed in Section 5.
2. Principle of SCA
2.1. Motivation of SCA
The SCA framework illustrates that, in the sparse domain (e.g., TF domain or frequency domain), the sparse coefficients of different monotone modals can cluster into the straight line. This phenomenon motivates that it is possible to separate the sparse coefficients based on the clustering features. Therefore, the clustering features in the SCA are the most important issue that should be well explored. Given a numerical example, where three harmonic sources are linearly mixed into two observations :
We first display the timeseries values of two observations by the scatter plot, which is shown in Figure 1(a). It can be seen that any useful information on three sources cannot be obtained in time domain. If we use STFT to transform the original signal into TF domain, then we plot the scatter of the real value of two observations, as shown in Figure 1(b). It can be obviously observed that there appear three clustering lines. Meanwhile, the directions of these straight lines correspond to the column vector of the mixing matrix . By estimating the clustering directions, we can obtain the mixing matrix, and then these sources can be separated by linear programming techniques or L1norm minimum algorithm. Therefore, it motivates us to first explore the reason that why the clustering features appear.
(a)
(b)
2.2. Lissajous Figure
Before analysis of the SCA, it is necessary to illustrate some essential background knowledge. In mathematics, a Lissajous figure is the graph of a system of parametric equations:which can describe the complex harmonic motion. This motion system was investigated in detail by Lissajous in 1857. The appearance of the figure is highly dependent on the ratio . For and , the figure is a line, whose direction is determined by , as shown in Figure 2(a). For and , the figure is an ellipse, whose shape is determined by and , as shown in Figure 2(b).
(a)
(b)
2.3. Sparse Transform Method
Selecting an appropriate sparse method is helpful for executing the SCA algorithm. The timefrequency analysis method has been accepted as the most effective method to achieve sparse in most of papers [10–15, 19–22]. Therefore, we select the STFT as the sparse method in this paper. Given the regular STFT aswhere is the moved window and is the truncated signal. Consider the modified STFT with additional phase shift than regular STFT as
According to Parseval’s theorem, the modified STFT can be written aswhere is the Fourier transform (FT) of the window function, , and is the FT of the signal. Consider a purely harmonic signal whose frequency is :
Considering that , the STFT of the harmonic signal can be written as
Equation (8) denotes us a first impression on how STFT works [23–27]. The STFT of the harmonic signal is constituted by a series of harmonic signals with the same frequency (which is consistent with the original signal) but different amplitudes (which is determined by both of the signal amplitude and the Fourier transform of window function). Given two harmonic signals with different frequencies:
Then, the corresponding STFTs can be written as
Based on the instantaneous mixture, two sources (9) are first mixed into two observations linearly:
According to the linearity property of STFT, the STFTs of two observations can be written as
Herein, it is needed to assume , which means that the frequency distance between two sources is larger than the frequency support of window function. Therefore, this assumption suggests that it is better to select the long window function, which makes sure that the frequency support is narrow enough. Then, we can have the following equations:
Equations (13) and (14) mean that there exists no overlap between two sources in the TF domain. Therefore, for , we can construct the following formulation:where denotes taking the real part. Then, for , we can construct the following formulation:
Equations (15) and (16) denote that, for the instantaneous mixture, when we plot the scatter of the real part of STFT of two observations, i.e., versus in each frequency bin, it is equal to plot a series of Lissajous figure of the special case (, ), i.e., straight lines. These lines have the same direction which is determined by the column vector of mixing matrix, e.g., or . This is the reason that why there can appear several clustering lines in the scatter plot (see Figure 1(b)), which can also be employed to estimate the mixing matrix.
2.4. Delay Mixture
For the signal of equation (9), we consider the delay mixture:where denotes time delay. According to equation (8), the STFT of two observations can be written as
Similar to equation (13), for , we have
Then, we can have the following equation:
And for ,
Then, we can have the following equation:
Equations (20) and (22) denote that, for the delay mixture, when we plot the scatter of the real part of STFT of two observations, i.e., versus in each frequency bin, it is equal to plot a series of Lissajous figure with the ellipse case (, ). The shape of these ellipses is determined by both amplitude and time delay in observations. We then utilize the numerical example to further illustrate above analysis. For the numerical example of equation (2), we consider the delay mixture, the scatter plot in the time domain and TF domain is displayed in Figures 3(a) and 3(b), respectively. It can be seen that, both in the time domain and the TF domain, the linear clustering features completely disappear, which lead to that we cannot obtain any useful information on the mixing matrix.
(a)
(b)
From above analysis, it can be known that the sparsity of the observation signals is highly sensitive to time delay. Considering the dynamic engineering in practice, we need to record the response signals by sensors located at different positions. Even if we can measure the signal at the same time, we cannot avoid the propagation delay of sources between all sensors. When the time delay of recorded signals heavily damages the sparsity of the signals, the SCA will fail to separate sources. This drawback of the SCA greatly restricts itself in practical application.
3. A Novel BSS Method
3.1. Theory of the Proposed Method
In this section, we develop an effective BSS method that can overcome the limitations in the SCA method. For a specified harmonic signal whose frequency is , we have the reconstruction formulation of STFT aswhere . Equation (23) denotes that if we can obtain the frequency band of each source in the TF representation, i.e., , we can reconstruct each source effectively. Therefore, it motivates us to first estimate the frequency band for each source. Considering from the viewpoint of signal energy, at the same time, the energy ratio of the same source recorded by two sensors should be determined, which depends on the source amplitude in recorded signals. Therefore, it inspires us to estimate the frequency band by calculating the TF energy of the signal. Based on the STFT of signals, we first define the formulation of frequency energy function as
For the signal of equation (18), we can have the following equation:
For the frequency band of , i.e., , according to equation (19), equation (25) can be written aswhere . Similarly, for the frequency band of , i.e., , equation (25) can be written as
Equations (26) and (27) denote that, for the frequency band of each source, the ratio of frequency energy is consistent, which is determined by the source amplitude in recorded signals. When we plot the scatter of the frequency energy data of two observations, i.e. versus , if , there will appear several clustering lines which can correspond to the frequency band of each source. Therefore, we can estimate the frequency band of each source according to these clustering features. For the delay mixture of numerical signals of equation (2), the scatter plot of frequency energy is shown as Figure 4(a). It can be seen that there are three obvious clustering lines.
(a)
(b)
In Figure 4(a), each point in the clustering lines corresponds to the frequency bin of source in TF representation. Then, it is required to develop an effective procedure to estimate these points. As observing the scatter plot, each clustering line has the maximum value located at the end of line (as red label “”), and other points have the same cosine distance with this end point. The point of maximum value can be detected by the peak data by summing all frequency energy data:
Therefore, we can first obtain the maximum value point, as shown in Figure 4(b) (as red label “”), and then calculate the cosine distance with it to estimate other points bywhere is the detected peak frequency and is a low value which is suggested to be empirically set as 0.004. If condition (29) is satisfied, these scatter points should belong to the same source. By repeating this procedure, we can obtain the frequency bands of all sources. After frequency band being estimated, we can reconstruct all sources by inversion of STFT. For the numerical signal, the three recovered sources are shown in Figure 5. It can be observed that three sources are well separated as monocomponents.
A known drawback of STFT reconstruction is the boundary effect, i.e., the amplitude in the beginning and end of the sources cannot be well recovered. In order to eliminate the boundary effect of the recovered sources, we introduce a padding line, aswhere is the moved window and denotes the discrete beginning and end time point of the analysed signal. The padding line of the numerical signal is shown in Figure 6. And then, each source is corrected through multiplying this padding line. The three sources that consider padding method are shown in Figure 7. It can be seen that the amplitudes of all sources are also recovered precisely.
Therefore, the proposed BSS method can be summarized as follows:(1)Utilize STFT to transform recorded signals into the TF domain(2)Obtain the frequency energy data according to equation (24)(3)Sum all frequency energy data and pick the peak data by the peak detection method(4)Calculate the cosine distance of other points with detected peak data to estimate the frequency band(5)Recover each source according to the estimated frequency band by inversion of STFT(6)Eliminate the boundary effect of each source by the padding method
3.2. Numerical Signal Analysis
In this subsection, we consider the performance of the proposed method in analysing a fivecomponent signal added with heavy Gaussian noise. This fivecomponent signal mixed into two observations is modelled aswhere the sampling frequency is 512 Hz and sampling time is 2 s. The waveforms of the two observations are plotted in Figure 8. In Figure 9, the timefrequency representations of two observations generated by STFT are displayed. Then, the scatter of the frequency energy data is calculated and plotted in Figure 10. It can be seen that there are five obvious clustering lines in Figure 10, which correspond to the frequency band of five sources. The eventual separated sources are plotted in Figure 11. It is shown that all sources are well separated into monocomponent signals.
(a)
(b)
(a)
(b)
4. Experimental Validation
In this section, we employ two experiments to validate the proposed method. These two experiments have been analysed in two published papers [15, 17]. Therefore, their analysed results can be the comparison references. To test the ability of the proposed method comparable to the SCA, we also display the results analysed by SCA.
4.1. Experiment 1
As shown in Figure 12, the structure is a uniform TC4 titaniumalloy column of , and three acceleration sensors are mounted on the beam which is excited by impact hammer. The sampling frequency is 2560 Hz, and we take 0.5 s data to analyse. The acceleration responses X1, X2, and X3 are shown in Figure 13. For comparative analysis, we list the identified results in [15], such as frequency and damping ratio, as shown in Table 1.

To show the processing procedure of the SCA method, we first utilize STFT to transform three recorded signals into the TF domain. Then, the scatter of original timeseries signals and the real part of the STFT results are plotted in Figures 14(a) and 14(b), respectively. It can be seen that the original signal in time domain cannot provide any useful information about the mixing matrix. However, in Figure 14(b), there obviously appear four clustering lines for the recorded signals, which should correspond to four column vectors of the mixing matrix. Then, we can separate each source using this mixing matrix, as shown in Figure 15. It can be seen that, through the SCA method, each source is well separated as monotone exponentially decaying sinusoid. Then, we utilize the monotone mode method to identify modal parameters, and the identification results are listed in Table 1. It is known that the SCA provides a similar modal identification result with [15].
(a)
(b)
Then, for the proposed method, according to equation (24), we first calculate the frequency energy function. The detected peaks of the sum of are shown in Figure 16(a), where “” denotes the detected peaks. It can be seen that six peaks are detected. Among them, two more sources are discovered than SCA, which is caused by unexpected interference. In Figure 16(b), we plot the scatter of frequency energy data, which also shows six obvious clustering lines. To facilitate comparison with SCA, we do not list the separated results of two interference sources. Therefore, four separated sources are shown in Figure 17. It is shown that each source is well recovered as the monomodal response. The estimated frequencies and damping are listed in Table 1, which show satisfactory accuracy [15] as well. It can be concluded that when the sparsity of signal in the TF domain is satisfactory, the monomodal sources can be well separated both by SCA and the proposed method.
(a)
(b)
4.2. Experiment 2
As shown in Figure 18, the structure is a uniform steel cantilever beam of , and five displacement sensors are settled upon the beam which is excited by impact hammer. The sampling frequency is 1600 Hz, and we take 1000 samples to analyse. The vibration responses X1–X5 are shown in Figure 19. Meanwhile, we first list the identified modal parameters by the method in reference [17] in Table 2.

To illustrate the ability of the proposed method applied in the underdetermined case, we only deal with three sensor data, i.e., X1, X2, and X5. The scatter plot of STFT results of the recorded signals is first displayed in Figure 20. It can be seen that, due to influence of time delay, there appear several ellipses, which make the SCA hard to estimate the accurate mixing matrix. The separated results by the SCA are also displayed in Figure 21. It can be observed that the recovered sources S3, S4, and S5 are unwillingly mixed with other sources.
For our proposed method, according to equation (24), we first calculate the frequency energy . And the scatter plot of frequency energy data is shown in Figure 22, where “” denotes the detected peaks. It can be seen that there appear five obvious clustering lines, which correspond to the frequency band of five sources. The sources separated by the proposed method are listed in Figure 23, which show all monomodal responses being recovered successfully. The identified parameters are listed in Table 2. Due to the recovered sources being close to harmonic signals, it cannot provide accurate estimation for damping, which is neglected in the identified parameters. From the identified parameters, it can be seen that the identified modal parameters show satisfactory accuracy with [17].
From above analysis, it can be known that, when the sparsity of the signal in the TF domain is heavily damaged by time delay, the SCA may fail to achieve accurate estimation of mixing matrix and recovery of sources. However, our proposed method can deal with it successfully. Furthermore, we consider the computational cost of the proposed method in dealing with the experimental data. Because it is helpful for judging whether the proposed method can be used in realtime applications or not. The tested computer configuration is as follows: Intel Core i76500 2.5 GHz, 8.0 GB of DDR3 RAM, Windows 10 OS, and MATLAB version R2016a. The required computation times in the first experiment and the second experiment are 0.91 s and 0.82 s, respectively. It can be concluded that the proposed methods can finish the processing within one second, which is in the acceptable range in realtime applications.
5. Conclusions
In this paper, we focus on a hot topic on application of underdetermined BSS in structural dynamic analysis. The SCA technique, as an effective underdetermined BSS method, has drawn many attentions in recent papers. However, all of these papers only pay attention to the application of SCA, but none of them illustrates its theorem in detail, which lead to some limitations of SCA being neglected. Therefore, we start with the classical TF method and harmonic signal, to illustrate the detailed theorem of SCA. Meanwhile, we point out its limitation in practical engineering. Furthermore, we propose a novel BSS method. By means of TF transform and estimation of frequency band, the proposed method can effectively deal with the problem existed in SCA.
The time delay between recorded signals is mainly determined by the propagation speed and distance. For the structures in experiments 1 and 2, they have the similar width and height, but distinct lengths. The structure in experiment 2 is obviously longer than the structure in experiment 1. Under the similar condition of propagation speed, the structure length has become the critical influence of the time delay. It is obvious that a longer structure will produce a larger time delay. Therefore, we suggest that the SCA technique should be applied to analyse small structures, which makes sure the time delay being in the acceptable range.
Data Availability
The data used to support the findings of this study are available from the corresponding author upon request.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
The author acknowledges Kaiping Yu and Kai Yang for their experiment data greatly. This work was supported by the Natural Science Foundation of China (grant no. 61901190) and University of Jinan Research Funds for Doctors (grant no. XBS160100391).
References
 R. Brincker, “Some elements of operational modal analysis,” Shock and Vibration, vol. 2014, Article ID 325839, 11 pages, 2014. View at: Publisher Site  Google Scholar
 L. Yu and H. Song, “Scaling mode shapes in outputonly structure by a masschangebased method,” Shock and Vibration, vol. 2017, 2017. View at: Publisher Site  Google Scholar
 F. Poncelet, G. Kerschen, J.C. Golinval, and D. Verhelst, “Outputonly modal analysis using blind source separation techniques,” Mechanical Systems and Signal Processing, vol. 21, no. 6, pp. 2335–2358, 2007. View at: Publisher Site  Google Scholar
 W. Zhou and D. Chelidze, “Blind source separation based vibration mode identification,” Mechanical Systems and Signal Processing, vol. 21, no. 8, pp. 3072–3087, 2007. View at: Publisher Site  Google Scholar
 S. Chauhan, R. Martell, R. J. Allemang, and D. L. Brown, “Application of independent component analysis and blind source separation techniques to operational modal analysis,” in Proceedings of the 25th IMAC, Orlando, FL, USA, January 2007. View at: Google Scholar
 S. I. McNiell and D. C. Zimmerman, “A framework for blind modal identification using joint approximate diagonalization,” Mechanical Systems and Signal Processing, vol. 22, no. 7, pp. 1526–1548, 2008. View at: Publisher Site  Google Scholar
 B. Hazra, A. J. Roffel, S. Narasimhan, and M. D. Pandey, “Modified crosscorrelation method for the blind identification of structures,” Journal of Engineering Mechanics, vol. 136, no. 7, pp. 889–897, 2010. View at: Publisher Site  Google Scholar
 J. Antoni and S. Braun, “Special issue: blind source separation,” Mechanical Systems and Signal Processing, vol. 19, no. 6, pp. 1163–1165, 2005. View at: Publisher Site  Google Scholar
 J. Antoni, “Blind separation of vibration components: principles and demonstrations,” Mechanical Systems and Signal Processing, vol. 19, no. 6, pp. 1166–1180, 2005. View at: Publisher Site  Google Scholar
 M. Zibulevsky and B. A. Pearlmutter, “Blind source separation by sparse decomposition in a signal dictionary,” Neural Computation, vol. 13, no. 4, pp. 863–882, 2001. View at: Publisher Site  Google Scholar
 P. Bofill and M. Zibulevsky, “Underdetermined blind source separation using sparse representations,” Signal Processing, vol. 81, no. 11, pp. 2353–2362, 2001. View at: Publisher Site  Google Scholar
 A. Sadhu, B. Hazra, S. Narasimhan et al., “Decentralized modal identification using sparse blind source separation,” Smart Material and Structures, vol. 20, no. 12, Article ID 125009, 2011. View at: Publisher Site  Google Scholar
 A. Sadhu, B. Hazra, and S. Narasimhan, “Decentralized modal identification of structures using parallel factor decomposition and sparse blind source separation,” Mechanical Systems and Signal Processing, vol. 41, no. 12, pp. 396–419, 2013. View at: Publisher Site  Google Scholar
 Y. Yang and S. Nagarajaiah, “Outputonly modal identification with limited sensors using sparse component analysis,” Journal of Sound and Vibration, vol. 332, no. 19, pp. 4741–4765, 2013. View at: Publisher Site  Google Scholar
 K. Yu, K. Yang, and Y. Bai, “Estimation of modal parameters using the sparse component analysis based underdetermined blind source separation,” Mechanical Systems and Signal Processing, vol. 45, no. 2, pp. 302–316, 2014. View at: Publisher Site  Google Scholar
 S. Qin, J. Guo, and C. Zhu, “Sparse component analysis using timefrequency representations for operational modal analysis,” Sensors, vol. 15, no. 3, pp. 6497–6519, 2015. View at: Publisher Site  Google Scholar
 K. Yang, K. Yu, and Q. Li, “Modal parameter extraction based on Hilbert transform and complex independent component analysis with reference,” Mechanical Systems and Signal Processing, vol. 40, no. 1, pp. 257–268, 2013. View at: Publisher Site  Google Scholar
 A. Sadhu, S. Narasimhan, and J. Antoni, “A review of outputonly structural mode identification literature employing blind source separation methods,” Mechanical Systems and Signal Processing, vol. 94, pp. 415–431, 2017. View at: Publisher Site  Google Scholar
 J. Zheng, H. Pan, T. Liu, and Q. Liu, “Extremepoint weighted mode decomposition,” Signal Processing, vol. 142, pp. 366–374, 2018. View at: Publisher Site  Google Scholar
 S. Wang, X. Chen, I. W. Selesnick, Y. Guo, C. Tong, and X. Zhang, “Matching synchrosqueezing transform: a useful tool for characterizing signals with fast varying instantaneous frequency and application to machine fault diagnosis,” Mechanical Systems and Signal Processing, vol. 100, pp. 242–288, 2018. View at: Publisher Site  Google Scholar
 S. Chen, X. Dong, G. Xing, Z. Peng, W. Zhang, and G. Meng, “Separation of overlapped nonstationary signals by ridge path regrouping and intrinsic chirp component decomposition,” IEEE Sensors Journal, vol. 17, no. 18, pp. 5994–6005, 2017. View at: Publisher Site  Google Scholar
 H. Ma, R. Song, X. Pang, and B. Wen, “Timevarying mesh stiffness calculation of cracked spur gears,” Engineering Failure Analysis, vol. 44, pp. 179–194, 2014. View at: Publisher Site  Google Scholar
 G. Yu, M. Yu, and C. Xu, “Synchroextracting transform,” IEEE Transactions on Industrial Electronics, vol. 64, no. 10, pp. 8042–8054, 2017. View at: Publisher Site  Google Scholar
 G. Yu, Z. Wang, and P. Zhao, “Multisynchrosqueezing transform,” IEEE Transactions on Industrial Electronics, vol. 66, no. 7, pp. 5441–5455, 2018. View at: Google Scholar
 G. Yu, Z. Wang, P. Zhao, and Z. Li, “Local maximum synchrosqueezing transform: an energyconcentrated timefrequency analysis tool,” Mechanical Systems and Signal Processing, vol. 117, pp. 537–552, 2019. View at: Publisher Site  Google Scholar
 K. Yu, H. Ma, J. Zeng et al., “An amplitude weak component detection technique based on normalized timefrequency coefficients and multisynchrosqueezing operation,” Measurement Science and Technology, vol. 30, no. 5, Article ID 055008, 2019. View at: Publisher Site  Google Scholar
 K. Yu, H. Ma, H. Han et al., “Second order multisynchrosqueezing transform for rubimpact detection of rotor systems,” Mechanism and Machine Theory, vol. 140, pp. 321–349, 2019. View at: Publisher Site  Google Scholar
Copyright
Copyright © 2019 Gang Yu. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.