各向异性介质多参数化走时反演敏感核（辐射模式）分析

原文参考：

(PDF) Sensitivity kernel analysis and inversion strategy for multi-parameter traveltime inversion in VTI media VTI介质初至波多参数走时反演敏感核分析及反演策略 (researchgate.net)

内容简介：

在阅读文章的基础上，很多朋友对此研究很感兴趣，想进一步深入学习，在此，提供了敏感核分析的代码，下面是Thomsen参数化的代码，其他参数化的代码一并打包上传。（好吧，发现不能打包上传，谁有办法，评论区告诉我！（需要代码的私信我就好））

思考：这种走时的敏感核与辐射模式的区别和联系是什么？（答案在文末）

说明：做该项研究当时用好几种办法推导过,结果都很相似,验证了我的推导是对的。

1.基于Fomel群速度近似的推导(推导较为繁琐，但是精度高)

clear all;

delta=0.1 ;
eps=0.2;
eta=(eps-delta)/(1+2*delta);
v=3000;
s0=1/v;
sh=1/v/sqrt(1+2*eps);
sn=1/v/sqrt(1+2*delta);
fi=0:1:360;
fi=fi/180*pi;

Ks0=zeros(1,361);
Ks01=zeros(1,361);
Ksn=zeros(1,361);
Ksn1=zeros(1,361);
Ksh=zeros(1,361);
Ksh1=zeros(1,361);
Keps=zeros(1,361);
Keps1=zeros(1,361);
Keta=zeros(1,361);
Keta1=zeros(1,361);
Kdel=zeros(1,361);
Kdel1=zeros(1,361);

Q=1+2*(eps-delta)/(1+2*delta);
A=s0*s0/(1+2*eps);
C=s0*s0;
E=A.*sin(fi).^2+C.*cos(fi).^2;
D=sqrt(E.*E+4.*(Q*Q-1)*A*C.*sin(fi).^2.*cos(fi).^2);
Sg=sqrt((1+2*Q)/2/(1+Q).*E+1/2/(1+Q).*D);

Ks0=1/2./Sg.*(E.*0+(1+2*Q)/2/(1+Q)*(2*s0/(1+2*eps).*sin(fi).^2+2*s0.*cos(fi).^2)...
    +D.*0+1/2/(1+Q)*1/2./D.*(2.*E.*(2*s0/(1+2*eps).*sin(fi).^2+2*s0.*cos(fi).^2)...
    +8*Q*0+4*(Q*Q-1)*2*s0/(1+2*eps)*C.*sin(fi).^2.*cos(fi).^2+4*(Q*Q-1)*A*2*s0...
    .*sin(fi).^2.*cos(fi).^2));
Keps=1/2./Sg.*(E.*(1/(1+Q)*2/(1+2*delta)-(1+2*Q)/2/(1+Q)^2*2/(1+2*delta))+(1+2*Q)/2/(1+Q)*...
    (-2*s0*s0/(1+2*eps)^2.*sin(fi).^2)...
    +D./(-2)/(1+Q)^2*2/(1+2*delta)+1/2/(1+Q)*1/2./D.*(2.*E.*(-2*s0*s0/(1+2*eps)^2.*sin(fi).^2)...
    +8*Q*2/(1+2*delta)*A*C.*sin(fi).^2.*cos(fi).^2+4*(Q*Q-1)*(-2)*s0*s0/(1+2*eps)^2*C.*sin(fi).^2.*cos(fi).^2+4*(Q*Q-1)*A*0*s0...
    .*sin(fi).^2.*cos(fi).^2));
Kdel=1/2./Sg.*(E.*(1/(1+Q)*(-2-4*eps)/(1+2*delta)^2-(1+2*Q)*2*(-2-4*eps)/(1+2*delta)^2/4/(1+Q)^2)+(1+2*Q)/2/(1+Q)*0.0...
    +D.*(-1)/2/(1+Q)^2*(-2-4*eps)/(1+2*delta)^2+1/2/(1+Q)*1/2./D.*(2.*E.*0.0...
    +8*Q*(-2-4*eps)/(1+2*delta)^2*A*C.*sin(fi).^2.*cos(fi).^2+4*(Q*Q-1)*0+4*(Q*Q-1)*A*0*s0...
    .*sin(fi).^2.*cos(fi).^2));



figure('units','normalized','position',[0,0,0.4,0.2])%设置成图框的大小
subplot(1,3,1);
polarplot(fi,Ks0);
pax=gca;
pax.ThetaZeroLocation =  'bottom';     % 将0度放在顶部 
title('Ks0');
subplot(1,3,2);
polarplot(fi,Keps);
pax=gca;
pax.ThetaZeroLocation =  'bottom';     % 将0度放在顶部 
title('Kepsilon');
subplot(1,3,3);
polarplot(fi,Kdel);
pax=gca;
pax.ThetaZeroLocation = 'bottom';     % 将0度放在顶部 
title('Kdelta');


% K1=K1-Ks0;
% K2=K2-Keps;
% K3=K3-Kdel;
% subplot(1,3,1);
% polarplot(fi,K1);
% subplot(1,3,2);
% polarplot(fi,K2);
% subplot(1,3,3);
% polarplot(fi,K3);

2.基于Sayers群速度近似的推导(推导相对容易，但是精度低)

% clear all;

delta=0.1;
eps=0.2;
v=3000;
fi=0:1:360;
fi=fi/180*pi;

A=zeros(1,361);

Ks=zeros(1,361);
Keps=zeros(1,361);
Kdel=zeros(1,361);
A=sqrt((1/(1+2*delta)-2/(1+2*eps)+1)*(cos(fi)).^2-(1/(1+2*delta)-1/(1+2*eps))*(cos(fi)).^4+1/(1+2*eps));
Ks=A;
Keps=-1*(sin(fi)).^4./((2*eps+1).^2*A*v);
Kdel=-(cos(fi).*sin(fi)).^2./((2*delta+1)^2.*A*v);
figure('units','normalized','position',[0,0,0.4,0.2])%设置成图框的大小
subplot(1,3,1);
polarplot(fi,Ks);
pax=gca;
pax.ThetaZeroLocation = 'bottom';     % 将0度放在顶部 
title('Ks0');
subplot(1,3,2);
polarplot(fi,Keps);
pax=gca;
pax.ThetaZeroLocation = 'bottom';     % 将0度放在顶部 
title('Kepsilon');
subplot(1,3,3);
polarplot(fi,Kdel);
pax=gca;
pax.ThetaZeroLocation = 'bottom';     % 将0度放在顶部 
title('Kdelta');

该敏感核的物理意义是射线在某特定群角方向透射时，对应参数对走时的影响强弱。

辐射模式是在地下某特定张角下，参数对反射能量影响的强弱。二者本质上其实差别不大。反射也可以理解为两次透射，所以，在VTI介质情况中，敏感核在α角度的强弱趋势与辐射模式在地下2α张角的趋势是一致的。

辐射模式参考文献：

A recipe for practical full-waveform inversion in anisotropic media: An analytical parameter resolution study | Request PDF (researchgate.net)