1482 lines
43 KiB
HLSL
1482 lines
43 KiB
HLSL
// Copyright Epic Games, Inc. All Rights Reserved.
|
|
|
|
/*
|
|
=============================================================================
|
|
ACES: Academy Color Encoding System
|
|
https://github.com/ampas/aces-dev/tree/v1.3
|
|
|
|
License Terms for Academy Color Encoding System Components
|
|
|
|
Academy Color Encoding System (ACES) software and tools are provided by the Academy under
|
|
the following terms and conditions: A worldwide, royalty-free, non-exclusive right to copy, modify, create
|
|
derivatives, and use, in source and binary forms, is hereby granted, subject to acceptance of this license.
|
|
|
|
Copyright © 2015 Academy of Motion Picture Arts and Sciences (A.M.P.A.S.). Portions contributed by
|
|
others as indicated. All rights reserved.
|
|
|
|
Performance of any of the aforementioned acts indicates acceptance to be bound by the following
|
|
terms and conditions:
|
|
|
|
* Copies of source code, in whole or in part, must retain the above copyright
|
|
notice, this list of conditions and the Disclaimer of Warranty.
|
|
* Use in binary form must retain the above copyright notice, this list of
|
|
conditions and the Disclaimer of Warranty in the documentation and/or other
|
|
materials provided with the distribution.
|
|
* Nothing in this license shall be deemed to grant any rights to trademarks,
|
|
copyrights, patents, trade secrets or any other intellectual property of
|
|
A.M.P.A.S. or any contributors, except as expressly stated herein.
|
|
* Neither the name "A.M.P.A.S." nor the name of any other contributors to this
|
|
software may be used to endorse or promote products derivative of or based on
|
|
this software without express prior written permission of A.M.P.A.S. or the
|
|
contributors, as appropriate.
|
|
|
|
This license shall be construed pursuant to the laws of the State of California,
|
|
and any disputes related thereto shall be subject to the jurisdiction of the courts therein.
|
|
|
|
Disclaimer of Warranty: THIS SOFTWARE IS PROVIDED BY A.M.P.A.S. AND CONTRIBUTORS "AS
|
|
IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
|
|
IMPLIED WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE, AND
|
|
NON-INFRINGEMENT ARE DISCLAIMED. IN NO EVENT SHALL A.M.P.A.S., OR ANY
|
|
CONTRIBUTORS OR DISTRIBUTORS, BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
|
|
SPECIAL, EXEMPLARY, RESITUTIONARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT
|
|
NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
|
|
DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY
|
|
OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
|
|
NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE,
|
|
EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
|
|
|
|
WITHOUT LIMITING THE GENERALITY OF THE FOREGOING, THE ACADEMY SPECIFICALLY
|
|
DISCLAIMS ANY REPRESENTATIONS OR WARRANTIES WHATSOEVER RELATED TO PATENT OR
|
|
OTHER INTELLECTUAL PROPERTY RIGHTS IN THE ACADEMY COLOR ENCODING SYSTEM, OR
|
|
APPLICATIONS THEREOF, HELD BY PARTIES OTHER THAN A.M.P.A.S.,WHETHER DISCLOSED
|
|
OR UNDISCLOSED.
|
|
=============================================================================
|
|
*/
|
|
|
|
/*
|
|
we do matrix/matrix and matrix/vectors multiplications are following textbook order, like in our current implementation in ACES1.0.
|
|
An example:
|
|
in ACES:
|
|
const float AP0_2_AP1_MAT[4][4] = mult_f44_f44( AP0_2_XYZ_MAT, XYZ_2_AP1_MAT);
|
|
const float RGB_out[3] = mult_f3_f44( RGB_in, AP0_2_AP1_MAT);
|
|
Equivalent math:
|
|
R_out = 1.4514393161 * R_in + -0.2365107469 * G_in + -0.2149285693 * B_in;
|
|
G_out = -0.0765537734 * R_in + 1.1762296998 * G_in + -0.0996759264 * B_in;
|
|
B_out = 0.0083161484 * R_in + -0.0060324498 * G_in + 0.9977163014 * B_in;
|
|
|
|
in HLSL:
|
|
static const float3x3 AP0_2_AP1_MAT = //mul( XYZ_2_AP1_MAT, AP0_2_XYZ_MAT );
|
|
{
|
|
1.4514393161, -0.2365107469, -0.2149285693,
|
|
-0.0765537734, 1.1762296998, -0.0996759264,
|
|
0.0083161484, -0.0060324498, 0.9977163014,
|
|
};
|
|
float3 RGB_out = mul(AP0_2_AP1_MAT, RGB_in.xyz);
|
|
*/
|
|
|
|
/////////////////////////////////////////////////////////////////////////////////////////
|
|
|
|
#pragma once
|
|
|
|
#include "ACESCommon.ush"
|
|
|
|
///////////////////////////////////////////////////////////////////////////////////////////
|
|
|
|
static const float HALF_MAX = 65504.0;
|
|
static const float HALF_POS_INF = 65535.0f;
|
|
|
|
static const float MIN_STOP_SDR = -6.5;
|
|
static const float MAX_STOP_SDR = 6.5;
|
|
|
|
static const float MIN_STOP_RRT = -15.;
|
|
static const float MAX_STOP_RRT = 18.;
|
|
|
|
static const float MIN_LUM_SDR = 0.02;
|
|
static const float MAX_LUM_SDR = 48.0;
|
|
|
|
static const float MIN_LUM_RRT = 0.0001;
|
|
static const float MAX_LUM_RRT = 10000.0;
|
|
|
|
float pow10(float x)
|
|
{
|
|
return pow(10, x);
|
|
}
|
|
|
|
// ------- Red modifier functions
|
|
float cubic_basis_shaper
|
|
(
|
|
float x,
|
|
float w // full base width of the shaper function (in degrees)
|
|
)
|
|
{
|
|
float M[4][4] = { { -1. / 6, 3. / 6, -3. / 6, 1. / 6 },
|
|
{ 3. / 6, -6. / 6, 3. / 6, 0. / 6 },
|
|
{ -3. / 6, 0. / 6, 3. / 6, 0. / 6 },
|
|
{ 1. / 6, 4. / 6, 1. / 6, 0. / 6 } };
|
|
|
|
float knots[5] = { -0.5 * w, -0.25 * w, 0, 0.25 * w, 0.5 * w };
|
|
|
|
float y = 0;
|
|
if ((x > knots[0]) && (x < knots[4])) {
|
|
float knot_coord = (x - knots[0]) * 4.0 / w;
|
|
int j = knot_coord;
|
|
float t = knot_coord - j;
|
|
|
|
float monomials[4] = { t * t * t, t * t, t, 1.0 };
|
|
|
|
// (if/else structure required for compatibility with CTL < v1.5.)
|
|
if (j == 3) {
|
|
y = monomials[0] * M[0][0] + monomials[1] * M[1][0] +
|
|
monomials[2] * M[2][0] + monomials[3] * M[3][0];
|
|
}
|
|
else if (j == 2) {
|
|
y = monomials[0] * M[0][1] + monomials[1] * M[1][1] +
|
|
monomials[2] * M[2][1] + monomials[3] * M[3][1];
|
|
}
|
|
else if (j == 1) {
|
|
y = monomials[0] * M[0][2] + monomials[1] * M[1][2] +
|
|
monomials[2] * M[2][2] + monomials[3] * M[3][2];
|
|
}
|
|
else if (j == 0) {
|
|
y = monomials[0] * M[0][3] + monomials[1] * M[1][3] +
|
|
monomials[2] * M[2][3] + monomials[3] * M[3][3];
|
|
}
|
|
else {
|
|
y = 0.0;
|
|
}
|
|
}
|
|
|
|
return y * 1.5;
|
|
}
|
|
|
|
// note: this matrix is not transposed unlike other matrix operations because instead of doing v_out = mul(M, v_in), we do v_out = mul(v_in, M)
|
|
// Textbook monomial to basis-function conversion matrix.
|
|
static const float3x3 M1 =
|
|
{
|
|
{ 0.5, -1.0, 0.5 },
|
|
{ -1.0, 1.0, 0.5 },
|
|
{ 0.5, 0.0, 0.0 }
|
|
};
|
|
|
|
|
|
/* ---- Functions to compress highlights ---- */
|
|
// allow for simulated white points without clipping
|
|
|
|
float roll_white_fwd(
|
|
float inValue, // color value to adjust (white scaled to around 1.0)
|
|
float new_wht, // white adjustment (e.g. 0.9 for 10% darkening)
|
|
float width // adjusted width (e.g. 0.25 for top quarter of the tone scale)
|
|
)
|
|
{
|
|
const float x0 = -1.0;
|
|
const float x1 = x0 + width;
|
|
const float y0 = -new_wht;
|
|
const float y1 = x1;
|
|
const float m1 = (x1 - x0);
|
|
const float a = y0 - y1 + m1;
|
|
const float b = 2 * (y1 - y0) - m1;
|
|
const float c = y0;
|
|
const float t = (-inValue - x0) / (x1 - x0);
|
|
float outValue = 0.0;
|
|
if (t < 0.0)
|
|
outValue = -(t * b + c);
|
|
else if (t > 1.0)
|
|
outValue = inValue;
|
|
else
|
|
outValue = -((t * a + b) * t + c);
|
|
return outValue;
|
|
}
|
|
|
|
float roll_white_rev(
|
|
float inValue, // color value to adjust (white scaled to around 1.0)
|
|
float new_wht, // white adjustment (e.g. 0.9 for 10% darkening)
|
|
float width // adjusted width (e.g. 0.25 for top quarter of the tone scale)
|
|
)
|
|
{
|
|
const float x0 = -1.0;
|
|
const float x1 = x0 + width;
|
|
const float y0 = -new_wht;
|
|
const float y1 = x1;
|
|
const float m1 = (x1 - x0);
|
|
const float a = y0 - y1 + m1;
|
|
const float b = 2. * (y1 - y0) - m1;
|
|
float c = y0;
|
|
float outValue = 0.0;
|
|
if (-inValue < y0)
|
|
outValue = -x0;
|
|
else if (-inValue > y1)
|
|
outValue = inValue;
|
|
else {
|
|
c = c + inValue;
|
|
const float discrim = sqrt(b * b - 4. * a * c);
|
|
const float t = (2. * c) / (-discrim - b);
|
|
outValue = -((t * (x1 - x0)) + x0);
|
|
}
|
|
return outValue;
|
|
}
|
|
|
|
float3 limit_to_primaries
|
|
(
|
|
float3 XYZ,
|
|
FColorSpace LIMITING_PRI
|
|
)
|
|
{
|
|
float3x3 XYZ_2_LIMITING_PRI_MAT = LIMITING_PRI.XYZtoRGB;
|
|
float3x3 LIMITING_PRI_2_XYZ_MAT = LIMITING_PRI.RGBtoXYZ;
|
|
|
|
// XYZ to limiting primaries
|
|
float3 rgb = mul(XYZ_2_LIMITING_PRI_MAT, XYZ);
|
|
|
|
// Clip any values outside the limiting primaries
|
|
float3 limitedRgb = clamp(rgb, 0.0.xxx, 1.0.xxx);
|
|
|
|
// Convert limited RGB to XYZ
|
|
return mul(LIMITING_PRI_2_XYZ_MAT, limitedRgb);
|
|
}
|
|
|
|
///+ TODO: check if valid
|
|
float interpolate1D(const float2 table[2], float value)
|
|
{
|
|
float t = saturate((value - table[0].x) / (table[1].x - table[0].x));
|
|
return lerp(table[0].y, table[1].y, t);
|
|
}
|
|
///-
|
|
|
|
float lookup_ACESmin(float minLum)
|
|
{
|
|
const float2 minTable[2] = { float2(log10(MIN_LUM_RRT), MIN_STOP_RRT),
|
|
float2(log10(MIN_LUM_SDR), MIN_STOP_SDR) };
|
|
|
|
return 0.18 * pow(2., interpolate1D(minTable, log10(minLum)));
|
|
}
|
|
|
|
float lookup_ACESmax(float maxLum)
|
|
{
|
|
const float2 maxTable[2] = { float2(log10(MAX_LUM_SDR), MAX_STOP_SDR),
|
|
float2(log10(MAX_LUM_RRT), MAX_STOP_RRT) };
|
|
|
|
return 0.18 * pow(2., interpolate1D(maxTable, log10(maxLum)));
|
|
}
|
|
|
|
struct TsPoint
|
|
{
|
|
float x; // ACES
|
|
float y; // luminance
|
|
float slope; //
|
|
};
|
|
|
|
struct TsParams
|
|
{
|
|
TsPoint Min;
|
|
TsPoint Mid;
|
|
TsPoint Max;
|
|
float coefsLow[6];
|
|
float coefsHigh[6];
|
|
};
|
|
|
|
void init_coefsLow(
|
|
TsPoint TsPointLow,
|
|
TsPoint TsPointMid,
|
|
out float coefsLow[5]
|
|
)
|
|
{
|
|
|
|
float knotIncLow = (log10(TsPointMid.x) - log10(TsPointLow.x)) / 3.;
|
|
// float halfKnotInc = (log10(TsPointMid.x) - log10(TsPointLow.x)) / 6.;
|
|
|
|
// Determine two lowest coefficients (straddling minPt)
|
|
coefsLow[0] = (TsPointLow.slope * (log10(TsPointLow.x) - 0.5 * knotIncLow)) + (log10(TsPointLow.y) - TsPointLow.slope * log10(TsPointLow.x));
|
|
coefsLow[1] = (TsPointLow.slope * (log10(TsPointLow.x) + 0.5 * knotIncLow)) + (log10(TsPointLow.y) - TsPointLow.slope * log10(TsPointLow.x));
|
|
// NOTE: if slope=0, then the above becomes just
|
|
// coefsLow[0] = log10(TsPointLow.y);
|
|
// coefsLow[1] = log10(TsPointLow.y);
|
|
// leaving it as a variable for now in case we decide we need non-zero slope extensions
|
|
|
|
// Determine two highest coefficients (straddling midPt)
|
|
coefsLow[3] = (TsPointMid.slope * (log10(TsPointMid.x) - 0.5 * knotIncLow)) + (log10(TsPointMid.y) - TsPointMid.slope * log10(TsPointMid.x));
|
|
coefsLow[4] = (TsPointMid.slope * (log10(TsPointMid.x) + 0.5 * knotIncLow)) + (log10(TsPointMid.y) - TsPointMid.slope * log10(TsPointMid.x));
|
|
|
|
// Middle coefficient (which defines the "sharpness of the bend") is linearly interpolated
|
|
float2 bendsLow[2] = { float2(MIN_STOP_RRT, 0.18),
|
|
float2(MIN_STOP_SDR, 0.35) };
|
|
float pctLow = interpolate1D(bendsLow, log2(TsPointLow.x / 0.18));
|
|
coefsLow[2] = log10(TsPointLow.y) + pctLow * (log10(TsPointMid.y) - log10(TsPointLow.y));
|
|
}
|
|
|
|
void init_coefsHigh(
|
|
TsPoint TsPointMid,
|
|
TsPoint TsPointMax,
|
|
out float coefsHigh[5]
|
|
)
|
|
{
|
|
|
|
float knotIncHigh = (log10(TsPointMax.x) - log10(TsPointMid.x)) / 3.;
|
|
// float halfKnotInc = (log10(TsPointMax.x) - log10(TsPointMid.x)) / 6.;
|
|
|
|
// Determine two lowest coefficients (straddling midPt)
|
|
coefsHigh[0] = (TsPointMid.slope * (log10(TsPointMid.x) - 0.5 * knotIncHigh)) + (log10(TsPointMid.y) - TsPointMid.slope * log10(TsPointMid.x));
|
|
coefsHigh[1] = (TsPointMid.slope * (log10(TsPointMid.x) + 0.5 * knotIncHigh)) + (log10(TsPointMid.y) - TsPointMid.slope * log10(TsPointMid.x));
|
|
|
|
// Determine two highest coefficients (straddling maxPt)
|
|
coefsHigh[3] = (TsPointMax.slope * (log10(TsPointMax.x) - 0.5 * knotIncHigh)) + (log10(TsPointMax.y) - TsPointMax.slope * log10(TsPointMax.x));
|
|
coefsHigh[4] = (TsPointMax.slope * (log10(TsPointMax.x) + 0.5 * knotIncHigh)) + (log10(TsPointMax.y) - TsPointMax.slope * log10(TsPointMax.x));
|
|
// NOTE: if slope=0, then the above becomes just
|
|
// coefsHigh[0] = log10(TsPointHigh.y);
|
|
// coefsHigh[1] = log10(TsPointHigh.y);
|
|
// leaving it as a variable for now in case we decide we need non-zero slope extensions
|
|
|
|
// Middle coefficient (which defines the "sharpness of the bend") is linearly interpolated
|
|
float2 bendsHigh[2] = { float2(MAX_STOP_SDR, 0.89) ,
|
|
float2(MAX_STOP_RRT, 0.90) };
|
|
float pctHigh = interpolate1D(bendsHigh, log2(TsPointMax.x / 0.18));
|
|
coefsHigh[2] = log10(TsPointMid.y) + pctHigh * (log10(TsPointMax.y) - log10(TsPointMid.y));
|
|
}
|
|
|
|
float shift(float inValue, float expShift)
|
|
{
|
|
return pow(2., (log2(inValue) - expShift));
|
|
}
|
|
|
|
TsParams init_TsParams(
|
|
float minLum,
|
|
float maxLum,
|
|
float expShift = 0
|
|
)
|
|
{
|
|
TsPoint MIN_PT = { lookup_ACESmin(minLum), minLum, 0.0 };
|
|
TsPoint MID_PT = { 0.18, 4.8, 1.55 };
|
|
TsPoint MAX_PT = { lookup_ACESmax(maxLum), maxLum, 0.0 };
|
|
float cLow[5];
|
|
init_coefsLow(MIN_PT, MID_PT, cLow);
|
|
float cHigh[5];
|
|
init_coefsHigh(MID_PT, MAX_PT, cHigh);
|
|
MIN_PT.x = shift(lookup_ACESmin(minLum), expShift);
|
|
MID_PT.x = shift(0.18, expShift);
|
|
MAX_PT.x = shift(lookup_ACESmax(maxLum), expShift);
|
|
|
|
TsParams P = {
|
|
{MIN_PT.x, MIN_PT.y, MIN_PT.slope},
|
|
{MID_PT.x, MID_PT.y, MID_PT.slope},
|
|
{MAX_PT.x, MAX_PT.y, MAX_PT.slope},
|
|
{cLow[0], cLow[1], cLow[2], cLow[3], cLow[4], cLow[4]},
|
|
{cHigh[0], cHigh[1], cHigh[2], cHigh[3], cHigh[4], cHigh[4]}
|
|
};
|
|
|
|
return P;
|
|
}
|
|
|
|
float ssts
|
|
(
|
|
const float x,
|
|
const TsParams C
|
|
)
|
|
{
|
|
const int N_KNOTS_LOW = 4;
|
|
const int N_KNOTS_HIGH = 4;
|
|
|
|
// Check for negatives or zero before taking the log. If negative or zero,
|
|
// set to HALF_MIN
|
|
float logx = log10(max(x, 1e-10));
|
|
float logy;
|
|
|
|
if (logx <= log10(C.Min.x)) {
|
|
|
|
logy = logx * C.Min.slope + (log10(C.Min.y) - C.Min.slope * log10(C.Min.x));
|
|
|
|
}
|
|
else if ((logx > log10(C.Min.x)) && (logx < log10(C.Mid.x))) {
|
|
|
|
float knot_coord = (N_KNOTS_LOW - 1) * (logx - log10(C.Min.x)) / (log10(C.Mid.x) - log10(C.Min.x));
|
|
int j = knot_coord;
|
|
float t = knot_coord - j;
|
|
|
|
float3 cf = { C.coefsLow[j], C.coefsLow[j + 1], C.coefsLow[j + 2] };
|
|
|
|
float3 monomials = { t * t, t, 1.0 };
|
|
logy = dot(monomials, mul(cf, M1));
|
|
|
|
}
|
|
else if ((logx >= log10(C.Mid.x)) && (logx < log10(C.Max.x))) {
|
|
|
|
float knot_coord = (N_KNOTS_HIGH - 1) * (logx - log10(C.Mid.x)) / (log10(C.Max.x) - log10(C.Mid.x));
|
|
int j = knot_coord;
|
|
float t = knot_coord - j;
|
|
|
|
float3 cf = { C.coefsHigh[j], C.coefsHigh[j + 1], C.coefsHigh[j + 2] };
|
|
|
|
float3 monomials = { t * t, t, 1.0 };
|
|
logy = dot(monomials, mul(cf, M1));
|
|
|
|
}
|
|
else { //if ( logIn >= log10(C.Max.x) ) {
|
|
|
|
logy = logx * C.Max.slope + (log10(C.Max.y) - C.Max.slope * log10(C.Max.x));
|
|
|
|
}
|
|
|
|
return pow10(logy);
|
|
|
|
}
|
|
|
|
float inv_ssts
|
|
(
|
|
const float y,
|
|
const TsParams C
|
|
)
|
|
{
|
|
const int N_KNOTS_LOW = 4;
|
|
const int N_KNOTS_HIGH = 4;
|
|
|
|
const float KNOT_INC_LOW = (log10(C.Mid.x) - log10(C.Min.x)) / (N_KNOTS_LOW - 1.);
|
|
const float KNOT_INC_HIGH = (log10(C.Max.x) - log10(C.Mid.x)) / (N_KNOTS_HIGH - 1.);
|
|
|
|
// KNOT_Y is luminance of the spline at each knot
|
|
float KNOT_Y_LOW[N_KNOTS_LOW];
|
|
///+warning: redefinition of 'i'
|
|
{
|
|
for (int i = 0; i < N_KNOTS_LOW; i = i + 1) {
|
|
KNOT_Y_LOW[i] = (C.coefsLow[i] + C.coefsLow[i + 1]) / 2.;
|
|
};
|
|
}
|
|
///-
|
|
|
|
float KNOT_Y_HIGH[N_KNOTS_HIGH];
|
|
///+ warning: redefinition of 'i'
|
|
{
|
|
for (int i = 0; i < N_KNOTS_HIGH; i = i + 1) {
|
|
KNOT_Y_HIGH[i] = (C.coefsHigh[i] + C.coefsHigh[i + 1]) / 2.;
|
|
};
|
|
}
|
|
///-
|
|
float logy = log10(max(y, 1e-10));
|
|
|
|
float logx;
|
|
if (logy <= log10(C.Min.y)) {
|
|
|
|
logx = log10(C.Min.x);
|
|
|
|
}
|
|
else if ((logy > log10(C.Min.y)) && (logy <= log10(C.Mid.y))) {
|
|
|
|
int j;
|
|
float3 cf;
|
|
if (logy > KNOT_Y_LOW[0] && logy <= KNOT_Y_LOW[1]) {
|
|
cf[0] = C.coefsLow[0]; cf[1] = C.coefsLow[1]; cf[2] = C.coefsLow[2]; j = 0;
|
|
}
|
|
else if (logy > KNOT_Y_LOW[1] && logy <= KNOT_Y_LOW[2]) {
|
|
cf[0] = C.coefsLow[1]; cf[1] = C.coefsLow[2]; cf[2] = C.coefsLow[3]; j = 1;
|
|
}
|
|
else if (logy > KNOT_Y_LOW[2] && logy <= KNOT_Y_LOW[3]) {
|
|
cf[0] = C.coefsLow[2]; cf[1] = C.coefsLow[3]; cf[2] = C.coefsLow[4]; j = 2;
|
|
}
|
|
|
|
const float3 tmp = mul(cf, M1);
|
|
|
|
float a = tmp[0];
|
|
float b = tmp[1];
|
|
float c = tmp[2];
|
|
c = c - logy;
|
|
|
|
const float d = sqrt(b * b - 4. * a * c);
|
|
|
|
const float t = (2. * c) / (-d - b);
|
|
|
|
logx = log10(C.Min.x) + (t + j) * KNOT_INC_LOW;
|
|
|
|
}
|
|
else if ((logy > log10(C.Mid.y)) && (logy < log10(C.Max.y))) {
|
|
|
|
int j;
|
|
float3 cf;
|
|
if (logy >= KNOT_Y_HIGH[0] && logy <= KNOT_Y_HIGH[1]) {
|
|
cf[0] = C.coefsHigh[0]; cf[1] = C.coefsHigh[1]; cf[2] = C.coefsHigh[2]; j = 0;
|
|
}
|
|
else if (logy > KNOT_Y_HIGH[1] && logy <= KNOT_Y_HIGH[2]) {
|
|
cf[0] = C.coefsHigh[1]; cf[1] = C.coefsHigh[2]; cf[2] = C.coefsHigh[3]; j = 1;
|
|
}
|
|
else if (logy > KNOT_Y_HIGH[2] && logy <= KNOT_Y_HIGH[3]) {
|
|
cf[0] = C.coefsHigh[2]; cf[1] = C.coefsHigh[3]; cf[2] = C.coefsHigh[4]; j = 2;
|
|
}
|
|
|
|
const float3 tmp = mul(cf, M1);
|
|
|
|
float a = tmp[0];
|
|
float b = tmp[1];
|
|
float c = tmp[2];
|
|
c = c - logy;
|
|
|
|
const float d = sqrt(b * b - 4. * a * c);
|
|
|
|
const float t = (2. * c) / (-d - b);
|
|
|
|
logx = log10(C.Mid.x) + (t + j) * KNOT_INC_HIGH;
|
|
|
|
}
|
|
else { //if ( logy >= log10(C.Max.y) ) {
|
|
|
|
logx = log10(C.Max.x);
|
|
|
|
}
|
|
|
|
return pow10(logx);
|
|
|
|
}
|
|
|
|
float3 ssts_f3
|
|
(
|
|
const float3 x,
|
|
const TsParams C
|
|
)
|
|
{
|
|
float3 outValue;
|
|
outValue[0] = ssts(x[0], C);
|
|
outValue[1] = ssts(x[1], C);
|
|
outValue[2] = ssts(x[2], C);
|
|
|
|
return outValue;
|
|
}
|
|
|
|
float3 inv_ssts_f3
|
|
(
|
|
const float3 x,
|
|
const TsParams C
|
|
)
|
|
{
|
|
float3 outValue;
|
|
outValue[0] = inv_ssts(x[0], C);
|
|
outValue[1] = inv_ssts(x[1], C);
|
|
outValue[2] = inv_ssts(x[2], C);
|
|
|
|
return outValue;
|
|
}
|
|
|
|
float moncurve_f(float x, float gamma, float offs)
|
|
{
|
|
// Forward monitor curve
|
|
float y;
|
|
const float fs = ((gamma - 1.0) / offs) * pow(offs * gamma / ((gamma - 1.0) * (1.0 + offs)), gamma);
|
|
const float xb = offs / (gamma - 1.0);
|
|
if (x >= xb)
|
|
y = pow((x + offs) / (1.0 + offs), gamma);
|
|
else
|
|
y = x * fs;
|
|
return y;
|
|
}
|
|
|
|
float moncurve_r(float y, float gamma, float offs)
|
|
{
|
|
// Reverse monitor curve
|
|
float x;
|
|
const float yb = pow(offs * gamma / ((gamma - 1.0) * (1.0 + offs)), gamma);
|
|
const float rs = pow((gamma - 1.0) / offs, gamma - 1.0) * pow((1.0 + offs) / gamma, gamma);
|
|
if (y >= yb)
|
|
x = (1.0 + offs) * pow(y, 1.0 / gamma) - offs;
|
|
else
|
|
x = y * rs;
|
|
return x;
|
|
}
|
|
|
|
float3 moncurve_f_f3(float3 x, float gamma, float offs)
|
|
{
|
|
float3 y;
|
|
y[0] = moncurve_f(x[0], gamma, offs);
|
|
y[1] = moncurve_f(x[1], gamma, offs);
|
|
y[2] = moncurve_f(x[2], gamma, offs);
|
|
return y;
|
|
}
|
|
|
|
float3 moncurve_r_f3(float3 y, float gamma, float offs)
|
|
{
|
|
float3 x;
|
|
x[0] = moncurve_r(y[0], gamma, offs);
|
|
x[1] = moncurve_r(y[1], gamma, offs);
|
|
x[2] = moncurve_r(y[2], gamma, offs);
|
|
return x;
|
|
}
|
|
|
|
float bt1886_f(float V, float gamma, float Lw, float Lb)
|
|
{
|
|
// The reference EOTF specified in Rec. ITU-R BT.1886
|
|
// L = a(max[(V+b),0])^g
|
|
float a = pow(pow(Lw, 1. / gamma) - pow(Lb, 1. / gamma), gamma);
|
|
float b = pow(Lb, 1. / gamma) / (pow(Lw, 1. / gamma) - pow(Lb, 1. / gamma));
|
|
float L = a * pow(max(V + b, 0.), gamma);
|
|
return L;
|
|
}
|
|
|
|
float bt1886_r(float L, float gamma, float Lw, float Lb)
|
|
{
|
|
// The reference EOTF specified in Rec. ITU-R BT.1886
|
|
// L = a(max[(V+b),0])^g
|
|
float a = pow(pow(Lw, 1. / gamma) - pow(Lb, 1. / gamma), gamma);
|
|
float b = pow(Lb, 1. / gamma) / (pow(Lw, 1. / gamma) - pow(Lb, 1. / gamma));
|
|
float V = pow(max(L / a, 0.), 1. / gamma) - b;
|
|
return V;
|
|
}
|
|
|
|
float3 bt1886_f_f3(float3 V, float gamma, float Lw, float Lb)
|
|
{
|
|
float3 L;
|
|
L[0] = bt1886_f(V[0], gamma, Lw, Lb);
|
|
L[1] = bt1886_f(V[1], gamma, Lw, Lb);
|
|
L[2] = bt1886_f(V[2], gamma, Lw, Lb);
|
|
return L;
|
|
}
|
|
|
|
float3 bt1886_r_f3(float3 L, float gamma, float Lw, float Lb)
|
|
{
|
|
float3 V;
|
|
V[0] = bt1886_r(L[0], gamma, Lw, Lb);
|
|
V[1] = bt1886_r(L[1], gamma, Lw, Lb);
|
|
V[2] = bt1886_r(L[2], gamma, Lw, Lb);
|
|
return V;
|
|
}
|
|
|
|
|
|
// SMPTE Range vs Full Range scaling formulas
|
|
float smpteRange_to_fullRange(float inValue)
|
|
{
|
|
const float REFBLACK = (64. / 1023.);
|
|
const float REFWHITE = (940. / 1023.);
|
|
|
|
return ((inValue - REFBLACK) / (REFWHITE - REFBLACK));
|
|
}
|
|
|
|
float fullRange_to_smpteRange(float inValue)
|
|
{
|
|
const float REFBLACK = (64. / 1023.);
|
|
const float REFWHITE = (940. / 1023.);
|
|
|
|
return (inValue * (REFWHITE - REFBLACK) + REFBLACK);
|
|
}
|
|
|
|
float3 smpteRange_to_fullRange_f3(float3 rgbIn)
|
|
{
|
|
float3 rgbOut;
|
|
rgbOut[0] = smpteRange_to_fullRange(rgbIn[0]);
|
|
rgbOut[1] = smpteRange_to_fullRange(rgbIn[1]);
|
|
rgbOut[2] = smpteRange_to_fullRange(rgbIn[2]);
|
|
|
|
return rgbOut;
|
|
}
|
|
|
|
float3 fullRange_to_smpteRange_f3(float3 rgbIn)
|
|
{
|
|
float3 rgbOut;
|
|
rgbOut[0] = fullRange_to_smpteRange(rgbIn[0]);
|
|
rgbOut[1] = fullRange_to_smpteRange(rgbIn[1]);
|
|
rgbOut[2] = fullRange_to_smpteRange(rgbIn[2]);
|
|
|
|
return rgbOut;
|
|
}
|
|
|
|
// Base functions from SMPTE ST 2084-2014
|
|
|
|
// Constants from SMPTE ST 2084-2014
|
|
static const float pq_m1 = 0.1593017578125; // ( 2610.0 / 4096.0 ) / 4.0;
|
|
static const float pq_m2 = 78.84375; // ( 2523.0 / 4096.0 ) * 128.0;
|
|
static const float pq_c1 = 0.8359375; // 3424.0 / 4096.0 or pq_c3 - pq_c2 + 1.0;
|
|
static const float pq_c2 = 18.8515625; // ( 2413.0 / 4096.0 ) * 32.0;
|
|
static const float pq_c3 = 18.6875; // ( 2392.0 / 4096.0 ) * 32.0;
|
|
|
|
static const float pq_C = 10000.0;
|
|
|
|
// Converts from the non-linear perceptually quantized space to linear cd/m^2
|
|
// Note that this is in float, and assumes normalization from 0 - 1
|
|
// (0 - pq_C for linear) and does not handle the integer coding in the Annex
|
|
// sections of SMPTE ST 2084-2014
|
|
float ST2084_2_Y(float N)
|
|
{
|
|
// Note that this does NOT handle any of the signal range
|
|
// considerations from 2084 - this assumes full range (0 - 1)
|
|
float Np = pow(N, 1.0 / pq_m2);
|
|
float L = Np - pq_c1;
|
|
if (L < 0.0)
|
|
L = 0.0;
|
|
L = L / (pq_c2 - pq_c3 * Np);
|
|
L = pow(L, 1.0 / pq_m1);
|
|
return L * pq_C; // returns cd/m^2
|
|
}
|
|
|
|
// Converts from linear cd/m^2 to the non-linear perceptually quantized space
|
|
// Note that this is in float, and assumes normalization from 0 - 1
|
|
// (0 - pq_C for linear) and does not handle the integer coding in the Annex
|
|
// sections of SMPTE ST 2084-2014
|
|
float Y_2_ST2084(float C)
|
|
//pq_r
|
|
{
|
|
// Note that this does NOT handle any of the signal range
|
|
// considerations from 2084 - this returns full range (0 - 1)
|
|
float L = C / pq_C;
|
|
float Lm = pow(L, pq_m1);
|
|
float N = (pq_c1 + pq_c2 * Lm) / (1.0 + pq_c3 * Lm);
|
|
N = pow(N, pq_m2);
|
|
return N;
|
|
}
|
|
|
|
float3 Y_2_ST2084_f3(float3 inValue)
|
|
{
|
|
// converts from linear cd/m^2 to PQ code values
|
|
|
|
float3 outValue;
|
|
outValue[0] = Y_2_ST2084(inValue[0]);
|
|
outValue[1] = Y_2_ST2084(inValue[1]);
|
|
outValue[2] = Y_2_ST2084(inValue[2]);
|
|
|
|
return outValue;
|
|
}
|
|
|
|
float3 ST2084_2_Y_f3(float3 inValue)
|
|
{
|
|
// converts from PQ code values to linear cd/m^2
|
|
|
|
float3 outValue;
|
|
outValue[0] = ST2084_2_Y(inValue[0]);
|
|
outValue[1] = ST2084_2_Y(inValue[1]);
|
|
outValue[2] = ST2084_2_Y(inValue[2]);
|
|
|
|
return outValue;
|
|
}
|
|
|
|
// Conversion of PQ signal to HLG, as detailed in Section 7 of ITU-R BT.2390-0
|
|
float3 ST2084_2_HLG_1000nits_f3(float3 PQ)
|
|
{
|
|
// ST.2084 EOTF (non-linear PQ to display light)
|
|
float3 displayLinear = ST2084_2_Y_f3(PQ);
|
|
|
|
// HLG Inverse EOTF (i.e. HLG inverse OOTF followed by the HLG OETF)
|
|
// HLG Inverse OOTF (display linear to scene linear)
|
|
float Y_d = 0.2627 * displayLinear[0] + 0.6780 * displayLinear[1] + 0.0593 * displayLinear[2];
|
|
const float L_w = 1000.;
|
|
const float L_b = 0.;
|
|
const float alpha = (L_w - L_b);
|
|
const float beta = L_b;
|
|
const float gamma = 1.2;
|
|
|
|
float3 sceneLinear;
|
|
if (Y_d == 0.) {
|
|
/* This case is to protect against pow(0,-N)=Inf error. The ITU document
|
|
does not offer a recommendation for this corner case. There may be a
|
|
better way to handle this, but for now, this works.
|
|
*/
|
|
sceneLinear[0] = 0.;
|
|
sceneLinear[1] = 0.;
|
|
sceneLinear[2] = 0.;
|
|
}
|
|
else {
|
|
sceneLinear[0] = pow((Y_d - beta) / alpha, (1. - gamma) / gamma) * ((displayLinear[0] - beta) / alpha);
|
|
sceneLinear[1] = pow((Y_d - beta) / alpha, (1. - gamma) / gamma) * ((displayLinear[1] - beta) / alpha);
|
|
sceneLinear[2] = pow((Y_d - beta) / alpha, (1. - gamma) / gamma) * ((displayLinear[2] - beta) / alpha);
|
|
}
|
|
|
|
// HLG OETF (scene linear to non-linear signal value)
|
|
const float a = 0.17883277;
|
|
const float b = 0.28466892; // 1.-4.*a;
|
|
const float c = 0.55991073; // 0.5-a*log(4.*a);
|
|
|
|
float3 HLG;
|
|
if (sceneLinear[0] <= 1. / 12) {
|
|
HLG[0] = sqrt(3. * sceneLinear[0]);
|
|
}
|
|
else {
|
|
HLG[0] = a * log(12. * sceneLinear[0] - b) + c;
|
|
}
|
|
if (sceneLinear[1] <= 1. / 12) {
|
|
HLG[1] = sqrt(3. * sceneLinear[1]);
|
|
}
|
|
else {
|
|
HLG[1] = a * log(12. * sceneLinear[1] - b) + c;
|
|
}
|
|
if (sceneLinear[2] <= 1. / 12) {
|
|
HLG[2] = sqrt(3. * sceneLinear[2]);
|
|
}
|
|
else {
|
|
HLG[2] = a * log(12. * sceneLinear[2] - b) + c;
|
|
}
|
|
|
|
return HLG;
|
|
}
|
|
|
|
|
|
// Conversion of HLG to PQ signal, as detailed in Section 7 of ITU-R BT.2390-0
|
|
float3 HLG_2_ST2084_1000nits_f3(float3 HLG)
|
|
{
|
|
const float a = 0.17883277;
|
|
const float b = 0.28466892; // 1.-4.*a;
|
|
const float c = 0.55991073; // 0.5-a*log(4.*a);
|
|
|
|
const float L_w = 1000.;
|
|
const float L_b = 0.;
|
|
const float alpha = (L_w - L_b);
|
|
const float beta = L_b;
|
|
const float gamma = 1.2;
|
|
|
|
// HLG EOTF (non-linear signal value to display linear)
|
|
// HLG to scene-linear
|
|
float3 sceneLinear;
|
|
if (HLG[0] >= 0. && HLG[0] <= 0.5) {
|
|
sceneLinear[0] = pow(HLG[0], 2.) / 3.;
|
|
}
|
|
else {
|
|
sceneLinear[0] = (exp((HLG[0] - c) / a) + b) / 12.;
|
|
}
|
|
if (HLG[1] >= 0. && HLG[1] <= 0.5) {
|
|
sceneLinear[1] = pow(HLG[1], 2.) / 3.;
|
|
}
|
|
else {
|
|
sceneLinear[1] = (exp((HLG[1] - c) / a) + b) / 12.;
|
|
}
|
|
if (HLG[2] >= 0. && HLG[2] <= 0.5) {
|
|
sceneLinear[2] = pow(HLG[2], 2.) / 3.;
|
|
}
|
|
else {
|
|
sceneLinear[2] = (exp((HLG[2] - c) / a) + b) / 12.;
|
|
}
|
|
|
|
float Y_s = 0.2627 * sceneLinear[0] + 0.6780 * sceneLinear[1] + 0.0593 * sceneLinear[2];
|
|
|
|
// Scene-linear to display-linear
|
|
float3 displayLinear;
|
|
displayLinear[0] = alpha * pow(Y_s, gamma - 1.) * sceneLinear[0] + beta;
|
|
displayLinear[1] = alpha * pow(Y_s, gamma - 1.) * sceneLinear[1] + beta;
|
|
displayLinear[2] = alpha * pow(Y_s, gamma - 1.) * sceneLinear[2] + beta;
|
|
|
|
// ST.2084 Inverse EOTF
|
|
float3 PQ = Y_2_ST2084_f3(displayLinear);
|
|
|
|
return PQ;
|
|
}
|
|
|
|
// Desaturation contants
|
|
static const float RRT_SAT_FACTOR = 0.96;
|
|
static const float ONE_MINUS_RRT_SAT_FACTOR = 0.04;
|
|
static const float3x3 RRT_SAT_MAT =
|
|
{
|
|
// {ONE_MINUS_RRT_SAT_FACTOR * AP1_RGB2Y.x + RRT_SAT_FACTOR, ONE_MINUS_RRT_SAT_FACTOR * AP1_RGB2Y.y, ONE_MINUS_RRT_SAT_FACTOR * AP1_RGB2Y.z},
|
|
// {ONE_MINUS_RRT_SAT_FACTOR * AP1_RGB2Y.x, ONE_MINUS_RRT_SAT_FACTOR * AP1_RGB2Y.y + RRT_SAT_FACTOR, ONE_MINUS_RRT_SAT_FACTOR * AP1_RGB2Y.z},
|
|
// {ONE_MINUS_RRT_SAT_FACTOR * AP1_RGB2Y.x, ONE_MINUS_RRT_SAT_FACTOR * AP1_RGB2Y.y, ONE_MINUS_RRT_SAT_FACTOR * AP1_RGB2Y.z + RRT_SAT_FACTOR},
|
|
|
|
{0.970889148672, 0.02696327063, 0.0021475807},
|
|
{0.010889148672, 0.98696327063, 0.0021475807},
|
|
{0.010889148672, 0.02696327063, 0.9621475807}
|
|
};
|
|
|
|
static const float3x3 RRT_SAT_MAT_INV =
|
|
{
|
|
{ 1.03032386 , -0.0280867405 , -0.00223706313 },
|
|
{ -0.0113428626 , 1.01357996 , -0.00223706337 },
|
|
{ -0.0113428626 , -0.0280867405 , 1.03942955 }
|
|
};
|
|
|
|
|
|
// "Glow" module constants
|
|
static const float RRT_GLOW_GAIN = 0.05;
|
|
static const float RRT_GLOW_MID = 0.08;
|
|
|
|
// Red modifier constants
|
|
static const float RRT_RED_SCALE = 0.82;
|
|
static const float RRT_RED_PIVOT = 0.03;
|
|
static const float RRT_RED_HUE = 0.;
|
|
static const float RRT_RED_WIDTH = 135.;
|
|
|
|
|
|
float3 rrt_sweeteners(float3 inValue)
|
|
{
|
|
float3 aces = inValue;
|
|
|
|
// --- Glow module --- //
|
|
float saturation = rgb_2_saturation(aces);
|
|
float ycIn = rgb_2_yc(aces);
|
|
float s = sigmoid_shaper((saturation - 0.4) / 0.2);
|
|
float addedGlow = 1 + glow_fwd(ycIn, RRT_GLOW_GAIN * s, RRT_GLOW_MID);
|
|
aces *= addedGlow;
|
|
|
|
// --- Red modifier --- //
|
|
float hue = rgb_2_hue(aces);
|
|
float centeredHue = center_hue(hue, RRT_RED_HUE);
|
|
float hueWeight = cubic_basis_shaper(centeredHue, RRT_RED_WIDTH);
|
|
|
|
aces.r += hueWeight * saturation * (RRT_RED_PIVOT - aces.r) * (1. - RRT_RED_SCALE);
|
|
|
|
// --- ACES to RGB rendering space --- //
|
|
aces = clamp(aces, 0, HALF_POS_INF);
|
|
float3 rgbPre = mul(AP0_2_AP1_MAT, aces);
|
|
rgbPre = clamp(rgbPre, 0, HALF_MAX);
|
|
|
|
// --- Global desaturation --- //
|
|
//rgbPre = mul(RRT_SAT_MAT, rgbPre);
|
|
rgbPre = lerp(dot(rgbPre, AP1_RGB2Y).xxx, rgbPre, RRT_SAT_FACTOR);
|
|
return rgbPre;
|
|
}
|
|
|
|
float3 inv_rrt_sweeteners(float3 inValue)
|
|
{
|
|
float3 rgbPost = inValue;
|
|
|
|
// --- Global desaturation --- //
|
|
rgbPost = mul(RRT_SAT_MAT_INV, rgbPost);
|
|
|
|
rgbPost = clamp(rgbPost, 0.0.xxx, HALF_MAX.xxx);
|
|
|
|
// --- RGB rendering space to ACES --- //
|
|
float3 aces = mul(AP1_2_AP0_MAT, rgbPost);
|
|
|
|
aces = clamp(aces, 0.0.xxx, HALF_MAX.xxx);
|
|
|
|
// --- Red modifier --- //
|
|
float hue = rgb_2_hue(aces);
|
|
float centeredHue = center_hue(hue, RRT_RED_HUE);
|
|
float hueWeight = cubic_basis_shaper(centeredHue, RRT_RED_WIDTH);
|
|
|
|
float minChan;
|
|
if (centeredHue < 0) { // min_f3(aces) = aces[1] (i.e. magenta-red)
|
|
minChan = aces[1];
|
|
}
|
|
else { // min_f3(aces) = aces[2] (i.e. yellow-red)
|
|
minChan = aces[2];
|
|
}
|
|
|
|
float a = hueWeight * (1. - RRT_RED_SCALE) - 1.;
|
|
float b = aces[0] - hueWeight * (RRT_RED_PIVOT + minChan) * (1. - RRT_RED_SCALE);
|
|
float c = hueWeight * RRT_RED_PIVOT * minChan * (1. - RRT_RED_SCALE);
|
|
|
|
aces[0] = (-b - sqrt(b * b - 4. * a * c)) / (2. * a);
|
|
|
|
// --- Glow module --- //
|
|
float saturation = rgb_2_saturation(aces);
|
|
float ycOut = rgb_2_yc(aces);
|
|
float s = sigmoid_shaper((saturation - 0.4) / 0.2);
|
|
float reducedGlow = 1. + glow_inv(ycOut, RRT_GLOW_GAIN * s, RRT_GLOW_MID);
|
|
|
|
aces *= reducedGlow;
|
|
return aces;
|
|
}
|
|
|
|
float Y_2_linCV(float Y, float Ymax, float Ymin)
|
|
{
|
|
return (Y - Ymin) / (Ymax - Ymin);
|
|
}
|
|
|
|
float linCV_2_Y(float linCV, float Ymax, float Ymin)
|
|
{
|
|
return linCV * (Ymax - Ymin) + Ymin;
|
|
}
|
|
|
|
float3 Y_2_linCV_f3(float3 Y, float Ymax, float Ymin)
|
|
{
|
|
float3 linCV;
|
|
linCV[0] = Y_2_linCV(Y[0], Ymax, Ymin);
|
|
linCV[1] = Y_2_linCV(Y[1], Ymax, Ymin);
|
|
linCV[2] = Y_2_linCV(Y[2], Ymax, Ymin);
|
|
return linCV;
|
|
}
|
|
|
|
float3 linCV_2_Y_f3(float3 linCV, float Ymax, float Ymin)
|
|
{
|
|
float3 Y;
|
|
Y[0] = linCV_2_Y(linCV[0], Ymax, Ymin);
|
|
Y[1] = linCV_2_Y(linCV[1], Ymax, Ymin);
|
|
Y[2] = linCV_2_Y(linCV[2], Ymax, Ymin);
|
|
return Y;
|
|
}
|
|
|
|
// Gamma compensation factor
|
|
static const float DIM_SURROUND_GAMMA = 0.9811;
|
|
|
|
float3 darkSurround_to_dimSurround(float3 linearCV)
|
|
{
|
|
float3 XYZ = mul(AP1_2_XYZ_MAT, linearCV);
|
|
|
|
float3 xyY = XYZ_2_xyY(XYZ);
|
|
xyY[2] = clamp(xyY[2], 0., HALF_POS_INF);
|
|
xyY[2] = pow(xyY[2], DIM_SURROUND_GAMMA);
|
|
XYZ = xyY_2_XYZ(xyY);
|
|
|
|
return mul(XYZ_2_AP1_MAT, XYZ);
|
|
}
|
|
|
|
float3 dimSurround_to_darkSurround(float3 linearCV)
|
|
{
|
|
float3 XYZ = mul(AP1_2_XYZ_MAT, linearCV);
|
|
|
|
float3 xyY = XYZ_2_xyY(XYZ);
|
|
xyY[2] = clamp(xyY[2], 0., HALF_POS_INF);
|
|
xyY[2] = pow(xyY[2], 1. / DIM_SURROUND_GAMMA);
|
|
XYZ = xyY_2_XYZ(xyY);
|
|
|
|
return mul(XYZ_2_AP1_MAT, XYZ);
|
|
}
|
|
|
|
float3 dark_to_dim(float3 XYZ)
|
|
{
|
|
float3 xyY = XYZ_2_xyY(XYZ);
|
|
xyY[2] = clamp(xyY[2], 0., HALF_POS_INF);
|
|
xyY[2] = pow(xyY[2], DIM_SURROUND_GAMMA);
|
|
return xyY_2_XYZ(xyY);
|
|
}
|
|
|
|
float3 dim_to_dark(float3 XYZ)
|
|
{
|
|
float3 xyY = XYZ_2_xyY(XYZ);
|
|
xyY[2] = clamp(xyY[2], 0., HALF_POS_INF);
|
|
xyY[2] = pow(xyY[2], 1. / DIM_SURROUND_GAMMA);
|
|
return xyY_2_XYZ(xyY);
|
|
}
|
|
|
|
float3 outputTransform
|
|
(
|
|
float3 inValue,
|
|
TsParams PARAMS,
|
|
FColorSpace DISPLAY_PRI,
|
|
FColorSpace LIMITING_PRI,
|
|
int EOTF,
|
|
int SURROUND,
|
|
bool STRETCH_BLACK = true,
|
|
bool D60_SIM = false,
|
|
bool LEGAL_RANGE = false
|
|
)
|
|
{
|
|
float3x3 XYZ_2_DISPLAY_PRI_MAT = DISPLAY_PRI.XYZtoRGB;
|
|
|
|
/*
|
|
NOTE: This is a bit of a hack - probably a more direct way to do this.
|
|
TODO: Fix in future version
|
|
*/
|
|
|
|
float Y_MIN = PARAMS.Min.y;
|
|
float Y_MAX = PARAMS.Max.y;
|
|
|
|
// RRT sweeteners
|
|
float3 rgbPre = rrt_sweeteners(inValue);
|
|
|
|
// Apply the tonescale independently in rendering-space RGB
|
|
float3 rgbPost = ssts_f3(rgbPre, PARAMS);
|
|
|
|
// At this point data encoded AP1, scaled absolute luminance (cd/m^2)
|
|
|
|
/* Scale absolute luminance to linear code value */
|
|
float3 linearCV = Y_2_linCV_f3(rgbPost, Y_MAX, Y_MIN);
|
|
|
|
// Rendering primaries to XYZ
|
|
float3 XYZ = mul(AP1_2_XYZ_MAT, linearCV);
|
|
|
|
// Apply gamma adjustment to compensate for dim surround
|
|
/*
|
|
NOTE: This is more or less a placeholder block and is largely inactive
|
|
in its current form. This section currently only applies for SDR, and
|
|
even then, only in very specific cases.
|
|
In the future it is fully intended for this module to be updated to
|
|
support surround compensation regardless of luminance dynamic range. */
|
|
/*
|
|
TOD0: Come up with new surround compensation algorithm, applicable
|
|
across all dynamic ranges and supporting dark/dim/normal surround.
|
|
*/
|
|
if (SURROUND == 0) { // Dark surround
|
|
/*
|
|
Current tone scale is designed for dark surround environment so no
|
|
adjustment is necessary.
|
|
*/
|
|
}
|
|
else if (SURROUND == 1) { // Dim surround
|
|
// INACTIVE for HDR and crudely implemented for SDR (see comment below)
|
|
if ((EOTF == 1) || (EOTF == 2) || (EOTF == 3)) {
|
|
/*
|
|
This uses a crude logical assumption that if the EOTF is BT.1886,
|
|
sRGB, or gamma 2.6 that the data is SDR and so the SDR gamma
|
|
compensation factor from v1.0 will apply.
|
|
*/
|
|
XYZ = dark_to_dim(XYZ); /*
|
|
This uses a local dark_to_dim function that is designed to take in
|
|
XYZ and return XYZ rather than AP1 as is currently in the functions
|
|
in 'ACESlib.ODT_Common.ctl' */
|
|
}
|
|
}
|
|
else if (SURROUND == 2) { // Normal surround
|
|
// INACTIVE - this does NOTHING
|
|
}
|
|
|
|
// Gamut limit to limiting primaries
|
|
// NOTE: Would be nice to just say
|
|
// if (LIMITING_PRI != DISPLAY_PRI)
|
|
// but you can't because Chromaticities do not work with bool comparison operator
|
|
// For now, limit no matter what.
|
|
XYZ = limit_to_primaries(XYZ, LIMITING_PRI);
|
|
|
|
// Apply CAT from ACES white point to assumed observer adapted white point
|
|
// TODO: Needs to expand from just supporting D60 sim to allow for any
|
|
// observer adapted white point.
|
|
if (D60_SIM == false) {
|
|
if (DISPLAY_PRI.White != WhiteStandard_D60) {
|
|
XYZ = mul(D60_2_D65_CAT, XYZ);
|
|
}
|
|
}
|
|
|
|
// CIE XYZ to display encoding primaries
|
|
linearCV = mul(XYZ_2_DISPLAY_PRI_MAT, XYZ);
|
|
|
|
// Scale to avoid clipping when device calibration is different from D60.
|
|
// To simulate D60, unequal code values are sent to the display.
|
|
// TODO: Needs to expand from just supporting D60 sim to allow for any
|
|
// observer adapted white point.
|
|
if (D60_SIM == true) {
|
|
/* TODO: The scale requires calling itself. Scale is same no matter the luminance.
|
|
Currently precalculated for D65, DCI. If DCI, roll_white_fwd is used also.
|
|
This needs a more complex algorithm to handle all cases.
|
|
*/
|
|
float SCALE = 1.0;
|
|
if (DISPLAY_PRI.White == WhiteStandard_D65) { // D65
|
|
SCALE = 0.96362;
|
|
}
|
|
else if (DISPLAY_PRI.White == WhiteStandard_DCI) { // DCI
|
|
linearCV[0] = roll_white_fwd(linearCV[0], 0.918, 0.5);
|
|
linearCV[1] = roll_white_fwd(linearCV[1], 0.918, 0.5);
|
|
linearCV[2] = roll_white_fwd(linearCV[2], 0.918, 0.5);
|
|
SCALE = 0.96;
|
|
}
|
|
linearCV *= SCALE;
|
|
}
|
|
|
|
|
|
// Clip values < 0 (i.e. projecting outside the display primaries)
|
|
// NOTE: P3 red and values close to it fall outside of Rec.2020 green-red
|
|
// boundary
|
|
linearCV = clamp(linearCV, 0.0.xxx, HALF_POS_INF.xxx);
|
|
|
|
// EOTF
|
|
// 0: ST-2084 (PQ)
|
|
// 1: BT.1886 (Rec.709/2020 settings)
|
|
// 2: sRGB (mon_curve w/ presets)
|
|
// moncurve_r with gamma of 2.4 and offset of 0.055 matches the EOTF found in IEC 61966-2-1:1999 (sRGB)
|
|
// 3: gamma 2.6
|
|
// 4: linear (no EOTF)
|
|
// 5: HLG
|
|
float3 outputCV;
|
|
if (EOTF == 0) { // ST-2084 (PQ)
|
|
// NOTE: This is a kludgy way of ensuring a PQ code value of 0. Ideally,
|
|
// luminance would map directly to code value, but colorists don't like
|
|
// that. Might just need the tonescale to go darker so that darkest
|
|
// values through the tone scale quantize to code value of 0.
|
|
if (STRETCH_BLACK == true) {
|
|
outputCV = Y_2_ST2084_f3(clamp(linCV_2_Y_f3(linearCV, Y_MAX, 0.0), 0.0, HALF_POS_INF));
|
|
}
|
|
else {
|
|
outputCV = Y_2_ST2084_f3(linCV_2_Y_f3(linearCV, Y_MAX, Y_MIN));
|
|
}
|
|
}
|
|
else if (EOTF == 1) { // BT.1886 (Rec.709/2020 settings)
|
|
outputCV = bt1886_r_f3(linearCV, 2.4, 1.0, 0.0);
|
|
}
|
|
else if (EOTF == 2) { // sRGB (mon_curve w/ presets)
|
|
outputCV = moncurve_r_f3(linearCV, 2.4, 0.055);
|
|
}
|
|
else if (EOTF == 3) { // gamma 2.6
|
|
outputCV = pow(linearCV, (1. / 2.6).xxx);
|
|
}
|
|
else if (EOTF == 4) { // linear
|
|
///+ have linear behave the same way at ST-2084 regarding STRETCH_BLACK option
|
|
if (STRETCH_BLACK == true) {
|
|
outputCV = clamp(linCV_2_Y_f3(linearCV, Y_MAX, 0.0), 0.0, HALF_POS_INF);
|
|
}
|
|
else {
|
|
outputCV = linCV_2_Y_f3(linearCV, Y_MAX, Y_MIN);
|
|
}
|
|
}
|
|
else if (EOTF == 5) { // HLG
|
|
// NOTE: HLG just maps ST.2084 output to HLG encoding.
|
|
// TODO: Restructure if/else tree to minimize this redundancy.
|
|
if (STRETCH_BLACK == true) {
|
|
outputCV = Y_2_ST2084_f3(clamp(linCV_2_Y_f3(linearCV, Y_MAX, 0.0), 0.0, HALF_POS_INF));
|
|
}
|
|
else {
|
|
outputCV = Y_2_ST2084_f3(linCV_2_Y_f3(linearCV, Y_MAX, Y_MIN));
|
|
}
|
|
outputCV = ST2084_2_HLG_1000nits_f3(outputCV);
|
|
}
|
|
|
|
if (LEGAL_RANGE == true) {
|
|
outputCV = fullRange_to_smpteRange_f3(outputCV);
|
|
}
|
|
|
|
return outputCV;
|
|
}
|
|
|
|
float3 outputTransform
|
|
(
|
|
float3 inValue,
|
|
float Y_MIN,
|
|
float Y_MID,
|
|
float Y_MAX,
|
|
FColorSpace DISPLAY_PRI,
|
|
FColorSpace LIMITING_PRI,
|
|
int EOTF,
|
|
int SURROUND,
|
|
bool STRETCH_BLACK = true,
|
|
bool D60_SIM = false,
|
|
bool LEGAL_RANGE = false
|
|
)
|
|
{
|
|
/*
|
|
NOTE: This is a bit of a hack - probably a more direct way to do this.
|
|
TODO: Fix in future version
|
|
*/
|
|
TsParams PARAMS_DEFAULT = init_TsParams(Y_MIN, Y_MAX);
|
|
float expShift = log2(inv_ssts(Y_MID, PARAMS_DEFAULT)) - log2(0.18);
|
|
TsParams PARAMS = init_TsParams(Y_MIN, Y_MAX, expShift);
|
|
|
|
return outputTransform(inValue, PARAMS, DISPLAY_PRI, LIMITING_PRI, EOTF, SURROUND, STRETCH_BLACK, D60_SIM, LEGAL_RANGE);
|
|
}
|
|
|
|
float3 invOutputTransform
|
|
(
|
|
float3 inValue,
|
|
float Y_MIN,
|
|
float Y_MID,
|
|
float Y_MAX,
|
|
FColorSpace DISPLAY_PRI,
|
|
FColorSpace LIMITING_PRI,
|
|
int EOTF,
|
|
int SURROUND,
|
|
bool STRETCH_BLACK = true,
|
|
bool D60_SIM = false,
|
|
bool LEGAL_RANGE = false
|
|
)
|
|
{
|
|
float3x3 DISPLAY_PRI_2_XYZ_MAT = DISPLAY_PRI.RGBtoXYZ;
|
|
|
|
/*
|
|
NOTE: This is a bit of a hack - probably a more direct way to do this.
|
|
TODO: Update in accordance with forward algorithm.
|
|
*/
|
|
TsParams PARAMS_DEFAULT = init_TsParams(Y_MIN, Y_MAX);
|
|
float expShift = log2(inv_ssts(Y_MID, PARAMS_DEFAULT)) - log2(0.18);
|
|
TsParams PARAMS = init_TsParams(Y_MIN, Y_MAX, expShift);
|
|
|
|
float3 outputCV = inValue;
|
|
|
|
if (LEGAL_RANGE == true) {
|
|
outputCV = smpteRange_to_fullRange_f3(outputCV);
|
|
}
|
|
|
|
// Inverse EOTF
|
|
// 0: ST-2084 (PQ)
|
|
// 1: BT.1886 (Rec.709/2020 settings)
|
|
// 2: sRGB (mon_curve w/ presets)
|
|
// moncurve_r with gamma of 2.4 and offset of 0.055 matches the EOTF found in IEC 61966-2-1:1999 (sRGB)
|
|
// 3: gamma 2.6
|
|
// 4: linear (no EOTF)
|
|
// 5: HLG
|
|
float3 linearCV;
|
|
if (EOTF == 0) { // ST-2084 (PQ)
|
|
if (STRETCH_BLACK == true) {
|
|
linearCV = Y_2_linCV_f3(ST2084_2_Y_f3(outputCV), Y_MAX, 0.);
|
|
}
|
|
else {
|
|
linearCV = Y_2_linCV_f3(ST2084_2_Y_f3(outputCV), Y_MAX, Y_MIN);
|
|
}
|
|
}
|
|
else if (EOTF == 1) { // BT.1886 (Rec.709/2020 settings)
|
|
linearCV = bt1886_f_f3(outputCV, 2.4, 1.0, 0.0);
|
|
}
|
|
else if (EOTF == 2) { // sRGB (mon_curve w/ presets)
|
|
linearCV = moncurve_f_f3(outputCV, 2.4, 0.055);
|
|
}
|
|
else if (EOTF == 3) { // gamma 2.6
|
|
linearCV = pow(outputCV, 2.6);
|
|
}
|
|
else if (EOTF == 4) { // linear
|
|
linearCV = Y_2_linCV_f3(outputCV, Y_MAX, Y_MIN);
|
|
}
|
|
else if (EOTF == 5) { // HLG
|
|
outputCV = HLG_2_ST2084_1000nits_f3(outputCV);
|
|
if (STRETCH_BLACK == true) {
|
|
linearCV = Y_2_linCV_f3(ST2084_2_Y_f3(outputCV), Y_MAX, 0.);
|
|
}
|
|
else {
|
|
linearCV = Y_2_linCV_f3(ST2084_2_Y_f3(outputCV), Y_MAX, Y_MIN);
|
|
}
|
|
}
|
|
|
|
// Un-scale
|
|
if (D60_SIM == true) {
|
|
/* TODO: The scale requires calling itself. Need an algorithm for this.
|
|
Scale is same no matter the luminance.
|
|
Currently using precalculated values for D65, DCI.
|
|
If DCI, roll_white_fwd is used also.
|
|
*/
|
|
float SCALE = 1.0;
|
|
if (DISPLAY_PRI.White == WhiteStandard_D65) { // D65
|
|
SCALE = 0.96362;
|
|
linearCV /= SCALE;
|
|
}
|
|
else if (DISPLAY_PRI.White == WhiteStandard_DCI) { // DCI
|
|
SCALE = 0.96;
|
|
linearCV[0] = roll_white_rev(linearCV[0] / SCALE, 0.918, 0.5);
|
|
linearCV[1] = roll_white_rev(linearCV[1] / SCALE, 0.918, 0.5);
|
|
linearCV[2] = roll_white_rev(linearCV[2] / SCALE, 0.918, 0.5);
|
|
}
|
|
|
|
}
|
|
|
|
// Encoding primaries to CIE XYZ
|
|
float3 XYZ = mul(DISPLAY_PRI_2_XYZ_MAT, linearCV);
|
|
|
|
// Undo CAT from assumed observer adapted white point to ACES white point
|
|
if (D60_SIM == false) {
|
|
if (DISPLAY_PRI.White != WhiteStandard_D60) {
|
|
XYZ = mul(D65_2_D60_CAT, XYZ);
|
|
}
|
|
}
|
|
|
|
// Apply gamma adjustment to compensate for dim surround
|
|
/*
|
|
NOTE: This is more or less a placeholder block and is largely inactive
|
|
in its current form. This section currently only applies for SDR, and
|
|
even then, only in very specific cases.
|
|
In the future it is fully intended for this module to be updated to
|
|
support surround compensation regardless of luminance dynamic range. */
|
|
/*
|
|
TOD0: Come up with new surround compensation algorithm, applicable
|
|
across all dynamic ranges and supporting dark/dim/normal surround.
|
|
*/
|
|
if (SURROUND == 0) { // Dark surround
|
|
/*
|
|
Current tone scale is designed for dark surround environment so no
|
|
adjustment is necessary.
|
|
*/
|
|
}
|
|
else if (SURROUND == 1) { // Dim surround
|
|
// INACTIVE for HDR and crudely implemented for SDR (see comment below)
|
|
if ((EOTF == 1) || (EOTF == 2) || (EOTF == 3)) {
|
|
/*
|
|
This uses a crude logical assumption that if the EOTF is BT.1886,
|
|
sRGB, or gamma 2.6 that the data is SDR and so the SDR gamma
|
|
compensation factor from v1.0 will apply.
|
|
*/
|
|
XYZ = dim_to_dark(XYZ); /*
|
|
This uses a local dim_to_dark function that is designed to take in
|
|
XYZ and return XYZ rather than AP1 as is currently in the functions
|
|
in 'ACESlib.ODT_Common.ctl' */
|
|
}
|
|
}
|
|
else if (SURROUND == 2) { // Normal surround
|
|
// INACTIVE - this does NOTHING
|
|
}
|
|
|
|
// XYZ to rendering primaries
|
|
linearCV = mul(XYZ_2_AP1_MAT, XYZ);
|
|
|
|
float3 rgbPost = linCV_2_Y_f3(linearCV, Y_MAX, Y_MIN);
|
|
|
|
// Apply the inverse tonescale independently in rendering-space RGB
|
|
float3 rgbPre = inv_ssts_f3(rgbPost, PARAMS);
|
|
|
|
// RRT sweeteners
|
|
float3 aces = inv_rrt_sweeteners(rgbPre);
|
|
|
|
return aces;
|
|
}
|
|
|
|
/* --- Gamut Compress Parameters --- */
|
|
// Distance from achromatic which will be compressed to the gamut boundary
|
|
// Values calculated to encompass the encoding gamuts of common digital cinema cameras
|
|
static const float LIM_CYAN = 1.147;
|
|
static const float LIM_MAGENTA = 1.264;
|
|
static const float LIM_YELLOW = 1.312;
|
|
|
|
// Percentage of the core gamut to protect
|
|
// Values calculated to protect all the colors of the ColorChecker Classic 24 as given by
|
|
// ISO 17321-1 and Ohta (1997)
|
|
static const float THR_CYAN = 0.815;
|
|
static const float THR_MAGENTA = 0.803;
|
|
static const float THR_YELLOW = 0.880;
|
|
|
|
// Aggressiveness of the compression curve
|
|
static const float PWR = 1.2;
|
|
|
|
|
|
|
|
// Calculate compressed distance
|
|
float compress(float dist, float lim, float thr, float pwr, bool invert)
|
|
{
|
|
float comprDist;
|
|
float scl;
|
|
float nd;
|
|
float p;
|
|
|
|
if (dist < thr) {
|
|
comprDist = dist; // No compression below threshold
|
|
}
|
|
else {
|
|
// Calculate scale factor for y = 1 intersect
|
|
scl = (lim - thr) / pow(pow((1.0 - thr) / (lim - thr), -pwr) - 1.0, 1.0 / pwr);
|
|
|
|
// Normalize distance outside threshold by scale factor
|
|
nd = (dist - thr) / scl;
|
|
p = pow(nd, pwr);
|
|
|
|
if (!invert) {
|
|
comprDist = thr + scl * nd / (pow(1.0 + p, 1.0 / pwr)); // Compress
|
|
}
|
|
else {
|
|
if (dist > (thr + scl)) {
|
|
comprDist = dist; // Avoid singularity
|
|
}
|
|
else {
|
|
comprDist = thr + scl * pow(-(p / (p - 1.0)), 1.0 / pwr); // Uncompress
|
|
}
|
|
}
|
|
}
|
|
|
|
return comprDist;
|
|
}
|
|
|
|
float max_f3(float3 a)
|
|
{
|
|
return max(a[0], max(a[1], a[2]));
|
|
}
|
|
|
|
float3 compressColor
|
|
(
|
|
const float3 ACES,
|
|
const bool invert = false
|
|
)
|
|
{
|
|
// Convert to ACEScg
|
|
float3 linAP1 = mul(AP0_2_AP1_MAT, ACES);
|
|
|
|
// Achromatic axis
|
|
float ach = max_f3(linAP1);
|
|
|
|
// Distance from the achromatic axis for each color component aka inverse RGB ratios
|
|
float3 comprLinAP1;
|
|
///+ added any(ACES < 0): don't try to compress colors outside the AP0 color space: clamping will happen anyway during RRT/ODT,
|
|
// and this allows to keep proper results when testing with a synthetic chart
|
|
if (ach < 1e-10f || any(ACES < 0)) {
|
|
///-
|
|
comprLinAP1 = linAP1;
|
|
}
|
|
else {
|
|
float3 dist;
|
|
dist[0] = (ach - linAP1[0]) / abs(ach);
|
|
dist[1] = (ach - linAP1[1]) / abs(ach);
|
|
dist[2] = (ach - linAP1[2]) / abs(ach);
|
|
|
|
float3 comprDist = {
|
|
compress(dist[0], LIM_CYAN, THR_CYAN, PWR, invert),
|
|
compress(dist[1], LIM_MAGENTA, THR_MAGENTA, PWR, invert),
|
|
compress(dist[2], LIM_YELLOW, THR_YELLOW, PWR, invert)
|
|
};
|
|
|
|
// Recalculate RGB from compressed distance and achromatic
|
|
comprLinAP1 = float3(
|
|
ach - comprDist[0] * abs(ach),
|
|
ach - comprDist[1] * abs(ach),
|
|
ach - comprDist[2] * abs(ach)
|
|
);
|
|
}
|
|
|
|
// Convert back to ACES2065-1
|
|
return mul(AP1_2_AP0_MAT, comprLinAP1);
|
|
}
|
|
|