datacreate.m

这是一个用MATLAB实现的磁法球体相关性自动反演的小程序

 function [X,Y,DeltaT]=datacreate(Xs0,Ys0,Zs0,number,space)
%numbe为采样点的个数
%space为采样的间隔
%%%地磁场参数赋初值%%%
M=1;                                    %磁化强度的值
I=45;                                    %磁倾角
A=45;                                    %磁偏角

%%%球体参数赋初值%%%
%Xs0=0;        
%Ys0=0;                                   %球体的位置
%Zs0=500;
u=4*pi*1.0e-7;                           %磁化率的值
v=100;                                   %球体的体积
m=v*M;                                   %球体的总磁矩


%%%幅度角度的转换%%%
I=I*pi/180;
A=A*pi/180;


%%一维数据采集正演计算%%%
i1=0;
for i=-space*number/2:space:space*number/2
     j1=0;    i1=i1+1;  
   for i=-space*number/2:space:space*number/2
    j1=j1+1;
    Xp=i;
	Yp=i;
        
    x1D(i1,j1)=Yp;                %%%%%注意设定坐标系与画图坐标系列不同
    y1D(i1,j1)=Xp;                %%%%%设定坐标系为正北向为X轴正东向为Y轴
    
	Rx=Xp-Xs0;
	Ry=Yp-Ys0;
    r=Zs0;

    R=sqrt(Rx*Rx+Ry*Ry+r*r);
    
    Hax=(u/4*pi)*m*((2*Rx*Rx-Ry*Ry-r*r)*cos(I)*cos(A)-3*r*Rx*sin(I)+3*Rx*Ry*cos(I)*sin(A))/(R*R*R*R*R);
    Hay=(u/4*pi)*m*((2*Ry*Ry-Rx*Rx-r*r)*cos(I)*cos(A)-3*r*Ry*sin(I)+3*Rx*Ry*cos(I)*sin(A))/(R*R*R*R*R);
    Za=(u/4*pi)*m*((2*r*r-Rx*Rx-Ry*Ry)*sin(I)-3*r*Rx*cos(I)*cos(A)-3*r*Ry*cos(I)*sin(A))/(R*R*R*R*R);

    Delta_T1D(i1,j1)=(Hax*cos(I)*cos(A)+Hay*cos(I)*sin(A)+Za*sin(I))*1.0e9;
   end
end
X=x1D;
DeltaT=Delta_T1D;
Y=y1D;

%%%画剖面图%%%
%xy=(-1000:10:1000);
%figure(1);
%plot(xy,Delta_T1D);