#!/usr/bin/env python3

# python3 -m pip install -U tifffile[all]
from tifffile import (imread, imwrite, TiffFile)

filename_input = "rgb-16bits.tif"
filename_output = "xyz-16bits.tif"


def extract_and_convert(filename_input):

	with TiffFile(filename_input) as tif:
		for page_index, page in enumerate(tif.pages):
			for line_index, line in enumerate(page.asarray()):
				for pixel_index, pixel in enumerate(line):

					print(f"==== Pixel {page_index}{line_index}{pixel_index} =====")

					print(f"Read 16 bits\t:\t", end="")
					print(f"Red = 0x{pixel[0]:04x}\t\tGreen = 0x{pixel[1]:04x}\t\tBlue = 0x{pixel[2]:04x}")

					"""
					to divide by the maximum value of 16 bits, you put your value inside [0...1] range
					where   0x0000 = 0.0000
							0xFF00 = 0.5000
							0xFFFF = 1.0000
					with this value, you can now move this value inside any space !
					if you want to put this value inside :
					8 bits space  =   <value> * 0xFF
					12 bits space =   <value> * 0xFFF
					16 bits space =   <value> * 0xFFFF
					32 bits space =   <value> * 0xFFFFFFFF
					"""

					# Map [0...1] range
					print(f"Into Map [0..1]\t:\t", end="")
					red   = (pixel[0] / 0xFFFF)
					green = (pixel[1] / 0xFFFF)
					blue  = (pixel[2] / 0xFFFF)
					print(f"Red = {red:<12f}\tGreen = {green:<12f}\tBlue = {blue:<12f}")

					# Gamma 2.6
					print(f"Gamma 2.6\t:\t", end="")
					red   = pow(red, 2.6)
					green = pow(green, 2.6)
					blue  = pow(blue, 2.6)
					print(f"Red = {red:<12f}\tGreen = {green:<12f}\tBlue = {blue:<12f}")

					# RGB -> XYZ
					# Matrix based SMPTE RP 431-2-2011 - D-Cinema Quality - Reference Projector and Environment
					print("RGB to XYZ\t:\t", end="")
					x = (red * 0.4451698156) + (green * 0.2771344092) + (blue * 0.1722826698)
					y = (red * 0.2094916779) + (green * 0.7215952542) + (blue * 0.0689130679)
					z = (red * 0.0000000000) + (green * 0.0470605601) + (blue * 0.9073553944)
					print(f" X  = {x:<12f}\t   Y  = {y:<12f}\t  Z  = {z:<12f}")

					# Gamma 1/2.6
					print("Gamma 1/2.6\t:\t", end="")
					x = pow(x, 1/2.6)
					y = pow(y, 1/2.6)
					z = pow(z, 1/2.6)
					print(f" X' = {x:<12f}\t   Y' = {y:<12f}\t  Z' = {z:<12f}")

					"""
					to multiply with 0xFFFF, you convert your [0...1] to 16 bits.
					if you multiply with 0xFF, you can couvert to 8 bits
					"""
					# Out of Map [0...1] range
					print(f"Out Map [0..1]\t:\t", end="")
					x = (x * 0xFFFF)
					y = (y * 0xFFFF)
					z = (z * 0xFFFF)
					print(f" X' = {x:<12f}\t   Y' = {y:<12f}\t  Z' = {z:<12f}")

					# Rounded
					x = int(round(x))
					y = int(round(y))
					z = int(round(z))
					print(f"Round+Int (dec)\t:\t X' = {x:<12f}\t   Y' = {y:<12f}\t  Z' = {z:<12f}")
					print(f"Round+Int (hex)\t:\t x' = 0x{x:04x}\t\t   Y' = 0x{y:04x}\t\t  Z' = 0x{z:04x}")

					# push 6 octets (2 octets for each components r,g,b) to bytes format
					# needed for yield
					yield x.to_bytes(2, byteorder='little') \
						+ y.to_bytes(2, byteorder='little') \
						+ z.to_bytes(2, byteorder='little')



# get informations about input file
image = imread(filename_input)

# shape = (width, height, number_of_components)
print("shape =", image.shape)
# dtype = size_of_component (ex: uint16)
print("dtype =", image.dtype)

print(f"\nwrite to {filename_output}... {image.shape} => {image.dtype}")
imwrite(
	filename_output,
	data=extract_and_convert(filename_input),    # iterate is better for memory
	photometric='rgb',   # don't care because XYZ have 3 components too
	shape=image.shape,
	dtype=image.dtype,
	metadata=None
)