SAR目标仿真,RD算法

距离多普勒算法(RDA)是在1976年至1978年为民用星载SAR提出的,它兼顾了成熟、简单、高效和精确等因素,至今仍是使用最广泛的成像算法。它分为三个步骤: 距离压缩(Range Compression, RC)、距离迁移校正(Range Cell Migration correction, RCMC)和方位压缩(Azimuth Compression, AC)。模拟仿真如下:

%============ SAR 三点回波模拟 正侧视RD算法 ======================
clc;
close all;
clear all;

c=3e8; %光速
f0=10e9; %雷达工作频率
lamd=c/f0; %雷达工作波长
Tp=1e-6; %
B=100e6; %
PRF=3600; %脉冲重复频率
fs=120e6; %
K=B/Tp; %

v=150; %等效雷达运动速度

x1=0;y1=0; %目标1位置
x2=0;y2=-50; %目标2位置
x3=0;y3=50; %目标3位置

Na=4096; % 方位向采样点数为4096
bet=18/180*pi; %
R0=1000;

Nr=256;
Tr=Nr/fs;
Ta=Na/PRF;

tr=2*R0/c+linspace(-Tr/2,Tr/2,Nr);
ta=linspace(-Ta/2,Ta/2,Na);
xr=v*ta;
yr=R0;
fr=fftshift(linspace(-fs/2,fs/2,Nr));
fa=linspace(-PRF/2,PRF/2,Na);
r=tr*c/2;
r_c=sqrt(xr.^2+yr^2);

%====================================产生回波===============================

rec_signal=zeros(Na,Nr);
for i=1:Na
puls_signal=zeros(1,Nr);
r1=sqrt((y1+yr).^2+(xr(i)-x1).^2);
cita1=atan(abs(x1-xr(i))/(yr+y1));
if (abs(cita1)<bet/2)
puls_signal=exp(j*(pi*K*(tr-2*r1/c).^2-2*pi*f0*2*r1/c)).*(abs(tr-2*r1/c)<Tp/2);
end

r2=sqrt((y2+yr).^2+(xr(i)-x2).^2);
cita2=atan(abs(x2-xr(i))/(yr+y2));
if (abs(cita2)<bet/2)
puls_signal=puls_signal+exp(j*(pi*K*(tr-2*r2/c).^2-2*pi*f0*2*r2/c)).*(abs(tr-2*r2/c)<Tp/2);%目标2回波
end

r3=sqrt((y3+yr).^2+(xr(i)-x3).^2);
cita3=atan(abs(x3-xr(i))/(yr+y3));
if (abs(cita3)<bet/2)
puls_signal=puls_signal+exp(j*(pi*K*(tr-2*r3/c).^2-2*pi*f0*2*r3/c)).*(abs(tr-2*r3/c)<Tp/2);%目标3回波
end
rec_signal(i,:)=puls_signal;
end
figure(1);
subplot(221);
imagesc(abs(rec_signal));
title(‘雷达原始信号1’);

%================================距离向脉压=================================

tran_signal=exp(1j*pi*fr.^2/K); %.*(abs(fr)<B/2); 这个范围限制最好要有,但为了方便改写C语言,不要范围图像也还可以,真正的改写要加上此范围
for i=1:Na
t = fft(rec_signal(i,:));
tt =t .*tran_signal;
rec_signal(i,:)=ifft(tt);
end
subplot(222);
imagesc(abs(rec_signal));
title(‘距离压缩后仿真图2’);
xlabel(‘距离向(采样点)’);
ylabel(‘方位向(采样点)’);

% %==============================插值法消除距离徙动===========================
signal_fft = fft(rec_signal,[],1);
rec_signal1=fftshift(signal_fft,1); %方位向变为频域
% rec_signal=fft(rec_signal,[],1);
% %方位向变为频域,只方位向FFT不行,需要FFT后进行fftshift操作

subplot(223);
imagesc(abs(rec_signal1));
title(‘方位向FFT后仿真图3’);
xlabel(‘距离向(采样点)’);
ylabel(‘方位向(采样点)’);

for i=1:Na
delt_r(i)=lamd.^2*R0*fa(i).^2/8/v.^2;
rec_signal(i,:)=interp1(r,rec_signal1(i,:),r+delt_r(i),’nearest’); %插值 分段线性插值 linear 球面插值spline 三次多项式插值pchip 临近插值nearest
end

subplot(224);
rec_signal_abs = abs(rec_signal);
imagesc(abs(rec_signal));
title(‘距离徙动校正后仿真图4’);
xlabel(‘距离向(采样点)’);
ylabel(‘方位向(采样点)’);

%================================方位向进行脉压===========================
%{
for ii=1:Nr
Ka=2*v.^2/lamd/r(ii);
puls_a=exp(-j*pi*fa.^2/Ka).*(abs(fa)<Ka*Ta/2);
rec_signal(ii,:)=ifft(rec_signal(ii,:).*puls_a);
end
%}
%以下为降低方位向旁瓣 18°可按照大斜视角做
for ii=1:Nr
Ka=sqrt( 1-(lamd.^2.*fa.^2) / ( 4*(v^2) ) );
theta = 2*pi.*tr(ii)*f0.*Ka;
puls_a=exp(j*2*pi.*tr(ii)*f0.*Ka);
rec_signal1(:,ii) = rec_signal(:,ii).*puls_a.’;
rec_signal2(:,ii)=ifft(rec_signal1(:,ii));
end
figure(2);
subplot(121);
imagesc(abs(rec_signal2));
title(‘目标图像5’);

% colormap(gray);
subplot(122);
rec_signal2_abs = abs(rec_signal2);
imagesc(abs(rec_signal2));
title(‘目标图像6’);
xlabel(‘距离向(采样点)’);
ylabel(‘方位向(采样点)’);

参考文献:https://blog.csdn.net/a1367666195/article/details/106614707

 

本文为原创文章,转载请注明出处!

超声检测技术

一  超声TOFD检测技术的起源和国外发展现状

TOFD(Time Of Flight Diffraction)技术是1972年国际原子能中心的哈韦尔(英国原子能权威人士-UKKAEA)提议下发展而来。TOFD最初的发展仅仅是作为定量工具,最初的想法是:使用常规技术探测到缺陷后使用TOFD进行精确的定量和监测在线设备裂纹的扩展(例如检测压力容器)。很多年以来TOFD一直在实验室里,各国做过大量实验直到八十年代才为业界所认同;在这些实验中,用事实证明了TOFD在可靠性和精度方面都是非常好的技术。

自上世纪 9O年代起,超声TOFD检测法在国外工业无损检测领域已得到广泛应用,欧、美、日均已推出相应的应用标准。1992年英国标准BS7706问世,1996年美国ASME标准将其列入规范案例2235和第ⅴ卷《无损检测》附录,2000年欧、日分别推出专用标准。用于不同壁厚的承压设备焊接接头的制造和在用检测。

二  超声TOFD检测技术国内发展现状

3.10  平行扫查(parallel scan)

探头运动方向与声束方向平行的扫查方式。

3.11  非平行扫查(non-parallel scan)

探头运动方向与声束方向垂直的扫查方式,一般指探头对称布置于焊缝中心线两侧沿焊缝长度方向(X轴)的扫查方式。

3.12  偏置非平行扫查(offset-scan)

偏移焊缝中心线一定距离的非平行扫查。

3.13  A扫描信号(A-scan signal)

超声波信号的波形显示图,通常水平轴表示声波的传播时间,垂直轴表示波幅。

3.14  TOFD图像(TOFD image)

TOFD数据的二维显示,是将扫查过程中采集的A扫描信号连续拼接而成,一个轴代表探头移动距离,另一个轴代表深度,一般用灰度表示A扫描信号的幅度。

四  超声TOFD检测原理

4.1.1  TOFD是一种利用超声波衍射现象、非基于波幅的自动超声检测方法。通常使用纵波斜探头,采用一发一收模式。

4.1.2  在工件无缺陷部位,发射超声脉冲后,首先到达接收探头的是直通波,然后是底面反射波。有缺陷存在时,在直通波和底面反射波之间,接收探头还会接收到缺陷处产生的衍射波。除上述波外,还有缺陷部位和底面因波型转换产生的横波,一般会迟于底面反射波到达接收探头。

本文为原创文章,转载请注明出处!

MFC实现三维图像绘制

绘制物体的面时考虑两个问题:

1. 如何绘制?

2. 如何消隐不可见面?

针对这两个问题,分别采取有效边表算法和Z-Buffer缓冲器算法。

有效边表:

#include “Pi3.h”
#include “T2.h”

class CAET
{
public:
CAET();
virtual ~CAET();
public:
double x;//当前扫描线与有效边交点的x坐标
int yMax;//边的最大y值
double k;//斜率的倒数(x的增量)
CPi3 ps;//边的起点
CPi3 pe;//边的终点
CAET *pNext;
CT2 ts;//起点纹理坐标
CT2 te;//终点纹理坐标
};
桶表:

#include “AET.h”
class CBucket
{
public:
CBucket();
virtual ~CBucket();
public:
int ScanLine; //扫描线
CAET *pET; //桶上的边表指针
CBucket *pNext;
};
具体绘制实现(这里涉及到物体的光照和纹理模型):

#include “Bucket.h”
#include “Vector3.h”
#include “Lighting.h”
#include “Pi3.h”
#include “Texture.h”
#include “T2.h”

class CFill
{
public:
CFill(void);
virtual ~CFill(void);
void SetPoint(CPi3 *p, int, CT2*);//类的初始化
void SetPoint(CPi3* p, int);//类的初始化
void CreateBucket();//创建桶
void CreateEdge();//边表
void AddEt(CAET *);//合并ET表
void EtOrder();//ET表排序
void Gouraud(CDC *, int,CLighting, CMaterial*,CP3, CTexture*,BOOL);//填充多边形
void Gouraud(CDC*, int, CLighting, CMaterial*, CP3);//填充多边形
CRGB Interpolation(double, double, double, CRGB, CRGB);//颜色线性插值
double Interpolation(double, double, double, double, double);//深度线性插值
CVector3 Interpolation(double, double, double, CVector3, CVector3);//法向量线性插值
CT2 Interpolation(double m, double m0, double m1, CT2 T0, CT2 T1);//纹理地址线性插值
void ClearMemory();//清理内存
void DeleteAETChain(CAET* pAET);//删除边表
CRGB GetTextureColor(double u, double v, CTexture* pTexture);//获得纹理颜色
void InitialDepthBuffer(int nWidth, int nHeight, double zDepth);//初始化深度缓冲器
public:
int PNum;//顶点个数
CT2* T;//纹理动态数组
CPi3 *P;//顶点坐标动态数组
CAET *pHeadE, *pCurrentE, *pEdge; //有效边表结点指针
CBucket *pHeadB, *pCurrentB; //桶表结点指针
double** zBuffer;//深度缓冲器
int nWidth, nHeight;//缓冲区宽度与高度
};

#include “pch.h”
#include “Fill.h”
#include “math.h”
#define ROUND(f) int(floor(f+0.5))//四舍五入宏定义

CFill::CFill(void)
{
PNum = 0;
P = NULL;
pEdge = NULL;
pHeadB = NULL;
pHeadE = NULL;
}

CFill::~CFill(void)
{
if (P != NULL)
{
delete[] P;
P = NULL;
}
for (int i = 0;i < nWidth;i++)
{
delete[] zBuffer[i];
zBuffer[i] = NULL;
}
if (zBuffer != NULL)
{
delete zBuffer;
zBuffer = NULL;
}
ClearMemory();
}

void CFill::SetPoint(CPi3 *p, int m, CT2* t)
{
if (P != NULL)
{
delete[] P;
P = NULL;
}
P = new CPi3[m];//创建一维动态数组
T = new CT2[m];
for (int i = 0; i < m; i++)
{
P[i] = p[i];
P[i].y = ROUND(P[i].y);
T[i] = t[i];
}
PNum = m;
}

void CFill::SetPoint(CPi3* p, int m)
{
if (P != NULL)
{
delete[] P;
P = NULL;
}
P = new CPi3[m];//创建一维动态数组
for (int i = 0; i < m; i++)
{
P[i] = p[i];
P[i].y = ROUND(P[i].y);
}
PNum = m;
}

void CFill::CreateBucket()//创建桶表
{
int yMin, yMax;
yMin = yMax = P[0].y;
for (int i = 0; i < PNum; i++)//查找多边形所覆盖的最小和最大扫描线
{
if (P[i].y < yMin)
{
yMin = P[i].y;//扫描线的最小值
}
if (P[i].y > yMax)
{
yMax = P[i].y;//扫描线的最大值
}
}
for (int y = yMin; y <= yMax; y++)
{
if (yMin == y)//如果是扫描线的最小值
{
pHeadB = new CBucket;//建立桶的头结点
pCurrentB = pHeadB;//pCurrentB为CBucket当前结点指针
pCurrentB->ScanLine = yMin;
pCurrentB->pET = NULL;//没有链接边表
pCurrentB->pNext = NULL;
}
else//其他扫描线
{
pCurrentB->pNext = new CBucket;//建立桶的其他结点
pCurrentB = pCurrentB->pNext;
pCurrentB->ScanLine = y;
pCurrentB->pET = NULL;
pCurrentB->pNext = NULL;
}
}
}

void CFill::CreateEdge()//创建边表
{
for (int i = 0; i < PNum; i++)
{
pCurrentB = pHeadB;
int j = (i + 1) % PNum;//边的第2个顶点,P[i]和P[j]点对构成边
if (P[i].y < P[j].y)//边的终点比起点高
{
pEdge = new CAET;
pEdge->x = P[i].x;//计算ET表的值
pEdge->yMax = P[j].y;
pEdge->k = (P[j].x – P[i].x) / (P[j].y – P[i].y);//代表1/k
pEdge->ps = P[i];//绑定顶点和颜色
pEdge->pe = P[j];
pEdge->ts= T[i];//绑定纹理
pEdge->te = T[j];
pEdge->pNext = NULL;
while (pCurrentB->ScanLine != P[i].y)//在桶内寻找当前边的yMin
{
pCurrentB = pCurrentB->pNext;//移到yMin所在的桶结点
}
}
if (P[j].y < P[i].y)//边的终点比起点低
{
pEdge = new CAET;
pEdge->x = P[j].x;
pEdge->yMax = P[i].y;
pEdge->k = (P[i].x – P[j].x) / (P[i].y – P[j].y);
pEdge->ps = P[i];
pEdge->pe = P[j];
pEdge->ts = T[i];//绑定纹理
pEdge->te = T[j];
pEdge->pNext = NULL;
while (pCurrentB->ScanLine != P[j].y)
{
pCurrentB = pCurrentB->pNext;
}
}
if (P[i].y != P[j].y)
{
pCurrentE = pCurrentB->pET;
if (pCurrentE == NULL)
{
pCurrentE = pEdge;
pCurrentB->pET = pCurrentE;
}
else
{
while (pCurrentE->pNext != NULL)
{
pCurrentE = pCurrentE->pNext;
}
pCurrentE->pNext = pEdge;
}
}
}
}

void CFill::Gouraud(CDC *pDC, int light_type,CLighting light, CMaterial* material, CP3 EyePoint, CTexture* pTexture,BOOL IfBump)//填充多边形
{
CAET *pT1 = NULL, *pT2 = NULL;
pHeadE = NULL;
for (pCurrentB = pHeadB; pCurrentB != NULL; pCurrentB = pCurrentB->pNext)
{
for (pCurrentE = pCurrentB->pET; pCurrentE != NULL; pCurrentE = pCurrentE->pNext)
{
pEdge = new CAET;
pEdge->x = pCurrentE->x;
pEdge->yMax = pCurrentE->yMax;
pEdge->k = pCurrentE->k;
pEdge->ps = pCurrentE->ps;
pEdge->pe = pCurrentE->pe;
pEdge->pNext = NULL;
pEdge->ts = pCurrentE->ts;
pEdge->te = pCurrentE->te;
pEdge->pNext = NULL;
AddEt(pEdge);
}
EtOrder();
pT1 = pHeadE;
if (pT1 == NULL)
{
return;
}
while (pCurrentB->ScanLine >= pT1->yMax)//下闭上开
{
CAET * pAETTEmp = pT1;
pT1 = pT1->pNext;
delete pAETTEmp;
pHeadE = pT1;
if (pHeadE == NULL)
return;
}
if (pT1->pNext != NULL)
{
pT2 = pT1;
pT1 = pT2->pNext;
}
while (pT1 != NULL)
{
if (pCurrentB->ScanLine >= pT1->yMax)//下闭上开
{
CAET* pAETTemp = pT1;
pT2->pNext = pT1->pNext;
pT1 = pT2->pNext;
delete pAETTemp;
}
else
{
pT2 = pT1;
pT1 = pT2->pNext;
}
}
CRGB Ca, Cb, Cf;//Ca、Cb代表边上任意点的颜色,Cf代表面上任意点的颜色
double za, zb, zf;
CVector3 va, vb, vf;
CT2 ta, tb, tf;//ta和tb代表边上任意点的纹理地址,tf代表面上任意点的纹理地址
ta = Interpolation(pCurrentB->ScanLine, pHeadE->ps.y, pHeadE->pe.y, pHeadE->ts, pHeadE->te);
tb = Interpolation(pCurrentB->ScanLine, pHeadE->pNext->ps.y, pHeadE->pNext->pe.y, pHeadE->pNext->ts, pHeadE->pNext->te);
Ca = Interpolation(pCurrentB->ScanLine, pHeadE->ps.y, pHeadE->pe.y, pHeadE->ps.c, pHeadE->pe.c);
Cb = Interpolation(pCurrentB->ScanLine, pHeadE->pNext->ps.y, pHeadE->pNext->pe.y, pHeadE->pNext->ps.c, pHeadE->pNext->pe.c);
za = Interpolation(pCurrentB->ScanLine, pHeadE->ps.y, pHeadE->pe.y, pHeadE->ps.z, pHeadE->pe.z);
zb = Interpolation(pCurrentB->ScanLine, pHeadE->pNext->ps.y, pHeadE->pNext->pe.y, pHeadE->pNext->ps.z, pHeadE->pNext->pe.z);
va = Interpolation(pCurrentB->ScanLine, pHeadE->ps.y, pHeadE->pe.y, pHeadE->ps.v, pHeadE->pe.v);
vb = Interpolation(pCurrentB->ScanLine, pHeadE->pNext->ps.y, pHeadE->pNext->pe.y, pHeadE->pNext->ps.v, pHeadE->pNext->pe.v);
BOOL Flag = FALSE;
double xb, xe;//扫描线和有效边相交区间的起点和终点坐标
for (pT1 = pHeadE; pT1 != NULL; pT1 = pT1->pNext)
{
if (Flag == FALSE)
{
xb = pT1->x;
Flag = TRUE;
}
else
{
xe = pT1->x;
for (double x = xb; x < xe; x++)//左闭右开
{
if (ROUND(x) + nWidth / 2 < nWidth && ROUND(x) + nWidth / 2 > 0 && pCurrentB->ScanLine + nHeight / 2 < nHeight && pCurrentB->ScanLine + nHeight / 2 > 0) {
zf = Interpolation(x, xb, xe, za, zb);
if (zf >= zBuffer[ROUND(x) + nWidth / 2][pCurrentB->ScanLine + nHeight / 2])//如果当前采样点的深度小于帧缓冲器中原采样点的深度)
{
tf = Interpolation(x, xb, xe, ta, tb);
tf.c = GetTextureColor(ROUND(tf.u), ROUND(tf.v), pTexture);
material->SetDiffuse(tf.c);//用纹理颜色作为材质的漫反射光反射率
material->SetAmbient(tf.c);//用纹理颜色作为材质的环境光反射率
if (light_type == 1) {
Cf = Interpolation(x, xb, xe, Ca, Cb);
}
else if (light_type == 2) {
vf = Interpolation(x, xb, xe, va, vb);
// 凹凸纹理
if (IfBump) {
CRGB frontU = GetTextureColor(ROUND(tf.u – 1), ROUND(tf.v), pTexture);
CRGB backU = GetTextureColor(ROUND(tf.u + 1), ROUND(tf.v), pTexture);
double Bu = (frontU.blue – backU.blue) / 2;
CRGB frontV = GetTextureColor(ROUND(tf.u), ROUND(tf.v – 1), pTexture);
CRGB backV = GetTextureColor(ROUND(tf.u), ROUND(tf.v + 1), pTexture);
double Bv = (frontV.blue – backV.blue) / 2;
CVector3 DNormal = CVector3(Bu, Bv, 0);
double BumpScale = 100.0;
vf = vf + BumpScale * DNormal;
vf = vf.Normalize();
}
Cf = light.Illuminate(EyePoint, CP3(ROUND(x), pCurrentB->ScanLine, zf), vf, material);
}
zBuffer[ROUND(x) + nWidth / 2][pCurrentB->ScanLine + nHeight / 2] = zf;//使用当前采样点的深度更新深度缓冲器
pDC->SetPixelV(ROUND(x), pCurrentB->ScanLine, RGB(Cf.red * 255, Cf.green * 255, Cf.blue * 255));
}
}
}
Flag = FALSE;
}
}
for (pT1 = pHeadE; pT1 != NULL; pT1 = pT1->pNext)//边的连续性
{
pT1->x = pT1->x + pT1->k;
}
}
}

void CFill::AddEt(CAET *pNewEdge)//合并ET表
{
CAET *pCE = pHeadE;
if (pCE == NULL)
{
pHeadE = pNewEdge;
pCE = pHeadE;
}
else
{
while (pCE->pNext != NULL)
{
pCE = pCE->pNext;
}
pCE->pNext = pNewEdge;
}
}
void CFill::EtOrder()//边表的冒泡排序算法
{
CAET *pT1 = NULL, *pT2 = NULL;
int Count = 1;
pT1 = pHeadE;
if (NULL == pT1)
{
return;
}
if (NULL == pT1->pNext)//如果该ET表没有再连ET表
{
return;//桶结点只有一条边,不需要排序
}
while (NULL != pT1->pNext)//统计结点的个数
{
Count++;
pT1 = pT1->pNext;
}
for (int i = 1; i < Count; i++)//冒泡排序
{
pT1 = pHeadE;
if (pT1->x > pT1->pNext->x)//按x由小到大排序
{
pT2 = pT1->pNext;
pT1->pNext = pT1->pNext->pNext;
pT2->pNext = pT1;
pHeadE = pT2;
}
else
{
if (pT1->x == pT1->pNext->x)
{
if (pT1->k > pT1->pNext->k)//按斜率由小到大排序
{
pT2 = pT1->pNext;
pT1->pNext = pT1->pNext->pNext;
pT2->pNext = pT1;
pHeadE = pT2;
}
}
}
pT1 = pHeadE;
while (pT1->pNext->pNext != NULL)
{
pT2 = pT1;
pT1 = pT1->pNext;
if (pT1->x > pT1->pNext->x)//按x由小到大排序
{
pT2->pNext = pT1->pNext;
pT1->pNext = pT1->pNext->pNext;
pT2->pNext->pNext = pT1;
pT1 = pT2->pNext;
}
else
{
if (pT1->x == pT1->pNext->x)
{
if (pT1->k > pT1->pNext->k)//按斜率由小到大排序
{
pT2->pNext = pT1->pNext;
pT1->pNext = pT1->pNext->pNext;
pT2->pNext->pNext = pT1;
pT1 = pT2->pNext;
}
}
}
}
}
}

CRGB CFill::GetTextureColor(double u, double v, CTexture* pTexture)//获得纹理颜色
{
/*检测图片的边界,防止越界*/
if (u < 0) u = 0; if (v < 0) v = 0;
if (u > pTexture->bmp.bmWidth – 1) u = pTexture->bmp.bmWidth – 1;
if (v > pTexture->bmp.bmHeight – 1) v = pTexture->bmp.bmHeight – 1;
/*查找对应纹理空间的颜色值*/
v = pTexture->bmp.bmHeight – 1 – v;
int position = ROUND(v * pTexture->bmp.bmWidthBytes + 4 * u);//颜色分量位置
COLORREF color = RGB(pTexture->image[position + 2], pTexture->image[position + 1], pTexture->image[position]);//获取颜色值
return CRGB(GetRValue(color) / 255.0, GetGValue(color) / 255.0, GetBValue(color) / 255.0);
}

CRGB CFill::Interpolation(double t, double t1, double t2, CRGB c1, CRGB c2)//颜色线性插值
{
CRGB c;
c = (t – t2) / (t1 – t2)*c1 + (t – t1) / (t2 – t1)*c2;
return c;
}

double CFill::Interpolation(double t, double t1, double t2, double z1,double z2)//线性插值
{
double z;
z = (t – t2) / (t1 – t2) * z1 + (t – t1) / (t2 – t1) * z2;
return z;
}

CVector3 CFill::Interpolation(double x, double xa, double xb, CVector3 va, CVector3 vb)//法向量线性插值
{
CVector3 v;
v = (x – xb) / (xa – xb) * va + (x – xa) / (xb – xa) * vb;
return v;
}

CT2 CFill::Interpolation(double m, double m0, double m1, CT2 T0, CT2 T1)//纹理地址线性插值
{
CT2 Texture;
Texture = (m1 – m) / (m1 – m0) * T0 + (m – m0) / (m1 – m0) * T1;
return Texture;
}

void CFill::ClearMemory()//安全删除所有桶和桶上面的边
{
DeleteAETChain(pHeadE);
CBucket *pBucket = pHeadB;
while (pBucket != NULL)// 针对每一个桶
{
CBucket * pBucketTemp = pBucket->pNext;
DeleteAETChain(pBucket->pET);
delete pBucket;
pBucket = pBucketTemp;
}
pHeadB = NULL;
pHeadE = NULL;
}

void CFill::DeleteAETChain(CAET* pAET)
{
while (pAET != NULL)
{
CAET* pAETTemp = pAET->pNext;
delete pAET;
pAET = pAETTemp;
}
}

void CFill::InitialDepthBuffer(int nWidth, int nHeight, double zDepth)//初始化深度缓冲
{
this->nWidth = nWidth, this->nHeight = nHeight;
zBuffer = new double* [nWidth];
for (int i = 0;i < nWidth;i++)
zBuffer[i] = new double[nHeight];
for (int i = 0;i < nWidth;i++)//初始化深度缓冲
for (int j = 0;j < nHeight;j++)
zBuffer[i][j] = zDepth;
}

————————————————
版权声明:本文为CSDN博主「仿生葡萄」的原创文章,遵循CC 4.0 BY-SA版权协议,转载请附上原文出处链接及本声明。
原文链接:https://blog.csdn.net/uchiha_rin/article/details/122132213

本文为原创文章,转载请注明出处!

小波探测算法

目前针对声波远探测测井异常数据的识别主要方法为波形幅度法,该方法采用波形幅度值和预设的临界值对各个深度点的波形进行坏道识别,认为某深度点的波形幅度值大于临界值,则确定该深度点的波形为坏道。事实上,声波远探测可以简化视为,随着时间增加,其波形幅度出现正负交替特征,叠加后的平均幅度(或者幅度绝对值的平均)变化并不能客观反映坏道特征。为此,本发明采用基于平均幅度的总体标准差方法,进行水平井声波远探测测井资料坏道信号识别,从而实现资料质量控制。

基于平均幅度的总体标准差计算公式为:

式中Xi为一个波列中某个时刻的数据,n为波列中数据个数,QC_Ave为平均幅度,QC_Std为总体标准差。

根据实际情况设置识别阈值σ,当时则视为坏道信号。

如图5在5973.5m、5979m、5994m等深度点平均幅度QC_Ave没有明显变化,依靠QC_Ave不能有效指示坏道信号;而QC_Std在上述井段出现明显的指状尖峰特征,指示存在波形异常,判断为坏道信号,根据波形和频率分析结果(图5中第三道与第四道)已证实上述深度的存在声波远探测测井数据异常。sheng

本文为原创文章,转载请注明出处!

OpendTect基础模块

为用户提供了一个开放式的平台。用户可在平台上进行二次开发,增加新的模块和应用程序。基础模块中,OpendTect提供了以下功能:

1、数据的输入与输出:OpendTect可支持SEGY、LAS、ASCII格式的地震、钻井、测井、层位等常用数据的输入与输出。

2、地震反射层位的追踪与编辑:通过基础模块,可以在OpendTect系统中开展2D、3D地震层位的对比追踪。同时可以对其它解释软件所解释的地震反射层位进行编辑。

3、常规地震属性计算:在基础模块中可以通过等时、沿层等方式对“三瞬”、均方根振幅、平均频率、平均能量等常规地震属性进行计算、分析。

4、特有地震属性计算:OpendTect针对断层、小型褶皱及烃类异常等开发了如相似性、中值曲率、吸收系数等特有的地震属性。这些特有的地震属性较之常规属性不仅更具针对性,应用效果也得到改善。

5、频谱分解:应用频谱分解针对不同频率段提取地震属性已经被广泛应用。OpendTect在基础模块中也提供了频谱分解技术。可将时深域转换为频率域,开展不同频段地震属性提取。

6、三维可视化:在OpendTect基础平台上,可对地震数据体、层位、属性计算结果等进行立体、透视、透明等三维可视化显示。

双平方根算子地震波偏移

  1. void main()
  2. {
  3.   int b1, b2, b3, nz, ny, nx, i,j,k,kk, nshot, ndim, is,order,n;
  4.   float sz, sy, sx, o1, o2, o3, dz, dy, dx, vel0, br1, br2, br3,vel1,vel2,vel3;
  5.   float *v,*inface1,*inface2;
  6.   dz=10;
  7.   dy=20;
  8.   dx=20;
  9.   nz=80;
  10.   ny=100;
  11.   nz=100;
  12.   b1=1;
  13.   b2=1;
  14.   b3=1;
  15.   order=2;
  16.   vel1=800;
  17.   vel2=2000;
  18.   vel3=4000;
  19.   vel0=1.0/vel1;
  20.   n=nx*ny*nz;
  21.   v=new float[n];
  22.   inface1=new float[nx*ny];
  23.   inface2=new float[nx*ny];
  24.   for(i=0;i<n;i++)
  25.   {
  26.   v[i]=1.0/vel1;
  27.   }
  28.   for(i=0;i<nx*ny;i++)
  29.   {
  30.           inface1[i]=20+i*nz;
  31.   }
  32.   for(j=0;j<ny;j++)
  33.   {
  34.   for(i=0;i<nx;i++)
  35.   {
  36.                   int kk=int(700-3*i/20);
  37.                   inface2[j*nx+i]=kk;
  38.   for(k=20;k<kk;k++)
  39.   {
  40.   v[(i*ny+j)*nz+k]=1.0/vel2;
  41.   }
  42.   for(k=kk;k<nz;k++)
  43.   {
  44.   v[(i*ny+j)*nz+k]=1.0/vel3;
  45.   }
  46.   }
  47.   }
  48.   ifstream vel(“D:double_v”);
  49.   for(i=0;i<n;i++)
  50.   {
  51.   vel<<v[i]<<“t”;
  52.   }
  53.   int ngx,ngy,nsx,nsy;
  54.   float dgy,dgx,dsx,dsy;

本文为原创文章,转载请注明出处!

开源ext2read代码走读之--“\\\\.\\PhysicalDrive0”意义

在ext2read中读取ext4文件系统的代码中,读取硬盘中的信息时,定义了以下的宏,那么这个宏是什么意思呢?

#define DEVICE “\\\\.\\PhysicalDrive0″是什么意思?

由于”\”是C/C+中转义符, “\\\\.\\”就相当于\\.\,那么以上的宏定义中的“\\\.\\PhysicalDrive0”就等价于“\\.\PhysicalDrive0”

在Windows中 \\.\ 前缀用于标识设备,其中的”.”表示本地计算机。
比如\\.\PhysicalDrive0表示本机的物理驱动器0(一般是主硬盘),
\\.\COM1表示本机的1号串行口。

本文为原创文章,转载请注明出处!

c 语言读取硬盘序列号

c++读取硬盘序列号:有时在将程序或者打包库提供给第三方时,为了防止传播泛滥,或者有一定的版权问题,需要绑定特定的计算机设备。此时就需要读取计算机的一些硬件资源(硬盘、cpu、bios等),来计算一个验证码,达到一机一码的目的。

软件查看硬盘序列号

借助diskgenius查看硬盘序列号,选中硬盘,即可看到在下方有序列号。不过貌似ssd和机械硬盘的序列号格式是不一样的

ssd: 12位序列号

机械硬盘: 8位序列号

c++实现读取硬盘序列号

参考资料

(https://msdn.microsoft.com/en-us/library/aa389273(v=vs.85)

代码实现

代码实现使用的是微软提供的wmi providers库, 读取不同硬件设备有不同的类,具体查看上面的msdn链接。

– 硬盘 [win32_physicalmemory class](https://msdn.microsoft.com/en-us/library/aa394346(v=vs.85)

– cpu [win32_processor class](https://msdn.microsoft.com/en-us/library/aa394373(v=vs.85)

– bios [win32_bios class](https://msdn.microsoft.com/en-us/library/aa394077(v=vs.85)

本文为原创文章,转载请注明出处!

IEEE与IBM格式

不同体系架构的计算机系统中, 浮点数在表
达和存储上存在较大差别。 在 SEGY 地震资料
中, 数据 可 采用 IEEE 或 者IBM 格式 浮 点 数 存
储 [ 1 ] 。 特定计算机硬件平台, 只 能处理与自 己 硬
件平台相同格式的 浮点数, 若浮点数格式不属于
目 标处理机能够处理的格式, 则需进行格式转换。
不同类型浮点数之间转换, 不仅涉及到阶码和尾
数的表达, 还涉及到浮点数存储的方式, 以及位和
位域的操作等问题。 为此, 本文设计了 一个函 数
库并提供 API 接口, 用户 可以简单直观地实现两
种类型浮点数之间的转换。

本文为原创文章,转载请注明出处!

光纤传感技术

  • 一、前言

光纤传感技术是20世纪70年代伴随光纤通信技术的发展而迅速发展起来的新型传感技术,国外一些发达国家对光纤传感技术的应用研究已取得丰硕成果,不少光纤传感系统已实用化,成为替代传统传感器的商品。

在油田的开发过程中,人们需要知道在产液或注水过程中有关井内流体的持性与状态的详细资料,这就要用到石油测井,其可靠性和准确性是至关重要的,而传统的电子基传感器无法在井下恶劣的环境诸如高温、高压、腐蚀、地磁地电干扰下工作。光纤传感器可以克服这些困难,其对电磁干扰不敏感而且能承受极端条件,包括高温、高压(几十兆帕以上)以及强烈的冲击与振动,可以高精度地测量井筒和井场环境参数,同时,光纤传感器具有分布式测量能力,可以测量被测量的空间分布,给出剖面信息。而且,光纤传感器横截面积小,外形短,在井筒中占据空间极小。

光纤传感器在地球物理测井领域取得了长足的进步,全世界各大石油生产公司、测井服务公司以及各种光纤传感器研发机构和企业都参加了研究、开发过程。为了开拓光纤传感器的应用领域,本文综述了光纤传感器在地球物理测井领域的研究与进展,希望其研究能够对进一步提高石油开发的水平作出贡献。

二、光纤传感器在测井上的研究进展

1、储层参数监测

(1)压力监测

由于开发方案的需要,对油藏压力的管理需要特别谨慎,这样做的目的是减少因在低于泡点压力的状态下开采所造成的原油损失,减少在注气过程中因油藏超压将原油挤入含水层所造成的原油损失。传统的井下压力监测采用的传感器主要有应变压力计和石英晶体压力计,应变式压力计受温度影响和滞后影响,而石英压力计会受到温度和压力急剧变化的影响。在压力监测时,这些传感器还涉及安装困难、长期稳定性差等问题。井下光纤传感器没有井下电子线路、易于安装、体积小、抗干扰能力强等优点,而这些正是井下监测所必需的。

美国CiDRA公司的在光纤压力监测研究方面处于前沿,他们的科研人员发现了布喇格光纤光栅传感器对压力的线性响应。已开发的传感器能够工作到175℃,200oC和稍高温度的产品正在开发,250℃是研发的下一个目标。不同温度和压力下的压力测量误差,在测试范围(0MPa~34.5MPa)内,均小于±6.89kPa,相当于电子测量系统的最好的水平。目前,CIDRA公司的光纤压力传感器的指标为:测程0~103MPa,过压极限129MPa,准确度±41.3kPa,分辨率2.06kPa,长期稳定性±34.5kPa/yr(连续保持150℃),工作温度范围25℃~175℃。1999年该公司在加利福尼亚的Baker油田进行了压力监测系统的试验,结果表明该系统具有非常高的精度,目前已经交付商业销售。2001年该公司的压力传感器在英国BP公司的几口井下安装,监测应力变化,结果表明其具有足够的可靠性。

美国斯伦贝谢油田服务公司Doll研究中心的TsutomuYamate等人对用布喇格光纤光栅传感器实行井下监测进行了长期的研究,他们研制成一种对温度不敏感的侧孔布喇格光纤光栅传感器,最高工作温度为300℃,最高测量压力82MPa,在最高测量压力下,对温度的灵敏度极小,可以适用于井下的压力监测。

(2)温度监测

分布式光纤温度传感器具有通过沿整个完井长度连续性采集温度资料来提供一条监测生产和油层的新途径的潜力。因为井的温度剖面的变化可以与其它地面采集的资料(流量、含水、井口压力等)以及裸眼测井曲线对比,从而为操作者提供有关出现在井下的变化的定性和定量信息。传统的测温工具只能在任何给定时间内测量某个点的温度,要测试全范围的温度,点式传感器只能在井中来回移动才能实现,不可避免地对井内环境平衡造成影响。光纤分布式温度传感器的优势在于光纤无须在检测区域内来回移动,能保证井内的温度平衡状态不受影响。而且由于光纤被置于毛细钢管内,因此凡毛细钢管能通达的地方都可进行光纤分布式温度传感器测试。

最广泛地应用于井下监测应用的光纤传感器之一就是喇曼反向散射分布式温度探测器,这种方法已经在测量井筒温度剖面(特别是在蒸汽驱井)中,得到了广泛的应用。分布式温度传感器要综合考虑测量的点数和连接器衰减,遇到的问题和解决方法为:

(1)光纤以及连接器对信号的衰减问题,解决的方法为尽量减少连接器的数目、采用布喇格光纤光栅传感器以及改进连接器的性能;

(2)井下安装时容易损坏,解决的方法为配备熟练工人、光纤传感器需要外部保护层、减小应力(包括射孔和温度引起的应力)。对于光纤分布式温度传感器系统,英国Sensa公司一直处于技术领先地位,有一系列产品问世,而且与各大石油公司合作,积极探索光纤分布式温度传感器在石油井下的应用。CiDRA公司也一直在研究光纤温度传感器,目前该公司的温度传感器技术指标为:测量范围0℃~175℃,准确度±1℃,分辨率0.1℃,长期稳定性±1℃/yr(150℃下连续使用)。

目前的光纤温度、压力传感器的最主要的缺点之一就是温度压力交叉敏感特性,如何消除或者利用这种交叉敏感特性是研究的热点。

(3)多相流监测

为了做好油藏监控和油田管理,最关键的环节是获得生产井和注水井稳定可信的总流量剖面和各相流体的持率。然而,大多数油井分层开采,每层含水量不同,而且有时流速较大,给利用常规生产测井设备测量和分析油井的生产状况带来了巨大的困难。液体在油管中的摩阻和从油藏中向井筒内的喷射使得压差密度仪器无法准确测量,电子探头更是无法探测到液体中的小油气泡。

光纤测量多相流有两种方法,第一种是美国斯伦贝谢公司的持气率光纤传感仪,该仪器能直接测量多相流中持气率。其四个光纤探头均匀地分布在井筒的横剖面中,其空间取向方位可用一个集成化的相对方位传感器准确测量,在气液混合物中,通过探头反射的光信号来确定持气率和泡沫数量(这二者与气体流量相关联)。此外,利用每个探头的测量值来建立一种井中气体流动的图像,这些图像资料特别适用于斜井和水平井,可以更好地了解多相流流型以及解释在倾斜条件下这些流型固有的相分离。最近,这种仪器已在世界各地成功地进行了测井实验。它提供的资料能直接测定和量化多相混合物中气体和液体,能准确诊断井眼问题,并有助于生产调整。仪器通过了三口井的现场测试。第二种是通过测量声速来确定两相混合流的相组分,因为混合流体的声速与各单相流体的声速和密度具有相关性,而这个相关性普遍存在于两相气/液和液/液混合流体系统中,同时也适用于多相混合流系统。

根据混合流体的声速确定各相流体的体积分数,就是测量流过流量计的各单相体积分数(即持率测量)。某一流体相持率是否等于该相流动体积分数,取决于该相相对于其它相是否存在严重的滑脱现象。对于不存在严重滑脱的油水两相混合流系统,可以用均匀流动模型进行分析;对于存在严重滑脱现象的流动状态,则必须应用更完善的滑脱模型来解释流量计测量的数据,才能准确地确定各相的流量。经流动循环实验表明:对于油水混合流体,流量计的长波长声速测量可以确定各相体积分数(即持率),而不受流动非均质性(如层状流动)的影响。

CiDRA公司挖掘了光纤传感器内在的优势,开发了井下光相多相流传感器。目前的样品只局限在测量准均匀流体:如油、水两相或油、水、气三相(气相体积份数小于20%)。为了考察这种新型的光纤多相流传感器在生产井中测量油/水/气三相的性能,CiDRA最近在一口测试井进行了实验。在测试井中混合了油、水和气体,混合物包括粘度为32API的油、7%矿化度的水和矿厂天然气(甲烷),测试温度100oF,压力<2.75MPa。在0%~100%含水率范围内,仪器测量误差小于±5%,精度满足要求。该流量计能够确定原油和盐水混合物中的持水率,在持水率全量程中其误差为±5%以内,满足生产要求。而且除了能够测量持水率之外,该仪器还测试了三相中气体的体积含量,只是测试中油水的比例已知。结果表明,该仪器能够求出以泡沫流流出型出现的液体中的气体体积百分数。

2、声波测量

与过去相比,勘探开发公司如今面临更大的风险和更复杂的钻井环境,因此获得准确的地层构造图和油藏机理具有重要意义。目前使用的地震测量方法,如拖曳等浮电缆检波器组、临时海底布放地震检波器和井下电缆布放地震检波器等,能提供目的产油区域的测量,但这些方法具有相对高的作业费用,不能下入井内或受环境条件的限制等,而且提供的图像不全面、不连续,分辨率不是很高,因此难于实现连续实时油藏动态监测。

基于光纤的井下地震检波器系统能够解决这些问题,它能提供整个油井寿命期间永久高分辨率四维油藏图像,极大方便了油藏管理。这种井下地震加速度检波器能接收地震波,并将其处理成地层和流体前缘图像。

永久井下光纤3分量地震测量具有高的灵敏度和方向性,能产生高精度空间图像,不仅能提供近井眼图像,而且能提供井眼周围地层图像,在某些情况下测量范围能达数千英尺。它在油井的整个寿命期间运行,能经受恶劣的环境条件(温度达175℃,压力达100MPa),且没有可移动部件和井下电子器件,被封装在直径2.5cm的保护外壳中,能经受强的冲击和振动,可安装到复杂的完井管柱及小的空间内。此外,该系统还具有动态范围大和信号频带宽的特点,其信号频带宽度为3Hz~800Hz,能记录从极低到极高频率的等效响应。

  • 3激光光纤核测井技术

激光技术和光纤技术可以用于研制井下传感器,用于在充有原油和泥浆等非透明流体的井中进行测井。对于激光光纤核传感器的研究在国外比较盛行,美国、德国、俄罗斯和比利时等国均有大量的有关研究论文。

激光光纤核传感器是在光纤通信和光纤传感器的基础上产生的,它利用了光致损耗和光致发光等物理效应,比常规核探测器具有更多的优越性,是典型的学科交叉。光纤核测井技术,实际上就是在特定的环境下的核探测技术,其典型的优点为:

(1)可以针对不同的核探测的能级范围,研制在该范围的敏感探头。

(2)因为应用了光致发光效应,可使探头位于千米的井下,而光电倍增管由传输光缆相连置于井上,远离了恶劣的井下环境(高温高压),从而延长其的使用寿命。

(3)光纤具有高速率、大容量传输能力,还能搭载其他井下仪器信号。

然而,激光光纤核探测器也有缺点,主要表现在耐高温和承受高压的保护涂层、传输光缆的机械强度以及耐辐射的传输光缆低衰减损耗。

三、结论与展望

从本文的分析可以看出,光纤传感器以其独特的优势,可以广泛应用于石油天然气井下的储层参数监测(包括温度、压力和多相流)、声波检测和激光光纤核测井之中,极大地丰富了石油和天然气公司对储层的了解,便于优化油气田开采和维护。值得一提的是,该系统能够及时得到注采的注水压力和温度,从而判断压力是否超标,从而预防由于压力超标导致的套管损坏,这是一个全新的领域,国内外关于此方面尚未有报道和介绍。

到目前为止,全世界各大石油生产和服务公司都投入了巨资来研究和开发光纤传感器在储层评价中的应用,还有相当多的光纤传感器研发机构也致力于这一新兴领域的工作。可以设想,下一代光纤传感器在克服自身的缺点和劣势以后,将大面积推广,能更有效地帮助实时了解油气开采动态的水平。各大油田公司能够充分利用这些有利的信息,实现和维持油田的最优化生产,从而使油藏达到最高的采收率。同时,由于因特网的飞速发展,光纤监测的井况参数可以及时传递,这使得石油行业相关的生产和服务公司能够更有效地分析和评价全世界的资产。

本文为原创文章,转载请注明出处!

声波反射成像

声反射成像测井是近年来发展起来的新型声波测井方法, 通
过探测从井旁裂缝或小构造反射回来的声场, 来研究井旁介质的构
造特征。 现有的处理方法大多都直接沿用地震勘探的信号处理方
法, 因此在应用于声反射信号处理时存在不少局限性。 本文主要研
究参数估计法和径向道滤波法, 通过这两种方法来增强来自井旁裂
缝的反射信号。
通过对基于声波阵列中不同振型的慢度不同而导出的参数估
计法进行研究, 发现此波场分离方法对共接收器阵列的分离结果不
是很理想。 在对反射波慢度进行分析后, 引入了发射器阵列数据,
提出了改进的波场分离方法。 对实测声波数据应用改进的参数估计
法, 成功的压制了井眼模式波。ng

本文为原创文章,转载请注明出处!