Problem 393932 · medium · Level 03 Linear Management & Searching

The Spectrum of a Motor Current

discrete Fourier transform · magnitude · phase · amplitude spectrum · cmath

A technician clamps a current probe on the cable of an electric motor and records N samples at fs samples per second. A healthy motor draws a steady current plus a few clean waves; a fault adds waves at new frequencies. The technician wants a short report of the waves the recording contains.

For every bin k, the DFT gives

X[k] = sum over n = 0 .. N-1 of  x[n] * exp(-2j * pi * k * n / N)

and it describes a wave at k * fs / N hertz. If the recording contains a * cos(2 * pi * k * n / N + p), the bin is X[k] = (a * N / 2) * exp(1j * p): its magnitude gives the amplitude a = 2 * |X[k]| / N and its angle gives the phase p, where the wave is in its cycle at the first sample. Bin 0 is special: X[0] is just the sum of the samples, so X[0].real / N is the offset (the steady current).

Write motor_spectrum(x, fs, threshold) that returns a tuple (offset, peaks):

  • offset is X[0].real / N,
  • peaks is a list with one tuple (frequency, amplitude, phase) for every bin k with 1 <= k < N / 2 whose amplitude 2 * |X[k]| / N is at least threshold, in order of k. The frequency is in hertz and the phase in degrees, between -180 and 180: math.degrees(cmath.phase(X[k])).

Press Run with plot([k * fs / N for k in range(N // 2)], [2 * abs(X[k]) / N for k in range(N // 2)], kind="stem") to see the spectrum once you have the list X.

Examples

Input:  x = [3, 1, -1, 1], fs = 4, threshold = 0.1
Output: (1.0, [(1.0, 2.0, 0.0)])
Explanation: the samples are 1 + 2 * cos(2 * pi * n / 4): an offset of 1 and a wave of
amplitude 2 at 1 Hz that starts at its peak (phase 0).

Input:  x = [0, 1, 0, -1], fs = 400, threshold = 0.1
Output: (0.0, [(100.0, 1.0, -90.0)])
Explanation: a sine is a cosine a quarter of a cycle late, so its phase is -90 degrees.

Constraints

  • 4 <= N <= 512; answers are compared with a tolerance of 1e-6
  • the tests keep every reported phase away from the jump between -180 and 180 degrees, and every amplitude well away from the threshold

Goals

  • Compute every useful bin of the DFT with two nested loops
  • Read a bin's magnitude as an amplitude and its angle as a phase
  • Report only the bins that stand above a noise threshold
Starting Python…