function [Fsmin,xb,yb,Rb]=ruidian(b,h,gama,phi0,c)
    %输入   b:坡度 h:坡高 gama:土的重度 phi0:土的内摩擦角 c:土的粘聚力
    %输出   Fsmin:最小安全系数 xb,yb:圆心坐标 Rb:滑弧半径
    pi=3.1416926;
    Fsmin=100.0;
    sum1=0;sum2=0;
    phi0=phi0*pi/180.0;
    alpha=atan(b);
    L=h/sin(alpha);
    m=L*cos(alpha);
    x0=m/2;
    y0=h/2;
    for theta=0:alpha/100:alpha
        for d=(0.25*L):(L/1000.0):(1.25*L)
            sum1=0;sum2=0;
            x=x0-cos(pi/2-alpha+theta)*d;
            y=y0+sin(pi/2-alpha+theta)*d;
            R=sqrt(x^2+y^2);
            x1=sqrt(R^2-(h-y)^2)+x;
            d1=x1/1000;
            for xd=0:d1:x1
                yd=y-sqrt(R^2-(xd-x)^2);
                beta=atan((xd-x)/(y-yd));
                n=d1/cos(beta);
                if xd<=m
                    y2=tan(alpha)*xd;
                    h1=abs(y2-yd);
                    W1=gama*h1*d1;
                    sum1=sum1+W1*cos(beta)*tan(phi0)+c*n;
                    sum2=sum2+W1*sin(beta);
                else
                    h1=abs(h-yd);
                    W1=gama*h1*d1;
                    sum1=sum1+W1*cos(beta)*tan(phi0)+c*n;
                    sum2=sum2+W1*sin(beta);
                end
            end
            Fs=sum1/sum2;
            if Fs < Fsmin

                Fsmin=Fs;
                xb=x;
                yb=y;
                Rb=R;
            end
        end
    end
end