diff --git a/marchenko3D/demo/marchenko3D/oneD/makeR.sh b/marchenko3D/demo/marchenko3D/oneD/makeR.sh
index e31cf348249a47faec8c10db9635619c0c2a7c69..5476c8335da4e2137cd88d54d99e98ec6a5f7419 100755
--- a/marchenko3D/demo/marchenko3D/oneD/makeR.sh
+++ b/marchenko3D/demo/marchenko3D/oneD/makeR.sh
@@ -1,5 +1,11 @@
#!/bin/bash
+echo "WARNING! The data size of the 3D reflection data is very large!!!"
+echo "When you are positive your machine has enough memory,"
+echo "delete the exit; statement from makeR.sh"
+
+exit;
+
#create the reflection data that was modeled in 1 1D medium
makeR1D file_in=shotx10y10.su file_out=reflx10y10.su verbose=2 nxrcv=201 nyrcv=61 nxsrc=201 nysrc=61
diff --git a/marchenko3D/marchenko3D.c b/marchenko3D/marchenko3D.c
index 1ed3108a5d1df175e587b4aa0277df95da7a7d8f..4ed77c7e5d14e15ef267d351cdb5312deddfbbd5 100644
--- a/marchenko3D/marchenko3D.c
+++ b/marchenko3D/marchenko3D.c
@@ -196,7 +196,7 @@ int main (int argc, char **argv)
float d1, d2, d3, f1, f2, f3, fxsb, fxse, fysb, fyse, ft, fx, fy, *xsyn, *ysyn, dxsrc, dysrc;
float *green, *f2p, *G_d, dt, dx, dy, dxs, dys, scl, mem;
float *f1plus, *f1min, *iRN, *Ni, *trace, *Gmin, *Gplus, *HomG;
- float scale, *tmpdata, tplmax;
+ float scale, *tmpdata, tplmax, tshift;
float *ixmask, *iymask, *ampscl, *Gd, *Image, dzim;
float grad2rad, px, py, src_anglex, src_angley, src_velox, src_veloy, *mutetest;
complex *Refl, *Fop;
@@ -393,6 +393,12 @@ int main (int argc, char **argv)
/* compute time shift for tilted plane waves */
if (plane_wave==1) {
+ grad2rad = 17.453292e-3;
+ px = sin(src_anglex*grad2rad)/src_velox;
+ py = sin(src_angley*grad2rad)/src_veloy;
+
+ tshift = fabs((nys-1)*dys*py) + fabs((nxs-1)*dxs*px);
+
for (j = 0; j < Nfoc; j++) {
itmin[j] = nt;
for (i=0; i<nys*nxs; i++) itmin[j] = MIN (itmin[j], muteW[j*nxs*nys+i]);
@@ -415,13 +421,13 @@ int main (int argc, char **argv)
for (j = ntapx; j < nxs-ntapx; j++)
tapersx[j] = 1.0;
for (j = nxs-ntapx; j < nxs; j++)
- tapersx[j] = (cos(PI*(j-(nxs-ntapx))/ntapx)+1)/2.0;
+ tapersx[j] = (cos(PI*(j+1-(nxs-ntapx))/ntapx)+1)/2.0;
for (j = 0; j < ntapy; j++)
tapersy[j] = (cos(PI*(j-ntapy)/ntapy)+1)/2.0;
for (j = ntapy; j < nys-ntapy; j++)
tapersy[j] = 1.0;
for (j = nys-ntapy; j < nys; j++)
- tapersy[j] = (cos(PI*(j-(nys-ntapy))/ntapy)+1)/2.0;
+ tapersy[j] = (cos(PI*(j+1-(nys-ntapy))/ntapy)+1)/2.0;
}
else {
for (j = 0; j < nxs; j++) tapersx[j] = 1.0;
@@ -883,6 +889,7 @@ int main (int argc, char **argv)
}
}
}
+ timeShift(Gmin, nts, npos*Nfoc, dt, tshift, fmin, fmax);
}
else {
applyMute3D(Gmin, muteW, smooth, 4, Nfoc, nxs, nys, nts, ixpos, iypos, npos, shift, tsynW);
diff --git a/raytime3d/Makefile b/raytime3d/Makefile
index 1125552ed09f052253f4cc22fc24ab3cc28949a7..47996e3c982f26a7e8a1c65c541f285cf2ce6f98 100644
--- a/raytime3d/Makefile
+++ b/raytime3d/Makefile
@@ -41,7 +41,8 @@ SRCC = $(PRG).c \
atopkge.c \
docpkge.c \
threadAffinity.c \
- getpars.c
+ getpars.c \
+ writeData3D.c
OBJC = $(SRCC:%.c=%.o)
diff --git a/raytime3d/raytime3d.c b/raytime3d/raytime3d.c
index f93bb963090d3480524bcd46d24524ae2b5f6f7e..a689661c756702dbe87911eb696e6f78aed07060 100644
--- a/raytime3d/raytime3d.c
+++ b/raytime3d/raytime3d.c
@@ -35,6 +35,8 @@ void vidale3d(float *slow0, float *time0, long nz, long nx, long ny, float h, lo
void src3D(float *time0, float *slow0, long nz, long nx, long ny, float h, float ox, float oy, float oz, long xs, long ys, long zs, long *cube);
+int writeData3D(FILE *fp, float *data, segy *hdrs, long n1, long n2);
+
/*********************** self documentation **********************/
char *sdoc[] = {
@@ -287,45 +289,45 @@ private (is,time0,ampl,nrx,nry,nrz,nr,cp_average,i,j,k,ix,iy,iz,hdrs,tmpdata,nwr
if (verbose) vmess("Writing src %li of %li sources",is+1,shot.n);
if (verbose>1) vmess("xsrc[%li]=%f ysrc[%li]=%f zsrc[%li]=%f",shot.x[is],shot.xs[is],shot.y[is],shot.ys[is],shot.z[is],shot.zs[is]);
- hdrs = (segy *) calloc(1,sizeof(segy));
- tmpdata = (float *)malloc((rec.nx)*sizeof(float));
+ hdrs = (segy *) calloc(rec.nz*rec.ny,sizeof(segy));
+ tmpdata = (float *)malloc((rec.n)*sizeof(float));
+
+ for (j = 0; j < rec.nz; j++) {
+ for (i = 0; i < rec.ny; i++) {
+ hdrs[j*rec.ny+i].fldr = is+1;
+ hdrs[j*rec.ny+i].tracl = i+1;
+ hdrs[j*rec.ny+i].tracf = i+1;
+ hdrs[j*rec.ny+i].scalco = -1000;
+ hdrs[j*rec.ny+i].scalel = -1000;
+ hdrs[j*rec.ny+i].sx = (long)(f2+(shot.x[is]-1)*d2)*1000;
+ hdrs[j*rec.ny+i].sy = (long)(f3+(shot.y[is]-1)*d3)*1000;
+ hdrs[j*rec.ny+i].sdepth = (long)(f1+(shot.z[is]-1)*d1)*1000;
+ hdrs[j*rec.ny+i].selev = -(long)(f1+(shot.z[is]-1)*d1)*1000;
+ hdrs[j*rec.ny+i].gy = (long)(mod.y0+rec.yr[j*rec.ny*rec.nx+i*rec.nx])*1000;
+ hdrs[j*rec.ny+i].gx = (long)(mod.z0+rec.zr[j*rec.ny*rec.nx])*1000;
+ hdrs[j*rec.ny+i].ns = rec.nx;
+ hdrs[j*rec.ny+i].ntr = rec.ny*shot.n;
+ hdrs[j*rec.ny+i].trwf = rec.ny*shot.n;
+ hdrs[j*rec.ny+i].f1 = mod.x0+rec.xr[0];
+ hdrs[j*rec.ny+i].f2 = mod.y0+rec.yr[0];
+ hdrs[j*rec.ny+i].dt = (long)(d2*1e6);
+ hdrs[j*rec.ny+i].d1 = rec.xr[1]-rec.xr[0];
+ hdrs[j*rec.ny+i].d2 = rec.yr[rec.nx]-rec.yr[0];
+
+ for (k = 0; k < rec.nx; k++) {
+ tmpdata[j*rec.ny*rec.nx+i*rec.nx+k] = time0[(rec.z[j*rec.ny*rec.nx+i*rec.nx+k]+1)*nxy+(rec.y[j*rec.ny*rec.nx+i*rec.nx+k]+1)*(nx+2)+rec.x[j*rec.ny*rec.nx+i*rec.nx+k]+1];
+ }
+ }
+ }
+
+ ret = writeData3D(fpt, &tmpdata[0], hdrs, rec.nx, rec.nz*rec.ny);
+ if (ret < 0 ) verr("error on writing output file.");
+
+ if (ray.geomspread==1) {
+ ret = writeData3D(fpt, &l[0], hdrs, rec.nx, rec.nz*rec.ny);
+ if (ret < 0 ) verr("error on writing output file.");
+ }
- for (i = 0; i < rec.ny; i++) {
- hdrs[0].fldr = is+1;
- hdrs[0].tracl = i+1;
- hdrs[0].tracf = i+1;
- hdrs[0].scalco = -1000;
- hdrs[0].scalel = -1000;
- hdrs[0].sx = (long)(f2+(shot.x[is]-1)*d2)*1000;
- hdrs[0].sy = (long)(f3+(shot.y[is]-1)*d3)*1000;
- hdrs[0].sdepth = (long)(f1+(shot.z[is]-1)*d1)*1000;
- hdrs[0].selev = -(long)(f1+(shot.z[is]-1)*d1)*1000;
- hdrs[0].gy = (long)(mod.y0+rec.yr[i*rec.nx])*1000;
- hdrs[0].ns = rec.nx;
- hdrs[0].ntr = rec.ny*shot.n;
- hdrs[0].trwf = rec.ny*shot.n;
- hdrs[0].f1 = mod.x0+rec.xr[0];
- hdrs[0].f2 = mod.y0+rec.yr[0];
- hdrs[0].dt = (long)(d2*1e6);
- hdrs[0].d1 = rec.xr[1]-rec.xr[0];
- hdrs[0].d2 = rec.yr[rec.nx]-rec.yr[0];
-
- for (k = 0; k < rec.nx; k++) {
- tmpdata[k] = time0[(rec.z[i*rec.nx+k]+1)*nxy+(rec.y[i*rec.nx+k]+1)*(nx+2)+rec.x[i*rec.nx+k]+1];
- }
-
- nwrite = fwrite( &hdrs[0], 1, TRCBYTES, fpt);
- assert(nwrite == TRCBYTES);
- nwrite = fwrite( tmpdata, sizeof(float), rec.nx, fpt );
- assert(nwrite == rec.nx);
-
- if (ray.geomspread==1) {
- nwrite = fwrite( &hdrs[0], 1, TRCBYTES, fpa);
- assert(nwrite == TRCBYTES);
- nwrite = fwrite( &l[i*rec.nx], sizeof(float), rec.nx, fpa );
- assert(nwrite == rec.nx);
- }
- }
writer++;
free(time0);
free(hdrs);
diff --git a/utils/HomG.c b/utils/HomG.c
index 676459e976e40a95a7351748cb1bd5750c7ca1cd..52cf1864bc34b440f49e4c60e3b92129a79d9a98 100755
--- a/utils/HomG.c
+++ b/utils/HomG.c
@@ -559,14 +559,16 @@ int main (int argc, char **argv)
nzvr = nyvr;
nyvr = ix;
tmp1 = (float *)calloc(nxvr*nyvr*nzvr,sizeof(float));
- for (it=0; it<ntr; it++) {
- for (iy=0; iy<nyvr*nxvr*nzvr; iy++) {
- tmp1[iy] = Ghom[it*nyvr*nxvr*nzvr+iy];
- }
- for (iy=0; iy<nyvr; iy++) {
- for (ix=0; ix<nxvr; ix++) {
- for (ir=0; ir<nzvr; ir++) {
- Ghom[it*nyvr*nxvr*nzvr+iy*nxvr*nzvr+ix*nzvr+ir] = tmp1[ir*nxvr*nyvr+ix*nyvr+iy];
+ for (isn=0; isn<nshots; isn++) {
+ for (it=0; it<ntr; it++) {
+ for (iy=0; iy<nyvr*nxvr*nzvr; iy++) {
+ tmp1[iy] = Ghom[isn*ntr*nyvr*nxvr*nzvr+it*nyvr*nxvr*nzvr+iy];
+ }
+ for (iy=0; iy<nyvr; iy++) {
+ for (ix=0; ix<nxvr; ix++) {
+ for (ir=0; ir<nzvr; ir++) {
+ Ghom[isn*ntr*nyvr*nxvr*nzvr+it*nyvr*nxvr*nzvr+iy*nxvr*nzvr+ix*nzvr+ir] = tmp1[ir*nxvr*nyvr+ix*nyvr+iy];
+ }
}
}
}
@@ -584,14 +586,16 @@ int main (int argc, char **argv)
nzvr = nxvr;
nxvr = ix;
tmp1 = (float *)calloc(nxvr*nyvr*nzvr,sizeof(float));
- for (it=0; it<ntr; it++) {
- for (iy=0; iy<nyvr*nxvr*nzvr; iy++) {
- tmp1[iy] = Ghom[it*nyvr*nxvr*nzvr+iy];
- }
- for (iy=0; iy<nyvr; iy++) {
- for (ix=0; ix<nxvr; ix++) {
- for (ir=0; ir<nzvr; ir++) {
- Ghom[it*nyvr*nxvr*nzvr+iy*nxvr*nzvr+ix*nzvr+ir] = tmp1[iy*nzvr*nxvr+ir*nxvr+ix];
+ for (isn=0; isn<nshots; isn++) {
+ for (it=0; it<ntr; it++) {
+ for (iy=0; iy<nyvr*nxvr*nzvr; iy++) {
+ tmp1[iy] = Ghom[isn*ntr*nyvr*nxvr*nzvr+it*nyvr*nxvr*nzvr+iy];
+ }
+ for (iy=0; iy<nyvr; iy++) {
+ for (ix=0; ix<nxvr; ix++) {
+ for (ir=0; ir<nzvr; ir++) {
+ Ghom[isn*ntr*nyvr*nxvr*nzvr+it*nyvr*nxvr*nzvr+iy*nxvr*nzvr+ix*nzvr+ir] = tmp1[iy*nzvr*nxvr+ir*nxvr+ix];
+ }
}
}
}
diff --git a/utils/combine_induced.c b/utils/combine_induced.c
index bdb37148f5f4fe00b96cff1597113798253685bc..b1c63a8679c6c1bbb2a847fd255df9cd1ab9c43a 100755
--- a/utils/combine_induced.c
+++ b/utils/combine_induced.c
@@ -41,21 +41,25 @@ char *sdoc[] = {
" ",
" Optional parameters: ",
" ",
-" file_out= ................ Filename of the output",
-" numb= .................... integer number of first file",
-" dnumb= ................... integer number of increment in files",
-" nzmax= ................... Maximum number of files read",
+" file_out=out.su .......... Filename of the output",
+" numb=1 ................... integer number of first file",
+" dnumb=1 .................. integer number of increment in files",
+" nzmax=0 .................. Maximum number of files read, 0 is no limit",
+" nshift=0 ................. Sample number shift",
+" shift=0.0 ................ Base shift in seconds for the data",
+" dtshift=0.0 .............. Shift per gather in seconds",
+" verbose=1 ................ Give information about the process (=1) or not (=0)",
NULL};
int main (int argc, char **argv)
{
FILE *fp_in, *fp_out;
char *fin, *fout, *ptr, fbegin[100], fend[100], fins[100], fin2[100], numb1[100];
- float *indata, *outdata, shift, dtshift, dt_time;
+ float *indata, *outdata, *loopdata, shift, dtshift;
float dt, dz, dy, dx, t0, x0, y0, z0, scl, dxrcv, dyrcv, dzrcv;
long nt, nz, ny, nx, nxyz, ntr, ix, iy, it, is, iz, pos, file_det, nxs, nys, nzs;
long numb, dnumb, ret, nzmax, verbose, nt_out, ishift, nshift, *sx, *sy, *sz;
- segy *hdr_in, *hdr_bin, *hdr_out;
+ segy *hdr_in, *hdr_loop, *hdr_out;
initargs(argc, argv);
requestdoc(1);
@@ -72,7 +76,6 @@ int main (int argc, char **argv)
if (!getparlong("nshift", &nshift)) nshift=0;
if (!getparfloat("shift", &shift)) shift=0.0;
if (!getparfloat("dtshift", &dtshift)) dtshift=0.0;
- if (!getparfloat("dt_time", &dt_time)) dt_time=0.004;
if (fin == NULL) verr("Incorrect downgoing input");
/*----------------------------------------------------------------------------*
@@ -140,15 +143,15 @@ int main (int argc, char **argv)
* Read in a single file to determine if the header values match
* and allocate the data
*----------------------------------------------------------------------------*/
- hdr_in = (segy *)calloc(nx*ny,sizeof(segy));
+ hdr_in = (segy *)calloc(nx*ny*nt,sizeof(segy));
indata = (float *)calloc(nxyz*nt,sizeof(float));
readSnapData3D(fin2, indata, hdr_in, nt, nx, ny, nz, 0, nx, 0, ny, 0, nz);
dt = ((float)hdr_in[0].dt)/1E6;
- sx = (long *)calloc(nx*ny,sizeof(float));
- sy = (long *)calloc(nx*ny,sizeof(float));
- sz = (long *)calloc(nx*ny,sizeof(float));
+ sx = (long *)calloc(nx*ny,sizeof(long));
+ sy = (long *)calloc(nx*ny,sizeof(long));
+ sz = (long *)calloc(nx*ny,sizeof(long));
for (ix = 0; ix < nx*ny; ix++) {
sx[ix] = hdr_in[0].sx;
@@ -165,14 +168,13 @@ int main (int argc, char **argv)
/*----------------------------------------------------------------------------*
* Parallel loop for reading in and combining the various shots
*----------------------------------------------------------------------------*/
-#pragma omp parallel default(shared) \
- private(indata,hdr_in,fins,fin2,fp_in,is,ix,iy,iz,it,ishift)
-{
- indata = (float *)calloc(nxyz*nt,sizeof(float));
- hdr_in = (segy *)calloc(nx*ny*nt,sizeof(segy));
+#pragma omp parallel for schedule(static,1) default(shared) \
+ private(loopdata,hdr_loop,fins,fin2,fp_in,is,ix,iy,iz,it,ishift)
-#pragma omp for
for (is = 0; is < nzs; is++) {
+ loopdata = (float *)calloc(nxyz*nt,sizeof(float));
+ hdr_loop = (segy *)calloc(nx*ny*nt,sizeof(segy));
+
if (verbose) vmess("Depth:%li out of %li",is+1,nzs);
sprintf(fins,"%li",is*dnumb+numb);
sprintf(fin2,"%s%s%s",fbegin,fins,fend);
@@ -181,10 +183,10 @@ int main (int argc, char **argv)
verr("Error opening file");
}
fclose(fp_in);
- readSnapData3D(fin2, indata, hdr_in, nt, nx, ny, nz, 0, nx, 0, ny, 0, nz);
- sx[is] = hdr_in[0].sx;
- sy[is] = hdr_in[0].sy;
- sz[is] = hdr_in[0].sdepth;
+ readSnapData3D(fin2, loopdata, hdr_loop, nt, nx, ny, nz, 0, nx, 0, ny, 0, nz);
+ sx[is] = hdr_loop[0].sx;
+ sy[is] = hdr_loop[0].sy;
+ sz[is] = hdr_loop[0].sdepth;
ishift = nshift*is;
if (verbose) vmess("Shifting %li timesteps for a total of %.3f seconds",ishift,shift+(dtshift*((float)iz)));
@@ -192,14 +194,14 @@ int main (int argc, char **argv)
for (iy = 0; iy < ny; iy++) {
for (ix = 0; ix < nx; ix++) {
for (iz = 0; iz < nz; iz++) {
- outdata[it*nxyz+iy*nx*nz+ix*nz+iz] += indata[(it-ishift+(nt/2))*nxyz+iy*nx*nz+ix*nz+iz];
+ outdata[it*nxyz+iy*nx*nz+ix*nz+iz] += loopdata[(it-ishift+(nt/2))*nxyz+iy*nx*nz+ix*nz+iz];
}
}
}
}
+ free(loopdata);free(hdr_loop);
}
- free(indata);free(hdr_in);
-}
+
/*----------------------------------------------------------------------------*
* Write out the data to file
@@ -232,7 +234,7 @@ int main (int argc, char **argv)
hdr_out[iy*nx+ix].offset = (hdr_out[iy*nx+ix].gx - hdr_out[iy*nx+ix].sx)/1000.0;
}
}
- ret = writeData3D(fp_out, &outdata[is*nxs*nt], hdr_out, nt, nxs);
+ ret = writeData3D(fp_out, &outdata[it*ny*nx*nz], hdr_out, nz, nx*ny);
if (ret < 0 ) verr("error on writing output file.");
}
diff --git a/utils/getFileInfo3D.c b/utils/getFileInfo3D.c
index 7bb65617eab8082096da5b7f2dd4edb788eb7272..f62dd4fccf91b791a1aa1743a29539b7b3056dc6 100644
--- a/utils/getFileInfo3D.c
+++ b/utils/getFileInfo3D.c
@@ -53,7 +53,7 @@ long getFileInfo3D(char *filename, long *n1, long *n2, long *n3, long *ngath, fl
else scl = hdr.scalco;
*n1 = hdr.ns;
- if ( (hdr.trid == 1) && (hdr.dt != 0) ) {
+ if ( (hdr.trid == 101) && (hdr.dt != 0) ) {
*d1 = ((float) hdr.dt)*1.e-6;
*f1 = ((float) hdr.delrt)/1000.;
}
diff --git a/utils/pwshift.c b/utils/pwshift.c
index afbaa72e674f0d5cbe7e3644c5925a279d9ff73a..74ccaa1ba8f5590eeed6059a59862436de08251a 100644
--- a/utils/pwshift.c
+++ b/utils/pwshift.c
@@ -64,7 +64,7 @@ int main (int argc, char **argv)
FILE *fp_gmin, *fp_ray, *fp_amp, *fp_out;
char *file_gmin, *file_ray, *file_amp, *file_out, fbr[100], fer[100], fba[100], fea[100], fins[100], fin2[100], numb1[100], *ptr;
float *gmin, *conv, *time, *amp, *image, fmin, fmax;
- float dt, dy, dx, dz, t0, y0, x0, scl, *shift;
+ float dt, dy, dx, dz, t0, y0, x0, scl, *shift, *block;
float dt_ray, dy_ray, dx_ray, t0_ray, y0_ray, x0_ray, scl_ray, px, py, src_velox, src_veloy, src_anglex, src_angley, grad2rad;
long nshots, nt, ny, nx, ntr, delay;
long nray, nt_ray, ny_ray, nx_ray, ntr_ray;
@@ -185,6 +185,8 @@ int main (int argc, char **argv)
gy = (long *)calloc(nray*nshots,sizeof(long));
gz = (long *)calloc(nray*nshots,sizeof(long));
+ block = (float *)calloc(nray*nshots*nt,sizeof(float));
+
readSnapData3D(file_gmin, gmin, hdr_gmin, nshots, nx, ny, nt, 0, nx, 0, ny, 0, nt);
if (verbose) vmess("Read in Gmin data");
@@ -279,7 +281,10 @@ int main (int argc, char **argv)
}
timeShift(&conv[0],nt,nx*ny,dt,&time[is*ny*nx],&[is*ny*nx],shift,fmin,fmax);
for (ix = 0; ix < ny*nx; ix++) {
- image[is*nshots+iy] += conv[ix*nt+delay]*dx*dy*dt;
+ image[is*nshots+iy] += conv[ix*nt+delay+2]*dx*dy*dt;
+ for (it=0; it<nt; it++) {
+ block[is*nshots*nt+iy*nt+it] += conv[ix*nt+it];
+ }
}
}
@@ -319,7 +324,39 @@ int main (int argc, char **argv)
fclose(fp_out);
- free(image); free(hdr_out);
+ fp_out = fopen("block.su", "w+");
+
+ if (nshots>1) dz = ((float)(gz[1] - gz[0]))/1000.0;
+ else dz = 1.0;
+
+ for (is = 0; is < nshots; is++) {
+ for (it = 0; it < nray; it++) {
+ hdr_out[it*nshots+is].fldr = it+1;
+ hdr_out[it*nshots+is].tracl = it+1;
+ hdr_out[it*nshots+is].tracf = it+1;
+ hdr_out[it*nshots+is].scalco = -1000;
+ hdr_out[it*nshots+is].scalel = -1000;
+ hdr_out[it*nshots+is].trid = 1;
+ hdr_out[it*nshots+is].ns = nt;
+ hdr_out[it*nshots+is].trwf = nray;
+ hdr_out[it*nshots+is].ntr = nray;
+ hdr_out[it*nshots+is].f1 = (((float)gz[0])/1000.0);
+ hdr_out[it*nshots+is].f2 = (((float)gx[0])/1000.0);
+ hdr_out[it*nshots+is].dt = ((int)(dt*1E6));
+ hdr_out[it*nshots+is].d1 = roundf(dz*1000.0)/1000.0;
+ hdr_out[it*nshots+is].d2 = roundf(dx*1000.0)/1000.0;
+ hdr_out[it*nshots+is].gx = gx[it*nshots+is];
+ hdr_out[it*nshots+is].gy = gy[it*nshots+is];
+ hdr_out[it*nshots+is].sdepth = gz[it*nshots+is];
+ }
+ }
+
+ ret = writeData3D(fp_out, &block[0], hdr_out, nt, nray*nshots);
+ if (ret < 0 ) verr("error on writing output file.");
+
+ fclose(fp_out);
+
+ free(image); free(hdr_out); free(block);
vmess("Wrote data");