function[re_pos1, re_pos2]=back_projection(ui_pos, P, re_plane_coeffs)
	% [re_pos1, re_pos2]=back_projection(ui_pos, P, re_plane_coeffs)
	%
	% Back projects the undistorted image space coordinates 
	% onto two planes in object space, defined by re_plane_coeffs
	% 
	% plane i is defined by
	%	re_plane_coeffs(i,1)*x + re_plane_coeffs(i,2)*y + re_plane_coeffs(i,3)*z + re_plane_coeffs(i,4) = 0
	%
	% re_pos1 is object space coordinate, assuming particle lies on plane 1
	% re_pos2 is object space coordinate, assuming particle lies on plane 2
	%
	%	Input:
	%		ui_pos				n_p x 2	undistorted image space coordinate of particle
	%		P						3 x 4		projection matrix from pinhole model
	%		re_plane_coeffs	2 x 4		coefficients of plane
	%	Output
	%		re_pos1				n_p x 3	object space position for plane 1
	%		re_pos2				n_p x 3	object space position for plane 2
	%
	
	%	Unit testing: 
	%	expect that we get ui1 = ui2 = ui_pos where
	%		ui1	= pinhole_model(re_pos1, model, false)
	%		ui2	= pinhole_model(re_pos1, model, false)

	%% input validation
	assert(size(ui_pos, 2) == 2);
	assert(size(P,1) == 3);
	assert(size(P,2) == 4);
	assert(size(re_plane_coeffs,1) == 2);
	assert(size(re_plane_coeffs,2) == 4);
	
	im_i		= ui_pos(:,1).';	% y1 coordinate of particle
	im_j		= ui_pos(:,2).';	% y2 coordinate of particle
	n_p		= size(ui_pos,1); % number of particles

	%% solve back projection in plane number 1
	% solve the back-projection of
	% lambda * [ui_pos;1] = P * [re_pos;1]
	% for the plane defined by  re_plane_coeffs(1,:) 
	%
	% solutions will satisfy
	% re_plane_coeffs(1,:) * [x;y;z;1] = 0
	% and the pinhole model projection of [x;y;z] will map to the (undistorted)
	% image coordinate ui_pos

	ij10		= [im_i; im_j; ones(1, n_p); zeros(1, n_p)];
	M1			= [P; re_plane_coeffs(1,:)];
	xtilde	= M1 \ ij10;
	re_pos1	= bsxfun(@rdivide, xtilde(1:3,:), xtilde(4,:));

	%% same thing in plane number 2
	M2			= [P; re_plane_coeffs(2,:)];
	xtilde	= M2 \ ij10;
	re_pos2	= bsxfun(@rdivide, xtilde(1:3,:), xtilde(4,:));
	
	%% transpose
	re_pos1	= re_pos1.';
	re_pos2	= re_pos2.';
end






