About Us(CRI)
Estimations are free. For more information,
please send a mail
-->here<--
Power Spectral Density computation (Spectral Analysis)
MicroJob Package Deal Service PD001A/B(Copyright Cygnus Research International (Apr 20, 2015))

Page 7 of PD001A/B User Guide

Table of Content
(1) Summary of PD001A/B
(1-1) What this Package Deal service, PD001A does
(1-2) What PD001B does
(1-3) Summary of selectable computational parameters
(1-3-1) Detrending data (Applicable to all cases)
(1-3-2) Confidence interval percentage (Applicable to all cases)
(1-3-3) Taper ratio of a Tukey window (Applicable to cases B1, B2 and B3 only)
(1-3-4) FDS bin-width (Applicable to all cases except for cases A1, B1 and C1)
(1-4) Computational procedure
(1-5) Your Data (Input data)
(1-6) Products
(1-6-1) Data file (Product file)
(1-6-2) Graphs (Figures) in PDF (PD001A only)
(1-7) Price and ordering procedure
(1-8) Confidentiality
(2) Detail of the product
(2-1) Product files
(2-1-1) Format of Product file
(2-1-2) Explanation of contents of product files
(2-1-2-1) Frequency
(2-1-2-2) Period
(2-1-2-3) PSD
(2-1-2-4) Confidence interval of PSD
(2-1-2-5) Percent of variance
(2-1-2-6) Amplitude (Amplitude spectrum)
(2-1-2-7) Phase
(2-1-3) How to use/interpret data in our product file
(2-1-3-1) Drawing graphs
(2-1-3-2) PSD and Confidence interval of PSD
(2-1-3-3) Estimating amplitude of variations at specific frequencies using amplitude spectrum
(Case 1) When the signal is sinusoid and the amplitude of the signal is constant or almost constant in time
(Case 2) When the signal itself is not a sinusoid but the composition of sinusoids of various frequencies
(Case 3) When the signal is sinusoid but the amplitude of the signal varies in time
(Case 3-1) Variation of amplitude itself is sinusoidal
(Example 1) Amplitude of sinusoidal oscillation of frequency f0 varies sinusoidally.  Frequency of amplitude variation is f1 and f1 is smaller than f0. Both f0 and f1 are constants.
(Example 2) Amplitude of sinusoidal oscillation of frequency f0 is summation of a constant and a sinusoidal variation of frequency f1.
(Case3-2) Amplitude of sinusoidal oscillation of frequency f0 varies but that variation does not resemble to a sinusoid.
(Case4) When the signal is periodic but the signal does not look like a sinusoid
(2-1-3-4) Obtain some idea how strong the variations within a specific frequency range is.
(2-1-3-5) Compute digitally filtered time series data
(2-1-3-6) Alternative method of amplitude estimation
(2-2) Graph
(3) How to prepare data for this package
(4) Guide for setting proper selectable parameters
(4-1) Detrending and removing average; Yes or No? Default setting is Yes.
(4-1-1) What is detrending? Why do we have to do that?
(4-1-2) Example of detrending
(4-1-3) When detrending becomes harmful
(4-2) Taper ratio of Tukey window; 0~100(%) Default setting is 10%
(4-2-1) Why we need to apply a window function
(4-2-2) What does window function do
(4-2-3) Why Hanning window is not necessarily a good solution
(4-2-4) Why we provide non-windowed and Tukey windowed result for PD001A
(4-2-5) The method to reduce dumping near the both ends of data while still using a Hanning window
(4-2-6) Implicit FDS
(4-2-7) Implicit FD'Smoothing' Really?
(4-2-8) Reduction of PSD and amplitude due to the window
(4-2-9) Does correction of reduction of PSD work well?
(4-2-10) Tukey window; how taper ratio affects its characteristics
(4-3) Bin-width of FDS; 1,3,5,7,... (Odd integer); Default setting is None, 3 and 7
(4-3-1) FDS works well for PSD
(4-3-2) FDS does NOT work well for amplitude spectrum
(4-3-3) Frequency resolution when we applied FDS
(4-3-4) Alternative method of FDS
(4-4) Percentage of confidence interval; any positive value larger than 0.0 and smaller than 100.0; Default is 95.0

Appendix 1; About numerical error generated by computers.
Appendix 2; Spectral leakage of power spectrum density function and window functions.
Appendix 3; Why we see a single peak instead of twin-peaks when data length is short.
Appendix 4; How the frequency of S1 affects the detectability of SR

(4-3) Bin-width of FDS; 1,3,5,7,... (Odd integer); Default setting is None, 3 and 7
The FDS we use is the simple un-weighted moving average (running mean) of PSD. For example, 5-bin FDS of PSD is

Moving average is used to smooth time series data quite often.

(4-3-1) FDS works well for PSD
Figure 4-30(a) shows how signal bin variance (PSD of signal bin multiplied by the frequency bandwidth of a bin) and data variance vary as we add more data when 3-bin FDS is applied. The data we used here contains single frequency sinusoid but nothing else and exactly the same as the one used for Figure 2-4. We divided data variance by 3 so that we can compare data variance directly with other values. If we compare this graph with Figure 4-24(b), we can see that FDS considerably reduced deviation of signal bin variance from data variance. Figure 4-30(b) shows the case when we made width of FDS larger. These graphs show that FDS works well to smooth PSD. Please note that the reduction of PSD is reasonable because frequency bandwidth of bins should be modified if we apply FDS as we will describe in (4-3-3). In any case, even though PSD of signal bin becomes smaller, equation (1) still holds because reduction of PSD of signal bin is compensated by increase of PSD of other nearby bins.

(4-3-2) FDS does NOT work well for amplitude spectrum
We showed what would happen if we apply FDS to amplitude spectrum in (2-1-3-3) (Figure 2-10, Figure 2-17). Let us explain why FDS does not work well for amplitude spectrum when window is not applied first. Figure 4-31(a) shows how amplitude of the signal bin varies when we applied 3-bin FDS as we add more data. Here, data contains only a single frequency sinusoid but nothing else as before. This graph is the same as Figure 2-10(a). To explain why the signal bin amplitude becomes larger than the correct value (black horizontal line; 1/3 of actual signal amplitude) when frequency of the signal bin does not match actual signal frequency, let us chose the case when data length is B in Figure 2-4. The data length B is indicated by short vertical solid black line in Figure 4-31(a). The expanded view of amplitude spectrum near the peak in this case is shown in Figure 4-31(b). The horizontal positions of plus (+) marks in this graph show the frequencies of bins and the horizontal position of vertical solid black line indicates signal frequency. The result of 3-bin FDS of the amplitude of Bin 2 is the simple 3-point average of amplitude of Bin 1, Bin 2 and Bin 3. Amplitudes of these bins are about 0.18, 0.90 and 0.31, respectively. Then, the average of these is (0.18+0.90+0.31)/3=0.46. Figure 4-31(c) shows reproduced time series of these bins (components of equation (3)). This graph indicates that variations of these components are not aligned but that is not really a surprise because frequencies and phases of these bins are all different. The problem of getting correct amplitude using FDS when frequency does not match, thus amplitude of Bin 1 and Bin 3 are not 0.0, is that FDS ignores the differences of frequencies and phases of bins. Now we add these three time series and let us call the result a synthesized time series. Figure 4-31(d) shows this synthesized time series (blue) and the data (purple). Synthesized time series contributes about 93% of data in term of variance. This graph shows that the amplitude of this synthesized time series is about the same as the amplitude of data, 1.0, but far smaller than the simple summation of amplitudes of three bins, 0.18+0.90+0.31=1.38.

On the other hand, when frequency of the signal bin matches actual signal frequency (data length A in Figure 2-4), the amplitudes of Bin 1 and Bin 3 become zero. Figure 4-32(a) shows time series of data in this case. Amplitude of Bin 2 is almost equal to the actual signal amplitude (1.0) and amplitudes of all other bins are zero as shown in Figure 4-32(b). Thus, simple 3-bin average is about 1/3 in this case.

Now, let us show why Hanning windowed amplitude spectrum becomes vastly larger than non-windowed amplitude spectrum when we applied FDS. For the sake of simplicity, let us chose data length to be data length A, so that frequency of the signal bin matches actual signal frequency and there are no spectral leakages. Figure 4-32(c) shows time series of data when we applied a Hanning window. If we applied a Hanning window, width of peak becomes wider as we described in (4-2-6) and now amplitudes of bin 1 and bin 3 are no longer zero as shown in Figure 4-32(d). The correction we described in (4-2-8) brings amplitude of bin 2 to the value almost equal to actual signal amplitude, 1.0, but it also amplifies amplitudes of Bin 1 and Bin 3 as shown in Figure 4-32(e). As a result, simple 3-bin average becomes about 2/3, which is about two times larger than the value when window is not applied. This is why Hanning windowed amplitude spectrum becomes vastly larger than non-windowed amplitude spectrum when we applied FDS.

Thus, in general, FDS does not produce meaningful result for amplitude spectrum and this is why we do not write amplitude in our product file when we applied FDS.

(4-3-3) Frequency resolution when we applied FDS
If we apply FDS, it is reasonable to think that the frequency resolution of resultant PSD should also be modified. Our n-bin FDS is a simple n-bin average and thus it is more natural to think that the frequency bandwidth of the resultant PSD is n times larger than the frequency bandwidth of original PSD. For example, in the case of 5-bin FDS shown in equation (25), it makes sense to pick up only Po(0) (uses Pi(0)~Pi(2)), Po(5) (uses Pi(3)~Pi(7)), Po(10) (uses Pi(8)~Pi(12)),..,Po(5xk; k=any positive integer),... but discard all others, and to make a frequency bandwidth 5 times larger. Alternative choice is picking up Po(2) (uses Pi(0)~Pi(4)), Po(7) (uses Pi(5)~Pi(9)), Po(12) (uses Pi(10)~pi(14)),..., Po(5xk-3; k=any positive integer),... but discarding all others. By doing so, PSD becomes smaller when we applied FDS but the frequency bandwidth becomes larger. This makes that the signal bin variance, PSD of signal bin multiplied by the frequency bandwidth, will not change by FDS, for example.

The primarily reason why we do not discard any of Po is that it would be convenient when we want to import our product file into applications for drawing graphs and/or doing further analyses. Another reason is that we can use intermediate Po such as Po(1) and Po(2) as interpolated Po between Po(1) and Po(6).

(4-3-4) Alternative method of FDS
Welch's method we described in (2-1-3-6) is also widely used to smooth PSD. At this moment, we do not have a plan to include that method in this package deal service for the reason that it might make ordering options too complicated. If you wish to use Welch's method instead of FDS we included in this package deal service, please contact us.

(4-4) Percentage of confidence interval; any positive value larger than 0.0 and smaller than 100.0; Default is 95.0
Our customer can choose any percentage value between 0.0 and 100.0 excluding 0.0 and 100.0. Commonly used percentage values are 80%, 90%, 95%, 99% and 99.5%, and we would recommend to choose one among these values unless you have a specific reason to choose otherwise.

Appendix 1; About numerical error generated by computers.
Numerical computations almost always produce some small computational errors in the process. This is because computers assign a fixed size of memory to each numbers they produce in the computational process and some numbers, such as the result of 1/3, are truncated so that they can be stored in a fixed sized memory. To demonstrate how this error is produced, let us assume we have a hypothetical computer, which always assigns a memory that can store 4 digits, this we call significant digit, for any number. Now let's compute (1/3-0.33)x100 from left to right with this computer. The result of 1/3 is 0.3333 or 3.333E-1. Now, subtract 0.33 from it. The value 0.33 is 0.3300 for our computer because significant digit is 4. The result of subtraction is 0.0033 or 3.3E-2. This will not be 3.333E-2 and now we have only two meaningful digits instead of four. Multiplying this by 100, we will get 0.3300 or 3.300E-1. Now let's do this computation manually. (1/3-0.33)x100=(1/3-33/100)x100=(100-99)/300x100=1/300x100=1/3. If we skip all the process in the middle but just compute the last of this, 1/3, with our hypothetical computer, result will be 0.3333 as before. Comparing this result with the previous result, 0.3300<->0.3333, we see that our computer generated computational error, 0.0033 during previous computation. In our actual computations we have 'about' 16 significant digits. The reason why we have to use a word 'about' is that computers do their business using binary numbers while we human use decimal numbers. What we would like to say here is that it is reasonable to think that the values of purple and red lines outside of the signal peak in Figure 2-5(a) are practically zero although actual numbers are not zero.

Appendix 2; Spectral leakage of power spectrum density function and window functions.
We described that spectral leakage occurs if we cannot achieve frequency match because signal bin alone cannot represent actual signal due to the difference of frequencies of signal bin and actual signal, and other bins are needed to hold equation (3) in (2-1-3-3). This is correct but, then, we might ask if using equation (3) is appropriate or not in that case. Our description in (2-1-3-3) merely indicates spectral leakage is not caused by computational problems.

In (4-2-1) we described that the 'cliff' at the boundary generated by the assumption of repetition of data (Figure 4-14(b)) is a cause of spectral leakage. The assumption we described there is consistent with the method we and many other people use to compute PSD and the explanation of the cause of spectral leakage is reasonable. However, it only gives us an idea how spectral leakage is generated but we cannot evaluate how large generated spectral leakage would be. Also we will have a problem explaining why we still have spectral leakage even after we applied a window function to remove 'cliffs'.

In this appendix we first introduce alternative way of looking at our limited length data. After that we show how spectral leakages are generated and how large they would be. Finally we describe how to reduce spectral leakages. The root cause of spectral leakages is the fact that the length of our data is limited.

(1) Alternative way of looking at our limited length data
Now let us introduce a new assumption. First, let us assume our data is a portion of infinitely long "Original data" x(n) shown in Figure A2-1(a).

We do not assume repetition of data described in (4-2-1).  Then we apply a window function w(n) shown in Figure A2-1(b) to this infinitely long data.  This window function is called a Boxcar window function due to its shape.  The result is "Our data 2" shown in Figure A2-1(c).  At this point all we need to know about "Original data" is the values between t=0 and T and we know them because our data covers that duration.

We can construct "Our data 2" in Figure A2-1(c) in a different way.  We start from our limited length data x'(n) shown in Figure A2-2(a).

This is the exact copy of the one we showed in Figure 4-14(a). Then we add infinitely long "Additional data" z(n) shown in Figure A2-2(b). The value of z(n) is zero for any n. This results that Fourier transform of z(n) is also zero for any frequency and, thus, adding z(n) to "Our data" x'(n) would not modify spectral characteristic of x'(n) at all. The result of adding z(n) to x'(n) is "Our data 2" y(n) shown in Figure A2-2(c). This "Our data 2" is exactly the same as the one shown in Figure A2-1(c). As we described before, we can add any number of data, value of which is exactly zero, to data we want to analyze to adjust data length so that we can achieve frequency match. We called it "zero padding". What we did here is the extreme form of zero padding. Since data length of "Our data 2" is infinite, frequency bandwidth of bins becomes zero (=1/infinite), and, thus, we cannot evaluate power spectrum "DENSITY" at discrete frequencies anymore.

The reason why we introduced these two methods to generate identical infinitely long data ("Our data 2") is that we will describe spectral leakage as a difference of spectral characteristic of infinitely long data shown in Figure A2-1(a) and "Our data 2" shown in Figure A2-1(c) (and Figure A2-2(c)). This approach makes sense only because spectral characteristic of "Our data 2" is the same as the spectral characteristic of "Our data" shown in Figure A2-2(a) because "Our data" is what we usually have, i.e. a limited length data. In other word, we consider limited length "Our data" as a result of the application of a Boxcar window function to an infinitely long data in this appendix. This leads us that if we study how changing the width of a Boxcar window function affect spectral of resultant windowed data (Figure A2-1), then we would know how changing the data length of "Our data" would affect spectral (Figure A2-2). Thus we will show how changing the width of a Boxcar window function affects spectral of resultant windowed data in this appendix.

Instead of using PSD we introduce similar quantity named PS for an infinitely long data as follows. First we extend the range of summation in equation (18) as

For a Boxcar windowed data, "Our data 2 (y)" in Figure A2-1(c), equation (A2-1) becomes

The transition from the one marked by blue underline to the one marked by red underline is possible because wm equals zero if m is less than 0 or greater than N-1 and equals 1 otherwise. Please note that rightmost part (red underline) of equation (A2-2) is the same as equation (18) except for the frequency. Since we can compute Yr and Yi at any frequency, we can compute them at frequencies where we computed PSD in section (2). In that case, the rightmost part of equation (A2-2) will be exactly the same as equation (18).

Using Yr and Yi we define a quantity PS for an infinitely long data as

This quantity PS is very much similar to PSD (equation (19)). One of the apparent differences is that there is no multiplication factor 2.0 in equation (A2-3). This is because we want to define PS both at positive and negative frequencies. Another apparent difference is that equation (A2-3) is not divided by T but this is simply because T is infinite for an infinitely long data. Aside from these apparent differences, the most important difference between PS and PSD is that PS is continuous in frequency and this becomes very convenient to understand why spectral leakage occurs.

Actual values of PS of "Our data 2 (y)" in figures A2-1(c) and A2-2(c) divided by T (equals width of a Boxcar window function) and multiplied by 2 at frequencies where we ordinarily compute PSD are exactly the same as PSD of "Our data (x')" in Figure A2-2(a). This is because equation (A2-2) will be exactly the same as equation (18) in that case as we described before.

(2) How spectral leakages are generated and how large they would be
In the case when we need to deal with real world data, we usually cannot compute equation (A2-1) because it is usually not possible to have an infinitely long data. However, in the case when we use hypothetical data we might be able to mathematically evaluate equation (A2-1). In the following we use hypothetical data we used in Case 1, (2-1-3-3), namely,

In this case Xr(f) and Xi(f) are

Equation (A2-5) shows that the Fourier transform of x(t) is not continuous but non-zero only at two specific frequencies. Using the result of this and the convolution theorem we described in (4-2-6), we can evaluate middle part (blue underline) of equation (A2-2) as follows.


Here W(f) is the Fourier transform (discrete) of a Boxcar window function,


k in the second line is an integer, greater equal 1 and less equal N/2.

The last part of equation (A2-6) (red underline) show that Fourier transform of Boxcar windowed x(t) is quite different from Fourier transform of x(t) shown by equation (A2-5). This difference is the spectral leakage. Decomposing W in equation (A2-6) into a real part (Wr) and an imaginary part (Wi) and substituting them into equation (A2-3), PS of the Boxcar windowed data is,

To show how the spectral leakage occurs by showing the convolution process graphically we have chosen two different cases. Figure A2-3(a) shows time series plot of data (Equation (A2-4)) we use.

Each data is shown as a small dot in this graph and straight blue lines connect dots. The first case is when window-width is 40 (data points). This case is equivalent to the data length A shown in Figure 2-4 in the sense that frequency match is achieved. From this graph we might get an impression that the second cycle is not completed in this case but the zero value point next to the end point of this case is actually the start point of the next cycle and does not belong to the second cycle. Another case is when window-width is 50 and this case is close to the data length C in Figure 2-4. We intentionally made these window-widths rather narrow so that following graphs are easy to see. In the following we set amplitude of x(t) (A in equation (A2-4)) 1.0. The sampling interval is 0.05, which results that fundamental frequency range is -10 to 10.

Before showing how the convolution process works we show PS for infinitely long but Boxcar windowed data similar to "Our data 2" in Figure A2-1(c) and PSD for limited length data similar to "Our data" in Figure A2-2(a) first. The red line in Figure A2-3(b) shows PS of Boxcar windowed infinitely long data x(t) for the case A. Since data length is infinite, we can evaluate continuous PS. We omitted PS in the negative frequency range because PS is symmetrical with respect to the zero frequency. The vertical blue arrow in this graph is PSD of limited length data multiplied by T and divided by 2 and we call it modified PSD. It is not an accidental that the tip of this arrow touches PS (red line). The modified PSD (blue arrow) and PS (red line) at the same frequency agree as we described before. Since this is the case when frequency match is achieved, PSDs of bins are zero except for a signal bin as we showed in Figure 2-5(a) (Purple line; "No window" case). Thus, there is only one blue arrow in this graph and the signal frequency shown by a short green solid line matches frequency of that arrow. The vertical dashed blue lines indicate frequencies of other bins. PS at those frequency is zero from the second line of equation (A2-7) because both f1 and f2 of equation (A2-6) become integer multiples of 1/T at those frequencies.

Figure A2-3(c) shows PS and modified PSD for the case B. Since PSD of all the bins are non-zero due to the spectral leakages, we omitted vertical dashed blue lines in this graph. The solid blue line (labeled as "Envelope") connecting tips of the arrows shows that the result shown in this graph is similar to the PSD shown in Figure 2-6(c) (Purple line). Again, modified PSD and PS at the same frequency agree. Figures A2-3(b) and (c) indicate that spectral leakage always occurs but we do not see it on graphs of PSD when we achieved frequency match because the effect of spectral leakage is zero at frequencies where we ordinarily compute PSD for that case. It is slightly difficult to see but the PS of the main lobe is not symmetrical with respect to the center of the main lobe. Due to this asymmetry the bin of the maximum amplitude may not be a signal bin when we cannot achieve frequency match as we described in (2-1-3-3).

Now let us show how the spectral leakage occurs by showing the convolution process graphically. We picked up bins marked by purple solid triangles in figures A2-3(b) and (c) for this purpose. The upper panel of Figure A2-4 shows the Fourier transform of an infinitely long data x and the lower panel of Figure A2-4 shows the Fourier transform of a Boxcar window function.

These are graphical presentations of equations (A2-5) and (A2-7), and the fundamental frequency ranges (horizontal width in the graph) of these two are the same because sampling interval of infinitely long data x(t) and Boxcar window function w(t) must be the same. There is no real part component in the upper graph of Figure A2-4. To calculate PS at frequency fT (purple solid triangle) using convolution, first thing we do is flipping lower graph left and right. This is already done so that frequency (horizontal axis) of the lower graph decreases from left to right unlike upper graph (red circle). The next thing we do is moving this lower graph horizontally so that the horizontal position of the center of lower graph marked by a black solid triangle matches the horizontal position of fT in the upper graph. The result of this manipulation is shown in Figure A2-5.

We omitted unnecessary part of graphs. Then, we draw lines vertically from the upper graph to the lower graph at frequencies f0 and -f0 (dashed lines H and L). Since the frequency difference between H and C is fT+f0=f1 and the frequency difference between L and C is fT-f0=f2, frequencies of these two dashed lines in the lower graph is f1 and f2 (Please see equation (A2-6)). If f1 goes beyond of the edge shown in Figure A2-4, we use periodicity of the Fourier transform as we show later,

As a fourth step, we read values of blue (real part) and red (imaginary part) lines at places where dashed lines H and L intersect those lines. Then, we calculate X multiplied by W as
XW(f1)={Xr(-f0)Wr(f1)-Xi(-f0)Wi(f1)}+i{Xi(-f0)Wr(f1)+Xr(-f0)Wi(f1)}
XW(f2)={Xr(f0)Wr(f2)-Xi(f0)Wi(f2)}+i{Xi(f0)Wr(f2)+Xr(f0)Wi(f2)}
After this we calculate sum of XW(f1) and XW(f2). This is how the convolution process is used to obtain the last part (red underline) of equation (A2-6). Figure A2-5 shows that this computation is the same as computation of weighted moving average (as we described before, in the strict sense it is not an average because the summation of weight is not equal to 1.0) of X with the weight W as we described in (4-2-6). As a final step, we decompose sum of XW(f1) and XW(f2) into a real part and an imaginary part and then substitute the result into equation (A2-3). This will give PS of a Boxcar windowed data at frequency fT,

In case of window-width A shown in Figure A2-5, Wr(f1), Wi(f1), Wr(f2) and Wi(f2) are all zero as we described before. Thus PS at fT is zero from equation (A2-8). Likewise, PS at frequencies where we ordinarily compute PSD except for a signal bin becomes zero and we do not see spectral leakages at those frequencies. However, if we compute PS at any other frequency, we will pick up a signal no matter how far away fT is from f0 (signal frequency). Figure A2-6 shows the case when window-width is B.

In this case, none of Wr(f1), Wi(f1), Wr(f2) and Wi(f2) is zero, and as a result, PSD becomes non-zero even if fT is far away from f0. This is how the spectral leakages occur.

Figure A2-7 shows convolution process for general cases when the Fourier transform of infinitely long data ((a):X) is more complicated than the one we have shown so far.

In this graph we show segment A covering fundamental frequency range and segment B which has a equal frequency width (horizontal width in this graph) of segment A. The contents of these two segments are identical due to the periodicity of X we described before. Let us assume we want to compute PS at frequency fT (purple triangle). We shift horizontal position of the graph of the Fourier transform of a Boxcar window function ((b):W) after flipping it left and right as before. For the convolution we chose segment E of X. The frequency width of segment E is the same as frequency width of segment A as well as frequency width of W. The segment E lacks the segment C of segment A at left but the segment D of segment B added at right compensates this missing segment. Again, due to the periodicity of Fourier transform contents of segments C and D are identical. In this way, the extent of convolution covers entire fundamental frequency range of X no matter where fT exists. We multiply segment E of X by W in the same way we did in figures A2-5 and A2-6; multiply two values at the same horizontal position together. The result of this multiplication is shown in Figure A2-7(c). Then we integrate this result to obtain Yr(fT) and Yi(fT) to compute PS at frequency fT using equation (A2-3). These graphs clearly show that the convolution is equivalent to the weighted moving average (in the strict sense it is not an average because the integral of weight is not equal to 1.0) of X with the weight W and the width of averaging is equal to the fundamental frequency range. As a result, PS at any frequency fT would pick up all the signals on X except for those on frequencies f that satisfy fT - f = +/- k/T (second line of equation (A2-7)).

(3) How to reduce spectral leakages
If we have an infinitely long data shown in Figure A2-1(a), we will not have spectral leakages, but of course, length of our real world data is always limited. Alternatively, if we can achieve a frequency match, we can remove spectral leakages but it is usually impossible to achieve a frequency match when we need to deal with real world data. Then, what we actually can do to reduce spectral leakages is proactively choosing a window function that would reduce spectral leakages instead of accepting an implicit Boxcar window function. As we have shown already, the magnitudes of spectral leakages would decrease if the magnitudes of the Fourier transform of a window function at frequencies other than zero become smaller. As an example of an alternative window function we use a Hanning window function here.

The Fourier transform (discrete) of a Hanning window function is

k in the third line is an integer, greater equal 2 and less equal N/2.

Figures A2-8(a) and (b) show Fourier transform of a Boxcar window function and a Hanning window function, respectively.

While we can see blue and red lines clearly at entire frequency range in Figure A2-8(a), it is difficult to see those lines beyond of +/- 2.0 in Figure A2-8(b). Figure A2-8(c) shows PS of these two windows. Please note that we used a logarithmic scaling for this graph unlike others. From this graph we can see that a Hanning window would considerable reduce spectral leakages and reduction rate would progressively increase as the difference between fT and f0 becomes larger. This graph also shows that, although applying a Hanning window removes 'cliff' described in (4-2-1), we will still have spectral leakages.

Another thing we can see in this graph is that while the PS of a Boxcar window function at frequencies marked by blue asterisks is zero, PS of a Hanning window function is not zero at those frequencies. Also, PS of a Hanning window function at zero frequency (black dashed vertical line) is smaller than the PS of a Boxcar window function at that frequency. These characteristics result the widening of a signal on the graph of PSD when we applied a Hanning window function as we described in (4-2-6).

These results suggest that it would be better if we apply a Hanning window before computing PSD. However, a Hanning window function and many other window functions considerably dump significant portion of data. In certain cases, applying these window functions may cause a problem as we described in (4-2-3).

Appendix 3; Why we see a single peak instead of twin-peaks when data length is short.
We described in the Example 1, (Case 3-1) that when the amplitude of single frequency sinusoid varies sinusoidally we would see a single peak if data length is short and we would see twin-peaks when data length is long. This behavior can be mathematically explained by the theorem called Convolution Theorem, which dictates that Fourier transform of the product of two functions (S0 and S1; first expression of equation (4)) is equal to the convolution of Fourier transform of individual function. As we described in (4-2-6) we compute PSD and amplitude spectrum using results of Fourier transform of data.

When time series is short relative to the period of S1, the maximum of the coefficient of Fourier transform of S1 is at zero frequency because frequency resolution is so large that f1 is included within zero frequency bin. Then, convolution would not modify frequency of S0 and the peak of PSD of our data will be at f0. At this stage, data is practically S0 except that amplitude is not constant, but is slowly varying. Please note that this data length is short relative to the period of S1 but it is fairly long relative to the period of SL and SH as shown in Figure 2-11.

As we increase data length, frequency resolution of the Fourier transform decreases and f1 leaves from zero frequency bin but moves toward higher frequency. Then, convolution of this result and the Fourier transform of S0 generates twin-peaks center of which is at f0. It is not our intention to write any educational material, so we will stop here. We will describe more about convolution in (4-2-6).

Appendix 4; How the frequency of S1 affects the detectability of SR
Let us assume that S1 in (4-1-3) is a sinusoid of amplitude of G and the frequency of it is the same as the frequency of k-th bin (frequency match is achieved). Let us also assume that SR is a sinusoid of amplitude of A and frequency of it is the same as the frequency of m-th bin. The latter assumption makes PSD of SR zero except for m-th bin. In this case the ratio of PSD of SR to the PSD due to added 'inclination' at the m-th bin is, from equation (17) and equation (14),

Now let us assume that the case when PSD due to added 'inclination' significantly affects the detectability of SR is the case when above ratio is less than 1.0. Using equation (A4-1) this condition can be expressed in term of amplitude ratio of SR (A) to S1 (G) as

Above equation shows that if frequencies of both S1 and SR are low (thus both k and m are small), it would be difficult to detect SR even if amplitude of SR is relatively large.

In case of the example we showed in Figure 4-10, k is 8. Substituting this value into (A4-2), the ratio A/G when the PSD of SR becomes equal to the PSD due to added 'inclination' is 0.076/m. This means that if amplitude of S1 (G) is 1.0, amplitude of SR (A) must be 0.076 for m=1 and 0.00076 for m=100 so that PSD of SR becomes equal to the PSD due to added 'inclination'. Since SR is so small compared with S1 in the example shown in Figure 4-10 as described before, the maximum m in Figure 4-10 is still about 100 times smaller to see SR in that figure.

On the other hand, k is 1600 for the example we showed in Figure 4-12. In this case, the ratio A/G when the PSD of SR becomes equal to the PSD due to added 'inclination' is 0.00038/m. If m is greater than about 120 (corresponding frequency is about 74) or so, PSD of SR becomes comparable to the PSD due to added 'inclination'. As a matter of fact we can see SR in the right half of Figure 4-12(d).

END of User Guide of PD001A (Copyright: Cygnus Research International, Apr. 21, 2014)

To read previous page, please click button below.