#!/usr/bin/env python3

# AD5933 measurement data analysis
#
# 2024-08-05 ch
#
import numpy as np
import matplotlib.pyplot as plt

from numpy import genfromtxt

calib = genfromtxt('calib.txt', delimiter=',')
data = genfromtxt('data.txt', delimiter=',')

#print(data)

f= data[1:,0]        
real= data[1:,1]
imag= data[1:,2]
realc= calib[1:,1]
imagc= calib[1:,2]

Rref=100000

phase=np.arctan2(real,imag)/np.pi*180
phasec=np.arctan2(realc,imagc)/np.pi*180

magnitudec=np.sqrt(realc*realc+imagc*imagc)
magnitude=np.sqrt(real*real+imag*imag)
gain=1/(Rref*magnitudec)
impedance=1/(gain*magnitude)

plt.subplot(221)
plt.suptitle('AD5933 measurement (Rref='+str(Rref)+"Ohm)", fontsize=16)
plt.title( "calib raw signals")
plt.plot(f, realc)
plt.plot(f, imagc)
plt.xlabel("f[Hz]");
plt.ylabel("[LSB]");
plt.legend(['real', 'imag'])
plt.grid()
plt.tight_layout()

plt.subplot(222)
plt.title( "raw signals")
plt.plot(f, real)
plt.plot(f, imag)
plt.xlabel("f[Hz]");
plt.ylabel("[LSB]");
plt.legend(['real', 'imag'])
plt.grid()
plt.tight_layout()

plt.subplot(223)
plt.title( "phase ")
plt.plot(f,phasec)
plt.plot(f,phase)
phaseCorrected=phase-phasec;
plt.plot(f,phaseCorrected)
plt.xlabel("f[Hz]");
plt.ylabel("phase[°]");
plt.legend(['calib', 'raw meas','result'])
plt.ylim(-200,200)
plt.grid()
plt.tight_layout()

plt.subplot(224)
plt.title( "impedance")
plt.plot(f, impedance)

plt.xlabel("f[Hz]");
plt.ylabel("|Z| [Ohm]");
#plt.legend(['real', 'imag'])
plt.grid()
plt.tight_layout()

plt.show()