-
Notifications
You must be signed in to change notification settings - Fork 34
/
Copy pathgrid3drcfs.m
111 lines (105 loc) · 3.75 KB
/
grid3drcfs.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
%GRID3DRCFS class to perform raytracing in 3D with the fast sweeping method
%
% Usage:
%
% Create and destroy instance of class
%
% g = grid3drcfs(par, nthreads)
% clear g
%
% Input for instantiation
% par: struct variable for grid definition with the following fields
% xmin: origin in X
% ymin: origin in Y
% zmin: origin in Z
% dx: cell size in X
% dy: cell size in Y
% dz: cell size in Z
% nx: number of cells in X
% ny: number of cells in Y
% nz: number of cells in Z
% nthreads: number of threads (optional, default = 1)
%
% Raytracing
% [tt] = g.raytrace(s, Tx, Rx, t0)
% [tt, rays] = g.raytrace(s, Tx, Rx, t0)
% [tt, rays, L] = g.raytrace(s, Tx, Rx, t0)
%
% Input
% g: grid instance
% s: slowness vector ( nSlowness by 1 )
% Tx: source coordinates, nTx by 3
% 1st column contains X coordinates, 2nd contains Y coordinates,
% 3rd contains Z coordinates
% Rx: receiver coordinates, nRx by 3
% 1st column contains X coordinates, 2nd contains Y coordinates,
% 3rd contains Z coordinates
% t0: source epoch, nTx by 1
% t0 is optional (0 if not given)
%
% *** IMPORTANT: Tx or Rx should _not_ lie on (or close to) an external
% face of the grid when rays are needed ***
% *** nTx must be equal to nRx, i.e. each row define one Tx-Rx pair ***
% *** nSlowness must equal g.nx*g.ny*g.nz ***
% *** Indexing of slowness values is done by "vectorizing" a 3D array,
% i.e. if slowness field s is of size (nx,ny,nz), enter s(:) as
% first argument
%
%
% Output
% tt: vector of traveltimes, nRx by 1
% rays: cell object containing the matrices of coordinates of the ray
% paths, nRx by 1. Each matrix is nPts by 3
% L: data kernel matrix (tt = L*s)
%
% -----------
%
% Bernard Giroux
% INRS-ETE
% 2015-03-04
classdef grid3drcfs < handle
properties (SetAccess = private, Hidden = true)
objectHandle; % Handle to the underlying C++ class instance
end
methods
% Constructor - Create a new C++ class instance
function this = grid3drcfs(varargin)
this.objectHandle = grid3drcfs_mex('new', varargin{:});
end
% Destructor - Destroy the C++ class instance
function delete(this)
grid3drcfs_mex('delete', this.objectHandle);
end
% setSlowness
function varargout = setSlowness(this, varargin)
[varargout{1:nargout}] = grid3drcfs_mex('setSlowness', this.objectHandle, varargin{:});
end
% raytrace
function varargout = raytrace(this, varargin)
[varargout{1:nargout}] = grid3drcfs_mex('raytrace', this.objectHandle, varargin{:});
end
% for saving in mat-files
function s = saveobj(obj)
s.xmin = grid2drcfs_mex('get_xmin', obj.objectHandle);
s.ymin = grid2drcfs_mex('get_ymin', obj.objectHandle);
s.zmin = grid2drcfs_mex('get_zmin', obj.objectHandle);
s.dx = grid2drcfs_mex('get_dx', obj.objectHandle);
s.dy = grid2drcfs_mex('get_dy', obj.objectHandle);
s.dz = grid2drcfs_mex('get_dz', obj.objectHandle);
s.nx = grid2drcfs_mex('get_nx', obj.objectHandle);
s.ny = grid2drcfs_mex('get_ny', obj.objectHandle);
s.nz = grid2drcfs_mex('get_nz', obj.objectHandle);
s.nthreads = grid2drcfs_mex('get_nthreads', obj.objectHandle);
end
end
methods(Static)
% for loading from mat-files
function obj = loadobj(s)
if isstruct(s)
obj = grid2drcfs(s, s.nthreads);
else
error('Wrong input arguments')
end
end
end
end