オーディオ考

オーディオに求めるものなど。

2025年10月 1日 (水)

SDRドングルでAM受信

Overview_20250930153901
写真1.ブレッドボードに組んだアップコンバータとドングルSDR


ワンセグチューナのUSBドングルを使って、パソコンを広帯域受信機(SDR)として使用する方法はよく知られていて、おおよそ30MHz~1.6GHzの範囲をオールモードでカバーすることができる。

構成としては写真2に示すようなAmazonでも売っているUSBのワンセグチューナーと、フリーのラジオ受信用ソフトSDR#、それにデバドラの置き換えユーティリティZadigを使って簡単に導入することができる(たとえばこのサイトで説明されている)。


Dsdt308sv
写真2.Amazonで1200円ほどの安価なSDRドングル  DS-DT308SV

ドングルに付属のロッドアンテナだけでは少々受信感度が低いので、感度よく受信するには大きめのアンテナをつなぐ必要がある。FM放送を受信する程度であれば、ロッドアンテナに1~2mほどのビニル線をつなぐだけでもわりと受信できる。

ところで、ちょっと残念なのはこのドングルは受信周波数の下限が24MHz程度なので、それ以下のAM放送や短波帯は受信できない。そこでよく紹介されているのは、ドングルを改造してトランスを介してQ信号をダイレクトに入力して、SDR#の受信設定をDirect Samplingにする方法と、もうひとつはアップコンバータを使用してターゲット周波数を高い周波数に変換して受信する方法だ。

今回は、後者のアップコンバートによる方法で、極力入手性のよい部品(っていうか部品箱に入ってた部品)でAMラジオの受信を試してみた。この方法ではドングルを改造しなくて済むので、工作難易度は低い。
アップコンバートは、任意の周波数を足して、ターゲットの周波数を高い周波数に変換する方法で、たとえば954kHz(TBS放送)に100MHzの発振周波数を足して100.954MHzとして受信する。
アップコンバータを構成する主な要素は、加算周波数の信号を生成する発振器と、加算周波数信号を足し込むDBM(ダブルバランスドミキサー)で、発振器はESP32からarduinoでSi5351(PLLオシレータ)を使用して100MHzを生成し、DBMは秋月で150円で入手できるNJM2594と、手持ちのフェライトコアとショットキーバリアダイオードを使って手作りしたものの2通りを試してみた。

100MHzを生成する部分はESP32+Si5351を使用したが、これは単に100MHzのオシレータを使ってもよい。
Si5351のモジュールは秋月でもAmazonでもaliでも売っているので、容易に入手可能で、サンプルプログラムも豊富にある。
今回は、ESP32をArduinoで使用し、こちらサンプルプログラムを使って100MHzを生成した。この中で、次のようにclk0を100MHzに設定した。

si5351.set_freq(10000000000ULL, SI5351_CLK0);

DBMにNJM2594を使用した回路は図1の通り。

Njm2594sch
図1.DBMにNJM2594を使用した回路


電源は、PCにUSB接続されたESP32のブレークアウトボードから5Vと3.3Vを供給している。
アンテナからの入力信号は、広帯域オペアンプNJU77701Fで高周波増幅してからDBMに入力している。
Si5351のみ電源が3.3V、NJM2594およびNJU77701の電源は5Vなので注意してほしい。
USBドングルチューナーのロッドアンテナはmcxコネクタになっているので、ここに回路の出力(J1)を接続する。
アンテナには、室内で2mほどのビニル線を接続した。

SDR#で受信している様子を図2に示す。

20250930_2594tbs_l
図2.DBMにNJM2594を使用したコンバータでの受信

放送内容が聴き取れたのは954kHzのTBSラジオのみだったが、いちおう目的は達成された。
実際にやってみると、Si5351による100MHzの発振周波数が実際には100.020MHzになっているため、TBSラジオの受信周波数は100.974MHzとなっている。録音は次の通り。

20250929TBS_NJM2594

ほかの局については、AFNはなにか放送されていることはわかるが、内容はまったく聴き取れない程度、NHK1,NHK2は僅かにSDR#のスペクトル表示がされているように見えるものの、音としてはまったく聞こえなかった。

次に、DBMにNJM2594のかわりに、トロイダルトランス2個とショットキーダイオードを使ったものを試してみた。回路を図3に示す。

Trans_dbmsch_20250930164701
図3.トロイダルトランスとショットキーダイオードのDBMを使用した回路

この回路はNJM2594の部分をトランスとダイオードを使ったDBMに置き換えている。
トロイダルコアは以前秋月で購入したTR-10-5-5ED使用したが、現在は販売終了となっているため、代替品はFT-37-43がよいと思われるが、入手困難な場合は、秋月で入手できるFT-50-43またはFT-50-61などでも使用可能と思われる。
今回はついでにAliで入手した正体不明の8x4xt3のトロイダルコアも試してみた。
ダイオードは部品箱に入っていたND412Gというショットキーを使用したが、これも秋月で入手可能な高周波用の1SS106などが代替可能だと思われる。
100MHz入力側のトランスT1、Mix出力側のトランスT2は、線を3本束ねて巻くトリファイラ巻きとして、回路図のように結線する。
巻線にはΦ0.2のホルマル線を使用し、巻き数は次の通り。


表1.トランスT1、T2の巻き数(トリファイラ巻き)
Trans_data_20250930164701

Aliの正体不明のトロイダルコアはμが低いかもしれないので、巻き数を若干増やし巻いた。根拠はまったくない……
トロイダルトランスを使用したDBMを写真3,写真4に示す。
また、コア単体の写真を写真5に示す。大きい方がTR-10-5-5ED。


Tr10trans
写真3.トロイダルコアTR-10-5-5EDを使用したDBM


Ali_trans_20250930164701
写真4.Aliで購入した正体不明のトロイダルコアを使用したDBM


Cores_20250930164701
写真5.トロイダルコア 大きい方がTR-10-5-5ED


それではまずはTR-10-5-5EDを使用したDBMによる受信の様子を図4に示す。

20250930_sdr2l
図4.TR-10-5-5EDを使用したDBMによる受信

NJM2594に比べて、受信感度がかなり上がっている。ノイズフロアも上がっているが、信号レベルがかなり上がっているため、聞こえる局が増えた。実際に聞いた感じでは、TBS(954kHz)ははっきり聞き取れるレベル、AFN(810kHz)も音楽や言葉が聞き取れる(英語は聞き取れない)。NHK第2(693kHz)は、何か受信できているのはわかるが聞き取れない。NHK第1(594kHz)、文化(1134kHz)、ニッポン(1242kHz)は画面上スペクトルが出ているが、聞こえなかった。
TBSの録音は次の通り。

20250930TBS_Troidal_DBM


次に、Aliで購入した正体不明のトロイダルコアを使用したもの。図5に受信の様子を示す。


20250930_sdr_ali2l
図5.Aliで購入したトロイダルコアを使用したDBMによる受信

これはTR-10-5-5EDとほぼ同等か、若干高感度という意外な結果となった(@_@)
TBSの録音はつぎのとおり。

20250930TBS_Troidal_Ali_DBM

以上の通り、TBSラジオは受信レベルが高く、完全に聞き取れる受信状態となった。ただ、ブーンというノイズが大きく、これはESP32に接続しているPCのUSB電源か、ESP32モジュールが由来の可能性が高い。電源に電池などを使い、発振器も単体のものを使用すれば改善する可能性がある。
アンテナは室内のビニル線2mほどという条件で、この程度の受信ができたので、満足な結果だったと思う。


今回の検証では上記のほかに、バーアンテナとバリコンによる同調回路、およびLCによるLPFも検討したが、いずれも感度が低下してしまい、期待した効果が得られなかった。
ほとんど秋月で入手可能な部品で構成できたので、中波受信に挑戦してみたいという方におすすめします。AM放送はあと3年ほどでなくなってしまうので、記念に実験してみてはいかがでしょうか?



| | | コメント (0)

2025年7月24日 (木)

カオスふたたび

Chaos
図1.夏によく合う涼しげなカオス


カオス現象について取り上げたのは2022年の記事で、あれからもう3年経ってしまった。
事の発端は、当時タモリ倶楽部で取り上げられた早稲田大学の「役に立たない機械コンテスト」に出品された、二重振り子を応用したメトロノームで、リズムがしっちゃかめっちゃかで実に役に立たなそうだったが、二重振り子はカオス現象を紹介するためによく使われる。

2022年当時、私はこのカオス現象を電子回路で再現させて、カオスの音を聴いてみたいと思ってあれこれ実験した。2つの発振回路が、相互に周波数を変更するような作用を行うように構成すれば再現できるのではないか……?そう思い何通りか実験してみたがどうにもういまく行かなかった。たいていは発振状態が安定してしまい、カオスと呼べる状態にはならなかったのだった。
そこでネットで、電子回路でカオス現象を再現する方法がないか調べてみたところ、chua回路なるものを発見した。これは1983年にchuaさんという技術者によって発明された回路だそうだ。さっそくLTSPICEでシミュレーションし、実際に回路も組み立てて検証したところ、カオス波形が得られた。カオス波形は、X信号とY信号のリサージュ波形として得られるものだった。

このカオス波形の音を聴いてみたいと思ったが、当時はXYのリサージュをうまいこと音声信号に変換する方法を思いつかず、苦肉の策でX+Y信号、X-Y信号の2種類を生成して聴いてみたが、ホワイトノイズのようにしかならず、結果として面白みに欠けるものだった。リサージュということなら位相差を演算すればいいのではないか?というところまでは考えたのだが、結局そこまではやらなかった。これらはすべて2022年の記事に記している。

ところで先日、トランジスタ技術の最新号を読んでいたら、アナログ乗算器を使った位相差演算回路の製作記事が載っていた。これを読んでいて、カオス信号の音声化がうまくいっていなかったことを思い出した。この位相差演算回路を使えば、カオス信号をもう少しそれらしく音声信号化できるのではないか……?
しかし記事ではアナログ乗算器というICを使っていて、経験のない分野なのでたとえばこれがLTSPICEでシミュレーションできるのか?など、わからない点も多く、またしても棚上げにしそうになったが、よく考えてみれば回路で実現するよりも、すべてプログラムでやってみた方が早く、しかもいろんな検証ができるのではないかと思い、pythonを使って検証してみることにした。

まずはpythonでchua回路のシミュレーションを考える。chua回路は次の微分方程式で与えられる。

Chua_function
Chua_constants_20250725070901

ただし、α、βは回路部品(抵抗、コンデンサ、インダクタ)に依存する定数、m0,m1は非線形素子の特性を決める定数。

pythonを使って、この微分方程式を数値的に解いてプロットしたグラフが、冒頭図1のカオスのグラフだ。実際にpythonで、この微分方程式を定義して、それを解く部分は次のとおり。(プログラムの全体は記事の下の方でリンクします。)


###Pythonコード抜粋

# Chua回路の微分方程式
def chua_circuit(t, state, alpha, beta, m0, m1):
x, y, z = state
# 非線形関数 h(x)
h = m1 * x + 0.5 * (m0 - m1) * (abs(x + 1) - abs(x - 1))
dxdt = alpha * (y - x - h)
dydt = x - y + z
dzdt = -beta * y
return [dxdt, dydt, dzdt]

# パラメータ設定
alpha = 15.6
beta = 28.0
m0 = -1.143
m1 = -0.714

# 初期条件と時間範囲
initial_state = [0.1, 0.0, 0.0]
t_span = (0, 100)
t_eval = np.linspace(t_span[0], t_span[1], 10000)

# 微分方程式を数値的に解く
solution = solve_ivp(chua_circuit, t_span, initial_state, t_eval=t_eval, args=(alpha, beta, m0, m1))

###ここまで

このコードでは、x,y,z空間上で、10000ポイント分の数値演算ができ、その3変数のうちのxとzをプロットしたものが図1だ。これは典型的なカオスのグラフとなっている。
実際のコードは次のとおり。

ダウンロード - 20250721chua_000.py

ダウンロード - variable_sine_tone_builder_wav.py

上の20250721chua_000.pyがchua回路のシミュレーションプログラム本体で、下の
variable_sine_tone_builder_wav.pyは、任意の正弦波を発生させてwavファイルに記録する関数を含むファイルだ。

使い方は、上の2つのファイルを同じホルダに配置し、20250721chua_000.pyを実行する。
すると、図1に示すカオスのグラフが表示され、同ホルダ内にカオスのグラフを音声化したものが生成される。
まずは理屈よりも先に、生成されたカオスの音を聴いていただこう(音量注意)。
20250805追記
91行目を次のように変更して、再生速度を2倍にしたパターンも録音してみた。こちらの方がスピーディでいいかもしれない。
ab.add_tone(freq=vFreq,level_db=-6,pan=vPan,dur_s=0.01,fade_ms=5)

chuaWav000.wav

20250805chuaWav2x.wav  (2倍速Ver.)

なんとなくそれらしいサウンドになったように思うが、どうだろうか?ヨガのBGM用に売れるかもしれない。
これは、図1のカオスグラフのxの変域に対して、約65Hz~3kHzの範囲の正弦波を割り付けている。
x<0の領域では、440/3|x|、x>=0の領域では440*3xとしている(20250721chua_000.pyの85、87行目)。
これは単にxの変域に周波数を割り付けているだけで、当初想定した位相差演算などはしていない。

それでは位相差を演算してそれを音に変換するとどうだろうか。

chuaWavHilbert000.wav

20250805chuaWav_hilbert2x.wav  (2倍速Ver.)


思ったほどおもしろい効果は得られず、個人的には上の単純割り付けの結果の方が良いように思う。
これはカオスの演算結果x,yに対してそれぞれをヒルベルト変換し、ポイント毎の瞬時位相を求めてその差をとることで、位相差を演算して、その結果を約70Hz~2.7kHzに割り付けている(20250721chua_000.pyの80,82行目)。

また上記のいずれに対しても、共通して、演算結果zの値に応じて、-1(左)~1(右)のパンポットを割り当てている。

以下、コードの簡単な説明。

まずはコード本体の20250721chua_000.pyにおいて、45行目で微分方程式を数値的に解いて、47行目で演算結果x,y,zを各10000個のリストデータとして得ている。61行目以降で、これらのデータのうちx,zをプロットして、図1のグラフを得ている。

次にデータの音声化は、ヒルベルト変換で位相差を使うか、単に変域のデータを使うかを72行目のhilbertフラグで設定する。hilbert=1ならヒルベルト変換による位相差、hilbert=0なら変域データ。

任意の正弦波を発生する関数はvariable_sine_tone_builder_wav.pyに含まれているので、20250721chua_000.pyの13行目以降でAudioBuiderライブラリをインポートしている。
使い方は、
ab = AudioBuilder(sample_rate=48000)
(74行目)として、AudioBuilderオブジェクトをabとして、48kHzサンプリングのインスタンス化をして、ab.add_tone(freq=vFreq,level_db=-6,pan=vPan,dur_s=0.02,fade_ms=10)
(91行目)で各ポイントデータに対してfreq(周波数)、level_db(信号レベル)、pan(パンポット-1~1)、s(信号長さ(s))、フェードイン/アウト時間(ms)をそれぞれ設定して、ポイントデータ毎の正弦波を追加していく。
すべてのデータ(10000個)に対して信号の割り付けができたら、
ab.write_wav("chuaWav000.wav", normalize="store_true", peak=0.99)
(95行目)で、信号ピークを0.99に正規化したオーディオ信号を"chuaWav000.wav"として書き出す。

以上の処理で、グラフと音声信号を生成している。


2022年の記事と比べると、格段にカオスっぽいサウンドが得られたと自負しているが、いかがでしょうか?

 

| | | コメント (0)

2025年7月 6日 (日)

FET入力無帰還A級パワーアンプ(基板頒布あり)

Pcb
写真1.新たに開発したFET入力無帰還パワーアンプ基板


無帰還A級25Wパワーアンプはこのブログのベストセラーで、多くの方に使っていただいています。ありがとうございます(^-^)

先月、無帰還ヘッドホンアンプを2機種発表しましたが、試聴会を実施したところ、FET入力とバイポーラ入力で好みが分かれることがわかりました。無帰還A級25Wパワーアンプはバイポーラのダイヤモンドバッファ入力としていますが、これをFET入力とした場合に音の傾向が変わり、聞く人の好みによって選択肢を提供できるのではないかと思い、今回新たに、FET入力タイプの無帰還A級アンプを開発しました。

今回は、入力回路をバイポーラのダイヤモンドバッファからJFETのゼロバイアスバッファに変更し、それに伴い、入力デバイスのバイアス電流のバランス調整が不要になったため、調整が簡素化でき、また従来回路では出力オフセット調整がゲインに影響を与えるため、部品のペア取りがシビアでしたが、今回はオフセット調整がゲインに影響しないため、再現性がより高くなりました。

今回開発したアンプ基板の写真を写真1に、回路図を図1に示します。

Pav100_schematic
図1.FET入力無帰還パワーアンプ回路図


今回は入力回路の変更と、それに伴う定数変更を行いましたが、保護回路は従来通りで、基板の寸法も同一なため、すでに従来の無帰還A級25Wパワーアンプをお使いの方は、放熱板に取り付けたパワートランジスタと熱結合デバイスはそのまま、基板のみ交換することで試していただくことができます。

今回当方で試作した結果は以下のとおりです。
シャーシおよび電源は、従来と同じものを使用しました。組み立てたアンプを写真2に示します。


Overview_20250706102601
写真2.FET入力無帰還A級パワーアンプ


写真2の状態で測定した諸特性は以下のとおりです。

図2、図3にひずみ率特性を示します。THD+Nはボトムで0.03%程度、全体でも10Wでおおむね0.2%以下でした。

Distl   Distr_20250706114501
図2.ひずみ率特性(Lch)                 図3.ひずみ率特性(Rch)



周波数特性は左右差がなかったので、代表して左チャンネルの特性を図4に示します。

Freql
図4.ゲイン周波数特性 -3dB特性はDC~約250kHz


出力インピーダンスの測定結果を表1に示します。測定は3.1Ω負荷によるON/OFF法で行いました。これも左右で差がありませんでした。

表1.出力インピーダンス実測値
Imp



プロテクトの動作電圧を表2に示します。これは従来の無帰還A級25Wパワーアンプとほぼ同等の結果です。

表2.プロテクト動作電圧
Prtct



以上により、全体の特性は従来の無帰還A級25Wパワーアンプとほぼ同等となりました。

この1週間ほど、実際に音楽を聴いていますが、ほんの僅かながら、従来アンプに比べて音質が柔らかく繊細になった印象があります。これは好みによって評価が分かれるところだと思いますが、製作者としては今回のアンプの音質はとても気に入っています。

【基板頒布情報】
生基板2枚ひと組(ステレオ分)と、製作説明書と回路図、部品表等の資料、製作例、LTSPICEのシミュレーションファイルをセットで2600円(税・送料込み)で頒布します。自力で部品収集、部品選別、調整できる方が対象です。

ご希望の方は表題に「FET入力無帰還電圧アンプ基板頒布」、本文にお名前、送付先郵便番号、ご住所、電話番号をお書きのうえ、

dj_higo_officialアットhigon.sakura.ne.jp
(アットを@に替えてお送りください)

までメールをお送りください。

代金の振込先のご案内メールをお送りします。入金が確認でき次第発送します。

回路図、部品表、製作マニュアル等

LTSPICEシミュレーションファイル



| | | コメント (0)

2025年5月12日 (月)

無帰還CR型フォノイコライザの製作(基板頒布あり)

20250512pcb2
写真1.今回製作した無帰還CR型RIAAイコライザ基板(片ch)


【20250622追記】残留ノイズ測定結果
入力ショート時の残留ノイズと、MCモード(ゲイン46dB@1kHz)時の入力換算圧音は次のとおりだった。
これはかなり満足な結果だと思う。

表A.残留ノイズおよび入力換算雑音
Zanryunoise


【これより本文】

以前から、いつか作ってみたいと考えていた無帰還CRフォノイコライザをとうとう作ったので報告する。

全体構成は、CR型RIAA回路の前段と後段にそれぞれ無帰還のフラットアンプを配置する一般的なもので、カートリッジはDL-103を使うこと前提とし、1kHzでのゲインを46dBとした。また、MMカートリッジも使用できるようにジャンパピンによって1kHzでのゲインを20dB下げて26dBとすることができ、入力抵抗(カートリッジに対する負荷)もジャンパピンで47kΩと220Ωのいずれかを選べるようにした。
回路図(片ch分)は次のPDFをダウンロードしてご覧ください。
ダウンロード - 20250512nnfb_eq.pdf

RIAA回路の前段と後段のフラットアンプは同じ回路で、無帰還A級25Wパワーアンプ無帰還電流駆動ヘッドホンアンプで実績のある回路方式を踏襲している。今回は、繊細なカートリッジの信号を受けるには経験的にバイポーラよりもFETの方が適していると考え、入力段をゼロバイアスFET入力とした。またJFETのコンプリメンタリペアは入手困難であるため、Nチャネルのみの構成とした。
入力回路で受けた信号はプッシュプルカレントミラーで電流増幅され、抵抗(前段回路ではR17)によってI/V変換されてから終段のプッシュプルバッファを経て出力される。

フラットアンプ単体の諸特性は以下のとおり。ゲインの周波数特性と10kHz方形波の応答は左右でほとんど差がないため、代表としてLchのみを示す。

20250512flatamp_fspec
図1.フラットアンプのゲイン周波数特性


20250512houeiha
図2.フラットアンプの10kHz方形波応答(上段:方形波入力、下段:2Vp-p出力)

図1に示すように、全帯域でフラットであり、ゲインは約35dB、-3dB周波数はおよそ200kHz。
図2の方形波の応答では、波形の歪みは見られない。

次に、Lch、Rchそれぞれのひずみ率について、100Hz、1kHz、10kHzにおいて、出力を約2Vp-pとした場合のひずみ率を表1に示す。

表1.周波数別ひずみ率
20250512hizumi_f

出力2Vp-pは実用域で、THD+Nがおおむね0.07%程度なので満足な結果だと考えられる。
次に1kHzの信号における、出力振幅vsひずみ率を図3(Lch)、図4(Rch)に示す。

20250512thdamp_l
図3.出力振幅vsひずみ率(Lch)


20250512thdamp_r
図4.出力振幅vsひずみ率(Rch)


以上のように、低振幅(低入力信号)のTHD+Nではノイズ成分が支配的だった。


次はいよいよCR型RIAA回路を含めた全体の特性について、以下に示す。

まずは本回路をLTSPICEでシミュレーションした結果を表2に示す。

表2.LTSPICEによるシミュレーション結果
20250512ltspice_sim

このように、シミュレーションでは全帯域にわたってRIAA理論値に対する誤差が0.04dB以内に収まっており、非常に優秀だ。
実際の回路での測定結果を表3に示す。


表3.実機でのRIAA誤差
20250512riaaerror

最悪値はLchの20kHzで+0.345dBとなった。シミュレーションよりも大幅に悪化したが、実用上は十分だと思われる。
これはデジタルオシロのカーソルを使用して信号のp-pを測定した。

回路について注記すべきことがいくつかある。
この回路ではカートリッジに対して+12Vから470kΩ(回路図のR1)で微小な直流電流を流している。これはこのようにバイアスをかけることで音質が向上するという説があるためで、もし効果がないと思われる場合は未実装とする。
今回開発した基板では、フィルムコンデンサをすべてWIMAのMKS2で統一した。終段のカップリングコンデンサ10μFは比較的入手しにくいため、5mmピッチとして電解コンでも代用できるようにした。電解コンを使用する場合はゼロ信号付近で音質に影響を与えやすいため、出力信号が0Vをまたがないようにオフセットが2V程度付加するように調整されたい。


トラ技で販売したSSDACのアクリルケースを流用して組み込みました(写真2)。

20250512nnfb_riaa_eq
写真2.無帰還CR型フォノイコライザ完成



本機の音質は、まったくクセがなくたいへん素直で、ノイズも低く、とても満足な結果となりました(^-^)


ご希望の方に基板を頒布します。
生基板2枚ひと組(ステレオ分)と、製作説明書と回路図、部品表等の資料、製作例、LTSPICEのシミュレーションファイルをセットで3000円(税・送料込み)で頒布します。自力で部品収集、部品選別、調整できる方が対象です。

ご希望の方は表題に「無帰還イコライザアンプ基板頒布」、本文にお名前、送付先郵便番号、ご住所、電話番号をお書きのうえ、

dj_higo_officialアットhigon.sakura.ne.jp
(アットを@に替えてお送りください)

までメールをお送りください。

代金の振込先のご案内メールをお送りします。入金が確認でき次第発送します。

製作マニュアルは以下をダウンロードしてご覧ください。
製作マニュアル

 

| | | コメント (0)

2025年2月15日 (土)

動作に謎のあるAMラジオ

20250215radiopic
写真1.再現してみたFET1石ラジオ


2025年になってはや1ヶ月経った。2028年まであと3年になってしまった……

技術関係のノートを見ていたら、何年か前にお遊びで実験したAMラジオの回路が出てきた。
これはごく普通にFETで一段高周波増幅をしてゲルマダイオードで検波してクリスタルイヤホンで聴くというものだった(図1)。
20250215amradio_20250215192901
図1.AMラジオ回路


C1はどこで買ったかわからないAM用のポリバリコン、L2はだいぶ前にたぶんaitendoで買ったバーアンテナ、ゲルマニウムダイオードはたぶん1N60同等品、FETは秋月で売っている高周波用のJFET、BF256Bだ。イヤホンはクリスタルイヤホンのフリをしたセラミックイヤホン。

実験現場は東京都中野区の鉄筋コンクリートのマンション3階の室内で、ゲルマラジオだと電灯線アンテナを使ってTBS(954kHz)だけがかろうじて聞こえる(何を言っているかは聞き取れない)程度だ。
それに対し、FET1石でこの回路を組むと、バーアンテナに50センチほどのアンテナ線をつけた状態で、
・TBSは内容が聴き取れる
・NHK第1、NHK第2、AFN、文化放送は内容までは聴き取れないが、受信できていることがわかる
という具合に、感度が大幅に改善する。
チューニングをしていると、局の近くでキュルキュルと再生式ラジオのような音がするので、どうも再生気味に動作しているような感じがするのだが、回路は再生がかかるようにはなっていない。ためしにFETのソースやドレインからバーアンテナの2次巻線にFB(フィードバック)するなど、明示的なFBをかけてあれこれいじってみたが、うまく動作しなかった。結局よくわからないままノートに書いてそのままになっていた。

今回ノートから出てきたこの回路をもう一度組んでみたが、やはり上に書いた受信状態は再現した。コツは電源に電池を使うこと。スイッチング電源やUSB電源などではノイズが乗ってラジオが聞こえないからだ。

なぜ今になってAMラジオを組んでみたのかというと、AM局廃止期限まであと3年しかなく、この回路ももう作ることもなくなるかもしれないと思ったからだ。備忘録としてブログに書いておこうというわけだ。
今ある国内のAMラジオ局は、秋田と北海道の3局、それにNHKとAFNを除いて、すべて2028年までに停止となる。ゲルマラジオが発端となって技術者になったものとしては少し寂しい気がする……


【20250218 追記】
回路をもう少し変更して追加実験を行った。今回検討した回路を図2に示す。20250218amradio_01
図2.追加実験を行った回路


今回試してみた回路は、FETのドレインからバーアンテナのアンテナコイルに明示的にFBをかけて、FB量をVR1で調整できるようにした。図中L2,L3がバーアンテナで、L2が同調用コイル、L3がアンテナコイル。

最初の図1の回路では、明示的なFBが無いにもかかわらず再生動作しているようだったが、再生を通り越して発振しているようでキュルキュル音がして、局によってはこの発振音がじゃまをして音声が聞き取れなかった。そこで今回はFB量を調整できるようにしてみた。

また図1の回路では検波回路にD1が付いていたが、これは無しにした方が若干聞き取りやすかった。

VR1で聴きやすいFB量を調整することによって、TBSだけでなくNHK第2も内容が聞き取れるようになり、AFNはかかっている曲がわかるようになった。NHK第1は音声が小さく聞き取りは難しかった。文化放送は微弱で入感があったり無かったり、ニッポン放送はまったく入感がなかった。

この回路では消費電流は10~15mA程度。今回は極力シンプルな回路で試したかったのでソース抵抗なども省いて、消費電流はかなり多めだと思うが、もし実用を考えるなら感度との兼ね合いを見ながらソースやドレインに抵抗を入れるのが現実的だと思う。

自作ラジオの醍醐味は、ゲルマラジオか、あるいは一石か二石のシンプルな回路でどれだけ高感度にできるかというところにあると思うが、中でも再生式や超再生式は非常に魅力のある分野で、奥が深く、のめり込むとキリがない。実用性ではスーパーヘテロダインにかなうものではないが、スーパーラジオはあまりにもあたりまえに高性能が得られるのでかえって面白みに欠けるような気がする。


| | | コメント (0)

2024年12月17日 (火)

RIAAフォノイコライザを FIR デジタルフィルタで実装する ④ 実践&聴き比べ編

【参考記事】
RIAAフォノイコライザをIIRデジタルフィルタで実装する ①準備編
RIAAフォノイコライザをIIRデジタルフィルタで実装する ②係数計算編
RIAAフォノイコライザをIIRデジタルフィルタで実装する ③シミュレーション編
RIAAフォノイコライザをIIRデジタルフィルタで実装する ④実装・評価編
RIAAフォノイコライザをIIRデジタルフィルタで実装する ⑤Octaveによる位相検証

RIAAフォノイコライザを FIR デジタルフィルタで実装する ① 予告編
RIAAフォノイコライザを FIR デジタルフィルタで実装する ② 係数計算編
RIAAフォノイコライザを FIR デジタルフィルタで実装する ③ シミュレーション編

前回までで、実際にダイレクトリッピングしたwavファイルを入力して、FIRおよびIIRフィルタでRIAAイコライジング処理をする具体的な方法を説明しました。

今回は私がもっている何枚かのアナログレコードからパソコンでダイレクトリッピングしたサンプルソースを紹介し、実際に聞き比べた結果をお知らせしたいと思います。

まずはダイレクトリッピングする際に使用した機材です。図1がそのブロック図、図2は使用したフラットアンプの回路です。

Rippingblockdiagram
図1.ダイレクトリッピングブロック図


Flatampsch
図2.フラットアンプ回路(電源は±12V)


実際のリッピング作業では、SoundForgePro11を使ってファイルフォーマットは96kHz24bitのwavで行い、今回の検証用に、同じSound ForgePro11を使って48kHz16bitにダウンサンプリングしました。
ダウンサンプリングした理由は、今回3種類のフィルタタイプによる比較を行うにあたって、どちらかというとロースペックにした方が差が出やすいのではないかと思ったからです。

今回試聴した曲を何曲か紹介します。以下はすべてダイレクトリッピングしてダウンサンプリングした48kHz16bitフォーマットですので、実際にOctaveでRIAAイコライザを通して聴くことができます。

かっこ内はアーティスト名と収録DISKです。

Minor Swing(Stephane Grappelli/Live in Sanfransisco)

Whisper Not(阿川泰子/Journey)

Autumn Leaves(Bill Evans/PORTRAIT IN JAZZ)

夏服のイブ(松田聖子/Seiko Town)

ロックンルージュ(松田聖子/Seiko Town)

Take Five(Dave Brubeck/TAKE FIVE(Single))

Mambo Mongo(Belmonte & His Afro Latin7/OLE)

NIGHT IN TUNISIA(HOT SALSA/Hot Salsa meets Swedish Jazz)


さて、実際に試聴してみた印象ですが、思ったほど違いはありませんでした。というかほぼ同じです。
ただ、ほんのわずかにですが、FIRの直線位相タイプは女性ボーカルのサ行やパーカッションのハイハット、トランペットの音、バイオリンの弦がこすれる音や、音が細く掠れる場面などが少し物足りないというか、リアリティに欠けるように聞こえる印象がしました。
FIRの非直線位相とIIRはほぼ違いがなく、FIRの直線位相に比べるとリアルな音、という印象です。
ピアノの音や、大編成の演奏の音はほとんど差を感じませんでした。

もっともこれはブラインドテストではなく、自分でどのEQでかけているかわかっていますから、先入観でそう感じた可能性は多いにあります。
機会があれば、少し広めの会場でたくさんの人に聴いてもらって比較をしてみたいと思います。

みなさんもぜひ試してみてくださいね(^-^)


【20250118追記】
2025年1月18日、秋葉原で行われた手作りアンプの会の三土会にて、この記事で比較した3種類のデジタルRIAAイコライザの試聴会を行いました。試聴方法は、イコライザタイプ(FIR非直線位相、IIR、FIR直線位相)に①~③の番号をつけて、タイプの異なる3曲について聴き比べを行い、あとで①~③のうちどれがよかったか、あるいは違いがわからなかったということをみなさんに伺いました。結果が出るまでどのタイプが何番かは伏せておきました。

参加人数は14~15名、会場は秋葉原の万世区民会館7階洋室Aで、この会場は会議室であるため音響的にはライブでした。

結果はとても意外で、FIR非直線位相(アナログ等価)を選んだ人が2名、IIR(アナログ等価)が2名、FIR直線位相が残り全員となりました(@_@)
また、違いがわからないという人はいませんでした。

今回の結果から、理論的に正しい音が心地よいとは限らないということです。

IIRタイプについてはナイキスト周波数に近づくほど誤差が大きくなり、このときゲイン誤差と位相誤差がトレードオフの関係にあるため、どちらを優先した方が良いか?という試聴会も過去に行ったことがありましたが、このときはゲイン誤差は耳で差がわかるが、位相誤差はまったくわからない、という結果になっています。
そうすると、人の耳は位相には鈍感らしいということがわかります。

何事もやってみないとわかりませんね。オーディオは奥が深い……


| | | コメント (0)

2024年12月16日 (月)

RIAAフォノイコライザを FIR デジタルフィルタで実装する ③ シミュレーション編

【参考記事】
RIAAフォノイコライザをIIRデジタルフィルタで実装する ①準備編
RIAAフォノイコライザをIIRデジタルフィルタで実装する ②係数計算編
RIAAフォノイコライザをIIRデジタルフィルタで実装する ③シミュレーション編
RIAAフォノイコライザをIIRデジタルフィルタで実装する ④実装・評価編
RIAAフォノイコライザをIIRデジタルフィルタで実装する ⑤Octaveによる位相検証

RIAAフォノイコライザを FIR デジタルフィルタで実装する ① 予告編
RIAAフォノイコライザを FIR デジタルフィルタで実装する ② 係数計算編

この記事はアナログレコード再生の際に使用するRIAAフォノイコライザをFIRフィルタで実装し、シミュレーションと実際の音出しを行って評価を行うことを目的とするものです。


前回までで、今回の目標である2種類のFIRフィルタの係数が計算できました。
2種類のFIRフィルタとはすなわち、直線位相タイプと非直線位相タイプを指します。
RIAA特性はアナログフィルタの時定数で定義されており、当然ながら位相も含めて録音時の逆特性で再生することが求められるため、RIAAイコライザをデジタルフィルタで再現する場合もアナログRIAAイコライザと等価な特性を目指すべきです。

しかしながらFIRデジタルフィルタでは、アナログでは実現できない直線位相フィルタが構成できるため、FIRフィルタで直線位相のRIAAイコライザを構成したらどのような音になるのか?アナログの特性に近似の、非直線位相のFIRフィルタや、IIRフィルタで構成した場合とで聴感上の違いがあるのか?という検証を行うことが今回の目的です。

それではさっそく始めましょう。

前回、直線位相FIRと非直線位相FIRの2種類のRIAAの係数を計算し、それぞれ8192タップの係数のファイルを得ました。
今回はこれらの係数ファイルを使って、特性のグラフ表示と、実際に音楽信号(wavファイル)を処理した出力(wavファイル)を得ることを目標とします。
今回もOctaveを使って進めます。

まずはFIRフィルタの係数ファイルと入力のwavファイルを与えると、結果をwavファイルで出力するプログラムです。

%%%%%ここから

pkg load signal;

% FIRフィルタの係数を読み込み
coeff = load('coeff_nonlin6db_48000.txt');%%%%%%%%%%%%  ①使用する係数ファイル

% 音声ファイルの読み込み
inputFile = "alice4816_direct_15s.wav";%%%%%%%%%%%%%  ②処理するwavファイル
[~,name,ext] = fileparts(inputFile);
inputFileName = [name,ext];
[input_signal, fs] = audioread(inputFile);%%%%%%%%%%%%% 入力信号とサンプリング周波数を取得

% フィルタリング
filtered_signal = filter(coeff, 1, input_signal);%%%%%%%%%%% フィルタ適用

% 出力ファイルを保存
outputFileName = sprintf('FIR_RIAA_%s',inputFileName);
audiowrite(outputFileName, filtered_signal, fs);

% 周波数特性を計算
[h, w] = freqz(coeff, 1, 4096, fs);%%%%%%%%%%%%%%%% 周波数特性を計算

% ゲイン(dB)と位相(degree)を計算
gain = 20 * log10(abs(h));%%%%%%%%%%%%%%%%%%% ゲイン(dB)
phase = angle(h) * (180 / pi);%%%%%%%%%%%%%%%%%%位相(degree)

% プロット
figure;

% ゲイン特性
subplot(2, 1, 1);
semilogx(w,gain);
grid;
title(sprintf('FIR RIAA EQ %dHz Sampling',fs));
xlabel('Frequency (Hz)');
ylabel('Gain (dB)');
xlim([10, 30000]);
ylim([-20, 30]);%%%%%%%%%% 必要に応じて調整

% 位相特性
subplot(2, 1, 2);
semilogx(w, phase);
grid;
xlabel('Frequency (Hz)');
ylabel('Phase (Deg)');
xlim([10, 30000]);

%%%%%ここまで

使い方は簡単で、まず上のコードをOctaveのエディタにコピペして、適当な名前でOctaveのワークスペースのディレクトリに保存します。
次に、使用する係数ファイルを①のところに記述します。前回のファイルはcoeff_lin6db_48000.txtとcoeff_nonlin6db_48000.txtの2種類でしたので、都度書き換えて使います。
次に、実際に処理したいwavファイルを②のところに指定します。
これらのファイルはすべてOctaveのワークスペースに配置します。
処理したいwavファイルが他の場所にある場合は、フルパスで指定します。形式は"c:/AAA/bbb/CCC.wav"のような書式で、バックスラッシュではなくスラッシュなので注意が必要です。フルパスを打ち込むのが面倒な場合は、Octaveのコマンドウィンドウを表示して、そこに目的のファイルをドラッグ&ドロップすればフルパスが表示されるので、それをコピペすると便利です。

あとは実行ボタンを押せば、特性が表示され、"FIR_RIAA_(入力ファイル名)"でRIAA処理されたファイルが出力されます。
今回は48kHzサンプリングでレコードからダイレクト録音したソースを使って検証しており、処理時間はだいたいwavファイルの再生時間と同じくらいかかります。つまり5分の音楽ファイルなら、処理に5分ほどかかりました。特性は図1のように表示されます。
実行ボタンを押したらすぐにコマンドウィンドウでエラーが出ていないか確認してください。エラーが出ているといくら待っても何も起こりません(^-^;

Fir_riaa_nonlinear_towav

図1.表示された非直線位相FIRの特性


直線位相タイプの係数も入れてみましょう。
リストの①のところにcoeff_lin6db_48000.txtと入力して実行します。図2に表示された特性を示します。


Fir_riaa_linear_towav
図2.表示された直線位相FIRの特性

ここまでで、処理されたwavファイルも出力されたと思います。出力wavファイルは同名の場合は上書きされてしまいますので、都度どこかに待避させてください。


次に、IIRフィルタで同じ処理を行ってみましょう。プログラムは次のとおりです。

%%%%%ここから

pkg load signal;

% フィルタ係数を定義
%{
%96kHz プリワーピングなし +6dB
#ゲインはb0のみで調整
b0 = 0.26;
b1 = [1 0];
b2 = [1 -0.96777105 0];
a1 = [1 0];
a2 = [1 -1.866859545 0.867284263];
%}

%
%48kHz プリワーピングなし +6dB
#ゲインはb0のみで調整
b0 = 0.5;
b1 = [1 0];
b2 = [1 -0.936564324 0];
a1 = [1 0];
a2 = [1 -1.749567588 0.751160265];
%}

% IIRフィルタ係数の計算
b = b0 * conv(b1, b2);
a = conv(a1, a2);

% 音声ファイルの読み込み
inputFile = "alice4816_direct_15s.wav";%%%%%%%%%③
[~,name,ext] = fileparts(inputFile);
inputFileName = [name,ext];
inputFileName
[input_signal, fs] = audioread(inputFile);

% フィルタリング
filtered_signal = filter(b, a, input_signal);

% 出力ファイルを保存
outputFileName = sprintf('IIR_RIAA_%s',inputFileName);
audiowrite(outputFileName, filtered_signal, fs);

% 周波数特性を計算
[h, w] = freqz(b, a, 4096, 'whole', fs);

% ゲイン(dB)と位相(degree)を計算
gain = 20 * log10(abs(h));%%%%%%%%%%% ゲイン (dB)
phase = angle(h) * (180 / pi);%%%%%%%%% 位相 (degree)

% プロット
figure;

%ゲイン特性
subplot(2,1,1);
semilogx(w,gain);
grid;
title(sprintf('IIR RIAA EQ %dHz Sampling',fs));
xlabel('Frequency (Hz)');
ylabel('Gain (dB)');
xlim([10 30000]);

%位相特性
subplot(2,1,2);
semilogx(w,phase);
grid;
xlabel('Frequency (Hz)');
ylabel('Phase (Deg)');
xlim([10 30000]);
ylim([-100, 0]);%%%%%%%%%%%%%%% 必要に応じて調整

%%%%%ここまで

このコードでは、96kHzサンプリングと48kHzサンプリングの2種類のフィルタ係数が書かれていて、96kHzはコメントアウトされ、48kHzサンプリングが使えるようになっています。IIRフィルタでの係数計算については
RIAAフォノイコライザをIIRデジタルフィルタで実装する ②係数計算編
RIAAフォノイコライザをIIRデジタルフィルタで実装する ③シミュレーション編
を参照してください。

処理したいwavファイルを③のところに指定して実行すれば、同じように"IIR_RIAA_(入力ファイル名)"でRIAA処理されたファイルが出力されます。
このIIRフィルタの処理は時間がほとんどかからず、即座に結果が出力されます。

図3に表示された特性を示します。


Iir_riaa_towav_20241214210201


図3.表示されたIIRフィルタによるRIAAの特性


これから何曲か処理して、聴き比べをしてみようと思います。はたして聴感的に差があるのか……!?


RIAAフォノイコライザを FIR デジタルフィルタで実装する ④ 実装&聴き比べ編へ続く

| | | コメント (0)

2024年12月14日 (土)

RIAAフォノイコライザを FIR デジタルフィルタで実装する ② 係数計算編

【参考記事】
RIAAフォノイコライザをIIRデジタルフィルタで実装する ①準備編
RIAAフォノイコライザをIIRデジタルフィルタで実装する ②係数計算編
RIAAフォノイコライザをIIRデジタルフィルタで実装する ③シミュレーション編
RIAAフォノイコライザをIIRデジタルフィルタで実装する ④実装・評価編
RIAAフォノイコライザをIIRデジタルフィルタで実装する ⑤Octaveによる位相検証

RIAAフォノイコライザを FIR デジタルフィルタで実装する ① 予告編


今回はフリーで使える数値解析用プログラミング言語Octaveを使って、アナログレコード再生用のRIAAイコライザのシミュレーションを行うために、まずはRIAAイコライザ回路をFIRフィルタで構成するためのフィルタ係数を求めます。

最終的な目標は、
・IIRフィルタでRIAAイコライザ回路(アナログ等価=非直線位相)
・FIRフィルタでRIAAイコライザ回路(アナログ等価=非直線位相)
・FIRフィルタでRIAAイコライザ回路(直線位相)
という3種類のシミュレーションと試聴を行い、それぞれの特性と実際の音の違いを検証することです。

ではさっそく、FIRフィルタでRIAAイコライザ回路を構成する作業を説明します。
デジタルフィルタとはなにか?などについては上に書いた過去記事やWikiなどを参照してください。

FIRフィルタは、フィルタ係数を対称形にすることで直線位相フィルタを構成できることが大きな特徴で、今回は直線位相と非直線位相(アナログに近似)の2種類のFIRフィルタ設計でRIAAイコライザ回路を構成します。

任意の周波数特性に対し、逆フーリエ変換を行うと、その周波数特性のインパルスを求めることができ、そのインパルスを係数とするFIRフィルタを構成すれば、もとの周波数特性をもつフィルタができます。

この作業を行う際に、元の周波数特性に複素共役対称性を持たせる(=負の周波数領域を追加して、対称性を持たせる)ことで、逆フーリエ変換の結果が対称性を持ち、直線位相のFIRフィルタを構成することができます。

逆に、元の周波数特性に複素共役対称性を適用せず、正の周波数領域のみで逆フーリエ変換を行った場合は、非直線位相(=アナログで構成するフィルタに近似)となります。

それでは実際の作業に入りましょう。作業環境はWindows11にGNU Octave 8.3.0です。

まずはFIRフィルタの特徴である、直線位相特性をもつRIAAイコライザのフィルタ係数生成です。つぎのコードをOctaveのエディタにコピペして実行すれば、RIAAイコライザを構成するFIRフィルタの係数ファイルが生成されます。

%%%%%ここから

pkg load signal;

% サンプリング周波数とフィルタ設定
fs = 48000;%%%%%%%%%%%%%%%  サンプリング周波数 (Hz)
N = 8192;%%%%%%%%%%%%%%%% FIR フィルタのタップ数
f = linspace(0, fs/2, N/2+1);%%%%%%%  周波数範囲(片側スペクトル)
current_dir=pwd;

% RIAA伝達関数の定義
T1 = 3.18e-3;% 時定数1
T2 = 318e-6;% 時定数2
T3 = 75e-6;% 時定数3
s = 2 * pi * 1i * f; % s = jω
Hs = (1 + s*T2) ./ ((1 + s*T1) .* (1 + s*T3));

% 時間領域のインパルス応答を逆フーリエ変換で取得
Hf_full = [Hs, conj(flip(Hs(2:end-1)))];%%%%%%%直線位相特性を得るため複素共役対称性を持たせる
h = real(ifft(Hf_full,N));%%%%%%%%%%%%%% 逆フーリエ変換でインパルス応答を取得
h = circshift(h, -N/2);%%%%%%%%%%%%%%%%ゼロ位相シフト
h_6db = 20*h;%%%%%%%%%%%%%%%%%%%%+6dB@1kHzに調整

% FIRフィルタの周波数特性
[Hfir, w] = freqz(h, 1, N, fs);

% gainとphase
gain = 20*log10(20*abs(Hfir));%%%%%ゲイン特性 (dB)
phase = angle(Hfir) * (180/pi);%%%%%位相特性 (degree)

% FIRフィルタの係数を保存
save_path=[pwd,sprintf("/coeff_lin6db_%d.txt",fs)]
save("-ascii", save_path, "h_6db");
#save(coeff_file, 'h_fir', '-ascii');

% ゲイン特性と位相特性のプロット
figure;
subplot(2,1,1);
semilogx(w, gain);
grid on;
xlim([10 24000]);
ylim([-30 30]);
xlabel('Frequency (Hz)');
ylabel('Gain (dB)');
title(sprintf('Linear FIR RIAA Response %dHzSampling',fs));

subplot(2,1,2);
semilogx(w, phase);
grid on;
xlim([10 24000]);
xlabel('Frequency (Hz)');
ylabel('Phase (Deg)');

% FIRフィルタ係数をファイルに保存
coeff_file = 'current_dir/coeff.txt';
disp(['FIR filter coefficients saved to ', coeff_file]);

%%%%%ここまで

サンプリング周波数fsとタップ数Nは任意の数値を入れられますが、今回はfs=48000,N=8192としています。
Nは2^nを選びます。
FIRフィルタの低域特性は分解能fs/Nで決まり、この場合は48000/8192=5.86 (Hz)
となるので、周波数帯域は5.86Hz~24kHzとなります。
ゲインは1kHzで+6dBになるように調整しています。

このコードの中で注目すべき箇所は、
Hf_full = [Hs, conj(flip(Hs(2:end-1)))];%%%%%%% 直線位相特性を得るため複素共役対称性を持たせる
h = real(ifft(Hf_full,N));%%%%%%%%%%%%%%  逆フーリエ変換でインパルス応答を取得

の部分です。直線位相特性を得るために、求める周波数特性に複素共役対称性をもたせて、周波数特性を正負対称にしています。
そして、ifft関数(逆フーリエ変換)を使って、この周波数特性のインパルスを得ています。

生成されたFIRフィルタの係数はcoeff_lin6db_48000.txtというファイルとして保存されます。これには8192個の係数が書かれています。
係数ファイルが生成されると同時に、この係数を使った場合のFIRフィルタの特性のグラフが表示されます(図1)。

Riaa_linear_gp

図1.FIRフィルタによるRIAAイコライザ(直線位相)のGain-Phase特性


次に、同じFIRですがこんどはアナログ回路に近い非直線位相のFIRフィルタ係数を生成します。

%%%%%ここから

pkg load signal;

% サンプリング周波数とフィルタ設定
fs = 48000;% サンプリング周波数 (Hz)
N = 8192;% FIR フィルタのタップ数
f = linspace(0, fs/2, N/2+1);% 周波数範囲(片側スペクトル)
current_dir=pwd;

% RIAA伝達関数の定義
T1 = 3.18e-3;%%%%%%%%%%%%%%%%% 時定数1
T2 = 318e-6;%%%%%%%%%%%%%%%%%% 時定数2
T3 = 75e-6;%%%%%%%%%%%%%%%%%%% 時定数3
s = 2 * pi * 1i * f; % s = jω
Hs = (1 + s*T2) ./ ((1 + s*T1) .* (1 + s*T3));%%%周波数特性

% 非直線位相を持つフィルタ設計のため、複素共役対称性を適用しない
%Hf_full = Hs;%%%%%%%%%%% 複素共役対称性をもたせず、そのまま使用

% 時間領域のインパルス応答を逆フーリエ変換で取得
h = real(ifft(Hs, N));%%%%%%%%% 時間領域の係数
h_6db = 36*h;%%%%%%%%%%%%+6dB@1kHzに調整

% FIRフィルタの周波数特性
[H, f_resp] = freqz(h, 1, N, fs);%%%%% 周波数特性の計算

% gainとphase
gain = 20*log10(36*abs(H));%%%%%%% ゲイン特性 (dB)
phase = angle(H) * (180/pi);%%%%%%% 位相特性 (degree)

% FIRフィルタの係数を保存
save_path=[pwd,sprintf("/coeff_nonlin6db_%d.txt",fs)]
save("-ascii", save_path, "h_6db");

% ゲイン特性と位相特性のプロット
figure;
subplot(2, 1, 1);
semilogx(f_resp, gain, 'b');
xlabel('Frequency (Hz)');
ylabel('Gain (dB)');
title(sprintf('non-Linear FIR RIAA Response %dHzSampling',fs));
xlim([10 24000]);
ylim([-30 30]);
grid on;

subplot(2, 1, 2);
semilogx(f_resp, phase, 'r');
xlabel('Frequency (Hz)');
ylabel('Phase (Deg)');
xlim([10 24000]);
ylim([-100 0]);
grid on;

%%%%%ここまで

使い方は上と同じですが、周波数特性に複素共役対称性を適用せず、正の周波数特性のみを使用し、非直線位相特性になっています。
得られる係数ファイルはcoeff_nonlin6db_48000.txt、特性図は図2のとおりです。

Riaa_nonlinear_gp_20241214180201

図2.FIRフィルタによるRIAAイコライザ(非直線位相)のGain-Phase特性


これでFIRの係数が計算できたので、次回はその係数をFIRフィルタに実装して、再度Gain-Phase特性を確認し、実際に音楽を再生します。


RIAAフォノイコライザを FIR デジタルフィルタで実装する ③ シミュレーション編へ続く





| | | コメント (3)

2024年12月10日 (火)

RIAAフォノイコライザを FIR デジタルフィルタで実装する ① 予告編

この記事は次の参考記事の続編です。
今回はRIAAフォノイコライザをFIRフィルタで実装するシミュレーションを行いました。

【参考記事】
RIAAフォノイコライザをIIRデジタルフィルタで実装する ①準備編
RIAAフォノイコライザをIIRデジタルフィルタで実装する ②係数計算編
RIAAフォノイコライザをIIRデジタルフィルタで実装する ③シミュレーション編
RIAAフォノイコライザをIIRデジタルフィルタで実装する ④実装・評価編
RIAAフォノイコライザをIIRデジタルフィルタで実装する ⑤Octaveによる位相検証

デジタルRIAAフォノイコライザ・サブ基板(基板頒布あり)
SSDAC128 and デジタルRIAAイコライザ デュアル基板(基板頒布あり)

上の参考記事で書きましたが、アナログレコードはRIAAカーブという周波数特性で録音されているので、正しい音で再生するにはRIAAイコライザを通す必要があります。従来はこれをアナログ回路で行っていましたが、IIRデジタルフィルタで構成することも可能で、IIRフィルタではゲイン特性、位相特性ともにアナログフィルタと等価な特性を持たせることができます。

ところで、デジタルフィルタにはIIRの他にFIRというタイプのものがあって、特徴は直線位相特性をもつことです。つまり全帯域において一律にサンプリング時間の遅れのみ発生し、アナログフィルタのように周波数による位相回りが起こらないという特徴があります。

※すみません、最初に記事を書いた時点で間違った思い込みがありました。
正しくは、「係数が対称なFIRフィルタは直線位相になる」であって、係数が対称でなければ直線位相にはなりません。
また直線位相特性が実現できることはFIRフィルタの重要な性質で、IIRでは実現できません。

過去記事ですでに述べたとおりRIAA規格によるレコード録音再生では、ゲイン特性に伴い位相特性も再現する必要があるため、当然ながら位相特性も含めてアナログ回路によるRIAAカーブと等価な特性を再現できるIIRフィルタで構成するのが妥当と考えられます。
ただ、正しい再生ができるかどうかとは別に、FIRフィルタでRIAAカーブを再現するとどのようになるのか、というのは興味があります。

FIRフィルタの係数はどのように求められるかというと、周波数特性を逆フーリエ変換したものがそのインパルス関数となるため、それをFIRフィルタの係数に使うことで実装できます。
ところが、検証に使っていたプログラミング言語”Octave”の逆フーリエ変換"ifft"の使い方がわからず、保留にしていました。

今回別の課題があって、Linux上でsoxやらbrutefirやらで遊んでいたときに、ふと、
「ひょっとしてChatGPT先生に訊いたら、FIRでRIAAイコライザを設計してくれるかも……!?」
と思い、実際にChatGPT先生に訊いてみたら、教えてもらえました\(^o^)/
なんと、Pythonのnumpyに逆フーリエ変換があって、それが使えるのです!

そういうわけで、ChatGPTに教えてもらったPythonのコードを手直ししたものを以下に示します。


############ 以下 コード

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

# サンプリング周波数
fs = 44100
nyquist = fs / 2

# RIAAカーブの時間定数
t1 = 3180e-6 # 3180 µs
t2 = 318e-6 # 318 µs
t3 = 75e-6 # 75 µs

# 周波数特性の計算
frequencies = np.linspace(10, nyquist, 4096) # 10Hzから始める
s = 2j * np.pi * frequencies # sはフーリエ演算子、jは虚数
riaa = (1 + t2 * s) / ((1 + t1 * s) * (1 + t3 * s)) #RIAAの伝達特性

# 振幅特性(対数スケールでのプロット用)
gain = 20 * np.log10(np.abs(riaa))
gain -= gain[np.argmin(np.abs(frequencies - 1000))] # 1kHzで0dBに正規化

# 位相特性(ラジアンを度に変換)
phase = np.angle(riaa, deg=True)

# FIRフィルタ設計(逆FFT)
N = 8192 # フィルタのタップ数
impulse_response = np.fft.ifft(riaa, n=N) #逆FFTを求める
fir_coeffs = np.real(impulse_response[:N])

# FIR係数を保存
np.savetxt("riaa_coefficients.dat", fir_coeffs)

# FIRフィルタの周波数応答を確認
w, h = freqz(fir_coeffs, worN=4096, fs=fs)
fir_gain = 20 * np.log10(np.abs(h))
fir_gain -= fir_gain[np.argmin(np.abs(w - 1000))] # 1kHzで0dBに正規化
fir_phase = np.angle(h, deg=True)

# グラフ描画
fig, ax1 = plt.subplots(figsize=(12, 8))

# 振幅特性(ゲイン)のプロット
ax1.set_title("RIAA EQ Frequency and Phase Response")
ax1.semilogx(frequencies, gain, label="Ideal RIAA Gain", linewidth=1.5)
ax1.semilogx(w, fir_gain, label="FIR Filter Gain", linewidth=1.5)
ax1.set_xlabel("Frequency (Hz) [Log Scale]")
ax1.set_ylabel("Gain (dB)")
ax1.set_xlim(10, 30000) # 横軸を10Hzから30kHzに設定
ax1.set_ylim(-30, 30) # 縦軸を30dBから-30dBに設定
ax1.grid(which="both", linestyle="--", linewidth=0.5)
ax1.axvline(1000, color="gray", linestyle="--", linewidth=0.8, label="1 kHz Reference")
ax1.legend(loc="upper left")

# 位相特性のプロット
ax2 = ax1.twinx() # 第二軸(位相特性用)
ax2.semilogx(frequencies, phase, label="Ideal RIAA Phase", color="tab:orange", linestyle="--", linewidth=1.5)
ax2.semilogx(w, fir_phase, label="FIR Filter Phase", color="tab:red", linestyle="--", linewidth=1.5)
ax2.set_ylabel("Phase (degrees)")
ax2.set_ylim(-180, 180) # 位相特性の範囲を-180度から180度に設定
ax2.legend(loc="lower right")

plt.show()


#############ここまで


この中でキモとなるのは、

s = 2j * np.pi * frequencies#フーリエ演算子s
riaa = (1 + t2 * s) / ((1 + t1 * s) * (1 + t3 * s))#RIAAの伝達関数

の、RIAA特性の伝達関数の部分と、

# FIRフィルタ設計(逆FFT)
N = 16384 # フィルタのタップ数
impulse_response = np.fft.ifft(riaa, n=N) #逆FFTを求める
fir_coeffs = np.real(impulse_response[:N])

# FIR係数を保存
np.savetxt("riaa_coefficients.dat", fir_coeffs)

ここの、逆FFTを行ってインパルス応答を求めて、FIR係数として保存している部分です。
このFIR係数を使えば実際のFIRフィルタに実装できます。

このプログラムを実行した結果は次のとおり。

20241210_fir_riaa_gp

図1.FIRで実装したRIAAイコライザの特性

図1は、青いプロットがRIAA特性の理論値、橙のプロットがFIRフィルタのシミュレーション結果です。
位相特性も含めてなかなかきれいに出ていますね。
低域で誤差が大きく出ているのは、サンプリング周波数(Fs)とタップ数(N)によって決まる分解能によるものです。この例では分解能Δfは

Δf=Fs/N=96000/8192=11.72(Hz)

となります。

そういうわけで、音が出せる状態になったらまた報告します(^-^)

RIAAフォノイコライザを FIR デジタルフィルタで実装する ② 係数計算編へ続く



| | | コメント (0)

2024年10月 3日 (木)

SSSRC( スーパーサンプリング・サンプリングレートコンバータ ) ソフトウェア公開

SSDACのPCアプリによる実装について、これまで次の記事で説明してきました。

SSDACをパソコン上のPythonで実装する ①

SSDACをパソコン上のPythonで実装する ②

これらの記事では、SSDACの目的や原理について説明し、
新たな発見として、スーパーサンプリングによってアップサンプリングしたwavデータは、一般的な市販品のオーバーサンプリングDACで再生しても新たにプリエコー、ポストエコーが発生しない、ということがわかりました。

ソフトウェアでスーパーサンプリングを行うことは、スプライン関数によって補間データを追加してファイルを生成していることから、スプライン関数補間によるサンプリングレートコンバートを行っているということになりますから、入力がCDフォーマットの44.1kHz16bitフォーマットの出力ファイルは、たとえば88.2kHz24bit(2倍スーパーサンプリング)、176.4kHz24bit(4倍スーパーサンプリング)という具合になります。

前回の記事では、88.2kHzサンプリング音源でプリエコー、ポストエコーが発生したとしても、その周波数は44.1kHzなので、聴感上の悪影響はほとんどないであろう事を予想していましたが、実際の検証では、スーパーサンプリングでアップサンプリングしたデジタル信号は、その後オーバーサンプリングDACで処理しても、プリエコー、ポストエコーが発生しない、という意外な結果を得ました。

また、前回までPC上のプログラミング言語Pythonを使ってスーパーサンプリングを実装し検証しましたが、Pythonはインタープリタ言語(逐次翻訳型言語)であることから非常に処理が遅く、15秒の音楽ファイルを処理するのに数時間かかるという具合で、およそ実用的ではありませんでした。

そこで今回は、wavファイルをスーパーサンプリング処理して同じくwavファイルで出力する、サンプリングレートコンバータ、名付けて、
スーパーサンプリング・サンプリングレートコンバータ(SSSRC)をC#コンパイラで実装しましたので、次の通り公開します。

今回作成したSSSRCソフトの仕様は次の通りです。

入力フォーマット : 44.1kHz 16bit (CDフォーマット)
出力フォーマット : 88.2kHz、176.4kHz、352.8kHz、705.6kHz(いずれも24bit)


●ダウンロードと使い方
このソフトはwindows10またはwindows11上で動作します。

ダウンロード - sssrc.zip

SSSRC.zipをダウンロードして任意の場所に解凍したら、SSSRC.exeをダブルクリックして実行します。
このとき、.NETのランタイムがインストールされていないと、次のようなエラーが出ることがあります。

Sssrcnet1
図1..NETランタイムがない場合のエラー

この場合は”Download it now”をクリックして、.NETランタイムをインストールしてから、あらためてSSSRC.exeを実行してください。


●使い方
図2はSSSRCの起動画面です。

Sssrc
図2.SSSRC起動画面


使い方は簡単で、希望するスーパーサンプリング倍率(Super Sampling Ratio)を選んで、変換したいwavファイル(44.1kHz16bitのみ)を、この画面上にドラッグ&ドロップすると、元のwavファイルと同じ場所に”SS_OUTPUT”ホルダが生成され、その中にアップサンプリングされたファイルが生成されます。複数ファイルを一括してドラッグ&ドロップすることができます。
変換時間は、私の古いパソコン(第5世代Corei7(2core)、メモリー16GB、windows11)で、4分(約40MB)のwavファイルを88.2kHzに変換するのに35秒、705.6kHzに変換するのに1分30秒かかりました。

実際にこのソフトで処理した信号の波形を見てみましょう。代表例として、1kHzの方形波と正弦波を4倍スーパーサンプリングでアップサンプリングした場合の波形を図3、図4に示します。

1korg
図3.元波形(44.1kHz16bitサンプリング)


1k4x
図4.4倍アップサンプリング波形(176.4kHz24bitサンプリング)

以上のように、スプライン関数による補間によりアップサンプリングしたファイルが出力されます。

ここで注意したいのは、出力ファイルの容量です。アップサンプリング倍率に加えて、ビット分解能が16bit→24bitと1.5倍となることから、
出力ファイルはアップサンプリング倍率x1.5倍の容量になります。
たとえば長さ4分で約40MBのwavファイルは、2倍アップサンプリングでは3倍の120MBに、16倍アップサンプリングでは24倍の960MBになりますので、ストレージを圧迫します。要注意です。

今後、今回のSSSRCで元ファイルをアップサンプリングして通常のDACで再生したものと、元ファイルをそのまま通常のDACで再生したもの、また、元ファイルをハードウェアのSSDACで再生したものの違いについて聴き比べと検証を行う予定です。


| | | コメント (0)