デルタシグマ変調(デルタシグマADC)のプログラム



https://www.google.com/search?q=%E3%83%87%E3%83%AB%E3%82%BF%E3%82%B7%E3%82%B0%E3%83%9E%E5%A4%89%E8%AA%BF&oq=%E3%81%A7%E3%82%8B%E3%81%9F%E3%81%97%E3%81%90&gs_lcrp=EgZjaHJvbWUqCQgBEAAYBBiABDIGCAAQRRg5MgkIARAAGAQYgAQyCwgCEEUYBBg7GIAEMgkIAxAAGAQYgAQyCwgEEAAYBBgKGIAEMgkIBRAAGAQYgAQyCQgGEAAYBBiABDIJCAcQABgEGIAEMgkICBAAGAQYgAQyCQgJEAAYBBiABNIBCDUxMjhqMGo3qAIIsAIB&sourceid=chrome&ie=UTF-8


```
import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import lfilter

# Parameters
fs = 100  # Sampling frequency
t = np.arange(0, 1, 1/fs)  # Time vector
Ain = np.sin(2 * np.pi * 5 * t)  # Input signal: 5 Hz sine wave
K = 10  # Gain of the integrator

# Simulate integrator output
integrator_output = np.cumsum(Ain) / fs * K

# Simulate quantizer output (adding quantization noise Q)
Q = np.random.normal(0, 0.1, len(t))  # Quantization noise
quantizer_output = integrator_output + Q

# Feedback signal
feedback_signal = Ain

# Output DOUT
DOUT = quantizer_output / K + feedback_signal

# Low-pass filter design
b, a = [1], [1, 0.5]  # Simple low-pass filter coefficients
DOUT_filtered = lfilter(b, a, DOUT)

# Plotting
plt.figure(figsize=(14, 8))

plt.subplot(4, 1, 1)
plt.plot(t, Ain, label='Input Signal $A_{in}$')
plt.title('Input Signal $A_{in}$')
plt.xlabel('Time [s]')
plt.ylabel('Amplitude')
plt.grid()

plt.subplot(4, 1, 2)
plt.plot(t, integrator_output, label='Integrator Output $K$')
plt.title('Integrator Output $K$')
plt.xlabel('Time [s]')
plt.ylabel('Amplitude')
plt.grid()

plt.subplot(4, 1, 3)
plt.plot(t, quantizer_output, label='Quantizer Output $Q$')
plt.title('Quantizer Output $Q$')
plt.xlabel('Time [s]')
plt.ylabel('Amplitude')
plt.grid()

plt.subplot(4, 1, 4)
plt.plot(t, DOUT, label='Output $D_{out}$')
plt.plot(t, DOUT_filtered, label='Filtered $D_{out}$', linestyle='dashed')
plt.title('Output $D_{out}$ and Filtered $D_{out}$')
plt.xlabel('Time [s]')
plt.ylabel('Amplitude')
plt.legend()
plt.grid()

plt.tight_layout()
plt.show()

```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3759109/64ee608d-c0f4-4cf6-d658-7fbffdf6e1b2.png)





```
import math
import numpy as np
import matplotlib.pyplot as plt

# Delta-Sigma modulator function
def delSig1st(x):
    v0 = 0
    y = 0
    z = []
    for i in range(len(x)):
        u = x[i] - y
        v = u + v0
        v0 = v
        if v >= 0:
            y = 1
        else:
            y = -1
        z.append(y)
    return z

# Signal generation
amp = 0.7
n = 256
omega = 2 * math.pi / 128 
x = np.arange(n)
y = amp * np.sin(omega * x)

# Modulation and spectrum calculation
z = delSig1st(y)

Fy = np.fft.fft(y)
spectY = np.abs(Fy)
spy = spectY[0: n//2]
spx = x[0: n//2]
Fz = np.fft.fft(z)
spectZ = np.abs(Fz)
spz = spectZ[0: n//2]

# Plotting and displaying
fig, axs = plt.subplots(2, 1, figsize=(10, 8))

# Plot time-domain signals
axs[0].plot(x, y, label='Input Signal')
axs[0].plot(x, z, label='Modulated Signal')
axs[0].set_title('Time Domain Signals')
axs[0].legend()
axs[0].grid()

# Plot frequency-domain signals
axs[1].plot(spx, spy, label='Input Signal Spectrum')
axs[1].plot(spx, spz, label='Modulated Signal Spectrum')
axs[1].set_title('Frequency Domain Signals')
axs[1].legend()
axs[1].grid()

plt.tight_layout()
plt.show()
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3759109/e1c3d9af-1382-9202-0eea-11636a9dbaf6.png)









```
def delSigNth(inp, z, N=1):
    dither = 0.0
    v = [0] * N  # Initialize N accumulators
    y = 0
    for i in range(len(inp)):
        u = inp[i] - y
        v[0] = u + v[0]
        for n in range(1, N):
            v[n] = v[n-1] + v[n]
        
        if v[-1] >= 0:
            y = 1
        else:
            y = -1
        
        z[i] = y + dither * random.random()

# Initialize parameters
nn = 1024 * (2 ** 8)  # Number of samples
amp = 0.7             # Amplitude
frq = 100.2           # Frequency
smfrq = 44100         # Sampling frequency
repn = 1              # Number of repetitions
N = 1                 # Order of Delta-Sigma Modulator

# Prepare data arrays
x = np.arange(nn)
ydat = np.zeros(nn)
z = np.zeros(nn)
sumdat = np.zeros(nn, dtype=complex)
logspect = np.zeros(nn)

# Generate Blackman window
w_black = np.blackman(nn)

# Generate signal and apply Nth-order delta-sigma modulation
for _ in range(repn):
    genSignal(amp, frq, smfrq, x, ydat)
    delSigNth(ydat, z, N)
    
    # Apply Blackman window
    z = z * w_black
    
    # Compute FFT and accumulate results
    sumdat += np.fft.fft(z) / repn

# Frequency axis
freq_axis = np.fft.fftfreq(nn, 1/smfrq)
half_nn = nn // 2  # Half the number of samples for positive frequencies

# Calculate magnitude spectrum in dB
spectZ = np.abs(sumdat) / nn
logspect = 20 * np.log10(spectZ + 1e-100)  # Avoid log(0) by adding a small number

# Plot the spectrum
plt.plot(freq_axis[:half_nn], logspect[:half_nn])
plt.xscale("log")
plt.ylim(-120, 0)
plt.xlabel("Frequency (Hz)")
plt.ylabel("Magnitude (dB)")
plt.title(f"Frequency Spectrum of {N}-order Delta-Sigma Modulated Signal")
plt.grid(True)
plt.show()
```
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3759109/12863e81-1a81-1537-7fac-e27c075d8e4d.png)




```
import numpy as np
import matplotlib.pyplot as plt

# Parameters
VREF = 2.0           # Reference voltage
num_samples = 1000   # Number of samples
frequency = 1.0      # Frequency of the sine wave
sampling_rate = 1000 # Sampling rate
t = np.arange(num_samples) / sampling_rate  # Time vector

# Generate a sine wave as input voltage
VIN = 0.5 * np.sin(2 * np.pi * frequency * t) + 1.0  # Sine wave centered at VREF/2

# Delta-Sigma modulator simulation
def delta_sigma_modulator(vin, vref):
    integrator = 0
    bitstream = []
    integrator_values = []
    for sample in vin:
        comparator = 1 if integrator > 0 else 0
        integrator += sample - comparator * vref
        bitstream.append(comparator)
        integrator_values.append(integrator)
    return bitstream, integrator_values

# Simulate the Delta-Sigma modulator
bitstream, integrator_values = delta_sigma_modulator(VIN, VREF)

# Plot the input sine wave
plt.figure(figsize=(12, 4))
plt.plot(t, VIN)
plt.title('Input Sine Wave')
plt.xlabel('Time [s]')
plt.ylabel('Amplitude [V]')
plt.grid(True)
plt.show()

# Plot the bitstream
plt.figure(figsize=(12, 4))
plt.step(t, bitstream, where='mid')
plt.ylim(-0.1, 1.1)
plt.title('Delta-Sigma Modulator Bitstream')
plt.xlabel('Sample Number')
plt.ylabel('Bit Value')
plt.grid(True)
plt.show()

# Plot the noise shaping (integrator values)
plt.figure(figsize=(12, 4))
plt.plot(t, integrator_values)
plt.title('Noise Shaping (Integrator Values)')
plt.xlabel('Time [s]')
plt.ylabel('Integrator Value')
plt.grid(True)
plt.show()

# Digital filter (simple moving average)
def digital_filter(bitstream, window_size):
    filtered_output = np.convolve(bitstream, np.ones(window_size) / window_size, mode='valid')
    return filtered_output

# Apply digital filter
window_size = 100
filtered_output = digital_filter(bitstream, window_size)

# Plot the filtered output
plt.figure(figsize=(12, 4))
plt.plot(filtered_output, marker='o')
plt.title('Filtered Output')
plt.xlabel('Sample Number')
plt.ylabel('Filtered Value')
plt.grid(True)
plt.show()

# Calculate the final ADC output
adc_output = np.mean(filtered_output) * VREF

```







![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3759109/770db806-69fc-354f-8d0e-f60bda8bb22f.png)







```
import numpy as np
import matplotlib.pyplot as plt
from scipy.fft import fft, fftfreq
from scipy.signal import butter, filtfilt

# パラメータ
n_bits = 1  # 量子化ビット数
n_points = 1000  # サンプル数
input_frequency = 1  # 入力信号の周波数（Hz）
output_range = 60000  # 出力範囲（nT）
sampling_rate = 100  # サンプリング周波数（Hz）
osr = 64  # オーバーサンプリング比

# サイン波入力信号の生成
time = np.arange(n_points) / sampling_rate
input_signal = np.sin(2 * np.pi * input_frequency * time) * output_range / 2

# Delta-Sigma変調
def delta_sigma_modulation(input_signal, output_range, n_points, n_bits):
    quantizer_output = []  # 量子化出力を格納するリスト
    integrator = 0  # 積分器の初期値
    integrator_output = []  # 積分器出力を格納するリスト
    for i in range(n_points):  # サンプル数分ループ
        integrator += input_signal[i]  # 積分器に入力信号を加算
        integrator_output.append(integrator)  # 積分器の出力を保存
        if integrator >= 0:  # 積分器の出力が0以上の場合
            quantizer_output.append(output_range)  # 出力範囲の正の値を量子化出力に追加
            integrator -= output_range  # 積分器から出力範囲の値を減算
        else:  # 積分器の出力が0未満の場合
            quantizer_output.append(-output_range)  # 出力範囲の負の値を量子化出力に追加
            integrator += output_range  # 積分器に出力範囲の値を加算
    return integrator_output, quantizer_output  # 積分器出力と量子化出力を返す

# Delta-Sigma変調信号の生成
integrator_output, quantizer_output = delta_sigma_modulation(input_signal, output_range, n_points, n_bits)

# ローパスフィルタの設計
def low_pass_filter(data, cutoff, fs, order=5):
    nyquist = 0.5 * fs  # ナイキスト周波数を計算
    normal_cutoff = cutoff / nyquist  # 正規化カットオフ周波数を計算
    b, a = butter(order, normal_cutoff, btype='low', analog=False)  # Butterworthフィルタの設計
    y = filtfilt(b, a, data)  # フィルタを適用
    return y  # フィルタ適用後のデータを返す

# 量子化出力にローパスフィルタを適用
filtered_output = low_pass_filter(quantizer_output, cutoff=5, fs=sampling_rate)

# プロット
plt.figure(figsize=(12, 18))

# 入力信号のプロット
plt.subplot(4, 1, 1)
plt.plot(time, input_signal, label='Input Sine Wave Signal')
plt.title('Input Signal')
plt.xlabel('Time [s]')
plt.ylabel('Amplitude')
plt.grid(True)
plt.legend()

# 積分器出力のプロット
plt.subplot(4, 1, 2)
plt.plot(time, integrator_output, label='Integrator Output')
plt.title('Integrator Output')
plt.xlabel('Time [s]')
plt.ylabel('Amplitude')
plt.grid(True)
plt.legend()

# 量子化出力のプロット
plt.subplot(4, 1, 3)
plt.plot(time, quantizer_output, label='Quantized Output (Delta-Sigma Modulated)')
plt.title('Quantized Output')
plt.xlabel('Time [s]')
plt.ylabel('Amplitude')
plt.grid(True)
plt.legend()

# フィルタ後出力のプロット
plt.subplot(4, 1, 4)
plt.plot(time, filtered_output, label='Filtered Output')
plt.title('Filtered Output')
plt.xlabel('Time [s]')
plt.ylabel('Amplitude')
plt.grid(True)
plt.legend()

plt.tight_layout()
plt.show()

```


![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3759109/865c6433-44e5-57a4-1ba9-1f1c98e9d6fc.png)


```

import numpy as np
import matplotlib.pyplot as plt

# Define the number of clock cycles
n_cycles = 100

# Generate a sample input signal u(n)
u = np.sin(np.linspace(0, 4 * np.pi, n_cycles))

# Initialize y(n) and v(n)
y = np.zeros(n_cycles)
v = np.zeros(n_cycles)

# Define the DSM equations
for n in range(1, n_cycles):
    y[n] = u[n] + y[n-1] - v[n-1]
    v[n] = np.round(y[n])  # Quantizer (simple rounding)

# Plot the signals
plt.figure(figsize=(12, 6))
plt.plot(u, label='u(n)', marker='o', linestyle='--')
plt.plot(y, label='y(n)', marker='x', linestyle='-.')
plt.plot(v, label='v(n)', marker='s', linestyle='-')
plt.xlabel('Clock Cycles (n)')
plt.ylabel('Amplitude')
plt.title('Delta-Sigma Modulation: Input, Intermediate, and Output Signals')
plt.legend()
plt.grid(True)
plt.show()
```


![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3759109/32d00bbc-b290-1a50-65bb-fccfa95adf6b.png)





```
import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import butter, lfilter

# パラメータ設定
fs = 1000  # サンプリング周波数 (Hz)
f = 5      # 入力サイン波の周波数 (Hz)
T = 1      # シミュレーション時間 (秒)
N = fs * T # サンプル数

# 時間軸の生成
t = np.linspace(0, T, N, endpoint=False)

# 入力信号（サイン波）の生成
input_signal = np.sin(2 * np.pi * f * t)

# デルタシグマ変調器の初期化
integrator = 0
quantizer_output = np.zeros(N)
feedback_signal = np.zeros(N)

# デルタシグマ変調器のシミュレーション
for i in range(N):
    integrator += input_signal[i] - feedback_signal[i-1]
    if integrator > 0:
        quantizer_output[i] = 1
    else:
        quantizer_output[i] = -1
    feedback_signal[i] = quantizer_output[i]

# ローパスフィルタの設計
def butter_lowpass(cutoff, fs, order=5):
    nyq = 0.5 * fs
    normal_cutoff = cutoff / nyq
    b, a = butter(order, normal_cutoff, btype='low', analog=False)
    return b, a

def lowpass_filter(data, cutoff, fs, order=5):
    b, a = butter_lowpass(cutoff, fs, order=order)
    y = lfilter(b, a, data)
    return y

# ローパスフィルタの適用
cutoff_frequency = 10  # カットオフ周波数 (Hz)
filtered_signal = lowpass_filter(quantizer_output, cutoff_frequency, fs)

# プロット
plt.figure(figsize=(12, 12))

# 入力信号のプロット
plt.subplot(4, 1, 1)
plt.plot(t, input_signal, label='Analog Input Signal (Sine Wave)')
plt.title('Analog Input Signal (Sine Wave)')
plt.xlabel('Time [s]')
plt.ylabel('Amplitude')
plt.grid(True)

# 積分器出力のプロット
plt.subplot(4, 1, 2)
plt.plot(t, np.cumsum(input_signal - feedback_signal), label='Integrator Output')
plt.title('Integrator Output')
plt.xlabel('Time [s]')
plt.ylabel('Amplitude')
plt.grid(True)

# デルタシグマ変調出力のプロット
plt.subplot(4, 1, 3)
plt.plot(t, quantizer_output, label='Delta-Sigma Modulated Signal')
plt.title('Delta-Sigma Modulated Signal')
plt.xlabel('Time [s]')
plt.ylabel('Output')
plt.grid(True)

# フィルタ後の信号のプロット
plt.subplot(4, 1, 4)
plt.plot(t, filtered_signal, label='Filtered Signal (Analog Output)')
plt.title('Filtered Signal (Analog Output)')
plt.xlabel('Time [s]')
plt.ylabel('Amplitude')
plt.grid(True)

plt.tight_layout()
plt.show()
```



```

import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import butter, lfilter

# Parameters
fs = 1000  # Sampling frequency
f = 5      # Frequency of the input signal
t = np.arange(0, 1, 1/fs)  # Time vector
input_signal = 0.5 * np.sin(2 * np.pi * f * t)  # Input signal

# Delta-Sigma Modulation
integrator_output = np.zeros(len(t))
quantizer_output = np.zeros(len(t))
feedback_signal = np.zeros(len(t))

for i in range(1, len(t)):
    integrator_output[i] = integrator_output[i-1] + input_signal[i] - feedback_signal[i-1]
    quantizer_output[i] = 1 if integrator_output[i] >= 0 else -1
    feedback_signal[i] = quantizer_output[i]

# Low-pass filter design
def butter_lowpass(cutoff, fs, order=5):
    nyquist = 0.5 * fs
    normal_cutoff = cutoff / nyquist
    b, a = butter(order, normal_cutoff, btype='low', analog=False)
    return b, a

def lowpass_filter(data, cutoff, fs, order=5):
    b, a = butter_lowpass(cutoff, fs, order=order)
    y = lfilter(b, a, data)
    return y

# Apply low-pass filter to integrator output
cutoff_frequency = 10  # Cutoff frequency of the LPF
filtered_output = lowpass_filter(integrator_output, cutoff_frequency, fs)

# Plotting
plt.figure(figsize=(12, 10))

plt.subplot(4, 1, 1)
plt.plot(t, input_signal, label='Input Signal')
plt.title('Delta-Sigma Modulation')
plt.legend()
plt.grid()

plt.subplot(4, 1, 2)
plt.plot(t, integrator_output, label='Integrator Output')
plt.legend()
plt.grid()

plt.subplot(4, 1, 3)
plt.plot(t, quantizer_output, label='Quantizer Output')
plt.legend()
plt.grid()

plt.subplot(4, 1, 4)
plt.plot(t, filtered_output, label='Filtered Output (LPF)')
plt.legend()
plt.grid()

plt.tight_layout()
plt.show()
```





```

import numpy as np
import matplotlib.pyplot as plt

# Function to compute the noise shaping characteristic
def noise_shaping(order, f):
    z = np.exp(2j * np.pi * f)
    H = (1 - z**-1)**order
    return np.abs(H)**2

# Frequency range
frequencies = np.linspace(0, 0.5, 1000)  # Normalized frequency (0 to 0.5)

# Plot noise shaping characteristics for different orders
plt.figure(figsize=(10, 6))
for order in range(1, 8):
    plt.plot(frequencies, 10 * np.log10(noise_shaping(order, frequencies)), label=f'Order {order}')

plt.xlabel('Normalized Frequency (f)')
plt.ylabel('Noise Power (dB)')
plt.title('Noise Shaping Characteristics of ΔΣ Modulators')
plt.legend()
plt.grid(True)
plt.ylim(-100, 60)
plt.show()



```

```
import numpy as np
import matplotlib.pyplot as plt

# Parameters for the input sine wave
amplitude = 0.8   # Amplitude of the sine wave
frequency = 1000  # Frequency of the sine wave in Hz
sampling_rate = 50000  # Sampling rate in Hz
duration = 0.005  # Duration in seconds

# Generate the time axis
t = np.linspace(0, duration, int(sampling_rate * duration), endpoint=False)

# Generate the input sine wave
input_signal = amplitude * np.sin(2 * np.pi * frequency * t)

# Delta-Sigma Modulation Parameters
integrator_output = 0
quantizer_output = np.zeros_like(t)
integrated_signal = np.zeros_like(t)

# Delta-Sigma Modulation (First-Order)
for i in range(1, len(t)):
    integrator_output += input_signal[i] - quantizer_output[i-1]  # Integrate the error
    integrated_signal[i] = integrator_output
    # 1-bit Quantizer
    if integrator_output > 0:
        quantizer_output[i] = 1
    else:
        quantizer_output[i] = -1

# Convert quantizer output back to analog using a simple low-pass filter
# Low-pass filter parameters
cutoff_freq = 2000  # Cutoff frequency in Hz
RC = 1 / (2 * np.pi * cutoff_freq)
alpha = 1 / (RC * sampling_rate)

# Apply low-pass filter to quantizer output
analog_output = np.zeros_like(t)
for i in range(1, len(t)):
    analog_output[i] = alpha * quantizer_output[i] + (1 - alpha) * analog_output[i-1]

# Plot the results
plt.figure(figsize=(14, 8))

# Plot input signal
plt.subplot(3, 1, 1)
plt.plot(t, input_signal, label='Input Sine Wave')
plt.title('Input Sine Wave')
plt.xlabel('Time [s]')
plt.ylabel('Amplitude')
plt.grid(True)
plt.legend()

# Plot 1-bit digital signal (quantizer output)
plt.subplot(3, 1, 2)
plt.step(t, quantizer_output, label='1-Bit Digital Signal (Quantizer Output)', where='mid')
plt.title('1-Bit Digital Signal (Delta-Sigma Modulation Output)')
plt.xlabel('Time [s]')
plt.ylabel('Level')
plt.grid(True)
plt.legend()

# Plot analog output (low-pass filter applied to digital signal)
plt.subplot(3, 1, 3)
plt.plot(t, analog_output, label='Reconstructed Analog Output')
plt.title('Reconstructed Analog Signal after Low-Pass Filter')
plt.xlabel('Time [s]')
plt.ylabel('Amplitude')
plt.grid(True)
plt.legend()

# Display the plots
plt.tight_layout()
plt.show()

```

```

import numpy as np
import matplotlib.pyplot as plt

# Parameters
fs = 1000  # Sampling frequency
f = 10     # Frequency of the input sine wave
T = 1      # Duration of the signal in seconds

# Time vector
t = np.linspace(0.0, T, T*fs, endpoint=False)

# Input sine wave (A)
input_signal = np.sin(2*np.pi*f*t)

# Delta-Sigma Modulation (B) - Assuming 5-bit quantization
quantized_signal = np.round(input_signal * 16) / 16  # 5-bit quantization

# Plotting
plt.figure(figsize=(10, 6))

plt.subplot(3, 1, 1)
plt.plot(t, input_signal, label='Input Sine Wave')
plt.title('Input Sine Wave')
plt.xlabel('Time')
plt.ylabel('Amplitude')
plt.grid(True)
plt.legend()

plt.subplot(3, 1, 2)
plt.plot(t, quantized_signal, label='5-bit DSM Output')
plt.title('5-bit DSM Output')
plt.xlabel('Time')
plt.ylabel('Quantized Amplitude')
plt.grid(True)
plt.legend()

plt.subplot(3, 1, 3)
plt.plot(t, input_signal - quantized_signal, label='Difference (A - B)')
plt.title('Difference between Input and DSM Output')
plt.xlabel('Time')
plt.ylabel('Amplitude Difference')
plt.grid(True)
plt.legend()

plt.tight_layout()
plt.show()

```




```
import numpy as np
import matplotlib.pyplot as plt
from scipy import fft, signal

# サンプリング周波数とシミュレーション時間
fs = 1000  # サンプリング周波数 (Hz)
T = 1.0    # シミュレーション時間 (秒)
t = np.arange(0, T, 1/fs)  # 時間ベクトル

# 入力信号 (例としてsin波を使用)
input_signal = np.sin(2 * np.pi * t)

# 量子化 (6ビット)
quantized_signal = np.round(input_signal * (2**5)) / (2**5)

# 量子化ノイズ
quantization_noise = input_signal - quantized_signal

# プロット
plt.figure(figsize=(12, 10))

# 入力信号、量子化された信号、量子化ノイズ
plt.subplot(3, 1, 1)
plt.plot(t, input_signal, label='Input Signal')
plt.plot(t, quantized_signal, label='Quantized Signal')
plt.plot(t, quantization_noise, label='Quantization Noise')
plt.title('Input Signal, Quantized Signal, Quantization Noise')
plt.xlabel('Time (s)')
plt.ylabel('Amplitude')
plt.legend()

# 量子化ノイズの周波数成分
plt.subplot(3, 2, 3)
freqs = fft.fftfreq(len(t), 1/fs)
fft_vals = fft.fft(quantization_noise)
plt.semilogy(freqs[:len(freqs)//2], np.abs(fft_vals[:len(freqs)//2]))
plt.title('Quantization Noise Frequency Components')
plt.xlabel('Frequency (Hz)')
plt.ylabel('Amplitude')

# ハイパスフィルタ処理
cutoff_hz = 50  # カットオフ周波数 (Hz)
b, a = signal.butter(4, cutoff_hz/(fs/2), 'high')
filtered_noise_hp = signal.lfilter(b, a, quantization_noise)

# ハイパスフィルタ処理後の周波数成分
plt.subplot(3, 2, 4)
fft_vals_hp = fft.fft(filtered_noise_hp)
plt.semilogy(freqs[:len(freqs)//2], np.abs(fft_vals_hp[:len(freqs)//2]))
plt.title('Highpass Filtered Quantization Noise Frequency Components')
plt.xlabel('Frequency (Hz)')
plt.ylabel('Amplitude')

# バンドパスフィルタ処理
bandpass_range = (20, 200)  # バンドパスの周波数帯域 (Hz)
b, a = signal.butter(4, [bandpass_range[0]/(fs/2), bandpass_range[1]/(fs/2)], 'bandpass')
filtered_noise_bp = signal.lfilter(b, a, quantization_noise)

# バンドパスフィルタ処理後の周波数成分
plt.subplot(3, 2, 5)
fft_vals_bp = fft.fft(filtered_noise_bp)
plt.semilogy(freqs[:len(freqs)//2], np.abs(fft_vals_bp[:len(freqs)//2]))
plt.title('Bandpass Filtered Quantization Noise Frequency Components')
plt.xlabel('Frequency (Hz)')
plt.ylabel('Amplitude')

plt.tight_layout()
plt.show()
```
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3759109/72681287-b484-aca5-8183-ed7685d26cdb.png)




```
import numpy as np
import matplotlib.pyplot as plt
from scipy import signal

# Parameters
fs = 1000  # Sampling frequency
T = 1  # Measurement time
t = np.linspace(0, T, int(fs * T), endpoint=False)  # Time axis

# Input sine wave (assuming a simple sine wave for demonstration)
input_signal = np.sin(2 * np.pi * t)

# Delta-Sigma modulation quantization (simplified to binary output here)
quantized_signal = np.zeros_like(input_signal)  # Initialize quantized signal

# Delta-Sigma modulator initial values
sigma1 = 0
sigma2 = 0

# Main loop: Delta-Sigma modulation
for i in range(len(input_signal)):
    # Quantizer
    if sigma2 > 0:
        quantized_signal[i] = 1
    else:
        quantized_signal[i] = 0
    
    # Update equations
    sigma1 = sigma1 + input_signal[i] - quantized_signal[i]
    sigma2 = sigma2 + sigma1 - quantized_signal[i]

# Designing a low-pass filter
order = 4
cutoff_freq = 10  # Hz
b, a = signal.butter(order, cutoff_freq / (fs / 2), 'low')

# Applying the filter to the quantized signal
filtered_signal = signal.filtfilt(b, a, quantized_signal)

# Plotting
plt.figure(figsize=(10, 8))

# Plot 1: Input signal
plt.subplot(4, 1, 1)
plt.plot(t, input_signal, label='Input Signal $x(t) = \sin(2\pi t)$')
plt.xlabel('Time [s]')
plt.ylabel('Amplitude')
plt.legend()

# Plot 2: Quantized signal (binary output)
plt.subplot(4, 1, 2)
plt.step(t, quantized_signal, label='Quantized Signal $y[n]$')
plt.xlabel('Time [s]')
plt.ylabel('Binary Output')
plt.legend()

# Plot 3: Filtered signal
plt.subplot(4, 1, 3)
plt.plot(t, filtered_signal, label='Filtered Signal')
plt.xlabel('Time [s]')
plt.ylabel('Amplitude')
plt.legend()

# Plot 4: FFT of the quantized signal
plt.subplot(4, 1, 4)
freq = np.fft.fftfreq(len(t), d=1/fs)
fft_quantized = np.fft.fft(quantized_signal)
plt.plot(freq, np.abs(fft_quantized), label='FFT of Quantized Signal')
plt.xlabel('Frequency [Hz]')
plt.ylabel('Magnitude')
plt.legend()

plt.tight_layout()
plt.show()
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3759109/4ac162f7-5393-4dff-0a03-545eb8510531.png)



```
import numpy as np
import matplotlib.pyplot as plt

# Parameters
fs = 1000  # Sampling frequency (Hz)
T = 1      # Duration of sampling (seconds)
f_analog = 10  # Analog sine wave frequency (Hz)
quantization_bits = 6  # Number of bits for quantization

# Time vector for analog signal
t_analog = np.linspace(0, T, fs*T, endpoint=False)
# Analog sine wave signal
analog_signal = np.sin(2 * np.pi * f_analog * t_analog)

# Sampling (creating digital signal)
digital_signal = analog_signal.copy()

# Digital filter (1/(1-z^(-1)))
filtered_signal = np.zeros_like(digital_signal)
for i in range(1, len(filtered_signal)):
    filtered_signal[i] = digital_signal[i] + digital_signal[i-1]

# Quantization (reduce to quantization_bits)
quantized_signal = np.round(filtered_signal * (2**(quantization_bits-1))) / (2**(quantization_bits-1))

# Feedback filter (z^(-1))
feedback_signal = np.zeros_like(quantized_signal)
for i in range(len(feedback_signal)-1):
    feedback_signal[i+1] = quantized_signal[i]

# Plotting
plt.figure(figsize=(10, 6))

plt.subplot(2, 2, 1)
plt.plot(t_analog, analog_signal)
plt.title('Analog Input Sine Wave')
plt.xlabel('Time (s)')
plt.ylabel('Amplitude')

plt.subplot(2, 2, 2)
plt.stem(t_analog, digital_signal, use_line_collection=True)
plt.title('Digital Sine Wave (Sampled)')
plt.xlabel('Time (s)')
plt.ylabel('Amplitude')

plt.subplot(2, 2, 3)
plt.plot(t_analog, filtered_signal)
plt.title('Digital Signal after Filter (1/(1-z^(-1)))')
plt.xlabel('Time (s)')
plt.ylabel('Amplitude')

plt.subplot(2, 2, 4)
plt.stem(t_analog, quantized_signal, use_line_collection=True)
plt.title('Quantized Digital Signal')
plt.xlabel('Time (s)')
plt.ylabel('Amplitude')

plt.tight_layout()
plt.show()

```


https://python-deltasigma.readthedocs.io/en/latest/#plotting-and-data-display-utilities



LTスパイスのシミュレーション
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3759109/266e42e8-ffaf-e21c-b1d2-53d1e1d22821.png)
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3759109/6df53686-6a29-9a13-4440-9efc0687f914.png)






```

import numpy as np
import matplotlib.pyplot as plt

# Parameters
t = np.linspace(0, 1, 1000)  # Time vector
f_sampling = 1000  # Sampling frequency
bits = 3  # Number of quantization bits

# Generate input signal (original sine wave)
input_signal = np.sin(2 * np.pi * t)

# Quantization process
quantized_signal = np.round(input_signal * (2**bits - 1)) / (2**bits - 1)

# Dithering (adding white noise)
dithered_signal = quantized_signal + np.random.uniform(-0.5/(2**bits), 0.5/(2**bits), size=len(t))

# Filters (example: 1 / (1 - z^(-1)))
filtered_signal = np.zeros_like(t)
for i in range(1, len(t)):
    filtered_signal[i] = dithered_signal[i] + dithered_signal[i-1]

# Plotting
plt.figure(figsize=(10, 6))

plt.subplot(411)
plt.plot(t, input_signal, label='Input Signal')
plt.xlabel('Time')
plt.ylabel('Amplitude')
plt.legend()

plt.subplot(412)
plt.plot(t, quantized_signal, label='Quantized Signal')
plt.xlabel('Time')
plt.ylabel('Amplitude')
plt.legend()

plt.subplot(413)
plt.plot(t, dithered_signal, label='Dithered Signal')
plt.xlabel('Time')
plt.ylabel('Amplitude')
plt.legend()

plt.subplot(414)
plt.plot(t, filtered_signal, label='Filtered Signal')
plt.xlabel('Time')
plt.ylabel('Amplitude')
plt.legend()

plt.tight_layout()
plt.show()
```


![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3759109/b7c81b04-8700-f90c-025e-d6155d48ecdc.png)


```

import numpy as np
import matplotlib.pyplot as plt

# パラメータ設定
N = 1024  # データ点数
t = np.arange(N)
fs = 1000  # サンプリング周波数 (Hz)
f_signal = 5  # 入力信号の周波数 (Hz)

# 入力信号 u(n) の生成
u = np.sin(2 * np.pi * f_signal * t / fs)

# 量子化誤差 e(n) の生成 (ランダム誤差)
e = np.random.normal(0, 0.1, N)

# 出力信号 v(n) の初期化
v = np.zeros(N)

# 2次デルタシグマ変調器の計算
for n in range(2, N):
    v[n] = u[n-1] + e[n] - e[n-1] - (e[n-1] - e[n-2])

# FFT の計算
fft_v = np.fft.fft(v)
fft_v = np.abs(fft_v[:N // 2]) * 2 / N  # 振幅を計算し、正規化
freqs = np.fft.fftfreq(N, 1 / fs)[:N // 2]  # 周波数軸

# プロット
plt.figure(figsize=(12, 6))

plt.subplot(2, 1, 1)
plt.plot(t, u, label='Input Signal u(n)')
plt.plot(t, v, label='Output Signal v(n)', alpha=0.75)
plt.xlabel('Time [samples]')
plt.ylabel('Amplitude')
plt.legend()
plt.title('Input and Output Time Series')

plt.subplot(2, 1, 2)
plt.plot(freqs, fft_v)
plt.xlabel('Frequency [Hz]')
plt.ylabel('Amplitude')
plt.title('FFT Spectrum of Output Signal v(n)')
plt.grid(True)

plt.tight_layout()
plt.show()

```


![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3759109/ca5d158b-5a2e-fe5f-cbd8-cfa5e9258346.png)

```

import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import butter, filtfilt

# アナログ入力信号の生成
fs = 1000  # サンプリング周波数
f = 5  # 入力信号の周波数
t = np.linspace(0, 1, fs)  # 時間軸
analog_signal = 0.5 * np.sin(2 * np.pi * f * t)  # 入力アナログ信号

# ΔΣ変調器の動作
for va in analog_signal:
    # 減算回路（Δ）の計算
    if integrator >= 0:
        # 積分器の出力が非負の場合
        delta = va - Vref
    else:
        # 積分器の出力が負の場合
        delta = va + Vref
    
    # 積分器の更新
    integrator += delta
    
    # 加算回路の出力（1ビットの生成）
    if integrator >= 0:
        bit = 1
    else:
        bit = 0
    
    # 出力ビットストリームにビットを追加
    output_bits.append(bit)
    
    # 積分器の出力を記録
    integrator_outputs.append(integrator)


# ローパスフィルタの適用
def butter_lowpass_filter(data, cutoff, fs, order=5):
    nyquist = 0.5 * fs
    normal_cutoff = cutoff / nyquist
    b, a = butter(order, normal_cutoff, btype='low', analog=False)
    y = filtfilt(b, a, data)
    return y

cutoff_frequency = 10  # カットオフ周波数
filtered_signal = butter_lowpass_filter(output_bits, cutoff_frequency, fs)

# プロット
fig, axs = plt.subplots(4, 1, figsize=(12, 10), sharex=True)

# アナログ入力信号のプロット
axs[0].plot(t, analog_signal, label='Analog Input Signal')
axs[0].set_title('Analog Input Signal')
axs[0].set_ylabel('Amplitude')
axs[0].legend()

# 積分器の出力のプロット
axs[1].plot(t, integrator_outputs, label='Integrator Output')
axs[1].set_title('Integrator Output')
axs[1].set_ylabel('Amplitude')
axs[1].legend()

# ΔΣ変調ビットストリームのプロット
axs[2].step(t, output_bits, where='mid', label='Delta-Sigma Modulated Bitstream')
axs[2].set_title('Delta-Sigma Modulated Bitstream')
axs[2].set_ylabel('Bit Value')
axs[2].legend()

# ローパスフィルタ後の信号のプロット
axs[3].plot(t, filtered_signal, label='Filtered Signal (Lowpass)')
axs[3].set_title('Filtered Signal (Lowpass)')
axs[3].set_xlabel('Time [s]')
axs[3].set_ylabel('Amplitude')
axs[3].legend()

plt.tight_layout()
plt.show()
```


![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3759109/15324229-c590-5cfb-6306-bfc778c888e9.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3759109/2a9db656-8687-00cd-68f9-b6f89afcbe0b.png)


```
import numpy as np
import matplotlib.pyplot as plt

# Parameters
Fs = 1000         # Sampling frequency (Hz)
Ts = 1 / Fs       # Sampling period
duration = 1.0    # Signal duration in seconds
f_input = 5       # Frequency of input signal (Hz)
amplitude = 1.0   # Amplitude of input signal

# Time vector
t = np.arange(0, duration, Ts)

# Input signal (sine wave)
input_signal = amplitude * np.sin(2 * np.pi * f_input * t)

# Delta-Sigma Modulator simulation
integrator_output = np.zeros_like(input_signal)
quantizer_output = np.zeros_like(input_signal)
feedback = np.zeros_like(input_signal)

# Initial states
integrator_state = 0.0
dac_output = 0.0

# Delta-Sigma Modulator Loop
for i in range(len(input_signal)):
    # Integrator: Accumulate input - feedback signal
    integrator_state += input_signal[i] - dac_output
    integrator_output[i] = integrator_state
    
    # Quantizer: Convert integrator output to 1-bit digital (1 or -1)
    if integrator_output[i] > 0:
        quantizer_output[i] = 1.0
    else:
        quantizer_output[i] = -1.0
    
    # DAC: Digital-to-Analog Conversion (1-bit DAC)
    dac_output = quantizer_output[i]
    
    # Feedback signal: Equal to DAC output
    feedback[i] = dac_output

# Plot the results
plt.figure(figsize=(10, 8))

# Plot input signal
plt.subplot(3, 1, 1)
plt.plot(t, input_signal, label='Input Signal')
plt.title('Input Signal')
plt.xlabel('Time [s]')
plt.ylabel('Amplitude')
plt.grid(True)


# Plot integrator output
plt.subplot(3, 1, 2)
plt.plot(t, integrator_output, label='Integrator Output', color='orange')
plt.title('Integrator Output')
plt.xlabel('Time [s]')
plt.ylabel('Amplitude')
plt.grid(True)

# Plot quantizer (1-bit output) signal
plt.subplot(3, 1, 3)
plt.step(t, quantizer_output, label='Quantizer Output (1-bit)', color='green')
plt.title('Quantizer Output (1-bit)')
plt.xlabel('Time [s]')
plt.ylabel('Digital Output')
plt.grid(True)

# Show the plots
plt.tight_layout()
plt.show()
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3759109/cf34a7fb-27fa-e270-564f-c0fef5d614cb.png)



https://lalsie.ist.hokudai.ac.jp/publication/dlcenter.php?fn=thesis/thesis_uchida.pdf

https://www.ssc.pe.titech.ac.jp/private/publications/2010/MathWorks_matsu_100602.pdf


```
import numpy as np
import matplotlib.pyplot as plt
from scipy.fft import fft
from scipy.signal.windows import hann  # Use Hann window for frequency analysis

# Parameters
Fs = 20000  # Sampling frequency
L = 10000   # Length of signal
Freq = 99   # Signal frequency

# Time and frequency vectors
t = np.arange(0, L) / Fs
f = Fs * np.arange(0, L // 2 + 1) / L

# Pure sine wave signal
pure = np.sin(2 * np.pi * Freq * t)

# 1-bit Delta-Sigma Modulator
dac = 0       # Initial DAC value
sigma_out = 0 # Initial integrator (sigma) register value
comp1_out = np.zeros(L)  # Output array for 1-bit modulator

# 1-bit Delta-Sigma loop
for i in range(L):
    delta_out = pure[i] - dac       # Difference between input and DAC output
    sigma_out += delta_out          # Integrator (accumulator)
    if sigma_out > 0:               # Comparator
        dac = 1
        comp1_out[i] = 1
    else:
        dac = -1
        comp1_out[i] = -1

# 2-bit Delta-Sigma Modulator
dac = 0       # Reset initial DAC value
sigma_out = 0 # Reset integrator value
comp2_out = np.zeros(L)  # Output array for 2-bit modulator

# 2-bit Delta-Sigma loop
for i in range(L):
    delta_out = pure[i] - dac       # Difference between input and DAC output
    sigma_out += delta_out          # Integrator (accumulator)
    if sigma_out > 1/3:
        dac = 1
        comp2_out[i] = 1
    elif 1/3 >= sigma_out > 0:
        dac = 1/3
        comp2_out[i] = 1/3
    elif 0 >= sigma_out > -1/3:
        dac = -1/3
        comp2_out[i] = -1/3
    else:
        dac = -1
        comp2_out[i] = -1

# Apply Hann window for FFT
win = hann(L)
pure_sp = np.abs(fft(pure * win)[:L//2+1])
pdm1_sp = np.abs(fft(comp1_out * win)[:L//2+1])
pdm2_sp = np.abs(fft(comp2_out * win)[:L//2+1])

# Convert to dB
pure_sp_dB = 20 * np.log10(pure_sp)
pdm1_sp_dB = 20 * np.log10(pdm1_sp)
pdm2_sp_dB = 20 * np.log10(pdm2_sp)

# Plot time-domain signals
plt.figure(figsize=(12, 10))

plt.subplot(3, 2, 1)
plt.plot(t, pure, label='Pure Signal')
plt.title("Pure Signal")
plt.xlabel("Time [s]")
plt.ylabel("Amplitude")
plt.grid(True)

plt.subplot(3, 2, 3)
plt.plot(t, comp1_out, label='1-bit Delta-Sigma Output', color='orange')
plt.title("1-bit Delta-Sigma Modulator Output")
plt.xlabel("Time [s]")
plt.ylabel("Amplitude")
plt.grid(True)

plt.subplot(3, 2, 5)
plt.plot(t, comp2_out, label='2-bit Delta-Sigma Output', color='green')
plt.title("2-bit Delta-Sigma Modulator Output")
plt.xlabel("Time [s]")
plt.ylabel("Amplitude")
plt.grid(True)

# Plot frequency-domain signals (Spectrum)
plt.subplot(1, 2, 2)
plt.plot(f, pure_sp_dB, label='Pure Signal')
plt.plot(f, pdm1_sp_dB, label='1-bit Delta-Sigma', color='orange')
plt.plot(f, pdm2_sp_dB, label='2-bit Delta-Sigma', color='green')
plt.title("Spectrum of Pure Signal and Delta-Sigma Modulator Outputs")
plt.xlabel("Frequency [Hz]")
plt.ylabel("Amplitude [dB]")
plt.legend(loc='lower right')
plt.grid(True)
plt.xlim(1, Fs/2)
plt.ylim(-80, max(pure_sp_dB))

# Show plots
plt.tight_layout()
plt.show()
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3759109/e1be4323-5bb4-9953-5454-82961c49a154.png)

```

import numpy as np
import matplotlib.pyplot as plt

# Parameters
n_samples = 100  # Number of samples
f = 0.05  # Frequency of the input sine wave

# Generate input sine wave
n = np.arange(n_samples)
x = np.sin(2 * np.pi * f * n)

# Initialize output array
y = np.zeros(n_samples)

# Implement the system
for i in range(1, n_samples):
    y[i] = x[i] - np.sign(y[i-1]) + y[i-1]

# Plot input and output signals
plt.figure(figsize=(12, 6))
plt.subplot(2, 1, 1)
plt.plot(n, x, label='Input x[n]')
plt.title("Input Signal")
plt.xlabel("n")
plt.ylabel("x[n]")
plt.grid()
plt.legend()

plt.subplot(2, 1, 2)
plt.plot(n, y, label='Output y[n]', color='orange')
plt.title("Output Signal")
plt.xlabel("n")
plt.ylabel("y[n]")
plt.grid()
plt.legend()

plt.tight_layout()
plt.show()

```

```


import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import lfilter
from scipy.fft import fft

# Parameters
fs = 26e6  # Sampling frequency
fin = 1e4  # Input frequency (10 kHz)
n_samples = 2**14  # Number of samples

# Generate input sine wave
t = np.arange(n_samples) / fs
x = 0.5 * np.sin(2 * np.pi * fin * t)

# First-order Delta-Sigma Modulator
def delta_sigma_1st_order(x):
    y = np.zeros_like(x)
    integrator = 0
    for i in range(len(x)):
        integrator += x[i] - y[i-1] if i > 0 else x[i]
        y[i] = 1 if integrator >= 0 else -1
    return y

# Second-order Delta-Sigma Modulator
def delta_sigma_2nd_order(x):
    y = np.zeros_like(x)
    integrator1 = 0
    integrator2 = 0
    for i in range(len(x)):
        integrator1 += x[i] - y[i-1] if i > 0 else x[i]
        integrator2 += integrator1 - y[i-1] if i > 0 else integrator1
        y[i] = 1 if integrator2 >= 0 else -1
    return y

# Apply 1st and 2nd order Delta-Sigma Modulators
y_1st = delta_sigma_1st_order(x)
y_2nd = delta_sigma_2nd_order(x)

# Compute FFT
def compute_fft(y, fs):
    n = len(y)
    yf = fft(y) / n
    xf = np.fft.fftfreq(n, 1 / fs)
    return xf[:n // 2], 20 * np.log10(np.abs(yf[:n // 2]))

xf, psd_1st = compute_fft(y_1st, fs)
_, psd_2nd = compute_fft(y_2nd, fs)

# Plot results
plt.figure(figsize=(12, 8))
plt.semilogx(xf, psd_1st, label="1st Order Noise Shaping", color='red')
plt.semilogx(xf, psd_2nd, label="2nd Order Noise Shaping", color='blue')
plt.xlabel("Frequency (Hz)")
plt.ylabel("Magnitude (dBFS)")
plt.title("1st and 2nd Order Delta-Sigma Modulator Noise Shaping")
plt.legend()
plt.grid(True)
plt.xlim([1e3, 1e7])
plt.ylim([-200, 0])
plt.show()
```



![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3759109/580cd960-5a55-8f93-11ba-58d1a72d7b5e.png)


```
import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import butter, lfilter, freqz

# Parameters for the input signal
fs = 1000  # Sampling frequency
t = np.linspace(0, 1, fs, endpoint=False)  # Time vector for 1 second
freq = 5  # Frequency of the input sine wave
input_signal = np.sin(2 * np.pi * freq * t)  # Generate a sine wave

# Function to perform delta-sigma modulation (DSM)
def delta_sigma_modulation(input_signal):
    """Simple Delta-Sigma Modulator"""
    integrator = 0
    output = np.zeros(len(input_signal))
    for i in range(len(input_signal)):
        integrator += input_signal[i]
        output[i] = 1 if integrator >= 0 else -1
        integrator -= output[i]
    return output

# Apply Delta-Sigma Modulation
dsm_output = delta_sigma_modulation(input_signal)

# Digital-to-Analog Conversion (DAC) using a low-pass filter
def butter_lowpass_filter(data, cutoff, fs, order=5):
    """Apply a Butterworth low-pass filter to the given data."""
    nyq = 0.5 * fs  # Nyquist Frequency
    normal_cutoff = cutoff / nyq  # Normalize cutoff frequency
    b, a = butter(order, normal_cutoff, btype='low', analog=False)
    y = lfilter(b, a, data)
    return y

# Parameters for low-pass filter
cutoff_frequency = 10  # Cutoff frequency for the low-pass filter (Hz)
output_signal = butter_lowpass_filter(dsm_output, cutoff_frequency, fs)

# FFT of input and output signals
input_fft = np.abs(np.fft.fft(input_signal)) / len(input_signal)
output_fft = np.abs(np.fft.fft(output_signal)) / len(output_signal)
xf = np.fft.fftfreq(len(input_signal), 1 / fs)  # Frequency bins

# Plotting the results
fig, axs = plt.subplots(3, 1, figsize=(10, 12))

# Input signal (time domain)
axs[0].plot(t, input_signal, label='Input Signal')
axs[0].set_title("Input Signal (Time Domain)")
axs[0].set_xlabel("Time [s]")
axs[0].set_ylabel("Amplitude")
axs[0].grid(True)

# Output signal after DSM and LPF (time domain)
axs[1].plot(t, output_signal, label='Output Signal (After DSM and LPF)', color='orange')
axs[1].set_title("Output Signal (Time Domain)")
axs[1].set_xlabel("Time [s]")
axs[1].set_ylabel("Amplitude")
axs[1].grid(True)

# FFT of input and output signals (frequency domain)
axs[2].plot(xf[:len(xf)//2], input_fft[:len(xf)//2], label='Input Signal FFT')
axs[2].plot(xf[:len(xf)//2], output_fft[:len(xf)//2], label='Output Signal FFT', color='orange')
axs[2].set_title("Frequency Spectrum of Input and Output Signals")
axs[2].set_xlabel("Frequency [Hz]")
axs[2].set_ylabel("Magnitude")
axs[2].grid(True)
axs[2].legend()

# Display plots
plt.tight_layout()
plt.show()

```

```

import numpy as np
import matplotlib.pyplot as plt

def delta_sigma_modulator_2bit(input_signal, n_levels=4):
    """
    2ビットのデルタシグマ変調器を実装する関数。

    Parameters:
    input_signal (numpy array): 入力信号（0〜1の範囲のアナログ信号）
    n_levels (int): 定量化レベルの数（デフォルトは4、2ビットのため）

    Returns:
    numpy array: 2ビットのデルタシグマ変調器の出力信号
    """
    # 入力信号の正規化（-1 〜 1の範囲にスケーリング）
    input_signal = 2 * (input_signal - np.min(input_signal)) / np.ptp(input_signal) - 1

    # 初期化
    integrator_output = 0       # 積分器の初期出力
    quantizer_output = 0        # 定量器の初期出力
    output_signal = []          # 変調器の出力信号
    quantizer_levels = np.linspace(-1, 1, n_levels)  # 定量化レベルの生成

    for sample in input_signal:
        # 積分器の出力を計算
        integrator_output += sample - quantizer_output

        # 定量器の出力（2ビット定量化）
        quantizer_output = quantizer_levels[np.argmin(np.abs(integrator_output - quantizer_levels))]
        
        # 変調器の出力を記録
        output_signal.append(quantizer_output)

    return np.array(output_signal)

# テスト用のサイン波入力信号を生成
fs = 1000  # サンプリング周波数
f = 5      # 信号の周波数
t = np.linspace(0, 1, fs, endpoint=False)
input_signal = 0.5 * (1 + np.sin(2 * np.pi * f * t))  # 0〜1の範囲のサイン波信号

# 2ビットデルタシグマ変調を適用
output_signal = delta_sigma_modulator_2bit(input_signal)

# 結果のプロット
plt.figure(figsize=(12, 6))

# 入力信号のプロット
plt.subplot(2, 1, 1)
plt.plot(t, input_signal, label='Input Signal')
plt.title('Input Signal')
plt.xlabel('Time [s]')
plt.ylabel('Amplitude')
plt.grid(True)
plt.legend()

# 2ビットデルタシグマ変調器の出力信号のプロット
plt.subplot(2, 1, 2)
plt.step(t, output_signal, where='mid', label='2-Bit Delta-Sigma Output', color='orange')
plt.title('2-Bit Delta-Sigma Modulator Output')
plt.xlabel('Time [s]')
plt.ylabel('Output (2-bit)')
plt.grid(True)
plt.legend()

# プロットを表示
plt.tight_layout()
plt.show()
```

```
import numpy as np
import matplotlib.pyplot as plt

# パラメータ設定
Fs = 1000  # サンプリング周波数
T = 1.0    # 信号長さ（秒）
t = np.linspace(0, T, int(Fs*T), endpoint=False)  # 時間軸
f_in = 5   # 入力信号の周波数

# 入力信号の生成 (サイン波)
input_signal = 0.5 * np.sin(2 * np.pi * f_in * t)

# 3ビット量子化器の設定（3ビット -> 2^3 = 8レベルの量子化）
quantizer_levels = 2**3
quantization_step = 1 / quantizer_levels  # 量子化ステップ幅
quantizer_min = -0.5  # 量子化範囲の最小値
quantizer_max = 0.5   # 量子化範囲の最大値

# デルタシグマ変調器の初期化
accumulator = 0
output_signal = np.zeros_like(input_signal)

# デルタシグマ変調器の実行
for i, sample in enumerate(input_signal):
    accumulator += sample  # 累積和を計算
    quantized_value = np.clip(accumulator, quantizer_min, quantizer_max)  # 量子化器への入力
    quantized_output = np.round(quantized_value / quantization_step) * quantization_step  # 量子化
    output_signal[i] = quantized_output  # 出力信号に量子化値を格納
    accumulator -= quantized_output  # 量子化誤差を累積和から差し引き

# 入力信号と出力信号のプロット
plt.figure(figsize=(12, 6))
plt.plot(t, input_signal, label='Input Signal (Analog)', color='blue')
plt.step(t, output_signal, label='Output Signal (3-bit Delta-Sigma)', color='orange', where='mid')
plt.xlabel('Time [s]')
plt.ylabel('Amplitude')
plt.title('3-bit Delta-Sigma Modulation')
plt.legend()
plt.grid()
plt.show()
```
```

import numpy as np
import matplotlib.pyplot as plt
import control as ctrl
from scipy import signal

# デルタシグマ変調器のパラメータ
order = 7  # ループフィルタの次数
fs = 1.0  # サンプリング周波数（正規化）
f0 = 0.1  # 入力信号の正規化周波数
A_in = 0.5  # 入力信号の振幅

# ループフィルタの設計（単純な7次積分器）
num = [1]
den = [1] + [0]*(order)
loop_filter = ctrl.TransferFunction(num, den)

# 量子化器（1ビット量子化器として近似）
def quantizer(x):
    return np.where(x >= 0, 1, -1)

# シミュレーション時間
T = 1000  # 総サンプル数
t = np.arange(T) / fs

# 入力信号
u = A_in * np.sin(2 * np.pi * f0 * t)

# シミュレーションの初期化
y = np.zeros(T)
v = np.zeros(T)
error = np.zeros(T)
integrator = 0.0

# デルタシグマ変調器のシミュレーション
for n in range(T):
    integrator += u[n] - v[n-1] if n > 0 else u[n]
    y[n] = integrator
    v[n] = quantizer(y[n])

# SNRの計算
def calculate_snr(signal, noise):
    power_signal = np.mean(signal**2)
    power_noise = np.mean(noise**2)
    return 10 * np.log10(power_signal / power_noise)

# 入力信号と出力信号のパワーを計算
signal_power = np.mean(u**2)
noise_power = np.mean((u - v)**2)
snr = calculate_snr(u, u - v)

print(f"SNR: {snr:.2f} dB")

# パワースペクトル密度のプロット
f, Pxx_den = signal.periodogram(v, fs)
plt.figure(figsize=(10, 6))
plt.semilogy(f, Pxx_den, label='Delta-Sigma Output')
plt.title('Power Spectral Density of Delta-Sigma Modulator Output')
plt.xlabel('Normalized Frequency')
plt.ylabel('Power/Frequency (dB)')
plt.legend()
plt.grid()
plt.show()

# 入力振幅に対するSNRのプロット
A_in_values = np.linspace(0.1, 1.0, 10)
snr_values = []

for A in A_in_values:
    u = A * np.sin(2 * np.pi * f0 * t)
    y = np.zeros(T)
    v = np.zeros(T)
    integrator = 0.0
    for n in range(T):
        integrator += u[n] - v[n-1] if n > 0 else u[n]
        y[n] = integrator
        v[n] = quantizer(y[n])
    noise = u - v
    snr_val = calculate_snr(u, noise)
    snr_values.append(snr_val)

plt.figure(figsize=(10, 6))
plt.plot(A_in_values, snr_values, marker='o')
plt.title('SNR vs Input Amplitude')
plt.xlabel('Input Amplitude (Normalized)')
plt.ylabel('SNR (dB)')
plt.grid()
plt.show()
```

```
import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import butter, filtfilt
from scipy.fftpack import fft

# Function for FFT and Power Spectrum Calculation
def calculate_fft(signal, fs):
    n = len(signal)
    freqs = np.fft.fftfreq(n, d=1/fs)
    fft_values = np.abs(fft(signal)) / n  # Normalize FFT magnitude
    power = np.abs(fft_values)**2  # Power spectrum
    power_db = 10 * np.log10(power)  # Power in decibels (dB)
    return freqs[:n // 2], power_db[:n // 2]  # Return positive frequencies only

# Low-pass filter function
def low_pass_filter(signal, cutoff, fs):
    nyquist = 0.5 * fs
    normal_cutoff = cutoff / nyquist
    b, a = butter(5, normal_cutoff, btype='low', analog=False)
    return filtfilt(b, a, signal)

# Second-order Delta-Sigma Modulator function
def second_order_delta_sigma_modulator(input_signal, oversampling_rate=64, a1=1, a2=1, b=2):
    n = len(input_signal)
    output_signal = np.zeros(n * oversampling_rate)
    integrator1 = 0  # First integrator state
    integrator2 = 0  # Second integrator state
    quantizer_output = 0  # Last quantized value
    
    # Oversample the input signal
    oversampled_input = np.repeat(input_signal, oversampling_rate)
    
    for i in range(n * oversampling_rate):
        # First integrator
        integrator1 += oversampled_input[i] - quantizer_output * a1
        # Second integrator
        integrator2 += integrator1 - quantizer_output * a2
        # Quantize the output (1-bit quantizer)
        if integrator2 >= 0:
            quantizer_output = 1
        else:
            quantizer_output = -1
        # Feedback loop
        output_signal[i] = quantizer_output
    
    return output_signal

# Parameters
fs = 1000  # Sampling frequency (Hz)
nn = 13  # Number of samples
A = 1  # Amplitude of input signal
Fin = 5  # Input frequency (Hz)
t = np.arange(0, 1, 1/fs)  # Time vector (1 second)
input_signal = 0.5 * np.sin(2 * np.pi * Fin * t)  # A 5 Hz sine wave input

# Perform second-order Delta-Sigma modulation
oversampling_rate = 64
delta_sigma_output = second_order_delta_sigma_modulator(input_signal, oversampling_rate)

# Apply low-pass filter (decimation process)
cutoff_frequency = 10  # Set the cutoff frequency of the LPF to 10 Hz
filtered_signal = low_pass_filter(delta_sigma_output, cutoff=cutoff_frequency, fs=fs * oversampling_rate)

# FFT and Power Spectrum Calculation
fft_freqs, power_db = calculate_fft(filtered_signal, fs * oversampling_rate)

# Plot input, modulated signal, filtered signal, and power spectrum
plt.figure(figsize=(12, 10))

# Original Input Signal
plt.subplot(4, 1, 1)
plt.plot(t, input_signal, label="Input Signal (5 Hz Sine)")
plt.title("Input Signal (5 Hz Sine)")
plt.grid(True)

# Delta-Sigma Modulated Signal
oversampled_time = np.linspace(0, 1, len(delta_sigma_output))
plt.subplot(4, 1, 2)
plt.step(oversampled_time, delta_sigma_output, where='mid', label="Delta-Sigma Modulated Signal (1-bit)")
plt.title("Second-Order Delta-Sigma Modulated Signal (1-bit)")
plt.grid(True)

# Filtered Signal after Low-Pass Filter
decimated_time = np.linspace(0, 1, len(filtered_signal))
plt.subplot(4, 1, 3)
plt.plot(decimated_time, filtered_signal, label="Filtered Signal (After Low-Pass Filter)")
plt.title("Filtered Signal (After Decimation)")
plt.grid(True)

# FFT of Filtered Signal - Power Spectrum
plt.subplot(4, 1, 4)
plt.plot(fft_freqs, power_db, label="Power Spectrum (dB)")
plt.title("Power Spectrum (dB)")
plt.grid(True)
plt.xlabel("Frequency (Hz)")
plt.ylabel("Power (dB)")

plt.tight_layout()
plt.show()
```