function robot = get_regressor(robot, opt)
% Create symbolic generilized coordiates, their first and second deriatives
ndof = robot.ndof;
q_sym = sym('q%d',[ndof+1,1],'real');
qd_sym = sym('qd%d',[ndof+1,1],'real');
q2d_sym = sym('q2d%d',[ndof+1,1],'real');
% init regressor
robot.regressor.m = sym('m%d',[ndof,1],'real');
robot.regressor.com = sym('com%d',[ndof,1],'real');
robot.regressor.I = sym('I%d',[ndof,1],'real');
robot.regressor.I_vec = inertiaMatrix2Vector(robot.regressor.I);
robot.regressor.mc = robot.regressor.m.*robot.regressor.com;
robot.regressor.pi = [robot.I_vec(:,i); robot.mc(:,i); robot.m(i)];
% init matrix
R = robot.R;
P = robot.t;
w = robot.vel.w ;
dw = robot.vel.dw ;
dv = robot.vel.dv ;
switch opt.LD_method
    case 'Direct'
        switch opt.KM_method
            case 'MDH'
                for i = 2:ndof+1
                    p_skew(:,:,i) = vec2skewSymMat(P(:,:,i));
                    w_skew(:,:,i) = vec2skewSymMat(w(:,i));
                    dw_skew(:,:,i) = vec2skewSymMat(dw(:,i));
                    dv_skew(:,:,i) = vec2skewSymMat(dv(:,i));
                    w_l(:,:,i) = vec2linearSymMat(w(:,i));
                    dw_l(:,:,i) = vec2linearSymMat(dw(:,i));
                    % size of matrix A is 6*10, need to -1
                    robot.regressor.A(:,:,i-1) = [dv(:,i),dw_skew(:,:,i)+w_skew(:,:,i)*w_skew(:,:,i),zeros(3,6); ...
                        zeros(3,1),-dv_skew(:,:,i),dw_l(:,:,i)+w_skew(:,:,i)*w_l(:,:,i)];
                end
                % construct matrix U, size: [6*n,10*n]
                % U_ = sym(zeros([6*ndof,10*ndof]));
                U_ = [];
                for i = 1:ndof
                    % tricky
                    for j = i:ndof
                        if(j == i)
                            TT = eye(6,6);
                            U_row = TT*robot.regressor.A(:,:,j);
                        else
                            TT = TT*Adjoint(RpToTrans(R(:,:,j),P(:,:,j)));
                            U_row = [U_row,TT*robot.regressor.A(:,:,j)];
                        end
                    end
                    U_ = [U_;zeros(6,(i-1)*10),U_row];
                end
                robot.regressor.U = U_;
                if(opt.debug)
                    sprintf('size of U_=%dx%d.',size(robot.regressor.U));
                end
                robot.regressor.K = [zeros(1,3),Z0]*;
        end
%     matlabFunction(Y_f,'File','autogen/standard_regressor_Two_bar',...
%                    'Vars',{q_sym,qd_sym,q2d_sym});
    case 'Lagrange'
        disp('TODO opt.LD_method Lagrange!')
        return;
    otherwise
        disp('Bad opt.KM_method!')
        return;
end