function coef=polyFitChebyshev(y, order)
% polyFitChebyshev returns the fitting coefficient of Chebyshev polynomial, starting from T_0(x)

%% Create Chebyshev polynomial, where each row is a polynomial
chebyPoly=zeros(order+1, order+1);
chebyPoly(1,end)=1;
chebyPoly(2,end-1:end)=[1 0];
for i=3:order+1
	chebyPoly(i, :)=2*[chebyPoly(i-1,2:end), 0]-chebyPoly(i-2,:);
end

%% Create A in Ax=b
pointCount=length(y);
x=linspace(-1, 1, pointCount)';
A=zeros(pointCount, order+1);
for i=1:order+1
	A(:,i)=polyval(chebyPoly(i,:), x);
end

coef=A\y(:);