前の記事でのAvisynthでの実験結果を基にして、以前の記事で書いた計算プログラムを修正してみました。
今回計算したのは、以下の２つの内容です。
１．8bit～10bitのYUVについて全ての組み合わせからRGB24の値を計算し、再現できる色数を調べる。
ただし、求められたR,G,Bの値のどれかが0～255の範囲に無い場合は
RGBの範囲外となる組み合わせだとみなし、その組み合わせは無効とする。
（以前の計算では飽和処理をしていた）
ついでに無効となる組み合わせの数も求める。
２．RGB24の16777216色それぞれを一度8bit～10bitのYUVに変換し、
更にそれをRGB24に戻して、再現できる色数を調べる。
ついでに入力と出力が同じ色になる色数も調べる。
■結果
●１について
※「除外数」というのは、RGB範囲外になるので無効とみなされたYUVの組み合わせの数
8bit-Rec601: 再現色数＝[2627304] 除外数＝[8510196]
8bit-PC.601: 再現色数＝[3964583] 除外数＝[12812633]
8bit-Rec709: 再現色数＝[2721568] 除外数＝[8415932]
8bit-PC.709: 再現色数＝[4106090] 除外数＝[12671126]
9bit-Rec601: 再現色数＝[15808872] 除外数＝[67484600]
9bit-PC.601: 再現色数＝[16711693] 除外数＝[102315476]
9bit-Rec709: 再現色数＝[16134243] 除外数＝[66732848]
9bit-PC.709: 再現色数＝[16777216] 除外数＝[101174033]
10bit-Rec601: 再現色数＝[16777216] 除外数＝[537495090]
10bit-PC.601: 再現色数＝[16777216] 除外数＝[817768098]
10bit-Rec709: 再現色数＝[16777216] 除外数＝[531483720]
10bit-PC.709: 再現色数＝[16777216] 除外数＝[808616695]
●２について
※元々RGBから変換したYUV値だからほとんど問題ないだろうということで、
最終的にRGB変換する際には飽和処理を行っている。
※「可逆色数」というのは、入力したRGB値と出力したRGB値が一致したものの数。
※パーセント表示は16777216色に対しての値で、小数第三位以下切捨て。
8bit-Rec601: 再現色数＝[2667904](15.90%) 可逆色数＝[2660636](15.85%)
8bit-PC.601: 再現色数＝[4007380](23.88%) 可逆色数＝[3996769](23.82%)
8bit-Rec709: 再現色数＝[2760514](16.45%) 可逆色数＝[2753782](16.41%)
8bit-PC.709: 再現色数＝[4146794](24.71%) 可逆色数＝[4141028](24.68%)
9bit-Rec601: 再現色数＝[13911017](82.91%) 可逆色数＝[13037997](77.71%)
9bit-PC.601: 再現色数＝[15102283](90.01%) 可逆色数＝[14802234](88.22%)
9bit-Rec709: 再現色数＝[13398699](79.86%) 可逆色数＝[12586505](75.02%)
9bit-PC.709: 再現色数＝[15037350](89.62%) 可逆色数＝[14428999](86.00%)
10bit-Rec601: 再現色数＝[16777142](99.99%) 可逆色数＝[16777142](99.99%)
10bit-PC.601: 再現色数＝[16777216](100.00%) 可逆色数＝[16777216](100.00%)
10bit-Rec709: 再現色数＝[16777216](100.00%) 可逆色数＝[16777216](100.00%)
10bit-PC.709: 再現色数＝[16777216](100.00%) 可逆色数＝[16777216](100.00%)
■結果について
あいかわらず、こんな計算でいいのかどうかはいまいち自信が持てませんが、
ちょっとバラつきがあるものの、8bitについてはだいたいにおいて
Avisynthでの実験結果に近い値が出ている感じですかね。
２については
●9bitではBT.601よりもBT.709のほうが再現色数が少ない。
●10bit-Rec601では16777216色にあと74色足りない。
という結果になっているのがちょっと気になります。
特に前者は何か計算間違えてるんじゃないかとちょっと不安・・・。
今回も以下にソースコードを置いておきますので、おかしなとこがあったら
コメントでツッコミをいただければ幸いです。
■ソースコード（クリックで展開、展開後ソース部分をダブルクリックでコピー可能）
●１のソースコード
/* yuvcolor_calc.cpp */
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#define CLIP_Y(Y) Y = (Y < 0) ? 0.0 : (Y > 1.0) ? 1.0 : Y
#define CLIP_C(C) C = (C < -0.5) ? -0.5 : (C > 0.5) ? 0.5 : C
#define SATURATE(X) X = (X < 0) ? 0 : (X > 255) ? 255 : X
/* RGB24の色数（256*256*256色） */
#define RGB_TABLE_SIZE 16777216
/* H.264のRound()がこんな感じなので */
int myRound(double x)
{
return (x >= 0) ? (int)floor(x + 0.5) : (int)(-floor(abs(x) + 0.5));
}
int main(void)
{
const struct {
char *name; // 色空間の名称
int bitDepth; // ビット深度
int y_min; // 輝度Yの最小値
int y_max; // 輝度Yの最大値
int uv_min; // 色差の最小値
int uv_max; // 色差の最大値
double r_cr; // YUV→RGB変換の係数１
double g_cb; // YUV→RGB変換の係数２
double g_cr; // YUV→RGB変換の係数３
double b_cb; // YUV→RGB変換の係数４
} matrix[] = {
{" 8bit-Rec601", 8, 16, 235, 16, 240, 1.402 , -0.344 , -0.714 , 1.772 },
{" 8bit-PC.601", 8, 0, 255, 0, 255, 1.402 , -0.344 , -0.714 , 1.772 },
{" 8bit-Rec709", 8, 16, 235, 16, 240, 1.5748, -0.1873, -0.4681, 1.8556},
{" 8bit-PC.709", 8, 0, 255, 0, 255, 1.5748, -0.1873, -0.4681, 1.8556},
{" 9bit-Rec601", 9, 32, 470, 32, 480, 1.402 , -0.344 , -0.714 , 1.772 },
{" 9bit-PC.601", 9, 0, 511, 0, 511, 1.402 , -0.344 , -0.714 , 1.772 },
{" 9bit-Rec709", 9, 32, 470, 32, 480, 1.5748, -0.1873, -0.4681, 1.8556},
{" 9bit-PC.709", 9, 0, 511, 0, 511, 1.5748, -0.1873, -0.4681, 1.8556},
{"10bit-Rec601", 10, 64, 940, 64, 960, 1.402 , -0.344 , -0.714 , 1.772 },
{"10bit-PC.601", 10, 0, 1023, 0, 1023, 1.402 , -0.344 , -0.714 , 1.772 },
{"10bit-Rec709", 10, 64, 940, 64, 960, 1.5748, -0.1873, -0.4681, 1.8556},
{"10bit-PC.709", 10, 0, 1023, 0, 1023, 1.5748, -0.1873, -0.4681, 1.8556},
{NULL}
};
bool *rgbTable = (bool *)malloc(RGB_TABLE_SIZE*sizeof(bool));
if(!rgbTable){
fprintf(stderr,"malloc(rgbTable) failed\n");
return -1;
}
double Y_a=0,Cb_a=0,Cr_a=0;
int r=0,g=0,b=0;
int i=0,count=0,exCount=0,num=0;
int rgbNumber;
for(int m=0; m<=11; m++){
printf("%s: ", matrix[m].name);
// rgbTableをクリア(falseにする)
for(i=0;i < RGB_TABLE_SIZE;i++){
rgbTable[i]=false;
}
count=0;
exCount=0;
for (int Y = matrix[m].y_min; Y <= matrix[m].y_max; Y++){
for (int Cb = matrix[m].uv_min; Cb <= matrix[m].uv_max; Cb++){
for (int Cr = matrix[m].uv_min; Cr <= matrix[m].uv_max; Cr++) {
// まずはデジタルYCbCr→アナログYCbCr変換
// (Y_a:0.0～1.0、Cb_a,Cr_a:-0.5～0.5)
// "アナログYCbCr"という表現はよろしくないけどスルーの方向で。
// bitDepthをnとすると、H.264で定義されてる量子化式はだいたいこんな感じなので、
// ここから逆算してY_a、Cb_a、Cr_aを求める。
// ●TVレンジ
// Y = ( 1 << (n-8) ) * (219 * Y_a + 16)
// Cb = ( 1 << (n-8) ) * (224 * Cb_a + 128)
// Cr = ( 1 << (n-8) ) * (224 * Cr_a + 128)
// ●フルレンジ
// Y = ( (1 << n) - 1) * Y_a
// Cb = ( (1 << n) - 1) * Cb_a + ( 1 << (n-1) )
// Cr = ( (1 << n) - 1) * Cr_a + ( 1 << (n-1) )
// minとmaxを用意しとけばTVレンジでもフルレンジでも
// １つの式でいけるので、そのように変形してみた。
Y_a = (Y - matrix[m].y_min) / (double)(matrix[m].y_max - matrix[m].y_min);
Cb_a = (Cb - ( 1 << (matrix[m].bitDepth - 1) ) ) / (double)(matrix[m].uv_max - matrix[m].uv_min);
Cr_a = (Cr - ( 1 << (matrix[m].bitDepth - 1) ) ) / (double)(matrix[m].uv_max - matrix[m].uv_min);
CLIP_Y(Y_a);
CLIP_C(Cb_a);
CLIP_C(Cr_a);
// 続いてアナログYCbCr→デジタルRGB変換
// 基本的には、まるも氏の
// http://www.marumo.ne.jp/db2002_5.htm#15
// の記事を参照。
// ただしこの記事は2002年のものであり、
// BT.709の係数はBT.709-1のものになっている。
// H.264ではBT.709-2以降の係数が使われているので、
// ここではBT.709-2の係数で計算している。
// YUV→RGBのアナログ変換式のみ載せておく。
// (R,G,B,Y:0.0～1.0、U,V:-0.5～0.5)
// ●BT.601
// R = Y ＋ 1.402 × V
// G = Y － 0.344 × U － 0.714 × V
// B = Y ＋ 1.772 × U
// ●BT.709-2
// R = Y ＋ 1.5748 × V
// G = Y － 0.1873 × U － 0.4681 × V
// B = Y ＋ 1.8556 × U
r = myRound(255.0 * (Y_a + matrix[m].r_cr * Cr_a));
g = myRound(255.0 * (Y_a + matrix[m].g_cb * Cb_a + matrix[m].g_cr * Cr_a));
b = myRound(255.0 * (Y_a + matrix[m].b_cb * Cb_a ));
// RGBの値が範囲外になるものは除外する。
if(r<0 || r>255 || g<0 || g>255 || b<0 || b>255){
exCount++;
continue;
}
/*
SATURATE(r);
SATURATE(g);
SATURATE(b);
*/
// 再現できた色はtrueにする
rgbNumber = ((r << 16) | (g << 8) | b);
rgbTable[rgbNumber] = true;
}
}
}
// trueになっている色の数をカウントして出力
num = 0;
for (i = 0; i < RGB_TABLE_SIZE; i++){
if (rgbTable[i]){
num++;
}
}
printf("再現色数＝[%d] 除外数＝[%d]\n", num, exCount);
}
free(rgbTable);
return 0;
}
●２のソースコード
/* yuvcolor_calc2.cpp */
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#define CLIP_Y(Y) Y = (Y < 0) ? 0.0 : (Y > 1.0) ? 1.0 : Y
#define CLIP_C(C) C = (C < -0.5) ? -0.5 : (C > 0.5) ? 0.5 : C
#define SATURATE(X) X = (X < 0) ? 0 : (X > 255) ? 255 : X
/* RGB24の色数（256*256*256色） */
#define RGB_TABLE_SIZE 16777216
/* H.264のRound()がこんな感じなので */
int myRound(double x)
{
return (x >= 0) ? (int)floor(x + 0.5) : (int)(-floor(abs(x) + 0.5));
}
int main(void)
{
const struct {
char *name; // 色空間の名称
int bitDepth; // ビット深度
int y_min; // 輝度Yの最小値
int y_max; // 輝度Yの最大値
int uv_min; // 色差の最小値
int uv_max; // 色差の最大値
double r_cr; // YUV→RGB変換のRを求める式のCrの係数
double g_cb; // YUV→RGB変換のGを求める式のCbの係数
double g_cr; // YUV→RGB変換のGを求める式のCrの係数
double b_cb; // YUV→RGB変換のBを求める式のCbの係数
double y_r; // RGB→YUV変換のYを求める式のRの係数
double y_g; // RGB→YUV変換のYを求める式のGの係数
double y_b; // RGB→YUV変換のYを求める式のBの係数
double cb_r; // RGB→YUV変換のCbを求める式のRの係数
double cb_g; // RGB→YUV変換のCbを求める式のGの係数
double cb_b; // RGB→YUV変換のCbを求める式のBの係数
double cr_r; // RGB→YUV変換のCrを求める式のRの係数
double cr_g; // RGB→YUV変換のCrを求める式のGの係数
double cr_b; // RGB→YUV変換のCrを求める式のBの係数
} matrix[] = {
{" 8bit-Rec601", 8, 16, 235, 16, 240, 1.402 , -0.344 , -0.714 , 1.772 ,
0.299 , 0.587 , 0.114 , -0.169 , -0.331 , 0.5 , 0.5 , -0.419 , -0.081},
{" 8bit-PC.601", 8, 0, 255, 0, 255, 1.402 , -0.344 , -0.714 , 1.772 ,
0.299 , 0.587 , 0.114 , -0.169 , -0.331 , 0.5 , 0.5 , -0.419 , -0.081},
{" 8bit-Rec709", 8, 16, 235, 16, 240, 1.5748, -0.1873, -0.4681, 1.8556 ,
0.2126 , 0.7152 , 0.0722 , -0.1146 , -0.3854 , 0.5 , 0.5 , -0.4542 , -0.0458},
{" 8bit-PC.709", 8, 0, 255, 0, 255, 1.5748, -0.1873, -0.4681, 1.8556 ,
0.2126 , 0.7152 , 0.0722 , -0.1146 , -0.3854 , 0.5 , 0.5 , -0.4542 , -0.0458},
{" 9bit-Rec601", 9, 32, 470, 32, 480, 1.402 , -0.344 , -0.714 , 1.772 ,
0.299 , 0.587 , 0.114 , -0.169 , -0.331 , 0.5 , 0.5 , -0.419 , -0.081},
{" 9bit-PC.601", 9, 0, 511, 0, 511, 1.402 , -0.344 , -0.714 , 1.772 ,
0.299 , 0.587 , 0.114 , -0.169 , -0.331 , 0.5 , 0.5 , -0.419 , -0.081},
{" 9bit-Rec709", 9, 32, 470, 32, 480, 1.5748, -0.1873, -0.4681, 1.8556 ,
0.2126 , 0.7152 , 0.0722 , -0.1146 , -0.3854 , 0.5 , 0.5 , -0.4542 , -0.0458},
{" 9bit-PC.709", 9, 0, 511, 0, 511, 1.5748, -0.1873, -0.4681, 1.8556 ,
0.2126 , 0.7152 , 0.0722 , -0.1146 , -0.3854 , 0.5 , 0.5 , -0.4542 , -0.0458},
{"10bit-Rec601", 10, 64, 940, 64, 960, 1.402 , -0.344 , -0.714 , 1.772 ,
0.299 , 0.587 , 0.114 , -0.169 , -0.331 , 0.5 , 0.5 , -0.419 , -0.081},
{"10bit-PC.601", 10, 0, 1023, 0, 1023, 1.402 , -0.344 , -0.714 , 1.772 ,
0.299 , 0.587 , 0.114 , -0.169 , -0.331 , 0.5 , 0.5 , -0.419 , -0.081},
{"10bit-Rec709", 10, 64, 940, 64, 960, 1.5748, -0.1873, -0.4681, 1.8556 ,
0.2126 , 0.7152 , 0.0722 , -0.1146 , -0.3854 , 0.5 , 0.5 , -0.4542 , -0.0458},
{"10bit-PC.709", 10, 0, 1023, 0, 1023, 1.5748, -0.1873, -0.4681, 1.8556 ,
0.2126 , 0.7152 , 0.0722 , -0.1146 , -0.3854 , 0.5 , 0.5 , -0.4542 , -0.0458},
{NULL}
};
bool *rgbTable = (bool *)malloc(RGB_TABLE_SIZE*sizeof(bool));
if(!rgbTable){
fprintf(stderr,"malloc(rgbTable) failed\n");
return -1;
}
double Y_a=0,Cb_a=0,Cr_a=0;
int Y=0,Cb=0,Cr=0;
int r=0,g=0,b=0;
int i=0,count=0,num=0,equalCount=0;
int rgbNumber;
for(int m=0; m<=11; m++){
printf("%s: ", matrix[m].name);
// rgbTableをクリア(falseにする)
for(i=0;i<RGB_TABLE_SIZE;i++){
rgbTable[i]=false;
}
count=0;
equalCount=0;
for (int r0 = 0; r0 <= 255; r0++){
for (int g0 = 0; g0 <= 255; g0++){
for (int b0 = 0; b0 <= 255; b0++) {
// まずはデジタルRGBからアナログYCbCr(Y_a:0.0～1.0、Cb_a,Cr_a:-0.5～0.5)を求める
// "アナログYCbCr"という表現はよろしくないけどスルーの方向で。
// 基本的には、まるも氏の
// http://www.marumo.ne.jp/db2002_5.htm#15
// の記事を参照。
// ただしこの記事は2002年のものであり、BT.709の係数はBT.709-1のものになっている。
// H.264ではBT.709-2以降の係数が使われているので、ここではBT.709-2の係数で計算している。
// RGB→YUVのアナログ変換式は以下のとおり。(R,G,B,Y:0.0～1.0、U,V:-0.5～0.5)
// ●BT.601
// Y = 0.299 × R ＋ 0.587 × G ＋ 0.114 × B
// U = －0.169 × R － 0.331 × G ＋ 0.500 × B
// V = 0.500 × R － 0.419 × G － 0.081 × B
// ●BT.709-2
// Y = 0.2126 × R ＋ 0.7152 × G ＋ 0.0722 × B
// U = －0.1146 × R － 0.3854 × G ＋ 0.5000 × B
// V = 0.5000 × R － 0.4542 × G － 0.0458 × B
//
Y_a = (matrix[m].y_r * r0 + matrix[m].y_g * g0 + matrix[m].y_b * b0)/255.0;
Cb_a = (matrix[m].cb_r * r0 + matrix[m].cb_g * g0 + matrix[m].cb_b * b0)/255.0;
Cr_a = (matrix[m].cr_r * r0 + matrix[m].cr_g * g0 + matrix[m].cr_b * b0)/255.0;
// あまり意味はないと思うけど一応飽和処理
CLIP_Y(Y_a);
CLIP_C(Cb_a);
CLIP_C(Cr_a);
// 次にアナログYCbCrを一度デジタルYCbCrに変換する
// bitDepthをnとすると、H.264で定義されてる量子化式はだいたいこんな感じ。
// ●TVレンジ
// Y = ( 1 << (n-8) ) * (219 * Y_a + 16)
// Cb = ( 1 << (n-8) ) * (224 * Cb_a + 128)
// Cr = ( 1 << (n-8) ) * (224 * Cr_a + 128)
// ●フルレンジ
// Y = ( (1<<n) - 1) * Y_a
// Cb = ( (1<<n) - 1) * Cb_a + ( 1 << (n-1) )
// Cr = ( (1<<n) - 1) * Cr_a + ( 1 << (n-1) )
//
// minとmaxを用意しとけばTVレンジでもフルレンジでも１つの式でいけるので変形してみた。
//
Y = myRound((matrix[m].y_max - matrix[m].y_min) * Y_a + matrix[m].y_min);
Cb = myRound((matrix[m].uv_max - matrix[m].uv_min) * Cb_a + ( 1 << (matrix[m].bitDepth - 1) ));
Cr = myRound((matrix[m].uv_max - matrix[m].uv_min) * Cr_a + ( 1 << (matrix[m].bitDepth - 1) ));
// あまり意味はないと思うけど一応飽和処理
Y = (Y < matrix[m].y_min) ? matrix[m].y_min : (Y > matrix[m].y_max) ? matrix[m].y_max : Y;
Cb = (Cb < matrix[m].uv_min) ? matrix[m].uv_min : (Cb > matrix[m].uv_max) ? matrix[m].uv_max : Cb;
Cr = (Cr < matrix[m].uv_min) ? matrix[m].uv_min : (Cr > matrix[m].uv_max) ? matrix[m].uv_max : Cr;
// 次に、デジタルYCbCrをあらためてアナログYCbCrに変換する
Y_a = (Y - matrix[m].y_min) / (double)(matrix[m].y_max - matrix[m].y_min);
Cb_a = (Cb - ( 1 << (matrix[m].bitDepth - 1) ) ) / (double)(matrix[m].uv_max - matrix[m].uv_min);
Cr_a = (Cr - ( 1 << (matrix[m].bitDepth - 1) ) ) / (double)(matrix[m].uv_max - matrix[m].uv_min);
// あまり意味はないと思うけど一応飽和処理
CLIP_Y(Y_a);
CLIP_C(Cb_a);
CLIP_C(Cr_a);
// 続いてアナログYCbCr→デジタルRGB変換
// 基本的には、まるも氏の
// http://www.marumo.ne.jp/db2002_5.htm#15
// の記事を参照。
// ただしこの記事は2002年のものであり、BT.709の係数はBT.709-1のものになっている。
// H.264ではBT.709-2以降の係数が使われているので、ここではBT.709-2の係数で計算している。
// YUV→RGBのアナログ変換式は以下のとおり。(R,G,B,Y:0.0～1.0、U,V:-0.5～0.5)
// ●BT.601
// R = Y ＋ 1.402 × V
// G = Y － 0.344 × U － 0.714 × V
// B = Y ＋ 1.772 × U
// ●BT.709-2
// R = Y ＋ 1.5748 × V
// G = Y － 0.1873 × U － 0.4681 × V
// B = Y ＋ 1.8556 × U
r = myRound(255.0 * (Y_a + matrix[m].r_cr * Cr_a));
g = myRound(255.0 * (Y_a + matrix[m].g_cb * Cb_a + matrix[m].g_cr * Cr_a));
b = myRound(255.0 * (Y_a + matrix[m].b_cb * Cb_a ));
// RGBの値が範囲外になるものは除外する。
if(r<0 || r>255 || g<0 || g>255 || b<0 || b>255) continue;
/*
// 飽和処理
SATURATE(r);
SATURATE(g);
SATURATE(b);
*/
// 入力したRGBと、変換後のRGBとが同じになるものをカウントする
if((r0==r) && (g0==g) && (b0==b)) equalCount++;
// 再現できた色はtrueにする
rgbNumber = ((r << 16) | (g << 8) | b);
rgbTable[rgbNumber] = true;
}
}
}
// trueになっている色の数をカウントして出力
num = 0;
for (i = 0; i < RGB_TABLE_SIZE; i++){
if (rgbTable[i]){
num++;
}
}
printf("再現色数＝[%d](%02.2f%%) 可逆色数＝[%d](%02.2f%%)\n", num, floor(num/(double)RGB_TABLE_SIZE*10000)/100.0, equalCount , floor(equalCount/(double)RGB_TABLE_SIZE*10000)/100.0);
}
free(rgbTable);
return 0;
}
最近のコメント