3D Land Seismic Survey

Explorasi seismic 3D merupakan teknologi pencitraan (imaging) bawah permukaan secara tiga dimensi.  Berbeda dengan seismic 2D yang mencitrakan point tertentu atau ‘titik’ maka seismic 3D adalah teknologi untuk mencitrakan ‘bidang’. Seismic 3D memiliki kelebihan untuk meng-eliminasi mis-tie dalam migrasi reflector miring, meningkatkan resolusi horizontal, dan memberikan citra yang lebih detail.

Berikut adalah terminologi yang sering digunakan dalam Explorasi Seismic 3D:

• Inline: garis-garis semu yang parallel dengan bentangan receiver.
• Crossline: garis semu yang tegak-lurus dengan Inline.
• CMP bin: kotak semu  di bawah permukaan dengan ukuran ½RI*½SI dimana RI adalah Interval receiver dan SI interval Source. CMP bin mengandung semua trace yang dimiliki oleh CMP yang sama.
• Patch: area dari reveiver yang merekam source yang sama.
• Swath: area dimana receiver merakam sumber-sumber tanpa adanya perpindahan crossline (crossline roll over).
• Salvo: sejumlah sumber tembakan yang direkam oleh patch yang sama.
• Fold:  banyaknya mid-point yang di-stack dalam CMP bin yang sama. Besaran Fold berbeda dari bin ke bin sejalan dengan perubahan offset dan azimuth serta berubah terhadap kedalaman sejalan dengan bertambahnya offset.  Fold=NS*NR*b², dimana NS dan NR  adalah banyaknya Source  dan Receiver dalam wilayah  tertentu dan b merupakan dimensi bin. Contoh jika per kilometer persegi terdapat 80 source dan 600 receiver dan dimensi bin 25m  maka Fold=80*600*25*25 m²/km²=30.
• Crossline Fold: setengah dari jumlah inline dalam satu patch. Jika dalam satu patch terdapat 8 inline maka Crossline Fold=8/2=4.
• Inline Fold: Fold/Crossline Fold. Untuk contoh kita 30/4=7.5. Dengan demikian Fold=Crossline Fold*Inline Fold=7.5*4=30.

Berikut adalah contoh untuk mendesign sebuah survey land 3D dengan kedalaman target=3000m, bin=25m dan Fold=30  dengan sistem split-spread (sumber di tengah). Dengan Interval lintasan receiver 400m:

• Receiver Interval dapat ditentukan dengan 2xbin=2x25=50m.
• Offset Maximum: katakanlah 90% dari kedalaman target, 3000mx90%=2700m.
• Jumlah masing-masing receiver pada setiap sisi split spread: 2700/50(receiver interval)=54 receiver.
• Total perekam setiap line setiap shot=2*54=108.
• Jumlah receiver yang harus diaktifkan jika hanya tersedia 900 receiver, 108* 8=864 receiver (untuk 1 patch). Maka kita dapat memiliki 8 lintasan receiver.
• Shot interval biasanya 2*bin=2*25=50m.
• Crossline fold=8(banyaknya line per patch)/2=4
• Inline Fold=30/4=7.5
• Shot line Interval (SI) dapat  ditentukan dengan NI=(Total perekam per line/2)*Receiver interval/SI.   7.5 =(108/2)*(50/SI). Jadi SI=360m.

Terdapat beberapa teknik shooting seismic 3D, diantaranya adalah Metoda Swath Shooting:

1. Lintasan-lintasan receiver dibentangkan secara parallel.
2. Sumber-sumber ledakan dipasang secara tegak lurus dengan lintasan receiver.
3. Sumber pertama diledakkan lalu dilakukan perekaman.
4. Sumber kedua-ketiga dst sampai ke-terakhir (dalam satu patch) diledakkan dengan perekaman dilakukan untuk masing-masing ledakan.
5. Serangkaian ledakan diatas disebut dengan Salvo-1.
6. Pindah ke source line berikutnya, lakukan hal yang sama sehingga diperoleh salvo-2, dst.
7. Beberapa salvo dilakukan sampai akhirnya sampai di ujung lintasan receiver sehingga diperoleh satu swath.
8. Roll-over sebesar setengah patch kearah crossline untuk memperoleh swath 2, dst sampai seluruh areal 3D.

Berikut Animasinya (klik untuk memperbesar) atau anda dapat mendownload versi resolusi tinggi di sini.

Inversi Matrix Singular-SVD (Part-2)

Pada bagian ini saya telah menjelaskan teori SVD serta aplikasinya dengan Matlab serta SVDLIBC. Pada bagian ini saya akan menjelaskan penggunaan redsvd serta menghitung matrix inverse dengan input hasil SVD dengan menggunakan kode C.

Saya telah berhasil men-test program redsvd dan kode C berikut pada PC untuk matrix berukuran 5000x10000, sementara out of memory jika saya  menggunakan Matlab.

Untuk dapat menggunakan redsvd anda harus memiliki platform Ubuntu 12.04.1 LTS. Setelah Ubuntu 12.04.1 LTS diinstall, anda perlu meng-install eigen3 melalui synaptic.

1. Install Synaptic
    sudo apt-get update
    sudo apt-get upgrade
    sudo apt-get install synaptic

2. Install eigen3
    Dash Home > search Synaptic  > Synaptic Package Manager
    Select libeigen3-dev dan libeigen3-doc lalu Apply

3. Install g++
    sudo apt-get update
    sudo apt-get install g++

4. Download redsvd di sini
    Setelah didownload (pada kasus ini saya mendownload redsvd-0.2.0.tar.bz2)
    tar xvjf redsvd-0.2.0.tar.bz2
    cd redsvd-0.2.0
    ./waf configure
    ./waf
    sudo ./waf install

5. Hitung SVD dengen redsvd
    redsvd -i file1 -o file1 -r 10 -f dense

-i file1 artinya memanggil file1 sebagai input, -o file1 mengoutputkan file1.U, file-1.S, file1.V  (V belum di-transpose) -r 10 artinya rank-10 atau mengoutputkan 10 non zero eigen terbesar -f dense artinya matrix input dense (bukan sparse). Dari pengalaman saya, nilai rank adalah nilai terbesar dari jumlah kolom atau baris dari matrix. Contoh matrix saya berukuran (5x4) saya pilih rank 5, matrix saya (3x7) saya pilih rank 7. Sebagai antisipasi, memilih rank yang besar tidak masalah sehingga tidak ada nilai eigen yang tidak di-outputkan.

Pada contoh ini saya memiliki file1 (matrix ini singular):
    1.0     0.0     0.0   
   -2.0     0.0     0.0
    4.0     6.0     1.0

6.  Untuk memperoleh inverse dari matrix file1, copy-lah kode dibawah ini lalu beri nama inverse.c

7. Jalankan dengan perintah gcc -o inverse inverse.c

8. ./inverse

9. Perhatikan untuk mengubah nama file-file input, jumlah baris dan kolom dari input dan nama output.

10. Berikut adalah inverse dari matrix file1
~/Documents$ cat output.txt
     0.200000   -0.400000   -0.000000
    -0.129730    0.259460    0.162162
    -0.021622    0.043243    0.027027


//This program calculates inverse of a matrix
//The input is U,S,V that are decomposed via SVD Decomposition
//Ainv=V*So*U' in which So=1/S
//by Agus Abdullah,PhD 2013.
#include<stdio.h>
#include <stdlib.h>

int main(void)
{
float **matrixU,*matrixS, **matrixV, **matrixI, **matrixD;
float **matrixSo,**matrixVSo, **matrixUT, **matrixVSoUT;
int row,col,i,j,nav1, nav2;
double eps;
//The matrix that you want to be inverted has row=??? and col=???
row=3;
col=3;
eps=0.0001;

//Read U
    matrixU=(float **)malloc(row*sizeof(float));
    for(i=0;i<row;i++)
    matrixU[i]=(float *)malloc(row*sizeof(float));

    FILE *inputU=fopen("file1.U","r");
    for(i=0;i<row;i++)
    {
        for(j=0;j<row;j++)
        {
    fscanf(inputU,"%f",&matrixU[i][j]);
        }
    }
    fclose(inputU);

printf("Step1 done... \n");

//Read V
    matrixV=(float **)malloc(col*sizeof(float));
    for(i=0;i<col;i++)
    matrixV[i]=(float *)malloc(row*sizeof(float));
    FILE *inputV=fopen("file1.V","r");
    for(i=0;i<col;i++)
    {
        for(j=0;j<row;j++)
        {
          fscanf(inputV,"%f",&matrixV[i][j]);
        }
    }
    fclose(inputV);
printf("Step2 done... \n");
//Read S
    matrixS=(float *)malloc(row*sizeof(float));
    FILE *inputS=fopen("file1.S","r");

    for(i=0;i<row;i++)
    {
    fscanf(inputS,"%f",&matrixS[i]);
    }
     fclose(inputS);

//Matrix Identity
    matrixI=(float **)malloc(row*sizeof(float));
    for(i=0;i<row;i++)
    matrixI[i]=(float *)malloc(row*sizeof(float));

    for(i=0;i<row;i++)
    {
        for(j=0;j<row;j++)
        {
          if(i==j)
        {           
        matrixI[i][j]=1.0;
         }
    else   
        {
        matrixI[i][j]=0.0;
        }          
         }
    }

printf("Step3 done... \n");
//Matrix S Diagonalization
    matrixD=(float **)malloc(row*sizeof(float));
    for(i=0;i<row;i++)
    matrixD[i]=(float *)malloc(row*sizeof(float));
    for(i=0;i<row;i++)
    {
        for(j=0;j<row;j++)
        {
    matrixD[i][j]=matrixS[i]*matrixI[i][j];
        }
    }

printf("Step4 done... \n");
//Matrix 1/S
    matrixSo=(float **)malloc(row*sizeof(float));
    for(i=0;i<row;i++)
    matrixSo[i]=(float *)malloc(row*sizeof(float));
    for(i=0;i<row;i++)
    {

        for(j=0;j<row;j++)
        if (matrixD[i][j] > eps)
        {
    matrixSo[i][j]=1/matrixD[i][j];
        }
        else
    matrixSo[i][j]=0.0;
    }

printf("Step5 done... \n");
//V*So --->V from redsvd is not transposed
    matrixVSo=(float **)malloc(col*sizeof(float));
    for(i=0;i<col;i++)
    matrixVSo[i]=(float *)malloc(row*sizeof(float));
    for(i=0;i<col;i++)
    {

        for(j=0;j<row;j++)
        {   matrixVSo[i][j]=0;
            for(nav1=0;nav1<row;nav1++)
            matrixVSo[i][j]+=matrixV[i][nav1]*matrixSo[nav1][j];
        }
    }

printf("Step6 done... \n");
//Compute Transpose of U
 matrixUT = malloc(row * sizeof(float));
 for (i = 0; i < row; i++)
 matrixUT[i] = malloc(row * sizeof(float));

   for (i = 0; i < row; i++)
   {
        for (j = 0; j < row; j++)
        matrixUT[j][i] = matrixU[i][j];
   }

printf("Step7 done... \n");
//V*So*U'
    matrixVSoUT=(float **)malloc(col*sizeof(float));
    for(i=0;i<col;i++)
    matrixVSoUT[i]=(float *)malloc(row*sizeof(float));
    for(i=0;i<col;i++)
    {
        for(j=0;j<row;j++)
        {   matrixVSoUT[i][j]=0;
            for(nav2=0;nav2<row;nav2++)
            matrixVSoUT[i][j]+=matrixVSo[i][nav2]*matrixUT[nav2][j];
        }
    }

printf("Step8 done... \n");

   FILE *f;
   f = fopen("output.txt","w");
   for(i=0; i<col; i++)
    {
        for(j=0; j<row; j++)
        {
        fprintf(f,"\t%.4f",matrixVSoUT[i][j]);
        }
        fprintf(f,"\n");
       }
   fclose(f);

printf("completed! \n");
//Free Memory
    for(i=0;i<row;i++)
    free(matrixI[i]);
    free(matrixI);
    for(i=0;i<row;i++)
    free(matrixU[i]);
    free(matrixU);
    for(i=0;i<col;i++)
    free(matrixV[i]);
    free(matrixV);
    for(i=0;i<row;i++)
    free(matrixD[i]);
    free(matrixD);
    free(matrixS);
    for(i=0;i<row;i++)
    free(matrixSo[i]);
    free(matrixSo);
    for(i=0;i<col;i++)
    free(matrixVSo[i]);
    free(matrixVSo);
    for (i=0; i< row; i++)
    free(matrixUT[i]);
    free(matrixUT);
    for(i=0;i<col;i++)
    free(matrixVSoUT[i]);
    free(matrixVSoUT);
}

C Language for Geoscientists (Part-7)

Kode C berikut adalah cara untuk mengkonversi format matrix. Malloc digunakan sehingga kode ini mampu menghandle data yang sangat besar.

Berikut adalah Kode format yang digunakan dalam C:
Code     Format
----    ------
%c     character
%d     signed integers
%i     signed integers
%e     scientific notation, with a lowercase "e"
%E     scientific notation, with a uppercase "E"
%f     floating point
%lf    double
%g     use %e or %f, whichever is shorter
%G     use %E or %f, whichever is shorter
%o     octal
%s     a string of characters
%u     unsigned integer
%x     unsigned hexadecimal, with lowercase letters
%X     unsigned hexadecimal, with uppercase letters
%p     a pointer
%%     a '%' sign


Untuk contoh ini, saya memiliki filein dengan format scientific:
1.7159866e-001  5.0353191e-002  8.8042207e-002 
5.8189807e-001  2.3275270e-001  6.8729395e-001

Dan Mengoutputkannya dengan floating point:
0.171599 0.050353 0.088042
0.581898 0.232753 0.687294

#include <stdio.h>
#include <stdlib.h>
float **matrixA;
int rowA,columnA,row,column;

void readA()
{
FILE *inputA=fopen("filein","r");
 for(row=0;row<rowA;row++)
    {
        for(column=0;column<columnA;column++)
        {
    fscanf(inputA,"%e",&matrixA[row][column]);
        }
    }
  fclose(inputA);
}

//Main
int main()
{
    //Define row and col
    rowA=2;columnA=3;
    matrixA=(float **)malloc(rowA*sizeof(float));
    for(row=0;row<rowA;row++)
    matrixA[row]=(float *)malloc(columnA*sizeof(float));
    readA();
   
   int i,j;
    for(i=0;i<rowA;++i)
        {
        for(j=0;j<columnA;++j)
        {
            printf("%f ",matrixA[i][j]);
         }
        printf("\n");
        }

    for(row=0;row<rowA;row++)
        free(matrixA[row]);
    free(matrixA);

}