1、GprMax中GprMax2D使用方法V1.3试验环境:操作系统:Windows7软件版本:MATLAB7.1&GprMaxV2.0参考资料:1UserGuideV2.pdf 一、GprMax2D软件使用1.1 直接运行.GprMaxV2.0Windows文件夹下GprMax2D.exe文件,会出现以下窗口(也能够在命令提醒符窗口输入命令运行) :1.2 输入文件名注意:要输入文件全路径;*.in文件只要出现任何语法错误或路径错误,软件全部会自动关闭,不会有任何错误提醒。出现以下画面(以自带例子文件bre1.in为例,*.in命令参考前面文章或1):运行完成会发觉.GprMaxV2.0Wind
2、ows文件夹下多了两个文件bre1.out、bre1.geo文件,复制到tools文件夹。 二、数据成像tools文件夹下有五个m文件:gprmax.m,gprmax2g.m, gprmax3g.m, gprmaxde.m , gprmaxso.m。这里只讲gprmax2g.m,gprmax.m这两个文件使用方法,其它三个以后有空再写。gprmax3g.m是处理GprMax3D几何数据;gprmaxde.m用来计算Debye公式(参考1)介电常数;gprmaxso.m用于计算激励函数。2.1 gprmax2g.m使用方法gprmax2g函数用于读取GprMax2D软件仿真探地雷达模型生成二进制
3、几何数据。gprmax2g函数原型:mesh,header,media = gprmax2g( filename )filename是.geo文件名;media: 存放介质类型,media.type;header: 存放模型几何参数;header.title: 模型名称;header.dx: 模型在X轴每次偏移大小(单位:m);header.dy: 模型在Y轴每次偏移大小(单位:m);header.dt: 最大许可时间步长(单位:秒);header.nx: 模型在X轴偏移次数;header.ny: 模型在Y轴偏移次数;例子:如输入文件*.in中定义: #domain: 2.5 0.65 #dx
4、_dy: 0.0025 0.0025那么:header.dx=0.0025; header.dy=0.0025;header.dt = 1/(c*sqrt(1/header.dx2+1/header.dy2); (其中c=,为光速,公式参考1);header.nx=2.5/0.0025=1000; header.ny=0.65/0.0025=260;mesh: 存放模型数据,M x N数组,其中M为Y轴方向Yee单元数目,N为X轴方向Yee单元数目;M=header.nx,N=header.ny;gprmax2g.m使用例子: filegeo = bre1.geo; meshdata,head
5、er,media=gprmax2g(filegeo); figure(1); MM,NN=size(meshdata); imagesc(1:NN)*header.dx,(1:MM)*header.dy,meshdata) axis(equal); xlabel(x(m); ylabel(y(m); 2.2 gprmax.m使用方法gprmax函数用于读取GprMax2D和GprMax3D软件仿真探地雷达模型生成二进制波形数据。gprmax函数原型:Header, Fields = gprmax( filename )filename是.out文件名;1) Header:该变量包含以下组员:H
6、eader.title: 模型名称;Header.iterations: 迭代次数;Header.dx: 在X轴每次偏移大小;Header.dy: 在Y轴每次偏移大小;Header.dt: 最大许可时间步长;Header.NSteps: 仿真次数;等于*.in文件中#analysis:命令第一个参数;例子:如输入文件*.in中定义: #domain: 2.5 0.65 #dx_dy: 0.0025 0.0025 #time_window: 12e-9 #analysis: 115 bre1.out b #tx: 0.0875 0.4525 MyLineSource 0.0 12e-9 #rx:
7、 0.1125 0.4525 #tx_steps: 0.02 0.0 #rx_steps: 0.02 0.0 #end_analysis: 那么,Header.iterations=ceil(Header.removed/Header.dt);Header.dx=0.0025;Header.dy=0.0025;Header.dt= 1/(c*sqrt(1/Header.dx2+1/Header.dy2);Header.NSteps=115;Header.TxStepX=0.02/Header.dx=8; Header.TxStepY=0.0/Header.dy=0;Header.RxStepX
8、=0.02/Header.dx=8; Header.RxStepY=0.0/Header.dy=0;Header.ntx: Header.nrx: Header.nrx_box: Header.tx=0.0875/Header.dx=35;Header.ty=0.4525/Header.dy=181;Header.source=MyLineSource;Header.delay=0;(等于#tx:命令第四个参数)Header.removed=12e-9;(等于#tx:命令第五个参数)Header.rx=0.1125/Header.dx=45;Header.ry=0.4525/Header.dy
9、=181; 2)Fields: 该变量包含以下组员:Fields.t: 每个波形时间轴。数组大小Header.iterations*1;Fields.ez: Z轴方向磁矢量数据。数组大小Header.iterations*1* Header.NSteps;Fields.hx: X轴方向电矢量数据。数组大小Header.iterations*1* Header.NSteps;Fields.hy: Y轴方向电矢量数据。数组大小Header.iterations*1* Header.NSteps;因为GprMax仿真是基于FDTD算法,以上三者关系为:具体可参考FDTD算法相关文件。 gprmax.
10、m使用例子: fileout = bre1.out; Header,Fields=gprmax(fileout); N=1:Header.NSteps; %移动次数 Position=Header.dx*Header.tx+(N-1)*(Header.dx*Header.TxStepX); %天线每次所在位置 Data(:,:)=Fields.ez(:,1,:); %转换数组格式 figure(2); imagesc(Position,Fields.t*1e9,Data);%画图 colorbar xlabel(Antena Position (m); ylabel(t(ns); GprMax
11、V2.0中GprMax2D输入文件命令试验环境:操作系统:Windows 7软件版本:MATLAB 7.1 & GprMaxV2.0 参考文件:GprMaxV2.0软件manual文件夹下UserGuideV2.pdf。GprMaxV2.0下载地址: 说明:翻译得不好,还望大家见谅,因为我也是边看边译。 3.1 GprMax2D命令通常注意事项为了描述GprMax2D/3D命令及其参数,我们作以下约定:f表示浮点数(如1.5或15e-1、1.5e1)i表示整数c表示字符str表示字符串file表示输入文件名 全部空间距离基础单位为米全部时间基础单位为秒全部频率参数基础单位是Hz 3.2 Gpr
12、Max2D 2.0版本共有32条命令:#title:#domain:#dx_dy:#time_step_stability_factor:#time_window:#messages: #number_of_media:#nips_number:#media_file:#geometry_file:#medium: #abc_type:#abc_order:#abc_stability_factors:#abc_optimization_angles:#abc_mixing_parameters:#pml_layers: #box:#cylinder:#x_segment:#y_segment
13、:#triangle: #analysis:#end_analysis:#tx:#rx:#rx_box:#snapshot:#tx_steps:#rx_steps: #line_source:#excitation_file: 以下能够运行于GprMax中命令在GprMax2D版本2.0中将不再支持:#scan:#csg:#extra_tx: 另外,部分命令旧参数规则亦发生了改变:#tx:#snapshot:3.3 GprMax2D命令参数为了愈加好地介绍GprMax2D命令,我们将它们分成四类:通常命令:包含用于指定大小和模型离散ABC相关命令: 许可定制和优化吸收边界条件介质和对象结构命令
14、: 用不一样参数来在模型中引入不一样介质和结构简单几何形状激励和输出命令:用来放置源代码和模型输出点 运行GprMax2D最低程度命令以下:#domain:#dx_dy:#time_window:最少一个#analysis:及和其对应结束命令#end_analysis:最少一个#tx:和#rx:,或#rx_box:命令 为了使#tx: 命令正确运行,同时需要一行新#line_source:命令 3.3.1通常命令#title:模型题目#title: strStr即是模型题目,必需是单行。 #domain:模型范围(单位:米)#domain: f1 f2f1和f2分别代表x和y轴上量度大小#do
15、main: 1.0 1.5表示1.0*1.5大小范围 #dx_dy:表示x和y轴上偏移量(如x、y)#dx_dy: f1 f2表示x轴偏移f1和y轴方向上偏移f2 #domain: 1.0 1.5#dx_dy: 0.1 0.1则模型单元数为10*15 #domain: 1.0 1.5#dx_dy: 0.01 0.01则模型单元数为100*150 最大许可时间步t和x、y约束关系以下:(3.1)C为光速,GprMax2D中计算t使用3.1式等号。#time step stability factor:经过该命令,你能够修改GprMax2D所计算t值,但必需满足3.1式要求。#time_windo
16、w:用于指定所需总模拟时间。语法:#time_window: f1或#time_window: i1总迭代次数和模拟时间窗口: (3.2)#number_of_media:但你需要使用大于10个介质时必需使用该命令,因为GprMax2D只初始化了10个介质使用空间。#number_of_media: i1i1大于10。 #media_file:使用该命令,你能够定义常见介质所描述结构参数文件路径和名称。#geometry_file:使用该命令,你能够定义模型几何信息二进制文件。这些信息能够用于创建模型图像和检验创建正确是否。#geometry_file: model.geo#messages:
17、使用该命令,你能够控制软件运行时在屏幕上输出信息。#messages:c1#nips_number:该命令仅当在GprMax2D需要时才用于你输入文件中。#nips_number:i14.1 GPRMAX3D命令通常注释大多数能够用在GPRMAX3D命令几乎和GPRMAX2D是一样。不过,有部分命令是GPRMAX3D独有。为了简明扼要,这里只叙述和GPRMAX2D不一样命令。基础空间和临时离散步骤分别为t和x,y,z。 4.2 命令清淡在版本2里面,GPRMAX3D一共有42个命令用在3D GPR建模中。她们是:#title:#domain:#dx_dy_dz:#time_step_stabi
18、lity_factor:#messages: #number_of_media:#nips_number:#media_file:#geometry_file:#medium: #abc_type:#abc_order:#abc_stability_factors:#abc_optimization_angles:#abc_mixing_parameters:#pml_layers: #box:#cylinder:#cylinder_new:#cylindrical_segment:#sphere:#plate:#edge:#triangle:#bowtie:#thin_wire: #anal
19、ysis:#end_analysis:#tx:#rx:#rx_box#snapshot:#tx_steps:#rx_steps: #huygens_surface:#hertzian_dipole:#voltage_source:#transmissio_line:#plane_wave:#excitation_file: 4.3 GPRMAX3D命令参数为了愈加好地介绍GprMax2D命令,我们将它们分成四类:1.通常命令:包含用于指定大小和模型离散2. ABC相关命令: 许可定制和优化吸收边界条件3.介质和对象结构命令: 用不一样参数来在模型中引入不一样介质和结构简单几何形状4.激励和输出
20、命令:用来放置源代码和模型输出点 4.3.1 通常命令#title:和GprMax2D使用方法形同#domain:用来指定模型尺寸(单位:米)#domain: f1 f2 f3f1,f2和f3分别是模型x,y,z方向上大小 #dx_dy_dz:指定x,y,z方向上增量#dx_dy_dz: f1 f2 f3f1是空间步x方向上增量,f2是空间步y方向上增量,f3是空间步z方向上增量。空间离散化控制最大许可时间步t和为了达成所需要仿真时间窗而提出处理方案。t和x,y,z之间关系是:(4.1)其中c为光速。在GPRMAX3D中,上式取等号。4.1中,一个小x,y,z值造成t得小值,这个t小值意味着为
21、了达成所给仿真时间而需要更多迭代次数。不过,需要指出更小x,y,z和t值会让模型更正确。 #time_step_stability_factor:使用方法和GPR2D相同。能够修改t值。#time_window:使用方法和GPR2D相同。#number_of_media:使用方法和GPR2D相同。#geometry_file:使用方法和GPR2D相同。#messages:使用方法和GPR2D相同。#nips_number:使用方法和GPR2D相同。仅当GPRMAX3D被要求用时,才必需用到她。 4.3.2 ABC 相关命令在GPRMAX 3D这些影响Higdon ABCs配置和性能命令和GPR
22、MAX 2D使用方法相同。不过,GPRMAX3D用了更有力PML ABC。#abc_type: pml#abc_type: higdon 能够设置ABC默认参数。 #pml_layers: i1i1是pml所占Yee单元数目。默认值为8。所占用Yee单元数越多,PML性能越好,不过所花计算资源也越多。而且,有两点要注意:a) PML是几何模型一部分。不过,PML层里域不参与计算而且用她们来计算就是错。所以,不要把源和接收器设置在这个区域。PML深度是用YEE 单元来度量而不是距离。b) 注意:目前PML仅用于非磁介质。所以,假如你模型需要用到磁性参数,你必需用Higdoon ABC而不是PML
23、。4.3.3 介质和对象结构命令在GPRMAX3D中,这些命令用来轻易引入模型中不一样介质和介质结构。在处理通常命令以后,GPRMAX3D建立一个初始化为free space(air)任意尺寸模型。 注意,free space和pec已经被定义在GPRMAX3D中,你无须再去定义这两种介质。所以,关键词free_space和pec能够直接用而不需要再定义参数。其它介质她们参数必需经过下面之一来设置:l #media_file:命令包含多种多样常常见介质定义(见附录A)l #medium: 输入文件中命令#medium:和GPRMAX2D中使用方法相同。而且,同一个介质文件能够不改变她们结构下既
24、用到GPRMAX2D又用到GPRMAX3D #thin_wire:引入细电线模型。细电线被用作介质标示符给#dege:对象结构命令。#thin_wire: f1 str1f1是电线半径,而且为了建立适宜细电线物理模型,必需必l小。Str1是介质ID。Thin wire被认为是pec。例:#thin_wire: 0.001 MyWire#box: f1 f2 f3 f4 f5 f6 str1f1 f2 f3是左下角坐标,f4 f5 f6是右上角坐标。Str1是#medium:定义标示符。 #plate:f1 f2 f3 f4 f5 f6 str1f1 f2 f3是金属板左下角坐标。f4 f5 f
25、6是金属板右上角坐标。str1 是介质标示符。 #triangle:三角行块。#triangle: f1 f2 f3 f4 f5 f6 f7 f8 f9 str1f1 f2 f3,f4 f5 f6 ,f7 f8 f9分别是三角形三个顶点坐标。str1是介质标示符。 #bowtie: 蝴蝶天线。由两个三角形块组成#bowtie:c1 c2 f1 f2 f3 f4 f5 str1C1是蝴蝶天线方向,x,y或z。c2是剩下方向。F1,f2,f3是天线馈电点x,y,z坐标。F4是天线元长度(完整蝴蝶天线长度二分之一)。f5是展开角。Str1是介质标示符。 #edge:一个天线。这个天线仅仅是YEE元边
26、缘,当建立电阻或细电线是能够用她。#edge:f1 f2 f3 f4 f5 f6 str1F1 f2 f3是边缘起始坐标,f4 f5 f6是边缘结束坐标。Str1是介质标示符。 #cylinder:有限维3D圆柱模型。#cylinder:c1 f1 f2 f3 f4 f5 str1C1是圆柱轴方向,能够是x,y或z。f1和f2分别是圆柱轴高低坐标。F3,f4是其它关键坐标用来表示两个圆柱圆形面中心。l X方向圆柱(f1,f3,f4)为(x1,y,z),(f2,f3,f4)为(x2,y,z)l y方向圆柱(f3,f1,f4)为(x,y1,z),(f3,f2,f4)为(x,y2,z)l z方向圆柱
27、(f3,f4,f1)为(x,y,z1),(f3,f4,f2)为(x,y,z2)f5是圆柱圆盘半径,str1是介质标示符。 #cylinder_new:有限维3D圆柱。和#cylinder:不一样是,它轴向能够是任意。#cylinder_new: f1 f2 f3 f4 f5 f6 f7 str1F1 f2 f3 是圆柱底面圆心坐标,f4 f5 f6是圆柱顶面圆心坐标。F7是圆柱底面半径。Str1是介质标示符。 #cylindrical_segment:圆柱一段。#cylindrical_segment:c1 f1 f2 f3 f4 f5 str1 c2 f6 f7C1是圆柱轴向,能够是x,y,
28、z。f1,f2是圆柱轴向上部和下部坐标。F3,f4是表示圆柱顶,底部其它两个坐标。F5是顶,底部半径。C2是片段改变方向。F6,f7是片段起始和结束点。#sphere:球体。#sphere: f1 f2 f3 f4 str1F1 f2 f3是球心坐标 f4是球半径。Str1是介质标示符4.3.4 激励和输出命令#excitation_file: 许可用户指定单个包含能够激励模型幅度值清单ASCII文件。这些值最少要和迭代次数相等。#excitation_file:str1Str1是ASCII文件名字。例:#excitation_file: mysource.dat #hertzian_dipo
29、le:定义最简单激励。#hertzian_dipole:f1 f2 str1 str2F1 f2分别是源波形幅度和频率。Str1是波形标示符。Str2是源标示符。例:#hertzian_dipole:1.0 600e6 ricker MyDipole和GPRMAX2D里面#line_source等价 #voltage_source:定义电压源。它引入一个电压器件位置,能够是一个硬源或一个内部集成电阻。#voltage_source:f1 f2 str1 f3 str2F1 f2是源波形振幅和频率。Str1是波形标示符。 F3是内部电阻R。例:#voltage_source:1.0 600e6
30、gaussian 50.0 MyVolt #transmission_line:定义能够刺激天线1D两线传输线参数。#transmission_line: f1 f2 str1 f3 f4 str2F1 f2是源波形振幅和频率,str1是激励类型。F3是传输线长度。F4是阻抗特征,str2是源标示符例:#transmission_line:1.0 600e6 gaussian 0.5 200.0 MyLine #plane_wave:描述平面波源#plane_wave:f1 f2 str1 f3 f4 f5 str2 str3F1 f2是源波形振幅和频率,str1是激励类型。F3 f4分别是0
31、到pi和0到2pi角度。 #huygens_surface:必需和平面波同时用。#huygens_surface:f1 f2 f3 f4 f5 f6 str1F1 f2 f3是huygens表面左下角坐标,f4 f5 f6是右上角坐标。Str1是要用到#plane_wave标示符例:#huygens_surface:0.2 0.2 0.2 0.8 0.8 0.8 MyHugens #tx:c1 f1 f2 f3 str1 f4 f5C1源极化方向,能够是x,y,z。f1 f2 f3是源坐标。Str1是源ID,f4是激励掩饰,f5是源移去时间。 #rx: f1 f2 f3F1 f2 f3是输出点
32、坐标 #rx_box:f1 f2 f3 f4 f5 f6 f7 f8 f9F1 f2 f3是输出体积左下坐标,f4 f5 f6是输出体积右上坐标,f7 f8 f9是定义每一方向输出点步数。最小值为f7为x,f8为y,f9为z。 #tx_steps: f1 f2 f3F1 f2 f3是x,y,y增量#rx_steps: f1 f2 f3同上 #snapshot:取得模型区域电磁场信息。#snapshot:i1 f1 f2 f3 f4 f5 f6 f7 f8 f9 f10 filed1 c1或#snapshot:i1 f1 f2 f3 f4 f5 f6 f7 f8 f9 i2 filed1 c1i
33、1是一个成为全局户激励点计数器。它决定哪一个snapshot被用激励源位置。她值在1和步数之间,这个步数定义为#analysis第一个参数f1 f2 f3是左下角坐标;f4 f5 f6是右上角坐标;f7f8f9是x,y,z方向上样本间隔。f10或i2是snapshot时间,或迭代次数。file1是存放snapshot文件名c1 能够是a或b表示snapshot格式是ASCII还是BINARYGprMax中GprMax3D使用方法很多网友问我GprMax 3D使用方法,没有那么多时间去一一叙述GprMax 3D原理。直接贴一个GprMax 3D例子上来,给大家参考吧。*注:这个例子也是很久之前参
34、考网上一个例子,起源应该是起源这里:(1)Second.in文件内容以下:#medium: 6.0 0.0 0.0 0.01 1.0 0.0 concrete-#domain: 0.6 1.2 0.65#dx_dy_dz: 0.01 0.01 0.01#time_window: 12e-9-#box: 0.0 0.0 0.0 0.6 1.2 0.5 concrete#box: 0.2 0.5 0.1 0.4 0.7 0.3 free_space-#hertzian_dipole: 1.0 900e6 ricker MyDipole#analysis: 21 second.out b#tx: x
35、 0.3 0.115 0.55 MyDipole 0.0 12e-9#rx: 0.3 0.365 0.55#tx_steps: 0.0 0.04 0.0#rx_steps: 0.0 0.04 0.0#end_analysis:-#geometry_file: second.geo#messages: y(2)GprMax 3D界面和运行方法跟GprMax 2D相同,这里不多说了。Matlab处理文件:name = second.geo;mesh,ID,header,media = gprmax3g(name);figurep1=patch(isosurface(mesh,1);set(p1,f
36、acecolor,1,0,0,edgecolor,none);p2=patch(isosurface(mesh,2);set(p2,facecolor,0,1,0,edgecolor,none);Header, Fields = gprmax(second.out);radargram(:,:) = Fields.ez(:,1,:);t=1e9*Fields.t;Y = Header.dy * Header.rx: Header.dy * Header.RxStepY :Header.dy * Header.rx + (Header.dy * Header.RxStepY * (Header.NSteps-1);figure;imagesc(Y,t,radargram);title(Ez);ylabel(Time (ns);xlabel(Scan length (m);figure;x,y=meshgrid(Y,t);surf(x,y,radargram,EdgeColor,none)view(1,1,0.5);title(Ez );ylabel(Time (m);xlabel(Scan length (m);处理结果图像:(见下页)(3)问题:geo文件显示不正确,没时间去找为何了,有空再弄吧。