%% MATLAB alapvető parancsok 
%Írta: Dr. Nagy Péter (pnagy@hds.bme.hu)
%Ha nem érthető, félreérthető valami, kérlek jelezd, hogy a hallgató társaidnak már könnyebb dolga legyen.
%A dokumentum elolvasása nem váltja ki a Matlab Onramp kurzus elvégzését (matlabacademy.mathworks.com) 
% A függvények leírása közel sem teljes körű. A cél, hogy képben legyünk a függvényekkel, amikre szükségünk lehet. A matlab help sokkal részletesebb leírást és példákat tartalmaz minden egyes esetben.
%Legutóbbi frissítés: 2019.07.24

%% Változó definiálása
%Egyszerű szám definiálása a következő módon történhet.
EgySzam=5;
%Ilyenkor az EgySzam nevű változó létrejön és 5-ös értéket vesz fel. A változó elnevezése bármi lehet, de nem kezdődhet számmal, nem tartalmazhat ékezetes vagy különleges karaktert. Illetve vannak
%előre definiált változók, függvények sem célszerű használni, mert azzal felülírjuk.  
%Például ilyen a pi és az i. 
AKorKeruletenekEsAtmereojenekHanyadosa=pi;
ImaginariusEgyseg=i;
%Nem javasolt pi=2; ezzel a pi- változót felül írnánk, ami elég komoly gondokat okozna a programok működésében
%A sorok végén a pontosvessző nem kötelező. Ha kitesszük, az eredmény nem írodik ki a konzolba. Ha nincs ott, akkor igen. Vegyük
%észre, hogy a változók a workspaceben is létrejönnek, ott le tudjuk ellenőrizni az értéküket (ha mátrix vagy vektor akkor kétszer rákkatintva) és típusukat.

%Egy változó lehet karakter (char) vagy szöveg (string) is, ezt aposztróffal jelöljük:
ValamiSzoveg='Szöveg';

%Vannak sor, oszlop vektorok és mátrixok is. A sorokat pontosvesszővel, az
%oszlopokat vesszővel vagy space-szel lehet elválasztani
OszlopVektor=[1;2; pi; 1.56];
SorVektor=[1 2 pi 1.56];
SorVektor2=[1,3, pi, 4];
EgysegMatrix=[1 0 0; 0 1 0; 0 0 1];

%Vektorok közti műveletnél nagyon fontos, hogy mit is szeretnénk csinálni.
%Ha vektort osztunk egy számmal az triviális:
SorVektorFele=SorVektor/2;
%Azonban gyakran előfordul, hogy két vektor elemeit szeretnénk elosztani.
%Ha nem az egész vektort szeretnék elosztani, hanem annak csak az elemeit, (nem vektor osztás, hanem elemenkénti osztás) 
%akkor ponttal tudjuk ezt jelezni a matlabnak

KetVektorHanyadosa=SorVektor/SorVektor2; %matematikai szabályok érvényesek, pl A*B esetén A oszlopainak és B sorainak a száma egyenlőnek kell lennie.
KetVektorElemeinekHanyadosa=SorVektor./SorVektor2; % két vektor vagy mátrix elemenkénti osztása. A két vektornak vagy mátrixank teljesen egyforma méretűnek kell lennie.



%% Alap függvények
%A függvények a változóktól eltérően nem tárolnak egy értéket, hanem azzal
%műveletet hajtanak végre. Ez a művelet egyszerűbb műveletek kombinációjából tevődik össze. Például a szinusz függvényt senki nem szeretné leprogramozni, hiába tanulta meg, hogyan lehetne sorfejtéssel előállítani.
%például a szinusz 5 fokot így számíthatjuk ki:
Szinusz5fok=sin(5/180*pi);
%Ennek a műveletnek egy bemenő változója volt a 5/180*pi és egy kimenő változója.
%A matlabban nagyon sok ehhez hasonló előre elkészített függvény, melynek
%leírását az F1 gomb leütésével hívhatjuk elő. Itt pontosan leírják, mik
%lehetnek a bemő és mik a kimenő változók. Például a sin() függvénybe egy
%vektort vagy mátrixot is írhatunk, akkor az a kimenet is ennek megfelelő
%lesz.
SokSzogSzinusza=sin([1 2 3 4 5]/180*pi);
SorVektorSzinusza=sin(SorVektor);

%Gyakori függvények:
%-Gyökovnás: sqrt()
%-Hatványozás: aAbediken=a^b; Ha a és b vektor, akkor elemenként célszerű a
%műveletet elvégezni a.^b
%-Szögfüggvények: sin(), cos(), tan(), asin(), acos(), atan()
%-Hiperbolikus függvények sinh(), cosh(), ...

%% Vektorok, mátrixok létrehozása
%Egyenlő osztásközű vektort a kettős pont segítségével hozhatunk létre
x0=1; %kezdő érték
dx=0.1; %különbség
x1=4; %maximum eddig az értékig állítjuk elő a vektort. Nem biztos, hogy ez az érték benne lesz.
EgyenloosztaskozuVektorKulonseggel=x0:dx:x1;

%Hasonlóan egyenközű vektort hozhatunk létre a linspace paranccsal
x0=1; %kezdő érték
n=20; %a létrehozandó vektor elemeinek száma
x1=8; %eddig az értékig állítjuk elő a vektort. Ez az érték is biztosan benne lesz.
EgyenloosztaskozuVektorDarabszammal=linspace(x0, x1, n);
%Ez egy jó példa több bemenetű (jelen esetben 3) és egy eredményt adó (a
%vektor) függvényre.
%Matlabban gyakori, hogy nem kell minden bemeneti változót megadni. Például
%az alábbi függvény is létrehozza a vektort az Matlab által alapnak
%tekinttet 100 darab elemmel. 
EgyenloosztaskozuVektor100darabbal=linspace(x0, x1); %Ez úgyanazt az eredményt adja, mint a linspace(x0, x1, 100);

%Az oszlop és sor vektorok egyszerűen transzponálással csinálhatúak
%egymásból, amit a ':
OszlopVektor2=SorVektor'; %a pontosság kedvvért ez az adjungálást jelenti, ami komplex számok esetén a képzetes rész előjelét is megfordítja. A tényleges transzponálást a .' paranccsal lehet előhívni.

% A következő függvények vektorok, mátrixok létrehozására alkamasak. Egy
% bemenet esetén egy négyzetes mátrix sor és egyszere oszlop számát adjuk
% meg
%Például 
ones(3)
%egy 3x3-as csupa egyesből álló mátrixot hoz létre.
% Ha két bemenetet adunk meg, akkor az első a sorok a második mindig az
% oszlopok számát jelöli
%Például
ones(1,3)
%egy 3 oszlopból álló sorvektort hoz létre.
%A függvények
%ones() csupa 1-esből álló vektort, mátrixot hoz létre
%zeros() csupa 0-ból álló vektort, mátrixot hoz létre
%rand() egyenlő eloszlású véletlenszámokból álló vektort, mátrixot hoz létre
%randn() normás eloszlású véletlenszámokból álló vektort, mátrixot hoz létre

%% Vektorok, mátrixok elemei
% Részletes leírás itt: https://www.mathworks.com/company/newsletters/articles/matrix-indexing-in-matlab.html
%Pár példa
v=rand(1, 5);
%vektor harmadik eleme 
v(3)
%a v harmadik és ötödik eleme
v([3,5])
%a v első 3 eleme
v(1:3)
%v utolsó 3 eleme
v(end-2:end); %vagy v(3:5)

A=rand(3,4);

% A mátrix első sorának  harmadik eleme
A(1, 3)

% A mátrix második sorának utolsó két eleme
A(2, (end-1):end)

%% Görbe (egyenes) illesztés
%lineáris görbe illesztés: https://www.mathworks.com/help/matlab/data_analysis/linear-regression.html
%további tutorialok: https://www.mathworks.com/help/curvefit/linear-and-nonlinear-regression.html





%% Diagram rajzolást segítő parancsok
x=linspace(0, 2*pi, 200); %x egy vektor, ami 0 és 2pi között egyenközűen felvett 200 db pontot tartamaz
f=sin(x); %f egy vektor, mely az x pontok szinuszát tartalmaza
f2=cos(x);

%Ábrázolás előtt érdemes egy üres rajzolási felületet létrehozni a 
figure();
%paranccsal. Ekkor egy új üres ablak jelenik meg. (ha többször lefutatjuk, akkor mindig egy új ablakot nyit meg.)
%Ábrázolni a plot parancsal tudunk. Alapvetően két vektort kell megadnunk, az egyik a függvény x a második az y koordinátákat jelöli.
plot(x, f);
% ha felül szeretnénk írni az ábrát, használjuk újra a plot parancsot.
% Például a koszinusz függvény kirajzolására
plot(x, f2);
%Ha egyszerre több függvényt szeretnénk ábrázolni, bekapcsolhatjuk, hogy az
%új ábra ne írja felül a régit. Ezt a
hold on
%paranccsal tehetjük meg. Ezek után lefutatva a 
plot(x, f);
%parancsot, a szinusz és a koszinusz függvény is látszik.

%tengelyfeliratokat a vizszintes tengelyre a 
xlabel('x')
%függvény segítségével, a függőlegesre pedig az
ylabel('f (x)')
%paranccsal hozhatunk létre. Mindkét függvény egy egyszerű stringet vár
%bementre ezért kell aposztrófokat használni.

%ha több görbét ábrázoltunk érdemes a jelölés egyértelműsíteni.
%Ezt a
legend('cos(x)', 'sin(x)')
%paranccsal tehetjük meg. A függvény annyi string bemenetet vár, ahány
%görbét már felrajzoltunk. A bemenetek sorrendje meg kell, hogy egyezzen.
%Mivel a törlés után először a koszinuszt és csak utána rajzoltuk fel a
%szinuszt, ezért első bemenet a 'cos(x)' és második a 'sin(x)'.

%Ábra kimentése
% az ábra mentéséhez először a rajzlapunk méretét kell fizikai méretre
% beálltani:
px=7.5;  %az ábra kívánt szélessége
py=5;     %az ábra kívánt magassága
% az aktuális "figure" méretének beállítása
set(gcf, 'PaperPositionMode', 'manual');
set(gcf, 'PaperUnits', 'centimeter');
set(gcf, 'PaperPosition', [0 0 px py]);
set(gcf, 'PaperSize', [px py]);

%Végül a print parancsal a kívánt fájlnévvel kimenhetjük az ábrát png ('-dpng')
%formátumba 600dpi ('-r600') felbontással
print('test_abra','-dpng', '-r600' )
%Ekkor a munkakönyvtárban megjelenik egy test_abra.png ami a digramunkat
%tartalmazza.

%A digramrajzolás megannyi további beállítási lehetőséget rejt, ezek nagy
%része a plot prancsnál a helpben bővebben le vannak leírva.
%% Jelfeldolgozás
%készítsünk az eddgiek alapján egy szinuszos jelet:
%az időlépések legyenek:
t=linspace(0, 20, 2048).'; % a .' jel a végén a sor vektort oszlop vektorrá transzponálja
% a választott körfrekvenciák
omega1=2;
omega2=5;
%választott jel amplitudók
A1=2;
A2=3;
tiszta_jel=A1*sin(omega1*t)+A2*cos(omega2*t);
%adjunk hozzá a jelhez egy 0 várható értékű, normális eloszlású zajt
A_zaj=0.1;
zaj=A_zaj*randn(size(t)); % a size(t) függvény visszaadja, az idővektor méretét, a randn() pedig ebben a méretben fog létrehozni normális eloszlású zajt

% vizsgált jelünk legyen a periodikus jel és a zaj összege.
vizsgalando_jel=tiszta_jel+zaj;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% szűrés
%Savitzky–Golay szűrés
%az előadáson is tanult egyik legismertebb szűrő a Savitzky–Golay szűrő
% ha nincs szükség a deriváltra, akkor az újabb Matlabban használhatjuk az
% sgolayfilt függvényt. Ennek a jelet, az illesztett polinom fokszámát és
% az ablak hosszát kell megadnunk
fokszam=3;
ablak_hossz=51; %csak páratlan lehet
szurt_jel1=sgolayfilt(vizsgalando_jel, fokszam, ablak_hossz);
figure();
%rajzoljuk ki az eredeti és a szűrt jelet:
plot(t, [vizsgalando_jel, szurt_jel1]);

%Savitzky–Golay szűrés
% Ha a deriváltra is szükségünk van nekünk kell elkészíteni a szűrő
% együtthatóit és azt "eltolva rászorozni" a jelre (bővebben előadáson). 
% A szűrőt a következőképp készíthetjük el:
[~,szuro] = sgolay(fokszam,ablak_hossz);
%A szűrő segítségével, így a szűrt jelet és annak deriváltjait is előállíthatjuk, (a derivált nem lehet nagyobb az illesztett görbe fokszámánál) 
% a szűrt jel előállatása:
p=0;
szurt_jel2 = conv(vizsgalando_jel, szuro(:,p+1), 'same'); %szűrt adat.  A con függvény eltolva rászorozza az jelünkre a szűrőt. Mivel szimmetrikus az ablak, a jel elején az ablaknak kellenének értékek korábbról, mint ami rendelkezésre áll. A 'same' opció azt jelenti, hogy ide 0 elemeket feltételez a program. A másik lehetőség a 'valid' ilyenkor, csak azokat az értékeket adja vissza, ahol a szűrő már tényleg "működik". Így viszont rövidebb lesz a jel
%rajzoljuk rá a korábbi ábránkra
hold on;
plot(t, szurt_jel2); % így a jel elejét és végét leszámítva, működik a szűrés

% a szűrt jel deriváltjának előállatása
p=1;
dt=t(2)-t(1); %időlépés vagy mintavételezés ideje
szurt_jel_derivaltja = conv(vizsgalando_jel, factorial(p)/(-dt)^p * szuro(:,p+1), 'same'); %szűrt adat


%periodikus jelek feldolgozásására
% aperiodikus jelekre a gyors Fourier-transzformációt szoktuk alkalmazni
%ehhez a tanszék írt egy függvényt, ami rögtön az amplitúdókat és a
%frekvenciákat adja vissza, ez a függvény az ext_real_fft, ami a tárgy
%honlpajáról letölthető. A függvénynek az aktuálisan használt mappában kell
%lennie, különben nem találja meg a matlab.
[freqs, amps, phases] = ext_real_fft( t, vizsgalando_jel);
%freqs a frekvenciák, amps az amplitudók és phases a fázisok vektora.
%plotoljuk ki ezeket is egy üres lapon
figure()
plot(freqs, amps)


% keresztkorreláció

%hozzuk létre a jel eltolt verzióját, egy másik zajfüggvényt ráadva
idokeses=0.3;
 vizsgalando_jel2=A1*sin(omega1*t-idokeses)+A2*cos(omega2*t-idokeses);
zajos_jel2=vizsgalando_jel2+0.1*randn(size(vizsgalando_jel2));

% a jelfeldolgozásban a kersztkorreláció alatt azt értjük, amit
% matematikában általában keresztkovarianciának hívnak. Hívjuk meg a matlab
% beépített függvényét, és kérjük, hogy az előadáson tanult képletnek
% megfelően normálja (0-1 közé) a kimenetet.

[R,lags] = xcov( vizsgalando_jel, vizsgalando_jel2, 'normalized');
%R hogy mennyire egyezik
%lags% mekkora idő eltolásnál. Mivel a Matlab egységnyi időben tud
%számolni, a lags-ot meg kell szorozni dt-vel hogy megkapjuk a tényleges
%időeltolásokat

figure()
plot(lags*dt, R) %kirajzolja a keresztkorreláció értékét az időeltolás függvényében
%ahol maximuma van ott egyezet legjobban a két jel

%autókorreláció
%hasonló a keresztkorrelációhoz, csak itt a jelet saját magával hasonlítjuk
%össze ismétlődéseket keresve

[R_auto1,lags_auto1] = xcov( vizsgalando_jel, vizsgalando_jel, 'normalized');
%vagy a második bemenet felesleges, ilyenkor automatikusan autokorrelációt
%feltételez a program
[R_auto2,lags_auto2] = xcov( vizsgalando_jel, 'normalized');


%% NI adatgyűjtő parancsok
devices = daq.getDevices; %kilistázza az elérhető eszközök neveit
%a devices egy strukturából álló vektor lesz, a n-edik eszköz nevét a devices(n).ID változó adja vissza.


s= daq.createSession('ni'); %ni (National Instruments) eszközökhöz készít egy munkafolyamatot (sessiont)
%Ehhez a munkafolyamathoz tudunk különböző mérési csatornákat adni.
%Egy munkafolyamathoz több mérőkártyából is választhatunk csatornákat. Akár több sessiont, munkafolyamatot is létrehozhatunk. Például célszerű a digitális és az analóg csatornákat külön munkafolymatban kezelni.

 ch1=addAnalogInputChannel(s,devices(1).ID,'ai0','voltage'); %egy analog mérőcsatornát adunk az s munkafolyamathoz. A függvény kimenete az új csatorna elérése. Ezt nem szükséges megadni, akkor is hozzáadódik a munkafolyamathoz. Azonban, ha a csatorna beállítsait módosítani szeretnénk akkor hasznos, ahogy azt majd később látjuk.
%A függyvény első bemenete a létrehozott munkafolyamat neve.
 % A függvény második bemenete a mérőeszköz neve. Alapból 'Dev0','Dev1'... hasonló neveik vannak a kártyáknak. Ha lefuttattuk a 
 %devices = daq.getDevices;  parancsot, akkor az első mérőkártya neve devices(1).ID.
 % A függvény 3-dik bemenete a mérőkártya csatornájának megnevezése. Az 'ai0' a 0-s csatornát jelöli. Az 'ai0','ai1','ai2', 'ai3' csatornák megadása esetén alapból differenciál mérést feltételez a program. Pl. 'ai0' esetén a 0-s és a 4-es bemenetek közti feszültséget méri majd a kártya. (Ez zajelnyomás szempontjából ideális. Hosszú vezetékek alkalmazása esetén mindenképp ilyen mérést kell használni.)
 % Az 'ai4',..,'ai7' beállítások estén 'SingleEnded' mérést feltételezünk, ami azt jelenti, hogy a földhöz képest mér feszültséget a program. Pl 'ai5' esetén az 5-ös bemenet és a föld közti feszültséget méri. (Ez zajelnyomás szempontjából nem ideális, azonban a differenciál méréshez képest nem kettő, hanem csak 1 csatornát használ.)
 %Ha meg szeretnénk változtatni, akkor a csatorna konfigurálását a 
 %ch1.TerminalConfig változó megváltoztatásával tehetjük.
 ch1.TerminalConfig = 'SingleEnded'
 %vagy
 ch1.TerminalConfig = 'Differential'
 % Az addAnalogInputChannel függvény negyedik bemenete a mérendő mennyiséget adja meg. Ha nem speciális kártyánk van, akkor ez általában a feszültség, 'voltage'.
 %Egy munkafolyamathoz több csatornát is hozzáadhatunk. Akár több eszközről is.
 
v_meres=inputSingleScan(s); %Az s munkafolyamat csatornáin hajt végre egyszeri mérést. Ez hasznos leeht teszteléshez. Illetve, ha ritkán akarunk mintát venni. A kimenet egy sor vektor, melyben az oszlopok száma megegyezik a csatornák számával. Az elsőként hozzáadott csatorna lesz az első, a másodikként a második és így tovább.

 %A munkafolymatnak két legfontosabb beállítása a mintavételezési frekvencia és a mérés időtartama.
 %A frekvenciát a munkafolyamat s.Rate változójának átírásával tehetjük meg. Pl  5000 Hz mintavételezést a
 s.Rate=5000; %pranccsal tudunk megadni.
 %A mérés hosszát pedig a s.DurationInSeconds változó beírásával. Pl 3 másodperces mérést az
 s.DurationInSeconds=3; %parancs segítségével adhatunk meg.

[M_data,timeStamps] = startForeground(s); % Ez az egyik legfontosabb parancs. Ha az "s" munkafolyamatot beálíltottuk, akkor ezzel a paranccsal indítjuk a mintavételezést. A mérés ideéig (s.DurationInSeconds) a program elfoglalt, semmilyen parancsot nem hajt végre. A mintavételezés befejeztével a mérési eredményeket az első kimenetén (M_data), a mérésekhez tartozó időt a timeStamps vektor tartalmazza. Az M_data egyes sorai tartoznak egy időpillanathoz (pl.: a 10-edik időpillanathoz tartozó mérési eredmény M_data(10,:)), az egyes oszlopok pedig az egyes mért mennyiségekhez (pl.: M_data(:, 1) a munkafolyamat első csatornájának mérési eredményi különböző időpillnatokban)  

%% Digitális csatorna
% A hőmérséklet mérésnél szükségünk van digitális csatorna használatára. Érdemes ezt egy új munkafolyamatot indítani a digitális csatornának. Például a laborban használt kártyák nem tudnak egyszerre a digitális és az analog csatornán is mintátvenni.
 %Például 
 s_digital = daq.createSession('ni');

%Digitális csatornát a addDigitalChannel paranccsal tudunk létrehozni. Ezzel led villogtatás, különböző kommunikációk és digitális szenzorok érzékelésére van leheteőség. A labor során ezeket nem használjuk, akit érdekel https://www.mathworks.com/help/daq/ref/adddigitalchannel.html oldalon tud utána olvasni.

%A
ch_counter=addCounterInputChannel(s_digital,devices(1).ID, 'ctr0', 'EdgeCount'); %paranccsal egy számláló tudunk inicializálni. Ez az s_digital munkafolyamathoz ad hozzá egy számlálót. A számolást a devices(1).ID mérőeszköz 'ctr0' csatornáján végzi. A laborban használt 6008-as kártyákon két bemenetet lehet használni számlálásra a 'ctr0'-t és 'ctr1'-et, azonban egyszerre csak az egyiket. 'EdgeCount' beállítás pedig azt jelenti, hogy a felfutó éleket (amikor a digitális jel 0-ról 1-re vált) számolja.

%A számláló aktuális értékét a már ismert 
impulzusokszama=inputSingleScan(s_digital);  %paranccsal tudjuk elérni.


 
 
 %A munkafolyamat n-edik csatornáját a
 removeChannel(s,n); %paranccsal törölhetjük az s munkafolyamatot. A csatornák sorszáma, a létrehozás sorrendjétől függ. A legkésőbb létrehozott csatorna lesz az utolsó. Ez a sorrend az adatok kiolvasásánál is számít.
%Természetesen ez korábbi mérési eredményen nem változtat, azonban egy új mérés elindítása után a törölt csatornán nem mér a program. 
 
 

 
 
daq.reset; %reseteli a munkafolyamaotot. Törli az eddigi sessiont, továbbá beolvassa újra az elérhető mérőeszközöket. Akkor lehet hasznos, ha futtattuk a daq.getDevices parancsot és utána csatlakoztattunk mérőkártyát a számítógéphez. A daq.getDevices nem érzékeli az újonnan (első lefutatása utáni) csatlakoztatott mérőkártyákat.

release(s); %az s munkafolyamat törlése. Ez akkor lehet hasznos, ha új munkafolyamatot szeretnénk indítani. (Az összes létrehozott csatornát törölni szeretnénk.) Továbbá, ha a munkafolymathoz bizonyos eszköz csatornáját hozzáadtuk azt az eszközt a matlab lefoglalja. Ez azt jelenti, hogy más program (pl az NI DAQmx vagy akár egy másik Matlab) nem férhet hozzá. A munkafolymat törlésével ezt a foglallást feloldjuk és a többi szoftver is hozzáfér. (Természetesen ez igaz minden szoftverre is. Pl nem lehet egyszerre Labviewból és Matlabból is használni a mérőkártyát. Ha nem tudjuk az "release"-lni az eszközt, célszerű a nem használt programot bezárni.)