圏9研究所 工作室

圏9研究所の開発情報資料など

STM32F411 BlackPill CMSIS-DSPによる楽器用チューナー（７）コード

BlackPill STM32

計測用コードのみ記載
１．コード
１）計測
（１）<fft_peak.h>

/**
 Project Name:FFT_F4

 File Name:fft_peak.h

 Description:
 USB speaker with isochronous transfer
 Generation Information :
 Device			:  STM32F411CEU6
 Version			:  01
 Created on		: 2023/12/20
 */
/**
 Copyright (c) 2023 luke24e-hbid

 This software is released under the MIT License.
 http://opensource.org/licenses/mit-license.php
 */

/* Define to prevent recursive inclusion -------------------------------------*/
#ifndef INC_FFT_PEAK_H_
#define INC_FFT_PEAK_H_

/* Includes ------------------------------------------------------------------*/
#include "main.h"
#define ARM_MATH_CM4
#include "arm_math.h"
#include "arm_const_structs.h"
#include "arm_common_tables.h"

/* Private defines -----------------------------------------------------------*/

/* Function prototype declaration ------------------------------------------*/
//	Read FFT Data
void read_data_sin(float32_t fin);

// FFT peak detect
float32_t rfft_spline(uint32_t nfft);

#endif /* INC_FFT_PEAK_H_ */
（２）<fft_peak.c>
/**
 Project Name:FFT_F4

 File Name:fft_peak.c

 Description:
 USB speaker with isochronous transfer
 Generation Information :
 Device			:  STM32F411CEU6
 Version			: 01
 Created on		: 2023/12/20
 */
/**
 Copyright (c) 2023 luke24e-hbid

 This software is released under the MIT License.
 http://opensource.org/licenses/mit-license.php
 */
/* Includes ------------------------------------------------------------------*/
// #include "main.h"
#include "fft_peak.h"
#include <stdlib.h>

/* Private macro -------------------------------------------------------------*/
// RFFT macro
#define FS1 16000                                 // sampling frequency 1 : 16000[Hz]
#define FT 4                                      // sampling period [Hz]
#define NS1 (FS1 / FT)                            // number of read sample data : 4000
#define NFFT 4096                                 // number of RFFT data
#define NPOWER (NFFT / 2)                         // number of RFFT power data
#define FRESO1 ((float32_t)FS1 / (float32_t)NFFT) // frequency resolution 1 : 16000/4096[Hz]

// measurement macro
#define FKMIN 55.0   // minimum frequency
#define FKMAX 880.0  // maximun frequency
#define PKMIN 200000 // minimum power
#define PKK 2        // minimum power rate

// Spline macro
#define NSPS 12         // number of spline data
#define NSPH 100        // number of spline interpolation
#define PSPH (NSPH / 2) //  magnification of spline interpolation

/* Private variables ---------------------------------------------------------*/
// global variable
int16_t data_sample[NFFT * 2]; // sampling data

/* Private user code ---------------------------------------------------------*/
//
// FFT peak detect
float32_t rfft_spline(uint32_t nfft)
{
  // RFFT parameter
  arm_rfft_fast_instance_f32 fft;
  float32_t fft_input[NFFT];
  float32_t fft_output[NFFT];
  float32_t fft_power[NFFT / 2];
  uint8_t ifftFlag = 0;
  float32_t maxValue;
  uint32_t maxIndex;

  // arm_spline_init_f32  parameter
  arm_spline_instance_f32 S;
  arm_spline_type type = ARM_SPLINE_NATURAL;
  float32_t x[NSPS];
  float32_t y[NSPS];
  uint32_t n = NSPS;
  float32_t coeffs[3 * (NSPS - 1)];
  float32_t tempBuffer[2 * (NSPS - 1)];

  // arm_spline_f32   parameter
  float32_t xq[NSPH];
  float32_t pDst[NSPH];
  uint32_t blockSize = NSPH;
  float32_t maxValue_s;
  uint32_t maxIndex_s;

  // measurement
  float32_t avg_power; // average of FFT power
  float32_t fm1;       // dominant frequency 1
  uint32_t fm1Index;
  uint32_t fm1Index_01;
  uint32_t fm1Index_02;
  uint32_t fm1_multi;
  float32_t fm1_01;
  float32_t fm1_02;
  float32_t xn;
  float32_t fs2; // sampling frequency 2
  float32_t freso2;
  uint16_t nk2;
  float32_t fm2;   // dominant frequency 2
  float32_t fmspl; // spline interpolation frequency

  uint16_t i;

  // Set FFT data 1
  for (i = 0; i < NS1; i++)
  {
    fft_input[i] = (float32_t)data_sample[i + nfft];
  }
  for (i = NS1; i < NFFT; i++)
  {
    fft_input[i] = 0;
  }

  // analyze signal 1
  arm_rfft_fast_init_f32(&fft, NFFT);
  arm_rfft_fast_f32(&fft, fft_input, fft_output, ifftFlag);
  arm_cmplx_mag_f32(fft_output, fft_power, NFFT / 2);
  // fill fft power with zeros
  fft_power[0] = 0.0;
  arm_mean_f32(fft_power, NPOWER, &avg_power);

  for (i = 0; i < NPOWER; i++)
  {
    if (i < (uint16_t)(FKMIN / FRESO1))
      fft_power[i] = 0;
    else if (i > (uint16_t)(FKMAX / FRESO1))
      fft_power[i] = 0;
    else
    {
      if (fft_power[i] < PKMIN)
        fft_power[i] = 0;
      else if (fft_power[i] < (avg_power * PKK))
        fft_power[i] = 0;
    }
  }
  // find dominant frequency 1-1
  arm_max_f32(fft_power, NPOWER, &maxValue, &maxIndex);
  fm1Index_01 = maxIndex;

  if (fm1Index_01 == 0)
  {
    return (0.0);
  }

  fm1_01 = (float32_t)fm1Index_01 * FRESO1;
  // fill fft power with zeros
  i = (uint16_t)fm1Index_01;
  do
  {
    fft_power[i++] = 0;
  } while (fft_power[i] != 0);
  i = (uint16_t)fm1Index_01;
  do
  {
    fft_power[i--] = 0;
  } while (fft_power[i] != 0);
  for (i = (uint16_t)(fm1Index_01 * 2 + fm1Index_01 / 2); i < NPOWER; i++)
  {
    fft_power[i] = 0;
  }

  // find dominant frequency 1-2
  arm_max_f32(fft_power, NPOWER, &maxValue, &maxIndex);
  fm1Index_02 = maxIndex;
  fm1_02 = (float32_t)fm1Index_02 * FRESO1;

  // find dominate frequency 1
  if (fm1Index_01 > fm1Index_02)
  {
    uint32_t t = fm1Index_01;
    fm1Index_01 = fm1Index_02;
    fm1Index_02 = t;
  }
  fm1_multi = 1;
  if (fm1Index_01 == 0)
    fm1Index = fm1Index_02;
  else
  {
    if (abs((int32_t)fm1Index_01 * 2 - (int32_t)fm1Index_02) < 2)
      fm1Index = fm1Index_01;
    else
    {
      fm1Index = fm1Index_01;
      fm1_multi = (uint32_t)((float32_t)fm1Index_01 / (float32_t)(fm1Index_02 - fm1Index_01) + 0.5);
    }
  }
  fm1 = (float32_t)fm1Index * FRESO1;

  // calculate sampling frequency 2
  xn = (log(fm1 * 16.0) - log(125.0)) / log(2.0) + 0.5;
  fs2 = 125 * (1 << (uint16_t)xn);
  freso2 = fs2 / NFFT;

  // Set FFT data 2
  nk2 = FS1 / (uint16_t)fs2;
  for (i = 0; i < (uint16_t)(fs2 / 4); i++)
  {
    fft_input[i] = 0;
    for (uint16_t j = 0; j < nk2; j++)
    {
      fft_input[i] = fft_input[i] + (float32_t)data_sample[i * nk2 + j + nfft];
    }
  }
  for (i = (uint16_t)(fs2 / 4); i < NFFT; i++)
  {
    fft_input[i] = 0;
  }
  // analyze signal 2
  arm_rfft_fast_init_f32(&fft, NFFT);
  arm_rfft_fast_f32(&fft, fft_input, fft_output, ifftFlag);
  arm_cmplx_mag_f32(fft_output, fft_power, NFFT / 2);
  // find dominant frequency 2
  maxValue = 0;
  maxIndex = 0;
  for (i = (fm1Index - 2) * nk2; i < (fm1Index + 2) * nk2; i++)
  {
    if (fft_power[i] > maxValue)
    {
      maxValue = fft_power[i];
      maxIndex = i;
    }
  }
  fm2 = (float32_t)maxIndex * freso2;

  // spline
  for (i = 0; i < NSPS; i++)
  {
    x[i] = (float32_t)(i + maxIndex - (NSPS / 2));
    y[i] = fft_power[i + maxIndex - (NSPS / 2)];
  }
  for (i = 0; i < NSPH; i++)
  {
    xq[i] = (float32_t)i / (float32_t)PSPH + (float32_t)(maxIndex - 1);
  }
  arm_spline_init_f32(&S, type, x, y, n, coeffs, tempBuffer);
  arm_spline_f32(&S, xq, pDst, blockSize);
  // find spline dominant frequency
  pDst[0] = 0.0;
  arm_max_f32(pDst, 100, &maxValue_s, &maxIndex_s);
  fmspl = xq[maxIndex_s] * freso2 / (float32_t)fm1_multi;

  // printf("%.3f\t%.0f\t%.3f\t%.3f\r\n", (double)fm1, fs2,(double)fm2, (double)fmspl);

  return (fmspl);
}
２）呼び出し例
　nfft : サンプリングデータ先頭
　fmspl : 計測値[Hz]
/* Private variables ---------------------------------------------------------*/
uint32_t g_nfft_adc = 0xffff;

/* Private user code ---------------------------------------------------------*/
// ADC DMA harf
void HAL_ADC_ConvHalfCpltCallback(ADC_HandleTypeDef *const hadc)
{
  g_nfft_adc = 0;
}
// ADC DMA full
void HAL_ADC_ConvCpltCallback(ADC_HandleTypeDef *const hadc)
{
  g_nfft_adc = 4000;
}

//	RFFT, Display
//
void display_tune(uint32_t nfft)
{
  float32_t fmspl;  // display frequency
  int16_t oct_cent; // octave cent
  int16_t cent1200; // cent
  int16_t note;     // note number
  int16_t cent;     // display cent +/- 50

  fmspl = rfft_spline(nfft);
  oct_cent = (uint16_t)(1200.0 * log(fmspl / 27.5) / log(2.0) + 50.5);
  cent1200 = oct_cent % 1200;
  note = (cent1200) / 100;
  cent = (cent1200) % 100 - 50;
}

//
//  main

extern int16_t data_sample[]; // sampling data
HAL_ADC_Start_DMA(&hadc1, (uint32_t *)data_sample, 8000);
HAL_TIM_Base_Start(&htim3);

while (1)
{
  if (g_nfft_adc != 0xffff)
  {
    display_tune(g_nfft_adc);
    g_nfft_adc = 0xffff;
  }
}

２．補足解説
１）コードについて
・コードは FFT計算とスプライン補間による基音周波数計測と基音周波数のセント値変換関数を記載
・呼び出し例としてAD変換DMA転送によるデータを計測を記載
・DMA転送はダブルバッファで構成
・表示用コードは音名表示に使用するフォントデータ量が大きいため記載していない

２）表示イメージ
（１）リソースの想定

リソース

デバイス

用途

I2C Master

OLED

音名、周波数等表示

LED 2port

2色LED 赤・黄緑

チューニング状態表示

　緑色：±5セント内

　橙色：±5セント外

　赤色：ADC入力レベル小

　　　　周波数測定範囲外

PWM 2ch

アナログメータ（電流計 両振れ）

セント値表示

 

（２）表示イメージ：YAMAHA TD-90
　音名、周波数等はOLEDで表示

jp.yamaha.com

 

おしまい

STM32F411 BlackPill CMSIS-DSPによる楽器用チューナー（６）定数、変数、処理

BlackPill STM32

１．定数一覧
１）FFT

シンボル

設定値

単位

内容

FS1

16000 = FSREAD / NK1

Hz

FFT1次サンプリング周波数

　マイク：160000　USB：48000/3

FT

4 (250msec)

Hz

1/サンプリング時間

NS1

4000 = FS1 / FT

 

読込サンプル数

NFFT

4096

 

FFTサンプル数

NPOWER

2048 = NFFT / 2

 

FFT Power サンプル数

FRESO1

FS1 / NFFT

Hz

FFT1次周波数ピッチ

FKMIN

55.0

Hz

計測下限周波数

FKMAX

880.0

Hz

計測上限周波数

PKMIN

1000000 * NK1

 

計測下限Power（無音レベル）

PKK

2

 

ピーク前後除去係数（Power平均値に掛ける値）

 

２）スプライン

シンボル

設定値

単位

内容

NSPS

12

 

サンプル数 FFTピーク±6相当

NSPH

100

 

補間数

PSPH

50 = NSPS / 2

 

補間倍率 = 1/補間ピッチ = NSPH / 2

 

２．変数一覧
１）FFT

変数名

G/L

内容

uint16_t data_sample[4096 * 2]

G

読込サンプルデータ 4000個 / 250msec

arm_rfft_fast_instance_f32 fft;

L

RFFT計算

float32_t fft_input[NFFT];

L

float32_t fft_output[NFFT];

L

float32_t fft_power[NFFT / 2];

L

uint8_t ifftFlag = 0;

L

float32_t maxValue;

L

arm_max_f32 引数

uint32_t maxIndex;

L

 

２）スプライン

変数名

G/L

内容

arm_spline_instance_f32 S;

L

スプライン係数

arm_spline_type type = ARM_SPLINE_NATURAL;

L

float32_t x[NSPS];

L

float32_t y[NSPS];

L

uint32_t n = NSPS;

L

float32_t coeffs[3 * (NSPS - 1)];

L

float32_t tempBuffer[2 * (NSPS - 1)];

L

float32_t xq[NSPH];

L

スプライン補間

float32_t pDst[NSPH];

L

uint32_t blockSize = NSPH;

L

float32_t maxValue_s;

L

arm_max_f32 引数

uint32_t maxIndex_S;

L

 

３）計測計算

処理

変数名

G/L

内容

FFT1次計算

float32_t avg_power

L

FFT power平均値

float32_t fm1_01, fm1_02

L

FFT1次powerピーク01, 02

uint32_t fm1Index_01, fm1Index_02

L

fm1_01, fm1_02の fft_powerインデックス

float32_t fm1

L

FFT1次計測周波数 fm1=基音*fm1_multi

uint32_t fm1Index

L

fm1の fft_powerインデックス

uint32_t fm1_multi

L

fm1基音の倍数

FFT２次準備

float32_t xn

L

サンプリング変換係数

float32_t fs2

L

FFT2次サンプリング周波数

float32_t freso2

L

FFT2次周波数ピッチ

uint32_t nk2

L

fs1 / fs2 サンプリング周波数比

FFT2次計算

float32_t fm2

L

FFT2次計測周波数

スプライン

float32_t fmspl

L

スプライン補間計測周波数

fmspl=スプライン補間値 / fm1_multi

セント値計算

int16_t oct_cent

L

オクターブセント値 1[oct] = 1200[cent]

int16_t cent1200

L

セント値：0 - 1199

int16_t note

L

音階値：0 - 11

int16_t cent

L

表示セント値：±50

 

３．計測手順
１）FFT １次

処理

手順

サンプルデータコピー

fft_input[0-3999] = sample_data[0-3999]

不足分ゼロ埋め

fft_input[4000-4095] = {0}

FFT計算

arm_rfft_fast_init_f32(&fft, NFFT);

arm_rfft_fast_f32(&fft, fft_input, fft_output, ifftFlag);

arm_cmplx_mag_f32(fft_output, fft_power, NFFT / 2);

0Hz値をゼロに

fft_power[0] = 0

fft_power 平均値算出

arm_mean_f32

ゼロ埋め

・計測範囲外　55Hz(55 / FRESO1)<, >880Hz(880 / FRESO1)

・fft_power < 平均値 * 2(PKK)

・fft_power < 無音レベル(PKMIN)

ピーク検出　1回目

 fm1_01

 fm1Index_01

arm_max_f32

ピーク前後 fft_power==0 までゼロ埋め

ピークIndex 2.5倍以上ゼロ埋め

ピーク検出　2回目

 fm1_02

 fm1Index_02

arm_max_f32

ピーク値1,2ソート

Index小さい順

計測周波数決定

 fm1

 fm1Index

fm1_multi

if power[fm1Index_01]==0　　　エラー

fm1_multi = 1

if power[fm1Index_02]==0　　　fm1Index=fm1Index_01

else

　if fm1Index_01 * 2 == fm1Index_02±1　　fm1Index=fm1Index_01

　else　

　　fm1Index = fm1Index_01

　　fm1_multi = fm1Index_01 / (fm1Index_02 - fm1Index_01)

fm1 =  (float32_t)fm1Index * FRESO1

 

２）FFT ２次

処理

手順

２次サンプリング周波数

と周波数比計算

　fs2, nk2, freso2

xn = (log(fm1 * 16) – log125)/log2

fs2 = 125 * 2 ^ xn

nk2 = fs1 / fs2

freso2 = fs2 / NFFT：4096;

サンプルデータコピー

　fft_input

nk2 個サンプル加算して１サンプル値とする

 データ数：fs2 / 4：FT

 残りゼロ埋め

FFT計算

FFT 1次と同じ

ピーク検出

　fm2

　fm2Index

fm1Index±2 相当の間のピーク fm2Index を探す

 (fm1Index – 2)*nk2 < 検索 < (fm1Index + 2)*nk2

 fm2 = (float32_t)fm2Index * freso2 / fm1_multi;

 

３）スプライン補間

処理

手順

スプライン係数計算データ

　x

　y

インデックス範囲：maxIndex ±6

データ数：12 = NSPS

スプライン補間 x値

　xq[]

補間区間数 = 100：NSPH

 fm2Index±１/2区間を補間倍率 50：PSPH = 100：NSPH/2で補間値算出

スプライン計算

arm_spline_init_f32(&S, type, x, y, n, coeffs, tempBuffer);

arm_spline_f32(&S, xq, pDst, blockSize);

ピーク検出

　fmspl

・arm_max 補間値yピーク値抽出

・計測周波数算出

　fmspl = xq[maxIndex_s] * freso2 / fm1_multi

 

４．表示計算
１）オクターブ値、セント値計算
（１）基本計算式
　基準周波数 = 27.5[Hz] A0
　オクターブセント：oct_cent = 1200 * log(fmspl / 27.5) / log(2)
　セント値        ：cent1200 = oct_cent % 1200
　音階値        ：note = cent1200 / 100
　表示セント値     ：cent = cent1200 % 100

（２）表示用計算式
　表示セント値範囲を±50 符号付にする

　基準周波数 = 27.5[Hz] A0
　オクターブセント：oct_cent = (int16_t)(1200 * log(fmspl / 27.5) / log(2) + 50.5)
　セント値        ：cent1200 = oct_cent % 1200
　音階値        ：note = cent1200 / 100
　表示セント値     ：cent = cent1200 % 100 - 50

５．補足：マイクアンプとFFT計算
１）マイクアンプ
（１）仕様
・購入先        秋月
・品名        高感度マイクアンプキット
・販売コード     105757
・種別        コンデンサーマイク
・電源電圧        2.5V〜5V
・周波数        50Hz〜6kHz
・出力        0V〜VCC
・DCオフセット    約VCC/2

（２）電源ノイズ対策
　VCC/2を仮想グランドとしていることもあり電源ノイズの影響が大きく未対策では計測レンジを小さくしなければならない
　マイクアンプの電源にはマイクロインダクター47mHを使ったパイ型フィルターを付けた

２）FFT計算データ処理
（１）DCオフセット対策
　DCオフセット相当量は0HzポイントのFFTパワー値に加算される
　対策
　・AD変換値は加工せずFFT計算し0Hzポイントのパワー値を無視する
　　　今回0Hzポイントは計測範囲外のためこの方法を採用
　・各AD変換値からDCオフセット相当量を引き算した後FFT計算する
　・各AD変換値からAD変換値の平均値を引き算した後FFT計算する

（２）サンプリング値不足対策
　サンプリング周波数とFFTサンプル数の不一致対策
　・不足するFFTサンプル値をゼロで埋める
　・ゼロで埋めた分は0HzポイントのFFTパワー値に加算されるため0Hzポイントのパワー値を無視する

続く

STM32F411 BlackPill CMSIS-DSPによる楽器用チューナー（４）CMSIS DSP実装

BlackPill STM32 開発環境

1．概要
　arm CMSIS-DSPライブラリの使用方法
　STM32CubeIDE内蔵のCMSIS-DSPはバージョンが古くスプライン補間関数が使えなかったので最新版を導入する

２．先生とお手本
１）先生
（１）CMSIS-DSP導入方法について

aptechlabs.com

ioloa.com

 

（２）スプライン補間について

qiita.com

qiita.com

https://www.rajgunesh.com/resources/downloads/numerical/cubicsplineinterpol.pdf

２）お手本
（１）FFT

klangstrom.dennisppaul.de

 

（２）スプライン補間

bbs.elecfans.com

 

３．STM32CubeIDEへの導入
１）ライブラリー入手
　STM32CubeIDE内蔵のCMSIS-DSPはバージョンが古いので最新版を導入する

（１）最新のライブラリー

github.com

　　.a .lib ファイルのある最新版CMSIS 5.7.0 をdownloadする

 

２）STM32CubeIDEへの導入
　cortexM4を搭載するSTM32F411での導入例

（１）ライブラリーファイルインポート
　ARM.CMSIS.5.7.0/CMSIS/DSP/Lib/ARM/arm_cortexM4lf_math.lib
　ARM.CMSIS.5.7.0/CMSIS/DSP/Lib/GCC/libarm_cortexM4lf_math.a

（２）インクルードファイルインポート
　ARM.CMSIS.5.7.0/CMSIS/DSP/Include
　　arm_math.h
　　arm_const_structs.h
　　arm_common_tables.h

（３）ライブラリー設定
　Librarys                arm_cortexM4lf_math
　Library search path    workspace を設定

（４）設定ハードコピー

 

4．コード例
１）三角関数

/* USER CODE BEGIN Includes */
#define ARM_MATH_CM4
#include "arm_math.h"
#include "arm_const_structs.h"
/* USER CODE END Includes */

float armsinVal = arm_sin_f32((float32_t)num*2*M_PI/255);
float armcosVal = arm_cos_f32((float32_t)num*2*M_PI/255);
２）Real FFT f32
（１）仕様
　データ型        float_f32 固定小数点32bit
　サンプル長    [32, 64, 128, ..., 4096] 
　入力バッファデータは保持されない
（２）コード例
　ピーク周波数の算出
　・fft_input データのRFFTを算出する
　・各ポイントのpower値を計算する
　・power最大値とそのインデックスを取り出す
/* Includes ------------------------------------------------------------------*/
#define ARM_MATH_CM4
#include "arm_math.h"
#include "arm_const_structs.h"
#include "arm_common_tables.h"

/* Private variables ---------------------------------------------------------*/
#define SAMPLE_BUFFER_LENGTH        4096
#define SAMPLE_BUFFER_LENGTH_HALF   (SAMPLE_BUFFER_LENGTH/2)
#define SAMPLING_RATE               16384

float fft_input[SAMPLE_BUFFER_LENGTH];
float fft_output[SAMPLE_BUFFER_LENGTH];
float fft_power[SAMPLE_BUFFER_LENGTH_HALF];
uint8_t ifftFlag = 0;

/* analyze signal */
arm_rfft_fast_instance_f32 fft;
arm_rfft_fast_init_f32(&fft, SAMPLE_BUFFER_LENGTH);
arm_rfft_fast_f32(&fft, fft_input, fft_output, ifftFlag);
arm_cmplx_mag_f32(fft_output, fft_power, SAMPLE_BUFFER_LENGTH_HALF);

/* find dominant frequency */
float32_t maxValue;
uint32_t maxIndex;

arm_max_f32(fft_power, SAMPLE_BUFFER_LENGTH_HALF, &maxValue, &maxIndex);
３）３次スプライン補間
　・12ポイントのデータから３次スプライン関数の係数を求める
　　ポイント終端は直線　ARM_SPLINE_NATURAL
　・xqに対する補間値pDstを算出する
　　パラメータ計算：12ポイント
　　補間配列：100セット
　　xq[]値を設定するコードは記載なし
/* Includes ------------------------------------------------------------------*/
#define ARM_MATH_CM4
#include "arm_math.h"
#include "arm_const_structs.h"
#include "arm_common_tables.h"

/* Private variables ---------------------------------------------------------*/
//arm_spline_init_f32  parameter
static arm_spline_instance_f32 S;
static arm_spline_type type = ARM_SPLINE_NATURAL;
static float32_t x[12];
static float32_t y[12];
static uint32_t n = 12;
static float32_t coeffs[3 * (12 - 1)];
static float32_t tempBuffer[2 * (12 - 1)];

//arm_spline_f32   parameter
static float32_t xq[100];
static float32_t pDst[100];
static uint32_t blockSize = 100;

/* analyze signal */
	arm_spline_init_f32(&S, type, x, y, n, coeffs, tempBuffer);
	arm_spline_f32(&S, xq, pDst, blockSize);

５．補足解説
１）STM32CubeIDEのCMSIS-DSPライブラリー
・Version 1.12.0にスプライン補間関数は実装されていない
・STM32CubeIDEの最新版に実装されているかどうかは未確認

２）CMSIS-DSP最新版
　先生のリンク先に記載されている通り
　・STM32CubeIDEプロジェクトへライブラリーをお手軽に導入するには .a .lib ファイルが必要
　・ライブラリーリンク先のCMSIS 5.7.0より新しいバージョンでは .a .lib ファイルを作る必要がある
　 ・ARM.CMSIS.5.7.0.pack をダウンロードしたら
　　右クリック　このアプリケーションで開く > アーカイブユーティリティを選択実行

２）RFFT処理時間
　実測値 5msec〜10msec程度（96MHz 4096ポイント）

３）スプライン補間で係数を求めるために必要なサンプル数
　先生として記載したリンク「時系列逐次スプライン補間の実践（後編）」の値をお手本にした
　　先頭から５，６点離れれば、以後その区分曲線はほとんど変化しない

続く

STM32F411 BlackPill CMSIS-DSPによる楽器用チューナー（３）トレースサポート実装

BlackPill STM32 開発環境

コーディング準備

１．トレースサポート実装
１）概要
　FFTの計算結果をモニターしたいので実装する
　デバッガのトレースサポート機能を使ってstdoutをデバッガのコンソールに接続してprintf関数等で状態を表示させる
　UARTを使って外部ターミナルに表示させることもできるがST LINKだけで動作するので実装してみる

２）実装方法
（１）お手本

lujji.github.io

deepbluembedded.com

moons.link

（２）手順
　詳細はお手本に詳しく記載されているので概要のみ
　・ mini ST-LINK V2 からSWO線を引き出しPB3に接続する
　・STM32CubeMXでトレース設定する
　・stdoutをリダイレクトしprintf関数等をコーディングする
　・デバッガのコンソールを設定する
　・デバッグ実行

２．補足解説
１）コード
　最後のお手本のコードがシンプル

２） mini ST-LINK V2 改造
　 難関の引き出し線はんだ付け手順について補足
　・ICピンピッチ0.5mmのためはんだ付けが難しい
　・導線　導線径がピンピッチの半分程度の単線を使う
　　　　　AWG30 ラッピングワイヤーとかUEW 0.26 ポリウレタン銅線とか
　・導線端全周漏れなくはんだ上げする
　　ポリウレタン銅線の被覆除去はムラが出ないようにはんだごての温度は高め
　・導線を基板に仮固定する
　　作業中に導線端がICピンから浮いたりズレたりしないように固定
　・フラックスを塗る
　・薄くハンダをのせたはんだごてを当てる
　　フラックスを塗っておけばブリッジしにくいので隣接するピンに当てても大丈夫
　　ブリッジした場合はこて先のハンダを除去しフラックスを塗り直して再度はんだごてを当てる
　・アルコールで洗浄する
　　筆とか歯ブラシを使うとやりやすい
　・ルーペではんだの状態を確認する
　・導線を接着剤等で基板に固定する

続く

STM32F411 BlackPill CMSIS-DSPによる楽器用チューナー（２）仕様・課題・処理

BlackPill STM32

目標仕様と課題及び処理方法について

１．目標仕様

項目

目標値

備考

用途

弦楽器の調弦

弦を使うピアノも含む

計測値

音名、セント値

オクターブ値は実用上必要ないので参考値とする

計測範囲

55[Hz]〜880[Hz]

下限がベース下限周波数より高いが下記により選定

・使用するマイクの帯域下限が50Hz

・倍音で計測しオクターブ値は参考

計測精度

±1[セント]

セント値=1200・log2(計測周波数/基準周波数)

計測周波数=基準周波数・2^(セント値/1200)

サンプリング周期

250[msec]

4回/秒 サンプリング設定と使いやすさの相関から選定

２．課題と対策
１）基音
（１）課題
　・基音より強い倍音が存在する場合がある
　・計測したい基音周波数が計測範囲外にある
（２）対策
　・特許公報H3-10118より
　　最も大きなスペクトル値よりも30dB程度の範囲内で小さなスペクトル値を示す周波数成分の周波数を基音として判定する
　・計測範囲外の基音については２次高調波と３次高調波の差または３次高調波と４次高調波の差を基音周波数とする

２）計測精度
（１）課題
　・CMSIS-DPS Real FFTライブラリーの最大ポイント数4096では分割ピッチが目標計測精度より荒い
　　分割ピッチ＝サンプリング周波数/4096
（２）対策
　・サンプリング周波数を基音に対して可変させて分割ピッチを最適化する
　・スプライン補間で誤差を少なくする

３．処理方法

処理

手順

詳細

データサンプリンング

計測：内蔵ADC

サンプリング周波数　16,000[Hz]

分解能 12bit

読込　DMA

デバッグ：

USB Speaker

サンプリング周波数　48,000[Hz] / 3

分解能 16bit / 16

読込　USB Speaker ISO プロジェクト流用

基音判定

CMSIS Real FFT １回目

R FFT計算

各ポイントのパワー計算

ピークパワーポイント選定

　ピーク平均値の２倍以下はピーク判定から除外

　大きい順に２個選定

if ピークパワーポイント最小値≒ピークパワーポイント間隔

　基音＝ピークパワーポイント最小値

　ピークパワーポイントの次数＝１

else

　基音＝ピークパワーポイント間隔

　スプライン補間で使用するピークパワーポイントの次数を計算

基音周波数計測

サンプリング周波数最適化

２回目サンプリング周波数 = 125 * 2 ^ xn

xn は２回目サンプリング周波数が基音周波数の概略16倍となる整数

CMSIS Real FFT ２回目

サンプリング周波数を変更

R FFT計算

各ポイントのパワー計算

CMSIS Spline Interpolation

基音またはR FFTピークパワーポイントを中央としてスプライン補間

補間したピークパワーポイントを抽出

基音周波数計算

基音周波数＝スプライン補間ピーク周波数/次数

４．解説
１）サンプリング周波数
（１）周波数選定
・FFTポイント数の倍数が処理しやすいがデバッグで使う48kHzとサンプリング周期をもとに選定

（２）サンプリング最適化周波数
・理論的には計測周波数の２倍以上と言われているが16倍にしている
・倍率を下げればFFT分割ピッチがより小さくなり計測精度が上がる
・現状実精度は目標値と同等レベルとなっているが音源の音程が揺れるため正確な検定は難しい

２）オクターブ値
・調弦する場合オクターブ値は必要ないので計測を省略している
・FFT１回目でピークポイントをより多く抽出して基音判定精度を上げればより正確にオクターブ値を計測できる

３）FFTパワー計測例
　高調波パワーが大きい場合
　Cello C弦　C1 : 65.406[Hz]

続く

STM32F411 BlackPill CMSIS-DSPによる楽器用チューナー（１）概要

BlackPill STM32

１．概要
　ARM社のCMSIS-DSPライブラリーを使ってSTM32F411 Cortec-M4で楽器用チューナーを作ってみる

　表示イメージ：YAMAHA TD-90
　　音名、周波数等はOLEDで表示

２．チューナー構成
１）構成と課題
　チューナーは周波数カウンターで出来ています
　楽器の音には多くの倍音成分が含まれていてシンプルな周波数カウンターでは聴感と異なる測定結果となる場合がある
　このため聴感上の基音となる周波数を識別して計測しなければならない

２）構成例
　特許公報から３種類抜粋
（１）特許公報　平3-42412
　発明の名称    ピッチ測定装置及び方法
　出願日        1982年12月13日　優先権主張　1981年12月14日 US330681
　概略        波形のゼロクロス極性と周期の相関により基音成分を抽出識別して周波数を計測する

（２）特許公報　6225818
　発明の名称    ピッチ情報生成装置、ピッチ情報生成方法、及びプログラム
　出願日        2014年4月30日
　概略        エンベロープ波形を生成して基音成分を抽出識別して周波数を計測する

（３）特許公報　H3-10118
　発明の名称    音符の表示装置
　出願日        1983年1月26日
　概略        FFTパワースペクトルで基音を抽出識別する

３．チューナー構成選定
　CMSIS-DSPライブラリーを使ってみることが目的でもあるためFFTで計測する

４．解説
・楽器用チューナーはシンプルな周波数計測とは異なり人間の聴感の基づく計測が必要なため難易度が高い
・従来の手法ではアナログ的なパターンマッチング法とデジタル的なFFT分析法がある
・40年以上前の技術で実用上スマホアプリで十分
・今回はCMSIS-DSPライブラリーの導入を主な目的とする
・特許公報は下のリンクから

www.j-platpat.inpit.go.jp

続く