はじめに
Kmeansの拡張でkを指定しないXmeansというのがありますが
本記事では画像の減色する際に類似度を指定して、色数を決めてみたいと思います。
普通にXmeansでやろうもんなら遅くて使いもんにならないだろうし、クラスタ数ではなく類似度で指定したいなと思ったので、考えた次第です。
アルゴリズム概要
- 画像を縮小します。
- kMin から kMax まで1ずつ増やして Kmeans 法を行います。
- Kmeans の前後の画像で類似度を計算し、それが類似度が指定されたものより大きければそれを採用します。
類似度について
$d({\bf x}, {\bf y}) / \sqrt{256^3\cdot 3}<\epsilon$を満たすとき、画素の色が一致するとします。
$d:{\mathbb R}^2\times{\mathbb R}^2\rightarrow{\mathbb R}$ はユークリッド距離関数とします。
画像の類似度は、一致する画素の数を全ての画素の数で割ったものとします。
実験
実際に動かしたものを載せます。
変数は
const [kMin, kMax, colorEpsilon, size, bitDepth, attempts] = [2, 16, 0.08, 32, 4, 2];
のように固定します。similarity(類似度)を変更したものが以下です。
処理時間はほぼ変わらないので記載していません。全体で100msで複数回のKmeansで20msぐらいです。
| 減色後の画像 | 類似度 | 色数 |
|---|---|---|
 |
0.8 | 4 |
 |
0.85 | 5 |
 |
0.9 | 5 |
 |
0.95 | 7 |
 |
0.99 | 12 |
類似度が大きくすれば、おおむね色数が増えていることが確認できました。
ソース
/**
* Determine the number of colors based on similarity
* @param {Mat} src source mat
* @param {Mat} dst destination mat
* @param {number} kMin cluster min count
* @param {number} kMax cluster max count
* @param {number} similarity similarity (ex. 0.95)
* @param {number} colorEpsilon color tolerance (ex. 0.08)
* @param {Size} dsize size in resizing
* @param {number} interpolation interpolation
* @param {number} bitDepth bit depth [1, 8]
* @param {number} attempts number of attempts
* @param {number} maxCount the maximum number of iterations/elements
* @param {number} epsilon the desired accuracy
* @returns {Mat} destination mat
*/
static fastXmeans(src, dst, kMin, kMax, similarity, colorEpsilon,
dsize, interpolation, bitDepth, attempts, maxCount = 1e4, epsilon = 1e-4) {
// Resize
const resized = new cv.Mat();
cv.resize(src, resized, dsize, 0, 0, interpolation);
// Initialize
dst = dst || new cv.Mat(src.size(), src.type());
const criteria = new cv.TermCriteria(cv.TermCriteria_EPS + cv.TermCriteria_MAX_ITER, maxCount, epsilon);
const flags = cv.KMEANS_PP_CENTERS;
// Create samples
const samples = new cv.Mat(resized.rows * resized.cols, 3, cv.CV_32F);
for(let y = 0; y < resized.rows; y += 1) {
for(let x = 0; x < resized.cols; x += 1) {
const idx = y * resized.cols + x;
for(let z = 0; z < 3; z += 1) {
samples.floatPtr(idx)[z] = resized.data[idx * 4 + z];
}
}
}
// Decide k(cluster count)
let startDate = new Date();
let x = -1;
const labelsArray = [], ui8CentersArray = [];
for(let k = kMin; k <= kMax; k += 1) {
// K-Means method
const [labels, centers] = [new cv.Mat(), new cv.Mat()];
cv.kmeans(samples, k, labels, criteria, attempts, flags, centers);
labelsArray.push(labels);
// Convert Float to UInt8
const floatToUint8 = x => Math.min(Math.max(Math.round(x), 0), 255);
const ui8Centers = new Uint8Array(k * 3);
for(let c = 0; c < k; c += 1) {
for(let z = 0; z < 3; z += 1) {
ui8Centers[c * 3 + z] = floatToUint8(centers.floatAt(c, z));
}
}
centers.delete();
ui8CentersArray.push(ui8Centers);
// Create kdst from labels and centers
const kdst = new cv.Mat(resized.size(), resized.type());
for(let y = 0; y < kdst.rows; y += 1) {
for(let x = 0; x < kdst.cols; x += 1) {
const p = (y * kdst.cols + x) * 4;
const c = labels.data[p];
for(let z = 0; z < 3; z += 1) {
kdst.data[p + z] = ui8Centers[c * 3 + z];
}
kdst.data[p + 3] = resized.data[p + 3];
}
}
// Compute sim
const sim = computeSimilarity(resized, kdst, colorEpsilon);
if(sim > similarity) {
x = k;
break;
}
}
if(x < 0) { x = kMax; } // not found.
const k = x;
const labels = labelsArray[x - kMin];
const ui8Centers = ui8CentersArray[x - kMin];
console.log('kmeans:', new Date() - startDate, 'ms');
// Create color map
const mem = 2 ** bitDepth;
const quant = (mem - 1) / (256 - 1);
const digit = Math.ceil(Math.log10(mem - 1));
const colorMap = {};
for(let y = 0; y < resized.rows; y += 1) {
for(let x = 0; x < resized.cols; x += 1) {
const p = (y * resized.cols + x) * 4;
const key = createColorKey(resized.data, p);
if(typeof colorMap[key] === 'undefined') {
colorMap[key] = labels.data[p];
}
}
}
// Release
resized.delete();
labelsArray.forEach(e => { e.delete(); });
// Create dst from ui8Centers
startDate = new Date();
for(let y = 0; y < dst.rows; y += 1) {
for(let x = 0; x < dst.cols; x += 1) {
const p = (y * dst.cols + x) * 4;
const key = createColorKey(src.data, p);
// Find nearest color
if(typeof colorMap[key] === 'undefined') {
let [minDist2, minC] = [Number.MAX_VALUE, -1];
for(let c = 0; c < k; c += 1) {
// Compute distance squared in color space
let dist2 = 0;
for(let z = 0; z < 3; z += 1) {
dist2 += (src.data[p + z] - ui8Centers[c * 3 + z]) ** 2;
}
if(dist2 < minDist2) {
minDist2 = dist2;
minC = c;
}
}
colorMap[key] = minC;
}
// Set color
const c = colorMap[key];
for(let z = 0; z < 3; z += 1) {
dst.data[p + z] = ui8Centers[c * 3 + z];
}
dst.data[p + 3] = src.data[p + 3];
}
}
console.log('after kmeans:', new Date() - startDate, 'ms');
return [ dst, k ];
/**
* Create color key
* @param {Uint8Array} data data
* @param {number} p pixel index
* @returns {string} color key
*/
function createColorKey(data, p) {
let key = '';
for(let z = 0; z < 3; z += 1) {
const c = Math.round(data[p + z] * quant); // Quantization
key += c.toString().padStart(digit, '0');
}
return key;
}
/**
* Compute image similarity
* @param {Mat} src source mat
* @param {Mat} dst destination mat
* @param {number} epsilon similarity epsilon ex. 0.08
* @returns {number} image similarity
*/
function computeSimilarity(src, dst, epsilon) {
const sqrt3 = 3 ** 0.5;
let diffCount = 0;
for(let p = 0; p < src.data.length; p += 4) {
const srcColor = getColor(src.data, p);
const dstColor = getColor(dst.data, p);
const diff = diffColor(srcColor, dstColor);
if(diff / (255 * sqrt3) > epsilon) {
diffCount++;
}
}
return 1 - diffCount / (src.cols * src.rows);
function diffColor(src, dst) {
return ((src[0] - dst[0]) ** 2 + (src[1] - dst[1]) ** 2 + (src[2] - dst[2]) ** 2) ** 0.5;
}
function getColor(data, p) {
return [ data[p], data[p + 1], data[p + 2] ];
}
}
}
関数の呼び出し方
// Create src and dst
const src = cv.imread('src-canvas');
const fx = fastXmeans;
const [kMin, kMax, similarity, colorEpsilon, size, bitDepth, attempts] = [2, 16, 0.9, 0.08, 32, 4, 2];
let [dst, k] = OpenCVWrapper.fastXmeans(src, null, kMin, kMax, similarity, colorEpsilon,
new cv.Size(size, size), cv.INTER_LINEAR, bitDepth, attempts);
// Show dst
cv.imshow('dst-canvas', dst);
// Remove src and dst
src.delete();
dst.delete();
最後に
関数名前にXmeansと付けていますが、最初にも言った通りXmeansは無関係なので、使うときは関数名を変えてください。