2018年6月11日月曜日

最小構成環境でSingleton

 Singleton パターン - Wikipedia

 組み込みの最小構成な(mallocが無い)環境ではwikipediaの例はうまく動かない。
 関数の中でstaticに宣言するとnewで確保されるようだ。

 解決方法としては、Singletonクラスのprivateでstaticに宣言し、適当な場所に実態を置いておく。
 こうすると、Singletonコンストラクタが呼ばれるのはグローバル変数の初期化時のはずなので、いつ初期化されるかわからない、という事態も避けることができるはず。
 まぁSingletonのコンストラクタにタイミングクリチカルな処理なんて書くな、という話なのだが。

***

 久しぶりにSTM32F1で遊んでる。なぜかUSART bootloader経由でプログラムを書き込むプログラムを作るところから始めてしまった。
 ウチの環境ST-Link(Nucleo割った)が使えないので、USART bootloaderが使えるとF1とかもそのまま使えるので便利。かな。
 STM32F0あたりだとPIC16F88より安い(チップ単体比)ので、ごく小規模なモノなら気兼ねなく使えるかも。

 STBee mini使うの久しぶりすぎてだいぶ手こずった。
 今回はprintfはおろか、USARTといった文字出力系の初期化もしていないので、バイナリサイズは極めて小さいが、デバッグはかなり手間取る。
 バイナリサイズが小さいのは転送速度が遅い(max115.2kbaud+オーバーヘッドがでかい)で相殺されてるので書き込み速度はUSB DFUを使った100k程度のバイナリと体感であまり変わらない。ビルドは早い気がする。

2018年6月6日水曜日

パルス圧縮

 ポカミスで変なバグ仕込んだ挙げ句に一晩熟成されてて全く原因わからず全波形ダンプしたので、ついでにパルス圧縮前と圧縮後の比較。

 全期間の波形

 10msec前後の拡大

 青が入力波形、赤がパルス圧縮後の波形。
 青は12bitADCの信号そのまま、赤はFFT後の値そのまま。横軸の単位はミリ秒。

 条件はパルス幅5msec、計測時間80msec、ディレイ2msec、Fc40kHz、ΔF5kHz、DAC500kSPS、ADC125kSPS、観測時30.3degC / 349.4m/sといった感じ。

 上の図では2msecのオフセットをしていないので、実時間に変換するには2msecを足す必要がある。

 0msecからの大きなピークは直接波を見ている。10msecのピークは天井からの反射波を見てる。
 パルス圧縮すると7msecにピークがある。実時間だと9msecあたりで、実際の天井の位置とおおよそ一致している。

 受信波がかなり乱れてるが、7msecから受信が始まり、パルス幅が5msecなので、12msecで受信が終わっている。心の目で見ると実際にそういうふうに受信されているように見える。
 パルス圧縮後の波形はおよそ0.3msec程度くらいに見える。ということは、実測の圧縮比は5msec/0.3msecで16くらいか。
 計算上は5msec*5kHzでパルス圧縮比25になる。だいぶ特性悪化してる。それでも半分までは落ちてないが。

***

 ゲート幅を狭くすれば意外と短時間で処理が終わる。
 受信素子を2個にすればモノパルス2次元追尾レーダーくらいは作れそうな気もする。うまく作れればサバゲの自動砲台とかも作れそう。
 ただ、探知距離が極めて短く、10m未満くらいが実用的な範囲。エアガンの射程を考えるともう少しレンジがほしい。直径15cmくらいのパラボラが有ればいいかな? パラボラの軸に合わせて送信素子を置いて、少しずらした2個の受信素子を置けば良さそう。

 捜索モードは素子1個を使ってゲート幅を広くし、レンジ内に十分な強度のエコーが有れば捕捉モードに移行して2個の素子で十分に受信できるようにし、その後に追尾モードに移行して2個の素子の強度の差が規定値になるようにする、というような挙動が必要なはず。

 砲台として使う場合は、ある程度の距離に入り、かつ接近率が0以上の場合は射撃を行い、接近率が0を下回った場合は射撃を停止する、という感じか。おそらく接近率が0未満の目標は無視して捜索モードに移行する、という処理になるはず。
 ということは、走って砲台に接近し、捕捉モードに入った時点で後退し捜索モードに移行させ、というサイクルを繰り返せば砲台の直近までは移動できそう。となると、ある程度の距離未満になった場合は接近率にかかわらず射撃する、という処理も必要かも。

 CIWSファランクスだと自前のレーダーで自分の弾丸の弾道を調べてフィードバックさせるけど、さすがに波長8mmの超音波でそれをやるのは無理な気がする。
 数百kHzの素子もあるから、それを使えばBB弾を捜索できるかもしれないけど、Cortex-M4じゃ無理そうだなぁ。

2018年6月2日土曜日

超音波距離計

 超音波距離計は飽きたと言ったな、あれは(ry


 FFTで作った相関器を試してみました。

 DAC800kSPS、ADC50kSPS、パルス幅10.24msec、f40kHz、ΔF5kHz、サンプリング期間8192ポイント、という感じのパラメータです。

 送信は十分に高い周波数で行っていますが、受信は40kHzの信号を50kSPSでサンプリングしています。明らかに元の波形は取れないわけですが、ちゃんと相関が取れています。崩れたなりの波形同士でちゃんと相関処理ができてるってことですかね。

 今まではDACとADCの周波数が同じだったので、DACで送信した波形をADCの相関に使用していました。そのあたりの変更がちょっと面倒でした。

 あとCMSIS DSP LibではFFTが4kポイントまでという制限、それからRAMが足りないという制限により、FFTは1kポイントで実行しています。


 上の図では1.5m先の天井付近に鋭いピークがあります。
 以前は天井で反射し、机で反射し、さらに天井で反射したようなタイミングのエコーがありましたが、今回はそれは見られません。ゴチャゴチャしてるどれかがそれなのかも。
 サンプリング周波数が低いので、そのあたりの悪影響が意外と大きいかもしれません。


 とりあえずパルス幅10.24msecで約140msecのサンプリングが可能でした(運が悪いとmallocでコケますが)。
 140msecだと20m程度の距離に相当しますから、1回のパルスで20m程度の幅のエコーを検出できます。以前は5m程度までだったので、大幅な改善です。とはいえ、これはFFTとかは関係なく、サンプリング周波数が下がったのが大きな要因だと思いますが。


 今回は全期間(140msec程度)をサンプリングした後にまとめて相関処理に入れていますが、なんとなく、計測時間よりも早く相関ができそうな気がします。ということは、ADCバッファの長さはせいぜい2048サンプル分有れば良いわけで、今は8192サンプル分を確保していますから、12KiB程度を開けることができます。
 その分に結果を保存しながらサンプリングしながら、というふうにすれば、かなり長い期間のサンプリングも可能な気がします。
 今の所、すべてのDAC/ADC転送が終わったのはDMAのカウンタをソフトウェアループで見て確認していますが、転送中に相関処理を行うならDMAの割り込み等を駆使する必要があります。そのあたりがかなり複雑になりそうですね。

2018年6月1日金曜日

STM32F4でFFT

 F4でFFTをやってみた。

 テストコード

bool flag(true);
uint8_t is_complex(0);
size_t length(0);
float SPS(0);
float freq(0);

const uint8_t ifft_flag(0);
const uint8_t do_bit_reverse(1);

flag = flag && (4 == sscanf(command,
                            "fft %hhu %u %f %f",
                            &is_complex,
                            &length,
                            &SPS,
                            &freq));

if (flag)
{
    printf("FFT %s %upoint %.1fSPS %.1fHz\n",
            is_complex ? "complex" : "real",
            length,
            SPS,
            freq);
}

if (is_complex)
{
    const arm_cfft_instance_f32 *S(0);

    if (flag)
    { // check length
        switch (length)
        {
        case 16:
            S = &arm_cfft_sR_f32_len16;
            break;
        case 32:
            S = &arm_cfft_sR_f32_len32;
            break;
        case 64:
            S = &arm_cfft_sR_f32_len64;
            break;
        case 128:
            S = &arm_cfft_sR_f32_len128;
            break;
        case 256:
            S = &arm_cfft_sR_f32_len256;
            break;
        case 512:
            S = &arm_cfft_sR_f32_len512;
            break;
        case 1024:
            S = &arm_cfft_sR_f32_len1024;
            break;
        case 2048:
            S = &arm_cfft_sR_f32_len2048;
            break;
        case 4096:
            S = &arm_cfft_sR_f32_len4096;
            break;
        default:
            flag = false;
            printf("length error\n");
            break;
        }
    }

    float *FFT_input(0);

    if (flag)
    { // alloc
        try
        {
            FFT_input = new float32_t[length * 2]();
        }
        catch (std::bad_alloc)
        {
            flag = false;
            printf("bad alloc\n");
        }
    }

    if (flag)
    { //generate wave
        for (size_t i(0); i < length; i++)
        {
            const float phase(i * pi * 2 / SPS * freq);
            FFT_input[i * 2 + 0] = cosf(phase);
            FFT_input[i * 2 + 1] = sinf(phase);
        }
    }

    float FFT_process_seconds(0);

    if (flag)
    { // FFT
        const uint32_t start(HAL_GetTick());
        arm_cfft_f32(S, FFT_input, ifft_flag, do_bit_reverse);
        const uint32_t stop(HAL_GetTick());

        FFT_process_seconds = (stop - start) * 0.001f;
    }

    if (flag)
    { // dump result
        printf(" %.3fsec\n\n re im\n", FFT_process_seconds);

        for (size_t i(0); i < length; i++)
        {
            const int j(i - length / 2);
            const int k(i < length / 2 ? (i + length / 2) : (i - length / 2));
            printf("%f %f %f\n",
                    j * SPS / length,
                    FFT_input[k * 2 + 0],
                    FFT_input[k * 2 + 1]);
        }
    }

    delete[] FFT_input;
    FFT_input = 0;
}
else
{
    arm_rfft_fast_instance_f32 S;

    if (flag)
    { // check length and initialize
        flag = (ARM_MATH_SUCCESS == arm_rfft_fast_init_f32(&S, length));

        if (!flag)
        {
            printf("init error\n");
        }
    }

    float *FFT_input(0);
    float *FFT_output(0);

    if (flag)
    { // alloc
        try
        {
            FFT_input = new float32_t[length]();
            FFT_output = new float32_t[length]();
        }
        catch (std::bad_alloc)
        {
            flag = false;
            printf("bad alloc\n");
        }
    }

    if (flag)
    { //generate wave
        for (size_t i(0); i < length; i++)
        {
            const float phase(i * pi * 2 / SPS * freq);
            FFT_input[i] = cosf(phase);
        }
    }

    float FFT_process_seconds(0);

    if (flag)
    { // FFT
        arm_rfft_fast_f32(&S, FFT_input, FFT_output, ifft_flag);
    }

    if (flag)
    { // dump result
        printf("%.3fsec\n\n re\n", FFT_process_seconds);

        for (size_t i(0); i < length; i++)
        {
            printf("%f %f\n",
                    i * SPS / length / 2,
                    FFT_output[i]);
        }
    }

    delete[] FFT_input;
    FFT_input = 0;

    delete[] FFT_output;
    FFT_output = 0;
}

 リンカでDrivers/CMSIS/Lib/GCC/libarm_cortexM4lf_math.aを読み込む必要がある。

 ARM_MATH_CM4(コアに合わせて選択)をdefine定義した後、arm_math.hとarm_const_structs.hをインクルードする。あとnewのstd::bad_allocをcatchするのに<new>もインクルードする必要がある。


 CFFTは複素数、RFFTは実数の処理で、それぞれ呼び方が違う。また引数のとり方も違い、CFFTはinputとoutputを同じ配列で行う。一方、RFFTはinputとoutputは別の配列を使う。
 CFFTは負の周波数も解析できる。一方、RFFTは負の周波数は扱えないが、CFFTの分解能の2倍くらいあると思う。1024ポイントでSPS1kHzとした場合、CFFTは1024ポイントで-0.5kHzから+0.5kHzの範囲を取る。RFFTは1024ポイントで0Hzから0.5kHzの範囲を取る。

 ちゃんと試してないけど、CFFTは16ポイントから4kポイントまで、RFFTは32ポイントから4kポイントまで解析できる。

 CFFT/RFFTを行う上記コードを含めるとバイナリサイズは256kくらいになった。この部分をコメントアウトすると128k位になる。FFT解析を行うには128kほどのFlashが必要になる。もちろん、FFTを行う上でメモリ領域も必要になる。


 const char *const commandに文字列を入れておけば、それを解析する。主にUART等からコマンドを受け取って実行する事を想定。もちろんsprintfで作ってもいいけど、そんな事するくらいなら直接呼べという話で。

 "fft 1 4096 100 -10" というコマンドを実行すると、複素数で100SPS/-10Hzの波形を生成しFFT解析する。
 "fft 0 4096 100 30" というコマンドを実行すると、実部のみで100SPS/30Hzの波形を生成しFFT解析する。
 それぞれの結果は以下の通り。



 CFFTは左端が-50Hz、右端が+50Hz、RFFTは左端が0Hz、右端が+50Hzになっている。どちらも結果は4096ポイント得られるので、RFFTのほうが倍の分解能になる。

 4096ポイントの場合、168MHzのCortex-M4FでCFFTが4msec未満、RFFTが1msec未満だった。CFFTは若干の時間がかかるが、RFFTは恐ろしく早い。
ただ、どちらも4kポイントだと32KiBのRAMが必要になる。自由に使えるヒープ領域は100KiB程度なので、あまり大規模なものを作ろうとするとちょっと心もとない。

***

 とりあえずFFTが動くようになったので、パルス圧縮の相関処理とかもできるわけだが、時間的・空間的な余裕がかなり少なそう。
 超音波距離計もちょっと飽きてきたし、気分転換に別のヤツで遊ぼうかな。
 FFT/IFFTが動くならOFDM変調とかできる。

2018年5月31日木曜日

CubeMXで新規プロジェクト まとめ

 前に書いたやつのまとめ。

 最低限、コンパイラが吐くエラーを自力で修正できる必要がある。あとMakefileをかなりいじるので、そのあたりの知識も。
 それとsyscalls.cというファイルが必要になる。どっかから拾ってきてください。

 Win10のWSLで実行することを前提にしています。まぁmakeとか一通り動いてarm-none-eabi-gccとかも一通り入ってればだいたい問題ない。

1 プロジェクトを作る

File > New Project (Ctrl-N)でNew Projectを表示する。
 使用するMCUを選択する。左側のMCU Filtersでコアやシリーズ、パッケージを選択したり、製品名を入力する。
 今回はSTBee F4miniを使用するので、STM32F405RGを選択する。
 MCUs List n itemsでアイテムをダブルクリックして終了。


 Project > Settings...でProject Settings画面を開く。
 Project Nameでプロジェクト名を入力する。
 Toolchain / IDEをMakefileにする。

 Code GeneratorタブでGenerate peripheral initialization as a pair of '.c/.h' files per peripheral、Set all free pins as analog (to optimize the power consumption)、Enable Full Assertの3個にチェックを入れる。

 OKでProject Settingsダイアログを閉じる。

2 初期設定

GPIOの設定を行う。

 左のツリーからRCCのHigh Speed Clock (HSE)とLow Speed Clocl (LSE)をCrystal/Ceramic Resonatorにする。

 左のツリーからUSB_OGT_FSのModeをDevice_Onlyにする。
 USB_DEVICEのClass For FS IPをCommunication Device Class (Virtual Port Com)にする。

 右のパッケージのイラストからPA0-WKUPをクリック、GPIO_Inputを選択。右クリックでEnter User Labelをクリックし、USER_SWと入力。同じようにPD2をGPIO_Outputにし、USER_LEDと入力。

3 クロックの設定

Clock Configurationタブを開く。
 Clock configuration Do you want to run automatic clock issues solver? というダイアログが出てくるので、Noで閉じる。
 左側のInput frequencyをボードの水晶の値である12にする。
 PLL Source MuxでHSEを選択、System Clock MuxでPLLCLKを選択する。
 右側のHCLK to AHB bus, core, memory and DMA (MHz)に168と入力し、Enterで反映する。少し待たされた後、赤枠のエラーが消えたりして設定が反映される。

4 その他の設定

Configurationタブを開く。
 左側のFREERTOSのEnableをチェックする。


 右側のFREERTOSアイコンをクリックし、FREERTOS Configurationダイアログを開く。
 Task and Queuesのタブを開く(デフォルトでこのタブのはず)。

 TasksのAddでNew Taskダイアログを開く。
 とりあえずLEDを点滅させるタスクをつくる。Task NameにLED_blink、PriorityにosPriorityLow、Stack Sizeに128、Enter FunctionにLED_blink_func、Code Generation OptionにAs externalを設定し、OKでダイアログを閉じる。
 コンソールのタスクも作っておく。Task Nameにconsole、PriorityにosPriorityBelowNormal、Stack Sizeに1024、Entry Functionにconsole_func、Code Generation OptionにAs externalを設定する。

 QueuesのAddでNew Queueダイアログを開く。
 Queue Nameにstdout_queue、Queue Sizeに256、Item Sizeにuint8_tとし、OKで閉じる。
 再びAddでダイアログを開き、stdin_queue、128、uint8_tとして追加。

 Config parametersタブでUSE_IDLE_HOOKとUSE_TRACE_FACILITYとUSE_STATS_FORMATTING_FUNCTIONSをEnableにする。

 OKでFREERTOS Configurationダイアログを閉じる。

5 コード生成

Project > Generate Code (Ctrl-Shift-G)をクリック。
 Warning: Code Generationと警告が出てくるが、Yesで続行する。

 しばらく待つとThe Code is successfulyとダイアログが出るので、Open Folderで出力フォルダを開くか、Closeでダイアログを閉じる。

6 コードの修正と追加

Makefileを開く。

 C_SOURCES =  \の下に Src/syscalls.c \ と Src/CCM.c \ を追加。

 # ASM sourcesの上に以下のコードを追加する。
"
CPP_SOURCES = \
Src/LED_blink.cpp \
Src/Console.cpp \
"

 CC = とかAS = とかの$(BINPATH)/$(PREFIX)を$(BINPATH)$(PREFIX)にする(スラッシュを消す)。
 CCの下にCPP = $(BINPATH)$(PREFIX)g++を追加する。

 LDFLAGS = にある-specs=nano.specsを削除する。

 CFLAGS += -MMD -MP -MF"$(@:%.o=%.d)" -MT"$(@:%.o=%.d)" を CFLAGS += -MMD -MP -MF"$(@:%.o=%.d)" に変更。

 # list of objectsの OBJECTS = $(addprefix $(BUILD_DIR)/,$(notdir $(C_SOURCES:.c=.o))) を OBJECTS = $(addprefix $(BUILD_DIR)/,$(C_SOURCES:.c=.o)) に変更。

 # list of ASM program objectsの上に以下のコードを追加。
"
OBJECTS += $(addprefix $(BUILD_DIR)/,$(CPP_SOURCES:.cpp=.o))
vpath %.cpp $(sort $(dir $(CPP_SOURCES)))
"

 vpath %.s $(sort $(dir $(ASM_SOURCES))) の下に以下のコードを追加。
"
DEPS = $(OBJECTS:%.o=%.d)
-include $(DEPS)
"

 OBJECTS += $(addprefix $(BUILD_DIR)/,$(notdir $(ASM_SOURCES:.s=.o))) を OBJECTS += $(addprefix $(BUILD_DIR)/,$(ASM_SOURCES:.s=.o)) に変更。

 $(BUILD_DIR)/%.o: %.c Makefile | $(BUILD_DIR) の下に mkdir -p $(dir $@); を追加。

 $(CC) -c $(CFLAGS) -Wa,-a,-ad,-alms=$(BUILD_DIR)/$(notdir $(<:.c=.lst)) $< -o $@ を $(CC) -c $(CFLAGS) -Wa,-a,-ad,-alms=$(BUILD_DIR)/$(<:.c=.lst) $< -o $@ に変更。

 $(BUILD_DIR)/%.o: %.s Makefile | $(BUILD_DIR)の上に以下のコードを追加。
"
$(BUILD_DIR)/%.o: %.cpp Makefile | $(BUILD_DIR)
mkdir -p $(dir $@);
$(CPP) -std=c++11 -c $(CFLAGS) -Wa,-a,-ad,-alms=$(BUILD_DIR)/$(<:.cpp=.lst) $< -o $@
"

 $(CC) $(OBJECTS) $(LDFLAGS) -o $@を$(CPP) $(OBJECTS) $(LDFLAGS) -o $@に変更する。


 STM32F405RGTx_FLASH.ldを開く。
 _Min_Stack_Size=の下に_heap_end = ORIGIN(RAM)+LENGTH(RAM)-4;を追加する。
 *(.ccmram*)の下に*CCM.oを追加する


 Inc/FreeRTOSConfig.hのUSER CODE END Definesの上に #define configAPPLICATION_ALLOCATED_HEAP 1 を追加する。


 SrcフォルダにCCM.cを追加。以下の行を追加。
"
#include <FreeRTOS.h>

#if( configAPPLICATION_ALLOCATED_HEAP == 1 )
uint8_t ucHeap[configTOTAL_HEAP_SIZE];
#endif
"


 SrcフォルダにLED_blink.cppを追加。以下の行を追加。
"
#include <cmsis_os.h>
#include <main.h>
#include <stm32f4xx_hal.h>

extern "C" void LED_blink_func(void const *argument)
{
    while (1)
    {
        HAL_GPIO_TogglePin(USER_LED_GPIO_Port, USER_LED_Pin);
        osDelay(HAL_GPIO_ReadPin(USER_SW_GPIO_Port, USER_SW_Pin) == GPIO_PIN_SET
                    ? 100
                    : 500);
    }
}
"


 SrcフォルダにConsole.cppを追加。以下のコードを追加。
"
#include <cmsis_os.h>
#include <stm32f4xx_hal.h>
#include <string.h>

extern "C" void console_func(void const *argument)
{
    while (1)
    {
        char str[100];
        scanf("%99[^\n]%*[^\n]", str);
        getchar();

        {
            char *const p(strchr(str, '\r'));

            if (p)
            {
                *p = '\0';
            }
        }

        if (!strcmp(str, "tasklist"))
        {
            char buff[10 * 45];
            osThreadList((uint8_t *)buff);
            printf(
                "Name            State   Priorty Stack   Num\n"
                "*******************************************\n"
                "%s"
                "*******************************************\n",
                buff);
        }
    }
}
"



 Srcにどこかから持ってきたsyscalls.cをコピーして入れておく。

 _read_rの中を以下のようにする。
"
uint8_t *p = (uint8_t *)ptr;
extern osMessageQId stdin_queueHandle;

while (len--)
{
    uint8_t data;
    xQueueReceive(stdin_queueHandle, &data, portMAX_DELAY);
    *p++ = data;

    if (data == '\n')
    {
        break;
    }
}

_ssize_t i = (uint32_t)p - (uint32_t)ptr;

return (i);
"

 _write_rの中を以下のようにする。
"
uint8_t *p = (uint8_t *)ptr;
int i = 0;
extern osMessageQId stdout_queueHandle;

if (__get_IPSR() == 0)
{
    for (i = 0; i < len; i++)
    {
        xQueueSend(stdout_queueHandle, p++, portMAX_DELAY);
    }
}
else
{
    for (i = 0; i < len; i++)
    {
        xQueueSendFromISR(stdout_queueHandle, p++, 0);
    }
}

return (i);
"


 Src/freertos.cを開く。
 #include <usbd_cdc_if.h> を追加する。

 vApplicationIdleHookに HAL_PWR_EnterSLEEPMode(PWR_MAINREGULATOR_ON, PWR_SLEEPENTRY_WFI); を追加する。

 StartDefaultTaskの無限ループの前に osDelay(500); を追加する。
 無限ループの中を以下のように書き換える。
"
static uint8_t buffer[32];
uint16_t count = 0;
BaseType_t dequeue_state = pdPASS;

while (dequeue_state == pdPASS && count < sizeof(buffer))
{
    dequeue_state = xQueueReceive(stdout_queueHandle, &buffer[count], 1);

    if (dequeue_state)
    {
        count++;
    }
}

if (count > 0)
{
    CDC_Transmit_FS(buffer, count);
}
osDelay(2);
"



 Src/usbd_cdc_if.cを開く。
 #include <cmsis_os.h> を追加。

 CDC_Control_FS関数のスイッチのケースCDC_GET_LINE_CODINGに以下のコードを追加。
"
*((uint32_t*)&pbuf[0]) = 115200;
pbuf[4] = 1;
pbuf[5] = 0;
pbuf[6] = 8;
"

 CDC_Receive_FS関数のreturnの前に以下のコードを追加。
"
extern osMessageQId stdin_queueHandle;
uint32_t i = 0;
for (i = 0; i < *Len; i++)
{
    xQueueSendFromISR(stdin_queueHandle, &Buf[i], 0);
}
"


 あとはmakeすればbuildフォルダにelfとかhexとかbinとかが作られてるので、任意のツールでマイコンに書き込めばOK。

 注意点として、このmakefileはヘッダの変更を認識しない。ヘッダファイルを書き換えた場合は、責任を持ってそれに依存するすべてのソースファイルをタッチするか、あるいはmake cleanでオブジェクトファイルを削除して再びビルドすること。

 Cubeが吐いたmakefileにはC_SOURCESに絶対パスを指定するファイルが含まれる場合がある。その場合は最初のスラッシュを削除し、相対パスで認識されるようにすること。

 未確認だが、割り込みの中でHALのタイムアウトを指定する関数を呼ぶ場合、うまく動かないと思う(タイムアウトしないはず)。対策はいろいろあるが、そもそも割り込みの中で何ミリ秒もポーリングするような処理はするな、ということで。
 ほぼ同じところが原因で、FreeRTOSのコンテキストスイッチを1kHz以外に設定するといろいろ不具合を起こすはず。


* WSLでGCCが古い問題
 Win10 WSLのUbuntuのaptで入れたarm-none-eabi-gccはバージョン4.9.3で更新が止まってる。新しいGCCを使いたい場合は、GNU Arm Embedded Toolchain | Downloads – Arm DeveloperからWindows ZIPをダウンロードしてきてスペースや全角文字を含まないディレクトリに解凍しておく。binフォルダにarm-none-eabi-gcc.exe等一式が入っているので、MakefileのBINPATHにこのフォルダのパスを入れておく。MakefileはWSLで実行するので、ディレクトリ名はC:\等ではなく/mnt/c/のようになる点に注意。そしてCC = $(BINPATH)$(PREFIX)gcc.exeのように拡張子を追加する。
 あとはmakeすれば新しいgccでビルドできる。

 C:\Devz\ARM\launchpad\bin\arm-none-eabi-gcc.exe のようなディレクトリ構成の場合、 BINPATHは /mnt/c/Devz/ARM/launchpad/bin/ のようになるはず。


更新:2018/06/15
 [update]依存関係を正常に処理できるように
 [add]BINPATHの例を追加

マッチドフィルタ/パルス圧縮/FFT相関




 上がマッチドフィルタ風、下がパルス圧縮風。
 受信信号と想定した正弦波は-1から+1の範囲を取る。ノイズと想定した疑似乱数は-1から+1の範囲を取る。

 マッチドフィルタ風は「ΔFが0なパルス圧縮」、つまり圧縮比が0のパルス圧縮と見ることができる。実際、上記のグラフはΔFの違いしか無い。

 マッチドフィルタでは近い範囲にエコーがあると判別ができない。パルス圧縮では十分に見分けることができる。

 レーダーの探知性能は送信出力とパルス幅の、いわゆる力積と呼ばれるもので、ざっくり言えばパルスの振幅とパルス幅の積がレーダー波のエネルギーになる。エネルギーが大きいほど遠距離の、あるいは小型な目標のエコーを拾うことができる(受信感度等は割愛)。
 この図ではマッチドフィルタもパルス圧縮も同じ振幅、同じパルス幅なので、どちらも探知性能としては同様。ただしパルス圧縮のほうが距離分解能が高い。また、パルス圧縮のほうが立ち上がり・立ち下がりが急峻であり、距離精度が改善する(目標のエコーが遠近方向に広がらずに表示される)。

***


 この図は、上がチャープ信号をFFT解析したもの、下はそれをIFFTに通したもの。

 FFT解析には下記ブログのソースコードを流用した。
 C#で画像処理 高速フーリエ変換(FFT)
 C#で画像処理 逆高速フーリエ変換(IFFT)

 FFT結果は凄まじいことになっているが、IFFTを通すとちゃんと元通りになる。

***


この図は、パルス圧縮の波形を単純な積で計算したものと、FFTを使って計算したものの比較。赤がfor2重ループで計算したもの、緑がFFTで計算したものだが、どちらも同じ結果になっている(最後の方が一致していないのはサンプル数の問題)。
 正弦波が±1なのに対し、ノイズが±5程度あるが、十分な強度としてパルスが見えており、狭い範囲にある2つのパルスを判別することができる。


 以下のページを参考にした
 離散フーリエ変換(3) - 畳み込み定理 - uhiaha888’s diary

 for2重ループは時間領域、FFTは周波数領域、というタイプなのかな?
 とある本によると、入力信号と基準波形をDFT計算し、その結果を乗算後に逆DFTに通す、というふうに書いてある。処理としては積をDFTに通せばいいらしい。
 それと、入力信号が1024サンプル、基準波形が128サンプル、といった場合、基準波形を1024サンプルになるようにゼロ埋めし、それらをFFT解析してから積をとり再びFFTに入れる、という処理になるらしい。
 また、次の入力信号の解析を行うときは、基準波形分の128サンプルの領域が重なるように処理する必要がありそうな気がする。とはいえ、これに関してはfor2重でも同じような問題が出てくるので、欠点と言えるほどではないかも。

 基準信号・入力信号ともにスペクトルの大部分が低周波領域に集中しているから、乗算を行うのは全体の2%程度で済んでいる(サンプリング周波数と目的の周波数の関係で変わるはず)。今回はサンプリング周波数とチャープ信号の中心周波数の比が200なので2%で済んでいるが、実際はこの比は20とかそれ未満の場合もあるだろうから、本来であればFFTのかなりの領域の積を取る必要があるはず。

 基本的にレーダーであれば基準信号がコロコロ変わることは少なく、近距離と遠距離で送信波を変える場合でも数種類で済むはず。ということは、基準波形のFFTは最初の1回行えば済むわけだから、データ解析のときは入力信号のFFTを行い、その積を再びFFTに入れる、という処理だけでいいはず。とはいえ、FFT2回分の計算量がどの程度になるのか。
 実際のレーダーでは50MサンプルくらいのFFT解析になるのかもしれない。とてつもない分量。

 非常にざっくりとした数値だけど、今回の16384ポイントの相関処理では、forだと220msec、FFTだと8msec程度で終わっている。かなり早い。

 入力信号、基準信号、出力波形は実部のみだが、入力信号と基準波形の結果、それとその積は虚部も必要になる。PCで処理するなら問題にならないデータ量だが、組み込み機器、特にワンチップマイコンで外付けRAMなしの場合は非常に厳しいと思う。
 必要なメモリは6x入力信号くらいか。入力信号が16384サンプルなら、必要なメモリは384KiB程度(float32の場合)。他に一時的なテーブルが必要になる。FFTのIN/OUTを16bitで計算するならその半分としても200KiB程度必要になる。
 STBeeF4miniだと100KiB程度しか使えないので、8192ポイント程度しかできないか? 8kポイントだと距離にして1.5m程度のゲート幅しか得られない。ただし8kポイントを取得するのに10msecかかるから、FFT2回が6msec程度で終わるなら結果を記録するメモリが足りる分だけ延々とサンプリングできる。


 まずはARMコアでFFTがどれくらい動くかを確かめる必要があるか。DSP_LibにFFTっぽいライブラリがあるけど、使ったことがないのでちゃんと動くかどうかをテストしなければ。それがちゃんと動けば、適当な値でFFTを通して処理時間を測って。十分に短い時間で処理できるならいろいろ遊べそうだ。
 FFTとかちゃんと動くなら、適当にパッケージの大きなF4を買ってきてSRAMを外においたほうが使いやすいかもしれない。
 別件でFFTがあると便利そうだと思ってる部分があるので、それも含めてFFTの評価は軽くやっておきたい。

2018年5月27日日曜日

独自float16

 マイコンで生成したint64_tの配列をPCに転送したいので、独自の浮動小数点を考えてみた。

 似たような規格として、既存の半精度浮動小数点数があるが、マイコンでのサポートがない(多分)、表現できる範囲が狭い、といった問題がある。
 前者に関しては、それっぽい値に変換するライブラリをかけば済む話だが、後者が問題で、IEEE 754の半精度浮動小数点数は10進でたかだか4.8桁程度しか表現できない。int64_tは10進で18.9桁ほどあるので、圧倒的に足りない。
 また、一般的な浮動小数点は1未満かつ-1より大きい値を表現できるが、あくまで整数の中間表現がほしいだけなので、1未満の分解能は必要ない。ということで、符号1bit、指数6bit、仮数9bit、の16bit浮動小数点を考えた次第。

 これにより、8バイトの64bit整数を2バイトで表現できるようになったので、4分の1のメモリサイズで済む。

 値が小さい時(<1024)は少し特殊な処理をしている。
 たぶん最小値(-2^63)を与えるとバグる。
 今回はそんな端っこまで使う予定はなかったので省いた。というかint32でも足りるんじゃないか、程度の範囲しか使っていないんだよなぁ。浮動仮数幅浮動小数点、なんてのもありかもしれない。

using original_float_16 = uint16_t;

original_float_16 convert_to_BIN(int64_t value)
{
    const bool is_negative(value < 0);

    if (is_negative)
    {
        value = -value;
    }

    uint16_t bit_data(0);

    if (value < 1024)
    {
        bit_data = value;
    }
    else
    {
        uint_fast8_t i(53);
        while (i > 0) // exit with break
        {
            if (value & ((uint64_t)1 << (i + 9)))
            {
                bit_data = ((i + 1) << 9) | ((value >> i) & 0x01FF);
                break;
            }
            i--;
        }
    }

    if (is_negative)
    {
        bit_data |= 0x8000;
    }

    return (bit_data);
}

int64_t convert_from_BIN(original_float_16 BIN)
{
    const uint_fast8_t sign((BIN >> 15) & 0x01);
    const uint_fast8_t exponent((BIN >> 9) & 0x3F);
    const uint_fast16_t fraction(BIN & 0x01FF);

    const int64_t result(exponent
                             ? (uint64_t)(fraction | 0x0200) << (exponent - 1)
                             : fraction);

    return (!sign ? result : -result);
}

 もうちょっとうまくかけそうな気もするけど、とりあえず動いてるので。

 fp→intは各ビットを取り出してビットシフトするだけなので比較的簡単。
 int→fpは最上位ビットの位置を探す必要があるのでループしてる。最上位ビットの位置を1クロックで返す命令とか無いものか。