$30 off During Our Annual Pro Sale. View Details »

逆FM音源

Fadis
February 08, 2020

 逆FM音源

与えられた楽器の音色に近いFM音源パラメータの探索を勾配法で解く方法を解説します
これは2020年2月8日に行われた カーネル/VM探検隊@関西 10回目 での発表資料です
サンプルコード: https://github.com/Fadis/ifm

Fadis

February 08, 2020
Tweet

More Decks by Fadis

Other Decks in Programming

Transcript

  1. ٯFMԻݯ
    NAOMASA MATSUBAYASHI
    INVERTED FREQUENCY MODULATION SYNTHESIZER
    https://github.com/Fadis/ifm
    αϯϓϧίʔυ

    View Slide

  2. Իָͷԋ૗͸೾ͷੜ੒
    εϐʔΧʔͷ
    ΞΠίϯ
    ిѹͷ೾ ۭؾͷ೾
    Ի֊
    Իྔ
    υ
    ίϯϐϡʔλͰԻָΛ૗ͰΔʹ͸
    Ի֊ʹରԠ͢Δप೾਺ͷ೾Λੜ੒͠ͳ͚Ε͹ͳΒͳ͍
    ?

    View Slide

  3. Իݯνοϓͷొ৔
    εϐʔΧʔͷ
    ΞΠίϯ
    ిѹͷ೾ ۭؾͷ೾
    Ի֊
    Իྔ
    υ
    ͜ͷੜ੒Λߦ͏ϋʔυ΢ΣΞ͕ొ৔ͨ͠

    View Slide

  4. υ
    Ϩ
    ϛ
    ϑΝ
    ι
    ԣ࣠: प೾਺ ॎ࣠: ࣌ؒ(ͦΕͧΕͷԻ֊ʹ͍ͭͯ໐Γ࢝Ί͔Β2ඵ෼)

    View Slide

  5. 65.406Hz 130.81Hz 196.22Hz 261.63Hz 327.03Hz 392.44Hz 457.84Hz 523.25Hz 588.66Hz
    υ
    Ϩ
    ϛ
    ϑΝ
    ι

    View Slide

  6. 65.406Hz 130.81Hz 196.22Hz 261.63Hz 327.03Hz 392.44Hz 457.84Hz 523.25Hz 588.66Hz
    υ
    Ϩ
    ϛ
    ϑΝ
    ι
    Ի֊Λָ࣋ͭثͷԻ͸
    ͋Δप೾਺ͷ೾ͱ
    ͦͷ੔਺ഒͷप೾਺ͷ೾ͷॏͶ߹Θͤ
    ͜ͷ͏ͪ࠷΋௿͍प೾਺ͷ೾Λ
    جԻͱݺͿ
    Ի֊্͕͕Δͱ
    جԻͷप೾਺্͕͕͍ͬͯ͘

    View Slide

  7. 65.406Hz 130.81Hz 196.22Hz 261.63Hz 327.03Hz 392.44Hz 457.84Hz 52
    υ
    Ϩ
    ϛ
    ϑΝ
    ι
    جԻͷ੔਺ഒͷ
    प೾਺Ͱग़Δ೾Λ
    ഒԻͱݺͿ
    جԻͷप೾਺্͕͕Δͱ
    ഒԻͷप೾਺΋্͕Δ

    View Slide

  8. άϥϯυϐΞϊ
    ΞϧταοΫε
    ϏϒϥϑΥϯ
    όΠΦϦϯ
    ϚϦϯό

    View Slide

  9. 261.63Hz 523.25Hz 784.88Hz 1046.5Hz 1308.1Hz 1569.8Hz 1831.4Hz 2093.0Hz 2354.6Hz
    άϥϯυϐΞϊ
    ΞϧταοΫε
    ϏϒϥϑΥϯ
    όΠΦϦϯ
    ϚϦϯό

    View Slide

  10. 261.63Hz 523.25Hz 784.88Hz 1046.5Hz 1308.1Hz 1569.8Hz 1831.4Hz 2093.0Hz 2354.6Hz
    άϥϯυϐΞϊ
    ΞϧταοΫε
    ϏϒϥϑΥϯ
    όΠΦϦϯ
    ϚϦϯό
    Ի͕4ΦΫλʔϒ໨ͷυʹͳΔ৚͕݅
    261.63HzͷجԻΛ࣋ͭ͜ͱͳͷͰ
    Ͳͷָث΋جԻʹ͸
    େ͖ͳҧ͍͸ݟΒΕͳ͍

    View Slide

  11. 1046.5Hz 1308.1Hz 1569.8Hz 1831.4Hz 2093.0Hz
    άϥϯυϐΞϊ
    ΞϧταοΫε
    ϏϒϥϑΥϯ
    όΠΦϦϯ
    ϚϦϯό
    2616.3Hz 2877.8Hz 3139.5Hz 3401.1Hz
    ഒԻͷҧ͍͕
    ָثͷԻ৭ͷҧ͍Λ
    ੜΈग़͢

    View Slide

  12. άϥϯυϐΞϊ
    ΞϧταοΫε
    ϏϒϥϑΥϯ
    όΠΦϦϯ
    ϚϦϯό

    View Slide

  13. άϥϯυϐΞϊ
    ԣ࣠: प೾਺ ॎ࣠: ࣌ؒ
    t
    ͋Δ࣌ࠁʹ͓͚Δશͯͷप೾਺ͷ੒෼ͷ૯࿨ΛऔΔͱ
    ͓͓Αͦͦͷ࣌ࠁͷԻͷৼ෯( Իྔ)͕ಘΒΕΔ

    View Slide

  14. Τϯϕϩʔϓ
    ָثͷԻͷৼ෯͸࣌ؒมԽ͢Δ
    ಛʹϐΞϊ΍ϚϦϯόͷΑ͏ͳ
    ࠷ॳʹՃ͑ͨྗ͚ͩͰ࠷ޙ·Ͱ໐Δָث͸
    ៉ྷͳݮਰৼಈͷۂઢΛඳ͘

    View Slide

  15. άϥϯυϐΞϊ
    ࠷ॳʹՃ͑ͨྗ͚ͩͰ࠷ޙ·Ͱ໐ΔָثͰ͸
    प೾਺ͷߴ͍ഒԻ΄Ͳૣ͘ফ͑Δ܏޲ʹ͋Δ
    όΠΦϦϯ
    ܧଓతʹྗ͕Ճ͑ΒΕΔָثͰ͸
    ໐Β͢ͷΛ΍ΊΔ·Ͱಉׂ͡߹ͷഒԻ͕ग़Δ܏޲ʹ͋Δ

    View Slide

  16. ཁٻ͞ΕͨԻ֊ʹରԠ͢ΔجԻͷαΠϯ೾ͱ
    ָثʹԠͨ͡ഒԻͷαΠϯ೾Λੜ੒͢Δ
    ͦΕΒͷৼ෯͕࣮ࡍͷָثͱಉ͡Α͏ʹ࣌ؒมԽ͍ͯ͠Ε͹
    ͦͷָثͰཁٻ͞ΕͨԻ֊ͷԻΛ໐Βͨ͠Α͏ʹฉ͑͜Δ
    ϐΞϊͷԻͰ
    υ
    ഒ཰ ৼ෯











    View Slide

  17. ϐΞϊͷԻͰ
    υ
    ໰୊఺
    αΠϯ೾ͷܭࢉ͸҆͘ͳ͍
    ഒ཰ ৼ෯











    View Slide

  18. प೾਺มௐ
    f (t) = sin (2πωc
    t + B sin (2πωm
    t))

    View Slide

  19. ωc
    = ωm
    2ωc
    = ωm
    3ωc
    = ωm
    4ωc
    = ωm
    f (t) = sin (2πωc
    t + B sin (2πωm
    t))

    View Slide

  20. ഒԻ੒෼ͷݮਰ
    f (t) = Ec
    (t) sin (2πωc
    t + Em
    (t) B sin (2πωm
    t))
    ৼ෯ͷ࣌ؒมԽΛ
    ίϯτϩʔϧ͢Δ
    ഒԻ੒෼ͷྔͷ࣌ؒมԽΛ
    ίϯτϩʔϧ͢Δ
    άϥϯυϐΞϊ

    View Slide

  21. '.Իݯ͸
    ΞίʔεςΟοΫͳָثͷ
    ಛ௃Λগͳ͍ύϥϝʔλͰදݱग़དྷΔ
    ѱ͘ͳ͍ϞσϧͰ͋Δ

    View Slide

  22. OΦϖϨʔλͷ'.Իݯͷܭࢉࣜ
    OΦϖϨʔλ਺
    &JΤϯϕϩʔϓ
    TJجԻͷप೾਺ʹର͢Δഒ཰
    ЋJKΦϖϨʔλK͕ΦϖϨʔλJΛมௐ͢Δڧ͞
    ЌJग़ྗʹؚ·ΕΔΦϖϨʔλJͷ৴߸ͷׂ߹
    ୈ12ճΧʔωϧ/VM୳ݕୂ Ҩ఻తFMԻݯ(2015೥)ΑΓ
    '.Իݯͷ໰୊఺

    View Slide

  23. Ͳ͜Λ͍ͬͨ͡Β
    ཉָ͍͠ثͷԻ৭ʹۙ෇͘ͷ͔
    ͬ͞ͺΓΘ͔ΒΜ
    OΦϖϨʔλͷ'.Իݯͷܭࢉࣜ
    ୈ12ճΧʔωϧ/VM୳ݕୂ Ҩ఻తFMԻݯ(2015೥)ΑΓ
    '.Իݯͷ໰୊఺
    ૂͬͨԻ৭Λ
    ࡞Δͷ͕
    ͱͯ΋೉͍͠

    View Slide

  24. ୈੈ୅
    ୈੈ୅
    ୈੈ୅
    Ҩ఻తΞϧΰϦζϜ
    ͜ͷૢ࡞Λ܁Γฦ͢͜ͱʹΑͬͯ
    ༩͑ΒΕͨ৚݅ʹΑ͘ద߹͢ΔݸମΛ
    ݟ͚ͭग़͢͜ͱ͕Ͱ͖Δ
    ୈ12ճΧʔωϧ/VM୳ݕୂ Ҩ఻తFMԻݯ(2015೥)ΑΓ
    ͦ͜Ͱ
    Ҩ఻తΞϧΰϦζϜͰ
    '.Իݯύϥϝʔλͷ
    ୳ࡧΛࢼΈͨ

    View Slide

  25. Ҩ఻తFMԻݯ
    ར఺
    &OEUP&OEͰԻ͔Β'.Իݯύϥϝʔλ͕ಘΒΕΔ
    ࣮༻ʹ଱͑Δ඼࣭ͷग़ྗ͕ಘΒΕΔ
    ܽ఺
    ηοτͷύϥϝʔλΛݟ͚ͭΔͷʹ(16Ͱ͔͔࣌ؒۙ͘Δ
    ෳ਺ͷ୳ࡧ݁Ռ͕ઢܗิؒͰ͖Δอূ͕ͳ͍

    View Slide

  26. ୈ8ճΧʔωϧ/VM୳ݕୂ@ؔ੢ χϡʔϥϧFMԻݯ(2017೥)ΑΓ
    χϡʔϥϧωοτϫʔΫͰ
    '.Իݯύϥϝʔλͷ
    ୳ࡧΛࢼΈͨ
    Caffe
    Φʔϓϯιʔεͳ
    σΟʔϓϥʔχϯάϑϨʔϜϫʔΫͷ1ͭ
    http://caffe.berkeleyvision.org/

    View Slide

  27. χϡʔϥϧFMԻݯ
    ར఺
    &OEUP&OEͰԻ͔Β'.Իݯύϥϝʔλ͕ಘΒΕΔ
    ֶशࡁΈͷωοτϫʔΫ͕͋Ε͹਺ඵͰύϥϝʔλΛಘΒΕΔ
    ܽ఺
    ग़ྗͷ඼࣭͕࣮༻ʹ଱͑ͳ͍
    ෳ਺ͷ୳ࡧ݁Ռ͕ઢܗิؒͰ͖Δอূ͕ͳ͍

    View Slide

  28. ωc
    = ωm
    2ωc
    = ωm
    3ωc
    = ωm
    4ωc
    = ωm
    5ωc
    = ωm
    f (t) = A cos (2πωc
    t + B sin (2πωm
    t))

    View Slide

  29. ഒԻ੒෼ͷग़ํʹ͸نଇੑ͕͋ΔΑ͏ʹݟ͑Δ
    ωc
    = ωm
    2ωc
    = ωm
    3ωc
    = ωm
    4ωc
    = ωm
    5ωc
    = ωm
    f (t) = A cos (2πωc
    t + B sin (2πωm
    t))

    View Slide

  30. Ͳ͜Λ͍ͬͨ͡Β
    ཉָ͍͠ثͷԻ৭ʹۙ෇͘ͷ͔
    ͬ͞ͺΓΘ͔ΒΜ
    OΦϖϨʔλͷ'.Իݯͷܭࢉࣜ
    ຊ౰ʹ?
    ୈ12ճΧʔωϧ/VM୳ݕୂ Ҩ఻తFMԻݯ(2015೥)ΑΓ

    View Slide

  31. sin (x)
    2πωm
    t
    +
    sin (x)
    ×
    B
    ×
    A
    ग़ྗ
    2πωc
    ×
    ×
    2ΦϖϨʔλ௚ྻͷ৔߹
    Τϯϕϩʔϓ͸
    ඍখͳ࣌ؒͷൣғͰ͸ఆ਺ͱݟ၏ͤΔͨΊ
    ͱ ʹؚ·ΕΔ
    A B
    ͋Δ࣌ࠁ෇ۙͰͷFMԻݯͷग़ྗΛߟ͑Δ

    View Slide

  32. ΩϟϦΞΛ ʹஔ͖׵͑Δ
    cos
    cos (x)
    2πωm
    t
    +
    sin (x)
    ×
    B
    ×
    A
    ग़ྗ
    2πωc
    ×
    ×
    f (t) = A cos (2πωc
    t + B sin (2πωm
    t))
    Ґ૬͕ ͣΕΔ͕ৼ෯ʹ͸Өڹ͕ͳ͍
    π
    2

    View Slide

  33. f (t) = A cos (2πωc
    t + B sin (2πωm
    t))
    ࡾ֯ؔ਺ͷՃ๏ఆཧ
    f (t) = A cos (2πωc
    t) cos (B sin (2πωm
    t))
    −A sin (2πωc
    t) sin (B sin (2πωm
    t))
    cos (a + b) = cos a cos b − sin a sin b

    View Slide

  34. f (t) = A cos (2πωc
    t) cos (B sin (2πωm
    t))
    −A sin (2πωc
    t) sin (B sin (2πωm
    t))
    cos (z sin θ) = J0
    (z) + 2


    k=1
    J2k
    (z) cos (2kθ)
    sin (z sin θ) = 2


    k=0
    J2k+1
    (z) sin ((2k + 1) θ)
    Olver, Frank W. NIST handbook of mathematical functions. Cambridge New York:
    Cambridge University Press NIST, 2010. p.226

    View Slide

  35. cos (z sin θ) = J0
    (z) + 2


    k=1
    J2k
    (z) cos (2kθ)
    sin (z sin θ) = 2


    k=0
    J2k+1
    (z) sin ((2k + 1) θ)
    Olver, Frank W. NIST handbook of mathematical functions. Cambridge New York:
    Cambridge University Press NIST, 2010. p.226
    f (t) = A cos (2πωc
    t)
    (
    J0
    (B) + 2


    k=1
    J2k
    (B) cos (2k (2πωm
    t)))
    −A sin (2πωc
    t)
    (
    2


    k=0
    J2k+1
    (B) sin ((2k + 1) (2πωm
    t)))

    View Slide

  36. f (t) = AJ0
    (B) cos (2πωc
    t)
    +2A


    k=1
    J2k
    (B) cos (2k (2πωm
    t)) cos (2πωc
    t)
    −2A

    ∑ J2k+1
    (B) sin ((2k + 1) (2πωm
    t)) sin (2πωc
    t)
    ࣜมܗͯ͠
    f (t) = A cos (2πωc
    t)
    (
    J0
    (B) + 2


    k=1
    J2k
    (B) cos (2k (2πωm
    t)))
    −A sin (2πωc
    t)
    (
    2


    k=0
    J2k+1
    (B) sin ((2k + 1) (2πωm
    t)))

    View Slide

  37. cos a cos b =
    1
    2 (cos (a − b) + cos (a + b))
    sin a sin b =
    1
    2 (cos (a + b) − cos (a − b))
    ࡾ֯ؔ਺ͷੵ࿨ެࣜ
    f (t) = AJ0
    (B) cos (2πωc
    t)
    +2A


    k=1
    J2k
    (B) cos (2k (2πωm
    t)) cos (2πωc
    t)
    −2A


    k=0
    J2k+1
    (B) sin ((2k + 1) (2πωm
    t)) sin (2πωc
    t)

    View Slide

  38. ࡾ֯ؔ਺ͷੵ࿨ެࣜ
    f (t) = AJ0
    (B) cos (2πωc
    t)
    +A


    k=1
    J2k
    (B) cos ((2k (2πωm) − 2πωc) t)
    +A


    k=1
    J2k
    (B) cos ((2k (2πωm) + 2πωc) t)
    −A


    k=0
    J2k+1
    (B) cos (((2k + 1) (2πωm) + 2πωc) t)
    +A


    k=0
    J2k+1
    (B) cos (((2k + 1) (2πωm) − 2πωc) t)
    2 ( ( ) ( ))

    View Slide

  39. +A∑
    k=1
    J2k
    (B) cos ((2k (2πωm) − 2πωc) t)
    +A


    k=1
    J2k
    (B) cos ((2k (2πωm) + 2πωc) t)
    −A


    k=0
    J2k+1
    (B) cos (((2k + 1) (2πωm) + 2πωc) t)
    +A


    k=0
    J2k+1
    (B) cos (((2k + 1) (2πωm) − 2πωc) t)
    ࣜΛ੔ཧ͢Δͱ
    f (t) = A


    n=−∞
    Jn
    (B) cos (2πt (nωm
    + ωc))

    View Slide

  40. f (t) = A


    n=−∞
    Jn
    (B) cos (2πt (nωm
    + ωc))
    ͜ͷৼ෯Ͱ ͜ͷप೾਺ͷ೾Λ
    શͯ଍͠߹ΘͤΔͱ
    FMԻݯ͕ग़ྗ͢Δ೾ͱ౳Ձ
    प೾਺มௐͷల։

    View Slide

  41. ͱ͜ΖͰ

    View Slide

  42. f (t) = A


    n=−∞
    Jn
    (B) cos (2πt (nωm
    + ωc))
    ͬͯԿ?
    J

    View Slide

  43. ୈҰछϕοηϧؔ਺
    Jn
    (x) = 2
    π
    ∫ π
    2
    0
    sin (x sin ω) sin nωdω (n = 2k + 1|k ∈ ℕ)
    Jn
    (x) = 2
    π
    ∫ π
    2
    0
    cos (x sin ω) cos nωdω (n = 2k|k ∈ ℕ)
    ͕ح਺ͷ࣌
    n
    ͕ۮ਺ͷ࣌
    n
    ஫ҙ: ͜͜Ͱࣗવ਺ ͸ ΛؚΉ
    ℕ 0

    View Slide

  44. ϙΠϯτ
    ࣍਺͕ߴ͍΄Ͳ
    ্ཱ͕ͪΓ͕؇΍͔ʹͳΔ
    ϙΠϯτ
    ͕େ͖͘ͳΔͱ
    पظؔ਺ͱݟ၏ͤΔ
    x

    View Slide

  45. ͕ෛͷ৔߹
    n

    View Slide

  46. f (t) = A


    n=−∞
    Jn
    (B) cos (2πt (nωm
    + ωc))
    ͕͜͜ෛͷ࣌
    प೾਺ ͷ೾͸Ґ૬͕ ͣΕͯग़ͯ͘Δ
    2πt (nωm
    + ωc) π
    ͔͠͠ਓؒͷௌ֮͸Իͷप೾਺ຖͷৼ෯͚ͩΛௌ͍͍ͯͯ
    Ґ૬ͷҧ͍Λ஌֮͢Δ͜ͱ͸Ͱ͖ͳ͍

    View Slide

  47. ࣮ࡍʹ஌֮͞ΕΔԻͷৼ෯͸
    ϕοηϧؔ਺ͷઈର஋ʹͳΔ

    View Slide

  48. f (t) = A


    n=−∞
    Jn
    (B) cos (2πt (nωm
    + ωc))
    ͷ৔߹
    ෛͷप೾਺ͷ೾Λ࣋ͭ͜ͱʹͳΔ
    nωm
    + ωc
    < 0
    ͳͷͰ
    ͋Δෛͷप೾਺ ͷ೾ؚ͕·Ε͍ͯΔ࣌
    ͦΕ͸प೾਺ ͷ೾ͱͯ͠؍ଌ͞ΕΔ
    cos (x) = cos (−x)
    ω
    |ω|
    ෛͷप೾਺੒෼

    View Slide

  49. 2ωc
    = ωm
    3ωc
    = ωm
    0 1 2 3 4 5 6 7 8 9
    −1 −2 −3 −4 −5 −6 −7 −8 −9 −10
    0 1 2 3 4 5 6
    −1 −2 −3 −4 −5 −6
    ͭ·ΓFMԻݯͷεϖΫτϧͷҙຯ͸͜͏
    2ωc
    2ωc 3ωc
    3ωc

    View Slide

  50. ྔࢠԽͯ͠
    ઢܗิؒͯ͠
    ૭ؔ਺Λ͔͚ͯ
    FFTׂ͖ͯͨ͠ʹ͸
    Α͘Ұக͍ͯ͠Δ
    ͷ৔߹ͷ
    ల։݁ՌͱFFTͷൺֱ
    3ωc
    = ωm

    View Slide

  51. ༩͑ΒΕͨεϖΫτϧʹ͍ۙ݁Ռ͕ಘΒΕΔ
    ͱ ͱ ͱ ΛٻΊ͍ͨ
    A B ωc
    ωm
    f (t) = A


    n=−∞
    Jn
    (B) cos (2πt (nωm
    + ωc))

    View Slide

  52. ͸ ͷ஋͔͠ฦ͞ͳ͍ҝ
    ͢΂ͯͷ ͷ૯࿨Ͱׂ͓͚ͬͯ͹ ͸ܭଌͨ͠ৼ෯ͦͷ΋ͷʹͳΔ
    cos −1 ≤ x ≤ 1
    J A
    f (t) = A
    1
    ∑∞
    n=−∞
    |Jn
    (B)|


    n=−∞
    Jn
    (B) cos (2πt (nωm
    + ωc))
    A (0.5) = 0.1

    View Slide

  53. ָثͷԻͷৼ෯͸جԻ෇ۙͰ࠷΋େ͖͘
    ԕ͍ഒԻ΄Ͳখ͘͞ͳΔ܏޲͕͋Δ
    جԻͷৼ෯Λ ʹͯ͠
    ԕ͍ഒԻ΄Ͳ࣍਺ͷߴ͍ Λ࢖͏Α͏ʹ͢Δͷ͕ಘࡦ
    J0
    J

    View Slide

  54. ָثͷԻͷৼ෯͸جԻ෇ۙͰ࠷΋େ͖͘
    ԕ͍ഒԻ΄Ͳখ͘͞ͳΔ܏޲͕͋Δ
    جԻͷৼ෯Λ ʹͯ͠
    ԕ͍ഒԻ΄Ͳ࣍਺ͷߴ͍ Λ࢖͏Α͏ʹ͢Δͷ͕ಘࡦ
    J0
    J
    ͭ·Γ ͸جԻͷप೾਺
    ͸ ͷ੔਺ഒͷ஋ʹͳΔ
    ωc
    ωm
    ωc

    View Slide

  55. Λେ͖͗͘͢͠Δͱप೾਺ͷߴ͍ഒԻ͕ՄௌҬΛ௒͑Δ
    αϯϓϦϯάप೾਺ʹରͯ͠ߴ͗͢ΔഒԻ
    ͸جԻͷ੔਺ഒͰͳ͍ͱ͜ΖʹϐʔΫΛ࡞Δ
    ωm
    261.63Hz 20kHz
    20ωc
    = ωm
    10ωc
    = ωm
    30ωc
    = ωm
    ਓؒͷௌ֮Ͱ
    ฉ͖औΕΔݶք
    ͜͏ͳΔͱԻ֊Λָ࣋ͭثͱͯ͠੒ཱ͠ͳ͍ҝ ͱ ͷൺ͸ҎԼͰे෼
    ωm
    ωc

    View Slide

  56. ͱ ͸ఆ·ͬͨ
    ͸ߴʑ௨ΓͳͷͰ
    ֤ Ͱ༩͑ΒΕͨεϖΫτϧʹ࠷΋ۙ͘ͳΔ ΛٻΊΕ͹ྑ͍
    A ωc
    ωm
    ωm
    B
    χϡʔϥϧωοτϫʔΫͱಉ͡ඇઢܗ࠷దԽ໰୊ͳͷͰ
    ޯ഑๏Ͱ߈ΊΔ

    View Slide

  57. جԻͷ੔਺ഒͷ೾ͷৼ෯ͷΈʹண໨͢Δ
    3ωc
    = ωm
    άϥϯυϐΞϊ
    J0
    (B) J1
    (B) J1
    (B) J2
    (B) J2
    (B) J3
    (B) J3
    (B) J4
    (B) J4
    (B) J5
    (B) J5
    (B) J6
    (B)
    D10
    = J3
    (B)2 − E2
    10
    L10
    = D10
    D3
    = 0
    L3
    = E3
    ੜ੒෺ ظ଴͢Δग़ྗ
    Λޡࠩٯ఻೻͢Δ͜ͱͰ ͷमਖ਼ͷྔΛܾఆ͢Δ
    Dn
    B

    View Slide

  58. جԻͷ੔਺ഒͷ೾ͷৼ෯ͷΈʹண໨͢Δ
    3ωc
    = ωm
    άϥϯυϐΞϊ
    J0
    (B) J1
    (B) J1
    (B) J2
    (B) J2
    (B) J3
    (B) J3
    (B) J4
    (B) J4
    (B) J5
    (B) J5
    (B) J6
    (B)
    D10
    = J3
    (B)2 − E2
    10
    L10
    = D10
    D3
    = 0
    L3
    = E3
    ΋ޙͰ࢖͏ͷͰٻΊ͓ͯ͘
    50

    n=1
    Ln
    ʹΑΒͣੜ͡Δޡࠩ
    B

    View Slide

  59. Adam
    α = 0.001
    β1
    = 0.9
    β2
    = 0.999
    ϵ = 10−8
    mt
    = β1
    mt−1
    + (1 − β1) gt
    vt
    = β2
    vt−1
    + (1 − β2) g2
    t
    ̂
    mt
    =
    mt
    1 − βt
    1
    ̂
    vt
    =
    vt
    1 − βt
    2
    θt
    = θt−1
    − α
    ̂
    mt
    ̂
    vt
    + ϵ
    χϡʔϥϧωοτϫʔΫͷֶशͰ
    Α͘༻͍ΒΕΔ࠷దԽΞϧΰϦζϜ
    ࣮૷͕؆୯ͳׂʹͦͦ͜͜ੑೳ͕ྑ͍
    ʹ͍ͭͯޡࠩͷޯ഑ ʹରͯ͠
    Λ ͚ͩมԽͤ͞Δ
    θ g
    θ −α
    ̂
    mt
    ̂
    vt
    + ϵ

    View Slide

  60. Jn
    (x) = 2
    π
    ∫ π
    2
    0
    sin (x sin ω) sin nωdω (n = 2k + 1|k ∈ ℕ)
    Jn
    (x) = 2
    π
    ∫ π
    2
    0
    cos (x sin ω) cos nωdω (n = 2k|k ∈ ℕ)
    ޡࠩٯ఻೻Λޮ཰Α͘ߦ͏ʹ͸ର৅ͷؔ਺Λඍ෼͢Δඞཁ͕͋Δ
    ͍ͭ͜͸ඍ෼Ͱ͖Δͷ͔

    View Slide

  61. dJn
    (x)
    dx
    =
    1
    2
    (Jn−1
    (x) − Jn+1
    (x))
    Good News!
    ୈछϕοηϧؔ਺ʹ͸ղੳతʹඍ෼Ͱ͖Δࣄ͕஌ΒΕ͍ͯΔ

    View Slide

  62. dJn
    (x)
    dx
    =
    1
    2
    (Jn−1
    (x) − Jn+1
    (x))
    Bad News!
    ୈछϕοηϧؔ਺ͷඍ෼ʹ͸ୈछϕοηϧؔ਺ؚ͕·Ε͍ͯΔ

    View Slide

  63. Jn
    (x) = 2
    π
    ∫ π
    2
    0
    sin (x sin ω) sin nωdω (n = 2k + 1|k ∈ ℕ)
    Jn
    (x) = 2
    π
    ∫ π
    2
    0
    cos (x sin ω) cos nωdω (n = 2k|k ∈ ℕ)
    ͜ͷ෦෼ʹ໘౗ͳੵ෼ؚ͕·Ε͍ͯΔ
    ୆ܗެࣜΛ࢖ͬͯ਺஋తʹղ͘͜ͱ͸Ͱ͖Δ͕
    ΛճٻΊΔͷʹେྔͷܭࢉ͕ඞཁʹͳΔҝ
    ద੾ͳ Λ୳͢൓෮ܭࢉͷதͰ΍Γͨ͘ͳ͍
    Jn
    (x)
    B
    ͳͥ໰୊͔

    View Slide

  64. ϕοηϧؔ਺ͷۙࣅ
    Tumakov, D.N. Lobachevskii J Math (2019) 40: 1725. https://doi.org/10.1134/S1995080219100287
    The Faster Methods for Computing Bessel Functions of the First Kind of an Integer Order
    with Application to Graphic Processors
    p0
    (x) =
    8
    x
    q0
    (x) = p (x)2
    t0
    (x) = x − 0.78539816
    J0
    (x) = 0.079577472p0
    (x) ((0.99999692 − 0.0010731477q0
    (x)) cos (t0
    (x)) + (0.015624982 + (−0.00014270786 + 0.0000059374342q0) q0
    (x)) p0
    (x) sin (t0
    (x)))
    p1
    (x) =
    0.63661977
    x
    q1
    (x) = p (x)2
    J1
    (x) = p (1 − (0.46263771 − 1.1771851q1
    (x)) q1
    (x)) cos (x − 2.3561945 (0.46263771 − 1.1771851q1
    (x)) p1
    (x))
    Jn
    (x) =
    2 (n − 1)
    x
    Jn−1
    (x) − Jn−2
    (x)

    View Slide

  65. ϕοηϧؔ਺ͷۙࣅ
    Tumakov, D.N. Lobachevskii J Math (2019) 40: 1725. https://doi.org/10.1134/S1995080219100287
    The Faster Methods for Computing Bessel Functions of the First Kind of an Integer Order
    with Application to Graphic Processors
    Jn
    (x) =
    2 (n − 1)
    x
    Jn−1
    (x) − Jn−2
    (x)
    ͜ͷۙࣅͷDPPMͳͱ͜Ζ
    ͭԼͷ࣍਺
    ͱ ͭԼͷ࣍਺
    ͕ط஌ͷ৔߹
    ΋ͷ͘͢͝؆୯ͳܭࢉͰ ͕ٻ·Δ
    Jn−1
    (x) Jn−2
    (x)
    Jn
    (x)
    Լͷ࣍਺͔Βܭࢉ͢Δͱര଎Ͱϕοηϧؔ਺ͷ஋͕ٻ·Δ

    View Slide

  66. J0
    (x)

    View Slide

  67. J1
    (x)

    View Slide

  68. J2
    (x)

    View Slide

  69. J6
    (x)

    View Slide

  70. J10
    (x)
    ۙ๣ͰਧͬඈͿ
    x = 0

    View Slide

  71. p0
    (x) =
    8
    x
    q0
    (x) = p (x)2
    t0
    (x) = x − 0.78539816
    J0
    (x) = 0.079577472p0
    (x) ((0.99999692 − 0.0010731477q0
    (x)) cos (t0
    (x)) + (0.015624982 + (−0.00014270786 + 0.0000059374342q0) q0
    (x)) p0
    (x) sin (t0
    (x)))
    p1
    (x) =
    0.63661977
    x
    q1
    (x) = p (x)2
    J1
    (x) = p (1 − (0.46263771 − 1.1771851q1
    (x)) q1
    (x)) cos (x − 2.3561945 (0.46263771 − 1.1771851q1
    (x)) p1
    (x))
    Jn
    (x) =
    2 (n − 1)
    x
    Jn−1
    (x) − Jn−2
    (x)
    ͜͜
    ͜͜
    ͜͜
    ֻ͕͔͍ͬͯΔҝ ۙ๣ͰਧͬඈͿఆΊʹ͋Δ
    1
    x
    x = 0

    View Slide

  72. Jn
    (x) = 2(n − 1)
    x
    Jn−1
    (x) − Jn−2
    (x) (x ≥ 2 (n − 1))
    Jn
    (x) = Jn−1
    (x) (x < 2 (n − 1))
    ݩͷ࿦จͰఏҊ͞Ε͍ͯΔճආํ๏
    ͩͬͨΒ1ͭԼͷ࣍਺ͷϕοηϧؔ਺Λͦͷ··࢖͏
    x < 2 (n − 1)

    View Slide

  73. J0
    (x)

    View Slide

  74. J1
    (x)

    View Slide

  75. J2
    (x)

    View Slide

  76. J6
    (x)

    View Slide

  77. J10
    (x)
    ਧͬඈͼʹ͘͘͸ͳΔ
    Ͱ΋ਧͬඈͿ
    Ͱͷ஋͕߹Θͳ͍
    x < 2 (n − 1)

    View Slide

  78. ద੾ͳ ͷ୳ࡧ͸
    ͜ͷลΓ͔ΒղΛ୳͢ҝ
    ͜ͷ··Ͱ͸͜ͷۙࣅ͸࢖͑ͳ͍
    B

    View Slide

  79. ϕοηϧؔ਺ͷۙࣅ
    ͔Β ·Ͱͷ஋Λଟ߲ࣜͰิؒ͢Δ
    x = 0 x = n
    Jn
    (x) = 2(n − 1)
    x
    Jn−1
    (x) − Jn−2
    (x) (x ≥ n)
    Jn
    (x) = an
    x4 + bn
    x3 (x < n)
    an
    = −
    Jn
    (n)
    n5
    +
    1
    n4
    dJn
    dx
    (n)
    bn
    =
    Jn
    (n)
    2n3

    1
    n2
    dJn
    dx
    (n)
    ͨͩ͠

    View Slide

  80. ଟ߲ࣜͷ܎਺Λܾఆ͢ΔͨΊʹ
    ͱ ͕ཁΔ
    Jn
    (n)
    dJn
    dx
    (n)
    ͜ͷͭͷ஋Λ ͋ͨΓ·Ͱ
    ࣄલʹٻΊ͓ͯ͘
    n = 60
    an
    = −
    Jn
    (n)
    n5
    +
    1
    n4
    dJn
    dx
    (n)
    bn
    =
    Jn
    (n)
    2n3

    1
    n2
    dJn
    dx
    (n)
    n
    0 1.000000 0.000000
    1 0.440051 0.383960
    2 0.351942 0.175971
    3 0.308550 0.073205
    4 0.279146 0.032248
    5 0.255199 0.011946
    6 0.245656 0.002787
    7 0.233618 -0.004337
    8 0.221478 -0.008807
    9 0.214564 -0.010998
    10 0.207680 -0.012935
    Jn
    (n) dJn
    dx
    (n)

    View Slide

  81. J0
    (x) = 1 −
    x2
    4
    +
    x4
    64

    x6
    2304
    (x < 1)
    J1
    (x) =
    x
    2

    x3
    16
    +
    x5
    384

    x7
    18432
    +
    x9
    1474560
    (x < 4)
    ͱ Ͱ ͕খ͍͞৔߹
    ୯७ͳςʔϥʔల։ʹ੾Γସ͑Δ
    J0
    J1
    x
    ͜ΕΒΛ࢖ͬͯϕοηϧؔ਺ͷۙࣅ͕
    ਧͬඈͿ෦෼ʹ֖Λ͢Δ

    View Slide

  82. J0
    (x)

    View Slide

  83. J1
    (x)

    View Slide

  84. J2
    (x)

    View Slide

  85. J6
    (x)

    View Slide

  86. ࣮༻ʹ଱͑Δۙࣅ
    J10
    (x)

    View Slide

  87. const auto t0 = std::chrono::steady_clock::now();
    for( unsigned int b_ = 0.f; b_ != 400; ++b_ ) {
    float b = b * 0.1f;
    for( unsigned int n = 0; n != 60; ++n ) {
    ifm::bessel_kind1( b, n );
    }
    }
    const auto t1 = std::chrono::steady_clock::now();
    for( unsigned int b_ = 0.f; b_ != 400; ++b_ ) {
    float b = b * 0.1f;
    float l1 = 0;
    float l2 = 0;
    for( unsigned int n = 0; n != 60; ++n ) {
    float y = ifm::bessel_kind1_approx_2019( b, n, l1, l2, pre );
    l2 = l1;
    l1 = y;
    }
    }
    const auto t2 = std::chrono::steady_clock::now();
    float e0 = float( std::chrono::duration_cast< std::chrono::microseconds >( t1 - t0 ).count() );
    float e1 = float( std::chrono::duration_cast< std::chrono::microseconds >( t2 - t1 ).count() );
    $ ./bessel_benchmark -b bessel2.mp
    ୆ܗެࣜ: 0.0484659Mbps(24000 samples in 495194microseconds)
    ۙࣅࣜ: 269.663Mbps(24000 samples in 89microseconds)
    1570෼ׂͷ
    ୆ܗެࣜʹΑΔܭࢉͱൺֱͯ͠
    5500ഒߴ଎
    Mega Bessel Per Second
    Intel(R) Core(TM) i5-6600 3.3GHz
    Linux 5.3.7 + gcc 8.3.0
    -O2 -march=nativeͰ1ίΞͷΈ࢖༻

    View Slide

  88. ϐΞϊͷ৔߹

    View Slide

  89. ϐΞϊͷ৔߹

    View Slide

  90. ϚϦϯόͷ৔߹

    View Slide

  91. ϚϦϯόͷ৔߹

    View Slide

  92. ϚϦϯόͷ৔߹

    View Slide

  93. ࣌ؒมԽ
    Τϯϕϩʔϓ͸؇΍͔ʹมԽ͢Δҝ
    ࣌ࠁ ʹ͍ͭͯٻΊͨ ͸
    ࣌ࠁ ͷ ΛٻΊΔࡍͷॳظ஋ͱͯ͠࢖͑Δ
    t B
    t + 1 B
    t
    t − 1
    t − 2
    t − 3 t + 1 t + 2 t + 3
    Bt
    Bt
    Bt+1
    Bt+2
    Bt−1
    Bt−2

    View Slide

  94. ࣌ؒมԽ
    Τϯϕϩʔϓ͸؇΍͔ʹมԽ͢Δҝ
    ࣌ࠁ ʹ͍ͭͯٻΊͨ ͸
    ࣌ࠁ ͷ ΛٻΊΔࡍͷॳظ஋ͱͯ͠࢖͑Δ
    t B
    t + 1 B
    t
    t − 1
    t − 2
    t − 3 t + 1 t + 2 t + 3
    Bt
    Bt
    Bt+1
    Bt+2
    Bt−1
    Bt−2
    ҰൠతͳFMԻݯ͸ Λ࣌ؒͰมԽͤ͞ΒΕͳ͍ͷͰ
    ࣌ࠁ Ͱ࠷΋௿͍ Λग़ͨ͠ Λͦͷ··Ҿ͖ܧ͙
    ωm
    t L ωm
    ωm
    ωm
    ωm
    ωm
    ωm
    ωm

    View Slide

  95. View Slide

  96. ΞίʔεςΟοΫͳָثͷৼ෯ͷมԽ͸
    छྨʹ෼ྨͰ͖Δ
    ܧଓతʹྗ͕Ճ͑ΒΕΔ
    ඇݮਰৼಈ

    ࠷ॳ͚ͩྗ͕Ճ͑ΒΕΔ
    ݮਰৼಈ

    View Slide

  97. Τϯϕϩʔϓਪఆ
    f (t) = (1 − β) e−αt + β
    ݮਰৼಈ߲ ඇݮਰৼಈ߲
    ೖྗ೾ܗ͔ΒٻΊͨΤϯϕϩʔϓ ͱ
    ͱ
    ͜ͷؔ਺ͷग़ྗͷ͕ࠩ࠷খʹͳΔΑ͏ʹ
    ͱ ΛٻΊΔ
    A (t) B (t)
    α β

    View Slide

  98. f (t) = (1 − β) e−αt + β
    ferror
    (t) = (fexpected
    (t) − f (t))
    2
    dferror
    (t)

    = 2(1 − β)te−αt
    ((1 − β)(−e−αt) − β + fexpected
    (t))
    dferror
    (t)

    = 2 (e−αt − 1) ((1 − β)(−e−αt) − β + fexpected
    (t))
    ͜Ε΋൓෮ͯ͠ ͱ Λमਖ਼͍͚ͯ͠͹ղ͚Δ
    α β

    View Slide

  99. ղ͚ͨ
    FFT͔ΒಘͨΤϯϕϩʔϓ͸
    ࣌ʑ૯࿨͕1Λ௒͑ΔͷͰ࠷େ஋͕1ʹͳΔΑ͏ʹ
    εέʔϧ͓ͯ͘͠

    View Slide

  100. Τϯϕϩʔϓਪఆ
    σΟέΠதؒϨϕϧ
    αεςΠϯϨϕϧ
    ΞλοΫதؒϨϕϧ
    ΞλοΫ1
    ΞλοΫ2
    σΟέΠ1
    σΟέΠ2
    αεςΠϯ
    ϦϦʔε
    Α͋͘Δ'.ԻݯͷΤϯϕϩʔϓ͸
    Λஈ֊ͷઢܗͳؔ਺Ͱۙࣅ͢Δ
    e−αt

    View Slide

  101. FYQͱEFDBZʹғ·Εͨ໘ੵͱ
    FYQͱEFBZʹғ·Εͨ໘ੵͱ
    FYQͱTVTUBJOʹғ·Εͨ໘ੵͷ
    ࿨͕࠷খͱͳΔ఺ ͱ ΛٻΊΔ
    n m
    n
    m

    View Slide

  102. Τϯϕϩʔϓਪఆ
    flinear
    (t, n, m) =
    (f (m) − f (n))
    m − n
    t +
    nf (m) − mf (n)
    n − m
    fconst
    (t, n, m) =
    nf (m) − mf (n)
    n − m
    ferror (n, m, l) =

    n
    0
    flinear
    (t,0,n) − f (t) dt
    +

    m
    n
    flinear
    (t, n, m) − f (t) dt
    +

    l
    m
    fconst (t, m, l) − f (t) dt
    ͱ ʹ͍ͭͯͷޯ഑ΛٻΊ
    ஋Λमਖ਼͍͚ͯ͠͹ྑ͍
    n m
    ݮਰৼಈ͢Δָثͷ৔߹
    Λ ʹ Λαϯϓϧ௕ݻఆ͢Δ
    fconst
    0 l

    View Slide

  103. ୈ12ճΧʔωϧ/VM୳ݕୂ Ҩ఻తFMԻݯ(2015೥)ΑΓ
    ຊ෺ͷϐΞϊͷԻ
    '.ԻݯͷϐΞϊͷԻ
    ܭࢉ࣌ؒ
    2࣌ؒ43෼42ඵ/keyframe
    ܭࢉ؀ڥ

    Intel Core i5-6600(Skylake) 3.3GHz 4ίΞ
    nvidia GeForceGTX1070 1.68GHz 1920CUDAίΞ
    FMԻݯͷ࢓༷
    4ΦϖϨʔλϑΟʔυόοΫ෇͖
    ιϑτ΢ΣΞFMԻݯ

    View Slide

  104. ࣮ࡍͷϐΞϊͷԻ χϡʔϥϧFMԻݯ Ҩ఻తFMԻݯ
    ϐΞϊͬΆ͘ฉ͑͜Δ෦෼΋͋Δ͕
    Ҩ఻తFMԻݯͱൺ΂Δͱ͠ΐͬͺ͍
    ϐΞϊ
    ݁Ռ
    ୈ8ճΧʔωϧ/VM୳ݕୂ@ؔ੢ χϡʔϥϧFMԻݯ(2017೥)ΑΓ
    ܭࢉ࣌ؒ
    3ඵ/keyframe
    ܭࢉ؀ڥ

    Intel Core i5-6600(Skylake) 3.3GHz 4ίΞ
    nvidia GeForceGTX1070 1.68GHz 1920CUDAίΞ
    FMԻݯͷ࢓༷
    8ΦϖϨʔλϑΟʔυόοΫ෇͖
    ιϑτ΢ΣΞFMԻݯ

    View Slide

  105. ຊ෺ͷϐΞϊ
    2ΦϖϨʔλFMԻݯϐΞϊ
    ܭࢉ࣌ؒ
    4෼1ඵ/keyframe
    ܭࢉ؀ڥ

    Intel Core i5-6600(Skylake) 3.3GHz 4ίΞ
    CPUͷΈ࢖༻
    FMԻݯͷ࢓༷
    2ΦϖϨʔλϑΟʔυόοΫͳ͠
    ιϑτ΢ΣΞFMԻݯ

    View Slide

  106. ຊ෺ͷϚϦϯό
    2ΦϖϨʔλFMԻݯϚϦϯό

    View Slide

  107. ຊ෺ͷόΠΦϦϯ
    2ΦϖϨʔλFMԻݯόΠΦϦϯ

    View Slide

  108. 3ΦϖϨʔλͷ৔߹
    sin (x)
    cos (x)
    2πω1
    t
    +
    ×
    E1
    ×
    E0
    ग़ྗ
    2πω0
    ×
    ×
    2πω2
    E2
    × sin (x)
    ×
    +
    f (t) = E0
    cos (2πω0
    t + E1
    sin (2πω1
    t) + E2
    sin (2πω2
    t))
    f (t) = E0


    n1
    =−∞


    n2
    =−∞
    Jn1
    (E1) Jn2
    (E2) cos (2πt (n1
    ω1
    + n2
    ω2
    + ω0))
    ҰԠղ͚Δ
    ΋ͬͱෳࡶͳFMԻݯʹ͍ͭͯղ͚ͳ͍ͷ?

    View Slide

  109. f (t) = E0


    n1
    =−∞


    n2
    =−∞
    Jn1
    (E1) Jn2
    (E2) cos (2πt (n1
    ω1
    + n2
    ω2
    + ω0))
    ͜Εͱ ͜Εͱ ͜Εͱ ͜ΕΛ
    ݟ͚ͭΔඞཁ͕͋Δ

    View Slide

  110. ͜ͷํ๏ͩͱඞཁͳࢼߦճ਺͕(op਺-1)৐Ͱ૿͑Δ

    View Slide

  111. ͜ͷ͋ͨΓ͕ͬͦ͝Γফ͍͑ͯΔͷ͸ԿͰ

    View Slide

  112. ഒԻ͕ظ଴͢Δप೾਺͔Β͔ᷮʹͣΕͯग़͍ͯΔ
    ޡࠩ͸ظ଴͢Δप೾਺ͷৼ෯͚ͩΛ࢖ͬͯٻΊ͍ͯΔͨΊ
    ͣΕͯग़͍ͯΔഒԻ͸ग़͍ͯͳ͍ͷͱಉ͡ѻ͍ʹͳΔ
    ΞίʔεςΟοΫͳָثͳΒ
    ͲΜͳʹ͖ͪΜͱௐ཯ͯ͋ͬͯ͠΋
    ߴप೾ͷഒԻͰ͸͍͘Β͔ͣΕΔ

    View Slide

  113. ·ͱΊ
    FMԻݯ͸ղ͚Δ
    ඇઢܗ࠷దԽ໰୊ͳͷͰ
    χϡʔϥϧωοτϫʔΫͷֶशͱಉ͡ख๏͕ޮ͘
    ͔͠͠3ΦϖϨʔλҎ্ʹରԠ͢Δʹ͸
    ΋͏গ͠޻෉͕ඞཁ

    View Slide