setpaths;
curdir = pwd;
expdir = [resultsdir '/rtm'];
if ~exist(expdir,'dir')
mkdir(expdir);
end
[v,o,d,n] = odnread([datadir '/marm_vp.odn']);
v = reshape(v,n);
v = [1500*ones(60,501); v(1:301,500 + [1:501])];
n = size(v);
[z,x] = odn2grid(o,d,n);
model.o = o;
model.d = d;
model.n = n;
model.nb = [50 50];
model.zsrc = 25;
model.xsrc = 0:100:5000;
model.zrec = 15;
model.xrec = 0:50:5000;
model.freq = [2.5:.5:20];
model.f0 = 10;
model.t0 = 0;
Q = speye(length(model.xsrc));
S = opKron(opSmooth(n(2),200),opSmooth(n(1),200));
m = 1e6./v(:).^2;
m0 = S*m;
dm = m - m0;
D = F(m,Q,model);
[D0,J0] = F(m0,Q,model);
dD = D - D0;
[od,dd,nd] = grid2odn(model.xrec,model.xsrc,model.freq);
Wf = oppKron2Lo(opDiag(1./model.freq),opDirac(nd(1)*nd(2)));
Wz = opKron(opDirac(n(2)),opDiag(sqrt(z)));
dmt = Wz*J0'*(Wf*dD);
odnwrite([expdir '/dmt.odn'],dmt,o,d,n);