【python】連続ウェーブレット変換の縦軸目盛のずれについて

Question

前提

pythonを勉強し始めた駆け出し学生です。
pythonで連続ウェーブレット変換をしようと試みています。
参考にさせていただいたサイトのコードにおいて、実行結果の縦軸目盛の表示にずれが生じています。
用意した波の周波数は2[Hz]、7[Hz]、13[Hz]であり、図のy軸の値に1ずつ足さなければplotしたい値（2、7、13）にならない状態です。
(参考サイト：https://camedphone.com/archives/729#comment-288)

実現したいこと

• 縦軸目盛の本来の値の表示をしたい

該当のソースコード

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

width = 6
 
# マザーウェーブレット：モルレーウェーブレット
def morlet(x, f, width):
    sf = f / width
    st = 1 / (2 * math.pi * sf)
    A = 1 / (st * math.sqrt(2 * math.pi))
    h = -np.power(x, 2) / (2 * st**2)
    co1 = 1j * 2 * math.pi * f * x
    return A * np.exp(co1) * np.exp(h)
 
# 連続ウェーブレット変換
def mycwt(Fs, data, fmax, wavelet_R=2):
    # Fs:           サンプリング周波数
    # data:         信号
    # wavelet_R:    マザーウェーブレットの長さ(秒)
    # fmax:         解析する最大周波数
 
    Ts = 1 / Fs     # サンプリング時間幅
    data_length = len(data) # 信号のサンプル数を取得
    
    # マザーウェーブレットの範囲
    wavelet_length = np.arange(-wavelet_R, wavelet_R, Ts)
 
    # 連続ウェーブレット変換後のデータを格納する配列の作成
    wn = np.zeros([fmax, data_length])
 
    # 連続ウェーブレット変換の実行
    for i in range(0, fmax):
        wn[i,:] = np.abs(np.convolve(data, morlet(wavelet_length, i+1, width), mode='same'))
        wn[i,:] = (2 * wn[i, :] / Fs)**2
    
    return wn
 
# 連続ウェーブレット変換後のカラーマップ作成関数
def cwt_plot(CWT, sample_time, fmax, fig_title):
    plt.imshow(CWT, cmap='jet', aspect='auto',vmax=abs(CWT).max(), vmin=-abs(CWT).max())  
    plt.title(fig_title)
    plt.xlabel("time[ms]")
    plt.ylabel("frequency[Hz]")
    plt.axis([0, len(sample_time), 0, fmax-1])
    plt.colorbar()
 
 
if __name__ == "__main__":
    # 信号作成 ----------------------------------------
    Fs = 1000 # サンプリング周波数
    Ts = 1 / Fs # 1ステップあたりの時間幅
    time_S = 5 # 信号は5秒分
    t_data = np.arange(0,time_S, Ts) # 5秒分の時間配列
    f1, f2, f3 = 2, 7, 13 # 2[Hz], 7[Hz], 13[Hz]の信号
    

    
    # signal作成 ----------------------------------------
    amp = np.linspace(0.5, 2, num = Fs* time_S) # 増幅用
 
    ss01 = np.sin(2*np.pi*f1*t_data)* amp
    ss02 = np.sin(2*np.pi*f2*t_data)
    ss03 = np.sin(2*np.pi*f3*t_data)* np.flipud(amp)
 
    signal = ss01 + ss02 + ss03 # 3つの信号の合成信号
 
    # 連続ウェーブレット変換 ----------------------------------------
    fmax=20 # 解析する最大周波数
    cwt_signal = mycwt(Fs=Fs, data=signal, fmax=fmax)
 

    # signalの図
    plt.figure() 
    plt.title("signal")
    plt.plot(t_data, signal)
    plt.xlabel("time[s]")
 

 
    # signal02を連続ウェーブレット変換した時のカラーマップの図
    plt.figure()
    fig_title02 = "cwt signal"
    cwt_plot(cwt_signal, t_data, fmax, fig_title02)
    
    plt.show()

試したこと

軸ラベルの修正で色々検索し、下記のようなコードを試したのですがエラーも出てしまい、pythonの理解もあまり深くなく、進展が無い状況です。

from typing import List, Tuple

def get_ticks_label_set(all_labels:List, num:int)-> Tuple[List, List]:
    length = len(all_labels)
    step = length//(num-1)
    pick_positions = np.arange(0, length, step)
    pick_labels = all_labels[::step]
    return pick_positions, pick_labels
x_positions, x_labels = get_ticks_label_set(t, num=10)
print(x_positions)
print(x_labels)
y_positions, y_labels = get_ticks_label_set(freqs[::-1], num=10)
print(y_positions)
print(y_labels)

python初心者のため色々なサイトを参考にしながらプログラムを改良しようと試みたのですが、うまくいかず困っています。
ご教授お願いいたします。

Comments

[jbpb0]

「def mycwt(...」の
wn[i,:] = np.abs(np.convolve(data, morlet(wavelet_length, i+1, width), mode='same'))
で「i+1」としてるので、たとえば「wn[0,:]」には「周波数1」の結果が入ります

そこのforループを
for i in range(1, fmax):
と変えて、「morlet(wavelet_length, i, width)」とすれば、たとえば「wn[1,:]」に「周波数1」の結果が入るようになり、縦軸と合うと思います

[loss]

ありがとうございます。
見事に解決することができました！
pythonに慣れてきたと思っていたのですが、まだまだでした･･･プログラミングでは0から数えるということを改めて認識することができました。
より能力の向上に努めようと思います。この度はありがとうございました。

Answer 1

「def mycwt(...」の

        wn[i,:] = np.abs(np.convolve(data, morlet(wavelet_length, i+1, width), mode='same'))

で「i+1」としてるので、たとえば「wn[0,:]」には「周波数1」の結果が入ります

　
そこのforループを

    for i in range(1, fmax):

と変えて、「morlet(wavelet_length, i, width)」とすれば、たとえば「wn[1,:]」に「周波数1」の結果が入るようになり、縦軸と合うと思います