#!/usr/bin/env python
# -*- coding: utf-8 -*-

from __future__ import division
import matrix

"""
matrix elements:
P_1:=(n'-n)/R_1 refraction force

vector elements:
(n⋅α
 x)

Refraction_1 = (
  1   -P_1
  0   1)

Translation_1 = (
  1      0
  D_12/n 1)

Reflection_1 = (
  1  2⋅n/R
  0 1)

"""

def point(angle, refraction_index, height):
	return matrix.vector(angle * refraction_index, height)

# note that radius is >0 if the light ray is coming to a convex   surface.
#           radius is <0                                 concave
def refraction_factor(refraction_index, refraction_index_2, radius):
	return (refraction_index_2 - refraction_index) / radius

def refract(refraction_factor):
	return matrix.matrix(2, 2, (1, -refraction_factor,
	                            0, 1))

def refract_P(matrix_1):
	return matrix_1[1][0] == 0

def translate(distance, refraction_index):
	return matrix.matrix(2, 2, (1, 0, distance / refraction_index, 1))

def translate_P(matrix_1):
	return matrix_1[0][1] == 0

def reflect(refraction_index, radius):
	return matrix.matrix(2, 2, (1, 2 * refraction_index / radius, 0, 1))

"""
lens equation:
	for refraction:
		n/D_01 + n'/D_12 = (n'-n)/R = P_1
	for reflection:
		f':=n'⋅R/(n'-n)=n'/refraction_factor
		f:=n⋅R/(n'-n)
		#n/D_01 + n'/D_02=n/f=n'/f'
		1/D_01+1/D_12=-2/R=1/f
lateral zoom:
	m_x=-n/n' ⋅ D_12/D_01
	m_α=?
	m_x⋅m_α=n/n'
for conjugated planes, the matrix is:
	(m_α⋅n/n'  M_12
	 0         m_x)
D being distances "on" the optical axis.
h planes:
	d_1=-P_2/P_L ⋅ n_1 ⋅ d_L/n_L
	d_2=-P_1/P_L ⋅ n_2 ⋅ d_L/n_L
"""

if __name__ == "__main__":
	print translate(10, 1)
	print refract(1)
	print translate(10, 1) * refract(1)

	assert(refract_P(refract(1)))
	assert(not refract_P(translate(10, 1)))
	assert(translate_P(translate(10, 1)))
	assert(not translate_P(refract(1)))