zad1

pkg load signal;
wf=wrecfile('ICP_sygnaly/046050jg.dta'); % otwarcie pliku
numer_segmentu = 3;

segment = wrec_dsgm(wf, numer_segmentu); %struktura zawierajaca naglowek->ciekawe informacje o segmencie

zbieranie = []; %inicjacja tablicy gdzie bede zbieral probki do liczenia mocy

moc_tablica = 0;
dlugosc_wektora_mocy = 1; % to nie jest dlugosc co ile mocy zostalo uwzglednionych

stopien_decymacji = round(wrec_freqsmp(wf)/0.25);
rzad_filtra = 660;

b = fir1(rzad_filtra, 0.9*1/stopien_decymacji, "low");


dlugosc_bloku = stopien_decymacji*2; % by decymacja za kazdym razem dawala 2 probki
stan1 = zeros(length(b)-1,1);

for(probka_startowa = 0:dlugosc_bloku:segment.nsmp-dlugosc_bloku-1)
  blok_probek = wrec_rdsig(wf,numer_segmentu,probka_startowa,dlugosc_bloku,{'#ICP'});
  %plot(blok_probek);

  %jak dlugi powinien byc wektor zer ktory bedzie stan1 zeby je najpierw zainicjowac
  %czy moze byc skalar zero, jedno zero jak
  [filtracja1, stan1] =  filter(b,1, blok_probek, stan1);
  wyjscie_decymacja1 = filtracja1(1:stopien_decymacji:length(filtracja1));  %komenda decimate od razu filtruje filtrem czebyszewa wiec nie musze wczesniej filtrowac
  %plot(wyjscie_1decymacja);

  %nie musze sprawdzac czy mam wystarczajaco duzo probek bo tak ustawilem dlugosc bloku by miec za kazdym razem wystarczajaco duzo probek



  %znowu nie musze sprawdzac czy uzbieralo mi sie wystarczajaco duzo probek bo tak ustawilem dlugosc bloku by za kazdym razem miec 4 probki



  zbieranie = [zbieranie; wyjscie_decymacja1];

  while(length(zbieranie) >= 60)  %while na wypadek gdyby zebralo sie wiecej niz 90 probek choc wiem, ze w moim przypadku nie jest to mozliwe
          transformata = abs(fft(zbieranie(1:60)));
          moc = transformata(1).^2+transformata(2).^2+transformata(3).^2+transformata(4).^2; %policzylem, ze w tych probkach znajduje sie moc ktora nas interesuje
          moc = moc./60./2;
          moc_tablica(dlugosc_wektora_mocy) = moc;

          dlugosc_wektora_mocy = dlugosc_wektora_mocy + 1;
         zbieranie = zbieranie(31:length(zbieranie));   %wyrzucam polowe probek ktore wykorzystalem
  endwhile

end

czas = (1:length(moc_tablica)).*2; % liczymy transformate z 60 probek ale liczymy ja co 30 probek a 30 probek odpowiada 30*4s = 120s = 2 min
plot(czas,moc_tablica);
title('fale B');
xlabel('czas [min]');
ylabel('moc');

wrec_close(wf);