logplus/logPlus/DepPairs.cpp

1236 lines
33 KiB
C++
Raw Normal View History

//////////////////////////////////////////////////////////////////
///
/// 测录井信息综合平台
///
/// 编制:王卫
///
/// 版权:中国石化集团石油工程技术研究院
///
/// 2010年10月-
///
// Depth.cpp : implementation file
//
#include "DepPairs.h"
#include "MemRdWt.h"
#include "math.h"
#include "BaseFun.h"
#include "geometryutils.h"
#include <qdir.h>
#include "qapplication.h"
#include <QMessageBox.h>
int GetDepPairs(float* curve1,float stdep1,float endep1,float rlev1,
float* curve2,float stdep2,float endep2,float rlev2,float cmins,
float winLength,float stepLength, float searchLength,
float *FirstDep,
float *SecondDep,
float *corc)
{
DepPairs depPairs(curve1,stdep1,endep1,rlev1,
curve2,stdep2,endep2,rlev2,cmins,
winLength,stepLength,searchLength);
depPairs.FirstDep=FirstDep;
depPairs.SecondDep=SecondDep;
depPairs.corc=corc;
float stdep=FirstDep[0];
float endep=SecondDep[0];
if(stdep==endep){
stdep=stdep1;
endep=endep1;
stdep=stdep1;
if(stdep>stdep2) stdep=stdep2;
endep=endep1;
if(endep<endep2) endep=endep2;
}
depPairs.StartDepth=stdep;
depPairs.EndDepth=endep;
return depPairs.AdmRun(curve1,curve2,NULL);
}
int getDepPairs(float* curve1,float stdep1,float endep1,float rlev1,
float* curve2,float stdep2,float endep2,float rlev2,
float *FirstDep,
float *SecondDep)
{
DepPairs depPairs(curve1,stdep1,endep1,rlev1,
curve2,stdep2,endep2,rlev2,NULL,
0,0,0);
depPairs.FirstDep=FirstDep;
depPairs.SecondDep=SecondDep;
depPairs.corc=NULL;
float stdep=FirstDep[0];
float endep=SecondDep[0];
if(stdep==endep){
stdep=stdep1;
endep=endep1;
stdep=stdep1;
if(stdep>stdep2) stdep=stdep2;
endep=endep1;
if(endep<endep2) endep=endep2;
}
depPairs.StartDepth=stdep;
depPairs.EndDepth=endep;
depPairs.ipar=1;
return depPairs.AdmRun(curve1,curve2,NULL);
}
int AutoCompZone(char *Filename,char *Filename1,char *CurveName1,char *CurveName2,
float *FirstDep,
float *SecondDep)
{
CMemRdWt well;
if(!well.Open(Filename,CMemRdWt::modeRead)) {
return 0;
};
CMemRdWt well1;
if(!well1.Open(Filename1,CMemRdWt::modeRead)) {
return 0;
}
int index=well.OpenCurve(CurveName1);
int index1=well1.OpenCurve(CurveName2);
if(index>-1&&index1>-1) {
Slf_CURVE curve,curve1;
well.GetCurveInfo(index,&curve);
int count=(curve.EndDepth-curve.StartDepth)/curve.DepLevel+1.5;
float *buffer=new float[count+1];
well.ReadCurveToFloatBuf(index,curve.StartDepth,count,buffer);
well1.GetCurveInfo(index1,&curve1);
int count1=(curve1.EndDepth-curve1.StartDepth)/curve1.DepLevel+1.5;
float *buffer1=new float[count1+1];
well1.ReadCurveToFloatBuf(index1,curve1.StartDepth,count1,buffer1);
int num=getDepPairs(buffer,curve.StartDepth,curve.EndDepth,curve.DepLevel,
buffer1,curve1.StartDepth,curve1.EndDepth,curve1.DepLevel,
FirstDep,
SecondDep
);
delete buffer;
delete buffer1;
return num;
}
return 0;
}
int AutoComp(char *Filename,char *Filename1,char *CurveName1,char *CurveName2,
float cmins,
float winLength,float stepLength, float searchLength,
float *FirstDep,
float *SecondDep,
float *corc)
{
CMemRdWt well;
if(!well.Open(Filename,CMemRdWt::modeRead)) {
return 0;
};
CMemRdWt well1;
if(!well1.Open(Filename1,CMemRdWt::modeRead)) {
return 0;
}
int index=well.OpenCurve(CurveName1);
int index1=well1.OpenCurve(CurveName2);
if(index>-1&&index1>-1) {
Slf_CURVE curve,curve1;
well.GetCurveInfo(index,&curve);
int count=(curve.EndDepth-curve.StartDepth)/curve.DepLevel+1.5;
float *buffer=new float[count+1];
well.ReadCurveToFloatBuf(index,curve.StartDepth,count,buffer);
well1.GetCurveInfo(index1,&curve1);
int count1=(curve1.EndDepth-curve1.StartDepth)/curve1.DepLevel+1.5;
float *buffer1=new float[count1+1];
well1.ReadCurveToFloatBuf(index1,curve1.StartDepth,count1,buffer1);
int num=GetDepPairs(buffer,curve.StartDepth,curve.EndDepth,curve.DepLevel,
buffer1,curve1.StartDepth,curve1.EndDepth,curve1.DepLevel,
cmins,
winLength,
stepLength,
searchLength,
FirstDep,
SecondDep,
corc
);
delete buffer;
delete buffer1;
return num;
}
return 0;
}
bool AutoCompute(char *Filename,char *Filename1,char *CurveName1,char *CurveName2)
{
CMemRdWt well;
if(!well.Open(Filename,CMemRdWt::modeRead)) {
return false;
};
CMemRdWt well1;
if(!well1.Open(Filename1)) {
return false;
}
int iIndex=well.OpenTable("OG_RESULT");
if(iIndex<0) {
return false;
}
int index=well.OpenCurve(CurveName1);
int index1=well1.OpenCurve(CurveName2);
if(index1>-1&&index1>-1) {
Slf_CURVE curve,curve1;
well.GetCurveInfo(index,&curve);
well1.GetCurveInfo(index1,&curve1);
int count=(curve.EndDepth-curve.StartDepth)/curve.DepLevel+1.5;
float *buffer=new float[count+1];
int count1=(curve1.EndDepth-curve1.StartDepth)/curve1.DepLevel+1.5;
float *buffer1=new float[count1+1];
well.ReadCurveToFloatBuf(index,curve.StartDepth,count,buffer);
well1.ReadCurveToFloatBuf(index1,curve1.StartDepth,count1,buffer1);
float cmins=0.9;
float winLength=4;
float stepLength=0.5;
float searchLength=8;
float* FirstDep=new float[5000];
float* SecondDep=new float[5000];
float* corc=new float[5000];
int num=GetDepPairs(buffer,curve.StartDepth,curve.EndDepth,curve.DepLevel,
buffer1,curve1.StartDepth,curve1.EndDepth,curve1.DepLevel,
cmins,
winLength,
stepLength,
searchLength,
FirstDep,
SecondDep,
corc
);
delete buffer;
delete buffer1;
char buf[80];
LAYER_DATA m_Result;
int iIndex1=well1.OpenTable("RESULT");
int len=well.GetTableRecordLength(iIndex);
if (iIndex1 < 0)
{
iIndex1=well1.Open_Set_Table("RESULT",0,35,
"NO,SDEP,EDEP,ZONE,RESULTNO,THICK,TT,MDEPTH1,MDEPTH2,MDEPTH3,MDEPTH4,MDEPTH5,MDEPTH6,MDEPTH7,MDEPTH8,MDEPTH9,MDEPTH10,OIL,POOROIL,OILWATER,WATEROIL,GAS,POORGAS,GASWATER,WATERGAS,DEST1,DEST2,DEST3,DEST4,DEST5,DEST6,DEST7,DEST8,DEST9,DEST10,SDEP1,EDEP1,SDEP2,EDEP2,SDEP3,EDEP3,SDEP4,EDEP4,SDEP5,EDEP5,SDEP6,EDEP6,SDEP7,EDEP7,SDEP8,EDEP8,SDEP9,EDEP9,SDEP10,EDEP10",
"4,4,4,8,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,48,48,48,48,48,48,48,48,48,48,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4",//字段长度
"1,4,4,6,1,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,6,6,6,6,6,6,6,6,6,6,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4",//字段类型
"0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0");//字段备注
}
int iIndex2=well.OpenTable("RESULT");
if (iIndex2 < 0)
{
iIndex2=well.Open_Set_Table("RESULT",0,35,
"NO,SDEP,EDEP,ZONE,RESULTNO,THICK,TT,MDEPTH1,MDEPTH2,MDEPTH3,MDEPTH4,MDEPTH5,MDEPTH6,MDEPTH7,MDEPTH8,MDEPTH9,MDEPTH10,OIL,POOROIL,OILWATER,WATEROIL,GAS,POORGAS,GASWATER,WATERGAS,DEST1,DEST2,DEST3,DEST4,DEST5,DEST6,DEST7,DEST8,DEST9,DEST10,SDEP1,EDEP1,SDEP2,EDEP2,SDEP3,EDEP3,SDEP4,EDEP4,SDEP5,EDEP5,SDEP6,EDEP6,SDEP7,EDEP7,SDEP8,EDEP8,SDEP9,EDEP9,SDEP10,EDEP10",
"4,4,4,8,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,48,48,48,48,48,48,48,48,48,48,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4",//字段长度
"1,4,4,6,1,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,6,6,6,6,6,6,6,6,6,6,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4",//字段类型
"0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0");//字段备注
}
count=well.GetTableRecordCount(iIndex);
well.SetTableRecordCount(iIndex2,count);
well1.SetTableRecordCount(iIndex1,count);
int j=0;
for(int k=0;k<count;k++) {
well.ReadTable(iIndex,k+1,&buf);
well.WriteTable(iIndex2,k+1,&buf);
memmove(&m_Result,buf,sizeof(LAYER_DATA));
if(m_Result.StartDepth<curve.StartDepth||m_Result.StartDepth>curve.EndDepth) {
continue;
}
if(m_Result.EndDepth<curve.StartDepth||m_Result.EndDepth>curve.EndDepth) {
continue;
}
well1.WriteTable(iIndex1,j+1,&buf);
j++;
}
well1.SetTableRecordCount(iIndex1,j);
well.CloseTable(iIndex);
well1.CloseTable(iIndex1);
well.CloseTable(iIndex2);
char strName[256];
memset(strName, 0, sizeof (strName));
strcat(strName,"OG_RESULT1");
well1.EShiftTableDepth(strName, num,FirstDep,SecondDep);
delete FirstDep;
delete SecondDep;
delete corc;
return true;
}
return false;
}
////////////////////////////////////////////////////////////////////////////
// DepPairs对象
// 显示边框
// 构造函数
DepPairs::DepPairs(float* curve1,float stdep1,float endep1,float rlev1,
float* curve2,float stdep2,float endep2,float rlev2,float cmins,
float winLength,float stepLength, float searchLength)
{
Curve1.stdep=stdep1;
Curve1.endep=endep1;
Curve1.rlev=rlev1;
Curve1.val=curve1;
Curve2.stdep=stdep2;
Curve2.endep=endep2;
Curve2.rlev=rlev2;
Curve2.val=curve2;
// num=i;
DepLevel=rlev1;
if(DepLevel>rlev2) DepLevel=rlev2;
StartDepth=stdep1;
if(StartDepth>stdep2) StartDepth=stdep2;
EndDepth=endep1;
if(EndDepth<endep2) EndDepth=endep2;
// int num=(EndDepth-StartDepth)/DepLevel+1.5;
m_nDepPairsNum=0;
wl=winLength;
sl=searchLength;
step=stepLength;
dr1=3.;
dr2=3.;
cmin=cmins;
elast=1.;
emp=0;
size=1;
type=0;
way=0;
ipar=0;
SetRef=0;
mCurrent=-1;
ShiftNum=0;
}
// 析构函数
DepPairs::~DepPairs()
{
}
// read depth pair from file
void DepPairs::qs(float *datx,float *daty,float *datc,int left,int right)
{
int i,j;
float x,y;
i=left;j=right;
x=datx[(left+right)/2];
do {
while(datx[i]<x&&i<right) i++;
while(datx[j]>x&&j>left ) j--;
if(i<=j) {
y=datx[i];datx[i]=datx[j];datx[j]=y;
y=daty[i];daty[i]=daty[j];daty[j]=y;
y=datc[i];datc[i]=datc[j];datc[j]=y;
i++;j--;
}
}while(i<=j);
if(left<j) qs(datx,daty,datc,left,j);
if(i<right) qs(datx,daty,datc,i,right);
}
void DepPairs::detscl(float sdeps,float edeps)
{
int k;
float tempdep;
tempdep = sdeps;
float xmin1=99999.;
float xmin2=99999.;
float xmax1=-99999.;
float xmax2=-99999.;
k=(int)((sdeps-sdeps)/ DepLevel+0.5);
float val=0;
while(tempdep<edeps)
{
val=Curve1.GetValue(tempdep);
if(val<xmin1) xmin1=val;
if(val>xmax1) xmax1=val;
val=Curve2.GetValue(tempdep);
if(val<xmin2) xmin2=val;
if(val>xmax2) xmax2=val;
k++;
tempdep=sdeps + k * DepLevel;
}
if(xmin1==xmax1)
{
xmin1=0.;
xmax1=1.;
}
if(xmin2==xmax2)
{
xmin2=0.;
xmax2=1.;
}
return;
}
float DepPairs::ampl(int n,float *x)
{
int i;
float xmin,xmax;
xmin=9999.;
xmax=-9999.;
for(i=0;i<n;i++) {
if(x[i]<xmin)xmin=x[i];
if(x[i]>xmax)xmax=x[i];
}
return xmax-xmin;
}
float DepPairs::corr(int n,float *x,float *y,int size,int emp)
{
float sx,sy,xave,yave,sqx,sqy,sxy,dp,dx,dy,c;
float r1,r2,amp1,amp2,ratio;
int i;
sx=0.;
sy=0.;
for(i=0;i<n;i++) {
sx=sx+x[i];
sy=sy+y[i];
}
xave=sx/n;
yave=sy/n;
sqx=0.;
sqy=0.;
for(i=0;i<n;i++) {
sqx=sqx+(x[i]-xave)*(x[i]-xave);
sqy=sqy+(y[i]-yave)*(y[i]-yave);
}
dx=sqrt(sqx/n);
dy=sqrt(sqy/n);
sxy=0.;
for(i=0;i<n;i++) {
sxy=sxy+(x[i]-xave)*(y[i]-yave);
}
dp=sxy/n;
c=0.;
if(dx*dy!=0.)c=dp/(dx*dy);
if(size==0.) {
r1=1.;
r2=1.;
if(dx!=0.)r1=dy/dx;
if(dy!=0.)r2=dx/dy;
if(r1>=r2)ratio=r1;
if(r1<r2)ratio=r2;
c=c/ratio;
}
if(emp!=0) {
amp1=ampl(n,x);
amp2=ampl(n,y);
ratio=amp1;
if(ratio>amp2)ratio=amp2;
c=c*ratio;
}
return c;
}
void DepPairs::cor(int iwl,int isl,float *x,float *y,float *c,int size,int emp)
{
int i,ix=isl,iy;
for(i=0;i<2*isl;i++) {
ix=isl;
iy=i;
c[i]=corr(iwl,&x[ix],&y[iy],size,emp);
}
}
double DepPairs::fun(DPoint &p1, DPoint &p2, DPoint &p3)
{
double s;
s= (p2.x-p1.x)*(p3.y-p1.y)+(p1.y-p2.y)*(p3.x-p1.x);
return s;
}
float DepPairs::act(float *x,int n,int way)
{
int i,i1,i2;
float sum,sq,ave,z,sqd;
sum=0.;
for(i=0;i<n;i++) sum=sum+x[i];
ave=sum/n;
sq=0.;
for(i=0;i<n;i++) sq=sq+(x[i]-ave)*(x[i]-ave);
sqd=sqrt(sq/n);
if(way==0) return sqd;
i1=0;
i2=0;
for(i=0;i<n;i++) {
if(i>=n/2&&x[i]>=ave)i1++;
if(i>=n/2&&x[i]<ave)i2++;
if(i<n/2&&x[i]>=ave)i2++;
if(i<n/2&&x[i]<ave)i1++;
}
z=sqd;
if(way==1&&i2>i1)z=-z;
if(way==2&&i2<=i1)z=-z;
return z;
}
void DepPairs::actvty(int iwl,int isl,float *dr,int way,float *x,float *c)
{
int n,n2,i;
n= (int)(*dr / DepLevel);
n2=n/2;
n=2*n2+1;
for(i=0;i<iwl+2*isl;i++) c[i]=act(&x[i],n,way);
for(i=0;i<iwl+2*isl;i++) x[i]=c[i];
}
void DepPairs::adjust(int n,int way,float *c)
{
int i;
for(i=0;i<n;i++) {
if(way==0)c[i]=abs(c[i]);
if(way==2)c[i]=-c[i];
}
}
int DepPairs::MMax(float depth,int isl,float *x,float *y,float *c,float cmins)
{
int i,iup,idown,n,idir1=0,idir2=0;
float xdep,ydep;
n=2*isl+1;
iup=1;
idown=-1;
if( c[0] < c[1] ) idir1 = iup;
if( c[0] >= c[1] ) idir1 = idown;
for(i=1;i<n-1;i++)
{
if( c[i] < c[i+1] )idir2=iup;
if( c[i] >= c[i+1] )idir2=idown;
if( idir1 == idown && idir2 == iup ) idir1=idir2;
if( idir1 == iup && idir2 == idown )
{
xdep=depth;
ydep=depth+(i-isl)*DepLevel;
if ( c[i] >= 0.5) {
fprintf(fp_temp,"%8.3f %8.3f %5.3f\n",ydep,xdep,c[i]);
// break;
}
idir1=idir2;
}
}
return 1;
}
int DepPairs::MMax(float ydep,float dep0,int isl,float *c,float cmins)
{
float xdep;
int ss=-1;
int n=1.5*isl+1;
float maxcor=-99999;
for(int i=0;i<n;i++)
{
if( c[i] >maxcor) {
maxcor=c[i];
ss=i;
}
}
if(ss>-1&&c[ss]>cmins)
{
xdep=dep0+ss*DepLevel;
fprintf(fp_temp,"%8.3f %8.3f %5.3f\n",xdep,ydep,c[ss]);
}
return 1;
}
int DepPairs::signum(float x)
{
int y;
y=0;
if(x<0.)y=-1;
if(x>0.)y=1;
return y;
}
int DepPairs::optim(int n,float *x1,float *y1,float *c1,float cmins,float elast)
{
float xmax=999999.;
int xnull=0;
float tfine=0,check,diff=0,fcrif=0,fine=0,fmin=0;
int i,iw,j,ip1,ip2,imin,nd,ids,nd1;
n=1;
float a,b,d;
int aa=0;
rewind(fp_temp);
while(!feof(fp_temp))
{
fscanf(fp_temp,"%f%f%f\n",&a,&b,&d);
if(ipar) aa++;
else if(d>=cmins) aa++;
}
if(aa<=0) return -1;
float *x=new float[aa+1000];
float *y=new float[aa+1000];
float *c=new float[aa+1000];
float *f=new float[aa+1000];
int *ncs=new int[aa+1000];
rewind(fp_temp);
while(!feof(fp_temp))
{
if(fscanf(fp_temp,"%f%f%f\n",&x[n],&y[n],&c[n])!=3)
{
break;
};
if(ipar) n++;
else if(c[n]>=cmins) n++;
}
fclose(fp_temp);
QString dir=GetLogdataPath();
dir+="\\temp";
CreateDir((char *)dir.toStdString().c_str());
dir+="\\admtemp";
fp_temp=fopen(dir.toStdString().c_str(),"w+");
if(ipar)
{
for(ip1=0;ip1<n;ip1++)
{
ncs[ip1]=0;
}
for(int i=1;i<n;i++)
{
for(j=i-1;j>0;j--)
{
if(x[i]<=x[j]){
ncs[i]=1;
break;
}
}
if(x[i]-x[i-1]>0&&x[i]-x[i-1]<0.5)ncs[i]=1;
if(!ncs[i])
{
for(j=i+1;j<n-1;j++)
{
if(x[i]>=x[j]){
ncs[i]=1;
break;
}
}
}
}
rewind(fp_temp);
j=0;
for(int i=1;i<n;i++){
if(ncs[i]) continue;
fprintf(fp_temp,"%f %f %f\n",x[i],y[i],c[i]);
j++;
}
}
else {
x[0]=xmax; y[0]=xmax; c[0]=0.;
x[n]=xmax; y[n]=xmax; c[n++]=0.;
f[0]=0; ncs[0]=xnull;
fine=elast;
for(ip1=1;ip1<n;ip1++)
{
f[ip1]=0;
ncs[ip1]=0;
}
for(ip1=1;ip1<n;ip1++)
{
if(x[ip1]>=xmax&&ip1!=n-1) goto l2100;
fmin=xmax;
tfine=0.;
for(iw=0;iw<ip1;iw++)
{
ip2=ip1-iw-1;
if(iw!=0&&x[ip2+1]<xmax) tfine=tfine+fine*c[ip2+1];
if(x[ip2]>=xmax&&ip2!=0) goto l1900;
if(tfine>fmin) goto l2000;
if(x[ip1]<xmax&&x[ip2]<xmax) goto l1500;
fcrif=f[ip2]+tfine;
goto l1800;
l1500:;
check=signum(x[ip1]-x[ip2])*signum(y[ip1]-y[ip2]);
if(check<=0.0) goto l1900;
diff=abs(( x[ip1]-x[ip2]) - (y[ip1]-y[ip2]));
fcrif=f[ip2]+tfine+diff;
l1800:;
if(fcrif>fmin) goto l1900;
fmin=fcrif;
imin=ip2;
l1900:;
}
l2000:;
f[ip1]=fmin;
ncs[ip1]=imin;
l2100:;
}
nd=0;
ids=n-1;
while(1)
{
ids=ncs[ids];
if(ids==0) break;
f[nd++]=ids;
}
for(i=0;i<nd;i++)
{
j=nd-i-1;
ncs[i]=f[j];
}
for(nd1=0,i=0;i<nd;i++)
{
if(x[ncs[i]]<xmax)
{
x[nd1]=x[ncs[i]];
y[nd1]=y[ncs[i]];
c[nd1++]=c[ncs[i]];
}
}
rewind(fp_temp);
fprintf(fp_temp,"%f %f %f\n",x[0],y[0],c[0]);
j=0;
float d1=0;
j++;
for(int i=1;i<nd1;i++) {
d=x[i]-x[i-1];
d1=y[j]-y[i-1];
if(d!=d1||i==nd1-1)
{
fprintf(fp_temp,"%f %f %f\n",x[i],y[i],c[i]);
j++;
}
}
}
delete x;
delete y;
delete c;
delete f;
delete ncs;
return j;
}
int DepPairs::GetDots(int n1,float sdeps,float edeps,float DepLevel,float *datx1,int n)
{
int count=(edeps-sdeps)/DepLevel+1.5;
DPoint*p1=new DPoint[count+1];
for(int i=0;i<count;i++)
{
if(datx1)
{
p1[i]=DPoint(sdeps+i*DepLevel,datx1[i]);
}
else p1[i]=DPoint(sdeps+i*DepLevel,Curve1.GetValue(sdeps+i*DepLevel));
}
int datn=GetGDots(count,p1,n);
delete p1;
return datn;
}
int DepPairs::GetGDots(int n1,DPoint *x,int n)
{
int k=0;
double s1, s2;
int *spos=new int[n1];
for(int i=1; i<n1-1; i++)
{
s1 = x[i].y-x[i-1].y;//fun(x[i-2], x[i-1], x[i]);
s2 = x[i].y-x[i+1].y;//fun(x[i-1], x[i], x[i+1]);
if(s1<0&&s2<0){
spos[k]=i;
k++;
}
}
int kk=0;
for(int j=0;j<n;j++)
{
if(kk)k=kk;
kk=0;
for(int i=1; i<k-1; i++)
{
s1 = x[spos[i]].y-x[spos[i-1]].y;//fun(x[i-2], x[i-1], x[i]);
s2 = x[spos[i]].y-x[spos[i+1]].y;//fun(x[i-1], x[i], x[i+1]);
if(s1>0&&s2<0||s1>0&&s2>0||s1<0&&s2>0){
continue;
}
spos[kk]=spos[i];
kk++;
}
}
if(n<=0) kk=k;
for(int i=0;i<kk-1;i++)
{
int ppos=-1;
float d=99999;
float dd;
int len=spos[i+1]-1-spos[i];
int isl=sl/DepLevel;
if(len>isl) len=isl;
for(int m=spos[i];m<spos[i]+len;m++)
{
dd=(x[m].y-x[m+1].y);
if(dd<d)
{
ppos=m;
d=dd;
}
}
if(ppos>0){
spos[i]=ppos-1;
}
}
for(int i=0; i<kk; i++)
{
fprintf(fp_temp,"%f %f %f\n",x[spos[i]].x,x[spos[i]].x,1);
}
delete spos;
return kk;
}
int DepPairs::paral(int n1,float *x1,float *y1,float *c1,int retype)
{
float a,b,d;
float *x=new float[n1+1];
float *y=new float[n1+1];
float *c=new float[n1+1];
float *spos=new float[n1+1];
memset(spos,0,n1*sizeof(float));
rewind(fp_temp);
int n=0;
while(!feof(fp_temp))
{
fscanf(fp_temp,"%f%f%f",&x[n],&y[n],&c[n]);
n++;
if(n>=n1) break;
}
float xmax=999999.,disp1,disp2,disp;
int i,k,j;
if(n<3) {
delete x;
delete y;
delete c;
delete spos;
return -1;
}
i=-1;
float v1,v2,v3;
fclose(fp_temp);
QString dir=GetLogdataPath();
dir+="\\temp";
CreateDir((char *)dir.toStdString().c_str());
dir+="\\admtemp";
fp_temp=fopen(dir.toStdString().c_str(),"w+");
while(1) {
i=i+1;
if(i>=n) break;
disp1=y[i]-x[i];
disp2=y[i+1]-x[i+1];
if(disp1==disp2) {
i+=2;
for(j=i;j<n;j++) {
disp=y[j]-x[j];
if(disp==disp1)x[j-1]=xmax;
if(disp!=disp1) break;
}
if(disp==disp1) break;
i=j-1;
}
}
for(k=0,i=0;i<n;i++) {
if(x[i]<xmax) {
fprintf(fp_temp,"%f %f %f\n",x[i],y[i],c[i]);
k++;
}
}
delete x;
delete y;
delete c;
delete spos;
return k;
}
void DepPairs::DepthSort(float *xxx,int number)
{
int i5,j5;
float tempf5;
for (i5=0; i5<number; i5++)
for (j5=i5; j5<number; j5++)
if (xxx[i5] > xxx[j5])
{
tempf5 = xxx[i5];
xxx[i5] = xxx[j5];
xxx[j5] = tempf5;
}
}
int DepPairs::GetCompare(int datn,float sdeps,float edeps,float* datx1,float* daty1,float* datc1,float cmins,float xelast)
{
float *x=new float[datn+10];
float *y=new float[datn+10];
float *c=new float[datn+10];
rewind(fp_temp);
float wl2,xdr2;
int num1,number,numt,i,iwl,isl,iwl2,istep,jsize,ixdr2,idr,idr2;
int n1,kno;
int cccc;
if(wl==0)wl=10.;
if(sl==0)sl=4.;
if(step==0)step=1;
if(dr1==0)dr1=3.;
if(dr2==0)dr2=3.;
if(type!=0 &&sl<dr2/2)sl=dr2/2;
iwl=(int)(wl/DepLevel+0.5);
iwl2=iwl/2;
iwl=2*iwl2+1;
isl=(int)(sl/DepLevel+0.5);
istep=(int)(step/DepLevel+0.5);
idr=(int)(dr2/DepLevel+0.5);
idr2=idr/2;
idr=2*idr2+1;
xdr2=idr2*DepLevel;
if( type == 0 ) xdr2=0;
ixdr2=idr2;
if( type == 0 ) ixdr2=0;
jsize=iwl+2*isl;
kno=1;
if(emp!=0||size==0) detscl(sdeps,edeps);
num1 = 0;
number = iwl + 2 * isl+2*ixdr2;
numt = number - istep;
wl2=iwl2*DepLevel;
int n=0;
int maxnum=-99999;
iwl=(int)(wl/DepLevel+0.5);
iwl2=iwl/2;
iwl=2*iwl2+1;
while(!feof(fp_temp))
{
fscanf(fp_temp,"%f%f%f",&x[n],&y[n],&c[n]);
if(n)iwl=(x[n]-x[n-1])/DepLevel;
if(maxnum<iwl) maxnum=iwl;
n++;
if(n>=datn) break;
}
fclose(fp_temp);
QString dir=GetLogdataPath();
dir+="\\temp";
CreateDir((char *)dir.toStdString().c_str());
dir+="\\admtemp";
fp_temp=fopen(dir.toStdString().c_str(),"w+");
number = maxnum + 2 * isl+2*ixdr2;
float *datx=new float[number+1];
float *daty=new float[number+1];
float *datc=new float[number+1];
cccc = 0;
for(int j=0;j<n-1;j++){
int k=isl;
for (i=0; i<number; i++) {
datx[i]=0.;
daty[i]=0.;
datc[i]=0.;
}
for(float dep=x[j];dep<x[j+1];dep+=DepLevel)
{
if(datx1)
{
int pos=(dep-sdeps)/DepLevel;
if(dep<sdeps||dep>edeps) datx[k++]=-9999.0;
else datx[k++]=datx1[pos];
}
else datx[k++]=Curve1.GetValue(dep);
}
k=0;
for(float dep=x[j]-isl*DepLevel/2;dep<x[j+1]+isl*DepLevel;dep+=DepLevel)
{
if(daty1)
{
int pos=(dep-sdeps)/DepLevel;
if(dep<sdeps||dep>edeps) daty[k++]=-9999.0;
else daty[k++]=daty1[pos];
}
else daty[k++]=Curve2.GetValue(dep);
}
iwl=(x[j+1]-x[j])/DepLevel;
cor(iwl,isl,datx,daty,datc,size,emp);
MMax(x[j],x[j]-isl*DepLevel/2,isl,datc,cmins);
}
datn=optim(0,datx1,daty1,datc1,cmins,xelast);
delete x;
delete y;
delete c;
return datn;
}
int DepPairs::joincurve(float sdeps,float edeps,float *datx1,float *daty1,float *datc1,
float wl,float sl,float step,float dr1,float dr2,
int type,int size,int way,int emp,int ipar,
float cmins,float xelast)
{
float depth,depth1,dept,deps;
float wl2,xdr2;
int num1,number,numt,i,iwl,isl,iwl2,istep,jsize,ixdr2,idr,idr2;
int n1,kno;
int cccc;
if(wl==0)wl=10.;
if(sl==0)sl=4.;
if(step==0)step=1;
if(dr1==0)dr1=3.;
if(dr2==0)dr2=3.;
if(type!=0 &&sl<dr2/2)sl=dr2/2;
iwl=(int)(wl/DepLevel+0.5);
iwl2=iwl/2;
iwl=2*iwl2+1;
isl=(int)(sl/DepLevel+0.5);
istep=(int)(step/DepLevel+0.5);
idr=(int)(dr2/DepLevel+0.5);
idr2=idr/2;
idr=2*idr2+1;
xdr2=idr2*DepLevel;
if( type == 0 ) xdr2=0;
ixdr2=idr2;
if( type == 0 ) ixdr2=0;
jsize=iwl+2*isl;
kno=1;
if(emp!=0||size==0) detscl(sdeps,edeps);
num1 = 0;
number = iwl + 2 * isl+2*ixdr2;
numt = number - istep;
wl2=iwl2*DepLevel;
depth=sdeps-step+wl2+sl+xdr2;
depth1 = depth + step + step/2;
dept = iwl2 + isl;
int ret=TRUE;
QString dir=GetLogdataPath();
dir+="\\temp";
CreateDir((char *)dir.toStdString().c_str());
dir+="\\admtemp";
fp_temp=fopen(dir.toStdString().c_str(),"w+");
if(fp_temp==NULL)
{
2026-02-03 14:40:58 +08:00
QMessageBox::warning(nullptr, "提示", "Not create file :admtemp !");
return -1;
}
// ::GetProgressBar()->setVisible(true);
// ::GetProgressBar()->setRange(sdeps,edeps);
// ::GetProgressBarHint()->setVisible(true);
char showname[200];
sprintf(showname,"自动产生深度对比线...");
int datn=0;
if(ipar){
// ::GetProgressBarHint()->setText("查找层界面....");
QApplication::processEvents();
datn=GetDots(n1,sdeps,edeps,DepLevel,datx1,1);
// ::GetProgressBar()->setValue(sdeps+(edeps-sdeps)*0.1);
QApplication::processEvents();
// ::GetProgressBarHint()->setText("开始地层对比...");
// ::GetProgressBar()->setValue(sdeps+(edeps-sdeps)*0.3);
QApplication::processEvents();
datn=GetCompare(datn,sdeps,edeps,datx1,daty1,datc1,cmins,xelast);
// ::GetProgressBar()->setValue(edeps);
QApplication::processEvents();
}
else{
float *datx=new float[number+1];
float *daty=new float[number+1];
float *datc=new float[number+1];
cccc = 0;
for (i=0; i<number; i++) {
datx[i]=0.;
daty[i]=0.;
datc[i]=0.;
}
while(1)
{
depth=depth+step;
// ::GetProgressBar()->setValue(depth);
// ::GetProgressBarHint()->setText(showname);
QApplication::processEvents();
if(depth>edeps-wl2-sl-xdr2) {
break;
}
if ( depth > depth1 )
{
if ( numt > 0 )
{
for(i=0;i<numt;i++)
{
datx[i] = datx[istep+i];
daty[i] = daty[istep+i];
}
num1 = numt;
}
}
for(i=num1;i<number;i++)
{
deps=depth+(i-dept)*DepLevel;
if(emp!=0||size==0)
{
datx[i]=(Curve1.GetValue(deps)-xmin1)/(xmax1-xmin1);
daty[i]=(Curve2.GetValue(deps)-xmin2)/(xmax2-xmin2);
}
else
{
datx[i]=Curve1.GetValue(deps);
daty[i]=Curve2.GetValue(deps);
}
}
if(type!=0)
{
actvty(iwl,isl,&dr1,way,datx,datc);
actvty(iwl,isl,&dr2,way,daty,datc);
}
cor(iwl,isl,datx,daty,datc,size,emp);
n1=iwl+isl*2;
if(type==0) adjust(n1,way,datc);
cccc ++;
MMax(depth,isl,datx,daty,datc,cmins);
}
fseek(fp_temp,0,0);
delete datx;
delete daty;
delete datc;
datn=optim(0,datx1,daty1,datc1,cmins,xelast);
if(datn>0)datn=paral(datn,datx1,daty1,datc1,ipar);
}
fclose(fp_temp);
// ::GetProgressBar()->setVisible(false);
// ::GetProgressBarHint()->setVisible(false);
QApplication::processEvents();
return datn;
}
int DepPairs::adm(float stem,float etem ,float cmins,
float winLength,float stepLength, float searchLength,
float *datx1,float *daty1,float *datc1
)
{
char str[120];
int i;
int datn;
cmin = cmins;
datn = joincurve(stem,etem,datx1,daty1,datc1,
winLength,searchLength,stepLength,dr1,dr2,type,size,
way,emp,ipar,cmins,elast);
if ( datn < 0 )
{
sprintf(str,"%f--%f层段相关性太差相关系数均小于截止值%f ",stem,etem,cmins);
QMessageBox::warning(nullptr, "提示", str);
}
else {
float *datx,*daty,*datc;
m_nDepPairsNum=0;
datx = new float[datn+100];
if ( datx == NULL )
{
QMessageBox::warning(nullptr, "提示", "自动校深:没有足够的内存.");
return m_nDepPairsNum;
}
daty = new float[datn+100];
if ( daty == NULL )
{
QMessageBox::warning(nullptr, "提示", "自动校深:没有足够的内存.");
delete datx;
return m_nDepPairsNum;
}
datc = new float[datn+100];
if ( datc == NULL )
{
QMessageBox::warning(nullptr, "提示", "自动校深:没有足够的内存.");
delete datx;
delete daty;
return m_nDepPairsNum;
}
QString dir=GetLogdataPath();
dir+="\\temp";
CreateDir((char *)dir.toStdString().c_str());
dir+="\\admtemp";
fp_temp=fopen(dir.toStdString().c_str(),"rt");
if(fp_temp==NULL)
{
2026-02-03 14:40:58 +08:00
QMessageBox::warning(nullptr, "提示", "Not create file :admtemp !");
delete [] datx;
delete [] daty;
delete [] datc;
QDir r;
r.remove(dir.toStdString().c_str());
return -1;
}
i=0;
while(!feof(fp_temp))
{
fscanf(fp_temp,"%f%f%f",&datx[i],&daty[i],&datc[i]);
i++;
if(i>=datn) break;
}
if(!i) {
QMessageBox::warning(nullptr, "提示", "自动对比文件读写错误!");
fclose(fp_temp);
delete [] datx;
delete [] daty;
delete [] datc;
QDir r;
r.remove(dir.toStdString().c_str());
return -1;
}
datn=i;
fclose(fp_temp);
qs(datx,daty,datc,0,datn-1);
if(last_edep!=stem)
{
float ddep=0;
if(datn) {
ddep=datx[0]-stem;
FirstDep[m_nDepPairsNum]=daty[0]-ddep;
}
else FirstDep[m_nDepPairsNum]=stem;
SecondDep[m_nDepPairsNum]=stem;
corc[m_nDepPairsNum]=1.;
m_nDepPairsNum++;
}
for(i=0;i<datn;i++)
{
FirstDep[m_nDepPairsNum]=daty[i];
SecondDep[m_nDepPairsNum]=datx[i];
if(corc) corc[m_nDepPairsNum]=datc[i];
m_nDepPairsNum++;
}
float ddep=0;
if(datn) {
ddep=etem-datx[datn-1];
FirstDep[m_nDepPairsNum]=daty[datn-1]+ddep;
}
else FirstDep[m_nDepPairsNum]=etem;
SecondDep[m_nDepPairsNum]=etem;
if(corc) corc[m_nDepPairsNum]=1;
m_nDepPairsNum++;
delete [] datx;
delete [] daty;
delete [] datc;
QDir r;
r.remove(dir.toStdString().c_str());
}
return m_nDepPairsNum;
}
int DepPairs::AdmRun(float *curve1,float *curve2,float *cc)
{
return adm(StartDepth,EndDepth ,cmin,wl,step,sl,curve1,curve2,cc);
}