{ "cells": [ { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import struct\n", "import random\n", "import itertools\n", "import datetime\n", "import multiprocessing\n", "from collections import defaultdict\n", "import json\n", "import traceback\n", "import glob\n", "\n", "from matplotlib import pyplot as plt\n", "import matplotlib\n", "from matplotlib import ticker\n", "import numpy as np\n", "from scipy import signal as sig\n", "from scipy import fftpack as fftpack\n", "import ipywidgets\n", "from IPython.display import set_matplotlib_formats\n", "import scipy\n", "\n", "from tqdm.notebook import tqdm\n", "import colorednoise\n", "import accupy\n", "\n", "np.set_printoptions(linewidth=240)" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "%matplotlib widget\n", "#%matplotlib inline\n", "#set_matplotlib_formats('png', 'pdf')\n", "#font = {'family' : 'normal',\n", "# 'weight' : 'normal',\n", "# 'size' : 6}\n", "#matplotlib.rc('font', **font)" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "sampling_rate = 10 # sp/s" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "# From https://github.com/mubeta06/python/blob/master/signal_processing/sp/gold.py\n", "preferred_pairs = {5:[[2],[1,2,3]], 6:[[5],[1,4,5]], 7:[[4],[4,5,6]],\n", " 8:[[1,2,3,6,7],[1,2,7]], 9:[[5],[3,5,6]], \n", " 10:[[2,5,9],[3,4,6,8,9]], 11:[[9],[3,6,9]]}\n", "\n", "def gen_gold(seq1, seq2):\n", " gold = [seq1, seq2]\n", " for shift in range(len(seq1)):\n", " gold.append(seq1 ^ np.roll(seq2, -shift))\n", " return gold\n", "\n", "def gold(n):\n", " n = int(n)\n", " \n", " min_key, max_key = min(preferred_pairs.keys()), max(preferred_pairs.keys())\n", " \n", " if n < min_key:\n", " raise KeyError('preferred pairs for %s bits unknown' % str(n))\n", " \n", " if n > max_key:\n", " t0, t1 = preferred_pairs[max_key]\n", " (seq0, _st0), (seq1, _st1) = sig.max_len_seq(max_key, taps=t0), sig.max_len_seq(max_key, taps=t1)\n", " raw_codes = gen_gold(seq0, seq1)\n", " # Hackety hack: This is super-experimental test code and guaranteed sub-optimal!\n", " \n", " raw_codes = np.array(raw_codes)\n", " out = raw_codes[:-1]\n", " for i in range(n - max_key):\n", " out = np.hstack([out[:out.shape[0]//2], out[out.shape[0]//2:]])\n", " return out\n", " \n", " else:\n", " t0, t1 = preferred_pairs[n]\n", " (seq0, _st0), (seq1, _st1) = sig.max_len_seq(n, taps=t0), sig.max_len_seq(n, taps=t1)\n", " return np.array(gen_gold(seq0, seq1))" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "def modulate(data, nbits=5, pad=True):\n", " # 0, 1 -> -1, 1\n", " mask = np.array(gold(nbits))*2 - 1\n", " \n", " sel = mask[data>>1]\n", " data_lsb_centered = ((data&1)*2 - 1)\n", "\n", " signal = (np.multiply(sel, np.tile(data_lsb_centered, (sel.shape[1], 1)).T).flatten() + 1) // 2\n", " if pad:\n", " print('mod: padding by', mask.shape[1])\n", " return np.hstack([ np.zeros(mask.shape[1]), signal, np.zeros(mask.shape[1]) ])\n", " else:\n", " return signal" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": [ "def correlate(sequence, nbits=5, decimation=1, mask_filter=lambda x: x, mask_max=None):\n", " codes = gold(nbits)\n", " mask = np.tile(np.array(codes)[:,:,np.newaxis]*2 - 1, (1, 1, decimation)).reshape((codes.shape[0], codes.shape[1] * decimation))\n", " mask = mask.astype(np.float64)\n", "\n", " # Our input signal has large DC bias. Remove DC bias to reduce numerical errors during correlation.\n", " sequence -= np.mean(sequence)\n", " \n", " #return np.array([np.correlate(sequence, row, mode='valid') for row in mask[:mask_max]])\n", " return np.array([sig.correlate(sequence, row, mode='valid') for row in mask[:mask_max]]) / len(sequence)" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "mean: 50.00341 len: 1946174\n" ] } ], "source": [ "with open('data/fmeas_export_ocxo_2day.bin', 'rb') as f:\n", " meas_data = np.copy(np.frombuffer(f.read(), dtype='float32'))\n", " print('mean:', np.mean(meas_data), 'len:', len(meas_data))\n", " meas_data -= np.mean(meas_data)\n", " meas_data = np.hstack([meas_data, meas_data[::-1]] * 50)\n", " def mains_noise_meas(n):\n", " last_valid = len(meas_data) - n\n", " start = random.randint(0, last_valid)\n", " return meas_data[start:start+n]" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "def load_noise_meas_params(capture_file):\n", " with open(capture_file, 'rb') as f:\n", " meas_data = np.copy(np.frombuffer(f.read(), dtype='float32'))\n", " print('mean:', np.mean(meas_data), 'len:', len(meas_data))\n", " meas_data -= np.mean(meas_data)\n", " return {'meas_data': meas_data}\n", "\n", "def mains_noise_measured(n, meas_data):\n", " last_valid = len(meas_data) - n\n", " start = np.random.randint(last_valid)\n", " return meas_data[start:start+n]\n", "\n", "def load_noise_synth_params(specfile):\n", " with open(specfile) as f:\n", " d = json.load(f)\n", " return {'spl_x': np.linspace(*d['x_spec']),\n", " 'spl_N': d['x_spec'][2],\n", " 'psd_spl': (d['t'], d['c'], d['k']) }\n", " \n", "def mains_noise_synthetic(n, psd_spl, spl_N, spl_x):\n", " noise = np.random.normal(size=spl_N) * 2\n", " spec = scipy.fftpack.fft(noise) **2\n", " spec *= np.exp(scipy.interpolate.splev(spl_x, psd_spl))\n", "\n", " spec **= 1/2\n", " \n", " renoise = scipy.fftpack.ifft(spec)\n", " return np.hstack([renoise[10000:], renoise[10000:][::-1]])[:n]\n", "\n", "noise_params = load_noise_synth_params('data/grid_freq_psd_spl_108pt.json')\n", "mains_noise = lambda n: mains_noise_synthetic(n, **noise_params)" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [], "source": [ "def generate_test_signal(duration, nbits=6, signal_amplitude=2.0e-3, decimation=10, seed=0, data=None):\n", " test_data = np.random.RandomState(seed=seed).randint(0, 2 * (2**nbits), duration) if data is None else data\n", " \n", " signal = np.repeat(modulate(test_data, nbits, pad=True) * 2.0 - 1, decimation) * signal_amplitude\n", " noise = mains_noise_meas(len(signal))\n", " \n", " return test_data, signal.astype(np.longdouble) + noise.astype(np.longdouble)" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([0.50806058, 0.49242794])" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "gold(12)[:2].mean(axis=1)" ] }, { "cell_type": "code", "execution_count": 20, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "amplitude @P_mod=10000.0kW: Δf=400.0µHz\n", "mod: padding by 1023\n", "[0 1 2 3] (49104,)\n" ] }, { "data": { "application/vnd.jupyter.widget-view+json": { "model_id": "70c8a0e9c03441fab580b5b6b219e065", "version_major": 2, "version_minor": 0 }, "text/plain": [ "Canvas(toolbar=Toolbar(toolitems=[('Home', 'Reset original view', 'home', 'home'), ('Back', 'Back to previous …" ] }, "metadata": {}, "output_type": "display_data" }, { "name": "stdout", "output_type": "stream", "text": [ "signal strength (off/on peak): +4.900e-05 | -6.406e-04\n", "signal strength (off/on peak): +7.699e-05 | +5.627e-04\n", "signal strength (off/on peak): +8.441e-05 | +4.206e-05\n", "signal strength (off/on peak): +8.178e-05 | -8.413e-05\n" ] } ], "source": [ "nbits = 10\n", "decimation = 8\n", "# working params: decimation=15 nbits=21 power=200kW\n", "# working params: decimation=8 nbits=10 power=10MW\n", "grid_frequency_characteristic = 25e9 # watt per hertz\n", "modulation_power = 10000e3 # watt\n", "data = np.array([0, 1, 2, 3])\n", "amplitude = modulation_power / grid_frequency_characteristic\n", "print(f'amplitude @P_mod={modulation_power/1e3}kW: Δf={amplitude*1e6}µHz')\n", "\n", "# Generate test data for this simulation run\n", "test_data, signal = generate_test_signal(duration=2, nbits=nbits, signal_amplitude=amplitude, decimation=decimation, seed=0, data=data)\n", "print(test_data, signal.shape)\n", "\n", "def correlate_test(sequence, nbits, decimation, mask_max):\n", " codes = gold(nbits)\n", " mask = np.tile(np.array(codes)[:,:,np.newaxis]*2 - 1, (1, 1, decimation)).reshape((codes.shape[0], codes.shape[1] * decimation))\n", " mask = mask.astype(np.float64)\n", "\n", " # Our input signal has large DC bias. Remove DC bias to reduce numerical errors during correlation.\n", " sequence -= np.mean(sequence)\n", " \n", " #return np.array([np.correlate(sequence, row, mode='valid') for row in mask[:mask_max]])\n", " return np.array([accupy.fsum(np.multiply(sequence, row)) for row in mask[:mask_max]]) / len(sequence)\n", "\n", "signal = signal.astype(np.float64)\n", "fig, axs = plt.subplots(len(data), figsize=(12, 5))\n", "for i, ax in enumerate(axs.flatten()):\n", " gen_per = len(data) + 2 # one period for each symbol, plus one period of padding at the beginning and end\n", " symbol_length = signal.shape[0]//gen_per\n", " l, r = int(max(0, (i + 1) * symbol_length - 1000)), int((i + 2) * symbol_length + 1000)\n", " signal_slice = signal[l:r]\n", " cor_an = correlate(signal_slice, nbits=nbits, decimation=decimation, mask_max=2)\n", " \n", " l, r = int((i + 1) * symbol_length), int((i + 2) * symbol_length)\n", " signal_slice = signal[l:r]\n", " #print('on peak [numpy]:', correlate(signal_slice, nbits=nbits, decimation=decimation, mask_max=2).flatten() / amplitude)\n", " #print('on peak [accupy]:', correlate_test(signal_slice, nbits=nbits, decimation=decimation, mask_max=2).flatten() / amplitude)\n", " \n", " l, r = int((i + 1 - 0.2) * symbol_length), int((i + 2 - 0.2) * symbol_length)\n", " signal_slice = signal[l:r]\n", " #print('off peak [numpy]:', correlate(signal_slice, nbits=nbits, decimation=decimation, mask_max=2).flatten() / amplitude)\n", " #print('off peak [accupy]:', correlate_test(signal_slice, nbits=nbits, decimation=decimation, mask_max=2).flatten() / amplitude)\n", " #ax.plot(cor_an.T)\n", "\n", " cwt_res = np.array([ sig.cwt(row, sig.ricker, [0.73 * decimation]).flatten() for row in cor_an ])\n", " print(f'signal strength (off/on peak): {np.sqrt(np.mean(cwt_res[0,cwt_res.shape[1]//6:cwt_res.shape[1]//3]**2)):+.3e} | {cwt_res[0,cwt_res.shape[1]//2]:+.3e}')\n", " ax.plot(cwt_res.T / amplitude)" ] }, { "cell_type": "code", "execution_count": 127, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "amplitude @P_mod=10.0kW: Δf=0.39999999999999997µHz\n", "mod: padding by 31\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ ":7: ComplexWarning: Casting complex values to real discards the imaginary part\n", " return test_data, signal.astype(np.longdouble) + noise.astype(np.longdouble)\n" ] }, { "name": "stdout", "output_type": "stream", "text": [ "mod: padding by 31\n", "mod: padding by 31\n", "mod: padding by 31\n", "mod: padding by 31\n", "mod: padding by 31\n", "mod: padding by 31\n", "mod: padding by 31\n", "mod: padding by 31\n", "mod: padding by 31\n", "mod: padding by 31\n", "mod: padding by 31\n", "mod: padding by 31\n", "mod: padding by 31\n", "mod: padding by 31\n", "mod: padding by 31\n", "mod: padding by 31\n", "mod: padding by 31\n", "mod: padding by 31\n", "mod: padding by 31\n", "[0.00031921745429844274216, 0.0003466253347083378983, 0.00044946868526334726654, 0.00043601785453840703518, 0.00027364766834844621423, 0.00034429714350871021226, 0.00029257480354291570764, 0.00043761373963434835483, 0.00034653346037204395042, 0.00041369145870726819648, 0.00031663491690404936546, 0.0003630671579533060333, 0.0002934121262870590337, 0.00036264881740159871053, 0.00033973548163424448722, 0.00045000162460528714934, 0.00036345615545028226934, 0.00034203880459996953548, 0.0004340491996267533479, 0.00027402160564997715637]\n", "mod: padding by 63\n", "mod: padding by 63\n", "mod: padding by 63\n", "mod: padding by 63\n", "mod: padding by 63\n", "mod: padding by 63\n", "mod: padding by 63\n", "mod: padding by 63\n", "mod: padding by 63\n", "mod: padding by 63\n", "mod: padding by 63\n", "mod: padding by 63\n", "mod: padding by 63\n", "mod: padding by 63\n", "mod: padding by 63\n", "mod: padding by 63\n", "mod: padding by 63\n", "mod: padding by 63\n", "mod: padding by 63\n", "mod: padding by 63\n", "[0.00043183462240684953605, 0.00049995218014044821913, 0.00034320159159449236953, 0.00047574377773946273916, 0.00042495372363872566055, 0.0005340590070278771213, 0.00043128485159190824515, 0.00034867807874549918169, 0.00039072607874113340072, 0.00034884834981040312958, 0.00042730792396368059213, 0.0005043675584057222778, 0.0003499736642140912821, 0.00021263777836169804578, 0.00023831087214244111956, 0.00030008253740133703718, 0.00021268569995328232606, 0.00054148296194613603103, 0.00040150545723533495933, 0.00028661636535881167316]\n", "mod: padding by 127\n", "mod: padding by 127\n", "mod: padding by 127\n", "mod: padding by 127\n", "mod: padding by 127\n", "mod: padding by 127\n", "mod: padding by 127\n", "mod: padding by 127\n", "mod: padding by 127\n", "mod: padding by 127\n", "mod: padding by 127\n", "mod: padding by 127\n", "mod: padding by 127\n", "mod: padding by 127\n", "mod: padding by 127\n", "mod: padding by 127\n", "mod: padding by 127\n", "mod: padding by 127\n", "mod: padding by 127\n", "mod: padding by 127\n", "[0.00026772134709063096758, 0.0004023169589146378664, 0.00052784935496986912216, 0.00033900150125753280712, 0.00022029107067617116654, 0.00041990310449284588383, 0.0004415131846159276476, 0.00040143078518334457215, 0.00031109331861264841798, 0.0002811855751776715795, 0.00044570966535767569383, 0.0003301271164337369956, 0.00033158808534416807198, 0.00026578644644711210834, 0.00034742527826464402753, 0.0005426992707322571906, 0.00044335106857497668567, 0.00032225907153724712904, 0.0002813426161001733043, 0.00041113085027627075947]\n", "mod: padding by 255\n", "mod: padding by 255\n", "mod: padding by 255\n", "mod: padding by 255\n", "mod: padding by 255\n", "mod: padding by 255\n", "mod: padding by 255\n", "mod: padding by 255\n", "mod: padding by 255\n", "mod: padding by 255\n", "mod: padding by 255\n", "mod: padding by 255\n", "mod: padding by 255\n", "mod: padding by 255\n", "mod: padding by 255\n", "mod: padding by 255\n", "mod: padding by 255\n", "mod: padding by 255\n", "mod: padding by 255\n", "mod: padding by 255\n", "[0.00019884404254742679324, 0.000242880777230771204, 0.00019589281771018505908, 0.00022349702779183853734, 0.00020954727853424042772, 0.00020473570231859194066, 0.00020169541990839266527, 0.00025024395700750437305, 0.0002751261591544219094, 0.00019445654879589896521, 0.00020805938128876005823, 0.00020078107591647889947, 0.00019010551607677558, 0.00023584019352316058905, 0.0002678021299787742911, 0.00017845149465711736266, 0.00016695541802082155936, 0.00016698843247300128494, 0.00014251122213695801984, 0.00028090032636289179365]\n", "mod: padding by 511\n", "mod: padding by 511\n", "mod: padding by 511\n", "mod: padding by 511\n", "mod: padding by 511\n", "mod: padding by 511\n", "mod: padding by 511\n", "mod: padding by 511\n", "mod: padding by 511\n", "mod: padding by 511\n", "mod: padding by 511\n", "mod: padding by 511\n", "mod: padding by 511\n", "mod: padding by 511\n", "mod: padding by 511\n", "mod: padding by 511\n", "mod: padding by 511\n", "mod: padding by 511\n", "mod: padding by 511\n", "mod: padding by 511\n", "[0.00019975970170629950027, 0.00013186823465058419483, 0.00023953255381256545738, 0.00017255565602586637237, 0.00018396445334660804652, 0.00014540194644455677518, 0.00019876358084066945335, 0.0001738424376447688047, 0.0001765109371402779291, 0.00018196424414407873752, 0.00023044940092346266965, 0.00019731842839493123688, 0.00023705358802857269644, 0.00015328320915330276564, 0.00017778635677371888602, 0.00021212946232256081042, 0.00026136191723243009158, 0.00016766188832612281848, 0.000192498527508538207, 0.00024108248581102643531]\n", "mod: padding by 1023\n", "mod: padding by 1023\n", "mod: padding by 1023\n", "mod: padding by 1023\n", "mod: padding by 1023\n", "mod: padding by 1023\n", "mod: padding by 1023\n", "mod: padding by 1023\n", "mod: padding by 1023\n", "mod: padding by 1023\n", "mod: padding by 1023\n", "mod: padding by 1023\n", "mod: padding by 1023\n", "mod: padding by 1023\n", "mod: padding by 1023\n", "mod: padding by 1023\n", "mod: padding by 1023\n", "mod: padding by 1023\n", "mod: padding by 1023\n", "mod: padding by 1023\n", "[0.00013530662801083781384, 0.0001497806916396048579, 0.00013664959261312633858, 0.00016589078984916737442, 0.00017011538645283205862, 0.00015336800379842749231, 0.0001694272184760911441, 0.00013875282881629243428, 0.000112969749633285926445, 0.00014931796422008114573, 0.0001372035855968592536, 0.00012312327359857258514, 0.00013905194521669023908, 0.00013901001422022015822, 0.00016339317685934381651, 0.00015711672887077335837, 0.00012437507736796685456, 0.00012540715390246491552, 0.00010657701913677854194, 0.00014676391008527084711]\n", "mod: padding by 2047\n", "mod: padding by 2047\n", "mod: padding by 2047\n", "mod: padding by 2047\n", "mod: padding by 2047\n", "mod: padding by 2047\n", "mod: padding by 2047\n", "mod: padding by 2047\n", "mod: padding by 2047\n", "mod: padding by 2047\n", "mod: padding by 2047\n", "mod: padding by 2047\n", "mod: padding by 2047\n", "mod: padding by 2047\n", "mod: padding by 2047\n", "mod: padding by 2047\n", "mod: padding by 2047\n", "mod: padding by 2047\n", "mod: padding by 2047\n", "mod: padding by 2047\n", "[9.397038360574742737e-05, 9.414892638810825941e-05, 8.643228267549556685e-05, 8.4662709794848760505e-05, 7.018825855548282147e-05, 9.3962524208973514975e-05, 8.850854906632382727e-05, 9.490061324129898401e-05, 9.9421482075907846084e-05, 9.87805754535881666e-05, 8.77129165764244642e-05, 8.627947590758716296e-05, 9.310894721386338135e-05, 0.00010087838856131649663, 7.6747263039331168876e-05, 8.2710693121826614583e-05, 0.000106444803775372911015, 9.682856374335765839e-05, 8.172981177790056443e-05, 9.233645079286106117e-05]\n", "mod: padding by 4094\n", "mod: padding by 4094\n", "mod: padding by 4094\n", "mod: padding by 4094\n", "mod: padding by 4094\n", "mod: padding by 4094\n", "mod: padding by 4094\n", "mod: padding by 4094\n", "mod: padding by 4094\n", "mod: padding by 4094\n", "mod: padding by 4094\n", "mod: padding by 4094\n", "mod: padding by 4094\n", "mod: padding by 4094\n", "mod: padding by 4094\n", "mod: padding by 4094\n", "mod: padding by 4094\n", "mod: padding by 4094\n", "mod: padding by 4094\n", "mod: padding by 4094\n", "[7.291179023159128505e-05, 6.921653169518604612e-05, 7.0212257862804909805e-05, 7.739735101427627361e-05, 7.05170185815625493e-05, 7.667140441838123699e-05, 6.684245257693558025e-05, 6.7274435477192927474e-05, 7.3140978179672748096e-05, 6.462540379490416508e-05, 7.8575870371873889974e-05, 7.0962264879147295656e-05, 6.887905626735870121e-05, 6.6741332451462875345e-05, 6.774053852476256235e-05, 8.4915545476866071356e-05, 6.755406828145308233e-05, 7.4150959965664274745e-05, 7.289475345040078811e-05, 7.417472792999099321e-05]\n", "mod: padding by 8188\n", "mod: padding by 8188\n", "mod: padding by 8188\n", "mod: padding by 8188\n", "mod: padding by 8188\n", "mod: padding by 8188\n", "mod: padding by 8188\n", "mod: padding by 8188\n", "mod: padding by 8188\n", "mod: padding by 8188\n", "mod: padding by 8188\n", "mod: padding by 8188\n", "mod: padding by 8188\n", "mod: padding by 8188\n", "mod: padding by 8188\n", "mod: padding by 8188\n", "mod: padding by 8188\n", "mod: padding by 8188\n", "mod: padding by 8188\n", "mod: padding by 8188\n", "[4.809424705445009225e-05, 5.3333793693925671324e-05, 5.6174142194133624293e-05, 5.070037019380478062e-05, 4.2797138004346089653e-05, 5.6184084783894627412e-05, 4.8948695230438478363e-05, 5.231607985651584901e-05, 4.8563366184299055167e-05, 5.0229443138648695e-05, 5.140676133340082623e-05, 5.1571059935388156693e-05, 4.799175351462962519e-05, 5.0143994275182861474e-05, 5.176115169674431516e-05, 5.4781771632211411004e-05, 5.1420164246213953733e-05, 5.584101664179133697e-05, 5.2387951249217037658e-05, 4.8454927642515831705e-05]\n", "mod: padding by 16376\n", "mod: padding by 16376\n", "mod: padding by 16376\n", "mod: padding by 16376\n", "mod: padding by 16376\n", "mod: padding by 16376\n", "mod: padding by 16376\n", "mod: padding by 16376\n", "mod: padding by 16376\n", "mod: padding by 16376\n", "mod: padding by 16376\n", "mod: padding by 16376\n", "mod: padding by 16376\n", "mod: padding by 16376\n", "mod: padding by 16376\n", "mod: padding by 16376\n", "mod: padding by 16376\n", "mod: padding by 16376\n", "mod: padding by 16376\n", "mod: padding by 16376\n", "[3.441001692841763238e-05, 3.9226152369629212293e-05, 3.5174358297006928404e-05, 3.858339666154567386e-05, 3.729609282512620058e-05, 3.8017607168678973304e-05, 3.6432437134237636557e-05, 3.5396590548610446987e-05, 4.2872778372966043538e-05, 3.881954270914244008e-05, 3.6638457035791679994e-05, 4.3987686908835178038e-05, 3.5810535332467513793e-05, 3.669588547824606978e-05, 3.9034162252655188724e-05, 3.686462573963708382e-05, 3.64018808798140893e-05, 4.140879713977473717e-05, 3.601891068378439981e-05, 4.2239933297790315508e-05]\n", "mod: padding by 32752\n", "mod: padding by 32752\n", "mod: padding by 32752\n", "mod: padding by 32752\n", "mod: padding by 32752\n", "mod: padding by 32752\n", "mod: padding by 32752\n", "mod: padding by 32752\n", "mod: padding by 32752\n", "mod: padding by 32752\n", "mod: padding by 32752\n", "mod: padding by 32752\n", "mod: padding by 32752\n", "mod: padding by 32752\n", "mod: padding by 32752\n", "mod: padding by 32752\n", "mod: padding by 32752\n", "mod: padding by 32752\n", "mod: padding by 32752\n", "mod: padding by 32752\n", "[3.0164568033183993325e-05, 3.340002097279428693e-05, 3.1692825762129764792e-05, 3.2341558069057274648e-05, 3.5104645491986034208e-05, 3.3804378516131733095e-05, 3.0847430711483434965e-05, 3.163598515247249634e-05, 3.1877253973163090688e-05, 3.145134169534852806e-05, 2.9469157225219549716e-05, 3.5022547350051354115e-05, 3.185664970143889292e-05, 3.2250533412830754613e-05, 3.2968988013916950676e-05, 3.1167159682343826453e-05, 3.0519324313846599056e-05, 3.6983652848933391153e-05, 3.4137279541673856477e-05, 3.4656056379915872357e-05]\n" ] } ], "source": [ "decimation = 10\n", "grid_frequency_characteristic = 25e9 # watt per hertz\n", "modulation_power = 10e3 # watt\n", "amplitude = modulation_power / grid_frequency_characteristic\n", "print(f'amplitude @P_mod={modulation_power/1e3}kW: Δf={amplitude*1e6}µHz')\n", "\n", "ps_off = []\n", "xs = np.array(list(range(5, 16)))\n", "for nbits in xs:\n", " data = np.array([1])\n", "\n", " def correlate_test(sequence, nbits, decimation, mask_max):\n", " codes = gold(nbits)\n", " mask = np.tile(np.array(codes)[:,:,np.newaxis]*2 - 1, (1, 1, decimation)).reshape((codes.shape[0], codes.shape[1] * decimation))\n", " mask = mask.astype(np.float64)\n", "\n", " # Our input signal has large DC bias. Remove DC bias to reduce numerical errors during correlation.\n", " sequence -= np.mean(sequence)\n", "\n", " return np.array([sig.correlate(sequence, row, mode='valid') for row in mask[:mask_max]])\n", " #return np.array([accupy.fsum(np.multiply(sequence, row)) for row in mask[:mask_max]]) / len(sequence)\n", "\n", " vals = []\n", " for i in range(20):\n", " test_data, signal = generate_test_signal(duration=2, nbits=nbits, signal_amplitude=amplitude, decimation=decimation, seed=random.randint(0, 0xffffffff), data=data)\n", " cor_an = correlate(signal, nbits=nbits, decimation=decimation, mask_max=2)\n", " vals.append(np.sqrt(np.mean(cor_an**2)))\n", " ps_off.append((np.mean(vals), np.std(vals)))\n", " print(vals)\n", "\n", "ps_off = np.array(ps_off)" ] }, { "cell_type": "code", "execution_count": 129, "metadata": {}, "outputs": [ { "data": { "application/vnd.jupyter.widget-view+json": { "model_id": "3dade90110af4097b123f3de508e843f", "version_major": 2, "version_minor": 0 }, "text/plain": [ "Canvas(toolbar=Toolbar(toolitems=[('Home', 'Reset original view', 'home', 'home'), ('Back', 'Back to previous …" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "" ] }, "execution_count": 129, "metadata": {}, "output_type": "execute_result" } ], "source": [ "fig, ax = plt.subplots(figsize=(12, 5))\n", "ax.errorbar(x=xs, y=ps_off[:,0], yerr=ps_off[:,1])" ] } ], "metadata": { "kernelspec": { "display_name": "ma-thesis-env", "language": "python", "name": "ma-thesis-env" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.9.2" } }, "nbformat": 4, "nbformat_minor": 4 }