Skip to content

Commit

Permalink
+add Base implementation of function BgrToLab (part 4: refactoring).
Browse files Browse the repository at this point in the history
  • Loading branch information
ermig1979 committed Feb 5, 2025
1 parent 83e524d commit 2073a5a
Show file tree
Hide file tree
Showing 2 changed files with 101 additions and 242 deletions.
237 changes: 47 additions & 190 deletions src/Simd/SimdBaseBgrToLab.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -26,217 +26,74 @@
#include "Simd/SimdMath.h"
#include "Simd/SimdLog.h"

#define SIMD_BGR_TO_LAB_OPENCV_COMPATIBILITY

namespace Simd
{
namespace Base
{
#define SIMD_BGR_TO_LAB_OPENCV_COMPATIBILITY
uint16_t LabGammaTab[LabGammaTabSize];
uint16_t LabCbrtTab[LabCbrtTabSize];
uint32_t LabCoeffsTab[9];

namespace LabV1
SIMD_INLINE float ApplyGamma(float x)
{
inline double FromRaw(uint64_t raw)
{
return *((double*)&raw);
}

SIMD_INLINE float mullAdd(float a, float b, float c)
{
return a * b + c;
}

#define CV_DESCALE(x,n) (((x) + (1 << ((n)-1))) >> (n))
const int xyz_shift = 12;

//const double D65[3] = { 0.950456, 1., 1.088754 };
const double D65[3] = { FromRaw(0x3fee6a22b3892ee8),
1.0 ,
FromRaw(0x3ff16b8950763a19) };

//static const double sRGB2XYZ_D65[] =
//{
// 0.412453, 0.357580, 0.180423,
// 0.212671, 0.715160, 0.072169,
// 0.019334, 0.119193, 0.950227
//};

const double sRGB2XYZ_D65[] =
{
FromRaw(0x3fda65a14488c60d),
FromRaw(0x3fd6e297396d0918),
FromRaw(0x3fc71819d2391d58),
FromRaw(0x3fcb38cda6e75ff6),
FromRaw(0x3fe6e297396d0918),
FromRaw(0x3fb279aae6c8f755),
FromRaw(0x3f93cc4ac6cdaf4b),
FromRaw(0x3fbe836eb4e98138),
FromRaw(0x3fee68427418d691)
};

enum { LAB_CBRT_TAB_SIZE = 1024, GAMMA_TAB_SIZE = 1024 };

static const float GammaTabScale((int)GAMMA_TAB_SIZE);

static uint16_t sRGBGammaTab_b[256];

#undef lab_shift
#define lab_shift xyz_shift
#define gamma_shift 3
#define lab_shift2 (lab_shift + gamma_shift)
#define LAB_CBRT_TAB_SIZE_B (256*3/2*(1<<gamma_shift))
static uint16_t LabCbrtTab_b[LAB_CBRT_TAB_SIZE_B];

//all constants should be presented through integers to keep bit-exactness
static const double gammaThreshold = double(809) / double(20000); // 0.04045
static const double gammaInvThreshold = double(7827) / double(2500000); // 0.0031308
static const double gammaLowScale = double(323) / double(25); // 12.92
static const double gammaPower = double(12) / double(5); // 2.4
static const double gammaXshift = double(11) / double(200); // 0.055
return x <= 0.04045f ? x * (1.f / 12.92f) : (float)std::pow((double)(x + 0.055) * (1. / 1.055), 2.4);
}

static const float lthresh = float(216) / float(24389); // 0.008856f = (6/29)^3
static const float lscale = float(841) / float(108); // 7.787f = (29/3)^3/(29*4)
static const float lbias = float(16) / float(116);
SIMD_INLINE float ApplyCbrt(float x)
{
static const float threshold = float(216) / float(24389);
static const float scale = float(841) / float(108);
static const float bias = float(16) / float(116);
return x < threshold ? (x * scale + bias) : cbrt(x);
}

static inline float applyGamma(float x)
static bool LabTabsInited = false;
void LabInitTabs()
{
if (LabTabsInited)
return;
const float intScale(255 * (1 << LAB_GAMMA_SHIFT));
for (int i = 0; i < LabGammaTabSize; i++)
{
//return x <= 0.04045f ? x*(1.f/12.92f) : (float)std::pow((double)(x + 0.055)*(1./1.055), 2.4);

double xd = x;
return (xd <= gammaThreshold ?
xd / gammaLowScale :
pow((xd + gammaXshift) / (1.0 + gammaXshift), gammaPower));
float x = float(i) / 255.0f;
LabGammaTab[i] = (uint16_t)(Round(intScale * ApplyGamma(x)));
}

static void initLabTabs()
const float tabScale(1.0 / (255.0 * (1 << LAB_GAMMA_SHIFT)));
const float shift2(1 << LAB_SHIFT2);
for (int i = 0; i < LabCbrtTabSize; i++)
{
static bool initialized = false;
if (!initialized)
{
float scale = 1.0f / float(GammaTabScale);
static const float intScale(255 * (1 << gamma_shift));
for (int i = 0; i < 256; i++)
{
float x = float(i) / 255.0f;
sRGBGammaTab_b[i] = (uint16_t)(Round(intScale * applyGamma(x)));
}

static const float cbTabScale(1.0 / (255.0 * (1 << gamma_shift)));
static const float lshift2(1 << lab_shift2);
for (int i = 0; i < LAB_CBRT_TAB_SIZE_B; i++)
{
float x = cbTabScale * float(i);
LabCbrtTab_b[i] = (uint16_t)(Round(lshift2 * (x < lthresh ? mullAdd(x, lscale, lbias) : cbrt(x))));
}
#if defined(SIMD_BGR_TO_LAB_OPENCV_COMPATIBILITY)
if(LabCbrtTab_b[324] == 17746)
LabCbrtTab_b[324] = 17745;
if (LabCbrtTab_b[49] == 9455)
LabCbrtTab_b[49] = 9454;
#endif
initialized = true;
}
float x = tabScale * float(i);
LabCbrtTab[i] = (uint16_t)Round(shift2 * ApplyCbrt(x));
}

struct RGB2Lab_b
{
typedef uint8_t channel_type;

RGB2Lab_b(int _srccn, int blueIdx)
: srccn(_srccn)
{
initLabTabs();

double whitePt[3];
for (int i = 0; i < 3; i++)
whitePt[i] = D65[i];

const double lshift(1 << lab_shift);
for (int i = 0; i < 3; i++)
{
double c[3];
for (int j = 0; j < 3; j++)
c[j] = sRGB2XYZ_D65[i * 3 + j];
coeffs[i * 3 + (blueIdx ^ 2)] = Round(lshift * c[0] / whitePt[i]);
coeffs[i * 3 + 1] = Round(lshift * c[1] / whitePt[i]);
coeffs[i * 3 + blueIdx] = Round(lshift * c[2] / whitePt[i]);

assert(coeffs[i * 3] >= 0 && coeffs[i * 3 + 1] >= 0 && coeffs[i * 3 + 2] >= 0 &&
coeffs[i * 3] + coeffs[i * 3 + 1] + coeffs[i * 3 + 2] < 2 * (1 << lab_shift));
}
}

void operator()(const uint8_t* src, uint8_t* dst, int n) const
{
const int Lscale = (116 * 255 + 50) / 100;
const int Lshift = -((16 * 255 * (1 << lab_shift2) + 50) / 100);
const uint16_t* tab = sRGBGammaTab_b;
int i, scn = srccn;
int C0 = coeffs[0], C1 = coeffs[1], C2 = coeffs[2],
C3 = coeffs[3], C4 = coeffs[4], C5 = coeffs[5],
C6 = coeffs[6], C7 = coeffs[7], C8 = coeffs[8];

for (int i = 0; i < n; i++, src += scn, dst += 3)
{
int R = tab[src[0]], G = tab[src[1]], B = tab[src[2]];
int fX = LabCbrtTab_b[CV_DESCALE(R * C0 + G * C1 + B * C2, lab_shift)];
int fY = LabCbrtTab_b[CV_DESCALE(R * C3 + G * C4 + B * C5, lab_shift)];
int fZ = LabCbrtTab_b[CV_DESCALE(R * C6 + G * C7 + B * C8, lab_shift)];

int L = CV_DESCALE(Lscale * fY + Lshift, lab_shift2);
int a = CV_DESCALE(500 * (fX - fY) + 128 * (1 << lab_shift2), lab_shift2);
int b = CV_DESCALE(200 * (fY - fZ) + 128 * (1 << lab_shift2), lab_shift2);

#if defined(_WIN32) && 0
if (L == 45 && a == 165)
{
std::cout << i << " old Lab = {" << L << ", " << a << ", " << b << "}";

int R = tab[src[0]], G = tab[src[1]], B = tab[src[2]];
int fX = LabCbrtTab_b[CV_DESCALE(R * C0 + G * C1 + B * C2, lab_shift)];
int fY = LabCbrtTab_b[CV_DESCALE(R * C3 + G * C4 + B * C5, lab_shift)] - 1;
int fZ = LabCbrtTab_b[CV_DESCALE(R * C6 + G * C7 + B * C8, lab_shift)];

int L = CV_DESCALE(Lscale * fY + Lshift, lab_shift2);
int a = CV_DESCALE(500 * (fX - fY) + 128 * (1 << lab_shift2), lab_shift2);
int b = CV_DESCALE(200 * (fY - fZ) + 128 * (1 << lab_shift2), lab_shift2);

std::cout << " new Lab = {" << L << ", " << a << ", " << b << "}";
std::cout << " index " << CV_DESCALE(R * C3 + G * C4 + B * C5, lab_shift);
std::cout << " value " << LabCbrtTab_b[CV_DESCALE(R * C3 + G * C4 + B * C5, lab_shift)];
std::cout << std::endl;
exit(0);
}
#if defined(SIMD_BGR_TO_LAB_OPENCV_COMPATIBILITY)
if (LabCbrtTab[324] == 17746)
LabCbrtTab[324] = 17745;
if (LabCbrtTab[49] == 9455)
LabCbrtTab[49] = 9454;
#endif
dst[0] = Base::RestrictRange(L);
dst[1] = Base::RestrictRange(a);
dst[2] = Base::RestrictRange(b);
}
}

int srccn;
int coeffs[9];
};
const double D65[3] = { 0.950456, 1., 1.088754 };
const double sRGB2XYZ_D65[9] = { 0.412453, 0.357580, 0.180423, 0.212671, 0.715160, 0.072169, 0.019334, 0.119193, 0.950227 };
const double shift(1 << LAB_SHIFT);
for (int i = 0; i < 3; i++)
for (int j = 0; j < 3; j++)
LabCoeffsTab[i * 3 + j] = Round(shift * sRGB2XYZ_D65[i * 3 + j] / D65[i]);
LabTabsInited = true;
}

//--------------------------------------------------------------------------------------------------

void BgrToLab(const uint8_t* bgr, size_t bgrStride, size_t width, size_t height, uint8_t* lab, size_t labStride)
{
#if 1
LabV1::RGB2Lab_b _lab(3, 0);
LabInitTabs();
for (size_t row = 0; row < height; ++row)
{
const uint8_t* pBgr = bgr + row * bgrStride;
const uint8_t* pBgr = bgr + row * bgrStride, * pEnd = pBgr + width * 3;
uint8_t* pLab = lab + row * labStride;
_lab(pBgr, pLab, width);
for (; pBgr < pEnd; pBgr += 3, pLab += 3)
RgbToLab(pBgr[2], pBgr[1], pBgr[0], pLab);
}

#else
for (size_t row = 0; row < height; ++row)
{
const uint8_t* pBgr = bgr + row * bgrStride;
uint8_t* pLab = lab + row * labStride;
for (const uint8_t* pBgrEnd = pBgr + width * 3; pBgr < pBgrEnd; pBgr += 3, pLab += 3)
BgrToLab(pBgr[0], pBgr[1], pBgr[2], pLab);
}
#endif
}
}
}
Expand Down
106 changes: 54 additions & 52 deletions src/Simd/SimdBgrToLab.h
Original file line number Diff line number Diff line change
Expand Up @@ -24,67 +24,69 @@
#ifndef __SimdBgrToLab_h__
#define __SimdBgrToLab_h__

#include "Simd/SimdDefs.h"
#include "Simd/SimdMath.h"

namespace Simd
{
namespace Base
{
namespace LabV0
const int LAB_GAMMA_SHIFT = 3;

const int LabGammaTabSize = 256;
extern uint16_t LabGammaTab[LabGammaTabSize];

const int LabCbrtTabSize = 256 * 3 / 2 * (1 << LAB_GAMMA_SHIFT);
extern uint16_t LabCbrtTab[LabCbrtTabSize];

const int LabCoeffsTabSize = 9;
extern uint32_t LabCoeffsTab[LabCoeffsTabSize];

void LabInitTabs();

const int LAB_SHIFT = 12;
const int LAB_ROUND = 1 << (LAB_SHIFT - 1);

const int LAB_SHIFT2 = LAB_SHIFT + LAB_GAMMA_SHIFT;
const int LAB_ROUND2 = 1 << (LAB_SHIFT2 - 1);

const int LAB_L_SCALE = (116 * 255 + 50) / 100;
const int LAB_L_SHIFT = -((16 * 255 * (1 << LAB_SHIFT2) + 50) / 100);

const int LAB_A_SCALE = 500;
const int LAB_B_SCALE = 200;
const int LAB_AB_SHIFT = 128 * (1 << LAB_SHIFT2);

SIMD_INLINE int LabDescale(int value)
{
return (value + LAB_ROUND) >> LAB_SHIFT;
}

SIMD_INLINE int LabDescale2(int value)
{
#define X(xyz) xyz.data[0]
#define Y(xyz) xyz.data[1]
#define Z(xyz) xyz.data[2]

typedef union
{
double data[3];
struct { double r, g, b; };
struct { double L, A, B; };
struct { double l, c, h; };
} DoubleTriplet;

const double eps = (6 * 6 * 6) / (29.0 * 29.0 * 29.0), kap = (29 * 29 * 29) / (3.0 * 3.0 * 3.0);

const DoubleTriplet xyzReferenceValues = { {0.95047, 1.0, 1.08883} };

SIMD_INLINE DoubleTriplet xyzFromRgb(DoubleTriplet rgb)
{
for (int i = 0; i < 3; ++i)
{
double v = rgb.data[i];
rgb.data[i] = (v > 0.04045) ? pow(((v + 0.055) / 1.055), 2.4) : (v / 12.92);
}
DoubleTriplet temp = { {
rgb.r * 0.4124564 + rgb.g * 0.3575761 + rgb.b * 0.1804375,
rgb.r * 0.2126729 + rgb.g * 0.7151522 + rgb.b * 0.0721750,
rgb.r * 0.0193339 + rgb.g * 0.1191920 + rgb.b * 0.9503041
} };
return temp;
}

SIMD_INLINE DoubleTriplet labFromXyz(DoubleTriplet xyz)
{
for (int i = 0; i < 3; ++i)
{
xyz.data[i] /= xyzReferenceValues.data[i];
double v = xyz.data[i];
xyz.data[i] = (v > eps) ? pow(v, (1 / 3.0)) : ((kap * v + 16) / 116.0);
}
DoubleTriplet temp = { {(116 * Y(xyz)) - 16, 500 * (X(xyz) - Y(xyz)), 200 * (Y(xyz) - Z(xyz))} };
return temp;
}

DoubleTriplet labFromRgb(DoubleTriplet rgb) { return labFromXyz(xyzFromRgb(rgb)); }
return (value + LAB_ROUND2) >> LAB_SHIFT2;
}

SIMD_INLINE void BgrToLab(int blue, int green, int red, uint8_t* lab)
SIMD_INLINE void RgbToLab(int red, int green, int blue, uint8_t *lab)
{
LabV0::DoubleTriplet _rgb{red / 255.0 , green / 255.0, blue / 255.0 };
LabV0::DoubleTriplet _lab = labFromRgb(_rgb);
lab[0] = int(_lab.data[0]);
lab[1] = int(_lab.data[1]);
lab[2] = int(_lab.data[2]);
int R = LabGammaTab[red];
int G = LabGammaTab[green];
int B = LabGammaTab[blue];

int iX = LabDescale(R * LabCoeffsTab[0] + G * LabCoeffsTab[1] + B * LabCoeffsTab[2]);
int iY = LabDescale(R * LabCoeffsTab[3] + G * LabCoeffsTab[4] + B * LabCoeffsTab[5]);
int iZ = LabDescale(R * LabCoeffsTab[6] + G * LabCoeffsTab[7] + B * LabCoeffsTab[8]);

int fX = LabCbrtTab[iX];
int fY = LabCbrtTab[iY];
int fZ = LabCbrtTab[iZ];

int L = LabDescale2(LAB_L_SCALE * fY + LAB_L_SHIFT);
int a = LabDescale2(LAB_A_SCALE * (fX - fY) + LAB_AB_SHIFT);
int b = LabDescale2(LAB_B_SCALE * (fY - fZ) + LAB_AB_SHIFT);

lab[0] = Base::RestrictRange(L);
lab[1] = Base::RestrictRange(a);
lab[2] = Base::RestrictRange(b);
}
}
}
Expand Down

0 comments on commit 2073a5a

Please sign in to comment.