直流电法正演模拟.doc

function T = Transfer_Fun(m)global rho;global h;N = size(rho,2);T = rho(N);for i = N-1:-1:1 A = 1-exp(-2*m*h(i); B = 1+exp(-2*m*h(i); T = rho(i)*(rho(i)*A+T*B)/(rho(i)*B+T*A);endT = T*m;function z1 = hankel1_47(fun,B)a = -3.05078187595e+00;s = 1.10599010095e-01;wt1 = 3-06 -9.73811660718e-06 1.64866227408e-05 -1.81501261160e-05 1.87556556369e-05 -1.46550406038e-05 1.53799733803e-05 -6.95628273934e-06 1.41881555665e-05 3.41445665537e-06 2-05 2.34962369042e-05 4.84340283290e-05 7.33732978590e-05 1.27703784430e-04 2.08120025730e-04 3.49803898913e-04 5.79107814687e-04 9.65887918451e-04 1.60401273703e-03 2.66903777685e-03 4.43111590040e-03 7.35631696247e-03 1.21782796293e-02 2.01097829218e-02 3.30096953061e-02 5.37143591532e-02 8.60516613299e-02 1.34267607144e-01 2.00125033067e-01 2.74027505792e-01 3-01 2.41655667461e-01 -5.40549161658e-02 -4.46912952135e-01 -1.92231885629e-01 5.52376753950e-01 -3.57429049025e-01 1.41510519002e-01 -4.61421935309e-02 1.48273761923e-02 -5.07479209193e-03 1.83829713749e-03 -6.67742804324e-04 2.21277518118e-04 -5.66248732755e-05 7.88229202853e-06;z1 = 0;for i = 1:size(wt1,1) lamda = (1/B)*10(a+(i-1)*s); z1 = z1+feval(fun,lamda)*wt1(i);endz1 = z1/B;function rho_s = DC1D_Forward(rho1,h1)global rho;global h;rho = rho1;h = h1;r = logspace(1,4,40);for i = 1:si