ktmig.c,v

seismic software,very useful

      "              (s2 is usually =(ncdppl-1)*dds, l2 is usually =0.0)        \n"
      "   (x3,y3) has the smallest s value and the largest l value              \n"
      "              (s3 is usually =0.0, l3 is usually =(nlines-1)*ddl         \n"
      "   where                                                                 \n"
      "         ncdppl is the number of traces per line in the master grid      \n"
      "         dds is the trace spacing (within a line) in the master grid     \n"
      "         nlines is the number of lines in the master grid                \n"
      "         ddl is the line spacing in the master grid                      \n"
      "   s is the coordinate (in m or ft --- NOT just an integer number)\n"
      "              of trace position within a 3d line                         \n"
      "         l is the coordinate (in m or ft --- NOT just an integer number)\n"
      "              of line position within a 3d survey                        \n"
      " 8. When doing velocity analysis, the output locations are               \n"
      "       the output inline locations (s):                          \n"
      " sstart      (ss1),      ss1+dds,      ...,       ss1+(nds-1)*dds,       \n"
      " sstart+ds   (ss2),      ss2+dds,      ...,       ss2+(nds-1)*dds        \n"
      " sstart+2*ds (ss3),      ss3+dds,      ...,       ss3+(nds-1)*dds        \n"
      "    ...                                                                  \n"
      " sstart+(ns-1)*ds (ssn), ssn+dds,      ...,       ssn+(nds-1)*dds      \n"
      "       the output crossline locations (l):                               \n"
      " lstart      (ll1),      ll1+ddl,      ...,       ll1+(ndl-1)*ddl,       \n"
      " lstart+dl   (ll2),      ll2+ddl,      ...,       ll2+(ndl-1)*ddl        \n"
      " lstart+2*dl (ll3),      ll3+ddl,      ...,       ll3+(ndl-1)*ddl        \n"
      "    ...                                                                  \n"
      " lstart+(nl-1)*ds (lln), lln+ddl,      ...,       lln+(ndl-1)*ddl      \n"
      " 9. The primary memory requirement is                                    \n"
      "            ns*nds*nlcore*nofo*ntau*4 (byte) + velocity grid size (byte) \n"
      "            + ns*nds*nlcore*nfft*4 (byte) \n"
      "    nlcore is set to nl*ndl.                                             \n"
      " 10. when nlines is less than 10, a 2.5-D rho filter is applied to       \n"
      "     compensate for 2D data input                                        \n"
      "\n"

      "AUTHOR:          Zhiming Li,             3/93                    \n";

void     migs(int nl, int ns, int nofo, float *mig, int ntau,
              FILE * imgfp, int ktrace, int ntl, float dldt,
              float *trace, float *tras, float dtl, float tminl, int nt,
              float *ss, float *sl, float *gs, float *gl, int *imute,
              float *fovt2, float *s, float *l, float *tm, int *iofs,
              float apers, float aperl, float *fold, float *w1,
              float *work1, float *work2, float *wsave, float *tracef,
              int *ifcut, int *ltaper, float ksmax, float klmax,
              float f0, float df, int nf, float ftaper, int nsave,
              int nfft, float dtau, float tau0, float angmax,
              int *indxw, int nindxw, int incdxw, float *resamp,
              int ires, int ntres, float *ss2, float *scs, float *scg,
              int ncpu, float *vr2, float *vi2, float *vq4, float *tau);

void     fhigint_(float *f0, float *df, int *nf, float *ftaper,
                  float *dt, int *ifcut, int *ltaper, float *wsave,
                  int *nsave, int *nfft);

void     radix_(int *nt, int *nfft);

void     bisear_(int *n, int *nnew, float *x, float *xnew, int *indx);

void     f2p5n_(float *data, int *nt, int *i3d, int *m, int *l);

void     vr22vi2_(float *vr2, float *vi2, int *ntau);
void     vi22vq4_(float *vi2, float *vq4, int *ntau);
void     vt2s2_(float *vr2, float *vq4, float *s2, float *tau,

                int *ntau);

void     pstm3d_(float *trace, int *nt, float *t0, float *dt, float *sx,
                 float *sy, float *gx, float *gy, int *imute,
                 float *fovt2, float *mig, int *nxy, int *ntm, float *x,
                 float *y, float *tm, float *aperx, float *apery,
                 float *fold, float *w1, float *f0, float *df, int *nf,
                 float *ftaper, int *ifcut, int *ltaper, float *tracef,
                 float *wsave, int *nsave, float *work1, float *work2,
                 int *nfft, float *kxmax, float *kymax, float *dtau,
                 float *tau0, float *angmax, int *indxw, int *nindxw,
                 int *incdx, int *np, float *resamp, int *ires,
                 int *ntres, float *s2, float *scs, float *scg);


segytrace tra, tr;
segybhdr bh;
segychdr ch;


main(int argc, char **argv)
{
   int      nlines, ncdppl;
   int      nt, ntau, it, ix, iy;
   FILE    *infp, *outfp, *velfp, *imgfp, *hdrfp, *hffp, *fldfp;
   float    x1, x2, x3, y1, y2, y3, s1, s2, s3, l1, l2, l3, sm, lm;
   char    *backupi = "", *backupo = "";
   char    *diskimg = "DISKIMAGE", *diskhdr = "DISKHEADER";
   char    *velfile, *datain, *dataout;
   char    *diskfld = "DISKFOLD", *hisfile = "KTMIG.HISFILE";
   int      ibacki, ibacko, isave, ihis = 0, ipre = 0;
   float    ds, dl, tmp;
   int      mtrace, ntrpre = 0, ktrace, itrace, jtrace, ntrace, nofo;
   int      ltrace;
   int     *imutel, cdppre, cdpnow, m1000, ifact, it0;
   int      vread = 0, ivdisk = 1;
   float   *velos;
   float   *mig, *fold, *tm, *tmig, *tgain, *s, *l, *tugain, *trace,

         *tras;
   float   *fovt2;
   float   *ofo;
   int     *off, itgain, nl, ns, ntras, ntv, is, il, io, nsv, nlv;
   int      nsread, nlread;
   long long mlimit;
   int      tmp_mlimit;
   int      nlcore, nlpos, nlmig, *iofs, itv;
   float    tpow, sstart, lstart, t0, dtv, t0v, dsv, dlv, s0v, l0v;
   float    tpowaf;
   float    xs, ys, xg, yg, dt, dtau, tau0, v, *vel;
   float   *ss, *sl, *gs, *gl, ms, ml, ofomin, dofo, obtol;
   char    *buf;
   int      ilv, isv, itmp, itmp2, ntl, itr, ir0, i;
   float    tminl, dtl, dldt, apers, aperl;
   float   *w1;
   float    f0, df, ftaper;
   int      nf, nfft, *ifcut, *ltaper, nsave;
   int      mm, ll;
   float   *work1, *work2, *tracef, *wsave;
   float    ddl, dds, ksmax, klmax;
   int      ndl, nds, ntmp, tapecntl;
   float    rmsmar, ampmax, amprms;
   char    *badtracefile = "BAD_TRACE_FILE";
   FILE    *btfp;
   float    dd, mx, my, ofs;
   int      ids, idl;
   int      cdp1, cdp2, cdp3;
   int      cdppds, cdppdds, cdppdl, cdppddl, cdpnum;
   int      nss, nll, iix, iiy, is0, il0;
   float    sinmin, sinmax, linmin, linmax;
   float    angmax;
   char    *jpfile;
   FILE    *jpfp;
   int      ntrend, nsout, nlout, isend, ilend;
   float    send, lend, dsout, dlout, dend, *disend;
   float    ps, pl;
   int      iofoin = 0;
   float    ofol, ofor;
   int      one = 1;
   float   *resamp;
   int      ires, ntres;
   int     *indxw, nindxw, incdxw;
   int      ntrdsk, ndskwt;

   float   *vr2, *vi2, *vq4, *tau, *ss2, *scs, *scg;
   int      i3d, ncpu, traceout;
   char    *envs;
   int      on2trace, ntotal;


   ghed     gh;
   int      n1, n2, n3, n4, n5, dtype, ierr, orient, gtype;
   float    o1, o2, o3, o4, o5, d1, d2, d3, d4, d5, ocdp2, oline3,
         dcdp2, dline3;
   float    gmin = 0., gmax, scale;

   long long lwork, ltmp;
   long long lpos;

   float    tr1, tr2, tr3, ln1, ln2, ln3, fline, ftrace;
   String   tracekey = "tracl", linekey = "tracr", trktype, lnktype;
   Value    trkval, lnkval;
   int      indxtrk, indxlnk;
   int      intline, inttrace, ikey;
   int      strace, jthtrace;
   int      mintrace, maxtrace;
   int      minline, maxline;
   int      sline, jthline;

   int      ispre, isnow, ilpre, ilnow;

   /* get parameters */
   initargs(argc, argv);
   askdoc(1);

   /* open input and output data sets */
   if (!getparstring("datain", &datain)) {
      infp = stdin;
   } else {
      infp = efopen(datain, "r");
   }
   if (!getparstring("dataout", &dataout)) {
      outfp = stdout;
   } else {
      outfp = efopen(dataout, "w");
   }
   fseek2g(infp, 0, 1);
   fseek2g(outfp, 0, 1);

   /* required parameters */
   if (!getparfloat("x1", &x1))
      err("must specify x1");
   if (!getparfloat("y1", &y1))
      err("must specify y1");
   if (!getparfloat("s1", &s1))
      err("must specify s1");
   if (!getparfloat("l1", &l1))
      err("must specify l1");
   if (!getparint("cdp1", &cdp1))
      err("must specify cdp1");
   if (!getparfloat("x2", &x2))
      err("must specify x2");
   if (!getparfloat("y2", &y2))
      err("must specify y2");
   if (!getparfloat("s2", &s2))
      err("must specify s2");
   if (!getparfloat("l2", &l2))
      err("must specify l2");
   if (!getparint("cdp2", &cdp2))
      err("must specify cdp2");
   if (!getparfloat("x3", &x3))
      err("must specify x3");
   if (!getparfloat("y3", &y3))
      err("must specify y3");
   if (!getparfloat("s3", &s3))
      err("must specify s3");
   if (!getparfloat("l3", &l3))
      err("must specify l3");
   if (!getparint("cdp3", &cdp3))
      err("must specify cdp3");
   if (!getparfloat("dds", &dds))
      err("Must specify dds!\n");
   if (!getparfloat("ddl", &ddl))
      ddl = 0.;
   if (!getparstring("velfile", &velfile))
      err("must specify velfile");

   /* optional parameters */
   if (!getparfloat("sstart", &sstart))
      sstart = s1;
   if (!getparfloat("lstart", &lstart))
      lstart = l1;
   if (!getparfloat("ds", &ds))
      ds = dds;
   if (!getparfloat("dl", &dl))
      dl = ddl;
   if (!getparint("ndl", &ndl))
      ndl = 1;
   if (ndl < 1)
      ndl = 1;
   if (!getparint("nds", &nds))
      nds = 1;
   if (nds < 1)
      nds = 1;
   if (!getparfloat("angmax", &angmax))
      angmax = 60.;

   if (ds < dds)
      err(" dds larger than ds ");
   if (dl < ddl)
      err(" ddl larger than dl ");

   if (dds > 0.) {
      tmp = (s2 - s1) / dds + 1.5;
      ncdppl = tmp;
   } else {
      ncdppl = 1;
   }

   if (ddl != 0.) {
      tmp = (l3 - l1) / ddl + 1.5;
      nlines = tmp;
   } else {
      nlines = 1;
   }
   if (!getparint("nl", &nlread))
      nlread = nlines;
   if (!getparint("ns", &nsread))
      nsread = ncdppl;
   nl = nlread;
   ns = nsread;
   if (dl == 0. && nl > 1)
      err("Must specify nonzero dl for migration");
   if (ds == 0. && ns > 1)
      err("Must specify nonzeor ds for migration");

   if (dds == 0.) {
      dds = 1.0;
      ds = dds;
   }
   if (ddl == 0.) {
      ddl = 1.0;
      dl = ddl;
   }

   if (ncdppl > 1) {
      tmp = (cdp2 - cdp1) / (ncdppl - 1.) + .5;
      cdppdds = tmp;
      tmp = (cdp2 - cdp1) / (ncdppl - 1.) + .5;
      cdppds = tmp;
   } else {
      cdppdds = nlines;
      cdppds = nlines;
   }
   if (cdppdds < 1)
      cdppdds = 1;
   if (cdppds < 1)
      cdppds = 1;
   if (nlines > 1) {
      tmp = (cdp3 - cdp1) / (nlines - 1.) + .5;
      cdppddl = tmp;
      tmp = (cdp3 - cdp1) / (nlines - 1.) + .5;
      cdppdl = tmp;
   } else {
      cdppddl = ncdppl;
      cdppdl = ncdppl;
   }
   if (cdppddl < 1)
      cdppddl = 1;
   if (cdppdl < 1)
      cdppdl = 1;

   if (cdppdds != 1 && cdppddl != 1) {
      warn("cdppdds=%d cdppddl=%d ", cdppdds, cdppddl);
      err("check cdp1, cdp2, cdp3 ");
   }
   if (cdppdds == 1) {
      cdpnum = 0;
   } else {
      cdpnum = 1;
   }


   tmp = ds / dds + 0.5;
   nss = tmp;
   tmp = dl / ddl + 0.5;
   nll = tmp;
   tmp = (sstart - s1) / dds + 0.5;
   is0 = tmp;
   tmp = (lstart - l1) / ddl + 0.5;
   il0 = tmp;


   if (!getparint("mtrace", &mtrace))
      mtrace = 1000;
   if (!getparint("ntras", &ntras))
      ntras = 2000000000;