色温度（Blackbody）に関する諸々。
計算式が間違っているかもしれないので要注意。
最初に
一時期流行りましたが、物理的に正しいかは？はなんとも言えません。
内部的にK（ケルビン）から計算された値をsRGB空間にリマップしていることが多く、
さらにPyro Blackbodyのように、トーンマップ処理が入っているものもあるので注意が必要です。
- 最初に
- HoudiniのBlackbody
- xyztorgbで変換
- Color Matrixで変換（Houdini）
- Color Matrixで変換（Nuke）
- プランキアン軌跡と各カラースペースのPrimaryで実装
- ACES AP0に変換してみた。
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で変換
//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);
Color Matrixで変換（Houdini）
//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のカラーマトリックスで計算する。
上記のものは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でもできる。
以下のようにノードを接続し、ConstantのRGBに、Blackbody関数のRGBを入力し、Color Matrixノードに変換のマトリックスを入力する。
値が微妙にずれる場合は小数点以下の桁数を増やす。
Pickerの値がHoudiniとNukeで一緒なのに、見た目が合わない場合は同じOCIOを選択し、ViewをRAWにする。
プランキアン軌跡と各カラースペースのPrimaryで実装
プランキアン軌跡（黒体軌跡）と呼ばれる黒体の色の軌跡の式（簡略版）から実装
White Pointも考慮したおかげかほぼほぼ同じ値に。
長すぎて2つのWrangleに別れました。
//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
//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をかける必要がある？？
