NMRチャートシミュレーターを作った

NMRの測定点が少ない場合にどのようなチャートになるのかを試してみたいと友人に言われたので、作製しました。

仕様

測定点数、ピークを表示する位置、分裂幅、を指定すると、NMRチャートのようなチャートを生成します。
ipywidgetsを利用しているため、Jupyter Notebook環境(Google Colaboratry)が必要です。

ソースコード

import numpy as np
import matplotlib.pyplot as plt
from ipywidgets import interact, FloatSlider, FloatLogSlider


def make_NMR_chart(range_max, N, peaktops):
    """
    Generate an NMR (Nuclear Magnetic Resonance) chart.

    Parameters:
    - range_max (float): The maximum range of the NMR chart.
    - N (int): The number of data points in the NMR chart.
    - peaktops (list): A list of peak frequencies for the NMR chart.

    Returns:
    - freq (ndarray): The frequency values of the NMR chart.
    - signal_positive (ndarray): The positive part of the NMR signal.

    """
    mesuring_time = N / range_max / 2
    t = np.linspace(0, mesuring_time, N, endpoint=False)
    signal = np.zeros(N)
    
    for peak in peaktops:
        signal += np.cos(2 * np.pi * peak * t)
    signal *= np.exp(-0.4 * t) 
    
    # Add Gaussian noise (noise intensity is specified by the second argument)
    signal += np.random.normal(0, 0.01, N)
    
    signal_positive = np.real(np.fft.fft(signal))[:N//2]
    freq = np.fft.fftfreq(N, 1 / range_max / 2)[:N//2]
    
    return freq, signal_positive - np.min(signal_positive) / np.max(signal_positive)

def plot_NMR_chart( N, peakshift, J1, J2):
    N = int(N)
    peaktops = np.array([(-J1-J2)/2,(-J1+J2)/2,(J1-J2)/2,(J1+J2)/2] ) + peakshift
    x, y = make_NMR_chart(10, N, peaktops)
    plt.plot(x, y)
    plt.show()

interact(plot_NMR_chart,
        N=FloatLogSlider(value=2048, min=5, max=14, base = 2, step=1),
        peakshift=FloatSlider(value=5, min=4, max=6, step=0.1),
        J1=FloatSlider(value=0, min=0, max=2, step=0.01),
        J2=FloatSlider(value=0, min=0, max=2, step=0.01),)

感想

1%のガウシアンノイズを載せているので、チャートはそれっぽい感じになっているかと思います。

実際のNMRに比べると、測定範囲は少ないのに測定点が非常に多いので、測定点数を少なくしておくことが重要そうです。

ギャラリー

128点、J1 = 0.6、J2 = 0.3。
ピークトップと谷の間に測定点がないため、最低限の測定点しかないことになる。

64点、J1 = 0.6、J2 = 0.3

化学シフトの中心をずらしてみる。本来4つのピークがあるはずだが、3つしか見えていない

512点、J1 = 0.6、J2 = 0.3
ノイズが少し強めに出ている。ピーク形状はシフトを移動させてもそれほど変化しない。