#! /usr/bin/perl -w # # A blackbody color tool hack # http://www.vendian.org/mncharity/dir3/blackbody/ # Mitchell Charity # 2001-Jun-22 # Thanks to http://www.fourmilab.ch/documents/specrend/ . # use POSIX qw(floor ceil pow exp abs); sub min { $_[0] <= $_[1] ? $_[0] : $_[1] } sub max { $_[0] >= $_[1] ? $_[0] : $_[1] } sub rint { my $f = floor($_[0]); ($_[0] - $f) < 0.5 ? $f : ($f + 1) } sub clamp {my($v,$l,$h) = @_; (($v)<($l) ? ($l) : ($v) > ($h) ? ($h) : ($v))}; sub feq { my($a,$b,$e)=@_; die if !defined $e; &abs($a-$b) <= $e } use Math::MatrixReal; use Math::MatrixReal::Ext1; #=== RGB Color Systems ======================================== $IlluminantD65 = [&xyz_from_xy( 0.312713, 0.329016 )]; # Daylight standard $IlluminantD50 = [&xyz_from_xy( 0.3457 , 0.3585 )]; # "Bright Incandescent Light" $IlluminantE = [&xyz_from_xy( 1/3 , 1/3 )]; # "Normalized Reference 5500 K" $IlluminantB = [&xyz_from_xy( 0.3840 , 0.3516 )]; # "Direct sunlight 4874 K" $IlluminantD65_another = [qw(0.3127268660 0.3290235126 0.3582496214)]; # http://www.srgb.com/PostScript/PSdescription.htm # CIE publication 15.2 for the 1 nm 2 degree observer. $IlluminantD65_sRGB = [ 0.3127, 0.3291 , 0.3583 ]; # http://www.srgb.com/sRGBstandard.pdf Adds up to 1.0001 !?! $IlluminantD65_sRGB_1 = [ 0.3127, 0.3290 , 0.3583 ]; $SystemsRGB = { # Same primaries as Rec709. name => "sRGB system (D65)", Rxyz => [&xyz_from_xy( 0.6400, 0.3300 )], Gxyz => [&xyz_from_xy( 0.3000, 0.6000 )], Bxyz => [&xyz_from_xy( 0.1500, 0.0600 )], Wxyz => $IlluminantD65_sRGB_1 }; $SystemRec709 = { name => "Rec709 system (D65)", Rxyz => [&xyz_from_xy( 0.640, 0.330 )], Gxyz => [&xyz_from_xy( 0.300, 0.600 )], Bxyz => [&xyz_from_xy( 0.150, 0.060 )], Wxyz => $IlluminantD65 }; $SystemEBU = { name => "EBU system (D65)", Rxyz => [&xyz_from_xy( 0.64, 0.33 )], Gxyz => [&xyz_from_xy( 0.29, 0.60 )], Bxyz => [&xyz_from_xy( 0.15, 0.06 )], Wxyz => $IlluminantD65 }; $SystemSMPTE = { name => "SMPTE-C RGB system (D65)", Rxyz => [&xyz_from_xy( 0.630, 0.340 )], Gxyz => [&xyz_from_xy( 0.310, 0.595 )], Bxyz => [&xyz_from_xy( 0.155, 0.070 )], Wxyz => $IlluminantD65 }; $SystemOldFAQ = { Rxyz => [0.64, 0.33, 0.03], Gxyz => [0.29, 0.60, 0.11], Bxyz => [0.15, 0.06, 0.79], Wxyz => [0.312713, 0.329016, 0.358271] }; sub useSystem { $System = $_[0]; $Map_XYZ_to_RGB = undef; } &useSystem($SystemsRGB); $suppress_nouse_warning = $SystemEBU; $suppress_nouse_warning = $SystemOldFAQ; $suppress_nouse_warning = $IlluminantD65_another; $suppress_nouse_warning = $IlluminantD65_sRGB; $suppress_nouse_warning = $Map_XYZ_to_RGB__sRGB; $OtherGamma = 2.3; sub gammacorrect_other { pow($_[0],(1.0/$OtherGamma)); } sub gammacorrect_45 { pow($_[0],0.45); } sub gammacorrect_simple { pow($_[0],(1.0/2.3)); } # 0.43 sub gammacorrect_Rec709 { my $x = shift; ((($x) > 0.018) ? (1.099 * pow($x,0.45) - 0.099) : (4.5 * $x)) } sub gammacorrect_sRGB { # curve ends up looking like 2.2 my $x = shift; ((($x) > 0.00304) ? (1.055 * pow($x,(1.0/2.4)) - 0.055) : (12.92 * $x)) } $GammaCorrect = \&gammacorrect_sRGB; # If you use this, don't forget to set the System also to Rec709. $Map_XYZ_to_RGB__ColorFAQ_Rec709 = &matrix_from_rows([ 3.240479 , -1.537150 , -0.498535 ], [-0.969256 , 1.875992 , 0.041556 ], [ 0.055648 , -0.204043 , 1.057311 ]); # http://www.w3.org/Graphics/Color/sRGB.html # http://www.srgb.com/sRGBstandard.pdf $Map_XYZ_to_RGB__sRGB = &matrix_from_rows([ 3.2410 , -1.5374 , -0.4986 ], [-0.9692 , 1.8760 , 0.0416 ], [ 0.0556 , -0.2040 , 1.0570 ]); # Deriving matrix from primaries may have some 0.0001 errors. # And depends on the version of D65 used. # $IlluminantD65_sRGB_1 matches. # $IlluminantD65 mumble has some errors. # $IlluminantD65_another gives _worse_ errors. $Map_XYZ_to_RGB = $Map_XYZ_to_RGB__sRGB; sub make_XYZ_RGB_conversion_matrix { my $test = 1; my @wxyY = (@{$System->{Wxyz}}[0,1],1.0); my $WXYZ = &matrix_from_cols([&XYZ_from_xyY(@wxyY)]); print "WXYZ = \n$WXYZ\n" if $test; my $rgbxyz = &matrix_from_cols($System->{Rxyz}, $System->{Gxyz}, $System->{Bxyz}); print "rgbxyz = \n$rgbxyz\n" if $test; my($dim,$rgbcoefs) = $rgbxyz->decompose_LR()->solve_LR($WXYZ); die if !defined $dim; print "rgbcoefs = \n$rgbcoefs\n" if $test; my $matrix_RGB_to_XYZ = &matrix_from_cols( $rgbxyz->column(1) * $rgbcoefs->element(1,1), $rgbxyz->column(2) * $rgbcoefs->element(2,1), $rgbxyz->column(3) * $rgbcoefs->element(3,1), ); print "M = \n$matrix_RGB_to_XYZ\n" if $test; my $matrix_XYZ_to_RGB = $matrix_RGB_to_XYZ->decompose_LR()->invert_LR(); print "Mi = \n$matrix_XYZ_to_RGB\n" if $test; $Map_XYZ_to_RGB = $matrix_XYZ_to_RGB; } #=== Color Matching Functions ======================================== foreach (`/bin/cat CMFs/ciexyzjv.txt`) { my($nm,$rbar,$gbar,$bbar) = split(/,\s*/o,$_); next if !defined $bbar; $ciexyzjv{$nm} = [$rbar,$gbar,$bbar]; } foreach (`/bin/cat CMFs/ciexyz31.txt`) { my($nm,$rbar,$gbar,$bbar) = split(/,\s*/o,$_); next if !defined $bbar; $ciexyz31{$nm} = [$rbar,$gbar,$bbar]; } foreach (`/bin/cat CMFs/ciexyz64.txt`) { my($nm,$rbar,$gbar,$bbar) = split(/,\s*/o,$_); next if !defined $bbar; $ciexyz64{&rint($nm)} = [$rbar,$gbar,$bbar]; } sub xyzBar_from_nm__helper { my($h,$nm) = @_; if(!defined $h->{$nm}) { return (0,0,0) if $nm < 360 || $nm > 830; my $dn = &floor($nm / 5) * 5; my $up = &ceil($nm / 5) * 5; my @dnval = @{$h->{$dn}}; my @upval = @{$h->{$up}}; my $mix = ($nm - $dn) / 5; my @val; for(my $i = 0; $i < 3; $i++) { $val[$i] = (((1 - $mix) * $dnval[$i]) + $mix * $upval[$i]); } return @val; } @{$h->{$nm}}; } sub xyzBar_from_nm__ciexyzjv { &xyzBar_from_nm__helper(\%ciexyzjv,shift); } sub xyzBar_from_nm__ciexyz31 { &xyzBar_from_nm__helper(\%ciexyz31,shift); } sub xyzBar_from_nm__ciexyz64 { &xyzBar_from_nm__helper(\%ciexyz64,shift); } $xyzBar_from_nm = \&xyzBar_from_nm__ciexyzjv; #=== Spectra ======================================== #fff1ed Rec709 Hamilton #fff1ed Rec709 CIE'31 original #fff1eb Rec709 CIE'64 #fff1eb xxx CIE'64 ; sRGB w Rec709 gamma #fff2e8 Rec709 CIE'31 Judd Vos #fff2e8 xxx CIE'31 Judd Vos; sRGB w Rec709 gamma #fff2ee sRGB CIE'64 #fff3ea sRGB CIE'31 Judd Vos #fff3ef sRGB CIE'31 original #dce5ff Rec709 CIE'31 Judd Vos; D50 #e0e8ff sRGB CIE'31 Judd Vos; D50 # Re Rec709 vs sRGB - its the gamma that makes the difference. sub XYZ_from_Sun { # Solar spectrum with no atmospheric filtering. # wehrli85.txt from http://rredc.nrel.gov/solar/standards/am0/ # wehrli85_.txt just has unix newlines. my($sumX,$sumY,$sumZ,$sumCount); my $nm_prev = 0; foreach (`/bin/cat OtherData/wehrli85_.txt`) { /^\s*(\d\S+)\s+(\S+)/o || next; my($nm,$P) = ($1,$2); # some of the data is every 1 nm, some every 2 nm. $P = $P * ($nm - $nm_prev); $nm_prev = $nm; next if $nm < 380 || $nm > 825; my($xBar,$yBar,$zBar) = $xyzBar_from_nm->($nm); my $X = $P * $xBar; my $Y = $P * $yBar; my $Z = $P * $zBar; $sumX += $X; $sumY += $Y; $sumZ += $Z; } $sumCount = 825 - 380; my $X = $sumX / $sumCount; my $Y = $sumY / $sumCount; my $Z = $sumZ / $sumCount; return ($X,$Y,$Z); } sub XYZ_from_Kurucz { my($file) = @_; my($sumX,$sumY,$sumZ,$sumCount); my $nm_prev = 380 -1; my @lines = `/bin/cat $file | /bin/gunzip`; while ($lines[0] !~ /^\s*\d/o) { shift(@lines); } foreach (@lines) { /^\s*(\d\S+)\s+(\S+)/o || die "failed $file line $_"; my($A,$P) = ($1,$2); my $nm = $A / 10; next if $nm < 380 || $nm > 825; $P = $P * ($nm - $nm_prev); $nm_prev = $nm; my($xBar,$yBar,$zBar) = $xyzBar_from_nm->($nm); my $X = $P * $xBar; my $Y = $P * $yBar; my $Z = $P * $zBar; $sumX += $X; $sumY += $Y; $sumZ += $Z; } $sumCount = 825 - 380; my $X = $sumX / $sumCount; my $Y = $sumY / $sumCount; my $Z = $sumZ / $sumCount; return ($X,$Y,$Z); } sub add_to_sums_nm { my($nm,$P) = @_; my($xBar,$yBar,$zBar) = $xyzBar_from_nm->($nm); my $X = $P * $xBar; my $Y = $P * $yBar; my $Z = $P * $zBar; $sumX += $X; $sumY += $Y; $sumZ += $Z; } sub XYZ_from_Silva { my($file) = @_; my $x = `/bin/cat $file`; $x =~ s/^[^\n]+\n *//; my $startA = 3505; my $datastep = 0.5; my $nm = ($startA / 10) - $datastep; local($sumX,$sumY,$sumZ,$sumCount); foreach (split(/\s+/o,$x)) { die "failed $file entry $_" if !/^\d/o; $nm += $datastep; my $P = $_ * $datastep; next if $nm < 380 || $nm > 825; &add_to_sums_nm($nm,$P); } $sumCount = 825 - 380; my $X = $sumX / $sumCount; my $Y = $sumY / $sumCount; my $Z = $sumZ / $sumCount; return ($X,$Y,$Z); } sub XYZ_from_Pickles { my($file) = @_; local($sumX,$sumY,$sumZ,$sumCount); my @lines = `/bin/cat $file | /bin/gunzip`; while ($lines[0] !~ /^\s*\d/o) { shift(@lines); } my $nm_prev = 380 - 0.5; foreach (@lines) { /^\s*(\d\S+)\s+(\S+)/o || die "failed $file line $_"; my($A,$P) = ($1,$2); my $nm = $A / 10; next if $nm < 380 || $nm > 825; $P = $P * ($nm - $nm_prev); $nm_prev = $nm; &add_to_sums_nm($nm,$P); } $sumCount = 825 - 380 + 0.5; my $X = $sumX / $sumCount; my $Y = $sumY / $sumCount; my $Z = $sumZ / $sumCount; return ($X,$Y,$Z); } sub lines_from_file { my($path) = @_; my $cmd = "/bin/cat $path"; $cmd .= "|gunzip" if $path =~ /\.gz/; $cmd .= ">deleteme; fdump columns=\"1,2\" rows=- pagewidth=255 showrow=no prhead=no page=no deleteme STDOUT;rm deleteme" if $path =~ /\.fits/; my @lines = `$cmd`; return @lines; } sub XYZ_from_Generic { my($file,$start_nm,$start_gap,$unitconvert) = @_; local($sumX,$sumY,$sumZ,$sumCount); my @lines = &lines_from_file($file); while (@lines && $lines[0] !~ /^\s*\d/o) { shift(@lines); } my $nm_prev; $nm_prev = $start_nm - $start_gap if defined $start_nm; foreach (@lines) { /^\s*(\d\S+)\s+(\S+)/o || do{ warn "skipping $file line $_"; next ; }; my($A,$P) = ($1,$2); my $nm = $A / 10; $nm_prev = $nm if !defined $nm_prev;#will lose entry if 1st entry is in range. next if $nm < 380 || $nm > 825; $width = ($nm - $nm_prev); $nm_prev = $nm; # print STDERR sprintf(" %.4f nm %.4e ",$nm,$P); # print STDERR "<$nm\t$P\n"; $P = $P * $width; # print STDERR sprintf(" -> %.4e (%.3f nm)",$P,$width); $P = &$unitconvert($nm,$P) if defined $unitconvert; # print STDERR sprintf(" -> %.4e\n",$P); # print STDERR ">$nm\t$P\n"; &add_to_sums_nm($nm,$P); } $start_gap = 0 if !defined $start_gap; $sumCount = 825 - 380 + $start_gap; my $X = $sumX / $sumCount; my $Y = $sumY / $sumCount; my $Z = $sumZ / $sumCount; return ($X,$Y,$Z); } sub energy_of_photon_at_nm { my($nm) = @_; my $h = 6.626176e-34; my $c = 2.997925e8; my $erg = 1.0e7; return ($h * $c / ($nm * 1.0e-9)) * $erg; } sub make_blackbody_sampler_for_K { my($K) = @_; # K -> nm -> W/m^2 per unit wavelength (ie, /dLambda) return sub { my($nm) = @_; my $wavelength = $nm * 1.0e-9; my $h = 6.626176e-34; my $c = 2.997925e8; my $k = 1.380662e-23; my $hc = $h * $c; my $constant1 = 2.0 * 3.1415927 * $c * $h * $c; my $constant2 = $hc / $k; my $P = ($constant1 * ( (1.0 / &pow($wavelength,5)) / ( &exp($constant2 / ($wavelength * $K)) - 1.0 ))); return $P; } } # Design note - really should do the integration over wavelength # inside the spectrum defining function. Else XYZ_from_spectrum # should be called XYZ_from_smooth_spectrum. sub XYZ_from_spectrum { my($spectrum,$nm_min,$nm_max) = @_; $nm_min = 380 if !defined $nm_min; $nm_max = 825 if !defined $nm_max; # Hardwired interval of 5 nm to match CMF tables. die if $nm_min % 5 != 0; die if $nm_max % 5 != 0; my($sumX,$sumY,$sumZ); for(my $nm = $nm_min; $nm <= $nm_max; $nm += 5) { my $P = $spectrum->($nm); $P = $P * 0.5 if $nm == $nm_min || $nm == $nm_max; my($xBar,$yBar,$zBar) = $xyzBar_from_nm->($nm); my $X = $P * $xBar; my $Y = $P * $yBar; my $Z = $P * $zBar; $sumX += $X; $sumY += $Y; $sumZ += $Z; } my $X = $sumX; my $Y = $sumY; my $Z = $sumZ; return ($X,$Y,$Z); } #=== Color conversion stuff ========================================- sub xy_from_xyx { my($x,$y,$z) = @_; return ($x,$y); } sub xyz_from_xy { my($x,$y) = @_; my $z = 1.0 - ($x + $y); return ($x,$y,$z); } sub xyY_from_XYZ { my($X,$Y,$Z) = @_; my $plusXYZ = $X + $Y + $Z; my $x = $X / $plusXYZ; my $y = $Y / $plusXYZ; return ($x,$y,$Y); } sub XYZ_from_xyY { my($x,$y,$Y) = @_; my $X = $x * $Y / $y; my $Z = (1.0 - $x - $y) * $Y / $y; return ($X,$Y,$Z); } sub RGBLinear_from_XYZ { my($X,$Y,$Z) = @_; if(defined $RGBLinear_from_XYZ__use_specrend_p && $RGBLinear_from_XYZ__use_specrend_p) { return &RGBLinear_from_XYZ__specrend(@_); } &make_XYZ_RGB_conversion_matrix if !defined $Map_XYZ_to_RGB; my $e = &matrix_from_cols([$X,$Y,$Z]); my $f = $Map_XYZ_to_RGB * $e; return &values_from_vector($f); } sub RGBLinearInGamut_from_RGBLinear__simple { my($R,$G,$B) = @_; return ($R,$G,$B) if $R >= 0 && $G >= 0 && $B >= 0; warn "clipping to gamut"; $R = &max(0,$R); $G = &max(0,$G); $B = &max(0,$B); return ($R,$G,$B); } sub RGBLinear_is_in_gamut_p { my($R,$G,$B) = @_; return $R >= 0 && $G >= 0 && $B >= 0; } sub RGBLinearInGamut_from_XYZ { my($X,$Y,$Z) = @_; my($R,$G,$B) = &RGBLinear_from_XYZ($X,$Y,$Z); return ($R,$G,$B) if &RGBLinear_is_in_gamut_p($R,$G,$B); die if $Y > 100; my($wx,$wy) = @{$System->{Wxyz}}; my($wX,$wY,$wZ) = &XYZ_from_xyY($wx,$wy,100.0); my($wR,$wG,$wB) = &RGBLinear_from_XYZ($wX,$wY,$wZ); # Algorithm taken from specrend. Unexamined!! # Find the primary with negative weight and calculate the # parameter of the point on the vector from the white point # to the original requested colour in RGB space. */ my $par; if ($R < $G && $R < $B) { $par = $wR / ($wR - $R); } elsif ($G < $R && $G < $B) { $par = $wG / ($wG - $G); } else { $par = $wB / ($wB - $B); } # Now finally calculate the gamut-constrained RGB weights. */ my $gR = $wR + $par * ($R - $wR); my $gG = $wG + $par * ($G - $wG); my $gB = $wB + $par * ($B - $wB); return ($gR,$gG,$gB); } sub RGBTick01_from_RGBLinear01 { my($R,$G,$B) = @_; die if $R < 0 || $R > 1; die if $G < 0 || $G > 1; die if $B < 0 || $B > 1; my $r = $GammaCorrect->($R); my $g = $GammaCorrect->($G); my $b = $GammaCorrect->($B); return ($r,$g,$b); } sub RGBTick0255_from_RGBTick01 { my($r,$g,$b) = @_; return (rint($r * 255), rint($g * 255), rint($b * 255)); } # convenience aggregates of above sub RGBLinear01_from_XYZ { my(@XYZ) = @_; my(@RGB0) = &RGBLinearInGamut_from_XYZ(@XYZ); my(@RGB01) = &normalize_values(@RGB0); return @RGB01; } sub RGBTick0255_from_RGBLinear01 { &RGBTick0255_from_RGBTick01(&RGBTick01_from_RGBLinear01(@_)); } sub RGBTick0255_from_xyY { my($x,$y,$Y) = @_; my(@XYZ) = &XYZ_from_xyY($x,$y,$Y); my(@RGB01) = &RGBLinear01_from_XYZ(@XYZ); return &RGBTick0255_from_RGBLinear01(@RGB01); } sub RGBTick0255_from_XYZ { my(@XYZ) = @_; my(@RGB01) = &RGBLinear01_from_XYZ(@XYZ); return &RGBTick0255_from_RGBLinear01(@RGB01); } sub RGBTick0255_from_blackbody_K { my($K) = @_; my $bbr = &make_blackbody_sampler_for_K($K); my(@XYZ) = &XYZ_from_spectrum($bbr,380,825); @XYZ = &normalize_values(@XYZ); my($x,$y) = &xyY_from_XYZ(@XYZ); my(@RGB) = &RGBLinear01_from_XYZ(@XYZ); my(@rgb) = &RGBTick0255_from_RGBLinear01(@RGB); return @rgb; } # specrend's algorithm sub RGBLinear_from_XYZ__specrend { my($X,$Y,$Z) = @_; my @r = @{$System->{Rxyz}}; my @g = @{$System->{Gxyz}}; my @b = @{$System->{Bxyz}}; my($x,$y) = &xyY_from_XYZ($X,$Y,$Z); my $z = 1-$x-$y; my @c = ($x,$y,$z); my(@ret, $d); $d = $r[0]*$g[1]*$b[2] - $g[0]*$r[1]*$b[2] - $r[0]*$b[1]*$g[2] + $b[0]*$r[1]*$g[2] + $g[0]*$b[1]*$r[2] - $b[0]*$g[1]*$r[2]; $ret[0] = (-$g[0]*$c[1]*$b[2] + $c[0]*$g[1]*$b[2] + $g[0]*$b[1]*$c[2] - $b[0]*$g[1]*$c[2] - $c[0]*$b[1]*$g[2] + $b[0]*$c[1]*$g[2]) / $d; $ret[1] = ($r[0]*$c[1]*$b[2] - $c[0]*$r[1]*$b[2] - $r[0]*$b[1]*$c[2] + $b[0]*$r[1]*$c[2] + $c[0]*$b[1]*$r[2] - $b[0]*$c[1]*$r[2]) / $d; $ret[2] = ($r[0]*$g[1]*$c[2] - $g[0]*$r[1]*$c[2] - $r[0]*$c[1]*$g[2] + $c[0]*$r[1]*$g[2] + $g[0]*$c[1]*$r[2] - $c[0]*$g[1]*$r[2]) / $d; return @ret; #if(0) { # Easier to be derivative than careful # $_specrend_calc = ' # d = xr*yg*zb - xg*yr*zb - xr*yb*zg + xb*yr*zg + xg*yb*zr - xb*yg*zr; # *r = (-xg*yc*zb + xc*yg*zb + xg*yb*zc - xb*yg*zc - xc*yb*zg + xb*yc*zg) / d; # *g = (xr*yc*zb - xc*yr*zb - xr*yb*zc + xb*yr*zc + xc*yb*zr - xb*yc*zr) / d; # *b = (xr*yg*zc - xg*yr*zc - xr*yc*zg + xc*yr*zg + xg*yc*zr - xc*yg*zr) / d; # '; # %__tmp = (x=>0,y=>1,z=>2,r=>0,g=>1,b=>2); # $_specrend_calc =~ s/([XYZ])([rgbc])/"\$$2\[$__tmp{$1}\]"/ge; # $_specrend_calc =~ s/d/\$d/g; # $_specrend_calc =~ s/\*([rgb])/"\$ret\[".$__tmp{$1}."\]"/ge; # print $_specrend_calc; #} } # new Luv stuff ########### UNTESTED! Ie, DOESNT WORK. sub XYZ_of_whitepoint { my($Wx,$Wy,$Wz) = @{$System->{Wxyz}}; return &XYZ_from_xyY($Wx,$Wy,100.0); } sub CMYFake_from_RGBLinear01 { my($R,$G,$B) = @_; my($C,$M,$Y) = (1-$R,1-$G,1-$B); return ($C,$M,$Y); } sub uvTick_from_XYZ { my($X,$Y,$Z) = @_; my $uT = (4 * $X / ($X + 15 * $Y + 3 * $Z)); my $vT = (9 * $Y / ($X + 15 * $Y + 3 * $Z)); return($uT,$vT); } sub L_from_XYZ { my($X,$Y,$Z) = @_; my($Xn,$Yn,$Zn) = &XYZ_of_whitepoint(); my $L = ((($Y/$Yn)>0.008856) ? (116 * &pow(($Y/$Yn),(1/3)) - 16) : 903.3*($Y/$Yn) ); #check-2nd src didnt mention this part. return $L; } sub LuvTick_from_XYZ { my(@XYZ) = @_; return (&L_from_XYZ(@XYZ), &uvTick_from_XYZ(@XYZ)); } sub Luv_from_XYZ { my(@XYZ) = @_; my $L = &L_from_XYZ(@XYZ); my(@XYZn) = &XYZ_of_whitepoint(); my($uT,$vT) = &uvTick_from_XYZ(@XYZ); my($uTn,$vTn) = &uvTick_from_XYZ(@XYZn); my $u = 13*($L)*($uT - $uTn); my $v = 13*($L)*($vT - $vTn); return ($L,$u,$v); } sub xy_from_uvTick { #http://www.cs.rit.edu/~ncs/color/t_convert.html#XYZ my($uT,$vT) = @_; my $x = 27 * $uT / ( 18 * $uT - 48 * $vT + 36 ); my $y = 12 * $vT / ( 18 * $uT - 48 * $vT + 36 ); return($x,$y); } #=== Helpers ========================================- sub matrix_bug_workaround { # Sigh. Workaround bug in matrix package. Worrisome. # 1E2 is ok, but 1e2 goes boom. my(@tuples) = @_; foreach $tuple (@tuples) { if(ref $tuple eq 'ARRAY') { foreach $value (@{$tuple}) { $value =~ tr/e/E/ if !ref $value; } } } @tuples; } sub matrix_from_cols { Math::MatrixReal::Ext1->new_from_cols([&matrix_bug_workaround(@_)]) } sub matrix_from_rows { Math::MatrixReal::Ext1->new_from_rows([&matrix_bug_workaround(@_)]) } sub values_from_vector { my($m) = @_; return ( $m->element(1,1), $m->element(2,1), $m->element(3,1) ); } sub normalize_values { my(@values) = @_; my $max = $values[0]; foreach (@values) { $max = $_ if $max < $_; } return (map { $_ / $max } @values); } sub countFromToNSamples { my($from,$to,$nsamples) = @_; my $step = ($to - $from) / $nsamples; map { ($_ * $step) + $from } (0 .. ($nsamples - 1)); } sub ppm_from_whPixelRow { my($w,$h,$pixelrow) = @_; my $ppm = "P3 \n# temp.ppm\n$w $h\n255\n"; my $scale = $w / @{$pixelrow}; my $line = ""; for(my $c=0; $c < $w; $c++) { my($r,$g,$b) = @{$pixelrow->[&clamp(rint($c/$scale),0,$#$pixelrow)]}; $line .= "$r $g $b\n"; } for(my $r=0; $r < $h; $r++) { $ppm .= $line; } $ppm; } #=== a simple image abstraction sub img_rgb_pack { my($r,$g,$b) = @_; die if $r < 0 || $r > 255; die if $g < 0 || $g > 255; die if $b < 0 || $b > 255; die if $r != rint($r); die if $g != rint($g); die if $b != rint($b); pack("CCC",$r,$g,$b); } sub img_rgb_unpack { my($s) = @_; unpack("CCC",$s); } sub img_make { my($w,$h) = @_; my(@cols,@arow); $arow[$w-1]=undef; $cols[$h-1]=\@arow; return \@cols; } sub img_pixel_set { my($img,$xcoor,$ycoor,$r,$g,$b) = @_; $img->[$ycoor]->[$xcoor] = &img_rgb_pack($r,$g,$b); } sub img_pixel_get { my($img,$xcoor,$ycoor) = @_; $img->[$ycoor]->[$xcoor]; } sub img_size { my($img) = @_; my $h = @{$img} +0; my $w = 0; foreach $row (@{$img}) { next if !defined $row; my $rw = @{$row} +0; $w = &max($w,$rw); } return ($w,$h); } sub img_to_ppm { my($img,$default_rgb) = @_; my($w,$h) = &img_size($img); $default_rgb = [0,0,0] if !defined $default_rgb; my $ppm = "P6 \n# temp.ppm\n$w $h\n255\n"; for(my $i = ($h - 1); $i >= 0; $i--) { for(my $j = 0; $j < $w; $j++) { my $pixel; $pixel = $img->[$i][$j] if defined $img->[$i]; my($r,$g,$b) = (defined $pixel ? &img_rgb_unpack($pixel) : @{$default_rgb}); $ppm .= pack("CCC",$r,$g,$b); } } $ppm; } sub img_from_gamut { my($step) = @_; $step = 0.01 if !defined $step; die if rint(0.1 / $step) * $step != 0.1; my $img = &img_make(1.0/$step,1.0/$step); my $xy2ij = sub { my($x,$y)=@_;return (floor($x/$step),floor($y/$step)); }; # grid for(my $flip = 0; $flip <= 1; $flip++) { for(my $y = 0.9; $y > 0.0; $y -= 0.1) { for(my $x = 1.0-$step; $x > 0.0; $x -= $step) { &img_pixel_set($img, $flip ? &$xy2ij($x,$y) : &$xy2ij($y,$x), 64,64,64); } } } # gamut triangle for(my $y = 1.0; $y > 0.0; $y -= $step) { for(my $x = 1.0; $x > 0.0; $x -= $step) { my $Y = 100.0; my(@XYZ) = &XYZ_from_xyY($x,$y,$Y); my(@RGB) = &RGBLinear_from_XYZ(@XYZ); next if !&RGBLinear_is_in_gamut_p(@RGB); my(@rgb0255) = &RGBTick0255_from_xyY($x,$y,$Y); &img_pixel_set($img,&$xy2ij($x,$y),@rgb0255); } } # blackbody curve for(my $K = 1000; $K <= 20000; $K += 100) { my $bbr = &make_blackbody_sampler_for_K($K); my(@XYZ) = &XYZ_from_spectrum($bbr,380,825); my($x,$y) = &xyY_from_XYZ(@XYZ); &img_pixel_set($img,&$xy2ij($x,$y),128,128,128); } return $img; } #=== more helpers sub bbr_table { my($maxval,$step) = @_; my $s = ""; for(my $K = 1000; $K <= $maxval; $K += $step) { my $bbr = &make_blackbody_sampler_for_K($K); my(@XYZ) = &XYZ_from_spectrum($bbr,380,825); @XYZ = &normalize_values(@XYZ); my($x,$y) = &xyY_from_XYZ(@XYZ); my(@RGB) = &RGBLinear01_from_XYZ(@XYZ); my(@rgb) = &RGBTick0255_from_RGBLinear01(@RGB); $s .= sprintf(" %5d K %.4f %.4f %.4f %.4f %.4f %.4f %.4f %.4f %3d %3d %3d #%02x%02x%02x\n", $K,$x,$y,@XYZ,@RGB,@rgb,@rgb); } $s; } #=== Output I needed for the webpages ======================================== sub want ($) { my $tf = @ARGV && $ARGV[0] eq $_[0]; shift(@ARGV) if $tf; $tf} if(want '-specrend') { $RGBLinear_from_XYZ__use_specrend_p = 1; &useSystem($SystemSMPTE); $xyzBar_from_nm = \&xyzBar_from_nm__ciexyz31; } if(want '-specrend_compare') { $RGBLinear_from_XYZ__use_specrend_p = 1; } if(want '-Rec709') { &useSystem($SystemRec709); $GammaCorrect = \&gammacorrect_Rec709; } if(want '-faqmatrix') { &useSystem($SystemRec709); $Map_XYZ_to_RGB = $Map_XYZ_to_RGB__ColorFAQ_Rec709; } if(want '-sRGB') { &useSystem($SystemsRGB); $GammaCorrect = \&gammacorrect_sRGB; } if(want '-sEBU') { &useSystem($SystemEBU); warn "Unknown gammacorrection for EBU"; } if(want '-sSMPTE') { &useSystem($SystemSMPTE); warn "Unknown gammacorrection for SMPTE"; } if(want '-CIEold') { $xyzBar_from_nm = \&xyzBar_from_nm__ciexyz31; } if(want '-CIE64') { $xyzBar_from_nm = \&xyzBar_from_nm__ciexyz64; } sub set_wp { warn "whitepoint specifier better be after system spec"; my %sys = %{$System}; $sys{Wxyz} = $_[0]; &useSystem(\%sys); } if(want '-D50') { &set_wp($IlluminantD50); } if(want '-E') { &set_wp($IlluminantE); } if(want '-B') { &set_wp($IlluminantB); } if(want '-a_weird_whitepoint') { &set_wp([0.37,0.33]); } if(want '-gamma') { $GammaCorrect = \&gammacorrect_other; $OtherGamma = shift; } if(want '-gamutT') { print &img_to_ppm(&img_from_gamut(0.02)); } if(want '-gamutQ') { print &img_to_ppm(&img_from_gamut(0.005)); } if(want '-table_specrend') { for(my $K = 1000; $K <= 10000; $K += 500) { my $bbr = &make_blackbody_sampler_for_K($K); my(@XYZ) = &XYZ_from_spectrum($bbr,380,825); my($x,$y) = &xyY_from_XYZ(@XYZ); my(@RGB) = &RGBLinear_from_XYZ(@XYZ); my(@RGB0) = &RGBLinearInGamut_from_XYZ(@XYZ); printf(" %5d K %.4f %.4f %.4f %.3f %.3f %.3f",$K,$x,$y,(1-$x-$y),@RGB0); print " (Approximation)" if !&RGBLinear_is_in_gamut_p(@RGB); print "\n"; } } if(want '-table1') { print &bbr_table(10000,500); } if(want '-lookuptable') { print &bbr_table(45000,100); } if(want '-scale0') { my($from,$to,$nsamp) = (1000,10000,100); my $step = ($to - $from) / $nsamp; my(@colors); for(my $K = $from; $K <= $to; $K += $step) { push(@colors,[&RGBTick0255_from_blackbody_K($K)]); } print &ppm_from_whPixelRow(300,50,\@colors); } if(want '-scale1') { my($from,$to) = (1000,30000); if(@ARGV && $ARGV[0] =~ /^\d+$/) {$to = shift;} my $step = 100; for(my $K = $from; $K <= $to; $K += $step) { push(@colors,[&RGBTick0255_from_blackbody_K($K)]); } print &ppm_from_whPixelRow(500,50,\@colors); } if(want '-K') { my $K = shift; my $bbr = &make_blackbody_sampler_for_K($K); my(@XYZ) = &XYZ_from_spectrum($bbr,380,825); @XYZ = &normalize_values(@XYZ); my($x,$y) = &xyY_from_XYZ(@XYZ); my(@RGB) = &RGBLinear01_from_XYZ(@XYZ); my(@rgb) = &RGBTick0255_from_RGBLinear01(@RGB); printf(" %5d K %.4f %.4f %.4f %.4f %.4f %.4f %.4f %.4f %3d %3d %3d #%02x%02x%02x\n",$K,$x,$y,@XYZ,@RGB,@rgb,@rgb); } if(want '-whitepoints') { sub wp { my($name,$Wxyz)=@_; my($x,$y) = @{$Wxyz}; my(@rgb) = &RGBTick0255_from_xyY($x,$y,100); printf(" %s %.4f %.4f %3d %3d %3d #%02x%02x%02x\n",$name,$x,$y,@rgb,@rgb); }; &wp('D65',$IlluminantD65); &wp('D50',$IlluminantD50); &wp('E',$IlluminantE); &wp('B',$IlluminantB); } if(want '-Sun') { my(@XYZ) = &XYZ_from_Sun(); my($x,$y) = &xyY_from_XYZ(@XYZ); my(@RGB) = &RGBLinear01_from_XYZ(@XYZ); my(@rgb) = &RGBTick0255_from_RGBLinear01(@RGB); # printf(" Sun %.4f %.4f %.4f %.4f %.4f %.4f %.4f %.4f %03d %03d %03d #%02x%02x%02x\n",$x,$y,@XYZ,@RGB,@rgb,@rgb); printf(" Sun %.5f %.5f %.4f %.4f %.4f %.4f %.4f %.4f %03d %03d %03d #%02x%02x%02x\n",$x,$y,@XYZ,@RGB,@rgb,@rgb); } sub dataline_from_spectra_XYZ { my(@XYZ) = @_; my($x,$y) = &xyY_from_XYZ(@XYZ); my(@RGB) = &RGBLinear01_from_XYZ(@XYZ); my(@rgb) = &RGBTick0255_from_RGBLinear01(@RGB); sprintf(" %7s %.4f %.4f %03d %03d %03d #%02x%02x%02x\n",$class,$x,$y,@rgb,@rgb); } sub dataline2_from_spectra_XYZ { my($XYZ,$class,$file,$description) = @_; my @XYZ = @{$XYZ}; my($x,$y) = &xyY_from_XYZ(@XYZ); my(@RGB) = &RGBLinear01_from_XYZ(@XYZ); my(@rgb) = &RGBTick0255_from_RGBLinear01(@RGB); sprintf(" %7s %.4f %.4f %03d %03d %03d #%02x%02x%02x \"%s\" /%s \n",$class,$x,$y,@rgb,@rgb,$description,$file); } sub dataline3 { my($XYZ,$class,$file,$description) = @_; my @XYZ = @{$XYZ}; my($x,$y) = &xyY_from_XYZ(@XYZ); my($uT,$vT) = &uvTick_from_XYZ(@XYZ); my(@RGB) = &RGBLinear01_from_XYZ(@XYZ); my(@rgb) = &RGBTick0255_from_RGBLinear01(@RGB); sprintf(" %25s %7s %.4f %.4f %.4f %.4f %5s %3s %.4f %.4f %.4f %03d %03d %03d #%02x%02x%02x Sp %7s \"%s\"\n", $file, "CIE31jv",$x,$y,$uT,$vT, "sRGB","D65", @RGB,@rgb,@rgb,$class,$description); } sub want_spectra { my($name,$f) = @_; return if ! want "-$name"; local $file; foreach $file (`/bin/ls Spectra/$name`) { chomp($file); local $path = "Spectra/$name/$file"; local $class = $file; &$f(); } } &want_spectra('Kurucz',sub { return if $file =~ /disable/i; $file =~ m!([^\./]+)\.dat(.gz)?$!o || die; $class = $1; print &dataline_from_spectra_XYZ(&XYZ_from_Kurucz($path)); }); &want_spectra('Silva',sub { return if $file =~ /html/; print &dataline_from_spectra_XYZ(&XYZ_from_Silva($path)); }); &want_spectra('Pickles',sub { return if $file !~ /\.dat/; $class =~ s/\.dat.*//; print &dataline_from_spectra_XYZ(&XYZ_from_Pickles($path)); }); sub flam_from_photlam { my($nm,$pd) = @_; return $pd * &energy_of_photon_at_nm($nm); } &want_spectra('GunnStryker',sub { if(!defined %GunnStryker_file_class_map) { foreach (&lines_from_file('Spectra/GunnStryker/gsspectype.fits.gz')) { /^(GS_\S+)\s+(\S+)/o || next; my($file,$class) = ($1,$2); $file =~ tr/A-Z/a-z/; $file =~ s/\..+//; $GunnStryker_file_class_map{$file} = $class; } } return if $file !~ /gs_.+\.fits/; $file =~ /^([^\.]+)/o || die; my $stem = $1; $class = $GunnStryker_file_class_map{$stem}; $class = '?' if !defined $class; print &dataline_from_spectra_XYZ(&XYZ_from_Generic($path,undef,undef, \&flam_from_photlam)); }); &want_spectra('FAST',sub { if(!defined %FAST_file_class_map) { foreach (&lines_from_file('Spectra/FAST/FILES')) { /^\s*(\S+)\s+(\S+)\s*(.*)/o || next; my($file,$class,$desc) = ($1,$2,$3); $FAST_file_class_map{$file} = [$class,$desc]; } } return if $file !~ /(.+\.gif)\.dat$/; $file = $1; die if !defined $FAST_file_class_map{$file}; local($class,$desc) = @{$FAST_file_class_map{$file}}; print &dataline2_from_spectra_XYZ([&XYZ_from_Generic($path,undef,0.8)], $class,$file,$desc); # print &dataline3([&XYZ_from_Generic($path,undef,0.8)], # $class,$file,$desc); }); if(want '-Misc') { my $report = sub { my($path,$start,$startstep) = @_; my(@XYZ) = &XYZ_from_Generic($path,$start,$startstep); my($x,$y) = &xyY_from_XYZ(@XYZ); my(@RGB) = &RGBLinear01_from_XYZ(@XYZ); my(@rgb) = &RGBTick0255_from_RGBLinear01(@RGB); printf(" %7s %.4f %.4f %03d %03d %03d #%02x%02x%02x $file\n",$class,$x,$y,@rgb,@rgb,$file); }; foreach $file (`/bin/ls Spectra/Misc`) { chomp($file); next if $file !~ /\.d([LM].+)/; local $class = $1; $class =~ s/\.gz//; $class .= "(V)"; &$report("Spectra/Misc/$file",400,0.5); } foreach $file (`/bin/ls Spectra/Misc`) { chomp($file); next if $file !~ /\.d([LM].+)/; local $class = $1; $class =~ s/\.gz//; $class .= "(V)"; &$report("Spectra/Misc/$file",400,0.5); } } if(want '-bbr_by_cmfs') { my $xy_from_K = sub { my($K) = @_; my($x,$y) = &xyY_from_XYZ(&XYZ_from_spectrum(&make_blackbody_sampler_for_K($K))); return ($x,$y); }; my $s = "# gnuplot file\n# tuples; 3 datasets: ciexyzjv ciexyz31 ciexyz64\n#\n"; map { $xyzBar_from_nm = $_; for(my $K = 1000; $K <= 30000;) { $s .= join(" ",&$xy_from_K($K))."\n"; $K += 250, next if $K < 2000; $K += 500, next if $K < 4000; $K += 1000, next if $K < 8000; $K += 5000, next if 1; } $s .= "\n\n"; } (\&xyzBar_from_nm__ciexyzjv, \&xyzBar_from_nm__ciexyz31, \&xyzBar_from_nm__ciexyz64); #plot [0:1] [0:1] 'deleteme.dat' index 0 #set key right bottom print STDERR <<'END_OF_SCRIPT'; set terminal jpeg size 400,400 set key right bottom set output 'deleteme0.jpg' set size ratio -1 plot [0.2:0.7] [0.2:0.45] \ 'deleteme.dat' index 0 title "CIE 1931 Judd Vos" with points, \ 'deleteme.dat' index 1 title "CIE 1931" with points lt 4, \ 'deleteme.dat' index 2 title "CIE 1964" with points lt 2 set output 'deleteme1.jpg' plot [0:1] [0:1] \ 'deleteme.dat' index 0 title "CIE 1931 Judd Vos" with lines, \ 'deleteme.dat' index 1 title "CIE 1931" with lines lt 4, \ 'deleteme.dat' index 2 title "CIE 1964" with lines lt 2 exit END_OF_SCRIPT print $s; } if(want '-datafile') { $xyzBar_from_nm = \&xyzBar_from_nm__ciexyz64; # &useSystem($SystemsRGB); # $GammaCorrect = \&gammacorrect_sRGB; print <<'END'; # Blackbody color datafile (bbr_color.txt) # Mitchell Charity # http://www.vendian.org/mncharity/dir3/blackbody/ # Version 2001-Jun-22 # History: # 2001-Jun-22 Switched to sRGB from Rec.709. # 2001-Jun-04 Initial release. # Fields: # K temperature (K) # CMF {" 2deg","10deg"} # either CIE 1931 2 degree CMFs with Judd Vos corrections # or CIE 1964 10 degree CMFs # x y chromaticity coordinates # P power in semi-arbitrary units # R G B {0-1}, normalized, mapped to gamut, logrithmic # (sRGB primaries and gamma correction) # r g b {0-255} # #rgb {00-ff} # END my $doit = sub { my(@XYZ) = &XYZ_from_spectrum($bbr,380,825); my($x,$y,$P) = &xyY_from_XYZ(@XYZ); my(@XYZn) = &normalize_values(@XYZ); my(@RGB) = &RGBLinear01_from_XYZ(@XYZn); my(@rgb) = &RGBTick0255_from_RGBLinear01(@RGB); printf(" %5d K %s %.4f %.4f %.3e %.4f %.4f %.4f %3d %3d %3d #%02x%02x%02x\n", $K,$cmfid,$x,$y,$P,@RGB,@rgb,@rgb); }; local $cmfid; for(local $K = 1000; $K <= 40000; $K += 100) { local $bbr = &make_blackbody_sampler_for_K($K); $xyzBar_from_nm = \&xyzBar_from_nm__ciexyzjv; $cmfid = " 2deg"; &$doit(); $xyzBar_from_nm = \&xyzBar_from_nm__ciexyz64; $cmfid = "10deg"; &$doit(); } print "# end\n"; print STDERR "dont forget to regenerate color version\n"; } if(want '-intense_scale') { my $baseline_Kmax = 10000; my $baseline_width = 500; my $Kmax = $baseline_Kmax; my $width = $baseline_width; my $step = 100; if(@ARGV && $ARGV[0] =~ /^\d+$/) { $Kmax = shift; die if $Kmax > $baseline_Kmax; $width = &rint($width * (($Kmax - 1000) / ($baseline_Kmax - 1000))); $step = &rint(($Kmax - 1000) / 90); # step=100 for maxK=10000 } my(@XYZarray,@RGBarray); my $maxY = -1; for(my $K = 1000; $K <= $Kmax; $K += 25) { my(@XYZ) = &XYZ_from_spectrum(&make_blackbody_sampler_for_K($K)); push(@XYZarray,[@XYZ]); $maxY = &max($XYZ[1],$maxY); } my $maxRGB = -1; foreach (@XYZarray) { my($x,$y,$Yold) = &xyY_from_XYZ(@{$_}); my $Y = $Yold / $maxY * 100; my(@XYZ) = &XYZ_from_xyY($x,$y,$Y); my(@RGB) = &RGBLinearInGamut_from_XYZ(@XYZ); push(@RGBarray,[@RGB]); $maxRGB = &max($maxRGB,&max($RGB[0],&max($RGB[1],$RGB[2]))); } my @colors; foreach (@RGBarray) { my(@RGB) = map { $_ / $maxRGB } @{$_}; my(@rgb) = &RGBTick0255_from_RGBLinear01(@RGB); push(@colors,[@rgb]); } print &ppm_from_whPixelRow($width,50,\@colors); } if(want '-test_against_sky_paper') { # A Physically-Based Nightsky Model # http://graphics.stanford.edu/~henrik/papers/nightsky/ # page 8. Columns mislabled? Look like u'v's. my $data = ' Temperature to Luv Conversion Temp. (K) u v 100000 0.18065 0.26589 50000 0.18132 0.26845 33333 0.18208 0.27118 25000 0.18293 0.27407 20000 0.18388 0.27708 16667 0.18494 0.28020 14286 0.18611 0.28340 12500 0.18739 0.28666 11111 0.18879 0.28995 10000 0.19031 0.29325 8000 0.19461 0.30139 6667 0.19960 0.30918 5714 0.20523 0.31645 5000 0.21140 0.32309 4444 0.21804 0.32906 4000 0.22507 0.33436 3636 0.23243 0.33901 3333 0.24005 0.34305 3077 0.24787 0.34653 2857 0.25585 0.34948 2667 0.26394 0.35198 2500 0.27210 0.35405 2353 0.28032 0.35575 2222 0.28854 0.35713 2105 0.29676 0;35822 2000 0.30496 0.35906 1905 0.31310 0.35968 1818 0.32119 0.36011 1739 0.32920 0.36038 1667 0.33713 0.36051 '; foreach $line (split(/\n/,$data)) { chomp($line); $line =~ /^\s*(\d\S+)\s+(\S+)\s+(\S+)/ || next; my($K,$ref_u,$ref_v) = ($1,$2,$3); my(@XYZ) = &XYZ_from_spectrum(&make_blackbody_sampler_for_K($K)); my($L,$uT,$vT) = &LuvTick_from_XYZ(@XYZ); printf("$line %6d K %8.5f %8.5f u'v' %.4f %.4f xy\n",$K,$uT,$vT, &xyY_from_XYZ(@XYZ)); } } warn "The color library ignored ".join(" ",@ARGV) if @ARGV; 1; __END__ Doables: Time for cleanup pass. Normalization and value magnitude (esp Y) arent being handled properly. Figure out what is going on with assorted D65 chromaticities. Further test Luv code. Why doesnt -test_against_sky_paper match? Check that liberties taken with XYZ from Mumble are acceptable. Add Lab? Crush out the remaining occasional 0.0001 xyz 0.001 RGB difference in specrend emulation. (wp precision, lambda range) Verify to-gamut mapping snarfed from specrend. Add 1K/5K ticks to bbr "spectra". Fix gamut img bug which leaves grid gap and missizes image. Add gamut img labels. Attach gamut img to associated bbr colors as single image. Combine gamut imgs into animation to clearly show color shift. Notes: