Документ взят из кэша поисковой машины. Адрес оригинального документа : http://lnfm1.sai.msu.ru/grav/english/lecture/filtering/Lab1.doc Дата изменения: Tue Jan 13 23:20:26 2009 Дата индексирования: Mon Oct 1 23:38:57 2012 Кодировка: Поисковые слова: cosmic rays

Theory of filtering and time series processing
L. Zotov
Laboratory work ?1

Purpose: Studying discrete Fourier transform with use of Matlab.

Matlab software - powerful tool for calculations and Geophysical
applications.
advantages - huge amount of functions and toolboxes, it's possible to write
programs and incorporate them into C++ programs
disadvantages - calculation speed, programs written in Matlab, perform
longer, than C or Fortran programs. They need optimization.

You can get free Matlab software and other programs for the 2008-2009
school year is OSU, filling License agreement and bringing empty DVD-R to
IT Service Desk in 025 Central Classroom
http://oit.osu.edu/site_license/slwin.html#matlab

If you are not familiar with Matlab, the best way is to run it and press F1
for Help. Tens of thousands of pages of information on its use, functions
and toolboxes will be at your disposal.


1) First part of the lab - generating the signal and it's fast Fourier
transformation.
Download program for the 1 lab at
http://lnfm1.sai.msu.ru/grav/english/lecture/filtering/
%---------------------------------------------------
clear;

N_signal=1024;
% generating two-sin signal
for (k=1:1:N_signal)
signal(k)=sin(2*pi/10*(k-1))+sin(2*pi/100*(k-1));
end;

plot(signal);

%fast Fourier transform calculation
Ftrns_signal=fft(signal);

%amplitude spectrum calculation
apl_spectrum=abs(Ftrns_signal);

plot(apl_spectrum);

















%-----------------------------------------------------
2) Second part of the lab - reading and processing the real signal
Read your data or download EOP C01 from
http://hpiers.obspm.fr/eop-pc/
http://hpiers.obspm.fr/iers/eop/eopc01/eopc01.1900-now
Remoove the header in the file.
% reading the signal from the file
cd C:\work\course\filtr\Lab1;

fin=fopen('eopc01.1900-now','rt');
A=fscanf(fin,'%f',[11 inf]);% A - array of data
fclose(fin);

%determining the size of the signal
l=size(A);
N=l(2);

%selecting the rows of the Array
Date=A(1,1:N);
X_pole=A(2,1:N);
Y_pole=A(4,1:N);
dt=Date(2)-Date(1);

plot3(X_pole,Y_pole,Date)

plot(X_pole);

%selecting the length of the part of the signal will be trasformed
N_ft=N;
%fast Fourier transform calculation
Ftrns_X=fft(X_pole,N_ft);


% frequency calculation
% for the half of the spectrum not counting 1 coefficient - an avarage
N_half=N_ft/2-1
freq=(2:N_half)/N_ft/dt;


% frequency calculation
% for the half of the spectrum not counting 1 coefficient - an avarage
% N_ft is odd or even - ?
N_half=floor(N_ft/2-1);
freq=(2:N_half)/N_ft/dt;

%amplitude spectrum calculation

apl_spectrum_X=abs(Ftrns_X)/N_ft;

% we plot spectrum multiplied by to to take into account energy in the
% second part
plot(freq, 2*apl_spectrum_X(2:N_half))
title('amlitude spectrum - module of Fourier-transformation ')
xlabel('frequency')


.