概要
はじめに
全3回に渡ってAC負荷アナライザを製作する手順を解説します。
第2回目の本記事では、dsPIC33AKを使用したマイコンソフトウェア設計について、具体的なプログラミング手法やポイントを取り上げます。
| 公開/ 変更日 | 特記 |
| 24/12/24 | 暫定版記事の公開。 完全版については、準備が整い次第、改めて公開いたします。 (開発の都合で、公開内容については修正される可能性が有りますのでご了承ください) |
| 24/12/26 | サンプル画像追加(fig.4-1)、字句修正。 |
| 24/12/28 | 皮相電力・有効電力の算出コード修正。 半田ごての計測例追加(fig.4-3-1)(fig.4-3-2) |
| 24/12/29 | ・各種演算処理呼び出し、ADC割り込み、リサンプリング処理追加 ・Pythonでの検算用コード追加 ・電気ストーブの計測例追加(fig.4-4-1)(fig.4-4-2) |
| 24/12/30 | ・各種演算処理呼び出し処理変更 ・各種演算実行時間比較表追加(fig.4-1-1) ・実行波形追加(fig.4-1-2) |
| 25/1/4 | ・字句修正 ・位相角算出の項目追加 |
| 25/1/5 | ・ブロック図修正(fig.2-3)(fig.2-4)(fig.2-5) ・字句修正 ・外部リンク追加 |
AC負荷アナライザは交流100Vを扱い感電のリスクを伴います。
この記事を参考にして生じた損害等については一切の責任を負いません。
関連リンク
内部リンク
外部リンク
本記事の内容は、以下のサイトの記事を参考にして製作を進めています。
開発環境
開発環境を示します。
| 項目 | 値 | リンク |
| 基板 | 自作基板 | |
| 統合開発環境 | MPLAB X IDE v6.20 | MPLAB® X IDE | Microchip Technology |
| コンパイラ | MPLAB XC DSC v3.10 | MPLAB® XC DSC Compiler | Microchip Technology |
全体ブロック図
以下に全体ブロック図を示します。

マイコンでは主に以下の3つの処理を行います。
・電圧・電流のサンプリング
・各種演算
・演算結果の送信



マイコンソフトウェア
1.検出部
1-1.アナログ電圧サンプリング
電圧は ADC1、電流は ADC2 を使ってシングルエンドモードで同時にサンプリングします。
当初は差動モードで計測する予定でしたが、後述の理由で変更しました。

サンプリング周期(fs)はアンチエイリアシングフィルタのカットオフ周波数である7.2kHzより遥かに高い 4.9152 MHz とし、ADC内蔵の 64倍オーバーサンプリング により、fs1 = 76.8 kHz ごとに結果を出力させます。なぜこのような中途半端な値にしたかと言うと、1サイクル(50Hzまたは60Hz)あたりのサンプル数をリサンプリング処理により256回 としたかったためです。これにより、後続の処理が 2のn乗 になるため、計算が効率的に行えるようになります。
これに関しては後述の”1-3.リサンプリング処理”にて説明します。
また、バッファは2サイクル分を確保し、Ping-Pongバッファのようにデータ入力とデータ処理をオーバーラップさせて効率的に処理を行っています。

シリコンエラッタ情報 DS80001139A では、ADCの差動モードにおいて演算結果に符号が付与されない問題が既知のエラッタとして公開されています。この問題については、ソフトウェアで対応する必要があるとされています。
しかし、公開されていない別のエラッタが存在するかも知れません(*1)。というのも、正入力と負入力の符号が変わるタイミングで、取得値に異常が生じる現象に遭遇したためです。
(*1)マイコンのエラッタではなく回路が原因で本現象が発生している可能性も有ります。
オシロスコープを用いてAD入力端子の波形を測定したところ、1.29Vを中心とした交流成分が確認されました。

この信号をシングルエンドモードで正入力のみ計測すると、綺麗なサイン波が得られます。一方で、差動モードで測定し、符号処理を加えると、ゼロクロス付近でグリッジ(瞬間的な不連続)が発生していることがわかりました。
パラメータをさまざまに調整して挙動を確認しましたが、問題を解消するには至りませんでした。そのため、今回はシングルエンドモードで計測することとしました。

1-2.サンプリングクロック生成
サンプリングクロックのトリガ(4.9152MHz)は、8MHzのFRCクロックをPLL1で600MHzに逓倍し、PLL VCOのFracDiv出力を使用して生成します。このクロックは、クロックジェネレーターのCLK Gen13を介して分周され出力されます。
直接クロックジェネレータの出力をサンプリングトリガとして使用できないため、CLC1を介しトリガにします。

CLC1の出力をRP52ピンから出力し、狙い通りの周波数が生成されていることを確認しました。

CLKGenのCLKxDIVレジスタは以下の通りINTDIV(15bit)とFRACDIV(9bit)部に分かれており、この機能により微細なクロック周波数が生成可能です。
具体的なINTDIVとFRACDIVの求め方は以下の通りです。
①以下の式が成り立つ最大のINTDIVを求める


- IF:入力クロック周波数(今回は600MHz)
- OF:出力周波数(今回は4.9152MHz)
②以下の式が成り立つFRACDIVを求める。




1-3.リサンプリング処理
4.9152MHzのクロックを使用した64倍オーバーサンプリングにより、データは76.8kHzでサンプリングされます。このデータを、商用周波数が60Hzの場合は5回に1回、50Hzの場合は6回に1回間引くことで、1サイクル(50Hzまたは60Hz)あたりのサンプル数を256回に調整します。
通常、サンプリング周波数を変更するマルチサンプリングレート変換処理には、アップサンプリング(サンプル数の増加)とダウンサンプリング(サンプル数の削減)の両方を使用します。これにより、例えば2.5倍といった非整数の間引き率も実現可能ですが、実際にはいくつかの考慮すべき点があります。
アップサンプリングを行う際、追加のサンプルを挿入する必要があるため、処理時間やメモリ使用量が増加するという問題があります。特に、非整数倍の間引き率を実現するためには、アップサンプリングの処理が複雑になり、計算量やメモリの消費が増える可能性があります。
そのため、今回は「5倍」と「6倍」といった整数倍の間引き率のみを使用し、アップサンプリングを避けることで、計算量やメモリ使用量を抑え、効率的にダウンサンプリングのみで処理を完結できるようにしています。
この際に単純に間引くと エイリアス の影響を受けて、信号が歪んでしまう可能性があります。そこで、通常はデシメーションフィルタ を使用して、信号の歪みを防ぎます。
この詳しい説明は、以下の東洋テクニカ様のサイトを参考にしてください。
AD変換の基礎 / 第2回 オーバーサンプリングによるアンチエイリアシング | 東陽テクニカ | “はかる”技術で未来を創る | オートモーティブ
サンプリングされた信号の周波数が実際の信号の周波数と異なって認識される現象
デシメーションフィルタ(Decimation filter)とは、信号処理において、サンプリングレートを下げるために使用されるフィルタです。デシメーションのプロセス自体は、入力信号を一定の割合で「間引く」ことで、信号のサンプル数を減らすことを指します。この過程で、データの量を減らし、計算量を削減する目的がありますが、そのためには適切なフィルタリングが必要です。
以下のグラフは、100Hzと450Hzの成分が含まれる信号を、サンプリング周波数 fs=1kHz で取得したデータです。この信号を250Hz(fs/4)のカットオフ周波数で帯域制限した後、間引き率2で間引いた場合と、帯域制限を行わずにそのまま間引いた場合の比較を示しています。
帯域制限せずにそのまま間引いた信号のスペクトルを見ると、50Hz(=250 -(450 – 250))の位置にエイリアスが出現していることが確認できます。
このような理由のため、フィルタで帯域制限を行ってから間引き処理を行います。
dsPICにはDSP命令を使用してフィルタリングと間引きを実行するFIRDecimate関数がありますので、通常はこれを利用します。

ただし、今回は、サンプリング周波数 fs = 4.9152 MHz で64倍のオーバーサンプリング を行い、76.8 kHz 間隔で結果を得ています。つまり現信号に対して 1000倍以上のレートになります。
例として、60Hzの基本周波数に加えて 3次、5次、7次、9次 の高調波がそれぞれ 20%、10%、5%、1% の割合で重畳したとして、これらの波形をLTSpiceでシミュレーションすると以下のように表されます。

Pythonで間引き前のフィルタ有無を比較できるソフトを作成し実行した結果、フィルタを適用しなくてもその影響はほぼ無視できるレベルであることが確認できました。したがって、今回はシンプルに間引いて処理しています。

■Pythonコード
import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import butter, lfilter
# Define parameters
fs = 76800 # Sampling frequency (76.8 kHz)
f_signal = 60 # Signal frequency (60 Hz)
n = 2 # Downsampling factor (adjust as needed)
fs_downsampled = fs // n # Downsampled frequency
duration = 0.1 # Signal duration in seconds (adjust as needed)
t = np.linspace(0, duration, int(fs * duration), endpoint=False) # Time vector
# Generate the signal with harmonic noise
harmonics = [3,5,7,9] # Harmonic multiples
amplitudes = [0.2,0.1,0.05,0.01] # Amplitudes of harmonics
signal = np.sin(2 * np.pi * f_signal * t) # Fundamental signal
for harmonic, amp in zip(harmonics, amplitudes):
signal += amp * np.sin(2 * np.pi * f_signal * harmonic * t)
# Define low-pass filter (LPF)
def butter_lowpass(cutoff, fs, order=4):
nyquist = 0.5 * fs
normal_cutoff = cutoff / nyquist
b, a = butter(order, normal_cutoff, btype='low', analog=False)
return b, a
def apply_lowpass_filter(data, cutoff, fs, order=4):
b, a = butter_lowpass(cutoff, fs, order=order)
return lfilter(b, a, data)
# Apply LPF with 19.2 kHz cutoff
cutoff_freq = fs / 4 # Cutoff frequency for LPF (19.2 kHz)
filtered_signal = apply_lowpass_filter(signal, cutoff_freq, fs)
# Downsample with LPF
signal_downsampled_with_lpf = filtered_signal[::n]
# Downsample without LPF
signal_downsampled_without_lpf = signal[::n]
# Create time vectors for downsampled signals
t_downsampled = np.linspace(0, duration, len(signal_downsampled_with_lpf), endpoint=False)
# Function to calculate FFT
def calculate_fft(data, fs):
n = len(data)
freq = np.fft.rfftfreq(n, d=1/fs)
fft_data = np.abs(np.fft.rfft(data)) / n
return freq, fft_data
# Calculate FFT for each signal
freq_original, fft_original = calculate_fft(signal, fs)
freq_lpf, fft_lpf = calculate_fft(signal, fs_downsampled)
freq_no_lpf, fft_no_lpf = calculate_fft(signal_downsampled_without_lpf, fs_downsampled)
def calculate_thd(freq, fft_data, f_signal, harmonic_count=9):
# 基本周波数のインデックス
fundamental_idx = np.argmin(np.abs(freq - f_signal))
fundamental_amplitude = fft_data[fundamental_idx]
# 高調波成分の振幅を取得
harmonic_amplitudes = []
for i in range(2, harmonic_count + 1):
harmonic_freq = f_signal * i
harmonic_idx = np.argmin(np.abs(freq - harmonic_freq))
harmonic_amplitudes.append(fft_data[harmonic_idx])
# THDを計算
harmonic_power_sum = np.sum(np.array(harmonic_amplitudes) ** 2)
thd = np.sqrt(harmonic_power_sum) / fundamental_amplitude * 100
return fundamental_amplitude, harmonic_amplitudes, thd
# THDの計算を実行
fundamental_amplitude, harmonic_amplitudes, thd = calculate_thd(freq_original, fft_original, f_signal)
# 結果表示
print(f"基本周波数の振幅: {fundamental_amplitude:.5f}")
print(f"高調波の振幅: {harmonic_amplitudes}")
print(f"全高調波歪 (THD): {thd:.2f}%")
# Plot the results
#plt.figure(figsize=(12, 8))
mb = 2000
# Original signal
plt.subplot(3, 1, 1)
plt.plot(t[::fs//mb], signal[::fs//mb], label="Original Signal", linewidth=1)
plt.scatter(t[::fs//mb], signal[::fs//mb], color='red', s=10, label="Sampling Points")
plt.title("Original Signal with Harmonic Noise")
plt.xlabel("Time (s)")
plt.ylabel("Amplitude")
plt.legend()
# Downsampled signal with LPF
plt.subplot(3, 1, 2)
plt.plot(t_downsampled[::fs//mb], signal_downsampled_with_lpf[::fs//mb], label="Downsampled with LPF", linewidth=1)
plt.scatter(t_downsampled[::fs//mb], signal_downsampled_with_lpf[::fs//mb], color='red', s=10, label="Sampling Points")
plt.title("Downsampled Signal with LPF")
plt.xlabel("Time (s)")
plt.ylabel("Amplitude")
plt.legend()
# Downsampled signal without LPF
plt.subplot(3, 1, 3)
plt.plot(t_downsampled[::fs//mb], signal_downsampled_without_lpf[::fs//mb], label="Downsampled without LPF", linewidth=1)
plt.scatter(t_downsampled[::fs//mb], signal_downsampled_without_lpf[::fs//mb], color='red', s=10, label="Sampling Points")
plt.title("Downsampled Signal without LPF")
plt.xlabel("Time (s)")
plt.ylabel("Amplitude")
plt.legend()
#plt.tight_layout()
#plt.show()
# Plot all signals in one figure
plt.figure(figsize=(12, 8))
# Plot time-domain signals
plt.subplot(2, 1, 1)
plt.plot(t, signal, label="Original Signal", linewidth=1)
plt.scatter(t[::fs//mb], signal[::fs//mb], color='red', s=10, label="Original Sampling Points")
plt.plot(t_downsampled, signal_downsampled_with_lpf, label="Downsampled with LPF", linewidth=1, linestyle='--')
plt.scatter(t_downsampled, signal_downsampled_with_lpf, color='green', s=10, label="LPF Sampling Points")
plt.plot(t_downsampled, signal_downsampled_without_lpf, label="Downsampled without LPF", linewidth=1, linestyle=':')
plt.scatter(t_downsampled, signal_downsampled_without_lpf, color='blue', s=10, label="No LPF Sampling Points")
plt.title("Comparison of Original and Downsampled Signals")
plt.xlabel("Time (s)")
plt.ylabel("Amplitude")
plt.legend()
# Plot FFT results
plt.subplot(2, 1, 2)
plt.plot(freq_original, fft_original, label="Original Signal FFT", linewidth=1)
#plt.plot(freq_lpf, fft_lpf, label="Downsampled with LPF FFT", linewidth=1, linestyle='--')
#plt.plot(freq_no_lpf, fft_no_lpf, label="Downsampled without LPF FFT", linewidth=1, linestyle=':')
plt.title("FFT Comparison")
plt.xlabel("Frequency (Hz)")
plt.ylabel("Amplitude")
plt.legend()
plt.tight_layout()
plt.show()
■ADC割り込み
/*----------------------------------------------------------------------------*/
/*
* @fn _AD1CH0Interrupt(void)
* @brief ADC1完了割り込み
* @detail
*/
/*----------------------------------------------------------------------------*/
void __attribute__((interrupt, no_auto_psv)) _AD1CH0Interrupt(void)
{
/*----------------------------------------------------------------------------*/
/* AD1とAD2側のADC完了待ち合わせ*/
/*----------------------------------------------------------------------------*/
while(AD1STATbits.CH0RDY == 0){;}
while(AD2STATbits.CH0RDY == 0){;}
/*----------------------------------------------------------------------------*/
/* リサンプリング処理*/
/*----------------------------------------------------------------------------*/
Performing_Resampling();
IFS4bits.AD1CH0IF = 0u;
IFS5bits.AD2CH0IF = 0u;
}
■リサンプリング処理
/*----------------------------------------------------------------------------*/
/*
* @fn Performing_Resampling(void)
* @brief サンプリングデータのリサンプリングを実行します
* @detail ADオーバーサンプリング(64倍)完了割り込みで呼ばれます。
*/
/*----------------------------------------------------------------------------*/
void Performing_Resampling(void)
{
static uint16_t FFT_Cnt = 0;
static fractional fra_Voltage = 0;
static fractional fra_Current = 0;
/*----------------------------------------------------------------------------*/
/* サンプリング*/
/*----------------------------------------------------------------------------*/
fra_Voltage = (AD1CH0DATA - VOLTAGE_OFFSET) * VOLTAGE_CALC_GAIN;
fra_Current = (AD2CH0DATA - CURRENT_OFFSET) * CURRENT_CALC_GAIN;
/*----------------------------------------------------------------------------*/
/* 間引きの実行*/
/*----------------------------------------------------------------------------*/
if (((g_Val.Freq.Status.Bits.is60Hz == 1u) &&
(g_Val.Resampling.u8_Decimate_Cnt > 3u)) ||
((g_Val.Freq.Status.Bits.is50Hz == 1u) &&
(g_Val.Resampling.u8_Decimate_Cnt > 4u)))
{
//実効値算出時の飽和を防止
g_Val.Voltage.fra_ResampleData[g_Val.Resampling.u16_Cnt +
g_Val.Resampling.u16_Buffer_Cnt] = fra_Voltage ;
g_Val.Current.fra_ResampleData[g_Val.Resampling.u16_Cnt +
g_Val.Resampling.u16_Buffer_Cnt] = fra_Current;
//2サンプルに1回FFT用のデータとして処理
if ( FFT_Cnt == 0u)
{
g_Val.Voltage.fra_FFT_Input[g_Val.Resampling.u16_FFT] =
fra_Voltage << 3u;
g_Val.Current.fra_FFT_Input[g_Val.Resampling.u16_FFT] =
fra_Current << 3u;
//カウンタ処理
if (g_Val.Resampling.u16_FFT < RESAMPLE_NUM)
{
g_Val.Resampling.u16_FFT ++;
}
FFT_Cnt = 1u;
}
else
{
FFT_Cnt = 0u;
}
//カウンタ処理
if (g_Val.Resampling.u16_Cnt < RESAMPLE_NUM)
{
g_Val.Resampling.u16_Cnt ++;
}
g_Val.Resampling.u8_Decimate_Cnt = 0u;
}
else
{
g_Val.Resampling.u8_Decimate_Cnt ++;
}
}
1-4.周波数検知部
周波数は、入力電圧をコンパレータで比較することによって算出されます。コンパレータの出力信号をSCCPのインプットキャプチャ機能を使用して、信号の周期を測定します。次に、この測定した周期を基に周波数変換を行い、商用周波数を求めます。
コンパレータの出力はリマッパブルピンを介し出力します。

この方式は一見良さそうに思えたものの、

立ち上がりエッジを拡大して観察するとチャタリングが発生(fig.3-1-4-3)し、インプットキャプチャで正確な値を測定することができませんでした。主な原因として、コンパレータは非常に高速な応答を示す一方で、入力信号が60Hzと非常に低速であることが挙げられます。
対策として、コンパレータのヒステリシスを45mV、デジタルフィルタを8サイクル分に設定することでグリッジの発生頻度を減らすことはできましたが、問題の根本的な解決には至りませんでした。

そこで以下の記事を参考にCLCモジュール2つとLPRCクロックを用い、デバウンス処理を行いました。
以下にデバウンス処理ブロック図を示します。電圧入力をコンパレータで比較し、CLC2とCLC3でデバウンス処理を実施します。このロジックはLPRCのクロックを2つのフリップフロップのクロック入力とし32kHz後に同一入力であれば、その信号を出力します。
出力信号の周期をSCCP2を用いて検出し、周波数へ変換します。

デバウンス処理を実行することで、チャタリングの除去ができました。

2.演算部
2-1.高速sqrt演算
SQRT(平方根計算)は処理負荷が高いことで知られていますが、世の中には精度よりも速度を重視した独自のSQRT算出アルゴリズムが存在します。今回は、32ビット浮動小数点数の平方根の逆数を高速に近似計算するアルゴリズムをご紹介します。
本ロジックは以下のサイトを参考にしました。完全にはアルゴリズムを理解していませんが、この方法を実際に動作させてみたところ、標準のsqrt関数に比べて約3倍高速でありながら、十分な精度を得ることができました。
/*----------------------------------------------------------------------------*/
/*
* @fn Calculate_Fast_Sqrt(float x)
* @brief 高速SQRT
* @param x 値
* @return xの平方根
* @detail Newton-Raphson法でSQRTを求める
* @detail 実行時間 約130usec
*/
/*----------------------------------------------------------------------------*/
float Calculate_Fast_Sqrt(float x)
{
/*----------------------------------------------------------------------------*/
/*ローカル変数の定義と初期化*/
/*----------------------------------------------------------------------------*/
long i;
float x2, y;
const float threehalfs = 1.5F;
/*----------------------------------------------------------------------------*/
/* 計算 */
/*----------------------------------------------------------------------------*/
x2 = x * 0.5F;
y = x;
i = * ( long * ) &y;
i = 0x5f3759df - ( i >> 1 );
y = * ( float * ) &i;
y = y * ( threehalfs - ( x2 * y * y ) );
y = y * ( threehalfs - ( x2 * y * y ) );
return (y * x);
}
2-2.実効値の算出
実効値は以下の式で求めます。
M は 256、u(t)とi(t)はデシメーションによって算出された各サンプルです。

dsPICにはDSP命令のVectorPower関数がありますので、これを利用します。
/*----------------------------------------------------------------------------*/
/*
* @fn Calculate_RMS(fractional *srcV,uint16_t Num)
* @brief 実効値の算出
* @param *Value 算出用バッファ
* @param Num 長さ
* @return 実効値
* @detail 実行時間 約2usec
*/
/*----------------------------------------------------------------------------*/
fractional Calculate_RMS(fractional *srcV, uint16_t Num )
{
/*----------------------------------------------------------------------------*/
/*ローカル変数の定義と初期化*/
/*----------------------------------------------------------------------------*/
fractional RMS = 0;
float Result;
/*----------------------------------------------------------------------------*/
/* 実効値の算出*/
/*----------------------------------------------------------------------------*/
if (Num > 0)
{
//要素のパワーの和を計算
Result = Fract2Float(VectorPower(Num, srcV));
//実効値を求める
//sqrtf (176cycles)
//RMS = Float2Fract(sqrtf(Result)) ;
//Calculate_Fast_Sqrt (61cycles)
RMS = Float2Fract((Calculate_Fast_Sqrt(Result))) ;
}
return RMS;
}
2-3.最大値の算出
最大値は VectorMax関数で求めます。
ただし、VectorMaxは参照渡しで関数内で引数の値を変更します。今回は生値を演算後送信する必要があるので、CPU処理にて求めます。
/*----------------------------------------------------------------------------*/
/*
* @fn Calculate_Max_Value(fractional *srcV,uint16_t Num)
* @brief 最大値の算出
* @param *Value 算出用バッファ
* @param Num 長さ
* @detail 実行時間 約7.74uusec
*/
/*----------------------------------------------------------------------------*/
fractional Calculate_Max_Value(fractional *srcV,uint16_t Num)
{
/*----------------------------------------------------------------------------*/
/*ローカル変数の定義と初期化*/
/*----------------------------------------------------------------------------*/
fractional Max_Value = 0;
int index;
/*----------------------------------------------------------------------------*/
/*最大値を求める*/
/*----------------------------------------------------------------------------*/
#if 0
if (Num > 0)
{
Max_Value = VectorMax(Num, srcV, &index) ;
}
#else
for (int i = 0 ; i < Num ; i ++)
{
if (*srcV > Max_Value){Max_Value = *srcV;}
*srcV ++;
}
#endif
return Max_Value;
}
2-4.DC成分の算出
DC成分は取得値の総和で求めます。

単純な1次元の総和を求めるDSP命令が無かったため、固定データ(=1)で埋められた配列を作成しVectorDotProduct関数を用いて算出しました。
/*----------------------------------------------------------------------------*/
/*
* @fn Calculate_DC(fractional *srcV,uint16_t Num)
* @brief DC値の算出
* @param *Value 算出用バッファ
* @param Num 長さ
* @return DC値
* @detail 実行時間 約1.5μsec
*/
/*----------------------------------------------------------------------------*/
fractional Calculate_DC(fractional *srcV, uint16_t Num )
{
/*----------------------------------------------------------------------------*/
/*ローカル変数の定義と初期化*/
/*----------------------------------------------------------------------------*/
int64_t DC = 0;
fractional Result = 0;
/*----------------------------------------------------------------------------*/
/* 実効値の算出*/
/*----------------------------------------------------------------------------*/
#if 1
if (Num > 0)
{
Result =( VectorDotProduct(Num,srcV,g_Val.fra_FixedData) / Num);
}
#else
if (Num > 0)
{
for (int i= 0 ;i < Num; i++)
{
DC += (*srcV ++);
}
Result = ((DC / Num) );
}
#endif
return Result;
}
2-5.有効電力の算出
有効電力は電圧と電流の各サンプルの積和で求めます。

dsPICにはDSP命令のVectorDotProduct関数がありますので、これを利用します。
/*----------------------------------------------------------------------------*/
/*
* @fn Calculate_Active_Power(fractional*Voltage,fractional *Current,uint16_t Num)
* @brief 有効電力の算出
* @param *Voltage 電圧
* @param *Current 電流
* @param Num 長さ
* @detail 実行時間 約2.78μsec
*/
/*----------------------------------------------------------------------------*/
fractional Calculate_Active_Power(fractional *Voltage,fractional *Current,uint16_t Num)
{
/*----------------------------------------------------------------------------*/
/*ローカル変数の定義と初期化*/
/*----------------------------------------------------------------------------*/
fractional Active_Power = 0;
/*----------------------------------------------------------------------------*/
/*有効電力の算出*/
/*----------------------------------------------------------------------------*/
if (Num > 0)
{
Active_Power = VectorDotProduct(Num,Voltage,Current) / Num << 3u;
}
return Active_Power;
}
2-6.皮相電力の算出
皮相電力は電圧実効値と電流実効値の積で求めます。

dsPICにはビルトインDSP関数の__builtin_mullss関数がありますので、これを利用します。
電圧実効値は2^18倍、電流実効値は2^22倍で表現しており、皮相電力は2^12倍で表現しますので、
演算結果を28bit(18+22-12)右シフトしています。
/*----------------------------------------------------------------------------*/
/*
* @fn Calculate_Apparent_Power(fractional Voltage_RMS,fractional Current_RMS)
* @brief 皮相電力の算出
* @param Voltage_RMS 電圧実効値
* @param Current_RMS 電流実効値
* @return 皮相電力
* @detail 実行時間 約10nsec
*/
/*----------------------------------------------------------------------------*/
fractional Calculate_Apparent_Power(fractional Voltage_RMS,fractional Current_RMS)
{
/*----------------------------------------------------------------------------*/
/*ローカル変数の定義と初期化*/
/*----------------------------------------------------------------------------*/
fractional Apparent_Power;
/*----------------------------------------------------------------------------*/
/*皮相電力の算出*/
/*----------------------------------------------------------------------------*/
Apparent_Power = __builtin_mulss(Voltage_RMS,Current_RMS) >> 28u;
return Apparent_Power;
}
2-7.無効電力の算出
無効電力は以下の式で求めます。

/*----------------------------------------------------------------------------*/
/*
* @fn Calculate_Reactive_Power(fractional Active_Power,fractional Apparent_Power)
* @brief 無効電力の算出
* @param Active_Power 有効電力
* @param Apparent_Power 皮相電力
* @return 無効電力
* @detail 実行時間 約610nsec
*/
/*----------------------------------------------------------------------------*/
fractional Calculate_Reactive_Power(fractional Active_Power,fractional Apparent_Power)
{
/*----------------------------------------------------------------------------*/
/*ローカル変数の定義と初期化*/
/*----------------------------------------------------------------------------*/
fractional ReActive_Power = 0;
float AcPow = Fract2Float(Active_Power);
float ApPow = Fract2Float(Apparent_Power);
/*----------------------------------------------------------------------------*/
/*無効電力の算出*/
/*----------------------------------------------------------------------------*/
if (AcPow < ApPow)
{
float ts = Calculate_Fast_Sqrt((ApPow * ApPow) - (AcPow * AcPow));
ReActive_Power = Float2Fract(ts);
}
return ReActive_Power;
}
2-8.力率
力率(PF)は有効電力と皮相電力の割合で求めます。

/*----------------------------------------------------------------------------*/
/*
* @fn Calculate_Power_Factor(fractional Active_Power,fractional Apparent_Power)
* @brief 力率の算出
* @param Active_Power 有効電力
* @param Apparent_Power 皮相電力
* @return 力率
* @detail 実行時間 約220nsec
*/
/*----------------------------------------------------------------------------*/
float Calculate_Power_Factor(fractional Active_Power,fractional Apparent_Power)
{
/*----------------------------------------------------------------------------*/
/*ローカル変数の定義と初期化*/
/*----------------------------------------------------------------------------*/
float Power_Factor = 0;
/*----------------------------------------------------------------------------*/
/*無効電力の算出*/
/*----------------------------------------------------------------------------*/
if (Apparent_Power > 0)
{
if (Active_Power < Apparent_Power)
{
Power_Factor = Fract2Float(Active_Power) /
Fract2Float( Apparent_Power);
}
else
{
Power_Factor = 1.0;
}
};
return Power_Factor;
}
2-9.高調波成分
高調波成分はFFTを使用し求めます。
詳しくは以下のリンクを参照してください。
dsPIC33AKシリーズ⑰ DSPの使い方2(FFT) – ぴくおの電子工作的な何かWP
15.36kHz(60Hzの場合)でリサンプリングされたデータを1/2に間引きを行い、FFTへの入力データとします。これによりFFTでは2サイクル分の解析を行うことになります。(fs =7.68kHz@60Hz)
■初期化関数
void Init_FFT(void)
{
/*------------------------------------------------------------------------*/
/*回転因子の初期化*/
/*------------------------------------------------------------------------*/
tFactor = TwidFactorInit(
FFT_POWEROF2,
(fractcomplex *)twiddleFactors,
0);
/*------------------------------------------------------------------------*/
/*窓関数の初期化*/
/*------------------------------------------------------------------------*/
HanningInit(FFT_POINT,fr_Window);
}
■実行関数
/*----------------------------------------------------------------------------*/
/*
* @fn Calculate_FFT(fractional *srcV,fractional *dstV)
* @brief FFTの実行
* @param *srcV ソースベクトル
* @param *dstV1 ディスティネーションベクトル(パワースペクトル)
* @detail 実行時間 約131usec
*/
/*----------------------------------------------------------------------------*/
void Calculate_FFT(fractional *srcV,fractional *dstV1)
{
/*----------------------------------------------------------------------------*/
/*ローカル変数の定義と初期化*/
/*----------------------------------------------------------------------------*/
int i;
/*------------------------------------------------------------------------*/
//①窓関数の適用(DSP関数)
*------------------------------------------------------------------------*/
VectorWindow(FFT_POINT,
(fractional *)fr_AfterWindow,
(fractional *)srcV,
(fractional *)fr_Window);
/*------------------------------------------------------------------------*/
//②入力データを複素数の実数に代入
/*------------------------------------------------------------------------*/
for (i = 0; i < FFT_POINT;i++)
{
frc_FftData[i].real = fr_AfterWindow[i];
frc_FftData[i].imag = 0;
}
/*------------------------------------------------------------------------*/
//③FFTの実行(DSP関数)
/*------------------------------------------------------------------------*/
FFTComplexIP(
FFT_POWEROF2,
(fractcomplex *)frc_FftData,
(fractcomplex *)twiddleFactors,
0);
/*------------------------------------------------------------------------*/
// ④ビット逆順ソート(DSP関数)
/*------------------------------------------------------------------------*/
BitReverseComplex( FFT_POWEROF2, frc_FftData );
/*------------------------------------------------------------------------*/
// ⑤複素数の絶対値自乗を取る(DSP関数)
/*------------------------------------------------------------------------*/
SquareMagnitudeCplx( FFT_POINT >>1 , frc_FftData, dstV1);
/*------------------------------------------------------------------------*/
// ⑥正規化
/*------------------------------------------------------------------------*/
for (i =0 ; i < (FFT_POINT >> 1); i++)
{
dstV1[i] = dstV1[i] << 5u;
}
}
2-10.THD(高調波歪み)
THD(Total Harmonic Distortion)は、全高調波歪みを意味し、信号の歪み具合を示す指標の一つです。特に電気・電子機器や音響機器などで、信号の純度や品質を評価する際に重要です。
THDは、ある信号に含まれる高調波成分の強度が、基準となる基本波(通常は最も低い周波数の成分)に対してどれだけの割合を占めるかを示します。高調波とは、基準となる周波数の整数倍の周波数を持つ成分です。例えば、基準となる周波数が50Hzの場合、100Hz、150Hz、200Hzなどが高調波になります。
THDは通常、次の式で定義されます:

ここで:
- H1:基本波(fundamental frequency)の振幅
- H2,H3,…:2次、3次、…といった高調波成分の振幅
dsPICのFFT結果は振幅ではなくパワーで出力される(SquareMagnitudeCplx関数)ため、THD(全高調波歪)は以下の式で定義されます。

ここで:
- P1:基本波(fundamental frequency)のパワー
- P2,P3,…:2次、3次、…といった高調波成分のパワー
/*----------------------------------------------------------------------------*/
/*
* @fn Calculate_THD(fractional *Value,uint16_t Num)
* @brief THDの算出
* @param *Value FFT結果バッファ
* @param Num FFT長さ
* @return 高調波歪
* @detail 実行時間 約2usec
*/
/*----------------------------------------------------------------------------*/
float Calculate_THD(fractional *FFTResult,uint16_t Num)
{
/*----------------------------------------------------------------------------*/
/*ローカル変数の定義と初期化*/
/*----------------------------------------------------------------------------*/
fractional fra_THD = 0;
float f4_THD = 0.0f;
float harmonic_power_sum = 0;
float fundamental_power = 0;
//基本周波数 * FFTポイント / サンプリング周波数
const int StepIndex = 60 * FFT_POINT / 7680;
/*----------------------------------------------------------------------------*/
/*2次以降の高調波を算出*/
/*----------------------------------------------------------------------------*/
for (int i = StepIndex * 2 ; i < Num / 2; i += StepIndex)
{
fra_THD += FFTResult[i] >> 5u;
}
harmonic_power_sum = Fract2Float(fra_THD);
/*----------------------------------------------------------------------------*/
/* 高調波歪の計算*/
/*----------------------------------------------------------------------------*/
fundamental_power = Fract2Float(FFTResult[StepIndex] >> 5u);
if (FFTResult[StepIndex] != 0)
{
f4_THD = Calculate_Fast_Sqrt(harmonic_power_sum / fundamental_power) * 100.0;
}
2-11.位相差
位相差を定義する方法としては以下のような方法が有ります。電圧や電流が純粋な正弦波であれば特に問題はありません。しかし、高調波が含まれる場合は、その定義が複雑になります。
1.基本波成分に基づく位相差
FFT(Fast Fourier Transform)などにより電圧、電流の基本波(50Hzまたは60Hz)成分を抽出し、その位相差から計算
2.全波形に基づく位相差
電圧、電流のゼロクロス点のタイミングから検出
3.電力を基にした定義
有効電力(P)、皮相電力(S)、無効電力(Q)、有効電力と皮相電力の位相角(Θ)の関係は以下となります。また力率はcosΘで表されます。つまり、力率から位相角を求める事は可能です。

ただし力率(Power Factor, PF)の値だけでは、進み位相か遅れ位相かを直接判断することはできません。力率は位相角のコサイン値で計算されるため、進み位相と遅れ位相で同じ値になるからです。
大手メーカの電力計では別途符号を求めている模様です。
今回は、1.基本波成分に基づく位相差について解説します。
複素平面では、複素数は実部 a と虚部 b を用いて表現されます。

振幅Aは原点から点zまでの距離を示し、位相角Θは複素数 zと実軸との間の角度を表します。

今回の解析では以下の条件を使用します:
・FFTポイント数:256ポイント
・サンプリング周波数:7.68kHz(60Hz信号の解析用)
これにより、周波数解像度は以下のように計算されます。

このため、FFT結果配列の2番目の要素(30Hz × 2 = 60Hz)が基本波に対応します。この要素の 実部 と 虚部 を用いて、基本波の位相角を計算します。
次に、FFTから求めた電圧ベクトルと電流ベクトルの位相角を比較することで、両ベクトルの角度差 Θ′を求めます。この角度差に基づき、進み位相か遅れ位相かを判断します。

/*----------------------------------------------------------------------------*/
/*
* @fn Calculate_Phase(fractcomplex *frc_Data)
* @brief 位相角の算出
* @param frc_Data FFT結果バッファ
* @return 位相
* @detail Time = 2usec
*/
/*----------------------------------------------------------------------------*/
float Calculate_Phase(fractcomplex *frc_Data)
{
float r,i;
frc_Data->real = frc_FftData[2].real;
frc_Data->imag = frc_FftData[2].imag;
r = Fract2Float(frc_Data->real);
i = Fract2Float(frc_Data->imag);
if (i != 0 )
{
return(atan( r/i));
}
else
{
return 0;
}
}
2-12.各種演算処理呼び出し
_CCT2Interrupt割り込みはゼロクロスの度に呼ばれます。
ここでサンプリングバッファの切り替えと各種演算処理を呼び出しています。
/*----------------------------------------------------------------------------*/
/*
* @fn _CCT2Interrupt(void)
* @brief インプットキャプチャ完了割り込み
* @detail
*/
/*----------------------------------------------------------------------------*/
void __attribute__((interrupt, no_auto_psv)) _CCT2Interrupt(void)
{
/*----------------------------------------------------------------------------*/
/* 各種演算処理*/
/*----------------------------------------------------------------------------*/
if (CCP2STATbits.ICBNE == 1u)
{
//周期値の読込
g_Val.Freq.u32_Period = CCP2BUF;
if (g_Val.Freq.u32_Period < 10000000)
{
//サンプリングバッファの入れ替え
if (g_Val.Resampling.u16_Buffer_Cnt == 0)
{
g_Val.Resampling.u16_Buffer_Cnt = 256;
}
else
{
g_Val.Resampling.u16_Buffer_Cnt = 0;
g_Val.Resampling.u16_FFT = 0;
}
g_Val.Status.Bits.PowerFailure = 0u;
//演算実効
Performing_Calculations();
}
}
IFS1bits.CCT2IF = 0u;
}
Performing_Calculations()関数は実効値や電力、FFTなど各種演算を呼び出します。
/*----------------------------------------------------------------------------*/
/*
* @fn Performing_Calculations(void)
* @brief 各種演算を実効します
* @detail ゼロクロス検出ごとに呼ばれます。
*/
/*----------------------------------------------------------------------------*/
void Performing_Calculations(void)
{
/*----------------------------------------------------------------------------*/
/* サンプリングカウンタのクリア */
/*----------------------------------------------------------------------------*/
g_Val.Resampling.u16_Cnt = 0;
/*----------------------------------------------------------------------------*/
/* 入力周波数の算出(1回/1周期)*/
/*----------------------------------------------------------------------------*/
g_Val.Freq.u32_Freq_Sum += Calculate_Input_Freqency(g_Val.Freq.u32_Period);
/*----------------------------------------------------------------------------*/
/* 電圧・電流実効値の算出(1回/1周期)*/
/*----------------------------------------------------------------------------*/
g_Val.Voltage.fra_RMS_Sum_X18 +=
Calculate_RMS(g_Val.Voltage.fra_ResampleData +
g_Val.Resampling.u16_Buffer_Cnt,RESAMPLE_NUM);
g_Val.Current.fra_RMS_Sum_X22 +=
Calculate_RMS(g_Val.Current.fra_ResampleData +
g_Val.Resampling.u16_Buffer_Cnt,RESAMPLE_NUM);
/*----------------------------------------------------------------------------*/
/* 電圧・電流DC値の算出(1回/1周期)*/
/*----------------------------------------------------------------------------*/
g_Val.Voltage.fra_DC_Sum_X18 +=
Calculate_DC(g_Val.Voltage.fra_ResampleData +
g_Val.Resampling.u16_Buffer_Cnt,RESAMPLE_NUM);
g_Val.Current.fra_DC_Sum_X22
Calculate_DC(g_Val.Current.fra_ResampleData +
g_Val.Resampling.u16_Buffer_Cnt,RESAMPLE_NUM);
/*----------------------------------------------------------------------------*/
/* 有効電力の算出(1回/1周期)*/
/*----------------------------------------------------------------------------*/
g_Val.Power.fra_Active_Power_Sum_X12 += Calculate_Active_Power(
g_Val.Voltage.fra_ResampleData + g_Val.Resampling.u16_Buffer_Cnt,
g_Val.Current.fra_ResampleData + g_Val.Resampling.u16_Buffer_Cnt,
RESAMPLE_NUM
);
/*----------------------------------------------------------------------------*/
/* 電圧・電流最大値の算出(1回/1周期)*/
/*----------------------------------------------------------------------------*/
g_Val.Voltage.fra_Peak_Sum_X18 +=
Calculate_Max_Value(g_Val.Voltage.fra_ResampleData +
g_Val.Resampling.u16_Buffer_Cnt,RESAMPLE_NUM);
g_Val.Current.fra_Peak_Sum_X22 +=
Calculate_Max_Value(g_Val.Current.fra_ResampleData +
g_Val.Resampling.u16_Buffer_Cnt,RESAMPLE_NUM);
/*----------------------------------------------------------------------------*/
/* 1回/N周期演算*/
/*----------------------------------------------------------------------------*/
if (( ++ g_Val.Average.u8_Cnt >= g_Val.Average.u8_Num) &&
(g_Val.Average.u8_Num > 0u))
{
/*------------------------------------------------------------------*/
/* 各平均値の算出*/
/*------------------------------------------------------------------*/
g_Val.Freq.u16_Freq_Ave = g_Val.Freq.u32_Freq_Sum / g_Val.Average.u8_Num;
g_Val.Voltage.fra_RMS_Ave_X18 = g_Val.Voltage.fra_RMS_Sum_X18 /
g_Val.Average.u8_Num ;
g_Val.Current.fra_RMS_Ave_X22 = g_Val.Current.fra_RMS_Sum_X22 /
g_Val.Average.u8_Num;
g_Val.Voltage.fra_DC_Ave_X18 = g_Val.Voltage.fra_DC_Sum_X18 /
g_Val.Average.u8_Num ;
g_Val.Current.fra_DC_Ave_X22 = g_Val.Current.fra_DC_Sum_X22 /
g_Val.Average.u8_Num;
g_Val.Voltage.fra_Peak_Ave_X18 = g_Val.Voltage.fra_Peak_Sum_X18 /
g_Val.Average.u8_Num;
g_Val.Current.fra_Peak_Ave_X22 = g_Val.Current.fra_Peak_Sum_X22 /
g_Val.Average.u8_Num;
g_Val.Power.fra_Active_Power_Ave_X12 =
g_Val.Power.fra_Active_Power_Sum_X12/ g_Val.Average.u8_Num;
/*------------------------------------------------------------------*/
/* 各積算値と平均化カウンタのクリア*/
/*------------------------------------------------------------------*/
g_Val.Average.u8_Cnt = 0;
g_Val.Freq.u32_Freq_Sum = 0u;
g_Val.Voltage.fra_RMS_Sum_X18 = 0u;
g_Val.Current.fra_RMS_Sum_X22 = 0u;
g_Val.Voltage.fra_DC_Sum_X18 = 0u;
g_Val.Current.fra_DC_Sum_X22 = 0u;
g_Val.Voltage.fra_Peak_Sum_X18 = 0u;
g_Val.Current.fra_Peak_Sum_X22 = 0u;
g_Val.Power.fra_Active_Power_Sum_X12 = 0u;
/*------------------------------------------------------------------*/
/* 周波数の判定*/
/*------------------------------------------------------------------*/
if ((g_Val.Freq.u16_Freq_Ave > 4500) &&
(g_Val.Freq.u16_Freq_Ave < 5500))
{g_Val.Freq.Status.Bits.is60Hz = 1u;}
else {g_Val.Freq.Status.Bits.is50Hz = 0u;}
if ((g_Val.Freq.u16_Freq_Ave > 5500) &&
(g_Val.Freq.u16_Freq_Ave < 6500))
{g_Val.Freq.Status.Bits.is60Hz = 1u;}
else {g_Val.Freq.Status.Bits.is60Hz = 0u;}
/*------------------------------------------------------------------*/
/* 皮相電力の算出*/
/*------------------------------------------------------------------*/
g_Val.Power.fra_Apparent_Power_Ave_X12 =
Calculate_Apparent_Power(
g_Val.Voltage.fra_RMS_Ave_X18,
g_Val.Current.fra_RMS_Ave_X22
);
/*------------------------------------------------------------------*/
/* 無効電力の算出*/
/*------------------------------------------------------------------*/
g_Val.Power.fra_Reactive_Power_Ave_X12 =
Calculate_Reactive_Power(
g_Val.Power.fra_Active_Power_Ave_X12,
g_Val.Power.fra_Apparent_Power_Ave_X12
);
/*------------------------------------------------------------------*/
/* 力率の算出*/
/*------------------------------------------------------------------*/
g_Val.Power.f4_Power_Factor_Ave =
Calculate_Power_Factor(
g_Val.Power.fra_Active_Power_Ave_X12,
g_Val.Power.fra_Apparent_Power_Ave_X12
);
/*------------------------------------------------------------------*/
/* 電圧・電流FFTの実行(1回/N周期)*/
/*------------------------------------------------------------------*/
Calculate_FFT(g_Val.Voltage.fra_FFT_Input,g_Val.Voltage.fra_FFT_Result);
g_Val.Voltage.f4_Phase_Ave = Calculate_Phase(&g_Val.Voltage.frc_FftData);
Calculate_FFT(g_Val.Current.fra_FFT_Input,g_Val.Current.fra_FFT_Result);
g_Val.Current.f4_Phase_Ave = Calculate_Phase(&g_Val.Current.frc_FftData);
/*------------------------------------------------------------------*/
/* 電圧・電流THDの算出*/
/*------------------------------------------------------------------*/
g_Val.Voltage.f4_THD_Ave =
Calculate_THD(g_Val.Voltage.fra_FFT_Result,RESAMPLE_NUM);
g_Val.Current.f4_THD_Ave =
Calculate_THD(g_Val.Current.fra_FFT_Result,RESAMPLE_NUM);
g_Val.Control.Bits.SendMinState = SEND_STATE_MINOR_MAKE_DATA;
g_Val.Control.Bits.SendMajState = SEND_STATE_MAJOR_NORMAL;
}
}
3.データ送信部
1.概要
各種演算した値や生のデータをPCに送信します。また平均化回数や送信の開始/停止信号をPCから受信します。
演算完了フラグが立つと、各種演算値 → 電圧瞬時値 → 電流瞬時値→ 電圧高調波→ 電流高調波データ の順にDMA転送を行います。全部で2kByte以上のデータを921kbpsで送るため、1回の送信には約25msec必要となります。
2.ブロック図

(fig.3-3-1)UART通信部ブロック図
3.通信プロトコル
| 項目 | 値 |
| ボーレート[bps] | 921600 |
| パリティ | 無し |
| ストップビット[bit] | 1 |
| フォーマット | バイナリ |
| 変換文字 | DLC(0x10) |
| エスケープバイト | ESC(0x1B) |
| デリミタ | CR(0x0C)+LF(0x0A) |
今回は通信量削減を目的とし、データはASCII形式ではなくバイナリ形式で送信します。しかし、バイナリデータにおいては、特定の値(例えば0x0Cや0x0A)がデータの区切り(デリミタ)として使われることがあるため、これらの値がデータに含まれると、区切りを正しく認識できなくなります。そこで、エスケープシーケンス処理を行うことで、この問題に対処します。
具体的な処理方法は以下の通りです:
- 送信データに特定の値が含まれる場合
送信データ内に0x0C(フォームフィード)、0x0A(ラインフィード)、または0x1B(エスケープ)が含まれている場合、それらの文字がデータの区切りとして誤認識されるのを防ぐため、エスケープ処理を行います。 - エスケープ処理の内容
該当する文字(0x0C、0x0A、0x1B)の前に、制御文字であるESC(0x1B)を付加し、その文字自体は0x10との排他的論理和(XOR)を取ります。これにより、データ内でこれらの文字が区切りとして使われることなく、正しく送受信できます。 - 受信側での処理
受信側では、データ列にESCが含まれている場合、その直後の文字と0x10を排他的論理和で処理します。これにより、元の文字(0x0C、0x0A、0x1B)が復元されます。
このエスケープシーケンス処理により、データ内で区切り文字を誤認識することなく、正確にデータを受信することができ、データの整合性を保つことができます。
//送信データ
0x3F,0x24,0x0C,0x55,0x32,0x87,CR(0x0C),LF(0x0A) // 変換前
↓
0x3F,0x24,ESC(0x1B),0x1C(0x0C XOR 0x10),0x55,0x32,0x87,CR(0x0C),LF(0x0A) // 変換後
4.ソースコード
/*----------------------------------------------------------------------------*/
/*
* @fn Uart1_Init_DMABuff(void)
* @brief DMA送信用初期化
* @detail
*/
/*----------------------------------------------------------------------------*/
void Uart1_Init_DMABuff(void)
{
DMACONbits.ON = 1U; //DMA有効
DMACONbits.PRIORITY = 0U; //Fixed
DMACONbits.SIDL = 0UL;
DMABUF = 0x00000000u;
DMALOW = 0x4000U;
DMAHIGH = 0x7FFFU;
DMA0SRC = (unsigned int)& s_Uart1.Tx.u8_Buff;
DMA0DST = (unsigned int)& U1TXREG;
DMA0CLR = 0x00000000U;
DMA0SET = 0x00000000U;
DMA0INV = 0x00000000U;
DMA0MSK = 0x00000000U;
DMA0PAT = 0x00000000U;
DMA0CH = 0x00000000U;
DMA0CHbits.SIZE = 0U; //8Bit
DMA0CHbits.TRMODE = 0U; //ワンショット
DMA0CHbits.DAMODE = 0U; //
DMA0CHbits.SAMODE = 1U; //Increment
DMA0SEL = 0x10U; //U1TXDone
DMA0STAT = 0x00000000U;
s_Uart1.Tx.u16_Head = 0;
}
/*----------------------------------------------------------------------------*/
/*
* @fn Uart1_Add_DMABuff(fractional val, char Length)
* @brief Uart送信データの作成
* @param val データソース
* @param type タイプ
* @return
* @detail
*/
/*----------------------------------------------------------------------------*/
void Uart1_Add_DMABuff(fractional val, char Length)
{
/*----------------------------------------------------------------------------*/
/* ローカル変数の定義と初期化*/
/*----------------------------------------------------------------------------*/
uint8_t Byte;
/*----------------------------------------------------------------------------*/
/* エスケープ処理の実効*/
/*----------------------------------------------------------------------------*/
for (int i = 3; i >= ( 4 - Length );i--)
{
Byte = ((val) >> (i * 8)) & 0xFF;
if ((Byte == ASCII_LF)||(Byte == ASCII_CR)||(Byte == ASCII_ESC))
{
s_Uart1.Tx.u8_Buff[s_Uart1.Tx.u16_Head] = ASCII_ESC;
if (s_Uart1.Tx.u16_Head < (UART1_TX_BUFF - 1u))
{
s_Uart1.Tx.u16_Head ++;
}
s_Uart1.Tx.u8_Buff[s_Uart1.Tx.u16_Head] = Byte ^ ASCII_DLE;
}
else
{
s_Uart1.Tx.u8_Buff[s_Uart1.Tx.u16_Head] = Byte;
}
if (s_Uart1.Tx.u16_Head < (UART1_TX_BUFF - 1u))
{
s_Uart1.Tx.u16_Head ++;
}
}
}
/*----------------------------------------------------------------------------*/
/*
* @fn Uart1_Add_DMABuff_Del(void)
* @brief デリミタの追加
* @detail
*/
/*----------------------------------------------------------------------------*/
void Uart1_Add_DMABuff_Del(void)
{
/*-----------------------------------------------------*/
/* デリミタの追加 */
/*-----------------------------------------------------*/
s_Uart1.Tx.u8_Buff[s_Uart1.Tx.u16_Head] = ASCII_CR;
if (s_Uart1.Tx.u16_Head < (UART1_TX_BUFF - 1u))
{
s_Uart1.Tx.u16_Head ++;
}
s_Uart1.Tx.u8_Buff[s_Uart1.Tx.u16_Head] = ASCII_LF;
if (s_Uart1.Tx.u16_Head < (UART1_TX_BUFF - 1u))
{
s_Uart1.Tx.u16_Head ++;
}
}
/*----------------------------------------------------------------------------*/
/*
* @fn Uart1_Make_NormalData(void)
* @brief 演算結果の作成と送信
* @detail
*/
/*----------------------------------------------------------------------------*/
void Uart1_Make_NormalData(void)
{
/*-----------------------------------------------------*/
/* 送信データの生成 */
/*-----------------------------------------------------*/
Uart1_Init_DMABuff();
Uart1_Add_DMABuff(((uint32_t)'D') << 24,1u);
//周波数(Freq)
Uart1_Add_DMABuff(((uint32_t)g_Val.Freq.u16_Freq_Ave) << 16,2u);
//電圧実効値(Urms)
Uart1_Add_DMABuff(g_Val.Voltage.fra_RMS_Ave_X18,4u);
//電圧ピーク(Upeak)
Uart1_Add_DMABuff(g_Val.Voltage.fra_Peak_Ave_X18,4u);
//電圧DC(Udc)
Uart1_Add_DMABuff(g_Val.Voltage.fra_DC_Ave_X18,4u);
//電圧THD(Uthd)
Uart1_Add_DMABuff(Float2Fract(g_Val.Voltage.f4_THD_Ave),4u);
//電圧位相(Uphase)
Uart1_Add_DMABuff(Float2Fract(g_Val.Voltage.f4_THD_Ave),4u);
//電流実効値(Irms)
Uart1_Add_DMABuff(g_Val.Current.fra_RMS_Ave_X22,4u);
//電流ピーク(Ipeak)
Uart1_Add_DMABuff(g_Val.Current.fra_Peak_Ave_X22,4u);
//電流DC(Idc)
Uart1_Add_DMABuff(g_Val.Current.fra_DC_Ave_X22,4u);
//電流THD(Ithd)
Uart1_Add_DMABuff(Float2Fract(g_Val.Current.f4_THD_Ave),4u);
//有効電力(P)
Uart1_Add_DMABuff(g_Val.Power.fra_Active_Power_Ave_X12,4u);
//皮相電力(S))
Uart1_Add_DMABuff(g_Val.Power.fra_Apparent_Power_Ave_X12,4u);
//無効電力(Q))
Uart1_Add_DMABuff(g_Val.Power.fra_Reactive_Power_Ave_X12,4u);
//積算電力(WP))
Uart1_Add_DMABuff(g_Val.Power.i64_Active_TotalPower_Ave_X5,4u);
//力率(λ)
Uart1_Add_DMABuff(Float2Fract(g_Val.Power.f4_Power_Factor_Ave),4u);
Uart1_Add_DMABuff_Del();
/*-----------------------------------------------------*/
/* DMA送信の開始 */
/*-----------------------------------------------------*/
DMA0CNT = s_Uart1.Tx.u16_Head;
DMA0CHbits.CHEN = 1U;
DMA0CHbits.CHREQ = 1U;
}
/*----------------------------------------------------------------------------*/
/*
* @fn Uart1_Make_VLiveData(void)
* @brief 電圧生データの作成と送信
* @detail
*/
/*----------------------------------------------------------------------------*/
void Uart1_Make_VLiveData(void)
{
/*-----------------------------------------------------*/
/* 送信データの生成 */
/*-----------------------------------------------------*/
Uart1_Init_DMABuff();
Uart1_Add_DMABuff(((uint32_t)'V') << 24,1u);
for (int i = 0;i < 256;i++)
{
Uart1_Add_DMABuff(g_Val.Voltage.fra_ResampleData[i] ,2u);
}
Uart1_Add_DMABuff_Del();
/*-----------------------------------------------------*/
/* DMA送信の開始 */
/*-----------------------------------------------------*/
DMA0CNT = s_Uart1.Tx.u16_Head;
DMA0CHbits.CHEN = 1U;
DMA0CHbits.CHREQ = 1U;
}
/*----------------------------------------------------------------------------*/
/*
* @fn Uart1_Make_ILiveData(void)
* @brief 電流生データの作成と送信
* @detail
*/
/*----------------------------------------------------------------------------*/
void Uart1_Make_ILiveData(void)
{
/*-----------------------------------------------------*/
/* 送信データの生成 */
/*-----------------------------------------------------*/
Uart1_Init_DMABuff();
Uart1_Add_DMABuff(((uint32_t)'I') << 24,1u);
for (int i = 0;i < 256;i++)
{
Uart1_Add_DMABuff(g_Val.Current.fra_ResampleData[i],2u);
}
Uart1_Add_DMABuff_Del();
/*-----------------------------------------------------*/
/* DMA送信の開始 */
/*-----------------------------------------------------*/
DMA0CNT = s_Uart1.Tx.u16_Head;
DMA0CHbits.CHEN = 1U;
DMA0CHbits.CHREQ = 1U;
}
/*----------------------------------------------------------------------------*/
/*
* @fn Uart1_Make_VFFTData(void)
* @brief 電圧FFTデータの作成と送信
* @detail
*/
/*----------------------------------------------------------------------------*/
void Uart1_Make_VFFTData(void)
{
/*-----------------------------------------------------*/
/* 送信データの生成 */
/*-----------------------------------------------------*/
Uart1_Init_DMABuff();
Uart1_Add_DMABuff(((uint32_t)'v') << 24,1u);
for (int i = 0;i < 128;i++)
{
Uart1_Add_DMABuff(g_Val.Voltage.fra_FFT_Result[i],4u);
}
Uart1_Add_DMABuff_Del();
/*-----------------------------------------------------*/
/* DMA送信の開始 */
/*-----------------------------------------------------*/
DMA0CNT = s_Uart1.Tx.u16_Head;
DMA0CHbits.CHEN = 1U;
DMA0CHbits.CHREQ = 1U;
}
/*----------------------------------------------------------------------------*/
/*
* @fn Uart1_Make_IFFTData(void)
* @brief 電流FFTデータの作成と送信
* @detail
*/
/*----------------------------------------------------------------------------*/
void Uart1_Make_IFFTData(void)
{
/*-----------------------------------------------------*/
/* 送信データの生成 */
/*-----------------------------------------------------*/
Uart1_Init_DMABuff();
Uart1_Add_DMABuff(((uint32_t)'i') << 24,1u);
for (int i = 0;i < 128;i++)
{
Uart1_Add_DMABuff(g_Val.Current.fra_FFT_Result[i],4u);
}
Uart1_Add_DMABuff_Del();
/*-----------------------------------------------------*/
/* データ送信 */
/*-----------------------------------------------------*/
DMA0CNT = s_Uart1.Tx.u16_Head;
DMA0CHbits.CHEN = 1U;
DMA0CHbits.CHREQ = 1U;
}
/*----------------------------------------------------------------------------*/
/*
* @fn Uart1_Make_TX_Data(void)
* @brief Uart送信データの作成
* @detail
*/
/*----------------------------------------------------------------------------*/
void Uart1_Idle_Process(void)
{
/*----------------------------------------------------------------------------*/
/* 文字列の生成 */
/*----------------------------------------------------------------------------*/
if (g_Val.Control.Bits.SendStart == 1u)
{
switch (g_Val.Control.Bits.SendMajState)
{
/*-----------------------------------------------------*/
/* 送信待機 */
/*-----------------------------------------------------*/
case SEND_STATE_MAJOR_WAIT:
break;
/*-----------------------------------------------------*/
/* 通常パラメータの送信 */
/*-----------------------------------------------------*/
case SEND_STATE_MAJOR_NORMAL:
switch (g_Val.Control.Bits.SendMinState)
{
case SEND_STATE_MINOR_MAKE_DATA:
Uart1_Make_NormalData();
DMA0CHbits.CHEN = 1U;
DMA0CHbits.CHREQ = 1U;
g_Val.Control.Bits.SendMinState = SEND_STATE_MINOR_WAIT_DONE;
break;
case SEND_STATE_MINOR_WAIT_DONE:
if (DMA0STATbits.DONE == 1u)
{
DMA0STATbits.DONE = 0u;
g_Val.Control.Bits.SendMinState = SEND_STATE_MINOR_MAKE_DATA;
g_Val.Control.Bits.SendMajState = SEND_STATE_MAJOR_VLIVE;
}
break;
default:
g_Val.Control.Bits.SendMinState = SEND_STATE_MINOR_MAKE_DATA;
break;
}
break;
/*-----------------------------------------------------*/
/* 生波形の送信 */
/*-----------------------------------------------------*/
case SEND_STATE_MAJOR_VLIVE:
switch (g_Val.Control.Bits.SendMinState)
{
case SEND_STATE_MINOR_MAKE_DATA:
Uart1_Make_VLiveData();
g_Val.Control.Bits.SendMinState = SEND_STATE_MINOR_WAIT_DONE;
DMA0CHbits.CHEN = 1U;
DMA0CHbits.CHREQ = 1U;
break;
case SEND_STATE_MINOR_WAIT_DONE:
if (DMA0STATbits.DONE == 1u)
{
DMA0STATbits.DONE = 0u;
g_Val.Control.Bits.SendMinState = SEND_STATE_MINOR_MAKE_DATA;
g_Val.Control.Bits.SendMajState = SEND_STATE_MAJOR_ILIVE;
}
break;
default:
g_Val.Control.Bits.SendMinState = SEND_STATE_MINOR_MAKE_DATA;
break;
}
break;
/*-----------------------------------------------------*/
/* 生波形の送信 */
/*----------------------------------------------------*/
case SEND_STATE_MAJOR_ILIVE:
switch (g_Val.Control.Bits.SendMinState)
{
case SEND_STATE_MINOR_MAKE_DATA:
Uart1_Make_ILiveData();
g_Val.Control.Bits.SendMinState = SEND_STATE_MINOR_WAIT_DONE;
DMA0CHbits.CHEN = 1U;
DMA0CHbits.CHREQ = 1U;
break;
case SEND_STATE_MINOR_WAIT_DONE:
if (DMA0STATbits.DONE == 1u)
{
DMA0STATbits.DONE = 0u;
g_Val.Control.Bits.SendMinState = SEND_STATE_MINOR_MAKE_DATA;
g_Val.Control.Bits.SendMajState = SEND_STATE_MAJOR_VFFT;
}
break;
default:
g_Val.Control.Bits.SendMinState = SEND_STATE_MINOR_MAKE_DATA;
break;
}
break;
/*-----------------------------------------------------*/
/* 電圧高調波の送信 */
/*-----------------------------------------------------*/
case SEND_STATE_MAJOR_VFFT:
switch (g_Val.Control.Bits.SendMinState)
{
case SEND_STATE_MINOR_MAKE_DATA:
Uart1_Make_VFFTData();
g_Val.Control.Bits.SendMinState = SEND_STATE_MINOR_WAIT_DONE;
DMA0CHbits.CHEN = 1U;
DMA0CHbits.CHREQ = 1U;
break;
case SEND_STATE_MINOR_WAIT_DONE:
if (DMA0STATbits.DONE == 1u)
{
DMA0STATbits.DONE = 0u;
g_Val.Control.Bits.SendMinState = SEND_STATE_MINOR_MAKE_DATA;
g_Val.Control.Bits.SendMajState = SEND_STATE_MAJOR_IFFT;
}
break;
default:
g_Val.Control.Bits.SendMinState = SEND_STATE_MINOR_MAKE_DATA;
break;
}
break;
/*-----------------------------------------------------*/
/* 電流高調波の送信 */
/*-----------------------------------------------------*/
case SEND_STATE_MAJOR_IFFT:
switch (g_Val.Control.Bits.SendMinState)
{
case SEND_STATE_MINOR_MAKE_DATA:
Uart1_Make_IFFTData();
g_Val.Control.Bits.SendMinState = SEND_STATE_MINOR_WAIT_DONE;
DMA0CHbits.CHEN = 1U;
DMA0CHbits.CHREQ = 1U;
break;
case SEND_STATE_MINOR_WAIT_DONE:
if (DMA0STATbits.DONE == 1u)
{
DMA0STATbits.DONE = 0u;
g_Val.Control.Bits.SendMinState = SEND_STATE_MINOR_MAKE_DATA;
g_Val.Control.Bits.SendMajState = SEND_STATE_MAJOR_WAIT;
}
break;
default:
g_Val.Control.Bits.SendMinState = SEND_STATE_MINOR_MAKE_DATA;
break;
}
break;
default:
g_Val.Control.Bits.SendMinState = SEND_STATE_MINOR_MAKE_DATA;
g_Val.Control.Bits.SendMajState = SEND_STATE_MAJOR_WAIT;
break;
}
}
else
{
g_Val.Control.Bits.SendMinState = 0u;
g_Val.Control.Bits.SendMajState = 0u;
}
}
/*----------------------------------------------------------------------------*/
/*
* @fn Uart1_Analyze_Rx_Data(void)
* @brief Uart1受信データの解析
* @detail
*/
/*----------------------------------------------------------------------------*/
void Uart1_Analyze_Rx_Data(void)
{
/*----------------------------------------------------------------------------*/
/* ローカル変数の定義と初期化*/
/*----------------------------------------------------------------------------*/
static uint8_t u8_Buff[16];
static uint8_t u8_Head = 0;
static int16_t u16_RxData = 0;
/*----------------------------------------------------------------------------*/
/* 受信データの解析*/
/*----------------------------------------------------------------------------*/
u16_RxData = Uart1_Deque_RX_Data();
if (u16_RxData != -1)
{
if (u16_RxData == ASCII_LF)
{
u8_Head = 0u;
switch (u8_Buff[0])
{
case 'N':
switch (u8_Buff[1])
{
case '0':g_Val.Average.u8_Num = 1u;break;
case '1':g_Val.Average.u8_Num = 2u;break;
case '2':g_Val.Average.u8_Num = 4u;break;
case '3':g_Val.Average.u8_Num = 8u;break;
case '4':g_Val.Average.u8_Num = 16u;break;
case '5':g_Val.Average.u8_Num = 32u;break;
case '6':g_Val.Average.u8_Num = 64u;break;
case '7':g_Val.Average.u8_Num = 128u;break;
default:break;
}
break;
case 'S':g_Val.Control.Bits.SendStart = 1u;break;
case 'E':g_Val.Control.Bits.SendStart = 0u;break;
case 'R':g_Val.Power.i64_Active_TotalPower_Ave_X5 = 0u;break;
default:break;
}
}
else if (u16_RxData == ASCII_CR)
{
u8_Buff[u8_Head] = 0;
u8_Head = 0u;
}
else
{
if (u8_Head < 15)
{
u8_Buff[u8_Head ++] = (uint8_t)u16_RxData;
}
}
}
}
ソースコード
コンフィグレーションファイル、クロック設定ファイルは以下のファイルをインクルードしてください。
■コンフィグレーションファイル (config.h)
■クロック設定ソースファイル(Clock_Driver.c)
■クロックヘッダーファイル(clock_driver.c)
■ソースコード全体
近日公開予定
結果
実行時間
平均化回数を1回に設定し、実行した際の各処理の実行時間を比較した結果を示します。
FFT処理が重いのは避けられない点ですが、ピーク値算出がDSPではなくCPU処理で行われているため、意外と時間がかかっています。
一括してピーク値を求めるのではなく、リサンプリングの度に計算を行う方法の方が、全体の実行時間を短縮できる可能性があります。
| 処理 | 関数 | 実行時間(μsec) | 処理回数 | 実行時間(μsec) |
| 周波数算出 | Calculate_Input_Freqency | 0.13 | 1 | 0.13 |
| 実効値算出 | Calculate_RMS | 1.89 | 2 | 3.78 |
| DC値算出 | Calculate_DC | 1.5 | 2 | 3 |
| ピーク算出 | Calculate_Max_Value | 7.74 | 2 | 15.48 |
| 有効電力算出 | Calculate_Active_Power | 2.78 | 1 | 2.78 |
| 皮相電力算出 | Calculate_Apparent_Power | 0.04 | 1 | 0.04 |
| 無効電力算出 | Calculate_Reactive_Power | 0.54 | 1 | 0.54 |
| 力率算出 | Calculate_Power_Factor | 0.2 | 1 | 0.2 |
| FFT実行 | Calculate_FFT | 132.97 | 2 | 265.94 |
| THD算出 | Calculate_THD | 2.11 | 2 | 4.22 |
| 位相算出 | Calculate_Phase | 0.68 | 2 | 1.36 |
| リサンプリング | Performing_Resampling | 0.1 | 256 | 25.6 |
実行波形
こちらは、平均化回数を2回に設定して動作させた際の波形です。ゼロクロスごとに演算が行われている様子が確認できます。具体的には、1回目のサンプリングでは演算に26μsec、2回目では324μsec、そしてその後のデータ送信には23.3msecを要しています。
この結果から、演算よりもデータ送信に多くの実行時間がかかっていることがわかります。

実行結果表示
市販の電力チェッカーで測定した値と比較したところ、多少の差異は見られたものの、検証する基準がない中で大きなずれは確認されなかったため、妥当と判断しました。
以下は、C#で作成したWindows用アプリケーションで、GrafanaライクなUIに仕上げています。
テストデータ
最初にダミーデータを使用して、正しく演算処理が行われているかを試行しました。

半田ごて
手持ちの温調式半田ごて(60W)で温度設定を変更し、計測を行いました。200℃に設定した場合、力率は73%で、450℃に設定すると電流波形が正弦波に近づき、力率は95%に向上しました。また、電流の高調波歪も約5%に低減され、安定した状態となりました。



電気ストーブ
電気ストーブでは弱・強設定とも、ほぼ歪が無い綺麗な電流で力率も99%を超過していました。



まとめ
以上で、dsPIC33AKを活用したAC負荷アナライザの設計と実装について解説しました。このプロジェクトは、さまざまな応用分野に展開可能です。次回の記事では、PCソフトウェアについて詳しく解説します。
記事についての注意点
本記事は慎重に内容を検討し正確さに努めておりますが、内容に誤りがあったとしても、この記事を参考にして生じた損害等については一切の責任を負いません。

コメント