Basic Laboratory Exercises

for:
Reflection Seismology using
Seismic Unix (CSM), Linux
and Perl
 (GEOL 4068)
Fall 2013
LSU
Baton Rouge
Juan M. Lorenzo
Contents
Acknowledgements 4

Introduction 5

Why do we need to study linux? 5

Why do we need OpenSource software? 5
Where do I get ssh? 5
Are you planning on doing any programming from home? 5
Where do I get Xming? 5
What do I do if I have a Macintosh computer? 5
How to run Xming: 6
Why do we need to know sh or Perl? 6
Linux 6
History of Linux 6
Q. What is a shell? 7
Q. What are the different shells?7
Q. Which one should I use? 7
Directory Structure of the Linux operating system 7
Seismic Unix 12

Introduction12
Perl (Hoffman, 2001) Practical Extraction and Report Language 21

Introduction to Perl 21

Q. Why use Perl? 21

When not to use Perl? 21
Can I use Perl to make simple, visually interactive programs? 21
Learning Perl on your own 21
Basic components of Perl 22

Input and Output 22

Documentation in Perl 24
Data Types 25
 Useful Operations 27
Incorporating SeismicUnix programs into Perl 30
MATLAB 32

Create a matrix of numbers 32

Sin function 32
Exercise 1 Simple 2D plotting 32
Exercise 2: Traveltime Equations33
Exercise 3: An ideal seismic wave signaturethe spike - 33
Exercise 4 Constant Phase and Linear Phase 34
Exercise 2- Matlab 34
% hyperbola 34
Matlab code for exercise 4 A study of the effects of constant phase and
linear phase on a seismic wavelet 35
CONSTANT PHASE SHIFTS 35
Matlab code for exercise 4 A study of the effects of constant phase and
linear phase on a seismic wavelet 36
Creating a shell script to log in automatically 36
Acknowledgements
These notes borrow from the Colorado School of Mines (Stockwell) for
S*nix, Universities of Indiana and Buffalo for linux and the University of
Florida for Matlab. Many students have also contributed to these notes:
Class of 2008: Erin Walden, Kody Kramer, Erin Elliott, Andrew Harrison,
Andrew Sampson, Ana Felix, JohnDAquin, Russell Crouch, Michael
Massengale, David Smolkin
Introduction
Why do we need to study linux?

Creative professional geophysicist and academics are able to explore new

ideas without constraints of black-box software.

Why do we need OpenSource software?

Scientifically, open source products can be verified independently by

anyone. Reproducibility is a core tenet of the scientific method. OpenSource
software replicates a scientific procedure.

Where do I get ssh?

Link to ssh: http://web.wm.edu/it/?&id=2948&svr=www

Are you planning on doing any programming from home?

Open SSH. Create a profile named lgc10.

Now go to File> Profile>Edit Profile. Edit the lgc10 profile.

On the Connections tab: Hostname field: lgc10.geol.lsu.edu) Username:

the user name given you in class by Dr. Lorenzo. Your password is of the
form XXXXXXX, where XX is a number given you by Dr. Lorenzo. The
password is case-sensitive. Save changes to your profile.

You can now connect to the lgc10 server using SSH.

Where do I get Xming?

Xming is the leading, free X Windows Server for Microsoft Windows.

For notes link to Xming: http://www.straightrunning.com/XmingNotes/

For download of X fonts, use Google, e.g:

http://sourceforge.net/projects/xming

For download of Xmin server, http://download.cnet.com/Xming/3000-

2094_4-10549058.html
What do I do if I have a Macintosh computer?

For Mac users connecting to lgc10@geol.lsu.edu, all you need to do is

going to the "Launch pad" on your desktop, open Utilities folder and double
click "Terminal". A terminal window will open. Type: ssh -X
yourname@lgc10.geol.lsu.edu. Then you are ready to go. (Courtesy of Chang
liu, fall 2013)

1) Download Cyberduck from internet (for Mac).

2) Open Cyberduck, click open connection.

3) Choose SFTP (Ssh File Transfer Protocol) option in the drop down
window on the very top.

4) Type "lgc10.geol.lsu.edu" in the "Server"

5) Use the login name and passwork Dr. Lorenzo gave you to login.

6) After putting all the information, click Connect. And you will see all the
files

7) Simply, Drag the file(s) into the folder you want it to be.

How to run Xming:

Making sure that youre still connected in SSH, run XLaunch to configure
Xming to connect to lgc10. Choose one window, then make sure that Start
no client is checked. Click Next>Next>Finish. Logout of
SSH(File>Disconnect) and then reconnect by selecting the lgc10 profile.

If you are having problems connecting, open the lgc10 profile in SSH and
go to Edit>Settings. Under the Tunneling option on the tree, make sure that
the Tunnel X11 Settings option is checked. Make sure to save your profile.

You now know you correctly edited the .login file if it reads DISPLAY:
undefined variable. If you get something with error in it, check to make
sure the setenv line is commented out.

Why do we need to know sh or Perl?

Shells are the basic sets of instructions for handling the operating system
and perl is a mature, widespread computer language ideal for file
manipulation. Perl can serve as a simple glue to make diverse pieces of
software talk to each other.

Linux
 The single-most advantage of linux is that the code is freely available so
many people around the world participate in its improvement continuously.
I first view Linux as a communal, philanthropic exercise which takes
advantage of the cooperative nature of our species. Linux is also a
collection of instructions in software that allow you to use the hardware in
your computer.

If well thought out, visually identifiable commands are friendlier if but

slower to use, (although especially tedious to write and computationally less
efficient). As part of linux there is a point and click WYSIWYG (What-you-
see-is-what-you-get)/GUI(Graphical-user-interface) to drive the same
instructions, visually.

History of Linux

Click here for a more comprehensive history of the subject by Ragib

Hasan at UIUC.

Linux was developed (for free) by Linus Torvald possibly inspired by at

least the GNU project (GNUs not Unix) , a software movement to provide
free and quality software

LINKS to sites that have important shell instructions:

Important Instructions in sh

Following: Linux Shells (Albing et al., 2007)

Q. What is a shell?

A shell is a convenient collection of command-line-instructions (actual

programs), written in a low-level language, such as C, which allow the user
to interact with files and the hardware and files. Shells have been around
since the start of the unix-type operating systems and have the advantage
that they interchangeable among different linux operating systems.
Although the instructions may have to be recompiled for each machine the
syntax remains constant and once learnt will last a career.

Example, ls.

ls stands for: list the contents of this directory

Q. What are the different shells?

sh: the original bourne-shell

csh: theC-shell
 The csh improves upon the sh because it introduces convenient
programming tools inherited from C

ksh: the k-shell

The commercial nature of this shell limited the growth of its popularity
from the start.

bash: the bourne-again-shell

The bash shell is ubiquitous among any linux-type operating systems you
might encounter. The bash shell inherits the advantages and experiences of
all prior shells.

Q. Which one should I use?

For this class the default is: csh

Directory Structure of the Linux operating system

In any operating system, linux programs and user directories are stored
in predictable locations. Exercise
 Logging in to your account

Type your login id, followed by your password

Running a Remote Session on "lgc10" and forwarding it to your local

machine CONTENTS

ssh X yourname@lgc10.geol.lsu.edu This command will redirect images

to the machine you are sitting at.

answer "yes" to the question involving "authenticity". You should only see
this question the first time you log on from each machine.

You should see a "prompt" such as

Lgc10:/home/yourname %

ls -l (see what's in your directory)

Changing System permissions and the stability of Linux

Every file and directory in linux has assigned codes which dictate the
degree of authority by each user of the computer to alter each file. There
are four types of user status on linux. First is the overall supreme
administrator known as root and who can do anything to any file on the
system. Next comes the specific original owner/user of each file. All users
can belong to one or several named groups of users. Finally anyone who is
not specified as belonging to your group or is not the supreme administrator
is considered belong to the outside world, or all other users. Within each
of the status levels: owner, group, world, binary codes or their letter
equivalents may be set to indicate whether a file may be only browsed
(read), modified (write), and/or executed as a program (executable).
Note that it is the files themselves that carry this important information with
them. The file permissions are consulted first to determine whether an
individual user has authority to manipulate the file in any way.

The purpose of this complex permission scheme is to provide an infinite

variety of protection schemes for the file systems but yet maintain an
unsinkable file system. In theory, and for much of practice, an individual
user will not be able to shut down the system; they will only be able to do
damage to themselves and not the files or others.

System permissions belonging to a file or directory can only be changed

by those users for whom files have had the proper permissions already
assigned. Initially it is root that sets all the first set of permissions for files
and directories when a user is given a space to work on the system. From
the first logon, the user has control of their assigned set of files and
directories.

If you want a file containing Perl code to become executable in the

system the creator of the file is required to change the appropriate
permission setting for that file. Following are the equivalent numeric codes
for the different types of permissions:

Read only - 4 Write only 2 Execute only - 1

Read and write 6 Write and execute 3 Read, write and execute
7 (add all three numbers together)

For example:

% ls l

My_perl_file r _ _ r _ _ r _ _

There are three spaces to explain the type access by user:

(read access), group (nothing/0) and everyone-else (dash/0),

respectively. The next three spaces show the same for the group to which
the user belongs and the final three for all other users.

In order to change permissions to allow the file to run as a program

enter the following:

chmod u+x

which only adds (+) the setting that gives only the owner (u)
executing privileges

Or, equivalently

chmod 600

In the numeric form the last two zeros mean that group and others
priviliges are nill. As you can see the numeric form can alter permissions for
all the three types of linux users at once.
 Here is a summary list of options used for setting file permissions and
understanding file types on the linux system

Abbreviation Stands Abbrev Stands for

of user status for iation of
file
permissio
n

u user r read

g group w write

o others x execute

a all

+ add

- remove

d director
y

l link

Examples:

Letters Numeri
symbols cal symbols

chmod chmod
u+rwx 700

chmod
u+rwx
chmod
chmod g+rw 761

chmod o+x

Q. Can I do any damage to another persons files?

 Yes, if the files belong to you. You can tell if you own the files by reading
the second column from the ls l instruction, which has the general form

drwxr-xr-x number of links your login name your group name

filesize(bytes) date etc.

Copy a file to your home directory

%cp /home/refseis13/pl/

use Control D or TAB to complete your file name

Copying files across the web securely using sftp

From your local machine type

%sftp loginname@remotemachinename

Once you are connected to the remote machine, the following basic
instructions will get you going:

help

get download a file over to the directory on the local machine

put upload a file to the remote machine

ls list CONTENTS of the remote machine

lls list directory CONTENTS of the local machine

pwd working directory name of the remote machine

lpwd working directory name of the local machine

(you can also type help once you are inside the remote machine)

Deleting files

%rm filename

Finding files

% locate filename

Renaming files

%mv filename
A free linux shell scripting tutorial:

http://www.freeos.com/guides/lsst/

Example 1 Example 2

The text The number 1 is assigned to

hello is assigned the variable named value. The
to the variable value of the variable is expressed
named output as $value
The value of the
variable is $1 is assinged value 2 from
expressed as the command line (outside the
$output The shell script). This number is the
variable name first value on the command line
can be any word. after the prog name

Arithmetic calculations are

carried out by a shell program
called expr.

%prog_name %prog_name 2

#! /bin/sh #!/bin/sh
echo "Enter the two numbers to
output=hello be added:"
read n1
echo $output read n2
answer=$(($n1+$n2))
echo $answer

Plotting your results

% gimp

Experiment capturing a screen dump, opening it and then printing it.

Repetitive tasks

for action in came. saw. conquered.

do

echo 'I ' $action

done
 The variable called action has three potential values. Each value is a
word that is sent to the screen using echo within the do done set of
instructions. The $ sign in front of action assigns its value to be sent to the
screen each time following the word I.

Creating an archive of directories and their contents

When it comes to collating all your directories and their contens into a
single manageable file that can keep a record of the directory structure use
the useful instruction called tar as follows:

%tar czvf tarred_file_name directory_to_archive

A file called tarred_file_name is created. Usually it is best to give your

tarred file a *.tar ending so you can automatically know what type of file it is
in future. In order to open up and generate all directory tree with all its
leaves (which are the files contained within ) use the following command:

% tar xzvf tarred_file_name

If you choose to get ONLY a LISTING of the contents of a tarred file

without rebuilding the directory tree and all its contents you can instead use
the following command:

% tar tzvf tarred_file_name >output_file or if you want to output the

listing to the screen use:

% tar tzvf tarred_file_name

Concatenating files

When you have one files you would like to append to another use the

cat file1 file2 > file3

Seismic Unix
Introduction

Examine a Seismic Data Set

%suxedit SH_geom_2s.su

>g 1 (this graphs the data)

Fourier Transform

>f1 (this graphs the strength of the frequency content at

trace #1)

>f24 (this graphs the strength of the frequency content at trace #24)
 >h (provides help to the user)

All data traces have a "header" that consists of descriptive variables, e.g.
length of the data set, date it was collected etc.

Display

% suximage < SH_geom_2s.su (The < or redirect symbol sends the data
set file into this program)

Bandpass Filtering

% sufilter <SH_geom_2s.su f=0,6,300,400 | suximage

(The | symbol or "pipe" streams the output of one program into the
mouth of the other program)

Exercise
Put all the above instructions into a script called "my_first.sh". Confirm
that this file runs correctly

Notch Filtering

% sufilter < SH_geom_2s.su f=3,6,40,50,60,70,80,180,200

amps=0.,1.,0.5,0,0.5,1.,1.,0

Notes: Verify your filter worked. Run suxedit and plot out the frequency
spectrum to examine whether a notch filter has been applied.

Application: To remove 50-60Hz electronic noise in data

Automatic gain control

In order to adjust for changes in signal strength in time along individual

traces apply the following:

% sugain < SH_geom_2s.su agc=1 wagc=0.05 | suxwigb title="AGC=1

WAGC=0.05s"

% sugain < SH_geom_2s.su agc=1 wagc=0.01 | suximage title="AGC=1

WAGC=0.01s"
 The "clip" value that appears is the amplitude number above which all
your traces are nulled out, i.e. they are assigned a white value, i.e. they are
lumped into a common meanngless value. As a test, why don't you run the
same instructions above, but this time include a clip command, say clip=1
and start varying the value of the clip by orders of magnitude. For example:

%sugain < SH_geom_2s.su agc=1 wagc=0.01 | suximage

title="AGC=1 WAGC=0.01s" clip=1

%sugain < SH_geom_2s.su agc=1 wagc=0.01 | suximage

title="AGC=1 WAGC=0.01s" clip=10

%sugain < SH_geom_2s.su agc=1 wagc=0.01 | suximage

title="AGC=1 WAGC=0.01s" clip=100

Exercise 2
Create a script that .

(1). reads a file SH_geom_2s su, (2) removes frequency below 120 Hz, (3)
applies automatic gain control to compensate for geometric spreading and
(4) plots it to the screen (5) hand in a hardcopy or e-mail me a *.gif file by
the next time we meet..

Iterative tests for filters

Create an iterative set of instructions that will allow you to test the data
set for the best set of filters. The ground roll is in the lower frequency range
(~<= 80Hz). Try at least the following four filter sets:

3,6,20,30 20,30, 60,70 70,80,300,400 80,120,500,600

for filter_values in $filter1 $filter2 $filter3 $filter4

do

sufilter < SH_geom_2s.su f=$filter_values amps=0,1,1,0 | sugain

wagc=0.1 a

gc=1 |suximage

echo 'hello' $filter_values

done

For this exercise you are expected to hand in two scripts and one image.
The first script will show the interactive tests you conducted for different
filters. The second script will show the final set of filters that best remove
the ground roll but keep the rest of the data. An image of the best-filtered
data set with the ground-roll minimized is what I expect to be handed in by
the next lab. Make sure you understand the accompanying linux script
exercise.

Killing bad traces

sukill <SH_geom_2s.su min=16 count=2 >SH_geom_2s_killed.su

min in the number of the first trace to kill and count is the number of
traces starting with min that will be deleted.

Compare the same file before and after traces 16, 17 and 18 have been
removed, e.g.

suximage < SH_geom_2s.su

suximage < SH_geom_2s _killed.su

Reordering traces

susort < SH_geom_2s.su tracf |suximage

tracf is the header value that is used to reorder the traces. In this case
tracf is the trace number. The negative sign implies that the reordered file
will have the traces ordered according to the decreasing value of tracf. So
if tracf = 1,2 ,3. 24 in the input file, tracf=24,23,22,21,20 in the output file

Here we reverse the order of the traces according to their sequential

trace number in the file (tracf) and traces 1 through 24 will be plotted in
reverse order, i.e. 24 through 1.

Cutting out a window of data

If only part of the data set needs to be used, for example only the first
half second we can use suwind

suwind <1001.su key=tracf min=1 max=24 tmin=0 tmax=0.5 |suximage

In this case we have selecting all traces who have values for tracf
between 1 and 24 and all samples between the time 0 s and half a second.

An example that shows how to kill traces, reorder and cut a window of
data from a certain data set
 While logged into your lgc10 account, go to directory ~/refseis13 and
copy over to your local directory the following file: LSU1_1999_TJHughes.tz

Untar this file.

This scripts that follow show the implementation of sukill, susort and
suwind

Examine the file Xamine.sh from the immediately previous section.

Identify susort, suwind, sukill. These S*nix programs are used to prepare a
pseudo-walkaway shotpoint gather for viewing.

Give me the reasons you think why susort was used, why suwind was
used and why sukill was used? Please give me one reason for each. You
will need to image the seismic data to see how the files look BEFORE they
are affected by suwind, susort, and sukill as well as AFTER. The differences
should allow you to see why each program was used and for which reasons..
This is a dropdead date and time with no extensions. You can answer in text
in three to four sentences only. But, you will have to view the data and
perform sukill, susort and suwind. You do not have to send me any images
you created. The reasons you give will show that you understand what
occurred.

How to locate the meaning of each header word

If you want to know what tracf means then type:

% sukeyword tracf

The output will appear on the screen explaining that tracf is the trace
number within the field record

How to change a header word value

Header values are changed according the following formula:

Header value = a + b * (i/j) + c(i/j)

For the following example a,b,c,I and j can be seen to represent:

a first value of each group of traces

b value of increment between traces in a shot gather

c increment in value between the first traces of adjacent shot gathers

 i trace number within the whole file ,e.g. 0,1,2,3,4,5,6,7,8 Note that I
starts at 0. We do not need to set I in the following example.

j number of traces to jump between shots

If b or c are equal to 0, then their products with (i/j) are also equal to 0
and there is no change to the patterc of the header value within adjacent
shot gathers.

For example:

key=offset (shot-receiver distance)

key=sx (x co-ordinate of shot position)

key=gx (x co-ordinate of geophone position)

X is increasing and positive

sx=0 gx= 33 36 39 SHOT #1

offset= 33 36 39

sx=3 gx= 36 39 42
SHOT #2

offset= 33 36 39

sx=6 gx= 39 42 45
SHOT #3

offset= 33 36 39

sushw < filename \

key=sx,offset,gx,fldr,tracf\

a=0,33,33,1001,1 \ # first value of each group of traces

 b=0,3,3,0,1 \ # increment between traces in a shot
gather

c=3,0,3,0,0 \ # increment between first traces of each shot

j=3,3,3,3,3 \ # number of traces to jump between shots

>filename_out

We can simplify the above into several steps:

Step 1: set the sx field of the first 3 traces to 0, the second set of 3 traces
to 3, the third set of 3 traces to 6; i.e. the shot stays at the same place for
whole shot gather and only increments when a new shot is taken (i.e. every
3 traces)

sushw < filename \

key=sx a=0 b=0 c=3 j=3

Step 2: set the offset field of the first shot (first set of 3 traces) to
33,36,39 , the second shot (next set of 3 traces) to 33,36,39, and thelast shot
(third set of 3 traces) to 33,36,39.

| sushw \

key=offset a=33 b=3 c=0 j=3

Step 3: set the X oordinate of the geophone position to 33, 36, 39 for the
first shot; to 36,39,42 for the second shot (next 3 traces), and to 39,42,45 for
the last shot (final 3 traces)

| sushw \

key=gx a=33 b=3 c=3 j=3

In a full script he above 3 steps together can look like:

 #!/bin/sh

set x

filename_in=1000.su

filename_out=1000_geom.su

sushw <$filename_in \

key=sx a=0 b=0 c=3 j=3 \

| sushw \

key=offset a=33 b=3 c=0 j=3 \

| sushw \

key=gx a=33 b=3 c=3 j=3 \

>$filename_out

or we can make a single call to sushw and place the variables together, in
its briefest form:

#!/bin/sh

set -x

filename_in=1000.su

filename_out=1000_geom.su

sushw <$filename_in \

a=0,33,33,1001,1 \

b=0,3,3,0,1 \

c=3,0,3,0,0 \

j=3,3,3,3,3 \

>$filename_out

How to calculate CMP/CDP in header

suchw imilar to sushw but where we use the header values to do the
math:
 value of key1 = (a + b* value of key2 + c * value of key3)/d

We use the formula that CMP = location of source + (source-receiver

offset )/2

suchw < filename_in \ #input file name

key1= cdp \ # output header word

key2= sx \ # first input header word - source x-coordinate

key3= offset \ #second input header word

a=202 \ # first CDP/CMP number * 2

b=2 \ # 2*CMP increment , e.g., 101, 102, 103

c=1 \ # multiplicand for key3

d=2 \ # divisor of all

< filename_out

An example for making CMP values in headers is available from lgc10 at

/home/refseis10/shell_exs/makecmp.sh

How to fix a data set with a variable time delay or a data set that has
false time breaks or how to cross-correlate two traces.

Cross correlation describes the similarity between two time series. For
us a trace consists of a series of amplitude values at regular intervals of time
or a time series. Mathematically, cross-correlation is like convolution, but
where none of the traces are reversed prior to the steps involving shifting,
multiplication and addition (See lecture PowerPoint Presentation entitled
XCor for cross-correlation and the PowerPoint presentation entitled
CMP for convolution, both hyperlinked fromthe main syllabus pertaining
to this class ).
 Lets start by assuming that the geology does not significantly change
from between two adjacent shots. Then, if for one shot gather, the recording
time accidentally starts at a different time with respeect to the shot going off
to that of another shot gather the true delay must be reset. Why?

Well, whereas delay keyword in the headers will have the same value the
data will be at the wrong time. We must change the delay header value so
that the data should appear at the correct time.

NO DATA
NO DATA ???
delay = NO DATA
T
delrt correction ???
0

Once the data is corrected for this wrong delay value then we must make
all the shot gathers have the same length in time starting at tmin=0 (shot
time) You will find however, that before you can do that the data you have
corrected to perhaps a later time now has missing data. What to do???

A worst-case scenario is that the seismograph started recording very late

after the shot went off and that you have irretrievably lost data. What to
do???

To see how this might be done copy to your directory, then modify
accordingly and run the following script that is located in
/home/refseis10/shell_exs/change_delay.sh.

A cross-correlation between traces can be used to estimate the

differences in the times between two identical events. You can see how this
might be done by looking at
/home/refseis10/shell_exs/study_CORRELATION

How to carry out Normal Moveout and Stacking

I recommend that once you have populated your header values for offset
and CDP you should sort the data before sending it to NMO.

susort <file_in >sorted_file_out cdp offset

 A brute stack can be obtained by first trying a constant-velocity stack,
say at 1500 m/s.

You can try various constant velocity stacks at different constant

velocities.

sunmo < sorted_file_in vnmo=1500 |sustack cdp |suximage clip=1

After you obtain an initial brute stack you are ready to start refining
many of your processing parameters. It is during this stage that your sunmo
can read the results of additional velocity analyses. More on that later

AN example of a script containing these instructions, among others, is

available from lgc10 at

/home/refseis10/shell_exs/study_NMO_STACK.sh
Perl (Hoffman, 2001) Practical Extraction and Report
Language

Introduction to Perl

Q. Why use Perl?

There are certainly better ways to write code, but here are my reasons
to use perl:

(1) It costs nothing, is mature and widely available

(2) Testing is quick; on-the-fly. Perl is an interpreted language which

means that code is translated into the machine language while it is running
one line at a time so that places where there are errors are easy to locate.

(3) Perl can easily incorporate shell programming scripts. Perl can be
used as a glue to organize a computational workplace. Perl can be used to
communicate between different modular command-line Open Source
programs.

(4) Perl can be used for more complicated programs that require setting
up functions or sub-routines that help keep complicated programs modular
and simple

(5) Handling text files and their content is carried out more easily than
with other programs

When not to use Perl?

When you want Perl to perform intensive numerical calculations.

When you know of an easier way that will save you time and frustration.

Can I use Perl to make simple, visually interactive programs?

Yes, even using well-known libraries such as GTK, Qt, and of course the
old, classsical Tk interface.
Learning Perl on your own

A great place to start is to use the online tutorials in linux. Use google to
find a Perl tutorial, e.g.: http://www.perl.com/pub/2000/10/begperl1.html

You can also consider subscribing to: http://www.perlmonks.org for free

help and Perl camaraderie.

Also, use Perl itself that comes with documentation. Check this out:

% info perl

perl Perl overview (this section)

perlintro Perl introduction for beginners

perltoc Perl documentation table of contents

Tutorials

perlreftut Perl references short introduction

perldsc Perl data structures intro

perllol Perl data structures: arrays of arrays

perlrequick Perl regular expressions quick start

perlretut Perl regular expressions tutorial

perlboot Perl OO tutorial for beginners

perltoot Perl OO tutorial, part 1

perltooc Perl OO tutorial, part 2

perlbot Perl OO tricks and examples

 perlstyle Perl style guide

Notes: Use the up and down arrow keys to move to the line you want to
select

Control C will get you out of any program

Basic components of Perl

Input and Output

Printing Hello World

In order to give you some courage to start working with this new
language, especially if you have not worked with one too extensively before,
lets consider writing one that is classical across most beginning tutorials
and that provides a stimulating output to the terminal.

#!/usr/bin/perl

#This is my first program in perl

print (Hello World\n\n);

In the above example there are at least five things to note.

(1) The first line denotes the location of the perl binary

(2) From now on all items that are output to the screen will be included
in parentheses and double inverted commas. Double-inverted commas
permit Perl to interpret the different items. For example some items are
read as text and others as special characters when needed. (Try out single
commas just to see what would happen). If you want to null the value of a
special character put a \ before it. For example \\n makes \n come out
just like the characters you see. (Try it out).
 (3) the \n is a shorthand code that means include a new line when the
rest of the text is written out. There is a new line before the start of writing
and there are two new lines after the start of writing.

(4) All lines except the first and the line commented out end with a ;
denoting the end of an instruction. Omission of the ; is a very common
mistake that we all make.

(5) The symbol #on the second line means that these words are
informational for the reader and will not be considered by Perl to be a
meaningful instruction.

Reading from and Writing to a file

If you want to read and write data to hard drive you must first tell the
system you are ready to access a part of the hard drive. This is done by
opening a FILEHANDLE or a file address. You must also provide a name.
The FILEHANDLE should be closed when you are done reading or writing to
the file.

Here is an example of opening a file:

#!/usr/bin/perl

open (FILE, filename) || die (cant open this file $!);

$i=0;

while ($read = <FILE>) {

$line[$i] = $read;

$i=$i + 1;

$imax = $i;

close (FILE);

for ($i=0;$i<$imax;$i=$i+1) {

print ($line[$i]);

}
 $! is a special operator indicating a system error has occurred.

<> is the line-reading operator which continues by itself until the end
of the file is encountered

Line reading continues as long as the value of the while statement is

true, i.e. as long as the content of the parentheses remains TRUE (=1).

Reading is quite straight forward except for the following:

(1) remember that lines of data may have invisible characters that you
may want to remove

(2) you can not read a file unless you know its internal makeup.

Here is an example of writing to a file:

#!/usr/bin/perl

$imax=3;

for ($i=0;$i<$imax;$i=$i+1) {

$line[$i] = $i;

print ($line[$i]);

open (FILE, > filename) || die (cant open this file $!);

for ($i=1; $i<3; $i=$i+1) {

print OUT $line[$i]

close (FILE);

Note that the only important difference between reading and writing is
that we have a redirect sign > before the filename.

Documentation in Perl
 There is another way of documenting perl programs that can later be
used to automatically generate a formatted description of the program to
newcomers. We call this using perlpod, which stands for perls plain old
documentation format, an html-like way of embedding documentation
within a perl script.

For example, here is a version of the same program above with a more
sophisticated and professional documented body. Make sure you leave a
space before the first line that starts with =

#! /usr/bin/perl

=pod

=head1 NAME

My first program

=head1 SYNOPSIS

perl Hello_World2.pl

This is my first program in perl

=head1 DESCRIPTION

Writes a few words to the in terminal

=head1 ARGUMENTS
None

=cut

print("Hello World\n\n");

=pod

=head1 AUTHOR

I am the author of this simple program

=head1 DATE

Sept-16-2013

=head2 TODO

also include lists of items

=cut

=pod

=head3 NOTES
 Although this is just my first program, I can use it as a template with

which to generate documentation in other programs that I write

=cut

print (Hello World\n\n);

Are there any advantges to perlpod ?

(that is, other than keeping notes on HOW the program works for the
next user?)

Yes, there are some advantages to using perlpod that outweigh the extra
time and thought required to place the comments inside your program. One
advantage is that it is relatively easy to convert your documentation (just the
documentation and not the rest of the program) into a different format, such
as PDF, or MSWord.

Data Types

Just as we saw in dealing with shell variables we distinguish between the

value stored on a hard drive and the name associated with that number.

A perl variable is a place to store the value, which is called the literal.

For example:

#!/usr/bin/perl

#This is my another program in perl

$number = 2;

$output_text = (Hello world);

Print (\n$output_text \n\n $number);

When writing out text, note that text consists of individual characters
strung together in a line, including minus signs, plus signs, spaces, tabs,
end-of-line-characters, etc. A string of characters is just that, a string. In
the example above we assign (Hello world) to the variable $output_text.
Lists of Variables (data) or Arrays (vaiable)
If you want to include various lines of texts it might be cleaner to break
up the text into different segments. In order to handle this we can create a
list of lines of text. The list consists of many scalar literals which are
assigned to ordered portions of the array.

#!/usr/bin/perl

#This is my third program in perl

$output_text[0] = (Hello world\n);

$output_text[1] = (I want to live\n);

$output_text[2] = (I want to flourish\n);

Print (\n@output_text \n);

List variables carry the @ sign at the beginning of their name and
will print out their whole content, as in the example above. The list is
ordered starting at 0 and not at 1.

Yes, you could also write the list with a different syntax:

#!/usr/bin/perl

#This is my third program in perl

@output_text = (Hello world\n,I want to live,I want to flourish);

 A list of variables is also known as an array and is identified with the @
symbol:

#!/usr/bin/perl -w

#PURPOSE: describe perl arrays

@output_text = ( Four score\n,and seven years ago\n,our fathers

landed\n);

print(@output_text[2]\n);

print($output_text[2]\n);

print(The number of values in the array is: [@output_text[$#output_text]

+1]\n)

print(The number of values in the array is: scalar(@output_text)\n\n);

Is there a difference between the two outputs?

There are a couple special arrays which will need later when we write
functions and perl programs that can interact with the user, that is they
require input from the user such as a number or a file name on the
command line : e.g.,

%perl sum.pl 1 2

The first variable is called @ARGV and keeps track of the order of the
values that follow the name of the program above (e.g., @ARGV[0], and
@ARGV[1]).

Another special variable @_ is needed to pass arrays to a subroutine (a

sub-program)

Scalars

Scalars are single-value data types. That is, only one value is assigned to
that variable and the value can be a string or a number. Scalars are
indicated by a $ sign at the beginning of the variable.
 There is one special variable in perl that is useful to know. Commonly
you will want to know the number of values your array. The length of your
array or the number of values in your array would be equal to the largest
index plus 1. For this purpose there is a special scalar variable in perl you
can use. This special variable has a literal value equal to the last index in
the array:

#!/usr/bin/perl -w

#PURPOSE: estimate array length

@output_text = (" Four score","and"," seven years ago","our fathers

landed");

$array_size = $#output_text + 1;

print("The number of values in the array is $array_size\n");

print("The last of value stored in the array is:\n);

print(\t\t@output_text[$#output_text]\n");

Note have inadvertently we have introduced, albeit briefly, how to carry

out some simple arithmetic from within perl.

Hashes

Hashes represent pairs of values and their names or keys. Because a

name can be a useful mnemonic for the associated value hashes are very
commonly employed to more easily keep track of lists of parameters. For
example, the following could be a hash for parameters in a seismic data set:

%seismic_data = ('sample_interval_s', 0.001, 'number_of_samples', 1001,

'first_time_s', 0);
 The syntax for retrieval of values from a hash is variable but there is one
form that is preferred because it is easier to read, e.g.,

print ($seismic_data {sample_interval_s}\n\n);

Useful Operations

For-loop/Do-loop in perl

Do-loops (herein for-loop) are a term inherited from Fortran (and bash).
In Perl there is a simple syntax to handle repetitive tasks that is very similar
to C and Fortran, and Matlab. After all, computers ARE supposed to be used
for doing repetitive tasks very fast. Here is how we do a loop:

#!usr/bin/perl

# NAME:

# PURPOSE: To show off for loops

$max = 10;

for ($i=0; $i<=$max; $i++) {

$output_number_array[$i] = $i+1;

for ($i=0; $i<=$max; $i++) {

print ("For index = $i \t value = \t $output_number_array[$i]\n ");

Inside the parentheses, after the for, there are three instructions. The
first instruction $i=0 provides the START of the loop. That is, the first
instruction is the first thing that is carried out in the loop. Remember this!

The second time the loop is run, the third instruction is carried out, i.e.
the $i value is updated by adding 1 to the previous value. At that point the
second instruction must be met for the calculations to enter the loop again.
If the second instruction is not me then the loop is exited and the $i
retains its previous value from the end of the last loop. To be safe, you can
examine the value of $i when the loop is exited.
 Note that we can work the index in reverse as well and that the values of
$i can increment by more than just 1 each time.

Perl operators

Various symbols exist in perl that are very similar to operators in other
programming languages. Operators can be of several types depending on
whether you are dealing with NUMBERS or CHARACTER STRINGS.

Arithmetic

+ addition

- subtraction

* multiplication

/ division

Numeric comparison

== equality

!= inequality

< less than

> greater than

<= less than or equal

>= greater than or equal

String comparison

eq equality

ne inequality

lt less than

gt greater than

le less than or equal

ge greater than or equal

 Boolean logic

&& (and)

and (also and)

or (or)

! (not)

Not (also not)

Miscellaneous

= assignment

. string concatenation

x string multiplication

.. range operator (creates a list of numbers)

Many operators can be combined with a "=" as follows:

$a += 1; # same as $a = $a + 1

$a -= 1; # same as $a = $a - 1

$a .= "\n"; # same as $a = $a . "\n";

Conditional if

An if statement allows perl to pass judgement on two variables. If the

judgement has a TRUE (1) outcome then the instructions inside the curly
braces are carried out, otherwise (FALSE ; =0) the perl language jumps to
the first line after the If statement.

An if statement in its shortest version looks as follows:

#!/usr/bin/perl
 $value[1] = 1.1;

$value[2] = 1.0;

if ($value[1] >= $value[2]) {

print (\You have entered the first set of instructions\n);

else {

print (\nYou have entered the second set of instructions\n);

How to execute System commands in Perl

All that you have learnt prior to perl regarding the linux OS and shell can
still be used within perl. Say, for example you wish to generate the following
working set of directories:

/home/loginID

/data /progs /images

/jpg /tiff

#!/usr/bin/perl

$HOME = (/login/loginID);

$DATA = $HOME

print(\nMaking directories @directory[1] \n);

system ( \\

mkdir p @directory[1] \\
 );

Incorporating SeismicUnix programs into Perl

Example 1

The following script incorporates suplane and sufilt into a perl

document:

We also use a for loop to evaluate a range of possible filtering

parameters. We maintain a 100 % open-pass width of 30 Hz, a decay of
3dB/octave (i.e. doubling of the base frequency) and a doubling of the
second filter value as in what follows:

f=3,6,36,72 f=6,12,42,84 f=12,24,54,108 f=24,48,78,156

Note that (1) we are doubling the second value in each list of filter
parameters: 6,12,24,48; (2) the gap (Hz) between the second and third
values is kept at 30 Hz.

#!usr/bin/perl

# NAME: filter_test.pl

# PURPOSE: To test a variety of bandpass filters

$max_number_ofcases = 4;

=pod

=head3 RULES of OPERATION

$f1 = 3 ;

$f2 = $f1 * 2;
$f3 = $f2 + 30;

$f4 = $f3 * 2;

=cut

=pod

Now we change the f values

=cut

for ($case_number = 1, $f1=3;

$case_number <= $max_number_ofcases;

$case_number++,$f1=$f1*2) {

$f2 = $f1 * 2;

$f3 = $f2 + 30;

$f4 = $f3 * 2;

$filter_parameter_array[$case_number] = ("$f1,$f2,$f3,$f4");

print("This is case number $case_number\n");

print("$filter_parameter_array[$case_number]\n\n");

=pod

Now we go ahead and apply our results

 =cut

$instructions = ("suplane | sufilter \\

f=$filter_parameter_array[1] \\

|suximage & ");

system($instructions);

print($instructions);

Example 2:

Xamine.pl is located in ~/LSU1_1999_TJHughes/seismics/pl/1999/Z,

sostart by copying the directory ~/LSU1_1999_TJHughes and all of its
contents (use cp r).

Project_Variables.pm will have to be modified to agree with the system

path to your personal directories.

Xamine.pl demonstrates both the traditional as well as the object-

oriented mode of incorporating seismic unix programs into Perl scripts.

The required file (sugain.pm) is a Perl package which is written in an

object-oriented version of Perl and pm indicates a package, or collection of
subroutines under a common program name, i.e., sugain.pm.

Each subroutine can be called to set one or more parameters at a time

for the gain functions. There are several advantages to hiding all these
options inside subroutines. First of all, any program can be tailored more to
the particular needs of the user. The user can be made to know less about
what the syntax should be, less about memorizing the parameter names and
more about performing the analysis.

For example, (a) the parameter names can be changed to suit the user
into a term that is more meaningful geophysically and more self-explanatory

(b) package name can also be changed when a new version (or instance)
is used inside a main perl program

(c) The user is NOT REQUIRED to use all the parameters. Parameters
that are not needed do not need to be called. In the seismic unix family of
programs all the parameters have set defaults which the user can not see
but can read about in the manual. As written now, the default values can be
changed by the user inside the package.

(d) if there are any subroutines that bear a similarity in their behavior to
other subroutines the common behavior can be factored out and shared
among the packages. Of course this requires more observant use of the
programming language and more planning ahead of what he different
subroutines do, but this is not too hard to do with Seismic Unix programs
because they are written to be independent of each other for the most part
and so do not share a lot of functionality.

Example 3

Modify Xamine.pl in order to kill bad traces. First you will have to
circumvent one problem:

The large number of files in the data set and the location of the different
traces that need to be killed.

Normally, on the command line you would enter:

sukill <SH_geom_2s.su min=16 count=2 >SH_geom_2s_killed.su

where min in the number of the first trace to kill and count is the number
of traces starting with min that will be deleted.

Incorporate this process into a perl script following the example of

Xamine.pl
Signal Processing a CDP data set using Perl

Killing bad traces

Recall from previous sections that if you open a shell window and enter
the following command you will delete two contiguous traces starting with
trace number 16.

sukill <SH_geom_2s.su min=16 count=2 >SH_geom_2s_killed.su

min in the number of the first trace to kill and count is the number of
traces starting with min that will be deleted.

Compare the same file before and after traces 16, 17 and 18 have been
removed, e.g.

suximage < SH_geom_2s.su

suximage < SH_geom_2s _killed.su

Now we can automate this process further for the case of a real data set
by using the following program to display and select the traces:

Select_tr_Sukill .pl

and the following program to apply the selected traces created by

Select_tr_Sukill.pl :

Sukill.pl

Reversing the polarity of bad traces

On occasions trace polarities are totally reversed, that is some traces in a
shot gather appear to have amplitudes that are equal but opposite in sign to
that of adjacent traces. Normally if we used only Seismic Unix modules for a
gather of 24 traces, we would employ:

suwind < input_filename key=tracef min=1 count=1 | suop op=neg >

temp1

suwind < input_filename key=tracef min=2 count=23 > temp2

cat temp1 temp2 > corrected_filename

If there are many individual files on which to perform this operation if

may more efficient to use the following program which allows a
simultaneous selection of many input file names as an arbitrary range of
traces whose amplitude will be effectively multiplied by -1:

Reverse_Polarity.pl

Frequency-wavenumber filters
Our data set may contain different types of noise: some random and some
coherent. Coherent noise is sometimes identifiable as linear because
unwanted data is aligned along one or more characteristic slopes and the
same noise is superimposed on good signal. In the T-X domain of shot gather
data, a slope has units of s/m or the inverse of velocity. Two-dimensional
integral transformations of these data map into the frequency-wavenumber
domain (f-k). At least one useful simplification occurs during the integral
transformation in that linear X-T data that share a common slope are
mapped on to a linear arrangement as well in the f-k domain. If the slope of
the noise is different enough to that of the good signal the noise and the
signal separate out into distinctly different regions of the f-k domain.

Try to analyze your data for noise, using the following:

Sudipfilt.pl

This collection of Seismic Unix flows will produce 6 plots across your
screen. Their size and location can be changed within the Perl program, so
you are encouraged to copy it over to the equivalent directory in your work
area.

The top row of three panels shows the data before it is f-k filtered and the
bottom row shows them after they have been filtered.

Currently the key dip-filter settings include two sets of four values each.
The inner two values mark the region of interest and the outer two mark the
limiting transition zoneThe units for f-k samples per trace, for ease of use,
rather than m/s. As when we estimated bandpass filter parameters we must
be careful not to choose outer dip-filter values that are too close numerically
to the inner ones or we will generate unwanted noise.

The inner region can be assigned for removal or not as can the region
outside the outer limits. As with the bandpass-filter parameters we do this
by applying weights of 1 through 0 where the following example shows a
zeroing effect to the inner-defined region:

@sudipfilter[1] = ( sudipfilt \\
dt=1 dx =1 \\
amps=1,0,0,1 \\
slopes=-11,-7,-4,-3 \\
);
You should try to change the amp values to the following and note the
differences:

@sudipfilter[1] = ( sudipfilt \\
 dt=1 dx =1 \\
amps=0,1,1,0 \\
slopes=-11,-7,-4,-3 \\
);

Modification of trace header geometry

The following notes are repeated from the corresponding section that
uses only shell scripts to handle seismic data. Hereon, we will use perl
scripts.

Header values are changed according the following formula:

Header value = first_val +

intra_gather_inc * (i % gather_size) +
inter_gather_inc (i %gather_size)

For the following example, first_val, intra_gather_inc,

inter_gather_inc, and gather_size represent:

first_val is first value of the header word for each group of traces or
trace gather.
intra_gather_inc increments the header word value between traces in, for
example, a shot gather.
inter_gather_inc increments the header word between adjacent (shot)
gathers
i trace number within the whole file: e.g., 0,1,2,3,4,5,6,7,8
Note that i starts at 0. We do not need to set i in the
following example.
gather_size number of traces to jump between shots

If intra_gather_inc or inter_gather_inc are equal to 0, then their

products with (i % gather_size) are also equal to 0 and there is no change
to the pattern of the header value within adjacent shot gathers. Note that
the symbol % means modular division. So that for example if the
gather_size = 24, and i = 25,
we would be dealing with the first trace of the next gather of of 24
traces that is:
25 % 24 = 1

A useful example, where we change the values of the headers whose

names are: offset, sx and gx

name=offset (shot-receiver distance)

name=sx (x co-ordinate of shot position)

name=gx (x co-ordinate of geophone position)

 X is increasing and positive, units are in cm. Consider that there are
24 traces per shot gather, although only three of the geophones are shown
here:

sx=0 gx= 450 750 1050

SHOT #1

offset= 450 750 1050

sx=3 gx= 750 1050 1350

SHOT #2

offset= 450 750 1050

sx=6 gx=1050 1350 1650

SHOT #3

offset=450 750 1050

In perl we write the following:

sushw->name(sx,offset,gx,fldr,tracf);

sushw->first_val(0,450,450,1001,1); # first value of each group of

traces

sushw->intra_gather_inc(0,300,300,0,1); # increment between first

traces of each shot gather

sushw->inter_gather_inc(300,0,300,1,0); # increment between traces in

a shot gather

sushw->gather_size(24,24,24,24,24); # number of traces to jump

between shots

We can simplify the above into several steps:

Step 1: set the sx field of the first 24 traces to 0, the second set of 24
traces to 3, the third set of 24 traces to 6; i.e. the shot stays at the same
place for whole shot gather and only increments when a new shot is taken
(i.e. every 24 traces)
 sushw->name(sx);
sushw->first_val(0);
sushw->intra_gather_inc(0);
sushw->inter_gather_inc(300);
sushw->gather_size(24);

Step 2: set the offset field of the first shot (first set of 24 traces) to
450,750,1050 , the second shot (next set of 24 traces) to 450,750,1050,
and the last shot (third set of 24 traces) 450,750,1050, Note that the
offsets DO NOT change between shot gathers.

sushw->name(offset);
sushw->first_val(450);
sushw->intra_gather_inc(300);
sushw->inter_gather_inc(0);
sushw->gather_size(24);

Step 3: set the X coordinate of the geophone position to 450, 750, 105,
for the first shot; to 750,1050,1350, for the second shot (next 24 traces),
and to 1050,1350,1650, ..for the last shot (final 24 traces)

sushw->name(gx);
sushw->first_val(450);
sushw->intra_gather_inc(300);
sushw->inter_gather_inc(300);
sushw->gather_size(24);

In a full script he above 3 steps together can look like:

#! /usr/bin/perl

# instantiation of programs

use SU;
use SeismicUnix qw ($in $out $on $go $to $suffix_ascii $off $suffix_su);
my $log = new message();
my $run = new flow();
my $setheader = new sushw();

# input and output files

my ($DATA_SEISMIC_SU) =
System_Variables::DATA_SEISMIC_SU();
my (@flow);
$file_in[1] = 'All_clean_kill';
$sufile_in[1] = $file_in[1].$suffix_su;
$inbound[1] = $DATA_SEISMIC_SU.'/'.$sufile_in[1];
 $outbound[1] = $DATA_SEISMIC_SU.'/'.
$file_in[1].'_geom'.$suffix_su;

# set up the program with the necessary variables

$setheader ->clear();
$setheader ->first_val(450,0,450,-100,1001,1);
$setheader ->intra_gather_inc(300,0,300,0,0,1);
$setheader ->inter_gather_inc(300,300,0,0,1,0);
$setheader ->gather_size(24,24,24,24,24,24);
$setheader ->name('gx','sx','offset','scalco','fldr','tracf');
$setheader[1] = $setheader ->Step();

# create a flow
@items = ($setheader[1],$in,$inbound[1],$out,$outbound[1]);
$flow[1] = $run->modules(\@items);

# run the flow

$run->flow(\$flow[1]);

# log process
print "$flow[1]\n";
#$log->file($flow[1]);

You probably have noted that the script modifies some other header
values, such as scalco, tracf, and fldr.

scalco=-100 is used to scale the units from meters to cm

fldr and tracf are used as additional counters within the gathers
(tracf=1,2,3,4 24) and between gathers (fldr=1001,1002,1003, etc.)
MATLAB
Create a matrix of numbers

a =[1 2 3 4 5]

a=

1 2 3 4 5

Sin function

>> sin(a)

ans =

0.841470984807897 0.909297426825682
0.141120008059867 -0.756802495307928 -0.958924274663138

Plot function

plot(ans)

Exercise 1 Simple 2D plotting

Make a plot showing at least 10 full sine waves within the plot.

What formula did you end up with? Complete by September 16, at 12.30

Exercise 2: Traveltime Equations

 Plot the traveltime equation for a hyperbola and the direct wave by
modifying the following code: Matlab Code

2
2h
T T
2
0
2

V
Plot for T= 0 to 1 seconds and X= 0 to 1000 m; Plot for V1=1500 m/s

Exercise 3: An ideal seismic wave signaturethe spike -

Seismic temporal resolution is defined, sometimes as the dominant

wavelength.

In reflection seismology, for each reflection interface, a scaled-down copy

of the incoming pulse is returned to the surface. A single interface is
mapped onto a broad pulse of 25(or less) to even 100 ms (or more) in
temporal width:

Plot the following:

 cos( T ) cos(2 T ) cos(3 T ) cos(4 T ) cos(5 T )
cos(6 T ) ....100 terms
Plot for T= - 2s to T = 2s, and amplitude 1 to 1. What happens as the
number or terms increases? (You should be generating a very narrow
function)

Conversely, in the ideal case, what is the ideal frequency content of a

narrow, high-resolution signal. Is it a signal with very high frequency
content or a signal that is narrow? Is it a signal that is both narrow and with
a very high-frequency content

You will need to generate matlab code to answer this exercise. E-mail me
the resulting code. E-mail me your answer to the question as well. Please e-
mail me this exercise by Friday, 23 September at 12.30 p.m. Hint: You will
need to learn how to use the forend construction so that you can
automatically add in 100 terms.

Exercise 4 Constant Phase and Linear Phase

Take the code you generated in exercise 3 and add a constant phase to
each of the frequency components and plot your results. Try adding a phase
value of 90, 180 and 270 degrees to each of the frequency components. Plot
each case. MATLAB CODE . Now add the subplot(2,2,2) case where phase
value = 360 degrees.

Take the code you generated in exercise 3 and add a phase to each of the
frequency component that is linearly dependent on frequency. Us the
following relation:

Phase = m * frequency where m has units of radians/Hz MATLAB CODE.

Finally, add one subplot (2,2,2) for the case where m = 40pi radians/100Hz.

You can use the links to my matlab code as an example.

MATLAB CODE

Exercise 2- Matlab

% hyperbola

% Author: Juan M. Lorenzo

% Date: Sept. 13 2005

 % Purpose: examine differences in Hyperbolic geometries

% as a function of thickness and velocity for a

% horizontal single layer case

xmin=-500;xmax=500;tmin=0;tmax=0.4

x=xmin:1:xmax; % meters

%first plot

%plot hyperobola

% traveltime equation for a hyperbola

V1=1500;

h = 50; %thickness of first layer in meters

for V1=1000:500:3000

t=sqrt(x.*x ./V1./V1 +1/V1/V1*4*h*h);

subplot(2,1,1)

plot(x,t)

axis ij;

hold on

end

Matlab code for exercise 4 A study of the effects of constant phase

and linear phase on a seismic wavelet

CONSTANT PHASE SHIFTS

% looking at constant phase

t= -pi/4:.001:pi/4;

clear Amplitude

phase = 0 ;

Amplitude = cos(2 .* pi .* 1 * t + phase);

for freq=2:100
 Amplitude = (Amplitude + cos (t .* 2 .* pi .* freq + phase));

end

subplot(2,2,1)

length(freq)

plot(t,Amplitude/100)

title('phase=0 deg')

phase = pi/2;

clear Amplitude

Amplitude = cos(2 .* pi .* 1 * t + phase);

for freq=2:100

Amplitude = (Amplitude + cos (t .* 2 .* pi .* freq + phase));

end

subplot(2,2,2)

plot(t,Amplitude/100)

title('phase=90 deg')

phase = -1 * pi ;

Amplitude = cos(2 .* pi .* 1 * t + phase);

for freq=2:100

Amplitude = (Amplitude + cos(t .* 2 .* pi .* freq + phase));

end

subplot(2,2,3)

plot(t,Amplitude/100)

title('phase=-180 deg')

Matlab code for exercise 4 A study of the effects of constant phase

and linear phase on a seismic wavelet
LINEAR PHASE SHIFTS

% looking at phase

t= -pi/4:.001:pi/4;

clear Amplitude

m = 20*pi/100

m = 0;

Amplitude = cos(2 .* pi .* 1 * t + m * freq);

for freq=2:100

Amplitude = (Amplitude + cos (t .* 2 .* pi .* freq + m * freq));

end

subplot(2,2,1)

length(freq)

plot(t,Amplitude/100)

title('phase slope = 0 rad/Hz')

m = pi/2 /100;

clear Amplitude

Amplitude = cos(2 .* pi .* 1 * t + m * freq);

for freq=2:100

Amplitude = (Amplitude + cos (t .* 2 .* pi .* freq + m * freq));

end

subplot(2,2,2)

plot(t,Amplitude/100)

title('phase = pi/2 rad/100 Hz ')

 m = 20 * pi /100 ;

Amplitude = cos(2 .* pi .* 1 * t + m * freq);

for freq=2:100

Amplitude = (Amplitude + cos(t .* 2 .* pi .* freq + m * freq));

end

subplot(2,2,3)

plot(t,Amplitude/100)

title('phase = 20pi rad/100 Hz')

Extras

Creating a shell script to log in automatically

% gedit a file called start.sh

(CONTENTS of start.sh: )

--leave xedit

%chmod 755 start.sh (make this file executable)