banjinyz.m

利用MATLAB工具箱

function [Di]=banjin(Rm)
Rm=1.222225;% 泥浆电阻率值
A=importdata('xing10-013.txt');%根据需要修改输入数据文件名
hh=1;%深度数据在数据文件xing10-013.txt（即矩阵A）中的列数
h=A(:,hh);
ss=3;%浅视电阻率在数据文件（即矩阵A）中的列数
Rlls=A(:,ss);
dd=2;%深视电阻率在数据文件xing10-013.txt（即矩阵A）中的列数
Rlld=A(:,dd);
oo=4;%冲洗带电阻率在数据文件xing10-013.txt（即矩阵A）中的列数
Rxo=A(:,oo);
DD=7;%井径在数据文件xing10-013.txt（即矩阵A）中的列数
Dh=A(:,DD);
Rxo=Rxo';
Rlls=Rlls';
Rlld=Rlld';
Sd=load('jiaozhengD.txt');%读入深视电阻率校正曲线函数的系数值
Ss=load('jiaozhengS.txt');%读入深视电阻率校正曲线函数的系数值
S=load('banjinXX.txt');%读入侵入半径校正曲线函数的系数值
x=log10(Rlld/Rm);
x0=log10(Rlls/Rm);
for t=1:9434%在数据文件xing10-013.xls中数据总行数
k=[6 8 10 12 14];
m=k-Dh(t);
if find(m==0)~=zeros(1,0)%当井径为[6 8 10 12 14]中的某个值时
    n=find(m==0);
    p1=Sd(:,n);
    p2=Ss(:,n);
    Rlldcc(t)=polyval(p1,x(t))*Rlld(t);%计算得到深视电阻率校正后井径8in下的深视电阻率
    Rllscc(t)=polyval(p2,x(t))*Rlls(t);%计算得到浅视电阻率校正后井径8in下的浅视电阻率 
else%当井径不为[6 8 10 12 14]中的任何值
  
    r1=find(m<0);
    [n11 n0]=size(r1);
   switch n0
       case{0}%当井径小于6in时
           n1=n0+1;
           p11=Sd(:,n1);
    p12=Sd(:,n1+1);
    y11=polyval(p11,x(t));
    y12=polyval(p12,x(t));
    y1=[y11 y12];
    x1=[k(n1) k(n1+1)];
    Rlldc(t)=interp1(x1,y1,Dh(t),'linear','extrap');
    Rlldcc(t)=Rlldc(t)*Rlld(t); %利用插值方法计算得到深视电阻率校正后井径8in下的深视电阻率
   
   p21=Ss(:,n1);
    p22=Ss(:,n1+1);
    y21=polyval(p21,x0(t));
    y22=polyval(p22,x0(t));
    y2=[y21 y22];
    Rllsc(t)=interp1(x1,y2,Dh(t),'linear','extrap');
    Rllscc(t)=Rllsc(t)*Rlls(t);%利用插值方法计算得到浅视电阻率校正后井径8in下的浅视电阻率
       case{5}%当井径大于14in时
           n1=n0-1;
          p11=Sd(:,n1);
    p12=Sd(:,n1+1);
    y11=polyval(p11,x(t));
    y12=polyval(p12,x(t));
    y1=[y11 y12];
    x1=[k(n1) k(n1+1)];
    Rlldc(t)=interp1(x1,y1,Dh(t),'linear','extrap');
    Rlldcc(t)=Rlldc(t)*Rlld(t);%利用插值方法计算得到深视电阻率校正后井径8in下的深视电阻率
    
   p21=Ss(:,n1);
    p22=Ss(:,n1+1);
    y21=polyval(p21,x0(t));
    y22=polyval(p22,x0(t));
    y2=[y21 y22];
    Rllsc(t)=interp1(x1,y2,Dh(t),'linear','extrap');
    Rllscc(t)=Rllsc(t)*Rlls(t);%利用插值方法计算得到浅视电阻率校正后井径8in下的浅视电阻  
       otherwise%当井径在[6 8 10 12 14]之间但不等于6或8或10或12或14
           n1=n0;
            p11=Sd(:,n1);
    p12=Sd(:,n1+1);
    y11=polyval(p11,x(t));
    y12=polyval(p12,x(t));
    y1=[y11 y12];
    x1=[k(n1) k(n1+1)];
    Rlldc(t)=interp1(x1,y1,Dh(t));
    Rlldcc(t)=Rlldc(t)*Rlld(t);%利用插值方法计算得到深视电阻率校正后井径8in下的深视电阻率
    
   p21=Ss(:,n1);
    p22=Ss(:,n1+1);
    y21=polyval(p21,x0(t));
    y22=polyval(p22,x0(t));
    y2=[y21 y22];
    Rllsc(t)=interp1(x1,y2,Dh(t));
    Rllscc(t)=Rllsc(t)*Rlls(t);%利用插值方法计算得到浅视电阻率校正后井径8in下的浅视电阻
   end
   
end
end
Xj=log10(Rlldcc./Rllscc);%转化为侵入半径校正曲线中的横坐标值
Yj=log10(Rlldcc./Rxo);%转化为侵入半径校正曲线中的纵坐标值
b=length(Xj);
for j=1:b
for i=1:8
    p=S(:,i)';
    Y(j,i)=polyval(p,Xj(j));
    p=zeros(1,3);
end
end
X=[20 30 40 50 60 80 100 120];
for i=1:b
    Di(i)=interp1(Y(i,:),X,Yj(i),'linear','extrap');%插值计算对应深度下的泥浆侵入半径

end

Di=Di*25.4/1000;%将存储的侵入半径从in转化为以m为单位
disp(Di);
h=h';

w=[h;Di];
fid=fopen('jieguo.xls','w');%将侵入半径数据存储到jieguo.xls文件中
fprintf(fid,'%12.5f %12.5f\n',w);
fclose(fid);