読者です 読者をやめる 読者になる 読者になる

zack_ff’s blog

ff14廃人がなにか書く

【python】辞書型の初期化

沢山のキーを持つ辞書を0で初期化しようと思った時た時にやり方がわからなかったのでメモ

listNmae = ["tmp0", "tmp1", "tmp2", "tmp3", "tmp4"]
dic = {name : 0 for name in listName}

これでいいみたい

【python】quaternionの処理で詰まったのでメモ

iphoneで取得したquaternionをpcで受け取って向いている方向をaxes3dで表示しようとした時にquaternionをどう処理していいかわからなくて詰まった。

調べているとtransforms3dを使うとquaternionからmatrixに変換できるらしい

コードで書くと

#初期位置
initPosition = [0, 1, 0]
#受け取ったquaternion
q = [w, x, y, z]
#matを求める
mat = transforms3d.quaternion.quat2mat(q)
#ドット積
result = numpy.dot(mat, initPosition)

これで初期位置からの変位が求められた
と思う…

【python】頭部伝達関数を重畳保留法を用いて畳み込んでみた

今回は重畳保留法(overlap_save)を使って頭部伝達関数を畳み込んでみました。

・頭部伝達関数

前回と同じ

・重畳保留法

以下のサイトを参考に書きました。
原理はわかっていません…
元気があったら勉強してみようと思います

www.softist.com

上記のサイトを参考にパラメータを決めていくと

N(fftブロック長) = 1024
L(フィルタ長) = 512
M(音源の切り出し区間) = 1024

MとNが同じなのでプログラムではMの部分にNを使いました。



ソースコード

今回も平面のみでやりました

プログラムの流れは以下のとおりです

・音源からM点抜き出す
・L点のhrtfの左側を0で埋めてN点の配列を作る
・この2つをfftして掛けあわせ、ifftする
・ifftして得られたデータの最初からL点を抜き出してストリームに流し込む

import numpy
import pyaudio
import scipy.io.wavfile as wav

def load_elev0hrtf():
    elev0Hrtf_L = {}
    elev0Hrtf_R = {}

    for i in range(72):
        str_i = str(i * 5)

        if len(str_i) < 2:
            str_i = "00" + str_i
        elif len(str_i) < 3:
            str_i = "0" + str_i

        fileName = "L0e" + str_i + "a.dat"
        filePath = "../hrtfs/elev0/" + fileName
        test = open(filePath, "r").read().split("\n")

        data = []

        for item in test:
            if item != '':
                data.append(float(item))

        elev0Hrtf_L[i] = data

    for i in range(72):
        str_i = str(i * 5)

        if len(str_i) < 2:
            str_i = "00" + str_i
        elif len(str_i) < 3:
            str_i = "0" + str_i

        fileName = "R0e" + str_i + "a.dat"
        filePath = "../hrtfs/elev0/" + fileName
        test = open(filePath, "r").read().split("\n")

        data = []

        for item in test:
            if item != '':
                data.append(float(item))

        elev0Hrtf_R[i] = data

    return elev0Hrtf_L, elev0Hrtf_R

def convolution(soundData, hrtf, N, L):
    tmpFilter = numpy.zeros(N)
    tmpFilter[L:] += hrtf

    spectrum = numpy.fft.fft(soundData, n=N)
    hrtfFft = numpy.fft.fft(tmpFilter, n=N)
    add = spectrum * hrtfFft
    result = numpy.fft.ifft(add, n=N).real
    return result[:L]

def play_elev0(soundData, hrtfL, hrtfR, L, N, streamObj, position):
    index = 0
    while(soundData[index:].size > L):
        resultData = numpy.empty((0, 2), dtype=numpy.int16)

        convL = convolution(soundData[index:index + N, 0], hrtfL[position], N, L)
        convR = convolution(soundData[index:index + N, 1], hrtfR[position], N, L)

        for i in range(convL.size):
            resultData = numpy.append(resultData, numpy.array([[int(convL[i]), int(convR[i])]], dtype=numpy.int16), axis=0)

        streamObj.write(bytes(resultData))
        index += L

    streamObj.close()

#wavファイルのサウンドデータをロード
soundDataPath = "ongaku.wav"
rate, soundData = wav.read(soundDataPath)

#HRTFのロード
hrtfL, hrtfR = load_elev0hrtf()

#PyAudioの初期設定
p = pyaudio.PyAudio()
stream = p.open(format = 8,
                channels = 2,
                rate = rate,
                output = True)
#パラメータ設定
L = 512
N = 1024
#音の方向(左回り)
position = 18

play_elev0(soundData, hrtfL, hrtfR, L, N, stream, position)
p.terminate()

【python】頭部伝達関数を重畳加算法を使って畳み込んでみた

頭部伝達関数(hrtf)を重畳加算法(overlap_add)を使って音源に畳み込んでみました。

・頭部伝達関数

音源位置から両耳の入り口までの伝達関数

今回は以下のサイトで公開されているhrtfを使用させていただきました

Head Related Transfer Function (HRTF) Database

・overlap_add

長い信号とfirフィルタの離散畳み込みを分割して効率的に処理する手法らしい

重畳加算法 - Wikipedia

ウィキを見てもらえばわかると思う

N = L + M - 1 = 2^p (N:fftブロック長、L:切り出し区間、M:フィルタ長)

と書かれている

使用するhrtfは512点なのでとりあえず

N = 1024

L = 513

M = 512

ウィキの図を見てもらうとオーバーラップ部分がM-1となっているので

overLap = 511

とした

 

・音源

フリーのwavファイル音源を使用しました

 

ソースコード

 とりあえず今回は水平面のhrtfのみ使用しました

import numpy
import pyaudio
import scipy.io.wavfile as scw

def load_elev0hrtf():
    elev0Hrtf_L = {}
    elev0Hrtf_R = {}

    for i in range(72):
        str_i = str(i * 5)

        if len(str_i) < 2:
            str_i = "00" + str_i
        elif len(str_i) < 3:
            str_i = "0" + str_i

        fileName = "L0e" + str_i + "a.dat"
        filePath = "../hrtfs/elev0/" + fileName
        test = open(filePath, "r").read().split("\n")

        data = []

        for item in test:
            if item != '':
                data.append(float(item))

        elev0Hrtf_L[i] = data

    for i in range(72):
        str_i = str(i * 5)

        if len(str_i) < 2:
            str_i = "00" + str_i
        elif len(str_i) < 3:
            str_i = "0" + str_i

        fileName = "R0e" + str_i + "a.dat"
        filePath = "../hrtfs/elev0/" + fileName
        test = open(filePath, "r").read().split("\n")

        data = []

        for item in test:
            if item != '':
                data.append(float(item))

        elev0Hrtf_R[i] = data

    return elev0Hrtf_L, elev0Hrtf_R

def convolution(data, hrtf, N, L):
    spectrum = numpy.fft.fft(data, n = N)
    hrtf_fft = numpy.fft.fft(hrtf, n = N)
    add = spectrum * hrtf_fft
    result = numpy.fft.ifft(add, n = N)
    return_data = result.real
    return return_data[:L], return_data[L:]

def play_elev0(sound_data, N, L, hrtfL, hrtfR, position, overLap, streamObj):
    index = 0
    overLap_L = numpy.zeros(overLap)
    overLap_R = numpy.zeros(overLap)

    while(sound_data[index:].size > L):
        result_data = numpy.empty((0, 2), dtype=numpy.int16)

        tmp_conv_L, add_L = convolution(sound_data[index:index + L, 0], hrtfL[position], N, L)
        tmp_conv_R, add_R = convolution(sound_data[index:index + L, 1], hrtfR[position], N, L)

        tmp_conv_L[:overLap] += overLap_L
        tmp_conv_R[:overLap] += overLap_R

        overLap_L = add_L
        overLap_R = add_R

        for i in range(tmp_conv_L.size):
            result_data = numpy.append(result_data, numpy.array([[int(tmp_conv_L[i]), int(tmp_conv_R[i])]], dtype=numpy.int16), axis=0)

        streamObj.write(bytes(result_data))
        index += L

    streamObj.close()


#初期設定など
soundDataPath = "./ongaku.wav"
rate, soundData = scw.read(soundDataPath)

p = pyaudio.PyAudio()
stream = p.open(format = 8,
                channels = 2,
                rate = rate,
                output = True)

hrtf_L, hrtf_R = load_elev0hrtf()
N = 1024
L = 513
overLap = 511
#音の再生する方向(0が正面で半時計回りで71まで)
position = 18

play_elev0(soundData, N, L, hrtf_L, hrtf_R, position, overLap, stream)
p.terminate()

hrtfをロードして畳み込んで再生といった流れです

間違いなどあればコメントください