Two-dimensional Cell Averaging Constant False Alarm Rate (CA-CFAR)の実装

0. 概要

CFARとは得られた信号からピーク値を抽出する閾値アルゴリズムのことである。例えば、以下のような図がCFARを説明するのによい。

https://jp.mathworks.com/help/phased/examples/constant-false-alarm-rate-cfar-detection.html

基本的にピーク値を抽出するのであれば、ある一定の閾値を設ければよいが、上記の図のような右肩上がりの（局在的変化する）信号情報の場合は固定値の閾値を設けたところで上手くピーク値を抽出することは出来ない。

そこで、閾値をAdaptiveに調整する必要があるということで生まれたのがCFARである。簡単にいうと、ある領域における背景ノイズ(偶然誤差や統計誤差)を計算して、それを閾値としている。最も単純なものがCA-CFARであり、他にもワイブルCFARなど様々な手法が提案されている。

ここまでは1次元の信号処理の話であるが、これを2次元に拡張したものもあり、それが2d-CA-CFARである。主に、信号を画像化したようなノイズのある画像から対象物を抽出するのに向いている。今回はこの2d-CA-CFARを実装していく。個人的にはめちゃくちゃよく使う方法なのに、なんでどこにも実装がないのか不思議である。

1. 理屈

理論は凄く簡単である。下のブロック図を参考にすると、すぐに理解ができる。

https://awrcorp.com/download/faq/english/docs/vss_system_blocks/cfar.htm
まず、あるWindowサイズを決定し、そのWindowsサイズに含まれる信号値の平均を計算する。
ここで重要なのは、ピーク値周辺のデータというのは、その周辺の信号値もピーク値に引っ張られて、全体的に高い値を示す傾向がある。このため、Guard Cellというものを使う。これは、ピーク値近傍の信号値は平均値に混ぜないといったものである。

なお、今回は2dなのである画素 (Pixel under test)に着目した時のGuard cell及びBackgroundは以下のようになる。

https://www.researchgate.net/figure/Sliding-windows-for-CFAR-algorithm_fig1_26548024

数式的に表現すると以下になる。0は不採用という意味。

2. 実装

実装は以下である。上手く動作しない場合はハイパーパラメータを調整するとよい。

# import
import sys
import cv2
import numpy as np
import random

# fields
GUARD_CELLS = 5
BG_CELLS = 10
ALPHA = 2
CFAR_UNITS = 1 + (GUARD_CELLS * 2) + (BG_CELLS * 2)
HALF_CFAR_UNITS = int(CFAR_UNITS/2) + 1
OUTPUT_IMG_DIR = "./"

# preparing
inputImg = cv2.imread("YOUR_IMAGE", 0)
estimateImg = np.zeros((inputImg.shape[0], inputImg.shape[1], 1), np.uint8)

# search
for i in range(inputImg.shape[0] - CFAR_UNITS):
  center_cell_x = i + BG_CELLS + GUARD_CELLS
  for j in range(inputImg.shape[1] - CFAR_UNITS):
    center_cell_y = j  + BG_CELLS + GUARD_CELLS
    average = 0
    for k in range(CFAR_UNITS):
      for l in range(CFAR_UNITS):
        if (k >= BG_CELLS) and (k < (CFAR_UNITS - BG_CELLS)) and (l >= BG_CELLS) and (l < (CFAR_UNITS - BG_CELLS)):
          continue
        average += inputImg[i + k, j + l]
    average /= (CFAR_UNITS * CFAR_UNITS) - ( ((GUARD_CELLS * 2) + 1) * ((GUARD_CELLS * 2) + 1) )

    if inputImg[center_cell_x, center_cell_y] > (average * ALPHA):
      estimateImg[center_cell_x, center_cell_y] = 255

# output
tmpName = OUTPUT_IMG_DIR + "est.png"
cv2.imwrite(tmpName, estimateImg)

Enjoy !