すあまの備忘録

誰得内容の自分のための非営利目的備忘録ブログ（筆者がわかっても内緒にしてください）

Houdiniの色温度の実装

houdini nuke color

色温度（Blackbody）に関する諸々。

計算式が間違っているかもしれないので要注意。

f:id:godofsuama:20200921170712g:plain

最初に

一時期流行りましたが、物理的に正しいかは？はなんとも言えません。

内部的にK（ケルビン）から計算された値をsRGB空間にリマップしていることが多く、

さらにPyro Blackbodyのように、トーンマップ処理が入っているものもあるので注意が必要です。 f:id:godofsuama:20200921153846p:plain

HoudiniのBlackbody

Houdiniではblackbody関数で得られるのはCIE表色法のXYZ三刺激値となっていて、

内部的にはそれをxyztorgb関数やctransform関数といったものでRGB値に変換しているような気がします。

xyztorgbはLinear sRGB(Rec.709)に変換し、ctransformはリニアNTSC空間となっています。

そうするとACES等でちゃんと使うには色域も変換する必要があるのでは？となります。

特に現状のHoudini(H18.0.499)ではOCIOのWorking Spaceが入力（値、画像）のカラースペースになるようなので手動で変換してあげる必要がありそうです。

xyztorgbで変換

f:id:godofsuama:20200921161049p:plain

//blackbody
v@temp = blackbody(chf("Temperature"),chf("Intensity"));

//set Color
v@Cd = xyztorgb(@temp);

こんな感じで記述するとLinear sRGB空間のRGB値で帰ってきます。これはWhite Pointの補正も込みな気がする。

その後好きなカラースペースにする場合は、ocio_transformを最後の行に追加して

v@Cd = ocio_transform("Utility - Linear - sRGB","変換したいカラースペース名",@Cd);

f:id:godofsuama:20200921161159p:plain

Color Matrixで変換（Houdini）

f:id:godofsuama:20200921161909p:plain

//blackbody
v@temp = blackbody(chf("Temperature"),chf("Intensity"));

//CIE XYZ to sRGB (no white balance)
3@mat3 = { {3.24096994, -1.53738318, -0.49861076}, 
{-0.96924364, 1.7859675, 0.04155506}, 
{0.05563008, -0.20397696, 1.05697151} };

//Multiply Matrix
float red = (@temp.x * getcomp(@mat3,0,0)) + (@temp.y * getcomp(@mat3,0,1)) + (@temp.z * getcomp(@mat3,0,2));
float green = (@temp.x * getcomp(@mat3,1,0)) + (@temp.y * getcomp(@mat3,1,1)) + (@temp.z * getcomp(@mat3,1,2));
float blue = (@temp.x * getcomp(@mat3,2,0)) + (@temp.y * getcomp(@mat3,2,1)) + (@temp.z * getcomp(@mat3,2,2));

//set Color
@Cd = set(red,green,blue);

こんな感じでCIE XYZ to sRGBのカラーマトリックスで計算する。

f:id:godofsuama:20200921163535p:plain

上記のものはWhite Pointの計算を入れていないのでちょっと値がずれるが概ね同じ。

変換のColor Matrixを見つけて入力すればblackbody関数のXYZから直接変換することができる。

RGBtoXYZしか見つからなかった場合は、invert(matrix)で変換してから計算してあげれば同じ結果になる。

Academy Color Encoding System - Wikipedia

Welcome to Bruce Lindbloom's Web Site

Color Matrixで変換（Nuke）

Color Matrixでの計算はNukeでもできる。

f:id:godofsuama:20200921163304p:plain

以下のようにノードを接続し、ConstantのRGBに、Blackbody関数のRGBを入力し、Color Matrixノードに変換のマトリックスを入力する。

f:id:godofsuama:20200921163706p:plain

値が微妙にずれる場合は小数点以下の桁数を増やす。

Pickerの値がHoudiniとNukeで一緒なのに、見た目が合わない場合は同じOCIOを選択し、ViewをRAWにする。

プランキアン軌跡と各カラースペースのPrimaryで実装

プランキアン軌跡（黒体軌跡）と呼ばれる黒体の色の軌跡の式（簡略版）から実装

f:id:godofsuama:20200921170028p:plain

White Pointも考慮したおかげかほぼほぼ同じ値に。

長すぎて2つのWrangleに別れました。

f:id:godofsuama:20200921165122p:plain

f:id:godofsuama:20200921165339p:plain

//Blackbody
float Temp = chf("Temperature");
float TempInt = chf("Intensity");

float tempX = 0;
float tempY = 0;
float tempZ = 0;

//x
if (Temp >= 1667 && Temp < 4000) {
    float xa = -0.2661239 * pow(10,9) / pow(Temp,3);
    float xb = -0.2343589 * pow(10,6) / pow(Temp,2);
    float xc = 0.8776956 * pow(10,3) / Temp;
    
    tempX = (xa + xb + xc) + 0.179910;
    
}else if (Temp >= 4000 && Temp <= 25000) {
    float xa = -3.0258469 * pow(10,9) / pow(Temp,3);
    float xb = 2.1070379 * pow(10,6) / pow(Temp,2);
    float xc = 0.2226347 * pow(10,3) / Temp;
    
    tempX = (xa + xb + xc) + 0.240390;
};

//y
if (Temp >= 1667 && Temp < 2222) {
    float ya = -1.1063814 * pow(tempX,3);
    float yb = -1.3481102 * pow(tempX,2);
    float yc = 2.18555823 * tempX;
    
    tempY = (ya + yb + yc) -0.20219683;
    
}else if (Temp >= 2222 && Temp < 4000) {
    float ya = -0.9549476 * pow(tempX,3);
    float yb = -1.3481102 * pow(tempX,2);
    float yc = 2.09137015 * tempX;
    
    tempY = (ya + yb + yc) -0.16748867;
    
}else if (Temp >= 4000 && Temp <= 25000) {
    float ya = 3.0817580 * pow(tempX,3);
    float yb = -5.8733867 * pow(tempX,2);
    float yc = 3.75112997 * tempX;
    
    tempY = (ya + yb + yc) -0.37001483;
};

float tempXn = (TempInt / tempY) * tempX;
float tempYn = (TempInt / tempY) * (1 - tempX - tempY);

v@Cd = set(tempXn,TempInt,tempYn);

PrimaryからColor Matrixを作成しXYZを変換するWrangle

Primaryは下記ページあたりから探す。

Welcome to Bruce Lindbloom's Web Site

White Pointは下記ページあたり

Standard illuminant - Wikipedia

D60はないのでACESのページのを参照

f:id:godofsuama:20200921165403p:plain

//white point
float wpnX = chf("wpx") / chf("wpy");
float wpnY = chf("wpy") / chf("wpy");
float wpnZ = chf("wpz") / chf("wpy");
vector wpn = set(wpnX,wpnY,wpnZ);

//sort
vector mX = set(chf("priRx"), chf("priGx"), chf("priBx"));
vector mY = set(chf("priRy"), chf("priGy"), chf("priBy"));
vector mZ = set(chf("priRz"), chf("priGz"), chf("priBz"));

//matrix
matrix3 m = set(mX,mY,mZ);
matrix3 mI = invert(m);
matrix3 mW = set(wpnX,wpnX,wpnX, wpnY,wpnY,wpnY, wpnZ,wpnZ, wpnZ);
matrix3 mW2 = mI * mW;

//
vector vecX = set(getcomp(mW2,0,0), 0,0);
vector vecY = set(0, getcomp(mW2,1,1) ,0);
vector vecZ = set(0 ,0 ,getcomp(mW2,2,2));

matrix3 mVec = set(vecX,vecY,vecZ);
3@rgb2xyz = m * mVec;
3@rgb2xyz = invert(@rgb2xyz);


//Color Matrix
float red = (@Cd.r * getcomp(@rgb2xyz,0,0)) + (@Cd.g * getcomp(@rgb2xyz,0,1)) + (@Cd.b * getcomp(@rgb2xyz,0,2));
float green = (@Cd.r * getcomp(@rgb2xyz,1,0)) + (@Cd.g * getcomp(@rgb2xyz,1,1)) + (@Cd.b * getcomp(@rgb2xyz,1,2));
float blue = (@Cd.r * getcomp(@rgb2xyz,2,0)) + (@Cd.g * getcomp(@rgb2xyz,2,1)) + (@Cd.b * getcomp(@rgb2xyz,2,2));

@Cd = set(red,green,blue);

ACES AP0に変換してみた。

Color MatrixはXYZtoACESで、 VEXの方はxyztorgbの後ocio_transformで、 formulaの方はAP0のPrimaryから変換。

どれも微妙にずれるのは精度の問題かもしれない。

もしくは相対色温度の補正値、1.4388/1.438をかける必要がある？？

f:id:godofsuama:20200921173238g:plain