8ビットBGRから64ビットCIE L*a*b*へ変換する

研究でCIE L*a*b*を使いたかったのだが、OpenCVのcvCvtColorだと

  • 8bit unsigned charだと数字が粗い
  • 32bit floatだと挙動がよくわからん(自分が悪いのかも)

ということで、自分で書いた。そしたら結構いい感じになったよ。(エントリー分けてそっちも書く)

とりあえず、適当だけどこんな感じ。

とにかく自分用メモ感がハンパないな。

/** 8bitのBGRを64bitのCIE L*a*b*に変換する */
EX_Lab64F toLab(EX_BGR8U bgr)
{
	// R,G,Bを[0,1]にスケーリング+ガンマ補正をかけてsRGBにする
	double sr = (double)bgr.R/255;
	if(sr<0.04045)
		sr /= 12.92;
	else
		sr = pow( (sr+0.055)/1.055, 2.4 );

	double sg = (double)bgr.G/255;
	if(sg<0.04045)
		sg /= 12.92;
	else
		sg = pow( (sg+0.055)/1.055, 2.4 );

	double sb = (double)bgr.B/255;
	if(sb<0.04045)
		sb /= 12.92;
	else
		sb = pow( (sb+0.055)/1.055, 2.4 );

	// XYZに変換
	double x = 0.412453*sr + 0.357580*sg + 0.180423*sb;
	double y = 0.212671*sr + 0.715160*sg + 0.072169*sb;
	double z = 0.019334*sr + 0.119193*sg + 0.950227*sb;

	x /= 0.950456;
	z /= 1.088754;

	// Labに変換
	EX_Lab64F lab;
	double fy;
	if(y>0.008856){
		fy = pow(y, 1./3);
		lab.L  = 116*fy;
	}else{
		fy = 7.787*y+16./116;
		lab.L  = 903.3*y;
	}
	lab.L -= 16;
	double fx = x > 0.008856 ? pow(x, 1./3) : 7.787*x+16./116;
	double fz = z > 0.008856 ? pow(z, 1./3) : 7.787*z+16./116;

	dst.a = 500*(fx-fy);
	dst.b = 200*(fy-fz);

	return lab;
}

/** 逆変換 */
EX_BGR8U toBGR( EX_Lab64F lab )
{
	double x,y,z;
	// XYZに変換
	if(lab.L < 7.9996){
		y = lab.L/903.3;
		x = 0.950456*(lab.a/3893.5+y);
		z = 1.088754*(y-lab.b/1557.4);
	}else{
		y = pow((lab.L+16)/116, 3.);
		x = 0.950456*pow( (lab.L+16)/116+lab.a/500, 3);
		z = 1.088754*pow( (lab.L+16)/116-lab.b/200, 3);
	}
	
	// sRGBに変換
	double r,g,b;
	r = 3.241*x -1.5374*y - 0.4986*z;
	g = -0.9692*x + 1.876*y + 0.0416*z;
	b = 0.0556*x - 0.2040*y + 1.057*z;

	// linear RGBに変換
	if( r > 0.00313080728 )
		r = pow(r, 10./24)*1.055 - 0.055;
	else
		r *= 12.92;

	if( g > 0.00313080728 )
		g = pow(g, 10./24)*1.055 - 0.055;
	else
		g *= 12.92;
	
	if( b > 0.00313080728 )
		b = pow(b, 10./24)*1.055 - 0.055;
	else
		b *= 12.92;
	


	EX_Lab8U bgr;
	bgr.L = cvRound(MAX(MIN(255*b, 255), 0));
	bgr.a = cvRound(MAX(MIN(255*g, 255), 0));
	bgr.b = cvRound(MAX(MIN(255*r, 255), 0));

	return bgr;
}