**Author:**- Mostafa Naghizadeh

1. Minimum Weighted Norm Interpolation (MWNI).

2. Multi-Step Auto-Regressive (MSAR) reconstruction.

3. The F-K domain interpolation (Gulunay)

This is the first version of the library and may not be as robust or well documented as it should be. Please keep track of bugs or missing/confusing instructions and report them to mostafan@ualberta.ca

2. Liu, B. and M. D. Sacchi (2004). Minimum weighted norm interpolation of seismic records. Geophysics 69 (6), 1560-1568.

3. Gulunay, N. (2003). Seismic trace interpolation in the fourier transform domain. Geophysics 68 (1), 355-369.

2. Unzip the package to your home directory.

>gunzip seisrec.tar.gz >tar -xf seisrec.tar

3. You need to add a line to the file .cshrc inside your home directory. Do the following steps:

> emacs .cshrc

setenv SEISREC $HOME/seisrec

> source .cshrc

4. Enter seisrec directory

> cd seisrec

5. Now you can configure the package using

./seisrec> configure -s

6. In order to compile and make executable file for examples, just type make in the command line:

./seisrec> make

1 #include <stdio.h> 2 #include "msar.h" 3 int main(void) { 4 int NDIM,ia,ind,iter_cg,iter_bl; unsigned long *nx,*nhan,nh,js,*npf; float fmin,fmax,*bmin,*bmax,*rzp,fmaxp; msar_plan ms_plan; mwni_plan mw_plan; FILE *fid; 5 NDIM=2; // Dimension of data nx=lvector(1,NDIM); // Vector for number of samples in each direction of filter rzp=vector(1,NDIM); bmin=vector(1,NDIM-1); bmax=vector(1,NDIM-1); nhan=lvector(1,NDIM-1); npf=lvector(1,NDIM-1); 6 // Number of whole and avialable traces and time smaples fid=fopen("data_count.dat","r"); fscanf(fid,"%lu",&nx[1]); fscanf(fid,"%lu",&nx[2]); fscanf(fid,"%lu",&nh); fclose(fid); 7 // Hanning dimension nhan[1]=3; 8 // zero-padding ratio rzp[1]=0; rzp[2]=1.5; 9 // band-limiting values fmin=0.035; fmax=0.07; bmin[1]=.45; bmax[1]=.55; 10 // type of MWNI ind=1; 11 // Number of iterations iter_cg=10; iter_bl=3; 12 // msar parameters js=8; fmaxp=.5; npf[1]=5; 13 // Initiating MSAR plan init_msar_plan(&mw_plan,&ms_plan, NDIM , nx, nh, rzp, fmin, fmax, bmin, bmax, ind, iter_cg, iter_bl, nhan, npf, fmaxp, js); 14 // Giving values of avialable traces fid=fopen("avail_data.dat","r"); for(ia=1;ia<=mw_plan.nh*mw_plan.n[1];ia++) { fscanf(fid,"%f",&mw_plan.data[ia]); } fclose(fid); 15 // Giving the location of avialable tarces fid=fopen("offsets.dat","r"); for(ia=1;ia<=mw_plan.nh;ia++) { fscanf(fid,"%lu",&mw_plan.h[ia]); } fclose(fid); 16 // Performing MSAR reconstruction perform_msar(&mw_plan, &ms_plan); 17 // Saving values of msar reconstructed traces fid=fopen("msar_rec_data.dat","w"); for(ia=1;ia<=mw_plan.nr*mw_plan.n[1];ia++) { fprintf(fid,"%8.5f \n",ms_plan.rdata[ia]); } fclose(fid); 18 // Saving values of mwni reconstructed traces fid=fopen("mwni_rec_data.dat","w"); for(ia=1;ia<=mw_plan.nr*mw_plan.n[1];ia++) { fprintf(fid,"%8.5f \n",mw_plan.datar[ia]); } fclose(fid); 19 // Saving Number of PFs extracted for each frequency fid=fopen("pfc.dat","w"); for(ia=1;ia<=ms_plan.ntz1/2+1;ia++) { fprintf(fid,"%lu \n",ms_plan.pfc[ia]); } fclose(fid); 20 // Finalizing MWNI paln finalize_msar_plan(&mw_plan,&ms_plan); 21 return 0; }

1. This line includes the C library stdio.h, which is required for standard input and output management of the data.

2. Here we include the library msar.h in our example in order to apply the MSAR reconstruction. It should be noticed that the library mwni.h is included inside msar.h and as a result it would be included in our example as well.

3. This is a standard C command to start a main C code.

4. These lines define the parameters which are going to be used inside this example. In C programming all the parameters must be defined at the top of the program before appearing inside the code. Besides the ordinary parameter formats of C such as int, float, unsigned long and FILE, we also define two other formats, msar_plan and mwni_plan, which will contain all the parameters needed for our job at the later stages.

5. In this line we define the dimension of the data. For this example we have a 2D seismic section as input so we have NDIM=2. In the following lines we allocate memory for the arrays we need to pass them as input to the algorithm. Further details of each of these parameters are given at below lines.

6. In this line, the number of samples in each direction (nx[1] (time) and nx[2](space)) and also the number of available traces (nh) are given to the code. These values are already saved in a file called “data_count.dat”.

7. This line defines the length of the Hanning window used for smoothing inside the MWNI. The reconstruction algorithm for 2D data (time-space) is applied for the elements of each single frequency. Since, in this example, we have only one space direction, the frequency slices are 1D. This means the Hanning window has to be 1D as well. For 3D data (time-space-space) the Hanning window would be a 2D window.

8. In order to have an artificial-free reconstruction in the boundaries of the data, it is important to zero-pad the sides of the data. This is because an abrupt cutting of events at the boundaries of the data causes severe artifacts in the Fourier domain. In this line user determines the ratio of zero-padding for each axis of data. Zeros are padded in both end of each axis with the indicated ratios. By setting rzp[1]=0 we do not zero-pad the time direction for this example. However for the spatial direction rzp[2]=1.5 asks to padd each end of the spatial direction with 1.5 size of the number of spatial sample (tarces). For example, if nx[2] is equal to 32 then the final size of spatial samples after zero-padding would be equal to 48+32+48=128 for this specific example. When the reconstruction is done, these zeros would be cut out from boundaries (Even though there are not zeros any more). This action would prevent any kind of artificial events appearance at the sides of the data.

9. This line determines the size of the low frequency band to be reconstructed by MWNI. We represent the frequency axis with the values in the range [0,.5]. Here we chose fmin (minmum frequency to be reconstructed using MWNI) equal to 0.035 and fmax (maximum frequency to be reconstructed using MWNI) equal to 0.07. It’s the user's duty to check the data and pick the proper values for fmin and fmax. The next step is to determine how much band-limited we want each frequency slice to be. This is done by the parameters bmin[1] and bmax[1]. Their value should be at the range (0,1). Zero means complete band-limiting (forcing all the wave-numbers to be zero) and one means no band-limiting (none of the wave-numbers are forced to be zero and all of the determined by the data). While we just define the band-limiting parameters for fmin and fmax, the band-limiting value for the frequencies in between is determined by a linear interpolation scheme. Notice that both bmin and bmax have only one element for 2D data. For the case of 3D data these parameters would have two elements each of them representing each spatial direction (for more see the 3D examples provided in example folder) and so for. Interested user are urged to try value higher than 1 for bmax and see how it effects the band-limiting for the frequencies in between.

10. In order to examine the different aspects of Fourier reconstruction, one can choose various type of applying MWNI by setting the parameter ind to 0, 1 and 2. Zero means Band-Limited Fourire Reconstruction (BLFR) without reweightening, one means MWNI (Band-limiting and weigthening) and two means Sparse Reconstruction (No band-limiting).

11. The optimization of MWNI method is done using Conjugate Gradient (CG) algorithm. The number of iterations for CG algorithm is determined by iter_cg which is set to 10 for this particular example. The CG routine has to be repeated inside an iterative routine to gain new weights for the data reconstruction purposes. The number of iterations for updating the weight function is determined by iter_bl which is fixed to be 3 for this example.

12. This line shows the parameters needed for prediction filetr estimation and MSAR method. The Maximum number of jumping steps for the MSAR algorithm is determined by js. While it is defined by user, but still it depends on the number of spatial sample and also the number of elements of prediction filter. In this case if js was bigger then the one is driven by data, it would be replaced by the new calculated js inside the code. On can also determine the maximum frequency (s)he wants to do the MSAR reconstruction by setting the value of fmaxp. Again this can not be bigger than 0.5. Also the length of prediction filters is determined by npf[1]=5 in this example. Notice that for 3D data we needed to define the length of Pf by setting another direction as npf[2].

13. In this line the above tuned parameters are fed to the subroutine inti_msar_plan. This subroutine will allocate memory and initial values for all the vectors and parameters needed for the MSAR reconstruction. For the MSAR plan we initiate one msar_plan (ms_plan) and one mwni_plan (mw_plan). These two plans will save all the parameters during application of the method and there is no need to define new parameter for each subroutine inside the code. This means by just passing the plans inside each subroutine all required information will be delivered to that specific subroutine.

14. In this line the values of available traces are given to algorithm. The available traces are already saved in the file "avail_data.dat" in a lexicographic order.

15. One also needs to give the locations of the available traces in the spatial direction. The locations of traces are also is saved a file (offsets.dat) in a lexicographic order (for more than one spatial dimension).

16. This line is the main part of applying the MSAR reconstruction. The subroutine called "perform_msar" will take the values of parameters initiated inside mw_paln and ms_plan and perform the MSAR reconstruction. The desired results are saved in the proper vectors located inside these plans.

17. In this line gives the result of MSAR reconstruction are saved in the file called "msar_rec_data.dat".

18. This line saves the result of MWNI part of data in the file "mwni_rec_data.dat".

19. In this line the number of extracted prediction filters for each frequency is saved in the file "pfc.dat".

20. Like all computationally intensive C programs, one should free all the vectors which have been allocated. In this line all the allocated memories for vectors inside the MSAR and MWNI plans are freed.

21. Who cares what this line does?

After running each example you can also make postscript plot and x plot presentation of the results by running (if that specific exampl has a file called PSplot.sh in its folder and you have Seismic Unix (SU) installed in your PC):

seisrec/your_example> PSplot.sh

In order to clean the results and end up with the source files that were in the folder before you run the example do the followings in your command line of x window:

seisrec/your_example> make clean seisrec/your_example> make