From 25891ce723f8bf2f5dbf9a8451e1a296cfe3c2e3 Mon Sep 17 00:00:00 2001 From: Kevin Fronczak Date: Tue, 25 Jun 2024 15:30:02 -0400 Subject: [PATCH 1/2] Update pyproject.toml --- pyproject.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pyproject.toml b/pyproject.toml index dd606f2..8851284 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -4,7 +4,7 @@ build-backend = "setuptools.build_meta" [project] name = "python-adc-eval" -version = "0.4.0b0" +version = "0.3.0" license = {text = "MIT"} description = "ADC Evaluation Library" readme = "README.rst" From 64ef0087d98e5d73c0abe2c038db4b2adc3f9501 Mon Sep 17 00:00:00 2001 From: Kevin Fronczak Date: Tue, 15 Apr 2025 14:01:19 -0400 Subject: [PATCH 2/2] Refactoring to allow addition of ADC models and simulation --- .github/workflows/tests.yml | 2 +- .gitignore | 1 + README.rst | 3 +- adc_eval/__init__.py | 4 - adc_eval/adcs/__init__.py | 1 + adc_eval/adcs/basic.py | 154 +++++++++++++++++++++ adc_eval/converters.py | 29 ---- adc_eval/eval/simulate.py | 49 +++++++ adc_eval/{ => eval}/spectrum.py | 118 +++++++++++----- adc_eval/filter.py | 229 +++++++++++++++++++++++++++++++ adc_eval/signals.py | 9 +- examples/basic_adc_simulation.py | 73 ++++++++++ requirements.txt | 2 + tests/test_calc_psd.py | 22 +-- tests/test_spectrum.py | 6 +- tests/test_spectrum_plotting.py | 2 +- tox.ini | 2 +- 17 files changed, 620 insertions(+), 86 deletions(-) create mode 100644 adc_eval/adcs/__init__.py create mode 100644 adc_eval/adcs/basic.py delete mode 100644 adc_eval/converters.py create mode 100644 adc_eval/eval/simulate.py rename adc_eval/{ => eval}/spectrum.py (72%) create mode 100644 adc_eval/filter.py create mode 100644 examples/basic_adc_simulation.py diff --git a/.github/workflows/tests.yml b/.github/workflows/tests.yml index 8b2ae6e..7b15d6d 100644 --- a/.github/workflows/tests.yml +++ b/.github/workflows/tests.yml @@ -14,7 +14,7 @@ jobs: matrix: platform: - ubuntu-latest - python-version: ['3.9', '3.10', '3.11', '3.12'] + python-version: ['3.10', '3.11', '3.12'] steps: - uses: actions/checkout@v3 diff --git a/.gitignore b/.gitignore index 93b30a0..c2ee334 100644 --- a/.gitignore +++ b/.gitignore @@ -5,4 +5,5 @@ python_adc_eval.egg-info dist .ruff_cache .coverage +.secret build diff --git a/README.rst b/README.rst index fa08338..d936b9c 100644 --- a/README.rst +++ b/README.rst @@ -42,7 +42,8 @@ Given an array of values representing the output of an ADC, the spectrum can be leak=, window=, no_plot=, - yaxis=<"power"/"fullscale"> + yaxis=<"power"/"fullscale"/"magnitude">, + single_sided= ) diff --git a/adc_eval/__init__.py b/adc_eval/__init__.py index 00d2b6a..b4593a3 100644 --- a/adc_eval/__init__.py +++ b/adc_eval/__init__.py @@ -1,5 +1 @@ """Initialization file for module.""" - -from . import spectrum -from . import converters -from . import signals diff --git a/adc_eval/adcs/__init__.py b/adc_eval/adcs/__init__.py new file mode 100644 index 0000000..b4593a3 --- /dev/null +++ b/adc_eval/adcs/__init__.py @@ -0,0 +1 @@ +"""Initialization file for module.""" diff --git a/adc_eval/adcs/basic.py b/adc_eval/adcs/basic.py new file mode 100644 index 0000000..01210fe --- /dev/null +++ b/adc_eval/adcs/basic.py @@ -0,0 +1,154 @@ +"""Basic ADC models.""" + +import numpy as np + + +def dac(samples, nbits=8, vref=1): + """Digital to analog converter.""" + quants = 2**nbits + dv = vref / quants + return samples * dv + + +class ADC: + """ + Generic ADC Class. + + ... + + Parameters + ---------- + nbits : int, default=8 + Number of bits for the ADC. + fs : int or float, default=1 + Sample rate for the ADC in Hz. + vref : int or float, default=1 + Reference level of the ADC in Volts ([0, +vref] conversion range). + seed : int, default=1 + Seed for random variable generation. + **kwargs + Extra arguments. + + Attributes + ------- + vin : float + Sets or returns the current input voltage level. Assumed +/-vref/2 input + vlsb : float + LSB voltage of the converter. vref/2^nbits + noise : float, default=0 + Sets or returns the stdev of the noise generated by the converter. + mismatch : float, default=0 + Sets or returns the stdev of the mismatch of the converter. + offset : tuple of float, default=(0, 0) + Sets the (mean, stdev) of the offset of the converter. + gain_error : tuple of float, default=(0, 0) + Sets the (mean, stdev) of the gain error of the converter. + distortion : list of float, default=[1] + Sets the harmonic distortion values with index=0 corresponding to HD1. + Example: For unity gain and only -30dB of HD3, input is [1, 0, 0.032] + dout : int + Digital output code for current vin value. + + Methods + ------- + run_step + + """ + + def __init__(self, nbits=8, fs=1, vref=1, seed=1, **kwargs): + """Initialization function for Generic ADC.""" + np.random.seed(seed) + self.nbits = nbits + self.fs = fs + self.vref = vref + self.seed = seed + self.err = {"noise": 0, "gain": 0, "dist": [1], "offset": 0, "mismatch": 0} + self.dbits = np.zeros(nbits) + self.dval = 0 + + @property + def vin(self): + """Return input value.""" + return self._vin + + @vin.setter + def vin(self, x): + """Set the input value.""" + x += self.vref / 2 + x = max(0, min(x, self.vref)) + self._vin = x + + @property + def vlsb(self): + """Return the LSB voltage.""" + return self.vref / 2**self.nbits + + @property + def noise(self): + """Return noise status.""" + return self.err["noise"] + + @noise.setter + def noise(self, stdev): + """Set noise stdev in Vrms.""" + self.err["noise"] = stdev + + @property + def mismatch(self): + """Return noise stdev.""" + print("WARNING: 'mismatch' feature not implemented for this class.") + return False + + @mismatch.setter + def mismatch(self, stdev): + """Set mismatch stdev.""" + print("WARNING: 'mismatch' feature not implemented for this class.") + pass + + @property + def offset(self): + """Return offset value.""" + return self.err["offset"] + + @offset.setter + def offset(self, values): + """Set offset mean and stdev.""" + self.err["offset"] = values[0] + values[1] * np.random.randn(1) + + @property + def gain_error(self): + """Return gain error status.""" + return self.err["gain"] + + @gain_error.setter + def gain_error(self, values): + """Set gain error mean and stdev.""" + self.err["gain"] = values[0] + values[1] * np.random.randn(1) + + @property + def distortion(self): + """Return distortion gains (1st-order indexed).""" + return self.err["dist"] + + @distortion.setter + def distortion(self, gains): + """Set distortion gains (1st-order indexed).""" + self.err["dist"] = gains + + @property + def dout(self): + """Return digital output code.""" + return int(self.dval) + + def run_step(self): + """Run a single ADC step.""" + vinx = self.vin + dval = int( + min(max(int((2**self.nbits) * vinx / self.vref), 0), 2**self.nbits - 1) + ) + bits = [int(x) for x in bin(dval)[2:]] + + while len(bits) < self.nbits: + bits.insert(0, 0) + self.dbits = bits + self.dval = dval diff --git a/adc_eval/converters.py b/adc_eval/converters.py deleted file mode 100644 index bea7480..0000000 --- a/adc_eval/converters.py +++ /dev/null @@ -1,29 +0,0 @@ -"""Analog <-> Digital converters behavioral models.""" - -import numpy as np - - -def analog2digital( - sig_f, sample_freq=1e6, sample_n=1024, sample_bits=8, vref=3.3, noisy_lsb=1 -): - """Analog to digital converter.""" - sample_quants = 2**sample_bits - sample_prd = 1 / sample_freq - t = np.arange(0, sample_n * sample_prd, sample_prd) - dv = vref / sample_quants - samples = np.rint(sig_f(t) / dv).astype(int) - # apply noise - if noisy_lsb: - noise = np.random.normal(0, 2 ** (noisy_lsb - 1), size=sample_n).astype(int) - samples += noise - # to be sure that samples fit the range: - samples[np.argwhere(samples >= sample_quants)] = sample_quants - 1 - samples[np.argwhere(samples < 0)] = 0 - return samples - - -def digital2analog(samples, sample_bits=8, vref=3.3): - """Digital to analog converter.""" - quants = 2**sample_bits - dv = vref / quants - return samples * dv diff --git a/adc_eval/eval/simulate.py b/adc_eval/eval/simulate.py new file mode 100644 index 0000000..7a4765b --- /dev/null +++ b/adc_eval/eval/simulate.py @@ -0,0 +1,49 @@ +"""Generic simulator class for adc evaluation.""" + +import numpy as np +from tqdm import tqdm + + +class Simulator: + """Class for handling simulation functions.""" + + def __init__(self, adc_obj, xarray): + """Initialize the simulator class.""" + self.dval = [] + self.adc = adc_obj + self.vin = self.calc_error(xarray) + + @property + def out(self): + """Return output value array.""" + return np.array(self.dval) + + def calc_error(self, vin): + """Using the adc obj, calculates global signal error.""" + vinx = vin + + # First calculate gain error + vinx *= (1 + self.adc.err["gain"]) * self.adc.err["dist"][0] + + # Next calculate the harmonic distortion + for index, gain in enumerate(self.adc.err["dist"]): + if index > 0: + vinx += gain * vin ** (index + 1) + + # Now add the offset + vinx += self.adc.err["offset"] + + # Now add random noise + vinx += self.adc.err["noise"] * np.random.randn(vin.size) + + return vinx + + def run(self): + with tqdm( + range(len(self.vin)), "RUNNING", unit=" samples", position=0, leave=True + ) as pbar: + for xval in self.vin: + self.adc.vin = xval + self.adc.run_step() + self.dval.append(self.adc.dout) + pbar.update() diff --git a/adc_eval/spectrum.py b/adc_eval/eval/spectrum.py similarity index 72% rename from adc_eval/spectrum.py rename to adc_eval/eval/spectrum.py index 3c4b35e..18ecd1a 100644 --- a/adc_eval/spectrum.py +++ b/adc_eval/eval/spectrum.py @@ -25,6 +25,7 @@ def enob(sndr, places=1): def sndr_sfdr(spectrum, freq, fs, nfft, leak, full_scale=0): """Get SNDR and SFDR.""" + # Zero the DC bin spectrum[0] = 0 bin_sig = np.argmax(spectrum) @@ -118,27 +119,35 @@ def find_harmonics(spectrum, freq, nfft, bin_sig, psig, harms=5, leak=20): return harm_stats -def calc_psd(data, fs, nfft=2**12, single_sided=False): +def calc_psd(data, fs, nfft=2**12): """Calculate the PSD using the Bartlett method.""" - nwindows = int(np.floor(len(data) / nfft)) + nwindows = max(1, int(np.floor(len(data) / nfft))) nfft = int(nfft) xs = data[0 : int(nwindows * nfft)] xt = xs.reshape(nwindows, nfft).T XF = abs(np.fft.fft(xt, nfft, axis=0) / nfft) ** 2 psd = np.mean(XF, axis=1) / (fs / nfft) # average the ffts and divide by bin width - freq = fs * np.linspace(0, 1, nfft) - if single_sided: - # First we double all the bins, then we halve the DC bin - psd = 2 * psd[0 : int(nfft / 2)] - psd[0] /= 2 - freq = freq[0 : int(nfft / 2)] - return (freq, psd) + psd += np.finfo(float).eps # Prevents zeros in the PSD + freq = np.fft.fftshift(np.fft.fftfreq(nfft, d=1 / fs)) + + # For single sided we double all the bins, then we halve the DC bin + psd_ss = 2 * psd[0 : int(nfft / 2)] + psd_ss[0] /= 2 + freq_ss = freq[int(nfft / 2) :] + + # Need to rotate DS PSD so 0Hz is in middle of graph + freq_ds = freq + psd_ds = np.concatenate([psd[int(nfft / 2) :], psd[0 : int(nfft / 2)]]) + return [freq_ss, psd_ss, freq_ds, psd_ds] -def get_spectrum(data, fs=1, nfft=2**12): + +def get_spectrum(data, fs=1, nfft=2**12, single_sided=True): """Get the power spectrum for an input signal.""" - (freq, psd) = calc_psd(np.array(data), fs=fs, nfft=nfft, single_sided=True) - return (freq, psd * fs / nfft) + (freq_ss, psd_ss, freq_ds, psd_ds) = calc_psd(np.array(data), fs=fs, nfft=nfft) + if single_sided: + return (freq_ss, psd_ss * fs / nfft) + return (freq_ds, psd_ds * fs / nfft) def plot_spectrum( @@ -151,6 +160,7 @@ def plot_spectrum( window="rectangular", no_plot=False, yaxis="power", + single_sided=True, ): """Plot Power Spectrum for input signal.""" wsize = data.size @@ -168,21 +178,38 @@ def plot_spectrum( "hanning": 1.633, }[window] - (freq, pwr) = get_spectrum(data * windows[window] * wscale, fs=fs, nfft=nfft) + (freq, pwr) = get_spectrum( + data * windows[window] * wscale, fs=fs, nfft=nfft, single_sided=single_sided + ) full_scale = dBW(dr**2 / 8) - scalar = 0 - yunits = "dB" - if yaxis.lower() == "fullscale": - scalar = full_scale - yunits = "dBFS" + yaxis_lut = { + "power": [0, "dB"], + "fullscale": [dBW(dr**2 / 8), "dBFS"], + "normalize": [max(dBW(pwr)), "dB Normalized"], + "magnitude": [0, "W"], + } + + lut_key = yaxis.lower() + scalar = yaxis_lut[lut_key][0] + yunits = yaxis_lut[lut_key][1] - pwr_dB = 10 * np.log10(pwr) - scalar + psd_out = 10 * np.log10(pwr) - scalar + if lut_key in ["magnitude"]: + psd_out = pwr + + f_ss = freq + psd_ss = pwr + if not single_sided: + # Get single-sided spectrum for SNDR and Harmonic stats + (f_ss, psd_ss) = get_spectrum( + data * windows[window] * wscale, fs=fs, nfft=nfft, single_sided=True + ) - sndr_stats = sndr_sfdr(pwr, freq, fs, nfft, leak=leak, full_scale=full_scale) + sndr_stats = sndr_sfdr(psd_ss, f_ss, fs, nfft, leak=leak, full_scale=full_scale) harm_stats = find_harmonics( - pwr, - freq, + psd_ss, + f_ss, nfft, sndr_stats["sig"]["bin"], sndr_stats["sig"]["power"], @@ -192,16 +219,18 @@ def plot_spectrum( stats = {**sndr_stats, **harm_stats} + xmin = 0 if single_sided else -fs / 2e6 + if not no_plot: plt_str = get_plot_string(stats, full_scale, fs, nfft, window) - fig, ax = plt.subplots(figsize=(15, 8)) - ax.plot(freq / 1e6, pwr_dB) + ax.plot(freq / 1e6, psd_out) ax.set_ylabel(f"Power Spectrum ({yunits})", fontsize=18) ax.set_xlabel("Frequency (MHz)", fontsize=16) ax.set_title("Output Power Spectrum", fontsize=16) - ax.set_xlim([0, fs / 2e6]) - ax.set_ylim([1.1 * min(pwr_dB), 0]) + ax.set_xlim([xmin, fs / 2e6]) + ax.set_ylim([1.1 * min(psd_out), 1]) + ax.annotate( plt_str, xy=(1, 1), @@ -217,11 +246,15 @@ def plot_spectrum( noise_dB = stats["noise"]["dBHz"] + full_scale # Add points for harmonics and largest spur + if not single_sided: + scalar += 3 for hindex in range(2, harmonics + 1): if stats["harm"][hindex]["dB"] > (noise_dB + 3): + fharm = stats["harm"][hindex]["freq"] + aharm = stats["harm"][hindex]["dB"] - scalar ax.plot( - stats["harm"][hindex]["freq"], - stats["harm"][hindex]["dB"] - scalar, + fharm, + aharm, marker="s", mec="r", ms=8, @@ -229,16 +262,33 @@ def plot_spectrum( mew=3, ) ax.text( - stats["harm"][hindex]["freq"], - stats["harm"][hindex]["dB"] - scalar + 3, + fharm, + aharm + 3, f"HD{hindex}", ha="center", weight="bold", ) + if not single_sided: + ax.plot( + -fharm, + aharm, + marker="s", + mec="r", + ms=8, + fillstyle="none", + mew=3, + ) + ax.text( + -fharm, + aharm + 3, + f"HD{hindex}", + ha="center", + weight="bold", + ) ax.tick_params(axis="both", which="major", labelsize=14) ax.grid() - return (pwr, stats) + return (freq, psd_out, stats) def get_plot_string(stats, full_scale, fs, nfft, window): @@ -284,9 +334,10 @@ def analyze( window="rectangular", no_plot=False, yaxis="fullscale", + single_sided=True, ): """Perform spectral analysis on input waveform.""" - (spectrum, stats) = plot_spectrum( + (freq, spectrum, stats) = plot_spectrum( data, fs=fs, nfft=nfft, @@ -296,6 +347,7 @@ def analyze( window=window, no_plot=no_plot, yaxis=yaxis, + single_sided=single_sided, ) - return (spectrum, stats) + return (freq, spectrum, stats) diff --git a/adc_eval/filter.py b/adc_eval/filter.py new file mode 100644 index 0000000..53e8c7b --- /dev/null +++ b/adc_eval/filter.py @@ -0,0 +1,229 @@ +"""Implements filters and decimation.""" + +import numpy as np +from scipy.signal import remez, freqz +import matplotlib.pyplot as plt +from adc_eval.eval import spectrum + + +class CICDecimate: + """ + Generic CIC Decimator Object. + + ... + + Parameters + ---------- + dec : int, default=2 + Output decimation factor. + order : int, default=1 + Filter order. + fs : int or float, default=1 + Sample rate for the filter in Hz. + + + Attributes + ---------- + gain : Gain normalization factor of CIC filter + out : Filtered and decimated output data + + + Methods + ------- + run + response + + """ + + def __init__(self, dec=2, order=1, fs=1): + """Initialize the CIC filter.""" + self.dec = dec + self.order = order + self.fs = fs + self.gain = self.dec**self.order + self._xout = None + self._xfilt = None + + @property + def out(self): + """Filtered and decimated output data.""" + return np.array(self._xout) + + def filt(self, xarray): + """CIC filtering routine.""" + yint = xarray + + # Integrate first + for _ in range(self.order): + yint = np.cumsum(yint) + + # Then comb, adding delays based on decimation factor + xcomb = yint + xcomb = np.insert(xcomb, 0, [0 for x in range(self.dec)]) + ycomb = xcomb + for _ in range(self.order): + ycomb = ycomb[self.dec :] - ycomb[0 : -self.dec] + xcomb = ycomb + xcomb = np.insert(xcomb, 0, [0 for x in range(self.dec)]) + ycomb = xcomb + + self._xfilt = ycomb / self.gain + + def decimate(self): + """decimation routine.""" + self._xout = self._xfilt[:: self.dec] + + def run(self, xarray): + """Runs filtering and decimation on input list.""" + self.filt(xarray) + self.decimate() + + def response(self, fft, no_plot=False): + """Plots the frequency response of the pre-decimated filter.""" + xin = np.zeros(fft) + xin[0] = 1 + self.filt(xin) + (freq, psd, stats) = spectrum.analyze( + self._xfilt * fft / np.sqrt(2), + fft, + fs=self.fs, + dr=1, + harmonics=0, + leak=1, + window="rectangular", + no_plot=no_plot, + yaxis="power", + single_sided=True, + ) + if not no_plot: + ax = plt.gca() + n1 = 0 + n2 = int(fft / (2 * self.dec)) + x = freq[n1:n2] / 1e6 + y1 = psd[n1:n2] - max(psd) + y2 = -2000 * np.ones(np.size(x)) + ax.plot(x, y2, alpha=0) + ax.plot(x, y1, alpha=0) + ax.fill_between(x, y1, y2, color="green", alpha=0.1) + ax.set_xticks(np.linspace(0, self.fs / 2e6, 9)) + return (freq, psd) + + +class FIRLowPass: + """ + Generic FIR Low Pass Filter. + + ... + + Parameters + ---------- + dec : int, optional + Output decimation rate. The default is 1. + fs : int or float, optional + Sample rate for the filter in Hz. The default is 1Hz. + bit_depth : int, optional + Bit depth to store coefficients. The default is 16b. + coeffs : list, optional + List of coefficients if pre-determined. + + + Attributes + ---------- + out : Filtered and decimated output data. + ntaps : Number of filter taps. + + + Methods + ------- + run + response + + """ + + def __init__(self, dec=1, fs=1, bit_depth=16, coeffs=None): + """Initialize the FIR LowPass Class.""" + self.coeffs = coeffs + self.dec = dec + self.fs = fs + self.bit_depth = bit_depth + self.ntaps = np.size(coeffs) if coeffs is None else 0 + self.yfilt = None + self._out = None + + @property + def out(self): + """Filtered and decimated output datat.""" + return np.array(self._out) + + def generate_taps(self, fbw, pbripple=1, stopatt=-60, deltaf=None): + """ + Generates FIR taps from key inputs. + + Parameters + ---------- + fbw : float + Bandwidth of the filter in Hz. + pbripple : float, optional + Acceptable passband ripple in percentage. The default is 1%. + stopatt : float, optional + Desired stop-band attenuation. The default is -60dB. + deltaf : float, optional + Desired transition band of the filter. The default is FS/100 if set to None. + + Returns + ------- + (ntaps, coeffs) : tuple + Minimum number of taps required, List of FIR tap coefficients. + + """ + if deltaf is None: + deltaf = self.fs / 100 + x1 = pbripple / 100 + x2 = 10 ** (stopatt / 20) + x3 = np.log10(1 / (10 * x1 * x2)) + _ntaps = 2 / 3 * x3 * self.fs / deltaf + + if self.ntaps > 0 and _ntaps > self.ntaps: + print( + f"WARNING: Required NTAPs calculated as {int(_ntaps)} but only {self.ntaps} were provided." + ) + elif self.ntaps == 0: + self.ntaps = int(_ntaps) + + _coeffs = remez( + self.ntaps, [0, fbw, fbw + deltaf, self.fs / 2], [1, 0], fs=self.fs + ) + self.coeffs = np.int32(_coeffs * 2**self.bit_depth).tolist() + return (self.ntaps, self.coeffs) + + def filt(self, xarray): + """Performs FIR filtering on input xarray.""" + _coeffs = np.array(self.coeffs) / 2**self.bit_depth + self.yfilt = np.convolve(_coeffs, xarray, mode="same") + + def decimate(self, xarray=None): + """Performs decimation on the input data.""" + if xarray is None: + xarray = self.yfilt + self._out = xarray[:: self.dec] + + def run(self, xarray): + """Runs FIR filtering and decimation on input xarray data.""" + self.filt(xarray) + self.decimate() + + def response(self, fft, no_plot=False): + """Plots the frequency response of the pre-decimated filter.""" + freq, mag = freqz(self.coeffs, [1], worN=fft, fs=self.fs) + yfft = spectrum.dBW(np.abs(mag)) + if not no_plot: + fig, ax = plt.subplots(figsize=(15, 8)) + ax.plot(freq / 1e6, yfft) + ax.grid(True) + ax.set_ylabel("Filter Magnitude Response (dB)", fontsize=18) + ax.set_xlabel("Frequency (MHz)", fontsize=16) + ax.set_title("FIR Low Pass Response", fontsize=16) + ax.set_xlim([0, self.fs / 2e6]) + ax.set_ylim([1.1 * min(yfft), 1]) + plt.show() + return (freq, yfft) diff --git a/adc_eval/signals.py b/adc_eval/signals.py index 0aeb71c..002e019 100644 --- a/adc_eval/signals.py +++ b/adc_eval/signals.py @@ -3,9 +3,14 @@ import numpy as np -def sin(t, peak=1.5, offset=1.65, freq=1e3, ph0=0): +def time(nsamp, fs=1): + """Create time array based on signal length and sample rate.""" + return 1 / fs * np.linspace(0, nsamp - 1, nsamp) + + +def sin(t, amp=0.5, offset=0.5, freq=1e3, ph0=0): """Generate a sine wave.""" - return offset + peak * np.sin(ph0 + 2 * np.pi * freq * t) + return offset + amp * np.sin(ph0 + 2 * np.pi * freq * t) def noise(t, mean=0, std=0.1): diff --git a/examples/basic_adc_simulation.py b/examples/basic_adc_simulation.py new file mode 100644 index 0000000..914615d --- /dev/null +++ b/examples/basic_adc_simulation.py @@ -0,0 +1,73 @@ +"""Runs a basic ADC simulation and plots the spectrum.""" +import numpy as np +import matplotlib.pyplot as plt +from scipy import signal +from adc_eval import signals +from adc_eval.adcs import basic +from adc_eval.eval import spectrum +from adc_eval.eval.simulate import Simulator + + +""" +Simulation Settings +""" +SEED = 42 +NBITS = 10 +FS = 200e6 +NLEN = 2**16 # Larger than NFFT to enable Bartlett method for PSD +NFFT = 2**12 +vref = 1 +fin_bin = NFFT / 4 - 31 +fin = fin_bin * FS/NFFT +vin_amp = 0.707 * vref / 2 + + +""" +VIN Generation +""" +t = signals.time(NLEN, FS) +vin = signals.sin(t, amp=vin_amp, offset=0, freq=fin) +vin += signals.sin(t, amp=vin_amp*0.2, offset=0, freq=(NFFT/2-15)*FS/NFFT) # Adds tone to show intermodulation + +""" +ADC Architecture Creation +""" +adc_dut = basic.ADC(nbits=NBITS, vref=vref, fs=FS, seed=SEED) + + +""" +Global ADC Error Settings +""" +adc_dut.noise = 0 # No internal noise generation +adc_dut.offset = (0, 0) # 10mV mean offset with no stdev +adc_dut.gain_err = (0, 0) # No internal gain error +adc_dut.distortion = [1, 0, 0.3] # HD3 only + + +""" +Run Simulation +""" +s = Simulator(adc_dut, vin) +s.run() + + +""" +Output Plotting +""" +(freq, ps, stats) = spectrum.analyze( + s.out, + NFFT, + fs=FS, + dr=2**NBITS, + harmonics=7, + leak=1, + window="rectangular", + no_plot=False, + yaxis="fullscale", + single_sided=True, +) +ax = plt.gca() +ax.set_title("ADC Spectrum") +ax.set_ylim([-100, 0]) +ax.set_yticks(np.linspace(-100, 0, 11)) +ax.set_xticks(np.linspace(0, FS/2e6, 9)) diff --git a/requirements.txt b/requirements.txt index 9a9169b..adfbc3c 100644 --- a/requirements.txt +++ b/requirements.txt @@ -1 +1,3 @@ matplotlib==3.9.0 +tqdm==4.67.1 +scipy==1.15.2 diff --git a/tests/test_calc_psd.py b/tests/test_calc_psd.py index 8aae2e8..1c5aa4f 100644 --- a/tests/test_calc_psd.py +++ b/tests/test_calc_psd.py @@ -3,7 +3,7 @@ import unittest import numpy as np from unittest import mock -from adc_eval import spectrum +from adc_eval.eval import spectrum class TestCalcPSD(unittest.TestCase): @@ -21,7 +21,7 @@ def test_calc_psd_randomized_dual(self): """Test calc_psd with random data.""" for i in range(0, 10): data = np.random.randn(self.nlen) - (freq, psd) = spectrum.calc_psd(data, 1, nfft=self.nfft, single_sided=False) + (_, _, freq, psd) = spectrum.calc_psd(data, 1, nfft=self.nfft) mean_val = np.mean(psd) self.assertTrue(self.bounds[0] <= mean_val <= self.bounds[1], msg=mean_val) @@ -29,7 +29,7 @@ def test_calc_psd_randomized_single(self): """Test calc_psd with random data and single-sided.""" for i in range(0, 10): data = np.random.randn(self.nlen) - (freq, psd) = spectrum.calc_psd(data, 1, nfft=self.nfft, single_sided=True) + (freq, psd, _, _) = spectrum.calc_psd(data, 1, nfft=self.nfft) mean_val = np.mean(psd) self.assertTrue( 2 * self.bounds[0] <= mean_val <= 2 * self.bounds[1], msg=mean_val @@ -38,7 +38,7 @@ def test_calc_psd_randomized_single(self): def test_calc_psd_zeros_dual(self): """Test calc_psd with zeros.""" data = np.zeros(self.nlen) - (freq, psd) = spectrum.calc_psd(data, 1, nfft=self.nfft, single_sided=False) + (_, _, freq, psd) = spectrum.calc_psd(data, 1, nfft=self.nfft) mean_val = np.mean(psd) self.assertTrue( self.bounds[0] - 1 <= mean_val <= self.bounds[1] - 1, msg=mean_val @@ -47,7 +47,7 @@ def test_calc_psd_zeros_dual(self): def test_calc_psd_zeros_single(self): """Test calc_psd with zeros and single-sided..""" data = np.zeros(self.nlen) - (freq, psd) = spectrum.calc_psd(data, 1, nfft=self.nfft, single_sided=True) + (freq, psd, _, _) = spectrum.calc_psd(data, 1, nfft=self.nfft) mean_val = np.mean(psd) self.assertTrue( self.bounds[0] - 1 <= mean_val <= self.bounds[1] - 1, msg=mean_val @@ -56,14 +56,14 @@ def test_calc_psd_zeros_single(self): def test_calc_psd_ones_dual(self): """Test calc_psd with ones.""" data = np.ones(self.nlen) - (freq, psd) = spectrum.calc_psd(data, 1, nfft=self.nfft, single_sided=False) + (_, _, freq, psd) = spectrum.calc_psd(data, 1, nfft=self.nfft) mean_val = np.mean(psd) self.assertTrue(self.bounds[0] <= mean_val <= self.bounds[1], msg=mean_val) def test_calc_psd_ones_single(self): """Test calc_psd with ones and single-sided.""" data = np.ones(self.nlen) - (freq, psd) = spectrum.calc_psd(data, 1, nfft=self.nfft, single_sided=True) + (freq, psd, _, _) = spectrum.calc_psd(data, 1, nfft=self.nfft) mean_val = np.mean(psd) self.assertTrue( 2 * self.bounds[0] <= mean_val <= 2 * self.bounds[1], msg=mean_val @@ -79,13 +79,13 @@ def test_calc_psd_two_sine_dual(self): a2 = 0.11 t = 1 / fs * np.linspace(0, self.nlen - 1, self.nlen) data = a1 * np.sin(2 * np.pi * f1 * t) + a2 * np.sin(2 * np.pi * f2 * t) - (freq, psd) = spectrum.calc_psd(data, fs, nfft=self.nfft, single_sided=False) + (_, _, freq, psd) = spectrum.calc_psd(data, fs, nfft=self.nfft) exp_peaks = [ round(a1**2 / 4 * self.nfft, 3), round(a2**2 / 4 * self.nfft, 3), ] - exp_f1 = [round(f1, 2), round(fs - f1, 2)] - exp_f2 = [round(f2, 2), round(fs - f2, 2)] + exp_f1 = [round(-f1, 2), round(f1, 2)] + exp_f2 = [round(-f2, 2), round(f2, 2)] peak1 = max(psd) ipeaks = np.where(psd >= peak1 * self.bounds[0])[0] @@ -114,7 +114,7 @@ def test_calc_psd_two_sine_single(self): a2 = 0.11 t = 1 / fs * np.linspace(0, self.nlen - 1, self.nlen) data = a1 * np.sin(2 * np.pi * f1 * t) + a2 * np.sin(2 * np.pi * f2 * t) - (freq, psd) = spectrum.calc_psd(data, fs, nfft=self.nfft, single_sided=True) + (freq, psd, _, _) = spectrum.calc_psd(data, fs, nfft=self.nfft) exp_peaks = [ round(a1**2 / 2 * self.nfft, 3), round(a2**2 / 2 * self.nfft, 3), diff --git a/tests/test_spectrum.py b/tests/test_spectrum.py index 625b314..3a46f34 100644 --- a/tests/test_spectrum.py +++ b/tests/test_spectrum.py @@ -3,7 +3,7 @@ import unittest import numpy as np from unittest import mock -from adc_eval import spectrum +from adc_eval.eval import spectrum class TestSpectrum(unittest.TestCase): @@ -47,7 +47,7 @@ def test_enob(self): for i in range(0, len(exp_val)): self.assertEqual(spectrum.enob(test_val, places=i), exp_val[i]) - @mock.patch("adc_eval.spectrum.calc_psd") + @mock.patch("adc_eval.eval.spectrum.calc_psd") def test_get_spectrum(self, mock_calc_psd): """Test that the get_spectrum method returns power spectrum.""" fs = 4 @@ -55,7 +55,7 @@ def test_get_spectrum(self, mock_calc_psd): data = np.array([1]) exp_spectrum = np.array([fs / nfft]) - mock_calc_psd.return_value = (None, data) + mock_calc_psd.return_value = (None, data, None, data) self.assertEqual( spectrum.get_spectrum(None, fs=fs, nfft=nfft), (None, exp_spectrum) diff --git a/tests/test_spectrum_plotting.py b/tests/test_spectrum_plotting.py index 03a0d9f..d3e5d2d 100644 --- a/tests/test_spectrum_plotting.py +++ b/tests/test_spectrum_plotting.py @@ -3,7 +3,7 @@ import unittest import numpy as np from unittest import mock -from adc_eval import spectrum +from adc_eval.eval import spectrum class TestSpectrumPlotting(unittest.TestCase): diff --git a/tox.ini b/tox.ini index fc9bbe1..f77285b 100644 --- a/tox.ini +++ b/tox.ini @@ -1,5 +1,5 @@ [tox] -envlist = lint,build,py39,py310,py311,py312 +envlist = lint,build,py310,py311,py312 skip_missing_interpreters = True skipsdist = True