# Creating a phase-velocity diagram with matlab

1. Apr 22, 2012

### Bugeye

1. The problem statement, all variables and given/known data
I'm given a time vector and y-values of a wave. The plot looks like this:

2. Relevant equations
I have been calculating the velocity using this equation:
$c=\frac{ω Δx}{Δ\Phi}$
Where $\Delta\Phi$ is the phase change between the phase that you calculate using the fourier transform and Δx is the space between the sensors that collected the data

3. The attempt at a solution
Code (Text):
function [ f velocities ] = dCurve( data, geophoneSpacing, freqLimit )
% Makes the dispersion curve given the raw geophone data, sampling frequency
% and geophone spacing. The resulting curve is plotted.

%Split up the data into useful variables.
%This can be ignored for the most part, this is just reformatting the data
%The results are the following variables:
%t-> the time vector for the equations
%sensorData-> a matrix that has the same number of columns as t and a row for each sensor's raw values
%sensorCount -> The # of sensors
%Fs -> the sampling frequency
datasize=size(data);
t=data(2:datasize(1));
sensorData=data(2:datasize(1),2:datasize(2))';
sensorCount=size(sensorData);
sensorCount=sensorCount(1);
Fs=1/(t(2)-t(1));

%Get the fourier transform
for i=1:sensorCount
[f, rawFft(i,:)]=dfft(Fs, sensorData(i,:));
end
absFft=abs(rawFft);
reFft=real(rawFft);
imFft=imag(rawFft);

%Get the phase between each sensor
for i=1:sensorCount
phi(i,:)=angle(rawFft(i,:));
end

%Get the phase change between each sensor
for i=1:(sensorCount-1) % Go through each pair of sensors. The pairs consist of sensors i and i+1
dPhi(i,:)=phi(i,:)-phi(i+1,:);
dPhi=dPhi-(2*pi*(dPhi>pi));
dPhi=dPhi+(2*pi*(dPhi<-pi));
end

%Calculate velocity
for i=1:(sensorCount-1)
w=2*pi()*f;
c(i,:)=(w.*geophoneSpacing)./(dPhi(i,:));
end

f=f(f<=freqLimit);
velocities=c(f<=freqLimit);

end

dfft is a function that creates a fft and a domain and returns it. The thing works fine if both of the equations are exactly the same except shifted a bit. For example, this works. Well, except for that dip at the end:
The initial functions:

The result after running them through the script

Unfortunately, it doesn't work as well as I would like. For example, those sine waves from above don't work. I get the following results:
Original:

Result:

I'm not sure if I'm using the wrong equations, if my code is wrong, or if I'm interpreting the graphs wrong. Can anyone help?

Last edited: Apr 22, 2012