[a,eckart_rot_matrix]=eckart(mass_weight_matrix,intervals)

function [a,eckart_rot_matrix]=eckart(mass_weight_matrix,intervals)

N=length(mass_weight_matrix(:,1,1));

trials=length(mass_weight_matrix(1,1,:));

mass_weight_comp=permute(mass_weight_matrix,[3 2 1]);

reshape_mass_weight_comp=reshape(mass_weight_comp,[trials,3*N]);

fft_mass_weight_comp=interpft(reshape_mass_weight_comp,trials*intervals);

reshape_mass_weight_matrix=reshape(fft_mass_weight_comp,[trials*intervals,3,N]);

mass_weight_matrix=permute(reshape_mass_weight_matrix,[3 2 1]);

[aligned_axes,a(1:N,1:3,1)]=align_inertia(mass_weight_matrix(1:N,1:3,1));

guess_lambda=aligned_axes';

eckart_rot_matrix(1:3,1:3,1)=guess_lambda;

for i=2:length(mass_weight_matrix(1,1,:))

    guess_matrix=mass_weight_matrix(1:N,1:3,i)*(guess_lambda');

    x=nleqn_jacob(@eckart_equations,[0;0;0],eps,eps,eps,guess_matrix,a(1:N,1:3,i-1));

    theta=x(1);
    phi=x(2);
    chi=x(3);

    lambda=direction_cosine(theta,phi,chi);

    guess_lambda=lambda*guess_lambda;

    eckart_rot_matrix(1:3,1:3,i)=guess_lambda;

    a(1:N,1:3,i)=guess_matrix*(lambda');

end

a=a(:,:,1:intervals:end);

eckart_rot_matrix=eckart_rot_matrix(:,:,1:intervals:end);