import numpy as np import matplotlib.pyplot as plot from scipy import signal, constants import argparse import os import datetime def rp_raw_to_numpy(rp_raw): # Convert raw buffer strings to numpy arrays buff_string = rp_raw.split(',') return np.array(list(map(float, buff_string))) RP_IP_ADDRS = { "creotech-1": "192.168.1.104", "creotech-2": "192.168.1.105", "mlabs-1": "rp-f05cc9", "mlabs-2": "rp-f0612e", } def main(): parser = argparse.ArgumentParser(description="Data plotting tool for Sayma DAC/TTL") parser.add_argument("dir", help="output directory", type=str) parser.add_argument("rps", metavar="NAME:CHANNEL", help="input of the RedPitaya at where data is collected; " "CHANNEL must be 1 or 2; " "NAME must be any of: " + " ".join(list(RP_IP_ADDRS.keys())), type=str, nargs=2) parser.add_argument("--log", help="path of the log file where the measurement record will be appended", type=str) parser.add_argument("--noplot", help="do not show data plot, which is blocking until the GUI is closed", action="store_true") args = parser.parse_args() # Must only compare 2 data name, channel = args.rps[0].split(':') y1_filename = 'rp_{}_y{}_raw.bin'.format(name, channel) name, channel = args.rps[1].split(':') y2_filename = 'rp_{}_y{}_raw.bin'.format(name, channel) if y1_filename == y2_filename: raise ValueError("Both files are the same.") with open(os.path.join(args.dir, y1_filename), 'rb') as f: y1_now_iso, y1_raw = [l.decode('utf-8') for l in f.readlines()] with open(os.path.join(args.dir, y2_filename), 'rb') as f: y2_now_iso, y2_raw = [l.decode('utf-8') for l in f.readlines()] if None in [y1_raw, y2_raw]: raise IOError("Raw RP string files cannot be opened.") if y1_now_iso != y2_now_iso: raise ValueError("Timestamps of raw RP files are not identical.") now_iso = y1_now_iso.rstrip() print("Reading raw RP data collected at {}.".format(now_iso)) y1 = rp_raw_to_numpy(y1_raw) y2 = rp_raw_to_numpy(y2_raw) # Define t as an array of timestamps (in seconds) for each sample # (Note that RedPitaya's oscilloscope has a sampling rate @ 125MHz) t = np.arange(y1.shape[0])/125e6 # Generate matrix Y by having arrays y1 and y2 as 2 rows y = np.c_[y1, y2].T # Element-wise multiply an array of cos(2pi*9e6*t) with each row in Y; # Then, downsample the array by 10 as Z z = signal.decimate(y*np.exp(1j*2*np.pi*9e6*t), q=10, ftype="fir", zero_phase=True)[:, 10:] # Downsample Z by 10 again. z = signal.decimate(z, q=10, ftype="fir", zero_phase=True)[:, 10:] # Element-wise multiply Z[0] with the conjugate of Z[1] to get the phase difference (i.e. angle(z0) - angle(z1)), and use the mean value. angle = np.angle(np.mean(z[0]*z[1].conj())) # Append the phase difference to the log file log = args.log if log is not None: with open(log, 'a') as f: f.write("{}\t{}\n".format(now_iso, angle)) print("Phase measurement record appended to log: {}".format(log)) # Print the phase difference print(angle) if not args.noplot: # Normalize y1 and y2 for plotting y1 /= abs(y1).max() y2 /= abs(y2).max() plot.plot(y1) plot.plot(y2) plot.show() if __name__ == "__main__": main()