1 / 36

What are the dominant frequencies?

What are the dominant frequencies?. Fourier transforms decompose a data sequence into a set of discrete spectral estimates – separate the variance of a time series as a function of frequency.

fritz-hester
Download Presentation

What are the dominant frequencies?

An Image/Link below is provided (as is) to download presentation Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author. Content is provided to you AS IS for your information and personal use only. Download presentation by click this link. While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server. During download, if you can't get a presentation, the file might be deleted by the publisher.

E N D

Presentation Transcript

  1. What are the dominant frequencies? Fourier transforms decompose a data sequence into a set of discrete spectral estimates – separate the variance of a time series as a function of frequency. A common use of Fourier transforms is to find the frequency components of a signal buried in a time domain signal.

  2. FAST FOURIER TRANSFORM (FFT) In practice, if the time series f(t) is not a power of 2, it should be padded with zeros

  3. What is the statistical significance of the peaks? Each spectral estimate has a confidence limit defined by a chi-squared distribution

  4. Spectral Analysis Approach 1. Remove mean and trend of time series 2. Pad series with zeroes to a power of 2 3. To reduce end effect (Gibbs’ phenomenon) use a window (Hanning, Hamming, Kaiser) to taper the series 4. Compute the Fourier transform of the series, multiplied times the window 5. Rescale Fourier transform by multiplying times 8/3 for the Hanning Window 6. Compute band-averages or block-segmented averages 7. Incorporate confidence intervals to spectral estimates

  5. 1. Remove mean and trend of time series (N = 1512) 2. Pad series with zeroes to a power of 2 (N = 2048) m Sea level at Mayport, FL July 1, 2007 (day “0” in the abscissa) to September 1, 2007 Raw data and Low-pass filtered data m High-pass filtered data

  6. Spectrum of raw data m2/cpd Spectrum of high-pass filtered data m2/cpd Cycles per day

  7. Hanning Window Hamming Window Value of the Window 3. To reduce end effect (Gibbs’ phenomenon) use a window (Hanning, Hamming, Kaiser) to taper the series Day from July 1, 2007

  8. Hanning Window Hamming Window Kaiser-Bessel, α = 2 Kaiser-Bessel, α = 3 Value of the Window Day from July 1, 2007 3. To reduce end effect (Gibbs’ phenomenon) use a window (Hanning, Hamming, Kaiser) to taper the series

  9. 4. Compute the Fourier transform of the series, multiplied times the window Raw series x Hanning Window (one to one) m Raw series x Hamming Window (one to one) m To reduce side-lobe effects Day from July 1, 2007

  10. 4. Compute the Fourier transform of the series, multiplied times the window High-pass series x Hanning Window (one to one) m High pass series x Hamming Window (one to one) m To reduce side-lobe effects Day from July 1, 2007

  11. High pass series x Kaiser-Bessel Window α=3 (one to one) m Day from July 1, 2007 4. Compute the Fourier transform of the series, multiplied times the window

  12. Windows reduce noise produced by side-lobe effects Noise reduction is effected at different frequencies with Hanning window m2/cpd Original from Raw Data with Hamming window m2/cpd Cycles per day

  13. with Hanning window m2/cpd with Hamming and Kaiser- Bessel (α=3) windows m2/cpd Cycles per day

  14. 5. Rescale Fourier transform by multiplying: times 8/3 for the Hanning Window times 2.5164 for the Hamming Window times ~8/3 for the Kaiser-Bessel (Depending on alpha)

  15. 6. Compute band-averages or block-segmented averages 7. Incorporate confidence intervals to spectral estimates Upper limit: Lower limit: 1-alpha is the confidence (or probability) nu are the degrees of freedom gamma is the ordinate reference value

  16. 0.995 0.990 0.975 0.950 0.900 0.750 0.500 0.250 0.100 0.050 0.025 0.010 0.005 7.88 6.63 5.02 3.84 2.71 1.32 0.45 0.10 0.02 0.00 0.00 0.00 0.00 10.60 9.21 7.38 5.99 4.61 2.77 1.39 0.58 0.21 0.10 0.05 0.02 0.01 12.84 11.34 9.35 7.81 6.25 4.11 2.37 1.21 0.58 0.35 0.22 0.11 0.07 14.86 13.28 11.14 9.49 7.78 5.39 3.36 1.92 1.06 0.71 0.48 0.30 0.21 16.75 15.09 12.83 11.07 9.24 6.63 4.35 2.67 1.61 1.15 0.83 0.55 0.41 18.55 16.81 14.45 12.59 10.64 7.84 5.35 3.45 2.20 1.64 1.24 0.87 0.68 20.28 18.48 16.01 14.07 12.02 9.04 6.35 4.25 2.83 2.17 1.69 1.24 0.99 21.95 20.09 17.53 15.51 13.36 10.22 7.34 5.07 3.49 2.73 2.18 1.65 1.34 23.59 21.67 19.02 16.92 14.68 11.39 8.34 5.90 4.17 3.33 2.70 2.09 1.73 25.19 23.21 20.48 18.31 15.99 12.55 9.34 6.74 4.87 3.94 3.25 2.56 2.16 26.76 24.72 21.92 19.68 17.28 13.70 10.34 7.58 5.58 4.57 3.82 3.05 2.60 28.30 26.22 23.34 21.03 18.55 14.85 11.34 8.44 6.30 5.23 4.40 3.57 3.07 29.82 27.69 24.74 22.36 19.81 15.98 12.34 9.30 7.04 5.89 5.01 4.11 3.57 31.32 29.14 26.12 23.68 21.06 17.12 13.34 10.17 7.79 6.57 5.63 4.66 4.07 32.80 30.58 27.49 25.00 22.31 18.25 14.34 11.04 8.55 7.26 6.26 5.23 4.60 34.27 32.00 28.85 26.30 23.54 19.37 15.34 11.91 9.31 7.96 6.91 5.81 5.14 35.72 33.41 30.19 27.59 24.77 20.49 16.34 12.79 10.09 8.67 7.56 6.41 5.70 37.16 34.81 31.53 28.87 25.99 21.60 17.34 13.68 10.86 9.39 8.23 7.01 6.26 38.58 36.19 32.85 30.14 27.20 22.72 18.34 14.56 11.65 10.12 8.91 7.63 6.84 40.00 37.57 34.17 31.41 28.41 23.83 19.34 15.45 12.44 10.85 9.59 8.26 7.43 41.40 38.93 35.48 32.67 29.62 24.93 20.34 16.34 13.24 11.59 10.28 8.90 8.03 42.80 40.29 36.78 33.92 30.81 26.04 21.34 17.24 14.04 12.34 10.98 9.54 8.64 44.18 41.64 38.08 35.17 32.01 27.14 22.34 18.14 14.85 13.09 11.69 10.20 9.26 45.56 42.98 39.36 36.42 33.20 28.24 23.34 19.04 15.66 13.85 12.40 10.86 9.89 46.93 44.31 40.65 37.65 34.38 29.34 24.34 19.94 16.47 14.61 13.12 11.52 10.52 48.29 45.64 41.92 38.89 35.56 30.43 25.34 20.84 17.29 15.38 13.84 12.20 11.16 49.64 46.96 43.19 40.11 36.74 31.53 26.34 21.75 18.11 16.15 14.57 12.88 11.81 50.99 48.28 44.46 41.34 37.92 32.62 27.34 22.66 18.94 16.93 15.31 13.56 12.46 52.34 49.59 45.72 42.56 39.09 33.71 28.34 23.57 19.77 17.71 16.05 14.26 13.12 53.67 50.89 46.98 43.77 40.26 34.80 29.34 24.48 20.60 18.49 16.79 14.95 13.79 55.00 52.19 48.23 44.99 41.42 35.89 30.34 25.39 21.43 19.28 17.54 15.66 14.46 56.33 53.49 49.48 46.19 42.58 36.97 31.34 26.30 22.27 20.07 18.29 16.36 15.13 57.65 54.78 50.73 47.40 43.75 38.06 32.34 27.22 23.11 20.87 19.05 17.07 15.82 58.96 56.06 51.97 48.60 44.90 39.14 33.34 28.14 23.95 21.66 19.81 17.79 16.50 60.27 57.34 53.20 49.80 46.06 40.22 34.34 29.05 24.80 22.47 20.57 18.51 17.19 Probability 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 Degrees of freedom

  17. N=1512 Includes low frequency

  18. N=1512 Excludes low frequency

  19. N=1512

  20. N=1512

  21. Least Squares Fit to Main Harmonics The observed flow u’ may be represented as the sum of M harmonics: u’ = u0 + ΣjM=1Aj sin (j t + j) For M = 1 harmonic (e.g. a diurnal or semidiurnal constituent): u’ = u0+ A1 sin (1t + 1) With the trigonometric identity: sin (A + B) = cosBsinA + cosAsinB u’ = u0 + a1 sin (1t ) + b1 cos (1t ) taking: a1= A1 cos 1 b1= A1 sin 1

  22. Using u’ = u0 + a1 sin (1t ) + b1 cos (1t ) Then: 2 = ΣN {u 2 - 2uu0 - 2ua1 sin (1t ) - 2ub1 cos (1t ) + u02 + 2u0a1 sin (1t ) + 2u0b1 cos (1t ) + 2a1 b1 sin (1t ) cos (1t ) + a12sin2 (1t ) + b12cos2 (1t ) } The squared errors between the observed current u and the harmonic representation may be expressed as 2 : 2 = ΣN [u - u’ ]2 = u 2 - 2uu’ + u’ 2 Then, to find the minimum distance between observed and theoretical values we need to minimize 2 with respect to u0 a1and b1, i.e., δ 2/ δu0 , δ 2/ δa1 , δ 2/ δb1 : δ2/ δu0 = ΣN {-2u +2u0 + 2a1 sin (1t ) + 2b1 cos (1t ) } = 0 δ2/ δa1 = ΣN { -2u sin (1t ) +2u0 sin (1t ) + 2b1 sin (1t ) cos (1t ) + 2a1 sin2(1t ) } = 0 δ2/ δb1 = ΣN {-2u cos (1t ) +2u0 cos (1t ) + 2a1 sin (1t ) cos (1t ) + 2b1 cos2(1t ) } = 0

  23. ΣN { u = u0 + a1 sin (1t ) + b1 cos (1t ) } ΣN { u sin (1t ) = u0 sin (1t ) + b1 sin (1t ) cos (1t ) + a1 sin2(1t ) } ΣN { u cos (1t ) = u0 cos (1t ) + a1 sin (1t ) cos (1t ) + b1 cos2(1t ) } ΣN u N ΣN sin (1t ) ΣN cos (1t ) u0 ΣN {-2u +2u0 + 2a1 sin (1t ) + 2b1 cos (1t ) } = 0 ΣN u sin (1t ) = ΣN sin (1t ) ΣN sin2(1t ) ΣN sin (1t ) cos (1t ) a1 ΣN u cos (1t ) ΣN cos (1t ) ΣN sin (1t ) cos (1t ) ΣN cos2(1t ) b1 ΣN {-2u sin (1t ) +2u0 sin (1t ) + 2b1 sin (1t ) cos (1t ) + 2a1 sin2(1t ) } = 0 ΣN { -2u cos (1t ) +2u0 cos (1t ) + 2a1 sin (1t ) cos (1t ) + 2b1 cos2(1t ) } = 0 X = A-1 B Rearranging: And in matrix form: B = A X

  24. Finally... The residual or mean is u0 The phase of constituent 1 is: 1 = atan( b1 / a1 ) The amplitude of constituent 1 is: A1 = ( b12+ a12 )½ Pay attention to the arc tangent function used. For example, in IDL you should use atan (b1,a1) and in MATLAB, you should use atan2

  25. Matrix A is then: N ΣN sin (1t ) ΣN cos (1t ) ΣN sin (2t ) ΣN cos (2t ) ΣN sin (1t ) ΣN sin2(1t ) ΣN sin (1t ) cos (1t ) ΣN sin (1t ) sin (2t ) ΣN sin (1t ) cos (2t ) ΣN cos (1t ) ΣN sin (1t ) cos (1t ) ΣN cos2(1t ) ΣN cos (1t ) sin (2t ) ΣN cos (1t ) cos (2t ) ΣN sin (2t ) ΣN sin (1t ) sin (2t ) ΣN cos (1t ) sin (2t ) ΣN sin2(2t ) ΣN sin (2t ) cos (2t ) ΣN cos (2t ) ΣN sin (1t ) cos (2t ) ΣN cos (1t ) cos (2t ) ΣN sin (2t ) cos (2t ) ΣN cos2 (2t ) Remember that: X = A-1 B and B = u0 a1 b1 a2 b2 ΣN u ΣN u sin (1t ) X = ΣN u cos (1t ) ΣN u sin (2t ) ΣN u cos (2t ) For M = 2 harmonics (e.g. diurnal and semidiurnal constituents): u’ = u0+ A1 sin (1t + 1) + A2 sin (2t + 2)

  26. Goodness of Fit: Σ [< uobs > - upred] 2 ------------------------------------- Σ [<uobs > - uobs] 2 Root mean square error: [1/N Σ (uobs - upred) 2] ½

  27. Fit with M2 only

  28. Fit with M2, K1

  29. Fit with M2, S2, K1

  30. Fit with M2, S2, K1, M4, M6

  31. Tidal Ellipse Parameters amplitude of the counter-clockwise rotary component phase of the clockwise rotary component phase of the counter-clockwise rotary component Ellipse Coordinates: ua, va, up, vpare the amplitudes and phases of the east-west and north-south components of velocity amplitude of the clockwise rotary component The characteristics of the tidal ellipses are: Major axis = M = Qcc + Qc minor axis = m = Qac - Qc ellipticity = m / M Phase = -0.5 (thetacc - thetac) Orientation = 0.5 (thetacc + thetac)

  32. M2 S2 K1

  33. Fit with M2 only (dotted) and M2 + M4 (continuous)

  34. Fit with M2 only (dotted) and M2 + M6 (continuous)

More Related