pythonnumpysignal-processingfftphase

Inaccurate phase returned by np.angle


'magnitude_50 Hz': 9.997827675356993, 'phase_50 HZ': -89.0677734968239, 
'magnitude_150 Hz': 4.990392258900833, 'phase_150 HZ': -57.231981462145704,

The magnitude returned is quite close to 10 and 5 respectively. But the phase is not 0 and 30 degrees.

I tried other methods like the math.atan2 and cmath.phase also and it provides similar results.

I would like to understand what is wrong with my phase calculation. My code is below.

def sine_wave(amplitude1: Union[int, float], amplitude2: Union[int, float], phase1: float, phase2: float, duration: Union[int, float],fund_freq_1: int, fund_freq_2: int, samp_freq: int) -> dict:
  
  # generating the time domain signal

  t = np.linspace(0, duration, int(samp_freq * duration))
  wave1 = amplitude1 * np.sin((2 * np.pi * fund_freq_1 * t)+phase1)
  wave2 = amplitude2 * np.sin((2 * np.pi * fund_freq_2 * t)+phase2)
  combined_wave = np.add(wave1, wave2)
  N = combined_wave.size
  T = 1/samp_freq

  # DFT
  f = np.fft.fftfreq(N, 1 / samp_freq)
  fft = np.fft.fft(combined_wave)

  index_one = np.where(np.isclose(f, fund_freq_1))
  magnitude_one = np.mean(np.abs(fft[index_one]) * (2 / N))
  phase_one = degrees(np.angle(fft[index_one]))
  # phase_one = atan2(fft[index_one].imag, fft[index_one].real)
  # phase_one = degrees(phase(fft[index_one]))

  index_two = np.where(np.isclose(f, fund_freq_2))
  magnitude_two = np.mean(np.abs(fft[index_two]) * (2 / N))
  phase_two = degrees(np.angle(fft[index_two]))
  # phase_two = atan2(fft[index_two].imag, fft[index_one].real)
  # phase_two = degrees(phase(fft[index_two]))

  return {'magnitude_{} Hz'.format(fund_freq_1): magnitude_one,
          'phase_{} HZ'.format(fund_freq_1): phase_one,
          'magnitude_{} Hz'.format(fund_freq_2): magnitude_two,
          'phase_{} HZ'.format(fund_freq_2): phase_two}

The code could be run like this

sine_wave(amplitude1=10, amplitude2=5, phase1=0, phase2=np.pi/6, duration=0.1, fund_freq_1=50, fund_freq_2=150, samp_freq=10000)

Solution

  • After performing the FFT the phase of the complex values correspond to the relative phase with a cosine. Since cos(x) has a 90 degrees phase difference with sin(x) you should expect your 0-degrees-phase sin to be detected with a phase of -90 degrees with respect to the corresponding cos at the same frequency. Similarly your 30-degrees-phase sin should be detected with a phase of -60 degrees. Your resulting values are indeed quite close.

    If you prefer to get the phase referenced to sin signals, then you may simply add 90 degrees to the result of np.angle:

    phase_one = degrees(np.angle(fft[index_one])) + 90
    phase_two = degrees(np.angle(fft[index_two])) + 90