I am doing ray tracing and I have 5th order polynomial. So I am solving this one as follows but the answer at the end "tCounter" are quite different from the truth (I have the data as reference and it is different totally). But I couldn't debug it. Could any one please help me?
The code is :
Poly = - (7240313003093287*X1^5)/590295810358705651712 + (6570582540807497*X1^4*X2)/576460752303423488 + (3195138777360751*X1^4)/576460752303423488 - (468053640978231*X1^3*X2^2)/36028797018963968 - (1058278253154999*X1^3*X2)/36028797018963968 - (8771588561102997*X1^3)/576460752303423488 - (1204292406378591*X1^2*X2^3)/36028797018963968 - (3734476854713205*X1^2*X2^2)/144115188075855872 + (8064171727270943*X1^2*X2)/576460752303423488 + (4358732576658615*X1^2)/288230376151711744 + (2443031512857249*X1*X2^4)/36028797018963968 + (2674690140041215*X1*X2^3)/36028797018963968 + (4992686117793699*X1*X2^2)/72057594037927936 + (4059700549521059*X1*X2)/288230376151711744 + (4188007860677011*X1)/36893488147419103232 - (7532606855015249*X2^5)/288230376151711744 - (8231536140053275*X2^4)/144115188075855872 - (5002126935193025*X2^3)/72057594037927936 - (8415527271549311*X2^2)/288230376151711744 - (2638814316081383*X2)/4503599627370496 + 4117132021192295/9007199254740992;
syms X1 X2 X3 t Rays1 Surface1 Rays2 Surface2 Surface3 Rays3 z;
Equation = Poly - X3;
tsun = matlabFunction(Equation, 'vars',[X2,X1,X3]);
Equation = tsun((Surface1+t*Rays1),(Surface2+t*Rays2),(Surface3+t*Rays3));
Answer = solve(Equation,t); % The answer is RootOf(.....,z) So I have to prepare the equation and get the roots of it in respect to z
a = char(Answer);
R = strrep(a,'RootOf(','');
R1 = strrep(R,', z)','');
b = sym(R1);
PolyCoeffs = coeffs(b,z);
tfun = matlabFunction(PolyCoeffs, 'vars',[Surface1,Surface2,Surface3,Rays1,Rays2,Rays3]);
tCounter = zeros(length(Rays),1);
NaNIndices = find(isnan(Surface(:,1))==1); % My data contains NaN and I am taking them out
tCounter(NaNIndices) = NaN;
NotNaNIndices = find(isnan(Surface(:,1))==0);
for i = NotNaNIndices'
Surface1New = Surface(i,1); %Surface is Nx3 matrix contains the data that I have to calculate the t according to it
Surface2New = Surface(i,2);
Surface3New = Surface(i,3);
Rays1New = Rays(i,1); %Rays is Nx3 matrix contains the data that I have to calculate the t according to it
Rays2New = Rays(i,2);
Rays3New = Rays(i,3);
P = tfun(Surface1New,Surface2New,Surface3New,Rays1New,Rays2New,Rays3New);
t = roots(P);
t(imag(t) ~= 0) = [];
t(t<0) = [];
t = min(t);
tCounter(i) = t;
end
Many thanks in advance
I think your mistake is here:
PolyCoeffs = coeffs(b,z);
PolyCoeffs = fliplr(PolyCoeffs);
You have to flip it because the function Coeffs give coefficient in the opposite order for example:
syms x
c = coeffs(16*x^2 + 19*x + 11)
c =
[ 11, 19, 16]
and function roots takes the other order for example:
The polynomial s3 – 6s2 – 72s – 27 is represented in MATLAB® software as
p = [1 -6 -72 -27]
The roots of this polynomial are returned in a column vector by
r = roots(p)
Good luck in your work.