%==================================================================
%
%	change_sphere_cart.m
%
%	Markus Junghoefer	[1995]
%   This software is protected by german copyright and international treaties.             
%   Copyright 2004 Markus Junghöfer & Peter Peyk. All Rights Reserved.                     
%                                                                                          
%   THIS SOFTWARE AND DOCUMENTATION IS PROVIDED "AS IS," AND COPYRIGHT HOLDERS MAKE        
%   NO REPRESENTATIONS OR WARRANTIES, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO,    
%   WARRANTIES OF MERCHANTABILITY OR FITNESS FOR ANY PARTICULAR PURPOSE OR THAT THE USE    
%   OF THE SOFTWARE OR DOCUMENTATION WILL NOT INFRINGE ANY THIRD PARTY PATENTS, COPYRIGHTS,
%   TRADEMARKS OR OTHER RIGHTS. COPYRIGHT HOLDERS WILL NOT BE LIABLE FOR ANY DIRECT,       
%   INDIRECT, SPECIAL OR CONSEQUENTIAL DAMAGES ARISING OUT OF ANY USE OF THE SOFTWARE      
%   OR DOCUMENTATION.                                                                      
%
%
%	If status == 1
%	This will change an (n x (theta,phi)) matrix in an (n x (x,y,z)) matrix 
%	If status == -1
%	This will change an (n x (x,y,z)) matrix in an (n x (theta,phi)) matrix 
%	if status == 2	%Besa to KN

%	Function definition

	function [X_out] = change_sphere_cart(X_in, ro, status)


%==================================================================
if nargin<3; status=[]; end
if nargin<2; ro=[]; end

if isempty(status)
	fprintf('\n');clc;
	fprintf('Change spherical in cartesian coordinates in "change_sphere_cart" !'); 
	pause(1); status=1;
end
if isempty(ro); ro=1; end

[m,n]=size(X_in);

if n==2
	
	X_in_tmp=ones(m,3);
	
	X_in_tmp(:,3)=ro.*ones(m,1);
	
	X_in_tmp(:,1:2)=X_in;
	
	X_in=X_in_tmp;	

end

if status == 1
 	
	%X_out=[x,y,z]
	
	X_out(:,1) = X_in(:,3) .* sin(X_in(:,1)) .*cos(X_in(:,2)); 				%x=r*sin(theta)*cos(phi)

 	X_out(:,2) = X_in(:,3) .* sin(X_in(:,1)) .*sin(X_in(:,2)); 				%y=r*sin(theta)*sin(phi)

	X_out(:,3) = X_in(:,3)	.* cos(X_in(:,1));								%z=r*cos(theta)

elseif status == -1															
	
	%X_out=[theta,phi,ro]
	
	X_out(:,3)= sqrt(X_in(:,1).^2 + X_in(:,2).^2 + X_in(:,3).^2);			%ro=sqrt(x^2+y^2+z^2)
	
	for i=1:m	
		
		X_out(i,1)=acos(X_in(i,3)/X_out(i,3));								%theta=acos(z/r) r darf niemals Null sein!	
		
		if sqrt(X_in(i,1)^2+X_in(i,2)^2)==0									%falls x und y Null sind wird phi zu Null definiert
		
			X_out(i,2)=0;
		
		else
			
			if X_in(i,2)>=0
			
				X_out(i,2)=acos(X_in(i,1)./sqrt(X_in(i,1)^2+X_in(i,2)^2));	%phi=acos(x/sqrt(x^2+y^2))
																			%Trennung in y>0 und y<0
			elseif X_in(i,2)<0
			
				X_out(i,2)=2*pi-acos(X_in(i,1)./sqrt(X_in(i,1)^2+X_in(i,2)^2));
		
			end
			
		end
	
	end

elseif status==2	%Besa to KN

	
	X_out=X_in;
	for i=1:m	
		if X_in(i,1)<0 & X_in(i,2)<0		%front left
			X_out(i,1)=-X_in(i,1);
			X_out(i,2)=180+X_in(i,2);
		elseif X_in(i,1)<0 & X_in(i,2)>=0	%back left
			X_out(i,1)=-X_in(i,1);
			X_out(i,2)=180+X_in(i,2);
		elseif X_in(i,1)>=0 & X_in(i,2)<0	%back right
			X_out(i,2)=360+X_in(i,2);
		end
	end
	X_out(:,1:2)=X_out(:,1:2).*pi./180;
	if ro~=1
		X_out(:,3)=ones(m,1).*ro;
	end
end

return;
%==================================================================
