
% ------------------------------------------------------------------------
% Function: gammaFit1d(x,y)
% ------------------------------------------------------------------------
% Aim:
% To fit the GSFE data
% ------------------------------------------------------------------------

function [fitresult, gof] = gammaFit1d(x, y)

    global fitcut

    if fitcut>=2

        R_0='';
        I_x='';
        fitfunc='';
        mmnn=1;

        for m=01:20

            if ( m <= fitcut )

                cosfit=subs(strrep('R_mn*cos(mm*2*pi*Ux)','mn',num2str(mmnn)),'mm',m);
                sinfit=subs(strrep('I_mn*sin(mm*2*pi*Ux)','mn',num2str(mmnn)),'mm',m);

                R_0=strcat(R_0,strrep('+R_mn','mn',num2str(mmnn)));

                if m~=1
                    I_x=strcat(I_x,strrep('+mn*I_mn','mn',num2str(mmnn)));
                end

                fitfunc=strcat(fitfunc,'+',char(cosfit),'+',char(sinfit));

                mmnn=mmnn+1;

            end

        end

        R_0=strcat('-(',R_0,')');
        I_x=strcat('-(',I_x,')','*');
        fitfunc=strcat(R_0,fitfunc);
        fitfunc=strrep(fitfunc,'+I_1*',I_x);

	elseif fitcut==1

        fitfunc='-R_1+R_1*cos(Ux*2*pi)';

	elseif fitcut==0

        fitfunc='-(R_1+R_2+R_3+R_4+R_5)+R_1*cos(Ux*2*pi)-(2*I_2+3*I_3+4*I_4+5*I_5)*sin(Ux*2*pi)+R_2*cos(2*Ux*2*pi)+I_2*sin(2*Ux*2*pi)+R_3*cos(3*Ux*2*pi)+I_3*sin(3*Ux*2*pi)+R_4*cos(4*Ux*2*pi)+I_4*sin(4*Ux*2*pi)+R_5*cos(5*Ux*2*pi)+I_5*sin(5*Ux*2*pi)'; 

    end

    %% Fit: 'Peierls-Nabarro model'.
    [xData, yData] = prepareCurveData( x, y );

    % Set up fittype and options.
    ft = fittype( fitfunc, 'independent', 'Ux', 'dependent', 'SFE' );
    opts = fitoptions( ft );
    opts.Algorithm = 'Levenberg-Marquardt';
    opts.DiffMinChange = 1e-20;
    opts.Display = 'Off';
    opts.Robust = 'LAR';
    opts.TolFun = 1e-12;
    opts.TolX = 1e-12;

    if fitcut>=2
        opts.StartPoint = rand(1,2*mmnn-3);
    elseif fitcut==1
        opts.StartPoint = rand(1,1);
    elseif fitcut==0
        opts.StartPoint = rand(1,9);
    end    

    % Fit model to data.
    [fitresult, gof] = fit( xData, yData, ft, opts );

    % if gof.rsquare<=0.9
    %     
    %     
    % end

    %% Plot fit with data.
    figure( 'Name', 'The resutl of fitting gamma surface' );
    h = plot( fitresult, xData, yData );
    legend( h, 'Input SFE data', 'SFE vs. Ux', 'Location', 'NorthEast' );
    % Title and Label axes
    title('The resutl of fitting gamma surface');
    xlabel('Displacement u/b');
    ylabel('SFE (J/m2)');
    grid on
end
