%% 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 = daqlist; %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('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.

 ch_1 = addinput(s, devices.DeviceID(1), "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 = daqlist;  parancsot, akkor az elsõ mérõkártya neve devices.DeviceID(1).
 % 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=read(s,1,"OutputFormat","Matrix"); %Az s munkafolyamat csatornáin hajt végre egyszeri mérést. Ez hasznos lehet 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 DurationInSeconds változó beírásával. Pl 3 másodperces mérést az
 DurationInSeconds=3; %parancs segítségével adhatunk meg.

[M_data,timeStamps] = read(s,seconds(DurationInSeconds),"OutputFormat","Matrix");  % Ez az egyik legfontosabb parancs. Ha az "s" munkafolyamatot beálíltottuk, akkor ezzel a paranccsal indítjuk a mintavételezést.A függvény második bemenete meghatározza a mérésünk hosszát. Két lehetõségünk van, ha számot adunk meg, akkor annyi pontot fog lemérni a beállított mintavételezési frekvenciával, ha a seconds() függvényt használjuk, akkor lehetõség van másodpercekben megadni a mérés hosszát.
% A program elfoglalt a beállított idõtartamig, 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("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 lehetõ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.

ch_counter=addinput(s_digital, devices.DeviceID(1),"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.DeviceID(1) 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=read(s_digital,1,"OutputFormat","Matrix");  %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. 
 

daqreset; %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.)
