function [Sx1,Sy1]=regression(Rx,Ry,p,Tn,k)
% input:
%     p:地震初至时间原始数据
%     Rx,Ry:检波器坐标组
%     Tn:道数
%     k:道集序号
%     email:2906310084@qq.com
%     2022-7-16
    peak=p((k-1)*Tn+1:k*Tn,:);
    t2=(peak(:,2)./10^6).^2;
    x0=ones(1,105);
    x1=Rx(k,:).^2+Ry(1,:).^2;
    x2=Rx(k,:);
    x3=Ry(k,:);
    b=[x0',x1',x2',x3'];%建立模型
    a=b\t2;
    a=a./a(2);
    Sx1=a(3)./(-2);
    Sy1=a(4)./(-2);
    scatter(Sx1,Sy1,'yO')%标记校正后炮点位置
    legend('','检波点','校正前炮点','校正后炮点')
end